simobserve.py

Go to the documentation of this file.

#
# This file was generated using xslt from its XML file
#
# Copyright 2009, Associated Universities Inc., Washington DC
#
import sys
import os
from  casac import *
import string
from taskinit import casalog
#from taskmanager import tm
import task_simobserve
def simobserve(project='sim', skymodel='', inbright='', indirection='', incell='', incenter='', inwidth='', complist='', compwidth='"8GHz"', setpointings=True, ptgfile='$project.ptg.txt', integration='10s', direction=[''], mapsize=['', ''], maptype='hexagonal', pointingspacing='', caldirection='', calflux='1Jy', obsmode='int', refdate='2014/01/01', hourangle='transit', totaltime='7200s', antennalist='', sdantlist='', sdant=0, thermalnoise='', user_pwv=1., t_ground=270., t_sky=260., tau0=0.1, seed=11111, leakage=0.0, graphics='both', verbose=False, overwrite=True):

        """mosaic simulation task
mosaic simulation task:

    This task simulates interferometric observations (currently 
    only ALMA can be done easily).  New functionality is actively 
    being added, so if you have changed versions of CASA, check 
    the inputs carefully.
    More information and examples are availible at 
     http://casaguides.nrao.edu/index.php?title=Simulating_Observations_in_CASA
    Please contact CASA experts with any questions, especially 
    about features noted below as *experimental*    
    -------------------------------
    project -- root filename for all output files.
    -------------------------------
    skymodel -- image model of the sky in FITS or CASA image, 
         with optional parameters that can be overridden
    inbright -- peak brightness in Jy/pixel, or "" for unchanged
       * NOTE: "unchanged" will take the numerical values in your image 
         and assume they are in Jy/pixel, even if it says some other unit 
         in the header. 
    indirection -- central direction, or "" for unchanged
    incell -- spatial pixel size, or "" for unchanged
    incenter -- frequency of center channel e.g. "89GHz", or "" for unchanged
    inwidth -- width of channels, or "" for unchanged - this should be a 
         string representing a quantity with units e.g. "10MHz"
       * NOTE: only works reliably with frequencies, not velocities
       * NOTE: it is not possible to change the number of spectral planes
         of the sky model, only to relabel them with different frequencies
         That kind of regridding can be accomplished with the CASA toolkit.
    -------------------------------
    complist -- component list model of the sky, added to or instead of skymodel
    compwidth -- bandwidth of components; if simulating from components only, 
         this defines the bandwidth of the MS and output images
    -------------------------------
    setpointings -- calculate a map of pointings, or if false, provide ptgfile
       * if graphics are on, display the pointings shown on the model image
    ptgfile -- a text file specifying directions in the following
         format, with optional integration times, e.g.
         #Epoch     RA          DEC      TIME(optional)
         J2000 23h59m28.10 -019d52m12.35 10.0
       * if the time column is not present in the file, it will use
         "integration" for all pointings.
       * NOTE: at this time the file should contain only science pointings:
         simobserve will observe these, then optionally the calibrator, 
         then the list of science pointings again, etc, until totaltime
         is used up. 
    integration --- Time interval for each integration e.g '10s'
       * NOTE: to simulate a "scan" longer than one integration, use 
         setpointings to generate a pointing file, and then edit the 
         file to increase the time at each point to be larger than 
         the parameter integration time.
    direction -- mosaic center direction e.g 'J2000 19h00m00 -40d00m00'
       * can optionally be a list of pointings
       * otherwise simobserve will pack mapsize according to maptype
    mapsize -- angular size of map 
       * set to "" to span the model image
    maptype -- hexagonal, square (rectangular raster), 
         "ALMA" for the same hex algorithm as the ALMA OT 
         or "ALMA2012" for the algorithm used in the Cycle 0 OT
    pointingspacing -- spacing in between beams e.g '1arcsec'
         or "0.25PB" to use 1/4 of the primary beam FWHM, or 
         "" will use Nyquist spacing, lambda/d/sqrt(3).
    -------------------------------
    obsmode -- observation mode to calculate visibilities from skymodel (which
         may have been modified above, (optionally) complist, and $ptgfile (which 
         may have been generated above)
       * this parameter takes two possible values:
       - interferometer (or i)
       - singledish (or s)
       * if graphics are on, display the array (like plotants), the uv 
         coverage, the synthesized (dirty) beam, and ephemeris information
       * if simulating from component list, you should specify
         compwidth, the desired bandwidth - 
         there is not currently a way to specify the spectrum of a component, 
         so simulations from a componentlist only will be continuum (1 chan)
    refdate -- date of simulated observation eg: '2014/05/21'
    hourangle -- hour angle of observation e.g. '-3h'
    totaltime --- total time of observation e.g '7200s' or if a number without
         units, interpreted as the number of times to repeat the map
    antennalist -- ascii file containing antenna positions.
         each row has x y z coordinates and antenna diameter; 
         header lines are required to specify the observatory name
         and coordinate system e.g. 
           # observatory=ALMA
           # coordsys=UTM
           # datum=WGS84
           # zone=19
       * standard arrays are found in your CASA data repository, 
         os.getenv("CASAPATH").split()[0]+"/data/alma/simmos/"        
       * if "", will not not produce an interferometric MS
       * a string of the form "alma;0.5arcsec" will be parsed into a full 12m ALMA
         configuration.  This only works for full ALMA and may fail to find the 
         standard configuration files on some systems - see casaguides.nrao.edu
    caldirection -- *NEW* an unresolved calibrator can be observed 
         interleaved with the science pointings.  This feature is 
         experimental, so please contact us with any questions.
       * The calibrator is implemented as a point source clean component 
         with this direction and flux=calflux
    sdant -- the index of the antenna in the list to use for total 
         power.  defaults to the first antenna on the list. 
    -------------------------------
    thermalnoise -- add thermal noise 
       * this parameter takes two possible values:
       - tsys-atm: J. Pardo's ATM library will be used to construct an
           atmospheric profile for the ALMA site:
           altitude 5000m, ground pressure 650mbar, relhum=20%, 
           a water layer of user_pwv at altitude of 2km,
           the sky brightness temperature returned by ATM, 
           and internally tabulated receiver temperatures
       - tsys-manual: instead of using the ATM model, specify the zenith 
           sky brightness and opacity manually.  Noise is added and then
           the visibility flux scale is referenced above the atmosphere. 
       * In either mode, noise is calculated using an antenna spillover 
           efficiency of 0.96, taper of 0.86, 
           surface accuracy of 25 and 300 microns for ALMA and EVLA
           respectively (using the Ruze formula for surface efficiency),
           correlator efficiencies of 0.95 and 0.91 for ALMA and EVLA, 
           receiver temperatures for ALMA of 
               17, 30, 37, 51, 65, 83,147,196,175,230 K interpolated between
               35, 75,110,145,185,230,345,409,675,867 GHz,
           for EVLA of 
               500, 70,  60,  55,  100, 130, 350 K interpolated between
               0.33,1.47,4.89,8.44,22.5,33.5,43.3 GHz,
           for SMA of
               67,  116, 134, 500 K interpolated between
               212.,310.,383.,660. GHz
       * These are only approximate numbers and do not take into account 
         performance at edges of receiver bands, neither are they guaranteed
         to reflect the most recent measurements.  Caveat emptor and use the 
         sm tool to add noise if you want more precise control.
    t_ground -- ground/spillover temperature in K 
    user_pwv -- precipitable water vapor if constructing an atmospheric model
    t_sky -- atmospheric temperature in K [for tsys-manual]
    tau0 -- zenith opacity at observing frequency [for tsys-manual]
       * see casaguides.nrao.edu for more information on noise, 
         in particular how to add a phase screen using the toolkit
    seed -- random number seed for noise generation
    -------------------------------
    leakage -- add cross polarization corruption of this fractional magnitude
    
    graphics -- view plots on the screen, saved to file, both, or neither
    verbose -- print extra information to the logger and terminal
    overwrite -- overwrite existing files in the project subdirectory



    -------------------------------
    How to specify a model image:
    -------------------------------
    * simobserve requires a CASA or fits image. If you merely have a grid of 
      numbers, you will need to write them out as fits or write a CASA script to
      read them in and use the ia tool to create an image and insert the data.
    
    * simobserve does NOT require a coordinate system in the header. If the
      coordinate information is incomplete, missing, or you would like to
      override it, set the appropriate "in" paremeters.  
      In the absence of other information, simobserve will assume that the axes of 
      your input correspond to RA, Dec, and (optionally) frequency and 
      (optionally) Stokes parameter.  

    * If you have a proper Coordinate System, simobserve will do its best to
      generate visibilities from that, and then create a synthesis image
      according to the specified user parameters.  Regridding the
      spectral dimension may not have complete flexibility yet.

    * You can manipulate an image header with the "imhead" task, or you can
      delve deeper with the ia and cs tools.  If you use the tools, you should
      be aware that a CoordinateSystem in CASA can exist independently of an
      Image.  Once the CoordinateSystem is detached from the image, it is the
      user's responsibility to do any manipulation e.g. axis reordering on
      both. Example:
        ia.open("myimage_filename")
        ia.summary() # see header as attached to the image
        csys=ia.coordsys() # detach the CoordinateSystem
        csys.summary() # examine it
        csys.setreferencepixel([100,100]) 
        arr=ia.getchunk() # get the data from the Image
        ia.done()
        csys.reorder([0,2,1])  # reorder the CoordinateSystem
        arr=arr.reorder([0,2,1])  # reorder the data         
        ia.fromshape(outfile="mynewimage_file",shape=[32,32,256],csys=csys.torecord(),overwrite=True)
            # make a new image, with the right shape and CoordinateSystem
        ia.putchunk(arr)   # put the data into the new image.
        csys.done()
        ia.done()

    -------------------------------
    Output produced: (not all will always exist, depending on input parameters)
    To support different runs with different arrays, the names have the
    configuration name from antennalist appended.
    -------------------------------

    project.[cfg].skymodel = 4d input sky model image (optionally) scaled
    project.[cfg].skymodel.flat.regrid.conv = input sky regridded to match the
        output image, and convolved with the output clean beam
    project.[cfg].skymodel.png = diagnostic figure of sky model with pointings

    project.[cfg].ptg.txt = list of mosaic pointings
    project.[cfg].quick.psf = psf calculated from uv coverage
    project.[cfg].ms = noise-free measurement set
    project.[cfg].noisy.ms = corrupted measurement set 
    project.[cfg].observe.png = diagnostic figure of uv coverage and 
         visibilities
            
    project.[cfg].simobserve.last = saved input parameters for simobserve task


        """
        if type(direction)==str: direction=[direction]
        if type(mapsize)==str: mapsize=[mapsize]

