Context Navigation

source: lmdz_wrf/trunk/tools/model_graphics.py @ 1165

Last change on this file since 1165 was 1160, checked in by lfita, 9 years ago
Seem to be fixed all the issues
File size: 179.1 KB

Line
1	# Python script to generate the files and the plots for a sensitivity study with multiple modle configurations and models
2	# L. Fita
3	# LMD, Jussieu/Palaisseau, France
4	# Script configuration get from ASCII file 'model_graphics.dat'
5
6	import numpy as np
7	import nc_var_tools as ncvar
8	import generic_tools as gen
9	import drawing_tools as drw
10	import time as tim
11	# To avoid errors like: '/bin/sh: 1: Syntax error: Bad fd number'
12	# Make sure that /bin/sh directs to bash and not to dash:
13	# http://stackoverflow.com/questions/15809060/sh-syntax-error-bad-fd-number
14	import subprocess as sub
15	import os
16
17	main = 'model_graphics.py'
18
19	errmsg = ncvar.errormsg
20	warnmsg = ncvar.warnmsg
21
22	def verify_configuration(ExpConf, debug):
23	""" Function to verify mandatory keys from experiment configuration
24	ExpConf= dictionary with the configuration
25	"""
26	fname = 'verify_configuration'
27
28	# List with the mandatory keys
29	mandatorykeys = ['pyHOME', 'cdoHOME', 'scratch', 'filescratch', 'figscratch', \
30	'diffscratch', 'figdiffscratch', 'figallmodexpscratch', 'addfiles', \
31	'addfigures', 'adddiffs', 'adddifffigures', 'addallmodexpfigures', 'debug', \
32	'models', 'modgraphchar', 'expgraphchar', 'ifold', 'ofold', 'warnmsg', \
33	'errmsg', 'titleoperations', 'kindfig', 'CFreftime', 'CFunitstime', 'mapval', \
34	'timekind', 'timefmt', 'timelabel' ]
35
36	# List with the optional keys
37	optionalkeys = ['specificvarplot', 'specificdiffopplot', 'specificdiffvarplot', \
38	'RefProj']
39
40	# Dictionary with the mandatory keys which are required as function of a given
41	# key from the dictionary
42	# if [key] must exist [keyB] = [values]
43	ifkeys = {'WRFexps': ['models','WRF'], 'LMDZexps': ['models','LMDZ'], \
44	'WRF_LMDZexps': ['models','WRF_LMDZ'], 'WRFheaders': ['models','WRF'], \
45	'WRF_LMDZheaders': ['models','WRF_LMDZ'], 'LMDZheaders': ['models','LMDZ']}
46
47	# Dictionary with the optional keys (at least one) which are required as function
48	# of a given key from the dictionary
49	# if [key] must exist at least one of the [keyA, keyB, keyC, ...] with values
50	ifopkeys = {'RefProj': ['CDOremapvar_dis bil', 'CDOremapvar_bic', \
51	'CDOremapvar_dis', 'CDOremapvar_nn', 'CDOremapvar_con', 'CDOremapvar_con2', \
52	'CDOremapvar_laf']}
53
54	for kdict in mandatorykeys:
55	if not ExpConf.has_key(kdict):
56	print errmsg
57	print ' ' + fname + "' configuration without required '" + kdict + "' !!"
58	quit(-1)
59	if debug:
60	print ExpConf[kdict]
61
62	for kdict in ifkeys.keys():
63	vals = ifkeys[kdict]
64	keyn = vals[0]
65	keyv = vals[1]
66	if type(ExpConf[keyn]) == type('A'):
67	if not ExpConf.has_key(kdict) and ExpConf[keyn] == keyv:
68	print errmsg
69	print ' ' + fname + "' configuration without '" + kdict + \
70	"' when it is required with configuration '" + keyn + '=' + keyv + \
71	"' !!"
72	quit(-1)
73	elif type(ExpConf[keyn]) == type(['A', 'B']):
74	keyvals = ExpConf[keyn]
75	if not ExpConf.has_key(kdict) and gen.searchInlist(keyvals, keyv):
76	print errmsg
77	print ' ' + fname + "' configuration without '" + kdict + \
78	"' when it is required with configuration '" + keyn + '= [..., ' + \
79	keyv + ", ...]' !!"
80	quit(-1)
81	else:
82	print errmsg
83	print ' ' +fname+ 'dictionary type ', type(ExpConf[keyn]), ' not ready!!'
84	quit(-1)
85	if debug:
86	print ExpConf[kdict]
87
88	for kdict in ifopkeys.keys():
89	if ExpConf.has_key(kdict):
90	opkeys = ifopkeys[kdict]
91	hasopkey = False
92	opdictvals = {}
93	for opkv in opkeys:
94	if ExpConf.has_key(opkv):
95	hasopkey = True
96	opdictvals[opkv] = ExpConf[opkv]
97
98	if not hasopkey:
99	print errmsg
100	print ' ' + fname + "' configuration wit '" + kdict + "' but " + \
101	"without any key value:", opkeys, "' when it is required !!"
102	quit(-1)
103
104	if debug:
105	print " Optional key '" + kdict + "' with values ________"
106	gen.printing_dictionary(opdictvals)
107
108	print fname + ': configuration seems to be fine'
109
110	return
111
112	def scratches(config):
113	""" Function to set-up if it is needed to start from the scratch
114	config: dictionary with the configuration
115	"""
116	fname = 'scratches'
117
118	if config['scratch'] == 'true':
119	scr = True
120	print warnmsg
121	print " " + fname + ": starting from the SCRATCH !!"
122	print " 10 seconds left!!"
123	filescr = True
124	figscr = True
125	difscr = True
126	figdifscr = True
127	moddifscr = True
128	figmoddifscr = True
129	figallmodexpscr = True
130	tim.sleep(10)
131	else:
132	scr = False
133	if config['filescratch'] == 'true':
134	filescr = True
135	print warnmsg
136	print " " + main + ": files starting from the SCRATCH !!"
137	print " 5 seconds left!!"
138	tim.sleep(5)
139	else:
140	filescr = False
141
142	if config['figscratch'] == 'true':
143	figscr = True
144	print warnmsg
145	print " " + main + ": figures starting from the SCRATCH !!"
146	print " 5 seconds left!!"
147	tim.sleep(5)
148	else:
149	figscr = False
150
151	if config['diffscratch'] == 'true':
152	difscr = True
153	print warnmsg
154	print " " + main + ": experiment differences starting from the SCRATCH !!"
155	print " 5 seconds left!!"
156	tim.sleep(5)
157	else:
158	difscr = False
159
160	if config['figdiffscratch'] == 'true':
161	figdifscr = True
162	print warnmsg
163	print " " + main + ": figures experiment differences starting from the SCRATCH !!"
164	print " 5 seconds left!!"
165	tim.sleep(5)
166	else:
167	figdifscr = False
168
169	if config['moddiffscratch'] == 'true':
170	moddifscr = True
171	print warnmsg
172	print " " + main + ": model differences starting from the SCRATCH !!"
173	print " 5 seconds left!!"
174	tim.sleep(5)
175	else:
176	moddifscr = False
177
178	if config['figmoddiffscratch'] == 'true':
179	figmoddifscr = True
180	print warnmsg
181	print " " + main + ": figures model differences starting from the SCRATCH !!"
182	print " 5 seconds left!!"
183	tim.sleep(5)
184	else:
185	figmoddifscr = False
186
187	if config['figallmodexpscratch'] == 'true':
188	figallmodexpscr = True
189	print warnmsg
190	print " " + main + ": figures all model-experiment starting from the SCRATCH !!"
191	print " 5 seconds left!!"
192	tim.sleep(5)
193	else:
194	figallmodexpscr = False
195
196	if config['addall'] == 'true':
197	print warnmsg
198	print " " + fname + ": adding new values everywhere !!"
199	addfils = True
200	addfigs = True
201	adddiffs = True
202	adddifffigs = True
203	addmoddiffs = True
204	addmoddifffigs = True
205	addallmodexpfigs = True
206	else:
207	if config['addfiles'] == 'true':
208	addfils = True
209	else:
210	addfils = False
211
212	if config['addfigures'] == 'true':
213	addfigs = True
214	else:
215	addfigs = False
216
217	if config['adddiffs'] == 'true':
218	adddiffs = True
219	else:
220	adddiffs = False
221
222	if config['adddifffigures'] == 'true':
223	adddifffigs = True
224	else:
225	adddifffigs = False
226
227	if config['addmoddiffs'] == 'true':
228	addmoddiffs = True
229	else:
230	addmoddiffs = False
231
232	if config['addmoddifffigures'] == 'true':
233	addmoddifffigs = True
234	else:
235	addmoddifffigs = False
236
237	if config['addallmodexpfigures'] == 'true':
238	addallmodexpfigs = True
239	else:
240	addallmodexpfigs = False
241
242	if config['debug'] == 'true':
243	debug = True
244	else:
245	debug = False
246
247	return scr, filescr, figscr, difscr, figdifscr, moddifscr, figmoddifscr, \
248	figallmodexpscr, addfils, addfigs, adddiffs, adddifffigs, addmoddiffs, \
249	addmoddifffigs, addallmodexpfigs, debug
250
251
252	class gc(object):
253	""" Class which holds all the graphical characteristics
254	[label]: label as it should appear in the graphics
255	[color]: specific color of the model
256	[linetype]: specific type of line of the model
257	[marker]: specific marker of the model
258	[sizes]: line width and point size for the model
259	[tmod]: time modification to apply to the model files (None for nothing)
260	'setorigin',[YYYYMMDDHHMISS]: re-set origin of times at [YYYYMMDDHHMISS]
261	"""
262	def __init__( self, label, color, linetype, marker, sizes, Tmod):
263	self.label = None
264	self.color = None
265	self.linetype = None
266	self.marker = None
267	self.sizes = None
268	self.tmod = None
269	if label is not None:
270	self.label = label
271	self.color = color
272	self.linetype = linetype
273	self.marker = marker
274	self.sizes = sizes
275	if Tmod != 'None':
276	self.tmod = Tmod
277
278	def graphical_characteristics(config, debug):
279	""" Function to provide the graphical characteristics of models and experiments
280	config= configuration of the experiment
281	"""
282	fname = 'graphical_characteristics'
283
284	modgraphchar = {}
285	modvals = config['modgraphchar'].split(':')
286	for modval in modvals:
287	modv = modval.split('\|')
288	modgraphchar[modv[0]] = gc(modv[1], modv[2], modv[3], modv[4], modv[5], modv[6])
289
290	expgraphchar = {}
291	expvals = config['expgraphchar'].split(':')
292	for expval in expvals:
293	expvs = expval.split('\|')
294	modexp = expvs[0] + '/' + expvs[1]
295	if not modgraphchar.has_key(expvs[0]):
296	print errmsg
297	print ' ' + fname + ": a model called '" +expvs[0]+ "' does not exist !!"
298	print ' existing models with graphic characteristics:', \
299	modgraphchar.keys()
300	quit(-1)
301
302	modvs = modgraphchar[expvs[0]]
303	expv = expvs[2:7]
304	if expv[1] == 'asmodel': expv[1] = modvs.color
305	if expv[2] == 'asmodel': expv[2] = modvs.linetype
306	if expv[3] == 'asmodel': expv[3] = modvs.marker
307	if expv[4] == 'asmodel': expv[4] = modvs.sizes
308
309	expgraphchar[modexp] = gc(expv[0], expv[1], expv[2], expv[3], expv[4], 'None')
310
311	# Running some tests
312	mods = config['models'].split(':')
313	for mod in mods:
314	if not modgraphchar.has_key(mod):
315	print errmsg
316	print ' ' + fname + ": model called '" + mod + "' does not have " + \
317	"graphical characteristics !!"
318	print ' existing models with graphic characteristics:', \
319	modgraphchar.keys()
320	quit(-1)
321
322	exps = config[mod + 'exps'].split(':')
323	for exp in exps:
324	modexp = mod + '/' + exp
325	if not expgraphchar.has_key(modexp):
326	print errmsg
327	print ' ' + fname + ": model/experiment called '" + modexp + \
328	"' does not have graphical characteristics !!"
329	print ' existing model/experiments with graphic ' + \
330	'characteristics:', expgraphchar.keys()
331	quit(-1)
332
333	if debug:
334	print ' ' + fname + ': model graphical characteristics _______'
335	for modn in modgraphchar.keys():
336	with gen.Capturing() as output:
337	gen.printing_class(modgraphchar[modn])
338	print "'" + modn + "'; " + ', '.join(output)
339	print ' ' + fname + ': experiment graphical characteristics _______'
340	for modexpn in expgraphchar.keys():
341	with gen.Capturing() as output:
342	gen.printing_class(expgraphchar[modexpn])
343	print "'" + modexpn + "'; " + ', '.join(output)
344
345	return modgraphchar, expgraphchar
346
347	def exp_headers(mod,config):
348	""" Function to provide the headers and the experiments of a given model
349	model= model
350	config= configuration of the experiment
351	"""
352	fname = 'exp_headers'
353
354	# No CASE in python!
355	# case ${mod} in ...
356	if mod == 'WRF':
357	expers = config['WRFexps'].split(':')
358	fhs = config['WRFheaders'].split(':')
359	elif mod == 'LMDZ':
360	expers = config['LMDZexps'].split(':')
361	fhs = config['LMDZheaders'].split(':')
362	elif mod == 'WRF_LMDZ':
363	expers = config['WRF_LMDZexps'].split(':')
364	fhs = config['WRF_LMDZheaders'].split(':')
365	else:
366	print errmsg
367	print " " + fname + ": model '" + mod + "' not ready!!"
368	quit(-1)
369
370	return expers, fhs
371
372	class VariableInf(object):
373	""" Class which holds all information of a given variable
374	name= CF name of the variable
375	header= header of the file from which it can be computed
376	model= file-header variable from which it can be directly computed
377	diag= file-header diagnostics from which it can be computed
378	"""
379	def __init__( self, name, fheader, model, diag):
380	self.name= None
381	self.fheader= None
382	self.model= None
383	self.diag= None
384	if name is not None:
385	self.name = name
386	self.fheader= fheader
387	self.model= model
388	self.diag= diag
389
390	class ModelInf(object):
391	""" Class which holds all information from a given model
392	name= Acronym of the model
393	model= Long description of the model
394	dim[x/y/z/t/s]n= name of the dimensions of the model
395	vdim[x/y/z/t/s]n= name of the variable-dimensions of the model
396	dimensions= list of the dimensions
397	vardimensions= list of the variable-dimensions
398	timemodif= time modification to apply to the model files ('None' for nothing)
399	'setorigin',[YYYYMMDDHHMISS]: re-set origin of times at [YYYYMMDDHHMISS]
400	"""
401	def __init__( self, name, model, dx, dy, dz, dt, ds, vdx, vdy, vdz, vdt, vds, \
402	tmod):
403	self.name= None
404	self.model= None
405	self.dimxn= None
406	self.dimyn= None
407	self.dimzn= None
408	self.dimtn= None
409	self.dimsn= None
410	self.vardxn= None
411	self.vardyn= None
412	self.vardzn= None
413	self.vardtn= None
414	self.vardsn= None
415	self.timemodif = None
416	if name is not None:
417	self.name = name
418	self.model= model
419	self.dimxn= dx
420	self.dimyn= dy
421	self.dimzn= dz
422	self.dimtn= dt
423	self.dimsn= ds
424	self.vardxn= vdx
425	self.vardyn= vdy
426	self.vardzn= vdz
427	self.vardtn= vdt
428	self.vardsn= vds
429	if tmod is not None:
430	self.timemodif = tmod
431	self.dimensions = [dx, dy, dz, dt, ds]
432	self.vardimensions = [vdx, vdy, vdz, vdt, vds]
433
434	def variable_compute(idir,var,ftests,db):
435	""" Function to retrieve the computation way of a given variable using a series of test files
436	iwdir= directory with the test files
437	var= name of the variable to test
438	filtests= dictionary with the test files for each header
439	"""
440	fname='variable_compute'
441
442	cancompute = True
443	for headerf in ftests.keys():
444	filen=ftests[headerf]
445	try:
446	with gen.Capturing() as output:
447	vmod, vdiag = ncvar.computevar_model(var, idir + '/' + filen)
448	except:
449	print errmsg
450	print 'ncvar.computevar_model(' + var + ', ' + idir + '/' + filen + ')'
451	for s1out in output: print s1out
452	quit(-1)
453
454	if vmod is None and vdiag is None:
455	cancompute = False
456	else:
457	cancompute = True
458	# Should be considered that variable can also be computed by both ways?
459	break
460
461	if not cancompute:
462	print warnmsg
463	print ' ' + fname + ": there is no way to compute '" + var + \
464	"' for model '" + mod
465	# Too extrict!
466	# quit(-1)
467
468	# Getting only the first possibility of model and diagnostic
469	if vmod is not None:
470	modV = vmod[0]
471	else:
472	modV = None
473	if vdiag is not None:
474	diagV = vdiag[0]
475	else:
476	diagV = None
477
478	varcomp = VariableInf(var, headerf, modV, diagV)
479
480	# a ';' list 'varcompute' it is created for each variable giving:
481	# [var]\|[vark]\|[headerf][varmod]\|[vardiag]
482	# This list will be used to compute a new file for each variable
483
484	if db:
485	print 'Variable information _______'
486	gen.printing_class(varcomp)
487
488	return varcomp
489
490	def get_operations_var(vc,db):
491	""" Function to provide the operations to make for each variable from 'VAR_' dictionary
492	vc= 'VAR_' dictionary, dictionary with all the parameters started with 'VAR_'
493	"""
494	fname = 'get_operations_var'
495
496	LH = len('VAR_')
497	iop = 1
498	# list of operations
499	doopers = {}
500	# operations by variable
501	operationsvar = {}
502	# individual operations by variable (not repeated)
503	indivoperationsvar = {}
504
505	# Variables with a calculation whic requires vertical interpolation to all the
506	# file (operations stating by 'VAR_pinterp')
507	LVp = len('VAR_pinterp')
508	varglobalp = []
509
510	for oper in vc.keys():
511	Loper = len(oper)
512	opn = oper[LH:Loper+1]
513	vns = vc[oper].split(':')
514	ops1 = opn.split('+')
515	doopers[opn] = ops1
516	for vn in vns:
517	if not operationsvar.has_key(vn):
518	operationsvar[vn] = [opn]
519	indivoperationsvar[vn] = ops1
520	else:
521	opers = operationsvar[vn]
522	opers.append(opn)
523	operationsvar[vn] = opers
524	opers1 = indivoperationsvar[vn]
525	for op1 in ops1:
526	if not gen.searchInlist(opers1, op1):
527	opers1.append(op1)
528	indivoperationsvar[vn] = opers1
529
530	if oper[0:LVp] == 'VAR_pinterp':
531	for vn in vns:
532	if not gen.searchInlist(varglobalp,vn): varglobalp.append(vn)
533
534	iop = iop + 1
535
536	if db:
537	print ' operations to make _______'
538	gen.printing_dictionary(doopers)
539	print ' operations by variable _______'
540	gen.printing_dictionary(operationsvar)
541	print ' individual operations by variable _______'
542	gen.printing_dictionary(indivoperationsvar)
543	print ' variables with the vertical interpolation for all the file _______'
544	print ' #', varglobalp
545
546	return doopers, operationsvar, indivoperationsvar, varglobalp
547
548	def pinterp_var(oper):
549	""" Function to retrieve characteristics of the vertical interpolation for the operation
550	oper= operation
551	Wether vertical interpolation is:
552	# 'global': for all file
553	# 'local': at a given step of the process
554	# 'none': no vertical interpolation
555	>>> pinterp_var('last+pinterp+xmean')
556	local
557	>>> pinterp_var('pinterp+tmean+xmean')
558	global
559	>>> pinterp_var('xmean')
560	none
561	"""
562	fname = 'pinterp_var'
563
564	pinn = 'pinterp'
565
566	Lpinterp = len(pinn)
567	if oper[0:Lpinterp] == pinn:
568	pkind = 'global'
569	return pkind
570
571	if oper.find(pinn) != -1:
572	pkind = 'local'
573	else:
574	pkind = 'none'
575
576	return pkind
577
578	def compvars_listconstruct(config, minf, Files, TestFiles, idir, odir, debug):
579	""" Function to construct the list of variables to compute
580	config= dictionary with the configuration of the execution
581	variables to compute are get with all values started by `VAR_' in 'module_graphics.dat', and then
582	providing a consecutive number of calculations separated by '+'
583	VAR_[calc1]+[calc2] = tas:wss
584	will compute first [calc1] and then [calc2] for 'tas' and 'wss'
585	minf= class with information about the model
586	Files= dictionary of files for heach header
587	TestFiles= dictionary of files to test how to compute the variable
588	odir= output directory
589	"""
590	fname='compvars_listconstruct'
591
592	varcomp = gen.get_specdictionary_HMT(config,H='VAR_')
593
594	if debug:
595	print ' variables to compute ________'
596	gen.printing_dictionary(varcomp)
597
598	# Getting operations by variable
599	opers, varoper, indivaroper, varglobp = get_operations_var(varcomp,debug)
600
601	strfiles = gen.dictKeysVals_stringList(Files)
602	strtestfiles = gen.dictKeysVals_stringList(TestFiles)
603
604	Svarcompute = ''
605
606	# Main dictionary with all the variables and their calculation
607	allvarcomp = {}
608
609	ivop = 0
610	for vn in varoper:
611	vcomp = variable_compute(idir,vn,TestFiles,False)
612	for op in varoper[vn]:
613	# Creation of a String as: [CFname]\|[operation]\|[fheader]\|[vmodel]\|[vdiag]\|[globalP]
614	# [CFname]: CF-name of the variable (must appear in 'variables_values.dat')
615	# [operation]: operation to compute
616	# [fheader]: header of file with the required variables
617	# [vmodel]: Variable from model output
618	# [vdiag]: Varible as diagnostic from different variables of the model output
619	# [globalP]: Wether vertical interpolation is:
620	# 'global': for all file
621	# 'local': at a given step of the process
622	# 'none': no vertical interpolation
623	globalP = pinterp_var(op)
624	if vcomp.model is not None:
625	if type(vcomp.model) == type(list([1, 2])):
626	Smodel = ':'.join(vcomp.model)
627	elif type(vcomp.model) == type('A'):
628	Smodel = vcomp.model
629	else:
630	Smodel = 'None'
631	if vcomp.diag is not None:
632	if type(vcomp.diag) == type(list([1, 2])):
633	Sdiag = ':'.join(vcomp.diag)
634	elif type(vcomp.diag) == type('A'):
635	Sdiag = vcomp.diag
636	else:
637	Sdiag = 'None'
638
639	# To get a combination of operations, all the precedent steps are kept
640	# thus, they are also available !
641	seqops = op.split('+')
642	for seqop in seqops:
643	if seqop == seqops[0]: seqopS=seqop
644	else: seqopS = seqopS + '+' + seqop
645	allvarcomp[vn+'_'+seqopS] = [vcomp.fheader, vcomp.model, vcomp.diag, \
646	globalP]
647
648	Svc = vcomp.name + '\|' + seqopS + '\|' + vcomp.fheader + '\|' + \
649	Smodel + '\|' + Sdiag + '\|' + globalP
650	if ivop == 0:
651	Svarcompute = Svc
652	else:
653	Svarcompute = Svarcompute + ',' + Svc
654
655	ivop = ivop + 1
656
657	if debug:
658	print ' Variables to compute _______'
659	gen.printing_dictionary(allvarcomp)
660
661	# Outwritting the varcompute to avoid next time (if it is not filescratch!)
662	objf = open(odir + '/varcompute.inf', 'w')
663	objf.write('files: ' + strfiles + '\n')
664	objf.write('testfiles: ' + strtestfiles + '\n')
665	objf.write('varcompute: ' + Svarcompute + '\n')
666	objf.write('itotv: ' + str(ivop) + '\n')
667	objf.close()
668
669	return allvarcomp, ivop
670
671	def read_varcomp_file(filen):
672	""" Function to read 'varcompute.inf' and reconstruct the dictionaries
673	"""
674	fname = 'read_varcomp_file'
675
676	allvarcomp = {}
677
678	objf = open(filen, 'r')
679
680	for line in objf:
681	vals = line.replace('\n','').split(' ')
682	if vals[0] == 'files:':
683	Files = gen.stringList_dictKeysVals(vals[1])
684	elif vals[0] == 'testfiles:':
685	TestFiles = gen.stringList_dictKeysVals(vals[1])
686	elif vals[0] == 'varcompute:':
687	for VvalsS in vals[1].split(','):
688
689	Vvals = VvalsS.split('\|')
690	if Vvals[3] == 'None':
691	mod = None
692	else:
693	mod = Vvals[3].split(':')
694	if Vvals[4] == 'None':
695	diag = None
696	else:
697	diag = Vvals[4].split(':')
698
699	allvarcomp[Vvals[0]+'_'+Vvals[1]] = [Vvals[2], mod, diag, Vvals[5]]
700
701	# # To get a combination of operations, all the precedent steps are kept
702	# # thus, they are also available !
703	# seqops = Vvals[1].split('+')
704	# for seqop in seqops:
705	# if seqop == seqops[0]: seqopS=seqop
706	# else: seqopS = seqopS + '+' + seqop
707	# allvarcomp[Vvals[0]+'_'+seqopS] = [Vvals[2], mod, diag, Vvals[5]]
708	elif vals[0] == 'itotv:':
709	ivop = int(vals[1])
710
711	objf.close()
712
713	return Files, TestFiles, allvarcomp, ivop
714
715	def pinterpS(gP):
716	""" For that variables which require vertical interpolation 'p' suffix to the file header is added
717	gP= kind of pinterp: 'global', 'local', 'none'
718	"""
719	fname = 'pinterpS'
720
721	if gP != 'none':
722	gPS = 'p'
723	else:
724	gPS = ''
725
726	return gPS
727
728	def compute_variable(minf, idir, usefiles, odir, cvar, gP, scr, pyH, Tref, Tunits, db):
729	""" Function to compute a variable
730	minf= class with the information of the model
731	idir= directory with the input files
732	usefiles= dictionary of files as dict([headerf]) = [file1], ..., [fileN]
733	odir= directory to write the output files
734	cvar= class with the information of the variable: 'name', 'fheader', 'varmod', 'vardiag'
735	gP= kind of vertical interpolation ('global', 'local', 'none')
736	scr= should it be done from the scratch?
737	pyH= location of the python HOME
738	Tref= CF time reference
739	Tunits= CF time units
740	"""
741	fname='compute_variable'
742
743	CFvarn=cvar.name
744	headerf = cvar.fheader
745	modvar = cvar.model
746	diagvar = cvar.diag
747
748	cfiles = usefiles[headerf]
749
750	# dimensions
751	dnx = minf.dimxn
752	dny = minf.dimyn
753	# var-dimensions
754	vdnx = minf.vardxn
755	vdny = minf.vardyn
756
757	# Computing separately and then joinging for all files
758	Ntotfiles = len(cfiles)
759	Nzeros = len(str(Ntotfiles))
760	NStot = str(Ntotfiles).zfill(Nzeros)
761
762	# For that variables which require vertical interpolation 'p' suffix to the
763	# file header is added
764	SgP = pinterpS(gP)
765
766	# Getting in working dir
767	os.chdir(odir)
768
769	# File to keep track of all operations
770	otrackf = open( odir + '/all_computevars.inf', 'a')
771
772	# Computing variable
773	ifile=1
774	for cf in cfiles:
775	ifS = str(ifile).zfill(Nzeros)
776	ifilen = odir + '/' + CFvarn + '_' + headerf + SgP + '_' + ifS + '-' + \
777	NStot + '.nc'
778	fileon = odir + '/' + CFvarn + '_' + headerf + SgP + '.nc'
779
780	if scr:
781	sout = sub.call('rm ' + ifilen + ' >& /dev/null', shell=True)
782	sout = sub.call('rm ' + fileon + ' >& /dev/null', shell=True)
783
784	if not os.path.isfile(ifilen) and not os.path.isfile(fileon):
785	# Since model direct values are retrieved from `variables_valules.dat' which was
786	# initially coincived as a way to only give variable attributes, range and color
787	# bars, if a variable has a diagnostic way to be computed, the later one will be
788	# preferred
789	if db:
790	print ' ' + fname + ": creation of variable file '" + CFvarn + \
791	"' in file '" + ifilen + "' ..."
792	print ' model variable:', modvar
793	print ' diagostic variable:', diagvar
794
795	if diagvar is None:
796	# model variable
797	if db: print ' '+fname+ ": variable direct from model '" +modvar+ "'"
798	values = modvar + ',0,-1,-1'
799	vs = modvar + ',' + vdnx + ',' + vdny + ',' + vdnz + ',' + vdnt
800
801	try:
802	with gen.Capturing() as output:
803	ncvar.DataSetSection_multivars(values, idir + '/' + cf, vs)
804	except:
805	print errmsg
806	print 'ncvar.DataSetSection_multivars('+values+', '+idir+'/'+cf+ \
807	', '+vs+')'
808	for s1out in output: print s1out
809	quit(-1)
810
811	newfile, loc = gen.search_sec_list(output,'succesfull')
812	ofile = newfile[0].split(' ')[7]
813	sub.call('mv ' + ofile + ' ' + ifilen, shell=True)
814
815	# Keeping track of the operations
816	pyins = pyH + "/nc_var.py -f " + idir + '/' + cf + \
817	" -o DataSetSection_multivars -v " + vs + " -S '" + values + "'"
818	otrackf.write('\n')
819	otrackf.write('# ' + CFvarn + " " + modvar + '\n')
820	otrackf.write('python ' + pyins + '\n')
821
822	# CF renaming of variable
823	ncvar.chvarname(CFvarn,ifilen,modvar)
824	else:
825	# diagnostic variable
826	if db: print ' '+fname+ ": variable as diagnostic '" +diagvar[0]+ "'"
827	dims = dnt+'@'+vdnt+','+dnz+'@'+vdnz+','+dny+'@'+vdny+','+dnx+'@'+vdnx
828	diagn = diagvar[0]
829	Ndiagvars = len(diagvar)
830	diagc = '@'.join(diagvar[1:])
831
832	values = '-f ' + idir + '/' + cf + " -d '" + dims + "' -v '" + \
833	diagn + '\|' + diagc + "'"
834	try:
835	with gen.Capturing() as output:
836	sout = sub.call('python ' +pyH+ '/diagnostics.py ' + values, \
837	shell=True)
838	except:
839	print errmsg
840	print 'python ' + pyH + '/diagnostics.py ' + values
841	for s1out in output: print s1out
842	quit(-1)
843	if db:
844	for s1out in output: print s1out
845
846	sout = sub.call(' mv diagnostics.nc ' + ifilen, shell=True)
847
848	# Keeping track of the operations
849	pyins = 'python ' + pyH + '/diagnostics.py ' + values
850	otrackf.write('\n')
851	otrackf.write('# ' + CFvarn + " " + diagn + '\n')
852	otrackf.write(pyins + '\n')
853
854	# Attaching necessary variables for the pressure interpolation
855	if gP != 'none':
856	if minf.name == 'WRF' or minf.name == 'WRF_LMDZ':
857	requiredinterpvars = ['P', 'PB', 'PSFC', 'PH', 'PHB', 'HGT', 'T',\
858	'QVAPOR', 'XLONG', 'XLAT', 'Times']
859	elif minf.name == 'LMDZ':
860	requiredinterpvars = ['pres', 'psol', 'geop', 'phis', 'temp', \
861	'ovap', 'lon', 'lat', 'time_counter']
862	else:
863	print errmsg
864	print fname + ": for model '" + minf.name + "' required " + \
865	" variables for vertical interpolation are not known !!"
866	quit(-1)
867
868	print " " + fname + ": adding variables:", requiredinterpvars, \
869	' to allow pressure interpolation'
870	for rqv in requiredinterpvars:
871	try:
872	with gen.Capturing() as output:
873	ncvar.fvaradd(idir+'/'+cf+','+rqv,ifilen)
874	except:
875	print errmsg
876	print 'fvaradd('+idir+'/'+cf+','+rqv+', '+ifilen+')'
877	for s1out in output: print s1out
878	quit(-1)
879	if db:
880	for s1out in output: print s1out
881
882	# CFification of files
883	if minf.name == 'WRF' or minf.name == 'WRF_LMDZ':
884	values = vdnx + ':' + vdny + ':'+ Tref + ':' + Tunits
885	try:
886	with gen.Capturing() as output:
887	ncvar.WRF_toCF(values, ifilen)
888	except:
889	print errmsg
890	print 'WRF_toCF('+ values +', ' + ifilen + ')'
891	for s1out in output: print s1out
892	quit(-1)
893	if db:
894	for s1out in output: print s1out
895	elif minf.name == 'LMDZ':
896	try:
897	with gen.Capturing() as output:
898	ncvar.LMDZ_toCF(Tunits + '!since!' + Tref, ifilen)
899	except:
900	print errmsg
901	print 'LMDZ_toCF('+ Tunits +'!since!'+ Tref + ', ' + ifilen + ')'
902	for s1out in output: print s1out
903	quit(-1)
904	if db:
905	for s1out in output: print s1out
906	else:
907	print errmsg
908	print ' ' + fname + ": no CFification for model '" +minf.name+ "' !!"
909	print " available ones: 'WRF', 'WRF_LMDZ', 'LMDZ'"
910	quit(-1)
911
912	ifile = ifile + 1
913
914	otrackf.close()
915
916	# Joining variable files
917	if not os.path.isfile(fileon):
918	if db:
919	print ' ' + fname + ": concatenating all variable files to create '" + \
920	fileon + "..."
921	try:
922	with gen.Capturing() as output:
923	ncvar.netcdf_fold_concatenation_HMT('./,time', CFvarn+'_'+headerf+ \
924	SgP+'_,-,.nc', 'all')
925	except:
926	print errmsg
927	print 'netcdf_fold_concatenation_HMT(./time, '+CFvarn+'_'+headerf+SgP+ \
928	'_,-,.nc, all)'
929	for s1out in output: print s1out
930	quit(-1)
931	if db:
932	for s1out in output: print s1out
933
934	sout = sub.call('mv netcdf_fold_concatenated_HMT.nc ' + fileon, shell=True)
935	if os.path.isfile(fileon):
936	sout = sub.call('rm '+CFvarn +'_'+ headerf+SgP + '_-.nc >& /dev/null',\
937	shell=True)
938
939	# Re-setting time of final concatenated file
940	if minf.timemodif is not None:
941	if db:
942	print ' ' + fname + ": Modifying times of the file by '" + \
943	minf.timemodif + "' !!"
944	try:
945	with gen.Capturing() as output:
946	ncvar.time_reset(minf.timemodif, fileon,'time')
947	except:
948	print errmsg
949	print 'time_reset(' + minf.timemodif + ', ' + fileon + ', time)'
950	for s1out in output: print s1out
951	quit(-1)
952	if db:
953	for s1out in output: print s1out
954
955	return
956
957	def compute_statistics(minf, config, idir, usefiles, odir, cvar, gP, Opers, scr, db):
958	""" Function to compute different statistics it will take previous steps if they
959	are availale
960	minf= class with the information of the model
961	config= dictionary with the configuration of the experiment
962	idir= directory with the input files
963	usefiles= ',' list of files to use [file1],...,[fileN]
964	odir= directory to write the output files
965	cvar= class with the information of the variable: 'name', 'fheader', 'varmod', 'vardiag'
966	gP= kind of vertical interpolation ('global', 'local', 'none')
967	Opers= kind of operation: (as possible multiple consecutive combination of operations separated by '+'
968	[calc1]+[calc2] will compute first [calc1] and then [calc2]
969	acc: temporal accumulated values
970	diff: differences between models
971	direct: no statistics
972	last: last temporal value
973	Lmean: latitudinal mean values
974	Lsec: latitudinal section (latitudinal value must be given, [var]@[lat])
975	lmean: longitudinal mean values
976	lsec: longitudinal section (longitudinal value must be given, [var]@[lat])
977	pinterp: pressure interpolation (to the given $plevels)
978	tmean: temporal mean values
979	tstd: temporal standard deviation values
980	tturb: Taylor's turbulence decomposition value (x - <x>) for time
981	tvar: temporal variance values
982	xmean: x-axis mean values
983	xvar: x-axis variance values
984	ymean: y-axis mean values
985	zsum: vertical aggregated values
986	scr= should it be done from the scratch?
987	"""
988	fname='compute_statistics'
989
990	# Getting a previous file name to continue with(for combinations)
991	CFvarn=cvar.name
992	headerf = cvar.fheader
993	modvar = cvar.model
994	diagvar = cvar.diag
995
996	cfiles = usefiles[headerf]
997
998	# dimensions
999	dnx = minf.dimxn
1000	dny = minf.dimyn
1001	# var-dimensions
1002	vdnx = minf.vardxn
1003	vdny = minf.vardyn
1004	# Model dimension-variable dictionary
1005	Modeldimvardict = {dnx: vdnx, dny: vdny, dnz: vdnz, dnt: vdnt}
1006
1007	# Some experiment configuration values
1008	# plevels= ':' separated list of pressures (in Pa) to use for the vertical interpolation
1009	plevels = config['plevels']
1010	# pyH= location of the python HOME
1011	pyH = config['pyHOME']
1012	# Tref= CF time reference
1013	Tref = config['CFreftime']
1014	# Tunits= CF time units
1015	Tunits = config['CFunitstime']
1016	# opsur = operation surnames
1017	opsur = gen.get_specdictionary_HMT(config, H='opsur_',M='',T='')
1018	Opsurs = {}
1019	for opk in opsur:
1020	opn = opk.split('_')[1]
1021	vals = opsur[opk]
1022	Opsurs[opn] = vals.split(':')
1023	# if db:
1024	# print ' ' + fname + ' operation surnames _______'
1025	# gen.printing_dictionary(Opsurs)
1026
1027	# Getting in working dir
1028	os.chdir(odir)
1029
1030	# File to keep track of all operations
1031	otrackf = open( odir + '/all_statsvars.inf', 'a')
1032
1033	# For that variables which require vertical interpolation 'p' suffix to the file
1034	# header is added.
1035	# After some operations some CF dimension-variables might not be kept
1036
1037	if gP != 'none':
1038	SgP = 'p'
1039	varnCFs = ['lon', 'lat', 'pres', 'time']
1040	else:
1041	SgP = ''
1042	varnCFs = ['lon', 'lat', 'time']
1043
1044	# CF dimension-variable dictionary
1045	CFdimvardict = {'lon': 'lon', 'lat': 'lat', 'pres': 'pres', 'time': 'time'}
1046
1047	# Input file
1048	ifilen = odir + '/' + CFvarn + '_' + headerf + SgP + '.nc'
1049
1050	# Computing variable statisitcs
1051	istats=0
1052
1053	# List of not computed statistics
1054	wrongstats = []
1055
1056	opers = Opers.split('+')
1057	Fopers = ''
1058	if db: print " computing statistics of variable:'", CFvarn, "' operation '" + \
1059	Opers + "'...'"
1060	for op in opers:
1061
1062	# File name from previous operations, and name of the variable within the
1063	# file (some operations change the name of the variable)
1064	if op == opers[0]:
1065	Fopers = op
1066	prevfile = ifilen
1067	vninF = CFvarn
1068	else:
1069	Fopers = Fopers + '_' + op
1070	prevfile = fileon
1071	fileon = odir + '/' + CFvarn + '_' + headerf + SgP + '_' + Fopers + '.nc'
1072
1073	# Adding required variables for the vertical interpolation in all operations
1074	SvarnCFs = ',' + ','.join(varnCFs)
1075
1076	if gP != 'none':
1077	if gP == 'local':
1078	CFvarnp = vninF + ',P,PB,PSFC,PH,PHB,HGT,T,QVAPOR,XLONG,XLAT,Times'+ \
1079	SvarnCFs
1080	else:
1081	CFvarnp = vninF + SvarnCFs
1082	else:
1083	CFvarnp = vninF + SvarnCFs
1084
1085	if scr:
1086	sout = sub.call('rm ' + fileon + ' >& /dev/null', shell=True)
1087
1088	hprevfile = prevfile[0:len(prevfile)-3]
1089
1090	if db: print ' ', op, "using '" + prevfile + "' "
1091	# Just produce file in case output file does not exist
1092	if not os.path.isfile(fileon):
1093	if db:
1094	print ' ' + fname + ": creation of file '" + fileon + "' ..."
1095	dofile = True
1096	else:
1097	dofile = False
1098
1099	# Variables to be kept in the final file
1100	varkeep = []
1101
1102	if op == 'acc':
1103	# temporal accumulated values
1104	print " " + fname + ": kind '" + op + "' not ready !!"
1105	wrongstats.append(CFvarn + '_' + opers)
1106	break
1107	elif op == 'direct':
1108	# no statistics
1109	if dofile:
1110	sout = sub.call('mv ' + prevfile + ' ' + fileon, shell=True)
1111	elif op == 'last':
1112	# last temporal value
1113	vals='time,-9,0,0'
1114	if dofile:
1115	try:
1116	with gen.Capturing() as output:
1117	ncvar.DataSetSection(vals,prevfile)
1118	except:
1119	print errmsg
1120	print 'DataSetSection('+vals+',', prevfile, ')'
1121	for s1out in output: print s1out
1122	quit(-1)
1123
1124	if db:
1125	for s1out in output: print s1out
1126
1127	sout = sub.call('mv ' + hprevfile + '_time_B-9-E0-I0.nc ' + fileon, \
1128	shell=True)
1129
1130	# Keeping the operations
1131	pyins=pyH + "/nc_var.py -o DataSetSection -S '" + vals + "' -f " + \
1132	prevfile
1133	otrackf.write("\n")
1134	otrackf.write("# " + CFvarn + " " + Fopers + "\n")
1135	otrackf.write(pyins + "\n")
1136
1137	elif op =='Lmean':
1138	# latitudinal mean values
1139	print " " + fname + ": kind '" + op + "' not ready !!"
1140	wrongstats.append(CFvarn + '_' + opers)
1141	break
1142	elif op =='Lsec':
1143	# latitudinal section (latitudinal value must be given, [var]@[lat])
1144	print " " + fname + ": kind '" + op + "' not ready !!"
1145	wrongstats.append(CFvarn + '_' + opers)
1146	break
1147	elif op =='lmean':
1148	# longitudinal mean values
1149	print " " + fname + ": kind '" + op + "' not ready !!"
1150	wrongstats.append(CFvarn + '_' + opers)
1151	break
1152	elif op =='lsec':
1153	# longitudinal section (longitudinal value must be given, [var]@[lon])
1154	print " " + fname + ": kind '" + op + "' not ready !!"
1155	wrongstats.append(CFvarn + '_' + opers)
1156	break
1157	elif op == 'pinterp':
1158	# pinterp: pressure interpolation (to the given $plevels)
1159	vals=plevels + ',1,1'
1160	if dofile:
1161	try:
1162	with gen.Capturing() as output:
1163	ncvar.pinterp(vals,prevfile,vninF)
1164	except:
1165	print errmsg
1166	print 'pinterp('+vals+',', prevfile, ','+vninF+')'
1167	for s1out in output: print s1out
1168	quit(-1)
1169
1170	if db:
1171	for s1out in output: print s1out
1172
1173	sout = sub.call('mv pinterp.nc ' + fileon, shell=True)
1174
1175	# Keeping the operations
1176	pyins=pyH + "/nc_var.py -o pinterp -S '" + vals + "' -f " + \
1177	prevfile + "-v " + vninF
1178	otrackf.write("\n")
1179	otrackf.write("# " + CFvarn + " " + Fopers + "\n")
1180	otrackf.write(pyins + "\n")
1181
1182	# adding CF lon,lat,time in WRF files
1183	if minf.name == 'WRF':
1184	values = vdnx + ':' + vdny + ':'+ Tref + ':' + Tunits
1185	ncvar.WRF_toCF(values, fileon)
1186
1187	# vertical interpolation variables are no more needed
1188	CFvarnp = vninF
1189	gP = 'none'
1190
1191	elif op == 'tmean':
1192	# temporal mean values
1193	vals='time\|-1,time,mean,' + ':'.join(varnCFs) + ':' + vdnz
1194	dims = gen.dictvar_listS(varnCFs,CFdimvardict,',','@') + dnz+'@'+vdnz
1195	if dofile:
1196	try:
1197	with gen.Capturing() as output:
1198	ncvar.file_oper_alongdims(vals,prevfile,CFvarnp)
1199	except:
1200	print errmsg
1201	print 'file_oper_alongdims('+vals+',', prevfile, ','+CFvarnp+')'
1202	for s1out in output: print s1out
1203	quit(-1)
1204
1205	if db:
1206	for s1out in output: print s1out
1207
1208	sout = sub.call('mv file_oper_alongdims_mean.nc '+fileon, shell=True)
1209
1210	# Keeping the operations
1211	pyins=pyH + "/nc_var.py -o file_oper_alongdims -S '" + vals + \
1212	"' -f " + prevfile + " -v " + CFvarnp
1213	otrackf.write("\n")
1214	otrackf.write("# " + CFvarn + " " + Fopers + "\n")
1215	otrackf.write(pyins + "\n")
1216
1217	# removing dimension variable-dimension 'time'
1218	varnCFs.remove('time')
1219
1220	varkeep.append('timestats')
1221
1222	elif op == 'tstd':
1223	# temporal standard deviation values
1224	vals='time\|-1,time,std,' + ':'.join(varnCFs) + ':' + vdnz
1225	dims = gen.dictvar_listS(varnCFs,CFdimvardict,',','@') + dnz+'@'+vdnz
1226	if dofile:
1227	try:
1228	with gen.Capturing() as output:
1229	ncvar.file_oper_alongdims(vals,prevfile,CFvarnp)
1230	except:
1231	print errmsg
1232	print 'file_oper_alongdims('+vals+',', prevfile, ','+CFvarnp+')'
1233	for s1out in output: print s1out
1234	quit(-1)
1235
1236	if db:
1237	for s1out in output: print s1out
1238
1239	sout = sub.call('mv file_oper_alongdims_std.nc '+fileon, shell=True)
1240
1241	# Keeping the operations
1242	pyins=pyH + "/nc_var.py -o file_oper_alongdims -S '" + vals + \
1243	"' -f " + prevfile + " -v " + CFvarnp
1244	otrackf.write("\n")
1245	otrackf.write("# " + CFvarn + " " + Fopers + "\n")
1246	otrackf.write(pyins + "\n")
1247
1248	# removing dimension variable-dimension 'time'
1249	varnCFs.remove('time')
1250
1251	varkeep.append('timestats')
1252
1253	elif op == 'tturb':
1254	# temporal turbulent values
1255	vals='time\|-1,time,turb,' + ':'.join(varnCFs) + ':' + vdnz
1256	dims = gen.dictvar_listS(varnCFs,CFdimvardict,',','@') + dnz+'@'+vdnz
1257	if dofile:
1258	try:
1259	with gen.Capturing() as output:
1260	ncvar.file_oper_alongdims(vals,prevfile,CFvarnp)
1261	except:
1262	print errmsg
1263	print 'file_oper_alongdims('+vals+',', prevfile, ','+CFvarnp+')'
1264	for s1out in output: print s1out
1265	quit(-1)
1266
1267	if db:
1268	for s1out in output: print s1out
1269
1270	sout = sub.call('mv file_oper_alongdims_turb.nc '+fileon, shell=True)
1271
1272	# Keeping the operations
1273	pyins=pyH + "/nc_var.py -o file_oper_alongdims -S '" + vals + \
1274	"' -f " + prevfile + " -v " + CFvarnp
1275	otrackf.write("\n")
1276	otrackf.write("# " + CFvarn + " " + Fopers + "\n")
1277	otrackf.write(pyins + "\n")
1278
1279	# removing dimension variable-dimension 'time'
1280	varkeep.append('timestats')
1281
1282	elif op == 'tvar':
1283	# temporal variance values (Taylor's turbulence)
1284	vals='time\|-1,time,var,' + ':'.join(varnCFs) + ':' + vdnz
1285	dims = gen.dictvar_listS(varnCFs,CFdimvardict,',','@') + dnz+'@'+vdnz
1286	if dofile:
1287	try:
1288	with gen.Capturing() as output:
1289	ncvar.file_oper_alongdims(vals,prevfile,CFvarnp)
1290	except:
1291	print errmsg
1292	print 'file_oper_alongdims('+vals+',', prevfile, ','+CFvarnp+')'
1293	for s1out in output: print s1out
1294	quit(-1)
1295
1296	if db:
1297	for s1out in output: print s1out
1298
1299	sout = sub.call('mv file_oper_alongdims_var.nc '+fileon, shell=True)
1300
1301	# Keeping the operations
1302	pyins=pyH + "/nc_var.py -o file_oper_alongdims -S '" + vals + \
1303	"' -f " + prevfile + " -v " + CFvarnp
1304	otrackf.write("\n")
1305	otrackf.write("# " + CFvarn + " " + Fopers + "\n")
1306	otrackf.write(pyins + "\n")
1307
1308	# removing dimension variable-dimension 'time'
1309	varnCFs.remove('time')
1310
1311	varkeep.append('timestats')
1312
1313	elif op == 'xmean':
1314	# x-axis mean values
1315	vals='lon\|-1,lon,mean,' + ':'.join(varnCFs) + ':' + vdnz
1316	dims = gen.dictvar_listS(varnCFs,CFdimvardict,',','@') + dnz+'@'+vdnz
1317	if dofile:
1318	try:
1319	with gen.Capturing() as output:
1320	ncvar.file_oper_alongdims(vals,prevfile,CFvarnp)
1321	except:
1322	print errmsg
1323	print 'file_oper_alongdims('+vals+',', prevfile, ','+CFvarnp+')'
1324	for s1out in output: print s1out
1325	quit(-1)
1326
1327	if db:
1328	for s1out in output: print s1out
1329
1330	sout = sub.call('mv file_oper_alongdims_mean.nc '+fileon, shell=True)
1331
1332	# Keeping the operations
1333	pyins=pyH + "/nc_var.py -o file_oper_alongdims -S '" + vals + \
1334	"' -f " + prevfile + " -v " + CFvarnp
1335	otrackf.write("\n")
1336	otrackf.write("# " + CFvarn + " " + Fopers + "\n")
1337	otrackf.write(pyins + "\n")
1338
1339	# removing dimension variable-dimension 'lon'
1340	varnCFs.remove('lon')
1341
1342	varkeep.append('lonstats')
1343
1344	elif op == 'xvar':
1345	# x-axis variance values
1346	vals='lon\|-1,lon,var,' + ':'.join(varnCFs) + ':' + vdnz
1347	dims = gen.dictvar_listS(varnCFs,CFdimvardict,',','@') + dnz+'@'+vdnz
1348	if dofile:
1349	try:
1350	with gen.Capturing() as output:
1351	ncvar.file_oper_alongdims(vals,prevfile,CFvarnp)
1352	except:
1353	print errmsg
1354	print 'file_oper_alongdims('+vals+',', prevfile, ','+CFvarnp+')'
1355	for s1out in output: print s1out
1356	quit(-1)
1357
1358	if db:
1359	for s1out in output: print s1out
1360
1361	sout = sub.call('mv file_oper_alongdims_var.nc '+fileon, shell=True)
1362
1363	# Keeping the operations
1364	pyins=pyH + "/nc_var.py -o file_oper_alongdims -S '" + vals + \
1365	"' -f " + prevfile + " -v " + CFvarnp
1366	otrackf.write("\n")
1367	otrackf.write("# " + CFvarn + " " + Fopers + "\n")
1368	otrackf.write(pyins + "\n")
1369
1370	# removing dimension variable-dimension 'lon'
1371	varnCFs.remove('lon')
1372
1373	varkeep.append('lonstats')
1374
1375	elif op == 'ymean':
1376	# y-axis mean values
1377	vals='lat\|-1,lat,mean,' + ':'.join(varnCFs) + ':' + vdnz
1378	dims = gen.dictvar_listS(varnCFs,CFdimvardict,',','@') + dnz+'@'+vdnz
1379	if dofile:
1380	try:
1381	with gen.Capturing() as output:
1382	ncvar.file_oper_alongdims(vals,prevfile,CFvarnp)
1383	except:
1384	print errmsg
1385	print 'file_oper_alongdims('+vals+',', prevfile, ','+CFvarnp+')'
1386	for s1out in output: print s1out
1387	quit(-1)
1388
1389	if db:
1390	for s1out in output: print s1out
1391
1392	sout = sub.call('mv file_oper_alongdims_mean.nc '+fileon, shell=True)
1393
1394	# Keeping the operations
1395	pyins=pyH + "/nc_var.py -o file_oper_alongdims -S '" + vals + \
1396	"' -f " + prevfile + " -v " + CFvarnp
1397	otrackf.write("\n")
1398	otrackf.write("# " + CFvarn + " " + Fopers + "\n")
1399	otrackf.write(pyins + "\n")
1400
1401	# removing dimension variable-dimension 'lat'
1402	varnCFs.remove('lat')
1403
1404	varkeep.append('latstats')
1405
1406	elif op == 'zsum':
1407	# vertical aggregated values
1408	print " " + fname + ": kind '" + op + "' not ready !!"
1409	wrongstats.append(CFvarn + '_' + opers)
1410	break
1411	else:
1412	print errmsg
1413	print ' ' + fname + ": operation '" + op + "' not ready !!"
1414	quit(-1)
1415
1416	# End of kind of operation
1417
1418	# Variable name in file (vninF) changed due to operation
1419	# but only if previous operation does not the same 'statistic'
1420	chvn = gen.dictionary_key_list(Opsurs, op)
1421	if chvn is not None:
1422	oldvninF = vninF
1423	vninF = vninF + chvn
1424	CFvarnp = CFvarnp.replace(oldvninF,vninF)
1425
1426	if dofile:
1427	if len(varkeep) > 0:
1428	varkeepS = ',' + ','.join(varkeep)
1429	else:
1430	varkeepS = ''
1431
1432	# Adding such CF dimension variables which might not be in the
1433	# post-operation file
1434	try:
1435	with gen.Capturing() as output:
1436	varinfile = ncvar.ivars(fileon)
1437	except:
1438	print errmsg
1439	print 'ivars('+fileon+')'
1440	for s1out in output: print s1out
1441	quit(-1)
1442
1443	for CFvn in varnCFs:
1444	if not gen.searchInlist(varinfile, CFvn):
1445	if db:
1446	print ' ' + fname + ": recupering CF variable '" + CFvn + "'"
1447	try:
1448	with gen.Capturing() as output:
1449	ncvar.fvaradd(prevfile+','+CFvn, fileon)
1450	except:
1451	print errmsg
1452	print 'fvaradd(' + prevfile + ',' + CFvn +', '+fileon+')'
1453	for s1out in output: print s1out
1454	quit(-1)
1455
1456	if db:
1457	for s1out in output: print s1out
1458
1459	totalvarkeeps = CFvarnp + ',' + ','.join(varnCFs) + varkeepS
1460	try:
1461	with gen.Capturing() as output:
1462	oclean = ncvar.cleaning_varsfile(totalvarkeeps,fileon)
1463	except:
1464	print errmsg
1465	print 'cleaning_varsfile('+totalvarkeeps+','+fileon+')'
1466	for s1out in output: print s1out
1467	quit(-1)
1468
1469	if db:
1470	for s1out in output: print s1out
1471
1472	# End of operations
1473	if db:
1474	print fname + ": successful calculation of '" + vninF + "' from '" + \
1475	CFvarn + "' '" + Opers + "' !!"
1476
1477	if len(wrongstats) > 1:
1478	print warnmsg
1479	print ' ' + fname + ": statisitcs not possible to compute:", wrongstats
1480
1481	return
1482
1483	def compute_vars(config, modinf, idir, odir, Files, allvarcomp, fscr, debug):
1484	""" Function to compute the variables
1485	config= Configuration of the experiment
1486	modinf= class with information about the model
1487	idir= input experiment folder
1488	odir= output experiment folder
1489	Files= dictionary with the header and the files' correspondence
1490	allvarcomp= dictionary with all the variables to compute and their information
1491	fscr= whether files should be done from the scratch or not
1492	"""
1493	fname = 'compute_vars'
1494
1495	for vopn in allvarcomp.keys():
1496	# variable & operation
1497	vn = vopn.split('_')[0]
1498	op = vopn.split('_')[1]
1499	vopnV = allvarcomp[vopn]
1500	fheader = vopnV[0]
1501	if type(vopnV[1]) == type([1, 2]) and len(vopnV[1]) == 1:
1502	vals = vopnV[1]
1503	model = vals[0]
1504	else:
1505	model = vopnV[1]
1506	diag = vopnV[2]
1507	globalP = vopnV[3]
1508
1509	if debug:
1510	print ' ' + vn + ': ' + op
1511
1512	# Computing CF variable
1513	if model is not None or diag is not None:
1514	vinf = VariableInf(vn,fheader,model,diag)
1515	if op.find('pinterp') != -1:
1516	SP = 'p'
1517	else:
1518	SP = ''
1519
1520	# Comppute variable
1521	finalfilen = odir + '/' + vn + '_' + fheader + SP + '.nc'
1522	if fscr:
1523	sout = sub.call('rm ' + finalfilen + ' >& /dev/null', shell=True)
1524
1525	if not os.path.isfile(finalfilen):
1526	compute_variable(modinf, idir, Files, odir, vinf, globalP, fscr, \
1527	config['pyHOME'], config['CFreftime'], config['CFunitstime'], debug)
1528
1529	# Compute variable statistics
1530	finalfilen = odir + '/' + vn + '_' + fheader + SP + '_' + \
1531	op.replace('+','_') + '.nc'
1532	if fscr:
1533	sout = sub.call('rm ' + finalfilen + ' >& /dev/null', shell=True)
1534
1535	if not os.path.isfile(finalfilen):
1536	compute_statistics(modinf, cnf, iwdir, Files, owdir, vinf, globalP, \
1537	op, fscr, debug)
1538
1539	else:
1540	print errmsg
1541	print ' ' + fname + ": neither 'model' or 'diag' variables for '" + vn \
1542	+ "' !!"
1543	quit(-1)
1544
1545	return
1546
1547	def get_plots_var(vc,pltk,db):
1548	""" Function to provide the plots to make for each variable from pltk ('DIRPLT_', or 'DIFFPLT_') dictionary
1549	vc= pltk dictionary, dictionary with all the parameters started with pltk
1550	pltk= kind of plots. Values started by 'DIRPLT_', or 'DIFFPLT_' in model configuration
1551	"""
1552	fname = 'get_plots_var'
1553
1554	LH = len(pltk)
1555	# list of plots
1556	doplots = {}
1557	# plots by variable
1558	plotsvar = {}
1559	# individual plots by variable (not repeated)
1560	indivplotsvar = {}
1561
1562	# Variables with a plot which requires vertical interpolation
1563	varplotp = []
1564
1565	ipl = 0
1566	for plot in vc.keys():
1567	Lplot = len(plot)
1568	pln = plot[LH:Lplot+1]
1569	varops = vc[plot].split(':')
1570	for vn in varops:
1571	VOn = vn.split('#')
1572	# Variables and operations in plot
1573	vplt = []
1574	oplt = []
1575	voplt = []
1576	for vop in VOn:
1577	vv = vop.split('\|')[0]
1578	oo = vop.split('\|')[1]
1579	vplt.append(vv)
1580	oplt.append(oo)
1581	if not gen.searchInlist(voplt,vop): voplt.append(vop)
1582	vpltS = '#'.join(vplt)
1583	opltS = '#'.join(oplt)
1584
1585	if not doplots.has_key(pln):
1586	doplots[pln] = [vn]
1587	else:
1588	plots = doplots[pln]
1589	plots.append(vn)
1590	doplots[pln] = plots
1591
1592	if not plotsvar.has_key(vn):
1593	plotsvar[vn] = [pln]
1594	else:
1595	plots = plotsvar[vn]
1596	plots.append(pln)
1597	plotsvar[vn] = plots
1598
1599	for VOv in voplt:
1600	if not indivplotsvar.has_key(VOv):
1601	indivplotsvar[VOv] = [pln]
1602	else:
1603	plots = indivplotsvar[VOv]
1604	plots.append(pln)
1605	indivplotsvar[VOv] = plots
1606
1607	if vc[plot].find('pinterp') != -1:
1608	if not gen.searchInlist(varplotp,vn): varplotp.append(pln)
1609
1610	ipl = ipl + 1
1611
1612	if db:
1613	print ' plots to make _______'
1614	gen.printing_dictionary(doplots)
1615	print ' plots by variable\|operation _______'
1616	gen.printing_dictionary(plotsvar)
1617	print ' individual plots by variable\|operation _______'
1618	gen.printing_dictionary(indivplotsvar)
1619	print ' plots with the vertical interpolation _______'
1620	print ' #', varplotp
1621
1622	return doplots, plotsvar, indivplotsvar, varplotp
1623
1624	def plots_listconstruct(config, pkind, finf, odir, debug):
1625	""" Function to create the list of plots to draw
1626	config= Configuration of the experiment
1627	pkind= kind of plots ('DIRPLT', 'DIFFPLT')
1628	finf= file with the instructions of the plots
1629	odir= output experiment folder
1630	"""
1631	fname='plots_listconstruct'
1632
1633	dirplot = gen.get_specdictionary_HMT(config, H=pkind+'_')
1634
1635	# No plots because there are not plots of this kind
1636	if dirplot is None:
1637	print warnmsg
1638	print ' ' + fname + ": no plots for type '" + pkind + "' !!"
1639	return {}, 0
1640
1641	if debug:
1642	print " '" + pkind + "' plots to draw ________"
1643	gen.printing_dictionary(dirplot)
1644
1645	# Getting plots by variable
1646	plots, varplot, indivarplot, plotp = get_plots_var(dirplot, pkind+'_', debug)
1647
1648	Svarplot = gen.dictKeysVals_stringList(plots,cV=':')
1649
1650	Nplots = 0
1651	for pl in Svarplot.split(','):
1652	Nplots = Nplots + len(pl.split(':'))
1653
1654	# Outwritting the varcompute to avoid next time (if it is not filescratch!)
1655	objf = open(finf, 'w')
1656	objf.write('plots: ' + Svarplot + '\n')
1657	objf.write('itotp: ' + str(Nplots) + '\n')
1658	objf.close()
1659
1660	return plots, Nplots
1661
1662	def read_plot_file(figfile):
1663	""" Function to read the file with the information about the plots to draw
1664	figfile= file with the information
1665	"""
1666	fname = 'read_plot_file'
1667
1668	if not os.path.isfile(figfile):
1669	print errormsg
1670	print ' ' + fname + ": figures file '" + figfile + "' does not exist !!"
1671	quit(-1)
1672
1673	objf = open(figfile, 'r')
1674	for line in objf:
1675	if line[0:1] != '#' and len(line) > 1:
1676	values = line.replace('\\n','').split(' ')
1677	if values[0] == 'plots:':
1678	if values[1] != 'none':
1679	plots = gen.stringList_dictKeysVals(values[1],cV=':')
1680	else:
1681	pots = {}
1682	elif values[0] == 'itotp:':
1683	Nplots = int(values[1])
1684
1685	objf.close()
1686
1687	return plots, Nplots
1688
1689	def varnoper(vn, oper, opsurdict):
1690	""" Function to provide name of the variable after operation as result of addition of the 'surnames'
1691	vn= variable name
1692	oper= '+' separated list of operations
1693	opsurdict= dictionary with the surname of operations
1694	>>> OpSurDict = {'mean': ['tmean', 'xmean', 'ymean']}
1695	>>> varnoper('tas', 'pinterp+tmean+xmean', OpSurDict)
1696	tasmeanmean
1697	"""
1698	fname = 'varnoper'
1699
1700	newvn = vn
1701
1702	for op in oper.split('+'):
1703	surname = gen.dictionary_key_list(opsurdict,op)
1704	if surname is not None:
1705	newvn = newvn + surname
1706
1707	return newvn
1708
1709	def draw_plot(kplot, CFvplot, fplot, vplot, dplot, pplot, finame, tfig, kfig, mapval,\
1710	tvals, expgraphc, od, pyH, fscr, db):
1711	""" Function to draw a plot
1712	kplot= kind of plot
1713	CFvplot= list of the CFnames of the variables
1714	fplot= list with the files in the plot
1715	vplot= list with the name of the variables in each file
1716	dplot= list of the dimensions in each file
1717	pplot= list of pictoric characteristics for each variable
1718	finame= name of the figure
1719	tfig= title of the figure
1720	kfig= kind of the figure
1721	mapval= value of the map
1722	tvals= list with time values
1723	tunits: units of time ('!' for spaces)
1724	timekind: kind of time ticks in the plot
1725	timefmt: format of time ticks in the plot
1726	timelabel: label of time-axis in the plot
1727	expgraphc= dictionary with expriment graphical characteristics
1728	od= output directory
1729	pyH= python HOME
1730	fscr= whether figure should be done from scratch
1731	"""
1732	fname = 'draw_plot'
1733
1734	if fscr:
1735	sout = sub.call('rm ' + finame + ' >& /dev/null', shell=True)
1736
1737	os.chdir(od)
1738
1739	# CF variable-dimensions
1740	CFvardims = {'lon':'lon', 'lat': 'lat', 'pres': 'pres', 'time': 'time'}
1741
1742	# Tracking file with all the figures
1743	trkobjf = open('all_figures.inf', 'a')
1744
1745	# Time-characteristics
1746	tunits = tvals[0]
1747	timekind = tvals[1]
1748	timefmt = tvals[2]
1749	timelabel = tvals[3]
1750
1751	# Graphical labels and configuration
1752
1753
1754	if not os.path.isfile(finame):
1755	if db:
1756	print " drawing '" + finame + "' ..."
1757	print ' plot kind:', kplot
1758	print ' CFvars:', CFvplot
1759	print ' files:', fplot
1760	print ' in file vars:', vplot
1761	print ' dims:', dplot
1762	print ' pictoric values:', pplot
1763	print ' figure name:', finame
1764	print ' title:', tfig.replace('!',' ')
1765	print ' kind:', kfig
1766	print ' map:', mapval
1767	print ' time units:', tunits
1768
1769	if kplot == 'diffmap2Dsfc':
1770	print " " + fname + ": kind of plot '" + kplot + "' not ready !!"
1771	quit(-1)
1772	elif kplot == 'diffmap2Dz':
1773	print " " + fname + ": kind of plot '" + kplot + "' not ready !!"
1774	quit(-1)
1775	elif kplot == 'map2Dsfc':
1776	print " " + fname + ": kind of plot '" + kplot + "' not ready !!"
1777	quit(-1)
1778	elif kplot == 'map3D':
1779	print " " + fname + ": kind of plot '" + kplot + "' not ready !!"
1780	quit(-1)
1781	elif kplot == '2lines':
1782	figtit = tfig.replace('!','\|')
1783
1784	lineAstdn = CFvplot[0]
1785	lineBstdn = CFvplot[1]
1786	figfs = ','.join(fplot)
1787
1788	lineA = pplot[0]
1789	Arange = str(lineA[0]) + ',' + str(lineA[1])
1790	# Values are changed from specific values in `model_graphics.dat'
1791	if lineA[3].find('%') == -1:
1792	Aline = lineA[2]
1793	Akind = lineA[3]
1794	Amark = lineA[4]
1795	Asize = lineA[5]
1796	else:
1797	Aline = 'red'
1798	Akind = '-'
1799	Amark = ','
1800	Asize = 2.
1801
1802	lineB = pplot[0]
1803	Brange = str(lineB[0]) + ',' + str(lineB[1])
1804	# Values are changed from specific values in `model_graphics.dat'
1805	if lineB[3].find('%') == -1:
1806	Bline = lineB[2]
1807	Bkind = lineB[3]
1808	Bmark = lineB[4]
1809	Bsize = lineB[5]
1810	else:
1811	Bline = 'blue'
1812	Bkind = '-'
1813	Bmark = ','
1814	Bsize = '2.'
1815
1816	# It is assumed that if the space variable is 'lon': is desired a
1817	# (lon, vals) plot if it is 'lat': then (vals, lat) plot
1818	if gen.searchInlist(dplot[0],'lon'):
1819	spacedim ='lon'
1820	axisvals = 'y'
1821	figmid = 'longitudinal'
1822	elif gen.searchInlist(dplot[0],'lat'):
1823	spacedim = 'lat'
1824	axisvals = 'x'
1825	figmid = 'meridional'
1826	elif gen.searchInlist(dplot[0],'pres'):
1827	spacedim = 'pres'
1828	axisvals = 'x'
1829	figmid = 'vertical'
1830	else:
1831	print errmsg
1832	print ' ' + fname + ": in '2lines' only ready for: 'lon', 'lat'," + \
1833	" 'pres' as common dimension !!"
1834	print ' dimensions in plot:', dplot[0]
1835	quit(-1)
1836
1837	figtit = '!'.join(tfig.split('!')[0:2]) + '!' + figmid + '!' + \
1838	'!'.join(tfig.split('!')[2:])
1839
1840	graphvals = spacedim + ':' + Arange + ':' + Brange + ':Extrs:' + \
1841	axisvals + ':' + ','.join(CFvplot) + ':' + Aline + ',' + Bline + \
1842	':2.,2.:' + Akind + ',' + Bkind + ':2.,2.:,:' + figtit + ':' + \
1843	spacedim + ':0\|auto:' + finame.replace('.'+kfig, '') + ':' + kfig + \
1844	':True'
1845
1846	fvarS = ','.join(vplot)
1847
1848	drwins = 'draw_2lines'
1849	plotins = 'python ' + pyH + '/drawing.py -f ' + figfs +' -o ' + drwins + \
1850	" -S '" + graphvals + "' -v " + fvarS
1851
1852	try:
1853	with gen.Capturing() as output:
1854	sout = sub.call(plotins, shell=True)
1855	except:
1856	print errmsg
1857	print drwins + '(' + graphvals + ',' + figfs + ',' + fvarS + ')'
1858	print sout
1859	for s1out in output: print s1out
1860	quit(-1)
1861
1862	# keeping all figures
1863	trkobjf.write('\n')
1864	trkobjf.write("# '" + kplot + "'; " + tfig.replace('!',' ') + '\n')
1865	trkobjf.write(plotins + '\n')
1866
1867	elif kplot == '2linesTime':
1868	figtit = tfig.replace('!','\|')
1869
1870	lineAstdn = CFvplot[0]
1871	lineBstdn = CFvplot[1]
1872	figfs = ','.join(fplot)
1873
1874	lineA = pplot[0]
1875	Arange = str(lineA[0]) + ',' + str(lineA[1])
1876	# Values are changed from specific values in `model_graphics.dat'
1877	if lineA[3].find('%') == -1:
1878	Aline = lineA[2]
1879	Akind = lineA[3]
1880	Amark = lineA[4]
1881	Asize = lineA[5]
1882	else:
1883	Aline = 'red'
1884	Akind = '-'
1885	Amark = ','
1886	Asize = '2.'
1887
1888	lineB = pplot[1]
1889	Brange = str(lineB[0]) + ',' + str(lineB[1])
1890	# Values are changed from specific values in `model_graphics.dat'
1891	if lineB[3].find('%') == -1:
1892	Bline = lineB[2]
1893	Bkind = lineB[3]
1894	Bmark = lineB[4]
1895	Bsize = lineB[5]
1896	else:
1897	Bline = 'blue'
1898	Bkind = '-'
1899	Bmark = ','
1900	Bsize = '2.'
1901
1902	timeaxis = 'x'
1903	figmid = 'evolution'
1904
1905	figtit = '!'.join(tfig.split('!')[0:2]) + '!' + figmid + '!' + \
1906	'!'.join(tfig.split('!')[2:])
1907
1908	graphvals = 'time:' + Arange + ':' + Brange + ':' + timekind + ';' + \
1909	timefmt + ':' + timeaxis + ':' + ','.join(CFvplot) + ':' + Aline + ',' \
1910	+ Bline + ':2.,2.:' + Akind + ',' + Bkind + ':2.,2.:,:' + figtit + ':' \
1911	+ timelabel + ':0\|auto:' + finame.replace('.'+kfig, '') + ':' + \
1912	kfig + ':True'
1913
1914	fvarS = ','.join(vplot)
1915
1916	drwins = 'draw_2lines_time'
1917	plotins = 'python ' + pyH + '/drawing.py -f ' + figfs +' -o ' + drwins + \
1918	" -S '" + graphvals + "' -v " + fvarS
1919
1920	try:
1921	with gen.Capturing() as output:
1922	sout = sub.call(plotins, shell=True)
1923	except:
1924	print errmsg
1925	print drwins + '(' + graphvals + ',' + figfs + ',' + fvarS + ')'
1926	print sout
1927	for s1out in output: print s1out
1928	quit(-1)
1929
1930	# keeping all figures
1931	trkobjf.write('\n')
1932	trkobjf.write("# '" + kplot + "'; " + tfig.replace('!',' ') + '\n')
1933	trkobjf.write(plotins + '\n')
1934
1935	elif kplot == 'Nlines':
1936	figtit = tfig.replace('!','\|')
1937
1938	linestdn = CFvplot[0]
1939	figfs = ','.join(fplot)
1940
1941	Nlines = len(fplot)
1942
1943	lablines = []
1944	for il in range(Nlines):
1945	linev = pplot[il]
1946	# Values are changed from specific values in `model_graphics.dat'
1947	if linev[3].find('%') == -1:
1948	line = linev[2]
1949	kind = linev[3]
1950	mark = linev[4]
1951	size = linev[5]
1952	if il == 0:
1953	vrange = str(linev[0]) + ',' + str(linev[1])
1954	kinds = kind
1955	lines = line
1956	marks = mark
1957	sizes = str(size)
1958	else:
1959	kinds = kinds + ',' + kind
1960	lines = lines + ',' + line
1961	marks = marks + '@' + mark
1962	sizes = sizes + ',' + str(size)
1963	else:
1964	lines = 'None'
1965	kinds = '-'
1966	marks = ','
1967	sizes = '2.'
1968
1969	secsf = fplot[il].split('/')
1970	Nsecsf = len(secsf)
1971	expvs = expgraphc[secsf[Nsecsf-3] + '/' + secsf[Nsecsf-2]]
1972	lablines.append(expvs.label)
1973
1974	# It is assumed that if the space variable is 'lon': is desired a
1975	# (lon, vals) plot if it is 'lat': then (vals, lat) plot
1976	if gen.searchInlist(dplot[0],'lon'):
1977	spacedim ='lon'
1978	sdlab ='lon ($degrees\ East$)'
1979	axisvals = 'y'
1980	figmid = 'longitudinal'
1981	elif gen.searchInlist(dplot[0],'lat'):
1982	spacedim = 'lat'
1983	sdlab = 'lat ($degrees\ North$)'
1984	axisvals = 'x'
1985	figmid = 'meridional'
1986	elif gen.searchInlist(dplot[0],'pres'):
1987	spacedim = 'pres'
1988	sdlab = 'pres ($Pa$)'
1989	axisvals = 'x'
1990	figmid = 'vertical'
1991	else:
1992	print errmsg
1993	print ' ' + fname + ": in 'Nlines' only ready for: 'lon', 'lat'," + \
1994	" 'pres' as common dimension !!"
1995	print ' dimensions in plot:', dplot[0]
1996	quit(-1)
1997
1998	figtit = '!'.join(tfig.split('!')[0:2]) + '!' + figmid + '!' + \
1999	'!'.join(tfig.split('!')[2:])
2000	figtit = figtit.replace('!',' ')
2001	leglabels = ','.join(lablines)
2002
2003	graphvals = spacedim + ':' + axisvals + ':' + sdlab + ':' + leglabels + \
2004	':' + CFvplot[0] +':'+ figtit + ':0\|auto:' + lines +':'+ kinds + ':' + \
2005	marks +':'+ sizes +':'+ sizes +':all:'+ finame.replace('.'+kfig, '') + \
2006	':' + kfig + ':True'
2007
2008	fvarS = vplot[0]
2009
2010	drwins = 'draw_lines'
2011	plotins = 'python ' + pyH + '/drawing.py -f ' + figfs +' -o ' + drwins + \
2012	" -S '" + graphvals + "' -v " + fvarS
2013
2014	try:
2015	with gen.Capturing() as output:
2016	sout = sub.call(plotins, shell=True)
2017	except:
2018	print errmsg
2019	print drwins + '(' + graphvals + ',' + figfs + ',' + fvarS + ')'
2020	print sout
2021	for s1out in output: print s1out
2022	quit(-1)
2023
2024	# keeping all figures
2025	trkobjf.write('\n')
2026	trkobjf.write("# '" + kplot + "'; " + tfig.replace('!',' ') + '\n')
2027	trkobjf.write(plotins + '\n')
2028
2029	elif kplot == 'Nlines_time':
2030	figtit = tfig.replace('!','\|')
2031
2032	linestdn = CFvplot[0]
2033	figfs = ','.join(fplot)
2034
2035	Nlines = len(fplot)
2036
2037	lablines = []
2038	for il in range(Nlines):
2039	linev = pplot[il]
2040	# Values are changed from specific values in `model_graphics.dat'
2041	if linev[3].find('%') == -1:
2042	line = linev[2]
2043	kind = linev[3]
2044	mark = linev[4]
2045	size = linev[5]
2046	if il == 0:
2047	vrange = str(linev[0]) + ',' + str(linev[1])
2048	kinds = kind
2049	lines = line
2050	marks = mark
2051	sizes = str(size)
2052	else:
2053	kinds = kinds + ',' + kind
2054	lines = lines + ',' + line
2055	marks = marks + '@' + mark
2056	sizes = sizes + ',' + str(size)
2057	else:
2058	lines = 'None'
2059	kinds = '-'
2060	marks = ','
2061	sizes = '2.'
2062
2063	secsf = fplot[il].split('/')
2064	Nsecsf = len(secsf)
2065	expvs = expgraphc[secsf[Nsecsf-3] + '/' + secsf[Nsecsf-2]]
2066	lablines.append(expvs.label)
2067
2068	valsaxis = 'y'
2069	figmid = 'evolution'
2070
2071	figtit = '!'.join(tfig.split('!')[0:2]) + '!' + figmid + '!' + \
2072	'!'.join(tfig.split('!')[2:])
2073	figtit = figtit.replace('!',' ')
2074	leglabels = ','.join(lablines)
2075
2076	graphvals = 'time;' + valsaxis + ';time;' + leglabels + ';' + \
2077	CFvplot[0] + ';' + figtit + ';None;time\|' + '\|'.join(tvals) + \
2078	';0\|auto;' + kfig + ';' + kinds + ';' + lines + ';' + marks + ';' + \
2079	sizes + ';' + sizes + ';all;-1;True'
2080
2081	fvarS = vplot[0]
2082
2083	drwins = 'draw_lines_time'
2084	plotins = 'python ' + pyH + '/drawing.py -f ' + figfs +' -o ' + drwins + \
2085	" -S '" + graphvals + "' -v " + fvarS
2086
2087	try:
2088	with gen.Capturing() as output:
2089	sout = sub.call(plotins, shell=True)
2090	except:
2091	print errmsg
2092	print drwins + '(' + graphvals + ',' + figfs + ',' + fvarS + ')'
2093	print sout
2094	for s1out in output: print s1out
2095	quit(-1)
2096
2097	sout = sub.call('mv lines_time.' + kfig + ' ' + finame, shell=True)
2098
2099	# keeping all figures
2100	trkobjf.write('\n')
2101	trkobjf.write("# '" + kplot + "'; " + tfig.replace('!',' ') + '\n')
2102	trkobjf.write(plotins + '\n')
2103
2104	elif kplot == 'shadcont2Dsfc':
2105	figtit = tfig.replace('!','\|')
2106	shdstdn = CFvplot[0]
2107	cntstdn = CFvplot[1]
2108	figfs = ','.join(fplot)
2109
2110	shad = pplot[0]
2111	srange = str(shad[0]) + ',' + str(shad[1])
2112	cbar = shad[2]
2113	cnt = pplot[1]
2114	crange = str(cnt[0]) + ',' + str(cnt[1])
2115	cfmt = cnt[3]
2116	ckind = cnt[4]
2117	cline = cnt[5]
2118
2119	graphvals = ','.join(CFvplot)
2120	graphvals = graphvals + ':lon\|-1,lat\|-1:lon\|-1,lat\|-1:lon:lat:auto:' + \
2121	cbar + ',auto,auto:' + ckind + ',' + cline + ':' + cfmt + ':' + \
2122	srange + ':' + crange + ',9:' + figtit + ':' + kfig + ':None:' + \
2123	mapval + ':True'
2124
2125	fvarS = ','.join(vplot)
2126
2127	drwins = 'draw_2D_shad_cont'
2128	plotins = "python " + pyH + "/drawing.py -f " + figfs + " -o " + drwins +\
2129	" -S '" + graphvals + "' -v " + fvarS
2130	try:
2131	with gen.Capturing() as output:
2132	sout = sub.call(plotins, shell=True)
2133	except:
2134	print errmsg
2135	print drwins + '(' + graphvals + ',' + figfs + ',' + fvarS + ')'
2136	print sout
2137	for s1out in output: print s1out
2138	quit(-1)
2139
2140	sout = sub.call('mv 2Dfields_shadow-contour.'+kfig+' '+finame, shell=True)
2141
2142	# keeping all figures
2143	trkobjf.write('\n')
2144	trkobjf.write("# '" + kplot + "'; " + tfig.replace('!',' ') + '\n')
2145	trkobjf.write(plotins + '\n')
2146
2147	elif kplot == 'shadconthovmsfc':
2148	figtit = tfig.replace('!','\|')
2149	shdstdn = CFvplot[0]
2150	cntstdn = CFvplot[1]
2151	figfs = ','.join(fplot)
2152
2153	shad = pplot[0]
2154	srange = str(shad[0]) + ',' + str(shad[1])
2155	cbar = shad[2]
2156	cnt = pplot[1]
2157	crange = str(cnt[0]) + ',' + str(cnt[1])
2158	cfmt = cnt[3]
2159	ckind = cnt[4]
2160	cline = cnt[5]
2161	# It is assumed that if the space variable is 'lon': is desired a
2162	# (lon, time) plot if it is 'lat': then (time, lat) plot
2163	if gen.searchInlist(dplot,'lon'):
2164	spacedim ='lon'
2165	figmid = 'longitudinal\|evolution\|of'
2166	figtit = '\|'.join(tfig.split('!')[0:2]) + '\|' + figmid + '\|' + \
2167	'\|'.join(tfig.split('!')[2:])
2168	reverse= 'None'
2169	dims= spacedim+'\|-1,time\|-1;'+spacedim+'\|-1,time\|-1;'+spacedim+';time'
2170	else:
2171	spacedim='lat'
2172	figmid = 'meridional\|evolution\|of'
2173	figtit = '\|'.join(tfig.split('!')[0:2]) + '\|' + figmid + '\|' + \
2174	'\|'.join(tfig.split('!')[2:])
2175	reverse='None'
2176	dims= spacedim+'\|-1,time\|-1;'+spacedim+'\|-1,time\|-1;time;'+spacedim
2177
2178	graphvals = ','.join(CFvplot)
2179	graphvals = graphvals + ';' + dims +';auto;' + cbar + ',auto,auto' + ';' \
2180	+ ckind + ',' + cline + ';' + cfmt + ';' + srange + ';' + crange + \
2181	',9;' + figtit + ';' + kfig + ';' + reverse + ';' + 'time\|' + tunits + \
2182	'\|'+ timekind + '\|' + timefmt + '\|' + timelabel + ';True'
2183
2184	fvarS = ','.join(vplot)
2185
2186	drwins = 'draw_2D_shad_cont_time'
2187	plotins = 'python ' + pyH + '/drawing.py -f ' + figfs +' -o ' + drwins + \
2188	" -S '" + graphvals + "' -v " + fvarS
2189
2190	try:
2191	with gen.Capturing() as output:
2192	sout = sub.call(plotins, shell=True)
2193	except:
2194	print errmsg
2195	print drwins + '(' + graphvals + ',' + figfs + ',' + fvarS + ')'
2196	print sout
2197	for s1out in output: print s1out
2198	quit(-1)
2199
2200	sout = sub.call('mv 2Dfields_shadow-contour.'+kfig+' '+finame, shell=True)
2201
2202	# keeping all figures
2203	trkobjf.write('\n')
2204	trkobjf.write("# '" + kplot + "'; " + tfig.replace('!',' ') + '\n')
2205	trkobjf.write(plotins + '\n')
2206
2207	elif kplot == 'shadcont2Dzsec':
2208	figtit = tfig.replace('!','\|')
2209	shdstdn = CFvplot[0]
2210	cntstdn = CFvplot[1]
2211	figfs = ','.join(fplot)
2212
2213	shad = pplot[0]
2214	srange = str(shad[0]) + ',' + str(shad[1])
2215	cbar = shad[2]
2216	cnt = pplot[1]
2217	crange = str(cnt[0]) + ',' + str(cnt[1])
2218	cfmt = cnt[3]
2219	ckind = cnt[4]
2220	cline = cnt[5]
2221
2222	# It is assumed that if the space variable is 'lon': is desired a
2223	# (lon, pres) plot it it is 'lat': then (pres, lat) plot
2224	if gen.searchInlist(dplot, 'lon'):
2225	spacedim='lon'
2226	figmid = 'long.\|vert.\|cross\|sec.\|of'
2227	figtit = '\|'.join(tfig.split('!')[0:2]) + '\|' + figmid + '\|' + \
2228	'\|'.join(tfig.split('!')[2:])
2229	dims= spacedim+'\|-1,pres\|-1:'+spacedim+'\|-1,pres\|-1:'+spacedim+':pres'
2230	else:
2231	spacedim='lat'
2232	figmid = 'mer.\|vert.\|cross\|sec.\|of'
2233	figtit = '\|'.join(tfig.split('!')[0:2]) + '\|' + figmid + '\|' + \
2234	'\|'.join(tfig.split('!')[2:])
2235	# NOT WORKING? dims= spacedim+'\|-1,pres\|-1:'+spacedim+'\|-1,pres\|-1:'+'pres:'+spacedim
2236	dims= spacedim+'\|-1,pres\|-1:'+spacedim+'\|-1,pres\|-1:'+spacedim+':pres'
2237
2238	reverse = 'flip@y'
2239
2240	graphvals = ','.join(CFvplot)
2241	graphvals = graphvals + ':' + dims + ':auto:' + cbar + ',auto,auto:' + \
2242	ckind + ',' + cline + ':' + cfmt + ':' + srange + ':' + crange + ',9:' \
2243	+ figtit + ':' + kfig + ':' + reverse + ':None:True'
2244
2245	fvarS = ','.join(vplot)
2246
2247	drwins = 'draw_2D_shad_cont'
2248	plotins = 'python ' + pyH + '/drawing.py -f ' + figfs + ' -o ' + drwins + \
2249	" -S '" + graphvals + "' -v " + fvarS
2250	try:
2251	with gen.Capturing() as output:
2252	sout0 = sub.call(plotins, shell=True)
2253	except:
2254	print errmsg
2255	print drwins + '(' + graphvals + ',' + figfs + ',' + fvarS + ')'
2256	for s1out in output: print s1out
2257	quit(-1)
2258
2259	sout = sub.call('mv 2Dfields_shadow-contour.'+kfig+' '+finame, shell=True)
2260
2261	# keeping all figures
2262	trkobjf.write('\n')
2263	trkobjf.write("# '" + kplot + "'; " + tfig.replace('!',' ') + '\n')
2264	trkobjf.write(plotins + '\n')
2265
2266	else:
2267	print errmsg
2268	print " " + fname + ": plot kind '" + kplot + "' not ready !!"
2269	quit(-1)
2270
2271	trkobjf.close()
2272
2273	return
2274
2275	def optitle(op,opexp):
2276	""" Function to creat the operation section in the title of a figure ('!' for spaces)
2277	op= '+' separated list of operations
2278	opexp= dictionary with the `explanation'(text in title) to appear for each operation
2279	>>> optitle('pinterp+xmean+tmean',{'pinterp':'pinterp', 'xmean':'xmean', 'tmean':'tmean'})
2280	$_{[pinterp\!xmean\!tmean]}$
2281	"""
2282	fname = 'optitle'
2283
2284	opvals = op.split('+')
2285	for op1 in opvals:
2286	if not opexp.has_key(op1):
2287	print errmsg
2288	print ' '+fname+": no explanation for operation '"+op1+"' !!"
2289	print ' provided:', opexp.keys()
2290	quit(-1)
2291	if op1 == opvals[0]:
2292	titop = '$_{[' + opexp[op1].replace(' ','!')
2293	if len(opvals) == 1: titop = titop + ']}$'
2294	elif op1 == opvals[len(opvals)-1]:
2295	titop = titop + '\!' + opexp[op1].replace(' ','!') + ']}$'
2296	else:
2297	titop = titop + '\!' + opexp[op1].replace(' ','!')
2298
2299	return titop
2300
2301	def create_figure_title(mod,exp,varops,explop):
2302	""" Function to create the title of a figure ('!' for spaces)
2303	mod: model name
2304	exp: experiment name
2305	varops: list with [var]\|[op] values from wich create the title
2306	explop: dictionary with the `explanation'(text in title) to appear for each operation
2307	>>> expops = {'pinterp':'pinterp', 'xmean':'xmean', 'tmean':'tmean', 'last':'last'}
2308	>>> create_figure_title('WRF','current',['ua\|pinterp+xmean+tmean', 'va\|pinterp+xmean+tmean'], expops)
2309	WRF!current!ua!&!va$_{[pinterp\!xmean\!tmean]}$
2310	>>> create_figure_title('WRF','current',['ua\|pinterp+xmean+tmean', 'ua\|pinterp+xmean+last'], expops)
2311	WRF!current!ua$_{[pinterp\!xmean\!tmean]}$!&!!$_{[pinterp\!xmean\!last]}$
2312	>>> create_figure_title('WRF','current',['ua\|pinterp+xmean+tmean', 'va\|pinterp+xmean+last'], expops)
2313	WRF!current!ua$_{[pinterp\!xmean\!tmean]}$!&!va$_{[pinterp\!xmean\!last]}$
2314	"""
2315	fname = 'create_figure_title'
2316
2317	titfig = mod + '!' + exp + '!'
2318
2319	varns = []
2320	opns = []
2321	for varop in varops:
2322	vn = varop.split('\|')[0]
2323	op = varop.split('\|')[1]
2324	varns.append(vn)
2325	opns.append(op)
2326	if varop == varops[0]:
2327	allvarequal = True
2328	allopequal = True
2329	varinfig = vn
2330	opinfig = op
2331	else:
2332	if vn != varinfig: allvarequal = False
2333	if op != opinfig: allopequal = False
2334
2335	Nvars = len(varns)
2336	Nops = len(opns)
2337
2338	if allvarequal:
2339	titfig = titfig + '!' + varns[0]
2340	if allopequal:
2341	opS = optitle(op,explop)
2342	titfig = titfig + opS
2343	else:
2344	for opn in opns:
2345	opS = optitle(opn,explop)
2346	if opn == opns[0]:
2347	titfig = titfig + opS
2348	elif opn == opns[Nops-1]:
2349	titfig = titfig + '!&!' + opS
2350	else:
2351	titfig = titfig + ',!' + opS
2352	else:
2353	if allopequal:
2354	opS = optitle(op,explop)
2355	for vn in varns:
2356	if vn == varns[0]:
2357	titfig = titfig + '!' + vn
2358	elif vn == varns[Nvars-1]:
2359	titfig = titfig + '!&!' + vn
2360	else:
2361	titfig = titfig + ',!' + vn
2362	titfig = titfig + opS
2363	else:
2364	for ivop in range(Nvars):
2365	vn = varns[ivop]
2366	op = opns[ivop]
2367	vnop = vn + optitle(op,explop)
2368
2369	if vn == varns[0]:
2370	titfig = titfig + '!' + vnop
2371	elif vn == varns[Nvars-1]:
2372	titfig = titfig + '!&!' + vnop
2373	else:
2374	titfig = titfig + ',!' + vnop
2375
2376	return titfig
2377
2378	def draw_plots(config, plots, mod, exp, odir, allvarcomp, figscr, debug):
2379	""" Function to draw all plots
2380	config= Configuration of the experiment
2381	plots= dictionary with the plots
2382	mod= model of the plots
2383	exp= experiment of the plots
2384	odir= output experiment folder
2385	allvarcomp= dictionary with all the variables to compute and their information
2386	figscr= whether figures should be done from the scratch or not
2387
2388	* Plot as
2389	{[kplot]} = [varn1]\|[op1]#[varn2]\|[op2]#[...[varnN]\|[opN]], ...
2390	[kplot] ___
2391	diffmap2Dsfc: 2D map of surface differences values of 1 variable
2392	diffmap2Dz: 2D map of 3D differences values of 1 variable
2393	map2Dsfc: 2D map of surface values of 1 variable
2394	map3D: 2D map of 3D values of 1 variable
2395	shadconthovmsfc: Hovmoeller diagrams of 2 variable at the surface in shadow and the other in contourn
2396	shadcont2Dsfc: 2D map of shadow (1st variable) and countour (2nd variable) [stvar1]#[stvar2]
2397	shadcont2Dzsec: 2D map of vertical section of 2 variables one in shadow and the other in contourn
2398	[varn]
2399	variable
2400	[op]
2401	'+' separated list of operations
2402	in figures with more than 1 variable, use '#' to separate the [varn]\|[op]
2403	"""
2404	fname = 'draw_plots'
2405
2406	os.chdir(odir)
2407
2408	# Dictionary with the operations with surnames for the operated variable
2409	opersurnames = {}
2410	opsur = gen.get_specdictionary_HMT(config, H='opsur_',M='',T='')
2411	for opsr in opsur.keys():
2412	opn = opsr.split('_')[1]
2413	vls = opsur[opsr].split(':')
2414	opersurnames[opn] = vls
2415
2416	# time values
2417	# Units time for the plots
2418	rd = config['CFreftime']
2419	tunits = config['CFunitstime'] + '!since!' + rd[0:4] + '-' + rd[4:6] + '-' + \
2420	rd[6:8] + '!' + rd[8:10] + ':' + rd[10:12] + ':' + rd[12:14]
2421	# time ticks kind
2422	tkind = config['timekind']
2423	# time ticks format
2424	tfmt = config['timefmt']
2425	# time axis label
2426	tlab = config['timelabel']
2427	timevals = [tunits, tkind, tfmt, tlab]
2428
2429	# Dictionary of plot specificities
2430	specplotkeyn = 'specificvarplot'
2431	plotspecifics = {}
2432	if config.has_key(specplotkeyn):
2433	# [minval]: minimum value
2434	# [maxval]: minimum value
2435	# [colorbar]: name of the colorbar (from matplotlib) to use
2436	# [cntformat]: format of the contour labels
2437	# [colorcnt]: color for the countor lines
2438	plotspecs = config[specplotkeyn].split(':')
2439	for pltspc in plotspecs:
2440	pltvls = pltspc.split('\|')
2441	vn = pltvls[0]
2442	op = pltvls[1]
2443	fn = pltvls[2]
2444	plotspecifics[fn + '_' + vn + '_' + op] = pltvls[3:]
2445	if debug:
2446	print 'Specific values for plots _______'
2447	gen.printing_dictionary(plotspecifics)
2448
2449	# Kind of figures
2450	kindfigure = config['kindfig']
2451
2452	# Graphical labels and configuration
2453	modgc, expgc = graphical_characteristics(config, debug)
2454	modv = modgc[mod]
2455	modlab = modv.label
2456	expv = expgc[mod+'/'+exp]
2457	explab = expv.label
2458
2459	# Map value
2460	mapvalue = config['mapval']
2461
2462	# pythone scripts HOME
2463	pyHOME = config['pyHOME']
2464
2465	# Title-text of operations
2466	opexplained = {}
2467	optits = config['titleoperations'].split(':')
2468	for optit in optits:
2469	opn = optit.split('\|')[0]
2470	opt = optit.split('\|')[1]
2471	opexplained[opn] = opt
2472	if debug:
2473	print 'Titles for operations _______'
2474	gen.printing_dictionary(opexplained)
2475
2476	for kplot in plots.keys():
2477	varsplt = plots[kplot]
2478	for varplt in varsplt:
2479	if debug:
2480	print " printing '" + kplot + "' var ':" + varplt + "'..."
2481	varops = varplt.split('#')
2482
2483	# CF variables in plot
2484	CFvarsplot = []
2485	# Files in plot
2486	filesplot = []
2487	# Variables in plot within the files
2488	varsplot = []
2489	# Dims in figure
2490	dimsplot = []
2491	# pictoric values in figure
2492	pictplot = []
2493	# Name of the figure
2494	figname = ''
2495	# Title of the figure
2496	titfigure = ''
2497
2498	ivp = 0
2499	for varop in varops:
2500	vn = varop.split('\|')[0]
2501	op = varop.split('\|')[1]
2502
2503	# CF variables in plot
2504	CFvarsplot.append(vn)
2505
2506	vnopS = vn + '_' + op
2507	if not allvarcomp.has_key(vnopS):
2508	print errmsg
2509	print ' ' + fname + ": no file for variable-operation '" + \
2510	vnopS + "' !!"
2511	vopvals = allvarcomp[vnopS]
2512	headf = vopvals[0]
2513	globalP = vopvals[3]
2514	gP = pinterpS(globalP)
2515
2516	filen = vn + '_' + headf + gP + '_' + op.replace('+','_') + '.nc'
2517	filesplot.append(filen)
2518	# Do we have processed the given variable?
2519	if not os.path.isfile(filen):
2520	print warnmsg
2521	print " " + fname + ": there is no file for variable '" + varop \
2522	+ "' skiping it !!"
2523	break
2524
2525	# Name of the variable inside the file
2526	vnsur = varnoper(vn, op, opersurnames)
2527	varsplot.append(vnsur)
2528
2529	# Dimensions in file
2530	try:
2531	with gen.Capturing() as output:
2532	dims = ncvar.idims(filen)
2533	except:
2534	print errmsg
2535	print 'ncvar.idims('+filen+')'
2536	for s1out in output: print s1out
2537	quit(-1)
2538
2539	dimsplot.append(dims)
2540
2541	# pictoric values for the figure
2542	Sfivaop = kplot + '_' + vn + '_' + op
2543	if plotspecifics.has_key(Sfivaop):
2544	pictvals = plotspecifics[Sfivaop]
2545	else:
2546	Vvals = gen.variables_values(vn)
2547	pictvals = [Vvals[2], Vvals[3], Vvals[6], '%g', 'fixc', 'black']
2548
2549	pictplot.append(pictvals)
2550
2551	# Header of the name of the figure
2552	if ivp == 0:
2553	figname = kplot +'_'+ vn + '_' + headf + '_' + op.replace('+','-')
2554	else:
2555	figname = figname + '-'+vn+'_'+headf+'_'+op.replace('+','-')
2556
2557	ivp = ivp + 1
2558	# End of variable-operation
2559	figname = figname + '.' + kindfigure
2560
2561	# Title of figure
2562	titfigure = create_figure_title(modlab, explab, varops, opexplained)
2563
2564	draw_plot(kplot, CFvarsplot, filesplot, varsplot, dimsplot, pictplot, \
2565	figname, titfigure, kindfigure, mapvalue, timevals, expgc, odir, \
2566	pyHOME, figscr, debug)
2567
2568	# End of variables-operations
2569
2570	# End of kind of plots
2571
2572	return
2573
2574	def get_differences_var(vc,kdiff,db):
2575	""" Function to provide the operations to make for each variable from kdiff ('DIFFOP_' or 'DIFFVAR_') dictionary
2576	vc= kdiff dictionary, dictionary with all the parameters started with kdiff
2577	kdiff= kind of differences 'op', 'var'
2578	"""
2579	fname = 'get_differences_var'
2580
2581	if kdiff == 'op':
2582	diffn = 'DIFFOP'
2583	elif kdiff == 'var':
2584	diffn = 'DIFFVAR'
2585	else:
2586	print errmsg
2587	print ' ' + fname + ": differences '" + kdiff + "' not ready !!"
2588	quit(-1)
2589	LD = len(diffn + '_')
2590
2591	iop = 1
2592	# list of operations by difference
2593	doopers = {}
2594	# operations by variable
2595	operationsvar = {}
2596	# individual operations by variable (not repeated)
2597	indivoperationsvar = {}
2598
2599	# Variables with a calculation whic requires vertical interpolation to all the
2600	# file (operations stating by [kdiff]'_pinterp')
2601	LVp = len(diffn + '_pinterp')
2602	varglobalp = []
2603
2604	for oper in vc.keys():
2605	Loper = len(oper)
2606	opn = oper[LD:Loper+1]
2607	vns = vc[oper].split(':')
2608	ops1 = opn.split('+')
2609	doopers[opn] = ops1
2610	for vn in vns:
2611	if not operationsvar.has_key(vn):
2612	operationsvar[vn] = [opn]
2613	indivoperationsvar[vn] = ops1
2614	else:
2615	opers = operationsvar[vn]
2616	opers.append(opn)
2617	operationsvar[vn] = opers
2618	opers1 = indivoperationsvar[vn]
2619	for op1 in ops1:
2620	if not gen.searchInlist(opers1, op1):
2621	opers1.append(op1)
2622	indivoperationsvar[vn] = opers1
2623
2624	if oper[0:LVp] == diffn + '_pinterp':
2625	for vn in vns:
2626	if not gen.searchInlist(varglobalp,vn): varglobalp.append(vn)
2627
2628	iop = iop + 1
2629
2630	if db:
2631	print ' operations to make _______'
2632	gen.printing_dictionary(doopers)
2633	print ' operations by variable _______'
2634	gen.printing_dictionary(operationsvar)
2635	print ' individual operations by variable _______'
2636	gen.printing_dictionary(indivoperationsvar)
2637	print ' variables with the vertical interpolation for all the file _______'
2638	print ' #', varglobalp
2639
2640	return doopers, operationsvar, indivoperationsvar, varglobalp
2641
2642	def diffvarop_listconstruct(config, mods, exps, reprj, odir, debug):
2643	""" Function to create the list of differences to compute and draw
2644	config= Configuration of the experiment
2645	mods= list with the models to use
2646	exps= list with the experiments to use
2647	reprj= list with whether re-rpojection needs to be done
2648	odir= output differences folder
2649	"""
2650	fname='diffvarop_listconstruct'
2651
2652	diffop = gen.get_specdictionary_HMT(config, H='DIFFOP_')
2653	diffvar = gen.get_specdictionary_HMT(config, H='DIFFVAR_')
2654
2655	if debug:
2656	if diffop is not None:
2657	print " 'op' differences ________"
2658	gen.printing_dictionary(diffop)
2659	if diffvar is not None:
2660	print " 'var' differences ________"
2661	gen.printing_dictionary(diffvar)
2662
2663	# Getting variable characteristics from each model/experiment
2664	idir1 = config['ofold'] + '/' + mods[0] + '/' + exps[0]
2665	idir2 = config['ofold'] + '/' + mods[1] + '/' + exps[1]
2666
2667	fvarcomp = 'varcompute.inf'
2668	if reprj[0]: fvc1 = 'reproj_' + fvarcomp
2669	else: fvc1 = fvarcomp
2670
2671	if reprj[1]: fvc2 = 'reproj_' + fvarcomp
2672	else: fvc2 = fvarcomp
2673
2674	files1, testfiles1, allcompvar1, Nvar1= read_varcomp_file(idir1 + '/' + fvc1)
2675	files2, testfiles2, allcompvar2, Nvar2= read_varcomp_file(idir2 + '/' + fvc2)
2676
2677	diffops = {}
2678	if diffop is not None:
2679	# Getting 'op' differences characteristics
2680	doops, opvar, indivopvar, varglobalp = get_differences_var(diffop,'op',debug)
2681
2682	Ndiffop = 0
2683	for vn in opvar.keys():
2684	vnops = opvar[vn]
2685	for op in vnops:
2686	vnop = vn + '_' + op
2687	if not allcompvar1.has_key(vnop):
2688	print warnmsg
2689	print ' ' + fname + ": DIFFOP variable+operation '" + vnop + \
2690	"' in '" + mods[0] + '/' + exps[0] + "' was not computed !!"
2691	print ' skipping it!'
2692	print ' variables computed:', allcompvar1
2693	break
2694	if not allcompvar2.has_key(vnop):
2695	print warnmsg
2696	print ' ' + fname + ": DIFFOP variable+operation '" + vnop + \
2697	"' in '" + mods[1] + '/' + exps[1] + "' was not computed !!"
2698	print ' skipping it!'
2699	print ' variables computed:', allcompvar2
2700	break
2701	vals1 = allcompvar1[vnop]
2702	vals2 = allcompvar2[vnop]
2703
2704	headerf1 = vals1[0]
2705	headerf2 = vals2[0]
2706
2707	diffops[vnop] = [mods[0],mods[1],exps[0],exps[1],headerf1,headerf2]
2708
2709	Ndiffop = Ndiffop + 1
2710	else:
2711	Ndiffop = 0
2712
2713	diffvars = {}
2714	if diffvar is not None:
2715	# Getting 'var' differences characteristics
2716	doops,opvar,indivopvar,varglobalp = get_differences_var(diffvar,'var',debug)
2717
2718	Ndiffvar = 0
2719	for vn in opvar.keys():
2720	vnops = opvar[vn]
2721	for op in vnops:
2722	vnop = vn + '_' + op
2723	if not allcompvar1.has_key(vnop):
2724	print warnmsg
2725	print ' ' + fname + ": DIFFVAR variable+operation '" + vnop + \
2726	"' in '" + mods[0] + '/' + exps[0] + "' was not computed !!"
2727	print ' skipping it!'
2728	print ' variables computed:', allcompvar1
2729	break
2730	if not allcompvar2.has_key(vnop):
2731	print warnmsg
2732	print ' ' + fname + ": DIFFVAR variable+operation '" + vnop + \
2733	"' in '" + mods[1] + '/' + exps[1] + "' was not computed !!"
2734	print ' skipping it!'
2735	print ' variables computed:', allcompvar2
2736	break
2737	vals1 = allcompvar1[vnop]
2738	vals2 = allcompvar2[vnop]
2739
2740	headerf1 = vals1[0]
2741	headerf2 = vals2[0]
2742
2743	diffvars[vnop] = [mods[0],mods[1],exps[0],exps[1],headerf1,headerf2]
2744
2745	Ndiffvar = Ndiffvar + 1
2746	else:
2747	Ndiffvar = 0
2748
2749	if diffop is not None:
2750	Sopdiffs = gen.dictKeysVals_stringList(diffops)
2751	else:
2752	Sopdiffs = 'none'
2753	# Outwritting the op diffferences file to avoid next time
2754	objf = open(odir + '/diffop.inf', 'w')
2755	objf.write('differences: ' + Sopdiffs + '\n')
2756	objf.write('Ndiff: ' + str(Ndiffop) + '\n')
2757	objf.close()
2758
2759	if diffvar is not None:
2760	Svardiffs = gen.dictKeysVals_stringList(diffvars)
2761	else:
2762	Svardiffs = 'none'
2763	# Outwritting the var diffferences file to avoid next time
2764	objf = open(odir + '/diffvar.inf', 'w')
2765	objf.write('differences: ' + Svardiffs + '\n')
2766	objf.write('Ndiff: ' + str(Ndiffvar) + '\n')
2767	objf.close()
2768
2769	return diffops, diffvars, Ndiffop, Ndiffvar
2770
2771	def read_diff_file(difffile):
2772	""" Function to read the file with the information about the differences
2773	difffile= file with the information
2774	"""
2775	fname = 'read_diff_file'
2776
2777	if not os.path.isfile(difffile):
2778	print errormsg
2779	print ' ' + fname + ": differences file '" + difffile + "' does not exist !!"
2780	quit(-1)
2781
2782	objf = open(difffile, 'r')
2783	for line in objf:
2784	if line[0:1] != '#' and len(line) > 1:
2785	values = line.replace('\n','').split(' ')
2786	if values[0] == 'differences:':
2787	if values[1] != 'none':
2788	diffs = gen.stringList_dictKeysVals(values[1])
2789	else:
2790	diffs = {}
2791	elif values[0] == 'Ndiff:':
2792	Ndiffs = int(values[1])
2793
2794	objf.close()
2795
2796	return diffs, Ndiffs
2797
2798	def compute_op_diffs(config, diffallop, odir, diffscr, debug):
2799	""" Function to compute operation differences
2800	config= experiment configuration
2801	diffallop= dictonary with the differences to compute
2802	odir= output folder
2803	diffscr= wether it should be done from the scratch
2804	"""
2805	fname = 'compute_op_diffs'
2806
2807	# output folder
2808	ofold = config['ofold']
2809
2810	# Home of the python scripts
2811	pyH = config['pyHOME']
2812
2813	# opsur = operation surnames
2814	opsur = gen.get_specdictionary_HMT(config, H='opsur_',M='',T='')
2815	Opsurs = {}
2816	for opk in opsur:
2817	opn = opk.split('_')[1]
2818	vals = opsur[opk]
2819	Opsurs[opn] = vals.split(':')
2820
2821	# Getting in working dir
2822	os.chdir(odir)
2823
2824	# CF dimensions
2825	CFdims = ['lon', 'lat', 'pres', 'time']
2826
2827	# Trakcing file
2828	trkobjf = open(odir + '/all_diffop.inf', 'a')
2829
2830	for vnop in diffallop.keys():
2831	vn = vnop.split('_')[0]
2832	op = vnop.split('_')[1]
2833	opS = op.replace('+','_')
2834
2835	# `p' header add related to 'pinterp'
2836	if opS.find('pinterp') != -1: SgP = 'p'
2837	else: SgP = ''
2838
2839	# Variable name within the file
2840	vninF = varnoper(vn, op, Opsurs)
2841
2842	vals = diffallop[vnop]
2843	mod1 = vals[0]
2844	mod2 = vals[1]
2845	exp1 = vals[2]
2846	exp2 = vals[3]
2847	headerf1 = vals[4]
2848	headerf2 = vals[5]
2849	modexpdiff = mod2 + '/' + exp2 + ' - ' + mod1 + '/' + exp1
2850	modexpdiffS = mod2 + '-' + exp2 + '_' + mod1 + '-' + exp1
2851
2852	ifile1 = ofold + '/' + mod1 + '/' + exp1 + '/' + vn + '_' + headerf1 + SgP + \
2853	'_' + opS + '.nc'
2854	ifile2 = ofold + '/' + mod2 + '/' + exp2 + '/' + vn + '_' + headerf2 + SgP + \
2855	'_' + opS + '.nc'
2856	difffilen = odir + '/' + vn + '_diffop_' + modexpdiffS + '_' + opS + '.nc'
2857
2858	if diffscr:
2859	sout = sub.call('rm ' + difffilen + ' >& /dev/null', shell=True)
2860
2861	if not os.path.isfile(difffilen):
2862	if debug:
2863	print " Computing operation difference '" + op + "' of '" + vn + \
2864	"' for: '" + modexpdiff
2865	print " var in file: '" + vninF + "'"
2866
2867	dfiles = 'add\|' + ifile2 + '\|' + vninF + ',sub\|' + ifile1 + '\|' + vninF
2868	values = '\|'.join(CFdims) + '@' + dfiles
2869	pyins = 'python ' + pyH + "/nc_var.py -o compute_opersvarsfiles -S '" + \
2870	values + "' -v " + vninF
2871
2872	try:
2873	with gen.Capturing() as output:
2874	ncvar.compute_opersvarsfiles(values,vninF)
2875	except:
2876	print errmsg
2877	print 'ncvar.compute_opersvarsfiles(' + values + ',' + vninF + ')'
2878	for s1out in output: print s1out
2879	quit(-1)
2880
2881	sout = sub.call('mv opersvarsfiles_'+vninF+'.nc ' + difffilen, shell=True)
2882
2883	# keeping all differences
2884	trkobjf.write('\n')
2885	trkobjf.write("#" + vn + ' ' + op + '\n')
2886	trkobjf.write( pyins + '\n')
2887
2888	trkobjf.close()
2889
2890	return
2891
2892	def compute_diff_statistics(Ecnf, filen, ofold, varn, Opers, scr, db):
2893	""" Function to compute statistics from var diff files
2894	Ecnf= dictionary with the configuration of the experiment
2895	filen= variable difference file to use
2896	ofold= directory to write the output files
2897	varn= name of the variable
2898	Opers= kind of operation: (as possible multiple consecutive combination of operations separated by '+'
2899	[calc1]+[calc2] will compute first [calc1] and then [calc2]
2900	acc: temporal accumulated values
2901	diff: differences between models
2902	direct: no statistics
2903	last: last temporal value
2904	Lmean: latitudinal mean values
2905	Lsec: latitudinal section (latitudinal value must be given, [var]@[lat])
2906	lmean: longitudinal mean values
2907	lsec: longitudinal section (longitudinal value must be given, [var]@[lat])
2908	pinterp: pressure interpolation (to the given $plevels)
2909	tmean: temporal mean values
2910	tstd: temporal standard deviation values
2911	tturb: Taylor's turbulence decomposition value (x - <x>) for time
2912	tvar: temporal variance values
2913	xmean: x-axis mean values
2914	xvar: x-axis variance values
2915	ymean: y-axis mean values
2916	zsum: vertical aggregated values
2917	scr= should it be done from the scratch?
2918	"""
2919	fname='compute_diff_statistics'
2920
2921	# Full header of file
2922	Lfilen = len(filen)
2923	Hfile = filen[0:Lfilen-3]
2924
2925	# CF dimension-variables name
2926	OpersS = '+'.join(Opers)
2927	if OpersS.find('pinterp') != -1 or filen.find('pinterp') != -1:
2928	varnCFs = ['lon', 'lat', 'pres', 'time']
2929	else:
2930	varnCFs = ['lon', 'lat', 'time']
2931
2932	# Variable-dimension CF dictionary
2933	CFdimvardict = {'lon': 'lon', 'lat': 'lat', 'pres': 'pres', 'time': 'time'}
2934
2935	# Python scripts HOME
2936	pyH = Ecnf['pyHOME']
2937
2938	# opsur = operation surnames
2939	opsur = gen.get_specdictionary_HMT(Ecnf, H='opsur_',M='',T='')
2940	Opsurs = {}
2941	for opk in opsur:
2942	opn = opk.split('_')[1]
2943	vals = opsur[opk]
2944	Opsurs[opn] = vals.split(':')
2945
2946	# Getting in working dir
2947	os.chdir(ofold)
2948
2949	# File to keep track of all the statistics of var differences
2950	otrackf = open(ofold + '/all_vardiffstatistics.inf', 'a')
2951
2952	if db:
2953	print " computing statistics of diff variable: '", varn ,"' operation '"+ \
2954	OpersS + "' using file '" + filen + " ...'"
2955
2956	for op in Opers:
2957	# File name from previous operations, and name of the variable within the
2958	# file (some operations change the name of the variable)
2959	if op == Opers[0]:
2960	Fopers = op
2961	prevfile = filen
2962	vninF = varn
2963	else:
2964	Fopers = Fopers + '_' + op
2965	prevfile = fileon
2966	fileon = Hfile + '_' + Fopers + '.nc'
2967
2968	SvarnCFs = ',' + ','.join(varnCFs)
2969	CFvarnp = vninF + SvarnCFs
2970
2971	if scr:
2972	sout = sub.call('rm ' + fileon + ' >& /dev/null', shell=True)
2973
2974	# Just produce file in case output file does not exist
2975	if not os.path.isfile(fileon):
2976	if db:
2977	print ' ' + fname + ": creation of file '" + fileon + "' ..."
2978	dofile = True
2979	else:
2980	dofile = False
2981
2982	# Variables to be kept in the final file
2983	varkeep = []
2984
2985	if op == 'acc':
2986	# temporal accumulated values
2987	print " " + fname + ": kind '" + op + "' not ready !!"
2988	wrongstats.append(varn + '_' + opers)
2989	break
2990	elif op == 'direct':
2991	# no statistics
2992	if dofile:
2993	sout = sub.call('mv ' + prevfile + ' ' + fileon, shell=True)
2994	elif op == 'last':
2995	# last temporal value
2996	vals='time,-9,0,0'
2997	if dofile:
2998	try:
2999	with gen.Capturing() as output:
3000	ncvar.DataSetSection(vals,prevfile)
3001	except:
3002	print errmsg
3003	print 'DataSetSection('+vals+',', prevfile, ')'
3004	for s1out in output: print s1out
3005	quit(-1)
3006
3007	if db:
3008	for s1out in output: print s1out
3009
3010	hprevfile = prevfile[0:len(prevfile)-3]
3011
3012	sout = sub.call('mv ' + hprevfile + '_time_B-9-E0-I0.nc ' + fileon, \
3013	shell=True)
3014
3015	# Keeping the operations
3016	pyins=pyH + "/nc_var.py -o DataSetSection -S '" + vals + "' -f " + \
3017	prevfile
3018	otrackf.write("\n")
3019	otrackf.write("# " + varn + " " + Fopers + "\n")
3020	otrackf.write(pyins + "\n")
3021
3022	elif op =='Lmean':
3023	# latitudinal mean values
3024	print " " + fname + ": kind '" + op + "' not ready !!"
3025	wrongstats.append(varn + '_' + opers)
3026	break
3027	elif op =='Lsec':
3028	# latitudinal section (latitudinal value must be given, [var]@[lat])
3029	print " " + fname + ": kind '" + op + "' not ready !!"
3030	wrongstats.append(varn + '_' + opers)
3031	break
3032	elif op =='lmean':
3033	# longitudinal mean values
3034	print " " + fname + ": kind '" + op + "' not ready !!"
3035	wrongstats.append(varn + '_' + opers)
3036	break
3037	elif op =='lsec':
3038	# longitudinal section (longitudinal value must be given, [var]@[lon])
3039	print " " + fname + ": kind '" + op + "' not ready !!"
3040	wrongstats.append(varn + '_' + opers)
3041	break
3042	elif op == 'tmean':
3043	# temporal mean values
3044	vals='time\|-1,time,mean,' + ':'.join(varnCFs)
3045	dims = gen.dictvar_listS(varnCFs,CFdimvardict,',','@')
3046	if dofile:
3047	try:
3048	with gen.Capturing() as output:
3049	ncvar.file_oper_alongdims(vals,prevfile,CFvarnp)
3050	except:
3051	print errmsg
3052	print 'file_oper_alongdims('+vals+',', prevfile, ','+CFvarnp+')'
3053	for s1out in output: print s1out
3054	quit(-1)
3055
3056	if db:
3057	for s1out in output: print s1out
3058
3059	sout = sub.call('mv file_oper_alongdims_mean.nc '+fileon, shell=True)
3060
3061	# Keeping the operations
3062	pyins=pyH + "/nc_var.py -o file_oper_alongdims -S '" + vals + \
3063	"' -f " + prevfile + " -v " + CFvarnp
3064	otrackf.write("\n")
3065	otrackf.write("# " + varn + " " + Fopers + "\n")
3066	otrackf.write(pyins + "\n")
3067
3068	# removing dimension variable-dimension 'time'
3069	varnCFs.remove('time')
3070
3071	varkeep.append('timestats')
3072
3073	elif op == 'tstd':
3074	# temporal standard deviation values
3075	vals='time\|-1,time,std,' + ':'.join(varnCFs) + ':' + vdnz
3076	dims = gen.dictvar_listS(varnCFs,CFdimvardict,',','@') + dnz+'@'+vdnz
3077	if dofile:
3078	try:
3079	with gen.Capturing() as output:
3080	ncvar.file_oper_alongdims(vals,prevfile,CFvarnp)
3081	except:
3082	print errmsg
3083	print 'file_oper_alongdims('+vals+',', prevfile, ','+CFvarnp+')'
3084	for s1out in output: print s1out
3085	quit(-1)
3086
3087	if db:
3088	for s1out in output: print s1out
3089
3090	sout = sub.call('mv file_oper_alongdims_std.nc '+fileon, shell=True)
3091
3092	# Keeping the operations
3093	pyins=pyH + "/nc_var.py -o file_oper_alongdims -S '" + vals + \
3094	"' -f " + prevfile + " -v " + CFvarnp
3095	otrackf.write("\n")
3096	otrackf.write("# " + varn + " " + Fopers + "\n")
3097	otrackf.write(pyins + "\n")
3098
3099	# removing dimension variable-dimension 'time'
3100	varnCFs.remove('time')
3101
3102	varkeep.append('timestats')
3103
3104	elif op == 'tturb':
3105	# temporal turbulent values
3106	vals='time\|-1,time,turb,' + ':'.join(varnCFs)
3107	dims = gen.dictvar_listS(varnCFs,CFdimvardict,',','@')
3108	if dofile:
3109	try:
3110	with gen.Capturing() as output:
3111	ncvar.file_oper_alongdims(vals,prevfile,CFvarnp)
3112	except:
3113	print errmsg
3114	print 'file_oper_alongdims('+vals+',', prevfile, ','+CFvarnp+')'
3115	for s1out in output: print s1out
3116	quit(-1)
3117
3118	if db:
3119	for s1out in output: print s1out
3120
3121	sout = sub.call('mv file_oper_alongdims_turb.nc '+fileon, shell=True)
3122
3123	# Keeping the operations
3124	pyins=pyH + "/nc_var.py -o file_oper_alongdims -S '" + vals + \
3125	"' -f " + prevfile + " -v " + CFvarnp
3126	otrackf.write("\n")
3127	otrackf.write("# " + varn + " " + Fopers + "\n")
3128	otrackf.write(pyins + "\n")
3129
3130	varkeep.append('timestats')
3131
3132	elif op == 'tvar':
3133	# temporal variance values (Taylor's turbulence)
3134	vals='time\|-1,time,var,' + ':'.join(varnCFs) + ':' + vdnz
3135	dims = gen.dictvar_listS(varnCFs,CFdimvardict,',','@') + dnz+'@'+vdnz
3136	if dofile:
3137	try:
3138	with gen.Capturing() as output:
3139	ncvar.file_oper_alongdims(vals,prevfile,CFvarnp)
3140	except:
3141	print errmsg
3142	print 'file_oper_alongdims('+vals+',', prevfile, ','+CFvarnp+')'
3143	for s1out in output: print s1out
3144	quit(-1)
3145
3146	if db:
3147	for s1out in output: print s1out
3148
3149	sout = sub.call('mv file_oper_alongdims_var.nc '+fileon, shell=True)
3150
3151	# Keeping the operations
3152	pyins=pyH + "/nc_var.py -o file_oper_alongdims -S '" + vals + \
3153	"' -f " + prevfile + " -v " + CFvarnp
3154	otrackf.write("\n")
3155	otrackf.write("# " + varn + " " + Fopers + "\n")
3156	otrackf.write(pyins + "\n")
3157
3158	# removing dimension variable-dimension 'time'
3159	varnCFs.remove('time')
3160
3161	varkeep.append('timestats')
3162
3163	elif op == 'xmean':
3164	# x-axis mean values
3165	vals='lon\|-1,lon,mean,' + ':'.join(varnCFs)
3166	dims = gen.dictvar_listS(varnCFs,CFdimvardict,',','@')
3167	if dofile:
3168	try:
3169	with gen.Capturing() as output:
3170	ncvar.file_oper_alongdims(vals,prevfile,CFvarnp)
3171	except:
3172	print errmsg
3173	print 'file_oper_alongdims('+vals+',', prevfile, ','+CFvarnp+')'
3174	for s1out in output: print s1out
3175	quit(-1)
3176
3177	if db:
3178	for s1out in output: print s1out
3179
3180	sout = sub.call('mv file_oper_alongdims_mean.nc '+fileon, shell=True)
3181
3182	# Keeping the operations
3183	pyins=pyH + "/nc_var.py -o file_oper_alongdims -S '" + vals + \
3184	"' -f " + prevfile + " -v " + CFvarnp
3185	otrackf.write("\n")
3186	otrackf.write("# " + varn + " " + Fopers + "\n")
3187	otrackf.write(pyins + "\n")
3188
3189	# removing dimension variable-dimension 'lon'
3190	varnCFs.remove('lon')
3191
3192	varkeep.append('lonstats')
3193
3194	elif op == 'xvar':
3195	# x-axis variance values
3196	vals='lon\|-1,lon,var,' + ':'.join(varnCFs)
3197	dims = gen.dictvar_listS(varnCFs,CFdimvardict,',','@')
3198	if dofile:
3199	try:
3200	with gen.Capturing() as output:
3201	ncvar.file_oper_alongdims(vals,prevfile,CFvarnp)
3202	except:
3203	print errmsg
3204	print 'file_oper_alongdims('+vals+',', prevfile, ','+CFvarnp+')'
3205	for s1out in output: print s1out
3206	quit(-1)
3207
3208	if db:
3209	for s1out in output: print s1out
3210
3211	sout = sub.call('mv file_oper_alongdims_var.nc '+fileon, shell=True)
3212
3213	# Keeping the operations
3214	pyins=pyH + "/nc_var.py -o file_oper_alongdims -S '" + vals + \
3215	"' -f " + prevfile + " -v " + CFvarnp
3216	otrackf.write("\n")
3217	otrackf.write("# " + varn + " " + Fopers + "\n")
3218	otrackf.write(pyins + "\n")
3219
3220	# removing dimension variable-dimension 'lon'
3221	varnCFs.remove('lon')
3222
3223	varkeep.append('lonstats')
3224
3225	elif op == 'ymean':
3226	# y-axis mean values
3227	vals='lat\|-1,lat,mean,' + ':'.join(varnCFs)
3228	dims = gen.dictvar_listS(varnCFs,CFdimvardict,',','@')
3229	if dofile:
3230	try:
3231	with gen.Capturing() as output:
3232	ncvar.file_oper_alongdims(vals,prevfile,CFvarnp)
3233	except:
3234	print errmsg
3235	print 'file_oper_alongdims('+vals+',', prevfile, ','+CFvarnp+')'
3236	for s1out in output: print s1out
3237	quit(-1)
3238
3239	if db:
3240	for s1out in output: print s1out
3241
3242	sout = sub.call('mv file_oper_alongdims_mean.nc '+fileon, shell=True)
3243
3244	# Keeping the operations
3245	pyins=pyH + "/nc_var.py -o file_oper_alongdims -S '" + vals + \
3246	"' -f " + prevfile + " -v " + CFvarnp
3247	otrackf.write("\n")
3248	otrackf.write("# " + varn + " " + Fopers + "\n")
3249	otrackf.write(pyins + "\n")
3250
3251	# removing dimension variable-dimension 'lat'
3252	varnCFs.remove('lat')
3253
3254	varkeep.append('latstats')
3255
3256	elif op == 'zsum':
3257	# vertical aggregated values
3258	print " " + fname + ": kind '" + op + "' not ready !!"
3259	wrongstats.append(varn + '_' + opers)
3260	break
3261	else:
3262	print errmsg
3263	print ' ' + fname + ": operation '" + op + "' not ready !!"
3264	quit(-1)
3265
3266	# End of kind of operation
3267
3268	# Variable name in file (vninF) changed due to operation
3269	# but only if previous operation does not the same 'statistic'
3270	chvn = gen.dictionary_key_list(Opsurs, op)
3271	if chvn is not None:
3272	oldvninF = vninF
3273	vninF = vninF + chvn
3274	CFvarnp = CFvarnp.replace(oldvninF,vninF)
3275
3276	if dofile:
3277	if len(varkeep) > 0:
3278	varkeepS = ',' + ','.join(varkeep)
3279	else:
3280	varkeepS = ''
3281
3282	# Adding such CF dimension variables which might not be in the
3283	# post-operation file
3284	try:
3285	with gen.Capturing() as output:
3286	varinfile = ncvar.ivars(fileon)
3287	except:
3288	print errmsg
3289	print 'ivars('+fileon+')'
3290	for s1out in output: print s1out
3291	quit(-1)
3292
3293	for CFvn in varnCFs:
3294	if not gen.searchInlist(varinfile, CFvn):
3295	if db:
3296	print ' ' + fname + ": recupering CF variable '" + CFvn + "'"
3297	try:
3298	with gen.Capturing() as output:
3299	ncvar.fvaradd(prevfile+','+CFvn, fileon)
3300	except:
3301	print errmsg
3302	print 'fvaradd(' + prevfile + ',' + CFvn +', ' + fileon + ')'
3303	for s1out in output: print s1out
3304	quit(-1)
3305
3306	if db:
3307	for s1out in output: print s1out
3308
3309	totalvarkeeps = CFvarnp + ',' + ','.join(varnCFs) + varkeepS
3310	try:
3311	with gen.Capturing() as output:
3312	oclean = ncvar.cleaning_varsfile(totalvarkeeps,fileon)
3313	except:
3314	print errmsg
3315	print 'cleaning_varsfile('+totalvarkeeps+','+fileon+')'
3316	for s1out in output: print s1out
3317	quit(-1)
3318
3319	if db:
3320	for s1out in output: print s1out
3321
3322	# End of operations
3323	otrackf.close()
3324
3325	return
3326
3327	def compute_var_diffs(config, diffallvar, odir, diffscr, debug):
3328	""" Function to compute variable differences
3329	config= experiment configuration
3330	diffallvar= dictonary with the differences to compute
3331	odir= output folder
3332	diffscr= wether it should be done from the scratch
3333	"""
3334	fname = 'compute_var_diffs'
3335
3336	# output folder
3337	ofold = config['ofold']
3338
3339	# Home of the python scripts
3340	pyH = config['pyHOME']
3341
3342	# opsur = operation surnames
3343	opsur = gen.get_specdictionary_HMT(config, H='opsur_',M='',T='')
3344	Opsurs = {}
3345	for opk in opsur:
3346	opn = opk.split('_')[1]
3347	vals = opsur[opk]
3348	Opsurs[opn] = vals.split(':')
3349
3350	# Getting in working dir
3351	os.chdir(odir)
3352
3353	# CF dimensions
3354	CFdims = ['lon', 'lat', 'pres', 'time']
3355
3356	# Trakcing file
3357	trkobjf = open(odir + '/all_diffvar.inf', 'a')
3358
3359	for vnop in diffallvar.keys():
3360	vn = vnop.split('_')[0]
3361	op = vnop.split('_')[1]
3362	opSfinal = op.replace('+','_')
3363
3364	lop = op.split('+')
3365	# `p' header add related to 'pinterp'
3366	if op.find('pinterp') != -1:
3367	print warnmsg
3368	print ' ' + fname + ': there is the vertical interpolation as operation'
3369	print ' variable differences will be computed from that files where '
3370	print ' vertical interpolation was just done'
3371	ipop = lop.index('pinterp')
3372	opS = '_' + '_'.join(lop[0:ipop+1])
3373	doops = lop[ipop+1:]
3374	SgP = 'p'
3375	vninF = varnoper(vn, '+'.join(lop[0:ipop+1]), Opsurs)
3376	print ' operations to perform:', doops
3377	else:
3378	SgP = ''
3379	opS = ''
3380	doops = lop
3381	vninF = vn
3382
3383	Sdoops = '+'.join(doops)
3384
3385	# Variable name within the file
3386
3387	vals = diffallvar[vnop]
3388	mod1 = vals[0]
3389	mod2 = vals[1]
3390	exp1 = vals[2]
3391	exp2 = vals[3]
3392	headerf1 = vals[4]
3393	headerf2 = vals[5]
3394	modexpdiff = mod2 + '/' + exp2 + ' - ' + mod1 + '/' + exp1
3395	modexpdiffS = mod2 + '-' + exp2 + '_' + mod1 + '-' + exp1
3396
3397	ifile1 = ofold + '/' + mod1 + '/' + exp1 + '/' + vn + '_' + headerf1 + SgP + \
3398	opS + '.nc'
3399	ifile2 = ofold + '/' + mod2 + '/' + exp2 + '/' + vn + '_' + headerf2 + SgP + \
3400	opS + '.nc'
3401	difffilen = odir + '/' + vn + '_diffvar_' + modexpdiffS + opS + '.nc'
3402	difffinalfilen = odir + '/' + vn + '_diffvar_' + modexpdiffS + opS + '_' + \
3403	Sdoops + '.nc'
3404
3405	if diffscr:
3406	sout = sub.call('rm ' + difffilen + ' >& /dev/null', shell=True)
3407	sout = sub.call('rm ' + difffinalfilen + ' >& /dev/null', shell=True)
3408
3409	if not os.path.isfile(difffinalfilen):
3410	if debug:
3411	print " Computing variable difference of '" + vn + "' for: '" + \
3412	modexpdiff + "' and performing operation '" + Sdoops + "'"
3413	print " var in file: '" + vninF + "'"
3414
3415	dfiles = 'add\|' + ifile2 + '\|' + vninF + ',sub\|' + ifile1 + '\|' + vninF
3416	values = '\|'.join(CFdims) + '@' + dfiles
3417	pyins = 'python ' + pyH + "/nc_var.py -o compute_opersvarsfiles -S '" + \
3418	values + "' -v " + vninF
3419
3420	try:
3421	with gen.Capturing() as output:
3422	ncvar.compute_opersvarsfiles(values,vninF)
3423	except:
3424	print errmsg
3425	print 'ncvar.compute_opersvarsfiles(' + values + ',' + vninF + ')'
3426	for s1out in output: print s1out
3427	quit(-1)
3428
3429	sout = sub.call('mv opersvarsfiles_'+vninF+'.nc ' + difffilen, shell=True)
3430
3431	# keeping all differences
3432	trkobjf.write('\n')
3433	trkobjf.write("#" + vn + ' ' + op + '\n')
3434	trkobjf.write( pyins + '\n')
3435
3436	# Computing statisitcs
3437	compute_diff_statistics(config, difffilen, odir, vninF, doops, diffscr, \
3438	debug)
3439
3440	trkobjf.close()
3441
3442	return
3443
3444
3445	def draw_diff_plots(config, plots, odir, allvarcomp, kdiff, figscr, debug):
3446	""" Function to draw all plots
3447	config= Configuration of the experiment
3448	plots= dictionary with the plots
3449	odir= output experiment folder
3450	allvarcomp= dictionary with all the variables to compute and their information
3451	kdiff= kind of differences:
3452	'diffop': plots from operation differences
3453	'diffvar': plots from variable differences
3454	figscr= whether figures should be done from the scratch or not
3455
3456	* Plot as
3457	{[kplot]} = [varn1]\|[op1]#[varn2]\|[op2]#[...[varnN]\|[opN]], ...
3458	[kplot] ___
3459	diffmap2Dsfc: 2D map of surface differences values of 1 variable
3460	diffmap2Dz: 2D map of 3D differences values of 1 variable
3461	map2Dsfc: 2D map of surface values of 1 variable
3462	map3D: 2D map of 3D values of 1 variable
3463	shadconthovmsfc: Hovmoeller diagrams of 2 variable at the surface in shadow and the other in contourn
3464	shadcont2Dsfc: 2D map of shadow (1st variable) and countour (2nd variable) [stvar1]#[stvar2]
3465	shadcont2Dzsec: 2D map of vertical section of 2 variables one in shadow and the other in contourn
3466	[varn]
3467	variable
3468	[op]
3469	'+' separated list of operations
3470	in figures with more than 1 variable, use '#' to separate the [varn]\|[op]
3471	"""
3472	fname = 'draw_diff_plots'
3473
3474	os.chdir(odir)
3475
3476	# Dictionary with the operations with surnames for the operated variable
3477	opersurnames = {}
3478	opsur = gen.get_specdictionary_HMT(config, H='opsur_',M='',T='')
3479	for opsr in opsur.keys():
3480	opn = opsr.split('_')[1]
3481	vls = opsur[opsr].split(':')
3482	opersurnames[opn] = vls
3483
3484	# time values
3485	# Units time for the plots
3486	rd = config['CFreftime']
3487	tunits = config['CFunitstime'] + '!since!' + rd[0:4] + '-' + rd[4:6] + '-' + \
3488	rd[6:8] + '!' + rd[8:10] + ':' + rd[10:12] + ':' + rd[12:14]
3489	# time ticks kind
3490	tkind = config['timekind']
3491	# time ticks format
3492	tfmt = config['timefmt']
3493	# time axis label
3494	tlab = config['timelabel']
3495	timevals = [tunits, tkind, tfmt, tlab]
3496
3497	if kdiff == 'diffop':
3498	specplotkeyn = 'specificdiffopplot'
3499	elif kdiff == 'diffvar':
3500	specplotkeyn = 'specificdiffvarplot'
3501	else:
3502	print errmsg
3503	print ' ' + fname + ": differences kind '" + kdiff + "' not ready !!"
3504	quit(-1)
3505
3506	# Dictionary of plot specificities
3507	plotspecifics = {}
3508	if config.has_key(specplotkeyn):
3509	# [minval]: minimum value
3510	# [maxval]: minimum value
3511	# [colorbar]: name of the colorbar (from matplotlib) to use
3512	# [cntformat]: format of the contour labels
3513	# [colorcnt]: color for the countor lines
3514	plotspecs = config[specplotkeyn].split(':')
3515	for pltspc in plotspecs:
3516	pltvls = pltspc.split('\|')
3517	vn = pltvls[0]
3518	op = pltvls[1]
3519	fn = pltvls[2]
3520	plotspecifics[fn + '_' + vn + '_' + op] = pltvls[3:]
3521	if debug:
3522	print 'Specific values for plots _______'
3523	gen.printing_dictionary(plotspecifics)
3524
3525	# Kind of figures
3526	kindfigure = config['kindfig']
3527
3528	# Graphical labels and configuration
3529	modgc, expgc = graphical_characteristics(config, debug)
3530
3531	# Map value
3532	mapvalue = config['mapval']
3533
3534	# pythone scripts HOME
3535	pyHOME = config['pyHOME']
3536
3537	# Title-text of operations
3538	opexplained = {}
3539	optits = config['titleoperations'].split(':')
3540	for optit in optits:
3541	opn = optit.split('\|')[0]
3542	opt = optit.split('\|')[1]
3543	opexplained[opn] = opt
3544	if debug:
3545	print 'Titles for operations _______'
3546	gen.printing_dictionary(opexplained)
3547
3548	for kplot in plots.keys():
3549	varsplt = plots[kplot]
3550	for varplt in varsplt:
3551	if debug:
3552	print " printing '" + kplot + "' var ':" + varplt + "'..."
3553	varops = varplt.split('#')
3554
3555	# CF variables in plot
3556	CFvarsplot = []
3557	# Files in plot
3558	filesplot = []
3559	# Variables in plot within the files
3560	varsplot = []
3561	# Dims in figure
3562	dimsplot = []
3563	# pictoric values in figure
3564	pictplot = []
3565	# Name of the figure
3566	figname = ''
3567	# Title of the figure
3568	titfigure = ''
3569
3570	ivp = 0
3571	for varop in varops:
3572	vn = varop.split('\|')[0]
3573	op = varop.split('\|')[1]
3574
3575	# CF variables in plot
3576	CFvarsplot.append(vn)
3577
3578	vnopS = vn + '_' + op
3579	if not allvarcomp.has_key(vnopS):
3580	print errmsg
3581	print ' ' + fname + ": no file for variable-operation '" + \
3582	vnopS + "' !!"
3583	vopvals = allvarcomp[vnopS]
3584	mod1 = vopvals[0]
3585	mod2 = vopvals[1]
3586	mod3 = vopvals[2]
3587	mod4 = vopvals[3]
3588	headf = vopvals[4]
3589
3590	modv = modgc[mod1]
3591	expv = expgc[mod1+'/'+exp1]
3592	modlab1 = modv.label
3593	explab1 = expv.label
3594	modv = modgc[mod2]
3595	expv = expgc[mod2+'/'+exp2]
3596	modlab2 = modv.label
3597	explab2 = expv.label
3598
3599	modexpdiff = explab2 + '-' + explab1
3600	modexpdiffS = mod2 + '-' + exp2 + '_' + mod1 + '-' + exp1
3601
3602	difffilen = odir + '/' + vn + '_' + kdiff + '_' + modexpdiffS + \
3603	'_' + op.replace('+','_') + '.nc'
3604
3605	filesplot.append(difffilen)
3606	# Do we have processed the given difference?
3607	if not os.path.isfile(difffilen):
3608	print warnmsg
3609	print " " + fname + ": there is no file for '" + kdiff + \
3610	"' difference '" + varop + "' skiping it !!"
3611	break
3612
3613	# Name of the variable inside the file
3614	vnsur = varnoper(vn, op, opersurnames)
3615	varsplot.append(vnsur)
3616
3617	# Dimensions in file
3618	try:
3619	with gen.Capturing() as output:
3620	dims = ncvar.idims(difffilen)
3621	except:
3622	print errmsg
3623	print 'ncvar.idims('+difffilen+')'
3624	for s1out in output: print s1out
3625	quit(-1)
3626
3627	dimsplot.append(dims)
3628
3629	# pictoric values for the figure
3630	Sfivaop = kplot + '_' + vn + '_' + op
3631	if plotspecifics.has_key(Sfivaop):
3632	pictvals = plotspecifics[Sfivaop]
3633	else:
3634	Vvals = gen.variables_values(vn)
3635	pictvals = [Vvals[2], Vvals[3], Vvals[6], '%g', 'fixsigc','black']
3636
3637	pictplot.append(pictvals)
3638
3639	# Header of the name of the figure
3640	if ivp == 0:
3641	figname = kplot + '_' + vn + '_' + kdiff + '_' + modexpdiffS + \
3642	'_' + op.replace('+','-')
3643	else:
3644	figname = figname + '-' + vn + '_' + op.replace('+','-')
3645
3646	# Title of the figure:
3647	if mod1 == mod2:
3648	modexptit = mod1 + '!' + exp2 + '-' + exp1
3649	else:
3650	modexptit = modexpdiff
3651
3652	ivp = ivp + 1
3653	# End of variable-operation
3654	figname = figname + '.' + kindfigure
3655
3656	# Title of figure
3657	titfigure = create_figure_title(modexptit, kdiff, varops, opexplained)
3658
3659	if len(titfigure) > 80:
3660	print warnmsg
3661	print ' ' + fname + ": figure title '" + titfigure.replace('!', ' ')+\
3662	"' larger than 80 characters (actually", len(titfigure), ') !!'
3663	print " simplifying it to '" + modexptit.replace('!',' ') + '!' + \
3664	kdiff + '!' + vnS + "'"
3665	titfigure = modexptit.replace(' ','!') + '!' + kdiff + '!' + vnS
3666
3667	draw_plot(kplot, CFvarsplot, filesplot, varsplot, dimsplot, pictplot, \
3668	figname, titfigure, kindfigure, mapvalue, timevals, expgc, odir, \
3669	pyHOME, figscr, debug)
3670
3671	# End of variables-operations
3672
3673	# End of kind of plots
3674
3675	return
3676
3677	def reproject_modexp(mod,exp,config,CDOops,fscr,fadd,debug):
3678	""" Function to re-project and re-compute files from a given model/experiment
3679	mod= model
3680	exp= experiment
3681	config= configuration
3682	CDOops= dictionary with the CDO operators to use for reproject and for which variable
3683	fscr= re-projected files to compute from the scratch
3684	fadd= re-projected files to be added
3685	"""
3686	fname = 'reproject_modexp'
3687
3688	odir = config['ofold'] + '/' + mod + '/' + exp
3689
3690	# NON re-projectable operations
3691	opNOreproj = config['NOreprojops'].split(':')
3692
3693	# Need to pass to analyze all the data?
3694	inffiles = ['reproj_varcompute.inf', 'all_reproj_computevars.inf']
3695	if fscr:
3696	for inff in inffiles:
3697	if debug: print " removing information file '" + inff + "' ..."
3698	ins = 'rm ' + odir + '/' + inff + ' >& /dev/null'
3699	sout = sub.call(ins, shell=True)
3700
3701	objf = open(odir+'/all_reproj_computevars.inf','w')
3702	objf.write("## Computation of reprojected variables \n")
3703	objf.close()
3704
3705	varcompf = odir + '/reproj_varcompute.inf'
3706	if fadd:
3707	sub.call('rm ' + varcompf +' >& /dev/null', shell=True)
3708
3709	if not os.path.isfile(varcompf):
3710	# Getting original computing information
3711	origvarcompf = odir + '/varcompute.inf'
3712	origfiles,origtestfiles,origallcompvar, Nvar = read_varcomp_file(origvarcompf)
3713
3714	Svarcompute = ''
3715	allcompvar = dict(origallcompvar)
3716	ivop = 0
3717	for vnop in origallcompvar.keys():
3718	vn = vnop.split('_')[0]
3719	op = vnop.split('_')[1]
3720	values = allcompvar[vnop]
3721	values[0] = 'reproj-' + values[0]
3722	allcompvar[vnop] = values
3723
3724	# Building strings from list values
3725	model = values[1]
3726	diag = values[2]
3727	if model is not None:
3728	if type(model) == type(list([1, 2])):
3729	Smodel = ':'.join(model)
3730	elif type(model) == type('A'):
3731	Smodel = model
3732	else:
3733	Smodel = 'None'
3734	if diag is not None:
3735	if type(diag) == type(list([1, 2])):
3736	Sdiag = ':'.join(diag)
3737	elif type(diag) == type('A'):
3738	Sdiag = diag
3739	else:
3740	Sdiag = 'None'
3741
3742	Svc = vn + '\|' + op + '\|' + values[0] + '\|' + Smodel + '\|' + Sdiag + \
3743	'\|' + values[3]
3744	if ivop == 0:
3745	Svarcompute = Svc
3746	else:
3747	Svarcompute = Svarcompute + ',' + Svc
3748	ivop = ivop + 1
3749
3750	origstrfiles = gen.dictKeysVals_stringList(origfiles)
3751	origstrtestfiles = gen.dictKeysVals_stringList(origtestfiles)
3752
3753	# Outwritting the varcompute to avoid next time (if it is not filescratch!)
3754	objf = open(odir + '/reproj_varcompute.inf', 'w')
3755	objf.write('files: ' + origstrfiles + '\n')
3756	objf.write('testfiles: ' + origstrtestfiles + '\n')
3757	objf.write('varcompute: ' + Svarcompute + '\n')
3758	objf.write('itotv: ' + str(ivop) + '\n')
3759	objf.close()
3760	else:
3761	print warnmsg
3762	print ' ' + main + ": getting variables to compute already from file !!"
3763	origfiles, origtestfiles, allcompvar, Nvar = read_varcomp_file(varcompf)
3764
3765	# End of avoiding to repeat all the experiment search
3766
3767	print " For experiment '" + exp + "' is required to re-project:", Nvar, \
3768	"variables"
3769
3770	if dbg:
3771	print 'Variables to re-project _______'
3772	gen.printing_dictionary(allcompvar)
3773
3774	# Re-projection variables
3775	print " Reprojecting variables ..."
3776	os.chdir(odir)
3777
3778	# Keeping track of reprojection
3779	trckobjf = open(odir + '/all_reproj_computevars.inf','a')
3780
3781	for vnop in allcompvar.keys():
3782	vn = vnop.split('_')[0]
3783	op = vnop.split('_')[1]
3784	values = allcompvar[vnop]
3785
3786	# Only are reprojected that files which do contain 'lon' and 'lat',
3787	# otherways, it is only re-projected the original file
3788	inNO, iopNO, iop = gen.list_coincidences(opNOreproj,op.split('+'))
3789
3790	cdoop = gen.dictionary_key_list(CDOops,vn).split('_')[1]
3791	if cdoop is None:
3792	print errmsg
3793	print ' ' + fname + ": no reprojecting CDO operator found for " + \
3794	" variable '" + vn + "' !!"
3795	print ' available values _______'
3796	gen.printing_dictionary(CDOops)
3797	quit(-1)
3798
3799	# Re-projecting also the 'original' (at least with 'pinterp') file for varsdif
3800	if op.find('pinterp') != -1:
3801	ifileorig = vn +'_'+ values[0].replace('reproj-','') + 'p_pinterp.nc'
3802	ofileorig = vn +'_'+ values[0] + 'p_pinterp.nc'
3803	ofilen = vn +'_'+ values[0] + 'p_pinterp.nc'
3804	else:
3805	ifileorig = vn +'_'+ values[0].replace('reproj-','') + '.nc'
3806	ofileorig = vn +'_'+ values[0] + '.nc'
3807
3808	ifileop = vn+'_'+values[0].replace('reproj-','')+'_'+op.replace('+','_')+'.nc'
3809	ofileop = vn + '_' + values[0] + '_' + op.replace('+','_') + '.nc'
3810
3811	if fscr:
3812	sout = sub.call('rm ' + odir+'/'+ofileorig + ' >& /dev/null', shell=True)
3813	sout = sub.call('rm ' + odir+'/'+ofileop + ' >& /dev/null', shell=True)
3814
3815	ifiles = [odir+'/'+ifileorig, odir+'/'+ifileop]
3816	ofiles = [odir+'/'+ofileorig, odir+'/'+ofileop]
3817
3818	for iif in range(2):
3819	ifile = ifiles[iif]
3820	ofile = ofiles[iif]
3821
3822	# Only projecting original file-variable if there are not 'lon' and 'lat'
3823	doreproj = True
3824	if iif == 1 and len(inNO) > 0:
3825	doreproj = False
3826	print warnmsg
3827	print ' ' + fname + ": operation '" + op + "' can not be " + \
3828	"reprojected due to the lack of 'lon' or 'lat' skipping it !!"
3829	print " NOTE: difference can only be computed as 'variable' " + \
3830	"difference"
3831	print " As operating along the differences from the original" + \
3832	" values from each model"
3833
3834	# Building up resources to compute re-projected statistics
3835	print ' ' + fname + ': computing from the original reprojected file'
3836	ModI = ModelInf(mod, mod, 'lon', 'lat', 'pres', 'time', 'depth', \
3837	'lon', 'lat', 'pres', 'time', 'depth', None)
3838	VarI = VariableInf(vn, values[0], values[2], values[3])
3839	usefs = {values[0]: ofileorig}
3840
3841	compute_statistics(ModI, config, config['ifold'], usefs, \
3842	config['ofold']+'/'+mod+'/'+exp,VarI,pinterp_var(op), op, fscr, dbg)
3843
3844	if not os.path.isfile(ofile) and doreproj:
3845	if debug and ofile == ofiles[0]:
3846	print " reprojecting '" + vn + "' '" + op + "' into '" + \
3847	ofile + "'"
3848
3849	ins = config['cdoHOME'] + '/cdo remap' + cdoop + ',' + \
3850	config['RefProjfile'] + ' ' + ifile + ' ' + ofile
3851	try:
3852	with gen.Capturing() as output:
3853	sout = sub.call(ins, shell=True)
3854	except:
3855	print errmsg
3856	print ' ' + fname + ': $ ' + ins
3857	for s1out in output: print s1out
3858	print sout
3859	quit(-1)
3860	if sout != 0:
3861	print errmsg
3862	print ' ' + fname + ': $ ' + ins
3863	sub.call(ins, shell=True)
3864	quit(-1)
3865
3866	if debug:
3867	for s1out in output: print s1out
3868	print sout
3869
3870	# Keeping track of re-prjection
3871	trckobjf.write('\n')
3872	if iif == 0:
3873	trckobjf.write('# ' + vn + '\n')
3874	else:
3875	trckobjf.write('# ' + vn + ' ' + op + '\n')
3876
3877	trckobjf.write(ins + '\n')
3878
3879	# CFing the reprojected file
3880	try:
3881	with gen.Capturing() as output:
3882	ncvar.CDO_toCF(ofile)
3883	except:
3884	print errmsg
3885	print 'ncvar.CDO_toCF(' + ofile + ')'
3886	for s1out in output: print s1out
3887	quit(-1)
3888
3889	if debug:
3890	for s1out in output: print s1out
3891
3892	trckobjf.close()
3893
3894	return
3895
3896	def allmodexps_listconstruct(config, modexps, odir, debug):
3897	""" Function to create the list of all model-experiments to draw
3898	config= Configuration of the experiment
3899	modexps= list with the models-experiments pairs to use
3900	odir= output differences folder
3901	"""
3902	fname='allmodexps_listconstruct'
3903
3904	allplts = gen.get_specdictionary_HMT(config, H='PLOTALLMODEXP_')
3905
3906	if debug:
3907	if allplts is not None:
3908	print " 'all mod-exps' plots ________"
3909	gen.printing_dictionary(allplts)
3910
3911	# Getting Model-Reference Projection
3912	ModComProj = cnf['RefProj'].split(':')
3913
3914	allmodexp = {}
3915	Nallmodexp = 0
3916	for modexp in modexps:
3917	if debug:
3918	print " " + fname + ": '" + modexp + "' ..."
3919	mod = modexp.split('/')[0]
3920	exp = modexp.split('/')[1]
3921
3922	# Getting variable characteristics from each model/experiment
3923	idir = config['ofold'] + '/' + modexp
3924
3925	reproj = not gen.searchInlist(ModComProj,mod)
3926
3927	fvarcomp = 'varcompute.inf'
3928	if reproj: fvc = 'reproj_' + fvarcomp
3929	else: fvc = fvarcomp
3930
3931	files, testfiles, allcompvar, Nvar= read_varcomp_file(idir + '/' + fvc)
3932
3933	if allplts is not None:
3934	# Getting all model-experiments plots characteristics
3935
3936	for allplt in allplts.keys():
3937	pltn = allplt.split('_')[1]
3938	if debug and modexp == modexps[0]:
3939	print ' ' + fname + ": plot '" + pltn + "' "
3940	varops = allplts[allplt].split(':')
3941	for varop in varops:
3942	vn = varop.split('\|')[0]
3943	op = varop.split('\|')[1]
3944	vnop = vn + '_' + op
3945	if debug and modexp == modexps[0]:
3946	print ' ' + fname + ": varop '" + vnop + "' "
3947	if not allcompvar.has_key(vnop):
3948	print warnmsg
3949	print ' ' + fname + ": variable+operation '" + vnop + \
3950	"' in '" + modexp + "' was not computed !!"
3951	print ' skipping it!'
3952	print ' computed values:', allcompvar.keys()
3953	break
3954
3955	vals = allcompvar[vnop]
3956
3957	headerf = vals[0]
3958	if modexp == modexps[0]:
3959	allmodexp[vnop] = [modexp,headerf]
3960	Nallmodexp = Nallmodexp + 1
3961	else:
3962	prevals = allmodexp[vnop]
3963	prevals = prevals + [modexp,headerf]
3964	allmodexp[vnop] = prevals
3965
3966	else:
3967	Nallmodexp = 0
3968
3969	if allplts is not None:
3970	Sallmodexps = gen.dictKeysVals_stringList(allmodexp)
3971	else:
3972	Sallmodexps = 'none'
3973
3974	# Outwritting the all model-experiment plots file to avoid next time
3975	objf = open(odir + '/allmodexp.inf', 'w')
3976	objf.write('allmodexps: ' + Sallmodexps + '\n')
3977	objf.write('Nallmodexp: ' + str(Nallmodexp) + '\n')
3978	objf.close()
3979
3980	return allmodexp, Nallmodexp
3981
3982	def allmodexps_read(filen):
3983	""" Function to read the file with the information about the all model-experiment plots
3984	filen= file with the information
3985	"""
3986	fname = 'allmodexps_read'
3987
3988	if not os.path.isfile(filen):
3989	print errormsg
3990	print ' ' + fname + ": all model-experiment file '" + filen + \
3991	"' does not exist !!"
3992	quit(-1)
3993
3994	objf = open(filen, 'r')
3995	for line in objf:
3996	if line[0:1] != '#' and len(line) > 1:
3997	values = line.replace('\n','').split(' ')
3998	if values[0] == 'allmodexps:':
3999	allplts = gen.stringList_dictKeysVals(values[1],cV=':')
4000	elif values[0] == 'Nallmodexp:':
4001	Nallmodexp = int(values[1])
4002
4003	objf.close()
4004
4005	return allplts, Nallmodexp
4006
4007	def draw_allmodexp_plots(config, allvarcomp, plots, modexp, odir, figscr, debug):
4008	""" Function to draw all model-experiment plots
4009	config= Configuration of the experiment
4010	allvarcomp = dictionary with the file headers for each variable
4011	plots= dictionary with the plots
4012	modexps= list with the names of the model-experiment pairs
4013	odir= output experiment folder
4014	figscr= whether figures should be done from the scratch or not
4015
4016	* Plot as
4017	{[kplot]} = [varn1]\|[op1]#[varn2]\|[op2]#[...[varnN]\|[opN]], ...
4018	[kplot] ___
4019	Nlines: Multiple lines with one line by model-experiment pairs
4020	[varn]
4021	variable
4022	[op]
4023	'+' separated list of operations
4024	in figures with more than 1 variable, use '#' to separate the [varn]\|[op]
4025	"""
4026	fname = 'draw_allmodexp_plots'
4027
4028	os.chdir(odir)
4029
4030	# Folder with the processed files
4031	infold = config['ofold']
4032
4033	# Dictionary with the operations with surnames for the operated variable
4034	opersurnames = {}
4035	opsur = gen.get_specdictionary_HMT(config, H='opsur_',M='',T='')
4036	for opsr in opsur.keys():
4037	opn = opsr.split('_')[1]
4038	vls = opsur[opsr].split(':')
4039	opersurnames[opn] = vls
4040
4041	# time values
4042	# Units time for the plots
4043	rd = config['CFreftime']
4044	tunits = config['CFunitstime'] + '!since!' + rd[0:4] + '-' + rd[4:6] + '-' + \
4045	rd[6:8] + '!' + rd[8:10] + ':' + rd[10:12] + ':' + rd[12:14]
4046	# time ticks kind
4047	tkind = config['timekind']
4048	# time ticks format
4049	tfmt = config['timefmt']
4050	# time axis label
4051	tlab = config['timelabel']
4052	timevals = [tunits, tkind, tfmt, tlab]
4053
4054	# Dictionary of plot specificities
4055	specplotkeyn = 'specificvarplot'
4056	plotspecifics = {}
4057	if config.has_key(specplotkeyn):
4058	# [minval]: minimum value
4059	# [maxval]: minimum value
4060	# [colorbar]: name of the colorbar (from matplotlib) to use
4061	# [cntformat]: format of the contour labels
4062	# [colorcnt]: color for the countor lines
4063	plotspecs = config[specplotkeyn].split(':')
4064	for pltspc in plotspecs:
4065	pltvls = pltspc.split('\|')
4066	vn = pltvls[0]
4067	op = pltvls[1]
4068	fn = pltvls[2]
4069	plotspecifics[fn + '_' + vn + '_' + op] = pltvls[3:]
4070	if debug:
4071	print 'Specific values for plots _______'
4072	gen.printing_dictionary(plotspecifics)
4073
4074	# Graphical labels and configuration
4075	modgc, expgc = graphical_characteristics(config, debug)
4076
4077	# Kind of figures
4078	kindfigure = config['kindfig']
4079
4080	# Map value
4081	mapvalue = config['mapval']
4082
4083	# pythone scripts HOME
4084	pyHOME = config['pyHOME']
4085
4086	# Title-text of operations
4087	opexplained = {}
4088	optits = config['titleoperations'].split(':')
4089	for optit in optits:
4090	opn = optit.split('\|')[0]
4091	opt = optit.split('\|')[1]
4092	opexplained[opn] = opt
4093	if debug:
4094	print 'Titles for operations _______'
4095	gen.printing_dictionary(opexplained)
4096
4097	for kplot in plots.keys():
4098	varsplt = plots[kplot]
4099	for varplt in varsplt:
4100	if debug:
4101	print " printing '" + kplot + "' var ':" + varplt + "'..."
4102	varops = varplt.split('#')
4103
4104	# CF variables in plot
4105	CFvarsplot = []
4106	# Files in plot
4107	filesplot = []
4108	# Variables in plot within the files
4109	varsplot = []
4110	# Dims in figure
4111	dimsplot = []
4112	# pictoric values in figure
4113	pictplot = []
4114	# Name of the figure
4115	figname = ''
4116	# Title of the figure
4117	titfigure = ''
4118
4119	ivp = 0
4120	for varop in varops:
4121	vn = varop.split('\|')[0]
4122	op = varop.split('\|')[1]
4123
4124	# CF variables in plot
4125	CFvarsplot.append(vn)
4126
4127	if op.find('pinterp') != -1: gP = 'p'
4128	else: gP = ''
4129
4130	vnopS = vn + '_' + op
4131	if not allvarcomp.has_key(vnopS):
4132	print errmsg
4133	print ' ' + fname + ": no file for variable-operation '" + \
4134	vnopS + "' !!"
4135	modexpvals = allvarcomp[vnopS]
4136	# modexpvalsS = allvarcomp[vnopS]
4137	# modexpvals = modexpvalsS[0].split('\|')
4138	Nmodexps = len(modexpvals)
4139
4140	for imodexp in range(0,Nmodexps,2):
4141	modexp = modexpvals[imodexp]
4142	headf = modexpvals[imodexp+1]
4143	if debug:
4144	print ' ' + fname + " model/experiment: '" + modexp + "'"
4145
4146	filen = infold + '/' + modexp + '/' + vn +'_'+ headf+gP + '_' + \
4147	op.replace('+','_') + '.nc'
4148
4149	filesplot.append(filen)
4150	# Do we have processed the given variable?
4151	if not os.path.isfile(filen):
4152	print warnmsg
4153	print " " + fname + ": there is no file for variable '" + \
4154	varop + "' skiping it !!"
4155	break
4156
4157	if imodexp == 0:
4158	# Name of the variable inside the file
4159	vnsur = varnoper(vn, op, opersurnames)
4160	varsplot.append(vnsur)
4161
4162	# Dimensions in file
4163	try:
4164	with gen.Capturing() as output:
4165	dims = ncvar.idims(filen)
4166	except:
4167	print errmsg
4168	print 'ncvar.idims('+filen+')'
4169	for s1out in output: print s1out
4170	quit(-1)
4171
4172	if imodexp == 0:
4173	# Dimensions in plot
4174	dimsplot.append(dims)
4175
4176	# Header of the name of the figure
4177	if ivp == 0:
4178	figname = kplot+'_'+vn+'_'+headf+'_'+ op.replace('+','-')
4179	else:
4180	figname = figname+'-'+vn+'_'+headf+'_'+op.replace('+','-')
4181
4182	# pictoric values for the figure
4183	Sfivaop = kplot + '_' + vn + '_' + op
4184	if plotspecifics.has_key(Sfivaop):
4185	pictvals = plotspecifics[Sfivaop]
4186	else:
4187	Vvals = gen.variables_values(vn)
4188	pictvals = [Vvals[2],Vvals[3],Vvals[6], '%g', 'fixc', 'black']
4189
4190	pictplot.append(pictvals)
4191
4192	ivp = ivp + 1
4193	# End of variable-operation
4194	figname = figname + '.' + kindfigure
4195
4196	# Title of figure
4197	titfigure = 'all!model-experiments!' + vn + optitle(op,opexplained)
4198
4199	draw_plot(kplot, CFvarsplot, filesplot, varsplot, dimsplot, pictplot, \
4200	figname, titfigure, kindfigure, mapvalue, timevals, expgc, odir, \
4201	pyHOME, figscr, debug)
4202
4203	# End of variables-operations
4204
4205	# End of kind of plots
4206
4207	return
4208
4209	# Files with information about the configuration of the script
4210	inffiles = ['varcompute.inf', 'all_computevars.inf', 'all_statsvars.inf']
4211
4212	####### #######
4213	## MAIN
4214	#######
4215	if not os.path.isfile('model_graphics.dat'):
4216	print errmsg
4217	print main + ": configuration file 'model_graphics.dat' does not exist!!"
4218	quit(-1)
4219
4220	# Getting configuration from external ASCII file 'model_graphics.dat'
4221	cnf = gen.get_configuration('model_graphics.dat', False)
4222
4223	verify_configuration(cnf, gen.Str_Bool(cnf['debug']))
4224
4225	# scratches
4226	scratch, filescratch, figscratch, diffscratch, figdiffscratch, moddiffscratch, \
4227	figmoddiffscratch, figallmodexpscratch, addfiles, addfigures, adddiffs, \
4228	adddifffigures, addmoddiffs, addmoddifffigures, addallmodexpfigures, dbg = \
4229	scratches(cnf)
4230
4231	if dbg:
4232	print 'Experiment configuration scratches _______'
4233	print ' scratch:', scratch, 'filescratch:', filescratch,'figscratch:', figscratch
4234	print ' diffscratch:', diffscratch, 'figdiffscratch:', figdiffscratch
4235	print ' moddiffscratch:', moddiffscratch, 'figmoddiffscratch:', figmoddiffscratch
4236	print 'Experiment configuration adding _______'
4237	print ' addfiles:', addfiles, 'addfigures:', addfigures
4238	print ' adddiffs:', adddiffs, 'adddifffigures:', adddifffigures
4239	print ' addmoddiffs:', addmoddiffs, 'addmoddifffigures:', addmoddifffigures
4240
4241	# Getting models
4242	mods = cnf['models'].split(':')
4243
4244	# Getting graphical characeristics
4245	modGC, expGC = graphical_characteristics(cnf, dbg)
4246
4247	# Models loop
4248	##
4249
4250	# dictionary with the experiments of each model
4251	modexps = {}
4252
4253	for mod in mods:
4254	print mod
4255	if scratch:
4256	if cnf['ofold'] != cnf['ifold']:
4257	sout = sub.call('rm -rf '+cnf['ofold']+'/'+mod+' > /dev/null', shell=True)
4258	else:
4259	print warnmsg
4260	print ' ' + main + ": input '" + cnf['ifold'] + "' and output '" + \
4261	cnf['ofold'] + "' are the same folder !!"
4262	print " Keeping output folder although it has 'scratch' start"
4263
4264	# Get experiments and headers of model
4265	exps, fheaders = exp_headers(mod,cnf)
4266
4267	modexps[mod] = exps
4268	# Characteristics of the model
4269	Modinf = ncvar.model_characteristics(mod,'None','False')
4270	dnx = Modinf.dimxn
4271	dny = Modinf.dimyn
4272	dnz = Modinf.dimzn
4273	dnt = Modinf.dimtn
4274	dns = Modinf.dimsn
4275	vdnx = Modinf.vardxn
4276	vdny = Modinf.vardyn
4277	vdnz = Modinf.vardzn
4278	vdnt = Modinf.vardtn
4279	vdns = Modinf.vardsn
4280
4281	modgraphv = modGC[mod]
4282	Modinf = ModelInf(mod, mod, dnx, dny, dnz, dnt, dns, vdnx, vdny, vdnz, vdnt, \
4283	vdns, modgraphv.tmod)
4284
4285	if dbg:
4286	print ' model characteristics _______'
4287	gen.printing_class(Modinf)
4288
4289	moddims = dnx + ',' + dny + ',' + dnz + ',' + dnt
4290	modvdims = vdnx + ',' + vdny + ',' + vdnz + ',' + vdnt
4291
4292	# Experiments loop
4293	##
4294	for exp in exps:
4295	print ' ' + exp + '...'
4296
4297	# input folder
4298	iwdir = cnf['ifold'] + '/' + mod + '/' + exp
4299
4300	# Does input folder exist?
4301	if not os.path.isdir(iwdir):
4302	if filescratch*addfiles:
4303	print errmsg
4304	print " " + main + ": folder '" + iwdir + "' does not exist !!"
4305	quit(-1)
4306	else:
4307	print warnmsg
4308	print " " + main + ": folder '" + iwdir + "' does not exist !!"
4309	print " but start of files from scratch 'filescratch':", \
4310	filescratch, "and there is not the add of any file 'addfiles':", \
4311	addfiles
4312	print ' so, keep moving'
4313
4314	owdir = cnf['ofold'] + '/' + mod + '/' + exp
4315	sout = sub.call('mkdir -p ' + owdir, shell=True)
4316
4317	# Need to pass to analyze all the data?
4318	if filescratch:
4319	for inff in inffiles:
4320	if dbg: print " removing information file '" + inff + "' ..."
4321	ins = 'rm ' + owdir + '/' + inff + ' >& /dev/null'
4322	sout = sub.call(ins, shell=True)
4323
4324	objf = open(owdir+'/all_computevars.inf','w')
4325	objf.write("## Computation of variables \n")
4326	objf.close()
4327	objf = open(owdir+'/all_statsvars.inf','w')
4328	objf.write("## Computation of statistics \n")
4329	objf.close()
4330
4331	varcompf = owdir + '/varcompute.inf'
4332	if addfiles:
4333	if dbg:
4334	print ' ' + main + ": adding variables to copmute removing " + \
4335	"file '" + varcompf + "' ..."
4336	sub.call('rm ' + varcompf +' >& /dev/null', shell=True)
4337
4338	if not os.path.isfile(varcompf):
4339	# Does input folder has header files?
4340	ih=1
4341	# Dictionary with the list of files for each headers
4342	files = {}
4343	# Dictionary with a file fromor each headers
4344	testfiles = {}
4345
4346	for fh in fheaders:
4347	if filescratch:
4348	ins = 'rm '+ owdir+'/_' + fh + '.nc >& /dev/null'
4349	sout = sub.call(ins, shell=True)
4350	files1h = gen.files_folder(iwdir,fh)
4351	if len(files1h) < 1:
4352	print errmsg
4353	print ' ' + main + ": folder '" + iwdir + "' does not " + \
4354	"contain files '" + fh + "*' !!"
4355	quit(-1)
4356	files[fh] = files1h
4357	testfiles[fh] = files1h[0]
4358
4359	if dbg:
4360	print ' Dictionary of files _______'
4361	gen.printing_dictionary(files)
4362
4363	print ' Dictionary of test-files _______'
4364	gen.printing_dictionary(testfiles)
4365
4366	allcompvar, Nvar = compvars_listconstruct(cnf, Modinf, files, testfiles, \
4367	iwdir, owdir, dbg)
4368	else:
4369	print warnmsg
4370	print ' ' + main + ": getting variables to compute already from file !!"
4371	files, testfiles, allcompvar, Nvar = read_varcomp_file(varcompf)
4372
4373	# End of avoiding to repeat all the experiment search
4374
4375	print " For experiment '" + exp + "' is required to compute:", Nvar, \
4376	"variables"
4377
4378	if dbg:
4379	print 'Variables to compute _______'
4380	gen.printing_dictionary(allcompvar)
4381
4382	# Computing variables
4383	print " Computing variables ..."
4384	compute_vars(cnf, Modinf, iwdir, owdir, files, allcompvar, filescratch, dbg)
4385
4386	# Figures of variables direct from files
4387	print " Plotting direct figures ..."
4388	dirfigf = owdir + '/directplotsdraw.inf'
4389	if figscratch:
4390	sout = sub.call('rm ' + dirfigf + ' >& /dev/null', shell=True)
4391
4392	objf = open(owdir+'/all_figures.inf','w')
4393	objf.write("## Drawing of all figures \n")
4394	objf.close()
4395
4396	if addfigures:
4397	if dbg:
4398	print ' ' + main + ": adding direct figures removing file '" + \
4399	dirfigf + "' ..."
4400	sout = sub.call('rm ' + dirfigf + ' >& /dev/null', shell=True)
4401
4402	if not os.path.isfile(dirfigf):
4403	listplots, Nplt = plots_listconstruct(cnf, 'DIRPLT', dirfigf, owdir, dbg)
4404	else:
4405	print warnmsg
4406	print ' ' + main + ": getting plots to draw already from file !!"
4407	listplots, Nplt = read_plot_file(dirfigf)
4408
4409	# End of avoiding to repeat all the plots search
4410
4411	print " For experiment '" + exp + "' is required to plot:", Nplt, "plots"
4412
4413	if dbg:
4414	print 'Plots to draw _______'
4415	gen.printing_dictionary(listplots)
4416
4417	draw_plots(cnf, listplots, mod, exp, owdir, allcompvar, figscratch, dbg)
4418
4419	# end of experiments loop
4420
4421	### ## #
4422	# Experiment differences
4423	## # ## #
4424
4425	# Experiments loop
4426	##
4427	print ' ***** ** * * ** * * ** * ** *****'
4428	print " ** '" + mod + "': Inter experiments differences "
4429	print ' ***** ** * * ** * * ** * ** *****'
4430	print ' experiments:', exps
4431	# There are two kind of differences:
4432	# DIFFop: differences between operations of each given variable.
4433	# Differences are computed directly from the last stage of the operation
4434	# DIFFvar: operations of differences between variables
4435	# First are computed the differences from the initial variable file
4436	# and then operations are made
4437	# NOTE: remember that: meanvar2 - meanvar1 = mean(var2 - var1)
4438	difffiles = ['diffop.inf', 'diffvar.inf']
4439
4440	Nexps = len(exps)
4441	for iexp1 in range(0,Nexps-1):
4442	exp1 = exps[iexp1]
4443	for iexp2 in range(iexp1+1,Nexps):
4444	exp2 = exps[iexp2]
4445	Sexps = exp2 + '-' + exp1
4446	print ' ' + Sexps + '...'
4447	owdir = cnf['ofold'] + '/' + mod + '/' + exp2 + '-' + exp1
4448	sout = sub.call('mkdir -p ' + owdir, shell=True)
4449
4450	# Removing files with list of differences if should be started from scratch or add differences
4451	for fdiff in difffiles:
4452	difff = owdir + '/' + fdiff
4453	if diffscratch:
4454	sub.call('rm -rf' + owdir +' >& /dev/null', shell=True)
4455	sub.call('rm ' + difff +' >& /dev/null', shell=True)
4456	objf = open(owdir+'/all_' + fdiff,'w')
4457	objf.write("## Computation and drawing of differences " + \
4458	"between '" + exp2 + "'-'" + exp1 + "'\n")
4459	objf.close()
4460	difff = owdir + '/all_vardiffstatistics.inf'
4461	sub.call('rm ' + difff +' >& /dev/null', shell=True)
4462	objf = open(owdir+'/all_' + fdiff,'w')
4463	objf.write("## Computation of all differences statistics " + \
4464	"between '" + exp2 + "'-'" + exp1 + "'\n")
4465	objf.close()
4466
4467	if adddiffs:
4468	if dbg:
4469	print ' ' + main + ": adding differences removing file '" + \
4470	difff + "' ..."
4471	sub.call('rm ' + difff +' >& /dev/null', shell=True)
4472
4473	# op and var differences
4474	diffvarcompf = owdir + '/' + difffiles[0]
4475	if not os.path.isfile(diffvarcompf):
4476	alldiffop, alldiffvar,Nopdiffs,Nvardiffs=diffvarop_listconstruct(cnf,\
4477	[mod, mod], [exp1, exp2], [False, False], owdir, dbg)
4478	else:
4479	print warnmsg
4480	print ' ' + main + ": getting 'op' and 'var'differences to " + \
4481	"compute already from file !!"
4482	alldiffop, Nopdiffops = read_diff_file(diffvarcompf)
4483	alldiffvar, Nopdiffvars = read_diff_file(owdir+'/'+difffiles[1])
4484
4485	# End of avoiding to repeat all the experiment search
4486
4487	print " For experiments '"+exp2+"'-'"+exp1+"' is required to compute:", \
4488	Nvar, "differences"
4489
4490	if dbg:
4491	print 'Differences to compute _______'
4492	gen.printing_dictionary(alldiffop)
4493	gen.printing_dictionary(alldiffvar)
4494
4495	# Computing differences
4496	##
4497	print " Computing operation differences ..."
4498	compute_op_diffs(cnf, alldiffop, owdir, diffscratch, dbg)
4499	print " Computing variable differences ..."
4500	compute_var_diffs(cnf, alldiffvar, owdir, diffscratch, dbg)
4501
4502	# Plotting operation differences
4503	##
4504	print " " + main + ": Plotting operation differences' figures ..."
4505	dirfigf = owdir + '/diffopplotsdraw.inf'
4506	if figscratch:
4507	sout = sub.call('rm ' + dirfigf + ' >& /dev/null', shell=True)
4508
4509	objf = open(owdir+'/all_diffopfigures.inf','w')
4510	objf.write("## Drawing of all operation difference figures \n")
4511	objf.close()
4512
4513	if adddifffigures:
4514	if dbg:
4515	print ' ' + main + ": adding differences' figures removing file '" +\
4516	dirfigf + "' ..."
4517
4518	sout = sub.call('rm ' + dirfigf + ' >& /dev/null', shell=True)
4519
4520	if not os.path.isfile(dirfigf):
4521	listplots, Nplt = plots_listconstruct(cnf, 'PLOTDIFFOP', dirfigf, owdir, \
4522	dbg)
4523	else:
4524	print warnmsg
4525	print ' ' + main + ": getting plots to draw already from file !!"
4526	listplots, Nplt = read_plot_file(dirfigf)
4527
4528	# End of avoiding to repeat all the plots search
4529
4530	print " For experiment 'operation' differences '" + Sexps + "' is " + \
4531	"required to plot:" , Nplt, "plots"
4532
4533	if dbg:
4534	print ' Plots to draw _______'
4535	gen.printing_dictionary(listplots)
4536
4537	draw_diff_plots(cnf, listplots, owdir, alldiffop, 'diffop', figdiffscratch, \
4538	dbg)
4539
4540	# Plotting variable differences
4541	##
4542	print " " + main + ": Plotting variable differences' figures ..."
4543	dirfigf = owdir + '/diffvarplotsdraw.inf'
4544	if figscratch:
4545	sout = sub.call('rm ' + dirfigf + ' >& /dev/null', shell=True)
4546
4547	objf = open(owdir+'/all_diffvarfigures.inf','w')
4548	objf.write("## Drawing of all variables difference figures \n")
4549	objf.close()
4550
4551	if adddifffigures:
4552	if dbg:
4553	print ' ' + main + ": adding differences' figures removing file '" +\
4554	dirfigf + "' ..."
4555	sout = sub.call('rm ' + dirfigf + ' >& /dev/null', shell=True)
4556
4557	if not os.path.isfile(dirfigf):
4558	listplots, Nplt = plots_listconstruct(cnf, 'PLOTDIFFVAR', dirfigf, owdir,\
4559	dbg)
4560	else:
4561	print warnmsg
4562	print ' ' + main + ": getting plots to draw already from file !!"
4563	listplots, Nplt = read_plot_file(dirfigf)
4564
4565	# End of avoiding to repeat all the plots search
4566
4567	print " For experiment 'variables' differences '" + Sexps + "' is " + \
4568	"required to plot:" , Nplt, "plots"
4569
4570	if dbg:
4571	print ' Plots to draw _______'
4572	gen.printing_dictionary(listplots)
4573
4574	draw_diff_plots(cnf, listplots, owdir, alldiffvar, 'diffvar', figdiffscratch,\
4575	dbg)
4576
4577	# end of mods loop
4578
4579	### ## #
4580	# Model differences
4581	### ## #
4582
4583	# Models already at the reference projection
4584	ModComProj = cnf['RefProj'].split(':')
4585	if dbg:
4586	print 'Models to use as reference for the projection:', ModComProj
4587
4588	# CDO operators to use for the remapping of each variable
4589	CDOremapvar0 = gen.get_specdictionary_HMT(cnf, H='CDOremapvar_')
4590	CDOremapvar = {}
4591	for Skey in CDOremapvar0.keys():
4592	CDOremapvar[Skey] = gen.str_list(CDOremapvar0[Skey],':')
4593
4594	if dbg:
4595	print 'Remaping operators to use _______'
4596	gen.printing_dictionary(CDOremapvar)
4597
4598	# Getting common base wights for grid point remapping with CDO
4599	frefgridn = 'RefProj.nc'
4600	for mod in mods:
4601	frefgriddes = cnf['ofold'] + '/' + mod + '/' + frefgridn
4602	if moddiffscratch:
4603	sout = sub.call('rm ' + frefgriddes + ' >& /dev/null', shell=True)
4604	if gen.searchInlist(ModComProj,mod) and not os.path.isfile(frefgriddes):
4605	# Looking for a file to use (all files are CF compilant, but we need to
4606	# be sure that at least they have 'lon', 'lat' variables)
4607	exps = modexps[mod]
4608	exp = exps[0]
4609	cdofdir = cnf['ofold'] + '/' + mod + '/' + exp
4610	cdofile = None
4611	for cdoop in CDOremapvar.keys():
4612	cdovars = CDOremapvar[cdoop]
4613	for cdovar in cdovars:
4614	cdofs= gen.files_folder_HMT(folder=cdofdir,head=cdovar+'_',tail='.nc')
4615	for cdof in cdofs:
4616	# Getting direct files as [var]_[fhead].nc
4617	if len(cdof.split('_')) == 2:
4618	cdofile = cdofdir + '/' + cdof
4619	cdovn = cdovar
4620	break
4621	if cdofile is None:
4622	print errmsg
4623	print ' ' + main + ": no proper file to get projection information " + \
4624	"from '" + cdofdir + "' has been found !!"
4625	quit(-1)
4626
4627	print " Creation of reference projection information file from '" + mod+ "'"
4628	if dbg:
4629	print " using file: '" + cdofile
4630
4631	# Getting only 1 time-step for disk space issues
4632	try:
4633	with gen.Capturing() as output:
4634	ncvar.DataSetSection('time,0,1,1', cdofile)
4635	except:
4636	print errmsg
4637	print 'ncvar.DataSetSection(time,0,1,1, ' + cdofile + ')'
4638	for s1out in output: print s1out
4639	quit(-1)
4640
4641	of = cdofile.split('.')[0] + '_time_B0-E1-I1.nc'
4642
4643	# Cleaning variables
4644	try:
4645	with gen.Capturing() as output:
4646	ncvar.selvar('lon@lon,lat@lat,time@time', of, cdovn)
4647	except:
4648	print errmsg
4649	print 'ncvar.selvar(lon@lon,lat@lat, ' + of + ', ' + cdovn + ')'
4650	for s1out in output: print s1out
4651	quit(-1)
4652	sout = sub.call('mv selvar_new.nc ' + frefgriddes, shell=True)
4653	sout = sub.call('rm ' + of + ' >& /dev/null', shell=True)
4654
4655	cnf['RefProjfile'] = frefgriddes
4656
4657	if dbg:
4658	for s1out in output: print s1out
4659	break
4660
4661	elif gen.searchInlist(ModComProj,mod) and os.path.isfile(frefgriddes):
4662	# Including reprojecting reference file in configuration dictionary
4663	cnf['RefProjfile'] = frefgriddes
4664
4665	Nmods = len(mods)
4666	for imod1 in range(0,Nmods-1):
4667	mod1 = mods[imod1]
4668	exps1 = modexps[mod1]
4669
4670	for imod2 in range(imod1+1,Nmods):
4671	mod2 = mods[imod2]
4672	exps2 = modexps[mod2]
4673
4674	Smods = mod2 + '-' + mod1
4675	print ' ***** ** * * ** * * ** * ** *****'
4676	print " ** '" + Smods + "': Inter models differences "
4677	print ' ***** ** * * ** * * ** * ** *****'
4678	print ' experiments mod1:', exps1
4679	print ' experiments mod2:', exps2
4680
4681	# Experiments loop
4682	##
4683	difffiles = ['diffop.inf', 'diffvar.inf']
4684
4685	Nexps1 = len(exps1)
4686	Nexps2 = len(exps2)
4687	for iexp1 in range(0,Nexps1):
4688	exp1 = exps1[iexp1]
4689	for iexp2 in range(0,Nexps2):
4690	exp2 = exps2[iexp2]
4691	Sexps = exp2 + '-' + exp1
4692	print ' ' + Sexps + '...'
4693	owdir = cnf['ofold'] + '/' + Smods + '/' + Sexps
4694	sout = sub.call('mkdir -p ' + owdir, shell=True)
4695
4696	# Getting the right projection in order to perform differences
4697	difmods = [mod1, mod2]
4698	difexps = [exp1, exp2]
4699	difreproj = [False, False]
4700
4701	for imod in range(2):
4702	mod = difmods[imod]
4703	exp = difexps[imod]
4704	if not gen.searchInlist(ModComProj,mod):
4705	print " Projecting files of model '" + mod + "' exp: '" + \
4706	exp + "' with file '" + cnf['RefProjfile'] + "'"
4707	reproject_modexp(mod, exp, cnf, CDOremapvar, moddiffscratch, \
4708	addmoddiffs, dbg)
4709	difreproj[imod] = True
4710
4711	# Removing files with list of differences if should be started from
4712	# scratch or add differences
4713	for fdiff in difffiles:
4714	difff = owdir + '/' + fdiff
4715	if moddiffscratch:
4716	sub.call('rm -rf' + owdir +' >& /dev/null', shell=True)
4717	sub.call('rm ' + difff +' >& /dev/null', shell=True)
4718	objf = open(owdir+'/all_' + fdiff,'w')
4719	objf.write("## Computation and drawing of differences " + \
4720	"between '" + mod2+'@'+exp2 + "'-'" + mod1+'@'+exp1 + "'\n")
4721	objf.close()
4722	difff = owdir + '/all_vardiffstatistics.inf'
4723	sub.call('rm ' + difff +' >& /dev/null', shell=True)
4724	objf = open(owdir+'/all_' + fdiff,'w')
4725	objf.write("## Computation of all differences statistics " + \
4726	"between '" + mod2+'@'+exp2 + "'-'" + mod1+'@'+exp1 + "'\n")
4727	objf.close()
4728
4729	if addmoddiffs:
4730	if dbg:
4731	print ' ' + main + ": adding model differences " + \
4732	"removing file '" + difff + "' ..."
4733	sub.call('rm ' + difff +' >& /dev/null', shell=True)
4734
4735	# op and var differences
4736	diffvarcompf = owdir + '/' + difffiles[0]
4737	if not os.path.isfile(owdir+'/'+difffiles[0]) or \
4738	not os.path.isfile(owdir+'/'+difffiles[1]):
4739	alldiffop, alldiffvar, Nopdiffs, Nvardiffs = \
4740	diffvarop_listconstruct(cnf, difmods, difexps, difreproj,owdir,\
4741	dbg)
4742	else:
4743	print warnmsg
4744	print ' ' + main + ": getting 'op' and 'var'differences to " + \
4745	"compute already from file !!"
4746	alldiffop, Nopdiffops = read_diff_file(diffvarcompf)
4747	alldiffvar, Nopdiffvars = read_diff_file(owdir+'/'+difffiles[1])
4748
4749	# End of avoiding to repeat all the experiment search
4750
4751	print " For experiments '" + mod2+'@'+exp2 + "'-'" + mod+'@'+exp1+\
4752	"' is required to compute:", Nvar, "differences"
4753
4754	if dbg:
4755	print ' op differences to compute _______'
4756	gen.printing_dictionary(alldiffop)
4757	print ' var differences to compute _______'
4758	gen.printing_dictionary(alldiffvar)
4759
4760	# Computing differences
4761	##
4762	print " Computing operation differences ..."
4763	compute_op_diffs(cnf, alldiffop, owdir, moddiffscratch, dbg)
4764	print " Computing variable differences ..."
4765	compute_var_diffs(cnf, alldiffvar, owdir, moddiffscratch, dbg)
4766
4767	# Plotting operation differences
4768	##
4769	print " " + main + ": Plotting operation differences' figures ..."
4770	dirfigf = owdir + '/diffopplotsdraw.inf'
4771	if figmoddiffscratch:
4772	sout = sub.call('rm ' + dirfigf + ' >& /dev/null', shell=True)
4773
4774	objf = open(owdir+'/all_diffopfigures.inf','w')
4775	objf.write("## Drawing of all operation difference figures \n")
4776	objf.close()
4777
4778	if addmoddifffigures:
4779	if dbg:
4780	print ' ' + main + ": adding model differences' figures " + \
4781	"removing file '" + dirfigf + "' ..."
4782	sout = sub.call('rm ' + dirfigf + ' >& /dev/null', shell=True)
4783
4784	if not os.path.isfile(dirfigf):
4785	listplots, Nplt = plots_listconstruct(cnf, 'PLOTDIFFOP', dirfigf,\
4786	owdir, dbg)
4787	else:
4788	print warnmsg
4789	print ' ' + main + ": getting plots to draw already from file !!"
4790	listplots, Nplt = read_plot_file(dirfigf)
4791
4792	# End of avoiding to repeat all the plots search
4793
4794	print " For experiment 'operation' differences '" + Sexps + "' is "+\
4795	"required to plot:" , Nplt, "plots"
4796
4797	if dbg:
4798	print 'Plots to draw _______'
4799	gen.printing_dictionary(listplots)
4800
4801	draw_diff_plots(cnf, listplots, owdir, alldiffop, 'diffop', \
4802	figmoddiffscratch, dbg)
4803
4804	# Plotting variable differences
4805	##
4806	print " " + main + ": Plotting variable differences' figures ..."
4807	dirfigf = owdir + '/diffvarplotsdraw.inf'
4808	if figscratch:
4809	sout = sub.call('rm ' + dirfigf + ' >& /dev/null', shell=True)
4810
4811	objf = open(owdir+'/all_diffvarfigures.inf','w')
4812	objf.write("## Drawing of all variables difference figures \n")
4813	objf.close()
4814
4815	if adddifffigures:
4816	if dbg:
4817	print ' ' + main + ": adding differences's figures " + \
4818	"removing file '" + dirfigf + "' ..."
4819	sout = sub.call('rm ' + dirfigf + ' >& /dev/null', shell=True)
4820
4821	if not os.path.isfile(dirfigf):
4822	listplots, Nplt = plots_listconstruct(cnf,'PLOTDIFFVAR',dirfigf, \
4823	owdir, dbg)
4824	else:
4825	print warnmsg
4826	print ' ' + main + ": getting plots to draw already from file !!"
4827	listplots, Nplt = read_plot_file(dirfigf)
4828
4829	# End of avoiding to repeat all the plots search
4830
4831	print " For experiment 'variables' differences '" + Sexps+ "' is "+ \
4832	"required to plot:" , Nplt, "plots"
4833
4834	if dbg:
4835	print ' Plots to draw _______'
4836	gen.printing_dictionary(listplots)
4837
4838	draw_diff_plots(cnf, listplots, owdir, alldiffvar, 'diffvar', \
4839	figdiffscratch, dbg)
4840
4841	# end of exp2 loop
4842	# end of exp1 loop
4843	# end of mod2 loop
4844	# end of mod1 loop
4845
4846	### ## #
4847	# All model/exps linear plots
4848	### ## #
4849
4850	allmodexps = []
4851	for mod in mods:
4852	for exp in modexps[mod]:
4853	allmodexps.append(mod + '/' + exp)
4854
4855	print ' ***** ** * * ** * * ** * ** *****'
4856	print ' ** Multi models, multi experiments plots '
4857	print ' ***** ** * * ** * * ** * ** *****'
4858	print ' Model/experiments:', allmodexps
4859
4860	owdir = cnf['ofold'] + '/allmodexps'
4861	sout = sub.call('mkdir -p ' + owdir, shell=True)
4862
4863	# Plotting all model-experiments lines
4864	##
4865	print " " + main + ": Plotting all model-experiments variable figures ..."
4866	allf = owdir + '/allmodexp.inf'
4867	dirfigf = owdir + '/allmodexpfigures.inf'
4868	allfigf = owdir+'/all_figures.inf'
4869	if figallmodexpscratch:
4870	sout = sub.call('rm ' + allf + ' >& /dev/null', shell=True)
4871	sout = sub.call('rm ' + dirfigf + ' >& /dev/null', shell=True)
4872	sout = sub.call('rm ' + allfigf + ' >& /dev/null', shell=True)
4873
4874	objf = open(allfigf,'w')
4875	objf.write("## Drawing of all variables figures for all model-experiments\n")
4876	objf.close()
4877
4878	if addallmodexpfigures:
4879	if dbg:
4880	print ' ' + main + ": adding model-experiment differences removing " + \
4881	" file '" + allfigf + "' ..."
4882	print ' ' + main + ": adding model-experiment differences' figures " + \
4883	"removing file '" + dirfigf + "' ..."
4884	sout = sub.call('rm ' + allfigf + ' >& /dev/null', shell=True)
4885	sout = sub.call('rm ' + dirfigf + ' >& /dev/null', shell=True)
4886
4887	if not os.path.isfile(dirfigf):
4888	listplots, Nplt = plots_listconstruct(cnf, 'PLOTALLMODEXP', dirfigf, owdir, dbg)
4889	else:
4890	print warnmsg
4891	print ' ' + main + ": getting plots to draw already from file !!"
4892	listplots, Nplt = read_plot_file(dirfigf)
4893
4894	# End of avoiding to repeat all the plots search
4895
4896	print " Is required to plot:", Nplt, " all model-experiment plots"
4897
4898	if dbg:
4899	print ' Plots to draw _______'
4900	gen.printing_dictionary(listplots)
4901
4902	if not os.path.isfile(allf):
4903	# Constructing the information of modexp / variable
4904	allmodexpvar, Nallmodexpvar = allmodexps_listconstruct(cnf, allmodexps, owdir, dbg)
4905	else:
4906	allmodexpvar, Nallmodexpvar = allmodexps_read(allf)
4907
4908	draw_allmodexp_plots(cnf, allmodexpvar, listplots, modexps, owdir, figallmodexpscratch, dbg)
4909
4910	print main + ': all files and figures have been properly done !!'

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

Download in other formats: