Context Navigation

source: trunk/GSASIIplot.py @ 986

Last change on this file since 986 was 986, checked in by vondreele, 10 years ago
remove getFont using freetype
Property svn:eol-style set to `native` Property svn:keywords set to `Date Author Revision URL Id`
File size: 159.2 KB

Line
1	# -- coding: utf-8 --
2	'''
3	GSASIIplot: plotting routines
4	===============================
5
6	'''
7	########### SVN repository information ###################
8	# $Date: 2013-07-10 09:28:58 +0000 (Wed, 10 Jul 2013) $
9	# $Author: vondreele $
10	# $Revision: 986 $
11	# $URL: trunk/GSASIIplot.py $
12	# $Id: GSASIIplot.py 986 2013-07-10 09:28:58Z vondreele $
13	########### SVN repository information ###################
14	import math
15	import time
16	import copy
17	import sys
18	import os.path
19	import numpy as np
20	import numpy.ma as ma
21	import numpy.linalg as nl
22	import wx
23	import wx.aui
24	import wx.glcanvas
25	import matplotlib as mpl
26	import mpl_toolkits.mplot3d.axes3d as mp3d
27	import GSASIIpath
28	GSASIIpath.SetVersionNumber("$Revision: 986 $")
29	import GSASIIgrid as G2gd
30	import GSASIIimage as G2img
31	import GSASIIpwd as G2pwd
32	import GSASIIIO as G2IO
33	import GSASIIpwdGUI as G2pdG
34	import GSASIIimgGUI as G2imG
35	import GSASIIphsGUI as G2phG
36	import GSASIIlattice as G2lat
37	import GSASIIspc as G2spc
38	import GSASIImath as G2mth
39	import pytexture as ptx
40	from OpenGL.GL import *
41	from OpenGL.GLU import *
42	from OpenGL.GLUT import *
43	from OpenGL.GLE import *
44	import gltext
45	from matplotlib.backends.backend_wx import _load_bitmap
46	from matplotlib.backends.backend_wxagg import FigureCanvasWxAgg as Canvas
47	from matplotlib.backends.backend_wxagg import NavigationToolbar2Wx as Toolbar
48
49	# useful degree trig functions
50	sind = lambda x: math.sin(x*math.pi/180.)
51	cosd = lambda x: math.cos(x*math.pi/180.)
52	tand = lambda x: math.tan(x*math.pi/180.)
53	asind = lambda x: 180.*math.asin(x)/math.pi
54	acosd = lambda x: 180.*math.acos(x)/math.pi
55	atan2d = lambda x,y: 180.*math.atan2(y,x)/math.pi
56	atand = lambda x: 180.*math.atan(x)/math.pi
57	# numpy versions
58	npsind = lambda x: np.sin(x*np.pi/180.)
59	npcosd = lambda x: np.cos(x*np.pi/180.)
60	nptand = lambda x: np.tan(x*np.pi/180.)
61	npacosd = lambda x: 180.*np.arccos(x)/np.pi
62	npasind = lambda x: 180.*np.arcsin(x)/np.pi
63	npatand = lambda x: 180.*np.arctan(x)/np.pi
64	npatan2d = lambda x,y: 180.*np.arctan2(x,y)/np.pi
65
66	class G2PlotMpl(wx.Panel):
67	'needs a doc string'
68	def __init__(self,parent,id=-1,dpi=None,**kwargs):
69	wx.Panel.__init__(self,parent,id=id,**kwargs)
70	mpl.rcParams['legend.fontsize'] = 10
71	self.figure = mpl.figure.Figure(dpi=dpi,figsize=(5,6))
72	self.canvas = Canvas(self,-1,self.figure)
73	self.toolbar = GSASIItoolbar(self.canvas)
74
75	self.toolbar.Realize()
76
77	sizer=wx.BoxSizer(wx.VERTICAL)
78	sizer.Add(self.canvas,1,wx.EXPAND)
79	sizer.Add(self.toolbar,0,wx.LEFT\|wx.EXPAND)
80	self.SetSizer(sizer)
81
82	class G2PlotOgl(wx.Panel):
83	'needs a doc string'
84	def __init__(self,parent,id=-1,dpi=None,**kwargs):
85	self.figure = wx.Panel.__init__(self,parent,id=id,**kwargs)
86	if 'win' in sys.platform: #Windows already double buffered
87	self.canvas = wx.glcanvas.GLCanvas(self,-1,**kwargs)
88	else: #fix from Jim Hester for X systems
89	attribs = (wx.glcanvas.WX_GL_DOUBLEBUFFER,)
90	self.canvas = wx.glcanvas.GLCanvas(self,-1,attribList=attribs,**kwargs)
91	self.camera = {}
92	sizer=wx.BoxSizer(wx.VERTICAL)
93	sizer.Add(self.canvas,1,wx.EXPAND)
94	self.SetSizer(sizer)
95
96	class G2Plot3D(wx.Panel):
97	'needs a doc string'
98	def __init__(self,parent,id=-1,dpi=None,**kwargs):
99	wx.Panel.__init__(self,parent,id=id,**kwargs)
100	self.figure = mpl.figure.Figure(dpi=dpi,figsize=(6,6))
101	self.canvas = Canvas(self,-1,self.figure)
102	self.toolbar = GSASIItoolbar(self.canvas)
103
104	self.toolbar.Realize()
105
106	sizer=wx.BoxSizer(wx.VERTICAL)
107	sizer.Add(self.canvas,1,wx.EXPAND)
108	sizer.Add(self.toolbar,0,wx.LEFT\|wx.EXPAND)
109	self.SetSizer(sizer)
110
111	class G2PlotNoteBook(wx.Panel):
112	'create a tabbed window for plotting'
113	def __init__(self,parent,id=-1):
114	wx.Panel.__init__(self,parent,id=id)
115	#so one can't delete a plot page!!
116	self.nb = wx.aui.AuiNotebook(self, \
117	style=wx.aui.AUI_NB_DEFAULT_STYLE ^ wx.aui.AUI_NB_CLOSE_ON_ACTIVE_TAB)
118	sizer = wx.BoxSizer()
119	sizer.Add(self.nb,1,wx.EXPAND)
120	self.SetSizer(sizer)
121	self.status = parent.CreateStatusBar()
122	self.status.SetFieldsCount(2)
123	self.status.SetStatusWidths([150,-1])
124	self.Bind(wx.aui.EVT_AUINOTEBOOK_PAGE_CHANGED, self.OnPageChanged)
125
126	self.plotList = []
127
128	def addMpl(self,name=""):
129	'Add a tabbed page with a matplotlib plot'
130	page = G2PlotMpl(self.nb)
131	self.nb.AddPage(page,name)
132
133	self.plotList.append(name)
134
135	return page.figure
136
137	def add3D(self,name=""):
138	'Add a tabbed page with a 3D plot'
139	page = G2Plot3D(self.nb)
140	self.nb.AddPage(page,name)
141
142	self.plotList.append(name)
143
144	return page.figure
145
146	def addOgl(self,name=""):
147	'Add a tabbed page with an openGL plot'
148	page = G2PlotOgl(self.nb)
149	self.nb.AddPage(page,name)
150
151	self.plotList.append(name)
152
153	return page.figure
154
155	def Delete(self,name):
156	'delete a tabbed page'
157	try:
158	item = self.plotList.index(name)
159	del self.plotList[item]
160	self.nb.DeletePage(item)
161	except ValueError: #no plot of this name - do nothing
162	return
163
164	def clear(self):
165	'clear all pages from plot window'
166	while self.nb.GetPageCount():
167	self.nb.DeletePage(0)
168	self.plotList = []
169	self.status.DestroyChildren()
170
171	def Rename(self,oldName,newName):
172	'rename a tab'
173	try:
174	item = self.plotList.index(oldName)
175	self.plotList[item] = newName
176	self.nb.SetPageText(item,newName)
177	except ValueError: #no plot of this name - do nothing
178	return
179
180	def OnPageChanged(self,event):
181	'respond to someone pressing a tab on the plot window'
182	if self.plotList:
183	self.status.SetStatusText('Better to select this from GSAS-II data tree',1)
184	self.status.DestroyChildren() #get rid of special stuff on status bar
185
186	class GSASIItoolbar(Toolbar):
187	'needs a doc string'
188	ON_MPL_HELP = wx.NewId()
189	ON_MPL_KEY = wx.NewId()
190	def __init__(self,plotCanvas):
191	Toolbar.__init__(self,plotCanvas)
192	POSITION_OF_CONFIGURE_SUBPLOTS_BTN = 6
193	self.DeleteToolByPos(POSITION_OF_CONFIGURE_SUBPLOTS_BTN)
194	parent = self.GetParent()
195	key = os.path.join(os.path.split(__file__)[0],'key.ico')
196	self.AddSimpleTool(self.ON_MPL_KEY,_load_bitmap(key),'Key press','Select key press')
197	wx.EVT_TOOL(self,self.ON_MPL_KEY,self.OnKey)
198	help = os.path.join(os.path.split(__file__)[0],'help.ico')
199	self.AddSimpleTool(self.ON_MPL_HELP,_load_bitmap(help),'Help on','Show help on')
200	wx.EVT_TOOL(self,self.ON_MPL_HELP,self.OnHelp)
201	def OnHelp(self,event):
202	'needs a doc string'
203	Page = self.GetParent().GetParent()
204	pageNo = Page.GetSelection()
205	bookmark = Page.GetPageText(pageNo)
206	bookmark = bookmark.strip(')').replace('(','_')
207	G2gd.ShowHelp(bookmark,self.TopLevelParent)
208	def OnKey(self,event):
209	'needs a doc string'
210	parent = self.GetParent()
211	if parent.Choice:
212	dlg = wx.SingleChoiceDialog(parent,'Select','Key press',list(parent.Choice))
213	if dlg.ShowModal() == wx.ID_OK:
214	sel = dlg.GetSelection()
215	event.key = parent.Choice[sel][0]
216	parent.keyPress(event)
217	dlg.Destroy()
218
219	################################################################################
220	##### PlotSngl
221	################################################################################
222
223	def PlotSngl(self,newPlot=False):
224	'''Single crystal structure factor plotting package - displays zone of reflections as rings proportional
225	to F, F**2, etc. as requested
226	'''
227	from matplotlib.patches import Circle
228	global HKL,HKLF
229
230	def OnSCMotion(event):
231	xpos = event.xdata
232	if xpos:
233	xpos = round(xpos) #avoid out of frame mouse position
234	ypos = round(event.ydata)
235	zpos = Data['Layer']
236	if '100' in Data['Zone']:
237	HKLtxt = '(%3d,%3d,%3d)'%(zpos,xpos,ypos)
238	elif '010' in Data['Zone']:
239	HKLtxt = '(%3d,%3d,%3d)'%(xpos,zpos,ypos)
240	elif '001' in Data['Zone']:
241	HKLtxt = '(%3d,%3d,%3d)'%(xpos,ypos,zpos)
242	Page.canvas.SetToolTipString(HKLtxt)
243	self.G2plotNB.status.SetFields(['HKL = '+HKLtxt,''])
244
245	def OnSCPick(event):
246	zpos = Data['Layer']
247	pos = event.artist.center
248	if '100' in Data['Zone']:
249	Page.canvas.SetToolTipString('(picked:(%3d,%3d,%3d))'%(zpos,pos[0],pos[1]))
250	hkl = np.array([zpos,pos[0],pos[1]])
251	elif '010' in Data['Zone']:
252	Page.canvas.SetToolTipString('(picked:(%3d,%3d,%3d))'%(pos[0],zpos,pos[1]))
253	hkl = np.array([pos[0],zpos,pos[1]])
254	elif '001' in Data['Zone']:
255	Page.canvas.SetToolTipString('(picked:(%3d,%3d,%3d))'%(pos[0],pos[1],zpos))
256	hkl = np.array([pos[0],pos[1],zpos])
257	h,k,l = hkl
258	hklf = HKLF[np.where(np.all(HKL-hkl == [0,0,0],axis=1))]
259	if len(hklf):
260	Fosq,sig,Fcsq = hklf[0]
261	HKLtxt = '(%3d,%3d,%3d %.2f %.3f %.2f %.2f)'%(h,k,l,Fosq,sig,Fcsq,(Fosq-Fcsq)/(scale*sig))
262	self.G2plotNB.status.SetFields(['','HKL, Fosq, sig, Fcsq, delFsq/sig = '+HKLtxt])
263
264	try:
265	plotNum = self.G2plotNB.plotList.index('Structure Factors')
266	Page = self.G2plotNB.nb.GetPage(plotNum)
267	if not newPlot:
268	Plot = Page.figure.gca() #get previous powder plot & get limits
269	xylim = Plot.get_xlim(),Plot.get_ylim()
270	Page.figure.clf()
271	Page.Choice = None
272	Plot = Page.figure.gca() #get a fresh plot after clf()
273	except ValueError:
274	Plot = self.G2plotNB.addMpl('Structure Factors').gca()
275	plotNum = self.G2plotNB.plotList.index('Structure Factors')
276	Page = self.G2plotNB.nb.GetPage(plotNum)
277	Page.canvas.mpl_connect('pick_event', OnSCPick)
278	Page.canvas.mpl_connect('motion_notify_event', OnSCMotion)
279	Page.Choice = None
280	Page.SetFocus()
281
282	Plot.set_aspect(aspect='equal')
283	HKLref = self.PatternTree.GetItemPyData(self.Sngl)[1]
284	Data = self.PatternTree.GetItemPyData(
285	G2gd.GetPatternTreeItemId(self,self.Sngl, 'HKL Plot Controls'))
286	Type = Data['Type']
287	scale = Data['Scale']
288	HKLmax = Data['HKLmax']
289	HKLmin = Data['HKLmin']
290	FosqMax = Data['FoMax']
291	FoMax = math.sqrt(FosqMax)
292	xlabel = ['k, h=','h, k=','h, l=']
293	ylabel = ['l','l','k']
294	zones = ['100','010','001']
295	pzone = [[1,2],[0,2],[0,1]]
296	izone = zones.index(Data['Zone'])
297	Plot.set_title(self.PatternTree.GetItemText(self.Sngl)[5:])
298	HKL = []
299	HKLF = []
300	for refl in HKLref:
301	H = np.array(refl[:3])
302	Fosq,sig,Fcsq = refl[5:8]
303	HKL.append(H)
304	HKLF.append([Fosq,sig,Fcsq])
305	if H[izone] == Data['Layer']:
306	A = 0
307	B = 0
308	if Type == 'Fosq':
309	A = scale*Fosq/FosqMax
310	B = scale*Fcsq/FosqMax
311	C = abs(A-B)
312	elif Type == 'Fo':
313	A = scale*math.sqrt(max(0,Fosq))/FoMax
314	B = scale*math.sqrt(max(0,Fcsq))/FoMax
315	C = abs(A-B)
316	elif Type == '\|DFsq\|/sig':
317	if sig > 0.:
318	A = scale*(Fosq-Fcsq)/(sig)
319	B = 0
320	elif Type == '\|DFsq\|>sig':
321	if sig > 0.:
322	A = (Fosq-Fcsq)/sig
323	if abs(A) < 1.0: A = 0
324	A *= scale
325	B = 0
326	elif Type == '\|DFsq\|>3sig':
327	if sig > 0.:
328	A = (Fosq-Fcsq)/sig
329	if abs(A) < 3.0: A = 0
330	A *= scale
331	B = 0
332	xy = (H[pzone[izone][0]],H[pzone[izone][1]])
333	if Type in ['\|DFsq\|/sig','\|DFsq\|>sig','\|DFsq\|>3sig']:
334	if A > 0.0:
335	Plot.add_artist(Circle(xy,radius=A,ec='g',fc='w',picker=3))
336	else:
337	Plot.add_artist(Circle(xy,radius=-A,ec='r',fc='w',picker=3))
338	else:
339	if A > 0.0 and A > B:
340	Plot.add_artist(Circle(xy,radius=A,ec='g',fc='w',picker=3))
341	if B:
342	Plot.add_artist(Circle(xy,radius=B,ec='b',fc='w'))
343	if A < B:
344	Plot.add_artist(Circle(xy,radius=A,ec='g',fc='w',picker=3))
345	radius = C
346	if radius > 0:
347	if A > B:
348	Plot.add_artist(Circle(xy,radius=radius,ec='g',fc='g'))
349	else:
350	Plot.add_artist(Circle(xy,radius=radius,ec='r',fc='r'))
351	HKL = np.array(HKL,dtype=np.int)
352	HKLF = np.array(HKLF)
353	Plot.set_xlabel(xlabel[izone]+str(Data['Layer']),fontsize=12)
354	Plot.set_ylabel(ylabel[izone],fontsize=12)
355	Plot.set_xlim((HKLmin[pzone[izone][0]],HKLmax[pzone[izone][0]]))
356	Plot.set_ylim((HKLmin[pzone[izone][1]],HKLmax[pzone[izone][1]]))
357	Page.canvas.draw()
358	# if not newPlot:
359	# Page.toolbar.push_current()
360	# Plot.set_xlim(xylim[0])
361	# Plot.set_ylim(xylim[1])
362	# xylim = []
363	# Page.toolbar.push_current()
364	# Page.toolbar.draw()
365	# else:
366	# Page.canvas.draw()
367
368	################################################################################
369	##### PlotPatterns
370	################################################################################
371
372	def PlotPatterns(G2frame,newPlot=False):
373	'''Powder pattern plotting package - displays single or multiple powder patterns as intensity vs
374	2-theta, q or TOF. Can display multiple patterns as "waterfall plots" or contour plots. Log I
375	plotting available.
376	'''
377	global HKL
378
379	# def OnKeyBox(event):
380	# if G2frame.G2plotNB.nb.GetSelection() == G2frame.G2plotNB.plotList.index('Powder Patterns'):
381	# event.key = cb.GetValue()[0]
382	# cb.SetValue(' key press')
383	# wx.CallAfter(OnPlotKeyPress,event)
384	#
385	def OnPlotKeyPress(event):
386	newPlot = False
387	if event.key == 'w':
388	if G2frame.Weight:
389	G2frame.Weight = False
390	else:
391	G2frame.Weight = True
392	G2frame.SinglePlot = True
393	newPlot = True
394	elif event.key == 'b':
395	if G2frame.SubBack:
396	G2frame.SubBack = False
397	else:
398	G2frame.SubBack = True
399	G2frame.SinglePlot = True
400	elif event.key == 'n':
401	if G2frame.Contour:
402	pass
403	else:
404	if G2frame.logPlot:
405	G2frame.logPlot = False
406	else:
407	G2frame.Offset[0] = 0
408	G2frame.logPlot = True
409	newPlot = True
410	elif event.key == 'u':
411	if G2frame.Contour:
412	G2frame.Cmax = min(1.0,G2frame.Cmax*1.2)
413	elif G2frame.logPlot:
414	pass
415	elif G2frame.Offset[0] < 100.:
416	G2frame.Offset[0] += 1.
417	elif event.key == 'd':
418	if G2frame.Contour:
419	G2frame.Cmax = max(0.0,G2frame.Cmax*0.8)
420	elif G2frame.logPlot:
421	pass
422	elif G2frame.Offset[0] > 0.:
423	G2frame.Offset[0] -= 1.
424	elif event.key == 'l':
425	G2frame.Offset[1] -= 1.
426	elif event.key == 'r':
427	G2frame.Offset[1] += 1.
428	elif event.key == 'o':
429	G2frame.Offset = [0,0]
430	elif event.key == 'c':
431	newPlot = True
432	if G2frame.Contour:
433	G2frame.Contour = False
434	else:
435	G2frame.Contour = True
436	G2frame.SinglePlot = False
437	G2frame.Offset = [0.,0.]
438	elif event.key == 'q':
439	newPlot = True
440	if G2frame.qPlot:
441	G2frame.qPlot = False
442	else:
443	G2frame.qPlot = True
444	elif event.key == 's':
445	if G2frame.Contour:
446	choice = [m for m in mpl.cm.datad.keys() if not m.endswith("_r")]
447	choice.sort()
448	dlg = wx.SingleChoiceDialog(G2frame,'Select','Color scheme',choice)
449	if dlg.ShowModal() == wx.ID_OK:
450	sel = dlg.GetSelection()
451	G2frame.ContourColor = choice[sel]
452	else:
453	G2frame.ContourColor = 'Paired'
454	dlg.Destroy()
455	else:
456	if G2frame.SinglePlot:
457	G2frame.SinglePlot = False
458	else:
459	G2frame.SinglePlot = True
460	newPlot = True
461	elif event.key == '+':
462	if G2frame.PickId:
463	G2frame.PickId = False
464	elif event.key == 'i': #for smoothing contour plot
465	choice = ['nearest','bilinear','bicubic','spline16','spline36','hanning',
466	'hamming','hermite','kaiser','quadric','catrom','gaussian','bessel',
467	'mitchell','sinc','lanczos']
468	dlg = wx.SingleChoiceDialog(G2frame,'Select','Interpolation',choice)
469	if dlg.ShowModal() == wx.ID_OK:
470	sel = dlg.GetSelection()
471	G2frame.Interpolate = choice[sel]
472	else:
473	G2frame.Interpolate = 'nearest'
474	dlg.Destroy()
475
476	PlotPatterns(G2frame,newPlot=newPlot)
477
478	def OnMotion(event):
479	xpos = event.xdata
480	if xpos: #avoid out of frame mouse position
481	ypos = event.ydata
482	Page.canvas.SetCursor(wx.CROSS_CURSOR)
483	try:
484	Parms,Parms2 = G2frame.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(G2frame,G2frame.PatternId, 'Instrument Parameters'))
485	if 'C' in Parms['Type'][0]:
486	wave = G2mth.getWave(Parms)
487	if G2frame.qPlot:
488	try:
489	xpos = 2.0asind(xposwave/(4*math.pi))
490	except ValueError: #avoid bad value in asin beyond upper limit
491	pass
492	dsp = 0.0
493	if abs(xpos) > 0.: #avoid possible singularity at beam center
494	dsp = wave/(2.*sind(abs(xpos)/2.0))
495	if G2frame.Contour:
496	G2frame.G2plotNB.status.SetStatusText('2-theta =%9.3f d =%9.5f pattern ID =%5d'%(xpos,dsp,int(ypos)),1)
497	else:
498	G2frame.G2plotNB.status.SetStatusText('2-theta =%9.3f d =%9.5f Intensity =%9.1f'%(xpos,dsp,ypos),1)
499	else: #TOF neutrons
500	dsp = 0.0
501	difC = Parms['difC'][1]
502	dsp = xpos/difC #rough approx.!
503	if G2frame.Contour:
504	G2frame.G2plotNB.status.SetStatusText('TOF =%9.3f d =%9.5f pattern ID =%5d'%(xpos,dsp,int(ypos)),1)
505	else:
506	G2frame.G2plotNB.status.SetStatusText('TOF =%9.3f d =%9.5f Intensity =%9.1f'%(xpos,dsp,ypos),1)
507	if G2frame.itemPicked:
508	Page.canvas.SetToolTipString('%9.3f'%(xpos))
509	if G2frame.PickId:
510	found = []
511	if G2frame.PatternTree.GetItemText(G2frame.PickId) in ['Index Peak List','Unit Cells List','Reflection Lists'] or \
512	'PWDR' in G2frame.PatternTree.GetItemText(PickId):
513	if len(HKL):
514	view = Page.toolbar._views.forward()[0][:2]
515	wid = view[1]-view[0]
516	found = HKL[np.where(np.fabs(HKL.T[5]-xpos) < 0.002*wid)]
517	if len(found):
518	h,k,l = found[0][:3]
519	Page.canvas.SetToolTipString('%d,%d,%d'%(int(h),int(k),int(l)))
520	else:
521	Page.canvas.SetToolTipString('')
522
523	except TypeError:
524	G2frame.G2plotNB.status.SetStatusText('Select PWDR powder pattern first',1)
525
526	def OnPick(event):
527	if G2frame.itemPicked is not None: return
528	PatternId = G2frame.PatternId
529	try:
530	Parms,Parms2 = G2frame.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(G2frame,G2frame.PatternId, 'Instrument Parameters'))
531	except TypeError:
532	return
533	if 'C' in Parms['Type'][0]:
534	wave = G2mth.getWave(Parms)
535	else:
536	difC = Parms['difC'][1]
537	PickId = G2frame.PickId
538	pick = event.artist
539	mouse = event.mouseevent
540	xpos = pick.get_xdata()
541	ypos = pick.get_ydata()
542	ind = event.ind
543	xy = list(zip(np.take(xpos,ind),np.take(ypos,ind))[0])
544	if G2frame.PatternTree.GetItemText(PickId) == 'Peak List':
545	if ind.all() != [0]: #picked a data point
546	data = G2frame.PatternTree.GetItemPyData(G2frame.PickId)
547	XY = G2mth.setPeakparms(Parms,Parms2,xy[0],xy[1]) #what happens for a q-plot???
548	data.append(XY)
549	G2pdG.UpdatePeakGrid(G2frame,data)
550	PlotPatterns(G2frame)
551	else: #picked a peak list line
552	G2frame.itemPicked = pick
553	elif G2frame.PatternTree.GetItemText(PickId) == 'Limits':
554	if ind.all() != [0]: #picked a data point
555	LimitId = G2gd.GetPatternTreeItemId(G2frame,PatternId, 'Limits')
556	data = G2frame.PatternTree.GetItemPyData(LimitId)
557	if 'C' in Parms['Type'][0]: #CW data - TOF later in an elif
558	if G2frame.qPlot: #qplot - convert back to 2-theta
559	xy[0] = 2.0asind(xy[0]wave/(4*math.pi))
560	if mouse.button==1:
561	data[1][0] = min(xy[0],data[1][1])
562	if mouse.button==3:
563	data[1][1] = max(xy[0],data[1][0])
564	G2frame.PatternTree.SetItemPyData(LimitId,data)
565	G2pdG.UpdateLimitsGrid(G2frame,data)
566	PlotPatterns(G2frame)
567	else: #picked a limit line
568	G2frame.itemPicked = pick
569	elif G2frame.PatternTree.GetItemText(PickId) == 'Reflection Lists' or \
570	'PWDR' in G2frame.PatternTree.GetItemText(PickId):
571	G2frame.itemPicked = pick
572	pick = str(pick)
573
574	def OnRelease(event):
575	if G2frame.itemPicked is None: return
576	Parms,Parms2 = G2frame.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(G2frame,G2frame.PatternId, 'Instrument Parameters'))
577	if 'C' in Parms['Type'][0]:
578	wave = G2mth.getWave(Parms)
579	else:
580	difC = Parms['difC'][1]
581	xpos = event.xdata
582	PickId = G2frame.PickId
583	if G2frame.PatternTree.GetItemText(PickId) in ['Peak List','Limits'] and xpos:
584	lines = []
585	for line in G2frame.Lines:
586	lines.append(line.get_xdata()[0])
587	# print G2frame.itemPicked.get_xdata()
588	lineNo = lines.index(G2frame.itemPicked.get_xdata()[0])
589	if lineNo in [0,1]:
590	LimitId = G2gd.GetPatternTreeItemId(G2frame,G2frame.PatternId, 'Limits')
591	data = G2frame.PatternTree.GetItemPyData(LimitId)
592	if G2frame.qPlot:
593	if 'C' in Parms['Type'][0]:
594	data[1][lineNo] = 2.0asind(wavexpos/(4*math.pi))
595	else:
596	data[1][lineNo] = 2math.piParms['difC'][1]/xpos
597	else:
598	data[1][lineNo] = xpos
599	G2frame.PatternTree.SetItemPyData(LimitId,data)
600	if G2frame.PatternTree.GetItemText(G2frame.PickId) == 'Limits':
601	G2pdG.UpdateLimitsGrid(G2frame,data)
602	else:
603	PeakId = G2gd.GetPatternTreeItemId(G2frame,G2frame.PatternId, 'Peak List')
604	data = G2frame.PatternTree.GetItemPyData(PeakId)
605	# print 'peaks',xpos
606	if event.button == 3:
607	del data[lineNo-2]
608	else:
609	if G2frame.qPlot:
610	data[lineNo-2][0] = 2.0asind(wavexpos/(4*math.pi))
611	else:
612	data[lineNo-2][0] = xpos
613	G2frame.PatternTree.SetItemPyData(PeakId,data)
614	G2pdG.UpdatePeakGrid(G2frame,data)
615	elif (G2frame.PatternTree.GetItemText(PickId) == 'Reflection Lists' or \
616	'PWDR' in G2frame.PatternTree.GetItemText(PickId)) and xpos:
617	Phases = G2frame.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(G2frame,PatternId,'Reflection Lists'))
618	pick = str(G2frame.itemPicked).split('(')[1].strip(')')
619	if 'line' not in pick: #avoid data points, etc.
620	num = Phases.keys().index(pick)
621	if num:
622	G2frame.refDelt = -(event.ydata-G2frame.refOffset)/(num*Ymax)
623	else: #1st row of refl ticks
624	G2frame.refOffset = event.ydata
625	PlotPatterns(G2frame)
626	G2frame.itemPicked = None
627
628	xylim = []
629	try:
630	plotNum = G2frame.G2plotNB.plotList.index('Powder Patterns')
631	Page = G2frame.G2plotNB.nb.GetPage(plotNum)
632	if not newPlot:
633	Plot = Page.figure.gca() #get previous powder plot & get limits
634	xylim = Plot.get_xlim(),Plot.get_ylim()
635	Page.figure.clf()
636	Plot = Page.figure.gca() #get a fresh plot after clf()
637	except ValueError:
638	newPlot = True
639	G2frame.Cmax = 1.0
640	Plot = G2frame.G2plotNB.addMpl('Powder Patterns').gca()
641	plotNum = G2frame.G2plotNB.plotList.index('Powder Patterns')
642	Page = G2frame.G2plotNB.nb.GetPage(plotNum)
643	Page.canvas.mpl_connect('key_press_event', OnPlotKeyPress)
644	Page.canvas.mpl_connect('motion_notify_event', OnMotion)
645	Page.canvas.mpl_connect('pick_event', OnPick)
646	Page.canvas.mpl_connect('button_release_event', OnRelease)
647	Page.SetFocus()
648	G2frame.G2plotNB.status.DestroyChildren()
649	if G2frame.Contour:
650	Page.Choice = (' key press','d: lower contour max','u: raise contour max',
651	'i: interpolation method','s: color scheme','c: contour off')
652	else:
653	if G2frame.logPlot:
654	Page.Choice = (' key press','n: log(I) off','l: offset left','r: offset right',
655	'c: contour on','q: toggle q plot','s: toggle single plot','+: no selection')
656	else:
657	Page.Choice = (' key press','l: offset left','r: offset right','d: offset down',
658	'u: offset up','o: reset offset','b: toggle subtract background','n: log(I) on','c: contour on',
659	'q: toggle q plot','s: toggle single plot','w: toggle divide by sig','+: no selection')
660	Page.keyPress = OnPlotKeyPress
661
662	PickId = G2frame.PickId
663	PatternId = G2frame.PatternId
664	colors=['b','g','r','c','m','k']
665	Lines = []
666	if G2frame.SinglePlot:
667	Pattern = G2frame.PatternTree.GetItemPyData(PatternId)
668	Pattern.append(G2frame.PatternTree.GetItemText(PatternId))
669	PlotList = [Pattern,]
670	Parms,Parms2 = G2frame.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(G2frame,
671	G2frame.PatternId, 'Instrument Parameters'))
672	ParmList = [Parms,]
673	Title = Pattern[-1]
674	else:
675	Title = os.path.split(G2frame.GSASprojectfile)[1]
676	PlotList = []
677	ParmList = []
678	item, cookie = G2frame.PatternTree.GetFirstChild(G2frame.root)
679	while item:
680	if 'PWDR' in G2frame.PatternTree.GetItemText(item):
681	Pattern = G2frame.PatternTree.GetItemPyData(item)
682	if len(Pattern) < 3: # put name on end if needed
683	Pattern.append(G2frame.PatternTree.GetItemText(item))
684	PlotList.append(Pattern)
685	ParmList.append(G2frame.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(G2frame,
686	item,'Instrument Parameters'))[0])
687	item, cookie = G2frame.PatternTree.GetNextChild(G2frame.root, cookie)
688	Ymax = 1.0
689	lenX = 0
690	if PickId:
691	if G2frame.PatternTree.GetItemText(PickId) in ['Reflection Lists']:
692	Phases = G2frame.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(G2frame,PatternId,'Reflection Lists'))
693	HKL = []
694	if Phases:
695	for peak in Phases[G2frame.RefList]:
696	HKL.append(peak[:6])
697	HKL = np.array(HKL)
698	else:
699	HKL = np.array(G2frame.HKL)
700	for Pattern in PlotList:
701	xye = Pattern[1]
702	Ymax = max(Ymax,max(xye[1]))
703	offset = G2frame.Offset[0]*Ymax/100.0
704	if G2frame.logPlot:
705	Title = 'log('+Title+')'
706	Plot.set_title(Title)
707	if G2frame.qPlot:
708	Plot.set_xlabel(r'$q, \AA^{-1}$',fontsize=14)
709	else:
710	if 'C' in ParmList[0]['Type'][0]:
711	Plot.set_xlabel(r'$\mathsf{2\theta}$',fontsize=14)
712	else:
713	Plot.set_xlabel(r'TOF, $\mathsf{\mu}$s',fontsize=14)
714	if G2frame.Weight:
715	Plot.set_ylabel(r'$\mathsf{I/\sigma(I)}$',fontsize=14)
716	else:
717	if 'C' in ParmList[0]['Type'][0]:
718	Plot.set_ylabel('Intensity',fontsize=14)
719	else:
720	Plot.set_ylabel('Normalized intensity',fontsize=14)
721	if G2frame.Contour:
722	ContourZ = []
723	ContourY = []
724	Nseq = 0
725	if len(PlotList) < 2:
726	G2frame.Contour = False
727	for N,Pattern in enumerate(PlotList):
728	Parms = ParmList[N]
729	if 'C' in Parms['Type'][0]:
730	wave = G2mth.getWave(Parms)
731	else:
732	difC = Parms['difC'][1]
733	ifpicked = False
734	LimitId = 0
735	xye = np.array(Pattern[1])
736	if PickId:
737	ifpicked = Pattern[2] == G2frame.PatternTree.GetItemText(PatternId)
738	LimitId = G2gd.GetPatternTreeItemId(G2frame,PatternId, 'Limits')
739	if G2frame.qPlot:
740	Id = G2gd.GetPatternTreeItemId(G2frame,G2frame.root, Pattern[2])
741	if 'C' in Parms['Type'][0]:
742	X = 4np.pinpsind((xye[0]-Parms['Zero'][1])/2.0)/wave
743	else:
744	X = 2np.piParms['difC'][1]/(xye[0]-Parms['Zero'][1])
745	else:
746	X = xye[0]-Parms['Zero'][1]
747	if not lenX:
748	lenX = len(X)
749	Y = xye[1]+offset*N
750	if LimitId:
751	limits = np.array(G2frame.PatternTree.GetItemPyData(LimitId))
752	if G2frame.qPlot:
753	if 'C' in Parms['Type'][0]:
754	limits = 4np.pinpsind(limits/2.0)/wave
755	else:
756	limits = 2np.pidifC/limits
757	Lines.append(Plot.axvline(limits[1][0],color='g',dashes=(5,5),picker=3.))
758	Lines.append(Plot.axvline(limits[1][1],color='r',dashes=(5,5),picker=3.))
759	if G2frame.Contour:
760
761	if lenX == len(X):
762	ContourY.append(N)
763	ContourZ.append(Y)
764	ContourX = X
765	Nseq += 1
766	Plot.set_ylabel('Data sequence',fontsize=12)
767	else:
768	X += G2frame.Offset[1].005N
769	if ifpicked:
770	Z = xye[3]+offset*N
771	W = xye[4]+offset*N
772	D = xye[5]-Ymax*G2frame.delOffset
773	if G2frame.logPlot:
774	Plot.semilogy(X,Y,colors[N%6]+'+',picker=3.,clip_on=False,nonposy='mask')
775	Plot.semilogy(X,Z,colors[(N+1)%6],picker=False,nonposy='mask')
776	Plot.semilogy(X,W,colors[(N+2)%6],picker=False,nonposy='mask')
777	elif G2frame.Weight:
778	DY = xye[1]*np.sqrt(xye[2])
779	DYmax = max(DY)
780	DZ = xye[3]*np.sqrt(xye[2])
781	DS = xye[5]np.sqrt(xye[2])-DYmaxG2frame.delOffset
782	Plot.plot(X,DY,colors[N%6]+'+',picker=3.,clip_on=False)
783	Plot.plot(X,DZ,colors[(N+1)%6],picker=False)
784	Plot.plot(X,DS,colors[(N+3)%6],picker=False)
785	Plot.axhline(0.,color=wx.BLACK)
786	else:
787	if G2frame.SubBack:
788	Plot.plot(X,Y-W,colors[N%6]+'+',picker=3.,clip_on=False)
789	Plot.plot(X,Z-W,colors[(N+1)%6],picker=False)
790	else:
791	Plot.plot(X,Y,colors[N%6]+'+',picker=3.,clip_on=False)
792	Plot.plot(X,Z,colors[(N+1)%6],picker=False)
793	Plot.plot(X,W,colors[(N+2)%6],picker=False)
794	Plot.plot(X,D,colors[(N+3)%6],picker=False)
795	Plot.axhline(0.,color=wx.BLACK)
796	Page.canvas.SetToolTipString('')
797	if G2frame.PatternTree.GetItemText(PickId) == 'Peak List':
798	tip = 'On data point: Pick peak - L or R MB. On line: L-move, R-delete'
799	Page.canvas.SetToolTipString(tip)
800	data = G2frame.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(G2frame,PatternId, 'Peak List'))
801	for item in data:
802	if G2frame.qPlot:
803	if 'C' in Parms['Type'][0]:
804	Lines.append(Plot.axvline(4math.pisind(item[0]/2.)/wave,color=colors[N%6],picker=2.))
805	else:
806	Lines.append(Plot.axvline(2math.pidifC/item[0],color=colors[N%6],picker=2.))
807	else:
808	Lines.append(Plot.axvline(item[0],color=colors[N%6],picker=2.))
809	if G2frame.PatternTree.GetItemText(PickId) == 'Limits':
810	tip = 'On data point: Lower limit - L MB; Upper limit - R MB. On limit: MB down to move'
811	Page.canvas.SetToolTipString(tip)
812	data = G2frame.LimitsTable.GetData()
813	else:
814	if G2frame.logPlot:
815	Plot.semilogy(X,Y,colors[N%6],picker=False,nonposy='clip')
816	else:
817	Plot.plot(X,Y,colors[N%6],picker=False)
818	if PickId and not G2frame.Contour:
819	Parms,Parms2 = G2frame.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(G2frame,PatternId, 'Instrument Parameters'))
820	if 'C' in Parms['Type'][0]:
821	wave = G2mth.getWave(Parms)
822	else:
823	difC = Parms['difC'][1]
824	if G2frame.PatternTree.GetItemText(PickId) in ['Index Peak List','Unit Cells List']:
825	peaks = np.array((G2frame.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(G2frame,PatternId, 'Index Peak List'))))
826	for peak in peaks:
827	if G2frame.qPlot:
828	Plot.axvline(4np.pisind(peak[0]/2.0)/wave,color='b')
829	else:
830	Plot.axvline(peak[0],color='b')
831	for hkl in G2frame.HKL:
832	if G2frame.qPlot:
833	Plot.axvline(4np.pisind(hkl[5]/2.0)/wave,color='r',dashes=(5,5))
834	else:
835	Plot.axvline(hkl[5],color='r',dashes=(5,5))
836	elif G2frame.PatternTree.GetItemText(PickId) in ['Reflection Lists'] or \
837	'PWDR' in G2frame.PatternTree.GetItemText(PickId):
838	refColors=['b','r','c','g','m','k']
839	Phases = G2frame.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(G2frame,PatternId,'Reflection Lists'))
840	for pId,phase in enumerate(Phases):
841	peaks = Phases[phase]
842	if not peaks:
843	continue
844	peak = np.array([[peak[4],peak[5]] for peak in peaks])
845	pos = G2frame.refOffset-pIdYmaxG2frame.refDelt*np.ones_like(peak)
846	if G2frame.qPlot:
847	Plot.plot(2*np.pi/peak.T[0],pos,refColors[pId%6]+'\|',mew=1,ms=8,picker=3.,label=phase)
848	else:
849	Plot.plot(peak.T[1],pos,refColors[pId%6]+'\|',mew=1,ms=8,picker=3.,label=phase)
850	if len(Phases):
851	handles,legends = Plot.get_legend_handles_labels() #got double entries in the legends for some reason
852	if handles:
853	Plot.legend(handles[::2],legends[::2],title='Phases',loc='best') #skip every other one
854
855	if G2frame.Contour:
856	acolor = mpl.cm.get_cmap(G2frame.ContourColor)
857	Img = Plot.imshow(ContourZ,cmap=acolor,vmin=0,vmax=Ymax*G2frame.Cmax,interpolation=G2frame.Interpolate,
858	extent=[ContourX[0],ContourX[-1],ContourY[0],ContourY[-1]],aspect='auto',origin='lower')
859	Page.figure.colorbar(Img)
860	else:
861	G2frame.Lines = Lines
862	if not newPlot:
863	Page.toolbar.push_current()
864	Plot.set_xlim(xylim[0])
865	Plot.set_ylim(xylim[1])
866	xylim = []
867	Page.toolbar.push_current()
868	Page.toolbar.draw()
869	else:
870	Page.canvas.draw()
871	# G2frame.Pwdr = True
872
873	################################################################################
874	##### PlotDeltSig
875	################################################################################
876
877	def PlotDeltSig(G2frame,kind):
878	'needs a doc string'
879	try:
880	plotNum = G2frame.G2plotNB.plotList.index('Error analysis')
881	Page = G2frame.G2plotNB.nb.GetPage(plotNum)
882	Page.figure.clf()
883	Plot = Page.figure.gca() #get a fresh plot after clf()
884	except ValueError:
885	newPlot = True
886	G2frame.Cmax = 1.0
887	Plot = G2frame.G2plotNB.addMpl('Error analysis').gca()
888	plotNum = G2frame.G2plotNB.plotList.index('Error analysis')
889	Page = G2frame.G2plotNB.nb.GetPage(plotNum)
890	Page.Choice = None
891	PatternId = G2frame.PatternId
892	Pattern = G2frame.PatternTree.GetItemPyData(PatternId)
893	Pattern.append(G2frame.PatternTree.GetItemText(PatternId))
894	wtFactor = Pattern[0]['wtFactor']
895	if kind == 'PWDR':
896	limits = G2frame.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(G2frame,PatternId, 'Limits'))[1]
897	xye = np.array(Pattern[1])
898	xmin = np.searchsorted(xye[0],limits[0])
899	xmax = np.searchsorted(xye[0],limits[1])
900	DS = xye[5][xmin:xmax]np.sqrt(wtFactorxye[2][xmin:xmax])
901	elif kind == 'HKLF':
902	refl = Pattern[1]
903	DS = []
904	for ref in refl:
905	if ref[6] > 0.:
906	DS.append((ref[5]-ref[7])/ref[6])
907	Page.SetFocus()
908	G2frame.G2plotNB.status.DestroyChildren()
909	DS.sort()
910	EDS = np.zeros_like(DS)
911	DX = np.linspace(0.,1.,num=len(DS),endpoint=True)
912	oldErr = np.seterr(invalid='ignore') #avoid problem at DX==0
913	T = np.sqrt(np.log(1.0/DX**2))
914	top = 2.515517+0.802853T+0.010328T**2
915	bot = 1.0+1.432788T+0.189269T*2+0.001308T**3
916	EDS = np.where(DX>0,-(T-top/bot),(T-top/bot))
917	low1 = np.searchsorted(EDS,-1.)
918	hi1 = np.searchsorted(EDS,1.)
919	slp,intcp = np.polyfit(EDS[low1:hi1],DS[low1:hi1],deg=1)
920	frac = 100.*(hi1-low1)/len(DS)
921	G2frame.G2plotNB.status.SetStatusText( \
922	'Over range -1. to 1. :'+' slope = %.3f, intercept = %.3f for %.2f%% of the fitted data'%(slp,intcp,frac),1)
923	Plot.set_title('Normal probability for '+Pattern[-1])
924	Plot.set_xlabel(r'expected $\mathsf{\Delta/\sigma}$',fontsize=14)
925	Plot.set_ylabel(r'observed $\mathsf{\Delta/\sigma}$',fontsize=14)
926	Plot.plot(EDS,DS,'r+',label='result')
927	Plot.plot([-2,2],[-2,2],'k',dashes=(5,5),label='ideal')
928	Plot.legend(loc='upper left')
929	np.seterr(invalid='warn')
930	Page.canvas.draw()
931
932	################################################################################
933	##### PlotISFG
934	################################################################################
935
936	def PlotISFG(G2frame,newPlot=False,type=''):
937	''' PLotting package for PDF analysis; displays I(q), S(q), F(q) and G(r) as single
938	or multiple plots with waterfall and contour plots as options
939	'''
940	if not type:
941	type = G2frame.G2plotNB.plotList[G2frame.G2plotNB.nb.GetSelection()]
942	if type not in ['I(Q)','S(Q)','F(Q)','G(R)']:
943	return
944	superMinusOne = unichr(0xaf)+unichr(0xb9)
945
946	def OnPlotKeyPress(event):
947	newPlot = False
948	if event.key == 'u':
949	if G2frame.Contour:
950	G2frame.Cmax = min(1.0,G2frame.Cmax*1.2)
951	elif G2frame.Offset[0] < 100.:
952	G2frame.Offset[0] += 1.
953	elif event.key == 'd':
954	if G2frame.Contour:
955	G2frame.Cmax = max(0.0,G2frame.Cmax*0.8)
956	elif G2frame.Offset[0] > 0.:
957	G2frame.Offset[0] -= 1.
958	elif event.key == 'l':
959	G2frame.Offset[1] -= 1.
960	elif event.key == 'r':
961	G2frame.Offset[1] += 1.
962	elif event.key == 'o':
963	G2frame.Offset = [0,0]
964	elif event.key == 'c':
965	newPlot = True
966	if G2frame.Contour:
967	G2frame.Contour = False
968	else:
969	G2frame.Contour = True
970	G2frame.SinglePlot = False
971	G2frame.Offset = [0.,0.]
972	elif event.key == 's':
973	if G2frame.Contour:
974	choice = [m for m in mpl.cm.datad.keys() if not m.endswith("_r")]
975	choice.sort()
976	dlg = wx.SingleChoiceDialog(G2frame,'Select','Color scheme',choice)
977	if dlg.ShowModal() == wx.ID_OK:
978	sel = dlg.GetSelection()
979	G2frame.ContourColor = choice[sel]
980	else:
981	G2frame.ContourColor = 'Paired'
982	dlg.Destroy()
983	else:
984	if G2frame.SinglePlot:
985	G2frame.SinglePlot = False
986	else:
987	G2frame.SinglePlot = True
988	elif event.key == 'i': #for smoothing contour plot
989	choice = ['nearest','bilinear','bicubic','spline16','spline36','hanning',
990	'hamming','hermite','kaiser','quadric','catrom','gaussian','bessel',
991	'mitchell','sinc','lanczos']
992	dlg = wx.SingleChoiceDialog(G2frame,'Select','Interpolation',choice)
993	if dlg.ShowModal() == wx.ID_OK:
994	sel = dlg.GetSelection()
995	G2frame.Interpolate = choice[sel]
996	else:
997	G2frame.Interpolate = 'nearest'
998	dlg.Destroy()
999	elif event.key == 't' and not G2frame.Contour:
1000	if G2frame.Legend:
1001	G2frame.Legend = False
1002	else:
1003	G2frame.Legend = True
1004
1005
1006	PlotISFG(G2frame,newPlot=newPlot,type=type)
1007
1008	def OnKeyBox(event):
1009	if G2frame.G2plotNB.nb.GetSelection() == G2frame.G2plotNB.plotList.index(type):
1010	event.key = cb.GetValue()[0]
1011	cb.SetValue(' key press')
1012	wx.CallAfter(OnPlotKeyPress,event)
1013
1014	def OnMotion(event):
1015	xpos = event.xdata
1016	if xpos: #avoid out of frame mouse position
1017	ypos = event.ydata
1018	Page.canvas.SetCursor(wx.CROSS_CURSOR)
1019	try:
1020	if G2frame.Contour:
1021	G2frame.G2plotNB.status.SetStatusText('R =%.3fA pattern ID =%5d'%(xpos,int(ypos)),1)
1022	else:
1023	G2frame.G2plotNB.status.SetStatusText('R =%.3fA %s =%.2f'%(xpos,type,ypos),1)
1024	except TypeError:
1025	G2frame.G2plotNB.status.SetStatusText('Select '+type+' pattern first',1)
1026
1027	xylim = []
1028	try:
1029	plotNum = G2frame.G2plotNB.plotList.index(type)
1030	Page = G2frame.G2plotNB.nb.GetPage(plotNum)
1031	if not newPlot:
1032	Plot = Page.figure.gca() #get previous plot & get limits
1033	xylim = Plot.get_xlim(),Plot.get_ylim()
1034	Page.figure.clf()
1035	Plot = Page.figure.gca()
1036	except ValueError:
1037	newPlot = True
1038	G2frame.Cmax = 1.0
1039	Plot = G2frame.G2plotNB.addMpl(type).gca()
1040	plotNum = G2frame.G2plotNB.plotList.index(type)
1041	Page = G2frame.G2plotNB.nb.GetPage(plotNum)
1042	Page.canvas.mpl_connect('key_press_event', OnPlotKeyPress)
1043	Page.canvas.mpl_connect('motion_notify_event', OnMotion)
1044
1045	Page.SetFocus()
1046	G2frame.G2plotNB.status.DestroyChildren()
1047	if G2frame.Contour:
1048	Page.Choice = (' key press','d: lower contour max','u: raise contour max',
1049	'i: interpolation method','s: color scheme','c: contour off')
1050	else:
1051	Page.Choice = (' key press','l: offset left','r: offset right','d: offset down','u: offset up',
1052	'o: reset offset','t: toggle legend','c: contour on','s: toggle single plot')
1053	Page.keyPress = OnPlotKeyPress
1054	PatternId = G2frame.PatternId
1055	PickId = G2frame.PickId
1056	Plot.set_title(type)
1057	if type == 'G(R)':
1058	Plot.set_xlabel(r'$R,\AA$',fontsize=14)
1059	else:
1060	Plot.set_xlabel(r'$Q,\AA$'+superMinusOne,fontsize=14)
1061	Plot.set_ylabel(r''+type,fontsize=14)
1062	colors=['b','g','r','c','m','k']
1063	name = G2frame.PatternTree.GetItemText(PatternId)[4:]
1064	Pattern = []
1065	if G2frame.SinglePlot:
1066	name = G2frame.PatternTree.GetItemText(PatternId)
1067	name = type+name[4:]
1068	Id = G2gd.GetPatternTreeItemId(G2frame,PatternId,name)
1069	Pattern = G2frame.PatternTree.GetItemPyData(Id)
1070	if Pattern:
1071	Pattern.append(name)
1072	PlotList = [Pattern,]
1073	else:
1074	PlotList = []
1075	item, cookie = G2frame.PatternTree.GetFirstChild(G2frame.root)
1076	while item:
1077	if 'PDF' in G2frame.PatternTree.GetItemText(item):
1078	name = type+G2frame.PatternTree.GetItemText(item)[4:]
1079	Id = G2gd.GetPatternTreeItemId(G2frame,item,name)
1080	Pattern = G2frame.PatternTree.GetItemPyData(Id)
1081	if Pattern:
1082	Pattern.append(name)
1083	PlotList.append(Pattern)
1084	item, cookie = G2frame.PatternTree.GetNextChild(G2frame.root, cookie)
1085	PDFdata = G2frame.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(G2frame,PatternId, 'PDF Controls'))
1086	numbDen = G2pwd.GetNumDensity(PDFdata['ElList'],PDFdata['Form Vol'])
1087	Xb = [0.,10.]
1088	Yb = [0.,-40.np.pinumbDen]
1089	Ymax = 1.0
1090	lenX = 0
1091	for Pattern in PlotList:
1092	xye = Pattern[1]
1093	Ymax = max(Ymax,max(xye[1]))
1094	offset = G2frame.Offset[0]*Ymax/100.0
1095	if G2frame.Contour:
1096	ContourZ = []
1097	ContourY = []
1098	Nseq = 0
1099	for N,Pattern in enumerate(PlotList):
1100	xye = Pattern[1]
1101	if PickId:
1102	ifpicked = Pattern[2] == G2frame.PatternTree.GetItemText(PatternId)
1103	X = xye[0]
1104	if not lenX:
1105	lenX = len(X)
1106	Y = xye[1]+offset*N
1107	if G2frame.Contour:
1108	if lenX == len(X):
1109	ContourY.append(N)
1110	ContourZ.append(Y)
1111	ContourX = X
1112	Nseq += 1
1113	Plot.set_ylabel('Data sequence',fontsize=12)
1114	else:
1115	X = xye[0]+G2frame.Offset[1].005N
1116	if ifpicked:
1117	Plot.plot(X,Y,colors[N%6]+'+',picker=3.,clip_on=False)
1118	Page.canvas.SetToolTipString('')
1119	else:
1120	if G2frame.Legend:
1121	Plot.plot(X,Y,colors[N%6],picker=False,label='Azm:'+Pattern[2].split('=')[1])
1122	else:
1123	Plot.plot(X,Y,colors[N%6],picker=False)
1124	if type == 'G(R)':
1125	Plot.plot(Xb,Yb,color='k',dashes=(5,5))
1126	elif type == 'F(Q)':
1127	Plot.axhline(0.,color=wx.BLACK)
1128	elif type == 'S(Q)':
1129	Plot.axhline(1.,color=wx.BLACK)
1130	if G2frame.Contour:
1131	acolor = mpl.cm.get_cmap(G2frame.ContourColor)
1132	Img = Plot.imshow(ContourZ,cmap=acolor,vmin=0,vmax=Ymax*G2frame.Cmax,interpolation=G2frame.Interpolate,
1133	extent=[ContourX[0],ContourX[-1],ContourY[0],ContourY[-1]],aspect='auto',origin='lower')
1134	Page.figure.colorbar(Img)
1135	elif G2frame.Legend:
1136	Plot.legend(loc='best')
1137	if not newPlot:
1138	Page.toolbar.push_current()
1139	Plot.set_xlim(xylim[0])
1140	Plot.set_ylim(xylim[1])
1141	xylim = []
1142	Page.toolbar.push_current()
1143	Page.toolbar.draw()
1144	else:
1145	Page.canvas.draw()
1146
1147	################################################################################
1148	##### PlotXY
1149	################################################################################
1150
1151	def PlotXY(G2frame,XY,newPlot=False,type=''):
1152	'''simple plot of xy data, used for diagnostic purposes
1153	'''
1154	def OnMotion(event):
1155	xpos = event.xdata
1156	if xpos: #avoid out of frame mouse position
1157	ypos = event.ydata
1158	Page.canvas.SetCursor(wx.CROSS_CURSOR)
1159	try:
1160	G2frame.G2plotNB.status.SetStatusText('X =%9.3f %s =%9.3f'%(xpos,type,ypos),1)
1161	except TypeError:
1162	G2frame.G2plotNB.status.SetStatusText('Select '+type+' pattern first',1)
1163
1164	try:
1165	plotNum = G2frame.G2plotNB.plotList.index(type)
1166	Page = G2frame.G2plotNB.nb.GetPage(plotNum)
1167	if not newPlot:
1168	Plot = Page.figure.gca()
1169	xylim = Plot.get_xlim(),Plot.get_ylim()
1170	Page.figure.clf()
1171	Plot = Page.figure.gca()
1172	except ValueError:
1173	newPlot = True
1174	Plot = G2frame.G2plotNB.addMpl(type).gca()
1175	plotNum = G2frame.G2plotNB.plotList.index(type)
1176	Page = G2frame.G2plotNB.nb.GetPage(plotNum)
1177	Page.canvas.mpl_connect('motion_notify_event', OnMotion)
1178
1179	Page.Choice = None
1180	Page.SetFocus()
1181	G2frame.G2plotNB.status.DestroyChildren()
1182	Plot.set_title(type)
1183	Plot.set_xlabel(r'X',fontsize=14)
1184	Plot.set_ylabel(r''+type,fontsize=14)
1185	colors=['b','g','r','c','m','k']
1186	Ymax = 1.0
1187	lenX = 0
1188	X,Y = XY[:2]
1189	Ymax = max(Ymax,max(Y))
1190	Plot.plot(X,Y,'k',picker=False)
1191	if not newPlot:
1192	Page.toolbar.push_current()
1193	Plot.set_xlim(xylim[0])
1194	Plot.set_ylim(xylim[1])
1195	xylim = []
1196	Page.toolbar.push_current()
1197	Page.toolbar.draw()
1198	else:
1199	Page.canvas.draw()
1200
1201	################################################################################
1202	##### PlotPowderLines
1203	################################################################################
1204
1205	def PlotPowderLines(G2frame):
1206	''' plotting of powder lines (i.e. no powder pattern) as sticks
1207	'''
1208	global HKL
1209
1210	def OnMotion(event):
1211	xpos = event.xdata
1212	if xpos: #avoid out of frame mouse position
1213	Page.canvas.SetCursor(wx.CROSS_CURSOR)
1214	G2frame.G2plotNB.status.SetFields(['','2-theta =%9.3f '%(xpos,)])
1215	if G2frame.PickId and G2frame.PatternTree.GetItemText(G2frame.PickId) in ['Index Peak List','Unit Cells List']:
1216	found = []
1217	if len(HKL):
1218	view = Page.toolbar._views.forward()[0][:2]
1219	wid = view[1]-view[0]
1220	found = HKL[np.where(np.fabs(HKL.T[5]-xpos) < 0.002*wid)]
1221	if len(found):
1222	h,k,l = found[0][:3]
1223	Page.canvas.SetToolTipString('%d,%d,%d'%(int(h),int(k),int(l)))
1224	else:
1225	Page.canvas.SetToolTipString('')
1226
1227	try:
1228	plotNum = G2frame.G2plotNB.plotList.index('Powder Lines')
1229	Page = G2frame.G2plotNB.nb.GetPage(plotNum)
1230	Page.figure.clf()
1231	Plot = Page.figure.gca()
1232	except ValueError:
1233	Plot = G2frame.G2plotNB.addMpl('Powder Lines').gca()
1234	plotNum = G2frame.G2plotNB.plotList.index('Powder Lines')
1235	Page = G2frame.G2plotNB.nb.GetPage(plotNum)
1236	Page.canvas.mpl_connect('motion_notify_event', OnMotion)
1237
1238	Page.Choice = None
1239	Page.SetFocus()
1240	Plot.set_title('Powder Pattern Lines')
1241	Plot.set_xlabel(r'$\mathsf{2\theta}$',fontsize=14)
1242	PickId = G2frame.PickId
1243	PatternId = G2frame.PatternId
1244	peaks = G2frame.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(G2frame,PatternId, 'Index Peak List'))
1245	for peak in peaks:
1246	Plot.axvline(peak[0],color='b')
1247	HKL = np.array(G2frame.HKL)
1248	for hkl in G2frame.HKL:
1249	Plot.axvline(hkl[5],color='r',dashes=(5,5))
1250	xmin = peaks[0][0]
1251	xmax = peaks[-1][0]
1252	delt = xmax-xmin
1253	xlim = [max(0,xmin-delt/20.),min(180.,xmax+delt/20.)]
1254	Plot.set_xlim(xlim)
1255	Page.canvas.draw()
1256	Page.toolbar.push_current()
1257
1258	################################################################################
1259	##### PlotPeakWidths
1260	################################################################################
1261
1262	def PlotPeakWidths(G2frame):
1263	''' Plotting of instrument broadening terms as function of 2-theta
1264	Seen when "Instrument Parameters" chosen from powder pattern data tree
1265	'''
1266	# sig = lambda Th,U,V,W: 1.17741math.sqrt(Utand(Th)*2+Vtand(Th)+W)*math.pi/18000.
1267	# gam = lambda Th,X,Y: (X/cosd(Th)+Ytand(Th))math.pi/18000.
1268	gamFW = lambda s,g: np.exp(np.log(s*5+2.69269s*4g+2.42843s3g*2+4.47163s*2g*3+0.07842sg4+g*5)/5.)
1269	# gamFW2 = lambda s,g: math.sqrt(s*2+(0.4654996g)*2)+.5345004g #Ubaldo Bafile - private communication
1270	PatternId = G2frame.PatternId
1271	limitID = G2gd.GetPatternTreeItemId(G2frame,PatternId, 'Limits')
1272	if limitID:
1273	limits = G2frame.PatternTree.GetItemPyData(limitID)
1274	else:
1275	return
1276	Parms,Parms2 = G2frame.PatternTree.GetItemPyData( \
1277	G2gd.GetPatternTreeItemId(G2frame,PatternId, 'Instrument Parameters'))
1278	if 'C' in Parms['Type'][0]:
1279	lam = G2mth.getWave(Parms)
1280	else:
1281	difC = Parms['difC'][0]
1282	peakID = G2gd.GetPatternTreeItemId(G2frame,PatternId, 'Peak List')
1283	if peakID:
1284	peaks = G2frame.PatternTree.GetItemPyData(peakID)
1285	else:
1286	peaks = []
1287
1288	try:
1289	plotNum = G2frame.G2plotNB.plotList.index('Peak Widths')
1290	Page = G2frame.G2plotNB.nb.GetPage(plotNum)
1291	Page.figure.clf()
1292	Plot = Page.figure.gca()
1293	except ValueError:
1294	Plot = G2frame.G2plotNB.addMpl('Peak Widths').gca()
1295	plotNum = G2frame.G2plotNB.plotList.index('Peak Widths')
1296	Page = G2frame.G2plotNB.nb.GetPage(plotNum)
1297	Page.Choice = None
1298	Page.SetFocus()
1299
1300	Page.canvas.SetToolTipString('')
1301	colors=['b','g','r','c','m','k']
1302	X = []
1303	Y = []
1304	Z = []
1305	W = []
1306	if 'C' in Parms['Type'][0]:
1307	Plot.set_title('Instrument and sample peak widths')
1308	Plot.set_xlabel(r'$q, \AA^{-1}$',fontsize=14)
1309	Plot.set_ylabel(r'$\Delta q/q, \Delta d/d$',fontsize=14)
1310	try:
1311	Xmin,Xmax = limits[1]
1312	X = np.linspace(Xmin,Xmax,num=101,endpoint=True)
1313	Q = 4.np.pinpsind(X/2.)/lam
1314	Z = np.ones_like(X)
1315	data = G2mth.setPeakparms(Parms,Parms2,X,Z)
1316	s = 1.17741np.sqrt(data[4])np.pi/18000.
1317	g = data[6]*np.pi/18000.
1318	G = gamFW(g,s)
1319	Y = s/nptand(X/2.)
1320	Z = g/nptand(X/2.)
1321	W = G/nptand(X/2.)
1322	Plot.plot(Q,Y,color='r',label='Gaussian')
1323	Plot.plot(Q,Z,color='g',label='Lorentzian')
1324	Plot.plot(Q,W,color='b',label='G+L')
1325
1326	fit = G2mth.setPeakparms(Parms,Parms2,X,Z,useFit=True)
1327	sf = 1.17741np.sqrt(fit[4])np.pi/18000.
1328	gf = fit[6]*np.pi/18000.
1329	Gf = gamFW(gf,sf)
1330	Yf = sf/nptand(X/2.)
1331	Zf = gf/nptand(X/2.)
1332	Wf = Gf/nptand(X/2.)
1333	Plot.plot(Q,Yf,color='r',dashes=(5,5),label='Gaussian fit')
1334	Plot.plot(Q,Zf,color='g',dashes=(5,5),label='Lorentzian fit')
1335	Plot.plot(Q,Wf,color='b',dashes=(5,5),label='G+L fit')
1336
1337	X = []
1338	Y = []
1339	Z = []
1340	W = []
1341	V = []
1342	for peak in peaks:
1343	X.append(4.0math.pisind(peak[0]/2.0)/lam)
1344	try:
1345	s = 1.17741math.sqrt(peak[4])math.pi/18000.
1346	except ValueError:
1347	s = 0.01
1348	g = peak[6]*math.pi/18000.
1349	G = gamFW(g,s)
1350	Y.append(s/tand(peak[0]/2.))
1351	Z.append(g/tand(peak[0]/2.))
1352	W.append(G/tand(peak[0]/2.))
1353	if len(peaks):
1354	Plot.plot(X,Y,'+',color='r',label='G peak')
1355	Plot.plot(X,Z,'+',color='g',label='L peak')
1356	Plot.plot(X,W,'+',color='b',label='G+L peak')
1357	Plot.legend(loc='best')
1358	Page.canvas.draw()
1359	except ValueError:
1360	print '**** ERROR - default U,V,W profile coefficients yield sqrt of negative value at 2theta =', \
1361	'%.3f'%(2*theta)
1362	G2frame.G2plotNB.Delete('Peak Widths')
1363	else:
1364	Plot.set_title('Instrument and sample peak coefficients')
1365	Plot.set_xlabel(r'$q, \AA^{-1}$',fontsize=14)
1366	Plot.set_ylabel(r'$\alpha, \beta, \Delta q/q, \Delta d/d$',fontsize=14)
1367	Xmin,Xmax = limits[1]
1368	T = np.linspace(Xmin,Xmax,num=101,endpoint=True)
1369	Z = np.ones_like(T)
1370	data = G2mth.setPeakparms(Parms,Parms2,T,Z)
1371	ds = T/difC
1372	Q = 2.*np.pi/ds
1373	A = data[4]
1374	B = data[6]
1375	S = 1.17741*np.sqrt(data[8])/T
1376	G = data[10]/T
1377	Plot.plot(Q,A,color='r',label='Alpha')
1378	Plot.plot(Q,B,color='g',label='Beta')
1379	Plot.plot(Q,S,color='b',label='Gaussian')
1380	Plot.plot(Q,G,color='m',label='Lorentzian')
1381
1382	fit = G2mth.setPeakparms(Parms,Parms2,T,Z)
1383	ds = T/difC
1384	Q = 2.*np.pi/ds
1385	Af = fit[4]
1386	Bf = fit[6]
1387	Sf = 1.17741*np.sqrt(fit[8])/T
1388	Gf = fit[10]/T
1389	Plot.plot(Q,Af,color='r',dashes=(5,5),label='Alpha fit')
1390	Plot.plot(Q,Bf,color='g',dashes=(5,5),label='Beta fit')
1391	Plot.plot(Q,Sf,color='b',dashes=(5,5),label='Gaussian fit')
1392	Plot.plot(Q,Gf,color='m',dashes=(5,5),label='Lorentzian fit')
1393
1394	T = []
1395	A = []
1396	B = []
1397	S = []
1398	G = []
1399	W = []
1400	Q = []
1401	V = []
1402	for peak in peaks:
1403	T.append(peak[0])
1404	A.append(peak[4])
1405	B.append(peak[6])
1406	Q.append(2.np.pidifC/peak[0])
1407	S.append(1.17741*np.sqrt(peak[8])/peak[0])
1408	G.append(peak[10]/peak[0])
1409
1410
1411	Plot.plot(Q,A,'+',color='r',label='Alpha peak')
1412	Plot.plot(Q,B,'+',color='g',label='Beta peak')
1413	Plot.plot(Q,S,'+',color='b',label='Gaussian peak')
1414	Plot.plot(Q,G,'+',color='m',label='Lorentzian peak')
1415	Plot.legend(loc='best')
1416	Page.canvas.draw()
1417
1418
1419	################################################################################
1420	##### PlotSizeStrainPO
1421	################################################################################
1422
1423	def PlotSizeStrainPO(G2frame,data,Start=False):
1424	'''Plot 3D mustrain/size/preferred orientation figure. In this instance data is for a phase
1425	'''
1426
1427	PatternId = G2frame.PatternId
1428	generalData = data['General']
1429	SGData = generalData['SGData']
1430	SGLaue = SGData['SGLaue']
1431	if Start: #initialize the spherical harmonics qlmn arrays
1432	ptx.pyqlmninit()
1433	Start = False
1434	# MuStrKeys = G2spc.MustrainNames(SGData)
1435	cell = generalData['Cell'][1:]
1436	A,B = G2lat.cell2AB(cell[:6])
1437	Vol = cell[6]
1438	useList = data['Histograms']
1439	phase = generalData['Name']
1440	plotType = generalData['Data plot type']
1441	plotDict = {'Mustrain':'Mustrain','Size':'Size','Preferred orientation':'Pref.Ori.'}
1442	for ptype in plotDict:
1443	G2frame.G2plotNB.Delete(ptype)
1444	if plotType in ['None']:
1445	return
1446
1447	for item in useList:
1448	if useList[item]['Show']:
1449	break
1450	else:
1451	return #nothing to show!!
1452
1453	numPlots = len(useList)
1454
1455	if plotType in ['Mustrain','Size']:
1456	Plot = mp3d.Axes3D(G2frame.G2plotNB.add3D(plotType))
1457	else:
1458	Plot = G2frame.G2plotNB.addMpl(plotType).gca()
1459	plotNum = G2frame.G2plotNB.plotList.index(plotType)
1460	Page = G2frame.G2plotNB.nb.GetPage(plotNum)
1461	Page.Choice = None
1462	Page.SetFocus()
1463	G2frame.G2plotNB.status.SetStatusText('',1)
1464	if not Page.IsShown():
1465	Page.Show()
1466
1467	for item in useList:
1468	if useList[item]['Show']:
1469	PHI = np.linspace(0.,360.,30,True)
1470	PSI = np.linspace(0.,180.,30,True)
1471	X = np.outer(npsind(PHI),npsind(PSI))
1472	Y = np.outer(npcosd(PHI),npsind(PSI))
1473	Z = np.outer(np.ones(np.size(PHI)),npcosd(PSI))
1474	try: #temp patch instead of 'mustrain' for old files with 'microstrain'
1475	coeff = useList[item][plotDict[plotType]]
1476	except KeyError:
1477	break
1478	if plotType in ['Mustrain','Size']:
1479	if coeff[0] == 'isotropic':
1480	X *= coeff[1][0]
1481	Y *= coeff[1][0]
1482	Z *= coeff[1][0]
1483	elif coeff[0] == 'uniaxial':
1484
1485	def uniaxCalc(xyz,iso,aniso,axes):
1486	Z = np.array(axes)
1487	cp = abs(np.dot(xyz,Z))
1488	sp = np.sqrt(1.-cp**2)
1489	R = isoaniso/np.sqrt((isocp)*2+(anisosp)**2)
1490	return R*xyz
1491
1492	iso,aniso = coeff[1][:2]
1493	axes = np.inner(A,np.array(coeff[3]))
1494	axes /= nl.norm(axes)
1495	Shkl = np.array(coeff[1])
1496	XYZ = np.dstack((X,Y,Z))
1497	XYZ = np.nan_to_num(np.apply_along_axis(uniaxCalc,2,XYZ,iso,aniso,axes))
1498	X,Y,Z = np.dsplit(XYZ,3)
1499	X = X[:,:,0]
1500	Y = Y[:,:,0]
1501	Z = Z[:,:,0]
1502
1503	elif coeff[0] == 'ellipsoidal':
1504
1505	def ellipseCalc(xyz,E,R):
1506	XYZ = xyz*E.T
1507	return np.inner(XYZ.T,R)
1508
1509	S6 = coeff[4]
1510	Sij = G2lat.U6toUij(S6)
1511	E,R = nl.eigh(Sij)
1512	XYZ = np.dstack((X,Y,Z))
1513	XYZ = np.nan_to_num(np.apply_along_axis(ellipseCalc,2,XYZ,E,R))
1514	X,Y,Z = np.dsplit(XYZ,3)
1515	X = X[:,:,0]
1516	Y = Y[:,:,0]
1517	Z = Z[:,:,0]
1518
1519	elif coeff[0] == 'generalized':
1520
1521	def genMustrain(xyz,SGData,A,Shkl):
1522	uvw = np.inner(A.T,xyz)
1523	Strm = np.array(G2spc.MustrainCoeff(uvw,SGData))
1524	sum = np.sum(np.multiply(Shkl,Strm))
1525	sum = np.where(sum > 0.01,sum,0.01)
1526	sum = np.sqrt(sum)*math.pi/.18 #centidegrees to radians!
1527	return sum*xyz
1528
1529	Shkl = np.array(coeff[4])
1530	if np.any(Shkl):
1531	XYZ = np.dstack((X,Y,Z))
1532	XYZ = np.nan_to_num(np.apply_along_axis(genMustrain,2,XYZ,SGData,A,Shkl))
1533	X,Y,Z = np.dsplit(XYZ,3)
1534	X = X[:,:,0]
1535	Y = Y[:,:,0]
1536	Z = Z[:,:,0]
1537
1538	if np.any(X) and np.any(Y) and np.any(Z):
1539	errFlags = np.seterr(all='ignore')
1540	Plot.plot_surface(X,Y,Z,rstride=1,cstride=1,color='g',linewidth=1)
1541	np.seterr(all='ignore')
1542	xyzlim = np.array([Plot.get_xlim3d(),Plot.get_ylim3d(),Plot.get_zlim3d()]).T
1543	XYZlim = [min(xyzlim[0]),max(xyzlim[1])]
1544	Plot.set_xlim3d(XYZlim)
1545	Plot.set_ylim3d(XYZlim)
1546	Plot.set_zlim3d(XYZlim)
1547	Plot.set_aspect('equal')
1548	if plotType == 'Size':
1549	Plot.set_title('Crystallite size for '+phase+'\n'+coeff[0]+' model')
1550	Plot.set_xlabel(r'X, $\mu$m')
1551	Plot.set_ylabel(r'Y, $\mu$m')
1552	Plot.set_zlabel(r'Z, $\mu$m')
1553	else:
1554	Plot.set_title(r'$\mu$strain for '+phase+'\n'+coeff[0]+' model')
1555	Plot.set_xlabel(r'X, $\mu$strain')
1556	Plot.set_ylabel(r'Y, $\mu$strain')
1557	Plot.set_zlabel(r'Z, $\mu$strain')
1558	else:
1559	h,k,l = generalData['POhkl']
1560	if coeff[0] == 'MD':
1561	print 'March-Dollase preferred orientation plot'
1562
1563	else:
1564	PH = np.array(generalData['POhkl'])
1565	phi,beta = G2lat.CrsAng(PH,cell[:6],SGData)
1566	SHCoef = {}
1567	for item in coeff[5]:
1568	L,N = eval(item.strip('C'))
1569	SHCoef['C%d,0,%d'%(L,N)] = coeff[5][item]
1570	ODFln = G2lat.Flnh(Start,SHCoef,phi,beta,SGData)
1571	X = np.linspace(0,90.0,26)
1572	Y = G2lat.polfcal(ODFln,'0',X,0.0)
1573	Plot.plot(X,Y,color='k',label=str(PH))
1574	Plot.legend(loc='best')
1575	Plot.set_title('Axial distribution for HKL='+str(PH)+' in '+phase)
1576	Plot.set_xlabel(r'$\psi$',fontsize=16)
1577	Plot.set_ylabel('MRD',fontsize=14)
1578	Page.canvas.draw()
1579
1580	################################################################################
1581	##### PlotTexture
1582	################################################################################
1583
1584	def PlotTexture(G2frame,data,Start=False):
1585	'''Pole figure, inverse pole figure, 3D pole distribution and 3D inverse pole distribution
1586	plotting.
1587	dict generalData contains all phase info needed which is in data
1588	'''
1589
1590	shModels = ['cylindrical','none','shear - 2/m','rolling - mmm']
1591	SamSym = dict(zip(shModels,['0','-1','2/m','mmm']))
1592	PatternId = G2frame.PatternId
1593	generalData = data['General']
1594	SGData = generalData['SGData']
1595	textureData = generalData['SH Texture']
1596	G2frame.G2plotNB.Delete('Texture')
1597	if not textureData['Order']:
1598	return #no plot!!
1599	SHData = generalData['SH Texture']
1600	SHCoef = SHData['SH Coeff'][1]
1601	cell = generalData['Cell'][1:7]
1602	Amat,Bmat = G2lat.cell2AB(cell)
1603	sq2 = 1.0/math.sqrt(2.0)
1604
1605	def rp2xyz(r,p):
1606	z = npcosd(r)
1607	xy = np.sqrt(1.-z**2)
1608	return xynpsind(p),xynpcosd(p),z
1609
1610	def OnMotion(event):
1611	SHData = data['General']['SH Texture']
1612	if event.xdata and event.ydata: #avoid out of frame errors
1613	xpos = event.xdata
1614	ypos = event.ydata
1615	if 'Inverse' in SHData['PlotType']:
1616	r = xpos2+ypos2
1617	if r <= 1.0:
1618	if 'equal' in G2frame.Projection:
1619	r,p = 2.npasind(np.sqrt(r)sq2),npatan2d(ypos,xpos)
1620	else:
1621	r,p = 2.*npatand(np.sqrt(r)),npatan2d(ypos,xpos)
1622	ipf = G2lat.invpolfcal(IODFln,SGData,np.array([r,]),np.array([p,]))
1623	xyz = np.inner(Amat.T,np.array([rp2xyz(r,p)]))
1624	y,x,z = list(xyz/np.max(np.abs(xyz)))
1625
1626	G2frame.G2plotNB.status.SetFields(['',
1627	'psi =%9.3f, beta =%9.3f, MRD =%9.3f hkl=%5.2f,%5.2f,%5.2f'%(r,p,ipf,x,y,z)])
1628
1629	elif 'Axial' in SHData['PlotType']:
1630	pass
1631
1632	else: #ordinary pole figure
1633	z = xpos2+ypos2
1634	if z <= 1.0:
1635	z = np.sqrt(z)
1636	if 'equal' in G2frame.Projection:
1637	r,p = 2.npasind(zsq2),npatan2d(ypos,xpos)
1638	else:
1639	r,p = 2.*npatand(z),npatan2d(ypos,xpos)
1640	pf = G2lat.polfcal(ODFln,SamSym[textureData['Model']],np.array([r,]),np.array([p,]))
1641	G2frame.G2plotNB.status.SetFields(['','phi =%9.3f, gam =%9.3f, MRD =%9.3f'%(r,p,pf)])
1642
1643	try:
1644	plotNum = G2frame.G2plotNB.plotList.index('Texture')
1645	Page = G2frame.G2plotNB.nb.GetPage(plotNum)
1646	Page.figure.clf()
1647	Plot = Page.figure.gca()
1648	if not Page.IsShown():
1649	Page.Show()
1650	except ValueError:
1651	Plot = G2frame.G2plotNB.addMpl('Texture').gca()
1652	plotNum = G2frame.G2plotNB.plotList.index('Texture')
1653	Page = G2frame.G2plotNB.nb.GetPage(plotNum)
1654	Page.canvas.mpl_connect('motion_notify_event', OnMotion)
1655
1656	Page.Choice = None
1657	Page.SetFocus()
1658	G2frame.G2plotNB.status.SetFields(['',''])
1659	PH = np.array(SHData['PFhkl'])
1660	phi,beta = G2lat.CrsAng(PH,cell,SGData)
1661	ODFln = G2lat.Flnh(Start,SHCoef,phi,beta,SGData)
1662	PX = np.array(SHData['PFxyz'])
1663	gam = atan2d(PX[0],PX[1])
1664	xy = math.sqrt(PX[0]2+PX[1]2)
1665	xyz = math.sqrt(PX[0]2+PX[1]2+PX[2]**2)
1666	psi = asind(xy/xyz)
1667	IODFln = G2lat.Glnh(Start,SHCoef,psi,gam,SamSym[textureData['Model']])
1668	if 'Axial' in SHData['PlotType']:
1669	X = np.linspace(0,90.0,26)
1670	Y = G2lat.polfcal(ODFln,SamSym[textureData['Model']],X,0.0)
1671	Plot.plot(X,Y,color='k',label=str(SHData['PFhkl']))
1672	Plot.legend(loc='best')
1673	Plot.set_title('Axial distribution for HKL='+str(SHData['PFhkl']))
1674	Plot.set_xlabel(r'$\psi$',fontsize=16)
1675	Plot.set_ylabel('MRD',fontsize=14)
1676
1677	else:
1678	npts = 201
1679	if 'Inverse' in SHData['PlotType']:
1680	X,Y = np.meshgrid(np.linspace(1.,-1.,npts),np.linspace(-1.,1.,npts))
1681	R,P = np.sqrt(X2+Y2).flatten(),npatan2d(X,Y).flatten()
1682	if 'equal' in G2frame.Projection:
1683	R = np.where(R <= 1.,2.npasind(Rsq2),0.0)
1684	else:
1685	R = np.where(R <= 1.,2.*npatand(R),0.0)
1686	Z = np.zeros_like(R)
1687	Z = G2lat.invpolfcal(IODFln,SGData,R,P)
1688	Z = np.reshape(Z,(npts,npts))
1689	CS = Plot.contour(Y,X,Z,aspect='equal')
1690	Plot.clabel(CS,fontsize=9,inline=1)
1691	try:
1692	Img = Plot.imshow(Z.T,aspect='equal',cmap=G2frame.ContourColor,extent=[-1,1,-1,1])
1693	except ValueError:
1694	pass
1695	Page.figure.colorbar(Img)
1696	Plot.set_title('Inverse pole figure for XYZ='+str(SHData['PFxyz']))
1697	Plot.set_xlabel(G2frame.Projection.capitalize()+' projection')
1698
1699	else:
1700	X,Y = np.meshgrid(np.linspace(1.,-1.,npts),np.linspace(-1.,1.,npts))
1701	R,P = np.sqrt(X2+Y2).flatten(),npatan2d(X,Y).flatten()
1702	if 'equal' in G2frame.Projection:
1703	R = np.where(R <= 1.,2.npasind(Rsq2),0.0)
1704	else:
1705	R = np.where(R <= 1.,2.*npatand(R),0.0)
1706	Z = np.zeros_like(R)
1707	Z = G2lat.polfcal(ODFln,SamSym[textureData['Model']],R,P)
1708	Z = np.reshape(Z,(npts,npts))
1709	try:
1710	CS = Plot.contour(Y,X,Z,aspect='equal')
1711	Plot.clabel(CS,fontsize=9,inline=1)
1712	except ValueError:
1713	pass
1714	Img = Plot.imshow(Z.T,aspect='equal',cmap=G2frame.ContourColor,extent=[-1,1,-1,1])
1715	Page.figure.colorbar(Img)
1716	Plot.set_title('Pole figure for HKL='+str(SHData['PFhkl']))
1717	Plot.set_xlabel(G2frame.Projection.capitalize()+' projection')
1718	Page.canvas.draw()
1719
1720	################################################################################
1721	##### PlotCovariance
1722	################################################################################
1723
1724	def PlotCovariance(G2frame,Data):
1725	'needs a doc string'
1726	if not Data:
1727	print 'No covariance matrix available'
1728	return
1729	varyList = Data['varyList']
1730	values = Data['variables']
1731	Xmax = len(varyList)
1732	covMatrix = Data['covMatrix']
1733	sig = np.sqrt(np.diag(covMatrix))
1734	xvar = np.outer(sig,np.ones_like(sig))
1735	covArray = np.divide(np.divide(covMatrix,xvar),xvar.T)
1736	title = ' for\n'+Data['title']
1737	newAtomDict = Data['newAtomDict']
1738
1739
1740	def OnPlotKeyPress(event):
1741	newPlot = False
1742	if event.key == 's':
1743	choice = [m for m in mpl.cm.datad.keys() if not m.endswith("_r")]
1744	choice.sort()
1745	dlg = wx.SingleChoiceDialog(G2frame,'Select','Color scheme',choice)
1746	if dlg.ShowModal() == wx.ID_OK:
1747	sel = dlg.GetSelection()
1748	G2frame.VcovColor = choice[sel]
1749	else:
1750	G2frame.VcovColor = 'RdYlGn'
1751	dlg.Destroy()
1752	PlotCovariance(G2frame,Data)
1753
1754	def OnMotion(event):
1755	#there is a problem here - reports wrong values
1756	if event.button:
1757	ytics = imgAx.get_yticks()
1758	ytics = np.where(ytics<len(varyList),ytics,-1)
1759	ylabs = [np.where(0<=i ,varyList[int(i)],' ') for i in ytics]
1760	imgAx.set_yticklabels(ylabs)
1761	if event.xdata and event.ydata: #avoid out of frame errors
1762	xpos = int(event.xdata+.5)
1763	ypos = int(event.ydata+.5)
1764	if -1 < xpos < len(varyList) and -1 < ypos < len(varyList):
1765	if xpos == ypos:
1766	value = values[xpos]
1767	name = varyList[xpos]
1768	if varyList[xpos] in newAtomDict:
1769	name,value = newAtomDict[name]
1770	msg = '%s value = %.4g, esd = %.4g'%(name,value,sig[xpos])
1771	else:
1772	msg = '%s - %s: %5.3f'%(varyList[xpos],varyList[ypos],covArray[xpos][ypos])
1773	Page.canvas.SetToolTipString(msg)
1774	G2frame.G2plotNB.status.SetFields(['',msg])
1775
1776	try:
1777	plotNum = G2frame.G2plotNB.plotList.index('Covariance')
1778	Page = G2frame.G2plotNB.nb.GetPage(plotNum)
1779	Page.figure.clf()
1780	Plot = Page.figure.gca()
1781	if not Page.IsShown():
1782	Page.Show()
1783	except ValueError:
1784	Plot = G2frame.G2plotNB.addMpl('Covariance').gca()
1785	plotNum = G2frame.G2plotNB.plotList.index('Covariance')
1786	Page = G2frame.G2plotNB.nb.GetPage(plotNum)
1787	Page.canvas.mpl_connect('motion_notify_event', OnMotion)
1788	Page.canvas.mpl_connect('key_press_event', OnPlotKeyPress)
1789	Page.Choice = ['s: to change colors']
1790	Page.keyPress = OnPlotKeyPress
1791	Page.SetFocus()
1792	G2frame.G2plotNB.status.SetFields(['',''])
1793	acolor = mpl.cm.get_cmap(G2frame.VcovColor)
1794	Img = Plot.imshow(covArray,aspect='equal',cmap=acolor,interpolation='nearest',origin='lower',
1795	vmin=-1.,vmax=1.)
1796	imgAx = Img.get_axes()
1797	ytics = imgAx.get_yticks()
1798	ylabs = [varyList[int(i)] for i in ytics[:-1]]
1799	imgAx.set_yticklabels(ylabs)
1800	colorBar = Page.figure.colorbar(Img)
1801	Plot.set_title('V-Cov matrix'+title)
1802	Plot.set_xlabel('Variable number')
1803	Plot.set_ylabel('Variable name')
1804	Page.canvas.draw()
1805
1806	################################################################################
1807	##### PlotTorsion
1808	################################################################################
1809
1810	def PlotTorsion(G2frame,phaseName,Torsion,TorName,Names=[],Angles=[],Coeff=[]):
1811	'needs a doc string'
1812
1813	global names
1814	names = Names
1815	sum = np.sum(Torsion)
1816	torsion = np.log(2*Torsion+1.)/sum
1817	tMin = np.min(torsion)
1818	tMax = np.max(torsion)
1819	torsion = 3.*(torsion-tMin)/(tMax-tMin)
1820	X = np.linspace(0.,360.,num=45)
1821
1822	def OnMotion(event):
1823	if event.xdata and event.ydata: #avoid out of frame errors
1824	xpos = event.xdata
1825	ypos = event.ydata
1826	msg = 'torsion,energy: %5.3f %5.3f'%(xpos,ypos)
1827	Page.canvas.SetToolTipString(msg)
1828	G2frame.G2plotNB.status.SetFields(['',msg])
1829	if len(Angles):
1830	fX = np.where(np.fabs(Angles-xpos)<1.0)[0]
1831	if len(fX):
1832	msg = 'atoms:'+names[fX[0]-1]
1833	Page.canvas.SetToolTipString(msg)
1834
1835	try:
1836	plotNum = G2frame.G2plotNB.plotList.index('Torsion')
1837	Page = G2frame.G2plotNB.nb.GetPage(plotNum)
1838	Page.figure.clf()
1839	Plot = Page.figure.gca()
1840	if not Page.IsShown():
1841	Page.Show()
1842	except ValueError:
1843	Plot = G2frame.G2plotNB.addMpl('Torsion').gca()
1844	plotNum = G2frame.G2plotNB.plotList.index('Torsion')
1845	Page = G2frame.G2plotNB.nb.GetPage(plotNum)
1846	Page.canvas.mpl_connect('motion_notify_event', OnMotion)
1847
1848	Page.SetFocus()
1849	G2frame.G2plotNB.status.SetFields(['',''])
1850	Plot.plot(X,torsion,'b+')
1851	if len(Coeff):
1852	X2 = np.linspace(0.,360.,45)
1853	Y2 = np.array([-G2mth.calcTorsionEnergy(x,Coeff)[1] for x in X2])
1854	Plot.plot(X2,Y2,'r')
1855	if len(Angles):
1856	Eval = np.array([-G2mth.calcTorsionEnergy(x,Coeff)[1] for x in Angles])
1857	Plot.plot(Angles,Eval,'ro')
1858	Plot.set_xlim((0.,360.))
1859	Plot.set_title('Torsion angles for '+TorName+' in '+phaseName)
1860	Plot.set_xlabel('angle',fontsize=16)
1861	Plot.set_ylabel('Energy',fontsize=16)
1862	Page.canvas.draw()
1863
1864	################################################################################
1865	##### PlotRama
1866	################################################################################
1867
1868	def PlotRama(G2frame,phaseName,Rama,RamaName,Names=[],PhiPsi=[],Coeff=[]):
1869	'needs a doc string'
1870
1871	global names
1872	names = Names
1873	rama = np.log(2*Rama+1.)
1874	ramaMax = np.max(rama)
1875	rama = np.reshape(rama,(45,45))
1876	global Phi,Psi
1877	Phi = []
1878	Psi = []
1879
1880	def OnPlotKeyPress(event):
1881	newPlot = False
1882	if event.key == 's':
1883	choice = [m for m in mpl.cm.datad.keys() if not m.endswith("_r")]
1884	choice.sort()
1885	dlg = wx.SingleChoiceDialog(G2frame,'Select','Color scheme',choice)
1886	if dlg.ShowModal() == wx.ID_OK:
1887	sel = dlg.GetSelection()
1888	G2frame.RamaColor = choice[sel]
1889	else:
1890	G2frame.RamaColor = 'RdYlGn'
1891	dlg.Destroy()
1892	PlotRama(G2frame,phaseName,Rama)
1893
1894	def OnMotion(event):
1895	if event.xdata and event.ydata: #avoid out of frame errors
1896	xpos = event.xdata
1897	ypos = event.ydata
1898	msg = 'phi/psi: %5.3f %5.3f'%(xpos,ypos)
1899	Page.canvas.SetToolTipString(msg)
1900	G2frame.G2plotNB.status.SetFields(['',msg])
1901	if len(Phi):
1902	fPhi = np.where(np.fabs(Phi-xpos)<1.0)[0]
1903	fPsi = np.where(np.fabs(Psi-ypos)<1.0)[0]
1904	if len(fPhi) and len(fPsi) and fPhi[0] == fPsi[0]:
1905	msg = 'atoms:'+names[fPhi[0]-1]
1906	Page.canvas.SetToolTipString(msg)
1907
1908	try:
1909	plotNum = G2frame.G2plotNB.plotList.index('Ramachandran')
1910	Page = G2frame.G2plotNB.nb.GetPage(plotNum)
1911	Page.figure.clf()
1912	Plot = Page.figure.gca()
1913	if not Page.IsShown():
1914	Page.Show()
1915	except ValueError:
1916	Plot = G2frame.G2plotNB.addMpl('Ramachandran').gca()
1917	plotNum = G2frame.G2plotNB.plotList.index('Ramachandran')
1918	Page = G2frame.G2plotNB.nb.GetPage(plotNum)
1919	Page.canvas.mpl_connect('motion_notify_event', OnMotion)
1920	Page.canvas.mpl_connect('key_press_event', OnPlotKeyPress)
1921
1922	Page.Choice = ['s: to change colors']
1923	Page.keyPress = OnPlotKeyPress
1924	Page.SetFocus()
1925	G2frame.G2plotNB.status.SetFields(['',''])
1926	acolor = mpl.cm.get_cmap(G2frame.RamaColor)
1927	if RamaName == 'All' or '-1' in RamaName:
1928	if len(Coeff):
1929	X,Y = np.meshgrid(np.linspace(-180.,180.,45),np.linspace(-180.,180.,45))
1930	Z = np.array([-G2mth.calcRamaEnergy(x,y,Coeff)[1] for x,y in zip(X.flatten(),Y.flatten())])
1931	Plot.contour(X,Y,np.reshape(Z,(45,45)))
1932	Img = Plot.imshow(rama,aspect='equal',cmap=acolor,interpolation='nearest',
1933	extent=[-180,180,-180,180],origin='lower')
1934	if len(PhiPsi):
1935	Phi,Psi = PhiPsi.T
1936	Phi = np.where(Phi>180.,Phi-360.,Phi)
1937	Psi = np.where(Psi>180.,Psi-360.,Psi)
1938	Plot.plot(Phi,Psi,'ro',picker=3)
1939	Plot.set_xlim((-180.,180.))
1940	Plot.set_ylim((-180.,180.))
1941	else:
1942	if len(Coeff):
1943	X,Y = np.meshgrid(np.linspace(0.,360.,45),np.linspace(0.,360.,45))
1944	Z = np.array([-G2mth.calcRamaEnergy(x,y,Coeff)[1] for x,y in zip(X.flatten(),Y.flatten())])
1945	Plot.contour(X,Y,np.reshape(Z,(45,45)))
1946	Img = Plot.imshow(rama,aspect='equal',cmap=acolor,interpolation='nearest',
1947	extent=[0,360,0,360],origin='lower')
1948	if len(PhiPsi):
1949	Phi,Psi = PhiPsi.T
1950	Plot.plot(Phi,Psi,'ro',picker=3)
1951	Plot.set_xlim((0.,360.))
1952	Plot.set_ylim((0.,360.))
1953	Plot.set_title('Ramachandran for '+RamaName+' in '+phaseName)
1954	Plot.set_xlabel(r'$\phi$',fontsize=16)
1955	Plot.set_ylabel(r'$\psi$',fontsize=16)
1956	colorBar = Page.figure.colorbar(Img)
1957	Page.canvas.draw()
1958
1959
1960	################################################################################
1961	##### PlotSeq
1962	################################################################################
1963
1964	def PlotSeq(G2frame,SeqData,SeqSig,SeqNames,sampleParm):
1965	'needs a doc string'
1966
1967	def OnKeyPress(event):
1968	if event.key == 's' and sampleParm:
1969	if G2frame.xAxis:
1970	G2frame.xAxis = False
1971	else:
1972	G2frame.xAxis = True
1973	Draw(False)
1974	try:
1975	plotNum = G2frame.G2plotNB.plotList.index('Sequential refinement')
1976	Page = G2frame.G2plotNB.nb.GetPage(plotNum)
1977	Page.figure.clf()
1978	Plot = Page.figure.gca()
1979	if not Page.IsShown():
1980	Page.Show()
1981	except ValueError:
1982	Plot = G2frame.G2plotNB.addMpl('Sequential refinement').gca()
1983	plotNum = G2frame.G2plotNB.plotList.index('Sequential refinement')
1984	Page = G2frame.G2plotNB.nb.GetPage(plotNum)
1985	Page.canvas.mpl_connect('key_press_event', OnKeyPress)
1986	G2frame.xAxis = False
1987	Page.Choice = ['s to toggle x-axis = sample environment parameter']
1988	Page.keyPress = OnKeyPress
1989
1990	def Draw(newPlot):
1991	Page.SetFocus()
1992	G2frame.G2plotNB.status.SetFields(['','press '])
1993	if len(SeqData):
1994	Plot.clear()
1995	if G2frame.xAxis:
1996	xName = sampleParm.keys()[0]
1997	X = sampleParm[xName]
1998	else:
1999	X = np.arange(0,len(SeqData[0]),1)
2000	xName = 'Data sequence number'
2001	for Y,sig,name in zip(SeqData,SeqSig,SeqNames):
2002	Plot.errorbar(X,Y,yerr=sig,label=name)
2003	Plot.legend(loc='best')
2004	Plot.set_ylabel('Parameter values')
2005	Plot.set_xlabel(xName)
2006	Page.canvas.draw()
2007	Draw(True)
2008
2009	################################################################################
2010	##### PlotExposedImage & PlotImage
2011	################################################################################
2012
2013	def PlotExposedImage(G2frame,newPlot=False,event=None):
2014	'''General access module for 2D image plotting
2015	'''
2016	plotNo = G2frame.G2plotNB.nb.GetSelection()
2017	if G2frame.G2plotNB.nb.GetPageText(plotNo) == '2D Powder Image':
2018	PlotImage(G2frame,newPlot,event,newImage=True)
2019	elif G2frame.G2plotNB.nb.GetPageText(plotNo) == '2D Integration':
2020	PlotIntegration(G2frame,newPlot,event)
2021
2022	def PlotImage(G2frame,newPlot=False,event=None,newImage=True):
2023	'''Plot of 2D detector images as contoured plot. Also plot calibration ellipses,
2024	masks, etc.
2025	'''
2026	from matplotlib.patches import Ellipse,Arc,Circle,Polygon
2027	import numpy.ma as ma
2028	Dsp = lambda tth,wave: wave/(2.*sind(tth/2.))
2029	global Data,Masks
2030	colors=['b','g','r','c','m','k']
2031	Data = G2frame.PatternTree.GetItemPyData(
2032	G2gd.GetPatternTreeItemId(G2frame,G2frame.Image, 'Image Controls'))
2033	Masks = G2frame.PatternTree.GetItemPyData(
2034	G2gd.GetPatternTreeItemId(G2frame,G2frame.Image, 'Masks'))
2035	try: #may be absent
2036	StrSta = G2frame.PatternTree.GetItemPyData(
2037	G2gd.GetPatternTreeItemId(G2frame,G2frame.Image, 'Stress/Strain'))
2038	except TypeError: #is missing
2039	StrSta = {}
2040
2041	def OnImMotion(event):
2042	Page.canvas.SetToolTipString('')
2043	sizexy = Data['size']
2044	if event.xdata and event.ydata and len(G2frame.ImageZ): #avoid out of frame errors
2045	Page.canvas.SetCursor(wx.CROSS_CURSOR)
2046	item = G2frame.itemPicked
2047	pixelSize = Data['pixelSize']
2048	scalex = 1000./pixelSize[0]
2049	scaley = 1000./pixelSize[1]
2050	if item and G2frame.PatternTree.GetItemText(G2frame.PickId) == 'Image Controls':
2051	if 'Text' in str(item):
2052	Page.canvas.SetToolTipString('%8.3f %8.3fmm'%(event.xdata,event.ydata))
2053	else:
2054	xcent,ycent = Data['center']
2055	xpos = event.xdata-xcent
2056	ypos = event.ydata-ycent
2057	tth,azm = G2img.GetTthAzm(event.xdata,event.ydata,Data)
2058	if 'line3' in str(item) or 'line4' in str(item) and not Data['fullIntegrate']:
2059	Page.canvas.SetToolTipString('%6d deg'%(azm))
2060	elif 'line1' in str(item) or 'line2' in str(item):
2061	Page.canvas.SetToolTipString('%8.3fdeg'%(tth))
2062	else:
2063	xpos = event.xdata
2064	ypos = event.ydata
2065	xpix = xpos*scalex
2066	ypix = ypos*scaley
2067	Int = 0
2068	if (0 <= xpix <= sizexy[0]) and (0 <= ypix <= sizexy[1]):
2069	Int = G2frame.ImageZ[ypix][xpix]
2070	tth,azm,dsp = G2img.GetTthAzmDsp(xpos,ypos,Data)
2071	Q = 2.*math.pi/dsp
2072	if G2frame.setPoly:
2073	G2frame.G2plotNB.status.SetFields(['','Polygon mask pick - LB next point, RB close polygon'])
2074	else:
2075	G2frame.G2plotNB.status.SetFields(\
2076	['','Detector 2-th =%9.3fdeg, dsp =%9.3fA, Q = %6.5fA-1, azm = %7.2fdeg, I = %6d'%(tth,dsp,Q,azm,Int)])
2077
2078	def OnImPlotKeyPress(event):
2079	try:
2080	PickName = G2frame.PatternTree.GetItemText(G2frame.PickId)
2081	except TypeError:
2082	return
2083	if PickName == 'Masks':
2084	Xpos = event.xdata
2085	if not Xpos: #got point out of frame
2086	return
2087	Ypos = event.ydata
2088	if event.key == 's':
2089	Masks['Points'].append([Xpos,Ypos,1.])
2090	elif event.key == 'r':
2091	tth = G2img.GetTth(Xpos,Ypos,Data)
2092	Masks['Rings'].append([tth,0.1])
2093	elif event.key == 'a':
2094	tth,azm = G2img.GetTthAzm(Xpos,Ypos,Data)
2095	azm = int(azm)
2096	Masks['Arcs'].append([tth,[azm-5,azm+5],0.1])
2097	elif event.key == 'p':
2098	G2frame.setPoly = True
2099	Masks['Polygons'].append([])
2100	G2frame.G2plotNB.status.SetFields(['','Polygon mask active - LB pick next point, RB close polygon'])
2101	G2imG.UpdateMasks(G2frame,Masks)
2102	elif PickName == 'Image Controls':
2103	if event.key == 'c':
2104	Xpos = event.xdata
2105	if not Xpos: #got point out of frame
2106	return
2107	Ypos = event.ydata
2108	dlg = wx.MessageDialog(G2frame,'Are you sure you want to change the center?',
2109	'Center change',style=wx.OK\|wx.CANCEL)
2110	try:
2111	if dlg.ShowModal() == wx.ID_OK:
2112	print 'move center to: ',Xpos,Ypos
2113	Data['center'] = [Xpos,Ypos]
2114	G2imG.UpdateImageControls(G2frame,Data,Masks)
2115	finally:
2116	dlg.Destroy()
2117	elif event.key == 'l':
2118	if G2frame.logPlot:
2119	G2frame.logPlot = False
2120	else:
2121	G2frame.logPlot = True
2122	PlotImage(G2frame,newImage=True)
2123
2124	def OnKeyBox(event):
2125	if G2frame.G2plotNB.nb.GetSelection() == G2frame.G2plotNB.plotList.index('2D Powder Image'):
2126	event.key = cb.GetValue()[0]
2127	cb.SetValue(' key press')
2128	if event.key in 'l':
2129	wx.CallAfter(OnImPlotKeyPress,event)
2130
2131	def OnImPick(event):
2132	if G2frame.PatternTree.GetItemText(G2frame.PickId) not in ['Image Controls','Masks']:
2133	return
2134	if G2frame.setPoly:
2135	polygon = Masks['Polygons'][-1]
2136	xpos,ypos = event.mouseevent.xdata,event.mouseevent.ydata
2137	if xpos and ypos: #point inside image
2138	if len(polygon) > 2 and event.mouseevent.button == 3:
2139	x0,y0 = polygon[0]
2140	polygon.append([x0,y0])
2141	G2frame.setPoly = False
2142	G2frame.G2plotNB.status.SetFields(['','Polygon closed - RB drag a vertex to change shape'])
2143	else:
2144	G2frame.G2plotNB.status.SetFields(['','New polygon point: %.1f,%.1f'%(xpos,ypos)])
2145	polygon.append([xpos,ypos])
2146	G2imG.UpdateMasks(G2frame,Masks)
2147	else:
2148	if G2frame.itemPicked is not None: return
2149	G2frame.itemPicked = event.artist
2150	G2frame.mousePicked = event.mouseevent
2151
2152	def OnImRelease(event):
2153	try:
2154	PickName = G2frame.PatternTree.GetItemText(G2frame.PickId)
2155	except TypeError:
2156	return
2157	if PickName not in ['Image Controls','Masks']:
2158	return
2159	pixelSize = Data['pixelSize']
2160	scalex = 1000./pixelSize[0]
2161	scaley = 1000./pixelSize[1]
2162	pixLimit = Data['pixLimit']
2163	if G2frame.itemPicked is None and PickName == 'Image Controls' and len(G2frame.ImageZ):
2164	# sizexy = Data['size']
2165	Xpos = event.xdata
2166	if not (Xpos and G2frame.ifGetRing): #got point out of frame
2167	return
2168	Ypos = event.ydata
2169	if Ypos and not Page.toolbar._active: #make sure zoom/pan not selected
2170	if event.button == 1:
2171	Xpix = Xpos*scalex
2172	Ypix = Ypos*scaley
2173	xpos,ypos,I,J = G2img.ImageLocalMax(G2frame.ImageZ,pixLimit,Xpix,Ypix)
2174	if I and J:
2175	xpos += .5 #shift to pixel center
2176	ypos += .5
2177	xpos /= scalex #convert to mm
2178	ypos /= scaley
2179	Data['ring'].append([xpos,ypos])
2180	elif event.button == 3:
2181	G2frame.dataFrame.GetStatusBar().SetStatusText('Calibrating...')
2182	if G2img.ImageCalibrate(G2frame,Data):
2183	G2frame.dataFrame.GetStatusBar().SetStatusText('Calibration successful - Show ring picks to check')
2184	print 'Calibration successful'
2185	else:
2186	G2frame.dataFrame.GetStatusBar().SetStatusText('Calibration failed - Show ring picks to diagnose')
2187	print 'Calibration failed'
2188	G2frame.ifGetRing = False
2189	G2imG.UpdateImageControls(G2frame,Data,Masks)
2190	return
2191	PlotImage(G2frame,newImage=False)
2192	return
2193	else:
2194	xpos = event.xdata
2195	if xpos: #avoid out of frame mouse position
2196	ypos = event.ydata
2197	if G2frame.ifGetRing: #delete a calibration ring pick
2198	xypos = [xpos,ypos]
2199	rings = Data['ring']
2200	for ring in rings:
2201	if np.allclose(ring,xypos,.01,0):
2202	rings.remove(ring)
2203	else:
2204	tth,azm,dsp = G2img.GetTthAzmDsp(xpos,ypos,Data)
2205	itemPicked = str(G2frame.itemPicked)
2206	if 'Line2D' in itemPicked and PickName == 'Image Controls':
2207	if 'line1' in itemPicked:
2208	Data['IOtth'][0] = max(tth,0.001)
2209	elif 'line2' in itemPicked:
2210	Data['IOtth'][1] = tth
2211	elif 'line3' in itemPicked:
2212	Data['LRazimuth'][0] = int(azm)
2213	if Data['fullIntegrate']:
2214	Data['LRazimuth'][1] = Data['LRazimuth'][0]+360
2215	elif 'line4' in itemPicked and not Data['fullIntegrate']:
2216	Data['LRazimuth'][1] = int(azm)
2217
2218	if Data['LRazimuth'][0] > Data['LRazimuth'][1]:
2219	Data['LRazimuth'][0] -= 360
2220
2221	azLim = np.array(Data['LRazimuth'])
2222	if np.any(azLim>360):
2223	azLim -= 360
2224	Data['LRazimuth'] = list(azLim)
2225
2226	if Data['IOtth'][0] > Data['IOtth'][1]:
2227	Data['IOtth'][0],Data['IOtth'][1] = Data['IOtth'][1],Data['IOtth'][0]
2228
2229	G2frame.InnerTth.SetValue("%8.2f" % (Data['IOtth'][0]))
2230	G2frame.OuterTth.SetValue("%8.2f" % (Data['IOtth'][1]))
2231	G2frame.Lazim.SetValue("%6d" % (Data['LRazimuth'][0]))
2232	G2frame.Razim.SetValue("%6d" % (Data['LRazimuth'][1]))
2233	elif 'Circle' in itemPicked and PickName == 'Masks':
2234	spots = Masks['Points']
2235	newPos = itemPicked.split(')')[0].split('(')[2].split(',')
2236	newPos = np.array([float(newPos[0]),float(newPos[1])])
2237	for spot in spots:
2238	if np.allclose(np.array([spot[:2]]),newPos):
2239	spot[:2] = xpos,ypos
2240	G2imG.UpdateMasks(G2frame,Masks)
2241	elif 'Line2D' in itemPicked and PickName == 'Masks':
2242	Obj = G2frame.itemPicked.findobj()
2243	rings = Masks['Rings']
2244	arcs = Masks['Arcs']
2245	polygons = Masks['Polygons']
2246	for ring in G2frame.ringList:
2247	if Obj == ring[0]:
2248	rN = ring[1]
2249	if ring[2] == 'o':
2250	rings[rN][0] = G2img.GetTth(xpos,ypos,Data)-rings[rN][1]/2.
2251	else:
2252	rings[rN][0] = G2img.GetTth(xpos,ypos,Data)+rings[rN][1]/2.
2253	for arc in G2frame.arcList:
2254	if Obj == arc[0]:
2255	aN = arc[1]
2256	if arc[2] == 'o':
2257	arcs[aN][0] = G2img.GetTth(xpos,ypos,Data)-arcs[aN][2]/2
2258	elif arc[2] == 'i':
2259	arcs[aN][0] = G2img.GetTth(xpos,ypos,Data)+arcs[aN][2]/2
2260	elif arc[2] == 'l':
2261	arcs[aN][1][0] = int(G2img.GetAzm(xpos,ypos,Data))
2262	else:
2263	arcs[aN][1][1] = int(G2img.GetAzm(xpos,ypos,Data))
2264	for poly in G2frame.polyList:
2265	if Obj == poly[0]:
2266	ind = G2frame.itemPicked.contains(G2frame.mousePicked)[1]['ind'][0]
2267	oldPos = np.array([G2frame.mousePicked.xdata,G2frame.mousePicked.ydata])
2268	pN = poly[1]
2269	for i,xy in enumerate(polygons[pN]):
2270	if np.allclose(np.array([xy]),oldPos,atol=1.0):
2271	polygons[pN][i] = xpos,ypos
2272	G2imG.UpdateMasks(G2frame,Masks)
2273	# else: #keep for future debugging
2274	# print str(G2frame.itemPicked),event.xdata,event.ydata,event.button
2275	PlotImage(G2frame,newImage=True)
2276	G2frame.itemPicked = None
2277
2278	try:
2279	plotNum = G2frame.G2plotNB.plotList.index('2D Powder Image')
2280	Page = G2frame.G2plotNB.nb.GetPage(plotNum)
2281	if not newPlot:
2282	Plot = Page.figure.gca() #get previous powder plot & get limits
2283	xylim = Plot.get_xlim(),Plot.get_ylim()
2284	if newImage:
2285	Page.figure.clf()
2286	Plot = Page.figure.gca() #get a fresh plot after clf()
2287
2288	except ValueError:
2289	Plot = G2frame.G2plotNB.addMpl('2D Powder Image').gca()
2290	plotNum = G2frame.G2plotNB.plotList.index('2D Powder Image')
2291	Page = G2frame.G2plotNB.nb.GetPage(plotNum)
2292	Page.canvas.mpl_connect('key_press_event', OnImPlotKeyPress)
2293	Page.canvas.mpl_connect('motion_notify_event', OnImMotion)
2294	Page.canvas.mpl_connect('pick_event', OnImPick)
2295	Page.canvas.mpl_connect('button_release_event', OnImRelease)
2296	xylim = []
2297	Page.Choice = None
2298	if not event: #event from GUI TextCtrl - don't want focus to change to plot!!!
2299	Page.SetFocus()
2300	Title = G2frame.PatternTree.GetItemText(G2frame.Image)[4:]
2301	G2frame.G2plotNB.status.DestroyChildren()
2302	if G2frame.logPlot:
2303	Title = 'log('+Title+')'
2304	Plot.set_title(Title)
2305	try:
2306	if G2frame.PatternTree.GetItemText(G2frame.PickId) in ['Image Controls',]:
2307	Page.Choice = (' key press','l: log(I) on',)
2308	if G2frame.logPlot:
2309	Page.Choice[1] = 'l: log(I) off'
2310	Page.keyPress = OnImPlotKeyPress
2311	except TypeError:
2312	pass
2313	size,imagefile = G2frame.PatternTree.GetItemPyData(G2frame.Image)
2314	if imagefile != G2frame.oldImagefile:
2315	imagefile = G2IO.CheckImageFile(G2frame,imagefile)
2316	if not imagefile:
2317	G2frame.G2plotNB.Delete('2D Powder Image')
2318	return
2319	G2frame.PatternTree.SetItemPyData(G2frame.Image,[size,imagefile])
2320	G2frame.ImageZ = G2IO.GetImageData(G2frame,imagefile,imageOnly=True)
2321	G2frame.oldImagefile = imagefile
2322
2323	imScale = 1
2324	if len(G2frame.ImageZ) > 1024:
2325	imScale = len(G2frame.ImageZ)/1024
2326	sizexy = Data['size']
2327	pixelSize = Data['pixelSize']
2328	scalex = 1000./pixelSize[0]
2329	scaley = 1000./pixelSize[1]
2330	Xmax = sizexy[0]*pixelSize[0]/1000.
2331	Ymax = sizexy[1]*pixelSize[1]/1000.
2332	xlim = (0,Xmax)
2333	ylim = (Ymax,0)
2334	Imin,Imax = Data['range'][1]
2335	acolor = mpl.cm.get_cmap(Data['color'])
2336	xcent,ycent = Data['center']
2337	Plot.set_xlabel('Image x-axis, mm',fontsize=12)
2338	Plot.set_ylabel('Image y-axis, mm',fontsize=12)
2339	#do threshold mask - "real" mask - others are just bondaries
2340	Zlim = Masks['Thresholds'][1]
2341	wx.BeginBusyCursor()
2342	try:
2343
2344	if newImage:
2345	MA = ma.masked_greater(ma.masked_less(G2frame.ImageZ,Zlim[0]),Zlim[1])
2346	MaskA = ma.getmaskarray(MA)
2347	A = G2img.ImageCompress(MA,imScale)
2348	AM = G2img.ImageCompress(MaskA,imScale)
2349	if G2frame.logPlot:
2350	A = np.where(A>0,np.log(A),0)
2351	AM = np.where(AM>0,np.log(AM),0)
2352	Imin,Imax = [np.amin(A),np.amax(A)]
2353	ImgM = Plot.imshow(AM,aspect='equal',cmap='Reds',
2354	interpolation='nearest',vmin=0,vmax=2,extent=[0,Xmax,Ymax,0])
2355	Img = Plot.imshow(A,aspect='equal',cmap=acolor,
2356	interpolation='nearest',vmin=Imin,vmax=Imax,extent=[0,Xmax,Ymax,0])
2357	if G2frame.setPoly:
2358	Img.set_picker(True)
2359
2360	Plot.plot(xcent,ycent,'x')
2361	#G2frame.PatternTree.GetItemText(item)
2362	if Data['showLines']:
2363	LRAzim = Data['LRazimuth'] #NB: integers
2364	Nazm = Data['outAzimuths']
2365	delAzm = float(LRAzim[1]-LRAzim[0])/Nazm
2366	AzmthOff = Data['azmthOff']
2367	IOtth = Data['IOtth']
2368	wave = Data['wavelength']
2369	dspI = wave/(2.0*sind(IOtth[0]/2.0))
2370	ellI = G2img.GetEllipse(dspI,Data) #=False if dsp didn't yield an ellipse (ugh! a parabola or a hyperbola)
2371	dspO = wave/(2.0*sind(IOtth[1]/2.0))
2372	ellO = G2img.GetEllipse(dspO,Data) #Ditto & more likely for outer ellipse
2373	Azm = np.array(range(LRAzim[0],LRAzim[1]+1))-AzmthOff
2374	if ellI:
2375	xyI = []
2376	for azm in Azm:
2377	xyI.append(G2img.GetDetectorXY(dspI,azm-90.,Data))
2378	xyI = np.array(xyI)
2379	arcxI,arcyI = xyI.T
2380	Plot.plot(arcxI,arcyI,picker=3)
2381	if ellO:
2382	xyO = []
2383	for azm in Azm:
2384	xyO.append(G2img.GetDetectorXY(dspO,azm-90.,Data))
2385	xyO = np.array(xyO)
2386	arcxO,arcyO = xyO.T
2387	Plot.plot(arcxO,arcyO,picker=3)
2388	if ellO and ellI:
2389	Plot.plot([arcxI[0],arcxO[0]],[arcyI[0],arcyO[0]],picker=3)
2390	Plot.plot([arcxI[-1],arcxO[-1]],[arcyI[-1],arcyO[-1]],picker=3)
2391	for i in range(Nazm):
2392	cake = LRAzim[0]+i*delAzm-AzmthOff
2393	if Data['centerAzm']:
2394	cake += delAzm/2.
2395	ind = np.searchsorted(Azm,cake)
2396	Plot.plot([arcxI[ind],arcxO[ind]],[arcyI[ind],arcyO[ind]],color='k',dashes=(5,5))
2397
2398	if G2frame.PatternTree.GetItemText(G2frame.PickId) in 'Image Controls':
2399	for xring,yring in Data['ring']:
2400	Plot.plot(xring,yring,'r+',picker=3)
2401	if Data['setRings']:
2402	# rings = np.concatenate((Data['rings']),axis=0)
2403	N = 0
2404	for ring in Data['rings']:
2405	xring,yring = np.array(ring).T[:2]
2406	Plot.plot(xring,yring,'+',color=colors[N%6])
2407	N += 1
2408	for ellipse in Data['ellipses']:
2409	cent,phi,[width,height],col = ellipse
2410	Plot.add_artist(Ellipse([cent[0],cent[1]],2width,2height,phi,ec=col,fc='none'))
2411	Plot.text(cent[0],cent[1],'+',color=col,ha='center',va='center')
2412	if G2frame.PatternTree.GetItemText(G2frame.PickId) in 'Stress/Strain':
2413	print 'plot stress/strain stuff'
2414	for ring in StrSta['d-zero']:
2415	xring,yring = ring['ImxyObs']
2416	Plot.plot(xring,yring,'r+')
2417	# for xring,yring in ring['ImxyCalc'].T:
2418	# Plot.add_artist(Polygon(ring['ImxyCalc'].T,ec='b',fc='none'))
2419	# Plot.plot(xring,yring)
2420	#masks - mask lines numbered after integration limit lines
2421	spots = Masks['Points']
2422	rings = Masks['Rings']
2423	arcs = Masks['Arcs']
2424	polygons = Masks['Polygons']
2425	for x,y,d in spots:
2426	Plot.add_artist(Circle((x,y),radius=d/2,fc='none',ec='r',picker=3))
2427	G2frame.ringList = []
2428	for iring,(tth,thick) in enumerate(rings):
2429	wave = Data['wavelength']
2430	x1,y1 = np.hsplit(np.array(G2img.makeIdealRing(G2img.GetEllipse(Dsp(tth+thick/2.,wave),Data))),2)
2431	x2,y2 = np.hsplit(np.array(G2img.makeIdealRing(G2img.GetEllipse(Dsp(tth-thick/2.,wave),Data))),2)
2432	G2frame.ringList.append([Plot.plot(x1,y1,'r',picker=3),iring,'o'])
2433	G2frame.ringList.append([Plot.plot(x2,y2,'r',picker=3),iring,'i'])
2434	G2frame.arcList = []
2435	for iarc,(tth,azm,thick) in enumerate(arcs):
2436	wave = Data['wavelength']
2437	x1,y1 = np.hsplit(np.array(G2img.makeIdealRing(G2img.GetEllipse(Dsp(tth+thick/2.,wave),Data),azm)),2)
2438	x2,y2 = np.hsplit(np.array(G2img.makeIdealRing(G2img.GetEllipse(Dsp(max(0.01,tth-thick/2.),wave),Data),azm)),2)
2439	G2frame.arcList.append([Plot.plot(x1,y1,'r',picker=3),iarc,'o'])
2440	G2frame.arcList.append([Plot.plot(x2,y2,'r',picker=3),iarc,'i'])
2441	G2frame.arcList.append([Plot.plot([x1[0],x2[0]],[y1[0],y2[0]],'r',picker=3),iarc,'l'])
2442	G2frame.arcList.append([Plot.plot([x1[-1],x2[-1]],[y1[-1],y2[-1]],'r',picker=3),iarc,'u'])
2443	G2frame.polyList = []
2444	for ipoly,polygon in enumerate(polygons):
2445	x,y = np.hsplit(np.array(polygon),2)
2446	G2frame.polyList.append([Plot.plot(x,y,'r+',picker=10),ipoly])
2447	Plot.plot(x,y,'r')
2448	if newImage:
2449	colorBar = Page.figure.colorbar(Img)
2450	Plot.set_xlim(xlim)
2451	Plot.set_ylim(ylim)
2452	Plot.invert_yaxis()
2453	if not newPlot and xylim:
2454	Page.toolbar.push_current()
2455	Plot.set_xlim(xylim[0])
2456	Plot.set_ylim(xylim[1])
2457	xylim = []
2458	Page.toolbar.push_current()
2459	Page.toolbar.draw()
2460	else:
2461	Page.canvas.draw()
2462	finally:
2463	wx.EndBusyCursor()
2464
2465	################################################################################
2466	##### PlotIntegration
2467	################################################################################
2468
2469	def PlotIntegration(G2frame,newPlot=False,event=None):
2470	'''Plot of 2D image after image integration with 2-theta and azimuth as coordinates
2471	'''
2472
2473	def OnMotion(event):
2474	Page.canvas.SetToolTipString('')
2475	Page.canvas.SetCursor(wx.CROSS_CURSOR)
2476	azm = event.ydata
2477	tth = event.xdata
2478	if azm and tth:
2479	G2frame.G2plotNB.status.SetFields(\
2480	['','Detector 2-th =%9.3fdeg, azm = %7.2fdeg'%(tth,azm)])
2481
2482	try:
2483	plotNum = G2frame.G2plotNB.plotList.index('2D Integration')
2484	Page = G2frame.G2plotNB.nb.GetPage(plotNum)
2485	if not newPlot:
2486	Plot = Page.figure.gca() #get previous plot & get limits
2487	xylim = Plot.get_xlim(),Plot.get_ylim()
2488	Page.figure.clf()
2489	Plot = Page.figure.gca() #get a fresh plot after clf()
2490
2491	except ValueError:
2492	Plot = G2frame.G2plotNB.addMpl('2D Integration').gca()
2493	plotNum = G2frame.G2plotNB.plotList.index('2D Integration')
2494	Page = G2frame.G2plotNB.nb.GetPage(plotNum)
2495	Page.canvas.mpl_connect('motion_notify_event', OnMotion)
2496	Page.views = False
2497	view = False
2498	Page.Choice = None
2499	if not event:
2500	Page.SetFocus()
2501
2502	Data = G2frame.PatternTree.GetItemPyData(
2503	G2gd.GetPatternTreeItemId(G2frame,G2frame.Image, 'Image Controls'))
2504	image = G2frame.Integrate[0]
2505	xsc = G2frame.Integrate[1]
2506	ysc = G2frame.Integrate[2]
2507	Imin,Imax = Data['range'][1]
2508	acolor = mpl.cm.get_cmap(Data['color'])
2509	Plot.set_title(G2frame.PatternTree.GetItemText(G2frame.Image)[4:])
2510	Plot.set_ylabel('azimuth',fontsize=12)
2511	Plot.set_xlabel('2-theta',fontsize=12)
2512	Img = Plot.imshow(image,cmap=acolor,vmin=Imin,vmax=Imax,interpolation='nearest', \
2513	extent=[ysc[0],ysc[-1],xsc[-1],xsc[0]],aspect='auto')
2514	colorBar = Page.figure.colorbar(Img)
2515	if Data['ellipses']:
2516	for ellipse in Data['ellipses']:
2517	ring = np.array(G2img.makeIdealRing(ellipse[:3])) #skip color
2518	x,y = np.hsplit(ring,2)
2519	tth,azm = G2img.GetTthAzm(x,y,Data)
2520	# azm = np.where(azm < 0.,azm+360,azm)
2521	Plot.plot(tth,azm,'b,')
2522	if not newPlot:
2523	Page.toolbar.push_current()
2524	Plot.set_xlim(xylim[0])
2525	Plot.set_ylim(xylim[1])
2526	xylim = []
2527	Page.toolbar.push_current()
2528	Page.toolbar.draw()
2529	else:
2530	Page.canvas.draw()
2531
2532	################################################################################
2533	##### PlotTRImage
2534	################################################################################
2535
2536	def PlotTRImage(G2frame,tax,tay,taz,newPlot=False):
2537	'''a test plot routine - not normally used
2538	'''
2539
2540	def OnMotion(event):
2541	Page.canvas.SetToolTipString('')
2542	Page.canvas.SetCursor(wx.CROSS_CURSOR)
2543	azm = event.xdata
2544	tth = event.ydata
2545	if azm and tth:
2546	G2frame.G2plotNB.status.SetFields(\
2547	['','Detector 2-th =%9.3fdeg, azm = %7.2fdeg'%(tth,azm)])
2548
2549	try:
2550	plotNum = G2frame.G2plotNB.plotList.index('2D Transformed Powder Image')
2551	Page = G2frame.G2plotNB.nb.GetPage(plotNum)
2552	if not newPlot:
2553	Plot = Page.figure.gca() #get previous plot & get limits
2554	xylim = Plot.get_xlim(),Plot.get_ylim()
2555	Page.figure.clf()
2556	Plot = Page.figure.gca() #get a fresh plot after clf()
2557
2558	except ValueError:
2559	Plot = G2frame.G2plotNB.addMpl('2D Transformed Powder Image').gca()
2560	plotNum = G2frame.G2plotNB.plotList.index('2D Transformed Powder Image')
2561	Page = G2frame.G2plotNB.nb.GetPage(plotNum)
2562	Page.canvas.mpl_connect('motion_notify_event', OnMotion)
2563	Page.views = False
2564	view = False
2565	Page.Choice = None
2566	Page.SetFocus()
2567
2568	Data = G2frame.PatternTree.GetItemPyData(
2569	G2gd.GetPatternTreeItemId(G2frame,G2frame.Image, 'Image Controls'))
2570	Imin,Imax = Data['range'][1]
2571	step = (Imax-Imin)/5.
2572	V = np.arange(Imin,Imax,step)
2573	acolor = mpl.cm.get_cmap(Data['color'])
2574	Plot.set_title(G2frame.PatternTree.GetItemText(G2frame.Image)[4:])
2575	Plot.set_xlabel('azimuth',fontsize=12)
2576	Plot.set_ylabel('2-theta',fontsize=12)
2577	Plot.contour(tax,tay,taz,V,cmap=acolor)
2578	if Data['showLines']:
2579	IOtth = Data['IOtth']
2580	if Data['fullIntegrate']:
2581	LRAzim = [-180,180]
2582	else:
2583	LRAzim = Data['LRazimuth'] #NB: integers
2584	Plot.plot([LRAzim[0],LRAzim[1]],[IOtth[0],IOtth[0]],picker=True)
2585	Plot.plot([LRAzim[0],LRAzim[1]],[IOtth[1],IOtth[1]],picker=True)
2586	if not Data['fullIntegrate']:
2587	Plot.plot([LRAzim[0],LRAzim[0]],[IOtth[0],IOtth[1]],picker=True)
2588	Plot.plot([LRAzim[1],LRAzim[1]],[IOtth[0],IOtth[1]],picker=True)
2589	if Data['setRings']:
2590	rings = np.concatenate((Data['rings']),axis=0)
2591	for xring,yring,dsp in rings:
2592	x,y = G2img.GetTthAzm(xring,yring,Data)
2593	Plot.plot(y,x,'r+')
2594	if Data['ellipses']:
2595	for ellipse in Data['ellipses']:
2596	ring = np.array(G2img.makeIdealRing(ellipse[:3])) #skip color
2597	x,y = np.hsplit(ring,2)
2598	tth,azm = G2img.GetTthAzm(x,y,Data)
2599	Plot.plot(azm,tth,'b,')
2600	if not newPlot:
2601	Page.toolbar.push_current()
2602	Plot.set_xlim(xylim[0])
2603	Plot.set_ylim(xylim[1])
2604	xylim = []
2605	Page.toolbar.push_current()
2606	Page.toolbar.draw()
2607	else:
2608	Page.canvas.draw()
2609
2610	################################################################################
2611	##### PlotStructure
2612	################################################################################
2613
2614	def PlotStructure(G2frame,data):
2615	'''Crystal structure plotting package. Can show structures as balls, sticks, lines,
2616	thermal motion ellipsoids and polyhedra
2617	'''
2618
2619	def FindPeaksBonds(XYZ):
2620	rFact = data['Drawing']['radiusFactor']
2621	Bonds = [[] for x in XYZ]
2622	for i,xyz in enumerate(XYZ):
2623	Dx = XYZ-xyz
2624	dist = np.sqrt(np.sum(np.inner(Dx,Amat)**2,axis=1))
2625	IndB = ma.nonzero(ma.masked_greater(dist,rFact*2.2))
2626	for j in IndB[0]:
2627	Bonds[i].append(Dx[j]/2.)
2628	Bonds[j].append(-Dx[j]/2.)
2629	return Bonds
2630
2631	ForthirdPI = 4.0*math.pi/3.0
2632	generalData = data['General']
2633	cell = generalData['Cell'][1:7]
2634	Vol = generalData['Cell'][7:8][0]
2635	Amat,Bmat = G2lat.cell2AB(cell) #Amat - crystal to cartesian, Bmat - inverse
2636	Gmat,gmat = G2lat.cell2Gmat(cell)
2637	A4mat = np.concatenate((np.concatenate((Amat,[[0],[0],[0]]),axis=1),[[0,0,0,1],]),axis=0)
2638	B4mat = np.concatenate((np.concatenate((Bmat,[[0],[0],[0]]),axis=1),[[0,0,0,1],]),axis=0)
2639	SGData = generalData['SGData']
2640	Mydir = generalData['Mydir']
2641	atomData = data['Atoms']
2642	mapPeaks = []
2643	drawingData = data['Drawing']
2644	if 'Map Peaks' in data:
2645	mapPeaks = np.array(data['Map Peaks'])
2646	peakMax = 100.
2647	if len(mapPeaks):
2648	peakMax = np.max(mapPeaks.T[0])
2649	resRBData = data['RBModels'].get('Residue',[])
2650	vecRBData = data['RBModels'].get('Vector',[])
2651	rbAtmDict = {}
2652	for rbObj in resRBData+vecRBData:
2653	exclList = ['X' for i in range(len(rbObj['Ids']))]
2654	rbAtmDict.update(dict(zip(rbObj['Ids'],exclList)))
2655	testRBObj = data.get('testRBObj',{})
2656	rbObj = testRBObj.get('rbObj',{})
2657	MCSA = data.get('MCSA',{})
2658	mcsaModels = MCSA.get('Models',[])
2659	if mcsaModels:
2660	XYZs,Types = G2mth.UpdateMCSAxyz(Bmat,MCSA)
2661	mcsaXYZ = []
2662	mcsaTypes = []
2663	for xyz,atyp in zip(XYZs,Types):
2664	for item in G2spc.GenAtom(xyz,SGData):
2665	mcsaXYZ.append(item[0])
2666	mcsaTypes.append(atyp)
2667	mcsaBonds = FindPeaksBonds(mcsaXYZ)
2668	# mcsaXYZ,atTypes = G2mth.UpdateMCSAxyz(Bmat,MCSA)
2669	# XYZeq = []
2670	# for xyz in mcsaXYZ:
2671	# XYZeq += G2spc.GenAtom(xyz,SGData)[0][1:] #skip self xyz
2672	#
2673	# mcsaBonds = FindPeaksBonds(mcsaXYZ)
2674	drawAtoms = drawingData.get('Atoms',[])
2675	mapData = {}
2676	flipData = {}
2677	rhoXYZ = []
2678	showBonds = False
2679	if 'Map' in generalData:
2680	mapData = generalData['Map']
2681	showBonds = mapData.get('Show bonds',False)
2682	if 'Flip' in generalData:
2683	flipData = generalData['Flip']
2684	flipData['mapRoll'] = [0,0,0]
2685	Wt = np.array([255,255,255])
2686	Rd = np.array([255,0,0])
2687	Gr = np.array([0,255,0])
2688	Bl = np.array([0,0,255])
2689	Or = np.array([255,128,0])
2690	uBox = np.array([[0,0,0],[1,0,0],[1,1,0],[0,1,0],[0,0,1],[1,0,1],[1,1,1],[0,1,1]])
2691	uEdges = np.array([
2692	[uBox[0],uBox[1]],[uBox[0],uBox[3]],[uBox[0],uBox[4]],[uBox[1],uBox[2]],
2693	[uBox[2],uBox[3]],[uBox[1],uBox[5]],[uBox[2],uBox[6]],[uBox[3],uBox[7]],
2694	[uBox[4],uBox[5]],[uBox[5],uBox[6]],[uBox[6],uBox[7]],[uBox[7],uBox[4]]])
2695	mD = 0.1
2696	mV = np.array([[[-mD,0,0],[mD,0,0]],[[0,-mD,0],[0,mD,0]],[[0,0,-mD],[0,0,mD]]])
2697	mapPeakVecs = np.inner(mV,Bmat)
2698
2699	uColors = [Rd,Gr,Bl,Wt, Wt,Wt,Wt,Wt, Wt,Wt,Wt,Wt]
2700	altDown = False
2701	shiftDown = False
2702	ctrlDown = False
2703
2704	def OnKeyBox(event):
2705	import Image
2706	# Draw() #make sure plot is fresh!!
2707	mode = cb.GetValue()
2708	if mode in ['jpeg','bmp','tiff',]:
2709	Fname = os.path.join(Mydir,generalData['Name']+'.'+mode)
2710	size = Page.canvas.GetSize()
2711	glPixelStorei(GL_UNPACK_ALIGNMENT, 1)
2712	if mode in ['jpeg',]:
2713	Pix = glReadPixels(0,0,size[0],size[1],GL_RGBA, GL_UNSIGNED_BYTE)
2714	im = Image.new("RGBA", (size[0],size[1]))
2715	else:
2716	Pix = glReadPixels(0,0,size[0],size[1],GL_RGB, GL_UNSIGNED_BYTE)
2717	im = Image.new("RGB", (size[0],size[1]))
2718	im.fromstring(Pix)
2719	im.save(Fname,mode)
2720	cb.SetValue(' save as/key:')
2721	G2frame.G2plotNB.status.SetStatusText('Drawing saved to: '+Fname,1)
2722	else:
2723	event.key = cb.GetValue()[0]
2724	cb.SetValue(' save as/key:')
2725	wx.CallAfter(OnKey,event)
2726
2727	def OnKey(event): #on key UP!!
2728	# Draw() #make sure plot is fresh!!
2729	try:
2730	keyCode = event.GetKeyCode()
2731	if keyCode > 255:
2732	keyCode = 0
2733	key = chr(keyCode)
2734	except AttributeError: #if from OnKeyBox above
2735	key = str(event.key).upper()
2736	indx = drawingData['selectedAtoms']
2737	cx,ct = drawingData['atomPtrs'][:2]
2738	if key in ['C']:
2739	drawingData['viewPoint'] = [[.5,.5,.5],[0,0]]
2740	drawingData['viewDir'] = [0,0,1]
2741	drawingData['oldxy'] = []
2742	V0 = np.array([0,0,1])
2743	V = np.inner(Amat,V0)
2744	V /= np.sqrt(np.sum(V**2))
2745	A = np.arccos(np.sum(V*V0))
2746	Q = G2mth.AV2Q(A,[0,1,0])
2747	drawingData['Quaternion'] = Q
2748	SetViewPointText(drawingData['viewPoint'][0])
2749	SetViewDirText(drawingData['viewDir'])
2750	Q = drawingData['Quaternion']
2751	G2frame.G2plotNB.status.SetStatusText('New quaternion: %.2f+, %.2fi+ ,%.2fj+, %.2fk'%(Q[0],Q[1],Q[2],Q[3]),1)
2752	elif key in ['N']:
2753	drawAtoms = drawingData['Atoms']
2754	pI = drawingData['viewPoint'][1]
2755	if indx:
2756	pI[0] = indx[pI[1]]
2757	Tx,Ty,Tz = drawAtoms[pI[0]][cx:cx+3]
2758	pI[1] += 1
2759	if pI[1] >= len(indx):
2760	pI[1] = 0
2761	else:
2762	Tx,Ty,Tz = drawAtoms[pI[0]][cx:cx+3]
2763	pI[0] += 1
2764	if pI[0] >= len(drawAtoms):
2765	pI[0] = 0
2766	drawingData['viewPoint'] = [[Tx,Ty,Tz],pI]
2767	SetViewPointText(drawingData['viewPoint'][0])
2768	G2frame.G2plotNB.status.SetStatusText('View point at atom '+drawAtoms[pI[0]][ct-1]+str(pI),1)
2769
2770	elif key in ['P']:
2771	drawAtoms = drawingData['Atoms']
2772	pI = drawingData['viewPoint'][1]
2773	if indx:
2774	pI[0] = indx[pI[1]]
2775	Tx,Ty,Tz = drawAtoms[pI[0]][cx:cx+3]
2776	pI[1] -= 1
2777	if pI[1] < 0:
2778	pI[1] = len(indx)-1
2779	else:
2780	Tx,Ty,Tz = drawAtoms[pI[0]][cx:cx+3]
2781	pI[0] -= 1
2782	if pI[0] < 0:
2783	pI[0] = len(drawAtoms)-1
2784	drawingData['viewPoint'] = [[Tx,Ty,Tz],pI]
2785	SetViewPointText(drawingData['viewPoint'][0])
2786	G2frame.G2plotNB.status.SetStatusText('View point at atom '+drawAtoms[pI[0]][ct-1]+str(pI),1)
2787	elif key in ['U','D','L','R'] and mapData['Flip'] == True:
2788	dirDict = {'U':[0,1],'D':[0,-1],'L':[-1,0],'R':[1,0]}
2789	SetMapRoll(dirDict[key])
2790	SetPeakRoll(dirDict[key])
2791	SetMapPeaksText(mapPeaks)
2792	Draw('key')
2793
2794	def GetTruePosition(xy,Add=False):
2795	View = glGetIntegerv(GL_VIEWPORT)
2796	Proj = glGetDoublev(GL_PROJECTION_MATRIX)
2797	Model = glGetDoublev(GL_MODELVIEW_MATRIX)
2798	Zmax = 1.
2799	if Add:
2800	Indx = GetSelectedAtoms()
2801	if G2frame.dataDisplay.GetPageText(getSelection()) == 'Map peaks':
2802	for i,peak in enumerate(mapPeaks):
2803	x,y,z = peak[1:4]
2804	X,Y,Z = gluProject(x,y,z,Model,Proj,View)
2805	XY = [int(X),int(View[3]-Y)]
2806	if np.allclose(xy,XY,atol=10) and Z < Zmax:
2807	Zmax = Z
2808	try:
2809	Indx.remove(i)
2810	ClearSelectedAtoms()
2811	for id in Indx:
2812	SetSelectedAtoms(id,Add)
2813	except:
2814	SetSelectedAtoms(i,Add)
2815	else:
2816	cx = drawingData['atomPtrs'][0]
2817	for i,atom in enumerate(drawAtoms):
2818	x,y,z = atom[cx:cx+3]
2819	X,Y,Z = gluProject(x,y,z,Model,Proj,View)
2820	XY = [int(X),int(View[3]-Y)]
2821	if np.allclose(xy,XY,atol=10) and Z < Zmax:
2822	Zmax = Z
2823	try:
2824	Indx.remove(i)
2825	ClearSelectedAtoms()
2826	for id in Indx:
2827	SetSelectedAtoms(id,Add)
2828	except:
2829	SetSelectedAtoms(i,Add)
2830
2831	def OnMouseDown(event):
2832	xy = event.GetPosition()
2833	if event.ShiftDown():
2834	if event.LeftIsDown():
2835	GetTruePosition(xy)
2836	elif event.RightIsDown():
2837	GetTruePosition(xy,True)
2838	else:
2839	drawingData['oldxy'] = list(xy)
2840
2841	def OnMouseMove(event):
2842	if event.ShiftDown(): #don't want any inadvertant moves when picking
2843	return
2844	newxy = event.GetPosition()
2845
2846	if event.Dragging():
2847	if event.AltDown() and rbObj:
2848	if event.LeftIsDown():
2849	SetRBRotation(newxy)
2850	Q = rbObj['Orient'][0]
2851	G2frame.G2plotNB.status.SetStatusText('New quaternion: %.2f+, %.2fi+ ,%.2fj+, %.2fk'%(Q[0],Q[1],Q[2],Q[3]),1)
2852	elif event.RightIsDown():
2853	SetRBTranslation(newxy)
2854	Tx,Ty,Tz = rbObj['Orig'][0]
2855	G2frame.G2plotNB.status.SetStatusText('New view point: %.4f, %.4f, %.4f'%(Tx,Ty,Tz),1)
2856	elif event.MiddleIsDown():
2857	SetRBRotationZ(newxy)
2858	Q = rbObj['Orient'][0]
2859	G2frame.G2plotNB.status.SetStatusText('New quaternion: %.2f+, %.2fi+ ,%.2fj+, %.2fk'%(Q[0],Q[1],Q[2],Q[3]),1)
2860	Draw('move')
2861	elif not event.ControlDown():
2862	if event.LeftIsDown():
2863	SetRotation(newxy)
2864	Q = drawingData['Quaternion']
2865	G2frame.G2plotNB.status.SetStatusText('New quaternion: %.2f+, %.2fi+ ,%.2fj+, %.2fk'%(Q[0],Q[1],Q[2],Q[3]),1)
2866	elif event.RightIsDown():
2867	SetTranslation(newxy)
2868	Tx,Ty,Tz = drawingData['viewPoint'][0]
2869	G2frame.G2plotNB.status.SetStatusText('New view point: %.4f, %.4f, %.4f'%(Tx,Ty,Tz),1)
2870	elif event.MiddleIsDown():
2871	SetRotationZ(newxy)
2872	Q = drawingData['Quaternion']
2873	G2frame.G2plotNB.status.SetStatusText('New quaternion: %.2f+, %.2fi+ ,%.2fj+, %.2fk'%(Q[0],Q[1],Q[2],Q[3]),1)
2874	Draw('move')
2875
2876
2877	def OnMouseWheel(event):
2878	if event.ShiftDown():
2879	return
2880	drawingData['cameraPos'] += event.GetWheelRotation()/24.
2881	drawingData['cameraPos'] = max(10,min(500,drawingData['cameraPos']))
2882	G2frame.G2plotNB.status.SetStatusText('New camera distance: %.2f'%(drawingData['cameraPos']),1)
2883	page = getSelection()
2884	if page:
2885	if G2frame.dataDisplay.GetPageText(page) == 'Draw Options':
2886	G2frame.dataDisplay.cameraPosTxt.SetLabel('Camera Position: '+'%.2f'%(drawingData['cameraPos']))
2887	G2frame.dataDisplay.cameraSlider.SetValue(drawingData['cameraPos'])
2888	Draw('wheel')
2889
2890	def getSelection():
2891	try:
2892	return G2frame.dataDisplay.GetSelection()
2893	except AttributeError:
2894	G2frame.G2plotNB.status.SetStatusText('Select this from Phase data window!')
2895	return 0
2896
2897	def SetViewPointText(VP):
2898	page = getSelection()
2899	if page:
2900	if G2frame.dataDisplay.GetPageText(page) == 'Draw Options':
2901	G2frame.dataDisplay.viewPoint.SetValue('%.3f %.3f %.3f'%(VP[0],VP[1],VP[2]))
2902
2903	def SetRBOrigText():
2904	page = getSelection()
2905	if page:
2906	if G2frame.dataDisplay.GetPageText(page) == 'RB Models':
2907	for i,sizer in enumerate(testRBObj['Sizers']['Xsizers']):
2908	sizer.SetValue('%8.5f'%(testRBObj['rbObj']['Orig'][0][i]))
2909
2910	def SetRBOrienText():
2911	page = getSelection()
2912	if page:
2913	if G2frame.dataDisplay.GetPageText(page) == 'RB Models':
2914	for i,sizer in enumerate(testRBObj['Sizers']['Osizers']):
2915	sizer.SetValue('%8.5f'%(testRBObj['rbObj']['Orient'][0][i]))
2916
2917	def SetViewDirText(VD):
2918	page = getSelection()
2919	if page:
2920	if G2frame.dataDisplay.GetPageText(page) == 'Draw Options':
2921	G2frame.dataDisplay.viewDir.SetValue('%.3f %.3f %.3f'%(VD[0],VD[1],VD[2]))
2922
2923	def SetMapPeaksText(mapPeaks):
2924	page = getSelection()
2925	if page:
2926	if G2frame.dataDisplay.GetPageText(page) == 'Map peaks':
2927	G2frame.MapPeaksTable.SetData(mapPeaks)
2928	panel = G2frame.dataDisplay.GetPage(page).GetChildren()
2929	names = [child.GetName() for child in panel]
2930	panel[names.index('grid window')].Refresh()
2931
2932	def ClearSelectedAtoms():
2933	page = getSelection()
2934	if page:
2935	if G2frame.dataDisplay.GetPageText(page) == 'Draw Atoms':
2936	G2frame.dataDisplay.GetPage(page).ClearSelection() #this is the Atoms grid in Draw Atoms
2937	elif G2frame.dataDisplay.GetPageText(page) == 'Map peaks':
2938	G2frame.dataDisplay.GetPage(page).ClearSelection() #this is the Atoms grid in Atoms
2939	elif G2frame.dataDisplay.GetPageText(page) == 'Atoms':
2940	G2frame.dataDisplay.GetPage(page).ClearSelection() #this is the Atoms grid in Atoms
2941
2942
2943	def SetSelectedAtoms(ind,Add=False):
2944	page = getSelection()
2945	if page:
2946	if G2frame.dataDisplay.GetPageText(page) == 'Draw Atoms':
2947	G2frame.dataDisplay.GetPage(page).SelectRow(ind,Add) #this is the Atoms grid in Draw Atoms
2948	elif G2frame.dataDisplay.GetPageText(page) == 'Map peaks':
2949	G2frame.dataDisplay.GetPage(page).SelectRow(ind,Add)
2950	elif G2frame.dataDisplay.GetPageText(page) == 'Atoms':
2951	Id = drawAtoms[ind][-3]
2952	for i,atom in enumerate(atomData):
2953	if atom[-1] == Id:
2954	G2frame.dataDisplay.GetPage(page).SelectRow(i) #this is the Atoms grid in Atoms
2955
2956	def GetSelectedAtoms():
2957	page = getSelection()
2958	Ind = []
2959	if page:
2960	if G2frame.dataDisplay.GetPageText(page) == 'Draw Atoms':
2961	Ind = G2frame.dataDisplay.GetPage(page).GetSelectedRows() #this is the Atoms grid in Draw Atoms
2962	elif G2frame.dataDisplay.GetPageText(page) == 'Map peaks':
2963	Ind = G2frame.dataDisplay.GetPage(page).GetSelectedRows()
2964	elif G2frame.dataDisplay.GetPageText(page) == 'Atoms':
2965	Ind = G2frame.dataDisplay.GetPage(page).GetSelectedRows() #this is the Atoms grid in Atoms
2966	return Ind
2967
2968	def SetBackground():
2969	R,G,B,A = Page.camera['backColor']
2970	glClearColor(R,G,B,A)
2971	glClear(GL_COLOR_BUFFER_BIT \| GL_DEPTH_BUFFER_BIT)
2972
2973	def SetLights():
2974	glEnable(GL_DEPTH_TEST)
2975	glShadeModel(GL_SMOOTH)
2976	glEnable(GL_LIGHTING)
2977	glEnable(GL_LIGHT0)
2978	glLightModeli(GL_LIGHT_MODEL_TWO_SIDE,0)
2979	glLightfv(GL_LIGHT0,GL_AMBIENT,[1,1,1,.8])
2980	glLightfv(GL_LIGHT0,GL_DIFFUSE,[1,1,1,1])
2981
2982	def GetRoll(newxy,rho):
2983	Q = drawingData['Quaternion']
2984	dxy = G2mth.prodQVQ(G2mth.invQ(Q),np.inner(Bmat,newxy+[0,]))
2985	dxy = np.array(dxy*rho.shape)
2986	roll = np.where(dxy>0.5,1,np.where(dxy<-.5,-1,0))
2987	return roll
2988
2989	def SetMapRoll(newxy):
2990	rho = mapData['rho']
2991	roll = GetRoll(newxy,rho)
2992	mapData['rho'] = np.roll(np.roll(np.roll(rho,roll[0],axis=0),roll[1],axis=1),roll[2],axis=2)
2993	drawingData['oldxy'] = list(newxy)
2994
2995	def SetPeakRoll(newxy):
2996	rho = mapData['rho']
2997	roll = GetRoll(newxy,rho)
2998	steps = 1./np.array(rho.shape)
2999	dxy = roll*steps
3000	for peak in mapPeaks:
3001	peak[1:4] += dxy
3002	peak[1:4] %= 1.
3003	peak[4] = np.sqrt(np.sum(np.inner(Amat,peak[1:4])**2))
3004
3005	def SetTranslation(newxy):
3006	#first get translation vector in screen coords.
3007	oldxy = drawingData['oldxy']
3008	if not len(oldxy): oldxy = list(newxy)
3009	dxy = newxy-oldxy
3010	drawingData['oldxy'] = list(newxy)
3011	V = np.array([-dxy[0],dxy[1],0.])
3012	#then transform to rotated crystal coordinates & apply to view point
3013	Q = drawingData['Quaternion']
3014	V = np.inner(Bmat,G2mth.prodQVQ(G2mth.invQ(Q),V))
3015	Tx,Ty,Tz = drawingData['viewPoint'][0]
3016	Tx += V[0]*0.01
3017	Ty += V[1]*0.01
3018	Tz += V[2]*0.01
3019	drawingData['viewPoint'][0] = Tx,Ty,Tz
3020	SetViewPointText([Tx,Ty,Tz])
3021
3022	def SetRBTranslation(newxy):
3023	#first get translation vector in screen coords.
3024	oldxy = drawingData['oldxy']
3025	if not len(oldxy): oldxy = list(newxy)
3026	dxy = newxy-oldxy
3027	drawingData['oldxy'] = list(newxy)
3028	V = np.array([-dxy[0],dxy[1],0.])
3029	#then transform to rotated crystal coordinates & apply to RB origin
3030	Q = drawingData['Quaternion']
3031	V = np.inner(Bmat,G2mth.prodQVQ(G2mth.invQ(Q),V))
3032	Tx,Ty,Tz = rbObj['Orig'][0]
3033	Tx -= V[0]*0.01
3034	Ty -= V[1]*0.01
3035	Tz -= V[2]*0.01
3036	rbObj['Orig'][0] = Tx,Ty,Tz
3037	SetRBOrigText()
3038
3039	def SetRotation(newxy):
3040	#first get rotation vector in screen coords. & angle increment
3041	oldxy = drawingData['oldxy']
3042	if not len(oldxy): oldxy = list(newxy)
3043	dxy = newxy-oldxy
3044	drawingData['oldxy'] = list(newxy)
3045	V = np.array([dxy[1],dxy[0],0.])
3046	A = 0.25np.sqrt(dxy[0]2+dxy[1]*2)
3047	# next transform vector back to xtal coordinates via inverse quaternion
3048	# & make new quaternion
3049	Q = drawingData['Quaternion']
3050	V = G2mth.prodQVQ(G2mth.invQ(Q),np.inner(Bmat,V))
3051	DQ = G2mth.AVdeg2Q(A,V)
3052	Q = G2mth.prodQQ(Q,DQ)
3053	drawingData['Quaternion'] = Q
3054	# finally get new view vector - last row of rotation matrix
3055	VD = np.inner(Bmat,G2mth.Q2Mat(Q)[2])
3056	VD /= np.sqrt(np.sum(VD**2))
3057	drawingData['viewDir'] = VD
3058	SetViewDirText(VD)
3059
3060	def SetRotationZ(newxy):
3061	#first get rotation vector (= view vector) in screen coords. & angle increment
3062	View = glGetIntegerv(GL_VIEWPORT)
3063	cent = [View[2]/2,View[3]/2]
3064	oldxy = drawingData['oldxy']
3065	if not len(oldxy): oldxy = list(newxy)
3066	dxy = newxy-oldxy
3067	drawingData['oldxy'] = list(newxy)
3068	V = drawingData['viewDir']
3069	A = [0,0]
3070	A[0] = dxy[1]*.25
3071	A[1] = dxy[0]*.25
3072	if newxy[0] > cent[0]:
3073	A[0] *= -1
3074	if newxy[1] < cent[1]:
3075	A[1] *= -1
3076	# next transform vector back to xtal coordinates & make new quaternion
3077	Q = drawingData['Quaternion']
3078	V = np.inner(Amat,V)
3079	Qx = G2mth.AVdeg2Q(A[0],V)
3080	Qy = G2mth.AVdeg2Q(A[1],V)
3081	Q = G2mth.prodQQ(Q,Qx)
3082	Q = G2mth.prodQQ(Q,Qy)
3083	drawingData['Quaternion'] = Q
3084
3085	def SetRBRotation(newxy):
3086	#first get rotation vector in screen coords. & angle increment
3087	oldxy = drawingData['oldxy']
3088	if not len(oldxy): oldxy = list(newxy)
3089	dxy = newxy-oldxy
3090	drawingData['oldxy'] = list(newxy)
3091	V = np.array([dxy[1],dxy[0],0.])
3092	A = 0.25np.sqrt(dxy[0]2+dxy[1]*2)
3093	# next transform vector back to xtal coordinates via inverse quaternion
3094	# & make new quaternion
3095	Q = rbObj['Orient'][0] #rotate RB to Cart
3096	QC = drawingData['Quaternion'] #rotate Cart to drawing
3097	V = G2mth.prodQVQ(G2mth.invQ(QC),V)
3098	V = G2mth.prodQVQ(G2mth.invQ(Q),V)
3099	DQ = G2mth.AVdeg2Q(A,V)
3100	Q = G2mth.prodQQ(Q,DQ)
3101	rbObj['Orient'][0] = Q
3102	SetRBOrienText()
3103
3104	def SetRBRotationZ(newxy):
3105	#first get rotation vector (= view vector) in screen coords. & angle increment
3106	View = glGetIntegerv(GL_VIEWPORT)
3107	cent = [View[2]/2,View[3]/2]
3108	oldxy = drawingData['oldxy']
3109	if not len(oldxy): oldxy = list(newxy)
3110	dxy = newxy-oldxy
3111	drawingData['oldxy'] = list(newxy)
3112	V = drawingData['viewDir']
3113	A = [0,0]
3114	A[0] = dxy[1]*.25
3115	A[1] = dxy[0]*.25
3116	if newxy[0] < cent[0]:
3117	A[0] *= -1
3118	if newxy[1] > cent[1]:
3119	A[1] *= -1
3120	# next transform vector back to RB coordinates & make new quaternion
3121	Q = rbObj['Orient'][0] #rotate RB to cart
3122	V = np.inner(Amat,V)
3123	V = -G2mth.prodQVQ(G2mth.invQ(Q),V)
3124	Qx = G2mth.AVdeg2Q(A[0],V)
3125	Qy = G2mth.AVdeg2Q(A[1],V)
3126	Q = G2mth.prodQQ(Q,Qx)
3127	Q = G2mth.prodQQ(Q,Qy)
3128	rbObj['Orient'][0] = Q
3129	SetRBOrienText()
3130
3131	def RenderBox():
3132	glEnable(GL_COLOR_MATERIAL)
3133	glLineWidth(2)
3134	glEnable(GL_BLEND)
3135	glBlendFunc(GL_SRC_ALPHA,GL_ONE_MINUS_SRC_ALPHA)
3136	glEnable(GL_LINE_SMOOTH)
3137	glBegin(GL_LINES)
3138	for line,color in zip(uEdges,uColors):
3139	glColor3ubv(color)
3140	glVertex3fv(line[0])
3141	glVertex3fv(line[1])
3142	glEnd()
3143	glColor4ubv([0,0,0,0])
3144	glDisable(GL_LINE_SMOOTH)
3145	glDisable(GL_BLEND)
3146	glDisable(GL_COLOR_MATERIAL)
3147
3148	def RenderUnitVectors(x,y,z):
3149	xyz = np.array([x,y,z])
3150	glEnable(GL_COLOR_MATERIAL)
3151	glLineWidth(1)
3152	glPushMatrix()
3153	glTranslate(x,y,z)
3154	glScalef(1/cell[0],1/cell[1],1/cell[2])
3155	glBegin(GL_LINES)
3156	for line,color in zip(uEdges,uColors)[:3]:
3157	glColor3ubv(color)
3158	glVertex3fv(-line[1]/2.)
3159	glVertex3fv(line[1]/2.)
3160	glEnd()
3161	glPopMatrix()
3162	glColor4ubv([0,0,0,0])
3163	glDisable(GL_COLOR_MATERIAL)
3164
3165	def RenderSphere(x,y,z,radius,color):
3166	glMaterialfv(GL_FRONT_AND_BACK,GL_DIFFUSE,color)
3167	glPushMatrix()
3168	glTranslate(x,y,z)
3169	glMultMatrixf(B4mat.T)
3170	q = gluNewQuadric()
3171	gluSphere(q,radius,20,10)
3172	glPopMatrix()
3173
3174	def RenderDots(XYZ,RC):
3175	glEnable(GL_COLOR_MATERIAL)
3176	XYZ = np.array(XYZ)
3177	glPushMatrix()
3178	for xyz,rc in zip(XYZ,RC):
3179	x,y,z = xyz
3180	r,c = rc
3181	glColor3ubv(c)
3182	glPointSize(r*50)
3183	glBegin(GL_POINTS)
3184	glVertex3fv(xyz)
3185	glEnd()
3186	glPopMatrix()
3187	glColor4ubv([0,0,0,0])
3188	glDisable(GL_COLOR_MATERIAL)
3189
3190	def RenderSmallSphere(x,y,z,radius,color):
3191	glMaterialfv(GL_FRONT_AND_BACK,GL_DIFFUSE,color)
3192	glPushMatrix()
3193	glTranslate(x,y,z)
3194	glMultMatrixf(B4mat.T)
3195	q = gluNewQuadric()
3196	gluSphere(q,radius,4,2)
3197	glPopMatrix()
3198
3199	def RenderEllipsoid(x,y,z,ellipseProb,E,R4,color):
3200	s1,s2,s3 = E
3201	glMaterialfv(GL_FRONT_AND_BACK,GL_DIFFUSE,color)
3202	glPushMatrix()
3203	glTranslate(x,y,z)
3204	glMultMatrixf(B4mat.T)
3205	glMultMatrixf(R4.T)
3206	glEnable(GL_NORMALIZE)
3207	glScale(s1,s2,s3)
3208	q = gluNewQuadric()
3209	gluSphere(q,ellipseProb,20,10)
3210	glDisable(GL_NORMALIZE)
3211	glPopMatrix()
3212
3213	def RenderBonds(x,y,z,Bonds,radius,color,slice=20):
3214	glMaterialfv(GL_FRONT_AND_BACK,GL_DIFFUSE,color)
3215	glPushMatrix()
3216	glTranslate(x,y,z)
3217	glMultMatrixf(B4mat.T)
3218	for bond in Bonds:
3219	glPushMatrix()
3220	Dx = np.inner(Amat,bond)
3221	Z = np.sqrt(np.sum(Dx**2))
3222	if Z:
3223	azm = atan2d(-Dx[1],-Dx[0])
3224	phi = acosd(Dx[2]/Z)
3225	glRotate(-azm,0,0,1)
3226	glRotate(phi,1,0,0)
3227	q = gluNewQuadric()
3228	gluCylinder(q,radius,radius,Z,slice,2)
3229	glPopMatrix()
3230	glPopMatrix()
3231
3232	def RenderLines(x,y,z,Bonds,color):
3233	glShadeModel(GL_FLAT)
3234	xyz = np.array([x,y,z])
3235	glEnable(GL_COLOR_MATERIAL)
3236	glLineWidth(1)
3237	glColor3fv(color)
3238	glPushMatrix()
3239	glBegin(GL_LINES)
3240	for bond in Bonds:
3241	glVertex3fv(xyz)
3242	glVertex3fv(xyz+bond)
3243	glEnd()
3244	glColor4ubv([0,0,0,0])
3245	glPopMatrix()
3246	glDisable(GL_COLOR_MATERIAL)
3247	glShadeModel(GL_SMOOTH)
3248
3249	def RenderPolyhedra(x,y,z,Faces,color):
3250	glShadeModel(GL_FLAT)
3251	glPushMatrix()
3252	glTranslate(x,y,z)
3253	glMaterialfv(GL_FRONT_AND_BACK,GL_DIFFUSE,color)
3254	glShadeModel(GL_SMOOTH)
3255	glMultMatrixf(B4mat.T)
3256	for face,norm in Faces:
3257	glPolygonMode(GL_FRONT_AND_BACK,GL_FILL)
3258	glFrontFace(GL_CW)
3259	glNormal3fv(norm)
3260	glBegin(GL_TRIANGLES)
3261	for vert in face:
3262	glVertex3fv(vert)
3263	glEnd()
3264	glPopMatrix()
3265	glShadeModel(GL_SMOOTH)
3266
3267	def RenderMapPeak(x,y,z,color,den):
3268	glShadeModel(GL_FLAT)
3269	xyz = np.array([x,y,z])
3270	glEnable(GL_COLOR_MATERIAL)
3271	glLineWidth(3)
3272	glColor3fv(color*den/255)
3273	glPushMatrix()
3274	glBegin(GL_LINES)
3275	for vec in mapPeakVecs:
3276	glVertex3fv(vec[0]+xyz)
3277	glVertex3fv(vec[1]+xyz)
3278	glEnd()
3279	glColor4ubv([0,0,0,0])
3280	glPopMatrix()
3281	glDisable(GL_COLOR_MATERIAL)
3282	glShadeModel(GL_SMOOTH)
3283
3284	def RenderBackbone(Backbone,BackboneColor,radius):
3285	glPushMatrix()
3286	glMultMatrixf(B4mat.T)
3287	glEnable(GL_COLOR_MATERIAL)
3288	glShadeModel(GL_SMOOTH)
3289	gleSetJoinStyle(TUBE_NORM_EDGE \| TUBE_JN_ANGLE \| TUBE_JN_CAP)
3290	glePolyCylinder(Backbone,BackboneColor,radius)
3291	glPopMatrix()
3292	glDisable(GL_COLOR_MATERIAL)
3293
3294	def RenderLabel(x,y,z,label,r,color):
3295	glPushMatrix()
3296	glTranslate(x,y,z)
3297	glMultMatrixf(B4mat.T)
3298	glDisable(GL_LIGHTING)
3299	glColor3fv(color)
3300	glRasterPos3f(0,0,0)
3301	if bool(glutBitmapCharacter): #seems to only exist in 32 bit Windows
3302	for c in list(label):
3303	glutBitmapCharacter(GLUT_BITMAP_8_BY_13,ord(c))
3304	else:
3305	text = gltext.TextElement(text=label,font=Font)
3306	text.draw_text(scale=0.025)
3307	# our_font.glPrint(0,0,label)
3308	glEnable(GL_LIGHTING)
3309	glPopMatrix()
3310
3311	def RenderMap(rho,rhoXYZ,indx,Rok):
3312	glShadeModel(GL_FLAT)
3313	cLevel = drawingData['contourLevel']
3314	XYZ = []
3315	RC = []
3316	for i,xyz in enumerate(rhoXYZ):
3317	if not Rok[i]:
3318	x,y,z = xyz
3319	I,J,K = indx[i]
3320	alpha = 1.0
3321	if cLevel < 1.:
3322	alpha = (abs(rho[I,J,K])/mapData['rhoMax']-cLevel)/(1.-cLevel)
3323	if rho[I,J,K] < 0.:
3324	XYZ.append(xyz)
3325	RC.append([0.1*alpha,Rd])
3326	else:
3327	XYZ.append(xyz)
3328	RC.append([0.1*alpha,Gr])
3329	RenderDots(XYZ,RC)
3330	glShadeModel(GL_SMOOTH)
3331
3332	def Draw(caller=''):
3333	#useful debug?
3334	# if caller:
3335	# print caller
3336	# end of useful debug
3337	mapData = generalData['Map']
3338	pageName = ''
3339	page = getSelection()
3340	if page:
3341	pageName = G2frame.dataDisplay.GetPageText(page)
3342	rhoXYZ = []
3343	if len(mapData['rho']):
3344	VP = np.array(drawingData['viewPoint'][0])-np.array([.5,.5,.5])
3345	contLevel = drawingData['contourLevel']*mapData['rhoMax']
3346	if 'delt-F' in mapData['MapType']:
3347	rho = ma.array(mapData['rho'],mask=(np.abs(mapData['rho'])<contLevel))
3348	else:
3349	rho = ma.array(mapData['rho'],mask=(mapData['rho']<contLevel))
3350	steps = 1./np.array(rho.shape)
3351	incre = np.where(VP>=0,VP%steps,VP%steps-steps)
3352	Vsteps = -np.array(VP/steps,dtype='i')
3353	rho = np.roll(np.roll(np.roll(rho,Vsteps[0],axis=0),Vsteps[1],axis=1),Vsteps[2],axis=2)
3354	indx = np.array(ma.nonzero(rho)).T
3355	rhoXYZ = indx*steps+VP-incre
3356	Nc = len(rhoXYZ)
3357	rcube = 2000.Vol/(ForthirdPINc)
3358	rmax = math.exp(math.log(rcube)/3.)**2
3359	radius = min(drawingData['mapSize']**2,rmax)
3360	view = np.array(drawingData['viewPoint'][0])
3361	Rok = np.sum(np.inner(Amat,rhoXYZ-view).T**2,axis=1)>radius
3362	Ind = GetSelectedAtoms()
3363	VS = np.array(Page.canvas.GetSize())
3364	aspect = float(VS[0])/float(VS[1])
3365	cPos = drawingData['cameraPos']
3366	Zclip = drawingData['Zclip']*cPos/200.
3367	Q = drawingData['Quaternion']
3368	Tx,Ty,Tz = drawingData['viewPoint'][0]
3369	cx,ct,cs,ci = drawingData['atomPtrs']
3370	bondR = drawingData['bondRadius']
3371	G,g = G2lat.cell2Gmat(cell)
3372	GS = G
3373	GS[0][1] = GS[1][0] = math.sqrt(GS[0][0]*GS[1][1])
3374	GS[0][2] = GS[2][0] = math.sqrt(GS[0][0]*GS[2][2])
3375	GS[1][2] = GS[2][1] = math.sqrt(GS[1][1]*GS[2][2])
3376	ellipseProb = G2lat.criticalEllipse(drawingData['ellipseProb']/100.)
3377
3378	SetBackground()
3379	glInitNames()
3380	glPushName(0)
3381
3382	glMatrixMode(GL_PROJECTION)
3383	glLoadIdentity()
3384	glViewport(0,0,VS[0],VS[1])
3385	gluPerspective(20.,aspect,cPos-Zclip,cPos+Zclip)
3386	gluLookAt(0,0,cPos,0,0,0,0,1,0)
3387	SetLights()
3388
3389	glMatrixMode(GL_MODELVIEW)
3390	glLoadIdentity()
3391	matRot = G2mth.Q2Mat(Q)
3392	matRot = np.concatenate((np.concatenate((matRot,[[0],[0],[0]]),axis=1),[[0,0,0,1],]),axis=0)
3393	glMultMatrixf(matRot.T)
3394	glMultMatrixf(A4mat.T)
3395	glTranslate(-Tx,-Ty,-Tz)
3396	if drawingData['unitCellBox']:
3397	RenderBox()
3398	if drawingData['showABC']:
3399	x,y,z = drawingData['viewPoint'][0]
3400	RenderUnitVectors(x,y,z)
3401	Backbones = {}
3402	BackboneColor = []
3403	time0 = time.time()
3404	# glEnable(GL_BLEND)
3405	# glBlendFunc(GL_SRC_ALPHA,GL_ONE_MINUS_SRC_ALPHA)
3406	for iat,atom in enumerate(drawingData['Atoms']):
3407	x,y,z = atom[cx:cx+3]
3408	Bonds = atom[-2]
3409	Faces = atom[-1]
3410	try:
3411	atNum = generalData['AtomTypes'].index(atom[ct])
3412	except ValueError:
3413	atNum = -1
3414	CL = atom[cs+2]
3415	atColor = np.array(CL)/255.
3416	if drawingData['showRigidBodies'] and atom[ci] in rbAtmDict:
3417	bndColor = Or
3418	else:
3419	bndColor = atColor
3420	if iat in Ind and G2frame.dataDisplay.GetPageText(getSelection()) != 'Map peaks':
3421	atColor = np.array(Gr)/255.
3422	# color += [.25,]
3423	radius = 0.5
3424	if atom[cs] != '':
3425	try:
3426	glLoadName(atom[-3])
3427	except: #problem with old files - missing code
3428	pass
3429	if 'balls' in atom[cs]:
3430	vdwScale = drawingData['vdwScale']
3431	ballScale = drawingData['ballScale']
3432	if atNum < 0:
3433	radius = 0.2
3434	elif 'H' == atom[ct]:
3435	if drawingData['showHydrogen']:
3436	if 'vdW' in atom[cs] and atNum >= 0:
3437	radius = vdwScale*generalData['vdWRadii'][atNum]
3438	else:
3439	radius = ballScale*drawingData['sizeH']
3440	else:
3441	radius = 0.0
3442	else:
3443	if 'vdW' in atom[cs]:
3444	radius = vdwScale*generalData['vdWRadii'][atNum]
3445	else:
3446	radius = ballScale*generalData['BondRadii'][atNum]
3447	RenderSphere(x,y,z,radius,atColor)
3448	if 'sticks' in atom[cs]:
3449	RenderBonds(x,y,z,Bonds,bondR,bndColor)
3450	elif 'ellipsoids' in atom[cs]:
3451	RenderBonds(x,y,z,Bonds,bondR,bndColor)
3452	if atom[cs+3] == 'A':
3453	Uij = atom[cs+5:cs+11]
3454	U = np.multiply(G2spc.Uij2U(Uij),GS)
3455	U = np.inner(Amat,np.inner(U,Amat).T)
3456	E,R = nl.eigh(U)
3457	R4 = np.concatenate((np.concatenate((R,[[0],[0],[0]]),axis=1),[[0,0,0,1],]),axis=0)
3458	E = np.sqrt(E)
3459	if atom[ct] == 'H' and not drawingData['showHydrogen']:
3460	pass
3461	else:
3462	RenderEllipsoid(x,y,z,ellipseProb,E,R4,atColor)
3463	else:
3464	if atom[ct] == 'H' and not drawingData['showHydrogen']:
3465	pass
3466	else:
3467	radius = ellipseProb*math.sqrt(abs(atom[cs+4]))
3468	RenderSphere(x,y,z,radius,atColor)
3469	elif 'lines' in atom[cs]:
3470	radius = 0.1
3471	RenderLines(x,y,z,Bonds,bndColor)
3472	# RenderBonds(x,y,z,Bonds,0.05,color,6)
3473	elif atom[cs] == 'sticks':
3474	radius = 0.1
3475	RenderBonds(x,y,z,Bonds,bondR,bndColor)
3476	elif atom[cs] == 'polyhedra':
3477	RenderPolyhedra(x,y,z,Faces,atColor)
3478	elif atom[cs] == 'backbone':
3479	if atom[ct-1].split()[0] in ['C','N']:
3480	if atom[2] not in Backbones:
3481	Backbones[atom[2]] = []
3482	Backbones[atom[2]].append(list(np.inner(Amat,np.array([x,y,z]))))
3483	BackboneColor.append(list(atColor))
3484
3485	if atom[cs+1] == 'type':
3486	RenderLabel(x,y,z,' '+atom[ct],radius,Gr)
3487	elif atom[cs+1] == 'name':
3488	RenderLabel(x,y,z,' '+atom[ct-1],radius,Gr)
3489	elif atom[cs+1] == 'number':
3490	RenderLabel(x,y,z,' '+str(iat),radius,Gr)
3491	elif atom[cs+1] == 'residue' and atom[ct-1] == 'CA':
3492	RenderLabel(x,y,z,' '+atom[ct-4],radius,Gr)
3493	elif atom[cs+1] == '1-letter' and atom[ct-1] == 'CA':
3494	RenderLabel(x,y,z,' '+atom[ct-3],radius,Gr)
3495	elif atom[cs+1] == 'chain' and atom[ct-1] == 'CA':
3496	RenderLabel(x,y,z,' '+atom[ct-2],radius,Gr)
3497	# glDisable(GL_BLEND)
3498	if len(rhoXYZ):
3499	RenderMap(rho,rhoXYZ,indx,Rok)
3500	if len(mapPeaks):
3501	XYZ = mapPeaks.T[1:4].T
3502	mapBonds = FindPeaksBonds(XYZ)
3503	for ind,[mag,x,y,z,d] in enumerate(mapPeaks):
3504	if ind in Ind and pageName == 'Map peaks':
3505	RenderMapPeak(x,y,z,Gr,1.0)
3506	else:
3507	RenderMapPeak(x,y,z,Wt,mag/peakMax)
3508	if showBonds:
3509	RenderLines(x,y,z,mapBonds[ind],Wt)
3510	if len(testRBObj) and pageName == 'RB Models':
3511	XYZ = G2mth.UpdateRBXYZ(Bmat,testRBObj['rbObj'],testRBObj['rbData'],testRBObj['rbType'])[0]
3512	rbBonds = FindPeaksBonds(XYZ)
3513	for ind,[x,y,z] in enumerate(XYZ):
3514	aType = testRBObj['rbAtTypes'][ind]
3515	name = ' '+aType+str(ind)
3516	color = np.array(testRBObj['AtInfo'][aType][1])
3517	RenderSphere(x,y,z,0.2,color/255.)
3518	RenderBonds(x,y,z,rbBonds[ind],0.03,Gr)
3519	RenderLabel(x,y,z,name,0.2,Or)
3520	if len(mcsaModels) > 1 and pageName == 'MC/SA': #skip the default MD entry
3521	for ind,[x,y,z] in enumerate(mcsaXYZ):
3522	aType = mcsaTypes[ind]
3523	name = ' '+aType+str(ind)
3524	color = np.array(MCSA['AtInfo'][aType][1])
3525	RenderSphere(x,y,z,0.2,color/255.)
3526	RenderBonds(x,y,z,mcsaBonds[ind],0.03,Gr)
3527	RenderLabel(x,y,z,name,0.2,Or)
3528	if Backbones:
3529	for chain in Backbones:
3530	Backbone = Backbones[chain]
3531	RenderBackbone(Backbone,BackboneColor,bondR)
3532	# print time.time()-time0
3533	Page.canvas.SwapBuffers()
3534
3535	def OnSize(event):
3536	Draw('size')
3537
3538	def OnFocus(event): #not needed?? Bind commented out below
3539	Draw('focus')
3540
3541	try:
3542	plotNum = G2frame.G2plotNB.plotList.index(generalData['Name'])
3543	Page = G2frame.G2plotNB.nb.GetPage(plotNum)
3544	except ValueError:
3545	Plot = G2frame.G2plotNB.addOgl(generalData['Name'])
3546	plotNum = G2frame.G2plotNB.plotList.index(generalData['Name'])
3547	Page = G2frame.G2plotNB.nb.GetPage(plotNum)
3548	Page.views = False
3549	view = False
3550	altDown = False
3551	Font = Page.GetFont()
3552	Page.SetFocus()
3553	Page.Choice = None
3554	if mapData['Flip']:
3555	choice = [' save as/key:','jpeg','tiff','bmp','u: roll up','d: roll down','l: roll left','r: roll right']
3556	else:
3557	choice = [' save as/key:','jpeg','tiff','bmp','c: center on 1/2,1/2,1/2','n: next','p: previous']
3558	cb = wx.ComboBox(G2frame.G2plotNB.status,style=wx.CB_DROPDOWN\|wx.CB_READONLY,choices=choice)
3559	cb.Bind(wx.EVT_COMBOBOX, OnKeyBox)
3560	cb.SetValue(' save as/key:')
3561	Page.canvas.Bind(wx.EVT_MOUSEWHEEL, OnMouseWheel)
3562	Page.canvas.Bind(wx.EVT_LEFT_DOWN, OnMouseDown)
3563	Page.canvas.Bind(wx.EVT_RIGHT_DOWN, OnMouseDown)
3564	Page.canvas.Bind(wx.EVT_MIDDLE_DOWN, OnMouseDown)
3565	Page.canvas.Bind(wx.EVT_KEY_UP, OnKey)
3566	Page.canvas.Bind(wx.EVT_MOTION, OnMouseMove)
3567	Page.canvas.Bind(wx.EVT_SIZE, OnSize)
3568	# Page.canvas.Bind(wx.EVT_SET_FOCUS, OnFocus)
3569	Page.camera['position'] = drawingData['cameraPos']
3570	Page.camera['viewPoint'] = np.inner(Amat,drawingData['viewPoint'][0])
3571	Page.camera['backColor'] = np.array(list(drawingData['backColor'])+[0,])/255.
3572	try:
3573	Page.canvas.SetCurrent()
3574	except:
3575	pass
3576	Draw('main')
3577
3578	################################################################################
3579	#### Plot Rigid Body
3580	################################################################################
3581
3582	def PlotRigidBody(G2frame,rbType,AtInfo,rbData,defaults):
3583	'''RB plotting package. Can show rigid body structures as balls & sticks
3584	'''
3585
3586	def FindBonds(XYZ):
3587	rbTypes = rbData['rbTypes']
3588	Radii = []
3589	for Atype in rbTypes:
3590	Radii.append(AtInfo[Atype][0])
3591	if Atype == 'H':
3592	Radii[-1] = 0.5
3593	Radii = np.array(Radii)
3594	Bonds = [[] for i in range(len(Radii))]
3595	for i,xyz in enumerate(XYZ):
3596	Dx = XYZ-xyz
3597	dist = np.sqrt(np.sum(Dx**2,axis=1))
3598	sumR = Radii[i]+Radii
3599	IndB = ma.nonzero(ma.masked_greater(dist-0.85*sumR,0.))
3600	for j in IndB[0]:
3601	Bonds[i].append(Dx[j]*Radii[i]/sumR[j])
3602	Bonds[j].append(-Dx[j]*Radii[j]/sumR[j])
3603	return Bonds
3604
3605	Wt = np.array([255,255,255])
3606	Rd = np.array([255,0,0])
3607	Gr = np.array([0,255,0])
3608	Bl = np.array([0,0,255])
3609	uBox = np.array([[0,0,0],[1,0,0],[0,1,0],[0,0,1]])
3610	uEdges = np.array([[uBox[0],uBox[1]],[uBox[0],uBox[2]],[uBox[0],uBox[3]]])
3611	uColors = [Rd,Gr,Bl]
3612	if rbType == 'Vector':
3613	atNames = [str(i)+':'+Ty for i,Ty in enumerate(rbData['rbTypes'])]
3614	XYZ = np.array([[0.,0.,0.] for Ty in rbData['rbTypes']])
3615	for imag,mag in enumerate(rbData['VectMag']):
3616	XYZ += mag*rbData['rbVect'][imag]
3617	Bonds = FindBonds(XYZ)
3618	elif rbType == 'Residue':
3619	# atNames = [str(i)+':'+Ty for i,Ty in enumerate(rbData['atNames'])]
3620	atNames = rbData['atNames']
3621	XYZ = np.copy(rbData['rbXYZ']) #don't mess with original!
3622	Seq = rbData['rbSeq']
3623	for ia,ib,ang,mv in Seq:
3624	va = XYZ[ia]-XYZ[ib]
3625	Q = G2mth.AVdeg2Q(ang,va)
3626	for im in mv:
3627	vb = XYZ[im]-XYZ[ib]
3628	vb = G2mth.prodQVQ(Q,vb)
3629	XYZ[im] = XYZ[ib]+vb
3630	Bonds = FindBonds(XYZ)
3631	elif rbType == 'Z-matrix':
3632	pass
3633
3634	# def SetRBOrigin():
3635	# page = getSelection()
3636	# if page:
3637	# if G2frame.dataDisplay.GetPageText(page) == 'Rigid bodies':
3638	# G2frame.MapPeaksTable.SetData(mapPeaks)
3639	# panel = G2frame.dataDisplay.GetPage(page).GetChildren()
3640	# names = [child.GetName() for child in panel]
3641	# panel[names.index('grid window')].Refresh()
3642
3643	def OnMouseDown(event):
3644	xy = event.GetPosition()
3645	defaults['oldxy'] = list(xy)
3646
3647	def OnMouseMove(event):
3648	newxy = event.GetPosition()
3649
3650	if event.Dragging():
3651	if event.LeftIsDown():
3652	SetRotation(newxy)
3653	Q = defaults['Quaternion']
3654	G2frame.G2plotNB.status.SetStatusText('New quaternion: %.2f+, %.2fi+ ,%.2fj+, %.2fk'%(Q[0],Q[1],Q[2],Q[3]),1)
3655	# elif event.RightIsDown():
3656	# SetRBOrigin(newxy)
3657	elif event.MiddleIsDown():
3658	SetRotationZ(newxy)
3659	Q = defaults['Quaternion']
3660	G2frame.G2plotNB.status.SetStatusText('New quaternion: %.2f+, %.2fi+ ,%.2fj+, %.2fk'%(Q[0],Q[1],Q[2],Q[3]),1)
3661	Draw('move')
3662
3663	def OnMouseWheel(event):
3664	defaults['cameraPos'] += event.GetWheelRotation()/24
3665	defaults['cameraPos'] = max(10,min(500,defaults['cameraPos']))
3666	G2frame.G2plotNB.status.SetStatusText('New camera distance: %.2f'%(defaults['cameraPos']),1)
3667	Draw('wheel')
3668
3669	def SetBackground():
3670	R,G,B,A = Page.camera['backColor']
3671	glClearColor(R,G,B,A)
3672	glClear(GL_COLOR_BUFFER_BIT \| GL_DEPTH_BUFFER_BIT)
3673
3674	def SetLights():
3675	glEnable(GL_DEPTH_TEST)
3676	glShadeModel(GL_FLAT)
3677	glEnable(GL_LIGHTING)
3678	glEnable(GL_LIGHT0)
3679	glLightModeli(GL_LIGHT_MODEL_TWO_SIDE,0)
3680	glLightfv(GL_LIGHT0,GL_AMBIENT,[1,1,1,.8])
3681	glLightfv(GL_LIGHT0,GL_DIFFUSE,[1,1,1,1])
3682
3683	# def SetRBOrigin(newxy):
3684	##first get translation vector in screen coords.
3685	# oldxy = defaults['oldxy']
3686	# if not len(oldxy): oldxy = list(newxy)
3687	# dxy = newxy-oldxy
3688	# defaults['oldxy'] = list(newxy)
3689	# V = np.array([dxy[0],-dxy[1],0.])/100.
3690	# Q = defaults['Quaternion']
3691	# V = G2mth.prodQVQ(G2mth.invQ(Q),V)
3692	# rbData['rbXYZ'] += V
3693	# PlotRigidBody(G2frame,rbType,AtInfo,rbData,defaults)
3694	#
3695	def SetRotation(newxy):
3696	#first get rotation vector in screen coords. & angle increment
3697	oldxy = defaults['oldxy']
3698	if not len(oldxy): oldxy = list(newxy)
3699	dxy = newxy-oldxy
3700	defaults['oldxy'] = list(newxy)
3701	V = np.array([dxy[1],dxy[0],0.])
3702	A = 0.25np.sqrt(dxy[0]2+dxy[1]*2)
3703	# next transform vector back to xtal coordinates via inverse quaternion
3704	# & make new quaternion
3705	Q = defaults['Quaternion']
3706	V = G2mth.prodQVQ(G2mth.invQ(Q),V)
3707	DQ = G2mth.AVdeg2Q(A,V)
3708	Q = G2mth.prodQQ(Q,DQ)
3709	defaults['Quaternion'] = Q
3710	# finally get new view vector - last row of rotation matrix
3711	VD = G2mth.Q2Mat(Q)[2]
3712	VD /= np.sqrt(np.sum(VD**2))
3713	defaults['viewDir'] = VD
3714
3715	def SetRotationZ(newxy):
3716	#first get rotation vector (= view vector) in screen coords. & angle increment
3717	View = glGetIntegerv(GL_VIEWPORT)
3718	cent = [View[2]/2,View[3]/2]
3719	oldxy = defaults['oldxy']
3720	if not len(oldxy): oldxy = list(newxy)
3721	dxy = newxy-oldxy
3722	defaults['oldxy'] = list(newxy)
3723	V = defaults['viewDir']
3724	A = [0,0]
3725	A[0] = dxy[1]*.25
3726	A[1] = dxy[0]*.25
3727	if newxy[0] > cent[0]:
3728	A[0] *= -1
3729	if newxy[1] < cent[1]:
3730	A[1] *= -1
3731	# next transform vector back to xtal coordinates & make new quaternion
3732	Q = defaults['Quaternion']
3733	Qx = G2mth.AVdeg2Q(A[0],V)
3734	Qy = G2mth.AVdeg2Q(A[1],V)
3735	Q = G2mth.prodQQ(Q,Qx)
3736	Q = G2mth.prodQQ(Q,Qy)
3737	defaults['Quaternion'] = Q
3738
3739	def RenderUnitVectors(x,y,z):
3740	xyz = np.array([x,y,z])
3741	glEnable(GL_COLOR_MATERIAL)
3742	glLineWidth(1)
3743	glPushMatrix()
3744	glTranslate(x,y,z)
3745	glBegin(GL_LINES)
3746	for line,color in zip(uEdges,uColors):
3747	glColor3ubv(color)
3748	glVertex3fv(-line[1])
3749	glVertex3fv(line[1])
3750	glEnd()
3751	glPopMatrix()
3752	glColor4ubv([0,0,0,0])
3753	glDisable(GL_COLOR_MATERIAL)
3754
3755	def RenderSphere(x,y,z,radius,color):
3756	glMaterialfv(GL_FRONT_AND_BACK,GL_DIFFUSE,color)
3757	glPushMatrix()
3758	glTranslate(x,y,z)
3759	q = gluNewQuadric()
3760	gluSphere(q,radius,20,10)
3761	glPopMatrix()
3762
3763	def RenderBonds(x,y,z,Bonds,radius,color,slice=20):
3764	glMaterialfv(GL_FRONT_AND_BACK,GL_DIFFUSE,color)
3765	glPushMatrix()
3766	glTranslate(x,y,z)
3767	for Dx in Bonds:
3768	glPushMatrix()
3769	Z = np.sqrt(np.sum(Dx**2))
3770	if Z:
3771	azm = atan2d(-Dx[1],-Dx[0])
3772	phi = acosd(Dx[2]/Z)
3773	glRotate(-azm,0,0,1)
3774	glRotate(phi,1,0,0)
3775	q = gluNewQuadric()
3776	gluCylinder(q,radius,radius,Z,slice,2)
3777	glPopMatrix()
3778	glPopMatrix()
3779
3780	def RenderLabel(x,y,z,label):
3781	glPushMatrix()
3782	glTranslate(x,y,z)
3783	glDisable(GL_LIGHTING)
3784	glColor3f(1.0,1.0,1.0)
3785	glRasterPos3f(0,0,0)
3786	if bool(glutBitmapCharacter):
3787	for c in list(label):
3788	glutBitmapCharacter(GLUT_BITMAP_8_BY_13,ord(c))
3789	else:
3790	text = gltext.TextElement(text=label,font=Font)
3791	text.draw_text(scale=0.025)
3792	glEnable(GL_LIGHTING)
3793	glPopMatrix()
3794
3795	def Draw(caller=''):
3796	#useful debug?
3797	# if caller:
3798	# print caller
3799	# end of useful debug
3800	cPos = defaults['cameraPos']
3801	VS = np.array(Page.canvas.GetSize())
3802	aspect = float(VS[0])/float(VS[1])
3803	Zclip = 500.0
3804	Q = defaults['Quaternion']
3805	SetBackground()
3806	glInitNames()
3807	glPushName(0)
3808
3809	glMatrixMode(GL_PROJECTION)
3810	glLoadIdentity()
3811	glViewport(0,0,VS[0],VS[1])
3812	gluPerspective(20.,aspect,0.1,500.)
3813	gluLookAt(0,0,cPos,0,0,0,0,1,0)
3814	SetLights()
3815
3816	glMatrixMode(GL_MODELVIEW)
3817	glLoadIdentity()
3818	matRot = G2mth.Q2Mat(Q)
3819	matRot = np.concatenate((np.concatenate((matRot,[[0],[0],[0]]),axis=1),[[0,0,0,1],]),axis=0)
3820	glMultMatrixf(matRot.T)
3821	RenderUnitVectors(0.,0.,0.)
3822	radius = 0.2
3823	for iat,atom in enumerate(XYZ):
3824	x,y,z = atom
3825	CL = AtInfo[rbData['rbTypes'][iat]][1]
3826	color = np.array(CL)/255.
3827	RenderSphere(x,y,z,radius,color)
3828	RenderBonds(x,y,z,Bonds[iat],0.05,color)
3829	RenderLabel(x,y,z,' '+atNames[iat])
3830	Page.canvas.SwapBuffers()
3831
3832	def OnSize(event):
3833	Draw('size')
3834
3835	try:
3836	plotNum = G2frame.G2plotNB.plotList.index('Rigid body')
3837	Page = G2frame.G2plotNB.nb.GetPage(plotNum)
3838	except ValueError:
3839	Plot = G2frame.G2plotNB.addOgl('Rigid body')
3840	plotNum = G2frame.G2plotNB.plotList.index('Rigid body')
3841	Page = G2frame.G2plotNB.nb.GetPage(plotNum)
3842	Page.views = False
3843	view = False
3844	altDown = False
3845	Page.SetFocus()
3846	Font = Page.GetFont()
3847	Page.canvas.Bind(wx.EVT_MOUSEWHEEL, OnMouseWheel)
3848	Page.canvas.Bind(wx.EVT_LEFT_DOWN, OnMouseDown)
3849	Page.canvas.Bind(wx.EVT_RIGHT_DOWN, OnMouseDown)
3850	Page.canvas.Bind(wx.EVT_MIDDLE_DOWN, OnMouseDown)
3851	Page.canvas.Bind(wx.EVT_MOTION, OnMouseMove)
3852	Page.canvas.Bind(wx.EVT_SIZE, OnSize)
3853	Page.camera['position'] = defaults['cameraPos']
3854	Page.camera['backColor'] = np.array([0,0,0,0])
3855	Page.canvas.SetCurrent()
3856	Draw('main')

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

Download in other formats: