J and D factors#
For a density \(\rho\), annihilation and decay factors are
These are NumPy postprocessing operations. Use physical posterior draws from either sampler; keep the same distance, aperture, cutoff and density convention when comparing results. The returned units are GeV² cm⁻⁵ and GeV cm⁻².
Spherical apertures#
The spherical DMModel.jfactor_cone method includes outer halo shells
projected into a finite cone through the truncated density.
jfactor_spherical_aperture integrates only to
min(dist_pc*sin(roi_deg), r_t_pc), omitting shells outside that sphere.
NFWModel.jfactor_small_angle_infinite_los uses the Evans et al. formula:
it caps the projected aperture at r_t_pc but integrates the untruncated
profile along the entire line of sight. These are different integrals.
The small-angle approximations enforce small_angle_limit_deg (default 1 degree)
by raising ValueError; the full cone does not use that bound. The
geometry derivation follows
Ullio & Valli (2016).
The spherical limit of AxisymmetricZhaoModel also supplies a D-factor
calculation for a spherical halo; there is no separate NumPy/SciPy spherical
DMModel.dfactor API.
"""Finite-cone factors for spherical and flattened halos."""
import numpy as np
from jeanspy.model import NFWModel
from jeanspy.axisymmetric import AxisymmetricZhaoModel
from jeanspy.axisymmetric_factors import jfactor, dfactor
# spherical-factors-start
halo = NFWModel(rs_pc=500., rhos_Msunpc3=.1, r_t_pc=5000.)
distance_pc, aperture_deg = 76000., .5
J_GeV2_cm_minus5 = halo.jfactor_cone(distance_pc, aperture_deg)
# The spherical limit of the spheroidal factor integrator also supplies D.
spherical = AxisymmetricZhaoModel(rhos_Msunpc3=.1, rs_pc=500., Q=1., alpha=1., beta=3., gamma=1., r_t_pc=5000.)
D_GeV_cm_minus2 = dfactor(spherical, distance_pc, aperture_deg,
n_mu=48, n_phi=32, n_radial=48)
assert J_GeV2_cm_minus5 > 0 and D_GeV_cm_minus2 > 0
print("log10 J, log10 D:", np.log10(J_GeV2_cm_minus5), np.log10(D_GeV_cm_minus2))
Axisymmetric finite-cone factors#
roi_deg is a circular cone half-angle in degrees, in [0,90). A finite
r_t_pc must be supplied. The observer must lie outside a sphere enclosing the
halo: dist_pc > r_t_pc * max(1,Q). J needs gamma < 1.5; an integrated central
aperture with a steeper annihilation cusp diverges and raises. No central core
or artificial radius floor regularizes this divergence. D is finite over the
supported density-slope domain. See the axisymmetric guide
for the physical parameters and viewing geometry.
The implementation evaluates the exact observer integral
Spheroidal volume coordinates have d³r = Q m² dm dmu dphi. For a unit
spheroidal direction with LOS component L and sky-plane length B, the cone
bounds m by m*(B*cos(theta)-L*sin(theta)) <= dist_pc*sin(theta) and the
observer distance is s² = dist_pc² + 2*dist_pc*m*L + m²*|e|².
This retains finite-distance geometry and contributions from outer shells
projected into the aperture. It does not replace the aperture by a sphere.
Radial quadrature is cusp-regularized analytically and split at rs; angles use
Gauss–Legendre and periodic azimuth quadrature. Expose accuracy with n_mu,
n_phi, n_radial (defaults 96,96,128). The cone/cutoff intersection can
require angular refinement. These settings are independent of the three Jeans
quadrature orders. J/D evaluation is a NumPy postprocessing operation for
chains from either backend, as in the spherical workflow.
For a flattened halo the factor depends on viewing inclination:
flattened = AxisymmetricZhaoModel(rhos_Msunpc3=.1, rs_pc=500., Q=.8, alpha=2., beta=4., gamma=.5, r_t_pc=5000.)
J_flat = jfactor(flattened, distance_pc, aperture_deg, inclination=1.1,
n_mu=48, n_phi=32, n_radial=48)
D_flat = dfactor(flattened, distance_pc, aperture_deg, inclination=1.1,
n_mu=48, n_phi=32, n_radial=48)
assert np.isfinite([J_flat, D_flat]).all() and min(J_flat, D_flat) > 0
# Refine these independent factor quadratures for the halo/aperture in use.