|
casa
$Rev:20696$
|
|
casa
$Rev:20696$
|
00001 from taskinit import * 00002 import shutil 00003 from parallel.parallel_task_helper import ParallelTaskHelper 00004 00005 00006 def fixplanets(vis, field, fixuvw=False, direction='', refant=0, reftime='first'): 00007 """ 00008 Fix FIELD, SOURCE, and UVW for given fields based on given direction or pointing 00009 table information 00010 00011 This task's main purpose is to correct observations which were performed 00012 with correct pointing and correlation but for which incorrect direction 00013 information was entered in the FIELD and SOURCE table of the MS. 00014 If you actually want to change the phase center of the visibilties in an MS, 00015 you should use task fixvis. 00016 00017 Input Parameters 00018 vis -- Name of the input visibility set 00019 00020 field -- field selection string 00021 00022 fixuvw -- recalc uvw coordinates? (default: False) 00023 00024 direction -- if set, don't use pointing table but set direction to this value 00025 example: 'J2000 19h30m00 -40d00m00', default= '' (use pointing table) 00026 00027 refant -- if using pointing table information, use it from this antenna 00028 default: 0 (antenna id 0) 00029 examples: 'DV06' (antenna with name DV06) 00030 3 (antenna id 3) 00031 00032 reftime -- if using pointing table information, use it from this timestamp 00033 default: 'first' 00034 examples: 'median' will use the median timestamp for the given field 00035 using only the unflagged maintable rows 00036 '2012/07/11/08:41:32' will use the given timestamp (must be 00037 within the observaton time) 00038 00039 Examples: 00040 00041 fixplanets('uid___A002_X1c6e54_X223.ms', 'Titan', True) 00042 will look up the pointing direction from antenna 0 for field 'Titan' in 00043 the POINTING table based on the first timestamp in the main table rows for 00044 this field, enter this direction in the FIELD and SOURCE tables, and then 00045 recalculate the UVW coordinates for this field. 00046 00047 fixplanets('uid___A002_X1c6e54_X223.ms', 'Titan', False, 'J2000 12h30m15 -02d12m00') 00048 will set the directions for field 'Titan' in the FIELD and SOURCE table to the 00049 given direction and not recalculate the UVW coordinates. 00050 (This can be useful for several purposes, among them preparing a concatenation 00051 of datasets. Only fields with the same direction will be recognised as identical. 00052 fixplanets can then be run again after the concatenation using parameters as in 00053 the first example above.) 00054 00055 """ 00056 00057 casalog.origin('fixplanets') 00058 00059 mst = mstool() 00060 tbt = tbtool() 00061 imt = None 00062 00063 try: 00064 fields = mst.msseltoindex(vis=vis,field=field)['field'] 00065 numfields = 0 00066 if(len(fields) == 0): 00067 casalog.post( "Field selection returned zero results.", 'WARN') 00068 return True 00069 00070 tbt.open(vis+"/FIELD") 00071 oldrefcol = [] 00072 if('PhaseDir_Ref' in tbt.colnames()): 00073 oldrefcol = tbt.getcol('PhaseDir_Ref') 00074 tbt.close() 00075 00076 fixplanetstemp = 'fixplanetstemp-'+os.path.basename(vis) # temp file name 00077 fixplanetstemp2 = 'fixplanetstemp2-'+os.path.basename(vis) # temp file name 00078 00079 for fld in fields: 00080 thenewra_rad = 0. 00081 thenewdec_rad = 0. 00082 thenewref = -1 00083 thenewrefstr = '' 00084 if(direction==''): # use information from the pointing table 00085 # find median timestamp for this field in the main table 00086 shutil.rmtree(fixplanetstemp, ignore_errors=True) 00087 thetime = 0 00088 if(reftime.lower()=='median'): 00089 tbt.open(vis) 00090 tttb = tbt.query('FIELD_ID=='+str(fld)+' AND FLAG_ROW==False', 00091 name=fixplanetstemp, columns='TIME') 00092 tttb.close() 00093 tttb = None 00094 tbt.close() 00095 tbt.open(fixplanetstemp) 00096 if(tbt.nrows()>0): 00097 thetime = tbt.getcell('TIME',tbt.nrows()/2) 00098 casalog.post( "MEDIAN TIME "+str(thetime), 'NORMAL') 00099 tbt.close() 00100 else: 00101 casalog.post( "No pointing table rows for field "+field, 'NORMAL') 00102 tbt.close() 00103 shutil.rmtree(fixplanetstemp, ignore_errors=True) 00104 return True 00105 elif(reftime.lower()=='first'): 00106 tbt.open(vis) 00107 tttb = tbt.query('FIELD_ID=='+str(fld), name=fixplanetstemp, columns='TIME') 00108 tttb.close() 00109 tttb = None 00110 tbt.close() 00111 tbt.open(fixplanetstemp) 00112 if(tbt.nrows()>0): 00113 thetime = tbt.getcell('TIME',0) 00114 casalog.post( "FIRST TIME "+str(thetime), 'NORMAL') 00115 tbt.close() 00116 else: 00117 casalog.post( "No pointing table rows for field "+field, 'NORMAL') 00118 tbt.close() 00119 shutil.rmtree(fixplanetstemp, ignore_errors=True) 00120 return True 00121 else: 00122 try: 00123 myqa = qa.quantity(reftime) 00124 thetime = qa.convert(myqa,'s')['value'] 00125 except Exception, instance: 00126 raise TypeError, "reftime parameter is not a valid date (e.g. YYYY/MM/DD/hh:mm:ss)" %reftime 00127 tbt.open(vis) 00128 tttb = tbt.query('FIELD_ID=='+str(fld), name=fixplanetstemp, columns='TIME') 00129 tttb.close() 00130 tttb = None 00131 tbt.close() 00132 tbt.open(fixplanetstemp) 00133 if(tbt.nrows()>0): 00134 thefirsttime = tbt.getcell('TIME',0) 00135 thelasttime = tbt.getcell('TIME',tbt.nrows()-1) 00136 tbt.close() 00137 else: 00138 casalog.post( "No pointing table rows for field "+field, 'NORMAL') 00139 tbt.close() 00140 shutil.rmtree(fixplanetstemp, ignore_errors=True) 00141 return True 00142 shutil.rmtree(fixplanetstemp, ignore_errors=True) 00143 if (thefirsttime<=thetime and thetime<=thelasttime): 00144 casalog.post( "GIVEN TIME "+reftime+" == "+str(thetime), 'NORMAL') 00145 else: 00146 casalog.post( "GIVEN TIME "+reftime+" == "+str(thetime)+" is not within the observation time (" 00147 +str(thefirsttime)+'s'+" - " 00148 +str(thelasttime)+'s'+")", 'SEVERE') 00149 raise TypeError 00150 shutil.rmtree(fixplanetstemp, ignore_errors=True) 00151 00152 # determine reference antenna 00153 antids = mst.msseltoindex(vis=vis,baseline=refant)['antenna1'] 00154 antid = -1 00155 if(len(antids) == 0): 00156 casalog.post( "Antenna selection returned zero results.", 'WARN') 00157 return False 00158 antid = int(antids[0]) 00159 casalog.post('Using antenna id '+str(antid)+' as reference antenna.', 'NORMAL') 00160 tbt.open(vis+'/ANTENNA') 00161 flgcol = tbt.getcol('FLAG_ROW') 00162 tbt.close() 00163 if(flgcol[antid]==True): 00164 casalog.post('Antenna id '+str(antid)+' is flagged. Please select a different one.', 'SEVERE') 00165 return False 00166 00167 # get direction for the timestamp from pointing table 00168 shutil.rmtree(fixplanetstemp2, ignore_errors=True) 00169 tbt.open(vis+'/POINTING') 00170 ttb = tbt.query('TRACKING==True AND NEARABS(TIME,'+str(thetime)+',INTERVAL/2.) AND ANTENNA_ID==' 00171 +str(antid), 00172 name=fixplanetstemp2) 00173 nr = ttb.nrows() 00174 ttb.close() 00175 ttb = None 00176 if(nr==0): 00177 shutil.rmtree(fixplanetstemp2, ignore_errors=True) 00178 ttb2 = tbt.query('TRACKING==True AND NEARABS(TIME,'+str(thetime)+',3.) AND ANTENNA_ID==' 00179 +str(antid), # search within 3 seconds 00180 name=fixplanetstemp2) 00181 nr = ttb2.nrows() 00182 ttb2.close() 00183 ttb2 = None 00184 if(nr==0): 00185 casalog.post( "Cannot find any POINTING table rows for antenna "+str(antid) 00186 +" with TRACKING==True within 3 seconds of TIME "+str(thetime), 'NORMAL') 00187 casalog.post( "Will try without requiring TRACKING==True ...", 'NORMAL') 00188 shutil.rmtree(fixplanetstemp2, ignore_errors=True) 00189 ttb3 = tbt.query('NEARABS(TIME,'+str(thetime)+',INTERVAL/2.) AND ANTENNA_ID==' 00190 +str(antid), 00191 name=fixplanetstemp2) 00192 nr = ttb3.nrows() 00193 ttb3.close() 00194 ttb3 = None 00195 if(nr==0): 00196 shutil.rmtree(fixplanetstemp2, ignore_errors=True) 00197 ttb4 = tbt.query('NEARABS(TIME,'+str(thetime)+',3.) AND ANTENNA_ID==' 00198 +str(antid), # search within 3 seconds 00199 name=fixplanetstemp2) 00200 nr = ttb4.nrows() 00201 ttb4.close() 00202 ttb4 = None 00203 if(nr==0): 00204 tbt.close() 00205 casalog.post( "Cannot find any POINTING table rows for antenna "+str(antid) 00206 +" within 3 seconds of TIME "+str(thetime), 'SEVERE') 00207 return False # give up 00208 tbt.close() 00209 tbt.open(fixplanetstemp2) 00210 thedir = tbt.getcell('DIRECTION',0) 00211 tbt.close() 00212 casalog.post( ' field id '+str(fld)+ ' AZ EL '+str(thedir[0])+" "+str(thedir[1]), 'NORMAL') 00213 thedirme = me.direction(rf='AZELGEO',v0=qa.quantity(thedir[0][0], 'rad'), 00214 v1=qa.quantity(thedir[1][0],'rad')) 00215 # convert AZEL to J2000 coordinates 00216 me.doframe(me.epoch(rf='UTC', v0=qa.quantity(thetime,'s'))) 00217 00218 tbt.open(vis+'/ANTENNA') 00219 thepos = tbt.getcell('POSITION',antid) 00220 theposref = tbt.getcolkeyword('POSITION', 'MEASINFO')['Ref'] 00221 theposunits = tbt.getcolkeyword('POSITION', 'QuantumUnits') 00222 tbt.close() 00223 casalog.post( "Ref. antenna position is "+str(thepos)+' ('+theposunits[0]+', '+theposunits[1]+', ' 00224 +theposunits[2]+')('+theposref+')', 'NORMAL') 00225 me.doframe(me.position(theposref, 00226 qa.quantity(thepos[0],theposunits[0]), 00227 qa.quantity(thepos[1],theposunits[1]), 00228 qa.quantity(thepos[2],theposunits[2])) 00229 ) 00230 thedirmemod = me.measure(v=thedirme, rf='J2000') 00231 #print thedirmemod 00232 thenewra_rad = thedirmemod['m0']['value'] 00233 thenewdec_rad = thedirmemod['m1']['value'] 00234 me.done() 00235 shutil.rmtree(fixplanetstemp2, ignore_errors=True) 00236 else: # direction is not an empty string, use this instead of the pointing table information 00237 if(type(direction)==str): 00238 dirstr = direction.split(' ') 00239 codes = [] 00240 if len(dirstr)==2: 00241 dirstr = ['J2000', dirstr[0], dirstr[1]] 00242 casalog.post('Assuming reference frame J2000 for parameter \'direction\'', 'NORMAL') 00243 elif not len(dirstr)==3: 00244 casalog.post('Incorrect format in parameter \'direction\'', 'SEVERE') 00245 return False 00246 try: 00247 thedir = me.direction(dirstr[0], dirstr[1], dirstr[2]) 00248 thenewra_rad = thedir['m0']['value'] 00249 thenewdec_rad = thedir['m1']['value'] 00250 except Exception, instance: 00251 casalog.post("*** Error \'%s\' when interpreting parameter \'direction\': " % (instance), 'SEVERE') 00252 return False 00253 try: 00254 tbt.open(vis+"/FIELD") 00255 thenewrefindex = tbt.getcolkeyword('PHASE_DIR', 'MEASINFO')['TabRefTypes'].tolist().index(dirstr[0]) 00256 thenewref = tbt.getcolkeyword('PHASE_DIR', 'MEASINFO')['TabRefCodes'][thenewrefindex] 00257 thenewrefstr = dirstr[0] 00258 tbt.close() 00259 except Exception, instance: 00260 casalog.post("PHASE_DIR is not a variable reference column. Will leave reference as is.", 'WARN') 00261 thenewref = -1 00262 if(0<thenewref and thenewref<32): 00263 casalog.post("*** Error when interpreting parameter \'direction\':\n presently only J2000 and solar system objects are supported.", 00264 'SEVERE') 00265 return False 00266 #endif 00267 00268 # modify FIELD table 00269 tbt.open(vis+'/FIELD', nomodify=False) 00270 numfields = tbt.nrows() 00271 theolddir = tbt.getcell('PHASE_DIR',fld) 00272 planetname = tbt.getcell('NAME',fld) 00273 casalog.post( 'object: '+planetname, 'NORMAL') 00274 casalog.post( 'old RA, DEC (rad) '+str(theolddir[0])+" " +str(theolddir[1]), 'NORMAL') 00275 casalog.post( 'new RA, DEC (rad) '+str(thenewra_rad)+" "+ str(thenewdec_rad), 'NORMAL') 00276 00277 pcol = tbt.getcol('PHASE_DIR') 00278 pcol[0][0][fld] = thenewra_rad 00279 pcol[1][0][fld] = thenewdec_rad 00280 tbt.putcol('PHASE_DIR',pcol) 00281 00282 pcol = tbt.getcol('DELAY_DIR') 00283 pcol[0][0][fld] = thenewra_rad 00284 pcol[1][0][fld] = thenewdec_rad 00285 tbt.putcol('DELAY_DIR',pcol) 00286 00287 pcol = tbt.getcol('REFERENCE_DIR') 00288 pcol[0][0][fld] = thenewra_rad 00289 pcol[1][0][fld] = thenewdec_rad 00290 tbt.putcol('REFERENCE_DIR',pcol) 00291 casalog.post("FIELD table PHASE_DIR, DELAY_DIR, and REFERENCE_DIR columns changed for field "+str(fld)+".", 'NORMAL') 00292 00293 if(thenewref!=-1): 00294 # modify reference of the direction columns permanently 00295 theoldref = -1 00296 theoldref2 = -1 00297 theoldref3 = -1 00298 try: 00299 pcol = tbt.getcol('PhaseDir_Ref') 00300 theoldref = pcol[fld] 00301 pcol[fld] = thenewref 00302 oldrefcol = pcol 00303 00304 pcol2 = tbt.getcol('DelayDir_Ref') 00305 theoldref2 = pcol2[fld] 00306 pcol2[fld] = thenewref 00307 00308 pcol3 = tbt.getcol('RefDir_Ref') 00309 theoldref3 = pcol3[fld] 00310 pcol3[fld] = thenewref 00311 00312 if not (theoldref==theoldref2 and theoldref==theoldref3): 00313 casalog.post("The three FIELD table direction reference frame entries for field "+str(fld) 00314 +" are not identical in the input data: " 00315 +str(theoldref)+", "+str(theoldref2)+", "+str(theoldref3) 00316 +". Will try to continue ...", 'WARN') 00317 theoldref=-1 00318 else: 00319 oldrefstrindex = tbt.getcolkeyword('PHASE_DIR', 'MEASINFO')['TabRefCodes'].tolist().index(theoldref) 00320 oldrefstr = tbt.getcolkeyword('PHASE_DIR', 'MEASINFO')['TabRefTypes'][oldrefstrindex] 00321 casalog.post("Original FIELD table direction reference frame entries for field "+str(fld) 00322 +" are all "+str(theoldref)+" ("+oldrefstr+")", 'NORMAL') 00323 tbt.putcol('PhaseDir_Ref', pcol) 00324 tbt.putcol('DelayDir_Ref', pcol2) 00325 tbt.putcol('RefDir_Ref', pcol3) 00326 00327 except Exception, instance: 00328 casalog.post("*** Error \'%s\' when writing reference frames in FIELD table: " % (instance), 00329 'SEVERE') 00330 return False 00331 if(theoldref != thenewref): 00332 casalog.post("FIELD table direction reference frame entries for field "+str(fld) 00333 +" set to "+str(thenewref)+" ("+thenewrefstr+")", 'NORMAL') 00334 else: 00335 casalog.post("FIELD table direction reference frame entries for field "+str(fld) 00336 +" unchanged.", 'NORMAL') 00337 00338 00339 if(fixuvw and (oldrefcol!=[]) and (thenewref>0)): 00340 # modify reference of phase dir for fixuvw 00341 pcol = tbt.getcol('PhaseDir_Ref') 00342 pcol[fld] = 0 # J2000 00343 tbt.putcol('PhaseDir_Ref', pcol) # J2000 00344 00345 tbt.close() 00346 00347 #modify SOURCE table 00348 tbt.open(vis+'/SOURCE', nomodify=False) 00349 sdir = tbt.getcol('DIRECTION') 00350 newsdir = sdir 00351 sname = tbt.getcol('NAME') 00352 for i in xrange(0,tbt.nrows()): 00353 if(sname[i]==planetname): 00354 #print 'i old dir ', i, " ", sdir[0][i], sdir[1][i] 00355 newsdir[0][i] = thenewra_rad 00356 newsdir[1][i] = thenewdec_rad 00357 #print ' new dir ', newsdir[0][i], newsdir[1][i] 00358 tbt.putcol('DIRECTION', newsdir) 00359 tbt.close() 00360 casalog.post("SOURCE table DIRECTION column changed.", 'NORMAL') 00361 00362 # end for 00363 00364 if(fixuvw): 00365 casalog.post("Fixing the UVW coordinates ...", 'NORMAL') 00366 00367 # similar to fixvis 00368 00369 fldids = [] 00370 for i in xrange(numfields): 00371 if (i in fields): 00372 fldids.append(i) 00373 00374 imt = imtool() 00375 imt.open(vis, usescratch=False) 00376 imt.calcuvw(fldids, refcode='J2000', reuse=False) 00377 imt.close() 00378 imt = None 00379 00380 if((oldrefcol!=[]) and (thenewref>0)): 00381 tbt.open(vis+'/FIELD', nomodify=False) 00382 tbt.putcol('PhaseDir_Ref', oldrefcol) 00383 tbt.close() 00384 00385 else: 00386 casalog.post("UVW coordinates not changed.", 'NORMAL') 00387 00388 if (ParallelTaskHelper.isParallelMS(vis)): 00389 casalog.post("Tidying up the MMS subtables ...", 'NORMAL') 00390 ParallelTaskHelper.restoreSubtableAgreement(vis) 00391 00392 mst = None 00393 tbt = None 00394 00395 return True 00396 00397 except Exception, instance: 00398 mst = None 00399 tbt = None 00400 imt = None 00401 casalog.post("*** Error \'%s\' " % (instance), 'SEVERE') 00402 return False 00403
1.8.0
