Changeset 184 for trunk/MESOSCALE/PLOT/PYTHON/mylib

Timestamp:

Jul 1, 2011, 5:14:38 AM (14 years ago)

Author:

aslmd

Message:

MESOSCALE: much improved python scripts for winds. lots of automatic settings.

File:

• trunk/MESOSCALE/PLOT/PYTHON/mylib/myplot.py (modified) (5 diffs)

trunk/MESOSCALE/PLOT/PYTHON/mylib/myplot.py

@@ -35 +35 @@
         return wlon,wlat
 
+def getproj (nc):
+    map_proj = getattr(nc, 'MAP_PROJ')
+    cen_lat  = getattr(nc, 'CEN_LAT')
+    if map_proj == 2:
+        if cen_lat > 10.:    
+            proj="npstere"
+            print "NP stereographic polar domain" 
+        else:            
+            proj="spstere"
+            print "SP stereographic polar domain"
+    elif map_proj == 1: 
+        print "lambert projection domain" 
+        proj="lcc"
+    elif map_proj == 3: 
+        print "mercator projection"
+        proj="merc"
+    return proj    
+
 def ptitle (name):
     from matplotlib.pyplot import title
@@ -43 +61 @@
     import numpy as np
     return [[np.min(lon2d),np.max(lon2d)],[np.min(lat2d),np.max(lat2d)]]
+
+def wrfinterv (lon2d,lat2d):
+    nx = len(lon2d[0,:])-1
+    ny = len(lon2d[:,0])-1
+    return [[lon2d[0,0],lon2d[nx,ny]],[lat2d[0,0],lat2d[nx,ny]]]
 
 def makeplotpngres (filename,res,pad_inches_value=0.25,folder='',disp=True):
@@ -58 +81 @@
     return
 
-def getcoord2d (nc,nlat='XLAT',nlon='XLONG'):
-        import numpy as np
-        lat = nc.variables[nlat][0,:,:]
-        lon = nc.variables[nlon][0,:,:]
-        if np.array(lat).ndim != 2:     [lon2d,lat2d] = np.meshgrid(lon,lat)
-        else:                           [lon2d,lat2d] = [lon,lat]
-        return lon2d,lat2d
+def dumpbdy (field):
+    nx = len(field[0,:])-1
+    ny = len(field[:,0])-1
+    return field[5:ny-5,5:nx-5]
+
+def getcoord2d (nc,nlat='XLAT',nlon='XLONG',is1d=False):
+    import numpy as np
+    if is1d:
+        lat = nc.variables[nlat][:]
+        lon = nc.variables[nlon][:]
+        [lon2d,lat2d] = np.meshgrid(lon,lat)
+    else:
+        lat = nc.variables[nlat][0,:,:]
+        lon = nc.variables[nlon][0,:,:]
+        [lon2d,lat2d] = [lon,lat]
+    return lon2d,lat2d
 
 def smooth (field, coeff):
@@ -94 +126 @@
         return improc
 
-def vectorfield (u, v, x, y, stride=3, scale=15., factor=250., color='black', csmooth=1):
-        ## scale regle la reference du vecteur
-        ## factor regle toutes les longueurs (dont la reference). l'AUGMENTER pour raccourcir les vecteurs.
-        import  matplotlib.pyplot               as plt
-        import  numpy                           as np
-        posx = np.max(x)*0.90 
-        posy = np.mean(y)
-        u = smooth(u,csmooth)
-        v = smooth(v,csmooth)
-        q = plt.quiver( x[::stride,::stride],\
-                        y[::stride,::stride],\
-                        u[::stride,::stride],\
-                        v[::stride,::stride],\
-                        angles='xy',color=color,\
-                        scale=factor,width=0.003 )
-        if color=='white':      kcolor='black'
-        elif color=='yellow':   kcolor=color
-        else:                   kcolor=color
-        p = plt.quiverkey(q,posx,posy,scale,\
-                        str(int(scale)),coordinates='data',color=kcolor)
-        return p
+def vectorfield (u, v, x, y, stride=3, scale=15., factor=250., color='black', csmooth=1, key=True):
+    ## scale regle la reference du vecteur
+    ## factor regle toutes les longueurs (dont la reference). l'AUGMENTER pour raccourcir les vecteurs.
+    import  matplotlib.pyplot               as plt
+    import  numpy                           as np
+    posx = np.max(x) + np.std(x) / 3.  ## pb pour les domaines globaux ...
+    posy = np.mean(y)
+    #posx = np.min(x)
+    #posy = np.max(x)
+    u = smooth(u,csmooth)
+    v = smooth(v,csmooth)
+    widthvec = 0.005 #0.003
+    q = plt.quiver( x[::stride,::stride],\
+                    y[::stride,::stride],\
+                    u[::stride,::stride],\
+                    v[::stride,::stride],\
+                    angles='xy',color=color,\
+                    scale=factor,width=widthvec )
+    if color=='white':     kcolor='black'
+    elif color=='yellow':  kcolor=color
+    else:                  kcolor=color
+    if key: p = plt.quiverkey(q,posx,posy,scale,\
+                   str(int(scale)),coordinates='data',color=kcolor)
+    return 
 
 def display (name):
