Skip to content

Commit

Permalink
Revert JP-2206: fix jump for miri 3 and 4 group ints #44
Browse files Browse the repository at this point in the history
This reverts commit cabd41e.
  • Loading branch information
jdavies-st committed Oct 14, 2021
1 parent 2f4ba67 commit 166dd54
Show file tree
Hide file tree
Showing 3 changed files with 36 additions and 109 deletions.
2 changes: 0 additions & 2 deletions src/stcal/jump/jump.py
Original file line number Diff line number Diff line change
Expand Up @@ -97,13 +97,11 @@ def detect_jumps(frames_per_group, data, gdq, pdq, err,
if len(wh_g[0] > 0):
pdq[wh_g] = np.bitwise_or(pdq[wh_g], dqflags["NO_GAIN_VALUE"])
pdq[wh_g] = np.bitwise_or(pdq[wh_g], dqflags["DO_NOT_USE"])
gain_2d[wh_g] = 0

wh_g = np.where(np.isnan(gain_2d))
if len(wh_g[0] > 0):
pdq[wh_g] = np.bitwise_or(pdq[wh_g], dqflags["NO_GAIN_VALUE"])
pdq[wh_g] = np.bitwise_or(pdq[wh_g], dqflags["DO_NOT_USE"])
gain_2d[wh_g] = 0.0

# Apply gain to the SCI, ERR, and readnoise arrays so they're in units
# of electrons
Expand Down
25 changes: 19 additions & 6 deletions src/stcal/jump/twopoint_difference.py
Original file line number Diff line number Diff line change
Expand Up @@ -122,7 +122,6 @@ def find_crs(data, group_dq, read_noise, normal_rej_thresh,

# Make all the first diffs for saturated groups be equal to
# 100,000 to put them above the good values in the sorted index
max_value_index = np.nanargmax(positive_first_diffs, axis=2)
first_diffs[np.isnan(first_diffs)] = 100000.

# Here we sort the 3D array along the last axis, which is the group
Expand Down Expand Up @@ -157,10 +156,25 @@ def find_crs(data, group_dq, read_noise, normal_rej_thresh,
sigma[:, :, np.newaxis]
ratio3d = np.reshape(ratio, (nrows, ncols, ndiffs))

# Get the row and column indices of pixels whose largest non-saturated ratio
# is above the threshold
# Get the group index for each pixel of the largest non-saturated
# group, assuming the indices are sorted. 2 is subtracted from ngroups
# because we are using differences and there is one less difference
# than the number of groups. This is a 2-D array.
max_value_index = ngroups - 2 - number_sat_groups

# Extract from the sorted group indices the index of the largest
# non-saturated group.
row, col = np.where(number_sat_groups >= 0)
max_ratio2d = np.reshape(ratio3d[row, col, max_value_index[row, col]], (nrows, ncols))
max_index1d = sort_index[row, col, max_value_index[row, col]]

# reshape to a 2-D array :
max_index1 = np.reshape(max_index1d, (nrows, ncols))
max_ratio2d = np.reshape(ratio3d[row, col, max_index1[row, col]],
(nrows, ncols))
max_index1d = sort_index[row, col, 1]
max_index2d = np.reshape(max_index1d, (nrows, ncols))
last_ratio = np.reshape(ratio3d[row, col, max_index2d[row, col]],
(nrows, ncols))

# Get the row and column indices of pixels whose largest non-saturated
# ratio is above the threshold, First search all the pixels that have
Expand All @@ -176,8 +190,7 @@ def find_crs(data, group_dq, read_noise, normal_rej_thresh,

# Finally, for pixels with only two good groups, compare the SNR of the
# last good group to the two diff threshold
row2cr, col2cr = np.where(np.logical_and(ndiffs - number_sat_groups == 2,
max_ratio2d > two_diff_rej_thresh))
row2cr, col2cr = np.where(last_ratio > two_diff_rej_thresh)
log.info(f'From highest outlier Two-point found {len(row4cr)} pixels \
with at least one CR and at least four groups')
log.info(f'From highest outlier Two-point found {len(row3cr)} pixels \
Expand Down
118 changes: 17 additions & 101 deletions tests/test_twopoint_difference.py
Original file line number Diff line number Diff line change
Expand Up @@ -10,7 +10,7 @@
@pytest.fixture(scope='function')
def setup_cube():

def _cube(ngroups, readnoise=10, nrows=204, ncols=204):
def _cube(ngroups, readnoise=10):
nints = 1
nrows = 204
ncols = 204
Expand Down Expand Up @@ -92,23 +92,15 @@ def test_5grps_cr2_negjumpflux(setup_cube):

def test_3grps_cr2_noflux(setup_cube):
ngroups = 3
data, gdq, nframes, read_noise, rej_threshold = setup_cube(ngroups, nrows=2, ncols=2)
data[0, 0, 0, 0] = 10.0
data[0, 1:3, 0, 0] = 1000
data[0, 0, 0, 1] = 100
data[0, 1, 0, 1] = 95
data[0, 2, 0, 1] = 94
data[0, 0, 1, 1] = 80
data[0, 1, 1, 1] = 100
data[0, 2, 1, 1] = 111
data[0, 0, 1, 0] = 90
data[0, 1, 1, 0] = 100
data[0, 2, 1, 0] = 81
data, gdq, nframes, read_noise, rej_threshold = setup_cube(ngroups)
data[0, 0, 100, 100] = 10.0
data[0, 1:4, 100, 100] = 1000
out_gdq, row_below_gdq, row_above_gdq = find_crs(data, gdq, read_noise, rej_threshold,
rej_threshold, rej_threshold, nframes,
False, 200, 10, DQFLAGS)
assert(4 == np.max(out_gdq)) # a CR was found
assert np.array_equal([0, 4, 0], out_gdq[0, :, 0, 0])
# assert(1,np.argmax(out_gdq[0, :, 100, 100])) # find the CR in the expected group
assert(np.array_equal([0, 4, 0], out_gdq[0, :, 100, 100]))


def test_4grps_cr2_noflux(setup_cube):
Expand Down Expand Up @@ -167,7 +159,7 @@ def test_5grps_twocrs_2nd_5th(setup_cube):
rej_threshold, rej_threshold, nframes,
False, 200, 10, DQFLAGS)
assert(4 == np.max(out_gdq)) # a CR was found
assert np.array_equal([0, 4, 0, 0, 4], out_gdq[0, :, 100, 100])
assert(np.array_equal([0, 4, 0, 0, 4], out_gdq[0, :, 100, 100]))


def test_5grps_twocrs_2nd_5thbig(setup_cube):
Expand Down Expand Up @@ -320,7 +312,7 @@ def test_5grps_cr2_negslope(setup_cube):

def test_6grps_1cr(setup_cube):
ngroups = 6
data, gdq, nframes, read_noise, rej_threshold = setup_cube(ngroups)
data, gdq, nframes, read_noise, rej_threshold = setup_cube(ngroups, readnoise=10)
nframes = 1
data[0, 0, 100, 100] = 0
data[0, 1, 100, 100] = 10
Expand All @@ -336,7 +328,7 @@ def test_6grps_1cr(setup_cube):

def test_7grps_1cr(setup_cube):
ngroups = 7
data, gdq, nframes, read_noise, rej_threshold = setup_cube(ngroups)
data, gdq, nframes, read_noise, rej_threshold = setup_cube(ngroups, readnoise=10)
nframes = 1
data[0, 0, 100, 100] = 0
data[0, 1, 100, 100] = 10
Expand All @@ -353,7 +345,7 @@ def test_7grps_1cr(setup_cube):

def test_8grps_1cr(setup_cube):
ngroups = 8
data, gdq, nframes, read_noise, rej_threshold = setup_cube(ngroups)
data, gdq, nframes, read_noise, rej_threshold = setup_cube(ngroups, readnoise=10)
nframes = 1
data[0, 0, 100, 100] = 0
data[0, 1, 100, 100] = 10
Expand All @@ -371,7 +363,7 @@ def test_8grps_1cr(setup_cube):

def test_9grps_1cr_1sat(setup_cube):
ngroups = 9
data, gdq, nframes, read_noise, rej_threshold = setup_cube(ngroups)
data, gdq, nframes, read_noise, rej_threshold = setup_cube(ngroups, readnoise=10)
nframes = 1
data[0, 0, 100, 100] = 0
data[0, 1, 100, 100] = 10
Expand All @@ -391,7 +383,7 @@ def test_9grps_1cr_1sat(setup_cube):

def test_10grps_1cr_2sat(setup_cube):
ngroups = 10
data, gdq, nframes, read_noise, rej_threshold = setup_cube(ngroups)
data, gdq, nframes, read_noise, rej_threshold = setup_cube(ngroups, readnoise=10)
nframes = 1
data[0, 0, 100, 100] = 0
data[0, 1, 100, 100] = 10
Expand All @@ -413,7 +405,7 @@ def test_10grps_1cr_2sat(setup_cube):

def test_11grps_1cr_3sat(setup_cube):
ngroups = 11
data, gdq, nframes, read_noise, rej_threshold = setup_cube(ngroups)
data, gdq, nframes, read_noise, rej_threshold = setup_cube(ngroups, readnoise=10)
nframes = 1
data[0, 0, 100, 100] = 0
data[0, 1, 100, 100] = 20
Expand All @@ -437,7 +429,7 @@ def test_11grps_1cr_3sat(setup_cube):

def test_11grps_0cr_3donotuse(setup_cube):
ngroups = 11
data, gdq, nframes, read_noise, rej_threshold = setup_cube(ngroups)
data, gdq, nframes, read_noise, rej_threshold = setup_cube(ngroups, readnoise=10)
nframes = 1
data[0, 0, 100, 100] = 0
data[0, 1, 100, 100] = 18
Expand All @@ -461,7 +453,7 @@ def test_11grps_0cr_3donotuse(setup_cube):

def test_5grps_nocr(setup_cube):
ngroups = 6
data, gdq, nframes, read_noise, rej_threshold = setup_cube(ngroups)
data, gdq, nframes, read_noise, rej_threshold = setup_cube(ngroups, readnoise=10)
nframes = 1
data[0, 0, 100, 100] = 0
data[0, 1, 100, 100] = 10
Expand All @@ -475,7 +467,7 @@ def test_5grps_nocr(setup_cube):

def test_6grps_nocr(setup_cube):
ngroups = 6
data, gdq, nframes, read_noise, rej_threshold = setup_cube(ngroups)
data, gdq, nframes, read_noise, rej_threshold = setup_cube(ngroups, readnoise=10)
nframes = 1
data[0, 0, 100, 100] = 0
data[0, 1, 100, 100] = 10
Expand Down Expand Up @@ -737,6 +729,7 @@ def test_median_with_saturation_odd_number_final_difference(setup_cube):

def test_first_last_group(setup_cube):
ngroups = 7
nframes = 1
data, gdq, nframes, read_noise, rej_threshold = setup_cube(ngroups, readnoise=25.0)

# set up the data so that if the first and last group are used in jump
Expand All @@ -763,80 +756,3 @@ def test_first_last_group(setup_cube):
assert outgdq[0, 0, 100, 100] == DQFLAGS['DO_NOT_USE']
assert outgdq[0, 6, 100, 100] == DQFLAGS['DO_NOT_USE']
assert outgdq[0, 3, 100, 100] == DQFLAGS['JUMP_DET']


def test_2group(setup_cube):
# test should not find a CR, can't do it with only one difference.
ngroups = 2
data, gdq, nframes, read_noise, rej_threshold = setup_cube(ngroups, readnoise=25.0, nrows=2, ncols=2)

data[0, 0, 0, 0] = 10000.0
# set groups 1,2 - to be around 30,000
data[0, 1, 0, 0] = 30000.0

outgdq, row_below_gdq, row_above_gdq = find_crs(data, gdq, read_noise, rej_threshold,
rej_threshold, rej_threshold, nframes,
False, 200, 10, DQFLAGS)
assert outgdq[0, 1, 0, 0] == 0
assert outgdq[0, 0, 0, 0] == 0


def test_4group(setup_cube):
ngroups = 4
data, gdq, nframes, read_noise, rej_threshold = setup_cube(ngroups, readnoise=25.0, nrows=2, ncols=2)

data[0, 0, 0, 0] = 10000.0
# set groups 1,2 - to be around 30,000
data[0, 1, 0, 0] = 30000.0
data[0, 2, 0, 0] = 30020.0
data[0, 3, 0, 0] = 30000.0

outgdq, row_below_gdq, row_above_gdq = find_crs(data, gdq, read_noise, rej_threshold,
rej_threshold, rej_threshold, nframes,
False, 200, 10, DQFLAGS)
assert outgdq[0, 1, 0, 0] == 4


def test_first_last_4group(setup_cube):
ngroups = 4
data, gdq, nframes, read_noise, rej_threshold = setup_cube(ngroups, readnoise=25.0, nrows=2, ncols=2)

# set up the data so that if the first and last group are used in jump
# detection it would cause a jump to be detected between group 0-1.
data[0, 0, 0, 0] = 10000.0
# set groups 1,2 - to be around 30,000
data[0, 1, 0, 0] = 30000.0
data[0, 2, 0, 0] = 30020.0
data[0, 3, 0, 0] = 30000.0
# treat as MIRI data with first and last flagged
gdq[0, 0, :, :] = DQFLAGS['DO_NOT_USE']
gdq[0, 3, :, :] = DQFLAGS['DO_NOT_USE']
outgdq, row_below_gdq, row_above_gdq = find_crs(data, gdq, read_noise, rej_threshold,
rej_threshold, rej_threshold, nframes,
False, 200, 10, DQFLAGS)

assert outgdq[0, 0, 0, 0] == DQFLAGS['DO_NOT_USE']
assert outgdq[0, 3, 0, 0] == DQFLAGS['DO_NOT_USE']
assert outgdq[0, 1, 0, 0] == 0


def test_first_last_3group(setup_cube):
ngroups = 3
data, gdq, nframes, read_noise, rej_threshold = setup_cube(ngroups, readnoise=25.0, nrows=2, ncols=2)

# set up the data so that if the first and last group are used in jump
# detection it would cause a jump to be detected between group 1-2
# and group 6-7. Add a jump between 3 and 4 just to make sure jump detection is working
# set group 1 to be 10,000
data[0, 0, 0, 0] = 10000.0
data[0, 1, 0, 0] = 10100.0
data[0, 2, 0, 0] = 30020.0

gdq[0, 2, 0, 0] = DQFLAGS['DO_NOT_USE'] # only flag the last group
outgdq, row_below_gdq, row_above_gdq = find_crs(data, gdq, read_noise, rej_threshold,
rej_threshold, rej_threshold, nframes,
False, 200, 10, DQFLAGS)

assert outgdq[0, 0, 0, 0] == 0
assert outgdq[0, 2, 0, 0] == DQFLAGS['DO_NOT_USE']
assert outgdq[0, 1, 0, 0] == 0

0 comments on commit 166dd54

Please sign in to comment.