public static double incompleteGamma(double a,double x){
  double ans, ax, c, r;
  if (x <= 0 || a <= 0) {
    return 0.0;
  }
  if (x > 1.0 && x > a) {
    return 1.0 - incompleteGammaComplement(a,x);
  }
  ax=a * Math.log(x) - x - lnGamma(a);
  if (ax < -MAXLOG) {
    return (0.0);
  }
  ax=Math.exp(ax);
  r=a;
  c=1.0;
  ans=1.0;
  do {
    r+=1.0;
    c*=x / r;
    ans+=c;
  }
 while (c / ans > MACHEP);
  return (ans * ax / a);
}

public static double incompleteGamma ( double a , double x ) { double ans , ax , c , r ; if ( x <= 0 || a <= 0 ) { return 0.0 ; } if ( x > 1.0 && x > a ) { return 1.0 - incompleteGammaComplement ( a , x ) ; } ax = a * Math . log ( x ) - x - lnGamma ( a ) ; if ( ax < - MAXLOG ) { return ( 0.0 ) ; } ax = Math . exp ( ax ) ; r = a ; c = 1.0 ; ans = 1.0 ; do { r += 1.0 ; c *= x / r ; ans += c ; } while ( c / ans > MACHEP ) ; return ( ans * ax / a ) ; }
returns the incomplete gamma function .
