Questions or issues? Please report them on the tutorials GitHub issues page.
Identification of the far 3 kpc arm in the Dame 2001 CO map¶
Authors¶
Summary¶
Simply reproduction of the key result is Dame & Thaddeus (2008) using spectral-cube.
The data are available here:
https://lweb.cfa.harvard.edu/rtdc/CO/download/COGAL_all_interp.fits.gz
In [1]:
%pip install spectral-cube dask numpy matplotlib astropy
%pip install --pre -U astroquery
Requirement already satisfied: spectral-cube in /opt/hostedtoolcache/Python/3.13.15/x64/lib/python3.13/site-packages (0.6.8.dev44+g88f163c7f) Requirement already satisfied: dask in /opt/hostedtoolcache/Python/3.13.15/x64/lib/python3.13/site-packages (2026.8.0) Requirement already satisfied: numpy in /opt/hostedtoolcache/Python/3.13.15/x64/lib/python3.13/site-packages (2.5.2) Requirement already satisfied: matplotlib in /opt/hostedtoolcache/Python/3.13.15/x64/lib/python3.13/site-packages (3.11.1) Requirement already satisfied: astropy in /opt/hostedtoolcache/Python/3.13.15/x64/lib/python3.13/site-packages (8.0.1) Requirement already satisfied: casa-formats-io>=0.1 in /opt/hostedtoolcache/Python/3.13.15/x64/lib/python3.13/site-packages (from spectral-cube) (0.3.0) Requirement already satisfied: joblib>=1.3 in /opt/hostedtoolcache/Python/3.13.15/x64/lib/python3.13/site-packages (from spectral-cube) (1.6.0) Requirement already satisfied: packaging>=19 in /opt/hostedtoolcache/Python/3.13.15/x64/lib/python3.13/site-packages (from spectral-cube) (26.3) Requirement already satisfied: radio-beam>=0.3.5 in /opt/hostedtoolcache/Python/3.13.15/x64/lib/python3.13/site-packages (from spectral-cube) (0.3.10) Requirement already satisfied: setuptools>=62.3.3 in /opt/hostedtoolcache/Python/3.13.15/x64/lib/python3.13/site-packages (from spectral-cube) (84.0.0) Requirement already satisfied: tqdm>=4.64 in /opt/hostedtoolcache/Python/3.13.15/x64/lib/python3.13/site-packages (from spectral-cube) (4.70.0) Requirement already satisfied: click>=8.1 in /opt/hostedtoolcache/Python/3.13.15/x64/lib/python3.13/site-packages (from dask) (8.5.0) Requirement already satisfied: cloudpickle>=3.0.0 in /opt/hostedtoolcache/Python/3.13.15/x64/lib/python3.13/site-packages (from dask) (3.1.2) Requirement already satisfied: fsspec>=2021.09.0 in /opt/hostedtoolcache/Python/3.13.15/x64/lib/python3.13/site-packages (from dask) (2026.7.0) Requirement already satisfied: partd>=1.4.0 in /opt/hostedtoolcache/Python/3.13.15/x64/lib/python3.13/site-packages (from dask) (1.4.2) Requirement already satisfied: pyyaml>=5.4.1 in /opt/hostedtoolcache/Python/3.13.15/x64/lib/python3.13/site-packages (from dask) (6.0.3) Requirement already satisfied: toolz>=0.12.0 in /opt/hostedtoolcache/Python/3.13.15/x64/lib/python3.13/site-packages (from dask) (1.1.0)
Requirement already satisfied: contourpy>=1.0.1 in /opt/hostedtoolcache/Python/3.13.15/x64/lib/python3.13/site-packages (from matplotlib) (1.3.3) Requirement already satisfied: cycler>=0.10 in /opt/hostedtoolcache/Python/3.13.15/x64/lib/python3.13/site-packages (from matplotlib) (0.12.1) Requirement already satisfied: fonttools>=4.28.2 in /opt/hostedtoolcache/Python/3.13.15/x64/lib/python3.13/site-packages (from matplotlib) (4.64.0) Requirement already satisfied: kiwisolver>=1.3.1 in /opt/hostedtoolcache/Python/3.13.15/x64/lib/python3.13/site-packages (from matplotlib) (1.5.1) Requirement already satisfied: pillow>=9 in /opt/hostedtoolcache/Python/3.13.15/x64/lib/python3.13/site-packages (from matplotlib) (12.3.0) Requirement already satisfied: pyparsing>=3 in /opt/hostedtoolcache/Python/3.13.15/x64/lib/python3.13/site-packages (from matplotlib) (3.3.2) Requirement already satisfied: python-dateutil>=2.7 in /opt/hostedtoolcache/Python/3.13.15/x64/lib/python3.13/site-packages (from matplotlib) (2.9.0.post0) Requirement already satisfied: astropy-iers-data>=0.2026.6.22.1.23.34 in /opt/hostedtoolcache/Python/3.13.15/x64/lib/python3.13/site-packages (from astropy) (0.2026.8.31.0.57.9) Requirement already satisfied: pyerfa>=2.0.1.3 in /opt/hostedtoolcache/Python/3.13.15/x64/lib/python3.13/site-packages (from astropy) (2.0.1.5) Requirement already satisfied: locket in /opt/hostedtoolcache/Python/3.13.15/x64/lib/python3.13/site-packages (from partd>=1.4.0->dask) (1.0.0) Requirement already satisfied: six>=1.5 in /opt/hostedtoolcache/Python/3.13.15/x64/lib/python3.13/site-packages (from python-dateutil>=2.7->matplotlib) (1.17.0) Requirement already satisfied: scipy>=1.8 in /opt/hostedtoolcache/Python/3.13.15/x64/lib/python3.13/site-packages (from radio-beam>=0.3.5->spectral-cube) (1.18.1)
Note: you may need to restart the kernel to use updated packages.
Requirement already satisfied: astroquery in /opt/hostedtoolcache/Python/3.13.15/x64/lib/python3.13/site-packages (0.4.12.dev10885)
Requirement already satisfied: numpy>=1.22 in /opt/hostedtoolcache/Python/3.13.15/x64/lib/python3.13/site-packages (from astroquery) (2.5.2) Requirement already satisfied: astropy>=5.0 in /opt/hostedtoolcache/Python/3.13.15/x64/lib/python3.13/site-packages (from astroquery) (8.0.1) Requirement already satisfied: requests>=2.26 in /opt/hostedtoolcache/Python/3.13.15/x64/lib/python3.13/site-packages (from astroquery) (2.34.2) Requirement already satisfied: beautifulsoup4>=4.10 in /opt/hostedtoolcache/Python/3.13.15/x64/lib/python3.13/site-packages (from astroquery) (4.15.0) Requirement already satisfied: html5lib>=0.999 in /opt/hostedtoolcache/Python/3.13.15/x64/lib/python3.13/site-packages (from astroquery) (1.1) Requirement already satisfied: keyring>=23.0 in /opt/hostedtoolcache/Python/3.13.15/x64/lib/python3.13/site-packages (from astroquery) (25.7.0) Requirement already satisfied: pyvo>=1.5 in /opt/hostedtoolcache/Python/3.13.15/x64/lib/python3.13/site-packages (from astroquery) (1.9.1) Requirement already satisfied: astropy-iers-data>=0.2026.6.22.1.23.34 in /opt/hostedtoolcache/Python/3.13.15/x64/lib/python3.13/site-packages (from astropy>=5.0->astroquery) (0.2026.8.31.0.57.9) Requirement already satisfied: packaging>=25.0 in /opt/hostedtoolcache/Python/3.13.15/x64/lib/python3.13/site-packages (from astropy>=5.0->astroquery) (26.3) Requirement already satisfied: pyerfa>=2.0.1.3 in /opt/hostedtoolcache/Python/3.13.15/x64/lib/python3.13/site-packages (from astropy>=5.0->astroquery) (2.0.1.5) Requirement already satisfied: PyYAML>=6.0.0 in /opt/hostedtoolcache/Python/3.13.15/x64/lib/python3.13/site-packages (from astropy>=5.0->astroquery) (6.0.3) Requirement already satisfied: soupsieve>=1.6.1 in /opt/hostedtoolcache/Python/3.13.15/x64/lib/python3.13/site-packages (from beautifulsoup4>=4.10->astroquery) (2.9.2) Requirement already satisfied: typing-extensions>=4.0.0 in /opt/hostedtoolcache/Python/3.13.15/x64/lib/python3.13/site-packages (from beautifulsoup4>=4.10->astroquery) (4.16.0) Requirement already satisfied: six>=1.9 in /opt/hostedtoolcache/Python/3.13.15/x64/lib/python3.13/site-packages (from html5lib>=0.999->astroquery) (1.17.0) Requirement already satisfied: webencodings in /opt/hostedtoolcache/Python/3.13.15/x64/lib/python3.13/site-packages (from html5lib>=0.999->astroquery) (0.6.1) Requirement already satisfied: SecretStorage>=3.2 in /opt/hostedtoolcache/Python/3.13.15/x64/lib/python3.13/site-packages (from keyring>=23.0->astroquery) (3.5.0) Requirement already satisfied: jeepney>=0.4.2 in /opt/hostedtoolcache/Python/3.13.15/x64/lib/python3.13/site-packages (from keyring>=23.0->astroquery) (0.9.0) Requirement already satisfied: jaraco.classes in /opt/hostedtoolcache/Python/3.13.15/x64/lib/python3.13/site-packages (from keyring>=23.0->astroquery) (3.4.0) Requirement already satisfied: jaraco.functools in /opt/hostedtoolcache/Python/3.13.15/x64/lib/python3.13/site-packages (from keyring>=23.0->astroquery) (4.6.0) Requirement already satisfied: jaraco.context in /opt/hostedtoolcache/Python/3.13.15/x64/lib/python3.13/site-packages (from keyring>=23.0->astroquery) (6.1.2) Requirement already satisfied: charset_normalizer<4,>=2 in /opt/hostedtoolcache/Python/3.13.15/x64/lib/python3.13/site-packages (from requests>=2.26->astroquery) (3.5.1) Requirement already satisfied: idna<4,>=2.5 in /opt/hostedtoolcache/Python/3.13.15/x64/lib/python3.13/site-packages (from requests>=2.26->astroquery) (3.19) Requirement already satisfied: urllib3<3,>=1.26 in /opt/hostedtoolcache/Python/3.13.15/x64/lib/python3.13/site-packages (from requests>=2.26->astroquery) (2.7.0) Requirement already satisfied: certifi>=2023.5.7 in /opt/hostedtoolcache/Python/3.13.15/x64/lib/python3.13/site-packages (from requests>=2.26->astroquery) (2026.7.22) Requirement already satisfied: cryptography>=2.0 in /opt/hostedtoolcache/Python/3.13.15/x64/lib/python3.13/site-packages (from SecretStorage>=3.2->keyring>=23.0->astroquery) (50.0.1) Requirement already satisfied: cffi>=2.0.0 in /opt/hostedtoolcache/Python/3.13.15/x64/lib/python3.13/site-packages (from cryptography>=2.0->SecretStorage>=3.2->keyring>=23.0->astroquery) (2.1.1) Requirement already satisfied: pycparser in /opt/hostedtoolcache/Python/3.13.15/x64/lib/python3.13/site-packages (from cffi>=2.0.0->cryptography>=2.0->SecretStorage>=3.2->keyring>=23.0->astroquery) (3.0) Requirement already satisfied: more-itertools in /opt/hostedtoolcache/Python/3.13.15/x64/lib/python3.13/site-packages (from jaraco.classes->keyring>=23.0->astroquery) (11.1.0)
Note: you may need to restart the kernel to use updated packages.
In [2]:
import numpy as np
from spectral_cube import SpectralCube
import radio_beam
from astropy import units as u
from astropy import convolution
from astropy.coordinates import SkyCoord
/opt/hostedtoolcache/Python/3.13.15/x64/lib/python3.13/site-packages/tqdm/auto.py:21: TqdmWarning: IProgress not found. Please update jupyter and ipywidgets. See https://ipywidgets.readthedocs.io/en/stable/user_install.html from .autonotebook import tqdm as notebook_tqdm
In [3]:
import pylab as pl
pl.rcParams['figure.facecolor'] = 'w'
pl.rcParams['image.origin'] = 'lower'
pl.rcParams['image.interpolation'] = 'none'
pl.rcParams['font.size'] = 16
In [4]:
# Enable dask for this example
from dask.diagnostics import ProgressBar
import dask
dask.config.set(scheduler='threads')
pbar = ProgressBar()
pbar.register()
We load up the data and append the 1/8 degree beam (there's no beam in the original data header):
In [5]:
from astropy.utils.data import download_file
datafile = download_file(
'https://lweb.cfa.harvard.edu/rtdc/CO/download/COGAL_all_interp.fits.gz',
cache=True, show_progress=True)
from spectral_cube import SpectralCube
cube = SpectralCube.read(datafile, use_dask=True)
cube = cube.with_beam(radio_beam.Beam(1/8*u.deg))
# cube = cube.with_spectral_unit(u.km / u.s)
In [6]:
cube
Out[6]:
DaskSpectralCube with shape=(493, 481, 2881) and unit=K and chunk size (322, 322, 322): n_x: 2881 type_x: GLON-CAR unit_x: deg range: nan deg: nan deg n_y: 481 type_y: GLAT-CAR unit_y: deg range: nan deg: nan deg n_s: 493 type_s: VOPT unit_s: m / s range: -319.800 m / s: 319.800 m / s
We then create a mask based on the best-fit linear arm model from their Figure 1:
(note that the header incorrectly specifies the velocity unit in m/s instead of km/s)
In [7]:
lmin,lmax = -12, 12.5
def v_of_l(ell):
return (56 + 4.0 * ell)*u.m/u.s
far3kpcarm_mask = [(cube.spectral_axis > v_of_l(l)-13*u.m/u.s) &
(cube.spectral_axis < v_of_l(l)+13*u.m/u.s)
for l in np.arange(lmax, lmin-0.01, -1/8.)]
far3kpcarm_mask = np.array(far3kpcarm_mask)
Find the corners to cut out of the image to match Fig 3:
In [8]:
xmin,ymin = cube.wcs.celestial.world_to_pixel(SkyCoord(lmin*u.deg,-2.5*u.deg,frame='galactic'))
xmax,ymax = cube.wcs.celestial.world_to_pixel(SkyCoord(lmax*u.deg,2.5*u.deg,frame='galactic'))
Calculate the moment map:
In [9]:
far3kpcarm = cube[:,int(ymin):int(ymax)+1,int(xmax):int(xmin)].with_mask(far3kpcarm_mask.T[:,None,:])
f3m0 = far3kpcarm.moment0()
[ ] | 0% Completed | 248.06 us
[########################################] | 100% Completed | 106.85 ms
In [10]:
pl.figure(figsize=(15,10))
ax = pl.subplot(projection=far3kpcarm.wcs.celestial)
pl.imshow(f3m0.value, vmax=5, vmin=0, origin='lower', cmap='gray_r', interpolation='none',)
cm = pl.cm.inferno.copy()
cm.set_bad((0,0,0,0))
pl.imshow(np.where(f3m0.value > 5, f3m0.value, np.nan), vmax=25, vmin=5, origin='lower', interpolation='none',
cmap=cm)
ax.coords[0].set_axislabel('Galactic Longitude')
ax.coords[1].set_axislabel('Galactic Latitude')
Compare to theirs:
