View source for CASA Contributed Script: vlacal-setjy

<div style="background-color: salmon">
<div style="background-color: salmon; margin: 20px">
<br>
This page is currently under development.
Please do not use this guide until this message has been removed.
<br><br>
</div>
</div>

__TOC__


== Introduction ==
Flux-density calibration for the VLA typically involves an observation of one of four primary calibrators: 3C286, 3C48, 3C138 and 3C147. Models for these sources are distributed as part of the CASA package and model visibilities are predicted using the task [https://casadocs.readthedocs.io/en/v6.4.0/api/tt/casatasks.imaging.setjy.html setjy]. The most current models available in 
[https://casadocs.readthedocs.io/en/v6.4.0/api/tt/casatasks.imaging.setjy.html setjy] are accessible using the flux-density standard [https://iopscience.iop.org/article/10.3847/1538-4365/aa6df9 Perley Butler 2017] (hereafter PB17), which is based on observations conducted on or before January 2016. [https://science.nrao.edu/facilities/vla/docs/manuals/oss/performance/fdscale See this link] for further discussion of VLA flux-density calibration.

Given the recent and significant variability in some primary calibrators (i.e., 3C48 and 3C138) NRAO staff are in the process of creating a new flux-density standard based on more recent observations.  Instead of packaging these new models with CASA, NRAO will host a database that is continually updated with the latest monitoring data (having an approximately monthly cadence). A new flux-density standard will be added to the task [https://casadocs.readthedocs.io/en/v6.4.0/api/tt/casatasks.imaging.setjy.html setjy]  to make use of these new models (estimated release date 2023). Until the new version of [https://casadocs.readthedocs.io/en/v6.4.0/api/tt/casatasks.imaging.setjy.html setjy] is released, this contributed script [http://www.aoc.nrao.edu/~jmarvil/vlacal-setjy.py vlacal-setjy.py] can be executed to obtain the newest models from the database.  This casaguide provides documentation and examples of how to use this script. Please ensure you are using a recent version of CASA to run the script (this guide was tested using CASA v6.4.0, and the script is known not to work in CASA versions <= 5.8).


== Obtaining the script ==

The script [http://www.aoc.nrao.edu/~jmarvil/vlacal-setjy.py vlacal-setjy.py is available for download here].  Please save a copy in your local working directory, e.g.,

<source lang="python">
# in CASA
!wget www.aoc.nrao.edu/~jmarvil/vlacal/vlacal-setjy.py
</source>

== Loading the functions ==

Next, you will want to run the script in CASA, i.e., 

<source lang="python">
# in CASA
execfile('vlacal-setjy.py')
</source>

This will import functions that you can now run from your CASA session.  You will need to run the script each time you start CASA if you want these functions to be available (or add the above to [https://casadocs.readthedocs.io/en/v6.4.0/api/configuration.html?highlight=startup.py#startup-py startup.py]). Note that there is no need to open the script in an editor or make any modifications to it.

Below, this guide describes in detail two functions and how to use them.  These functions represent two different approaches depending on the state of your data.  If you already calibrated and/or imaged your data then the first function may offer a more simple solution to correcting the flux-density scale.  Alternatively, if you are starting with uncalibrated visibilities then the second function is the recommended approach. You only need to follow one of these methods to obtain accurate results.


== Calculating a flux-density ratio ==

One function provided by this script will calculate the ratio of the PB17 flux to the new NRAO flux standard. This ratio may be of interest to users who have already calibrated their data using the PB17 flux standard.  In this case, the ratio could be used to scale the visibility amplitude or the image intensity. The function to calculate this ratio is called ''vlacal_setjy_calculate_ratio''. The following functional arguments will need to be provided:
* '''source''': The string name of the VLA primary calibrator, i.e., one of '3C286', '3C48', '3C138' or '3C147'
* '''band''': The single-letter string name of the VLA receiver band, i.e., one of 'P', 'L', 'S', 'C', 'X', 'U', 'K', 'A' or 'Q'
* '''date''': The date to determine the flux-density scale (typically the date your observation took place), formatted as 'YYYY/MM/DD'

For example, for source 3C48 at Ka-band on June 1, 2019 you would run:
<source lang="python">
# in CASA
vlacal_setjy_calculate_ratio( source='3C48', band='A', date='2019/06/01' )
</source>

which will produce the following output:
<pre style="background-color:lightgrey;">
Obtained flux of 1.074 Jy on mjd 58635.0 for source 3C48 at reference frequency 33.25 GHz
This interpolated flux differs from Perley-Butler 2017 by 32.65%
</pre>

Based on this output, one could scale their results by a factor of 1.3265 to account for the estimated flux of 3C48 at their time of observation. Images can be scaled with the task [https://casadocs.readthedocs.io/en/latest/api/tt/casatasks.analysis.immath.html immath].  For example, to apply the above result of 32.65% to an image, run [https://casadocs.readthedocs.io/en/latest/api/tt/casatasks.analysis.immath.html immath] using ''mode='evalexpr' '' and ''expr='1.3265*IM0' ''.   Visibilities can be scaled by using the task [https://casadocs.readthedocs.io/en/latest/api/tt/casatasks.calibration.gencal.html gencal] to create a time-independent amplitude calibration table.  For example, to apply the above result of 32.65% to a MS, use [https://casadocs.readthedocs.io/en/latest/api/tt/casatasks.calibration.gencal.html gencal] with ''caltype='amp''' and ''parameter=[1.3265]'' to create a new calibration table, then apply this table with the task [https://casadocs.readthedocs.io/en/latest/api/tt/casatasks.calibration.applycal.html applycal].

In this example, we would recommend the following line be included in publications using these data:

<blockquote>"These data were calibrated using 3C48 assuming a flux density of 1.074 Jy at 33.25 GHz, as obtained using the NRAO VLA flux-density standard."</blockquote>


== Recalibrating your visibilities ==

=== The setjy replacement function ===
Another function provided by this script will fill the MODEL_DATA column of your measurement set.  This is intended to replace the running of the task [https://casadocs.readthedocs.io/en/v6.4.0/api/tt/casatasks.imaging.setjy.html setjy] in the VLA calibration procedure. The setjy replacement function is called ''vlacal_setjy_predict_model_data''. The following functional arguments will need to be provided:
* '''vis''': The string name of the measurement set to be processed
* '''field''': The string field ID(s) of the primary calibrator to be processed, e.g., '0'.  All selected fields need to be for the same source. 
* '''source''': The string name of the VLA primary calibrator corresponding to the field parameter, i.e., one of '3C286', '3C48', '3C138' or '3C147'
* '''spw''': The CASA spectral window selection string to be processed. All SPWs need to be from the same receiver band.  
* '''band''': The single-letter string name of the VLA receiver band corresponding to the spw parameter, i.e., one of 'P', 'L', 'S', 'C', 'X', 'U', 'K', 'A' or 'Q'
* '''date''': The date to determine the flux-density scale (typically the date your observation took place), formatted as'YYYY/MM/DD'

If your measurement set contains multiple primary flux calibrators and/or multiple receiver bands then you will need to run the above function once per source/band combination, identifying the appropriate field-source and spw-band relationships for each combination. 


=== Example function usage ===

In this example, we will use the following measurement set: [http://www.aoc.nrao.edu/~jmarvil/vlacal_example.ms.tar.gz vlacal_example.ms.tar.gz] (82 MB)

We will download and untar the example measurement set, then run the task [https://casadocs.readthedocs.io/en/v6.4.0/api/tt/casatasks.information.listobs.html listobs] to view a summary of its contents.

<source lang="python">
# in CASA
listobs( vis='vlacal_example.ms', verbose=True )
</source>


A portion of the [https://casadocs.readthedocs.io/en/v6.4.0/api/tt/casatasks.information.listobs.html listobs] output is shown below, as it appears in the log file and the logger window:

<pre style="background-color: #fffacd;">
================================================================================
           MeasurementSet Name:  vlacal_example.ms      MS Version 2
================================================================================
   Observer: Lorant Sjouwerman     Project: uid://evla/pdb/36543887  
Observation: EVLA
Computing scan and subscan properties...
Data records: 84721       Total elapsed time = 518 seconds
   Observed from   30-Nov-2019/08:53:26.0   to   30-Nov-2019/09:02:04.0 (UTC)
   
   ObservationID = 0         ArrayID = 0
  Date        Timerange (UTC)          Scan  FldId FieldName             nRows     SpwIds   Average Interval(s)    ScanIntent
  30-Nov-2019/08:53:26.0 - 08:53:50.0    27      0 0137+331=3C48            18977  [0,1,2,3,4,5,6,7]  [2, 2, 2, 2, 2, 2, 2, 2] [CALIBRATE_AMPLI#UNSPECIFIED,CALIBRATE_PHASE#UNSPECIFIED]
              08:53:54.0 - 08:54:20.0    28      0 0137+331=3C48            19160  [8,9,10,11,12,13,14,15]  [2, 2, 2, 2, 2, 2, 2, 2] [CALIBRATE_AMPLI#UNSPECIFIED,CALIBRATE_PHASE#UNSPECIFIED]
              09:01:08.0 - 09:01:34.0    38      1 J0319+4130               23096  [0,1,2,3,4,5,6,7]  [2, 2, 2, 2, 2, 2, 2, 2] [CALIBRATE_AMPLI#UNSPECIFIED,CALIBRATE_PHASE#UNSPECIFIED]
              09:01:38.0 - 09:02:04.0    39      1 J0319+4130               23488  [8,9,10,11,12,13,14,15]  [2, 2, 2, 2, 2, 2, 2, 2] [CALIBRATE_AMPLI#UNSPECIFIED,CALIBRATE_PHASE#UNSPECIFIED]
           (nRows = Total number of rows per scan) 
Fields: 2
  ID   Code Name                RA               Decl           Epoch   SrcId      nRows
  0    NONE 0137+331=3C48       01:37:41.299431 +33.09.35.13299 J2000   0          38137
  1    NONE J0319+4130          03:19:48.160102 +41.30.42.10305 J2000   1          46584
Spectral Windows:  (16 unique spectral windows and 1 unique polarization setups)
  SpwID  Name           #Chans   Frame   Ch0(MHz)  ChanWid(kHz)  TotBW(kHz) CtrFreq(MHz) BBC Num  Corrs  
  0      EVLA_Q#A0C0#68     64   TOPO   39988.000      2000.000    128000.0  40051.0000       12  RR  LL
  1      EVLA_Q#A0C0#69     64   TOPO   40116.000      2000.000    128000.0  40179.0000       12  RR  LL
  2      EVLA_Q#A0C0#70     64   TOPO   40244.000      2000.000    128000.0  40307.0000       12  RR  LL
  3      EVLA_Q#A0C0#71     64   TOPO   40372.000      2000.000    128000.0  40435.0000       12  RR  LL
  4      EVLA_Q#A0C0#72     64   TOPO   40500.000      2000.000    128000.0  40563.0000       12  RR  LL
  5      EVLA_Q#A0C0#73     64   TOPO   40628.000      2000.000    128000.0  40691.0000       12  RR  LL
  6      EVLA_Q#A0C0#74     64   TOPO   40756.000      2000.000    128000.0  40819.0000       12  RR  LL
  7      EVLA_Q#A0C0#75     64   TOPO   40884.000      2000.000    128000.0  40947.0000       12  RR  LL
  8      EVLA_K#A0C0#84     64   TOPO   20188.000      2000.000    128000.0  20251.0000       12  RR  LL
  9      EVLA_K#A0C0#85     64   TOPO   20316.000      2000.000    128000.0  20379.0000       12  RR  LL
  10     EVLA_K#A0C0#86     64   TOPO   20444.000      2000.000    128000.0  20507.0000       12  RR  LL
  11     EVLA_K#A0C0#87     64   TOPO   20572.000      2000.000    128000.0  20635.0000       12  RR  LL
  12     EVLA_K#A0C0#88     64   TOPO   20700.000      2000.000    128000.0  20763.0000       12  RR  LL
  13     EVLA_K#A0C0#89     64   TOPO   20828.000      2000.000    128000.0  20891.0000       12  RR  LL
  14     EVLA_K#A0C0#90     64   TOPO   20956.000      2000.000    128000.0  21019.0000       12  RR  LL
  15     EVLA_K#A0C0#91     64   TOPO   21084.000      2000.000    128000.0  21147.0000       12  RR  LL

</pre>

From the [https://casadocs.readthedocs.io/en/v6.4.0/api/tt/casatasks.information.listobs.html listobs] output we can see that field ID 0 is the VLA primary flux calibrator 3C48, so we will set field='0' and source='3C48'. We also see that 3C48 was observed in two different receiver bands, i.e., K and Q. This means we will need to run the function 'vlacal_setjy_predict_model_data' twice, first with band='Q' and spw='0~7', then with band='K' and spw='8~15'.  Finally, we can read the observation date from the listobs output and format it as required by the function, i.e., date='2019/11/30'.

We are now ready to run the function to predict the model visibilities for field ID 0:

<source lang="python">
# in CASA
vlacal_setjy_predict_model_data( vis='vlacal_example.ms', field='0', source='3C48', spw='0~7', band='Q', date='2019/11/30' )
vlacal_setjy_predict_model_data( vis='vlacal_example.ms', field='0', source='3C48', spw='8~15', band='K', date='2019/11/30' )

</source>

which will produce the following output:
<pre style="background-color:lightgrey;">
Interpolated 3C48 flux on mjd 58817 of 0.854 Jy at reference frequency 45.0 GHz
Filling MODEL_DATA column using task_ft for data selection: field='0', spw='0~7'

Interpolated 3C48 flux on mjd 58817 of 1.469 Jy at reference frequency 22.25 GHz
Filling MODEL_DATA column using task_ft for data selection: field='0', spw='8~15'
</pre>

The MODEL_DATA column has been filled using the most current monitoring data from the new NRAO VLA flux-density standard. The standard calibration procedure may now continue, having run this function instead of [https://casadocs.readthedocs.io/en/v6.4.0/api/tt/casatasks.imaging.setjy.html setjy]. We recommend the following line be included in publications using these data:

<blockquote>"These data were calibrated using 3C48 assuming a flux density of 0.854 Jy at 45 GHz and 1.469 Jy at 22.25 GHz, as obtained using the NRAO VLA flux-density standard."</blockquote>


Questions about this tutorial? Please contact the [https://science.nrao.edu/observing/helpdesk NRAO Helpdesk].

<div style="background-color: #CCFF99;">
Last checked on CASA Version 6.4.0.
</div>