source: lmdz_wrf/trunk/tools/drawing.py @ 1071

Last change on this file since 1071 was 1071, checked in by lfita, 8 years ago

Getting on `draw_2lines'

File size: 166.6 KB
Line 
1import numpy as np
2import os
3from netCDF4 import Dataset as NetCDFFile
4import drawing_tools as drw
5import generic_tools as gen
6from optparse import OptionParser
7import sys
8from cStringIO import StringIO
9
10## e.g. # drawing.py -f /media/data2/etudes/WRF_LMDZ/WL_HyMeX/IIphase/medic950116/wlmdza/wrfout/wrfout_d01_1995-01-13_00:00:00 -o create_movie -S 'draw_2D_shad#Time@WRFTimes@10@95@191@1#tas:East_West|-1,North_South|-1,Time|2:longitude:latitude:Summer:270.,300.:tas|at|t=0:pdf:None:None' -v T2
11## e.g. # drawing.py -f wrfout_d01_1980-03-01_00\:00\:00_Time_B0-E48-I1.nc -o draw_2D_shad -S 'tas:East_West|-1,North_South|-1,Time|2:longitude:latitude:Summer:270.,300.:tas|at|t=0:pdf:None:cyl,i' -v T2
12## e.g. # drawing.py -f ~/etudes/domains/MEDCORDEX/geo_em.d01.nc -o draw_2D_shad_cont -S 'landmask,height:Time|0:Time|0:XLONG_M:XLAT_M:rainbow:fixc,k:%.2f:0,1:0.,3000.,10:landmask & height:pdf:False:lcc,i' -v LANDMASK,HGT_M
13## e.g. # drawing.py -f ~/etudes/domains/MEDCORDEX/geo_em.d01.nc -o draw_2D_shad_cont -S 'height,landmask:Time|0:Time|0:XLONG_M:XLAT_M:terrain:fixc,k:None:0.,3000.:0,1,10:MEDCORDEX height & landmask:pdf:False:lcc,i' -v HGT_M,LANDMASK
14## e.g. # drawing.py -o draw_2D_shad_line -f 'mean_dtcon-pluc-pres_lat.nc,mean_dtcon-pluc-pres_lat.nc' -S 'dtcon,prc:bottom_top|-1,south_north|-1:latmean:presmean:seismic,k:-5.,5.:monthly|dtcon|&|prc:pdf:flip@y:None:True' -v 'dtconmean,prcmean'
15## e.g. # drawing.py -f 'geo_em.d02.nc' -o draw_topo_geogrid -S '0.,3000.:None:FF_3dom d02:png:cyl,i'
16## e.g. # drawing.py -o draw_topo_geogrid_boxes -S '0.,3000.:None:FF domains:pdf:lcc,i:d01,d02,d03:0' -f geo_em.d01.nc,geo_em.d02.nc,geo_em.d03.nc
17## e.g. # drawing.py -o draw_trajectories -f 'WRF/control/trajectory.dat@-1@/home/lluis/etudes/domains/WL_HyMeX_HighRes_C/geo_em.d03.nc#XLONG_M#XLAT_M,WRF_LMDZ/wlmdza/trajectory.dat@-1@/home/lluis/etudes/domains/WL_HyMeX_C/geo_em.d01.nc#XLONG_M#XLAT_M,WRF_LMDZ/wlmdzb/trajectory.dat@-1@/home/lluis/etudes/domains/WL_HyMeX_C/geo_em.d01.nc#XLONG_M#XLAT_M,WRF_LMDZ/wlmdzb_ii/trajectory.dat@-1@/home/lluis/etudes/domains/WL_HyMeX_C/geo_em.d01.nc#XLONG_M#XLAT_M' -S '$WRF_{CRM}$,$LMDZ_{AR4.0}$,$LMDZ_{NPv3.1}$,$LMDZ_{NPv3.1b}$|None|medicane trajectories|pdf|cyl,i' -v obs/trajectory.dat,satellite,-1
18## e.g. # drawing.py -o draw_vals_trajectories -f WRF_LMDZ/wlmdza/tevolboxtraj_T2.nc,WRF_LMDZ/wlmdzb/tevolboxtraj_T2.nc,WRF/control/tevolboxtraj_T2.nc -S 'mean:-1:$WRF_{CRM}$,$LMDZ_{AR4.0}$,$LMDZ_{NPv3.1}$@4:tas:time|($[DD]^[HH]$):exct,6,h:$%d^{%H}$:trajectory|following|mean:pdf' -v T2
19## e.g. # drawing.py -o draw_2D_shad_time -f 'netcdf_concatenated.nc' -S 'dtcon:Time|-1,bottom_top|-1:presmean:time:seismic:-3.e-6,3.e-6:monthly|dtcon:pdf:transpose:time|hours!since!1949-12-01|exct,2,d|%d|date!([DD])|x:True' -v 'dtconmean'
20## e.g. # drawing.py -o variable_values -S PSFC
21## e.g. # drawing.py -o draw_timeSeries -f wrfout_d01_1979-12-01_00:00:00_bottom_top_B6-E6-I1_south_north_B3-E3-I1_west_east_B26-E26-I1.nc -S 'dt_con:time|($[DD]^{[HH]}$):exct,12,h:$%d^{%H}$:time|evolution|at|-1|6|3|26:1:pdf' -v 'LDQCON,time'
22## e.g. # drawing.py -f wrfout_d01_1979-12-01_00:00:00 -o draw_Neighbourghood_evol -S 'q:Time|-1|Times,bottom_top|6|ZNU,south_north|3|XLAT,west_east|26|XLONG:south_north,west_east:5:auto:time|($[DD]^{[HH]}$),time|($[DD]^{[HH]}$):exct,2,h|exct,1,d:$%d^{%H}$,$%d^{%H}$:5|pts|neighbourghood|temporal|evolution:0.0,0.004:BuPu:pdf:True' -v QVAPOR
23## e.g. # drawing.py -o draw_lines_time -f wrfout_d01_1980-03-01_00:00:00_Time_B0-E48-I1_south_north_B15-E15-I1_west_east_B15-E15-I1.nc -S 'time;y;time ([DD]${[HH]}$);file1;tas;evolution;time|hours!since!1949-12-01_00:00:00|exct,12,h|%d$^{%H}$;pdf' -v T2
24## e.g. # drawing.py -o draw_barbs -f ERAI_pl199501_131-132.nc -S 'X|lon|lon|-1,Y|lat|lat|-1,Z|lev|lev|4,T|time|time|0:auto,auto,auto:wind,ms-1:cyl,c:ERA-Interim|winds|at|1000|hPa|on|1996|January|1st|at|00|UTC:pdf:ERAI_pl199501_131-132' -v var131,var132
25## e.g. # ~/etudes/WRF_LMDZ/svn/LMDZ_WRF/tools/drawing.py -o draw_points -S 'tslist.dat,#,3,2,1:SuperStorm|sfc|stations:cyl,i:legend:auto:None:0:png:stations_loc' -f $HOME'/etudes/domains/WRFdynamicoSENS_SuperStorm/geo_em.d02.nc,XLONG_M,XLAT_M,HGT_M,Time|0,height,0.,3000.,terrain,m'
26## e.g. # drawing.py -o draw_points -S 'tslist.dat,#,3,2,1:SuperStorm|sfc|stations:cyl,i:labelled,8,black:auto:None:0:png:stations_loc' -f $HOME'/etudes/domains/WRFdynamicoSENS_SuperStorm/geo_em_west_east_B25-E180-I1_south_north_B160-E262-I1.nc,XLONG_M,XLAT_M,HGT_M,Time|0,height,0.,1500.,terrain,m'
27## e.g. # drawing.py -o draw_ptZvals -f geo_v2_2012102123_RR1.nc -S 'pracc:lon,lat:o:80:2,42,7,47,:values!of!values:Blues:cyl,l:pdf' -v pr
28## e.g. # drawing.py -f carteveg5km.nc -o draw_points_lonlat -S 'longitude:latitude:pdf:points!veget|type:green:.:0.5:None:0:legend'
29## e.g. # drawing.py -o draw_vectors -f wrfout_d01_2001-11-11_00:00:00 -S 'T|Time|Times|2,Y|south_north|XLAT|-1,X|west_east|XLONG|-1:3@3,wind@rainbow,9:10m wind,ms-1:cyl,l:WRF 10 m winds:pdf:winds' -v U10,V10
30## e.g. # drawing.py -o draw_basins -f routing.py -S 'Y|y|nav_lat|-1,X|x|nav_lon|-1:1@1,rainbow,9:basins,-:cyl,l:ORCDHIEE river-basins:pdf:basins_named' -v nav_lon,nav_lat,trip,basins
31## e.g. # drawing.py -o draw_river_desc -f diver_desc.nc -S 'Y|lat|lat|-1,X|lon|lon|-1:red,green:Blues:cyl,l:ORCDHIEE rivers:pdf:0:or_rivers -v Amazon
32## e.g. # drawing.py -o draw_vertical_levels -f wrfout_d01_2001-11-11_00:00:00 -S 'true:false:wrfout!vertical!levels!(standard!40):png:4' -v WRFz
33## e.g. $ drawing.py -o draw_subbasin -f Caceres_subbasin.nc -S 'Caceres:None:cyl,l:2,True:Caceres:pdf:0:Caceres_subbasin'
34
35main = 'drawing.py'
36
37errormsg = 'ERROR -- error -- ERROR -- error'
38warnmsg = 'WARNING -- waring -- WARNING -- warning'
39fillValue=1.e20
40
41namegraphics = ['create_movie', 'draw_2D_shad', 'draw_2D_shad_time',                 \
42  'draw_2D_shad_cont', 'draw_2D_shad_cont_time', 'draw_2D_shad_line',                \
43  'draw_2D_shad_line_time', 'draw_barbs', 'draw_basins',                             \
44  'draw_lines', 'draw_lines_time', 'draw_Neighbourghood_evol',                       \
45  'draw_points', 'draw_points_lonlat',                                               \
46  'draw_ptZvals', 'draw_river_desc', 'draw_subbasin', 'draw_timeSeries',             \
47  'draw_topo_geogrid',                                                               \
48  'draw_topo_geogrid_boxes', 'draw_trajectories', 'draw_vals_trajectories',          \
49  'draw_vectors',  'draw_vertical_levels', 'list_graphics', 'variable_values']
50
51def draw_2D_shad(ncfile, values, varn):
52    """ plotting a fields with shading
53    draw_2D_shad(ncfile, values, varn)
54      ncfile= file to use
55      values=[vnamefs]:[dimvals]:[dimxvn]:[dimyvn]:[colorbar]:[sminv],[smaxv]:[figt]:
56       [kindfig]:[reverse]:[mapv]:[close]
57        [vnamefs]: Name in the figure of the variable to be shaded
58        [dimvals]: ',' list of [dimname]|[value] telling at which dimension of the
59          variable a given value is required (-1, all the length)
60        [dimx/yvn]: name of the variables with the values of the final dimensions (x,y)
61        [colorbar]: name of the color bar
62        [smin/axv]: minimum and maximum value for the shading or:
63          'Srange': for full range
64          'Saroundmean@val': for mean-xtrm,mean+xtrm where xtrm = np.min(mean-min@val,max@val-mean)
65          'Saroundminmax@val': for min*val,max*val
66          'Saroundpercentile@val': for median-xtrm,median+xtrm where xtrm = np.min(median-percentile_(val),
67            percentile_(100-val)-median)
68          'Smean@val': for -xtrm,xtrm where xtrm = np.min(mean-min*@val,max*@val-mean)
69          'Smedian@val': for -xtrm,xtrm where xtrm = np.min(median-min@val,max@val-median)
70          'Spercentile@val': for -xtrm,xtrm where xtrm = np.min(median-percentile_(val),
71             percentile_(100-val)-median)
72        [figt]: title of the figure ('|' for spaces)
73        [kindfig]: kind of figure
74        [reverse]: Transformation of the values
75          * 'transpose': reverse the axes (x-->y, y-->x)
76          * 'flip'@[x/y]: flip the axis x or y
77        [mapv]: map characteristics: [proj],[res]
78          see full documentation: http://matplotlib.org/basemap/
79          [proj]: projection
80            * 'cyl', cilindric
81            * 'lcc', lamvbert conformal
82          [res]: resolution:
83            * 'c', crude
84            * 'l', low
85            * 'i', intermediate
86            * 'h', high
87            * 'f', full
88      valules= 'rh:z|-1,x|-1:z|-1,x|-1:lat:pressure:BuPu:0.,100.:rh:pdf:flip@y:None'
89      varn= [varsn] name of the variable to plot with shading
90    """
91
92    fname = 'draw_2D_shad'
93    if values == 'h':
94        print fname + '_____________________________________________________________'
95        print draw_2D_shad.__doc__
96        quit()
97
98    expectargs = '[vnamefs]:[dimvals]:[dimxvn]:[dimyvn]:[colorbar]:[sminv],[smaxv]'+\
99      ':[figt]:[kindfig]:[reverse]:[mapv]:[close]'
100 
101    drw.check_arguments(fname,values,expectargs,':')
102
103    vnamesfig = values.split(':')[0]
104    dimvals= values.split(':')[1].replace('|',':')
105    vdimxn = values.split(':')[2]
106    vdimyn = values.split(':')[3]
107    colbarn = values.split(':')[4]
108    shadminmax = values.split(':')[5]
109    figtitle = values.split(':')[6].replace('|',' ')
110    figkind = values.split(':')[7]
111    revals = values.split(':')[8]
112    mapvalue = values.split(':')[9]
113#    varn = values.split(':')[10]
114
115    ncfiles = ncfile
116   
117    if not os.path.isfile(ncfiles):
118        print errormsg
119        print '  ' + fname + ': shading file "' + ncfiles + '" does not exist !!'
120        quit(-1)   
121
122    objsf = NetCDFFile(ncfiles, 'r')
123   
124    varns = varn.split(',')[0]
125
126    if  not objsf.variables.has_key(varns):
127        print errormsg
128        print '  ' + fname + ': shading file "' + ncfiles +                          \
129          '" does not have variable "' +  varns + '" !!'
130        quit(-1)
131
132# Variables' values
133    objvars = objsf.variables[varns]
134
135    valshad, dimsshad = drw.slice_variable(objvars, dimvals.replace(',','|'))
136
137# Dimensions names
138##    print fname + ' obj dimnames: ', objvars.dimensions, dimvals, len(dimvals.split(','))
139##    dimnamesv = []
140##    for idd in range(len(objvars.dimensions)):
141##        cutdim = False
142##        for idc in range(len(dimvals.split(','))):
143##            dimcutn = dimvals.split(',')[idc].split(':')[0]
144##            print objvars.dimensions[idd], dimcutn
145##            if objvars.dimensions[idd] == dimcutn:
146##                cutdim = True
147##                break
148##        if not cutdim: dimnamesv.append(objvars.dimensions[idd])
149    dimnamesv = [vdimxn, vdimyn]
150
151    if drw.searchInlist(objvars.ncattrs(),'units'):
152        varunits = objvars.getncattr('units')
153    else:
154        print warnmsg
155        print '  ' + fname + ": variable '" + varn + "' without units!!"
156        varunits = '-'
157
158    if  not objsf.variables.has_key(vdimxn):
159        print errormsg
160        print '  ' + fname + ': shading file "' + ncfiles +                          \
161          '" does not have dimension variable "' +  vdimxn + '" !!'
162        quit(-1)
163    if  not objsf.variables.has_key(vdimyn):
164        print errormsg
165        print '  ' + fname + ': shading file "' + ncfiles +                          \
166          '" does not have dimension variable "' +  vdimyn + '" !!'
167        quit(-1)
168
169    objdimx = objsf.variables[vdimxn]
170    objdimy = objsf.variables[vdimyn]
171    if drw.searchInlist(objdimx.ncattrs(),'units'):
172        odimxu = objdimx.getncattr('units')
173    else:
174        print warnmsg
175        print '  ' + fname + ": variable dimension '" + vdimxn + "' without units!!"
176        odimxu = '-'
177
178    if drw.searchInlist(objdimy.ncattrs(),'units'):
179        odimyu = objdimy.getncattr('units')
180    else:
181        print warnmsg
182        print '  ' + fname + ": variable dimension '" + vdimyn + "' without units!!"
183        odimyu = '-'
184
185    odimxv, odimyv = drw.dxdy_lonlatDIMS(objdimx[:], objdimy[:], objdimx.dimensions,     \
186      objdimy.dimensions, dimvals.replace(':','|').split(','))
187
188
189#    if len(objdimx.shape) <= 2:
190##        odimxv = objdimx[valshad.shape]
191##        odimyv = objdimy[valshad.shape]
192#        odimxv = objdimx[:]
193#        odimyv = objdimy[:]
194
195#    elif len(objdimx.shape) == 3:
196##        dimcut = [0, slice(0,valshad.shape[0]), slice(0,valshad.shape[1])]
197##        odimxv = objdimx[tuple(dimcut)]
198##        odimyv = objdimy[tuple(dimcut)]
199#        odimxv = objdimx[0,:]
200#        odimyv = objdimy[0,:]
201#    else:
202#        print errormsg
203#        print '  ' + fname + ': shape of dimension variable:', objdimx.shape,        \
204#          ' not ready!!'
205#        quit(-1)
206
207    shading_nx = []
208    if shadminmax.split(',')[0][0:1] != 'S':
209            shading_nx.append(np.float(shadminmax.split(',')[0]))
210    else:
211        shading_nx.append(shadminmax.split(',')[0])
212
213    if shadminmax.split(',')[1][0:1] != 'S':
214        shading_nx.append(np.float(shadminmax.split(',')[1]))
215    else:
216        shading_nx.append(shadminmax.split(',')[1])
217
218    if mapvalue == 'None': mapvalue = None
219
220    drw.plot_2D_shadow(valshad, vnamesfig, odimxv, odimyv, odimxu, odimyu, dimnamesv,\
221      colbarn, shading_nx, varunits, figtitle, figkind, revals, mapvalue, True)
222
223    return
224
225def draw_2D_shad_time(ncfile, values, varn):
226    """ plotting a fields with shading with time values
227    draw_2D_shad(ncfile, values, varn)
228      ncfile= file to use
229      values=[vnamefs]~[dimvals]~[dimxvn]~[dimyvn]~[colorbar]~[sminv],[smaxv]~[figt]~
230       [kindfig]~[reverse]~[timevals]~[close]
231        [vnamefs]: Name in the figure of the variable to be shaded
232        [dimvals]: ',' list of [dimname]|[value] telling at which dimension of the
233          variable a given value is required (-1, all the length, [beg]@[end] for an interval)
234        [dimx/yvn]: name of the variables with the values of the final dimensions (x,y)
235        [colorbar]: name of the color bar
236        [smin/axv]: minimum and maximum value for the shading or:
237          'Srange': for full range
238          'Saroundmean@val': for mean-xtrm,mean+xtrm where xtrm = np.min(mean-min@val,max@val-mean)
239          'Saroundminmax@val': for min*val,max*val
240          'Saroundpercentile@val': for median-xtrm,median+xtrm where xtrm = np.min(median-percentile_(val),
241            percentile_(100-val)-median)
242          'Smean@val': for -xtrm,xtrm where xtrm = np.min(mean-min*@val,max*@val-mean)
243          'Smedian@val': for -xtrm,xtrm where xtrm = np.min(median-min@val,max@val-median)
244          'Spercentile@val': for -xtrm,xtrm where xtrm = np.min(median-percentile_(val),
245             percentile_(100-val)-median)
246        [figt]: title of the figure ('|' for spaces)
247        [kindfig]: kind of figure
248        [reverse]: Transformation of the values
249          * 'transpose': reverse the axes (x-->y, y-->x)
250          * 'flip'@[x/y]: flip the axis x or y
251        [timevals]: [timen]|[units]|[kind]|[tfmt]|[label]|[timeaxis] time labels characteristics
252           [timen]; name of the time variable
253           [units]; units string according to CF conventions ([tunits] since
254             [YYYY]-[MM]-[DD] [[HH]:[MI]:[SS]], '!' for spaces)
255           [kind]; kind of output
256             'Nval': according to a given number of values as 'Nval',[Nval]
257             'exct': according to an exact time unit as 'exct',[tunit];
258               tunit= [Nunits],[tu]; [tu]= 'c': centuries, 'y': year, 'm': month,
259                'w': week, 'd': day, 'h': hour, 'i': minute, 's': second,
260                'l': milisecond
261           [tfmt]; desired format
262           [label]; label at the graph ('!' for spaces)
263        [close]: should figure be closed (finished)
264      values='dtcon~Time|-1,bottom_top|-1~presmean~time~seismic~-3.e-6,3.e-6~monthly|'
265        'dtcon~pdf~transpose~time|hours!since!1949-12-01|exct,2,d|%d|date!([DD])~True
266      varn= [varsn] name of the variable to plot with shading
267    """
268    fname = 'draw_2D_shad_time'
269
270    if values == 'h':
271        print fname + '_____________________________________________________________'
272        print draw_2D_shad_time.__doc__
273        quit()
274
275    farguments = '[vnamefs]~[dimvals]~[dimxvn]~[dimyvn]~[colorbar]~[sminv],[smaxv]~'+\
276      '[figt]~[kindfig]~[reverse]~[timevals]~[close]'
277    drw.check_arguments(fname,values,farguments,'~')
278
279    vnamesfig = values.split('~')[0]
280    dimvals= values.split('~')[1].replace('|',':')
281    vdimxn = values.split('~')[2]
282    vdimyn = values.split('~')[3]
283    colbarn = values.split('~')[4]
284    shadminmax = values.split('~')[5]
285    figtitle = values.split('~')[6].replace('|',' ')
286    figkind = values.split('~')[7]
287    revals = values.split('~')[8]
288    timevals = values.split('~')[9]
289    close = values.split('~')[10]
290
291    ncfiles = ncfile
292   
293    if not os.path.isfile(ncfiles):
294        print errormsg
295        print '  ' + fname + ': shading file "' + ncfiles + '" does not exist !!'
296        quit(-1)   
297
298    objsf = NetCDFFile(ncfiles, 'r')
299   
300    varns = varn.split(',')[0]
301
302    if  not objsf.variables.has_key(varns):
303        print errormsg
304        print '  ' + fname + ': shading file "' + ncfiles +                          \
305          '" does not have variable "' +  varns + '" !!'
306        quit(-1)
307
308# Variables' values
309    objvars = objsf.variables[varns]
310
311    valshad, dimsshad = drw.slice_variable(objvars, dimvals.replace(',','|'))
312
313    dimnamesv = [vdimxn, vdimyn]
314
315    varunits = objvars.getncattr('units')
316
317    if  not objsf.variables.has_key(vdimxn):
318        print errormsg
319        print '  ' + fname + ': shading file "' + ncfiles +                          \
320          '" does not have dimension variable "' +  vdimxn + '" !!'
321        quit(-1)
322    if  not objsf.variables.has_key(vdimyn):
323        print errormsg
324        print '  ' + fname + ': shading file "' + ncfiles +                          \
325          '" does not have dimension variable "' +  vdimyn + '" !!'
326        quit(-1)
327
328    objdimx = objsf.variables[vdimxn]
329    objdimy = objsf.variables[vdimyn]
330    odimxu = objdimx.getncattr('units')
331    odimyu = objdimy.getncattr('units')
332
333    if len(objdimx.shape) <= 2:
334        odimxv0 = objdimx[:]
335        odimyv0 = objdimy[:]
336
337    elif len(objdimx.shape) == 3:
338        odimxv0 = objdimx[0,:]
339        odimyv0 = objdimy[0,:]
340    else:
341        print errormsg
342        print '  ' + fname + ': shape of dimension variable:', objdimx.shape,        \
343          ' not ready!!'
344        quit(-1)
345
346    timename = timevals.split('|')[0]
347    timeunit = timevals.split('|')[1].replace('!',' ')
348    timekind = timevals.split('|')[2]
349    timefmt = timevals.split('|')[3]
350    timelabel = timevals.split('|')[4].replace('!',' ')
351
352# Dimensional values
353    odxv, dimsdxv = drw.slice_variable(objsf.variables[vdimxn], dimvals.replace(',','|'))
354    odyv, dimsdyv = drw.slice_variable(objsf.variables[vdimyn], dimvals.replace(',','|'))
355
356    if vdimxn == timename:
357        odimxv = objsf.variables[vdimxn][:]
358        odimxu = timelabel
359        timeaxis = 'x'
360        odimyv = objsf.variables[vdimyn]
361        odimyu = odimyv.getncattr('units')
362        timepos, timelabels = drw.CFtimes_plot(odxv, timeunit, timekind, timefmt)
363    elif vdimyn == timename:
364        odimyv = objsf.variables[vdimyn]
365        odimyu = timelabel
366        timeaxis = 'y'
367        odimxv = objsf.variables[vdimxn]
368        odimxu = odimxv.getncattr('units')
369        timepos, timelabels = drw.CFtimes_plot(odyv, timeunit, timekind, timefmt)
370    else:
371        print errormsg
372        print '  ' + fname + ": time variable '" + timename + "' not found!!"
373        quit(-1)
374
375    shading_nx = []
376    if shadminmax.split(',')[0][0:1] != 'S':
377        shading_nx.append(np.float(shadminmax.split(',')[0]))
378    else:
379        shading_nx.append(shadminmax.split(',')[0])
380
381    if shadminmax.split(',')[1][0:1] != 'S':
382        shading_nx.append(np.float(shadminmax.split(',')[1]))
383    else:
384        shading_nx.append(shadminmax.split(',')[1])
385
386    closeval = drw.Str_Bool(close)
387
388    drw.plot_2D_shadow_time(valshad, vnamesfig, odxv, odyv, odimxu, odimyu,          \
389      dimnamesv, colbarn, shading_nx, varunits, figtitle, figkind, revals, timeaxis, \
390      timepos, timelabels, closeval)
391
392    return
393
394def draw_2D_shad_cont(ncfile, values, varn):
395    """ plotting two fields, one with shading and the other with contour lines
396    draw_2D_shad_cont(ncfile, values, varn)
397      ncfile= [ncfilevars],[ncfilevarc] files to use (one value, same file)
398      values=[vnamefs]:[dimvals]:[dimvalc]:[dimxvn]:[dimyvn]:[colorbar]:[ckind]:[clabfmt]:[sminv],[smaxv]:[sminc],[smaxv],[Nlev]:[figt]:[kindfig]:[reverse]:[mapv]
399        [vnamefs],[vnamefc]: Name in the figure of the shaded and the contour variables
400        [dimvals/c]: list of [dimname]|[value] telling at which dimension of the
401          variable a given value is required (no dimension name, all the length)
402        [dimx/yvn]: names of the variables with the values of the dimensions for the plot
403        [colorbar]: name of the color bar
404        [ckind]: kind of contours
405          'cmap': as it gets from colorbar
406          'fixc,[colname]': fixed color [colname], all stright lines
407          'fixsigc,[colname]': fixed color [colname], >0 stright, <0 dashed  line
408        [clabfmt]: format of the labels in the contour (None, also possible)
409        [smin/axv]: minimum and maximum value for the shading
410        [sminc]:[smaxv]:[Nlev]: minimum, maximum and number of values for the contour
411        [figt]: title of the figure ('|' for spaces)
412        [kindfig]: kind of figure
413        [reverse]: does the values be transposed? 'True/False',
414        [mapv]: map characteristics: [proj],[res]
415          see full documentation: http://matplotlib.org/basemap/
416          [proj]: projection
417            * 'cyl', cilindric
418            * 'lcc', lamvbert conformal
419          [res]: resolution:
420            * 'c', crude
421            * 'l', low
422            * 'i', intermediate
423            * 'h', high
424            * 'f', full
425      valules= 'rh,ta:z|-1,x|-1:z|-1,x|-1:lat:pressure:BuPu:fixsigc,black:%d:0.,100.:195.,305.,7:Meridonal|average|of|rh|&|ta:pdf:flip@y:None'
426      varn= [varsn],[varcn] name of the variable to plot with shading variable with contour
427    """
428
429    fname = 'draw_2D_shad_cont'
430    if values == 'h':
431        print fname + '_____________________________________________________________'
432        print draw_2D_shad_cont.__doc__
433        quit()
434
435    expectargs = '[vnamefs]:[dimvals]:[dimvalc]:[dimxvn]:[dimyvn]:[colorbar]:'
436    expectargs = expectargs + '[ckind]:[clabfmt]:[sminv],[smaxv]:[sminc],[smaxv],'
437    expectargs = expectargs + '[Nlev]:[figt]:[kindfig]:[reverse]:[mapv]'
438 
439    drw.check_arguments(fname,values,expectargs,':')
440
441    vnamesfig = values.split(':')[0].split(',')
442    dimvals= values.split(':')[1].replace('|',':')
443    dimvalc= values.split(':')[2].replace('|',':')
444    vdimxn = values.split(':')[3]
445    vdimyn = values.split(':')[4]
446    colbarn = values.split(':')[5]
447    countkind = values.split(':')[6]
448    countlabelfmt = values.split(':')[7]
449    shadminmax = values.split(':')[8]
450    contlevels = values.split(':')[9]
451    figtitle = values.split(':')[10].replace('|',' ')
452    figkind = values.split(':')[11]
453    revals = values.split(':')[12]
454    mapvalue = values.split(':')[13]
455
456    if2filenames = ncfile.find(',')
457
458    if if2filenames != -1:
459        ncfiles = ncfile.split(',')[0]
460        ncfilec = ncfile.split(',')[1]
461    else:
462        ncfiles = ncfile
463        ncfilec = ncfile
464
465    if not os.path.isfile(ncfiles):
466        print errormsg
467        print '  ' + fname + ': shading file "' + ncfiles + '" does not exist !!'
468        quit(-1)   
469
470    if not os.path.isfile(ncfilec):
471        print errormsg
472        print '  ' + fname + ': contour file "' + ncfilec + '" does not exist !!'
473        quit(-1)   
474
475    objsf = NetCDFFile(ncfiles, 'r')
476    objcf = NetCDFFile(ncfilec, 'r')
477   
478    varns = varn.split(',')[0]
479    varnc = varn.split(',')[1]
480
481    if  not objsf.variables.has_key(varns):
482        print errormsg
483        print '  ' + fname + ': shading file "' + ncfiles +                          \
484          '" does not have variable "' +  varns + '" !!'
485        quit(-1)
486
487    if  not objcf.variables.has_key(varnc):
488        print errormsg
489        print '  ' + fname + ': contour file "' + ncfilec +                          \
490          '" does not have variable "' +  varnc + '" !!'
491        quit(-1)
492
493# Variables' values
494    objvars = objsf.variables[varns]
495    objvarc = objcf.variables[varnc]
496
497    if len(objvars.shape) != len(objvarc.shape):
498        print errormsg
499        print '  ' + fname + ': shading variable "' + varns + '" has a shape: ',     \
500          objvars.shape,  'different than contour variable "' +  varnc + '": ',      \
501          objvarc.shape,' !!!'
502        quit(-1)
503
504    for idim in range(len(objvars.shape)):
505        if objvars.shape[idim] != objvarc.shape[idim]:
506            print errormsg
507            print '  ' + fname + ': shading variable "' + varns + '" has a shape: ', \
508              objvars.shape,  'different than contour variable "' +  varnc + '": ',  \
509              objvarc.shape,' !!!'
510            quit(-1)
511
512    valshad, dimsshad = drw.slice_variable(objvars, dimvals.replace(',','|'))
513    valcont, dimscont = drw.slice_variable(objvarc, dimvalc.replace(',','|'))
514
515# Dimensions names
516##    print fname + ' obj dimnames: ', objvars.dimensions, dimvals, len(dimvals.split(','))
517##    dimnamesv = []
518##    for idd in range(len(objvars.dimensions)):
519##        cutdim = False
520##        for idc in range(len(dimvals.split(','))):
521##            dimcutn = dimvals.split(',')[idc].split(':')[0]
522##            print objvars.dimensions[idd], dimcutn
523##            if objvars.dimensions[idd] == dimcutn:
524##                cutdim = True
525##                break
526##        if not cutdim: dimnamesv.append(objvars.dimensions[idd])
527    dimnamesv = [vdimxn, vdimyn]
528
529    varunits = []
530    varunits.append(objvars.getncattr('units'))
531    varunits.append(objvarc.getncattr('units'))
532
533    if  not objsf.variables.has_key(vdimxn):
534        print errormsg
535        print '  ' + fname + ': shading file "' + ncfiles +                          \
536          '" does not have dimension variable "' +  vdimxn + '" !!'
537        quit(-1)
538    if  not objsf.variables.has_key(vdimyn):
539        print errormsg
540        print '  ' + fname + ': shading file "' + ncfiles +                          \
541          '" does not have dimension variable "' +  vdimyn + '" !!'
542        quit(-1)
543
544    objdimx = objsf.variables[vdimxn]
545    objdimy = objsf.variables[vdimyn]
546    odimxu = objdimx.getncattr('units')
547    odimyu = objdimy.getncattr('units')
548
549# Getting only that dimensions with coincident names
550    odimxv, odimyv = drw.dxdy_lonlatDIMS(objdimx[:], objdimy[:], objdimx.dimensions,     \
551      objdimy.dimensions, dimvals.replace(':','|').split(','))
552
553#    dimnvx = objdimx.dimensions
554#    cutslice = []
555#    for idimn in objdimx.dimensions:
556#        found = False
557#        for dimsn in dimsshad:
558#            if idimn == dimsn:
559#                cutslice.append(slice(0,len(objsf.dimensions[idimn])))
560#                found = True
561#        if not found: cutslice.append(0)
562#
563#    odimxv = objdimx[tuple(cutslice)]
564#
565#    dimnvy = objdimy.dimensions
566#    cutslice = []
567#    for idimn in objdimy.dimensions:
568#        found = False
569#        for dimsn in dimsshad:
570#            if idimn == dimsn:
571#                cutslice.append(slice(0,len(objsf.dimensions[idimn])))
572#                found = True
573#        if not found: cutslice.append(0)
574#
575#    odimyv = objdimy[tuple(cutslice)]
576
577#    if len(objdimx.shape) <= 2:
578#        odimxv = objdimx[:]
579#        odimyv = objdimy[:]
580#    elif len(objdimx.shape) == 3:
581#        odimxv = objdimx[0,:]
582#        odimyv = objdimy[0,:]
583#    else:
584#        print errormsg
585#        print '  ' + fname + ': shape of dimension variable:', objdimx.shape,        \
586#          ' not ready!!'
587#        quit(-1)
588
589    if countlabelfmt == 'None': 
590        countlfmt = None
591    else:
592        countlfmt = countlabelfmt
593
594    shading_nx = np.zeros((2), dtype=np.float)
595    shading_nx[0] = np.float(shadminmax.split(',')[0])
596    shading_nx[1] = np.float(shadminmax.split(',')[1])
597
598    clevmin = np.float(contlevels.split(',')[0])
599    clevmax = np.float(contlevels.split(',')[1])
600    Nclevels = int(contlevels.split(',')[2])
601
602    levels_cont = drw.pretty_int(clevmin, clevmax, Nclevels)
603
604    if len(levels_cont) <= 1: 
605        print warnmsg
606        print '  ' + fname + ': wrong contour levels:', levels_cont,' !!'
607        del(levels_cont)
608        levels_cont = np.zeros((Nclevels), dtype=np.float)
609        levels_cont = np.arange(7)*(clevmax - clevmin)/(Nclevels-1)
610        print '    generating default ones: ',levels_cont
611
612    if mapvalue == 'None': mapvalue = None
613    if revals == 'None': revals = None
614
615    drw.plot_2D_shadow_contour(valshad, valcont, vnamesfig, odimxv, odimyv, odimxu,  \
616      odimyu, dimnamesv, colbarn, countkind, countlfmt, shading_nx, levels_cont,     \
617      varunits, figtitle, figkind, revals, mapvalue)
618
619    return
620
621def draw_2D_shad_cont_time(ncfile, values, varn):
622    """ plotting two fields, one with shading and the other with contour lines being
623    one of the dimensions of time characteristics
624    draw_2D_shad_cont(ncfile, values, varn)
625      ncfile= [ncfilevars],[ncfilevarc] files to use (one value, same file)
626      values=[vnamefs];[dimvals];[dimvalc];[dimxvn];[dimyvn];[colorbar];[ckind];[clabfmt];[sminv],[smaxv];[sminc],[smaxv],[Nlev];[figt];[kindfig];[reverse];[timevals]
627        [vnamefs],[vnamefc]: Name in the figure of the shaded and the contour variables
628        [dimvals/c]: list of [dimname]|[value] telling at which dimension of the
629          variable a given value is required (no dimension name, all the length)
630        [dimxvn]: name of the variables with the values of the dimension of the x-axis
631        [dimyvn]: name of the variables with the values of the dimension of the y-axis
632        [colorbar]: name of the color bar
633        [ckind]: kind of contours
634          'cmap': as it gets from colorbar
635          'fixc,[colname]': fixed color [colname], all stright lines
636          'fixsigc,[colname]': fixed color [colname], >0 stright, <0 dashed  line
637        [clabfmt]: format of the labels in the contour (None, also possible)
638        [smin/axv]: minimum and maximum value for the shading
639        [sminc]:[smaxv]:[Nlev]: minimum, maximum and number of values for the contour
640        [figt]: title of the figure ('|' for spaces)
641        [kindfig]: kind of figure
642        [reverse]: modification to the dimensions:
643          'transposed': transpose matrices
644          'flip',[x/y]: flip only the dimension [x] or [y]
645        [timevals]: [timen]|[units]|[kind]|[tfmt]|[label] time labels characteristics
646           [timen]; name of the time variable
647           [units]; units string according to CF conventions ([tunits] since
648             [YYYY]-[MM]-[DD] [[HH]:[MI]:[SS]], '!' for spaces)
649           [kind]; kind of output
650             'Nval': according to a given number of values as 'Nval',[Nval]
651             'exct': according to an exact time unit as 'exct',[tunit];
652               tunit= [Nunits],[tu]; [tu]= 'c': centuries, 'y': year, 'm': month,
653                'w': week, 'd': day, 'h': hour, 'i': minute, 's': second,
654                'l': milisecond
655           [tfmt]; desired format
656           [label]; label at the graph ('!' for spaces)
657      valules= 'rh,ta;z|-1,x|-1;z|-1,x|-1;lat;pressure;BuPu;fixsigc,black;%d;0.,100.;195.,305.,7;Meridonal|average|of|rh|&|ta;pdf;flip@y;time!hours!since!1949/12/01|exct,5d|%d|date!([DD])'
658      varn= [varsn],[varcn] name of the variable to plot with shading variable with contour
659    """
660
661    fname = 'draw_2D_shad_cont_time'
662    if values == 'h':
663        print fname + '_____________________________________________________________'
664        print draw_2D_shad_cont_time.__doc__
665        quit()
666
667    expectargs = '[vnamefs];[dimvals];[dimvalc];[dimxvn];[dimyvn];[colorbar];' +     \
668      '[ckind];[clabfmt];[sminv],[smaxv];[sminc],[smaxv],[Nlev];[figt];[kindfig];' + \
669      '[reverse];[timevals]'
670 
671    drw.check_arguments(fname,values,expectargs,';')
672
673    vnamesfig = values.split(';')[0].split(',')
674    dimvals= values.split(';')[1].replace('|',':')
675    dimvalc= values.split(';')[2].replace('|',':')
676    vdimxn = values.split(';')[3]
677    vdimyn = values.split(';')[4]
678    colbarn = values.split(';')[5]
679    countkind = values.split(';')[6]
680    countlabelfmt = values.split(';')[7]
681    shadminmax = values.split(';')[8]
682    contlevels = values.split(';')[9]
683    figtitle = values.split(';')[10].replace('|',' ')
684    figkind = values.split(';')[11]
685    revals = values.split(';')[12]
686    timevals = values.split(';')[13]
687
688    if2filenames = ncfile.find(',')
689
690    if if2filenames != -1:
691        ncfiles = ncfile.split(',')[0]
692        ncfilec = ncfile.split(',')[1]
693    else:
694        ncfiles = ncfile
695        ncfilec = ncfile
696
697    if not os.path.isfile(ncfiles):
698        print errormsg
699        print '  ' + fname + ': shading file "' + ncfiles + '" does not exist !!'
700        quit(-1)   
701
702    if not os.path.isfile(ncfilec):
703        print errormsg
704        print '  ' + fname + ': contour file "' + ncfilec + '" does not exist !!'
705        quit(-1)   
706
707    objsf = NetCDFFile(ncfiles, 'r')
708    objcf = NetCDFFile(ncfilec, 'r')
709   
710    varns = varn.split(',')[0]
711    varnc = varn.split(',')[1]
712
713    if  not objsf.variables.has_key(varns):
714        print errormsg
715        print '  ' + fname + ': shading file "' + ncfiles +                          \
716          '" does not have variable "' +  varns + '" !!'
717        quit(-1)
718
719    if  not objcf.variables.has_key(varnc):
720        print errormsg
721        print '  ' + fname + ': contour file "' + ncfilec +                          \
722          '" does not have variable "' +  varnc + '" !!'
723        quit(-1)
724
725# Variables' values
726    objvars = objsf.variables[varns]
727    objvarc = objcf.variables[varnc]
728
729    if len(objvars.shape) != len(objvarc.shape):
730        print errormsg
731        print '  ' + fname + ': shading variable "' + varns + '" has a shape: ',     \
732          objvars.shape,  'different than contour variable "' +  varnc + '": ',      \
733          objvarc.shape,' !!!'
734        quit(-1)
735
736    for idim in range(len(objvars.shape)):
737        if objvars.shape[idim] != objvarc.shape[idim]:
738            print errormsg
739            print '  ' + fname + ': shading variable "' + varns + '" has a shape: ', \
740              objvars.shape,  'different than contour variable "' +  varnc + '": ',  \
741              objvarc.shape,' !!!'
742            quit(-1)
743
744    valshad, dimsshad = drw.slice_variable(objvars, dimvals.replace(',','|'))
745    valcont, dimscont = drw.slice_variable(objvarc, dimvalc.replace(',','|'))
746
747    dimnamesv = [vdimxn, vdimyn]
748
749    varunits = []
750    varunits.append(objvars.getncattr('units'))
751    varunits.append(objvarc.getncattr('units'))
752
753    if  not objsf.variables.has_key(vdimxn):
754        print errormsg
755        print '  ' + fname + ': shading file "' + ncfiles +                          \
756          '" does not have dimension variable "' +  vdimxn + '" !!'
757        quit(-1)
758    if  not objsf.variables.has_key(vdimyn):
759        print errormsg
760        print '  ' + fname + ': shading file "' + ncfiles +                          \
761          '" does not have dimension variable "' +  vdimyn + '" !!'
762        quit(-1)
763
764    timename = timevals.split('|')[0]
765    timeunit = timevals.split('|')[1].replace('!',' ')
766    timekind = timevals.split('|')[2]
767    timefmt = timevals.split('|')[3]
768    timelabel = timevals.split('|')[4].replace('!',' ')
769
770    if vdimxn == timename:
771        timevals = objsf.variables[vdimxn][:]
772        timedims = objsf.variables[vdimxn].dimensions
773        dimt = 'x'
774        ovalaxis = objsf.variables[vdimyn]
775        ovalu = ovalaxis.getncattr('units')
776    elif vdimyn == timename:
777        timevals = objsf.variables[vdimyn][:]
778        timedims = objsf.variables[vdimyn].dimensions
779        dimt = 'y'
780        ovalaxis = objsf.variables[vdimxn]
781        ovalu = ovalaxis.getncattr('units')
782    else:
783        print errormsg
784        print '  ' + fname + ": time variable '" + timename + "' not found!!"
785        quit(-1)
786
787    timepos, timelabels = drw.CFtimes_plot(timevals, timeunit, timekind, timefmt)
788
789# Getting only that dimensions with coincident names
790    dimnvx = ovalaxis.dimensions
791
792    cutslice = []
793    for idimn in dimnvx:
794        found = False
795        for dimsn in dimsshad:
796            if idimn == dimsn:
797                cutslice.append(slice(0,len(objsf.dimensions[idimn])))
798                found = True
799        if not found: cutslice.append(slice(0,len(objsf.dimensions[idimn])))
800
801    ovalaxisv = ovalaxis[tuple(cutslice)]
802
803    if countlabelfmt == 'None': 
804        countlfmt = None
805    else:
806        countlfmt = countlabelfmt
807
808    shading_nx = np.zeros((2), dtype=np.float)
809    shading_nx[0] = np.float(shadminmax.split(',')[0])
810    shading_nx[1] = np.float(shadminmax.split(',')[1])
811
812    clevmin = np.float(contlevels.split(',')[0])
813    clevmax = np.float(contlevels.split(',')[1])
814    Nclevels = int(contlevels.split(',')[2])
815
816    levels_cont = drw.pretty_int(clevmin, clevmax, Nclevels)
817
818    if len(levels_cont) <= 1: 
819        print warnmsg
820        print '  ' + fname + ': wrong contour levels:', levels_cont,' !!'
821        del(levels_cont)
822        levels_cont = np.zeros((Nclevels), dtype=np.float)
823        levels_cont = np.arange(7)*(clevmax - clevmin)/(Nclevels-1)
824        print '    generating default ones: ',levels_cont
825
826    if revals == 'None': revals = None
827
828    drw.plot_2D_shadow_contour_time(valshad, valcont, vnamesfig, ovalaxisv,         \
829      timevals, timepos, timelabels, ovalu, timelabel, dimt, dimnamesv, colbarn,    \
830      countkind, countlfmt, shading_nx, levels_cont, varunits, figtitle, figkind,   \
831      revals)
832
833    return
834
835def draw_2D_shad_line(ncfile, values, varn):
836    """ plotting a fields with shading and another with line
837    draw_2D_shad_line(ncfile, values, varn)
838      ncfile= [ncfiles],[ncfilel] file to use for the shading and for the line
839      values=[vnamefs],[vnamefl]:[dimvals]:[dimxvn]:[dimyvn]:[colorbar],[colline]:[sminv],[smaxv]:[figt]:
840       [kindfig]:[reverse]:[mapv]:[close]
841        [vnamefs]: Name in the figure of the variable to be shaded
842        [vnamefl]: Name in the figure of the variable to be lined
843        [dimvals]: ',' list of [dimname]|[value] telling at which dimension of the
844          variable a given value is required (-1, all the length)
845        [dimx/yvn]: name of the variables with the values of the final dimensions (x,y)
846        [colorbar]: name of the color bar
847        [colline]: name of the color for the line
848        [smin/axv]: minimum and maximum value for the shading
849        [figt]: title of the figure ('|' for spaces)
850        [kindfig]: kind of figure
851        [reverse]: Transformation of the values
852          * 'transpose': reverse the axes (x-->y, y-->x)
853          * 'flip'@[x/y]: flip the axis x or y
854        [mapv]: map characteristics: [proj],[res]
855          see full documentation: http://matplotlib.org/basemap/
856          [proj]: projection
857            * 'cyl', cilindric
858            * 'lcc', lamvbert conformal
859          [res]: resolution:
860            * 'c', crude
861            * 'l', low
862            * 'i', intermediate
863            * 'h', high
864            * 'f', full
865      valules= 'rh:z|-1,x|-1:z|-1,x|-1:lat:pressure:BuPu:0.,100.:rh:pdf:flip@y:None'
866      varn= [varsn],[varnl] name of the variable to plot with shading and with line
867    """
868
869    fname = 'draw_2D_shad_line'
870    if values == 'h':
871        print fname + '_____________________________________________________________'
872        print draw_2D_shad_line.__doc__
873        quit()
874
875    farguments = '[vnamefs],[vnamefl]:[dimvals]:[dimxvn]:[dimyvn]:' +                \
876      '[colorbar],[colline]:[sminv],[smaxv]:[figt]:[kindfig]:[reverse]:' +           \
877      '[mapv]:[close]'
878    drw.check_arguments(fname,values,farguments,':')
879
880    vnamesfig = values.split(':')[0].split(',')[0]
881    dimvals= values.split(':')[1].replace('|',':')
882    vdimxn = values.split(':')[2]
883    vdimyn = values.split(':')[3]
884    colbarn = values.split(':')[4].split(',')[0]
885    shadminmax = values.split(':')[5]
886    figtitle = values.split(':')[6].replace('|',' ')
887    figkind = values.split(':')[7]
888    revals = values.split(':')[8]
889    mapvalue = values.split(':')[9]
890#    varn = values.split(':')[10]
891
892    ncfiles = ncfile.split(',')[0]
893   
894    if not os.path.isfile(ncfiles):
895        print errormsg
896        print '  ' + fname + ': shading file "' + ncfiles + '" does not exist !!'
897        quit(-1)   
898
899    objsf = NetCDFFile(ncfiles, 'r')
900   
901    varns = varn.split(',')[0]
902
903    if  not objsf.variables.has_key(varns):
904        print errormsg
905        print '  ' + fname + ': shading file "' + ncfiles +                          \
906          '" does not have variable "' +  varns + '" !!'
907        quit(-1)
908
909# Variables' values
910    objvars = objsf.variables[varns]
911
912    valshad, dimsshad = drw.slice_variable(objvars, dimvals.replace(',','|'))
913
914# Dimensions names
915##    print fname + ' obj dimnames: ', objvars.dimensions, dimvals, len(dimvals.split(','))
916##    dimnamesv = []
917##    for idd in range(len(objvars.dimensions)):
918##        cutdim = False
919##        for idc in range(len(dimvals.split(','))):
920##            dimcutn = dimvals.split(',')[idc].split(':')[0]
921##            print objvars.dimensions[idd], dimcutn
922##            if objvars.dimensions[idd] == dimcutn:
923##                cutdim = True
924##                break
925##        if not cutdim: dimnamesv.append(objvars.dimensions[idd])
926    dimnamesv = [vdimxn, vdimyn]
927
928    varunits = objvars.getncattr('units')
929
930    if  not objsf.variables.has_key(vdimxn):
931        print errormsg
932        print '  ' + fname + ': shading file "' + ncfiles +                          \
933          '" does not have dimension variable "' +  vdimxn + '" !!'
934        quit(-1)
935    if  not objsf.variables.has_key(vdimyn):
936        print errormsg
937        print '  ' + fname + ': shading file "' + ncfiles +                          \
938          '" does not have dimension variable "' +  vdimyn + '" !!'
939        quit(-1)
940
941    objdimx = objsf.variables[vdimxn]
942    objdimy = objsf.variables[vdimyn]
943    odimxu = objdimx.getncattr('units')
944    odimyu = objdimy.getncattr('units')
945
946    if len(objdimx.shape) <= 2:
947#        odimxv = objdimx[valshad.shape]
948#        odimyv = objdimy[valshad.shape]
949        odimxv = objdimx[:]
950        odimyv = objdimy[:]
951
952    elif len(objdimx.shape) == 3:
953#        dimcut = [0, slice(0,valshad.shape[0]), slice(0,valshad.shape[1])]
954#        odimxv = objdimx[tuple(dimcut)]
955#        odimyv = objdimy[tuple(dimcut)]
956        odimxv = objdimx[0,:]
957        odimyv = objdimy[0,:]
958    else:
959        print errormsg
960        print '  ' + fname + ': shape of dimension variable:', objdimx.shape,        \
961          ' not ready!!'
962        quit(-1)
963
964    shading_nx = np.zeros((2), dtype=np.float)
965    shading_nx[0] = np.float(shadminmax.split(',')[0])
966    shading_nx[1] = np.float(shadminmax.split(',')[1])
967
968    if mapvalue == 'None': mapvalue = None
969
970# line plot
971##
972    ncfilel = ncfile.split(',')[1]
973    vnamelfig = values.split(':')[0].split(',')[1]
974    varnl = varn.split(',')[1]
975    colline = values.split(':')[4].split(',')[1]
976
977    objlf = NetCDFFile(ncfilel,'r')
978    objlvar = objlf.variables[varnl]
979
980    linevals = objlvar[:]
981    varlunits = objlvar.units
982
983    drw.plot_2D_shadow_line(valshad, linevals, vnamesfig, vnamelfig, odimxv, odimyv, \
984      odimxu, odimyu, dimnamesv, colbarn, colline, shading_nx, varunits, varlunits,  \
985      figtitle, figkind, revals, mapvalue, True)
986
987    objsf.close()
988    objlf.close()
989
990    return
991
992def draw_2D_shad_line_time(ncfile, values, varn):
993    """ plotting a fields with shading and a line with time values
994    draw_2D_shad_line(ncfile, values, varn)
995      ncfile= [ncfiles],[ncfilel] files to use to draw with shading and the line
996      values= [vnamefs],[vanemefl]:[dimvals]:[dimxvn]:[dimyvn]:[colorbar]:[sminv],[smaxv]:[figt]:
997       [kindfig]:[reverse]:[timevals]:[close]
998        [vnamefs]: Name in the figure of the variable to be shaded
999        [vnamefl]: Name in the figure of the variable to be lined
1000        [dimvals]: ',' list of [dimname]|[value] telling at which dimension of the
1001          variable a given value is required (-1, all the length)
1002        [dimx/yvn]: name of the variables with the values of the final dimensions (x,y)
1003        [colorbar]: name of the color bar
1004        [smin/axv]: minimum and maximum value for the shading
1005        [figt]: title of the figure ('|' for spaces)
1006        [kindfig]: kind of figure
1007        [reverse]: Transformation of the values
1008          * 'transpose': reverse the axes (x-->y, y-->x)
1009          * 'flip'@[x/y]: flip the axis x or y
1010        [timevals]: [timen]|[units]|[kind]|[tfmt]|[label]|[timeaxis] time labels characteristics
1011           [timen]; name of the time variable
1012           [units]; units string according to CF conventions ([tunits] since
1013             [YYYY]-[MM]-[DD] [[HH]:[MI]:[SS]], '!' for spaces)
1014           [kind]; kind of output
1015             'Nval': according to a given number of values as 'Nval',[Nval]
1016             'exct': according to an exact time unit as 'exct',[tunit];
1017               tunit= [Nunits],[tu]; [tu]= 'c': centuries, 'y': year, 'm': month,
1018                'w': week, 'd': day, 'h': hour, 'i': minute, 's': second,
1019                'l': milisecond
1020           [tfmt]; desired format
1021           [label]; label at the graph ('!' for spaces)
1022        [close]: should figure be closed (finished)
1023      values='dtcon,prc:Time|-1,bottom_top|-1:presmean:time:seismic:-3.e-6,3.e-6:monthly|'
1024        'dtcon:pdf:transpose:time|hours!since!1949-12-01|exct,2,d|%d|date!([DD])|x:True
1025      varn= [varsn].[varln] name of the variable to plot with shading and to plot with line
1026    """
1027    fname = 'draw_2D_shad_line_time'
1028    if values == 'h':
1029        print fname + '_____________________________________________________________'
1030        print draw_2D_shad__line_time.__doc__
1031        quit()
1032
1033    farguments = '[vnamefs],[vanemefl]:[dimvals]:[dimxvn]:[dimyvn]:' +               \
1034      '[colorbar]:[sminv],[smaxv]:[figt]:[kindfig]:[reverse]:[timevals]:[close]'
1035    drw.check_arguments(fname,values,farguments,':')
1036
1037    vnamesfig = values.split(':')[0].split(',')[0]
1038    dimvals= values.split(':')[1].replace('|',':')
1039    vdimxn = values.split(':')[2]
1040    vdimyn = values.split(':')[3]
1041    colbarn = values.split(':')[4]
1042    shadminmax = values.split(':')[5]
1043    figtitle = values.split(':')[6].replace('|',' ')
1044    figkind = values.split(':')[7]
1045    revals = values.split(':')[8]
1046    timevals = values.split(':')[9]
1047    close = values.split(':')[10]
1048
1049    ncfiles = ncfile.split(',')[0]
1050   
1051    if not os.path.isfile(ncfiles):
1052        print errormsg
1053        print '  ' + fname + ': shading file "' + ncfiles + '" does not exist !!'
1054        quit(-1)   
1055
1056    objsf = NetCDFFile(ncfiles, 'r')
1057   
1058    varns = varn.split(',')[0]
1059
1060    if  not objsf.variables.has_key(varns):
1061        print errormsg
1062        print '  ' + fname + ': shading file "' + ncfiles +                          \
1063          '" does not have variable "' +  varns + '" !!'
1064        quit(-1)
1065
1066# Variables' values
1067    objvars = objsf.variables[varns]
1068
1069    valshad, dimsshad = drw.slice_variable(objvars, dimvals.replace(',','|'))
1070
1071    dimnamesv = [vdimxn, vdimyn]
1072
1073    varunits = objvars.getncattr('units')
1074
1075    if  not objsf.variables.has_key(vdimxn):
1076        print errormsg
1077        print '  ' + fname + ': shading file "' + ncfiles +                          \
1078          '" does not have dimension variable "' +  vdimxn + '" !!'
1079        quit(-1)
1080    if  not objsf.variables.has_key(vdimyn):
1081        print errormsg
1082        print '  ' + fname + ': shading file "' + ncfiles +                          \
1083          '" does not have dimension variable "' +  vdimyn + '" !!'
1084        quit(-1)
1085
1086    objdimx = objsf.variables[vdimxn]
1087    objdimy = objsf.variables[vdimyn]
1088    odimxu = objdimx.getncattr('units')
1089    odimyu = objdimy.getncattr('units')
1090
1091    if len(objdimx.shape) <= 2:
1092        odimxv = objdimx[:]
1093        odimyv = objdimy[:]
1094
1095    elif len(objdimx.shape) == 3:
1096        odimxv = objdimx[0,:]
1097        odimyv = objdimy[0,:]
1098    else:
1099        print errormsg
1100        print '  ' + fname + ': shape of dimension variable:', objdimx.shape,        \
1101          ' not ready!!'
1102        quit(-1)
1103
1104    timename = timevals.split('|')[0]
1105    timeunit = timevals.split('|')[1].replace('!',' ')
1106    timekind = timevals.split('|')[2]
1107    timefmt = timevals.split('|')[3]
1108    timelabel = timevals.split('|')[4].replace('!',' ')
1109
1110    if vdimxn == timename:
1111        odimxv = objsf.variables[vdimxn][:]
1112        odimxu = timelabel
1113        timeaxis = 'x'
1114        odimyv = objsf.variables[vdimyn]
1115        odimyu = odimyv.getncattr('units')
1116        timepos, timelabels = drw.CFtimes_plot(odimxv, timeunit, timekind, timefmt)
1117    elif vdimyn == timename:
1118        odimyv = objsf.variables[vdimyn][:]
1119        odimyu = timelabel
1120        timeaxis = 'y'
1121        odimxv = objsf.variables[vdimxn]
1122        odimxu = odimxv.getncattr('units')
1123        timepos, timelabels = drw.CFtimes_plot(odimyv, timeunit, timekind, timefmt)
1124    else:
1125        print errormsg
1126        print '  ' + fname + ": time variable '" + timename + "' not found!!"
1127        quit(-1)
1128
1129    shading_nx = np.zeros((2), dtype=np.float)
1130    shading_nx[0] = np.float(shadminmax.split(',')[0])
1131    shading_nx[1] = np.float(shadminmax.split(',')[1])
1132
1133    closeval = drw.Str_Bool(close)
1134
1135    drw.plot_2D_shadow_time(valshad, vnamesfig, odimxv, odimyv, odimxu, odimyu,      \
1136      dimnamesv, colbarn, shading_nx, varunits, figtitle, figkind, revals, timeaxis, \
1137      timepos, timelabels, False)
1138
1139# Line values
1140##
1141    ncfilel = ncfile.split(',')[1]
1142
1143    vnamelfig = values.split(':')[0].split(',')[1]
1144    varnl = varn.split(',')[1]
1145
1146    objlf = NetCDFFile(ncfilel,'r')
1147    objlvar = objlf.variables[varnl]
1148
1149    linevals = objlvar[:]
1150    if reva0 == 'tranpose':
1151        plt.plot (linevals, odimxv, '-', color='k')
1152    else:
1153        plt.plot (odimxv, linevals, '-', color='k')
1154
1155    objsf.close()
1156    objsl.close()
1157
1158    return
1159
1160def draw_barbs(ncfile, values, varns):
1161    """ Function to plot wind barbs
1162      values= [dimname]|[vardimname]|[value]:[vecvals]:[windlabs]:[mapvalues]:
1163        [gtit]:[kindfig]:[figuren]
1164        'X/Y/Z/T'|[dimname]|[vardimname]|[value]: ',', list for each basic dimension '|' separated of:
1165          [dimname]: name of the dimension in the file
1166          [vardimname]: name of the variable with the values for the dimension in the file
1167          [value]: which value of the given dimension (-1, all; [ibeg]@[iend], i-range beginning, end)
1168          No value takes all the range of the dimension
1169        [vecvals]= [frequency],[color],[length]
1170          [frequency]: [xfreq]@[yfreq] frequency of values allong each axis ('None', all grid points;
1171            'auto', computed automatically to have 20 vectors along each axis)
1172          [color]: color of the vectors ('auto', for 'red')
1173          [length]: length of the wind barbs ('auto', for 9)
1174        [windlabs]= [windname],[windunits]
1175          [windname]: name of the wind variable in the graph
1176          [windunits]: units of the wind variable in the graph ('None', for the value in the file)
1177        [mapvalues]= map characteristics: [proj],[res]
1178          see full documentation: http://matplotlib.org/basemap/
1179          [proj]: projection
1180            * 'cyl', cilindric
1181            * 'lcc', lambert conformal
1182          [res]: resolution:
1183            * 'c', crude
1184            * 'l', low
1185            * 'i', intermediate
1186            * 'h', high
1187            * 'f', full
1188        gtit= title of the graph ('|', for spaces)
1189        kindfig= kind of figure
1190        figuren= name of the figure
1191      ncfile= file to use
1192      varns= [uwind],[ywind] ',' list of the name of the variables with the u-wind,y-wind component
1193    """
1194    fname = 'draw_barbs'
1195
1196    if values == 'h':
1197        print fname + '_____________________________________________________________'
1198        print draw_barbs.__doc__
1199        quit()
1200
1201    expectargs = '[X/Y/Z/T]|[dimname]|[vardimname]|[value]:[vecvals]:[windlabs]:' +  \
1202      '[mapvalues]:[gtit]:[kindfig]:[figuren]'
1203 
1204    drw.check_arguments(fname,values,expectargs,':')
1205
1206    dimvals = values.split(':')[0]
1207    vecvals = values.split(':')[1]
1208    windlabels = values.split(':')[2]
1209    mapvalues = values.split(':')[3]
1210    gtit = values.split(':')[4]
1211    kindfig = values.split(':')[5]
1212    figuren = values.split(':')[6]
1213
1214    of = NetCDFFile(ncfile,'r')
1215
1216    dims = {}
1217    for dimv in dimvals.split(','):
1218        dns = dimv.split('|')
1219        dims[dns[0]] = [dns[1], dns[2], dns[3]]
1220
1221    varNs = []
1222    for dn in dims.keys():
1223        if dn == 'X':
1224            varNs.append(dims[dn][1])
1225            dimx = len(of.dimensions[dims[dn][0]])
1226        elif dn == 'Y':
1227            varNs.append(dims[dn][1])
1228            dimy = len(of.dimensions[dims[dn][0]])
1229
1230    ivar = 0
1231    for wvar in varns.split(','):
1232        if not drw.searchInlist(of.variables.keys(), wvar):
1233            print errormsg
1234            print '  ' + fname + ": file does not have variable '" + wvar + "' !!"
1235            quit(-1)
1236        if ivar == 0:
1237            varNs.append(wvar)
1238        else:
1239            varNs.append(wvar)
1240
1241    ivar = 0
1242    for varN in varNs:
1243        varslice = []
1244
1245        ovarN = of.variables[varN]
1246        vard = ovarN.dimensions
1247        for vdn in vard:
1248            found = False
1249            for dd in dims.keys():
1250                if dims[dd][0] == vdn:
1251                    if dims[dd][2].find('@') != -1:
1252                        rvals = dims[dd][2].split('@')
1253                        varslice.append(slice(int(rvals[0]), int(rvals[1])))
1254                    elif dims[dd][2] == '-1':
1255                        varslice.append(slice(0,len(of.dimensions[dims[dd][0]])))
1256                    else:
1257                        varslice.append(int(dims[dd][2]))
1258
1259                    found = True
1260                    break
1261            if not found:
1262                varslice.append(slice(0,len(of.dimensions[dims[dd][0]])))
1263
1264        if varN == dims['X'][1]:
1265            lonvals0 = np.squeeze(ovarN[tuple(varslice)])
1266        elif varN == dims['Y'][1]:
1267            latvals0 = np.squeeze(ovarN[tuple(varslice)])
1268        elif ivar == 2:
1269            uwvals = np.squeeze(np.array(ovarN[tuple(varslice)]))
1270        elif ivar == 3:
1271            vwvals = np.squeeze(ovarN[tuple(varslice)])           
1272
1273        ivar = ivar + 1
1274
1275#    print 'Final shapes:',lonvals0.shape,':',latvals0.shape,':',uwvals.shape,':',
1276#      vwvals.shape
1277
1278    if len(uwvals.shape) != 2 or len(vwvals.shape) != 2:
1279        print errormsg
1280        print '  ' + fname + ': wrong size of the wind fields! they must be ' +      \
1281          '2-dimensional!'
1282        print '    u-winds shape:',uwvals.shape,'dims:',of.variables[varNs[2]]
1283        print '    v-winds shape:',vwvals.shape,'dims:',of.variables[varNs[3]]
1284        print '      provide more values for their dimensions!!'
1285        quit(-1)
1286
1287    if len(lonvals0.shape) == 1:
1288        lonvals, latvals = np.meshgrid(lonvals0, latvals0)
1289    else:
1290        lonvals = lonvals0
1291        latvals = latvals0
1292
1293# Vecor values
1294    if vecvals.split(',')[0] == 'None':
1295        freqv = None
1296    else:
1297        freqv = vecvals.split(',')[0] 
1298    colorv = vecvals.split(',')[1]
1299    lengthv = vecvals.split(',')[2]
1300
1301# Vector labels
1302    windname = windlabels.split(',')[0]
1303    windunits = windlabels.split(',')[1]
1304
1305    drw.plot_barbs(lonvals, latvals, uwvals, vwvals, freqv, colorv, lengthv, 
1306      windname, windunits, mapvalues, gtit, kindfig, figuren)
1307
1308    return
1309 
1310def draw_topo_geogrid(ncfile, values):
1311    """ plotting geo_em.d[nn].nc topography from WPS files
1312    draw_topo_geogrid(ncfile, values)
1313      ncfile= geo_em.d[nn].nc file to use
1314      values= [minTopo],[maxTopo]:[lonlatL]:[title]:[graphic_kind]:[mapvalues]
1315        [min/max]Topo: minimum and maximum values of topography to draw
1316        lonlatL: limits of longitudes and latitudes [lonmin, latmin, lonmax, latmax] or None
1317        title: title of the graph ('!' for spaces)
1318        graphic_kind: kind of figure (jpg, pdf, png)
1319        mapvalues: map characteristics [proj],[res]
1320          see full documentation: http://matplotlib.org/basemap/
1321          [proj]: projection
1322            * 'cyl', cilindric
1323            * 'lcc', lambert conformal
1324          [res]: resolution:
1325            * 'c', crude
1326            * 'l', low
1327            * 'i', intermediate
1328            * 'h', high
1329            * 'f', full
1330    """
1331    fname = 'draw_topo_geogrid'
1332
1333    if values == 'h':
1334        print fname + '_____________________________________________________________'
1335        print draw_topo_geogrid.__doc__
1336        quit()
1337
1338    expectargs = '[minTopo],[maxTopo]:[lonlatL]:[title]:[graphic_kind]:[mapvalues]'
1339 
1340    drw.check_arguments(fname,values,expectargs,':')
1341
1342    mintopo = np.float(values.split(':')[0].split(',')[0])
1343    maxtopo = np.float(values.split(':')[0].split(',')[1])
1344
1345    lonlatLS = values.split(':')[1]
1346    lonlatLv = lonlatLS.split(',')[0]
1347
1348    if lonlatLv == 'None':
1349        lonlatL = None
1350    else:
1351        lonlatL = np.zeros((4), dtype=np.float)
1352        lonlatL[0] = np.float(lonlatLS.split(',')[0])
1353        lonlatL[1] = np.float(lonlatLS.split(',')[1])
1354        lonlatL[2] = np.float(lonlatLS.split(',')[2])
1355        lonlatL[3] = np.float(lonlatLS.split(',')[3])
1356
1357    grtit = values.split(':')[2].replace('!',' ')
1358    kindfig = values.split(':')[3]
1359    mapvalues = values.split(':')[4]
1360
1361    if not os.path.isfile(ncfile):
1362        print errormsg
1363        print '  ' + fname + ': domain file "' + ncfile + '" does not exist !!'
1364        quit(-1)   
1365
1366    objdomf = NetCDFFile(ncfile, 'r')
1367   
1368    objhgt = objdomf.variables['HGT_M']
1369    objlon = objdomf.variables['XLONG_M']
1370    objlat = objdomf.variables['XLAT_M']
1371
1372    topography = objhgt[0,:,:]
1373
1374    drw.plot_topo_geogrid(topography, objlon, objlat, mintopo, maxtopo, lonlatL,     \
1375      grtit, kindfig, mapvalues, True)
1376
1377    objdomf.close()
1378
1379    return
1380
1381def draw_topo_geogrid_boxes(ncfiles, values):
1382    """ plotting different geo_em.d[nn].nc topography from WPS files
1383    draw_topo_geogrid_boxes(ncfile, values)
1384      ncfiles= ',' list of geo_em.d[nn].nc files to use (fisrt as topographyc reference)
1385      values= [minTopo],[maxTopo]:[lonlatL]:[title]:[graphic_kind]:[mapvalues]:[labels]:[legloc]
1386        [min/max]Topo: minimum and maximum values of topography to draw
1387        lonlatL: limits of longitudes and latitudes [lonmin, latmin, lonmax, latmax] or None
1388        title: title of the graph ('!' for spaces)
1389        graphic_kind: kind of figure (jpg, pdf, png)
1390        mapvalues: map characteristics [proj],[res]
1391          see full documentation: http://matplotlib.org/basemap/
1392          [proj]: projection
1393            * 'cyl', cilindric
1394            * 'lcc', lambert conformal
1395          [res]: resolution:
1396            * 'c', crude
1397            * 'l', low
1398            * 'i', intermediate
1399            * 'h', high
1400            * 'f', full
1401        legloc= location of the legend (0, autmoatic)
1402          1: 'upper right', 2: 'upper left', 3: 'lower left', 4: 'lower right',
1403          5: 'right', 6: 'center left', 7: 'center right', 8: 'lower center',
1404          9: 'upper center', 10: 'center'
1405        labels= labels to write in the graph
1406    """
1407#    import matplotlib as mpl
1408#    mpl.use('Agg')
1409    import matplotlib.pyplot as plt
1410
1411    fname = 'draw_topo_geogrid_boxes'
1412
1413    if values == 'h':
1414        print fname + '_____________________________________________________________'
1415        print draw_topo_geogrid_boxes.__doc__
1416        quit()
1417
1418    mintopo = np.float(values.split(':')[0].split(',')[0])
1419    maxtopo = np.float(values.split(':')[0].split(',')[1])
1420
1421    lonlatLS = values.split(':')[1]
1422    lonlatLv = lonlatLS.split(',')[0]
1423
1424    if lonlatLv == 'None':
1425        lonlatL = None
1426    else:
1427        lonlatL = np.zeros((4), dtype=np.float)
1428        lonlatL[0] = np.float(lonlatLS.split(',')[0])
1429        lonlatL[1] = np.float(lonlatLS.split(',')[1])
1430        lonlatL[2] = np.float(lonlatLS.split(',')[2])
1431        lonlatL[3] = np.float(lonlatLS.split(',')[3])
1432
1433    grtit = values.split(':')[2].replace('!', ' ')
1434    kindfig = values.split(':')[3]
1435    mapvalues = values.split(':')[4]
1436    labels = values.split(':')[5]
1437    legloc = int(values.split(':')[6])
1438
1439    ncfile = ncfiles.split(',')[0]
1440    if not os.path.isfile(ncfile):
1441        print errormsg
1442        print '  ' + fname + ': domain file "' + ncfile + '" does not exist !!'
1443        quit(-1)   
1444
1445    objdomf = NetCDFFile(ncfile, 'r')
1446   
1447    objhgt = objdomf.variables['HGT_M']
1448    objlon0 = objdomf.variables['XLONG_M']
1449    objlat0 = objdomf.variables['XLAT_M']
1450
1451    topography = objhgt[0,:,:]
1452
1453    Nfiles = len(ncfiles.split(','))
1454    boxlabels = labels.split(',')
1455
1456    Xboxlines = []
1457    Yboxlines = []
1458
1459    for ifile in range(Nfiles):
1460        ncfile = ncfiles.split(',')[ifile]
1461#        print ifile, ncfile
1462        if not os.path.isfile(ncfile):
1463            print errormsg
1464            print '  ' + fname + ': domain file "' + ncfile + '" does not exist !!'
1465            quit(-1)   
1466
1467        objdomfi = NetCDFFile(ncfile, 'r')
1468   
1469        objlon = objdomfi.variables['XLONG_M']
1470        objlat = objdomfi.variables['XLAT_M']
1471
1472        dx = objlon.shape[2]
1473        dy = objlon.shape[1]
1474
1475        Xboxlines.append(objlon[0,0,:])
1476        Yboxlines.append(objlat[0,0,:])
1477        Xboxlines.append(objlon[0,dy-1,:])
1478        Yboxlines.append(objlat[0,dy-1,:])
1479        Xboxlines.append(objlon[0,:,0])
1480        Yboxlines.append(objlat[0,:,0])
1481        Xboxlines.append(objlon[0,:,dx-1])
1482        Yboxlines.append(objlat[0,:,dx-1])
1483
1484        objdomfi.close()
1485
1486    drw.plot_topo_geogrid_boxes(topography, Xboxlines, Yboxlines, boxlabels,         \
1487      objlon0, objlat0, mintopo, maxtopo, lonlatL, grtit, kindfig, mapvalues, legloc,\
1488      True)
1489
1490    objdomf.close()
1491
1492    return
1493
1494def movievalslice(origslice, dimmovien, framenum):
1495    """ Function to provide variable slice according to a geneation of a movie
1496    movievals(origslice, dimmovien, framenum)
1497      [origslice]= slice original as [dimname1]|[val1],[...,[dimnameN]|[valN]]
1498        ([val] = -1, full length)
1499      [dimmovien]= name of the dimension to produce the movie
1500      [framenum]= value of the frame to substitue in [origslice] as
1501        [dimmovien]|[framenum]
1502    >>> movievalslice('East_West|-1,North_South|-1,Time|2','Time',0)
1503    East_West|-1,North_South|-1,Time|0
1504    """
1505
1506    fname = 'movievalslice'
1507
1508    if origslice == 'h':
1509        print fname + '_____________________________________________________________'
1510        print movievalslice.__doc__
1511        quit()
1512   
1513    dims = origslice.split(',')
1514
1515    movieslice = ''
1516    idim = 0
1517
1518    for dimn in dims:
1519        dn = dimn.split('|')[0]
1520        if dn == dimmovien:
1521            movieslice = movieslice + dn + '|' + str(framenum)
1522        else:
1523            movieslice = movieslice + dimn
1524        if idim < len(dims)-1: movieslice = movieslice + ','
1525
1526        idim = idim + 1
1527
1528    return movieslice
1529
1530class Capturing(list):
1531    """ Class to capture function output as a list
1532    from: http://stackoverflow.com/questions/16571150/how-to-capture-stdout-output-from-a-python-function-call
1533    """
1534#    from cStringIO import StringIO
1535
1536    def __enter__(self):
1537        self._stdout = sys.stdout
1538        sys.stdout = self._stringio = StringIO()
1539        return self
1540    def __exit__(self, *args):
1541        self.extend(self._stringio.getvalue().splitlines())
1542        sys.stdout = self._stdout
1543
1544def create_movie(netcdfile, values, variable):
1545    """ Function to create a movie assuming ImageMagick installed!
1546      values= [graph]#[movie_dimension]#[graph_values]
1547        [graph]: which graphic
1548        [movie_dimension]: [dimnmovie]@[dimvmovie]@[moviedelay]@[interval]
1549          [dimnmovie]; name of the dimension from which make the movie
1550          [dimvmovie]; name of the variable with the values of the dimension
1551          [moviedelay]; delay between frames
1552          [interval]; [beg]@[end]@[freq] or -1 (all)
1553        [graph_values]: values to generate the graphic
1554      netcdfile= netCDF file
1555      variable= variable to use (when applicable)
1556    """ 
1557    fname = 'create_movie'
1558
1559    if values == 'h':
1560        print fname + '_____________________________________________________________'
1561        print create_movie.__doc__
1562        quit()
1563
1564    graph = values.split('#')[0]
1565    movie_dim = values.split('#')[1]
1566    graph_vals = values.split('#')[2]
1567
1568    ncobj = NetCDFFile(netcdfile, 'r')
1569
1570# Movie dimension
1571##
1572    dimnmovie = movie_dim.split('@')[0]
1573    dimvmovie = movie_dim.split('@')[1]
1574    moviedelay = movie_dim.split('@')[2]
1575    moviebeg = int(movie_dim.split('@')[3])
1576
1577    if not drw.searchInlist(ncobj.dimensions.keys(),dimnmovie):
1578        print errormsg
1579        print '  ' + fname + ": file '" + netcdfile + "' has not dimension '" +      \
1580          dimnmovie + "' !!!"
1581        quit(-1)
1582
1583    objdmovie = ncobj.dimensions[dimnmovie]
1584    dmovie = len(objdmovie)
1585    if moviebeg != -1:
1586        moviend = int(movie_dim.split('@')[4])
1587        moviefreq = int(movie_dim.split('@')[5])
1588    else:
1589        moviebeg = 0
1590        moviend = dmovie
1591        moviefreq = 1
1592
1593    if dimvmovie == 'WRFTimes':
1594        objvdmovie = ncobj.variables['Times']
1595        vdmovieunits = ''
1596        valsdmovie = []
1597        for it in range(objvdmovie.shape[0]):
1598            valsdmovie.append(drw.datetimeStr_conversion(objvdmovie[it,:],           \
1599              'WRFdatetime', 'Y/m/d H-M-S'))
1600    elif dimvmovie == 'CFtime':
1601        objvdmovie = ncobj.variables['time']
1602        vdmovieunits = ''
1603        print objvdmovie.units
1604        valsdmovie0 = drw.netCDFdatetime_realdatetime(objvdmovie.units, 'standard',  \
1605          objvdmovie[:])
1606        valsdmovie = []
1607        for it in range(objvdmovie.shape[0]):
1608            valsdmovie.append(drw.datetimeStr_conversion(valsdmovie0[it,:],          \
1609              'matYmdHMS', 'Y/m/d H-M-S'))
1610    else:
1611        if  not drw.searchInlist(ncobj.variables.keys(),dimvmovie):
1612            print errormsg
1613            print '  ' + fname + ": file '" + netcdfile + "' has not variable '" +   \
1614              dimvmovie + "' !!!"
1615            quit(-1)
1616        vdmovieunits = objvdmovie.getncattr('units')
1617        objvdmovie = ncobj.variables[dimvmovie]
1618        if len(objvdmovie.shape) == 1:
1619            vasldmovie = objvdmovie[:]
1620        else:
1621            print errormsg
1622            print '  ' + fname + ': shape', objvdmovie.shape, 'of variable with ' +  \
1623              'dimension movie values not ready!!!'
1624            quit(-1)
1625
1626    ncobj.close()
1627    os.system('rm frame_*.png > /dev/null')
1628
1629# graphic
1630##
1631    if graph == 'draw_2D_shad':
1632        graphvals = graph_vals.split(':')
1633
1634        for iframe in range(moviebeg,moviend,moviefreq):
1635            iframeS = str(iframe).zfill(4)
1636
1637            drw.percendone((iframe-moviebeg)/moviefreq,(moviend-moviebeg)/moviefreq, \
1638              5, 'frames')
1639            titgraph = dimnmovie + '|=|' + str(valsdmovie[iframe]) + '|' +           \
1640              vdmovieunits
1641
1642            graphvals[1] = movievalslice(graphvals[1],dimnmovie,iframe)
1643            graphvals[6] = titgraph
1644            graphvals[7] = 'png'
1645
1646            graphv = drw.numVector_String(graphvals, ":")
1647
1648            with Capturing() as output:
1649                draw_2D_shad(netcdfile, graphv, variable)
1650
1651            os.system('mv 2Dfields_shadow.png frame_' + iframeS + '.png')
1652    else:
1653        print errormsg
1654        print '  ' + fname + ": graphic '" +  graph + "' not defined !!!"
1655        quit(-1)
1656
1657    os.system('convert -delay ' + moviedelay + ' -loop 0 frame_*.png create_movie.gif')
1658
1659    print "Succesfuly creation of movie file 'create_movie.gif' !!!"
1660
1661    return
1662
1663def draw_lines(ncfilens, values, varname):
1664    """ Function to draw different lines at the same time from different files
1665    draw_lines(ncfilens, values, varname):
1666      ncfilens= [filen] ',' separated list of netCDF files
1667      values= [dimvname]:[valuesaxis]:[dimtit]:[leglabels]:[vtit]:[title]:[locleg]:[colns]:[lines]
1668       [points]:[lwdths]:[psizes]:[freqv]:[figname][graphk]
1669        [dimvname]: ',' list of names of the variable with he values of the common dimension
1670        [valuesaxis]: which axis will be used for the values ('x', or 'y')
1671        [dimtit]: title for the common dimension
1672        [leglabels]: ',' separated list of names for the legend
1673        [vartit]: name of the variable in the graph
1674        [title]: title of the plot ('|' for spaces)
1675        [locleg]: location of the legend (0, autmoatic)
1676          1: 'upper right', 2: 'upper left', 3: 'lower left', 4: 'lower right',
1677          5: 'right', 6: 'center left', 7: 'center right', 8: 'lower center',
1678          9: 'upper center', 10: 'center'
1679        [colns]= ',' list of color names ('None' for automatic, single value for all the same)
1680        [lines]= ',' list of style of lines ('None' for automatic, single value for all the same)
1681        [points]= ',' list of style of points ('None' for automatic, single value for all the same)
1682        [lwdths]= ',' list of withs of lines ('None' for automatic, single value for all the same)
1683        [psizes]= ',' list of size of points ('None' for automatic, single value for all the same)
1684        [freqv]= frequency of values ('all' for all values)
1685        [figname]= name of the figure
1686        [graphk]: kind of the graphic
1687      varname= variable to plot
1688      values= 'XLAT:x:latitude:32x32:$wss^{*}$:wss Taylor's turbulence term:pdf'
1689    """
1690
1691    fname = 'draw_lines'
1692
1693    if values == 'h':
1694        print fname + '_____________________________________________________________'
1695        print draw_lines.__doc__
1696        quit()
1697
1698    expectargs = '[dimvname]:[valuesaxis]:[dimtit]:[leglabels]:[vtit]:[title]:' +    \
1699      '[locleg]:[colns]:[lines]:[points]:[lwdths]:[psizes]:[freqv]:[graphk]'
1700    drw.check_arguments(fname,values,expectargs,':')
1701
1702    ncfiles = ncfilens.split(',')
1703    dimvnames = values.split(':')[0]
1704    valuesaxis = values.split(':')[1]
1705    dimtit = values.split(':')[2]
1706    leglabels = values.split(':')[3].replace('_','\_')
1707    vartit = values.split(':')[4]
1708    title = values.split(':')[5].replace('|',' ')
1709    locleg = values.split(':')[6]
1710    colns = str_list(values.split(':')[7], ',')
1711    lines = str_list(values.split(':')[8], ',')
1712    points = str_list(values.split(':')[9], ',')
1713    lwdths = str_list(values.split(':')[10], ',')
1714    psizes = str_list(values.split(':')[11], ',')
1715    freqv0 = values.split(':')[12]
1716    figname = values.split(':')[13]
1717    graphk = values.split(':')[14]
1718
1719    Nfiles = len(ncfiles)
1720
1721# Getting trajectotries
1722##
1723
1724    varvalues = []
1725    dimvalues = []
1726
1727    print '  ' + fname
1728    ifn = 0
1729    for ifile in ncfiles:
1730        filen = ifile.split('@')[0]
1731
1732        print '    filen:',filen
1733
1734        if not os.path.isfile(filen):
1735            print errormsg
1736            print '  ' + fname + ": netCDF file '" + filen + "' does not exist !!"
1737            quit(-1)
1738
1739        objfile = NetCDFFile(filen, 'r')
1740
1741        if dimvnames.find(',') != -1:
1742            dimvname = dimvnames.split(',')
1743        else:
1744            dimvname = [dimvnames]
1745   
1746        found = False
1747        for dvn in dimvname:
1748            if objfile.variables.has_key(dvn):
1749                found = True
1750                break
1751   
1752        if not found:
1753            print errormsg
1754            print '  ' + fname + ": netCDF file '" + filen +                         \
1755              "' does not have variables '" + dimvnames + "' !!"
1756            quit(-1)
1757
1758        if not objfile.variables.has_key(varname):
1759            print errormsg
1760            print '  ' + fname + ": netCDF file '" + filen +                         \
1761              "' does not have variable '" + varname + "' !!"
1762            quit(-1)
1763
1764        vvobj = objfile.variables[varname]
1765        if len(vvobj.shape) != 1:
1766            print errormsg
1767            print '  ' + fname + ': wrong shape:',vvobj.shape," of variable '" +     \
1768              varname +  "' !!"
1769            quit(-1)
1770
1771        for dimvn in dimvname:
1772            if drw.searchInlist(objfile.variables, dimvn):
1773                vdobj = objfile.variables[dimvn]
1774                if len(vdobj.shape) != 1:
1775                    print errormsg
1776                    print '  ' + fname + ': wrong shape:',vdobj.shape,               \
1777                      " of variable '" + dimvn +  "' !!"
1778                    quit(-1)
1779                break
1780
1781        varvalues.append(vvobj[:])
1782        dimvalues.append(vdobj[:])
1783
1784        if ifn == 0:
1785            varunits = vvobj.units
1786
1787        objfile.close()
1788
1789        ifn = ifn + 1
1790
1791    if freqv0 == 'all':
1792        freqv = None
1793    else:
1794        freqv = int(freqv0)
1795
1796    drw.plot_lines(dimvalues, varvalues, valuesaxis, dimtit, leglabels.split(','),   \
1797      vartit, varunits, title, locleg, colns, lines, points, lwdths, psizes, freqv,  \
1798      figname,graphk)
1799
1800    return
1801
1802def draw_lines_time(ncfilens, values, varname0):
1803    """ Function to draw different lines at the same time from different files with times
1804    draw_lines_time(ncfilens, values, varname):
1805      ncfilens= [filen] ',' separated list of netCDF files
1806      values= [dimvname];[valuesaxis];[dimtit];[leglabels];[vtit];[title];[rangevals];[timevals];
1807        [legvals];[graphk];[collines];[points];[linewidths];[pointsizes];[pointfreq];[period]
1808        [dimvname]: ',' list of names of the variables with he values of the common dimension
1809        [valuesaxis]: which axis will be used for the values ('x', or 'y')
1810        [dimtit]: title for the common dimension ('|' for spaces)
1811        [leglabels]: ',' separated list of names for the legend ('None', no legend '!' for spaces)
1812        [vartit]: name of the variable in the graph
1813        [title]: title of the plot ('|' for spaces)
1814        [rangevals]: Range of the axis with the values ('None' for 'auto','auto')
1815          [vmin],[vmax]: minimum and maximum values where [vmNN] can also be:
1816            'auto': the computed minimumm or maximum of the values 
1817        [timevals]: [timen]|[units]|[kind]|[tfmt] time labels characteristics
1818           [timen]; name of the time variable
1819           [units]; units string according to CF conventions ([tunits] since
1820             [YYYY]-[MM]-[DD] [[HH]:[MI]:[SS]], '!' for spaces)
1821           [kind]; kind of output
1822             'Nval': according to a given number of values as 'Nval',[Nval]
1823             'exct': according to an exact time unit as 'exct',[tunit];
1824               tunit= [Nunits],[tu]; [tu]= 'c': centuries, 'y': year, 'm': month,
1825                'w': week, 'd': day, 'h': hour, 'i': minute, 's': second,
1826                'l': milisecond
1827           [tfmt]; desired format
1828        [legvals]=[locleg]|[fontsize]:
1829          [locleg]: location of the legend (0, autmoatic)
1830            1: 'upper right', 2: 'upper left', 3: 'lower left', 4: 'lower right',
1831            5: 'right', 6: 'center left', 7: 'center right', 8: 'lower center',
1832            9: 'upper center', 10: 'center'
1833          [fontsize]: font size for the legend (auto for 12)
1834        [graphk]: kind of the graphic
1835        [lines]: ',' list of type of lines, None for automatic, single value all the same
1836        [collines]: ',' list of colors for the lines, None for automatic, single
1837          value all the same
1838        [points]: ',' list of type of points for the lines, None for automatic, single
1839          value all the same
1840        [linewidths]: ',' list of widths for the lines, None for automatic, single
1841          value all the same
1842        [pointsizes]: ',' list of widths for the lines, None for automatic, single
1843          value all the same
1844        [pointfreq]: frequency of point plotting, 'all' for all time steps
1845        [period]: which period to plot
1846          '-1': all period
1847          [beg],[end]: beginning and end of the period in reference time-units of first file
1848      varname0= ',' list of variable names to plot (assuming only 1 variable per file)
1849      values= 'time;y;time ([DD]${[HH]}$);32x32;$wss^{*}$;wss Taylor's turbulence term;time|hours!since!1949-12-01_00:00:00;exct,12,h|%d$^{%H}$;2;pdf'
1850    """
1851
1852    fname = 'draw_lines_time'
1853
1854    if values == 'h':
1855        print fname + '_____________________________________________________________'
1856        print draw_lines_time.__doc__
1857        quit()
1858
1859    expectargs = '[dimvname];[valuesaxis];[dimtit];[leglabels];[vtit];[title];'
1860    expectargs = expectargs + '[rangevals];[timevals];[legvals];[graphk];[lines];'
1861    expectargs = expectargs + '[collines];[points];[linewidths];[pointsizes];'
1862    expectargs = expectargs + '[pointfreq];[period]'
1863    drw.check_arguments(fname,values,expectargs,';')
1864
1865    ncfiles = ncfilens.split(',')
1866    dimvname0 = values.split(';')[0]
1867    valuesaxis = values.split(';')[1]
1868    dimtit = values.split(';')[2].replace('|',' ')
1869    leglabels = values.split(';')[3].replace('_','\_').replace('!',' ')
1870    vartit = values.split(';')[4]
1871    title = values.split(';')[5].replace('|',' ')
1872    rangevals = values.split(';')[6]
1873    timevals = values.split(';')[7]
1874    legvalues = values.split(';')[8]
1875    graphk = values.split(';')[9]
1876    lines0 = values.split(';')[10]
1877    collines0 = values.split(';')[11]
1878    points0 = values.split(';')[12]
1879    linewidths0 = values.split(';')[13]
1880    pointsizes0 = values.split(';')[14]
1881    pointfreq0 = values.split(';')[15]
1882    period = values.split(';')[16]
1883
1884    Nfiles = len(ncfiles)
1885
1886# Multiple variable-dimension names?
1887    if dimvname0.find(',') != -1:
1888        dimvname = dimvname0.split(',')
1889    else:
1890        dimvname = [dimvname0]
1891
1892# Multiple variables?
1893    if varname0.find(',') != -1:
1894        varname = varname0.split(',')
1895    else:
1896        varname = [varname0]
1897
1898# Multiple lines types?
1899    if lines0.find(',') != -1:
1900        lines = lines0.split(',')
1901    elif lines0 == 'None':
1902        lines = None
1903    else:
1904        lines = []
1905        for il in range(Nfiles):
1906            lines.append(lines0)
1907
1908# Multiple color names?
1909    if collines0.find(',') != -1:
1910        collines = collines0.split(',')
1911    elif collines0 == 'None':
1912        collines = None
1913    else:
1914        collines = []
1915        for ip in range(Nfiles):
1916            collines.append(collines0)
1917
1918# Multiple point types?
1919    if points0.find(',') != -1:
1920        if len(points0) == 1:
1921            points = []
1922            for ip in range(Nfiles):
1923                points.append(points0)
1924        else:
1925            points = points0.split(',')
1926    elif points0 == 'None':
1927        points = None
1928    else:
1929        points = []
1930        for ip in range(Nfiles):
1931            points.append(points0)
1932
1933# Multiple line sizes?
1934    if linewidths0.find(',') != -1:
1935        linewidths = []
1936        Nlines = len(linewidths0.split(','))
1937        for il in range(Nlines):
1938          linewidths.append(np.float(linewidths0.split(',')[il]))
1939    elif linewidths0 == 'None':
1940        linewidths = None
1941    else:
1942        linewidths = [np.float(linewidths0)]
1943
1944# Multiple point sizes?
1945    if pointsizes0.find(',') != -1:
1946        pointsizes = []
1947        Npts = len(pointsizes0.split(','))
1948        for ip in Npts:
1949          pointsizes.append(np.float(pointsizes0.split(',')[ip]))
1950    elif pointsizes0 == 'None':
1951        pointsizes = None
1952    else:
1953        pointsizes = [np.float(pointsizes0)]
1954
1955    timename = timevals.split('|')[0]
1956    timeunit = timevals.split('|')[1].replace('!',' ')
1957    timekind = timevals.split('|')[2]
1958    timefmt = timevals.split('|')[3]
1959
1960    if rangevals == 'None':
1961        valmin = 'auto'
1962        valmax = 'auto'
1963    else:
1964        valmin = rangevals.split(',')[0]
1965        valmax = rangevals.split(',')[1]
1966        if valmin != 'auto': valmin = np.float(valmin)
1967        if valmax != 'auto': valmax = np.float(valmax)
1968
1969    locleg = int(legvalues.split('|')[0])
1970    if legvalues.split('|')[1] == 'auto':
1971        legfontsize = 12
1972    else:
1973        legfontsize = int(legvalues.split('|')[1])
1974
1975# Getting values
1976##
1977    varvalues = []
1978    dimvalues = []
1979    timvalues = []
1980    timvals0 = timvalues
1981
1982    print '  ' + fname
1983    ifn = 0
1984    mintval = 1.e20
1985    maxtval = -1.e20
1986
1987    for ifile in ncfiles:
1988        filen = ifile.split('@')[0]
1989
1990        print '    filen:',filen
1991
1992        if not os.path.isfile(filen):
1993            print errormsg
1994            print '  ' + fname + ": netCDF file '" + filen + "' does not exist !!"
1995            quit(-1)
1996
1997        objfile = NetCDFFile(filen, 'r')
1998
1999        founddvar = False
2000        for dvar in dimvname:
2001            if objfile.variables.has_key(dvar):
2002                founddvar = True
2003                vdobj = objfile.variables[dvar]
2004                if len(vdobj.shape) != 1:
2005                    print errormsg
2006                    print '  ' + fname + ': wrong shape:',vdobj.shape," of " +       \
2007                      "variable '" + dvar +  "' !!"
2008                    quit(-1)
2009                break
2010        if not founddvar:
2011            print errormsg
2012            print '  ' + fname + ": netCDF file '" + filen +                         \
2013            "' has any variable '", dimvname, "' !!"
2014            quit(-1)
2015
2016        foundvar = False
2017        for var in varname:
2018            if objfile.variables.has_key(var):
2019                foundvar = True
2020                vvobj = objfile.variables[var]
2021                if len(vvobj.shape) != 1:
2022                    print errormsg
2023                    print '  ' + fname + ': wrong shape:',vvobj.shape," of " +       \
2024                      "variable '" + var +  "' !!"
2025                    quit(-1)
2026
2027                break
2028        if not foundvar:
2029            print errormsg
2030            print '  ' + fname + ": netCDF file '" + filen +                         \
2031              "' has any variable '", varname, "' !!"
2032            quit(-1)
2033        if vdobj.units.find('month') != 1:
2034            print warnmsg
2035            print '  ' + fname + ": tranforming time units from 'months' to 'days'!!"
2036            timevals0, tunits0 = gen.CFmonthU_daysU(vdobj[:], vdobj.units) 
2037        else:
2038            timevals0 = vdobj[:]
2039            tunits0 = str(vdobj.units) 
2040
2041# Getting period
2042        if ifn > 0: 
2043# Referring all times to the same reference time!
2044            reftvals = drw.coincident_CFtimes(timevals0, timeunit, tunits0)
2045        else:
2046            reftvals = timevals0
2047
2048        dimt = len(vdobj[:])
2049
2050        if period == '-1':
2051            varvalues.append(vvobj[:])
2052            dimvalues.append(reftvals)
2053            mindvals = np.min(reftvals)
2054            maxdvals = np.max(reftvals)
2055        else:
2056            ibeg=-1
2057            iend=-1
2058            tbeg = np.float(period.split(',')[0])
2059            tend = np.float(period.split(',')[1])
2060
2061            for it in range(dimt-1):
2062                if reftvals[it] <= tbeg and reftvals[it+1] > tbeg: ibeg = it
2063                if reftvals[it] <= tend and reftvals[it+1] > tend: iend = it + 1
2064                if ibeg != -1 and iend != -1: break
2065
2066            if ibeg == -1 and iend == -1:
2067                print warnmsg
2068                print '  ' + fname + ': Period:',tbeg,',',tend,'not found!!'
2069                print '    ibeg:',ibeg,'iend:',iend
2070                print '    period in file:',np.min(reftvals), np.max(reftvals)
2071                print '    getting all the period in file !!!'
2072                ibeg = 0
2073                iend = dimt
2074            elif iend == -1:
2075                iend = dimt
2076                print warnmsg
2077                print '  ' + fname + ': end of Period:',tbeg,',',tend,'not found!!'
2078                print '    getting last available time instead'
2079                print '    ibeg:',ibeg,'iend:',iend
2080                print '    period in file:',np.min(reftvals), np.max(reftvals)
2081            elif ibeg == -1:
2082                ibeg = 0
2083                print warnmsg
2084                print '  ' + fname + ': beginning of Period:',tbeg,',',tend,         \
2085                  'not found!!'
2086                print '    getting first available time instead'
2087                print '    ibeg:',ibeg,'iend:',iend
2088                print '    period in file:',np.min(reftvals), np.max(reftvals)
2089
2090            varvalues.append(vvobj[ibeg:iend])
2091            dimvalues.append(reftvals[ibeg:iend])
2092            mindvals = np.min(reftvals[ibeg:iend])
2093            maxdvals = np.max(reftvals[ibeg:iend])
2094
2095            dimt = iend - ibeg
2096
2097        if mindvals < mintval: mintval = mindvals
2098        if maxdvals > maxtval: maxtval = maxdvals
2099        print '  ' + fname + ": file '" + filen + "' period:", mindvals, '->', maxdvals
2100
2101        if ifn == 0:
2102            varunits = drw.units_lunits(vvobj.units)
2103
2104        objfile.close()
2105
2106        ifn = ifn + 1
2107
2108# Times
2109
2110    dtvals = (maxtval - mintval)/dimt
2111#    dti = mintval-dtvals/2.
2112#    dte = maxtval+dtvals/2.
2113    dti = mintval
2114    dte = maxtval
2115    tvals = np.arange(dti, dte, dtvals)
2116
2117    dtiS = drw.datetimeStr_conversion(str(dti) + ',' + timeunit, 'cfTime',           \
2118      'Y/m/d H-M-S')
2119    dteS = drw.datetimeStr_conversion(str(dte) + ',' + timeunit, 'cfTime',           \
2120      'Y/m/d H-M-S')
2121
2122    print '  ' + fname + ': plotting from: ' + dtiS + ' to ' + dteS
2123
2124    timepos, timelabels = drw.CFtimes_plot(tvals, timeunit, timekind, timefmt)
2125
2126#    print 'Lluis min/max tval +/- dtval/2:', mintval-dtvals/2., maxtval+dtvals/2.,'dt:', len(tvals)
2127#    for it in range(len(timepos)):
2128#        print timepos[it], timelabels[it]
2129
2130    if leglabels != 'None':
2131        legvals = leglabels.split(',')
2132    else:
2133        legvals = None
2134
2135    if pointfreq0 == 'all':
2136        pointfreq = None
2137    else:
2138        pointfreq = int(pointfreq0)
2139
2140    drw.plot_lines_time(dimvalues, varvalues, valuesaxis, dimtit, legvals, vartit,   \
2141      varunits, timepos, timelabels, title, locleg, legfontsize, graphk, valmin,     \
2142      valmax, lines, collines, points, linewidths, pointsizes, pointfreq)
2143
2144    return
2145
2146def draw_Neighbourghood_evol(filen, values, variable):
2147    """ Function to draw the temporal evolution of a neighbourghood around a point
2148    draw_Neighbourghood_evol(filen, values, variable)
2149      filen= netCDF file name
2150      values= [gvarname]:[dimsval]:[neigdims]:[Nneig]:[Ncol]:[timetits]:[tkinds]:
2151       [timefmts]:[gtitle]:[shadxtrms]:[cbar]:[gkind]:[ofile]
2152        [dimsval]: [dimn1]|[val1]|[dimv1],...,[dimnN]|[valN]|[dimvN] dimension names, values to get
2153          (-1, for all; no name/value pair given full length) and variable with values of the dimension
2154          NOTE: when dimsval[X,Y] == neigdims[X,Y], valX,valY --> valX,valY-Nneig/2, valX,valY+Nneig/2
2155        [neigdims]: [dimnX],[dimnY] dimensions mnames along which the neigbourghood should be defined
2156        [Nneig]: Number of grid points of the full side of the box (odd value)
2157        [Ncol]: Number of columns ('auto': square final plot)
2158        [gvarname]: name of the variable to appear in the graph
2159        [timetits]: [titX],[titY] titles of the axes ('|' for spaces)
2160        [tkinds]: [tkindX]|[tkindY] kinds of time to appear in the graph
2161          'Nval': according to a given number of values as 'Nval',[Nval]
2162          'exct': according to an exact time unit as 'exct',[tunit];
2163            tunit= [Nunits],[tu]; [tu]= 'c': centuries, 'y': year, 'm': month,
2164              'w': week, 'd': day, 'h': hour, 'i': minute, 's': second,
2165              'l': milisecond
2166        [timefmts]: [tfmtX],[tfmtY] format of the time labels
2167        [gtitle]: title of the graphic ('|' for spaces)
2168        [shadxtrms]: Extremes for the shading
2169        [cbar]: colorbar to use
2170        [gkind]: kind of graphical output
2171        [ofile]: True/False whether the netcdf with data should be created or not
2172      variable= name of the variable
2173      values = 'q:Time|-1|Times,bottom_top|6|ZNU,south_north|3|XLAT,west_east|26|XLONG:south_north,west_east:5:auto:time|($[DD]^{[HH]}$),time|($[DD]^{[HH]}$):exct,2,h|exct,1,d:$%d^{%H}$,$%d^{%H}$:5|pts|neighbourghood|temporal|evolution:0.0,0.004:BuPu:pdf:True'
2174    """ 
2175
2176    fname = 'draw_Neighbourghood_evol'
2177
2178    if values == 'h':
2179        print fname + '_____________________________________________________________'
2180        print draw_Neighbourghood_evol.__doc__
2181        quit()
2182
2183    expectargs = '[gvarname]:[dimsval]:[neigdims]:[Nneig]:[Ncol]:' +                 \
2184      '[timetits]:[tkinds]:[timefmts]:[gtitle]:[shadxtrms]:[cbar]:[gkind]:[ofile]'
2185 
2186    drw.check_arguments(fname,values,expectargs,':')
2187
2188    gvarname = values.split(':')[0]
2189    dimsval = values.split(':')[1].split(',')
2190    neigdims = values.split(':')[2].split(',')
2191    Nneig = int(values.split(':')[3])
2192    Ncol0 = values.split(':')[4]
2193    timetits = values.split(':')[5].split(',')
2194    timekinds = values.split(':')[6].split('|')
2195    timefmts = values.split(':')[7].split(',')
2196    gtitle = values.split(':')[8].replace('|',' ')
2197    shadxtrms = values.split(':')[9].split(',')
2198    cbar = values.split(':')[10]
2199    gkind = values.split(':')[11]
2200    ofile = values.split(':')[12]
2201
2202    if Ncol0 != 'auto': 
2203        Ncol = int(Ncol0)
2204    else:
2205        Ncol = Ncol0
2206
2207    timetits[0] = timetits[0].replace('|',' ')
2208    timetits[1] = timetits[1].replace('|',' ')
2209
2210    if np.mod(Nneig,2) == 0:
2211        print errormsg
2212        print '  ' + fname + ": an odd value for 'Nneig':",Nneig,'is required !!!'
2213        quit(-1)
2214
2215    Nneig2 = int(Nneig/2)
2216
2217# Values to slice the variable
2218    dimvslice = {}
2219    dimvvalues = {}
2220    for dimvs in dimsval:
2221        dimn = dimvs.split('|')[0]
2222        dimv = int(dimvs.split('|')[1])
2223        dimnv = dimvs.split('|')[2]
2224
2225        dimvvalues[dimn] = dimnv
2226        dimvslice[dimn] = dimv
2227
2228    ncobj = NetCDFFile(filen, 'r')
2229
2230    varobj = ncobj.variables[variable]
2231
2232    slicevar = []
2233    newdimn = []
2234    newdimsvar = {}
2235
2236    for dimn in varobj.dimensions:
2237        if not drw.searchInlist(dimvslice.keys(), dimn):
2238            dimsize = len(ncobj.dimensions[dimn])
2239            slicevar.append(slice(0, dimsize+1))
2240            newdimn.append(dimn)
2241            newdimsvar[dimn] = dimseize
2242
2243        for dimslicen in dimvslice.keys():
2244            if dimn == dimslicen:
2245                if dimvslice[dimn] != -1:
2246                    if drw.searchInlist(neigdims, dimn):
2247                        slicevar.append(slice(dimvslice[dimn]-Nneig2,                \
2248                          dimvslice[dimn]+Nneig2+1))
2249                        newdimn.append(dimn)
2250                        newdimsvar[dimn] = Nneig
2251                        break
2252                    else:
2253                        slicevar.append(slice(dimvslice[dimn], dimvslice[dimn]+1))
2254                        break
2255                else:
2256                    dimsize = len(ncobj.dimensions[dimn])
2257                    slicevar.append(slice(0, dimsize+1))
2258                    newdimn.append(dimn)
2259                    newdimsvar[dimn] = dimsize
2260                    break
2261 
2262    varv = varobj[tuple(slicevar)]
2263
2264    if len(newdimn) != 3:
2265        print errormsg
2266        print '  ' + fname + ': sliced variable with shape=', varv.shape,            \
2267          ' must have three dimensions',len(varv.shape),'given !!'
2268        quit(-1)
2269
2270    newdims = []
2271    for nwdims in newdimn:
2272        newdims.append(newdimsvar[nwdims])
2273
2274# The dimension which is not in the neighbourhood dimensions must be time!
2275    for dim1 in newdimn:
2276        if not drw.searchInlist(neigdims, dim1):
2277            dimt = newdimsvar[dim1]
2278            dimtime = dim1
2279
2280    if Ncol == 'auto':
2281        dimtsqx = int(np.sqrt(dimt)) + 1
2282        dimtsqy = int(np.sqrt(dimt)) + 1
2283    else:
2284        dimtsqx = int(Ncol)
2285        dimtsqy = dimt/dimtsqx + 1
2286
2287    neighbourghood = np.ones((dimtsqy*Nneig,dimtsqx*Nneig), dtype=np.float)*fillValue
2288
2289    for it in range(dimt):
2290        ity = int(it/dimtsqx)
2291        itx = it-ity*dimtsqx
2292
2293        itty = (dimtsqy - ity - 1)*Nneig + Nneig2
2294        ittx = itx*Nneig + Nneig2
2295
2296        neighbourghood[itty-Nneig2:itty+Nneig2+1,ittx-Nneig2:ittx+Nneig2+1]=         \
2297          varv[it,::-1,:]
2298
2299    variablevals = drw.variables_values(variable)
2300    if drw.searchInlist(varobj.ncattrs(), 'units'):
2301        vunits = varobj.units
2302    else:
2303        vunits = variablevals[5]
2304
2305# Time values at the X/Y axes
2306    if ncobj.variables[dimvvalues[dimtime]].dtype == '|S1':
2307        print '    ' + fname + ': WRF time variable!'
2308        refdate = '19491201000000'
2309        tunitsval = 'hours'
2310        dimtvalues = np.zeros((dimt), dtype=np.float)
2311        tvals = ncobj.variables[dimvvalues[dimtime]]
2312        yrref=refdate[0:4]
2313        monref=refdate[4:6]
2314        dayref=refdate[6:8]
2315        horref=refdate[8:10]
2316        minref=refdate[10:12]
2317        secref=refdate[12:14]
2318
2319        refdateS = yrref + '/' + monref + '/' + dayref + '_' + horref + ':' +        \
2320          minref + ':' + secref
2321        tunits = tunitsval + ' since ' + refdateS
2322        for it in range(dimt):
2323            wrfdates = drw.datetimeStr_conversion(tvals[it,:],'WRFdatetime', 'matYmdHMS')
2324            dimtvalues[it] = drw.realdatetime1_CFcompilant(wrfdates, refdate, tunitsval)
2325    else:
2326        dimtvalues = ncobj.variables[dimvvalues[dimtime]][:]
2327        tunits = ncobj.variables[newdimsvar[dimtime]].units
2328
2329    dimxv = dimtvalues[0:dimtsqx]
2330    dimyv = dimtvalues[0:dimt:dimtsqx]
2331
2332    dimn = ['time','time']
2333
2334    if ofile == 'True':
2335        ofilen = 'Neighbourghood_evol.nc'
2336        newnc = NetCDFFile(ofilen, 'w')
2337# Dimensions
2338        newdim = newnc.createDimension('time',None)
2339        newdim = newnc.createDimension('y',dimtsqy*Nneig)
2340        newdim = newnc.createDimension('x',dimtsqx*Nneig)
2341# Dimension values
2342        newvar = newnc.createVariable('time','f8',('time'))
2343        newvar[:] = np.arange(dimt)
2344        newattr = drw.basicvardef(newvar, 'time','time',tunits)
2345# Neighbourhghood variable
2346        newvar = newnc.createVariable(variable + 'neigevol', 'f4', ('y','x'),        \
2347          fill_value=fillValue)
2348        newvar[:] = neighbourghood
2349
2350        newnc.sync()
2351        newnc.close()
2352        print fname + ": Successfull generation of file '" + ofilen + "' !!"
2353
2354# Time ticks
2355    timeposX, timelabelsX = drw.CFtimes_plot(dimxv, tunits, timekinds[0], timefmts[0])
2356    timeposY, timelabelsY = drw.CFtimes_plot(dimyv, tunits, timekinds[1], timefmts[1])
2357
2358    timepos = [timeposX[0:len(timeposX)], timeposY[len(timeposY):0:-1]]
2359    timelabels = [timelabelsX[0:len(timeposX)], timelabelsY[0:len(timeposY)]]
2360
2361    for i in range(2):
2362        if shadxtrms[i][0:1] != 'S':
2363            shadxtrms[i] = np.float(shadxtrms[i])
2364
2365    drw.plot_Neighbourghood_evol(neighbourghood, dimxv, dimyv, gvarname, timetits,   \
2366      timepos, timelabels, cbar, Nneig, shadxtrms, vunits, gtitle, gkind, True)
2367
2368def draw_points(filen, values):
2369    """ Function to plot a series of points
2370      [values]= [ptasciifile]:[gtit]:[mapvalues]:[kindfigure]:[pointcolor]:[pointlabels]:
2371        [locleg]:[figureko]:[figuren]
2372        [ptasciifile]:[file],[comchar],[collon],[collat],[lab]
2373          [file]: column ASCII file with the location of the points
2374          [comchar]: '|' list of characters for commentaries
2375          [collon]: number of column with the longitude of the points
2376          [collat]: number of column with the latitude of the points
2377          [collab]: number of column with the labels of the points ('None', and will get
2378            the values from the [pointlabels] variable
2379        [gtit]: title of the figure ('|' for spaces)
2380        [mapvalues]: drawing coastaline ([proj],[res]) or None
2381          [proj]: projection
2382             * 'cyl', cilindric
2383             * 'lcc', lambert conformal
2384          [res]: resolution:
2385             * 'c', crude
2386             * 'l', low
2387             * 'i', intermediate
2388             * 'h', high
2389             * 'f', full
2390        [kindfigure]: kind of figure
2391          'legend': only points in the map with the legend with the names
2392          'labelled',[txtsize],[txtcol]: points with the names and size, color of text
2393        [pointcolor]: color for the points ('auto' for "red")
2394        [pointlabels]: ',' of labels [only used if [collab]='None'] ('None' for no labels)
2395        [locleg]: location of the legend (0, autmoatic)
2396          1: 'upper right', 2: 'upper left', 3: 'lower left', 4: 'lower right',
2397          5: 'right', 6: 'center left', 7: 'center right', 8: 'lower center',
2398          9: 'upper center', 10: 'center'
2399        [figureko]: kind of the output file (pdf, png, ...)
2400        [figuren]: name of the figure
2401      [filen]= [ncfile],[lonvarn],[latvarn][,[varn],[dimvals],[vargn],[min],[max],[cbar],[varu]]
2402        [ncfile]: netCDF to use to geolocalize the points
2403        [lonvarn]: name of the variable with the longitudes
2404        [latvarn]: name of the variable with the latitudes
2405        [varn]: optional variable to add staff into the graph
2406        [dimval]: '@' list of [dimn]|[dimval] to get the values for [varn]
2407          [dimn]: name of the dimension
2408          [dimval]: value of the dimension (no value all range)
2409        [vargn]: name of the variable in the graph
2410        [min]: minimum value for the extra variable
2411        [max]: maximum value for the extra variable
2412        [cbar]: color bar
2413        [varu]: units of the variable
2414    """
2415    fname = 'draw_points'
2416
2417    if values == 'h':
2418        print fname + '_____________________________________________________________'
2419        print draw_points.__doc__
2420        quit()
2421
2422    expectargs = '[ptasciifile]:[gtit]:[mapvalues]:[kindfigure]:[pointcolor]:' +     \
2423      '[pointlabels]:[locleg]:[figurek]:[figuren]'
2424 
2425    drw.check_arguments(fname,values,expectargs,':')
2426
2427    ptasciifile = values.split(':')[0]
2428    gtit = values.split(':')[1]
2429    mapvalues = values.split(':')[2]
2430    kindfigure = values.split(':')[3]
2431    pointcolor = values.split(':')[4]
2432    pointlabels = values.split(':')[5]
2433    locleg = int(values.split(':')[6])
2434    figureko = values.split(':')[7]
2435    figuren = values.split(':')[8]
2436
2437# Getting station locations
2438##
2439    filev = ptasciifile.split(',')[0]
2440    comchar = ptasciifile.split(',')[1].split('|')
2441    collon = int(ptasciifile.split(',')[2])
2442    collat = int(ptasciifile.split(',')[3])
2443    collab = ptasciifile.split(',')[4]
2444
2445    if not os.path.isfile(filev):
2446        print errormsg
2447        print '  ' + fname + ": file '" + filev + "' does not exist!!"
2448        quit(-1)
2449
2450# Getting points position and labels
2451    oascii = open(filev, 'r')
2452    xptval = []
2453    yptval = []
2454    if collab != 'None':
2455        ptlabels = []
2456        for line in oascii:
2457            if not drw.searchInlist(comchar, line[0:1]):
2458                linevals = drw.reduce_spaces(line)
2459                xptval.append(np.float(linevals[collon].replace('\n','')))
2460                yptval.append(np.float(linevals[collat].replace('\n','')))
2461                ptlabels.append(linevals[int(collab)].replace('\n',''))
2462    else:
2463        ptlabels = None
2464        for line in oascii:
2465            if  not drw.searchInlist(comchar, line[0:1]):
2466                linevals = drw.reduce_spaces(line)
2467                xptval.append(np.float(linevals[collon].replace('\n','')))
2468                yptval.append(np.float(linevals[collat].replace('\n','')))
2469
2470    oascii.close()
2471
2472    if pointlabels != 'None' and collab == 'None':
2473        ptlabels = pointlabels.split(',')
2474
2475# Getting localization of the points
2476##
2477    filev = filen.split(',')
2478    Nvals = len(filev)
2479
2480    ncfile = filev[0]
2481    lonvarn = filev[1]
2482    latvarn = filev[2]
2483    varn = None
2484    varextrav = None
2485    if Nvals == 10:
2486        varn = filev[3]
2487        dimvals = filev[4]
2488        varextrav = [filev[5], np.float(filev[6]), np.float(filev[7]), filev[8],     \
2489          filev[9]]
2490   
2491    if not os.path.isfile(ncfile):
2492        print errormsg
2493        print '  ' + fname + ": file '" + ncfile + "' does not exist!!"
2494        quit(-1)
2495
2496    onc = NetCDFFile(ncfile, 'r')
2497   
2498    lonv, latv = drw.lonlat2D(onc.variables[lonvarn], onc.variables[latvarn])
2499
2500    if varn is not None:
2501        objextra = onc.variables[varn]
2502        vard = objextra.dimensions
2503        dd = {}
2504        for dn in dimvals.split('@'):
2505            ddn = dn.split('|')[0]
2506            ddv = dn.split('|')[1]
2507            dd[ddn] = ddv
2508        slicevar = []
2509        for dv in vard:
2510            found= False
2511            for dn in dd.keys():
2512                if dn == dv:
2513                    slicevar.append(int(dd[dn]))
2514                    found = True
2515                    break
2516            if not found:
2517                slicevar.append(slice(0,len(onc.dimensions[dv])))
2518
2519        varextra = np.squeeze(objextra[tuple(slicevar)])
2520
2521    if mapvalues == 'None':
2522        mapV = None
2523    else:
2524        mapV = mapvalues
2525
2526    drw.plot_points(xptval, yptval, lonv, latv, varextra, varextrav, gtit, mapV,     \
2527      kindfigure, pointcolor, ptlabels, locleg, figureko, figuren)
2528
2529    onc.close()
2530
2531    return
2532
2533def draw_points_lonlat(filen, values):
2534    """ Function to plot a series of lon/lat points
2535     filen= name of the file
2536     values= [lonvarname]:[latvarname]:[gkind]:[gtit]:[ptcolor]:[pttype]:[ptsize]:[labels]:[locleg]:[figureK]
2537       [lonvarname]: name of the variable longitude
2538       [latvarname]: name of the variable latitude
2539       [gkind]: kind of graphical output
2540       [gtit]: graphic title '!' for spaces
2541       [ptcolor]: color of the points ('auto', for "red")
2542       [pttype]: type of point
2543       [ptsize]: size of point
2544       [labels]: ',' list of labels to use
2545       [locleg]: location of the legend (0, automatic)
2546         1: 'upper right', 2: 'upper left', 3: 'lower left', 4: 'lower right',
2547         5: 'right', 6: 'center left', 7: 'center right', 8: 'lower center',
2548         9: 'upper center', 10: 'center'
2549       [figureK]= kind of figure
2550         'legend': only points in the map with the legend with the names
2551         'labelled',[txtsize],[txtcol]: points with the names and size, color of text
2552    """
2553    fname = 'draw_points_lonlat'
2554
2555    if values == 'h':
2556        print fname + '_____________________________________________________________'
2557        print draw_points_lonlat.__doc__
2558        quit()
2559
2560    expectargs = '[lonvarname]:[latvarname]:[gkind]:[gtit]:[ptcolor]:[pttype]:' +    \
2561      '[ptsize]:[labels]:[locleg]:[figureK]'
2562 
2563    drw.check_arguments(fname,values,expectargs,':')
2564
2565    lonname = values.split(':')[0]
2566    latname = values.split(':')[1]
2567    kindfigure = values.split(':')[2]
2568    gtit = values.split(':')[3].replace('!',' ')
2569    pointcolor = values.split(':')[4]
2570    pointtype = values.split(':')[5]
2571    pointsize = np.float(values.split(':')[6])
2572    labelsv = values.split(':')[7]
2573    loclegend = values.split(':')[8]
2574    figureK = values.split(':')[9]
2575
2576    fname = 'points_lonlat'
2577
2578    onc = NetCDFFile(filen, 'r')
2579    if not onc.variables.has_key(lonname):
2580        print errormsg
2581        print fname + ": file '" + filen + "' does not have longitudes '" + lonname +\
2582          "' !!"
2583        quit(-1)
2584    if not onc.variables.has_key(lonname):
2585        print errormsg
2586        print fname + ": file '" + filen + "' does not have longitudes '" + lonname +\
2587          "' !!"
2588        quit(-1)
2589   
2590    olon = onc.variables[lonname]
2591    olat = onc.variables[latname]
2592
2593    Ndimlon = len(olon.shape)
2594    if Ndimlon == 1:
2595        dx = olon.shape[0]
2596        dy = olat.shape[0]
2597        if dx == dy: 
2598            lonvals = olon[:]
2599            latvals = olat[:]
2600        else: 
2601            lonvals0 = np.zeros((dy,dx), dtype=np.float)
2602            latvals0 = np.zeros((dy,dx), dtype=np.float)
2603            for iL in range(dy):
2604                lonvals0[iL,:] = olon[:]
2605            for il in range(dx):
2606                latvals0[:,il] = olat[:]
2607            lonvals = lonvals0.flatten()
2608            latvals = latvals0.flatten()
2609
2610    elif Ndimlon == 2:
2611        lonvals = olon[:].flatten()
2612        latvals = olat[:].flatten()
2613    elif Ndimlon == 3:
2614        lonvals = olon[1,:,:].flatten()
2615        latvals = olat[1,:,:].flatten()
2616# Playing for Anna
2617#        lonvals = olon[:].flatten()
2618#        latvals = olat[:].flatten()
2619    elif Ndimlon == 4:
2620        lonvals = olon[1,0,:,:].flatten()
2621        latvals = olat[1,0,:,:].flatten()
2622    else:
2623        print errormsg
2624        print '  ' + fname + ': longitude size:',len(olon),' not ready!!'
2625        quit(-1)
2626
2627    if labelsv == 'None':
2628        labels = None
2629    else:
2630        labels = labelsv.split(',')
2631
2632    drw.plot_list_points(lonvals, latvals, lonname, latname, gtit, figureK, pointcolor, pointtype,    \
2633      pointsize, labels, loclegend, kindfigure, fname)
2634
2635    onc.close()
2636
2637    return
2638
2639def draw_timeSeries(filen, values, variables):
2640    """ Function to draw a time-series
2641    draw_timeSeries(filen, values, variable):
2642      filen= name of the file
2643      values= [gvarname]:[timetit]:[tkind]:[timefmt]:[title]:[locleg]:[gkind]:[colorlines]:[pointtype]:[pointfreq]
2644      [gvarname]: name of the variable to appear in the graph
2645      [timetit]: title of the time axis (assumed x-axis, '|' for spaces)
2646      [tkind]: kind of time to appear in the graph (assumed x-axis)
2647        'Nval': according to a given number of values as 'Nval',[Nval]
2648        'exct': according to an exact time unit as 'exct',[tunit];
2649          tunit= [Nunits],[tu]; [tu]= 'c': centuries, 'y': year, 'm': month,
2650            'w': week, 'd': day, 'h': hour, 'i': minute, 's': second,
2651            'l': milisecond
2652      [timefmt]: format of the time labels
2653      [title]: title of the graphic ('|' for spaces)
2654      [locleg]: location of the legend (0, automatic)
2655        1: 'upper right', 2: 'upper left', 3: 'lower left', 4: 'lower right',
2656        5: 'right', 6: 'center left', 7: 'center right', 8: 'lower center',
2657        9: 'upper center', 10: 'center'
2658      [gkind]: kind of graphical output
2659      [colorlines]: ',' list of colors for the lines, None for automatic, single
2660          value all the same
2661      [pointtype]: ',' list of type of points for the lines, None for automatic, single
2662          value all the same
2663      [pointfreq]: frequency of point plotting, 'all' for all time steps
2664      variables= [varname],[timename] names of variable and variable with times
2665      draw_timeSeries('wrfout_d01_1979-12-01_00:00:00_bottom_top_B6-E6-I1_south_north_B3-E3-I1_west_east_B26-E26-I1.nc', 'dt_con:time|($[DD]^{[HH]}$):exct,12,h:$%d^{%H}$:time|evolution|at|-1|6|3|26:1:pdf', 'LDQCON,time')
2666    """
2667
2668    fname = 'draw_timeSeries'
2669
2670    if values == 'h':
2671        print fname + '_____________________________________________________________'
2672        print draw_timeSeries.__doc__
2673        quit()
2674
2675    expectargs = '[gvarname]:[timetit]:[tkind]:[timefmt]:[title]:' +                 \
2676      '[locleg]:[gkind]:[colorlines]:[pointtype]:[pointfreq]'
2677 
2678    drw.check_arguments(fname,values,expectargs,':')
2679
2680    gvarname = values.split(':')[0]
2681    timetit = values.split(':')[1].replace('|',' ')
2682    tkind = values.split(':')[2]
2683    timefmt = values.split(':')[3]
2684    title = values.split(':')[4].replace('|',' ')
2685    locleg = int(values.split(':')[5])
2686    gkind = values.split(':')[6]
2687    colorlines = values.split(':')[7]
2688    pointtype = values.split(':')[8]
2689    pointfreq0 = values.split(':')[9]
2690   
2691    ncobj = NetCDFFile(filen, 'r')
2692
2693    variable = variables.split(',')[0]
2694    timevar = variables.split(',')[1]
2695
2696    if not ncobj.variables.has_key(variable):
2697        print errormsg
2698        print '  ' + fname + ": file '" +  filen + "' does not have variable '" +    \
2699          variable + "' !!"
2700        quit(-1)
2701
2702    if not ncobj.variables.has_key(timevar):
2703        print errormsg
2704        print '  ' + fname + ": file '" +  filen + "' does not have variable time '" \
2705          + timevar + "' !!"
2706        quit(-1)
2707
2708    varobj = ncobj.variables[variable]
2709    timeobj = ncobj.variables[timevar]
2710
2711    dimt = len(timeobj[:])
2712    varvals = np.zeros((2,dimt), dtype=np.float)
2713
2714    gunits = varobj.getncattr('units')
2715    tunits = timeobj.getncattr('units')
2716
2717    varvals[0,:], valpot, newgunits, Spot = drw.pot_values(varobj[:].flatten(), gunits)
2718    varvals[1,:] = timeobj[:]
2719
2720    tseriesvals = []
2721    tseriesvals.append(varvals)
2722
2723    if colorlines == 'None': 
2724        collines = None
2725    else:
2726        collines = colorlines.split(',')
2727    if pointtype == 'None': 
2728        pttype = None
2729    else:
2730        pttype = pointtype.split(',')
2731
2732    if pointfreq0 == 'all':
2733        pointfreq = None
2734    else:
2735        pointfreq = int(pointfreq0)
2736
2737    drw.plot_TimeSeries(tseriesvals, Spot + drw.units_lunits(gunits), tunits,        \
2738      'TimeSeries', gvarname, timetit, tkind, timefmt, title,                        \
2739      gvarname.replace('_','\_'), locleg, gkind, collines, pttype, pointfreq)
2740
2741    return
2742
2743#draw_timeSeries('wrfout_d01_1979-12-01_00:00:00_bottom_top_B6-E6-I1_south_north_B3-E3-I1_west_east_B26-E26-I1.nc', 'dt_con:time|($[DD]^{[HH]}$):exct,12,h:$%d^{%H}$:time|evolution|at|-1|6|3|26:1:pdf:None:None', 'LDQCON,time')
2744
2745def draw_trajectories(trjfilens, values, observations):
2746    """ Function to draw different trajectories at the same time
2747    draw_trajectories(trjfilens, values, observations):
2748      trjfilens= [filen]@[Tint]@[map] ',' separated list of files with trajectories,
2749         time intervals and reference maps (first one will be used to plot)
2750        [filen]: name of the file to use (lines with '#', not readed) as:
2751          [t-step] [x] [y]
2752        [Tint]: interval of time as [Tbeg]@[Tend] or -1 for all the interval
2753        [map]: [file]#[lonname]#[latname]
2754          [file]; with the [lon],[lat] matrices
2755          [lonname],[latname]; names of the longitudes and latitudes variables
2756      values=[leglabels]|[lonlatlims]|[title]|[graphk]|[mapkind]
2757        [leglabels]: ',' separated list of names for the legend
2758        [lonlatlims]: ',' list of limits of the map [lonmin, latmin, lonmax, latmax] or None
2759        [title]: title of the plot ('!' for spaces)
2760        [graphk]: kind of the graphic
2761        [mapkind]: drawing coastaline ([proj],[res]) or None
2762          [proj]: projection
2763             * 'cyl', cilindric
2764             * 'lcc', lambert conformal
2765          [res]: resolution:
2766             * 'c', crude
2767             * 'l', low
2768             * 'i', intermediate
2769             * 'h', high
2770             * 'f', full
2771      obsevations= [obsfile],[obsname],[Tint],[null]
2772        [obsfile]: name fo the File with the observations as [t-step] [lat] [lon]
2773        [obsname]: name of the observations in the graph
2774        [Tint]: interval of time as [Tbeg]@[Tend] or -1 for all the interval
2775        [null]: null value for the observed trajectory
2776    """
2777
2778    fname = 'draw_trajectories'
2779
2780    if values == 'h':
2781        print fname + '_____________________________________________________________'
2782        print draw_trajectories.__doc__
2783        quit()
2784
2785    expectargs = '[leglabels]|[lonlatlims]|[title]|[graphk]|[mapkind]'
2786 
2787    drw.check_arguments(fname,values,expectargs,'|')
2788
2789    trjfiles = trjfilens.split(',')
2790    leglabels = values.split('|')[0]
2791    lonlatlims = values.split('|')[1]
2792    title = values.split('|')[2].replace('!',' ')
2793    graphk = values.split('|')[3]
2794    mapkind = values.split('|')[4]
2795
2796    Nfiles = len(trjfiles)
2797
2798# Getting trajectotries
2799##
2800
2801    lontrjvalues = []
2802    lattrjvalues = []
2803
2804    print '  ' + fname
2805    ifn = 0
2806    for ifile in trjfiles:
2807        filen = ifile.split('@')[0]
2808        Tint = ifile.split('@')[1]
2809
2810        print '    trajectory:',filen
2811
2812        if Tint != '-1':
2813            Tbeg = Tint
2814            Tend = ifile.split('@')[2]
2815            mapv = ifile.split('@')[3]
2816        else:
2817            mapv = ifile.split('@')[2]
2818
2819        if not os.path.isfile(filen):
2820            print errormsg
2821            print '  ' + fname + ": trajectory file '" + filen + "' does not exist !!"
2822            quit(-1)
2823
2824# Charging longitude and latitude values
2825##
2826        lonvals, latvals = drw.lonlat_values(mapv.split('#')[0], mapv.split('#')[1], \
2827          mapv.split('#')[2])
2828
2829        if ifn == 0: mapref = mapv
2830        ifn = ifn + 1
2831
2832        objfile = open(filen, 'r')
2833        trjtimev = []
2834        trjxv = []
2835        trjyv = []
2836
2837        for line in objfile:
2838            if line[0:1] != '#':
2839                trjtimev.append(int(line.split(' ')[0]))
2840                trjxv.append(int(line.split(' ')[1]))
2841                trjyv.append(int(line.split(' ')[2]))
2842
2843        objfile.close()
2844
2845        if Tint != '-1':
2846            lontrjvalues.append(lonvals[trjyv[Tint:Tend+1], trjxv[Tint:Tend+1]])
2847            lattrjvalues.append(latvals[trjyv[Tint:Tend+1], trjxv[Tint:Tend+1]])
2848        else:
2849            lontrjvalues.append(lonvals[trjyv[:], trjxv[:]])
2850            lattrjvalues.append(latvals[trjyv[:], trjxv[:]])
2851
2852# lonlatlimits
2853##
2854
2855    if lonlatlims == 'None':
2856        lonlatlimsv = None
2857    else:
2858        lonlatlimsv = np.zeros((4), dtype=np.float)
2859        lonlatlimsv[0] = np.float(lonlatlims.split(',')[0])
2860        lonlatlimsv[1] = np.float(lonlatlims.split(',')[1])
2861        lonlatlimsv[2] = np.float(lonlatlims.split(',')[2])
2862        lonlatlimsv[3] = np.float(lonlatlims.split(',')[3])
2863
2864# lon/lat objects
2865##
2866    objnc = NetCDFFile(mapref.split('#')[0])
2867    lonobj = objnc.variables[mapref.split('#')[1]]
2868    latobj = objnc.variables[mapref.split('#')[2]]
2869
2870# map
2871##
2872    if mapkind == 'None':
2873        mapkindv = None
2874    else:
2875        mapkindv = mapkind
2876
2877    if observations is None:
2878        obsname = None
2879    else:
2880        obsfile = observations.split(',')[0]
2881        obsname = observations.split(',')[1]
2882        Tint = observations.split(',')[2]
2883        null = np.float(observations.split(',')[3])
2884        print '    observational trajectory:',obsfile
2885
2886        if not os.path.isfile(obsfile):
2887            print errormsg
2888            print '  ' + fname + ": observations trajectory file '" + obsfile +      \
2889              "' does not exist !!"
2890            quit(-1)
2891
2892        objfile = open(obsfile, 'r')
2893        obstrjtimev = []
2894        obstrjxv = []
2895        obstrjyv = []
2896
2897        for line in objfile:
2898            if line[0:1] != '#':
2899                lon = np.float(line.split(' ')[2])
2900                lat = np.float(line.split(' ')[1])
2901                if not lon == null and not lat == null:
2902                    obstrjtimev.append(int(line.split(' ')[0]))
2903                    obstrjxv.append(lon)
2904                    obstrjyv.append(lat)
2905                else:
2906                    obstrjtimev.append(int(line.split(' ')[0]))
2907                    obstrjxv.append(None)
2908                    obstrjyv.append(None)
2909
2910        objfile.close()
2911
2912        if Tint != '-1':
2913            Tint = int(observations.split(',')[2].split('@')[0])
2914            Tend = int(observations.split(',')[2].split('@')[1])
2915            lontrjvalues.append(obstrjxv[Tint:Tend+1])
2916            lattrjvalues.append(obstrjyv[Tint:Tend+1])
2917        else:
2918            lontrjvalues.append(obstrjxv[:])
2919            lattrjvalues.append(obstrjyv[:])
2920
2921    drw.plot_Trajectories(lontrjvalues, lattrjvalues, leglabels.split(','),          \
2922      lonobj, latobj, lonlatlimsv, title, graphk, mapkindv, obsname)
2923
2924    objnc.close()
2925
2926    return
2927
2928def draw_vals_trajectories(ncfile, values, variable):
2929    """ Function to draw values from the outputs from 'compute_tevolboxtraj'
2930    draw_vals_trajectories(ncfile, values, variable)
2931    ncfile= [ncfile] ',' list of files to use
2932    values= [statisticskind]:[Tint]:[labels]@[locleg]:[gvarname]:[timetit]:[tkind]:[timefmt]:[title]:[gkind]
2933      [statisticskind]=[statistics][kind]
2934        [statistics]: which statistics to use, from: 'center', 'min', 'max', 'mean',
2935        'mean2', 'stdev'
2936        [kind]: 'box', 'circle' statistics taking the values from a box or a circle
2937        'trj': value following the trajectory
2938      [Tint]: [Tbeg]@[Tend] or None, interval of time to plot or -1 for all the times
2939      [labels]: ',' separated list of labels for the legend
2940      [locleg]: location of the legend (0, automatic)
2941        1: 'upper right', 2: 'upper left', 3: 'lower left', 4: 'lower right',
2942        5: 'right', 6: 'center left', 7: 'center right', 8: 'lower center',
2943        9: 'upper center', 10: 'center'
2944      [gvarname]: name of the variable to appear in the graph
2945      [timetit]: title of the time axis (assumed x-axis, '|' for spaces)
2946      [tkind]: kind of time to appear in the graph (assumed x-axis)
2947        'Nval': according to a given number of values as 'Nval',[Nval]
2948        'exct': according to an exact time unit as 'exct',[tunit];
2949          tunit= [Nunits],[tu]; [tu]= 'c': centuries, 'y': year, 'm': month,
2950            'w': week, 'd': day, 'h': hour, 'i': minute, 's': second,
2951            'l': milisecond
2952      [timefmt]: format of the time labels
2953      [title]: title of the graphic ('|' for spaces)
2954      [gkind]: kind of graphical output
2955    variable= variable to use
2956    """
2957
2958    fname = 'draw_vals_trajectories'
2959
2960    if values == 'h':
2961        print fname + '_____________________________________________________________'
2962        print draw_vals_trajectories.__doc__
2963        quit()
2964
2965    sims = ncfile.split(',')
2966
2967    if len(values.split(':')) != 9:
2968        print errormsg
2969        print '  ' + fname  + ': wrong number of values!', len(values.split(':')),   \
2970          'given 9 needed!!'
2971        print '    ',values.split(':')
2972        quit(-1)
2973
2974    statisticskind = values.split(':')[0]
2975    Tint = values.split(':')[1]
2976    labels = values.split(':')[2]
2977    gvarname = values.split(':')[3]
2978    timetit = values.split(':')[4].replace('|',' ')
2979    tkind = values.split(':')[5]
2980    timefmt = values.split(':')[6]
2981    title = values.split(':')[7].replace('|',' ')
2982    gkind = values.split(':')[8]
2983
2984    leglabels = labels.split('@')[0].split(',')
2985    locleg = int(labels.split('@')[1])
2986
2987    Nsims = len(sims)
2988
2989    if Tint != '-1':
2990        tini = np.float(Tint.split('@')[0])
2991        tend = np.float(Tint.split('@')[1])
2992    else:
2993        tini = -1.
2994        tend = -1.
2995
2996    vartimetrjv = []
2997
2998    print '  ' + fname
2999    for trjfile in sims:
3000        print '    ' + trjfile + ' ...'
3001        if not os.path.isfile(trjfile):
3002            print errormsg
3003            print '  ' + fname + ": trajectory file: '" + trjfile +                  \
3004              "' does not exist !!"
3005            quit(-1)
3006
3007        trjobj = NetCDFFile(trjfile, 'r')
3008        otim = trjobj.variables['time']
3009        if not trjobj.variables.has_key(statisticskind + '_' + variable):
3010            print errormsg
3011            print '  ' + fname + ": file '" + trjfile + "' does not have variable '"+\
3012              statisticskind + '_' + variable + "' !!"
3013            quit(-1)
3014        ovar = trjobj.variables[statisticskind + '_' + variable]
3015        dimt = otim.shape[0]
3016
3017        if trjfile == sims[0]:
3018            gunits = ovar.getncattr('units')
3019            lname = ovar.getncattr('long_name')
3020            tunits = otim.getncattr('units')
3021
3022        if tini != -1:
3023            tiniid = -1
3024            tendid = -1       
3025            for itv in range(dimt):
3026                if otim[itv] <= tini and otim[itv+1] >= tini: tiniid = itv
3027                if otim[itv] <= tend and otim[itv+1] >= tend: tendid = itv
3028
3029            if tiniid == -1 or tendid == -1:
3030                print errormsg
3031                print '  ' + main + ' time interval ', tini,',',tend,' not found: ',     \
3032                  tendid, ',', tiniid, ' !!'
3033                print '    data interval [',otim[0], otim[dimt-1],']'
3034                quit(-1)
3035            dimt = tendid - tiniid + 1
3036
3037        else:
3038            dimt = otim.shape[0]
3039
3040        valsv = np.zeros((2,dimt), dtype=np.float)
3041# Checking for time consistency
3042        if otim.getncattr('units') != tunits:
3043            print warnmsg
3044            print '  ' + fname + ': different time units in the plot!!'
3045            newtimes = drw.coincident_CFtimes(otim[:], tunits, otim.getncattr('units'))
3046        else:
3047            newtimes = otim[:]
3048
3049        if tini == -1:
3050            valsv[1,:] = newtimes[:]
3051            valsv[0,:] = ovar[:]
3052        else:
3053            valsv[1,:] = newtimes[tiniid:tendid+1]
3054            valsv[0,:] = ovar[tiniid:tendid+1]
3055
3056        vartimetrjv.append(valsv)
3057        trjobj.close()
3058
3059    drw.plot_TimeSeries(vartimetrjv, drw.units_lunits(gunits), tunits,               \
3060      'val_trajectories_' + statisticskind, gvarname, timetit, tkind, timefmt, title,\
3061      leglabels, locleg, gkind)
3062
3063def variable_values(values):
3064    """ Function to give back values for a given variable
3065      values= [varname] name of the variable
3066    """
3067
3068    fname = 'variable_values'
3069
3070    values = drw.variables_values(values)
3071
3072    print fname,'values:',values
3073    print fname,'variable_name:',values[0]
3074    print fname,'standard_name:',values[1]
3075    print fname,'min,max:',str(values[2]) + ',' + str(values[3])
3076    print fname,'long_name:',values[4]
3077    print fname,'units:',values[5]
3078    print fname,'shad_colors:',values[6]
3079    print fname,'all_values:',drw.numVector_String(values,',')
3080
3081    return
3082
3083def draw_ptZvals(ncfile, values, variable):
3084    """ Function to plot a given list of points and values
3085      ncfile= netCDF file to use
3086      values= [fvname]:[XYvar]:[pointype]:[pointsize]:[graphlimits]:[nxtype]:
3087        [figuretitle]:[colorbar]:[mapvalue]:[kindfig]
3088        fvname: name of the variable in the graph
3089        XYvar: [lon],[lat] variable names
3090        ptype: type of the point
3091        ptsize: size of the point
3092        graphlimits: minLON,minLAT,maxLON,maxLAT limits of the graph 'None' for the full size
3093        nxtype: minimum and maximum type
3094          'auto': values taken from the extrems of the data
3095          [min],[max]: given minimum and maximum values
3096        figtitle: title of the figure
3097        cbar: color bar
3098        mapv: map characteristics: [proj],[res]
3099          see full documentation: http://matplotlib.org/basemap/
3100          [proj]: projection
3101            * 'cyl', cilindric
3102            * 'lcc', lambert-conformal
3103          [res]: resolution:
3104            * 'c', crude
3105            * 'l', low
3106            * 'i', intermediate
3107            * 'h', high
3108            * 'f', full
3109        kfig: kind of figure
3110      variable= name of the variable to plot
3111    """
3112    fname = 'draw_ptZvals'
3113    import numpy.ma as ma
3114
3115    if values == 'h':
3116        print fname + '_____________________________________________________________'
3117        print draw_ptZvals.__doc__
3118        quit()
3119
3120    expectargs = '[fvname]:[XYvar]:[pointype]:[pointsize]:[graphlmits]:[nxtype]:' +  \
3121      '[figuretit]:[colorbar]:[mapvalue]:[kindfig]'
3122 
3123    drw.check_arguments(fname,values,expectargs,':')
3124
3125    fvname = values.split(':')[0]
3126    XYvar = values.split(':')[1]
3127    pointype = values.split(':')[2]
3128    pointsize = int(values.split(':')[3])
3129    graphlimits = values.split(':')[4]
3130    nxtype = values.split(':')[5]
3131    figuretitle = values.split(':')[6].replace('!',' ')
3132    colorbar = values.split(':')[7]
3133    mapvalue = values.split(':')[8]
3134    kindfig = values.split(':')[9]
3135
3136    onc = NetCDFFile(ncfile, 'r')
3137   
3138    if not onc.variables.has_key(variable):
3139        print errormsg
3140        print '  ' + fname + ": file '" + ncfile + "' does not have variable '" +    \
3141          variable + "' !!"
3142        quit(-1)
3143
3144# points
3145    lonvarn = XYvar.split(',')[0]
3146    latvarn = XYvar.split(',')[1]
3147
3148    if not onc.variables.has_key(lonvarn):
3149        print errormsg
3150        print '  ' + fname + ": file '" + ncfile + "' does not have longitude " +    \
3151          "variable '" + lonvarn + "' !!"
3152        quit(-1)
3153   
3154    if not onc.variables.has_key(latvarn):
3155        print errormsg
3156        print '  ' + fname + ": file '" + ncfile + "' does not have latitude " +    \
3157          "variable '" + latvarn + "' !!"
3158        quit(-1)
3159
3160    olonvar = onc.variables[lonvarn]
3161    olatvar = onc.variables[latvarn]
3162    ovarvar = onc.variables[variable]
3163
3164    Lpts = len(olonvar[:].flatten())
3165
3166    pointvalues = ma.masked_array(np.zeros((Lpts,3), dtype=np.float))
3167    pointvalues[:,0] = olonvar[:].flatten()
3168    pointvalues[:,1] = olatvar[:].flatten()
3169    pointvalues[:,2] = ovarvar[:].flatten()
3170
3171    varattrs = ovarvar.ncattrs()
3172    if drw.searchInlist(varattrs, 'units'):
3173        fvunits = ovarvar.getncattr('units')
3174    else:
3175        fvunits = drw.variables_values(variable)[5]
3176
3177# map value
3178    if mapvalue == 'None': mapvalue = None
3179
3180# Graph limits
3181    if graphlimits != 'None':
3182        graphlts = np.zeros((4), dtype=np.float)
3183        for il in range(4): graphlts[il] = np.float(graphlimits.split(',')[il])
3184        pointvalues[:,0] = ma.masked_outside(pointvalues[:,0], graphlts[0],          \
3185          graphlts[2])
3186        pointvalues[:,1] = ma.masked_outside(pointvalues[:,1], graphlts[3],          \
3187          graphlts[2])
3188
3189#        for ip in range(Lpts): 
3190#            if pointvalues[ip,0] < graphlts[0] or pointvalues[ip,0] > graphlts[2]    \
3191#              or pointvalues[ip,1] < graphlts[1] or pointvalues[ip,1] > graphlts[3]:
3192#                print ip,pointvalues[ip,0:2], graphlts
3193#                pointvalues[ip,2] = None
3194    else:
3195        graphlts = None
3196
3197    drw.plot_ptZvals(fvname,fvunits,pointvalues,pointype,pointsize,graphlts, nxtype, \
3198      figuretitle,colorbar,mapvalue,kindfig)
3199
3200    return
3201
3202#draw_ptZvals('OBSnetcdf.nc', 'pracc:lon,lat:o:80:2,42,7,47,:values!of!values:Blues:cyl,l:pdf', 'pr')
3203
3204def draw_vectors(ncfile, values, varns):
3205    """ Function to plot wind vectors
3206      values= [dimname]|[vardimname]|[value]:[vecvals]:[windlabs]:[mapvalues]:
3207        [gtit]:[kindfig]:[figuren]
3208        'X/Y/Z/T'|[dimname]|[vardimname]|[value]: ',', list for each basic dimension '|' separated of:
3209          [dimname]: name of the dimension in the file
3210          [vardimname]: name of the variable with the values for the dimension in the file
3211          [value]: which value of the given dimension (-1, all; [ibeg]@[iend], i-range beginning, end)
3212          No value takes all the range of the dimension
3213        [vecvals]= [frequency],[color],[length]
3214          [frequency]: [xfreq]@[yfreq] frequency of values allong each axis ('None', all grid points;
3215            'auto', computed automatically to have 20 vectors along each axis)
3216          [color]: color of the vectors
3217            'singlecol'@[colorn]: all the vectors same color ('auto': for 'red')
3218            'wind'@[colorbar]: color of the vectors according to wind speed sqrt(u^2+v^2) and given [colorbar]
3219              all vectors the same length
3220            '3rdvar'@[colorbar]@[varn]@[units]: color of the vectors according to a 3rd variable (to be added at -v) and given [colorbar]
3221              all vectors the same length
3222          [length]: length of the wind vectors ('auto', for 9)
3223        [windlabs]= [windname],[windunits]
3224          [windname]: name of the wind variable in the graph
3225          [windunits]: units of the wind variable in the graph ('None', for the value in the file)
3226        [mapvalues]= map characteristics: [proj],[res]
3227          see full documentation: http://matplotlib.org/basemap/
3228          [proj]: projection
3229            * 'cyl', cilindric
3230            * 'lcc', lambert conformal
3231          [res]: resolution:
3232            * 'c', crude
3233            * 'l', low
3234            * 'i', intermediate
3235            * 'h', high
3236            * 'f', full
3237        gtit= title of the graph ('|', for spaces)
3238        kindfig= kind of figure
3239        figuren= name of the figure
3240      ncfile= file to use
3241      varns= [uwind],[ywind] ',' list of the name of the variables with the u-wind,y-wind component
3242    """
3243    fname = 'draw_vectors'
3244
3245    if values == 'h':
3246        print fname + '_____________________________________________________________'
3247        print draw_vectors.__doc__
3248        quit()
3249
3250    expectargs = '[X/Y/Z/T]|[dimname]|[vardimname]|[value]:[vecvals]:[windlabs]:' +  \
3251      '[mapvalues]:[gtit]:[kindfig]:[figuren]'
3252 
3253    drw.check_arguments(fname,values,expectargs,':')
3254
3255    dimvals = values.split(':')[0]
3256    vecvals = values.split(':')[1]
3257    windlabels = values.split(':')[2]
3258    mapvalues = values.split(':')[3]
3259    gtit = values.split(':')[4]
3260    kindfig = values.split(':')[5]
3261    figuren = values.split(':')[6]
3262
3263    of = NetCDFFile(ncfile,'r')
3264
3265    dims = {}
3266    for dimv in dimvals.split(','):
3267        dns = dimv.split('|')
3268        dims[dns[0]] = [dns[1], dns[2], dns[3]]
3269
3270    varNs = []
3271    for dn in dims.keys():
3272        if dn == 'X':
3273            varNs.append(dims[dn][1])
3274            dimx = len(of.dimensions[dims[dn][0]])
3275        elif dn == 'Y':
3276            varNs.append(dims[dn][1])
3277            dimy = len(of.dimensions[dims[dn][0]])
3278
3279    ivar = 0
3280    for wvar in varns.split(','):
3281        if not drw.searchInlist(of.variables.keys(), wvar):
3282            print errormsg
3283            print '  ' + fname + ": file does not have variable '" + wvar + "' !!"
3284            quit(-1)
3285        if ivar == 0:
3286            varNs.append(wvar)
3287        else:
3288            varNs.append(wvar)
3289
3290    ivar = 0
3291    for varN in varNs:
3292        varslice = []
3293
3294        ovarN = of.variables[varN]
3295        vard = ovarN.dimensions
3296        for vdn in vard:
3297            found = False
3298            for dd in dims.keys():
3299                if dims[dd][0] == vdn:
3300                    if dims[dd][2].find('@') != -1:
3301                        rvals = dims[dd][2].split('@')
3302                        varslice.append(slice(int(rvals[0]), int(rvals[1])))
3303                    elif dims[dd][2] == '-1':
3304                        varslice.append(slice(0,len(of.dimensions[dims[dd][0]])))
3305                    else:
3306                        varslice.append(int(dims[dd][2]))
3307
3308                    found = True
3309                    break
3310            if not found:
3311                varslice.append(slice(0,len(of.dimensions[dims[dd][0]])))
3312
3313        if varN == dims['X'][1]:
3314            lonvals0 = np.squeeze(ovarN[tuple(varslice)])
3315        elif varN == dims['Y'][1]:
3316            latvals0 = np.squeeze(ovarN[tuple(varslice)])
3317        elif ivar == 2:
3318            uwvals = np.squeeze(np.array(ovarN[tuple(varslice)]))
3319        elif ivar == 3:
3320            vwvals = np.squeeze(ovarN[tuple(varslice)])
3321
3322        ivar = ivar + 1
3323
3324#    print 'Final shapes:',lonvals0.shape,':',latvals0.shape,':',uwvals.shape,':',
3325#      vwvals.shape
3326
3327    if len(uwvals.shape) != 2 or len(vwvals.shape) != 2:
3328        print errormsg
3329        print '  ' + fname + ': wrong size of the wind fields! they must be ' +      \
3330          '2-dimensional!'
3331        print '    u-winds shape:',uwvals.shape,'dims:',of.variables[varNs[2]]
3332        print '    v-winds shape:',vwvals.shape,'dims:',of.variables[varNs[3]]
3333        print '      provide more values for their dimensions!!'
3334        quit(-1)
3335
3336    if len(lonvals0.shape) == 1:
3337        lonvals, latvals = np.meshgrid(lonvals0, latvals0)
3338    else:
3339        lonvals = lonvals0
3340        latvals = latvals0
3341
3342# Vector values
3343    if vecvals.split(',')[0] == 'None':
3344        freqv = None
3345    else:
3346        freqv = vecvals.split(',')[0] 
3347    colorvals = vecvals.split(',')[1]
3348    coln = colorvals.split('@')[0]
3349    colv = colorvals.split('@')[1]
3350    if coln == 'singlecol':
3351        colorv = colv
3352    elif coln == 'wind':
3353        colorv = np.sqrt(uwvals**2 + vwvals**2)
3354    elif coln == '3rdvar':
3355        if len(varn.split(',')) != 3:
3356            print errormsg
3357            print '  ' + fname + ": color of vectors should be according to '" +     \
3358              coln + "' but a third varibale is not provided !!"
3359            quit(-1)
3360        ocolvec = of.variables[varNs[4]]
3361        colorv = ocolvec[:]
3362        stdvn, lonvn, unitsvn = drw.var_3desc(ocolvec)
3363        colorvals = colorvals + '@' + stdvn + '@' + unitsvn
3364    else:
3365        print errormsg
3366        print '  ' + fname + ": color type '" + coln + "' not ready !!"
3367        quit(-1)
3368
3369    lengthv = vecvals.split(',')[2]
3370
3371# Vector labels
3372    windname = windlabels.split(',')[0]
3373    windunits = windlabels.split(',')[1]
3374
3375    drw.plot_vector(lonvals, latvals, uwvals, vwvals, freqv, colorvals, colorv,      \
3376      lengthv, windname, windunits, mapvalues, gtit, kindfig, figuren)
3377
3378    of.close()
3379
3380    return
3381
3382def draw_basins(ncfile, values, varns):
3383    """ Function to plot river basins with their discharge vector and basins id (from 'routing.nc')
3384      values= [lonlatbox]:[mapres]:[cbarname]:[xtrmbasin]:[mapdraw]:[veclength]:[freq]:
3385        [ifreq]:[plotcountry]:[basinidn]:[gtit]:[kindfig]:[figuren]
3386        [lonlatbox]= [lonSW],[lonNE],[latSW],[latNE] coordinates of the lon/lat box
3387        [mapres]= resolution of the mapping information
3388        [cbarname]= colorbar name for the colors
3389        [xtrmbasin]= [minbasin],[maxbasin] minimum and maximum basin numbers
3390        [mapdraw]= whether to draw the map (and project the data) or not ('True/False')
3391        [veclength]= length of the vectors of discharge at each grid cell
3392        [freq]= frequency of values allong each axis (None, all grid points;
3393          'auto', computed automatically to have 20 vectors along each axis)
3394        [plotcountry]= whether country lines should be plotted or not ('True/False')
3395        [plotbasinid]= whether id of the basins should be plotted or not ('True/False')
3396        [gtit]= title of the graph ('|', for spaces)
3397        [kindfig]= kind of figure
3398        [figuren]= name of the figure
3399      ncfile= file to use
3400    """
3401    fname = 'draw_basins'
3402
3403    if values == 'h':
3404        print fname + '_____________________________________________________________'
3405        print draw_vectors.__doc__
3406        quit()
3407
3408    expectargs = '[lonlatbox]:[mapres]:[cbarname]:[xtrmbasin]:[mapdraw]:' +          \
3409      '[veclength]:[freq]:[plotcountry]:[basinidn]:[gtit]:[kindfig]:[figuren]'
3410 
3411    drw.check_arguments(fname,values,expectargs,':')
3412
3413    varn='basins'
3414    lonname = 'nav_lon'
3415    latname = 'nav_lat'
3416    flowname = 'trip'
3417
3418    lonlims =[]
3419    latlims =[]
3420
3421    lonlims.append(np.float(values.split(':')[0].split(',')[0]))
3422    lonlims.append(np.float(values.split(':')[0].split(',')[1]))
3423    latlims.append(np.float(values.split(':')[0].split(',')[2]))
3424    latlims.append(np.float(values.split(':')[0].split(',')[3]))
3425    map_res = values.split(':')[1]
3426    cbarname = values.split(':')[2]
3427    vtit = 'basins'
3428    minbasin = np.int(values.split(':')[3].split(',')[0])
3429    maxbasin = np.int(values.split(':')[3].split(',')[1])
3430    mapdraw = gen.Str_Bool(values.split(':')[4])
3431    veclength = np.float(values.split(':')[5])
3432    freq0 = values.split(':')[6]
3433    plotcountry = gen.Str_Bool(values.split(':')[7])
3434    plotbasinid = gen.Str_Bool(values.split(':')[8])
3435    gtit = values.split(':')[9].replace('|',' ')
3436    kindfig = values.split(':')[10]
3437    figuren = values.split(':')[11]
3438
3439    if freq0 == 'None': freq = None
3440
3441    ofile = NetCDFFile(ncfile, 'r')
3442
3443    obasins = ofile.variables[varn]
3444    olon = ofile.variables[lonname]
3445    olat = ofile.variables[latname]
3446    oflow = ofile.variables[flowname]
3447
3448    lons = olon[:]
3449    lats = olat[:]
3450
3451    lon, lat = drw.lonlat2D(lons, lats)
3452
3453    nlon = lonlims[0]
3454    xlon = lonlims[1]
3455    nlat = latlims[0]
3456    xlat = latlims[1]
3457
3458    imin, imax, jmin, jmax = gen.ijlonlat(lon, lat, nlon, xlon, nlat, xlat)
3459
3460    drw.plot_basins(lon[jmin:jmax,imin:imax], lat[jmin:jmax,imin:imax],              \
3461      oflow[jmin:jmax,imin:imax], freq, cbarname+'@basin@-',                         \
3462      obasins[jmin:jmax,imin:imax], veclength, minbasin, maxbasin, 'outflow', '-',   \
3463      'cyl,'+map_res, plotcountry, plotbasinid, gtit, kindfig, figuren)
3464
3465    ofile.close()
3466
3467    return
3468
3469def draw_basinsold(ncfile, values, varns):
3470    """ Function to plot wind basins
3471      values= [dimname]|[vardimname]|[value]:[vecvals]:[windlabs]:[mapvalues]:
3472        [gtit]:[kindfig]:[figuren]
3473        'X/Y/Z/T'|[dimname]|[vardimname]|[value]: ',', list for each basic dimension '|' separated of:
3474          [dimname]: name of the dimension in the file
3475          [vardimname]: name of the variable with the values for the dimension in the file
3476          [value]: which value of the given dimension (-1, all; [ibeg]@[iend], i-range beginning, end)
3477          No value takes all the range of the dimension
3478        [vecvals]= [frequency],[color],[length]
3479          [frequency]: [xfreq]@[yfreq] frequency of values allong each axis ('None', all grid points;
3480            'auto', computed automatically to have 20 vectors along each axis)
3481          [color]: [colorbar]@[varn]@[units]: color of the vectors according to a 3rd variable (to be added at -v) and given [colorbar]
3482              all vectors the same length
3483          [length]: length of the wind vectors ('auto', for 9)
3484        [windlabs]= [windname],[windunits]
3485          [windname]: name of the wind variable in the graph
3486          [windunits]: units of the wind variable in the graph ('None', for the value in the file)
3487        [mapvalues]= map characteristics: [proj],[res]
3488          see full documentation: http://matplotlib.org/basemap/
3489          [proj]: projection
3490            * 'cyl', cilindric
3491            * 'lcc', lambert conformal
3492          [res]: resolution:
3493            * 'c', crude
3494            * 'l', low
3495            * 'i', intermediate
3496            * 'h', high
3497            * 'f', full
3498        gtit= title of the graph ('|', for spaces)
3499        kindfig= kind of figure
3500        figuren= name of the figure
3501      ncfile= file to use
3502      varns= [lon],[lat],[outflow],[basinID] ',' list of the name of the variables with the lon,lat, the outflow and the basin ID
3503    """
3504    fname = 'draw_basins'
3505
3506    if values == 'h':
3507        print fname + '_____________________________________________________________'
3508        print draw_vectors.__doc__
3509        quit()
3510
3511    expectargs = '[X/Y/Z/T]|[dimname]|[vardimname]|[value]:[vecvals]:[windlabs]:' +  \
3512      '[mapvalues]:[gtit]:[kindfig]:[figuren]'
3513 
3514    drw.check_arguments(fname,values,expectargs,':')
3515
3516    dimvals = values.split(':')[0]
3517    vecvals = values.split(':')[1]
3518    windlabels = values.split(':')[2]
3519    mapvalues = values.split(':')[3]
3520    gtit = values.split(':')[4]
3521    kindfig = values.split(':')[5]
3522    figuren = values.split(':')[6]
3523
3524    of = NetCDFFile(ncfile,'r')
3525
3526    dims = {}
3527    for dimv in dimvals.split(','):
3528        dns = dimv.split('|')
3529        dims[dns[0]] = [dns[1], dns[2], dns[3]]
3530
3531    varNs = []
3532    for dn in dims.keys():
3533        if dn == 'X':
3534            varNs.append(dims[dn][1])
3535            dimx = len(of.dimensions[dims[dn][0]])
3536        elif dn == 'Y':
3537            varNs.append(dims[dn][1])
3538            dimy = len(of.dimensions[dims[dn][0]])
3539
3540    ivar = 0
3541    for wvar in varns.split(','):
3542        if not drw.searchInlist(of.variables.keys(), wvar):
3543            print errormsg
3544            print '  ' + fname + ": file does not have variable '" + wvar + "' !!"
3545            quit(-1)
3546        if ivar == 0:
3547            varNs.append(wvar)
3548        else:
3549            varNs.append(wvar)
3550
3551    ivar = 0
3552    for varN in varNs:
3553        varslice = []
3554
3555        ovarN = of.variables[varN]
3556        vard = ovarN.dimensions
3557        for vdn in vard:
3558            found = False
3559            for dd in dims.keys():
3560                if dims[dd][0] == vdn:
3561                    if dims[dd][2].find('@') != -1:
3562                        rvals = dims[dd][2].split('@')
3563                        varslice.append(slice(int(rvals[0]), int(rvals[1])))
3564                    elif dims[dd][2] == '-1':
3565                        varslice.append(slice(0,len(of.dimensions[dims[dd][0]])))
3566                    else:
3567                        varslice.append(int(dims[dd][2]))
3568
3569                    found = True
3570                    break
3571            if not found:
3572                varslice.append(slice(0,len(of.dimensions[dims[dd][0]])))
3573
3574        if varN == dims['X'][1]:
3575            lonvals0 = np.squeeze(ovarN[tuple(varslice)])
3576        elif varN == dims['Y'][1]:
3577            latvals0 = np.squeeze(ovarN[tuple(varslice)])
3578
3579        ivar = ivar + 1
3580
3581    if len(lonvals0.shape) == 1:
3582        lonvals, latvals = np.meshgrid(lonvals0, latvals0)
3583    else:
3584        lonvals = lonvals0
3585        latvals = latvals0
3586
3587# Vector values
3588    if vecvals.split(',')[0] == 'None':
3589        freqv = None
3590    else:
3591        freqv = vecvals.split(',')[0] 
3592
3593    colorvals = vecvals.split(',')[1]
3594    if len(varn.split(',')) != 3:
3595        print errormsg
3596        print '  ' + fname + ": color of vectors should be according to '" +         \
3597          coln + "' but a third varibale is not provided !!"
3598        quit(-1)
3599
3600    ocolvec = of.variables[varNs[3]]
3601    colorv = ocolvec[:]
3602    stdvn, lonvn, unitsvn = drw.var_3desc(ocolvec)
3603    colorvals = colorvals + '@' + stdvn + '@' + unitsvn
3604
3605    lengthv = vecvals.split(',')[2]
3606
3607# Vector labels
3608    windname = windlabels.split(',')[0]
3609    windunits = windlabels.split(',')[1]
3610
3611# Vector angles
3612    oflow = ofile.variables[varNs[2]]
3613    angle = (oflow[:] - 1)*np.pi/4
3614    xflow = np.where(oflow[:] < 9, np.float(lengthv)*np.sin(angle), 0.)
3615    yflow = np.where(oflow[:] < 9, np.float(lengthv)*np.cos(angle), 0.)
3616
3617    drw.plot_basins(lonvals, latvals, xflow, yflow, freqv, colorvals, colorv,      \
3618      lengthv, windname, windunits, mapvalues, gtit, kindfig, figuren)
3619
3620    of.close()
3621
3622    return
3623
3624def draw_river_desc(ncfile, values, riverns):
3625    """ Function to plot rivers' description from ORCHIDEE's routing scheme
3626      values= [dimname]|[vardimname]|[value]:[basinvals]:[upstreamvals]:[mapvalues]:
3627        [gtit]:[kindfig]:[legloc]:[figuren]
3628        'X/Y'|[dimname]|[vardimname]|[value]: ',', list for each basic dimension '|' separated of:
3629          [dimname]: name of the dimension in the file
3630          [vardimname]: name of the variable with the values for the dimension in the file
3631          [value]: which value of the given dimension (-1, all; [ibeg]@[iend], i-range beginning, end)
3632          No value takes all the range of the dimension
3633        [basinsvals]= [colorline]
3634          [basincolor]: ',' list of colors of the line to use to mark the basins contours (single value also possible)
3635        [upstreamvals]= [upstreamvarn],[colorbar]
3636          [upstreamcolor]: colorbar to use to plot the basins upstream values
3637        [mapvalues]= map characteristics: [proj],[res]
3638          see full documentation: http://matplotlib.org/basemap/
3639          [proj]: projection
3640            * 'cyl', cilindric
3641            * 'lcc', lambert conformal
3642          [res]: resolution:
3643            * 'c', crude
3644            * 'l', low
3645            * 'i', intermediate
3646            * 'h', high
3647            * 'f', full
3648        gtit= title of the graph ('|', for spaces)
3649        kindfig= kind of figure
3650        legloc= location of the legend (0, automatic)
3651          1: 'upper right', 2: 'upper left', 3: 'lower left', 4: 'lower right',
3652          5: 'right', 6: 'center left', 7: 'center right', 8: 'lower center',
3653          9: 'upper center', 10: 'center'      kfig= kind of figure
3654        figuren= name of the figure
3655      ncfile= file to use
3656      riverns= ',' list of the name of the rivers to plot
3657    """
3658    import numpy.ma as ma
3659    fname = 'draw_river_desc'
3660
3661    if values == 'h':
3662        print fname + '_____________________________________________________________'
3663        print draw_river_desc.__doc__
3664        quit()
3665
3666    expectargs = '[X/Y/Z/T]|[dimname]|[vardimname]|[value]:[basinvals]:' +           \
3667      '[upstreamvals]:[mapvalues]:[gtit]:[kindfig]:[legloc]:[figuren]'
3668 
3669    drw.check_arguments(fname,values,expectargs,':')
3670
3671    dimvals = values.split(':')[0]
3672    basinvals = values.split(':')[1]
3673    upstreamvals = values.split(':')[2]
3674    mapvals = values.split(':')[3]
3675    gtit = values.split(':')[4]
3676    kindfig = values.split(':')[5]
3677    legloc = int(values.split(':')[6])
3678    figuren = values.split(':')[7]
3679
3680    basincol = basinvals
3681    if basincol.find(',') != 1:
3682        basincolor = basincol.split(',')
3683    else:
3684        basincolor = [basincol]
3685
3686    upstreamcolor = upstreamvals
3687
3688    of = NetCDFFile(ncfile,'r')
3689
3690    dims = {}
3691    for dimv in dimvals.split(','):
3692        dns = dimv.split('|')
3693        dims[dns[0]] = [dns[1], dns[2], dns[3]]
3694
3695    varNs = []
3696    for dn in dims.keys():
3697        if dn == 'X':
3698            varNs.append(dims[dn][1])
3699            dimx = len(of.dimensions[dims[dn][0]])
3700        elif dn == 'Y':
3701            varNs.append(dims[dn][1])
3702            dimy = len(of.dimensions[dims[dn][0]])
3703
3704    if riverns.find(',') != -1:
3705        riverns = riverns.split(',')
3706    else:
3707        riverns = [riverns]
3708
3709    rivers = []
3710    riversubbasins = {}
3711    riversupstream = {}
3712    riversoutflow = {}
3713    for rivern in riverns:
3714        print rivern
3715
3716# subBasins
3717        basinvar = rivern + '_coding'
3718        if not drw.searchInlist(of.variables.keys(), basinvar):
3719            print errormsg
3720            print '  ' + fname + ": file does not have variable '" + basinvar + "' !!"
3721            quit(-1)
3722        rivers.append(rivern)
3723        obasin = of.variables[basinvar]
3724        riversubbasins[rivern] = obasin[:]
3725        if rivern == riverns[0]:
3726            finalmask = obasin[:].mask
3727        else:
3728            finalmask = finalmask * obasin[:].mask
3729
3730# upstream
3731        upstreamvar = rivern + '_upstream'
3732        if not drw.searchInlist(of.variables.keys(), upstreamvar):
3733            print errormsg
3734            print '  ' + fname + ": file does not have variable '" + upstreamvar + "' !!"
3735            quit(-1)
3736        oupstream = of.variables[upstreamvar]
3737        riversupstream[rivern] = oupstream[:]
3738        if rivern == riverns[0]:
3739            uunits = oupstream.getncattr('units')
3740
3741# River metadata
3742        fracvar = rivern + '_frac'
3743        if not drw.searchInlist(of.variables.keys(), fracvar):
3744            print errormsg
3745            print '  ' + fname + ": file does not have variable '" + fracvar + "' !!"
3746            quit(-1)
3747        ofrac = of.variables[fracvar]
3748        riversoutflow[rivern] = [ofrac.getncattr('Longitude_of_outflow_point'),      \
3749          ofrac.getncattr('Latitude_of_outflow_point')]
3750
3751    ivar = 0
3752    for varN in varNs:
3753        varslice = []
3754
3755        ovarN = of.variables[varN]
3756        vard = ovarN.dimensions
3757        for vdn in vard:
3758            found = False
3759            for dd in dims.keys():
3760                if dims[dd][0] == vdn:
3761                    if dims[dd][2].find('@') != -1:
3762                        rvals = dims[dd][2].split('@')
3763                        varslice.append(slice(int(rvals[0]), int(rvals[1])))
3764                    elif dims[dd][2] == '-1':
3765                        varslice.append(slice(0,len(of.dimensions[dims[dd][0]])))
3766                    else:
3767                        varslice.append(int(dims[dd][2]))
3768
3769                    found = True
3770                    break
3771            if not found:
3772                varslice.append(slice(0,len(of.dimensions[dims[dd][0]])))
3773
3774        if varN == dims['X'][1]:
3775            lonvals0 = np.squeeze(ovarN[tuple(varslice)])
3776        elif varN == dims['Y'][1]:
3777            latvals0 = np.squeeze(ovarN[tuple(varslice)])
3778
3779        ivar = ivar + 1
3780
3781    if len(lonvals0.shape) == 1:
3782        lonvals, latvals = np.meshgrid(lonvals0, latvals0)
3783    else:
3784        lonvals = lonvals0
3785        latvals = latvals0
3786
3787# Masking only the lon,lat with rivers
3788    malonvals = ma.masked_array(lonvals, mask=finalmask)
3789    malatvals = ma.masked_array(latvals, mask=finalmask)
3790
3791    if mapvals == 'None':
3792        mapvalues = None
3793    else:
3794        mapvalues = mapvals
3795
3796    drw.plot_river_desc(malonvals, malatvals, rivers, riversubbasins, riversupstream, riversoutflow,  \
3797      basincolor, upstreamcolor, uunits, mapvalues, gtit, kindfig, legloc, figuren)
3798
3799    of.close()
3800
3801def draw_vertical_levels(ncfile, values, varn):
3802    """ plotting vertical levels distribution
3803    draw_topo_geogrid(ncfile, values, varn)
3804      ncfile= file to use
3805      values= [zlogs]:[plogs]:[title]:[graphic_kind]:[legloc]
3806        zlogs= zlog,dzlog
3807        zlog: to use logarithmic scale on the height axis ('true/false')
3808        dzlog: to use logarithmic scale on the difference of height between levels axis ('true/false')
3809        plogs= plog,dplog
3810        plog: to use logarithmic scale on the height axis ('true/false')
3811        dplog: to use logarithmic scale on the difference of height between levels axis ('true/false')
3812        title: title of the graph ('!' for spaces)
3813        graphic_kind: kind of figure (jpg, pdf, png)
3814        legloc= location of the legend (0, automatic)
3815          1: 'upper right', 2: 'upper left', 3: 'lower left', 4: 'lower right',
3816          5: 'right', 6: 'center left', 7: 'center right', 8: 'lower center',
3817          9: 'upper center', 10: 'center'      kfig= kind of figure
3818      varn= [varnheight],[varnpres]
3819        varnheight: name of the variable with the height of the vertical levels
3820          'WRFz': for WRF z-levels (computed as (PH + PHB)/g, from a PHB(0,i,j) = 0)
3821        varnpres: name of the variable with the pressure of the vertical levels ('None', for no pressure plot)
3822          'WRFp': for WRF p-levels (computed as P + PB, from a PHB(0,i,j) = 0)
3823    """
3824    fname = 'draw_vertical_levels'
3825
3826    if values == 'h':
3827        print fname + '_____________________________________________________________'
3828        print draw_vertical_levels.__doc__
3829        quit()
3830
3831    expectargs = '[zlogs]:[plogs]:[title]:[graphic_kind]:[legloc]'
3832 
3833    drw.check_arguments(fname,values,expectargs,':')
3834
3835    zlog = values.split(':')[0].split(',')[0]
3836    dzlog = values.split(':')[0].split(',')[1]
3837    plog = values.split(':')[1].split(',')[0]
3838    dplog = values.split(':')[1].split(',')[1]
3839    title = values.split(':')[2].replace('!',' ')
3840    kindfig = values.split(':')[3]
3841    legloc = int(values.split(':')[4])
3842
3843    if varn.find(',') == -1:
3844        varnheight = varn
3845        varnpres = None
3846        pvals = None
3847        print warnmsg
3848        print '  ' + fname + ': assuming no pressure variable!!'
3849    else:
3850        varnheight = varn.split(',')[0]
3851        varnpres = varn.split(',')[1]
3852        if varnpres == 'None': 
3853            varnpres = None
3854            pvals = None
3855
3856    if not os.path.isfile(ncfile):
3857        print errormsg
3858        print '  ' + fname + ': file "' + ncfile + '" does not exist !!'
3859        quit(-1)   
3860
3861    objf = NetCDFFile(ncfile, 'r')
3862
3863    if varnheight == 'WRFz': 
3864        if not gen.searchInlist(objf.variables,'PH'):
3865            print errormsg
3866            print '  ' + fname + ": WRF file '" + ncfile + "' does not have " +      \
3867              "variable 'PH' !!"
3868            quit(-1)
3869        if not gen.searchInlist(objf.variables,'PHB'):
3870            print errormsg
3871            print '  ' + fname + ": WRF file '" + ncfile + "' does not have " +      \
3872              "variable 'PHB' !!"
3873            quit(-1)
3874
3875        objph = objf.variables['PH']
3876        objphb = objf.variables['PHB']
3877        geop = objph[:] + objphb[:]
3878
3879        ijz0 = gen.index_mat(geop[0,], 0.)
3880        zvals = geop[0, :, ijz0[0], ijz0[1]] / 9.8
3881    else:
3882        if not gen.searchInlist(objf.variables, varnheight):
3883            print errormsg
3884            print '  ' + fname + ": file '" + ncfile + "' does not have height " +   \
3885              " variable '" + varnheight + "' !!"
3886            quit(-1)
3887        objvar = objf.variables[varn]
3888        if len(objvar.shape) == 4:
3889            print warnmsg
3890            print '  ' + fname + ": assuming that height variable '" + varnheight +  \
3891              "' with shape: dt,dz,dy,dx. Tacking first time-step"
3892
3893            ijz0 = gen.index_mat(objvar[0,0,], 0.)
3894            zvals = objvar[0, :, ijz0[0], ijz0[1]]
3895        elif len(objvar.shape) == 3:
3896            print warnmsg
3897            print '  ' + fname + ": assuming that height variable '" + varnheight +  \
3898              "' with shape: dz,dy,dx"
3899
3900            ijz0 = gen.index_mat(objvar[0,], 0.)
3901            zvals = objvar[:, ijz0[0], ijz0[1]]
3902       
3903        elif len(objvar.shape) == 2:
3904            print warnmsg
3905            print '  ' + fname + ": assuming that height variable '" + varnheight +  \
3906              "' with shape: dz,dyx"
3907
3908            ijz0 = gen.index_mat(objvar[0,], 0.)
3909            zvals = objvar[:, ijz0[0]]
3910        else:
3911            zvals = objvar[:]
3912
3913# Pressure
3914    if varnpres is not None:
3915        if varnpres == 'WRFp': 
3916            if not gen.searchInlist(objf.variables,'P'):
3917                print errormsg
3918                print '  ' + fname + ": WRF file '" + ncfile + "' does not have " +      \
3919                  "variable 'P' !!"
3920                quit(-1)
3921            if not gen.searchInlist(objf.variables,'PB'):
3922                print errormsg
3923                print '  ' + fname + ": WRF file '" + ncfile + "' does not have " +      \
3924                  "variable 'PB' !!"
3925                quit(-1)
3926
3927            objph = objf.variables['P']
3928            objphb = objf.variables['PB']
3929            pres = objph[:] + objphb[:]
3930
3931            pvals = pres[0, :, ijz0[0], ijz0[1]]
3932        else:
3933            if not gen.searchInlist(objf.variables, varnpres):
3934                print errormsg
3935                print '  ' + fname + ": file '" + ncfile + "' does not have pressure " + \
3936                  " variable '" + varnpres + "' !!"
3937                quit(-1)
3938            objvar = objf.variables[varnpres]
3939            if len(objvar.shape) == 4:
3940                print warnmsg
3941                print '  ' + fname + ": assuming that pressure variable '" + varnpres +  \
3942                  "' with shape: dt,dz,dy,dx. Tacking first time-step"
3943   
3944                pvals = objvar[0, :, ijz0[0], ijz0[1]]
3945            elif len(objvar.shape) == 3:
3946                print warnmsg
3947                print '  ' + fname + ": assuming that pressure variable '" + varnpres +  \
3948                  "' with shape: dz,dy,dx"
3949   
3950                pvals = objvar[:, ijz0[0], ijz0[1]]
3951           
3952            elif len(objvar.shape) == 2:
3953                print warnmsg
3954                print '  ' + fname + ": assuming that pressure variable '" + varnpres +  \
3955                  "' with shape: dz,dyx"
3956   
3957                pvals = objvar[:, ijz0[0]]
3958            else:
3959                pvals = objvar[:]
3960
3961# Logarithmic axes
3962    if zlog == 'true':
3963        zlogv = True
3964    elif zlog == 'false':
3965        zlogv = False
3966    else:
3967        print errormsg
3968        print '  ' + fname + ": wrong value for zlog: '" + zlog + "' !!"
3969        print "    must be either: 'true' or 'false'"
3970        quit(-1)
3971
3972    if dzlog == 'true':
3973        dzlogv = True
3974    elif dzlog == 'false':
3975        dzlogv = False
3976    else:
3977        print errormsg
3978        print '  ' + fname + ": wrong value for dzlog: '" + dzlog + "' !!"
3979        print "    must be either: 'true' or 'false'"
3980        quit(-1)
3981
3982    if pvals is not None:
3983        if plog == 'true':
3984            plogv = True
3985        elif plog == 'false':
3986            plogv = False
3987        else:
3988            print errormsg
3989            print '  ' + fname + ": wrong value for plog: '" + plog + "' !!"
3990            print "    must be either: 'true' or 'false'"
3991            quit(-1)
3992        if dplog == 'true':
3993            dplogv = True
3994        elif dplog == 'false':
3995            dplogv = False
3996        else:
3997            print errormsg
3998            print '  ' + fname + ": wrong value for dplog: '" + dplog + "' !!"
3999            print "    must be either: 'true' or 'false'"
4000            quit(-1)
4001
4002    drw.plot_vertical_lev(zvals, pvals, zlogv, dzlogv, plogv, dplogv, title, kindfig, legloc)
4003
4004    objf.close()
4005
4006    return
4007
4008def draw_subbasin(ncfile, values):
4009    """ Function to plot subbasin from 'routnig.nc' ORCDHIEE
4010      ncfile= file to use produced with nc_var.py#subbasin function
4011      values= [subasiname]:[rangecolors]:[mapv]:[basinlinewidth]:[drawsubid]:[gtit]:[figkind]:[legloc]:[figurename]
4012        [subasiname]= name of the subbasin ('!' for spaces)
4013        [rcolor]= '@', list of 'r|g|b' 1-based colors (as much as first level sub-flow). 'None' for automatic
4014        [mapv]= map characteristics: [proj],[res]
4015          see full documentation: http://matplotlib.org/basemap/
4016          [proj]: projection
4017            * 'cyl', cilindric
4018            * 'lcc', lambert conformal
4019          [res]: resolution:
4020            * 'c', crude
4021            * 'l', low
4022            * 'i', intermediate
4023            * 'h', high
4024            * 'f', full
4025        [basinlinewidth]= with of the line to draw the basin
4026        [drawsubid]= wehther sub-flow ids should be plot or not
4027        [graphtit]= title of the graph ('|', for spaces)
4028        [lloc]= location of the legend (0, automatic)
4029          1: 'upper right', 2: 'upper left', 3: 'lower left', 4: 'lower right',
4030          5: 'right', 6: 'center left', 7: 'center right', 8: 'lower center',
4031          9: 'upper center', 10: 'center'      kfig= kind of figure
4032        [figname]= name of the figure
4033    """
4034    fname = 'draw_subbasin'
4035
4036    if values == 'h':
4037        print fname + '_____________________________________________________________'
4038        print draw_subbasin.__doc__
4039        quit()
4040
4041    expectargs = '[subasiname]:[rangecolors]:[mapv]:[basinlinewidth]:[drawsubid]:' + \
4042      '[gtit]:[figkind]:[legloc]:[figurename]'
4043 
4044    drw.check_arguments(fname,values,expectargs,':')
4045
4046    subbasiname = values.split(':')[0].replace('!',' ')
4047    rangecolors = values.split(':')[1]
4048    mapv = values.split(':')[2]
4049    basinlinewidth = np.float(values.split(':')[3])
4050    drawsubid = gen.Str_Bool(values.split(':')[4])
4051    gtit = values.split(':')[5].replace('!',' ')
4052    figkind = values.split(':')[6]
4053    legloc = int(values.split(':')[7])
4054    figurename = values.split(':')[8]
4055
4056    if not os.path.isfile(ncfile):
4057        print errormsg
4058        print '  ' + fname + ': file "' + ncfile + '" does not exist !!'
4059        quit(-1)   
4060
4061    objf = NetCDFFile(ncfile, 'r')
4062
4063    searchvars = ['lon', 'lat', 'lonsubflow', 'latsubflow', 'outsubflow']
4064    for searchvar in searchvars: 
4065        if not gen.searchInlist(objf.variables,searchvar):
4066            print errormsg
4067            print '  ' + fname + ": WRF file '" + ncfile + "' does not have " +      \
4068              "variable '" + searchvar + "' !!"
4069            quit(-1)
4070   
4071# lon,lat
4072    olon = objf.variables['lon']
4073    olat = objf.variables['lat']
4074    lon = olon[:]
4075    lat = olat[:]
4076
4077# sub-flow names
4078    osubnames = objf.variables['subflow']
4079    subnames = drw.get_str_nc(osubnames, osubnames.shape[1])
4080
4081# sub-flow lat, lon
4082    latlonsub = {}
4083    outflowsub = {}
4084    osublon = objf.variables['lonsubflow']
4085    osublat = objf.variables['latsubflow']
4086    oNsubflow = objf.variables['Nsubflow']
4087    ooutsubflow = objf.variables['outsubflow']
4088    Nsubflow = oNsubflow[:]
4089    isub = 0
4090    for Ssub in subnames:
4091        sublatlon = []
4092        suboutflow = []
4093        for igrid in range(Nsubflow[isub]):
4094            sublatlon.append([osublat[isub,igrid], osublon[isub,igrid]])
4095            suboutflow.append(ooutsubflow[isub,igrid])
4096        latlonsub[Ssub] = sublatlon
4097        outflowsub[Ssub] = suboutflow
4098        isub = isub + 1
4099
4100# colors
4101    if rangecolors == 'None':
4102        rangecols = None
4103    else:
4104        cols = rangecolors.split('@')
4105        Ncols = len(cols)
4106        rangecols = []
4107        for icol in range(Ncols):
4108            cval = cols[icol].split('|')
4109            rangecols.append([np.float(cval[0]),np.float(cval[1]),np.float(cval[2])])
4110
4111    drw.plot_subbasin(subbasiname, lon, lat, subnames, latlonsub, outflowsub,        \
4112      rangecols, mapv, basinlinewidth, drawsubid, gtit, figkind, legloc, figurename)
4113
4114    objf.close()
4115
4116    return
4117
4118def draw_2lines(values,ncfiles,varnames):
4119    """ Fucntion to plot two lines in different axes (x/x2 or y/y2)
4120      values= [varnA]:[varnB]:[commonvardim]:[varangeA]:[varangeB]:[varangeaxis]:[axisvals]:[figvarns]:[colors]:
4121       [widths]:[styles]:[sizemarks]:[marks]:[graphtitle]:[labelaxis]:[lloc]:[figname]:[figkind]
4122        [varnA]: name of the variable A in fileA
4123        [varnB]: name of the variable B in fileB
4124        [commonvardim]: name of the common variable-dimension
4125        [varangeA]: ',' separated list of range (min,max) for A values ('None', automatic; 'Extrs' from values extremes)
4126        [varangeB]: ',' separated list of range (min,max) for B values ('None', automatic; 'Extrs' from values extremes)
4127        [varangeaxis]: ',' separated list of range (min,max) for common axis values ('None', automatic; 'Extrs' from
4128          values extremes)
4129        [axisvals]: which is the axis to plot the values ('x' or 'y')
4130        [figvarns]: ',' separated list of names of the variables in the  plot
4131        [colors]: ',' list with color names of the lines for the variables ('None', automatic)
4132        [widths]: ',' list with widths of the lines for the variables ('None', automatic)
4133        [styles]: ',' list with the styles of the lines ('None', automatic)
4134        [sizemarks]: ',' list with the size of the markers of the lines ('None', automatic)
4135        [marks]: ',' list with the markers of the lines ('None', automatic)
4136        [graphtitle]: title of the figure ('!' for spaces)
4137        [labelaxis]: label in the figure of the common axis ('!' for spaces)
4138        [lloc]: location of the legend (0, automatic)
4139          1: 'upper right', 2: 'upper left', 3: 'lower left', 4: 'lower right',
4140          5: 'right', 6: 'center left', 7: 'center right', 8: 'lower center',
4141          9: 'upper center', 10: 'center'      kfig= kind of figure
4142        [figname]: name of the figure
4143        [figkind]: kind of figure
4144      ncfiles= ',' separated list of files to use
4145      varnames=  ',' separated list of variables names in the files to plot
4146    """
4147    function = 'draw_2lines'
4148
4149    if values == 'h':
4150        print fname + '_____________________________________________________________'
4151        print draw_2lines.__doc__
4152        quit()
4153
4154    expectargs = '[varnA]:[varnB]:[commonvardim]:[varangeA]:[varangeB]:' +           \
4155     '[varangeaxis]:[axisvals]:[figvarns]:[colors]:[widths]:[styles]:[sizemarks]:' + \
4156     '[marks]:[graphtitle]:[labelaxis]:[lloc]:[figname]:[figkind]'
4157 
4158    drw.check_arguments(fname,values,expectargs,':')
4159
4160    varns = [values.split(':')[0], values.split(':')[1]]
4161    commonvardim = values.split(':')[2]
4162    varangeA0 = values.split(':')[3]
4163    varangeB0 = values.split(':')[4]
4164    varangeaxis0 = values.split(':')[5]
4165    axisvals = values.split(':')[6]
4166    figvarns = values.split(':')[7].split(',')
4167    colors = gen.str_list(values.split(':')[8],',')
4168    widths = gen.str_list(values.split(':')[9],',')
4169    styles = gen.str_list(values.split(':')[10],',')
4170    sizemarks = gen.str_list(values.split(':')[11],',')
4171    marks = gen.str_list(values.split(':')[12],',')
4172    graphtitle = values.split(':')[13].replace('!',' ')
4173    labelaxis = values.split(':')[14].replace('!',' ')
4174    lloc = np.int(values.split(':')[15])
4175    figname = values.split(':')[16]
4176    figkind = values.split(':')[17]
4177
4178    files = ncfiles.split(',')
4179    invarns = varnames.split(',')
4180
4181    varunits = []
4182
4183    # Values line A
4184    if not os.path.isfile(files[0]):
4185        print errormsg
4186        print '  ' + fname + ": file '" + files[0] + "' does not exist !!"
4187        quit(-1)
4188    oncA = NetCDFFile(files[0], 'r')
4189
4190    if not gen.searchInfile(oncA.variables.kesy(), invarns[0]):
4191        print errormsg
4192        print '  ' + fname + ": A file '" + files[0] + "' does not have variable '" +\
4193          invarns[0] + "' !!"
4194        quit(-1)
4195    objvA = oncA.variables[invarns[0]]
4196    varvalsA = objvA[:]
4197    varangeA = []
4198
4199    if gen.searchInfile(objvA.ncattrs(), 'units'):
4200        varunits.append(objva.getncattr('units'))
4201    else:
4202        valsA = gen.variables_values(invarns[0])
4203        varunits.append(valsA[5])
4204        if varangeA0 == 'None':
4205            varangeA = [valsA[2], valsA[3]]
4206        elif varangeA0 == 'Extrs':
4207            varangeA = [np.min(varvalsA), np.max(varvalsA)]
4208        else:
4209            for iv in range(2): varangeA[iv] = np.float(varangeA0[iv])
4210
4211    if not gen.searchInfile(oncA.variables.kesy(), commvardim):
4212        print errormsg
4213        print '  ' + fname + ": A file '" + files[0] + "' does not have common " +   \
4214          "dimvar '" + commonvardim + "' !!"
4215        quit(-1)
4216    objvd = oncA.variables[commonvardim]
4217    varvalsaxisA = objvd[:]   
4218
4219    oncA.close()
4220
4221    # Values line B
4222    if not os.path.isfile(files[1]):
4223        print errormsg
4224        print '  ' + fname + ": file '" + files[1] + "' does not exist !!"
4225        quit(-1)
4226    oncB = NetCDFFile(files[1], 'r')
4227
4228    if not gen.searchInfile(oncB.variables.keys(), invarns[1]):
4229        print errormsg
4230        print '  ' + fname + ": B file '" + files[1] + "' does not have variable '" +\
4231          invarns[1] + "' !!"
4232        quit(-1)
4233    objvB = oncB.variables[invarns[1]]
4234    varvalsB = objvB[:]
4235    varangeB = []
4236
4237    if gen.searchInfile(objvB.ncattrs(), 'units'):
4238        varunits.append(objvB.getncattr('units'))
4239    else:
4240        valsB = gen.variables_values(invarns[1])
4241        varunits.append(valsB[5])
4242        if varangeB0 == 'None':
4243            varangeB = [valsB[2], valsB[3]]
4244        elif varangeB0 == 'Extrs':
4245            varangeA = [np.min(varvalsA), np.max(varvalsA)]
4246        else:
4247            for iv in range(2): varangeB[iv] = np.float(varangeB0[iv])
4248
4249    if not gen.searchInfile(oncB.variables.keys(), commvardim):
4250        print errormsg
4251        print '  ' + fname + ": B file '" + files[1] + "' does not have common " +   \
4252          "dimvar '" + commonvardim + "' !!"
4253        quit(-1)
4254    objvd = oncB.variables[commonvardim]
4255    varvalsaxisB = objvd[:]
4256
4257    # Range of the axis
4258    varangeaxis = []
4259    if gen.searchInfile(objvd.ncattrs(), 'units'):
4260        dimvarunits.append(objvd.getncattr('units'))
4261    else:
4262        valsVD = gen.variables_values(commonvardim)
4263        varunits.append(valsVD[5])
4264        if varangeaxis0 == 'None':
4265            varangeaxis = [valsVD[2], valsVD[3]]
4266        elif varangeaxis0 == 'Extrs':
4267            varangeaxis[0] = np.min([np.min(varvalsaxisA), np.min(varvalsaxisB)])
4268            varangeaxis[1] = np.max([np.max(varvalsaxisA), np.max(varvalsaxisB)])
4269        else:
4270            for iv in range(2): varangeaxis0[iv] = np.float(varangeaxis[iv])
4271   
4272    oncB.close()
4273
4274    plot_2lines(varvalsA, varvalsB, varvalsaxisA, varvalsaxisB, varangeA, varrangeB, \
4275      varangeaxis, axisvals, varns, varunits, colors, widths, styles, sizemarks,     \
4276      marks, graphtitle, labelaxis, lloc, figname, figkind):
4277
4278#quit()
4279
4280####### ###### ##### #### ### ## #
4281
4282ngraphics = "'" + drw.numVector_String(namegraphics, "', '") + "'"
4283 
4284### Options
4285##string_operation="operation to make: " + '\n' + " out, output values -S inidim1,[inidim2,...]:enddim1,[enddim2,...]"
4286string_operation="""operation to make:
4287  draw_topo_geogrid, draws topography from a WPS geo_em.d[nn].nc: -S [minTopo],[maxTopo]:[SW_lon],[SW_lat],[NE_lon],[NE_lat]:[title]:[graphic_kind]:[projection],[res_coastline]
4288  draw_2D_shad_cont, draws two 2D fields, first with shading second with contour lines: -v [varns],[varnc] -S [vnamefs],[vnamefc],[dimxvn],[dimyvn],[colorbar],[ckind],[clabfmt],[sminv]:[smaxv],[sminc]:[smaxv]:[Nlev],[figt],[kindfig],[reverse]
4289    [ckind]:
4290      'cmap': as it gets from colorbar
4291      'fixc,[colname]': fixed color [colname], all stright lines
4292      'fixsignc,[colname]': fixed color [colname], >0 stright, <0 dashed  line
4293"""
4294
4295#print string_operation
4296
4297parser = OptionParser()
4298parser.add_option("-f", "--netCDF_file", dest="ncfile", 
4299                  help="file to use", metavar="FILE")
4300parser.add_option("-o", "--operation", type='choice', dest="operation", 
4301       choices=namegraphics, 
4302                  help="operation to make: " + ngraphics, metavar="OPER")
4303parser.add_option("-S", "--valueS", dest="values", 
4304                  help="[WHEN APPLICABLE] values to use according to the operation", metavar="VALUES")
4305parser.add_option("-v", "--variable", dest="varname",
4306                  help="[WHEN APPLICABLE] variable to check", metavar="VAR")
4307
4308(opts, args) = parser.parse_args()
4309
4310#######    #######
4311## MAIN
4312    #######
4313
4314# Not checking file operation
4315Notcheckingfile = ['draw_2D_shad_cont', 'draw_2D_shad_cont_time',                    \
4316  'draw_2D_shad_line', 'draw_2D_shad_line_time', 'draw_lines', 'draw_lines_time',    \
4317  'draw_points', 'draw_topo_geogrid_boxes', 'draw_trajectories',                     \
4318  'draw_vals_trajectories', 'variable_values']
4319
4320####### ###### ##### #### ### ## #
4321errormsg='ERROR -- error -- ERROR -- error'
4322
4323varn=opts.varname
4324oper=opts.operation
4325
4326if opts.ncfile is not None and not os.path.isfile(opts.ncfile) and                   \
4327  not drw.searchInlist(Notcheckingfile, oper):
4328    print errormsg
4329    print '  ' + main + ': File ' + opts.ncfile + ' does not exist !!'
4330    quit(-1)
4331
4332if oper == 'create_movie':
4333    create_movie(opts.ncfile, opts.values, opts.varname)
4334elif oper == 'draw_2D_shad':
4335    draw_2D_shad(opts.ncfile, opts.values, opts.varname)
4336elif oper == 'draw_2D_shad_time':
4337    draw_2D_shad_time(opts.ncfile, opts.values, opts.varname)
4338elif oper == 'draw_2D_shad_cont':
4339    draw_2D_shad_cont(opts.ncfile, opts.values, opts.varname)
4340elif oper == 'draw_2D_shad_cont_time':
4341    draw_2D_shad_cont_time(opts.ncfile, opts.values, opts.varname)
4342elif oper == 'draw_2D_shad_line':
4343    draw_2D_shad_line(opts.ncfile, opts.values, opts.varname)
4344elif oper == 'draw_2D_shad_line_time':
4345    draw_2D_shad_line_time(opts.ncfile, opts.values, opts.varname)
4346elif oper == 'draw_barbs':
4347    draw_barbs(opts.ncfile, opts.values, opts.varname)
4348elif oper == 'draw_basins':
4349    draw_basins(opts.ncfile, opts.values, opts.varname)
4350elif oper == 'draw_Neighbourghood_evol':
4351    draw_Neighbourghood_evol(opts.ncfile, opts.values, opts.varname)
4352elif oper == 'draw_lines':
4353    draw_lines(opts.ncfile, opts.values, opts.varname)
4354elif oper == 'draw_lines_time':
4355    draw_lines_time(opts.ncfile, opts.values, opts.varname)
4356elif oper == 'draw_points':
4357    draw_points(opts.ncfile, opts.values)
4358elif oper == 'draw_points_lonlat':
4359    draw_points_lonlat(opts.ncfile, opts.values)
4360elif oper == 'draw_ptZvals':
4361    draw_ptZvals(opts.ncfile, opts.values, opts.varname)
4362elif oper == 'draw_river_desc':
4363    draw_river_desc(opts.ncfile, opts.values, opts.varname)
4364elif oper == 'draw_subbasin':
4365    draw_subbasin(opts.ncfile, opts.values)
4366elif oper == 'draw_timeSeries':
4367    draw_timeSeries(opts.ncfile, opts.values, opts.varname)
4368elif oper == 'draw_topo_geogrid':
4369    draw_topo_geogrid(opts.ncfile, opts.values)
4370elif oper == 'draw_topo_geogrid_boxes':
4371    draw_topo_geogrid_boxes(opts.ncfile, opts.values)
4372elif oper == 'draw_trajectories':
4373    draw_trajectories(opts.ncfile, opts.values, opts.varname)
4374elif oper == 'draw_vals_trajectories':
4375    draw_vals_trajectories(opts.ncfile, opts.values, opts.varname)
4376elif oper == 'draw_vectors':
4377    draw_vectors(opts.ncfile, opts.values, opts.varname)
4378elif oper == 'draw_vertical_levels':
4379    draw_vertical_levels(opts.ncfile, opts.values, opts.varname)
4380elif oper == 'list_graphics':
4381# From: http://www.diveintopython.net/power_of_introspection/all_together.html
4382    import drawing as myself
4383    object = myself
4384    for opern in namegraphics:
4385        if  opern != 'list_graphics': 
4386            print opern + '_______ ______ _____ ____ ___ __ _'
4387            print getattr(object, opern).__doc__
4388elif oper == 'variable_values':
4389    variable_values(opts.values)
4390else:
4391    print errormsg
4392    print '  ' + main + ": the graphic '" + oper + "' is not ready !!"
4393    print errormsg
4394    quit()
Note: See TracBrowser for help on using the repository browser.