source: trunk/GSASIIplot.py @ 453

Last change on this file since 453 was 453, checked in by vondreele, 10 years ago

refactor GSASIIgrid
new Hessian based least squares
new GSASIImath.py
work on focus issues

  • Property svn:keywords set to Date Author Revision URL Id
File size: 115.6 KB
Line 
1#GSASII plotting routines
2########### SVN repository information ###################
3# $Date: 2012-01-19 22:56:20 +0000 (Thu, 19 Jan 2012) $
4# $Author: vondreele $
5# $Revision: 453 $
6# $URL: trunk/GSASIIplot.py $
7# $Id: GSASIIplot.py 453 2012-01-19 22:56:20Z vondreele $
8########### SVN repository information ###################
9import math
10import time
11import copy
12import os.path
13import numpy as np
14import numpy.linalg as nl
15import wx
16import wx.aui
17import wx.glcanvas
18import matplotlib as mpl
19import mpl_toolkits.mplot3d.axes3d as mp3d
20import GSASIIpath
21import GSASIIgrid as G2gd
22import GSASIIimage as G2img
23import GSASIIpwd as G2pwd
24import GSASIIIO as G2IO
25import GSASIIpwdGUI as G2pdG
26import GSASIIimgGUI as G2imG
27import GSASIIphsGUI as G2phG
28import GSASIIlattice as G2lat
29import GSASIIspc as G2spc
30import pytexture as ptx
31from  OpenGL.GL import *
32from OpenGL.GLU import *
33from OpenGL.GLUT import *
34from OpenGL.GLE import *
35from matplotlib.backends.backend_wx import _load_bitmap
36from matplotlib.backends.backend_wxagg import FigureCanvasWxAgg as Canvas
37from matplotlib.backends.backend_wxagg import NavigationToolbar2Wx as Toolbar
38
39# useful degree trig functions
40sind = lambda x: math.sin(x*math.pi/180.)
41cosd = lambda x: math.cos(x*math.pi/180.)
42tand = lambda x: math.tan(x*math.pi/180.)
43asind = lambda x: 180.*math.asin(x)/math.pi
44acosd = lambda x: 180.*math.acos(x)/math.pi
45atan2d = lambda x,y: 180.*math.atan2(y,x)/math.pi
46atand = lambda x: 180.*math.atan(x)/math.pi
47# numpy versions
48npsind = lambda x: np.sin(x*np.pi/180.)
49npcosd = lambda x: np.cos(x*np.pi/180.)
50npacosd = lambda x: 180.*np.arccos(x)/np.pi
51npasind = lambda x: 180.*np.arcsin(x)/np.pi
52npatand = lambda x: 180.*np.arctan(x)/np.pi
53npatan2d = lambda x,y: 180.*np.arctan2(x,y)/np.pi
54   
55class G2PlotMpl(wx.Panel):   
56    def __init__(self,parent,id=-1,dpi=None,**kwargs):
57        wx.Panel.__init__(self,parent,id=id,**kwargs)
58        mpl.rcParams['legend.fontsize'] = 10
59        self.figure = mpl.figure.Figure(dpi=dpi,figsize=(5,7))
60        self.canvas = Canvas(self,-1,self.figure)
61        self.toolbar = GSASIItoolbar(self.canvas)
62
63        self.toolbar.Realize()
64       
65        sizer=wx.BoxSizer(wx.VERTICAL)
66        sizer.Add(self.canvas,1,wx.EXPAND)
67        sizer.Add(self.toolbar,0,wx.LEFT|wx.EXPAND)
68        self.SetSizer(sizer)
69       
70class G2PlotOgl(wx.Panel):
71    def __init__(self,parent,id=-1,dpi=None,**kwargs):
72        self.figure = wx.Panel.__init__(self,parent,id=id,**kwargs)
73        self.canvas = wx.glcanvas.GLCanvas(self,-1,**kwargs)
74        self.camera = {}
75        sizer=wx.BoxSizer(wx.VERTICAL)
76        sizer.Add(self.canvas,1,wx.EXPAND)
77        self.SetSizer(sizer)
78       
79class G2Plot3D(wx.Panel):
80    def __init__(self,parent,id=-1,dpi=None,**kwargs):
81        wx.Panel.__init__(self,parent,id=id,**kwargs)
82        self.figure = mpl.figure.Figure(dpi=dpi,figsize=(6,6))
83        self.canvas = Canvas(self,-1,self.figure)
84        self.toolbar = GSASIItoolbar(self.canvas)
85
86        self.toolbar.Realize()
87       
88        sizer=wx.BoxSizer(wx.VERTICAL)
89        sizer.Add(self.canvas,1,wx.EXPAND)
90        sizer.Add(self.toolbar,0,wx.LEFT|wx.EXPAND)
91        self.SetSizer(sizer)
92                             
93class G2PlotNoteBook(wx.Panel):
94    def __init__(self,parent,id=-1):
95        wx.Panel.__init__(self,parent,id=id)
96        #so one can't delete a plot page!!
97        self.nb = wx.aui.AuiNotebook(self, \
98            style=wx.aui.AUI_NB_DEFAULT_STYLE ^ wx.aui.AUI_NB_CLOSE_ON_ACTIVE_TAB)
99        sizer = wx.BoxSizer()
100        sizer.Add(self.nb,1,wx.EXPAND)
101        self.SetSizer(sizer)
102        self.status = parent.CreateStatusBar()
103        self.status.SetFieldsCount(2)
104        self.status.SetStatusWidths([150,-1])
105        self.Bind(wx.aui.EVT_AUINOTEBOOK_PAGE_CHANGED, self.OnPageChanged)
106       
107        self.plotList = []
108           
109    def addMpl(self,name=""):
110        page = G2PlotMpl(self.nb)
111        self.nb.AddPage(page,name)
112       
113        self.plotList.append(name)
114       
115        return page.figure
116       
117    def add3D(self,name=""):
118        page = G2Plot3D(self.nb)
119        self.nb.AddPage(page,name)
120       
121        self.plotList.append(name)
122       
123        return page.figure
124       
125    def addOgl(self,name=""):
126        page = G2PlotOgl(self.nb)
127        self.nb.AddPage(page,name)
128       
129        self.plotList.append(name)
130       
131        return page.figure
132       
133    def Delete(self,name):
134        try:
135            item = self.plotList.index(name)
136            del self.plotList[item]
137            self.nb.DeletePage(item)
138        except ValueError:          #no plot of this name - do nothing
139            return     
140               
141    def clear(self):
142        while self.nb.GetPageCount():
143            self.nb.DeletePage(0)
144        self.plotList = []
145        self.status.DestroyChildren()
146       
147    def Rename(self,oldName,newName):
148        try:
149            item = self.plotList.index(oldName)
150            self.plotList[item] = newName
151            self.nb.SetPageText(item,newName)
152        except ValueError:          #no plot of this name - do nothing
153            return     
154       
155    def OnPageChanged(self,event):       
156        if self.plotList:
157            self.status.SetStatusText('Better to select this from GSAS-II data tree',1)
158        self.status.DestroyChildren()                           #get rid of special stuff on status bar
159       
160class GSASIItoolbar(Toolbar):
161    ON_MPL_HELP = wx.NewId()
162    def __init__(self,plotCanvas):
163        Toolbar.__init__(self,plotCanvas)
164        POSITION_OF_CONFIGURE_SUBPLOTS_BTN = 6
165        self.DeleteToolByPos(POSITION_OF_CONFIGURE_SUBPLOTS_BTN)
166        help = os.path.join(os.path.split(sys.argv[0])[0],'help.ico')
167        self.AddSimpleTool(self.ON_MPL_HELP,_load_bitmap(help),'Help on','Show help on')
168        wx.EVT_TOOL(self,self.ON_MPL_HELP,self.OnHelp)
169    def OnHelp(self,event):
170        Page = self.GetParent().GetParent()
171        pageNo = Page.GetSelection()
172        print 'Matplotlib help on '+Page.GetPageText(pageNo)
173        G2gd.ShowHelp(Page.GetPageText(pageNo),self.TopLevelParent)
174           
175def PlotSngl(self,newPlot=False):
176    '''Single crystal structure factor plotting package - displays zone of reflections as rings proportional
177        to F, F**2, etc. as requested
178    '''
179    from matplotlib.patches import Circle
180    global HKL,HKLF
181
182    def OnSCMotion(event):
183        xpos = event.xdata
184        if xpos:
185            xpos = round(xpos)                                        #avoid out of frame mouse position
186            ypos = round(event.ydata)
187            zpos = Data['Layer']
188            if '100' in Data['Zone']:
189                HKLtxt = '(%3d,%3d,%3d)'%(zpos,xpos,ypos)
190            elif '010' in Data['Zone']:
191                HKLtxt = '(%3d,%3d,%3d)'%(xpos,zpos,ypos)
192            elif '001' in Data['Zone']:
193                HKLtxt = '(%3d,%3d,%3d)'%(xpos,ypos,zpos)
194            Page.canvas.SetToolTipString(HKLtxt)
195            self.G2plotNB.status.SetFields(['HKL = '+HKLtxt,''])
196               
197    def OnSCPick(event):
198        zpos = Data['Layer']
199        pos = event.artist.center
200        if '100' in Data['Zone']:
201            Page.canvas.SetToolTipString('(picked:(%3d,%3d,%3d))'%(zpos,pos[0],pos[1]))
202            hkl = np.array([zpos,pos[0],pos[1]])
203        elif '010' in Data['Zone']:
204            Page.canvas.SetToolTipString('(picked:(%3d,%3d,%3d))'%(pos[0],zpos,pos[1]))
205            hkl = np.array([pos[0],zpos,pos[1]])
206        elif '001' in Data['Zone']:
207            Page.canvas.SetToolTipString('(picked:(%3d,%3d,%3d))'%(pos[0],pos[1],zpos))
208            hkl = np.array([pos[0],pos[1],zpos])
209        h,k,l = hkl
210        hklf = HKLF[np.where(np.all(HKL-hkl == [0,0,0],axis=1))]
211        if len(hklf):
212            Fosq,sig,Fcsq = hklf[0]
213            HKLtxt = '(%3d,%3d,%3d %.2f %.3f %.2f %.2f)'%(h,k,l,Fosq,sig,Fcsq,(Fosq-Fcsq)/(scale*sig))
214            self.G2plotNB.status.SetFields(['','HKL, Fosq, sig, Fcsq, delFsq/sig = '+HKLtxt])
215                                 
216    def OnSCKeyPress(event):
217        print event.key
218
219    try:
220        plotNum = self.G2plotNB.plotList.index('Structure Factors')
221        Page = self.G2plotNB.nb.GetPage(plotNum)
222        if not newPlot:
223            Plot = Page.figure.gca()          #get previous powder plot & get limits
224            xylim = Plot.get_xlim(),Plot.get_ylim()
225        Page.figure.clf()
226        Plot = Page.figure.gca()          #get a fresh plot after clf()
227    except ValueError:
228        Plot = self.G2plotNB.addMpl('Structure Factors').gca()
229        plotNum = self.G2plotNB.plotList.index('Structure Factors')
230        Page = self.G2plotNB.nb.GetPage(plotNum)
231#        Page.canvas.mpl_connect('key_press_event', OnSCKeyPress)
232        Page.canvas.mpl_connect('pick_event', OnSCPick)
233        Page.canvas.mpl_connect('motion_notify_event', OnSCMotion)
234    Page.SetFocus()
235   
236    Plot.set_aspect(aspect='equal')
237    HKLref = self.PatternTree.GetItemPyData(self.Sngl)
238    Data = self.PatternTree.GetItemPyData( \
239        G2gd.GetPatternTreeItemId(self,self.Sngl, 'HKL Plot Controls'))
240    Type = Data['Type']           
241    scale = Data['Scale']
242    HKLmax = Data['HKLmax']
243    HKLmin = Data['HKLmin']
244    FosqMax = Data['FoMax']
245    FoMax = math.sqrt(FosqMax)
246    ifFc = Data['ifFc']
247    xlabel = ['k, h=','h, k=','h, l=']
248    ylabel = ['l','l','k']
249    zones = ['100','010','001']
250    pzone = [[1,2],[0,2],[0,1]]
251    izone = zones.index(Data['Zone'])
252    Plot.set_title(self.PatternTree.GetItemText(self.Sngl)[5:])
253    HKL = []
254    HKLF = []
255    for H,Fosq,sig,Fcsq,x,x,x in HKLref:
256        HKL.append(H)
257        HKLF.append([Fosq,sig,Fcsq])
258        if H[izone] == Data['Layer']:
259            B = 0
260            if Type == 'Fosq':
261                A = scale*Fosq/FosqMax
262                B = scale*Fcsq/FosqMax
263                C = abs(A-B)
264            elif Type == 'Fo':
265                A = scale*math.sqrt(max(0,Fosq))/FoMax
266                B = scale*math.sqrt(max(0,Fcsq))/FoMax
267                C = abs(A-B)
268            elif Type == '|DFsq|/sig':
269                A = abs(Fosq-Fcsq)/(scale*sig)
270            elif Type == '|DFsq|>sig':
271                A = abs(Fosq-Fcsq)/(scale*sig)
272                if A < 1.0: A = 0                   
273            elif Type == '|DFsq|>3sig':
274                A = abs(Fosq-Fcsq)/(scale*sig)
275                if A < 3.0: A = 0                   
276            xy = (H[pzone[izone][0]],H[pzone[izone][1]])
277            if A > 0.0:
278                Plot.add_artist(Circle(xy,radius=A,ec='g',fc='w',picker=3))
279            if B:
280                Plot.add_artist(Circle(xy,radius=B,ec='b',fc='w'))
281                radius = C
282                if radius > 0:
283                    if A > B:
284                        Plot.add_artist(Circle(xy,radius=radius,ec='g',fc='g'))
285                    else:                   
286                        Plot.add_artist(Circle(xy,radius=radius,ec='r',fc='r'))
287    HKL = np.array(HKL,dtype=np.int)
288    HKLF = np.array(HKLF)
289    Plot.set_xlabel(xlabel[izone]+str(Data['Layer']),fontsize=12)
290    Plot.set_ylabel(ylabel[izone],fontsize=12)
291    Plot.set_xlim((HKLmin[pzone[izone][0]],HKLmax[pzone[izone][0]]))
292    Plot.set_ylim((HKLmin[pzone[izone][1]],HKLmax[pzone[izone][1]]))
293    if not newPlot:
294        Page.toolbar.push_current()
295        Plot.set_xlim(xylim[0])
296        Plot.set_ylim(xylim[1])
297        xylim = []
298        Page.toolbar.push_current()
299        Page.toolbar.draw()
300    else:
301        Page.canvas.draw()
302       
303def PlotPatterns(self,newPlot=False):
304    '''Powder pattern plotting package - displays single or multiple powder patterns as intensity vs
305    2-theta or q (future TOF). Can display multiple patterns as "waterfall plots" or contour plots. Log I
306    plotting available.
307    '''
308    global HKL
309   
310    def OnKeyBox(event):
311        if self.G2plotNB.nb.GetSelection() == self.G2plotNB.plotList.index('Powder Patterns'):
312            event.key = cb.GetValue()[0]
313            cb.SetValue(' key press')
314            OnPlotKeyPress(event)
315                       
316    def OnPlotKeyPress(event):
317        newPlot = False
318        if event.key == 'w':
319            if self.Weight:
320                self.Weight = False
321            else:
322                self.Weight = True
323            print 'plot weighting:',self.Weight
324        elif event.key == 'n':
325            if self.Contour:
326                pass
327            else:
328                if self.logPlot:
329                    self.logPlot = False
330                else:
331                    self.Offset[0] = 0
332                    self.logPlot = True
333        elif event.key == 'u':
334            if self.Contour:
335                self.Cmax = min(1.0,self.Cmax*1.2)
336            elif self.logPlot:
337                pass
338            elif self.Offset[0] < 100.:
339                self.Offset[0] += 1.
340        elif event.key == 'd':
341            if self.Contour:
342                self.Cmax = max(0.0,self.Cmax*0.8)
343            elif self.logPlot:
344                pass
345            elif self.Offset[0] > 0.:
346                self.Offset[0] -= 1.
347        elif event.key == 'l':
348            self.Offset[1] -= 1.
349        elif event.key == 'r':
350            self.Offset[1] += 1.
351        elif event.key == 'o':
352            self.Offset = [0,0]
353        elif event.key == 'c':
354            newPlot = True
355            if self.Contour:
356                self.Contour = False
357            else:
358                self.Contour = True
359                self.SinglePlot = False
360                self.Offset = [0.,0.]
361        elif event.key == 'q':
362            newPlot = True
363            if self.qPlot:
364                self.qPlot = False
365            else:
366                self.qPlot = True
367        elif event.key == 's':
368            if self.Contour:
369                choice = [m for m in mpl.cm.datad.keys() if not m.endswith("_r")]
370                choice.sort()
371                dlg = wx.SingleChoiceDialog(self,'Select','Color scheme',choice)
372                if dlg.ShowModal() == wx.ID_OK:
373                    sel = dlg.GetSelection()
374                    self.ContourColor = choice[sel]
375                else:
376                    self.ContourColor = 'Paired'
377                dlg.Destroy()
378            else:               
379                if self.SinglePlot:
380                    self.SinglePlot = False
381                else:
382                    self.SinglePlot = True
383        elif event.key == '+':
384            if self.PickId:
385                self.PickId = False
386        elif event.key == 'i':                  #for smoothing contour plot
387            choice = ['nearest','bilinear','bicubic','spline16','spline36','hanning',
388               'hamming','hermite','kaiser','quadric','catrom','gaussian','bessel',
389               'mitchell','sinc','lanczos']
390            dlg = wx.SingleChoiceDialog(self,'Select','Interpolation',choice)
391            if dlg.ShowModal() == wx.ID_OK:
392                sel = dlg.GetSelection()
393                self.Interpolate = choice[sel]
394            else:
395                self.Interpolate = 'nearest'
396            dlg.Destroy()
397           
398        PlotPatterns(self,newPlot=newPlot)
399       
400    def OnMotion(event):
401        xpos = event.xdata
402        if xpos:                                        #avoid out of frame mouse position
403            ypos = event.ydata
404            Page.canvas.SetCursor(wx.CROSS_CURSOR)
405            try:
406                Values,Names = self.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(self,self.PatternId, 'Instrument Parameters'))[1::2]
407                Parms = dict(zip(Names,Values))
408                try:
409                    wave = Parms['Lam']
410                except KeyError:
411                    wave = Parms['Lam1']
412                if self.qPlot:
413                    xpos = 2.0*asind(xpos*wave/(4*math.pi))
414                dsp = 0.0
415                if abs(xpos) > 0.:                  #avoid possible singularity at beam center
416                    dsp = wave/(2.*sind(abs(xpos)/2.0))
417                if self.Contour:
418                    self.G2plotNB.status.SetStatusText('2-theta =%9.3f d =%9.5f pattern ID =%5d'%(xpos,dsp,int(ypos)),1)
419                else:
420                    self.G2plotNB.status.SetStatusText('2-theta =%9.3f d =%9.5f Intensity =%9.1f'%(xpos,dsp,ypos),1)
421                if self.itemPicked:
422                    Page.canvas.SetToolTipString('%9.3f'%(xpos))
423                if self.PickId:
424                    found = []
425                    if self.PatternTree.GetItemText(self.PickId) in ['Index Peak List','Unit Cells List','Reflection Lists']:
426                        if len(HKL):
427                            view = Page.toolbar._views.forward()[0][:2]
428                            wid = view[1]-view[0]
429                            found = HKL[np.where(np.fabs(HKL.T[5]-xpos) < 0.002*wid)]
430                        if len(found):
431                            h,k,l = found[0][:3] 
432                            Page.canvas.SetToolTipString('%d,%d,%d'%(int(h),int(k),int(l)))
433                        else:
434                            Page.canvas.SetToolTipString('')
435
436            except TypeError:
437                self.G2plotNB.status.SetStatusText('Select PWDR powder pattern first',1)
438                                                   
439    def OnPick(event):
440        if self.itemPicked is not None: return
441        PatternId = self.PatternId
442        try:
443            Values,Names = self.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(self,self.PatternId, 'Instrument Parameters'))[1::2]
444        except TypeError:
445            return
446        Parms = dict(zip(Names,Values))
447        try:
448            wave = Parms['Lam']
449        except KeyError:
450            wave = Parms['Lam1']
451        PickId = self.PickId
452        pick = event.artist
453        mouse = event.mouseevent       
454        xpos = pick.get_xdata()
455        ypos = pick.get_ydata()
456        ind = event.ind
457        xy = list(zip(np.take(xpos,ind),np.take(ypos,ind))[0])
458        if self.PatternTree.GetItemText(PickId) == 'Peak List':
459            if ind.all() != [0]:                                    #picked a data point
460                if 'C' in Parms['Type']:                            #CW data - TOF later in an elif
461                    ins = [Parms[x] for x in ['U','V','W','X','Y']]
462                    if self.qPlot:                              #qplot - convert back to 2-theta
463                        xy[0] = 2.0*asind(xy[0]*wave/(4*math.pi))
464                    sig = ins[0]*tand(xy[0]/2.0)**2+ins[1]*tand(xy[0]/2.0)+ins[2]
465                    gam = ins[3]/cosd(xy[0]/2.0)+ins[4]*tand(xy[0]/2.0)           
466                    data = self.PatternTree.GetItemPyData(self.PickId)
467                    XY = [xy[0],0, xy[1],1, sig,0, gam,0]       #default refine intensity 1st
468                data.append(XY)
469                G2pdG.UpdatePeakGrid(self,data)
470                PlotPatterns(self)
471            else:                                                   #picked a peak list line
472                self.itemPicked = pick
473        elif self.PatternTree.GetItemText(PickId) == 'Limits':
474            if ind.all() != [0]:                                    #picked a data point
475                LimitId = G2gd.GetPatternTreeItemId(self,PatternId, 'Limits')
476                data = self.PatternTree.GetItemPyData(LimitId)
477                if 'C' in Parms['Type']:                            #CW data - TOF later in an elif
478                    if self.qPlot:                              #qplot - convert back to 2-theta
479                        xy[0] = 2.0*asind(xy[0]*wave/(4*math.pi))
480                if mouse.button==1:
481                    data[1][0] = min(xy[0],data[1][1])
482                if mouse.button==3:
483                    data[1][1] = max(xy[0],data[1][0])
484                self.PatternTree.SetItemPyData(LimitId,data)
485                G2pdG.UpdateLimitsGrid(self,data)
486                PlotPatterns(self)
487            else:                                                   #picked a limit line
488                self.itemPicked = pick
489        elif self.PatternTree.GetItemText(PickId) == 'Reflection Lists':
490            self.itemPicked = pick
491            pick = str(pick)
492       
493    def OnRelease(event):
494        if self.itemPicked is None: return
495        Values,Names = self.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(self,self.PatternId, 'Instrument Parameters'))[1::2]
496        Parms = dict(zip(Names,Values))
497        try:
498            wave = Parms['Lam']
499        except KeyError:
500            wave = Parms['Lam1']
501        xpos = event.xdata
502        PickId = self.PickId
503        if self.PatternTree.GetItemText(PickId) in ['Peak List','Limits'] and xpos:
504            lines = []
505            for line in self.Lines: lines.append(line.get_xdata()[0])
506            lineNo = lines.index(self.itemPicked.get_xdata()[0])
507            if  lineNo in [0,1]:
508                LimitId = G2gd.GetPatternTreeItemId(self,self.PatternId, 'Limits')
509                data = self.PatternTree.GetItemPyData(LimitId)
510                if self.qPlot:
511                    data[1][lineNo] = 2.0*asind(wave*xpos/(4*math.pi))
512                else:
513                    data[1][lineNo] = xpos
514                self.PatternTree.SetItemPyData(LimitId,data)
515                if self.PatternTree.GetItemText(self.PickId) == 'Limits':
516                    G2pdG.UpdateLimitsGrid(self,data)
517            else:
518                PeakId = G2gd.GetPatternTreeItemId(self,self.PatternId, 'Peak List')
519                data = self.PatternTree.GetItemPyData(PeakId)
520#                print 'peaks',xpos
521                if event.button == 3:
522                    del data[lineNo-2]
523                else:
524                    if self.qPlot:
525                        data[lineNo-2][0] = 2.0*asind(wave*xpos/(4*math.pi))
526                    else:
527                        data[lineNo-2][0] = xpos
528                self.PatternTree.SetItemPyData(PeakId,data)
529                G2pdG.UpdatePeakGrid(self,data)
530        elif self.PatternTree.GetItemText(PickId) == 'Reflection Lists' and xpos:
531            Phases = self.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(self,PatternId,'Reflection Lists'))
532            pick = str(self.itemPicked).split('(')[1].strip(')')
533            if 'Line' not in pick:       #avoid data points, etc.
534                num = Phases.keys().index(pick)
535                if num:
536                    self.refDelt = -(event.ydata-self.refOffset)/(num*Ymax)
537                else:       #1st row of refl ticks
538                    self.refOffset = event.ydata
539        PlotPatterns(self)
540        self.itemPicked = None   
541
542    xylim = []
543    try:
544        plotNum = self.G2plotNB.plotList.index('Powder Patterns')
545        Page = self.G2plotNB.nb.GetPage(plotNum)
546        if not newPlot:
547            Plot = Page.figure.gca()          #get previous powder plot & get limits
548            xylim = Plot.get_xlim(),Plot.get_ylim()
549        Page.figure.clf()
550        Plot = Page.figure.gca()          #get a fresh plot after clf()
551    except ValueError:
552        newPlot = True
553        self.Cmax = 1.0
554        Plot = self.G2plotNB.addMpl('Powder Patterns').gca()
555        plotNum = self.G2plotNB.plotList.index('Powder Patterns')
556        Page = self.G2plotNB.nb.GetPage(plotNum)
557        Page.canvas.mpl_connect('key_press_event', OnPlotKeyPress)
558        Page.canvas.mpl_connect('motion_notify_event', OnMotion)
559        Page.canvas.mpl_connect('pick_event', OnPick)
560        Page.canvas.mpl_connect('button_release_event', OnRelease)
561    Page.SetFocus()
562    self.G2plotNB.status.DestroyChildren()
563    if self.Contour:
564        Choice = (' key press','d: lower contour max','u: raise contour max',
565            'i: interpolation method','s: color scheme','c: contour off')
566    else:
567        if self.logPlot:
568            Choice = (' key press','n: log(I) off','l: offset left','r: offset right',
569                'c: contour on','q: toggle q plot','s: toggle single plot','+: no selection')
570        else:
571            Choice = (' key press','l: offset left','r: offset right','d: offset down',
572                'u: offset up','o: reset offset','n: log(I) on','c: contour on',
573                'q: toggle q plot','s: toggle single plot','+: no selection')
574    cb = wx.ComboBox(self.G2plotNB.status,style=wx.CB_DROPDOWN|wx.CB_READONLY,
575        choices=Choice)
576    cb.Bind(wx.EVT_COMBOBOX, OnKeyBox)
577    cb.SetValue(' key press')
578   
579    PickId = self.PickId
580    PatternId = self.PatternId
581    colors=['b','g','r','c','m','k']
582    Lines = []
583    if self.SinglePlot:
584        Pattern = self.PatternTree.GetItemPyData(PatternId)
585        Pattern.append(self.PatternTree.GetItemText(PatternId))
586        PlotList = [Pattern,]
587        ParmList = [self.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(self,
588            self.PatternId, 'Instrument Parameters'))[1],]
589    else:       
590        PlotList = []
591        ParmList = []
592        item, cookie = self.PatternTree.GetFirstChild(self.root)
593        while item:
594            if 'PWDR' in self.PatternTree.GetItemText(item):
595                Pattern = self.PatternTree.GetItemPyData(item)
596                if len(Pattern) < 3:                    # put name on end if needed
597                    Pattern.append(self.PatternTree.GetItemText(item))
598                PlotList.append(Pattern)
599                ParmList.append(self.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(self,
600                    item,'Instrument Parameters'))[1])
601            item, cookie = self.PatternTree.GetNextChild(self.root, cookie)               
602    Ymax = 1.0
603    lenX = 0
604    if self.PatternTree.GetItemText(PickId) in ['Reflection Lists']:
605        Phases = self.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(self,PatternId,'Reflection Lists'))
606        HKL = []
607        if Phases:
608            for peak in Phases[self.RefList]:
609                HKL.append(peak[:6])
610            HKL = np.array(HKL)
611    else:
612        HKL = np.array(self.HKL)
613    for Pattern in PlotList:
614        xye = Pattern[1]
615        Ymax = max(Ymax,max(xye[1]))
616    offset = self.Offset[0]*Ymax/100.0
617    Title = 'Powder Patterns: '+os.path.split(self.GSASprojectfile)[1]
618    if self.logPlot:
619        Title = 'log('+Title+')'
620    Plot.set_title(Title)
621    if self.qPlot:
622        Plot.set_xlabel(r'$q, \AA^{-1}$',fontsize=14)
623    else:       
624        Plot.set_xlabel(r'$\mathsf{2\theta}$',fontsize=14)
625    Plot.set_ylabel('Intensity',fontsize=12)
626    if self.Contour:
627        ContourZ = []
628        ContourY = []
629        Nseq = 0
630    if len(PlotList) < 2:
631        self.Contour = False
632    for N,Pattern in enumerate(PlotList):
633        Parms = ParmList[N]
634        ifpicked = False
635        LimitId = 0
636        xye = np.array(Pattern[1])
637        if PickId:
638            ifpicked = Pattern[2] == self.PatternTree.GetItemText(PatternId)
639            LimitId = G2gd.GetPatternTreeItemId(self,PatternId, 'Limits')
640        if self.qPlot:
641            Id = G2gd.GetPatternTreeItemId(self,self.root, Pattern[2])
642            X = 4*np.pi*npsind(xye[0]/2.0)/Parms[1]
643        else:
644            X = xye[0]
645        if not lenX:
646            lenX = len(X)           
647        Y = xye[1]+offset*N
648        if LimitId:
649            limits = np.array(self.PatternTree.GetItemPyData(LimitId))
650            if self.qPlot:
651                limits = 4*np.pi*npsind(limits/2.0)/Parms[1]
652            Lines.append(Plot.axvline(limits[1][0],color='g',dashes=(5,5),picker=3.))   
653            Lines.append(Plot.axvline(limits[1][1],color='r',dashes=(5,5),picker=3.))                   
654        if self.Contour:
655           
656            if lenX == len(X):
657                ContourY.append(N)
658                ContourZ.append(Y)
659                ContourX = X
660                Nseq += 1
661                Plot.set_ylabel('Data sequence',fontsize=12)
662        else:
663            X += self.Offset[1]*.005*N
664            if ifpicked:
665                Z = xye[3]+offset*N
666                W = xye[4]+offset*N
667                D = xye[5]-Ymax*self.delOffset
668                if self.Weight:
669                    W2 = np.sqrt(xye[2])
670                    D *= W2-Ymax*self.delOffset
671                if self.logPlot:
672                    Plot.semilogy(X,Y,colors[N%6]+'+',picker=3.,clip_on=False,nonposy='mask')
673                    Plot.semilogy(X,Z,colors[(N+1)%6],picker=False,nonposy='mask')
674                    Plot.semilogy(X,W,colors[(N+2)%6],picker=False,nonposy='mask')
675                else:
676                    Plot.plot(X,Y,colors[N%6]+'+',picker=3.,clip_on=False)
677                    Plot.plot(X,Z,colors[(N+1)%6],picker=False)
678                    Plot.plot(X,W,colors[(N+2)%6],picker=False)
679                    Plot.plot(X,D,colors[(N+3)%6],picker=False)
680                    Plot.axhline(0.,color=wx.BLACK)
681                Page.canvas.SetToolTipString('')
682                if self.PatternTree.GetItemText(PickId) == 'Peak List':
683                    tip = 'On data point: Pick peak - L or R MB. On line: L-move, R-delete'
684                    Page.canvas.SetToolTipString(tip)
685                    data = self.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(self,PatternId, 'Peak List'))
686                    for item in data:
687                        if self.qPlot:
688                            Lines.append(Plot.axvline(4*math.pi*sind(item[0]/2.)/Parms[1],color=colors[N%6],picker=2.))
689                        else:
690                            Lines.append(Plot.axvline(item[0],color=colors[N%6],picker=2.))
691                if self.PatternTree.GetItemText(PickId) == 'Limits':
692                    tip = 'On data point: Lower limit - L MB; Upper limit - R MB. On limit: MB down to move'
693                    Page.canvas.SetToolTipString(tip)
694                    data = self.LimitsTable.GetData()
695            else:
696                if self.logPlot:
697                    Plot.semilogy(X,Y,colors[N%6],picker=False,nonposy='clip')
698                else:
699                    Plot.plot(X,Y,colors[N%6],picker=False)
700    if PickId and not self.Contour:
701        Values,Names = self.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(self,PatternId, 'Instrument Parameters'))[1::2]
702        Parms = dict(zip(Names,Values))
703        try:
704            wave = Parms['Lam']
705        except KeyError:
706            wave = Parms['Lam1']
707        if self.PatternTree.GetItemText(PickId) in ['Index Peak List','Unit Cells List']:
708            peaks = np.array((self.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(self,PatternId, 'Index Peak List'))))
709            for peak in peaks:
710                if self.qPlot:
711                    Plot.axvline(4*np.pi*sind(peak[0]/2.0)/wave,color='b')
712                else:
713                    Plot.axvline(peak[0],color='b')
714            for hkl in self.HKL:
715                if self.qPlot:
716                    Plot.axvline(4*np.pi*sind(hkl[5]/2.0)/wave,color='r',dashes=(5,5))
717                else:
718                    Plot.axvline(hkl[5],color='r',dashes=(5,5))
719        elif self.PatternTree.GetItemText(PickId) in ['Reflection Lists']:
720            Phases = self.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(self,PatternId,'Reflection Lists'))
721            for pId,phase in enumerate(Phases):
722                peaks = Phases[phase]
723                peak = np.array([[peak[4],peak[5]] for peak in peaks])
724                pos = self.refOffset-pId*Ymax*self.refDelt*np.ones_like(peak)
725                if self.qPlot:
726                    Plot.plot(2*np.pi/peak.T[0],pos,colors[pId%6]+'|',mew=1,ms=8,picker=3.,label=phase)
727                else:
728                    Plot.plot(peak.T[1],pos,colors[pId%6]+'|',mew=1,ms=8,picker=3.,label=phase)
729            handles,legends = Plot.get_legend_handles_labels()  #got double entries in the legends for some reason
730            Plot.legend(handles[::2],legends[::2],title='Phases',loc='best')    #skip every other one
731           
732    if self.Contour:
733        acolor = mpl.cm.get_cmap(self.ContourColor)
734        Img = Plot.imshow(ContourZ,cmap=acolor,vmin=0,vmax=Ymax*self.Cmax,interpolation=self.Interpolate, 
735            extent=[ContourX[0],ContourX[-1],ContourY[0],ContourY[-1]],aspect='auto',origin='lower')
736        Page.figure.colorbar(Img)
737    else:
738        self.Lines = Lines
739    if not newPlot:
740        Page.toolbar.push_current()
741        Plot.set_xlim(xylim[0])
742        Plot.set_ylim(xylim[1])
743        xylim = []
744        Page.toolbar.push_current()
745        Page.toolbar.draw()
746    else:
747        Page.canvas.draw()
748    self.Pwdr = True
749    if self.oldFocus:
750        self.SetFocus()
751   
752def PlotISFG(self,newPlot=False,type=''):
753    ''' PLotting package for PDF analysis; displays I(q), S(q), F(q) and G(r) as single
754    or multiple plots with waterfall and contour plots as options
755    '''
756    if not type:
757        type = self.G2plotNB.plotList[self.G2plotNB.nb.GetSelection()]
758    if type not in ['I(Q)','S(Q)','F(Q)','G(R)']:
759        return
760    superMinusOne = unichr(0xaf)+unichr(0xb9)
761   
762    def OnPlotKeyPress(event):
763        newPlot = False
764        if event.key == 'u':
765            if self.Contour:
766                self.Cmax = min(1.0,self.Cmax*1.2)
767            elif self.Offset[0] < 100.:
768                self.Offset[0] += 1.
769        elif event.key == 'd':
770            if self.Contour:
771                self.Cmax = max(0.0,self.Cmax*0.8)
772            elif self.Offset[0] > 0.:
773                self.Offset[0] -= 1.
774        elif event.key == 'l':
775            self.Offset[1] -= 1.
776        elif event.key == 'r':
777            self.Offset[1] += 1.
778        elif event.key == 'o':
779            self.Offset = [0,0]
780        elif event.key == 'c':
781            newPlot = True
782            if self.Contour:
783                self.Contour = False
784            else:
785                self.Contour = True
786                self.SinglePlot = False
787                self.Offset = [0.,0.]
788        elif event.key == 's':
789            if self.Contour:
790                choice = [m for m in mpl.cm.datad.keys() if not m.endswith("_r")]
791                choice.sort()
792                dlg = wx.SingleChoiceDialog(self,'Select','Color scheme',choice)
793                if dlg.ShowModal() == wx.ID_OK:
794                    sel = dlg.GetSelection()
795                    self.ContourColor = choice[sel]
796                else:
797                    self.ContourColor = 'Paired'
798                dlg.Destroy()
799            else:               
800                if self.SinglePlot:
801                    self.SinglePlot = False
802                else:
803                    self.SinglePlot = True
804        elif event.key == 'i':                  #for smoothing contour plot
805            choice = ['nearest','bilinear','bicubic','spline16','spline36','hanning',
806               'hamming','hermite','kaiser','quadric','catrom','gaussian','bessel',
807               'mitchell','sinc','lanczos']
808            dlg = wx.SingleChoiceDialog(self,'Select','Interpolation',choice)
809            if dlg.ShowModal() == wx.ID_OK:
810                sel = dlg.GetSelection()
811                self.Interpolate = choice[sel]
812            else:
813                self.Interpolate = 'nearest'
814            dlg.Destroy()
815        elif event.key == 't' and not self.Contour:
816            if self.Legend:
817                self.Legend = False
818            else:
819                self.Legend = True
820           
821           
822        PlotISFG(self,newPlot=newPlot,type=type)
823       
824    def OnKeyBox(event):
825        if self.G2plotNB.nb.GetSelection() == self.G2plotNB.plotList.index(type):
826            event.key = cb.GetValue()[0]
827            cb.SetValue(' key press')
828            OnPlotKeyPress(event)
829                       
830    def OnMotion(event):
831        xpos = event.xdata
832        if xpos:                                        #avoid out of frame mouse position
833            ypos = event.ydata
834            Page.canvas.SetCursor(wx.CROSS_CURSOR)
835            try:
836                if self.Contour:
837                    self.G2plotNB.status.SetStatusText('R =%.3fA pattern ID =%5d'%(xpos,int(ypos)),1)
838                else:
839                    self.G2plotNB.status.SetStatusText('R =%.3fA %s =%.2f'%(xpos,type,ypos),1)                   
840            except TypeError:
841                self.G2plotNB.status.SetStatusText('Select '+type+' pattern first',1)
842   
843    xylim = []
844    try:
845        plotNum = self.G2plotNB.plotList.index(type)
846        Page = self.G2plotNB.nb.GetPage(plotNum)
847        if not newPlot:
848            Plot = Page.figure.gca()          #get previous plot & get limits
849            xylim = Plot.get_xlim(),Plot.get_ylim()
850        Page.figure.clf()
851        Plot = Page.figure.gca()
852    except ValueError:
853        newPlot = True
854        self.Cmax = 1.0
855        Plot = self.G2plotNB.addMpl(type).gca()
856        plotNum = self.G2plotNB.plotList.index(type)
857        Page = self.G2plotNB.nb.GetPage(plotNum)
858        Page.canvas.mpl_connect('key_press_event', OnPlotKeyPress)
859        Page.canvas.mpl_connect('motion_notify_event', OnMotion)
860   
861    Page.SetFocus()
862    self.G2plotNB.status.DestroyChildren()
863    if self.Contour:
864        Choice = (' key press','d: lower contour max','u: raise contour max',
865            'i: interpolation method','s: color scheme','c: contour off')
866    else:
867        Choice = (' key press','l: offset left','r: offset right','d: offset down','u: offset up',
868            'o: reset offset','t: toggle legend','c: contour on','s: toggle single plot')
869    cb = wx.ComboBox(self.G2plotNB.status,style=wx.CB_DROPDOWN|wx.CB_READONLY,
870        choices=Choice)
871    cb.Bind(wx.EVT_COMBOBOX, OnKeyBox)
872    cb.SetValue(' key press')
873    PatternId = self.PatternId
874    PickId = self.PickId
875    Plot.set_title(type)
876    if type == 'G(R)':
877        Plot.set_xlabel(r'$R,\AA$',fontsize=14)
878    else:
879        Plot.set_xlabel(r'$Q,\AA$'+superMinusOne,fontsize=14)
880    Plot.set_ylabel(r''+type,fontsize=14)
881    colors=['b','g','r','c','m','k']
882    name = self.PatternTree.GetItemText(PatternId)[4:]
883    Pattern = []   
884    if self.SinglePlot:
885        name = self.PatternTree.GetItemText(PatternId)
886        name = type+name[4:]
887        Id = G2gd.GetPatternTreeItemId(self,PatternId,name)
888        Pattern = self.PatternTree.GetItemPyData(Id)
889        if Pattern:
890            Pattern.append(name)
891        PlotList = [Pattern,]
892    else:
893        PlotList = []
894        item, cookie = self.PatternTree.GetFirstChild(self.root)
895        while item:
896            if 'PDF' in self.PatternTree.GetItemText(item):
897                name = type+self.PatternTree.GetItemText(item)[4:]
898                Id = G2gd.GetPatternTreeItemId(self,item,name)
899                Pattern = self.PatternTree.GetItemPyData(Id)
900                if Pattern:
901                    Pattern.append(name)
902                    PlotList.append(Pattern)
903            item, cookie = self.PatternTree.GetNextChild(self.root, cookie)
904    PDFdata = self.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(self,PatternId, 'PDF Controls'))
905    numbDen = G2pwd.GetNumDensity(PDFdata['ElList'],PDFdata['Form Vol'])
906    Xb = [0.,10.]
907    Yb = [0.,-40.*np.pi*numbDen]
908    Ymax = 1.0
909    lenX = 0
910    for Pattern in PlotList:
911        xye = Pattern[1]
912        Ymax = max(Ymax,max(xye[1]))
913    offset = self.Offset[0]*Ymax/100.0
914    if self.Contour:
915        ContourZ = []
916        ContourY = []
917        Nseq = 0
918    for N,Pattern in enumerate(PlotList):
919        xye = Pattern[1]
920        if PickId:
921            ifpicked = Pattern[2] == self.PatternTree.GetItemText(PatternId)
922        X = xye[0]
923        if not lenX:
924            lenX = len(X)           
925        Y = xye[1]+offset*N
926        if self.Contour:
927            if lenX == len(X):
928                ContourY.append(N)
929                ContourZ.append(Y)
930                ContourX = X
931                Nseq += 1
932                Plot.set_ylabel('Data sequence',fontsize=12)
933        else:
934            X = xye[0]+self.Offset[1]*.005*N
935            if ifpicked:
936                Plot.plot(X,Y,colors[N%6]+'+',picker=3.,clip_on=False)
937                Page.canvas.SetToolTipString('')
938            else:
939                if self.Legend:
940                    Plot.plot(X,Y,colors[N%6],picker=False,label='Azm:'+Pattern[2].split('=')[1])
941                else:
942                    Plot.plot(X,Y,colors[N%6],picker=False)
943            if type == 'G(R)':
944                Plot.plot(Xb,Yb,color='k',dashes=(5,5))
945            elif type == 'F(Q)':
946                Plot.axhline(0.,color=wx.BLACK)
947            elif type == 'S(Q)':
948                Plot.axhline(1.,color=wx.BLACK)
949    if self.Contour:
950        acolor = mpl.cm.get_cmap(self.ContourColor)
951        Img = Plot.imshow(ContourZ,cmap=acolor,vmin=0,vmax=Ymax*self.Cmax,interpolation=self.Interpolate, 
952            extent=[ContourX[0],ContourX[-1],ContourY[0],ContourY[-1]],aspect='auto',origin='lower')
953        Page.figure.colorbar(Img)
954    elif self.Legend:
955        Plot.legend(loc='best')
956    if not newPlot:
957        Page.toolbar.push_current()
958        Plot.set_xlim(xylim[0])
959        Plot.set_ylim(xylim[1])
960        xylim = []
961        Page.toolbar.push_current()
962        Page.toolbar.draw()
963    else:
964        Page.canvas.draw()
965       
966def PlotXY(self,XY,newPlot=False,type=''):
967    '''simple plot of xy data, used for diagnostic purposes
968    '''
969    def OnMotion(event):
970        xpos = event.xdata
971        if xpos:                                        #avoid out of frame mouse position
972            ypos = event.ydata
973            Page.canvas.SetCursor(wx.CROSS_CURSOR)
974            try:
975                self.G2plotNB.status.SetStatusText('X =%9.3f %s =%9.3f'%(xpos,type,ypos),1)                   
976            except TypeError:
977                self.G2plotNB.status.SetStatusText('Select '+type+' pattern first',1)
978
979    try:
980        plotNum = self.G2plotNB.plotList.index(type)
981        Page = self.G2plotNB.nb.GetPage(plotNum)
982        if not newPlot:
983            Plot = Page.figure.gca()
984            xylim = Plot.get_xlim(),Plot.get_ylim()
985        Page.figure.clf()
986        Plot = Page.figure.gca()
987    except ValueError:
988        newPlot = True
989        Plot = self.G2plotNB.addMpl(type).gca()
990        plotNum = self.G2plotNB.plotList.index(type)
991        Page = self.G2plotNB.nb.GetPage(plotNum)
992        Page.canvas.mpl_connect('motion_notify_event', OnMotion)
993   
994    Page.SetFocus()
995    self.G2plotNB.status.DestroyChildren()
996    Plot.set_title(type)
997    Plot.set_xlabel(r'X',fontsize=14)
998    Plot.set_ylabel(r''+type,fontsize=14)
999    colors=['b','g','r','c','m','k']
1000    Ymax = 1.0
1001    lenX = 0
1002    X,Y = XY[:2]
1003    Ymax = max(Ymax,max(Y))
1004    Plot.plot(X,Y,'k',picker=False)
1005    if not newPlot:
1006        Page.toolbar.push_current()
1007        Plot.set_xlim(xylim[0])
1008        Plot.set_ylim(xylim[1])
1009        xylim = []
1010        Page.toolbar.push_current()
1011        Page.toolbar.draw()
1012    else:
1013        Page.canvas.draw()
1014
1015def PlotPowderLines(self):
1016    ''' plotting of powder lines (i.e. no powder pattern) as sticks
1017    '''
1018    global HKL
1019
1020    def OnMotion(event):
1021        xpos = event.xdata
1022        if xpos:                                        #avoid out of frame mouse position
1023            Page.canvas.SetCursor(wx.CROSS_CURSOR)
1024            self.G2plotNB.status.SetFields(['','2-theta =%9.3f '%(xpos,)])
1025            if self.PickId and self.PatternTree.GetItemText(self.PickId) in ['Index Peak List','Unit Cells List']:
1026                found = []
1027                if len(HKL):
1028                    view = Page.toolbar._views.forward()[0][:2]
1029                    wid = view[1]-view[0]
1030                    found = HKL[np.where(np.fabs(HKL.T[5]-xpos) < 0.002*wid)]
1031                if len(found):
1032                    h,k,l = found[0][:3] 
1033                    Page.canvas.SetToolTipString('%d,%d,%d'%(int(h),int(k),int(l)))
1034                else:
1035                    Page.canvas.SetToolTipString('')
1036
1037    try:
1038        plotNum = self.G2plotNB.plotList.index('Powder Lines')
1039        Page = self.G2plotNB.nb.GetPage(plotNum)
1040        Page.figure.clf()
1041        Plot = Page.figure.gca()
1042    except ValueError:
1043        Plot = self.G2plotNB.addMpl('Powder Lines').gca()
1044        plotNum = self.G2plotNB.plotList.index('Powder Lines')
1045        Page = self.G2plotNB.nb.GetPage(plotNum)
1046        Page.canvas.mpl_connect('motion_notify_event', OnMotion)
1047       
1048    Page.SetFocus()
1049    Plot.set_title('Powder Pattern Lines')
1050    Plot.set_xlabel(r'$\mathsf{2\theta}$',fontsize=14)
1051    PickId = self.PickId
1052    PatternId = self.PatternId
1053    peaks = self.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(self,PatternId, 'Index Peak List'))
1054    for peak in peaks:
1055        Plot.axvline(peak[0],color='b')
1056    HKL = np.array(self.HKL)
1057    for hkl in self.HKL:
1058        Plot.axvline(hkl[5],color='r',dashes=(5,5))
1059    xmin = peaks[0][0]
1060    xmax = peaks[-1][0]
1061    delt = xmax-xmin
1062    xlim = [max(0,xmin-delt/20.),min(180.,xmax+delt/20.)]
1063    Plot.set_xlim(xlim)
1064    Page.canvas.draw()
1065    Page.toolbar.push_current()
1066
1067def PlotPeakWidths(self):
1068    ''' Plotting of instrument broadening terms as function of 2-theta (future TOF)
1069    Seen when "Instrument Parameters" chosen from powder pattern data tree
1070    '''
1071    PatternId = self.PatternId
1072    limitID = G2gd.GetPatternTreeItemId(self,PatternId, 'Limits')
1073    if limitID:
1074        limits = self.PatternTree.GetItemPyData(limitID)
1075    else:
1076        return
1077    instParms = self.PatternTree.GetItemPyData( \
1078        G2gd.GetPatternTreeItemId(self,PatternId, 'Instrument Parameters'))
1079    if instParms[0][0] in ['PXC','PNC']:
1080        lam = instParms[1][1]
1081        if len(instParms[1]) == 13:
1082            GU,GV,GW,LX,LY = instParms[0][6:11]
1083        else:
1084            GU,GV,GW,LX,LY = instParms[0][4:9]
1085    peakID = G2gd.GetPatternTreeItemId(self,PatternId, 'Peak List')
1086    if peakID:
1087        peaks = self.PatternTree.GetItemPyData(peakID)
1088    else:
1089        peaks = []
1090   
1091    try:
1092        plotNum = self.G2plotNB.plotList.index('Peak Widths')
1093        Page = self.G2plotNB.nb.GetPage(plotNum)
1094        Page.figure.clf()
1095        Plot = Page.figure.gca()
1096    except ValueError:
1097        Plot = self.G2plotNB.addMpl('Peak Widths').gca()
1098        plotNum = self.G2plotNB.plotList.index('Peak Widths')
1099        Page = self.G2plotNB.nb.GetPage(plotNum)
1100    Page.SetFocus()
1101   
1102    Page.canvas.SetToolTipString('')
1103    colors=['b','g','r','c','m','k']
1104    Xmin,Xmax = limits[1]
1105    Xmin = min(0.5,max(Xmin,1))
1106    Xmin /= 2
1107    Xmax /= 2
1108    nPts = 100
1109    delt = (Xmax-Xmin)/nPts
1110    thetas = []
1111    for i in range(nPts):
1112        thetas.append(Xmin+i*delt)
1113    X = []
1114    Y = []
1115    Z = []
1116    W = []
1117    sig = lambda Th,U,V,W: 1.17741*math.sqrt(U*tand(Th)**2+V*tand(Th)+W)*math.pi/18000.
1118    gam = lambda Th,X,Y: (X/cosd(Th)+Y*tand(Th))*math.pi/18000.
1119    gamFW = lambda s,g: math.exp(math.log(s**5+2.69269*s**4*g+2.42843*s**3*g**2+4.47163*s**2*g**3+0.07842*s*g**4+g**5)/5.)
1120#    gamFW2 = lambda s,g: math.sqrt(s**2+(0.4654996*g)**2)+.5345004*g  #Ubaldo Bafile - private communication
1121    for theta in thetas:
1122        X.append(4.0*math.pi*sind(theta)/lam)              #q
1123        s = sig(theta,GU,GV,GW)
1124        g = gam(theta,LX,LY)
1125        G = gamFW(g,s)
1126        Y.append(s/tand(theta))
1127        Z.append(g/tand(theta))
1128        W.append(G/tand(theta))
1129    Plot.set_title('Instrument and sample peak widths')
1130    Plot.set_ylabel(r'$\Delta q/q, \Delta d/d$',fontsize=14)
1131    Plot.set_xlabel(r'$q, \AA^{-1}$',fontsize=14)
1132    Plot.plot(X,Y,color='r',label='Gaussian')
1133    Plot.plot(X,Z,color='g',label='Lorentzian')
1134    Plot.plot(X,W,color='b',label='G+L')
1135    X = []
1136    Y = []
1137    Z = []
1138    W = []
1139    V = []
1140    for peak in peaks:
1141        X.append(4.0*math.pi*sind(peak[0]/2.0)/lam)
1142        try:
1143            s = 1.17741*math.sqrt(peak[4])*math.pi/18000.
1144        except ValueError:
1145            s = 0.01
1146        g = peak[6]*math.pi/18000.
1147        G = gamFW(g,s)
1148        Y.append(s/tand(peak[0]/2.))
1149        Z.append(g/tand(peak[0]/2.))
1150        W.append(G/tand(peak[0]/2.))
1151    Plot.plot(X,Y,'+',color='r',label='G peak')
1152    Plot.plot(X,Z,'+',color='g',label='L peak')
1153    Plot.plot(X,W,'+',color='b',label='G+L peak')
1154    Plot.legend(loc='best')
1155    Page.canvas.draw()
1156   
1157def PlotSizeStrainPO(self,data,Start=False):
1158    '''Plot 3D mustrain/size/preferred orientation figure. In this instance data is for a phase
1159    '''
1160   
1161    PatternId = self.PatternId
1162    generalData = data['General']
1163    SGData = generalData['SGData']
1164    SGLaue = SGData['SGLaue']
1165    if Start:                   #initialize the spherical harmonics qlmn arrays
1166        ptx.pyqlmninit()
1167        Start = False
1168    MuStrKeys = G2spc.MustrainNames(SGData)
1169    cell = generalData['Cell'][1:]
1170    A,B = G2lat.cell2AB(cell[:6])
1171    Vol = cell[6]
1172    useList = data['Histograms']
1173    phase = generalData['Name']
1174    plotType = generalData['Data plot type']
1175    plotDict = {'Mustrain':'Mustrain','Size':'Size','Preferred orientation':'Pref.Ori.'}
1176    for ptype in plotDict:
1177        self.G2plotNB.Delete(ptype)       
1178
1179    for item in useList:
1180        if useList[item]['Show']:
1181            break
1182    else:
1183        return            #nothing to show!!
1184   
1185    numPlots = len(useList)
1186
1187    if plotType in ['Mustrain','Size']:
1188        Plot = mp3d.Axes3D(self.G2plotNB.add3D(plotType))
1189    else:
1190        Plot = self.G2plotNB.addMpl(plotType).gca()       
1191    plotNum = self.G2plotNB.plotList.index(plotType)
1192    Page = self.G2plotNB.nb.GetPage(plotNum)
1193    Page.SetFocus()
1194    self.G2plotNB.status.SetStatusText('Adjust frame size to get desired aspect ratio',1)
1195    if not Page.IsShown():
1196        Page.Show()
1197   
1198    for item in useList:
1199        if useList[item]['Show']:
1200            PHI = np.linspace(0.,360.,30,True)
1201            PSI = np.linspace(0.,180.,30,True)
1202            X = np.outer(npsind(PHI),npsind(PSI))
1203            Y = np.outer(npcosd(PHI),npsind(PSI))
1204            Z = np.outer(np.ones(np.size(PHI)),npcosd(PSI))
1205            try:        #temp patch instead of 'mustrain' for old files with 'microstrain'
1206                coeff = useList[item][plotDict[plotType]]
1207            except KeyError:
1208                break
1209            if plotType in ['Mustrain','Size']:
1210                if coeff[0] == 'isotropic':
1211                    X *= coeff[1][0]
1212                    Y *= coeff[1][0]
1213                    Z *= coeff[1][0]                               
1214                elif coeff[0] == 'uniaxial':
1215                   
1216                    def uniaxCalc(xyz,iso,aniso,axes):
1217                        Z = np.array(axes)
1218                        cp = abs(np.dot(xyz,Z))
1219                        sp = np.sqrt(1.-cp**2)
1220                        R = iso*aniso/np.sqrt((iso*cp)**2+(aniso*sp)**2)
1221                        return R*xyz
1222                       
1223                    iso,aniso = coeff[1][:2]
1224                    axes = np.inner(A,np.array(coeff[3]))
1225                    axes /= nl.norm(axes)
1226                    Shkl = np.array(coeff[1])
1227                    XYZ = np.dstack((X,Y,Z))
1228                    XYZ = np.nan_to_num(np.apply_along_axis(uniaxCalc,2,XYZ,iso,aniso,axes))
1229                    X,Y,Z = np.dsplit(XYZ,3)
1230                    X = X[:,:,0]
1231                    Y = Y[:,:,0]
1232                    Z = Z[:,:,0]
1233               
1234                elif coeff[0] == 'ellipsoidal':
1235                   
1236                    def ellipseCalc(xyz,E,R):
1237                        XYZ = xyz*E.T
1238                        return np.inner(XYZ.T,R)
1239                       
1240                    S6 = coeff[4]
1241                    Sij = G2lat.U6toUij(S6)
1242                    E,R = nl.eigh(Sij)
1243                    XYZ = np.dstack((X,Y,Z))
1244                    XYZ = np.nan_to_num(np.apply_along_axis(ellipseCalc,2,XYZ,E,R))
1245                    X,Y,Z = np.dsplit(XYZ,3)
1246                    X = X[:,:,0]
1247                    Y = Y[:,:,0]
1248                    Z = Z[:,:,0]
1249                   
1250                elif coeff[0] == 'generalized':
1251                   
1252                    def genMustrain(xyz,SGData,A,Shkl):
1253                        uvw = np.inner(A.T,xyz)
1254                        Strm = np.array(G2spc.MustrainCoeff(uvw,SGData))
1255                        sum = np.sum(np.multiply(Shkl,Strm))
1256                        sum = np.where(sum > 0.01,sum,0.01)
1257                        sum = np.sqrt(sum)*math.pi/0.018      #centidegrees to radians!
1258                        return sum*xyz
1259                       
1260                    Shkl = np.array(coeff[4])
1261                    if np.any(Shkl):
1262                        XYZ = np.dstack((X,Y,Z))
1263                        XYZ = np.nan_to_num(np.apply_along_axis(genMustrain,2,XYZ,SGData,A,Shkl))
1264                        X,Y,Z = np.dsplit(XYZ,3)
1265                        X = X[:,:,0]
1266                        Y = Y[:,:,0]
1267                        Z = Z[:,:,0]
1268                           
1269                if np.any(X) and np.any(Y) and np.any(Z):
1270                    errFlags = np.seterr(all='ignore')
1271                    Plot.plot_surface(X,Y,Z,rstride=1,cstride=1,color='g',linewidth=1)
1272                    np.seterr(all='ignore')
1273                    xyzlim = np.array([Plot.get_xlim3d(),Plot.get_ylim3d(),Plot.get_zlim3d()]).T
1274                    XYZlim = [min(xyzlim[0]),max(xyzlim[1])]
1275                    Plot.set_xlim3d(XYZlim)
1276                    Plot.set_ylim3d(XYZlim)
1277                    Plot.set_zlim3d(XYZlim)
1278                    Plot.set_aspect('equal')
1279                if plotType == 'Size':
1280                    Plot.set_title('Crystallite size for '+phase+'\n'+coeff[0]+' model')
1281                    Plot.set_xlabel(r'X, $\mu$m')
1282                    Plot.set_ylabel(r'Y, $\mu$m')
1283                    Plot.set_zlabel(r'Z, $\mu$m')
1284                else:   
1285                    Plot.set_title(r'$\mu$strain for '+phase+'\n'+coeff[0]+' model')
1286                    Plot.set_xlabel(r'X, $\mu$strain')
1287                    Plot.set_ylabel(r'Y, $\mu$strain')
1288                    Plot.set_zlabel(r'Z, $\mu$strain')
1289            else:
1290                h,k,l = generalData['POhkl']
1291                if coeff[0] == 'MD':
1292                    print 'March-Dollase preferred orientation plot'
1293               
1294                else:
1295                    PH = np.array(generalData['POhkl'])
1296                    phi,beta = G2lat.CrsAng(PH,cell[:6],SGData)
1297                    SHCoef = {}
1298                    for item in coeff[5]:
1299                        L,N = eval(item.strip('C'))
1300                        SHCoef['C%d,0,%d'%(L,N)] = coeff[5][item]                       
1301                    ODFln = G2lat.Flnh(Start,SHCoef,phi,beta,SGData)
1302                    X = np.linspace(0,90.0,26)
1303                    Y = G2lat.polfcal(ODFln,'0',X,0.0)
1304                    Plot.plot(X,Y,color='k',label=str(PH))
1305                    Plot.legend(loc='best')
1306                    Plot.set_title('Axial distribution for HKL='+str(PH)+' in '+phase)
1307                    Plot.set_xlabel(r'$\psi$',fontsize=16)
1308                    Plot.set_ylabel('MRD',fontsize=14)
1309    Page.canvas.draw()
1310   
1311def PlotTexture(self,data,newPlot=False,Start=False):
1312    '''Pole figure, inverse pole figure, 3D pole distribution and 3D inverse pole distribution
1313    plotting.
1314    dict generalData contains all phase info needed which is in data
1315    '''
1316
1317    shModels = ['cylindrical','none','shear - 2/m','rolling - mmm']
1318    SamSym = dict(zip(shModels,['0','-1','2/m','mmm']))
1319    PatternId = self.PatternId
1320    generalData = data['General']
1321    SGData = generalData['SGData']
1322    textureData = generalData['SH Texture']
1323    self.G2plotNB.Delete('Texture')
1324    if not textureData['Order']:
1325        return                  #no plot!!
1326    SHData = generalData['SH Texture']
1327    SHCoef = SHData['SH Coeff'][1]
1328    cell = generalData['Cell'][1:7]
1329    Amat,Bmat = G2lat.cell2AB(cell)
1330    sq2 = 1.0/math.sqrt(2.0)
1331   
1332    def rp2xyz(r,p):
1333        z = npcosd(r)
1334        xy = np.sqrt(1.-z**2)
1335        return xy*npsind(p),xy*npcosd(p),z
1336           
1337    def OnMotion(event):
1338        SHData = data['General']['SH Texture']
1339        if event.xdata and event.ydata:                 #avoid out of frame errors
1340            xpos = event.xdata
1341            ypos = event.ydata
1342            if 'Inverse' in SHData['PlotType']:
1343                r = xpos**2+ypos**2
1344                if r <= 1.0:
1345                    if 'equal' in self.Projection: 
1346                        r,p = 2.*npasind(np.sqrt(r)*sq2),npatan2d(ypos,xpos)
1347                    else:
1348                        r,p = 2.*npatand(np.sqrt(r)),npatan2d(ypos,xpos)
1349                    ipf = G2lat.invpolfcal(IODFln,SGData,np.array([r,]),np.array([p,]))
1350                    xyz = np.inner(Amat.T,np.array([rp2xyz(r,p)]))
1351                    y,x,z = list(xyz/np.max(np.abs(xyz)))
1352                   
1353                    self.G2plotNB.status.SetFields(['',
1354                        'psi =%9.3f, beta =%9.3f, MRD =%9.3f hkl=%5.2f,%5.2f,%5.2f'%(r,p,ipf,x,y,z)])
1355                                   
1356            elif 'Axial' in SHData['PlotType']:
1357                pass
1358               
1359            else:                       #ordinary pole figure
1360                z = xpos**2+ypos**2
1361                if z <= 1.0:
1362                    z = np.sqrt(z)
1363                    if 'equal' in self.Projection: 
1364                        r,p = 2.*npasind(z*sq2),npatan2d(ypos,xpos)
1365                    else:
1366                        r,p = 2.*npatand(z),npatan2d(ypos,xpos)
1367                    pf = G2lat.polfcal(ODFln,SamSym[textureData['Model']],np.array([r,]),np.array([p,]))
1368                    self.G2plotNB.status.SetFields(['','phi =%9.3f, gam =%9.3f, MRD =%9.3f'%(r,p,pf)])
1369   
1370    if self.Projection == '3D display' and 'Axial'  not in SHData['PlotType']:               
1371        Plot = mp3d.Axes3D(self.G2plotNB.add3D('Texture'))
1372    else:
1373        Plot = self.G2plotNB.addMpl('Texture').gca()
1374    plotNum = self.G2plotNB.plotList.index('Texture')
1375    Page = self.G2plotNB.nb.GetPage(plotNum)
1376    if self.Projection == '3D display':
1377        pass
1378    else:
1379        Page.canvas.mpl_connect('motion_notify_event', OnMotion)
1380
1381    Page.SetFocus()
1382    self.G2plotNB.status.SetFields(['',''])   
1383    PH = np.array(SHData['PFhkl'])
1384    phi,beta = G2lat.CrsAng(PH,cell,SGData)
1385    ODFln = G2lat.Flnh(Start,SHCoef,phi,beta,SGData)
1386    PX = np.array(SHData['PFxyz'])
1387    gam = atan2d(PX[0],PX[1])
1388    xy = math.sqrt(PX[0]**2+PX[1]**2)
1389    xyz = math.sqrt(PX[0]**2+PX[1]**2+PX[2]**2)
1390    psi = asind(xy/xyz)
1391    IODFln = G2lat.Glnh(Start,SHCoef,psi,gam,SamSym[textureData['Model']])
1392    if 'Axial' in SHData['PlotType']:
1393        X = np.linspace(0,90.0,26)
1394        Y = G2lat.polfcal(ODFln,SamSym[textureData['Model']],X,0.0)
1395        Plot.plot(X,Y,color='k',label=str(SHData['PFhkl']))
1396        Plot.legend(loc='best')
1397        Plot.set_title('Axial distribution for HKL='+str(SHData['PFhkl']))
1398        Plot.set_xlabel(r'$\psi$',fontsize=16)
1399        Plot.set_ylabel('MRD',fontsize=14)
1400       
1401    else:       
1402        npts = 201
1403        if 'Inverse' in SHData['PlotType']:
1404            X,Y = np.meshgrid(np.linspace(1.,-1.,npts),np.linspace(-1.,1.,npts))
1405            R,P = np.sqrt(X**2+Y**2).flatten(),npatan2d(X,Y).flatten()
1406            if 'equal' in self.Projection:
1407                R = np.where(R <= 1.,2.*npasind(R*sq2),0.0)
1408            else:
1409                R = np.where(R <= 1.,2.*npatand(R),0.0)
1410            Z = np.zeros_like(R)
1411            Z = G2lat.invpolfcal(IODFln,SGData,R,P)
1412            Z = np.reshape(Z,(npts,npts))
1413            CS = Plot.contour(Y,X,Z,aspect='equal')
1414            Plot.clabel(CS,fontsize=9,inline=1)
1415            try:
1416                Img = Plot.imshow(Z.T,aspect='equal',cmap=self.ContourColor,extent=[-1,1,-1,1])
1417            except ValueError:
1418                pass
1419            if newPlot:
1420#                Page.figure.colorbar(Img)    #colorbar fails - crashes gsasii
1421                newPlot = False
1422            Plot.set_title('Inverse pole figure for XYZ='+str(SHData['PFxyz']))
1423            Plot.set_xlabel(self.Projection.capitalize()+' projection')
1424                       
1425        else:
1426            X,Y = np.meshgrid(np.linspace(1.,-1.,npts),np.linspace(-1.,1.,npts))
1427            R,P = np.sqrt(X**2+Y**2).flatten(),npatan2d(X,Y).flatten()
1428            if 'equal' in self.Projection:
1429                R = np.where(R <= 1.,2.*npasind(R*sq2),0.0)
1430            else:
1431                R = np.where(R <= 1.,2.*npatand(R),0.0)
1432            Z = np.zeros_like(R)
1433            Z = G2lat.polfcal(ODFln,SamSym[textureData['Model']],R,P)
1434            Z = np.reshape(Z,(npts,npts))
1435            try:
1436                CS = Plot.contour(Y,X,Z,aspect='equal')
1437                Plot.clabel(CS,fontsize=9,inline=1)
1438            except ValueError:
1439                pass
1440            Img = Plot.imshow(Z.T,aspect='equal',cmap=self.ContourColor,extent=[-1,1,-1,1])
1441            if newPlot:
1442#                Page.figure.colorbar(Img)    #colorbar fails - crashes gsasii
1443                newPlot = False
1444            Plot.set_title('Pole figure for HKL='+str(SHData['PFhkl']))
1445            Plot.set_xlabel(self.Projection.capitalize()+' projection')
1446    Page.canvas.draw()
1447
1448def PlotCovariance(self,Data={}):
1449    if not Data:
1450        Data = self.PatternTree.GetItemPyData(
1451            G2gd.GetPatternTreeItemId(self,self.root, 'Covariance'))
1452    if not Data:
1453        print 'No covariance matrix available'
1454        return
1455    varyList = Data['varyList']
1456    values = Data['variables']
1457    Xmax = len(varyList)
1458    covMatrix = Data['covMatrix']
1459    sig = np.sqrt(np.diag(covMatrix))
1460    xvar = np.outer(sig,np.ones_like(sig))
1461    covArray = np.divide(np.divide(covMatrix,xvar),xvar.T)
1462    title = ' for\n'+Data['title']
1463    newAtomDict = Data['newAtomDict']
1464
1465    def OnPlotKeyPress(event):
1466        newPlot = False
1467        if event.key == 's':
1468            choice = [m for m in mpl.cm.datad.keys() if not m.endswith("_r")]
1469            choice.sort()
1470            dlg = wx.SingleChoiceDialog(self,'Select','Color scheme',choice)
1471            if dlg.ShowModal() == wx.ID_OK:
1472                sel = dlg.GetSelection()
1473                self.VcovColor = choice[sel]
1474            else:
1475                self.VcovColor = 'RdYlGn'
1476            dlg.Destroy()
1477        PlotCovariance(self)
1478
1479    def OnMotion(event):
1480        if event.button:
1481            ytics = imgAx.get_yticks()
1482            ytics = np.where(ytics<len(varyList),ytics,-1)
1483            ylabs = [np.where(0<=i ,varyList[int(i)],' ') for i in ytics]
1484            imgAx.set_yticklabels(ylabs)
1485           
1486        if event.xdata and event.ydata:                 #avoid out of frame errors
1487            xpos = int(event.xdata+.5)
1488            ypos = int(event.ydata+.5)
1489            if -1 < xpos < len(varyList) and -1 < ypos < len(varyList):
1490                if xpos == ypos:
1491                    value = values[xpos]
1492                    name = varyList[xpos]
1493                    if varyList[xpos] in newAtomDict:
1494                        name,value = newAtomDict[name]                       
1495                    msg = '%s value = %.4g, esd = %.4g'%(name,value,sig[xpos])
1496                else:
1497                    msg = '%s - %s: %5.3f'%(varyList[xpos],varyList[ypos],covArray[xpos][ypos])
1498                Page.canvas.SetToolTipString(msg)
1499                self.G2plotNB.status.SetFields(['Key: s to change colors',msg])
1500    try:
1501        plotNum = self.G2plotNB.plotList.index('Covariance')
1502        Page = self.G2plotNB.nb.GetPage(plotNum)
1503        Page.figure.clf()
1504        Plot = Page.figure.gca()
1505        if not Page.IsShown():
1506            Page.Show()
1507    except ValueError:
1508        Plot = self.G2plotNB.addMpl('Covariance').gca()
1509        plotNum = self.G2plotNB.plotList.index('Covariance')
1510        Page = self.G2plotNB.nb.GetPage(plotNum)
1511        Page.canvas.mpl_connect('motion_notify_event', OnMotion)
1512        Page.canvas.mpl_connect('key_press_event', OnPlotKeyPress)
1513
1514    Page.SetFocus()
1515    self.G2plotNB.status.SetFields(['',''])   
1516    acolor = mpl.cm.get_cmap(self.VcovColor)
1517    Img = Plot.imshow(covArray,aspect='equal',cmap=acolor,interpolation='nearest',origin='lower')
1518    imgAx = Img.get_axes()
1519    ytics = imgAx.get_yticks()
1520    ylabs = [varyList[int(i)] for i in ytics[:-1]]
1521    imgAx.set_yticklabels(ylabs)
1522    colorBar = Page.figure.colorbar(Img)
1523    Plot.set_title('V-Cov matrix'+title)
1524    Plot.set_xlabel('Variable number')
1525    Plot.set_ylabel('Variable name')
1526    Page.canvas.draw()
1527   
1528def PlotSeq(self,SeqData,SeqSig,SeqNames,sampleParm):
1529   
1530    def OnKeyPress(event):
1531        if event.key == 's' and sampleParm:
1532            if self.xAxis:
1533                self.xAxis = False
1534            else:
1535                self.xAxis = True
1536            Draw(False)
1537    try:
1538        plotNum = self.G2plotNB.plotList.index('Sequential refinement')
1539        Page = self.G2plotNB.nb.GetPage(plotNum)
1540        Page.figure.clf()
1541        Plot = Page.figure.gca()
1542        if not Page.IsShown():
1543            Page.Show()
1544    except ValueError:
1545        Plot = self.G2plotNB.addMpl('Sequential refinement').gca()
1546        plotNum = self.G2plotNB.plotList.index('Sequential refinement')
1547        Page = self.G2plotNB.nb.GetPage(plotNum)
1548        Page.canvas.mpl_connect('key_press_event', OnKeyPress)
1549        self.xAxis = False
1550       
1551    def Draw(newPlot):
1552        Page.SetFocus()
1553        self.G2plotNB.status.SetFields(['','press s to toggle x-axis = sample environment parameter'])
1554        if len(SeqData):
1555            Plot.clear()
1556            if self.xAxis:   
1557                xName = sampleParm.keys()[0]
1558                X = sampleParm[xName]
1559            else:
1560                X = np.arange(0,len(SeqData[0]),1)
1561                xName = 'Data sequence number'
1562            for Y,sig,name in zip(SeqData,SeqSig,SeqNames):
1563                Plot.errorbar(X,Y,yerr=sig,label=name)       
1564            Plot.legend(loc='best')
1565            Plot.set_ylabel('Parameter values')
1566            Plot.set_xlabel(xName)
1567            Page.canvas.draw()           
1568    Draw(True)
1569           
1570def PlotExposedImage(self,newPlot=False,event=None):
1571    '''General access module for 2D image plotting
1572    '''
1573    plotNo = self.G2plotNB.nb.GetSelection()
1574    if self.G2plotNB.nb.GetPageText(plotNo) == '2D Powder Image':
1575        PlotImage(self,newPlot,event,newImage=True)
1576    elif self.G2plotNB.nb.GetPageText(plotNo) == '2D Integration':
1577        PlotIntegration(self,newPlot,event)
1578
1579def PlotImage(self,newPlot=False,event=None,newImage=True):
1580    '''Plot of 2D detector images as contoured plot. Also plot calibration ellipses,
1581    masks, etc.
1582    '''
1583    from matplotlib.patches import Ellipse,Arc,Circle,Polygon
1584    import numpy.ma as ma
1585    Dsp = lambda tth,wave: wave/(2.*sind(tth/2.))
1586    global Data,Masks
1587    colors=['b','g','r','c','m','k']
1588    Data = self.PatternTree.GetItemPyData(
1589        G2gd.GetPatternTreeItemId(self,self.Image, 'Image Controls'))
1590    Masks = self.PatternTree.GetItemPyData(
1591        G2gd.GetPatternTreeItemId(self,self.Image, 'Masks'))
1592
1593    def OnImMotion(event):
1594        Page.canvas.SetToolTipString('')
1595        sizexy = Data['size']
1596        if event.xdata and event.ydata:                 #avoid out of frame errors
1597            Page.canvas.SetCursor(wx.CROSS_CURSOR)
1598            item = self.itemPicked
1599            pixelSize = Data['pixelSize']
1600            scalex = 1000./pixelSize[0]
1601            scaley = 1000./pixelSize[1]
1602            if item and self.PatternTree.GetItemText(self.PickId) == 'Image Controls':
1603                if 'Text' in str(item):
1604                    Page.canvas.SetToolTipString('%8.3f %8.3fmm'%(event.xdata,event.ydata))
1605                else:
1606                    xcent,ycent = Data['center']
1607                    xpos = event.xdata-xcent
1608                    ypos = event.ydata-ycent
1609                    tth,azm = G2img.GetTthAzm(event.xdata,event.ydata,Data)
1610                    if 'line3' in  str(item) or 'line4' in str(item) and not Data['fullIntegrate']:
1611                        Page.canvas.SetToolTipString('%6d deg'%(azm))
1612                    elif 'line1' in  str(item) or 'line2' in str(item):
1613                        Page.canvas.SetToolTipString('%8.3fdeg'%(tth))                           
1614            else:
1615                xpos = event.xdata
1616                ypos = event.ydata
1617                xpix = xpos*scalex
1618                ypix = ypos*scaley
1619                Int = 0
1620                if (0 <= xpix <= sizexy[0]) and (0 <= ypix <= sizexy[1]):
1621                    Int = self.ImageZ[ypix][xpix]
1622                tth,azm,dsp = G2img.GetTthAzmDsp(xpos,ypos,Data)
1623                Q = 2.*math.pi/dsp
1624                if self.setPoly:
1625                    self.G2plotNB.status.SetFields(['','Polygon mask pick - LB next point, RB close polygon'])
1626                else:
1627                    self.G2plotNB.status.SetFields(\
1628                        ['','Detector 2-th =%9.3fdeg, dsp =%9.3fA, Q = %6.5fA-1, azm = %7.2fdeg, I = %6d'%(tth,dsp,Q,azm,Int)])
1629
1630    def OnImPlotKeyPress(event):
1631        try:
1632            PickName = self.PatternTree.GetItemText(self.PickId)
1633        except TypeError:
1634            return
1635        if PickName == 'Masks':
1636            Xpos = event.xdata
1637            if not Xpos:            #got point out of frame
1638                return
1639            Ypos = event.ydata
1640            if event.key == 's':
1641                Masks['Points'].append([Xpos,Ypos,1.])
1642            elif event.key == 'r':
1643                tth = G2img.GetTth(Xpos,Ypos,Data)
1644                Masks['Rings'].append([tth,0.1])
1645            elif event.key == 'a':
1646                tth,azm = G2img.GetTthAzm(Xpos,Ypos,Data)
1647                azm = int(azm)               
1648                Masks['Arcs'].append([tth,[azm-5,azm+5],0.1])
1649            elif event.key == 'p':
1650                self.setPoly = True
1651                Masks['Polygons'].append([])
1652                self.G2plotNB.status.SetFields(['','Polygon mask active - LB pick next point, RB close polygon'])
1653            G2imG.UpdateMasks(self,Masks)
1654        elif PickName == 'Image Controls':
1655            if event.key == 'c':
1656                Xpos = event.xdata
1657                if not Xpos:            #got point out of frame
1658                    return
1659                Ypos = event.ydata
1660                dlg = wx.MessageDialog(self,'Are you sure you want to change the center?',
1661                    'Center change',style=wx.OK|wx.CANCEL)
1662                try:
1663                    if dlg.ShowModal() == wx.ID_OK:
1664                        print 'move center to: ',Xpos,Ypos
1665                        Data['center'] = [Xpos,Ypos]
1666                        G2imG.UpdateImageControls(self,Data,Masks)
1667                finally:
1668                    dlg.Destroy()
1669            elif event.key == 'l':
1670                if self.logPlot:
1671                    self.logPlot = False
1672                else:
1673                    self.logPlot = True
1674        PlotImage(self,newImage=True)
1675           
1676    def OnKeyBox(event):
1677        if self.G2plotNB.nb.GetSelection() == self.G2plotNB.plotList.index('2D Powder Image'):
1678            event.key = cb.GetValue()[0]
1679            cb.SetValue(' key press')
1680            if event.key in 'l':
1681                OnImPlotKeyPress(event)
1682                       
1683    def OnImPick(event):
1684        if self.PatternTree.GetItemText(self.PickId) not in ['Image Controls','Masks']:
1685            return
1686        if self.setPoly:
1687            polygon = Masks['Polygons'][-1]
1688            xpos,ypos = event.mouseevent.xdata,event.mouseevent.ydata
1689            if xpos and ypos:                       #point inside image
1690                if len(polygon) > 2 and event.mouseevent.button == 3:
1691                    x0,y0 = polygon[0]
1692                    polygon.append([x0,y0])
1693                    self.setPoly = False
1694                    self.G2plotNB.status.SetFields(['','Polygon closed - RB drag a vertex to change shape'])
1695                else:
1696                    self.G2plotNB.status.SetFields(['','New polygon point: %.1f,%.1f'%(xpos,ypos)])
1697                    polygon.append([xpos,ypos])
1698                G2imG.UpdateMasks(self,Masks)
1699        else:
1700            if self.itemPicked is not None: return
1701            self.itemPicked = event.artist
1702            self.mousePicked = event.mouseevent
1703       
1704    def OnImRelease(event):
1705        try:
1706            PickName = self.PatternTree.GetItemText(self.PickId)
1707        except TypeError:
1708            return
1709        if PickName not in ['Image Controls','Masks']:
1710            return
1711        pixelSize = Data['pixelSize']
1712        scalex = 1000./pixelSize[0]
1713        scaley = 1000./pixelSize[1]
1714        pixLimit = Data['pixLimit']
1715        if self.itemPicked is None and PickName == 'Image Controls':
1716#            sizexy = Data['size']
1717            Xpos = event.xdata
1718            if not (Xpos and self.ifGetRing):                   #got point out of frame
1719                return
1720            Ypos = event.ydata
1721            if Ypos and not Page.toolbar._active:         #make sure zoom/pan not selected
1722                if event.button == 1:
1723                    Xpix = Xpos*scalex
1724                    Ypix = Ypos*scaley
1725                    xpos,ypos,I,J = G2img.ImageLocalMax(self.ImageZ,pixLimit,Xpix,Ypix)
1726                    if I and J:
1727                        xpos += .5                              #shift to pixel center
1728                        ypos += .5
1729                        xpos /= scalex                          #convert to mm
1730                        ypos /= scaley
1731                        Data['ring'].append([xpos,ypos])
1732                elif event.button == 3:
1733                    self.dataFrame.GetStatusBar().SetStatusText('Calibrating...')
1734                    if G2img.ImageCalibrate(self,Data):
1735                        self.dataFrame.GetStatusBar().SetStatusText('Calibration successful - Show ring picks to check')
1736                        print 'Calibration successful'
1737                    else:
1738                        self.dataFrame.GetStatusBar().SetStatusText('Calibration failed - Show ring picks to diagnose')
1739                        print 'Calibration failed'
1740                    self.ifGetRing = False
1741                    G2imG.UpdateImageControls(self,Data,Masks)
1742                    return
1743                PlotImage(self,newImage=False)
1744            return
1745        else:
1746            xpos = event.xdata
1747            if xpos:                                        #avoid out of frame mouse position
1748                ypos = event.ydata
1749                if self.ifGetRing:                          #delete a calibration ring pick
1750                    xypos = [xpos,ypos]
1751                    rings = Data['ring']
1752                    for ring in rings:
1753                        if np.allclose(ring,xypos,.01,0):
1754                            rings.remove(ring)                                                                       
1755                else:
1756                    tth,azm,dsp = G2img.GetTthAzmDsp(xpos,ypos,Data)
1757                    itemPicked = str(self.itemPicked)
1758                    if 'Line2D' in itemPicked and PickName == 'Image Controls':
1759                        if 'line1' in itemPicked:
1760                            Data['IOtth'][0] = max(tth,0.001)
1761                        elif 'line2' in itemPicked:
1762                            Data['IOtth'][1] = tth
1763                        elif 'line3' in itemPicked:
1764                            Data['LRazimuth'][0] = int(azm)
1765                            if Data['fullIntegrate']:
1766                                Data['LRazimuth'][1] = Data['LRazimuth'][0]+360
1767                        elif 'line4' in itemPicked and not Data['fullIntegrate']:
1768                            Data['LRazimuth'][1] = int(azm)
1769                           
1770                        if Data['LRazimuth'][0] > Data['LRazimuth'][1]:
1771                            Data['LRazimuth'][0] -= 360
1772                           
1773                        azLim = np.array(Data['LRazimuth'])   
1774                        if np.any(azLim>360):
1775                            azLim -= 360
1776                            Data['LRazimuth'] = list(azLim)
1777                           
1778                        if  Data['IOtth'][0] > Data['IOtth'][1]:
1779                            Data['IOtth'][0],Data['IOtth'][1] = Data['IOtth'][1],Data['IOtth'][0]
1780                           
1781                        self.InnerTth.SetValue("%8.2f" % (Data['IOtth'][0]))
1782                        self.OuterTth.SetValue("%8.2f" % (Data['IOtth'][1]))
1783                        self.Lazim.SetValue("%6d" % (Data['LRazimuth'][0]))
1784                        self.Razim.SetValue("%6d" % (Data['LRazimuth'][1]))
1785                    elif 'Circle' in itemPicked and PickName == 'Masks':
1786                        spots = Masks['Points']
1787                        newPos = itemPicked.split(')')[0].split('(')[2].split(',')
1788                        newPos = np.array([float(newPos[0]),float(newPos[1])])
1789                        for spot in spots:
1790                            if np.allclose(np.array([spot[:2]]),newPos):
1791                                spot[:2] = xpos,ypos
1792                        G2imG.UpdateMasks(self,Masks)
1793                    elif 'Line2D' in itemPicked and PickName == 'Masks':
1794                        Obj = self.itemPicked.findobj()
1795                        rings = Masks['Rings']
1796                        arcs = Masks['Arcs']
1797                        polygons = Masks['Polygons']
1798                        for ring in self.ringList:
1799                            if Obj == ring[0]:
1800                                rN = ring[1]
1801                                if ring[2] == 'o':
1802                                    rings[rN][0] = G2img.GetTth(xpos,ypos,Data)-rings[rN][1]/2.
1803                                else:
1804                                    rings[rN][0] = G2img.GetTth(xpos,ypos,Data)+rings[rN][1]/2.
1805                        for arc in self.arcList:
1806                            if Obj == arc[0]:
1807                                aN = arc[1]
1808                                if arc[2] == 'o':
1809                                    arcs[aN][0] = G2img.GetTth(xpos,ypos,Data)-arcs[aN][2]/2
1810                                elif arc[2] == 'i':
1811                                    arcs[aN][0] = G2img.GetTth(xpos,ypos,Data)+arcs[aN][2]/2
1812                                elif arc[2] == 'l':
1813                                    arcs[aN][1][0] = int(G2img.GetAzm(xpos,ypos,Data))
1814                                else:
1815                                    arcs[aN][1][1] = int(G2img.GetAzm(xpos,ypos,Data))
1816                        for poly in self.polyList:
1817                            if Obj == poly[0]:
1818                                ind = self.itemPicked.contains(self.mousePicked)[1]['ind'][0]
1819                                oldPos = np.array([self.mousePicked.xdata,self.mousePicked.ydata])
1820                                pN = poly[1]
1821                                for i,xy in enumerate(polygons[pN]):
1822                                    if np.allclose(np.array([xy]),oldPos,atol=1.0):
1823                                        polygons[pN][i] = xpos,ypos
1824                        G2imG.UpdateMasks(self,Masks)
1825#                    else:                  #keep for future debugging
1826#                        print str(self.itemPicked),event.xdata,event.ydata,event.button
1827                PlotImage(self,newImage=True)
1828            self.itemPicked = None
1829           
1830    try:
1831        plotNum = self.G2plotNB.plotList.index('2D Powder Image')
1832        Page = self.G2plotNB.nb.GetPage(plotNum)
1833        if not newPlot:
1834            Plot = Page.figure.gca()          #get previous powder plot & get limits
1835            xylim = Plot.get_xlim(),Plot.get_ylim()
1836        if newImage:
1837            Page.figure.clf()
1838            Plot = Page.figure.gca()          #get a fresh plot after clf()
1839       
1840    except ValueError:
1841        Plot = self.G2plotNB.addMpl('2D Powder Image').gca()
1842        plotNum = self.G2plotNB.plotList.index('2D Powder Image')
1843        Page = self.G2plotNB.nb.GetPage(plotNum)
1844        Page.canvas.mpl_connect('key_press_event', OnImPlotKeyPress)
1845        Page.canvas.mpl_connect('motion_notify_event', OnImMotion)
1846        Page.canvas.mpl_connect('pick_event', OnImPick)
1847        Page.canvas.mpl_connect('button_release_event', OnImRelease)
1848        xylim = []
1849    if not event:                       #event from GUI TextCtrl - don't want focus to change to plot!!!
1850        Page.SetFocus()
1851    Title = self.PatternTree.GetItemText(self.Image)[4:]
1852    self.G2plotNB.status.DestroyChildren()
1853    if self.logPlot:
1854        Title = 'log('+Title+')'
1855    Plot.set_title(Title)
1856    try:
1857        if self.PatternTree.GetItemText(self.PickId) in ['Image Controls',]:
1858            Choice = (' key press','l: log(I) on',)
1859            if self.logPlot:
1860                Choice[1] = 'l: log(I) off'
1861            cb = wx.ComboBox(self.G2plotNB.status,style=wx.CB_DROPDOWN|wx.CB_READONLY,
1862                choices=Choice)
1863            cb.Bind(wx.EVT_COMBOBOX, OnKeyBox)
1864            cb.SetValue(' key press')
1865    except TypeError:
1866        pass
1867    size,imagefile = self.PatternTree.GetItemPyData(self.Image)
1868    if imagefile != self.oldImagefile:
1869        imagefile = G2IO.CheckImageFile(self,imagefile)
1870        if not imagefile:
1871            self.G2plotNB.Delete('2D Powder Image')
1872            return
1873        self.PatternTree.SetItemPyData(self.Image,[size,imagefile])
1874        self.ImageZ = G2IO.GetImageData(self,imagefile,imageOnly=True)
1875#        print self.ImageZ.shape,self.ImageZ.size,Data['size'] #might be useful debugging line
1876        self.oldImagefile = imagefile
1877
1878    imScale = 1
1879    if len(self.ImageZ) > 1024:
1880        imScale = len(self.ImageZ)/1024
1881    sizexy = Data['size']
1882    pixelSize = Data['pixelSize']
1883    scalex = 1000./pixelSize[0]
1884    scaley = 1000./pixelSize[1]
1885    Xmax = sizexy[0]*pixelSize[0]/1000.
1886    Ymax = sizexy[1]*pixelSize[1]/1000.
1887    xlim = (0,Xmax)
1888    ylim = (Ymax,0)
1889    Imin,Imax = Data['range'][1]
1890    acolor = mpl.cm.get_cmap(Data['color'])
1891    xcent,ycent = Data['center']
1892    Plot.set_xlabel('Image x-axis, mm',fontsize=12)
1893    Plot.set_ylabel('Image y-axis, mm',fontsize=12)
1894    #do threshold mask - "real" mask - others are just bondaries
1895    Zlim = Masks['Thresholds'][1]
1896    wx.BeginBusyCursor()
1897    try:
1898           
1899        if newImage:                   
1900            MA = ma.masked_greater(ma.masked_less(self.ImageZ,Zlim[0]),Zlim[1])
1901            MaskA = ma.getmaskarray(MA)
1902            A = G2img.ImageCompress(MA,imScale)
1903            AM = G2img.ImageCompress(MaskA,imScale)
1904            if self.logPlot:
1905                A = np.where(A>0,np.log(A),0)
1906                AM = np.where(AM>0,np.log(AM),0)
1907                Imin,Imax = [np.amin(A),np.amax(A)]
1908            ImgM = Plot.imshow(AM,aspect='equal',cmap='Reds',
1909                interpolation='nearest',vmin=0,vmax=2,extent=[0,Xmax,Ymax,0])
1910            Img = Plot.imshow(A,aspect='equal',cmap=acolor,
1911                interpolation='nearest',vmin=Imin,vmax=Imax,extent=[0,Xmax,Ymax,0])
1912            if self.setPoly:
1913                Img.set_picker(True)
1914   
1915        Plot.plot(xcent,ycent,'x')
1916        if Data['showLines']:
1917            LRAzim = Data['LRazimuth']                  #NB: integers
1918            Nazm = Data['outAzimuths']
1919            delAzm = float(LRAzim[1]-LRAzim[0])/Nazm
1920            AzmthOff = Data['azmthOff']
1921            IOtth = Data['IOtth']
1922            wave = Data['wavelength']
1923            dspI = wave/(2.0*sind(IOtth[0]/2.0))
1924            ellI = G2img.GetEllipse(dspI,Data)           #=False if dsp didn't yield an ellipse (ugh! a parabola or a hyperbola)
1925            dspO = wave/(2.0*sind(IOtth[1]/2.0))
1926            ellO = G2img.GetEllipse(dspO,Data)           #Ditto & more likely for outer ellipse
1927            Azm = np.array(range(LRAzim[0],LRAzim[1]+1))-AzmthOff
1928            if ellI:
1929                xyI = []
1930                for azm in Azm:
1931                    xyI.append(G2img.GetDetectorXY(dspI,azm-90.,Data))
1932                xyI = np.array(xyI)
1933                arcxI,arcyI = xyI.T
1934                Plot.plot(arcxI,arcyI,picker=3)
1935            if ellO:
1936                xyO = []
1937                for azm in Azm:
1938                    xyO.append(G2img.GetDetectorXY(dspO,azm-90.,Data))
1939                xyO = np.array(xyO)
1940                arcxO,arcyO = xyO.T
1941                Plot.plot(arcxO,arcyO,picker=3)
1942            if ellO and ellI:
1943                Plot.plot([arcxI[0],arcxO[0]],[arcyI[0],arcyO[0]],picker=3)
1944                Plot.plot([arcxI[-1],arcxO[-1]],[arcyI[-1],arcyO[-1]],picker=3)
1945            for i in range(Nazm):
1946                cake = LRAzim[0]+i*delAzm-AzmthOff
1947                ind = np.searchsorted(Azm,cake)
1948                Plot.plot([arcxI[ind],arcxO[ind]],[arcyI[ind],arcyO[ind]],color='k',dashes=(5,5))
1949                   
1950        for xring,yring in Data['ring']:
1951            Plot.plot(xring,yring,'r+',picker=3)
1952        if Data['setRings']:
1953#            rings = np.concatenate((Data['rings']),axis=0)
1954            N = 0
1955            for ring in Data['rings']:
1956                xring,yring = np.array(ring).T[:2]
1957                Plot.plot(xring,yring,'+',color=colors[N%6])
1958                N += 1           
1959        for ellipse in Data['ellipses']:
1960            cent,phi,[width,height],col = ellipse
1961            Plot.add_artist(Ellipse([cent[0],cent[1]],2*width,2*height,phi,ec=col,fc='none'))
1962            Plot.text(cent[0],cent[1],'+',color=col,ha='center',va='center')
1963        #masks - mask lines numbered after integration limit lines
1964        spots = Masks['Points']
1965        rings = Masks['Rings']
1966        arcs = Masks['Arcs']
1967        polygons = Masks['Polygons']
1968        for x,y,d in spots:
1969            Plot.add_artist(Circle((x,y),radius=d/2,fc='none',ec='r',picker=3))
1970        self.ringList = []
1971        for iring,(tth,thick) in enumerate(rings):
1972            wave = Data['wavelength']
1973            x1,y1 = np.hsplit(np.array(G2img.makeIdealRing(G2img.GetEllipse(Dsp(tth+thick/2.,wave),Data))),2)
1974            x2,y2 = np.hsplit(np.array(G2img.makeIdealRing(G2img.GetEllipse(Dsp(tth-thick/2.,wave),Data))),2)
1975            self.ringList.append([Plot.plot(x1,y1,'r',picker=3),iring,'o'])           
1976            self.ringList.append([Plot.plot(x2,y2,'r',picker=3),iring,'i'])
1977        self.arcList = []
1978        for iarc,(tth,azm,thick) in enumerate(arcs):           
1979            wave = Data['wavelength']
1980            x1,y1 = np.hsplit(np.array(G2img.makeIdealRing(G2img.GetEllipse(Dsp(tth+thick/2.,wave),Data),azm)),2)
1981            x2,y2 = np.hsplit(np.array(G2img.makeIdealRing(G2img.GetEllipse(Dsp(max(0.01,tth-thick/2.),wave),Data),azm)),2)
1982            self.arcList.append([Plot.plot(x1,y1,'r',picker=3),iarc,'o'])           
1983            self.arcList.append([Plot.plot(x2,y2,'r',picker=3),iarc,'i'])
1984            self.arcList.append([Plot.plot([x1[0],x2[0]],[y1[0],y2[0]],'r',picker=3),iarc,'l'])
1985            self.arcList.append([Plot.plot([x1[-1],x2[-1]],[y1[-1],y2[-1]],'r',picker=3),iarc,'u'])
1986        self.polyList = []
1987        for ipoly,polygon in enumerate(polygons):
1988            x,y = np.hsplit(np.array(polygon),2)
1989            self.polyList.append([Plot.plot(x,y,'r+',picker=10),ipoly])
1990            Plot.plot(x,y,'r')           
1991        if newImage:
1992            colorBar = Page.figure.colorbar(Img)
1993        Plot.set_xlim(xlim)
1994        Plot.set_ylim(ylim)
1995        if not newPlot and xylim:
1996            Page.toolbar.push_current()
1997            Plot.set_xlim(xylim[0])
1998            Plot.set_ylim(xylim[1])
1999            xylim = []
2000            Page.toolbar.push_current()
2001            Page.toolbar.draw()
2002        else:
2003            Page.canvas.draw()
2004    finally:
2005        wx.EndBusyCursor()
2006       
2007def PlotIntegration(self,newPlot=False,event=None):
2008    '''Plot of 2D image after image integration with 2-theta and azimuth as coordinates
2009    '''
2010           
2011    def OnMotion(event):
2012        Page.canvas.SetToolTipString('')
2013        Page.canvas.SetCursor(wx.CROSS_CURSOR)
2014        azm = event.ydata
2015        tth = event.xdata
2016        if azm and tth:
2017            self.G2plotNB.status.SetFields(\
2018                ['','Detector 2-th =%9.3fdeg, azm = %7.2fdeg'%(tth,azm)])
2019                               
2020    try:
2021        plotNum = self.G2plotNB.plotList.index('2D Integration')
2022        Page = self.G2plotNB.nb.GetPage(plotNum)
2023        if not newPlot:
2024            Plot = Page.figure.gca()          #get previous plot & get limits
2025            xylim = Plot.get_xlim(),Plot.get_ylim()
2026        Page.figure.clf()
2027        Plot = Page.figure.gca()          #get a fresh plot after clf()
2028       
2029    except ValueError:
2030        Plot = self.G2plotNB.addMpl('2D Integration').gca()
2031        plotNum = self.G2plotNB.plotList.index('2D Integration')
2032        Page = self.G2plotNB.nb.GetPage(plotNum)
2033        Page.canvas.mpl_connect('motion_notify_event', OnMotion)
2034        Page.views = False
2035        view = False
2036    if not event:
2037        Page.SetFocus()
2038       
2039    Data = self.PatternTree.GetItemPyData(
2040        G2gd.GetPatternTreeItemId(self,self.Image, 'Image Controls'))
2041    image = self.Integrate[0]
2042    xsc = self.Integrate[1]
2043    ysc = self.Integrate[2]
2044    Imin,Imax = Data['range'][1]
2045    acolor = mpl.cm.get_cmap(Data['color'])
2046    Plot.set_title(self.PatternTree.GetItemText(self.Image)[4:])
2047    Plot.set_ylabel('azimuth',fontsize=12)
2048    Plot.set_xlabel('2-theta',fontsize=12)
2049    Img = Plot.imshow(image,cmap=acolor,vmin=Imin,vmax=Imax,interpolation='nearest', \
2050        extent=[ysc[0],ysc[-1],xsc[-1],xsc[0]],aspect='auto')
2051    colorBar = Page.figure.colorbar(Img)
2052    if Data['ellipses']:           
2053        for ellipse in Data['ellipses']:
2054            ring = np.array(G2img.makeIdealRing(ellipse[:3])) #skip color
2055            x,y = np.hsplit(ring,2)
2056            tth,azm = G2img.GetTthAzm(x,y,Data)
2057#            azm = np.where(azm < 0.,azm+360,azm)
2058            Plot.plot(tth,azm,'b,')
2059    if not newPlot:
2060        Page.toolbar.push_current()
2061        Plot.set_xlim(xylim[0])
2062        Plot.set_ylim(xylim[1])
2063        xylim = []
2064        Page.toolbar.push_current()
2065        Page.toolbar.draw()
2066    else:
2067        Page.canvas.draw()
2068               
2069def PlotTRImage(self,tax,tay,taz,newPlot=False):
2070    '''a test plot routine - not normally used
2071    ''' 
2072           
2073    def OnMotion(event):
2074        Page.canvas.SetToolTipString('')
2075        Page.canvas.SetCursor(wx.CROSS_CURSOR)
2076        azm = event.xdata
2077        tth = event.ydata
2078        if azm and tth:
2079            self.G2plotNB.status.SetFields(\
2080                ['','Detector 2-th =%9.3fdeg, azm = %7.2fdeg'%(tth,azm)])
2081                               
2082    try:
2083        plotNum = self.G2plotNB.plotList.index('2D Transformed Powder Image')
2084        Page = self.G2plotNB.nb.GetPage(plotNum)
2085        if not newPlot:
2086            Plot = Page.figure.gca()          #get previous plot & get limits
2087            xylim = Plot.get_xlim(),Plot.get_ylim()
2088        Page.figure.clf()
2089        Plot = Page.figure.gca()          #get a fresh plot after clf()
2090       
2091    except ValueError:
2092        Plot = self.G2plotNB.addMpl('2D Transformed Powder Image').gca()
2093        plotNum = self.G2plotNB.plotList.index('2D Transformed Powder Image')
2094        Page = self.G2plotNB.nb.GetPage(plotNum)
2095        Page.canvas.mpl_connect('motion_notify_event', OnMotion)
2096        Page.views = False
2097        view = False
2098    Page.SetFocus()
2099       
2100    Data = self.PatternTree.GetItemPyData(
2101        G2gd.GetPatternTreeItemId(self,self.Image, 'Image Controls'))
2102    Imin,Imax = Data['range'][1]
2103    step = (Imax-Imin)/5.
2104    V = np.arange(Imin,Imax,step)
2105    acolor = mpl.cm.get_cmap(Data['color'])
2106    Plot.set_title(self.PatternTree.GetItemText(self.Image)[4:])
2107    Plot.set_xlabel('azimuth',fontsize=12)
2108    Plot.set_ylabel('2-theta',fontsize=12)
2109    Plot.contour(tax,tay,taz,V,cmap=acolor)
2110    if Data['showLines']:
2111        IOtth = Data['IOtth']
2112        if Data['fullIntegrate']:
2113            LRAzim = [-180,180]
2114        else:
2115            LRAzim = Data['LRazimuth']                  #NB: integers
2116        Plot.plot([LRAzim[0],LRAzim[1]],[IOtth[0],IOtth[0]],picker=True)
2117        Plot.plot([LRAzim[0],LRAzim[1]],[IOtth[1],IOtth[1]],picker=True)
2118        if not Data['fullIntegrate']:
2119            Plot.plot([LRAzim[0],LRAzim[0]],[IOtth[0],IOtth[1]],picker=True)
2120            Plot.plot([LRAzim[1],LRAzim[1]],[IOtth[0],IOtth[1]],picker=True)
2121    if Data['setRings']:
2122        rings = np.concatenate((Data['rings']),axis=0)
2123        for xring,yring,dsp in rings:
2124            x,y = G2img.GetTthAzm(xring,yring,Data)
2125            Plot.plot(y,x,'r+')           
2126    if Data['ellipses']:           
2127        for ellipse in Data['ellipses']:
2128            ring = np.array(G2img.makeIdealRing(ellipse[:3])) #skip color
2129            x,y = np.hsplit(ring,2)
2130            tth,azm = G2img.GetTthAzm(x,y,Data)
2131            Plot.plot(azm,tth,'b,')
2132    if not newPlot:
2133        Page.toolbar.push_current()
2134        Plot.set_xlim(xylim[0])
2135        Plot.set_ylim(xylim[1])
2136        xylim = []
2137        Page.toolbar.push_current()
2138        Page.toolbar.draw()
2139    else:
2140        Page.canvas.draw()
2141       
2142def PlotStructure(self,data):
2143    '''Crystal structure plotting package. Can show structures as balls, sticks, lines,
2144    thermal motion ellipsoids and polyhedra
2145    '''
2146    generalData = data['General']
2147    cell = generalData['Cell'][1:7]
2148    Amat,Bmat = G2lat.cell2AB(cell)         #Amat - crystal to cartesian, Bmat - inverse
2149    A4mat = np.concatenate((np.concatenate((Amat,[[0],[0],[0]]),axis=1),[[0,0,0,1],]),axis=0)
2150    B4mat = np.concatenate((np.concatenate((Bmat,[[0],[0],[0]]),axis=1),[[0,0,0,1],]),axis=0)
2151    Mydir = generalData['Mydir']
2152    atomData = data['Atoms']
2153    drawingData = data['Drawing']
2154    drawAtoms = drawingData['Atoms']
2155    cx,ct,cs = drawingData['atomPtrs']
2156    Wt = [255,255,255]
2157    Rd = [255,0,0]
2158    Gr = [0,255,0]
2159    Bl = [0,0,255]
2160    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]])
2161    uEdges = np.array([
2162        [uBox[0],uBox[1]],[uBox[0],uBox[3]],[uBox[0],uBox[4]],[uBox[1],uBox[2]], 
2163        [uBox[2],uBox[3]],[uBox[1],uBox[5]],[uBox[2],uBox[6]],[uBox[3],uBox[7]], 
2164        [uBox[4],uBox[5]],[uBox[5],uBox[6]],[uBox[6],uBox[7]],[uBox[7],uBox[4]]])
2165    uColors = [Rd,Gr,Bl,Wt, Wt,Wt,Wt,Wt, Wt,Wt,Wt,Wt]
2166    altDown = False
2167    shiftDown = False
2168    ctrlDown = False
2169   
2170    def OnKeyBox(event):
2171        import Image
2172        Draw()                          #make sure plot is fresh!!
2173        mode = cb.GetValue()
2174        Fname = Mydir+'\\'+generalData['Name']+'.'+mode
2175        size = Page.canvas.GetSize()
2176        glPixelStorei(GL_UNPACK_ALIGNMENT, 1)
2177        if mode in ['jpeg',]:
2178            Pix = glReadPixels(0,0,size[0],size[1],GL_RGBA, GL_UNSIGNED_BYTE)
2179            im = Image.new("RGBA", (size[0],size[1]))
2180        else:
2181            Pix = glReadPixels(0,0,size[0],size[1],GL_RGB, GL_UNSIGNED_BYTE)
2182            im = Image.new("RGB", (size[0],size[1]))
2183        im.fromstring(Pix)
2184        im.save(Fname,mode)
2185        cb.SetValue(' Save as:')
2186        self.G2plotNB.status.SetStatusText('Drawing saved to: '+Fname,1)
2187   
2188    def GetTruePosition(xy):
2189        View = glGetIntegerv(GL_VIEWPORT)
2190        Proj = glGetDoublev(GL_PROJECTION_MATRIX)
2191        Model = glGetDoublev(GL_MODELVIEW_MATRIX)
2192        Zmax = 1.
2193        for i,atom in enumerate(drawAtoms):
2194            x,y,z = atom[cx:cx+3]
2195            X,Y,Z = gluProject(x,y,z,Model,Proj,View)
2196            XY = [int(X),int(View[3]-Y)]
2197            if np.allclose(xy,XY,atol=10) and Z < Zmax:
2198                Zmax = Z
2199                SetSelectedAtoms(i)
2200                   
2201    def OnMouseDown(event):
2202        xy = event.GetPosition()
2203        if event.ShiftDown():
2204            GetTruePosition(xy)
2205        else:
2206            drawingData['Rotation'][3] = xy
2207        Draw()
2208       
2209    def OnMouseMove(event):
2210        newxy = event.GetPosition()
2211        page = getSelection()
2212        if event.ControlDown() and drawingData['showABC']:
2213            if event.LeftIsDown():
2214                SetTestRot(newxy)
2215            elif event.RightIsDown():
2216                SetTestPos(newxy)
2217            elif event.MiddleIsDown():
2218                SetTestRotZ(newxy)
2219            x,y,z = drawingData['testPos'][0]
2220            self.G2plotNB.status.SetStatusText('moving test point %.4f,%.4f,%.4f'%(x,y,z),1)
2221               
2222               
2223        if event.Dragging() and not event.ControlDown():
2224            if event.LeftIsDown():
2225                SetRotation(newxy)
2226                angX,angY,angZ = drawingData['Rotation'][:3]
2227                self.G2plotNB.status.SetStatusText('New rotation: %.2f, %.2f ,%.2f'%(angX,angY,angZ),1)
2228            elif event.RightIsDown():
2229                SetTranslation(newxy)
2230                Tx,Ty,Tz = drawingData['viewPoint'][0]
2231                self.G2plotNB.status.SetStatusText('New view point: %.4f, %.4f, %.4f'%(Tx,Ty,Tz),1)
2232            elif event.MiddleIsDown():
2233                SetRotationZ(newxy)
2234                angX,angY,angZ = drawingData['Rotation'][:3]
2235                self.G2plotNB.status.SetStatusText('New rotation: %.2f, %.2f, %.2f'%(angX,angY,angZ),1)
2236        Draw()
2237       
2238    def OnMouseWheel(event):
2239        drawingData['cameraPos'] += event.GetWheelRotation()/24
2240        drawingData['cameraPos'] = max(10,min(500,drawingData['cameraPos']))
2241        self.G2plotNB.status.SetStatusText('New camera distance: %.2f'%(drawingData['cameraPos']),1)
2242        page = getSelection()
2243        if page:
2244            if self.dataDisplay.GetPageText(page) == 'Draw Options':
2245                panel = self.dataDisplay.GetPage(page).GetChildren()[0].GetChildren()
2246                names = [child.GetName() for child in panel]
2247                panel[names.index('cameraPos')].SetLabel('Camera Position: '+'%.2f'%(drawingData['cameraPos']))
2248                panel[names.index('cameraSlider')].SetValue(drawingData['cameraPos'])
2249        Draw()
2250       
2251    def getSelection():
2252        try:
2253            return self.dataDisplay.GetSelection()
2254        except AttributeError:
2255            print self.dataDisplay.GetLabel()
2256            self.G2plotNB.status.SetStatusText('Select this from Phase data window!')
2257            return 0
2258           
2259    def SetViewPointText(VP):
2260        page = getSelection()
2261        if page:
2262            if self.dataDisplay.GetPageText(page) == 'Draw Options':
2263                panel = self.dataDisplay.GetPage(page).GetChildren()[0].GetChildren()
2264                names = [child.GetName() for child in panel]
2265                panel[names.index('viewPoint')].SetValue('%.3f, %.3f, %.3f'%(VP[0],VP[1],VP[2]))
2266           
2267    def ClearSelectedAtoms():
2268        page = getSelection()
2269        if page:
2270            if self.dataDisplay.GetPageText(page) == 'Draw Atoms':
2271                self.dataDisplay.GetPage(page).ClearSelection()      #this is the Atoms grid in Draw Atoms
2272            elif self.dataDisplay.GetPageText(page) == 'Atoms':
2273                self.dataDisplay.GetPage(page).ClearSelection()      #this is the Atoms grid in Atoms
2274                   
2275    def SetSelectedAtoms(ind):
2276        page = getSelection()
2277        if page:
2278            if self.dataDisplay.GetPageText(page) == 'Draw Atoms':
2279                self.dataDisplay.GetPage(page).SelectRow(ind)      #this is the Atoms grid in Draw Atoms
2280            elif self.dataDisplay.GetPageText(page) == 'Atoms':
2281                Id = drawAtoms[ind][-2]
2282                for i,atom in enumerate(atomData):
2283                    if atom[-1] == Id:
2284                        self.dataDisplay.GetPage(page).SelectRow(i)      #this is the Atoms grid in Atoms
2285                 
2286    def GetSelectedAtoms():
2287        page = getSelection()
2288        Ind = []
2289        if page:
2290            if self.dataDisplay.GetPageText(page) == 'Draw Atoms':
2291                Ind = self.dataDisplay.GetPage(page).GetSelectedRows()      #this is the Atoms grid in Draw Atoms
2292            elif self.dataDisplay.GetPageText(page) == 'Atoms':
2293                Ind = self.dataDisplay.GetPage(page).GetSelectedRows()      #this is the Atoms grid in Atoms
2294        return Ind
2295                                       
2296    def OnKey(event):           #on key UP!!
2297        keyCode = event.GetKeyCode()
2298        if keyCode > 255:
2299            keyCode = 0
2300        key,xyz = chr(keyCode),event.GetPosition()
2301        indx = drawingData['selectedAtoms']
2302        if key in ['c','C']:
2303            drawingData['viewPoint'] = [[.5,.5,.5],[0,0]]
2304            drawingData['testPos'] = [[-.1,-.1,-.1],[0.0,0.0,0.0],[0,0]]
2305            drawingData['Rotation'] = [0.0,0.0,0.0,[]]
2306            SetViewPointText(drawingData['viewPoint'][0])
2307        elif key in ['n','N']:
2308            drawAtoms = drawingData['Atoms']
2309            pI = drawingData['viewPoint'][1]
2310            if indx:
2311                pI[0] = indx[pI[1]]
2312                Tx,Ty,Tz = drawAtoms[pI[0]][cx:cx+3]
2313                pI[1] += 1
2314                if pI[1] >= len(indx):
2315                    pI[1] = 0
2316            else:
2317                Tx,Ty,Tz = drawAtoms[pI[0]][cx:cx+3]               
2318                pI[0] += 1
2319                if pI[0] >= len(drawAtoms):
2320                    pI[0] = 0
2321            drawingData['viewPoint'] = [[Tx,Ty,Tz],pI]
2322            SetViewPointText(drawingData['viewPoint'][0])
2323               
2324        elif key in ['p','P']:
2325            drawAtoms = drawingData['Atoms']
2326            pI = drawingData['viewPoint'][1]
2327            if indx:
2328                pI[0] = indx[pI[1]]
2329                Tx,Ty,Tz = drawAtoms[pI[0]][cx:cx+3]
2330                pI[1] -= 1
2331                if pI[1] < 0:
2332                    pI[1] = len(indx)-1
2333            else:
2334                Tx,Ty,Tz = drawAtoms[pI[0]][cx:cx+3]               
2335                pI[0] -= 1
2336                if pI[0] < 0:
2337                    pI[0] = len(drawAtoms)-1
2338            drawingData['viewPoint'] = [[Tx,Ty,Tz],pI]
2339            SetViewPointText(drawingData['viewPoint'][0])           
2340        Draw()
2341           
2342    def SetBackground():
2343        R,G,B,A = Page.camera['backColor']
2344        glClearColor(R,G,B,A)
2345        glClear(GL_COLOR_BUFFER_BIT | GL_DEPTH_BUFFER_BIT)
2346       
2347    def SetLights():
2348        glEnable(GL_DEPTH_TEST)
2349        glShadeModel(GL_SMOOTH)
2350        glEnable(GL_LIGHTING)
2351        glEnable(GL_LIGHT0)
2352        glLightModeli(GL_LIGHT_MODEL_TWO_SIDE,0)
2353        glLightfv(GL_LIGHT0,GL_AMBIENT,[1,1,1,.8])
2354        glLightfv(GL_LIGHT0,GL_DIFFUSE,[1,1,1,1])
2355       
2356    def SetTranslation(newxy):
2357        Tx,Ty,Tz = drawingData['viewPoint'][0]
2358        anglex,angley,anglez,oldxy = drawingData['Rotation']
2359        Rx = G2lat.rotdMat(anglex,0)
2360        Ry = G2lat.rotdMat(angley,1)
2361        Rz = G2lat.rotdMat(anglez,2)
2362        dxy = list(newxy-oldxy)+[0,]
2363        dxy = np.inner(Bmat,np.inner(Rz,np.inner(Ry,np.inner(Rx,dxy))))
2364        Tx -= dxy[0]*0.01
2365        Ty += dxy[1]*0.01
2366        Tz -= dxy[2]*0.01
2367        drawingData['Rotation'][3] = newxy
2368        drawingData['viewPoint'][0] =  Tx,Ty,Tz
2369        SetViewPointText([Tx,Ty,Tz])
2370       
2371    def SetTestPos(newxy):
2372        Tx,Ty,Tz = drawingData['testPos'][0]
2373        anglex,angley,anglez,oldxy = drawingData['Rotation']
2374        Rx = G2lat.rotdMat(anglex,0)
2375        Ry = G2lat.rotdMat(angley,1)
2376        Rz = G2lat.rotdMat(anglez,2)
2377        dxy = list(newxy-oldxy)+[0,]
2378        dxy = np.inner(Rz,np.inner(Ry,np.inner(Rx,dxy)))
2379        Tx += dxy[0]*0.001
2380        Ty -= dxy[1]*0.001
2381        Tz += dxy[2]*0.001
2382        drawingData['Rotation'][3] = newxy
2383        drawingData['testPos'][0] =  Tx,Ty,Tz
2384       
2385    def SetTestRot(newxy):
2386        Txyz = np.array(drawingData['testPos'][0])
2387        oldxy = drawingData['testPos'][2]
2388        Ax,Ay,Az = drawingData['testPos'][1]
2389        Vxyz = np.array(drawingData['viewPoint'][0])
2390        Dxyz = np.inner(Amat,Txyz-Vxyz)
2391        dxy = list(newxy-oldxy)+[0,]
2392        Ax += dxy[1]*0.01
2393        Ay += dxy[0]*0.01
2394        Rx = G2lat.rotdMat(Ax,0)
2395        Ry = G2lat.rotdMat(Ay,1)
2396        Dxyz = np.inner(Ry,np.inner(Rx,Dxyz))       
2397        Dxyz = np.inner(Bmat,Dxyz)+Vxyz
2398        drawingData['testPos'][1] = [Ax,Ay,Az]
2399        drawingData['testPos'][2] = newxy
2400        drawingData['testPos'][0] = Dxyz
2401       
2402    def SetTestRotZ(newxy):
2403        Txyz = np.array(drawingData['testPos'][0])
2404        oldxy = drawingData['testPos'][2]
2405        Ax,Ay,Az = drawingData['testPos'][1]
2406        Vxyz = np.array(drawingData['viewPoint'][0])
2407        Dxyz = np.inner(Amat,Txyz-Vxyz)       
2408        dxy = list(newxy-oldxy)+[0,]
2409        Az += (dxy[0]+dxy[1])*.01
2410        Rz = G2lat.rotdMat(Az,2)
2411        Dxyz = np.inner(Rz,Dxyz)       
2412        Dxyz = np.inner(Bmat,Dxyz)+Vxyz
2413        drawingData['testPos'][1] = [Ax,Ay,Az]
2414        drawingData['testPos'][2] = newxy
2415        drawingData['testPos'][0] = Dxyz
2416                             
2417    def SetRotation(newxy):       
2418        anglex,angley,anglez,oldxy = drawingData['Rotation']
2419        dxy = newxy-oldxy
2420        anglex += dxy[1]*.25
2421        angley += dxy[0]*.25
2422        oldxy = newxy
2423        drawingData['Rotation'] = [anglex,angley,anglez,oldxy]
2424       
2425    def SetRotationZ(newxy):                       
2426        anglex,angley,anglez,oldxy = drawingData['Rotation']
2427        dxy = newxy-oldxy
2428        anglez += (dxy[0]+dxy[1])*.25
2429        oldxy = newxy
2430        drawingData['Rotation'] = [anglex,angley,anglez,oldxy]
2431       
2432    def RenderBox():
2433        glEnable(GL_COLOR_MATERIAL)
2434        glLineWidth(2)
2435        glEnable(GL_BLEND)
2436        glBlendFunc(GL_SRC_ALPHA,GL_ONE_MINUS_SRC_ALPHA)
2437        glEnable(GL_LINE_SMOOTH)
2438        glBegin(GL_LINES)
2439        for line,color in zip(uEdges,uColors):
2440            glColor3ubv(color)
2441            glVertex3fv(line[0])
2442            glVertex3fv(line[1])
2443        glEnd()
2444        glColor4ubv([0,0,0,0])
2445        glDisable(GL_LINE_SMOOTH)
2446        glDisable(GL_BLEND)
2447        glDisable(GL_COLOR_MATERIAL)
2448       
2449    def RenderUnitVectors(x,y,z):
2450        xyz = np.array([x,y,z])
2451        glEnable(GL_COLOR_MATERIAL)
2452        glLineWidth(1)
2453        glPushMatrix()
2454        glTranslate(x,y,z)
2455        glScalef(1/cell[0],1/cell[1],1/cell[2])
2456        glBegin(GL_LINES)
2457        for line,color in zip(uEdges,uColors)[:3]:
2458            glColor3ubv(color)
2459            glVertex3fv(line[0])
2460            glVertex3fv(line[1])
2461        glEnd()
2462        glPopMatrix()
2463        glColor4ubv([0,0,0,0])
2464        glDisable(GL_COLOR_MATERIAL)
2465               
2466    def RenderSphere(x,y,z,radius,color):
2467        glMaterialfv(GL_FRONT_AND_BACK,GL_DIFFUSE,color)
2468        glPushMatrix()
2469        glTranslate(x,y,z)
2470        glMultMatrixf(B4mat.T)
2471        q = gluNewQuadric()
2472        gluSphere(q,radius,20,10)
2473        glPopMatrix()
2474       
2475    def RenderEllipsoid(x,y,z,ellipseProb,E,R4,color):
2476        s1,s2,s3 = E
2477        glMaterialfv(GL_FRONT_AND_BACK,GL_DIFFUSE,color)
2478        glPushMatrix()
2479        glTranslate(x,y,z)
2480        glMultMatrixf(B4mat.T)
2481        glMultMatrixf(R4.T)
2482        glEnable(GL_NORMALIZE)
2483        glScale(s1,s2,s3)
2484        q = gluNewQuadric()
2485        gluSphere(q,ellipseProb,20,10)
2486        glDisable(GL_NORMALIZE)
2487        glPopMatrix()
2488       
2489    def RenderBonds(x,y,z,Bonds,radius,color,slice=20):
2490        glMaterialfv(GL_FRONT_AND_BACK,GL_DIFFUSE,color)
2491        glPushMatrix()
2492        glTranslate(x,y,z)
2493        glMultMatrixf(B4mat.T)
2494        for bond in Bonds:
2495            glPushMatrix()
2496            Dx = np.inner(Amat,bond)
2497            Z = np.sqrt(np.sum(Dx**2))
2498            azm = atan2d(-Dx[1],-Dx[0])
2499            phi = acosd(Dx[2]/Z)
2500            glRotate(-azm,0,0,1)
2501            glRotate(phi,1,0,0)
2502            q = gluNewQuadric()
2503            gluCylinder(q,radius,radius,Z,slice,2)
2504            glPopMatrix()           
2505        glPopMatrix()
2506               
2507    def RenderLines(x,y,z,Bonds,color):
2508        xyz = np.array([x,y,z])
2509        glEnable(GL_COLOR_MATERIAL)
2510        glLineWidth(1)
2511        glColor3fv(color)
2512        glPushMatrix()
2513        glBegin(GL_LINES)
2514        for bond in Bonds:
2515            glVertex3fv(xyz)
2516            glVertex3fv(xyz+bond)
2517        glEnd()
2518        glColor4ubv([0,0,0,0])
2519        glPopMatrix()
2520        glDisable(GL_COLOR_MATERIAL)
2521       
2522    def RenderPolyhedra(x,y,z,Faces,color):
2523        glPushMatrix()
2524        glTranslate(x,y,z)
2525        glMaterialfv(GL_FRONT_AND_BACK,GL_DIFFUSE,color)
2526        glShadeModel(GL_SMOOTH)
2527        glMultMatrixf(B4mat.T)
2528        for face,norm in Faces:
2529            glPolygonMode(GL_FRONT_AND_BACK,GL_FILL)
2530            glFrontFace(GL_CW)
2531            glNormal3fv(norm)
2532            glBegin(GL_TRIANGLES)
2533            for vert in face:
2534                glVertex3fv(vert)
2535            glEnd()
2536        glPopMatrix()
2537       
2538    def RenderBackbone(Backbone,BackboneColor,radius):
2539        glPushMatrix()
2540        glMultMatrixf(B4mat.T)
2541        glEnable(GL_COLOR_MATERIAL)
2542        glShadeModel(GL_SMOOTH)
2543        gleSetJoinStyle(TUBE_NORM_EDGE | TUBE_JN_ANGLE | TUBE_JN_CAP)
2544        glePolyCylinder(Backbone,BackboneColor,radius)
2545        glPopMatrix()       
2546        glDisable(GL_COLOR_MATERIAL)
2547       
2548    def RenderLabel(x,y,z,label,r):       
2549        glPushMatrix()
2550        glTranslate(x,y,z)
2551        glMultMatrixf(B4mat.T)
2552        glDisable(GL_LIGHTING)
2553        glColor3f(0,1.,0)
2554        glRasterPos3f(r,r,r)
2555        for c in list(label):
2556            glutBitmapCharacter(GLUT_BITMAP_8_BY_13,ord(c))
2557        glEnable(GL_LIGHTING)
2558        glPopMatrix()
2559                           
2560    def Draw():
2561        import numpy.linalg as nl
2562        Ind = GetSelectedAtoms()
2563        VS = np.array(Page.canvas.GetSize())
2564        aspect = float(VS[0])/float(VS[1])
2565        cPos = drawingData['cameraPos']
2566        Zclip = drawingData['Zclip']*cPos/200.
2567        anglex,angley,anglez = drawingData['Rotation'][:3]
2568        Tx,Ty,Tz = drawingData['viewPoint'][0]
2569        cx,ct,cs = drawingData['atomPtrs']
2570        bondR = drawingData['bondRadius']
2571        G,g = G2lat.cell2Gmat(cell)
2572        GS = G
2573        GS[0][1] = GS[1][0] = math.sqrt(GS[0][0]*GS[1][1])
2574        GS[0][2] = GS[2][0] = math.sqrt(GS[0][0]*GS[2][2])
2575        GS[1][2] = GS[2][1] = math.sqrt(GS[1][1]*GS[2][2])
2576        ellipseProb = G2lat.criticalEllipse(drawingData['ellipseProb']/100.)
2577       
2578        SetBackground()
2579        glInitNames()
2580        glPushName(0)
2581       
2582        glMatrixMode(GL_PROJECTION)
2583        glLoadIdentity()
2584        glViewport(0,0,VS[0],VS[1])
2585        gluPerspective(20.,aspect,cPos-Zclip,cPos+Zclip)
2586        gluLookAt(0,0,cPos,0,0,0,0,1,0)
2587        SetLights()           
2588           
2589        glMatrixMode(GL_MODELVIEW)
2590        glLoadIdentity()
2591        glRotate(anglez,0,0,1)
2592        glRotate(anglex,cosd(anglez),-sind(anglez),0)
2593        glRotate(angley,sind(anglez),cosd(anglez),0)
2594        glMultMatrixf(A4mat.T)
2595        glTranslate(-Tx,-Ty,-Tz)
2596        if drawingData['unitCellBox']:
2597            RenderBox()
2598        if drawingData['showABC']:
2599            x,y,z = drawingData['testPos'][0]
2600#            if altDown:
2601#                self.G2plotNB.status.SetStatusText('moving test point %.4f,%.4f,%.4f'%(x,y,z),1)
2602#            else:
2603#                self.G2plotNB.status.SetStatusText('test point %.4f,%.4f,%.4f'%(x,y,z),1)           
2604            RenderUnitVectors(x,y,z)
2605        Backbone = []
2606        BackboneColor = []
2607        time0 = time.time()
2608        for iat,atom in enumerate(drawingData['Atoms']):
2609            x,y,z = atom[cx:cx+3]
2610            Bonds = atom[-2]
2611            Faces = atom[-1]
2612            try:
2613                atNum = generalData['AtomTypes'].index(atom[ct])
2614            except ValueError:
2615                atNum = -1
2616            CL = atom[cs+2]
2617            color = np.array(CL)/255.
2618            if iat in Ind:
2619                color = np.array(Gr)/255.
2620            radius = 0.5
2621            if atom[cs] != '':
2622                glLoadName(atom[-3])                   
2623            if 'balls' in atom[cs]:
2624                vdwScale = drawingData['vdwScale']
2625                ballScale = drawingData['ballScale']
2626                if atNum < 0:
2627                    radius = 0.2
2628                elif 'H' == atom[ct]:
2629                    if drawingData['showHydrogen']:
2630                        if 'vdW' in atom[cs] and atNum >= 0:
2631                            radius = vdwScale*generalData['vdWRadii'][atNum]
2632                        else:
2633                            radius = ballScale*drawingData['sizeH']
2634                    else:
2635                        radius = 0.0
2636                else:
2637                    if 'vdW' in atom[cs]:
2638                        radius = vdwScale*generalData['vdWRadii'][atNum]
2639                    else:
2640                        radius = ballScale*generalData['BondRadii'][atNum]
2641                RenderSphere(x,y,z,radius,color)
2642                if 'sticks' in atom[cs]:
2643                    RenderBonds(x,y,z,Bonds,bondR,color)
2644            elif 'ellipsoids' in atom[cs]:
2645                RenderBonds(x,y,z,Bonds,bondR,color)
2646                if atom[cs+3] == 'A':                   
2647                    Uij = atom[cs+5:cs+11]
2648                    U = np.multiply(G2spc.Uij2U(Uij),GS)
2649                    U = np.inner(Amat,np.inner(U,Amat).T)
2650                    E,R = nl.eigh(U)
2651                    R4 = np.concatenate((np.concatenate((R,[[0],[0],[0]]),axis=1),[[0,0,0,1],]),axis=0)
2652                    E = np.sqrt(E)
2653                    if atom[ct] == 'H' and not drawingData['showHydrogen']:
2654                        pass
2655                    else:
2656                        RenderEllipsoid(x,y,z,ellipseProb,E,R4,color)                   
2657                else:
2658                    if atom[ct] == 'H' and not drawingData['showHydrogen']:
2659                        pass
2660                    else:
2661                        radius = ellipseProb*math.sqrt(abs(atom[cs+4]))
2662                        RenderSphere(x,y,z,radius,color)
2663            elif 'lines' in atom[cs]:
2664                radius = 0.1
2665                RenderLines(x,y,z,Bonds,color)
2666#                RenderBonds(x,y,z,Bonds,0.05,color,6)
2667            elif atom[cs] == 'sticks':
2668                radius = 0.1
2669                RenderBonds(x,y,z,Bonds,bondR,color)
2670            elif atom[cs] == 'polyhedra':
2671                RenderPolyhedra(x,y,z,Faces,color)
2672            elif atom[cs] == 'backbone':
2673                if atom[ct-1].split()[0] in ['C','N']:
2674                    Backbone.append(list(np.inner(Amat,np.array([x,y,z]))))
2675                    BackboneColor.append(list(color))
2676                   
2677            if atom[cs+1] == 'type':
2678                RenderLabel(x,y,z,atom[ct],radius)
2679            elif atom[cs+1] == 'name':
2680                RenderLabel(x,y,z,atom[ct-1],radius)
2681            elif atom[cs+1] == 'number':
2682                RenderLabel(x,y,z,str(iat+1),radius)
2683            elif atom[cs+1] == 'residue' and atom[ct-1] == 'CA':
2684                RenderLabel(x,y,z,atom[ct-4],radius)
2685            elif atom[cs+1] == '1-letter' and atom[ct-1] == 'CA':
2686                RenderLabel(x,y,z,atom[ct-3],radius)
2687            elif atom[cs+1] == 'chain' and atom[ct-1] == 'CA':
2688                RenderLabel(x,y,z,atom[ct-2],radius)
2689        if Backbone:
2690            RenderBackbone(Backbone,BackboneColor,bondR)
2691#        print time.time()-time0
2692        Page.canvas.SwapBuffers()
2693       
2694    def OnSize(event):
2695        Draw()
2696       
2697    def OnFocus(event):
2698        Draw()
2699       
2700    try:
2701        plotNum = self.G2plotNB.plotList.index(generalData['Name'])
2702        Page = self.G2plotNB.nb.GetPage(plotNum)       
2703    except ValueError:
2704        Plot = self.G2plotNB.addOgl(generalData['Name'])
2705        plotNum = self.G2plotNB.plotList.index(generalData['Name'])
2706        Page = self.G2plotNB.nb.GetPage(plotNum)
2707        Page.views = False
2708        view = False
2709        altDown = False
2710    Page.SetFocus()
2711    cb = wx.ComboBox(self.G2plotNB.status,style=wx.CB_DROPDOWN|wx.CB_READONLY,
2712        choices=(' save as:','jpeg','tiff','bmp'))
2713    cb.Bind(wx.EVT_COMBOBOX, OnKeyBox)
2714    cb.SetValue(' save as:')
2715    Page.canvas.Bind(wx.EVT_MOUSEWHEEL, OnMouseWheel)
2716    Page.canvas.Bind(wx.EVT_LEFT_DOWN, OnMouseDown)
2717    Page.canvas.Bind(wx.EVT_RIGHT_DOWN, OnMouseDown)
2718    Page.canvas.Bind(wx.EVT_MIDDLE_DOWN, OnMouseDown)
2719    Page.canvas.Bind(wx.EVT_KEY_UP, OnKey)
2720    Page.canvas.Bind(wx.EVT_MOTION, OnMouseMove)
2721    Page.canvas.Bind(wx.EVT_SIZE, OnSize)
2722    Page.canvas.Bind(wx.EVT_SET_FOCUS, OnFocus)
2723    Page.camera['position'] = drawingData['cameraPos']
2724    Page.camera['viewPoint'] = np.inner(Amat,drawingData['viewPoint'][0])
2725    Page.camera['backColor'] = np.array(list(drawingData['backColor'])+[0,])/255.
2726    Page.canvas.SetCurrent()
2727    Draw()
2728       
Note: See TracBrowser for help on using the repository browser.