Context Navigation

source: trunk/UTIL/PYTHON/planetoplot_v2/ppplot.py @ 1062

Last change on this file since 1062 was 1062, checked in by aslmd, 11 years ago
planetoplot_v2. added trans back showcb attributes to ppclass. simplified example scripts. added the web home page example.
File size: 35.8 KB

Rev	Line
[910]	1	###############################################
	2	## PLANETOPLOT ##
	3	## --> PPPLOT ##
	4	###############################################
	5	## Author: Aymeric Spiga. 02-03/2013 ##
	6	###############################################
	7	# python built-in librairies
	8	import time as timelib
	9	# added librairies
	10	import numpy as np
	11	import matplotlib.pyplot as mpl
	12	from matplotlib.cm import get_cmap
	13	from mpl_toolkits.basemap import Basemap
[936]	14	from matplotlib.ticker import FormatStrFormatter,MaxNLocator
[910]	15	# personal librairies
	16	import ppcompute
	17	###############################################
	18
	19	#################################
	20	# global variables and settings #
	21	#################################
	22
	23	# matplotlib settings
	24	# http://matplotlib.org/users/customizing.html
	25	# -------------------------------
	26	mpl.rcParams['font.family'] = "serif"
[947]	27	mpl.rcParams['axes.color_cycle'] = "b,r,g,k"
[928]	28	mpl.rcParams['contour.negative_linestyle'] = "dashed" # ou "solid"
[910]	29	mpl.rcParams['verbose.level'] = "silent"
	30	mpl.rcParams['lines.linewidth'] = 1.5
	31	mpl.rcParams['lines.markersize'] = 10
	32	mpl.rcParams['xtick.major.pad'] = 10
	33	mpl.rcParams['ytick.major.pad'] = 10
	34
	35	# global variables
	36	# -------------------------------
	37	# - where settings files are located
	38	# (None means planetoplot_v2 in PYTHONPATH)
	39	whereset = None
	40	# - some good default settings.
	41	# (bounds)
	42	how_many_sigma = 3.0
	43	# (contours)
	44	ccol = 'black'
	45	cline = 0.55
[1062]	46	cline = 0.8
[910]	47	# (vectors)
	48	widthvec = 0.002
	49	reducevec = 30.
	50	# (colorbar)
	51	zeorientation="vertical"
	52	zefrac = 0.05
	53	# (save figures)
	54	pad_inches_value=0.25
[928]	55	# (negative mode)
	56	def_negative = False
[910]	57	###############################################
	58
[928]	59	### settings for 'negative-like' mode
	60	if def_negative:
	61	mpl.rc('figure', facecolor='k', edgecolor='k')
	62	mpl.rcParams['text.color'] = 'w'
	63	mpl.rc('axes',labelcolor='w',facecolor='k',edgecolor='w')
	64	mpl.rcParams['xtick.color'] = 'w'
	65	mpl.rcParams['ytick.color'] = 'w'
	66	mpl.rcParams['grid.color'] = 'w'
	67	mpl.rc('savefig',facecolor='k',edgecolor='k')
	68
[910]	69	##########################
	70	# executed when imported #
	71	##########################
	72	###########################################
	73	# we load user-defined automatic settings #
	74	###########################################
	75	# initialize the warning variable about file not present...
	76	files_not_present = ""
	77	# ... and the whereset variable
	78	whereset = ppcompute.findset(whereset)
	79
	80	# - variable settings
	81	# -------------------------------
	82	zefile = "set_var.txt"
[933]	83	vf = {} ; vc = {} ; vl = {} ; vu = {}
[910]	84	try:
	85	f = open(whereset+zefile, 'r')
	86	for line in f:
	87	if "#" in line: pass
	88	else:
[933]	89	var, format, colorb, label, units = line.strip().split(';')
[910]	90	ind = var.strip()
	91	vf[ind] = format.strip()
	92	vc[ind] = colorb.strip()
	93	vl[ind] = label.strip()
