Determining the velocity scale for SMA spectra, and getting data into CASA
Unlike ALMA and VLA, the SMA does not observe in the observatory's topocentric frame, but rather employs doppler tracking of a single specified line. It can be said to be in a fixed frame that matches "topocentric at source transit", which can be translated to closely approximate LSRK for that one spectral line. As such, the raw data cannot be directly read into CASA, as it is neither topocentric nor LSRK. While SMA data can be read in as LSRK and subsequently imaged in that single line, the uvdata data cannot be correctly converted to a proper topocentric frame inside CASA. Thus, other lines cannot be correctly imaged with a velocity scale that is not shifted and smeared. This page describes the only method that I have found to convert SMA data to proper topocentric frame, along with other attemtped methods that did not work.
1. SWARM correlator (via MIR = recommended method)
- If you are only interested in a single line that was doppler-tracked online:
- UVFITS files written by MIR encode the frequency frame as 257 (which means LSRK RADIO)
- Using data from project 2018B-A020, I found that the default of importuvfits is LSRK, so there should be no need to run au.setMeasFreqRef to change it
- If there are multiple sources in the project, MIR has a function (uti_doppler_fix) to convert the data scan-by-scan to track a different source. Then another set of UVFITS files can be created.
- If you are interested in more than one line, then you can convert the data to TOPO in MIR and output one set of UVFITS files that can be used to image any source correctly in CASA.
- Load the data into MIR, one receiver band at a time
- apply_tsys as normal
- Run uti_to_topo_complete.pro (which uses the IDL INTERPOLATE function and my sub function transit_time.pro)
- It also requires ASCII listings of several variables during the track, which can be obtained from the SMAOC engineering database, as summarized in this screenshot
- Run autofits one-source-at-a-time to write out all the sideband/spw combinations for that source as individual FITS files
- Run toddTools.importsma (which uses MIRtoCASA.py to convert these to measurement sets)
- Run toddTools.concatSMAdataset to produce a single ms which is properly delineated by scan
- In previous version: uti_to_topo.pro, I had to then run toddTools.shiftSpectralAxis on each spw, using the velocity = IFLO_VRADIAL(transit) - V_LSRK + IFLO_EROTATE(transit) (-34.09967 - -15.9 + -0.3602725) = -18.560 km/s (relative to the center frequency of each spw), but this is now done in uti_to_topo_complete (by calling transit_time.pro)
- Background info (email from Todd Hunter to Mark Gurwell on Oct 4, 2019):
- After spending another day on it, we have worked out the exact correction needed, and now the SMA spectra in CASA align extremely well with the ALMA spectra (to better than 0.05 km/s). We checked the method on all 3 lines still visible in the source at the time of the ALMA observations (1 in LSB and 2 in USB in different windows). In summary, in addition to removing the time-variable part of the LO, and the difference between the TOPO velocity at transit (IFLO_VRADIAL) and the V_LSRK that was dopplerTracked (-15.9 km/s), we must also remove the value of IFLO_EROTATE at transit (-0.3602725 km/s). This puts the data into a true TOPO frame, fixed throughout the observation. CASA then handles all the terms needed for tclean (or cvel) to convert to LSRK (including the earth orbital motion term) just as it does for ALMA and VLA data taken in this frame. We are proceeding with calibration and imaging in CASA.
2. SWARM correlator (via Taco's sma2casa.py)
- In CASA 5.6, importfitsidi is labeling the data as frame GEO by default (this has been true ever since CASA 5.3 due to CAS-10790), which is close to TOPO, i.e. differing only by the diurnal correction.
- Plotting the subsequent raw uvspectra (in this case the 199.574 GHz line) has the spectral lines appearing at the expected TOPO frequency and TOPO velocity (-34.09967), as computed by NRAO Dopset for Mauna Kea on LST day 65367 (see dopset2 in Todd's uap directory).
- Transforming the data using plotms to LSRK frame does put the line at the expected LSRK velocity (-15.9 kms), but causes odd features, likely due to inappropriate removal of the diurnal correction, because SMA already removed it during observing (although it was wrongly done prior to April 2019, Haystack was being used as the reference position instead of Mauna Kea)
- The correct solution: change smaImportFix.py to use the new specframe parameter of importfitsidi (introduced in CASA 5.3) to 'LSRK' (or use au.setMeasFreqRef(vis, newcode=1) to relabel a previously-imported ms), and then either:
- Apply a velocity shift of -dsm_as_iflo_vradial_d/1000 + V_LSRK (in this case: +34.09967 - 15.90 = +18.20 km/s) to the resulting cubes
- or adjust the spw frequency scale in the uvdata (one correction per spw) prior to calibration and imaging, using tt.shiftSpectralAxis(vis, spw, shift_kms), which produces the correct-looking spectrum below:
- Potential trap to avoid: one can alternatively set spectframe='TOPO' in importfitsidi, which should recover the behavior in CASA version < 5.3. If you do this, and then let plotms transform to LSRK, you get this promising result.
- The x-axis becomes finely sampled as the topocentric source velocity changes during the track, which is similar to the behavior of ALMA data which is always observed in TOPO:
- However, this is not correct! because SMA data were not observed in TOPO and should not behave this way. Indeed, this line profile does not match the profile observed by ALMA, nor the profile from the same dataset in MIR/IDL (exported in the LSRK frame via UVFITS to CASA).
- As an aside, the regular pattern of low points in all of these SMA profiles above is due to low coherence in the early part of the SMA track (scans < 18), and disappear if you exclude that timerange.
- I also checked that the high frequency receiver data behaves the same with sma2casa.py (in this case the 240 Rx USB with the 247.840 GHz maser line, where I have also avoided plotting the noisy antenna SMA7):
- However, any line that was not the doppler-tracked line will still have a residual doppler shift that has not been corrected by this method, which is time-dependent and will smear the line profile at the few 0.1 km/s level.
3. SWARM correlator (via Jun-Hui's SWARM2CASA)
- No experience with this code yet.
Table of results on the internal wiki
- The SMA doppler tracks a source velocity using the LSRK system. The call to slalib is this one.
- The SMA frequency scale is topocentric at source transit.
- If you load the data in miriad, then immediately write it out as uvfits, then load it into CASA for calibration and imaging, then you are fine. CASA will interpret the spectral axis correctly and lines will appear at their rest frequency (for the LSR velocity that was doppler-tracked).
- However, if you load and calibrate the data in miriad, export the spectra using frequency, and read it into another package (AIPS or CASA), I've found that you need to shift by the total radial velocity being tracked at transit if you want the lines to appear at their rest frequency (for the LSR velocity that was doppler-tracked). Fortunately, this radial velocity is written to the engineering database, but only reliably since 2006. These data are older than that, but the value can be computed with NRAO dopset.
- For each SMA track, you need to query the database for the velocity being tracked at source transit. The command to issue for recent data (>2006) is:
- obscon.sma.hawaii.edu% databaseValue -m hal9000:m5 dsm_as_iflo_vradial_d -s "2010-02-23 15:10:00"
- Older dates will simply report the first stored value, which is 10793.4
- Alternatively, you can use the engineering data GUI and select DSM_AS_IFLO_VRADIAL_D
- NRAO dopset takes the siderealDay as input (the online version for the VLA site is here):
- For Mauna Kea, just run dopset2 (available in uap area on Todd's workstation in Charlottesville).
- siderealDay = floor(JD*(1.002737909350795)-2400000.5) + LST/24.
- = floor(JD*(1.002737909350795)-2400000.5) + 0.7222
- average shift for LSB for 2004+2005 may = 20.62 MHz
- Examining spectra after conversion to class
- Pointing center in oldest data: 17:20:53.387 -35:47:01.63
- SMA1: 53.44, 57.9: offset from pointing center = +0.053s, +3.73" = (+0.645", +3.73")
- SMA2: 53.2, 59.6: offset from pointing center = -0.187s, +2.03" = (-2.28", +2.03")
- 2010 I(N) phase center: 17:20:55.0 -d -35:45:07.0
- Frequency shift
- The shift necessary for the combination of 2004-05-27 with 2005-05-01 appears to be the value for 2004-05-27 (i.e. the first dataset in the combination) which is +17.28.
- The equivalent value for USB is +18.07.
- These values can be used in class like so: modify freq 218324.723+17.28, at which point the line will show up at the correct LSR velocity. About -5 to -6 for SMA1 and -8 to -9 for SMA2.
- DBCON of May 12 + May 27, 2004: 216.678460 and 216.67327789, diff = 5.182 MHz, nearly consistent with 23.19-17.28=5.91 MHz.
- DBCON of 2004 + 2005: 216.67327789 and 216.6814054, diff = 8.128 MHz, consistent with 25.42-17.28=8.14MHz
- The famous half-channel labelling error
- Taco's description in SMA log #14506 (2007Nov26): "...all of our data, since day one, has had 1/2 of the chunks labelled with frequencies and velocities 1/2 channel step too high, and the other half have been labelled with frequencies and velocities 1/2 channel step too low. The sign of the shift depends on the sideband used in the 2nd correlator mixing (the third mixing, if you count the mixing at the SIS junction). So chunks s01 and s02 were shifted one way, and s03 and s04 were shifted the other direction. The new version of dataCatcher has been copied into the default area, so data taken from now on should not need this correction."
- If you use miriad smalod version equal or newer than 27-Nov-2007, then it takes care of the problem automatically.
- Mark has worked out a sequence of MIR commands which will correct the labelling problem in old data sets.
- Charlie Qi's relevant log entry #14628 (2007Dec19): Charlie put Mark and Taco's fix into a mir program called uti_fqhdr_fix.pro. For all data taken before November 26th 2007, uti_fahdr_fix can be used to fix this frequency header problem. For data taken after November 26th 2007, the program won't do anything to the data. To use it: IDL> uti_fqhdr_fix. It can be used at any step during calibration. But be careful not to use it more than once.