############################################################################# # # # A VLA (B-config) 43 GHz continuum fast-switching observation of M82 # # # # Initial by GM 2009-05-27 (Beta 2.4) init minimally annotated version # # Revised by GM 2009-06-01 added opacity=0.067, gaincurve=T # # added scan=based phase solution # # use weighting='briggs' in clean # # Revised by GM 2009-06-02 added more plotting, split m82 # # removed 'briggs' weighting # # # # This calibration example is a good demonstration of calibration and # # and imaging of a reasonably phase-coherent 43 GHz VLA continuum # # fast-switching observation of M82. The observation is 2h50m # # long, uses 3.33s integrations, has 2 spws, and the dataset has # # ~1.5M visibilities (time/spw/baseline points). The phases are # # interestingly variable, but quite tolerable and easily calibrated. # # Although reference pointing was done every few 90 min or so, the # # deterioration of the pointing is evident in the data and calibration # # amplitudes as a function of time after each pointing scan (the pointing # # scans are not filled, but the gaps where they occur are readily evident). # # # # The resulting image reveals ~dozen or so discrete sources in M82 (SNRs # # presumably); in the B-array, the large-scale continuum structure # # of the galaxy is filtered out (see the D-config M82 example). The # # comparison of pre-calibration, dirty and cleaned images is a # # compelling demonstration of the important of calibration and # # deconvolution. # # # # Required for this script (in the directory from which it is run): # # - this script # # - The AT303_C050423.xp1 VLA archive file # # - The 3C147_Q.im model image from the data repository # # # # Wishlist: # # - Add some cleanup and other niceties # # - Pending fluxscale improvements, do a selective fluxscale # # to avoid using solutions from poorly pointed data to set # # the f.d. scale. # # - Improve the imaging with robust weighting? taper? clean boxes? # # - Add some image analysis (peak, rms, etc.) # # - Make a regression out of this reduction # # # ############################################################################# # some file and table names to use below arch='AT303_C050423.xp1' fdmodel='3C147_Q.im' msname='m82B.ms' phasecal1='m82B.intphase.G' ampcal='m82B.scanamp.G' phasecal2='m82B.scanphase.G' fluxcal=ampcal+'flx' # Fill from VLA archive file #---------------------------------------- importvla(archivefiles=arch,vis=msname) # Scan listing, initial plots #---------------------------------------- listobs(msname) plotxy(vis=msname,xaxis='time',yaxis='amp',selectdata=T,correlation="RR,LL") plotxy(vis=msname,xaxis='time',yaxis='phase',selectdata=T,correlation="RR,LL") plotxy(vis=msname,xaxis='time',yaxis='phase',selectdata=T,correlation="RR,LL",antenna='VA01 & VA02',plotsymbol='o') # Basic commanded flagging #---------------------------------------- # quack (to lop off first few integrations in each # scan; all "low" integrations are already flagged online, # but there may be one or two integrations that # have only a few unflagged baselines; NB: the quackinterval # is a bit longer than what one would use in AIPS because # casa fills a bit more of the online-flagged data) flagdata(vis=msname,mode='quack',quackinterval=14) # zoom in on this plot to see the effect of quack plotxy(vis=msname,xaxis='time',yaxis='amp',selectdata=T,correlation="RR,LL") # Pre-calibration Images #---------------------------------------- # uncal'd calibrator clean(vis=msname,imagename='m82B.uncal_cal_dirty', imsize=[512,512], cell=['0.03arcsec','0.03arcsec'], niter=0, field='2') # uncal'd m82 clean(vis=msname,imagename='m82B.uncal_dirty', imsize=[2048,2048], cell=['0.03arcsec','0.03arcsec'], niter=0, field='3') # view viewer('m82B.uncal_cal_dirty') # Calibration! #---------------------------------------- The calibration steps: 1. Set the f.d. calibrator model (the gain calibrator already set to 1) 2. Solve for per-integration phase 3. Solve for per-scan amplitude, using the per-int phase 4. Solve for per-scan phase, using amplitude (discard per-int phase) Notes: 1. Using opacity=0.067 2. Using VLA gain curve 3. Reference antenna (phase) VA08 (near array center) # Set f.d. model: setjy(vis=msname,field='0542+498',modimage=fdmodel) # f.d. calibrator, model and data side-by-side plotxy(vis=msname,xaxis='uvdist',yaxis='amp',datacolumn='model',field='0',selectdata=T,correlation="RR,LL",subplot=121,plotrange=[0,0,0,1]) plotxy(vis=msname,xaxis='uvdist',yaxis='amp',datacolumn='data',field='0',selectdata=T,correlation="RR,LL",subplot=122,plotrange=[0,0,0,1]) # gain calibrator, model and data side-by-side plotxy(vis=msname,xaxis='uvdist',yaxis='amp',datacolumn='model',field='1,2',selectdata=T,correlation="RR,LL",subplot=121,plotrange=[0,0,0,1]) plotxy(vis=msname,xaxis='uvdist',yaxis='amp',datacolumn='data',field='1,2',selectdata=T,correlation="RR,LL",subplot=122,plotrange=[0,0,0,1]) # gaincal (fast phase-only) gaincal(vis=msname,caltable=phasecal1, field='0,1,2', calmode='p',solint='int', refant='VA08', opacity=0.067, gaincurve=T) # plot the solution plotcal(caltable=phasecal1,yaxis='phase') plotcal(caltable=phasecal1,yaxis='snr') plotcal(caltable=phasecal1,yaxis='phase',subplot=313,iteration='antenna') # gaincal (slower amp-only, using phase from above) gaincal(vis=msname,caltable=ampcal, field='0,1,2', calmode='a',solint='inf', refant='VA08', gaintable=phasecal1,interp='nearest', opacity=0.067, gaincurve=T) # plot the solution plotcal(caltable=ampcal,yaxis='amp') plotcal(caltable=ampcal,yaxis='snr') plotcal(caltable=ampcal,yaxis='amp',subplot=313,iteration='antenna',plotrange=[0,0,0,0.45]) # A more modest timescale phase solution for better SNR, phase-tracking gaincal(vis=msname,caltable=phasecal2, field='0,1,2', calmode='p',solint='50s', refant='VA08', gaintable=ampcal,interp='nearest', opacity=0.067, gaincurve=T) # plot the solution plotcal(caltable=phasecal2,yaxis='phase') plotcal(caltable=phasecal2,yaxis='snr') plotcal(caltable=phasecal2,yaxis='phase',subplot=313,iteration='antenna',plotrange=[0,0,-180,180]) # fluxscale to transfer proper scale to phase-calibrator # (this isn't quite ideal since it uses all of the phase-calibrator # scans and all antennas, and many of these solutions are poor due # to deteriorating pointing....) fluxscale(vis=msname, caltable=ampcal,fluxtable=fluxcal, reference='0542+498') # plot the solutions, before and after plotcal(caltable=ampcal,yaxis='amp',subplot=121,plotrange=[0,0,0,0.45]) plotcal(caltable=fluxcal,yaxis='amp',subplot=122,plotrange=[0,0,0,0.45]) # accum to show interpolation effects in plocal # (only for illustration) # accumulate phase on calibrators accum(vis=msname, incrtable=phasecal2, caltable='m82B.accum1', field='0,1,2',calfield='0,1,2',interp='nearest') # ...and onto M82 accum(vis=msname,tablein='m82B.accum1', incrtable=phasecal2, caltable='m82B.accum2', field='3',calfield='1,2',interp='linear') # accumulate amplitude on calibrators accum(vis=msname,tablein='m82B.accum2', incrtable=fluxcal, caltable='m82B.accum3', field='0,1,2',calfield='0,1,2',interp='nearest') # ...and onto M82 accum(vis=msname,tablein='m82B.accum3', incrtable=fluxcal, caltable='m82B.accum4', field='3',calfield='1,2',interp='linear') # plot accum'd solutions, zoom in on pieces to see what's happening plotcal(caltable='m82B.accum4',yaxis='amp',antenna='VA02',plotrange=[0,0,0,0.45]) plotcal(caltable='m82B.accum4',yaxis='phase',antenna='VA02') # apply the calibration to each source # (alternatively could use m82B.accum4 from above) # calibrators w/ nearest interpolation applycal(vis=msname, field='0,1,2', gaintable=[phasecal2,fluxcal], interp=['nearest','nearest'], opacity=0.067, gaincurve=T) # target with linear interpolation applycal(vis=msname, field='3', gaintable=[phasecal2,fluxcal], gainfield=[[1,2],[1,2]], interp=['linear','linear'], opacity=0.067, gaincurve=T) # f.d. calibrator, model and corrected data side-by-side # AMP plotxy(vis=msname,xaxis='uvdist',yaxis='amp',datacolumn='data',field='0',selectdata=T,correlation="RR,LL",subplot=131,plotrange=[0,0,0,1.2],timebin='500') plotxy(vis=msname,xaxis='uvdist',yaxis='amp',datacolumn='model',field='0',selectdata=T,correlation="RR,LL",subplot=132,plotrange=[0,0,0,1.2],timebin='500') plotxy(vis=msname,xaxis='uvdist',yaxis='amp',datacolumn='corrected',field='0',selectdata=T,correlation="RR,LL",subplot=133,plotrange=[0,0,0,1.2],timebin='500') # PHASE plotxy(vis=msname,xaxis='uvdist',yaxis='phase',datacolumn='data',field='0',selectdata=T,correlation="RR,LL",subplot=131,plotrange=[0,0,-180,180],timebin='500') plotxy(vis=msname,xaxis='uvdist',yaxis='phase',datacolumn='model',field='0',selectdata=T,correlation="RR,LL",subplot=132,plotrange=[0,0,-180,180],timebin='500') plotxy(vis=msname,xaxis='uvdist',yaxis='phase',datacolumn='corrected',field='0',selectdata=T,correlation="RR,LL",subplot=133,plotrange=[0,0,-180,180],timebin='500') # gain calibrator, model and corrected data side-by-side # AMP plotxy(vis=msname,xaxis='uvdist',yaxis='amp',datacolumn='data',field='1,2',selectdata=T,correlation="RR,LL",subplot=131,plotrange=[0,0,0,1.5],timebin='100') plotxy(vis=msname,xaxis='uvdist',yaxis='amp',datacolumn='model',field='1,2',selectdata=T,correlation="RR,LL",subplot=132,plotrange=[0,0,0,1.5],timebin='100') plotxy(vis=msname,xaxis='uvdist',yaxis='amp',datacolumn='corrected',field='1,2',selectdata=T,correlation="RR,LL",subplot=133,plotrange=[0,0,0,1.5],timebin='100') # PHASE plotxy(vis=msname,xaxis='uvdist',yaxis='phase',datacolumn='data',field='1,2',selectdata=T,correlation="RR,LL",subplot=131,plotrange=[0,0,-180,180],timebin='100') plotxy(vis=msname,xaxis='uvdist',yaxis='phase',datacolumn='model',field='1,2',selectdata=T,correlation="RR,LL",subplot=132,plotrange=[0,0,-180,180],timebin='100') plotxy(vis=msname,xaxis='uvdist',yaxis='phase',datacolumn='corrected',field='1,2',selectdata=T,correlation="RR,LL",subplot=133,plotrange=[0,0,-180,180],timebin='100') # clean to make dirty image: # calibrator clean(vis=msname,imagename='m82B.cal_dirty', imsize=[512,512], cell=['0.03arcsec','0.03arcsec'], niter=0, field='2') # m82 clean(vis=msname,imagename='m82B.dirty', imsize=[2048,2048], cell=['0.03arcsec','0.03arcsec'], niter=0, field='3') # split out m82 for careful imaging... split(vis=msname,outputvis='m82B_m82only.ms',field='3'); # (see Crystal's demo for imaging/deconvolution...) # end