ACA Simulation CASA 6.5.4: Difference between revisions
| (14 intermediate revisions by the same user not shown) | |||
| Line 6: | Line 6: | ||
We can simulate ALMA observations that use the main array of 12 m antennas augmented by the ACA and Total Power antennas by generating simulated observations for each component separately, and then combining the Measurement Sets in {{simanalyze_6.5.4}}. This technique is general and can be used to simulate observations using multiple 12 m array configurations, as well. Total power observations can be simulated either in an independent run of {{simobserve_6.5.4}}, or along with an interferometric simulation. Note that if you simulate total power and an interferometric observation simultaneously with {{simobserve_6.5.4}}, they must have the same set of pointing centers and the same integration and total time, which is probably not realistic. (For example it is generally recommended to observe a larger area by 1/2 primary beam in total power mode to combine with a 12 m ALMA mosaic). | We can simulate ALMA observations that use the main array of 12 m antennas augmented by the ACA and Total Power antennas by generating simulated observations for each component separately, and then combining the Measurement Sets in {{simanalyze_6.5.4}}. This technique is general and can be used to simulate observations using multiple 12 m array configurations, as well. Total power observations can be simulated either in an independent run of {{simobserve_6.5.4}}, or along with an interferometric simulation. Note that if you simulate total power and an interferometric observation simultaneously with {{simobserve_6.5.4}}, they must have the same set of pointing centers and the same integration and total time, which is probably not realistic. (For example it is generally recommended to observe a larger area by 1/2 primary beam in total power mode to combine with a 12 m ALMA mosaic). | ||
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]]<br> | ||
Note that this guide is written for interactive CASA, so the script won't run correctly with ''"casa -c script.py"'' or ''"execfile('script.py')"'' in this case. However, it will run by copy and pasting the commands into terminal. | |||
<source lang="python"> | |||
## interactive CASA ## | |||
# default("taskname") | |||
# parameter1 = value1 | |||
# parameter2 = value2 | |||
# taskname() | |||
## scripted CASA ## | |||
# taskname(parameter1=value1, parameter2=value2) | |||
</source> | |||
One advantage of interactive CASA is that python global variables will carry forward into following tasks, so they don't need to be redefined every time. Be aware that any defined parameter will override the default value for any of the following tasks that contain that parameter, unless the task is initiated with ''default("taskname")''. | |||
=====Set simobserve as current task===== | =====Set simobserve as current task===== | ||
Reset all parameters to default, and then set the project name to ''m51c'' | Reset all parameters to default, and then set the project name to ''m51c''. | ||
<source lang="python"> | <source lang="python"> | ||
# In CASA | |||
# Set simobserve to default parameters | # Set simobserve to default parameters | ||
default("simobserve") | default("simobserve") | ||
| Line 18: | Line 35: | ||
=====Specify sky model image===== | =====Specify sky model image===== | ||
We'll use an | We'll use an H-alpha image of M51 as a model of the sky, for this example. The ''curl'' command will copy the file from a URL and rename it. | ||
<source lang="python"> | <source lang="python"> | ||
# Model sky = | # In CASA | ||
# Model sky = H-alpha image of M51 | |||
import os | import os | ||
os.system('curl https://casaguides.nrao.edu/images/3/3f/M51ha.fits.txt -f -o M51ha.fits') | os.system('curl https://casaguides.nrao.edu/images/3/3f/M51ha.fits.txt -f -o M51ha.fits') | ||
| Line 32: | Line 51: | ||
<source lang="python"> | <source lang="python"> | ||
# In CASA | |||
# Set model image parameters: | # Set model image parameters: | ||
indirection = "B1950 23h59m59.96s -34d59m59.50s" | indirection = "B1950 23h59m59.96s -34d59m59.50s" | ||
| Line 55: | Line 76: | ||
<source lang="python"> | <source lang="python"> | ||
# In CASA | |||
# have simobserve calculate mosaic pointing locations: | # have simobserve calculate mosaic pointing locations: | ||
setpointings = True | setpointings = True | ||
| Line 64: | Line 87: | ||
=====Calculate Visibilities===== | =====Calculate Visibilities===== | ||
When you provide a desired resolution in the '''antennalist''' parameter, {{simobserve_6.5.4}} will determine what array configuration to use. | When you provide a desired resolution in the '''antennalist''' parameter, {{simobserve_6.5.4}} will determine what array configuration to use. | ||
<source lang="python"> | <source lang="python"> | ||
# In CASA | |||
obsmode = "int" | obsmode = "int" | ||
antennalist = "ALMA;0.5arcsec" | antennalist = "ALMA;0.5arcsec" | ||
| Line 74: | Line 100: | ||
<source lang="python"> | <source lang="python"> | ||
# In CASA | |||
graphics = "both" | graphics = "both" | ||
simobserve() | simobserve() | ||
| Line 88: | Line 116: | ||
We need to change the mapping parameters to specify a square region a bit larger than the interferometric map. | We need to change the mapping parameters to specify a square region a bit larger than the interferometric map. | ||
<source lang="python"> | <source lang="python"> | ||
# In CASA | |||
integration = "10s" | integration = "10s" | ||
mapsize = "1.3arcmin" | mapsize = "1.3arcmin" | ||
| Line 96: | Line 126: | ||
<source lang="python"> | <source lang="python"> | ||
# In CASA | |||
obsmode = "sd" | obsmode = "sd" | ||
sdantlist = "aca.tp.cfg" | sdantlist = "aca.tp.cfg" | ||
| Line 106: | Line 138: | ||
<source lang="python"> | <source lang="python"> | ||
# In CASA | |||
simobserve() | simobserve() | ||
</source> | </source> | ||
| Line 116: | Line 150: | ||
<source lang="python"> | <source lang="python"> | ||
# In CASA | |||
integration = "10s" | integration = "10s" | ||
mapsize = "1arcmin" | mapsize = "1arcmin" | ||
| Line 122: | Line 158: | ||
</source> | </source> | ||
We can specify an integral number of times to repeat the mosaic by setting | We can specify an integral number of times to repeat the mosaic by setting '''totaltime''' to an integer string without units. | ||
<source lang="python"> | <source lang="python"> | ||
# In CASA | |||
obsmode = "int" | obsmode = "int" | ||
refdate = "2012/12/02" | refdate = "2012/12/02" | ||
| Line 134: | Line 172: | ||
<source lang="python"> | <source lang="python"> | ||
# In CASA | |||
simobserve() | simobserve() | ||
</source> | </source> | ||
== Deconvolve the visibilities back into an image == | == Deconvolve the visibilities back into an image == | ||
Next we use {{simanalyze_6.5.4}} to combine the three measurement sets and create a single image. | Next we use {{simanalyze_6.5.4}} to combine the three measurement sets and create a single image. | ||
<font color="red"> | |||
WARNING: There is a [https://casadocs.readthedocs.io/en/latest/notebooks/introduction.html#simobserve-/-simanalyze known issue] with {{simanalyze_6.5.4}} in CASA 6. This section is validated for CASA 5.6. | |||
</font> | |||
There are many ways to do this, and you may wish to discuss options with scientists at your ARC. In this example we will use the total power image as a model when deconvolving the ACA image, and then use the result as a model when deconvolving the 12m interferometric image. This method tends to give low weight to the large spatial scales, but is simple to illustrate. | There are many ways to do this, and you may wish to discuss options with scientists at your ARC. In this example we will use the total power image as a model when deconvolving the ACA image, and then use the result as a model when deconvolving the 12m interferometric image. This method tends to give low weight to the large spatial scales, but is simple to illustrate. | ||
It's possible to get better results if one used multiscale clean in the | It's possible to get better results if one used multiscale clean in the {{tclean_6.5.4}} task (again using the lower resolution image as a model when deconvolving the higher resolution one). An alternative would be to create an image independently from each dataset, and then use the CASA {{feather_6.5.4}} task to combine them entirely in the image plane. | ||
If given a total power and interferometric measurement set, {{simanalyze_6.5.4}} will automatically create the total power image, | If given a total power and interferometric measurement set, {{simanalyze_6.5.4}} will automatically create the total power image, | ||
| Line 150: | Line 193: | ||
<div style="background-color: #E0FFFF; padding: 10px; font-family: 'Consolas', monospace;> | <div style="background-color: #E0FFFF; padding: 10px; font-family: 'Consolas', monospace;> | ||
If you decide to concatenate the two measurement sets and image all visibilities simultaneously, it is critical that the relative weights be set between the two different interferometric arrays. Simulated data has weights=1, since the thermal noise is generated uniformly per baseline. However, in reality the 7m baselines have lower sensitivity than the 12m baselines, and their weights must be decreased by that sensitivity ratio. | If you decide to concatenate the two measurement sets and image all visibilities simultaneously, it is critical that the relative weights be set between the two different interferometric arrays. Simulated data has weights=1, since the thermal noise is generated uniformly per baseline. However, in reality the 7m baselines have lower sensitivity than the 12m baselines, and their weights must be decreased by that sensitivity ratio. {{simalma_6.5.4}} uses the '''visweightscale''' parameter of {{concat_6.5.4}} to apply that lower weight of (7./12)**2 to the 7m visibilities. If you wish to combine data manually, you must do this step yourself (see the guide: [[simalma_CASA_6.5.4]]). | ||
</div> | </div> | ||
| Line 157: | Line 200: | ||
<source lang="python"> | <source lang="python"> | ||
# In CASA | |||
default("simanalyze") | default("simanalyze") | ||
project = "m51c" | project = "m51c" | ||
| Line 174: | Line 219: | ||
<source lang="python"> | <source lang="python"> | ||
# In CASA | |||
default("simanalyze") | default("simanalyze") | ||
project = "m51c" | project = "m51c" | ||
Latest revision as of 17:33, 27 March 2024
Last checked on CASA Version 6.5.4
Overview: ALMA 12m + ACA (7m) + Total Power
We can simulate ALMA observations that use the main array of 12 m antennas augmented by the ACA and Total Power antennas by generating simulated observations for each component separately, and then combining the Measurement Sets in simanalyze. This technique is general and can be used to simulate observations using multiple 12 m array configurations, as well. Total power observations can be simulated either in an independent run of simobserve, or along with an interferometric simulation. Note that if you simulate total power and an interferometric observation simultaneously with simobserve, they must have the same set of pointing centers and the same integration and total time, which is probably not realistic. (For example it is generally recommended to observe a larger area by 1/2 primary beam in total power mode to combine with a 12 m ALMA mosaic).
To create a script of the Python code on this page, see: Extracting scripts from these tutorials
Note that this guide is written for interactive CASA, so the script won't run correctly with "casa -c script.py" or "execfile('script.py')" in this case. However, it will run by copy and pasting the commands into terminal.
## interactive CASA ##
# default("taskname")
# parameter1 = value1
# parameter2 = value2
# taskname()
## scripted CASA ##
# taskname(parameter1=value1, parameter2=value2)
One advantage of interactive CASA is that python global variables will carry forward into following tasks, so they don't need to be redefined every time. Be aware that any defined parameter will override the default value for any of the following tasks that contain that parameter, unless the task is initiated with default("taskname").
Set simobserve as current task
Reset all parameters to default, and then set the project name to m51c.
# In CASA
# Set simobserve to default parameters
default("simobserve")
# Our project name will be m51c, and all simulation products will be placed in a subdirectory m51c/
project = "m51c"
Specify sky model image
We'll use an H-alpha image of M51 as a model of the sky, for this example. The curl command will copy the file from a URL and rename it.
# In CASA
# Model sky = H-alpha image of M51
import os
os.system('curl https://casaguides.nrao.edu/images/3/3f/M51ha.fits.txt -f -o M51ha.fits')
skymodel = "M51ha.fits"
Note that simobserve will not modify your original input image. Rather, it will make a copy: m51c/m51c.skymodel
We will override most of the parameters in the Halpha FITS image to make the image more suitable to a sub-millimeter ALMA observation.
# In CASA
# Set model image parameters:
indirection = "B1950 23h59m59.96s -34d59m59.50s"
incell = "0.1arcsec"
inbright = "0.004"
incenter = "330.076GHz"
inwidth = "50MHz"
These parameters are plausible for observing a sub-mm emission line in a galaxy:
- place the source in the southern hemisphere with the indirection parameter,
- set the pixel size to 0.1arcsec, to simulate an observation of a galaxy that is smaller in angular size. (M51 itself would require a quite large mosaic, and in any case we'd like the input model pixels to be significantly smaller than the synthesized beam.)
- set the peak brightness to 0.004 Jy/pixel
- set the frequency to 330GHz
- since it's a 2D image we'll set the single "channel" width to be 50MHz.
Simulate 12m interferometric observation
We'll begin by simulating the observation as seen by the main 12 m ALMA array. We'll have simobserve calculate a hexagonal mosaic of pointings.
The interface for simobserve provides an integration parameter, which is the dwell time at each pointing in the mosaic. We'll set that to 10s. A real observation might dwell on each pointing for ~20s and record 3 sets of visibilities (i.e. 6s integration time) for each pointing. simobserve does not currently offer a capability to reproduce this behaviour exactly, but the difference in UV coverage should be very minor between the real observation and the simulated one. You could try setting integration to ~20s, which would more closely match the dwell time for some real observations.
# In CASA
# have simobserve calculate mosaic pointing locations:
setpointings = True
integration = "10s"
mapsize = "1arcmin"
maptype = "hex"
pointingspacing = "9arcsec" # this could also be specified in units of the primary beam, e.g. "0.5PB"
Calculate Visibilities
When you provide a desired resolution in the antennalist parameter, simobserve will determine what array configuration to use.
# In CASA
obsmode = "int"
antennalist = "ALMA;0.5arcsec"
totaltime = "3600s"
Run simobserve, displaying graphics to the screen and also writing graphics output to files that can be found in the project directory, m51c.
# In CASA
graphics = "both"
simobserve()
Simulate 12m total power observation
Next we'll simulate a total power raster map of the same area, on a square grid. CASA simulation tools can not simulate true on-the-fly mapping (with smearing on timescales smaller than an integration time), but a square grid with a short integration time will provide a very accurate approximation.
By virtue of CASA's global parameters, we already have project and image world coordinate system parameters set correctly.
We need to change the mapping parameters to specify a square region a bit larger than the interferometric map.
# In CASA
integration = "10s"
mapsize = "1.3arcmin"
maptype = "square"
We'll observe on a different day. This doesn't really matter, but if you choose to simulate two different 12m ALMA configurations and combine them, you will have issues concatenating the datasets if they are simulated on the same day with the same antenna names. So, it's a good habit to change the day.
# In CASA
obsmode = "sd"
sdantlist = "aca.tp.cfg"
sdant = 0
refdate = "2012/12/01"
totaltime = "2h"
Run simobserve, displaying graphics to screen and to files.
# In CASA
simobserve()
Simulate 7m ACA observation
Next we'll add an ACA mosaic, with its larger primary beam.
# In CASA
integration = "10s"
mapsize = "1arcmin"
maptype = "hex"
pointingspacing = "15arcsec"
We can specify an integral number of times to repeat the mosaic by setting totaltime to an integer string without units.
# In CASA
obsmode = "int"
refdate = "2012/12/02"
antennalist = "aca.i.cfg"
totaltime = "3"
Run simobserve, displaying graphics to screen and to files.
# In CASA
simobserve()
Deconvolve the visibilities back into an image
Next we use simanalyze to combine the three measurement sets and create a single image.
WARNING: There is a known issue with simanalyze in CASA 6. This section is validated for CASA 5.6.
There are many ways to do this, and you may wish to discuss options with scientists at your ARC. In this example we will use the total power image as a model when deconvolving the ACA image, and then use the result as a model when deconvolving the 12m interferometric image. This method tends to give low weight to the large spatial scales, but is simple to illustrate.
It's possible to get better results if one used multiscale clean in the tclean task (again using the lower resolution image as a model when deconvolving the higher resolution one). An alternative would be to create an image independently from each dataset, and then use the CASA feather task to combine them entirely in the image plane.
If given a total power and interferometric measurement set, simanalyze will automatically create the total power image, then use it as a model and deconvolve the interferometric image.
If you decide to concatenate the two measurement sets and image all visibilities simultaneously, it is critical that the relative weights be set between the two different interferometric arrays. Simulated data has weights=1, since the thermal noise is generated uniformly per baseline. However, in reality the 7m baselines have lower sensitivity than the 12m baselines, and their weights must be decreased by that sensitivity ratio. simalma uses the visweightscale parameter of concat to apply that lower weight of (7./12)**2 to the 7m visibilities. If you wish to combine data manually, you must do this step yourself (see the guide: simalma_CASA_6.5.4).
First image total power and ACA with total power as a model
# In CASA
default("simanalyze")
project = "m51c"
vis = '$project.aca.i.ms,$project.aca.tp.sd.ms'
imsize = [512,512]
cell = '0.2arcsec'
analyze = True
showpsf = False
showresidual = False
showconvolved = True
simanalyze()
Next add the 12m interferometric data
Here we explicitly have to set the modelimage to the result from the previous run.
# In CASA
default("simanalyze")
project = "m51c"
vis = '$project.ALMA_0.5arcsec.ms'
imsize = [512,512]
cell = '0.1arcsec'
modelimage = "$project.aca.i.image"
analyze = True
showpsf = False
showresidual = False
showconvolved = True
simanalyze()
