Changeset 405 for trunk/UTIL/PYTHON/myplot.py

Timestamp:

Nov 21, 2011, 6:27:42 PM (13 years ago)

Author:

aslmd

Message:

GRAPHICS: 1. handling of anomaly in profiles as well. 2. fixed xmin ymin etc... that were not working with 1D plot. 3. fixed localtime for mesoscale files, still a bit of work though.

File:

• trunk/UTIL/PYTHON/myplot.py (modified) (10 diffs)

trunk/UTIL/PYTHON/myplot.py

@@ -67 +67 @@
 
 ## Author: AS + TN + AC
-def reducefield (input,d4=None,d3=None,d2=None,d1=None,yint=False,alt=None):
+def reducefield (input,d4=None,d3=None,d2=None,d1=None,yint=False,alt=None,anomaly=False):
     ### we do it the reverse way to be compliant with netcdf "t z y x" or "t y x" or "y x"
     ### it would be actually better to name d4 d3 d2 d1 as t z y x
+    ### ... note, anomaly is only computed over d1 and d2 for the moment
     import numpy as np
     from mymath import max,mean
+    csmooth = 12 ## a fair amount of grid points
     dimension = np.array(input).ndim
     shape = np.array(input).shape
     #print 'd1,d2,d3,d4: ',d1,d2,d3,d4
     print 'IN  REDUCEFIELD dim,shape: ',dimension,shape
+    if anomaly: print 'ANOMALY ANOMALY'
     output = input
     error = False
@@ -110 +113 @@
              output = mean(input[d4,:,:,:],axis=0)
              output = reduce_zaxis(output[d3,:,:],ax=0,yint=yint,vert=alt,indice=d3)
+             if anomaly: output = 100. * ((output / smooth(output,csmooth)) - 1.) 
              output = mean(output[d2,:],axis=0)
              output = mean(output[d1],axis=0)
@@ -115 +119 @@
              output = mean(input[d4,:,:,:],axis=0)
              output = reduce_zaxis(output[d3,:,:],ax=0,yint=yint,vert=alt,indice=d3)
+             if anomaly: output = 100. * ((output / smooth(output,csmooth)) - 1.) 
              output = mean(output[d2,:],axis=0)
         elif d4 is not None and d3 is not None and d1 is not None: 
              output = mean(input[d4,:,:,:],axis=0)
              output = reduce_zaxis(output[d3,:,:],ax=0,yint=yint,vert=alt,indice=d3)
+             if anomaly: output = 100. * ((output / smooth(output,csmooth)) - 1.) 
              output = mean(output[:,d1],axis=1)
         elif d4 is not None and d2 is not None and d1 is not None: 
              output = mean(input[d4,:,:,:],axis=0)
+             if anomaly:
+                 for k in range(output.shape[0]):  output[k,:,:] = 100. * ((output[k,:,:] / smooth(output[k,:,:],csmooth)) - 1.)
              output = mean(output[:,d2,:],axis=1)
              output = mean(output[:,d1],axis=1)
+             #noperturb = smooth1d(output,window_len=7)
+             #lenlen = len(output) ; output = output[1:lenlen-7] ; yeye = noperturb[4:lenlen-4]
+             #plot(output) ; plot(yeye) ; show() ; plot(output-yeye) ; show()
         elif d3 is not None and d2 is not None and d1 is not None: 
              output = reduce_zaxis(input[:,d3,:,:],ax=1,yint=yint,vert=alt,indice=d3) 
+             if anomaly:
+                 for k in range(output.shape[0]):  output[k,:,:] = 100. * ((output[k,:,:] / smooth(output[k,:,:],csmooth)) - 1.)
              output = mean(output[:,d2,:],axis=1)
              output = mean(output[:,d1],axis=1) 
@@ -131 +144 @@
              output = mean(input[d4,:,:,:],axis=0) 
              output = reduce_zaxis(output[d3,:,:],ax=0,yint=yint,vert=alt,indice=d3)
+             if anomaly: output = 100. * ((output / smooth(output,csmooth)) - 1.) 
         elif d4 is not None and d2 is not None:  
              output = mean(input[d4,:,:,:],axis=0) 
