radio-astro-tools / spectral-cube

Library for reading and analyzing astrophysical spectral data cubes
http://spectral-cube.rtfd.org
BSD 3-Clause "New" or "Revised" License
95 stars 61 forks source link

Handling of "K_Tb" or "K (Tb)" BUNIT #862

Closed ashill closed 1 year ago

ashill commented 1 year ago

FITS files from the DRAO Synthesis Telescope, including the Canadian Galactic Plane Survey (https://www.cadc-ccda.hia-iha.nrc-cnrc.gc.ca/en/cgps/) have BUNIT set to K (Tb) or K_Tb. This means Kelvins of brightness temperature, but it's physically just K. SpectralCube.read doesn't know what to make of that unit. This line u.add_enabled_units(u.def_unit(['K (Tb)', 'K_Tb'], represents=u.K)) makes both work. Does it make sense to add that line either to convert_bunit in cube_utils.py or somewhere in astropy.units?

I originally hacked cube_utils.py so that convert_bunit checks for either of those strings in the BUNIT header field and interprets them as K, but the add_enabled_units approach seems cleaner. (convert_bunit already has a hack to check for jy/beam; this would just be an additional check.)

keflavich commented 1 year ago

I don't think this is something we should integrate directly into spectral-cube, but we should integrate it into the docs. Could you post an error message from trying to load with these wrong units? Then we can add to the docs, "If you encounter this message, do this:

from astropy import units as u
from spectral_cube import SpectralCube
u.add_enabled_units(u.def_unit(['K (Tb)', 'K_Tb'], represents=u.K))
cube = SpectralCube.read('DRAO_Cube.fits')

(aside question to self/others: is there a context manager for enabled units?)

ashill commented 1 year ago

It's the error message on line 510 of cube_utils.py:

warnings.warn("Could not parse unit {0}".format(bunit),
    AstropyUserWarning)

so the displayed error message is

WARNING: Could not parse unit K_Tb [spectral_cube.cube_utils]
ashill commented 1 year ago

Could change the warning message with this diff in cube_utils.py. (I can do a pull request if needed.)

         except ValueError:
-            warnings.warn("Could not parse unit {0}".format(bunit),
+            warnings.warn("Could not parse unit {0}. "
+                    "If you know the correct unit, try "
+                    "u.add_enabled_units(u.def_unit(['{0}'], represents=u.correct_unit))".format(bunit),
                           AstropyUserWarning)
keflavich commented 1 year ago

@ashill that would be great!

ashill commented 1 year ago

https://github.com/radio-astro-tools/spectral-cube/pull/863

e-koch commented 1 year ago

Added in #863