Context Navigation

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

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

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

Download in other formats: