PPdisk simdata (CASA 3.1): Difference between revisions

From CASA Guides
Jump to navigationJump to search
Line 8: Line 8:


<br>
<br>
*Explanation of the script:
====Explanation of the script====


Set simdata2 as current task and reset all parameters
=====Set simdata2 as current task and reset all parameters=====
   default("simdata2")   
   default("simdata2")   
Specify sky model image
=====Specify sky model image=====
   modelimage        =  "input50pc_672GHz.fits"   
   modelimage        =  "input50pc_672GHz.fits"   
Image coordinate system can be verified
=====Image coordinate system can be verified=====
   imhead("input50pc_672GHz.fits")   
   imhead("input50pc_672GHz.fits")   
Image center can be identified
=====Image center can be identified=====
   # ia.open("input50pc_672GHz.fits")
   # ia.open("input50pc_672GHz.fits")
   # ia.shape()
   # ia.shape()
Line 27: Line 27:
   # '-022.59.59.602743'
   # '-022.59.59.602743'
   # ia.done()
   # ia.done()
Brightness scale can be viewed with 'imstat' task
=====Brightness scale can be viewed with 'imstat' task=====
   # imstat("input50pc_672GHz.fits")
   # imstat("input50pc_672GHz.fits")
   # ...
   # ...
Line 33: Line 33:
   # ...
   # ...
   # that's 0.0652 mJy/pixel.   
   # that's 0.0652 mJy/pixel.   
Let's leave the brightness of the image as it is.
=====Let's leave the brightness of the image as it is=====
   inbright          =  "unchanged"   
   inbright          =  "unchanged"   
Let's call our project psim2
=====Let's call our project psim2=====
   project            =  "psim2"   
   project            =  "psim2"   
We'll leave the sky model the way it is -- simdata2 will create psim2.skymodel CASA image since this  
=====We'll leave the sky model the way it is -- simdata2 will create psim2.skymodel CASA image since this  
   model is a fits file, and most but not all of CASA routines can operate directly on fits
   model is a fits file, and most but not all of CASA routines can operate directly on fits=====
   modifymodel        =  False   
   modifymodel        =  False   
   skymodel          =  "input50pc_672GHz.fits"   
   skymodel          =  "input50pc_672GHz.fits"   
We need to decide where to point the telescope.  The image is 2/3 arcsec in size, so we only need  
=====We need to decide where to point the telescope.  The image is 2/3 arcsec in size, so we only need  
   one pointing.  We could put that in a text file ourself, or let simdata2 create the ascii
   one pointing.  We could put that in a text file ourself, or let simdata2 create the ascii
   pointing file for us.  
   pointing file for us.=====
   setpointings      =  True
   setpointings      =  True
   direction          =  "J2000 18h00m00.031s -22d59m59.6s"
   direction          =  "J2000 18h00m00.031s -22d59m59.6s"
   mapsize            =  "0.76arcsec"   
   mapsize            =  "0.76arcsec"   
The default pointingspacing is fine - we'll only fit one pointing in the small mapsize
=====The default pointingspacing is fine - we'll only fit one pointing in the small mapsize
the default calculation maptype hexagonal is ok too since only one will fit anyway.
the default calculation maptype hexagonal is ok too since only one will fit anyway.=====
    
    
We do want to calculate visibilities in a measurement set:
=====We do want to calculate visibilities in a measurement set:
let's do a 20 min snapshot observation using out20 configuration:
let's do a 20 min snapshot observation using out20 configuration:=====
   predict            =  True
   predict            =  True
   totaltime          =  "1200s"   
   totaltime          =  "1200s"   
Use appropriate antenna configurations based on desired angular resolution
=====Use appropriate antenna configurations based on desired angular resolution=====
   repodir=os.getenv("CASAPATH").split(' ')[0]
   repodir=os.getenv("CASAPATH").split(' ')[0]
   antennalist        =  repodir+"/data/alma/simmos/alma.out20.cfg"   
   antennalist        =  repodir+"/data/alma/simmos/alma.out20.cfg"   
Deconvolve the visibilities back into an image
=====Deconvolve the visibilities back into an image=====
   image              =  True
   image              =  True
   vis                =  "$project.ms"
   vis                =  "$project.ms"
   imsize            =  [192, 192]   
   imsize            =  [192, 192]   
