Context Navigation

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

Last change on this file since 1263 was 1263, checked in by lfita, 8 years ago
Adding the use of `dimxyfmt' on 'Nlines'
File size: 204.0 KB

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

Note: See TracBrowser for help on using the repository browser.

Download in other formats: