Repository navigation
WCS.world_to_array_index_values ignores units of inputs #11173
Description
Description
When equivalent inputs are given to WCS.world_to_array_index_values but in different units, different array indices are returned. Presumably this is because the unit information is being ignored and it is being assumed that the inputs are in the same unit as stored in the WCS.
Expected behavior
Equivalent inputs in different units produce same output.
Actual behavior
Steps to Reproduce
>>> import astropy.units as u
>>> import numpy as np
>>> from astropy.wcs import WCS
>>> wcs = WCS()
>>> wcs.wcs.ctype = 'HPLN-TAN', 'HPLT-TAN'
>>> wcs.wcs.crval = 0.00089530541880571, 0.00038493926472938695
>>> wcs.wcs.crpix = 512.5, 512.5
>>> wcs.wcs.pc = np.array([[ 0.99999706, 0.00242302], [-0.00242302, 0.99999706]])
>>> wcs.wcs.cdelt = 0.0006674422222222221, 0.0006674422222222221
>>> wcs.wcs.cunit = ['deg', 'deg']
>>> # Define inputs
>>> a = [-900., -900.] * u.arcsec
>>> b = [0., -200.] *u.arcsec
>>> # The following two commands should return the same result, but don't.
>>> wcs.world_to_array_index_values(a.to(u.deg), b.to(u.deg))
(array([511, 428]), array([136, 136]))
>>> wcs.world_to_array_index_values(a.to(u.arcsec), b.to(u.arcsec))
(array([-9223372036854775808, 31755]),
array([-9223372036854775808, 510]))
Activity
pllim commentedon Dec 17, 2020
Are you sure it takes Quantity? I can't tell from the API doc. On https://docs.astropy.org/en/latest/api/astropy.wcs.wcsapi.BaseWCSWrapper.html , it says
world_to_array_index_valuesis the same asworld_to_pixel_values, and thenworld_to_pixel_valuessays it "takesworld_n_dimscalars or arrays as input in units given byworld_axis_units," which is a separate attribute. 🤷DanRyanIrish commentedon Dec 17, 2020
Ah. So you're suggesting this is deliberate behavior, not a bug? Is there a deliberate reason behind that? Without knowing the history and context I would think it should take a quantity and convert to the relevant unit. If there is a reason it shouldn't however, then perhaps an error or warning should be raised if the user enters a quantity. Otherwise it's very confusing for a user who understandably might expect the function to convert the quantity before the operation.
pllim commentedon Dec 17, 2020
I will have to defer to the authors of APE 14 -- @astrofrog @eteq @Cadair @keflavich
keflavich commentedon Dec 17, 2020
I suspect there is a bug here, in that this low-level API function being called (
world_to_array_index_values) should not accept quantities at all. Low-level WCS should (probably) only accept floats and should explicitly declare that the floats should havedegunits. Or they should do the unit conversion, but that to me seems to be in scope for the high-level API.@DanRyanIrish you may want to use
world_to_array_index, leaving off the_values, to use the high-level API instead.DanRyanIrish commentedon Dec 17, 2020
Hi @keflavich. Thanks for your thoughts. My specific use case in within ndcube and the
_valuesmethod is the correct one to use. So it's a matter of deciding where the unit handling should be done, in astropy or ndcube. Whichever direction astropy wants to go (perform unit conversion or explicitly not accept Quantities) is fine so long as there's clarity. That way ndcube can proceed with confidence.keflavich commentedon Dec 17, 2020
Because the low-level WCS is the FITS WCS here, I suspect it can only handle values in
degunits, so I suggest going with that forndcubedevelopment. However, we still should leave this issue open, as the_valuesmethod needs to either raise an error or properly treat quantities.Cadair commentedon Dec 17, 2020
Surely this is only true for
astropy.wcs.WCSand not a general APE 14 low level WCS.