Changeset 942

Timestamp:

Apr 25, 2013, 10:51:36 PM (12 years ago)

Author:

aslmd

Message:

UTIL PYTHON planetoplot_v2. handling time for mesoscale files. added -E
to pp_reload to loop on several colormaps to try. corrected a bug for 1D

plots with linestyle and colors. added a functionality to easily add te

xt (or crosses, dots, etc...) on a map through a text file add_text.txt

Location:

trunk/UTIL/PYTHON/planetoplot_v2

Files:

• pp.py (modified) (1 diff)
• pp_reload.py (modified) (2 diffs)
• ppclass.py (modified) (3 diffs)
• ppcompute.py (modified) (1 diff)
• ppplot.py (modified) (6 diffs)
• set_area.txt (modified) (1 diff)
• set_var.txt (modified) (1 diff)

trunk/UTIL/PYTHON/planetoplot_v2/pp.py

@@ -53 +53 @@
 parser.add_option('-o','--output',action='store',dest='filename',type="string",default="myplot",help="name of output files")
 parser.add_option('-d','--directory',action='store',dest='folder',type="string",default="./",help="directory of output files")
-parser.add_option('-s','--changetime',action='store',dest='changetime',type="string",default=None,help="transformation on time axis : [None] | mars_sol2ls")
+parser.add_option('-s','--changetime',action='store',dest='changetime',type="string",default=None,\
+                  help="transformation on time axis : [None] | mars_sol2ls | mars_meso_ls | mars_meso_sol | mars_meso_utc | mars_meso_lt ")
 # plot --> upper case
 # -- generic

trunk/UTIL/PYTHON/planetoplot_v2/pp_reload.py

@@ -14 +14 @@
 parser.add_option('-P','--proj',action='store',dest='proj',type="string",default=None,help='[2D] Define a new map projection')
 parser.add_option('-C','--colorb',action='store',dest='colorb',type="string",default=None,help="[2D] Define a new colormap")
+parser.add_option('-E','--explore',action='store_true',dest='explore',default=False,help="[2D] Try many colormaps")
 (opt,args) = parser.parse_args()
 # -----------------------------------------------------------------------
+clb = ["Greys","Blues","YlOrRd",\
+       "jet","spectral","hot",\
+       "RdBu","RdYlBu","Paired",\
+       "gist_ncar","gist_rainbow","gist_stern"]
+# -----------------------------------------------------------------------
+
 for files in args:
+
     yeah = pp()
     yeah.defineplot(loadfile=files)
     yeah.out = opt.out
+
     if opt.proj is not None:
       for plot in yeah.p:
@@ -32 +41 @@
       for plot in yeah.p:
         plot.title = opt.title
-    yeah.makeplot()
 
+    if opt.explore:
+      for cm in clb:
+       for plot in yeah.p:
+         plot.colorb = cm
+         plot.title = cm
+       yeah.makeplot()
+    else:         
+      yeah.makeplot()
+

trunk/UTIL/PYTHON/planetoplot_v2/ppclass.py

@@ -538 +538 @@
               obj.file = self.file[i]
               obj.var = self.var[j]
+              # get methods
+              obj.method_x = mx ; obj.method_y = my
+              obj.method_z = mz ; obj.method_t = mt
               # indicate the computation method
               obj.compute = self.compute
@@ -547 +550 @@
               obj.changetime = self.changetime
               obj.performtimechange()
-              ### get methods 
-              obj.method_x = mx ; obj.method_y = my
-              obj.method_z = mz ; obj.method_t = mt           
               ### get index
               obj.getindextime(dalist=st,ind=t,stride=self.stridet)
@@ -1206 +1206 @@
     # --------------------------
     def performtimechange(self):
-        if self.changetime is not None \
-          and self.name_t != "t grid points":
+        if self.changetime is not None:
             if self.verbose: print "**** OK. Converting time axis:",self.changetime
+            ### option added by T. Navarro 
             if self.changetime == "mars_sol2ls": 
                 self.field_t = ppcompute.mars_sol2ls(self.field_t)
+            ### options added by A. Spiga
+            elif "mars_meso" in self.changetime:
+                if 'Times' not in self.f.variables.keys():
+                    if self.verbose: print "!! WARNING !! Variable Times not in file. Cannot proceed to change of time axis."
+                else:
+                    # get the array of strings describing dates
+                    dates = self.f.variables['Times']
+                    dates.set_auto_maskandscale(False) # necessary to solve the api Times bug!
+                    # get ls sol utc from those strings
+                    ls, sol, utc = ppcompute.mars_date(dates[:])
+                    # populate self.field_t with the right output from mars_date
+                    if self.changetime == "mars_meso_ls": 
+                        self.field_t = ls 
+                        self.name_t = "Ls"
+                    elif self.changetime == "mars_meso_sol": 
+                        self.field_t = sol 
+                        self.name_t = "sol"
+                    elif self.changetime == "mars_meso_utc" \
+                        and ( self.changetime == "mars_meso_lt" \
+                              and not hasattr(self.f,'CEN_LON') ): 
+                        self.field_t = ppcompute.timecorrect(utc)
+                        self.name_t = "utc"
+                        if self.method_t == "fixed": 
+                            self.field_t = self.field_t % 24 # so that the user is not mistaken!
+                    elif self.changetime == "mars_meso_lt":
+                        self.field_t = ppcompute.timecorrect(utc) + getattr(self.f,'CEN_LON') / 15.
+                        self.field_t = ppcompute.timecorrect(self.field_t)
+                        self.name_t = "local time (center of domain)"
+                        if self.method_t == "fixed": 
+                            self.field_t = self.field_t % 24 # so that the user is not mistaken!
             else:
                 print "!! WARNING !! This time change is not implemented. Nothing is done."
+            if self.verbose: print "**** OK. new t axis values [%5.1f,%5.1f]" % (self.field_t.min(),self.field_t.max())
 
     # get list of index to be retrieved for time axis

trunk/UTIL/PYTHON/planetoplot_v2/ppcompute.py

@@ -231 +231 @@
       return solout
 
+# mars_date
+# author A. Spiga
+# ------------------------------
+# get Ls, sol, utc from a string with format yyyy-mm-dd_hh:00:00
+# -- argument timechar is a vector of such strings indicating dates
+# -- example: timechar = nc.variables['Times'][:] in mesoscale files
+# NB: uses mars_sol2ls function above
+# ------------------------------
+def mars_date(timechar):
+    # some preliminary information
+    days_in_month = [61, 66, 66, 65, 60, 54, 50, 46, 47, 47, 51, 56]
+    plus_in_month = [ 0, 61,127,193,258,318,372,422,468,515,562,613]
+    # get utc and sol from strings
+    utc = [] ; sol = []
+    for zetime in timechar:
+        dautc = float(zetime[11]+zetime[12]) + float(zetime[14]+zetime[15])/37.
+        dasol = dautc / 24.
+        dasol = dasol + plus_in_month[ int(zetime[5]+zetime[6])-1 ] + int(zetime[8]+zetime[9])
+        dasol = dasol - 1 ##les sols GCM commencent a 0
+        utc.append(dautc)
+        sol.append(dasol)
+    sol = np.array(sol)
+    utc = np.array(utc)
+    # get ls from sol
+    ls = mars_sol2ls(sol)
+    return ls, sol, utc
+
+# timecorrect
+# author A. Spiga
+# -----------------------------
+# ensure time axis is monotonic
+# correct negative values
+# -----------------------------
+def timecorrect(time):
+    for ind in range(len(time)-1):
+        if time[ind] < 0.: time[ind] = time[ind] + 24.
+        if time[ind+1] < time[ind]: time[ind+1] = time[ind+1] + 24.
+    return time
+

trunk/UTIL/PYTHON/planetoplot_v2/ppplot.py

@@ -220 +220 @@
     amax = ppcompute.max(field)
     if np.abs(amin) < 1.e-15: 
