cometspec.rates

Rate-matrix construction, level-population solver, g-factors, spectrum synthesis, and line-spread-function builders.

Radiative rate matrix, steady-state solver, g-factors, spectrum synthesis, and line-spread-function helpers.

Routines

cometspec.rates.build_rate_matrix_nbar(lines, *, include_stim_emission=False, verbose=True, A_col='A_ul', upper_id_col='upper_id', lower_id_col='lower_id', g_upper_col='g_upper', g_lower_col='g_lower')[source]

Build the radiative rate matrix from normalized transitions.

Parameters:

  • lines (astropy.table.Table) – Transition table with solar incident irradiance quantities attached.

  • include_stim_emission (bool, optional, default False) – Include stimulated emission in downward rates.

  • verbose (bool, optional, default True) – Print diagnostic information.

  • A_col (str, optional, default "A_ul") – Einstein A column name.

  • upper_id_col (str, optional, default "upper_id") – Upper-level ID column name.

  • lower_id_col (str, optional, default "lower_id") – Lower-level ID column name.

  • g_upper_col (str, optional, default "g_upper") – Upper-level degeneracy column name.

  • g_lower_col (str, optional, default "g_lower") – Lower-level degeneracy column name.

Returns:

tuple[numpy.ndarray, dict[int, str], astropy.table.Table] – Rate matrix, a dictionary that maps each row/column index of M back to its physical level identifier (the string from upper_id_col / lower_id_col), and the input line lists with new columns (__nu_Hz, __R_lu, __R_ul, __upper_idx, __lower_idx).

cometspec.rates.solve_with_normalization(M, *, verbose=True)[source]

Solve the steady-state system Mn = 0 with sum(n) = 1.

Parameters:

  • M (numpy.ndarray) – Rate matrix.

  • verbose (bool, optional, default True) – Print solver diagnostics.

Returns:

numpy.ndarray – Normalized non-negative level populations.

Raises:

RuntimeError – If no positive normalized solution can be formed.

Notes

If np.linalg.lstsq fails to converge (e.g. because M contains NaN or Inf from extreme parameter values), a LinAlgError is caught and an array of NaN is returned. In the MCMC context this propagates to a non-finite model flux so the walker step is rejected (log-likelihood → -inf) rather than crashing. In the direct-modelling context best_model will contain NaN, signalling that the chosen parameters are unphysical.

cometspec.rates.solve_with_normalization_fast(M, A_work, b_work)[source]

Fast buffer-reusing version of solve_with_normalization().

Parameters:
Returns:

numpy.ndarray – Normalized non-negative level populations.

Raises:

RuntimeError – If no positive normalized solution can be formed.

Notes

The first row of M is overwritten with ones (and b[0] = 1) to impose sum(n) = 1 directly, so the resulting square system is non-singular for any well-posed rate matrix and is solved with scipy.linalg.solve() (LU, with overwrite_a/overwrite_b and check_finite=False to skip the input scan and reuse the work buffers in-place — several times faster than the SVD-based numpy.linalg.lstsq() used by solve_with_normalization()). If the matrix is singular (e.g. NaN/Inf collision rates at extreme parameter values), the raised LinAlgError is caught and an array of NaN is returned. In the MCMC context this propagates to a non-finite model flux so the walker step is rejected (log-likelihood -> -inf) rather than crashing. In the direct-modelling context best_model will contain NaN, signalling that the chosen parameters are unphysical.

cometspec.rates.g_factors(lines_with_rates, n, *, A_col='A_ul')[source]

Compute per-line and total photon/energy g-factors.

The g-factors are in units of photons/s/molecule for g_phot and erg/s/molecule for g_energy.

Parameters:

  • lines_with_rates (astropy.table.Table) – Transition table with cached rate-matrix indices.

  • n (numpy.ndarray) – Level populations.

  • A_col (str, optional, default "A_ul") – Einstein A column name.

Returns:

tuple[numpy.ndarray, numpy.ndarray, float, float] – (g_phot, g_energy, sum_g_phot, sum_g_energy).

cometspec.rates.g_factors_fast_from_cache(*, ui, A_ul, hnu, n, out_g_ph=None, out_g_en=None)[source]

Fast g-factor evaluation using precomputed arrays.

Parameters:
Returns:

tuple[numpy.ndarray, numpy.ndarray] – (g_ph, g_en) arrays.

cometspec.rates.synth_spectrum_from_lines(df_lines, N_col_cm2, *, g_line_energy=None, g_line_phot=None, fwhm_A=0.02, dlam_A=0.05, lam_min=None, lam_max=None, lam_col='Wave_vac_AA', Omega_sr=None, grid=None, lsf=None, v_shift_kms=0.0, dlam_shift_A=0.0)[source]

Build a synthetic emission spectrum from line g-factors.

Parameters:

  • df_lines (astropy.table.Table) – Transition table. Required

  • N_col_cm2 (float) – Column density in cm^-2. Required

  • g_line_energy (numpy.ndarray, optional, default None) – Per-line energy g-factors. If not provided g_line_phot must be used to compute line intensities.

  • g_line_phot (numpy.ndarray, optional, default None) – Per-line photon g-factors. If not provided, g_line_energy must be used to compute line intensities.

  • fwhm_A (float, optional, default 0.02) – Gaussian FWHM in Angstrom if no custom LSF is provided.

  • dlam_A (float, optional, default 0.05) – Wavelength step for an auto-generated grid.

  • lam_min (float, optional, default None) – Optional minimum wavelength (after wavelength shifts) for an auto-generated grid. The minimum grid value will be at least 5 FWHM below the minimum line wavelength.

  • lam_max (float, optional, default None) – Optional maximum wavelength (after wavelength shifts) for an auto-generated grid. The maximum grid value will be at least 5 FWHM above the maximum line wavelength.

  • lam_col (str, optional, default "Wave_vac_AA" If None, it will be set to "Wave_vac_AA") – Wavelength column in df_lines.

  • Omega_sr (float, optional, default None) – Optional solid angle scaling in sr.

  • grid (numpy.ndarray, optional, default None) – Optional output wavelength grid.

  • lsf (Callable[[numpy.ndarray], numpy.ndarray], optional, default None) – Optional custom line-spread function.

  • v_shift_kms (float, optional, default 0.0) – Velocity shift in km/s for the synthetic spectrum lines.

  • dlam_shift_A (float, optional, default 0.0) – Additive wavelength shift in Angstrom for the synthetic spectrum lines.

Returns:

tuple[numpy.ndarray, numpy.ndarray] – (wavelength_grid, synthetic_flux).

Raises:

ValueError – If required inputs are missing.

cometspec.rates.make_lsf(params, mode)[source]

Create a line-spread function callable from parameter values.

  • If LSF mode == “2Gauss”, the parameters are sigma1, sigma2, and ratio.

  • If LSF mode == “Gauss”, the parameter is sigma.

  • If LSF mode == “Gauss_Lorentz”, the parameters are sigma_G, fwhm_L, and ratio.

  • If LSF mode == “Lorentz”, the parameter is fwhm_L.

Parameters:

  • params (dict[str, float]) – LSF parameter dictionary.

  • mode (str) – LSF mode: 2Gauss, Gauss, Gauss_Lorentz, or Lorentz.

Returns:

Callable[[numpy.ndarray], numpy.ndarray] – LSF callable. Returned as a small picklable class instance (not a closure), so it can be sent to spawn-based multiprocessing workers when used as the lsf argument to cometspec.mcmc.mcmc_fitting().

Raises:

ValueError – If mode is invalid.