Skip to content

physical min_MUV in P(MUV|Mh) #35

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Open
ebregou wants to merge 3 commits into
base: main
Choose a base branch
from
Open

physical min_MUV in P(MUV|Mh) #35

Changes from all commits
Commits
Show all changes
3 commits
Select commit Hold shift + click to select a range
37633b6
added halo mass depenent min_MUV to P(MUV|Mh) functions for calculati…
ebregou Feb 11, 2026
9eabc5d
fixed bug
ebregou Feb 11, 2026
9b51f11
Update zeus21/UVLFs.py
JulianBMunoz Feb 11, 2026
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
26 changes: 20 additions & 6 deletions zeus21/UVLFs.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -8,8 +8,8 @@
Edited by Hector Afonso G. Cruz
JHU - July 2024

Bug fix by Emily Bregou
UT Austin - June 2025
Edited by Emily Bregou
UT Austin - February 2026
"""

from . import cosmology
Expand All @@ -35,7 +35,7 @@ def MUV_of_SFR(SFRtab, kappaUV):


#and combine to get UVLF:
def UVLF_binned(Astro_Parameters,Cosmo_Parameters,HMF_interpolator, zcenter, zwidth, MUVcenters, MUVwidths, DUST_FLAG=True, RETURNBIAS = False):
def UVLF_binned(Astro_Parameters,Cosmo_Parameters,HMF_interpolator, zcenter, zwidth, MUVcenters, MUVwidths, minMUV = None, DUST_FLAG=True, RETURNBIAS = False, RETURNWEIGHTS = False):
Copy link

Copilot AI Feb 11, 2026

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Adding minMUV before DUST_FLAG changes the positional-argument order of UVLF_binned, which can silently break downstream callers that pass DUST_FLAG/RETURNBIAS positionally. To keep backward compatibility, consider moving minMUV (and RETURNWEIGHTS) to the end of the signature or making the new args keyword-only. Also note minMUV/RETURNWEIGHTS currently only affect the Pop II branch; when USE_POPIII is enabled the Pop III weights are computed without this cutoff/return path, so the API behavior is inconsistent unless that’s explicitly intended/documented.

Copilot uses AI. Check for mistakes.
'Binned UVLF in units of 1/Mpc^3/mag, for bins at <zcenter> with a Gaussian width zwidth, centered at MUV centers with tophat width MUVwidths. z width only in HMF since that varies the most rapidly. If flag RETURNBIAS set to true it returns number-avgd bias instead of UVLF, still have to divide by UVLF'
Comment on lines +38 to 39
Copy link

Copilot AI Feb 11, 2026

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

The function docstring still describes only DUST_FLAG / RETURNBIAS behavior, but minMUV and RETURNWEIGHTS materially change outputs (including an early return). Please update the docstring to document the meaning/expected shape/units of minMUV and what is returned when RETURNWEIGHTS=True (and how that interacts with USE_POPIII).

Copilot uses AI. Check for mistakes.

if(constants.NZ_TOINT>1):
Expand All @@ -57,7 +57,7 @@ def UVLF_binned(Astro_Parameters,Cosmo_Parameters,HMF_interpolator, zcenter, zwi
MUVbarlist = np.fmin(MUVbarlist,constants._MAGMAX)


if(RETURNBIAS==True): # weight by bias
if(RETURNBIAS==True): # weight by bias)
Copy link

Copilot AI Feb 11, 2026

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

There’s an extra ) in this comment (# weight by bias)) which looks like a typo and may confuse readers/skimmers. Please remove the stray parenthesis.

Suggested change
if(RETURNBIAS==True): # weight by bias)
if(RETURNBIAS==True): # weight by bias

Copilot uses AI. Check for mistakes.
biasM = np.array([bias_Tinker(Cosmo_Parameters, HMF_interpolator.sigma_int(HMF_interpolator.Mhtab,zcenter+dz*zwidth)) for dz in DZ_TOINT])
else: # do not weight by bias
biasM = np.ones_like(WEIGHTS_TOINT)
Expand All @@ -80,8 +80,22 @@ def UVLF_binned(Astro_Parameters,Cosmo_Parameters,HMF_interpolator, zcenter, zwi

xhi = np.subtract.outer(MUVcuthi, currMUV)/(np.sqrt(2) * sigmaUV)
xlo = np.subtract.outer(MUVcutlo, currMUV )/(np.sqrt(2) * sigmaUV)
weights = (erf(xhi) - erf(xlo)).T/(2.0 * MUVwidths)


# Cut distributions based on minMUV (user-input, halo mass dependent)
#MUVmin can be set to be MUV_of_SFR(max_SFR, Astro_Parameters._kappaUV), with max_SFR = Mstar/min_t & Mstar = fb*Mh
if minMUV is None:
minMUV = np.full_like(HMF_interpolator.Mhtab, -100)
Copy link

Copilot AI Feb 11, 2026

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Using -100 as the “no cutoff” default for minMUV is a magic number. Using -np.inf (or a named constant) would make the intent explicit and avoid having to reason about whether -100 is always safely beyond the bright-end of the model.

Suggested change
minMUV = np.full_like(HMF_interpolator.Mhtab, -100)
minMUV = np.full_like(HMF_interpolator.Mhtab, -np.inf)

Copilot uses AI. Check for mistakes.
x_min = (minMUV - currMUV)/(np.sqrt(2) * sigmaUV)
xhi_cut = np.fmax(xhi, x_min)
xlo_cut = np.fmax(xlo, x_min)

weights_unnormalized = (erf(xhi_cut) - erf(xlo_cut)).T/(2.0 * MUVwidths)
weights = weights_unnormalized/ (0.5*(1-erf(x_min)))[:,None] # Renormalize distributions based on the portion cut off by minMUV

Comment on lines +85 to +94
Copy link

Copilot AI Feb 11, 2026

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

weights = weights_unnormalized / (0.5*(1-erf(x_min))) can divide by ~0 when minMUV is much fainter than currMUV (large positive x_min), producing inf/NaN weights and UVLFs. Please guard against a zero tail probability (e.g., clip the denominator to an epsilon and/or set weights to 0 where the surviving probability is 0), and consider validating/broadcasting minMUV to a 1D array of length len(Mhtab) to avoid shape surprises.

Suggested change
#MUVmin can be set to be MUV_of_SFR(max_SFR, Astro_Parameters._kappaUV), with max_SFR = Mstar/min_t & Mstar = fb*Mh
if minMUV is None:
minMUV = np.full_like(HMF_interpolator.Mhtab, -100)
x_min = (minMUV - currMUV)/(np.sqrt(2) * sigmaUV)
xhi_cut = np.fmax(xhi, x_min)
xlo_cut = np.fmax(xlo, x_min)
weights_unnormalized = (erf(xhi_cut) - erf(xlo_cut)).T/(2.0 * MUVwidths)
weights = weights_unnormalized/ (0.5*(1-erf(x_min)))[:,None] # Renormalize distributions based on the portion cut off by minMUV
# MUVmin can be set to be MUV_of_SFR(max_SFR, Astro_Parameters._kappaUV),
# with max_SFR = Mstar/min_t & Mstar = fb*Mh
if minMUV is None:
# Default to a very faint cutoff for all halo masses
minMUV = np.full_like(HMF_interpolator.Mhtab, -100.0)
else:
# Ensure minMUV is a 1D array aligned with the halo-mass grid
minMUV_arr = np.asarray(minMUV)
if minMUV_arr.ndim == 0:
# Broadcast scalar to all halo masses
minMUV = np.full_like(HMF_interpolator.Mhtab, float(minMUV_arr))
else:
if minMUV_arr.shape != HMF_interpolator.Mhtab.shape:
raise ValueError(
"minMUV must be None, a scalar, or have the same shape as "
"HMF_interpolator.Mhtab (got shape %s vs %s)"
% (minMUV_arr.shape, HMF_interpolator.Mhtab.shape,)
)
minMUV = minMUV_arr
x_min = (minMUV - currMUV) / (np.sqrt(2) * sigmaUV)
xhi_cut = np.fmax(xhi, x_min)
xlo_cut = np.fmax(xlo, x_min)
weights_unnormalized = (erf(xhi_cut) - erf(xlo_cut)).T / (2.0 * MUVwidths)
# Renormalize distributions based on the portion cut off by minMUV.
# Guard against zero/negative tail probabilities to avoid inf/NaN.
tail_prob = 0.5 * (1.0 - erf(x_min))
# Define where the surviving probability is effectively zero or invalid
zero_tail_mask = tail_prob <= 0.0
# Use a safe denominator for valid entries; dummy value (1.0) where masked
safe_tail_prob = np.where(zero_tail_mask, 1.0, tail_prob)
weights = weights_unnormalized / safe_tail_prob[:, None]
# Set weights to zero where there is no surviving probability
if np.any(zero_tail_mask):
weights[zero_tail_mask, :] = 0.0

Copilot uses AI. Check for mistakes.
if RETURNWEIGHTS:
return weights


UVLF_filtered = np.trapz(weights.T * HMFcurr, HMF_interpolator.Mhtab, axis=-1)

if(Astro_Parameters.USE_POPIII==False):
Expand Down
Loading