diff --git a/xrtpy/response/effective_area.py b/xrtpy/response/effective_area.py index 4796fa8a7..24f7eb060 100644 --- a/xrtpy/response/effective_area.py +++ b/xrtpy/response/effective_area.py @@ -11,6 +11,7 @@ import astropy.time import numpy as np +import scipy.interpolate import scipy.io import sunpy.io.special import sunpy.time @@ -523,29 +524,34 @@ def _CCD_contamination_transmission(self): return np.array([abs(transmittance[i] ** 2) for i in range(4000)]) @property - def channel_wavelength(self): + def wavelength(self): """Array of wavelengths for every X-ray channel in Angstroms (Å).""" - return Channel(self.name).wavelength + _wave = self._channel.wavelength.to_value("AA") + delta_wave = 0.01 + return np.arange(_wave[0], _wave[-1], delta_wave) * u.Angstrom @property def channel_geometry_aperture_area(self): """XRT flight model geometry aperture area.""" - return Channel(self.name).geometry.geometry_aperture_area + return self._channel.geometry.geometry_aperture_area @property def channel_transmission(self): """XRT channel transmission.""" - return Channel(self.name).transmission + return np.interp( + self.wavelength, self._channel.wavelength, self._channel.transmission + ) + + def _contamination_interpolator(self, x, y): + return np.interp(self.wavelength.to_value("Angstrom"), x, y) @property def _interpolated_CCD_contamination_transmission(self): """Interpolate filter contam transmission to the wavelength.""" - CCD_contam_transmission = np.interp( - self.channel_wavelength.to_value("AA"), + return self._contamination_interpolator( self.n_DEHP_wavelength, self._CCD_contamination_transmission, ) - return CCD_contam_transmission @cached_property def _filter_contamination_transmission(self): @@ -587,12 +593,10 @@ def _filter_contamination_transmission(self): @property def _interpolated_filter_contamination_transmission(self): """Interpolate filter contam transmission to the wavelength.""" - Filter_contam_transmission = np.interp( - self.channel_wavelength.to_value("AA"), + return self._contamination_interpolator( self.n_DEHP_wavelength, self._filter_contamination_transmission, ) - return Filter_contam_transmission @u.quantity_input def effective_area(self) -> u.cm**2: diff --git a/xrtpy/response/temperature_from_filter_ratio.py b/xrtpy/response/temperature_from_filter_ratio.py index a9664966a..9ef52c019 100644 --- a/xrtpy/response/temperature_from_filter_ratio.py +++ b/xrtpy/response/temperature_from_filter_ratio.py @@ -530,7 +530,7 @@ def calculate_TE_errors(map1, map2, T_e, EM, model_ratio, tresp1, tresp2, Trange Narukage's K factor for image 2 """ - wvl = tresp1.channel_wavelength + wvl = tresp1.wavelength eVe = tresp1.ev_per_electron gain = tresp1.ccd_gain_right # (h*c/lambda) * 1/(eV per electron) * 1/gain diff --git a/xrtpy/response/tests/test_effective_area.py b/xrtpy/response/tests/test_effective_area.py index 962fa6ad6..917270c65 100644 --- a/xrtpy/response/tests/test_effective_area.py +++ b/xrtpy/response/tests/test_effective_area.py @@ -105,6 +105,11 @@ def get_IDL_data_files(): return sorted(filter_data_files) +# NOTE: This is marked as xfail because the IDL results that this test compares against +# are incorrect due to the use of quadratic interpolation in the contamination curves +# which leads to ringing near the edges in the contamination curve. +# See https://github.com/HinodeXRT/xrtpy/pull/284#issuecomment-2334503108 +@pytest.mark.xfail @pytest.mark.parametrize("filename", get_IDL_data_files()) def test_effective_area_compare_idl(filename): with Path.open(filename) as f: @@ -112,11 +117,16 @@ def test_effective_area_compare_idl(filename): filter_obs_date = " ".join(f.readline().split()[1:]) # NOTE: Annoyingly the date strings use "Sept" instead of "Sep" for "September" filter_obs_date = filter_obs_date.replace("Sept", "Sep") - IDL_effective_area = np.loadtxt(filename, skiprows=3)[:, 1] * u.cm**2 + IDL_data = np.loadtxt(filename, skiprows=3) + IDL_wavelength = IDL_data[:, 0] * u.AA + IDL_effective_area = IDL_data[:, 1] * u.cm**2 instance = EffectiveAreaFundamental(filter_name, filter_obs_date) actual_effective_area = instance.effective_area() + IDL_effective_area = np.interp( + instance.wavelength, IDL_wavelength, IDL_effective_area + ) assert u.allclose( actual_effective_area, IDL_effective_area, - atol=1e-2 * u.cm**2, + rtol=1e-6, )