Changeset 57 for CedaCdoScripts
- Timestamp:
- 23/09/13 17:14:25 (7 years ago)
- Location:
- CedaCdoScripts/Misc
- Files:
-
- 3 added
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
CedaCdoScripts/Misc/cdo_basic.py
r54 r57 8 8 debug = True 9 9 cleanUpAsUGo = True 10 args = getopt.getopt( sys.argv[1:], 'i:m:e:r:b:x:o:', ['onlyVars=','scratch=','omitVars=' ] )10 args = getopt.getopt( sys.argv[1:], 'i:m:e:r:b:x:o:', ['onlyVars=','scratch=','omitVars=','version=','log='] ) 11 11 ad = { 'i':'inst', 'm':'model', 'e':'ensem', 'r':'realm', 'b':'batch', 'x':'expt', 'o':'outDir' } 12 12 13 knownops = ['am',' bm','mclim5','cdd','tr']13 knownops = ['am','sm','bm','mclim5','cdd','tr', 'fd'] 14 14 optag = args[-1][0] 15 15 assert optag in knownops, 'First argument [%s] should be operation in: %s' % (optag, str(knmownops)) … … 22 22 fullYears = True 23 23 24 print args25 24 edict = {} 26 25 for a in args[0]: 27 26 edict[a[0][1:]] = a[1] 27 28 def writeat( fh, tag, str,log): 29 if len(str) + len(tag) + 2 > 1024 : 30 log.warn( 'TRUNCATING global attribute %s: %s' % (tag,str) ) 31 strl = 1022 - len(tag) 32 str = str[:strl] 33 fh.write( '%s %s\n' % (tag,str) ) 28 34 29 35 fullYearNeeded = False … … 50 56 ##batch = 'b011_' 51 57 ##onlyVars = ['tas'] 58 project = 'cmip5' 52 59 ensem = edict['e'] 53 60 inst = edict['i'] … … 56 63 experiment = edict['x'] 57 64 batch = edict['b'] + '_' 65 processing_batch = edict['b'] 58 66 opdir = edict['o'] 59 67 onlyVars = None … … 64 72 ## qcDir = '%s%s/%s/%s/%s/' % (qcBaseDir,product,inst, model, expt) 65 73 66 tstring = '%s%2.2i%2.2i-%2.2i%2.2i%2.2i' % time.gmtime()[0:6]67 74 vstring = '%s%2.2i%2.2i' % time.gmtime()[0:3] 68 75 doctag = optag 69 76 if optag in ['cdd']: 70 77 doctag = 'precip_indices' 71 if optag in ['tr' ]:78 if optag in ['tr','fd']: 72 79 doctag = 'temperature_indices' 73 80 doc_filename = 'ceda_cmip5_processing_%s.pdf' % doctag 74 log_filename = 'ceda_cmip5_processing_log_%s_%s.txt' % (optag,tstring) 81 tstring = '%s%2.2i%2.2i-%2.2i%2.2i%2.2i' % time.gmtime()[0:6] 82 log_filename = edict.get( '-log', 'ceda_cmip5_processing_log_%s_%s.txt' % (optag,tstring) ) 75 83 batching='byFile' 76 84 ofreq = 'yr' … … 79 87 batching = '5yr' 80 88 monthsPerBatch = 60 81 elif optag in ['cdd','tr' ]:89 elif optag in ['cdd','tr','fd']: 82 90 batching = 'yr' 83 91 ifreq = 'day' … … 90 98 91 99 scriptName = 'cdo_basic.py' 92 scriptVersion = 0. 8100 scriptVersion = 0.9 93 101 log = logging.getLogger( __name__ ) 94 102 fh = logging.FileHandler( log_filename ) … … 108 116 omitVars = string.split( edict['-omitVars'], ',' ) 109 117 log.info( 'Omitting variables: %s' % edict['-omitVars'] ) 118 if edict.has_key( '-version' ): 119 version = edict['-version'] 120 log.info( 'Processing version: %s' % edict['-version'] ) 121 else: 122 version = None 110 123 ## log.info( 'CDO version information:' ) 111 124 ## os.popen( 'cdo -V > .cdo_version' ).readlines() … … 123 136 if not os.path.isdir(tmpdir): 124 137 os.mkdir( tmpdir ) 125 print 'Creating temporary directory: %s ' % tmpdir138 log.info( 'Creating temporary directory: %s ' % tmpdir ) 126 139 127 140 if tmpdir[-1] != '/': … … 133 146 assert os.path.isfile( '%scdo_basic_lock' % tmpdir ), 'Failed to make lock file in work directory %s' % tmpdir 134 147 135 136 148 plurality = 'single' 137 149 mdbase = "http://badc.nerc.ac.uk/browse/badc/cmip5/metadata/processing/" … … 139 151 rcout = open( 'cdo_stdout.txt', 'w' ) 140 152 rcerr = open( 'cdo_stderr.txt', 'w' ) 141 def run_cmd( cmd, comment='', strict=True ):153 def run_cmd( cmd, comment='', strict=True, append=None, out=None ): 142 154 if type(cmd) == type('x'): 143 155 args = map( lambda x: string.strip(x, "'"), string.split( cmd )) … … 145 157 args = cmd 146 158 log.info( cmd ) 147 proc = subprocess.Popen( args, stdout=rcout, stderr=rcerr, env={ 'LD_LIBRARY_PATH':'/usr/local/lib'} ) 159 if append != None: 160 fho = open( append, 'a' ) 161 clo = True 162 elif out != None: 163 fho = open( out, 'w' ) 164 clo = True 165 else: 166 fho = rcout 167 clo = False 168 proc = subprocess.Popen( args, stdout=fho, stderr=rcerr, env={ 'LD_LIBRARY_PATH':'/usr/local/lib'} ) 148 169 std = proc.communicate() 149 170 rc = proc.returncode … … 156 177 elif rc != 0: 157 178 print 'Failed (%s) [rc=%s] executing %s' % (comment,rc,str(cmd)) 179 if clo: 180 fho.close() 158 181 return rc 159 182 160 183 dir_tmpl_s = '/disks/fuse/cmip5/output1/MOHC/HadGEM2-ES/%(expt)s/%(ifreq)s/atmos/Amon/' 184 dir_tmpl_s = '/badc/cmip5/data/cmip5/output1/MOHC/HadGEM2-ES/%(expt)s/%(ifreq)s/atmos/Amon/' 161 185 if ifreq == 'mon': 162 186 mip = {'atmos':'Amon', 'ocean':'Omon', 'land':'Lmon' }[realm] … … 166 190 raise 'Should not be here' 167 191 168 dir_tmpl_r = '/ disks/fuse/cmip5/output1/%(inst)s/%(model)s/%(expt)s/%(ifreq)s/%(realm)s/%(mip)s/%(ensem)s/'169 dir_tmpl_e = '/ disks/fuse/cmip5/output1/%(inst)s/%(model)s/%(expt)s/%(ifreq)s/%(realm)s/%(mip)s/%(ensem)s/latest/'192 dir_tmpl_r = '/badc/cmip5/data/cmip5/output1/%(inst)s/%(model)s/%(expt)s/%(ifreq)s/%(realm)s/%(mip)s/%(ensem)s/' 193 dir_tmpl_e = '/badc/cmip5/data/cmip5/output1/%(inst)s/%(model)s/%(expt)s/%(ifreq)s/%(realm)s/%(mip)s/%(ensem)s/latest/' 170 194 dir_tmpl = dir_tmpl_e + '%(var)s/' 171 195 bldir_tmpl = '/home/users/prototype/cmip_processing/cmip5_bl/output1_%(inst)s_%(model)s_%(expt)s/' … … 180 204 if len(b2) == 0: 181 205 return None 182 return b2[-1] 183 184 def get_version_dir(): 206 var = string.split( b2[-1], '_' )[0] 207 return var 208 209 def get_version_dir(version=None): 185 210 186 211 pdict = { 'inst':inst, 'model':model, 'expt':expt, 'ensem':ensem, 'ifreq':ifreq, 'realm':realm, 'mip':mip } 187 212 rdir = dir_tmpl_r % pdict 213 if version != None: 214 this_v = rdir + version 215 assert os.path.isdir( this_v ), 'Version directory %s not found' % this_v 216 return this_v 217 188 218 vlist = glob.glob( rdir + 'v*' ) 189 219 assert len(vlist) > 0, 'No version directories found in %s' % rdir … … 200 230 ## edir = dir_tmpl % { 'inst':inst, 'model':model, 'expt':expt, 'var':var, 'ensem':ensem, 'ifreq':ifreq, 'realm':realm, 'mip':mip } 201 231 assert os.path.isdir( edir), 'Error in get_expt_fl: %s is not a directory' % edir 232 log.info( 'edir:: %s' % edir ) 202 233 fl = glob.glob( edir + '*.nc' ) 203 234 assert len(fl) > 0, 'No netcdf files found in %s' % edir … … 211 242 if os.path.isfile( bldir + blfile ): 212 243 blvars = map( get_var_from_bl_record, open( bldir + blfile ).readlines() ) 213 if var in blvars: 214 print 'skipping ',var 215 log.info( 'Skipping black-listed ADS expt=%s, var=%s' % (expt,var) ) 216 return fm 244 if var in blvars: 245 log.info( 'Skipping black-listed ADS expt=%s, var=%s' % (expt,var) ) 246 return fm 247 else: 248 log.info( 'Var %s not in black list [%s] %s' % (var, bldir + blfile, str(blvars)) ) 249 else: 250 log.info( 'Black list file %s not present' % (bldir+blfile) ) 251 else: 252 log.info( 'Black list directory %s not present' % bldir ) 217 253 for f in fl: 218 254 fn = string.split(f, '/')[-1] … … 268 304 onlyVars = ['pr'] 269 305 omitVars = None 270 if optag in ['tr' ]:306 if optag in ['tr','fd']: 271 307 onlyVars = ['tasmin'] 272 308 omitVars = None … … 274 310 for expt in expts: 275 311 ##dir = dir_tmpl_e % locals() 276 dir = get_version_dir( )312 dir = get_version_dir( version=version) 277 313 vars = glob.glob( '%s/*' % dir ) 278 print expt, dir, vars279 314 thisvarl = map( lambda x: string.split(x,'/')[-1], vars ) 280 315 thisvarl.sort() … … 312 347 lon1,lon2,lat1,lat2 = (60.,120.,30.,60.) 313 348 cdoop, ofreq,opName, proc_description = { 'am':('divdpy -yearsum -muldpm', 'yr','Annual mean','Annual mean of monthly mean data'), \ 349 'sm':('divdpy -seassum -muldpm', 'yr','Seasonal mean','Annual seasonal of monthly mean data'), \ 314 350 'bm':('fldmean -sellonlatbox,%s,%s,%s,%s' % (lon1,lon2,lat1,lat2),ifreq,'Box mean','Spatial mean over a latitude-longitude box'), \ 315 351 'mclim5':('ymonavg', '5yr','5-year monthly climatology','Monthly climatology for 5 year period: 5 year mean for each month'), 316 352 'cdd':('eca_cdd','yr','Continuous dry days','Maximum of continuous dry days in a year'), 317 'tr':('eca_tr','yr','Tropical nights','Count of tropical nights per year') }[optag] 353 'tr':('eca_tr','yr','Tropical nights','Count of tropical nights per year'), 354 'fd':('eca_fd','yr','Frost days','Count of frost days (Tmin below zero) per year') }[optag] 318 355 ## opdir = '%s%s' % (batch,optag) 319 356 if go and not os.path.isdir( opdir ): … … 327 364 k=-1 328 365 print 'INFO0080: ',k, str(k), str( [1,2,3] ) 329 for expt,var,fm in ye: 366 ##for expt,var,fm in ye: 367 def run_evf( expt,var,fm ): 330 368 kk = 0 331 dir = get_version_dir( )369 dir = get_version_dir( version=version) 332 370 ##edir = dir_tmpl % locals() 333 371 edir = '%s/%s/' % (dir,var) … … 367 405 ifiles[key] += '%s,' % file[0] 368 406 407 startstart = None 369 408 nbfb = 0 370 409 nbfbb = 0 … … 395 434 nodir = False 396 435 if nbfb > 0: 397 print 'FILE DATE BOUNDARIES DO NOT MATCH ASSUMPTIONS: %s, %s' % (expt,var)398 print 'merging files'399 print nbfb, nbfbb400 436 log.info( 'File date boundaries do not match end of year (%s,%s): merging: %s, %s' % (nbfb, nbfbb, expt,var) ) 401 437 if nbfbb > 0: … … 441 477 cmd1 = None 442 478 cmd = '/usr/local/2/bin/cdo -s %s %s/%s %s' % (cdoop, tmpdir, '.tmp1.nc', ofile) 443 print cmd479 log.info( 'cdo command: %s\n' % cmd ) 444 480 ## raise 'bad file boundary' 445 481 … … 471 507 nf = len(fm) 472 508 if smem[0][0] != 1: 473 print 'Omitting first file: %s' % str(fm[0])509 log.info( 'Omitting first file: %s' % str(fm[0]) ) 474 510 i0 = 1 475 511 nf = len(fm)-1 … … 514 550 ofiles[key] += ' %s' % (ofile) 515 551 if nbfb > 0: 516 print 'FILE DATE BOUNDARIES DO NOT MATCH ASSUMPTIONS: %s, %s' % (expt,var)552 log.info( 'FILE DATE BOUNDARIES DO NOT MATCH ASSUMPTIONS: %s, %s' % (expt,var) ) 517 553 elif batching in ['5yr','yr']: 518 554 startstart = fm[0][1] … … 566 602 567 603 for k in range(len(tfm)-1): 568 print 'INFO0100: ',k, str(k), str( [1,2,3])569 print 'INFO0101: ',k, tfm[k], tfm[k+1]604 log.info( 'INFO0100: ',k, str(k), str( [1,2,3] ) ) 605 log.info( 'INFO0101: ',k, tfm[k], tfm[k+1] ) 570 606 a = str(tfm[k]) 571 607 b = str(tfm[k+1]) … … 587 623 if debug: 588 624 log.info( 'Trying to log file maximum value' ) 589 cmd 0 = 'echo %s%s >> /home/users/prototype/cmip_processing/wrk/cdd1/debug.txt' % (edir,tfm[0][0])590 cmd = '/usr/local/2/bin/cdo info %s 2>&1 >> /home/users/prototype/cmip_processing/wrk/cdd1/debug.txt' % ofile591 rc = run_cmd( cmdo )592 rc = run_cmd( cmd )625 cmdo = 'echo %s%s' % (edir,tfm[0][0]) 626 cmd = '/usr/local/2/bin/cdo info %s' % ofile 627 rc = run_cmd( cmdo, append='/home/users/prototype/cmip_processing/wrk/cdd1/debug.txt' ) 628 rc = run_cmd( cmd, append='/home/users/prototype/cmip_processing/wrk/cdd1/debug.txt' ) 593 629 else: 594 630 print 'single: %s' % cmd … … 635 671 else: 636 672 raise 'Code missing' 637 638 attfile = '%s/temp_attribute_file.txt' % tmpdir 639 log.info( 'Merging files' ) 640 if go: 641 for expt in expts: 673 return tid, crd, table_id, realization, initialization_method, physics_version, startstart 674 675 676 for expt,var,fm in ye: 677 tid, crd, table_id, realization, initialization_method, physics_version, startstart = run_evf( expt,var,fm ) 678 679 def proc_expt(expt,version,optag): 680 attfile = '%s/temp_attribute_file.txt' % tmpdir 642 681 varl = varld[expt] 643 682 ## for var in varld[expt]: 644 683 for var in varl: 684 if optag in ['am', 'sm']: 685 if var == 'tasmax': 686 cmval = "time: maximum within days time: mean over days" 687 elif var == 'tasmin': 688 cmval = "time: minimum within days time: mean over days" 689 else: 690 cmval = "time: mean" 691 else: 692 assert True, 'Need to set cell_methods for this operator' 693 645 694 key = '%s_%s' % (expt,var) 646 695 if key in badkeys: 647 print 'Skipping %s because or earlier error' % key696 log.info( 'Skipping %s because or earlier error' % key ) 648 697 elif string.strip( ofiles[key] ) == '': 649 print 'Skipping %s because no output files registered' % key698 log.info( 'Skipping %s because no output files registered' % key ) 650 699 log.warning( 'WARN004: unexplained missing output file for key %s' % (key) ) 651 700 else: 652 print 'Trying %s' % key653 dir = get_version_dir( )701 log.info( 'Trying %s' % key ) 702 dir = get_version_dir( version=version) 654 703 ## vdir = dir_tmpl % locals() 655 704 vdir = '%s/%s/' % (dir,var) … … 667 716 for k in range(36): 668 717 tid += cx[random.randint(0,len(cx)-1 )] 718 this_mip = string.split( table_id )[1] 719 version = string.split( string.strip( dir, '/' ), '/' )[-1] 720 esgf_dsid = '%s.%s.%s.%s.%s.%s.%s.%s.%s.%s' % (project,product,inst,expt,model,ifreq,realm,this_mip,ensem,version) 721 ##cmip5.output1.MIROC.MIROC5.sstClim.monClim.atmos.Amon.r1i1p1.v20111129.xml 669 722 oo = open( attfile, 'w' ) 670 723 kkey = '%s_%s' % (expt,var) … … 681 734 oo.write( '%s %s\n' % ('operation_qualifiers ', ' ' ) ) 682 735 oo.write( '%s %s\n' % ('input_frequency ', ifreq ) ) 736 oo.write( '%s %s\n' % ('input_esgf_dsid ', esgf_dsid ) ) 683 737 oo.write( '%s %s\n' % ('model_id ', model ) ) 684 738 oo.write( '%s %s\n' % ('experiment_id ', expt ) ) … … 694 748 oo.write( '%s %s\n' % ('input_dir ', vdir ) ) 695 749 oo.write( '%s %s\n' % ('plurality ', plurality ) ) 750 oo.write( '%s %s\n' % ('proc_code_repository_url ', 'http://proj.badc.rl.ac.uk/svn/exarch/CedaCdoScripts' ) 751 oo.write( '%s %s\n' % ('proc_code_home_url ', 'http://proj.badc.rl.ac.uk/exarch/wiki/PackageRBPS' ) 696 752 oo.write( '%s %s\n' % ('proc_method_url ', '%s/documentation/%s/%s' % (mdbase, plurality, doc_filename) ) ) 697 753 oo.write( '%s %s\n' % ('proc_log_url ', '%s/logfiles/%s/%s/%s' % (mdbase, plurality, optag, log_filename) ) ) … … 715 771 attr += ',' + dat 716 772 attr += ']' 717 oo.write( '%s %s\n' % ('input_files ', attr))773 writeat(oo, 'input_files ', attr, log ) 718 774 else: 719 oo.write( '%s %s\n' % ('input_files ', ifiles[kkey] ))775 writeat(oo, 'input_files ', ifiles[kkey], log ) 720 776 721 777 if len( tids[kkey] ) > 511: … … 726 782 oo.write( '%s %s\n' % ('input_creation_date ', crds[kkey]) ) 727 783 oo.write( '%s %s\n' % ('tracking_id ', tid) ) 728 oo.write( '%s %s; %s\n' % ('history ', hists[kkey], '%s: %s: %s' % (creation_date, 'CDO script', proc_description)))784 writeat( oo, 'history ', '%s; %s' % (string.replace(hists[kkey],'\n',';;'), '%s: %s: %s' % (creation_date, 'CDO script', proc_description)), log ) 729 785 oo.close() 730 786 tmpofilename = tmpdir + '.cdo_basic_temp.nc' … … 745 801 assert os.path.isfile( string.strip(ofiles[kkey]) ), 'File "%s" not found [2], kkey = %s' % (ofiles[kkey], kkey) 746 802 cmd = '/usr/local/2/bin/cdo -s setgatts,%s %s %s/%s' % (attfile,ofiles[kkey],opdir,ofilename) 747 print cmd803 log.info( 'Command: %s\n' % cmd ) 748 804 rc = run_cmd( cmd ) 749 805 750 log.info( 'Script completed' ) 751 rcout.close() 752 rcerr.close() 753 fh.close() 754 for k in eeproc.keys(): 755 sh[procid][k] = eeproc[k] 756 print sh.keys() 757 sh.close() 806 ###cmd = '/usr/local/2/bin/ncatted -a cell_methods,%s,o,c,%s %s' % (var,cmval,ofilename) 807 cmd = ['/usr/local/2/bin/ncatted','-a','cell_methods,%s,o,c,%s' % (var,cmval),'%s/%s' % (opdir,ofilename)] 808 rc = run_cmd( cmd ) 809 810 def close_all(): 811 log.info( 'Script completed' ) 812 fh.close() 813 rcout.close() 814 rcerr.close() 815 for k in eeproc.keys(): 816 sh[procid][k] = eeproc[k] 817 sh.close() 818 819 if go: 820 log.info( 'Merging files' ) 821 for expt in expts: 822 proc_expt(expt,version,optag) 823 824 close_all() 758 825 759 826 os.popen( 'rm %scdo_basic_lock' % tmpdir ).readlines()
Note: See TracChangeset
for help on using the changeset viewer.