Data Formats

hod_mod consumes three types of input data:

  • HDF5 — primary format for galaxy surveys (BGS/LS10, mocks) produced by the companion sum_stat package; stores the full covariance matrix and cosmological metadata.

  • CSV + JSON — paper benchmark datasets bundled in data/{paper_name}/; plain CSV files with a metadata.json sidecar describing cosmology and column meanings.

  • Benchmark-observables JSON tree — the multi-wavelength compilation of Sensitivity benchmark: the existing observables the model must reproduce in $HOD_MOD_DATA_DIR/benchmark_observables/; one self-describing JSON per (reference, observable, sample) with provenance and uncertainties (see the section below).

All spatial quantities in sum_stat are stored in Mpc (h-free). SumStatReader converts to Mpc/h automatically using the H0 attribute embedded in each file, so all arrays returned by the reader are already in the h-unit system required by hod_mod.

—

sum_stat HDF5 Schema

Single-statistic files (e.g. a w_p-only measurement) use one of the following top-level group structures. Joint files nest these same groups together.

w_p (projected correlation function)

wp/
  sep_centres       (N_rp,)  float64   — projected separation bin centres [Mpc]
  xi                (N_rp,)  float64   — w_p values [Mpc]   (name "xi" is historical)
  cov               (N_rp, N_rp) float64 — covariance matrix [Mpc²]
  bin_edges         (N_rp+1,) float64  — bin edges [Mpc]
  cosmology/
    H0              scalar  — Hubble constant H₀ [km/s/Mpc]
    Om0             scalar  — Ω_m(z=0)
    Ob0             scalar  — Ω_b(z=0)
    Ok0             scalar  — Ω_k(z=0)
  attrs:
    pi_max_Mpc      float   — line-of-sight integration limit [Mpc]
    estimator       str     — "landy-szalay" | "hamilton"
    survey          str
    n_gal           int     — galaxy count in sample

SMF (stellar mass function)

smf/
  log10mstar_centres (N_m,) float64  — log₁₀(M*/M⊙) bin centres
  phi               (N_m,)  float64  — Φ(M*) [Mpc⁻³ dex⁻¹]
  phi_err           (N_m,)  float64  — 1σ uncertainty [Mpc⁻³ dex⁻¹]
  cov               (N_m, N_m) float64
  bin_edges         (N_m+1,) float64
  cosmology/        — same sub-group as above

ESD (excess surface density / weak-lensing ΔΣ)

esd/
  rp_centres        (N_R,)  float64  — projected radius [Mpc]
  delta_sigma       (N_R,)  float64  — ΔΣ(R) [M⊙/pc²]
  cov               (N_R, N_R) float64
  cosmology/

Joint files

A joint file combines the three statistics above and adds a pre-computed joint covariance block:

smf/{sample_key}/           — SMF group (structure as above)
twopcf/{sample_key}/        — w_p group, named "wp_{sample_key}" is also common
esd/{sample_key}/           — ESD group
joint_covariance/
  data_vector       (N_tot,) float64   — [phi | wp | delta_sigma]
  cov               (N_tot, N_tot)     — full joint covariance
  err_jackknife     (N_tot,)            — √diag(cov)
  mstar_centres     (N_smf,)            — same as smf/log10mstar_centres
  rp_centres        (N_wp,)             — same as twopcf/sep_centres
  attrs:
    n_bins_smf      int
    n_bins_wp       int
    n_bins_ds       int

File naming convention:

{SURVEY}_VLIM_ANY_Mstar{MSTAR_LO}-{MSTAR_HI}_z{Z_MIN}-{Z_MAX}-{STAT}.h5

Examples:
  LS10_VLIM_ANY_Mstar10.5-12.0_z0.05-0.18-wp-pimax100-sys-comb.h5
  MOCK_VLIM_ANY_Mstar11.39_z0.05-0.35-wp-pimax100.h5

—

Unit Conversion

SumStatReader.from_hdf5() reads h = H0/100 from the embedded cosmology group and applies the following conversions before returning arrays:

Quantity

sum_stat unit

hod_mod unit

Conversion

\(r_p\) (separation)

Mpc

Mpc/h

\(r_p^{h} = r_p \times h\)

\(w_p\) (correlation)

Mpc

Mpc/h

\(w_p^{h} = w_p \times h\)

Cov(\(w_p\))

Mpc²

(Mpc/h)²

\(C^{h} = C \times h^2\)

\(\Phi(M_*)\) (SMF)

