CMS 3D CMS Logo

CMSSW_4_4_3_patch1/src/CommonTools/Statistics/src/GammaContinuedFraction.cc

Go to the documentation of this file.
00001 #include "CommonTools/Statistics/src/GammaContinuedFraction.h"
00002 #include <cmath>
00003 #include <iostream>
00004 
00005 #define ITMAX 100        // maximum allowed number of iterations
00006 #define EPS 3.0e-7       // relative accuracy
00007 #define FPMIN 1.0e-30    // number near the smallest representable floating-point number
00008 
00009 float
00010 GammaContinuedFraction( float a, float x )
00011 {
00012   int i;
00013   float an,del;
00014 
00015   /* Set up for evaluating continued fraction by modified Lentz's method (par.5.2
00016      in Numerical Recipes in C) with b_0 = 0 */
00017   double b = x+1.0-a;
00018   double c = 1.0/FPMIN;
00019   double d = 1.0/b;
00020   double h = d;
00021   for (i=1;i<=ITMAX;i++) {
00022     an = -i*(i-a);
00023     b += 2.0;
00024     d=an*d+b;
00025     if (fabs(d) < FPMIN) d=FPMIN;
00026     c=b+an/c;
00027     if (fabs(c) < FPMIN) c=FPMIN;
00028     d=1.0/d;
00029     del=d*c;
00030     h *= del;
00031     if (fabs(del-1.0) < EPS) break;
00032   }
00033   if( i > ITMAX ) std::cerr << "GammaContinuedFraction::a too large, "
00034                        << "ITMAX too small" << std::endl;
00035   return h;
00036 }
00037 #undef ITMAX
00038 #undef EPS
00039 #undef FPMIN
00040 /* (C) Copr. 1986-92 Numerical Recipes Software B2.. */