Move gaussian_noise() to utils.c and synergise it with poisson_noise()

1 files changed, 11 insertions, 12 deletions

diff --git a/src/utils.c b/src/utils.c
index 7e61f250..bfa5c7d4 100644
--- a/src/utils.c
+++ b/src/utils.c
@@ -120,24 +120,23 @@ void progress_bar(int val, int total, const char *text)
 }
 
 
-static int fake_poisson_noise(double expected)
+double gaussian_noise(double expected, double stddev)
 {
-	double x1, x2, w;
-	double noise, rf;
-
-	do {
+	double x1, x2, noise;
 
-		x1 = 2.0 * ((double)random()/RAND_MAX) - 1.0;
-		x2 = 2.0 * ((double)random()/RAND_MAX) - 1.0;
-		w = pow(x1, 2.0) + pow(x2, 2.0);
+	/* A uniformly distributed random number between 0 and 1 */
+	x1 = ((double)random()/RAND_MAX);
+	x2 = ((double)random()/RAND_MAX);
 
-	} while ( w >= 1.0 );
+	noise = sqrt(-2.0*log(x1)) * cos(2.0*M_PI*x2);
 
-	w = sqrt((-2.0*log(w))/w);
-	noise = w * x1;
+	return expected + noise*stddev;
+}
 
-	rf = expected + noise*sqrt(expected);
 
+static int fake_poisson_noise(double expected)
+{
+	double rf = gaussian_noise(expected, sqrt(expected));
 	return (int)rf;
 }