Context Navigation

source: trunk/GSASIIplot.py @ 1379

Last change on this file since 1379 was 1379, checked in by toby, 9 years ago
more wx2.9: partial fix for GL, fails until rotated, etc.
Property svn:eol-style set to `native` Property svn:keywords set to `Date Author Revision URL Id`
File size: 189.0 KB

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