source: trunk/GSASIIplot.py @ 466

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

fix for PE tif file reading
invert y-axis on images

  • Property svn:keywords set to Date Author Revision URL Id
File size: 116.1 KB
Line 
1#GSASII plotting routines
2########### SVN repository information ###################
3# $Date: 2012-02-01 19:21:44 +0000 (Wed, 01 Feb 2012) $
4# $Author: vondreele $
5# $Revision: 466 $
6# $URL: trunk/GSASIIplot.py $
7# $Id: GSASIIplot.py 466 2012-02-01 19:21:44Z 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(__file__)[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   
750def PlotISFG(self,newPlot=False,type=''):
751    ''' PLotting package for PDF analysis; displays I(q), S(q), F(q) and G(r) as single
752    or multiple plots with waterfall and contour plots as options
753    '''
754    if not type:
755        type = self.G2plotNB.plotList[self.G2plotNB.nb.GetSelection()]
756    if type not in ['I(Q)','S(Q)','F(Q)','G(R)']:
757        return
758    superMinusOne = unichr(0xaf)+unichr(0xb9)
759   
760    def OnPlotKeyPress(event):
761        newPlot = False
762        if event.key == 'u':
763            if self.Contour:
764                self.Cmax = min(1.0,self.Cmax*1.2)
765            elif self.Offset[0] < 100.:
766                self.Offset[0] += 1.
767        elif event.key == 'd':
768            if self.Contour:
769                self.Cmax = max(0.0,self.Cmax*0.8)
770            elif self.Offset[0] > 0.:
771                self.Offset[0] -= 1.
772        elif event.key == 'l':
773            self.Offset[1] -= 1.
774        elif event.key == 'r':
775            self.Offset[1] += 1.
776        elif event.key == 'o':
777            self.Offset = [0,0]
778        elif event.key == 'c':
779            newPlot = True
780            if self.Contour:
781                self.Contour = False
782            else:
783                self.Contour = True
784                self.SinglePlot = False
785                self.Offset = [0.,0.]
786        elif event.key == 's':
787            if self.Contour:
788                choice = [m for m in mpl.cm.datad.keys() if not m.endswith("_r")]
789                choice.sort()
790                dlg = wx.SingleChoiceDialog(self,'Select','Color scheme',choice)
791                if dlg.ShowModal() == wx.ID_OK:
792                    sel = dlg.GetSelection()
793                    self.ContourColor = choice[sel]
794                else:
795                    self.ContourColor = 'Paired'
796                dlg.Destroy()
797            else:               
798                if self.SinglePlot:
799                    self.SinglePlot = False
800                else:
801                    self.SinglePlot = True
802        elif event.key == 'i':                  #for smoothing contour plot
803            choice = ['nearest','bilinear','bicubic','spline16','spline36','hanning',
804               'hamming','hermite','kaiser','quadric','catrom','gaussian','bessel',
805               'mitchell','sinc','lanczos']
806            dlg = wx.SingleChoiceDialog(self,'Select','Interpolation',choice)
807            if dlg.ShowModal() == wx.ID_OK:
808                sel = dlg.GetSelection()
809                self.Interpolate = choice[sel]
810            else:
811                self.Interpolate = 'nearest'
812            dlg.Destroy()
813        elif event.key == 't' and not self.Contour:
814            if self.Legend:
815                self.Legend = False
816            else:
817                self.Legend = True
818           
819           
820        PlotISFG(self,newPlot=newPlot,type=type)
821       
822    def OnKeyBox(event):
823        if self.G2plotNB.nb.GetSelection() == self.G2plotNB.plotList.index(type):
824            event.key = cb.GetValue()[0]
825            cb.SetValue(' key press')
826            OnPlotKeyPress(event)
827                       
828    def OnMotion(event):
829        xpos = event.xdata
830        if xpos:                                        #avoid out of frame mouse position
831            ypos = event.ydata
832            Page.canvas.SetCursor(wx.CROSS_CURSOR)
833            try:
834                if self.Contour:
835                    self.G2plotNB.status.SetStatusText('R =%.3fA pattern ID =%5d'%(xpos,int(ypos)),1)
836                else:
837                    self.G2plotNB.status.SetStatusText('R =%.3fA %s =%.2f'%(xpos,type,ypos),1)                   
838            except TypeError:
839                self.G2plotNB.status.SetStatusText('Select '+type+' pattern first',1)
840   
841    xylim = []
842    try:
843        plotNum = self.G2plotNB.plotList.index(type)
844        Page = self.G2plotNB.nb.GetPage(plotNum)
845        if not newPlot:
846            Plot = Page.figure.gca()          #get previous plot & get limits
847            xylim = Plot.get_xlim(),Plot.get_ylim()
848        Page.figure.clf()
849        Plot = Page.figure.gca()
850    except ValueError:
851        newPlot = True
852        self.Cmax = 1.0
853        Plot = self.G2plotNB.addMpl(type).gca()
854        plotNum = self.G2plotNB.plotList.index(type)
855        Page = self.G2plotNB.nb.GetPage(plotNum)
856        Page.canvas.mpl_connect('key_press_event', OnPlotKeyPress)
857        Page.canvas.mpl_connect('motion_notify_event', OnMotion)
858   
859    Page.SetFocus()
860    self.G2plotNB.status.DestroyChildren()
861    if self.Contour:
862        Choice = (' key press','d: lower contour max','u: raise contour max',
863            'i: interpolation method','s: color scheme','c: contour off')
864    else:
865        Choice = (' key press','l: offset left','r: offset right','d: offset down','u: offset up',
866            'o: reset offset','t: toggle legend','c: contour on','s: toggle single plot')
867    cb = wx.ComboBox(self.G2plotNB.status,style=wx.CB_DROPDOWN|wx.CB_READONLY,
868        choices=Choice)
869    cb.Bind(wx.EVT_COMBOBOX, OnKeyBox)
870    cb.SetValue(' key press')
871    PatternId = self.PatternId
872    PickId = self.PickId
873    Plot.set_title(type)
874    if type == 'G(R)':
875        Plot.set_xlabel(r'$R,\AA$',fontsize=14)
876    else:
877        Plot.set_xlabel(r'$Q,\AA$'+superMinusOne,fontsize=14)
878    Plot.set_ylabel(r''+type,fontsize=14)
879    colors=['b','g','r','c','m','k']
880    name = self.PatternTree.GetItemText(PatternId)[4:]
881    Pattern = []   
882    if self.SinglePlot:
883        name = self.PatternTree.GetItemText(PatternId)
884        name = type+name[4:]
885        Id = G2gd.GetPatternTreeItemId(self,PatternId,name)
886        Pattern = self.PatternTree.GetItemPyData(Id)
887        if Pattern:
888            Pattern.append(name)
889        PlotList = [Pattern,]
890    else:
891        PlotList = []
892        item, cookie = self.PatternTree.GetFirstChild(self.root)
893        while item:
894            if 'PDF' in self.PatternTree.GetItemText(item):
895                name = type+self.PatternTree.GetItemText(item)[4:]
896                Id = G2gd.GetPatternTreeItemId(self,item,name)
897                Pattern = self.PatternTree.GetItemPyData(Id)
898                if Pattern:
899                    Pattern.append(name)
900                    PlotList.append(Pattern)
901            item, cookie = self.PatternTree.GetNextChild(self.root, cookie)
902    PDFdata = self.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(self,PatternId, 'PDF Controls'))
903    numbDen = G2pwd.GetNumDensity(PDFdata['ElList'],PDFdata['Form Vol'])
904    Xb = [0.,10.]
905    Yb = [0.,-40.*np.pi*numbDen]
906    Ymax = 1.0
907    lenX = 0
908    for Pattern in PlotList:
909        xye = Pattern[1]
910        Ymax = max(Ymax,max(xye[1]))
911    offset = self.Offset[0]*Ymax/100.0
912    if self.Contour:
913        ContourZ = []
914        ContourY = []
915        Nseq = 0
916    for N,Pattern in enumerate(PlotList):
917        xye = Pattern[1]
918        if PickId:
919            ifpicked = Pattern[2] == self.PatternTree.GetItemText(PatternId)
920        X = xye[0]
921        if not lenX:
922            lenX = len(X)           
923        Y = xye[1]+offset*N
924        if self.Contour:
925            if lenX == len(X):
926                ContourY.append(N)
927                ContourZ.append(Y)
928                ContourX = X
929                Nseq += 1
930                Plot.set_ylabel('Data sequence',fontsize=12)
931        else:
932            X = xye[0]+self.Offset[1]*.005*N
933            if ifpicked:
934                Plot.plot(X,Y,colors[N%6]+'+',picker=3.,clip_on=False)
935                Page.canvas.SetToolTipString('')
936            else:
937                if self.Legend:
938                    Plot.plot(X,Y,colors[N%6],picker=False,label='Azm:'+Pattern[2].split('=')[1])
939                else:
940                    Plot.plot(X,Y,colors[N%6],picker=False)
941            if type == 'G(R)':
942                Plot.plot(Xb,Yb,color='k',dashes=(5,5))
943            elif type == 'F(Q)':
944                Plot.axhline(0.,color=wx.BLACK)
945            elif type == 'S(Q)':
946                Plot.axhline(1.,color=wx.BLACK)
947    if self.Contour:
948        acolor = mpl.cm.get_cmap(self.ContourColor)
949        Img = Plot.imshow(ContourZ,cmap=acolor,vmin=0,vmax=Ymax*self.Cmax,interpolation=self.Interpolate, 
950            extent=[ContourX[0],ContourX[-1],ContourY[0],ContourY[-1]],aspect='auto',origin='lower')
951        Page.figure.colorbar(Img)
952    elif self.Legend:
953        Plot.legend(loc='best')
954    if not newPlot:
955        Page.toolbar.push_current()
956        Plot.set_xlim(xylim[0])
957        Plot.set_ylim(xylim[1])
958        xylim = []
959        Page.toolbar.push_current()
960        Page.toolbar.draw()
961    else:
962        Page.canvas.draw()
963       
964def PlotXY(self,XY,newPlot=False,type=''):
965    '''simple plot of xy data, used for diagnostic purposes
966    '''
967    def OnMotion(event):
968        xpos = event.xdata
969        if xpos:                                        #avoid out of frame mouse position
970            ypos = event.ydata
971            Page.canvas.SetCursor(wx.CROSS_CURSOR)
972            try:
973                self.G2plotNB.status.SetStatusText('X =%9.3f %s =%9.3f'%(xpos,type,ypos),1)                   
974            except TypeError:
975                self.G2plotNB.status.SetStatusText('Select '+type+' pattern first',1)
976
977    try:
978        plotNum = self.G2plotNB.plotList.index(type)
979        Page = self.G2plotNB.nb.GetPage(plotNum)
980        if not newPlot:
981            Plot = Page.figure.gca()
982            xylim = Plot.get_xlim(),Plot.get_ylim()
983        Page.figure.clf()
984        Plot = Page.figure.gca()
985    except ValueError:
986        newPlot = True
987        Plot = self.G2plotNB.addMpl(type).gca()
988        plotNum = self.G2plotNB.plotList.index(type)
989        Page = self.G2plotNB.nb.GetPage(plotNum)
990        Page.canvas.mpl_connect('motion_notify_event', OnMotion)
991   
992    Page.SetFocus()
993    self.G2plotNB.status.DestroyChildren()
994    Plot.set_title(type)
995    Plot.set_xlabel(r'X',fontsize=14)
996    Plot.set_ylabel(r''+type,fontsize=14)
997    colors=['b','g','r','c','m','k']
998    Ymax = 1.0
999    lenX = 0
1000    X,Y = XY[:2]
1001    Ymax = max(Ymax,max(Y))
1002    Plot.plot(X,Y,'k',picker=False)
1003    if not newPlot:
1004        Page.toolbar.push_current()
1005        Plot.set_xlim(xylim[0])
1006        Plot.set_ylim(xylim[1])
1007        xylim = []
1008        Page.toolbar.push_current()
1009        Page.toolbar.draw()
1010    else:
1011        Page.canvas.draw()
1012
1013def PlotPowderLines(self):
1014    ''' plotting of powder lines (i.e. no powder pattern) as sticks
1015    '''
1016    global HKL
1017
1018    def OnMotion(event):
1019        xpos = event.xdata
1020        if xpos:                                        #avoid out of frame mouse position
1021            Page.canvas.SetCursor(wx.CROSS_CURSOR)
1022            self.G2plotNB.status.SetFields(['','2-theta =%9.3f '%(xpos,)])
1023            if self.PickId and self.PatternTree.GetItemText(self.PickId) in ['Index Peak List','Unit Cells List']:
1024                found = []
1025                if len(HKL):
1026                    view = Page.toolbar._views.forward()[0][:2]
1027                    wid = view[1]-view[0]
1028                    found = HKL[np.where(np.fabs(HKL.T[5]-xpos) < 0.002*wid)]
1029                if len(found):
1030                    h,k,l = found[0][:3] 
1031                    Page.canvas.SetToolTipString('%d,%d,%d'%(int(h),int(k),int(l)))
1032                else:
1033                    Page.canvas.SetToolTipString('')
1034
1035    try:
1036        plotNum = self.G2plotNB.plotList.index('Powder Lines')
1037        Page = self.G2plotNB.nb.GetPage(plotNum)
1038        Page.figure.clf()
1039        Plot = Page.figure.gca()
1040    except ValueError:
1041        Plot = self.G2plotNB.addMpl('Powder Lines').gca()
1042        plotNum = self.G2plotNB.plotList.index('Powder Lines')
1043        Page = self.G2plotNB.nb.GetPage(plotNum)
1044        Page.canvas.mpl_connect('motion_notify_event', OnMotion)
1045       
1046    Page.SetFocus()
1047    Plot.set_title('Powder Pattern Lines')
1048    Plot.set_xlabel(r'$\mathsf{2\theta}$',fontsize=14)
1049    PickId = self.PickId
1050    PatternId = self.PatternId
1051    peaks = self.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(self,PatternId, 'Index Peak List'))
1052    for peak in peaks:
1053        Plot.axvline(peak[0],color='b')
1054    HKL = np.array(self.HKL)
1055    for hkl in self.HKL:
1056        Plot.axvline(hkl[5],color='r',dashes=(5,5))
1057    xmin = peaks[0][0]
1058    xmax = peaks[-1][0]
1059    delt = xmax-xmin
1060    xlim = [max(0,xmin-delt/20.),min(180.,xmax+delt/20.)]
1061    Plot.set_xlim(xlim)
1062    Page.canvas.draw()
1063    Page.toolbar.push_current()
1064
1065def PlotPeakWidths(self):
1066    ''' Plotting of instrument broadening terms as function of 2-theta (future TOF)
1067    Seen when "Instrument Parameters" chosen from powder pattern data tree
1068    '''
1069    PatternId = self.PatternId
1070    limitID = G2gd.GetPatternTreeItemId(self,PatternId, 'Limits')
1071    if limitID:
1072        limits = self.PatternTree.GetItemPyData(limitID)
1073    else:
1074        return
1075    instParms = self.PatternTree.GetItemPyData( \
1076        G2gd.GetPatternTreeItemId(self,PatternId, 'Instrument Parameters'))
1077    if instParms[0][0] in ['PXC','PNC']:
1078        lam = instParms[1][1]
1079        if len(instParms[1]) == 13:
1080            GU,GV,GW,LX,LY = instParms[0][6:11]
1081        else:
1082            GU,GV,GW,LX,LY = instParms[0][4:9]
1083    peakID = G2gd.GetPatternTreeItemId(self,PatternId, 'Peak List')
1084    if peakID:
1085        peaks = self.PatternTree.GetItemPyData(peakID)
1086    else:
1087        peaks = []
1088   
1089    try:
1090        plotNum = self.G2plotNB.plotList.index('Peak Widths')
1091        Page = self.G2plotNB.nb.GetPage(plotNum)
1092        Page.figure.clf()
1093        Plot = Page.figure.gca()
1094    except ValueError:
1095        Plot = self.G2plotNB.addMpl('Peak Widths').gca()
1096        plotNum = self.G2plotNB.plotList.index('Peak Widths')
1097        Page = self.G2plotNB.nb.GetPage(plotNum)
1098    Page.SetFocus()
1099   
1100    Page.canvas.SetToolTipString('')
1101    colors=['b','g','r','c','m','k']
1102    Xmin,Xmax = limits[1]
1103    Xmin = min(0.5,max(Xmin,1))
1104    Xmin /= 2
1105    Xmax /= 2
1106    nPts = 100
1107    delt = (Xmax-Xmin)/nPts
1108    thetas = []
1109    for i in range(nPts):
1110        thetas.append(Xmin+i*delt)
1111    X = []
1112    Y = []
1113    Z = []
1114    W = []
1115    sig = lambda Th,U,V,W: 1.17741*math.sqrt(U*tand(Th)**2+V*tand(Th)+W)*math.pi/18000.
1116    gam = lambda Th,X,Y: (X/cosd(Th)+Y*tand(Th))*math.pi/18000.
1117    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.)
1118#    gamFW2 = lambda s,g: math.sqrt(s**2+(0.4654996*g)**2)+.5345004*g  #Ubaldo Bafile - private communication
1119    for theta in thetas:
1120        X.append(4.0*math.pi*sind(theta)/lam)              #q
1121        s = sig(theta,GU,GV,GW)
1122        g = gam(theta,LX,LY)
1123        G = gamFW(g,s)
1124        Y.append(s/tand(theta))
1125        Z.append(g/tand(theta))
1126        W.append(G/tand(theta))
1127    Plot.set_title('Instrument and sample peak widths')
1128    Plot.set_ylabel(r'$\Delta q/q, \Delta d/d$',fontsize=14)
1129    Plot.set_xlabel(r'$q, \AA^{-1}$',fontsize=14)
1130    Plot.plot(X,Y,color='r',label='Gaussian')
1131    Plot.plot(X,Z,color='g',label='Lorentzian')
1132    Plot.plot(X,W,color='b',label='G+L')
1133    X = []
1134    Y = []
1135    Z = []
1136    W = []
1137    V = []
1138    for peak in peaks:
1139        X.append(4.0*math.pi*sind(peak[0]/2.0)/lam)
1140        try:
1141            s = 1.17741*math.sqrt(peak[4])*math.pi/18000.
1142        except ValueError:
1143            s = 0.01
1144        g = peak[6]*math.pi/18000.
1145        G = gamFW(g,s)
1146        Y.append(s/tand(peak[0]/2.))
1147        Z.append(g/tand(peak[0]/2.))
1148        W.append(G/tand(peak[0]/2.))
1149    Plot.plot(X,Y,'+',color='r',label='G peak')
1150    Plot.plot(X,Z,'+',color='g',label='L peak')
1151    Plot.plot(X,W,'+',color='b',label='G+L peak')
1152    Plot.legend(loc='best')
1153    Page.canvas.draw()
1154   
1155def PlotSizeStrainPO(self,data,Start=False):
1156    '''Plot 3D mustrain/size/preferred orientation figure. In this instance data is for a phase
1157    '''
1158   
1159    PatternId = self.PatternId
1160    generalData = data['General']
1161    SGData = generalData['SGData']
1162    SGLaue = SGData['SGLaue']
1163    if Start:                   #initialize the spherical harmonics qlmn arrays
1164        ptx.pyqlmninit()
1165        Start = False
1166    MuStrKeys = G2spc.MustrainNames(SGData)
1167    cell = generalData['Cell'][1:]
1168    A,B = G2lat.cell2AB(cell[:6])
1169    Vol = cell[6]
1170    useList = data['Histograms']
1171    phase = generalData['Name']
1172    plotType = generalData['Data plot type']
1173    plotDict = {'Mustrain':'Mustrain','Size':'Size','Preferred orientation':'Pref.Ori.'}
1174    for ptype in plotDict:
1175        self.G2plotNB.Delete(ptype)       
1176
1177    for item in useList:
1178        if useList[item]['Show']:
1179            break
1180    else:
1181        return            #nothing to show!!
1182   
1183    numPlots = len(useList)
1184
1185    if plotType in ['Mustrain','Size']:
1186        Plot = mp3d.Axes3D(self.G2plotNB.add3D(plotType))
1187    else:
1188        Plot = self.G2plotNB.addMpl(plotType).gca()       
1189    plotNum = self.G2plotNB.plotList.index(plotType)
1190    Page = self.G2plotNB.nb.GetPage(plotNum)
1191    Page.SetFocus()
1192    self.G2plotNB.status.SetStatusText('Adjust frame size to get desired aspect ratio',1)
1193    if not Page.IsShown():
1194        Page.Show()
1195   
1196    for item in useList:
1197        if useList[item]['Show']:
1198            PHI = np.linspace(0.,360.,30,True)
1199            PSI = np.linspace(0.,180.,30,True)
1200            X = np.outer(npsind(PHI),npsind(PSI))
1201            Y = np.outer(npcosd(PHI),npsind(PSI))
1202            Z = np.outer(np.ones(np.size(PHI)),npcosd(PSI))
1203            try:        #temp patch instead of 'mustrain' for old files with 'microstrain'
1204                coeff = useList[item][plotDict[plotType]]
1205            except KeyError:
1206                break
1207            if plotType in ['Mustrain','Size']:
1208                if coeff[0] == 'isotropic':
1209                    X *= coeff[1][0]
1210                    Y *= coeff[1][0]
1211                    Z *= coeff[1][0]                               
1212                elif coeff[0] == 'uniaxial':
1213                   
1214                    def uniaxCalc(xyz,iso,aniso,axes):
1215                        Z = np.array(axes)
1216                        cp = abs(np.dot(xyz,Z))
1217                        sp = np.sqrt(1.-cp**2)
1218                        R = iso*aniso/np.sqrt((iso*cp)**2+(aniso*sp)**2)
1219                        return R*xyz
1220                       
1221                    iso,aniso = coeff[1][:2]
1222                    axes = np.inner(A,np.array(coeff[3]))
1223                    axes /= nl.norm(axes)
1224                    Shkl = np.array(coeff[1])
1225                    XYZ = np.dstack((X,Y,Z))
1226                    XYZ = np.nan_to_num(np.apply_along_axis(uniaxCalc,2,XYZ,iso,aniso,axes))
1227                    X,Y,Z = np.dsplit(XYZ,3)
1228                    X = X[:,:,0]
1229                    Y = Y[:,:,0]
1230                    Z = Z[:,:,0]
1231               
1232                elif coeff[0] == 'ellipsoidal':
1233                   
1234                    def ellipseCalc(xyz,E,R):
1235                        XYZ = xyz*E.T
1236                        return np.inner(XYZ.T,R)
1237                       
1238                    S6 = coeff[4]
1239                    Sij = G2lat.U6toUij(S6)
1240                    E,R = nl.eigh(Sij)
1241                    XYZ = np.dstack((X,Y,Z))
1242                    XYZ = np.nan_to_num(np.apply_along_axis(ellipseCalc,2,XYZ,E,R))
1243                    X,Y,Z = np.dsplit(XYZ,3)
1244                    X = X[:,:,0]
1245                    Y = Y[:,:,0]
1246                    Z = Z[:,:,0]
1247                   
1248                elif coeff[0] == 'generalized':
1249                   
1250                    def genMustrain(xyz,SGData,A,Shkl):
1251                        uvw = np.inner(A.T,xyz)
1252                        Strm = np.array(G2spc.MustrainCoeff(uvw,SGData))
1253                        sum = np.sum(np.multiply(Shkl,Strm))
1254                        sum = np.where(sum > 0.01,sum,0.01)
1255                        sum = np.sqrt(sum)*math.pi/0.018      #centidegrees to radians!
1256                        return sum*xyz
1257                       
1258                    Shkl = np.array(coeff[4])
1259                    if np.any(Shkl):
1260                        XYZ = np.dstack((X,Y,Z))
1261                        XYZ = np.nan_to_num(np.apply_along_axis(genMustrain,2,XYZ,SGData,A,Shkl))
1262                        X,Y,Z = np.dsplit(XYZ,3)
1263                        X = X[:,:,0]
1264                        Y = Y[:,:,0]
1265                        Z = Z[:,:,0]
1266                           
1267                if np.any(X) and np.any(Y) and np.any(Z):
1268                    errFlags = np.seterr(all='ignore')
1269                    Plot.plot_surface(X,Y,Z,rstride=1,cstride=1,color='g',linewidth=1)
1270                    np.seterr(all='ignore')
1271                    xyzlim = np.array([Plot.get_xlim3d(),Plot.get_ylim3d(),Plot.get_zlim3d()]).T
1272                    XYZlim = [min(xyzlim[0]),max(xyzlim[1])]
1273                    Plot.set_xlim3d(XYZlim)
1274                    Plot.set_ylim3d(XYZlim)
1275                    Plot.set_zlim3d(XYZlim)
1276                    Plot.set_aspect('equal')
1277                if plotType == 'Size':
1278                    Plot.set_title('Crystallite size for '+phase+'\n'+coeff[0]+' model')
1279                    Plot.set_xlabel(r'X, $\mu$m')
1280                    Plot.set_ylabel(r'Y, $\mu$m')
1281                    Plot.set_zlabel(r'Z, $\mu$m')
1282                else:   
1283                    Plot.set_title(r'$\mu$strain for '+phase+'\n'+coeff[0]+' model')
1284                    Plot.set_xlabel(r'X, $\mu$strain')
1285                    Plot.set_ylabel(r'Y, $\mu$strain')
1286                    Plot.set_zlabel(r'Z, $\mu$strain')
1287            else:
1288                h,k,l = generalData['POhkl']
1289                if coeff[0] == 'MD':
1290                    print 'March-Dollase preferred orientation plot'
1291               
1292                else:
1293                    PH = np.array(generalData['POhkl'])
1294                    phi,beta = G2lat.CrsAng(PH,cell[:6],SGData)
1295                    SHCoef = {}
1296                    for item in coeff[5]:
1297                        L,N = eval(item.strip('C'))
1298                        SHCoef['C%d,0,%d'%(L,N)] = coeff[5][item]                       
1299                    ODFln = G2lat.Flnh(Start,SHCoef,phi,beta,SGData)
1300                    X = np.linspace(0,90.0,26)
1301                    Y = G2lat.polfcal(ODFln,'0',X,0.0)
1302                    Plot.plot(X,Y,color='k',label=str(PH))
1303                    Plot.legend(loc='best')
1304                    Plot.set_title('Axial distribution for HKL='+str(PH)+' in '+phase)
1305                    Plot.set_xlabel(r'$\psi$',fontsize=16)
1306                    Plot.set_ylabel('MRD',fontsize=14)
1307    Page.canvas.draw()
1308   
1309def PlotTexture(self,data,newPlot=False,Start=False):
1310    '''Pole figure, inverse pole figure, 3D pole distribution and 3D inverse pole distribution
1311    plotting.
1312    dict generalData contains all phase info needed which is in data
1313    '''
1314
1315    shModels = ['cylindrical','none','shear - 2/m','rolling - mmm']
1316    SamSym = dict(zip(shModels,['0','-1','2/m','mmm']))
1317    PatternId = self.PatternId
1318    generalData = data['General']
1319    SGData = generalData['SGData']
1320    textureData = generalData['SH Texture']
1321    self.G2plotNB.Delete('Texture')
1322    if not textureData['Order']:
1323        return                  #no plot!!
1324    SHData = generalData['SH Texture']
1325    SHCoef = SHData['SH Coeff'][1]
1326    cell = generalData['Cell'][1:7]
1327    Amat,Bmat = G2lat.cell2AB(cell)
1328    sq2 = 1.0/math.sqrt(2.0)
1329   
1330    def rp2xyz(r,p):
1331        z = npcosd(r)
1332        xy = np.sqrt(1.-z**2)
1333        return xy*npsind(p),xy*npcosd(p),z
1334           
1335    def OnMotion(event):
1336        SHData = data['General']['SH Texture']
1337        if event.xdata and event.ydata:                 #avoid out of frame errors
1338            xpos = event.xdata
1339            ypos = event.ydata
1340            if 'Inverse' in SHData['PlotType']:
1341                r = xpos**2+ypos**2
1342                if r <= 1.0:
1343                    if 'equal' in self.Projection: 
1344                        r,p = 2.*npasind(np.sqrt(r)*sq2),npatan2d(ypos,xpos)
1345                    else:
1346                        r,p = 2.*npatand(np.sqrt(r)),npatan2d(ypos,xpos)
1347                    ipf = G2lat.invpolfcal(IODFln,SGData,np.array([r,]),np.array([p,]))
1348                    xyz = np.inner(Amat.T,np.array([rp2xyz(r,p)]))
1349                    y,x,z = list(xyz/np.max(np.abs(xyz)))
1350                   
1351                    self.G2plotNB.status.SetFields(['',
1352                        'psi =%9.3f, beta =%9.3f, MRD =%9.3f hkl=%5.2f,%5.2f,%5.2f'%(r,p,ipf,x,y,z)])
1353                                   
1354            elif 'Axial' in SHData['PlotType']:
1355                pass
1356               
1357            else:                       #ordinary pole figure
1358                z = xpos**2+ypos**2
1359                if z <= 1.0:
1360                    z = np.sqrt(z)
1361                    if 'equal' in self.Projection: 
1362                        r,p = 2.*npasind(z*sq2),npatan2d(ypos,xpos)
1363                    else:
1364                        r,p = 2.*npatand(z),npatan2d(ypos,xpos)
1365                    pf = G2lat.polfcal(ODFln,SamSym[textureData['Model']],np.array([r,]),np.array([p,]))
1366                    self.G2plotNB.status.SetFields(['','phi =%9.3f, gam =%9.3f, MRD =%9.3f'%(r,p,pf)])
1367   
1368    if self.Projection == '3D display' and 'Axial'  not in SHData['PlotType']:               
1369        Plot = mp3d.Axes3D(self.G2plotNB.add3D('Texture'))
1370    else:
1371        Plot = self.G2plotNB.addMpl('Texture').gca()
1372    plotNum = self.G2plotNB.plotList.index('Texture')
1373    Page = self.G2plotNB.nb.GetPage(plotNum)
1374    if self.Projection == '3D display':
1375        pass
1376    else:
1377        Page.canvas.mpl_connect('motion_notify_event', OnMotion)
1378
1379    Page.SetFocus()
1380    self.G2plotNB.status.SetFields(['',''])   
1381    PH = np.array(SHData['PFhkl'])
1382    phi,beta = G2lat.CrsAng(PH,cell,SGData)
1383    ODFln = G2lat.Flnh(Start,SHCoef,phi,beta,SGData)
1384    PX = np.array(SHData['PFxyz'])
1385    gam = atan2d(PX[0],PX[1])
1386    xy = math.sqrt(PX[0]**2+PX[1]**2)
1387    xyz = math.sqrt(PX[0]**2+PX[1]**2+PX[2]**2)
1388    psi = asind(xy/xyz)
1389    IODFln = G2lat.Glnh(Start,SHCoef,psi,gam,SamSym[textureData['Model']])
1390    if 'Axial' in SHData['PlotType']:
1391        X = np.linspace(0,90.0,26)
1392        Y = G2lat.polfcal(ODFln,SamSym[textureData['Model']],X,0.0)
1393        Plot.plot(X,Y,color='k',label=str(SHData['PFhkl']))
1394        Plot.legend(loc='best')
1395        Plot.set_title('Axial distribution for HKL='+str(SHData['PFhkl']))
1396        Plot.set_xlabel(r'$\psi$',fontsize=16)
1397        Plot.set_ylabel('MRD',fontsize=14)
1398       
1399    else:       
1400        npts = 201
1401        if 'Inverse' in SHData['PlotType']:
1402            X,Y = np.meshgrid(np.linspace(1.,-1.,npts),np.linspace(-1.,1.,npts))
1403            R,P = np.sqrt(X**2+Y**2).flatten(),npatan2d(X,Y).flatten()
1404            if 'equal' in self.Projection:
1405                R = np.where(R <= 1.,2.*npasind(R*sq2),0.0)
1406            else:
1407                R = np.where(R <= 1.,2.*npatand(R),0.0)
1408            Z = np.zeros_like(R)
1409            Z = G2lat.invpolfcal(IODFln,SGData,R,P)
1410            Z = np.reshape(Z,(npts,npts))
1411            CS = Plot.contour(Y,X,Z,aspect='equal')
1412            Plot.clabel(CS,fontsize=9,inline=1)
1413            try:
1414                Img = Plot.imshow(Z.T,aspect='equal',cmap=self.ContourColor,extent=[-1,1,-1,1])
1415            except ValueError:
1416                pass
1417            if newPlot:
1418#                Page.figure.colorbar(Img)    #colorbar fails - crashes gsasii
1419                newPlot = False
1420            Plot.set_title('Inverse pole figure for XYZ='+str(SHData['PFxyz']))
1421            Plot.set_xlabel(self.Projection.capitalize()+' projection')
1422                       
1423        else:
1424            X,Y = np.meshgrid(np.linspace(1.,-1.,npts),np.linspace(-1.,1.,npts))
1425            R,P = np.sqrt(X**2+Y**2).flatten(),npatan2d(X,Y).flatten()
1426            if 'equal' in self.Projection:
1427                R = np.where(R <= 1.,2.*npasind(R*sq2),0.0)
1428            else:
1429                R = np.where(R <= 1.,2.*npatand(R),0.0)
1430            Z = np.zeros_like(R)
1431            Z = G2lat.polfcal(ODFln,SamSym[textureData['Model']],R,P)
1432            Z = np.reshape(Z,(npts,npts))
1433            try:
1434                CS = Plot.contour(Y,X,Z,aspect='equal')
1435                Plot.clabel(CS,fontsize=9,inline=1)
1436            except ValueError:
1437                pass
1438            Img = Plot.imshow(Z.T,aspect='equal',cmap=self.ContourColor,extent=[-1,1,-1,1])
1439            if newPlot:
1440#                Page.figure.colorbar(Img)    #colorbar fails - crashes gsasii
1441                newPlot = False
1442            Plot.set_title('Pole figure for HKL='+str(SHData['PFhkl']))
1443            Plot.set_xlabel(self.Projection.capitalize()+' projection')
1444    Page.canvas.draw()
1445
1446def PlotCovariance(self,Data={}):
1447    if not Data:
1448        Data = self.PatternTree.GetItemPyData(
1449            G2gd.GetPatternTreeItemId(self,self.root, 'Covariance'))
1450    if not Data:
1451        print 'No covariance matrix available'
1452        return
1453    varyList = Data['varyList']
1454    values = Data['variables']
1455    Xmax = len(varyList)
1456    covMatrix = Data['covMatrix']
1457    sig = np.sqrt(np.diag(covMatrix))
1458    xvar = np.outer(sig,np.ones_like(sig))
1459    covArray = np.divide(np.divide(covMatrix,xvar),xvar.T)
1460    title = ' for\n'+Data['title']
1461    newAtomDict = Data['newAtomDict']
1462
1463    def OnPlotKeyPress(event):
1464        newPlot = False
1465        if event.key == 's':
1466            choice = [m for m in mpl.cm.datad.keys() if not m.endswith("_r")]
1467            choice.sort()
1468            dlg = wx.SingleChoiceDialog(self,'Select','Color scheme',choice)
1469            if dlg.ShowModal() == wx.ID_OK:
1470                sel = dlg.GetSelection()
1471                self.VcovColor = choice[sel]
1472            else:
1473                self.VcovColor = 'RdYlGn'
1474            dlg.Destroy()
1475        PlotCovariance(self)
1476
1477    def OnMotion(event):
1478        if event.button:
1479            ytics = imgAx.get_yticks()
1480            ytics = np.where(ytics<len(varyList),ytics,-1)
1481            ylabs = [np.where(0<=i ,varyList[int(i)],' ') for i in ytics]
1482            imgAx.set_yticklabels(ylabs)
1483           
1484        if event.xdata and event.ydata:                 #avoid out of frame errors
1485            xpos = int(event.xdata+.5)
1486            ypos = int(event.ydata+.5)
1487            if -1 < xpos < len(varyList) and -1 < ypos < len(varyList):
1488                if xpos == ypos:
1489                    value = values[xpos]
1490                    name = varyList[xpos]
1491                    if varyList[xpos] in newAtomDict:
1492                        name,value = newAtomDict[name]                       
1493                    msg = '%s value = %.4g, esd = %.4g'%(name,value,sig[xpos])
1494                else:
1495                    msg = '%s - %s: %5.3f'%(varyList[xpos],varyList[ypos],covArray[xpos][ypos])
1496                Page.canvas.SetToolTipString(msg)
1497                self.G2plotNB.status.SetFields(['Key: s to change colors',msg])
1498    try:
1499        plotNum = self.G2plotNB.plotList.index('Covariance')
1500        Page = self.G2plotNB.nb.GetPage(plotNum)
1501        Page.figure.clf()
1502        Plot = Page.figure.gca()
1503        if not Page.IsShown():
1504            Page.Show()
1505    except ValueError:
1506        Plot = self.G2plotNB.addMpl('Covariance').gca()
1507        plotNum = self.G2plotNB.plotList.index('Covariance')
1508        Page = self.G2plotNB.nb.GetPage(plotNum)
1509        Page.canvas.mpl_connect('motion_notify_event', OnMotion)
1510        Page.canvas.mpl_connect('key_press_event', OnPlotKeyPress)
1511
1512    Page.SetFocus()
1513    self.G2plotNB.status.SetFields(['',''])   
1514    acolor = mpl.cm.get_cmap(self.VcovColor)
1515    Img = Plot.imshow(covArray,aspect='equal',cmap=acolor,interpolation='nearest',origin='lower')
1516    imgAx = Img.get_axes()
1517    ytics = imgAx.get_yticks()
1518    ylabs = [varyList[int(i)] for i in ytics[:-1]]
1519    imgAx.set_yticklabels(ylabs)
1520    colorBar = Page.figure.colorbar(Img)
1521    Plot.set_title('V-Cov matrix'+title)
1522    Plot.set_xlabel('Variable number')
1523    Plot.set_ylabel('Variable name')
1524    Page.canvas.draw()
1525   
1526def PlotSeq(self,SeqData,SeqSig,SeqNames,sampleParm):
1527   
1528    def OnKeyPress(event):
1529        if event.key == 's' and sampleParm:
1530            if self.xAxis:
1531                self.xAxis = False
1532            else:
1533                self.xAxis = True
1534            Draw(False)
1535    try:
1536        plotNum = self.G2plotNB.plotList.index('Sequential refinement')
1537        Page = self.G2plotNB.nb.GetPage(plotNum)
1538        Page.figure.clf()
1539        Plot = Page.figure.gca()
1540        if not Page.IsShown():
1541            Page.Show()
1542    except ValueError:
1543        Plot = self.G2plotNB.addMpl('Sequential refinement').gca()
1544        plotNum = self.G2plotNB.plotList.index('Sequential refinement')
1545        Page = self.G2plotNB.nb.GetPage(plotNum)
1546        Page.canvas.mpl_connect('key_press_event', OnKeyPress)
1547        self.xAxis = False
1548       
1549    def Draw(newPlot):
1550        Page.SetFocus()
1551        self.G2plotNB.status.SetFields(['','press s to toggle x-axis = sample environment parameter'])
1552        if len(SeqData):
1553            Plot.clear()
1554            if self.xAxis:   
1555                xName = sampleParm.keys()[0]
1556                X = sampleParm[xName]
1557            else:
1558                X = np.arange(0,len(SeqData[0]),1)
1559                xName = 'Data sequence number'
1560            for Y,sig,name in zip(SeqData,SeqSig,SeqNames):
1561                Plot.errorbar(X,Y,yerr=sig,label=name)       
1562            Plot.legend(loc='best')
1563            Plot.set_ylabel('Parameter values')
1564            Plot.set_xlabel(xName)
1565            Page.canvas.draw()           
1566    Draw(True)
1567           
1568def PlotExposedImage(self,newPlot=False,event=None):
1569    '''General access module for 2D image plotting
1570    '''
1571    plotNo = self.G2plotNB.nb.GetSelection()
1572    if self.G2plotNB.nb.GetPageText(plotNo) == '2D Powder Image':
1573        PlotImage(self,newPlot,event,newImage=True)
1574    elif self.G2plotNB.nb.GetPageText(plotNo) == '2D Integration':
1575        PlotIntegration(self,newPlot,event)
1576
1577def PlotImage(self,newPlot=False,event=None,newImage=True):
1578    '''Plot of 2D detector images as contoured plot. Also plot calibration ellipses,
1579    masks, etc.
1580    '''
1581    from matplotlib.patches import Ellipse,Arc,Circle,Polygon
1582    import numpy.ma as ma
1583    Dsp = lambda tth,wave: wave/(2.*sind(tth/2.))
1584    global Data,Masks
1585    colors=['b','g','r','c','m','k']
1586    Data = self.PatternTree.GetItemPyData(
1587        G2gd.GetPatternTreeItemId(self,self.Image, 'Image Controls'))
1588    Masks = self.PatternTree.GetItemPyData(
1589        G2gd.GetPatternTreeItemId(self,self.Image, 'Masks'))
1590
1591    def OnImMotion(event):
1592        Page.canvas.SetToolTipString('')
1593        sizexy = Data['size']
1594        if event.xdata and event.ydata:                 #avoid out of frame errors
1595            Page.canvas.SetCursor(wx.CROSS_CURSOR)
1596            item = self.itemPicked
1597            pixelSize = Data['pixelSize']
1598            scalex = 1000./pixelSize[0]
1599            scaley = 1000./pixelSize[1]
1600            if item and self.PatternTree.GetItemText(self.PickId) == 'Image Controls':
1601                if 'Text' in str(item):
1602                    Page.canvas.SetToolTipString('%8.3f %8.3fmm'%(event.xdata,event.ydata))
1603                else:
1604                    xcent,ycent = Data['center']
1605                    xpos = event.xdata-xcent
1606                    ypos = event.ydata-ycent
1607                    tth,azm = G2img.GetTthAzm(event.xdata,event.ydata,Data)
1608                    if 'line3' in  str(item) or 'line4' in str(item) and not Data['fullIntegrate']:
1609                        Page.canvas.SetToolTipString('%6d deg'%(azm))
1610                    elif 'line1' in  str(item) or 'line2' in str(item):
1611                        Page.canvas.SetToolTipString('%8.3fdeg'%(tth))                           
1612            else:
1613                xpos = event.xdata
1614                ypos = event.ydata
1615                xpix = xpos*scalex
1616                ypix = ypos*scaley
1617                Int = 0
1618                if (0 <= xpix <= sizexy[0]) and (0 <= ypix <= sizexy[1]):
1619                    Int = self.ImageZ[ypix][xpix]
1620                tth,azm,dsp = G2img.GetTthAzmDsp(xpos,ypos,Data)
1621                Q = 2.*math.pi/dsp
1622                if self.setPoly:
1623                    self.G2plotNB.status.SetFields(['','Polygon mask pick - LB next point, RB close polygon'])
1624                else:
1625                    self.G2plotNB.status.SetFields(\
1626                        ['','Detector 2-th =%9.3fdeg, dsp =%9.3fA, Q = %6.5fA-1, azm = %7.2fdeg, I = %6d'%(tth,dsp,Q,azm,Int)])
1627
1628    def OnImPlotKeyPress(event):
1629        try:
1630            PickName = self.PatternTree.GetItemText(self.PickId)
1631        except TypeError:
1632            return
1633        if PickName == 'Masks':
1634            Xpos = event.xdata
1635            if not Xpos:            #got point out of frame
1636                return
1637            Ypos = event.ydata
1638            if event.key == 's':
1639                Masks['Points'].append([Xpos,Ypos,1.])
1640            elif event.key == 'r':
1641                tth = G2img.GetTth(Xpos,Ypos,Data)
1642                Masks['Rings'].append([tth,0.1])
1643            elif event.key == 'a':
1644                tth,azm = G2img.GetTthAzm(Xpos,Ypos,Data)
1645                azm = int(azm)               
1646                Masks['Arcs'].append([tth,[azm-5,azm+5],0.1])
1647            elif event.key == 'p':
1648                self.setPoly = True
1649                Masks['Polygons'].append([])
1650                self.G2plotNB.status.SetFields(['','Polygon mask active - LB pick next point, RB close polygon'])
1651            G2imG.UpdateMasks(self,Masks)
1652        elif PickName == 'Image Controls':
1653            if event.key == 'c':
1654                Xpos = event.xdata
1655                if not Xpos:            #got point out of frame
1656                    return
1657                Ypos = event.ydata
1658                dlg = wx.MessageDialog(self,'Are you sure you want to change the center?',
1659                    'Center change',style=wx.OK|wx.CANCEL)
1660                try:
1661                    if dlg.ShowModal() == wx.ID_OK:
1662                        print 'move center to: ',Xpos,Ypos
1663                        Data['center'] = [Xpos,Ypos]
1664                        G2imG.UpdateImageControls(self,Data,Masks)
1665                finally:
1666                    dlg.Destroy()
1667            elif event.key == 'l':
1668                if self.logPlot:
1669                    self.logPlot = False
1670                else:
1671                    self.logPlot = True
1672        PlotImage(self,newImage=True)
1673           
1674    def OnKeyBox(event):
1675        if self.G2plotNB.nb.GetSelection() == self.G2plotNB.plotList.index('2D Powder Image'):
1676            event.key = cb.GetValue()[0]
1677            cb.SetValue(' key press')
1678            if event.key in 'l':
1679                OnImPlotKeyPress(event)
1680                       
1681    def OnImPick(event):
1682        if self.PatternTree.GetItemText(self.PickId) not in ['Image Controls','Masks']:
1683            return
1684        if self.setPoly:
1685            polygon = Masks['Polygons'][-1]
1686            xpos,ypos = event.mouseevent.xdata,event.mouseevent.ydata
1687            if xpos and ypos:                       #point inside image
1688                if len(polygon) > 2 and event.mouseevent.button == 3:
1689                    x0,y0 = polygon[0]
1690                    polygon.append([x0,y0])
1691                    self.setPoly = False
1692                    self.G2plotNB.status.SetFields(['','Polygon closed - RB drag a vertex to change shape'])
1693                else:
1694                    self.G2plotNB.status.SetFields(['','New polygon point: %.1f,%.1f'%(xpos,ypos)])
1695                    polygon.append([xpos,ypos])
1696                G2imG.UpdateMasks(self,Masks)
1697        else:
1698            if self.itemPicked is not None: return
1699            self.itemPicked = event.artist
1700            self.mousePicked = event.mouseevent
1701       
1702    def OnImRelease(event):
1703        try:
1704            PickName = self.PatternTree.GetItemText(self.PickId)
1705        except TypeError:
1706            return
1707        if PickName not in ['Image Controls','Masks']:
1708            return
1709        pixelSize = Data['pixelSize']
1710        scalex = 1000./pixelSize[0]
1711        scaley = 1000./pixelSize[1]
1712        pixLimit = Data['pixLimit']
1713        if self.itemPicked is None and PickName == 'Image Controls':
1714#            sizexy = Data['size']
1715            Xpos = event.xdata
1716            if not (Xpos and self.ifGetRing):                   #got point out of frame
1717                return
1718            Ypos = event.ydata
1719            if Ypos and not Page.toolbar._active:         #make sure zoom/pan not selected
1720                if event.button == 1:
1721                    Xpix = Xpos*scalex
1722                    Ypix = Ypos*scaley
1723                    xpos,ypos,I,J = G2img.ImageLocalMax(self.ImageZ,pixLimit,Xpix,Ypix)
1724                    if I and J:
1725                        xpos += .5                              #shift to pixel center
1726                        ypos += .5
1727                        xpos /= scalex                          #convert to mm
1728                        ypos /= scaley
1729                        Data['ring'].append([xpos,ypos])
1730                elif event.button == 3:
1731                    self.dataFrame.GetStatusBar().SetStatusText('Calibrating...')
1732                    if G2img.ImageCalibrate(self,Data):
1733                        self.dataFrame.GetStatusBar().SetStatusText('Calibration successful - Show ring picks to check')
1734                        print 'Calibration successful'
1735                    else:
1736                        self.dataFrame.GetStatusBar().SetStatusText('Calibration failed - Show ring picks to diagnose')
1737                        print 'Calibration failed'
1738                    self.ifGetRing = False
1739                    G2imG.UpdateImageControls(self,Data,Masks)
1740                    return
1741                PlotImage(self,newImage=False)
1742            return
1743        else:
1744            xpos = event.xdata
1745            if xpos:                                        #avoid out of frame mouse position
1746                ypos = event.ydata
1747                if self.ifGetRing:                          #delete a calibration ring pick
1748                    xypos = [xpos,ypos]
1749                    rings = Data['ring']
1750                    for ring in rings:
1751                        if np.allclose(ring,xypos,.01,0):
1752                            rings.remove(ring)                                                                       
1753                else:
1754                    tth,azm,dsp = G2img.GetTthAzmDsp(xpos,ypos,Data)
1755                    itemPicked = str(self.itemPicked)
1756                    if 'Line2D' in itemPicked and PickName == 'Image Controls':
1757                        if 'line1' in itemPicked:
1758                            Data['IOtth'][0] = max(tth,0.001)
1759                        elif 'line2' in itemPicked:
1760                            Data['IOtth'][1] = tth
1761                        elif 'line3' in itemPicked:
1762                            Data['LRazimuth'][0] = int(azm)
1763                            if Data['fullIntegrate']:
1764                                Data['LRazimuth'][1] = Data['LRazimuth'][0]+360
1765                        elif 'line4' in itemPicked and not Data['fullIntegrate']:
1766                            Data['LRazimuth'][1] = int(azm)
1767                           
1768                        if Data['LRazimuth'][0] > Data['LRazimuth'][1]:
1769                            Data['LRazimuth'][0] -= 360
1770                           
1771                        azLim = np.array(Data['LRazimuth'])   
1772                        if np.any(azLim>360):
1773                            azLim -= 360
1774                            Data['LRazimuth'] = list(azLim)
1775                           
1776                        if  Data['IOtth'][0] > Data['IOtth'][1]:
1777                            Data['IOtth'][0],Data['IOtth'][1] = Data['IOtth'][1],Data['IOtth'][0]
1778                           
1779                        self.InnerTth.SetValue("%8.2f" % (Data['IOtth'][0]))
1780                        self.OuterTth.SetValue("%8.2f" % (Data['IOtth'][1]))
1781                        self.Lazim.SetValue("%6d" % (Data['LRazimuth'][0]))
1782                        self.Razim.SetValue("%6d" % (Data['LRazimuth'][1]))
1783                    elif 'Circle' in itemPicked and PickName == 'Masks':
1784                        spots = Masks['Points']
1785                        newPos = itemPicked.split(')')[0].split('(')[2].split(',')
1786                        newPos = np.array([float(newPos[0]),float(newPos[1])])
1787                        for spot in spots:
1788                            if np.allclose(np.array([spot[:2]]),newPos):
1789                                spot[:2] = xpos,ypos
1790                        G2imG.UpdateMasks(self,Masks)
1791                    elif 'Line2D' in itemPicked and PickName == 'Masks':
1792                        Obj = self.itemPicked.findobj()
1793                        rings = Masks['Rings']
1794                        arcs = Masks['Arcs']
1795                        polygons = Masks['Polygons']
1796                        for ring in self.ringList:
1797                            if Obj == ring[0]:
1798                                rN = ring[1]
1799                                if ring[2] == 'o':
1800                                    rings[rN][0] = G2img.GetTth(xpos,ypos,Data)-rings[rN][1]/2.
1801                                else:
1802                                    rings[rN][0] = G2img.GetTth(xpos,ypos,Data)+rings[rN][1]/2.
1803                        for arc in self.arcList:
1804                            if Obj == arc[0]:
1805                                aN = arc[1]
1806                                if arc[2] == 'o':
1807                                    arcs[aN][0] = G2img.GetTth(xpos,ypos,Data)-arcs[aN][2]/2
1808                                elif arc[2] == 'i':
1809                                    arcs[aN][0] = G2img.GetTth(xpos,ypos,Data)+arcs[aN][2]/2
1810                                elif arc[2] == 'l':
1811                                    arcs[aN][1][0] = int(G2img.GetAzm(xpos,ypos,Data))
1812                                else:
1813                                    arcs[aN][1][1] = int(G2img.GetAzm(xpos,ypos,Data))
1814                        for poly in self.polyList:
1815                            if Obj == poly[0]:
1816                                ind = self.itemPicked.contains(self.mousePicked)[1]['ind'][0]
1817                                oldPos = np.array([self.mousePicked.xdata,self.mousePicked.ydata])
1818                                pN = poly[1]
1819                                for i,xy in enumerate(polygons[pN]):
1820                                    if np.allclose(np.array([xy]),oldPos,atol=1.0):
1821                                        polygons[pN][i] = xpos,ypos
1822                        G2imG.UpdateMasks(self,Masks)
1823#                    else:                  #keep for future debugging
1824#                        print str(self.itemPicked),event.xdata,event.ydata,event.button
1825                PlotImage(self,newImage=True)
1826            self.itemPicked = None
1827           
1828    try:
1829        plotNum = self.G2plotNB.plotList.index('2D Powder Image')
1830        Page = self.G2plotNB.nb.GetPage(plotNum)
1831        if not newPlot:
1832            Plot = Page.figure.gca()          #get previous powder plot & get limits
1833            xylim = Plot.get_xlim(),Plot.get_ylim()
1834        if newImage:
1835            Page.figure.clf()
1836            Plot = Page.figure.gca()          #get a fresh plot after clf()
1837       
1838    except ValueError:
1839        Plot = self.G2plotNB.addMpl('2D Powder Image').gca()
1840        plotNum = self.G2plotNB.plotList.index('2D Powder Image')
1841        Page = self.G2plotNB.nb.GetPage(plotNum)
1842        Page.canvas.mpl_connect('key_press_event', OnImPlotKeyPress)
1843        Page.canvas.mpl_connect('motion_notify_event', OnImMotion)
1844        Page.canvas.mpl_connect('pick_event', OnImPick)
1845        Page.canvas.mpl_connect('button_release_event', OnImRelease)
1846        xylim = []
1847    if not event:                       #event from GUI TextCtrl - don't want focus to change to plot!!!
1848        Page.SetFocus()
1849    Title = self.PatternTree.GetItemText(self.Image)[4:]
1850    self.G2plotNB.status.DestroyChildren()
1851    if self.logPlot:
1852        Title = 'log('+Title+')'
1853    Plot.set_title(Title)
1854    try:
1855        if self.PatternTree.GetItemText(self.PickId) in ['Image Controls',]:
1856            Choice = (' key press','l: log(I) on',)
1857            if self.logPlot:
1858                Choice[1] = 'l: log(I) off'
1859            cb = wx.ComboBox(self.G2plotNB.status,style=wx.CB_DROPDOWN|wx.CB_READONLY,
1860                choices=Choice)
1861            cb.Bind(wx.EVT_COMBOBOX, OnKeyBox)
1862            cb.SetValue(' key press')
1863    except TypeError:
1864        pass
1865    size,imagefile = self.PatternTree.GetItemPyData(self.Image)
1866    if imagefile != self.oldImagefile:
1867        imagefile = G2IO.CheckImageFile(self,imagefile)
1868        if not imagefile:
1869            self.G2plotNB.Delete('2D Powder Image')
1870            return
1871        self.PatternTree.SetItemPyData(self.Image,[size,imagefile])
1872        self.ImageZ = G2IO.GetImageData(self,imagefile,imageOnly=True)
1873#        print self.ImageZ.shape,self.ImageZ.size,Data['size'] #might be useful debugging line
1874        self.oldImagefile = imagefile
1875
1876    imScale = 1
1877    if len(self.ImageZ) > 1024:
1878        imScale = len(self.ImageZ)/1024
1879    sizexy = Data['size']
1880    pixelSize = Data['pixelSize']
1881    scalex = 1000./pixelSize[0]
1882    scaley = 1000./pixelSize[1]
1883    Xmax = sizexy[0]*pixelSize[0]/1000.
1884    Ymax = sizexy[1]*pixelSize[1]/1000.
1885    xlim = (0,Xmax)
1886    ylim = (Ymax,0)
1887    Imin,Imax = Data['range'][1]
1888    acolor = mpl.cm.get_cmap(Data['color'])
1889    xcent,ycent = Data['center']
1890    Plot.set_xlabel('Image x-axis, mm',fontsize=12)
1891    Plot.set_ylabel('Image y-axis, mm',fontsize=12)
1892    #do threshold mask - "real" mask - others are just bondaries
1893    Zlim = Masks['Thresholds'][1]
1894    wx.BeginBusyCursor()
1895    try:
1896           
1897        if newImage:                   
1898            MA = ma.masked_greater(ma.masked_less(self.ImageZ,Zlim[0]),Zlim[1])
1899            MaskA = ma.getmaskarray(MA)
1900            A = G2img.ImageCompress(MA,imScale)
1901            AM = G2img.ImageCompress(MaskA,imScale)
1902            if self.logPlot:
1903                A = np.where(A>0,np.log(A),0)
1904                AM = np.where(AM>0,np.log(AM),0)
1905                Imin,Imax = [np.amin(A),np.amax(A)]
1906            ImgM = Plot.imshow(AM,aspect='equal',cmap='Reds',
1907                interpolation='nearest',vmin=0,vmax=2,extent=[0,Xmax,Ymax,0])
1908            Img = Plot.imshow(A,aspect='equal',cmap=acolor,
1909                interpolation='nearest',vmin=Imin,vmax=Imax,extent=[0,Xmax,Ymax,0])
1910            if self.setPoly:
1911                Img.set_picker(True)
1912   
1913        Plot.plot(xcent,ycent,'x')
1914        if Data['showLines']:
1915            LRAzim = Data['LRazimuth']                  #NB: integers
1916            Nazm = Data['outAzimuths']
1917            delAzm = float(LRAzim[1]-LRAzim[0])/Nazm
1918            AzmthOff = Data['azmthOff']
1919            IOtth = Data['IOtth']
1920            wave = Data['wavelength']
1921            dspI = wave/(2.0*sind(IOtth[0]/2.0))
1922            ellI = G2img.GetEllipse(dspI,Data)           #=False if dsp didn't yield an ellipse (ugh! a parabola or a hyperbola)
1923            dspO = wave/(2.0*sind(IOtth[1]/2.0))
1924            ellO = G2img.GetEllipse(dspO,Data)           #Ditto & more likely for outer ellipse
1925            Azm = np.array(range(LRAzim[0],LRAzim[1]+1))-AzmthOff
1926            if ellI:
1927                xyI = []
1928                for azm in Azm:
1929                    xyI.append(G2img.GetDetectorXY(dspI,azm-90.,Data))
1930                xyI = np.array(xyI)
1931                arcxI,arcyI = xyI.T
1932                Plot.plot(arcxI,arcyI,picker=3)
1933            if ellO:
1934                xyO = []
1935                for azm in Azm:
1936                    xyO.append(G2img.GetDetectorXY(dspO,azm-90.,Data))
1937                xyO = np.array(xyO)
1938                arcxO,arcyO = xyO.T
1939                Plot.plot(arcxO,arcyO,picker=3)
1940            if ellO and ellI:
1941                Plot.plot([arcxI[0],arcxO[0]],[arcyI[0],arcyO[0]],picker=3)
1942                Plot.plot([arcxI[-1],arcxO[-1]],[arcyI[-1],arcyO[-1]],picker=3)
1943            for i in range(Nazm):
1944                cake = LRAzim[0]+i*delAzm-AzmthOff
1945                ind = np.searchsorted(Azm,cake)
1946                Plot.plot([arcxI[ind],arcxO[ind]],[arcyI[ind],arcyO[ind]],color='k',dashes=(5,5))
1947                   
1948        for xring,yring in Data['ring']:
1949            Plot.plot(xring,yring,'r+',picker=3)
1950        if Data['setRings']:
1951#            rings = np.concatenate((Data['rings']),axis=0)
1952            N = 0
1953            for ring in Data['rings']:
1954                xring,yring = np.array(ring).T[:2]
1955                Plot.plot(xring,yring,'+',color=colors[N%6])
1956                N += 1           
1957        for ellipse in Data['ellipses']:
1958            cent,phi,[width,height],col = ellipse
1959            Plot.add_artist(Ellipse([cent[0],cent[1]],2*width,2*height,phi,ec=col,fc='none'))
1960            Plot.text(cent[0],cent[1],'+',color=col,ha='center',va='center')
1961        #masks - mask lines numbered after integration limit lines
1962        spots = Masks['Points']
1963        rings = Masks['Rings']
1964        arcs = Masks['Arcs']
1965        polygons = Masks['Polygons']
1966        for x,y,d in spots:
1967            Plot.add_artist(Circle((x,y),radius=d/2,fc='none',ec='r',picker=3))
1968        self.ringList = []
1969        for iring,(tth,thick) in enumerate(rings):
1970            wave = Data['wavelength']
1971            x1,y1 = np.hsplit(np.array(G2img.makeIdealRing(G2img.GetEllipse(Dsp(tth+thick/2.,wave),Data))),2)
1972            x2,y2 = np.hsplit(np.array(G2img.makeIdealRing(G2img.GetEllipse(Dsp(tth-thick/2.,wave),Data))),2)
1973            self.ringList.append([Plot.plot(x1,y1,'r',picker=3),iring,'o'])           
1974            self.ringList.append([Plot.plot(x2,y2,'r',picker=3),iring,'i'])
1975        self.arcList = []
1976        for iarc,(tth,azm,thick) in enumerate(arcs):           
1977            wave = Data['wavelength']
1978            x1,y1 = np.hsplit(np.array(G2img.makeIdealRing(G2img.GetEllipse(Dsp(tth+thick/2.,wave),Data),azm)),2)
1979            x2,y2 = np.hsplit(np.array(G2img.makeIdealRing(G2img.GetEllipse(Dsp(max(0.01,tth-thick/2.),wave),Data),azm)),2)
1980            self.arcList.append([Plot.plot(x1,y1,'r',picker=3),iarc,'o'])           
1981            self.arcList.append([Plot.plot(x2,y2,'r',picker=3),iarc,'i'])
1982            self.arcList.append([Plot.plot([x1[0],x2[0]],[y1[0],y2[0]],'r',picker=3),iarc,'l'])
1983            self.arcList.append([Plot.plot([x1[-1],x2[-1]],[y1[-1],y2[-1]],'r',picker=3),iarc,'u'])
1984        self.polyList = []
1985        for ipoly,polygon in enumerate(polygons):
1986            x,y = np.hsplit(np.array(polygon),2)
1987            self.polyList.append([Plot.plot(x,y,'r+',picker=10),ipoly])
1988            Plot.plot(x,y,'r')           
1989        if newImage:
1990            colorBar = Page.figure.colorbar(Img)
1991        Plot.set_xlim(xlim)
1992        Plot.set_ylim(ylim)
1993        Plot.invert_yaxis()
1994        if not newPlot and xylim:
1995            Page.toolbar.push_current()
1996            Plot.set_xlim(xylim[0])
1997            Plot.set_ylim(xylim[1])
1998            xylim = []
1999            Page.toolbar.push_current()
2000            Page.toolbar.draw()
2001        else:
2002            Page.canvas.draw()
2003    finally:
2004        wx.EndBusyCursor()
2005       
2006def PlotIntegration(self,newPlot=False,event=None):
2007    '''Plot of 2D image after image integration with 2-theta and azimuth as coordinates
2008    '''
2009           
2010    def OnMotion(event):
2011        Page.canvas.SetToolTipString('')
2012        Page.canvas.SetCursor(wx.CROSS_CURSOR)
2013        azm = event.ydata
2014        tth = event.xdata
2015        if azm and tth:
2016            self.G2plotNB.status.SetFields(\
2017                ['','Detector 2-th =%9.3fdeg, azm = %7.2fdeg'%(tth,azm)])
2018                               
2019    try:
2020        plotNum = self.G2plotNB.plotList.index('2D Integration')
2021        Page = self.G2plotNB.nb.GetPage(plotNum)
2022        if not newPlot:
2023            Plot = Page.figure.gca()          #get previous plot & get limits
2024            xylim = Plot.get_xlim(),Plot.get_ylim()
2025        Page.figure.clf()
2026        Plot = Page.figure.gca()          #get a fresh plot after clf()
2027       
2028    except ValueError:
2029        Plot = self.G2plotNB.addMpl('2D Integration').gca()
2030        plotNum = self.G2plotNB.plotList.index('2D Integration')
2031        Page = self.G2plotNB.nb.GetPage(plotNum)
2032        Page.canvas.mpl_connect('motion_notify_event', OnMotion)
2033        Page.views = False
2034        view = False
2035    if not event:
2036        Page.SetFocus()
2037       
2038    Data = self.PatternTree.GetItemPyData(
2039        G2gd.GetPatternTreeItemId(self,self.Image, 'Image Controls'))
2040    image = self.Integrate[0]
2041    xsc = self.Integrate[1]
2042    ysc = self.Integrate[2]
2043    Imin,Imax = Data['range'][1]
2044    acolor = mpl.cm.get_cmap(Data['color'])
2045    Plot.set_title(self.PatternTree.GetItemText(self.Image)[4:])
2046    Plot.set_ylabel('azimuth',fontsize=12)
2047    Plot.set_xlabel('2-theta',fontsize=12)
2048    Img = Plot.imshow(image,cmap=acolor,vmin=Imin,vmax=Imax,interpolation='nearest', \
2049        extent=[ysc[0],ysc[-1],xsc[-1],xsc[0]],aspect='auto')
2050    colorBar = Page.figure.colorbar(Img)
2051    if Data['ellipses']:           
2052        for ellipse in Data['ellipses']:
2053            ring = np.array(G2img.makeIdealRing(ellipse[:3])) #skip color
2054            x,y = np.hsplit(ring,2)
2055            tth,azm = G2img.GetTthAzm(x,y,Data)
2056#            azm = np.where(azm < 0.,azm+360,azm)
2057            Plot.plot(tth,azm,'b,')
2058    if not newPlot:
2059        Page.toolbar.push_current()
2060        Plot.set_xlim(xylim[0])
2061        Plot.set_ylim(xylim[1])
2062        xylim = []
2063        Page.toolbar.push_current()
2064        Page.toolbar.draw()
2065    else:
2066        Page.canvas.draw()
2067               
2068def PlotTRImage(self,tax,tay,taz,newPlot=False):
2069    '''a test plot routine - not normally used
2070    ''' 
2071           
2072    def OnMotion(event):
2073        Page.canvas.SetToolTipString('')
2074        Page.canvas.SetCursor(wx.CROSS_CURSOR)
2075        azm = event.xdata
2076        tth = event.ydata
2077        if azm and tth:
2078            self.G2plotNB.status.SetFields(\
2079                ['','Detector 2-th =%9.3fdeg, azm = %7.2fdeg'%(tth,azm)])
2080                               
2081    try:
2082        plotNum = self.G2plotNB.plotList.index('2D Transformed Powder Image')
2083        Page = self.G2plotNB.nb.GetPage(plotNum)
2084        if not newPlot:
2085            Plot = Page.figure.gca()          #get previous plot & get limits
2086            xylim = Plot.get_xlim(),Plot.get_ylim()
2087        Page.figure.clf()
2088        Plot = Page.figure.gca()          #get a fresh plot after clf()
2089       
2090    except ValueError:
2091        Plot = self.G2plotNB.addMpl('2D Transformed Powder Image').gca()
2092        plotNum = self.G2plotNB.plotList.index('2D Transformed Powder Image')
2093        Page = self.G2plotNB.nb.GetPage(plotNum)
2094        Page.canvas.mpl_connect('motion_notify_event', OnMotion)
2095        Page.views = False
2096        view = False
2097    Page.SetFocus()
2098       
2099    Data = self.PatternTree.GetItemPyData(
2100        G2gd.GetPatternTreeItemId(self,self.Image, 'Image Controls'))
2101    Imin,Imax = Data['range'][1]
2102    step = (Imax-Imin)/5.
2103    V = np.arange(Imin,Imax,step)
2104    acolor = mpl.cm.get_cmap(Data['color'])
2105    Plot.set_title(self.PatternTree.GetItemText(self.Image)[4:])
2106    Plot.set_xlabel('azimuth',fontsize=12)
2107    Plot.set_ylabel('2-theta',fontsize=12)
2108    Plot.contour(tax,tay,taz,V,cmap=acolor)
2109    if Data['showLines']:
2110        IOtth = Data['IOtth']
2111        if Data['fullIntegrate']:
2112            LRAzim = [-180,180]
2113        else:
2114            LRAzim = Data['LRazimuth']                  #NB: integers
2115        Plot.plot([LRAzim[0],LRAzim[1]],[IOtth[0],IOtth[0]],picker=True)
2116        Plot.plot([LRAzim[0],LRAzim[1]],[IOtth[1],IOtth[1]],picker=True)
2117        if not Data['fullIntegrate']:
2118            Plot.plot([LRAzim[0],LRAzim[0]],[IOtth[0],IOtth[1]],picker=True)
2119            Plot.plot([LRAzim[1],LRAzim[1]],[IOtth[0],IOtth[1]],picker=True)
2120    if Data['setRings']:
2121        rings = np.concatenate((Data['rings']),axis=0)
2122        for xring,yring,dsp in rings:
2123            x,y = G2img.GetTthAzm(xring,yring,Data)
2124            Plot.plot(y,x,'r+')           
2125    if Data['ellipses']:           
2126        for ellipse in Data['ellipses']:
2127            ring = np.array(G2img.makeIdealRing(ellipse[:3])) #skip color
2128            x,y = np.hsplit(ring,2)
2129            tth,azm = G2img.GetTthAzm(x,y,Data)
2130            Plot.plot(azm,tth,'b,')
2131    if not newPlot:
2132        Page.toolbar.push_current()
2133        Plot.set_xlim(xylim[0])
2134        Plot.set_ylim(xylim[1])
2135        xylim = []
2136        Page.toolbar.push_current()
2137        Page.toolbar.draw()
2138    else:
2139        Page.canvas.draw()
2140       
2141def PlotStructure(self,data):
2142    '''Crystal structure plotting package. Can show structures as balls, sticks, lines,
2143    thermal motion ellipsoids and polyhedra
2144    '''
2145    generalData = data['General']
2146    cell = generalData['Cell'][1:7]
2147    Amat,Bmat = G2lat.cell2AB(cell)         #Amat - crystal to cartesian, Bmat - inverse
2148    A4mat = np.concatenate((np.concatenate((Amat,[[0],[0],[0]]),axis=1),[[0,0,0,1],]),axis=0)
2149    B4mat = np.concatenate((np.concatenate((Bmat,[[0],[0],[0]]),axis=1),[[0,0,0,1],]),axis=0)
2150    Mydir = generalData['Mydir']
2151    atomData = data['Atoms']
2152    drawingData = data['Drawing']
2153    drawAtoms = drawingData['Atoms']
2154    cx,ct,cs = drawingData['atomPtrs']
2155    Wt = [255,255,255]
2156    Rd = [255,0,0]
2157    Gr = [0,255,0]
2158    Bl = [0,0,255]
2159    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]])
2160    uEdges = np.array([
2161        [uBox[0],uBox[1]],[uBox[0],uBox[3]],[uBox[0],uBox[4]],[uBox[1],uBox[2]], 
2162        [uBox[2],uBox[3]],[uBox[1],uBox[5]],[uBox[2],uBox[6]],[uBox[3],uBox[7]], 
2163        [uBox[4],uBox[5]],[uBox[5],uBox[6]],[uBox[6],uBox[7]],[uBox[7],uBox[4]]])
2164    uColors = [Rd,Gr,Bl,Wt, Wt,Wt,Wt,Wt, Wt,Wt,Wt,Wt]
2165    altDown = False
2166    shiftDown = False
2167    ctrlDown = False
2168   
2169    def OnKeyBox(event):
2170        import Image
2171        Draw()                          #make sure plot is fresh!!
2172        mode = cb.GetValue()
2173        Fname = Mydir+'\\'+generalData['Name']+'.'+mode
2174        size = Page.canvas.GetSize()
2175        glPixelStorei(GL_UNPACK_ALIGNMENT, 1)
2176        if mode in ['jpeg',]:
2177            Pix = glReadPixels(0,0,size[0],size[1],GL_RGBA, GL_UNSIGNED_BYTE)
2178            im = Image.new("RGBA", (size[0],size[1]))
2179        else:
2180            Pix = glReadPixels(0,0,size[0],size[1],GL_RGB, GL_UNSIGNED_BYTE)
2181            im = Image.new("RGB", (size[0],size[1]))
2182        im.fromstring(Pix)
2183        im.save(Fname,mode)
2184        cb.SetValue(' Save as:')
2185        self.G2plotNB.status.SetStatusText('Drawing saved to: '+Fname,1)
2186   
2187    def GetTruePosition(xy,Add=False):
2188        View = glGetIntegerv(GL_VIEWPORT)
2189        Proj = glGetDoublev(GL_PROJECTION_MATRIX)
2190        Model = glGetDoublev(GL_MODELVIEW_MATRIX)
2191        Zmax = 1.
2192        if Add:
2193            Indx = GetSelectedAtoms()
2194        for i,atom in enumerate(drawAtoms):
2195            x,y,z = atom[cx:cx+3]
2196            X,Y,Z = gluProject(x,y,z,Model,Proj,View)
2197            XY = [int(X),int(View[3]-Y)]
2198            if np.allclose(xy,XY,atol=10) and Z < Zmax:
2199                Zmax = Z
2200                try:
2201                    Indx.remove(i)
2202                    ClearSelectedAtoms()
2203                    for id in Indx:
2204                        SetSelectedAtoms(id,Add)
2205                except:
2206                    SetSelectedAtoms(i,Add)
2207                   
2208    def OnMouseDown(event):
2209        xy = event.GetPosition()
2210        if event.ShiftDown():
2211            if event.LeftIsDown():
2212                GetTruePosition(xy)
2213            elif event.RightIsDown():
2214                GetTruePosition(xy,True)
2215        else:
2216            drawingData['Rotation'][3] = xy
2217        Draw()
2218       
2219    def OnMouseMove(event):
2220        if event.ShiftDown():
2221            return       
2222        newxy = event.GetPosition()
2223        page = getSelection()
2224        if event.ControlDown() and drawingData['showABC']:
2225            if event.LeftIsDown():
2226                SetTestRot(newxy)
2227            elif event.RightIsDown():
2228                SetTestPos(newxy)
2229            elif event.MiddleIsDown():
2230                SetTestRotZ(newxy)
2231            x,y,z = drawingData['testPos'][0]
2232            self.G2plotNB.status.SetStatusText('moving test point %.4f,%.4f,%.4f'%(x,y,z),1)
2233                               
2234        if event.Dragging() and not event.ControlDown():
2235            if event.LeftIsDown():
2236                SetRotation(newxy)
2237                angX,angY,angZ = drawingData['Rotation'][:3]
2238                self.G2plotNB.status.SetStatusText('New rotation: %.2f, %.2f ,%.2f'%(angX,angY,angZ),1)
2239            elif event.RightIsDown():
2240                SetTranslation(newxy)
2241                Tx,Ty,Tz = drawingData['viewPoint'][0]
2242                self.G2plotNB.status.SetStatusText('New view point: %.4f, %.4f, %.4f'%(Tx,Ty,Tz),1)
2243            elif event.MiddleIsDown():
2244                SetRotationZ(newxy)
2245                angX,angY,angZ = drawingData['Rotation'][:3]
2246                self.G2plotNB.status.SetStatusText('New rotation: %.2f, %.2f, %.2f'%(angX,angY,angZ),1)
2247        Draw()
2248       
2249    def OnMouseWheel(event):
2250        if event.ShiftDown():
2251            return
2252        drawingData['cameraPos'] += event.GetWheelRotation()/24
2253        drawingData['cameraPos'] = max(10,min(500,drawingData['cameraPos']))
2254        self.G2plotNB.status.SetStatusText('New camera distance: %.2f'%(drawingData['cameraPos']),1)
2255        page = getSelection()
2256        if page:
2257            if self.dataDisplay.GetPageText(page) == 'Draw Options':
2258                panel = self.dataDisplay.GetPage(page).GetChildren()[0].GetChildren()
2259                names = [child.GetName() for child in panel]
2260                panel[names.index('cameraPos')].SetLabel('Camera Position: '+'%.2f'%(drawingData['cameraPos']))
2261                panel[names.index('cameraSlider')].SetValue(drawingData['cameraPos'])
2262        Draw()
2263       
2264    def getSelection():
2265        try:
2266            return self.dataDisplay.GetSelection()
2267        except AttributeError:
2268            print self.dataDisplay.GetLabel()
2269            self.G2plotNB.status.SetStatusText('Select this from Phase data window!')
2270            return 0
2271           
2272    def SetViewPointText(VP):
2273        page = getSelection()
2274        if page:
2275            if self.dataDisplay.GetPageText(page) == 'Draw Options':
2276                panel = self.dataDisplay.GetPage(page).GetChildren()[0].GetChildren()
2277                names = [child.GetName() for child in panel]
2278                panel[names.index('viewPoint')].SetValue('%.3f, %.3f, %.3f'%(VP[0],VP[1],VP[2]))
2279           
2280    def ClearSelectedAtoms():
2281        page = getSelection()
2282        if page:
2283            if self.dataDisplay.GetPageText(page) == 'Draw Atoms':
2284                self.dataDisplay.GetPage(page).ClearSelection()      #this is the Atoms grid in Draw Atoms
2285            elif self.dataDisplay.GetPageText(page) == 'Atoms':
2286                self.dataDisplay.GetPage(page).ClearSelection()      #this is the Atoms grid in Atoms
2287                   
2288    def SetSelectedAtoms(ind,Add=False):
2289        page = getSelection()
2290        if page:
2291            if self.dataDisplay.GetPageText(page) == 'Draw Atoms':
2292                self.dataDisplay.GetPage(page).SelectRow(ind,Add)      #this is the Atoms grid in Draw Atoms
2293            elif self.dataDisplay.GetPageText(page) == 'Atoms':
2294                Id = drawAtoms[ind][-2]
2295                for i,atom in enumerate(atomData):
2296                    if atom[-1] == Id:
2297                        self.dataDisplay.GetPage(page).SelectRow(i)      #this is the Atoms grid in Atoms
2298                 
2299    def GetSelectedAtoms():
2300        page = getSelection()
2301        Ind = []
2302        if page:
2303            if self.dataDisplay.GetPageText(page) == 'Draw Atoms':
2304                Ind = self.dataDisplay.GetPage(page).GetSelectedRows()      #this is the Atoms grid in Draw Atoms
2305            elif self.dataDisplay.GetPageText(page) == 'Atoms':
2306                Ind = self.dataDisplay.GetPage(page).GetSelectedRows()      #this is the Atoms grid in Atoms
2307        return Ind
2308                                       
2309    def OnKey(event):           #on key UP!!
2310        keyCode = event.GetKeyCode()
2311        if keyCode > 255:
2312            keyCode = 0
2313        key,xyz = chr(keyCode),event.GetPosition()
2314        indx = drawingData['selectedAtoms']
2315        if key in ['c','C']:
2316            drawingData['viewPoint'] = [[.5,.5,.5],[0,0]]
2317            drawingData['testPos'] = [[-.1,-.1,-.1],[0.0,0.0,0.0],[0,0]]
2318            drawingData['Rotation'] = [0.0,0.0,0.0,[]]
2319            SetViewPointText(drawingData['viewPoint'][0])
2320        elif key in ['n','N']:
2321            drawAtoms = drawingData['Atoms']
2322            pI = drawingData['viewPoint'][1]
2323            if indx:
2324                pI[0] = indx[pI[1]]
2325                Tx,Ty,Tz = drawAtoms[pI[0]][cx:cx+3]
2326                pI[1] += 1
2327                if pI[1] >= len(indx):
2328                    pI[1] = 0
2329            else:
2330                Tx,Ty,Tz = drawAtoms[pI[0]][cx:cx+3]               
2331                pI[0] += 1
2332                if pI[0] >= len(drawAtoms):
2333                    pI[0] = 0
2334            drawingData['viewPoint'] = [[Tx,Ty,Tz],pI]
2335            SetViewPointText(drawingData['viewPoint'][0])
2336               
2337        elif key in ['p','P']:
2338            drawAtoms = drawingData['Atoms']
2339            pI = drawingData['viewPoint'][1]
2340            if indx:
2341                pI[0] = indx[pI[1]]
2342                Tx,Ty,Tz = drawAtoms[pI[0]][cx:cx+3]
2343                pI[1] -= 1
2344                if pI[1] < 0:
2345                    pI[1] = len(indx)-1
2346            else:
2347                Tx,Ty,Tz = drawAtoms[pI[0]][cx:cx+3]               
2348                pI[0] -= 1
2349                if pI[0] < 0:
2350                    pI[0] = len(drawAtoms)-1
2351            drawingData['viewPoint'] = [[Tx,Ty,Tz],pI]
2352            SetViewPointText(drawingData['viewPoint'][0])           
2353        Draw()
2354           
2355    def SetBackground():
2356        R,G,B,A = Page.camera['backColor']
2357        glClearColor(R,G,B,A)
2358        glClear(GL_COLOR_BUFFER_BIT | GL_DEPTH_BUFFER_BIT)
2359       
2360    def SetLights():
2361        glEnable(GL_DEPTH_TEST)
2362        glShadeModel(GL_SMOOTH)
2363        glEnable(GL_LIGHTING)
2364        glEnable(GL_LIGHT0)
2365        glLightModeli(GL_LIGHT_MODEL_TWO_SIDE,0)
2366        glLightfv(GL_LIGHT0,GL_AMBIENT,[1,1,1,.8])
2367        glLightfv(GL_LIGHT0,GL_DIFFUSE,[1,1,1,1])
2368       
2369    def SetTranslation(newxy):
2370        Tx,Ty,Tz = drawingData['viewPoint'][0]
2371        anglex,angley,anglez,oldxy = drawingData['Rotation']
2372        Rx = G2lat.rotdMat(anglex,0)
2373        Ry = G2lat.rotdMat(angley,1)
2374        Rz = G2lat.rotdMat(anglez,2)
2375        dxy = list(newxy-oldxy)+[0,]
2376        dxy = np.inner(Bmat,np.inner(Rz,np.inner(Ry,np.inner(Rx,dxy))))
2377        Tx -= dxy[0]*0.01
2378        Ty += dxy[1]*0.01
2379        Tz -= dxy[2]*0.01
2380        drawingData['Rotation'][3] = newxy
2381        drawingData['viewPoint'][0] =  Tx,Ty,Tz
2382        SetViewPointText([Tx,Ty,Tz])
2383       
2384    def SetTestPos(newxy):
2385        Tx,Ty,Tz = drawingData['testPos'][0]
2386        anglex,angley,anglez,oldxy = drawingData['Rotation']
2387        Rx = G2lat.rotdMat(anglex,0)
2388        Ry = G2lat.rotdMat(angley,1)
2389        Rz = G2lat.rotdMat(anglez,2)
2390        dxy = list(newxy-oldxy)+[0,]
2391        dxy = np.inner(Rz,np.inner(Ry,np.inner(Rx,dxy)))
2392        Tx += dxy[0]*0.001
2393        Ty -= dxy[1]*0.001
2394        Tz += dxy[2]*0.001
2395        drawingData['Rotation'][3] = newxy
2396        drawingData['testPos'][0] =  Tx,Ty,Tz
2397       
2398    def SetTestRot(newxy):
2399        Txyz = np.array(drawingData['testPos'][0])
2400        oldxy = drawingData['testPos'][2]
2401        Ax,Ay,Az = drawingData['testPos'][1]
2402        Vxyz = np.array(drawingData['viewPoint'][0])
2403        Dxyz = np.inner(Amat,Txyz-Vxyz)
2404        dxy = list(newxy-oldxy)+[0,]
2405        Ax += dxy[1]*0.01
2406        Ay += dxy[0]*0.01
2407        Rx = G2lat.rotdMat(Ax,0)
2408        Ry = G2lat.rotdMat(Ay,1)
2409        Dxyz = np.inner(Ry,np.inner(Rx,Dxyz))       
2410        Dxyz = np.inner(Bmat,Dxyz)+Vxyz
2411        drawingData['testPos'][1] = [Ax,Ay,Az]
2412        drawingData['testPos'][2] = newxy
2413        drawingData['testPos'][0] = Dxyz
2414       
2415    def SetTestRotZ(newxy):
2416        Txyz = np.array(drawingData['testPos'][0])
2417        oldxy = drawingData['testPos'][2]
2418        Ax,Ay,Az = drawingData['testPos'][1]
2419        Vxyz = np.array(drawingData['viewPoint'][0])
2420        Dxyz = np.inner(Amat,Txyz-Vxyz)       
2421        dxy = list(newxy-oldxy)+[0,]
2422        Az += (dxy[0]+dxy[1])*.01
2423        Rz = G2lat.rotdMat(Az,2)
2424        Dxyz = np.inner(Rz,Dxyz)       
2425        Dxyz = np.inner(Bmat,Dxyz)+Vxyz
2426        drawingData['testPos'][1] = [Ax,Ay,Az]
2427        drawingData['testPos'][2] = newxy
2428        drawingData['testPos'][0] = Dxyz
2429                             
2430    def SetRotation(newxy):       
2431        anglex,angley,anglez,oldxy = drawingData['Rotation']
2432        dxy = newxy-oldxy
2433        anglex += dxy[1]*.25
2434        angley += dxy[0]*.25
2435        oldxy = newxy
2436        drawingData['Rotation'] = [anglex,angley,anglez,oldxy]
2437       
2438    def SetRotationZ(newxy):                       
2439        anglex,angley,anglez,oldxy = drawingData['Rotation']
2440        dxy = newxy-oldxy
2441        anglez += (dxy[0]+dxy[1])*.25
2442        oldxy = newxy
2443        drawingData['Rotation'] = [anglex,angley,anglez,oldxy]
2444       
2445    def RenderBox():
2446        glEnable(GL_COLOR_MATERIAL)
2447        glLineWidth(2)
2448        glEnable(GL_BLEND)
2449        glBlendFunc(GL_SRC_ALPHA,GL_ONE_MINUS_SRC_ALPHA)
2450        glEnable(GL_LINE_SMOOTH)
2451        glBegin(GL_LINES)
2452        for line,color in zip(uEdges,uColors):
2453            glColor3ubv(color)
2454            glVertex3fv(line[0])
2455            glVertex3fv(line[1])
2456        glEnd()
2457        glColor4ubv([0,0,0,0])
2458        glDisable(GL_LINE_SMOOTH)
2459        glDisable(GL_BLEND)
2460        glDisable(GL_COLOR_MATERIAL)
2461       
2462    def RenderUnitVectors(x,y,z):
2463        xyz = np.array([x,y,z])
2464        glEnable(GL_COLOR_MATERIAL)
2465        glLineWidth(1)
2466        glPushMatrix()
2467        glTranslate(x,y,z)
2468        glScalef(1/cell[0],1/cell[1],1/cell[2])
2469        glBegin(GL_LINES)
2470        for line,color in zip(uEdges,uColors)[:3]:
2471            glColor3ubv(color)
2472            glVertex3fv(line[0])
2473            glVertex3fv(line[1])
2474        glEnd()
2475        glPopMatrix()
2476        glColor4ubv([0,0,0,0])
2477        glDisable(GL_COLOR_MATERIAL)
2478               
2479    def RenderSphere(x,y,z,radius,color):
2480        glMaterialfv(GL_FRONT_AND_BACK,GL_DIFFUSE,color)
2481        glPushMatrix()
2482        glTranslate(x,y,z)
2483        glMultMatrixf(B4mat.T)
2484        q = gluNewQuadric()
2485        gluSphere(q,radius,20,10)
2486        glPopMatrix()
2487       
2488    def RenderEllipsoid(x,y,z,ellipseProb,E,R4,color):
2489        s1,s2,s3 = E
2490        glMaterialfv(GL_FRONT_AND_BACK,GL_DIFFUSE,color)
2491        glPushMatrix()
2492        glTranslate(x,y,z)
2493        glMultMatrixf(B4mat.T)
2494        glMultMatrixf(R4.T)
2495        glEnable(GL_NORMALIZE)
2496        glScale(s1,s2,s3)
2497        q = gluNewQuadric()
2498        gluSphere(q,ellipseProb,20,10)
2499        glDisable(GL_NORMALIZE)
2500        glPopMatrix()
2501       
2502    def RenderBonds(x,y,z,Bonds,radius,color,slice=20):
2503        glMaterialfv(GL_FRONT_AND_BACK,GL_DIFFUSE,color)
2504        glPushMatrix()
2505        glTranslate(x,y,z)
2506        glMultMatrixf(B4mat.T)
2507        for bond in Bonds:
2508            glPushMatrix()
2509            Dx = np.inner(Amat,bond)
2510            Z = np.sqrt(np.sum(Dx**2))
2511            azm = atan2d(-Dx[1],-Dx[0])
2512            phi = acosd(Dx[2]/Z)
2513            glRotate(-azm,0,0,1)
2514            glRotate(phi,1,0,0)
2515            q = gluNewQuadric()
2516            gluCylinder(q,radius,radius,Z,slice,2)
2517            glPopMatrix()           
2518        glPopMatrix()
2519               
2520    def RenderLines(x,y,z,Bonds,color):
2521        xyz = np.array([x,y,z])
2522        glEnable(GL_COLOR_MATERIAL)
2523        glLineWidth(1)
2524        glColor3fv(color)
2525        glPushMatrix()
2526        glBegin(GL_LINES)
2527        for bond in Bonds:
2528            glVertex3fv(xyz)
2529            glVertex3fv(xyz+bond)
2530        glEnd()
2531        glColor4ubv([0,0,0,0])
2532        glPopMatrix()
2533        glDisable(GL_COLOR_MATERIAL)
2534       
2535    def RenderPolyhedra(x,y,z,Faces,color):
2536        glPushMatrix()
2537        glTranslate(x,y,z)
2538        glMaterialfv(GL_FRONT_AND_BACK,GL_DIFFUSE,color)
2539        glShadeModel(GL_SMOOTH)
2540        glMultMatrixf(B4mat.T)
2541        for face,norm in Faces:
2542            glPolygonMode(GL_FRONT_AND_BACK,GL_FILL)
2543            glFrontFace(GL_CW)
2544            glNormal3fv(norm)
2545            glBegin(GL_TRIANGLES)
2546            for vert in face:
2547                glVertex3fv(vert)
2548            glEnd()
2549        glPopMatrix()
2550       
2551    def RenderBackbone(Backbone,BackboneColor,radius):
2552        glPushMatrix()
2553        glMultMatrixf(B4mat.T)
2554        glEnable(GL_COLOR_MATERIAL)
2555        glShadeModel(GL_SMOOTH)
2556        gleSetJoinStyle(TUBE_NORM_EDGE | TUBE_JN_ANGLE | TUBE_JN_CAP)
2557        glePolyCylinder(Backbone,BackboneColor,radius)
2558        glPopMatrix()       
2559        glDisable(GL_COLOR_MATERIAL)
2560       
2561    def RenderLabel(x,y,z,label,r):       
2562        glPushMatrix()
2563        glTranslate(x,y,z)
2564        glMultMatrixf(B4mat.T)
2565        glDisable(GL_LIGHTING)
2566        glColor3f(0,1.,0)
2567        glRasterPos3f(r,r,r)
2568        for c in list(label):
2569            glutBitmapCharacter(GLUT_BITMAP_8_BY_13,ord(c))
2570        glEnable(GL_LIGHTING)
2571        glPopMatrix()
2572                           
2573    def Draw():
2574        import numpy.linalg as nl
2575        Ind = GetSelectedAtoms()
2576        VS = np.array(Page.canvas.GetSize())
2577        aspect = float(VS[0])/float(VS[1])
2578        cPos = drawingData['cameraPos']
2579        Zclip = drawingData['Zclip']*cPos/200.
2580        anglex,angley,anglez = drawingData['Rotation'][:3]
2581        Tx,Ty,Tz = drawingData['viewPoint'][0]
2582        cx,ct,cs = drawingData['atomPtrs']
2583        bondR = drawingData['bondRadius']
2584        G,g = G2lat.cell2Gmat(cell)
2585        GS = G
2586        GS[0][1] = GS[1][0] = math.sqrt(GS[0][0]*GS[1][1])
2587        GS[0][2] = GS[2][0] = math.sqrt(GS[0][0]*GS[2][2])
2588        GS[1][2] = GS[2][1] = math.sqrt(GS[1][1]*GS[2][2])
2589        ellipseProb = G2lat.criticalEllipse(drawingData['ellipseProb']/100.)
2590       
2591        SetBackground()
2592        glInitNames()
2593        glPushName(0)
2594       
2595        glMatrixMode(GL_PROJECTION)
2596        glLoadIdentity()
2597        glViewport(0,0,VS[0],VS[1])
2598        gluPerspective(20.,aspect,cPos-Zclip,cPos+Zclip)
2599        gluLookAt(0,0,cPos,0,0,0,0,1,0)
2600        SetLights()           
2601           
2602        glMatrixMode(GL_MODELVIEW)
2603        glLoadIdentity()
2604        glRotate(anglez,0,0,1)
2605        glRotate(anglex,cosd(anglez),-sind(anglez),0)
2606        glRotate(angley,sind(anglez),cosd(anglez),0)
2607        glMultMatrixf(A4mat.T)
2608        glTranslate(-Tx,-Ty,-Tz)
2609        if drawingData['unitCellBox']:
2610            RenderBox()
2611        if drawingData['showABC']:
2612            x,y,z = drawingData['testPos'][0]
2613#            if altDown:
2614#                self.G2plotNB.status.SetStatusText('moving test point %.4f,%.4f,%.4f'%(x,y,z),1)
2615#            else:
2616#                self.G2plotNB.status.SetStatusText('test point %.4f,%.4f,%.4f'%(x,y,z),1)           
2617            RenderUnitVectors(x,y,z)
2618        Backbone = []
2619        BackboneColor = []
2620        time0 = time.time()
2621        for iat,atom in enumerate(drawingData['Atoms']):
2622            x,y,z = atom[cx:cx+3]
2623            Bonds = atom[-2]
2624            Faces = atom[-1]
2625            try:
2626                atNum = generalData['AtomTypes'].index(atom[ct])
2627            except ValueError:
2628                atNum = -1
2629            CL = atom[cs+2]
2630            color = np.array(CL)/255.
2631            if iat in Ind:
2632                color = np.array(Gr)/255.
2633            radius = 0.5
2634            if atom[cs] != '':
2635                glLoadName(atom[-3])                   
2636            if 'balls' in atom[cs]:
2637                vdwScale = drawingData['vdwScale']
2638                ballScale = drawingData['ballScale']
2639                if atNum < 0:
2640                    radius = 0.2
2641                elif 'H' == atom[ct]:
2642                    if drawingData['showHydrogen']:
2643                        if 'vdW' in atom[cs] and atNum >= 0:
2644                            radius = vdwScale*generalData['vdWRadii'][atNum]
2645                        else:
2646                            radius = ballScale*drawingData['sizeH']
2647                    else:
2648                        radius = 0.0
2649                else:
2650                    if 'vdW' in atom[cs]:
2651                        radius = vdwScale*generalData['vdWRadii'][atNum]
2652                    else:
2653                        radius = ballScale*generalData['BondRadii'][atNum]
2654                RenderSphere(x,y,z,radius,color)
2655                if 'sticks' in atom[cs]:
2656                    RenderBonds(x,y,z,Bonds,bondR,color)
2657            elif 'ellipsoids' in atom[cs]:
2658                RenderBonds(x,y,z,Bonds,bondR,color)
2659                if atom[cs+3] == 'A':                   
2660                    Uij = atom[cs+5:cs+11]
2661                    U = np.multiply(G2spc.Uij2U(Uij),GS)
2662                    U = np.inner(Amat,np.inner(U,Amat).T)
2663                    E,R = nl.eigh(U)
2664                    R4 = np.concatenate((np.concatenate((R,[[0],[0],[0]]),axis=1),[[0,0,0,1],]),axis=0)
2665                    E = np.sqrt(E)
2666                    if atom[ct] == 'H' and not drawingData['showHydrogen']:
2667                        pass
2668                    else:
2669                        RenderEllipsoid(x,y,z,ellipseProb,E,R4,color)                   
2670                else:
2671                    if atom[ct] == 'H' and not drawingData['showHydrogen']:
2672                        pass
2673                    else:
2674                        radius = ellipseProb*math.sqrt(abs(atom[cs+4]))
2675                        RenderSphere(x,y,z,radius,color)
2676            elif 'lines' in atom[cs]:
2677                radius = 0.1
2678                RenderLines(x,y,z,Bonds,color)
2679#                RenderBonds(x,y,z,Bonds,0.05,color,6)
2680            elif atom[cs] == 'sticks':
2681                radius = 0.1
2682                RenderBonds(x,y,z,Bonds,bondR,color)
2683            elif atom[cs] == 'polyhedra':
2684                RenderPolyhedra(x,y,z,Faces,color)
2685            elif atom[cs] == 'backbone':
2686                if atom[ct-1].split()[0] in ['C','N']:
2687                    Backbone.append(list(np.inner(Amat,np.array([x,y,z]))))
2688                    BackboneColor.append(list(color))
2689                   
2690            if atom[cs+1] == 'type':
2691                RenderLabel(x,y,z,atom[ct],radius)
2692            elif atom[cs+1] == 'name':
2693                RenderLabel(x,y,z,atom[ct-1],radius)
2694            elif atom[cs+1] == 'number':
2695                RenderLabel(x,y,z,str(iat+1),radius)
2696            elif atom[cs+1] == 'residue' and atom[ct-1] == 'CA':
2697                RenderLabel(x,y,z,atom[ct-4],radius)
2698            elif atom[cs+1] == '1-letter' and atom[ct-1] == 'CA':
2699                RenderLabel(x,y,z,atom[ct-3],radius)
2700            elif atom[cs+1] == 'chain' and atom[ct-1] == 'CA':
2701                RenderLabel(x,y,z,atom[ct-2],radius)
2702        if Backbone:
2703            RenderBackbone(Backbone,BackboneColor,bondR)
2704#        print time.time()-time0
2705        Page.canvas.SwapBuffers()
2706       
2707    def OnSize(event):
2708        Draw()
2709       
2710    def OnFocus(event):
2711        Draw()
2712       
2713    try:
2714        plotNum = self.G2plotNB.plotList.index(generalData['Name'])
2715        Page = self.G2plotNB.nb.GetPage(plotNum)       
2716    except ValueError:
2717        Plot = self.G2plotNB.addOgl(generalData['Name'])
2718        plotNum = self.G2plotNB.plotList.index(generalData['Name'])
2719        Page = self.G2plotNB.nb.GetPage(plotNum)
2720        Page.views = False
2721        view = False
2722        altDown = False
2723    Page.SetFocus()
2724    cb = wx.ComboBox(self.G2plotNB.status,style=wx.CB_DROPDOWN|wx.CB_READONLY,
2725        choices=(' save as:','jpeg','tiff','bmp'))
2726    cb.Bind(wx.EVT_COMBOBOX, OnKeyBox)
2727    cb.SetValue(' save as:')
2728    Page.canvas.Bind(wx.EVT_MOUSEWHEEL, OnMouseWheel)
2729    Page.canvas.Bind(wx.EVT_LEFT_DOWN, OnMouseDown)
2730    Page.canvas.Bind(wx.EVT_RIGHT_DOWN, OnMouseDown)
2731    Page.canvas.Bind(wx.EVT_MIDDLE_DOWN, OnMouseDown)
2732    Page.canvas.Bind(wx.EVT_KEY_UP, OnKey)
2733    Page.canvas.Bind(wx.EVT_MOTION, OnMouseMove)
2734    Page.canvas.Bind(wx.EVT_SIZE, OnSize)
2735    Page.canvas.Bind(wx.EVT_SET_FOCUS, OnFocus)
2736    Page.camera['position'] = drawingData['cameraPos']
2737    Page.camera['viewPoint'] = np.inner(Amat,drawingData['viewPoint'][0])
2738    Page.camera['backColor'] = np.array(list(drawingData['backColor'])+[0,])/255.
2739    Page.canvas.SetCurrent()
2740    Draw()
2741       
Note: See TracBrowser for help on using the repository browser.