-        zevmin = 0.
         cmin = 0.
     else:
@@ -419 +418 @@
         if self.lstyle is not None and self.color is not None:
             mpl.plot(x,y,self.color+self.lstyle,marker=self.marker,label=self.label)
+        elif self.color is not None:
+            mpl.plot(x,y,color=self.color,marker=self.marker,label=self.label)
+        elif self.lstyle is not None:
+            mpl.plot(x,y,self.lstyle,marker=self.marker,label=self.label)
         else:
             mpl.plot(x,y,marker=self.marker,label=self.label)
@@ -672 +675 @@
             ## draw colorbar. settings are different with projections. or if not mapmode.
             if not self.mapmode: orientation=zeorientation ; frac = 0.075 ; pad = 0.03 ; lu = 0.5
-            elif self.proj in ['moll']: orientation="horizontal" ; frac = 0.075 ; pad = 0.03 ; lu = 1.0
+            elif self.proj in ['moll']: orientation="horizontal" ; frac = 0.08 ; pad = 0.03 ; lu = 1.0
             elif self.proj in ['cyl']: orientation="vertical" ; frac = 0.023 ; pad = 0.03 ; lu = 0.5
             else: orientation = zeorientation ; frac = zefrac ; pad = 0.03 ; lu = 0.5
@@ -690 +693 @@
         ### not expecting NaN in self.addvecx and self.addvecy. masked arrays is just enough.
         stride = 3
+        #stride = 1
         if self.addvecx is not None and self.addvecy is not None and self.mapmode:
                 ### for metwinds only ???
@@ -696 +700 @@
                 # reference vector is scaled
                 zescale = ppcompute.mean(np.sqrt(self.addvecx*self.addvecx+self.addvecy*self.addvecy))
+                #zescale = 25.
                 # make vector field
                 q = m.quiver( x[::stride,::stride],y[::stride,::stride],\
@@ -704 +709 @@
                 #keyh = 1.025 ; keyv = 1.05 # upper right corner over colorbar
                 keyh = 0.98 ; keyv = 1.06
+                keyh = 0.98 ; keyv = 1.07
                 #keyh = -0.03 ; keyv = 1.08 # upper left corner
                 p = mpl.quiverkey(q,keyh,keyv,\
                                   zescale,str(int(zescale)),\
                                   color='black',labelpos='S',labelsep = 0.07)
+        ############################################################################################
+        ### TEXT. ANYWHERE. add_text.txt should be present with lines x ; y ; text ; color
+        ############################################################################################
+        try:
+            f = open("add_text.txt", 'r')
+            for line in f:
+              if "#" in line: pass
+              else:
+                  userx, usery, usert, userc = line.strip().split(';')
+                  userc = userc.strip()
+                  usert = usert.strip()
+                  userx = float(userx.strip())
+                  usery = float(usery.strip())
+                  if self.mapmode: userx,usery = m(userx,usery)
+                  mpl.text(userx,usery,usert,\
+                           color = userc,\
+                           horizontalalignment='center',\
+                           verticalalignment='center')
+            f.close()
+        except IOError:
+            pass

trunk/UTIL/PYTHON/planetoplot_v2/set_area.txt

@@ -24 +24 @@
 Xanadu ; -40 ; 20 ; 40 ; 120
 Hyperboreae ; 80 ; 87 ; -70 ; -10
-VM ; -15 ; 0 ; -90 ; -60
+VM ; -18 ; 0 ; -95 ; -40
 VM_alt ; -15 ; 0 ; 270 ; 300

trunk/UTIL/PYTHON/planetoplot_v2/set_var.txt

@@ -58 +58 @@
 ISR ; %.0f ; gist_heat ; Incoming Solar Radiation ; W m$^{-2}$
 SWDOWNZ ; %.0f ; gist_heat ; $\Phi^{SW}\downarrow$ ; W m$^{-2}$
+ZMAX_TH ; %.1f ; spectral ; PBL depth ; km