Context Navigation

source: trunk/GSASIIplot.py @ 1205

Last change on this file since 1205 was 1205, checked in by vondreele, 9 years ago
Add models to SASD data tree Add error bar plotting to SASD data Add I*Q^{4 option to log SASD data got rid of the if, else blocks for all key driven toggles in G2plot}
Property svn:eol-style set to `native` Property svn:keywords set to `Date Author Revision URL Id`
File size: 173.1 KB

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

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

Download in other formats: