diff --git a/examples/des_y1_3x2pt/des_y1_3x2pt.py b/examples/des_y1_3x2pt/des_y1_3x2pt.py index d51e79b24..fb10ffa43 100755 --- a/examples/des_y1_3x2pt/des_y1_3x2pt.py +++ b/examples/des_y1_3x2pt/des_y1_3x2pt.py @@ -7,6 +7,8 @@ import firecrown.likelihood.number_counts as nc from firecrown.likelihood.two_point import TwoPoint from firecrown.likelihood.gaussian import ConstGaussian +from firecrown.modeling_tools import ModelingTools +from firecrown.ccl_factory import CCLFactory # The likelihood used for DES Y1 3x2pt analysis is a Gaussian likelihood, which @@ -111,6 +113,8 @@ def build_likelihood(_): # file and the sources their respective dndz. likelihood.read(sacc_data) + modeling_tools = ModelingTools(ccl_factory=CCLFactory(require_nonlinear_pk=True)) + # This script will be loaded by the appropriated connector. The framework # will call the factory function that should return a Likelihood instance. - return likelihood + return likelihood, modeling_tools diff --git a/examples/des_y1_3x2pt/des_y1_3x2pt_PT.py b/examples/des_y1_3x2pt/des_y1_3x2pt_PT.py index 8e5076394..9431b9886 100755 --- a/examples/des_y1_3x2pt/des_y1_3x2pt_PT.py +++ b/examples/des_y1_3x2pt/des_y1_3x2pt_PT.py @@ -20,6 +20,8 @@ from firecrown.parameters import ParamsMap from firecrown.modeling_tools import ModelingTools from firecrown.likelihood.likelihood import Likelihood +from firecrown.ccl_factory import CCLFactory +from firecrown.updatable import get_default_params_map saccfile = os.path.expanduser( @@ -137,7 +139,9 @@ def build_likelihood(_) -> tuple[Likelihood, ModelingTools]: nk_per_decade=20, ) - modeling_tools = ModelingTools(pt_calculator=pt_calculator) + modeling_tools = ModelingTools( + pt_calculator=pt_calculator, ccl_factory=CCLFactory(require_nonlinear_pk=True) + ) # Note that the ordering of the statistics is relevant, because the data # vector, theory vector, and covariance matrix will be organized to # follow the order used here. @@ -173,9 +177,11 @@ def run_likelihood() -> None: z, nz = src0_tracer.z, src0_tracer.nz lens_z, lens_nz = lens0_tracer.z, lens0_tracer.nz - # Define a ccl.Cosmology object using default parameters - ccl_cosmo = ccl.CosmologyVanillaLCDM() - ccl_cosmo.compute_nonlin_power() + # Setting cosmological parameters + cosmo_params = get_default_params_map(tools) + tools.update(cosmo_params) + tools.prepare() + ccl_cosmo = tools.get_ccl_cosmology() cs = CclSetup() c_1, c_d, c_2 = pyccl.nl_pt.translate_IA_norm( @@ -223,10 +229,6 @@ def run_likelihood() -> None: # Apply the systematics parameters likelihood.update(systematics_params) - tools.update(systematics_params) - - # Prepare the cosmology object - tools.prepare(ccl_cosmo) # Compute the log-likelihood, using the ccl.Cosmology object as the input log_like = likelihood.compute_loglike(tools) diff --git a/examples/des_y1_3x2pt/des_y1_3x2pt_default_factory.ini b/examples/des_y1_3x2pt/des_y1_3x2pt_default_factory.ini new file mode 100644 index 000000000..8f6a211c1 --- /dev/null +++ b/examples/des_y1_3x2pt/des_y1_3x2pt_default_factory.ini @@ -0,0 +1,60 @@ +[runtime] +sampler = test +root = ${PWD} + +[DEFAULT] +fatal_errors = T + +[output] +filename = output/des_y1_3x2pt_samples.txt +format = text +verbosity = 0 + +[pipeline] +modules = consistency camb firecrown_likelihood +values = ${FIRECROWN_DIR}/examples/des_y1_3x2pt/des_y1_3x2pt_values.ini +likelihoods = firecrown +quiet = T +debug = T +timing = T +extra_output = TwoPoint/NumberCountsScale_lens0 TwoPoint/NumberCountsScale_lens1 TwoPoint/NumberCountsScale_lens2 TwoPoint/NumberCountsScale_lens3 TwoPoint/NumberCountsScale_lens4 + +[consistency] +file = ${CSL_DIR}/utility/consistency/consistency_interface.py + +[camb] +file = ${CSL_DIR}/boltzmann/camb/camb_interface.py + +mode = all +lmax = 2500 +feedback = 0 +zmin = 0.0 +zmax = 4.0 +nz = 100 +kmin = 1e-4 +kmax = 50.0 +nk = 1000 + +[firecrown_likelihood] +;; Fix this to use an environment variable to find the files. +;; Set FIRECROWN_DIR to the base of the firecrown installation (or build, if you haven't +;; installed it) +file = ${FIRECROWN_DIR}/firecrown/connector/cosmosis/likelihood.py +likelihood_source = firecrown.likelihood.factories.build_two_point_likelihood +likelihood_config = ${FIRECROWN_DIR}/examples/des_y1_3x2pt/des_y1_3x2pt_experiment.yaml +;; Connector settings +require_nonlinear_pk = True +sampling_parameters_sections = firecrown_two_point + +[test] +fatal_errors = T +save_dir = des_y1_3x2pt_output + +[metropolis] +samples = 1000 +nsteps = 1 + +[emcee] +walkers = 64 +samples = 400 +nsteps = 10 diff --git a/examples/des_y1_3x2pt/des_y1_3x2pt_experiment.yaml b/examples/des_y1_3x2pt/des_y1_3x2pt_experiment.yaml new file mode 100644 index 000000000..298317a39 --- /dev/null +++ b/examples/des_y1_3x2pt/des_y1_3x2pt_experiment.yaml @@ -0,0 +1,29 @@ +--- + +# This file contains all information about the DES Y1 3x2pt experiment. + +# The 'data_source:' field points to the SACC file containing the DES Y1 3x2pt data vector +# and covariance matrix, which is used for analysis. + +# The 'two_point_factory:' field points to the factory that will create the TwoPoint +# objects. These objects represent the chosen theoretical models for each data point +# found in the SACC file. + +data_source: + sacc_data_file: des_y1_3x2pt_sacc_data.fits + +# The two point statistics are defined by the TwoPoint objects. The TwoPoint statistics +# are created using the factories defined in this file. +two_point_factory: + correlation_space: real + number_counts_factory: + global_systematics: [] + per_bin_systematics: + - type: PhotoZShiftFactory + weak_lensing_factory: + global_systematics: + - alphag: 1 + type: LinearAlignmentSystematicFactory + per_bin_systematics: + - type: MultiplicativeShearBiasFactory + - type: PhotoZShiftFactory diff --git a/examples/des_y1_3x2pt/des_y1_cosmic_shear_TATT.py b/examples/des_y1_3x2pt/des_y1_cosmic_shear_TATT.py index 640726aae..647268bdc 100644 --- a/examples/des_y1_3x2pt/des_y1_cosmic_shear_TATT.py +++ b/examples/des_y1_3x2pt/des_y1_cosmic_shear_TATT.py @@ -11,7 +11,9 @@ from firecrown.parameters import ParamsMap from firecrown.modeling_tools import ModelingTools from firecrown.likelihood.likelihood import Likelihood - +from firecrown.ccl_factory import CCLFactory +from firecrown.updatable import get_default_params_map +from firecrown.metadata_types import TracerNames, TRACER_NAMES_TOTAL SACCFILE = os.path.expanduser( os.path.expandvars( @@ -38,7 +40,9 @@ def build_likelihood(_) -> tuple[Likelihood, ModelingTools]: nk_per_decade=20, ) - modeling_tools = ModelingTools(pt_calculator=pt_calculator) + modeling_tools = ModelingTools( + pt_calculator=pt_calculator, ccl_factory=CCLFactory(require_nonlinear_pk=True) + ) likelihood = ConstGaussian(statistics=list(stats.values())) # Read the two-point data from the sacc file @@ -102,14 +106,33 @@ def run_likelihood() -> None: src0_tracer = sacc_data.get_tracer("src0") z, nz = src0_tracer.z, src0_tracer.nz - # Define a ccl.Cosmology object using default parameters - ccl_cosmo = ccl.CosmologyVanillaLCDM() - ccl_cosmo.compute_nonlin_power() - # Bare CCL setup a_1 = 1.0 a_2 = 0.5 a_d = 0.5 + # Set the parameters for our systematics + systematics_params = ParamsMap( + { + "ia_a_1": a_1, + "ia_a_2": a_2, + "ia_a_d": a_d, + "src0_delta_z": 0.000, + "src1_delta_z": 0.003, + "src2_delta_z": -0.001, + "src3_delta_z": 0.002, + } + ) + # Prepare the cosmology object + params = ParamsMap(get_default_params_map(tools) | systematics_params) + + # Apply the systematics parameters + likelihood.update(params) + + # Prepare the cosmology object + tools.update(params) + tools.prepare() + ccl_cosmo = tools.get_ccl_cosmology() + c_1, c_d, c_2 = pyccl.nl_pt.translate_IA_norm( ccl_cosmo, z=z, a1=a_1, a1delta=a_d, a2=a_2, Om_m2_for_c2=False ) @@ -131,25 +154,6 @@ def run_likelihood() -> None: pk_im = ptc.get_biased_pk2d(tracer1=ptt_i, tracer2=ptt_m) pk_ii = ptc.get_biased_pk2d(tracer1=ptt_i, tracer2=ptt_i) - # Set the parameters for our systematics - systematics_params = ParamsMap( - { - "ia_a_1": a_1, - "ia_a_2": a_2, - "ia_a_d": a_d, - "src0_delta_z": 0.000, - "src1_delta_z": 0.003, - "src2_delta_z": -0.001, - "src3_delta_z": 0.002, - } - ) - - # Apply the systematics parameters - likelihood.update(systematics_params) - - # Prepare the cosmology object - tools.prepare(ccl_cosmo) - # Compute the log-likelihood, using the ccl.Cosmology object as the input log_like = likelihood.compute_loglike(tools) @@ -180,11 +184,11 @@ def make_plot(ccl_cosmo, nz, pk_ii, pk_im, two_point_0, z): import numpy as np # pylint: disable-msg=import-outside-toplevel import matplotlib.pyplot as plt # pylint: disable-msg=import-outside-toplevel - ells = two_point_0.ells - cells_gg = two_point_0.cells[("shear", "shear")] - cells_gi = two_point_0.cells[("shear", "intrinsic_pt")] - cells_ii = two_point_0.cells[("intrinsic_pt", "intrinsic_pt")] - cells_total = two_point_0.cells["total"] + ells = two_point_0.ells_for_xi + cells_gg = two_point_0.cells[TracerNames("shear", "shear")] + cells_gi = two_point_0.cells[TracerNames("shear", "intrinsic_pt")] + cells_ii = two_point_0.cells[TracerNames("intrinsic_pt", "intrinsic_pt")] + cells_total = two_point_0.cells[TRACER_NAMES_TOTAL] # pylint: enable=no-member # Code that computes effect from IA using that Pk2D object t_lens = ccl.WeakLensingTracer(ccl_cosmo, dndz=(z, nz)) diff --git a/examples/des_y1_3x2pt/des_y1_cosmic_shear_pk_modifier.py b/examples/des_y1_3x2pt/des_y1_cosmic_shear_pk_modifier.py index 6166b1fcf..40620eb6e 100644 --- a/examples/des_y1_3x2pt/des_y1_cosmic_shear_pk_modifier.py +++ b/examples/des_y1_3x2pt/des_y1_cosmic_shear_pk_modifier.py @@ -12,10 +12,12 @@ import firecrown.likelihood.weak_lensing as wl from firecrown.likelihood.two_point import TwoPoint from firecrown.likelihood.gaussian import ConstGaussian -from firecrown.parameters import ParamsMap, create +from firecrown.parameters import ParamsMap, register_new_updatable_parameter from firecrown.modeling_tools import ModelingTools, PowerspectrumModifier from firecrown.likelihood.likelihood import Likelihood - +from firecrown.ccl_factory import CCLFactory +from firecrown.updatable import get_default_params_map +from firecrown.metadata_types import TracerNames SACCFILE = os.path.expanduser( os.path.expandvars( @@ -36,7 +38,7 @@ def __init__(self, pk_to_modify: str = "delta_matter:delta_matter"): super().__init__() self.pk_to_modify = pk_to_modify self.vD19 = pyccl.baryons.BaryonsvanDaalen19() - self.f_bar = create() + self.f_bar = register_new_updatable_parameter(default_value=0.5) def compute_p_of_k_z(self, tools: ModelingTools) -> pyccl.Pk2D: """Compute the 3D power spectrum P(k, z).""" @@ -56,7 +58,9 @@ def build_likelihood(_) -> tuple[Likelihood, ModelingTools]: # Define the power spectrum modification and add it to the ModelingTools pk_modifier = vanDaalen19Baryonfication(pk_to_modify="delta_matter:delta_matter") - modeling_tools = ModelingTools(pk_modifiers=[pk_modifier]) + modeling_tools = ModelingTools( + pk_modifiers=[pk_modifier], ccl_factory=CCLFactory(require_nonlinear_pk=True) + ) # Create the likelihood from the statistics likelihood = ConstGaussian(statistics=list(stats.values())) @@ -123,18 +127,7 @@ def run_likelihood() -> None: src0_tracer = sacc_data.get_tracer("src0") z, nz = src0_tracer.z, src0_tracer.nz - # Define a ccl.Cosmology object using default parameters - ccl_cosmo = ccl.CosmologyVanillaLCDM() - ccl_cosmo.compute_nonlin_power() - f_bar = 0.5 - - # Calculate the barynic effects directly with CCL - vD19 = pyccl.BaryonsvanDaalen19(fbar=f_bar) - pk_baryons = vD19.include_baryonic_effects( - cosmo=ccl_cosmo, pk=ccl_cosmo.get_nonlin_power() - ) - # Set the parameters for our systematics systematics_params = ParamsMap( { @@ -145,13 +138,22 @@ def run_likelihood() -> None: "src3_delta_z": 0.002, } ) + # Prepare the cosmology object + params = ParamsMap(get_default_params_map(tools) | systematics_params) - # Apply the systematics parameters - likelihood.update(systematics_params) + tools.update(params) + tools.prepare() - # Prepare the cosmology object - tools.update(systematics_params) - tools.prepare(ccl_cosmo) + ccl_cosmo = tools.get_ccl_cosmology() + + # Calculate the barynic effects directly with CCL + vD19 = pyccl.BaryonsvanDaalen19(fbar=f_bar) + pk_baryons = vD19.include_baryonic_effects( + cosmo=ccl_cosmo, pk=ccl_cosmo.get_nonlin_power() + ) + + # Apply the systematics parameters + likelihood.update(params) # Compute the log-likelihood, using the ccl.Cosmology object as the input log_like = likelihood.compute_loglike(tools) @@ -166,7 +168,7 @@ def run_likelihood() -> None: # Predict CCL Cl wl_tracer = ccl.WeakLensingTracer(ccl_cosmo, dndz=(z, nz)) - ell = two_point_0.ells + ell = two_point_0.ells_for_xi cl_dm = ccl.angular_cl( cosmo=ccl_cosmo, tracer1=wl_tracer, @@ -191,7 +193,7 @@ def make_plot(ell, cl_dm, cl_baryons, two_point_0): """Create and show a diagnostic plot.""" import matplotlib.pyplot as plt # pylint: disable-msg=import-outside-toplevel - cl_firecrown = two_point_0.cells[("shear", "shear")] + cl_firecrown = two_point_0.cells[TracerNames("shear", "shear")] plt.plot(ell, cl_firecrown / cl_dm, label="firecrown w/ baryons") plt.plot(ell, cl_baryons / cl_dm, ls="--", label="CCL w/ baryons") diff --git a/firecrown/ccl_factory.py b/firecrown/ccl_factory.py new file mode 100644 index 000000000..14f03e003 --- /dev/null +++ b/firecrown/ccl_factory.py @@ -0,0 +1,233 @@ +"""This module contains the CCLFactory class. + +The CCLFactory class is a factory class that creates instances of the +`pyccl.Cosmology` class. +""" + +from typing import Annotated +from enum import Enum, auto + +# To be moved to the import from typing when migrating to Python 3.11 +from typing_extensions import NotRequired, TypedDict + +import numpy as np +import numpy.typing as npt +from pydantic import ( + BaseModel, + ConfigDict, + BeforeValidator, + SerializerFunctionWrapHandler, + SerializationInfo, + Field, + field_serializer, + model_serializer, +) + +import pyccl +from pyccl.neutrinos import NeutrinoMassSplits + +from firecrown.updatable import Updatable +from firecrown.parameters import register_new_updatable_parameter +from firecrown.utils import YAMLSerializable + +PowerSpec = TypedDict( + "PowerSpec", + { + "a": npt.NDArray[np.float64], + "k": npt.NDArray[np.float64], + "delta_matter:delta_matter": npt.NDArray[np.float64], + }, +) + +Background = TypedDict( + "Background", + { + "a": npt.NDArray[np.float64], + "chi": npt.NDArray[np.float64], + "h_over_h0": npt.NDArray[np.float64], + }, +) + +CCLCalculatorArgs = TypedDict( + "CCLCalculatorArgs", + { + "background": Background, + "pk_linear": NotRequired[PowerSpec], + "pk_nonlin": NotRequired[PowerSpec], + }, +) + + +class PoweSpecAmplitudeParameter(YAMLSerializable, str, Enum): + """This class defines the two-point correlation space. + + The two-point correlation space can be either real or harmonic. The real space + corresponds measurements in terms of angular separation, while the harmonic space + corresponds to measurements in terms of spherical harmonics decomposition. + """ + + @staticmethod + def _generate_next_value_(name, _start, _count, _last_values): + return name.lower() + + AS = auto() + SIGMA8 = auto() + + +def _validade_amplitude_parameter(value): + if isinstance(value, str): + try: + return PoweSpecAmplitudeParameter( + value.lower() + ) # Convert from string to Enum + except ValueError as exc: + raise ValueError( + f"Invalid value for PoweSpecAmplitudeParameter: {value}" + ) from exc + return value + + +def _validate_neutrino_mass_splits(value): + if isinstance(value, str): + try: + return NeutrinoMassSplits(value.lower()) # Convert from string to Enum + except ValueError as exc: + raise ValueError(f"Invalid value for NeutrinoMassSplits: {value}") from exc + return value + + +class CCLFactory(Updatable, BaseModel): + """Factory class for creating instances of the `pyccl.Cosmology` class.""" + + model_config = ConfigDict(extra="allow") + require_nonlinear_pk: Annotated[bool, Field(frozen=True)] = False + amplitude_parameter: Annotated[ + PoweSpecAmplitudeParameter, + BeforeValidator(_validade_amplitude_parameter), + Field(frozen=True), + ] = PoweSpecAmplitudeParameter.SIGMA8 + mass_split: Annotated[ + NeutrinoMassSplits, + BeforeValidator(_validate_neutrino_mass_splits), + Field(frozen=True), + ] = NeutrinoMassSplits.NORMAL + + def __init__(self, **data): + """Initialize the CCLFactory object.""" + parameter_prefix = parameter_prefix = data.pop("parameter_prefix", None) + BaseModel.__init__(self, **data) + Updatable.__init__(self, parameter_prefix=parameter_prefix) + + self._ccl_cosmo: None | pyccl.Cosmology = None + + ccl_cosmo = pyccl.CosmologyVanillaLCDM() + + self.Omega_c = register_new_updatable_parameter( + default_value=ccl_cosmo["Omega_c"] + ) + self.Omega_b = register_new_updatable_parameter( + default_value=ccl_cosmo["Omega_b"] + ) + self.h = register_new_updatable_parameter(default_value=ccl_cosmo["h"]) + self.n_s = register_new_updatable_parameter(default_value=ccl_cosmo["n_s"]) + self.Omega_k = register_new_updatable_parameter( + default_value=ccl_cosmo["Omega_k"] + ) + self.Neff = register_new_updatable_parameter(default_value=ccl_cosmo["Neff"]) + self.m_nu = register_new_updatable_parameter(default_value=ccl_cosmo["m_nu"]) + self.w0 = register_new_updatable_parameter(default_value=ccl_cosmo["w0"]) + self.wa = register_new_updatable_parameter(default_value=ccl_cosmo["wa"]) + self.T_CMB = register_new_updatable_parameter(default_value=ccl_cosmo["T_CMB"]) + + match self.amplitude_parameter: + case PoweSpecAmplitudeParameter.AS: + # VanillaLCDM has does not have A_s, so we need to add it + self.A_s = register_new_updatable_parameter(default_value=2.1e-9) + case PoweSpecAmplitudeParameter.SIGMA8: + assert ccl_cosmo["sigma8"] is not None + self.sigma8 = register_new_updatable_parameter( + default_value=ccl_cosmo["sigma8"] + ) + + @model_serializer(mode="wrap") + def serialize_model(self, nxt: SerializerFunctionWrapHandler, _: SerializationInfo): + """Serialize the CCLFactory object.""" + model_dump = nxt(self) + exclude_params = [param.name for param in self._sampler_parameters] + list( + self._internal_parameters.keys() + ) + + return {k: v for k, v in model_dump.items() if k not in exclude_params} + + @field_serializer("amplitude_parameter") + @classmethod + def serialize_amplitude_parameter(cls, value: PoweSpecAmplitudeParameter) -> str: + """Serialize the amplitude parameter.""" + return value.name + + @field_serializer("mass_split") + @classmethod + def serialize_mass_split(cls, value: NeutrinoMassSplits) -> str: + """Serialize the mass split parameter.""" + return value.name + + def model_post_init(self, __context) -> None: + """Initialize the WeakLensingFactory object.""" + + def create( + self, calculator_args: CCLCalculatorArgs | None = None + ) -> pyccl.Cosmology: + """Create a `pyccl.Cosmology` object.""" + if not self.is_updated(): + raise ValueError("Parameters have not been updated yet.") + + if self._ccl_cosmo is not None: + raise ValueError("CCLFactory object has already been created.") + + # pylint: disable=duplicate-code + ccl_args = { + "Omega_c": self.Omega_c, + "Omega_b": self.Omega_b, + "h": self.h, + "n_s": self.n_s, + "Omega_k": self.Omega_k, + "Neff": self.Neff, + "m_nu": self.m_nu, + "w0": self.w0, + "wa": self.wa, + "T_CMB": self.T_CMB, + "mass_split": self.mass_split.value, + } + # pylint: enable=duplicate-code + match self.amplitude_parameter: + case PoweSpecAmplitudeParameter.AS: + ccl_args["A_s"] = self.A_s + case PoweSpecAmplitudeParameter.SIGMA8: + ccl_args["sigma8"] = self.sigma8 + + assert ("A_s" in ccl_args) or ("sigma8" in ccl_args) + + if calculator_args is not None: + ccl_args.update(calculator_args) + if ("pk_nonlin" not in ccl_args) and self.require_nonlinear_pk: + ccl_args["nonlinear_model"] = "halofit" + else: + ccl_args["nonlinear_model"] = None + + return pyccl.CosmologyCalculator(**ccl_args) + + if self.require_nonlinear_pk: + ccl_args["matter_power_spectrum"] = "halofit" + + self._ccl_cosmo = pyccl.Cosmology(**ccl_args) + return self._ccl_cosmo + + def _reset(self) -> None: + """Reset the CCLFactory object.""" + self._ccl_cosmo = None + + def get(self) -> pyccl.Cosmology: + """Return the `pyccl.Cosmology` object.""" + if self._ccl_cosmo is None: + raise ValueError("CCLFactory object has not been created yet.") + return self._ccl_cosmo diff --git a/firecrown/connector/cobaya/ccl.py b/firecrown/connector/cobaya/ccl.py index 34a9e7668..7d0f6cc1b 100644 --- a/firecrown/connector/cobaya/ccl.py +++ b/firecrown/connector/cobaya/ccl.py @@ -7,11 +7,11 @@ import numpy as np import numpy.typing as npt -import pyccl from cobaya.theory import Theory from firecrown.connector.mapping import mapping_builder, MappingCAMB +from firecrown.ccl_factory import CCLCalculatorArgs class CCLConnector(Theory): @@ -110,9 +110,7 @@ def must_provide(self, **requirements) -> None: This version does nothing. """ - def calculate( - self, state: dict[str, float], want_derived=True, **params_values - ) -> None: + def calculate(self, state: dict, want_derived=True, **params_values) -> None: """Calculate the current cosmology, and set state["pyccl"] to the result. :param state: The state dictionary to update. @@ -120,10 +118,9 @@ def calculate( :param params_values: The values of the parameters to use. """ self.map.set_params_from_camb(**params_values) - pyccl_params_values = self.map.asdict() - # This is the dictionary appropriate for CCL creation + # This is the dictionary appropriate for CCL creation chi_arr = self.provider.get_comoving_radial_distance(self.z_bg) hoh0_arr = self.provider.get_Hubble(self.z_bg) / self.map.get_H0() k, z, pk = self.provider.get_Pk_grid() @@ -135,21 +132,23 @@ def calculate( self.a_Pk = self.map.redshift_to_scale_factor(z) pk_a = self.map.redshift_to_scale_factor_p_k(pk) - cosmo = pyccl.CosmologyCalculator( - **pyccl_params_values, - background={"a": self.a_bg, "chi": chi_arr, "h_over_h0": hoh0_arr}, - pk_linear={ - "a": self.a_Pk, - "k": k, - "delta_matter:delta_matter": pk_a, - }, - nonlinear_model="halofit", - ) - state["pyccl"] = cosmo - - def get_pyccl(self) -> pyccl.Cosmology: - """Return the current cosmology. - - :return: The current cosmology + pyccl_args: CCLCalculatorArgs = { + "background": {"a": self.a_bg, "chi": chi_arr, "h_over_h0": hoh0_arr}, + "pk_linear": {"a": self.a_Pk, "k": k, "delta_matter:delta_matter": pk_a}, + } + state["pyccl_args"] = pyccl_args + state["pyccl_params"] = pyccl_params_values + + def get_pyccl_args(self) -> CCLCalculatorArgs: + """Return the current CCL arguments. + + :return: The current CCL arguments + """ + return self.current_state["pyccl_args"] + + def get_pyccl_params(self) -> dict[str, float]: + """Return the current cosmological parameters. + + :return: The current cosmological parameters. """ - return self.current_state["pyccl"] + return self.current_state["pyccl_params"] diff --git a/firecrown/connector/cobaya/likelihood.py b/firecrown/connector/cobaya/likelihood.py index 16c05d985..54dcc3042 100644 --- a/firecrown/connector/cobaya/likelihood.py +++ b/firecrown/connector/cobaya/likelihood.py @@ -12,6 +12,7 @@ from firecrown.likelihood.likelihood import load_likelihood, NamedParameters from firecrown.parameters import ParamsMap +from firecrown.ccl_factory import PoweSpecAmplitudeParameter class LikelihoodConnector(Likelihood): @@ -81,21 +82,29 @@ def get_requirements( """Returns a dictionary. Returns a dictionary with keys corresponding the contained likelihood's - required parameter, plus "pyccl". All values are None. + required parameter, plus "pyccl_args" and "pyccl_params". All values are None. Required by Cobaya. :return: a dictionary """ likelihood_requires: dict[ str, None | dict[str, npt.NDArray[np.float64]] | dict[str, object] - ] = {"pyccl": None} + ] = {"pyccl_args": None, "pyccl_params": None} required_params = ( self.likelihood.required_parameters() + self.tools.required_parameters() ) + # Cosmological parameters differ from Cobaya's, so we need to remove them. + required_params -= self.tools.ccl_factory.required_parameters() for param_name in required_params.get_params_names(): likelihood_requires[param_name] = None + if ( + self.tools.ccl_factory.amplitude_parameter + == PoweSpecAmplitudeParameter.SIGMA8 + ): + likelihood_requires["sigma8"] = None + return likelihood_requires def must_provide(self, **requirements) -> None: @@ -110,18 +119,21 @@ def logp(self, **params_values) -> float: Required by Cobaya. :params values: The values of the parameters to use. """ - pyccl = self.provider.get_pyccl() + pyccl_args = self.provider.get_pyccl_args() + pyccl_params = self.provider.get_pyccl_params() - self.likelihood.update(ParamsMap(params_values)) - self.tools.update(ParamsMap(params_values)) - self.tools.prepare(pyccl) + derived = params_values.pop("_derived", {}) + params = ParamsMap(params_values | pyccl_params) + self.likelihood.update(params) + self.tools.update(params) + self.tools.prepare(calculator_args=pyccl_args) loglike = self.likelihood.compute_loglike(self.tools) derived_params_collection = self.likelihood.get_derived_parameters() assert derived_params_collection is not None for section, name, val in derived_params_collection: - params_values["_derived"][f"{section}__{name}"] = val + derived[f"{section}__{name}"] = val self.likelihood.reset() self.tools.reset() diff --git a/firecrown/connector/cosmosis/likelihood.py b/firecrown/connector/cosmosis/likelihood.py index c7446df3c..9163ceafa 100644 --- a/firecrown/connector/cosmosis/likelihood.py +++ b/firecrown/connector/cosmosis/likelihood.py @@ -11,7 +11,6 @@ import cosmosis.datablock from cosmosis.datablock import option_section from cosmosis.datablock import names as section_names -import pyccl as ccl from firecrown.connector.mapping import mapping_builder, MappingCosmoSIS from firecrown.likelihood.gaussfamily import GaussFamily @@ -52,17 +51,16 @@ def __init__(self, config: cosmosis.datablock) -> None: if likelihood_source == "": likelihood_source = config[option_section, "firecrown_config"] - require_nonlinear_pk = config.get_bool( - option_section, "require_nonlinear_pk", False - ) - build_parameters = extract_section(config, option_section) - sections = config.get_string(option_section, "sampling_parameters_sections", "") - sections = sections.split() + sections_str: str = config.get_string( + option_section, "sampling_parameters_sections", "" + ) + assert isinstance(sections_str, str) + sections = sections_str.split() self.firecrown_module_name = option_section - self.sampling_sections = sections + self.sampling_sections: list[str] = sections self.likelihood: Likelihood try: self.likelihood, self.tools = load_likelihood( @@ -75,9 +73,7 @@ def __init__(self, config: cosmosis.datablock) -> None: raise # We have to do some extra type-fiddling here because mapping_builder # has a declared return type of the base class. - new_mapping = mapping_builder( - input_style="CosmoSIS", require_nonlinear_pk=require_nonlinear_pk - ) + new_mapping = mapping_builder(input_style="CosmoSIS") assert isinstance(new_mapping, MappingCosmoSIS) self.map = new_mapping @@ -88,6 +84,8 @@ def __init__(self, config: cosmosis.datablock) -> None: required_parameters = ( self.likelihood.required_parameters() + self.tools.required_parameters() ) + required_parameters -= self.tools.ccl_factory.required_parameters() + if len(required_parameters) != 0: msg = ( f"The configured likelihood has required " @@ -111,10 +109,7 @@ def execute(self, sample: cosmosis.datablock) -> int: sample, "cosmological_parameters" ) self.map.set_params_from_cosmosis(cosmological_params) - - ccl_cosmo = ccl.CosmologyCalculator( - **self.map.asdict(), **self.map.calculate_ccl_args(sample) - ) + ccl_args = self.map.calculate_ccl_args(sample) # TODO: Future development will need to capture elements that get put into the # datablock. This probably will be in a different "physics module" and not in @@ -122,9 +117,11 @@ def execute(self, sample: cosmosis.datablock) -> int: # calculations. e.g., data_vector/firecrown_theory data_vector/firecrown_data firecrown_params = self.calculate_firecrown_params(sample) + firecrown_params = ParamsMap(firecrown_params | self.map.asdict()) + firecrown_params.use_lower_case_keys(True) self.update_likelihood_and_tools(firecrown_params) - self.tools.prepare(ccl_cosmo) + self.tools.prepare(calculator_args=ccl_args) loglike = self.likelihood.compute_loglike(self.tools) derived_params_collection = self.likelihood.get_derived_parameters() diff --git a/firecrown/connector/mapping.py b/firecrown/connector/mapping.py index ab1a61e0f..dc2ff5a3f 100644 --- a/firecrown/connector/mapping.py +++ b/firecrown/connector/mapping.py @@ -9,7 +9,7 @@ import typing import warnings from abc import ABC -from typing import Any, Type, final +from typing import Type, final import cosmosis.datablock import numpy as np @@ -18,6 +18,7 @@ from firecrown.descriptors import TypeFloat, TypeString from firecrown.likelihood.likelihood import NamedParameters +from firecrown.ccl_factory import CCLCalculatorArgs, PowerSpec, Background def build_ccl_background_dict( @@ -25,7 +26,7 @@ def build_ccl_background_dict( a: npt.NDArray[np.float64], chi: npt.NDArray[np.float64], h_over_h0: npt.NDArray[np.float64], -) -> dict[str, npt.NDArray[np.float64]]: +) -> Background: """Builds the CCL dictionary of background quantities. :param a: The scale factor array @@ -66,13 +67,8 @@ class Mapping(ABC): wa = TypeFloat() T_CMB = TypeFloat() - def __init__(self, *, require_nonlinear_pk: bool = False) -> None: - """Initialize the Mapping object. - - :param require_nonlinear_pk: Whether the mapping requires the - non-linear power spectrum - """ - self.require_nonlinear_pk = require_nonlinear_pk + def __init__(self) -> None: + """Initialize the Mapping object.""" self.m_nu: float | list[float] | None = None def get_params_names(self) -> list[str]: @@ -139,7 +135,6 @@ def set_params( Omega_k: float, Neff: float, m_nu: float | list[float], - m_nu_type: str, w0: float, wa: float, T_CMB: float, @@ -158,7 +153,6 @@ def set_params( :param Omega_k: curvature of the universe :param Neff: effective number of relativistic neutrino species :param m_nu: effective mass of neutrinos - :param m_nu_type: type of massive neutrinos :param w0: constant of the CPL parameterization of the dark energy equation of state :param wa: linear coefficient of the CPL parameterization of the @@ -186,7 +180,6 @@ def set_params( self.Omega_g = None self.Neff = Neff self.m_nu = m_nu - self.m_nu_type = m_nu_type self.w0 = w0 self.wa = wa self.T_CMB = T_CMB @@ -220,7 +213,7 @@ def redshift_to_scale_factor_p_k( p_k_out = np.flipud(p_k) return p_k_out - def asdict(self) -> dict[str, None | float | list[float] | str]: + def asdict(self) -> dict[str, float | list[float]]: """Return a dictionary containing the cosmological constants. :return: the dictionary, containing keys: @@ -234,29 +227,34 @@ def asdict(self) -> dict[str, None | float | list[float] | str]: - ``Omega_k``: curvature of the universe - ``Neff``: effective number of relativistic neutrino species - ``m_nu``: effective mass of neutrinos - - ``m_nu_type``: type of massive neutrinos - ``w0``: constant of the CPL parameterization of the dark energy equation of state - ``wa``: linear coefficient of the CPL parameterization of the dark energy equation of state - ``T_CMB``: cosmic microwave background temperature today """ - return { + cosmo_dict: dict[str, float | list[float]] = { "Omega_c": self.Omega_c, "Omega_b": self.Omega_b, "h": self.h, - "A_s": self.A_s, - "sigma8": self.sigma8, "n_s": self.n_s, "Omega_k": self.Omega_k, - "Omega_g": self.Omega_g, "Neff": self.Neff, - "m_nu": self.m_nu, - "mass_split": self.m_nu_type, "w0": self.w0, "wa": self.wa, "T_CMB": self.T_CMB, } + if self.A_s is not None: + cosmo_dict["A_s"] = self.A_s + if self.sigma8 is not None: + cosmo_dict["sigma8"] = self.sigma8 + # Currently we do not support Omega_g + # if self.Omega_g is not None: + # cosmo_dict["Omega_g"] = self.Omega_g + if self.m_nu is not None: + cosmo_dict["m_nu"] = self.m_nu + + return cosmo_dict def get_H0(self) -> float: """Return the value of H0. @@ -346,7 +344,6 @@ def set_params_from_cosmosis(self, cosmosis_params: NamedParameters) -> None: Neff = delta_neff + 3.046 omega_nu = cosmosis_params.get_float("omega_nu") m_nu = omega_nu * h * h * 93.14 - m_nu_type = "normal" w0 = cosmosis_params.get_float("w") wa = cosmosis_params.get_float("wa") @@ -360,7 +357,6 @@ def set_params_from_cosmosis(self, cosmosis_params: NamedParameters) -> None: Omega_k=Omega_k, Neff=Neff, m_nu=m_nu, - m_nu_type=m_nu_type, w0=w0, wa=-wa, # Is this minus sign here correct? T_CMB=2.7255, @@ -368,13 +364,14 @@ def set_params_from_cosmosis(self, cosmosis_params: NamedParameters) -> None: ) # pylint: enable=duplicate-code - def calculate_ccl_args(self, sample: cosmosis.datablock) -> dict[str, Any]: + def calculate_ccl_args(self, sample: cosmosis.datablock) -> CCLCalculatorArgs: """Calculate the arguments necessary for CCL for this sample. :param sample: the datablock for the current sample :return: the arguments required by CCL """ - ccl_args: dict[str, Any] = {} + pk_linear: None | PowerSpec = None + pk_nonlin: None | PowerSpec = None if sample.has_section("matter_power_lin"): k = self.transform_k_h_to_k(sample["matter_power_lin", "k_h"]) z_mpl = sample["matter_power_lin", "z"] @@ -382,7 +379,7 @@ def calculate_ccl_args(self, sample: cosmosis.datablock) -> dict[str, Any]: p_k = self.transform_p_k_h3_to_p_k(sample["matter_power_lin", "p_k"]) p_k = self.redshift_to_scale_factor_p_k(p_k) - ccl_args["pk_linear"] = { + pk_linear = { "a": scale_mpl, "k": k, "delta_matter:delta_matter": p_k, @@ -395,15 +392,11 @@ def calculate_ccl_args(self, sample: cosmosis.datablock) -> dict[str, Any]: p_k = self.transform_p_k_h3_to_p_k(sample["matter_power_nl", "p_k"]) p_k = self.redshift_to_scale_factor_p_k(p_k) - ccl_args["pk_nonlin"] = { + pk_nonlin = { "a": scale_mpl, "k": k, "delta_matter:delta_matter": p_k, } - elif self.require_nonlinear_pk: - ccl_args["nonlinear_model"] = "halofit" - else: - ccl_args["nonlinear_model"] = None # TODO: We should have several configurable modes for this module. # In all cases, an exception will be raised (causing a program @@ -437,9 +430,15 @@ def calculate_ccl_args(self, sample: cosmosis.datablock) -> dict[str, Any]: # h_over_h0 = np.flip(sample["distances", "h"]) * hubble_radius_today h_over_h0 = self.transform_h_to_h_over_h0(sample["distances", "h"]) - ccl_args["background"] = build_ccl_background_dict( + background: Background = build_ccl_background_dict( a=scale_distances, chi=chi, h_over_h0=h_over_h0 ) + ccl_args: CCLCalculatorArgs = {"background": background} + + if pk_linear is not None: + ccl_args.update({"pk_linear": pk_linear}) + if pk_nonlin is not None: + ccl_args.update({"pk_nonlin": pk_nonlin}) return ccl_args @@ -487,7 +486,6 @@ def set_params_from_camb(self, **params_values) -> None: m_nu = params_values["mnu"] Omega_k0 = params_values["omk"] - m_nu_type = "normal" h0 = H0 / 100.0 h02 = h0 * h0 Omega_b0 = ombh2 / h02 @@ -513,7 +511,6 @@ def set_params_from_camb(self, **params_values) -> None: sigma8=None, A_s=As, m_nu=m_nu, - m_nu_type=m_nu_type, w0=w, wa=wa, Neff=Neff, diff --git a/firecrown/connector/numcosmo/numcosmo.py b/firecrown/connector/numcosmo/numcosmo.py index 0f88d04cd..0e39f934e 100644 --- a/firecrown/connector/numcosmo/numcosmo.py +++ b/firecrown/connector/numcosmo/numcosmo.py @@ -4,9 +4,9 @@ be used without an installation of NumCosmo. """ -from typing import Any +import warnings + import numpy as np -import pyccl as ccl from numcosmo_py import Nc, Ncm, GObject, var_dict_to_dict, dict_to_var_dict @@ -17,6 +17,12 @@ from firecrown.parameters import ParamsMap from firecrown.connector.mapping import Mapping, build_ccl_background_dict from firecrown.modeling_tools import ModelingTools +from firecrown.ccl_factory import ( + CCLCalculatorArgs, + PowerSpec, + CCLFactory, + PoweSpecAmplitudeParameter, +) def get_hiprim(hi_cosmo: Nc.HICosmo) -> Nc.HIPrimPowerLaw: @@ -48,27 +54,33 @@ class MappingNumCosmo(GObject.Object): def __init__( self, - require_nonlinear_pk: bool = False, + require_nonlinear_pk: None | bool = None, p_ml: None | Nc.PowspecML = None, p_mnl: None | Nc.PowspecMNL = None, dist: None | Nc.Distance = None, ) -> None: """Initialize a MappingNumCosmo object. - :param require_nonlinear_pk: whether to require a nonlinear power spectrum :param p_ml: optional PowspecML object :param p_mnl: optional PowspecMNL object :param dist: optional Distance object """ - super().__init__( # type: ignore - require_nonlinear_pk=require_nonlinear_pk, p_ml=p_ml, p_mnl=p_mnl, dist=dist - ) + super().__init__(p_ml=p_ml, p_mnl=p_mnl, dist=dist) # type: ignore self.mapping: Mapping self._mapping_name: str self._p_ml: None | Nc.PowspecML self._p_mnl: None | Nc.PowspecMNL self._dist: Nc.Distance + if require_nonlinear_pk is not None: + warnings.warn( + "The require_nonlinear_pk argument is deprecated and will be removed " + "in future versions. This configuration is now handled by the " + "likelihood factory function.", + DeprecationWarning, + stacklevel=2, + ) + if not hasattr(self, "_p_ml"): self._p_ml = None @@ -101,28 +113,6 @@ def _set_mapping_name(self, value: str) -> None: setter=_set_mapping_name, ) - def _get_require_nonlinear_pk(self) -> bool: - """Return whether nonlinear power spectra are required. - - :return: whether nonlinear power spectra are required - """ - return self.mapping.require_nonlinear_pk - - def _set_require_nonlinear_pk(self, value: bool) -> None: - """Set whether nonlinear power spectra are required. - - :param value: whether nonlinear power spectra are required - """ - self.mapping.require_nonlinear_pk = value - - require_nonlinear_pk = GObject.Property( - type=bool, - default=False, - flags=GObject.ParamFlags.READWRITE, - getter=_get_require_nonlinear_pk, - setter=_set_require_nonlinear_pk, - ) - def _get_p_ml(self) -> None | Nc.PowspecML: """Return the NumCosmo PowspecML object. @@ -184,7 +174,7 @@ def _set_dist(self, value: Nc.Distance) -> None: ) def set_params_from_numcosmo( - self, mset: Ncm.MSet + self, mset: Ncm.MSet, ccl_factory: CCLFactory ) -> None: # pylint: disable-msg=too-many-locals """Set the parameters of the contained Mapping object. @@ -207,10 +197,7 @@ def set_params_from_numcosmo( T_gamma0 = hi_cosmo.T_gamma0() m_nu: float | list[float] = 0.0 - if hi_cosmo.NMassNu() == 0: - m_nu_type = "normal" - else: - m_nu_type = "list" + if hi_cosmo.NMassNu() > 0: assert hi_cosmo.NMassNu() <= 3 m_nu = [hi_cosmo.MassNuInfo(i)[0] for i in range(hi_cosmo.NMassNu())] @@ -224,19 +211,31 @@ def set_params_from_numcosmo( case _: raise ValueError(f"NumCosmo object {type(hi_cosmo)} not supported.") - hiprim = get_hiprim(hi_cosmo) + A_s = None + sigma8 = None + match ccl_factory.amplitude_parameter: + case PoweSpecAmplitudeParameter.SIGMA8: + if self._p_ml is None: + raise ValueError( + "PowspecML object must be provided when using sigma8." + ) + sigma8 = self._p_ml.sigma_tophat_R(hi_cosmo, 1.0e-7, 0.0, 8.0 / h) + case PoweSpecAmplitudeParameter.AS: + A_s = get_hiprim(hi_cosmo).SA_Ampl() + + assert (A_s is not None) or (sigma8 is not None) # pylint: disable=duplicate-code self.mapping.set_params( Omega_c=Omega_c, Omega_b=Omega_b, h=h, - A_s=hiprim.SA_Ampl(), - n_s=hiprim.props.n_SA, + A_s=A_s, + sigma8=sigma8, + n_s=get_hiprim(hi_cosmo).props.n_SA, Omega_k=Omega_k, Neff=Neff, m_nu=m_nu, - m_nu_type=m_nu_type, w0=w0, wa=wa, T_CMB=T_gamma0, @@ -245,13 +244,14 @@ def set_params_from_numcosmo( def calculate_ccl_args( # pylint: disable-msg=too-many-locals self, mset: Ncm.MSet - ) -> dict[str, Any]: + ) -> CCLCalculatorArgs: """Calculate the arguments necessary for CCL for this sample. :param mset: the NumCosmo MSet object from which to get the parameters :return: a dictionary of the arguments required by CCL """ - ccl_args: dict[str, Any] = {} + pk_linear: None | PowerSpec = None + pk_nonlin: None | PowerSpec = None hi_cosmo = mset.peek(Nc.HICosmo.id()) assert isinstance(hi_cosmo, Nc.HICosmo) @@ -265,8 +265,7 @@ def calculate_ccl_args( # pylint: disable-msg=too-many-locals np.array(p_m_spline.peek_zm().dup_array()).reshape(len(k), len(z)) ) p_k = self.mapping.redshift_to_scale_factor_p_k(p_k) - - ccl_args["pk_linear"] = { + pk_linear = { "a": scale, "k": k, "delta_matter:delta_matter": p_k, @@ -282,16 +281,11 @@ def calculate_ccl_args( # pylint: disable-msg=too-many-locals np.array(p_mnl_spline.peek_zm().dup_array()).reshape(len(k), len(z)) ) p_mnl = self.mapping.redshift_to_scale_factor_p_k(p_mnl) - - ccl_args["pk_nonlin"] = { + pk_nonlin = { "a": scale_mpnl, "k": k, "delta_matter:delta_matter": p_mnl, } - elif self.mapping.require_nonlinear_pk: - ccl_args["nonlinear_model"] = "halofit" - else: - ccl_args["nonlinear_model"] = None d_spline = self._dist.comoving_distance_spline.peek_spline() z_dist = np.array(d_spline.get_xv().dup_array()) @@ -309,15 +303,16 @@ def calculate_ccl_args( # pylint: disable-msg=too-many-locals chi = chi[a_unique_indices] h_over_h0 = h_over_h0[a_unique_indices] - ccl_args["background"] = build_ccl_background_dict( - a=scale_distances, chi=chi, h_over_h0=h_over_h0 - ) - - ccl_args["background"] = { - "a": scale_distances, - "chi": chi, - "h_over_h0": h_over_h0, + ccl_args: CCLCalculatorArgs = { + "background": build_ccl_background_dict( + a=scale_distances, chi=chi, h_over_h0=h_over_h0 + ) } + if pk_linear: + ccl_args["pk_linear"] = pk_linear + if pk_nonlin: + ccl_args["pk_nonlin"] = pk_nonlin + return ccl_args def create_params_map(self, model_list: list[str], mset: Ncm.MSet) -> ParamsMap: @@ -341,15 +336,22 @@ def create_params_map(self, model_list: list[str], mset: Ncm.MSet) -> ParamsMap: model_dict = { param: model.param_get_by_name(param) for param in param_names } - shared_keys = set(model_dict).intersection(params_map) - if len(shared_keys) > 0: - raise RuntimeError( - f"The following keys `{shared_keys}` appear " - f"in more than one model used by the " - f"module {model_list}." - ) - params_map = ParamsMap({**params_map, **model_dict}) + params_map = self._update_params_map(model_list, params_map, model_dict) + params_map = self._update_params_map( + model_list, params_map, self.mapping.asdict() + ) + + return params_map + + def _update_params_map(self, model_list, params_map, model_dict): + shared_keys = set(model_dict).intersection(params_map) + if len(shared_keys) > 0: + raise RuntimeError( + f"The following keys `{shared_keys}` appear in more than one model " + f"used by the module {model_list} or cosmological parameters." + ) + params_map = ParamsMap({**params_map, **model_dict}) return params_map @@ -370,7 +372,6 @@ def __init__(self): super().__init__() self.likelihood: Likelihood self.tools: ModelingTools - self.ccl_cosmo: None | ccl.Cosmology = None self._model_list: list[str] self._nc_mapping: MappingNumCosmo self._likelihood_source: None | str = None @@ -566,16 +567,13 @@ def do_prepare( # pylint: disable-msg=arguments-differ self.likelihood.reset() self.tools.reset() - self._nc_mapping.set_params_from_numcosmo(mset) + self._nc_mapping.set_params_from_numcosmo(mset, self.tools.ccl_factory) ccl_args = self._nc_mapping.calculate_ccl_args(mset) - self.ccl_cosmo = ccl.CosmologyCalculator( - **self._nc_mapping.mapping.asdict(), **ccl_args - ) params_map = self._nc_mapping.create_params_map(self.model_list, mset) self.likelihood.update(params_map) self.tools.update(params_map) - self.tools.prepare(self.ccl_cosmo) + self.tools.prepare(calculator_args=ccl_args) def do_m2lnL_val(self, _) -> float: # pylint: disable-msg=arguments-differ """Implements the virtual method `m2lnL`. @@ -615,7 +613,6 @@ def __init__(self) -> None: super().__init__() self.likelihood: ConstGaussian self.tools: ModelingTools - self.ccl_cosmo: ccl.Cosmology self.dof: int self.len: int self._model_list: list[str] @@ -675,7 +672,6 @@ def _configure_object(self) -> None: assert nrows == ncols self.set_size(nrows) - self.ccl_cosmo = None self.dof = nrows self.len = nrows self.peek_cov().set_from_array( # pylint: disable-msg=no-member @@ -843,17 +839,13 @@ def do_prepare( # pylint: disable-msg=arguments-differ self.likelihood.reset() self.tools.reset() - self._nc_mapping.set_params_from_numcosmo(mset) + self._nc_mapping.set_params_from_numcosmo(mset, self.tools.ccl_factory) ccl_args = self._nc_mapping.calculate_ccl_args(mset) - - self.ccl_cosmo = ccl.CosmologyCalculator( - **self._nc_mapping.mapping.asdict(), **ccl_args - ) params_map = self._nc_mapping.create_params_map(self._model_list, mset) self.likelihood.update(params_map) self.tools.update(params_map) - self.tools.prepare(self.ccl_cosmo) + self.tools.prepare(calculator_args=ccl_args) # pylint: disable-next=arguments-differ def do_mean_func(self, _, vp) -> None: diff --git a/firecrown/likelihood/factories.py b/firecrown/likelihood/factories.py new file mode 100644 index 000000000..a00510a9a --- /dev/null +++ b/firecrown/likelihood/factories.py @@ -0,0 +1,224 @@ +""" +Factory functions for creating likelihoods from SACC files. + +This module provides factory functions to create likelihood objects by combining a SACC +file and a set of statistic factories. Users can define their own custom statistic +factories for advanced use cases or rely on the generic factory functions provided here +for simpler scenarios. + +For straightforward contexts where all data in the SACC file is utilized, the generic +factories simplify the process. The user only needs to supply the SACC file and specify +which statistic factories to use, and the likelihood factory will handle the creation of +the likelihood object, assembling the necessary components automatically. + +These functions are particularly useful when the full set of statistics present in a +SACC file is being used without the need for complex customization. +""" + +from typing import Annotated +from enum import Enum, auto +from pathlib import Path + +import yaml +from pydantic import BaseModel, ConfigDict, BeforeValidator + +import sacc +from firecrown.likelihood.likelihood import Likelihood, NamedParameters +from firecrown.likelihood.gaussian import ConstGaussian +from firecrown.likelihood.weak_lensing import WeakLensingFactory +from firecrown.likelihood.number_counts import NumberCountsFactory +from firecrown.likelihood.two_point import TwoPoint +from firecrown.data_functions import ( + extract_all_real_data, + extract_all_harmonic_data, + check_two_point_consistence_real, + check_two_point_consistence_harmonic, +) +from firecrown.modeling_tools import ModelingTools +from firecrown.utils import YAMLSerializable + + +class TwoPointCorrelationSpace(YAMLSerializable, str, Enum): + """This class defines the two-point correlation space. + + The two-point correlation space can be either real or harmonic. The real space + corresponds measurements in terms of angular separation, while the harmonic space + corresponds to measurements in terms of spherical harmonics decomposition. + """ + + @staticmethod + def _generate_next_value_(name, _start, _count, _last_values): + return name.lower() + + REAL = auto() + HARMONIC = auto() + + +def _validate_correlation_space(value): + if isinstance(value, str): + try: + return TwoPointCorrelationSpace(value) # Convert from string to Enum + except ValueError as exc: + raise ValueError( + f"Invalid value for TwoPointCorrelationSpace: {value}" + ) from exc + return value + + +class TwoPointFactory(BaseModel): + """Factory class for WeakLensing objects.""" + + model_config = ConfigDict(extra="forbid", frozen=True) + + correlation_space: Annotated[ + TwoPointCorrelationSpace, BeforeValidator(_validate_correlation_space) + ] + weak_lensing_factory: WeakLensingFactory + number_counts_factory: NumberCountsFactory + + def model_post_init(self, __context) -> None: + """Initialize the WeakLensingFactory object.""" + + +class DataSourceSacc(BaseModel): + """Model for the data source in a likelihood configuration.""" + + sacc_data_file: str + + def model_post_init(self, __context) -> None: + """Initialize the DataSourceSacc object.""" + sacc_data_file = Path(self.sacc_data_file) + if not sacc_data_file.exists(): + raise FileNotFoundError(f"File {sacc_data_file} does not exist") + + def get_sacc_data(self) -> sacc.Sacc: + """Load the SACC data file.""" + return sacc.Sacc.load_fits(self.sacc_data_file) + + +class TwoPointExperiment(BaseModel): + """Model for the two-point experiment in a likelihood configuration.""" + + two_point_factory: TwoPointFactory + data_source: DataSourceSacc + + def model_post_init(self, __context) -> None: + """Initialize the TwoPointExperiment object.""" + + +def build_two_point_likelihood( + build_parameters: NamedParameters, +) -> tuple[Likelihood, ModelingTools]: + """ + Build a likelihood object for two-point statistics from a SACC file. + + This function creates a likelihood object for two-point statistics using a SACC file + and a set of statistic factories. The user must provide the SACC file and specify + which statistic factories to use. The likelihood object is created by combining the + SACC file with the specified statistic factories. + + :param build_parameters: A NamedParameters object containing the following + parameters: + - sacc_file: The SACC file containing the data. + - statistic_factories: A YAML file containing the statistic factories to use. + """ + likelihood_config_file = build_parameters.get_string("likelihood_config") + + with open(likelihood_config_file, "r", encoding="utf-8") as f: + likelihood_config = yaml.safe_load(f) + + if likelihood_config is None: + raise ValueError("No likelihood config found.") + + exp = TwoPointExperiment.model_validate(likelihood_config, strict=True) + modeling_tools = ModelingTools() + + # Load the SACC file + sacc_data = exp.data_source.get_sacc_data() + + match exp.two_point_factory.correlation_space: + case TwoPointCorrelationSpace.REAL: + likelihood = _build_two_point_likelihood_real( + sacc_data, + exp.two_point_factory.weak_lensing_factory, + exp.two_point_factory.number_counts_factory, + ) + case TwoPointCorrelationSpace.HARMONIC: + likelihood = _build_two_point_likelihood_harmonic( + sacc_data, + exp.two_point_factory.weak_lensing_factory, + exp.two_point_factory.number_counts_factory, + ) + + return likelihood, modeling_tools + + +def _build_two_point_likelihood_harmonic( + sacc_data: sacc.Sacc, + wl_factory: WeakLensingFactory, + nc_factory: NumberCountsFactory, +): + """ + Build a likelihood object for two-point statistics in harmonic space. + + This function creates a likelihood object for two-point statistics in harmonic space + using a SACC file and a set of statistic factories. The user must provide the SACC + file and specify which statistic factories to use. The likelihood object is created + by combining the SACC file with the specified statistic factories. + + :param sacc_data: The SACC file containing the data. + :param wl_factory: The weak lensing statistic factory. + :param nc_factory: The number counts statistic factory. + + :return: A likelihood object for two-point statistics in harmonic space. + """ + tpms = extract_all_harmonic_data(sacc_data) + if len(tpms) == 0: + raise ValueError( + "No two-point measurements in harmonic space found in the SACC file." + ) + + check_two_point_consistence_harmonic(tpms) + + two_points = TwoPoint.from_measurement( + tpms, wl_factory=wl_factory, nc_factory=nc_factory + ) + + likelihood = ConstGaussian.create_ready(two_points, sacc_data.covariance.dense) + + return likelihood + + +def _build_two_point_likelihood_real( + sacc_data: sacc.Sacc, + wl_factory: WeakLensingFactory, + nc_factory: NumberCountsFactory, +): + """ + Build a likelihood object for two-point statistics in real space. + + This function creates a likelihood object for two-point statistics in real space + using a SACC file and a set of statistic factories. The user must provide the SACC + file and specify which statistic factories to use. The likelihood object is created + by combining the SACC file with the specified statistic factories. + + :param sacc_data: The SACC file containing the data. + :param wl_factory: The weak lensing statistic factory. + :param nc_factory: The number counts statistic factory. + + :return: A likelihood object for two-point statistics in real space. + """ + tpms = extract_all_real_data(sacc_data) + if len(tpms) == 0: + raise ValueError( + "No two-point measurements in real space found in the SACC file." + ) + check_two_point_consistence_real(tpms) + + two_points = TwoPoint.from_measurement( + tpms, wl_factory=wl_factory, nc_factory=nc_factory + ) + + likelihood = ConstGaussian.create_ready(two_points, sacc_data.covariance.dense) + + return likelihood diff --git a/firecrown/likelihood/likelihood.py b/firecrown/likelihood/likelihood.py index 2aa751552..9059f3bad 100644 --- a/firecrown/likelihood/likelihood.py +++ b/firecrown/likelihood/likelihood.py @@ -321,7 +321,9 @@ def convert_to_basic_dict( def load_likelihood_from_module_type( - module: types.ModuleType, build_parameters: NamedParameters + module: types.ModuleType, + build_parameters: NamedParameters, + build_likelihood_name: str = "build_likelihood", ) -> tuple[Likelihood, ModelingTools]: """Loads a likelihood from a module type. @@ -337,28 +339,29 @@ def load_likelihood_from_module_type( function parameters :return : a tuple of the likelihood and the modeling tools """ - if not hasattr(module, "build_likelihood"): + if not hasattr(module, build_likelihood_name): if not hasattr(module, "likelihood"): raise AttributeError( f"Firecrown initialization module {module.__name__} in " f"{module.__file__} does not define " - f"a `build_likelihood` factory function." + f"a `{build_likelihood_name}` factory function." ) warnings.warn( - "The use of a likelihood variable in Firecrown's initialization " - "module is deprecated. Any parameters passed to the likelihood " - "will be ignored. The module should define a `build_likelihood` " - "factory function.", + f"The use of a likelihood variable in Firecrown's initialization " + f"module is deprecated. Any parameters passed to the likelihood " + f"will be ignored. The module should define the `{build_likelihood_name}` " + f"factory function.", category=DeprecationWarning, ) likelihood = module.likelihood tools = ModelingTools() else: - if not callable(module.build_likelihood): + build_likelihood = getattr(module, build_likelihood_name) + if not callable(build_likelihood): raise TypeError( - "The factory function `build_likelihood` must be a callable." + f"The factory function `{build_likelihood_name}` must be a callable." ) - build_return = module.build_likelihood(build_parameters) + build_return = build_likelihood(build_parameters) if isinstance(build_return, tuple): likelihood, tools = build_return else: @@ -443,13 +446,26 @@ def load_likelihood_from_module( :return : a tuple of the likelihood and the modeling tools """ try: - mod = importlib.import_module(module) + # Try importing the entire string as a module first + try: + mod = importlib.import_module(module) + func = "build_likelihood" + except ImportError as sub_exc: + # If it fails, split and try importing as module.function + if "." not in module: + raise sub_exc + module_name, func = module.rsplit(".", 1) + mod = importlib.import_module(module_name) except ImportError as exc: raise ValueError( - f"Unrecognized Firecrown initialization module {module}." + f"Unrecognized Firecrown initialization module '{module}'. " + f"The module must be either a module_name or a module_name.func " + f"where func is the factory function." ) from exc - return load_likelihood_from_module_type(mod, build_parameters) + return load_likelihood_from_module_type( + mod, build_parameters, build_likelihood_name=func + ) def load_likelihood( diff --git a/firecrown/modeling_tools.py b/firecrown/modeling_tools.py index 35e947967..9a53eb6a1 100644 --- a/firecrown/modeling_tools.py +++ b/firecrown/modeling_tools.py @@ -15,6 +15,7 @@ from firecrown.models.cluster.abundance import ClusterAbundance from firecrown.updatable import Updatable, UpdatableCollection +from firecrown.ccl_factory import CCLFactory, CCLCalculatorArgs class ModelingTools(Updatable): @@ -30,6 +31,7 @@ def __init__( pt_calculator: None | pyccl.nl_pt.EulerianPTCalculator = None, pk_modifiers: None | Collection[PowerspectrumModifier] = None, cluster_abundance: None | ClusterAbundance = None, + ccl_factory: None | CCLFactory = None, ): super().__init__() self.ccl_cosmo: None | pyccl.Cosmology = None @@ -39,6 +41,7 @@ def __init__( self.powerspectra: dict[str, pyccl.Pk2D] = {} self._prepared: bool = False self.cluster_abundance = cluster_abundance + self.ccl_factory = CCLFactory() if ccl_factory is None else ccl_factory def add_pk(self, name: str, powerspectrum: pyccl.Pk2D) -> None: """Add a :python:`pyccl.Pk2D` to the table of power spectra.""" @@ -71,7 +74,7 @@ def has_pk(self, name: str) -> bool: return False return True - def prepare(self, ccl_cosmo: pyccl.Cosmology) -> None: + def prepare(self, *, calculator_args: None | CCLCalculatorArgs = None) -> None: """Prepare the Cosmology for use in likelihoods. This method will prepare the ModelingTools for use in likelihoods. This @@ -88,16 +91,17 @@ def prepare(self, ccl_cosmo: pyccl.Cosmology) -> None: if self.ccl_cosmo is not None: raise RuntimeError("Cosmology has already been set") - self.ccl_cosmo = ccl_cosmo + + self.ccl_cosmo = self.ccl_factory.create(calculator_args) if self.pt_calculator is not None: - self.pt_calculator.update_ingredients(ccl_cosmo) + self.pt_calculator.update_ingredients(self.ccl_cosmo) for pkm in self.pk_modifiers: self.add_pk(name=pkm.name, powerspectrum=pkm.compute_p_of_k_z(tools=self)) if self.cluster_abundance is not None: - self.cluster_abundance.update_ingredients(ccl_cosmo) + self.cluster_abundance.update_ingredients(self.ccl_cosmo) self._prepared = True diff --git a/firecrown/parameters.py b/firecrown/parameters.py index 81b6d9f6d..79f017af7 100644 --- a/firecrown/parameters.py +++ b/firecrown/parameters.py @@ -35,6 +35,28 @@ def parameter_get_full_name(prefix: None | str, param: str) -> str: return param +def _validade_params_map_value(name: str, value: float | list[float]) -> None: + """Check if the value is a float or a list of floats. + + Raises a TypeError if the value is not a float or a list of floats. + + :param name: name of the parameter + :param value: value to be checked + """ + if not isinstance(value, (float, list)): + raise TypeError( + f"Value for parameter {name} is not a float or a list of floats: " + f"{type(value)}" + ) + + if isinstance(value, list): + if not all(isinstance(v, float) for v in value): + raise TypeError( + f"Value for parameter {name} is not a float or a list of floats: " + f"{type(value)}" + ) + + class ParamsMap(dict[str, float]): """A specialized dict in which all keys are strings and values are floats. @@ -44,6 +66,9 @@ class ParamsMap(dict[str, float]): def __init__(self, *args, **kwargs) -> None: super().__init__(*args, **kwargs) + for name, value in self.items(): + _validade_params_map_value(name, value) + self.lower_case: bool = False def use_lower_case_keys(self, enable: bool) -> None: @@ -62,11 +87,14 @@ def get_from_full_name(self, full_name: str) -> float: Raises a KeyError if the parameter is not found. """ - if self.lower_case: - full_name = full_name.lower() - if full_name in self.keys(): return self[full_name] + + if self.lower_case: + full_name_lower = full_name.lower() + if full_name_lower in self.keys(): + return self[full_name_lower] + raise KeyError(f"Key {full_name} not found.") def get_from_prefix_param(self, prefix: None | str, param: str) -> float: @@ -109,6 +137,14 @@ def __add__(self, other: RequiredParameters) -> RequiredParameters: """ return RequiredParameters(self.params_set | other.params_set) + def __sub__(self, other: RequiredParameters) -> RequiredParameters: + """Return a new RequiredParameters with the names in self but not in other. + + Note that this function returns a new object that does not share state + with either argument to the subtraction operator. + """ + return RequiredParameters(self.params_set - other.params_set) + def __eq__(self, other: object): """Compare two RequiredParameters objects for equality. diff --git a/firecrown/updatable.py b/firecrown/updatable.py index ff11af301..ed7136f6e 100644 --- a/firecrown/updatable.py +++ b/firecrown/updatable.py @@ -46,9 +46,9 @@ def __init__(self, parameter: str) -> None: self.parameter = parameter msg = ( f"The parameter `{parameter}` is required to update " - "something in this likelihood.\nIt should have been supplied" - "by the sampling framework.\nThe object being updated was:\n" - "{context}\n" + f"something in this likelihood.\nIt should have been supplied " + f"by the sampling framework.\nThe object being updated was:\n" + f"{self}\n" ) super().__init__(msg) @@ -351,9 +351,9 @@ def __init__(self, iterable: None | Iterable[T] = None) -> None: super().__init__(iterable) self._updated: bool = False for item in self: - if not isinstance(item, Updatable): + if not isinstance(item, (Updatable | UpdatableCollection)): raise TypeError( - "All the items in an UpdatableCollection must be updatable" + f"All the items in an UpdatableCollection must be updatable {item}" ) @final @@ -435,3 +435,25 @@ def __setitem__(self, key, value): ) super().__setitem__(key, cast(T, value)) + + +def get_default_params(*args: Updatable) -> dict[str, float]: + """Get a ParamsMap with the default values of all parameters in the updatables. + + :param args: updatables to get the default parameters from + :return: a ParamsMap with the default values of all parameters + """ + updatable_collection = UpdatableCollection(args) + required_parameters = updatable_collection.required_parameters() + default_parameters = required_parameters.get_default_values() + + return default_parameters + + +def get_default_params_map(*args: Updatable) -> ParamsMap: + """Get a ParamsMap with the default values of all parameters in the updatables. + + :param args: updatables to get the default parameters from + :return: a ParamsMap with the default values of all parameters + """ + return ParamsMap(get_default_params(*args)) diff --git a/tests/conftest.py b/tests/conftest.py index df9789cc4..d135c1f47 100644 --- a/tests/conftest.py +++ b/tests/conftest.py @@ -11,8 +11,8 @@ import numpy as np import numpy.typing as npt -import pyccl +from firecrown.updatable import get_default_params_map from firecrown.utils import upper_triangle_indices from firecrown.likelihood.statistic import TrivialStatistic from firecrown.parameters import ParamsMap @@ -119,7 +119,6 @@ def fixture_mapping_cosmosis() -> MappingCosmoSIS: Omega_k=0.0, Neff=3.046, m_nu=0.0, - m_nu_type="normal", w0=-1.0, wa=0.0, T_CMB=2.7255, @@ -156,8 +155,9 @@ def fixture_tools_with_vanilla_cosmology(): """Return a ModelingTools object containing the LCDM cosmology from pyccl.""" result = ModelingTools() - result.update(ParamsMap()) - result.prepare(pyccl.CosmologyVanillaLCDM()) + params = get_default_params_map(result) + result.update(params) + result.prepare() return result diff --git a/tests/connector/cobaya/test_model_likelihood.py b/tests/connector/cobaya/test_model_likelihood.py index c6446544d..609add67e 100644 --- a/tests/connector/cobaya/test_model_likelihood.py +++ b/tests/connector/cobaya/test_model_likelihood.py @@ -1,7 +1,6 @@ """Unit tests for the cobaya Mapping connector.""" import pytest -import pyccl as ccl from cobaya.model import get_model, Model from cobaya.log import LoggedError from firecrown.connector.cobaya.ccl import CCLConnector @@ -69,7 +68,7 @@ def test_cobaya_ccl_model(fiducial_params): "likelihood": { "test_lk": { "external": lambda _self=None: 0.0, - "requires": {"pyccl": None}, + "requires": {"pyccl_args": None, "pyccl_params": None}, } }, "theory": { @@ -82,20 +81,21 @@ def test_cobaya_ccl_model(fiducial_params): assert isinstance(model_fiducial, Model) model_fiducial.logposterior({}) - cosmo = model_fiducial.provider.get_pyccl() - assert isinstance(cosmo, ccl.Cosmology) + cosmo_args = model_fiducial.provider.get_pyccl_args() + cosmo_params = model_fiducial.provider.get_pyccl_params() + assert isinstance(cosmo_args, dict) h = fiducial_params["H0"] / 100.0 - assert cosmo["H0"] == pytest.approx(fiducial_params["H0"], rel=1.0e-5) - assert cosmo["Omega_c"] == pytest.approx( + assert cosmo_params["h"] * 100.0 == pytest.approx(fiducial_params["H0"], rel=1.0e-5) + assert cosmo_params["Omega_c"] == pytest.approx( fiducial_params["omch2"] / h**2, rel=1.0e-5 ) - assert cosmo["Omega_b"] == pytest.approx( + assert cosmo_params["Omega_b"] == pytest.approx( fiducial_params["ombh2"] / h**2, rel=1.0e-5 ) - assert cosmo["Omega_k"] == pytest.approx(0.0, rel=1.0e-5) - assert cosmo["A_s"] == pytest.approx(fiducial_params["As"], rel=1.0e-5) - assert cosmo["n_s"] == pytest.approx(fiducial_params["ns"], rel=1.0e-5) + assert cosmo_params["Omega_k"] == pytest.approx(0.0, rel=1.0e-5) + assert cosmo_params["A_s"] == pytest.approx(fiducial_params["As"], rel=1.0e-5) + assert cosmo_params["n_s"] == pytest.approx(fiducial_params["ns"], rel=1.0e-5) # The following test fails because of we are using the default # neutrino hierarchy, which is normal, while CAMB depends on the # parameter which we do not have access to. diff --git a/tests/connector/cosmosis/test_cosmosis_module.py b/tests/connector/cosmosis/test_cosmosis_module.py index 43de9f902..dc4394b05 100644 --- a/tests/connector/cosmosis/test_cosmosis_module.py +++ b/tests/connector/cosmosis/test_cosmosis_module.py @@ -555,29 +555,3 @@ def test_mapping_cosmosis_pk_nonlin(mapping_cosmosis): mapping_cosmosis.transform_p_k_h3_to_p_k(matter_power_nl_p_k) ), ) - - -def test_mapping_cosmosis_pk_nonlin_nonlinear_model(mapping_cosmosis): - block = DataBlock() - - block.put_double_array_1d("distances", "d_m", np.geomspace(0.1, 10.0, 100)) - block.put_double_array_1d("distances", "z", np.linspace(0.0, 2.0, 10)) - block.put_double_array_1d("distances", "h", np.geomspace(0.1, 10.0, 100)) - - mapping_cosmosis.require_nonlinear_pk = True - ccl_args = mapping_cosmosis.calculate_ccl_args(block) - - assert "background" in ccl_args - assert "pk_linear" not in ccl_args - assert "pk_nonlin" not in ccl_args - assert "nonlinear_model" in ccl_args - assert ccl_args["nonlinear_model"] - - mapping_cosmosis.require_nonlinear_pk = False - ccl_args = mapping_cosmosis.calculate_ccl_args(block) - - assert "background" in ccl_args - assert "pk_linear" not in ccl_args - assert "pk_nonlin" not in ccl_args - assert "nonlinear_model" in ccl_args - assert not ccl_args["nonlinear_model"] diff --git a/tests/connector/numcosmo/test_numcosmo_connector.py b/tests/connector/numcosmo/test_numcosmo_connector.py index cc8ca58d1..eae4a4e0d 100644 --- a/tests/connector/numcosmo/test_numcosmo_connector.py +++ b/tests/connector/numcosmo/test_numcosmo_connector.py @@ -18,10 +18,7 @@ @pytest.fixture(name="factory_plain") def fixture_factory_plain(): """Create a NumCosmoFactory instance.""" - map_cosmo = MappingNumCosmo( - require_nonlinear_pk=True, - dist=Nc.Distance.new(6.0), - ) + map_cosmo = MappingNumCosmo(dist=Nc.Distance.new(6.0)) build_parameters = NamedParameters() return NumCosmoFactory( "tests/likelihood/lkdir/lkscript.py", @@ -34,10 +31,7 @@ def fixture_factory_plain(): @pytest.fixture(name="factory_const_gauss") def fixture_factory_const_gauss(): """Create a NumCosmoFactory instance.""" - map_cosmo = MappingNumCosmo( - require_nonlinear_pk=True, - dist=Nc.Distance.new(6.0), - ) + map_cosmo = MappingNumCosmo(dist=Nc.Distance.new(6.0)) build_parameters = NamedParameters() return NumCosmoFactory( "tests/likelihood/gauss_family/lkscript_const_gaussian.py", diff --git a/tests/connector/numcosmo/test_numcosmo_mapping.py b/tests/connector/numcosmo/test_numcosmo_mapping.py index a116ede25..13640d058 100644 --- a/tests/connector/numcosmo/test_numcosmo_mapping.py +++ b/tests/connector/numcosmo/test_numcosmo_mapping.py @@ -14,6 +14,7 @@ MappingNumCosmo, NumCosmoFactory, ) +from firecrown.ccl_factory import CCLFactory, PoweSpecAmplitudeParameter Ncm.cfg_init() @@ -24,11 +25,7 @@ def test_numcosmo_mapping_create_params_map_non_existing_model(): cosmo = Nc.HICosmoDEXcdm() - map_cosmo = MappingNumCosmo( - require_nonlinear_pk=True, - dist=Nc.Distance.new(6.0), - ) - + map_cosmo = MappingNumCosmo(dist=Nc.Distance.new(6.0)) mset = Ncm.MSet() mset.set(cosmo) @@ -45,11 +42,7 @@ def test_numcosmo_mapping_create_params_map_absent_model(): cosmo = Nc.HICosmoDEXcdm() - map_cosmo = MappingNumCosmo( - require_nonlinear_pk=True, - dist=Nc.Distance.new(6.0), - ) - + map_cosmo = MappingNumCosmo(dist=Nc.Distance.new(6.0)) mset = Ncm.MSet() mset.set(cosmo) @@ -65,11 +58,7 @@ def test_numcosmo_mapping_create_params_map_two_models_sharing_parameters(): with an existing type but not present in the model set.""" cosmo = Nc.HICosmoDEXcdm() - - map_cosmo = MappingNumCosmo( - require_nonlinear_pk=True, - dist=Nc.Distance.new(6.0), - ) + map_cosmo = MappingNumCosmo(dist=Nc.Distance.new(6.0)) mset = Ncm.MSet() mset.set(cosmo) @@ -153,16 +142,12 @@ def test_numcosmo_mapping_unsupported(): cosmo = Nc.HICosmoDEJbp() - map_cosmo = MappingNumCosmo( - require_nonlinear_pk=True, - dist=Nc.Distance.new(6.0), - ) - + map_cosmo = MappingNumCosmo(dist=Nc.Distance.new(6.0)) mset = Ncm.MSet() mset.set(cosmo) with pytest.raises(ValueError, match="NumCosmo object .* not supported."): - map_cosmo.set_params_from_numcosmo(mset) + map_cosmo.set_params_from_numcosmo(mset, CCLFactory()) def test_numcosmo_mapping_missing_hiprim(): @@ -170,18 +155,16 @@ def test_numcosmo_mapping_missing_hiprim(): cosmo = Nc.HICosmoDECpl() - map_cosmo = MappingNumCosmo( - require_nonlinear_pk=True, - dist=Nc.Distance.new(6.0), - ) - + map_cosmo = MappingNumCosmo(dist=Nc.Distance.new(6.0)) mset = Ncm.MSet() mset.set(cosmo) with pytest.raises( ValueError, match="NumCosmo object must include a HIPrim object." ): - map_cosmo.set_params_from_numcosmo(mset) + map_cosmo.set_params_from_numcosmo( + mset, CCLFactory(amplitude_parameter=PoweSpecAmplitudeParameter.AS) + ) def test_numcosmo_mapping_invalid_hiprim(): @@ -191,64 +174,16 @@ def test_numcosmo_mapping_invalid_hiprim(): prim = Nc.HIPrimAtan() cosmo.add_submodel(prim) - map_cosmo = MappingNumCosmo( - require_nonlinear_pk=True, - dist=Nc.Distance.new(6.0), - ) - + map_cosmo = MappingNumCosmo(dist=Nc.Distance.new(6.0)) mset = Ncm.MSet() mset.set(cosmo) with pytest.raises( ValueError, match="NumCosmo HIPrim object type .* not supported." ): - map_cosmo.set_params_from_numcosmo(mset) - - -def test_numcosmo_mapping_no_p_mnl_require_nonlinear_pk(): - """Test the NumCosmo mapping connector with a model without p_mnl but - with require_nonlinear_pk=True.""" - - cosmo = Nc.HICosmoDECpl() - prim = Nc.HIPrimPowerLaw() - cosmo.add_submodel(prim) - - map_cosmo = MappingNumCosmo( - require_nonlinear_pk=True, - dist=Nc.Distance.new(6.0), - ) - - mset = Ncm.MSet() - mset.set(cosmo) - - map_cosmo.set_params_from_numcosmo(mset) - - ccl_args = map_cosmo.calculate_ccl_args(mset) - - assert ccl_args["nonlinear_model"] == "halofit" - - -def test_numcosmo_mapping_no_p_mnl(): - """Test the NumCosmo mapping connector with a model without p_mnl and - require_nonlinear_pk=False.""" - - cosmo = Nc.HICosmoDECpl() - prim = Nc.HIPrimPowerLaw() - cosmo.add_submodel(prim) - - map_cosmo = MappingNumCosmo( - require_nonlinear_pk=False, - dist=Nc.Distance.new(6.0), - ) - - mset = Ncm.MSet() - mset.set(cosmo) - - map_cosmo.set_params_from_numcosmo(mset) - - ccl_args = map_cosmo.calculate_ccl_args(mset) - - assert ccl_args["nonlinear_model"] is None + map_cosmo.set_params_from_numcosmo( + mset, CCLFactory(amplitude_parameter=PoweSpecAmplitudeParameter.AS) + ) @pytest.mark.parametrize( @@ -261,7 +196,6 @@ def test_numcosmo_mapping(numcosmo_cosmo_fixture, request): cosmo = numcosmo_cosmo["cosmo"] map_cosmo = MappingNumCosmo( - require_nonlinear_pk=True, p_ml=numcosmo_cosmo["p_ml"], p_mnl=numcosmo_cosmo["p_mnl"], dist=numcosmo_cosmo["dist"], @@ -270,7 +204,7 @@ def test_numcosmo_mapping(numcosmo_cosmo_fixture, request): mset = Ncm.MSet() mset.set(cosmo) - map_cosmo.set_params_from_numcosmo(mset) + map_cosmo.set_params_from_numcosmo(mset, CCLFactory()) ccl_args = map_cosmo.calculate_ccl_args(mset) ccl_cosmo = ccl.CosmologyCalculator(**map_cosmo.mapping.asdict(), **ccl_args) @@ -290,7 +224,6 @@ def test_numcosmo_serialize_mapping(numcosmo_cosmo_fixture, request): cosmo = numcosmo_cosmo["cosmo"] map_cosmo = MappingNumCosmo( - require_nonlinear_pk=True, p_ml=numcosmo_cosmo["p_ml"], p_mnl=numcosmo_cosmo["p_mnl"], dist=numcosmo_cosmo["dist"], @@ -307,8 +240,8 @@ def test_numcosmo_serialize_mapping(numcosmo_cosmo_fixture, request): mset = Ncm.MSet() mset.set(cosmo) - map_cosmo.set_params_from_numcosmo(mset) - map_cosmo_dup.set_params_from_numcosmo(mset) + map_cosmo.set_params_from_numcosmo(mset, CCLFactory()) + map_cosmo_dup.set_params_from_numcosmo(mset, CCLFactory()) if map_cosmo_dup.p_ml is None: assert map_cosmo_dup.p_ml is None @@ -343,7 +276,6 @@ def test_numcosmo_data( cosmo = numcosmo_cosmo["cosmo"] map_cosmo = MappingNumCosmo( - require_nonlinear_pk=True, p_ml=numcosmo_cosmo["p_ml"], p_mnl=numcosmo_cosmo["p_mnl"], dist=numcosmo_cosmo["dist"], @@ -395,7 +327,6 @@ def test_numcosmo_gauss_cov( cosmo = numcosmo_cosmo["cosmo"] map_cosmo = MappingNumCosmo( - require_nonlinear_pk=True, p_ml=numcosmo_cosmo["p_ml"], p_mnl=numcosmo_cosmo["p_mnl"], dist=numcosmo_cosmo["dist"], @@ -451,7 +382,6 @@ def test_numcosmo_serialize_likelihood( numcosmo_cosmo = request.getfixturevalue(numcosmo_cosmo_fixture) map_cosmo = MappingNumCosmo( - require_nonlinear_pk=True, p_ml=numcosmo_cosmo["p_ml"], p_mnl=numcosmo_cosmo["p_mnl"], dist=numcosmo_cosmo["dist"], @@ -490,3 +420,27 @@ def test_numcosmo_serialize_likelihood( assert id(fc_data_dup.model_list) != id(fc_data.model_list) assert isinstance(fc_data_dup.model_list, list) assert fc_data_dup.model_list == fc_data.model_list + + +def test_numcosmo_mapping_deprecated_require_nonlinear_pk(): + """Test the MappingNumCosmo deprecated require_nonlinear_pk.""" + + with pytest.deprecated_call(): + _ = MappingNumCosmo(require_nonlinear_pk=True) + + +def test_numcosmo_mapping_sigma8_missing_pk(): + """Test the MappingNumCosmo with sigma8 as a parameter but missing pk.""" + + cosmo = Nc.HICosmoDEXcdm() + + map_cosmo = MappingNumCosmo(dist=Nc.Distance.new(6.0)) + mset = Ncm.MSet() + mset.set(cosmo) + + with pytest.raises( + ValueError, match="PowspecML object must be provided when using sigma8." + ): + map_cosmo.set_params_from_numcosmo( + mset, CCLFactory(amplitude_parameter=PoweSpecAmplitudeParameter.SIGMA8) + ) diff --git a/tests/connector/test_mapping.py b/tests/connector/test_mapping.py index 3fef86935..27824c2b8 100644 --- a/tests/connector/test_mapping.py +++ b/tests/connector/test_mapping.py @@ -53,7 +53,6 @@ def test_conversion_from_cosmosis_camb(): assert p.Omega_g is None assert p.Neff == pytest.approx(3.046) assert p.m_nu == pytest.approx(0.3015443336635814) - assert p.m_nu_type == "normal" # Currently the only option assert p.w0 == cosmosis_params["w"] assert p.wa == cosmosis_params["wa"] assert p.T_CMB == 2.7255 # currently hard-wired @@ -135,7 +134,6 @@ def test_sigma8_and_A_s(): "Omega_k": 0.0, "Neff": 3.046, "m_nu": 0.0, - "m_nu_type": "normal", "w0": -1.0, "wa": 0.0, "T_CMB": 2.7255, diff --git a/tests/integration/test_des_y1_3x2pt.py b/tests/integration/test_des_y1_3x2pt.py index 75c3a7a89..99fee1ab3 100644 --- a/tests/integration/test_des_y1_3x2pt.py +++ b/tests/integration/test_des_y1_3x2pt.py @@ -15,6 +15,7 @@ def test_des_y1_3x2pt_cosmosis(): cd examples/des_y1_3x2pt cosmosis des_y1_3x2pt.ini cosmosis des_y1_3x2pt_PT.ini + cosmosis des_y1_3x2pt_default_factory.ini """, ], capture_output=True, diff --git a/tests/likelihood/gauss_family/lkscript_const_gaussian.py b/tests/likelihood/gauss_family/lkscript_const_gaussian.py index c58be271b..e8a378933 100644 --- a/tests/likelihood/gauss_family/lkscript_const_gaussian.py +++ b/tests/likelihood/gauss_family/lkscript_const_gaussian.py @@ -7,6 +7,8 @@ from firecrown.likelihood.gaussian import ConstGaussian from firecrown.likelihood.supernova import Supernova +from firecrown.modeling_tools import ModelingTools +from firecrown.ccl_factory import CCLFactory, PoweSpecAmplitudeParameter def build_likelihood(_): @@ -24,4 +26,9 @@ def build_likelihood(_): sacc_data.add_data_point("supernova_distance_mu", tracer_tuple, -3.0, z=0.3) sacc_data.add_covariance([4.0, 9.0, 16.0]) likelihood.read(sacc_data) - return likelihood + + tools = ModelingTools( + ccl_factory=CCLFactory(amplitude_parameter=PoweSpecAmplitudeParameter.SIGMA8) + ) + + return likelihood, tools diff --git a/tests/likelihood/gauss_family/lkscript_two_point.py b/tests/likelihood/gauss_family/lkscript_two_point.py index b8ec49e34..fff9029dc 100644 --- a/tests/likelihood/gauss_family/lkscript_two_point.py +++ b/tests/likelihood/gauss_family/lkscript_two_point.py @@ -12,6 +12,8 @@ from firecrown.likelihood.number_counts import ( NumberCounts, ) +from firecrown.modeling_tools import ModelingTools +from firecrown.ccl_factory import CCLFactory, PoweSpecAmplitudeParameter def build_likelihood(params: NamedParameters): @@ -42,4 +44,9 @@ def build_likelihood(params: NamedParameters): likelihood = ConstGaussian(statistics=[two_point]) likelihood.read(sacc_data) - return likelihood + + tools = ModelingTools( + ccl_factory=CCLFactory(amplitude_parameter=PoweSpecAmplitudeParameter.SIGMA8) + ) + + return likelihood, tools diff --git a/tests/likelihood/gauss_family/statistic/test_binned_cluster_number_counts.py b/tests/likelihood/gauss_family/statistic/test_binned_cluster_number_counts.py index 9597320e1..58e5f2dcb 100644 --- a/tests/likelihood/gauss_family/statistic/test_binned_cluster_number_counts.py +++ b/tests/likelihood/gauss_family/statistic/test_binned_cluster_number_counts.py @@ -4,18 +4,18 @@ import sacc import pytest import pyccl +from firecrown.updatable import get_default_params_map from firecrown.models.cluster.recipes.cluster_recipe import ClusterRecipe from firecrown.likelihood.source import SourceSystematic from firecrown.modeling_tools import ModelingTools from firecrown.models.cluster.properties import ClusterProperty -from firecrown.parameters import ParamsMap from firecrown.models.cluster.abundance import ClusterAbundance from firecrown.likelihood.binned_cluster_number_counts import ( BinnedClusterNumberCounts, ) -def test_create_binned_number_counts(): +def test_create_binned_number_counts() -> None: recipe = Mock(spec=ClusterRecipe) bnc = BinnedClusterNumberCounts(ClusterProperty.NONE, "Test", recipe) assert bnc is not None @@ -36,7 +36,7 @@ def test_create_binned_number_counts(): assert bnc.systematics == systematics -def test_get_data_vector(): +def test_get_data_vector() -> None: recipe = Mock(spec=ClusterRecipe) bnc = BinnedClusterNumberCounts(ClusterProperty.NONE, "Test", recipe) dv = bnc.get_data_vector() @@ -44,7 +44,7 @@ def test_get_data_vector(): assert len(dv) == 0 -def test_read_throws_if_no_property(cluster_sacc_data: sacc.Sacc): +def test_read_throws_if_no_property(cluster_sacc_data: sacc.Sacc) -> None: recipe = Mock(spec=ClusterRecipe) bnc = BinnedClusterNumberCounts(ClusterProperty.NONE, "my_survey", recipe) @@ -55,7 +55,7 @@ def test_read_throws_if_no_property(cluster_sacc_data: sacc.Sacc): bnc.read(cluster_sacc_data) -def test_read_single_property(cluster_sacc_data: sacc.Sacc): +def test_read_single_property(cluster_sacc_data: sacc.Sacc) -> None: recipe = Mock(spec=ClusterRecipe) bnc = BinnedClusterNumberCounts(ClusterProperty.COUNTS, "my_survey", recipe) @@ -75,7 +75,7 @@ def test_read_single_property(cluster_sacc_data: sacc.Sacc): assert len(bnc.sacc_indices) == 2 -def test_read_multiple_properties(cluster_sacc_data: sacc.Sacc): +def test_read_multiple_properties(cluster_sacc_data: sacc.Sacc) -> None: recipe = Mock(spec=ClusterRecipe) bnc = BinnedClusterNumberCounts( (ClusterProperty.COUNTS | ClusterProperty.MASS), "my_survey", recipe @@ -88,18 +88,16 @@ def test_read_multiple_properties(cluster_sacc_data: sacc.Sacc): assert len(bnc.sacc_indices) == 4 -def test_compute_theory_vector(cluster_sacc_data: sacc.Sacc): +def test_compute_theory_vector(cluster_sacc_data: sacc.Sacc) -> None: recipe = Mock(spec=ClusterRecipe) recipe.evaluate_theory_prediction.return_value = 1.0 tools = ModelingTools() hmf = pyccl.halos.MassFuncBocquet16() - cosmo = pyccl.cosmology.CosmologyVanillaLCDM() - params = ParamsMap() - tools.cluster_abundance = ClusterAbundance(13, 17, 0, 2, hmf) + params = get_default_params_map(tools) tools.update(params) - tools.prepare(cosmo) + tools.prepare() bnc = BinnedClusterNumberCounts(ClusterProperty.COUNTS, "my_survey", recipe) bnc.read(cluster_sacc_data) diff --git a/tests/likelihood/gauss_family/statistic/test_supernova.py b/tests/likelihood/gauss_family/statistic/test_supernova.py index 4a433bfc8..c3932b736 100644 --- a/tests/likelihood/gauss_family/statistic/test_supernova.py +++ b/tests/likelihood/gauss_family/statistic/test_supernova.py @@ -4,11 +4,10 @@ import pytest import sacc -import pyccl +from firecrown.updatable import get_default_params_map from firecrown.likelihood.supernova import Supernova from firecrown.modeling_tools import ModelingTools -from firecrown.parameters import ParamsMap @pytest.fixture(name="minimal_stat") @@ -48,7 +47,7 @@ def fixture_sacc_data() -> sacc.Sacc: def test_missing_sacc_tracer_fails_read( minimal_stat: Supernova, missing_sacc_tracer: sacc.Sacc -): +) -> None: with pytest.raises( ValueError, match="The SACC file does not contain the MiscTracer sn_fake_sample", @@ -58,7 +57,7 @@ def test_missing_sacc_tracer_fails_read( def test_wrong_tracer_type_fails_read( minimal_stat: Supernova, wrong_tracer_type: sacc.Sacc -): +) -> None: with pytest.raises( ValueError, match=f"The SACC tracer {minimal_stat.sacc_tracer} is not a MiscTracer", @@ -66,7 +65,7 @@ def test_wrong_tracer_type_fails_read( minimal_stat.read(wrong_tracer_type) -def test_read_works(minimal_stat: Supernova, good_sacc_data: sacc.Sacc): +def test_read_works(minimal_stat: Supernova, good_sacc_data: sacc.Sacc) -> None: """After read() is called, we should be able to get the statistic's :class:`DataVector` and also should be able to call @@ -78,9 +77,10 @@ def test_read_works(minimal_stat: Supernova, good_sacc_data: sacc.Sacc): assert data_vector[0] == 16.95 tools = ModelingTools() - tools.update(ParamsMap()) - tools.prepare(pyccl.CosmologyVanillaLCDM()) - params = ParamsMap({"sn_fake_sample_M": 1.1}) + params = get_default_params_map(tools) + params.update({"sn_fake_sample_M": 1.1}) + tools.update(params) + tools.prepare() minimal_stat.update(params) theory_vector = minimal_stat.compute_theory_vector(tools) assert len(theory_vector) == 1 diff --git a/tests/likelihood/gauss_family/statistic/test_two_point.py b/tests/likelihood/gauss_family/statistic/test_two_point.py index 2f58450ab..3eaba3858 100644 --- a/tests/likelihood/gauss_family/statistic/test_two_point.py +++ b/tests/likelihood/gauss_family/statistic/test_two_point.py @@ -8,8 +8,7 @@ from numpy.testing import assert_allclose import pytest -import pyccl - +from firecrown.updatable import get_default_params_map from firecrown.modeling_tools import ModelingTools from firecrown.parameters import ParamsMap @@ -194,20 +193,21 @@ def test_two_point_src0_src0_window(sacc_galaxy_cells_src0_src0_window) -> None: statistic.read(sacc_data) tools = ModelingTools() - tools.update(ParamsMap()) - tools.prepare(pyccl.CosmologyVanillaLCDM()) + params = get_default_params_map(tools) + tools.update(params) + tools.prepare() assert statistic.window is not None statistic.reset() - statistic.update(ParamsMap()) - tools.update(ParamsMap()) + statistic.update(params) + tools.update(params) result1 = statistic.compute_theory_vector(tools) assert all(np.isfinite(result1)) statistic.reset() - statistic.update(ParamsMap()) - tools.update(ParamsMap()) + statistic.update(params) + tools.update(params) result2 = statistic.compute_theory_vector(tools) assert np.array_equal(result1, result2) @@ -223,20 +223,21 @@ def test_two_point_src0_src0_no_window(sacc_galaxy_cells_src0_src0_no_window) -> statistic.read(sacc_data) tools = ModelingTools() - tools.update(ParamsMap()) - tools.prepare(pyccl.CosmologyVanillaLCDM()) + params = get_default_params_map(tools) + tools.update(params) + tools.prepare() assert statistic.window is None statistic.reset() - statistic.update(ParamsMap()) - tools.update(ParamsMap()) + statistic.update(params) + tools.update(params) result1 = statistic.compute_theory_vector(tools) assert all(np.isfinite(result1)) statistic.reset() - statistic.update(ParamsMap()) - tools.update(ParamsMap()) + statistic.update(params) + tools.update(params) result2 = statistic.compute_theory_vector(tools) assert np.array_equal(result1, result2) @@ -360,8 +361,9 @@ def test_two_point_src0_src0_cuts(sacc_galaxy_cells_src0_src0) -> None: statistic.read(sacc_data) tools = ModelingTools() - tools.update(ParamsMap()) - tools.prepare(pyccl.CosmologyVanillaLCDM()) + params = get_default_params_map(tools) + tools.update(params) + tools.prepare() assert statistic.window is None assert statistic.ells is not None @@ -369,7 +371,7 @@ def test_two_point_src0_src0_cuts(sacc_galaxy_cells_src0_src0) -> None: assert all(statistic.ells >= 50) assert all(statistic.ells <= 200) - statistic.update(ParamsMap()) + statistic.update(params) statistic.compute_theory_vector(tools) @@ -383,10 +385,11 @@ def test_two_point_lens0_lens0_cuts(sacc_galaxy_xis_lens0_lens0) -> None: ) statistic.read(sacc_data) - param_map = ParamsMap({"lens0_bias": 1.0}) tools = ModelingTools() - tools.update(param_map) - tools.prepare(pyccl.CosmologyVanillaLCDM()) + params = get_default_params_map(tools) + params.update({"lens0_bias": 1.0}) + tools.update(params) + tools.prepare() assert statistic.window is None assert statistic.thetas is not None @@ -394,7 +397,7 @@ def test_two_point_lens0_lens0_cuts(sacc_galaxy_xis_lens0_lens0) -> None: assert all(statistic.thetas >= 0.1) assert all(statistic.thetas <= 0.5) - statistic.update(param_map) + statistic.update(params) statistic.compute_theory_vector(tools) @@ -411,10 +414,11 @@ def test_two_point_lens0_lens0_config(sacc_galaxy_xis_lens0_lens0) -> None: ) statistic.read(sacc_data) - param_map = ParamsMap({"lens0_bias": 1.0}) tools = ModelingTools() - tools.update(param_map) - tools.prepare(pyccl.CosmologyVanillaLCDM()) + params = get_default_params_map(tools) + params.update({"lens0_bias": 1.0}) + tools.update(params) + tools.prepare() assert statistic.window is None assert statistic.thetas is not None @@ -425,7 +429,7 @@ def test_two_point_lens0_lens0_config(sacc_galaxy_xis_lens0_lens0) -> None: # on how many unique ells we get from the log-binning. assert len(statistic.ells_for_xi) == 175 - statistic.update(param_map) + statistic.update(params) statistic.compute_theory_vector(tools) @@ -637,13 +641,15 @@ def test_from_measurement_compute_theory_vector_window( assert isinstance(two_point_with_window, TwoPoint) assert two_point_with_window.ready - req_params = two_point_with_window.required_parameters() + tools = ModelingTools() + req_params = ( + two_point_with_window.required_parameters() + tools.required_parameters() + ) default_values = req_params.get_default_values() params = ParamsMap(default_values) - tools = ModelingTools() tools.update(params) - tools.prepare(pyccl.CosmologyVanillaLCDM()) + tools.prepare() two_point_with_window.update(params) prediction = two_point_with_window.compute_theory_vector(tools) @@ -661,14 +667,15 @@ def test_from_measurement_compute_theory_vector_window_check( assert isinstance(two_point_without_window, TwoPoint) assert two_point_without_window.ready - req_params = two_point_with_window.required_parameters() + tools = ModelingTools() + req_params = ( + two_point_with_window.required_parameters() + tools.required_parameters() + ) default_values = req_params.get_default_values() params = ParamsMap(default_values) - tools = ModelingTools() tools.update(params) - tools.prepare(pyccl.CosmologyVanillaLCDM()) - + tools.prepare() two_point_with_window.update(params) two_point_without_window.update(params) diff --git a/tests/likelihood/lkdir/lk_derived_parameter.py b/tests/likelihood/lkdir/lk_derived_parameter.py index 9b490f89e..78d90e20e 100644 --- a/tests/likelihood/lkdir/lk_derived_parameter.py +++ b/tests/likelihood/lkdir/lk_derived_parameter.py @@ -4,9 +4,14 @@ """ from firecrown.likelihood.likelihood import NamedParameters +from firecrown.modeling_tools import ModelingTools +from firecrown.ccl_factory import CCLFactory, PoweSpecAmplitudeParameter from . import lkmodule def build_likelihood(_: NamedParameters): """Return a DerivedParameterLikelihood object.""" - return lkmodule.derived_parameter_likelihood() + tools = ModelingTools( + ccl_factory=CCLFactory(amplitude_parameter=PoweSpecAmplitudeParameter.SIGMA8) + ) + return lkmodule.derived_parameter_likelihood(), tools diff --git a/tests/likelihood/lkdir/lk_needing_param.py b/tests/likelihood/lkdir/lk_needing_param.py index e228cfdd4..28e9f745d 100644 --- a/tests/likelihood/lkdir/lk_needing_param.py +++ b/tests/likelihood/lkdir/lk_needing_param.py @@ -4,9 +4,15 @@ """ from firecrown.likelihood.likelihood import NamedParameters +from firecrown.modeling_tools import ModelingTools +from firecrown.ccl_factory import CCLFactory, PoweSpecAmplitudeParameter from . import lkmodule def build_likelihood(params: NamedParameters): """Return a ParameterizedLikelihood object.""" - return lkmodule.parameterized_likelihood(params) + tools = ModelingTools( + ccl_factory=CCLFactory(amplitude_parameter=PoweSpecAmplitudeParameter.SIGMA8) + ) + + return lkmodule.parameterized_likelihood(params), tools diff --git a/tests/likelihood/lkdir/lk_sampler_parameter.py b/tests/likelihood/lkdir/lk_sampler_parameter.py index 7cf9a2d37..f89bba054 100644 --- a/tests/likelihood/lkdir/lk_sampler_parameter.py +++ b/tests/likelihood/lkdir/lk_sampler_parameter.py @@ -5,9 +5,14 @@ """ from firecrown.likelihood.likelihood import NamedParameters +from firecrown.modeling_tools import ModelingTools +from firecrown.ccl_factory import CCLFactory, PoweSpecAmplitudeParameter from . import lkmodule def build_likelihood(params: NamedParameters): """Return a SamplerParameterLikelihood object.""" - return lkmodule.sampler_parameter_likelihood(params) + tools = ModelingTools( + ccl_factory=CCLFactory(amplitude_parameter=PoweSpecAmplitudeParameter.SIGMA8) + ) + return lkmodule.sampler_parameter_likelihood(params), tools diff --git a/tests/likelihood/lkdir/lkscript.py b/tests/likelihood/lkdir/lkscript.py index 011abc77d..7478afa3c 100644 --- a/tests/likelihood/lkdir/lkscript.py +++ b/tests/likelihood/lkdir/lkscript.py @@ -3,11 +3,14 @@ """ from firecrown.modeling_tools import ModelingTools +from firecrown.ccl_factory import CCLFactory, PoweSpecAmplitudeParameter from . import lkmodule def build_likelihood(_): """Return an EmptyLikelihood object.""" - tools = ModelingTools() + tools = ModelingTools( + ccl_factory=CCLFactory(amplitude_parameter=PoweSpecAmplitudeParameter.SIGMA8) + ) tools.test_attribute = "test" # type: ignore[unused-ignore] return (lkmodule.empty_likelihood(), tools) diff --git a/tests/likelihood/test_factories.py b/tests/likelihood/test_factories.py new file mode 100644 index 000000000..45c15da1c --- /dev/null +++ b/tests/likelihood/test_factories.py @@ -0,0 +1,315 @@ +"""Tests for the module firecrown.likelihood.factories. +""" + +import re +from pathlib import Path +import pytest + +import sacc + +from firecrown.likelihood.factories import ( + build_two_point_likelihood, + DataSourceSacc, + TwoPointCorrelationSpace, + TwoPointExperiment, + TwoPointFactory, +) +from firecrown.likelihood.weak_lensing import WeakLensingFactory +from firecrown.likelihood.number_counts import NumberCountsFactory +from firecrown.likelihood.likelihood import Likelihood, NamedParameters +from firecrown.modeling_tools import ModelingTools + + +def test_two_point_factory_dict() -> None: + two_point_factory_dict = { + "correlation_space": TwoPointCorrelationSpace.HARMONIC, + "weak_lensing_factory": {"per_bin_systematics": [], "global_systematics": []}, + "number_counts_factory": {"per_bin_systematics": [], "global_systematics": []}, + } + two_point_factory = TwoPointFactory.model_validate(two_point_factory_dict) + assert isinstance(two_point_factory, TwoPointFactory) + assert isinstance(two_point_factory.weak_lensing_factory, WeakLensingFactory) + assert isinstance(two_point_factory.number_counts_factory, NumberCountsFactory) + assert two_point_factory.correlation_space == TwoPointCorrelationSpace.HARMONIC + assert two_point_factory.weak_lensing_factory.per_bin_systematics == [] + assert two_point_factory.weak_lensing_factory.global_systematics == [] + assert two_point_factory.number_counts_factory.per_bin_systematics == [] + assert two_point_factory.number_counts_factory.global_systematics == [] + + +def test_two_point_factor_direct() -> None: + two_point_factory = TwoPointFactory( + correlation_space=TwoPointCorrelationSpace.HARMONIC, + weak_lensing_factory=WeakLensingFactory( + per_bin_systematics=[], global_systematics=[] + ), + number_counts_factory=NumberCountsFactory( + per_bin_systematics=[], global_systematics=[] + ), + ) + + assert two_point_factory.correlation_space == TwoPointCorrelationSpace.HARMONIC + assert two_point_factory.weak_lensing_factory.per_bin_systematics == [] + assert two_point_factory.weak_lensing_factory.global_systematics == [] + assert two_point_factory.number_counts_factory.per_bin_systematics == [] + assert two_point_factory.number_counts_factory.global_systematics == [] + + +def test_two_point_factor_invalid_correlation_space_type() -> None: + two_point_factory_dict = { + "correlation_space": 1.2, + "weak_lensing_factory": {"per_bin_systematics": [], "global_systematics": []}, + "number_counts_factory": {"per_bin_systematics": [], "global_systematics": []}, + } + with pytest.raises( + ValueError, + match=re.compile( + ".*validation error for TwoPointFactory.*correlation_space.*", re.DOTALL + ), + ): + TwoPointFactory.model_validate(two_point_factory_dict) + + +def test_two_point_factor_invalid_correlation_space_option() -> None: + two_point_factory_dict = { + "correlation_space": "invalid", + "weak_lensing_factory": {"per_bin_systematics": [], "global_systematics": []}, + "number_counts_factory": {"per_bin_systematics": [], "global_systematics": []}, + } + with pytest.raises( + ValueError, + match=re.compile( + ".*validation error for TwoPointFactory.*correlation_space.*", re.DOTALL + ), + ): + TwoPointFactory.model_validate(two_point_factory_dict) + + +def test_data_source_sacc_dict() -> None: + data_source_sacc_dict = {"sacc_data_file": "tests/bug_398.sacc.gz"} + data_source_sacc = DataSourceSacc.model_validate(data_source_sacc_dict) + assert isinstance(data_source_sacc, DataSourceSacc) + assert data_source_sacc.sacc_data_file == "tests/bug_398.sacc.gz" + + +def test_data_source_sacc_direct() -> None: + data_source_sacc = DataSourceSacc(sacc_data_file="tests/bug_398.sacc.gz") + assert data_source_sacc.sacc_data_file == "tests/bug_398.sacc.gz" + + +def test_data_source_sacc_invalid_file() -> None: + data_source_sacc_dict = {"sacc_data_file": "tests/file_not_found.sacc.gz"} + with pytest.raises( + FileNotFoundError, match=".*File tests/file_not_found.sacc.gz does not exist.*" + ): + DataSourceSacc.model_validate(data_source_sacc_dict) + + +def test_data_source_sacc_get_sacc_data() -> None: + data_source_sacc_dict = {"sacc_data_file": "tests/bug_398.sacc.gz"} + data_source_sacc = DataSourceSacc.model_validate(data_source_sacc_dict) + assert isinstance(data_source_sacc, DataSourceSacc) + assert data_source_sacc.sacc_data_file == "tests/bug_398.sacc.gz" + sacc_data = data_source_sacc.get_sacc_data() + assert sacc_data is not None + assert isinstance(sacc_data, sacc.Sacc) + + +def test_two_point_experiment_dict() -> None: + two_point_experiment_dict = { + "two_point_factory": { + "correlation_space": TwoPointCorrelationSpace.HARMONIC, + "weak_lensing_factory": { + "per_bin_systematics": [], + "global_systematics": [], + }, + "number_counts_factory": { + "per_bin_systematics": [], + "global_systematics": [], + }, + }, + "data_source": {"sacc_data_file": "tests/bug_398.sacc.gz"}, + } + two_point_experiment = TwoPointExperiment.model_validate(two_point_experiment_dict) + assert isinstance(two_point_experiment, TwoPointExperiment) + assert isinstance(two_point_experiment.two_point_factory, TwoPointFactory) + assert isinstance(two_point_experiment.data_source, DataSourceSacc) + assert ( + two_point_experiment.two_point_factory.correlation_space + == TwoPointCorrelationSpace.HARMONIC + ) + assert ( + two_point_experiment.two_point_factory.weak_lensing_factory.per_bin_systematics + == [] + ) + assert ( + two_point_experiment.two_point_factory.weak_lensing_factory.global_systematics + == [] + ) + assert ( + two_point_experiment.two_point_factory.number_counts_factory.per_bin_systematics + == [] + ) + assert ( + two_point_experiment.two_point_factory.number_counts_factory.global_systematics + == [] + ) + assert two_point_experiment.data_source.sacc_data_file == "tests/bug_398.sacc.gz" + + +def test_two_point_experiment_direct() -> None: + two_point_experiment = TwoPointExperiment( + two_point_factory=TwoPointFactory( + correlation_space=TwoPointCorrelationSpace.HARMONIC, + weak_lensing_factory=WeakLensingFactory( + per_bin_systematics=[], global_systematics=[] + ), + number_counts_factory=NumberCountsFactory( + per_bin_systematics=[], global_systematics=[] + ), + ), + data_source=DataSourceSacc(sacc_data_file="tests/bug_398.sacc.gz"), + ) + assert isinstance(two_point_experiment, TwoPointExperiment) + assert isinstance(two_point_experiment.two_point_factory, TwoPointFactory) + assert isinstance(two_point_experiment.data_source, DataSourceSacc) + assert ( + two_point_experiment.two_point_factory.correlation_space + == TwoPointCorrelationSpace.HARMONIC + ) + assert ( + two_point_experiment.two_point_factory.weak_lensing_factory.per_bin_systematics + == [] + ) + assert ( + two_point_experiment.two_point_factory.weak_lensing_factory.global_systematics + == [] + ) + assert ( + two_point_experiment.two_point_factory.number_counts_factory.per_bin_systematics + == [] + ) + assert ( + two_point_experiment.two_point_factory.number_counts_factory.global_systematics + == [] + ) + assert two_point_experiment.data_source.sacc_data_file == "tests/bug_398.sacc.gz" + + +def test_build_two_point_likelihood_real(tmp_path: Path) -> None: + tmp_experiment_file = tmp_path / "experiment.yaml" + tmp_experiment_file.write_text( + """ +two_point_factory: + correlation_space: real + weak_lensing_factory: + per_bin_systematics: [] + global_systematics: [] + number_counts_factory: + per_bin_systematics: [] + global_systematics: [] +data_source: + sacc_data_file: examples/des_y1_3x2pt/des_y1_3x2pt_sacc_data.fits +""" + ) + + build_parameters = NamedParameters({"likelihood_config": str(tmp_experiment_file)}) + likelihood, tools = build_two_point_likelihood(build_parameters) + assert isinstance(likelihood, Likelihood) + assert isinstance(tools, ModelingTools) + + +def test_build_two_point_likelihood_harmonic(tmp_path: Path) -> None: + tmp_experiment_file = tmp_path / "experiment.yaml" + tmp_experiment_file.write_text( + """ +two_point_factory: + correlation_space: harmonic + weak_lensing_factory: + per_bin_systematics: [] + global_systematics: [] + number_counts_factory: + per_bin_systematics: [] + global_systematics: [] +data_source: + sacc_data_file: tests/bug_398.sacc.gz +""" + ) + + build_parameters = NamedParameters({"likelihood_config": str(tmp_experiment_file)}) + likelihood, tools = build_two_point_likelihood(build_parameters) + assert isinstance(likelihood, Likelihood) + assert isinstance(tools, ModelingTools) + + +def test_build_two_point_likelihood_real_no_real_data(tmp_path: Path) -> None: + tmp_experiment_file = tmp_path / "experiment.yaml" + tmp_experiment_file.write_text( + """ +two_point_factory: + correlation_space: real + weak_lensing_factory: + per_bin_systematics: [] + global_systematics: [] + number_counts_factory: + per_bin_systematics: [] + global_systematics: [] +data_source: + sacc_data_file: tests/bug_398.sacc.gz +""" + ) + + build_parameters = NamedParameters({"likelihood_config": str(tmp_experiment_file)}) + with pytest.raises( + ValueError, + match=re.compile( + "No two-point measurements in real space found in the SACC file.", re.DOTALL + ), + ): + _ = build_two_point_likelihood(build_parameters) + + +def test_build_two_point_likelihood_harmonic_no_harmonic_data(tmp_path: Path) -> None: + tmp_experiment_file = tmp_path / "experiment.yaml" + tmp_experiment_file.write_text( + """ +two_point_factory: + correlation_space: harmonic + weak_lensing_factory: + per_bin_systematics: [] + global_systematics: [] + number_counts_factory: + per_bin_systematics: [] + global_systematics: [] +data_source: + sacc_data_file: examples/des_y1_3x2pt/des_y1_3x2pt_sacc_data.fits +""" + ) + + build_parameters = NamedParameters({"likelihood_config": str(tmp_experiment_file)}) + with pytest.raises( + ValueError, + match=re.compile( + "No two-point measurements in harmonic space found in the SACC file.", + re.DOTALL, + ), + ): + _ = build_two_point_likelihood(build_parameters) + + +def test_build_two_point_likelihood_empty_likelihood_config(tmp_path: Path) -> None: + tmp_experiment_file = tmp_path / "experiment.yaml" + tmp_experiment_file.write_text("") + + build_parameters = NamedParameters({"likelihood_config": str(tmp_experiment_file)}) + with pytest.raises(ValueError, match="No likelihood config found."): + _ = build_two_point_likelihood(build_parameters) + + +def test_build_two_point_likelihood_invalid_likelihood_config(tmp_path: Path) -> None: + tmp_experiment_file = tmp_path / "experiment.yaml" + tmp_experiment_file.write_text("I'm not a valid YAML file.") + + build_parameters = NamedParameters({"likelihood_config": str(tmp_experiment_file)}) + with pytest.raises(ValueError, match=".*validation error for TwoPointExperiment.*"): + _ = build_two_point_likelihood(build_parameters) diff --git a/tests/likelihood/test_likelihood.py b/tests/likelihood/test_likelihood.py index f811d6064..92e10291c 100644 --- a/tests/likelihood/test_likelihood.py +++ b/tests/likelihood/test_likelihood.py @@ -124,7 +124,9 @@ def test_load_likelihood_from_module(): def test_load_likelihood_from_module_not_in_path(): with pytest.raises( ValueError, - match="Unrecognized Firecrown initialization module lkscript_not_in_path.", + match="Unrecognized Firecrown initialization module 'lkscript_not_in_path'. " + "The module must be either a module_name or a module_name.func where func " + "is the factory function.", ): _ = load_likelihood_from_module("lkscript_not_in_path", NamedParameters()) @@ -136,3 +138,16 @@ def test_load_likelihood_not_in_path(): "lkscript_not_in_path.", ): _ = load_likelihood("lkscript_not_in_path", NamedParameters()) + + +def test_load_likelihood_from_module_function_missing_likelihood_config(): + dir_path = os.path.dirname(os.path.realpath(__file__)) + module_path = os.path.join(dir_path, "lkdir") + + sys.path.append(module_path) + + with pytest.raises(KeyError, match="likelihood_config"): + _ = load_likelihood_from_module( + "firecrown.likelihood.factories.build_two_point_likelihood", + NamedParameters(), + ) diff --git a/tests/test_bug398.py b/tests/test_bug398.py index bf99b6fbc..266abbcca 100644 --- a/tests/test_bug398.py +++ b/tests/test_bug398.py @@ -2,16 +2,16 @@ import os import sacc -import pyccl from numpy.testing import assert_allclose +from firecrown.updatable import get_default_params_map import firecrown.likelihood.weak_lensing as wl from firecrown.likelihood.two_point import TwoPoint from firecrown.likelihood.gaussian import ConstGaussian from firecrown.modeling_tools import ModelingTools from firecrown.likelihood.likelihood import Likelihood, NamedParameters -from firecrown.parameters import ParamsMap +from firecrown.ccl_factory import CCLFactory, PoweSpecAmplitudeParameter def build_likelihood( @@ -42,7 +42,9 @@ def build_likelihood( sacc_data_type="galaxy_density_cl", ) - modeling_tools = ModelingTools() + modeling_tools = ModelingTools( + ccl_factory=CCLFactory(amplitude_parameter=PoweSpecAmplitudeParameter.SIGMA8) + ) likelihood = ConstGaussian(statistics=[src2_src2, lens0_lens0, lens0_src2]) likelihood.read(sacc_data) @@ -119,7 +121,7 @@ def build_likelihood( ] -def test_broken_window_function(): +def test_broken_window_function() -> None: likelihood, modeling_tools = build_likelihood( build_parameters=NamedParameters({"sacc_data": SACC_FILE}) ) @@ -127,13 +129,13 @@ def test_broken_window_function(): assert modeling_tools is not None -def test_eval_cl_window_src2_src2(): - tools = ModelingTools() - cosmo = pyccl.CosmologyVanillaLCDM() - params = ParamsMap() - +def test_eval_cl_window_src2_src2() -> None: + tools = ModelingTools( + ccl_factory=CCLFactory(amplitude_parameter=PoweSpecAmplitudeParameter.SIGMA8) + ) + params = get_default_params_map(tools) tools.update(params) - tools.prepare(cosmo) + tools.prepare() sacc_data = sacc.Sacc.load_fits(SACC_FILE) src2 = wl.WeakLensing(sacc_tracer="src2") @@ -151,13 +153,13 @@ def test_eval_cl_window_src2_src2(): assert_allclose(theory_vector, SRC2_SRC2_CL_VANILLA_LCDM, rtol=1e-8) -def test_eval_cl_window_lens0_src2(): - tools = ModelingTools() - cosmo = pyccl.CosmologyVanillaLCDM() - params = ParamsMap() - +def test_eval_cl_window_lens0_src2() -> None: + tools = ModelingTools( + ccl_factory=CCLFactory(amplitude_parameter=PoweSpecAmplitudeParameter.SIGMA8) + ) + params = get_default_params_map(tools) tools.update(params) - tools.prepare(cosmo) + tools.prepare() sacc_data = sacc.Sacc.load_fits(SACC_FILE) src2 = wl.WeakLensing(sacc_tracer="src2") @@ -176,13 +178,13 @@ def test_eval_cl_window_lens0_src2(): assert_allclose(theory_vector, LENS0_SRC2_CL_VANILLA_LCDM, rtol=1e-8) -def test_eval_cl_window_lens0_lens0(): - tools = ModelingTools() - cosmo = pyccl.CosmologyVanillaLCDM() - params = ParamsMap() - +def test_eval_cl_window_lens0_lens0() -> None: + tools = ModelingTools( + ccl_factory=CCLFactory(amplitude_parameter=PoweSpecAmplitudeParameter.SIGMA8) + ) + params = get_default_params_map(tools) tools.update(params) - tools.prepare(cosmo) + tools.prepare() sacc_data = sacc.Sacc.load_fits(SACC_FILE) lens0 = wl.WeakLensing(sacc_tracer="lens0") @@ -200,13 +202,13 @@ def test_eval_cl_window_lens0_lens0(): assert_allclose(theory_vector, LENS0_LENS0_CL_VANILLA_LCDM, rtol=1e-8) -def test_compute_likelihood_src0_src0(): - tools = ModelingTools() - cosmo = pyccl.CosmologyVanillaLCDM() - params = ParamsMap() - +def test_compute_likelihood_src0_src0() -> None: + tools = ModelingTools( + ccl_factory=CCLFactory(amplitude_parameter=PoweSpecAmplitudeParameter.SIGMA8) + ) + params = get_default_params_map(tools) tools.update(params) - tools.prepare(cosmo) + tools.prepare() sacc_data = sacc.Sacc.load_fits(SACC_FILE) src0 = wl.WeakLensing(sacc_tracer="src0") diff --git a/tests/test_bug_413.py b/tests/test_bug_413.py index 1091a7fcc..16bdfcb65 100644 --- a/tests/test_bug_413.py +++ b/tests/test_bug_413.py @@ -4,12 +4,12 @@ import pyccl import sacc +from firecrown.updatable import get_default_params_map from firecrown.likelihood.number_counts import ( NumberCounts, PTNonLinearBiasSystematic, ) from firecrown.likelihood.two_point import TwoPoint -from firecrown.parameters import ParamsMap from firecrown.modeling_tools import ModelingTools @@ -42,17 +42,6 @@ def make_twopoint_with_optional_systematics( ) statistic.read(sacc_data) # Note the two sources have identical parameters. - param_map = ParamsMap( - { - "lens0_bias": 1.1, - "lens0_b_2": 1.05, - "lens0_b_s": 0.99, - "lens1_bias": 1.1, - "lens1_b_2": 1.05, - "lens1_b_s": 0.99, - } - ) - statistic.update(param_map) tools = ModelingTools( pt_calculator=pyccl.nl_pt.EulerianPTCalculator( with_NC=True, @@ -62,8 +51,21 @@ def make_twopoint_with_optional_systematics( nk_per_decade=20, ) ) - tools.update(param_map) - tools.prepare(pyccl.CosmologyVanillaLCDM()) + params = get_default_params_map(tools) + params.update( + { + "lens0_bias": 1.1, + "lens0_b_2": 1.05, + "lens0_b_s": 0.99, + "lens1_bias": 1.1, + "lens1_b_2": 1.05, + "lens1_b_s": 0.99, + } + ) + + tools.update(params) + tools.prepare() + statistic.update(params) _ = a.get_tracers(tools) _ = b.get_tracers(tools) assert a.tracers[0].has_pt == first_source_has_systematic diff --git a/tests/test_ccl_factory.py b/tests/test_ccl_factory.py new file mode 100644 index 000000000..3e12d6140 --- /dev/null +++ b/tests/test_ccl_factory.py @@ -0,0 +1,367 @@ +"""Test the CCLFactory object.""" + +import numpy as np +import pytest +import pyccl +from pyccl.neutrinos import NeutrinoMassSplits + +from firecrown.ccl_factory import ( + CCLFactory, + PoweSpecAmplitudeParameter, + CCLCalculatorArgs, +) +from firecrown.updatable import get_default_params_map +from firecrown.parameters import ParamsMap +from firecrown.utils import base_model_from_yaml, base_model_to_yaml + + +@pytest.fixture(name="amplitude_parameter", params=list(PoweSpecAmplitudeParameter)) +def fixture_amplitude_parameter(request): + return request.param + + +@pytest.fixture(name="neutrino_mass_splits", params=list(NeutrinoMassSplits)) +def fixture_neutrino_mass_splits(request): + return request.param + + +@pytest.fixture( + name="require_nonlinear_pk", + params=[True, False], + ids=["require_nonlinear_pk", "no_require_nonlinear_pk"], +) +def fixture_require_nonlinear_pk(request): + return request.param + + +Z_ARRAY = np.linspace(0.0, 5.0, 100) +A_ARRAY = 1.0 / (1.0 + np.flip(Z_ARRAY)) +K_ARRAY = np.geomspace(1.0e-5, 10.0, 100) +A_GRID, K_GRID = np.meshgrid(A_ARRAY, K_ARRAY, indexing="ij") + +CHI_ARRAY = np.linspace(100.0, 0.0, 100) +H_OVER_H0_ARRAY = np.linspace(1.0, 100.0, 100) +# Simple power spectrum model for testing +PK_ARRAY = ( + 2.0e-9 + * (K_GRID / 0.05) ** (0.96) + * (np.log(1 + 2.34 * K_GRID / 0.30) / (2.34 * K_GRID / 0.30)) + * (1 / (1 + (K_GRID / (0.1 * 0.05)) ** (3 * 0.96)) * A_GRID) +) + +BACKGROUND = { + "background": {"a": A_ARRAY, "chi": CHI_ARRAY, "h_over_h0": H_OVER_H0_ARRAY} +} + +PK_LINEAR = { + "pk_linear": {"k": K_ARRAY, "a": A_ARRAY, "delta_matter:delta_matter": PK_ARRAY} +} + +PK_NONLIN = { + "pk_nonlin": { + "k": np.linspace(0.1, 1.0, 100), + "a": np.linspace(0.1, 1.0, 100), + "delta_matter:delta_matter": PK_ARRAY, + } +} + + +@pytest.fixture( + name="calculator_args", + params=[BACKGROUND, BACKGROUND | PK_LINEAR, BACKGROUND | PK_LINEAR | PK_NONLIN], + ids=["background", "background_linear_pk", "background_linear_pk_nonlinear_pk"], +) +def fixture_calculator_args(request) -> CCLCalculatorArgs: + return request.param + + +def test_ccl_factory_simple( + amplitude_parameter: PoweSpecAmplitudeParameter, + neutrino_mass_splits: NeutrinoMassSplits, + require_nonlinear_pk: bool, +) -> None: + ccl_factory = CCLFactory( + amplitude_parameter=amplitude_parameter, + mass_split=neutrino_mass_splits, + require_nonlinear_pk=require_nonlinear_pk, + ) + + assert ccl_factory is not None + assert ccl_factory.amplitude_parameter == amplitude_parameter + assert ccl_factory.mass_split == neutrino_mass_splits + assert ccl_factory.require_nonlinear_pk == require_nonlinear_pk + + default_params = get_default_params_map(ccl_factory) + + ccl_factory.update(default_params) + + cosmo = ccl_factory.create() + + assert cosmo is not None + assert isinstance(cosmo, pyccl.Cosmology) + + +def test_ccl_factory_ccl_args( + amplitude_parameter: PoweSpecAmplitudeParameter, + neutrino_mass_splits: NeutrinoMassSplits, + require_nonlinear_pk: bool, + calculator_args: CCLCalculatorArgs, +) -> None: + ccl_factory = CCLFactory( + amplitude_parameter=amplitude_parameter, + mass_split=neutrino_mass_splits, + require_nonlinear_pk=require_nonlinear_pk, + ) + + if require_nonlinear_pk and "pk_linear" not in calculator_args: + pytest.skip( + "Nonlinear power spectrum requested but " + "linear power spectrum not provided." + ) + + assert ccl_factory is not None + assert ccl_factory.amplitude_parameter == amplitude_parameter + assert ccl_factory.mass_split == neutrino_mass_splits + assert ccl_factory.require_nonlinear_pk == require_nonlinear_pk + + default_params = get_default_params_map(ccl_factory) + + ccl_factory.update(default_params) + + cosmo = ccl_factory.create(calculator_args) + + assert cosmo is not None + assert isinstance(cosmo, pyccl.Cosmology) + + +def test_ccl_factory_update() -> None: + ccl_factory = CCLFactory() + + assert ccl_factory is not None + + default_params = get_default_params_map(ccl_factory) + + ccl_factory.update(default_params) + + cosmo = ccl_factory.create() + + assert cosmo is not None + assert isinstance(cosmo, pyccl.Cosmology) + + new_params = ParamsMap(default_params.copy()) + new_params["Omega_c"] = 0.1 + + ccl_factory.reset() + ccl_factory.update(new_params) + + cosmo = ccl_factory.create() + + assert cosmo is not None + assert isinstance(cosmo, pyccl.Cosmology) + + +def test_ccl_factory_amplitude_parameter( + amplitude_parameter: PoweSpecAmplitudeParameter, +) -> None: + ccl_factory = CCLFactory(amplitude_parameter=amplitude_parameter) + + assert ccl_factory is not None + + default_params = get_default_params_map(ccl_factory) + + ccl_factory.update(default_params) + + cosmo = ccl_factory.create() + + assert cosmo is not None + assert isinstance(cosmo, pyccl.Cosmology) + + +def test_ccl_factory_neutrino_mass_splits( + neutrino_mass_splits: NeutrinoMassSplits, +) -> None: + ccl_factory = CCLFactory(neutrino_mass_splits=neutrino_mass_splits) + + assert ccl_factory is not None + + default_params = get_default_params_map(ccl_factory) + + ccl_factory.update(default_params) + + cosmo = ccl_factory.create() + + assert cosmo is not None + assert isinstance(cosmo, pyccl.Cosmology) + + +def test_ccl_factory_get() -> None: + ccl_factory = CCLFactory() + + assert ccl_factory is not None + + default_params = get_default_params_map(ccl_factory) + + ccl_factory.update(default_params) + + cosmo = ccl_factory.create() + + assert cosmo is not None + assert isinstance(cosmo, pyccl.Cosmology) + + cosmo2 = ccl_factory.get() + + assert cosmo2 is not None + assert isinstance(cosmo2, pyccl.Cosmology) + + assert cosmo is cosmo2 + + +def test_ccl_factory_get_not_created() -> None: + ccl_factory = CCLFactory() + + with pytest.raises(ValueError, match="CCLFactory object has not been created yet."): + ccl_factory.get() + + +def test_ccl_factory_create_twice() -> None: + ccl_factory = CCLFactory() + + assert ccl_factory is not None + + default_params = get_default_params_map(ccl_factory) + + ccl_factory.update(default_params) + + cosmo = ccl_factory.create() + + assert cosmo is not None + assert isinstance(cosmo, pyccl.Cosmology) + + with pytest.raises(ValueError, match="CCLFactory object has already been created."): + ccl_factory.create() + + +def test_ccl_factory_create_not_updated() -> None: + ccl_factory = CCLFactory() + + assert ccl_factory is not None + + with pytest.raises(ValueError, match="Parameters have not been updated yet."): + ccl_factory.create() + + +def test_ccl_factory_tofrom_yaml() -> None: + ccl_factory = CCLFactory() + + assert ccl_factory is not None + + default_params = get_default_params_map(ccl_factory) + + ccl_factory.update(default_params) + + cosmo = ccl_factory.create() + + assert cosmo is not None + assert isinstance(cosmo, pyccl.Cosmology) + + yaml_str = base_model_to_yaml(ccl_factory) + + ccl_factory2 = base_model_from_yaml(CCLFactory, yaml_str) + + assert ccl_factory2 is not None + ccl_factory2.update(default_params) + + cosmo2 = ccl_factory2.create() + + assert cosmo2 is not None + assert isinstance(cosmo2, pyccl.Cosmology) + + assert cosmo == cosmo2 + + +def test_ccl_factory_tofrom_yaml_all_options( + amplitude_parameter: PoweSpecAmplitudeParameter, + neutrino_mass_splits: NeutrinoMassSplits, + require_nonlinear_pk: bool, +) -> None: + ccl_factory = CCLFactory( + amplitude_parameter=amplitude_parameter, + mass_split=neutrino_mass_splits, + require_nonlinear_pk=require_nonlinear_pk, + ) + + assert ccl_factory is not None + assert ccl_factory.amplitude_parameter == amplitude_parameter + assert ccl_factory.mass_split == neutrino_mass_splits + assert ccl_factory.require_nonlinear_pk == require_nonlinear_pk + + default_params = get_default_params_map(ccl_factory) + + ccl_factory.update(default_params) + + cosmo = ccl_factory.create() + + assert cosmo is not None + assert isinstance(cosmo, pyccl.Cosmology) + + yaml_str = base_model_to_yaml(ccl_factory) + + ccl_factory2 = base_model_from_yaml(CCLFactory, yaml_str) + + assert ccl_factory2 is not None + assert ccl_factory2.amplitude_parameter == amplitude_parameter + assert ccl_factory2.mass_split == neutrino_mass_splits + assert ccl_factory2.require_nonlinear_pk == require_nonlinear_pk + ccl_factory2.update(default_params) + + cosmo2 = ccl_factory2.create() + + assert cosmo2 is not None + assert isinstance(cosmo2, pyccl.Cosmology) + + assert cosmo == cosmo2 + + +def test_ccl_factory_invalid_amplitude_parameter() -> None: + with pytest.raises( + ValueError, + match=".*Invalid value for PoweSpecAmplitudeParameter: Im not a valid value.*", + ): + CCLFactory(amplitude_parameter="Im not a valid value") + + +def test_ccl_factory_invalid_mass_splits() -> None: + with pytest.raises( + ValueError, + match=".*Invalid value for NeutrinoMassSplits: Im not a valid value.*", + ): + CCLFactory(mass_split="Im not a valid value") + + +def test_ccl_factory_from_dict() -> None: + ccl_factory_dict = { + "amplitude_parameter": PoweSpecAmplitudeParameter.SIGMA8, + "mass_split": NeutrinoMassSplits.EQUAL, + "require_nonlinear_pk": True, + } + + ccl_factory = CCLFactory.model_validate(ccl_factory_dict) + + assert ccl_factory is not None + assert ccl_factory.amplitude_parameter == PoweSpecAmplitudeParameter.SIGMA8 + assert ccl_factory.mass_split == NeutrinoMassSplits.EQUAL + assert ccl_factory.require_nonlinear_pk is True + + +def test_ccl_factory_from_dict_wrong_type() -> None: + ccl_factory_dict = { + "amplitude_parameter": 0.32, + "mass_split": NeutrinoMassSplits.EQUAL, + "require_nonlinear_pk": True, + } + + with pytest.raises( + ValueError, + match=".*Input should be.*", + ): + CCLFactory.model_validate(ccl_factory_dict) diff --git a/tests/test_modeling_tools.py b/tests/test_modeling_tools.py index 4ec3959e6..d546e8d78 100644 --- a/tests/test_modeling_tools.py +++ b/tests/test_modeling_tools.py @@ -4,8 +4,8 @@ import pytest import pyccl +from firecrown.updatable import get_default_params_map from firecrown.modeling_tools import ModelingTools, PowerspectrumModifier -from firecrown.parameters import ParamsMap @pytest.fixture(name="dummy_powerspectrum") @@ -15,7 +15,7 @@ def make_dummy_powerspectrum() -> pyccl.Pk2D: return pyccl.Pk2D.__new__(pyccl.Pk2D) -def test_default_constructed_state(): +def test_default_constructed_state() -> None: tools = ModelingTools() # Default constructed state is pretty barren... assert tools.ccl_cosmo is None @@ -24,23 +24,23 @@ def test_default_constructed_state(): assert len(tools.powerspectra) == 0 -def test_default_constructed_no_tools(): +def test_default_constructed_no_tools() -> None: tools = ModelingTools() with pytest.raises(RuntimeError): _ = tools.get_pk("nonesuch") -def test_adding_pk_and_getting(dummy_powerspectrum: pyccl.Pk2D): +def test_adding_pk_and_getting(dummy_powerspectrum: pyccl.Pk2D) -> None: tools = ModelingTools() tools.add_pk("silly", dummy_powerspectrum) - cosmo = pyccl.CosmologyVanillaLCDM() - tools.update(ParamsMap()) - tools.prepare(cosmo) + params = get_default_params_map(tools) + tools.update(params) + tools.prepare() assert tools.get_pk("silly") == dummy_powerspectrum -def test_no_adding_pk_twice(dummy_powerspectrum: pyccl.Pk2D): +def test_no_adding_pk_twice(dummy_powerspectrum: pyccl.Pk2D) -> None: tools = ModelingTools() tools.add_pk("silly", dummy_powerspectrum) assert tools.powerspectra["silly"] == dummy_powerspectrum @@ -48,38 +48,38 @@ def test_no_adding_pk_twice(dummy_powerspectrum: pyccl.Pk2D): tools.add_pk("silly", dummy_powerspectrum) -def test_modeling_tool_prepare_without_update(): +def test_modeling_tool_prepare_without_update() -> None: tools = ModelingTools() - cosmo = pyccl.CosmologyVanillaLCDM() with pytest.raises(RuntimeError, match="ModelingTools has not been updated."): - tools.prepare(cosmo) + tools.prepare() -def test_modeling_tool_preparing_twice(): +def test_modeling_tool_preparing_twice() -> None: tools = ModelingTools() - cosmo = pyccl.CosmologyVanillaLCDM() - tools.update(ParamsMap()) - tools.prepare(cosmo) + params = get_default_params_map(tools) + tools.update(params) + tools.prepare() with pytest.raises(RuntimeError, match="ModelingTools has already been prepared"): - tools.prepare(cosmo) + tools.prepare() -def test_modeling_tool_wrongly_setting_ccl_cosmo(): +def test_modeling_tool_wrongly_setting_ccl_cosmo() -> None: tools = ModelingTools() cosmo = pyccl.CosmologyVanillaLCDM() - tools.update(ParamsMap()) + params = get_default_params_map(tools) + tools.update(params) tools.ccl_cosmo = cosmo with pytest.raises(RuntimeError, match="Cosmology has already been set"): - tools.prepare(cosmo) + tools.prepare() -def test_modeling_tools_get_cosmo_without_setting(): +def test_modeling_tools_get_cosmo_without_setting() -> None: tools = ModelingTools() with pytest.raises(RuntimeError, match="Cosmology has not been set"): tools.get_ccl_cosmology() -def test_modeling_tools_get_pt_calculator_without_setting(): +def test_modeling_tools_get_pt_calculator_without_setting() -> None: tools = ModelingTools() with pytest.raises(RuntimeError, match="A PT calculator has not been set"): tools.get_pt_calculator() @@ -98,11 +98,11 @@ def compute_p_of_k_z(self, tools: ModelingTools) -> pyccl.Pk2D: return pyccl.Pk2D.__new__(pyccl.Pk2D) -def test_modeling_tools_add_pk_modifiers(): +def test_modeling_tools_add_pk_modifiers() -> None: tools = ModelingTools(pk_modifiers=[DummyPowerspectrumModifier()]) - cosmo = pyccl.CosmologyVanillaLCDM() - tools.update(ParamsMap()) - tools.prepare(cosmo) + params = get_default_params_map(tools) + tools.update(params) + tools.prepare() assert tools.get_pk("dummy:dummy") is not None assert isinstance(tools.get_pk("dummy:dummy"), pyccl.Pk2D) diff --git a/tests/test_parameters.py b/tests/test_parameters.py index 5041ee5f6..5982bdd56 100644 --- a/tests/test_parameters.py +++ b/tests/test_parameters.py @@ -103,7 +103,7 @@ def test_get_params_names_does_not_allow_mutation(): def test_params_map(): - my_params = ParamsMap({"a": 1}) + my_params = ParamsMap({"a": 1.0}) x = my_params.get_from_prefix_param(None, "a") assert x == 1 with pytest.raises(KeyError): @@ -112,6 +112,20 @@ def test_params_map(): _ = my_params.get_from_prefix_param(None, "no_such_name") +def test_params_map_wrong_type(): + with pytest.raises( + TypeError, match="Value for parameter a is not a float or a list of floats.*" + ): + _ = ParamsMap({"a": "not a float or a list of floats"}) + + +def test_params_map_wrong_type_list(): + with pytest.raises( + TypeError, match="Value for parameter a is not a float or a list of floats.*" + ): + _ = ParamsMap({"a": ["not a float or a list of floats"]}) + + def test_parameter_get_full_name_reject_empty_name(): with pytest.raises(ValueError): _ = parameter_get_full_name(None, "") diff --git a/tests/test_pt_systematics.py b/tests/test_pt_systematics.py index 727005c45..8857c84f6 100644 --- a/tests/test_pt_systematics.py +++ b/tests/test_pt_systematics.py @@ -12,6 +12,7 @@ import pyccl.nl_pt as pt import sacc +from firecrown.updatable import get_default_params_map import firecrown.likelihood.weak_lensing as wl import firecrown.likelihood.number_counts as nc from firecrown.likelihood.two_point import ( @@ -21,7 +22,7 @@ ) from firecrown.likelihood.gaussian import ConstGaussian from firecrown.modeling_tools import ModelingTools -from firecrown.parameters import ParamsMap +from firecrown.ccl_factory import CCLFactory, PoweSpecAmplitudeParameter @pytest.fixture(name="weak_lensing_source") @@ -85,9 +86,13 @@ def test_pt_systematics(weak_lensing_source, number_counts_source, sacc_data): nk_per_decade=4, cosmo=ccl_cosmo, ) - modeling_tools = ModelingTools(pt_calculator=pt_calculator) - modeling_tools.update(ParamsMap()) - modeling_tools.prepare(ccl_cosmo) + modeling_tools = ModelingTools( + pt_calculator=pt_calculator, + ccl_factory=CCLFactory(amplitude_parameter=PoweSpecAmplitudeParameter.SIGMA8), + ) + params = get_default_params_map(modeling_tools) + modeling_tools.update(params) + modeling_tools.prepare() # Bare CCL setup a_1 = 1.0 @@ -117,7 +122,7 @@ def test_pt_systematics(weak_lensing_source, number_counts_source, sacc_data): pk_gg = pt_calculator.get_biased_pk2d(tracer1=ptt_g, tracer2=ptt_g) # Set the parameters for our systematics - systematics_params = ParamsMap( + params.update( { "ia_a_1": a_1, "ia_a_2": a_2, @@ -132,7 +137,7 @@ def test_pt_systematics(weak_lensing_source, number_counts_source, sacc_data): ) # Apply the systematics parameters - likelihood.update(systematics_params) + likelihood.update(params) # Make things faster by only using a couple of ells for s in likelihood.statistics: @@ -219,7 +224,7 @@ def test_pt_systematics(weak_lensing_source, number_counts_source, sacc_data): cl_cs_theory = cl_GG + 2 * cl_GI + cl_II cl_gg_theory = cl_gg + 2 * cl_gm + cl_mm - print("IDS: ", id(s0), id(s1)) + # print("IDS: ", id(s0), id(s1)) assert np.allclose(cells_GG, cells_GG_m, atol=0, rtol=1e-127) @@ -282,9 +287,13 @@ def test_pt_mixed_systematics(sacc_data): nk_per_decade=4, cosmo=ccl_cosmo, ) - modeling_tools = ModelingTools(pt_calculator=pt_calculator) - modeling_tools.update(ParamsMap()) - modeling_tools.prepare(ccl_cosmo) + modeling_tools = ModelingTools( + pt_calculator=pt_calculator, + ccl_factory=CCLFactory(amplitude_parameter=PoweSpecAmplitudeParameter.SIGMA8), + ) + params = get_default_params_map(modeling_tools) + modeling_tools.update(params) + modeling_tools.prepare() # Bare CCL setup a_1 = 1.0 @@ -305,7 +314,7 @@ def test_pt_mixed_systematics(sacc_data): pk_mi = pt_calculator.get_biased_pk2d(tracer1=ptt_m, tracer2=ptt_i) # Set the parameters for our systematics - systematics_params = ParamsMap( + params.update( { "ia_a_1": a_1, "ia_a_2": a_2, @@ -316,7 +325,7 @@ def test_pt_mixed_systematics(sacc_data): ) # Apply the systematics parameters - likelihood.update(systematics_params) + likelihood.update(params) # Make things faster by only using a couple of ells for s in likelihood.statistics: diff --git a/tutorial/inferred_zdist.qmd b/tutorial/inferred_zdist.qmd index 42eb7380b..dabb27ca9 100644 --- a/tutorial/inferred_zdist.qmd +++ b/tutorial/inferred_zdist.qmd @@ -121,16 +121,16 @@ Markdown(f"```yaml\n{default_values_yaml}\n```") Now, we must configure the cosmology and prepare the two-point statistics for analysis: ```{python} -import pyccl - from firecrown.modeling_tools import ModelingTools +from firecrown.ccl_factory import CCLFactory +from firecrown.updatable import get_default_params_map +from firecrown.parameters import ParamsMap -ccl_cosmo = pyccl.CosmologyVanillaLCDM() -ccl_cosmo.compute_nonlin_power() +tools = ModelingTools(ccl_factory=CCLFactory(require_nonlinear_pk=True)) +params = get_default_params_map(tools, two_point_lens0_lens0) -tools = ModelingTools() tools.update(params) -tools.prepare(ccl_cosmo) +tools.prepare() two_point_lens0_lens0.update(params) theory_vector = two_point_lens0_lens0.compute_theory_vector(tools) diff --git a/tutorial/two_point_factories.qmd b/tutorial/two_point_factories.qmd index bfe8eb88c..ef3b6c90a 100644 --- a/tutorial/two_point_factories.qmd +++ b/tutorial/two_point_factories.qmd @@ -190,16 +190,17 @@ Markdown(f"```yaml\n{default_values_yaml}\n```") Finally, we can prepare both likelihoods and compare the loglike values. ```{python} -import pyccl - from firecrown.modeling_tools import ModelingTools +from firecrown.ccl_factory import CCLFactory +from firecrown.updatable import get_default_params_map +from firecrown.parameters import ParamsMap -ccl_cosmo = pyccl.CosmologyVanillaLCDM() -ccl_cosmo.compute_nonlin_power() +tools = ModelingTools(ccl_factory=CCLFactory(require_nonlinear_pk=True)) +params = get_default_params_map(tools, likelihood) tools = ModelingTools() tools.update(params) -tools.prepare(ccl_cosmo) +tools.prepare() likelihood.update(params) likelihood_ready.update(params) diff --git a/tutorial/two_point_generators.qmd b/tutorial/two_point_generators.qmd index dd38577d4..1fa32ba48 100644 --- a/tutorial/two_point_generators.qmd +++ b/tutorial/two_point_generators.qmd @@ -148,14 +148,18 @@ all_two_point_functions = tp.TwoPoint.from_metadata( ## Setting Up the Two-Point Statistics Setting up the two-point statistics requires a cosmology and a set of parameters. -The parameters necessary in our analysis depend on the two-point statistic measurements and the systematics we are considering. -In the code below, we extract the required parameters for the two-point statistics and create a `ParamsMap` object with the parameters' default values: +The parameters necessary in our analysis depend on the cosmology, two-point statistic measurements and the systematics we are considering. +The `ModelingTools` class is used to manage the cosmology and general computation objects. +In the code below, we extract the required parameters for the cosmology and the two-point statistics and create a `ParamsMap` object with the parameters' default values: ```{python} +from firecrown.modeling_tools import ModelingTools +from firecrown.ccl_factory import CCLFactory +from firecrown.updatable import get_default_params from firecrown.parameters import ParamsMap -req_params = all_two_point_functions.required_parameters() -default_values = req_params.get_default_values() +tools = ModelingTools(ccl_factory=CCLFactory(require_nonlinear_pk=True)) +default_values = get_default_params(tools, all_two_point_functions) params = ParamsMap(default_values) ``` @@ -167,7 +171,7 @@ These default values are: # | code-fold: true import yaml -default_values_yaml = yaml.dump(default_values, default_flow_style=False) +default_values_yaml = yaml.dump(dict(params), default_flow_style=False) Markdown(f"```yaml\n{default_values_yaml}\n```") ``` @@ -175,16 +179,8 @@ Markdown(f"```yaml\n{default_values_yaml}\n```") Lastly, we must configure the cosmology and prepare the two-point statistics for analysis: ```{python} -import pyccl - -from firecrown.modeling_tools import ModelingTools - -ccl_cosmo = pyccl.CosmologyVanillaLCDM() -ccl_cosmo.compute_nonlin_power() - -tools = ModelingTools() tools.update(params) -tools.prepare(ccl_cosmo) +tools.prepare() all_two_point_functions.update(params) ```