Context Navigation

source: trunk/GSASIIplot.py @ 936

Last change on this file since 936 was 936, checked in by vondreele, 10 years ago
more MC/SA stuff further drawing optimization
Property svn:eol-style set to `native` Property svn:keywords set to `Date Author Revision URL Id`
File size: 157.5 KB

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

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

Download in other formats: