From 3f36e423871039ab295f5ea95c56af49ef5067ce Mon Sep 17 00:00:00 2001 From: Anne Philipp <anne.philipp@univie.ac.at> Date: Tue, 4 Dec 2018 22:46:22 +0100 Subject: [PATCH] outsourced some checks from CONTROL class to checks module; translated namelist generation by genshi templating; corrected a bug in grib2 conversion --- source/python/classes/ControlFile.py | 97 +++++------------- source/python/classes/EcFlexpart.py | 132 ++++++++++++++----------- source/python/mods/checks.py | 112 ++++++++++++++++++--- source/python/mods/get_mars_data.py | 3 +- source/python/mods/prepare_flexpart.py | 3 +- 5 files changed, 200 insertions(+), 147 deletions(-) diff --git a/source/python/classes/ControlFile.py b/source/python/classes/ControlFile.py index 6f2db6b..61e7228 100644 --- a/source/python/classes/ControlFile.py +++ b/source/python/classes/ControlFile.py @@ -60,7 +60,7 @@ import inspect sys.path.append('../') import _config from mods.tools import my_error, silent_remove -from mods.checks import check_grid_area +from mods.checks import check_grid, check_area, check_levels # ------------------------------------------------------------------------------ # CLASS @@ -344,27 +344,29 @@ class ControlFile(object): print('Basetime = ' + str(self.basetime)) sys.exit(1) - # assure consistency of levelist and level - if not self.levelist and not self.level: - print('Warning: neither levelist nor level \ - specified in CONTROL file') - sys.exit(1) - elif not self.levelist and self.level: - self.levelist = '1/to/' + self.level - elif (self.levelist and not self.level) or \ - (self.levelist[-1] != self.level[-1]): - self.level = self.levelist.split('/')[-1] - else: - pass + self.levelist, self.level = check_levels(self.levelist, self.level) + + # # assure consistency of levelist and level + # if not self.levelist and not self.level: + # print('Warning: neither levelist nor level \ + # specified in CONTROL file') + # sys.exit(1) + # elif not self.levelist and self.level: + # self.levelist = '1/to/' + self.level + # elif (self.levelist and not self.level) or \ + # (self.levelist[-1] != self.level[-1]): + # self.level = self.levelist.split('/')[-1] + # else: + # pass - # check if max level is a valid level - if int(self.level) not in _config.MAX_LEVEL_LIST: - print('ERROR: ') - print('LEVEL must be the maximum level of a specified ' - 'level list from ECMWF, e.g.') - print(_config.MAX_LEVEL_LIST) - print('Check parameter "LEVEL" or the max level of "LEVELIST"!') - sys.exit(1) + # # check if max level is a valid level + # if int(self.level) not in _config.MAX_LEVEL_LIST: + # print('ERROR: ') + # print('LEVEL must be the maximum level of a specified ' + # 'level list from ECMWF, e.g.') + # print(_config.MAX_LEVEL_LIST) + # print('Check parameter "LEVEL" or the max level of "LEVELIST"!') + # sys.exit(1) # prepare step list if "/" signs are found if '/' in self.step: @@ -476,39 +478,12 @@ class ControlFile(object): self.accmaxstep='12' - self.grid, self.area = check_grid_area(self.grid, self.area, - self.upper, self.lower, - self.left, self.right) + self.grid = check_grid(self.grid) + self.area = check_area(self.grid, self.area, self.upper, self.lower, + self.left, self.right) - # convert grid and area components to correct format and input - #if 'N' in self.grid: # Gaussian output grid - # self.area = 'G' - # else: - # if '/' in self.grid: - # gridx, gridy = self.grid.split('/') - # if gridx == gridy: - # self.grid = gridx - - - # # check on grid format - # if float(self.grid) / 100. >= 0.5: - # # grid is defined in 1/1000 degrees; old format - # self.grid = '{}/{}'.format(float(self.grid) / 1000., - # float(self.grid) / 1000.) - # self.area = '{}/{}/{}/{}'.format(float(self.upper) / 1000., - # float(self.left) / 1000., - # float(self.lower) / 1000., - # float(self.right) / 1000.) - # elif float(self.grid) / 100. < 0.5: - # # grid is defined in normal degree; new format - # self.grid = '{}/{}'.format(float(self.grid), float(self.grid)) - # self.area = '{}/{}/{}/{}'.format(float(self.upper), - # float(self.left), - # float(self.lower), - # float(self.right)) - return def to_list(self): @@ -554,23 +529,3 @@ class ControlFile(object): return sorted(l) - def check_ppid(self, ppid): - '''Sets the current PPID. - - Parameters - ---------- - ppid : :obj:`int` or :obj:`None` - Contains the ppid number provided by the command line parameter - of is None otherwise. - - Return - ------ - - ''' - - if not ppid: - self.ppid = str(os.getppid()) - else: - self.ppid = ppid - - return diff --git a/source/python/classes/EcFlexpart.py b/source/python/classes/EcFlexpart.py index 5da789a..9e872d3 100644 --- a/source/python/classes/EcFlexpart.py +++ b/source/python/classes/EcFlexpart.py @@ -146,12 +146,6 @@ class EcFlexpart(object): self.dtime = c.dtime i = 0 if fluxes and c.maxstep <= 24: - # no forecast beyond one day is needed! - # Thus, prepare flux data manually as usual - # with only forecast fields with start times at 00/12 - # (but without 00/12 fields since these are - # the initialisation times of the flux fields - # and therefore are zero all the time) self.types[str(c.acctype)] = {'times': str(c.acctime), 'steps': '{}/to/{}/by/{}'.format( c.dtime, c.accmaxstep, c.dtime)} @@ -522,8 +516,6 @@ class EcFlexpart(object): retr_param_dict['area'] = self.area retr_param_dict['gaussian'] = self.gaussian - if pk == 'OG__SL': - pass if pk == 'OG_OROLSM__SL' and not oro: oro = True # in CERA20C (class EP) there is no stream "OPER"! @@ -657,50 +649,72 @@ class EcFlexpart(object): from genshi.template.text import NewTextTemplate from genshi.template import TemplateLoader - - loader = TemplateLoader(_config.PATH_TEMPLATES, auto_reload=False) - namelist_template = loader.load(_config.TEMPFILE_NAMELIST, - cls=NewTextTemplate) - - self.inputdir = c.inputdir - area = np.asarray(self.area.split('/')).astype(float) - grid = np.asarray(self.grid.split('/')).astype(float) - - if area[1] > area[3]: - area[1] -= 360 - maxl = int((area[3] - area[1]) / grid[1]) + 1 - maxb = int((area[0] - area[2]) / grid[0]) + 1 - - stream = namelist_template.generate( - maxl = str(maxl), - maxb = str(maxb), - mlevel = str(self.level), - mlevelist = str(self.levelist), - mnauf = str(self.resol), - metapar = '77', - rlo0 = str(area[1]), - rlo1 = str(area[3]), - rla0 = str(area[2]), - rla1 = str(area[0]), - momega = str(c.omega), - momegadiff = str(c.omegadiff), - mgauss = str(c.gauss), - msmooth = str(c.smooth), - meta = str(c.eta), - metadiff = str(c.etadiff), - mdpdeta = str(c.dpdeta) - ) - - namelistfile = os.path.join(self.inputdir, _config.FILE_NAMELIST) - - with open(namelistfile, 'w') as f: - f.write(stream.render('text')) + from genshi.template.eval import UndefinedError + + try: + loader = TemplateLoader(_config.PATH_TEMPLATES, auto_reload=False) + namelist_template = loader.load(_config.TEMPFILE_NAMELIST, + cls=NewTextTemplate) + + self.inputdir = c.inputdir + area = np.asarray(self.area.split('/')).astype(float) + grid = np.asarray(self.grid.split('/')).astype(float) + + if area[1] > area[3]: + area[1] -= 360 + maxl = int((area[3] - area[1]) / grid[1]) + 1 + maxb = int((area[0] - area[2]) / grid[0]) + 1 + + stream = namelist_template.generate( + maxl = str(maxl), + maxb = str(maxb), + mlevel = str(self.level), + mlevelist = str(self.levelist), + mnauf = str(self.resol), + metapar = '77', + rlo0 = str(area[1]), + rlo1 = str(area[3]), + rla0 = str(area[2]), + rla1 = str(area[0]), + momega = str(c.omega), + momegadiff = str(c.omegadiff), + mgauss = str(c.gauss), + msmooth = str(c.smooth), + meta = str(c.eta), + metadiff = str(c.etadiff), + mdpdeta = str(c.dpdeta) + ) + except UndefinedError as e: + print('... ERROR ' + str(e)) + + sys.exit('\n... error occured while trying to generate namelist ' + + _config.TEMPFILE_NAMELIST) + except OSError as e: + print('... ERROR CODE: ' + str(e.errno)) + print('... ERROR MESSAGE:\n \t ' + str(e.strerror)) + + sys.exit('\n... error occured while trying to generate template ' + + _config.TEMPFILE_NAMELIST) + + try: + namelistfile = os.path.join(self.inputdir, _config.FILE_NAMELIST) + + with open(namelistfile, 'w') as f: + f.write(stream.render('text')) + except OSError as e: + print('... ERROR CODE: ' + str(e.errno)) + print('... ERROR MESSAGE:\n \t ' + str(e.strerror)) + + sys.exit('\n... error occured while trying to write ' + + namelistfile) return def deacc_fluxes(self, inputfiles, c): - '''Goes through all flux fields in ordered time and de-accumulate + '''De-accumulate and disaggregate flux data. + + Goes through all flux fields in ordered time and de-accumulate the fields. Afterwards the fields are disaggregated in time. Different versions of disaggregation is provided for rainfall data (darain, modified linear) and the surface fluxes and @@ -987,14 +1001,14 @@ class EcFlexpart(object): # skip the rest of the for loop and start with next timestep/product if not gid: continue - +#============================================================================================ # remove old fort.* files and open new ones # they are just valid for a single product for k, f in fdict.iteritems(): fortfile = os.path.join(c.inputdir, 'fort.' + k) silent_remove(fortfile) fdict[k] = open(fortfile, 'w') - +#============================================================================================ # create correct timestamp from the three time informations cdate = str(codes_get(gid, 'date')) ctime = '{:0>2}'.format(codes_get(gid, 'time')/100) @@ -1018,7 +1032,7 @@ class EcFlexpart(object): fwrf = open(os.path.join(c.outputdir, 'WRF' + cdate + '.' + ctime + '.000.grb2'), 'w') olddate = cdate[:] - +#============================================================================================ # savedfields remembers which fields were already used. savedfields = [] # sum of cloud liquid and ice water content @@ -1082,10 +1096,10 @@ class EcFlexpart(object): codes_release(gid) gid = codes_new_from_index(iid) - +#============================================================================================ for f in fdict.values(): f.close() - +#============================================================================================ # call for Fortran program to convert e.g. reduced_gg grids to # regular_ll and calculate detadot/dp pwd = os.getcwd() @@ -1096,7 +1110,7 @@ class EcFlexpart(object): print('Check parameters CLASS, TYPE, STREAM, START_DATE\n') my_error(c.mailfail, 'fort.21 is empty while parameter eta \ is set to 1 in CONTROL file') - +#============================================================================================ # write out all output to log file before starting fortran programm sys.stdout.flush() @@ -1104,7 +1118,7 @@ class EcFlexpart(object): p = subprocess.check_call([os.path.join( c.exedir, _config.FORTRAN_EXECUTABLE)], shell=True) os.chdir(pwd) - +#============================================================================================ # create name of final output file, e.g. EN13040500 (ENYYMMDDHH) if c.maxstep > 12: suffix = cdate[2:8] + '.' + ctime + '.' + cstep @@ -1114,7 +1128,7 @@ class EcFlexpart(object): print("outputfile = " + fnout) # collect for final processing self.outputfilelist.append(os.path.basename(fnout)) - +#============================================================================================ # create outputfile and copy all data from intermediate files # to the outputfile (final GRIB input files for FLEXPART) orolsm = os.path.basename(glob.glob(c.inputdir + @@ -1134,7 +1148,7 @@ class EcFlexpart(object): with open(os.path.join(c.outputdir, 'OMEGA'), 'wb') as fout: shutil.copyfileobj(open(os.path.join(c.inputdir, 'fort.25'), 'rb'), fout) - +#============================================================================================ if c.wrf: fwrf.close() @@ -1144,7 +1158,9 @@ class EcFlexpart(object): def process_output(self, c): - '''The grib files are postprocessed depending on the selection in + '''Postprocessing of FLEXPART input files. + + The grib files are postprocessed depending on the selection in CONTROL file. The resulting files are moved to the output directory if its not equal to the input directory. The following modifications might be done if @@ -1179,7 +1195,7 @@ class EcFlexpart(object): ofile = os.path.join(self.inputdir, ofile) if c.format.lower() == 'grib2': - p = subprocess.check_call(['codes_set', '-s', 'edition=2,', + p = subprocess.check_call(['grib_set', '-s', 'edition=2,', 'productDefinitionTemplateNumber=8', ofile, ofile + '_2']) p = subprocess.check_call(['mv', ofile + '_2', ofile]) diff --git a/source/python/mods/checks.py b/source/python/mods/checks.py index 0abe1e8..d0f9498 100644 --- a/source/python/mods/checks.py +++ b/source/python/mods/checks.py @@ -25,27 +25,16 @@ # MODULES # ------------------------------------------------------------------------------ - +import _config # ------------------------------------------------------------------------------ # FUNCTIONS # ------------------------------------------------------------------------------ -def check_grid_area(grid, area, upper, lower, left , right): - ''' - - - ''' - # if area was provided - # decompose area into its 4 components - if area: - components = area.split('/') - upper, left, lower, right = components - - if 'N' in grid: # Gaussian output grid - area = 'G' - return grid, area +def check_grid(grid): + if 'N' in grid: + return grid if '/' in grid: gridx, gridy = grid.split('/') if gridx == gridy: @@ -62,6 +51,22 @@ def check_grid_area(grid, area, upper, lower, left , right): # grid is defined in normal degree; new format grid = '{}/{}'.format(float(grid), float(grid)) + return grid + +def check_area(grid, area, upper, lower, left , right): + ''' + + + ''' + if 'N' in grid: # Gaussian output grid + area = 'G' + return area + + # if area was provided decompose area into its 4 components + if area: + components = area.split('/') + upper, left, lower, right = components + # determine area format if (float(upper) / 1000. >= 0.05 and float(lower) / 1000. >= 0.05 and @@ -85,4 +90,79 @@ def check_grid_area(grid, area, upper, lower, left , right): raise ValueError('The area components have different ' 'formats: %s ' (area)) - return grid, area \ No newline at end of file + return area + +def check_levels(levelist, level): + ''' + + Parameters + ---------- + par : :obj:`` + ... + + Return + ------ + + ''' + # assure consistency of levelist and level + if not levelist and not level: + raise ValueError('ERROR: neither levelist nor level ' + 'specified in CONTROL file') + elif not levelist and level: + levelist = '1/to/' + level + elif (levelist and not level) or \ + (levelist[-1] != level[-1]): + level = levelist.split('/')[-1] + else: + pass + + # check if max level is a valid level + if int(level) not in _config.MAX_LEVEL_LIST: + raise ValueError('ERROR: \n' + 'LEVEL must be the maximum level of a specified ' + 'level list from ECMWF, e.g. {} \n' + 'Check parameter "LEVEL" or the max level of ' + '"LEVELIST"!'.format(str(_config.MAX_LEVEL_LIST))) + + return levelist, level + + +def check_ppid(c, ppid): + '''Sets the current PPID. + + Parameters + ---------- + c : :obj:`ControlFile` + Contains all the parameters of CONTROL file and + command line. + + ppid : :obj:`int` or :obj:`None` + Contains the ppid number provided by the command line parameter + of is None otherwise. + + Return + ------ + + ''' + + if not ppid: + c.ppid = str(os.getppid()) + else: + c.ppid = ppid + + return + + +def check_(): + ''' + + Parameters + ---------- + par : :obj:`` + ... + + Return + ------ + + ''' + return diff --git a/source/python/mods/get_mars_data.py b/source/python/mods/get_mars_data.py index c184906..c7c112c 100755 --- a/source/python/mods/get_mars_data.py +++ b/source/python/mods/get_mars_data.py @@ -265,7 +265,8 @@ def mk_dates(c, fluxes): return start, end, chunk def remove_old(pattern, inputdir): - '''Deletes old retrieval files matching the pattern. + '''Deletes old retrieval files from current input directory + matching the pattern. Parameters ---------- diff --git a/source/python/mods/prepare_flexpart.py b/source/python/mods/prepare_flexpart.py index 14e234a..0b95534 100755 --- a/source/python/mods/prepare_flexpart.py +++ b/source/python/mods/prepare_flexpart.py @@ -62,6 +62,7 @@ import socket sys.path.append(os.path.dirname(os.path.abspath( inspect.getfile(inspect.currentframe()))) + '/../') import _config +from checks import check_ppid from classes.UioFiles import UioFiles from classes.ControlFile import ControlFile from tools import clean_up, get_cmdline_arguments, read_ecenv, make_dir @@ -125,7 +126,7 @@ def prepare_flexpart(ppid, c): ------ ''' - c.check_ppid(ppid) + check_ppid(c, ppid) c.ecapi = ecapi -- GitLab