[933]	94	vu[ind] = units.strip()
[910]	95	f.close()
	96	except IOError:
	97	files_not_present = files_not_present + zefile + " "
	98
	99	## - file settings
	100	## -------------------------------
	101	#zefile = "set_file.txt"
	102	#prefix_t = {} ; suffix_t = {} ; name_t = {} ; proj_t = {} ; vecx_t = {} ; vecy_t = {}
	103	#try:
	104	# f = open(whereset+zefile, 'r')
	105	# for line in f:
	106	# if "#" in line: pass
	107	# else:
	108	# prefix, suffix, name, proj, vecx, vecy = line.strip().split(';')
	109	# ind = name.strip()
	110	# prefix_t[ind] = prefix.strip()
	111	# suffix_t[ind] = suffix.strip()
	112	# #name_t[ind] = name.strip()
	113	# proj_t[ind] = proj.strip()
	114	# vecx_t[ind] = vecx.strip()
	115	# vecy_t[ind] = vecy.strip()
	116	# f.close()
	117	#except IOError:
	118	# files_not_present = files_not_present + zefile + " "
	119
	120	# - multiplot settings
	121	# -------------------------------
	122	zefile = "set_multiplot.txt"
	123	subv_t = {} ; subh_t = {} ; wspace_t = {} ; hspace_t = {} ; font_t = {}
	124	try:
	125	f = open(whereset+zefile, 'r')
	126	for line in f:
	127	if "#" in line: pass
	128	else:
	129	num, subv, subh, wspace, hspace, font = line.strip().split(';')
	130	ind = int(num.strip())
	131	subv_t[ind] = int(subv.strip())
	132	subh_t[ind] = int(subh.strip())
	133	wspace_t[ind] = float(wspace.strip())
	134	hspace_t[ind] = float(hspace.strip())
	135	font_t[ind] = float(font.strip())
	136	f.close()
	137	except IOError:
	138	files_not_present = files_not_present + zefile + " "
	139
	140	# - background settings
	141	# -------------------------------
	142	zefile = "set_back.txt"
	143	back = {}
	144	try:
	145	f = open(whereset+zefile, 'r')
	146	for line in f:
	147	if "#" in line: pass
	148	else:
	149	name, link = line.strip().split(';')
	150	ind = name.strip()
	151	back[ind] = link.strip()
	152	f.close()
	153	except IOError:
	154	files_not_present = files_not_present + zefile + " "
	155
	156	# - area settings
	157	# -------------------------------
	158	zefile = "set_area.txt"
	159	area = {}
	160	try:
	161	f = open(whereset+zefile, 'r')
	162	for line in f:
	163	if "#" in line: pass
	164	else:
	165	name, wlat1, wlat2, wlon1, wlon2 = line.strip().split(';')
	166	area[name.strip()] = [[float(wlon1.strip()),float(wlon2.strip())],\
	167	[float(wlat1.strip()),float(wlat2.strip())]]
	168	f.close()
	169	except IOError:
	170	files_not_present = files_not_present + zefile + " "
	171	# A note to the user about missing files
	172	if files_not_present != "":
	173	print "warning: files "+files_not_present+" not in "+whereset+" ; those presets will be missing"
	174
	175	# TBD: should change vector color with colormaps
	176	#"gist_heat": "white",\
	177	#"hot": "white",\
	178	#"gray": "red",\
	179
	180	####################
	181	# useful functions #
	182	####################
	183
[1007]	184	# continuity with matplotlib
	185	def close():
	186	mpl.close()
	187
	188	# a function for predefined figure sizes
	189	def figuref(x=16,y=9):
[1029]	190	fig = mpl.figure(figsize=(x,y))
	191	return fig
[1007]	192
[977]	193	# a function to change color lines according to a color map (idea by A. Pottier)
	194	# ... two optional arguments: color maps + a number telling how much intervals are needed
	195	def rainbow(cb="jet",num=8):
	196	ax = mpl.gca()
	197	pal = mpl.cm.get_cmap(name=cb)
	198	ax._get_lines.set_color_cycle([pal(i) for i in np.linspace(0,0.9,num)])
	199
[910]	200	# a function to define subplot
	201	# ... user can change settings in set_multiplot.txt read above
	202	# -------------------------------
	203	def definesubplot(numplot, fig):
	204	try:
	205	mpl.rcParams['font.size'] = font_t[numplot]
	206	except:
	207	mpl.rcParams['font.size'] = 18
	208	try:
	209	fig.subplots_adjust(wspace = wspace_t[numplot], hspace = hspace_t[numplot])
	210	subv, subh = subv_t[numplot],subh_t[numplot]
	211	except:
[1029]	212	print "!! WARNING !! (ignore if superpose mode) no preset found from set_multiplot.txt, or this setting file was not found."
[910]	213	subv = 1 ; subh = numplot
	214	return subv,subh
	215
	216	# a function to calculate automatically bounds (or simply prescribe those)
	217	# -------------------------------
	218	def calculate_bounds(field,vmin=None,vmax=None):
	219	# prescribed cases first
	220	zevmin = vmin
	221	zevmax = vmax
	222	# computed cases
	223	if zevmin is None or zevmax is None:
	224	# select values
	225	ind = np.where(field < 9e+35)
	226	fieldcalc = field[ ind ] # field must be a numpy array
	227	# calculate stdev and mean
	228	dev = np.std(fieldcalc)*how_many_sigma
	229	damean = ppcompute.mean(fieldcalc)
	230	# fill min/max if needed
	231	if vmin is None: zevmin = damean - dev
	232	if vmax is None: zevmax = damean + dev
	233	# special case: negative values with stddev while field is positive
[933]	234	if zevmin < 0. and ppcompute.min(fieldcalc) >= 0.: zevmin = 0.
[920]	235	# check that bounds are not too tight given the field
	236	amin = ppcompute.min(field)
	237	amax = ppcompute.max(field)
[936]	238	if np.abs(amin) < 1.e-15:
	239	cmin = 0.
	240	else:
	241	cmin = 100.*np.abs((amin - zevmin)/amin)
[920]	242	cmax = 100.*np.abs((amax - zevmax)/amax)
	243	if cmin > 150. or cmax > 150.:
	244	print "!! WARNING !! Bounds are a bit too tight. Might need to reconsider those."
	245	print "!! WARNING !! --> actual",amin,amax,"adopted",zevmin,zevmax
