|
casa
$Rev:20696$
|
|
casa
$Rev:20696$
|
00001 ############################################### 00002 # # 00003 # Regression/Benchmarking Script for G192 # 00004 # # 00005 ############################################### 00006 00007 import time 00008 import os 00009 00010 os.system('rm -rf g192*.ms g192_a* g192.*.im') 00011 00012 startTime = time.time() 00013 startProc = time.clock() 00014 00015 print '--Import--' 00016 datapath=os.environ.get('CASAPATH').split()[0] +'/data/regression/ATST1/G192/' 00017 00018 default('importvla') 00019 importvla(archivefiles=[datapath + 'AS758_C030425.xp1',datapath+'AS758_C030425.xp2',datapath+'AS758_C030425.xp3',datapath+'AS758_C030426.xp4',datapath+'AS758_C030426.xp5'], 00020 vis='g192_a.ms',bandname='K',frequencytol=10000000.0) 00021 importtime = time.time() 00022 print '--Observation summary--' 00023 default('listobs') 00024 listobs(vis='g192_a.ms') 00025 #listtime = time.time() 00026 print '--Flag auto-correlations--' 00027 default('tflagdata') 00028 tflagdata(vis='g192_a.ms',mode='manual',autocorr=True) 00029 flagtime = time.time() 00030 print '--Setjy--' 00031 default('setjy') 00032 setjy(vis='g192_a.ms',field='4',scalebychan=False,standard='Perley-Taylor 99') #set flux density for 1331+305 (3C286) 00033 setjytime = time.time() 00034 print '--Gaincal--' 00035 #Select data for gain calibrators and drop outer channesl that may bias the solution 00036 #correct for specified zenith opacity and VLA gain curves (default) 00037 default('gaincal') 00038 gaincal(vis='g192_a.ms',caltable='g192_a.gcal', 00039 field='0,2,3,4',spw='0:3~117', gaintype='G', 00040 opacity=0.062,solint='inf',combine='',refant='VA05') 00041 gaintime = time.time() 00042 print '--Bandpass--' 00043 #Select bandpass calibrator. Arrange to solve for a single solution over the entire observation 00044 default('bandpass') 00045 bandpass(vis='g192_a.ms',caltable='g192_a.bcal', 00046 field='3', 00047 opacity=0.062,gaintable='g192_a.gcal',gainfield='3',interp='nearest', 00048 solint='inf',combine='scan', 00049 refant='VA05') 00050 bptime = time.time() 00051 print '--Fluxscale--' 00052 #Transfer the flux density scale from the flux calibrator to the gain calibrators 00053 #Solutions are written to a table (g192_a.fluxcal) on disk. 00054 default('fluxscale') 00055 fluxscale(vis='g192_a.ms',caltable='g192_a.gcal',fluxtable='g192_a.fluxcal', 00056 reference=['1331+305'],transfer=['0530+135','05309+13319']) 00057 fstime = time.time() 00058 00059 print '--Correct--' 00060 #Apply calibration solutions to the data 00061 default('applycal') 00062 applycal(vis='g192_a.ms', 00063 field='0,1,2', 00064 opacity=0.062, 00065 gaintable=['g192_a.fluxcal','g192_a.bcal'],gainfield='0,2'); 00066 correcttime = time.time() 00067 00068 print '--Split (Cal/src data)--' 00069 #split out the data of interest into a new visibility data set 00070 default('split') 00071 split(vis='g192_a.ms',outputvis='g192_cal.split.ms', 00072 # field=4,spw=0,nchan=100,start=9,step=1,datacolumn='CORRECTED_DATA') 00073 field='4',spw='0:9~108',datacolumn='corrected') 00074 splitcaltime = time.time() 00075 default('split') 00076 split(vis='g192_a.ms',outputvis='g192_src.split.ms', 00077 # field=1,spw=0,nchan=100,start=9,step=1,datacolumn='CORRECTED_DATA') 00078 field='1',spw='0:9~108',datacolumn='corrected') 00079 splitsrctime = time.time() 00080 00081 print '--Flag bad time range--' 00082 #flag data in the specified time range for the source and spw 00083 tflagdata(vis="g192_src.split.ms", field="0", spw="0", 00084 timerange="2003/04/26/02:45:00.0~2003/04/26/02:46:30.0") 00085 flagsrctime=time.time() 00086 00087 print '--Clean src line--' 00088 #image the source 00089 default('clean') 00090 clean(vis='g192_src.split.ms',imagename='g192_a2',mode='channel', 00091 psfmode='hogbom',niter=5000,gain=0.1,threshold=15.,mask='', 00092 nchan=40,start=34,width=1,field='0',spw='0', 00093 imsize=[512,512],cell=[.5,.5],weighting='natural') 00094 cleantime = time.time() 00095 00096 print '--Write FITS--' 00097 #export the data to fits 00098 ia.open(infile='g192_a2.image') 00099 ia.tofits(outfile='g192_a2.fits') 00100 ia.close() 00101 writefitstime=time.time() 00102 00103 print '--Contsub (image plane)--' 00104 #do image plane continuum subtraction; channels specified will be used to make a continuum image: g192_cont.im 00105 ia.open(infile='g192_a2.image') 00106 myim=ia.continuumsub(outline='g192.line.im',outcont='g192.cont.im', 00107 channels=range(0,1),fitorder=0) 00108 x=myim.statistics(list=False) 00109 thistest_con=x['rms'][0] 00110 ia.close() 00111 myim.close() 00112 contsubtime=time.time() 00113 00114 endProc = time.clock() 00115 endTime = time.time() 00116 00117 # Regression 00118 00119 test_name_cal = 'G192--Calibrater maximum amplitude test' 00120 test_name_src = 'G192--Source maximum amplitude test' 00121 #test_name_con = 'G192--Continuum subtraction rms test' 00122 test_name_immax = 'G192--Image maximum test' 00123 test_name_imrms = 'G192--Image rms test' 00124 00125 ms.open('g192_cal.split.ms') 00126 thistest_cal=max(ms.range(["amplitude"]).get('amplitude')) 00127 ms.close() 00128 ms.open('g192_src.split.ms') 00129 thistest_src=max(ms.range(["amplitude"]).get('amplitude')) 00130 ms.close() 00131 ia.open('g192_a2.image') 00132 # ia.statistics returns dictionary with 'return','statsout' 00133 # get the second value in the dictionary (statsout) 00134 statistics=ia.statistics() 00135 # note thistest_immax will be a list with one value 00136 thistest_immax=statistics['max'] 00137 # note thistest_imrms will be a list with one value 00138 thistest_imrms=statistics['rms'] 00139 00140 calmax=2.7573 00141 srcmax= 25.116 00142 contsubrms= 0.00283678 00143 immax=0.026332 00144 imrms= 0.0020570 00145 #data set size = 634.9 MB - 5 VLA archive (xp) files 00146 00147 diff_cal=abs((calmax-thistest_cal)/calmax) 00148 diff_src=abs((srcmax-thistest_src)/srcmax) 00149 diff_con=abs((contsubrms-thistest_con)/contsubrms) 00150 diff_immax=abs((immax-thistest_immax[0])/immax) 00151 diff_imrms=abs((imrms-thistest_imrms[0])/imrms) 00152 00153 import datetime 00154 datestring=datetime.datetime.isoformat(datetime.datetime.today()) 00155 outfile='g192.'+datestring+'.log' 00156 logfile=open(outfile,'w') 00157 00158 print >>logfile,'' 00159 print >>logfile,'********** Data Summary *********' 00160 print >>logfile,'* Observer: unavailable Project: AS758 *' 00161 print >>logfile,'* Observation: VLA(27 antennas) *' 00162 print >>logfile,'* Telescope Observation Date Observer Project *' 00163 print >>logfile,'* VLA [ 4.55803e+09, 4.55803e+09]unavailable AS758 *' 00164 print >>logfile,'* VLA [ 4.55803e+09, 4.55803e+09]unavailable AS758 *' 00165 print >>logfile,'* VLA [ 4.55803e+09, 4.55803e+09]unavailable AS758 *' 00166 print >>logfile,'* VLA [ 4.55803e+09, 4.55804e+09]unavailable AS758 *' 00167 print >>logfile,'* VLA [ 4.55804e+09, 4.55804e+09]unavailable AS758 *' 00168 print >>logfile,'* Data records: 1200015 Total integration time = 19347.5 seconds *' 00169 print >>logfile,'* Observed from 25-Apr-2003/22:03:38 to 26-Apr-2003/03:26:05 *' 00170 print >>logfile,'* Fields: 6 *' 00171 print >>logfile,'* ID Name Right Ascension Declination Epoch *' 00172 print >>logfile,'* 0 0530+135 05:30:56.42 +13.31.55.15 J2000 (gaincal) *' 00173 print >>logfile,'* 1 05582+16320 05:58:13.53 +16.31.58.29 J2000 (target) *' 00174 print >>logfile,'* 2 05309+13319 05:30:56.42 +13.31.55.15 J2000 (gaincal) *' 00175 print >>logfile,'* 3 0319+415 03:19:48.16 +41.30.42.10 J2000 (bandpass) *' 00176 print >>logfile,'* 4 1331+305 13:31:08.29 +30.30.32.96 J2000 (fluxcal) *' 00177 print >>logfile,'* 6 KTIP 21:20:00.00 +60.00.00.00 J2000 *' 00178 print >>logfile,'* Data descriptions: 3 (3 spectral windows and 2 polarization setups) *' 00179 print >>logfile,'* ID #Chans Frame Ch1(MHz) Resoln(kHz) TotBW(kHz) Ref(MHz) Corrs *' 00180 print >>logfile,'* 0 127 LSRK 23692.5072 24.4170056 3100.95971 23694.045 RR *' 00181 print >>logfile,'* 2 1 TOPO 22485.1 50000 50000 22485.1 RR RL LR LL ' 00182 print >>logfile,'* 3 1 TOPO 22435.1 50000 50000 22435.1 RR RL LR LL ' 00183 print >>logfile,'*********************************' 00184 print >>logfile,'' 00185 print >>logfile,'********** Regression ***********' 00186 print >>logfile,'* *' 00187 00188 regstate=True 00189 if (diff_cal < 0.05): 00190 print >>logfile,'* Passed cal max amplitude test *' 00191 else: 00192 regstate=False 00193 print >>logfile,'* Failed cal max amplitude test *' 00194 print >>logfile,' Cal max amp '+str(thistest_cal)+' ('+str(calmax)+')' 00195 00196 if (diff_src < 0.05): 00197 print >>logfile,'* Passed src max amplitude test *' 00198 else: 00199 regstate=False 00200 print >>logfile,'* Failed src max amplitude test *' 00201 print >>logfile,' Src max amp '+str(thistest_src)+' ('+str(srcmax)+')' 00202 00203 if (diff_con < 0.05): 00204 print >>logfile,'* Passed contsub rms test *' 00205 else: 00206 regstate=False 00207 print >>logfile,'* Failed contsub rms test *' 00208 print >>logfile,' Contsub rms '+str(thistest_con)+' ('+str(contsubrms)+')' 00209 00210 if (diff_immax < 0.05): 00211 print >>logfile,'* Passed image max test *' 00212 else: 00213 regstate=False 00214 print >>logfile,'* Failed image max test *' 00215 print >>logfile,' Image max '+str(thistest_immax)+' ('+str(immax)+')' 00216 00217 if (diff_imrms < 0.05): 00218 print >>logfile,'* Passed image rms test *' 00219 else: 00220 regstate=False 00221 print >>logfile,'* Failed image rms test *' 00222 print >>logfile,' Image rms '+str(thistest_imrms)+' ('+str(imrms)+')' 00223 00224 00225 if (regstate): 00226 print >>logfile,'---' 00227 print >>logfile,'Passed Regression test for G192' 00228 print >>logfile,'---' 00229 else: 00230 print >>logfile,'----FAILED Regression test for G192' 00231 print >>logfile,'*********************************' 00232 00233 print >>logfile,'' 00234 print >>logfile,'' 00235 print >>logfile,'********* Benchmarking *****************' 00236 print >>logfile,'* *' 00237 print >>logfile,'Total wall clock time was: '+str(endTime - startTime) 00238 print >>logfile,'Total CPU time was: '+str(endProc - startProc) 00239 print >>logfile,'Processing rate MB/s was: '+str(634.9/(endTime - startTime)) 00240 print >>logfile,'* Breakdown: *' 00241 print >>logfile,'* import time was: '+str(importtime-startTime) 00242 print >>logfile,'* tflagdata time was: '+str(flagtime-importtime) 00243 print >>logfile,'* setjy time was: '+str(setjytime-flagtime) 00244 print >>logfile,'* gaincal time was: '+str(gaintime-setjytime) 00245 print >>logfile,'* bandpass time was: '+str(bptime-gaintime) 00246 print >>logfile,'* fluxscale time was: '+str(fstime-bptime) 00247 print >>logfile,'* applycal time was: '+str(correcttime-fstime) 00248 print >>logfile,'* split-cal time was: '+str(splitcaltime-correcttime) 00249 print >>logfile,'* split-src time was: '+str(splitsrctime-splitcaltime) 00250 print >>logfile,'* flag-src time was: '+str(flagsrctime-splitsrctime) 00251 print >>logfile,'* clean-src time was: '+str(cleantime-flagsrctime) 00252 #print '* exportfits tiem was: ',exportfitstime-cleantime 00253 #print '* contsub time was: ',contsubtime-exportfitstime 00254 print >>logfile,'*****************************************' 00255 print >>logfile,'basho (test cpu) time was: 1500 seconds' 00256 00257 logfile.close()
1.8.0
