From 659836598711f31c499c9687d3426dd3f620acdd Mon Sep 17 00:00:00 2001
From: Yuanyuan Zhang <ynzhang@umich.edu>
Date: Wed, 17 Jul 2024 16:23:38 +0200
Subject: [PATCH] adding a configuration of the pipeline to run for the SDSS
 data vector (N+M)

---
 ...DSS_counts_mean_mass_redshift_richness.ini | 58 +++++++++++++++++++
 .../cluster_SDSS_redshift_richness.py         | 58 +++++++++++++++++++
 2 files changed, 116 insertions(+)
 create mode 100644 examples/cluster_number_counts/cluster_SDSS_counts_mean_mass_redshift_richness.ini
 create mode 100644 examples/cluster_number_counts/cluster_SDSS_redshift_richness.py

diff --git a/examples/cluster_number_counts/cluster_SDSS_counts_mean_mass_redshift_richness.ini b/examples/cluster_number_counts/cluster_SDSS_counts_mean_mass_redshift_richness.ini
new file mode 100644
index 00000000..edfcd655
--- /dev/null
+++ b/examples/cluster_number_counts/cluster_SDSS_counts_mean_mass_redshift_richness.ini
@@ -0,0 +1,58 @@
+[runtime]
+sampler = test
+root = ${PWD}
+
+[default]
+fatal_errors = T
+
+[output]
+filename = output_rp/number_counts_samples.txt
+format = text
+verbosity = 0
+
+[pipeline]
+modules = consistency camb firecrown_likelihood
+values = ${FIRECROWN_DIR}/examples/cluster_number_counts/cluster_richness_values.ini
+likelihoods = firecrown
+quiet = T
+debug = T
+timing = T
+
+[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 = 1.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_DIR}/examples/cluster_number_counts/cluster_SDSS_redshift_richness.py
+sampling_parameters_sections = firecrown_number_counts
+use_cluster_counts = True
+use_mean_log_mass = True
+
+[test]
+fatal_errors = T
+save_dir = output_counts_mean_mass
+
+[metropolis]
+samples = 1000
+nsteps = 1
+
+[emcee]
+walkers = 20
+samples = 4000
+nsteps = 10
diff --git a/examples/cluster_number_counts/cluster_SDSS_redshift_richness.py b/examples/cluster_number_counts/cluster_SDSS_redshift_richness.py
new file mode 100644
index 00000000..67e9607b
--- /dev/null
+++ b/examples/cluster_number_counts/cluster_SDSS_redshift_richness.py
@@ -0,0 +1,58 @@
+"""Likelihood factory function for cluster number counts."""
+
+import os
+
+import pyccl as ccl
+import sacc
+
+from firecrown.likelihood.gaussian import ConstGaussian
+from firecrown.likelihood.binned_cluster_number_counts import (
+    BinnedClusterNumberCounts,
+)
+from firecrown.likelihood.likelihood import Likelihood, NamedParameters
+from firecrown.modeling_tools import ModelingTools
+from firecrown.models.cluster.abundance import ClusterAbundance
+from firecrown.models.cluster.properties import ClusterProperty
+from firecrown.models.cluster.recipes.murata_binned_spec_z import (
+    MurataBinnedSpecZRecipe,
+)
+
+
+def get_cluster_abundance() -> ClusterAbundance:
+    """Creates and returns a ClusterAbundance object."""
+    hmf = ccl.halos.MassFuncBocquet16()
+    min_mass, max_mass = 13.0, 16.0
+    min_z, max_z = 0.2, 0.8
+    cluster_abundance = ClusterAbundance(min_mass, max_mass, min_z, max_z, hmf)
+
+    return cluster_abundance
+
+
+def build_likelihood(
+    build_parameters: NamedParameters,
+) -> tuple[Likelihood, ModelingTools]:
+    """Builds the likelihood for Firecrown."""
+    # Pull params for the likelihood from build_parameters
+    average_on = ClusterProperty.NONE
+    if build_parameters.get_bool("use_cluster_counts", True):
+        average_on |= ClusterProperty.COUNTS
+    if build_parameters.get_bool("use_mean_log_mass", False):
+        average_on |= ClusterProperty.MASS
+
+    survey_name = "SDSSCluster_redshift_richness"
+    likelihood = ConstGaussian(
+        [BinnedClusterNumberCounts(average_on, survey_name, MurataBinnedSpecZRecipe())]
+    )
+
+    # Read in sacc data
+    sacc_file_nm = "SDSSTEST_redshift_richness_sacc_data.fits"
+    sacc_path = os.path.expanduser(
+        os.path.expandvars("${FIRECROWN_DIR}/examples/cluster_number_counts/")
+    )
+    sacc_data = sacc.Sacc.load_fits(os.path.join(sacc_path, sacc_file_nm))
+    likelihood.read(sacc_data)
+
+    cluster_abundance = get_cluster_abundance()
+    modeling_tools = ModelingTools(cluster_abundance=cluster_abundance)
+
+    return likelihood, modeling_tools