error parsing regions with dimensions in degrees, radians, physical/image pixels

[cristobal-sifon]

I created this in the aplpy site (https://github.com/aplpy/aplpy/issues/381) but it seems it may be more relevant here.

using aplpy.FITSFigure.show_regions I can load a region file with the following content just fine:

fk5
circle(10.0,-5.0,10")

but when the units correspond to a letter (i.e., 'd' for degrees, 'r' for radians, 'p' for physical pixels or 'i' for image pixels), the command cannot be parsed. Take a file with the following content

fk5
circle(10.0,-5.0,0.02d)

This file loads without a problem in ds9, but when passed to aplpy.FITSFigure.show_regions it issues the following warning:

/Users/cristobal/anaconda3/lib/python3.6/site-packages/pyregion/ds9_region_parser.py:107: UserWarning: Failed to parse : circle(10.0,-5.0,0.02d)

If this happens with all the entries in my reg file, then the following exception is raised:

---------------------------------------------------------------------------
ValueError                                Traceback (most recent call last)
<ipython-input-9-e910de502d80> in <module>()
      6 regfile = mask.regions()
      7 #img.show_regions(regfile)
----> 8 img.show_regions('test.reg')

/Users/cristobal/anaconda3/envs/pyraf/lib/python2.7/site-packages/aplpy/regions.pyc in show_regions(self, region_file, layer, **kwargs)

/Users/cristobal/anaconda3/envs/pyraf/lib/python2.7/site-packages/aplpy/decorators.pyc in _auto_refresh(f, *args, **kwargs)
     23     mydata.nesting = getattr(mydata, 'nesting', 0) + 1
     24     try:
---> 25         return f(*args, **kwargs)
     26     finally:
     27         mydata.nesting -= 1

/Users/cristobal/anaconda3/envs/pyraf/lib/python2.7/site-packages/aplpy/regions.pyc in show_regions(self, region_file, layer, **kwargs)
     63         """
     64 
---> 65         PC, TC = ds9(region_file, wcs.WCS(self._header).sub([wcs.WCSSUB_CELESTIAL]), **kwargs)
     66 
     67         # ffpc = self._ax1.add_collection(PC)

/Users/cristobal/anaconda3/envs/pyraf/lib/python2.7/site-packages/aplpy/regions.pyc in ds9(region_file, header, zorder, **kwargs)
     94     # read region file
     95     if isinstance(region_file, six.string_types):
---> 96         rr = pyregion.open(region_file)
     97     elif isinstance(region_file, pyregion.ShapeList):
     98         rr = region_file

/Users/cristobal/anaconda3/envs/pyraf/lib/python2.7/site-packages/pyregion/core.pyc in open(fname)
    264     """
    265     region_string = _builtin_open(fname).read()
--> 266     return parse(region_string)
    267 
    268 

/Users/cristobal/anaconda3/envs/pyraf/lib/python2.7/site-packages/pyregion/core.pyc in parse(region_string)
    246     sss2 = _check_wcs(sss1)
    247 
--> 248     shape_list, comment_list = rp.filter_shape2(sss2)
    249     return ShapeList(shape_list, comment_list=comment_list)
    250 

ValueError: need more than 0 values to unpack

N.B. all this is with pyregion 1.2.0 re: https://github.com/aplpy/aplpy/issues/369

[Cadair]

We have some functionality in ndcube which can do this kind of thing for any ape 14 WCS.

[keflavich]

could you link it here?

[Cadair]

https://github.com/sunpy/ndcube/blob/main/ndcube/wcs/wrappers/resampled_wcs.py#L7

I couldn't dig it up on my phone :see_no_evil:

These are things which were planning on being upstreamed to astropy but were instead taken into ndcube to incubate them for a bit.

[keflavich]

Thanks. That looks useful, but reproject doesn't support APE14 now, so we can't use it directly.

[keflavich]

Uh, rather, reproject does support APE14... https://github.com/astropy/reproject/pull/166... but... my example doesn't.

[Cadair]

Why doesn't your example?

[keflavich]

from ndcube.wcs.wrappers import ResampledLowLevelWCS
from spectral_cube import SpectralCube
import reproject


# read whatever random cube I have laying around
cube = SpectralCube.read('G000.00+00.00_H2CO_2pol.fits')[:25,:20,:21]
ww = ResampledLowLevelWCS(cube.wcs, [2,2,1])

reproject.reproject_interp(cube[0].hdu, ww)

gives

Traceback (most recent call last):
  File "<ipython-input-24-b22d4859483b>", line 1, in <module>
    reproject.reproject_interp(cube[0].hdu, ww)
  File "/home/adam/repos/astropy/astropy/utils/decorators.py", line 536, in wrapper
    return function(*args, **kwargs)
  File "/home/adam/anaconda3/envs/python3.9/lib/python3.9/site-packages/reproject/interpolation/high_level.py", line 79, in reproject_interp
    wcs_out, shape_out = parse_output_projection(output_projection, shape_out=shape_out,
  File "/home/adam/anaconda3/envs/python3.9/lib/python3.9/site-packages/reproject/utils.py", line 129, in parse_output_projection
    raise TypeError('output_projection should either be a Header, a WCS '
TypeError: output_projection should either be a Header, a WCS object, or a filename

[keflavich]

It should work: https://github.com/astropy/reproject/blob/43b0d8a4a5641cfbe6adbc3b1f2d7598f2fd5930/reproject/utils.py#L117 but doesn't.

[keflavich]

Oh. reproject wants a high-level WCS, but we only have a low-level WCS. What do we do to wrap this to get a high-level WCS?

[keflavich]

aha, got it:

from ndcube.wcs.wrappers import ResampledLowLevelWCS
from spectral_cube import SpectralCube
import reproject

reproject.reproject_interp(cube[0].hdu, ww)
from astropy.wcs.wcsapi import HighLevelWCSWrapper

# read whatever random cube I have laying around
cube = SpectralCube.read('G000.00+00.00_H2CO_2pol.fits')[:25,:20,:21]
ww = ResampledLowLevelWCS(cube.wcs.celesital, 2)

reproject.reproject_interp(cube[0].hdu, HighLevelWCSWrapper(ww), shape_out=(10,10))

[keflavich]

spectral-cube's reproject only supports FITS WCS right now though

[keflavich]

and ResampledLowLevelWCS doesn't support serialization right now (no to_header methhod)