#
#    The following is work around to avoid a bug with current python translation
#
        mytmp = {}

        mytmp['project'] = project
        mytmp['skymodel'] = skymodel
        mytmp['inbright'] = inbright
        mytmp['indirection'] = indirection
        mytmp['incell'] = incell
        mytmp['incenter'] = incenter
        mytmp['inwidth'] = inwidth
        mytmp['complist'] = complist
        mytmp['compwidth'] = compwidth
        mytmp['setpointings'] = setpointings
        mytmp['ptgfile'] = ptgfile
        mytmp['integration'] = integration
        mytmp['direction'] = direction
        mytmp['mapsize'] = mapsize
        mytmp['maptype'] = maptype
        mytmp['pointingspacing'] = pointingspacing
        mytmp['caldirection'] = caldirection
        mytmp['calflux'] = calflux
        mytmp['obsmode'] = obsmode
        mytmp['refdate'] = refdate
        mytmp['hourangle'] = hourangle
        mytmp['totaltime'] = totaltime
        mytmp['antennalist'] = antennalist
        mytmp['sdantlist'] = sdantlist
        mytmp['sdant'] = sdant
        mytmp['thermalnoise'] = thermalnoise
        mytmp['user_pwv'] = user_pwv
        mytmp['t_ground'] = t_ground
        mytmp['t_sky'] = t_sky
        mytmp['tau0'] = tau0
        mytmp['seed'] = seed
        mytmp['leakage'] = leakage
        mytmp['graphics'] = graphics
        mytmp['verbose'] = verbose
        mytmp['overwrite'] = overwrite
        pathname='file:///'+os.environ.get('CASAPATH').split()[0]+'/share/xml/'
        trec = casac.utils().torecord(pathname+'simobserve.xml')

        casalog.origin('simobserve')
        if trec.has_key('simobserve') and casac.utils().verify(mytmp, trec['simobserve']) :
            result = task_simobserve.simobserve(project, skymodel, inbright, indirection, incell, incenter, inwidth, complist, compwidth, setpointings, ptgfile, integration, direction, mapsize, maptype, pointingspacing, caldirection, calflux, obsmode, refdate, hourangle, totaltime, antennalist, sdantlist, sdant, thermalnoise, user_pwv, t_ground, t_sky, tau0, seed, leakage, graphics, verbose, overwrite)

        else :
          result = False
        return result