source: trunk/GSASIIplot.py @ 432

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

remove test prints from GSASIIIO
new Gmat2AB routine
work on the DData GUI bug - doesn't crash but bad GUI shown
work on ellipsoidal crystallites - all now work.

  • Property svn:keywords set to Date Author Revision URL Id
File size: 115.3 KB
Line 
1#GSASII plotting routines
2########### SVN repository information ###################
3# $Date: 2011-12-05 22:00:00 +0000 (Mon, 05 Dec 2011) $
4# $Author: vondreele $
5# $Revision: 432 $
6# $URL: trunk/GSASIIplot.py $
7# $Id: GSASIIplot.py 432 2011-12-05 22:00:00Z vondreele $
8########### SVN repository information ###################
9import math
10import time
11import copy
12import os.path
13import numpy as np
14import numpy.linalg as nl
15import wx
16import wx.aui
17import wx.glcanvas
18import matplotlib as mpl
19import mpl_toolkits.mplot3d.axes3d as mp3d
20import GSASIIpath
21import GSASIIgrid as G2gd
22import GSASIIimage as G2img
23import GSASIIpwd as G2pwd
24import GSASIIIO as G2IO
25import GSASIIpwdGUI as G2pdG
26import GSASIIimgGUI as G2imG
27import GSASIIphsGUI as G2phG
28import GSASIIlattice as G2lat
29import GSASIIspc as G2spc
30import pytexture as ptx
31from  OpenGL.GL import *
32from OpenGL.GLU import *
33from OpenGL.GLUT import *
34from OpenGL.GLE import *
35from matplotlib.backends.backend_wx import _load_bitmap
36from matplotlib.backends.backend_wxagg import FigureCanvasWxAgg as Canvas
37from matplotlib.backends.backend_wxagg import NavigationToolbar2Wx as Toolbar
38
39# useful degree trig functions
40sind = lambda x: math.sin(x*math.pi/180.)
41cosd = lambda x: math.cos(x*math.pi/180.)
42tand = lambda x: math.tan(x*math.pi/180.)
43asind = lambda x: 180.*math.asin(x)/math.pi
44acosd = lambda x: 180.*math.acos(x)/math.pi
45atan2d = lambda x,y: 180.*math.atan2(y,x)/math.pi
46atand = lambda x: 180.*math.atan(x)/math.pi
47# numpy versions
48npsind = lambda x: np.sin(x*np.pi/180.)
49npcosd = lambda x: np.cos(x*np.pi/180.)
50npacosd = lambda x: 180.*np.arccos(x)/np.pi
51npasind = lambda x: 180.*np.arcsin(x)/np.pi
52npatand = lambda x: 180.*np.arctan(x)/np.pi
53npatan2d = lambda x,y: 180.*np.arctan2(x,y)/np.pi
54   
55class G2PlotMpl(wx.Panel):   
56    def __init__(self,parent,id=-1,dpi=None,**kwargs):
57        wx.Panel.__init__(self,parent,id=id,**kwargs)
58        mpl.rcParams['legend.fontsize'] = 10
59        self.figure = mpl.figure.Figure(dpi=dpi,figsize=(5,7))
60        self.canvas = Canvas(self,-1,self.figure)
61        self.toolbar = GSASIItoolbar(self.canvas)
62
63        self.toolbar.Realize()
64       
65        sizer=wx.BoxSizer(wx.VERTICAL)
66        sizer.Add(self.canvas,1,wx.EXPAND)
67        sizer.Add(self.toolbar,0,wx.LEFT|wx.EXPAND)
68        self.SetSizer(sizer)
69       
70class G2PlotOgl(wx.Panel):
71    def __init__(self,parent,id=-1,dpi=None,**kwargs):
72        self.figure = wx.Panel.__init__(self,parent,id=id,**kwargs)
73        self.canvas = wx.glcanvas.GLCanvas(self,-1,**kwargs)
74        self.camera = {}
75        sizer=wx.BoxSizer(wx.VERTICAL)
76        sizer.Add(self.canvas,1,wx.EXPAND)
77        self.SetSizer(sizer)
78       
79class G2Plot3D(wx.Panel):
80    def __init__(self,parent,id=-1,dpi=None,**kwargs):
81        wx.Panel.__init__(self,parent,id=id,**kwargs)
82        self.figure = mpl.figure.Figure(dpi=dpi,figsize=(6,6))
83        self.canvas = Canvas(self,-1,self.figure)
84        self.toolbar = GSASIItoolbar(self.canvas)
85
86        self.toolbar.Realize()
87       
88        sizer=wx.BoxSizer(wx.VERTICAL)
89        sizer.Add(self.canvas,1,wx.EXPAND)
90        sizer.Add(self.toolbar,0,wx.LEFT|wx.EXPAND)
91        self.SetSizer(sizer)
92                             
93class G2PlotNoteBook(wx.Panel):
94    def __init__(self,parent,id=-1):
95        wx.Panel.__init__(self,parent,id=id)
96        #so one can't delete a plot page!!
97        self.nb = wx.aui.AuiNotebook(self, \
98            style=wx.aui.AUI_NB_DEFAULT_STYLE ^ wx.aui.AUI_NB_CLOSE_ON_ACTIVE_TAB)
99        sizer = wx.BoxSizer()
100        sizer.Add(self.nb,1,wx.EXPAND)
101        self.SetSizer(sizer)
102        self.status = parent.CreateStatusBar()
103        self.status.SetFieldsCount(2)
104        self.status.SetStatusWidths([150,-1])
105        self.Bind(wx.aui.EVT_AUINOTEBOOK_PAGE_CHANGED, self.OnPageChanged)
106       
107        self.plotList = []
108           
109    def addMpl(self,name=""):
110        page = G2PlotMpl(self.nb)
111        self.nb.AddPage(page,name)
112       
113        self.plotList.append(name)
114       
115        return page.figure
116       
117    def add3D(self,name=""):
118        page = G2Plot3D(self.nb)
119        self.nb.AddPage(page,name)
120       
121        self.plotList.append(name)
122       
123        return page.figure
124       
125    def addOgl(self,name=""):
126        page = G2PlotOgl(self.nb)
127        self.nb.AddPage(page,name)
128       
129        self.plotList.append(name)
130       
131        return page.figure
132       
133    def Delete(self,name):
134        try:
135            item = self.plotList.index(name)
136            del self.plotList[item]
137            self.nb.DeletePage(item)
138        except ValueError:          #no plot of this name - do nothing
139            return     
140               
141    def clear(self):
142        while self.nb.GetPageCount():
143            self.nb.DeletePage(0)
144        self.plotList = []
145        self.status.DestroyChildren()
146       
147    def Rename(self,oldName,newName):
148        try:
149            item = self.plotList.index(oldName)
150            self.plotList[item] = newName
151            self.nb.SetPageText(item,newName)
152        except ValueError:          #no plot of this name - do nothing
153            return     
154       
155    def OnPageChanged(self,event):       
156        if self.plotList:
157            self.status.SetStatusText('Better to select this from GSAS-II data tree',1)
158        self.status.DestroyChildren()                           #get rid of special stuff on status bar
159       
160class GSASIItoolbar(Toolbar):
161    ON_MPL_HELP = wx.NewId()
162    def __init__(self,plotCanvas):
163        Toolbar.__init__(self,plotCanvas)
164        POSITION_OF_CONFIGURE_SUBPLOTS_BTN = 6
165        self.DeleteToolByPos(POSITION_OF_CONFIGURE_SUBPLOTS_BTN)
166        help = os.path.join(os.path.split(sys.argv[0])[0],'help.ico')
167        self.AddSimpleTool(self.ON_MPL_HELP,_load_bitmap(help),'Help on','Show help on')
168        wx.EVT_TOOL(self,self.ON_MPL_HELP,self.OnHelp)
169    def OnHelp(self,event):
170        Page = self.GetParent().GetParent()
171        pageNo = Page.GetSelection()
172       
173        print 'Matplotlib help on '+Page.GetPageText(pageNo)
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            coeff = useList[item][plotDict[plotType]]
1204            if plotType in ['Mustrain','Size']:
1205                if coeff[0] == 'isotropic':
1206                    X *= coeff[1][0]
1207                    Y *= coeff[1][0]
1208                    Z *= coeff[1][0]                               
1209                elif coeff[0] == 'uniaxial':
1210                   
1211                    def uniaxCalc(xyz,iso,aniso,axes):
1212                        Z = np.array(axes)
1213                        cp = abs(np.dot(xyz,Z))
1214                        sp = np.sqrt(1.-cp**2)
1215                        R = iso*aniso/np.sqrt((iso*cp)**2+(aniso*sp)**2)
1216                        return R*xyz
1217                       
1218                    iso,aniso = coeff[1][:2]
1219                    axes = np.inner(A,np.array(coeff[3]))
1220                    axes /= nl.norm(axes)
1221                    Shkl = np.array(coeff[1])
1222                    XYZ = np.dstack((X,Y,Z))
1223                    XYZ = np.nan_to_num(np.apply_along_axis(uniaxCalc,2,XYZ,iso,aniso,axes))
1224                    X,Y,Z = np.dsplit(XYZ,3)
1225                    X = X[:,:,0]
1226                    Y = Y[:,:,0]
1227                    Z = Z[:,:,0]
1228               
1229                elif coeff[0] == 'ellipsoidal':
1230                   
1231                    def ellipseCalc(xyz,E,R):
1232                        XYZ = xyz*E.T
1233                        return np.inner(XYZ.T,R)
1234                       
1235                    S6 = coeff[4]
1236                    Sij = G2lat.U6toUij(S6)
1237                    E,R = nl.eigh(Sij)
1238                    XYZ = np.dstack((X,Y,Z))
1239                    XYZ = np.nan_to_num(np.apply_along_axis(ellipseCalc,2,XYZ,E,R))
1240                    X,Y,Z = np.dsplit(XYZ,3)
1241                    X = X[:,:,0]
1242                    Y = Y[:,:,0]
1243                    Z = Z[:,:,0]
1244                   
1245                elif coeff[0] == 'generalized':
1246                   
1247                    def genMustrain(xyz,SGData,A,Shkl):
1248                        uvw = np.inner(A.T,xyz)
1249                        Strm = np.array(G2spc.MustrainCoeff(uvw,SGData))
1250                        sum = np.sum(np.multiply(Shkl,Strm))
1251                        sum = np.where(sum > 0.01,sum,0.01)
1252                        sum = np.sqrt(sum)*math.pi/0.018      #centidegrees to radians!
1253                        return sum*xyz
1254                       
1255                    Shkl = np.array(coeff[4])
1256                    if np.any(Shkl):
1257                        XYZ = np.dstack((X,Y,Z))
1258                        XYZ = np.nan_to_num(np.apply_along_axis(genMustrain,2,XYZ,SGData,A,Shkl))
1259                        X,Y,Z = np.dsplit(XYZ,3)
1260                        X = X[:,:,0]
1261                        Y = Y[:,:,0]
1262                        Z = Z[:,:,0]
1263                           
1264                if np.any(X) and np.any(Y) and np.any(Z):
1265                    errFlags = np.seterr(all='ignore')
1266                    Plot.plot_surface(X,Y,Z,rstride=1,cstride=1,color='g',linewidth=1)
1267                    np.seterr(all='ignore')
1268                    xyzlim = np.array([Plot.get_xlim3d(),Plot.get_ylim3d(),Plot.get_zlim3d()]).T
1269                    XYZlim = [min(xyzlim[0]),max(xyzlim[1])]
1270                    Plot.set_xlim3d(XYZlim)
1271                    Plot.set_ylim3d(XYZlim)
1272                    Plot.set_zlim3d(XYZlim)
1273                    Plot.set_aspect('equal')
1274                if plotType == 'Size':
1275                    Plot.set_title('Crystallite size for '+phase+'\n'+coeff[0]+' model')
1276                    Plot.set_xlabel(r'X, $\mu$m')
1277                    Plot.set_ylabel(r'Y, $\mu$m')
1278                    Plot.set_zlabel(r'Z, $\mu$m')
1279                else:   
1280                    Plot.set_title(r'$\mu$strain for '+phase+'\n'+coeff[0]+' model')
1281                    Plot.set_xlabel(r'X, $\mu$strain')
1282                    Plot.set_ylabel(r'Y, $\mu$strain')
1283                    Plot.set_zlabel(r'Z, $\mu$strain')
1284            else:
1285                h,k,l = generalData['POhkl']
1286                if coeff[0] == 'MD':
1287                    print 'March-Dollase preferred orientation plot'
1288               
1289                else:
1290                    PH = np.array(generalData['POhkl'])
1291                    phi,beta = G2lat.CrsAng(PH,cell[:6],SGData)
1292                    SHCoef = {}
1293                    for item in coeff[5]:
1294                        L,N = eval(item.strip('C'))
1295                        SHCoef['C%d,0,%d'%(L,N)] = coeff[5][item]                       
1296                    ODFln = G2lat.Flnh(Start,SHCoef,phi,beta,SGData)
1297                    X = np.linspace(0,90.0,26)
1298                    Y = G2lat.polfcal(ODFln,'0',X,0.0)
1299                    Plot.plot(X,Y,color='k',label=str(PH))
1300                    Plot.legend(loc='best')
1301                    Plot.set_title('Axial distribution for HKL='+str(PH)+' in '+phase)
1302                    Plot.set_xlabel(r'$\psi$',fontsize=16)
1303                    Plot.set_ylabel('MRD',fontsize=14)
1304    Page.canvas.draw()
1305   
1306def PlotTexture(self,data,newPlot=False,Start=False):
1307    '''Pole figure, inverse pole figure, 3D pole distribution and 3D inverse pole distribution
1308    plotting.
1309    dict generalData contains all phase info needed which is in data
1310    '''
1311
1312    shModels = ['cylindrical','none','shear - 2/m','rolling - mmm']
1313    SamSym = dict(zip(shModels,['0','-1','2/m','mmm']))
1314    PatternId = self.PatternId
1315    generalData = data['General']
1316    SGData = generalData['SGData']
1317    textureData = generalData['SH Texture']
1318    self.G2plotNB.Delete('Texture')
1319    if not textureData['Order']:
1320        return                  #no plot!!
1321    SHData = generalData['SH Texture']
1322    SHCoef = SHData['SH Coeff'][1]
1323    cell = generalData['Cell'][1:7]
1324    Amat,Bmat = G2lat.cell2AB(cell)
1325    sq2 = 1.0/math.sqrt(2.0)
1326   
1327    def rp2xyz(r,p):
1328        z = npcosd(r)
1329        xy = np.sqrt(1.-z**2)
1330        return xy*npsind(p),xy*npcosd(p),z
1331           
1332    def OnMotion(event):
1333        SHData = data['General']['SH Texture']
1334        if event.xdata and event.ydata:                 #avoid out of frame errors
1335            xpos = event.xdata
1336            ypos = event.ydata
1337            if 'Inverse' in SHData['PlotType']:
1338                r = xpos**2+ypos**2
1339                if r <= 1.0:
1340                    if 'equal' in self.Projection: 
1341                        r,p = 2.*npasind(np.sqrt(r)*sq2),npatan2d(ypos,xpos)
1342                    else:
1343                        r,p = 2.*npatand(np.sqrt(r)),npatan2d(ypos,xpos)
1344                    ipf = G2lat.invpolfcal(IODFln,SGData,np.array([r,]),np.array([p,]))
1345                    xyz = np.inner(Amat.T,np.array([rp2xyz(r,p)]))
1346                    y,x,z = list(xyz/np.max(np.abs(xyz)))
1347                   
1348                    self.G2plotNB.status.SetFields(['',
1349                        'psi =%9.3f, beta =%9.3f, MRD =%9.3f hkl=%5.2f,%5.2f,%5.2f'%(r,p,ipf,x,y,z)])
1350                                   
1351            elif 'Axial' in SHData['PlotType']:
1352                pass
1353               
1354            else:                       #ordinary pole figure
1355                z = xpos**2+ypos**2
1356                if z <= 1.0:
1357                    z = np.sqrt(z)
1358                    if 'equal' in self.Projection: 
1359                        r,p = 2.*npasind(z*sq2),npatan2d(ypos,xpos)
1360                    else:
1361                        r,p = 2.*npatand(z),npatan2d(ypos,xpos)
1362                    pf = G2lat.polfcal(ODFln,SamSym[textureData['Model']],np.array([r,]),np.array([p,]))
1363                    self.G2plotNB.status.SetFields(['','phi =%9.3f, gam =%9.3f, MRD =%9.3f'%(r,p,pf)])
1364   
1365    if self.Projection == '3D display' and 'Axial'  not in SHData['PlotType']:               
1366        Plot = mp3d.Axes3D(self.G2plotNB.add3D('Texture'))
1367    else:
1368        Plot = self.G2plotNB.addMpl('Texture').gca()
1369    plotNum = self.G2plotNB.plotList.index('Texture')
1370    Page = self.G2plotNB.nb.GetPage(plotNum)
1371    if self.Projection == '3D display':
1372        pass
1373    else:
1374        Page.canvas.mpl_connect('motion_notify_event', OnMotion)
1375
1376    Page.SetFocus()
1377    self.G2plotNB.status.SetFields(['',''])   
1378    PH = np.array(SHData['PFhkl'])
1379    phi,beta = G2lat.CrsAng(PH,cell,SGData)
1380    ODFln = G2lat.Flnh(Start,SHCoef,phi,beta,SGData)
1381    PX = np.array(SHData['PFxyz'])
1382    gam = atan2d(PX[0],PX[1])
1383    xy = math.sqrt(PX[0]**2+PX[1]**2)
1384    xyz = math.sqrt(PX[0]**2+PX[1]**2+PX[2]**2)
1385    psi = asind(xy/xyz)
1386    IODFln = G2lat.Glnh(Start,SHCoef,psi,gam,SamSym[textureData['Model']])
1387    if 'Axial' in SHData['PlotType']:
1388        X = np.linspace(0,90.0,26)
1389        Y = G2lat.polfcal(ODFln,SamSym[textureData['Model']],X,0.0)
1390        Plot.plot(X,Y,color='k',label=str(SHData['PFhkl']))
1391        Plot.legend(loc='best')
1392        Plot.set_title('Axial distribution for HKL='+str(SHData['PFhkl']))
1393        Plot.set_xlabel(r'$\psi$',fontsize=16)
1394        Plot.set_ylabel('MRD',fontsize=14)
1395       
1396    else:       
1397        npts = 201
1398        if 'Inverse' in SHData['PlotType']:
1399            X,Y = np.meshgrid(np.linspace(1.,-1.,npts),np.linspace(-1.,1.,npts))
1400            R,P = np.sqrt(X**2+Y**2).flatten(),npatan2d(X,Y).flatten()
1401            if 'equal' in self.Projection:
1402                R = np.where(R <= 1.,2.*npasind(R*sq2),0.0)
1403            else:
1404                R = np.where(R <= 1.,2.*npatand(R),0.0)
1405            Z = np.zeros_like(R)
1406            Z = G2lat.invpolfcal(IODFln,SGData,R,P)
1407            Z = np.reshape(Z,(npts,npts))
1408            CS = Plot.contour(Y,X,Z,aspect='equal')
1409            Plot.clabel(CS,fontsize=9,inline=1)
1410            try:
1411                Img = Plot.imshow(Z.T,aspect='equal',cmap=self.ContourColor,extent=[-1,1,-1,1])
1412            except ValueError:
1413                pass
1414            if newPlot:
1415#                Page.figure.colorbar(Img)    #colorbar fails - crashes gsasii
1416                newPlot = False
1417            Plot.set_title('Inverse pole figure for XYZ='+str(SHData['PFxyz']))
1418            Plot.set_xlabel(self.Projection.capitalize()+' projection')
1419                       
1420        else:
1421            X,Y = np.meshgrid(np.linspace(1.,-1.,npts),np.linspace(-1.,1.,npts))
1422            R,P = np.sqrt(X**2+Y**2).flatten(),npatan2d(X,Y).flatten()
1423            if 'equal' in self.Projection:
1424                R = np.where(R <= 1.,2.*npasind(R*sq2),0.0)
1425            else:
1426                R = np.where(R <= 1.,2.*npatand(R),0.0)
1427            Z = np.zeros_like(R)
1428            Z = G2lat.polfcal(ODFln,SamSym[textureData['Model']],R,P)
1429            Z = np.reshape(Z,(npts,npts))
1430            try:
1431                CS = Plot.contour(Y,X,Z,aspect='equal')
1432                Plot.clabel(CS,fontsize=9,inline=1)
1433            except ValueError:
1434                pass
1435            Img = Plot.imshow(Z.T,aspect='equal',cmap=self.ContourColor,extent=[-1,1,-1,1])
1436            if newPlot:
1437#                Page.figure.colorbar(Img)    #colorbar fails - crashes gsasii
1438                newPlot = False
1439            Plot.set_title('Pole figure for HKL='+str(SHData['PFhkl']))
1440            Plot.set_xlabel(self.Projection.capitalize()+' projection')
1441    Page.canvas.draw()
1442
1443def PlotCovariance(self,Data={}):
1444    if not Data:
1445        Data = self.PatternTree.GetItemPyData(
1446            G2gd.GetPatternTreeItemId(self,self.root, 'Covariance'))
1447    if not Data:
1448        print 'No covariance matrix available'
1449        return
1450    varyList = Data['varyList']
1451    values = Data['variables']
1452    Xmax = len(varyList)
1453    covMatrix = Data['covMatrix']
1454    sig = np.sqrt(np.diag(covMatrix))
1455    xvar = np.outer(sig,np.ones_like(sig))
1456    covArray = np.divide(np.divide(covMatrix,xvar),xvar.T)
1457    title = ' for\n'+Data['title']
1458    newAtomDict = Data['newAtomDict']
1459
1460    def OnPlotKeyPress(event):
1461        newPlot = False
1462        if event.key == 's':
1463            choice = [m for m in mpl.cm.datad.keys() if not m.endswith("_r")]
1464            choice.sort()
1465            dlg = wx.SingleChoiceDialog(self,'Select','Color scheme',choice)
1466            if dlg.ShowModal() == wx.ID_OK:
1467                sel = dlg.GetSelection()
1468                self.VcovColor = choice[sel]
1469            else:
1470                self.VcovColor = 'RdYlGn'
1471            dlg.Destroy()
1472        PlotCovariance(self)
1473
1474    def OnMotion(event):
1475        if event.button:
1476            ytics = imgAx.get_yticks()
1477            ytics = np.where(ytics<len(varyList),ytics,-1)
1478            ylabs = [np.where(0<=i ,varyList[int(i)],' ') for i in ytics]
1479            imgAx.set_yticklabels(ylabs)
1480           
1481        if event.xdata and event.ydata:                 #avoid out of frame errors
1482            xpos = int(event.xdata+.5)
1483            ypos = int(event.ydata+.5)
1484            if -1 < xpos < len(varyList) and -1 < ypos < len(varyList):
1485                if xpos == ypos:
1486                    value = values[xpos]
1487                    name = varyList[xpos]
1488                    if varyList[xpos] in newAtomDict:
1489                        name,value = newAtomDict[name]                       
1490                    msg = '%s value = %.4g, esd = %.4g'%(name,value,sig[xpos])
1491                else:
1492                    msg = '%s - %s: %5.3f'%(varyList[xpos],varyList[ypos],covArray[xpos][ypos])
1493                Page.canvas.SetToolTipString(msg)
1494                self.G2plotNB.status.SetFields(['Key: s to change colors',msg])
1495    try:
1496        plotNum = self.G2plotNB.plotList.index('Covariance')
1497        Page = self.G2plotNB.nb.GetPage(plotNum)
1498        Page.figure.clf()
1499        Plot = Page.figure.gca()
1500        if not Page.IsShown():
1501            Page.Show()
1502    except ValueError:
1503        Plot = self.G2plotNB.addMpl('Covariance').gca()
1504        plotNum = self.G2plotNB.plotList.index('Covariance')
1505        Page = self.G2plotNB.nb.GetPage(plotNum)
1506        Page.canvas.mpl_connect('motion_notify_event', OnMotion)
1507        Page.canvas.mpl_connect('key_press_event', OnPlotKeyPress)
1508
1509    Page.SetFocus()
1510    self.G2plotNB.status.SetFields(['',''])   
1511    acolor = mpl.cm.get_cmap(self.VcovColor)
1512    Img = Plot.imshow(covArray,aspect='equal',cmap=acolor,interpolation='nearest',origin='lower')
1513    imgAx = Img.get_axes()
1514    ytics = imgAx.get_yticks()
1515    ylabs = [varyList[int(i)] for i in ytics[:-1]]
1516    imgAx.set_yticklabels(ylabs)
1517    colorBar = Page.figure.colorbar(Img)
1518    Plot.set_title('V-Cov matrix'+title)
1519    Plot.set_xlabel('Variable number')
1520    Plot.set_ylabel('Variable name')
1521    Page.canvas.draw()
1522   
1523def PlotSeq(self,SeqData,SeqSig,SeqNames,sampleParm):
1524   
1525    def OnKeyPress(event):
1526        if event.key == 's' and sampleParm:
1527            if self.xAxis:
1528                self.xAxis = False
1529            else:
1530                self.xAxis = True
1531            Draw(False)
1532    try:
1533        plotNum = self.G2plotNB.plotList.index('Sequential refinement')
1534        Page = self.G2plotNB.nb.GetPage(plotNum)
1535        Page.figure.clf()
1536        Plot = Page.figure.gca()
1537        if not Page.IsShown():
1538            Page.Show()
1539    except ValueError:
1540        Plot = self.G2plotNB.addMpl('Sequential refinement').gca()
1541        plotNum = self.G2plotNB.plotList.index('Sequential refinement')
1542        Page = self.G2plotNB.nb.GetPage(plotNum)
1543        Page.canvas.mpl_connect('key_press_event', OnKeyPress)
1544        self.xAxis = False
1545       
1546    def Draw(newPlot):
1547        Page.SetFocus()
1548        self.G2plotNB.status.SetFields(['','press s to toggle x-axis = sample environment parameter'])
1549        if len(SeqData):
1550            Plot.clear()
1551            if self.xAxis:   
1552                xName = sampleParm.keys()[0]
1553                X = sampleParm[xName]
1554            else:
1555                X = np.arange(0,len(SeqData[0]),1)
1556                xName = 'Data sequence number'
1557            for Y,sig,name in zip(SeqData,SeqSig,SeqNames):
1558                Plot.errorbar(X,Y,yerr=sig,label=name)       
1559            Plot.legend(loc='best')
1560            Plot.set_ylabel('Parameter values')
1561            Plot.set_xlabel(xName)
1562            Page.canvas.draw()           
1563    Draw(True)
1564           
1565def PlotExposedImage(self,newPlot=False,event=None):
1566    '''General access module for 2D image plotting
1567    '''
1568    plotNo = self.G2plotNB.nb.GetSelection()
1569    if self.G2plotNB.nb.GetPageText(plotNo) == '2D Powder Image':
1570        PlotImage(self,newPlot,event,newImage=True)
1571    elif self.G2plotNB.nb.GetPageText(plotNo) == '2D Integration':
1572        PlotIntegration(self,newPlot,event)
1573
1574def PlotImage(self,newPlot=False,event=None,newImage=True):
1575    '''Plot of 2D detector images as contoured plot. Also plot calibration ellipses,
1576    masks, etc.
1577    '''
1578    from matplotlib.patches import Ellipse,Arc,Circle,Polygon
1579    import numpy.ma as ma
1580    Dsp = lambda tth,wave: wave/(2.*sind(tth/2.))
1581    global Data,Masks
1582    colors=['b','g','r','c','m','k']
1583    Data = self.PatternTree.GetItemPyData(
1584        G2gd.GetPatternTreeItemId(self,self.Image, 'Image Controls'))
1585    Masks = self.PatternTree.GetItemPyData(
1586        G2gd.GetPatternTreeItemId(self,self.Image, 'Masks'))
1587
1588    def OnImMotion(event):
1589        Page.canvas.SetToolTipString('')
1590        sizexy = Data['size']
1591        if event.xdata and event.ydata:                 #avoid out of frame errors
1592            Page.canvas.SetCursor(wx.CROSS_CURSOR)
1593            item = self.itemPicked
1594            pixelSize = Data['pixelSize']
1595            scalex = 1000./pixelSize[0]
1596            scaley = 1000./pixelSize[1]
1597            if item and self.PatternTree.GetItemText(self.PickId) == 'Image Controls':
1598                if 'Text' in str(item):
1599                    Page.canvas.SetToolTipString('%8.3f %8.3fmm'%(event.xdata,event.ydata))
1600                else:
1601                    xcent,ycent = Data['center']
1602                    xpos = event.xdata-xcent
1603                    ypos = event.ydata-ycent
1604                    tth,azm = G2img.GetTthAzm(event.xdata,event.ydata,Data)
1605                    if 'line3' in  str(item) or 'line4' in str(item) and not Data['fullIntegrate']:
1606                        Page.canvas.SetToolTipString('%6d deg'%(azm))
1607                    elif 'line1' in  str(item) or 'line2' in str(item):
1608                        Page.canvas.SetToolTipString('%8.3fdeg'%(tth))                           
1609            else:
1610                xpos = event.xdata
1611                ypos = event.ydata
1612                xpix = xpos*scalex
1613                ypix = ypos*scaley
1614                Int = 0
1615                if (0 <= xpix <= sizexy[0]) and (0 <= ypix <= sizexy[1]):
1616                    Int = self.ImageZ[ypix][xpix]
1617                tth,azm,dsp = G2img.GetTthAzmDsp(xpos,ypos,Data)
1618                Q = 2.*math.pi/dsp
1619                if self.setPoly:
1620                    self.G2plotNB.status.SetFields(['','Polygon mask pick - LB next point, RB close polygon'])
1621                else:
1622                    self.G2plotNB.status.SetFields(\
1623                        ['','Detector 2-th =%9.3fdeg, dsp =%9.3fA, Q = %6.5fA-1, azm = %7.2fdeg, I = %6d'%(tth,dsp,Q,azm,Int)])
1624
1625    def OnImPlotKeyPress(event):
1626        try:
1627            PickName = self.PatternTree.GetItemText(self.PickId)
1628        except TypeError:
1629            return
1630        if PickName == 'Masks':
1631            Xpos = event.xdata
1632            if not Xpos:            #got point out of frame
1633                return
1634            Ypos = event.ydata
1635            if event.key == 's':
1636                Masks['Points'].append([Xpos,Ypos,1.])
1637            elif event.key == 'r':
1638                tth = G2img.GetTth(Xpos,Ypos,Data)
1639                Masks['Rings'].append([tth,0.1])
1640            elif event.key == 'a':
1641                tth,azm = G2img.GetTthAzm(Xpos,Ypos,Data)
1642                azm = int(azm)               
1643                Masks['Arcs'].append([tth,[azm-5,azm+5],0.1])
1644            elif event.key == 'p':
1645                self.setPoly = True
1646                Masks['Polygons'].append([])
1647                self.G2plotNB.status.SetFields(['','Polygon mask active - LB pick next point, RB close polygon'])
1648            G2imG.UpdateMasks(self,Masks)
1649        elif PickName == 'Image Controls':
1650            if event.key == 'c':
1651                Xpos = event.xdata
1652                if not Xpos:            #got point out of frame
1653                    return
1654                Ypos = event.ydata
1655                dlg = wx.MessageDialog(self,'Are you sure you want to change the center?',
1656                    'Center change',style=wx.OK|wx.CANCEL)
1657                try:
1658                    if dlg.ShowModal() == wx.ID_OK:
1659                        print 'move center to: ',Xpos,Ypos
1660                        Data['center'] = [Xpos,Ypos]
1661                        G2imG.UpdateImageControls(self,Data,Masks)
1662                finally:
1663                    dlg.Destroy()
1664            elif event.key == 'l':
1665                if self.logPlot:
1666                    self.logPlot = False
1667                else:
1668                    self.logPlot = True
1669        PlotImage(self,newImage=True)
1670           
1671    def OnKeyBox(event):
1672        if self.G2plotNB.nb.GetSelection() == self.G2plotNB.plotList.index('2D Powder Image'):
1673            event.key = cb.GetValue()[0]
1674            cb.SetValue(' key press')
1675            if event.key in 'l':
1676                OnImPlotKeyPress(event)
1677                       
1678    def OnImPick(event):
1679        if self.PatternTree.GetItemText(self.PickId) not in ['Image Controls','Masks']:
1680            return
1681        if self.setPoly:
1682            polygon = Masks['Polygons'][-1]
1683            xpos,ypos = event.mouseevent.xdata,event.mouseevent.ydata
1684            if xpos and ypos:                       #point inside image
1685                if len(polygon) > 2 and event.mouseevent.button == 3:
1686                    x0,y0 = polygon[0]
1687                    polygon.append([x0,y0])
1688                    self.setPoly = False
1689                    self.G2plotNB.status.SetFields(['','Polygon closed - RB drag a vertex to change shape'])
1690                else:
1691                    self.G2plotNB.status.SetFields(['','New polygon point: %.1f,%.1f'%(xpos,ypos)])
1692                    polygon.append([xpos,ypos])
1693                G2imG.UpdateMasks(self,Masks)
1694        else:
1695            if self.itemPicked is not None: return
1696            self.itemPicked = event.artist
1697            self.mousePicked = event.mouseevent
1698       
1699    def OnImRelease(event):
1700        try:
1701            PickName = self.PatternTree.GetItemText(self.PickId)
1702        except TypeError:
1703            return
1704        if PickName not in ['Image Controls','Masks']:
1705            return
1706        pixelSize = Data['pixelSize']
1707        scalex = 1000./pixelSize[0]
1708        scaley = 1000./pixelSize[1]
1709        pixLimit = Data['pixLimit']
1710        if self.itemPicked is None and PickName == 'Image Controls':
1711#            sizexy = Data['size']
1712            Xpos = event.xdata
1713            if not (Xpos and self.ifGetRing):                   #got point out of frame
1714                return
1715            Ypos = event.ydata
1716            if Ypos and not Page.toolbar._active:         #make sure zoom/pan not selected
1717                if event.button == 1:
1718                    Xpix = Xpos*scalex
1719                    Ypix = Ypos*scaley
1720                    xpos,ypos,I,J = G2img.ImageLocalMax(self.ImageZ,pixLimit,Xpix,Ypix)
1721                    if I and J:
1722                        xpos += .5                              #shift to pixel center
1723                        ypos += .5
1724                        xpos /= scalex                          #convert to mm
1725                        ypos /= scaley
1726                        Data['ring'].append([xpos,ypos])
1727                elif event.button == 3:
1728                    self.dataFrame.GetStatusBar().SetStatusText('Calibrating...')
1729                    if G2img.ImageCalibrate(self,Data):
1730                        self.dataFrame.GetStatusBar().SetStatusText('Calibration successful - Show ring picks to check')
1731                        print 'Calibration successful'
1732                    else:
1733                        self.dataFrame.GetStatusBar().SetStatusText('Calibration failed - Show ring picks to diagnose')
1734                        print 'Calibration failed'
1735                    self.ifGetRing = False
1736                    G2imG.UpdateImageControls(self,Data,Masks)
1737                    return
1738                PlotImage(self,newImage=False)
1739            return
1740        else:
1741            xpos = event.xdata
1742            if xpos:                                        #avoid out of frame mouse position
1743                ypos = event.ydata
1744                if self.ifGetRing:                          #delete a calibration ring pick
1745                    xypos = [xpos,ypos]
1746                    rings = Data['ring']
1747                    for ring in rings:
1748                        if np.allclose(ring,xypos,.01,0):
1749                            rings.remove(ring)                                                                       
1750                else:
1751                    tth,azm,dsp = G2img.GetTthAzmDsp(xpos,ypos,Data)
1752                    itemPicked = str(self.itemPicked)
1753                    if 'Line2D' in itemPicked and PickName == 'Image Controls':
1754                        if 'line1' in itemPicked:
1755                            Data['IOtth'][0] = max(tth,0.001)
1756                        elif 'line2' in itemPicked:
1757                            Data['IOtth'][1] = tth
1758                        elif 'line3' in itemPicked:
1759                            Data['LRazimuth'][0] = int(azm)
1760                            if Data['fullIntegrate']:
1761                                Data['LRazimuth'][1] = Data['LRazimuth'][0]+360
1762                        elif 'line4' in itemPicked and not Data['fullIntegrate']:
1763                            Data['LRazimuth'][1] = int(azm)
1764                           
1765                        if Data['LRazimuth'][0] > Data['LRazimuth'][1]:
1766                            Data['LRazimuth'][0] -= 360
1767                           
1768                        azLim = np.array(Data['LRazimuth'])   
1769                        if np.any(azLim>360):
1770                            azLim -= 360
1771                            Data['LRazimuth'] = list(azLim)
1772                           
1773                        if  Data['IOtth'][0] > Data['IOtth'][1]:
1774                            Data['IOtth'][0],Data['IOtth'][1] = Data['IOtth'][1],Data['IOtth'][0]
1775                           
1776                        self.InnerTth.SetValue("%8.2f" % (Data['IOtth'][0]))
1777                        self.OuterTth.SetValue("%8.2f" % (Data['IOtth'][1]))
1778                        self.Lazim.SetValue("%6d" % (Data['LRazimuth'][0]))
1779                        self.Razim.SetValue("%6d" % (Data['LRazimuth'][1]))
1780                    elif 'Circle' in itemPicked and PickName == 'Masks':
1781                        spots = Masks['Points']
1782                        newPos = itemPicked.split(')')[0].split('(')[2].split(',')
1783                        newPos = np.array([float(newPos[0]),float(newPos[1])])
1784                        for spot in spots:
1785                            if np.allclose(np.array([spot[:2]]),newPos):
1786                                spot[:2] = xpos,ypos
1787                        G2imG.UpdateMasks(self,Masks)
1788                    elif 'Line2D' in itemPicked and PickName == 'Masks':
1789                        Obj = self.itemPicked.findobj()
1790                        rings = Masks['Rings']
1791                        arcs = Masks['Arcs']
1792                        polygons = Masks['Polygons']
1793                        for ring in self.ringList:
1794                            if Obj == ring[0]:
1795                                rN = ring[1]
1796                                if ring[2] == 'o':
1797                                    rings[rN][0] = G2img.GetTth(xpos,ypos,Data)-rings[rN][1]/2.
1798                                else:
1799                                    rings[rN][0] = G2img.GetTth(xpos,ypos,Data)+rings[rN][1]/2.
1800                        for arc in self.arcList:
1801                            if Obj == arc[0]:
1802                                aN = arc[1]
1803                                if arc[2] == 'o':
1804                                    arcs[aN][0] = G2img.GetTth(xpos,ypos,Data)-arcs[aN][2]/2
1805                                elif arc[2] == 'i':
1806                                    arcs[aN][0] = G2img.GetTth(xpos,ypos,Data)+arcs[aN][2]/2
1807                                elif arc[2] == 'l':
1808                                    arcs[aN][1][0] = int(G2img.GetAzm(xpos,ypos,Data))
1809                                else:
1810                                    arcs[aN][1][1] = int(G2img.GetAzm(xpos,ypos,Data))
1811                        for poly in self.polyList:
1812                            if Obj == poly[0]:
1813                                ind = self.itemPicked.contains(self.mousePicked)[1]['ind'][0]
1814                                oldPos = np.array([self.mousePicked.xdata,self.mousePicked.ydata])
1815                                pN = poly[1]
1816                                for i,xy in enumerate(polygons[pN]):
1817                                    if np.allclose(np.array([xy]),oldPos,atol=1.0):
1818                                        polygons[pN][i] = xpos,ypos
1819                        G2imG.UpdateMasks(self,Masks)
1820#                    else:                  #keep for future debugging
1821#                        print str(self.itemPicked),event.xdata,event.ydata,event.button
1822                PlotImage(self,newImage=True)
1823            self.itemPicked = None
1824           
1825    try:
1826        plotNum = self.G2plotNB.plotList.index('2D Powder Image')
1827        Page = self.G2plotNB.nb.GetPage(plotNum)
1828        if not newPlot:
1829            Plot = Page.figure.gca()          #get previous powder plot & get limits
1830            xylim = Plot.get_xlim(),Plot.get_ylim()
1831        if newImage:
1832            Page.figure.clf()
1833            Plot = Page.figure.gca()          #get a fresh plot after clf()
1834       
1835    except ValueError:
1836        Plot = self.G2plotNB.addMpl('2D Powder Image').gca()
1837        plotNum = self.G2plotNB.plotList.index('2D Powder Image')
1838        Page = self.G2plotNB.nb.GetPage(plotNum)
1839        Page.canvas.mpl_connect('key_press_event', OnImPlotKeyPress)
1840        Page.canvas.mpl_connect('motion_notify_event', OnImMotion)
1841        Page.canvas.mpl_connect('pick_event', OnImPick)
1842        Page.canvas.mpl_connect('button_release_event', OnImRelease)
1843        xylim = []
1844    if not event:                       #event from GUI TextCtrl - don't want focus to change to plot!!!
1845        Page.SetFocus()
1846    Title = self.PatternTree.GetItemText(self.Image)[4:]
1847    self.G2plotNB.status.DestroyChildren()
1848    if self.logPlot:
1849        Title = 'log('+Title+')'
1850    Plot.set_title(Title)
1851    try:
1852        if self.PatternTree.GetItemText(self.PickId) in ['Image Controls',]:
1853            Choice = (' key press','l: log(I) on',)
1854            if self.logPlot:
1855                Choice[1] = 'l: log(I) off'
1856            cb = wx.ComboBox(self.G2plotNB.status,style=wx.CB_DROPDOWN|wx.CB_READONLY,
1857                choices=Choice)
1858            cb.Bind(wx.EVT_COMBOBOX, OnKeyBox)
1859            cb.SetValue(' key press')
1860    except TypeError:
1861        pass
1862    size,imagefile = self.PatternTree.GetItemPyData(self.Image)
1863    if imagefile != self.oldImagefile:
1864        imagefile = G2IO.CheckImageFile(self,imagefile)
1865        if not imagefile:
1866            self.G2plotNB.Delete('2D Powder Image')
1867            return
1868        self.PatternTree.SetItemPyData(self.Image,[size,imagefile])
1869        self.ImageZ = G2IO.GetImageData(self,imagefile,imageOnly=True)
1870#        print self.ImageZ.shape,self.ImageZ.size,Data['size'] #might be useful debugging line
1871        self.oldImagefile = imagefile
1872
1873    imScale = 1
1874    if len(self.ImageZ) > 1024:
1875        imScale = len(self.ImageZ)/1024
1876    sizexy = Data['size']
1877    pixelSize = Data['pixelSize']
1878    scalex = 1000./pixelSize[0]
1879    scaley = 1000./pixelSize[1]
1880    Xmax = sizexy[0]*pixelSize[0]/1000.
1881    Ymax = sizexy[1]*pixelSize[1]/1000.
1882    xlim = (0,Xmax)
1883    ylim = (Ymax,0)
1884    Imin,Imax = Data['range'][1]
1885    acolor = mpl.cm.get_cmap(Data['color'])
1886    xcent,ycent = Data['center']
1887    Plot.set_xlabel('Image x-axis, mm',fontsize=12)
1888    Plot.set_ylabel('Image y-axis, mm',fontsize=12)
1889    #do threshold mask - "real" mask - others are just bondaries
1890    Zlim = Masks['Thresholds'][1]
1891    wx.BeginBusyCursor()
1892    try:
1893           
1894        if newImage:                   
1895            MA = ma.masked_greater(ma.masked_less(self.ImageZ,Zlim[0]),Zlim[1])
1896            MaskA = ma.getmaskarray(MA)
1897            A = G2img.ImageCompress(MA,imScale)
1898            AM = G2img.ImageCompress(MaskA,imScale)
1899            if self.logPlot:
1900                A = np.where(A>0,np.log(A),0)
1901                AM = np.where(AM>0,np.log(AM),0)
1902                Imin,Imax = [np.amin(A),np.amax(A)]
1903            ImgM = Plot.imshow(AM,aspect='equal',cmap='Reds',
1904                interpolation='nearest',vmin=0,vmax=2,extent=[0,Xmax,Xmax,0])
1905            Img = Plot.imshow(A,aspect='equal',cmap=acolor,
1906                interpolation='nearest',vmin=Imin,vmax=Imax,extent=[0,Xmax,Ymax,0])
1907            if self.setPoly:
1908                Img.set_picker(True)
1909   
1910        Plot.plot(xcent,ycent,'x')
1911        if Data['showLines']:
1912            LRAzim = Data['LRazimuth']                  #NB: integers
1913            Nazm = Data['outAzimuths']
1914            delAzm = float(LRAzim[1]-LRAzim[0])/Nazm
1915            AzmthOff = Data['azmthOff']
1916            IOtth = Data['IOtth']
1917            wave = Data['wavelength']
1918            dspI = wave/(2.0*sind(IOtth[0]/2.0))
1919            ellI = G2img.GetEllipse(dspI,Data)           #=False if dsp didn't yield an ellipse (ugh! a parabola or a hyperbola)
1920            dspO = wave/(2.0*sind(IOtth[1]/2.0))
1921            ellO = G2img.GetEllipse(dspO,Data)           #Ditto & more likely for outer ellipse
1922            Azm = np.array(range(LRAzim[0],LRAzim[1]+1))-AzmthOff
1923            if ellI:
1924                xyI = []
1925                for azm in Azm:
1926                    xyI.append(G2img.GetDetectorXY(dspI,azm-90.,Data))
1927                xyI = np.array(xyI)
1928                arcxI,arcyI = xyI.T
1929                Plot.plot(arcxI,arcyI,picker=3)
1930            if ellO:
1931                xyO = []
1932                for azm in Azm:
1933                    xyO.append(G2img.GetDetectorXY(dspO,azm-90.,Data))
1934                xyO = np.array(xyO)
1935                arcxO,arcyO = xyO.T
1936                Plot.plot(arcxO,arcyO,picker=3)
1937            if ellO and ellI:
1938                Plot.plot([arcxI[0],arcxO[0]],[arcyI[0],arcyO[0]],picker=3)
1939                Plot.plot([arcxI[-1],arcxO[-1]],[arcyI[-1],arcyO[-1]],picker=3)
1940            for i in range(Nazm):
1941                cake = LRAzim[0]+i*delAzm-AzmthOff
1942                ind = np.searchsorted(Azm,cake)
1943                Plot.plot([arcxI[ind],arcxO[ind]],[arcyI[ind],arcyO[ind]],color='k',dashes=(5,5))
1944                   
1945        for xring,yring in Data['ring']:
1946            Plot.plot(xring,yring,'r+',picker=3)
1947        if Data['setRings']:
1948#            rings = np.concatenate((Data['rings']),axis=0)
1949            N = 0
1950            for ring in Data['rings']:
1951                xring,yring = np.array(ring).T[:2]
1952                Plot.plot(xring,yring,'+',color=colors[N%6])
1953                N += 1           
1954        for ellipse in Data['ellipses']:
1955            cent,phi,[width,height],col = ellipse
1956            Plot.add_artist(Ellipse([cent[0],cent[1]],2*width,2*height,phi,ec=col,fc='none'))
1957            Plot.text(cent[0],cent[1],'+',color=col,ha='center',va='center')
1958        #masks - mask lines numbered after integration limit lines
1959        spots = Masks['Points']
1960        rings = Masks['Rings']
1961        arcs = Masks['Arcs']
1962        polygons = Masks['Polygons']
1963        for x,y,d in spots:
1964            Plot.add_artist(Circle((x,y),radius=d/2,fc='none',ec='r',picker=3))
1965        self.ringList = []
1966        for iring,(tth,thick) in enumerate(rings):
1967            wave = Data['wavelength']
1968            x1,y1 = np.hsplit(np.array(G2img.makeIdealRing(G2img.GetEllipse(Dsp(tth+thick/2.,wave),Data))),2)
1969            x2,y2 = np.hsplit(np.array(G2img.makeIdealRing(G2img.GetEllipse(Dsp(tth-thick/2.,wave),Data))),2)
1970            self.ringList.append([Plot.plot(x1,y1,'r',picker=3),iring,'o'])           
1971            self.ringList.append([Plot.plot(x2,y2,'r',picker=3),iring,'i'])
1972        self.arcList = []
1973        for iarc,(tth,azm,thick) in enumerate(arcs):           
1974            wave = Data['wavelength']
1975            x1,y1 = np.hsplit(np.array(G2img.makeIdealRing(G2img.GetEllipse(Dsp(tth+thick/2.,wave),Data),azm)),2)
1976            x2,y2 = np.hsplit(np.array(G2img.makeIdealRing(G2img.GetEllipse(Dsp(max(0.01,tth-thick/2.),wave),Data),azm)),2)
1977            self.arcList.append([Plot.plot(x1,y1,'r',picker=3),iarc,'o'])           
1978            self.arcList.append([Plot.plot(x2,y2,'r',picker=3),iarc,'i'])
1979            self.arcList.append([Plot.plot([x1[0],x2[0]],[y1[0],y2[0]],'r',picker=3),iarc,'l'])
1980            self.arcList.append([Plot.plot([x1[-1],x2[-1]],[y1[-1],y2[-1]],'r',picker=3),iarc,'u'])
1981        self.polyList = []
1982        for ipoly,polygon in enumerate(polygons):
1983            x,y = np.hsplit(np.array(polygon),2)
1984            self.polyList.append([Plot.plot(x,y,'r+',picker=10),ipoly])
1985            Plot.plot(x,y,'r')           
1986        if newImage:
1987            colorBar = Page.figure.colorbar(Img)
1988        Plot.set_xlim(xlim)
1989        Plot.set_ylim(ylim)
1990        if not newPlot and xylim:
1991            Page.toolbar.push_current()
1992            Plot.set_xlim(xylim[0])
1993            Plot.set_ylim(xylim[1])
1994            xylim = []
1995            Page.toolbar.push_current()
1996            Page.toolbar.draw()
1997        else:
1998            Page.canvas.draw()
1999    finally:
2000        wx.EndBusyCursor()
2001       
2002def PlotIntegration(self,newPlot=False,event=None):
2003    '''Plot of 2D image after image integration with 2-theta and azimuth as coordinates
2004    '''
2005           
2006    def OnMotion(event):
2007        Page.canvas.SetToolTipString('')
2008        Page.canvas.SetCursor(wx.CROSS_CURSOR)
2009        azm = event.ydata
2010        tth = event.xdata
2011        if azm and tth:
2012            self.G2plotNB.status.SetFields(\
2013                ['','Detector 2-th =%9.3fdeg, azm = %7.2fdeg'%(tth,azm)])
2014                               
2015    try:
2016        plotNum = self.G2plotNB.plotList.index('2D Integration')
2017        Page = self.G2plotNB.nb.GetPage(plotNum)
2018        if not newPlot:
2019            Plot = Page.figure.gca()          #get previous plot & get limits
2020            xylim = Plot.get_xlim(),Plot.get_ylim()
2021        Page.figure.clf()
2022        Plot = Page.figure.gca()          #get a fresh plot after clf()
2023       
2024    except ValueError:
2025        Plot = self.G2plotNB.addMpl('2D Integration').gca()
2026        plotNum = self.G2plotNB.plotList.index('2D Integration')
2027        Page = self.G2plotNB.nb.GetPage(plotNum)
2028        Page.canvas.mpl_connect('motion_notify_event', OnMotion)
2029        Page.views = False
2030        view = False
2031    if not event:
2032        Page.SetFocus()
2033       
2034    Data = self.PatternTree.GetItemPyData(
2035        G2gd.GetPatternTreeItemId(self,self.Image, 'Image Controls'))
2036    image = self.Integrate[0]
2037    xsc = self.Integrate[1]
2038    ysc = self.Integrate[2]
2039    Imin,Imax = Data['range'][1]
2040    acolor = mpl.cm.get_cmap(Data['color'])
2041    Plot.set_title(self.PatternTree.GetItemText(self.Image)[4:])
2042    Plot.set_ylabel('azimuth',fontsize=12)
2043    Plot.set_xlabel('2-theta',fontsize=12)
2044    Img = Plot.imshow(image,cmap=acolor,vmin=Imin,vmax=Imax,interpolation='nearest', \
2045        extent=[ysc[0],ysc[-1],xsc[-1],xsc[0]],aspect='auto')
2046    colorBar = Page.figure.colorbar(Img)
2047    if Data['ellipses']:           
2048        for ellipse in Data['ellipses']:
2049            ring = np.array(G2img.makeIdealRing(ellipse[:3])) #skip color
2050            x,y = np.hsplit(ring,2)
2051            tth,azm = G2img.GetTthAzm(x,y,Data)
2052#            azm = np.where(azm < 0.,azm+360,azm)
2053            Plot.plot(tth,azm,'b,')
2054    if not newPlot:
2055        Page.toolbar.push_current()
2056        Plot.set_xlim(xylim[0])
2057        Plot.set_ylim(xylim[1])
2058        xylim = []
2059        Page.toolbar.push_current()
2060        Page.toolbar.draw()
2061    else:
2062        Page.canvas.draw()
2063               
2064def PlotTRImage(self,tax,tay,taz,newPlot=False):
2065    '''a test plot routine - not normally used
2066    ''' 
2067           
2068    def OnMotion(event):
2069        Page.canvas.SetToolTipString('')
2070        Page.canvas.SetCursor(wx.CROSS_CURSOR)
2071        azm = event.xdata
2072        tth = event.ydata
2073        if azm and tth:
2074            self.G2plotNB.status.SetFields(\
2075                ['','Detector 2-th =%9.3fdeg, azm = %7.2fdeg'%(tth,azm)])
2076                               
2077    try:
2078        plotNum = self.G2plotNB.plotList.index('2D Transformed Powder Image')
2079        Page = self.G2plotNB.nb.GetPage(plotNum)
2080        if not newPlot:
2081            Plot = Page.figure.gca()          #get previous plot & get limits
2082            xylim = Plot.get_xlim(),Plot.get_ylim()
2083        Page.figure.clf()
2084        Plot = Page.figure.gca()          #get a fresh plot after clf()
2085       
2086    except ValueError:
2087        Plot = self.G2plotNB.addMpl('2D Transformed Powder Image').gca()
2088        plotNum = self.G2plotNB.plotList.index('2D Transformed Powder Image')
2089        Page = self.G2plotNB.nb.GetPage(plotNum)
2090        Page.canvas.mpl_connect('motion_notify_event', OnMotion)
2091        Page.views = False
2092        view = False
2093    Page.SetFocus()
2094       
2095    Data = self.PatternTree.GetItemPyData(
2096        G2gd.GetPatternTreeItemId(self,self.Image, 'Image Controls'))
2097    Imin,Imax = Data['range'][1]
2098    step = (Imax-Imin)/5.
2099    V = np.arange(Imin,Imax,step)
2100    acolor = mpl.cm.get_cmap(Data['color'])
2101    Plot.set_title(self.PatternTree.GetItemText(self.Image)[4:])
2102    Plot.set_xlabel('azimuth',fontsize=12)
2103    Plot.set_ylabel('2-theta',fontsize=12)
2104    Plot.contour(tax,tay,taz,V,cmap=acolor)
2105    if Data['showLines']:
2106        IOtth = Data['IOtth']
2107        if Data['fullIntegrate']:
2108            LRAzim = [-180,180]
2109        else:
2110            LRAzim = Data['LRazimuth']                  #NB: integers
2111        Plot.plot([LRAzim[0],LRAzim[1]],[IOtth[0],IOtth[0]],picker=True)
2112        Plot.plot([LRAzim[0],LRAzim[1]],[IOtth[1],IOtth[1]],picker=True)
2113        if not Data['fullIntegrate']:
2114            Plot.plot([LRAzim[0],LRAzim[0]],[IOtth[0],IOtth[1]],picker=True)
2115            Plot.plot([LRAzim[1],LRAzim[1]],[IOtth[0],IOtth[1]],picker=True)
2116    if Data['setRings']:
2117        rings = np.concatenate((Data['rings']),axis=0)
2118        for xring,yring,dsp in rings:
2119            x,y = G2img.GetTthAzm(xring,yring,Data)
2120            Plot.plot(y,x,'r+')           
2121    if Data['ellipses']:           
2122        for ellipse in Data['ellipses']:
2123            ring = np.array(G2img.makeIdealRing(ellipse[:3])) #skip color
2124            x,y = np.hsplit(ring,2)
2125            tth,azm = G2img.GetTthAzm(x,y,Data)
2126            Plot.plot(azm,tth,'b,')
2127    if not newPlot:
2128        Page.toolbar.push_current()
2129        Plot.set_xlim(xylim[0])
2130        Plot.set_ylim(xylim[1])
2131        xylim = []
2132        Page.toolbar.push_current()
2133        Page.toolbar.draw()
2134    else:
2135        Page.canvas.draw()
2136       
2137def PlotStructure(self,data):
2138    '''Crystal structure plotting package. Can show structures as balls, sticks, lines,
2139    thermal motion ellipsoids and polyhedra
2140    '''
2141    generalData = data['General']
2142    cell = generalData['Cell'][1:7]
2143    Amat,Bmat = G2lat.cell2AB(cell)         #Amat - crystal to cartesian, Bmat - inverse
2144    A4mat = np.concatenate((np.concatenate((Amat,[[0],[0],[0]]),axis=1),[[0,0,0,1],]),axis=0)
2145    B4mat = np.concatenate((np.concatenate((Bmat,[[0],[0],[0]]),axis=1),[[0,0,0,1],]),axis=0)
2146    Mydir = generalData['Mydir']
2147    atomData = data['Atoms']
2148    drawingData = data['Drawing']
2149    drawAtoms = drawingData['Atoms']
2150    cx,ct,cs = drawingData['atomPtrs']
2151    Wt = [255,255,255]
2152    Rd = [255,0,0]
2153    Gr = [0,255,0]
2154    Bl = [0,0,255]
2155    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]])
2156    uEdges = np.array([
2157        [uBox[0],uBox[1]],[uBox[0],uBox[3]],[uBox[0],uBox[4]],[uBox[1],uBox[2]], 
2158        [uBox[2],uBox[3]],[uBox[1],uBox[5]],[uBox[2],uBox[6]],[uBox[3],uBox[7]], 
2159        [uBox[4],uBox[5]],[uBox[5],uBox[6]],[uBox[6],uBox[7]],[uBox[7],uBox[4]]])
2160    uColors = [Rd,Gr,Bl,Wt, Wt,Wt,Wt,Wt, Wt,Wt,Wt,Wt]
2161    altDown = False
2162    shiftDown = False
2163    ctrlDown = False
2164   
2165    def OnKeyBox(event):
2166        import Image
2167        Draw()                          #make sure plot is fresh!!
2168        mode = cb.GetValue()
2169        Fname = Mydir+'\\'+generalData['Name']+'.'+mode
2170        size = Page.canvas.GetSize()
2171        glPixelStorei(GL_UNPACK_ALIGNMENT, 1)
2172        if mode in ['jpeg',]:
2173            Pix = glReadPixels(0,0,size[0],size[1],GL_RGBA, GL_UNSIGNED_BYTE)
2174            im = Image.new("RGBA", (size[0],size[1]))
2175        else:
2176            Pix = glReadPixels(0,0,size[0],size[1],GL_RGB, GL_UNSIGNED_BYTE)
2177            im = Image.new("RGB", (size[0],size[1]))
2178        im.fromstring(Pix)
2179        im.save(Fname,mode)
2180        cb.SetValue(' Save as:')
2181        self.G2plotNB.status.SetStatusText('Drawing saved to: '+Fname,1)
2182   
2183    def GetTruePosition(xy):
2184        View = glGetIntegerv(GL_VIEWPORT)
2185        Proj = glGetDoublev(GL_PROJECTION_MATRIX)
2186        Model = glGetDoublev(GL_MODELVIEW_MATRIX)
2187        Zmax = 1.
2188        for i,atom in enumerate(drawAtoms):
2189            x,y,z = atom[cx:cx+3]
2190            X,Y,Z = gluProject(x,y,z,Model,Proj,View)
2191            XY = [int(X),int(View[3]-Y)]
2192            if np.allclose(xy,XY,atol=10) and Z < Zmax:
2193                Zmax = Z
2194                SetSelectedAtoms(i)
2195                   
2196    def OnMouseDown(event):
2197        xy = event.GetPosition()
2198        if event.ShiftDown():
2199            GetTruePosition(xy)
2200        else:
2201            drawingData['Rotation'][3] = xy
2202        Draw()
2203       
2204    def OnMouseMove(event):
2205        newxy = event.GetPosition()
2206        page = getSelection()
2207        if event.ControlDown() and drawingData['showABC']:
2208            if event.LeftIsDown():
2209                SetTestRot(newxy)
2210            elif event.RightIsDown():
2211                SetTestPos(newxy)
2212            elif event.MiddleIsDown():
2213                SetTestRotZ(newxy)
2214            x,y,z = drawingData['testPos'][0]
2215            self.G2plotNB.status.SetStatusText('moving test point %.4f,%.4f,%.4f'%(x,y,z),1)
2216               
2217               
2218        if event.Dragging() and not event.ControlDown():
2219            if event.LeftIsDown():
2220                SetRotation(newxy)
2221                angX,angY,angZ = drawingData['Rotation'][:3]
2222                self.G2plotNB.status.SetStatusText('New rotation: %.2f, %.2f ,%.2f'%(angX,angY,angZ),1)
2223            elif event.RightIsDown():
2224                SetTranslation(newxy)
2225                Tx,Ty,Tz = drawingData['viewPoint'][0]
2226                self.G2plotNB.status.SetStatusText('New view point: %.4f, %.4f, %.4f'%(Tx,Ty,Tz),1)
2227            elif event.MiddleIsDown():
2228                SetRotationZ(newxy)
2229                angX,angY,angZ = drawingData['Rotation'][:3]
2230                self.G2plotNB.status.SetStatusText('New rotation: %.2f, %.2f, %.2f'%(angX,angY,angZ),1)
2231        Draw()
2232       
2233    def OnMouseWheel(event):
2234        drawingData['cameraPos'] += event.GetWheelRotation()/24
2235        drawingData['cameraPos'] = max(10,min(500,drawingData['cameraPos']))
2236        self.G2plotNB.status.SetStatusText('New camera distance: %.2f'%(drawingData['cameraPos']),1)
2237        page = getSelection()
2238        if page:
2239            if self.dataDisplay.GetPageText(page) == 'Draw Options':
2240                panel = self.dataDisplay.GetPage(page).GetChildren()[0].GetChildren()
2241                names = [child.GetName() for child in panel]
2242                panel[names.index('cameraPos')].SetLabel('Camera Position: '+'%.2f'%(drawingData['cameraPos']))
2243                panel[names.index('cameraSlider')].SetValue(drawingData['cameraPos'])
2244        Draw()
2245       
2246    def getSelection():
2247        try:
2248            return self.dataDisplay.GetSelection()
2249        except AttributeError:
2250            print self.dataDisplay.GetLabel()
2251            self.G2plotNB.status.SetStatusText('Select this from Phase data window!')
2252            return 0
2253           
2254    def SetViewPointText(VP):
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('viewPoint')].SetValue('%.3f, %.3f, %.3f'%(VP[0],VP[1],VP[2]))
2261           
2262    def ClearSelectedAtoms():
2263        page = getSelection()
2264        if page:
2265            if self.dataDisplay.GetPageText(page) == 'Draw Atoms':
2266                self.dataDisplay.GetPage(page).ClearSelection()      #this is the Atoms grid in Draw Atoms
2267            elif self.dataDisplay.GetPageText(page) == 'Atoms':
2268                self.dataDisplay.GetPage(page).ClearSelection()      #this is the Atoms grid in Atoms
2269                   
2270    def SetSelectedAtoms(ind):
2271        page = getSelection()
2272        if page:
2273            if self.dataDisplay.GetPageText(page) == 'Draw Atoms':
2274                self.dataDisplay.GetPage(page).SelectRow(ind)      #this is the Atoms grid in Draw Atoms
2275            elif self.dataDisplay.GetPageText(page) == 'Atoms':
2276                Id = drawAtoms[ind][-2]
2277                for i,atom in enumerate(atomData):
2278                    if atom[-1] == Id:
2279                        self.dataDisplay.GetPage(page).SelectRow(i)      #this is the Atoms grid in Atoms
2280                 
2281    def GetSelectedAtoms():
2282        page = getSelection()
2283        Ind = []
2284        if page:
2285            if self.dataDisplay.GetPageText(page) == 'Draw Atoms':
2286                Ind = self.dataDisplay.GetPage(page).GetSelectedRows()      #this is the Atoms grid in Draw Atoms
2287            elif self.dataDisplay.GetPageText(page) == 'Atoms':
2288                Ind = self.dataDisplay.GetPage(page).GetSelectedRows()      #this is the Atoms grid in Atoms
2289        return Ind
2290                                       
2291    def OnKey(event):           #on key UP!!
2292        keyCode = event.GetKeyCode()
2293        if keyCode > 255:
2294            keyCode = 0
2295        key,xyz = chr(keyCode),event.GetPosition()
2296        indx = drawingData['selectedAtoms']
2297        if key in ['c','C']:
2298            drawingData['viewPoint'] = [[.5,.5,.5],[0,0]]
2299            drawingData['testPos'] = [[-.1,-.1,-.1],[0.0,0.0,0.0],[0,0]]
2300            drawingData['Rotation'] = [0.0,0.0,0.0,[]]
2301            SetViewPointText(drawingData['viewPoint'][0])
2302        elif key in ['n','N']:
2303            drawAtoms = drawingData['Atoms']
2304            pI = drawingData['viewPoint'][1]
2305            if indx:
2306                pI[0] = indx[pI[1]]
2307                Tx,Ty,Tz = drawAtoms[pI[0]][cx:cx+3]
2308                pI[1] += 1
2309                if pI[1] >= len(indx):
2310                    pI[1] = 0
2311            else:
2312                Tx,Ty,Tz = drawAtoms[pI[0]][cx:cx+3]               
2313                pI[0] += 1
2314                if pI[0] >= len(drawAtoms):
2315                    pI[0] = 0
2316            drawingData['viewPoint'] = [[Tx,Ty,Tz],pI]
2317            SetViewPointText(drawingData['viewPoint'][0])
2318               
2319        elif key in ['p','P']:
2320            drawAtoms = drawingData['Atoms']
2321            pI = drawingData['viewPoint'][1]
2322            if indx:
2323                pI[0] = indx[pI[1]]
2324                Tx,Ty,Tz = drawAtoms[pI[0]][cx:cx+3]
2325                pI[1] -= 1
2326                if pI[1] < 0:
2327                    pI[1] = len(indx)-1
2328            else:
2329                Tx,Ty,Tz = drawAtoms[pI[0]][cx:cx+3]               
2330                pI[0] -= 1
2331                if pI[0] < 0:
2332                    pI[0] = len(drawAtoms)-1
2333            drawingData['viewPoint'] = [[Tx,Ty,Tz],pI]
2334            SetViewPointText(drawingData['viewPoint'][0])           
2335        Draw()
2336           
2337    def SetBackground():
2338        R,G,B,A = Page.camera['backColor']
2339        glClearColor(R,G,B,A)
2340        glClear(GL_COLOR_BUFFER_BIT | GL_DEPTH_BUFFER_BIT)
2341       
2342    def SetLights():
2343        glEnable(GL_DEPTH_TEST)
2344        glShadeModel(GL_SMOOTH)
2345        glEnable(GL_LIGHTING)
2346        glEnable(GL_LIGHT0)
2347        glLightModeli(GL_LIGHT_MODEL_TWO_SIDE,0)
2348        glLightfv(GL_LIGHT0,GL_AMBIENT,[1,1,1,.8])
2349        glLightfv(GL_LIGHT0,GL_DIFFUSE,[1,1,1,1])
2350       
2351    def SetTranslation(newxy):
2352        Tx,Ty,Tz = drawingData['viewPoint'][0]
2353        anglex,angley,anglez,oldxy = drawingData['Rotation']
2354        Rx = G2lat.rotdMat(anglex,0)
2355        Ry = G2lat.rotdMat(angley,1)
2356        Rz = G2lat.rotdMat(anglez,2)
2357        dxy = list(newxy-oldxy)+[0,]
2358        dxy = np.inner(Bmat,np.inner(Rz,np.inner(Ry,np.inner(Rx,dxy))))
2359        Tx -= dxy[0]*0.01
2360        Ty += dxy[1]*0.01
2361        Tz -= dxy[2]*0.01
2362        drawingData['Rotation'][3] = newxy
2363        drawingData['viewPoint'][0] =  Tx,Ty,Tz
2364        SetViewPointText([Tx,Ty,Tz])
2365       
2366    def SetTestPos(newxy):
2367        Tx,Ty,Tz = drawingData['testPos'][0]
2368        anglex,angley,anglez,oldxy = drawingData['Rotation']
2369        Rx = G2lat.rotdMat(anglex,0)
2370        Ry = G2lat.rotdMat(angley,1)
2371        Rz = G2lat.rotdMat(anglez,2)
2372        dxy = list(newxy-oldxy)+[0,]
2373        dxy = np.inner(Rz,np.inner(Ry,np.inner(Rx,dxy)))
2374        Tx += dxy[0]*0.001
2375        Ty -= dxy[1]*0.001
2376        Tz += dxy[2]*0.001
2377        drawingData['Rotation'][3] = newxy
2378        drawingData['testPos'][0] =  Tx,Ty,Tz
2379       
2380    def SetTestRot(newxy):
2381        Txyz = np.array(drawingData['testPos'][0])
2382        oldxy = drawingData['testPos'][2]
2383        Ax,Ay,Az = drawingData['testPos'][1]
2384        Vxyz = np.array(drawingData['viewPoint'][0])
2385        Dxyz = np.inner(Amat,Txyz-Vxyz)
2386        dxy = list(newxy-oldxy)+[0,]
2387        Ax += dxy[1]*0.01
2388        Ay += dxy[0]*0.01
2389        Rx = G2lat.rotdMat(Ax,0)
2390        Ry = G2lat.rotdMat(Ay,1)
2391        Dxyz = np.inner(Ry,np.inner(Rx,Dxyz))       
2392        Dxyz = np.inner(Bmat,Dxyz)+Vxyz
2393        drawingData['testPos'][1] = [Ax,Ay,Az]
2394        drawingData['testPos'][2] = newxy
2395        drawingData['testPos'][0] = Dxyz
2396       
2397    def SetTestRotZ(newxy):
2398        Txyz = np.array(drawingData['testPos'][0])
2399        oldxy = drawingData['testPos'][2]
2400        Ax,Ay,Az = drawingData['testPos'][1]
2401        Vxyz = np.array(drawingData['viewPoint'][0])
2402        Dxyz = np.inner(Amat,Txyz-Vxyz)       
2403        dxy = list(newxy-oldxy)+[0,]
2404        Az += (dxy[0]+dxy[1])*.01
2405        Rz = G2lat.rotdMat(Az,2)
2406        Dxyz = np.inner(Rz,Dxyz)       
2407        Dxyz = np.inner(Bmat,Dxyz)+Vxyz
2408        drawingData['testPos'][1] = [Ax,Ay,Az]
2409        drawingData['testPos'][2] = newxy
2410        drawingData['testPos'][0] = Dxyz
2411                             
2412    def SetRotation(newxy):       
2413        anglex,angley,anglez,oldxy = drawingData['Rotation']
2414        dxy = newxy-oldxy
2415        anglex += dxy[1]*.25
2416        angley += dxy[0]*.25
2417        oldxy = newxy
2418        drawingData['Rotation'] = [anglex,angley,anglez,oldxy]
2419       
2420    def SetRotationZ(newxy):                       
2421        anglex,angley,anglez,oldxy = drawingData['Rotation']
2422        dxy = newxy-oldxy
2423        anglez += (dxy[0]+dxy[1])*.25
2424        oldxy = newxy
2425        drawingData['Rotation'] = [anglex,angley,anglez,oldxy]
2426       
2427    def RenderBox():
2428        glEnable(GL_COLOR_MATERIAL)
2429        glLineWidth(2)
2430        glEnable(GL_BLEND)
2431        glBlendFunc(GL_SRC_ALPHA,GL_ONE_MINUS_SRC_ALPHA)
2432        glEnable(GL_LINE_SMOOTH)
2433        glBegin(GL_LINES)
2434        for line,color in zip(uEdges,uColors):
2435            glColor3ubv(color)
2436            glVertex3fv(line[0])
2437            glVertex3fv(line[1])
2438        glEnd()
2439        glColor4ubv([0,0,0,0])
2440        glDisable(GL_LINE_SMOOTH)
2441        glDisable(GL_BLEND)
2442        glDisable(GL_COLOR_MATERIAL)
2443       
2444    def RenderUnitVectors(x,y,z):
2445        xyz = np.array([x,y,z])
2446        glEnable(GL_COLOR_MATERIAL)
2447        glLineWidth(1)
2448        glPushMatrix()
2449        glTranslate(x,y,z)
2450        glScalef(1/cell[0],1/cell[1],1/cell[2])
2451        glBegin(GL_LINES)
2452        for line,color in zip(uEdges,uColors)[:3]:
2453            glColor3ubv(color)
2454            glVertex3fv(line[0])
2455            glVertex3fv(line[1])
2456        glEnd()
2457        glPopMatrix()
2458        glColor4ubv([0,0,0,0])
2459        glDisable(GL_COLOR_MATERIAL)
2460               
2461    def RenderSphere(x,y,z,radius,color):
2462        glMaterialfv(GL_FRONT_AND_BACK,GL_DIFFUSE,color)
2463        glPushMatrix()
2464        glTranslate(x,y,z)
2465        glMultMatrixf(B4mat.T)
2466        q = gluNewQuadric()
2467        gluSphere(q,radius,20,10)
2468        glPopMatrix()
2469       
2470    def RenderEllipsoid(x,y,z,ellipseProb,E,R4,color):
2471        s1,s2,s3 = E
2472        glMaterialfv(GL_FRONT_AND_BACK,GL_DIFFUSE,color)
2473        glPushMatrix()
2474        glTranslate(x,y,z)
2475        glMultMatrixf(B4mat.T)
2476        glMultMatrixf(R4.T)
2477        glEnable(GL_NORMALIZE)
2478        glScale(s1,s2,s3)
2479        q = gluNewQuadric()
2480        gluSphere(q,ellipseProb,20,10)
2481        glDisable(GL_NORMALIZE)
2482        glPopMatrix()
2483       
2484    def RenderBonds(x,y,z,Bonds,radius,color,slice=20):
2485        glMaterialfv(GL_FRONT_AND_BACK,GL_DIFFUSE,color)
2486        glPushMatrix()
2487        glTranslate(x,y,z)
2488        glMultMatrixf(B4mat.T)
2489        for bond in Bonds:
2490            glPushMatrix()
2491            Dx = np.inner(Amat,bond)
2492            Z = np.sqrt(np.sum(Dx**2))
2493            azm = atan2d(-Dx[1],-Dx[0])
2494            phi = acosd(Dx[2]/Z)
2495            glRotate(-azm,0,0,1)
2496            glRotate(phi,1,0,0)
2497            q = gluNewQuadric()
2498            gluCylinder(q,radius,radius,Z,slice,2)
2499            glPopMatrix()           
2500        glPopMatrix()
2501               
2502    def RenderLines(x,y,z,Bonds,color):
2503        xyz = np.array([x,y,z])
2504        glEnable(GL_COLOR_MATERIAL)
2505        glLineWidth(1)
2506        glColor3fv(color)
2507        glPushMatrix()
2508        glBegin(GL_LINES)
2509        for bond in Bonds:
2510            glVertex3fv(xyz)
2511            glVertex3fv(xyz+bond)
2512        glEnd()
2513        glColor4ubv([0,0,0,0])
2514        glPopMatrix()
2515        glDisable(GL_COLOR_MATERIAL)
2516       
2517    def RenderPolyhedra(x,y,z,Faces,color):
2518        glPushMatrix()
2519        glTranslate(x,y,z)
2520        glMaterialfv(GL_FRONT_AND_BACK,GL_DIFFUSE,color)
2521        glShadeModel(GL_SMOOTH)
2522        glMultMatrixf(B4mat.T)
2523        for face,norm in Faces:
2524            glPolygonMode(GL_FRONT_AND_BACK,GL_FILL)
2525            glFrontFace(GL_CW)
2526            glNormal3fv(norm)
2527            glBegin(GL_TRIANGLES)
2528            for vert in face:
2529                glVertex3fv(vert)
2530            glEnd()
2531        glPopMatrix()
2532       
2533    def RenderBackbone(Backbone,BackboneColor,radius):
2534        glPushMatrix()
2535        glMultMatrixf(B4mat.T)
2536        glEnable(GL_COLOR_MATERIAL)
2537        glShadeModel(GL_SMOOTH)
2538        gleSetJoinStyle(TUBE_NORM_EDGE | TUBE_JN_ANGLE | TUBE_JN_CAP)
2539        glePolyCylinder(Backbone,BackboneColor,radius)
2540        glPopMatrix()       
2541        glDisable(GL_COLOR_MATERIAL)
2542       
2543    def RenderLabel(x,y,z,label,r):       
2544        glPushMatrix()
2545        glTranslate(x,y,z)
2546        glMultMatrixf(B4mat.T)
2547        glDisable(GL_LIGHTING)
2548        glColor3f(0,1.,0)
2549        glRasterPos3f(r,r,r)
2550        for c in list(label):
2551            glutBitmapCharacter(GLUT_BITMAP_8_BY_13,ord(c))
2552        glEnable(GL_LIGHTING)
2553        glPopMatrix()
2554                           
2555    def Draw():
2556        import numpy.linalg as nl
2557        Ind = GetSelectedAtoms()
2558        VS = np.array(Page.canvas.GetSize())
2559        aspect = float(VS[0])/float(VS[1])
2560        cPos = drawingData['cameraPos']
2561        Zclip = drawingData['Zclip']*cPos/200.
2562        anglex,angley,anglez = drawingData['Rotation'][:3]
2563        Tx,Ty,Tz = drawingData['viewPoint'][0]
2564        cx,ct,cs = drawingData['atomPtrs']
2565        bondR = drawingData['bondRadius']
2566        G,g = G2lat.cell2Gmat(cell)
2567        GS = G
2568        GS[0][1] = GS[1][0] = math.sqrt(GS[0][0]*GS[1][1])
2569        GS[0][2] = GS[2][0] = math.sqrt(GS[0][0]*GS[2][2])
2570        GS[1][2] = GS[2][1] = math.sqrt(GS[1][1]*GS[2][2])
2571        ellipseProb = G2lat.criticalEllipse(drawingData['ellipseProb']/100.)
2572       
2573        SetBackground()
2574        glInitNames()
2575        glPushName(0)
2576       
2577        glMatrixMode(GL_PROJECTION)
2578        glLoadIdentity()
2579        glViewport(0,0,VS[0],VS[1])
2580        gluPerspective(20.,aspect,cPos-Zclip,cPos+Zclip)
2581        gluLookAt(0,0,cPos,0,0,0,0,1,0)
2582        SetLights()           
2583           
2584        glMatrixMode(GL_MODELVIEW)
2585        glLoadIdentity()
2586        glRotate(anglez,0,0,1)
2587        glRotate(anglex,cosd(anglez),-sind(anglez),0)
2588        glRotate(angley,sind(anglez),cosd(anglez),0)
2589        glMultMatrixf(A4mat.T)
2590        glTranslate(-Tx,-Ty,-Tz)
2591        if drawingData['unitCellBox']:
2592            RenderBox()
2593        if drawingData['showABC']:
2594            x,y,z = drawingData['testPos'][0]
2595#            if altDown:
2596#                self.G2plotNB.status.SetStatusText('moving test point %.4f,%.4f,%.4f'%(x,y,z),1)
2597#            else:
2598#                self.G2plotNB.status.SetStatusText('test point %.4f,%.4f,%.4f'%(x,y,z),1)           
2599            RenderUnitVectors(x,y,z)
2600        Backbone = []
2601        BackboneColor = []
2602        time0 = time.time()
2603        for iat,atom in enumerate(drawingData['Atoms']):
2604            x,y,z = atom[cx:cx+3]
2605            Bonds = atom[-2]
2606            Faces = atom[-1]
2607            try:
2608                atNum = generalData['AtomTypes'].index(atom[ct])
2609            except ValueError:
2610                atNum = -1
2611            CL = atom[cs+2]
2612            color = np.array(CL)/255.
2613            if iat in Ind:
2614                color = np.array(Gr)/255.
2615            radius = 0.5
2616            if atom[cs] != '':
2617                glLoadName(atom[-3])                   
2618            if 'balls' in atom[cs]:
2619                vdwScale = drawingData['vdwScale']
2620                ballScale = drawingData['ballScale']
2621                if atNum < 0:
2622                    radius = 0.2
2623                elif 'H' == atom[ct]:
2624                    if drawingData['showHydrogen']:
2625                        if 'vdW' in atom[cs] and atNum >= 0:
2626                            radius = vdwScale*generalData['vdWRadii'][atNum]
2627                        else:
2628                            radius = ballScale*drawingData['sizeH']
2629                    else:
2630                        radius = 0.0
2631                else:
2632                    if 'vdW' in atom[cs]:
2633                        radius = vdwScale*generalData['vdWRadii'][atNum]
2634                    else:
2635                        radius = ballScale*generalData['BondRadii'][atNum]
2636                RenderSphere(x,y,z,radius,color)
2637                if 'sticks' in atom[cs]:
2638                    RenderBonds(x,y,z,Bonds,bondR,color)
2639            elif 'ellipsoids' in atom[cs]:
2640                RenderBonds(x,y,z,Bonds,bondR,color)
2641                if atom[cs+3] == 'A':                   
2642                    Uij = atom[cs+5:cs+11]
2643                    U = np.multiply(G2spc.Uij2U(Uij),GS)
2644                    U = np.inner(Amat,np.inner(U,Amat).T)
2645                    E,R = nl.eigh(U)
2646                    R4 = np.concatenate((np.concatenate((R,[[0],[0],[0]]),axis=1),[[0,0,0,1],]),axis=0)
2647                    E = np.sqrt(E)
2648                    if atom[ct] == 'H' and not drawingData['showHydrogen']:
2649                        pass
2650                    else:
2651                        RenderEllipsoid(x,y,z,ellipseProb,E,R4,color)                   
2652                else:
2653                    if atom[ct] == 'H' and not drawingData['showHydrogen']:
2654                        pass
2655                    else:
2656                        radius = ellipseProb*math.sqrt(abs(atom[cs+4]))
2657                        RenderSphere(x,y,z,radius,color)
2658            elif 'lines' in atom[cs]:
2659                radius = 0.1
2660                RenderLines(x,y,z,Bonds,color)
2661#                RenderBonds(x,y,z,Bonds,0.05,color,6)
2662            elif atom[cs] == 'sticks':
2663                radius = 0.1
2664                RenderBonds(x,y,z,Bonds,bondR,color)
2665            elif atom[cs] == 'polyhedra':
2666                RenderPolyhedra(x,y,z,Faces,color)
2667            elif atom[cs] == 'backbone':
2668                if atom[ct-1].split()[0] in ['C','N']:
2669                    Backbone.append(list(np.inner(Amat,np.array([x,y,z]))))
2670                    BackboneColor.append(list(color))
2671                   
2672            if atom[cs+1] == 'type':
2673                RenderLabel(x,y,z,atom[ct],radius)
2674            elif atom[cs+1] == 'name':
2675                RenderLabel(x,y,z,atom[ct-1],radius)
2676            elif atom[cs+1] == 'number':
2677                RenderLabel(x,y,z,str(iat+1),radius)
2678            elif atom[cs+1] == 'residue' and atom[ct-1] == 'CA':
2679                RenderLabel(x,y,z,atom[ct-4],radius)
2680            elif atom[cs+1] == '1-letter' and atom[ct-1] == 'CA':
2681                RenderLabel(x,y,z,atom[ct-3],radius)
2682            elif atom[cs+1] == 'chain' and atom[ct-1] == 'CA':
2683                RenderLabel(x,y,z,atom[ct-2],radius)
2684        if Backbone:
2685            RenderBackbone(Backbone,BackboneColor,bondR)
2686#        print time.time()-time0
2687        Page.canvas.SwapBuffers()
2688       
2689    def OnSize(event):
2690        Draw()
2691       
2692    def OnFocus(event):
2693        Draw()
2694       
2695    try:
2696        plotNum = self.G2plotNB.plotList.index(generalData['Name'])
2697        Page = self.G2plotNB.nb.GetPage(plotNum)       
2698    except ValueError:
2699        Plot = self.G2plotNB.addOgl(generalData['Name'])
2700        plotNum = self.G2plotNB.plotList.index(generalData['Name'])
2701        Page = self.G2plotNB.nb.GetPage(plotNum)
2702        Page.views = False
2703        view = False
2704        altDown = False
2705    Page.SetFocus()
2706    cb = wx.ComboBox(self.G2plotNB.status,style=wx.CB_DROPDOWN|wx.CB_READONLY,
2707        choices=(' save as:','jpeg','tiff','bmp'))
2708    cb.Bind(wx.EVT_COMBOBOX, OnKeyBox)
2709    cb.SetValue(' save as:')
2710    Page.canvas.Bind(wx.EVT_MOUSEWHEEL, OnMouseWheel)
2711    Page.canvas.Bind(wx.EVT_LEFT_DOWN, OnMouseDown)
2712    Page.canvas.Bind(wx.EVT_RIGHT_DOWN, OnMouseDown)
2713    Page.canvas.Bind(wx.EVT_MIDDLE_DOWN, OnMouseDown)
2714    Page.canvas.Bind(wx.EVT_KEY_UP, OnKey)
2715    Page.canvas.Bind(wx.EVT_MOTION, OnMouseMove)
2716    Page.canvas.Bind(wx.EVT_SIZE, OnSize)
2717    Page.canvas.Bind(wx.EVT_SET_FOCUS, OnFocus)
2718    Page.camera['position'] = drawingData['cameraPos']
2719    Page.camera['viewPoint'] = np.inner(Amat,drawingData['viewPoint'][0])
2720    Page.camera['backColor'] = np.array(list(drawingData['backColor'])+[0,])/255.
2721    Page.canvas.SetCurrent()
2722    Draw()
2723       
Note: See TracBrowser for help on using the repository browser.