diff --git a/python/CONTROL_OD.glob.025 b/python/CONTROL_OD.glob.025 new file mode 100644 index 0000000000000000000000000000000000000000..fc96137ab014cd7519112fb2f4177d9144151400 --- /dev/null +++ b/python/CONTROL_OD.glob.025 @@ -0,0 +1,34 @@ +DAY1 +DAY2 +DTIME 3 +M_TYPE AN FC FC FC FC FC FC FC FC FC FC FC AN FC FC FC FC FC FC FC FC FC FC FC +M_TIME 00 00 00 00 00 00 00 00 00 00 00 00 12 12 12 12 12 12 12 12 12 12 12 12 +M_STEP 00 01 02 03 04 05 06 07 08 09 10 11 00 01 02 03 04 05 06 07 08 09 10 11 +M_CLASS OD +M_STREAM OPER +M_NUMBER OFF +M_EXPVER +M_GRID 250 +M_LEFT -179750 +M_LOWER -90000 +M_UPPER 90000 +M_RIGHT 180000 +M_LEVEL 137 +M_RESOL 799 +M_GAUSS 0 +M_ACCURACY 24 +M_OMEGA 0 +M_OMEGADIFF 0 +M_ETA 1 +M_ETADIFF 0 +M_DPDETA 1 +M_SMOOTH 0 +M_FORMAT GRIB2 +M_ADDPAR /186/187/188/235/139/39 +PREFIX ENH +ECSTORAGE 0 +ECTRANS 1 +ECFSDIR ectmp:/${USER}/econdemand/ +MAILOPS ${USER} +MAILFAIL ${USER} +EOF diff --git a/python/CONTROL_ZAMG_SYNTH b/python/CONTROL_ZAMG_SYNTH new file mode 100644 index 0000000000000000000000000000000000000000..c930b2faeef880da33b357dba9534112aa537f75 --- /dev/null +++ b/python/CONTROL_ZAMG_SYNTH @@ -0,0 +1,31 @@ +DAY1 +DAY2 +DTIME 3 +M_TYPE PF PF PF PF PF PF PF PF +M_TIME 00 00 00 00 12 12 12 12 +M_STEP 00 03 06 09 00 03 06 09 +M_STREAM ENFO +PREFIX EN +M_UPPER 90000 +M_LOWER -90000 +M_LEFT -179500 +M_RIGHT 180000 +M_CLASS OD +M_NUMBER 1/TO/50 +M_GRID 500 +M_RESOL 319 +M_SMOOTH 0 +M_GAUSS 1 +M_ETA 1 +M_ETAPAR 77 +M_DPDETA 1 +M_LEVEL 91 +M_LEVELIST 1/TO/91 +M_ADDPAR 186/187/188/235/139/39/134 +M_FORMAT GRIB1 +ECSTORAGE 0 +ECTRANS 0 +DEBUG 1 +ECFSDIR ectmp:/${USER}/ecops +MAILOPS ${USER} +MAILFAIL ${USER} diff --git a/python/CONTROL_ZAMG_SYNTH_ENAN b/python/CONTROL_ZAMG_SYNTH_ENAN new file mode 100644 index 0000000000000000000000000000000000000000..1f7e74ca2a3c8339e5a3037f8da2d3e4b81afa0c --- /dev/null +++ b/python/CONTROL_ZAMG_SYNTH_ENAN @@ -0,0 +1,33 @@ +DAY1 +DAY2 +DTIME 3 +M_TYPE AN FC AN FC AN FC AN FC +M_TIME 00 18 06 06 12 06 18 18 +M_STEP 00 09 00 03 00 09 00 03 +M_ACCTIME 06/18 +M_ACCMAXSTEP 12 +M_STREAM ELDA +PREFIX EN +M_UPPER 90000 +M_LOWER -90000 +M_LEFT -179500 +M_RIGHT 180000 +M_CLASS OD +M_NUMBER 0/to/25 +M_GRID 500 +M_RESOL 319 +M_SMOOTH 0 +M_GAUSS 0 +M_ETA 1 +M_ETAPAR 77 +M_DPDETA 1 +M_LEVEL 137 +M_LEVELIST 1/TO/137 +M_ADDPAR 186/187/188/235/139/39/134 +M_FORMAT GRIB1 +ECSTORAGE 1 +ECTRANS 0 +DEBUG 0 +ECFSDIR ectmp:/${USER}/ecops +MAILOPS ${USER} +MAILFAIL ${USER} diff --git a/python/CONTROL_ZAMG_SYNTH_ENFO b/python/CONTROL_ZAMG_SYNTH_ENFO new file mode 100644 index 0000000000000000000000000000000000000000..ee0b7535897c467ca658ed4d353d6340001fb8a4 --- /dev/null +++ b/python/CONTROL_ZAMG_SYNTH_ENFO @@ -0,0 +1,31 @@ +DAY1 +DAY2 +DTIME 3 +M_TYPE PF PF PF PF PF PF PF PF +M_TIME 00 00 00 00 12 12 12 12 +M_STEP 00 03 06 09 00 03 06 09 +M_STREAM ENFO +PREFIX OP +M_UPPER 90000 +M_LOWER -90000 +M_LEFT -179500 +M_RIGHT 180000 +M_CLASS OD +M_NUMBER 1/TO/50 +M_GRID 500 +M_RESOL 319 +M_SMOOTH 0 +M_GAUSS 1 +M_ETA 0 +M_ETAPAR 77 +M_DPDETA 1 +M_LEVEL 91 +M_LEVELIST 1/TO/91 +M_ADDPAR 186/187/188/235/139/39/134 +M_FORMAT GRIB1 +ECSTORAGE 1 +ECTRANS 0 +DEBUG 0 +ECFSDIR ectmp:/${USER}/ecops +MAILOPS ${USER} +MAILFAIL ${USER} diff --git a/python/CONTROL_ZAMG_SYNTH_ENFOCF b/python/CONTROL_ZAMG_SYNTH_ENFOCF new file mode 100644 index 0000000000000000000000000000000000000000..2f8c377b2d9612b82f482ba871e277a56c54293a --- /dev/null +++ b/python/CONTROL_ZAMG_SYNTH_ENFOCF @@ -0,0 +1,31 @@ +DAY1 +DAY2 +DTIME 3 +M_TYPE CF CF CF CF CF CF CF CF +M_TIME 00 00 00 00 12 12 12 12 +M_STEP 00 03 06 09 00 03 06 09 +M_STREAM ENFO +PREFIX OP +M_UPPER 90000 +M_LOWER -90000 +M_LEFT -179500 +M_RIGHT 180000 +M_CLASS OD +M_NUMBER 0 +M_GRID 500 +M_RESOL 319 +M_SMOOTH 0 +M_GAUSS 1 +M_ETA 0 +M_ETAPAR 77 +M_DPDETA 1 +M_LEVEL 91 +M_LEVELIST 1/TO/91 +M_ADDPAR 186/187/188/235/139/39/134 +M_FORMAT GRIB1 +ECSTORAGE 1 +ECTRANS 0 +DEBUG 0 +ECFSDIR ectmp:/${USER}/ecops +MAILOPS ${USER} +MAILFAIL ${USER} diff --git a/python/FlexpartTools.py b/python/FlexpartTools.py index 41a806b89e7016502672bbf39355b30a41f441f7..345595415a22a2aa725cea6f7c33bd86facedef3 100644 --- a/python/FlexpartTools.py +++ b/python/FlexpartTools.py @@ -48,6 +48,7 @@ except ImportError: ecapi=False from gribapi import * from GribTools import GribTools +from opposite import opposite def interpret_args_and_control(*args,**kwargs): @@ -411,6 +412,32 @@ def darain(alist): return sfeld +def add_previousday_ifneeded(attrs): + + if attrs['type']=='FC' and 'acc' not in attrs['target']: + steps=attrs['step'].split('/') + times=attrs['time'].split('/') + addpday=False + #for t in times: + #for s in steps: + if int(times[0])+int(steps[0])>23: + addpday=True + if(addpday): + dates=attrs['date'].split('/') + start_date=dates[0] + syear=int(start_date[:4]) + smonth=int(start_date[4:6]) + sday=int(start_date[6:]) + start = datetime.date( year = syear, month = smonth, day = sday ) + startm1=start- datetime.timedelta(days=1) + + attrs['date']='/'.join([startm1.strftime("%Y%m%d")]+dates[1:]) + print('CHANGED FC start date to '+startm1.strftime("%Y%m%d")+'to accomodate TIME='+times[0]+', STEP='+steps[0]) + + return + + + class Control: 'class containing the information of the flex_extract control file' @@ -703,6 +730,9 @@ class MARSretrieval: mclass = attrs.get('marsclass') del attrs['marsclass'] attrs['class'] = mclass + + add_previousday_ifneeded(attrs) + target = attrs.get('target') if not int(self.public): del attrs['target'] @@ -768,7 +798,7 @@ class EIFlexpart: self.mreq_count = 0 self.types=dict() try: - if c.maxstep>len(c.type): # Pure forecast mode + if c.maxstep>24: #len(c.type): # Pure forecast mode c.type=[c.type[0]] # AP changed this index from 1 to 0 c.step=['{:0>3}'.format(int(c.step[0]))] c.time=[c.time[0]] @@ -803,7 +833,7 @@ class EIFlexpart: btlist=[13,14,15,16,17,18,19,20,21,22,23,0] - if ((ty.upper() == 'AN' and mod(int(c.time[i]),int(c.dtime))==0) or \ + if ((ty.upper() == 'AN' and mod(int(c.time[i]),int(c.dtime))==0 and int(c.step[i])==0) or \ (ty.upper() != 'AN' and mod(int(c.step[i]),int(c.dtime))==0 and \ mod(int(c.step[i]),int(c.dtime))==0) ) and \ (int(c.time[i]) in btlist or c.maxstep>24): @@ -1205,6 +1235,19 @@ class EIFlexpart: for ofile in self.outputfilelist: p=subprocess.check_call(['grib_set','-s','edition=2,productDefinitionTemplateNumber=8',ofile,ofile+'_2']) p=subprocess.check_call(['mv',ofile+'_2',ofile]) + if c.debug==0: + inputfiles=glob.glob('*.grb') + for grb in inputfiles: + try: + os.remove(grb) + except: + pass + if c.stream=='ELDA': + opposite(self.inputdir+'/'+c.prefix) + for i in range(len(self.outputfilelist)): + if self.outputfilelist[i][-4:]!='N000' : + j=int(self.outputfilelist[i][-3:]) + self.outputfilelist.append(self.outputfilelist[i][:-3]+'{:0>3}'.format(j+25)) if int(c.ectrans)==1 and c.ecapi==False: for ofile in self.outputfilelist: @@ -1275,7 +1318,11 @@ class EIFlexpart: table128=init128(c.flexextractdir+'/grib_templates/ecmwf_grib1_table_128') wrfpars=toparamId('sp/mslp/skt/2t/10u/10v/2d/z/lsm/sst/ci/sd/stl1/stl2/stl3/stl4/swvl1/swvl2/swvl3/swvl4',table128) # index_keys=["date","time","stepRange"] - index_keys=["date","time","step"] + if '/' in c.number: + index_keys=["number","date","time","step"] + else: + index_keys=["date","time","step"] + indexfile=c.inputdir+"/date_time_stepRange.idx" silentremove(indexfile) grib=GribTools(inputfiles.files) @@ -1332,6 +1379,12 @@ class EIFlexpart: if timestamp<slimit or timestamp>elimit: continue + else: + if c.maxstep<24: + if cdateH<c.start_date+'00': + continue + if cdateH>c.end_date+'23': + continue @@ -1431,6 +1484,12 @@ class EIFlexpart: suffix=cdate[2:8]+'.{:0>2}'.format(time/100)+'.{:0>3}'.format(step) else: suffix=cdateH[2:10] + try: + numberindex=index_keys.index('number') + if len(index_vals[numberindex])>1: + suffix=suffix+'.N{:0>3}'.format(int(prod[numberindex])) + except: + pass fnout+=suffix print "outputfile = " + fnout @@ -1459,11 +1518,16 @@ class EIFlexpart: grib_index_release(iid) + return + def deacc_fluxes(self, inputfiles, c): table128=init128(c.flexextractdir+'/grib_templates/ecmwf_grib1_table_128') pars=toparamId(self.params['OG_acc_SL'][0],table128) - index_keys=["date","time","step"] + if '/' in c.number: + index_keys=["number","date","time","step"] + else: + index_keys=["date","time","step"] indexfile=c.inputdir+"/date_time_stepRange.idx" silentremove(indexfile) grib=GribTools(inputfiles.files) @@ -1474,7 +1538,7 @@ class EIFlexpart: index_vals = [] for key in index_keys: key_vals = grib_index_get(iid,key) - + print(key_vals) l=[] for k in key_vals: l.append(int(k)) @@ -1520,16 +1584,24 @@ class EIFlexpart: break fnout=c.inputdir+'/' + numbersuffix='' + try: + numberindex=index_keys.index('number') + if len(index_vals[numberindex])>1: + numbersuffix='.N{:0>3}'.format(int(prod[numberindex])) + except: + pass + if c.maxstep>12: - fnout+='flux'+sdate.strftime('%Y%m%d')+'.{:0>2}'.format(time/100)+'.{:0>3}'.format(step-2*int(c.dtime)) - gnout=c.inputdir+'/flux'+sdate.strftime('%Y%m%d')+'.{:0>2}'.format(time/100)+'.{:0>3}'.format(step-int(c.dtime)) - hnout=c.inputdir+'/flux'+sdate.strftime('%Y%m%d')+'.{:0>2}'.format(time/100)+'.{:0>3}'.format(step) + fnout+='flux'+sdate.strftime('%Y%m%d')+'.{:0>2}'.format(time/100)+'.{:0>3}'.format(step-2*int(c.dtime))+numbersuffix + gnout=c.inputdir+'/flux'+sdate.strftime('%Y%m%d')+'.{:0>2}'.format(time/100)+'.{:0>3}'.format(step-int(c.dtime))+numbersuffix + hnout=c.inputdir+'/flux'+sdate.strftime('%Y%m%d')+'.{:0>2}'.format(time/100)+'.{:0>3}'.format(step)+numbersuffix g=open(gnout,'w') h=open(hnout,'w') else: - fnout+='flux'+fdate.strftime('%Y%m%d%H') - gnout=c.inputdir+'/flux'+(fdate+datetime.timedelta(hours=int(c.dtime))).strftime('%Y%m%d%H') - hnout=c.inputdir+'/flux'+sdates.strftime('%Y%m%d%H') + fnout+='flux'+fdate.strftime('%Y%m%d%H')+numbersuffix + gnout=c.inputdir+'/flux'+(fdate+datetime.timedelta(hours=int(c.dtime))).strftime('%Y%m%d%H')+numbersuffix + hnout=c.inputdir+'/flux'+sdates.strftime('%Y%m%d%H')+numbersuffix g=open(gnout,'w') h=open(hnout,'w') print "outputfile = " + fnout @@ -1555,7 +1627,7 @@ class EIFlexpart: values=(reshape(values,(nj,ni))).flatten()/fak vdp.append(values[:]) # save the accumulated values - if c.marsclass.upper() == 'EA' or \ + if c.marsclass.upper() in ('EA') or \ step<=int(c.dtime): svdp.append(values[:]/int(c.dtime)) else: diff --git a/python/job.temp b/python/job.temp index 2e03859f783ffd1f1cd504d3404f732fc0a4dee6..8afd7c9f3db08e869aba0403efddd8cd13e178fd 100644 --- a/python/job.temp +++ b/python/job.temp @@ -10,7 +10,7 @@ #SBATCH --output=flex_ecmwf.%j.out #SBATCH --error=flex_ecmwf.%j.out #SBATCH --mail-type=FAIL -#SBATCH --time=12:00:00 +#SBATCH --time=24:00:00 ## CRAY specific batch requests ##PBS -N flex_ecmwf diff --git a/python/job.temp.o b/python/job.temp.o index c3209640b214ac32469239c0e6ab7f2ef78ba11f..010404c60aa7647de06678a38a5863578596dc85 100644 --- a/python/job.temp.o +++ b/python/job.temp.o @@ -10,7 +10,7 @@ #SBATCH --output=flex_ecmwf.%j.out #SBATCH --error=flex_ecmwf.%j.out #SBATCH --mail-type=FAIL -#SBATCH --time=12:00:00 +#SBATCH --time=24:00:00 ## CRAY specific batch requests ##PBS -N flex_ecmwf diff --git a/python/opposite.py b/python/opposite.py new file mode 100644 index 0000000000000000000000000000000000000000..979398940adf61e201a7376e0ed51c688a9a02ab --- /dev/null +++ b/python/opposite.py @@ -0,0 +1,53 @@ +import os,sys,glob +from gribapi import * +import numpy + +def opposite(prefix): + fclist=glob.glob(prefix+'*'+'N000') + + for fcn in fclist: + nprefix=fcn.split('N000')[0] + f=open(nprefix+'N{:0>3}'.format(0),'rb') + fcvalues=[] + while 1: + try: + fid=grib_new_from_file(f) + if fid is not None: + fcvalues.append(grib_get_array(fid,'values')) + grib_release(fid) + else: + break + except: + break + for i in range(1,26): + try: + g=open(nprefix+'N{:0>3}'.format(i),'rb') + h=open(nprefix+'N{:0>3}'.format(i+25),'wb') + j=0 + while 1: + try: + gid=grib_new_from_file(g) + if gid is not None: + values=grib_get_array(gid,'values') + grib_set_array(gid,'values',values-2*(values-fcvalues[j])) + grib_set(gid, 'number', i+25) + grib_write(gid,h) + grib_release(gid) + j+=1 + else: + break + except: + break + g.close() + h.close() + print('wrote '+nprefix+'N{:0>3}'.format(i+25)) + except: + break + f.close() + return + +if __name__ == "__main__": + opposite(os.path.expandvars('EN')) + + + diff --git a/python/prepareFLEXPART.py b/python/prepareFLEXPART.py index c8df93d33dfec6870985ea9ea843056c22d97d53..cec924d298f0bf7e1a39fae9b94e6cda24a31791 100755 --- a/python/prepareFLEXPART.py +++ b/python/prepareFLEXPART.py @@ -36,6 +36,7 @@ from UIOTools import UIOFiles #from string import strip from GribTools import GribTools from FlexpartTools import EIFlexpart, Control,interpret_args_and_control, cleanup +from opposite import opposite hostname=socket.gethostname() ecapi= 'ecmwf' not in hostname