-        from    os                              import system
-        system("display "+name+" > /dev/null 2> /dev/null &")
-        return name
+    from os import system
+    system("display "+name+" > /dev/null 2> /dev/null &")
+    return name
 
 def findstep (wlon):
-    steplon = int((wlon[1]-wlon[0])/3.)
-    step = 60.
-    if steplon < 60.: step = 30.
-    if steplon < 30.: step = 15.
-    if steplon < 15.: step = 10.
-    if steplon < 10.: step = 5.
-    if steplon < 5.:  step = 1.
+    steplon = int((wlon[1]-wlon[0])/4.)  #3
+    step = 120.
+    while step > steplon and step > 15. :       step = step / 2.
+    if step <= 15.:
+        while step > steplon and step > 5.  :   step = step - 5.
+    if step <= 5.:
+        while step > steplon and step > 1.  :   step = step - 1.
+    if step <= 1.:
+        step = 1. 
     return step
 
@@ -137 +174 @@
     meanlon = 0.5*(wlon[0]+wlon[1])
     meanlat = 0.5*(wlat[0]+wlat[1])
+    if   wlat[0] >= 80.:   blat =  40. 
+    elif wlat[0] <= -80.:  blat = -40. 
+    else:                  blat = wlat[0]
     h = 2000.
-    if   char == "cyl":     m = Basemap(projection='cyl',llcrnrlat=wlat[0],urcrnrlat=wlat[1],llcrnrlon=wlon[0],urcrnrlon=wlon[1])
-    elif char == "moll":    m = Basemap(projection='moll',lon_0=meanlon)
-    elif char == "ortho":   m = Basemap(projection='ortho',lon_0=meanlon,lat_0=meanlat)
-    elif char == "lcc":     m = Basemap(projection='lcc',lat_1=meanlat,lat_0=meanlat,lon_0=meanlon,\
+    radius = 3397200
+    if   char == "cyl":     m = Basemap(rsphere=radius,projection='cyl',\
                               llcrnrlat=wlat[0],urcrnrlat=wlat[1],llcrnrlon=wlon[0],urcrnrlon=wlon[1])
-    elif char == "npstere": m = Basemap(projection='npstere', boundinglat=wlat[0], lon_0=0.)
-    elif char == "spstere": m = Basemap(projection='spstere', boundinglat=wlat[0], lon_0=0.)
-    elif char == "nsper":   m = Basemap(projection='nsper',lon_0=meanlon,lat_0=meanlat,satellite_height=h*1000.)
-    fontsizemer = int(mpl.rcParams['font.size']*2./3.)
-    if char in ["cyl","lcc"]:   step = findstep(wlon)
-    else:                       step = 10.
+    elif char == "moll":    m = Basemap(rsphere=radius,projection='moll',lon_0=meanlon)
+    elif char == "ortho":   m = Basemap(rsphere=radius,projection='ortho',lon_0=meanlon,lat_0=meanlat)
+    elif char == "lcc":     m = Basemap(rsphere=radius,projection='lcc',lat_1=meanlat,lat_0=meanlat,lon_0=meanlon,\
+                              llcrnrlat=wlat[0],urcrnrlat=wlat[1],llcrnrlon=wlon[0],urcrnrlon=wlon[1])
+    elif char == "npstere": m = Basemap(rsphere=radius,projection='npstere', boundinglat=blat, lon_0=0.)
+    elif char == "spstere": m = Basemap(rsphere=radius,projection='spstere', boundinglat=blat, lon_0=0.)
+    elif char == "nsper":   m = Basemap(rsphere=radius,projection='nsper',lon_0=meanlon,lat_0=meanlat,satellite_height=h*1000.)
+    elif char == "merc":    m = Basemap(rsphere=radius,projection='merc',lat_ts=0.,\
+                              llcrnrlat=wlat[0],urcrnrlat=wlat[1],llcrnrlon=wlon[0],urcrnrlon=wlon[1])
+    fontsizemer = int(mpl.rcParams['font.size']*3./4.)
+    if char in ["cyl","lcc","merc"]:   step = findstep(wlon)
+    else:                              step = 10.
     m.drawmeridians(np.r_[-180.:180.:step*2.], labels=[0,0,0,1], color='grey', fontsize=fontsizemer)
     m.drawparallels(np.r_[-90.:90.:step], labels=[1,0,0,0], color='grey', fontsize=fontsizemer)