dysh for GBTIDL Users#

Some key differences#

  • No longer tied to a single input and output file, can load as many SDFITS files as RAM allows

  • No longer tied to a finite number of data containers, you can have as many ScanBlock/Spectrum objects as RAM allows

  • Python allows for well-trackable local and global variables

Rough equivalents#

In the tables below, it is assumed you have executed the following commands in Python:

from dysh.fits import GBTFITSLoad, GBTOnline, GBTOffline

If you are running from the dysh shell, then the above commands are already available.

File I/O and Metadata Operations#

GBTIDL

dysh

Notes

online

sdf = GBTOnline(); sdf = GBTOnline([projID])

Monitor the currently online project or a specific project projID. Load the SDFITS into sdf.

offline

sdf = GBTOffline([projID])

Load data from project/session projID into sdf.

filein, [SDFITS file]

sdfits = GBTFITSLoad( [SDFITSfile] )

This will directly load a SDFITS file into python.

dirin, [SDFITS dir]

sdfits = GBTFITSLoad( [SDFITS dir] )

This will load each SDFITS file in a directory into python, in the case it is a project directory where each file is a different VEGAS bank. (A, B, C, etc.). The implementation is identical to the single SDFITS file load.

fileout

There is no direct translation; in dysh, simply specify a filename with the Spectrum.write or ScanBlock.write functions.

keep

spec.write( [SDFITSfile], format=[fmt] )

Write the spectrum spec to the output SDFITSfile specified, using the format fmt.

summary

sdfits.summary()

In dysh, you can also input a list of scan numbers you would like to see, or modify which SDFITS columns are printed.

list

sdfits.summary( verbose=True )

Using the verbose=True option returns similar results to GBTIDL’s list. It is possible to provide a scan number as scan=x.

header, scan_info

spec.meta

Return meta information about the spectrum spec. In GBTIDL, this returns a text printout. In dysh, this returns a dictionary with all SDFITS keywords and values.

Calibration and Data Retrieval#

GBTIDL

dysh

explanation

select

sdf.select( [col1]=[val1], [col2]=[val2] )

Select subsets of the GBTFITSLoad object sdf based on the constraints on col1 and col2 given by val1 and val2

flag

sdf.flag( [col1]=[val1], [col2]=[val2] )

Create flagging rules of the GBTFITSLoad object sdf based on the constraints on col1 and col2 given by val1 and val2. See flagging examples given in the Recipes.

gettp, [scan_num]

sdfits.gettp( scan=[scan_num] )

Get a total power scan. Returns a ScanBlock.

getps, [scan_num]

sdfits.getps( scan=[scan_num] )

Retrieve and calibrate position-switched data. Returns a ScanBlock.

getfs, [scan_num]

sdfits.getfs( scan=[scan_num] )

Retrieve and calibrate frequency-switched data. Returns a ScanBlock.

getnod, [scan_num]

sdfits.getnod( scan=[scan_num] )

Retrieve and calibrate nodding data. Returns a ScanBlock.

show

spec.plot()

Plot the spectrum spec. Calibration routines in GBTIDL do this automatically.

Spectrum Operations#

GBTIDL

dysh

Notes

baseline, bshape, bsubtract, nfit

spec.baseline( [deg], exclude=[ex_reg],remove=True/False )

Compute and optionally remove a baseline with degree deg and exclusion regions ex_reg.

gsmooth, [int]

spec.smooth( “gauss”, [int] )

Convolve spec with a Gaussian kernel. Default behavior is to decimate, which isn’t true in GBTIDL.

hanning

spec.smooth( “hann” )

Convolve spec with a Hanning window. Default behavior is to decimate, which isn’t true in GBTIDL.

boxcar, [int]

spec.smooth( “box”, [int] )

Convolve spec with a boxcar window. This matches the GBTIDL results only when integer int is an odd number. Default behavior is to decimate, which isn’t true in GBTIDL

add/subtract/divide/scale

spec + other, spec - other, spec/other, spec*other

Perform element-wise math between spectrum spec and other vector/scalar other.

bias, [bias]

spec + bias

Add scalar bias to spectrum spec.

accum

average_spectra(list of spectra)

average_spectra has to be imported first (from dysh.spectra import average_spectrum). Alternatively, spec.average(list of other spectra) will average spec with the spectra in the list.

ave

scanblock.timeaverage(); scan.timeaverage()

Average a set of scans inside ScanBlock scanblock or ScanBase scan.

stats

spec.stats()

Print statistics of spectrum spec.

fshift, vshift, xshift

spec.find_shift(other)

Find the shift between spec and other.

gshift

spec.shift(shift); spec.align_to(other)

