Context Navigation

ppclass.py @ 923

Last change on this file since 923 was 923, checked in by aslmd, 12 years ago

UTIL PYTHON planetoplot_v2

PPCLASS AND PPPLOT

added global plot settings in pp() objects
self is returned for each method to allow e.g. one-line get + attribution
added labeling of 1D plots
fine-tuning of plot appearance (e.g. vector key)
corrected a problem with plotin (same object can be referred to all along)
fixed meanarea for 1D requests

PP.PY

no more -f for pp.py (files are simply given as arguments!)
added missing options to pp.py
nice informative header

PP_RELOAD.PY

pp_reload.py can now change colorbars and 1D markers

EXAMPLES

general update of examples + few additions
added a nice example tide with mixed 2D+1D plots

File size: 78.5 KB

Line
1	###############################################
2	## PLANETOPLOT ##
3	## --> PPCLASS ##
4	## A generic and versatile Python module ##
5	## ... to read netCDF files and plot ##
6	###############################################
7	## Author: Aymeric Spiga. 02-03/2013 ##
8	###############################################
9	# python built-in librairies
10	import os
11	import time as timelib
12	import pickle
13	# added librairies
14	import numpy as np
15	import netCDF4
16	import matplotlib.pyplot as mpl
17	# personal librairies
18	import ppplot
19	import ppcompute
20	###############################################
21
22	###################################
23	#### HEADER ##
24	#### ... executed when imported ##
25	###################################
26	# where settings files are located...
27	whereset = None
28	whereset = ppcompute.findset(whereset)
29	# ... we load user-defined automatic settings from set_ppclass.txt
30	zefile = "set_ppclass.txt"
31	glob_listx = [] ; glob_listy = [] ; glob_listz = [] ; glob_listt = []
32	glob_listarea = []
33	try:
34	f = open(whereset+zefile, 'r') ; lines = f.readlines()
35	for stuff in lines[5].strip().split(';'): glob_listx.append(stuff)
36	for stuff in lines[8].strip().split(';'): glob_listy.append(stuff)
37	for stuff in lines[11].strip().split(';'): glob_listz.append(stuff)
38	for stuff in lines[14].strip().split(';'): glob_listt.append(stuff)
39	for stuff in lines[17].strip().split(';'): glob_listarea.append(stuff)
40	except IOError:
41	print "warning: "+zefile+" not in "+whereset+" ; no presets."
42
43	##################################
44	#### USEFUL GENERIC FUNCTIONS ####
45	##################################
46
47	# inspect variables and dimensions in a netCDF file
48	def inspect(filename):
49	print "**** INSPECT FILE",filename
50	test = netCDF4.Dataset(filename)
51	print "**** VARIABLES: ",test.variables.keys()
52	for dim in test.dimensions.keys():
53	output = "**** DIMENSION: "+str(dim)+" "+str(len(test.dimensions[dim]))
54	try: output = output + " ----> "+str(test.variables[dim][0])+" "+str(test.variables[dim][-1])
55	except: pass
56	print output ; output = ""
57
58	# check a tab and exit if wrong. if just one string make it a list.
59	# (if allownumber, convert this into a string).
60	def checktab(tab,mess="",allownone=False,allownumber=False):
61	if tab is None:
62	if not allownone: print "pp.define: no "+mess ; exit()
63	else: pass
64	else:
65	if not isinstance(tab, list):
66	if isinstance(tab, str):
67	tab = [tab]
68	elif (isinstance(tab, int) or isinstance(tab, float)) and allownumber:
69	tab = [str(tab)]
70	else:
71	print "pp.define: "+mess+" should be either a string or a list of strings!" ; exit()
72	elif isinstance(tab, list):
73	if isinstance(tab[0],str):
74	pass
75	elif (isinstance(tab[0], int) or isinstance(tab[0], float)) and allownumber:
76	for iii in range(len(tab)): tab[iii] = str(tab[iii])
77	else:
78	print "pp.define: "+mess+" should be either a string or a list of strings!" ; exit()
79	return tab
80
81	# determine which method is to be applied to a given dimension
82	def findmethod(tab):
83	if tab is None: output = "free"
84	elif tab[0,0] != tab[0,1]: output = "comp"
85	else: output = "fixed"
86	return output
87
88	# read what is given by the user (version of T. Navarro)
89	def readslices(saxis):
90	if saxis == None:
91	zesaxis = None
92	else:
93	zesaxis = np.empty((len(saxis),2))
94	for i in range(len(saxis)):
95	a = separatenames(saxis[i])
96	if len(a) == 1:
97	zesaxis[i,:] = float(a[0])
98	else:
99	zesaxis[i,0] = float(a[0])
100	zesaxis[i,1] = float(a[1])
101	return zesaxis
102	# (needed by readslices)
103	# look for comas in the input name to separate different names (files, variables,etc ..)
104	def separatenames (name):
105	if name is None: names = None
106	else:
107	names = [] ; stop = 0 ; currentname = name
108	while stop == 0:
109	indexvir = currentname.find(',')
110	if indexvir == -1: stop = 1 ; name1 = currentname
111	else: name1 = currentname[0:indexvir]
112	names = np.concatenate((names,[name1]))
113	currentname = currentname[indexvir+1:len(currentname)]
114	return names
115
116	#######################
117	### THE MAIN OBJECT ###
118	#######################
119	class pp():
120
121	# print out a help string when help is invoked on the object
122	def __repr__(self):
123	whatprint = 'pp object. \"help(pp)\" for more information\n'
124	return whatprint
125
126	# default settings
127	# -- user can define settings by two methods.
128	# -- 1. yeah = pp(file="file.nc")
129	# -- 2. yeah = pp() ; yeah.file = "file.nc"
130	def __init__(self,file=None,var="notset",\
131	filegoal=None,vargoal=None,\
132	x=None,y=None,z=None,t=None,\
133	stridex=1,stridey=1,\
134	stridez=1,stridet=1,\
135	compute="mean",\
136	verbose=False,noproj=False,\
137	superpose=False,\
138	plotin=None,\
139	forcedimplot=-1,\
140	out="gui",\
141	filename="myplot",\
142	folder="./",\
143	xlabel=None,ylabel=None,\
144	xcoeff=None,ycoeff=None,\
145	proj=None,\
146	vmin=None,vmax=None,\
147	div=None,\
148	colorb=None,\
149	lstyle=None,\
150	marker=None,\
151	color=None,\
152	label=None,\
153	title=None):
154	self.request = None
155	self.nfin = 0 ; self.nvin = 0
156	self.nplotx = None ; self.nploty = None
157	self.nplotz = None ; self.nplott = None
158	self.status = "init"
159	self.fig = None ; self.subv = None ; self.subh = None
160	self.n = 0 ; self.howmanyplots = 0
161	self.nplot = 0
162	self.p = None
163	self.customplot = False
164	## what could be defined by the user
165	self.file = file
166	self.var = var
167	self.filegoal = filegoal
168	self.vargoal = vargoal
169	self.x = x ; self.y = y ## if None, free dimension
170	self.z = z ; self.t = t ## if None, free dimension
171	self.stridex = stridex ; self.stridey = stridey
172	self.stridez = stridez ; self.stridet = stridet
173	self.compute = compute
174	self.verbose = verbose
175	self.noproj = noproj
176	self.plotin = plotin
177	self.superpose = superpose
178	self.forcedimplot = forcedimplot
179	self.out = out
180	self.filename = filename
181	self.folder = folder
182	## here are user-defined plot settings
183	## -- if not None, valid on all plots in the pp() objects
184	self.xlabel = xlabel ; self.xcoeff = xcoeff
185	self.ylabel = ylabel ; self.ycoeff = ycoeff
186	self.proj = proj
187	self.vmin = vmin ; self.vmax = vmax
188	self.div = div
189	self.colorb = colorb
190	self.lstyle = lstyle
191	self.marker = marker
192	self.color = color
193	self.label = label
194	self.title = title
195
196	# print status
197	def printstatus(self):
198	if self.filename == "THIS_IS_A_CLONE":
199	pass
200	else:
201	print "**** PPCLASS. Done step: " + self.status
202
203	#####################################################
204	# EMULATE OPERATORS + - * / ** << FOR PP() OBJECTS #
205	#####################################################
206
207	# define the operation <<
208	# ... e.g. obj2 << obj1
209	# ... means: get init for pp object obj2 from another pp object obj1
210	# ... (this should solve the affectation trap obj2 = obj1)
211	def __lshift__(self,other):
212	if other.__class__.__name__ == "pp":
213	self.file = other.file
214	self.var = other.var
215	self.filegoal = other.filegoal
216	self.vargoal = other.vargoal
217	self.x = other.x ; self.y = other.y ## if None, free dimension
218	self.z = other.z ; self.t = other.t ## if None, free dimension
219	self.stridex = other.stridex ; self.stridey = other.stridey
220	self.stridez = other.stridez ; self.stridet = other.stridet
221	self.verbose = other.verbose
222	self.noproj = other.noproj
223	self.plotin = other.plotin
224	self.superpose = other.superpose
225	self.forcedimplot = other.forcedimplot
226	self.out = other.out
227	self.filename = other.filename
228	self.folder = other.folder
229	self.xlabel = other.xlabel ; self.xcoeff = other.xcoeff
230	self.ylabel = other.ylabel ; self.ycoeff = other.ycoeff
231	self.proj = other.proj
232	self.vmin = other.vmin ; self.vmax = other.vmax
233	self.div = other.div
234	self.colorb = other.colorb
235	self.lstyle = other.lstyle
236	self.marker = other.marker
237	self.color = other.color
238	self.label = other.label
239	self.title = other.title
240	else:
241	print "!! ERROR !! argument must be a pp object." ; exit()
242
243	# check the compatibility of two objects for operations
244	# --> if other is a pp class, test sizes and return isnum = False
245	# --> if other is an int or a float, return isnum = True
246	# --> otherwise, just print an error and exit
247	def checktwo(self,other):
248	if other.__class__.__name__ == "pp":
249	isnum = False
250	if self.status in ["init","defined"] or other.status in ["init","define"]:
251	print "!! ERROR !! Please use .retrieve to get fields for plots with one of your pp operands." ; exit()
252	if self.nfin != other.nfin or \
253	self.nvin != other.nvin or \
254	self.nplott != other.nplott or \
255	self.nplotz != other.nploty or \
256	self.nploty != other.nploty or \
257	self.nplotx != other.nplotx :
258	print "!! ERROR !! The two operands do not have the same number of files, variables, t z y x requests."
259	exit()
260	elif isinstance(other,int) or isinstance(other,float):
261	isnum = True
262	else:
263	print "!! ERROR !! The operand is neither a pp class nor an integer or a float." ; exit()
264	return isnum
265
266	# define a selective copy of a pp() object for operations
267	# ... copy.copy() is not conservative (still acts like a pointer)
268	# ... copy.deepcopy() does not work with netCDF objects
269	# so what is done here is a copy of everything except
270	# (to avoid sharing with self and therefore modifying self through operations)
271	# - request attribute of pp() object
272	# - field attribute of the onerequest() objects
273	def selective_copy(self):
274	if self.status in ["init","defined"]:
275	print "!! ERROR !! Please use .retrieve to get fields for the object you want to copy from." ; exit()
276	the_clone = pp()
277	for k, v in vars(self).items():
278	if k != "request":
279	setattr(the_clone,k,v)
280	the_clone.verbose = False
281	the_clone.filename = "THIS_IS_A_CLONE" # trick to avoid additional outputs
282	the_clone.define()
283	for i in range(self.nfin):
284	for j in range(self.nvin):
285	for t in range(self.nplott):
286	for z in range(self.nplotz):
287	for y in range(self.nploty):
288	for x in range(self.nplotx):
289	obj_ref = self.request[i][j][t][z][y][x]
290	obj = the_clone.request[i][j][t][z][y][x]
291	for k, v in vars(obj_ref).items():
292	if k != "field":
293	setattr(obj,k,v)
294	the_clone.status = "retrieved"
295	the_clone.filename = self.filename
296	return the_clone
297
298	# define the operation + on two objects. or with an int/float.
299	# ... with selective_copy the self object is not modified.
300	def __add__(self,other):
301	isnum = self.checktwo(other)
302	the_clone = self.selective_copy()
303	for i in range(self.nfin):
304	for j in range(self.nvin):
305	for t in range(self.nplott):
306	for z in range(self.nplotz):
307	for y in range(self.nploty):
308	for x in range(self.nplotx):
309	obj = the_clone.request[i][j][t][z][y][x]
310	obj_ref = self.request[i][j][t][z][y][x]
311	if not isnum:
312	ope = other.request[i][j][t][z][y][x].field
313	if obj_ref.field.shape != ope.shape:
314	print "!! ERROR !! The two fields for operation do not have the same shape.",self.field.shape,other.field.shape
315	exit()
316	else:
317	ope = other
318	goal = self.vargoal[j] + self.filegoal[i]
319	if ("vector" in goal) or ("contour" in goal):
320	if self.verbose: print "!! WARNING !! No operation was made on contours and vectors. This can be debatted actually."
321	obj.field = obj_ref.field
322	else:
323	obj.field = obj_ref.field + ope
324	return the_clone
325
326	# define the operation - on two objects. or with an int/float.
327	# ... with selective_copy the self object is not modified.
328	def __sub__(self,other):
329	isnum = self.checktwo(other)
330	the_clone = self.selective_copy()
331	for i in range(self.nfin):
332	for j in range(self.nvin):
333	for t in range(self.nplott):
334	for z in range(self.nplotz):
335	for y in range(self.nploty):
336	for x in range(self.nplotx):
337	obj = the_clone.request[i][j][t][z][y][x]
338	obj_ref = self.request[i][j][t][z][y][x]
339	if not isnum:
340	ope = other.request[i][j][t][z][y][x].field
341	if obj_ref.field.shape != ope.shape:
342	print "!! ERROR !! The two fields for operation do not have the same shape.",self.field.shape,other.field.shape
343	exit()
344	else:
345	ope = other
346	goal = self.vargoal[j] + self.filegoal[i]
347	if ("vector" in goal) or ("contour" in goal):
348	if self.verbose: print "!! WARNING !! No operation was made on contours and vectors. This can be debatted actually."
349	obj.field = obj_ref.field
350	else:
351	obj.field = obj_ref.field - ope
352	return the_clone
353
354	# define the operation * on two objects. or with an int/float.
355	# ... with selective_copy the self object is not modified.
356	def __mul__(self,other):
357	isnum = self.checktwo(other)
358	the_clone = self.selective_copy()
359	for i in range(self.nfin):
360	for j in range(self.nvin):
361	for t in range(self.nplott):
362	for z in range(self.nplotz):
363	for y in range(self.nploty):
364	for x in range(self.nplotx):
365	obj = the_clone.request[i][j][t][z][y][x]
366	obj_ref = self.request[i][j][t][z][y][x]
367	if not isnum:
368	ope = other.request[i][j][t][z][y][x].field
369	if obj_ref.field.shape != ope.shape:
370	print "!! ERROR !! The two fields for operation do not have the same shape.",self.field.shape,other.field.shape
371	exit()
372	else:
373	ope = other
374	goal = self.vargoal[j] + self.filegoal[i]
375	if ("vector" in goal) or ("contour" in goal):
376	if self.verbose: print "!! WARNING !! No operation was made on contours and vectors. This can be debatted actually."
377	obj.field = obj_ref.field
378	else:
379	obj.field = obj_ref.field * ope
380	return the_clone
381
382	# define the operation / on two objects. or with an int/float.
383	# ... with selective_copy the self object is not modified.
384	def __div__(self,other):
385	isnum = self.checktwo(other)
386	the_clone = self.selective_copy()
387	for i in range(self.nfin):
388	for j in range(self.nvin):
389	for t in range(self.nplott):
390	for z in range(self.nplotz):
391	for y in range(self.nploty):
392	for x in range(self.nplotx):
393	obj = the_clone.request[i][j][t][z][y][x]
394	obj_ref = self.request[i][j][t][z][y][x]
395	if not isnum:
396	ope = other.request[i][j][t][z][y][x].field
397	if obj_ref.field.shape != ope.shape:
398	print "!! ERROR !! The two fields for operation do not have the same shape.",self.field.shape,other.field.shape
399	exit()
400	else:
401	ope = other
402	goal = self.vargoal[j] + self.filegoal[i]
403	if ("vector" in goal) or ("contour" in goal):
404	if self.verbose: print "!! WARNING !! No operation was made on contours and vectors. This can be debatted actually."
405	obj.field = obj_ref.field
406	else:
407	obj.field = obj_ref.field / ope
408	return the_clone
409
410	# define the reverse operation float/int + object
411	def __radd__(self,other):
412	isnum = self.checktwo(other)
413	if not isnum: print "!! ERROR !! Operand should be a number" ; exit()
414	return self.__add__(other)
415
416	# define the reverse operation float/int - object
417	def __rsub__(self,other):
418	isnum = self.checktwo(other)
419	if not isnum: print "!! ERROR !! Operand should be a number" ; exit()
420	return self.__sub__(other)
421
422	# define the reverse operation float/int * object
423	def __rmul__(self,other):
424	isnum = self.checktwo(other)
425	if not isnum: print "!! ERROR !! Operand should be a number" ; exit()
426	return self.__mul__(other)
427
428	# define the reverse operation float/int / object
429	def __rdiv__(self,other):
430	isnum = self.checktwo(other)
431	if not isnum: print "!! ERROR !! Operand should be a number" ; exit()
432	return self.__div__(other)
433
434	# define the operation ** on one object.
435	# ... with selective_copy the self object is not modified.
436	def __pow__(self,num):
437	the_clone = self.selective_copy()
438	if isinstance(num,int) or isinstance(num,float):
439	for i in range(self.nfin):
440	for j in range(self.nvin):
441	for t in range(self.nplott):
442	for z in range(self.nplotz):
443	for y in range(self.nploty):
444	for x in range(self.nplotx):
445	obj = the_clone.request[i][j][t][z][y][x]
446	obj_ref = self.request[i][j][t][z][y][x]
447	goal = self.vargoal[j] + self.filegoal[i]
448	if ("vector" in goal) or ("contour" in goal):
449	if self.verbose: print "!! WARNING !! No operation was made on contours and vectors. This can be debatted actually."
450	obj.field = obj_ref.field
451	else:
452	obj.field = obj_ref.field ** num
453	else:
454	print "!! ERROR !! To define a power, either an int or a float is needed." ; exit()
455	return the_clone
456
457	### TBD: reverse power? for exponentials?
458
459	##############################################################################################
460	# define method
461	# ---------
462	# ... (file and var are either one string or a vector of strings)
463	# ... the goal of define is to define a 2D array of onerequest() objects (see class below)
464	# given the number of file, var, x, y, z, t asked by the user
465	# ... objectives for file or var are given through filegoal and vargoal
466	# --> possible values: main contour vector
467	# ---------
468	# ... then onerequest() objects are being defined more precisely
469	# by getting index_x index_y index_z index_t
470	# and setting method_x method_y method_z method_t to either
471	# - "free" for free dimensions (plot dimensions)
472	# - "comp" for averages, max, min
473	# - "fixed" for fixed dimensions (possibly several i.e. multislice)
474	##############################################################################################
475	def define(self):
476	self.printstatus()
477	# initial check and get dimensions
478	self.file = checktab(self.file,mess="file")
479	self.nfin = len(self.file)
480	if self.verbose:
481	for i in range(self.nfin): inspect(self.file[i])
482	self.var = checktab(self.var,mess="var")
483	self.nvin = len(self.var)
484	# check goal tabs for files and variables
485	# ... default is to plot everything
486	if self.filegoal is None: self.filegoal = ["main"]*self.nfin
487	if self.vargoal is None: self.vargoal = ["main"]*self.nvin
488	self.filegoal = checktab(self.filegoal, mess="filegoal")
489	self.vargoal = checktab(self.vargoal, mess="vargoal")
490	if len(self.filegoal) != self.nfin: print "!! ERROR !! filegoal must be the same size as file." ; exit()
491	if len(self.vargoal) != self.nvin: print "!! ERROR !! vargoal must be the same size as var." ; exit()
492	# variables: initial check
493	self.x = checktab(self.x,mess="x",allownone=True,allownumber=True)
494	self.y = checktab(self.y,mess="y",allownone=True,allownumber=True)
495	self.z = checktab(self.z,mess="z",allownone=True,allownumber=True)
496	self.t = checktab(self.t,mess="t",allownone=True,allownumber=True)
497	# for the moment not var- nor file- dependent.
498	# but this could be the case.
499	sx = readslices(self.x) ; sy = readslices(self.y)
500	sz = readslices(self.z) ; st = readslices(self.t)
501	# get methods
502	mx = findmethod(sx) ; my = findmethod(sy)
503	mz = findmethod(sz) ; mt = findmethod(st)
504	# get number of plots to be done
505	if mx == "fixed": self.nplotx = sx.size/2
506	else: self.nplotx = 1
507	if my == "fixed": self.nploty = sy.size/2
508	else: self.nploty = 1
509	if mz == "fixed": self.nplotz = sz.size/2
510	else: self.nplotz = 1
511	if mt == "fixed": self.nplott = st.size/2
512	else: self.nplott = 1
513	if self.verbose: print "**** OK. Plots over x,y,z,t -->",self.nplotx,self.nploty,self.nplotz,self.nplott
514	# create the list of onerequest() objects
515	self.request = [[[[[[ \
516	onerequest() \
517	for x in range(self.nplotx)] for y in range(self.nploty)] \
518	for z in range(self.nplotz)] for t in range(self.nplott)] \
519	for j in range(self.nvin)] for i in range(self.nfin)]
520	# loop on onerequest() objects
521	for i in range(self.nfin):
522	for j in range(self.nvin):
523	for t in range(self.nplott):
524	for z in range(self.nplotz):
525	for y in range(self.nploty):
526	for x in range(self.nplotx):
527	obj = self.request[i][j][t][z][y][x]
528	# fill in names for files and variables
529	obj.verbose = self.verbose
530	obj.file = self.file[i]
531	obj.var = self.var[j]
532	# indicate the computation method
533	obj.compute = self.compute
534	# open the files (the same file might be opened several times but this is cheap)
535	obj.openfile()
536	### get x,y,z,t dimensions from file
537	obj.getdim()
538	### get methods
539	obj.method_x = mx ; obj.method_y = my
540	obj.method_z = mz ; obj.method_t = mt
541	### get index
542	obj.getindextime(dalist=st,ind=t,stride=self.stridet)
543	obj.getindexvert(dalist=sz,ind=z,stride=self.stridez)
544	obj.getindexhori(dalistx=sx,dalisty=sy,indx=x,indy=y,stridex=self.stridex,stridey=self.stridey)
545	# change status
546	self.status = "defined"
547	return self
548
549	##############################################################################################
550	# retrieve method
551	# --> for each element onerequest() in the array, get field .var from .f file
552	# --> see below the onerequest() class:
553	# - only get what is needed for computing and plotting
554	# - averages etc... are computed here
555	# --> RESULT: each onerequest() object has now its attribute .field filled
556	# --> if one wants to perform operations on fields, this should be done after retrieve()
557	##############################################################################################
558	def retrieve(self):
559	self.printstatus()
560	# check if things were done OK before
561	if self.status != "defined": print "!! ERROR !! Please use .define() to define your pp object." ; exit()
562	## first get fields
563	## ... only what is needed is extracted from the files
564	## ... and computations are performed
565	for i in range(self.nfin):
566	for j in range(self.nvin):
567	for t in range(self.nplott):
568	for z in range(self.nplotz):
569	for y in range(self.nploty):
570	for x in range(self.nplotx):
571	obj = self.request[i][j][t][z][y][x]
572	obj.getfield()
573	obj.computations()
574	# change status
575	self.status = "retrieved"
576	return self
577
578	##########################################################
579	# get: a shortcut method for the define + retrieve chain #
580	##########################################################
581	def get(self):
582	self.define()
583	self.retrieve()
584	return self
585
586	########################################
587	# smooth: smooth the field in 1D or 2D #
588	########################################
589	## TBD: smooth not OK with masked array in the end of retrieve()
590	def smooth(self,window):
591	if self.verbose:
592	print "!! WARNING !! Performing a smoothing with a window size",window
593	print "!! WARNING !! To come back to unsmoothed file, use .get() again"
594	for i in range(self.nfin):
595	for j in range(self.nvin):
596	for t in range(self.nplott):
597	for z in range(self.nplotz):
598	for y in range(self.nploty):
599	for x in range(self.nplotx):
600	obj = self.request[i][j][t][z][y][x]
601	if obj.field.ndim == 1:
602	print "!! ERROR !! 1D smoothing not supported yet because reduces array sizes."
603	exit()
604	# TBD TBD TBD
605	#obj.field = ppcompute.smooth1d(obj.field,window=window)
606	elif obj.field.ndim == 2:
607	obj.field = ppcompute.smooth2d(obj.field,window=window)
608
609	##############################################################################################
610	# defineplot method
611	# --> defineplot first defines arrays of plot objects and set each of them
612	# ... simple looping except cases where goal is not main (e.g. contour or vector)
613	# --> principle: each onerequest() object gives birth to a subplot
614	# --> defineplot vs. makeplot: defining plot and actually plotting it are clearly separated
615	# --> THE KEY OUPUT OF defineplot IS AN ARRAY self.p OF PLOT OBJECTS
616	# optional arguments
617	# --> extraplot: to indicate a number of plots to be added afterwards (use self.plotin)
618	# --> loadfile: to use self.p from a previously saved file
619	##############################################################################################
620	def defineplot(self,extraplot=0,loadfile=None):
621	# -----------------------------------------------------
622	# LOAD MODE: load a self.p object. count plots from it.
623	# -----------------------------------------------------
624	if loadfile is not None:
625	try: filehandler = open(loadfile, 'r') ; self.p = pickle.load(filehandler)
626	except IOError: print "!! ERROR !! Cannot find object file to load." ; exit()
627	self.status = "definedplot" ; self.plotin = None
628	self.nplot = len(self.p) ; self.howmanyplots = self.nplot
629	return
630	# -----------------------------------------------------
631	# REGULAR MODE
632	# -----------------------------------------------------
633	self.printstatus()
634	# check if things were done OK before
635	if self.status in ["init","defined"]:
636	print "!! ERROR !! Please use .retrieve() to get fields for plots with your pp object." ; exit()
637	# check self.plotin (an existing fig on which to add plots afterwards)
638	if self.plotin.__class__.__name__ == "pp":
639	if self.plotin.fig is None:
640	self.plotin = None # this is an additional security in case
641	# a pp object is given without figure opened yet.
642	elif self.plotin is not None:
643	print "!! ERROR !! plotin argument must be a pp object." ; exit()
644	# initialize the array of subplot objects
645	# either something new or attributes coming from plotin object
646	if self.plotin is None: self.p = [ ]
647	else: self.p = self.plotin.p
648	# create an array of subplot objects
649	# ... in theory the order of looping can be changed without any harm
650	# ... the only important thing is to keep i,j,t,z,y,x resp. for file,var,t,z,y,x
651	count = 0
652	for i in range(self.nfin):
653	if self.filegoal[i] == "main":
654	for j in range(self.nvin):
655	if self.vargoal[j] == "main":
656	for t in range(self.nplott):
657	for z in range(self.nplotz):
658	for y in range(self.nploty):
659	for x in range(self.nplotx):
660	# look at dimension and append the right plot object
661	obj = self.request[i][j][t][z][y][x]
662	dp = obj.dimplot
663	if dp == 1 or self.forcedimplot == 1: plobj = ppplot.plot1d()
664	elif dp == 2 or self.forcedimplot == 2: plobj = ppplot.plot2d()
665	elif dp == 0: print "**** OK. VALUES VALUES VALUES",obj.field
666	else: print "!! ERROR !! 3D or 4D plots not supported" ; exit()
667	# load abscissa and ordinate in obj
668	obj.definecoord()
669	# we start to define things here before appending
670	# (convenient: could be overridden by the user before makeplot)
671	# ... the if loop is necessary so that we can loop above on the dp=0 case
672	if dp in [1,2]:
673	# and define what to do in plobj
674	plobj.invert = obj.invert_axes
675	plobj.swap = obj.swap_axes
676	# axis labels
677	plobj.xlabel = obj.absclab ; plobj.ylabel = obj.ordilab
678	# superpose or not (this is mostly for saving purpose)
679	plobj.superpose = self.superpose
680	# get title, colormaps, labels, etc.. from var
681	plobj.var = obj.var
682	plobj.define_from_var()
683	# generic 1D/2D: load field and coord in plot object
684	plobj.field = obj.field # field to be plotted
685	plobj.absc = obj.absc # abscissa (or longitude)
686	plobj.ordi = obj.ordi # ordinate (or latitude)
687	# -- useless in 1D but not used anyway
688	# specific 1D plot stuff
689	if dp == 1:
690	# -- a default label
691	plobj.label = ""
692	if self.nfin > 1: plobj.label = plobj.label + " file #"+str(i+1)
693	if self.nvin > 1: plobj.label = plobj.label + " var #"+str(j+1)
694	if self.nplott > 1: plobj.label = plobj.label + " t #"+str(t+1)
695	if self.nplotz > 1: plobj.label = plobj.label + " z #"+str(z+1)
696	if self.nploty > 1: plobj.label = plobj.label + " y #"+str(y+1)
697	if self.nplotx > 1: plobj.label = plobj.label + " x #"+str(x+1)
698	# specific 2d plot stuff
699	if dp == 2:
700	# -- light grey background for missing values
701	if type(plobj.field).__name__ in 'MaskedArray': plobj.axisbg = '0.75'
702	# -- define if it is a map or a plot
703	plobj.mapmode = ( obj.method_x+obj.method_y == "freefree" \
704	and "grid points" not in obj.name_x \
705	and not self.noproj )
706	# possible user-defined plot settings shared by all plots
707	if self.div is not None: plobj.div = self.div
708	if self.xlabel is not None: plobj.xlabel = self.xlabel
709	if self.xcoeff is not None: plobj.xcoeff = self.xcoeff
710	if self.ylabel is not None: plobj.ylabel = self.ylabel
711	if self.ycoeff is not None: plobj.ycoeff = self.ycoeff
712	if self.title is not None: plobj.title = self.title
713	# -- 1D specific
714	if dp == 1:
715	if self.lstyle is not None: plobj.lstyle = self.lstyle
716	if self.marker is not None: plobj.marker = self.marker
717	if self.color is not None: plobj.color = self.color
718	if self.label is not None: plobj.label = self.label
719	# -- 2D specific
720	elif dp == 2:
721	if self.proj is not None and not self.noproj: plobj.proj = self.proj
722	if self.vmin is not None: plobj.vmin = self.vmin
723	if self.vmax is not None: plobj.vmax = self.vmax
724	if self.colorb is not None: plobj.colorb = self.colorb
725	# finally append plot object
726	self.p.append(plobj)
727	count = count + 1
728	# self.nplot is number of plot to be defined in this call to defineplot()
729	# (because of self.plotin this might less than length of self.p)
730	self.nplot = count
731	# --- superimposed contours and vectors ---
732	# we have to start another loop because we need forward information
733	# TBD: there is probably a more flexible way to do that
734	count = 0
735	for i in range(self.nfin):
736	for j in range(self.nvin):
737	for t in range(self.nplott):
738	for z in range(self.nplotz):
739	for y in range(self.nploty):
740	for x in range(self.nplotx):
741	goal = self.vargoal[j] + self.filegoal[i]
742	obj = self.request[i][j][t][z][y][x]
743	if "mainmain" in goal and obj.dimplot == 2:
744	# the plot object we consider in the loop
745	pl = self.p[count]
746	# -- see if there is a contour requested...
747	# (we use try because we might be at the end of the list)
748	found = 0
749	try: condvar = self.vargoal[j+1]
750	except: condvar = "itisok"
751	try: condfile = self.filegoal[i+1]
752	except: condfile = "itisok"
753	# ... get contour
754	##########################################
755	# NB: contour is expected to be right after main otherwise it is not displayed
756	##########################################
757	if condvar == "contour":
758	plobj.addcontour = self.request[i][j+1][t][z][y][x].field ; found += 1
759	if condfile == "contour":
760	plobj.addcontour = self.request[i+1][j][t][z][y][x].field ; found += 1
761	# see if there is a vector requested...
762	# (we use try because we might be at the end of the list)
763	try: condvar = self.vargoal[j+found+1]+self.vargoal[j+found+2]
764	except: condvar = "itisok"
765	try: condfile = self.filegoal[i+found+1]+self.filegoal[i+found+2]
766	except: condfile = "itisok"
767	# ... get vector and go directly to the next iteration
768	# (in some cases we would do this twice but this is cheap)
769	if "vector" in condvar:
770	plobj.addvecx = self.request[i][j+found+1][t][z][y][x].field
771	plobj.addvecy = self.request[i][j+found+2][t][z][y][x].field
772	if "vector" in condfile:
773	plobj.addvecx = self.request[i+found+1][j][t][z][y][x].field
774	plobj.addvecy = self.request[i+found+2][j][t][z][y][x].field
775	count = count + 1
776	# COUNT PLOTS. if 0 just exit.
777	# self.howmanyplots is self.nplot + possible extraplots
778	self.howmanyplots = self.nplot + extraplot
779	if self.howmanyplots > 0:
780	if self.verbose: print "**** OK. expect %i plots" % (self.howmanyplots)
781	else:
782	exit() # because this means that we only had 0D values !
783	# final status
784	self.status = "definedplot"
785	return self
786
787	##############################################################################################
788	# makeplot method
789	# --> after defineplot and before makeplot, user-defined plot settings can be easily given
790	# simply by modifying the attributes of each elements of self.p
791	# --> to change only one plot setting, no need to call defineplot again before makeplot
792	# --> in the end, the array self.p of plot objects is saved for easy and convenient replotting
793	# --> see practical examples in the folder 'examples'
794	##############################################################################################
795	def makeplot(self):
796	self.printstatus()
797	# a few initial operations
798	# ------------------------
799	if "definedplot" not in self.status:
800	print "!! ERROR !! Please use .defineplot() before .makeplot() can be used with your pp object." ; exit()
801	# open a figure and define subplots
802	# ---------------------------------
803	if self.plotin is None:
804	# start from scratch
805	self.fig = mpl.figure(figsize=(16,8))
806	self.subv,self.subh = ppplot.definesubplot(self.howmanyplots,self.fig)
807	self.n = 0
808	## adapted space for labels etc
809	## ... except for ortho because there is no label anyway
810	self.customplot = self.p[0].field.ndim == 2 \
811	and self.p[0].mapmode == True \
812	and self.p[0].proj not in ["ortho"]
813	if self.customplot:
814	margin = 0.07
815	self.fig.subplots_adjust(left=margin,right=1-margin,bottom=margin,top=1-margin)
816	else:
817	# start from an existing figure.
818	# extraplot must have been set in the call to the previous figure.
819	self.fig = self.plotin.fig
820	self.subv,self.subh = self.plotin.subv,self.plotin.subh
821	self.n = self.plotin.n
822	self.howmanyplots = self.plotin.howmanyplots
823	self.customplot = self.plotin.customplot
824	# LOOP on all subplots
825	# NB: cannot use 'for pl in self.p' if self.plotin not None
826	# --------------------
827	for count in range(self.nplot):
828	# the plot object we consider in the loop
829	pl = self.p[self.n]
830	# before making the plot, create a subplot. the first one is numbered 1 not 0.
831	# ... if pl.superpose, we use only one and only figure
832	# ... (and we have to be careful with not doing things several times)
833	if pl.superpose:
834	if self.n == 0:
835	self.fig.add_subplot(1,1,1,axisbg=pl.axisbg) # define one subplot (still needed for user-defined font sizes)
836	sav = pl.xlabel,pl.ylabel,pl.xcoeff,pl.ycoeff,pl.title,pl.swaplab # save titles and labels
837	else:
838	pl.invert = False ; pl.lstyle = None # don't invert again axis
839	# set saved titles and labels
840	if self.plotin is None:
841	pl.xlabel,pl.ylabel,pl.xcoeff,pl.ycoeff,pl.title,pl.swaplab = sav
842	else:
843	prev_plot = self.plotin.p[self.n-1]
844	pl.xlabel = prev_plot.xlabel
845	pl.ylabel = prev_plot.ylabel
846	pl.xcoeff = prev_plot.xcoeff
847	pl.ycoeff = prev_plot.ycoeff
848	pl.title = prev_plot.title
849	pl.swaplab = prev_plot.swaplab
850	else:
851	self.fig.add_subplot(self.subv,self.subh,self.n+1,axisbg=pl.axisbg)
852	if self.verbose: print "**** Done subplot %i / %i " %( self.n+1,self.howmanyplots )
853	# finally make the plot
854	pl.make()
855	# increment plot count (and propagate this in plotin)
856	self.n = self.n+1
857	if self.plotin is not None: self.plotin.n = self.n
858	# once completed show the plot (cannot show intermediate plotin)
859	# ... added a fix (customplot=True) for the label problem in basemap
860	print "**** PPCLASS. Done step: makeplot"
861	if (self.n == self.howmanyplots):
862	ppplot.save(mode=self.out,filename=self.filename,folder=self.folder,custom=self.customplot)
863	mpl.close()
864	# SAVE A PICKLE FILE WITH THE self.p ARRAY OF OBJECTS
865	if self.verbose: print "**** Saving session in "+self.filename + ".ppobj"
866	savfile = self.folder + "/" + self.filename + ".ppobj"
867	try:
868	filehandler = open(savfile, 'w')
869	pickle.dump(self.p, filehandler)
870	except IOError:
871	print "!! WARNING !! Saved object file not written. Probably do not have permission to write here."
872	return self
873
874	###########################################################
875	# plot: a shortcut method for the defineplot + plot chain #
876	###########################################################
877	def plot(self,extraplot=0):
878	self.defineplot(extraplot=extraplot)
879	self.makeplot()
880	return self
881
882	#######################################################
883	# getplot: a shortcut method for the get + plot chain #
884	#######################################################
885	def getplot(self,extraplot=0):
886	self.get()
887	self.plot(extraplot=extraplot)
888	return self
889
890	###################################################################
891	# getdefineplot: a shortcut method for the get + defineplot chain #
892	###################################################################
893	def getdefineplot(self,extraplot=0):
894	self.get()
895	self.defineplot(extraplot=extraplot)
896	return self
897
898	##############################################################
899	# f: operation on two pp objects being on status 'definedplot'
900	# this allows for one field being function of another one
901	# e.g. u.f(v) means u will be displayed as a function of v
902	# ... no need to do defineplot after u.f(v), makeplot directly
903	##############################################################
904	def f(self,other):
905	# preamble: for this operation to work, defineplot() must have been done
906	if self.status != "definedplot":
907	if self.verbose: print "!! WARNING !! performing defineplot on operand"
908	self.defineplot()
909	if other.status != "definedplot":
910	if self.verbose: print "!! WARNING !! performing defineplot on operand"
911	other.defineplot()
912	# check total number of plots
913	if self.howmanyplots != other.howmanyplots:
914	print "!! ERROR !! The two operands do not have the same number of subplots."
915	exit()
916	# and now operation.
917	count = 0
918	while count < self.howmanyplots:
919	sobj = self.p[count] ; oobj = other.p[count]
920	if sobj.field.ndim !=1 or oobj.field.ndim !=1:
921	if self.verbose: print "!! WARNING !! Flattening arrays because more than one-dimensional."
922	sobj.field = np.ravel(sobj.field)
923	oobj.field = np.ravel(oobj.field)
924	sobj.absc = oobj.field
925	sobj.xlabel = oobj.ylabel
926	if sobj.absc.size > sobj.field.size:
927	if self.verbose:
928	print "!! WARNING !! Trying to define y=f(x) with x and y not at the same size.",sobj.absc.size,sobj.field.size
929	print "!! WARNING !! Modifying x to fit y size but please check."
930	sobj.absc = sobj.absc[0:sobj.field.size]
931	count = count + 1
932	return self
933
934	###########################################################
935	# copyopt: get options from e.g. a parser
936	# ... allow for simple scripting and user-defined settings
937	# ... must be called between defineplot and makeplot
938	# REQUIRED: attributes of opt must be the same as in the pp object
939	###########################################################
940	def getopt(self,opt):
941	# -- if only one, or less than the number of plots --> we take the first one
942	# -- if as many as number of plots --> OK, each plot has its own setting
943	# (except a few cases such as trans)
944	for iii in range(self.howmanyplots):
945	###
946	try: self.p[iii].trans = opt.trans
947	except: pass
948	###
949	try: self.p[iii].div = opt.div
950	except: pass
951	###
952	try: self.p[iii].logy = opt.logy
953	except: pass
954	###
955	try: self.p[iii].colorb = opt.colorb[iii]
956	except:
957	try: self.p[iii].colorb = opt.colorb[0]
958	except: pass
959	###
960	try: self.p[iii].title = opt.title[iii]
961	except:
962	try: self.p[iii].title = opt.title[0]
963	except: pass
964	###
965	try: self.p[iii].xlabel = opt.xlabel[iii]
966	except:
967	try: self.p[iii].xlabel = opt.xlabel[0]
968	except: pass
969	###
970	try: self.p[iii].ylabel = opt.ylabel[iii]
971	except:
972	try: self.p[iii].ylabel = opt.ylabel[0]
973	except: pass
974	###
975	try: self.p[iii].lstyle = opt.lstyle[iii]
976	except:
977	try: self.p[iii].lstyle = opt.lstyle[0]
978	except: pass
979	###
980	try: self.p[iii].color = opt.color[iii]
981	except:
982	try: self.p[iii].color = opt.color[0]
983	except: pass
984	###
985	try: self.p[iii].marker = opt.marker[iii]
986	except:
987	try: self.p[iii].marker = opt.marker[0]
988	except: pass
989	###
990	try: self.p[iii].label = opt.label[iii]
991	except:
992	try: self.p[iii].label = opt.label[0]
993	except: pass
994	###
995	try: self.p[iii].proj = opt.proj[iii]
996	except:
997	try: self.p[iii].proj = opt.proj[0]
998	except: pass
999	###
1000	try: self.p[iii].back = opt.back[iii]
1001	except:
1002	try: self.p[iii].back = opt.back[0]
1003	except: pass
1004	###
1005	try: self.p[iii].area = opt.area[iii]
1006	except:
1007	try: self.p[iii].area = opt.area[0]
1008	except: pass
1009	###
1010	try: self.p[iii].blon = opt.blon[iii]
1011	except:
1012	try: self.p[iii].blon = opt.blon[0]
1013	except: pass
1014	###
1015	try: self.p[iii].blat = opt.blat[iii]
1016	except:
1017	try: self.p[iii].blat = opt.blat[0]
1018	except: pass
1019	###
1020	try: self.p[iii].vmin = opt.vmin[iii]
1021	except:
1022	try: self.p[iii].vmin = opt.vmin[0]
1023	except: pass
1024	###
1025	try: self.p[iii].vmax = opt.vmax[iii]
1026	except:
1027	try: self.p[iii].vmax = opt.vmax[0]
1028	except: pass
1029
1030	##########################################################
1031	### THE ONEREQUEST SUBOBJECT TO PP (ON WHICH IT LOOPS) ###
1032	##########################################################
1033	class onerequest():
1034
1035	# default settings. mostly initialized to diagnose problem, except dimplot, nplot, verbose, swap_axes, invert_axes
1036	# -------------------------------
1037	def __init__(self):
1038	self.file = '!! file: I am not set, damned !!'
1039	self.f = None
1040	self.dim = None
1041	self.var = '!! var: I am not set, damned !!'
1042	self.index_x = [] ; self.index_y = [] ; self.index_z = [] ; self.index_t = []
1043	self.index_x2d = [] ; self.index_y2d = []
1044	self.method_x = '!! method_x: I am not set, damned !!'
1045	self.method_y = '!! method_y: I am not set, damned !!'
1046	self.method_z = '!! method_z: I am not set, damned !!'
1047	self.method_t = '!! method_t: I am not set, damned !!'
1048	self.field = None
1049	self.name_x = None ; self.name_y = None ; self.name_z = None ; self.name_t = None
1050	self.dim_x = None ; self.dim_y = None ; self.dim_z = None ; self.dim_t = None
1051	self.field_x = None ; self.field_y = None ; self.field_z = None ; self.field_t = None
1052	self.dimplot = 0
1053	self.nplot = 1
1054	self.absc = None ; self.ordi = None ; self.absclab = None ; self.ordilab = None
1055	self.verbose = True
1056	self.swap_axes = False ; self.invert_axes = False
1057	self.compute = None
1058
1059	# open a file. for now it is netcdf. TBD for other formats.
1060	# check that self.var is inside.
1061	# -------------------------------
1062	def openfile(self):
1063	if not os.path.exists(self.file): print '!! ERROR !! I could not find the following file: '+self.file ; exit()
1064	if not os.path.isfile(self.file): print '!! ERROR !! This does not appear to be a file: '+self.file ; exit()
1065	self.f = netCDF4.Dataset(self.file)
1066	if self.verbose: print "**** OK. Opened file "+self.file
1067	if self.var not in self.f.variables.keys():
1068	print '!! ERROR !! File '+self.file+' does not contain variable: '+self.var
1069	print '..... try instead with ',self.f.variables.keys() ; exit()
1070
1071	# copy attributes from another existing object
1072	# --------------------------------------------
1073	def copy(self,source):
1074	for k, v in vars(source).items():
1075	setattr(self,k,v)
1076
1077	# get x,y,z,t dimensions from NETCDF file
1078	# TBD: user could request for a specific altitude dimension
1079	# TBD: staggered variables could request specific dimensions
1080	# -------------------------------
1081	def getdim(self):
1082	# GET SIZES OF EACH DIMENSION
1083	if self.verbose: print "**** OK. Found variable "+self.var
1084	shape = self.f.variables[self.var].shape
1085	self.dim = len(shape)
1086	if self.dim == 1:
1087	if self.verbose: print "**** OK. 1D field. I assume this varies with time."
1088	self.dim_x = 1 ; self.dim_y = 1 ; self.dim_z = 1 ; self.dim_t = shape[0]
1089	elif self.dim == 2:
1090	if self.verbose: print "**** OK. 2D field. I assume this is not-time-varying lat-lon map."
1091	self.dim_x = shape[1] ; self.dim_y = shape[0] ; self.dim_z = 1 ; self.dim_t = 1
1092	elif self.dim == 3:
1093	if self.verbose: print "**** OK. 3D field. I assume this is time-varying lat-lon map."
1094	self.dim_x = shape[2] ; self.dim_y = shape[1] ; self.dim_z = 1 ; self.dim_t = shape[0]
1095	elif self.dim == 4:
1096	if self.verbose: print "**** OK. 4D field."
1097	self.dim_x = shape[3] ; self.dim_y = shape[2] ; self.dim_z = shape[1] ; self.dim_t = shape[0]
1098	# LONGITUDE. Try preset fields. If not present set grid points axis.
1099	self.name_x = "nothing"
1100	for c in glob_listx:
1101	if c in self.f.variables.keys():
1102	self.name_x = c
1103	if self.name_x == "nothing":
1104	self.field_x = np.array(range(self.dim_x))
1105	self.name_x = "x grid points"
1106	else:
1107	self.field_x = self.f.variables[self.name_x]
1108	# LATITUDE. Try preset fields. If not present set grid points axis.
1109	self.name_y = "nothing"
1110	for c in glob_listy:
1111	if c in self.f.variables.keys():
1112	self.name_y = c
1113	if self.name_y == "nothing":
1114	self.field_y = np.array(range(self.dim_y))
1115	self.name_y = "y grid points"
1116	else:
1117	self.field_y = self.f.variables[self.name_y]
1118	# ensure that lon and lat are 2D fields
1119	# 1. simple 1D case (not time-varying)
1120	if len(self.field_x.shape)*len(self.field_y.shape) == 1:
1121	if self.verbose: print "**** OK. recasting lon and lat as 2D fields."
1122	[self.field_x,self.field_y] = np.meshgrid(self.field_x,self.field_y)
1123	# 2. complex 3D case (time-varying, actually just copied over time axis)
1124	elif len(self.field_x.shape)*len(self.field_y.shape) == 9:
1125	if self.verbose: print "**** OK. reducing lon and lat as 2D fields. get rid of time."
1126	self.field_x = self.field_x[0,:,:]
1127	self.field_y = self.field_y[0,:,:]
1128	# if xy axis are apparently undefined, set 2D grid points axis.
1129	if "grid points" not in self.name_x:
1130	if self.field_x.all() == self.field_x[0,0]:
1131	print "!! WARNING !! xy axis look undefined. creating non-dummy ones."
1132	self.field_x = np.array(range(self.dim_x)) ; self.name_x = "x grid points"
1133	self.field_y = np.array(range(self.dim_y)) ; self.name_y = "y grid points"
1134	[self.field_x,self.field_y] = np.meshgrid(self.field_x,self.field_y)
1135	if self.dim_x > 1:
1136	if self.verbose: print "**** OK. x axis %4.0f values [%5.1f,%5.1f]" % (self.dim_x,self.field_x.min(),self.field_x.max())
1137	if self.dim_y > 1:
1138	if self.verbose: print "**** OK. y axis %4.0f values [%5.1f,%5.1f]" % (self.dim_y,self.field_y.min(),self.field_y.max())
1139	# ALTITUDE. Try preset fields. If not present set grid points axis.
1140	# WARNING: how do we do if several are available?
1141	self.name_z = "nothing"
1142	for c in glob_listz:
1143	if c in self.f.variables.keys():
1144	self.name_z = c
1145	if self.name_z == "nothing":
1146	self.field_z = np.array(range(self.dim_z))
1147	self.name_z = "z grid points"
1148	else:
1149	self.field_z = self.f.variables[self.name_z][:] # specify dimension
1150	# TBD: have to check that this is not a 3D field
1151	if self.dim_z > 1:
1152	if self.verbose: print "**** OK. z axis %4.0f values [%5.1f,%5.1f]" % (self.dim_z,self.field_z.min(),self.field_z.max())
1153	# TIME. Try preset fields.
1154	self.name_t = "nothing"
1155	for c in glob_listt:
1156	if c in self.f.dimensions.keys():
1157	self.name_t = c
1158	try:
1159	# speed up: only get first value, last one.
1160	dafirst = self.f.variables[self.name_t][0]
1161	dalast = self.f.variables[self.name_t][self.dim_t-1]
1162	self.field_t = np.linspace(dafirst,dalast,num=self.dim_t)
1163	if dafirst == dalast: self.field_t = np.array([dafirst])
1164	except:
1165	# ... or if a problem encountered, define a simple time axis
1166	self.field_t = np.array(range(self.dim_t))
1167	self.name_t = "t grid points"
1168	if self.dim_t > 1:
1169	if self.verbose: print "**** OK. t axis %4.0f values [%5.1f,%5.1f]" % (self.dim_t,self.field_t.min(),self.field_t.max())
1170
1171	# get list of index to be retrieved for time axis
1172	### TBD: il faudrait ne prendre que les indices qui correspondent a l interieur d un plot (dans all)
1173	# -------------------------------
1174	def getindextime(self,dalist=None,ind=None,stride=1):
1175	if self.method_t == "free":
1176	self.index_t = np.arange(0,self.dim_t,stride)
1177	if self.dim_t > 1:
1178	self.dimplot = self.dimplot + 1
1179	if self.verbose: print "**** OK. t values. all."
1180	else:
1181	self.method_t = "fixed"
1182	if self.verbose: print "**** OK. no t dimension."
1183	elif self.method_t == "comp":
1184	start = np.argmin( np.abs( self.field_t - dalist[ind][0] ) )
1185	stop = np.argmin( np.abs( self.field_t - dalist[ind][1] ) )
1186	self.index_t = np.arange(start,stop,stride)
1187	if self.verbose: print "**** OK. t values. comp over interval ",self.field_t[start],self.field_t[stop]," nvalues=",self.index_t.size
1188	elif self.method_t == "fixed":
1189	self.index_t.append( np.argmin( np.abs( self.field_t - dalist[ind][0] ) ))
1190	if self.verbose: print "**** OK. t values",self.field_t[self.index_t]
1191	else:
1192	print "!! ERROR !! method "+self.method_t+" not supported"
1193	self.index_t = np.array(self.index_t)
1194
1195	# get list of index to be retrieved for vertical axis
1196	### TBD: il faudrait ne prendre que les indices qui correspondent a l interieur d un plot (dans all)
1197	# -------------------------------
1198	def getindexvert(self,dalist=None,ind=None,stride=1):
1199	if self.method_z == "free":
1200	self.index_z = np.arange(0,self.dim_z,stride)
1201	if self.dim_z > 1:
1202	self.dimplot = self.dimplot + 1
1203	if self.verbose: print "**** OK. z values. all."
1204	else:
1205	self.method_z = "fixed"
1206	if self.verbose: print "**** OK. no z dimension."
1207	elif self.method_z == "comp":
1208	start = np.argmin( np.abs( self.field_z - dalist[ind][0] ) )
1209	stop = np.argmin( np.abs( self.field_z - dalist[ind][1] ) )
1210	self.index_z = np.arange(start,stop,stride)
1211	if self.verbose: print "**** OK. z values. comp over interval",self.field_z[start],self.field_z[stop]," nvalues=",self.index_z.size
1212	elif self.method_z == "fixed":
1213	self.index_z.append( np.argmin( np.abs( self.field_z - dalist[ind][0] ) ))
1214	if self.verbose: print "**** OK. z values",self.field_z[self.index_z]
1215	else:
1216	if self.verbose: print "!! ERROR !! method "+self.method_z+" not supported"
1217	self.index_z = np.array(self.index_z)
1218
1219	# get list of index to be retrieved for horizontal grid
1220	# --> index_x and index_y are slices to be retrieved from NETCDF files
1221	# --> index_x2d and index_y2d are the actual (x,y) coordinates corresponding to each relevant point
1222	# [this is slightly more complicated because 2D arrays for lat-lon projection possibly irregular]
1223	# NB: to append index we use lists (the most convenient) then we convert into a numpy.array
1224	### TBD: il faudrait ne prendre que les indices qui correspondent a l interieur d un plot (dans all)
1225	# -------------------------------
1226	def getindexhori(self,dalistx=None,dalisty=None,indx=None,indy=None,stridex=1,stridey=1):
1227	## get what is the method over x and y axis
1228	test = self.method_x+self.method_y
1229	## CASE 0, EASY CASES:
1230	## - LAT IS FREE (we do here what must be done whatever LON case is)
1231	## - LON IS FREE (we do here what must be done whatever LAT case is)
1232	## - LAT IS COMP AND LON IS FREE
1233	## - LON IS COMP AND LAT IS FREE
1234	if self.method_x == "free" or test in ["compfree","compcomp"]:
1235	self.index_x = range(0,self.dim_x,stridex)
1236	if self.dim_x > 1:
1237	if self.method_x == "free": self.dimplot = self.dimplot + 1
1238	if self.verbose: print "**** OK. x values. all."
1239	else:
1240	self.method_x = "fixed"
1241	if self.verbose: print "**** OK. no x dimension."
1242	if self.method_y == "free" or test in ["freecomp","compcomp"]:
1243	self.index_y = range(0,self.dim_y,stridey)
1244	if self.dim_y > 1:
1245	if self.method_y == "free": self.dimplot = self.dimplot + 1
1246	if self.verbose: print "**** OK. y values. all."
1247	else:
1248	self.method_y = "fixed"
1249	if self.verbose: print "**** OK. no y dimension."
1250	## CASE 0 above, this is just for continuity.
1251	if self.method_x in ["free","comp"] and self.method_y in ["free","comp"]:
1252	self.index_x2d = self.index_x
1253	self.index_y2d = self.index_y
1254	## AND NOW THE LITTLE BIT MORE COMPLICATED CASES
1255	## CASE 1 LAT AND LON ARE FIXED
1256	elif test == "fixedfixed":
1257	idy,idx = np.unravel_index( np.argmin( ( self.field_x - dalistx[indx][0])2 + (self.field_y - dalisty[indy][0])2 ), self.field_x.shape )
1258	#TBD: pb with staggered coord
1259	if idx not in self.index_x: self.index_x.append(idx)
1260	if idy not in self.index_y: self.index_y.append(idy)
1261	self.index_x2d.append(idx)
1262	self.index_y2d.append(idy)
1263	## CASE 2 LON IS FIXED BUT NOT LAT
1264	elif test in ["fixedfree","fixedcomp"]:
1265	# find where are requested x values for each y on the free dimension
1266	# NB: this does not work for non-bijective cases e.g. polar stereographic
1267	for iy in range(self.dim_y):
1268	idx = np.argmin( np.abs( self.field_x[iy,:] - dalistx[indx][0] ) )
1269	# if comp is requested we select only indexes which yield values between requested min and max
1270	storeval = (self.method_y == "comp") and (self.field_y[iy,idx] > dalisty[indy][0]) and (self.field_y[iy,idx] < dalisty[indy][1])
1271	storeval = storeval or (self.method_y == "free")
1272	if storeval:
1273	if idx not in self.index_x: self.index_x.append(idx)
1274	if iy not in self.index_y and self.method_y == "comp": self.index_y.append(iy)
1275	if idx not in self.index_x2d or iy not in self.index_y2d:
1276	self.index_x2d.append(idx)
1277	self.index_y2d.append(iy)
1278	## CASE 3 LAT IS FIXED BUT NOT LON
1279	elif test in ["freefixed","compfixed"]:
1280	# find where are requested y values for each x on the free dimension
1281	# NB: this does not work for non-bijective cases e.g. polar stereographic
1282	for ix in range(self.dim_x):
1283	idy = np.argmin( np.abs( self.field_y[:,ix] - dalisty[indy][0] ) )
1284	# if comp is requested we select only indexes which yield values between requested min and max
1285	storeval = (self.method_x == "comp") and (self.field_x[idy,ix] > dalistx[indx][0]) and (self.field_x[idy,ix] < dalistx[indx][1])
1286	storeval = storeval or (self.method_x == "free")
1287	if storeval:
1288	if idy not in self.index_y: self.index_y.append(idy)
1289	if ix not in self.index_x and self.method_x == "comp": self.index_x.append(ix)
1290	if ix not in self.index_x2d or idy not in self.index_y2d:
1291	self.index_x2d.append(ix)
1292	self.index_y2d.append(idy)
1293	## check index tab
1294	if len(self.index_x) == 0 or len(self.index_y) == 0:
1295	print "!! ERROR !! no indices found. check prescribed latitudes or longitudes" ; exit()
1296	## ensure the array is a numpy array for getfield to work
1297	self.index_x = np.array(self.index_x)
1298	self.index_y = np.array(self.index_y)
1299	self.index_x2d = np.array(self.index_x2d)
1300	self.index_y2d = np.array(self.index_y2d)
1301	### print extrema
1302	printx = self.field_x[np.ix_(self.index_y2d, self.index_x2d)]
1303	printy = self.field_y[np.ix_(self.index_y2d, self.index_x2d)]
1304	if self.verbose:
1305	print "**** OK. x values (min,max).", printx.min(),printx.max()
1306	print "**** OK. y values (min,max).", printy.min(),printy.max()
1307
1308	# get the field from the NETCDF file and perform averages
1309	# -------------------------------
1310	def getfield(self):
1311	## first tell what is to be done
1312	if self.dimplot > 2: print "**** !! ERROR !! "+str(self.dimplot)+"D plots not supported!" ; exit()
1313	elif self.dimplot == 0 and self.verbose: print "**** OK. 0D value requested."
1314	elif self.dimplot == 1 and self.verbose: print "**** OK. 1D plot requested."
1315	elif self.verbose: print "**** OK. 2D section requested."
1316	# well, now get field from netcdf file
1317	# part below is necessary otherwise there is an index error below
1318	if self.index_x.size == 1: self.index_x = self.index_x[0]
1319	if self.index_y.size == 1: self.index_y = self.index_y[0]
1320	if self.index_z.size == 1: self.index_z = self.index_z[0]
1321	if self.index_t.size == 1: self.index_t = self.index_t[0]
1322	# then retrieve what is requested by user
1323	# each self.dim case corresponds to tests in the beginning of getdim.
1324	time0 = timelib.time()
1325	if self.verbose: print "**** OK. I am getting values from files. Please wait."
1326	if self.dim == 1:
1327	nt = self.index_t.size ; nz = 1 ; ny = 1 ; nx = 1
1328	self.field = self.f.variables[self.var][self.index_t]
1329	elif self.dim == 2:
1330	nt = 1 ; nz = 1 ; ny = self.index_y.size ; nx = self.index_x.size
1331	self.field = self.f.variables[self.var][self.index_y,self.index_x]
1332	elif self.dim == 3:
1333	nt = self.index_t.size ; nz = 1 ; ny = self.index_y.size ; nx = self.index_x.size
1334	self.field = self.f.variables[self.var][self.index_t,self.index_y,self.index_x]
1335	# this is far faster than retrieving each term with a loop
1336	elif self.dim == 4:
1337	nt = self.index_t.size ; nz = self.index_z.size ; ny = self.index_y.size ; nx = self.index_x.size
1338	self.field = self.f.variables[self.var][self.index_t,self.index_z,self.index_y,self.index_x]
1339	else:
1340	print "!! ERROR !! field would have more than four dimensions ?" ; exit()
1341	# NB: ... always 4D array but possibly with "size 1" dimensions
1342	# ... if one dimension is missing because 1D 2D or 3D requests, make it appear again
1343	self.field = np.reshape(self.field,(nt,nz,ny,nx))
1344	if self.verbose: print "**** OK. I got %7.1e values. This took me %6.4f seconds" % (nxnynz*nt,timelib.time() - time0)
1345	if self.verbose: print "**** OK. I got var "+self.var+" with shape",self.field.shape
1346	# reduce coordinates to useful points
1347	# ... TBD: this should be ordered in the case of non-regular projections
1348	if self.method_x in ["free","comp"] and self.method_y in ["free","comp"]:
1349	# we need 2D coordinates (free) or we get broadcast problem (comp) so we use np.ix
1350	self.field_x = self.field_x[np.ix_(self.index_y2d, self.index_x2d)]
1351	self.field_y = self.field_y[np.ix_(self.index_y2d, self.index_x2d)]
1352	else:
1353	# we are OK with 1D coordinates
1354	self.field_x = self.field_x[self.index_y2d, self.index_x2d]
1355	self.field_y = self.field_y[self.index_y2d, self.index_x2d]
1356	self.field_z = self.field_z[self.index_z]
1357	self.field_t = self.field_t[self.index_t]
1358	# now have to obtain the new indexes which correspond to the extracted self.field
1359	# for it to work with unique index, ensure that any index_* is a numpy array
1360	if not isinstance(self.index_x, np.ndarray): self.index_x = np.array([self.index_x])
1361	if not isinstance(self.index_y, np.ndarray): self.index_y = np.array([self.index_y])
1362	if not isinstance(self.index_z, np.ndarray): self.index_z = np.array([self.index_z])
1363	if not isinstance(self.index_t, np.ndarray): self.index_t = np.array([self.index_t])
1364	for val in self.index_x: self.index_x2d[np.where(self.index_x2d == val)] = np.where(self.index_x == val)[0]
1365	for val in self.index_y: self.index_y2d[np.where(self.index_y2d == val)] = np.where(self.index_y == val)[0]
1366	for val in self.index_z: self.index_z [np.where(self.index_z == val)] = np.where(self.index_z == val)[0]
1367	for val in self.index_t: self.index_t [np.where(self.index_t == val)] = np.where(self.index_t == val)[0]
1368	##### VERY EXPENSIVE
1369	## recast self.field with 2D horizontal arrays because we might have extracted
1370	## more than what is to be actually plot or computed, in particular for comps on 2D lat,lon coordinates
1371	#self.field = self.field[np.ix_(self.index_t,self.index_z,self.index_y2d,self.index_x2d)]
1372	#(nt,nz,ny,nx) = self.field.shape
1373	# extract relevant horizontal points
1374	# TBD: is compfree OK with computing on irregular grid?
1375	test = self.method_x + self.method_y
1376	if test in ["fixedfixed","freefree"]:
1377	pass
1378	elif test in ["fixedfree","fixedcomp"] or test in ["freefixed","compfixed"]:
1379	time0 = timelib.time()
1380	# prepare the loop on all relevant horizontal points
1381	if self.method_x in ["comp","free"]:
1382	nnn = self.index_x2d.shape[0] ; what_I_am_supposed_to_do = "keepx"
1383	elif self.method_y in ["comp","free"]:
1384	nnn = self.index_y2d.shape[0] ; what_I_am_supposed_to_do = "keepy"
1385	# LOOP to extract only useful values over horizontal dimensions
1386	# only take diagonal terms, do not loop on all self.index_x2d*self.index_y2d
1387	# ... this method is fast enough, perhaps there is a faster way though
1388	# ... (for sure the method with np.diag is much slower)
1389	for iii in range(nnn):
1390	ix = self.index_x2d[iii] ; iy = self.index_y2d[iii]
1391	for iz in self.index_z:
1392	for it in self.index_t:
1393	if what_I_am_supposed_to_do == "keepx": self.field[it,iz,0,ix] = self.field[it,iz,iy,ix]
1394	elif what_I_am_supposed_to_do == "keepy": self.field[it,iz,iy,0] = self.field[it,iz,iy,ix]
1395	if self.verbose: print "**** OK. I got to pick the right values for your request. This took me %6.4f seconds" % (timelib.time() - time0)
1396	# we only keep the one value that was modified on the dimension which is not free
1397	if what_I_am_supposed_to_do == "keepx": self.field = self.field[:,:,0,:] ; ny = 1 ; self.field = np.reshape(self.field,(nt,nz,ny,nx))
1398	elif what_I_am_supposed_to_do == "keepy": self.field = self.field[:,:,:,0] ; nx = 1 ; self.field = np.reshape(self.field,(nt,nz,ny,nx))
1399	# make a mask in case there are non-NaN missing values. (what about NaN missing values?)
1400	# ... this is important for computations below (see ppcompute)
1401	masked = np.ma.masked_where(np.abs(self.field) > 1e25,self.field)
1402	if masked.mask.any() == True:
1403	if self.verbose: print "!! WARNING !! Values over +-1e25 are considered missing values."
1404	self.field = masked
1405	self.field.set_fill_value([np.NaN])
1406
1407	# perform computations
1408	# -------------------------------
1409	# available: mean, max, min, meanarea
1410	# TB: integrals? for derivatives, define a function self.dx()
1411	def computations(self):
1412	nt,nz,ny,nx = self.field.shape
1413	# treat the case of mean on fields normalized with grid mesh area
1414	# ... this is done in the .area() method.
1415	# after that self.field contains field*area/totarea
1416	if "area" in self.compute:
1417	if "comp" in self.method_x+self.method_y:
1418	self.area()
1419	else:
1420	if self.verbose: print "!! WARNING !! No area accounted for (computing on t and/or z axis)."
1421	# now ready to compute [TBD: we would like to have e.g. mean over x,y and min over t ??]
1422	if self.method_t == "comp":
1423	if self.verbose: print "**** OK. Computing over t axis."
1424	if "mean" in self.compute: self.field = ppcompute.mean(self.field,axis=0)
1425	elif self.compute == "min": self.field = ppcompute.min(self.field,axis=0)
1426	elif self.compute == "max": self.field = ppcompute.max(self.field,axis=0)
1427	else: print "!! ERROR !! operation not supported." ; exit()
1428	nt = 1 ; self.field = np.reshape(self.field,(nt,nz,ny,nx))
1429	if self.method_z == "comp":
1430	if self.verbose: print "**** OK. Computing over z axis."
1431	if "mean" in self.compute: self.field = ppcompute.mean(self.field,axis=1)
1432	elif self.compute == "min": self.field = ppcompute.min(self.field,axis=1)
1433	elif self.compute == "max": self.field = ppcompute.max(self.field,axis=1)
1434	else: print "!! ERROR !! operation not supported." ; exit()
1435	nz = 1 ; self.field = np.reshape(self.field,(nt,nz,ny,nx))
1436	if self.method_y == "comp":
1437	if self.verbose: print "**** OK. Computing over y axis."
1438	if self.compute == "mean": self.field = ppcompute.mean(self.field,axis=2)
1439	elif self.compute == "min": self.field = ppcompute.min(self.field,axis=2)
1440	elif self.compute == "max": self.field = ppcompute.max(self.field,axis=2)
1441	elif self.compute == "meanarea": self.field = ppcompute.sum(self.field,axis=2)
1442	else: print "!! ERROR !! operation not supported." ; exit()
1443	ny = 1 ; self.field = np.reshape(self.field,(nt,nz,ny,nx))
1444	if self.field_x.ndim == 2: self.field_x = self.field_x[0,:] # TBD: this is OK for regular grid but not for irregular
1445	if self.method_x == "comp":
1446	if self.verbose: print "**** OK. Computing over x axis."
1447	if self.compute == "mean": self.field = ppcompute.mean(self.field,axis=3)
1448	elif self.compute == "min": self.field = ppcompute.min(self.field,axis=3)
1449	elif self.compute == "max": self.field = ppcompute.max(self.field,axis=3)
1450	elif self.compute == "meanarea": self.field = ppcompute.sum(self.field,axis=3)
1451	else: print "!! ERROR !! operation not supported." ; exit()
1452	nx = 1 ; self.field = np.reshape(self.field,(nt,nz,ny,nx))
1453	if self.field_y.ndim == 2: self.field_y = self.field_y[:,0] # TBD: this is OK for regular grid but not for irregular
1454	# remove all dimensions with size 1 to prepare plot (and check the resulting dimension with dimplot)
1455	self.field = np.squeeze(self.field)
1456	if self.field.ndim != self.dimplot:
1457	print "!! ERROR !! Problem: self.field is different than plot dimensions", self.field.ndim, self.dimplot ; exit()
1458	if self.verbose:
1459	print "**** OK. Final shape for "+self.var+" after averaging and squeezing",self.field.shape
1460
1461	# get areas for computations and ponderate self.field by area/totarea
1462	# -------------------------------------------------------------------
1463	def area(self):
1464	if self.verbose: print "**** OK. Get area array for computations."
1465	# create a request object for area
1466	# ... and copy known attributes from self
1467	aire = onerequest()
1468	aire.copy(self)
1469	# get area field name
1470	aire.var = "nothing"
1471	for c in glob_listarea:
1472	if c in aire.f.variables.keys():
1473	aire.var = c
1474	# do not try to calculate areas automatically
1475	if aire.var == "nothing":
1476	print "!! ERROR !! area variable not found... needs to be added in set_ppclass.txt?"
1477	exit()
1478	# define area request dimensions
1479	aire.getdim()
1480	# ensure this is a 2D horizontal request and define indexes
1481	# ... areas are not supposed to vary with time and height
1482	aire.method_x = "free" ; aire.method_y = "free"
1483	aire.getindexhori() ; aire.dimplot = 2
1484	aire.method_z = "fixed" ; aire.field_z = np.array([0]) ; aire.index_z = np.array([0])
1485	aire.method_t = "fixed" ; aire.field_t = np.array([0]) ; aire.index_t = np.array([0])
1486	# read the 2D area array in netCDF file
1487	aire.getfield()
1488	aire.field = np.squeeze(aire.field)
1489	# reduce with self horizontal indexes
1490	if "fixed" in self.method_x+self.method_y:
1491	aire.field = aire.field[self.index_y,self.index_x]
1492	# calculate total area
1493	# ... 2D comp is easy. 1D comp is a bit less easy but simple array manipulation.
1494	if "free" in self.method_x+self.method_y:
1495	if self.method_x == "free":
1496	totarea = ppcompute.sum(aire.field,axis=0)
1497	totarea = np.reshape(totarea,(1,totarea.size))
1498	totarea = np.tile(totarea,(1,self.index_x))
1499	elif self.method_y == "free":
1500	totarea = ppcompute.sum(aire.field,axis=1)
1501	totarea = np.reshape(totarea,(totarea.size,1))
1502	totarea = np.tile(totarea,(1,self.index_x.size))
1503	elif self.method_x == "comp" and self.method_y == "comp":
1504	totarea = ppcompute.sum(ppcompute.sum(aire.field,axis=1),axis=0)
1505	else:
1506	if self.verbose: print "!! WARNING !! Not account for areas. Only averaging over z and/or t axis."
1507	# normalize by total area
1508	print "**** OK. I can now normalize field by areas."
1509	aire.field = aire.field / totarea
1510	# tile area array over self t and z axis so that area field could be multiplied with self.field
1511	aire.field = np.tile(aire.field,(self.index_t.size,self.index_z.size,1,1))
1512	if self.field.shape != aire.field.shape:
1513	print "!! ERROR !! Problem in area(). Check array shapes." ; exit()
1514	else:
1515	self.field = self.field*aire.field
1516
1517	# define coordinates for plot
1518	# -------------------------------
1519	def definecoord(self):
1520	I_got_abs = False ; I_got_ord = False
1521	# here is the thing. time is usually taken as an abscissa so we start with time.
1522	if self.method_t == "free":
1523	self.absc = self.field_t ; self.absclab = self.name_t
1524	I_got_abs = True
1525	# then we usually have x as an abscissa.
1526	if self.method_x == "free":
1527	if I_got_abs:
1528	self.ordi = self.field_x ; self.ordilab = self.name_x
1529	I_got_ord = True
1530	else:
1531	self.absc = self.field_x ; self.absclab = self.name_x
1532	I_got_abs = True
1533	# ... or we have y
1534	if self.method_y == "free":
1535	if I_got_abs:
1536	self.ordi = self.field_y ; self.ordilab = self.name_y
1537	I_got_ord = True
1538	else:
1539	self.absc = self.field_y ; self.absclab = self.name_y
1540	I_got_abs = True
1541	# ... and we end with z because it is usually not an abscissa (profiles).
1542	if self.method_z == "free":
1543	if self.field_z[0] > self.field_z[1]:
1544	self.invert_axes = True # the axis will be turned upside-down
1545	if I_got_abs:
1546	self.ordi = self.field_z ; self.ordilab = self.name_z
1547	I_got_ord = True
1548	else:
1549	self.absc = self.field_z ; self.absclab = self.name_z
1550	I_got_abs = True
1551	self.swap_axes = True # says that altitude is not supposed to remain as an abscissa
1552	if I_got_abs and self.verbose: print "**** OK. abscissa:",self.absclab, self.absc.shape
1553	if I_got_ord and self.verbose: print "**** OK. ordinate:",self.ordilab, self.ordi.shape

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

Context Navigation

source: trunk/UTIL/PYTHON/planetoplot_v2/ppclass.py @ 923

Download in other formats: