jwspecabund.dust
Dust correction for emission-line fluxes.
Implements Salim+18/Noll+09 (Calzetti base + UV bump + slope deviation) and Cardelli+89 Milky Way extinction, plus A_V derivation from the Balmer decrement and Lyα escape fraction computation.
Functions
|
Cardelli, Clayton & Mathis (1989) MW extinction curve. |
|
Derive A_V from a single, explicitly-chosen Balmer decrement. |
|
Derive A_V from a Balmer decrement. |
|
Derive A_V from every available Balmer decrement. |
|
Compute Lyα escape fraction from dust-corrected Balmer lines. |
|
Compute Lyα escape fraction with MC error propagation. |
|
Apply dust correction to a set of emission-line fluxes. |
|
Salim+18/Noll+09 modified Calzetti attenuation curve. |
- jwspecabund.dust.cardelli_extinction(wave_A, Av, Rv=3.1)[source]
Cardelli, Clayton & Mathis (1989) MW extinction curve.
- jwspecabund.dust.compute_Av_balmer_pair(fluxes, errors, pair, law='salim', **kwargs)[source]
Derive A_V from a single, explicitly-chosen Balmer decrement.
Like
compute_Av_multi_balmer()but restricted to exactly one(numerator, denominator)pair — useful when you want to force the dust correction onto a high-SNR pair (e.g.("Ha", "HBETA")) and ignore lower-SNR Balmer lines.- Parameters:
fluxes (dict) – Observed (not dust-corrected) line fluxes and 1σ errors.
errors (dict) – Observed (not dust-corrected) line fluxes and 1σ errors.
pair (tuple of str) –
(numerator, denominator)line names from the Balmer ladder ("Ha","HBETA","HGAMMA","HDELTA","H9","H10"), e.g.("Ha", "HBETA")for Hα/Hβ.law (str) –
"salim"(default) or"cardelli".**kwargs – Extra arguments to the attenuation law (Rv, delta, B).
- Returns:
Same shape as
compute_Av_multi_balmer():"Av","Av_err","individual","n_lines"(1 if usable, else 0),"anchor"(the denominator line).- Return type:
- Raises:
ValueError – If either line name is not in the Balmer ladder.
- jwspecabund.dust.compute_Av_from_balmer(flux_num, flux_den, flux_num_err, flux_den_err, law='salim', intrinsic_ratio=0.468, wave_num_A=4341.68, wave_den_A=4862.68, **kwargs)[source]
Derive A_V from a Balmer decrement.
By default uses Hgamma/Hbeta, but any pair of Balmer lines can be used by specifying the wavelengths and intrinsic ratio.
- Parameters:
flux_num (float) – Observed flux of the numerator line (e.g. Hgamma).
flux_den (float) – Observed flux of the denominator line (e.g. Hbeta).
flux_num_err (float) – Error on numerator flux.
flux_den_err (float) – Error on denominator flux.
law (str) –
"salim"(default) or"cardelli".intrinsic_ratio (float) – Intrinsic flux ratio (default 0.468 for Hgamma/Hbeta).
wave_num_A (float) – Rest wavelength of numerator line in Angstroms (default Hgamma).
wave_den_A (float) – Rest wavelength of denominator line in Angstroms (default Hbeta).
**kwargs – Extra arguments to the attenuation law.
- Returns:
(Av, Av_err).- Return type:
- jwspecabund.dust.compute_Av_multi_balmer(fluxes, errors, law='salim', snr_min=3.0, anchor='HBETA', **kwargs)[source]
Derive A_V from every available Balmer decrement.
Ratios each detected Balmer line bluer than the chosen anchor against it, and returns individual and weighted-average A_V values. Only lines at shorter wavelength than the anchor are used: anchor Hβ therefore uses Hγ/Hβ, Hδ/Hβ, H9/Hβ, H10/Hβ (Hα is excluded), while anchor Hα uses Hβ/Hα, Hγ/Hα, Hδ/Hα, H9/Hα, H10/Hα (all of them).
- Parameters:
fluxes (dict) – Emission-line fluxes (observed, not dust-corrected).
errors (dict) – Corresponding 1σ errors.
law (str) –
"salim"(default) or"cardelli".snr_min (float) – Minimum SNR for a Balmer line to be included (default 3.0).
anchor (str) – Reference (denominator) line for the decrement.
"HBETA"(default) uses only bluer lines (Hγ/Hβ, Hδ/Hβ, H9/Hβ, H10/Hβ), excluding Hα."Ha"uses every other Balmer line (Hβ/Hα, Hγ/Hα, Hδ/Hα, H9/Hα, H10/Hα).**kwargs – Extra arguments to the attenuation law (Rv, delta, B).
- Returns:
Keys:
"Av"(weighted mean),"Av_err"(error on mean),"individual"(list of dicts with per-line results),"n_lines"(number of lines used),"anchor"(anchor name).- Return type:
- jwspecabund.dust.compute_lya_escape_fraction(lya_flux, lya_flux_err, fluxes_corr, errors_corr, snr_min=3.0)[source]
Compute Lyα escape fraction from dust-corrected Balmer lines.
Each available Balmer line predicts an intrinsic Lyα flux via Case B recombination (T=10⁴ K, n_e=100 cm⁻³). The escape fraction is the ratio of observed Lyα to intrinsic Lyα.
When multiple Balmer lines are available, f_esc is computed as the inverse-variance weighted mean of the individual estimates.
The observed Lyα flux should NOT be dust-corrected — the escape fraction encapsulates both dust attenuation and resonant scattering.
- Parameters:
lya_flux (float) – Observed (not dust-corrected) Lyα flux.
lya_flux_err (float) – 1σ error on the observed Lyα flux.
fluxes_corr (dict) – Dust-corrected emission-line fluxes keyed by line name.
errors_corr (dict) – Dust-corrected 1σ errors keyed by line name.
snr_min (float) – Minimum SNR for a Balmer line to be included (default 3.0).
- Returns:
Keys:
"f_esc"(weighted mean),"f_esc_err","individual"(list of per-line dicts withline,f_esc,f_esc_err,lya_intrinsic,lya_intrinsic_err),"n_lines"(number used).- Return type:
- jwspecabund.dust.compute_lya_escape_fraction_mc(lya_flux, lya_flux_err, fluxes_corr, errors_corr, Av, Av_err, Av_prior='gaussian', dust_law='salim', snr_min=3.0, n_mc=1000, rng=None, **dust_kwargs)[source]
Compute Lyα escape fraction with MC error propagation.
Draws
n_mcsamples from the observed Lyα flux, dust-corrected Balmer fluxes, and A_V, recomputing the dust correction at each draw. Returns the median and 16th/84th percentile uncertainties.- Parameters:
lya_flux (float) – Observed (not dust-corrected) Lyα flux.
lya_flux_err (float) – 1σ error on Lyα flux.
fluxes_corr (dict) – Dust-corrected Balmer fluxes (used as central values; the MC re-applies dust correction from observed fluxes internally).
errors_corr (dict) – Dust-corrected 1σ errors.
Av (float) – Central A_V value.
Av_err (float) – 1σ error on A_V (Gaussian) or half-width (uniform).
Av_prior (str) – Prior shape for A_V:
"gaussian"(default, draws from N(Av, Av_err) clipped at 0) or"uniform"(draws from U(max(Av − Av_err, 0), Av + Av_err)).dust_law (str) –
"salim"or"cardelli".snr_min (float) – Minimum SNR for a Balmer line to be included (default 3.0).
n_mc (int) – Number of Monte Carlo draws (default 1000).
rng (np.random.Generator or None) – Random number generator (for reproducibility).
**dust_kwargs – Extra arguments to the dust law (Rv, delta, B).
- Returns:
Keys:
"f_esc"(median),"f_esc_err"(lo, hi) as 16th/84th percentile half-widths,"f_esc_posterior"(1D array of MC samples),"individual"(per-line results),"n_lines".- Return type:
- jwspecabund.dust.dust_correct_fluxes(line_fluxes, Av, law='salim', **kwargs)[source]
Apply dust correction to a set of emission-line fluxes.
- Parameters:
- Returns:
{line_name: (corrected_flux, corrected_flux_err)}.- Return type: