Context Navigation

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

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

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

Download in other formats: