Changeset 1845 in lmdz_wrf for trunk/tools/nc_var_tools.py

Timestamp:

Mar 22, 2018, 9:30:10 PM (7 years ago)

Author:

lfita

Message:

Introducing into the flow of:

• `CFmorzization': Function to provide a CF-compilation version of a variable within a file

File:

• trunk/tools/nc_var_tools.py (modified) (9 diffs)

trunk/tools/nc_var_tools.py

@@ -30 +30 @@
 # basicvardef: Function to give the basic attributes to a variable (std_name, long_name, units)
 # CDO_toCF: Function to pass a CDO output file to CF-conventions
+# CFmorzization: Function to provide a CF-compilation version of a variable within a file
 # changevartype: Function to change the type of a variable (when possible)
 # chdimname: Changing the name of the dimension
@@ -22676 +22677 @@
         - Christensen, Gutowski, Nikulin, and Legutke; 2013. CORDEX Archive Design 
           (http://cordex.org/publications/report-and-document-archives/)
-        - Taylor, K., and C. Doutriaux, 2011. “CMIP5 Model Output Requirements: File 
-          Contents and Format, Data Structure and Meta Data” 
+        - Taylor, K., and C. Doutriaux, 2011. CMIP5 Model Output Requirements: File 
+          Contents and Format, Data Structure and Meta Data
           (http://cmip-pcmdi.llnl.gov/cmip5/output_req.html#metadata)
 
-        values=[dimvarns]:[globattrfile]
+        values=[dimvarns]:[globattrfile]:[projectfile]
           [dimvarns]: ',' separated list for identification of axes, dimensions and 
             variable-dimensions as '[AXIS]|[dimn]|[vardimn]'
@@ -22693 +22694 @@
               institution  [value]
               model_id [value] (CORDEX starts with <institute_id>-)
-              experiment_id     [value] (same CV as CMIP5 with "evaluation" in addition)
+              experiment_id [value] (same CV as CMIP5 with "evaluation" in addition)
               experiment [value] (long version of "experiment_id")
               contact [value] ([Name], [email] of contact person)
-              product   [value] (single value: "output")
+              product [value] (single value: "output")
               Conventions [value] (later than CF1.4 if present / single value: "CF-1.4")
-              creation_date     [YYYY-MM-DDTHH:MM:SSZ] (format according to UNIDATA conventions, same in CMIP5 & CORDEX)
-              frequency [value] (same CV as CMIP5 with "sem" in addition; 3hr highest freq. in CORDEX; no "yr")
+              creation_date [YYYY-MM-DDTHH:MM:SSZ] (format according to UNIDATA conventions, same in CMIP5 & CORDEX)
+              frequency [value] (same CV as CMIP5 with "sem" in addition; 3hr highest freq. in CORDEX; no "yr")
               driving_model_id [value] (CMIP5 institute_id, <CMIP5 institute_id>-<member in CMIP5 CV list of model_id>)
               driving_model_ensemble_member [value] (CMIP5 CV, e.g. "r1i1p1")
@@ -22706 +22707 @@
                 members of the CMIP5 list of forcing agents are not necessarily forcing agents of RCM (besides implicitly
                 in driving_experimnt)
-              rcm_version_ID [value] (<free text string>; valid characters only indicates model modifs during the project
+              rcm_version_ID [value] (<free text string>; valid characters only indicates model modifs during the project
                 e.g.: parameterizations, small upgrades)
               project_id [value] (single value)
-              CORDEX_domain     [value] ((e.g.: "AFR-44", "AFR-44i", includes resolution acronym)
+              CORDEX_domain [value] ((e.g.: "AFR-44", "AFR-44i", includes resolution acronym)
               tracking_id [value] (not required or explained in ADD; it is, however, strongly recommended to include it 
                 as prescribed CMIP5.)
-              basetime  [value] (reference time "1949-12-01T00:00:00Z" formatted acording to UNIDATA conventions)
-              calendar  [value] (CF conventions as in driving model)
+              basetime [value] (reference time "1949-12-01T00:00:00Z" formatted acording to UNIDATA conventions)
+              calendar [value] (CF conventions as in driving model)
               grid [value] (according to CF conventions
+          [projectfile]: ASCII file with all the CF-mandatory projection information. A line for each attribute. Name of 
+            the variable to be used to refer the projection inside the file will be the same as global attribute `grid'
+            See more information from: http://cfconventions.org/Data/cf-conventions/cf-conventions-1.7/build/ch05s06.html
+             (...) (specific projection attributes and values)
+             [name] [value] [kind] ('I': integer, 'R': float, 'D': double float, 'S': string)
         ncfile= netCDF file to use
         variable= ',' list of variables to CF transform ('all' for all variables). A separated file will be created
