From a3ef55c318e269547eacce97dfea43ca8ff22e65 Mon Sep 17 00:00:00 2001
From: Guilhem Lavaux <guilhem.lavaux@gmail.com>
Date: Fri, 30 Nov 2012 17:16:45 -0500
Subject: [PATCH] Multiply the power spectrum by the volume to have proper
 units

---
 src/fourier/details/euclidian_spectrum_1d.hpp | 9 +++++++--
 1 file changed, 7 insertions(+), 2 deletions(-)

diff --git a/src/fourier/details/euclidian_spectrum_1d.hpp b/src/fourier/details/euclidian_spectrum_1d.hpp
index 291470d..6e2c21e 100644
--- a/src/fourier/details/euclidian_spectrum_1d.hpp
+++ b/src/fourier/details/euclidian_spectrum_1d.hpp
@@ -69,17 +69,22 @@ namespace CosmoTool
     std::complex<T> *d = rand_map.data();
     long idx;
     const DimArray& dims = rand_map.getDims();
+    const std::vector<double>& delta_k = rand_map.get_delta_k();
     long plane_size;
     bool alleven = rand_map.allDimensionsEven();
+    double V = 1;
+
+    for (int p = 0; p < delta_k.size(); p++)
+      V *= (2*M_PI/delta_k[p]);
 
     for (long p = 1; p < rand_map.size(); p++)
       {
-	double A_k = std::sqrt(0.5*f(rand_map.get_K(p)));
+	double A_k = std::sqrt(0.5*V*f(rand_map.get_K(p)));
 	d[p] = std::complex<T>(gsl_ran_gaussian(rng, A_k),
 			       gsl_ran_gaussian(rng, A_k));
       }
     // Generate the mean value
-    d[0] = std::complex<T>(gsl_ran_gaussian(rng, std::sqrt(f(0))), 0);
+    d[0] = std::complex<T>(gsl_ran_gaussian(rng, std::sqrt(V*f(0))), 0);
 
     if (!rand_map.firstDimensionEven())
       return;