Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

JP-3793: Fix numbers of flagged groups #338

Open
wants to merge 5 commits into
base: main
Choose a base branch
from
Open
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
1 change: 1 addition & 0 deletions changes/338.bugfix.rst
Original file line number Diff line number Diff line change
@@ -0,0 +1 @@
Fix a bug in the jump step with counting the number of usable groups.
39 changes: 27 additions & 12 deletions src/stcal/jump/twopoint_difference.py
Original file line number Diff line number Diff line change
Expand Up @@ -225,24 +225,39 @@
num_flagged_grps = 0
# determine the number of groups with all pixels set to DO_NOT_USE
ngrps = dat.shape[1]
max_flagged_grps = 0
total_flagged_grps = 0
for integ in range(nints):
num_flagged_grps = 0
for grp in range(dat.shape[1]):
if np.all(np.bitwise_and(gdq[integ, grp, :, :], dnu_flag)):
num_flagged_grps += 1
if num_flagged_grps > max_flagged_grps:
max_flagged_grps = num_flagged_grps
total_flagged_grps += num_flagged_grps
if only_use_ints:
total_sigclip_groups = nints
else:
total_sigclip_groups = nints * (ngrps - num_flagged_grps)
total_groups = nints * (ngrps - num_flagged_grps)
total_diffs = nints * (ngrps - 1 - num_flagged_grps)
total_usable_diffs = total_diffs - num_flagged_grps
if ((ngrps < minimum_groups and only_use_ints and nints < minimum_sigclip_groups) or
(not only_use_ints and nints * ngrps < minimum_sigclip_groups and
total_groups < minimum_groups)):
total_sigclip_groups = nints * ngrps - num_flagged_grps

Check warning on line 241 in src/stcal/jump/twopoint_difference.py

View check run for this annotation

Codecov / codecov/patch

src/stcal/jump/twopoint_difference.py#L241

Added line #L241 was not covered by tests

min_usable_groups = ngrps - max_flagged_grps
total_groups = nints * ngrps - total_flagged_grps
min_usable_diffs = min_usable_groups - 1
sig_clip_grps_fails = False
total_noise_min_grps_fails = False

# Determine whether there are enough usable groups the two sigma clip options
if (only_use_ints and nints < minimum_sigclip_groups) \
or \
(not only_use_ints and total_sigclip_groups < minimum_sigclip_groups):
sig_clip_grps_fails = True
if min_usable_groups < minimum_groups:
total_noise_min_grps_fails = True

if total_noise_min_grps_fails and sig_clip_grps_fails:
log.info("Jump Step was skipped because exposure has less than the minimum number of usable groups")
dummy = np.zeros((dataa.shape[1] - 1, dataa.shape[2], dataa.shape[3]),
dtype=np.float32)
return gdq, row_below_gdq, row_above_gdq, 0, dummy
dummy = np.zeros((ngroups - 1, nrows, ncols), dtype=np.float32)
return gdq, row_below_gdq, row_above_gdq, -99, dummy
else:
# set 'saturated' or 'do not use' pixels to nan in data
dat[gdq & (dnu_flag | sat_flag) != 0] = np.nan
Expand Down Expand Up @@ -301,8 +316,8 @@

warnings.resetwarnings()
else: # There are not enough groups for sigma clipping

if total_usable_diffs >= min_diffs_single_pass:
if min_usable_diffs >= min_diffs_single_pass:
# There are enough diffs in all ints to look for more than one jump

# compute 'ratio' for each group. this is the value that will be
# compared to 'threshold' to classify jumps. subtract the median of
Expand Down
83 changes: 83 additions & 0 deletions tests/test_jump_twopoint_difference.py
Original file line number Diff line number Diff line change
Expand Up @@ -22,6 +22,89 @@ def _cube(ngroups, nints=1, nrows=204, ncols=204, readnoise=10):

return _cube