Specify number of iteration of cleaning task with proper threshold and weighting
=====Specify number of iteration of cleaning task with proper threshold and weighting=====
   niter              =  10000
   niter              =  10000
   threshold          =  "1e-7Jy"
   threshold          =  "1e-7Jy"
   weighting          =  "natural"     
   weighting          =  "natural"     
We'd like to calculate a difference and fidelity image, and see some diagnostics:
=====We'd like to calculate a difference and fidelity image, and see some diagnostics:=====
   analyze            =  True   
   analyze            =  True   
And see the array but not the UV coverage:
=====And see the array but not the UV coverage:=====
   showarray          =  True
   showarray          =  True
   showuv            =  False   
   showuv            =  False   
Plot both to the screen and the png files, and giving us lots of messages:
=====Plot both to the screen and the png files, and giving us lots of messages:=====
   graphics          =  "both"
   graphics          =  "both"
   verbose            =  True
   verbose            =  True
   overwrite = True
   overwrite = True


*To run the script, type  
===To run the script, type ====


   CASA<> execfile("Ppdisk.simdata2.txt")
   CASA<> execfile("Ppdisk.simdata2.txt")

Revision as of 18:34, 8 June 2010

Simulating Observations in CASA

Protoplanetary disk


Explanation of the script

Set simdata2 as current task and reset all parameters
 default("simdata2")   
Specify sky model image
 modelimage         =  "input50pc_672GHz.fits"  
Image coordinate system can be verified
 imhead("input50pc_672GHz.fits")  
Image center can be identified
 # ia.open("input50pc_672GHz.fits")
 # ia.shape()
 # [257L, 257L, 1L, 1L]
 # ia.toworld([128.5,128.5])
 # {'numeric': array([  4.71239120e+00,  -4.01423802e-01,   1.00000000e+00,  6.72000001e+11])}
 # qa.formxxx("4.71239120rad",format='hms',prec=5)
 # '18:00:00.03052'
 # qa.formxxx("-0.401423802rad",format='dms',prec=5)
 # '-022.59.59.602743'
 # ia.done()
Brightness scale can be viewed with 'imstat' task
 # imstat("input50pc_672GHz.fits")
 # ...
 #  'max': array([  6.52469971e-05]),
 # ...
 # that's 0.0652 mJy/pixel.   
Let's leave the brightness of the image as it is
 inbright           =  "unchanged"  
Let's call our project psim2
 project            =  "psim2"  

=====We'll leave the sky model the way it is -- simdata2 will create psim2.skymodel CASA image since this

 model is a fits file, and most but not all of CASA routines can operate directly on fits=====
 modifymodel        =  False  
 skymodel           =  "input50pc_672GHz.fits"   

=====We need to decide where to point the telescope. The image is 2/3 arcsec in size, so we only need

 one pointing.  We could put that in a text file ourself, or let simdata2 create the ascii
 pointing file for us.===== 
 setpointings       =  True
 direction          =  "J2000 18h00m00.031s -22d59m59.6s"
 mapsize            =  "0.76arcsec"  

=====The default pointingspacing is fine - we'll only fit one pointing in the small mapsize the default calculation maptype hexagonal is ok too since only one will fit anyway.=====

=====We do want to calculate visibilities in a measurement set: let's do a 20 min snapshot observation using out20 configuration:=====

 predict            =  True
 totaltime          =  "1200s"  
Use appropriate antenna configurations based on desired angular resolution
 repodir=os.getenv("CASAPATH").split(' ')[0]
 antennalist        =  repodir+"/data/alma/simmos/alma.out20.cfg"  
Deconvolve the visibilities back into an image
 image              =  True
 vis                =  "$project.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:
 showarray          =  True
 showuv             =  False  
Plot both to the screen and the png files, and giving us lots of messages:
 graphics           =  "both"
 verbose            =  True
 overwrite = True

To run the script, type =

 CASA<> execfile("Ppdisk.simdata2.txt")
 CASA<> go simdata2


  • Output results:


Input:
Psim2.skymodel.png
Predict:
Psim2.predict.png
Image:
Psim2.image.png
Analyze:
Psim2.analysis.png