Shift spec by shift;Find the shift between spec and other and apply. Equivalent to spec.shift(spec.find_shift(other)).

write_ascii, filename

spec.write(filename, format=”ascii.basic”)

Write the contents of spec to a text file in ascii format with minimum header information.

Side-by-side Examples#

The following are examples in GBTIDL and dysh that produce equivalent results. The dysh examples assume they are being run from the dysh shell, so that the modules imported on startup are available (e.g., from astropy import units as u).

OTF Mapping#

This example is based on the observations under TGBT17A_506_11, the same observations used for the on-the-fly data reduction section of the users guide.

sdfits = GBTFITSLoad("TGBT17A_506_11.raw.vegas")
scan_block = sdfits.getsigref(scan=list(range(14,27)), ref=27, fdnum=0, ifnum=0, plnum=0)
scan_block.write("dysh.fits")
filein,"TGBT17A_506_11.raw.vegas"
fileout,"gbtidl.fits"
for s=14,26,1 do begin &$
    getsigref,s,27,/avgref,/keepints &$
    keep &$
endfor

Extragalactic 21 cm Observations using Position Switching#

This example was borrowed from GBTdocs position switched tutorial. A more thorough dysh example using the same data can be found in example reduction of an HI survey. The GBTIDL version of this example can only be run with a display.

sdfits = GBTFITSLoad("/home/astro-util/HIsurvey/Session02")

tp0 = sdfits.gettp(scan=299, ifnum=0, plnum=0, fdnum=0).timeaverage()
tsys0 = tp0.meta["TSYS"]
tp1 = sdfits.gettp(scan=299, ifnum=0, plnum=1, fdnum=0).timeaverage()
tsys1 = tp1.meta["TSYS"]

sigref0a = sdfits.getsigref(scan=[296], ref=295, ifnum=0, fdnum=0, plnum=0, t_sys=tsys0, units="flux", zenith_opacity=0.08).timeaverage()
sigref0b = sdfits.getsigref(scan=[298], ref=297, ifnum=0, fdnum=0, plnum=0, t_sys=tsys0, units="flux", zenith_opacity=0.08).timeaverage()
sigref0 = sigref0a.average(sigref0b)
sigref1a = sdfits.getsigref(scan=[296], ref=295, ifnum=0, fdnum=0, plnum=1, t_sys=tsys1, units="flux", zenith_opacity=0.08).timeaverage()
sigref1b = sdfits.getsigref(scan=[298], ref=297, ifnum=0, fdnum=0, plnum=1, t_sys=tsys1, units="flux", zenith_opacity=0.08).timeaverage()
sigref1 = sigref1a.average(sigref1b)
sigref = sigref0.average(sigref1)

sigref_smo = sigref.smooth(kernel="gauss", width=100, decimate=0)

region = [[1.402*u.GHz, 1.4045*u.GHz],
          [1.40506*u.GHz, 1.4054*u.GHz],
          [1.4072*u.GHz, 1.4115*u.GHz]]

sigref_smo.baseline(1, model="poly", include=region, remove=True)

kms = u.km/u.s
stats_b = sigref_smo[2000*kms:2500*kms].stats()
stats_r = sigref_smo[3500*kms:4000*kms].stats()

rms = (stats_b["rms"] + stats_r["rms"])/2.

cog = sigref_smo.cog(bchan=100, echan=200)
dirin,'/home/astro-util/HIsurvey/Session02'

freeze      ; keep from plotting on the screen, to speed up the processing

gettp, 299, plnum=0, /quiet
tsys0 = !g.s.tsys
gettp, 299, plnum=1, /quiet
tsys1 = !g.s.tsys

for i=295,297,2 do begin &$
  for p=0,1 do begin &$
    if (p eq 0) then tsys=tsys0[0] else tsys=tsys1[0] &$
    getsigref, i+1, i, plnum=p, tsys=tsys, unit='Jy', /quiet &$
    accum &$
  endfor &$
endfor
ave

gsmooth, 100, /decimate

setxunit, 'GHz'
setx, 1.401, 1.412

show

region = [1.402, 1.4045, 1.40506, 1.4054, 1.4072, 1.4115]
reg1 = xtochan(region)
Nregion, reg1
nfit, 1
bshape
baseline

velo
show
stats, 2000, 2500, ret=mystats
rms1 = mystats.rms
stats, 3500, 4000, ret=mystats
rms2 = mystats.rms

rms = (rms1 + rms2) / 2

gmeasure, 1, 0.5, brange=2815, erange=3142, rms=rms
gmeasure, 1, 0.2, brange=2815, erange=3142, rms=rms