JP-3739: Update NRM geometry, affine distortion use (#8974)

changes/8974.ami.rst

@@ -0,0 +1 @@
+Use mask and pupil geometry constants from NRM reference file, and apply affine distortion from commissioning to LG model as default.

docs/jwst/ami_analyze/description.rst

@@ -47,11 +47,14 @@ other options:
 :--chooseholes: If not None, fit only certain fringes e.g. ['B4','B5','B6','C2']
                 (default=None)
 
-:--affine2d: ASDF file containing user-defined affine parameters (default=None)
+:--affine2d: ASDF file containing user-defined affine parameters (default='commissioning')
 
 :--run_bpfix: Run Fourier bad pixel fix on cropped data (default=True)
 
 
+Note that the `affine2d` default argument is a special case; 'commissioning' is currently the only string other than an ASDF filename that is accepted. If `None` is passed, it will perform a rotation search (least-squares fit to a PSF model) and use that for the affine transform.
+
+
 Creating ASDF files
 ^^^^^^^^^^^^^^^^^^^
 The optional arguments `bandpass` and `affine2d` must be written to `ASDF <https://asdf-standard.readthedocs.io/>`_ 

jwst/ami/ami_analyze.py

@@ -1,4 +1,3 @@
-#  Module for applying the LG-PLUS algorithm to an AMI exposure
 import logging
 import numpy as np
 import copy
@@ -30,7 +29,7 @@ def apply_LG_plus(
     run_bpfix,
 ):
     """
-    Applies the image plane algorithm to an AMI image
+    Applies the image plane algorithm (LG-PLUS) to an AMI exposure
 
     Parameters
     ----------
@@ -126,15 +125,16 @@ def apply_LG_plus(
     # Throughput (combined filter and source spectrum) calculated here
     bandpass = utils.handle_bandpass(bandpass, throughput_model)
 
-    rotsearch_d = np.append(
+    if affine2d is None:
+        log.info("Searching for best-fit affine transform")
+        rotsearch_d = np.append(
         np.arange(
             rotsearch_parameters[0], rotsearch_parameters[1], rotsearch_parameters[2]
         ),
         rotsearch_parameters[1],
-    )
+        )
 
-    log.info(f"Initial values to use for rotation search: {rotsearch_d}")
-    if affine2d is None:
+        log.info(f"Initial values to use for rotation search: {rotsearch_d}")
         # affine2d object, can be overridden by user input affine.
         # do rotation search on uncropped median image (assuming rotation constant over exposure)
         # replace remaining NaNs in median image with median of surrounding 8 (non-NaN) pixels
@@ -159,6 +159,7 @@ def apply_LG_plus(
 
         affine2d = find_rotation(
             meddata,
+            nrm_model,
             psf_offset,
             rotsearch_d,
             mx,
@@ -173,6 +174,15 @@ def apply_LG_plus(
             oversample,
             holeshape,
         )
+        log.info(f'Found rotation: {affine2d.rotradccw:.4f} rad ({np.rad2deg(affine2d.rotradccw):.4f} deg)')
+        # the affine2d returned here has only rotation...
+        # to use rotation and scaling/shear, do some matrix multiplication here??
+
+    log.info('Using affine transform with parameters:')
+    log.info(f'\tmx={affine2d.mx:.6f}\tmy={affine2d.my:.6f}')
+    log.info(f'\tsx={affine2d.sx:.6f}\tsy={affine2d.sy:.6f}')
+    log.info(f'\txo={affine2d.xo:.6f}\tyo={affine2d.yo:.6f}')
+    log.info(f'\trotradccw={affine2d.rotradccw}')
 
     niriss = instrument_data.NIRISS(filt,
                                     nrm_model,

jwst/ami/ami_analyze_step.py

@@ -9,6 +9,13 @@
 
 __all__ = ["AmiAnalyzeStep"]
 
+# Affine parameters from commissioning
+MX_COMM = 9.92820e-01
+MY_COMM = 9.98540e-01
+SX_COMM = 6.18605e-03
+SY_COMM = -7.27008e-03
+XO_COMM = 0.
+YO_COMM = 0.
 
 class AmiAnalyzeStep(Step):
     """Performs analysis of an AMI mode exposure by applying the LG algorithm."""
@@ -24,7 +31,7 @@ class AmiAnalyzeStep(Step):
         usebp = boolean(default=True) # If True, exclude pixels marked DO_NOT_USE from fringe fitting
         firstfew = integer(default=None) # If not None, process only the first few integrations
         chooseholes = string(default=None) # If not None, fit only certain fringes e.g. ['B4','B5','B6','C2']
-        affine2d = string(default=None) # ASDF file containing user-defined affine parameters
+        affine2d = string(default='commissioning') # ASDF file containing user-defined affine parameters OR 'commssioning'
         run_bpfix = boolean(default=True) # Run Fourier bad pixel fix on cropped data
     """
 
@@ -158,6 +165,10 @@ def process(self, input):
         psf_offset = [float(a) for a in self.psf_offset.split()]
         rotsearch_parameters = [float(a) for a in self.rotation_search.split()]
 
+        # handle command-line None input interpreted as string
+        if str(affine2d).lower() == 'none':
+            affine2d = None
+
         self.log.info(f"Oversampling factor = {oversample}")
         self.log.info(f"Initial rotation guess = {rotate} deg")
         self.log.info(f"Initial values to use for psf offset = {psf_offset}")
@@ -180,11 +191,19 @@ def process(self, input):
             # If there's a user-defined bandpass or affine, handle it
             if bandpass is not None:
                 bandpass = self.override_bandpass()
-
             if affine2d is not None:
-                # if it is None, handled in apply_LG_plus
-                affine2d = self.override_affine2d()
-
+                if affine2d == 'commissioning':
+                    affine2d = utils.Affine2d(mx=MX_COMM,
+                        my=MY_COMM,
+                        sx=SX_COMM,
+                        sy=SY_COMM,
+                        xo=XO_COMM,
+                        yo=YO_COMM,
+                        name='commissioning')
+                    self.log.info("Using affine parameters from commissioning.")
+                else:
+                    affine2d = self.override_affine2d()
+            # and if it is None, rotation search done in apply_LG_plus
 
             # Get the name of the NRM reference file to use
             nrm_reffile = self.get_reference_file(input_model, 'nrm')

jwst/ami/find_affine2d_parameters.py

@@ -1,7 +1,3 @@
-#
-#  Module for calculation of the optimal rotation for an image plane
-#
-
 import logging
 import numpy as np
 
@@ -14,8 +10,9 @@
 
 def create_afflist_rot(rotdegs):
     """
-    Create a list of affine objects with various rotations to use in order to
-    go through and find which fits an image plane data best.
+    Create a list of affine objects with various rotations.
+
+    Find which affine rotation fits an image plane data best.
 
     Parameters
     ----------
@@ -35,7 +32,7 @@ def create_afflist_rot(rotdegs):
     return alist
 
 
-def find_rotation(imagedata, psf_offset, rotdegs, mx, my, sx, sy, xo, yo,
+def find_rotation(imagedata, nrm_model, psf_offset, rotdegs, mx, my, sx, sy, xo, yo,
                   pixel, npix, bandpass, over, holeshape):
     """
     Create an affine2d object using the known rotation and scale.
@@ -45,6 +42,9 @@ def find_rotation(imagedata, psf_offset, rotdegs, mx, my, sx, sy, xo, yo,
     imagedata: 2D float array
         image data
 
+    nrm_model: NRMModel datamodel
+        datamodel containing mask geometry information
+
     psf_offset: 2D float array
         offset from image center in detector pixels
 
@@ -98,7 +98,7 @@ def find_rotation(imagedata, psf_offset, rotdegs, mx, my, sx, sy, xo, yo,
     crosscorr_rots = []
 
     for (rot, aff) in zip(rotdegs, affine2d_list):
-        jw = lg_model.NrmModel(mask='jwst_g7s6c', holeshape=holeshape, over=over, affine2d=aff)
+        jw = lg_model.LgModel(nrm_model, mask='jwst_ami', holeshape=holeshape, over=over, affine2d=aff)
 
         jw.set_pixelscale(pixel)
         # psf_offset in data coords & pixels.  Does it get rotated?  Second order errors poss.

jwst/ami/instrument_data.py

@@ -1,12 +1,7 @@
-#
-#  Module for defining data format, wavelength info, an mask geometry for these
-#   instrument: NIRISS AMI
-#
-
 import logging
 import numpy as np
 
-from .mask_definitions import NRM_mask_definitions
+from .mask_definition_ami import NRMDefinition
 from . import utils
 from . import bp_fix
 from stdatamodels.jwst.datamodels import dqflags
@@ -32,12 +27,18 @@ def __init__(self,
                  run_bpfix=True
                  ):
         """
-        Initialize NIRISS class
+        Initialize NIRISS class for NIRISS/AMI instrument.
+
+        Module for defining all instrument characteristics including data format, 
+        wavelength info, and mask geometry.
 
         Parameters
         ----------
         filt: string
             filter name
+        
+        nrm_model: NRMModel datamodel
+            datamodel containing mask geometry information
 
         chooseholes: list
             None, or e.g. ['B2', 'B4', 'B5', 'B6'] for a four-hole mask
@@ -91,12 +92,12 @@ def __init__(self,
         # only one NRM on JWST:
         self.telname = "JWST"
         self.instrument = "NIRISS"
-        self.arrname = "jwst_g7s6c"
-        self.holeshape = "hex"
-        self.mask = NRM_mask_definitions(
+        self.arrname = "jwst_ami"
+        self.holeshape = 'hex'
+        self.mask = NRMDefinition(
+            self.nrm_model,
             maskname=self.arrname,
-            chooseholes=self.chooseholes,
-            holeshape=self.holeshape,
+            chooseholes=self.chooseholes
         )
 
         # save affine deformation of pupil object or create a no-deformation object.

jwst/ami/leastsqnrm.py

@@ -490,7 +490,7 @@ def weighted_operations(img, model, dqm=None):
     Aw = flatmodel * weights[:, np.newaxis]
     bw = flatimg * weights
     # resids are pixel value residuals, flattened to 1d vector
-    x, rss, rank, singvals = np.linalg.lstsq(Aw, bw)
+    x, rss, rank, singvals = np.linalg.lstsq(Aw, bw, rcond=None)
 
     # actual residuals in image:
     res = flatimg - np.dot(flatmodel, x)