Spectralindex: Difference between revisions

From CASA Guides
Jump to navigationJump to search
(Created page with "'''Return to Analysis Utilities''' This page documents the '''spectralindex''' function of Python module analysisUtils. This function is designed to ...")
 
 
(39 intermediate revisions by the same user not shown)
Line 1: Line 1:
'''Return to [[Analysis Utilities]]'''
'''Return to [[Analysis Utilities]]'''


This page documents the '''spectralindex''' function of Python module [[Analysis Utilities|analysisUtils]].
This page documents the '''spectralindex''' function of the linfit class in the Python module [[Analysis Utilities|analysisUtils]].


This function is designed to fit a power-law spectral index to the results output from casa's fluxscale task. Currently, it also requires the output from the listobs task to determine the center frequencies of each spectral window. It runs a brief Monte-Carlo series of fits to determine the uncertainty on the fitted slope on the basis of the error bars on each flux density. Finally, it produces a plot. Note: spectralIndex is a synonym for spectralindex.  
This function is designed to fit a power-law spectral index to the results output from casa's [[fluxscale]] task (using scipy.optimize.leastsq). Currently, it also requires a [[listobs]]-generated ascii file as the first argument in order to determine the center frequencies of each spectral window. (To avoid the need for this file in the future, we have submitted an [https://bugs.nrao.edu/browse/CAS-3770 improvement request] for [[fluxscale]] to also print these frequencies.)  By default, it estimates the uncertainties
of the fitted slope and amplitude from the square root of the diagonals of the covariance matrix.  Alternatively, there is an option to run a Monte-Carlo series of fits to empirically determine the uncertainties. Finally, it produces a plot summarizing the results. By default, it will use all spectral windows, but it will exclude those with an "INSUFFICIENT DATA" status from fluxscale (which can arise in high-frequency EVLA data if you have not split off the pointing spws).  


Usage:
'''Usage:'''


spectralIndex(filename, yfilename, source='', verbose=False, maxpoints=0, trials=2000, spw='', help=False)  
<tt>au.spectralindex(filename, yfilename="", source="", verbose=False, maxpoints=0, trials=0, spw="", plotdir="", labelspw=False, referenceFrame='TOPO', plaintext=False, lineNumbers=None, columns=None, yscale=1.0, plotunits='linear') </tt>


     filename: contains a listobs output file
     filename: contains a listobs output file
     yfilename: contains a fluxscale output file
     yfilename: contains a fluxscale output file
     source: sourcename to choose from the (possibly) multi-source fluxscale 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)
     maxpoints: the maximum number of spws to select for the fit starting from the first one (0=no max)
     trials: number of Monte-Carlo fits to run to estimate the fit uncertainties
     trials: if > 0, use a Monte-Carlo technique estimate the fit uncertainties,
     spw: the spws to use, e.g. ''=all, '1~3,5,6~8'=[1,2,3,5,6,7,8]
      otherwise, use the sqrt(covarMatrix) 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 (in case you don't have permission in the pwd)
    labelspw: draw spw numeric labels adjacent to each point
    referenceFrame: the frequency reference frame in the listobs file, if not 'TOPO'
    plaintext: if True, then read the freqs, flux densities and uncertainties from a single plain text file,
      specified by filename, with one measurement per line with a format of
      "frequency_in_GHz flux_density uncertainty"
    lineNumbers: which lines to read, where 1 is the first line in the file (default=all)
    columns: which columns to read for freq, flux, error (starting at 1) (default=[1,2,3])
    yscale: factor by which to scale the y values after reading them from the file
    plotunits: 'linear' or 'log' (only used if trials==0)
 
==Examples==
==Examples==
<source lang="Python">
# In CASA
CASA <6>: au.spectralindex('13A-287.sb24170282.eb24176850.56505.0217034838.ms.listobs','fluxdensities.log',labelspw=T)
</source>
<pre>
Skipping spw 0 due to INSUFFICIENT DATA result from fluxscale.
Skipping spw 1 due to INSUFFICIENT DATA result from fluxscale.
Read 64 x values for J1215-1731 (avg=44700.000)
Read 64 y values for J1215-1731 (avg=0.998)
Error-weighted fit: Slope: -0.564+-0.086  Flux D. @ 47.240GHz: 0.966+-0.002
  Un-weighted fit: Slope: -0.584        Flux D. @ 47.240GHz: 0.965
Plot saved in /tmp/fluxdensities.log.J1215-1731.png
</pre>
[[File:fluxdensities.log.J1215-1731.png]]
==File formats==
The listobs file should be a file generated by the listfile option of {{listobs}} (or at least the excerpt containing the spw information). It is keying off the referenceFrame string (default = "TOPO"), so there can be any number of leading or trailing comment characters on each line, and any number of unrelated lines of text at the top or bottom.  It also works if the spw name is present, or not present (casa < 4.1).  Here is an example of the old-style listobs format (which still works):
<pre>
SpwID  #Chans Frame Ch1(MHz)    ChanWid(kHz)  TotBW(kHz)  Corrs
  0        128 TOPO  39926      1000          128000      RR  LL
  1        128 TOPO  40054      1000          128000      RR  LL
  2        128 TOPO  40182      1000          128000      RR  LL
  3        128 TOPO  40310      1000          128000      RR  LL
  4        128 TOPO  40438      1000          128000      RR  LL
  5        128 TOPO  40566      1000          128000      RR  LL
  6        128 TOPO  40694      1000          128000      RR  LL
  7        128 TOPO  40822      1000          128000      RR  LL
  8        128 TOPO  42248      1000          128000      RR  LL
  9        128 TOPO  42376      1000          128000      RR  LL
  10        128 TOPO  42504      1000          128000      RR  LL
  11        128 TOPO  42632      1000          128000      RR  LL
  12        128 TOPO  42760      1000          128000      RR  LL
  13        128 TOPO  42888      1000          128000      RR  LL
  14        128 TOPO  43016      1000          128000      RR  LL
  15        128 TOPO  43144      1000          128000      RR  LL
</pre>
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.
<pre>
# Found reference field(s): J1331+3030
# Found transfer field(s):  J1717-3342 J1924-2914
# Flux density for J1717-3342 in SpW=0 is: 1.81919 +/- 0.0271038 (SNR = 67.1192, N= 40)
# Flux density for J1717-3342 in SpW=1 is: 1.79095 +/- 0.0272306 (SNR = 65.7698, N= 40)
# Flux density for J1717-3342 in SpW=2 is: 1.7854 +/- 0.0267479 (SNR = 66.7493, N= 40)
# Flux density for J1717-3342 in SpW=3 is: 1.77486 +/- 0.0272127 (SNR = 65.222, N= 40)
# Flux density for J1717-3342 in SpW=4 is: 1.77261 +/- 0.0278025 (SNR = 63.7573, N= 40)
# Flux density for J1717-3342 in SpW=5 is: 1.77175 +/- 0.0280525 (SNR = 63.1585, N= 40)
# Flux density for J1717-3342 in SpW=6 is: 1.76183 +/- 0.0287407 (SNR = 61.3007, N= 40)
# Flux density for J1717-3342 in SpW=7 is: 1.76226 +/- 0.0288519 (SNR = 61.0795, N= 40)
# Flux density for J1717-3342 in SpW=8 is: 1.7605 +/- 0.0303092 (SNR = 58.0849, N= 40)
# Flux density for J1717-3342 in SpW=9 is: 1.76177 +/- 0.0301251 (SNR = 58.4819, N= 40)
# Flux density for J1717-3342 in SpW=10 is: 1.75611 +/- 0.0296369 (SNR = 59.2543, N= 40)
# Flux density for J1717-3342 in SpW=11 is: 1.74256 +/- 0.0300252 (SNR = 58.0368, N= 40)
# Flux density for J1717-3342 in SpW=12 is: 1.7471 +/- 0.0304127 (SNR = 57.4463, N= 40)
# Flux density for J1717-3342 in SpW=13 is: 1.75201 +/- 0.0304666 (SNR = 57.5059, N= 40)
# Flux density for J1717-3342 in SpW=14 is: 1.74798 +/- 0.0303736 (SNR = 57.5495, N= 40)
# Flux density for J1717-3342 in SpW=15 is: 1.74099 +/- 0.0306128 (SNR = 56.8714, N= 40)
# Flux density for J1924-2914 in SpW=0 is: 13.5637 +/- 0.0392228 (SNR = 345.812, N= 40)
# Flux density for J1924-2914 in SpW=1 is: 13.393 +/- 0.0400481 (SNR = 334.424, N= 40)
# Flux density for J1924-2914 in SpW=2 is: 13.3472 +/- 0.0388534 (SNR = 343.528, N= 40)
# Flux density for J1924-2914 in SpW=3 is: 13.3187 +/- 0.0386627 (SNR = 344.483, N= 40)
# Flux density for J1924-2914 in SpW=4 is: 13.3088 +/- 0.0386251 (SNR = 344.564, N= 40)
# Flux density for J1924-2914 in SpW=5 is: 13.2596 +/- 0.0392739 (SNR = 337.62, N= 40)
# Flux density for J1924-2914 in SpW=6 is: 13.2172 +/- 0.0413275 (SNR = 319.816, N= 40)
# Flux density for J1924-2914 in SpW=7 is: 13.2503 +/- 0.0422232 (SNR = 313.815, N= 40)
# Flux density for J1924-2914 in SpW=8 is: 13.1581 +/- 0.0434202 (SNR = 303.042, N= 40)
# Flux density for J1924-2914 in SpW=9 is: 13.1242 +/- 0.0443572 (SNR = 295.875, N= 40)
# Flux density for J1924-2914 in SpW=10 is: 13.0839 +/- 0.0432023 (SNR = 302.853, N= 40)
# Flux density for J1924-2914 in SpW=11 is: 13.0526 +/- 0.0424958 (SNR = 307.15, N= 40)
# Flux density for J1924-2914 in SpW=12 is: 13.0662 +/- 0.0436696 (SNR = 299.206, N= 40)
# Flux density for J1924-2914 in SpW=13 is: 13.057 +/- 0.0440618 (SNR = 296.335, N= 40)
# Flux density for J1924-2914 in SpW=14 is: 13.0354 +/- 0.0437052 (SNR = 298.257, N= 40)
# Flux density for J1924-2914 in SpW=15 is: 13.0241 +/- 0.0444637 (SNR = 292.915, N= 40)
</pre>

Latest revision as of 12:56, 8 August 2013

Return to Analysis Utilities

This page documents the spectralindex function of the linfit class in the Python module analysisUtils.

This function is designed to fit a power-law spectral index to the results output from casa's fluxscale task (using scipy.optimize.leastsq). Currently, it also requires a listobs-generated ascii file as the first argument in order to determine the center frequencies of each spectral window. (To avoid the need for this file in the future, we have submitted an improvement request for fluxscale to also print these frequencies.) By default, it estimates the uncertainties of the fitted slope and amplitude from the square root of the diagonals of the covariance matrix. Alternatively, there is an option to run a Monte-Carlo series of fits to empirically determine the uncertainties. Finally, it produces a plot summarizing the results. By default, it will use all spectral windows, but it will exclude those with an "INSUFFICIENT DATA" status from fluxscale (which can arise in high-frequency EVLA data if you have not split off the pointing spws).

Usage:

au.spectralindex(filename, yfilename="", source="", verbose=False, maxpoints=0, trials=0, spw="", plotdir="", labelspw=False, referenceFrame='TOPO', plaintext=False, lineNumbers=None, columns=None, yscale=1.0, plotunits='linear')

   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 starting from the first one (0=no max)
   trials: if > 0, use a Monte-Carlo technique estimate the fit uncertainties,
      otherwise, use the sqrt(covarMatrix) 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 (in case you don't have permission in the pwd)
   labelspw: draw spw numeric labels adjacent to each point
   referenceFrame: the frequency reference frame in the listobs file, if not 'TOPO'
   plaintext: if True, then read the freqs, flux densities and uncertainties from a single plain text file, 
      specified by filename, with one measurement per line with a format of 
      "frequency_in_GHz flux_density uncertainty"
   lineNumbers: which lines to read, where 1 is the first line in the file (default=all)
   columns: which columns to read for freq, flux, error (starting at 1) (default=[1,2,3])
   yscale: factor by which to scale the y values after reading them from the file
   plotunits: 'linear' or 'log' (only used if trials==0)

Examples

# In CASA
CASA <6>: au.spectralindex('13A-287.sb24170282.eb24176850.56505.0217034838.ms.listobs','fluxdensities.log',labelspw=T)
Skipping spw 0 due to INSUFFICIENT DATA result from fluxscale.
Skipping spw 1 due to INSUFFICIENT DATA result from fluxscale.
Read 64 x values for J1215-1731 (avg=44700.000)
Read 64 y values for J1215-1731 (avg=0.998)
Error-weighted fit: Slope: -0.564+-0.086  Flux D. @ 47.240GHz: 0.966+-0.002
   Un-weighted fit: Slope: -0.584         Flux D. @ 47.240GHz: 0.965
Plot saved in /tmp/fluxdensities.log.J1215-1731.png

Fluxdensities.log.J1215-1731.png

File formats

The listobs file should be a file generated by the listfile option of listobs (or at least the excerpt containing the spw information). It is keying off the referenceFrame string (default = "TOPO"), so there can be any number of leading or trailing comment characters on each line, and any number of unrelated lines of text at the top or bottom. It also works if the spw name is present, or not present (casa < 4.1). Here is an example of the old-style listobs format (which still works):

 SpwID  #Chans Frame Ch1(MHz)    ChanWid(kHz)  TotBW(kHz)  Corrs
  0         128 TOPO  39926       1000          128000      RR  LL
  1         128 TOPO  40054       1000          128000      RR  LL
  2         128 TOPO  40182       1000          128000      RR  LL
  3         128 TOPO  40310       1000          128000      RR  LL
  4         128 TOPO  40438       1000          128000      RR  LL
  5         128 TOPO  40566       1000          128000      RR  LL
  6         128 TOPO  40694       1000          128000      RR  LL
  7         128 TOPO  40822       1000          128000      RR  LL
  8         128 TOPO  42248       1000          128000      RR  LL
  9         128 TOPO  42376       1000          128000      RR  LL
  10        128 TOPO  42504       1000          128000      RR  LL
  11        128 TOPO  42632       1000          128000      RR  LL
  12        128 TOPO  42760       1000          128000      RR  LL
  13        128 TOPO  42888       1000          128000      RR  LL
  14        128 TOPO  43016       1000          128000      RR  LL
  15        128 TOPO  43144       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.

# Found reference field(s): J1331+3030
# Found transfer field(s):  J1717-3342 J1924-2914
# Flux density for J1717-3342 in SpW=0 is: 1.81919 +/- 0.0271038 (SNR = 67.1192, N= 40)
# Flux density for J1717-3342 in SpW=1 is: 1.79095 +/- 0.0272306 (SNR = 65.7698, N= 40)
# Flux density for J1717-3342 in SpW=2 is: 1.7854 +/- 0.0267479 (SNR = 66.7493, N= 40)
# Flux density for J1717-3342 in SpW=3 is: 1.77486 +/- 0.0272127 (SNR = 65.222, N= 40)
# Flux density for J1717-3342 in SpW=4 is: 1.77261 +/- 0.0278025 (SNR = 63.7573, N= 40)
# Flux density for J1717-3342 in SpW=5 is: 1.77175 +/- 0.0280525 (SNR = 63.1585, N= 40)
# Flux density for J1717-3342 in SpW=6 is: 1.76183 +/- 0.0287407 (SNR = 61.3007, N= 40)
# Flux density for J1717-3342 in SpW=7 is: 1.76226 +/- 0.0288519 (SNR = 61.0795, N= 40)
# Flux density for J1717-3342 in SpW=8 is: 1.7605 +/- 0.0303092 (SNR = 58.0849, N= 40)
# Flux density for J1717-3342 in SpW=9 is: 1.76177 +/- 0.0301251 (SNR = 58.4819, N= 40)
# Flux density for J1717-3342 in SpW=10 is: 1.75611 +/- 0.0296369 (SNR = 59.2543, N= 40)
# Flux density for J1717-3342 in SpW=11 is: 1.74256 +/- 0.0300252 (SNR = 58.0368, N= 40)
# Flux density for J1717-3342 in SpW=12 is: 1.7471 +/- 0.0304127 (SNR = 57.4463, N= 40)
# Flux density for J1717-3342 in SpW=13 is: 1.75201 +/- 0.0304666 (SNR = 57.5059, N= 40)
# Flux density for J1717-3342 in SpW=14 is: 1.74798 +/- 0.0303736 (SNR = 57.5495, N= 40)
# Flux density for J1717-3342 in SpW=15 is: 1.74099 +/- 0.0306128 (SNR = 56.8714, N= 40)
# Flux density for J1924-2914 in SpW=0 is: 13.5637 +/- 0.0392228 (SNR = 345.812, N= 40)
# Flux density for J1924-2914 in SpW=1 is: 13.393 +/- 0.0400481 (SNR = 334.424, N= 40)
# Flux density for J1924-2914 in SpW=2 is: 13.3472 +/- 0.0388534 (SNR = 343.528, N= 40)
# Flux density for J1924-2914 in SpW=3 is: 13.3187 +/- 0.0386627 (SNR = 344.483, N= 40)
# Flux density for J1924-2914 in SpW=4 is: 13.3088 +/- 0.0386251 (SNR = 344.564, N= 40)
# Flux density for J1924-2914 in SpW=5 is: 13.2596 +/- 0.0392739 (SNR = 337.62, N= 40)
# Flux density for J1924-2914 in SpW=6 is: 13.2172 +/- 0.0413275 (SNR = 319.816, N= 40)
# Flux density for J1924-2914 in SpW=7 is: 13.2503 +/- 0.0422232 (SNR = 313.815, N= 40)
# Flux density for J1924-2914 in SpW=8 is: 13.1581 +/- 0.0434202 (SNR = 303.042, N= 40)
# Flux density for J1924-2914 in SpW=9 is: 13.1242 +/- 0.0443572 (SNR = 295.875, N= 40)
# Flux density for J1924-2914 in SpW=10 is: 13.0839 +/- 0.0432023 (SNR = 302.853, N= 40)
# Flux density for J1924-2914 in SpW=11 is: 13.0526 +/- 0.0424958 (SNR = 307.15, N= 40)
# Flux density for J1924-2914 in SpW=12 is: 13.0662 +/- 0.0436696 (SNR = 299.206, N= 40)
# Flux density for J1924-2914 in SpW=13 is: 13.057 +/- 0.0440618 (SNR = 296.335, N= 40)
# Flux density for J1924-2914 in SpW=14 is: 13.0354 +/- 0.0437052 (SNR = 298.257, N= 40)
# Flux density for J1924-2914 in SpW=15 is: 13.0241 +/- 0.0444637 (SNR = 292.915, N= 40)