Skip to content

Reproject unable to handle CompImageHDU #538

@SpacialTree

Description

@SpacialTree

I was trying to reproject some images to mosaic them together, and encountered an issue where trying to reproject files packed into .fits.fz results in this error: ValueError: mmap length is greater than file size

Code to reproduce:

hdu1 = fits.open('file1.fits.fz')[1]
hdu2 = fits.open('file2.fits.fz')[1]

hdu_list = [hdu1, hdu2]
wcs_out, shape_out = find_optimal_celestial_wcs(hdu_list)

hdu_reproj, footprint = reproject_interp(hdu1, wcs_out, shape_out=shape_out)

Using a tuple of the data and header instead of the HDU works:

hdu_reproj, footprint = reproject_interp((hdu1.data, hdu1.header), wcs_out, shape_out=shape_out)

@keflavich thinks it could be a bug with how hdu_to_numpy_memmap handles CompImageHDU

Full traceback:

---------------------------------------------------------------------------
ValueError                                Traceback (most recent call last)
Cell In[23], [line 1](vscode-notebook-cell:?execution_count=23&line=1)
----> hdu_reproj, footprint = reproject_interp(hdu, wcs_out, shape_out=shape_out)

File /red/adamginsburg/repos/reproject/reproject/interpolation/high_level.py:126, in reproject_interp(input_data, output_projection, shape_out, hdu_in, order, roundtrip_coords, output_array, output_footprint, return_footprint, block_size, parallel, return_type)
     17 def reproject_interp(
     18     input_data,
     19     output_projection,
   (...)
     29     return_type=None,
     30 ):
     31     """
     32     Reproject data to a new projection using interpolation (this is typically
     33     the fastest way to reproject an image).
   (...)
    123         indicate valid values.
    124     """
--> [126](https://vscode-remote+tunnel-002bc0703a-002ds9ufhpc.vscode-resource.vscode-cdn.net/red/adamginsburg/repos/reproject/reproject/interpolation/high_level.py:126)     array_in, wcs_in = parse_input_data(input_data, hdu_in=hdu_in)
    127     wcs_out, shape_out = parse_output_projection(
    128         output_projection, shape_in=array_in.shape, shape_out=shape_out, output_array=output_array
    129     )
    131     if isinstance(order, str):

File /red/adamginsburg/repos/reproject/reproject/utils.py:128, in parse_input_data(input_data, hdu_in, source_hdul)
    126     return parse_input_data(input_data[hdu_in], source_hdul=input_data)
    127 elif isinstance(input_data, PrimaryHDU | ImageHDU | CompImageHDU):
--> [128](https://vscode-remote+tunnel-002bc0703a-002ds9ufhpc.vscode-resource.vscode-cdn.net/red/adamginsburg/repos/reproject/reproject/utils.py:128)     return (hdu_to_numpy_memmap(input_data), WCS(input_data.header, fobj=source_hdul))
    129 elif isinstance(input_data, tuple) and isinstance(input_data[0], np.ndarray | da.core.Array):
    130     if isinstance(input_data[1], Header):

File /red/adamginsburg/repos/reproject/reproject/utils.py:95, in hdu_to_numpy_memmap(hdu)
     86 if (
     87     hdu.header.get("BSCALE", 1) != 1
     88     or hdu.header.get("BZERO", 0) != 0
   (...)
     91     or hdu.fileinfo()["file"].compression is not None
     92 ):
     93     return hdu.data
---> [95](https://vscode-remote+tunnel-002bc0703a-002ds9ufhpc.vscode-resource.vscode-cdn.net/red/adamginsburg/repos/reproject/reproject/utils.py:95) return np.memmap(
     96     hdu.fileinfo()["file"].name,
     97     mode="r",
     98     dtype=hdu.data.dtype.newbyteorder(">"),
     99     shape=hdu.data.shape,
    100     offset=hdu.fileinfo()["datLoc"],
    101 )

File /orange/adamginsburg/miniconda3/envs/python312/lib/python3.12/site-packages/numpy/core/memmap.py:268, in memmap.__new__(subtype, filename, dtype, mode, offset, shape, order)
    266 bytes -= start
    267 array_offset = offset - start
--> [268](https://vscode-remote+tunnel-002bc0703a-002ds9ufhpc.vscode-resource.vscode-cdn.net/orange/adamginsburg/miniconda3/envs/python312/lib/python3.12/site-packages/numpy/core/memmap.py:268) mm = mmap.mmap(fid.fileno(), bytes, access=acc, offset=start)
    270 self = ndarray.__new__(subtype, shape, dtype=descr, buffer=mm,
    271                        offset=array_offset, order=order)
    272 self._mmap = mm

ValueError: mmap length is greater than file size

Metadata

Metadata

Assignees

No one assigned

    Labels

    No labels
    No labels

    Type

    No type

    Projects

    No projects

    Milestone

    No milestone

    Relationships

    None yet

    Development

    No branches or pull requests

    Issue actions