source: trunk/UTIL/PYTHON/pp.py @ 647

Last change on this file since 647 was 647, checked in by tnavarro, 13 years ago

Corrected bug in reducefield for masked arrays with grid area. Possibility to see stream function for a lat/vert slice with --stream option. Output of both data and axis with save txt option in 1D. Small bug corrected: title is the file name for multiple plots with multiple plots. Cleanup in myplot.py
  • Property svn:executable set to *
File size: 9.4 KB
Line 
1#!/usr/bin/env python
2
3### A. Spiga + T. Navarro + A. Colaitis
4
5###########################################################################################
6###########################################################################################
7### What is below relate to running the file as a command line executable (very convenient)
8if __name__ == "__main__":
9    import sys
10    from optparse import OptionParser    ### to be replaced by argparse
11    from api_wrapper import api_onelevel
12    from gcm_transformations import call_zrecast
13    from netCDF4 import Dataset
14    from myplot import getlschar, separatenames, readslices, adjust_length, whatkindfile, errormess
15    from os import system
16    from planetoplot import planetoplot
17    from myscript import getparseroptions
18    import glob
19    import numpy as np
20
21
22    #############################
23    ### Get options and variables
24    parser = OptionParser() ; getparseroptions(parser) ; (opt,args) = parser.parse_args()
25    if opt.file is None:                                errormess("I want to eat one file at least ! Use winds.py -f name_of_my_file. Or type winds.py -h")
26    if opt.var is None and opt.anomaly is True:         errormess("Cannot ask to compute anomaly if no variable is set")
27    if opt.fref is not None and opt.operat is None:     errormess("you must specify an operation when using a reference file")
28    if opt.operat in ["+","-"] and opt.fref is None:    errormess("you must specifiy a reference file when using inter-file operations")
29    if opt.fref is not None and opt.operat is not None and opt.itp is not None: interpref=True
30    else:   interpref=False
31    if opt.rate is not None:      opt.save = "avi"
32    elif opt.save == "avi":       opt.rate = 8   ## this is a default value for -S avi
33    if opt.save == "html":        opt.rate = -1  ## this is convenient because everything is done in planetoplot with mrate
34
35    #############################
36    ### Get infos about slices
37    zeslat  = readslices(opt.slat) ; zeslon  = readslices(opt.slon) ; zesvert = readslices(opt.svert) ; zestime = readslices(opt.stime)
38    reffile = opt.fref
39    zexaxis = [opt.xmin,opt.xmax] ; zeyaxis=[opt.ymin,opt.ymax]
40
41    #############################
42    ### Catch multiple files
43    if "*" in opt.file[0] or "?" in opt.file[0]: 
44        yeah = glob.glob(opt.file[0]) ; yeah.sort()
45        opt.file[0] = yeah[0] 
46        for file in yeah[1:]: opt.file[0] = opt.file[0] + "," + file
47
48    #############################
49    ### 1. LOOP ON FILE LISTS TO BE PUT IN DIFFERENT FIGURES
50    for i in range(len(opt.file)):
51
52      zefiles = separatenames(opt.file[i])
53
54      typefile = whatkindfile(Dataset(zefiles[0])) ; stralt = None
55      if typefile in ["meso"]:         
56          [lschar,zehour,zehourin] = getlschar ( zefiles[0] )
57          if opt.var is None:  opt.var = ["HGT"] ; opt.clb = "nobar"
58      elif typefile in ["geo"]:
59          [lschar,zehour,zehourin] = ["",0,0]
60          if opt.var is None:  opt.var = ["HGT_M"] ; opt.clb = "nobar"
61      else:                                       
62          [lschar,zehour,zehourin] = ["",0,0]
63          if opt.var is None:  opt.var = ["phisinit"] ; opt.clb = "nobar"
64
65      if opt.vmin is not None : zevmin  = opt.vmin[min(i,len(opt.vmin)-1)]
66      else:                     zevmin = None
67      if opt.vmax is not None : zevmax  = opt.vmax[min(i,len(opt.vmax)-1)]
68      else:                     zevmax = None
69      #print "vmin, zevmin", opt.vmin, zevmin ; print "vmax, zevmax", opt.vmax, zevmax
70
71      #############################
72      ### 2. LOOP ON VAR LISTS TO BE PUT IN DIFFERENT FIGURES
73      for j in range(len(opt.var)):
74
75        zevars = separatenames(opt.var[j])
76
77        inputnvert = separatenames(opt.lvl)
78        if np.array(inputnvert).size == 1:
79            zelevel = float(inputnvert[0])
80            ze_interp_levels = [-9999.]
81        elif np.array(inputnvert).size > 2:
82            zelevel = -99.
83            start = float(inputnvert[0])
84            stop = float(inputnvert[1])
85            if np.array(inputnvert).size == 2:  numsample = 20
86            else:                               numsample = float(inputnvert[2])
87            if stop > start:   
88               ## altitude coordinates
89               ze_interp_levels = np.linspace(start,stop,numsample)
90            else:
91               ## pressure coordinates
92               ze_interp_levels = np.logspace(np.log10(start),np.log10(stop),numsample)
93
94        #########################################################
95        if opt.itp is not None:
96         if opt.itp > 0:
97          #####
98          ##### MESOSCALE : written by AS
99          #####
100          if typefile in ["meso"]:
101            if zelevel == 0. and opt.itp == 4:  zelevel = 0.010
102            ### winds or no winds
103            if opt.winds            :  zefields = 'uvmet'
104            elif opt.var[j] in ['deltat','DELTAT'] : zefields = 'tk,TSURF'
105            else                    :  zefields = ''
106            ### var or no var
107            if zefields == ''       :  zefields = opt.var[j] 
108            else                    :  zefields = zefields + "," + opt.var[j]
109            if opt.var2 is not None :  zefields = zefields + "," + opt.var2 
110            ### call fortran routines
111            for fff in range(len(zefiles)):
112                newname = api_onelevel (  path_to_input   = '', \
113                                               input_name      = zefiles[fff], \
114                                               fields          = zefields, \
115                                               interp_method   = opt.itp, \
116                                               interp_level    = ze_interp_levels, \
117                                               onelevel        = zelevel, \
118                                               nocall          = opt.nocall )
119                if fff == 0: zetab = newname
120                else:        zetab = np.append(zetab,newname)
121            if interpref:
122                reffile = api_onelevel (  path_to_input   = '', \
123                                               input_name      = opt.fref, \
124                                               fields          = zefields, \
125                                               interp_method   = opt.itp, \
126                                               interp_level    = ze_interp_levels, \
127                                               onelevel        = zelevel, \
128                                               nocall          = opt.nocall )
129            zefiles = zetab #; print zefiles
130            zelevel = 0 ## so that zelevel could play again the role of nvert
131          #####
132          ##### GCM : written by AC
133          #####
134          elif typefile == "gcm":
135            inputvar = zevars
136            if opt.var2 is not None : inputvar = np.append(inputvar,opt.var2)
137            interpolated_files=""
138            interpolated_files=call_zrecast(interp_mode=opt.itp,\
139                    input_name=zefiles,\
140                    fields=inputvar,\
141                    limites = ze_interp_levels,\
142                    predefined=opt.intas)
143
144            zefiles=interpolated_files
145            if interpref:
146               interpolated_ref=""
147               interpolated_ref=call_zrecast(interp_mode=opt.itp,\
148                    input_name=[opt.fref],\
149                    fields=zevars,\
150                    predefined=opt.intas)
151
152               reffile=interpolated_ref[0]
153          else:
154            print "type not supported"
155            exit()
156
157        #############
158        ### Main call
159        name = planetoplot (zefiles,level=int(zelevel),vertmode=opt.itp,\
160                proj=opt.proj,back=opt.back,target=opt.tgt,stride=opt.ste,var=zevars,\
161                colorb=opt.clb,winds=opt.winds,\
162                addchar=lschar,interv=[zehour,zehourin],vmin=zevmin,vmax=zevmax,\
163                tile=opt.tile,zoom=opt.zoom,display=opt.display,\
164                hole=opt.hole,save=opt.save,\
165                anomaly=opt.anomaly,var2=opt.var2,ndiv=opt.ndiv,\
166                mult=opt.mult,zetitle=separatenames(opt.zetitle),\
167                slon=zeslon,slat=zeslat,svert=zesvert,stime=zestime,\
168                outputname=opt.out,resolution=opt.res,\
169                ope=opt.operat,fileref=reffile,minop=opt.vminope,maxop=opt.vmaxope,titleref=opt.titref,\
170                invert_y=opt.inverty,xaxis=zexaxis,yaxis=zeyaxis,ylog=opt.logy,yintegral=opt.column,\
171                blat=opt.blat,blon=opt.blon,tsat=opt.tsat,flagnolow=opt.nolow,\
172                mrate=opt.rate,mquality=opt.quality,trans=opt.trans,zarea=opt.area,axtime=opt.axtime,\
173                redope=opt.redope,seevar=opt.seevar,xlab=opt.xlab,ylab=opt.ylab,lbls=separatenames(opt.labels),\
174                lstyle=separatenames(opt.linestyle),cross=readslices(opt.mark),streamflag=opt.stream)
175        print 'DONE: '+name
176        system("rm -f to_be_erased")
177 
178    #########################################################
179    ### Generate a .sh file with the used command saved in it
180    command = "" 
181    for arg in sys.argv: command = command + arg + ' '
182    #if typefile not in ["meso","mesoapi"]: name = 'pycommand'
183    if opt.save == "gui":    name = 'pycommand'
184    elif opt.save == "avi":  system("mv -f movie*.avi "+name+".avi")
185    elif opt.save == "html": system("cat $PYTHONPATH/header.html > anim.html ; cat zepics >> anim.html ; cat $PYTHONPATH/body.html >> anim.html ; rm -rf zepics "+name+" ; mkdir "+name+" ; mv anim.html image*png "+name) 
186    f = open(name+'.sh', 'w')
187    f.write(command)
188
189    #print "********** OPTIONS: ", opt
190    print "********************************************************** END"
Note: See TracBrowser for help on using the repository browser.