From 903cc3c2f05febab2a9dc9938836ac30cbd9317b Mon Sep 17 00:00:00 2001 From: David Law Date: Fri, 12 Apr 2024 09:12:16 -0400 Subject: [PATCH 1/3] JP-3491: Modify MIRI LRS WCS mapping code (#8411) Co-authored-by: Howard Bushouse --- CHANGES.rst | 7 +++ jwst/assign_wcs/miri.py | 103 ++++++++++++++++++---------------------- 2 files changed, 53 insertions(+), 57 deletions(-) diff --git a/CHANGES.rst b/CHANGES.rst index e0c94003e0..c589dccd29 100644 --- a/CHANGES.rst +++ b/CHANGES.rst @@ -6,6 +6,13 @@ ami - Replaced deprecated ``np.mat()`` with ``np.asmatrix()``. [#8415] +assign_wcs +---------- + +- Change MIRI LRS WCS code to handle the tilted trace via a centroid shift as a function + of pixel row rather than a rotation of the pixel coordinates. The practical impact is + to ensure that iso-lambda is along pixel rows after this change. [#8411] + associations ------------ diff --git a/jwst/assign_wcs/miri.py b/jwst/assign_wcs/miri.py index 4727fea77b..f8f4372687 100644 --- a/jwst/assign_wcs/miri.py +++ b/jwst/assign_wcs/miri.py @@ -289,57 +289,29 @@ def lrs_distortion(input_model, reference_files): # Now deal with the fact that the spectral trace isn't perfectly up and down along detector. # This information is contained in the xcenter/ycenter values in the CDP table, but we'll handle it - # as a simple rotation using a linear fit to this relation provided by the CDP. - - z = np.polyfit(xcen, ycen, 1) - slope = 1. / z[0] - traceangle = np.arctan(slope) * 180. / np.pi # trace angle in degrees - rot = models.Rotation2D(traceangle) # Rotation model - - # Now include this rotation in our overall transform - # First shift to a frame relative to the trace zeropoint, then apply the rotation - # to correct for the curved trace. End in a rotated frame relative to zero at the reference point - # and where yrot is aligned with the spectral trace) - xysubtoxyrot = models.Shift(-zero_point[0]) & models.Shift(-zero_point[1]) | rot - - # Next shift back to the subarray frame, and then map to v2v3 - xyrottov2v3 = models.Shift(zero_point[0]) & models.Shift(zero_point[1]) | det_to_v2v3 - - # The two models together - xysubtov2v3 = xysubtoxyrot | xyrottov2v3 - - # Work out the spectral component of the transform - # First compute the reference trace in the rotated-Y frame - xcenrot, ycenrot = rot(xcen, ycen) - # The input table of wavelengths isn't perfect, and the delta-wavelength - # steps show some unphysical behaviour - # Therefore fit with a spline for the ycenrot->wavelength transform - # Reverse vectors so that yinv is increasing (needed for spline fitting function) - yrev = ycenrot[::-1] - wrev = wavetab[::-1] + # as a simple x shift using a linear fit to this relation provided by the CDP. + # First convert the values in CDP table to subarray x/y + xcen_subarray = xcen + zero_point[0] + ycen_subarray = ycen + zero_point[1] + + # Fit for X shift as a function of Y # Spline fit with enforced smoothness - spl = UnivariateSpline(yrev, wrev, s=0.002) - # Evaluate the fit at the rotated-y reference points - wavereference = spl(yrev) - # wavereference now contains the wavelengths corresponding to regularly-sampled ycenrot, create the model - wavemodel = models.Tabular1D(lookup_table=wavereference, points=yrev, name='waveref', + spl = UnivariateSpline(ycen_subarray[::-1], xcen_subarray[::-1] - zero_point[0], s=0.002) + # Evaluate the fit at the y reference points + xshiftref = spl(ycen_subarray) + # This function will give slit dX as a function of Y subarray pixel value + dxmodel = models.Tabular1D(lookup_table=xshiftref, points=ycen_subarray, name='xshiftref', bounds_error=False, fill_value=np.nan) - # Now construct the inverse spectral transform. - # First we need to create a spline going from wavereference -> ycenrot - spl2 = UnivariateSpline(wavereference[::-1], ycenrot, s=0.002) - # Make a uniform grid of wavelength points from min to max, sampled according - # to the minimum delta in the input table - dw = np.amin(np.absolute(np.diff(wavereference))) - wmin = np.amin(wavereference) - wmax = np.amax(wavereference) - wgrid = np.arange(wmin, wmax, dw) - # Evaluate the rotated y locations of the grid - ygrid = spl2(wgrid) - # ygrid now contains the rotated y pixel locations corresponding to - # regularly-sampled wavelengths, create the model - wavemodel.inverse = models.Tabular1D(lookup_table=ygrid, points=wgrid, name='waverefinv', - bounds_error=False, fill_value=np.nan) + # Fit for the wavelength as a function of Y + # Reverse the vectors so that yinv is increasing (needed for spline fitting function) + # Spline fit with enforced smoothness + spl = UnivariateSpline(ycen_subarray[::-1], wavetab[::-1], s=0.002) + # Evaluate the fit at the y reference points + wavereference = spl(ycen_subarray) + # This model will now give the wavelength corresponding to a given Y subarray pixel value + wavemodel = models.Tabular1D(lookup_table=wavereference, points=ycen_subarray, name='waveref', + bounds_error=False, fill_value=np.nan) # Wavelength barycentric correction try: @@ -352,24 +324,41 @@ def lrs_distortion(input_model, reference_files): wavemodel = wavemodel | velocity_corr log.info("Applied Barycentric velocity correction : {}".format(velocity_corr[1].amplitude.value)) + # What is the effective slit X as a function of subarray x,y? + xmodel = models.Mapping([0], n_inputs=2) - (models.Mapping([1], n_inputs=2) | dxmodel) + # What is the effective Y as a function of subarray x,y? + ymodel = models.Mapping([1], n_inputs=2) + # What is the effective XY as a function of subarray x,y? + xymodel = models.Mapping((0, 1, 0, 1)) | xmodel & ymodel + # Define a shift by the reference point and immediately back again + # This doesn't do anything effectively, but it stores the reference point for later use in pathloss + reftransform = models.Shift(-zero_point[0]) & models.Shift(-zero_point[1]) | models.Shift(+zero_point[0]) & models.Shift(+zero_point[1]) + # Put the transforms together + xytov2v3 = reftransform | xymodel | det_to_v2v3 + # Construct the full distortion model (xsub,ysub -> v2,v3,wavelength) - lrs_wav_model = xysubtoxyrot | models.Mapping([1], n_inputs=2) | wavemodel - dettotel = models.Mapping((0, 1, 0, 1)) | xysubtov2v3 & lrs_wav_model + lrs_wav_model = models.Mapping([1], n_inputs=2) | wavemodel + dettotel = models.Mapping((0, 1, 0, 1)) | xytov2v3 & lrs_wav_model # Construct the inverse distortion model (v2,v3,wavelength -> xsub,ysub) - # Transform to get xrot from v2,v3 - v2v3toxrot = subarray_dist.inverse | xysubtoxyrot | models.Mapping([0], n_inputs=2) - # wavemodel.inverse gives yrot from wavelength - # v2,v3,lambda -> xrot,yrot - xform1 = v2v3toxrot & wavemodel.inverse - dettotel.inverse = xform1 | xysubtoxyrot.inverse + # Go from v2,v3 to slit-x + v2v3_to_xdet = det_to_v2v3.inverse | models.Mapping([0], n_inputs=2) + # Go from lambda to real y + lam_to_y = wavemodel.inverse + # Go from slit-x and real y to real-x + backwards = models.Mapping([0], n_inputs=2) + (models.Mapping([1], n_inputs=2) | dxmodel) + # Go from v2,v3,lam to real x + aa = v2v3_to_xdet & lam_to_y | backwards + # Go from v2,v3,lam to real y + bb = models.Mapping([2], n_inputs=3) | lam_to_y + # Go from v2,v3,lam, to real x,y + dettotel.inverse = models.Mapping((0, 1, 2, 0, 1, 2)) | aa & bb # Bounding box is the subarray bounding box, because we're assuming subarray coordinates passed in dettotel.bounding_box = bb_sub[::-1] return dettotel - def ifu(input_model, reference_files): """ The MIRI MRS WCS pipeline. From 6f412d72301e671974cb5847d4fbe2b00cbf99dd Mon Sep 17 00:00:00 2001 From: Melanie Clarke Date: Fri, 12 Apr 2024 10:04:15 -0400 Subject: [PATCH 2/3] JP-3543: Allow alternate resampling weight type in outlier detection (#8290) Co-authored-by: Howard Bushouse --- CHANGES.rst | 3 +++ jwst/outlier_detection/outlier_detection.py | 4 ++-- jwst/outlier_detection/outlier_detection_scaled.py | 2 +- jwst/outlier_detection/outlier_detection_spec.py | 2 +- jwst/outlier_detection/outlier_detection_step.py | 3 ++- 5 files changed, 9 insertions(+), 5 deletions(-) diff --git a/CHANGES.rst b/CHANGES.rst index c589dccd29..b912f10ceb 100644 --- a/CHANGES.rst +++ b/CHANGES.rst @@ -47,6 +47,9 @@ outlier_detection - Add association id to ``outlier_i2d`` intermediate filenames. [#8418] +- Pass the ``weight_type`` parameter to all resampling function calls so that + the default weighting can be overridden by the input step parameter. [#8290] + pipeline -------- diff --git a/jwst/outlier_detection/outlier_detection.py b/jwst/outlier_detection/outlier_detection.py index fc29873094..ef584132ec 100644 --- a/jwst/outlier_detection/outlier_detection.py +++ b/jwst/outlier_detection/outlier_detection.py @@ -113,7 +113,7 @@ def _convert_inputs(self): dq=self.inputs.dq[i]) image.meta = self.inputs.meta image.wht = build_driz_weight(image, - weight_type='ivm', + weight_type=self.outlierpars['weight_type'], good_bits=bits) self.input_models.append(image) self.converted = True @@ -199,7 +199,7 @@ def do_detection(self): for i in range(len(self.input_models)): drizzled_models[i].wht = build_driz_weight( self.input_models[i], - weight_type='ivm', + weight_type=pars['weight_type'], good_bits=pars['good_bits']) # Initialize intermediate products used in the outlier detection diff --git a/jwst/outlier_detection/outlier_detection_scaled.py b/jwst/outlier_detection/outlier_detection_scaled.py index a7252bd9af..e3b8c2cfcf 100644 --- a/jwst/outlier_detection/outlier_detection_scaled.py +++ b/jwst/outlier_detection/outlier_detection_scaled.py @@ -102,7 +102,7 @@ def do_detection(self): for image in self.input_models: image.wht = resample_utils.build_driz_weight( image, - weight_type='ivm', + weight_type=pars['weight_type'], good_bits=pars['good_bits'] ) diff --git a/jwst/outlier_detection/outlier_detection_spec.py b/jwst/outlier_detection/outlier_detection_spec.py index 1cc60db1c7..4adf2b4b06 100644 --- a/jwst/outlier_detection/outlier_detection_spec.py +++ b/jwst/outlier_detection/outlier_detection_spec.py @@ -88,7 +88,7 @@ def do_detection(self): for i in range(len(self.input_models)): drizzled_models[i].wht = resample_utils.build_driz_weight( self.input_models[i], - weight_type='ivm', + weight_type=pars['weight_type'], good_bits=pars['good_bits']) # Initialize intermediate products used in the outlier detection diff --git a/jwst/outlier_detection/outlier_detection_step.py b/jwst/outlier_detection/outlier_detection_step.py index 077fdcf782..2c90435045 100644 --- a/jwst/outlier_detection/outlier_detection_step.py +++ b/jwst/outlier_detection/outlier_detection_step.py @@ -102,7 +102,8 @@ def process(self, input_data): # Setup outlier detection parameters pars = { - 'weight_type': self.weight_type, + 'weight_type': self.weight_type, # for calling the resample step + 'wht_type': self.weight_type, # for calling the resample class directly 'pixfrac': self.pixfrac, 'kernel': self.kernel, 'fillval': self.fillval, From 9429dc65e509bec07b7912c5e6f228257d1cfacf Mon Sep 17 00:00:00 2001 From: Brett Graham Date: Mon, 15 Apr 2024 08:42:19 -0400 Subject: [PATCH 3/3] Remove configuration for non-existant schema tests (#8366) Co-authored-by: Howard Bushouse --- pyproject.toml | 3 --- 1 file changed, 3 deletions(-) diff --git a/pyproject.toml b/pyproject.toml index 833faa2d64..acccccce70 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -223,9 +223,6 @@ norecursedirs = [ "build", "venv", ] -asdf_schema_tests_enabled = true -asdf_schema_validate_default = false -asdf_schema_root = "jwst/transforms/resources/schemas jwst/datamodels/schemas" junit_family = "xunit2" inputs_root = "jwst-pipeline" results_root = "jwst-pipeline-results"