JP-1885 Changes for when there are only three and four good groups fo…

…r jump detection (#5915)

* Improve Variable Names and Comments

In preperation for the code walkthrough I updated the comments and renamed a few variables to be more descriptive.
  I also updated a couple of the test cases.

* added ignore for profile information and anything starting with “.”

* Update .gitignore

* Create twopoint_difference.py

* CR Neighbor flagging

Initial Changes to add flagging neihbors of CRs

* remove conflicting changes

Partial changes on wrong branch removed.

* update jump.py

Try to fix conflict

* Remove Pycharm file

* Remove Pycharm File

* Delete workspace.xml

* Delete vcs.xml

* Delete other.xml

* Delete modules.xml

* Delete misc.xml

* Delete Project_Default.xml

* Delete codeStyleConfig.xml

* Delete coverage file

* Delete coverage file

* Delete coverage file

* Update twopoint_difference.py

Changes for when there are two and three difference for single pixels.

* Changes for 3 and 4 Group Integrations

Updates to the code to now use different thresholds for 3 and 4 groups integrations

* Update jump.py

Mostly changest to the interfaces to pass the new 3 and 4 group thresholds to twopoint_difference.

* Update twopoint_difference.py

Changest for to add the capability to work with 3 and 4 group integrations. This involved now tracking 3 different ways to flag the outlier differences.
Also, removed extraneous DEL commands that really don't free up memory.

* Update test_detect_jumps.py

No new tests were added. Just changes for the new interface to detect_jumps.

* Update test_twopoint_difference.py

Several new test cases were added to handle the new 3 and 4 group integrations. The tests are for both the full array search for the first jump and the followup pixel-by-pixel search for an other jumps.

* update tests

Added a test with varying saturation and CRs for each pixel.

* Update twopoint_difference.py

Updated comments to correctly reflect the current version of the code.

* Update jump_step.py

change default thresholds to be higher for the smaller number of groups.

* Update jump.py

Removed a line break to help git not have a conflict.

* Update test_twopoint_difference.py

add space after asset to match master.

* Update test_twopoint_difference.py

Trying to deal with changes in master to assert statements to avoid having to maunally resolve everything.

* Update test_twopoint_difference.py

more attempts to help git find less conflicts

* Update test_detect_jumps.py

fixing white space issues

* Update test_twopoint_difference.py

add some more white space to try to match master

* Update test_twopoint_difference.py

yet more white space

* Update test_twopoint_difference.py

again with the white space

* Update test_twopoint_difference.py

now change ] , to ],

* Update test_twopoint_difference.py

* Update test_twopoint_difference.py

more matching white space

* Update test_twopoint_difference.py

white space updates number 10

* Update jump.py

Added back missing comma at end of line 120

* Update test_twopoint_difference.py

added back expected fail to 4grps two CR case. Should still fail.

* Update test_twopoint_difference.py

Updated calling parameters to find_crs for 5grps_twcrs_2nd_5thbig. Missed the line when resolving merge conflict.

* Update test_twopoint_difference.py

added back in expected fail for 4grps two CRs.

* Pep8 and variable names

Pep8 and variable name updates

* make array smaller

* respond to comments on white space

* fix spelling

* PEP 8 updates

* Update test_twopoint_difference.py

Fix change in parameter name to find_median routines, comment XFAIL tests and rework and XFAIL to a normal expected pass version.

* Update test_detect_jumps.py

Fixed change to calling parameter name.

* Update twopoint_difference.py

removed unused variable

* Update test_jump_step.py

Changed test for minimum number of groups to be 3 instead of 5.

* Update CHANGES.rst

Add jump step changes

* Update test_jump_step.py

Remove redundant test.

* Clean up two_point difference

Update comments to explain parameters better. Also, create a new local variable for tracking the index of the largest difference to use.

* Update test_twopoint_difference.py

Remove two tests that had no asserts and didn't test anything.

* Fixes for PEP8

Fixed white space issues

* Update description.rst

Update read the docs to match changes for three and four group integrations.

* Update CHANGES.rst

Moved jump changes to be in 1.1.1 and added PR to description

CHANGES.rst

