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

Move common resample code to stcal #8986

Open
wants to merge 21 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/8986.resample.rst
Original file line number Diff line number Diff line change
@@ -0,0 +1 @@
Moved common resample code to ``stcal`` (see https://github.com/spacetelescope/stcal/pull/320). Adjusted code in ``jwst.resample`` to depend on the common code from ``stcal.resample``.
2 changes: 2 additions & 0 deletions docs/conf.py
Original file line number Diff line number Diff line change
Expand Up @@ -92,6 +92,8 @@ def check_sphinx_version(expected_version):
'gwcs': ('https://gwcs.readthedocs.io/en/stable/', None),
'stdatamodels': ('https://stdatamodels.readthedocs.io/en/latest/', None),
'stcal': ('https://stcal.readthedocs.io/en/latest/', None),
'drizzle': ('https://drizzlepac.readthedocs.io/en/latest/', None),
'tweakwcs': ('https://tweakwcs.readthedocs.io/en/latest/', None),
}

if sys.version_info[0] == 2:
Expand Down
6 changes: 3 additions & 3 deletions docs/jwst/outlier_detection/outlier_detection_spec.rst
Original file line number Diff line number Diff line change
Expand Up @@ -21,10 +21,10 @@ This routine performs identical operations to the imaging mode, with the followi
alongside the median datamodel (in the ``err`` extension) if ``save_intermediate_results``
is `True`.

#. Resampling is handled by a different class, :py:class:`~jwst.resample.resample_spec.ResampleSpecData`
instead of :py:class:`~jwst.resample.resample.ResampleData`.
#. Resampling is handled by a different class, :py:class:`~jwst.resample.resample_spec.ResampleSpec`
instead of :py:class:`~jwst.resample.resample.ResampleImage`.

#. The resampled images are written out to disk with suffix "outlier_s2d" instead of
#. The resampled images are written out to disk with suffix "outlier_s2d" instead of
"outlier_i2d" if the ``save_intermediate_results`` parameter is set to `True`.

#. The ``in_memory`` parameter has no effect, and all operations are performed in memory.
Expand Down
12 changes: 6 additions & 6 deletions docs/jwst/resample/main.rst
Original file line number Diff line number Diff line change
Expand Up @@ -98,16 +98,16 @@ this information and hence it can keep track of only 32 input images.
First bit corresponds to the first input image, second bit corresponds to the
second input image, and so on. If the number of input images is larger than 32,
then it is necessary to have multiple context images ("planes") to hold
information about all input images
information about all input images,
with the first plane encoding which of the first 32 images contributed
to the output data pixel, second plane representing next 32 input images
to the output data pixel, the second plane representing next 32 input images
(number 33-64), etc. For this reason, context array is a 3D array of the type
`numpy.int32` and shape ``(np, ny, nx)`` where ``nx`` and ``ny``
are dimensions of image's data. ``np`` is the number of "planes" equal to
are dimensions of the image data. ``np`` is the number of "planes" computed as
``(number of input images - 1) // 32 + 1``. If a bit at position ``k`` in a
pixel with coordinates ``(p, y, x)`` is 0 then input image number
pixel with coordinates ``(p, y, x)`` is 0, then input image number
``32 * p + k`` (0-indexed) did not contribute to the output data pixel
with array coordinates ``(y, x)`` and if that bit is 1 then input image number
with array coordinates ``(y, x)`` and if that bit is 1, then input image number
``32 * p + k`` did contribute to the pixel ``(y, x)`` in the resampled image.

As an example, let's assume we have 8 input images. Then, when ``'CON'`` pixel
Expand Down Expand Up @@ -143,7 +143,7 @@ context array ``con``, one can do something like this:
>>> np.flatnonzero([v & (1 << k) for v in con[:, y, x] for k in range(32)])

For convenience, this functionality was implemented in the
:py:func:`~jwst.resample.resample_utils.decode_context` function.
:py:func:`~drizzle.utils.decode_context` function.


References
Expand Down
4 changes: 2 additions & 2 deletions docs/jwst/resample/resample.rst
Original file line number Diff line number Diff line change
@@ -1,6 +1,6 @@
.. resample_:
Python Interface to Drizzle: ResampleData()
===========================================
Python Interface to Drizzle: ResampleImage()
============================================

