JP-2256: Apply median filter to IVM weight

[mcara]

Resolves JP-2256

This PR addresses the issue of lower flux in saturated sources in resampled images when resampling is done using the IVM weighting. By applying median filtering to the weight mask we will be able to replace low weights of saturated pixels with a median value of the weight from neighbors (suggested by @schlafly).

Checklist for maintainers

• [x] added entry in CHANGES.rst within the relevant release section
• [ ] updated or added relevant tests
• [ ] updated relevant documentation
• [x] added relevant milestone
• [x] added relevant label(s)
• [x] ran regression tests, post a link to the Jenkins job below. How to run regression tests on a PR
• [ ] Make sure the JIRA ticket is resolved properly

[schlafly]

Thanks Mihai. I like this, but I think we only want to touch the pixels that are partially saturated. i.e., some pseudocode might be:

tmp_weight = orig_weight.copy()
msaturated = (dq & PARTIALLY_SATURATED) != 0
tmp_weight[masaturated | (orig_weight == 0)] = np.nan
tmp_weight = astropy.convolution.interpolate_replace_nans(tmp_weight, Gaussian2DKernel(x_stddev=2, y_stddev=2))
interp_weight = orig_weight.copy()
interp_weight[msaturated] = new_weight[msaturated]

i.e., roughly: make a temporary weight array, replacing bad pixels and partially saturated pixels with NaN. Interpolate over those pixels to replace them. Then replace the values of partially saturated pixels with the ones from the interpolated array.

The little bit of extra care around weight zero pixels is intended to make sure that known bad pixel weights don't bleed into good neighbors. interpolate_replace_nans does a Gaussian convolution weighting rather than the median I originally proposed; I don't think the choice of statistic there matters much, and I didn't quickly find a routine that has the NaN filtering I wanted that used medians.

[mcara]

@schlafly Yes, I started implementing alternative approach a little earlier. Will make a new commit soon.

[codecov[bot]]

Codecov Report

Attention: Patch coverage is 91.66667% with 3 lines in your changes are missing coverage. Please review.

Project coverage is 75.67%. Comparing base (2fb073e) to head (d4c79b3). Report is 17 commits behind head on master.

Files	Patch %	Lines
jwst/resample/resample_utils.py	91.66%	3 Missing :warning:

Additional details and impacted files

@@            Coverage Diff             @@
##           master    #7563      +/-   ##
==========================================
+ Coverage   75.31%   75.67%   +0.36%     
==========================================
  Files         474      476       +2     
  Lines       38965    39465     +500     
==========================================
+ Hits        29345    29866     +521     
+ Misses       9620     9599      -21     

Flag	Coverage Δ		*Carryforward flag
nightly	77.63% <ø> (+0.29%)	:arrow_up:	Carriedforward from 2fb073e

*This pull request uses carry forward flags. Click here to find out more.

:umbrella: View full report in Codecov by Sentry.
:loudspeaker: Have feedback on the report? Share it here.

[mcara]

Regression tests are running here: https://plwishmaster.stsci.edu:8081/job/RT/job/JWST-Developers-Pull-Requests/678/

I expect to see differences in pixel values for tests that use IVM weighting.

[mcara]

@schlafly @hbushouse The weighting types 'ivm-med' and 'ivm-smed' allow specifying window size for the median, for example: 'ivm-med5', or 'ivm-med7', ... Would that be of interest? (I thought it might, depending on how many adjacent pixels could saturate at once) If so, the current spec specification would not allow specifying window size and it would need to be changed from option to string. An alternative would be to introduce another weight_wndsize or so parameter to the resample step.

[schlafly]

Sorry for dropping this. I think in general going out 1 or 2 FWHM is going to be what people want. The PSF usually falls off fast enough that you should always be able to get a normal pixel within an FWHM of a partially saturated pixel? But the FWHM size in pixels is variable enough that it's good to make this configurable as you have done. I don't have a sense for what options would best be of use here, though.

[mcara]

@hbushouse @nden What is your preference for an option to specify median filter size? see https://github.com/spacetelescope/jwst/pull/7563#issuecomment-1530095360

[hbushouse]

@hbushouse @nden What is your preference for an option to specify median filter size? see #7563 (comment)

Given that the whole purpose of this scheme is to reset/rebalance the weighting for pixels that are only partially saturated, I agree that it's probably only necessary to go out 1 or 2 FWHM to find pixels that have no groups saturated. The rationale being that if the core pixel(s) is only partially saturated, the PSF drops off fast enough in the neighboring pixels that the next closest 1 or 2 pix probably aren't saturated at all. I'm OK with adding a param to set this, with a suitable default in the spec. Of course the definition of "suitable" will highly filter (wavelength) dependent, because of the variation in FWHM with wavelength. To properly optimize for each instrument and filter will likely require parameter reference files in CRDS. But that would be up to the teams to create and deliver.

[melanieclarke]

@mcara - I'm going to close this PR, since this solution did not attract sufficient buy-in, and since it is now significantly out of date.