static double lgamma()

in commons-numbers-gamma/src/main/java/org/apache/commons/numbers/gamma/BoostGamma.java [757:815]
44 lines of code
12 McCabe index (conditional complexity)

    static double lgamma(double z, int[] sign) {
        double result;
        int sresult = 1;
        if (z <= -ROOT_EPSILON) {
            // reflection formula:
            if (Math.rint(z) == z) {
                // Pole error
                return Double.NaN;
            }

            double t = sinpx(z);
            z = -z;
            if (t < 0) {
                t = -t;
            } else {
                sresult = -sresult;
            }

            // This summation can have large magnitudes with opposite signs.
            // Use an extended precision sum to reduce cancellation.
            result = Sum.of(-lgamma(z)).add(-Math.log(t)).add(LOG_PI).getAsDouble();
        } else if (z < ROOT_EPSILON) {
            if (z == 0) {
                // Pole error
                return Double.NaN;
            }
            if (4 * Math.abs(z) < EPSILON) {
                result = -Math.log(Math.abs(z));
            } else {
                result = Math.log(Math.abs(1 / z - EULER));
            }
            if (z < 0) {
                sresult = -1;
            }
        } else if (z < 15) {
            result = lgammaSmall(z, z - 1, z - 2);
        // The z > 3 condition is always true
        //} else if (z > 3 && z < 100) {
        } else if (z < 100) {
            // taking the log of tgamma reduces the error, no danger of overflow here:
            result = Math.log(tgamma(z));
        } else {
            // regular evaluation:
            final double zgh = z + Lanczos.GMH;
            result = Math.log(zgh) - 1;
            result *= z - 0.5f;
            //
            // Only add on the lanczos sum part if we're going to need it:
            //
            if (result * EPSILON < 20) {
                result += Math.log(Lanczos.lanczosSumExpGScaled(z));
            }
        }

        if (nonZeroLength(sign)) {
            sign[0] = sresult;
        }
        return result;
    }