def test_sigclip_not_enough_groups_ng10(setup_cube):
ngroups = 10
data, gdq, nframes, read_noise, rej_threshold = setup_cube(ngroups, nints=8, nrows=2, ncols=2, readnoise=8)
out_gdq, row_below_gdq, rows_above_gdq, total_crs, stddev = find_crs(
data, gdq, read_noise, rej_threshold, rej_threshold, rej_threshold,
nframes, False, 200, 10, DQFLAGS,
minimum_sigclip_groups=10, minimum_groups=5
)
assert total_crs == -99

def test_sigclip_not_enough_groups_ng5(setup_cube):
ngroups = 5
data, gdq, nframes, read_noise, rej_threshold = setup_cube(ngroups, nints=8, nrows=2, ncols=2, readnoise=8)
out_gdq, row_below_gdq, rows_above_gdq, total_crs, stddev = find_crs(
data, gdq, read_noise, rej_threshold, rej_threshold, rej_threshold,
nframes, False, 200, 10, DQFLAGS,
minimum_sigclip_groups=11, minimum_groups=6
)
assert total_crs == -99

def test_sigclip_not_enough_groups_ng12(setup_cube):
ngroups = 12

data, gdq, nframes, read_noise, rej_threshold = setup_cube(ngroups, nints=8, nrows=2, ncols=2, readnoise=8)
data[0, 0:2, 0:, 0] = 1
out_gdq, row_below_gdq, rows_above_gdq, total_crs, stddev = find_crs(
data, gdq, read_noise, rej_threshold, rej_threshold, rej_threshold,
nframes, False, 200, 10, DQFLAGS,
minimum_sigclip_groups=11, minimum_groups=16
)
assert total_crs == -99

def test_sigclip_with_num_small_groups_ni190(setup_cube):
ngroups = 12
nints = 190

data, gdq, nframes, read_noise, rej_threshold = setup_cube(ngroups, nints=nints, nrows=2, ncols=2, readnoise=8)
data[0, 0:2, 0:, 0] = 1
out_gdq, row_below_gdq, rows_above_gdq, total_crs, stddev = find_crs(
data, gdq, read_noise, rej_threshold, rej_threshold, rej_threshold,
nframes, False, 200, 10, DQFLAGS,
minimum_sigclip_groups=11, minimum_groups=16
)
assert total_crs != -99

def test_nosigclip_with_enough_groups_ni1(setup_cube):
ngroups = 12
nints = 1

data, gdq, nframes, read_noise, rej_threshold = setup_cube(ngroups, nints=nints, nrows=2, ncols=2, readnoise=8)
data[0, 0:2, 0:, 0] = 1
out_gdq, row_below_gdq, rows_above_gdq, total_crs, stddev = find_crs(
data, gdq, read_noise, rej_threshold, rej_threshold, rej_threshold,
nframes, False, 200, 10, DQFLAGS,
minimum_sigclip_groups=11, minimum_groups=5
)
assert total_crs != -99

def test_nosigclip_with_num_small_groups_ni1(setup_cube):
ngroups = 4
nints = 1

data, gdq, nframes, read_noise, rej_threshold = setup_cube(ngroups, nints=nints, nrows=2, ncols=2, readnoise=8)
data[0, 0:2, 0:, 0] = 1
out_gdq, row_below_gdq, rows_above_gdq, total_crs, stddev = find_crs(
data, gdq, read_noise, rej_threshold, rej_threshold, rej_threshold,
nframes, False, 200, 10, DQFLAGS,
minimum_sigclip_groups=11, minimum_groups=16
)
assert total_crs == -99

def test_nosigclip_with_num_small_groups_ni100(setup_cube):
ngroups = 4
nints = 100

data, gdq, nframes, read_noise, rej_threshold = setup_cube(ngroups, nints=nints, nrows=2, ncols=2, readnoise=8)
data[0, 0:2, 0:, 0] = 1
out_gdq, row_below_gdq, rows_above_gdq, total_crs, stddev = find_crs(
data, gdq, read_noise, rej_threshold, rej_threshold, rej_threshold,
nframes, False, 200, 10, DQFLAGS,
minimum_sigclip_groups=100, minimum_groups=16
)
assert total_crs != -99

def test_varying_groups(setup_cube):
ngroups = 5
Expand Down
Loading