diff --git a/.github/workflows/ci.yml b/.github/workflows/ci.yml index 046e051..0dd8905 100644 --- a/.github/workflows/ci.yml +++ b/.github/workflows/ci.yml @@ -22,7 +22,6 @@ jobs: envs: | - linux: check-style - linux: check-security - - linux: py310-xdist - linux: py311-xdist - linux: py313-xdist @@ -30,7 +29,6 @@ jobs: - windows: py312-xdist - linux: py312-xdist-cov coverage: codecov - - linux: py312-xdist-devdeps secrets: CODECOV_TOKEN: ${{ secrets.CODECOV_TOKEN }} diff --git a/CHANGES.rst b/CHANGES.rst index 2f9978c..fdfc8c6 100644 --- a/CHANGES.rst +++ b/CHANGES.rst @@ -4,6 +4,15 @@ Release Notes ============= +2.0.1 (Unreleased) +================== + +- Update ``utils.calc_pixmap`` code to be ready for upcoming changes in GWCS + due to which inverse WCS transformations will respect bounding box by + allowing the caller of ``utils.calc_pixmap`` to disable the bounding box(es) + on both or one of the input WCS when computing the pixel map. [#164] + + 2.0.0 (2024-10-23) ================== diff --git a/drizzle/tests/test_utils.py b/drizzle/tests/test_utils.py index 5e263de..0aa182b 100644 --- a/drizzle/tests/test_utils.py +++ b/drizzle/tests/test_utils.py @@ -102,6 +102,63 @@ def test_translated_map(wcs_type): assert_almost_equal(pixmap[2:, 2:], ok_pixmap[2:, 2:], decimal=5) +def test_disable_gwcs_bbox(): + """ + Map a pixel array to a translated version ofitself. + """ + first_wcs = wcs_from_file( + "j8bt06nyq_sip_flt.fits", + ext=1, + wcs_type="gwcs" + ) + second_wcs = wcs_from_file( + "j8bt06nyq_sip_flt.fits", + ext=1, + crpix_shift=(-2, -2), # shift loaded WCS by adding this to CRPIX + wcs_type="gwcs" + ) + + ok_pixmap = np.indices(first_wcs.array_shape, dtype='float64') - 2.0 + ok_pixmap = ok_pixmap.transpose() + + # Mapping an array to a translated array + + # disable both bounding boxes: + pixmap = calc_pixmap(first_wcs, second_wcs, disable_bbox="both") + assert_almost_equal(pixmap[2:, 2:], ok_pixmap[2:, 2:], decimal=5) + assert np.all(np.isfinite(pixmap[:2, :2])) + assert np.all(np.isfinite(pixmap[-2:, -2:])) + # check bbox was restored + assert first_wcs.bounding_box is not None + assert second_wcs.bounding_box is not None + + # disable "from" bounding box: + pixmap = calc_pixmap(second_wcs, first_wcs, disable_bbox="from") + assert_almost_equal(pixmap[:-2, :-2], ok_pixmap[:-2, :-2] + 4.0, decimal=5) + assert np.all(np.logical_not(np.isfinite(pixmap[-2:, -2:]))) + # check bbox was restored + assert first_wcs.bounding_box is not None + assert second_wcs.bounding_box is not None + + # disable "to" bounding boxes: + pixmap = calc_pixmap(first_wcs, second_wcs, disable_bbox="to") + assert_almost_equal(pixmap[2:, 2:], ok_pixmap[2:, 2:], decimal=5) + assert np.all(np.isfinite(pixmap[:2, :2])) + assert np.all(pixmap[:2, :2] < 0.0) + assert np.all(np.isfinite(pixmap[-2:, -2:])) + # check bbox was restored + assert first_wcs.bounding_box is not None + assert second_wcs.bounding_box is not None + + # enable all bounding boxes: + pixmap = calc_pixmap(first_wcs, second_wcs, disable_bbox="none") + assert_almost_equal(pixmap[2:, 2:], ok_pixmap[2:, 2:], decimal=5) + assert np.all(np.logical_not(np.isfinite(pixmap[:2, :2]))) + # check bbox was restored + assert first_wcs.bounding_box is not None + assert second_wcs.bounding_box is not None + + def test_estimate_pixel_scale_ratio(): w = wcs_from_file("j8bt06nyq_flt.fits", ext=1) pscale = estimate_pixel_scale_ratio(w, w, w.wcs.crpix, (0, 0)) diff --git a/drizzle/utils.py b/drizzle/utils.py index d0d4ebf..93535eb 100644 --- a/drizzle/utils.py +++ b/drizzle/utils.py @@ -7,7 +7,7 @@ _DEG2RAD = math.pi / 180.0 -def calc_pixmap(wcs_from, wcs_to, shape=None): +def calc_pixmap(wcs_from, wcs_to, shape=None, disable_bbox="to"): """ Calculate a discretized on a grid mapping between the pixels of two images using provided WCS of the original ("from") image and the destination ("to") @@ -35,6 +35,14 @@ def calc_pixmap(wcs_from, wcs_to, shape=None): ``numpy.ndarray`` order. When provided, it takes precedence over the ``wcs_from.array_shape`` property. + disable_bbox : {"to", "from", "both", "none"}, optional + Indicates whether to use or not to use the bounding box of either + (both) ``wcs_from`` or (and) ``wcs_to`` when computing pixel map. When + ``disable_bbox`` is "none", pixel coordinates outside of the bounding + box are set to `NaN` only if ``wcs_from`` or (and) ``wcs_to`` sets + world coordinates to NaN when input pixel coordinates are outside of + the bounding box. + Returns ------- pixmap : numpy.ndarray @@ -57,27 +65,39 @@ def calc_pixmap(wcs_from, wcs_to, shape=None): If ``bounding_box`` is not available, a `ValueError` will be raised. """ + if (bbox_from := getattr(wcs_from, "bounding_box", None)) is not None: + try: + # to avoid dependency on astropy just to check whether + # the bounding box is an instance of + # modeling.bounding_box.ModelBoundingBox, we try to + # directly use and bounding_box(order='F') and if it fails, + # fall back to converting the bounding box to a tuple + # (of intervals): + bbox_from = bbox_from.bounding_box(order='F') + except AttributeError: + bbox_from = tuple(bbox_from) + + if (bbox_to := getattr(wcs_to, "bounding_box", None)) is not None: + try: + # to avoid dependency on astropy just to check whether + # the bounding box is an instance of + # modeling.bounding_box.ModelBoundingBox, we try to + # directly use and bounding_box(order='F') and if it fails, + # fall back to converting the bounding box to a tuple + # (of intervals): + bbox_to = bbox_to.bounding_box(order='F') + except AttributeError: + bbox_to = tuple(bbox_to) + if shape is None: shape = wcs_from.array_shape - if shape is None: - if (bbox := getattr(wcs_from, "bounding_box", None)) is not None: - try: - # to avoid dependency on astropy just to check whether - # the bounding box is an instance of - # modeling.bounding_box.ModelBoundingBox, we try to - # directly use and bounding_box(order='F') and if it fails, - # fall back to converting the bounding box to a tuple - # (of intervals): - bbox = bbox.bounding_box(order='F') - except AttributeError: - bbox = tuple(bbox) - - if (nd := np.ndim(bbox)) == 1: - bbox = (bbox, ) - if nd > 1: - shape = tuple( - int(math.ceil(lim[1] + 0.5)) for lim in bbox[::-1] - ) + if shape is None and bbox_from is not None: + if (nd := np.ndim(bbox_from)) == 1: + bbox_from = (bbox_from, ) + if nd > 1: + shape = tuple( + int(math.ceil(lim[1] + 0.5)) for lim in bbox_from[::-1] + ) if shape is None: raise ValueError( @@ -85,7 +105,22 @@ def calc_pixmap(wcs_from, wcs_to, shape=None): ) y, x = np.indices(shape, dtype=np.float64) - x, y = wcs_to.world_to_pixel_values(*wcs_from.pixel_to_world_values(x, y)) + + # temporarily disable the bounding box for the "from" WCS: + if disable_bbox in ["from", "both"] and bbox_from is not None: + wcs_from.bounding_box = None + if disable_bbox in ["to", "both"] and bbox_to is not None: + wcs_to.bounding_box = None + try: + x, y = wcs_to.world_to_pixel_values( + *wcs_from.pixel_to_world_values(x, y) + ) + finally: + if bbox_from is not None: + wcs_from.bounding_box = bbox_from + if bbox_to is not None: + wcs_to.bounding_box = bbox_to + pixmap = np.dstack([x, y]) return pixmap