#include <math.h>
#include "nrutil.h"

double gammln(double xx)
/** ${\tt xx} > 0$ のとき $\ln \Gamma({\tt xx})$ を返す。**/
{
    /** 内部での計算は倍精度で行う。
	  5桁精度でよいならこの工夫は不要。**/
	double x,y,tmp,ser;
	static double cof[6]={76.18009172947146,-86.50532032941677,
		24.01409824083091,-1.231739572450155,
		0.1208650973866179e-2,-0.5395239384953e-5};
	int j;

	y=x=xx;
	tmp=x+5.5;
	tmp -= (x+0.5)*log(tmp);
	ser=1.000000000190015;
	for (j=0;j<=5;j++) ser += cof[j]/++y;
	return -tmp+log(2.5066282746310005*ser/x);
}