Merge pull request #148 from joezuntz/des-planck-update

DES + planck and related update

examples/des-y3-planck.ini

@@ -0,0 +1,282 @@
+[runtime]
+sampler = test
+verbosity = standard
+
+; This parameter is used several times in this file, so is
+; put in the DEFAULT section and is referenced below as %(2PT_FILE)s
+[DEFAULT]
+2PT_FILE = likelihood/des-y3/2pt_NG_final_2ptunblind_02_24_21_wnz_covupdate.v2.fits
+
+[pipeline]
+modules =  consistency  bbn_consistency
+           camb_planck  fast_pt
+           fits_nz  lens_photoz_width  lens_photoz_bias  source_photoz_bias
+           IA  pk_to_cl_gg  pk_to_cl
+           add_magnification  add_intrinsic
+           2pt_shear 2pt_gal  2pt_gal_shear
+           shear_m_bias   add_point_mass
+           2pt_like  shear_ratio_like planck_lite
+
+timing=F
+debug=F
+priors = examples/des-y3-priors.ini
+values = examples/des-y3-values.ini
+extra_output = cosmological_parameters/sigma_8 cosmological_parameters/sigma_12 data_vector/2pt_chi2
+
+
+; It's worth switching this to T when sampling using multinest, polychord,
+; or other samplers that can take advantage of differences in calculation speeds between
+; different parameters.
+fast_slow = F
+first_fast_module = shear_m_bias
+; For some use cases this might be faster:
+;first_fast_module=lens_photoz_width
+
+
+; This is the sampler used for the DES-Y3 releases.
+[polychord]
+base_dir = output/y3-polychord-checkpoints
+polychord_outfile_root = y3
+resume = F
+feedback = 3
+fast_fraction = 0.1
+
+;Minimum settings
+live_points = 250
+num_repeats = 30
+tolerance = 0.1
+
+;Settings for paper runs
+; live_points = 500
+; num_repeats=60
+; tolerance=0.01
+; boost_posteriors=10.0
+
+[output]
+filename = output/des-y3.txt
+format=text
+
+; you can set this to T if you don't want 
+privacy = F
+
+[emcee]
+walkers = 80
+samples = 10000
+nsteps = 5
+
+
+; This sampler makes a 1D likelihood slice in each dimension
+; It's handy for doing a thorough comparison to other codes
+[star]
+nsample_dimension = 20
+
+[test]
+save_dir=output/des-y3
+fatal_errors=T
+
+
+[consistency]
+file = utility/consistency/consistency_interface.py
+
+
+[camb_planck]
+file = boltzmann/camb/camb_interface.py
+mode=all
+lmax=4500
+max_eta_k=500000
+feedback=0         ;amount of output to print
+AccuracyBoost=1.1 ;CAMB accuracy boost parameter
+do_tensors=T
+do_lensing=T
+NonLinear = both
+halofit_version = takahashi
+zmin_background = 0.
+zmax_background = 4.
+nz_background = 401
+kmax = 50.0
+kmax_extrapolate = 500.0
+nk=700
+
+
+
+[bbn_consistency]
+file = utility/bbn_consistency/bbn_consistency.py
+
+[fits_nz]
+file = number_density/load_nz_fits/load_nz_fits.py
+nz_file = %(2PT_FILE)s
+data_sets = lens source
+prefix_section = T
+prefix_extension = T
+
+[fits_nz_lens]
+file = number_density/load_nz_fits/load_nz_fits.py
+nz_file = %(2PT_FILE)s
+data_sets = lens
+prefix_section = T
+prefix_extension = T
+
+[lens_photoz_width]
+file = number_density/photoz_width/photoz_width.py
+mode = stretch
+sample = nz_lens
+bias_section = lens_photoz_errors
+interpolation = linear
+
+[lens_photoz_bias]
+file = number_density/photoz_bias/photoz_bias.py
+mode = additive
+sample = nz_lens
+bias_section = lens_photoz_errors
+interpolation = linear
+
+; hyperrank and source_photoz_bias are exclusive
+[hyperrank]
+file = number_density/nz_multirank/nz_multirank.py
+nz_file = %(2PT_FILE)s
+data_set = source
+dimensions = 3
+bin_ranks= 1 2 4 
+
+[source_photoz_bias]
+file = number_density/photoz_bias/photoz_bias.py
+mode = additive
+sample = nz_source
+bias_section = wl_photoz_errors
+interpolation = linear
+
+[fast_pt]
+file = structure/fast_pt/fast_pt_interface.py
+do_ia = T
+k_res_fac = 0.5
+verbose = F
+
+[IA]
+file = intrinsic_alignments/tatt/tatt_interface.py
+sub_lowk = F
+do_galaxy_intrinsic = F
+ia_model = tatt
+
+[pk_to_cl_gg]
+file = structure/projection/project_2d.py
+lingal-lingal = lens-lens
+do_exact = lingal-lingal
+do_rsd = True
+ell_min_linspaced = 1
+ell_max_linspaced = 4
+n_ell_linspaced = 5
+ell_min_logspaced = 5.
+ell_max_logspaced = 5.e5
+n_ell_logspaced = 80
+limber_ell_start = 200
+ell_max_logspaced=1.e5
+auto_only=lingal-lingal
+sig_over_dchi_exact = 3.5
+
+[pk_to_cl]
+file = structure/projection/project_2d.py
+ell_min_logspaced = 0.1
+ell_max_logspaced = 5.0e5
+n_ell_logspaced = 100 
+shear-shear = source-source
+shear-intrinsic = source-source
+intrinsic-intrinsic = source-source
+intrinsicb-intrinsicb=source-source
+lingal-shear = lens-source
+lingal-intrinsic = lens-source
+lingal-magnification = lens-lens
+magnification-shear = lens-source
+magnification-magnification = lens-lens
+magnification-intrinsic = lens-source 
+verbose = F
+get_kernel_peaks = F
+sig_over_dchi = 20. 
+shear_kernel_dchi = 10. 
+
+[add_magnification]
+file = structure/magnification/add_magnification.py
+include_intrinsic=T
+
+[add_intrinsic]
+file=shear/add_intrinsic/add_intrinsic.py
+shear-shear=T
+position-shear=T
+perbin=F
+
+[add_eb]
+file = intrinsic_alignments/add_b_mode/add_b_mode_cl.py
+
+[2pt_shear]
+file = shear/cl_to_xi_fullsky/cl_to_xi_interface.py
+ell_max = 40000
+xi_type = EB
+theta_file=%(2PT_FILE)s
+bin_avg = T
+; these get
+input_section_name = shear_cl  shear_cl_bb
+output_section_name = shear_xi_plus  shear_xi_minus
+
+[2pt_gal]
+file = shear/cl_to_xi_fullsky/cl_to_xi_interface.py
+ell_max = 40000
+xi_type='00'
+theta_file=%(2PT_FILE)s
+bin_avg = T
+
+[2pt_gal_shear]
+file = shear/cl_to_xi_fullsky/cl_to_xi_interface.py
+ell_max = 40000
+xi_type='02'
+theta_file=%(2PT_FILE)s
+bin_avg = T
+
+[shear_m_bias]
+file = shear/shear_bias/shear_m_bias.py
+m_per_bin = True
+; Despite the parameter name, this can operate on xi as well as C_ell.
+cl_section = shear_xi_plus shear_xi_minus
+cross_section = galaxy_shear_xi
+verbose = F
+
+[add_point_mass]
+file=shear/point_mass/add_gammat_point_mass.py
+add_togammat = False
+use_fiducial = True
+sigcrit_inv_section = sigma_crit_inv_lens_source
+
+[2pt_like]
+file = likelihood/2pt/2pt_point_mass/2pt_point_mass.py
+do_pm_marg = True
+do_pm_sigcritinv = True
+sigma_a = 10000.0
+no_det_fac = False
+include_norm = True
+data_file = %(2PT_FILE)s
+data_sets = xip xim gammat wtheta
+make_covariance=F
+covmat_name=COVMAT
+
+; we put these in a separate file because they are long
+%include examples/des-y3-scale-cuts.ini
+
+
+; The small-scale galaxy galaxy-lensing correlations have uncertain
+; enough modelling that we can't use them directly, but the ratio between
+; sets of correlations can be used, since it only depends on geometry
+[shear_ratio_like]
+file = likelihood/des-y3/shear_ratio/shear_ratio_likelihood.py
+data_file = likelihood/des-y3/shear_ratio/2pt_NG_final_2ptunblind_02_24_21_wnz_covupdate_sr.pkl
+theta_min_1 = 9.0  6.0  4.5  2.5  2.5
+theta_min_2 = 9.0  6.0  4.5  2.5  2.5
+theta_min_3 = 2.5  2.5  4.5  2.5  2.5 
+theta_max = 26.83313651 17.63634989 13.61215672 11.32891161 10.01217238
+include_norm = T
+
+
+[planck_lite]
+file = likelihood/planck_py/planck_py_interface.py 
+use_low_ell_bins = T
+spectra = TTTEEE
+year = 2018
+
+

