astropy / reproject

Python-based Astronomical image reprojection :milky_way: - maintainer @astrofrog
https://reproject.readthedocs.io
BSD 3-Clause "New" or "Revised" License
107 stars 66 forks source link

Missing non-celestial reprojection routines #278

Open keflavich opened 2 years ago

keflavich commented 2 years ago

The documentation claims "Non-celestial images/data are not currently supported but will be added in future releases." but there's no open issue about this that I can find.

We need to support reprojection between different position-velocity slices.

Cadair commented 2 years ago

In general the interpolation method should be able to handle any reprojection which doesn't change the axis type?

keflavich commented 2 years ago

yeah, I had to do some awful hackery in this example: https://github.com/CentralMolecularZone/DataSetVisualizations/blob/main/FullCMZVisualizations.ipynb

(look down where there are big code blocks)

ramseykarim commented 1 year ago

Hi! I know this is a bit of an older issue, but since it's still open, I wanted to chime in with a related issue I'm having (and please let me know if it's better off in a separate issue).

I have been able to reproject position-velocity images (typically created with pvextractor) onto each other in the past. If s1 and s2 are two PV slices from different spectral lines created with pvextractor over the same on-sky path, I would use the line

s2_reproject = reproject_interp((s2.data, s2.header), WCS(s1.header), shape_out=s1.data.shape, return_footprint=False)

That had worked in the past, and unfortunately I don't have a snapshot of all my reproject, pvextractor, astropy, and numpy versions when that worked. I have since updated everything via pip as of April 2023. Now to reproject PV slices, I have to set their rest frequencies to be equal or else I get all NaNs. I see in @keflavich's example that they might be doing something like that too.

s2.header['RESTFRQ'] = s1.header['RESTFRQ']

Something must have changed (in astropy coordinates I guess?) because my old code that produced good reprojected PV diagrams before my updates now does not work unless I apply the rest frequency fix.

Is there a good reason that the rest frequencies should be the same? All the use cases I can imagine for my work would involve different rest frequencies, i.e. overlaying PV cuts of different spectral lines. If the headers both describe the spectral axes in velocity, then I would expect reproject to use that with no further questions (besides checking/adjusting radio vs optical convention I suppose).

Thanks, really appreciate all the work on this package!