casa
$Rev:20696$
|
00001 ########################################################################3 00002 # task_importgmrt.py 00003 # 00004 # 00005 # Copyright (C) 2009 00006 # Associated Universities, Inc. Washington DC, USA. 00007 # 00008 # This script is free software; you can redistribute it and/or modify it 00009 # under the terms of the GNU Library General Public License as published by 00010 # the Free Software Foundation; either version 2 of the License, or (at your 00011 # option) any later version. 00012 # 00013 # This library is distributed in the hope that it will be useful, but WITHOUT 00014 # ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or 00015 # FITNESS FOR A PARTICULAR PURPOSE. See the GNU Library General Public 00016 # License for more details. 00017 # 00018 # You should have received a copy of the GNU Library General Public License 00019 # along with this library; if not, write to the Free Software Foundation, 00020 # Inc., 675 Massachusetts Ave, Cambridge, MA 02139, USA. 00021 # 00022 # Correspondence concerning AIPS++ should be adressed as follows: 00023 # Internet email: aips2-request@nrao.edu. 00024 # Postal address: AIPS++ Project Office 00025 # National Radio Astronomy Observatory 00026 # 520 Edgemont Road 00027 # Charlottesville, VA 22903-2475 USA 00028 # 00029 # <author> 00030 # Shannon Jaeger (University of Calgary) 00031 # </author> 00032 # 00033 # <summary> 00034 # </summary> 00035 # 00036 # <reviewed reviwer="" date="" tests="" demos=""> 00037 # </reviewed 00038 # 00039 # <etymology> 00040 # importgmrt stands for import GMRT UV FITS data. 00041 # </etymology> 00042 # 00043 # <synopsis> 00044 # task_importgmrt.py is a Python script providing an easy to use task 00045 # for importing GMRT FITS data files into the CASA system. A measurement 00046 # set is produced from the GMRT FITS file. 00047 # 00048 # The given flag files are read and data is flagged according to 00049 # the information found in these files. It is expected that the 00050 # FLAG files are in the old AIPS TVFLAG file format. 00051 # 00052 # </synopsis> 00053 # 00054 # <example> 00055 # <srcblock> 00056 00057 # </srblock> 00058 # </example> 00059 # 00060 # <motivation> 00061 # To provide a user-friendly method for importing GMRT FITS files. 00062 # </motivation> 00063 # 00064 # <todo> 00065 # </todo> 00066 00067 import numpy 00068 import os 00069 from importuvfits import * 00070 from taskinit import * 00071 #from fg import * 00072 00073 def importgmrt( fitsfile, flagfile, vis ): 00074 retValue=False 00075 00076 # First check if the output file exists. If it 00077 # does then we abort. CASA policy prevents files 00078 # from being over-written. 00079 # 00080 # Also if the vis name given is the empty string 00081 # we use the default name, input fitsfile with 00082 # "fits" removed 00083 if ( len( vis ) < 1 ): 00084 vis = '' 00085 fits_list=fitsfile.split( '.' ) 00086 last = len(fits_list) 00087 if ( fits_list[last-1].lower() == 'fits' or \ 00088 fits_list[last-1].lower() == 'fit' ): 00089 last = last -1 00090 for i in range( last ): 00091 if ( len ( vis ) < 1 ): 00092 vis = fits_list[i] 00093 else: 00094 vis = vis + '.'+fits_list[i] 00095 vis=vis+'.MS' 00096 00097 casalog.post( "The vis paramter is empty, consequently the" \ 00098 +" the default visibilty file\nname will be used, "\ 00099 + vis, 'WARN' ) 00100 if ( len( vis ) > 0 and os.path.exists( vis ) ): 00101 casalog.post( 'Visibility file, '+vis+\ 00102 ' exists. import GMRT can not proceed, please\n'+\ 00103 'remove it or set vis to a different value.', 'SEVERE' ) 00104 return retValue 00105 00106 00107 # Loop through the list of flag files and make sure each one 00108 # of them exists 00109 if isinstance( flagfile, str ): 00110 if ( len( flagfile ) > 0 ): 00111 flagfile=[flagfile] 00112 else: 00113 flagfile=[] 00114 00115 ok = True 00116 for i in range( len( flagfile ) ): 00117 if ( not os.path.exists( flagfile[i] ) ): 00118 casalog.post( 'Unable to locate FLAG file '+ flagfile[i], 00119 'SEVERE' ) 00120 ok = False 00121 00122 if ( not ok ): 00123 return retValue 00124 00125 00126 # Ok, we've done our preliminary checks, now let's get 00127 # down to business and import our fitsfile. 00128 try: 00129 casalog.post( 'Starting import ...', 'NORMAL' ) 00130 importuvfits( fitsfile, vis ) 00131 except Exception, instance: 00132 casalog.post( str(instance), 'SEVERE' ) 00133 return retValue 00134 00135 # mytb, myms, myfg = gentools(['tb', 'ms', 'fg']) 00136 mytb, myms = gentools(['tb', 'ms']) 00137 aflocal = casac.agentflagger() 00138 00139 # Write history 00140 try: 00141 param_names = importgmrt.func_code.co_varnames[:importgmrt.func_code.co_argcount] 00142 param_vals = [eval(p) for p in param_names] 00143 ok &= write_history(myms, vis, 'importgmrt', param_names, 00144 param_vals, casalog) 00145 except Exception, instance: 00146 casalog.post("*** Error \'%s\' updating HISTORY" % (instance), 00147 'WARN') 00148 00149 # CASA's importuvfits doesn't understand the GMRT antenna 00150 # names very well and/or the FITS files tend to put the 00151 # antenna names in a place that CASA doesn't look for them. 00152 # Luckily the information seems to be in the "station" 00153 # information. So let's fix up the Name column. 00154 casalog.post( 'Correcting GMRT Antenna names.', 'NORMAL1' ) 00155 try: 00156 mytb.open( vis+'/ANTENNA', nomodify=False ) 00157 stations = mytb.getcol('STATION') 00158 names = [] 00159 for idx in range( 0, len( stations ) ): 00160 names.append( stations[idx].split(':')[0] ) 00161 mytb.putcol( 'NAME', numpy.array( names ) ) 00162 mytb.done() 00163 except Exception, instance: 00164 casalog.post( 'Unable to properly name the antennas, but continuing') 00165 casalog.post( str(instance), 'WARN' ) 00166 00167 # If we don't have a flagfile then we are done! 00168 # Yippee awe eh! 00169 if ( len( flagfile ) < 1 ): 00170 return True 00171 00172 # We've imported let's find out the observation start 00173 # and end times, this will be needed for flagging. 00174 # We use it to find the different days the observation 00175 # was done on. 00176 startObs = '' 00177 endObs = '' 00178 try: 00179 myms.open( vis ) 00180 trange=myms.range( 'time' )['time'] 00181 myms.done() 00182 00183 # Now have the observation time range in a numpy array. 00184 startObs = qa.time( str( trange[0] )+'s', prec=8, form='ymd' )[0] 00185 endObs = qa.time( str( trange[1] )+'s', prec=8, form='ymd' )[0] 00186 00187 except Exception, instance: 00188 casalog.post( 'Unable to find obaservation start/en times', 'SEVERE' ) 00189 casalog.post( str(instance), 'SEVERE' ) 00190 return retValue 00191 00192 days=[] 00193 startYY = startObs.split('/')[0] 00194 startMM = startObs.split('/')[1] 00195 startDD = startObs.split('/')[2] 00196 00197 endYY = endObs.split('/')[0] 00198 endMM =endObs.split('/')[1] 00199 endDD =endObs.split('/')[2] 00200 00201 for year in range( int( startYY ), int( endYY ) + 1 ): 00202 for month in range( int( startMM ), int( endMM ) + 1 ): 00203 for day in range( int( startDD ), int( endDD ) + 1 ): 00204 days.append( str(year) + '/' + str(month) + '/' + str(day) + '/' ) 00205 casalog.post( "DAYS: "+str(days), 'DEBUG1') 00206 00207 00208 # Read through the flag file(s) and assemble it into 00209 # a data structure. The hope is to minimize the 00210 # number of calls to "flagdata" to make the flagging 00211 # fast! 00212 # 00213 # The expected file format is what is the flag file 00214 # format used by TVFLAG/UVFLAG in AIPS, which is as follows: 00215 # 1. If a line starts with '!' it's a comment 00216 # 2. 00217 # 00218 casalog.post( 'Starting the flagging ...', 'NORMAL' ) 00219 00220 # First lets save the flag information as we got it from the 00221 # flag file. 00222 # myfg.open( vis ) 00223 # myfg.saveflagversion( 'none', 'No flagging performed yet' ) 00224 # myfg.done() 00225 aflocal.open( vis ) 00226 aflocal.saveflagversion( 'none', 'No flagging performed yet' ) 00227 aflocal.done() 00228 00229 for file in flagfile: 00230 casalog.post( 'Reading flag file '+file, 'NORMAL2' ) 00231 FLAG_FILE = open( file, 'r' ) 00232 line = FLAG_FILE.readline() 00233 ant_names = [] 00234 00235 while ( len( line ) > 0 ): 00236 casalog.post( 'Read from flag file: '+str(line), 'DEBUG1' ) 00237 00238 # Default antenna list and time range 00239 antennas = [] 00240 baselines = [] 00241 timerange = [] 00242 00243 # Skip comment lines, and the end of file. 00244 if ( len(line) < 1 or ( line[0] == '!' and line[2:6] != 'PANT' ) ): 00245 line = FLAG_FILE.readline() 00246 continue 00247 00248 # Store the antenna name list 00249 if ( line[0:6] == '! PANT' ): 00250 ant_names= line.split(':')[1].lstrip(' ').split(' ') 00251 line = FLAG_FILE.readline() 00252 continue 00253 00254 # First divide lines up on spaces 00255 casalog.post( 'LINE READ: '+line, 'DEBUG1' ) 00256 step1=line.split(' ') 00257 if ( len( step1 ) < 2 ): 00258 # We want a line with something on it. 00259 line = FLAG_FILE.readline() 00260 continue 00261 00262 # Parse each bit of the string looking for the 00263 # antenna and timrange information. 00264 for i in range( len( step1 ) ): 00265 step2=step1[i].split('=') 00266 00267 if ( len( step2 ) != 2 ): 00268 # We want something with an equals sign in it! 00269 continue 00270 casalog.post( 'keyword/value: '+str(step2), 'DEBUG1' ) 00271 00272 00273 if ( step2[0].upper().startswith( 'ANT' ) ): 00274 # Should be a comma separated list which 00275 # is acceptable to flag data 00276 #antennas=list(step2[1].split(',') 00277 antennas=step2[1] 00278 casalog.post( "antennas: "+antennas, 'DEBUG1') 00279 00280 00281 if ( step2[0].upper().startswith( 'BASEL' ) ): 00282 # Should be a comma separated list which 00283 # is acceptable to flag data 00284 baselines=step2[1] 00285 casalog.post( "baselines: "+baselines, 'DEBUG1') 00286 00287 if ( step2[0].upper() == 'TIMERANGE' ): 00288 # Time in the FLAG file is in the format 00289 # start-day,start-hour,start-min,star-sec,\ 00290 # end-day,end-hour,end-min,end-sec 00291 # 00292 # CASA expects data to be of the form 00293 # YYYY/MM/DD/HH:MM:SS~YYYY/MM/DD/HH:MM:SS 00294 # http://casa.nrao.edu/Memos/msselection/node8.html 00295 # 00296 # We use the day # to index into our days list and 00297 # append the time on the end. But day # starts at 00298 # 1 and indexing starts at 0 so we need to adjust. 00299 times = step2[1].split(',') 00300 casalog.post( "time: "+str(times), 'DEBUG1') 00301 if ( int(times[0]) < 1 and int(times[4]) > len(days) ): 00302 # All data for this antenna 00303 timerange = '' 00304 elif ( int(times[0]) < 1 ): 00305 # We want all time before the end time 00306 timerange = '<' + days[int(times[4])-1] \ 00307 + times[5] + ':' + times[6] + ':' +times[7] 00308 elif ( int(times[4]) > len(days) ): 00309 # We want all times after the start time 00310 timerange = '>' + days[int(times[0])-1] \ 00311 + times[1]+':'+times[2]+':'+times[3] 00312 else: 00313 # We ahve a fully specified time range 00314 timerange = days[int(times[0])-1] \ 00315 + times[1]+':'+times[2]+':'+times[3]\ 00316 +'~'\ 00317 +days[int(times[4])-1] \ 00318 +times[5]+':'+times[6]+':'+times[7] 00319 00320 # Construct the antenna query string 00321 # 00322 # Verify the syntax of the ANTENNA & BASELINE keywords, 00323 # if ANTENNA will always have only one. I suspect we 00324 # with to loop through the baselines only here. 00325 # 00326 # OLD CODE FOR creating ANT SELECTION STRING 00327 # 00328 #antStr='' 00329 #selectAntStr='' 00330 #if ( len( antennas ) > 0 and len( baselines ) > 0 ): 00331 # # We are selecting baselines. 00332 # for ant1 in antennas: 00333 # for ant2 in baselines: 00334 # if ( len( antStr ) < 1 ): 00335 # antStr=str(int(ant1)-1)+'&'+str(int(ant2)-1) 00336 # selectAntStr=antStr 00337 #elif ( len( antennas ) > 0 ): 00338 # antStr=str(int(antennas)-1) 00339 # # A hack to make sure we always have some unflagged data. 00340 # # CASA seg faults if the myfg.setdata() results in no data. 00341 # selectAntStr+=antStr+','+antennas 00342 00343 00344 antStr='' 00345 #if ( len( antennas ) > 0 and len( baselines ) > 0 ): 00346 # # We are selecting baselines. 00347 # for ant1 in antennas: 00348 # for ant2 in baselines: 00349 # if ( len( antStr ) < 1 ): 00350 # antStr=str(int(ant1)-1)+'&'+str(int(ant2)-1) 00351 # selectAntStr=antStr 00352 #elif ( len( antennas ) > 0 ): 00353 if ( len( antennas ) > 0 ): 00354 antStr = ant_names[ int(antennas)-1 ].strip( ' ' ) 00355 00356 try: 00357 # casalog.post( "flagdata( "+vis+", mode='manualflag', antenna='" 00358 # +antStr+"', timerange='"+timerange+"' )", 'NORMAL') 00359 casalog.post( "flagdata( "+vis+", mode='manual', antenna='" 00360 +antStr+"', timerange='"+timerange+"' )", 'NORMAL') 00361 #flagdata( vis, mode='manualflag', selectdata=True, \ 00362 # antenna=antStr, timerange=timerange, \ 00363 # flagbackup=False ) 00364 # myfg.open( vis ) 00365 # #print "myfg.setdata( time='",timerange,"' )" 00366 # myfg.setdata( time=timerange ) 00367 # #print "myfg.setmanualflags( baseline='",antStr,"', time='",timerange,"' )" 00368 # myfg.setmanualflags( baseline=antStr, time=timerange, unflag=False, clipexpr='' ) 00369 # myfg.run() 00370 # myfg.done() 00371 aflocal.open(vis) 00372 aflocal.selectdata(timerange=timerange) 00373 aflocal.parsemanualparameters(antenna=antStr, timerange=timerange) 00374 aflocal.init() 00375 aflocal.run(True, True) 00376 aflocal.done() 00377 except Exception, instance: 00378 casalog.post( 'Unable to flag data from flag file '+file\ 00379 +'.\nAntennas='+antennas+' and timerange='\ 00380 +timerange, 'WARN' ) 00381 casalog.post( str(instance), 'SEVERE' ) 00382 return retValue 00383 00384 line = FLAG_FILE.readline() 00385 00386 FLAG_FILE.close() 00387 00388 # Save a Flag version so we can revert back to it if we wish 00389 # myfg.open( vis ) 00390 # myfg.saveflagversion( 'import', 'Flagged the data from the GMRT flag file' ) 00391 # myfg.done() 00392 aflocal.open( vis ) 00393 aflocal.saveflagversion( 'import', 'Flagged the data from the GMRT flag file' ) 00394 aflocal.done() 00395 00396 retValue = True 00397 return retValue 00398 00399 00400 00401 00402