From 31aa51321204602164e042d310629e8371919d46 Mon Sep 17 00:00:00 2001 From: "niko.sarcevic" Date: Tue, 12 May 2026 00:16:22 -0400 Subject: [PATCH 1/9] Add conditional luminosity function utilities --- src/lfkit/api/lumfunc.py | 197 +++++++- .../photometry/conditional_lf_integrals.py | 155 +++++++ src/lfkit/photometry/conditional_lf_models.py | 427 ++++++++++++++++++ src/lfkit/utils/types.py | 4 + 4 files changed, 782 insertions(+), 1 deletion(-) create mode 100644 src/lfkit/photometry/conditional_lf_integrals.py create mode 100644 src/lfkit/photometry/conditional_lf_models.py diff --git a/src/lfkit/api/lumfunc.py b/src/lfkit/api/lumfunc.py index 43712ad..6e1b6c1 100644 --- a/src/lfkit/api/lumfunc.py +++ b/src/lfkit/api/lumfunc.py @@ -29,6 +29,14 @@ schechter_evolving_from_m, schechter_double_from_m, ) +from lfkit.photometry.conditional_lf_models import ( + central_lognormal_conditional_lf, + central_satellite_conditional_lf, + conditional_schechter, + conditional_schechter_double, + conditional_schechter_evolving, + satellite_modified_schechter_conditional_lf, +) from lfkit.photometry.magnitudes import ( absolute_magnitude, absolute_magnitude_from_luminosity_distance, @@ -67,6 +75,7 @@ FloatInput, ParameterModel, ParameterValue, + ConditionalParameter, ) if TYPE_CHECKING: @@ -417,7 +426,7 @@ def phi( Args: absolute_mag: Absolute magnitude values where the LF is evaluated. - z: Redshift values. Required for evolving Schechter models. + z: Redshift values. Required for evolving and conditional models. Returns: Luminosity-function values evaluated at the input magnitudes. @@ -444,6 +453,66 @@ def phi( **self.parameters_dict, ) + if self.model == "conditional_schechter": + if z is None: + raise ValueError("z is required for a conditional luminosity function.") + + return conditional_schechter( + np.asarray(absolute_mag, dtype=float), + np.asarray(z, dtype=float), + **self.parameters_dict, + ) + + if self.model == "conditional_evolving_schechter": + if z is None: + raise ValueError("z is required for a conditional luminosity function.") + + return conditional_schechter_evolving( + np.asarray(absolute_mag, dtype=float), + np.asarray(z, dtype=float), + **self.parameters_dict, + ) + + if self.model == "conditional_double_schechter": + if z is None: + raise ValueError("z is required for a conditional luminosity function.") + + return conditional_schechter_double( + np.asarray(absolute_mag, dtype=float), + np.asarray(z, dtype=float), + **self.parameters_dict, + ) + + if self.model == "central_lognormal_conditional": + if z is None: + raise ValueError("z is required for a conditional luminosity function.") + + return central_lognormal_conditional_lf( + np.asarray(absolute_mag, dtype=float), + np.asarray(z, dtype=float), + **self.parameters_dict, + ) + + if self.model == "satellite_modified_schechter_conditional": + if z is None: + raise ValueError("z is required for a conditional luminosity function.") + + return satellite_modified_schechter_conditional_lf( + np.asarray(absolute_mag, dtype=float), + np.asarray(z, dtype=float), + **self.parameters_dict, + ) + + if self.model == "central_satellite_conditional": + if z is None: + raise ValueError("z is required for a conditional luminosity function.") + + return central_satellite_conditional_lf( + np.asarray(absolute_mag, dtype=float), + np.asarray(z, dtype=float), + **self.parameters_dict, + ) + raise ValueError(f"Unsupported luminosity-function model '{self.model}'.") def phi_from_m( @@ -561,6 +630,132 @@ def absolute_magnitude_limit( e_correction=e_corr, ) + @classmethod + def conditional_schechter( + cls, + *, + phi_star: ConditionalParameter, + m_star: ConditionalParameter, + alpha: ConditionalParameter, + ) -> "LuminosityFunction": + """Create a conditional Schechter luminosity function.""" + return cls( + model="conditional_schechter", + parameters={ + "phi_star": phi_star, + "m_star": m_star, + "alpha": alpha, + }, + ) + + @classmethod + def conditional_evolving_schechter( + cls, + *, + phi_model: str = "linear_p", + phi_kwargs: Mapping[str, ParameterValue] | None = None, + m_star_model: str = "linear_q", + m_star_kwargs: Mapping[str, ParameterValue] | None = None, + alpha_model: str = "constant", + alpha_kwargs: Mapping[str, ParameterValue] | None = None, + ) -> "LuminosityFunction": + """Create a conditional Schechter LF using LF parameter models.""" + return cls( + model="conditional_evolving_schechter", + parameters={ + "phi_model": phi_model, + "phi_kwargs": {} if phi_kwargs is None else dict(phi_kwargs), + "m_star_model": m_star_model, + "m_star_kwargs": {} if m_star_kwargs is None else dict(m_star_kwargs), + "alpha_model": alpha_model, + "alpha_kwargs": {} if alpha_kwargs is None else dict(alpha_kwargs), + }, + ) + + @classmethod + def conditional_double_schechter( + cls, + *, + phi_star: ConditionalParameter, + m_star: ConditionalParameter, + alpha: float, + beta: float, + m_transition: ConditionalParameter, + ) -> "LuminosityFunction": + """Create a conditional double-power-law Schechter LF.""" + return cls( + model="conditional_double_schechter", + parameters={ + "phi_star": phi_star, + "m_star": m_star, + "alpha": alpha, + "beta": beta, + "m_transition": m_transition, + }, + ) + + @classmethod + def central_lognormal_conditional( + cls, + *, + mean_absolute_mag: ConditionalParameter, + sigma_log_luminosity: ConditionalParameter, + amplitude: ConditionalParameter = 1.0, + ) -> "LuminosityFunction": + """Create a central-galaxy lognormal conditional LF.""" + return cls( + model="central_lognormal_conditional", + parameters={ + "mean_absolute_mag": mean_absolute_mag, + "sigma_log_luminosity": sigma_log_luminosity, + "amplitude": amplitude, + }, + ) + + @classmethod + def satellite_modified_schechter_conditional( + cls, + *, + phi_star: ConditionalParameter, + m_star: ConditionalParameter, + alpha: ConditionalParameter, + ) -> "LuminosityFunction": + """Create a satellite modified-Schechter conditional LF.""" + return cls( + model="satellite_modified_schechter_conditional", + parameters={ + "phi_star": phi_star, + "m_star": m_star, + "alpha": alpha, + }, + ) + + @classmethod + def central_satellite_conditional( + cls, + *, + central_mean_absolute_mag: ConditionalParameter, + central_sigma_log_luminosity: ConditionalParameter, + satellite_phi_star: ConditionalParameter, + satellite_alpha: ConditionalParameter, + central_amplitude: ConditionalParameter = 1.0, + satellite_m_star: ConditionalParameter | None = None, + satellite_luminosity_fraction: ConditionalParameter = 0.562, + ) -> "LuminosityFunction": + """Create a central-plus-satellite conditional LF.""" + return cls( + model="central_satellite_conditional", + parameters={ + "central_mean_absolute_mag": central_mean_absolute_mag, + "central_sigma_log_luminosity": central_sigma_log_luminosity, + "satellite_phi_star": satellite_phi_star, + "satellite_alpha": satellite_alpha, + "central_amplitude": central_amplitude, + "satellite_m_star": satellite_m_star, + "satellite_luminosity_fraction": satellite_luminosity_fraction, + }, + ) + @staticmethod def available_parameter_models() -> dict[str, list[str]]: """Return available LF parameter evolution models.""" diff --git a/src/lfkit/photometry/conditional_lf_integrals.py b/src/lfkit/photometry/conditional_lf_integrals.py new file mode 100644 index 0000000..29098d4 --- /dev/null +++ b/src/lfkit/photometry/conditional_lf_integrals.py @@ -0,0 +1,155 @@ +"""Conditional luminosity function integration utilities. + +This module provides small numerical helpers for conditional luminosity +functions of the form ``Phi(M | x)``, where ``M`` is absolute magnitude and +``x`` is an external conditioning variable. + +The conditioning variable is intentionally generic. It may represent halo mass, +environment, galaxy type, richness, stellar mass, or any other quantity. This +module does not implement HOD or halo-model machinery. + +The goal is to support conditional luminosity function evaluation and +integration while keeping halo-model calculations outside LFKit. +""" + +from __future__ import annotations + +from collections.abc import Callable + +import numpy as np + +from lfkit.utils.types import FloatArray, FloatInput +from lfkit.utils.validators import validate_array + +__all__ = [ + "evaluate_conditional_lf", + "integrate_conditional_lf", + "integrate_weighted_conditional_lf", +] + + +def evaluate_conditional_lf( + absolute_mag: FloatInput, + condition: FloatInput, + conditional_lf: Callable[[FloatArray, FloatArray], FloatArray], +) -> FloatArray: + """Evaluate a conditional luminosity function. + + Args: + absolute_mag: Absolute magnitude values. + condition: Values of the conditioning variable. + conditional_lf: Callable returning ``Phi(M | x)``. The callable must + accept absolute magnitude and condition arrays. + + Returns: + Conditional luminosity function evaluated at the requested absolute + magnitude and condition values. + + Raises: + ValueError: If the inputs or evaluated conditional luminosity function + contain non-finite values, or if the evaluated conditional + luminosity function contains negative values. + """ + absolute_mag_arr = validate_array(absolute_mag, name="absolute_mag") + condition_arr = validate_array(condition, name="condition") + + phi = np.asarray(conditional_lf(absolute_mag_arr, condition_arr), dtype=float) + + if not np.all(np.isfinite(phi)): + raise ValueError("conditional_lf returned NaN or infinite values.") + + if np.any(phi < 0.0): + raise ValueError( + "conditional_lf returned negative values, which are not allowed." + ) + + return np.asarray(phi, dtype=np.float64) + + +def integrate_conditional_lf( + absolute_mag: FloatInput, + condition: FloatInput, + conditional_lf: Callable[[FloatArray, FloatArray], FloatArray], + *, + axis: int = -1, +) -> FloatArray: + """Integrate a conditional luminosity function over absolute magnitude. + + This computes ``integral Phi(M | x) dM``, where ``x`` is the conditioning + variable. + + Args: + absolute_mag: Absolute magnitude grid. + condition: Values of the conditioning variable. + conditional_lf: Callable returning ``Phi(M | x)``. + axis: Axis corresponding to the absolute magnitude grid. + + Returns: + Conditional luminosity function integrated over absolute magnitude. + + Raises: + ValueError: If the inputs or evaluated conditional luminosity function + contain invalid values. + """ + absolute_mag_arr = validate_array(absolute_mag, name="absolute_mag") + + phi = evaluate_conditional_lf( + absolute_mag=absolute_mag_arr, + condition=condition, + conditional_lf=conditional_lf, + ) + + return np.asarray( + np.trapezoid(phi, x=absolute_mag_arr, axis=axis), + dtype=np.float64, + ) + + +def integrate_weighted_conditional_lf( + absolute_mag: FloatInput, + condition: FloatInput, + conditional_lf: Callable[[FloatArray, FloatArray], FloatArray], + weight: Callable[[FloatArray, FloatArray], FloatArray], + *, + axis: int = -1, +) -> FloatArray: + """Integrate a weighted conditional luminosity function. + + This computes ``integral w(M, x) Phi(M | x) dM``, where ``x`` is the + conditioning variable. + + Args: + absolute_mag: Absolute magnitude grid. + condition: Values of the conditioning variable. + conditional_lf: Callable returning ``Phi(M | x)``. + weight: Callable returning weights ``w(M, x)``. + axis: Axis corresponding to the absolute magnitude grid. + + Returns: + Weighted conditional luminosity function integrated over absolute + magnitude. + + Raises: + ValueError: If the inputs, evaluated conditional luminosity function, + or weights contain invalid values. + """ + absolute_mag_arr = validate_array(absolute_mag, name="absolute_mag") + condition_arr = validate_array(condition, name="condition") + + phi = evaluate_conditional_lf( + absolute_mag=absolute_mag_arr, + condition=condition_arr, + conditional_lf=conditional_lf, + ) + + weight_arr = np.asarray(weight(absolute_mag_arr, condition_arr), dtype=float) + + if not np.all(np.isfinite(weight_arr)): + raise ValueError("weight returned NaN or infinite values.") + + weighted_phi = np.asarray(weight_arr * phi, dtype=float) + + return np.asarray( + np.trapezoid(weighted_phi, x=absolute_mag_arr, axis=axis), + dtype=np.float64, + ) diff --git a/src/lfkit/photometry/conditional_lf_models.py b/src/lfkit/photometry/conditional_lf_models.py new file mode 100644 index 0000000..c255809 --- /dev/null +++ b/src/lfkit/photometry/conditional_lf_models.py @@ -0,0 +1,427 @@ +"""Conditional luminosity-function model utilities. + +This module provides conditional wrappers around existing LFKit luminosity +function models. + +A conditional luminosity function has the form ``Phi(M | x)``, where ``M`` is +absolute magnitude and ``x`` is an external conditioning variable. The +conditioning variable is intentionally generic. It may represent redshift, +halo mass, environment, galaxy type, richness, stellar mass, or any other +quantity. + +This module does not implement HOD or halo-model machinery. +""" + +from __future__ import annotations + +from collections.abc import Mapping +from typing import cast + +import numpy as np + +from lfkit.photometry.lf_parameter_models import evaluate_lf_parameters +from lfkit.photometry.luminosities import ( + luminosity_ratio, + magnitude_difference_from_luminosity_ratio, +) +from lfkit.photometry.luminosity_function import schechter, schechter_double +from lfkit.utils.types import ( + ConditionalParameter, + FloatArray, + FloatInput, + ParameterValue, +) +from lfkit.utils.validators import validate_array + + +__all__ = [ + "conditional_schechter", + "conditional_schechter_evolving", + "conditional_schechter_double", + "central_lognormal_conditional_lf", + "satellite_modified_schechter_conditional_lf", + "central_satellite_conditional_lf", +] + + +def conditional_schechter( + absolute_mag: FloatInput, + condition: FloatInput, + *, + phi_star: ConditionalParameter, + m_star: ConditionalParameter, + alpha: ConditionalParameter, +) -> FloatArray: + """Evaluate a conditional Schechter luminosity function. + + Args: + absolute_mag: Absolute magnitude value(s). + condition: Values of the conditioning variable. + phi_star: Schechter normalization. May be scalar, array-like, or + callable of ``condition``. + m_star: Characteristic absolute magnitude. May be scalar, array-like, + or callable of ``condition``. + alpha: Faint-end slope. May be scalar, array-like, or callable of + ``condition``. + + Returns: + Conditional Schechter luminosity-function values. + """ + condition_arr = validate_array(condition, name="condition") + + return schechter( + absolute_mag, + phi_star=_evaluate_conditional_parameter( + phi_star, + condition_arr, + name="phi_star", + ), + m_star=_evaluate_conditional_parameter( + m_star, + condition_arr, + name="m_star", + ), + alpha=_evaluate_conditional_parameter( + alpha, + condition_arr, + name="alpha", + ), + ) + + +def conditional_schechter_evolving( + absolute_mag: FloatInput, + condition: FloatInput, + *, + phi_model: str = "linear_p", + phi_kwargs: Mapping[str, ParameterValue] | None = None, + m_star_model: str = "linear_q", + m_star_kwargs: Mapping[str, ParameterValue] | None = None, + alpha_model: str = "constant", + alpha_kwargs: Mapping[str, ParameterValue] | None = None, +) -> FloatArray: + """Evaluate a conditional Schechter LF using LFKit parameter models. + + This is the conditional-LF analogue of ``schechter_evolving``. The + conditioning variable is passed to LFKit's registered parameter models. + + Args: + absolute_mag: Absolute magnitude value(s). + condition: Values of the conditioning variable. + phi_model: Evolution/condition model for ``phi_star``. + phi_kwargs: Keyword arguments passed to the selected ``phi_star`` model. + m_star_model: Evolution/condition model for ``M_star``. + m_star_kwargs: Keyword arguments passed to the selected ``M_star`` model. + alpha_model: Evolution/condition model for ``alpha``. + alpha_kwargs: Keyword arguments passed to the selected ``alpha`` model. + + Returns: + Conditional Schechter luminosity-function values. + + Raises: + ValueError: If an unsupported parameter model is requested. + """ + condition_arr = validate_array(condition, name="condition") + + phi_star, m_star, alpha = evaluate_lf_parameters( + condition_arr, + phi_model=phi_model, + phi_kwargs=phi_kwargs, + m_star_model=m_star_model, + m_star_kwargs=m_star_kwargs, + alpha_model=alpha_model, + alpha_kwargs=alpha_kwargs, + ) + + return schechter( + absolute_mag, + phi_star=phi_star, + m_star=m_star, + alpha=alpha, + ) + + +def conditional_schechter_double( + absolute_mag: FloatInput, + condition: FloatInput, + *, + phi_star: ConditionalParameter, + m_star: ConditionalParameter, + alpha: float, + beta: float, + m_transition: ConditionalParameter, +) -> FloatArray: + """Evaluate a conditional double-power-law Schechter luminosity function. + + Args: + absolute_mag: Absolute magnitude value(s). + condition: Values of the conditioning variable. + phi_star: Overall normalization. May be scalar, array-like, or + callable of ``condition``. + m_star: Characteristic absolute magnitude. May be scalar, array-like, + or callable of ``condition``. + alpha: Bright/intermediate faint-end slope parameter. + beta: Additional faint-end slope modifier. + m_transition: Transition magnitude. May be scalar, array-like, or + callable of ``condition``. + + Returns: + Conditional double-power-law Schechter luminosity-function values. + """ + condition_arr = validate_array(condition, name="condition") + + return schechter_double( + absolute_mag, + phi_star=_evaluate_conditional_parameter( + phi_star, + condition_arr, + name="phi_star", + ), + m_star=_evaluate_conditional_parameter( + m_star, + condition_arr, + name="m_star", + ), + alpha=alpha, + beta=beta, + m_transition=_evaluate_conditional_parameter( + m_transition, + condition_arr, + name="m_transition", + ), + ) + + +def central_lognormal_conditional_lf( + absolute_mag: FloatInput, + condition: FloatInput, + *, + mean_absolute_mag: ConditionalParameter, + sigma_log_luminosity: ConditionalParameter, + amplitude: ConditionalParameter = 1.0, +) -> FloatArray: + """Evaluate a central-galaxy lognormal conditional LF in magnitudes. + + Args: + absolute_mag: Absolute magnitude value(s). + condition: Values of the conditioning variable. + mean_absolute_mag: Mean central-galaxy absolute magnitude. May be + scalar, array-like, or callable of ``condition``. + sigma_log_luminosity: Scatter in ``log10(L)`` at fixed condition. May + be scalar, array-like, or callable of ``condition``. + amplitude: Non-negative amplitude of the central component. May be + scalar, array-like, or callable of ``condition``. + + Returns: + Central-galaxy conditional luminosity-function values. + """ + absolute_mag_arr = validate_array(absolute_mag, name="absolute_mag") + condition_arr = validate_array(condition, name="condition") + + mean_absolute_mag_arr = _evaluate_conditional_parameter( + mean_absolute_mag, + condition_arr, + name="mean_absolute_mag", + ) + sigma_log_luminosity_arr = _evaluate_conditional_parameter( + sigma_log_luminosity, + condition_arr, + name="sigma_log_luminosity", + ) + amplitude_arr = _evaluate_conditional_parameter( + amplitude, + condition_arr, + name="amplitude", + ) + + if np.any(sigma_log_luminosity_arr <= 0.0): + raise ValueError("sigma_log_luminosity must be positive.") + + if np.any(amplitude_arr < 0.0): + raise ValueError("amplitude must be non-negative.") + + delta_log_luminosity = -0.4 * (absolute_mag_arr - mean_absolute_mag_arr) + + phi = ( + amplitude_arr + * 0.4 + / (np.sqrt(2.0 * np.pi) * sigma_log_luminosity_arr) + * np.exp(-0.5 * (delta_log_luminosity / sigma_log_luminosity_arr) ** 2.0) + ) + + return _validate_lf_output(phi, name="central_lognormal_conditional_lf") + + +def satellite_modified_schechter_conditional_lf( + absolute_mag: FloatInput, + condition: FloatInput, + *, + phi_star: ConditionalParameter, + m_star: ConditionalParameter, + alpha: ConditionalParameter, +) -> FloatArray: + """Evaluate a satellite modified-Schechter conditional LF. + + This is the Cacciato-style satellite form with ``exp[-(L/L_star)^2]`` + instead of the standard Schechter ``exp[-L/L_star]``. + + Args: + absolute_mag: Absolute magnitude value(s). + condition: Values of the conditioning variable. + phi_star: Satellite normalization. May be scalar, array-like, or + callable of ``condition``. + m_star: Satellite characteristic absolute magnitude. May be scalar, + array-like, or callable of ``condition``. + alpha: Satellite faint-end slope. May be scalar, array-like, or + callable of ``condition``. + + Returns: + Satellite conditional luminosity-function values. + """ + absolute_mag_arr = validate_array(absolute_mag, name="absolute_mag") + condition_arr = validate_array(condition, name="condition") + + phi_star_arr = _evaluate_conditional_parameter( + phi_star, + condition_arr, + name="phi_star", + ) + m_star_arr = _evaluate_conditional_parameter( + m_star, + condition_arr, + name="m_star", + ) + alpha_arr = _evaluate_conditional_parameter( + alpha, + condition_arr, + name="alpha", + ) + + if np.any(phi_star_arr < 0.0): + raise ValueError("phi_star must be non-negative.") + + x = luminosity_ratio(absolute_mag_arr, m_star_arr) + + phi = 0.4 * np.log(10.0) * phi_star_arr * x ** (alpha_arr + 1.0) * np.exp(-(x**2.0)) + + return _validate_lf_output( + phi, + name="satellite_modified_schechter_conditional_lf", + ) + + +def central_satellite_conditional_lf( + absolute_mag: FloatInput, + condition: FloatInput, + *, + central_mean_absolute_mag: ConditionalParameter, + central_sigma_log_luminosity: ConditionalParameter, + satellite_phi_star: ConditionalParameter, + satellite_alpha: ConditionalParameter, + central_amplitude: ConditionalParameter = 1.0, + satellite_m_star: ConditionalParameter | None = None, + satellite_luminosity_fraction: ConditionalParameter = 0.562, +) -> FloatArray: + """Evaluate the sum of central and satellite conditional LF components. + + Args: + absolute_mag: Absolute magnitude value(s). + condition: Values of the conditioning variable. + central_mean_absolute_mag: Mean central-galaxy absolute magnitude. May + be scalar, array-like, or callable of ``condition``. + central_sigma_log_luminosity: Scatter in ``log10(L)`` for the central + component. May be scalar, array-like, or callable of ``condition``. + satellite_phi_star: Satellite normalization. May be scalar, array-like, + or callable of ``condition``. + satellite_alpha: Satellite faint-end slope. May be scalar, array-like, + or callable of ``condition``. + central_amplitude: Non-negative amplitude of the central component. May + be scalar, array-like, or callable of ``condition``. + satellite_m_star: Satellite characteristic absolute magnitude. If + omitted, it is derived from ``central_mean_absolute_mag`` and + ``satellite_luminosity_fraction``. + satellite_luminosity_fraction: Ratio ``L_sat_star / L_c`` used only + when ``satellite_m_star`` is omitted. + + Returns: + Total central-plus-satellite conditional luminosity-function values. + """ + condition_arr = validate_array(condition, name="condition") + + central_mean_absolute_mag_arr = _evaluate_conditional_parameter( + central_mean_absolute_mag, + condition_arr, + name="central_mean_absolute_mag", + ) + + central_phi = central_lognormal_conditional_lf( + absolute_mag=absolute_mag, + condition=condition_arr, + mean_absolute_mag=central_mean_absolute_mag_arr, + sigma_log_luminosity=central_sigma_log_luminosity, + amplitude=central_amplitude, + ) + + if satellite_m_star is None: + satellite_luminosity_fraction_arr = _evaluate_conditional_parameter( + satellite_luminosity_fraction, + condition_arr, + name="satellite_luminosity_fraction", + ) + + if np.any(satellite_luminosity_fraction_arr <= 0.0): + raise ValueError("satellite_luminosity_fraction must be positive.") + + satellite_m_star_arr = central_mean_absolute_mag_arr + ( + magnitude_difference_from_luminosity_ratio( + satellite_luminosity_fraction_arr + ) + ) + else: + satellite_m_star_arr = _evaluate_conditional_parameter( + satellite_m_star, + condition_arr, + name="satellite_m_star", + ) + + satellite_phi = satellite_modified_schechter_conditional_lf( + absolute_mag=absolute_mag, + condition=condition_arr, + phi_star=satellite_phi_star, + m_star=satellite_m_star_arr, + alpha=satellite_alpha, + ) + + return _validate_lf_output( + central_phi + satellite_phi, + name="central_satellite_conditional_lf", + ) + + +def _evaluate_conditional_parameter( + parameter: ConditionalParameter, + condition: FloatArray, + *, + name: str, +) -> FloatArray: + """Evaluate a scalar, array-like, or callable conditional parameter.""" + if callable(parameter): + values = parameter(condition) + else: + values = cast(ParameterValue, parameter) + + return validate_array(values, name=name) + + +def _validate_lf_output( + phi: FloatInput, + *, + name: str, +) -> FloatArray: + """Validate luminosity-function model output.""" + phi_arr = validate_array(phi, name=name) + + if np.any(phi_arr < 0.0): + raise ValueError(f"{name} returned negative values, which are not allowed.") + + return np.asarray(phi_arr, dtype=np.float64) diff --git a/src/lfkit/utils/types.py b/src/lfkit/utils/types.py index cb426d7..d961be1 100644 --- a/src/lfkit/utils/types.py +++ b/src/lfkit/utils/types.py @@ -16,6 +16,9 @@ Cosmology: TypeAlias = Any LuminosityFunction: TypeAlias = Callable[[FloatArray, FloatArray], FloatArray] +ConditionalParameter: TypeAlias = ( + ParameterValue | Callable[[FloatArray], ParameterValue] +) __all__ = [ "Cosmology", @@ -24,4 +27,5 @@ "LuminosityFunction", "ParameterModel", "ParameterValue", + "ConditionalParameter", ] From 041b946494b146f2759a5e4a786ffd42c2dac455 Mon Sep 17 00:00:00 2001 From: "niko.sarcevic" Date: Tue, 12 May 2026 00:16:28 -0400 Subject: [PATCH 2/9] Test conditional luminosity function utilities --- .../benchmarks/test_cacciato_clf_reference.py | 244 +++++++++ ...est_photometry_conditional_lf_integrals.py | 350 ++++++++++++ .../test_photometry_conditional_lf_models.py | 509 ++++++++++++++++++ 3 files changed, 1103 insertions(+) create mode 100644 tests/benchmarks/test_cacciato_clf_reference.py create mode 100644 tests/test_photometry_conditional_lf_integrals.py create mode 100644 tests/test_photometry_conditional_lf_models.py diff --git a/tests/benchmarks/test_cacciato_clf_reference.py b/tests/benchmarks/test_cacciato_clf_reference.py new file mode 100644 index 0000000..be34087 --- /dev/null +++ b/tests/benchmarks/test_cacciato_clf_reference.py @@ -0,0 +1,244 @@ +"""Reference benchmarks for Cacciato-style conditional luminosity functions. + +These tests compare LFKit's Cacciato-style CLF helpers against explicit +reference formulae from the Cacciato et al. conditional luminosity-function +parametrization. + +They are intentionally written as reference/benchmark tests rather than normal +unit tests. They check that LFKit preserves the expected central lognormal, +satellite modified-Schechter, and central-plus-satellite behaviour. +""" + +import numpy as np +import pytest + +from lfkit.photometry.conditional_lf_models import ( + central_lognormal_conditional_lf, + central_satellite_conditional_lf, + satellite_modified_schechter_conditional_lf, +) +from lfkit.photometry.luminosities import magnitude_difference_from_luminosity_ratio + + +pytestmark = pytest.mark.benchmark + + +def _central_magnitude_from_halo_mass( + log_halo_mass: np.ndarray, + *, + log_m1: float, + log_l0: float, + gamma1: float, + gamma2: float, + m_reference: float = 0.0, +) -> np.ndarray: + """Return central absolute magnitude from the Cacciato Lc(M) relation.""" + halo_mass = 10.0**log_halo_mass + m1 = 10.0**log_m1 + + luminosity_c = 10.0**log_l0 * (halo_mass / m1) ** gamma1 + luminosity_c /= (1.0 + halo_mass / m1) ** (gamma1 - gamma2) + + return m_reference - 2.5 * np.log10(luminosity_c) + + +def _satellite_phi_star_from_halo_mass( + log_halo_mass: np.ndarray, + *, + b0: float, + b1: float, + b2: float, +) -> np.ndarray: + """Return satellite normalization from the Cacciato quadratic relation.""" + log_m12 = log_halo_mass - 12.0 + + return 10.0 ** (b0 + b1 * log_m12 + b2 * log_m12**2.0) + + +def _reference_central_lognormal_magnitude_clf( + absolute_mag: np.ndarray, + mean_absolute_mag: np.ndarray, + sigma_log_luminosity: float, + amplitude: float = 1.0, +) -> np.ndarray: + """Return the reference central lognormal CLF in magnitude units.""" + delta_log_luminosity = -0.4 * (absolute_mag - mean_absolute_mag) + + return ( + amplitude + * 0.4 + / (np.sqrt(2.0 * np.pi) * sigma_log_luminosity) + * np.exp(-0.5 * (delta_log_luminosity / sigma_log_luminosity) ** 2.0) + ) + + +def _reference_satellite_modified_schechter_magnitude_clf( + absolute_mag: np.ndarray, + m_star: np.ndarray, + phi_star: np.ndarray, + alpha: float, +) -> np.ndarray: + """Return the reference satellite modified-Schechter CLF in magnitudes.""" + luminosity_ratio = 10.0 ** (-0.4 * (absolute_mag - m_star)) + + return ( + 0.4 + * np.log(10.0) + * phi_star + * luminosity_ratio ** (alpha + 1.0) + * np.exp(-(luminosity_ratio**2.0)) + ) + + +def test_cacciato_central_lognormal_matches_reference_formula() -> None: + """Tests central Cacciato-style CLF against the explicit formula.""" + + absolute_mag = np.array([-23.0, -22.0, -21.0, -20.0]) + log_halo_mass = np.array([11.5, 12.0, 12.5, 13.0]) + + mean_absolute_mag = _central_magnitude_from_halo_mass( + log_halo_mass, + log_m1=11.24, + log_l0=9.95, + gamma1=3.18, + gamma2=0.245, + ) + + result = central_lognormal_conditional_lf( + absolute_mag=absolute_mag, + condition=log_halo_mass, + mean_absolute_mag=mean_absolute_mag, + sigma_log_luminosity=0.157, + amplitude=1.0, + ) + + expected = _reference_central_lognormal_magnitude_clf( + absolute_mag=absolute_mag, + mean_absolute_mag=mean_absolute_mag, + sigma_log_luminosity=0.157, + amplitude=1.0, + ) + + np.testing.assert_allclose(result, expected, rtol=1.0e-12, atol=0.0) + + +def test_cacciato_satellite_modified_schechter_matches_reference_formula() -> None: + """Tests satellite Cacciato-style CLF against the explicit formula.""" + + absolute_mag = np.array([-23.0, -22.0, -21.0, -20.0]) + log_halo_mass = np.array([11.5, 12.0, 12.5, 13.0]) + + central_mag = _central_magnitude_from_halo_mass( + log_halo_mass, + log_m1=11.24, + log_l0=9.95, + gamma1=3.18, + gamma2=0.245, + ) + satellite_m_star = central_mag + magnitude_difference_from_luminosity_ratio(0.562) + + phi_star = _satellite_phi_star_from_halo_mass( + log_halo_mass, + b0=-1.17, + b1=1.53, + b2=-0.217, + ) + + result = satellite_modified_schechter_conditional_lf( + absolute_mag=absolute_mag, + condition=log_halo_mass, + phi_star=phi_star, + m_star=satellite_m_star, + alpha=-1.18, + ) + + expected = _reference_satellite_modified_schechter_magnitude_clf( + absolute_mag=absolute_mag, + m_star=satellite_m_star, + phi_star=phi_star, + alpha=-1.18, + ) + + np.testing.assert_allclose(result, expected, rtol=1.0e-12, atol=0.0) + + +def test_cacciato_central_satellite_matches_sum_of_components() -> None: + """Tests that the Cacciato-style total CLF equals central plus satellite.""" + + absolute_mag = np.array([-23.0, -22.0, -21.0, -20.0]) + log_halo_mass = np.array([11.5, 12.0, 12.5, 13.0]) + + central_mag = _central_magnitude_from_halo_mass( + log_halo_mass, + log_m1=11.24, + log_l0=9.95, + gamma1=3.18, + gamma2=0.245, + ) + phi_star = _satellite_phi_star_from_halo_mass( + log_halo_mass, + b0=-1.17, + b1=1.53, + b2=-0.217, + ) + + result = central_satellite_conditional_lf( + absolute_mag=absolute_mag, + condition=log_halo_mass, + central_mean_absolute_mag=central_mag, + central_sigma_log_luminosity=0.157, + satellite_phi_star=phi_star, + satellite_alpha=-1.18, + central_amplitude=1.0, + satellite_m_star=None, + satellite_luminosity_fraction=0.562, + ) + + satellite_m_star = central_mag + magnitude_difference_from_luminosity_ratio(0.562) + + central = _reference_central_lognormal_magnitude_clf( + absolute_mag=absolute_mag, + mean_absolute_mag=central_mag, + sigma_log_luminosity=0.157, + amplitude=1.0, + ) + satellite = _reference_satellite_modified_schechter_magnitude_clf( + absolute_mag=absolute_mag, + m_star=satellite_m_star, + phi_star=phi_star, + alpha=-1.18, + ) + + expected = central + satellite + + np.testing.assert_allclose(result, expected, rtol=1.0e-12, atol=0.0) + + +def test_cacciato_satellite_m_star_is_fainter_than_central_for_fraction_below_one() -> None: + """Tests the Cacciato convention Ls_star = 0.562 Lc in magnitude space.""" + + central_mag = np.array([-22.0, -21.0, -20.0]) + + satellite_m_star = central_mag + magnitude_difference_from_luminosity_ratio(0.562) + + assert np.all(satellite_m_star > central_mag) + np.testing.assert_allclose( + satellite_m_star - central_mag, + magnitude_difference_from_luminosity_ratio(0.562), + ) + + +def test_cacciato_satellite_phi_star_peaks_near_reference_mass_range() -> None: + """Tests that the quadratic satellite normalization remains positive.""" + + log_halo_mass = np.linspace(11.0, 14.5, 32) + + phi_star = _satellite_phi_star_from_halo_mass( + log_halo_mass, + b0=-1.17, + b1=1.53, + b2=-0.217, + ) + + assert np.all(np.isfinite(phi_star)) + assert np.all(phi_star > 0.0) \ No newline at end of file diff --git a/tests/test_photometry_conditional_lf_integrals.py b/tests/test_photometry_conditional_lf_integrals.py new file mode 100644 index 0000000..e2187b2 --- /dev/null +++ b/tests/test_photometry_conditional_lf_integrals.py @@ -0,0 +1,350 @@ +"""Unit tests for ``lfkit.photometry.conditional_lf_integrals.py``.""" + +import numpy as np +import pytest + +from lfkit.photometry.conditional_lf_integrals import ( + evaluate_conditional_lf, + integrate_conditional_lf, + integrate_weighted_conditional_lf, +) + + +def test_evaluate_conditional_lf_accepts_scalar_inputs() -> None: + """Tests that scalar inputs are evaluated as float arrays.""" + + def conditional_lf(absolute_mag, condition): + return condition * np.exp(-0.1 * absolute_mag) + + result = evaluate_conditional_lf( + absolute_mag=-20.0, + condition=2.0, + conditional_lf=conditional_lf, + ) + + expected = 2.0 * np.exp(2.0) + + assert isinstance(result, np.ndarray) + assert result.dtype == np.float64 + assert result.shape == () + assert result == pytest.approx(expected) + + +def test_evaluate_conditional_lf_accepts_array_inputs() -> None: + """Tests that array inputs are evaluated element-wise.""" + + absolute_mag = np.array([-22.0, -21.0, -20.0]) + condition = np.array([1.0, 2.0, 3.0]) + + def conditional_lf(absolute_mag, condition): + return condition * (absolute_mag + 23.0) + + result = evaluate_conditional_lf( + absolute_mag=absolute_mag, + condition=condition, + conditional_lf=conditional_lf, + ) + + expected = np.array([1.0, 4.0, 9.0]) + + np.testing.assert_allclose(result, expected) + assert result.dtype == np.float64 + + +def test_evaluate_conditional_lf_accepts_broadcastable_inputs() -> None: + """Tests that broadcastable magnitude and condition arrays are supported.""" + + absolute_mag = np.array([-22.0, -21.0, -20.0]) + condition = np.array([[1.0], [2.0]]) + + def conditional_lf(absolute_mag, condition): + return condition * (absolute_mag + 23.0) + + result = evaluate_conditional_lf( + absolute_mag=absolute_mag, + condition=condition, + conditional_lf=conditional_lf, + ) + + expected = np.array( + [ + [1.0, 2.0, 3.0], + [2.0, 4.0, 6.0], + ] + ) + + np.testing.assert_allclose(result, expected) + assert result.dtype == np.float64 + + +def test_evaluate_conditional_lf_rejects_non_finite_absolute_mag() -> None: + """Tests that non-finite absolute magnitudes are rejected.""" + + def conditional_lf(absolute_mag, condition): + return np.ones_like(absolute_mag) + + with pytest.raises(ValueError, match="absolute_mag contains NaN or infinite values."): + evaluate_conditional_lf( + absolute_mag=[-22.0, np.nan, -20.0], + condition=1.0, + conditional_lf=conditional_lf, + ) + + +def test_evaluate_conditional_lf_rejects_non_finite_condition() -> None: + """Tests that non-finite condition values are rejected.""" + + def conditional_lf(absolute_mag, condition): + return np.ones_like(absolute_mag) + + with pytest.raises(ValueError, match="condition contains NaN or infinite values."): + evaluate_conditional_lf( + absolute_mag=[-22.0, -21.0, -20.0], + condition=np.inf, + conditional_lf=conditional_lf, + ) + + +def test_evaluate_conditional_lf_rejects_non_finite_result() -> None: + """Tests that non-finite conditional luminosity values are rejected.""" + + def conditional_lf(absolute_mag, condition): + return np.array([1.0, np.nan, 3.0]) + + with pytest.raises( + ValueError, + match="conditional_lf returned NaN or infinite values.", + ): + evaluate_conditional_lf( + absolute_mag=[-22.0, -21.0, -20.0], + condition=1.0, + conditional_lf=conditional_lf, + ) + + +def test_evaluate_conditional_lf_rejects_negative_result() -> None: + """Tests that negative conditional luminosity values are rejected.""" + + def conditional_lf(absolute_mag, condition): + return np.array([1.0, -2.0, 3.0]) + + with pytest.raises( + ValueError, + match="conditional_lf returned negative values, which are not allowed.", + ): + evaluate_conditional_lf( + absolute_mag=[-22.0, -21.0, -20.0], + condition=1.0, + conditional_lf=conditional_lf, + ) + + +def test_integrate_conditional_lf_integrates_over_last_axis() -> None: + """Tests conditional luminosity-function integration over magnitude.""" + + absolute_mag = np.array([-22.0, -21.0, -20.0]) + condition = np.array([[1.0], [2.0]]) + + def conditional_lf(absolute_mag, condition): + return condition * (absolute_mag + 23.0) + + result = integrate_conditional_lf( + absolute_mag=absolute_mag, + condition=condition, + conditional_lf=conditional_lf, + ) + + expected = np.trapezoid( + np.array( + [ + [1.0, 2.0, 3.0], + [2.0, 4.0, 6.0], + ] + ), + x=absolute_mag, + axis=-1, + ) + + np.testing.assert_allclose(result, expected) + assert result.dtype == np.float64 + + +def test_integrate_conditional_lf_integrates_over_requested_axis() -> None: + """Tests conditional luminosity-function integration over a custom axis.""" + + absolute_mag = np.array([-22.0, -21.0, -20.0]) + condition = np.array([[1.0, 2.0]]) + + def conditional_lf(absolute_mag, condition): + return (absolute_mag[:, None] + 23.0) * condition + + result = integrate_conditional_lf( + absolute_mag=absolute_mag, + condition=condition, + conditional_lf=conditional_lf, + axis=0, + ) + + phi = np.array( + [ + [1.0, 2.0], + [2.0, 4.0], + [3.0, 6.0], + ] + ) + expected = np.trapezoid(phi, x=absolute_mag, axis=0) + + np.testing.assert_allclose(result, expected) + assert result.dtype == np.float64 + + +def test_integrate_conditional_lf_rejects_invalid_conditional_lf() -> None: + """Tests that integration propagates invalid conditional LF errors.""" + + def conditional_lf(absolute_mag, condition): + return np.array([1.0, -1.0, 3.0]) + + with pytest.raises( + ValueError, + match="conditional_lf returned negative values, which are not allowed.", + ): + integrate_conditional_lf( + absolute_mag=[-22.0, -21.0, -20.0], + condition=1.0, + conditional_lf=conditional_lf, + ) + + +def test_integrate_weighted_conditional_lf_integrates_weighted_values() -> None: + """Tests weighted conditional luminosity-function integration.""" + + absolute_mag = np.array([-22.0, -21.0, -20.0]) + condition = np.array([[1.0], [2.0]]) + + def conditional_lf(absolute_mag, condition): + return condition * (absolute_mag + 23.0) + + def weight(absolute_mag, condition): + return absolute_mag + 24.0 + + result = integrate_weighted_conditional_lf( + absolute_mag=absolute_mag, + condition=condition, + conditional_lf=conditional_lf, + weight=weight, + ) + + phi = np.array( + [ + [1.0, 2.0, 3.0], + [2.0, 4.0, 6.0], + ] + ) + weight_arr = np.array([2.0, 3.0, 4.0]) + expected = np.trapezoid(phi * weight_arr, x=absolute_mag, axis=-1) + + np.testing.assert_allclose(result, expected) + assert result.dtype == np.float64 + + +def test_integrate_weighted_conditional_lf_integrates_over_requested_axis() -> None: + """Tests weighted conditional LF integration over a custom axis.""" + + absolute_mag = np.array([-22.0, -21.0, -20.0]) + condition = np.array([[1.0, 2.0]]) + + def conditional_lf(absolute_mag, condition): + return (absolute_mag[:, None] + 23.0) * condition + + def weight(absolute_mag, condition): + return absolute_mag[:, None] + 24.0 + + result = integrate_weighted_conditional_lf( + absolute_mag=absolute_mag, + condition=condition, + conditional_lf=conditional_lf, + weight=weight, + axis=0, + ) + + phi = np.array( + [ + [1.0, 2.0], + [2.0, 4.0], + [3.0, 6.0], + ] + ) + weight_arr = np.array( + [ + [2.0], + [3.0], + [4.0], + ] + ) + expected = np.trapezoid(phi * weight_arr, x=absolute_mag, axis=0) + + np.testing.assert_allclose(result, expected) + assert result.dtype == np.float64 + + +def test_integrate_weighted_conditional_lf_rejects_non_finite_weight() -> None: + """Tests that non-finite weights are rejected.""" + + def conditional_lf(absolute_mag, condition): + return np.ones_like(absolute_mag) + + def weight(absolute_mag, condition): + return np.array([1.0, np.inf, 3.0]) + + with pytest.raises(ValueError, match="weight returned NaN or infinite values."): + integrate_weighted_conditional_lf( + absolute_mag=[-22.0, -21.0, -20.0], + condition=1.0, + conditional_lf=conditional_lf, + weight=weight, + ) + + +def test_integrate_weighted_conditional_lf_allows_negative_weight() -> None: + """Tests that negative weights are allowed.""" + + absolute_mag = np.array([-22.0, -21.0, -20.0]) + + def conditional_lf(absolute_mag, condition): + return np.ones_like(absolute_mag) + + def weight(absolute_mag, condition): + return np.array([1.0, -2.0, 3.0]) + + result = integrate_weighted_conditional_lf( + absolute_mag=absolute_mag, + condition=1.0, + conditional_lf=conditional_lf, + weight=weight, + ) + + expected = np.trapezoid(np.array([1.0, -2.0, 3.0]), x=absolute_mag) + + assert result == pytest.approx(expected) + assert result.dtype == np.float64 + + +def test_integrate_weighted_conditional_lf_rejects_invalid_conditional_lf() -> None: + """Tests that weighted integration propagates invalid conditional LF errors.""" + + def conditional_lf(absolute_mag, condition): + return np.array([1.0, -1.0, 3.0]) + + def weight(absolute_mag, condition): + return np.ones_like(absolute_mag) + + with pytest.raises( + ValueError, + match="conditional_lf returned negative values, which are not allowed.", + ): + integrate_weighted_conditional_lf( + absolute_mag=[-22.0, -21.0, -20.0], + condition=1.0, + conditional_lf=conditional_lf, + weight=weight, + ) diff --git a/tests/test_photometry_conditional_lf_models.py b/tests/test_photometry_conditional_lf_models.py new file mode 100644 index 0000000..f83e5eb --- /dev/null +++ b/tests/test_photometry_conditional_lf_models.py @@ -0,0 +1,509 @@ +"""Unit tests for ``lfkit.photometry.conditional_lf_models.py``.""" + +import numpy as np +import pytest + +from lfkit.photometry.conditional_lf_models import ( + central_lognormal_conditional_lf, + central_satellite_conditional_lf, + conditional_schechter, + conditional_schechter_double, + conditional_schechter_evolving, + satellite_modified_schechter_conditional_lf, +) +from lfkit.photometry.luminosities import ( + luminosity_ratio, + magnitude_difference_from_luminosity_ratio, +) +from lfkit.photometry.luminosity_function import schechter, schechter_double + + +def test_conditional_schechter_matches_schechter_for_scalar_parameters() -> None: + """Tests that the conditional Schechter wrapper matches Schechter.""" + + absolute_mag = np.array([-22.0, -21.0, -20.0]) + condition = np.array([0.1, 0.2, 0.3]) + + result = conditional_schechter( + absolute_mag=absolute_mag, + condition=condition, + phi_star=1.0e-3, + m_star=-21.0, + alpha=-1.1, + ) + + expected = schechter( + absolute_mag, + phi_star=1.0e-3, + m_star=-21.0, + alpha=-1.1, + ) + + np.testing.assert_allclose(result, expected) + assert result.dtype == np.float64 + + +def test_conditional_schechter_accepts_callable_parameters() -> None: + """Tests that conditional Schechter parameters can be callables.""" + + absolute_mag = np.array([-22.0, -21.0, -20.0]) + condition = np.array([0.0, 1.0, 2.0]) + + result = conditional_schechter( + absolute_mag=absolute_mag, + condition=condition, + phi_star=lambda x: 1.0e-3 * (1.0 + x), + m_star=lambda x: -21.0 - 0.1 * x, + alpha=lambda x: -1.0 - 0.05 * x, + ) + + expected = schechter( + absolute_mag, + phi_star=np.array([1.0e-3, 2.0e-3, 3.0e-3]), + m_star=np.array([-21.0, -21.1, -21.2]), + alpha=np.array([-1.0, -1.05, -1.1]), + ) + + np.testing.assert_allclose(result, expected) + assert result.dtype == np.float64 + + +def test_conditional_schechter_rejects_non_finite_condition() -> None: + """Tests that non-finite condition values are rejected.""" + + with pytest.raises(ValueError, match="condition contains NaN or infinite values."): + conditional_schechter( + absolute_mag=[-22.0, -21.0, -20.0], + condition=[0.0, np.nan, 1.0], + phi_star=1.0e-3, + m_star=-21.0, + alpha=-1.1, + ) + + +def test_conditional_schechter_rejects_non_finite_callable_parameter() -> None: + """Tests that non-finite callable parameter values are rejected.""" + + with pytest.raises(ValueError, match="phi_star contains NaN or infinite values."): + conditional_schechter( + absolute_mag=[-22.0, -21.0, -20.0], + condition=[0.0, 1.0, 2.0], + phi_star=lambda x: np.array([1.0e-3, np.nan, 2.0e-3]), + m_star=-21.0, + alpha=-1.1, + ) + + +def test_conditional_schechter_evolving_matches_explicit_parameter_models() -> None: + """Tests the conditional evolving Schechter wrapper with simple models.""" + + absolute_mag = np.array([-22.0, -21.0, -20.0]) + condition = np.array([0.0, 1.0, 2.0]) + + result = conditional_schechter_evolving( + absolute_mag=absolute_mag, + condition=condition, + phi_model="constant", + phi_kwargs={"phi_star": 1.0e-3}, + m_star_model="constant", + m_star_kwargs={"m_star": -21.0}, + alpha_model="constant", + alpha_kwargs={"alpha": -1.1}, + ) + + expected = schechter( + absolute_mag, + phi_star=1.0e-3, + m_star=-21.0, + alpha=-1.1, + ) + + np.testing.assert_allclose(result, expected) + assert result.dtype == np.float64 + + +def test_conditional_schechter_evolving_rejects_unknown_model() -> None: + """Tests that unknown LF parameter models are rejected.""" + + with pytest.raises(ValueError): + conditional_schechter_evolving( + absolute_mag=[-22.0, -21.0, -20.0], + condition=[0.0, 1.0, 2.0], + phi_model="not_a_model", + phi_kwargs={"value": 1.0e-3}, + m_star_model="constant", + m_star_kwargs={"value": -21.0}, + alpha_model="constant", + alpha_kwargs={"value": -1.1}, + ) + + +def test_conditional_schechter_double_matches_double_schechter() -> None: + """Tests that the conditional double-Schechter wrapper matches the model.""" + + absolute_mag = np.array([-22.0, -21.0, -20.0]) + condition = np.array([0.0, 1.0, 2.0]) + + result = conditional_schechter_double( + absolute_mag=absolute_mag, + condition=condition, + phi_star=1.0e-3, + m_star=-21.0, + alpha=-1.0, + beta=-0.5, + m_transition=-19.5, + ) + + expected = schechter_double( + absolute_mag, + phi_star=1.0e-3, + m_star=-21.0, + alpha=-1.0, + beta=-0.5, + m_transition=-19.5, + ) + + np.testing.assert_allclose(result, expected) + assert result.dtype == np.float64 + + +def test_conditional_schechter_double_accepts_callable_parameters() -> None: + """Tests callable parameters for the conditional double-Schechter model.""" + + absolute_mag = np.array([-22.0, -21.0, -20.0]) + condition = np.array([0.0, 1.0, 2.0]) + + result = conditional_schechter_double( + absolute_mag=absolute_mag, + condition=condition, + phi_star=lambda x: 1.0e-3 * (1.0 + x), + m_star=lambda x: -21.0 - 0.1 * x, + alpha=-1.0, + beta=-0.5, + m_transition=lambda x: -19.5 - 0.2 * x, + ) + + expected = schechter_double( + absolute_mag, + phi_star=np.array([1.0e-3, 2.0e-3, 3.0e-3]), + m_star=np.array([-21.0, -21.1, -21.2]), + alpha=-1.0, + beta=-0.5, + m_transition=np.array([-19.5, -19.7, -19.9]), + ) + + np.testing.assert_allclose(result, expected) + assert result.dtype == np.float64 + + +def test_central_lognormal_conditional_lf_matches_expected_formula() -> None: + """Tests the central lognormal conditional LF formula.""" + + absolute_mag = np.array([-22.0, -21.0, -20.0]) + condition = np.array([0.0, 1.0, 2.0]) + mean_absolute_mag = np.array([-21.0, -21.0, -21.0]) + sigma_log_luminosity = np.array([0.2, 0.2, 0.2]) + amplitude = np.array([1.0, 2.0, 3.0]) + + result = central_lognormal_conditional_lf( + absolute_mag=absolute_mag, + condition=condition, + mean_absolute_mag=mean_absolute_mag, + sigma_log_luminosity=sigma_log_luminosity, + amplitude=amplitude, + ) + + delta_log_luminosity = -0.4 * (absolute_mag - mean_absolute_mag) + expected = ( + amplitude + * 0.4 + / (np.sqrt(2.0 * np.pi) * sigma_log_luminosity) + * np.exp(-0.5 * (delta_log_luminosity / sigma_log_luminosity) ** 2.0) + ) + + np.testing.assert_allclose(result, expected) + assert result.dtype == np.float64 + + +def test_central_lognormal_conditional_lf_accepts_callable_parameters() -> None: + """Tests callable parameters for the central lognormal conditional LF.""" + + absolute_mag = np.array([-22.0, -21.0, -20.0]) + condition = np.array([0.0, 1.0, 2.0]) + + result = central_lognormal_conditional_lf( + absolute_mag=absolute_mag, + condition=condition, + mean_absolute_mag=lambda x: -21.0 - 0.1 * x, + sigma_log_luminosity=lambda x: 0.2 + 0.01 * x, + amplitude=lambda x: 1.0 + x, + ) + + mean_absolute_mag = np.array([-21.0, -21.1, -21.2]) + sigma_log_luminosity = np.array([0.2, 0.21, 0.22]) + amplitude = np.array([1.0, 2.0, 3.0]) + + delta_log_luminosity = -0.4 * (absolute_mag - mean_absolute_mag) + expected = ( + amplitude + * 0.4 + / (np.sqrt(2.0 * np.pi) * sigma_log_luminosity) + * np.exp(-0.5 * (delta_log_luminosity / sigma_log_luminosity) ** 2.0) + ) + + np.testing.assert_allclose(result, expected) + assert result.dtype == np.float64 + + +def test_central_lognormal_conditional_lf_rejects_zero_sigma() -> None: + """Tests that zero central scatter is rejected.""" + + with pytest.raises(ValueError, match="sigma_log_luminosity must be positive."): + central_lognormal_conditional_lf( + absolute_mag=[-22.0, -21.0, -20.0], + condition=[0.0, 1.0, 2.0], + mean_absolute_mag=-21.0, + sigma_log_luminosity=0.0, + amplitude=1.0, + ) + + +def test_central_lognormal_conditional_lf_rejects_negative_sigma() -> None: + """Tests that negative central scatter is rejected.""" + + with pytest.raises(ValueError, match="sigma_log_luminosity must be positive."): + central_lognormal_conditional_lf( + absolute_mag=[-22.0, -21.0, -20.0], + condition=[0.0, 1.0, 2.0], + mean_absolute_mag=-21.0, + sigma_log_luminosity=-0.2, + amplitude=1.0, + ) + + +def test_central_lognormal_conditional_lf_rejects_negative_amplitude() -> None: + """Tests that negative central amplitude is rejected.""" + + with pytest.raises(ValueError, match="amplitude must be non-negative."): + central_lognormal_conditional_lf( + absolute_mag=[-22.0, -21.0, -20.0], + condition=[0.0, 1.0, 2.0], + mean_absolute_mag=-21.0, + sigma_log_luminosity=0.2, + amplitude=-1.0, + ) + + +def test_satellite_modified_schechter_conditional_lf_matches_expected_formula() -> None: + """Tests the satellite modified-Schechter conditional LF formula.""" + + absolute_mag = np.array([-22.0, -21.0, -20.0]) + condition = np.array([0.0, 1.0, 2.0]) + phi_star = np.array([1.0e-3, 2.0e-3, 3.0e-3]) + m_star = np.array([-21.0, -21.1, -21.2]) + alpha = np.array([-1.0, -1.1, -1.2]) + + result = satellite_modified_schechter_conditional_lf( + absolute_mag=absolute_mag, + condition=condition, + phi_star=phi_star, + m_star=m_star, + alpha=alpha, + ) + + x = luminosity_ratio(absolute_mag, m_star) + expected = 0.4 * np.log(10.0) * phi_star * x ** (alpha + 1.0) * np.exp(-(x**2.0)) + + np.testing.assert_allclose(result, expected) + assert result.dtype == np.float64 + + +def test_satellite_modified_schechter_conditional_lf_accepts_callable_parameters() -> None: + """Tests callable parameters for the satellite modified-Schechter LF.""" + + absolute_mag = np.array([-22.0, -21.0, -20.0]) + condition = np.array([0.0, 1.0, 2.0]) + + result = satellite_modified_schechter_conditional_lf( + absolute_mag=absolute_mag, + condition=condition, + phi_star=lambda x: 1.0e-3 * (1.0 + x), + m_star=lambda x: -21.0 - 0.1 * x, + alpha=lambda x: -1.0 - 0.05 * x, + ) + + phi_star = np.array([1.0e-3, 2.0e-3, 3.0e-3]) + m_star = np.array([-21.0, -21.1, -21.2]) + alpha = np.array([-1.0, -1.05, -1.1]) + + x = luminosity_ratio(absolute_mag, m_star) + expected = 0.4 * np.log(10.0) * phi_star * x ** (alpha + 1.0) * np.exp(-(x**2.0)) + + np.testing.assert_allclose(result, expected) + assert result.dtype == np.float64 + + +def test_satellite_modified_schechter_conditional_lf_rejects_negative_phi_star() -> None: + """Tests that negative satellite normalization is rejected.""" + + with pytest.raises(ValueError, match="phi_star must be non-negative."): + satellite_modified_schechter_conditional_lf( + absolute_mag=[-22.0, -21.0, -20.0], + condition=[0.0, 1.0, 2.0], + phi_star=-1.0e-3, + m_star=-21.0, + alpha=-1.1, + ) + + +def test_central_satellite_conditional_lf_equals_sum_with_explicit_satellite_m_star() -> None: + """Tests that the total conditional LF equals central plus satellite parts.""" + + absolute_mag = np.array([-22.0, -21.0, -20.0]) + condition = np.array([0.0, 1.0, 2.0]) + + result = central_satellite_conditional_lf( + absolute_mag=absolute_mag, + condition=condition, + central_mean_absolute_mag=-21.0, + central_sigma_log_luminosity=0.2, + satellite_phi_star=1.0e-3, + satellite_alpha=-1.1, + central_amplitude=1.0, + satellite_m_star=-20.5, + ) + + central = central_lognormal_conditional_lf( + absolute_mag=absolute_mag, + condition=condition, + mean_absolute_mag=-21.0, + sigma_log_luminosity=0.2, + amplitude=1.0, + ) + satellite = satellite_modified_schechter_conditional_lf( + absolute_mag=absolute_mag, + condition=condition, + phi_star=1.0e-3, + m_star=-20.5, + alpha=-1.1, + ) + expected = central + satellite + + np.testing.assert_allclose(result, expected) + assert result.dtype == np.float64 + + +def test_central_satellite_conditional_lf_derives_satellite_m_star() -> None: + """Tests that satellite M-star is derived from the luminosity fraction.""" + + absolute_mag = np.array([-22.0, -21.0, -20.0]) + condition = np.array([0.0, 1.0, 2.0]) + central_mean_absolute_mag = np.array([-21.0, -21.1, -21.2]) + satellite_luminosity_fraction = np.array([0.5, 0.6, 0.7]) + + result = central_satellite_conditional_lf( + absolute_mag=absolute_mag, + condition=condition, + central_mean_absolute_mag=central_mean_absolute_mag, + central_sigma_log_luminosity=0.2, + satellite_phi_star=1.0e-3, + satellite_alpha=-1.1, + central_amplitude=1.0, + satellite_m_star=None, + satellite_luminosity_fraction=satellite_luminosity_fraction, + ) + + satellite_m_star = central_mean_absolute_mag + ( + magnitude_difference_from_luminosity_ratio(satellite_luminosity_fraction) + ) + + central = central_lognormal_conditional_lf( + absolute_mag=absolute_mag, + condition=condition, + mean_absolute_mag=central_mean_absolute_mag, + sigma_log_luminosity=0.2, + amplitude=1.0, + ) + satellite = satellite_modified_schechter_conditional_lf( + absolute_mag=absolute_mag, + condition=condition, + phi_star=1.0e-3, + m_star=satellite_m_star, + alpha=-1.1, + ) + expected = central + satellite + + np.testing.assert_allclose(result, expected) + assert result.dtype == np.float64 + + +def test_central_satellite_conditional_lf_rejects_zero_luminosity_fraction() -> None: + """Tests that zero satellite luminosity fraction is rejected.""" + + with pytest.raises( + ValueError, + match="satellite_luminosity_fraction must be positive.", + ): + central_satellite_conditional_lf( + absolute_mag=[-22.0, -21.0, -20.0], + condition=[0.0, 1.0, 2.0], + central_mean_absolute_mag=-21.0, + central_sigma_log_luminosity=0.2, + satellite_phi_star=1.0e-3, + satellite_alpha=-1.1, + central_amplitude=1.0, + satellite_m_star=None, + satellite_luminosity_fraction=0.0, + ) + + +def test_central_satellite_conditional_lf_rejects_negative_luminosity_fraction() -> None: + """Tests that negative satellite luminosity fraction is rejected.""" + + with pytest.raises( + ValueError, + match="satellite_luminosity_fraction must be positive.", + ): + central_satellite_conditional_lf( + absolute_mag=[-22.0, -21.0, -20.0], + condition=[0.0, 1.0, 2.0], + central_mean_absolute_mag=-21.0, + central_sigma_log_luminosity=0.2, + satellite_phi_star=1.0e-3, + satellite_alpha=-1.1, + central_amplitude=1.0, + satellite_m_star=None, + satellite_luminosity_fraction=-0.5, + ) + + +def test_central_satellite_conditional_lf_propagates_invalid_central_component() -> None: + """Tests that invalid central-component parameters are propagated.""" + + with pytest.raises(ValueError, match="central_sigma_log_luminosity must be positive.|sigma_log_luminosity must be positive."): + central_satellite_conditional_lf( + absolute_mag=[-22.0, -21.0, -20.0], + condition=[0.0, 1.0, 2.0], + central_mean_absolute_mag=-21.0, + central_sigma_log_luminosity=0.0, + satellite_phi_star=1.0e-3, + satellite_alpha=-1.1, + central_amplitude=1.0, + satellite_m_star=-20.5, + ) + + +def test_central_satellite_conditional_lf_propagates_invalid_satellite_component() -> None: + """Tests that invalid satellite-component parameters are propagated.""" + + with pytest.raises(ValueError, match="phi_star must be non-negative."): + central_satellite_conditional_lf( + absolute_mag=[-22.0, -21.0, -20.0], + condition=[0.0, 1.0, 2.0], + central_mean_absolute_mag=-21.0, + central_sigma_log_luminosity=0.2, + satellite_phi_star=-1.0e-3, + satellite_alpha=-1.1, + central_amplitude=1.0, + satellite_m_star=-20.5, + ) From 201135a2f2d496005720fb7dc2a4bfd0a7e01e14 Mon Sep 17 00:00:00 2001 From: "niko.sarcevic" Date: Tue, 12 May 2026 00:16:41 -0400 Subject: [PATCH 3/9] Document conditional luminosity function examples --- ...it.photometry.conditional_lf_integrals.rst | 7 + ...lfkit.photometry.conditional_lf_models.rst | 7 + docs/api/lfkit.photometry.rst | 2 + ...nditional_luminosity_function_examples.rst | 788 ++++++++++++++++++ docs/examples/index.rst | 16 +- 5 files changed, 819 insertions(+), 1 deletion(-) create mode 100644 docs/api/lfkit.photometry.conditional_lf_integrals.rst create mode 100644 docs/api/lfkit.photometry.conditional_lf_models.rst create mode 100644 docs/examples/conditional_luminosity_function_examples.rst diff --git a/docs/api/lfkit.photometry.conditional_lf_integrals.rst b/docs/api/lfkit.photometry.conditional_lf_integrals.rst new file mode 100644 index 0000000..1912c31 --- /dev/null +++ b/docs/api/lfkit.photometry.conditional_lf_integrals.rst @@ -0,0 +1,7 @@ +lfkit.photometry.conditional\_lf\_integrals module +================================================== + +.. automodule:: lfkit.photometry.conditional_lf_integrals + :members: + :show-inheritance: + :undoc-members: diff --git a/docs/api/lfkit.photometry.conditional_lf_models.rst b/docs/api/lfkit.photometry.conditional_lf_models.rst new file mode 100644 index 0000000..13f14a0 --- /dev/null +++ b/docs/api/lfkit.photometry.conditional_lf_models.rst @@ -0,0 +1,7 @@ +lfkit.photometry.conditional\_lf\_models module +=============================================== + +.. automodule:: lfkit.photometry.conditional_lf_models + :members: + :show-inheritance: + :undoc-members: diff --git a/docs/api/lfkit.photometry.rst b/docs/api/lfkit.photometry.rst index 254d341..c7c6fae 100644 --- a/docs/api/lfkit.photometry.rst +++ b/docs/api/lfkit.photometry.rst @@ -8,6 +8,8 @@ Submodules :maxdepth: 2 lfkit.photometry.catalog_completeness + lfkit.photometry.conditional_lf_integrals + lfkit.photometry.conditional_lf_models lfkit.photometry.lf_integrals lfkit.photometry.lf_parameter_models lfkit.photometry.lf_redshift_density diff --git a/docs/examples/conditional_luminosity_function_examples.rst b/docs/examples/conditional_luminosity_function_examples.rst new file mode 100644 index 0000000..62f1ba5 --- /dev/null +++ b/docs/examples/conditional_luminosity_function_examples.rst @@ -0,0 +1,788 @@ +.. |lfkitlogo| image:: /_static/logos/lfkit_logo-icon.png + :alt: LFKit logo + :width: 50px + +|lfkitlogo| Conditional luminosity-function examples +==================================================== + +This page shows how to use :class:`lfkit.LuminosityFunction` to construct and +evaluate conditional luminosity functions. + +A conditional luminosity function has the form :math:`\Phi(M \mid x)`, where +:math:`M` is absolute magnitude and :math:`x` is an external conditioning +variable. In these examples, the public :class:`lfkit.LuminosityFunction` +interface uses redshift as the conditioning variable through ``phi(M, z)``. + +Conditional luminosity functions are useful when the luminosity function +parameters depend on another quantity, such as redshift, halo mass, galaxy +type, environment, richness, or stellar mass. The examples below focus on +redshift-dependent behaviour because it fits naturally into the main LFKit API. + +The examples include conditional Schechter models and central/satellite +components. The central component is represented by a narrow lognormal term, +while the satellite component is represented by a modified Schechter-like term. + +All examples below are executable via ``.. plot::``. + +The number-density units follow the normalization of the luminosity function. +For example, if the amplitudes are supplied in comoving +:math:`{\rm Mpc}^{-3}\,{\rm mag}^{-1}`, then :math:`\Phi(M \mid z)` has units of +:math:`{\rm Mpc}^{-3}\,{\rm mag}^{-1}`. + + +Conditional Schechter luminosity function +----------------------------------------- + +A conditional Schechter luminosity function allows one or more Schechter +parameters to depend on the conditioning variable. + +This example makes the normalization and characteristic magnitude depend on +redshift. The faint-end slope is kept fixed. The result is still a Schechter +luminosity function at each redshift, but the curve changes as the conditioning +variable changes. + +This is useful when a single fixed luminosity function is too restrictive but a +full tabulated model is not needed. + +.. plot:: + :include-source: True + :width: 520 + + import numpy as np + import matplotlib.pyplot as plt + import cmasher as cmr + + from lfkit import LuminosityFunction + + LABEL_SIZE = 15 + TICK_SIZE = 13 + TITLE_SIZE = 17 + LEGEND_SIZE = 15 + + colors_blue = cmr.take_cmap_colors("cmr.guppy", 3, cmap_range=(0.72, 0.96)) + + absolute_mag = np.linspace(-24.0, -14.0, 500) + redshifts = [0.1, 0.6, 1.1] + + lf = LuminosityFunction.conditional_schechter( + phi_star=lambda z: 1.0e-3 * (1.0 + z) ** 0.8, + m_star=lambda z: -20.5 - 0.7 * (z - 0.1), + alpha=-1.1, + ) + + fig, ax = plt.subplots(figsize=(7.0, 5.0)) + + for z_value, color in zip(redshifts, colors_blue): + z = np.full_like(absolute_mag, z_value) + phi = lf.phi(absolute_mag, z) + + ax.plot( + absolute_mag, + phi, + lw=3, + color=color, + label=rf"$z={z_value}$", + ) + + ax.set_yscale("log") + ax.invert_xaxis() + ax.set_xlabel("Absolute magnitude $M$", fontsize=LABEL_SIZE) + ax.set_ylabel(r"$\Phi(M \mid z)$ [$\mathrm{Mpc}^{-3}\,\mathrm{mag}^{-1}$]", fontsize=LABEL_SIZE) + ax.set_title("Conditional Schechter luminosity function", fontsize=TITLE_SIZE) + ax.tick_params(axis="both", labelsize=TICK_SIZE) + ax.legend(frameon=True, fontsize=LEGEND_SIZE, loc="best") + plt.tight_layout() + + +Conditional Schechter surface +----------------------------- + +The same conditional Schechter model can be visualized across the full +magnitude-redshift plane. + +The filled colour scale shows :math:`\log_{10}\Phi(M \mid z)`. The white +contours mark constant :math:`\log_{10}\Phi(M \mid z)` levels at -5, -4, -3, +and -2. These contours make it easier to see where equal-abundance regions sit +as both magnitude and redshift change. + +This plot is useful for checking that the conditional model varies smoothly +across the range where it will be used. + +.. plot:: + :include-source: True + :width: 560 + + import numpy as np + import matplotlib.pyplot as plt + import cmasher as cmr + + from lfkit import LuminosityFunction + + LABEL_SIZE = 15 + TICK_SIZE = 13 + TITLE_SIZE = 17 + LEGEND_SIZE = 15 + + absolute_mag = np.linspace(-24.0, -16.0, 220) + redshift = np.linspace(0.0, 1.5, 180) + + mag_grid, z_grid = np.meshgrid(absolute_mag, redshift) + + lf = LuminosityFunction.conditional_schechter( + phi_star=lambda z: 1.0e-3 * (1.0 + z) ** 0.8, + m_star=lambda z: -20.5 - 0.7 * (z - 0.1), + alpha=-1.1, + ) + + phi = lf.phi(mag_grid, z_grid) + log_phi = np.log10(phi) + + fig, ax = plt.subplots(figsize=(7.2, 5.0)) + + mesh = ax.pcolormesh( + absolute_mag, + redshift, + log_phi, + shading="auto", + cmap="cmr.guppy", + ) + + contour_levels = [-5.0, -4.0, -3.0, -2.0] + contours = ax.contour( + absolute_mag, + redshift, + log_phi, + levels=contour_levels, + colors="white", + linewidths=1.2, + ) + ax.clabel(contours, inline=True, fontsize=TICK_SIZE, fmt=r"$10^{%.0f}$") + + ax.invert_xaxis() + ax.set_xlabel("Absolute magnitude $M$", fontsize=LABEL_SIZE) + ax.set_ylabel("Redshift $z$", fontsize=LABEL_SIZE) + ax.set_title("Conditional Schechter LF surface", fontsize=TITLE_SIZE) + ax.tick_params(axis="both", labelsize=TICK_SIZE) + + cbar = fig.colorbar(mesh, ax=ax) + cbar.set_label( + r"$\log_{10}\Phi(M \mid z)$ [$\log_{10}(\mathrm{Mpc}^{-3}\,\mathrm{mag}^{-1})$]", + fontsize=LABEL_SIZE, + ) + cbar.ax.tick_params(labelsize=TICK_SIZE) + + plt.tight_layout() + + +Conditional evolving Schechter model +------------------------------------ + +A conditional evolving Schechter model uses LFKit's registered parameter models +to define how :math:`\phi_\star`, :math:`M_\star`, and :math:`\alpha` depend on +the conditioning variable. + +This example is similar to the standard evolving Schechter interface, but it is +included here because it can be used in the same conditional-LF workflow. The +conditioning variable is passed through ``phi(M, z)``. + +This is useful when the desired parameter evolution already matches one of +LFKit's registered parameter models. + +.. plot:: + :include-source: True + :width: 520 + + import numpy as np + import matplotlib.pyplot as plt + import cmasher as cmr + + from lfkit import LuminosityFunction + + LABEL_SIZE = 15 + TICK_SIZE = 13 + TITLE_SIZE = 17 + LEGEND_SIZE = 15 + + colors_red = cmr.take_cmap_colors("cmr.guppy", 3, cmap_range=(0.03, 0.26)) + + absolute_mag = np.linspace(-24.0, -14.0, 500) + redshifts = [0.1, 0.6, 1.1] + + lf = LuminosityFunction.conditional_evolving_schechter( + phi_model="linear_p", + phi_kwargs={"phi_0_star": 1.0e-3, "p": 0.7}, + m_star_model="linear_q", + m_star_kwargs={"m_0_star": -20.5, "q": 0.8, "z_ref": 0.1}, + alpha_model="constant", + alpha_kwargs={"alpha": -1.1}, + ) + + fig, ax = plt.subplots(figsize=(7.0, 5.0)) + + for z_value, color in zip(redshifts, colors_red): + z = np.full_like(absolute_mag, z_value) + phi = lf.phi(absolute_mag, z) + + ax.plot( + absolute_mag, + phi, + lw=3, + color=color, + label=rf"$z={z_value}$", + ) + + ax.set_yscale("log") + ax.invert_xaxis() + ax.set_xlabel("Absolute magnitude $M$", fontsize=LABEL_SIZE) + ax.set_ylabel(r"$\Phi(M \mid z)$ [$\mathrm{Mpc}^{-3}\,\mathrm{mag}^{-1}$]", fontsize=LABEL_SIZE) + ax.set_title("Conditional evolving Schechter model", fontsize=TITLE_SIZE) + ax.tick_params(axis="both", labelsize=TICK_SIZE) + ax.legend(frameon=True, fontsize=LEGEND_SIZE, loc="best") + plt.tight_layout() + + +Central lognormal conditional luminosity function +------------------------------------------------- + +A central-galaxy conditional luminosity function can be represented by a narrow +lognormal component in luminosity, written here in magnitude space. + +This example shows a central component whose mean absolute magnitude becomes +brighter with redshift. The scatter is kept fixed. The peak therefore shifts in +absolute magnitude while retaining a similar width. + +This type of component is useful in central/satellite decompositions where the +central galaxy population is concentrated around a characteristic luminosity at +fixed condition. + +.. plot:: + :include-source: True + :width: 520 + + import numpy as np + import matplotlib.pyplot as plt + import cmasher as cmr + + from lfkit import LuminosityFunction + + LABEL_SIZE = 15 + TICK_SIZE = 13 + TITLE_SIZE = 17 + LEGEND_SIZE = 15 + + colors_blue = cmr.take_cmap_colors("cmr.guppy", 3, cmap_range=(0.72, 0.96)) + + absolute_mag = np.linspace(-24.0, -16.0, 500) + redshifts = [0.1, 0.6, 1.1] + + lf = LuminosityFunction.central_lognormal_conditional( + mean_absolute_mag=lambda z: -20.8 - 0.6 * (z - 0.1), + sigma_log_luminosity=0.18, + amplitude=lambda z: 8.0e-4 * (1.0 + z) ** 0.4, + ) + + fig, ax = plt.subplots(figsize=(7.0, 5.0)) + + for z_value, color in zip(redshifts, colors_blue): + z = np.full_like(absolute_mag, z_value) + phi = lf.phi(absolute_mag, z) + + ax.plot( + absolute_mag, + phi, + lw=3, + color=color, + label=rf"$z={z_value}$", + ) + + ax.set_yscale("log") + ax.invert_xaxis() + ax.set_xlabel("Absolute magnitude $M$", fontsize=LABEL_SIZE) + ax.set_ylabel(r"$\Phi_{\rm cen}(M \mid z)$ [$\mathrm{Mpc}^{-3}\,\mathrm{mag}^{-1}$]", fontsize=LABEL_SIZE) + ax.set_title("Central lognormal conditional LF", fontsize=TITLE_SIZE) + ax.tick_params(axis="both", labelsize=TICK_SIZE) + ax.legend(frameon=True, fontsize=LEGEND_SIZE, loc="best") + plt.tight_layout() + + +Satellite modified-Schechter conditional luminosity function +------------------------------------------------------------ + +A satellite conditional luminosity function can be represented by a modified +Schechter-like component. + +This example uses a satellite model with an exponential cutoff +:math:`\exp[-(L/L_\star)^2]`. Compared with the central lognormal component, +the satellite component is broader and contributes more strongly across a wider +range of faint magnitudes. + +This is useful for modeling satellite populations separately from central +galaxies. + +.. plot:: + :include-source: True + :width: 520 + + import numpy as np + import matplotlib.pyplot as plt + import cmasher as cmr + + from lfkit import LuminosityFunction + + LABEL_SIZE = 15 + TICK_SIZE = 13 + TITLE_SIZE = 17 + LEGEND_SIZE = 15 + + colors_red = cmr.take_cmap_colors("cmr.guppy", 3, cmap_range=(0.03, 0.26)) + + absolute_mag = np.linspace(-24.0, -14.0, 500) + redshifts = [0.1, 0.6, 1.1] + + lf = LuminosityFunction.satellite_modified_schechter_conditional( + phi_star=lambda z: 1.2e-3 * (1.0 + z) ** 0.5, + m_star=lambda z: -19.9 - 0.5 * (z - 0.1), + alpha=lambda z: -1.05 - 0.10 * z, + ) + + fig, ax = plt.subplots(figsize=(7.0, 5.0)) + + for z_value, color in zip(redshifts, colors_red): + z = np.full_like(absolute_mag, z_value) + phi = lf.phi(absolute_mag, z) + + ax.plot( + absolute_mag, + phi, + lw=3, + color=color, + label=rf"$z={z_value}$", + ) + + ax.set_yscale("log") + ax.invert_xaxis() + ax.set_xlabel("Absolute magnitude $M$", fontsize=LABEL_SIZE) + ax.set_ylabel(r"$\Phi_{\rm sat}(M \mid z)$ [$\mathrm{Mpc}^{-3}\,\mathrm{mag}^{-1}$]", fontsize=LABEL_SIZE) + ax.set_title("Satellite modified-Schechter conditional LF", fontsize=TITLE_SIZE) + ax.tick_params(axis="both", labelsize=TICK_SIZE) + ax.legend(frameon=True, fontsize=LEGEND_SIZE, loc="best") + plt.tight_layout() + + +Central and satellite components +-------------------------------- + +The central and satellite conditional luminosity functions can be combined into +a single central-plus-satellite model. + +This plot separates the central component, the satellite component, and their +sum at a fixed redshift. The central component is narrow and localized around +the mean central magnitude. The satellite component is broader and dominates +over a wider faint-magnitude range. + +This decomposition is useful when checking which part of the galaxy population +drives the total luminosity function. + +.. plot:: + :include-source: True + :width: 520 + + import numpy as np + import matplotlib.pyplot as plt + import cmasher as cmr + + from lfkit import LuminosityFunction + + LABEL_SIZE = 15 + TICK_SIZE = 13 + TITLE_SIZE = 17 + LEGEND_SIZE = 15 + + colors_blue = cmr.take_cmap_colors("cmr.guppy", 3, cmap_range=(0.72, 0.96)) + colors_red = cmr.take_cmap_colors("cmr.guppy", 3, cmap_range=(0.03, 0.26)) + c_mid = 0.5 * (np.array(colors_blue[1]) + np.array(colors_red[1])) + + absolute_mag = np.linspace(-24.0, -14.0, 500) + z_value = 0.6 + z = np.full_like(absolute_mag, z_value) + + central = LuminosityFunction.central_lognormal_conditional( + mean_absolute_mag=lambda z: -20.8 - 0.6 * (z - 0.1), + sigma_log_luminosity=0.18, + amplitude=lambda z: 8.0e-4 * (1.0 + z) ** 0.4, + ) + + satellite = LuminosityFunction.satellite_modified_schechter_conditional( + phi_star=lambda z: 1.2e-3 * (1.0 + z) ** 0.5, + m_star=lambda z: -19.9 - 0.5 * (z - 0.1), + alpha=lambda z: -1.05 - 0.10 * z, + ) + + total = LuminosityFunction.central_satellite_conditional( + central_mean_absolute_mag=lambda z: -20.8 - 0.6 * (z - 0.1), + central_sigma_log_luminosity=0.18, + central_amplitude=lambda z: 8.0e-4 * (1.0 + z) ** 0.4, + satellite_phi_star=lambda z: 1.2e-3 * (1.0 + z) ** 0.5, + satellite_alpha=lambda z: -1.05 - 0.10 * z, + satellite_m_star=lambda z: -19.9 - 0.5 * (z - 0.1), + ) + + central_phi = central.phi(absolute_mag, z) + satellite_phi = satellite.phi(absolute_mag, z) + total_phi = total.phi(absolute_mag, z) + + fig, ax = plt.subplots(figsize=(7.0, 5.0)) + + ax.plot( + absolute_mag, + central_phi, + lw=3, + color=colors_blue[1], + label="Central", + ) + ax.plot( + absolute_mag, + satellite_phi, + lw=3, + color=colors_red[1], + label="Satellite", + ) + ax.plot( + absolute_mag, + total_phi, + lw=3, + color=c_mid, + label="Central + satellite", + ) + + ax.set_yscale("log") + ax.invert_xaxis() + ax.set_xlabel("Absolute magnitude $M$", fontsize=LABEL_SIZE) + ax.set_ylabel(r"$\Phi(M \mid z)$ [$\mathrm{Mpc}^{-3}\,\mathrm{mag}^{-1}$]", fontsize=LABEL_SIZE) + ax.set_title(r"Central and satellite components at $z=0.6$", fontsize=TITLE_SIZE) + ax.tick_params(axis="both", labelsize=TICK_SIZE) + ax.legend(frameon=True, fontsize=LEGEND_SIZE, loc="best") + plt.tight_layout() + + +Central-plus-satellite evolution +-------------------------------- + +The central-plus-satellite model can also be evaluated at several redshifts. + +This plot shows how the total conditional luminosity function changes when the +central and satellite parameters vary with redshift. The model combines a +narrow central component with a broader satellite component at each redshift. + +This is useful when the total luminosity function is needed, but the model +still keeps a physically interpretable central/satellite split. + +.. plot:: + :include-source: True + :width: 520 + + import numpy as np + import matplotlib.pyplot as plt + import cmasher as cmr + + from lfkit import LuminosityFunction + + LABEL_SIZE = 15 + TICK_SIZE = 13 + TITLE_SIZE = 17 + LEGEND_SIZE = 15 + + colors_blue = cmr.take_cmap_colors("cmr.guppy", 3, cmap_range=(0.72, 0.96)) + + absolute_mag = np.linspace(-24.0, -14.0, 500) + redshifts = [0.1, 0.6, 1.1] + + lf = LuminosityFunction.central_satellite_conditional( + central_mean_absolute_mag=lambda z: -20.8 - 0.6 * (z - 0.1), + central_sigma_log_luminosity=0.18, + central_amplitude=lambda z: 8.0e-4 * (1.0 + z) ** 0.4, + satellite_phi_star=lambda z: 1.2e-3 * (1.0 + z) ** 0.5, + satellite_alpha=lambda z: -1.05 - 0.10 * z, + satellite_m_star=lambda z: -19.9 - 0.5 * (z - 0.1), + ) + + fig, ax = plt.subplots(figsize=(7.0, 5.0)) + + for z_value, color in zip(redshifts, colors_blue): + z = np.full_like(absolute_mag, z_value) + phi = lf.phi(absolute_mag, z) + + ax.plot( + absolute_mag, + phi, + lw=3, + color=color, + label=rf"$z={z_value}$", + ) + + ax.set_yscale("log") + ax.invert_xaxis() + ax.set_xlabel("Absolute magnitude $M$", fontsize=LABEL_SIZE) + ax.set_ylabel(r"$\Phi_{\rm cen+sat}(M \mid z)$ [$\mathrm{Mpc}^{-3}\,\mathrm{mag}^{-1}$]", fontsize=LABEL_SIZE) + ax.set_title("Central-plus-satellite conditional LF", fontsize=TITLE_SIZE) + ax.tick_params(axis="both", labelsize=TICK_SIZE) + ax.legend(frameon=True, fontsize=LEGEND_SIZE, loc="best") + plt.tight_layout() + + +Integrated conditional number density +------------------------------------- + +A conditional luminosity function can be integrated over absolute magnitude at +each value of the conditioning variable. + +This example integrates the central component, satellite component, and total +central-plus-satellite luminosity function over a fixed absolute-magnitude +range. The result shows how the selected number density changes with redshift. + +This type of calculation is useful when the conditional luminosity function is +used as an ingredient in redshift-distribution or abundance calculations. + +.. plot:: + :include-source: True + :width: 520 + + import numpy as np + import matplotlib.pyplot as plt + import cmasher as cmr + + from lfkit import LuminosityFunction + + LABEL_SIZE = 15 + TICK_SIZE = 13 + TITLE_SIZE = 17 + LEGEND_SIZE = 15 + + colors_blue = cmr.take_cmap_colors("cmr.guppy", 3, cmap_range=(0.72, 0.96)) + colors_red = cmr.take_cmap_colors("cmr.guppy", 3, cmap_range=(0.03, 0.26)) + c_mid = 0.5 * (np.array(colors_blue[1]) + np.array(colors_red[1])) + + redshift = np.linspace(0.05, 1.5, 180) + + central = LuminosityFunction.central_lognormal_conditional( + mean_absolute_mag=lambda z: -20.8 - 0.6 * (z - 0.1), + sigma_log_luminosity=0.18, + amplitude=lambda z: 8.0e-4 * (1.0 + z) ** 0.4, + ) + + satellite = LuminosityFunction.satellite_modified_schechter_conditional( + phi_star=lambda z: 1.2e-3 * (1.0 + z) ** 0.5, + m_star=lambda z: -19.9 - 0.5 * (z - 0.1), + alpha=lambda z: -1.05 - 0.10 * z, + ) + + total = LuminosityFunction.central_satellite_conditional( + central_mean_absolute_mag=lambda z: -20.8 - 0.6 * (z - 0.1), + central_sigma_log_luminosity=0.18, + central_amplitude=lambda z: 8.0e-4 * (1.0 + z) ** 0.4, + satellite_phi_star=lambda z: 1.2e-3 * (1.0 + z) ** 0.5, + satellite_alpha=lambda z: -1.05 - 0.10 * z, + satellite_m_star=lambda z: -19.9 - 0.5 * (z - 0.1), + ) + + n_central = central.integrated_number_density( + redshift, + m_bright=-24.0, + m_faint=-14.0, + n_m=800, + ) + n_satellite = satellite.integrated_number_density( + redshift, + m_bright=-24.0, + m_faint=-14.0, + n_m=800, + ) + n_total = total.integrated_number_density( + redshift, + m_bright=-24.0, + m_faint=-14.0, + n_m=800, + ) + + fig, ax = plt.subplots(figsize=(7.0, 5.0)) + + ax.plot(redshift, n_central, lw=3, color=colors_blue[1], label="Central") + ax.plot(redshift, n_satellite, lw=3, color=colors_red[1], label="Satellite") + ax.plot(redshift, n_total, lw=3, color=c_mid, label="Central + satellite") + + ax.set_yscale("log") + ax.set_xlabel("Redshift $z$", fontsize=LABEL_SIZE) + ax.set_ylabel(r"Integrated number density [$\mathrm{Mpc}^{-3}$]", fontsize=LABEL_SIZE) + ax.set_title(r"Integrated conditional LF over $-24 \leq M \leq -14$", fontsize=TITLE_SIZE) + ax.tick_params(axis="both", labelsize=TICK_SIZE) + ax.legend(frameon=True, fontsize=LEGEND_SIZE, loc="best") + plt.tight_layout() + + +Central fraction +---------------- + +The relative contribution of central and satellite galaxies can be summarized +as a fraction of the integrated central-plus-satellite luminosity function. + +This example computes the central fraction over a fixed absolute-magnitude +range. The fraction changes with redshift because the central and satellite +components have different conditional parameter dependence. + +This is a compact diagnostic for checking whether the selected population is +central-dominated, satellite-dominated, or mixed. + +.. plot:: + :include-source: True + :width: 520 + + import numpy as np + import matplotlib.pyplot as plt + import cmasher as cmr + + from lfkit import LuminosityFunction + + LABEL_SIZE = 15 + TICK_SIZE = 13 + TITLE_SIZE = 17 + LEGEND_SIZE = 15 + + colors_blue = cmr.take_cmap_colors("cmr.guppy", 3, cmap_range=(0.72, 0.96)) + colors_red = cmr.take_cmap_colors("cmr.guppy", 3, cmap_range=(0.03, 0.26)) + + redshift = np.linspace(0.05, 1.5, 180) + + central = LuminosityFunction.central_lognormal_conditional( + mean_absolute_mag=lambda z: -20.8 - 0.6 * (z - 0.1), + sigma_log_luminosity=0.18, + amplitude=lambda z: 8.0e-4 * (1.0 + z) ** 0.4, + ) + + satellite = LuminosityFunction.satellite_modified_schechter_conditional( + phi_star=lambda z: 1.2e-3 * (1.0 + z) ** 0.5, + m_star=lambda z: -19.9 - 0.5 * (z - 0.1), + alpha=lambda z: -1.05 - 0.10 * z, + ) + + n_central = central.integrated_number_density( + redshift, + m_bright=-24.0, + m_faint=-14.0, + n_m=800, + ) + n_satellite = satellite.integrated_number_density( + redshift, + m_bright=-24.0, + m_faint=-14.0, + n_m=800, + ) + + central_fraction = n_central / (n_central + n_satellite) + satellite_fraction = n_satellite / (n_central + n_satellite) + + fig, ax = plt.subplots(figsize=(7.0, 5.0)) + + ax.plot( + redshift, + central_fraction, + lw=3, + color=colors_blue[1], + label="Central fraction", + ) + ax.plot( + redshift, + satellite_fraction, + lw=3, + color=colors_red[1], + label="Satellite fraction", + ) + + ax.set_ylim(-0.05, 1.05) + ax.set_xlabel("Redshift $z$", fontsize=LABEL_SIZE) + ax.set_ylabel("Fraction of integrated LF", fontsize=LABEL_SIZE) + ax.set_title(r"Central and satellite fractions over $-24 \leq M \leq -14$", fontsize=TITLE_SIZE) + ax.tick_params(axis="both", labelsize=TICK_SIZE) + ax.legend(frameon=True, fontsize=LEGEND_SIZE, loc="center right") + plt.tight_layout() + + +Central-plus-satellite surface +------------------------------ + +The full central-plus-satellite conditional luminosity function can be shown as +a two-dimensional surface. + +The filled colour scale shows :math:`\log_{10}\Phi_{\rm cen+sat}(M \mid z)`. +The white contours mark constant :math:`\log_{10}\Phi_{\rm cen+sat}(M \mid z)` +levels at -5, -4, -3, and -2. + +This view is useful for checking whether the central peak, satellite tail, and +redshift dependence combine smoothly over the full magnitude-redshift range. + +.. plot:: + :include-source: True + :width: 560 + + import numpy as np + import matplotlib.pyplot as plt + import cmasher as cmr + + from lfkit import LuminosityFunction + + LABEL_SIZE = 15 + TICK_SIZE = 13 + TITLE_SIZE = 17 + LEGEND_SIZE = 15 + + absolute_mag = np.linspace(-24.0, -14.0, 220) + redshift = np.linspace(0.0, 1.5, 180) + + mag_grid, z_grid = np.meshgrid(absolute_mag, redshift) + + lf = LuminosityFunction.central_satellite_conditional( + central_mean_absolute_mag=lambda z: -20.8 - 0.6 * (z - 0.1), + central_sigma_log_luminosity=0.18, + central_amplitude=lambda z: 8.0e-4 * (1.0 + z) ** 0.4, + satellite_phi_star=lambda z: 1.2e-3 * (1.0 + z) ** 0.5, + satellite_alpha=lambda z: -1.05 - 0.10 * z, + satellite_m_star=lambda z: -19.9 - 0.5 * (z - 0.1), + ) + + phi = lf.phi(mag_grid, z_grid) + log_phi = np.log10(phi) + + fig, ax = plt.subplots(figsize=(7.2, 5.0)) + + mesh = ax.pcolormesh( + absolute_mag, + redshift, + log_phi, + shading="auto", + cmap="cmr.guppy", + ) + + contour_levels = [-5.0, -4.0, -3.0, -2.0] + contours = ax.contour( + absolute_mag, + redshift, + log_phi, + levels=contour_levels, + colors="white", + linewidths=1.2, + ) + ax.clabel(contours, inline=True, fontsize=TICK_SIZE, fmt=r"$10^{%.0f}$") + + ax.invert_xaxis() + ax.set_xlabel("Absolute magnitude $M$", fontsize=LABEL_SIZE) + ax.set_ylabel("Redshift $z$", fontsize=LABEL_SIZE) + ax.set_title("Central-plus-satellite conditional LF surface", fontsize=TITLE_SIZE) + ax.tick_params(axis="both", labelsize=TICK_SIZE) + + cbar = fig.colorbar(mesh, ax=ax) + cbar.set_label( + r"$\log_{10}\Phi_{\rm cen+sat}(M \mid z)$ [$\log_{10}(\mathrm{Mpc}^{-3}\,\mathrm{mag}^{-1})$]", + fontsize=LABEL_SIZE, + ) + cbar.ax.tick_params(labelsize=TICK_SIZE) + + plt.tight_layout() diff --git a/docs/examples/index.rst b/docs/examples/index.rst index 6237843..02ce404 100644 --- a/docs/examples/index.rst +++ b/docs/examples/index.rst @@ -49,6 +49,19 @@ Working, executable examples for using LFKit. +++ *API:* LuminosityFunction + .. grid-item-card:: + :link: conditional_luminosity_function_examples + :link-type: doc + :shadow: md + + **Conditional luminosity-function examples** + ^^^ + Build conditional luminosity functions, including redshift-dependent + Schechter models and central/satellite components. + + +++ + *API:* LuminosityFunction + .. grid-item-card:: :link: catalog_completeness_examples :link-type: doc @@ -69,4 +82,5 @@ Working, executable examples for using LFKit. kcorrect_examples poggianti_examples luminosity_function_examples - catalog_completeness_examples \ No newline at end of file + conditional_luminosity_function_examples + catalog_completeness_examples From a441f3e01cc9c192ebabb7874cf9ac3613deef48 Mon Sep 17 00:00:00 2001 From: "niko.sarcevic" Date: Tue, 12 May 2026 00:16:47 -0400 Subject: [PATCH 4/9] Fix Poggianti97 downloader tests --- src/lfkit/utils/download_poggianti97_data.py | 6 ++++-- tests/test_utils_download_pogg97_data.py | 4 ++-- 2 files changed, 6 insertions(+), 4 deletions(-) diff --git a/src/lfkit/utils/download_poggianti97_data.py b/src/lfkit/utils/download_poggianti97_data.py index 83dcd7f..b374101 100644 --- a/src/lfkit/utils/download_poggianti97_data.py +++ b/src/lfkit/utils/download_poggianti97_data.py @@ -28,8 +28,10 @@ CATALOG = "J/A+AS/122/399" -# Package data directory where CSVs are stored -PKG_DATADIR = Path("src/lfkit/data/poggianti1997") +# Package data directory where CSVs are stored. +# Resolve relative to this file so running tests/scripts from repo root +# cannot accidentally write into the current working directory. +PKG_DATADIR = Path(__file__).resolve().parents[1] / "data" / "poggianti1997" def main() -> None: diff --git a/tests/test_utils_download_pogg97_data.py b/tests/test_utils_download_pogg97_data.py index adad66a..65446a9 100644 --- a/tests/test_utils_download_pogg97_data.py +++ b/tests/test_utils_download_pogg97_data.py @@ -63,7 +63,7 @@ def test_main_writes_csv_files(tmp_path, monkeypatch): assert len(written) == 2 -def test_main_uses_correct_catalog(monkeypatch): +def test_main_uses_correct_catalog(tmp_path, monkeypatch): """Tests that main calls Vizier.get_catalogs with the expected catalog identifier.""" called = {} @@ -72,7 +72,7 @@ def fake_get_catalogs(catalog): return {} monkeypatch.setattr(dl.Vizier, "get_catalogs", fake_get_catalogs) - monkeypatch.setattr(dl, "PKG_DATADIR", Path("unused")) + monkeypatch.setattr(dl, "PKG_DATADIR", tmp_path / "unused") dl.main() From bf6daacc640cc320f707f4201f02c4dcd6fb9a30 Mon Sep 17 00:00:00 2001 From: "niko.sarcevic" Date: Tue, 12 May 2026 00:16:51 -0400 Subject: [PATCH 5/9] Update project test configuration --- pyproject.toml | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/pyproject.toml b/pyproject.toml index 81c04de..966659f 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -65,9 +65,10 @@ exclude = ["src/lfkit/_version.py"] [tool.pytest.ini_options] testpaths = ["tests"] -addopts = "-ra -m 'not slow'" +addopts = "-ra -m 'not slow and not benchmark'" markers = [ "slow: marks tests as slow (deselect with '-m \"not slow\"')", + "benchmark: reference benchmark tests that are not part of the core unit-test suite", ] [tool.tox] From 6b0eade1ad9765d9456feeefa4eee1637890cdda Mon Sep 17 00:00:00 2001 From: "niko.sarcevic" Date: Wed, 13 May 2026 23:15:09 -0400 Subject: [PATCH 6/9] updated exampels in docs: typos, added more examples --- docs/about/lf_overview.rst | 2 +- .../catalog_completeness_examples.rst | 4 +- ...nditional_luminosity_function_examples.rst | 1031 ++++++++++++----- docs/examples/index.rst | 2 +- docs/index.rst | 4 +- .../benchmarks/test_cacciato_hod_reference.py | 0 6 files changed, 756 insertions(+), 287 deletions(-) create mode 100644 tests/benchmarks/test_cacciato_hod_reference.py diff --git a/docs/about/lf_overview.rst b/docs/about/lf_overview.rst index 53b50eb..7adf3c9 100644 --- a/docs/about/lf_overview.rst +++ b/docs/about/lf_overview.rst @@ -303,7 +303,7 @@ and normalize the result over the supplied redshift grid. ) This is useful for constructing LF-dependent redshift trends. -Without a volume weight, the result is the luminosity-function +Without a volume weight, the result is the luminosity function selection factor only, not a full survey redshift distribution. A volume-weighted trend can be computed by passing a callable: diff --git a/docs/examples/catalog_completeness_examples.rst b/docs/examples/catalog_completeness_examples.rst index 48e9e47..307a96e 100644 --- a/docs/examples/catalog_completeness_examples.rst +++ b/docs/examples/catalog_completeness_examples.rst @@ -2,8 +2,8 @@ :alt: LFKit logo :width: 50px -|lfkitlogo| Catalog completeness examples -========================================= +|lfkitlogo| Catalog completeness +================================ This page shows how to use :class:`lfkit.LuminosityFunction` to compute observed and missing galaxy number densities for a magnitude-limited catalog. diff --git a/docs/examples/conditional_luminosity_function_examples.rst b/docs/examples/conditional_luminosity_function_examples.rst index 62f1ba5..f01be1c 100644 --- a/docs/examples/conditional_luminosity_function_examples.rst +++ b/docs/examples/conditional_luminosity_function_examples.rst @@ -2,47 +2,46 @@ :alt: LFKit logo :width: 50px -|lfkitlogo| Conditional luminosity-function examples -==================================================== +|lfkitlogo| Conditional luminosity functions +============================================ -This page shows how to use :class:`lfkit.LuminosityFunction` to construct and -evaluate conditional luminosity functions. +This page shows how to evaluate conditional luminosity functions with LFKit. A conditional luminosity function has the form :math:`\Phi(M \mid x)`, where :math:`M` is absolute magnitude and :math:`x` is an external conditioning -variable. In these examples, the public :class:`lfkit.LuminosityFunction` -interface uses redshift as the conditioning variable through ``phi(M, z)``. +variable. The conditioning variable is generic: it can represent redshift, halo +mass, environment, galaxy type, richness, stellar mass, or another quantity. -Conditional luminosity functions are useful when the luminosity function -parameters depend on another quantity, such as redshift, halo mass, galaxy -type, environment, richness, or stellar mass. The examples below focus on -redshift-dependent behaviour because it fits naturally into the main LFKit API. +The examples below use redshift as the conditioning variable because it is a +natural choice for luminosity function evolution. The same functions can be +used with any other conditioning variable by replacing ``z`` with the desired +quantity. -The examples include conditional Schechter models and central/satellite -components. The central component is represented by a narrow lognormal term, -while the satellite component is represented by a modified Schechter-like term. +The examples include: -All examples below are executable via ``.. plot::``. +* a conditional Schechter luminosity function, +* a conditional Schechter model using LFKit parameter models, +* a lognormal component, +* a modified Schechter-like component, +* a two-component lognormal plus modified-Schechter model, +* integrated number densities and component fractions. -The number-density units follow the normalization of the luminosity function. -For example, if the amplitudes are supplied in comoving -:math:`{\rm Mpc}^{-3}\,{\rm mag}^{-1}`, then :math:`\Phi(M \mid z)` has units of -:math:`{\rm Mpc}^{-3}\,{\rm mag}^{-1}`. +The number-density units follow the normalization supplied to the luminosity +function. For example, if the amplitudes are supplied in +:math:`{\rm Mpc}^{-3}\,{\rm mag}^{-1}`, then :math:`\Phi(M \mid z)` has units +of :math:`{\rm Mpc}^{-3}\,{\rm mag}^{-1}`. Conditional Schechter luminosity function ----------------------------------------- -A conditional Schechter luminosity function allows one or more Schechter -parameters to depend on the conditioning variable. +A conditional Schechter luminosity function allows the Schechter parameters to +depend on an external variable. This example makes the normalization and characteristic magnitude depend on -redshift. The faint-end slope is kept fixed. The result is still a Schechter -luminosity function at each redshift, but the curve changes as the conditioning -variable changes. - -This is useful when a single fixed luminosity function is too restrictive but a -full tabulated model is not needed. +redshift. The faint-end slope is kept fixed. At each redshift, the model is +still a Schechter luminosity function, but the curve evolves smoothly with the +conditioning variable. .. plot:: :include-source: True @@ -52,29 +51,32 @@ full tabulated model is not needed. import matplotlib.pyplot as plt import cmasher as cmr - from lfkit import LuminosityFunction + from lfkit.photometry.conditional_lf_models import ( + conditional_schechter, + ) LABEL_SIZE = 15 TICK_SIZE = 13 TITLE_SIZE = 17 LEGEND_SIZE = 15 - colors_blue = cmr.take_cmap_colors("cmr.guppy", 3, cmap_range=(0.72, 0.96)) + colors = cmr.take_cmap_colors("cmr.guppy", 3, cmap_range=(0.72, 0.96)) absolute_mag = np.linspace(-24.0, -14.0, 500) redshifts = [0.1, 0.6, 1.1] - lf = LuminosityFunction.conditional_schechter( - phi_star=lambda z: 1.0e-3 * (1.0 + z) ** 0.8, - m_star=lambda z: -20.5 - 0.7 * (z - 0.1), - alpha=-1.1, - ) - fig, ax = plt.subplots(figsize=(7.0, 5.0)) - for z_value, color in zip(redshifts, colors_blue): + for z_value, color in zip(redshifts, colors): z = np.full_like(absolute_mag, z_value) - phi = lf.phi(absolute_mag, z) + + phi = conditional_schechter( + absolute_mag, + z, + phi_star=lambda z: 1.0e-3 * (1.0 + z) ** 0.8, + m_star=lambda z: -20.5 - 0.7 * (z - 0.1), + alpha=-1.1, + ) ax.plot( absolute_mag, @@ -87,7 +89,10 @@ full tabulated model is not needed. ax.set_yscale("log") ax.invert_xaxis() ax.set_xlabel("Absolute magnitude $M$", fontsize=LABEL_SIZE) - ax.set_ylabel(r"$\Phi(M \mid z)$ [$\mathrm{Mpc}^{-3}\,\mathrm{mag}^{-1}$]", fontsize=LABEL_SIZE) + ax.set_ylabel( + r"$\Phi(M \mid z)$ [$\mathrm{Mpc}^{-3}\,\mathrm{mag}^{-1}$]", + fontsize=LABEL_SIZE, + ) ax.set_title("Conditional Schechter luminosity function", fontsize=TITLE_SIZE) ax.tick_params(axis="both", labelsize=TICK_SIZE) ax.legend(frameon=True, fontsize=LEGEND_SIZE, loc="best") @@ -97,16 +102,11 @@ full tabulated model is not needed. Conditional Schechter surface ----------------------------- -The same conditional Schechter model can be visualized across the full -magnitude-redshift plane. - -The filled colour scale shows :math:`\log_{10}\Phi(M \mid z)`. The white -contours mark constant :math:`\log_{10}\Phi(M \mid z)` levels at -5, -4, -3, -and -2. These contours make it easier to see where equal-abundance regions sit -as both magnitude and redshift change. +The same model can be shown across the full magnitude-redshift plane. -This plot is useful for checking that the conditional model varies smoothly -across the range where it will be used. +The filled colour scale shows :math:`\log_{10}\Phi(M \mid z)`. The contours +mark constant abundance levels. This is a useful diagnostic for checking that +the conditional model behaves smoothly across the region where it will be used. .. plot:: :include-source: True @@ -116,25 +116,27 @@ across the range where it will be used. import matplotlib.pyplot as plt import cmasher as cmr - from lfkit import LuminosityFunction + from lfkit.photometry.conditional_lf_models import ( + conditional_schechter, + ) LABEL_SIZE = 15 TICK_SIZE = 13 TITLE_SIZE = 17 - LEGEND_SIZE = 15 absolute_mag = np.linspace(-24.0, -16.0, 220) redshift = np.linspace(0.0, 1.5, 180) mag_grid, z_grid = np.meshgrid(absolute_mag, redshift) - lf = LuminosityFunction.conditional_schechter( + phi = conditional_schechter( + mag_grid, + z_grid, phi_star=lambda z: 1.0e-3 * (1.0 + z) ** 0.8, m_star=lambda z: -20.5 - 0.7 * (z - 0.1), alpha=-1.1, ) - phi = lf.phi(mag_grid, z_grid) log_phi = np.log10(phi) fig, ax = plt.subplots(figsize=(7.2, 5.0)) @@ -166,7 +168,8 @@ across the range where it will be used. cbar = fig.colorbar(mesh, ax=ax) cbar.set_label( - r"$\log_{10}\Phi(M \mid z)$ [$\log_{10}(\mathrm{Mpc}^{-3}\,\mathrm{mag}^{-1})$]", + r"$\log_{10}\Phi(M \mid z)$ " + r"[$\log_{10}(\mathrm{Mpc}^{-3}\,\mathrm{mag}^{-1})$]", fontsize=LABEL_SIZE, ) cbar.ax.tick_params(labelsize=TICK_SIZE) @@ -174,19 +177,15 @@ across the range where it will be used. plt.tight_layout() -Conditional evolving Schechter model ------------------------------------- - -A conditional evolving Schechter model uses LFKit's registered parameter models -to define how :math:`\phi_\star`, :math:`M_\star`, and :math:`\alpha` depend on -the conditioning variable. +Conditional Schechter model with LFKit parameter models +------------------------------------------------------- -This example is similar to the standard evolving Schechter interface, but it is -included here because it can be used in the same conditional-LF workflow. The -conditioning variable is passed through ``phi(M, z)``. +LFKit can also evaluate conditional Schechter models using its registered +parameter models. This is useful when the desired dependence follows one of the +standard LFKit parameterizations. -This is useful when the desired parameter evolution already matches one of -LFKit's registered parameter models. +Here, the normalization and characteristic magnitude evolve with the +conditioning variable, while the faint-end slope is constant. .. plot:: :include-source: True @@ -196,32 +195,35 @@ LFKit's registered parameter models. import matplotlib.pyplot as plt import cmasher as cmr - from lfkit import LuminosityFunction + from lfkit.photometry.conditional_lf_models import ( + conditional_schechter_evolving, + ) LABEL_SIZE = 15 TICK_SIZE = 13 TITLE_SIZE = 17 LEGEND_SIZE = 15 - colors_red = cmr.take_cmap_colors("cmr.guppy", 3, cmap_range=(0.03, 0.26)) + colors = cmr.take_cmap_colors("cmr.guppy", 3, cmap_range=(0.03, 0.26)) absolute_mag = np.linspace(-24.0, -14.0, 500) redshifts = [0.1, 0.6, 1.1] - lf = LuminosityFunction.conditional_evolving_schechter( - phi_model="linear_p", - phi_kwargs={"phi_0_star": 1.0e-3, "p": 0.7}, - m_star_model="linear_q", - m_star_kwargs={"m_0_star": -20.5, "q": 0.8, "z_ref": 0.1}, - alpha_model="constant", - alpha_kwargs={"alpha": -1.1}, - ) - fig, ax = plt.subplots(figsize=(7.0, 5.0)) - for z_value, color in zip(redshifts, colors_red): + for z_value, color in zip(redshifts, colors): z = np.full_like(absolute_mag, z_value) - phi = lf.phi(absolute_mag, z) + + phi = conditional_schechter_evolving( + absolute_mag, + z, + phi_model="linear_p", + phi_kwargs={"phi_0_star": 1.0e-3, "p": 0.7}, + m_star_model="linear_q", + m_star_kwargs={"m_0_star": -20.5, "q": 0.8, "z_ref": 0.1}, + alpha_model="constant", + alpha_kwargs={"alpha": -1.1}, + ) ax.plot( absolute_mag, @@ -234,26 +236,25 @@ LFKit's registered parameter models. ax.set_yscale("log") ax.invert_xaxis() ax.set_xlabel("Absolute magnitude $M$", fontsize=LABEL_SIZE) - ax.set_ylabel(r"$\Phi(M \mid z)$ [$\mathrm{Mpc}^{-3}\,\mathrm{mag}^{-1}$]", fontsize=LABEL_SIZE) + ax.set_ylabel( + r"$\Phi(M \mid z)$ [$\mathrm{Mpc}^{-3}\,\mathrm{mag}^{-1}$]", + fontsize=LABEL_SIZE, + ) ax.set_title("Conditional evolving Schechter model", fontsize=TITLE_SIZE) ax.tick_params(axis="both", labelsize=TICK_SIZE) ax.legend(frameon=True, fontsize=LEGEND_SIZE, loc="best") plt.tight_layout() -Central lognormal conditional luminosity function -------------------------------------------------- - -A central-galaxy conditional luminosity function can be represented by a narrow -lognormal component in luminosity, written here in magnitude space. +Lognormal conditional component +------------------------------- -This example shows a central component whose mean absolute magnitude becomes -brighter with redshift. The scatter is kept fixed. The peak therefore shifts in -absolute magnitude while retaining a similar width. +A narrow lognormal component can represent a population concentrated around a +characteristic luminosity at fixed condition. -This type of component is useful in central/satellite decompositions where the -central galaxy population is concentrated around a characteristic luminosity at -fixed condition. +This example uses a mean absolute magnitude that becomes brighter with +redshift. The scatter is fixed, so the peak shifts while retaining a similar +width. .. plot:: :include-source: True @@ -263,29 +264,32 @@ fixed condition. import matplotlib.pyplot as plt import cmasher as cmr - from lfkit import LuminosityFunction + from lfkit.photometry.conditional_lf_models import ( + lognormal_conditional_lf, + ) LABEL_SIZE = 15 TICK_SIZE = 13 TITLE_SIZE = 17 LEGEND_SIZE = 15 - colors_blue = cmr.take_cmap_colors("cmr.guppy", 3, cmap_range=(0.72, 0.96)) + colors = cmr.take_cmap_colors("cmr.guppy", 3, cmap_range=(0.72, 0.96)) absolute_mag = np.linspace(-24.0, -16.0, 500) redshifts = [0.1, 0.6, 1.1] - lf = LuminosityFunction.central_lognormal_conditional( - mean_absolute_mag=lambda z: -20.8 - 0.6 * (z - 0.1), - sigma_log_luminosity=0.18, - amplitude=lambda z: 8.0e-4 * (1.0 + z) ** 0.4, - ) - fig, ax = plt.subplots(figsize=(7.0, 5.0)) - for z_value, color in zip(redshifts, colors_blue): + for z_value, color in zip(redshifts, colors): z = np.full_like(absolute_mag, z_value) - phi = lf.phi(absolute_mag, z) + + phi = lognormal_conditional_lf( + absolute_mag, + z, + mean_absolute_mag=lambda z: -20.8 - 0.6 * (z - 0.1), + sigma_log_luminosity=0.18, + amplitude=lambda z: 8.0e-4 * (1.0 + z) ** 0.4, + ) ax.plot( absolute_mag, @@ -298,27 +302,24 @@ fixed condition. ax.set_yscale("log") ax.invert_xaxis() ax.set_xlabel("Absolute magnitude $M$", fontsize=LABEL_SIZE) - ax.set_ylabel(r"$\Phi_{\rm cen}(M \mid z)$ [$\mathrm{Mpc}^{-3}\,\mathrm{mag}^{-1}$]", fontsize=LABEL_SIZE) - ax.set_title("Central lognormal conditional LF", fontsize=TITLE_SIZE) + ax.set_ylabel( + r"$\Phi_{\rm lognormal}(M \mid z)$ " + r"[$\mathrm{Mpc}^{-3}\,\mathrm{mag}^{-1}$]", + fontsize=LABEL_SIZE, + ) + ax.set_title("Lognormal conditional LF component", fontsize=TITLE_SIZE) ax.tick_params(axis="both", labelsize=TICK_SIZE) ax.legend(frameon=True, fontsize=LEGEND_SIZE, loc="best") plt.tight_layout() -Satellite modified-Schechter conditional luminosity function ------------------------------------------------------------- +Modified Schechter conditional component +---------------------------------------- -A satellite conditional luminosity function can be represented by a modified -Schechter-like component. - -This example uses a satellite model with an exponential cutoff -:math:`\exp[-(L/L_\star)^2]`. Compared with the central lognormal component, -the satellite component is broader and contributes more strongly across a wider +The modified Schechter component uses a squared exponential cutoff in luminosity +ratio. It is broader than the lognormal component and contributes over a wider range of faint magnitudes. -This is useful for modeling satellite populations separately from central -galaxies. - .. plot:: :include-source: True :width: 520 @@ -327,29 +328,32 @@ galaxies. import matplotlib.pyplot as plt import cmasher as cmr - from lfkit import LuminosityFunction + from lfkit.photometry.conditional_lf_models import ( + modified_schechter_conditional_lf, + ) LABEL_SIZE = 15 TICK_SIZE = 13 TITLE_SIZE = 17 LEGEND_SIZE = 15 - colors_red = cmr.take_cmap_colors("cmr.guppy", 3, cmap_range=(0.03, 0.26)) + colors = cmr.take_cmap_colors("cmr.guppy", 3, cmap_range=(0.03, 0.26)) absolute_mag = np.linspace(-24.0, -14.0, 500) redshifts = [0.1, 0.6, 1.1] - lf = LuminosityFunction.satellite_modified_schechter_conditional( - phi_star=lambda z: 1.2e-3 * (1.0 + z) ** 0.5, - m_star=lambda z: -19.9 - 0.5 * (z - 0.1), - alpha=lambda z: -1.05 - 0.10 * z, - ) - fig, ax = plt.subplots(figsize=(7.0, 5.0)) - for z_value, color in zip(redshifts, colors_red): + for z_value, color in zip(redshifts, colors): z = np.full_like(absolute_mag, z_value) - phi = lf.phi(absolute_mag, z) + + phi = modified_schechter_conditional_lf( + absolute_mag, + z, + phi_star=lambda z: 1.2e-3 * (1.0 + z) ** 0.5, + m_star=lambda z: -19.9 - 0.5 * (z - 0.1), + alpha=lambda z: -1.05 - 0.10 * z, + ) ax.plot( absolute_mag, @@ -362,26 +366,123 @@ galaxies. ax.set_yscale("log") ax.invert_xaxis() ax.set_xlabel("Absolute magnitude $M$", fontsize=LABEL_SIZE) - ax.set_ylabel(r"$\Phi_{\rm sat}(M \mid z)$ [$\mathrm{Mpc}^{-3}\,\mathrm{mag}^{-1}$]", fontsize=LABEL_SIZE) - ax.set_title("Satellite modified-Schechter conditional LF", fontsize=TITLE_SIZE) + ax.set_ylabel( + r"$\Phi_{\rm modSch}(M \mid z)$ " + r"[$\mathrm{Mpc}^{-3}\,\mathrm{mag}^{-1}$]", + fontsize=LABEL_SIZE, + ) + ax.set_title("Modified Schechter conditional LF component", fontsize=TITLE_SIZE) ax.tick_params(axis="both", labelsize=TICK_SIZE) ax.legend(frameon=True, fontsize=LEGEND_SIZE, loc="best") plt.tight_layout() -Central and satellite components --------------------------------- +Standard, modified, and lognormal component shapes +-------------------------------------------------- + +It is useful to compare the component shapes at fixed condition. The standard +Schechter form has the usual exponential cutoff in luminosity ratio. The +modified Schechter component uses a squared exponential cutoff, making the +bright-end suppression sharper. The lognormal component is localized around a +mean luminosity and is useful for narrow populations. + +.. plot:: + :include-source: True + :width: 520 + + import numpy as np + import matplotlib.pyplot as plt + import cmasher as cmr + + from lfkit.photometry.conditional_lf_models import ( + conditional_schechter, + lognormal_conditional_lf, + modified_schechter_conditional_lf, + ) -The central and satellite conditional luminosity functions can be combined into -a single central-plus-satellite model. + LABEL_SIZE = 15 + TICK_SIZE = 13 + TITLE_SIZE = 17 + LEGEND_SIZE = 15 -This plot separates the central component, the satellite component, and their -sum at a fixed redshift. The central component is narrow and localized around -the mean central magnitude. The satellite component is broader and dominates -over a wider faint-magnitude range. + colors = cmr.take_cmap_colors("cmr.guppy", 3, cmap_range=(0.08, 0.92)) -This decomposition is useful when checking which part of the galaxy population -drives the total luminosity function. + absolute_mag = np.linspace(-24.0, -14.0, 600) + z_value = 0.6 + z = np.full_like(absolute_mag, z_value) + + phi_schechter = conditional_schechter( + absolute_mag, + z, + phi_star=1.0e-3, + m_star=-20.5, + alpha=-1.1, + ) + + phi_modified = modified_schechter_conditional_lf( + absolute_mag, + z, + phi_star=1.0e-3, + m_star=-20.5, + alpha=-1.1, + ) + + phi_lognormal = lognormal_conditional_lf( + absolute_mag, + z, + mean_absolute_mag=-20.5, + sigma_log_luminosity=0.20, + amplitude=1.0e-3, + ) + + fig, ax = plt.subplots(figsize=(7.0, 5.0)) + + ax.plot( + absolute_mag, + phi_schechter, + lw=3, + color=colors[0], + label="Standard Schechter", + ) + ax.plot( + absolute_mag, + phi_modified, + lw=3, + color=colors[1], + label="Modified Schechter", + ) + ax.plot( + absolute_mag, + phi_lognormal, + lw=3, + color=colors[2], + label="Lognormal", + ) + + ax.set_yscale("log") + ax.invert_xaxis() + ax.set_xlabel("Absolute magnitude $M$", fontsize=LABEL_SIZE) + ax.set_ylabel( + r"$\Phi(M \mid z=0.6)$ " + r"[$\mathrm{Mpc}^{-3}\,\mathrm{mag}^{-1}$]", + fontsize=LABEL_SIZE, + ) + ax.set_title("Conditional LF component shapes", fontsize=TITLE_SIZE) + ax.tick_params(axis="both", labelsize=TICK_SIZE) + ax.legend(frameon=True, fontsize=LEGEND_SIZE, loc="best") + + plt.tight_layout() + + +Two-component conditional luminosity function +--------------------------------------------- + +The lognormal and modified Schechter components can be combined into a +two-component conditional luminosity function. + +This plot separates the lognormal component, the modified Schechter component, +and their sum at a fixed redshift. This is a useful way to check which component +dominates different magnitude ranges. .. plot:: :include-source: True @@ -391,7 +492,11 @@ drives the total luminosity function. import matplotlib.pyplot as plt import cmasher as cmr - from lfkit import LuminosityFunction + from lfkit.photometry.conditional_lf_models import ( + lognormal_conditional_lf, + modified_schechter_conditional_lf, + two_component_conditional_lf, + ) LABEL_SIZE = 15 TICK_SIZE = 13 @@ -400,82 +505,82 @@ drives the total luminosity function. colors_blue = cmr.take_cmap_colors("cmr.guppy", 3, cmap_range=(0.72, 0.96)) colors_red = cmr.take_cmap_colors("cmr.guppy", 3, cmap_range=(0.03, 0.26)) - c_mid = 0.5 * (np.array(colors_blue[1]) + np.array(colors_red[1])) + total_color = 0.5 * (np.array(colors_blue[1]) + np.array(colors_red[1])) absolute_mag = np.linspace(-24.0, -14.0, 500) z_value = 0.6 z = np.full_like(absolute_mag, z_value) - central = LuminosityFunction.central_lognormal_conditional( + lognormal_phi = lognormal_conditional_lf( + absolute_mag, + z, mean_absolute_mag=lambda z: -20.8 - 0.6 * (z - 0.1), sigma_log_luminosity=0.18, amplitude=lambda z: 8.0e-4 * (1.0 + z) ** 0.4, ) - satellite = LuminosityFunction.satellite_modified_schechter_conditional( + modified_phi = modified_schechter_conditional_lf( + absolute_mag, + z, phi_star=lambda z: 1.2e-3 * (1.0 + z) ** 0.5, m_star=lambda z: -19.9 - 0.5 * (z - 0.1), alpha=lambda z: -1.05 - 0.10 * z, ) - total = LuminosityFunction.central_satellite_conditional( - central_mean_absolute_mag=lambda z: -20.8 - 0.6 * (z - 0.1), - central_sigma_log_luminosity=0.18, - central_amplitude=lambda z: 8.0e-4 * (1.0 + z) ** 0.4, - satellite_phi_star=lambda z: 1.2e-3 * (1.0 + z) ** 0.5, - satellite_alpha=lambda z: -1.05 - 0.10 * z, - satellite_m_star=lambda z: -19.9 - 0.5 * (z - 0.1), + total_phi = two_component_conditional_lf( + absolute_mag, + z, + lognormal_mean_absolute_mag=lambda z: -20.8 - 0.6 * (z - 0.1), + lognormal_sigma_log_luminosity=0.18, + lognormal_amplitude=lambda z: 8.0e-4 * (1.0 + z) ** 0.4, + modified_phi_star=lambda z: 1.2e-3 * (1.0 + z) ** 0.5, + modified_m_star=lambda z: -19.9 - 0.5 * (z - 0.1), + modified_alpha=lambda z: -1.05 - 0.10 * z, ) - central_phi = central.phi(absolute_mag, z) - satellite_phi = satellite.phi(absolute_mag, z) - total_phi = total.phi(absolute_mag, z) - fig, ax = plt.subplots(figsize=(7.0, 5.0)) ax.plot( absolute_mag, - central_phi, + lognormal_phi, lw=3, color=colors_blue[1], - label="Central", + label="Lognormal component", ) ax.plot( absolute_mag, - satellite_phi, + modified_phi, lw=3, color=colors_red[1], - label="Satellite", + label="Modified Schechter component", ) ax.plot( absolute_mag, total_phi, lw=3, - color=c_mid, - label="Central + satellite", + color=total_color, + label="Two-component total", ) ax.set_yscale("log") ax.invert_xaxis() ax.set_xlabel("Absolute magnitude $M$", fontsize=LABEL_SIZE) - ax.set_ylabel(r"$\Phi(M \mid z)$ [$\mathrm{Mpc}^{-3}\,\mathrm{mag}^{-1}$]", fontsize=LABEL_SIZE) - ax.set_title(r"Central and satellite components at $z=0.6$", fontsize=TITLE_SIZE) + ax.set_ylabel( + r"$\Phi(M \mid z)$ [$\mathrm{Mpc}^{-3}\,\mathrm{mag}^{-1}$]", + fontsize=LABEL_SIZE, + ) + ax.set_title(r"Two-component conditional LF at $z=0.6$", fontsize=TITLE_SIZE) ax.tick_params(axis="both", labelsize=TICK_SIZE) ax.legend(frameon=True, fontsize=LEGEND_SIZE, loc="best") plt.tight_layout() -Central-plus-satellite evolution --------------------------------- - -The central-plus-satellite model can also be evaluated at several redshifts. - -This plot shows how the total conditional luminosity function changes when the -central and satellite parameters vary with redshift. The model combines a -narrow central component with a broader satellite component at each redshift. +Two-component evolution +----------------------- -This is useful when the total luminosity function is needed, but the model -still keeps a physically interpretable central/satellite split. +The two-component conditional luminosity function can be evaluated across +several redshifts. This example shows how the full model changes when both +components depend on the conditioning variable. .. plot:: :include-source: True @@ -485,32 +590,35 @@ still keeps a physically interpretable central/satellite split. import matplotlib.pyplot as plt import cmasher as cmr - from lfkit import LuminosityFunction + from lfkit.photometry.conditional_lf_models import ( + two_component_conditional_lf, + ) LABEL_SIZE = 15 TICK_SIZE = 13 TITLE_SIZE = 17 LEGEND_SIZE = 15 - colors_blue = cmr.take_cmap_colors("cmr.guppy", 3, cmap_range=(0.72, 0.96)) + colors = cmr.take_cmap_colors("cmr.guppy", 3, cmap_range=(0.72, 0.96)) absolute_mag = np.linspace(-24.0, -14.0, 500) redshifts = [0.1, 0.6, 1.1] - lf = LuminosityFunction.central_satellite_conditional( - central_mean_absolute_mag=lambda z: -20.8 - 0.6 * (z - 0.1), - central_sigma_log_luminosity=0.18, - central_amplitude=lambda z: 8.0e-4 * (1.0 + z) ** 0.4, - satellite_phi_star=lambda z: 1.2e-3 * (1.0 + z) ** 0.5, - satellite_alpha=lambda z: -1.05 - 0.10 * z, - satellite_m_star=lambda z: -19.9 - 0.5 * (z - 0.1), - ) - fig, ax = plt.subplots(figsize=(7.0, 5.0)) - for z_value, color in zip(redshifts, colors_blue): + for z_value, color in zip(redshifts, colors): z = np.full_like(absolute_mag, z_value) - phi = lf.phi(absolute_mag, z) + + phi = two_component_conditional_lf( + absolute_mag, + z, + lognormal_mean_absolute_mag=lambda z: -20.8 - 0.6 * (z - 0.1), + lognormal_sigma_log_luminosity=0.18, + lognormal_amplitude=lambda z: 8.0e-4 * (1.0 + z) ** 0.4, + modified_phi_star=lambda z: 1.2e-3 * (1.0 + z) ** 0.5, + modified_m_star=lambda z: -19.9 - 0.5 * (z - 0.1), + modified_alpha=lambda z: -1.05 - 0.10 * z, + ) ax.plot( absolute_mag, @@ -523,8 +631,12 @@ still keeps a physically interpretable central/satellite split. ax.set_yscale("log") ax.invert_xaxis() ax.set_xlabel("Absolute magnitude $M$", fontsize=LABEL_SIZE) - ax.set_ylabel(r"$\Phi_{\rm cen+sat}(M \mid z)$ [$\mathrm{Mpc}^{-3}\,\mathrm{mag}^{-1}$]", fontsize=LABEL_SIZE) - ax.set_title("Central-plus-satellite conditional LF", fontsize=TITLE_SIZE) + ax.set_ylabel( + r"$\Phi_{\rm total}(M \mid z)$ " + r"[$\mathrm{Mpc}^{-3}\,\mathrm{mag}^{-1}$]", + fontsize=LABEL_SIZE, + ) + ax.set_title("Two-component conditional LF", fontsize=TITLE_SIZE) ax.tick_params(axis="both", labelsize=TICK_SIZE) ax.legend(frameon=True, fontsize=LEGEND_SIZE, loc="best") plt.tight_layout() @@ -536,12 +648,10 @@ Integrated conditional number density A conditional luminosity function can be integrated over absolute magnitude at each value of the conditioning variable. -This example integrates the central component, satellite component, and total -central-plus-satellite luminosity function over a fixed absolute-magnitude -range. The result shows how the selected number density changes with redshift. - -This type of calculation is useful when the conditional luminosity function is -used as an ingredient in redshift-distribution or abundance calculations. +This example uses LFKit's conditional luminosity-function integration helper to +integrate the lognormal component, the modified Schechter component, and the +two-component total over a fixed absolute-magnitude range. The result shows how +the selected number density changes with redshift. .. plot:: :include-source: True @@ -551,7 +661,14 @@ used as an ingredient in redshift-distribution or abundance calculations. import matplotlib.pyplot as plt import cmasher as cmr - from lfkit import LuminosityFunction + from lfkit.photometry.conditional_lf_integrals import ( + integrate_conditional_luminosity_function, + ) + from lfkit.photometry.conditional_lf_models import ( + lognormal_conditional_lf, + modified_schechter_conditional_lf, + two_component_conditional_lf, + ) LABEL_SIZE = 15 TICK_SIZE = 13 @@ -560,56 +677,79 @@ used as an ingredient in redshift-distribution or abundance calculations. colors_blue = cmr.take_cmap_colors("cmr.guppy", 3, cmap_range=(0.72, 0.96)) colors_red = cmr.take_cmap_colors("cmr.guppy", 3, cmap_range=(0.03, 0.26)) - c_mid = 0.5 * (np.array(colors_blue[1]) + np.array(colors_red[1])) + total_color = 0.5 * (np.array(colors_blue[1]) + np.array(colors_red[1])) redshift = np.linspace(0.05, 1.5, 180) + absolute_mag = np.linspace(-24.0, -14.0, 800) - central = LuminosityFunction.central_lognormal_conditional( - mean_absolute_mag=lambda z: -20.8 - 0.6 * (z - 0.1), - sigma_log_luminosity=0.18, - amplitude=lambda z: 8.0e-4 * (1.0 + z) ** 0.4, + _, z_grid = np.meshgrid(absolute_mag, redshift) + + n_lognormal = integrate_conditional_luminosity_function( + absolute_mag=absolute_mag, + condition=z_grid, + conditional_lf=lambda absolute_mag, condition: lognormal_conditional_lf( + absolute_mag, + condition, + mean_absolute_mag=lambda z: -20.8 - 0.6 * (z - 0.1), + sigma_log_luminosity=0.18, + amplitude=lambda z: 8.0e-4 * (1.0 + z) ** 0.4, + ), + axis=1, ) - satellite = LuminosityFunction.satellite_modified_schechter_conditional( - phi_star=lambda z: 1.2e-3 * (1.0 + z) ** 0.5, - m_star=lambda z: -19.9 - 0.5 * (z - 0.1), - alpha=lambda z: -1.05 - 0.10 * z, + n_modified = integrate_conditional_luminosity_function( + absolute_mag=absolute_mag, + condition=z_grid, + conditional_lf=lambda absolute_mag, condition: modified_schechter_conditional_lf( + absolute_mag, + condition, + phi_star=lambda z: 1.2e-3 * (1.0 + z) ** 0.5, + m_star=lambda z: -19.9 - 0.5 * (z - 0.1), + alpha=lambda z: -1.05 - 0.10 * z, + ), + axis=1, ) - total = LuminosityFunction.central_satellite_conditional( - central_mean_absolute_mag=lambda z: -20.8 - 0.6 * (z - 0.1), - central_sigma_log_luminosity=0.18, - central_amplitude=lambda z: 8.0e-4 * (1.0 + z) ** 0.4, - satellite_phi_star=lambda z: 1.2e-3 * (1.0 + z) ** 0.5, - satellite_alpha=lambda z: -1.05 - 0.10 * z, - satellite_m_star=lambda z: -19.9 - 0.5 * (z - 0.1), + n_total = integrate_conditional_luminosity_function( + absolute_mag=absolute_mag, + condition=z_grid, + conditional_lf=lambda absolute_mag, condition: two_component_conditional_lf( + absolute_mag, + condition, + lognormal_mean_absolute_mag=lambda z: -20.8 - 0.6 * (z - 0.1), + lognormal_sigma_log_luminosity=0.18, + lognormal_amplitude=lambda z: 8.0e-4 * (1.0 + z) ** 0.4, + modified_phi_star=lambda z: 1.2e-3 * (1.0 + z) ** 0.5, + modified_m_star=lambda z: -19.9 - 0.5 * (z - 0.1), + modified_alpha=lambda z: -1.05 - 0.10 * z, + ), + axis=1, ) - n_central = central.integrated_number_density( + fig, ax = plt.subplots(figsize=(7.0, 5.0)) + + ax.plot( redshift, - m_bright=-24.0, - m_faint=-14.0, - n_m=800, + n_lognormal, + lw=3, + color=colors_blue[1], + label="Lognormal component", ) - n_satellite = satellite.integrated_number_density( + ax.plot( redshift, - m_bright=-24.0, - m_faint=-14.0, - n_m=800, + n_modified, + lw=3, + color=colors_red[1], + label="Modified Schechter component", ) - n_total = total.integrated_number_density( + ax.plot( redshift, - m_bright=-24.0, - m_faint=-14.0, - n_m=800, + n_total, + lw=3, + color=total_color, + label="Two-component total", ) - fig, ax = plt.subplots(figsize=(7.0, 5.0)) - - ax.plot(redshift, n_central, lw=3, color=colors_blue[1], label="Central") - ax.plot(redshift, n_satellite, lw=3, color=colors_red[1], label="Satellite") - ax.plot(redshift, n_total, lw=3, color=c_mid, label="Central + satellite") - ax.set_yscale("log") ax.set_xlabel("Redshift $z$", fontsize=LABEL_SIZE) ax.set_ylabel(r"Integrated number density [$\mathrm{Mpc}^{-3}$]", fontsize=LABEL_SIZE) @@ -619,18 +759,16 @@ used as an ingredient in redshift-distribution or abundance calculations. plt.tight_layout() -Central fraction ----------------- +Component fractions +------------------- -The relative contribution of central and satellite galaxies can be summarized -as a fraction of the integrated central-plus-satellite luminosity function. +The relative contribution of each component can be summarized as a fraction of +the integrated two-component luminosity function. -This example computes the central fraction over a fixed absolute-magnitude -range. The fraction changes with redshift because the central and satellite -components have different conditional parameter dependence. - -This is a compact diagnostic for checking whether the selected population is -central-dominated, satellite-dominated, or mixed. +This example uses LFKit's conditional luminosity-function integration helper to +compute the integrated lognormal and modified Schechter components. This is a +compact diagnostic for checking whether the selected population is dominated by +the lognormal component, the modified Schechter component, or a mixture of both. .. plot:: :include-source: True @@ -640,7 +778,13 @@ central-dominated, satellite-dominated, or mixed. import matplotlib.pyplot as plt import cmasher as cmr - from lfkit import LuminosityFunction + from lfkit.photometry.conditional_lf_integrals import ( + integrate_conditional_luminosity_function, + ) + from lfkit.photometry.conditional_lf_models import ( + lognormal_conditional_lf, + modified_schechter_conditional_lf, + ) LABEL_SIZE = 15 TICK_SIZE = 13 @@ -651,73 +795,77 @@ central-dominated, satellite-dominated, or mixed. colors_red = cmr.take_cmap_colors("cmr.guppy", 3, cmap_range=(0.03, 0.26)) redshift = np.linspace(0.05, 1.5, 180) + absolute_mag = np.linspace(-24.0, -14.0, 800) - central = LuminosityFunction.central_lognormal_conditional( - mean_absolute_mag=lambda z: -20.8 - 0.6 * (z - 0.1), - sigma_log_luminosity=0.18, - amplitude=lambda z: 8.0e-4 * (1.0 + z) ** 0.4, - ) + _, z_grid = np.meshgrid(absolute_mag, redshift) - satellite = LuminosityFunction.satellite_modified_schechter_conditional( - phi_star=lambda z: 1.2e-3 * (1.0 + z) ** 0.5, - m_star=lambda z: -19.9 - 0.5 * (z - 0.1), - alpha=lambda z: -1.05 - 0.10 * z, + n_lognormal = integrate_conditional_luminosity_function( + absolute_mag=absolute_mag, + condition=z_grid, + conditional_lf=lambda absolute_mag, condition: lognormal_conditional_lf( + absolute_mag, + condition, + mean_absolute_mag=lambda z: -20.8 - 0.6 * (z - 0.1), + sigma_log_luminosity=0.18, + amplitude=lambda z: 8.0e-4 * (1.0 + z) ** 0.4, + ), + axis=1, ) - n_central = central.integrated_number_density( - redshift, - m_bright=-24.0, - m_faint=-14.0, - n_m=800, - ) - n_satellite = satellite.integrated_number_density( - redshift, - m_bright=-24.0, - m_faint=-14.0, - n_m=800, + n_modified = integrate_conditional_luminosity_function( + absolute_mag=absolute_mag, + condition=z_grid, + conditional_lf=lambda absolute_mag, condition: modified_schechter_conditional_lf( + absolute_mag, + condition, + phi_star=lambda z: 1.2e-3 * (1.0 + z) ** 0.5, + m_star=lambda z: -19.9 - 0.5 * (z - 0.1), + alpha=lambda z: -1.05 - 0.10 * z, + ), + axis=1, ) - central_fraction = n_central / (n_central + n_satellite) - satellite_fraction = n_satellite / (n_central + n_satellite) + n_total = n_lognormal + n_modified + + lognormal_fraction = n_lognormal / n_total + modified_fraction = n_modified / n_total fig, ax = plt.subplots(figsize=(7.0, 5.0)) ax.plot( redshift, - central_fraction, + lognormal_fraction, lw=3, color=colors_blue[1], - label="Central fraction", + label="Lognormal fraction", ) ax.plot( redshift, - satellite_fraction, + modified_fraction, lw=3, color=colors_red[1], - label="Satellite fraction", + label="Modified Schechter fraction", ) ax.set_ylim(-0.05, 1.05) ax.set_xlabel("Redshift $z$", fontsize=LABEL_SIZE) ax.set_ylabel("Fraction of integrated LF", fontsize=LABEL_SIZE) - ax.set_title(r"Central and satellite fractions over $-24 \leq M \leq -14$", fontsize=TITLE_SIZE) + ax.set_title(r"Component fractions over $-24 \leq M \leq -14$", fontsize=TITLE_SIZE) ax.tick_params(axis="both", labelsize=TICK_SIZE) ax.legend(frameon=True, fontsize=LEGEND_SIZE, loc="center right") plt.tight_layout() -Central-plus-satellite surface ------------------------------- - -The full central-plus-satellite conditional luminosity function can be shown as -a two-dimensional surface. +Two-component LF surface +----------------------------------------- -The filled colour scale shows :math:`\log_{10}\Phi_{\rm cen+sat}(M \mid z)`. -The white contours mark constant :math:`\log_{10}\Phi_{\rm cen+sat}(M \mid z)` -levels at -5, -4, -3, and -2. +The full two-component conditional luminosity function can be shown as a +surface in the magnitude-redshift plane. -This view is useful for checking whether the central peak, satellite tail, and -redshift dependence combine smoothly over the full magnitude-redshift range. +The filled colour scale shows :math:`\log_{10}\Phi_{\rm total}(M \mid z)`. +The contours mark constant abundance levels. This view is useful for checking +whether the narrow component, broad component, and redshift dependence combine +smoothly. .. plot:: :include-source: True @@ -727,28 +875,30 @@ redshift dependence combine smoothly over the full magnitude-redshift range. import matplotlib.pyplot as plt import cmasher as cmr - from lfkit import LuminosityFunction + from lfkit.photometry.conditional_lf_models import ( + two_component_conditional_lf, + ) LABEL_SIZE = 15 TICK_SIZE = 13 TITLE_SIZE = 17 - LEGEND_SIZE = 15 absolute_mag = np.linspace(-24.0, -14.0, 220) redshift = np.linspace(0.0, 1.5, 180) mag_grid, z_grid = np.meshgrid(absolute_mag, redshift) - lf = LuminosityFunction.central_satellite_conditional( - central_mean_absolute_mag=lambda z: -20.8 - 0.6 * (z - 0.1), - central_sigma_log_luminosity=0.18, - central_amplitude=lambda z: 8.0e-4 * (1.0 + z) ** 0.4, - satellite_phi_star=lambda z: 1.2e-3 * (1.0 + z) ** 0.5, - satellite_alpha=lambda z: -1.05 - 0.10 * z, - satellite_m_star=lambda z: -19.9 - 0.5 * (z - 0.1), + phi = two_component_conditional_lf( + mag_grid, + z_grid, + lognormal_mean_absolute_mag=lambda z: -20.8 - 0.6 * (z - 0.1), + lognormal_sigma_log_luminosity=0.18, + lognormal_amplitude=lambda z: 8.0e-4 * (1.0 + z) ** 0.4, + modified_phi_star=lambda z: 1.2e-3 * (1.0 + z) ** 0.5, + modified_m_star=lambda z: -19.9 - 0.5 * (z - 0.1), + modified_alpha=lambda z: -1.05 - 0.10 * z, ) - phi = lf.phi(mag_grid, z_grid) log_phi = np.log10(phi) fig, ax = plt.subplots(figsize=(7.2, 5.0)) @@ -775,14 +925,333 @@ redshift dependence combine smoothly over the full magnitude-redshift range. ax.invert_xaxis() ax.set_xlabel("Absolute magnitude $M$", fontsize=LABEL_SIZE) ax.set_ylabel("Redshift $z$", fontsize=LABEL_SIZE) - ax.set_title("Central-plus-satellite conditional LF surface", fontsize=TITLE_SIZE) + ax.set_title("Two-component conditional LF surface", fontsize=TITLE_SIZE) ax.tick_params(axis="both", labelsize=TICK_SIZE) cbar = fig.colorbar(mesh, ax=ax) cbar.set_label( - r"$\log_{10}\Phi_{\rm cen+sat}(M \mid z)$ [$\log_{10}(\mathrm{Mpc}^{-3}\,\mathrm{mag}^{-1})$]", + r"$\log_{10}\Phi_{\rm total}(M \mid z)$ " + r"[$\log_{10}(\mathrm{Mpc}^{-3}\,\mathrm{mag}^{-1})$]", fontsize=LABEL_SIZE, ) cbar.ax.tick_params(labelsize=TICK_SIZE) plt.tight_layout() + + +Halo-mass conditional luminosity function +----------------------------------------- + +The conditioning variable does not need to be redshift. In halo-model +applications, a conditional luminosity function is often written as +:math:`\Phi(M \mid M_h)`, where :math:`M_h` is halo mass. + +This example uses log halo mass as the conditioning variable and lets the +lognormal mean magnitude become brighter in more massive halos. + +.. plot:: + :include-source: True + :width: 520 + + import numpy as np + import matplotlib.pyplot as plt + import cmasher as cmr + + from lfkit.photometry.conditional_lf_models import lognormal_conditional_lf + + LABEL_SIZE = 15 + TICK_SIZE = 13 + TITLE_SIZE = 17 + LEGEND_SIZE = 15 + + colors = cmr.take_cmap_colors("cmr.guppy", 4, cmap_range=(0.08, 0.92)) + + absolute_mag = np.linspace(-24.0, -16.0, 600) + log_halo_masses = [11.5, 12.0, 12.5, 13.0] + + fig, ax = plt.subplots(figsize=(7.0, 5.0)) + + for log_mh, color in zip(log_halo_masses, colors): + condition = np.full_like(absolute_mag, log_mh) + + phi = lognormal_conditional_lf( + absolute_mag, + condition, + mean_absolute_mag=lambda log_mh: -20.0 - 0.8 * (log_mh - 12.0), + sigma_log_luminosity=0.18, + amplitude=lambda log_mh: 5.0e-4 * 10.0 ** (0.3 * (log_mh - 12.0)), + ) + + ax.plot( + absolute_mag, + phi, + lw=3, + color=color, + label=rf"$\log_{{10}} M_h={log_mh}$", + ) + + ax.set_yscale("log") + ax.invert_xaxis() + ax.set_xlabel("Absolute magnitude $M$", fontsize=LABEL_SIZE) + ax.set_ylabel( + r"$\Phi(M \mid M_h)$ [$\mathrm{Mpc}^{-3}\,\mathrm{mag}^{-1}$]", + fontsize=LABEL_SIZE, + ) + ax.set_title("Halo-mass conditional lognormal LF", fontsize=TITLE_SIZE) + ax.tick_params(axis="both", labelsize=TICK_SIZE) + ax.legend(frameon=True, fontsize=LEGEND_SIZE, loc="best") + + plt.tight_layout() + + +Mean magnitude from a conditional luminosity function +----------------------------------------------------- + +Weighted integrals can be used to compute summary statistics of a conditional +luminosity function. For example, the mean absolute magnitude at fixed condition +is + +:math:`\langle M \rangle(x) = \int M \Phi(M \mid x)\,dM / \int \Phi(M \mid x)\,dM`. + +.. plot:: + :include-source: True + :width: 520 + + import numpy as np + import matplotlib.pyplot as plt + import cmasher as cmr + + from lfkit.photometry.conditional_lf_integrals import ( + integrate_conditional_luminosity_function, + integrate_weighted_conditional_luminosity_function, + ) + from lfkit.photometry.conditional_lf_models import two_component_conditional_lf + + LABEL_SIZE = 15 + TICK_SIZE = 13 + TITLE_SIZE = 17 + + redshift = np.linspace(0.05, 1.5, 180) + absolute_mag = np.linspace(-24.0, -14.0, 800) + + mag_grid, z_grid = np.meshgrid(absolute_mag, redshift) + + number_density = integrate_conditional_luminosity_function( + absolute_mag=absolute_mag, + condition=z_grid, + conditional_lf=lambda absolute_mag, condition: two_component_conditional_lf( + absolute_mag, + condition, + lognormal_mean_absolute_mag=lambda z: -20.8 - 0.6 * (z - 0.1), + lognormal_sigma_log_luminosity=0.18, + lognormal_amplitude=lambda z: 8.0e-4 * (1.0 + z) ** 0.4, + modified_phi_star=lambda z: 1.2e-3 * (1.0 + z) ** 0.5, + modified_m_star=lambda z: -19.9 - 0.5 * (z - 0.1), + modified_alpha=lambda z: -1.05 - 0.10 * z, + ), + axis=1, + ) + + weighted_magnitude = integrate_weighted_conditional_luminosity_function( + absolute_mag=absolute_mag, + condition=z_grid, + conditional_lf=lambda absolute_mag, condition: two_component_conditional_lf( + absolute_mag, + condition, + lognormal_mean_absolute_mag=lambda z: -20.8 - 0.6 * (z - 0.1), + lognormal_sigma_log_luminosity=0.18, + lognormal_amplitude=lambda z: 8.0e-4 * (1.0 + z) ** 0.4, + modified_phi_star=lambda z: 1.2e-3 * (1.0 + z) ** 0.5, + modified_m_star=lambda z: -19.9 - 0.5 * (z - 0.1), + modified_alpha=lambda z: -1.05 - 0.10 * z, + ), + weight=lambda absolute_mag, condition: absolute_mag, + axis=1, + ) + + mean_magnitude = weighted_magnitude / number_density + + fig, ax = plt.subplots(figsize=(7.0, 5.0)) + + ax.plot( + redshift, + mean_magnitude, + lw=3, + color=cmr.take_cmap_colors("cmr.guppy", 1, cmap_range=(0.7, 0.9))[0], + ) + + ax.invert_yaxis() + ax.set_xlabel("Redshift $z$", fontsize=LABEL_SIZE) + ax.set_ylabel(r"Mean absolute magnitude $\langle M \rangle$", fontsize=LABEL_SIZE) + ax.set_title("Mean magnitude of the conditional LF", fontsize=TITLE_SIZE) + ax.tick_params(axis="both", labelsize=TICK_SIZE) + + plt.tight_layout() + + +Selection-limited conditional number density +-------------------------------------------- + +Instead of integrating over a fixed absolute-magnitude range by hand, LFKit can +integrate a luminosity function callable over finite magnitude bounds. This is +useful for survey-like selections where only galaxies brighter than a limiting +absolute magnitude contribute to the selected sample. + +Here, the limiting absolute magnitude becomes brighter with redshift. The +example compares the full number density over a fixed magnitude range with the +number density brighter than the redshift-dependent limit. + +.. plot:: + :include-source: True + :width: 520 + + import numpy as np + import matplotlib.pyplot as plt + import cmasher as cmr + + from lfkit.photometry.conditional_lf_models import two_component_conditional_lf + from lfkit.photometry.lf_integrals import integrated_number_density + + LABEL_SIZE = 15 + TICK_SIZE = 13 + TITLE_SIZE = 17 + LEGEND_SIZE = 15 + + redshift = np.linspace(0.05, 1.5, 180) + + def lf(absolute_mag, z): + return two_component_conditional_lf( + absolute_mag, + z, + lognormal_mean_absolute_mag=lambda z: -20.8 - 0.6 * (z - 0.1), + lognormal_sigma_log_luminosity=0.18, + lognormal_amplitude=lambda z: 8.0e-4 * (1.0 + z) ** 0.4, + modified_phi_star=lambda z: 1.2e-3 * (1.0 + z) ** 0.5, + modified_m_star=lambda z: -19.9 - 0.5 * (z - 0.1), + modified_alpha=lambda z: -1.05 - 0.10 * z, + ) + + limiting_mag = -18.5 - 1.2 * redshift + + n_total = integrated_number_density( + redshift, + lf, + m_bright=-24.0, + m_faint=-14.0, + n_m=800, + ) + + n_selected = integrated_number_density( + redshift, + lf, + m_bright=-24.0, + m_faint=limiting_mag, + n_m=800, + ) + + colors = cmr.take_cmap_colors("cmr.guppy", 3, cmap_range=(0.12, 0.9)) + + fig, ax = plt.subplots(figsize=(7.0, 5.0)) + + ax.plot( + redshift, + n_total, + lw=3, + color=colors[0], + label="Full magnitude range", + ) + ax.plot( + redshift, + n_selected, + lw=3, + color=colors[2], + label="Brighter than limit", + ) + + ax.set_yscale("log") + ax.set_xlabel("Redshift $z$", fontsize=LABEL_SIZE) + ax.set_ylabel( + r"Integrated number density [$\mathrm{Mpc}^{-3}$]", + fontsize=LABEL_SIZE, + ) + ax.set_title("Selection-limited conditional number density", fontsize=TITLE_SIZE) + ax.tick_params(axis="both", labelsize=TICK_SIZE) + ax.legend(frameon=True, fontsize=LEGEND_SIZE, loc="best") + + plt.tight_layout() + + +Selection fraction +------------------ + +The selected fraction is the ratio between the number density brighter than the +redshift-dependent limiting magnitude and the number density over the full +reference magnitude range. + +.. plot:: + :include-source: True + :width: 520 + + import numpy as np + import matplotlib.pyplot as plt + import cmasher as cmr + + from lfkit.photometry.conditional_lf_models import two_component_conditional_lf + from lfkit.photometry.lf_integrals import integrated_number_density + + LABEL_SIZE = 15 + TICK_SIZE = 13 + TITLE_SIZE = 17 + + redshift = np.linspace(0.05, 1.5, 180) + + def lf(absolute_mag, z): + return two_component_conditional_lf( + absolute_mag, + z, + lognormal_mean_absolute_mag=lambda z: -20.8 - 0.6 * (z - 0.1), + lognormal_sigma_log_luminosity=0.18, + lognormal_amplitude=lambda z: 8.0e-4 * (1.0 + z) ** 0.4, + modified_phi_star=lambda z: 1.2e-3 * (1.0 + z) ** 0.5, + modified_m_star=lambda z: -19.9 - 0.5 * (z - 0.1), + modified_alpha=lambda z: -1.05 - 0.10 * z, + ) + + limiting_mag = -18.5 - 1.2 * redshift + + n_total = integrated_number_density( + redshift, + lf, + m_bright=-24.0, + m_faint=-14.0, + n_m=800, + ) + + n_selected = integrated_number_density( + redshift, + lf, + m_bright=-24.0, + m_faint=limiting_mag, + n_m=800, + ) + + selected_fraction = n_selected / n_total + + color = cmr.take_cmap_colors("cmr.guppy", 1, cmap_range=(0.72, 0.9))[0] + + fig, ax = plt.subplots(figsize=(7.0, 5.0)) + + ax.plot( + redshift, + selected_fraction, + lw=3, + color=color, + ) + + ax.set_ylim(-0.05, 1.05) + ax.set_xlabel("Redshift $z$", fontsize=LABEL_SIZE) + ax.set_ylabel("Selected fraction", fontsize=LABEL_SIZE) + ax.set_title("Fraction brighter than the limiting magnitude", fontsize=TITLE_SIZE) + ax.tick_params(axis="both", labelsize=TICK_SIZE) + + plt.tight_layout() diff --git a/docs/examples/index.rst b/docs/examples/index.rst index 02ce404..01ac140 100644 --- a/docs/examples/index.rst +++ b/docs/examples/index.rst @@ -54,7 +54,7 @@ Working, executable examples for using LFKit. :link-type: doc :shadow: md - **Conditional luminosity-function examples** + **Conditional luminosity function examples** ^^^ Build conditional luminosity functions, including redshift-dependent Schechter models and central/satellite components. diff --git a/docs/index.rst b/docs/index.rst index c361bee..d08b85a 100644 --- a/docs/index.rst +++ b/docs/index.rst @@ -12,10 +12,10 @@ It provides a clean interface for building theoretical luminosity functions, including redshift-dependent parameter models, and for connecting apparent magnitude limits to observable and missing galaxy number densities. -LFKit is designed to be science-use-case agnostic: the same luminosity-function +LFKit is designed to be science-use-case agnostic: the same luminosity function machinery can be used in photometric-redshift modelling, intrinsic-alignment modelling, cluster science, GW-cosmology catalog completeness, or any other -analysis that needs luminosity-function-based number densities. +analysis that needs luminosity function-based number densities. Getting started --------------- diff --git a/tests/benchmarks/test_cacciato_hod_reference.py b/tests/benchmarks/test_cacciato_hod_reference.py new file mode 100644 index 0000000..e69de29 From c8992093c48b34279f41e2562f08ed88d339784e Mon Sep 17 00:00:00 2001 From: "niko.sarcevic" Date: Wed, 13 May 2026 23:15:27 -0400 Subject: [PATCH 7/9] updated tests to match new function names --- .../benchmarks/test_cacciato_clf_reference.py | 34 +-- .../benchmarks/test_cacciato_hod_reference.py | 245 ++++++++++++++++++ tests/test_api_lumfunc.py | 6 +- ...est_photometry_conditional_lf_integrals.py | 56 ++-- .../test_photometry_conditional_lf_models.py | 206 +++++++-------- tests/test_photometry_lf_integrals.py | 8 +- 6 files changed, 402 insertions(+), 153 deletions(-) diff --git a/tests/benchmarks/test_cacciato_clf_reference.py b/tests/benchmarks/test_cacciato_clf_reference.py index be34087..d527ab9 100644 --- a/tests/benchmarks/test_cacciato_clf_reference.py +++ b/tests/benchmarks/test_cacciato_clf_reference.py @@ -1,8 +1,8 @@ """Reference benchmarks for Cacciato-style conditional luminosity functions. -These tests compare LFKit's Cacciato-style CLF helpers against explicit -reference formulae from the Cacciato et al. conditional luminosity-function -parametrization. +These tests compare LFKit's generic conditional luminosity function helpers +against explicit reference formulae from the Cacciato et al. conditional +luminosity function parametrization. They are intentionally written as reference/benchmark tests rather than normal unit tests. They check that LFKit preserves the expected central lognormal, @@ -13,9 +13,9 @@ import pytest from lfkit.photometry.conditional_lf_models import ( - central_lognormal_conditional_lf, - central_satellite_conditional_lf, - satellite_modified_schechter_conditional_lf, + lognormal_conditional_lf, + modified_schechter_conditional_lf, + two_component_conditional_lf, ) from lfkit.photometry.luminosities import magnitude_difference_from_luminosity_ratio @@ -104,7 +104,7 @@ def test_cacciato_central_lognormal_matches_reference_formula() -> None: gamma2=0.245, ) - result = central_lognormal_conditional_lf( + result = lognormal_conditional_lf( absolute_mag=absolute_mag, condition=log_halo_mass, mean_absolute_mag=mean_absolute_mag, @@ -144,7 +144,7 @@ def test_cacciato_satellite_modified_schechter_matches_reference_formula() -> No b2=-0.217, ) - result = satellite_modified_schechter_conditional_lf( + result = modified_schechter_conditional_lf( absolute_mag=absolute_mag, condition=log_halo_mass, phi_star=phi_star, @@ -182,16 +182,16 @@ def test_cacciato_central_satellite_matches_sum_of_components() -> None: b2=-0.217, ) - result = central_satellite_conditional_lf( + result = two_component_conditional_lf( absolute_mag=absolute_mag, condition=log_halo_mass, - central_mean_absolute_mag=central_mag, - central_sigma_log_luminosity=0.157, - satellite_phi_star=phi_star, - satellite_alpha=-1.18, - central_amplitude=1.0, - satellite_m_star=None, - satellite_luminosity_fraction=0.562, + lognormal_mean_absolute_mag=central_mag, + lognormal_sigma_log_luminosity=0.157, + modified_phi_star=phi_star, + modified_alpha=-1.18, + lognormal_amplitude=1.0, + modified_m_star=None, + modified_luminosity_fraction=0.562, ) satellite_m_star = central_mag + magnitude_difference_from_luminosity_ratio(0.562) @@ -241,4 +241,4 @@ def test_cacciato_satellite_phi_star_peaks_near_reference_mass_range() -> None: ) assert np.all(np.isfinite(phi_star)) - assert np.all(phi_star > 0.0) \ No newline at end of file + assert np.all(phi_star > 0.0) diff --git a/tests/benchmarks/test_cacciato_hod_reference.py b/tests/benchmarks/test_cacciato_hod_reference.py index e69de29..dc5df65 100644 --- a/tests/benchmarks/test_cacciato_hod_reference.py +++ b/tests/benchmarks/test_cacciato_hod_reference.py @@ -0,0 +1,245 @@ +"""Reference benchmark for Cacciato-style HOD occupation integrals. + +This benchmark compares the original luminosity-space Cacciato conditional +luminosity functions against LFKit's magnitude-space conditional luminosity +function helpers. The goal is to check that integrating the central and +satellite CLFs over a luminosity window gives the same halo occupation numbers +when expressed through LFKit. +""" + +from __future__ import annotations + +import numpy as np +import pytest + +from lfkit.photometry.conditional_lf_models import ( + lognormal_conditional_lf, + modified_schechter_conditional_lf, + two_component_conditional_lf, +) +from lfkit.photometry.luminosities import ( + magnitude_difference_from_luminosity_ratio, +) + + +pytestmark = pytest.mark.benchmark + + +LOGE = 0.4342944819 + + +def cacciato_lc( + halo_mass: np.ndarray, + *, + log_l0: float = 9.95, + log_m1: float = 11.24, + gamma1: float = 3.18, + gamma2: float = 0.245, +) -> np.ndarray: + """Return the Cacciato central luminosity relation.""" + + l0 = 10.0**log_l0 + m1 = 10.0**log_m1 + + return l0 * (halo_mass / m1) ** gamma1 / ( + 1.0 + (halo_mass / m1) ** (gamma1 - gamma2) + ) + + +def cacciato_phi_star_satellite( + halo_mass: np.ndarray, + *, + b0: float = -1.17, + b1: float = 1.53 * 0.739, + b2: float = -0.217 * 0.739**2, +) -> np.ndarray: + """Return the Cacciato satellite normalization.""" + + log_m12 = np.log10(halo_mass / 1.0e12) + + return 10.0 ** (b0 + b1 * log_m12 + b2 * log_m12**2.0) + + +def cacciato_lstar_satellite( + halo_mass: np.ndarray, + *, + luminosity_fraction: float = 0.562, +) -> np.ndarray: + """Return the Cacciato satellite characteristic luminosity.""" + + return luminosity_fraction * cacciato_lc(halo_mass) + + +def magnitude_from_luminosity_ratio( + luminosity: np.ndarray, +) -> np.ndarray: + """Return magnitude for luminosity relative to unit reference luminosity.""" + + return magnitude_difference_from_luminosity_ratio(luminosity) + + +def central_luminosity_clf_reference( + luminosity: np.ndarray, + halo_mass: np.ndarray, + *, + sigma_c: float = 0.157, +) -> np.ndarray: + """Return the Cacciato central CLF in luminosity units.""" + + luminosity_grid = np.atleast_1d(luminosity).reshape(1, -1) + halo_mass_grid = np.atleast_1d(halo_mass).reshape(-1, 1) + + lc = cacciato_lc(halo_mass_grid) + + return ( + LOGE + / (np.sqrt(2.0 * np.pi) * sigma_c) + * np.exp( + -( + np.log10(luminosity_grid) + - np.log10(lc) + ) + ** 2.0 + / (2.0 * sigma_c**2.0) + ) + / luminosity_grid + ) + + +def satellite_luminosity_clf_reference( + luminosity: np.ndarray, + halo_mass: np.ndarray, + *, + alpha_s: float = -1.18, +) -> np.ndarray: + """Return the Cacciato satellite CLF in luminosity units.""" + + luminosity_grid = np.atleast_1d(luminosity).reshape(1, -1) + halo_mass_grid = np.atleast_1d(halo_mass).reshape(-1, 1) + + l_star = cacciato_lstar_satellite(halo_mass_grid) + phi_star = cacciato_phi_star_satellite(halo_mass_grid) + + return ( + phi_star + / luminosity_grid + * (luminosity_grid / l_star) ** (alpha_s + 1.0) + * np.exp(-((luminosity_grid / l_star) ** 2.0)) + ) + + +def test_lfkit_cacciato_central_occupation_matches_luminosity_reference() -> None: + """Compare central occupation integrals from luminosity and magnitude CLFs.""" + + halo_mass = np.geomspace(1.0e8, 1.0e16, 128) + luminosity = np.geomspace(1.5e7, 5.6e9, 2048) + + reference_clf = central_luminosity_clf_reference( + luminosity=luminosity, + halo_mass=halo_mass, + ) + reference_nc = np.trapezoid(reference_clf, x=luminosity, axis=1) + + magnitude = magnitude_from_luminosity_ratio(luminosity) + magnitude_grid = magnitude[::-1] + + central_mag = magnitude_from_luminosity_ratio(cacciato_lc(halo_mass)) + + lfkit_clf = lognormal_conditional_lf( + absolute_mag=magnitude_grid.reshape(1, -1), + condition=np.log10(halo_mass).reshape(-1, 1), + mean_absolute_mag=central_mag.reshape(-1, 1), + sigma_log_luminosity=0.157, + amplitude=1.0, + ) + lfkit_nc = np.trapezoid(lfkit_clf, x=magnitude_grid, axis=1) + + np.testing.assert_allclose( + lfkit_nc, + reference_nc, + rtol=5.0e-4, + atol=1.0e-8, + ) + + +def test_lfkit_cacciato_satellite_occupation_matches_luminosity_reference() -> None: + """Compare satellite occupation integrals from luminosity and magnitude CLFs.""" + + halo_mass = np.geomspace(1.0e8, 1.0e16, 128) + luminosity = np.geomspace(1.5e7, 5.6e9, 2048) + + reference_clf = satellite_luminosity_clf_reference( + luminosity=luminosity, + halo_mass=halo_mass, + ) + reference_ns = np.trapezoid(reference_clf, x=luminosity, axis=1) + + magnitude = magnitude_from_luminosity_ratio(luminosity) + magnitude_grid = magnitude[::-1] + + central_mag = magnitude_from_luminosity_ratio(cacciato_lc(halo_mass)) + satellite_m_star = central_mag + magnitude_difference_from_luminosity_ratio(0.562) + phi_star = cacciato_phi_star_satellite(halo_mass) + + lfkit_clf = modified_schechter_conditional_lf( + absolute_mag=magnitude_grid.reshape(1, -1), + condition=np.log10(halo_mass).reshape(-1, 1), + phi_star=phi_star.reshape(-1, 1), + m_star=satellite_m_star.reshape(-1, 1), + alpha=-1.18, + ) + lfkit_ns = np.trapezoid(lfkit_clf, x=magnitude_grid, axis=1) + + np.testing.assert_allclose( + lfkit_ns, + reference_ns, + rtol=5.0e-4, + atol=1.0e-8, + ) + + +def test_lfkit_cacciato_total_occupation_matches_luminosity_reference() -> None: + """Compare total occupation integrals from luminosity and magnitude CLFs.""" + + halo_mass = np.geomspace(1.0e8, 1.0e16, 128) + luminosity = np.geomspace(1.5e7, 5.6e9, 2048) + + reference_central = central_luminosity_clf_reference( + luminosity=luminosity, + halo_mass=halo_mass, + ) + reference_satellite = satellite_luminosity_clf_reference( + luminosity=luminosity, + halo_mass=halo_mass, + ) + reference_total = np.trapezoid( + reference_central + reference_satellite, + x=luminosity, + axis=1, + ) + + magnitude = magnitude_from_luminosity_ratio(luminosity) + magnitude_grid = magnitude[::-1] + + central_mag = magnitude_from_luminosity_ratio(cacciato_lc(halo_mass)) + phi_star = cacciato_phi_star_satellite(halo_mass) + + lfkit_clf = two_component_conditional_lf( + absolute_mag=magnitude_grid.reshape(1, -1), + condition=np.log10(halo_mass).reshape(-1, 1), + lognormal_mean_absolute_mag=central_mag.reshape(-1, 1), + lognormal_sigma_log_luminosity=0.157, + modified_phi_star=phi_star.reshape(-1, 1), + modified_alpha=-1.18, + lognormal_amplitude=1.0, + modified_m_star=None, + modified_luminosity_fraction=0.562, + ) + lfkit_total = np.trapezoid(lfkit_clf, x=magnitude_grid, axis=1) + + np.testing.assert_allclose( + lfkit_total, + reference_total, + rtol=5.0e-4, + atol=1.0e-8, + ) diff --git a/tests/test_api_lumfunc.py b/tests/test_api_lumfunc.py index 36e8f4a..0abcd37 100644 --- a/tests/test_api_lumfunc.py +++ b/tests/test_api_lumfunc.py @@ -170,10 +170,10 @@ def fake_schechter_double( def test_phi_raises_for_unsupported_model(): - """Tests that unsupported luminosity-function models fail clearly.""" + """Tests that unsupported luminosity function models fail clearly.""" lf = LuminosityFunction(model="bad_model", parameters={}) - with pytest.raises(ValueError, match="Unsupported luminosity-function model"): + with pytest.raises(ValueError, match="Unsupported luminosity function model"): lf.phi(np.array([-21.0, -20.0])) @@ -855,7 +855,7 @@ def test_phi_from_m_raises_for_unsupported_model(): """Tests that phi_from_m rejects unsupported LF models.""" lf = LuminosityFunction(model="bad_model", parameters={}) - with pytest.raises(ValueError, match="Unsupported luminosity-function model"): + with pytest.raises(ValueError, match="Unsupported luminosity function model"): lf.phi_from_m(object(), np.array([0.1]), np.array([22.0])) diff --git a/tests/test_photometry_conditional_lf_integrals.py b/tests/test_photometry_conditional_lf_integrals.py index e2187b2..e6f36c2 100644 --- a/tests/test_photometry_conditional_lf_integrals.py +++ b/tests/test_photometry_conditional_lf_integrals.py @@ -4,19 +4,19 @@ import pytest from lfkit.photometry.conditional_lf_integrals import ( - evaluate_conditional_lf, - integrate_conditional_lf, - integrate_weighted_conditional_lf, + evaluate_conditional_luminosity_function, + integrate_conditional_luminosity_function, + integrate_weighted_conditional_luminosity_function, ) -def test_evaluate_conditional_lf_accepts_scalar_inputs() -> None: +def test_evaluate_conditional_luminosity_function_accepts_scalar_inputs() -> None: """Tests that scalar inputs are evaluated as float arrays.""" def conditional_lf(absolute_mag, condition): return condition * np.exp(-0.1 * absolute_mag) - result = evaluate_conditional_lf( + result = evaluate_conditional_luminosity_function( absolute_mag=-20.0, condition=2.0, conditional_lf=conditional_lf, @@ -30,7 +30,7 @@ def conditional_lf(absolute_mag, condition): assert result == pytest.approx(expected) -def test_evaluate_conditional_lf_accepts_array_inputs() -> None: +def test_evaluate_conditional_luminosity_function_accepts_array_inputs() -> None: """Tests that array inputs are evaluated element-wise.""" absolute_mag = np.array([-22.0, -21.0, -20.0]) @@ -39,7 +39,7 @@ def test_evaluate_conditional_lf_accepts_array_inputs() -> None: def conditional_lf(absolute_mag, condition): return condition * (absolute_mag + 23.0) - result = evaluate_conditional_lf( + result = evaluate_conditional_luminosity_function( absolute_mag=absolute_mag, condition=condition, conditional_lf=conditional_lf, @@ -51,7 +51,7 @@ def conditional_lf(absolute_mag, condition): assert result.dtype == np.float64 -def test_evaluate_conditional_lf_accepts_broadcastable_inputs() -> None: +def test_evaluate_conditional_luminosity_function_accepts_broadcastable_inputs() -> None: """Tests that broadcastable magnitude and condition arrays are supported.""" absolute_mag = np.array([-22.0, -21.0, -20.0]) @@ -60,7 +60,7 @@ def test_evaluate_conditional_lf_accepts_broadcastable_inputs() -> None: def conditional_lf(absolute_mag, condition): return condition * (absolute_mag + 23.0) - result = evaluate_conditional_lf( + result = evaluate_conditional_luminosity_function( absolute_mag=absolute_mag, condition=condition, conditional_lf=conditional_lf, @@ -77,35 +77,35 @@ def conditional_lf(absolute_mag, condition): assert result.dtype == np.float64 -def test_evaluate_conditional_lf_rejects_non_finite_absolute_mag() -> None: +def test_evaluate_conditional_luminosity_function_rejects_non_finite_absolute_mag() -> None: """Tests that non-finite absolute magnitudes are rejected.""" def conditional_lf(absolute_mag, condition): return np.ones_like(absolute_mag) with pytest.raises(ValueError, match="absolute_mag contains NaN or infinite values."): - evaluate_conditional_lf( + evaluate_conditional_luminosity_function( absolute_mag=[-22.0, np.nan, -20.0], condition=1.0, conditional_lf=conditional_lf, ) -def test_evaluate_conditional_lf_rejects_non_finite_condition() -> None: +def test_evaluate_conditional_luminosity_function_rejects_non_finite_condition() -> None: """Tests that non-finite condition values are rejected.""" def conditional_lf(absolute_mag, condition): return np.ones_like(absolute_mag) with pytest.raises(ValueError, match="condition contains NaN or infinite values."): - evaluate_conditional_lf( + evaluate_conditional_luminosity_function( absolute_mag=[-22.0, -21.0, -20.0], condition=np.inf, conditional_lf=conditional_lf, ) -def test_evaluate_conditional_lf_rejects_non_finite_result() -> None: +def test_evaluate_conditional_luminosity_function_rejects_non_finite_result() -> None: """Tests that non-finite conditional luminosity values are rejected.""" def conditional_lf(absolute_mag, condition): @@ -115,14 +115,14 @@ def conditional_lf(absolute_mag, condition): ValueError, match="conditional_lf returned NaN or infinite values.", ): - evaluate_conditional_lf( + evaluate_conditional_luminosity_function( absolute_mag=[-22.0, -21.0, -20.0], condition=1.0, conditional_lf=conditional_lf, ) -def test_evaluate_conditional_lf_rejects_negative_result() -> None: +def test_evaluate_conditional_luminosity_function_rejects_negative_result() -> None: """Tests that negative conditional luminosity values are rejected.""" def conditional_lf(absolute_mag, condition): @@ -132,7 +132,7 @@ def conditional_lf(absolute_mag, condition): ValueError, match="conditional_lf returned negative values, which are not allowed.", ): - evaluate_conditional_lf( + evaluate_conditional_luminosity_function( absolute_mag=[-22.0, -21.0, -20.0], condition=1.0, conditional_lf=conditional_lf, @@ -140,7 +140,7 @@ def conditional_lf(absolute_mag, condition): def test_integrate_conditional_lf_integrates_over_last_axis() -> None: - """Tests conditional luminosity-function integration over magnitude.""" + """Tests conditional luminosity function integration over magnitude.""" absolute_mag = np.array([-22.0, -21.0, -20.0]) condition = np.array([[1.0], [2.0]]) @@ -148,7 +148,7 @@ def test_integrate_conditional_lf_integrates_over_last_axis() -> None: def conditional_lf(absolute_mag, condition): return condition * (absolute_mag + 23.0) - result = integrate_conditional_lf( + result = integrate_conditional_luminosity_function( absolute_mag=absolute_mag, condition=condition, conditional_lf=conditional_lf, @@ -170,7 +170,7 @@ def conditional_lf(absolute_mag, condition): def test_integrate_conditional_lf_integrates_over_requested_axis() -> None: - """Tests conditional luminosity-function integration over a custom axis.""" + """Tests conditional luminosity function integration over a custom axis.""" absolute_mag = np.array([-22.0, -21.0, -20.0]) condition = np.array([[1.0, 2.0]]) @@ -178,7 +178,7 @@ def test_integrate_conditional_lf_integrates_over_requested_axis() -> None: def conditional_lf(absolute_mag, condition): return (absolute_mag[:, None] + 23.0) * condition - result = integrate_conditional_lf( + result = integrate_conditional_luminosity_function( absolute_mag=absolute_mag, condition=condition, conditional_lf=conditional_lf, @@ -208,7 +208,7 @@ def conditional_lf(absolute_mag, condition): ValueError, match="conditional_lf returned negative values, which are not allowed.", ): - integrate_conditional_lf( + integrate_conditional_luminosity_function( absolute_mag=[-22.0, -21.0, -20.0], condition=1.0, conditional_lf=conditional_lf, @@ -216,7 +216,7 @@ def conditional_lf(absolute_mag, condition): def test_integrate_weighted_conditional_lf_integrates_weighted_values() -> None: - """Tests weighted conditional luminosity-function integration.""" + """Tests weighted conditional luminosity function integration.""" absolute_mag = np.array([-22.0, -21.0, -20.0]) condition = np.array([[1.0], [2.0]]) @@ -227,7 +227,7 @@ def conditional_lf(absolute_mag, condition): def weight(absolute_mag, condition): return absolute_mag + 24.0 - result = integrate_weighted_conditional_lf( + result = integrate_weighted_conditional_luminosity_function( absolute_mag=absolute_mag, condition=condition, conditional_lf=conditional_lf, @@ -259,7 +259,7 @@ def conditional_lf(absolute_mag, condition): def weight(absolute_mag, condition): return absolute_mag[:, None] + 24.0 - result = integrate_weighted_conditional_lf( + result = integrate_weighted_conditional_luminosity_function( absolute_mag=absolute_mag, condition=condition, conditional_lf=conditional_lf, @@ -297,7 +297,7 @@ def weight(absolute_mag, condition): return np.array([1.0, np.inf, 3.0]) with pytest.raises(ValueError, match="weight returned NaN or infinite values."): - integrate_weighted_conditional_lf( + integrate_weighted_conditional_luminosity_function( absolute_mag=[-22.0, -21.0, -20.0], condition=1.0, conditional_lf=conditional_lf, @@ -316,7 +316,7 @@ def conditional_lf(absolute_mag, condition): def weight(absolute_mag, condition): return np.array([1.0, -2.0, 3.0]) - result = integrate_weighted_conditional_lf( + result = integrate_weighted_conditional_luminosity_function( absolute_mag=absolute_mag, condition=1.0, conditional_lf=conditional_lf, @@ -342,7 +342,7 @@ def weight(absolute_mag, condition): ValueError, match="conditional_lf returned negative values, which are not allowed.", ): - integrate_weighted_conditional_lf( + integrate_weighted_conditional_luminosity_function( absolute_mag=[-22.0, -21.0, -20.0], condition=1.0, conditional_lf=conditional_lf, diff --git a/tests/test_photometry_conditional_lf_models.py b/tests/test_photometry_conditional_lf_models.py index f83e5eb..c4dd2d0 100644 --- a/tests/test_photometry_conditional_lf_models.py +++ b/tests/test_photometry_conditional_lf_models.py @@ -4,12 +4,12 @@ import pytest from lfkit.photometry.conditional_lf_models import ( - central_lognormal_conditional_lf, - central_satellite_conditional_lf, conditional_schechter, conditional_schechter_double, conditional_schechter_evolving, - satellite_modified_schechter_conditional_lf, + lognormal_conditional_lf, + modified_schechter_conditional_lf, + two_component_conditional_lf, ) from lfkit.photometry.luminosities import ( luminosity_ratio, @@ -196,8 +196,8 @@ def test_conditional_schechter_double_accepts_callable_parameters() -> None: assert result.dtype == np.float64 -def test_central_lognormal_conditional_lf_matches_expected_formula() -> None: - """Tests the central lognormal conditional LF formula.""" +def test_lognormal_conditional_lf_matches_expected_formula() -> None: + """Tests the lognormal conditional LF formula.""" absolute_mag = np.array([-22.0, -21.0, -20.0]) condition = np.array([0.0, 1.0, 2.0]) @@ -205,7 +205,7 @@ def test_central_lognormal_conditional_lf_matches_expected_formula() -> None: sigma_log_luminosity = np.array([0.2, 0.2, 0.2]) amplitude = np.array([1.0, 2.0, 3.0]) - result = central_lognormal_conditional_lf( + result = lognormal_conditional_lf( absolute_mag=absolute_mag, condition=condition, mean_absolute_mag=mean_absolute_mag, @@ -225,13 +225,13 @@ def test_central_lognormal_conditional_lf_matches_expected_formula() -> None: assert result.dtype == np.float64 -def test_central_lognormal_conditional_lf_accepts_callable_parameters() -> None: - """Tests callable parameters for the central lognormal conditional LF.""" +def test_lognormal_conditional_lf_accepts_callable_parameters() -> None: + """Tests callable parameters for the lognormal conditional LF.""" absolute_mag = np.array([-22.0, -21.0, -20.0]) condition = np.array([0.0, 1.0, 2.0]) - result = central_lognormal_conditional_lf( + result = lognormal_conditional_lf( absolute_mag=absolute_mag, condition=condition, mean_absolute_mag=lambda x: -21.0 - 0.1 * x, @@ -255,11 +255,11 @@ def test_central_lognormal_conditional_lf_accepts_callable_parameters() -> None: assert result.dtype == np.float64 -def test_central_lognormal_conditional_lf_rejects_zero_sigma() -> None: - """Tests that zero central scatter is rejected.""" +def test_lognormal_conditional_lf_rejects_zero_sigma() -> None: + """Tests that zero lognormal scatter is rejected.""" with pytest.raises(ValueError, match="sigma_log_luminosity must be positive."): - central_lognormal_conditional_lf( + lognormal_conditional_lf( absolute_mag=[-22.0, -21.0, -20.0], condition=[0.0, 1.0, 2.0], mean_absolute_mag=-21.0, @@ -268,11 +268,11 @@ def test_central_lognormal_conditional_lf_rejects_zero_sigma() -> None: ) -def test_central_lognormal_conditional_lf_rejects_negative_sigma() -> None: - """Tests that negative central scatter is rejected.""" +def test_lognormal_conditional_lf_rejects_negative_sigma() -> None: + """Tests that negative lognormal scatter is rejected.""" with pytest.raises(ValueError, match="sigma_log_luminosity must be positive."): - central_lognormal_conditional_lf( + lognormal_conditional_lf( absolute_mag=[-22.0, -21.0, -20.0], condition=[0.0, 1.0, 2.0], mean_absolute_mag=-21.0, @@ -281,11 +281,11 @@ def test_central_lognormal_conditional_lf_rejects_negative_sigma() -> None: ) -def test_central_lognormal_conditional_lf_rejects_negative_amplitude() -> None: - """Tests that negative central amplitude is rejected.""" +def test_lognormal_conditional_lf_rejects_negative_amplitude() -> None: + """Tests that negative lognormal amplitude is rejected.""" with pytest.raises(ValueError, match="amplitude must be non-negative."): - central_lognormal_conditional_lf( + lognormal_conditional_lf( absolute_mag=[-22.0, -21.0, -20.0], condition=[0.0, 1.0, 2.0], mean_absolute_mag=-21.0, @@ -294,8 +294,8 @@ def test_central_lognormal_conditional_lf_rejects_negative_amplitude() -> None: ) -def test_satellite_modified_schechter_conditional_lf_matches_expected_formula() -> None: - """Tests the satellite modified-Schechter conditional LF formula.""" +def test_modified_schechter_conditional_lf_matches_expected_formula() -> None: + """Tests the modified-Schechter conditional LF formula.""" absolute_mag = np.array([-22.0, -21.0, -20.0]) condition = np.array([0.0, 1.0, 2.0]) @@ -303,7 +303,7 @@ def test_satellite_modified_schechter_conditional_lf_matches_expected_formula() m_star = np.array([-21.0, -21.1, -21.2]) alpha = np.array([-1.0, -1.1, -1.2]) - result = satellite_modified_schechter_conditional_lf( + result = modified_schechter_conditional_lf( absolute_mag=absolute_mag, condition=condition, phi_star=phi_star, @@ -312,19 +312,21 @@ def test_satellite_modified_schechter_conditional_lf_matches_expected_formula() ) x = luminosity_ratio(absolute_mag, m_star) - expected = 0.4 * np.log(10.0) * phi_star * x ** (alpha + 1.0) * np.exp(-(x**2.0)) + expected = 0.4 * np.log(10.0) * phi_star * x ** (alpha + 1.0) * np.exp( + -(x**2.0) + ) np.testing.assert_allclose(result, expected) assert result.dtype == np.float64 -def test_satellite_modified_schechter_conditional_lf_accepts_callable_parameters() -> None: - """Tests callable parameters for the satellite modified-Schechter LF.""" +def test_modified_schechter_conditional_lf_accepts_callable_parameters() -> None: + """Tests callable parameters for the modified-Schechter conditional LF.""" absolute_mag = np.array([-22.0, -21.0, -20.0]) condition = np.array([0.0, 1.0, 2.0]) - result = satellite_modified_schechter_conditional_lf( + result = modified_schechter_conditional_lf( absolute_mag=absolute_mag, condition=condition, phi_star=lambda x: 1.0e-3 * (1.0 + x), @@ -337,17 +339,19 @@ def test_satellite_modified_schechter_conditional_lf_accepts_callable_parameters alpha = np.array([-1.0, -1.05, -1.1]) x = luminosity_ratio(absolute_mag, m_star) - expected = 0.4 * np.log(10.0) * phi_star * x ** (alpha + 1.0) * np.exp(-(x**2.0)) + expected = 0.4 * np.log(10.0) * phi_star * x ** (alpha + 1.0) * np.exp( + -(x**2.0) + ) np.testing.assert_allclose(result, expected) assert result.dtype == np.float64 -def test_satellite_modified_schechter_conditional_lf_rejects_negative_phi_star() -> None: - """Tests that negative satellite normalization is rejected.""" +def test_modified_schechter_conditional_lf_rejects_negative_phi_star() -> None: + """Tests that negative modified-Schechter normalization is rejected.""" with pytest.raises(ValueError, match="phi_star must be non-negative."): - satellite_modified_schechter_conditional_lf( + modified_schechter_conditional_lf( absolute_mag=[-22.0, -21.0, -20.0], condition=[0.0, 1.0, 2.0], phi_star=-1.0e-3, @@ -356,154 +360,154 @@ def test_satellite_modified_schechter_conditional_lf_rejects_negative_phi_star() ) -def test_central_satellite_conditional_lf_equals_sum_with_explicit_satellite_m_star() -> None: - """Tests that the total conditional LF equals central plus satellite parts.""" +def test_two_component_conditional_lf_equals_sum_with_explicit_modified_m_star() -> None: + """Tests that the two-component LF equals lognormal plus modified parts.""" absolute_mag = np.array([-22.0, -21.0, -20.0]) condition = np.array([0.0, 1.0, 2.0]) - result = central_satellite_conditional_lf( + result = two_component_conditional_lf( absolute_mag=absolute_mag, condition=condition, - central_mean_absolute_mag=-21.0, - central_sigma_log_luminosity=0.2, - satellite_phi_star=1.0e-3, - satellite_alpha=-1.1, - central_amplitude=1.0, - satellite_m_star=-20.5, + lognormal_mean_absolute_mag=-21.0, + lognormal_sigma_log_luminosity=0.2, + modified_phi_star=1.0e-3, + modified_alpha=-1.1, + lognormal_amplitude=1.0, + modified_m_star=-20.5, ) - central = central_lognormal_conditional_lf( + lognormal = lognormal_conditional_lf( absolute_mag=absolute_mag, condition=condition, mean_absolute_mag=-21.0, sigma_log_luminosity=0.2, amplitude=1.0, ) - satellite = satellite_modified_schechter_conditional_lf( + modified = modified_schechter_conditional_lf( absolute_mag=absolute_mag, condition=condition, phi_star=1.0e-3, m_star=-20.5, alpha=-1.1, ) - expected = central + satellite + expected = lognormal + modified np.testing.assert_allclose(result, expected) assert result.dtype == np.float64 -def test_central_satellite_conditional_lf_derives_satellite_m_star() -> None: - """Tests that satellite M-star is derived from the luminosity fraction.""" +def test_two_component_conditional_lf_derives_modified_m_star() -> None: + """Tests that modified M-star is derived from the luminosity fraction.""" absolute_mag = np.array([-22.0, -21.0, -20.0]) condition = np.array([0.0, 1.0, 2.0]) - central_mean_absolute_mag = np.array([-21.0, -21.1, -21.2]) - satellite_luminosity_fraction = np.array([0.5, 0.6, 0.7]) + lognormal_mean_absolute_mag = np.array([-21.0, -21.1, -21.2]) + modified_luminosity_fraction = np.array([0.5, 0.6, 0.7]) - result = central_satellite_conditional_lf( + result = two_component_conditional_lf( absolute_mag=absolute_mag, condition=condition, - central_mean_absolute_mag=central_mean_absolute_mag, - central_sigma_log_luminosity=0.2, - satellite_phi_star=1.0e-3, - satellite_alpha=-1.1, - central_amplitude=1.0, - satellite_m_star=None, - satellite_luminosity_fraction=satellite_luminosity_fraction, + lognormal_mean_absolute_mag=lognormal_mean_absolute_mag, + lognormal_sigma_log_luminosity=0.2, + modified_phi_star=1.0e-3, + modified_alpha=-1.1, + lognormal_amplitude=1.0, + modified_m_star=None, + modified_luminosity_fraction=modified_luminosity_fraction, ) - satellite_m_star = central_mean_absolute_mag + ( - magnitude_difference_from_luminosity_ratio(satellite_luminosity_fraction) + modified_m_star = lognormal_mean_absolute_mag + ( + magnitude_difference_from_luminosity_ratio(modified_luminosity_fraction) ) - central = central_lognormal_conditional_lf( + lognormal = lognormal_conditional_lf( absolute_mag=absolute_mag, condition=condition, - mean_absolute_mag=central_mean_absolute_mag, + mean_absolute_mag=lognormal_mean_absolute_mag, sigma_log_luminosity=0.2, amplitude=1.0, ) - satellite = satellite_modified_schechter_conditional_lf( + modified = modified_schechter_conditional_lf( absolute_mag=absolute_mag, condition=condition, phi_star=1.0e-3, - m_star=satellite_m_star, + m_star=modified_m_star, alpha=-1.1, ) - expected = central + satellite + expected = lognormal + modified np.testing.assert_allclose(result, expected) assert result.dtype == np.float64 -def test_central_satellite_conditional_lf_rejects_zero_luminosity_fraction() -> None: - """Tests that zero satellite luminosity fraction is rejected.""" +def test_two_component_conditional_lf_rejects_zero_luminosity_fraction() -> None: + """Tests that zero modified luminosity fraction is rejected.""" with pytest.raises( ValueError, - match="satellite_luminosity_fraction must be positive.", + match="modified_luminosity_fraction must be positive.", ): - central_satellite_conditional_lf( + two_component_conditional_lf( absolute_mag=[-22.0, -21.0, -20.0], condition=[0.0, 1.0, 2.0], - central_mean_absolute_mag=-21.0, - central_sigma_log_luminosity=0.2, - satellite_phi_star=1.0e-3, - satellite_alpha=-1.1, - central_amplitude=1.0, - satellite_m_star=None, - satellite_luminosity_fraction=0.0, + lognormal_mean_absolute_mag=-21.0, + lognormal_sigma_log_luminosity=0.2, + modified_phi_star=1.0e-3, + modified_alpha=-1.1, + lognormal_amplitude=1.0, + modified_m_star=None, + modified_luminosity_fraction=0.0, ) -def test_central_satellite_conditional_lf_rejects_negative_luminosity_fraction() -> None: - """Tests that negative satellite luminosity fraction is rejected.""" +def test_two_component_conditional_lf_rejects_negative_luminosity_fraction() -> None: + """Tests that negative modified luminosity fraction is rejected.""" with pytest.raises( ValueError, - match="satellite_luminosity_fraction must be positive.", + match="modified_luminosity_fraction must be positive.", ): - central_satellite_conditional_lf( + two_component_conditional_lf( absolute_mag=[-22.0, -21.0, -20.0], condition=[0.0, 1.0, 2.0], - central_mean_absolute_mag=-21.0, - central_sigma_log_luminosity=0.2, - satellite_phi_star=1.0e-3, - satellite_alpha=-1.1, - central_amplitude=1.0, - satellite_m_star=None, - satellite_luminosity_fraction=-0.5, + lognormal_mean_absolute_mag=-21.0, + lognormal_sigma_log_luminosity=0.2, + modified_phi_star=1.0e-3, + modified_alpha=-1.1, + lognormal_amplitude=1.0, + modified_m_star=None, + modified_luminosity_fraction=-0.5, ) -def test_central_satellite_conditional_lf_propagates_invalid_central_component() -> None: - """Tests that invalid central-component parameters are propagated.""" +def test_two_component_conditional_lf_propagates_invalid_lognormal_component() -> None: + """Tests that invalid lognormal-component parameters are propagated.""" - with pytest.raises(ValueError, match="central_sigma_log_luminosity must be positive.|sigma_log_luminosity must be positive."): - central_satellite_conditional_lf( + with pytest.raises(ValueError, match="sigma_log_luminosity must be positive."): + two_component_conditional_lf( absolute_mag=[-22.0, -21.0, -20.0], condition=[0.0, 1.0, 2.0], - central_mean_absolute_mag=-21.0, - central_sigma_log_luminosity=0.0, - satellite_phi_star=1.0e-3, - satellite_alpha=-1.1, - central_amplitude=1.0, - satellite_m_star=-20.5, + lognormal_mean_absolute_mag=-21.0, + lognormal_sigma_log_luminosity=0.0, + modified_phi_star=1.0e-3, + modified_alpha=-1.1, + lognormal_amplitude=1.0, + modified_m_star=-20.5, ) -def test_central_satellite_conditional_lf_propagates_invalid_satellite_component() -> None: - """Tests that invalid satellite-component parameters are propagated.""" +def test_two_component_conditional_lf_propagates_invalid_modified_component() -> None: + """Tests that invalid modified-component parameters are propagated.""" with pytest.raises(ValueError, match="phi_star must be non-negative."): - central_satellite_conditional_lf( + two_component_conditional_lf( absolute_mag=[-22.0, -21.0, -20.0], condition=[0.0, 1.0, 2.0], - central_mean_absolute_mag=-21.0, - central_sigma_log_luminosity=0.2, - satellite_phi_star=-1.0e-3, - satellite_alpha=-1.1, - central_amplitude=1.0, - satellite_m_star=-20.5, + lognormal_mean_absolute_mag=-21.0, + lognormal_sigma_log_luminosity=0.2, + modified_phi_star=-1.0e-3, + modified_alpha=-1.1, + lognormal_amplitude=1.0, + modified_m_star=-20.5, ) diff --git a/tests/test_photometry_lf_integrals.py b/tests/test_photometry_lf_integrals.py index 888f7d7..acbe731 100644 --- a/tests/test_photometry_lf_integrals.py +++ b/tests/test_photometry_lf_integrals.py @@ -93,7 +93,7 @@ def test_integrated_number_density_accepts_scalar_redshift() -> None: def test_integrated_number_density_accepts_broadcastable_scalar_lf_output() -> None: - """Tests that scalar luminosity-function outputs are broadcast.""" + """Tests that scalar luminosity function outputs are broadcast.""" def scalar_lf(m_abs: np.ndarray, z: np.ndarray) -> float: """Return a scalar LF value.""" @@ -156,7 +156,7 @@ def test_integrated_number_density_rejects_nonfinite_magnitude_upper_bound() -> def test_integrated_number_density_rejects_nonfinite_lf_values() -> None: - """Tests that non-finite luminosity-function values are rejected.""" + """Tests that non-finite luminosity function values are rejected.""" def bad_lf(m_abs: np.ndarray, z: np.ndarray) -> np.ndarray: """Return non-finite LF values.""" @@ -172,7 +172,7 @@ def bad_lf(m_abs: np.ndarray, z: np.ndarray) -> np.ndarray: def test_integrated_number_density_rejects_negative_lf_values() -> None: - """Tests that negative luminosity-function values are rejected.""" + """Tests that negative luminosity function values are rejected.""" def bad_lf(m_abs: np.ndarray, z: np.ndarray) -> np.ndarray: """Return negative LF values.""" @@ -188,7 +188,7 @@ def bad_lf(m_abs: np.ndarray, z: np.ndarray) -> np.ndarray: def test_integrated_number_density_rejects_unbroadcastable_lf_values() -> None: - """Tests that unbroadcastable luminosity-function outputs are rejected.""" + """Tests that unbroadcastable luminosity function outputs are rejected.""" def bad_lf(m_abs: np.ndarray, z: np.ndarray) -> np.ndarray: """Return LF values with an invalid shape.""" From c84fd8d3dfac4c04b6709f440cc0c56baacf97a2 Mon Sep 17 00:00:00 2001 From: "niko.sarcevic" Date: Wed, 13 May 2026 23:15:54 -0400 Subject: [PATCH 8/9] udpated function names to be more general --- src/lfkit/api/lumfunc.py | 38 ++--- src/lfkit/photometry/catalog_completeness.py | 6 +- .../photometry/conditional_lf_integrals.py | 59 +++---- src/lfkit/photometry/conditional_lf_models.py | 159 +++++++++--------- src/lfkit/photometry/lf_integrals.py | 8 +- src/lfkit/photometry/lf_parameter_models.py | 12 +- src/lfkit/photometry/lf_redshift_density.py | 2 +- src/lfkit/photometry/luminosity_function.py | 14 +- 8 files changed, 146 insertions(+), 152 deletions(-) diff --git a/src/lfkit/api/lumfunc.py b/src/lfkit/api/lumfunc.py index 6e1b6c1..8f35741 100644 --- a/src/lfkit/api/lumfunc.py +++ b/src/lfkit/api/lumfunc.py @@ -1,4 +1,4 @@ -r"""Public luminosity-function interface. +r"""Public luminosity function interface. This module provides the user-facing :class:`LuminosityFunction` API for evaluating luminosity functions in absolute- or apparent-magnitude space. @@ -30,12 +30,12 @@ schechter_double_from_m, ) from lfkit.photometry.conditional_lf_models import ( - central_lognormal_conditional_lf, - central_satellite_conditional_lf, conditional_schechter, conditional_schechter_double, conditional_schechter_evolving, - satellite_modified_schechter_conditional_lf, + lognormal_conditional_lf, + modified_schechter_conditional_lf, + two_component_conditional_lf, ) from lfkit.photometry.magnitudes import ( absolute_magnitude, @@ -88,7 +88,7 @@ class LuminosityFunction: - """User-facing wrapper for luminosity-function evaluation.""" + """User-facing wrapper for luminosity function evaluation.""" def __init__( self, @@ -97,10 +97,10 @@ def __init__( parameters: Mapping[str, object], meta: Mapping[str, object] | None = None, ) -> None: - """Store a luminosity-function model and its parameters. + """Store a luminosity function model and its parameters. Args: - model: Name of the luminosity-function model. + model: Name of the luminosity function model. parameters: Model parameters passed to the underlying LF function. meta: Optional metadata describing the LF source or calibration. """ @@ -487,7 +487,7 @@ def phi( if z is None: raise ValueError("z is required for a conditional luminosity function.") - return central_lognormal_conditional_lf( + return lognormal_conditional_lf( np.asarray(absolute_mag, dtype=float), np.asarray(z, dtype=float), **self.parameters_dict, @@ -497,7 +497,7 @@ def phi( if z is None: raise ValueError("z is required for a conditional luminosity function.") - return satellite_modified_schechter_conditional_lf( + return modified_schechter_conditional_lf( np.asarray(absolute_mag, dtype=float), np.asarray(z, dtype=float), **self.parameters_dict, @@ -507,13 +507,13 @@ def phi( if z is None: raise ValueError("z is required for a conditional luminosity function.") - return central_satellite_conditional_lf( + return two_component_conditional_lf( np.asarray(absolute_mag, dtype=float), np.asarray(z, dtype=float), **self.parameters_dict, ) - raise ValueError(f"Unsupported luminosity-function model '{self.model}'.") + raise ValueError(f"Unsupported luminosity function model '{self.model}'.") def phi_from_m( self, @@ -575,7 +575,7 @@ def phi_from_m( **self.parameters_dict, ) - raise ValueError(f"Unsupported luminosity-function model '{self.model}'.") + raise ValueError(f"Unsupported luminosity function model '{self.model}'.") def parameters( self, @@ -746,13 +746,13 @@ def central_satellite_conditional( return cls( model="central_satellite_conditional", parameters={ - "central_mean_absolute_mag": central_mean_absolute_mag, - "central_sigma_log_luminosity": central_sigma_log_luminosity, - "satellite_phi_star": satellite_phi_star, - "satellite_alpha": satellite_alpha, - "central_amplitude": central_amplitude, - "satellite_m_star": satellite_m_star, - "satellite_luminosity_fraction": satellite_luminosity_fraction, + "lognormal_mean_absolute_mag": central_mean_absolute_mag, + "lognormal_sigma_log_luminosity": central_sigma_log_luminosity, + "lognormal_amplitude": central_amplitude, + "modified_phi_star": satellite_phi_star, + "modified_alpha": satellite_alpha, + "modified_m_star": satellite_m_star, + "modified_luminosity_fraction": satellite_luminosity_fraction, }, ) diff --git a/src/lfkit/photometry/catalog_completeness.py b/src/lfkit/photometry/catalog_completeness.py index 5dc0ff1..1acece3 100644 --- a/src/lfkit/photometry/catalog_completeness.py +++ b/src/lfkit/photometry/catalog_completeness.py @@ -1,4 +1,4 @@ -r"""Catalog completeness utilities for luminosity-function models. +r"""Catalog completeness utilities for luminosity function models. This module provides helpers for estimating the observed and missing galaxy population implied by a magnitude-limited catalog. These functions are useful @@ -9,13 +9,13 @@ limit and call the generic LF integration helpers to return number densities or fractions. -The core API accepts a luminosity-function callable with signature +The core API accepts a luminosity function callable with signature lf(absolute_mag, z) where ``absolute_mag`` and ``z`` are NumPy arrays that can be broadcast together. This keeps the completeness machinery independent of any specific -luminosity-function parameterization. +luminosity function parameterization. """ from __future__ import annotations diff --git a/src/lfkit/photometry/conditional_lf_integrals.py b/src/lfkit/photometry/conditional_lf_integrals.py index 29098d4..b605719 100644 --- a/src/lfkit/photometry/conditional_lf_integrals.py +++ b/src/lfkit/photometry/conditional_lf_integrals.py @@ -1,8 +1,8 @@ """Conditional luminosity function integration utilities. -This module provides small numerical helpers for conditional luminosity -functions of the form ``Phi(M | x)``, where ``M`` is absolute magnitude and -``x`` is an external conditioning variable. +This module provides numerical helpers for conditional luminosity functions of +the form ``Phi(M | x)``, where ``M`` is absolute magnitude and ``x`` is an +external conditioning variable. The conditioning variable is intentionally generic. It may represent halo mass, environment, galaxy type, richness, stellar mass, or any other quantity. This @@ -22,13 +22,13 @@ from lfkit.utils.validators import validate_array __all__ = [ - "evaluate_conditional_lf", - "integrate_conditional_lf", - "integrate_weighted_conditional_lf", + "evaluate_conditional_luminosity_function", + "integrate_conditional_luminosity_function", + "integrate_weighted_conditional_luminosity_function", ] -def evaluate_conditional_lf( +def evaluate_conditional_luminosity_function( absolute_mag: FloatInput, condition: FloatInput, conditional_lf: Callable[[FloatArray, FloatArray], FloatArray], @@ -38,17 +38,17 @@ def evaluate_conditional_lf( Args: absolute_mag: Absolute magnitude values. condition: Values of the conditioning variable. - conditional_lf: Callable returning ``Phi(M | x)``. The callable must - accept absolute magnitude and condition arrays. + conditional_lf: Callable returning ``Phi(M | x)`` for absolute + magnitude ``M`` and condition ``x``. Returns: - Conditional luminosity function evaluated at the requested absolute - magnitude and condition values. + Conditional luminosity function values evaluated at the requested + absolute magnitudes and conditioning values. Raises: - ValueError: If the inputs or evaluated conditional luminosity function - contain non-finite values, or if the evaluated conditional - luminosity function contains negative values. + ValueError: If the inputs contain non-finite values, or if the + evaluated conditional luminosity function contains non-finite or + negative values. """ absolute_mag_arr = validate_array(absolute_mag, name="absolute_mag") condition_arr = validate_array(condition, name="condition") @@ -66,7 +66,7 @@ def evaluate_conditional_lf( return np.asarray(phi, dtype=np.float64) -def integrate_conditional_lf( +def integrate_conditional_luminosity_function( absolute_mag: FloatInput, condition: FloatInput, conditional_lf: Callable[[FloatArray, FloatArray], FloatArray], @@ -75,25 +75,23 @@ def integrate_conditional_lf( ) -> FloatArray: """Integrate a conditional luminosity function over absolute magnitude. - This computes ``integral Phi(M | x) dM``, where ``x`` is the conditioning - variable. - Args: absolute_mag: Absolute magnitude grid. condition: Values of the conditioning variable. - conditional_lf: Callable returning ``Phi(M | x)``. + conditional_lf: Callable returning ``Phi(M | x)`` for absolute + magnitude ``M`` and condition ``x``. axis: Axis corresponding to the absolute magnitude grid. Returns: Conditional luminosity function integrated over absolute magnitude. Raises: - ValueError: If the inputs or evaluated conditional luminosity function - contain invalid values. + ValueError: If the inputs contain invalid values, or if the evaluated + conditional luminosity function contains invalid values. """ absolute_mag_arr = validate_array(absolute_mag, name="absolute_mag") - phi = evaluate_conditional_lf( + phi = evaluate_conditional_luminosity_function( absolute_mag=absolute_mag_arr, condition=condition, conditional_lf=conditional_lf, @@ -105,7 +103,7 @@ def integrate_conditional_lf( ) -def integrate_weighted_conditional_lf( +def integrate_weighted_conditional_luminosity_function( absolute_mag: FloatInput, condition: FloatInput, conditional_lf: Callable[[FloatArray, FloatArray], FloatArray], @@ -115,14 +113,13 @@ def integrate_weighted_conditional_lf( ) -> FloatArray: """Integrate a weighted conditional luminosity function. - This computes ``integral w(M, x) Phi(M | x) dM``, where ``x`` is the - conditioning variable. - Args: absolute_mag: Absolute magnitude grid. condition: Values of the conditioning variable. - conditional_lf: Callable returning ``Phi(M | x)``. - weight: Callable returning weights ``w(M, x)``. + conditional_lf: Callable returning ``Phi(M | x)`` for absolute + magnitude ``M`` and condition ``x``. + weight: Callable returning weights ``w(M, x)`` for absolute magnitude + ``M`` and condition ``x``. axis: Axis corresponding to the absolute magnitude grid. Returns: @@ -130,13 +127,13 @@ def integrate_weighted_conditional_lf( magnitude. Raises: - ValueError: If the inputs, evaluated conditional luminosity function, - or weights contain invalid values. + ValueError: If the inputs contain invalid values, or if the evaluated + conditional luminosity function or weights contain invalid values. """ absolute_mag_arr = validate_array(absolute_mag, name="absolute_mag") condition_arr = validate_array(condition, name="condition") - phi = evaluate_conditional_lf( + phi = evaluate_conditional_luminosity_function( absolute_mag=absolute_mag_arr, condition=condition_arr, conditional_lf=conditional_lf, diff --git a/src/lfkit/photometry/conditional_lf_models.py b/src/lfkit/photometry/conditional_lf_models.py index c255809..2cc1eaa 100644 --- a/src/lfkit/photometry/conditional_lf_models.py +++ b/src/lfkit/photometry/conditional_lf_models.py @@ -1,4 +1,4 @@ -"""Conditional luminosity-function model utilities. +"""Conditional luminosity function model utilities. This module provides conditional wrappers around existing LFKit luminosity function models. @@ -38,9 +38,9 @@ "conditional_schechter", "conditional_schechter_evolving", "conditional_schechter_double", - "central_lognormal_conditional_lf", - "satellite_modified_schechter_conditional_lf", - "central_satellite_conditional_lf", + "lognormal_conditional_lf", + "modified_schechter_conditional_lf", + "two_component_conditional_lf", ] @@ -65,7 +65,7 @@ def conditional_schechter( ``condition``. Returns: - Conditional Schechter luminosity-function values. + Conditional Schechter luminosity function values. """ condition_arr = validate_array(condition, name="condition") @@ -102,7 +102,7 @@ def conditional_schechter_evolving( ) -> FloatArray: """Evaluate a conditional Schechter LF using LFKit parameter models. - This is the conditional-LF analogue of ``schechter_evolving``. The + This is the conditional LF analogue of ``schechter_evolving``. The conditioning variable is passed to LFKit's registered parameter models. Args: @@ -116,7 +116,7 @@ def conditional_schechter_evolving( alpha_kwargs: Keyword arguments passed to the selected ``alpha`` model. Returns: - Conditional Schechter luminosity-function values. + Conditional Schechter luminosity function values. Raises: ValueError: If an unsupported parameter model is requested. @@ -166,7 +166,7 @@ def conditional_schechter_double( callable of ``condition``. Returns: - Conditional double-power-law Schechter luminosity-function values. + Conditional double-power-law Schechter luminosity function values. """ condition_arr = validate_array(condition, name="condition") @@ -192,7 +192,7 @@ def conditional_schechter_double( ) -def central_lognormal_conditional_lf( +def lognormal_conditional_lf( absolute_mag: FloatInput, condition: FloatInput, *, @@ -200,20 +200,20 @@ def central_lognormal_conditional_lf( sigma_log_luminosity: ConditionalParameter, amplitude: ConditionalParameter = 1.0, ) -> FloatArray: - """Evaluate a central-galaxy lognormal conditional LF in magnitudes. + """Evaluate a lognormal conditional luminosity function in magnitudes. Args: absolute_mag: Absolute magnitude value(s). condition: Values of the conditioning variable. - mean_absolute_mag: Mean central-galaxy absolute magnitude. May be - scalar, array-like, or callable of ``condition``. - sigma_log_luminosity: Scatter in ``log10(L)`` at fixed condition. May - be scalar, array-like, or callable of ``condition``. - amplitude: Non-negative amplitude of the central component. May be - scalar, array-like, or callable of ``condition``. + mean_absolute_mag: Mean absolute magnitude. + May be scalar, array-like, or callable of ``condition``. + sigma_log_luminosity: Scatter in ``log10(L)`` at fixed condition. + May be scalar, array-like, or callable of ``condition``. + amplitude: Non-negative amplitude of the component. + May be scalar, array-like, or callable of ``condition``. Returns: - Central-galaxy conditional luminosity-function values. + Lognormal conditional luminosity function values. """ absolute_mag_arr = validate_array(absolute_mag, name="absolute_mag") condition_arr = validate_array(condition, name="condition") @@ -249,10 +249,10 @@ def central_lognormal_conditional_lf( * np.exp(-0.5 * (delta_log_luminosity / sigma_log_luminosity_arr) ** 2.0) ) - return _validate_lf_output(phi, name="central_lognormal_conditional_lf") + return _validate_lf_output(phi, name="lognormal_conditional_lf") -def satellite_modified_schechter_conditional_lf( +def modified_schechter_conditional_lf( absolute_mag: FloatInput, condition: FloatInput, *, @@ -260,23 +260,20 @@ def satellite_modified_schechter_conditional_lf( m_star: ConditionalParameter, alpha: ConditionalParameter, ) -> FloatArray: - """Evaluate a satellite modified-Schechter conditional LF. + """Evaluate a modified Schechter conditional luminosity function. - This is the Cacciato-style satellite form with ``exp[-(L/L_star)^2]`` - instead of the standard Schechter ``exp[-L/L_star]``. + This uses a squared exponential cutoff in luminosity ratio instead of the + standard Schechter exponential cutoff. Args: absolute_mag: Absolute magnitude value(s). condition: Values of the conditioning variable. - phi_star: Satellite normalization. May be scalar, array-like, or - callable of ``condition``. - m_star: Satellite characteristic absolute magnitude. May be scalar, - array-like, or callable of ``condition``. - alpha: Satellite faint-end slope. May be scalar, array-like, or - callable of ``condition``. + phi_star: Component normalization. + m_star: Characteristic absolute magnitude. + alpha: Faint-end slope. Returns: - Satellite conditional luminosity-function values. + Modified Schechter conditional luminosity function values. """ absolute_mag_arr = validate_array(absolute_mag, name="absolute_mag") condition_arr = validate_array(condition, name="condition") @@ -306,95 +303,95 @@ def satellite_modified_schechter_conditional_lf( return _validate_lf_output( phi, - name="satellite_modified_schechter_conditional_lf", + name="modified_schechter_conditional_lf", ) -def central_satellite_conditional_lf( +def two_component_conditional_lf( absolute_mag: FloatInput, condition: FloatInput, *, - central_mean_absolute_mag: ConditionalParameter, - central_sigma_log_luminosity: ConditionalParameter, - satellite_phi_star: ConditionalParameter, - satellite_alpha: ConditionalParameter, - central_amplitude: ConditionalParameter = 1.0, - satellite_m_star: ConditionalParameter | None = None, - satellite_luminosity_fraction: ConditionalParameter = 0.562, + lognormal_mean_absolute_mag: ConditionalParameter, + lognormal_sigma_log_luminosity: ConditionalParameter, + modified_phi_star: ConditionalParameter, + modified_alpha: ConditionalParameter, + lognormal_amplitude: ConditionalParameter = 1.0, + modified_m_star: ConditionalParameter | None = None, + modified_luminosity_fraction: ConditionalParameter = 0.562, ) -> FloatArray: - """Evaluate the sum of central and satellite conditional LF components. + """Evaluate the sum of lognormal and modified Schechter components. Args: absolute_mag: Absolute magnitude value(s). condition: Values of the conditioning variable. - central_mean_absolute_mag: Mean central-galaxy absolute magnitude. May - be scalar, array-like, or callable of ``condition``. - central_sigma_log_luminosity: Scatter in ``log10(L)`` for the central + lognormal_mean_absolute_mag: Mean absolute magnitude of the lognormal component. May be scalar, array-like, or callable of ``condition``. - satellite_phi_star: Satellite normalization. May be scalar, array-like, - or callable of ``condition``. - satellite_alpha: Satellite faint-end slope. May be scalar, array-like, - or callable of ``condition``. - central_amplitude: Non-negative amplitude of the central component. May - be scalar, array-like, or callable of ``condition``. - satellite_m_star: Satellite characteristic absolute magnitude. If - omitted, it is derived from ``central_mean_absolute_mag`` and - ``satellite_luminosity_fraction``. - satellite_luminosity_fraction: Ratio ``L_sat_star / L_c`` used only - when ``satellite_m_star`` is omitted. + lognormal_sigma_log_luminosity: Scatter in ``log10(L)`` for the + lognormal component. May be scalar, array-like, or callable of + ``condition``. + modified_phi_star: Normalization of the modified Schechter component. + May be scalar, array-like, or callable of ``condition``. + modified_alpha: Faint-end slope of the modified Schechter component. + May be scalar, array-like, or callable of ``condition``. + lognormal_amplitude: Non-negative amplitude of the lognormal component. + May be scalar, array-like, or callable of ``condition``. + modified_m_star: Characteristic absolute magnitude of the modified + Schechter component. If omitted, it is derived from + ``lognormal_mean_absolute_mag`` and ``modified_luminosity_fraction``. + modified_luminosity_fraction: Ratio used to derive the modified + Schechter characteristic luminosity from the lognormal mean + luminosity when ``modified_m_star`` is omitted. Returns: - Total central-plus-satellite conditional luminosity-function values. + Combined conditional luminosity function values. """ condition_arr = validate_array(condition, name="condition") - central_mean_absolute_mag_arr = _evaluate_conditional_parameter( - central_mean_absolute_mag, + lognormal_mean_absolute_mag_arr = _evaluate_conditional_parameter( + lognormal_mean_absolute_mag, condition_arr, - name="central_mean_absolute_mag", + name="lognormal_mean_absolute_mag", ) - central_phi = central_lognormal_conditional_lf( + lognormal_phi = lognormal_conditional_lf( absolute_mag=absolute_mag, condition=condition_arr, - mean_absolute_mag=central_mean_absolute_mag_arr, - sigma_log_luminosity=central_sigma_log_luminosity, - amplitude=central_amplitude, + mean_absolute_mag=lognormal_mean_absolute_mag_arr, + sigma_log_luminosity=lognormal_sigma_log_luminosity, + amplitude=lognormal_amplitude, ) - if satellite_m_star is None: - satellite_luminosity_fraction_arr = _evaluate_conditional_parameter( - satellite_luminosity_fraction, + if modified_m_star is None: + modified_luminosity_fraction_arr = _evaluate_conditional_parameter( + modified_luminosity_fraction, condition_arr, - name="satellite_luminosity_fraction", + name="modified_luminosity_fraction", ) - if np.any(satellite_luminosity_fraction_arr <= 0.0): - raise ValueError("satellite_luminosity_fraction must be positive.") + if np.any(modified_luminosity_fraction_arr <= 0.0): + raise ValueError("modified_luminosity_fraction must be positive.") - satellite_m_star_arr = central_mean_absolute_mag_arr + ( - magnitude_difference_from_luminosity_ratio( - satellite_luminosity_fraction_arr - ) + modified_m_star_arr = lognormal_mean_absolute_mag_arr + ( + magnitude_difference_from_luminosity_ratio(modified_luminosity_fraction_arr) ) else: - satellite_m_star_arr = _evaluate_conditional_parameter( - satellite_m_star, + modified_m_star_arr = _evaluate_conditional_parameter( + modified_m_star, condition_arr, - name="satellite_m_star", + name="modified_m_star", ) - satellite_phi = satellite_modified_schechter_conditional_lf( + modified_phi = modified_schechter_conditional_lf( absolute_mag=absolute_mag, condition=condition_arr, - phi_star=satellite_phi_star, - m_star=satellite_m_star_arr, - alpha=satellite_alpha, + phi_star=modified_phi_star, + m_star=modified_m_star_arr, + alpha=modified_alpha, ) return _validate_lf_output( - central_phi + satellite_phi, - name="central_satellite_conditional_lf", + lognormal_phi + modified_phi, + name="two_component_conditional_lf", ) @@ -418,7 +415,7 @@ def _validate_lf_output( *, name: str, ) -> FloatArray: - """Validate luminosity-function model output.""" + """Validate luminosity function model output.""" phi_arr = validate_array(phi, name=name) if np.any(phi_arr < 0.0): diff --git a/src/lfkit/photometry/lf_integrals.py b/src/lfkit/photometry/lf_integrals.py index b045cde..81567f1 100644 --- a/src/lfkit/photometry/lf_integrals.py +++ b/src/lfkit/photometry/lf_integrals.py @@ -1,15 +1,15 @@ r"""Luminosity-function integration utilities. -This module provides generic numerical integrals of luminosity-function +This module provides generic numerical integrals of luminosity function callables over finite absolute-magnitude ranges. -The core API accepts a luminosity-function callable with signature +The core API accepts a luminosity function callable with signature lf(absolute_mag, z) where ``absolute_mag`` and ``z`` are NumPy arrays that can be broadcast together. This keeps the integration machinery independent of any specific -luminosity-function parameterization, catalog selection, or cosmology backend. +luminosity function parameterization, catalog selection, or cosmology backend. These helpers are intentionally generic. Catalog completeness, LF-dependent redshift densities, luminosity-density calculations, and selection-weighted @@ -83,7 +83,7 @@ def lf_weighted_integral( weight_fn: Callable[[FloatArray, FloatArray], FloatArray], n_m: int = 512, ) -> FloatArray: - r"""Return a weighted luminosity-function integral. + r"""Return a weighted luminosity function integral. This computes diff --git a/src/lfkit/photometry/lf_parameter_models.py b/src/lfkit/photometry/lf_parameter_models.py index b947423..4489786 100644 --- a/src/lfkit/photometry/lf_parameter_models.py +++ b/src/lfkit/photometry/lf_parameter_models.py @@ -1,10 +1,10 @@ r"""Luminosity-function parameter evolution models for LFKit. This module provides helper functions for evaluating redshift-dependent -luminosity-function parameters such as ``phi_star(z)``, ``M_star(z)``, +luminosity function parameters such as ``phi_star(z)``, ``M_star(z)``, and ``alpha(z)``. -These helpers are used by the main luminosity-function evaluators but do +These helpers are used by the main luminosity function evaluators but do not evaluate the luminosity function themselves. Built-in options include constant evolution and simple linearized forms @@ -201,7 +201,7 @@ def get_parameter_model( *, model_kind: str, ) -> ParameterModel: - r"""Return a registered luminosity-function parameter model. + r"""Return a registered luminosity function parameter model. Args: model_name: Name of the requested model. @@ -225,7 +225,7 @@ def get_parameter_model( def available_lf_parameter_models() -> dict[str, list[str]]: - """Return available luminosity-function parameter evolution models.""" + """Return available luminosity function parameter evolution models.""" return { "phi_star": sorted(PHI_STAR_MODELS), "m_star": sorted(M_STAR_MODELS), @@ -241,7 +241,7 @@ def _register_parameter_model( model_kind: str, overwrite: bool = False, ) -> None: - """Register a luminosity-function parameter evolution model.""" + """Register a luminosity function parameter evolution model.""" if not name: raise ValueError(f"{model_kind} model name cannot be empty.") @@ -315,7 +315,7 @@ def evaluate_lf_parameters( alpha_model: str = "constant", alpha_kwargs: Mapping[str, ParameterValue] | None = None, ) -> tuple[FloatArray, FloatArray, FloatArray]: - r"""Evaluate evolving luminosity-function parameters at redshift ``z``. + r"""Evaluate evolving luminosity function parameters at redshift ``z``. Args: z: Redshift value or array-like of redshift values. diff --git a/src/lfkit/photometry/lf_redshift_density.py b/src/lfkit/photometry/lf_redshift_density.py index 9a619b1..759bf92 100644 --- a/src/lfkit/photometry/lf_redshift_density.py +++ b/src/lfkit/photometry/lf_redshift_density.py @@ -1,6 +1,6 @@ r"""Luminosity-function redshift-density utilities. -This module provides helpers for converting a luminosity-function callable into +This module provides helpers for converting a luminosity function callable into an LF-integrated redshift-density curve. The core operation is diff --git a/src/lfkit/photometry/luminosity_function.py b/src/lfkit/photometry/luminosity_function.py index c4b5f78..2565573 100644 --- a/src/lfkit/photometry/luminosity_function.py +++ b/src/lfkit/photometry/luminosity_function.py @@ -1,7 +1,7 @@ """Luminosity function utilities for LFKit. This module provides simple standalone functions for evaluating -common galaxy luminosity-function parameterization. +common galaxy luminosity function parameterization. All luminosity functions in this module are defined in rest-frame absolute magnitude space. Functions with names ending in ``_from_m`` are convenience @@ -101,7 +101,7 @@ def schechter( alpha: Faint-end slope. Can be a scalar or array-like. Returns: - NumPy array of luminosity-function values evaluated at + NumPy array of luminosity function values evaluated at ``absolute_mag``. """ x = luminosity_ratio(absolute_mag, m_star) @@ -163,7 +163,7 @@ def schechter_evolving( ``alpha`` evolution model. Returns: - NumPy array of luminosity-function values evaluated at + NumPy array of luminosity function values evaluated at ``absolute_mag`` and ``z``. Raises: @@ -237,7 +237,7 @@ def schechter_double( m_transition: Transition magnitude M_t corresponding to L_t. Returns: - NumPy array of luminosity-function values evaluated at + NumPy array of luminosity function values evaluated at ``absolute_mag``. """ absolute_mag = validate_array(absolute_mag, name="absolute_mag") @@ -324,7 +324,7 @@ def schechter_from_m( e_correction: Optional evolution-correction term(s). Returns: - NumPy array of luminosity-function values corresponding to the + NumPy array of luminosity function values corresponding to the supplied apparent magnitudes. """ z = validate_array(z, name="z") @@ -405,7 +405,7 @@ def schechter_evolving_from_m( e_correction: Optional evolution-correction term(s). Returns: - NumPy array of luminosity-function values corresponding to the + NumPy array of luminosity function values corresponding to the supplied apparent magnitudes. """ z = validate_array(z, name="z") @@ -505,7 +505,7 @@ def schechter_double_from_m( e_correction: Optional evolution-correction term(s). Returns: - NumPy array of luminosity-function values corresponding to the + NumPy array of luminosity function values corresponding to the supplied apparent magnitudes. """ z = validate_array(z, name="z") From ff74e3c5a87f3b8f0fc17332c0e2603c5cf42ebe Mon Sep 17 00:00:00 2001 From: "niko.sarcevic" Date: Wed, 13 May 2026 23:21:01 -0400 Subject: [PATCH 9/9] added cacciato ref to bench fiels --- tests/benchmarks/test_cacciato_clf_reference.py | 2 ++ tests/benchmarks/test_cacciato_hod_reference.py | 2 ++ 2 files changed, 4 insertions(+) diff --git a/tests/benchmarks/test_cacciato_clf_reference.py b/tests/benchmarks/test_cacciato_clf_reference.py index d527ab9..1d4fae8 100644 --- a/tests/benchmarks/test_cacciato_clf_reference.py +++ b/tests/benchmarks/test_cacciato_clf_reference.py @@ -7,6 +7,8 @@ They are intentionally written as reference/benchmark tests rather than normal unit tests. They check that LFKit preserves the expected central lognormal, satellite modified-Schechter, and central-plus-satellite behaviour. + +For more information, see https://arxiv.org/abs/1207.0503. """ import numpy as np diff --git a/tests/benchmarks/test_cacciato_hod_reference.py b/tests/benchmarks/test_cacciato_hod_reference.py index dc5df65..d2fbce8 100644 --- a/tests/benchmarks/test_cacciato_hod_reference.py +++ b/tests/benchmarks/test_cacciato_hod_reference.py @@ -5,6 +5,8 @@ function helpers. The goal is to check that integrating the central and satellite CLFs over a luminosity window gives the same halo occupation numbers when expressed through LFKit. + +For more information, see https://arxiv.org/abs/1207.0503. """ from __future__ import annotations