htswanalysis/MACS/lib/gsl/gsl-1.11/histogram/test1d.c

   1 /* histogram/test.c
   2  *
   3  * Copyright (C) 1996, 1997, 1998, 1999, 2000, 2007 Brian Gough
   4  *
   5  * This program is free software; you can redistribute it and/or modify
   6  * it under the terms of the GNU General Public License as published by
   7  * the Free Software Foundation; either version 3 of the License, or (at
   8  * your option) any later version.
   9  *
  10  * This program is distributed in the hope that it will be useful, but
  11  * WITHOUT ANY WARRANTY; without even the implied warranty of
  12  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  13  * General Public License for more details.
  14  *
  15  * You should have received a copy of the GNU General Public License
  16  * along with this program; if not, write to the Free Software
  17  * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA.
  18  */
  19
  20 #include <config.h>
  21 #include <stdlib.h>
  22 #include <stdio.h>
  23 #include <gsl/gsl_histogram.h>
  24 #include <gsl/gsl_test.h>
  25 #include <gsl/gsl_ieee_utils.h>
  26
  27 #define N 397
  28 #define NR 10
  29
  30 void
  31 test1d (void)
  32 {
  33   double xr[NR + 1] =
  34   {0.0, 1.0, 2.0, 3.0, 4.0, 5.0, 6.0, 7.0, 8.0, 9.0, 10.0};
  35
  36   gsl_histogram *h, *h1, *hr, *g;
  37   size_t i, j;
  38
  39   gsl_ieee_env_setup ();
  40
  41   h = gsl_histogram_calloc (N);
  42   h1 = gsl_histogram_calloc (N);
  43   g = gsl_histogram_calloc (N);
  44
  45   gsl_test (h->range == 0, "gsl_histogram_alloc returns valid range pointer");
  46   gsl_test (h->bin == 0, "gsl_histogram_alloc returns valid bin pointer");
  47   gsl_test (h->n != N, "gsl_histogram_alloc returns valid size");
  48
  49
  50   hr = gsl_histogram_calloc_range (NR, xr);
  51
  52   gsl_test (hr->range == 0, "gsl_histogram_calloc_range returns valid range pointer");
  53   gsl_test (hr->bin == 0, "gsl_histogram_calloc_range returns valid bin pointer");
  54   gsl_test (hr->n != NR, "gsl_histogram_calloc_range returns valid size");
  55
  56   {
  57     int status = 0;
  58     for (i = 0; i <= NR; i++)
  59       {
  60         if (hr->range[i] != xr[i])
  61           {
  62             status = 1;
  63           }
  64       };
  65
  66     gsl_test (status, "gsl_histogram_calloc_range creates range");
  67   }
  68
  69   for (i = 0; i <= NR; i++)
  70     {
  71       hr->range[i] = 0.0;
  72     }
  73
  74   {
  75     int status = gsl_histogram_set_ranges (hr, xr, NR+1);
  76
  77     for (i = 0; i <= NR; i++)
  78       {
  79         if (hr->range[i] != xr[i])
  80           {
  81             status = 1;
  82           }
  83       };
  84
  85     gsl_test (status, "gsl_histogram_set_range sets range");
  86   }
  87
  88
  89   for (i = 0; i < N; i++)
  90     {
  91       gsl_histogram_accumulate (h, (double) i, (double) i);
  92     };
  93
  94   {
  95     int status = 0;
  96
  97     for (i = 0; i < N; i++)
  98       {
  99         if (h->bin[i] != (double) i)
 100           {
 101             status = 1;
 102           }
 103       };
 104
 105     gsl_test (status, "gsl_histogram_accumulate writes into array");
 106   }
 107
 108
 109   {
 110     int status = 0;
 111
 112     for (i = 0; i < N; i++)
 113       {
 114         if (gsl_histogram_get (h, i) != i)
 115           status = 1;
 116       };
 117     gsl_test (status, "gsl_histogram_get reads from array");
 118   }
 119
 120   for (i = 0; i <= N; i++)
 121     {
 122       h1->range[i] = 100.0 + i;
 123     }
 124
 125   gsl_histogram_memcpy (h1, h);
 126
 127   {
 128     int status = 0;
 129     for (i = 0; i <= N; i++)
 130       {
 131         if (h1->range[i] != h->range[i])
 132           status = 1;
 133       };
 134     gsl_test (status, "gsl_histogram_memcpy copies bin ranges");
 135   }
 136
 137   {
 138     int status = 0;
 139     for (i = 0; i < N; i++)
 140       {
 141         if (gsl_histogram_get (h1, i) != gsl_histogram_get (h, i))
 142           status = 1;
 143       };
 144     gsl_test (status, "gsl_histogram_memcpy copies bin values");
 145   }
 146
 147   gsl_histogram_free (h1);
 148
 149   h1 = gsl_histogram_clone (h);
 150
 151   {
 152     int status = 0;
 153     for (i = 0; i <= N; i++)
 154       {
 155         if (h1->range[i] != h->range[i])
 156           status = 1;
 157       };
 158     gsl_test (status, "gsl_histogram_clone copies bin ranges");
 159   }
 160
 161   {
 162     int status = 0;
 163     for (i = 0; i < N; i++)
 164       {
 165         if (gsl_histogram_get (h1, i) != gsl_histogram_get (h, i))
 166           status = 1;
 167       };
 168     gsl_test (status, "gsl_histogram_clone copies bin values");
 169   }
 170
 171   gsl_histogram_reset (h);
 172
 173   {
 174     int status = 0;
 175
 176     for (i = 0; i < N; i++)
 177       {
 178         if (h->bin[i] != 0)
 179           status = 1;
 180       }
 181     gsl_test (status, "gsl_histogram_reset zeros array");
 182   }
 183
 184
 185   {
 186     int status = 0;
 187
 188     for (i = 0; i < N; i++)
 189       {
 190         gsl_histogram_increment (h, (double) i);
 191
 192         for (j = 0; j <= i; j++)
 193           {
 194             if (h->bin[j] != 1)
 195               {
 196                 status = 1;
 197               }
 198           }
 199
 200         for (j = i + 1; j < N; j++)
 201           {
 202             if (h->bin[j] != 0)
 203               {
 204                 status = 1;
 205               }
 206           }
 207       }
 208
 209     gsl_test (status, "gsl_histogram_increment increases bin value");
 210   }
 211
 212   {
 213     int status = 0;
 214     for (i = 0; i < N; i++)
 215       {
 216         double x0 = 0, x1 = 0;
 217
 218         gsl_histogram_get_range (h, i, &x0, &x1);
 219
 220         if (x0 != i || x1 != i + 1)
 221           {
 222             status = 1;
 223           }
 224       }
 225     gsl_test (status, "gsl_histogram_getbinrange returns bin range");
 226   }
 227
 228   {
 229     int status = 0;
 230     if (gsl_histogram_max (h) != N)
 231       status = 1;
 232     gsl_test (status, "gsl_histogram_max returns maximum");
 233   }
 234
 235   {
 236     int status = 0;
 237     if (gsl_histogram_min (h) != 0)
 238       status = 1;
 239     gsl_test (status, "gsl_histogram_min returns minimum");
 240   }
 241
 242   {
 243     int status = 0;
 244     if (gsl_histogram_bins (h) != N)
 245       status = 1;
 246     gsl_test (status, "gsl_histogram_bins returns number of bins");
 247   }
 248
 249   h->bin[2] = 123456.0;
 250   h->bin[4] = -654321;
 251
 252   {
 253     double max = gsl_histogram_max_val (h);
 254     gsl_test (max != 123456.0, "gsl_histogram_max_val finds maximum value");
 255   }
 256
 257   {
 258     double min = gsl_histogram_min_val (h);
 259     gsl_test (min != -654321.0, "gsl_histogram_min_val finds minimum value");
 260   }
 261
 262   {
 263     size_t imax = gsl_histogram_max_bin (h);
 264     gsl_test (imax != 2, "gsl_histogram_max_bin finds maximum value bin");
 265   }
 266
 267   {
 268     size_t imin = gsl_histogram_min_bin (h);
 269     gsl_test (imin != 4, "gsl_histogram_min_bin find minimum value bin");
 270   }
 271
 272   for (i = 0; i < N; i++)
 273     {
 274       h->bin[i] = i + 27;
 275       g->bin[i] = (i + 27) * (i + 1);
 276     }
 277
 278   {
 279     double sum=gsl_histogram_sum (h);
 280     gsl_test(sum != N*27+((N-1)*N)/2, "gsl_histogram_sum sums all bin values");
 281   }
 282
 283   gsl_histogram_memcpy (h1, g);
 284   gsl_histogram_add (h1, h);
 285
 286   {
 287     int status = 0;
 288     for (i = 0; i < N; i++)
 289       {
 290         if (h1->bin[i] != g->bin[i] + h->bin[i])
 291           status = 1;
 292       }
 293     gsl_test (status, "gsl_histogram_add histogram addition");
 294   }
 295
 296   gsl_histogram_memcpy (h1, g);
 297   gsl_histogram_sub (h1, h);
 298
 299   {
 300     int status = 0;
 301     for (i = 0; i < N; i++)
 302       {
 303         if (h1->bin[i] != g->bin[i] - h->bin[i])
 304           status = 1;
 305       }
 306     gsl_test (status, "gsl_histogram_sub histogram subtraction");
 307   }
 308
 309
 310   gsl_histogram_memcpy (h1, g);
 311   gsl_histogram_mul (h1, h);
 312
 313   {
 314     int status = 0;
 315     for (i = 0; i < N; i++)
 316       {
 317         if (h1->bin[i] != g->bin[i] * h->bin[i])
 318           status = 1;
 319       }
 320     gsl_test (status, "gsl_histogram_mul histogram multiplication");
 321   }
 322
 323
 324   gsl_histogram_memcpy (h1, g);
 325   gsl_histogram_div (h1, h);
 326
 327   {
 328     int status = 0;
 329     for (i = 0; i < N; i++)
 330       {
 331         if (h1->bin[i] != g->bin[i] / h->bin[i])
 332           status = 1;
 333       }
 334     gsl_test (status, "gsl_histogram_div histogram division");
 335   }
 336
 337   gsl_histogram_memcpy (h1, g);
 338   gsl_histogram_scale (h1, 0.5);
 339
 340   {
 341     int status = 0;
 342     for (i = 0; i < N; i++)
 343       {
 344         if (h1->bin[i] != 0.5 * g->bin[i])
 345           status = 1;
 346       }
 347     gsl_test (status, "gsl_histogram_scale histogram scaling");
 348   }
 349
 350   gsl_histogram_memcpy (h1, g);
 351   gsl_histogram_shift (h1, 0.25);
 352
 353   {
 354     int status = 0;
 355     for (i = 0; i < N; i++)
 356       {
 357         if (h1->bin[i] != 0.25 + g->bin[i])
 358           status = 1;
 359       }
 360     gsl_test (status, "gsl_histogram_shift histogram shift");
 361   }
 362
 363
 364   gsl_histogram_free (h);       /* free whatever is in h */
 365
 366   h = gsl_histogram_calloc_uniform (N, 0.0, 1.0);
 367
 368   gsl_test (h->range == 0,
 369             "gsl_histogram_calloc_uniform returns valid range pointer");
 370   gsl_test (h->bin == 0,
 371             "gsl_histogram_calloc_uniform returns valid bin pointer");
 372   gsl_test (h->n != N,
 373             "gsl_histogram_calloc_uniform returns valid size");
 374
 375   gsl_histogram_accumulate (h, 0.0, 1.0);
 376   gsl_histogram_accumulate (h, 0.1, 2.0);
 377   gsl_histogram_accumulate (h, 0.2, 3.0);
 378   gsl_histogram_accumulate (h, 0.3, 4.0);
 379
 380   {
 381     size_t i1, i2, i3, i4;
 382     double expected;
 383     int status = gsl_histogram_find (h, 0.0, &i1);
 384     status = gsl_histogram_find (h, 0.1, &i2);
 385     status = gsl_histogram_find (h, 0.2, &i3);
 386     status = gsl_histogram_find (h, 0.3, &i4);
 387
 388     for (i = 0; i < N; i++)
 389       {
 390         if (i == i1)
 391           {
 392             expected = 1.0;
 393           }
 394         else if (i == i2)
 395           {
 396             expected = 2.0;
 397           }
 398         else if (i == i3)
 399           {
 400             expected = 3.0;
 401           }
 402         else if (i == i4)
 403           {
 404             expected = 4.0;
 405           }
 406         else
 407           {
 408             expected = 0.0;
 409           }
 410
 411         if (h->bin[i] != expected)
 412           {
 413             status = 1;
 414           }
 415
 416       }
 417     gsl_test (status, "gsl_histogram_find returns index");
 418   }
 419
 420
 421   {
 422     FILE *f = fopen ("test.txt", "w");
 423     gsl_histogram_fprintf (f, h, "%.19e", "%.19e");
 424     fclose (f);
 425   }
 426
 427   {
 428     FILE *f = fopen ("test.txt", "r");
 429     gsl_histogram *hh = gsl_histogram_calloc (N);
 430     int status = 0;
 431
 432     gsl_histogram_fscanf (f, hh);
 433
 434     for (i = 0; i < N; i++)
 435       {
 436         if (h->range[i] != hh->range[i])
 437           status = 1;
 438         if (h->bin[i] != hh->bin[i])
 439           status = 1;
 440       }
 441     if (h->range[N] != hh->range[N])
 442       status = 1;
 443
 444     gsl_test (status, "gsl_histogram_fprintf and fscanf");
 445
 446     gsl_histogram_free (hh);
 447     fclose (f);
 448   }
 449
 450   {
 451     FILE *f = fopen ("test.dat", "wb");
 452     gsl_histogram_fwrite (f, h);
 453     fclose (f);
 454   }
 455
 456   {
 457     FILE *f = fopen ("test.dat", "rb");
 458     gsl_histogram *hh = gsl_histogram_calloc (N);
 459     int status = 0;
 460
 461     gsl_histogram_fread (f, hh);
 462
 463     for (i = 0; i < N; i++)
 464       {
 465         if (h->range[i] != hh->range[i])
 466           status = 1;
 467         if (h->bin[i] != hh->bin[i])
 468           status = 1;
 469       }
 470     if (h->range[N] != hh->range[N])
 471       status = 1;
 472
 473     gsl_test (status, "gsl_histogram_fwrite and fread");
 474
 475     gsl_histogram_free (hh);
 476     fclose (f);
 477   }
 478
 479   gsl_histogram_free (h);
 480   gsl_histogram_free (g);
 481   gsl_histogram_free (h1);
 482   gsl_histogram_free (hr);
 483 }