@@ -22724 +22730 @@
     if values == 'h':
         print fname + '_____________________________________________________________'
-        print CFmorzizarion.__doc__
+        print CFmorzization.__doc__
         quit()
 
-    expectargs = '[dimvarns]:[globattrfilen]'
+    expectargs = '[dimvarns]:[globattrfilen]:[projectfile]'
     gen.check_arguments(fname, values, expectargs, ':')
 
     dimvarns = values.split(':')[0].split(',')
     globalattrfilen = values.split(':')[1]
+    projfilen = values.split(':')[2]
 
     # CF Mandatory global attributes
@@ -22770 +22777 @@
             quit(-1)
 
+    # Getting projection
+    if not os.path.isfile(projfilen):
+        print errormsg
+        print '  '+fname+ "': ASCII file with projection information '" + projfilen+ \
+          + "' does not exist !!"
+        quit(-1)
+
+    pinf = {}
+    opinff = open(projfilen, 'r')
+    for line in opinff:
+        if line[0:1] != '#' and len(line) > 1:
+            linevals = line.replace('\n','').replace('\t','').split(' ')
+            pinf[linevals[0]] = [linevals[1].replace('!',' '), linevals[2]]
+    opinff.close()
+
     # Getting variables to process
     onc = NetCDFFile(ncfile,'r')
@@ -22787 +22809 @@
     # VariableName: from `variables_values.dat'
     # Domain: CORDEX_domain
-    # GCMModelName:     driving_model_id
+    # GCMModelName: driving_model_id
     # CMIP5ExperimentName: driving_experiment_name
     # CMIP5EnsembleMember: driving_model_ensemble_member
@@ -22797 +22819 @@
     #   CMIP5EnsembleMember_RCMModelName_RCMVersionID_Frequency_StartTime-EndTime.nc
 
+    # Expected axes
+    axes = ['X', 'Y', 'Z', 'T']
+
+    oncvars = onc.variables.keys()
+    oncvars.sort()
+
+    # Axes information
+    CFdimvals = {}
+    filedimvals = {}
+    CFaxisvardimvals = {}
+    for axn in axisinf.keys():
+        axisv = axisinf[axn]
+        dimn = axisv[0]
+        vardimn = axisv[1]
+
+        # Axis values from file
+        if vardimn != 'WRFtime':
+            if not onc.variables.has_key(vardimn):
+                print errormsg
+                print '  '+fname+"': netcdf file '"+ ncfile + "' does not contain "+ \
+                  "variable '" + vardimn + "' !!"
+                print '    available ones:', oncvars
+                onc.close()
+                quit(-1)
+            oaxisv = onc.variables[vardimn]
+        else:
+            print infmsg
+            print '    ' + fname + ": creation of variable 'time' from WRF 'Times' !!"
+            odimvar = onc.variables['Times']
+            timewrfv = odimvar[:]
+            refT = gen.datetimeStr_conversion(gattr['basetime'], 'Y-m-dTH:M:SZ',     \
+              'YmdHMS')
+            tvals, urefvals = compute_WRFtime(timewrfv, refdate=refT,                \
+              tunitsval='minutes')
+            oaxisv = tvals
+     
+
+        # CF values of axis
+        CFdimvalues = gen.CFcorValues(axn)
+        print '  ' + axn + " ..."
+        print '    dimension in file:', dimn, ' variable in file:', vardimn
+        print '    CF-values _______'
+        print '    dimn:', CFdimvalues['dimn']
+        print '    variable-dimension. standard_name:', CFdimvalues['stdn'],         \
+          'long_name:', CFdimvalues['longname'], 'units', CFdimvalues['longname'],   \
+          'max, allowed rank:', CFdimvalues['maxrank']
+        print '    others:'
+        for ivn in CFdimvalues.keys():
+            if ivn != 'dimn' and ivn != 'stdn' and ivn != 'longname' and             \
+              ivn != 'units' and ivn != 'maxrank':
+                print ivn, ':', CFdimvalues[ivn]
+
+        # File values of the axis and assigning dimension length from rank in file
+        CFvardimvalues = []
+        if len(oaxisv.shape) == CFdimvalues['maxrank']:
+            axisv = oaxisv[:]
+            if CFdimvalues['length'] != -1:CFdimvalues['length'] = oaxisv.shape[0]
+            CFvardimvalues.append[dimn]
+        elif len(oaxisv.shape) == CFdimvalues['maxrank'] + 1:
+            print warnmsg
+            print '  '+ fname + ": dimension variable '" + vardimn + "' for axis '"+ \
+              axn + "' has a rank 1 unit higher that the maximum CF allowed:",       \
+              CFdimvalues['maxrank']
+            print '    assuming extra-dimension due to temporal values. ' +          \
+              'dimensions of variable:', oaxisv.dimensions
+            dimTn = axisinf['T'][0]
+            varslice = []
+            for dn in oaxisv.dimensions:
+                if dn == dimTn: varslice.append(0)
+                else: 
+                    CFvardimvalues.append(dn)
+                    varslice.append(slice(0,len(onc.dimensions[dn])))
+                    if dn == dimn and CFdimvalues['length'] != -1: 
+                        CFdimvalues['length'] = len(onc.dimensions[dn])
+            axsiv = oaxisv[tuple(varslice)]
+        else:
+            print errormsg
+            print '  ' + fname + ": axis-dimension '" + axn, ' of var: ', vardimn +  \
+              "' has maximum allowed rank of '", CFdimvalues['maxrank'],             \
+              ' while variable:', len(oaxisv.shape), ' and file-dimensions:',        \
+              oaxisv.dimensions, '!!'
+            print '    I do not know what to do'
+            onc.close()
+            quit(-1)
+        if CFdimvalues['length'] == -1: CFdimvalues['length'] = None
+        CFdimvals[axn] = CFdimvalues
+        filedimvals[axn] = axsiv[:]
+        CFaxisvardimvals[axn] = CFvardimvalues
+
+        # Some extra axis-stuff...
+        if axn == 'T':
+            dimt = len(axisv)
+            CFitime = axisv[0]
+            CFetime = axisv[dimt]
+            tunits = oaxisv.units
+            tcal = oaxisv.calendar
+            Sinit=datetimeStr_conversion(str(CFitime),'cftime,'+tunits,'Y-m-dTH:M:SZ')
+            Sendt=datetimeStr_conversion(str(CFetime),'cftime,'+tunits,'Y-m-dTH:M:SZ')
+            secfreq = CFtime_freq(axisv[0], axisv[1], tunits)
+            # Labelling frequency
+            if np.mod(secfreq, 100 * 365 * 24 * 3600.) == 0.:
+                Sfreq = str(int(secfreq/(100 * 365 * 24 * 3600.))) + 'C'
+            elif np.mod(secfreq, 365 * 24 * 3600.) == 0.: 
+                Sfreq = str(int(secfreq/(365 * 24 * 3600.))) + 'Y'
+            elif np.mod(secfreq, 7 * 24 * 3600.) == 0.: 
+                Sfreq = str(int(secfreq/(7 * 24 * 3600.))) + 'W'
+            elif np.mod(secfreq, 24 * 3600.) == 0.:
+                Sfreq = str(int(secfreq/(24.*3600.))) + 'D'
+            elif np.mod(secfreq,3600.) == 0.: 
+                Sfreq = str(int(secfreq/3600.)) + 'H'
+            elif np.mod(secfreq,60.) == 0.:
+                Sfreq = str(int(secfreq/60.)) + 'M'
+            else: Sfreq = 'UNKNOWN'
+
+    # Processing variables in file
+    ##
     for vn in varns:
         # getting CF information