[910]	246	return zevmin, zevmax
	247	#### treat vmin = vmax for continuity
	248	#if vmin == vmax: zevmin = damean - dev ; zevmax = damean + dev
	249
	250	# a function to solve the problem with blank bounds !
	251	# -------------------------------
	252	def bounds(what_I_plot,zevmin,zevmax,miss=9e+35):
	253	small_enough = 1.e-7
	254	if zevmin < 0: what_I_plot[ what_I_plot < zevmin(1.-small_enough) ] = zevmin(1.-small_enough)
	255	else: what_I_plot[ what_I_plot < zevmin(1.+small_enough) ] = zevmin(1.+small_enough)
	256	what_I_plot[ what_I_plot > miss ] = -miss
	257	what_I_plot[ what_I_plot > zevmax ] = zevmax*(1.-small_enough)
	258	return what_I_plot
	259
[991]	260	# a function to change labels with modulo
	261	# ---------------------------------------
	262	def labelmodulo(ax,mod):
	263	mpl.draw()
	264	strtab = []
	265	for tick in ax.get_xaxis().get_ticklabels():
[1027]	266	onetick = tick.get_text()
	267	if len(onetick) > 0:
	268	num = float(onetick)
	269	strtab.append(num % mod)
[991]	270	ax.get_xaxis().set_ticklabels(strtab)
	271	return ax
	272
[990]	273	# a function to output an ascii file
	274	# ----------------------------------
	275	def writeascii (field=None,absc=None,name=None):
[1002]	276	field = np.array(field)
	277	absc = np.array(absc)
[990]	278	if name is None:
	279	name = "prof"
	280	for ttt in timelib.gmtime():
	281	name = name + "_" + str(ttt)
	282	if field is None or absc is None:
	283	print "!! WARNING !! Not printing the file, incorrect field or absc."
	284	else:
	285	if field.ndim == 1:
	286	myfile = open(name, 'w')
	287	for ix in range(len(absc)):
	288	myfile.write("%15.5e%15.5e\n" % ( absc[ix], field[ix] ))
	289	myfile.close()
	290	else:
	291	print "!! WARNING !! Not printing the file, 2D fields not supported yet."
	292	return
	293
[910]	294	# a generic function to show (GUI) or save a plot (PNG,EPS,PDF,...)
	295	# -------------------------------
[917]	296	def save(mode="gui",filename="plot",folder="./",includedate=True,res=150,custom=False):
[933]	297	if mode != "nothing":
	298	# a few settings
	299	possiblesave = ['eps','ps','svg','png','jpg','pdf']
	300	# now the main instructions
	301	if mode == "gui":
	302	mpl.show()
	303	elif mode in possiblesave:
	304	## name of plot
	305	name = folder+'/'+filename
	306	if includedate:
	307	for ttt in timelib.gmtime():
	308	name = name + "_" + str(ttt)
	309	name = name +"."+mode
	310	## save file
	311	print "**** Saving file in "+mode+" format... Please wait."
	312	if not custom:
	313	# ... regular plots
	314	mpl.savefig(name,dpi=res,pad_inches=pad_inches_value,bbox_inches='tight')
	315	else:
	316	# ... mapping mode, adapted space for labels etc...
	317	mpl.savefig(name,dpi=res)
	318	else:
	319	print "!! ERROR !! File format not supported. Supported: ",possiblesave
[910]	320
[1050]	321	# a necessary addition for people used to matplotlib
	322	def show():
	323	mpl.show()
	324
[910]	325	##################################
	326	# a generic class to make a plot #
	327	##################################
	328	class plot():
	329
	330	# print out a help string when help is invoked on the object
	331	# -------------------------------
	332	def __repr__(self):
	333	whatprint = 'plot object. \"help(plot)\" for more information\n'
	334	return whatprint
	335
	336	# default settings
	337	# -- user can define settings by two methods.
	338	# -- 1. yeah = plot2d(title="foo")
	339	# -- 2. yeah = pp() ; yeah.title = "foo"
	340	# -------------------------------
	341	def __init__(self,\
	342	var=None,\
[1050]	343	f=None,\
	344	x=None,\
[910]	345	xlabel="",\
	346	ylabel="",\
	347	div=20,\
	348	logx=False,\
	349	logy=False,\
	350	swap=False,\
	351	swaplab=True,\
	352	invert=False,\
	353	xcoeff=1.,\
	354	ycoeff=1.,\
[933]	355	fmt=None,\
[1050]	356	colorbar="jet",\
[933]	357	units="",\
[991]	358	modx=None,\
[1004]	359	xmin=None,\
	360	ymin=None,\
	361	xmax=None,\
	362	ymax=None,\
[1029]	363	nxticks=10,\
	364	nyticks=10,\
[910]	365	title=""):
	366	## what could be defined by the user
	367	self.var = var
[1050]	368	self.f = f
	369	self.x = x
[910]	370	self.xlabel = xlabel
	371	self.ylabel = ylabel
	372	self.title = title
	373	self.div = div
	374	self.logx = logx
	375	self.logy = logy
	376	self.swap = swap
	377	self.swaplab = swaplab # NB: swaplab only used if swap=True
	378	self.invert = invert
	379	self.xcoeff = xcoeff
	380	self.ycoeff = ycoeff
