spectralindex
Return to directory of Todd's CASA extensions
This function is designed to fit a power-law spectral index to the results
output from casa's fluxscale task, OR from a simple list of frequencies,
flux densities and uncertainties. In the former use case, currently it also requires the output
from the listobs task to determine the center frequencies of each
spectral window. It has two options for determining the uncertainty of
the fitted parameters: using the covariance matrix (obtained by setting
trials=0) or running a Monte Carlo simulation (default, with trials=2000). The
latter appears to be more accurate in most cases. It produces a plot summarizing the results.
Finally, if you have only one flux density measurement and a spectral index, you can
use this function to ask what the flux density would be at a different frequency.
Usage:
au.spectralindex(self, filename='', yfilename='', source='', verbose=False, maxpoints=0,
trials=2000, spw='', plotdir='', labelspw=False, referenceFrame='TOPO',
plaintext=False, lineNumbers=None,columns=None, yscale=1.0,
plotunits='linear', freqs=[], fluxes=[], errors=[], plotfile='',
showplot=True, xaxis=[], spectralIndex=None, silent=False,
axisEqual=True)
- filename: contains a listobs output file
- yfilename: contains a fluxscale output file
- source: sourcename to choose from the (possibly) multi-source fluxscale file
- maxpoints: the maximum number of spws to select for the fit (0=no max)
- trials: if > 0, use a Monte-Carlo technique estimate the fit uncertainties, otherwise, use the sqrt(covar_matrix) from scipy.optimize.leastsq (default) There is a minimum number of 100 trials, and ~1000 is recommended.
- spw: the spws to use, e.g. ''=all, '1~3,5,6~8'=[1,2,3,5,6,7,8]
- plotdir: the directory in which to write the plotfile
- labelspw: draw spw numeric labels adjacent to each point
- referenceFrame: the frequency reference frame, if not 'TOPO'
- plaintext: if True, then read the freqs, flux densities and uncertainties from a single plain text file, specified by filename
- lineNumbers: which lines to read, where 1 is the first line in the file
- columns: which columns to read for freq, flux, error (starting at 1)
- yscale: factor to scale the y values after reading from the file
- plotunits: 'linear' or 'log' (only used if trials==0)
- freqs,fluxes,errors: alternative approach to using filename, just pass 3 lists
- xaxis: a list of points for which to compute model y values
- spectralindex: If you have only one flux measurement, you can specify it along with a spectral index (via this parameter) and it will compute expected flux densities at the frequencies specified by the xaxis parameter.
- axisEqual: if True, then when plotting, set pl.axis('equal')
The
listobs file should look like this. It is keying off the string TOPO, so there can be any number
of leading comment characters on each line, and any number of unrelated lines of text at the top.
It can also have an additional column (the spw name) between the spw and the #Chans, as this
is the new format of listobs output.
SpwID #Chans Frame Ch1(MHz) ChanWid(kHz) TotBW(kHz) Corrs
0 128 TOPO 42738 1000 128000 RR LL
1 128 TOPO 42866 1000 128000 RR LL
2 128 TOPO 42994 1000 128000 RR LL
3 128 TOPO 43122 1000 128000 RR LL
4 128 TOPO 43250 1000 128000 RR LL
5 128 TOPO 43378 1000 128000 RR LL
6 128 TOPO 43506 1000 128000 RR LL
7 128 TOPO 43634 1000 128000 RR LL
8 128 TOPO 44488 1000 128000 RR LL
9 128 TOPO 44616 1000 128000 RR LL
10 128 TOPO 44744 1000 128000 RR LL
11 128 TOPO 44872 1000 128000 RR LL
12 128 TOPO 45000 1000 128000 RR LL
13 128 TOPO 45128 1000 128000 RR LL
14 128 TOPO 45256 1000 128000 RR LL
15 128 TOPO 45384 1000 128000 RR LL
Similarly, the
fluxscale file should look like this. It can contain the results for
multiple sources. It is keying off of the string 'SpW', so there can be any number
of leading comment characters on each line, and any number of unrelated lines of text at the top.
# Opening MS: 10C186_B_Q_d1_3s.ms for calibration.
# Found reference field(s): J1331+3030
# Found transfer field(s): J1717-3342 J1924-2914
# Flux density for J1717-3342 in SpW=0 is: 2.33891 ± 0.0620232 (SNR = 37.7102, N= 44)
# Flux density for J1717-3342 in SpW=1 is: 2.33168 ± 0.0608224 (SNR = 38.3358, N= 44)
# Flux density for J1717-3342 in SpW=2 is: 2.33139 ± 0.0642111 (SNR = 36.3082, N= 44)
# Flux density for J1717-3342 in SpW=3 is: 2.33639 ± 0.0697579 (SNR = 33.4929, N= 44)
# Flux density for J1717-3342 in SpW=4 is: 2.31025 ± 0.0631668 (SNR = 36.5738, N= 44)
# Flux density for J1717-3342 in SpW=5 is: 2.32929 ± 0.0667446 (SNR = 34.8986, N= 44)
# Flux density for J1717-3342 in SpW=6 is: 2.35698 ± 0.0698231 (SNR = 33.7564, N= 44)
# Flux density for J1717-3342 in SpW=7 is: 2.32728 ± 0.0668059 (SNR = 34.8364, N= 44)
# Flux density for J1717-3342 in SpW=8 is: 2.23395 ± 0.0394774 (SNR = 56.588, N= 46)
# Flux density for J1717-3342 in SpW=9 is: 2.23898 ± 0.0404523 (SNR = 55.3487, N= 46)
# Flux density for J1717-3342 in SpW=10 is: 2.20212 ± 0.0408321 (SNR = 53.9311, N= 46)
# Flux density for J1717-3342 in SpW=11 is: 2.18275 ± 0.0392774 (SNR = 55.5727, N= 46)
# Flux density for J1717-3342 in SpW=12 is: 2.17859 ± 0.0388576 (SNR = 56.0659, N= 46)
# Flux density for J1717-3342 in SpW=13 is: 2.16182 ± 0.0405544 (SNR = 53.3066, N= 46)
# Flux density for J1717-3342 in SpW=14 is: 2.16303 ± 0.0413542 (SNR = 52.305, N= 46)
# Flux density for J1717-3342 in SpW=15 is: 2.19681 ± 0.0400311 (SNR = 54.8777, N= 46)
# Flux density for J1924-2914 in SpW=0 is: 14.6731 ± 0.138529 (SNR = 105.921, N= 44)
# Flux density for J1924-2914 in SpW=1 is: 14.5996 ± 0.136761 (SNR = 106.753, N= 44)
# Flux density for J1924-2914 in SpW=2 is: 14.6242 ± 0.146673 (SNR = 99.7061, N= 44)
# Flux density for J1924-2914 in SpW=3 is: 14.6802 ± 0.158394 (SNR = 92.6816, N= 44)
# Flux density for J1924-2914 in SpW=4 is: 14.5259 ± 0.141069 (SNR = 102.97, N= 44)
# Flux density for J1924-2914 in SpW=5 is: 14.6075 ± 0.151812 (SNR = 96.2206, N= 44)
# Flux density for J1924-2914 in SpW=6 is: 14.7261 ± 0.157892 (SNR = 93.2666, N= 44)
# Flux density for J1924-2914 in SpW=7 is: 14.5783 ± 0.14938 (SNR = 97.5919, N= 44)
# Flux density for J1924-2914 in SpW=8 is: 13.9264 ± 0.0748849 (SNR = 185.971, N= 46)
# Flux density for J1924-2914 in SpW=9 is: 13.9368 ± 0.0751502 (SNR = 185.452, N= 46)
# Flux density for J1924-2914 in SpW=10 is: 13.7826 ± 0.0785357 (SNR = 175.494, N= 46)
# Flux density for J1924-2914 in SpW=11 is: 13.7178 ± 0.0778214 (SNR = 176.273, N= 46)
# Flux density for J1924-2914 in SpW=12 is: 13.6833 ± 0.0748109 (SNR = 182.905, N= 46)
# Flux density for J1924-2914 in SpW=13 is: 13.6308 ± 0.0772038 (SNR = 176.555, N= 46)
# Flux density for J1924-2914 in SpW=14 is: 13.6372 ± 0.0792322 (SNR = 172.117, N= 46)
# Flux density for J1924-2914 in SpW=15 is: 13.7463 ± 0.0735882 (SNR = 186.8, N= 46)
Examples:
1. Using a single plain text file:
CASA <4>: au.spectralindex('plain.txt', plaintext=True)
Error-weighted fit: Slope: 2.250+-0.109 Flux D. @ 115.000GHz: 0.233+-0.005
Un-weighted fit: Slope: 2.234 Flux D. @ 115.000GHz: 0.239
Plot saved in ./plain.txt.png
2. Using a listobs output file and a fluxscale output file:
CASA <5>: au.spectralIndex('listobs.txt','fluxscale_d2_new.txt','J17')
Read 16 x values for J1717-3342 (avg=44125.000)
Read 16 y values for J1717-3342 (avg=2.264)
Error-weighted fit: Slope: -1.454+-0.640 Flux D. @ 42.802GHz: 2.362+-0.014
Un-weighted fit: Slope: -1.453 Flux D. @ 42.802GHz: 2.364
Plot saved in ./fluxscale_d2_new.txt.J1717-3342.png
3. Using a list of flux density measurements, and finding the value at another frequency:
CASA <9>: au.spectralindex(freqs=[22.3,43.3,222], fluxes=[.0003,.001,.044], errors=[.00001,.0001,.001], xaxis=[115])
Value at 22.300000 = 0.000295 +- 0.000022 (scaled MAD = 0.000022)
Error-weighted fit: Slope: 2.177+-0.222 Flux D. @ 22.300GHz: 0.295+-0.022 mJy
Un-weighted fit: Slope: 2.198 Flux D. @ 22.300GHz: 0.270 mJy
Value at 115.000000 = 0.010475 +- 0.000451 (scaled MAD = 0.000448)
Out[9]:
{'covar': array([[ 0.00043185, -0.00038152],
[-0.00038152, 0.0015951 ]]),
'intercept': -2.5972443931713598,
'interceptUncertainty': 0.50968131460806199,
'meanOfLogX': 1.7770485779507217,
'spectralIndex': 2.1767855043259088,
'spectralIndexUncertainty': 0.22226760634713091,
'yaxis': [0.010476089728640246],
'yaxisUncertainty': [0.00045094643771614152]}
--
ToddHunter - 2012-03-07