JP-3848 MIRI LRS s_region and resample WCS #9193

jemorrison · 2025-02-14T15:18:47Z

This PR replaces: #9104

NOT TO BE MERGED until change to MIRI LRS specwcs reference file is made:spacetelescope/stdatamodels#393

This PR corrects the s_region for MIRI LRS output from assign_wcs and resample. In addition an error was found in the WCS of resampled MIRI LRS data.
In assign_wcs the s_region is determined using the V2,V3 corners of the slit. These values are provided by a new
reference file.

Tasks

news fragment change types...

changes/<PR#>.general.rst: infrastructure or miscellaneous change
changes/<PR#>.docs.rst
changes/<PR#>.stpipe.rst
changes/<PR#>.datamodels.rst
changes/<PR#>.scripts.rst
changes/<PR#>.set_telescope_pointing.rst
changes/<PR#>.pipeline.rst

stage 1

changes/<PR#>.group_scale.rst
changes/<PR#>.dq_init.rst
changes/<PR#>.emicorr.rst
changes/<PR#>.saturation.rst
changes/<PR#>.ipc.rst
changes/<PR#>.firstframe.rst
changes/<PR#>.lastframe.rst
changes/<PR#>.reset.rst
changes/<PR#>.superbias.rst
changes/<PR#>.refpix.rst
changes/<PR#>.linearity.rst
changes/<PR#>.rscd.rst
changes/<PR#>.persistence.rst
changes/<PR#>.dark_current.rst
changes/<PR#>.charge_migration.rst
changes/<PR#>.jump.rst
changes/<PR#>.clean_flicker_noise.rst
changes/<PR#>.ramp_fitting.rst
changes/<PR#>.gain_scale.rst

stage 2

changes/<PR#>.assign_wcs.rst
changes/<PR#>.badpix_selfcal.rst
changes/<PR#>.msaflagopen.rst
changes/<PR#>.nsclean.rst
changes/<PR#>.imprint.rst
changes/<PR#>.background.rst
changes/<PR#>.extract_2d.rst
changes/<PR#>.master_background.rst
changes/<PR#>.wavecorr.rst
changes/<PR#>.srctype.rst
changes/<PR#>.straylight.rst
changes/<PR#>.wfss_contam.rst
changes/<PR#>.flatfield.rst
changes/<PR#>.fringe.rst
changes/<PR#>.pathloss.rst
changes/<PR#>.barshadow.rst
changes/<PR#>.photom.rst
changes/<PR#>.pixel_replace.rst
changes/<PR#>.resample_spec.rst
changes/<PR#>.residual_fringe.rst
changes/<PR#>.cube_build.rst
changes/<PR#>.extract_1d.rst
changes/<PR#>.resample.rst

stage 3

changes/<PR#>.assign_mtwcs.rst
changes/<PR#>.mrs_imatch.rst
changes/<PR#>.tweakreg.rst
changes/<PR#>.skymatch.rst
changes/<PR#>.exp_to_source.rst
changes/<PR#>.outlier_detection.rst
changes/<PR#>.tso_photometry.rst
changes/<PR#>.stack_refs.rst
changes/<PR#>.align_refs.rst
changes/<PR#>.klip.rst
changes/<PR#>.spectral_leak.rst
changes/<PR#>.source_catalog.rst
changes/<PR#>.combine_1d.rst
changes/<PR#>.ami.rst

other

changes/<PR#>.wfs_combine.rst
changes/<PR#>.white_light.rst
changes/<PR#>.cube_skymatch.rst
changes/<PR#>.engdb_tools.rst
changes/<PR#>.guider_cds.rst

codecov · 2025-02-14T15:29:30Z

Codecov Report

Attention: Patch coverage is 99.15254% with 1 line in your changes missing coverage. Please review.

Project coverage is 73.72%. Comparing base (16aa8b0) to head (7c4c40f).
Report is 3 commits behind head on main.

Files with missing lines	Patch %	Lines
jwst/resample/resample_spec.py	94.11%	1 Missing ⚠️

Additional details and impacted files

@@            Coverage Diff             @@
##             main    #9193      +/-   ##
==========================================
+ Coverage   73.66%   73.72%   +0.05%     
==========================================
  Files         368      368              
  Lines       36392    36464      +72     
==========================================
+ Hits        26808    26882      +74     
+ Misses       9584     9582       -2

☔ View full report in Codecov by Sentry.
📢 Have feedback on the report? Share it here.

emolter · 2025-02-14T16:03:33Z

pyproject.toml

@@ -34,7 +34,8 @@ dependencies = [
    "scikit-image>=0.20.0",
    "scipy>=1.14.1",
    "spherical-geometry>=1.2.22",
-    "stdatamodels @ git+https://github.com/spacetelescope/stdatamodels.git@main",
+    #"stdatamodels @ git+https://github.com/spacetelescope/stdatamodels.git@main",
+    "stdatamodels @ stdatamodels-2.2.1.dev20+gf85885d",


Suggested change

"stdatamodels @ stdatamodels-2.2.1.dev20+gf85885d",

"stdatamodels @ git+https://github.com/jemorrison/stdatamodels.git@MIRI_LRS_Specwcs",

drlaw1558 · 2025-02-17T14:21:27Z

jwst/resample/resample_spec.py


        # single model: use size of x_tan_array
        # to be consistent with method before
        if len(input_models) == 1:
-            x_size = int(np.ceil(xstop))
+            ny = int(np.ceil(xstop))


I think we want to remove this IF statement; see general comments.

Made this change and it does make the s_regions look better

drlaw1558 · 2025-02-17T14:28:48Z

@jemorrison This is looking much better, but a few thoughts:

The across-slice positioning now looks good, both for S_REGIONs and for individual pixels. The along-slice is still a little weird though; from poking at this I think it's because the rotation that's being applied can change what the 'edge' pixels are slightly between cal and s2d frames (actual SOURCE positions look ok). One thing that I think we should fix though- in the individual s2d files it looks like the right edge of the resampled array is getting cut off; if I remove the ny = int(np.ceil(xstop)) and corresponding IF statement (see comment above) though then this looks better.
For some reason the pipeline log says again that there are NaNs in the bounding box during the assign_wcs step and that the S_REGIONS are not being updated. As a result I can't tell whether the new S_REGION code there is doing the right thing or if it's simply not changing the S_REGION from the rate files again. I'd thought that the latest pipeline was running this recalculation though, which is what kicked off much of this work in the first place?
I've got the new reference file ready to go as soon as the datamodels PR is merged.

jemorrison · 2025-02-17T14:40:34Z

@drlaw1558 are you supplying the new reference file when you run assign_wcs (spec2) in the override method ?

drlaw1558 · 2025-02-17T14:48:11Z

@jemorrison Hm, I thought I had been but it looks like the new ref file wasn't getting picked up properly. With that fixed the assign_wcs indeed runs and results look good.

drlaw1558 · 2025-02-17T15:30:28Z

LG2M now:

drlaw1558

Looks good from the science side.

emolter · 2025-02-17T15:40:10Z

changes/9193.resample.rst

@@ -0,0 +1 @@
+Fix MIRI LRS s_region and WCS  in resample_spec


Suggested change

Fix MIRI LRS s_region and WCS in resample_spec

Fix MIRI LRS s_region and WCS in resample_spec

nitpick, but I think there's an accidental double-space here

changes/9193.resample.rst

emolter · 2025-02-17T15:44:23Z

jwst/assign_wcs/assign_wcs.py

@@ -72,8 +73,13 @@ def load_wcs(input_model, reference_files={}, nrs_slit_y_range=None):

        if output_model.meta.exposure.type.lower() not in exclude_types:
            imaging_types = IMAGING_TYPES.copy()
-            imaging_types.update(['mir_lrs-fixedslit', 'mir_lrs-slitless'])
-            if output_model.meta.exposure.type.lower() in imaging_types:
+            imaging_lrs_types = ['mir_lrs-fixedslit'] # keep it separate


Suggested change

imaging_lrs_types = ['mir_lrs-fixedslit'] # keep it separate

imaging_lrs_types = ['mir_lrs-fixedslit']

once the PR has been merged, this comment seems like it might be more confusing than helpful

emolter · 2025-02-17T15:45:51Z

jwst/assign_wcs/assign_wcs.py

+                # uses slits corners in V2, V3 that are read in from the
+                # lrs specwcs reference file


can this information be placed in the docstring of update_s_region_lrs()?

Do you mean something different than:

Determine S_REGION using V2,V3 of the slit corners from reference file

no that's fine if it's already there, I just didn't find the comment here to be particularly helpful, in the sense that someone who wanted to know how the s_region for LRS was computed would look in the docstring of the function dedicated to that. i.e. IMO the comment here could be removed

emolter · 2025-02-17T15:47:29Z

jwst/assign_wcs/miri.py

+    # remove commented out how after testing
+    #if input_model.meta.exposure.type.lower() == 'mir_lrs-fixedslit':
+    #    bb_sub = (bb_sub[0], (dxmodel.points[0].min(), dxmodel.points[0].max()))


don't forget to remove this before merge

thanks I had forgotten

emolter · 2025-02-17T15:50:24Z

jwst/assign_wcs/util.py

+                          [a[3], b[3]]])
+
+    update_s_region_keyword(model, footprint)
+    print('in  assign_wcs util.py',model.meta.wcsinfo.s_region)


