Changeset ff99eae in flex_extract.git


Ignore:
Timestamp:
Jun 1, 2018, 8:34:59 PM (6 years ago)
Author:
Anne Philipp <anne.philipp@…>
Branches:
master, ctbto, dev
Children:
e1228f3
Parents:
ccab809
Message:

completed application of pep8 style guide and pylint investigations. added documentation almost everywhere

Files:
1 added
7 edited
8 moved

Legend:

Unmodified
Added
Removed
  • python/ControlFile.py

    r812283d rff99eae  
    22# -*- coding: utf-8 -*-
    33#************************************************************************
    4 # TODO AP
     4# ToDo AP
    55# - write a test class
    66#************************************************************************
     
    3737#    - __init__
    3838#    - __str__
    39 #    - tolist
     39#    - to_list
     40#
     41# @Class Attributes:
     42#    - start_date
     43#    - end_date
     44#    - accuracy
     45#    - omega
     46#    - cwc
     47#    - omegadiff
     48#    - etadiff
     49#    - level
     50#    - levelist
     51#    - step
     52#    - maxstep
     53#    - prefix
     54#    - makefile
     55#    - basetime
     56#    - date_chunk
     57#    - grib2flexpart
     58#    - exedir
     59#    - flexpart_root_scripts
     60#    - ecmwfdatadir
    4061#
    4162#*******************************************************************************
     
    4667import os
    4768import inspect
     69
    4870# software specific module from flex_extract
    49 import Tools
     71from tools import get_list_as_string, my_error
    5072
    5173# ------------------------------------------------------------------------------
    5274# CLASS
    5375# ------------------------------------------------------------------------------
    54 class ControlFile:
     76class ControlFile(object):
    5577    '''
    5678    Class containing the information of the flex_extract CONTROL file.
     
    114136                        for d in dd:
    115137                            data.append(d)
    116                     pass
    117138                if len(data) == 2:
    118139                    if '$' in data[1]:
     
    126147                                data[1] = data[1][:i] + var + data[1][k+1:]
    127148                            else:
    128                                 Tools.myerror(None,
    129                                               'Could not find variable ' +
    130                                               data[1][j+1:k] +
    131                                               ' while reading ' +
    132                                               filename)
     149                                my_error(None, 'Could not find variable ' +
     150                                         data[1][j+1:k] + ' while reading ' +
     151                                         filename)
    133152                        setattr(self, data[0].lower() + '_expanded', data[1])
    134153                    else:
     
    160179        if not hasattr(self, 'levelist'):
    161180            if not hasattr(self, 'level'):
    162                 print('Warning: neither levelist nor level \
    163                        specified in CONTROL file')
     181                print 'Warning: neither levelist nor level \
     182                       specified in CONTROL file'
    164183            else:
    165184                self.levelist = '1/to/' + self.level
     
    224243        return ', '.join("%s: %s" % item for item in attrs.items())
    225244
    226     def tolist(self):
     245    def to_list(self):
    227246        '''
    228247        @Description:
     
    255274                pass
    256275            else:
    257                 if type(item[1]) is list:
     276                if isinstance(item[1], list):
    258277                    stot = ''
    259278                    for s in item[1]:
     
    265284
    266285        return sorted(l)
     286
     287    # def to_dict(self):
     288        # '''
     289
     290        # '''
     291        # parameters_dict = vars(self)
     292
     293        # # remove unneeded parameter
     294        # parameters_dict.pop('_expanded', None)
     295        # parameters_dict.pop('exedir', None)
     296        # parameters_dict.pop('flexpart_root_scripts', None)
     297        # parameters_dict.pop('ecmwfdatadir', None)
     298
     299        # parameters_dict_str = {}
     300        # for key, value in parameters_dict.iteritems():
     301            # if isinstance(value, list):
     302                # parameters_dict_str[str(key)] = get_list_as_string(value, ' ')
     303            # else:
     304                # parameters_dict_str[str(key)] = str(value)
     305
     306        # return parameters_dict_str
  • python/EcFlexpart.py

    r991df6a rff99eae  
    22# -*- coding: utf-8 -*-
    33#************************************************************************
    4 # TODO AP
     4# ToDo AP
    55# - specifiy file header documentation
    66# - add class description in header information
     
    2121#        - extended with class Control
    2222#        - removed functions mkdir_p, daterange, years_between, months_between
    23 #        - added functions darain, dapoly, toparamId, init128, normalexit,
    24 #          myerror, cleanup, install_args_and_control,
     23#        - added functions darain, dapoly, to_param_id, init128, normal_exit,
     24#          my_error, clean_up, install_args_and_control,
    2525#          interpret_args_and_control,
    2626#        - removed function __del__ in class EIFLexpart
     
    4040#        - applied PEP8 style guide
    4141#        - added documentation
    42 #        - removed function getFlexpartTime in class ECFlexpart
     42#        - removed function getFlexpartTime in class EcFlexpart
    4343#        - outsourced class ControlFile
    4444#        - outsourced class MarsRetrieval
    45 #        - changed class name from EIFlexpart to ECFlexpart
     45#        - changed class name from EIFlexpart to EcFlexpart
    4646#        - applied minor code changes (style)
    4747#        - removed "dead code" , e.g. retrieval of Q since it is not needed
     
    7070#    - deacc_fluxes
    7171#
     72# @Class Attributes:
     73#    - dtime
     74#    - basetime
     75#    - server
     76#    - marsclass
     77#    - stream
     78#    - resol
     79#    - accuracy
     80#    - number
     81#    - expver
     82#    - glevelist
     83#    - area
     84#    - grid
     85#    - level
     86#    - levelist
     87#    - types
     88#    - dates
     89#    - area
     90#    - gaussian
     91#    - params
     92#    - inputdir
     93#    - outputfilelist
     94#
    7295#*******************************************************************************
    73 
     96#pylint: disable=unsupported-assignment-operation
     97# this is disabled because its an error in pylint for this specific case
     98#pylint: disable=consider-using-enumerate
     99# this is not useful in this case
    74100# ------------------------------------------------------------------------------
    75101# MODULES
     
    78104import shutil
    79105import os
    80 import sys
    81 import inspect
    82106import glob
    83 import datetime
    84 from numpy import *
    85 from argparse import ArgumentParser, ArgumentDefaultsHelpFormatter
    86 ecapi = True
    87 try:
    88     import ecmwfapi
    89 except ImportError:
    90     ecapi = False
    91 from gribapi import *
     107from datetime import datetime, timedelta
     108import numpy as np
     109from gribapi import grib_set, grib_index_select, grib_new_from_index, grib_get,\
     110                    grib_write, grib_get_values, grib_set_values, grib_release,\
     111                    grib_index_release, grib_index_get
    92112
    93113# software specific classes and modules from flex_extract
    94 from GribTools import GribTools
    95 from Tools import init128, toparamId, silentremove, product
    96 from ControlFile import ControlFile
    97 from MARSretrieval import MARSretrieval
    98 import Disagg
     114from Gribtools import Gribtools
     115from tools import init128, to_param_id, silent_remove, product, my_error
     116from MarsRetrieval import MarsRetrieval
     117import disaggregation
    99118
    100119# ------------------------------------------------------------------------------
    101120# CLASS
    102121# ------------------------------------------------------------------------------
    103 class ECFlexpart:
     122class EcFlexpart(object):
    104123    '''
    105124    Class to retrieve FLEXPART specific ECMWF data.
     
    111130        '''
    112131        @Description:
    113             Creates an object/instance of ECFlexpart with the
     132            Creates an object/instance of EcFlexpart with the
    114133            associated settings of its attributes for the retrieval.
    115134
    116135        @Input:
    117             self: instance of ECFlexpart
     136            self: instance of EcFlexpart
    118137                The current object of the class.
    119138
     
    142161        # different mars types for retrieving data for flexpart
    143162        self.types = dict()
    144         try:
    145             if c.maxstep > len(c.type):    # Pure forecast mode
    146                 c.type = [c.type[1]]
    147                 c.step = ['{:0>3}'.format(int(c.step[0]))]
    148                 c.time = [c.time[0]]
    149                 for i in range(1, c.maxstep + 1):
    150                     c.type.append(c.type[0])
    151                     c.step.append('{:0>3}'.format(i))
    152                     c.time.append(c.time[0])
    153         except:
    154             pass
     163
     164        if c.maxstep > len(c.type):    # Pure forecast mode
     165            c.type = [c.type[1]]
     166            c.step = ['{:0>3}'.format(int(c.step[0]))]
     167            c.time = [c.time[0]]
     168            for i in range(1, c.maxstep + 1):
     169                c.type.append(c.type[0])
     170                c.step.append('{:0>3}'.format(i))
     171                c.time.append(c.time[0])
    155172
    156173        self.inputdir = c.inputdir
     
    175192                    btlist = [13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 0]
    176193
    177                 if mod(i, int(c.dtime)) == 0 and \
    178                     (c.maxstep > 24 or i in btlist):
     194                if i % int(c.dtime) == 0 and c.maxstep > 24 or i in btlist:
    179195
    180196                    if ty not in self.types.keys():
     
    182198
    183199                    if ti not in self.types[ty]['times']:
    184                         if len(self.types[ty]['times']) > 0:
     200                        if self.types[ty]['times']:
    185201                            self.types[ty]['times'] += '/'
    186202                        self.types[ty]['times'] += ti
    187203
    188204                    if st not in self.types[ty]['steps']:
    189                         if len(self.types[ty]['steps']) > 0:
     205                        if self.types[ty]['steps']:
    190206                            self.types[ty]['steps'] += '/'
    191207                        self.types[ty]['steps'] += st
    192208                i += 1
    193209
    194         # Different grids need different retrievals
    195         # SH = Spherical Harmonics, GG = Gaussian Grid,
    196         # OG = Output Grid, ML = MultiLevel, SL = SingleLevel
    197         self.params = {'SH__ML': '', 'SH__SL': '',
    198                        'GG__ML': '', 'GG__SL': '',
    199                        'OG__ML': '', 'OG__SL': '',
    200                        'OG_OROLSM_SL': '', 'OG_acc_SL': ''}
     210
    201211        self.marsclass = c.marsclass
    202212        self.stream = c.stream
     
    205215        self.accuracy = c.accuracy
    206216        self.level = c.level
    207         try:
     217
     218        if c.levelist:
    208219            self.levelist = c.levelist
    209         except:
     220        else:
    210221            self.levelist = '1/to/' + c.level
    211222
     
    213224        self.glevelist = '1/to/' + c.level
    214225
    215         try:
     226        if c.gaussian:
    216227            self.gaussian = c.gaussian
    217         except:
     228        else:
    218229            self.gaussian = ''
    219230
    220         try:
     231        if c.expver:
    221232            self.expver = c.expver
    222         except:
     233        else:
    223234            self.expver = '1'
    224235
    225         try:
     236        if c.number:
    226237            self.number = c.number
    227         except:
     238        else:
    228239            self.number = '0'
    229240
     
    250261        # 4) Download also data for WRF
    251262
     263
     264        # Different grids need different retrievals
     265        # SH = Spherical Harmonics, GG = Gaussian Grid,
     266        # OG = Output Grid, ML = MultiLevel, SL = SingleLevel
     267        self.params = {'SH__ML': '', 'SH__SL': '',
     268                       'GG__ML': '', 'GG__SL': '',
     269                       'OG__ML': '', 'OG__SL': '',
     270                       'OG_OROLSM_SL': '', 'OG_acc_SL': ''}
     271
    252272        if fluxes is False:
    253273            self.params['SH__SL'] = ['LNSP', 'ML', '1', 'OFF']
     
    255275            self.params['OG__SL'] = ["141/151/164/165/166/167/168/129/172", \
    256276                                     'SFC', '1', self.grid]
    257             if len(c.addpar) > 0:
     277            if c.addpar:
    258278                if c.addpar[0] == '/':
    259279                    c.addpar = c.addpar[1:]
     
    277297                self.params['SH__ML'] = ['U/V/D', 'ML', self.glevelist, 'OFF']
    278298            else:
    279                 print('Warning: This is a very costly parameter combination, \
    280                        use only for debugging!')
     299                print 'Warning: This is a very costly parameter combination, \
     300                       use only for debugging!'
    281301                self.params['GG__SL'] = ['Q', 'ML', '1', \
    282302                                         '{}'.format((int(self.resol) + 1) / 2)]
     
    287307                self.params['OG__ML'][0] += '/W'
    288308
    289             try:
    290                 # add cloud water content if necessary
    291                 if c.cwc == '1':
    292                     self.params['OG__ML'][0] += '/CLWC/CIWC'
    293             except:
    294                 pass
    295 
    296             try:
    297                 # add vorticity and geopotential height for WRF if necessary
    298                 if c.wrf == '1':
    299                     self.params['OG__ML'][0] += '/Z/VO'
    300                     if '/D' not in self.params['OG__ML'][0]:
    301                         self.params['OG__ML'][0] += '/D'
    302                     #wrf_sfc = 'sp/msl/skt/2t/10u/10v/2d/z/lsm/sst/ci/sd/stl1/ /
    303                     #           stl2/stl3/stl4/swvl1/swvl2/swvl3/swvl4'.upper()
    304                     wrf_sfc = '134/235/167/165/166/168/129/172/34/31/141/ \
    305                                139/170/183/236/39/40/41/42'.upper()
    306                     lwrt_sfc = wrf_sfc.split('/')
    307                     for par in lwrt_sfc:
    308                         if par not in self.params['OG__SL'][0]:
    309                             self.params['OG__SL'][0] += '/' + par
    310             except:
    311                 pass
     309            # add cloud water content if necessary
     310            if c.cwc == '1':
     311                self.params['OG__ML'][0] += '/CLWC/CIWC'
     312
     313            # add vorticity and geopotential height for WRF if necessary
     314            if c.wrf == '1':
     315                self.params['OG__ML'][0] += '/Z/VO'
     316                if '/D' not in self.params['OG__ML'][0]:
     317                    self.params['OG__ML'][0] += '/D'
     318                #wrf_sfc = 'sp/msl/skt/2t/10u/10v/2d/z/lsm/sst/ci/sd/stl1/ /
     319                #           stl2/stl3/stl4/swvl1/swvl2/swvl3/swvl4'.upper()
     320                wrf_sfc = '134/235/167/165/166/168/129/172/34/31/141/ \
     321                           139/170/183/236/39/40/41/42'.upper()
     322                lwrt_sfc = wrf_sfc.split('/')
     323                for par in lwrt_sfc:
     324                    if par not in self.params['OG__SL'][0]:
     325                        self.params['OG__SL'][0] += '/' + par
     326
    312327        else:
    313328            self.params['OG_acc_SL'] = ["LSP/CP/SSHF/EWSS/NSSS/SSR", \
     
    328343
    329344        @Input:
    330             self: instance of ECFlexpart
     345            self: instance of EcFlexpart
    331346                The current object of the class.
    332347
     
    352367
    353368        self.inputdir = c.inputdir
    354         area = asarray(self.area.split('/')).astype(float)
    355         grid = asarray(self.grid.split('/')).astype(float)
     369        area = np.asarray(self.area.split('/')).astype(float)
     370        grid = np.asarray(self.grid.split('/')).astype(float)
    356371
    357372        if area[1] > area[3]:
    358373            area[1] -= 360
    359         zyk = abs((area[3] - area[1] - 360.) + grid[1]) < 1.e-6
    360374        maxl = int((area[3] - area[1]) / grid[1]) + 1
    361375        maxb = int((area[0] - area[2]) / grid[0]) + 1
     
    364378            f.write('&NAMGEN\n')
    365379            f.write(',\n  '.join(['maxl = ' + str(maxl), 'maxb = ' + str(maxb),
    366                     'mlevel = ' + self.level,
    367                     'mlevelist = ' + '"' + self.levelist + '"',
    368                     'mnauf = ' + self.resol, 'metapar = ' + '77',
    369                     'rlo0 = ' + str(area[1]), 'rlo1 = ' + str(area[3]),
    370                     'rla0 = ' + str(area[2]), 'rla1 = ' + str(area[0]),
    371                     'momega = ' + c.omega, 'momegadiff = ' + c.omegadiff,
    372                     'mgauss = ' + c.gauss, 'msmooth = ' + c.smooth,
    373                     'meta = ' + c.eta, 'metadiff = ' + c.etadiff,
    374                     'mdpdeta = ' + c.dpdeta]))
     380                                  'mlevel = ' + self.level,
     381                                  'mlevelist = ' + '"' + self.levelist + '"',
     382                                  'mnauf = ' + self.resol,
     383                                  'metapar = ' + '77',
     384                                  'rlo0 = ' + str(area[1]),
     385                                  'rlo1 = ' + str(area[3]),
     386                                  'rla0 = ' + str(area[2]),
     387                                  'rla1 = ' + str(area[0]),
     388                                  'momega = ' + c.omega,
     389                                  'momegadiff = ' + c.omegadiff,
     390                                  'mgauss = ' + c.gauss,
     391                                  'msmooth = ' + c.smooth,
     392                                  'meta = ' + c.eta,
     393                                  'metadiff = ' + c.etadiff,
     394                                  'mdpdeta = ' + c.dpdeta]))
    375395
    376396            f.write('\n/\n')
     
    386406
    387407        @Input:
    388             self: instance of ECFlexpart
     408            self: instance of EcFlexpart
    389409                The current object of the class.
    390410
     
    449469    # ------  on demand path  --------------------------------------------------
    450470                if self.basetime is None:
    451                     MR = MARSretrieval(self.server,
    452                             marsclass=self.marsclass, stream=mfstream,
    453                             type=mftype, levtype=pv[1], levelist=pv[2],
    454                             resol=self.resol, gaussian=gaussian,
    455                             accuracy=self.accuracy, grid=pv[3],
    456                             target=mftarget, area=area, date=mfdate,
    457                             time=mftime, number=self.number, step=mfstep,
    458                             expver=self.expver, param=pv[0])
    459 
    460                     MR.displayInfo()
    461                     MR.dataRetrieve()
     471                    MR = MarsRetrieval(self.server,
     472                                       marsclass=self.marsclass,
     473                                       stream=mfstream,
     474                                       type=mftype,
     475                                       levtype=pv[1],
     476                                       levelist=pv[2],
     477                                       resol=self.resol,
     478                                       gaussian=gaussian,
     479                                       accuracy=self.accuracy,
     480                                       grid=pv[3],
     481                                       target=mftarget,
     482                                       area=area,
     483                                       date=mfdate,
     484                                       time=mftime,
     485                                       number=self.number,
     486                                       step=mfstep,
     487                                       expver=self.expver,
     488                                       param=pv[0])
     489
     490                    MR.display_info()
     491                    MR.data_retrieve()
    462492    # ------  operational path  ------------------------------------------------
    463493                else:
     
    475505                        tm1 = -1
    476506
    477                     maxtime = datetime.datetime.strptime(
    478                                 mfdate.split('/')[-1] + mftime.split('/')[tm1],
    479                                 '%Y%m%d%H') + datetime.timedelta(
    480                                 hours=int(mfstep.split('/')[sm1]))
    481                     elimit = datetime.datetime.strptime(
    482                                 mfdate.split('/')[-1] +
    483                                 self.basetime, '%Y%m%d%H')
     507                    maxdate = datetime.strptime(mfdate.split('/')[-1] +
     508                                                mftime.split('/')[tm1],
     509                                                '%Y%m%d%H')
     510                    istep = int(mfstep.split('/')[sm1])
     511                    maxtime = maxdate + timedelta(hours=istep)
     512
     513                    elimit = datetime.strptime(mfdate.split('/')[-1] +
     514                                               self.basetime, '%Y%m%d%H')
    484515
    485516                    if self.basetime == '12':
     
    492523                        # if maxtime-elimit<12h reduce step for last time
    493524                        # A split of the MARS job into 2 is likely necessary.
    494                             maxtime = elimit - datetime.timedelta(hours=24)
    495                             mfdate = '/'.join(('/'.join(mfdate.split('/')[:-1]),
    496                                                 datetime.datetime.strftime(
    497                                                 maxtime, '%Y%m%d')))
    498 
    499                             MR = MARSretrieval(self.server,
    500                                             marsclass=self.marsclass,
    501                                             stream=self.stream, type=mftype,
    502                                             levtype=pv[1], levelist=pv[2],
    503                                             resol=self.resol, gaussian=gaussian,
    504                                             accuracy=self.accuracy, grid=pv[3],
    505                                             target=mftarget, area=area,
    506                                             date=mfdate, time=mftime,
    507                                             number=self.number, step=mfstep,
    508                                             expver=self.expver, param=pv[0])
    509 
    510                             MR.displayInfo()
    511                             MR.dataRetrieve()
    512 
    513                             maxtime = elimit - datetime.timedelta(hours=12)
    514                             mfdate = datetime.datetime.strftime(maxtime,
    515                                                                 '%Y%m%d')
     525                            maxtime = elimit - timedelta(hours=24)
     526                            mfdate = '/'.join(['/'.join(mfdate.split('/')[:-1]),
     527                                               datetime.strftime(maxtime,
     528                                                                 '%Y%m%d')])
     529
     530                            MR = MarsRetrieval(self.server,
     531                                               marsclass=self.marsclass,
     532                                               stream=self.stream,
     533                                               type=mftype,
     534                                               levtype=pv[1],
     535                                               levelist=pv[2],
     536                                               resol=self.resol,
     537                                               gaussian=gaussian,
     538                                               accuracy=self.accuracy,
     539                                               grid=pv[3],
     540                                               target=mftarget,
     541                                               area=area,
     542                                               date=mfdate,
     543                                               time=mftime,
     544                                               number=self.number,
     545                                               step=mfstep,
     546                                               expver=self.expver,
     547                                               param=pv[0])
     548
     549                            MR.display_info()
     550                            MR.data_retrieve()
     551
     552                            maxtime = elimit - timedelta(hours=12)
     553                            mfdate = datetime.strftime(maxtime, '%Y%m%d')
    516554                            mftime = '00'
    517555                            mftarget = self.inputdir + "/" + ftype + pk + \
     
    519557                                       '.' + str(os.getpid()) + ".grb"
    520558
    521                             MR = MARSretrieval(self.server,
    522                                             marsclass=self.marsclass,
    523                                             stream=self.stream, type=mftype,
    524                                             levtype=pv[1], levelist=pv[2],
    525                                             resol=self.resol, gaussian=gaussian,
    526                                             accuracy=self.accuracy, grid=pv[3],
    527                                             target=mftarget, area=area,
    528                                             date=mfdate, time=mftime,
    529                                             number=self.number, step=mfstep,
    530                                             expver=self.expver, param=pv[0])
    531 
    532                             MR.displayInfo()
    533                             MR.dataRetrieve()
     559                            MR = MarsRetrieval(self.server,
     560                                               marsclass=self.marsclass,
     561                                               stream=self.stream,
     562                                               type=mftype,
     563                                               levtype=pv[1],
     564                                               levelist=pv[2],
     565                                               resol=self.resol,
     566                                               gaussian=gaussian,
     567                                               accuracy=self.accuracy,
     568                                               grid=pv[3],
     569                                               target=mftarget,
     570                                               area=area,
     571                                               date=mfdate,
     572                                               time=mftime,
     573                                               number=self.number,
     574                                               step=mfstep,
     575                                               expver=self.expver,
     576                                               param=pv[0])
     577
     578                            MR.display_info()
     579                            MR.data_retrieve()
    534580                        # --------------  non flux data ------------------------
    535581                        else:
    536                             MR = MARSretrieval(self.server,
    537                                             marsclass=self.marsclass,
    538                                             stream=self.stream, type=mftype,
    539                                             levtype=pv[1], levelist=pv[2],
    540                                             resol=self.resol, gaussian=gaussian,
    541                                             accuracy=self.accuracy, grid=pv[3],
    542                                             target=mftarget, area=area,
    543                                             date=mfdate, time=mftime,
    544                                             number=self.number, step=mfstep,
    545                                             expver=self.expver, param=pv[0])
    546 
    547                             MR.displayInfo()
    548                             MR.dataRetrieve()
     582                            MR = MarsRetrieval(self.server,
     583                                               marsclass=self.marsclass,
     584                                               stream=self.stream,
     585                                               type=mftype,
     586                                               levtype=pv[1],
     587                                               levelist=pv[2],
     588                                               resol=self.resol,
     589                                               gaussian=gaussian,
     590                                               accuracy=self.accuracy,
     591                                               grid=pv[3],
     592                                               target=mftarget,
     593                                               area=area,
     594                                               date=mfdate,
     595                                               time=mftime,
     596                                               number=self.number,
     597                                               step=mfstep,
     598                                               expver=self.expver,
     599                                               param=pv[0])
     600
     601                            MR.display_info()
     602                            MR.data_retrieve()
    549603                    else: # basetime == 0 ??? #AP
    550604
    551                         maxtime = elimit - datetime.timedelta(hours=24)
    552                         mfdate = datetime.datetime.strftime(maxtime,'%Y%m%d')
     605                        maxtime = elimit - timedelta(hours=24)
     606                        mfdate = datetime.strftime(maxtime, '%Y%m%d')
    553607                        mftimesave = ''.join(mftime)
    554608
     
    556610                            times = mftime.split('/')
    557611                            while ((int(times[0]) +
    558                                    int(mfstep.split('/')[0]) <= 12) and
    559                                   (pk != 'OG_OROLSM__SL') and 'acc' not in pk):
     612                                    int(mfstep.split('/')[0]) <= 12) and
     613                                   (pk != 'OG_OROLSM__SL') and 'acc' not in pk):
    560614                                times = times[1:]
    561615                            if len(times) > 1:
     
    564618                                mftime = times[0]
    565619
    566                         MR = MARSretrieval(self.server,
    567                                         marsclass=self.marsclass,
    568                                         stream=self.stream, type=mftype,
    569                                         levtype=pv[1], levelist=pv[2],
    570                                         resol=self.resol, gaussian=gaussian,
    571                                         accuracy=self.accuracy, grid=pv[3],
    572                                         target=mftarget, area=area,
    573                                         date=mfdate, time=mftime,
    574                                         number=self.number, step=mfstep,
    575                                         expver=self.expver, param=pv[0])
    576 
    577                         MR.displayInfo()
    578                         MR.dataRetrieve()
     620                        MR = MarsRetrieval(self.server,
     621                                           marsclass=self.marsclass,
     622                                           stream=self.stream,
     623                                           type=mftype,
     624                                           levtype=pv[1],
     625                                           levelist=pv[2],
     626                                           resol=self.resol,
     627                                           gaussian=gaussian,
     628                                           accuracy=self.accuracy,
     629                                           grid=pv[3],
     630                                           target=mftarget,
     631                                           area=area,
     632                                           date=mfdate,
     633                                           time=mftime,
     634                                           number=self.number,
     635                                           step=mfstep,
     636                                           expver=self.expver,
     637                                           param=pv[0])
     638
     639                        MR.display_info()
     640                        MR.data_retrieve()
    579641
    580642                        if (int(mftimesave.split('/')[0]) == 0 and
    581                             int(mfstep.split('/')[0]) == 0 and
    582                             pk != 'OG_OROLSM__SL'):
    583                             mfdate = datetime.datetime.strftime(elimit,'%Y%m%d')
     643                                int(mfstep.split('/')[0]) == 0 and
     644                                pk != 'OG_OROLSM__SL'):
     645
     646                            mfdate = datetime.strftime(elimit, '%Y%m%d')
    584647                            mftime = '00'
    585648                            mfstep = '000'
     
    588651                                       '.' + str(os.getpid()) + ".grb"
    589652
    590                             MR = MARSretrieval(self.server,
    591                                         marsclass=self.marsclass,
    592                                         stream=self.stream, type=mftype,
    593                                         levtype=pv[1], levelist=pv[2],
    594                                         resol=self.resol, gaussian=gaussian,
    595                                         accuracy=self.accuracy, grid=pv[3],
    596                                         target=mftarget, area=area,
    597                                         date=mfdate, time=mftime,
    598                                         number=self.number, step=mfstep,
    599                                         expver=self.expver, param=pv[0])
    600 
    601                             MR.displayInfo()
    602                             MR.dataRetrieve()
    603 
    604         print("MARS retrieve done... ")
     653                            MR = MarsRetrieval(self.server,
     654                                               marsclass=self.marsclass,
     655                                               stream=self.stream,
     656                                               type=mftype,
     657                                               levtype=pv[1],
     658                                               levelist=pv[2],
     659                                               resol=self.resol,
     660                                               gaussian=gaussian,
     661                                               accuracy=self.accuracy,
     662                                               grid=pv[3],
     663                                               target=mftarget,
     664                                               area=area,
     665                                               date=mfdate,
     666                                               time=mftime,
     667                                               number=self.number,
     668                                               step=mfstep,
     669                                               expver=self.expver,
     670                                               param=pv[0])
     671
     672                            MR.display_info()
     673                            MR.data_retrieve()
     674
     675        print "MARS retrieve done... "
    605676
    606677        return
     
    621692
    622693        @Input:
    623             self: instance of ECFlexpart
     694            self: instance of EcFlexpart
    624695                The current object of the class.
    625696
     
    642713        '''
    643714
    644         print('\n\nPostprocessing:\n Format: {}\n'.format(c.format))
     715        print '\n\nPostprocessing:\n Format: {}\n'.format(c.format)
    645716
    646717        if c.ecapi is False:
     
    652723                c.destination = os.getenv('DESTINATION')
    653724            print('ectrans: {}\n gateway: {}\n destination: {}\n '
    654                     .format(c.ectrans, c.gateway, c.destination))
    655 
    656         print('Output filelist: \n')
    657         print(self.outputfilelist)
     725                  .format(c.ectrans, c.gateway, c.destination))
     726
     727        print 'Output filelist: \n'
     728        print self.outputfilelist
    658729
    659730        if c.format.lower() == 'grib2':
    660731            for ofile in self.outputfilelist:
    661732                p = subprocess.check_call(['grib_set', '-s', 'edition=2, \
    662                                             productDefinitionTemplateNumber=8',
    663                                             ofile, ofile + '_2'])
     733                                           productDefinitionTemplateNumber=8',
     734                                           ofile, ofile + '_2'])
    664735                p = subprocess.check_call(['mv', ofile + '_2', ofile])
    665736
     
    667738            for ofile in self.outputfilelist:
    668739                p = subprocess.check_call(['ectrans', '-overwrite', '-gateway',
    669                                             c.gateway, '-remote', c.destination,
    670                                             '-source', ofile])
     740                                           c.gateway, '-remote', c.destination,
     741                                           '-source', ofile])
    671742                print('ectrans:', p)
    672743
     
    674745            for ofile in self.outputfilelist:
    675746                p = subprocess.check_call(['ecp', '-o', ofile,
    676                                             os.path.expandvars(c.ecfsdir)])
     747                                           os.path.expandvars(c.ecfsdir)])
    677748
    678749        if c.outputdir != c.inputdir:
     
    692763                if '.' in fname[-1]:
    693764                    l = fname[-1].split('.')
    694                     timestamp = datetime.datetime.strptime(l[0][-6:] + l[1],
    695                                                            '%y%m%d%H')
    696                     timestamp += datetime.timedelta(hours=int(l[2]))
    697                     cdate = datetime.datetime.strftime(timestamp, '%Y%m%d')
    698                     chms = datetime.datetime.strftime(timestamp, '%H%M%S')
     765                    timestamp = datetime.strptime(l[0][-6:] + l[1],
     766                                                  '%y%m%d%H')
     767                    timestamp += timedelta(hours=int(l[2]))
     768                    cdate = datetime.strftime(timestamp, '%Y%m%d')
     769                    chms = datetime.strftime(timestamp, '%H%M%S')
    699770                else:
    700771                    cdate = '20' + fname[-1][-8:-2]
     
    708779            # generate pathnames file
    709780            pwd = os.path.abspath(c.outputdir)
    710             with open(pwd + '/pathnames','w') as f:
     781            with open(pwd + '/pathnames', 'w') as f:
    711782                f.write(pwd + '/Options/\n')
    712783                f.write(pwd + '/\n')
     
    720791
    721792            # read template COMMAND file
    722             with open(os.path.expandvars(
    723                      os.path.expanduser(c.flexpart_root_scripts)) +
    724                      '/../Options/COMMAND', 'r') as f:
     793            with open(os.path.expandvars(os.path.expanduser(
     794                c.flexpart_root_scripts)) + '/../Options/COMMAND', 'r') as f:
    725795                lflist = f.read().split('\n')
    726796
     
    747817            # afterwards switch back to the working dir
    748818            os.chdir(c.outputdir)
    749             p = subprocess.check_call([os.path.expandvars(
    750                         os.path.expanduser(c.flexpart_root_scripts)) +
    751                         '/../FLEXPART_PROGRAM/grib2flexpart',
    752                         'useAvailable', '.'])
     819            p = subprocess.check_call([
     820                os.path.expandvars(os.path.expanduser(c.flexpart_root_scripts))
     821                + '/../FLEXPART_PROGRAM/grib2flexpart', 'useAvailable', '.'])
    753822            os.chdir(pwd)
    754823
     
    771840
    772841        @Input:
    773             self: instance of ECFlexpart
     842            self: instance of EcFlexpart
    774843                The current object of the class.
    775844
    776             inputfiles: instance of UIOFiles
     845            inputfiles: instance of UioFiles
    777846                Contains a list of files.
    778847
     
    796865        table128 = init128(c.ecmwfdatadir +
    797866                           '/grib_templates/ecmwf_grib1_table_128')
    798         wrfpars = toparamId('sp/mslp/skt/2t/10u/10v/2d/z/lsm/sst/ci/sd/\
     867        wrfpars = to_param_id('sp/mslp/skt/2t/10u/10v/2d/z/lsm/sst/ci/sd/\
    799868                            stl1/stl2/stl3/stl4/swvl1/swvl2/swvl3/swvl4',
    800869                            table128)
     
    802871        index_keys = ["date", "time", "step"]
    803872        indexfile = c.inputdir + "/date_time_stepRange.idx"
    804         silentremove(indexfile)
    805         grib = GribTools(inputfiles.files)
     873        silent_remove(indexfile)
     874        grib = Gribtools(inputfiles.files)
    806875        # creates new index file
    807876        iid = grib.index(index_keys=index_keys, index_file=indexfile)
     
    811880        for key in index_keys:
    812881            index_vals.append(grib_index_get(iid, key))
    813             print(index_vals[-1])
     882            print index_vals[-1]
    814883            # index_vals looks for example like:
    815884            # index_vals[0]: ('20171106', '20171107', '20171108') ; date
     
    842911                # remove old fort.* files and open new ones
    843912                for k, f in fdict.iteritems():
    844                     silentremove(c.inputdir + "/fort." + k)
     913                    silent_remove(c.inputdir + "/fort." + k)
    845914                    fdict[k] = open(c.inputdir + '/fort.' + k, 'w')
    846915
    847916                cdate = str(grib_get(gid, 'date'))
    848917                time = grib_get(gid, 'time')
    849                 type = grib_get(gid, 'type')
    850918                step = grib_get(gid, 'step')
    851919                # create correct timestamp from the three time informations
    852920                # date, time, step
    853                 timestamp = datetime.datetime.strptime(
    854                                 cdate + '{:0>2}'.format(time/100), '%Y%m%d%H')
    855                 timestamp += datetime.timedelta(hours=int(step))
    856 
    857                 cdateH = datetime.datetime.strftime(timestamp, '%Y%m%d%H')
    858                 chms = datetime.datetime.strftime(timestamp, '%H%M%S')
     921                timestamp = datetime.strptime(cdate + '{:0>2}'.format(time/100),
     922                                              '%Y%m%d%H')
     923                timestamp += timedelta(hours=int(step))
     924
     925                cdateH = datetime.strftime(timestamp, '%Y%m%d%H')
    859926
    860927                if c.basetime is not None:
    861                     slimit = datetime.datetime.strptime(
    862                                 c.start_date + '00', '%Y%m%d%H')
     928                    slimit = datetime.strptime(c.start_date + '00', '%Y%m%d%H')
    863929                    bt = '23'
    864930                    if c.basetime == '00':
    865931                        bt = '00'
    866                         slimit = datetime.datetime.strptime(
    867                                     c.end_date + bt, '%Y%m%d%H') - \
    868                                     datetime.timedelta(hours=12-int(c.dtime))
     932                        slimit = datetime.strptime(c.end_date + bt, '%Y%m%d%H')\
     933                            - timedelta(hours=12-int(c.dtime))
    869934                    if c.basetime == '12':
    870935                        bt = '12'
    871                         slimit = datetime.datetime.strptime(
    872                                     c.end_date + bt, '%Y%m%d%H') - \
    873                                  datetime.timedelta(hours=12-int(c.dtime))
    874 
    875                     elimit = datetime.datetime.strptime(
    876                                 c.end_date + bt, '%Y%m%d%H')
     936                        slimit = datetime.strptime(c.end_date + bt, '%Y%m%d%H')\
     937                            - timedelta(hours=12-int(c.dtime))
     938
     939                    elimit = datetime.strptime(c.end_date + bt, '%Y%m%d%H')
    877940
    878941                    if timestamp < slimit or timestamp > elimit:
     
    901964                paramId = grib_get(gid, 'paramId')
    902965                gridtype = grib_get(gid, 'gridType')
    903                 datatype = grib_get(gid, 'dataType')
    904966                levtype = grib_get(gid, 'typeOfLevel')
    905967                if paramId == 133 and gridtype == 'reduced_gg':
     
    9391001                        savedfields.append(paramId)
    9401002                    else:
    941                         print('duplicate ' + str(paramId) + ' not written')
     1003                        print 'duplicate ' + str(paramId) + ' not written'
    9421004
    9431005                try:
     
    9631025                os.chdir(c.inputdir)
    9641026                if os.stat('fort.21').st_size == 0 and int(c.eta) == 1:
    965                     print('Parameter 77 (etadot) is missing, most likely it is \
    966                            not available for this type or date/time\n')
    967                     print('Check parameters CLASS, TYPE, STREAM, START_DATE\n')
    968                     myerror(c, 'fort.21 is empty while parameter eta is set \
     1027                    print 'Parameter 77 (etadot) is missing, most likely it is \
     1028                           not available for this type or date/time\n'
     1029                    print 'Check parameters CLASS, TYPE, STREAM, START_DATE\n'
     1030                    my_error(c, 'fort.21 is empty while parameter eta is set \
    9691031                                to 1 in CONTROL file')
    9701032
    9711033                # create the corresponding output file fort.15
    9721034                # (generated by CONVERT2) + fort.16 (paramId 167 and 168)
    973                 p = subprocess.check_call([os.path.expandvars(
    974                         os.path.expanduser(c.exedir)) + '/CONVERT2'], shell=True)
     1035                p = subprocess.check_call(
     1036                    [os.path.expandvars(os.path.expanduser(c.exedir)) +
     1037                     '/CONVERT2'], shell=True)
    9751038                os.chdir(pwd)
    9761039
     
    9831046                    suffix = cdateH[2:10]
    9841047                fnout += suffix
    985                 print("outputfile = " + fnout)
     1048                print "outputfile = " + fnout
    9861049                self.outputfilelist.append(fnout) # needed for final processing
    9871050
     
    10051068                            open(c.inputdir + '/fort.25', 'rb'), fout)
    10061069
    1007         try:
    1008             if c.wrf == '1':
    1009                 fwrf.close()
    1010         except:
    1011             pass
     1070        if c.wrf == '1':
     1071            fwrf.close()
    10121072
    10131073        grib_index_release(iid)
     
    10251085
    10261086        @Input:
    1027             self: instance of ECFlexpart
     1087            self: instance of EcFlexpart
    10281088                The current object of the class.
    10291089
    1030             inputfiles: instance of UIOFiles
     1090            inputfiles: instance of UioFiles
    10311091                Contains a list of files.
    10321092
     
    10501110        table128 = init128(c.ecmwfdatadir +
    10511111                           '/grib_templates/ecmwf_grib1_table_128')
    1052         pars = toparamId(self.params['OG_acc_SL'][0], table128)
     1112        pars = to_param_id(self.params['OG_acc_SL'][0], table128)
    10531113        index_keys = ["date", "time", "step"]
    10541114        indexfile = c.inputdir + "/date_time_stepRange.idx"
    1055         silentremove(indexfile)
    1056         grib = GribTools(inputfiles.files)
     1115        silent_remove(indexfile)
     1116        grib = Gribtools(inputfiles.files)
    10571117        # creates new index file
    10581118        iid = grib.index(index_keys=index_keys, index_file=indexfile)
     
    10611121        index_vals = []
    10621122        for key in index_keys:
    1063             key_vals = grib_index_get(iid,key)
    1064             print(key_vals)
     1123            key_vals = grib_index_get(iid, key)
     1124            print key_vals
    10651125            # have to sort the steps for disaggregation,
    10661126            # therefore convert to int first
     
    10831143            stepsdict[str(p)] = []
    10841144
     1145        print 'maxstep: ', c.maxstep
     1146
    10851147        for prod in product(*index_vals):
    10861148            # e.g. prod = ('20170505', '0', '12')
     
    10921154
    10931155            gid = grib_new_from_index(iid)
    1094             # do convert2 program if gid at this time is not None,
    1095             # therefore save in hid
    1096             hid = gid
    10971156            if gid is not None:
    10981157                cdate = grib_get(gid, 'date')
    10991158                time = grib_get(gid, 'time')
    1100                 type = grib_get(gid, 'type')
    11011159                step = grib_get(gid, 'step')
    11021160                # date+time+step-2*dtime
    11031161                # (since interpolated value valid for step-2*dtime)
    1104                 sdate = datetime.datetime(year=cdate/10000,
    1105                                           month=mod(cdate, 10000)/100,
    1106                                           day=mod(cdate, 100),
    1107                                           hour=time/100)
    1108                 fdate = sdate + datetime.timedelta(
    1109                                     hours=step-2*int(c.dtime))
    1110                 sdates = sdate + datetime.timedelta(hours=step)
     1162                sdate = datetime(year=cdate/10000,
     1163                                 month=(cdate % 10000)/100,
     1164                                 day=(cdate % 100),
     1165                                 hour=time/100)
     1166                fdate = sdate + timedelta(hours=step-2*int(c.dtime))
     1167                sdates = sdate + timedelta(hours=step)
     1168                elimit = None
    11111169            else:
    11121170                break
     
    11261184            else:
    11271185                fnout = c.inputdir + '/flux' + fdate.strftime('%Y%m%d%H')
    1128                 gnout = c.inputdir + '/flux' + (fdate+datetime.timedelta(
    1129                     hours = int(c.dtime))).strftime('%Y%m%d%H')
     1186                gnout = c.inputdir + '/flux' + (fdate +
     1187                                                timedelta(hours=int(c.dtime))
     1188                                               ).strftime('%Y%m%d%H')
    11301189                hnout = c.inputdir + '/flux' + sdates.strftime('%Y%m%d%H')
    11311190                g = open(gnout, 'w')
    11321191                h = open(hnout, 'w')
    11331192
    1134             print("outputfile = " + fnout)
     1193            print "outputfile = " + fnout
    11351194            f = open(fnout, 'w')
    11361195
     
    11561215                        fak = 3600.
    11571216
    1158                     values = (reshape(values, (nj, ni))).flatten() / fak
     1217                    values = (np.reshape(values, (nj, ni))).flatten() / fak
    11591218                    vdp.append(values[:])  # save the accumulated values
    11601219                    if step <= int(c.dtime):
     
    11641223
    11651224                    print(cparamId, atime, step, len(values),
    1166                           values[0], std(values))
     1225                          values[0], np.std(values))
    11671226                    # save the 1/3-hourly or specific values
    11681227                    # svdp.append(values[:])
     
    11721231                        if len(svdp) > 3:
    11731232                            if cparamId == '142' or cparamId == '143':
    1174                                 values = Disagg.darain(svdp)
     1233                                values = disaggregation.darain(svdp)
    11751234                            else:
    1176                                 values = Disagg.dapoly(svdp)
     1235                                values = disaggregation.dapoly(svdp)
    11771236
    11781237                            if not (step == c.maxstep and c.maxstep > 12 \
     
    11971256
    11981257                        if c.basetime is not None:
    1199                             elimit = datetime.datetime.strptime(c.end_date +
    1200                                                                 c.basetime,
    1201                                                                 '%Y%m%d%H')
     1258                            elimit = datetime.strptime(c.end_date +
     1259                                                       c.basetime, '%Y%m%d%H')
    12021260                        else:
    1203                             elimit = sdate + datetime.timedelta(2*int(c.dtime))
     1261                            elimit = sdate + timedelta(2*int(c.dtime))
    12041262
    12051263                        # squeeze out information of last two steps contained
    12061264                        # in svdp
    12071265                        # if step+int(c.dtime) == c.maxstep and c.maxstep>12
    1208                         # or sdates+datetime.timedelta(hours = int(c.dtime))
     1266                        # or sdates+timedelta(hours = int(c.dtime))
    12091267                        # >= elimit:
    12101268                        # Note that svdp[0] has not been popped in this case
    12111269
    1212                         if (step == c.maxstep and c.maxstep > 12
    1213                             or sdates == elimit):
     1270                        if step == c.maxstep and c.maxstep > 12 or \
     1271                           sdates == elimit:
    12141272
    12151273                            values = svdp[3]
    12161274                            grib_set_values(gid, values)
    12171275                            grib_set(gid, 'step', 0)
    1218                             truedatetime = fdate + datetime.timedelta(
    1219                                      hours=2*int(c.dtime))
     1276                            truedatetime = fdate + timedelta(hours=
     1277                                                             2*int(c.dtime))
    12201278                            grib_set(gid, 'time', truedatetime.hour * 100)
    12211279                            grib_set(gid, 'date', truedatetime.year * 10000 +
     
    12261284                            #values = (svdp[1]+svdp[2])/2.
    12271285                            if cparamId == '142' or cparamId == '143':
    1228                                 values = Disagg.darain(list(reversed(svdp)))
     1286                                values = disaggregation.darain(list(reversed(svdp)))
    12291287                            else:
    1230                                 values = Disagg.dapoly(list(reversed(svdp)))
    1231 
    1232                             grib_set(gid, 'step',0)
    1233                             truedatetime = fdate + datetime.timedelta(
    1234                                      hours=int(c.dtime))
     1288                                values = disaggregation.dapoly(list(reversed(svdp)))
     1289
     1290                            grib_set(gid, 'step', 0)
     1291                            truedatetime = fdate + timedelta(hours=int(c.dtime))
    12351292                            grib_set(gid, 'time', truedatetime.hour * 100)
    12361293                            grib_set(gid, 'date', truedatetime.year * 10000 +
     
    12501307        grib_index_release(iid)
    12511308
     1309        exit()
    12521310        return
  • python/GribTools.py

    rccab809 rff99eae  
    11#!/usr/bin/env python
    22# -*- coding: utf-8 -*-
    3 #************************************************************************
    4 # TODO AP
    5 # - GribTools name möglicherweise etwas verwirrend.
    6 # - change self.filename in self.filenames!!!
    7 # - bis auf --init-- und index wird keine Funktion verwendet!?
    8 #************************************************************************
    93#*******************************************************************************
    104# @Author: Anne Fouilloux (University of Oslo)
     
    3630# @Class Content:
    3731#    - __init__
    38 #    - getkeys
    39 #    - setkeys
     32#    - get_keys
     33#    - set_keys
    4034#    - copy
    4135#    - index
     36#
     37# @Class Attributes:
     38#    - filenames
    4239#
    4340#*******************************************************************************
     
    4744# ------------------------------------------------------------------------------
    4845import os
    49 from gribapi import *
     46from gribapi import grib_new_from_file, grib_is_defined, grib_get, \
     47                    grib_release, grib_set, grib_write, grib_index_read, \
     48                    grib_index_new_from_file, grib_index_add_file,  \
     49                    grib_index_write
    5050
    5151# ------------------------------------------------------------------------------
    5252# CLASS
    5353# ------------------------------------------------------------------------------
    54 class GribTools:
     54class Gribtools(object):
    5555    '''
    5656    Class for GRIB utilities (new methods) based on GRIB API
     
    6262        '''
    6363        @Description:
    64             Initialise an object of GribTools and assign a list
     64            Initialise an object of Gribtools and assign a list
    6565            of filenames.
    6666
     
    7373        '''
    7474
    75         self.filename = filenames
     75        self.filenames = filenames
    7676
    7777        return
    7878
    7979
    80     def getkeys(self, keynames, wherekeynames=[], wherekeyvalues=[]):
     80    def get_keys(self, keynames, wherekeynames=[], wherekeyvalues=[]):
    8181        '''
    8282        @Description:
     
    8888                List of keynames.
    8989
    90             wherekeynames: list of ???, optional
     90            wherekeynames: list of strings, optional
    9191                Default value is an empty list.
    9292
    93             wherekeyvalues: list of ???, optional
     93            wherekeyvalues: list of strings, optional
    9494                Default value is an empty list.
    9595
     
    9999        '''
    100100
    101         fileid = open(self.filename, 'r')
     101        fileid = open(self.filenames, 'r')
    102102
    103103        return_list = []
     
    136136
    137137
    138     def setkeys(self, fromfile, keynames, keyvalues, wherekeynames=[],
    139                 wherekeyvalues=[], strict=False, filemode='w'):
     138    def set_keys(self, fromfile, keynames, keyvalues, wherekeynames=[],
     139                 wherekeyvalues=[], strict=False, filemode='w'):
    140140        '''
    141141        @Description:
     
    150150                Filename of the input file to read the grib messages from.
    151151
    152             keynames: list of ???
     152            keynames: list of strings
    153153                List of keynames. Default is an empty list.
    154154
    155             keyvalues: list of ???
     155            keyvalues: list of strings
    156156                List of keynames. Default is an empty list.
    157157
    158             wherekeynames: list of ???, optional
     158            wherekeynames: list of strings, optional
    159159                Default value is an empty list.
    160160
    161             wherekeyvalues: list of ???, optional
     161            wherekeyvalues: list of strings, optional
    162162                Default value is an empty list.
    163163
     
    174174
    175175        '''
    176         fout = open(self.filename, filemode)
     176        fout = open(self.filenames, filemode)
    177177        fin = open(fromfile)
    178178
     
    196196                i += 1
    197197
    198 #AP is it secured that the order of keynames is equal to keyvalues?
    199198            if select:
    200199                i = 0
     
    203202                    i += 1
    204203
    205 #AP this is a redundant code section
    206 # delete the if/else :
    207 #
    208 #           grib_write(gid_in, fout)
    209 #
    210             if strict:
    211                 if select:
    212                     grib_write(gid_in, fout)
    213             else:
    214                 grib_write(gid_in, fout)
    215 #AP end
     204            grib_write(gid_in, fout)
     205
    216206            grib_release(gid_in)
    217207
     
    238228                function. Default is True.
    239229
    240             keynames: list of ???, optional
     230            keynames: list of strings, optional
    241231                List of keynames. Default is an empty list.
    242232
    243             keyvalues: list of ???, optional
     233            keyvalues: list of strings, optional
    244234                List of keynames. Default is an empty list.
    245235
     
    252242
    253243        fin = open(filename_in)
    254         fout = open(self.filename, filemode)
     244        fout = open(self.filenames, filemode)
    255245
    256246        while 1:
     
    307297                Grib index id.
    308298        '''
    309         print("... index will be done")
    310         self.iid = None
    311 
    312         if (os.path.exists(index_file)):
    313             self.iid = grib_index_read(index_file)
    314             print("Use existing index file: %s " % (index_file))
     299        print "... index will be done"
     300        iid = None
     301
     302        if os.path.exists(index_file):
     303            iid = grib_index_read(index_file)
     304            print "Use existing index file: %s " % (index_file)
    315305        else:
    316             for file in self.filename:
    317                 print("Inputfile: %s " % (file))
    318                 if self.iid is None:
    319                     self.iid = grib_index_new_from_file(file, index_keys)
     306            for filename in self.filenames:
     307                print "Inputfile: %s " % (filename)
     308                if iid is None:
     309                    iid = grib_index_new_from_file(filename, index_keys)
    320310                else:
    321311                    print 'in else zweig'
    322                     grib_index_add_file(self.iid, file)
    323 
    324             if self.iid is not None:
    325                 grib_index_write(self.iid, index_file)
    326 
    327         print('... index done')
    328 
    329         return self.iid
    330 
    331 
    332 
    333 
    334 
     312                    grib_index_add_file(iid, filename)
     313
     314            if iid is not None:
     315                grib_index_write(iid, index_file)
     316
     317        print '... index done'
     318
     319        return iid
  • python/MarsRetrieval.py

    r991df6a rff99eae  
    11#!/usr/bin/env python
    22# -*- coding: utf-8 -*-
    3 #************************************************************************
    4 # TODO AP
    5 # -
    6 #************************************************************************
    73#*******************************************************************************
    84# @Author: Anne Fouilloux (University of Oslo)
     
    139#
    1410#   November 2015 - Leopold Haimberger (University of Vienna):
    15 #        - optimized displayInfo
    16 #        - optimized dataRetrieve and seperate between python and shell
     11#        - optimized display_info
     12#        - optimized data_retrieve and seperate between python and shell
    1713#          script call
    1814#
     
    3733# @Class Content:
    3834#    - __init__
    39 #    - displayInfo
    40 #    - dataRetrieve
     35#    - display_info
     36#    - data_retrieve
     37#
     38# @Class Attributes:
     39#    - server
     40#    - marsclass
     41#    - dtype
     42#    - levtype
     43#    - levelist
     44#    - repres
     45#    - date
     46#    - resol
     47#    - stream
     48#    - area
     49#    - time
     50#    - step
     51#    - expver
     52#    - number
     53#    - accuracy
     54#    - grid
     55#    - gaussian
     56#    - target
     57#    - param
    4158#
    4259#*******************************************************************************
     
    4865import os
    4966
    50 ecapi = True
    51 try:
    52     import ecmwfapi
    53 except ImportError:
    54     ecapi = False
    55 
    5667# ------------------------------------------------------------------------------
    5768# CLASS
    5869# ------------------------------------------------------------------------------
    59 class MARSretrieval:
     70class MarsRetrieval(object):
    6071    '''
    6172    Class for submitting MARS retrievals.
     
    6879    '''
    6980
    70     def __init__(self, server, marsclass = "ei", type = "", levtype = "",
    71                  levelist = "", repres = "", date = "", resol = "", stream = "",
    72                  area = "", time = "", step = "", expver = "1", number = "",
    73                  accuracy = "", grid = "", gaussian = "", target = "",
    74                  param = ""):
     81    def __init__(self, server, marsclass="ei", dtype="", levtype="",
     82                 levelist="", repres="", date="", resol="", stream="",
     83                 area="", time="", step="", expver="1", number="",
     84                 accuracy="", grid="", gaussian="", target="",
     85                 param=""):
    7586        '''
    7687        @Description:
    77             Initialises the instance of the MARSretrieval class and
     88            Initialises the instance of the MarsRetrieval class and
    7889            defines and assigns a set of the necessary retrieval parameters
    7990            for the FLEXPART input data.
     
    8495
    8596        @Input:
    86             self: instance of MARSretrieval
     97            self: instance of MarsRetrieval
    8798                For description see class documentation.
    8899
     
    96107                Default is the ERA-Interim dataset "ei".
    97108
    98             type: string, optional
     109            dtype: string, optional
    99110                Determines the type of fields to be retrieved.
    100111                Selects between observations, images or fields.
     
    275286        self.server = server
    276287        self.marsclass = marsclass
    277         self.type = type
     288        self.dtype = dtype
    278289        self.levtype = levtype
    279290        self.levelist = levelist
     
    296307
    297308
    298     def displayInfo(self):
     309    def display_info(self):
    299310        '''
    300311        @Description:
     
    302313
    303314        @Input:
    304             self: instance of MARSretrieval
     315            self: instance of MarsRetrieval
    305316                For description see class documentation.
    306317
     
    314325        # with their corresponding values
    315326        for item in attrs.items():
    316             if item[0] in ('server'):
     327            if item[0] in 'server':
    317328                pass
    318329            else:
    319                 print(item[0] + ': ' + str(item[1]))
     330                print item[0] + ': ' + str(item[1])
    320331
    321332        return
    322333
    323     def dataRetrieve(self):
     334    def data_retrieve(self):
    324335        '''
    325336        @Description:
     
    330341
    331342        @Input:
    332             self: instance of MARSretrieval
     343            self: instance of MarsRetrieval
    333344                For description see class documentation.
    334345
     
    344355        s = 'ret'
    345356        for k, v in attrs.iteritems():
    346             if k in ('server'):
     357            if k in 'server':
    347358                continue
    348359            if k == 'marsclass':
     
    360371                self.server.execute(s, target)
    361372            except:
    362                 print('MARS Request failed, \
    363                     have you already registered at apps.ecmwf.int?')
     373                print 'MARS Request failed, \
     374                    have you already registered at apps.ecmwf.int?'
    364375                raise IOError
    365376            if os.stat(target).st_size == 0:
    366                 print('MARS Request returned no data - please check request')
     377                print 'MARS Request returned no data - please check request'
    367378                raise IOError
    368379        # MARS request via extra process in shell
     
    373384                                 stderr=subprocess.PIPE, bufsize=1)
    374385            pout = p.communicate(input=s)[0]
    375             print(pout.decode())
     386            print pout.decode()
    376387
    377388            if 'Some errors reported' in pout.decode():
    378                 print('MARS Request failed - please check request')
     389                print 'MARS Request failed - please check request'
    379390                raise IOError
    380391
    381392            if os.stat(target).st_size == 0:
    382                 print('MARS Request returned no data - please check request')
     393                print 'MARS Request returned no data - please check request'
    383394                raise IOError
    384395
  • python/UioFiles.py

    r991df6a rff99eae  
    11#!/usr/bin/env python
    22# -*- coding: utf-8 -*-
    3 #************************************************************************
    4 # TODO AP
    5 # - checken welche regelmässigen methoden auf diese Files noch angewendet werden
    6 # und dann hier implementieren
    7 # cleanup hier rein
    8 #************************************************************************
    93#*******************************************************************************
    104# @Author: Anne Fouilloux (University of Oslo)
     
    159#
    1610#    November 2015 - Leopold Haimberger (University of Vienna):
    17 #        - modified method listFiles to work with glob instead of listdir
    18 #        - added pattern search in method listFiles
     11#        - modified method list_files to work with glob instead of listdir
     12#        - added pattern search in method list_files
    1913#
    2014#    February 2018 - Anne Philipp (University of Vienna):
    2115#        - applied PEP8 style guide
    2216#        - added documentation
    23 #        - optimisation of method listFiles since it didn't work correctly
     17#        - optimisation of method list_files since it didn't work correctly
    2418#          for sub directories
    25 #        - additional speed up of method listFiles
     19#        - additional speed up of method list_files
    2620#        - modified the class so that it is initiated with a pattern instead
    2721#          of suffixes. Gives more precision in selection of files.
     
    4034# @Class Content:
    4135#    - __init__
    42 #    - listFiles
    43 #    - deleteFiles
     36#    - list_files
     37#    - delete_files
     38#
     39# @Class Attributes:
     40#    - pattern
     41#    - files
    4442#
    4543#*******************************************************************************
     
    4947# ------------------------------------------------------------------------------
    5048import os
    51 import glob
    5249import fnmatch
    53 import time
    5450
    5551# software specific module from flex_extract
    56 import profiling
    57 from Tools import silentremove
     52#import profiling
     53from tools import silent_remove
    5854
    5955# ------------------------------------------------------------------------------
     
    6157# ------------------------------------------------------------------------------
    6258
    63 class UIOFiles:
     59class UioFiles(object):
    6460    '''
    6561    Class to manipulate files. At initialisation it has the attribute
     
    7672
    7773        @Input:
    78             self: instance of UIOFiles
     74            self: instance of UioFiles
    7975                Description see class documentation.
    8076
     
    8783
    8884        self.pattern = pattern
     85        self.files = None
    8986
    9087        return
    9188
    9289    #@profiling.timefn
    93     def listFiles(self, path, callid=0):
     90    def list_files(self, path, callid=0):
    9491        '''
    9592        @Description:
     
    9895
    9996        @Input:
    100             self: instance of UIOFiles
     97            self: instance of UioFiles
    10198                Description see class documentation.
    10299
     
    132129        if subdirs:
    133130            for subdir in subdirs:
    134                 self.listFiles(os.path.join(path, subdir), callid=1)
     131                self.list_files(os.path.join(path, subdir), callid=1)
    135132
    136133        return
    137134
    138     def deleteFiles(self):
     135    def delete_files(self):
    139136        '''
    140137        @Description:
     
    142139
    143140        @Input:
    144             self: instance of UIOFiles
     141            self: instance of UioFiles
    145142                Description see class documentation.
    146143
     
    149146        '''
    150147
    151         for f in self.files:
    152             silentremove(f)
     148        for old_file in self.files:
     149            silent_remove(old_file)
    153150
    154151        return
  • python/disaggregation.py

    r991df6a rff99eae  
    22# -*- coding: utf-8 -*-
    33#************************************************************************
    4 # TODO AP
     4# ToDo AP
    55# - check alist of size 4 ?
    66# - write a test, IMPORTANT
     
    2020#        - added structured documentation
    2121#        - outsourced the disaggregation functions dapoly and darain
    22 #          to a new module named Disagg
     22#          to a new module named disaggregation
    2323#
    2424# @License:
     
    2929#
    3030# @Module Description:
    31 #    Disaggregation of deaccumulated flux data from an ECMWF model FG field.
     31#    disaggregationregation of deaccumulated flux data from an ECMWF model FG field.
    3232#    Initially the flux data to be concerned are:
    3333#    - large-scale precipitation
     
    7070        using a cubic polynomial solution which conserves the integrals
    7171        of the fluxes within each timespan.
    72         Disaggregation is done for 4 accumluated timespans which generates
     72        disaggregationregation is done for 4 accumluated timespans which generates
    7373        a new, disaggregated value which is output at the central point
    7474        of the 4 accumulation timespans. This new point is used for linear
     
    111111        field using a modified linear solution which conserves the integrals
    112112        of the fluxes within each timespan.
    113         Disaggregation is done for 4 accumluated timespans which generates
     113        disaggregationregation is done for 4 accumluated timespans which generates
    114114        a new, disaggregated value which is output at the central point
    115115        of the 4 accumulation timespans. This new point is used for linear
     
    144144
    145145    return nfield
    146 
    147 
    148 
    149 
    150 
    151 
    152 
    153 
    154 
    155 
  • python/get_mars_data.py

    r991df6a rff99eae  
    11#!/usr/bin/env python
    22# -*- coding: utf-8 -*-
    3 #************************************************************************
    4 # TODO AP
    5 # - add function docstrings!!!!
    6 #************************************************************************
    73#*******************************************************************************
    84# @Author: Anne Fouilloux (University of Oslo)
     
    139#
    1410#    November 2015 - Leopold Haimberger (University of Vienna):
    15 #        - moved the getEIdata program into a function "getMARSdata"
     11#        - moved the getEIdata program into a function "get_mars_data"
    1612#        - moved the AgurmentParser into a seperate function
    1713#        - adatpted the function for the use in flex_extract
    18 #        - renamed file to getMARSdata
     14#        - renamed file to get_mars_data
    1915#
    2016#    February 2018 - Anne Philipp (University of Vienna):
     
    4238# @Program Content:
    4339#    - main
    44 #    - getMARSdata
     40#    - get_mars_data
    4541#
    4642#*******************************************************************************
     
    5450import inspect
    5551try:
    56     ecapi=True
     52    ecapi = True
    5753    import ecmwfapi
    5854except ImportError:
    59     ecapi=False
     55    ecapi = False
     56
     57# software specific classes and modules from flex_extract
     58from tools import my_error, normal_exit, interpret_args_and_control
     59from EcFlexpart import EcFlexpart
     60from UioFiles import UioFiles
    6061
    6162# add path to pythonpath so that python finds its buddies
    62 localpythonpath = os.path.dirname(os.path.abspath(
     63LOCAL_PYTHON_PATH = os.path.dirname(os.path.abspath(
    6364    inspect.getfile(inspect.currentframe())))
    64 if localpythonpath not in sys.path:
    65     sys.path.append(localpythonpath)
    66 
    67 # software specific classes and modules from flex_extract
    68 from ControlFile import ControlFile
    69 from Tools import myerror, normalexit, \
    70                           interpret_args_and_control
    71 from ECFlexpart import ECFlexpart
    72 from UIOFiles import UIOFiles
     65if LOCAL_PYTHON_PATH not in sys.path:
     66    sys.path.append(LOCAL_PYTHON_PATH)
    7367
    7468# ------------------------------------------------------------------------------
     
    7872    '''
    7973    @Description:
    80         If getMARSdata is called from command line, this function controls
     74        If get_mars_data is called from command line, this function controls
    8175        the program flow and calls the argumentparser function and
    82         the getMARSdata function for retrieving EC data.
     76        the get_mars_data function for retrieving EC data.
    8377
    8478    @Input:
     
    8983    '''
    9084    args, c = interpret_args_and_control()
    91     getMARSdata(args, c)
    92     normalexit(c)
     85    get_mars_data(c)
     86    normal_exit(c)
    9387
    9488    return
    9589
    96 def getMARSdata(args, c):
     90def get_mars_data(c):
    9791    '''
    9892    @Description:
     
    10397
    10498    @Input:
    105         args: instance of ArgumentParser
    106             Contains the commandline arguments from script/program call.
    107 
    10899        c: instance of class ControlFile
    109100            Contains all the parameters of CONTROL file, which are e.g.:
     
    126117        os.makedirs(c.inputdir)
    127118
    128     print("Retrieving EC data!")
    129     print("start date %s " % (c.start_date))
    130     print("end date %s " % (c.end_date))
     119    print "Retrieving EC data!"
     120    print "start date %s " % (c.start_date)
     121    print "end date %s " % (c.end_date)
    131122
    132123    if ecapi:
     
    160151    # --------------  flux data ------------------------------------------------
    161152    print 'removing old flux content of ' + c.inputdir
    162     tobecleaned = UIOFiles('*_acc_*.' + str(os.getppid()) + '.*.grb')
    163     tobecleaned.listFiles(c.inputdir)
    164     tobecleaned.deleteFiles()
     153    tobecleaned = UioFiles('*_acc_*.' + str(os.getppid()) + '.*.grb')
     154    tobecleaned.list_files(c.inputdir)
     155    tobecleaned.delete_files()
    165156
    166157    # if forecast for maximum one day (upto 23h) are to be retrieved,
     
    172163        while day < endp1:
    173164            # retrieve MARS data for the whole period
    174             flexpart = ECFlexpart(c, fluxes=True)
     165            flexpart = EcFlexpart(c, fluxes=True)
    175166            tmpday = day + datechunk - datetime.timedelta(days=1)
    176167            if tmpday < endp1:
     
    186177                flexpart.retrieve(server, dates, c.inputdir)
    187178            except IOError:
    188                 myerror(c, 'MARS request failed')
     179                my_error(c, 'MARS request failed')
    189180
    190181            day += datechunk
     
    199190        while day <= end:
    200191            # retrieve MARS data for the whole period
    201             flexpart = ECFlexpart(c, fluxes=True)
     192            flexpart = EcFlexpart(c, fluxes=True)
    202193            tmpday = day + datechunk - datetime.timedelta(days=1)
    203194            if tmpday < end:
     
    213204                flexpart.retrieve(server, dates, c.inputdir)
    214205            except IOError:
    215                 myerror(c, 'MARS request failed')
     206                my_error(c, 'MARS request failed')
    216207
    217208            day += datechunk
     
    219210    # --------------  non flux data --------------------------------------------
    220211    print 'removing old non flux content of ' + c.inputdir
    221     tobecleaned = UIOFiles('*__*.' + str(os.getppid()) + '.*.grb')
    222     tobecleaned.listFiles(c.inputdir)
    223     tobecleaned.deleteFiles()
     212    tobecleaned = UioFiles('*__*.' + str(os.getppid()) + '.*.grb')
     213    tobecleaned.list_files(c.inputdir)
     214    tobecleaned.delete_files()
    224215
    225216    day = start
    226217    while day <= end:
    227             # retrieve all non flux MARS data for the whole period
    228             flexpart = ECFlexpart(c, fluxes=False)
    229             tmpday = day + datechunk - datetime.timedelta(days=1)
    230             if tmpday < end:
    231                 dates = day.strftime("%Y%m%d") + "/to/" + \
    232                         tmpday.strftime("%Y%m%d")
    233             else:
    234                 dates = day.strftime("%Y%m%d") + "/to/" + \
    235                         end.strftime("%Y%m%d")
    236 
    237             print "retrieve " + dates + " in dir " + c.inputdir
    238 
    239             try:
    240                 flexpart.retrieve(server, dates, c.inputdir)
    241             except IOError:
    242                 myerror(c, 'MARS request failed')
    243 
    244             day += datechunk
     218        # retrieve all non flux MARS data for the whole period
     219        flexpart = EcFlexpart(c, fluxes=False)
     220        tmpday = day + datechunk - datetime.timedelta(days=1)
     221        if tmpday < end:
     222            dates = day.strftime("%Y%m%d") + "/to/" + \
     223                    tmpday.strftime("%Y%m%d")
     224        else:
     225            dates = day.strftime("%Y%m%d") + "/to/" + \
     226                    end.strftime("%Y%m%d")
     227
     228        print "retrieve " + dates + " in dir " + c.inputdir
     229
     230        try:
     231            flexpart.retrieve(server, dates, c.inputdir)
     232        except IOError:
     233            my_error(c, 'MARS request failed')
     234
     235        day += datechunk
    245236
    246237    return
     
    248239if __name__ == "__main__":
    249240    main()
    250 
  • python/install.py

    r991df6a rff99eae  
    22# -*- coding: utf-8 -*-
    33#************************************************************************
    4 # TODO AP
    5 # - localpythonpath should not be set in module load section!
     4# ToDo AP
    65# - create a class Installation and divide installation in 3 subdefs for
    76#   ecgate, local and cca seperatly
     
    1817#        - applied PEP8 style guide
    1918#        - added documentation
     19#        - moved install_args_and_control in here
    2020#
    2121# @License:
     
    4545# MODULES
    4646# ------------------------------------------------------------------------------
    47 import datetime
    4847import os
    4948import sys
     
    5150import subprocess
    5251import inspect
    53 from argparse import ArgumentParser,ArgumentDefaultsHelpFormatter
    54 
    55 # add path to pythonpath so that python finds its buddies
    56 localpythonpath = os.path.dirname(os.path.abspath(
    57     inspect.getfile(inspect.currentframe())))
    58 if localpythonpath not in sys.path:
    59     sys.path.append(localpythonpath)
     52from argparse import ArgumentParser, ArgumentDefaultsHelpFormatter
    6053
    6154# software specific classes and modules from flex_extract
    6255from ControlFile import ControlFile
     56
     57# add path to pythonpath so that python finds its buddies
     58LOCAL_PYTHON_PATH = os.path.dirname(os.path.abspath(
     59    inspect.getfile(inspect.currentframe())))
     60if LOCAL_PYTHON_PATH not in sys.path:
     61    sys.path.append(LOCAL_PYTHON_PATH)
    6362
    6463# ------------------------------------------------------------------------------
     
    7877    '''
    7978
    80     os.chdir(localpythonpath)
     79    os.chdir(LOCAL_PYTHON_PATH)
    8180    args, c = install_args_and_control()
    8281
     
    8483        install_via_gateway(c, args.install_target)
    8584    else:
    86         print('Please specify installation target (local|ecgate|cca)')
    87         print('use -h or --help for help')
     85        print 'Please specify installation target (local|ecgate|cca)'
     86        print 'use -h or --help for help'
    8887
    8988    sys.exit()
     
    152151    try:
    153152        c = ControlFile(args.controlfile)
    154     except:
    155         print('Could not read CONTROL file "' + args.controlfile + '"')
    156         print('Either it does not exist or its syntax is wrong.')
    157         print('Try "' + sys.argv[0].split('/')[-1] +
    158               ' -h" to print usage information')
     153    except IOError:
     154        print 'Could not read CONTROL file "' + args.controlfile + '"'
     155        print 'Either it does not exist or its syntax is wrong.'
     156        print 'Try "' + sys.argv[0].split('/')[-1] + \
     157              ' -h" to print usage information'
    159158        exit(1)
    160159
    161160    if args.install_target != 'local':
    162         if (args.ecgid is None or args.ecuid is None or args.gateway is None
    163             or args.destination is None):
    164             print('Please enter your ECMWF user id and group id as well as \
     161        if args.ecgid is None or args.ecuid is None or args.gateway is None \
     162           or args.destination is None:
     163            print 'Please enter your ECMWF user id and group id as well as \
    165164                   the \nname of the local gateway and the ectrans \
    166                    destination ')
    167             print('with command line options --ecuid --ecgid \
    168                    --gateway --destination')
    169             print('Try "' + sys.argv[0].split('/')[-1] +
    170                   ' -h" to print usage information')
    171             print('Please consult ecaccess documentation or ECMWF user support \
    172                    for further details')
     165                   destination '
     166            print 'with command line options --ecuid --ecgid \
     167                   --gateway --destination'
     168            print 'Try "' + sys.argv[0].split('/')[-1] + \
     169                  ' -h" to print usage information'
     170            print 'Please consult ecaccess documentation or ECMWF user support \
     171                   for further details'
    173172            sys.exit(1)
    174173        else:
     
    178177            c.destination = args.destination
    179178
    180     try:
     179    if args.makefile:
    181180        c.makefile = args.makefile
    182     except:
    183         pass
    184181
    185182    if args.install_target == 'local':
     
    240237                            c.flexpart_root_scripts
    241238                else:
    242                     data='export FLEXPART_ROOT_SCRIPTS=$HOME'
     239                    data = 'export FLEXPART_ROOT_SCRIPTS=$HOME'
    243240            if target.lower() != 'local':
    244241                if '--workdir' in data:
     
    292289        fo.close()
    293290
    294 
    295 
    296291    if target.lower() == 'local':
    297292        # compile CONVERT2
    298293        if c.flexpart_root_scripts is None or c.flexpart_root_scripts == '../':
    299             print('Warning: FLEXPART_ROOT_SCRIPTS has not been specified')
    300             print('Only CONVERT2 will be compiled in ' + ecd + '/../src')
     294            print 'Warning: FLEXPART_ROOT_SCRIPTS has not been specified'
     295            print 'Only CONVERT2 will be compiled in ' + ecd + '/../src'
    301296        else:
    302297            c.flexpart_root_scripts = os.path.expandvars(os.path.expanduser(
    303                                         c.flexpart_root_scripts))
     298                c.flexpart_root_scripts))
    304299            if os.path.abspath(ecd) != os.path.abspath(c.flexpart_root_scripts):
    305300                os.chdir('/')
     
    311306                try:
    312307                    os.makedirs(c.flexpart_root_scripts + '/ECMWFDATA7.1')
    313                 except:
     308                finally:
    314309                    pass
    315310                os.chdir(c.flexpart_root_scripts + '/ECMWFDATA7.1')
     
    329324            p = subprocess.check_call(['rm'] + flist)
    330325        try:
    331             print(('Using makefile: ' + makefile))
     326            print 'Using makefile: ' + makefile
    332327            p = subprocess.check_call(['make', '-f', makefile])
    333             p = subprocess.check_call(['ls', '-l','CONVERT2'])
    334         except:
    335             print('compile failed - please edit ' + makefile +
    336                   ' or try another Makefile in the src directory.')
    337             print('most likely GRIB_API_INCLUDE_DIR, GRIB_API_LIB '
    338                     'and EMOSLIB must be adapted.')
    339             print('Available Makefiles:')
    340             print(glob.glob('Makefile*'))
    341 
     328            p = subprocess.check_call(['ls', '-l', 'CONVERT2'])
     329        except subprocess.CalledProcessError as e:
     330            print 'compile failed with the following error:'
     331            print e.output
     332            print 'please edit ' + makefile + \
     333                  ' or try another Makefile in the src directory.'
     334            print 'most likely GRIB_API_INCLUDE_DIR, GRIB_API_LIB  \
     335                   and EMOSLIB must be adapted.'
     336            print 'Available Makefiles:'
     337            print glob.glob('Makefile*')
    342338    elif target.lower() == 'ecgate':
    343339        os.chdir('/')
     
    352348                                       'ecgate:/home/ms/' + c.ecgid + '/' +
    353349                                       c.ecuid + '/ECMWFDATA7.1.tar'])
    354         except:
    355             print('ecaccess-file-put failed! Probably the eccert key has expired.')
     350        except subprocess.CalledProcessError as e:
     351            print 'ecaccess-file-put failed! \
     352                   Probably the eccert key has expired.'
    356353            exit(1)
    357         p = subprocess.check_call(['ecaccess-job-submit',
    358                                    '-queueName',
    359                                    target,
    360                                    ecd + 'python/compilejob.ksh'])
    361         print('compilejob.ksh has been submitted to ecgate for '
    362                 'installation in ' + c.ec_flexpart_root_scripts +
    363                 '/ECMWFDATA7.1')
    364         print('You should get an email with subject flexcompile within '
    365                 'the next few minutes')
     354
     355        try:
     356            p = subprocess.check_call(['ecaccess-job-submit',
     357                                       '-queueName',
     358                                       target,
     359                                       ecd + 'python/compilejob.ksh'])
     360            print 'compilejob.ksh has been submitted to ecgate for  \
     361                   installation in ' + c.ec_flexpart_root_scripts + \
     362                   '/ECMWFDATA7.1'
     363            print 'You should get an email with subject flexcompile within  \
     364                   the next few minutes'
     365        except subprocess.CalledProcessError as e:
     366            print 'ecaccess-job-submit failed!'
     367            exit(1)
    366368
    367369    elif target.lower() == 'cca':
     
    377379                                       'cca:/home/ms/' + c.ecgid + '/' +
    378380                                       c.ecuid + '/ECMWFDATA7.1.tar'])
    379         except:
    380             print('ecaccess-file-put failed! '
    381                     'Probably the eccert key has expired.')
     381        except subprocess.CalledProcessError as e:
     382            print 'ecaccess-file-put failed! \
     383                   Probably the eccert key has expired.'
    382384            exit(1)
    383385
    384         p=subprocess.check_call(['ecaccess-job-submit',
    385                                 '-queueName',
    386                                 target,
    387                                 ecd + 'python/compilejob.ksh'])
    388         print('compilejob.ksh has been submitted to cca for installation in ' +
    389               c.ec_flexpart_root_scripts + '/ECMWFDATA7.1')
    390         print('You should get an email with subject flexcompile '
    391                 'within the next few minutes')
     386        try:
     387            p = subprocess.check_call(['ecaccess-job-submit',
     388                                       '-queueName',
     389                                       target,
     390                                       ecd + 'python/compilejob.ksh'])
     391            print 'compilejob.ksh has been submitted to cca for installation in ' +\
     392                  c.ec_flexpart_root_scripts + '/ECMWFDATA7.1'
     393            print 'You should get an email with subject flexcompile \
     394                   within the next few minutes'
     395        except subprocess.CalledProcessError as e:
     396            print 'ecaccess-job-submit failed!'
     397            exit(1)
    392398
    393399    else:
    394         print('ERROR: unknown installation target ', target)
    395         print('Valid targets: ecgate, cca, local')
     400        print 'ERROR: unknown installation target ', target
     401        print 'Valid targets: ecgate, cca, local'
    396402
    397403    return
  • python/plot_retrieved.py

    rccab809 rff99eae  
    22# -*- coding: utf-8 -*-
    33#************************************************************************
    4 # TODO AP
     4# ToDo AP
    55# - documentation der Funktionen
    66# - docu der progam functionality
     
    1919#        - created function main and moved the two function calls for
    2020#          arguments and plotting into it
    21 #        - added function getBasics to extract the boundary conditions
     21#        - added function get_basics to extract the boundary conditions
    2222#          of the data fields from the first grib file it gets.
    2323#
     
    3333# @Program Content:
    3434#    - main
    35 #    - getBasics
    36 #    - plotRetrieved
    37 #    - plotTS
    38 #    - plotMap
    39 #    - getPlotArgs
     35#    - get_basics
     36#    - plot_retrieved
     37#    - plot_timeseries
     38#    - plot_map
     39#    - get_plot_args
    4040#
    4141#*******************************************************************************
     
    4949import inspect
    5050import sys
    51 import glob
    5251from argparse import ArgumentParser, ArgumentDefaultsHelpFormatter
    53 
    54 from matplotlib.pylab import *
    55 import matplotlib.patches as mpatches
    56 from mpl_toolkits.basemap import Basemap, addcyclic
    57 import matplotlib.colors as mcolors
    58 from matplotlib.font_manager import FontProperties
    59 from matplotlib.patches import Polygon
    60 import matplotlib.cm as cmx
    61 import matplotlib.colors as colors
    62 #from rasotools.utils import stats
     52import matplotlib
     53import matplotlib.pyplot as plt
     54from mpl_toolkits.basemap import Basemap
     55from eccodes import codes_grib_new_from_file, codes_get, codes_release, \
     56                    codes_get_values
     57import numpy as np
     58
     59# software specific classes and modules from flex_extract
     60from ControlFile import ControlFile
     61from UioFiles import UioFiles
     62
     63# add path to pythonpath so that python finds its buddies
     64LOCAL_PYTHON_PATH = os.path.dirname(os.path.abspath(
     65    inspect.getfile(inspect.currentframe())))
     66if LOCAL_PYTHON_PATH not in sys.path:
     67    sys.path.append(LOCAL_PYTHON_PATH)
    6368
    6469font = {'family': 'monospace', 'size': 12}
    6570matplotlib.rcParams['xtick.major.pad'] = '20'
    66 
    6771matplotlib.rc('font', **font)
    68 
    69 from eccodes import *
    70 import numpy as np
    71 
    72 # add path to pythonpath so that python finds its buddies
    73 localpythonpath = os.path.dirname(os.path.abspath(
    74     inspect.getfile(inspect.currentframe())))
    75 if localpythonpath not in sys.path:
    76     sys.path.append(localpythonpath)
    77 
    78 # software specific classes and modules from flex_extract
    79 from Tools import silentremove
    80 from ControlFile import ControlFile
    81 #from GribTools import GribTools
    82 from UIOFiles import UIOFiles
    83 
    8472# ------------------------------------------------------------------------------
    8573# FUNCTION
     
    8876    '''
    8977    @Description:
    90         If plotRetrieved is called from command line, this function controls
     78        If plot_retrieved is called from command line, this function controls
    9179        the program flow and calls the argumentparser function and
    92         the plotRetrieved function for plotting the retrieved GRIB data.
     80        the plot_retrieved function for plotting the retrieved GRIB data.
    9381
    9482    @Input:
     
    9886        <nothing>
    9987    '''
    100     args, c = getPlotArgs()
    101     plotRetrieved(args, c)
     88    args, c = get_plot_args()
     89    plot_retrieved(c)
    10290
    10391    return
    10492
    105 def getBasics(ifile, verb=False):
     93def get_basics(ifile, verb=False):
    10694    """
    10795    @Description:
     
    113101        ifile: string
    114102            Contains the full absolute path to the ECMWF grib file.
     103
    115104        verb (opt): bool
    116105            Is True if there should be extra output in verbose mode.
     
    127116            'iDirectionIncrementInDegrees'
    128117    """
    129     from eccodes import *
    130118
    131119    data = {}
     
    141129
    142130        # information needed from grib message
    143         keys = [
    144                 'Ni',
     131        keys = ['Ni',
    145132                'Nj',
    146133                'latitudeOfFirstGridPointInDegrees',
     
    149136                'longitudeOfLastGridPointInDegrees',
    150137                'jDirectionIncrementInDegrees',
    151                 'iDirectionIncrementInDegrees'
    152                ]
    153 
    154         if verb: print('\nInformations are: ')
     138                'iDirectionIncrementInDegrees']
     139
     140        if verb:
     141            print '\nInformations are: '
    155142        for key in keys:
    156143            # Get the value of the key in a grib message.
    157             data[key] = codes_get(gid,key)
    158             if verb: print "%s = %s" % (key,data[key])
    159         if verb: print '\n'
     144            data[key] = codes_get(gid, key)
     145            if verb:
     146                print "%s = %s" % (key, data[key])
     147        if verb:
     148            print '\n'
    160149
    161150        # Free the memory for the message referred as gribid.
     
    164153    return data
    165154
    166 def getFilesPerDate(files, datelist):
     155def get_files_per_date(files, datelist):
    167156    '''
    168157    @Description:
     
    171160
    172161    @Input:
    173         files: instance of UIOFiles
     162        files: instance of UioFiles
    174163            For description see class documentation.
    175164            It contains the attribute "files" which is a list of pathes
     
    185174
    186175    filelist = []
    187     for file in files:
    188         fdate = file[-8:]
    189         ddate = datetime.datetime.strptime(fdate, '%y%m%d%H')
     176    for filename in files:
     177        filedate = filename[-8:]
     178        ddate = datetime.datetime.strptime(filedate, '%y%m%d%H')
    190179        if ddate in datelist:
    191             filelist.append(file)
     180            filelist.append(filename)
    192181
    193182    return filelist
    194183
    195 def plotRetrieved(args, c):
     184def plot_retrieved(c):
    196185    '''
    197186    @Description:
    198187        Reads GRIB data from a specified time period, a list of levels
    199188        and a specified list of parameter.
    200 
    201     @Input:
    202         args: instance of ArgumentParser
    203             Contains the commandline arguments from script/program call.
    204 
    205         c: instance of class ControlFile
    206             Contains all necessary information of a CONTROL file. The parameters
    207             are: DAY1, DAY2, DTIME, MAXSTEP, TYPE, TIME, STEP, CLASS, STREAM,
    208             NUMBER, EXPVER, GRID, LEFT, LOWER, UPPER, RIGHT, LEVEL, LEVELIST,
    209             RESOL, GAUSS, ACCURACY, OMEGA, OMEGADIFF, ETA, ETADIFF, DPDETA,
    210             SMOOTH, FORMAT, ADDPAR, WRF, CWC, PREFIX, ECSTORAGE, ECTRANS,
    211             ECFSDIR, MAILOPS, MAILFAIL, GRIB2FLEXPART, DEBUG, INPUTDIR,
    212             OUTPUTDIR, FLEXPART_ROOT_SCRIPTS
    213             For more information about format and content of the parameter see
    214             documentation.
    215 
    216     @Return:
    217         <nothing>
    218     '''
    219     start = datetime.datetime.strptime(c.start_date, '%Y%m%d%H')
    220     end = datetime.datetime.strptime(c.end_date, '%Y%m%d%H')
    221 
    222     # create datelist between start and end date
    223     datelist = [start] # initialise datelist with first date
    224     run_date = start
    225     while run_date < end:
    226         run_date += datetime.timedelta(hours=int(c.dtime))
    227         datelist.append(run_date)
    228 
    229     print 'datelist: ', datelist
    230 
    231     c.paramIds = asarray(c.paramIds, dtype='int')
    232     c.levels = asarray(c.levels, dtype='int')
    233     c.area = asarray(c.area)
    234 
    235     # index_keys = ["date", "time", "step"]
    236     # indexfile = c.inputdir + "/date_time_stepRange.idx"
    237     # silentremove(indexfile)
    238     # grib = GribTools(inputfiles.files)
    239     # # creates new index file
    240     # iid = grib.index(index_keys=index_keys, index_file=indexfile)
    241 
    242     # # read values of index keys
    243     # index_vals = []
    244     # for key in index_keys:
    245         # index_vals.append(grib_index_get(iid, key))
    246         # print(index_vals[-1])
    247         # # index_vals looks for example like:
    248         # # index_vals[0]: ('20171106', '20171107', '20171108') ; date
    249         # # index_vals[1]: ('0', '1200', '1800', '600') ; time
    250         # # index_vals[2]: ('0', '12', '3', '6', '9') ; stepRange
    251 
    252 
    253 
    254 
    255     #index_keys = ["date", "time", "step", "paramId"]
    256     #indexfile = c.inputdir + "/date_time_stepRange.idx"
    257     #silentremove(indexfile)
    258 
    259     files = UIOFiles(c.prefix+'*')
    260     files.listFiles(c.inputdir)
    261     ifiles = getFilesPerDate(files.files, datelist)
    262     ifiles.sort()
    263 
    264     gdict = getBasics(ifiles[0], verb=False)
    265 
    266     fdict = dict()
    267     fmeta = dict()
    268     fstamp = dict()
    269     for p in c.paramIds:
    270         for l in c.levels:
    271             key = '{:0>3}_{:0>3}'.format(p, l)
    272             fdict[key] = []
    273             fmeta[key] = []
    274             fstamp[key] = []
    275 
    276     for file in ifiles:
    277         f = open(file)
    278         print( "Opening file for reading data --- %s" % file)
    279         fdate = datetime.datetime.strptime(file[-8:], "%y%m%d%H")
    280 
    281         # Load in memory a grib message from a file.
    282         gid = codes_grib_new_from_file(f)
    283         while(gid is not None):
    284             gtype = codes_get(gid, 'type')
    285             paramId = codes_get(gid, 'paramId')
    286             parameterName = codes_get(gid, 'parameterName')
    287             level = codes_get(gid, 'level')
    288 
    289             if paramId in c.paramIds and level in c.levels:
    290                 key = '{:0>3}_{:0>3}'.format(paramId, level)
    291                 print 'key: ', key
    292                 if len(fstamp[key]) != 0 :
    293                     for i in range(len(fstamp[key])):
    294                         if fdate < fstamp[key][i]:
    295                             fstamp[key].insert(i, fdate)
    296                             fmeta[key].insert(i, [paramId, parameterName, gtype,
    297                                                 fdate, level])
    298                             fdict[key].insert(i, flipud(reshape(
    299                                             codes_get_values(gid),
    300                                             [gdict['Nj'], gdict['Ni']])))
    301                             break
    302                         elif fdate > fstamp[key][i] and i == len(fstamp[key])-1:
    303                             fstamp[key].append(fdate)
    304                             fmeta[key].append([paramId, parameterName, gtype,
    305                                                 fdate, level])
    306                             fdict[key].append(flipud(reshape(
    307                                             codes_get_values(gid),
    308                                             [gdict['Nj'], gdict['Ni']])))
    309                             break
    310                         elif fdate > fstamp[key][i] and i != len(fstamp[key])-1 \
    311                              and fdate < fstamp[key][i+1]:
    312                             fstamp[key].insert(i, fdate)
    313                             fmeta[key].insert(i, [paramId, parameterName, gtype,
    314                                                     fdate, level])
    315                             fdict[key].insert(i, flipud(reshape(
    316                                             codes_get_values(gid),
    317                                             [gdict['Nj'], gdict['Ni']])))
    318                             break
    319                         else:
    320                             pass
    321                 else:
    322                     fstamp[key].append(fdate)
    323                     fmeta[key].append((paramId, parameterName, gtype,
    324                                        fdate, level))
    325                     fdict[key].append(flipud(reshape(
    326                         codes_get_values(gid), [gdict['Nj'], gdict['Ni']])))
    327 
    328             codes_release(gid)
    329 
    330             # Load in memory a grib message from a file.
    331             gid = codes_grib_new_from_file(f)
    332 
    333         f.close()
    334     #print 'fstamp: ', fstamp
    335     #exit()
    336 
    337     for k in fdict.keys():
    338         print 'fmeta: ', len(fmeta),fmeta
    339         fml = fmeta[k]
    340         fdl = fdict[k]
    341         print 'fm1: ', len(fml), fml
    342         #print 'fd1: ', fdl
    343         #print zip(fdl, fml)
    344         for fd, fm in zip(fdl, fml):
    345             print fm
    346             ftitle = fm[1] + ' {} '.format(fm[-1]) + \
    347                 datetime.datetime.strftime(fm[3], '%Y%m%d%H') #+ ' ' + stats(fd)
    348             pname = '_'.join(fm[1].split()) + '_{}_'.format(fm[-1]) + \
    349                 datetime.datetime.strftime(fm[3], '%Y%m%d%H')
    350             plotMap(c, fd, fm, gdict, ftitle, pname, 'png')
    351 
    352     for k in fdict.keys():
    353         fml = fmeta[k]
    354         fdl = fdict[k]
    355         fsl = fstamp[k]
    356         if fdl:
    357             fm = fml[0]
    358             fd = fdl[0]
    359             ftitle = fm[1] + ' {} '.format(fm[-1]) + \
    360                 datetime.datetime.strftime(fm[3], '%Y%m%d%H') #+ ' ' + stats(fd)
    361             pname = '_'.join(fm[1].split()) + '_{}_'.format(fm[-1]) + \
    362                 datetime.datetime.strftime(fm[3], '%Y%m%d%H')
    363             lat = -20
    364             lon = 20
    365             plotTS(c, fdl, fml, fsl, lat, lon,
    366                            gdict, ftitle, pname, 'png')
    367 
    368     return
    369 
    370 def plotTS(c, flist, fmetalist, ftimestamps, lat, lon,
    371                    gdict, ftitle, filename, fending, show=False):
    372     '''
    373     @Description:
    374 
    375     @Input:
    376         c:
    377 
    378         flist:
    379             The actual data values to be plotted from the grib messages.
    380 
    381         fmetalist: list of strings
    382             Contains some meta date for the data field to be plotted:
    383             parameter id, parameter Name, grid type, datetime, level
    384 
    385         ftimestamps: list of datetime
    386             Contains the time stamps.
    387 
    388         lat:
    389 
    390         lon:
    391 
    392         gdict:
    393 
    394         ftitle: string
    395             The title of the timeseries.
    396 
    397         filename: string
    398             The time series is stored in a file with this name.
    399 
    400         fending: string
    401             Contains the type of plot, e.g. pdf or png
    402 
    403         show: boolean
    404             Decides if the plot is shown after plotting or hidden.
    405 
    406     @Return:
    407         <nothing>
    408     '''
    409     print 'plotting timeseries'
    410 
    411     t1 = time.time()
    412 
    413     llx = gdict['longitudeOfFirstGridPointInDegrees']
    414     if llx > 180. :
    415         llx -= 360.
    416     lly = gdict['latitudeOfLastGridPointInDegrees']
    417     dxout = gdict['iDirectionIncrementInDegrees']
    418     dyout = gdict['jDirectionIncrementInDegrees']
    419     urx = gdict['longitudeOfLastGridPointInDegrees']
    420     ury = gdict['latitudeOfFirstGridPointInDegrees']
    421     numxgrid = gdict['Ni']
    422     numygrid = gdict['Nj']
    423 
    424     farr = asarray(flist)
    425     #(time, lat, lon)
    426 
    427     lonindex = linspace(llx, urx, numxgrid)
    428     latindex = linspace(lly, ury, numygrid)
    429     #print len(lonindex), len(latindex), farr.shape
    430 
    431     #latindex = (lat + 90) * 180 / (gdict['Nj'] - 1)
    432     #lonindex = (lon + 179) * 360 / gdict['Ni']
    433     #print latindex, lonindex
    434 
    435 
    436     ts = farr[:, 0, 0]#latindex[0], lonindex[0]]
    437 
    438     fig = plt.figure(figsize=(12, 6.7))
    439 
    440     plt.plot(ftimestamps, ts)
    441     plt.title(ftitle)
    442 
    443     plt.savefig(c.outputdir+'/'+filename+'_TS.'+fending, facecolor=fig.get_facecolor(), edgecolor='none',format=fending)
    444     print 'created ', c.outputdir + '/' + filename
    445     if show == True:
    446         plt.show()
    447     fig.clf()
    448     plt.close(fig)
    449 
    450     print time.time() - t1, 's'
    451 
    452     return
    453 
    454 def plotMap(c, flist, fmetalist, gdict, ftitle, filename, fending, show=False):
    455     '''
    456     @Description:
    457189
    458190    @Input:
     
    468200            documentation.
    469201
    470         flist
    471 
    472         fmetalist
    473 
    474         gdict: dict
    475             Contains basic informations of the ECMWF grib files, e.g.
    476             'Ni', 'Nj', 'latitudeOfFirstGridPointInDegrees',
    477             'longitudeOfFirstGridPointInDegrees',
    478             'latitudeOfLastGridPointInDegrees',
    479             'longitudeOfLastGridPointInDegrees',
    480             'jDirectionIncrementInDegrees',
    481             'iDirectionIncrementInDegrees'
    482 
    483         ftitle: string
    484             The titel of the plot.
    485 
    486         filename: string
    487             The plot is stored in a file with this name.
    488 
    489         fending: string
    490             Contains the type of plot, e.g. pdf or png
    491 
    492         show: boolean
    493             Decides if the plot is shown after plotting or hidden.
    494 
    495202    @Return:
    496203        <nothing>
    497204    '''
    498     print 'plotting map'
    499 
    500     t1 = time.time()
    501 
    502     fig = plt.figure(figsize=(12, 6.7))
    503     #mbaxes = fig.add_axes([0.05, 0.15, 0.8, 0.7])
    504 
    505     llx = gdict['longitudeOfFirstGridPointInDegrees'] #- 360
    506     if llx > 180. :
    507         llx -= 360.
    508     lly = gdict['latitudeOfLastGridPointInDegrees']
    509     dxout = gdict['iDirectionIncrementInDegrees']
    510     dyout = gdict['jDirectionIncrementInDegrees']
    511     urx = gdict['longitudeOfLastGridPointInDegrees']
    512     ury = gdict['latitudeOfFirstGridPointInDegrees']
    513     numxgrid = gdict['Ni']
    514     numygrid = gdict['Nj']
    515 
    516     m = Basemap(projection='cyl', llcrnrlon=llx, llcrnrlat=lly,
    517                 urcrnrlon=urx, urcrnrlat=ury,resolution='i')
    518 
    519     lw = 0.5
    520     m.drawmapboundary()
    521     x = linspace(llx, urx, numxgrid)
    522     y = linspace(lly, ury, numygrid)
    523 
    524     xx, yy = m(*meshgrid(x, y))
    525 
    526     #s = m.contourf(xx, yy, flist)
    527 
    528     s = plt.imshow(flist.T,
    529                     extent=(llx, urx, lly, ury),
    530                     alpha=1.0,
    531                     interpolation='nearest'
    532                     #vmin=vn,
    533                     #vmax=vx,
    534                     #cmap=my_cmap,
    535                     #levels=levels,
    536                     #cmap=my_cmap,
    537                     #norm=LogNorm(vn,vx)
    538                     )
    539 
    540     title(ftitle, y=1.08)
    541     cb = m.colorbar(s, location="right", pad="10%")
    542     #cb.set_label('Contribution per cell (ng m$^{-3}$)',size=14)
    543 
    544     thickline = np.arange(lly,ury+1,10.)
    545     thinline  = np.arange(lly,ury+1,5.)
    546     m.drawparallels(thickline,color='gray',dashes=[1,1],linewidth=0.5,labels=[1,1,1,1], xoffset=1.) # draw parallels
    547     m.drawparallels(np.setdiff1d(thinline,thickline),color='lightgray',dashes=[1,1],linewidth=0.5,labels=[0,0,0,0]) # draw parallels
    548 
    549     thickline = np.arange(llx,urx+1,10.)
    550     thinline  = np.arange(llx,urx+1,5.)
    551     m.drawmeridians(thickline,color='gray',dashes=[1,1],linewidth=0.5,labels=[1,1,1,1],yoffset=1.) # draw meridians
    552     m.drawmeridians(np.setdiff1d(thinline,thickline),color='lightgray',dashes=[1,1],linewidth=0.5,labels=[0,0,0,0]) # draw meridians
    553 
    554     m.drawcoastlines()
    555     m.drawcountries()
    556 
    557     plt.savefig(c.outputdir+'/'+filename+'_MAP.'+fending, facecolor=fig.get_facecolor(), edgecolor='none',format=fending)
    558     print 'created ', c.outputdir + '/' + filename
    559     if show == True:
    560         plt.show()
    561     fig.clf()
    562     plt.close(fig)
    563 
    564     print time.time() - t1, 's'
     205    start = datetime.datetime.strptime(c.start_date, '%Y%m%d%H')
     206    end = datetime.datetime.strptime(c.end_date, '%Y%m%d%H')
     207
     208    # create datelist between start and end date
     209    datelist = [start] # initialise datelist with first date
     210    run_date = start
     211    while run_date < end:
     212        run_date += datetime.timedelta(hours=int(c.dtime))
     213        datelist.append(run_date)
     214
     215    print 'datelist: ', datelist
     216
     217    c.paramIds = np.asarray(c.paramIds, dtype='int')
     218    c.levels = np.asarray(c.levels, dtype='int')
     219    c.area = np.asarray(c.area)
     220
     221    files = UioFiles(c.prefix+'*')
     222    files.list_files(c.inputdir)
     223    ifiles = get_files_per_date(files.files, datelist)
     224    ifiles.sort()
     225
     226    gdict = get_basics(ifiles[0], verb=False)
     227
     228    fdict = dict()
     229    fmeta = dict()
     230    fstamp = dict()
     231    for p in c.paramIds:
     232        for l in c.levels:
     233            key = '{:0>3}_{:0>3}'.format(p, l)
     234            fdict[key] = []
     235            fmeta[key] = []
     236            fstamp[key] = []
     237
     238    for filename in ifiles:
     239        f = open(filename)
     240        print "Opening file for reading data --- %s" % filename
     241        fdate = datetime.datetime.strptime(filename[-8:], "%y%m%d%H")
     242
     243        # Load in memory a grib message from a file.
     244        gid = codes_grib_new_from_file(f)
     245        while gid is not None:
     246            gtype = codes_get(gid, 'type')
     247            paramId = codes_get(gid, 'paramId')
     248            parameterName = codes_get(gid, 'parameterName')
     249            level = codes_get(gid, 'level')
     250
     251            if paramId in c.paramIds and level in c.levels:
     252                key = '{:0>3}_{:0>3}'.format(paramId, level)
     253                print 'key: ', key
     254                if fstamp[key]:
     255                    for i in range(len(fstamp[key])):
     256                        if fdate < fstamp[key][i]:
     257                            fstamp[key].insert(i, fdate)
     258                            fmeta[key].insert(i, [paramId, parameterName, gtype,
     259                                                  fdate, level])
     260                            fdict[key].insert(i, np.flipud(np.reshape(
     261                                codes_get_values(gid),
     262                                [gdict['Nj'], gdict['Ni']])))
     263                            break
     264                        elif fdate > fstamp[key][i] and i == len(fstamp[key])-1:
     265                            fstamp[key].append(fdate)
     266                            fmeta[key].append([paramId, parameterName, gtype,
     267                                               fdate, level])
     268                            fdict[key].append(np.flipud(np.reshape(
     269                                codes_get_values(gid),
     270                                [gdict['Nj'], gdict['Ni']])))
     271                            break
     272                        elif fdate > fstamp[key][i] and i != len(fstamp[key])-1 \
     273                             and fdate < fstamp[key][i+1]:
     274                            fstamp[key].insert(i, fdate)
     275                            fmeta[key].insert(i, [paramId, parameterName, gtype,
     276                                                  fdate, level])
     277                            fdict[key].insert(i, np.flipud(np.reshape(
     278                                codes_get_values(gid),
     279                                [gdict['Nj'], gdict['Ni']])))
     280                            break
     281                        else:
     282                            pass
     283                else:
     284                    fstamp[key].append(fdate)
     285                    fmeta[key].append((paramId, parameterName, gtype,
     286                                       fdate, level))
     287                    fdict[key].append(np.flipud(np.reshape(
     288                        codes_get_values(gid), [gdict['Nj'], gdict['Ni']])))
     289
     290            codes_release(gid)
     291
     292            # Load in memory a grib message from a file.
     293            gid = codes_grib_new_from_file(f)
     294
     295        f.close()
     296
     297    for k in fdict.iterkeys():
     298        print 'fmeta: ', len(fmeta), fmeta
     299        fml = fmeta[k]
     300        fdl = fdict[k]
     301        print 'fm1: ', len(fml), fml
     302        for fd, fm in zip(fdl, fml):
     303            print fm
     304            ftitle = fm[1] + ' {} '.format(fm[-1]) + \
     305                datetime.datetime.strftime(fm[3], '%Y%m%d%H')
     306            pname = '_'.join(fm[1].split()) + '_{}_'.format(fm[-1]) + \
     307                datetime.datetime.strftime(fm[3], '%Y%m%d%H')
     308            plot_map(c, fd, fm, gdict, ftitle, pname, 'png')
     309
     310    for k in fdict.iterkeys():
     311        fml = fmeta[k]
     312        fdl = fdict[k]
     313        fsl = fstamp[k]
     314        if fdl:
     315            fm = fml[0]
     316            fd = fdl[0]
     317            ftitle = fm[1] + ' {} '.format(fm[-1]) + \
     318                datetime.datetime.strftime(fm[3], '%Y%m%d%H')
     319            pname = '_'.join(fm[1].split()) + '_{}_'.format(fm[-1]) + \
     320                datetime.datetime.strftime(fm[3], '%Y%m%d%H')
     321            lat = -20.
     322            lon = 20.
     323            plot_timeseries(c, fdl, fml, fsl, lat, lon, gdict,
     324                            ftitle, pname, 'png')
    565325
    566326    return
    567327
    568 def getPlotArgs():
     328def plot_timeseries(c, flist, fmetalist, ftimestamps, lat, lon,
     329                    gdict, ftitle, filename, fending, show=False):
    569330    '''
    570331    @Description:
    571         Assigns the command line arguments and reads CONTROL file
    572         content. Apply default values for non mentioned arguments.
     332        Creates a timeseries plot for a given lat/lon position.
    573333
    574334    @Input:
    575         <nothing>
    576 
    577     @Return:
    578         args: instance of ArgumentParser
    579             Contains the commandline arguments from script/program call.
    580 
    581335        c: instance of class ControlFile
    582336            Contains all necessary information of a CONTROL file. The parameters
     
    589343            For more information about format and content of the parameter see
    590344            documentation.
     345
     346        flist: numpy array, 2d
     347            The actual data values to be plotted from the grib messages.
     348
     349        fmetalist: list of strings
     350            Contains some meta date for the data field to be plotted:
     351            parameter id, parameter Name, grid type, datetime, level
     352
     353        ftimestamps: list of datetime
     354            Contains the time stamps.
     355
     356        lat: float
     357            The latitude for which the timeseries should be plotted.
     358
     359        lon: float
     360            The longitude for which the timeseries should be plotted.
     361
     362        gdict: dict
     363            Contains basic informations of the ECMWF grib files, e.g.
     364            'Ni', 'Nj', 'latitudeOfFirstGridPointInDegrees',
     365            'longitudeOfFirstGridPointInDegrees',
     366            'latitudeOfLastGridPointInDegrees',
     367            'longitudeOfLastGridPointInDegrees',
     368            'jDirectionIncrementInDegrees',
     369            'iDirectionIncrementInDegrees'
     370
     371        ftitle: string
     372            The title of the timeseries.
     373
     374        filename: string
     375            The time series is stored in a file with this name.
     376
     377        fending: string
     378            Contains the type of plot, e.g. pdf or png
     379
     380        show: boolean
     381            Decides if the plot is shown after plotting or hidden.
     382
     383    @Return:
     384        <nothing>
     385    '''
     386    print 'plotting timeseries'
     387
     388    t1 = time.time()
     389
     390    #llx = gdict['longitudeOfFirstGridPointInDegrees']
     391    #if llx > 180. :
     392    #    llx -= 360.
     393    #lly = gdict['latitudeOfLastGridPointInDegrees']
     394    #dxout = gdict['iDirectionIncrementInDegrees']
     395    #dyout = gdict['jDirectionIncrementInDegrees']
     396    #urx = gdict['longitudeOfLastGridPointInDegrees']
     397    #ury = gdict['latitudeOfFirstGridPointInDegrees']
     398    #numxgrid = gdict['Ni']
     399    #numygrid = gdict['Nj']
     400
     401    farr = np.asarray(flist)
     402    #(time, lat, lon)
     403
     404    #lonindex = linspace(llx, urx, numxgrid)
     405    #latindex = linspace(lly, ury, numygrid)
     406
     407    ts = farr[:, 0, 0]
     408
     409    fig = plt.figure(figsize=(12, 6.7))
     410
     411    plt.plot(ftimestamps, ts)
     412    plt.title(ftitle)
     413
     414    plt.savefig(c.outputdir + '/' + filename + '_TS.' + fending,
     415                facecolor=fig.get_facecolor(),
     416                edgecolor='none',
     417                format=fending)
     418    print 'created ', c.outputdir + '/' + filename
     419    if show:
     420        plt.show()
     421    fig.clf()
     422    plt.close(fig)
     423
     424    print time.time() - t1, 's'
     425
     426    return
     427
     428def plot_map(c, flist, fmetalist, gdict, ftitle, filename, fending, show=False):
     429    '''
     430    @Description:
     431        Creates a basemap plot with imshow for a given data field.
     432
     433    @Input:
     434        c: instance of class ControlFile
     435            Contains all necessary information of a CONTROL file. The parameters
     436            are: DAY1, DAY2, DTIME, MAXSTEP, TYPE, TIME, STEP, CLASS, STREAM,
     437            NUMBER, EXPVER, GRID, LEFT, LOWER, UPPER, RIGHT, LEVEL, LEVELIST,
     438            RESOL, GAUSS, ACCURACY, OMEGA, OMEGADIFF, ETA, ETADIFF, DPDETA,
     439            SMOOTH, FORMAT, ADDPAR, WRF, CWC, PREFIX, ECSTORAGE, ECTRANS,
     440            ECFSDIR, MAILOPS, MAILFAIL, GRIB2FLEXPART, DEBUG, INPUTDIR,
     441            OUTPUTDIR, FLEXPART_ROOT_SCRIPTS
     442            For more information about format and content of the parameter see
     443            documentation.
     444
     445        flist: numpy array, 2d
     446            The actual data values to be plotted from the grib messages.
     447
     448        fmetalist: list of strings
     449            Contains some meta date for the data field to be plotted:
     450            parameter id, parameter Name, grid type, datetime, level
     451
     452        gdict: dict
     453            Contains basic informations of the ECMWF grib files, e.g.
     454            'Ni', 'Nj', 'latitudeOfFirstGridPointInDegrees',
     455            'longitudeOfFirstGridPointInDegrees',
     456            'latitudeOfLastGridPointInDegrees',
     457            'longitudeOfLastGridPointInDegrees',
     458            'jDirectionIncrementInDegrees',
     459            'iDirectionIncrementInDegrees'
     460
     461        ftitle: string
     462            The titel of the plot.
     463
     464        filename: string
     465            The plot is stored in a file with this name.
     466
     467        fending: string
     468            Contains the type of plot, e.g. pdf or png
     469
     470        show: boolean
     471            Decides if the plot is shown after plotting or hidden.
     472
     473    @Return:
     474        <nothing>
     475    '''
     476    print 'plotting map'
     477
     478    t1 = time.time()
     479
     480    fig = plt.figure(figsize=(12, 6.7))
     481    #mbaxes = fig.add_axes([0.05, 0.15, 0.8, 0.7])
     482
     483    llx = gdict['longitudeOfFirstGridPointInDegrees'] #- 360
     484    if llx > 180.:
     485        llx -= 360.
     486    lly = gdict['latitudeOfLastGridPointInDegrees']
     487    #dxout = gdict['iDirectionIncrementInDegrees']
     488    #dyout = gdict['jDirectionIncrementInDegrees']
     489    urx = gdict['longitudeOfLastGridPointInDegrees']
     490    ury = gdict['latitudeOfFirstGridPointInDegrees']
     491    #numxgrid = gdict['Ni']
     492    #numygrid = gdict['Nj']
     493
     494    m = Basemap(projection='cyl', llcrnrlon=llx, llcrnrlat=lly,
     495                urcrnrlon=urx, urcrnrlat=ury, resolution='i')
     496
     497    #lw = 0.5
     498    m.drawmapboundary()
     499    #x = linspace(llx, urx, numxgrid)
     500    #y = linspace(lly, ury, numygrid)
     501
     502    #xx, yy = m(*meshgrid(x, y))
     503
     504    #s = m.contourf(xx, yy, flist)
     505
     506    s = plt.imshow(flist.T,
     507                   extent=(llx, urx, lly, ury),
     508                   alpha=1.0,
     509                   interpolation='nearest'
     510                   #vmin=vn,
     511                   #vmax=vx,
     512                   #cmap=my_cmap,
     513                   #levels=levels,
     514                   #cmap=my_cmap,
     515                   #norm=LogNorm(vn,vx)
     516                  )
     517
     518    plt.title(ftitle, y=1.08)
     519    cb = m.colorbar(s, location="right", pad="10%")
     520    cb.set_label('label', size=14)
     521
     522    thickline = np.arange(lly, ury+1, 10.)
     523    thinline = np.arange(lly, ury+1, 5.)
     524    m.drawparallels(thickline,
     525                    color='gray',
     526                    dashes=[1, 1],
     527                    linewidth=0.5,
     528                    labels=[1, 1, 1, 1],
     529                    xoffset=1.)
     530    m.drawparallels(np.setdiff1d(thinline, thickline),
     531                    color='lightgray',
     532                    dashes=[1, 1],
     533                    linewidth=0.5,
     534                    labels=[0, 0, 0, 0])
     535
     536    thickline = np.arange(llx, urx+1, 10.)
     537    thinline = np.arange(llx, urx+1, 5.)
     538    m.drawmeridians(thickline,
     539                    color='gray',
     540                    dashes=[1, 1],
     541                    linewidth=0.5,
     542                    labels=[1, 1, 1, 1],
     543                    yoffset=1.)
     544    m.drawmeridians(np.setdiff1d(thinline, thickline),
     545                    color='lightgray',
     546                    dashes=[1, 1],
     547                    linewidth=0.5,
     548                    labels=[0, 0, 0, 0])
     549
     550    m.drawcoastlines()
     551    m.drawcountries()
     552
     553    plt.savefig(c.outputdir + '/' + filename + '_MAP.' + fending,
     554                facecolor=fig.get_facecolor(),
     555                edgecolor='none',
     556                format=fending)
     557    print 'created ', c.outputdir + '/' + filename
     558    if show:
     559        plt.show()
     560    fig.clf()
     561    plt.close(fig)
     562
     563    print time.time() - t1, 's'
     564
     565    return
     566
     567def get_plot_args():
     568    '''
     569    @Description:
     570        Assigns the command line arguments and reads CONTROL file
     571        content. Apply default values for non mentioned arguments.
     572
     573    @Input:
     574        <nothing>
     575
     576    @Return:
     577        args: instance of ArgumentParser
     578            Contains the commandline arguments from script/program call.
     579
     580        c: instance of class ControlFile
     581            Contains all necessary information of a CONTROL file. The parameters
     582            are: DAY1, DAY2, DTIME, MAXSTEP, TYPE, TIME, STEP, CLASS, STREAM,
     583            NUMBER, EXPVER, GRID, LEFT, LOWER, UPPER, RIGHT, LEVEL, LEVELIST,
     584            RESOL, GAUSS, ACCURACY, OMEGA, OMEGADIFF, ETA, ETADIFF, DPDETA,
     585            SMOOTH, FORMAT, ADDPAR, WRF, CWC, PREFIX, ECSTORAGE, ECTRANS,
     586            ECFSDIR, MAILOPS, MAILFAIL, GRIB2FLEXPART, DEBUG, INPUTDIR,
     587            OUTPUTDIR, FLEXPART_ROOT_SCRIPTS
     588            For more information about format and content of the parameter see
     589            documentation.
    591590    '''
    592591    parser = ArgumentParser(description='Plot retrieved GRIB data from ' + \
     
    597596    parser.add_argument("--start_date", dest="start_date",
    598597                        help="start date YYYYMMDD")
    599     parser.add_argument( "--end_date", dest="end_date",
    600                          help="end_date YYYYMMDD")
     598    parser.add_argument("--end_date", dest="end_date",
     599                        help="end_date YYYYMMDD")
    601600
    602601    parser.add_argument("--start_step", dest="start_step",
    603602                        help="start step in hours")
    604     parser.add_argument( "--end_step", dest="end_step",
    605                          help="end step in hours")
     603    parser.add_argument("--end_step", dest="end_step",
     604                        help="end step in hours")
    606605
    607606# some arguments that override the default in the CONTROL file
     
    635634    except IOError:
    636635        try:
    637             c = ControlFile(localpythonpath + args.controlfile)
    638 
    639         except:
     636            c = ControlFile(LOCAL_PYTHON_PATH + args.controlfile)
     637        except IOError:
    640638            print 'Could not read CONTROL file "' + args.controlfile + '"'
    641639            print 'Either it does not exist or its syntax is wrong.'
     
    682680if __name__ == "__main__":
    683681    main()
    684 
  • python/prepare_flexpart.py

    r991df6a rff99eae  
    22# -*- coding: utf-8 -*-
    33#************************************************************************
    4 # TODO AP
    5 # - wieso cleanup in main wenn es in prepareflexpart bereits abgefragt wurde?
    6 #   doppelt gemoppelt?
     4# ToDo AP
    75# - wieso start=startm1 wenn basetime = 0 ?  wenn die fluxes nicht mehr
    86#   relevant sind? verstehe ich nicht
     
    2927#        - added documentation
    3028#        - minor changes in programming style for consistence
    31 #        - BUG: removed call of cleanup-Function after call of
     29#        - BUG: removed call of clean_up-Function after call of
    3230#               prepareFlexpart in main since it is already called in
    3331#               prepareFlexpart at the end!
    3432#        - created function main and moved the two function calls for
    35 #          arguments and prepareFLEXPART into it
     33#          arguments and prepare_flexpart into it
    3634#
    3735# @License:
     
    4442#    This program prepares the final version of the grib files which are
    4543#    then used by FLEXPART. It converts the bunch of grib files extracted
    46 #    via getMARSdata by doing for example the necessary conversion to get
     44#    via get_mars_data by doing for example the necessary conversion to get
    4745#    consistent grids or the disaggregation of flux data. Finally, the
    4846#    program combines the data fields in files per available hour with the
     
    5250# @Program Content:
    5351#    - main
    54 #    - prepareFLEXPART
     52#    - prepare_flexpart
    5553#
    5654#*******************************************************************************
     
    5957# MODULES
    6058# ------------------------------------------------------------------------------
    61 import shutil
    6259import datetime
    63 #import time
    6460import os
    6561import inspect
    6662import sys
    6763import socket
    68 from argparse import ArgumentParser, ArgumentDefaultsHelpFormatter
    6964
    70 hostname = socket.gethostname()
    71 ecapi = 'ecmwf' not in hostname
     65# software specific classes and modules from flex_extract
     66from UioFiles import UioFiles
     67from tools import interpret_args_and_control, clean_up
     68from EcFlexpart import EcFlexpart
     69
     70ecapi = 'ecmwf' not in socket.gethostname()
    7271try:
    7372    if ecapi:
     
    7776
    7877# add path to pythonpath so that python finds its buddies
    79 localpythonpath = os.path.dirname(os.path.abspath(
     78LOCAL_PYTHON_PATH = os.path.dirname(os.path.abspath(
    8079    inspect.getfile(inspect.currentframe())))
    81 if localpythonpath not in sys.path:
    82     sys.path.append(localpythonpath)
     80if LOCAL_PYTHON_PATH not in sys.path:
     81    sys.path.append(LOCAL_PYTHON_PATH)
    8382
    84 # software specific classes and modules from flex_extract
    85 from UIOFiles import UIOFiles
    86 from Tools import interpret_args_and_control, cleanup
    87 from ECFlexpart import ECFlexpart
    8883# ------------------------------------------------------------------------------
    8984# FUNCTION
     
    9287    '''
    9388    @Description:
    94         If prepareFLEXPART is called from command line, this function controls
     89        If prepare_flexpart is called from command line, this function controls
    9590        the program flow and calls the argumentparser function and
    96         the prepareFLEXPART function for preparation of GRIB data for FLEXPART.
     91        the prepare_flexpart function for preparation of GRIB data for FLEXPART.
    9792
    9893    @Input:
     
    10398    '''
    10499    args, c = interpret_args_and_control()
    105     prepareFLEXPART(args, c)
     100    prepare_flexpart(args, c)
    106101
    107102    return
    108103
    109 def prepareFLEXPART(args, c):
     104def prepare_flexpart(args, c):
    110105    '''
    111106    @Description:
    112         Lists all grib files retrieved from MARS with getMARSdata and
     107        Lists all grib files retrieved from MARS with get_mars_data and
    113108        uses prepares data for the use in FLEXPART. Specific data fields
    114109        are converted to a different grid and the flux data are going to be
     
    157152    # one day after the end date is needed
    158153    startm1 = start - datetime.timedelta(days=1)
    159     endp1 = end + datetime.timedelta(days=1)
     154#    endp1 = end + datetime.timedelta(days=1)
    160155
    161156    # get all files with flux data to be deaccumulated
    162     inputfiles = UIOFiles('*OG_acc_SL*.' + c.ppid + '.*')
    163     inputfiles.listFiles(c.inputdir)
     157    inputfiles = UioFiles('*OG_acc_SL*.' + c.ppid + '.*')
     158    inputfiles.list_files(c.inputdir)
    164159
    165160    # create output dir if necessary
     
    168163
    169164    # deaccumulate the flux data
    170     flexpart = ECFlexpart(c, fluxes=True)
     165    flexpart = EcFlexpart(c, fluxes=True)
    171166    flexpart.write_namelist(c, 'fort.4')
    172167    flexpart.deacc_fluxes(inputfiles, c)
    173168
    174     print('Prepare ' + start.strftime("%Y%m%d") +
    175           "/to/" + end.strftime("%Y%m%d"))
     169    print 'Prepare ' + start.strftime("%Y%m%d") + \
     170          "/to/" + end.strftime("%Y%m%d")
    176171
    177172    # get a list of all files from the root inputdir
    178     inputfiles = UIOFiles('????__??.*' + c.ppid + '.*')
    179     inputfiles.listFiles(c.inputdir)
     173    inputfiles = UioFiles('????__??.*' + c.ppid + '.*')
     174    inputfiles.list_files(c.inputdir)
    180175
    181176    # produce FLEXPART-ready GRIB files and
     
    185180        start = startm1
    186181
    187     flexpart = ECFlexpart(c, fluxes=False)
     182    flexpart = EcFlexpart(c, fluxes=False)
    188183    flexpart.create(inputfiles, c)
    189184    flexpart.process_output(c)
     
    192187    # otherwise delete temporary files
    193188    if int(c.debug) != 0:
    194         print('Temporary files left intact')
     189        print 'Temporary files left intact'
    195190    else:
    196         cleanup(c)
     191        clean_up(c)
    197192
    198193    return
  • python/profiling.py

    r991df6a rff99eae  
    22# -*- coding: utf-8 -*-
    33#************************************************************************
    4 # TODO AP
     4# ToDo AP
    55# - check of license of book content
    66#************************************************************************
     
    6666        result = fn(*args, **kwargs)
    6767        t2 = time.time()
    68         print("@timefn:" + fn.func_name + " took " + str(t2 - t1) + " seconds")
     68        print "@timefn:" + fn.func_name + " took " + str(t2 - t1) + " seconds"
    6969
    7070        return result
  • python/submit.py

    r991df6a rff99eae  
    11#!/usr/bin/env python
    22# -*- coding: utf-8 -*-
    3 #************************************************************************
    4 # TODO AP
    5 #
    6 # - Change History ist nicht angepasst ans File! Original geben lassen
    7 # - dead code ? what to do?
    8 # - seperate operational and reanlysis for clarification
    9 # - divide in two submits , ondemand und operational
    10 # -
    11 #************************************************************************
    12 
    133#*******************************************************************************
    144# @Author: Anne Fouilloux (University of Oslo)
     
    3727#    program flow.
    3828#    If it is supposed to work locally then it works through the necessary
    39 #    functions getMARSdata and prepareFlexpart. Otherwise it prepares
     29#    functions get_mars_data and prepareFlexpart. Otherwise it prepares
    4030#    a shell job script which will do the necessary work on the
    4131#    ECMWF server and is submitted via ecaccess-job-submit.
     
    5242import os
    5343import sys
    54 import glob
    5544import subprocess
    5645import inspect
    57 # add path to pythonpath so that python finds its buddies
    58 localpythonpath = os.path.dirname(os.path.abspath(
    59     inspect.getfile(inspect.currentframe())))
    60 if localpythonpath not in sys.path:
    61     sys.path.append(localpythonpath)
    6246
    6347# software specific classes and modules from flex_extract
    64 from Tools import interpret_args_and_control, normalexit
    65 from getMARSdata import getMARSdata
    66 from prepareFLEXPART import prepareFLEXPART
     48from tools import interpret_args_and_control, normal_exit
     49from get_mars_data import get_mars_data
     50from prepare_flexpart import prepare_flexpart
     51
     52# add path to pythonpath so that python finds its buddies
     53LOCAL_PYTHON_PATH = os.path.dirname(os.path.abspath(
     54    inspect.getfile(inspect.currentframe())))
     55if LOCAL_PYTHON_PATH not in sys.path:
     56    sys.path.append(LOCAL_PYTHON_PATH)
    6757
    6858# ------------------------------------------------------------------------------
     
    8474        <nothing>
    8575    '''
    86     calledfromdir = os.getcwd()
     76
     77    called_from_dir = os.getcwd()
    8778    args, c = interpret_args_and_control()
     79
    8880    # on local side
    8981    if args.queue is None:
    9082        if c.inputdir[0] != '/':
    91             c.inputdir = os.path.join(calledfromdir, c.inputdir)
     83            c.inputdir = os.path.join(called_from_dir, c.inputdir)
    9284        if c.outputdir[0] != '/':
    93             c.outputdir = os.path.join(calledfromdir, c.outputdir)
    94         getMARSdata(args, c)
    95         prepareFLEXPART(args, c)
    96         normalexit(c)
     85            c.outputdir = os.path.join(called_from_dir, c.outputdir)
     86        get_mars_data(args, c)
     87        prepare_flexpart(args, c)
     88        normal_exit(c)
    9789    # on ECMWF server
    9890    else:
     
    139131
    140132    # put all parameters of ControlFile instance into a list
    141     clist = c.tolist()
     133    clist = c.to_list() # ondemand
    142134    colist = []  # operational
    143135    mt = 0
    144136
    145 #AP wieso 2 for loops?
    146137    for elem in clist:
    147138        if 'maxstep' in elem:
     
    152143            elem = 'start_date ' + '${MSJ_YEAR}${MSJ_MONTH}${MSJ_DAY}'
    153144        if 'end_date' in elem:
    154 #AP Fehler?! Muss end_date heissen
    155             elem = 'start_date ' + '${MSJ_YEAR}${MSJ_MONTH}${MSJ_DAY}'
     145            elem = 'end_date ' + '${MSJ_YEAR}${MSJ_MONTH}${MSJ_DAY}'
    156146        if 'base_time' in elem:
    157147            elem = 'base_time ' + '${MSJ_BASETIME}'
     
    173163        p = subprocess.check_call(['ecaccess-job-submit', '-queueName',
    174164                                   queue, 'job.ksh'])
    175     except:
    176         print('ecaccess-job-submit failed, probably eccert has expired')
     165    except subprocess.CalledProcessError as e:
     166        print 'ecaccess-job-submit failed!'
     167        print 'Error Message: '
     168        print e.output
    177169        exit(1)
    178170
    179     print('You should get an email with subject flex.hostname.pid')
     171    print 'You should get an email with subject flex.hostname.pid'
    180172
    181173    return
  • python/test_suite.py

    r991df6a rff99eae  
    22# -*- coding: utf-8 -*-
    33#************************************************************************
    4 # TODO AP
    5 #
     4# ToDo AP
    65# - provide more tests
    76# - provide more documentation
    8 # -
    97#************************************************************************
    108
     
    2826# @Program Functionality:
    2927#    This script triggers the ECMWFDATA test suite. Call with
    30 #    testsuite.py [test group]
     28#    test_suite.py [test group]
    3129#
    3230# @Program Content:
     
    4644# ------------------------------------------------------------------------------
    4745try:
    48     taskfile = open('testsuite.json')
    49 except:
    50     print 'could not open suite definition file testsuite.json'
     46    taskfile = open('test_suite.json')
     47except IOError:
     48    print 'could not open suite definition file test_suite.json'
    5149    exit()
    5250
     
    7371    try:
    7472        tk, tv = g, tasks[g]
    75     except:
    76         continue
     73    finally:
     74        pass
    7775    garglist = []
    7876    for ttk, ttv in tv.iteritems():
     
    8078            if ttk != 'script':
    8179                garglist.append('--' + ttk)
    82                 if '$' == ttv[0]:
     80                if ttv[0] == '$':
    8381                    garglist.append(os.path.expandvars(ttv))
    8482                else:
     
    8987            for tttk, tttv in ttv.iteritems():
    9088                if isinstance(tttv, basestring):
    91                         arglist.append('--' + tttk)
    92                         if '$' in tttv[0]:
    93                             arglist.append(os.path.expandvars(tttv))
    94                         else:
    95                             arglist.append(tttv)
     89                    arglist.append('--' + tttk)
     90                    if '$' in tttv[0]:
     91                        arglist.append(os.path.expandvars(tttv))
     92                    else:
     93                        arglist.append(tttv)
    9694            print 'Command: ', ' '.join([tv['script']] + garglist + arglist)
    9795            o = '../test/' + tk + '_' + ttk + '_' + '_'.join(ttv.keys())
     
    10199                p = subprocess.check_call([tv['script']] + garglist + arglist,
    102100                                          stdout=f, stderr=f)
    103             except:
     101            except subprocess.CalledProcessError as e:
    104102                f.write('\nFAILED\n')
    105103                print 'FAILED'
     
    111109print str(jobcounter-jobfailed) + ' successful, ' + str(jobfailed) + ' failed'
    112110print 'If tasks have been submitted via ECACCESS please check emails'
    113 
    114 #print tasks
  • python/tools.py

    r991df6a rff99eae  
    22# -*- coding: utf-8 -*-
    33#************************************************************************
    4 # TODO AP
    5 # - check myerror
    6 # - check normalexit
    7 # - check getListAsString
     4# ToDo AP
     5# - check my_error
     6# - check normal_exit
     7# - check get_list_as_string
    88# - seperate args and control interpretation
    99#************************************************************************
     
    1515# @Change History:
    1616#    October 2014 - Anne Fouilloux (University of Oslo)
    17 #        - created functions silentremove and product (taken from ECMWF)
     17#        - created functions silent_remove and product (taken from ECMWF)
    1818#
    1919#    November 2015 - Leopold Haimberger (University of Vienna)
    20 #        - created functions: interpret_args_and_control, cleanup
    21 #          myerror, normalexit, init128, toparamId
     20#        - created functions: interpret_args_and_control, clean_up
     21#          my_error, normal_exit, init128, to_param_id
    2222#
    2323#    April 2018 - Anne Philipp (University of Vienna):
    2424#        - applied PEP8 style guide
    2525#        - added documentation
    26 #        - moved all functions from file FlexpartTools to this file Tools
    27 #        - added function getListAsString
     26#        - moved all functions from file Flexparttools to this file tools
     27#        - added function get_list_as_string
    2828#
    2929# @License:
     
    3939# @Module Content:
    4040#    - interpret_args_and_control
    41 #    - cleanup
    42 #    - myerror
    43 #    - normalexit
     41#    - clean_up
     42#    - my_error
     43#    - normal_exit
    4444#    - product
    45 #    - silentremove
     45#    - silent_remove
    4646#    - init128
    47 #    - toparamId
    48 #    - getListAsString
     47#    - to_param_id
     48#    - get_list_as_string
    4949#
    5050#*******************************************************************************
     
    5757import sys
    5858import glob
     59import inspect
     60import subprocess
    5961import traceback
    60 from numpy import *
    61 from gribapi import *
    6262from argparse import ArgumentParser, ArgumentDefaultsHelpFormatter
     63import numpy as np
    6364
    6465# software specific class from flex_extract
     
    123124    parser.add_argument("--flexpart_root_scripts", dest="flexpart_root_scripts",
    124125                        help="FLEXPART root directory (to find grib2flexpart \
    125                         and COMMAND file)\n\ Normally ECMWFDATA resides in \
     126                        and COMMAND file)\n Normally ECMWFDATA resides in \
    126127                        the scripts directory of the FLEXPART distribution")
    127128
    128     # this is only used by prepareFLEXPART.py to rerun a postprocessing step
     129    # this is only used by prepare_flexpart.py to rerun a postprocessing step
    129130    parser.add_argument("--ppid", dest="ppid",
    130131                        help="Specify parent process id for \
    131                         rerun of prepareFLEXPART")
     132                        rerun of prepare_flexpart")
    132133
    133134    # arguments for job submission to ECMWF, only needed by submit.py
     
    152153    except IOError:
    153154        try:
    154             c = ControlFile(localpythonpath + args.controlfile)
    155         except:
    156             print('Could not read CONTROL file "' + args.controlfile + '"')
    157             print('Either it does not exist or its syntax is wrong.')
    158             print('Try "' + sys.argv[0].split('/')[-1] +
    159                   ' -h" to print usage information')
     155            LOCAL_PYTHON_PATH = os.path.dirname(os.path.abspath(
     156                inspect.getfile(inspect.currentframe())))
     157            c = ControlFile(LOCAL_PYTHON_PATH + args.controlfile)
     158        except IOError:
     159            print 'Could not read CONTROL file "' + args.controlfile + '"'
     160            print 'Either it does not exist or its syntax is wrong.'
     161            print 'Try "' + sys.argv[0].split('/')[-1] + \
     162                  ' -h" to print usage information'
    160163            exit(1)
    161164
    162165    # check for having at least a starting date
    163166    if  args.start_date is None and getattr(c, 'start_date') is None:
    164         print('start_date specified neither in command line nor \
    165                in CONTROL file ' + args.controlfile)
    166         print('Try "' + sys.argv[0].split('/')[-1] +
    167               ' -h" to print usage information')
     167        print 'start_date specified neither in command line nor \
     168               in CONTROL file ' + args.controlfile
     169        print 'Try "' + sys.argv[0].split('/')[-1] + \
     170              ' -h" to print usage information'
    168171        exit(1)
    169172
     
    206209        l = args.area.split('/')
    207210        if afloat:
    208             for i in range(len(l)):
    209                 l[i] = str(int(float(l[i]) * 1000))
     211            for i, item in enumerate(l):
     212                item = str(int(float(item) * 1000))
    210213        c.upper, c.left, c.lower, c.right = l
    211214
     
    218221        if 'to' in args.step.lower():
    219222            if 'by' in args.step.lower():
    220                 ilist = arange(int(l[0]), int(l[2]) + 1, int(l[4]))
     223                ilist = np.arange(int(l[0]), int(l[2]) + 1, int(l[4]))
    221224                c.step = ['{:0>3}'.format(i) for i in ilist]
    222225            else:
    223                 myerror(None, args.step + ':\n' +
    224                         'please use "by" as well if "to" is used')
     226                my_error(None, args.step + ':\n' +
     227                         'please use "by" as well if "to" is used')
    225228        else:
    226229            c.step = l
     
    239242
    240243
    241 def cleanup(c):
     244def clean_up(c):
    242245    '''
    243246    @Description:
     
    263266    '''
    264267
    265     print("cleanup")
     268    print "clean_up"
    266269
    267270    cleanlist = glob.glob(c.inputdir + "/*")
    268271    for cl in cleanlist:
    269272        if c.prefix not in cl:
    270             silentremove(cl)
     273            silent_remove(cl)
    271274        if c.ecapi is False and (c.ectrans == '1' or c.ecstorage == '1'):
    272             silentremove(cl)
    273 
    274     print("Done")
     275            silent_remove(cl)
     276
     277    print "Done"
    275278
    276279    return
    277280
    278281
    279 def myerror(c, message='ERROR'):
     282def my_error(c, message='ERROR'):
    280283    '''
    281284    @Description:
     
    303306        <nothing>
    304307    '''
    305     # uncomment if user wants email notification directly from python
    306     #try:
    307         #target = c.mailfail
    308     #except AttributeError:
    309         #target = os.getenv('USER')
    310 
    311     #if(type(target) is not list):
    312         #target = [target]
    313 
    314     print(message)
    315 
    316     # uncomment if user wants email notification directly from python
    317     #for t in target:
    318     #p = subprocess.Popen(['mail','-s ECMWFDATA v7.0 ERROR', os.path.expandvars(t)],
    319     #                     stdin = subprocess.PIPE, stdout = subprocess.PIPE,
    320     #                     stderr = subprocess.PIPE, bufsize = 1)
    321     #tr = '\n'.join(traceback.format_stack())
    322     #pout = p.communicate(input = message+'\n\n'+tr)[0]
    323     #print 'Email sent to '+os.path.expandvars(t) # +' '+pout.decode()
     308
     309    print message
     310
     311    # comment if user does not want email notification directly from python
     312    try:
     313        target = []
     314        target.extend(c.mailfail)
     315    except AttributeError:
     316        target = []
     317        target.extend(os.getenv('USER'))
     318
     319    for t in target:
     320        p = subprocess.Popen(['mail', '-s ECMWFDATA v7.0 ERROR',
     321                              os.path.expandvars(t)],
     322                             stdin=subprocess.PIPE,
     323                             stdout=subprocess.PIPE,
     324                             stderr=subprocess.PIPE,
     325                             bufsize=1)
     326        tr = '\n'.join(traceback.format_stack())
     327        pout = p.communicate(input=message + '\n\n' + tr)[0]
     328        print 'Email sent to ' + os.path.expandvars(t) + ' ' + pout.decode()
    324329
    325330    exit(1)
     
    328333
    329334
    330 def normalexit(c, message='Done!'):
     335def normal_exit(c, message='Done!'):
    331336    '''
    332337    @Description:
     
    354359
    355360    '''
    356     # Uncomment if user wants notification directly from python
    357     #try:
    358         #target = c.mailops
    359         #if(type(target) is not list):
    360             #target = [target]
    361         #for t in target:
    362             #p = subprocess.Popen(['mail','-s ECMWFDATA v7.0 normal exit',
    363             #                      os.path.expandvars(t)],
    364             #                      stdin = subprocess.PIPE,
    365             #                      stdout = subprocess.PIPE,
    366             #                      stderr = subprocess.PIPE, bufsize = 1)
    367             #pout = p.communicate(input = message+'\n\n')[0]
    368             #print pout.decode()
    369     #except:
    370         #pass
    371 
    372     print(message)
     361    print message
     362
     363    # comment if user does not want notification directly from python
     364    try:
     365        target = []
     366        target.extend(c.mailops)
     367        for t in target:
     368            p = subprocess.Popen(['mail', '-s ECMWFDATA v7.0 normal exit',
     369                                  os.path.expandvars(t)],
     370                                 stdin=subprocess.PIPE,
     371                                 stdout=subprocess.PIPE,
     372                                 stderr=subprocess.PIPE,
     373                                 bufsize=1)
     374            pout = p.communicate(input=message+'\n\n')[0]
     375            print 'Email sent to ' + os.path.expandvars(t) + ' ' + pout.decode()
     376    finally:
     377        pass
    373378
    374379    return
     
    411416
    412417
    413 def silentremove(filename):
     418def silent_remove(filename):
    414419    '''
    415420    @Description:
     
    435440
    436441
    437 def init128(fn):
     442def init128(filepath):
    438443    '''
    439444    @Description:
     
    441446
    442447    @Input:
    443         fn: string
     448        filepath: string
    444449            Path to file of ECMWF grib table number 128.
    445450
     
    451456    '''
    452457    table128 = dict()
    453     with open(fn) as f:
     458    with open(filepath) as f:
    454459        fdata = f.read().split('\n')
    455460    for data in fdata:
     
    460465
    461466
    462 def toparamId(pars, table):
     467def to_param_id(pars, table):
    463468    '''
    464469    @Description:
     
    486491    ipar = []
    487492    for par in cpar:
    488         found = False
    489493        for k, v in table.iteritems():
    490494            if par == k or par == v:
    491495                ipar.append(int(k))
    492                 found = True
    493496                break
    494         if found is False:
    495             print('Warning: par ' + par + ' not found in table 128')
     497        else:
     498            print 'Warning: par ' + par + ' not found in table 128'
    496499
    497500    return ipar
    498501
    499 def getListAsString(listObj):
     502def get_list_as_string(list_obj, concatenate_sign=', '):
    500503    '''
    501504    @Description:
     
    503506
    504507    @Input:
    505         listObj: list
     508        list_obj: list
    506509            A list with arbitrary content.
    507510
    508     @Return:
    509         strOfList: string
     511        concatenate_sign: string, optional
     512            A string which is used to concatenate the single
     513            list elements. Default value is ", ".
     514
     515    @Return:
     516        str_of_list: string
    510517            The content of the list as a single string.
    511518    '''
    512519
    513     strOfList = ", ".join( str(l) for l in listObj)
    514 
    515     return strOfList
     520    str_of_list = concatenate_sign.join(str(l) for l in list_obj)
     521
     522    return str_of_list
Note: See TracChangeset for help on using the changeset viewer.
hosted by ZAMG