From 8851681dc1440b419bdb3be88d2f244d027877d9 Mon Sep 17 00:00:00 2001 From: MARCHAND MANON Date: Tue, 10 Sep 2024 17:04:02 +0200 Subject: [PATCH] feat: add approximate option for faster MOC generation fix: return an empty MOC when no image data was found --- CHANGELOG.md | 7 +++++- python/mocpy/moc/moc.py | 51 ++++++++++++++++++++++++++++++++--------- 2 files changed, 46 insertions(+), 12 deletions(-) diff --git a/CHANGELOG.md b/CHANGELOG.md index 8905885..374cb12 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -28,10 +28,15 @@ way more precise than the given WCS [#166] Small cones/boxes is faster for non-overlapping cones/boxes. * `MOC.from_fits_images` can now loop through the HDUList to only keep images with the parameter `hdu_index` set to -1 [#110] +* `MOC.from_fits_image` now has an 'approximate' option that returns a rough approximation +of the footprint of the image data from the corners of a square deduced from its WCS and +does not apply any mask. ### Fixed -* fix healpix order corresponding to 1 pixel on the image calculation in `from_fits_image` [#169] +* fix healpix order corresponding to 1 pixel on the image calculation in `MOC.from_fits_image` [#169] +* `MOC.from_fits_images` will return an empty MOC and emit a warning if there are no images in the +FITS file instead of returning an error. ## [0.16.2] diff --git a/python/mocpy/moc/moc.py b/python/mocpy/moc/moc.py index c51f286..c717869 100644 --- a/python/mocpy/moc/moc.py +++ b/python/mocpy/moc/moc.py @@ -499,7 +499,7 @@ def get_boundaries(self, order=None): return Boundaries.get(self, order) @classmethod - def from_fits_image(cls, hdu, max_norder, mask=None): + def from_fits_image(cls, hdu, max_norder, mask=None, approximate=False): """ Create a `~mocpy.moc.MOC` from an image stored as a FITS file. @@ -530,6 +530,17 @@ def from_fits_image(cls, hdu, max_norder, mask=None): # Compute a WCS from the header of the image w = wcs.WCS(header) + corners = w.calc_footprint(header) + + if approximate: + if np.isfinite(corners).all(): + sky_corners = SkyCoord(corners[:, 0], corners[:, 1], unit=u.deg) + return MOC.from_polygon_skycoord(sky_corners, max_depth=max_norder) + raise ValueError( + "Corners of at least one of the images cannot be " + "calculated with its WCS, the 'approximate' method " + "cannot be used for this image." + ) if mask is None: data = hdu.data @@ -563,8 +574,6 @@ def from_fits_image(cls, hdu, max_norder, mask=None): # Compute the deepest HEALPix order containing at least one 1 pixel of the image # We want the order so that area_hpx_cell >= area_img_pixel # <=> 4pi / (12 * 2^(2*order)) in [steradians] >= area_img_pixel in [steradians] - corners = w.calc_footprint(header) - healpix_order_computed = True if np.isfinite(corners).all(): sky_corners = SkyCoord(corners[:, 0], corners[:, 1], unit=u.deg) @@ -606,7 +615,7 @@ def from_fits_image(cls, hdu, max_norder, mask=None): return moc # noqa: RET504 @classmethod - def from_fits_images(cls, path_l, max_norder, hdu_index=0): + def from_fits_images(cls, path_l, max_norder, hdu_index=0, approximate=False): """ Load a MOC from a set of FITS file images. @@ -620,6 +629,10 @@ def from_fits_images(cls, path_l, max_norder, hdu_index=0): Index of the the HDUs containing the image in each FITS file (default = 0) If set to -1, all the HUD will be taken in account, and only the ones corresponding to images will be kept. + approximate : bool, optional + A faster but less precise way to build the MOC out of the images. This does + not mask the boolean values, and will approximate each image as a polygon + defined by the footprint deduced from the WCS. Default is False. Returns ------- @@ -635,21 +648,37 @@ def from_fits_images(cls, path_l, max_norder, hdu_index=0): with fits.open(filename) as hdul: for hdu in hdul: if ( - isinstance(hdu, (fits.ImageHDU, fits.PrimaryHDU)) - and len(hdu.data.shape) == 2 + isinstance( + hdu, (fits.ImageHDU, fits.PrimaryHDU, fits.CompImageHDU) + ) + and hdu.header + and hdu.header["NAXIS"] == 2 ): - mocs.append(MOC.from_fits_image(hdu, max_norder)) + mocs.append( # noqa: PERF401 + MOC.from_fits_image( + hdu, + max_norder, + approximate=approximate, + ) + ) else: for filename in path_l: with fits.open(filename) as hdul: - mocs.append(MOC.from_fits_image(hdu=hdul[hdu_index], - max_norder=max_norder,)) + mocs.append( + MOC.from_fits_image( + hdul[hdu_index], max_norder, approximate=approximate + ) + ) + if len(mocs) == 0: + warnings.warn( + "No image HDU found, returning an empty MOC.", UserWarning, stacklevel=2 + ) + return MOC.new_empty(max_depth=max_norder) if len(mocs) == 1: return mocs[0] - return mocs[0].union(*mocs[1:]) # this is the fastest way to do multi union - + return mocs[0].union(*mocs[1:]) # this is the fastest way to do multi union @classmethod def from_vizier_table(cls, table_id, max_depth=None, nside=None):