%pip install -r https://github.com/radio-astro-tools/tutorials/raw/master/requirements.txt
%pip install --pre -U astroquery

Collecting spectral-cube (from -r https://github.com/radio-astro-tools/tutorials/raw/master/requirements.txt (line 14))
  Cloning https://github.com/radio-astro-tools/spectral-cube/ to /tmp/pip-install-uif3qivd/spectral-cube_1149857666194f28810b3004ed2042bf
  Running command git clone --filter=blob:none --quiet https://github.com/radio-astro-tools/spectral-cube/ /tmp/pip-install-uif3qivd/spectral-cube_1149857666194f28810b3004ed2042bf

  Resolved https://github.com/radio-astro-tools/spectral-cube/ to commit 8ae6705d7e64739cbb08904efc660dbed98874ab

  Installing build dependencies ... -

 \

 |

 /

 done

  Getting requirements to build wheel ... -

 done

  Preparing metadata (pyproject.toml) ... -

import pylab as pl
pl.rcParams['figure.facecolor'] = 'w'
pl.rcParams['image.origin'] = 'lower'
pl.rcParams['image.interpolation'] = 'none'
pl.rcParams['font.size'] = 16

from dask.diagnostics import ProgressBar
import dask
dask.config.set(scheduler='threads')
pbar = ProgressBar()
pbar.register()

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

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)

cube

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

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)

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'))

far3kpcarm = cube[:,int(ymin):int(ymax)+1,int(xmax):int(xmin)].with_mask(far3kpcarm_mask.T[:,None,:])
f3m0 = far3kpcarm.moment0()

[                                        ] | 0% Completed | 206.11 us

[########################################] | 100% Completed | 109.39 ms

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')

Identification of the far 3 kpc arm in the Dame 2001 CO map¶