From 9e34063799c84ed577fab13e18fd1986ff993d24 Mon Sep 17 00:00:00 2001 From: Nabil Freij Date: Tue, 20 Oct 2026 09:21:25 -0700 Subject: [PATCH 1/2] Fold the temperature response on the spectra wavelengths Interpolating the spectra onto the channel grid lost or multiplied the flux of lines narrower than the channel sampling. --- changelog/231.bugfix.rst | 1 + sunkit_instruments/response/tests/test_channel.py | 13 +++++++++++++ sunkit_instruments/response/thermal.py | 13 +++++++------ 3 files changed, 21 insertions(+), 6 deletions(-) create mode 100644 changelog/231.bugfix.rst diff --git a/changelog/231.bugfix.rst b/changelog/231.bugfix.rst new file mode 100644 index 00000000..4ccc53d4 --- /dev/null +++ b/changelog/231.bugfix.rst @@ -0,0 +1 @@ +`~sunkit_instruments.response.SourceSpectra.temperature_response` now interpolates the wavelength response onto the wavelengths of the spectra, instead of the spectra onto the channel wavelengths, so that spectra sampled more finely than the channel keep their flux. diff --git a/sunkit_instruments/response/tests/test_channel.py b/sunkit_instruments/response/tests/test_channel.py index 69a95f80..3bff52d2 100644 --- a/sunkit_instruments/response/tests/test_channel.py +++ b/sunkit_instruments/response/tests/test_channel.py @@ -90,3 +90,16 @@ def test_temperature_response(fake_channel, fake_spectra, obstime): temp_response = fake_spectra.temperature_response(fake_channel, obstime=obstime) assert isinstance(temp_response, u.Quantity) assert temp_response.shape == fake_spectra.temperature.shape + + +@pytest.mark.parametrize("offset", [0, 0.5]) +def test_temperature_response_line(fake_channel, offset): + # A unit-integral line on (0) or between (0.5) channel wavelengths gives the response at the line + line = fake_channel.wavelength[50] + offset * np.diff(fake_channel.wavelength)[50] + delta = 0.02 * u.AA + wavelength = line + np.arange(-100, 101) * delta + data = np.zeros((1, wavelength.size)) * u.Unit("photon cm3 s-1 sr-1 Angstrom-1") + data[0, 100] = u.Quantity(1, "photon cm3 s-1 sr-1") / delta + spectra = SourceSpectra([1e6] * u.K, wavelength, data) + expected = np.interp(line, fake_channel.wavelength, fake_channel.wavelength_response()) + assert u.allclose(spectra.temperature_response(fake_channel), expected * u.Quantity(1, "photon cm3 s-1 sr-1")) diff --git a/sunkit_instruments/response/thermal.py b/sunkit_instruments/response/thermal.py index 7de91ee5..209697bd 100644 --- a/sunkit_instruments/response/thermal.py +++ b/sunkit_instruments/response/thermal.py @@ -165,13 +165,14 @@ def temperature_response( attrs={'unit': 'Angstrom'}),], attrs={'unit': wave_response.unit.to_string()} ) - spec_interp = self._da.interp( - wavelength=da_response.wavelength, + # Integrate on the spectra wavelengths so that narrow lines keep their flux + response_interp = da_response.interp( + wavelength=self._da.wavelength, kwargs={'bounds_error': False, 'fill_value': 0.0}, ) - spec_interp_weighted = spec_interp * da_response - temp_response = spec_interp_weighted.integrate(coord='wavelength') - final_unit = (u.Unit(spec_interp.unit) + spec_weighted = self._da * response_interp + temp_response = spec_weighted.integrate(coord='wavelength') + final_unit = (u.Unit(self._da.unit) * u.Unit(da_response.attrs['unit']) - * u.Unit(spec_interp_weighted.wavelength.unit)) + * u.Unit(self._da.wavelength.unit)) return u.Quantity(temp_response.data, final_unit) From 2e071d5b267be288a15ec1234ca709355e1c4138 Mon Sep 17 00:00:00 2001 From: Nabil Freij Date: Tue, 20 Oct 2026 10:21:13 -0700 Subject: [PATCH 2/2] Integrate in ascending wavelength order --- sunkit_instruments/response/tests/test_channel.py | 5 +++-- sunkit_instruments/response/thermal.py | 2 +- 2 files changed, 4 insertions(+), 3 deletions(-) diff --git a/sunkit_instruments/response/tests/test_channel.py b/sunkit_instruments/response/tests/test_channel.py index 3bff52d2..9dc3aaa7 100644 --- a/sunkit_instruments/response/tests/test_channel.py +++ b/sunkit_instruments/response/tests/test_channel.py @@ -100,6 +100,7 @@ def test_temperature_response_line(fake_channel, offset): wavelength = line + np.arange(-100, 101) * delta data = np.zeros((1, wavelength.size)) * u.Unit("photon cm3 s-1 sr-1 Angstrom-1") data[0, 100] = u.Quantity(1, "photon cm3 s-1 sr-1") / delta - spectra = SourceSpectra([1e6] * u.K, wavelength, data) expected = np.interp(line, fake_channel.wavelength, fake_channel.wavelength_response()) - assert u.allclose(spectra.temperature_response(fake_channel), expected * u.Quantity(1, "photon cm3 s-1 sr-1")) + for order in [slice(None), slice(None, None, -1)]: + spectra = SourceSpectra([1e6] * u.K, wavelength[order], data[:, order]) + assert u.allclose(spectra.temperature_response(fake_channel), expected * u.Quantity(1, "photon cm3 s-1 sr-1")) diff --git a/sunkit_instruments/response/thermal.py b/sunkit_instruments/response/thermal.py index 209697bd..8e8e0b89 100644 --- a/sunkit_instruments/response/thermal.py +++ b/sunkit_instruments/response/thermal.py @@ -170,7 +170,7 @@ def temperature_response( wavelength=self._da.wavelength, kwargs={'bounds_error': False, 'fill_value': 0.0}, ) - spec_weighted = self._da * response_interp + spec_weighted = (self._da * response_interp).sortby('wavelength') temp_response = spec_weighted.integrate(coord='wavelength') final_unit = (u.Unit(self._da.unit) * u.Unit(da_response.attrs['unit'])