From b74ecc066e8cb8ed1221e5469a7c2c27133959b8 Mon Sep 17 00:00:00 2001 From: Will Barnes Date: Tue, 8 Oct 2024 15:48:52 -0400 Subject: [PATCH] free-free fixes (#333) --- fiasco/collections.py | 2 +- fiasco/gaunt.py | 9 ++++----- 2 files changed, 5 insertions(+), 6 deletions(-) diff --git a/fiasco/collections.py b/fiasco/collections.py index cc1f099d..41ca94fd 100644 --- a/fiasco/collections.py +++ b/fiasco/collections.py @@ -115,8 +115,8 @@ def free_free(self, wavelength: u.angstrom): free_free = u.Quantity(np.zeros(self.temperature.shape + wavelength.shape), 'erg cm^3 s^-1 Angstrom^-1') for ion in self: + ff = ion.free_free(wavelength) try: - ff = ion.free_free(wavelength) abundance = ion.abundance ioneq = ion.ioneq except MissingDatasetException as e: diff --git a/fiasco/gaunt.py b/fiasco/gaunt.py index cf918234..60b753ed 100644 --- a/fiasco/gaunt.py +++ b/fiasco/gaunt.py @@ -124,15 +124,14 @@ def _free_free_itoh(self, temperature: u.K, wavelength: u.angstrom, atomic_numbe lower_u = const.h * const.c / const.k_B / tmp upper_u = 1. / 2.5 * (np.log10(lower_u) + 1.5) t = 1. / 1.25 * (log10_temperature - 7.25) - itoh_coefficients = self._itoh['a'][atomic_number-1] + itoh_coefficients = self._itoh['a'][self._itoh['Z']==atomic_number].squeeze() # calculate Gaunt factor gf = u.Quantity(np.zeros(upper_u.shape)) - for j in range(11): - for i in range(11): + for i in range(itoh_coefficients.shape[0]): + for j in range(itoh_coefficients.shape[1]): gf += (itoh_coefficients[i, j] * (t**i))[:, np.newaxis] * (upper_u**j) # apply NaNs where Itoh approximation is not valid - gf = np.where(np.logical_and(np.log10(lower_u) >= -4., np.log10(lower_u) <= 1.0), - gf, np.nan) + gf = np.where(np.logical_and(np.log10(lower_u) >= -4., np.log10(lower_u) <= 1.0), gf, np.nan) gf[np.where(np.logical_or(log10_temperature <= 6.0, log10_temperature >= 8.5)), :] = np.nan return gf