[933]	381	self.fmt = fmt
	382	self.units = units
[1050]	383	self.colorbar = colorbar
[991]	384	self.modx = modx
[1004]	385	self.xmin = xmin
	386	self.ymin = ymin
	387	self.xmax = xmax
	388	self.ymax = ymax
[1029]	389	self.nxticks = nxticks
	390	self.nyticks = nyticks
[910]	391	## other useful arguments
	392	## ... not used here in ppplot but need to be attached to plot object
	393	self.axisbg = "white"
	394	self.superpose = False
	395
	396	# check
	397	# -------------------------------
	398	def check(self):
[1050]	399	if self.f is None: print "!! ERROR !! Please define a field to be plotted" ; exit()
	400	else: self.f = np.array(self.f) # ensure this is a numpy array
[910]	401
	402	# define_from_var
	403	# ... this uses settings in set_var.txt
	404	# -------------------------------
	405	def define_from_var(self):
	406	if self.var is not None:
	407	if self.var.upper() in vl.keys():
	408	self.title = vl[self.var.upper()]
[933]	409	self.units = vu[self.var.upper()]
[910]	410
	411	# make
	412	# this is generic to all plots
	413	# -------------------------------
	414	def make(self):
	415	self.check()
[928]	416	# labels, title, etc...
[910]	417	mpl.xlabel(self.xlabel)
	418	mpl.ylabel(self.ylabel)
	419	if self.swap:
	420	if self.swaplab:
	421	mpl.xlabel(self.ylabel)
	422	mpl.ylabel(self.xlabel)
	423	mpl.title(self.title)
	424	# if masked array, set masked values to filled values (e.g. np.nan) for plotting purposes
[1050]	425	if type(self.f).__name__ in 'MaskedArray':
	426	self.f[self.f.mask] = self.f.fill_value
[910]	427
	428	################################
	429	# a subclass to make a 1D plot #
	430	################################
	431	class plot1d(plot):
	432
	433	# print out a help string when help is invoked on the object
	434	# -------------------------------
	435	def __repr__(self):
	436	whatprint = 'plot1d object. \"help(plot1d)\" for more information\n'
	437	return whatprint
	438
	439	# default settings
	440	# -- user can define settings by two methods.
	441	# -- 1. yeah = plot1d(title="foo")
	442	# -- 2. yeah = pp() ; yeah.title = "foo"
	443	# -------------------------------
	444	def __init__(self,\
[1050]	445	linestyle=None,\
[947]	446	color=None,\
[923]	447	marker='x',\
[1050]	448	legend=None):
[910]	449	## get initialization from parent class
	450	plot.__init__(self)
	451	## what could be defined by the user
[1050]	452	self.linestyle = linestyle
[910]	453	self.color = color
	454	self.marker = marker
[1050]	455	self.legend = legend
[910]	456
	457	# define_from_var
	458	# ... this uses settings in set_var.txt
	459	# -------------------------------
	460	def define_from_var(self):
	461	# get what is done in the parent class
	462	plot.define_from_var(self)
	463	# add specific stuff
	464	if self.var is not None:
	465	if self.var.upper() in vl.keys():
[943]	466	self.ylabel = vl[self.var.upper()] + " (" + vu[self.var.upper()] + ")"
[910]	467	self.title = ""
[933]	468	self.fmt = vf[self.var.upper()]
[910]	469
	470	# make
	471	# -------------------------------
	472	def make(self):
	473	# get what is done in the parent class
	474	plot.make(self)
[933]	475	if self.fmt is None: self.fmt = '%.0f'
[910]	476	# add specific stuff
	477	mpl.grid(color='grey')
[1050]	478	if self.linestyle == "": self.linestyle = " " # to allow for no line at all with ""
[977]	479	# set dummy x axis if not defined
[1050]	480	if self.x is None:
	481	self.x = np.array(range(self.f.shape[0]))
[977]	482	print "!! WARNING !! dummy coordinates on x axis"
[1002]	483	else:
[1050]	484	self.x = np.array(self.x) # ensure this is a numpy array
[910]	485	# swapping if requested
[1050]	486	if self.swap: x = self.f ; y = self.x
	487	else: x = self.x ; y = self.f
[910]	488	# coefficients on axis
	489	x=xself.xcoeff ; y=yself.ycoeff
	490	# check axis
	491	if x.size != y.size:
	492	print "!! ERROR !! x and y sizes don't match on 1D plot.", x.size, y.size
	493	exit()
	494	# make the 1D plot
	495	# either request linestyle or let matplotlib decide
[1050]	496	if self.linestyle is not None and self.color is not None:
	497	mpl.plot(x,y,self.color+self.linestyle,marker=self.marker,label=self.legend)
[942]	498	elif self.color is not None:
[1050]	499	mpl.plot(x,y,color=self.color,marker=self.marker,label=self.legend)
	500	elif self.linestyle is not None:
	501	mpl.plot(x,y,linestyle=self.linestyle,marker=self.marker,label=self.legend)
[910]	502	else:
[1050]	503	mpl.plot(x,y,marker=self.marker,label=self.legend)
[1027]	504	# AXES
	505	ax = mpl.gca()
