Skip to content

Commit

Permalink
Add fillval option to gwcs_blot
Browse files Browse the repository at this point in the history
  • Loading branch information
melanieclarke committed Sep 23, 2024
1 parent 1cf374b commit 09c4eed
Show file tree
Hide file tree
Showing 2 changed files with 31 additions and 4 deletions.
11 changes: 7 additions & 4 deletions src/stcal/outlier_detection/utils.py
Original file line number Diff line number Diff line change
Expand Up @@ -223,7 +223,7 @@ def flag_resampled_crs(
return mask1_smoothed & mask2


def gwcs_blot(median_data, median_wcs, blot_shape, blot_wcs, pix_ratio):
def gwcs_blot(median_data, median_wcs, blot_shape, blot_wcs, pix_ratio, fillval=0.0):
"""
Resample the median data to recreate an input image based on
the blot wcs.
Expand All @@ -236,7 +236,7 @@ def gwcs_blot(median_data, median_wcs, blot_shape, blot_wcs, pix_ratio):
median_wcs : gwcs.wcs.WCS
The wcs for the median data.
blot_shape : list of int
blot_shape : tuple of int
The target blot data shape.
blot_wcs : gwcs.wcs.WCS
Expand All @@ -245,6 +245,9 @@ def gwcs_blot(median_data, median_wcs, blot_shape, blot_wcs, pix_ratio):
pix_ratio : float
Pixel ratio.
fillval : float, optional
Fill value for missing data.
Returns
-------
blotted : numpy.ndarray
Expand All @@ -259,15 +262,15 @@ def gwcs_blot(median_data, median_wcs, blot_shape, blot_wcs, pix_ratio):
log.debug("Sci shape: {}".format(blot_shape))
log.info('Blotting {} <-- {}'.format(blot_shape, median_data.shape))

outsci = np.zeros(blot_shape, dtype=np.float32)
outsci = np.full(blot_shape, fillval, dtype=np.float32)

# Currently tblot cannot handle nans in the pixmap, so we need to give some
# other value. -1 is not optimal and may have side effects. But this is
# what we've been doing up until now, so more investigation is needed
# before a change is made. Preferably, fix tblot in drizzle.
pixmap[np.isnan(pixmap)] = -1
tblot(median_data, pixmap, outsci, scale=pix_ratio, kscale=1.0,
interp='linear', exptime=1.0, misval=0.0, sinscl=1.0)
interp='linear', exptime=1.0, misval=fillval, sinscl=1.0)

return outsci

Expand Down
24 changes: 24 additions & 0 deletions tests/outlier_detection/test_utils.py
Original file line number Diff line number Diff line change
Expand Up @@ -118,6 +118,30 @@ def test_gwcs_blot():
np.testing.assert_equal(blotted, median_data[:blot_shape[0], :blot_shape[1]])


@pytest.mark.parametrize('fillval', [0.0, np.nan])
def test_gwcs_blot_fillval(fillval):
# set up a very simple wcs that scales by 1x
output_frame = gwcs.Frame2D(name="world")
forward_transform = models.Scale(1) & models.Scale(1)

median_shape = (10, 10)
median_data = np.arange(100, dtype=np.float32).reshape((10, 10))
median_wcs = gwcs.WCS(forward_transform, output_frame=output_frame)
blot_shape = (20, 20)
blot_wcs = gwcs.WCS(forward_transform, output_frame=output_frame)
pix_ratio = 1.0

blotted = gwcs_blot(median_data, median_wcs, blot_shape, blot_wcs,
pix_ratio, fillval=fillval)

# since the blot data is larger and the wcs are equivalent the blot
# will contain the median data + some fill values
assert blotted.shape == blot_shape
np.testing.assert_equal(blotted[:median_shape[0], :median_shape[1]], median_data)
np.testing.assert_equal(blotted[median_shape[0]:, :], fillval)
np.testing.assert_equal(blotted[:, median_shape[1]:], fillval)


def test_calc_gwcs_pixmap():
# generate 2 wcses with different scales
output_frame = gwcs.Frame2D(name="world")
Expand Down

0 comments on commit 09c4eed

Please sign in to comment.