Function get_normal_pair

contrib/tablefunc/tablefunc.c:282–313

* get_normal_pair() * Assigns normally distributed (Gaussian) values to a pair of provided * parameters, with mean 0, standard deviation 1. * * This routine implements Algorithm P (Polar method for normal deviates) * from Knuth's _The_Art_of_Computer_Programming_, Volume 2, 3rd ed., pages * 122-126. Knuth cites his source as "The polar method", G. E. P. Box, M. E. * Muller, and G. Marsaglia

 *
 */
static void
get_normal_pair(float8 *x1, float8 *x2)
{
	float8		u1,
				u2,
				v1,
				v2,
				s;

	do
	{
		u1 = (float8) random() / (float8) MAX_RANDOM_VALUE;
		u2 = (float8) random() / (float8) MAX_RANDOM_VALUE;

		v1 = (2.0 * u1) - 1.0;
		v2 = (2.0 * u2) - 1.0;

		s = v1 * v1 + v2 * v2;
	} while (s >= 1.0);

	if (s == 0)
	{
		*x1 = 0;
		*x2 = 0;
	}
	else
	{
		s = sqrt((-2.0 * log(s)) / s);
		*x1 = v1 * s;
		*x2 = v2 * s;
	}
}

/*
 * crosstab - create a crosstab of rowids and values columns from a