dfreproject.reproject ===================== .. py:module:: dfreproject.reproject Attributes ---------- .. autoapisummary:: dfreproject.reproject.EPSILON Classes ------- .. autoapisummary:: dfreproject.reproject.Reproject Functions --------- .. autoapisummary:: dfreproject.reproject.atan2d dfreproject.reproject.sincosd dfreproject.reproject.interpolate_image dfreproject.reproject.calculate_reprojection Module Contents --------------- .. py:data:: EPSILON :value: 1e-10 .. py:function:: atan2d(y, x) PyTorch implementation of WCSLib's atan2d function .. py:function:: sincosd(angle_deg) PyTorch implementation of WCSLib's sincosd function .. py:function:: interpolate_image(source_image: torch.Tensor, grid: torch.Tensor, interpolation_mode: str) -> torch.Tensor JIT-compiled image interpolation using grid_sample .. py:class:: Reproject(source_hdus: List[astropy.io.fits.PrimaryHDU], target_wcs: astropy.wcs.WCS, shape_out: Tuple[int, int], device=None) .. py:attribute:: batch_source_images .. py:attribute:: batch_source_wcs_params .. py:attribute:: target_wcs_params .. py:attribute:: target_grid .. py:method:: calculate_skyCoords(x_grid=None, y_grid=None) -> Tuple[torch.Tensor, torch.Tensor] Calculate sky coordinates. There are four primary steps: 1. Apply shift 2. Apply SIP distortion 3. Apply CD matrix 4. Apply transformation to celestial coordinates uing the Gnomonic Projection These steps use the target wcs parameters :param x_grid: Batch of x-coordinates. If None, uses target grid x-coordinates :type x_grid: torch.Tensor, optional :param y_grid: Batch of y-coordinates. If None, uses target grid y-coordinates :type y_grid: torch.Tensor, optional :returns: Batched RA and Dec coordinates :rtype: tuple .. py:method:: calculate_sourceCoords() Calculate source image pixel coordinates corresponding to each target image pixel. This function repeats the same steps in self.calculate_skyCoords() except in the opposite order and with the source coordinate wcs. :returns: Batch of source image pixel coordinates :rtype: torch.Tensor .. py:method:: interpolate_source_image(interpolation_mode='bilinear') -> torch.Tensor Interpolate the source image at the calculated source coordinates with flux conservation. This method performs the actual pixel resampling needed for dfreproject while preserving the total flux (photometric accuracy). It implements a footprint-based approach similar to that used in reproject_interp from the Astropy package. The method uses a combined tensor approach for computational efficiency, performing both image resampling and footprint tracking in a single operation. Total flux is conserved locally (via footprint correction). :param interpolation_mode: The interpolation mode to use when sampling the source image. Options include: - 'nearest' : Nearest neighbor interpolation (no interpolation) - 'bilinear' : Bilinear interpolation (default) - 'bicubic' : Bicubic interpolation These correspond to the modes available in torch.nn.functional.grid_sample. :type interpolation_mode: str, default 'bilinear' :returns: The reprojected image with the same shape as the target image. Pixel values are interpolated from the source image according to the WCS transformation with flux conservation preserved. :rtype: torch.Tensor .. rubric:: Notes This implementation uses a two-step flux conservation approach: 1. Local flux conservation: The image and a "ones" tensor are interpolated together, and the interpolated image is divided by the interpolated ones tensor (footprint) to correct for any flux spreading during interpolation. Areas in the target image that map outside the source image boundaries will be filled with zeros (using 'zeros' padding_mode). This method is particularly suitable for high-precision photometry with extended sources, as it properly preserves both the background noise characteristics and the flux distribution of sources. .. py:function:: calculate_reprojection(source_hdus: Union[astropy.io.fits.PrimaryHDU, List[astropy.io.fits.PrimaryHDU]], target_wcs: astropy.wcs.WCS, shape_out: Tuple[int, int], order='nearest', device=None) Reproject an astronomical image from a source WCS to a target WCS. This high-level function provides a convenient interface for image dfreproject, handling all the necessary steps: WCS extraction, tensor creation, and interpolation. It converts FITS HDU objects to the internal representation, performs the dfreproject, and returns the resulting image as a PyTorch tensor. :param source_hdus: The source image HDU list containing the image data to be reprojected and its associated WCS information in the header. :type source_hdus: Union[PrimaryHDU, List[PrimaryHDU]] :param target_wcs: WCS information fopr the target. :type target_wcs: WCS :param shape_out: Shape of the resampled array :type shape_out: Tuple[int, int] :param order: The interpolation method to use when resampling the source image. Options: - 'nearest' : Nearest neighbor interpolation (fastest, default) - 'bilinear' : Bilinear interpolation (good balance of speed/quality) - 'bicubic' : Bicubic interpolation (highest quality, slowest) :type order: str, default 'nearest' :param device: Device to use for computations. Defaults to GPU if available otherwise uses CPU :type device: str :returns: The reprojected image as a numpy ndarray with the same shape as target_hdu.data. :rtype: numpy.ndarray .. rubric:: Notes This function automatically: - Detects and uses GPU acceleration if available - Handles byte order conversion for tensor creation - Converts data to float32 for processing - Creates WCSHeader objects from FITS headers To save the result as a FITS file, you will need to convert the tensor back to a NumPy array and create a new FITS HDU with the target WCS header. .. rubric:: Examples >>> from astropy.io import fits >>> from astropy.wcs import WCS >>> import torch >>> from dfreproject.reproject import calculate_reprojection >>> >>> # Open source and target images >>> source_hdu = fits.open('source_image.fits')[0] >>> target_hdu = fits.open('target_grid.fits')[0] >>> target_wcs = WCS(target_hdu.header) >>> >>> # Perform dfreproject with bilinear interpolation >>> reprojected = calculate_reprojection( ... source_hdus=source_hdu, ... target_wcs=target_wcs, ... shape_out=target_hdu[0].data.shape, ... order='bilinear' ... ) >>> >>> # Convert back to NumPy and save as FITS >>> reprojected_np = reprojected.cpu().numpy() >>> output_hdu = fits.PrimaryHDU(data=reprojected_np, header=target_hdu.header) >>> output_hdu.writeto('reprojected_image.fits', overwrite=True)