integrated example:simbad + splatalogue + alma + spectral_cube

Author: keflavich

docs/gallery.rst

@@ -107,3 +107,48 @@ Finding the mass of a specific planet:
         >>> print findvalue( kepler68b, 'mass')
         0.02105109
 
+Example 6
++++++++++
+
+Grab some data from ALMA, then analyze it using the Spectral Cube package after
+identifying some spectral lines in the data.
+
+.. code-block:: python
+
+   from astroquery.alma import Alma
+   from astroquery.splatalogue import Splatalogue
+   from astroquery.simbad import Simbad
+   from astropy import units as u
+   from astropy import constants
+   from spectral_cube import SpectralCube
+
+   m83table = Alma.query_object('M83', public=True)
+   m83urls = Alma.stage_data(m83table['Asdm_uid'])
+   m83files = Alma.download_and_extract_files(m83urls['URL'])
+   # Sometimes there can be duplicates
+   m83files = list(set(m83files))
+   
+   Simbad.add_votable_fields('rvel')
+   m83simbad = Simbad.query_object('M83')
+   rvel = m83simbad['RVel_Rvel'][0]*u.Unit(m83simbad['RVel_Rvel'].unit)
+
+   for fn in m83files:
+       if 'line' in fn:
+           cube = SpectralCube.read(fn)
+           # Convert frequencies to their rest frequencies
+           frange = u.Quantity([cube.spectral_axis.min(),
+                                cube.spectral_axis.max()]) * (1+rvel/constants.c)
+
+           # Query the top 20 most common species in the frequency range of the
+           # cube with an upper energy state <= 50K
+           lines = Splatalogue.query_lines(frange[0], frange[1], top20='top20',
+                                           energy_max=50, energy_type='eu_k',
+                                           only_NRAO_recommended=True)
+           lines.pprint()
+
+           vcube = cube.with_spectral_unit(u.km/u.s,
+                                           rest_value=lines['Freq-GHz'][0]*u.GHz,
+                                           velocity_convention='radio')
+           # Ugly, I know, but it's a clever python trick...
+           linename = "{Species}{Resolved QNs}".format(**dict(zip(row.colnames,row.data)))
+           vcube.write('M83_ALMA_{linename}.fits'.format(linename=linename))