[910]	506	# make log axes and/or invert ordinate
	507	# ... this must be after plot so that axis bounds are well-defined
	508	# ... also inverting must be after making the thing logarithmic
[977]	509	if self.logx: mpl.xscale("log") # not mpl.semilogx() because excludes log on y
	510	if self.logy: mpl.yscale("log") # not mpl.semilogy() because excludes log on x
[1027]	511	if self.invert: ax.set_ylim(ax.get_ylim()[::-1])
	512	if self.xmin is not None and self.xmax is not None:
	513	if self.xmin > self.xmax:
	514	ax.set_xlim(ax.get_xlim()[::-1])
	515	self.xmin,self.xmax = self.xmax,self.xmin
[923]	516	# add a label for line(s)
[1050]	517	if self.legend is not None:
	518	if self.legend != "":
[923]	519	mpl.legend(loc="best",fancybox=True)
[933]	520	# format labels
[936]	521	if self.swap: ax.xaxis.set_major_formatter(FormatStrFormatter(self.fmt))
	522	else: ax.yaxis.set_major_formatter(FormatStrFormatter(self.fmt))
[933]	523	# plot limits: ensure that no curve would be outside the window
	524	# x-axis
	525	x1, x2 = ax.get_xbound()
	526	xmin = ppcompute.min(x)
	527	xmax = ppcompute.max(x)
	528	if xmin < x1: x1 = xmin
	529	if xmax > x2: x2 = xmax
[977]	530	if self.xmin is not None: x1 = self.xmin
	531	if self.xmax is not None: x2 = self.xmax
[933]	532	ax.set_xbound(lower=x1,upper=x2)
	533	# y-axis
	534	y1, y2 = ax.get_ybound()
	535	ymin = ppcompute.min(y)
	536	ymax = ppcompute.max(y)
	537	if ymin < y1: y1 = ymin
	538	if ymax > y2: y2 = ymax
[977]	539	if self.ymin is not None: y1 = self.ymin
	540	if self.ymax is not None: y2 = self.ymax
[933]	541	ax.set_ybound(lower=y1,upper=y2)
[936]	542	## set with .div the number of ticks. (is it better than automatic?)
[981]	543	if not self.logx:
[1029]	544	ax.xaxis.set_major_locator(MaxNLocator(self.div/2)) #TBD: with nxticks?
[981]	545	else:
	546	print "!! WARNING. in logx mode, ticks are set automatically."
[991]	547	## specific modulo labels
	548	if self.modx is not None:
	549	ax = labelmodulo(ax,self.modx)
[910]	550
[1050]	551	# makeshow = make + show
	552	# -------------------------------
	553	def makeshow(self):
	554	self.make()
	555	mpl.show()
	556
[910]	557	################################
	558	# a subclass to make a 2D plot #
	559	################################
	560	class plot2d(plot):
	561
	562	# print out a help string when help is invoked on the object
	563	# -------------------------------
	564	def __repr__(self):
	565	whatprint = 'plot2d object. \"help(plot2d)\" for more information\n'
	566	return whatprint
	567
	568	# default settings
	569	# -- user can define settings by two methods.
	570	# -- 1. yeah = plot2d(title="foo")
	571	# -- 2. yeah = pp() ; yeah.title = "foo"
	572	# -------------------------------
	573	def __init__(self,\
[1050]	574	y=None,\
[910]	575	mapmode=False,\
	576	proj="cyl",\
	577	back=None,\
	578	trans=1.0,\
[1050]	579	vx=None,\
	580	vy=None,\
	581	c=None,\
[910]	582	blon=None,\
	583	blat=None,\
	584	area=None,\
	585	vmin=None,\
	586	vmax=None,\
[964]	587	showcb=True,\
[972]	588	wscale=None,\
[1050]	589	svx=1,\
	590	svy=1,\
[1040]	591	leftcorrect=False,\
[910]	592	colorvec="black"):
	593	## get initialization from parent class
	594	plot.__init__(self)
	595	## what could be defined by the user
[1050]	596	self.y = y
[910]	597	self.mapmode = mapmode
	598	self.proj = proj
	599	self.back = back
	600	self.trans = trans
[1050]	601	self.vx = vx
	602	self.vy = vy
[910]	603	self.colorvec = colorvec
[1050]	604	self.c = c
[910]	605	self.blon = blon ; self.blat = blat
	606	self.area = area
	607	self.vmin = vmin ; self.vmax = vmax
[964]	608	self.showcb = showcb
[972]	609	self.wscale = wscale
[1050]	610	self.svx = svx
	611	self.svy = svy
[1040]	612	self.leftcorrect = leftcorrect
[910]	613
	614	# define_from_var
	615	# ... this uses settings in set_var.txt
	616	# -------------------------------
	617	def define_from_var(self):
	618	# get what is done in the parent class
	619	plot.define_from_var(self)
	620	# add specific stuff
	621	if self.var is not None:
	622	if self.var.upper() in vl.keys():
[1050]	623	self.colorbar = vc[self.var.upper()]
[910]	624	self.fmt = vf[self.var.upper()]
	625
	626	# make
	627	# -------------------------------
	628	def make(self):
	629	# get what is done in the parent class...
	630	plot.make(self)