Mpc⁻³ dex⁻¹

(Mpc/h)⁻³ dex⁻¹

\(\Phi^{h} = \Phi / h^3\)

Cov(\(\Phi\))

Mpc⁻⁶

(Mpc/h)⁻⁶

\(C_\Phi^{h} = C_\Phi / h^6\)

\(\Delta\Sigma(R)\) (ESD)

\(M_\odot/\mathrm{pc}^2\)

\(M_\odot/\mathrm{pc}^2\)

(invariant — pc absorbs the h)

\(\log_{10}(M_*/M_\odot)\)

dimensionless

dimensionless

(no change needed)

Warning

The xi dataset in the wp/ HDF5 group stores the projected correlation function \(w_p(r_p)\), not the 3D correlation function \(\xi(r)\). This naming is historical (TreeCorr uses xi as the generic correlation variable).

—

Paper Benchmark Data (CSV + JSON)

All published benchmark datasets are stored under data/{paper_name}/ as plain CSV files with an accompanying metadata.json.

Directory layout

data/
  guo2018_sdss/
    metadata.json
    wp_mstar10_lowz.csv
    ...
  leauthaud2012_cosmos/
    metadata.json
    ds_photo_z2_thresh106.csv
    wp_photo_z2_thresh106.csv
    ...

wp CSV (projected correlation function)

rp_hMpc           — projected separation r_p  [h⁻¹ Mpc]
wp_hMpc           — w_p(r_p)                  [h⁻¹ Mpc]
wp_err_hMpc       — 1σ uncertainty             [h⁻¹ Mpc]

Lines beginning with # are comments (header / provenance notes) and are ignored by the reader.

ΔΣ CSV (excess surface density)

R_hMpc            — projected radius R         [h⁻¹ Mpc]
ds_Msun_h_pc2     — ΔΣ(R)                      [M⊙ h pc⁻²]
ds_err_Msun_h_pc2 — 1σ uncertainty             [M⊙ h pc⁻²]

metadata.json

One JSON sidecar per dataset. Fields:

paper             — citation string (e.g. "Guo et al. 2018")
arxiv             — arXiv identifier
survey            — survey name (e.g. "SDSS BOSS LOWZ")
sample            — description of the galaxy sample
z_eff             — effective redshift
pi_max_hMpc       — line-of-sight integration limit [h⁻¹ Mpc]  (wp only)
cosmology         — dict: Omega_m, h, sigma8, n_s, Omega_b
observable        — "wp" | "wp+ds" | …
columns_wp        — dict mapping column names to descriptions
columns_ds        — dict mapping column names to descriptions
status            — "ready" | "NEEDS_DATA" |
                    "NOT_APPLICABLE_FOR_PROJECTED_BENCHMARKS"
published_params  — dict of best-fit HOD parameters from the paper (optional)
published_param_errors — uncertainties on published_params (optional)
notes             — free-text remarks

—

Benchmark-Observables JSON Tree

The multi-wavelength benchmark compilation of Sensitivity benchmark: the existing observables the model must reproduce is materialised as a third format: a self-describing JSON tree in the data repository, $HOD_MOD_DATA_DIR/benchmark_observables/, with one file per (bibliographic reference, observable, sample):

benchmark_observables/
  README.md            — schema + operator workflow
  index.json           — every file with provenance + extraction flag
  <wavelength>/<tracer>/<RefKey>__<observable>[__<sample>].json

Each file carries the reference (citation + arXiv/DOI links), the sample definition and measurement cosmology, per-column units, a provenance block (observed | observed_derived_fit | simulated | placeholder, with a needs_operator_extraction flag naming the published table still to digitise), and the data arrays with uncertainties. The observed entries are ingested from the data/{paper_name}/ CSVs above (their metadata.json is the authoritative source for arXiv/DOI) and from the in-package X-ray band measurements; the simulated entries are forward-model fiducials with the forecast noise, standing in until the operator extracts the published table. Full layout, schema table and workflow: The benchmark data tree. Regenerate with:

python -m hod_mod.scripts.data.make_benchmark_observables

—

Reading Data in Python

from hod_mod.data_io.sum_stat_reader import SumStatReader

# ── HDF5 single w_p file ──────────────────────────────────────────────────
reader = SumStatReader.from_hdf5(
    "LS10_VLIM_ANY_Mstar10.5-12.0_z0.05-0.18-wp-pimax100-sys-comb.h5"
)
d = reader.wp()
# d["rp"]    shape (N_rp,)      Mpc/h
# d["wp"]    shape (N_rp,)      Mpc/h
# d["cov"]   shape (N_rp, N_rp) (Mpc/h)²
# d["pi_max"]                   float, Mpc/h