.. automodapi:: jwst.resample.resample
4 changes: 2 additions & 2 deletions docs/jwst/resample_spec/resample_spec.rst
Original file line number Diff line number Diff line change
@@ -1,6 +1,6 @@
.. resample_spec_:
Python Interface to Drizzle: ResampleSpecData()
===============================================
Python Interface to Drizzle: ResampleSpec()
===========================================

.. automodapi:: jwst.resample.resample_spec
8 changes: 5 additions & 3 deletions jwst/outlier_detection/imaging.py
Original file line number Diff line number Diff line change
Expand Up @@ -96,15 +96,17 @@ def detect_outliers(
return input_models

if resample_data:
resamp = resample.ResampleData(
resamp = resample.ResampleImage(
input_models,
single=True,
blendheaders=False,
wht_type=weight_type,
weight_type=weight_type,
pixfrac=pixfrac,
kernel=kernel,
fillval=fillval,
good_bits=good_bits,
enable_ctx=False,
enable_var=False,
compute_err=None,
Comment on lines +107 to +109
Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Are these ever set to non-default values in JWST? They are hard-coded here, and I don't think (correct me if I'm wrong) they are ever part of the **kwargs built by get_drizpars(). If that's the case, I suggest removing them as arguments to JWST's version of ResampleImage and hard-coding them at the location(s) where ResampleImage calls stcal's Resample.

Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Yes, they are:

if self.single:
resamp = resample_spec.ResampleSpec(
container,
enable_var=False,
compute_err="driz_err",
**self.drizpars
)
drizzled_library = resamp.resample_many_to_many()
else:
resamp = resample_spec.ResampleSpec(
container,
enable_var=True,
compute_err="from_var",
**self.drizpars
)

Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

gotcha, thanks.

)
median_data, median_wcs = median_with_resampling(
input_models,
Expand Down
6 changes: 3 additions & 3 deletions jwst/outlier_detection/spec.py
Original file line number Diff line number Diff line change
Expand Up @@ -98,15 +98,15 @@ def detect_outliers(
if resample_data is True:
# Start by creating resampled/mosaic images for
# each group of exposures
resamp = resample_spec.ResampleSpecData(
resamp = resample_spec.ResampleSpec(
input_models,
single=True,
blendheaders=False,
wht_type=weight_type,
weight_type=weight_type,
pixfrac=pixfrac,
kernel=kernel,
fillval=fillval,
good_bits=good_bits,
compute_err="driz_err",
)

median_data, median_wcs, median_err = median_with_resampling(
Expand Down
12 changes: 6 additions & 6 deletions jwst/outlier_detection/tests/test_outlier_detection.py
Original file line number Diff line number Diff line change
Expand Up @@ -13,7 +13,7 @@
from jwst.outlier_detection.utils import _flag_resampled_model_crs
from jwst.resample.tests.test_resample_step import miri_rate_model
from jwst.outlier_detection.utils import median_with_resampling, median_without_resampling
from jwst.resample.resample import ResampleData
from jwst.resample.resample import ResampleImage

OUTLIER_DO_NOT_USE = np.bitwise_or(
datamodels.dqflags.pixel["DO_NOT_USE"], datamodels.dqflags.pixel["OUTLIER"]
Expand Down Expand Up @@ -670,18 +670,18 @@ def test_drizzle_and_median_with_resample(three_sci_as_asn, tmp_cwd):

def make_resamp(input_models):
"""All defaults are same as what is run by default by outlier detection"""
in_memory = not input_models._on_disk
resamp = ResampleData(
resamp = ResampleImage(
input_models,
output="",
single=True,
blendheaders=False,
wht_type="ivm",
weight_type="ivm",
pixfrac=1.0,
kernel="square",
fillval="INDEF",
good_bits="~DO_NOT_USE",
in_memory=in_memory,
asn_id="test",
enable_var=False,
enable_ctx=False,
compute_err="driz_err",
)
return resamp
91 changes: 50 additions & 41 deletions jwst/outlier_detection/utils.py
Original file line number Diff line number Diff line change
Expand Up @@ -6,7 +6,7 @@

from jwst.lib.pipe_utils import match_nans_and_flags
from jwst.resample.resample import compute_image_pixel_area
from jwst.resample.resample_utils import build_driz_weight
from stcal.resample.utils import build_driz_weight
from stcal.outlier_detection.utils import (
compute_weight_threshold,
gwcs_blot,
Expand Down Expand Up @@ -121,7 +121,12 @@ def median_without_resampling(
drizzled_err = drizzled_model.err.copy()
else:
drizzled_err = None
weight = build_driz_weight(drizzled_model, weight_type=weight_type, good_bits=good_bits)
weight = build_driz_weight(
drizzled_model,
weight_type=weight_type,
good_bits=good_bits,
flag_name_map=datamodels.dqflags.pixel,
)
if i == 0:
median_wcs = copy.deepcopy(drizzled_model.meta.wcs)
input_shape = (ngroups,) + drizzled_data.shape
Expand Down Expand Up @@ -185,7 +190,7 @@ def median_with_resampling(
----------
input_models : ModelLibrary
The input datamodels.
resamp : resample.resample.ResampleData object
resamp : resample.resample.ResampleImage object
The controlling object for the resampling process.
maskpt : float
The weight threshold for masking out low weight pixels.
Expand All @@ -208,63 +213,67 @@ def median_with_resampling(
The median data array.
median_wcs : gwcs.WCS
A WCS corresponding to the median data.
median_error : np.ndarray, optional
A median error estimate, returned only if `return_error` is True.
median_error : np.ndarray, None, optional
A median error estimate, returned only if `return_error` is `True`.
If ``resamp.compute_err`` is not set to "driz_err", `None` will be
returned.
"""
if not resamp.single:
raise ValueError("median_with_resampling should only be used for resample_many_to_many")

in_memory = not input_models.on_disk
indices_by_group = list(input_models.group_indices.values())
ngroups = len(indices_by_group)
median_err = None

eval_med_err = False
if return_error and resamp.compute_err == "driz_err":
eval_med_err = True
log.warning(
"Returning median_error has been disabled since input "
"'resamp' object does not have 'compute_err' attribute set to "
"'driz_err'."
)

if save_intermediate_results:
# create an empty image model for the median data
median_model = datamodels.ImageModel(None)

with input_models:
for i, indices in enumerate(indices_by_group):
drizzled_model = resamp.resample_group(
input_models, indices, compute_error=return_error
)
for i, indices in enumerate(indices_by_group):
drizzled_model = resamp.resample_group(indices)

if save_intermediate_results:
# write the drizzled model to file
_fileio.save_drizzled(drizzled_model, make_output_path)
if save_intermediate_results:
# write the drizzled model to file
_fileio.save_drizzled(drizzled_model, make_output_path)

if i == 0:
median_wcs = resamp.output_wcs
input_shape = (ngroups,) + drizzled_model.data.shape
dtype = drizzled_model.data.dtype
computer = MedianComputer(input_shape, in_memory, buffer_size, dtype)
if return_error:
err_computer = MedianComputer(input_shape, in_memory, buffer_size, dtype)
else:
err_computer = None
if save_intermediate_results:
# update median model's meta with meta from the first model:
median_model.update(drizzled_model)
median_model.meta.wcs = median_wcs

weight_threshold = compute_weight_threshold(drizzled_model.wht, maskpt)
drizzled_model.data[drizzled_model.wht < weight_threshold] = np.nan
computer.append(drizzled_model.data, i)
if return_error:
drizzled_model.err[drizzled_model.wht < weight_threshold] = np.nan
err_computer.append(drizzled_model.err, i)
del drizzled_model
if i == 0:
median_wcs = resamp.output_wcs
input_shape = (ngroups,) + drizzled_model.data.shape
dtype = drizzled_model.data.dtype
computer = MedianComputer(input_shape, in_memory, buffer_size, dtype)
if eval_med_err:
err_computer = MedianComputer(input_shape, in_memory, buffer_size, dtype)
else:
err_computer = None
if save_intermediate_results:
# update median model's meta with meta from the first model:
median_model.update(drizzled_model)
median_model.meta.wcs = median_wcs

weight_threshold = compute_weight_threshold(drizzled_model.wht, maskpt)
drizzled_model.data[drizzled_model.wht < weight_threshold] = np.nan
computer.append(drizzled_model.data, i)
if eval_med_err:
drizzled_model.err[drizzled_model.wht < weight_threshold] = np.nan
err_computer.append(drizzled_model.err, i)
del drizzled_model

# Perform median combination on set of drizzled mosaics
median_data = computer.evaluate()
if return_error:
if eval_med_err:
median_err = err_computer.evaluate()
else:
median_err = None

if save_intermediate_results:
# Save median model to fits
median_model.data = median_data
if return_error:
if eval_med_err:
median_model.err = median_err
# drizzled model already contains asn_id
make_output_path = partial(make_output_path, asn_id=None)
Expand Down
Loading