JP-3547: Flux conservation for spectral resampling

[melanieclarke]

Resolves JP-3547 Resolves JP-3544 Resolves JP-3328

Closes #8297 Closes #8293 Closes #7790

Several fixes associated with spectral flux conservation through resampling. Also added missing documentation for spectral resampling modules.

For NIRSpec:

• create the output WCS by resampling linearly in spatial offset from center rather than in slit units. This makes the output spatial area the same as the input by default.
• fix the output WCS for moving targets in spec3

For MIRI:

• reverse the sense of the scaling when pixel_scale_ratio is provided, so that output_scale = input_scale / pixel_scale_ratio, as documented.
• fix a bug in array sizing when pixel_scale_ratio is not 1.0
• do not apply pixel_scale_ratio to the wavelength array

For both:

• implement support for a directly specified pixel_scale
• scale the output pixel area keywords when pixel_scale or pixel_scale_ratio are set
• if the input data is flux density instead of surface brightness, additionally scale the flux in the input images before resampling if pixel_scale or pixel_scale_ratio are set.

Checklist for PR authors (skip items if you don't have permissions or they are not applicable)

• [x] added entry in CHANGES.rst within the relevant release section
• [x] updated or added relevant tests
• [x] updated relevant documentation
• [x] added relevant milestone
• [x] added relevant label(s)
• [x] ran regression tests, post a link to the Jenkins job below. How to run regression tests on a PR
• [x] All comments are resolved
• [x] Make sure the JIRA ticket is resolved properly

[codecov[bot]]

Codecov Report

Attention: Patch coverage is 95.06726% with 11 lines in your changes missing coverage. Please review.

Project coverage is 60.47%. Comparing base (ebd929c) to head (3e3d718). Report is 2 commits behind head on master.

Files	Patch %	Lines
jwst/resample/resample_spec.py	94.56%	10 Missing :warning:
jwst/resample/resample.py	95.83%	1 Missing :warning:

Additional details and impacted files

```diff @@ Coverage Diff @@ ## master #8596 +/- ## ========================================== + Coverage 60.16% 60.47% +0.30% ========================================== Files 370 369 -1 Lines 38666 38408 -258 ========================================== - Hits 23264 23226 -38 + Misses 15402 15182 -220 ```

:umbrella: View full report in Codecov by Sentry.
:loudspeaker: Have feedback on the report? Share it here.

[melanieclarke]

Started regression tests here: https://plwishmaster.stsci.edu:8081/job/RT/job/JWST-Developers-Pull-Requests/1580/

[melanieclarke]

Regression tests showed a problem with data centering for some spec3 data sets (e.g. in s2d products for test_nirspec_mos_fs_spec3). I pushed a fix and started another run here: https://plwishmaster.stsci.edu:8081/job/RT/job/JWST-Developers-Pull-Requests/1589/

[melanieclarke]

Restarted regression tests here: https://plwishmaster.stsci.edu:8081/job/RT/job/JWST-Developers-Pull-Requests/1606/

There were some unrelated failures for reference file changes in the last run.

[melanieclarke]

I think this is ready for review now. It would be helpful if @nden or @mcara could review from the WCS perspective. @stscirij, could you check on the documentation updates, and see if the note about updating the extraction aperture if the pixel parameters change is sufficient?

[melanieclarke]

Regression test failures are all unrelated (NIRSpec IFU reference file change), or else downstream of NIRSpec MOS or FS resampling.

Inspecting the data on a local run, resampled products show expected differences:

• slight differences in output array shape
• small changes to flux values
• small changes to rectification details (e.g. which rows the source peak falls on, which pixels are NaN due to slit bounding boxes)

Extracted spectra also have slightly different fluxes, as expected: the overall flux level is scaled for most default reductions, due to flux conservation changes. The magnitude of the change depends on the data, but can be up to about 5%. Extraction aperture centering values (EXTRYSTR, EXTRYSTP) have changed from before, but look appropriate to the new resampled arrays.

[melanieclarke]

Restarting regression tests here, to pick up the changes to the moving target spec3 test and the change to a helper function for iscale computation: https://plwishmaster.stsci.edu:8081/job/RT/job/JWST-Developers-Pull-Requests/1614/

[mcara]

Is there anything for spectral data that is equivalent to output_model.meta.photometry.pixelarea_steradians. Shouldn't this be adjusted accordingly if pixel_ratio != 0?

[melanieclarke]

Are changes to build_nirspec_output_wcs() required to fix the flux issue related to pixel_ratio and flux? If not, it would have been easier (for me/reviewer) to review two separate PRs that deal with distinct issues.

Yes, sorry, I know it's a big PR, but the changes to the output WCS were necessary for flux conservation for NIRSpec. The two issues are interrelated: the fixes I implemented for pixel_ratio only work if the default behavior when pixel_ratio=1 preserves input spatial area. The existing version of the spectral WCS did not: it resampled linearly in slit units instead of sky units, so it was implicitly changing the sky area without accounting for it in either the area keywords or the flux values.

[melanieclarke]

Is there anything for spectral data that is equivalent to output_model.meta.photometry.pixelarea_steradians. Shouldn't this be adjusted accordingly if pixel_ratio != 0?

Yes, that's handled at the end of the init function for ResampleSpecData.

[mcara]

Is there anything for spectral data that is equivalent to output_model.meta.photometry.pixelarea_steradians. Shouldn't this be adjusted accordingly if pixel_ratio != 0?

Yes, that's handled at the end of the init function for ResampleSpecData.

~I don't see it.~

Yes, has been added to this PR but it was not in the original (on master) code.

[mcara]

I think I finally figured out the math that justifies why the scaling must be sqrt(pix_tatio) and now I have a better understanding of this. However, according to this newly acquired understanding of the problem, it seems to me that something similar must be done for the spectral coordinate too, unless it is guaranteed that all input images are on the same wavelength grid and this will be the grid used for the output image. For example, if one resamples two non-NIRSPEC input images whose spectral grid is (for example):

image 1: 500.0, 501.0, 502.0, 503.0, 504.0, ...
image 2: 500.5, 501.5, 502.5, 503.5, 504.5, ...

then, the resampled image' grid along the spectral axis will be, thanks to build_interpolated_output_wcs():

output image spectral axis grid: 500.0, 500.5, 501.0, 501.5, 502.0, 502.5, 503.0, 503.5, 504.0, ...

Also see https://github.com/spacetelescope/jwst/blob/ce95ab6babbea9dd5e400b76f3fd3826a72b2921/jwst/resample/resample_spec.py#L91-L99)

So, if data are not in units of spectral density (i.e., MJ/m), the pixel ratio along the spectral axis (0.5 in the above example) should also be taken into account when computing iscale.

At first glance NIRSPEC data are not affected by this unless user-provided output_wcs has a finer-sampled spectral grid. For the case of self.input_models[0].meta.instrument.name != "NIRSPEC", I do not know what the units of the data are. So please verify that they are in density, then nothing needs to be done and if not - then a correction needs to be applied to the iscale.

---

One more issue: reported pixel area as computed/reported in https://github.com/spacetelescope/jwst/blob/6b6f4f2188bb18e3d32c105f59bdb7136e9d488f/jwst/resample/resample_spec.py#L192-L195 may be incorrect. In the resampled image, the grid along the spatial axis is a grid along the slit and it is a 1D thing. There is no "area" here at first glance but in reality, each "pixel" along the spectral axis was integrated in the direction perpendicular to the slit. So, I wonder whether self.blank_output.meta.photometry.pixelarea_steradians should be multiplied by slit width. Of course this is not needed if nominal_area = self.input_models[0].meta.photometry.pixelarea_steradians already takes into account slit width.

[melanieclarke]

So please verify that they are in density, then nothing needs to be done and if not - then a correction needs to be applied to the iscale.

The units are in density, in the spectral axis. Jy is a flux density unit (10^−26 W m−2 Hz−1). Only the spatial area changes need to be accounted for for flux conservation.

[melanieclarke]

So, I wonder whether self.blank_output.meta.photometry.pixelarea_steradians should be multiplied by slit width. Of course this is not needed if nominal_area = self.input_models[0].meta.photometry.pixelarea_steradians already takes into account slit width.

I wrestled with that too, developing this. The nominal_area, set in the photom step, takes into account the slit width, which is not otherwise recorded anywhere in the data model (that I could find).

I initially thought we might need to apply a direct update to the nominal area based on the WCS, the way you did for imaging, but eventually realized it wasn't possible for this reason. The nominal area includes the slit width, the WCS does not. That's why my changes to the resample_spec init function are structured the way they are: pixel_ratio can relatively correct the nominal area, for one of the spatial axes, but we can't compute a whole new spatial area from the WCS.

[melanieclarke]

Thanks so much for your thoughtful review @mcara! I really appreciate you taking the time to work it through with me.