Context Navigation

source: trunk/GSASIIplot.py @ 1004

Last change on this file since 1004 was 1004, checked in by toby, 10 years ago
recapture focus on plot window so keyboard commands work
Property svn:eol-style set to `native` Property svn:keywords set to `Date Author Revision URL Id`
File size: 160.3 KB

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

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

Download in other formats: