Adds SRM model with unit handling by jajmitchell · Pull Request #227 · sunpy/sunkit-spex

jajmitchell · 2025-08-04T14:31:26Z

SRM model handles units

Adds in the SRM model developed in the PR https://github.com/sunpy/sunkit-spex/pull/163.

This enables physical models to be piped through the SRM model with correct unit handling and then fit to data in count space.

SRM model handles units

Adds changelog

Fixes test_models

Import requires modules in test_models

Fixes examples

Tidy Comments

Merge main

Removes astropy fitting simulated data example

Pre-commit

Renames input and output units as user facing

Fixes Tests

Copilot

Pull request overview

This PR introduces a unit-aware SRM (MatrixModel) so that photon-space models can be propagated through an instrument response into count space with correct unit handling, and updates the simulated-data example and docs accordingly.

Changes:

Replace the old scalar MatrixModel with an SRM-style MatrixModel that knows about input/output axes, bin widths, and spectral units, and update tests to exercise the new behavior.
Rework the fitting_simulated_data example to use astropy.units, SpectralAxis, and Spectrum, piping a physical photon model through the SRM into count space and fitting in count space.
Tidy up documentation examples (RST link targets, math formatting) and pin the Sphinx version for docs builds, plus add a changelog entry.

Reviewed changes

Copilot reviewed 10 out of 10 changed files in this pull request and generated 5 comments.

Show a summary per file

File	Description
`sunkit_spex/models/instrument_response.py`	Replaces the previous minimal `MatrixModel` with a unit-aware SRM implementation that uses input/output spectral axes, bin widths, and a conversion parameter `c` to map photon-space spectra into count-space, and declares input/return units for Astropy modeling.
`sunkit_spex/tests/test_models.py`	Extends the model tests to construct `MatrixModel` with explicit axes and units and to verify that composition with `StraightLineModel` produces the expected result.
`examples/fitting_simulated_data.py`	Refactors the simulated-data fitting example to use `astropy.units`, `quantity_support`, `SpectralAxis`, and `Spectrum`, and to fit a count-space model built from a photon-space model piped through the SRM; updates plotting and summary table to work with quantities.
`examples/fitting_attenuated_RHESSI_spectra.py`	Adjusts the narrative docstring to use an RST-style reference name for the Knuth & Glesener 2020 citation and adds the corresponding hyperlink target.
`examples/fitting_RHESSI_spectra.py`	Cleans up the comparison table comments and converts various values to inline LaTeX math for consistency, while keeping the example’s numerical results unchanged.
`examples/fitting_NuSTAR_spectra-simultaneously.py`	Updates the intro docstring to use an RST reference for Cooper et al. 2021 and adds the link target; leaves the fitting workflow intact.
`examples/fitting_NuSTAR_spectra-glesener2020.py`	Updates documentation-style comparison tables in comments to use inline LaTeX math and an RST reference for Glesener et al. 2020.
`examples/fitting_NuSTAR_spectra-duncan2021.py`	Switches the introduction and comparison comments to use RST references and LaTeX math, and slightly reflows long explanatory comments.
`examples/fitting_attenuated_RHESSI_spectra.py`	Introduces an RST hyperlink target for the Knuth & Glesener 2020 reference in the docstring.
`pyproject.toml`	Pins the Sphinx dependency for docs builds to `<9.0` to avoid incompatibilities with future major releases.
`changelog/227.feature.rst`	Adds a changelog entry announcing that `MatrixModel` now has unit handling.

💡 Add Copilot custom instructions for smarter, more guided reviews. Learn how to get started.

sunkit_spex/models/instrument_response.py

examples/fitting_NuSTAR_spectra-glesener2020.py

examples/fitting_RHESSI_spectra.py

Copilot · 2026-01-30T12:46:30Z

sunkit_spex/models/instrument_response.py

+    def evaluate(self, x, c):
        # Requires input must have a specific dimensionality
-        return model_y @ self.matrix
+
+        input_widths = np.diff(self.input_axis)
+        output_widths = np.diff(self.output_axis)
+
+        flux = ((x * input_widths) @ self.matrix * c) / output_widths
+
+        if hasattr(c, "unit"):
+            return flux
+        return flux.value


The new unit-handling logic in MatrixModel.evaluate (using input_axis/output_axis widths, matrix multiplication, and unit-aware c) is only tested for dimensionless units in test_MatrixModel, and there is no test that covers a fully unitful case like the fitting_simulated_data example (keV, ph/keV, ct/keV). To guard against regressions in the SRM unit handling, it would be valuable to add a test that constructs a MatrixModel with physical units on the axes and spectra, asserts the returned quantity has the expected ct keV^-1 units, and checks the numerical result against a known value.

i agree with the Copilot note

I wrote a similar comment later then found this. I agree too.

Co-authored-by: Copilot <175728472+Copilot@users.noreply.github.com>

jajmitchell · 2026-01-30T15:25:19Z

I have read through the co-pilot suggestions, I think it makes sense to merge as I am going to work immediately on the fitter which will require a different approach of unit handling and testing anyway

settwi · 2026-01-30T18:13:34Z

hey @jajmitchell i'm taking a look through this, i started a review and will request a couple minor changes but i think it looks good overall. i might just push some commits up fixing typos if that is ok

settwi · 2026-01-30T18:14:29Z

does f38afa1 really remove the example?

settwi · 2026-01-30T18:08:32Z

examples/fitting_simulated_data.py


-ph_mod_4astropyfit = StraightLineModel(**guess_cont) + GaussianModel(**guess_line)
-count_model_4astropyfit = (ph_mod_4fit | srm_model) + GaussianModel(**guess_gauss)
+# guess_cont = {"slope": -0.5 * u.ph / u.keV**2, "intercept": 80 * u.ph / u.keV}


@jajmitchell i'm taking a first pass through this PR. i think we need to remove commented blocks if they aren't part of the actual example

settwi · 2026-01-30T18:09:17Z

sunkit_spex/models/instrument_response.py

+
+    c = Parameter(fixed=True)
+
+    def __init__(self, matrix, input_axis, output_axis, input_spec_units, output_spec_units, c):


i think this could benefit from type annotations

settwi · 2026-01-30T18:09:30Z

sunkit_spex/models/instrument_response.py

+        self.input_axis = input_axis
+        self.output_axis = output_axis
+        self.matrix = matrix
+        super().__init__(c)


what is "c" ? i'm surprised Ruff didn't complain about this

I think, I think, c is a stand in for when the matrix model will have response parameters for gain slopes and offsets. The other reason was because Astropy complains when models don't have a fittable parameter. @jajmitchell can correct me on this.

Yea it's the way get the units + model to work with astropy

settwi · 2026-01-30T18:10:13Z

sunkit_spex/models/instrument_response.py

+    _input_units_allow_dimensionless = True
+
+    def evaluate(self, x, c):
        # Requires input must have a specific dimensionality


# Requires input to have a specific dimensionality

sunkit_spex/models/instrument_response.py

settwi · 2026-01-30T19:06:58Z

examples/fitting_simulated_data.py

+    sim_count_model + (2 * np_rand.random(sim_count_model.size)) * np.sqrt(sim_count_model.value) * u.ct / u.keV
+)
+
+obs_spec = Spectrum(sim_count_model_wn.reshape(-1), spectral_axis=ph_energies)


@jajmitchell this is not used after plotting, especially now since the astropy model is commented out.

i think this example needs to get reworked from the ground up to remove stray bits like this. if i was a user and saw this object created then ignored, i would be rather confused. i was actually just confused looking through this example

settwi · 2026-01-30T19:19:49Z

sunkit_spex/models/instrument_response.py

+    def evaluate(self, x, c):
        # Requires input must have a specific dimensionality
-        return model_y @ self.matrix
+
+        input_widths = np.diff(self.input_axis)
+        output_widths = np.diff(self.output_axis)
+
+        flux = ((x * input_widths) @ self.matrix * c) / output_widths
+
+        if hasattr(c, "unit"):
+            return flux
+        return flux.value


i agree with the Copilot note

settwi · 2026-01-30T19:27:54Z

examples/fitting_simulated_data.py

+    plt.plot(ph_energies_centers, (ph_model | srm_model)(ph_energies), label="photon model features")
+    plt.plot(ph_energies_centers, GaussianModel(**sim_gauss)(ph_energies), label="gaussian feature")
+    plt.plot(ph_energies_centers, sim_count_model, label="total sim. spectrum")
+    plt.plot(obs_spec._spectral_axis, obs_spec.data, label="total sim. spectrum + noise", lw=0.5)


These should be plotted using plt.stairs to show that the X-ray spectra are binned

I would go further and say any models should use plt.stairs but the final data plt.plot should use plt.errorbar.

settwi · 2026-01-30T19:32:14Z

examples/fitting_simulated_data.py

-true_vals = np.array(tuple(sim_cont.values())[-2:] + tuple(sim_line.values())[-3:] + tuple(sim_gauss.values())[-3:])
-guess_vals = np.array(
-    tuple(guess_cont.values())[-2:] + tuple(guess_line.values())[-3:] + tuple(guess_gauss.values())[-3:]
+    tuple(sim_cont) + tuple(f"{p}1" for p in tuple(sim_line)) + ("C",) + tuple(f"{p}2" for p in tuple(sim_gauss))


The scipy method does not respect the fixed=True status on the photon/ct scaling parameter

- pass units rather than dict to matrix model - make conversion factor transparent to user (input/output units) - fix quantity_support units in simulated data example

KriSun95 · 2026-01-30T19:57:17Z

examples/fitting_simulated_data.py

+    plt.plot(ph_energies_centers, (ph_model | srm_model)(ph_energies), label="photon model features")
+    plt.plot(ph_energies_centers, GaussianModel(**sim_gauss)(ph_energies), label="gaussian feature")
+    plt.plot(ph_energies_centers, sim_count_model, label="total sim. spectrum")
+    plt.plot(obs_spec._spectral_axis, obs_spec.data, label="total sim. spectrum + noise", lw=0.5)


I would go further and say any models should use plt.stairs but the final data plt.plot should use plt.errorbar.

KriSun95 · 2026-01-30T19:59:20Z

sunkit_spex/models/instrument_response.py

+    def evaluate(self, x, c):
        # Requires input must have a specific dimensionality
-        return model_y @ self.matrix
+
+        input_widths = np.diff(self.input_axis)
+        output_widths = np.diff(self.output_axis)
+
+        flux = ((x * input_widths) @ self.matrix * c) / output_widths
+
+        if hasattr(c, "unit"):
+            return flux
+        return flux.value


I wrote a similar comment later then found this. I agree too.

KriSun95 · 2026-01-30T20:02:19Z

sunkit_spex/models/instrument_response.py

+        self.input_axis = input_axis
+        self.output_axis = output_axis
+        self.matrix = matrix
+        super().__init__(c)


I think, I think, c is a stand in for when the matrix model will have response parameters for gain slopes and offsets. The other reason was because Astropy complains when models don't have a fittable parameter. @jajmitchell can correct me on this.

KriSun95 · 2026-01-30T20:31:51Z

sunkit_spex/models/instrument_response.py

+        input_widths = np.diff(self.input_axis)
+        output_widths = np.diff(self.output_axis)
+
+        flux = ((x * input_widths) @ self.matrix * c) / output_widths


Just looking at the maths here, I think the /output_widths here is a typo. This would likely move the returned flux away from the data units, would it not?

jajmitchell · 2026-02-04T14:50:30Z

@settwi @KriSun95, thanks for the comments, I will make some changes based on the suggestions, however perhaps at this stage its best to remove the fitting_simulated_data.py example, as currently the Astropy models can't be used for fitting so this will all dramatically change, unless you think we should include an example fitting simulated data with the legacy code? If so, that should probably be added in a separate pr and the current example removed in this one.

settwi · 2026-02-04T18:35:07Z

@settwi @KriSun95, thanks for the comments, I will make some changes based on the suggestions, however perhaps at this stage its best to remove the fitting_simulated_data.py example, as currently the Astropy models can't be used for fitting so this will all dramatically change, unless you think we should include an example fitting simulated data with the legacy code? If so, that should probably be added in a separate pr and the current example removed in this one.

Sure, good idea.

Feel free to revert the matrix commit I made btw. I’m not sure if it lines up with the intention of the “c” parameter. Making it transparent (you don’t need to pass it explicitly to the constructor) seemed like a better approach to me

jajmitchell · 2026-02-13T12:39:48Z

@settwi @samaloney @KriSun95 would it make sense to sort this and merge it, or given that the unit handling will completely change within the next week or two with the implementation of the fitter, should we close this PR and add the SRM in the fitting PR, along with the changes to the unit handling of the other models?

samaloney · 2026-02-13T17:12:14Z

I don't have a strong opinion either way, but prob lean towards merging. Also there's now bunch of changes to fix the doc builds so if we don't merge this I'll cherry pick then to a separate PR.

KriSun95 · 2026-02-13T17:14:57Z

Sure! I'm happy to vote to merge then decide on dimension handling stuff later.

When it comes to the fitting_simulated_data.py, can we keep it and just show the fit with Scipy minimize? I know the Astropy fit doesn't work but I still think it's useful to see a full fit from start to finish with what the newer code is capable of just now.

A simulated/custom data fit example for the legacy code already exists in the Fitting Custom Spectra example.

jajmitchell and others added 13 commits August 4, 2025 16:24

SRM model handles units

7a2eec3

SRM model handles units

Adds changelog

947b14f

Adds changelog

Fixes test_models

cda6f81

Fixes test_models

Import requires modules in test_models

5c9eb2f

Import requires modules in test_models

Fixes examples

ea38918

Fixes examples

Tidy comments

0c3357c

Tidy Comments

Merge branch 'main' into srm-model

75a4223

Merge main

Removes astropy fitting simulated data example

f38afa1

Removes astropy fitting simulated data example

Pre-commit

a3005b8

Pre-commit

Renames input and output units as user facing

bcca308

Renames input and output units as user facing

Fixes Tests

fa249b2

Fixes Tests

Fix examples for updates to sphin and docutils

f5e75b9

Pin sphinx<9 so doc build works.

f65c176

samaloney requested a review from Copilot January 30, 2026 12:40

Copilot started reviewing on behalf of samaloney January 30, 2026 12:40 View session

Copilot AI reviewed Jan 30, 2026

View reviewed changes

samaloney and others added 3 commits January 30, 2026 14:20

Disable parallel doc build

8afdd79

Update examples/fitting_NuSTAR_spectra-glesener2020.py

b016d88

Co-authored-by: Copilot <175728472+Copilot@users.noreply.github.com>

Update examples/fitting_RHESSI_spectra.py

8c42844

Co-authored-by: Copilot <175728472+Copilot@users.noreply.github.com>

jajmitchell requested review from KriSun95 and samaloney January 30, 2026 15:25

samaloney approved these changes Jan 30, 2026

View reviewed changes

jajmitchell requested a review from settwi January 30, 2026 15:44

sunpy deleted a comment from Copilot AI Jan 30, 2026

settwi reviewed Jan 30, 2026

View reviewed changes

updates to matrix model

d336fed

- pass units rather than dict to matrix model - make conversion factor transparent to user (input/output units) - fix quantity_support units in simulated data example

KriSun95 requested changes Jan 30, 2026

View reviewed changes

KriSun95 reviewed Jan 30, 2026

View reviewed changes


		c = Parameter(fixed=True)

		def __init__(self, matrix, input_axis, output_axis, input_spec_units, output_spec_units, c):

Comments

Conversation

jajmitchell commented Aug 4, 2025

SRM model handles units

Uh oh!

Copilot AI left a comment

Choose a reason for hiding this comment

Pull request overview

Reviewed changes

Uh oh!

Uh oh!

Uh oh!

Uh oh!

Copilot AI Jan 30, 2026

Choose a reason for hiding this comment

Uh oh!

settwi Jan 30, 2026 • edited Loading Uh oh! There was an error while loading. Please reload this page.

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

jajmitchell commented Jan 30, 2026

Uh oh!

settwi commented Jan 30, 2026

Uh oh!

settwi commented Jan 30, 2026

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Uh oh!

Choose a reason for hiding this comment

Uh oh!

settwi Jan 30, 2026 • edited Loading Uh oh! There was an error while loading. Please reload this page.

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

jajmitchell commented Feb 4, 2026 • edited Loading Uh oh! There was an error while loading. Please reload this page.

Uh oh!

Uh oh!

settwi commented Feb 4, 2026

Uh oh!

jajmitchell commented Feb 13, 2026

Uh oh!

samaloney commented Feb 13, 2026

Uh oh!

KriSun95 commented Feb 13, 2026 • edited Loading Uh oh! There was an error while loading. Please reload this page.

Uh oh!

Uh oh!

Reviewers

Assignees

Labels

Projects

Milestone

settwi Jan 30, 2026 •

edited

Loading

settwi Jan 30, 2026 •

edited

Loading

jajmitchell commented Feb 4, 2026 •

edited

Loading

KriSun95 commented Feb 13, 2026 •

edited

Loading