[977]	631	if self.fmt is None: self.fmt = "%.1e"
[910]	632	# ... then add specific stuff
	633	############################################################################################
	634	### PRE-SETTINGS
	635	############################################################################################
[960]	636	# set dummy xy axis if not defined
[1050]	637	if self.x is None:
	638	self.x = np.array(range(self.f.shape[0]))
[960]	639	self.mapmode = False
	640	print "!! WARNING !! dummy coordinates on x axis"
[1050]	641	if self.y is None:
	642	self.y = np.array(range(self.f.shape[1]))
[960]	643	self.mapmode = False
	644	print "!! WARNING !! dummy coordinates on y axis"
[977]	645	# transposing if necessary
[1050]	646	shape = self.f.shape
[977]	647	if shape[0] != shape[1]:
[1050]	648	if len(self.x) == shape[0] and len(self.y) == shape[1]:
[977]	649	print "!! WARNING !! Transposing axes"
[1050]	650	self.f = np.transpose(self.f)
[910]	651	# bound field
[1050]	652	zevmin, zevmax = calculate_bounds(self.f,vmin=self.vmin,vmax=self.vmax)
	653	what_I_plot = bounds(self.f,zevmin,zevmax)
[910]	654	# define contour field levels. define color palette
	655	ticks = self.div + 1
	656	zelevels = np.linspace(zevmin,zevmax,ticks)
[1050]	657	palette = get_cmap(name=self.colorbar)
[910]	658	# do the same thing for possible contourline entries
[1050]	659	if self.c is not None:
[910]	660	# if masked array, set masked values to filled values (e.g. np.nan) for plotting purposes
[1050]	661	if type(self.c).__name__ in 'MaskedArray':
	662	self.c[self.c.mask] = self.c.fill_value
	663	zevminc, zevmaxc = calculate_bounds(self.c)
	664	what_I_contour = bounds(self.c,zevminc,zevmaxc)
[910]	665	ticks = self.div + 1
	666	zelevelsc = np.linspace(zevminc,zevmaxc,ticks)
	667	############################################################################################
	668	### MAIN PLOT
	669	### NB: contour lines are done before contour shades otherwise colorar error
	670	############################################################################################
	671	if not self.mapmode:
	672	## A SIMPLE 2D PLOT
	673	###################
	674	# swapping if requested
[1050]	675	if self.swap: x = self.y ; y = self.x
	676	else: x = self.x ; y = self.y
[910]	677	# coefficients on axis
	678	x=xself.xcoeff ; y=yself.ycoeff
	679	# make shaded and line contours
[1050]	680	if self.c is not None:
[943]	681	objC = mpl.contour(x, y, what_I_contour, \
[910]	682	zelevelsc, colors = ccol, linewidths = cline)
[1007]	683	#mpl.clabel(objC, inline=1, fontsize="small",\
	684	# inline_spacing=1,fmt="%.0f")
[910]	685	mpl.contourf(x, y, \
[1050]	686	self.f, \
[910]	687	zelevels, cmap=palette)
[1007]	688	#mpl.pcolor(x,y,\
[1050]	689	# self.f, \
[1007]	690	# cmap=palette)
[910]	691	# make log axes and/or invert ordinate
[1004]	692	ax = mpl.gca()
[910]	693	if self.logx: mpl.semilogx()
	694	if self.logy: mpl.semilogy()
[1004]	695	if self.invert: ax.set_ylim(ax.get_ylim()[::-1])
[1027]	696	if self.xmin is not None and self.xmax is not None:
	697	if self.xmin > self.xmax:
	698	ax.set_xlim(ax.get_xlim()[::-1])
	699	self.xmin,self.xmax = self.xmax,self.xmin
[1004]	700	if self.xmin is not None: ax.set_xbound(lower=self.xmin)
	701	if self.xmax is not None: ax.set_xbound(upper=self.xmax)
	702	if self.ymin is not None: ax.set_ybound(lower=self.ymin)
	703	if self.ymax is not None: ax.set_ybound(upper=self.ymax)
[1029]	704	# set the number of ticks
	705	if not self.logx:
	706	ax.xaxis.set_major_locator(MaxNLocator(self.nxticks))
	707	else:
	708	print "!! WARNING. in logx mode, ticks are set automatically."
[1027]	709	## specific modulo labels
	710	if self.modx is not None:
	711	ax = labelmodulo(ax,self.modx)
[910]	712	else:
	713	## A 2D MAP USING PROJECTIONS (basemap)
	714	#######################################
	715	mpl.xlabel("") ; mpl.ylabel("")
	716	# additional security in case self.proj is None here
	717	# ... we set cylindrical projection (the simplest one)
	718	if self.proj is None: self.proj = "cyl"
	719	# get lon and lat in 2D version.
	720	# (but first ensure we do have 2D coordinates)
[1050]	721	if self.x.ndim == 1: [self.x,self.y] = np.meshgrid(self.x,self.y)
	722	elif self.x.ndim > 2: print "!! ERROR !! lon and lat arrays must be 1D or 2D"
[910]	723	# get lat lon intervals and associated settings
[1050]	724	wlon = [np.min(self.x),np.max(self.x)]
	725	wlat = [np.min(self.y),np.max(self.y)]
[910]	726	# -- area presets are in set_area.txt
	727	if self.area is not None:
	728	if self.area in area.keys():
	729	wlon, wlat = area[self.area]
	730	# -- settings for meridians and parallels
[936]	731	steplon = int(abs(wlon[1]-wlon[0])/6.)
	732	steplat = int(abs(wlat[1]-wlat[0])/3.)
	733	#mertab = np.r_[wlon[0]:wlon[1]:steplon] ; merlab = [0,0,0,1]
	734	#partab = np.r_[wlat[0]:wlat[1]:steplat] ; parlab = [1,0,0,0]
[972]	735	mertab = np.r_[-180.:180.:steplon] ; merlab = [0,0,0,1] #-360:360.
[936]	736	partab = np.r_[-90.:90.:steplat] ; parlab = [1,0,0,0]
[910]	737	format = '%.1f'
	738	# -- center of domain and bounding lats
	739	lon_0 = 0.5*(wlon[0]+wlon[1])
	740	lat_0 = 0.5*(wlat[0]+wlat[1])
	741	# some tests, bug fixes, and good-looking settings
	742	# ... cyl is good for global and regional
	743	if self.proj == "cyl":
[917]	744	format = '%.0f'
[910]	745	# ... global projections
	746	elif self.proj in ["ortho","moll","robin"]:
	747	wlat[0] = None ; wlat[1] = None ; wlon[0] = None ; wlon[1] = None
	748	steplon = 30. ; steplat = 30.
[1062]	749	if self.proj in ["moll"]: steplon = 60.
	750	if self.proj in ["robin"]: steplon = 90.
[910]	751	mertab = np.r_[-360.:360.:steplon]
	752	partab = np.r_[-90.:90.:steplat]
	753	if self.proj == "ortho":
	754	merlab = [0,0,0,0] ; parlab = [0,0,0,0]
	755	# in ortho projection, blon and blat can be used to set map center
	756	if self.blon is not None: lon_0 = self.blon
	757	if self.blat is not None: lat_0 = self.blat
	758	elif self.proj == "moll":
	759	merlab = [0,0,0,0]
	760	format = '%.0f'
	761	# ... regional projections
	762	elif self.proj in ["lcc","laea","merc"]:
	763	if self.proj in ["lcc","laea"] and wlat[0] == -wlat[1]:
	764	print "!! ERROR !! with Lambert lat1 must be different than lat2" ; exit()
	765	if wlat[0] < -80. and wlat[1] > 80.:
	766	print "!! ERROR !! set an area (not global)" ; exit()
	767	format = '%.0f'
	768	elif self.proj in ["npstere","spstere"]:
	769	# in polar projections, blat gives the bounding lat
	770	# if not set, set something reasonable
	771	if self.blat is None: self.blat = 60.
	772	# help the user who forgets self.blat would better be negative in spstere
	773	# (this actually serves for the default setting just above)
	774	if self.proj == "spstere" and self.blat > 0: self.blat = -self.blat
[932]	775	# labels
	776	mertab = np.r_[-360.:360.:15.]
	777	partab = np.r_[-90.:90.:5.]
[910]	778	# ... unsupported projections
	779	else:
	780	print "!! ERROR !! unsupported projection. supported: "+\
	781	"cyl, npstere, spstere, ortho, moll, robin, lcc, laea, merc"
	782	# finally define projection
	783	m = Basemap(projection=self.proj,\
	784	lat_0=lat_0,lon_0=lon_0,\
	785	boundinglat=self.blat,\
	786	llcrnrlat=wlat[0],urcrnrlat=wlat[1],\
	787	llcrnrlon=wlon[0],urcrnrlon=wlon[1])
[1040]	788	# in some case need to translated to the left for colorbar + labels
	789	# TBD: break stuff. a better solution should be found.
	790	if self.leftcorrect:
[1039]	791	ax = mpl.gca()
	792	pos = ax.get_position().bounds
	793	newpos = [0.,pos[1],pos[2],pos[3]]
	794	ax.set_position(newpos)
[910]	795	# draw meridians and parallels
	796	ft = int(mpl.rcParams['font.size']*3./4.)
[932]	797	zelatmax = 85.
	798	m.drawmeridians(mertab,labels=merlab,color='grey',linewidth=0.75,fontsize=ft,fmt=format,latmax=zelatmax)
	799	m.drawparallels(partab,labels=parlab,color='grey',linewidth=0.75,fontsize=ft,fmt=format,latmax=zelatmax)
[910]	800	# define background (see set_back.txt)
	801	if self.back is not None:
	802	if self.back in back.keys():
	803	print "**** info: loading a background, please wait.",self.back
[1029]	804	if self.back not in ["coast","sea"]:
	805	try: m.warpimage(back[self.back],scale=0.75)
	806	except: print "!! ERROR !! no background image could be loaded. probably not connected to the internet?"
	807	elif self.back == "coast":
	808	m.drawcoastlines()
	809	elif self.back == "sea":
	810	m.drawlsmask(land_color='white',ocean_color='aqua')
[910]	811	else:
	812	print "!! ERROR !! requested background not defined. change name or fill in set_back.txt" ; exit()
	813	# define x and y given the projection
[1050]	814	x, y = m(self.x, self.y)
[910]	815	# contour field. first line contour then shaded contour.
[1050]	816	if self.c is not None:
[943]	817	objC2 = m.contour(x, y, what_I_contour, \
[910]	818	zelevelsc, colors = ccol, linewidths = cline)
[1054]	819	#mpl.clabel(objC2, inline=1, fontsize=10,manual=True,fmt='-%2.0f$^{\circ}$C',colors='r')
[1062]	820	m.contourf(x, y, what_I_plot, zelevels, cmap = palette, alpha = self.trans, antialiased=True)
[910]	821	############################################################################################
	822	### COLORBAR
	823	############################################################################################
[964]	824	if self.trans > 0. and self.showcb:
[910]	825	## draw colorbar. settings are different with projections. or if not mapmode.
[933]	826	if not self.mapmode: orientation=zeorientation ; frac = 0.075 ; pad = 0.03 ; lu = 0.5
[942]	827	elif self.proj in ['moll']: orientation="horizontal" ; frac = 0.08 ; pad = 0.03 ; lu = 1.0
[933]	828	elif self.proj in ['cyl']: orientation="vertical" ; frac = 0.023 ; pad = 0.03 ; lu = 0.5
	829	else: orientation = zeorientation ; frac = zefrac ; pad = 0.03 ; lu = 0.5
[910]	830	zelevpal = np.linspace(zevmin,zevmax,num=min([ticks/2+1,21]))
[917]	831	zecb = mpl.colorbar(fraction=frac,pad=pad,\
[910]	832	format=self.fmt,orientation=orientation,\
	833	ticks=zelevpal,\
[933]	834	extend='neither',spacing='proportional')
[910]	835	if zeorientation == "horizontal": zecb.ax.set_xlabel(self.title) ; self.title = ""
[933]	836	# colorbar title --> units
	837	if self.units not in ["dimless",""]:
	838	zecb.ax.set_title("["+self.units+"]",fontsize=3.*mpl.rcParams['font.size']/4.,x=lu,y=1.025)
	839
[910]	840	############################################################################################
	841	### VECTORS. must be after the colorbar. we could also leave possibility for streamlines.
	842	############################################################################################
[1050]	843	### not expecting NaN in self.vx and self.vy. masked arrays is just enough.
	844	if self.vx is not None and self.vy is not None:
[982]	845	# vectors on map projection or simple 2D mapping
[1002]	846	if self.mapmode:
[1050]	847	[vecx,vecy] = m.rotate_vector(self.vx,self.vy,self.x,self.y) # for metwinds only ?
[1002]	848	else:
[1050]	849	vecx,vecy = self.vx,self.vy
[1002]	850	if x.ndim < 2 and y.ndim < 2: x,y = np.meshgrid(x,y)
[910]	851	# reference vector is scaled
[1002]	852	if self.wscale is None:
[1050]	853	self.wscale = ppcompute.mean(np.sqrt(self.vxself.vx+self.vyself.vy))
[910]	854	# make vector field
[981]	855	if self.mapmode:
[1050]	856	q = m.quiver( x[::self.svy,::self.svx],y[::self.svy,::self.svx],\
	857	vecx[::self.svy,::self.svx],vecy[::self.svy,::self.svx],\
[981]	858	angles='xy',color=self.colorvec,pivot='tail',\
	859	scale=self.wscale*reducevec,width=widthvec )
	860	else:
[1050]	861	q = mpl.quiver( x[::self.svy,::self.svx],y[::self.svy,::self.svx],\
	862	vecx[::self.svy,::self.svx],vecy[::self.svy,::self.svx],\
[981]	863	angles='xy',color=self.colorvec,pivot='tail',\
	864	scale=self.wscale*reducevec,width=widthvec )
[930]	865	# make vector key.
[933]	866	#keyh = 1.025 ; keyv = 1.05 # upper right corner over colorbar
[972]	867	keyh = 0.97 ; keyv = 1.05
[930]	868	#keyh = -0.03 ; keyv = 1.08 # upper left corner
[923]	869	p = mpl.quiverkey(q,keyh,keyv,\
[972]	870	self.wscale,str(int(self.wscale)),\
	871	fontproperties={'size': 'small'},\
[910]	872	color='black',labelpos='S',labelsep = 0.07)
[942]	873	############################################################################################
	874	### TEXT. ANYWHERE. add_text.txt should be present with lines x ; y ; text ; color
	875	############################################################################################
	876	try:
	877	f = open("add_text.txt", 'r')
	878	for line in f:
	879	if "#" in line: pass
	880	else:
	881	userx, usery, usert, userc = line.strip().split(';')
	882	userc = userc.strip()
	883	usert = usert.strip()
	884	userx = float(userx.strip())
	885	usery = float(usery.strip())
	886	if self.mapmode: userx,usery = m(userx,usery)
	887	mpl.text(userx,usery,usert,\
	888	color = userc,\
	889	horizontalalignment='center',\
	890	verticalalignment='center')
	891	f.close()
	892	except IOError:
	893	pass
[1050]	894
	895	# makeshow = make + show
	896	# -------------------------------
	897	def makeshow(self):
	898	self.make()
	899	mpl.show()

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

Download in other formats: