Difference between revisions of "Protoplanetary Disk Simulation (CASA 3.4)"

From CASA Guides
Line 2: Line 2:
[[Category: Simulations]]
[[Category: Simulations]]
* '''This is an advanced simulation tutorial.  New users are recommended to begin with the [[Simulation Guide for New Users (CASA 3.4)]].'''
* '''This guide is applicable to CASA version 3.4.  For older versions of CASA see [[PPdisk_simdata_(CASA_3.3)]].'''
* '''This guide is applicable to CASA version 3.4.  For older versions of CASA see [[PPdisk_simdata_(CASA_3.3)]].'''
* '''To create a script of the Python code on this page see [[Extracting scripts from these tutorials]].'''
* '''To create a script of the Python code on this page see [[Extracting scripts from these tutorials]].'''

Revision as of 16:03, 11 June 2012

Simulating Observations in CASA


For this CASA Guide we will use a protoplanetary disk model from S. Wolf. If you use this FITS data for anything more than learning CASA, please cite Wolf & D'Angelo 2005.

Script with Explanation

Set simobserve as the current task and reset all parameters.


Review the image coordinate system using task imhead.

# This reports image header parameters in the Log Messages window

We now us the ia (image analysis) and qa (units and quantities) tools from the CASA Toolkit to find the image center. In comparison to tasks, tools are a more advanced way of manipulating data in CASA. You can learn more about tools using the tool reference manual.

When data are being manipulated with tools the data file must be explicitly opened and closed.


Next, get the right ascension and declination of the image center. We get the number of pixels along each axis using ia.shape. Then, we get the RA and Dec values for the center pixel using ia.toworld.

axesLength = ia.shape()
# Divide the first two elements of axesLength by 2.
center_pixel = [ x / 2.0 for x in axesLength[:2] ]
# Feed center_pixel to ia.toworld and and save the RA and Dec to ra_radians and dec_radians
(ra_radians, dec_radians) = ia.toworld( center_pixel )['numeric'][:2]

Use the qa tool to convert the image center from radians to sexagesimal coordinates.

ra_hms  = qa.formxxx(str(ra_radians)+"rad",format='hms',prec=5)
dec_dms = qa.formxxx(str(dec_radians)+"rad",format='dms',prec=5)
Brightness scale can be viewed with 'imstat' task
# Default parameters are adequate for this
# ...
#  'max': array([  6.52469971e-05]),
# ...
# that's 0.0652 mJy/pixel.

Let's call our project "psim2". This defines the root prefix for any output files from simobserve.

project = "psim2"

We'll set skymodel to the FITS file downloaded, above, and leave all skymodel subparameters at their default values. simobserve will create CASA image psim2.skymodel.

skymodel = "input50pc_672GHz.fits"

We will specify the sky position for the center of the observation and set the map size to the size of the model image. Since the model image is 2/3 arcseconds across, we should only need one pointing. In this case, pointingspacing and maptype can be left at their default values.

setpointings       =  True
direction          =  "J2000 18h00m00.031s -22d59m59.6s"
mapsize            =  "0.76arcsec"

We do want to simulate an interferometric observation, so we set obsmode accordingly. We'll set totaltime to a 20-minute snapshot observation.

obsmode            =  "int"
totaltime          =  "1200s"

We want to use the appropriate antenna configuration for the desired angular resolution. Configuration 20, alma.out20.cfg, is the largest "compact" configuration.

antennalist        =  "alma.out20.cfg"

Now run simobserve.

Deconvolve the visibilities back into an image
default ("simanalyze")
project            =  "psim2"
image              =  True
# Prior image to use in clean
modelimage         =  "input50pc_672GHz.fits"
vis                =  project+".alma.out20.ms"
imsize             =  [192, 192]
Specify number of iteration of cleaning task with proper threshold and weighting
niter              =  10000
threshold          =  "1e-7Jy"
weighting          =  "natural"
We'd like to calculate a difference and fidelity image, and see some diagnostics:
analyze            =  True
And see the array but not the UV coverage:
showuv             =  False
showresidual       =  True  
showconvolved      =  True
Plot both to the screen and the png files with lots of messages:
graphics           =  "both"
verbose            =  True
overwrite = True

Run simdata

# This commands CASA to execute simanalyze

  • Output results:


Simulating Observations in CASA

Last checked on CASA Version 3.4.0.