|
casa
$Rev:20696$
|
|
casa
$Rev:20696$
|
00001 ########################################################################## 00002 # listcal_regression.py 00003 # 00004 # Regression test for task listcal. 00005 # 00006 # Current number of test scenarios: 2 00007 # 00008 # The testing is performed by comparing the listcal output at runtime 00009 # with "standards" contained in the data repository. Metadata (all 00010 # output that is not a floating point number) is compared character by 00011 # character; any difference will cause the test to fail. Data (all 00012 # floating point numbers in output) are assumed to be sequential 00013 # amplitude-phase pairs and are required to be equal to within a 00014 # minimum precision value specified below in this script. 00015 # 00016 # Data files are built from scratch using functions defined in this script. 00017 00018 import time 00019 import listing as lt 00020 import regression_utility as tstutl 00021 00022 startTime = time.time() 00023 00024 print "BEGIN: listcal_regression.py" 00025 00026 testPassed = 0 00027 testFailed = 0 00028 pathName = os.environ.get('CASAPATH').split()[0] 00029 localData = pathName + '/data/regression/ngc4826/' 00030 automate = true # set to false for testing or debugging 00031 regressionDir = 'listcal_regression' 00032 if (not os.path.exists(regressionDir)): os.mkdir(regressionDir) 00033 00034 if(automate): 00035 print "Running in automated mode." 00036 print " - All MS data will be rebuilt from scratch." 00037 print " - All test files will be removed after testing." 00038 else: 00039 print "Running in non-automated mode!" 00040 00041 # For testing: 00042 # casalog.filter('DEBUG2') 00043 00044 #============================================================================= 00045 # METHOD: load_ngc4826 00046 # 00047 # Load, edit, and calibrate NGC4826 data. 00048 # 00049 def load_ngc4826(prefix,msname,caltable): 00050 # Clear out previous run results 00051 os.system('rm -rf '+prefix+'ngc4826.tutorial.*') 00052 ########################################################################## 00053 # Import and concatenate sources 00054 # Data Description: 00055 # N4826 - BIMA SONG Data # 00056 # This data is from the BIMA Survey of Nearby Galaxies (BIMA SONG) # 00057 # Helfer, Thornley, Regan, et al., 2003, ApJS, 145, 259 # 00058 # 16apr98 00059 # source=ngc4826 00060 # phasecal=1310+323 00061 # fluxcal=3c273, Flux = 23 Jy on 16apr98 00062 # passcal= none - data were observed with online bandpass correction. 00063 # NOTE: This data has been filled into MIRIAD, line-length correction 00064 # done, and then exported as separate files for each source. 00065 # 3c273 was not line length corrected since it was observed 00066 # for such a short amount of time that it did not need it. 00067 # From miriad: source Vlsr = 408; delta V is 20 km/s 00068 # NOTE: This data contains correlations of only one polarization, 'YY'. 00069 # The antennas contain 'X' and 'Y' feeds, but the 'X' data was not 00070 # correlated. 00071 ########################################################################## 00072 # [ This section derived from the 2008 summer school tutorial. ] 00073 # [ http://casa.nrao.edu/Tutorial/SIworkshop2008/Scripts/ngc4826_tutorial.py ] 00074 # USB spectral windows written separately by miriad for 16apr98 00075 # Assumes these are in sub-directory called "fitsfiles" of working directory 00076 print '--Importuvfits (16apr98)--' 00077 default('importuvfits') 00078 print "Starting from the uvfits files exported by miriad" 00079 print "The USB spectral windows were written separately by miriad for 16apr98" 00080 pathName = os.environ.get('CASAPATH').split()[0] 00081 localData = pathName + '/data/regression/ngc4826/' 00082 importuvfits(fitsfile= localData + 'fitsfiles/3c273.fits5', vis=prefix+'ngc4826.tutorial.3c273.5.ms') 00083 importuvfits(fitsfile= localData + 'fitsfiles/3c273.fits6', vis=prefix+'ngc4826.tutorial.3c273.6.ms') 00084 importuvfits(fitsfile= localData + 'fitsfiles/3c273.fits7', vis=prefix+'ngc4826.tutorial.3c273.7.ms') 00085 importuvfits(fitsfile= localData + 'fitsfiles/3c273.fits8', vis=prefix+'ngc4826.tutorial.3c273.8.ms') 00086 importuvfits(fitsfile= localData + 'fitsfiles/1310+323.ll.fits9', vis=prefix+'ngc4826.tutorial.1310+323.ll.9.ms') 00087 importuvfits(fitsfile= localData + 'fitsfiles/1310+323.ll.fits10', vis=prefix+'ngc4826.tutorial.1310+323.ll.10.ms') 00088 importuvfits(fitsfile= localData + 'fitsfiles/1310+323.ll.fits11', vis=prefix+'ngc4826.tutorial.1310+323.ll.11.ms') 00089 importuvfits(fitsfile= localData + 'fitsfiles/1310+323.ll.fits12', vis=prefix+'ngc4826.tutorial.1310+323.ll.12.ms') 00090 importuvfits(fitsfile= localData + 'fitsfiles/1310+323.ll.fits13', vis=prefix+'ngc4826.tutorial.1310+323.ll.13.ms') 00091 importuvfits(fitsfile= localData + 'fitsfiles/1310+323.ll.fits14', vis=prefix+'ngc4826.tutorial.1310+323.ll.14.ms') 00092 importuvfits(fitsfile= localData + 'fitsfiles/1310+323.ll.fits15', vis=prefix+'ngc4826.tutorial.1310+323.ll.15.ms') 00093 importuvfits(fitsfile= localData + 'fitsfiles/1310+323.ll.fits16', vis=prefix+'ngc4826.tutorial.1310+323.ll.16.ms') 00094 importuvfits(fitsfile= localData + 'fitsfiles/ngc4826.ll.fits5', vis=prefix+'ngc4826.tutorial.ngc4826.ll.5.ms') 00095 importuvfits(fitsfile= localData + 'fitsfiles/ngc4826.ll.fits6', vis=prefix+'ngc4826.tutorial.ngc4826.ll.6.ms') 00096 importuvfits(fitsfile= localData + 'fitsfiles/ngc4826.ll.fits7', vis=prefix+'ngc4826.tutorial.ngc4826.ll.7.ms') 00097 importuvfits(fitsfile= localData + 'fitsfiles/ngc4826.ll.fits8', vis=prefix+'ngc4826.tutorial.ngc4826.ll.8.ms') 00098 ########################################################################## 00099 print '--Concat--' 00100 default('concat') 00101 concat(vis=[prefix+'ngc4826.tutorial.3c273.5.ms', 00102 prefix+'ngc4826.tutorial.3c273.6.ms', 00103 prefix+'ngc4826.tutorial.3c273.7.ms', 00104 prefix+'ngc4826.tutorial.3c273.8.ms', 00105 prefix+'ngc4826.tutorial.1310+323.ll.9.ms', 00106 prefix+'ngc4826.tutorial.1310+323.ll.10.ms', 00107 prefix+'ngc4826.tutorial.1310+323.ll.11.ms', 00108 prefix+'ngc4826.tutorial.1310+323.ll.12.ms', 00109 prefix+'ngc4826.tutorial.1310+323.ll.13.ms', 00110 prefix+'ngc4826.tutorial.1310+323.ll.14.ms', 00111 prefix+'ngc4826.tutorial.1310+323.ll.15.ms', 00112 prefix+'ngc4826.tutorial.1310+323.ll.16.ms', 00113 prefix+'ngc4826.tutorial.ngc4826.ll.5.ms', 00114 prefix+'ngc4826.tutorial.ngc4826.ll.6.ms', 00115 prefix+'ngc4826.tutorial.ngc4826.ll.7.ms', 00116 prefix+'ngc4826.tutorial.ngc4826.ll.8.ms'], 00117 concatvis=msname, 00118 freqtol="",dirtol="1arcsec",async=False) 00119 ########################################################################## 00120 # Fix up the MS (temporary, changes to importfits underway) 00121 print '--Fixing up spw rest frequencies in MS--' 00122 vis=msname 00123 tb.open(vis+'/SOURCE',nomodify=false) 00124 spwid=tb.getcol('SPECTRAL_WINDOW_ID') 00125 #spwid.setfield(-1,int) 00126 # 64bit imported from ngc4826_tutorial_regression 00127 spwid.setfield(-1,'int32') 00128 tb.putcol('SPECTRAL_WINDOW_ID',spwid) 00129 tb.close() 00130 # This ensures that the rest freq will be found for all spws. 00131 ########################################################################## 00132 # 16 APR Calibration 00133 ########################################################################## 00134 print '--Clearcal--' 00135 print 'Create scratch columns and initialize in '+msname 00136 # Force create/initialize of scratch columns 00137 # NOTE: plotxy will not run properly without this step. 00138 clearcal(vis=msname) 00139 # But this data is relatively clean, and flagging will not improve results. 00140 ########################################################################## 00141 # Flag end channels 00142 print '--Flagdata--' 00143 default('tflagdata') 00144 print "" 00145 print "Flagging edge channels in all spw" 00146 print " 0~3:0~1;62~63 , 4~11:0~1;30~31, 12~15:0~1;62~63 " 00147 print "" 00148 tflagdata(vis=msname, mode='manual', 00149 spw='0~3:0;1;62;63,4~11:0;1;30;31,12~15:0;1;62;63') 00150 # Flag correlator glitch 00151 print "" 00152 print "Flagging bad correlator field 8 antenna 3&9 spw 15 all channels" 00153 print " timerange 1998/04/16/06:19:00.0~1998/04/16/06:20:00.0" 00154 print "" 00155 tflagdata(vis=msname, mode='manual', field='8', spw='15', antenna='3&9', timerange='1998/04/16/06:19:00.0~1998/04/16/06:20:00.0') 00156 print "Completed pre-calibration flagging" 00157 ########################################################################## 00158 # Use Flagmanager to save a copy of the flags so far 00159 print '--Flagmanager--' 00160 default('flagmanager') 00161 print "Now will use flagmanager to save a copy of the flags we just made" 00162 print "These are named myflags" 00163 flagmanager(vis=msname,mode='save',versionname='myflags', 00164 comment='My flags',merge='replace') 00165 ########################################################################## 00166 # CALIBRATION 00167 ########################################################################## 00168 # Bandpasses are very flat because of observing mode used (online bandpass 00169 # correction) so bandpass calibration is unecessary for these data. 00170 ########################################################################## 00171 # Derive gain calibration solutions. 00172 # We will use VLA-like G (per-scan) calibration: 00173 ########################################################################## 00174 # Set the flux density of 3C273 to 23 Jy 00175 print '--Setjy (3C273)--' 00176 default('setjy') 00177 setjy(vis=msname,field='0',fluxdensity=[23.0,0.,0.,0.],spw='0~3',scalebychan=False,standard='Perley-Taylor 99') 00178 # Not really necessary to set spw but you get lots of warning messages if 00179 # you don't 00180 ########################################################################## 00181 # Gain calibration 00182 print '--Gaincal--' 00183 default('gaincal') 00184 # This should be combining all spw for the two calibrators for single 00185 # scan-based solutions 00186 print 'Gain calibration for fields 0,1 and spw 0~11' 00187 print 'Using solint=inf combining over spw' 00188 print 'Output table ngc4826.tutorial.16apr98.gcal' 00189 gaincal(vis=msname, caltable=caltable, 00190 field='0,1', spw='0~11', gaintype='G', minsnr=2.0, 00191 refant='ANT5', gaincurve=False, opacity=0.0, 00192 solint='inf', combine='spw') 00193 #============================================================================= 00194 00195 #============================================================================= 00196 # 00197 # The use of this calibration function has not yet been implemented. 00198 # It can be used to generate and MS and calibration tables to test listcal. 00199 # 00200 def load_jupiter6cm(prefix,msname,caltable): 00201 ###################################################################### 00202 # Use Case Script for Jupiter 6cm VLA # 00203 # Trimmed down from Use Case jupiter6cm_usecase.py # 00204 # Updated STM 2008-05-15 (Beta Patch 2.0) # 00205 # Updated STM 2008-06-11 (Beta Patch 2.0) # 00206 # Updated STM 2008-06-12 (Beta Patch 2.0) for summer school demo # 00207 # This is a VLA 6cm dataset that was observed in 1999 to set the # 00208 # flux scale for calibration of the VLA. Included in the program # 00209 # were observations of the planets, including Jupiter. # 00210 # This is D-configuration data, with resolution of around 14" # 00211 # Includes polarization imaging and analysis # 00212 ###################################################################### 00213 import time 00214 import os 00215 #===================================================================== 00216 # This script has some interactive commands: scriptmode = True 00217 # if you are running it and want it to stop during interactive parts. 00218 scriptmode = False 00219 ## #===================================================================== 00220 ## # Set up some useful variables - these will be set during the script 00221 ## # also, but if you want to restart the script in the middle here 00222 ## # they are in one place: 00223 ## 00224 ## # This will prefix all output file names 00225 ## prefix='jupiter6cm.demo' 00226 ## 00227 ## # Clean up old files 00228 ## os.system('rm -rf '+prefix+'*') 00229 ## 00230 ## # This is the output MS file name 00231 ## msfile = prefix + '.ms' 00232 ## 00233 ## # 00234 #===================================================================== 00235 # Calibration variables 00236 # 00237 # Use same prefix as rest of script 00238 calprefix = prefix 00239 00240 # spectral windows to process 00241 usespw = '' 00242 usespwlist = ['0','1'] 00243 00244 # prior calibration to apply 00245 usegaincurve = True 00246 gainopacity = 0.0 00247 00248 # reference antenna 11 (11=VLA:N1) 00249 calrefant = '11' 00250 00251 gtable = calprefix + '.gcal' 00252 ftable = calprefix + '.fluxscale' 00253 atable = calprefix + '.accum' 00254 00255 # 00256 #===================================================================== 00257 # Polarization calibration setup 00258 # 00259 dopolcal = True 00260 00261 ptable = calprefix + '.pcal' 00262 xtable = calprefix + '.polx' 00263 00264 # Pol leakage calibrator 00265 poldfield = '0137+331' 00266 00267 # Pol angle calibrator 00268 polxfield = '1331+305' 00269 # At Cband the fractional polarization of this source is 0.112 and 00270 # the R-L PhaseDiff = 66deg (EVPA = 33deg) 00271 polxfpol = 0.112 00272 polxrlpd_deg = 66.0 00273 # Dictionary of IPOL in the spw 00274 polxipol = {'0' : 7.462, 00275 '1' : 7.510} 00276 00277 # Make Stokes lists for setjy 00278 polxiquv = {} 00279 for spw in ['0','1']: 00280 ipol = polxipol[spw] 00281 fpol = polxfpol 00282 ppol = ipol*fpol 00283 rlpd = polxrlpd_deg*pi/180.0 00284 qpol = ppol*cos(rlpd) 00285 upol = ppol*sin(rlpd) 00286 polxiquv[spw] = [ipol,qpol,upol,0.0] 00287 00288 # 00289 # Split output setup 00290 # 00291 srcname = 'JUPITER' 00292 srcsplitms = calprefix + '.' + srcname + '.split.ms' 00293 calname = '0137+331' 00294 calsplitms = calprefix + '.' + calname + '.split.ms' 00295 00296 # 00297 #===================================================================== 00298 # 00299 # Intensity imaging parameters 00300 # 00301 # Same prefix for this imaging demo output 00302 # 00303 imprefix = prefix 00304 00305 # This is D-config VLA 6cm (4.85GHz) obs 00306 # Check the observational status summary 00307 # Primary beam FWHM = 45'/f_GHz = 557" 00308 # Synthesized beam FWHM = 14" 00309 # RMS in 10min (600s) = 0.06 mJy (thats now, but close enough) 00310 00311 # Set the output image size and cell size (arcsec) 00312 # 4" will give 3.5x oversampling 00313 clncell = [4.,4.] 00314 00315 # 280 pix will cover to 2xPrimaryBeam 00316 # clean will say to use 288 (a composite integer) for efficiency 00317 clnalg = 'clark' 00318 clnmode = '' 00319 # For Cotton-Schwab use 00320 clnmode = 'csclean' 00321 clnimsize = [288,288] 00322 00323 # iterations 00324 clniter = 10000 00325 00326 # Also set flux residual threshold (0.04 mJy) 00327 # From our listobs: 00328 # Total integration time = 85133.2 seconds 00329 # With rms of 0.06 mJy in 600s ==> rms = 0.005 mJy 00330 # Set to 10x thermal rms 00331 clnthreshold=0.05 00332 00333 # 00334 # Filenames 00335 # 00336 imname1 = imprefix + '.clean1' 00337 clnimage1 = imname1+'.image' 00338 clnmodel1 = imname1+'.model' 00339 clnresid1 = imname1+'.residual' 00340 clnmask1 = imname1+'.clean_interactive.mask' 00341 00342 imname2 = imprefix + '.clean2' 00343 clnimage2 = imname2+'.image' 00344 clnmodel2 = imname2+'.model' 00345 clnresid2 = imname2+'.residual' 00346 clnmask2 = imname2+'.clean_interactive.mask' 00347 00348 imname3 = imprefix + '.clean3' 00349 clnimage3 = imname3+'.image' 00350 clnmodel3 = imname3+'.model' 00351 clnresid3 = imname3+'.residual' 00352 clnmask3 = imname3+'.clean_interactive.mask' 00353 00354 # 00355 # Selfcal parameters 00356 # 00357 # reference antenna 11 (11=VLA:N1) 00358 calrefant = '11' 00359 00360 # 00361 # Filenames 00362 # 00363 selfcaltab1 = imprefix + '.selfcal1.gtable' 00364 00365 selfcaltab2 = imprefix + '.selfcal2.gtable' 00366 smoothcaltab2 = imprefix + '.smoothcal2.gtable' 00367 00368 # 00369 #===================================================================== 00370 # 00371 # Polarization imaging parameters 00372 # 00373 # New prefix for polarization imaging output 00374 # 00375 polprefix = prefix + '.polimg' 00376 00377 # Set up clean slightly differently 00378 polclnalg = 'hogbom' 00379 polclnmode = 'csclean' 00380 00381 polimname = polprefix + '.clean' 00382 polimage = polimname+'.image' 00383 polmodel = polimname+'.model' 00384 polresid = polimname+'.residual' 00385 polmask = polimname+'.clean_interactive.mask' 00386 00387 # 00388 # Other files 00389 # 00390 ipolimage = polimage+'.I' 00391 qpolimage = polimage+'.Q' 00392 upolimage = polimage+'.U' 00393 00394 poliimage = polimage+'.poli' 00395 polaimage = polimage+'.pola' 00396 00397 # 00398 #===================================================================== 00399 # Start processing 00400 #===================================================================== 00401 # Get to path to the CASA home and stip off the name 00402 pathname=os.environ.get('CASAPATH').split()[0] 00403 # This is where the UVFITS data should be 00404 fitsdata=pathname+'/data/regression/jupiter6cm/jupiter6cm.fits' 00405 #===================================================================== 00406 # Data Import and List 00407 #===================================================================== 00408 # Import the data from FITS to MS 00409 print '--Import--' 00410 # Safest to start from task defaults 00411 default('importuvfits') 00412 print "Use importuvfits to read UVFITS and make an MS" 00413 # Set up the MS filename and save as new global variable 00414 msfile = prefix + '.ms' 00415 print "MS will be called "+msfile 00416 # Use task importuvfits 00417 fitsfile = fitsdata 00418 vis = msfile 00419 importuvfits() 00420 #===================================================================== 00421 # Data Examination and Flagging 00422 # REMOVED: All flagging was interactive. Could be replaced with 00423 # automatic flagging. 00424 #===================================================================== 00425 # Calibration 00426 #===================================================================== 00427 # Set the fluxes of the primary calibrator(s) 00428 # 00429 print '--Setjy--' 00430 default('setjy') 00431 00432 print "Use setjy to set flux of 1331+305 (3C286)" 00433 00434 vis = msfile 00435 00436 # 00437 # 1331+305 = 3C286 is our primary calibrator 00438 field = '1331+305' 00439 00440 scalebychan=False 00441 00442 # Setjy knows about this source so we dont need anything more 00443 standard='Perley-Taylor 99' 00444 00445 setjy() 00446 00447 # 00448 # You should see something like this in the logger and casapy.log file: 00449 # 00450 # 1331+305 spwid= 0 [I=7.462, Q=0, U=0, V=0] Jy, (Perley-Taylor 99) 00451 # 1331+305 spwid= 1 [I=7.51, Q=0, U=0, V=0] Jy, (Perley-Taylor 99) 00452 # 00453 00454 print "Look in logger for the fluxes (should be 7.462 and 7.510 Jy)" 00455 00456 # 00457 #===================================================================== 00458 # 00459 # Initial gain calibration 00460 # 00461 print '--Gaincal--' 00462 default('gaincal') 00463 00464 print "Solve for antenna gains on 1331+305 and 0137+331" 00465 print "We have 2 single-channel continuum spw" 00466 print "Do not want bandpass calibration" 00467 00468 vis = msfile 00469 00470 # set the name for the output gain caltable 00471 caltable = gtable 00472 00473 print "Output gain cal table will be "+gtable 00474 00475 # Gain calibrators are 1331+305 and 0137+331 (FIELD_ID 7 and 0) 00476 # We have 2 IFs (SPW 0,1) with one channel each 00477 00478 # selection is via the field and spw strings 00479 field = '1331+305,0137+331' 00480 spw = '' 00481 00482 # a-priori calibration application 00483 gaincurve = usegaincurve 00484 opacity = gainopacity 00485 00486 # scan-based G solutions for both amplitude and phase 00487 gaintype = 'G' 00488 calmode = 'ap' 00489 00490 # one solution per scan 00491 solint = 'inf' 00492 combine = '' 00493 00494 # do not apply parallactic angle correction (yet) 00495 parang = False 00496 00497 # reference antenna 00498 refant = calrefant 00499 00500 # minimum SNR 3 00501 minsnr = 3 00502 00503 gaincal() 00504 00505 # 00506 #===================================================================== 00507 # 00508 # Bootstrap flux scale 00509 # 00510 print '--Fluxscale--' 00511 default('fluxscale') 00512 00513 print "Use fluxscale to rescale gain table to make new one" 00514 00515 vis = msfile 00516 00517 # set the name for the output rescaled caltable 00518 fluxtable = ftable 00519 00520 print "Output scaled gain cal table is "+ftable 00521 00522 # point to our first gain cal table 00523 caltable = gtable 00524 00525 # we will be using 1331+305 (the source we did setjy on) as 00526 # our flux standard reference 00527 reference = '1331+305' 00528 00529 # we want to transfer the flux to our other gain cal source 0137+331 00530 # to bring its gain amplitues in line with the absolute scale 00531 transfer = '0137+331' 00532 00533 fluxscale() 00534 00535 # You should see in the logger something like: 00536 #Flux density for 0137+331 in SpW=0 is: 00537 # 5.42575 +/- 0.00285011 (SNR = 1903.7, nAnt= 27) 00538 #Flux density for 0137+331 in SpW=1 is: 00539 # 5.46569 +/- 0.00301326 (SNR = 1813.88, nAnt= 27) 00540 00541 # 00542 #--------------------------------------------------------------------- 00543 # Plot calibration 00544 # 00545 print '--PlotCal--' 00546 default('plotcal') 00547 00548 showgui = True 00549 00550 caltable = ftable 00551 multiplot = True 00552 yaxis = 'amp' 00553 00554 showgui = True 00555 00556 plotcal() 00557 00558 print "" 00559 print "-------------------------------------------------" 00560 print "Plotcal" 00561 print "Looking at amplitude in cal-table "+caltable 00562 00563 # Pause script if you are running in scriptmode 00564 if scriptmode: 00565 user_check=raw_input('Return to continue script\n') 00566 00567 # 00568 # Now go back and plot to file 00569 # 00570 showgui = False 00571 00572 yaxis = 'amp' 00573 00574 figfile = caltable + '.plotcal.amp.png' 00575 print "Plotting calibration to file "+figfile 00576 #saveinputs('plotcal',caltable.plotcal.amp.saved') 00577 plotcal() 00578 00579 yaxis = 'phase' 00580 00581 figfile = caltable + '.plotcal.phase.png' 00582 print "Plotting calibration to file "+figfile 00583 #saveinputs('plotcal',caltable.plotcal.phase.saved') 00584 plotcal() 00585 00586 # 00587 #===================================================================== 00588 # Polarization Calibration 00589 #===================================================================== 00590 # 00591 if (dopolcal): 00592 print '--Polcal (D)--' 00593 default('polcal') 00594 00595 print "Solve for polarization leakage on 0137+331" 00596 print "Pretend it has unknown polarization" 00597 00598 vis = msfile 00599 00600 # Start with the un-fluxscaled gain table 00601 gaintable = gtable 00602 00603 # use settings from gaincal 00604 gaincurve = usegaincurve 00605 opacity = gainopacity 00606 00607 # Output table 00608 caltable = ptable 00609 00610 # Use a 3C48 tracked through a range of PA 00611 field = '0137+331' 00612 spw = '' 00613 00614 # No need for further selection 00615 selectdata=False 00616 00617 # Polcal mode (D+QU = unknown pol for D) 00618 poltype = 'D+QU' 00619 00620 # One solution for entire dataset 00621 solint = 'inf' 00622 combine = 'scan' 00623 00624 # reference antenna 00625 refant = calrefant 00626 00627 # minimum SNR 3 00628 minsnr = 3 00629 00630 #saveinputs('polcal',calprefix+'.polcal.saved') 00631 polcal() 00632 00633 #===================================================================== 00634 # 00635 # List polcal solutions 00636 # 00637 print '--Listcal (PolD)--' 00638 00639 listfile = caltable + '.list' 00640 00641 print "Listing calibration to file "+listfile 00642 00643 listcal() 00644 00645 #===================================================================== 00646 # 00647 # Plot polcal solutions 00648 # 00649 print '--Plotcal (PolD)--' 00650 00651 iteration = '' 00652 showgui = False 00653 00654 xaxis = 'real' 00655 yaxis = 'imag' 00656 figfile = caltable + '.plotcal.reim.png' 00657 print "Plotting calibration to file "+figfile 00658 #saveinputs('plotcal',caltable+'.plotcal.reim.saved') 00659 plotcal() 00660 00661 xaxis = 'antenna' 00662 yaxis = 'amp' 00663 figfile = caltable + '.plotcal.antamp.png' 00664 print "Plotting calibration to file "+figfile 00665 #saveinputs('plotcal',caltable+'.plotcal.antamp.saved') 00666 plotcal() 00667 00668 xaxis = 'antenna' 00669 yaxis = 'phase' 00670 figfile = caltable + '.plotcal.antphase.png' 00671 print "Plotting calibration to file "+figfile 00672 #saveinputs('plotcal',caltable+'.plotcal.antphase.saved') 00673 plotcal() 00674 00675 xaxis = 'antenna' 00676 yaxis = 'snr' 00677 figfile = caltable + '.plotcal.antsnr.png' 00678 print "Plotting calibration to file "+figfile 00679 #saveinputs('plotcal',caltable+'.plotcal.antsnr.saved') 00680 plotcal() 00681 00682 #===================================================================== 00683 # Do Chi (X) pol angle calibration 00684 #===================================================================== 00685 # First set the model 00686 print '--Setjy--' 00687 default('setjy') 00688 00689 vis = msfile 00690 00691 print "Use setjy to set IQU fluxes of "+polxfield 00692 field = polxfield 00693 00694 scalebychan=False 00695 00696 for spw in usespwlist: 00697 fluxdensity = polxiquv[spw] 00698 00699 #saveinputs('setjy',calprefix+'.setjy.polspw.'+spw+'.saved') 00700 setjy() 00701 00702 # 00703 # Polarization (X-term) calibration 00704 # 00705 print '--PolCal (X)--' 00706 default('polcal') 00707 00708 print "Polarization R-L Phase Calibration (linear approx)" 00709 00710 vis = msfile 00711 00712 # Start with the G and D tables 00713 gaintable = [gtable,ptable] 00714 00715 # use settings from gaincal 00716 gaincurve = usegaincurve 00717 opacity = gainopacity 00718 00719 # Output table 00720 caltable = xtable 00721 00722 # previously set with setjy 00723 field = polxfield 00724 spw = '' 00725 00726 selectdata=False 00727 00728 # Solve for Chi 00729 poltype = 'X' 00730 solint = 'inf' 00731 combine = 'scan' 00732 00733 # reference antenna 00734 refant = calrefant 00735 00736 # minimum SNR 3 00737 minsnr = 3 00738 00739 #saveinputs('polcal',calprefix+'.polcal.X.saved') 00740 polcal() 00741 #============================================================================= 00742 00743 # # # # TESTS BEGIN HERE # # # # # # # # # # # # # # # # # # # # # # # # # # # 00744 00745 testNum = 0 00746 00747 ########################################################################## 00748 # TEST1 00749 # Using ngc4826 data. This data contains multiple spectral windows with 00750 # variable numbers of channels; only one correlation, YY. 00751 # 00752 00753 testNum += 1 00754 print "" 00755 print "* TEST " + str(testNum) + ": default test; using ngc4826 tutorial data and calibration" 00756 print """ 00757 - This data has multiple spws with variable numbers of channels. 00758 - This data has only one correlation, YY 00759 > Using default input values where possible 00760 """ 00761 00762 prefix = regressionDir+'/test'+str(testNum)+'/' 00763 msname = prefix+"ngc4826.tutorial.ms" 00764 caltableName = prefix+"ngc4826.tutorial.16apr98.gcal" 00765 outputFilename = prefix+'listcal.ngc4826.default.out' 00766 standardFileName = localData + 'listcal.default.out' 00767 00768 # Use existing data or load data from scratch? 00769 if (not lt.resetData([msname,caltableName], automate)): 00770 print "Using preexisting data." 00771 lt.removeOut(outputFilename) 00772 else: 00773 print "Building data from scratch." 00774 tstutl.maketestdir(prefix) # create test dir, overwrite preexisting 00775 load_ngc4826(prefix,msname,caltableName) # Build data from scratch 00776 00777 # Setup listcal input and run 00778 default(listcal) 00779 vis = msname # Name of input visibility file 00780 caltable = caltableName # Input calibration table to list 00781 field = '' # Field name or index; ''==>all 00782 antenna = '' # Antenna name or index; ''==>all; antenna='3' 00783 spw = '' # Spectral window and channel: ''==>all; 00784 # spw='5:0~10' 00785 listfile = outputFilename # Disk file to write output: ''==>to terminal 00786 pagerows = 50 # Rows per page 00787 async = False # If true the taskname must be started using 00788 go(listcal) 00789 00790 # Remove first line of listcal output (contains hard-coded path) 00791 lt.listcalFix(outputFilename) 00792 00793 compareFilename = prefix + 'compare' 00794 if (lt.runTests(outputFilename,standardFileName,'1.000',compareFilename)): 00795 print "Passed listcal output test" 00796 testPassed +=1 00797 else: 00798 print "FAILED listcal output test" 00799 testFailed +=1 00800 00801 ########################################################################## 00802 # TEST - NGC4826 00803 # Using ngc4826 data from above. 00804 # 00805 testNum += 1 00806 print "" 00807 print "* TEST " + str(testNum) + ": using ngc4826 tutorial data and calibration." 00808 print """ 00809 - Using same data as above 00810 > Using all non-default values where possible 00811 00812 """ 00813 00814 prefix = regressionDir+'/test'+str(testNum)+'/' 00815 tstutl.maketestdir(prefix) 00816 outputFilename = prefix+'listcal.ngc4826.nondefault.out' 00817 standardFileName = localData + 'listcal.nondefault.out' 00818 compareFilename = prefix + 'compare' 00819 00820 default(listcal) 00821 vis = msname # Name of input visibility file 00822 caltable = caltableName # Input calibration table to list 00823 field = '1310+323' # Field name or index; ''==>all 00824 antenna = '3~5,10' # Antenna name or index; ''==>all; antenna='3' 00825 spw = '0' # Spectral window and channel: ''==>all; spw='5:0~10' 00826 listfile = outputFilename # Disk file to write output: ''==>to terminal 00827 pagerows = 9 # Rows per page 00828 async = False # If true the taskname must be started using 00829 go(listcal) 00830 00831 # Remove first line of listcal output (contains hard-coded path) 00832 lt.listcalFix(outputFilename) 00833 00834 if (lt.runTests(outputFilename,standardFileName,'1.000',compareFilename)): 00835 print "Passed listcal output test" 00836 testPassed +=1 00837 else: 00838 print "FAILED listcal output test" 00839 testFailed +=1 00840 00841 ########################################################################## 00842 # Test complete, summarize. 00843 # 00844 00845 print "" 00846 print "* listcal regression test complete" 00847 print "SUMMARY:" 00848 print " number of tests PASSED: " + str(testPassed) 00849 print " number of tests FAILED: " + str(testFailed) 00850 print "" 00851 00852 if testFailed > 0: 00853 print '' 00854 print 'Regression FAILED' 00855 print '' 00856 else: 00857 print '' 00858 print 'Regression PASSED' 00859 print '' 00860 00861 # # If running in automated mode, remove all data files. 00862 # if (automate): 00863 # print "Removing all listcal_regression files..." 00864 # os.system('rm -f '+regressionDir) # remove all listcal output 00865 00866 print "END: listcal_regression.py" 00867 00868 if (testFailed > 0): 00869 regstate=False 00870 else: 00871 regstate=True 00872 00873 endTime = time.time()
1.8.0
