Context Navigation

source: trunk/GSASIIplot.py @ 985

Last change on this file since 985 was 985, checked in by vondreele, 10 years ago
remove glFreeType & add gltext - works much better implement gltext - works in 64 bit python/opengl needs a bit of work for model rotation
Property svn:eol-style set to `native` Property svn:keywords set to `Date Author Revision URL Id`
File size: 159.6 KB

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

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

Download in other formats: