Converting Jy/sr to K for model data cubes without specified beam

[SmirnGreg]

Hello!

I am getting, to my understanding, a false alarm from SpectralCube.to().

I have modeled data in Jy/sr units with rest frequency in the file header (RESTFREQ) and frequency axis in Hz. The data have no beam as it was never smoothed. Now, I want to convert Jy/sr to K, which seems to me as a trivial problem knowing the frequencies of the channel.

from spectral_cube import SpectralCube
cube = SpectralCube.read("13CO J=3-2_image.fits")
cube

#SpectralCube with shape=(31, 100, 100) and unit=Jy / sr:
# n_x:    100  type_x: RA---TAN  unit_x: deg    range:    71.719863 deg:   71.722334 deg
# n_y:    100  type_y: DEC--TAN  unit_y: deg    range:    16.998856 deg:   17.001219 deg
# n_s:     31  type_s: FREQ      unit_s: Hz     range: 220396846199.799 Hz:220400522044.321 Hz

"""This crashes"""
cube.to(u.K)

# Traceback (most recent call last):
#  File "/home/smirnov/anaconda3/envs/latest/lib/python3.8/site-packages/IPython/core/interactiveshell.py", line 3437, in run_code
#    exec(code_obj, self.user_global_ns, self.user_ns)
#  File "<ipython-input-6-284474804884>", line 1, in <module>
#    cube.to(u.K)
#  File "/home/smirnov/anaconda3/envs/latest/lib/python3.8/site-packages/spectral_cube/utils.py", line 49, in wrapper
#    return function(self, *args, **kwargs)
#  File "/home/smirnov/anaconda3/envs/latest/lib/python3.8/site-packages/spectral_cube/spectral_cube.py", line 2484, in to
#    factor = cube_utils.bunit_converters(self, unit, equivalencies=equivalencies)
#  File "/home/smirnov/anaconda3/envs/latest/lib/python3.8/site-packages/spectral_cube/cube_utils.py", line 663, in bunit_converters
#    raise ValueError("To convert cubes with Jy/beam units, "
#ValueError: To convert cubes with Jy/beam units, the cube needs to have a beam defined.

"""u.Quantity can make this conversion"""
cube.filled().to(u.K, equivalencies=u.brightness_temperature(cube.header["RESTFREQ"] * u.Hz))
#<Quantity [[[0., 0., 0., ..., 0., 0., 0.],
#            [0., 0., 0., ..., 0., 0., 0.],
#            [0., 0., 0., ..., 0., 0., 0.],
#            ...,
#            [0., 0., 0., ..., 0., 0., 0.],
#            [0., 0., 0., ..., 0., 0., 0.],
#            [0., 0., 0., ..., 0., 0., 0.]],
#           [[0., 0., 0., ..., 0., 0., 0.],
#            [0., 0., 0., ..., 0., 0., 0.],
#            [0., 0., 0., ..., 0., 0., 0.],
#            ...,
#            [0., 0., 0., ..., 0., 0., 0.],
#            [0., 0., 0., ..., 0., 0., 0.],
#            [0., 0., 0., ..., 0., 0., 0.]],
#           [[0., 0., 0., ..., 0., 0., 0.],
#            [0., 0., 0., ..., 0., 0., 0.],
#            [0., 0., 0., ..., 0., 0., 0.],
#            ...,
#            [0., 0., 0., ..., 0., 0., 0.],
#            [0., 0., 0., ..., 0., 0., 0.],
#            [0., 0., 0., ..., 0., 0., 0.]],
#           ...,
#           [[0., 0., 0., ..., 0., 0., 0.],
#            [0., 0., 0., ..., 0., 0., 0.],
#            [0., 0., 0., ..., 0., 0., 0.],
#            ...,
#     ...

"""This does not work too"""
cube.to(u.K, equivalencies=u.brightness_temperature(cube.header["RESTFREQ"] * u.Hz))

# Traceback (most recent call last):
#  File "/home/smirnov/anaconda3/envs/latest/lib/python3.8/site-packages/IPython/core/interactiveshell.py", line 3437, in run_code
#    exec(code_obj, self.user_global_ns, self.user_ns)
#  File "<ipython-input-8-b4bbb9513f9e>", line 1, in <module>
#    cube.to(u.K, equivalencies=u.brightness_temperature(cube.header["RESTFREQ"] * u.Hz))
#  File "/home/smirnov/anaconda3/envs/latest/lib/python3.8/site-packages/spectral_cube/utils.py", line 49, in wrapper
#    return function(self, *args, **kwargs)
#  File "/home/smirnov/anaconda3/envs/latest/lib/python3.8/site-packages/spectral_cube/spectral_cube.py", line 2484, in to
#    factor = cube_utils.bunit_converters(self, unit, equivalencies=equivalencies)
#  File "/home/smirnov/anaconda3/envs/latest/lib/python3.8/site-packages/spectral_cube/cube_utils.py", line 663, in bunit_converters
#    raise ValueError("To convert cubes with Jy/beam units, "
# ValueError: To convert cubes with Jy/beam units, the cube needs to have a beam defined.

spectral_cube.__version__
Out[12]: '0.6.0'
astropy.__version__
Out[14]: '5.0.1'

13CO J=3-2_image.zip

The check appears here:

https://github.com/radio-astro-tools/spectral-cube/blob/7d0b9c7b1ef9259fa698af954a1087b24ddf5acf/spectral_cube/cube_utils.py#L661-L664

Can or has_perangarea just be omitted here? If the beam is involved neither in the original cube unit nor in the target unit, would anyone need it for the transformation?

Thanks for checking this.

[SmirnGreg]

Specifying different beams does not change the final value (as expected)

cube1 = cube.with_beam(Beam(1e-10*u.arcsec)).to(u.K)
cube2 = cube.with_beam(Beam(1e-3*u.arcsec)).to(u.K)

np.nanmax((cube1/cube2).filled().value)
#Out[31]: 1.0

np.nanmin((cube1/cube2).filled().value)
#Out[32]: 1.0

[keflavich]

OK, right, the problem is that a cube with Jy/sr units can't convert to K, but it should be able to. I think this is an easy fix. Thanks for identifying it.

[keflavich]

and I think your fix is right; that looks like a logical error. Would you be willing to submit a PR removing or has_perangarea from that statement and add a test of a Jy/sr->K conversion?