remove print

emolter · 2025-02-17T15:50:44Z

jwst/assign_wcs/util.py

+        log.info("The s_region will not be updated")
+
+    lam = 7.0 # wavelength does not matter for s region assign a value in
+              # wavelength of MIRI LRS


keep on a single line if possible

emolter · 2025-02-17T16:00:08Z

jwst/resample/resample_spec.py


        # turn the size into a numpy shape in (y, x) order
        output_wcs.array_shape = output_array_size[::-1]
        output_wcs.pixel_shape = output_array_size
        bounding_box = wcs_bbox_from_shape(output_array_size[::-1])
        output_wcs.bounding_box = bounding_box

-        return output_wcs
+        s_region = find_miri_lrs_sregion(sregion_list,output_wcs)


is MIRI LRS the only mode that uses the build_interpolated_output_wcs function? It looks like this is being applied to any data that ends up here. If it applies to more than just MIRI LRS, I suggest changing the name of this function to be more generic. If build_interpolated_output_wcs only applies to MIRI LRS, then I might suggest adding a line to the docstring of that function specifying that.

@melanieclarke I wanted to double check that build_interpolated_output_wcs - is only applied to MIRI LRS data.

It is currently used only for MIRI LRS FS but it's written to be more generally applicable, to any instrument other than NIRSpec. There is internal handling for dispersion directions that MIRI doesn't use, for example. It might be worthwhile to keep the implementation general for these new changes, in case we wanted to use it to try resampling for other spectral modes.

emolter · 2025-02-17T16:17:46Z

jwst/resample/resample_spec.py

@@ -1011,3 +1025,85 @@ def compute_spectral_pixel_scale(wcs, fiducial=None, disp_axis=1):

    pixel_scale = compute_scale(wcs, fiducial, disp_axis=disp_axis)
    return float(pixel_scale)
+
+def find_miri_lrs_sregion(sregion_list,wcs):


This appears to only need sregion_list[0] and not all the s_regions from all the models. Is there a way to leverage that fact to move the call to this function closer to the end of the step, so that we don't need to add a _s_region attribute to ResampleSpec?

I can check if I'm == 0 and only set the s_region to the one case. Is that what you would like ?

I guess it doesn't really matter, I just thought it might be cleaner if it received only a single reference model.

emolter · 2025-02-17T16:18:47Z

jwst/resample/resample_spec.py

I'm a bit concerned about the increasing spaghettification of the code with respect to how the s_region is handled for different modes, as follows:

This PR adds an _s_region attribute to the ResampleSpec class but only updates it for MIRI LRS data. Then the output gets the s_region applied to it inside the update_output_model method. This lives within an if _s_region is not None to apply it only to the appropriate data.