@@ -22803 +22941 @@
         stdvarn = varinf[1]
         longvarn = varinf[4].replace('|', ' ')
-
+        utsvarn = varinf[5]
+
+        # getting variable from file
         if not onc.variables.has_key(vn):
             print errormsg
             print '  '+fname+ "': netcdf file '" + ncfile + "' does not contain " +  \
               "variable '" + vn + "' !!"
+            print '    available ones:', oncvars
             onc.close()
             quit(-1)
-
+        ovar = onc.variables[vn]
+        # Looking for dimensions of variable
+        varaxes = []
+        CFvardims = []
+        for axn in axes:
+            axisv = axisinf[axn]
+            dimn = axisv[0]
+            if gen.searchInfile(ovar.dimensions, dimn): 
+                varaxes.append[axn]
+                CFvardims.append(CFdimvals[axn]['dimn'])
+                
         filen = gattr['institude_id'] + '_' + cfvarn + '_' + gattr['CORDEX_domain'] +\
           '_' + gattr['driving_model_id'] + '_' + gattr['driving_experiment_name'] + \
           '_' + gattr['driving_model_ensemble_member'] + '_' + gattr['model_id'] +   \
-          '_' + gattr['rcm_version_id'] + '_' + Sfrequency + '_' StartTime-EndTime + \
+          '_' + gattr['rcm_version_id'] + '_' + Sfreq + '_' + Sinit + '-' + Sendt +  \
           '.nc'
+
+        onewnc = NetCDFFile(filen,'w')
+        # Definition only of that coordinates from the variable
+        for axn in varaxes:
+            dvals = CFdimvals[axn]
+            newdim = onewnc.createDimension(dvals['dimn'], dvals['legnth'])
+
+        # Definition only of that variables-coordinates from the axes in the variable
+        for axn in varaxes:
+            dvals = CFdimvals[axn]
+            newvar = onewnc.createVariable(dvals['dimn'], 'f8',                      \
+              tuple(CFaxesdimvarvals[axn]))
+            basicvardef(newvar, dvals['stdn'], dvals['longname'], dvals['units'])
+            for ivn in dvals.keys():
+                if ivn != 'dimn' and ivn != 'stdn' and ivn != 'longname' and         \
+                  ivn != 'units' and ivn != 'maxrank' and ivn != 'length':
+                    set_attribute(newvar,ivn,dvals[ivn])
+
+        # Variable
+        newvar=onewnc.createVariable(cfvarn, 'f4', tuple(CFvardims),                 \
+          fill_Value=gen.fillValueR)
+        basicvardef(newvar, stdvarn, longvarn, utsvarn)
+        set_attribute(newvar, 'coordinates', ' '.join(CFvardims))
+        # Setting values, but taking into account pre-existing masked values
+        varv = ovar[:]
+        if type(varv) == type(gen.mamat):
+            imaskv = varv.fill_value
+            varv = np.where(varv == fill_value, gen.fillValueR, varv)
+
+        newvar[:] = varv[:]
+        onewnc.sync()
+
+        # Global attributes
+        for attrn in gattrN:
+            set_attribute(onewnc, attrn, gattr[attrn])
+
+        onewnc.sync()
+        add_global_PyNCplot(onewnc, 'PyNCplot', fname, '1.0')
+        onewnc.sync()
+
+        print "Successfull output of file '" + filen + "' !!"
+        
 
     onc.close()