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

Timestamp:

Jan 28, 2016, 2:48:41 PM (9 years ago)

Author:

lfita

Message:

Fixing right computation of `Npoints' in 'lonlatProj'

File:

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

trunk/tools/nc_var_tools.py

@@ -17719 +17719 @@
      dyres= resolution on the SN axis [m]
      project= projection
-       'lonlat': regular lon/lat projection
+       'lonlat': standard lon/lat projection tacking number of points from [dx/dy]res at the Equator
+       'lonlat_dxdyFix': projection with a fixed grid distance
      wfile= Should an output file be written (True/False)
     """
@@ -17730 +17731 @@
     RadDeg = 180./np.pi
 
-    Npts = 0
+    Npts = 1
     if project == 'lonlat':
+        Nlon = int(2. * np.pi * EarthR / dxres) + 1
+        Nlat = int(np.pi * EarthR / dxres) + 1
+
+        lonProj = np.zeros((Nlon), dtype=np.float)
+        latProj = np.zeros((Nlat), dtype=np.float)
+
+        lonProj = -180. + np.arange(Nlon) * 2. * np.pi * RadDeg / (Nlon - 1)
+        latProj = -90. + np.arange(Nlat) * np.pi * RadDeg / (Nlat - 1)
+
+        dimx = Nlon
+        dimy = Nlat
+        Npts = dimx*dimy
+        Npoints = np.arange((Npts)).reshape(dimy,dimx)
+
+    elif project == 'lonlat_dxdyFix':
         lonProj = []
         latProj = []
@@ -17746 +17762 @@
             while lon <= 2.*np.pi:
 #                print lonarc, latarc, lon, lat, [lon*RadDeg, lat*RadDeg]
-                lonProj.append(lon*RadDeg) 
+                lonProj.append(-180. + lon*RadDeg) 
                 latProj.append(lat*RadDeg) 
                 lon = lon + lonarc
@@ -17752 +17768 @@
         lonProj.append(0.)
         latProj.append(90.)
-
+        Npts = Npts + 1
+
+        dimx = Npts
+        dimy = Npts
+        Npoints = np.arange(Npts)
     else:
         print errormsg
@@ -17764 +17784 @@
 
 # Dimensions
-        Npoints = len(lonProj)
-        newdim = onewnc.createDimension('Npts', Npoints)
+        newdim = onewnc.createDimension('x', dimx)
+        newdim = onewnc.createDimension('y', dimy)
+        newdim = onewnc.createDimension('Npts', Npts)
         newdim = onewnc.createDimension('Ncoord', 2)
 
 # Variables
-        newvar = onewnc.createVariable('lon', 'f8', ('Npts'))
+        newvar = onewnc.createVariable('lon', 'f8', ('x'))
         basicvardef(newvar, 'longitude', 'Longitude', 'degress_East')
-        newvar[:] = lonProj[:]
+        newvar[:] = lonProj
         newvar.setncattr('axis', 'X')
         newvar.setncattr('_CoordinateAxisType', 'Lon')   
 
-        newvar = onewnc.createVariable('lat', 'f8', ('Npts'))
+        newvar = onewnc.createVariable('lat', 'f8', ('y'))
         basicvardef(newvar, 'latitude', 'Latitude', 'degress_North')
-        newvar[:] = latProj[:]
+        newvar[:] = latProj
         newvar.setncattr('axis', 'Y')
         newvar.setncattr('_CoordinateAxisType', 'Lat')   
 
-        newvar = onewnc.createVariable('points', 'i8', ('Npts'))
+        if len(Npoints.shape) == 2:
+            newvar = onewnc.createVariable('points', 'i8', ('y', 'x'))
+        else:
+            newvar = onewnc.createVariable('points', 'i8', ('Npts'))
+
         basicvardef(newvar, 'points', '# point of the projection', '-')
-        newvar[:] = np.int64(np.arange(Npoints))
+        newvar[:] = Npoints
         newvar.setncattr('coordinates',  'lon lat')
+        print Npoints.shape
     
 
@@ -17801 +17827 @@
         print fname + "Successfull written of '" + ofile + "' !!'"
 
-    return lonlatProj
+    return lonProj, latProj
 
 #lonlatProj(100000.,100000.,'lonlat',True)