|
casa
$Rev:20696$
|
|
casa
$Rev:20696$
|
00001 ############################### 00002 # 00003 # IRC+10216 HC3N Reduction Script 00004 # using new sd tasks 00005 # Nod data 00006 # 00007 # tasks used 00008 # sdlist 00009 # sdreduce (new) 00010 # sdplot 00011 # sdstat 00012 # sdsave 00013 ############################### 00014 import time 00015 import os 00016 00017 os.system('rm -rf IRC+10216_rawACSmod IRC+10216_rawACSmod_cal IRC+10216_rawACSmod_cal_sm IRC+10216_rawACSmod_cal_sm_bs irc_hc3n_reducedSCAN0_CYCLE0_BEAM0_IF0.txt irc_hc3n_reduced.eps irc_cs_fit.txt') 00018 00019 #enable/disable plotting 00020 doplot = False 00021 00022 casapath=os.environ['CASAPATH'].split()[0] 00023 datapath=casapath+'/data/regression/ATST5/IRC+10216/IRC+10216_rawACSmod' 00024 copystring='cp -r '+datapath+' .' 00025 os.system(copystring) 00026 00027 startTime=time.time() 00028 startProc=time.clock() 00029 00030 # Project: AGBT06A_018_01 00031 #Observation: GBT(1 antennas) 00032 # 00033 #Tue Jan 30 20:24:29 2007 NORMAL ms::summary: 00034 #Data records: 4992 Total integration time = 4244.15 seconds 00035 # Observed from 10:31:01 to 11:41:45 00036 # 00037 #Tue Jan 30 20:24:29 2007 NORMAL ms::summary: 00038 #Fields: 3 00039 # ID Name Right Ascension Declination Epoch 00040 # 0 IRC+10216 09:47:57.38 +13.16.43.70 J2000 00041 # 1 IRC+10216 09:47:57.38 +13.16.43.70 J2000 00042 # 2 IRC+10216 09:47:57.38 +13.16.43.70 J2000 00043 # 00044 #Tue Jan 30 20:24:29 2007 NORMAL ms::summary: 00045 #Spectral Windows: (13 unique spectral windows and 1 unique polarization setups) 00046 # SpwID #Chans Frame Ch1(MHz) Resoln(kHz) TotBW(kHz) Ref(MHz) Corrs 00047 # 0 4096 LSRK 48263.5858 12.2067945 49999.0303 48288.5853 RR LL H2CO 00048 # 1 4096 LSRK 45899.1486 12.2067945 49999.0303 45924.1481 RR LL H213CO 00049 # 2 4096 LSRK 43857.1942 12.2067945 49999.0303 43882.1937 RR LL 13CS 00050 # 3 4096 LSRK 48970.0031 12.2067945 49999.0303 48995.0026 RR LL CS 00051 # 7 4096 LSRK 46226.6452 12.2067949 49999.0319 46251.6447 RR LL 00052 # 17 4096 LSRK 45469.0638 12.2067955 49999.0343 45494.0633 RR LL HC3N 00053 # 18 4096 LSRK 45280.5214 12.2067955 49999.0343 45305.5209 RR LL HCC13CN 00054 # 19 4096 LSRK 44054.8331 12.2067955 49999.0343 44079.8326 RR LL CH3OH 00055 # 20 4096 LSRK 44146.1063 12.2067955 49999.0343 44171.1058 RR LL HCCC15N 00056 # 27 4096 LSRK 43941.6784 12.206796 49999.0365 43966.6779 RR LL HN13CO 00057 # 28 4096 LSRK 42625.1528 12.206796 49999.0365 42650.1523 RR LL H15NCO 00058 # 29 4096 LSRK 41574.856 12.206796 49999.0365 41599.8556 RR LL HNC18O 00059 # 30 4096 LSRK 43402.4488 12.206796 49999.0365 43427.4483 RR LL SiO 00060 00061 # Scans: 230-231,233-236,253-256 Setup 1 H2CO et al 00062 # Scans: 237-240,249-252 Setup 2 HC3N et al 00063 # Scans: 241-248 Setup 3 SiO et al 00064 00065 00066 asap_init() #load ASAP module 00067 00068 #scan numbers (zero-based) as compared to GBTIDL 00069 #changes made to get to IRC+10216_rawACSmod 00070 # -- merge spectral windows with tolerance 00071 00072 # summary 00073 #default(sdlist) 00074 #infile='IRC+10216_rawACSmod' 00075 #sdlist() 00076 00077 if doplot: 00078 localplotlevel = 1 00079 else: 00080 localplotlevel = 0 00081 00082 # calibartion,averaging, smoothing, and baseline removal 00083 # calibrate nod scans for CS line (IF=3) 00084 default(sdreduce) 00085 infile='IRC+10216_rawACSmod' 00086 fluxunit='K' 00087 calmode='nod' 00088 scanlist=[236,237,238,239,248,249,250,251] 00089 iflist=[17] 00090 average=True 00091 scanaverage=False 00092 timeaverage=True # average in time 00093 tweight='tintsys' # weighted by iteg time and Tsys for time averaging 00094 polaverage=True # average polarization 00095 pweight='tsys' # weighted by Tsys for polarization averaging 00096 tau=0.09 # do opacity correction 00097 # do boxcar smoothing with channel width=5 00098 kernel='boxcar' 00099 kwidth=5 00100 # output 00101 #fit and remove baselines 00102 maskmode='auto' 00103 thresh=5 00104 avg_limit=4 00105 blfunc='poly' 00106 order=2 00107 overwrite=True 00108 plotlevel=localplotlevel 00109 sdreduce() 00110 localoutfile=infile+'_cal' 00111 00112 #plotting the reslut 00113 #plot the spectrum and save to a postscript file 00114 if doplot: 00115 default(sdplot) 00116 infile=localoutfile 00117 specunit='GHz' 00118 outfile='irc_hc3n_reduced.eps' 00119 #sd.plotter.set_histogram(hist=True) # draw spectrum using histogram # histogram 00120 #sd.plotter.axhline(color='r',linewidth=2) # zline # zline 00121 sdplot() 00122 else: 00123 print "Plotting the result is skipped." 00124 00125 # statistics 00126 default(sdstat) 00127 # select line free regions to get rms 00128 infile=localoutfile 00129 masklist=[200,1500] 00130 xstat=sdstat() 00131 rms=xstat['rms'] 00132 #rms= 00133 # 00134 # select the line region 00135 masklist=[1800,2400] 00136 xstat=sdstat() 00137 xstat 00138 max=xstat['max'] 00139 sum=xstat['sum'] 00140 median=xstat['median'] 00141 mean=xstat['mean'] 00142 00143 # Save the spectrum 00144 # in different formats 00145 default(sdsave) 00146 infile=localoutfile 00147 outfile='irc_hc3n_reduced' 00148 outform='ASCII' 00149 overwrite=True 00150 sdsave() 00151 #outfile='irc_hc3n_reduced.ms' 00152 #outform='MS2' 00153 #sdsave() 00154 00155 # 00156 endProc = time.clock() 00157 endTime = time.time() 00158 00159 # --- end of irc cs script 00160 #irc_max=1.827 00161 #irc_rms=0.022 00162 #irc_sum=474.123 00163 # Regression values of CASA 2.3(#6654)+ASAP 2.2.0(#1448) 00164 # on 64bit REL5.2 (2008/12/01) 00165 irc_max=1.827 00166 irc_rms=0.02213 00167 irc_sum=474.1 00168 diff_max = abs((irc_max-max)/irc_max) 00169 diff_rms = abs((irc_rms-rms)/irc_rms) 00170 diff_sum = abs((irc_sum-sum)/irc_sum) 00171 00172 import datetime 00173 datestring=datetime.datetime.isoformat(datetime.datetime.today()) 00174 outfile='irc.hc3n.task'+datestring+'.log' 00175 logfile=open(outfile,'w') 00176 00177 print >>logfile,'' 00178 print >>logfile,'********** Regression ***********' 00179 print >>logfile,'* *' 00180 if (diff_max < 0.05): print >>logfile,'* Passed spectrum max test ' 00181 print >>logfile,'* Spectrum max '+str(max) 00182 if (diff_rms < 0.05): print >>logfile,'* Passed spectrum rms test ' 00183 print >>logfile,'* Spectrum rms '+str(rms) 00184 if (diff_sum < 0.05): print >>logfile,'* Passed spectrum (line) sum test' 00185 print >>logfile,'* Line integral '+str(sum) 00186 if ((diff_max<0.05) & (diff_rms<0.05) & (diff_sum<0.05)): 00187 regstate=True 00188 print '' 00189 print 'Regression PASSED' 00190 print '' 00191 print >>logfile,'---' 00192 print >>logfile,'Passed Regression test for IRC-HC3N' 00193 print >>logfile,'---' 00194 else: 00195 regstate=False 00196 print '' 00197 print 'Regression FAILED' 00198 print '' 00199 print >>logfile,'----FAILED Regression test for IRC-HC3N' 00200 print >>logfile,'*********************************' 00201 00202 print >>logfile,'' 00203 print >>logfile,'' 00204 print >>logfile,'********* Benchmarking *****************' 00205 print >>logfile,'* *' 00206 print >>logfile,'Total wall clock time was: '+str(endTime - startTime) 00207 print >>logfile,'Total CPU time was: '+str(endProc - startProc) 00208 #print 'Processing rate MB/s was: ', 35.1/(endTime - startTime) 00209 00210 logfile.close()
1.8.0
