JP-1860: Fix ramp_fit counting of segments #5653
Conversation
Codecov Report
@@ Coverage Diff @@
## master #5653 +/- ##
==========================================
+ Coverage 74.46% 74.68% +0.22%
==========================================
Files 416 416
Lines 38109 38632 +523
Branches 5101 5101
==========================================
+ Hits 28378 28854 +476
- Misses 9731 9778 +47
*This pull request uses carry forward flags. Click here to find out more.
Continue to review full report at Codecov.
|
jwst/ramp_fitting/ramp_fit.py
Outdated
| gdq_cr = np.bitwise_and( gdq[nint, :, :, :], dqflags.group['JUMP_DET']) | ||
| temp_max_cr = int((gdq_cr.sum(axis=0)).max()/dqflags.group['JUMP_DET']) | ||
| max_cr = max( max_cr, temp_max_cr ) | ||
| gdq_cr = np.bitwise_and(gdq[nint], dqflags.group['JUMP_DET'] | dqflags.group['DO_NOT_USE']) |
There was a problem hiding this comment.
The expression:
dqflags.group['JUMP_DET'] | dqflags.group['DO_NOT_USE']always evaluates to 5. So the above is really:
gdq_cr = np.bitwise_and(gdq[nint], 5)which will pick out all the 1, 4 and 5 values in gdq. For example:
In [24]: gdq
Out[24]: array([0, 1, 4, 5, 8])
In [25]: np.bitwise_and(gdq, 5)
Out[25]: array([0, 1, 4, 5, 0])
But then when you sum the resulting array above, instead of getting 3, the 4 and 5 are going to give you a much higher count, 10 in this case. When really all you needed was 4 segments for 3 CRs.
Instead of summing you may want to count them using np.count_nonzero. Summing was fine when the values were all 4 and that could be corrected by dividing by 4, but not now.
There was a problem hiding this comment.
Also, please write a quick unit test for this so we don't break it again. 😄
There was a problem hiding this comment.
Hmmm, could've sworn I had put in an appropriate expression that basically returned an array of 1's and 0's (equivalent to True/False), such that the .sum() really did just count up the number of hits, as opposed to the sum of all the flag values. Will look at that again.
Unit tests: to quote the Beatles, "yeah, yeah, yeah, ..."
There was a problem hiding this comment.
Also also, I would pick either:
np.bitwise_and(gdq[nint], np.bitwise_or(dqflags.group['JUMP_DET'], dqflags.group['DO_NOT_USE']))or
gdq[nint] & (dqflags.group['JUMP_DET'] | dqflags.group['DO_NOT_USE'])not mix up the function calls and operators (which point to the function calls).
Also also also, might be nice to put something like
DO_NOT_USE = dqflags.group['DO_NOT_USE']globally in this module as we've done elsewhere to make these sorts of things a bit more clear and concise:
gdq[nint] & (JUMP_DET | DO_NOT_USE)
| # ramps, to use as a surrogate for the number of segments along the ramps | ||
| # Note that we only care about flags that are NOT in the first or last groups, | ||
| # because exclusion of a first or last group won't result in an additional segment. | ||
| max_cr = np.count_nonzero(np.bitwise_and(gdq[:, 1:-1], JUMP_DET | DO_NOT_USE), axis=1).max() |
There was a problem hiding this comment.
note that max_cr is bumped up by 1 later in line 2956, to account for the fact that n flags implies n+1 segments.
There was a problem hiding this comment.
Minor comment: Is there a reason not to bump now? Or is the later use a conceptually different value, segments, as opposed to "maximum rejects".
There was a problem hiding this comment.
No really good reason I suppose, although they are a bit different conceptually. This line simply computes the max number "interrupts" in any ramp, and then that number of interrupts turns into n+1 segments for fitting. Besides, the line is so long already, I didn't want to add yet another operation on the end of it!
| # ramps, to use as a surrogate for the number of segments along the ramps | ||
| # Note that we only care about flags that are NOT in the first or last groups, | ||
| # because exclusion of a first or last group won't result in an additional segment. | ||
| max_cr = np.count_nonzero(np.bitwise_and(gdq[:, 1:-1], JUMP_DET | DO_NOT_USE), axis=1).max() |
There was a problem hiding this comment.
Minor comment: Is there a reason not to bump now? Or is the later use a conceptually different value, segments, as opposed to "maximum rejects".
|
Errors in CI appear to be due to access problems with crds_cache. Trying again. |
Update the function in ramp_fitting that attempts to estimate the max number of segments that will be needed to fit any pixel, given the number of CR and DO_NOT_USE flags that break up each pixel's ramp. Previously the function only looked for CR (JUMP_DET) flags, but a recent update to the saturation step allows for adding DO_NOT_USE flags to groups with values below the A-to-D floor (< 0). So now the ramp for a given pixel can be broken into multiple segments due to both JUMP_DET and DO_NOT_USE flags, so we need to check for both. Checking for only JUMP_DET was causing the data arrays to be created too small to hold data from all segments, leading to an "index out of bounds" error when trying to stuff data into the array.
Fixes #5630 / JP-1860