likelihood/des-y3/shear_ratio/shear_ratio_likelihood.py

@@ -6,8 +6,10 @@
 
 def get_ratio_from_gammat(gammat1, gammat2, inv_cov):
     #Given two gammats, calculate the ratio
-    s2 = (1./float(np.matrix(np.ones(len(gammat1)))*np.matrix(inv_cov)*np.matrix(np.ones(len(gammat1))).T))
-    ratio = s2*float(np.matrix(gammat1/gammat2)*np.matrix(inv_cov)*np.matrix(np.ones(len(gammat1))).T)
+    n1 = len(gammat1)
+    unit1 = np.ones(n1)
+    s2 = (1./float(unit1 @ inv_cov @ unit1.T))
+    ratio = s2*float(((gammat1/gammat2) @ inv_cov @ unit1.T))
     return ratio
 
 def radians_to_arcmin(r):

tests/test_cosmosis_standard_library.py

@@ -87,6 +87,12 @@ def test_des_y3(capsys):
     check_likelihood(capsys, "6043.23", "6043.34", "6043.37", "6043.33")
     check_no_camb_warnings(capsys)
 
+def test_des_y3_plus_planck(capsys):
+    run_cosmosis("examples/des-y3-planck.ini")
+    check_likelihood(capsys, "5679.6", "5679.7")
+    check_no_camb_warnings(capsys)
+
+
 def test_des_y3_class(capsys):
     run_cosmosis("examples/des-y3-class.ini")
     # class is not consistent across systems to the level needed here