docs: function gammaln( x )  {
  	var isNegative;
  	var nadj;
  	var flg;
  	var p3;
  	var p2;
  	var p1;
  	var p;
  	var q;
  	var t;
  	var w;
  	var y;
  	var z;
  	var r;
  
  	// Special cases: NaN
  	if ( x !== x ) {
  		return x;
  	}
  	// Special case: +-infinity
  	if ( x === PINF || x === NINF ) {
  		return x;
  	}
  	// Special case: 0
  	if ( x === 0 ) {
  		return PINF;
  	}
  	if ( x < 0 ) {
  		isNegative = true;
  		x = -x;
  	} else {
  		isNegative = false;
  	}
  	// If |x| < 2**-70, return -ln(|x|)
  	if ( x < TINY ) {
  		return -ln( x );
  	}
  	if ( isNegative ) {
  		// If |x| >= 2**52, must be -integer
  		if ( x >= TWO52 ) {
  			return PINF;
  		}
  		t = sinpi( x );
  		if ( t === 0 ) {
  			return PINF;
  		}
  		nadj = ln( PI / abs( t*x ) );
  	}
  	// If x equals 1 or 2, return 0
  	if ( x === 1 || x === 2 ) {
  		return 0;
  	}
  	// If x < 2, use lgamma(x) = lgamma(x+1) - log(x)
  	if ( x < 2 ) {
  		if ( x <= 0.9 ) {
  			r = -ln( x );
  
  			// 0.7316 <= x <=  0.9
  			if ( x >= ( YMIN - 1 + 0.27 ) ) {
  				y = 1 - x;
  				flg = 0;
  			}
  			// 0.2316 <= x < 0.7316
  			else if ( x >= (YMIN - 1 - 0.27) ) {
  				y = x - (TC - 1);
  				flg = 1;
  			}
  			// 0 < x < 0.2316
  			else {
  				y = x;
  				flg = 2;
  			}
  		} else {
  			r = 0;
  
  			// 1.7316 <= x < 2
  			if ( x >= (YMIN + 0.27) ) {
  				y = 2 - x;
  				flg = 0;
  			}
  			// 1.2316 <= x < 1.7316
  			else if ( x >= (YMIN - 0.27) ) {
  				y = x - TC;
  				flg = 1;
  			}
  			// 0.9 < x < 1.2316
  			else {
  				y = x - 1;
  				flg = 2;
  			}
  		}
  		switch ( flg ) {
  		case 0:
  			z = y * y;
  			p1 = A1C + z*polyvalA1( z );
  			p2 = z * (A2C + z*polyvalA2( z ));
  			p = y*p1 + p2;
  			r += ( p - 0.5 * y );
  			break;
  		case 1:
  			z = y * y;
  			w = z * y;
  			p1 = T1C + w*polyvalT1( w );
  			p2 = T2C + w*polyvalT2( w );
  			p3 = T3C + w*polyvalT3( w );
  			p = z*p1 - (TT - w*(p2+y*p3));
  			r += ( TF + p );
  			break;
  		case 2:
  			p1 = y * (UC + y*polyvalU( y ));
  			p2 = VC + y*polyvalV( y );
  			r += ( -0.5 * y + p1/p2 );
  			break;
  		}
  	}
  	// 2 <= x < 8
  	else if ( x < 8 ) {
  		flg = trunc( x );
  		y = x - flg;
  		p = y * (SC + y*polyvalS( y ));
  		q = RC + y*polyvalR( y );
  		r = 0.5*y + p/q;
  		z = 1.0; // gammaln(1+s) = ln(s) + gammaln(s)
  		switch ( flg ) {
  		case 7:
  			z *= (y + 6);
  			/* falls through */
  		case 6:
  			z *= (y + 5);
  			/* falls through */
  		case 5:
  			z *= (y + 4);
  			/* falls through */
  		case 4:
  			z *= (y + 3);
  			/* falls through */
  		case 3:
  			z *= ( y + 2 );
  			r += ln( z );
  		}
  	}
  	// 8 <= x < 2**58
  	else if ( x < TWO58 ) {
  		t = ln( x );
  		z = 1 / x;
  		y = z * z;
  		w = WC + z*polyvalW( y );
  		r = (x-0.5)*(t-1) + w;
  	}
  	// 2**58 <= x <= Inf
  	else {
  		r = x * ( ln(x)-1 );
  	}
  	if ( isNegative ) {
  		r = nadj - r;
  	}
  	return r;
  }