source: trunk/UTIL/PYTHON/gcm.py @ 361

Last change on this file since 361 was 361, checked in by acolaitis, 14 years ago

PYTHON.

Added zrecast capabilites for gcm.py. Interp modes AGL and P are working. Adding more is not difficult.
One should make sure that his zrecast.e is up to date and in his $PATH.

Note: we could also use f2py to compile a python module... However, most people do have zrecast in their path.

Note2: When using inter-file comparisons, such as field differences, one should make sure to use zrecast (or -i with gcm.py) to
make sure that operations are done at the same altitudes ! (levels in different simulations a priori do not represent same altitudes)

Note3: zrecast has been modified in previous commit to work with teta instead of temp, if needed.
  • Property svn:executable set to *
File size: 13.1 KB
Line 
1#!/usr/bin/env python
2
3### T. Navarro + A. Spiga + 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 zrecast_wrapper import call_zrecast
13    from netCDF4 import Dataset
14    from myplot import getlschar, separatenames, readslices, adjust_length
15    from os import system
16    from planetoplot import planetoplot
17   
18    import numpy as np
19
20    #############################
21    ### Get options and variables
22    parser = OptionParser()
23    parser.add_option('-f', '--file',   action='append',dest='namefile', type="string",  default=None,  help='[NEEDED] name of WRF file (append). Plot files separated by comas in the same figure')
24    parser.add_option('-l', '--level',  action='store',dest='nvert',     type="string",  default="0",   help='level or start,stop,step (def=0)(-i 2: p,mbar)(-i 3,4: z,km)')
25    parser.add_option('-p', '--proj',   action='store',dest='proj',      type="string",  default=None,  help='projection')
26    parser.add_option('-b', '--back',   action='store',dest='back',      type="string",  default=None,  help='background image (def: None)')
27    parser.add_option('-t', '--target', action='store',dest='target',    type="string",  default=None,  help='destination folder')
28    parser.add_option('-s', '--stride', action='store',dest='stride',    type="int",     default=3,     help='stride vectors (def=3)')
29    parser.add_option('-v', '--var',    action='append',dest='var',      type="string",  default=None,  help='variable color-shaded (append)')
30    parser.add_option('-n', '--num',    action='store',dest='numplot',   type="int",     default=None,  help='plot number (def=2)(<0: plot LT -*numplot*)')
31    parser.add_option('-i', '--interp', action='store',dest='interp',    type="int",     default=None,  help='interpolation (2: p, 3: z-amr, 4:z-als)')
32    parser.add_option('-c', '--color',  action='store',dest='colorb',    type="string",  default="def", help='change colormap (nobar: no colorbar)')
33    parser.add_option('-x', '--no-vect',action='store_false',dest='winds',               default=True,  help='no wind vectors')
34    parser.add_option('-m', '--min',    action='append',dest='vmin',     type="float",   default=None,  help='bounding minimum value (append)')   
35    parser.add_option('-M', '--max',    action='append',dest='vmax',     type="float",   default=None,  help='bounding maximum value (append)') 
36    parser.add_option('-T', '--tiled',  action='store_true',dest='tile',                 default=False, help='draw a tiled plot (no blank zone)')
37    parser.add_option('-z', '--zoom',   action='store',dest='zoom',      type="float",   default=None,  help='zoom factor in %')
38    parser.add_option('-N', '--no-api', action='store_true',dest='nocall',               default=False, help='do not recreate api file')
39    parser.add_option('-d', '--display',action='store_false',dest='display',             default=True,  help='do not pop up created images')
40    parser.add_option('-e', '--itstep', action='store',dest='itstep',    type="int",     default=None,  help='stride time (def=4)')
41    parser.add_option('-H', '--hole',   action='store_true',dest='hole',                 default=False, help='holes above max and below min')
42    parser.add_option('-S', '--save',   action='store',dest='save',      type="string",  default="gui", help='save mode (png,eps,svg,pdf or gui)(def=gui)')
43    parser.add_option('-a', '--anomaly',action='store_true',dest='anomaly',              default=False, help='compute and plot relative anomaly in %')
44    parser.add_option('-w', '--with',   action='store',dest='var2',      type="string",  default=None,  help='variable contoured')
45    parser.add_option('--div',          action='store',dest='ndiv',      type="int",     default=10,    help='number of divisions in colorbar (def: 10)')
46    parser.add_option('-F', '--first',  action='store',dest='first',     type="int",     default=1,     help='first subscript to plot (def: 1)')
47    parser.add_option('--mult',         action='store',dest='mult',      type="float",   default=1.,    help='a multiplicative factor to plotted field')
48    parser.add_option('--title',        action='store',dest='zetitle',   type="string",  default="fill",help='customize the whole title')
49    #parser.add_option('-V', action='store', dest='comb',        type="float",   default=None,  help='a defined combination of variables')
50
51    parser.add_option('-O','--output',  action='store',dest='output',    type="string",  default=None,  help='Output file name')
52    parser.add_option('--res',          action='store',dest='res',       type="float",   default=200.,  help='Resolution for png outputs. --save png must be specified. (def=200.)')
53    parser.add_option('--operation',    action='store',dest='operation', type="string",  default=None,  help='Basic operation to perform between input file and a specified reference file. Available: addition "+", soustraction "-". Operation is done as: file1 +or- fileref. The reference file must be specified using -fref option')
54    parser.add_option('--fref',         action='store',dest='fref',     type="string",  default=None,  help='Reference namefile for the --operation option.')
55    parser.add_option('--mope',         action='store',dest='vminope', type="float",   default=0.,  help='bounding minimum value for inter-file operation')
56    parser.add_option('--Mope',         action='store',dest='vmaxope', type="float",   default=0.,  help='bounding maximum value for inter-file operation')
57    parser.add_option('--titleref',     action='store',dest='titref',  type="string",  default="fill",  help='title for the reference file. Default: title of figure (1)')
58
59    ############# T.N. changes
60    #parser.add_option('-o','--operation',action='store',dest='operation',type="string",   default=None ,help='matrix of operations between files (for now see code, sorry)')
61    parser.add_option('--lat',          action='append',dest='slat',type="string",   default=None, help='slices along latitude. One value, or two values separated by comas for averaging')
62    parser.add_option('--lon',          action='append',dest='slon', type="string",   default=None, help='slices along longitude. One value, or two values separated by comas for averaging')
63    parser.add_option('--vert',         action='append',dest='svert',type="string",   default=None, help='slices along vertical axis. One value, or two values separated by comas for averaging') # quelles coordonnees ?
64    parser.add_option('--time',         action='append',dest='stime',type="string",   default=None, help='slices along time. One value, or two values separated by comas for averaging') # quelles coordonnees ?
65
66    (opt,args) = parser.parse_args()
67    if opt.namefile is None: 
68        print "I want to eat one file at least ! Use winds.py -f name_of_my_file. Or type winds.py -h"
69        exit()
70    if opt.var is None and opt.anomaly is True:
71        print "Cannot ask to compute anomaly if no variable is set"
72        exit()
73    print "Options:", opt
74
75#    listvar = ''
76#    if opt.var is None:
77#        zerange = [-999999]
78#    else:
79#        zelen = len(opt.var)
80#        zerange = range(zelen)
81#        #if zelen == 1: listvar = opt.var[0] + ','
82#        #else         :
83#        for jjj in zerange: listvar += opt.var[jjj] + ','
84#        listvar = listvar[0:len(listvar)-1]
85#        vmintab = adjust_length (opt.vmin, zelen)
86#        vmaxtab = adjust_length (opt.vmax, zelen)
87
88
89    ################################
90    ### General check
91 
92    if opt.fref is not None:
93       if opt.operation is None:
94          print "you must specify an operation when using a reference file"
95          exit()
96    if opt.operation is not None:
97       if opt.fref is None:
98          print "you must specifiy a reference file when using inter-file operations"
99          exit()
100
101    if opt.fref is not None:
102       if opt.operation is not None:
103          if opt.interp is not None:
104             interpref=True
105
106    ################################
107
108
109    print "namefile, length", opt.namefile, len(opt.namefile)
110
111    zeslat  = readslices(opt.slat)
112    zeslon  = readslices(opt.slon)
113    zesvert = readslices(opt.svert)
114    zestime = readslices(opt.stime)
115    print "slat,zeslat", opt.slat, zeslat
116    print "slon,zeslon", opt.slon, zeslon
117    print "svert,zesvert", opt.svert, zesvert
118    print "stime,zestime", opt.stime, zestime
119
120    for i in range(len(opt.namefile)):
121      for j in range(len(opt.var)):
122
123        zenamefiles = separatenames(opt.namefile[i])
124        print "zenamefiles", zenamefiles
125       
126        if opt.vmin is not None : zevmin  = opt.vmin[min(i,len(opt.vmin)-1)]
127        else: zevmin = None
128        if opt.vmax is not None : zevmax  = opt.vmax[min(i,len(opt.vmax)-1)]
129        else: zevmax = None
130        print "vmin, zevmin", opt.vmin, zevmin
131        print "vmax, zevmax", opt.vmax, zevmax
132       
133        zevar = separatenames(opt.var[j])
134        zevars = zevar[j]
135        zevar = zevar[0]
136        print "var, zevar", opt.var, zevar
137       
138        #checkcoherence(len(zenamefiles),len(opt.slat),len(opt.slon),len(opt.stime))
139       
140        zefile = zenamefiles[0]
141             
142        #zelevel = opt.nvert   
143        stralt = None
144        #[lschar,zehour,zehourin] = getlschar ( zefile )  ## getlschar from wrfout (or simply return "" if another file)
145        [lschar,zehour,zehourin] = ["",0,0]  ## dummy
146   
147        inputnvert = separatenames(opt.nvert)
148        if np.array(inputnvert).size == 1:
149            zelevel = float(inputnvert[0])
150            ze_interp_levels = [-9999.]
151        else:
152            zelevel = -99.
153            ze_interp_levels = np.linspace(float(inputnvert[0]),float(inputnvert[1]),float(inputnvert[2]))
154        print 'level: ', zelevel
155        print 'interp_levels: ',ze_interp_levels
156
157        #####################################################
158        ### Call Fortran routines for vertical interpolations  -->  zrecast for GCM ?       
159#        if opt.interp is not None:
160#            if zelevel == 0. and opt.interp == 4:  zelevel = 0.010
161#            ### winds or no winds
162#            if opt.winds            :  zefields = 'uvmet'
163#            else                    :  zefields = ''
164#            ### var or no var
165#            #if opt.var is None      :  pass
166#            if zefields == ''       :  zefields = listvar
167#            else                    :  zefields = zefields + "," + listvar
168#            if opt.var2 is not None : zefields = zefields + "," + opt.var2 
169#            print zefields
170#            zefile = api_onelevel (  path_to_input   = '', \
171#                                     input_name      = zefile, \
172#                                     fields          = zefields, \
173#                                     interp_method   = opt.interp, \
174#                                     interp_level    = ze_interp_levels, \
175#                                     onelevel        = zelevel, \
176#                                     nocall          = opt.nocall )
177#            print zefile
178#            zelevel = 0 ## so that zelevel could play again the role of nvert
179
180# A.C.  #####################################################
181        ### Call Fortran routines for vertical interpolations  -->  zrecast
182
183        if opt.interp is not None:
184           interpolated_files=""
185           interpolated_files=call_zrecast(interp_mode=opt.interp,\
186                    input_name=zenamefiles,\
187                    fields=zevars)
188
189           zenamefiles=interpolated_files
190        if interpref:
191           interpolated_ref=""
192           interpolated_ref=call_zrecast(interp_mode=opt.interp,\
193                    input_name=[opt.fref],\
194                    fields=zevars)
195
196           reffile=interpolated_ref[0]
197        else:
198           reffile=opt.fref
199
200        #############
201        ### Main call
202        name = planetoplot (zenamefiles,level=int(zelevel),vertmode=opt.interp,\
203                proj=opt.proj,back=opt.back,target=opt.target,stride=opt.stride,var=zevar,\
204                numplot=opt.numplot,colorb=opt.colorb,winds=opt.winds,\
205                addchar=lschar,interv=[zehour,zehourin],vmin=zevmin,vmax=zevmax,\
206                tile=opt.tile,zoom=opt.zoom,display=opt.display,\
207                itstep=opt.itstep,hole=opt.hole,save=opt.save,\
208                anomaly=opt.anomaly,var2=opt.var2,ndiv=opt.ndiv,first=opt.first,\
209                mult=opt.mult,zetitle=opt.zetitle,\
210                slon=zeslon,slat=zeslat,svert=zesvert,stime=zestime,\
211                outputname=opt.output,resolution=opt.res,\
212                ope=opt.operation,fileref=reffile,minop=opt.vminope,maxop=opt.vmaxope,titleref=opt.titref)
213        print 'Done: '+name
214        system("rm -f to_be_erased")
215 
216    #########################################################
217    ### Generate a .sh file with the used command saved in it
218    command = ""
219    for arg in sys.argv: command = command + arg + ' '
220    name = 'pycommand'
221    f = open(name+'.sh', 'w')
222    f.write(command)
Note: See TracBrowser for help on using the repository browser.