source: trunk/GSASIIplot.py @ 431

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

3D plots for mustrain & size now have correct aspect ratio
title includes model name

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