# ── HDF5 joint file ───────────────────────────────────────────────────────
joint_reader = SumStatReader.from_hdf5("joint_stats.h5")
j = joint_reader.joint()
# j["data_vector"]    shape (N_smf + N_wp + N_ds,)
# j["cov"]            shape (N_tot, N_tot)
# j["n_bins_smf"], j["n_bins_wp"], j["n_bins_ds"]

# ── Paper benchmark CSV + metadata.json ───────────────────────────────────
import json
import pandas as pd

meta = json.load(open("data/guo2018_sdss/metadata.json"))
df   = pd.read_csv("data/guo2018_sdss/wp_mstar10_lowz.csv", comment="#")
# df.columns: rp_hMpc, wp_hMpc, wp_err_hMpc

# ── Benchmark-observables JSON tree ───────────────────────────────────────
import os

root = os.path.join(os.environ["HOD_MOD_DATA_DIR"], "benchmark_observables")
d = json.load(open(os.path.join(
    root, "optical/galaxies/ZuMandelbaum2015__wp__10p2_10p6.json")))
# d["reference"]["citation"], d["reference"]["arxiv"]
# d["provenance"]["type"]  ("observed" here)
# d["data"]["rp_hMpc"], d["data"]["wp_hMpc"], d["data"]["wp_err_hMpc"]

idx = json.load(open(os.path.join(root, "index.json")))
todo = [k for k, v in idx.items() if v["needs_operator_extraction"]]
h = meta["cosmology"]["h"]

Read summary statistics produced by the sum_stat package.

The sum_stat package stores measurements in two formats:

  • HDF5 — LS10/BGS two-point functions and stellar mass functions. See Data Formats for the full schema.

  • FITS binary tables — GAMA and COSMOS stellar mass / luminosity functions.

All angular and projected distances are stored in physical Mpc by sum_stat. This module converts everything to h-units (Mpc/h) on the way out, using the Hubble constant stored in the file’s cosmology/H0 sub-group.

Conversion rules

  • Distances: r_h    = r_Mpc * h

  • Volumes: V_h    = V_Mpc * h^3

  • Number densities: phi_h  = phi_Mpc / h^3

  • Covariances scale as the square of the primary quantity.

Examples

Read a single-stat projected correlation function file:

reader = SumStatReader.from_hdf5(
    "/path/to/sum_stat/data/twopcf/"
    "LS10_VLIM_ANY_Mstar10.0-12.0_z0.05-0.18-wp-pimax100-sys-comb.h5"
)
data = reader.wp()
rp   = data["rp"]     # (N,) array, Mpc/h
wp   = data["wp"]     # (N,) array, Mpc/h
cov  = data["cov"]    # (N, N) array, (Mpc/h)^2

Read a joint SMF+wp+ESD file:

reader = SumStatReader.from_hdf5(
    "/path/to/sum_stat/data/BGS_Mstar10.00/"
    "joint_smf_wprp_deltasigma-sys-comb.h5"
)
jt = reader.joint()
# jt["data_vector"], jt["cov"] — full joint data vector and covariance
# jt["n_bins_smf"], jt["n_bins_wp"], jt["n_bins_ds"]

Read a GAMA stellar mass function:

reader = SumStatReader.from_fits(
    "/path/to/sum_stat/data/GAMA/gama_smf_z0.060_0.100.fits"
)
data = reader.smf()

References

class hod_mod.data_io.sum_stat_reader.SumStatReader(path: str, fmt: str, _cache: dict)[source]

Bases: object

Unified reader for sum_stat HDF5 and FITS measurement files.

Do not instantiate directly; use the class methods:

  • from_hdf5() — for HDF5 files from the twopcf, lf_smf, mocks, or BGS_Mstar* directories.

  • from_fits() — for FITS files from the GAMA or COSMOS directories.

attrs() → dict[source]

File-level attributes (creation date, version).

esd() → dict[source]

Return the excess surface mass density ΔΣ(R).

Returns:

  • dict with keys

  • * rp — projected separation bin centres, Mpc/h

  • * delta_sigma — ΔΣ(R) in M_sun pc⁻²

  • * cov — covariance matrix, (M_sun pc⁻²)²

  • * attrs — raw HDF5 group attributes

