-
Notifications
You must be signed in to change notification settings - Fork 32
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-3593: Saturation step should increase flagging for group 3 saturation #259
Conversation
Codecov ReportAll modified and coverable lines are covered by tests ✅
Additional details and impacted files@@ Coverage Diff @@
## main #259 +/- ##
==========================================
+ Coverage 83.79% 83.81% +0.01%
==========================================
Files 35 35
Lines 6998 7006 +8
==========================================
+ Hits 5864 5872 +8
Misses 1134 1134 ☔ View full report in Codecov by Sentry. |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
The "resolves" links in the description of the PR should be changed to JP-3593, instead of JP-nnnn, so it links directly to the JP ticket. The RCAL-nnnn should be removed.
Other than that, everything looks good to me.
Is there a RCAL ticket? The romancal tests show some failures: @schlafly is there someone from roman who can look at these changes? |
Thanks, yes, romancal will definitely need to be changed to accommodate this. "number of frames per group" isn't really a good concept for Roman since it can vary from group to group. I'm having trouble immediately understanding how group 3 is important here. It's not immediately clear to me that we would want this feature in Roman. Can you provide a little more context as to what the purpose is? If we wanted to not change the API, would an option be setting nframes = None by default, and only doing the new block if nframes is not None? |
Sorry, now reading the linked https://jira.stsci.edu/browse/JP-3593, I understand the motivation here. If I have understood correctly, the issue is that the ramp actually is saturating in group 3 due to a very large CR. In Roman we address this via this dilution factor approach here: |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I was under the impression that this new flagging of group 2 (when group 3 is saturated) would only occur when nframes>1. It appears to me that nframes is only ever used to compute the values in "read_pattern", but is not used anywhere in this code to trigger the new algorithm. So does that mean it's going to get applied to all exposures, including those with nframes=1?
src/stcal/saturation/saturation.py
Outdated
mask &= np.where(gdq[ints, 2, :, :] & saturated, True, False) | ||
|
||
# Flag the 2nd group for the pixels passing that gauntlet | ||
gdq[ints, 1][mask] |= saturated |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Does the "|" (or) operator really work as desired in this case? Because we track flags in individual bits, don't we need to use "np.bitwise_or()" here?
Re only triggering when nframes > 1, the reason this works is that when nframes = 1, len(read_pattern[1]) = 1, so the bit on line 93 becomes just a cut against the normal saturation level. As you say, it runs even when nframes = 1, but won't flag anything extra in that case. You could code around it to avoid the checks entirely if wanted. |
Thanks for explanation. Makes sense now. Just wasn't obvious from a cursory look at the code. |
src/stcal/saturation/saturation.py
Outdated
@@ -85,6 +84,20 @@ def flag_saturated_pixels( | |||
plane = data[ints, group, :, :] | |||
|
|||
if read_pattern is not None: | |||
if group == 2: |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I don't think this code will work here; it's looking for where SATURATED is set in groupdq for group==2, but the flags aren't passed from flag_temp into groupdq until after this code block executes.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
flag_temp? I don't see that used anywhere in here. Where is it?
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Sorry, copy-paste issue from my review of the jwst PR. Same general idea though- gdq for the third group isn't set until the code block after all the (if group==2) stuff, so it can never fulfill the criterion of identifying things where saturation was set in the third group.
src/stcal/saturation/saturation.py
Outdated
mask = data[ints, 0, ...] / np.mean(read_pattern[0]) * read_pattern[2][-1] < sat_thresh | ||
|
||
# Identify groups with suspiciously large values in the second group | ||
mask &= data[ints, 1, ...] > sat_thresh / len(read_pattern[1]) |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I think we also need something to grow this mask by n_pix_grow_sat using the adjacent_pixels routine. For pixels bordering pixels that saturate in group 3, the saturation flag is set too. Those pixels might not be intrinsically bright enough to trigger the 'suspiciously large' criterion, but will likewise be affected by the CR and restricted to only fitting a ramp to groups 1+2.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
@drlaw1558 thank you for the comments! I'm committing the changes soon. Can I please ask you to take look again when they are in?
src/stcal/saturation/saturation.py
Outdated
|
||
# now, flag any pixels that border saturated pixels | ||
if n_pix_grow_sat > 0: | ||
mask = adjacent_pixels(mask, saturated, n_pix_grow_sat) |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I don't think we can call adjacent_pixels directly; it takes an array with dq flag values, and mask is an array of True/False.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
yes, I thought I changed this but now it should be ok. @drlaw1558 can you please verify?
mask &= np.where(gdq[ints, 2, :, :] & saturated, True, False) | ||
|
||
# Flag the 2nd group for the pixels passing that gauntlet | ||
gdq[ints, 1][mask] |= saturated |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I think we need to use a temporary array here; if we set SATURATED for things passing the gauntlet above in the gdq array, and then grow the gdq array, then we'll be growing many pixels that we didn't intend to which didn't pass the gauntlet. I.e., perhaps define a zero-valued gdq_temp, set SATURATED for things that pass the gauntlet, grow gdq_temp, then add gdq_temp into gdq[ints,1, ...]
closing this PR as a new branch is dealing with the issues found for NIRISS and NIRCam data https://github.com/drlaw1558/stcal/tree/jp3593_take2 |
Resolves JP-3593
Closes spacetelescope/jwst#8413
This PR adds a check for saturation bias in group 2 for nframes > 1.
Checklist
CHANGES.rst
(either inBug Fixes
orChanges to API
)