@@ -109,6 +109,13 @@ general
 - Update data products, ``calwebb_image3``, and ``source_catalog`` docs to include
   information about the segmentation map product [#5949]
 
+jump
+-----------------
+
+- Update the step to detect jumps in three and four group integations [#5915].
+- Change the default S/N ratio for not flagging neighbors to be a higher value to
+  better reflect the correct IPC.
+
 lib
 ---
 
@@ -955,6 +962,7 @@ fringe
 
 - Update the fringe step to handle 3D inputs for MIRI MRS TSO mode. [#5202]
 
+
 master_background
 -----------------
 

docs/jwst/jump/description.rst

@@ -13,7 +13,7 @@ or residual non-linearities or noise can lead to the false detection of jumps in
 due to departure from linearity.
 
 The ``jump`` step will automatically skip execution if the input data contain fewer
-than 5 groups per integration, because the baseline algorthim requires four first
+than 3 groups per integration, because the baseline algorthim requires two first
 differences to work.
 
 Algorithm
@@ -47,7 +47,12 @@ The two-point difference method is applied to each integration as follows:
 * If a jump is found in a given pixel, iterate the above steps with the
   jump-impacted group excluded, looking for additional lower-level jumps
   that still exceed the rejection threshold.
-* Stop iterating on a given pixel when no new jumps are found
+* Stop iterating on a given pixel when no new jumps are found or only one
+  difference remains.
+* If the there are only three differences (four groups), the standard median
+  is used rather than the clipped median.
+* If there are only two differences (three groups), the smallest one is compared to the larger
+  one and if the larger one is above a threshold, it is flagged as a jump.
 
 Note that any ramp values flagged as SATURATED in the input GROUPDQ array
 are not used in any of the above calculations and hence will never be

jwst/jump/jump.py

@@ -12,7 +12,7 @@
 
 
 def detect_jumps(input_model, gain_model, readnoise_model,
-                 rejection_threshold, max_cores,
+                 rejection_thresh, three_grp_thresh, four_grp_thresh, max_cores,
                  max_jump_to_flag_neighbors, min_jump_to_flag_neighbors,
                  flag_4_neighbors):
     """
@@ -100,21 +100,24 @@ def detect_jumps(input_model, gain_model, readnoise_model,
     slices = []
     # Slice up data, gdq, readnoise_2d into slices
     # Each element of slices is a tuple of
-    # (data, gdq, readnoise_2d, rejection_threshold, nframes)
+    # (data, gdq, readnoise_2d, rejection_thresh, three_grp_thresh, four_grp_thresh, nframes)
     for i in range(numslices - 1):
         slices.insert(i, (data[:, :, i * yincrement:(i + 1) * yincrement, :],
                           gdq[:, :, i * yincrement:(i + 1) * yincrement, :],
                           readnoise_2d[i * yincrement:(i + 1) * yincrement, :],
-                          rejection_threshold, frames_per_group, flag_4_neighbors,
+                          rejection_thresh, three_grp_thresh, four_grp_thresh,
+                          frames_per_group, flag_4_neighbors,
                           max_jump_to_flag_neighbors, min_jump_to_flag_neighbors))
     # last slice get the rest
     slices.insert(numslices - 1, (data[:, :, (numslices - 1) * yincrement:nrows, :],
                                   gdq[:, :, (numslices - 1) * yincrement:nrows, :],
                                   readnoise_2d[(numslices - 1) * yincrement:nrows, :],
-                                  rejection_threshold, frames_per_group, flag_4_neighbors,
+                                  rejection_thresh, three_grp_thresh, four_grp_thresh,
+                                  frames_per_group, flag_4_neighbors,
                                   max_jump_to_flag_neighbors, min_jump_to_flag_neighbors))
     if numslices == 1:
-        gdq, row_below_dq, row_above_dq = twopt.find_crs(data, gdq, readnoise_2d, rejection_threshold,
+        gdq, row_below_dq, row_above_dq = twopt.find_crs(data, gdq, readnoise_2d, rejection_thresh,
+                                                         three_grp_thresh, four_grp_thresh,
                                                          frames_per_group, flag_4_neighbors,
                                                          max_jump_to_flag_neighbors,
                                                          min_jump_to_flag_neighbors)

jwst/jump/jump_step.py

@@ -16,9 +16,11 @@ class JumpStep(Step):
 
     spec = """
         rejection_threshold = float(default=4.0,min=0) # CR sigma rejection threshold
+        three_group_rejection_threshold = float(default=6.0,min=0) # CR sigma rejection threshold
+        four_group_rejection_threshold = float(default=5.0,min=0) # CR sigma rejection threshold
         maximum_cores = option('none', 'quarter', 'half', 'all', default='none') # max number of processes to create
         flag_4_neighbors = boolean(default=True) # flag the four perpendicular neighbors of each CR
-        max_jump_to_flag_neighbors = float(default=200) # maximum jump sigma that will trigger neighbor flagging
+        max_jump_to_flag_neighbors = float(default=1000) # maximum jump sigma that will trigger neighbor flagging
         min_jump_to_flag_neighbors = float(default=10) # minimum jump sigma that will trigger neighbor flagging
     """
 
@@ -30,15 +32,17 @@ def process(self, input):
             tstart = time.time()
             # Check for an input model with NGROUPS<=2
             ngroups = input_model.data.shape[1]
-            if ngroups <= 4:
-                self.log.warning('Can not apply jump detection when NGROUPS<=4;')
+            if ngroups <= 2:
+                self.log.warning('Cannot apply jump detection when NGROUPS<=2;')
                 self.log.warning('Jump step will be skipped')
                 result = input_model.copy()
                 result.meta.cal_step.jump = 'SKIPPED'
                 return result
 
             # Retrieve the parameter values
             rej_thresh = self.rejection_threshold
+            three_grp_rej_thresh = self.three_group_rejection_threshold
+            four_grp_rej_thresh = self.four_group_rejection_threshold
             max_cores = self.maximum_cores
             max_jump_to_flag_neighbors = self.max_jump_to_flag_neighbors
             min_jump_to_flag_neighbors = self.min_jump_to_flag_neighbors
@@ -62,10 +66,9 @@ def process(self, input):
 
             # Call the jump detection routine
             result = detect_jumps(input_model, gain_model, readnoise_model,
-                                  rej_thresh, max_cores,
+                                  rej_thresh, three_grp_rej_thresh, four_grp_rej_thresh, max_cores,
                                   max_jump_to_flag_neighbors, min_jump_to_flag_neighbors,
                                   flag_4_neighbors)
-
             gain_model.close()
             readnoise_model.close()
             tstop = time.time()

jwst/jump/tests/test_detect_jumps.py

@@ -21,7 +21,7 @@ def test_nocrs_noflux(setup_inputs):
     All pixel values are zero. So slope should be zero
     """
     model, rnoise, gain = setup_inputs(ngroups=5)
-    out_model = detect_jumps(model, gain, rnoise, 4.0, 1, 200, 4, True)
+    out_model = detect_jumps(model, gain, rnoise, 4.0, 5.0, 6.0, 1, 200, 4, True)
     assert np.max(out_model.groupdq) == GOOD
 
 
@@ -33,7 +33,7 @@ def test_nocrs_noflux_badgain_pixel(setup_inputs):
     model, rnoise, gain = setup_inputs(ngroups=5, nrows=20, ncols=20)
     gain.data[7, 7] = -10  # bad gain
     gain.data[17, 17] = np.nan  # bad gain
-    out_model = detect_jumps(model, gain, rnoise, 4.0, 1, 200, 4, True)
+    out_model = detect_jumps(model, gain, rnoise, 4.0, 5.0, 6.0, 1, 200, 4, True)
 
     # 2 bits are set for each pixel, so use bitwise_and to check is set
     assert np.bitwise_and(out_model.pixeldq[7, 7], NO_GAIN_VALUE) == NO_GAIN_VALUE
@@ -48,7 +48,7 @@ def test_nocrs_noflux_subarray(setup_inputs):
     extracted from the full frame versions.
     """
     model, rnoise, gain = setup_inputs(ngroups=5, subarray=True)
-    out_model = detect_jumps(model, gain, rnoise, 4.0, 1, 200, 4, True)
+    out_model = detect_jumps(model, gain, rnoise, 4.0, 5.0, 6.0, 1, 200, 4, True)
     assert np.max(out_model.groupdq) == GOOD
 
 
@@ -74,7 +74,7 @@ def test_onecr_10_groups_neighbors_flagged(setup_inputs):
     model.data[0, 7, 5, 5] = 160.0
     model.data[0, 8, 5, 5] = 170.0
     model.data[0, 9, 5, 5] = 180.0
-    out_model = detect_jumps(model, gain, rnoise, 4.0, 1, 200, 4, True)
+    out_model = detect_jumps(model, gain, rnoise, 4.0, 5.0, 6.0, 1, 200, 4, True)
 
     assert np.max(out_model.groupdq[0, 5, 5, 5]) == JUMP_DET
     assert out_model.groupdq[0, 5, 5, 6] == JUMP_DET
@@ -109,7 +109,7 @@ def test_nocr_100_groups_nframes1(setup_inputs):
     model.data[0, 9, 5, 5] = 68.0
     for i in range(10, 100):
         model.data[0, i, 5, 5] = i * 5
-    out_model = detect_jumps(model, gain, rnoise, 4.0, 1, 200, 4, True)
+    out_model = detect_jumps(model, gain, rnoise, 4.0, 5.0, 6.0, 1, 200, 4, True)
     assert np.max(out_model.groupdq) == GOOD
 
 
@@ -146,7 +146,7 @@ def test_twoints_onecr_each_10_groups_neighbors_flagged(setup_inputs):
     model.data[1, 7, 15, 5] = 160.0
     model.data[1, 8, 15, 5] = 170.0
     model.data[1, 9, 15, 5] = 180.0
-    out_model = detect_jumps(model, gain, rnoise, 4.0, 1, 200, 4, True)
+    out_model = detect_jumps(model, gain, rnoise, 4.0, 5.0, 6.0, 1, 200, 4, True)
 
     assert np.max(out_model.groupdq[0, 5, 5, 5]) == JUMP_DET
     assert out_model.groupdq[0, 5, 5, 6] == JUMP_DET
@@ -222,7 +222,8 @@ def test_multiple_neighbor_jumps_firstlastbad(setup_inputs):
     model.groupdq[0, -1, :, :] = 1
 
     # run jump detection
-    out_model = detect_jumps(model, gain, rnoise, rejection_threshold=200.0,
+    out_model = detect_jumps(model, gain, rnoise, rejection_thresh=200.0, three_grp_thresh=200,
+                             four_grp_thresh=200,
                              max_cores=None, max_jump_to_flag_neighbors=200,
                              min_jump_to_flag_neighbors=10, flag_4_neighbors=True)
 
@@ -281,7 +282,8 @@ def test_nirspec_saturated_pix(setup_inputs):
     model.groupdq[0, :, 4, 4] = [0, 0, 2, 2, 2, 2, 2]
 
     # run jump detection
-    out_model = detect_jumps(model, gain, rnoise, rejection_threshold=200.0,
+    out_model = detect_jumps(model, gain, rnoise, rejection_thresh=200.0, three_grp_thresh=200,
+                             four_grp_thresh=200,
                              max_cores=None, max_jump_to_flag_neighbors=200,
                              min_jump_to_flag_neighbors=10, flag_4_neighbors=True)
 
@@ -334,7 +336,7 @@ def test_flagging_of_CRs_across_slice_boundaries(setup_inputs):
         model.data[1, 7, yincrement, 25] = 160.0
         model.data[1, 8, yincrement, 25] = 170.0
         model.data[1, 9, yincrement, 25] = 180.0
-        out_model = detect_jumps(model, gain, rnoise, 4.0, max_cores, 200, 4, True)
+        out_model = detect_jumps(model, gain, rnoise, 4.0, 5.0, 6.0,  max_cores, 200, 4, True)
 
         # check that the neighbors of the CR on the last row were flagged
         assert out_model.groupdq[0, 5, yincrement - 1, 5] == JUMP_DET
@@ -385,7 +387,7 @@ def test_twoints_onecr_10_groups_neighbors_flagged_multi(setup_inputs):
     model.data[1, 7, 15, 5] = 160.0
     model.data[1, 8, 15, 5] = 170.0
     model.data[1, 9, 15, 5] = 180.0
-    out_model = detect_jumps(model, gain, rnoise, 4.0, 'half', 200, 4, True)
+    out_model = detect_jumps(model, gain, rnoise, 4.0, 5.0, 6.0,  'half', 200, 4, True)
 
     assert np.max(out_model.groupdq[0, 5, 5, 5]) == JUMP_DET
     assert out_model.groupdq[0, 5, 5, 6] == JUMP_DET
@@ -409,7 +411,7 @@ def test_every_pixel_CR_neighbors_flagged(setup_inputs):
     ingain = 200
     inreadnoise = 7.0
     ngroups = 10
-    model, rnoise, gain = setup_inputs(ngroups=ngroups, gain=ingain,
+    model, rnoise, gain = setup_inputs(ngroups=ngroups, gain=ingain, nrows=100, ncols=100,
                                        readnoise=inreadnoise, deltatime=grouptime)
 
     # two segments perfect fit, second segment has twice the slope
@@ -423,7 +425,7 @@ def test_every_pixel_CR_neighbors_flagged(setup_inputs):
     model.data[0, 7, :, :] = 160.0
     model.data[0, 8, :, :] = 170.0
     model.data[0, 9, :, :] = 180.0
-    out_model = detect_jumps(model, gain, rnoise, 4.0, 'half', 200, 4, True)
+    out_model = detect_jumps(model, gain, rnoise, 4.0, 5.0, 6.0, 'half', 200, 4, True)
 
     assert np.max(out_model.groupdq[0, 5, 5, 5]) == JUMP_DET
     assert out_model.groupdq[0, 5, 5, 6] == JUMP_DET
@@ -491,7 +493,7 @@ def test_crs_on_edge_with_neighbor_flagging(setup_inputs):
     model.data[0, 8, 15, -1] = 170.0
     model.data[0, 9, 15, -1] = 180.0
 
-    out_model = detect_jumps(model, gain, rnoise, 4.0, 1, 200, 10, True)
+    out_model = detect_jumps(model, gain, rnoise, 4.0, 5.0, 6.0, 1, 200, 10, True)
 
     # flag CR and three neighbors of first row CR
     assert out_model.groupdq[0, 5, 0, 15] == JUMP_DET
@@ -538,7 +540,7 @@ def test_onecr_10_groups(setup_inputs):
     model.data[0, 7, 5, 5] = 160.0
     model.data[0, 8, 5, 5] = 170.0
     model.data[0, 9, 5, 5] = 180.0
-    out_model = detect_jumps(model, gain, rnoise, 4.0, 1, 200, 10, False)
+    out_model = detect_jumps(model, gain, rnoise, 4.0, 5.0, 6.0, 1, 200, 10, False)
 
     assert out_model.groupdq[0, 5, 5, 5] == JUMP_DET
     assert out_model.groupdq[0, 5, 4, 5] == GOOD
@@ -577,7 +579,7 @@ def test_onecr_10_groups_fullarray(setup_inputs):
     model.data[0, 7, 5, 10] = 400
     model.data[0, 8, 5, 10] = 410
     model.data[0, 9, 5, 10] = 420
-    out_model = detect_jumps(model, gain, rnoise, 4.0, 1, 200, 10, False)
+    out_model = detect_jumps(model, gain, rnoise, 4.0, 5.0, 6.0, 1, 200, 10, False)
 
     # The jump is in group 5 for columns 0-9
     assert_array_equal(out_model.groupdq[0, 5, 5, 0:10], JUMP_DET)
@@ -611,7 +613,7 @@ def test_onecr_50_groups(setup_inputs):
     model.data[0, 9, 5, 5] = 180.0
     model.data[0, 10:30, 5, 5] = np.arange(190, 290, 5)
     model.data[0, 30:50, 5, 5] = np.arange(500, 600, 5)
-    out_model = detect_jumps(model, gain, rnoise, 4.0, 1, 200, 10, False)
+    out_model = detect_jumps(model, gain, rnoise, 4.0, 5.0, 6.0, 1, 200, 10, False)
 
     # CR in group 5
     assert out_model.groupdq[0, 5, 5, 5] == JUMP_DET
@@ -649,7 +651,7 @@ def test_single_CR_neighbor_flag(setup_inputs):
     model.data[0, 9, 3, 3] = 180.0
 
     # Flag neighbors
-    out_model = detect_jumps(model, gain, rnoise, 4.0, 1, 200, 4, True)
+    out_model = detect_jumps(model, gain, rnoise, 4.0, 5.0, 6.0, 1, 200, 4, True)
 
     assert np.max(out_model.groupdq[0, 5, 3, 3]) == JUMP_DET
     assert out_model.groupdq[0, 5, 3, 4] == JUMP_DET
@@ -658,7 +660,7 @@ def test_single_CR_neighbor_flag(setup_inputs):
     assert out_model.groupdq[0, 5, 4, 3] == JUMP_DET
 
     # Do not flag neighbors
-    out_model = detect_jumps(model, gain, rnoise, 4.0, 1, 200, 4, False)
+    out_model = detect_jumps(model, gain, rnoise, 4.0, 5.0, 6.0, 1, 200, 4, False)
 
     assert np.max(out_model.groupdq[0, 5, 3, 3]) == JUMP_DET
     assert out_model.groupdq[0, 5, 3, 4] == GOOD
@@ -693,11 +695,11 @@ def test_proc(setup_inputs):
     model.data[0, 8, 2, 3] = 170.0
     model.data[0, 9, 2, 3] = 180.0
 
-    out_model_a = detect_jumps(model, gain, rnoise, 4.0, None, 200, 4, True)
-    out_model_b = detect_jumps(model, gain, rnoise, 4.0, 'half', 200, 4, True)
+    out_model_a = detect_jumps(model, gain, rnoise, 4.0, 5.0, 6.0, None, 200, 4, True)
+    out_model_b = detect_jumps(model, gain, rnoise, 4.0, 5.0, 6.0, 'half', 200, 4, True)
     assert_array_equal(out_model_a.groupdq, out_model_b.groupdq)
 
-    out_model_c = detect_jumps(model, gain, rnoise, 4.0, 'all', 200, 4, True)
+    out_model_c = detect_jumps(model, gain, rnoise, 4.0, 5.0, 6.0, 'all', 200, 4, True)
     assert_array_equal(out_model_a.groupdq, out_model_c.groupdq)
 
 
@@ -756,7 +758,7 @@ def test_adjacent_CRs(setup_inputs):
     model.data[0, 8, y, x] = 160.0
     model.data[0, 9, y, x] = 190.0
 
-    out_model = detect_jumps(model, gain, rnoise, 4.0, 'half', 200, 4, True)
+    out_model = detect_jumps(model, gain, rnoise, 4.0, 5.0, 6.0, 'half', 200, 4, True)
 
     # 1st CR (centered at x=2, y=3)
     assert out_model.groupdq[0, 5, 2, 2] == JUMP_DET

jwst/jump/tests/test_jump_step.py

@@ -253,29 +253,12 @@ def test_two_group_integration(generate_miri_reffiles, max_cores, setup_inputs):
     assert(out_model.meta.cal_step.jump == 'SKIPPED')
 
 
-def test_four_group_integration(generate_miri_reffiles, setup_inputs):
+def test_three_group_integration(generate_miri_reffiles, setup_inputs):
     override_gain, override_readnoise = generate_miri_reffiles
     grouptime = 3.0
     ingain = 6
     inreadnoise = 7
-    ngroups = 4
-    xsize = 103
-    ysize = 102
-    model1, gdq, rnModel, pixdq, err, gain = setup_inputs(ngroups=ngroups,
-                                                          nrows=ysize, ncols=xsize,
-                                                          gain=ingain, readnoise=inreadnoise,
-                                                          deltatime=grouptime)
-    out_model = JumpStep.call(model1, override_gain=override_gain,
-                              override_readnoise=override_readnoise, maximum_cores='none')
-    assert(out_model.meta.cal_step.jump == 'SKIPPED')
-
-
-def test_five_group_integration(generate_miri_reffiles, setup_inputs):
-    override_gain, override_readnoise = generate_miri_reffiles
-    grouptime = 3.0
-    ingain = 6
-    inreadnoise = 7
-    ngroups = 5
+    ngroups = 3
     xsize = 103
     ysize = 102
     model1, gdq, rnModel, pixdq, err, gain = setup_inputs(ngroups=ngroups,