classmethod from_fits(path: str) → SumStatReader[source]

Open a GAMA or COSMOS FITS SMF/LF file.

The FITS file must have an HDU-1 binary table with at least the columns log10mstar, phi, and phi_err. Units are assumed to be Mpc (not Mpc/h); no h-conversion is applied because GAMA and COSMOS files do not embed cosmological parameters. Pass the \(h\) value explicitly when calling smf() if needed.

Parameters:

path (str) – Absolute or relative path to the FITS file.

Returns:

SumStatReader

classmethod from_hdf5(path: str) → SumStatReader[source]

Open a sum_stat HDF5 file.

Parameters:

path (str) – Absolute or relative path to the HDF5 file.

Returns:

SumStatReader

h() → float | None[source]

Hubble constant h = H0/100 embedded in the file (HDF5 only).

joint() → dict[source]

Return the full joint data vector and covariance matrix.

The joint data vector has layout [φ_SMF (n_smf), w_p (n_wp), ΔΣ (n_ds)].

Returns:

  • dict with keys

  • * data_vector — (n_total,) concatenated data vector

  • * cov — (n_total, n_total) joint covariance

  • * err_jackknife — sqrt(diag(cov))

  • * mstar_centres — log10(M*/M_sun) bin centres for SMF section

  • * rp_centres — r_p bin centres [Mpc/h] for wp and ΔΣ sections

  • * n_bins_smf, n_bins_wp, n_bins_ds — section lengths

  • * attrs — raw HDF5 group attributes

joint_bgs(probes: tuple = ('wp', 'esd_hsc', 'esd_des')) → dict[source]

Extract a multi-probe sub-data-vector from the new BGS joint format.

Selects the requested probes from the full 286-element joint data vector and returns the corresponding sub-block of the joint covariance matrix, with h-unit conversion applied:

  • wp — w_p multiplied by h (Mpc → Mpc/h); covariance × h²

  • esd_* — ΔΣ left as M_sun pc⁻² (invariant); covariance unchanged

  • Cross-terms — scaled by h¹ (WP axis) × 1 (ESD axis) = h

Parameters:

probes (tuple of str) – Any subset of ('smf', 'wp', 'esd_hsc', 'esd_des', 'esd_kids', 'wtheta', 'knn'). Default: ('wp', 'esd_hsc', 'esd_des').

Returns:

  • dict with keys

  • * data_vector — (n_sel,) sub-vector in h-units

  • * cov — (n_sel, n_sel) sub-block covariance in h-units

  • * rp_wp — (30,) r_p bin centres for WP [Mpc/h] (if wp requested)

  • * rp_esd — (30,) r_p bin centres for ESD [Mpc/h] (if any esd requested)

  • * slices_out — {probe: slice} mapping probe name to its position – in the returned data_vector

  • * h — embedded Hubble constant h = H0/100

  • * attrs — raw HDF5 group attributes

list_groups() → list[source]

List available statistic types in this file.

number_density() → dict[source]

Return the galaxy number density n of the sample.

Returns:

  • dict with keys

  • * n — number density in h³ Mpc⁻³

  • * n_err — uncertainty in h³ Mpc⁻³

  • * cov — (1, 1) variance in (h³ Mpc⁻³)²

  • * estimator — 'sum(w_i / Vmax_i)'

  • * attrs — raw HDF5 group attributes

smf(h: float | None = None) → dict[source]

Return the stellar mass function Φ(M*).

For FITS files (GAMA, COSMOS) the covariance matrix is not stored; cov will be a diagonal matrix built from phi_err.

Parameters:

h (float, optional) – Hubble constant h = H0/100 for Mpc → Mpc/h conversion of number densities. Only used for FITS files; HDF5 files carry h internally.

Returns:

  • dict with keys

  • * log10mstar — log10(M*/M_sun) bin centres

  • * phi — Φ(M*) in h³ Mpc⁻³ dex⁻¹

  • * phi_err — uncertainty

  • * cov — covariance matrix, (h³ Mpc⁻³ dex⁻¹)²

wp() → dict[source]

Return the projected correlation function w_p(r_p).

Returns:

  • dict with keys

  • * rp — projected separation bin centres, Mpc/h

  • * wp — projected correlation function, Mpc/h

  • * cov — covariance matrix, (Mpc/h)²

  • * pi_max — line-of-sight integration limit, Mpc/h

  • * estimator — 'landy-szalay' or similar

  • * attrs — raw HDF5 group attributes