@@ -248 +262 @@
     from numpy import array
     from string import rstrip
-    namefile = rstrip( rstrip( namefile, chars="_z"), chars="_zabg")
+    namefile = rstrip( rstrip( rstrip( namefile, chars="_z"), chars="_zabg"), chars="_p")
     #### we assume that wrfout is next to wrfout_z and wrfout_zabg
     nc  = Dataset(namefile)
@@ -394 +408 @@
         [lon2d,lat2d] = [lon,lat]
     return lon2d,lat2d
+
+## FROM COOKBOOK http://www.scipy.org/Cookbook/SignalSmooth
+def smooth1d(x,window_len=11,window='hanning'):
+    import numpy
+    """smooth the data using a window with requested size.
+    This method is based on the convolution of a scaled window with the signal.
+    The signal is prepared by introducing reflected copies of the signal 
+    (with the window size) in both ends so that transient parts are minimized
+    in the begining and end part of the output signal.
+    input:
+        x: the input signal 
+        window_len: the dimension of the smoothing window; should be an odd integer
+        window: the type of window from 'flat', 'hanning', 'hamming', 'bartlett', 'blackman'
+            flat window will produce a moving average smoothing.
+    output:
+        the smoothed signal
+    example:
+    t=linspace(-2,2,0.1)
+    x=sin(t)+randn(len(t))*0.1
+    y=smooth(x)
+    see also: 
+    numpy.hanning, numpy.hamming, numpy.bartlett, numpy.blackman, numpy.convolve
+    scipy.signal.lfilter
+    TODO: the window parameter could be the window itself if an array instead of a string   
+    """
+    x = numpy.array(x)
+    if x.ndim != 1:
+        raise ValueError, "smooth only accepts 1 dimension arrays."
+    if x.size < window_len:
+        raise ValueError, "Input vector needs to be bigger than window size."
+    if window_len<3:
+        return x
+    if not window in ['flat', 'hanning', 'hamming', 'bartlett', 'blackman']:
+        raise ValueError, "Window is on of 'flat', 'hanning', 'hamming', 'bartlett', 'blackman'"
+    s=numpy.r_[x[window_len-1:0:-1],x,x[-1:-window_len:-1]]
+    #print(len(s))
+    if window == 'flat': #moving average
+        w=numpy.ones(window_len,'d')
+    else:
+        w=eval('numpy.'+window+'(window_len)')
+    y=numpy.convolve(w/w.sum(),s,mode='valid')
+    return y
 
 ## Author: AS
@@ -607 +663 @@
 def fmtvar (whichvar="def"):
     fmtvar    =     { \
-             "tk":           "%.0f",\
+             "TK":           "%.0f",\
              "T_NADIR_DAY":  "%.0f",\
              "T_NADIR_NIT":  "%.0f",\
              "TEMP_DAY":     "%.0f",\
              "TEMP_NIGHT":   "%.0f",\
-             "tpot":         "%.0f",\
+             "TPOT":         "%.0f",\
              "TSURF":        "%.0f",\
              "def":          "%.1e",\
@@ -623 +679 @@
              "VMR_ICE":      "%.1e",\
              "MTOT":         "%.1f",\
-             "anomaly":      "%.1f",\
+             "ANOMALY":      "%.1f",\
              "W":            "%.1f",\
              "WMAX_TH":      "%.1f",\
              "QSURFICE":     "%.0f",\
-             "Um":           "%.0f",\
+             "UM":           "%.0f",\
              "ALBBARE":      "%.2f",\
              "TAU":          "%.1f",\
@@ -650 +706 @@
              "def":          "spectral",\
              "HGT":          "spectral",\
-             "tk":           "gist_heat",\
+             "TK":           "gist_heat",\
              "TSURF":        "RdBu_r",\
              "QH2O":         "PuBu",\
@@ -666 +722 @@
              "W":            "jet",\
              "WMAX_TH":      "spectral",\
-             "anomaly":      "RdBu_r",\
+             "ANOMALY":      "RdBu_r",\
              "QSURFICE":     "hot_r",\
              "ALBBARE":      "spectral",\