From 9e79eff8dbef8fdef65465d6be973e672f625702 Mon Sep 17 00:00:00 2001 From: Maria Pena-Guerrero Date: Tue, 21 May 2024 15:57:23 -0400 Subject: [PATCH 01/10] check for saturation bias in group 2 for nframes > 1 --- jwst/saturation/saturation.py | 16 +++++++++++++++- 1 file changed, 15 insertions(+), 1 deletion(-) diff --git a/jwst/saturation/saturation.py b/jwst/saturation/saturation.py index fe81f8d370..fed75b12b5 100644 --- a/jwst/saturation/saturation.py +++ b/jwst/saturation/saturation.py @@ -49,6 +49,7 @@ def flag_saturation(input_model, ref_model, n_pix_grow_sat): """ data = input_model.data + nframes = input_model.meta.exposure.nframes # Create the output model as a copy of the input output_model = input_model.copy() @@ -69,7 +70,7 @@ def flag_saturation(input_model, ref_model, n_pix_grow_sat): ref_sub_model.close() gdq_new, pdq_new, zframe = flag_saturated_pixels( - data, gdq, pdq, sat_thresh, sat_dq, ATOD_LIMIT, dqflags.pixel, + data, gdq, pdq, sat_thresh, sat_dq, ATOD_LIMIT, dqflags.pixel, nframes, n_pix_grow_sat=n_pix_grow_sat, zframe=zframe) # Save the flags in the output GROUPDQ array @@ -117,6 +118,7 @@ def irs2_flag_saturation(input_model, ref_model, n_pix_grow_sat): nints = data.shape[0] ngroups = data.shape[1] detector = input_model.meta.instrument.detector + nframes = input_model.meta.exposure.nframes # create a mask of the appropriate size irs2_mask = x_irs2.make_mask(input_model) @@ -160,6 +162,18 @@ def irs2_flag_saturation(input_model, ref_model, n_pix_grow_sat): irs2_mask, detector) # check for saturation flag_temp = np.where(sci_temp >= sat_thresh, SATURATED, 0) + if group == 2 and nframes > 1: + # Look for pixels for which SATURATION is set in group 3. If the count + # value in group 1 is less than 1/10 of the value that would be required + # to trigger SATURATION (i.e., it isn't obviously a very bright source), + # then also set the SATURATION flag for this pixel in group 2. + saturation_in_goup3 = np.where(flag_temp[:, 2, ...] == SATURATED, True, False) + if np.any(saturation_in_goup3): + nints, nrows, ncols = np.shape(saturation_in_goup3) + for ni in range(nints): + satgp3mask = saturation_in_goup3[ni] + flag_temp[ni, 1][satgp3mask] = np.where(sci_temp[ni, 0][satgp3mask] < ATOD_LIMIT/10., + SATURATED, flag_temp[ni, 1][satgp3mask]) # check for A/D floor flaglow_temp = np.where(sci_temp <= 0, AD_FLOOR | DONOTUSE, 0) From f2ba9c8019aabca2773e56ba9776d78b2769c13c Mon Sep 17 00:00:00 2001 From: Maria Pena-Guerrero Date: Tue, 21 May 2024 16:12:51 -0400 Subject: [PATCH 02/10] added PR number --- CHANGES.rst | 5 +++++ 1 file changed, 5 insertions(+) diff --git a/CHANGES.rst b/CHANGES.rst index c23e2bf3a4..5e8cb63ea6 100644 --- a/CHANGES.rst +++ b/CHANGES.rst @@ -117,6 +117,11 @@ pipeline - Removed unused ``scale_detection`` argument from ``calwebb_tso3`` pipeline. [#8438] +saturation +---------- + +- Adds a check for saturation bias in group 2 for nframes > 1. [#8499] + ramp_fitting ------------ From d4658b39f7d25e0ccf7ad6a4e13b60f52c50cfc3 Mon Sep 17 00:00:00 2001 From: Maria Pena-Guerrero Date: Wed, 22 May 2024 11:49:38 -0400 Subject: [PATCH 03/10] fixed tests and code --- jwst/saturation/saturation.py | 5 +++-- jwst/saturation/tests/test_saturation.py | 7 +++++-- 2 files changed, 8 insertions(+), 4 deletions(-) diff --git a/jwst/saturation/saturation.py b/jwst/saturation/saturation.py index fed75b12b5..86a7e7b0c0 100644 --- a/jwst/saturation/saturation.py +++ b/jwst/saturation/saturation.py @@ -172,8 +172,9 @@ def irs2_flag_saturation(input_model, ref_model, n_pix_grow_sat): nints, nrows, ncols = np.shape(saturation_in_goup3) for ni in range(nints): satgp3mask = saturation_in_goup3[ni] - flag_temp[ni, 1][satgp3mask] = np.where(sci_temp[ni, 0][satgp3mask] < ATOD_LIMIT/10., - SATURATED, flag_temp[ni, 1][satgp3mask]) + gp1mask = np.where(sci_temp[ni, 0][satgp3mask] < ATOD_LIMIT/10., True, False) + if np.any(gp1mask): + flag_temp[ni, 1, ...] = np.where(gp1mask, SATURATED, flag_temp[ni, 1, ...]) # check for A/D floor flaglow_temp = np.where(sci_temp <= 0, AD_FLOOR | DONOTUSE, 0) diff --git a/jwst/saturation/tests/test_saturation.py b/jwst/saturation/tests/test_saturation.py index be74671006..2206f7e4f4 100644 --- a/jwst/saturation/tests/test_saturation.py +++ b/jwst/saturation/tests/test_saturation.py @@ -28,8 +28,8 @@ def test_basic_saturation_flagging(setup_nrc_cube): # Add ramp values up to the saturation limit data.data[0, 0, 5, 5] = 0 data.data[0, 1, 5, 5] = 20000 - data.data[0, 2, 5, 5] = 40000 - data.data[0, 3, 5, 5] = 60000 # Signal reaches saturation limit + data.data[0, 2, 5, 5] = 60000 # Signal reaches saturation limit + data.data[0, 3, 5, 5] = 60000 data.data[0, 4, 5, 5] = 62000 # Set saturation value in the saturation model @@ -386,6 +386,7 @@ def _cube(ngroups, nrows, ncols): data_model.meta.subarray.xsize = ncols data_model.meta.subarray.ysize = nrows data_model.meta.exposure.ngroups = ngroups + data_model.meta.exposure.nframes = 3 data_model.meta.instrument.name = 'NIRCAM' data_model.meta.instrument.detector = 'NRCA1' data_model.meta.observation.date = '2017-10-01' @@ -432,6 +433,7 @@ def _cube(xstart, ystart, ngroups, nrows, ncols): data_model.meta.subarray.xsize = ncols data_model.meta.subarray.ystart = ystart data_model.meta.subarray.ysize = nrows + data_model.meta.exposure.nframes = 3 # create a saturation model for the saturation step saturation_model = SaturationModel((1032, 1024)) @@ -477,6 +479,7 @@ def _cube(): data_model.meta.exposure.nrs_normal = 16 data_model.meta.exposure.nrs_reference = 4 data_model.meta.exposure.readpatt = 'NRSIRS2RAPID' + data_model.meta.exposure.nframes = 3 # create a saturation model for the saturation step saturation_model = SaturationModel((2048, 2048)) From f8e394cf07c7c9889143e8756ac90fe97ca77aa9 Mon Sep 17 00:00:00 2001 From: Maria Pena-Guerrero Date: Wed, 22 May 2024 19:00:17 -0400 Subject: [PATCH 04/10] avoid breaking Roman code --- jwst/saturation/saturation.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/jwst/saturation/saturation.py b/jwst/saturation/saturation.py index 86a7e7b0c0..c897e90c1f 100644 --- a/jwst/saturation/saturation.py +++ b/jwst/saturation/saturation.py @@ -70,8 +70,8 @@ def flag_saturation(input_model, ref_model, n_pix_grow_sat): ref_sub_model.close() gdq_new, pdq_new, zframe = flag_saturated_pixels( - data, gdq, pdq, sat_thresh, sat_dq, ATOD_LIMIT, dqflags.pixel, nframes, - n_pix_grow_sat=n_pix_grow_sat, zframe=zframe) + data, gdq, pdq, sat_thresh, sat_dq, ATOD_LIMIT, dqflags.pixel, + n_pix_grow_sat=n_pix_grow_sat, zframe=zframe, nframes=nframes) # Save the flags in the output GROUPDQ array output_model.groupdq = gdq_new From 0455eb1743c16849991a0d1f2af218c229b44482 Mon Sep 17 00:00:00 2001 From: Maria Pena-Guerrero Date: Tue, 28 May 2024 09:22:26 -0400 Subject: [PATCH 05/10] added changes to implement read_pattern solution --- jwst/saturation/saturation.py | 8 ++++++++ 1 file changed, 8 insertions(+) diff --git a/jwst/saturation/saturation.py b/jwst/saturation/saturation.py index c897e90c1f..e91e39bb83 100644 --- a/jwst/saturation/saturation.py +++ b/jwst/saturation/saturation.py @@ -49,7 +49,9 @@ def flag_saturation(input_model, ref_model, n_pix_grow_sat): """ data = input_model.data + ngroups = input_model.meta.exposure.ngroups nframes = input_model.meta.exposure.nframes + readpatt = input_model.meta.exposure.readpatt # Create the output model as a copy of the input output_model = input_model.copy() @@ -69,6 +71,12 @@ def flag_saturation(input_model, ref_model, n_pix_grow_sat): sat_dq = ref_sub_model.dq.copy() ref_sub_model.close() + #print('\n ******** using read_pattern = ', readpatt) + #read_pattern = [[x + 1 + groupstart * nframes for x in range(nframes)] for groupstart in range(ngroups)] + #gdq_new, pdq_new, zframe = flag_saturated_pixels( + # data, gdq, pdq, sat_thresh, sat_dq, ATOD_LIMIT, dqflags.pixel, + # n_pix_grow_sat=n_pix_grow_sat, read_pattern=read_pattern, zframe=zframe) + print('\n ******** using nframes = ', nframes) gdq_new, pdq_new, zframe = flag_saturated_pixels( data, gdq, pdq, sat_thresh, sat_dq, ATOD_LIMIT, dqflags.pixel, n_pix_grow_sat=n_pix_grow_sat, zframe=zframe, nframes=nframes) From 85a7bab3c4ead8d90df3c82cd902514a1804436a Mon Sep 17 00:00:00 2001 From: Maria Pena-Guerrero Date: Fri, 31 May 2024 17:21:10 -0400 Subject: [PATCH 06/10] removing previous code and adding modification to read_pattern code --- jwst/saturation/saturation.py | 42 ++++++++++++------------ jwst/saturation/tests/test_saturation.py | 2 ++ 2 files changed, 23 insertions(+), 21 deletions(-) diff --git a/jwst/saturation/saturation.py b/jwst/saturation/saturation.py index e91e39bb83..42ed2c6a95 100644 --- a/jwst/saturation/saturation.py +++ b/jwst/saturation/saturation.py @@ -51,7 +51,6 @@ def flag_saturation(input_model, ref_model, n_pix_grow_sat): data = input_model.data ngroups = input_model.meta.exposure.ngroups nframes = input_model.meta.exposure.nframes - readpatt = input_model.meta.exposure.readpatt # Create the output model as a copy of the input output_model = input_model.copy() @@ -71,15 +70,10 @@ def flag_saturation(input_model, ref_model, n_pix_grow_sat): sat_dq = ref_sub_model.dq.copy() ref_sub_model.close() - #print('\n ******** using read_pattern = ', readpatt) - #read_pattern = [[x + 1 + groupstart * nframes for x in range(nframes)] for groupstart in range(ngroups)] - #gdq_new, pdq_new, zframe = flag_saturated_pixels( - # data, gdq, pdq, sat_thresh, sat_dq, ATOD_LIMIT, dqflags.pixel, - # n_pix_grow_sat=n_pix_grow_sat, read_pattern=read_pattern, zframe=zframe) - print('\n ******** using nframes = ', nframes) + read_pattern = [[x + 1 + groupstart * nframes for x in range(nframes)] for groupstart in range(ngroups)] gdq_new, pdq_new, zframe = flag_saturated_pixels( data, gdq, pdq, sat_thresh, sat_dq, ATOD_LIMIT, dqflags.pixel, - n_pix_grow_sat=n_pix_grow_sat, zframe=zframe, nframes=nframes) + n_pix_grow_sat=n_pix_grow_sat, read_pattern=read_pattern, zframe=zframe) # Save the flags in the output GROUPDQ array output_model.groupdq = gdq_new @@ -127,6 +121,7 @@ def irs2_flag_saturation(input_model, ref_model, n_pix_grow_sat): ngroups = data.shape[1] detector = input_model.meta.instrument.detector nframes = input_model.meta.exposure.nframes + read_pattern = [[x + 1 + groupstart * nframes for x in range(nframes)] for groupstart in range(ngroups)] # create a mask of the appropriate size irs2_mask = x_irs2.make_mask(input_model) @@ -170,19 +165,24 @@ def irs2_flag_saturation(input_model, ref_model, n_pix_grow_sat): irs2_mask, detector) # check for saturation flag_temp = np.where(sci_temp >= sat_thresh, SATURATED, 0) - if group == 2 and nframes > 1: - # Look for pixels for which SATURATION is set in group 3. If the count - # value in group 1 is less than 1/10 of the value that would be required - # to trigger SATURATION (i.e., it isn't obviously a very bright source), - # then also set the SATURATION flag for this pixel in group 2. - saturation_in_goup3 = np.where(flag_temp[:, 2, ...] == SATURATED, True, False) - if np.any(saturation_in_goup3): - nints, nrows, ncols = np.shape(saturation_in_goup3) - for ni in range(nints): - satgp3mask = saturation_in_goup3[ni] - gp1mask = np.where(sci_temp[ni, 0][satgp3mask] < ATOD_LIMIT/10., True, False) - if np.any(gp1mask): - flag_temp[ni, 1, ...] = np.where(gp1mask, SATURATED, flag_temp[ni, 1, ...]) + if group == 2: + # Identify groups which we wouldn't expect to saturate by the third group, + # on the basis of the first group + mask = sci_temp[ints, 0, ...] / np.mean(read_pattern[0]) * read_pattern[2][-1] < sat_thresh + + # Identify groups with suspiciously large values in the second group + mask &= sci_temp[ints, 1, ...] > sat_thresh / len(read_pattern[1]) + + # Identify groups that are saturated in the third group + dq_temp = x_irs2.from_irs2(groupdq[ints, 2, :, :], irs2_mask, detector) + gp3mask = np.where(dq_temp & SATURATED, True, False) + mask &= gp3mask + + # Flag the 2nd group for the pixels passing that gauntlet + dq_temp = x_irs2.from_irs2(groupdq[ints, 1, :, :], irs2_mask, detector) + dq_temp[mask] |= SATURATED + x_irs2.to_irs2(flagarray, dq_temp, irs2_mask, detector) + # check for A/D floor flaglow_temp = np.where(sci_temp <= 0, AD_FLOOR | DONOTUSE, 0) diff --git a/jwst/saturation/tests/test_saturation.py b/jwst/saturation/tests/test_saturation.py index 2206f7e4f4..c39997107b 100644 --- a/jwst/saturation/tests/test_saturation.py +++ b/jwst/saturation/tests/test_saturation.py @@ -433,6 +433,7 @@ def _cube(xstart, ystart, ngroups, nrows, ncols): data_model.meta.subarray.xsize = ncols data_model.meta.subarray.ystart = ystart data_model.meta.subarray.ysize = nrows + data_model.meta.exposure.ngroups = ngroups data_model.meta.exposure.nframes = 3 # create a saturation model for the saturation step @@ -480,6 +481,7 @@ def _cube(): data_model.meta.exposure.nrs_reference = 4 data_model.meta.exposure.readpatt = 'NRSIRS2RAPID' data_model.meta.exposure.nframes = 3 + data_model.meta.exposure.ngroups = 5 # create a saturation model for the saturation step saturation_model = SaturationModel((2048, 2048)) From 17677b03fa8cc4a844262b4d5b27f6808de5d1ec Mon Sep 17 00:00:00 2001 From: Maria Pena-Guerrero Date: Fri, 7 Jun 2024 15:27:51 -0400 Subject: [PATCH 07/10] fixing calling the right array and index accessing --- jwst/saturation/saturation.py | 13 +++++++++---- 1 file changed, 9 insertions(+), 4 deletions(-) diff --git a/jwst/saturation/saturation.py b/jwst/saturation/saturation.py index 42ed2c6a95..422f7ea50e 100644 --- a/jwst/saturation/saturation.py +++ b/jwst/saturation/saturation.py @@ -168,16 +168,21 @@ def irs2_flag_saturation(input_model, ref_model, n_pix_grow_sat): if group == 2: # Identify groups which we wouldn't expect to saturate by the third group, # on the basis of the first group - mask = sci_temp[ints, 0, ...] / np.mean(read_pattern[0]) * read_pattern[2][-1] < sat_thresh + scigp1 = x_irs2.from_irs2(data[ints, 0, :, :], irs2_mask, detector) + mask = scigp1 / np.mean(read_pattern[0]) * read_pattern[2][-1] < sat_thresh # Identify groups with suspiciously large values in the second group - mask &= sci_temp[ints, 1, ...] > sat_thresh / len(read_pattern[1]) + scigp2 = x_irs2.from_irs2(data[ints, 1, :, :], irs2_mask, detector) + mask &= scigp2 > sat_thresh / len(read_pattern[1]) # Identify groups that are saturated in the third group - dq_temp = x_irs2.from_irs2(groupdq[ints, 2, :, :], irs2_mask, detector) - gp3mask = np.where(dq_temp & SATURATED, True, False) + gp3mask = np.where(flag_temp & SATURATED, True, False) mask &= gp3mask + # now, flag any pixels that border saturated pixels + if n_pix_grow_sat > 0: + mask = adjacency_sat(mask, SATURATED, n_pix_grow_sat) + # Flag the 2nd group for the pixels passing that gauntlet dq_temp = x_irs2.from_irs2(groupdq[ints, 1, :, :], irs2_mask, detector) dq_temp[mask] |= SATURATED From b63260390ef71cd356507a0b212aa92fb7977420 Mon Sep 17 00:00:00 2001 From: Maria Pena-Guerrero Date: Mon, 17 Jun 2024 11:45:28 -0400 Subject: [PATCH 08/10] addressing comments --- jwst/saturation/saturation.py | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/jwst/saturation/saturation.py b/jwst/saturation/saturation.py index 422f7ea50e..ecac585b24 100644 --- a/jwst/saturation/saturation.py +++ b/jwst/saturation/saturation.py @@ -179,13 +179,13 @@ def irs2_flag_saturation(input_model, ref_model, n_pix_grow_sat): gp3mask = np.where(flag_temp & SATURATED, True, False) mask &= gp3mask - # now, flag any pixels that border saturated pixels - if n_pix_grow_sat > 0: - mask = adjacency_sat(mask, SATURATED, n_pix_grow_sat) - # Flag the 2nd group for the pixels passing that gauntlet dq_temp = x_irs2.from_irs2(groupdq[ints, 1, :, :], irs2_mask, detector) dq_temp[mask] |= SATURATED + # flag any pixels that border saturated pixels + if n_pix_grow_sat > 0: + dq_temp = adjacency_sat(dq_temp, SATURATED, n_pix_grow_sat) + # set the flags in flagarray for group 2, i.e. index 1 x_irs2.to_irs2(flagarray, dq_temp, irs2_mask, detector) # check for A/D floor From 850aa05aa9b022c2a149cb595ddc2837f81aadf0 Mon Sep 17 00:00:00 2001 From: Maria Pena-Guerrero Date: Mon, 17 Jun 2024 12:48:08 -0400 Subject: [PATCH 09/10] addressing comments --- jwst/saturation/saturation.py | 2 ++ 1 file changed, 2 insertions(+) diff --git a/jwst/saturation/saturation.py b/jwst/saturation/saturation.py index ecac585b24..d200c8917c 100644 --- a/jwst/saturation/saturation.py +++ b/jwst/saturation/saturation.py @@ -187,6 +187,8 @@ def irs2_flag_saturation(input_model, ref_model, n_pix_grow_sat): dq_temp = adjacency_sat(dq_temp, SATURATED, n_pix_grow_sat) # set the flags in flagarray for group 2, i.e. index 1 x_irs2.to_irs2(flagarray, dq_temp, irs2_mask, detector) + np.bitwise_or(groupdq[ints, 1, ...], flagarray, + groupdq[ints, 1, ...]) # check for A/D floor flaglow_temp = np.where(sci_temp <= 0, AD_FLOOR | DONOTUSE, 0) From 9b01f262897db27756230a87eedb532a254dee80 Mon Sep 17 00:00:00 2001 From: Maria Pena-Guerrero Date: Mon, 17 Jun 2024 12:49:04 -0400 Subject: [PATCH 10/10] addressing comments --- jwst/saturation/saturation.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/jwst/saturation/saturation.py b/jwst/saturation/saturation.py index d200c8917c..cd1c9663b0 100644 --- a/jwst/saturation/saturation.py +++ b/jwst/saturation/saturation.py @@ -185,7 +185,7 @@ def irs2_flag_saturation(input_model, ref_model, n_pix_grow_sat): # flag any pixels that border saturated pixels if n_pix_grow_sat > 0: dq_temp = adjacency_sat(dq_temp, SATURATED, n_pix_grow_sat) - # set the flags in flagarray for group 2, i.e. index 1 + # set the flags in dq array for group 2, i.e. index 1 x_irs2.to_irs2(flagarray, dq_temp, irs2_mask, detector) np.bitwise_or(groupdq[ints, 1, ...], flagarray, groupdq[ints, 1, ...])