Other spectroscopic data have their s_region computed on this line basically right at the end of ResampleSpecStep in _process_slit() and not within any kind of update_output_model or similar function.

Imaging data get their s_region on this line inside the finalize method, which is also accessible to, and called by, ResampleSpec. This, too, needs an if statement around it to ensure it only gets applied to imaging data.

Is there any way to put these three in the same place so the code is easier to understand?

I know they're all computed slightly differently, but is there really no way to compute the LRS s_region closer to the end of the step? I noticed that the computation really only needs the zeroth s_region, not an entire list, perhaps that helps achieve this?

Well I need this first input model s_region to determine the slit width and I need to output wcs also.
Would you prefer that I move the 'find_miri_lrs_region' to resample_spec_step ?

I think yes, see https://github.com/spacetelescope/jwst/pull/9193/files#r1958823896

melanieclarke · 2025-02-17T22:15:51Z

@jemorrison - there are some unit test failures that need attention when the code structure has settled. In particular, it looks like some of the tests I wrote to make sure flux is conserved are failing with the updates to the MIRI WCS.

emolter · 2025-02-17T22:29:49Z

jwst/resample/resample_spec_step.py

        self.update_slit_metadata(result)
-        update_s_region_spectral(result)
+        if result.meta.exposure.type.lower() != 'mir_lrs-fixedslit':
+            update_s_region_spectral(result)



Per my other comment, I assume the reason you put this if statement here is because the LRS data in question also follow this path. Could you put the call to find_miri_lrs_sregion() here? It appears to me that all the necessary inputs are already here, and it would make it so all spectral data get their s_region updated in the same place. Plus you wouldn't have to have ResampleSpec._s_region

e.g.,

if result.meta.exposure.type.lower() == 'mir_lrs-fixedslit': result.meta.wcsinfo.s_region = find_miri_lrs_sregion(input_models[0], result.meta.wcs) else: update_s_region_spectral(result)

…cwcs model

emolter

This looks much cleaner now from an overall code design standpoint - this is exactly the flow I had in mind. I made a few additional small comments, but overall this looks great and I'm happy for it to be merged

emolter · 2025-02-25T17:46:16Z

jwst/assign_wcs/miri.py

+        zero_point = refmodel.meta.x_ref_slitless - 1, \
+            refmodel.meta.y_ref_slitless - 1


Does this fit on one line?

emolter · 2025-02-25T17:47:17Z

jwst/assign_wcs/util.py

+
+    Returns
+    -------
+    None
+        s_region for model is updated in place .
+


Suggested change

Returns

-------

None

s_region for model is updated in place .

if you think it's useful to specify the update is done in-place, it can be put in the "extended summary"

emolter · 2025-02-25T17:48:18Z

jwst/assign_wcs/util.py

+        log.info("The V2,V3 coordinates of the MIRI LRS-Fixed slit contains NaN values.")
+        log.info("The s_region will not be updated")
+
+    lam = 7.0 # wavelength does not matter for s region assign a value in range of LRS


Suggested change

lam = 7.0 # wavelength does not matter for s region assign a value in range of LRS

lam = 7.0 # wavelength does not matter for s_region so just assign a value in range of LRS

emolter · 2025-02-25T17:51:40Z

jwst/resample/resample_utils.py

+        s_region for the resample data.
+    """
+    # use the first sregion to set the width of the slit
+


Suggested change

emolter · 2025-02-25T17:52:55Z

jwst/resample/tests/test_resample_step.py

@@ -482,6 +489,7 @@ def test_pixel_scale_ratio_spec_miri(miri_cal, ratio, units):
    assert result1.data.shape[0] == result2.data.shape[0]

    # spatial dimension is scaled
+    print(result1.data.shape[1], result2.data.shape[1])


don't forget to remove print statement before merge

emolter · 2025-02-25T17:54:51Z

jwst/resample/tests/test_resample_step.py

I know you already have an s_region test in the regtest, but would it be straightforward to write a unit test for just the find_miri_lrs_sregion function on its own?

jemorrison requested review from a team as code owners February 14, 2025 15:18

github-actions bot added assign_wcs resample MIRI spec2 pipeline regression_testing labels Feb 14, 2025

jemorrison added this to the Build 11.3 milestone Feb 14, 2025

github-actions bot added the installation label Feb 14, 2025

jemorrison mentioned this pull request Feb 14, 2025

JP-3848 Update S_REGION for MIRI LRS fixed slit data #9104

Closed

10 tasks

emolter reviewed Feb 14, 2025

View reviewed changes

jemorrison requested review from nden, drdavella and drlaw1558 and removed request for drdavella February 14, 2025 18:10

drlaw1558 reviewed Feb 17, 2025

View reviewed changes

jemorrison force-pushed the JP-3848_MIRI_LRS_Second branch from 1eae471 to 1743b60 Compare February 17, 2025 15:19

drlaw1558 approved these changes Feb 17, 2025

View reviewed changes

emolter reviewed Feb 17, 2025

View reviewed changes

braingram mentioned this pull request Feb 19, 2025

Add miri lrs spec wcs model spacetelescope/stdatamodels#393

Open

5 tasks

github-actions bot added the testing label Feb 20, 2025

jemorrison added 2 commits February 25, 2025 10:26

updates to MIRI LRS s_region and resample WCS

4937b83

changed pyproject.toml to point to stdatamodels with new MIRI LRS spe…

936df7c

…cwcs model

jemorrison added 8 commits February 25, 2025 10:26

added change log

273dc22

update location of stdatamodels to use

fa6924c

fix ny of s2d region

4415b41

update from review

882361c

update to where resampled s_region determined

466045f

working on unit tests

f81e583

fix unit tests

eb221c8

remove comment

7c4c40f

jemorrison force-pushed the JP-3848_MIRI_LRS_Second branch from 3280f74 to 7c4c40f Compare February 25, 2025 17:26

emolter reviewed Feb 25, 2025

View reviewed changes

	"stdatamodels @ stdatamodels-2.2.1.dev20+gf85885d",
	"stdatamodels @ git+https://github.com/jemorrison/stdatamodels.git@MIRI_LRS_Specwcs",

	Fix MIRI LRS s_region and WCS in resample_spec
	Fix MIRI LRS s_region and WCS in resample_spec

	imaging_lrs_types = ['mir_lrs-fixedslit'] # keep it separate
	imaging_lrs_types = ['mir_lrs-fixedslit']

		# uses slits corners in V2, V3 that are read in from the
		# lrs specwcs reference file

		zero_point = refmodel.meta.x_ref_slitless - 1, \
		refmodel.meta.y_ref_slitless - 1

	lam = 7.0 # wavelength does not matter for s region assign a value in range of LRS
	lam = 7.0 # wavelength does not matter for s_region so just assign a value in range of LRS

JP-3848 MIRI LRS s_region and resample WCS #9193

Are you sure you want to change the base?

JP-3848 MIRI LRS s_region and resample WCS #9193

Conversation

jemorrison commented Feb 14, 2025 • edited Loading

Tasks

stage 1

stage 2

stage 3

other

codecov bot commented Feb 14, 2025 • edited Loading

Codecov Report

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

drlaw1558 commented Feb 17, 2025 • edited Loading

jemorrison commented Feb 17, 2025

drlaw1558 commented Feb 17, 2025

drlaw1558 commented Feb 17, 2025 • edited Loading

drlaw1558 left a comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

emolter Feb 17, 2025 • edited Loading

Choose a reason for hiding this comment

melanieclarke commented Feb 17, 2025

emolter Feb 17, 2025 • edited Loading

Choose a reason for hiding this comment

emolter left a comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

jemorrison commented Feb 14, 2025 •

edited

Loading

codecov bot commented Feb 14, 2025 •

edited

Loading

drlaw1558 commented Feb 17, 2025 •

edited

Loading

drlaw1558 commented Feb 17, 2025 •

edited

Loading

emolter Feb 17, 2025 •

edited

Loading

emolter Feb 17, 2025 •

edited

Loading