Changeset 1871 in lmdz_wrf

Timestamp:

Mar 27, 2018, 3:22:39 PM (7 years ago)

Author:

lfita

Message:

Working version of `CFmorzization' with extrafiles

File:

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

trunk/tools/nc_var_tools.py

@@ -21726 +21726 @@
 
                 newvar = newnc.createVariable(vn, nctype(ovar.dtype), tuple(dimvar))
-                print '    Lluis shapes newvar:', newvar.shape, 'ovar:', ovar.shape
                 newvar[:] = ovar[:]
 
@@ -22911 +22910 @@
     return od
 
+def getvar_listonc(varn, fns, listoc):
+    """ Function to get a variable object from a list of netCDF file object
+      varn= name of the variable to be found
+      fns= list with the name of the files
+      listoc= list of NetCDF file objects
+    """
+    fname = 'getvar_listonc'
+
+    allfs = {}
+    for ioc in range(len(listoc)):
+        oc = listoc[ioc]
+        # Seems that ObjNetcdf.filepath() requires netcdf >= 4.1.2....
+        #  Thus avoiding to use the functionality
+        #allfs[oc.filepath()] = ...
+        allfs[fns[ioc]] = [oc.dimensions, oc.variables.keys()]
+
+    # Looping for all files
+    found = False
+    for ioc in range(len(listoc)):
+        oc = listoc[ioc]
+        # Seems that ObjNetcdf.filepath() requires netcdf >= 4.1.2....
+        #  Thus avoiding to use the functionality
+        oncvars = allfs[fns[ioc]]
+        oncvars.sort()
+        if oc.variables.has_key(varn):
+            ov = oc.variables[varn]
+            found = True
+            break
+
+    if not found:
+        print errormsg
+        print '  ' + fname + "': any netcdf file '" + ','.join(fns) + "' contain " + \
+          "variable '" + varn + "' !!"
+        print '    available ones _______'
+        ioc = 0
+        for fn in fns:
+            print fn, ':', allfs[fn][1]
+            listoc[ioc].close()
+            ioc = ioc + 1
+        quit(-1)
+
+    return ov
+
 def projection2D_information(PROJn, AXinf, projinf, dx, dy, Fns, oncs):
     """ Function to get the basic information of the 2D-projection: new dimensions, variables, ...
@@ -22938 +22980 @@
       'fileYvals': 'variable with the reference values for the projection ' +        \
         'along Y axis (e.g.: values for rlat)'}
+    nonCFattr = list(projfileexpl.keys())
     CFprojweb = 'http://cfconventions.org/Data/cf-conventions/cf-conventions-1.7/' + \
       'build/ch05s06.html'
@@ -22957 +23000 @@
           'proj_units']
         # Removing the non-CF attributes
-        nonCFattr = ['x_resolution', 'y_resolution', 'proj_units']
         for atn in proj_attr:
             if not projinf.has_key(atn):
@@ -23008 +23050 @@
           "grid_north_pole_longitude", 'fileXvals', 'fileYvals']
         # Removing the non-CF attributes
-        nonCFattr = ['fileXvals', 'fileYvals']
         for atn in proj_attr:
             if not projinf.has_key(atn):
@@ -23053 +23094 @@
         oYvar = ovar[0]
         if len(oYvar.shape) == 3:
-            newvarv['rlat'] = oYvar[0,0,:]
+            newvarv['rlat'] = oYvar[0,:,0]
         elif len(oYvar.shape) == 2:
-            newvarv['rlat'] = oYvar[0,:]
+            newvarv['rlat'] = oYvar[:,0]
         elif len(oYvar.shape) == 1:
             newvarv['rlat'] = oYvar[:]
@@ -23073 +23114 @@
           "standard_parallel", 'x_resolution', 'y_resolution', 'proj_units']
         # Removing the non-CF attributes
-        nonCFattr = ['x_resolution', 'y_resolution', 'proj_units']
         for atn in proj_attr:
             if not projinf.has_key(atn):
@@ -23087 +23127 @@
                 quit(-1)
 
-        newdim = {}
-        newvar = {}
-        newvarv = {}
-
-def getvar_listonc(varn, fns, listoc):
-    """ Function to get a variable object from a list of netCDF file object
-      varn= name of the variable to be found
-      fns= list with the name of the files
-      listoc= list of NetCDF file objects
-    """
-    fname = 'getvar_listonc'
-
-    allfs = {}
-    for ioc in range(len(listoc)):
-        oc = listoc[ioc]
-        # Seems that ObjNetcdf.filepath() requires netcdf >= 4.1.2....
-        #  Thus avoiding to use the functionality
-        #allfs[oc.filepath()] = ...
-        allfs[fns[ioc]] = [oc.dimensions, oc.variables.keys()]
-
-    # Looping for all files
-    found = False
-    for ioc in range(len(listoc)):
-        oc = listoc[ioc]
-        # Seems that ObjNetcdf.filepath() requires netcdf >= 4.1.2....
-        #  Thus avoiding to use the functionality
-        oncvars = allfs[fns[ioc]]
-        oncvars.sort()
-        if oc.variables.has_key(varn):
-            ov = oc.variables[varn]
-            found = True
-            break
-
-    if not found:
-        print errormsg
-        print '  ' + fname + "': any netcdf file '" + ','.join(fns) + "' contain " + \
-          "variable '" + varn + "' !!"
-        print '    available ones _______'
-        ioc = 0
-        for fn in fns:
-            print fn, ':', allfs[fn][1]
-            listoc[ioc].close()
-            ioc = ioc + 1
-        quit(-1)
-
-    return ov
-
-def projection2D_information(PROJn, AXinf, projinf, dx, dy, Fns, oncs):
-    """ Function to get the basic information of the 2D-projection: new dimensions, variables, ...
-      Reproducing: class Projection from wrfncxnj.py
-        PROJn= name of the projection
-        AXinf= dictionary with the informatino of the axes 'X', 'Y', 'Z', 'T'
-        projinf= dictionary with all the information of the projection
-        Fns= list with the name of the files
-        oncs= list of netCDF objects from which retrieve the information
-      Returns:
-        newdim: dictionary new dimensions to be added {dimn1: dimlength1, ..., dimnN: dimlengthN} 
-        newvar: dictionary new variables to be added {varn1: [stdn1, longn1, units1], ..., varnM: [stdnM, longnM, unitsM]}
-        newvarv: dictionary values new variables to be added {varn1: [values,...], ..., varnM: [values,...]}
-        nonCFattr: list of attributes from PROJn non-CF and non to be added to the file
-    """
-    fname = 'projection2D_information'
-
-    availproj = ['Lambert_Conformal', 'Rotated_Pole/Rotate_Latituded_Longitude',     \
-      'Mercator']
-
-    # Explanation of certain attributes in projection file
-    projfileexpl = {'x_resolution': 'Resolution in file alon x axis',                \
-      'y_resolution': 'Resolution in file alon y axis',                              \
-      'proj_units': 'units of the projection (e.g.: km, m, degrees,...)',            \
-      'fileXvals': 'variable with the reference values for the projection ' +        \
-        'along X axis (e.g.: values for rlon)',                                      \
-      'fileYvals': 'variable with the reference values for the projection ' +        \
-        'along Y axis (e.g.: values for rlat)'}
-    CFprojweb = 'http://cfconventions.org/Data/cf-conventions/cf-conventions-1.7/' + \
-      'build/ch05s06.html'
-
-    isrotated = False
-    if PROJn.find('rotated') != -1 or PROJn.find('Rotated') != -1:
-        isrotated = True
-        projpiec = PROJn.split('_')
-        if PROJn.find('rotated') != -1:
-            projpiec.remove('rotated')
-        else:
-            projpiec.remove('Rotated')
-        projn = '_'.join(projpiec)
-
-    if PROJn == "Lambert_Conformal":
-        proj_attr = ["grid_mapping_name", "cone_type", "northern_parallel",          \
-          "southern_parallel", "longitude_of_central_meridian",                      \
-          "latitude_of_projection_origin", 'x_resolution', 'y_resolution',           \
-          'proj_units']
-        # Removing the non-CF attributes
-        nonCFattr = ['x_resolution', 'y_resolution', 'proj_units']
-        for atn in proj_attr:
-            if not projinf.has_key(atn):
-                print errormsg
-                print '  ' + fname + ": to process projection '" + PROJn +           \
-                  "' is required to provide argument '" + atn + "' in the " +        \
-                  "'[projectfile]' !!"
-                if projfileexpl.has_key(atn): 
-                    print '    ', atn, ':', projfileexpl[atn]
-                else:
-                    print '    ', atn, ': CF-attribute, please visit: ', CFprojweb
-                print '    values provided:', projinf
-                quit(-1)
         newdim = {}
         newvar = {}
@@ -23227 +23161 @@
         AXinf['Y'] = axv
 
-    elif PROJn == "Rotated_Pole" or PROJn == 'Rotated_Latitude_Longitude':
-        proj_attr = ["grid_mapping_name", "grid_north_pole_latitude",                \
-          "grid_north_pole_longitude", 'fileXvals', 'fileYvals']
-        # Removing the non-CF attributes
-        nonCFattr = ['fileXvals', 'fileYvals']
-        for atn in proj_attr:
-            if not projinf.has_key(atn):
-                print errormsg
-                print '  ' + fname + ": to process projection '" + PROJn +           \
-                  "' is required to provide argument '" + atn + "' in the " +        \
-                  "'[projectfile]' !!"
-                if projfileexpl.has_key(atn): 
-                    print '    ', atn, ':', projfileexpl[atn]
-                else:
-                    print '    ', atn, ': CF-attribute, please visit: ', CFprojweb
-                print '    values provided:', projinf
-                quit(-1)
-        newdim = {}
-        newvar = {}
-        newvarv = {}
-
-        newdim['rlon'] = dx
-        newvar['rlon'] = ['grid_longitude', 'longitude in rotated pole grid',        \
-          'degrees']
-
-        ovar = getvar_listonc(projinf['fileXvals'][0], Fns, oncs)
-        oXvar = ovar[0]
-        if len(oXvar.shape) == 3:
-            newvarv['rlon'] = oXvar[0,0,:]
-        elif len(oXvar.shape) == 2:
-            newvarv['rlon'] = oXvar[0,:]
-        elif len(oXvar.shape) == 1:
-            newvarv['rlon'] = oXvar[:]
-        else:
-            print errormsg
-            print '  ' + fname + ": rank ", oXvar.shape, " of values for X ref " +   \
-              "not ready !!"
-            quit(-1)
-        # updating values on axes
-        axv = AXinf['X']
-        axv['dimn'] = 'rlon'
-        axv['vdimn'] = 'lon'
-        AXinf['X'] = axv
-
-        newdim['rlat'] = dy
-        newvar['rlat'] = ['grid_latitude', 'latitude in rotated pole grid', 'degrees']
-        ovar = getvar_listonc(projinf['fileYvals'][0], Fns, oncs)
-        oYvar = ovar[0]
-        if len(oYvar.shape) == 3:
-            newvarv['rlat'] = oYvar[0,0,:]
-        elif len(oYvar.shape) == 2:
-            newvarv['rlat'] = oYvar[0,:]
-        elif len(oYvar.shape) == 1:
-            newvarv['rlat'] = oYvar[:]
-        else:
-            print errormsg
-            print '  ' + fname + ": rank ", oYvar.shape, " of values for Y ref " +   \
-              "not ready !!"
-            quit(-1)
-        # updating values on axes
-        axv = AXinf['Y']
-        axv['dimn'] = 'rlat'
-        axv['vdimn'] = 'lat'
-        AXinf['Y'] = axv
-
-    elif PROJn == "Mercator":
-        proj_attr = ["grid_mapping_name", "longitude_of_projection_origin",          \
-          "standard_parallel", 'x_resolution', 'y_resolution', 'proj_units']
-        # Removing the non-CF attributes
-        nonCFattr = ['x_resolution', 'y_resolution', 'proj_units']
-        for atn in proj_attr:
-            if not projinf.has_key(atn):
-                print errormsg
-                print '  ' + fname + ": to process projection '" + PROJn +           \
-                  "' is required to provide argument '" + atn + "' in the " +        \
-                  "'[projectfile]' !!"
-                if projfileexpl.has_key(atn): 
-                    print '    ', atn, ':', projfileexpl[atn]
-                else:
-                    print '    ', atn, ': CF-attribute, please visit: ', CFprojweb
-                print '    values provided:', projinf
-                quit(-1)
-
-        newdim = {}
-        newvar = {}
-        newvarv = {}
-
-        newdim['x'] = dx
-        newvar['x'] = ['projection_x_coordinate', 'x coordinate of projection', 'km']
-        Xres = np.float(projinf['x_resolution'][0])
-        if projinf['proj_units'][0] == 'km':
-            scale = 1.
-        elif projinf['proj_units'][0] == 'm':
-            scale = 1000.
-        else:
-            print errormsg
-            print '  ' + fname + ": Resolution units '" + projinf['proj_units'][0] + \
-              "' not ready !!"
-            print '    available ones:', ['km', 'm']
-            quit(-1)
-        newvarv['x'] = (np.arange(1,dx+1)-dx/2)*Xres/scale
-        # updating values on axes
-        axv = AXinf['X']
-        axv['dimn'] = 'x'
-        axv['vdimn'] = 'lon'
-        AXinf['X'] = axv
-
-        newdim['y'] = dy
-        newvar['y'] = ['projection_y_coordinate', 'y coordinate of projection', 'km']
-        Yres = np.float(projinf['y_resolution'][0])
-        newvarv['y'] = (np.arange(1,dy+1)-dy/2)*Yres/scale
-        # updating values on axes
-        axv = AXinf['Y']
-        axv['dimn'] = 'y'
-        axv['vdimn'] = 'lat'
-        AXinf['Y'] = axv
-
     else:
         print errormsg
@@ -23362 +23179 @@
     """
     fname = 'CFvars'
-    inCFvarfile = 'CFvariables.dat'
-    availactions = ['new_z_dim']
+
+    availactions = ['new_p_dim','new_z_dim']
+
+    folder = os.path.dirname(os.path.realpath(__file__))
+    inCFvarfile = folder + '/' + 'CFvariables.dat'
 
     oCFf = open(inCFvarfile, 'r')
@@ -23670 +23490 @@
     allfileNs = {}
     for ionc in range(len(allonc)):
-        onc = allonc[ionc]
-        allfileNs[afNs[ionc]] = [onc.dimensions.keys(), onc.variables.keys()]
+        oc = allonc[ionc]
+        allfileNs[afNs[ionc]] = [oc.dimensions.keys(), oc.variables.keys()]
 
     # Creation of new file name
@@ -23837 +23657 @@
         addCFdimvals = {}
         for newdn in newDim.keys():
+            print infmsg
+            print '      ' + fname + ": adding dimension '" + newdn + "' ..."
             addCFdimvals[newdn] = {'dimn': newdn, 'vdimn': newdn,                    \
               'stdn': newVar[newdn][0], 'longname': newVar[newdn][1],                \
               'units': newVar[newdn][2], 'maxrank': 1,                               \
-              'length': filedimvals['X'].shape[1]}
+              'length': newDim[newdn]}
             filedimvals[newdn] = newVARv[newdn]
+            gen.printing_dictionary(addCFdimvals[newdn])
 
     print 'Values for variables-axes from file ________'
@@ -23895 +23718 @@
         utsvarn = varinf[5]
 
+        oncvars = onc.variables.keys()
+        oncvars.sort()
+
         # getting variable from file
         if not onc.variables.has_key(vn):
@@ -23952 +23778 @@
                 if not gen.searchInlist(onewnc.variables.keys(), dvals['dimn']):
                     ddvals = addCFdimvals[dvals['dimn']]
-                    newvar = onewnc.createVariable(dvals['dimn'], 'f8',              \
-                      tuple(dvals['dimn']))
+                    newvar = onewnc.createVariable(ddvals['vdimn'], 'f8',            \
+                      tuple([ddvals['dimn']]))
                     basicvardef(newvar, ddvals['stdn'], ddvals['longname'],          \
                       ddvals['units'])
-                    newvar[:] = filedimvals[dvals['dimn']]
+                    newvar[:] = filedimvals[ddvals['dimn']]
                 if not gen.searchInlist(onewnc.variables.keys(), dvals['vdimn']):
                     newvar = onewnc.createVariable(dvals['vdimn'], 'f4',             \
@@ -23967 +23793 @@
                           ivn != 'units' and ivn != 'maxrank' and ivn != 'length':
                             set_attribute(newvar,ivn,dvals[ivn])
+            onewnc.sync()
 
         # Variable