Template cosmo

py/picca/bin/picca_cf.py

@@ -274,6 +274,13 @@ def main(cmdargs=None):
         help=("Equation of state of dark energy of fiducial " "LambdaCDM cosmology"),
     )
 
+    parser.add_argument(
+        "--distance-template",
+        type=str,
+        required=False,
+        help=("Template cosmology to use in distance calculations. "),
+    )
+
     parser.add_argument(
         "--no-project",
         action="store_true",
@@ -380,6 +387,7 @@ def main(cmdargs=None):
         Ok=args.fid_Ok,
         wl=args.fid_wl,
         blinding=blinding,
+        template=args.distance_template
     )
 
     t0 = time.time()

py/picca/bin/picca_dmat.py

@@ -282,6 +282,13 @@ def main(cmdargs=None):
         help="Equation of state of dark energy of fiducial LambdaCDM cosmology",
     )
 
+    parser.add_argument(
+        "--distance-template",
+        type=str,
+        required=False,
+        help=("Template cosmology to use in distance calculations. "),
+    )
+
     parser.add_argument(
         "--remove-same-half-plate-close-pairs",
         action="store_true",
@@ -403,6 +410,7 @@ def main(cmdargs=None):
         Ok=args.fid_Ok,
         wl=args.fid_wl,
         blinding=blinding,
+        template=args.distance_template
     )
 
     t0 = time.time()

py/picca/bin/picca_xcf.py

@@ -279,6 +279,13 @@ def main(cmdargs=None):
         help="Equation of state of dark energy of fiducial LambdaCDM cosmology",
     )
 
+    parser.add_argument(
+        "--distance-template",
+        type=str,
+        required=False,
+        help=("Template cosmology to use in distance calculations. "),
+    )
+
     parser.add_argument(
         "--no-project",
         action="store_true",
@@ -380,6 +387,7 @@ def main(cmdargs=None):
         Ok=args.fid_Ok,
         wl=args.fid_wl,
         blinding=blinding,
+        template=args.distance_template
     )
 
     t0 = time.time()

py/picca/bin/picca_xdmat.py

@@ -296,6 +296,13 @@ def main(cmdargs=None):
         help="Equation of state of dark energy of fiducial LambdaCDM cosmology",
     )
 
+    parser.add_argument(
+        "--distance-template",
+        type=str,
+        required=False,
+        help=("Template cosmology to use in distance calculations. "),
+    )
+
     parser.add_argument(
         "--rej",
         type=float,
@@ -392,6 +399,7 @@ def main(cmdargs=None):
         Ok=args.fid_Ok,
         wl=args.fid_wl,
         blinding=blinding,
+        template=args.distance_template
     )
 
     t0 = time.time()

py/picca/constants.py

@@ -8,6 +8,7 @@
 from scipy import interpolate
 from scipy.constants import speed_of_light as speed_light
 from picca.utils import userprint
+import camb
 
 # TODO: this constant is unused. They should be removed at some point
 BOSS_LAMBDA_MIN = 3600. # [Angstrom]
@@ -155,7 +156,8 @@ def get_dist_v(self, z):
         """
         raise NotImplementedError("Function should be specified at run-time")
 
-    def __init__(self,Om,Ok=0.,Or=0.,wl=-1.,blinding=False,verbose=True):
+    def __init__(self,Om,Ok=0.,Or=0.,wl=-1.,blinding=False
+                 ,verbose=True, template=None):
         """Initializes the methods for this instance
 
         Args:
@@ -171,6 +173,10 @@ def __init__(self,Om,Ok=0.,Or=0.,wl=-1.,blinding=False,verbose=True):
                 Hubble constant at redshift 0 (in km/s/Mpc)
         """
 
+        #initialise template (if not None)
+        self._use_template = False
+        self._init_template(template)
+
         # WARNING: This is introduced due to historical issues in how this class
         # is coded. Using H0=100 implies that we are returning the distances
         # in Mpc/h instead of Mpc. This class should be fixed at some point to
@@ -188,7 +194,10 @@ def __init__(self,Om,Ok=0.,Or=0.,wl=-1.,blinding=False,verbose=True):
         z_max = 10.
         delta_z = z_max/num_bins
         z = np.arange(num_bins, dtype=float)*delta_z
-        hubble = H0*np.sqrt(Ol*(1. + z)**(3.*(1. + wl)) +
+        if self._use_template:
+            hubble = self._camb_results.hubble_parameter(z) / self._h_camb
+        else:
+            hubble = H0*np.sqrt(Ol*(1. + z)**(3.*(1. + wl)) +
                             Ok*(1. + z)**2 +
                             Om*(1. + z)**3 +
                             Or*(1. + z)**4)
@@ -223,6 +232,13 @@ def __init__(self,Om,Ok=0.,Or=0.,wl=-1.,blinding=False,verbose=True):
                           1./3.)
         self.get_dist_v = interpolate.interp1d(z, dist_v)
 
+    def _init_template(self, template):
+        if template is not None:
+            pars = camb.read_ini(template)
+            self._camb_results = camb.get_results(pars)
+            self._h_camb = pars.H0 / 100
+            self._use_template = True
+
 ### Absorber names and wavelengths [Angstrom]
 ABSORBER_IGM = {
     "Halpha"      : 6562.8,