Impala Documentation

template<class DstArrayT, class SrcArrayT>
double Impala::Core::Feature::IntWeibullNgbPnLoop< DstArrayT, SrcArrayT >::DiGamma ( double x ) const [inline]
Definition at line 277 of file IntWeibullNgbPnLoop.h.
Referenced by Impala::Core::Feature::IntWeibullNgbPnLoop< DstArrayT, SrcArrayT >::Done().
00278     {
00279         /* force into the interval 1..3 */
00280         if( x < 0.0 )
00281             return DiGamma(1.0-x)+M_PI/tan(M_PI*(1.0-x)) ;  /* refection formula */
00282         else if( x < 1.0 )
00283             return DiGamma(1.0+x)-1.0/x ;
00284         else if ( x == 1.0)
00285             return -MMGAMMA ;
00286         else if ( x == 2.0)
00287             return 1.0-MMGAMMA ;
00288         else if ( x == 3.0)
00289             return 1.5-MMGAMMA ;
00290         else if ( x > 3.0)
00291             /* duplication formula */
00292             return 0.5*(DiGamma(x/2.0)+DiGamma((x+1.0)/2.0))+M_LN2 ;
00293         else
00294         {
00295             static double Kncoe[] =
00296                 { .30459198558715155634315638246624251,
00297                   .72037977439182833573548891941219706, -.12454959243861367729528855995001087,
00298                   .27769457331927827002810119567456810e-1, -.67762371439822456447373550186163070e-2,
00299                   .17238755142247705209823876688592170e-2, -.44817699064252933515310345718960928e-3,
00300                   .11793660000155572716272710617753373e-3, -.31253894280980134452125172274246963e-4,
00301                   .83173997012173283398932708991137488e-5, -.22191427643780045431149221890172210e-5,
00302                   .59302266729329346291029599913617915e-6, -.15863051191470655433559920279603632e-6,
00303                   .42459203983193603241777510648681429e-7, -.11369129616951114238848106591780146e-7,
00304                   .304502217295931698401459168423403510e-8, -.81568455080753152802915013641723686e-9,
00305                   .21852324749975455125936715817306383e-9, -.58546491441689515680751900276454407e-10,
00306                   .15686348450871204869813586459513648e-10, -.42029496273143231373796179302482033e-11,
00307                   .11261435719264907097227520956710754e-11, -.30174353636860279765375177200637590e-12,
00308                   .80850955256389526647406571868193768e-13, -.21663779809421233144009565199997351e-13,
00309                   .58047634271339391495076374966835526e-14, -.15553767189204733561108869588173845e-14,
00310                   .41676108598040807753707828039353330e-15, -.11167065064221317094734023242188463e-15 } ;
00311 
00312             double Tn_1 = 1.0 ; /* T_{n-1}(x), started at n=1 */
00313             double Tn = x-2.0 ; /* T_{n}(x) , started at n=1 */
00314             double resul = Kncoe[0] + Kncoe[1]*Tn ;
00315 
00316             x -= 2.0 ;
00317 
00318             for(int n = 2 ; n < sizeof(Kncoe)/sizeof(double) ;n++)
00319             {
00320                 const double Tn1 = 2.0 * x * Tn - Tn_1 ;    /* Chebyshev recursion, Eq. 22.7.4 Abramowitz-Stegun */
00321                 resul += Kncoe[n]*Tn1 ;
00322                 Tn_1 = Tn ;
00323                 Tn = Tn1 ;
00324             }
00325             return resul ;
00326         }
00327     } //  DiGamma(double x)