LCOV - code coverage report
Current view: top level - libavutil/tests - lfg.c (source / functions) Hit Total Coverage
Test: coverage.info Lines: 62 66 93.9 %
Date: 2017-10-18 21:45:51 Functions: 2 2 100.0 %

          Line data    Source code
       1             : /*
       2             :  * This file is part of FFmpeg.
       3             :  *
       4             :  * FFmpeg is free software; you can redistribute it and/or
       5             :  * modify it under the terms of the GNU Lesser General Public
       6             :  * License as published by the Free Software Foundation; either
       7             :  * version 2.1 of the License, or (at your option) any later version.
       8             :  *
       9             :  * FFmpeg is distributed in the hope that it will be useful,
      10             :  * but WITHOUT ANY WARRANTY; without even the implied warranty of
      11             :  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
      12             :  * Lesser General Public License for more details.
      13             :  *
      14             :  * You should have received a copy of the GNU Lesser General Public
      15             :  * License along with FFmpeg; if not, write to the Free Software
      16             :  * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
      17             :  */
      18             : 
      19             : #include "libavutil/log.h"
      20             : #include "libavutil/timer.h"
      21             : #include "libavutil/lfg.h"
      22             : 
      23             : static const double Z_TABLE[31][10] = {
      24             :     {0.5000,  0.5040,  0.5080,  0.5120,  0.5160,  0.5199,  0.5239,  0.5279,  0.5319,  0.5359},
      25             :     {0.5398,  0.5438,  0.5478,  0.5517,  0.5557,  0.5596,  0.5636,  0.5675,  0.5714,  0.5753},
      26             :     {0.5793,  0.5832,  0.5871,  0.5910,  0.5948,  0.5987,  0.6026,  0.6064,  0.6103,  0.6141},
      27             :     {0.6179,  0.6217,  0.6255,  0.6293,  0.6331,  0.6368,  0.6406,  0.6443,  0.6480,  0.6517},
      28             :     {0.6554,  0.6591,  0.6628,  0.6664,  0.6700,  0.6736,  0.6772,  0.6808,  0.6844,  0.6879},
      29             :     {0.6915,  0.6950,  0.6985,  0.7019,  0.7054,  0.7088,  0.7123,  0.7157,  0.7190,  0.7224},
      30             :     {0.7257,  0.7291,  0.7324,  0.7357,  0.7389,  0.7422,  0.7454,  0.7486,  0.7517,  0.7549},
      31             :     {0.7580,  0.7611,  0.7642,  0.7673,  0.7704,  0.7734,  0.7764,  0.7794,  0.7823,  0.7852},
      32             :     {0.7881,  0.7910,  0.7939,  0.7967,  0.7995,  0.8023,  0.8051,  0.8078,  0.8106,  0.8133},
      33             :     {0.8159,  0.8186,  0.8212,  0.8238,  0.8264,  0.8289,  0.8315,  0.8340,  0.8365,  0.8389},
      34             :     {0.8413,  0.8438,  0.8461,  0.8485,  0.8508,  0.8531,  0.8554,  0.8577,  0.8599,  0.8621},
      35             :     {0.8643,  0.8665,  0.8686,  0.8708,  0.8729,  0.8749,  0.8770,  0.8790,  0.8810,  0.8830},
      36             :     {0.8849,  0.8869,  0.8888,  0.8907,  0.8925,  0.8944,  0.8962,  0.8980,  0.8997,  0.9015},
      37             :     {0.9032,  0.9049,  0.9066,  0.9082,  0.9099,  0.9115,  0.9131,  0.9147,  0.9162,  0.9177},
      38             :     {0.9192,  0.9207,  0.9222,  0.9236,  0.9251,  0.9265,  0.9279,  0.9292,  0.9306,  0.9319},
      39             :     {0.9332,  0.9345,  0.9357,  0.9370,  0.9382,  0.9394,  0.9406,  0.9418,  0.9429,  0.9441},
      40             :     {0.9452,  0.9463,  0.9474,  0.9484,  0.9495,  0.9505,  0.9515,  0.9525,  0.9535,  0.9545},
      41             :     {0.9554,  0.9564,  0.9573,  0.9582,  0.9591,  0.9599,  0.9608,  0.9616,  0.9625,  0.9633},
      42             :     {0.9641,  0.9649,  0.9656,  0.9664,  0.9671,  0.9678,  0.9686,  0.9693,  0.9699,  0.9706},
      43             :     {0.9713,  0.9719,  0.9726,  0.9732,  0.9738,  0.9744,  0.9750,  0.9756,  0.9761,  0.9767},
      44             :     {0.9772,  0.9778,  0.9783,  0.9788,  0.9793,  0.9798,  0.9803,  0.9808,  0.9812,  0.9817},
      45             :     {0.9821,  0.9826,  0.9830,  0.9834,  0.9838,  0.9842,  0.9846,  0.9850,  0.9854,  0.9857},
      46             :     {0.9861,  0.9864,  0.9868,  0.9871,  0.9875,  0.9878,  0.9881,  0.9884,  0.9887,  0.9890},
      47             :     {0.9893,  0.9896,  0.9898,  0.9901,  0.9904,  0.9906,  0.9909,  0.9911,  0.9913,  0.9916},
      48             :     {0.9918,  0.9920,  0.9922,  0.9925,  0.9927,  0.9929,  0.9931,  0.9932,  0.9934,  0.9936},
      49             :     {0.9938,  0.9940,  0.9941,  0.9943,  0.9945,  0.9946,  0.9948,  0.9949,  0.9951,  0.9952},
      50             :     {0.9953,  0.9955,  0.9956,  0.9957,  0.9959,  0.9960,  0.9961,  0.9962,  0.9963,  0.9964},
      51             :     {0.9965,  0.9966,  0.9967,  0.9968,  0.9969,  0.9970,  0.9971,  0.9972,  0.9973,  0.9974},
      52             :     {0.9974,  0.9975,  0.9976,  0.9977,  0.9977,  0.9978,  0.9979,  0.9979,  0.9980,  0.9981},
      53             :     {0.9981,  0.9982,  0.9982,  0.9983,  0.9984,  0.9984,  0.9985,  0.9985,  0.9986,  0.9986},
      54             :     {0.9987,  0.9987,  0.9987,  0.9988,  0.9988,  0.9989,  0.9989,  0.9989,  0.9990,  0.9990} };
      55             : 
      56             : // Inverse cumulative distribution function
      57        1998 : static double inv_cdf(double u)
      58             : {
      59             :     const double a[4] = { 2.50662823884,
      60             :                          -18.61500062529,
      61             :                           41.39119773534,
      62             :                          -25.44106049637};
      63             : 
      64             :     const double b[4] = {-8.47351093090,
      65             :                           23.08336743743,
      66             :                          -21.06224101826,
      67             :                            3.13082909833};
      68             : 
      69             :     const double c[9] = {0.3374754822726147,
      70             :                           0.9761690190917186,
      71             :                           0.1607979714918209,
      72             :                           0.0276438810333863,
      73             :                           0.0038405729373609,
      74             :                           0.0003951896511919,
      75             :                           0.0000321767881768,
      76             :                           0.0000002888167364,
      77             :                           0.0000003960315187};
      78             : 
      79             :     double r;
      80        1998 :     double x = u - 0.5;
      81             : 
      82             :     // Beasley-Springer
      83        1998 :      if (fabs(x) < 0.42) {
      84             : 
      85        1680 :         double y = x * x;
      86        3360 :         r        = x * (((a[3]*y+a[2])*y+a[1])*y+a[0]) /
      87        1680 :                         ((((b[3]*y+b[2])*y+b[1])*y+b[0])*y+1.0);
      88             :     }
      89             :     else {// Moro
      90         318 :         r = u;
      91         318 :         if (x > 0.0)
      92         159 :             r = 1.0 - u;
      93         318 :         r = log(-log(r));
      94         636 :         r = c[0] + r*(c[1]+r*(c[2]+r*(c[3]+r*(c[4]+r*(c[5]+r*(c[6]+
      95         318 :                                                         r*(c[7]+r*c[8])))))));
      96         318 :         if (x < 0.0)
      97         159 :             r = -r;
      98             :     }
      99             : 
     100        1998 :     return r;
     101             : }
     102           1 : int main(void)
     103             : {
     104           1 :     int x = 0;
     105             :     int i, j;
     106             :     AVLFG state;
     107           1 :     av_lfg_init(&state, 0xdeadbeef);
     108       10001 :     for (j = 0; j < 10000; j++) {
     109     6250000 :         for (i = 0; i < 624; i++) {
     110             :             //av_log(NULL, AV_LOG_ERROR, "%X\n", av_lfg_get(&state));
     111     6240000 :             x += av_lfg_get(&state);
     112             :         }
     113             :     }
     114           1 :     av_log(NULL, AV_LOG_ERROR, "final value:%X\n", x);
     115             : 
     116             :     /* BMG usage example */
     117             :     {
     118           1 :         double mean   = 1000;
     119           1 :         double stddev = 53;
     120           1 :         double samp_mean = 0.0, samp_stddev = 0.0, QH = 0;
     121           1 :         double Z, p_value = -1, tot_samp = 1000;
     122           1 :         double *PRN_arr = av_malloc_array(tot_samp, sizeof(double));
     123             : 
     124           1 :         if (!PRN_arr) {
     125           0 :             fprintf(stderr, "failed to allocate memory!\n");
     126           0 :             return 1;
     127             :         }
     128             : 
     129           1 :         av_lfg_init(&state, 42);
     130         501 :         for (i = 0; i < tot_samp; i += 2) {
     131             :             double bmg_out[2];
     132         500 :             av_bmg_get(&state, bmg_out);
     133         500 :             PRN_arr[i  ] = bmg_out[0] * stddev + mean;
     134         500 :             PRN_arr[i+1] = bmg_out[1] * stddev + mean;
     135         500 :             samp_mean   += PRN_arr[i] + PRN_arr[i+1];
     136         500 :             samp_stddev += PRN_arr[i] * PRN_arr[i] + PRN_arr[i+1] * PRN_arr[i+1];
     137        1000 :             printf("PRN%d : %f\n"
     138             :                    "PRN%d : %f\n",
     139        1000 :                    i, PRN_arr[i], i+1, PRN_arr[i+1]);
     140             :         }
     141           1 :         samp_mean /= tot_samp;
     142           1 :         samp_stddev /= (tot_samp - 1);
     143           1 :         samp_stddev -= (tot_samp * 1.0 / (tot_samp - 1))*samp_mean*samp_mean;
     144           1 :         samp_stddev = sqrt(samp_stddev);
     145           1 :         Z = (mean - samp_mean) / (stddev / sqrt(tot_samp));
     146             :         {
     147           1 :             int x, y, a, b, flag = 0;
     148             : 
     149           1 :             if (Z < 0.0) {
     150           1 :                 flag = !flag;
     151           1 :                 Z = Z * -1.0;
     152             :             }
     153             : 
     154           1 :             a = (int)(Z * 100);
     155           1 :             b = ((int)Z * 100);
     156           1 :             x = Z * 10;
     157           1 :             y = (b > 0) ? a % b : a;
     158           1 :             y = y % 10;
     159           1 :             if (x > 30 || y > 9) {
     160           0 :                 av_log(NULL, AV_LOG_INFO, "error: out of bounds! tried to access"
     161             :                                           "Z_TABLE[%d][%d]\n", x, y);
     162           0 :                 goto SKIP;
     163             :             }
     164           1 :             p_value = flag ? 1 - Z_TABLE[x][y] : Z_TABLE[x][y];
     165             :         }
     166             : 
     167        1001 : SKIP:   for (i = 0; i < tot_samp; ++i) {
     168             : 
     169        1000 :             if ( i < (tot_samp - 1)) {
     170             :                 double H_diff;
     171         999 :                 H_diff  = inv_cdf((i + 2.0 - (3.0/8.0)) / (tot_samp + (1.0/4.0)));
     172         999 :                 H_diff -= inv_cdf((i + 1.0 - (3.0/8.0)) / (tot_samp + (1.0/4.0)));
     173             : 
     174         999 :                 QH     += ((PRN_arr[i + 1] - PRN_arr[i]) / H_diff);
     175             :             }
     176             :         }
     177           1 :         QH = 1.0 - QH / ((tot_samp - 1.0) * samp_stddev);
     178             : 
     179           1 :         printf("sample mean  : %f\n"
     180             :                "true mean    : %f\n"
     181             :                "sample stddev: %f\n"
     182             :                "true stddev  : %f\n"
     183             :                "z-score      : %f\n"
     184             :                "p-value      : %f\n"
     185             :                "QH[normality]: %f\n",
     186             :                samp_mean, mean, samp_stddev, stddev, Z, p_value, QH);
     187             : 
     188           1 :         av_freep(&PRN_arr);
     189             :     }
     190           1 :     return 0;
     191             : }

Generated by: LCOV version 1.13