source: trunk/GSASIIplot.py @ 244

Last change on this file since 244 was 244, checked in by vondreele, 11 years ago

fix to reading GEsum data - made "memory error"s
improve read of missing image dialog
GSASIIplot.py now has a commented image size print line

  • Property svn:keywords set to Date Author Revision URL Id
File size: 79.3 KB
Line 
1#GSASII plotting routines
2########### SVN repository information ###################
3# $Date: 2011-02-12 16:26:33 +0000 (Sat, 12 Feb 2011) $
4# $Author: vondreele $
5# $Revision: 244 $
6# $URL: trunk/GSASIIplot.py $
7# $Id: GSASIIplot.py 244 2011-02-12 16:26:33Z 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 GSASIIIO as G2IO
24import GSASIIpwdGUI as G2pdG
25import GSASIIimgGUI as G2imG
26import GSASIIphsGUI as G2phG
27import GSASIIlattice as G2lat
28import GSASIIspc as G2spc
29from  OpenGL.GL import *
30from OpenGL.GLU import *
31from OpenGL.GLUT import *
32from OpenGL.GLE import *
33from matplotlib.backends.backend_wxagg import FigureCanvasWxAgg as Canvas
34from matplotlib.backends.backend_wxagg import NavigationToolbar2Wx as Toolbar
35
36# useful degree trig functions
37sind = lambda x: math.sin(x*math.pi/180.)
38cosd = lambda x: math.cos(x*math.pi/180.)
39tand = lambda x: math.tan(x*math.pi/180.)
40asind = lambda x: 180.*math.asin(x)/math.pi
41acosd = lambda x: 180.*math.acos(x)/math.pi
42atan2d = lambda x,y: 180.*math.atan2(y,x)/math.pi
43atand = lambda x: 180.*math.atan(x)/math.pi
44# numpy versions
45npsind = lambda x: np.sin(x*np.pi/180.)
46npcosd = lambda x: np.cos(x*np.pi/180.)
47npacosd = lambda x: 180.*np.arccos(x)/np.pi
48   
49class G2PlotMpl(wx.Panel):   
50    def __init__(self,parent,id=-1,dpi=None,**kwargs):
51        wx.Panel.__init__(self,parent,id=id,**kwargs)
52        self.figure = mpl.figure.Figure(dpi=dpi,figsize=(5,7))
53        self.canvas = Canvas(self,-1,self.figure)
54        self.toolbar = Toolbar(self.canvas)
55
56        self.toolbar.Realize()
57       
58        sizer=wx.BoxSizer(wx.VERTICAL)
59        sizer.Add(self.canvas,1,wx.EXPAND)
60        sizer.Add(self.toolbar,0,wx.LEFT|wx.EXPAND)
61        self.SetSizer(sizer)
62       
63class G2PlotOgl(wx.Panel):
64    def __init__(self,parent,id=-1,dpi=None,**kwargs):
65        self.figure = wx.Panel.__init__(self,parent,id=id,**kwargs)
66        self.canvas = wx.glcanvas.GLCanvas(self,-1,**kwargs)
67        self.camera = {}
68        sizer=wx.BoxSizer(wx.VERTICAL)
69        sizer.Add(self.canvas,1,wx.EXPAND)
70        self.SetSizer(sizer)
71       
72class G2Plot3D(wx.Panel):
73    def __init__(self,parent,id=-1,dpi=None,**kwargs):
74        wx.Panel.__init__(self,parent,id=id,**kwargs)
75        self.figure = mpl.figure.Figure(dpi=dpi,figsize=(6,6))
76        self.canvas = Canvas(self,-1,self.figure)
77        self.toolbar = Toolbar(self.canvas)
78
79        self.toolbar.Realize()
80       
81        sizer=wx.BoxSizer(wx.VERTICAL)
82        sizer.Add(self.canvas,1,wx.EXPAND)
83        sizer.Add(self.toolbar,0,wx.LEFT|wx.EXPAND)
84        self.SetSizer(sizer)
85                             
86class G2PlotNoteBook(wx.Panel):
87    def __init__(self,parent,id=-1):
88        wx.Panel.__init__(self,parent,id=id)
89        #so one can't delete a plot page!!
90        self.nb = wx.aui.AuiNotebook(self, \
91            style=wx.aui.AUI_NB_DEFAULT_STYLE ^ wx.aui.AUI_NB_CLOSE_ON_ACTIVE_TAB)
92        sizer = wx.BoxSizer()
93        sizer.Add(self.nb,1,wx.EXPAND)
94        self.SetSizer(sizer)
95        self.status = parent.CreateStatusBar()
96        self.status.SetFieldsCount(2)
97        self.status.SetStatusWidths([125,-1])
98        self.Bind(wx.aui.EVT_AUINOTEBOOK_PAGE_CHANGED, self.OnPageChanged)
99       
100        self.plotList = []
101           
102    def addMpl(self,name=""):
103        page = G2PlotMpl(self.nb)
104        self.nb.AddPage(page,name)
105       
106        self.plotList.append(name)
107       
108        return page.figure
109       
110    def add3D(self,name=""):
111        page = G2Plot3D(self.nb)
112        self.nb.AddPage(page,name)
113       
114        self.plotList.append(name)
115       
116        return page.figure
117       
118    def addOgl(self,name=""):
119        page = G2PlotOgl(self.nb)
120        self.nb.AddPage(page,name)
121       
122        self.plotList.append(name)
123       
124        return page.figure
125       
126    def Delete(self,name):
127        try:
128            item = self.plotList.index(name)
129            del self.plotList[item]
130            self.nb.DeletePage(item)
131        except ValueError:          #no plot of this name - do nothing
132            return     
133               
134    def clear(self):
135        while self.nb.GetPageCount():
136            self.nb.DeletePage(0)
137        self.plotList = []
138        self.status.DestroyChildren()
139       
140    def Rename(self,oldName,newName):
141        try:
142            item = self.plotList.index(oldName)
143            self.plotList[item] = newName
144            self.nb.SetPageText(item,newName)
145        except ValueError:          #no plot of this name - do nothing
146            return     
147       
148    def OnPageChanged(self,event):       
149        if self.plotList:
150            self.status.SetStatusText('Better to select this from GSAS-II data tree',1)
151        self.status.DestroyChildren()                           #get rid of special stuff on status bar
152       
153def PlotSngl(self,newPlot=False):
154    from matplotlib.patches import Circle
155    global HKL,HKLF
156
157    def OnSCMotion(event):
158        xpos = event.xdata
159        if xpos:
160            xpos = round(xpos)                                        #avoid out of frame mouse position
161            ypos = round(event.ydata)
162            zpos = Data['Layer']
163            if '100' in Data['Zone']:
164                HKLtxt = '(%3d,%3d,%3d)'%(zpos,xpos,ypos)
165            elif '010' in Data['Zone']:
166                HKLtxt = '(%3d,%3d,%3d)'%(xpos,zpos,ypos)
167            elif '001' in Data['Zone']:
168                HKLtxt = '(%3d,%3d,%3d)'%(xpos,ypos,zpos)
169            Page.canvas.SetToolTipString(HKLtxt)
170            self.G2plotNB.status.SetFields(['HKL = '+HKLtxt,''])
171               
172    def OnSCPick(event):
173        zpos = Data['Layer']
174        pos = event.artist.center
175        if '100' in Data['Zone']:
176            Page.canvas.SetToolTipString('(picked:(%3d,%3d,%3d))'%(zpos,pos[0],pos[1]))
177            hkl = np.array([zpos,pos[0],pos[1]])
178        elif '010' in Data['Zone']:
179            Page.canvas.SetToolTipString('(picked:(%3d,%3d,%3d))'%(pos[0],zpos,pos[1]))
180            hkl = np.array([pos[0],zpos,pos[1]])
181        elif '001' in Data['Zone']:
182            Page.canvas.SetToolTipString('(picked:(%3d,%3d,%3d))'%(pos[0],pos[1],zpos))
183            hkl = np.array([pos[0],pos[1],zpos])
184        h,k,l = hkl
185        hklf = HKLF[np.where(np.all(HKL-hkl == [0,0,0],axis=1))]
186        if len(hklf):
187            Fosq,sig,Fcsq = hklf[0]
188            HKLtxt = '(%3d,%3d,%3d %.2f %.3f %.2f %.2f)'%(h,k,l,Fosq,sig,Fcsq,(Fosq-Fcsq)/(scale*sig))
189            self.G2plotNB.status.SetFields(['','HKL, Fosq, sig, Fcsq, delFsq/sig = '+HKLtxt])
190                                 
191    def OnSCKeyPress(event):
192        print event.key
193
194    try:
195        plotNum = self.G2plotNB.plotList.index('Structure Factors')
196        Page = self.G2plotNB.nb.GetPage(plotNum)
197        if not newPlot:
198            Plot = Page.figure.gca()          #get previous powder plot & get limits
199            xylim = Plot.get_xlim(),Plot.get_ylim()
200        Page.figure.clf()
201        Plot = Page.figure.gca()          #get a fresh plot after clf()
202    except ValueError,error:
203        Plot = self.G2plotNB.addMpl('Structure Factors').gca()
204        plotNum = self.G2plotNB.plotList.index('Structure Factors')
205        Page = self.G2plotNB.nb.GetPage(plotNum)
206#        Page.canvas.mpl_connect('key_press_event', OnSCKeyPress)
207        Page.canvas.mpl_connect('pick_event', OnSCPick)
208        Page.canvas.mpl_connect('motion_notify_event', OnSCMotion)
209    Page.SetFocus()
210   
211    Plot.set_aspect(aspect='equal')
212    HKLref = self.PatternTree.GetItemPyData(self.Sngl)
213    Data = self.PatternTree.GetItemPyData( \
214        G2gd.GetPatternTreeItemId(self,self.Sngl, 'HKL Plot Controls'))
215    Type = Data['Type']           
216    scale = Data['Scale']
217    HKLmax = Data['HKLmax']
218    HKLmin = Data['HKLmin']
219    FosqMax = Data['FoMax']
220    FoMax = math.sqrt(FosqMax)
221    ifFc = Data['ifFc']
222    xlabel = ['k, h=','h, k=','h, l=']
223    ylabel = ['l','l','k']
224    zones = ['100','010','001']
225    pzone = [[1,2],[0,2],[0,1]]
226    izone = zones.index(Data['Zone'])
227    Plot.set_title(self.PatternTree.GetItemText(self.Sngl)[5:])
228    HKL = []
229    HKLF = []
230    for H,Fosq,sig,Fcsq,x,x,x in HKLref:
231        HKL.append(H)
232        HKLF.append([Fosq,sig,Fcsq])
233        if H[izone] == Data['Layer']:
234            B = 0
235            if Type == 'Fosq':
236                A = scale*Fosq/FosqMax
237                B = scale*Fcsq/FosqMax
238                C = abs(A-B)
239            elif Type == 'Fo':
240                A = scale*math.sqrt(max(0,Fosq))/FoMax
241                B = scale*math.sqrt(max(0,Fcsq))/FoMax
242                C = abs(A-B)
243            elif Type == '|DFsq|/sig':
244                A = abs(Fosq-Fcsq)/(scale*sig)
245            elif Type == '|DFsq|>sig':
246                A = abs(Fosq-Fcsq)/(scale*sig)
247                if A < 1.0: A = 0                   
248            elif Type == '|DFsq|>3sig':
249                A = abs(Fosq-Fcsq)/(scale*sig)
250                if A < 3.0: A = 0                   
251            xy = (H[pzone[izone][0]],H[pzone[izone][1]])
252            if A > 0.0:
253                Plot.add_artist(Circle(xy,radius=A,ec='g',fc='w',picker=3))
254            if B:
255                Plot.add_artist(Circle(xy,radius=B,ec='b',fc='w'))
256                radius = C
257                if radius > 0:
258                    if A > B:
259                        Plot.add_artist(Circle(xy,radius=radius,ec='g',fc='g'))
260                    else:                   
261                        Plot.add_artist(Circle(xy,radius=radius,ec='r',fc='r'))
262    HKL = np.array(HKL,dtype=np.int)
263    HKLF = np.array(HKLF)
264    Plot.set_xlabel(xlabel[izone]+str(Data['Layer']),fontsize=12)
265    Plot.set_ylabel(ylabel[izone],fontsize=12)
266    Plot.set_xlim((HKLmin[pzone[izone][0]],HKLmax[pzone[izone][0]]))
267    Plot.set_ylim((HKLmin[pzone[izone][1]],HKLmax[pzone[izone][1]]))
268    if not newPlot:
269        Page.toolbar.push_current()
270        Plot.set_xlim(xylim[0])
271        Plot.set_ylim(xylim[1])
272        xylim = []
273        Page.toolbar.push_current()
274        Page.toolbar.draw()
275    else:
276        Page.canvas.draw()
277       
278def PlotPatterns(self,newPlot=False):
279    global HKL
280   
281    def OnPick(event):
282        if self.itemPicked is not None: return
283        PatternId = self.PatternId
284        PickId = self.PickId
285        pick = event.artist
286        mouse = event.mouseevent
287        xpos = pick.get_xdata()
288        ypos = pick.get_ydata()
289        ind = event.ind
290        xy = zip(np.take(xpos,ind),np.take(ypos,ind))[0]
291        if self.PatternTree.GetItemText(PickId) == 'Peak List':
292            if ind.all() != [0]:                                    #picked a data point
293                inst = self.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(self,PatternId, 'Instrument Parameters'))
294                if len(inst[1]) == 11:
295                    ins = inst[1][4:10]
296                else:
297                    ins = inst[1][6:12]   
298                sig = ins[0]*tand(xy[0]/2.0)**2+ins[1]*tand(xy[0]/2.0)+ins[2]
299                gam = ins[3]/cosd(xy[0]/2.0)+ins[4]*tand(xy[0]/2.0)           
300                data = self.PatternTree.GetItemPyData(self.PickId)
301                XY = [xy[0],0, xy[1],1, sig,0, gam,0]       #default refine intensity 1st   
302                data.append(XY)
303                G2pdG.UpdatePeakGrid(self,data)
304                PlotPatterns(self)
305            else:                                                   #picked a peak list line
306                self.itemPicked = pick
307        elif self.PatternTree.GetItemText(PickId) == 'Limits':
308            if ind.all() != [0]:                                    #picked a data point
309                LimitId = G2gd.GetPatternTreeItemId(self,PatternId, 'Limits')
310                data = self.PatternTree.GetItemPyData(LimitId)
311                if mouse.button==1:
312                    data[1][0] = min(xy[0],data[1][1])
313                if mouse.button==3:
314                    data[1][1] = max(xy[0],data[1][0])
315                self.PatternTree.SetItemPyData(LimitId,data)
316                G2pdG.UpdateLimitsGrid(self,data)
317                PlotPatterns(self)
318            else:                                                   #picked a limit line
319                self.itemPicked = pick
320       
321    def OnPlotKeyPress(event):
322        newPlot = False
323        if event.key == 'w':
324            if self.Weight:
325                self.Weight = False
326            else:
327                self.Weight = True
328            print 'plot weighting:',self.Weight
329        elif event.key == 'l':
330            if self.Contour:
331                pass
332            else:
333                if self.logPlot:
334                    self.logPlot = False
335                else:
336                    self.Offset = 0
337                    self.logPlot = True
338        elif event.key == 'u':
339            if self.Contour:
340                self.Cmax = min(1.0,self.Cmax*1.2)
341            elif self.logPlot:
342                pass
343            elif self.Offset < 100.:
344                self.Offset += 1.
345        elif event.key == 'd':
346            if self.Contour:
347                self.Cmax = max(0.0,self.Cmax*0.8)
348            elif self.logPlot:
349                pass
350            elif self.Offset > 0.:
351                self.Offset -= 1.
352        elif event.key == 'c':
353            newPlot = True
354            if self.Contour:
355                self.Contour = False
356            else:
357                self.Contour = True
358                self.SinglePlot = False
359                self.Offset = 0
360        elif event.key == 's':
361            if self.Contour:
362                choice = [m for m in mpl.cm.datad.keys() if not m.endswith("_r")]
363                dlg = wx.SingleChoiceDialog(self,'Select','Color scheme',choice)
364                if dlg.ShowModal() == wx.ID_OK:
365                    sel = dlg.GetSelection()
366                    self.ContourColor = choice[sel]
367                else:
368                    self.ContourColor = 'Paired'
369                dlg.Destroy()
370            else:               
371                if self.SinglePlot:
372                    self.SinglePlot = False
373                else:
374                    self.SinglePlot = True
375        elif event.key == '+':
376            if self.PickId:
377                self.PickId = False
378        elif event.key == 'i':
379            choice = ['nearest','bilinear','bicubic','spline16','spline36','hanning',
380               'hamming','hermite','kaiser','quadric','catrom','gaussian','bessel',
381               'mitchell','sinc','lanczos']
382            dlg = wx.SingleChoiceDialog(self,'Select','Interpolation',choice)
383            if dlg.ShowModal() == wx.ID_OK:
384                sel = dlg.GetSelection()
385                self.Interpolate = choice[sel]
386            else:
387                self.Interpolate = 'nearest'
388            dlg.Destroy()
389           
390        PlotPatterns(self,newPlot=newPlot)
391       
392    def OnKeyBox(event):
393        if self.G2plotNB.nb.GetSelection() == self.G2plotNB.plotList.index('Powder Patterns'):
394            event.key = cb.GetValue()[0]
395            cb.SetValue(' key press')
396            OnPlotKeyPress(event)
397                       
398    def OnMotion(event):
399        xpos = event.xdata
400        if xpos:                                        #avoid out of frame mouse position
401            ypos = event.ydata
402            Page.canvas.SetCursor(wx.CROSS_CURSOR)
403            try:
404                wave = self.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(self, \
405                    self.PatternId, 'Instrument Parameters'))[0][1]
406                dsp = 0.0
407                if abs(xpos) > 0.:                  #avoid possible singularity at beam center
408                    dsp = wave/(2.*sind(abs(xpos)/2.0))
409                if self.Contour:
410                    self.G2plotNB.status.SetStatusText('2-theta =%9.3f d =%9.5f pattern ID =%5d'%(xpos,dsp,int(ypos)),1)
411                else:
412                    self.G2plotNB.status.SetStatusText('2-theta =%9.3f d =%9.5f Intensity =%9.1f'%(xpos,dsp,ypos),1)
413                if self.itemPicked:
414                    Page.canvas.SetToolTipString('%9.3f'%(xpos))
415                if self.PickId and self.PatternTree.GetItemText(self.PickId) in ['Index Peak List','Unit Cells List']:
416                    found = []
417                    if len(HKL):
418                        view = Page.toolbar._views.forward()[0][:2]
419                        wid = view[1]-view[0]
420                        found = HKL[np.where(np.fabs(HKL.T[5]-xpos) < 0.002*wid)]
421                    if len(found):
422                        h,k,l = found[0][:3] 
423                        Page.canvas.SetToolTipString('%d,%d,%d'%(int(h),int(k),int(l)))
424                    else:
425                        Page.canvas.SetToolTipString('')
426
427            except TypeError:
428                self.G2plotNB.status.SetStatusText('Select PWDR powder pattern first',1)
429                                                   
430    def OnRelease(event):
431        if self.itemPicked is None: return
432        xpos = event.xdata
433        if xpos:                                        #avoid out of frame mouse position
434            lines = []
435            for line in self.Lines: lines.append(line.get_xdata()[0])
436            lineNo = lines.index(self.itemPicked.get_xdata()[0])
437            if  lineNo in [0,1]:
438                LimitId = G2gd.GetPatternTreeItemId(self,self.PatternId, 'Limits')
439                data = self.PatternTree.GetItemPyData(LimitId)
440#                print 'limits',xpos
441                data[1][lineNo] = xpos
442                self.PatternTree.SetItemPyData(LimitId,data)
443                if self.PatternTree.GetItemText(self.PickId) == 'Limits':
444                    G2pdG.UpdateLimitsGrid(self,data)
445            else:
446                PeakId = G2gd.GetPatternTreeItemId(self,self.PatternId, 'Peak List')
447                data = self.PatternTree.GetItemPyData(PeakId)
448#                print 'peaks',xpos
449                data[lineNo-2][0] = xpos
450                self.PatternTree.SetItemPyData(PeakId,data)
451                G2pdG.UpdatePeakGrid(self,data)
452        PlotPatterns(self)
453        self.itemPicked = None   
454
455    xylim = []
456    try:
457        plotNum = self.G2plotNB.plotList.index('Powder Patterns')
458        Page = self.G2plotNB.nb.GetPage(plotNum)
459        if not newPlot:
460            Plot = Page.figure.gca()          #get previous powder plot & get limits
461            xylim = Plot.get_xlim(),Plot.get_ylim()
462        Page.figure.clf()
463        Plot = Page.figure.gca()          #get a fresh plot after clf()
464    except ValueError,error:
465        newPlot = True
466        self.Cmax = 1.0
467        Plot = self.G2plotNB.addMpl('Powder Patterns').gca()
468        plotNum = self.G2plotNB.plotList.index('Powder Patterns')
469        Page = self.G2plotNB.nb.GetPage(plotNum)
470        Page.canvas.mpl_connect('key_press_event', OnPlotKeyPress)
471        Page.canvas.mpl_connect('motion_notify_event', OnMotion)
472        Page.canvas.mpl_connect('pick_event', OnPick)
473        Page.canvas.mpl_connect('button_release_event', OnRelease)
474    Page.SetFocus()
475    self.G2plotNB.status.DestroyChildren()
476    if self.Contour:
477        Choice = (' key press','d: lower contour max','u: raise contour max',
478            'i: interpolation method','s: color scheme','c: contour off')
479    else:
480        if self.logPlot:
481            Choice = (' key press','l: log(I) off',
482                'c: contour on','s: toggle single plot','+: no selection')
483        else:
484            Choice = (' key press','d: offset down','u: offset up','l: log(I) on',
485                'c: contour on','s: toggle single plot','+: no selection')
486    cb = wx.ComboBox(self.G2plotNB.status,style=wx.CB_DROPDOWN|wx.CB_READONLY,
487        choices=Choice)
488    cb.Bind(wx.EVT_COMBOBOX, OnKeyBox)
489    cb.SetValue(' key press')
490   
491    PickId = self.PickId
492    PatternId = self.PatternId
493    colors=['b','g','r','c','m','k']
494    Lines = []
495    if self.SinglePlot:
496        Pattern = self.PatternTree.GetItemPyData(self.PatternId)
497        Pattern.append(self.PatternTree.GetItemText(self.PatternId))
498        PlotList = [Pattern,]
499    else:       
500        PlotList = []
501        item, cookie = self.PatternTree.GetFirstChild(self.root)
502        while item:
503            if 'PWDR' in self.PatternTree.GetItemText(item):
504                Pattern = self.PatternTree.GetItemPyData(item)
505                if len(Pattern) < 3:                    # put name on end if needed
506                    Pattern.append(self.PatternTree.GetItemText(item))
507                PlotList.append(Pattern)
508            item, cookie = self.PatternTree.GetNextChild(self.root, cookie)               
509    Ymax = 1.0
510    lenX = 0
511    HKL = np.array(self.HKL)
512    for Pattern in PlotList:
513        xye = Pattern[1]
514        Ymax = max(Ymax,max(xye[1]))
515    offset = self.Offset*Ymax/100.0
516    Plot.set_title('Powder Patterns: '+os.path.split(self.GSASprojectfile)[1])
517    Plot.set_xlabel(r'$\mathsf{2\theta}$',fontsize=14)
518    Plot.set_ylabel('Intensity',fontsize=12)
519    if self.Contour:
520        ContourZ = []
521        ContourY = []
522        Nseq = 0
523    for N,Pattern in enumerate(PlotList):
524        ifpicked = False
525        LimitId = 0
526        xye = np.array(Pattern[1])
527        if PickId:
528            ifpicked = Pattern[2] == self.PatternTree.GetItemText(PatternId)
529            LimitId = G2gd.GetPatternTreeItemId(self,PatternId, 'Limits')
530        X = xye[0]
531        if not lenX:
532            lenX = len(X)           
533        Y = xye[1]+offset*N
534        if LimitId:
535            limits = self.PatternTree.GetItemPyData(LimitId)
536            Lines.append(Plot.axvline(limits[1][0],color='g',dashes=(5,5),picker=3.))   
537            Lines.append(Plot.axvline(limits[1][1],color='r',dashes=(5,5),picker=3.))                   
538        if self.Contour:
539           
540            if lenX == len(X):
541                ContourY.append(N)
542                ContourZ.append(Y)
543                ContourX = X
544                Nseq += 1
545                Plot.set_ylabel('Data sequence',fontsize=12)
546        else:
547            if ifpicked:
548                Z = xye[3]+offset*N
549                W = xye[4]+offset*N
550                D = xye[5]+offset*N-Ymax*.02
551                if self.Weight:
552                    W2 = np.sqrt(xye[2])
553                    D *= W2-Ymax*.02
554                if self.logPlot:
555                    Plot.semilogy(X,Y,colors[N%6]+'+',picker=3.,clip_on=False,nonposy='mask')
556                    Plot.semilogy(X,Z,colors[(N+1)%6],picker=False,nonposy='mask')
557                    Plot.semilogy(X,W,colors[(N+2)%6],picker=False,nonposy='mask')
558                else:
559                    Plot.plot(X,Y,colors[N%6]+'+',picker=3.,clip_on=False)
560                    Plot.plot(X,Z,colors[(N+1)%6],picker=False)
561                    Plot.plot(X,W,colors[(N+2)%6],picker=False)
562                    Plot.plot(X,D,colors[(N+3)%6],picker=False)
563                    Plot.axhline(0.,color=wx.BLACK)
564                Page.canvas.SetToolTipString('')
565                if self.PatternTree.GetItemText(PickId) == 'Peak List':
566                    tip = 'On data point: Pick peak - L or R MB.On line: MB down to move'
567                    Page.canvas.SetToolTipString(tip)
568                    data = self.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(self,PatternId, 'Peak List'))
569                    for item in data:
570                        Lines.append(Plot.axvline(item[0],color=colors[N%6],picker=2.))
571                if self.PatternTree.GetItemText(PickId) == 'Limits':
572                    tip = 'On data point: Lower limit - L MB; Upper limit - R MB. On limit: MB down to move'
573                    Page.canvas.SetToolTipString(tip)
574                    data = self.LimitsTable.GetData()
575            else:
576                if self.logPlot:
577                    Plot.semilogy(X,Y,colors[N%6],picker=False,nonposy='clip')
578                else:
579                    Plot.plot(X,Y,colors[N%6],picker=False)
580    if PickId and self.PatternTree.GetItemText(PickId) in ['Index Peak List','Unit Cells List']:
581        peaks = np.array((self.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(self,PatternId, 'Index Peak List'))))
582        for peak in peaks:
583            Plot.axvline(peak[0],color='b')
584        for hkl in self.HKL:
585            Plot.axvline(hkl[5],color='r',dashes=(5,5))
586    if self.Contour:
587        acolor = mpl.cm.get_cmap(self.ContourColor)
588        Img = Plot.imshow(ContourZ,cmap=acolor,vmin=0,vmax=Ymax*self.Cmax,interpolation=self.Interpolate, 
589            extent=[ContourX[0],ContourX[-1],ContourY[0],ContourY[-1]],aspect='auto',origin='lower')
590        Page.figure.colorbar(Img)
591    else:
592        self.Lines = Lines
593    if not newPlot:
594        Page.toolbar.push_current()
595        Plot.set_xlim(xylim[0])
596        Plot.set_ylim(xylim[1])
597        xylim = []
598        Page.toolbar.push_current()
599        Page.toolbar.draw()
600    else:
601        Page.canvas.draw()
602    self.Pwdr = True
603
604def PlotPowderLines(self):
605    global HKL
606
607    def OnMotion(event):
608        xpos = event.xdata
609        if xpos:                                        #avoid out of frame mouse position
610            Page.canvas.SetCursor(wx.CROSS_CURSOR)
611            self.G2plotNB.status.SetFields(['','2-theta =%9.3f '%(xpos,)])
612            if self.PickId and self.PatternTree.GetItemText(self.PickId) in ['Index Peak List','Unit Cells List']:
613                found = []
614                if len(HKL):
615                    view = Page.toolbar._views.forward()[0][:2]
616                    wid = view[1]-view[0]
617                    found = HKL[np.where(np.fabs(HKL.T[5]-xpos) < 0.002*wid)]
618                if len(found):
619                    h,k,l = found[0][:3] 
620                    Page.canvas.SetToolTipString('%d,%d,%d'%(int(h),int(k),int(l)))
621                else:
622                    Page.canvas.SetToolTipString('')
623
624    try:
625        plotNum = self.G2plotNB.plotList.index('Powder Lines')
626        Page = self.G2plotNB.nb.GetPage(plotNum)
627        Page.figure.clf()
628        Plot = Page.figure.gca()
629    except ValueError,error:
630        Plot = self.G2plotNB.addMpl('Powder Lines').gca()
631        plotNum = self.G2plotNB.plotList.index('Powder Lines')
632        Page = self.G2plotNB.nb.GetPage(plotNum)
633        Page.canvas.mpl_connect('motion_notify_event', OnMotion)
634       
635    Page.SetFocus()
636    Plot.set_title('Powder Pattern Lines')
637    Plot.set_xlabel(r'$\mathsf{2\theta}$',fontsize=14)
638    PickId = self.PickId
639    PatternId = self.PatternId
640    peaks = self.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(self,PatternId, 'Index Peak List'))
641    for peak in peaks:
642        Plot.axvline(peak[0],color='b')
643    HKL = np.array(self.HKL)
644    for hkl in self.HKL:
645        Plot.axvline(hkl[5],color='r',dashes=(5,5))
646    xmin = peaks[0][0]
647    xmax = peaks[-1][0]
648    delt = xmax-xmin
649    xlim = [max(0,xmin-delt/20.),min(180.,xmax+delt/20.)]
650    Plot.set_xlim(xlim)
651    Page.canvas.draw()
652
653def PlotPeakWidths(self):
654    PatternId = self.PatternId
655    limitID = G2gd.GetPatternTreeItemId(self,PatternId, 'Limits')
656    if limitID:
657        limits = self.PatternTree.GetItemPyData(limitID)
658    else:
659        return
660    instParms = self.PatternTree.GetItemPyData( \
661        G2gd.GetPatternTreeItemId(self,PatternId, 'Instrument Parameters'))
662    if instParms[0][0] == 'PXC':
663        lam = instParms[1][1]
664        if len(instParms[1]) == 13:
665            GU,GV,GW,LX,LY = instParms[0][6:11]
666        else:
667            GU,GV,GW,LX,LY = instParms[0][4:9]
668    peakID = G2gd.GetPatternTreeItemId(self,PatternId, 'Peak List')
669    if peakID:
670        peaks = self.PatternTree.GetItemPyData(peakID)
671    else:
672        peaks = []
673   
674    try:
675        plotNum = self.G2plotNB.plotList.index('Peak Widths')
676        Page = self.G2plotNB.nb.GetPage(plotNum)
677        Page.figure.clf()
678        Plot = Page.figure.gca()
679    except ValueError,error:
680        Plot = self.G2plotNB.addMpl('Peak Widths').gca()
681        plotNum = self.G2plotNB.plotList.index('Peak Widths')
682        Page = self.G2plotNB.nb.GetPage(plotNum)
683    Page.SetFocus()
684   
685    Page.canvas.SetToolTipString('')
686    colors=['b','g','r','c','m','k']
687    Xmin,Xmax = limits[1]
688    Xmin = min(0.5,max(Xmin,1))
689    Xmin /= 2
690    Xmax /= 2
691    nPts = 100
692    delt = (Xmax-Xmin)/nPts
693    thetas = []
694    for i in range(nPts):
695        thetas.append(Xmin+i*delt)
696    X = []
697    Y = []
698    Z = []
699    W = []
700    V = []
701    sig = lambda Th,U,V,W: 1.17741*math.sqrt(U*tand(Th)**2+V*tand(Th)+W)*math.pi/18000.
702    gam = lambda Th,X,Y: (X/cosd(Th)+Y*tand(Th))*math.pi/18000.
703    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.)
704    gamFW2 = lambda s,g: math.sqrt(s**2+(0.4654996*g)**2)+.5345004*#Ubaldo Bafile - private communication
705    for theta in thetas:
706        X.append(4.0*math.pi*sind(theta)/lam)              #q
707        s = sig(theta,GU,GV,GW)
708        g = gam(theta,LX,LY)
709        G = gamFW(g,s)
710        H = gamFW2(g,s)
711        Y.append(s/tand(theta))
712        Z.append(g/tand(theta))
713        W.append(G/tand(theta))
714        V.append(H/tand(theta))
715    Plot.set_title('Instrument and sample peak widths')
716    Plot.set_ylabel(r'$\Delta q/q, \Delta d/d$',fontsize=14)
717    Plot.set_xlabel(r'$q, \AA^{-1}$',fontsize=14)
718    Plot.plot(X,Y,color='r',label='Gaussian')
719    Plot.plot(X,Z,color='g',label='Lorentzian')
720    Plot.plot(X,W,color='b',label='G+L')
721    Plot.plot(X,V,color='k',label='G+L2')
722    X = []
723    Y = []
724    Z = []
725    W = []
726    V = []
727    for peak in peaks:
728        X.append(4.0*math.pi*sind(peak[0]/2.0)/lam)
729        s = 1.17741*math.sqrt(peak[4])*math.pi/18000.
730        g = peak[6]*math.pi/18000.
731        G = gamFW(g,s)
732        H = gamFW2(g,s)
733        Y.append(s/tand(peak[0]/2.))
734        Z.append(g/tand(peak[0]/2.))
735        W.append(G/tand(peak[0]/2.))
736        V.append(H/tand(peak[0]/2.))
737    Plot.plot(X,Y,'+',color='r',label='G peak')
738    Plot.plot(X,Z,'+',color='g',label='L peak')
739    Plot.plot(X,W,'+',color='b',label='G+L peak')
740    Plot.plot(X,V,'+',color='k',label='G+L2 peak')
741    Plot.legend(loc='best')
742    Page.canvas.draw()
743   
744def PlotStrain(self,data):
745    # in this instance data is for a phase
746    PatternId = self.PatternId
747    generalData = data['General']
748    SGData = generalData['SGData']
749    MuStrKeys = G2spc.MustrainNames(SGData)
750    cell = generalData['Cell'][1:]
751    A,B = G2lat.cell2AB(cell[:6])
752    Vol = cell[6]
753    useList = data['Histograms']
754    numPlots = len(useList)
755   
756    try:
757        plotNum = self.G2plotNB.plotList.index('Microstrain')
758        Page = self.G2plotNB.nb.GetPage(plotNum)
759        Page.figure.clf()
760        Plot = mp3d.Axes3D(Page.figure)
761    except ValueError,error:
762        Plot = mp3d.Axes3D(self.G2plotNB.add3D('Microstrain'))
763        plotNum = self.G2plotNB.plotList.index('Microstrain')
764        Page = self.G2plotNB.nb.GetPage(plotNum)
765    Page.SetFocus()
766    self.G2plotNB.status.SetStatusText('Adjust frame size to get desired aspect ratio',1)
767   
768    for item in useList:
769        if useList[item]['Show']:
770            muStrain = useList[item]['Mustrain']
771            PHI = np.linspace(0.,360.,30,True)
772            PSI = np.linspace(0.,180.,30,True)
773            X = np.outer(npsind(PHI),npsind(PSI))
774            Y = np.outer(npcosd(PHI),npsind(PSI))
775            Z = np.outer(np.ones(np.size(PHI)),npcosd(PSI))
776            if muStrain[0] == 'isotropic':
777                muiso = muStrain[1][0]*math.pi/0.018
778                X *= muiso
779                Y *= muiso
780                Z *= muiso               
781               
782            elif muStrain[0] == 'uniaxial':
783                def uniaxMustrain(xyz,muiso,muaniso,axes):
784                    cp = abs(np.dot(xyz,axes))
785                    S = muiso+muaniso*cp
786                    return S*xyz*math.pi/0.018
787                muiso,muaniso = muStrain[1][:2]
788                axes = np.inner(A,np.array(muStrain[3]))
789                axes /= nl.norm(axes)
790                Shkl = np.array(muStrain[1])
791                Shape = X.shape[0]
792                XYZ = np.dstack((X,Y,Z))
793                XYZ = np.nan_to_num(np.apply_along_axis(uniaxMustrain,2,XYZ,muiso,muaniso,axes))
794                X,Y,Z = np.dsplit(XYZ,3)
795                X = X[:,:,0]
796                Y = Y[:,:,0]
797                Z = Z[:,:,0]
798               
799            elif muStrain[0] == 'generalized':
800                def genMustrain(xyz,SGData,A,Shkl):
801                    uvw = np.inner(A.T,xyz)
802                    Strm = np.array(G2spc.MustrainCoeff(uvw,SGData))
803                    sum = np.sqrt(np.sum(np.multiply(Shkl,Strm)))*math.pi/0.018
804                    return sum*xyz
805                Shkl = np.array(muStrain[4])
806                if np.any(Shkl):
807                    Shape = X.shape[0]
808                    XYZ = np.dstack((X,Y,Z))
809                    XYZ = np.nan_to_num(np.apply_along_axis(genMustrain,2,XYZ,SGData,A,Shkl))
810                    X,Y,Z = np.dsplit(XYZ,3)
811                    X = X[:,:,0]
812                    Y = Y[:,:,0]
813                    Z = Z[:,:,0]
814                   
815            if np.any(X) and np.any(Y) and np.any(Z):
816                Plot.plot_surface(X,Y,Z,rstride=1,cstride=1,color='g')
817               
818            Plot.set_xlabel('X')
819            Plot.set_ylabel('Y')
820            Plot.set_zlabel('Z')
821    Page.canvas.draw()
822           
823def PlotExposedImage(self,newPlot=False,event=None):
824    plotNo = self.G2plotNB.nb.GetSelection()
825    if self.G2plotNB.nb.GetPageText(plotNo) == '2D Powder Image':
826        PlotImage(self,newPlot,event)
827    elif self.G2plotNB.nb.GetPageText(plotNo) == '2D Integration':
828        PlotIntegration(self,newPlot,event)
829
830def PlotImage(self,newPlot=False,event=None):
831    from matplotlib.patches import Ellipse,Arc,Circle,Polygon
832    import numpy.ma as ma
833    Dsp = lambda tth,wave: wave/(2.*sind(tth/2.))
834    global Data,Masks
835    Data = self.PatternTree.GetItemPyData(
836        G2gd.GetPatternTreeItemId(self,self.Image, 'Image Controls'))
837    Masks = self.PatternTree.GetItemPyData(
838        G2gd.GetPatternTreeItemId(self,self.Image, 'Masks'))
839
840    def OnImMotion(event):
841        Page.canvas.SetToolTipString('')
842        sizexy = Data['size']
843        if event.xdata and event.ydata:                 #avoid out of frame errors
844            Page.canvas.SetCursor(wx.CROSS_CURSOR)
845            item = self.itemPicked
846            pixelSize = Data['pixelSize']
847            scalex = 1000./pixelSize[0]
848            scaley = 1000./pixelSize[1]
849            if item and self.PatternTree.GetItemText(self.PickId) == 'Image Controls':
850                if 'Text' in str(item):
851                    Page.canvas.SetToolTipString('%8.3f %8.3fmm'%(event.xdata,event.ydata))
852                else:
853                    xcent,ycent = Data['center']
854                    xpos = event.xdata-xcent
855                    ypos = event.ydata-ycent
856                    if 'line3' in  str(item) or 'line4' in str(item) and not Data['fullIntegrate']:
857                        ang = int(atan2d(xpos,ypos))
858                        Page.canvas.SetToolTipString('%6d deg'%(ang))
859                    elif 'line1' in  str(item) or 'line2' in str(item):
860                        tth = G2img.GetTth(event.xdata,event.ydata,Data)
861                        Page.canvas.SetToolTipString('%8.3fdeg'%(tth))                           
862            else:
863                xpos = event.xdata
864                ypos = event.ydata
865                xpix = xpos*scalex
866                ypix = ypos*scaley
867                Int = 0
868                if (0 <= xpix <= sizexy[0]) and (0 <= ypix <= sizexy[1]):
869                    Int = self.ImageZ[ypix][xpix]
870                tth,azm,dsp = G2img.GetTthAzmDsp(xpos,ypos,Data)
871                Q = 2.*math.pi/dsp
872                if self.setPoly:
873                    self.G2plotNB.status.SetFields(['','Polygon mask pick - LB next point, RB close polygon'])
874                else:
875                    self.G2plotNB.status.SetFields(\
876                        ['','Detector 2-th =%9.3fdeg, dsp =%9.3fA, Q = %6.5fA-1, azm = %7.2fdeg, I = %6d'%(tth,dsp,Q,azm,Int)])
877
878    def OnImPlotKeyPress(event):
879        if self.PatternTree.GetItemText(self.PickId) == 'Masks':
880            Xpos = event.xdata
881            if not Xpos:            #got point out of frame
882                return
883            Ypos = event.ydata
884            if event.key == 's':
885                print 'spot mask @ ',Xpos,Ypos
886                Masks['Points'].append([Xpos,Ypos,1.])
887            elif event.key == 'r':
888                tth = G2img.GetTth(Xpos,Ypos,Data)
889                print 'ring mask @ ',Xpos,Ypos,tth
890                Masks['Rings'].append([tth,0.1])
891            elif event.key == 'a':
892                tth,azm = G2img.GetTthAzm(Xpos,Ypos,Data)
893                azm = int(azm)               
894                print 'arc mask @ ', Xpos,Ypos
895                Masks['Arcs'].append([tth,[azm-5,azm+5],0.1])
896            elif event.key == 'p':
897                self.setPoly = True
898                Masks['Polygons'].append([])
899                self.G2plotNB.status.SetFields(['','Polygon mask active - LB pick next point, RB close polygon'])
900            G2imG.UpdateMasks(self,Masks)
901        else:
902            Xpos = event.xdata
903            if not Xpos:            #got point out of frame
904                return
905            Ypos = event.ydata
906            if event.key == 'c':
907                print 'move center to: ',Xpos,Ypos
908                Data['center'] = [Xpos,Ypos]
909                G2imG.UpdateImageControls(self,Data,Masks)
910        PlotImage(self)
911           
912    def OnImPick(event):
913        if self.PatternTree.GetItemText(self.PickId) not in ['Image Controls','Masks']:
914            return
915        if self.setPoly:
916            polygon = Masks['Polygons'][-1]
917            xpos,ypos = event.mouseevent.xdata,event.mouseevent.ydata
918            if xpos and ypos:                       #point inside image
919                if len(polygon) > 2 and event.mouseevent.button == 3:
920                    x0,y0 = polygon[0]
921                    polygon.append([x0,y0])
922                    self.setPoly = False
923                    self.G2plotNB.status.SetFields(['','Polygon closed - RB drag a vertex to change shape'])
924                else:
925                    self.G2plotNB.status.SetFields(['','New polygon point: %.1f,%.1f'%(xpos,ypos)])
926                    polygon.append([xpos,ypos])
927                G2imG.UpdateMasks(self,Masks)
928        else:
929            if self.itemPicked is not None: return
930            self.itemPicked = event.artist
931            self.mousePicked = event.mouseevent
932       
933    def OnImRelease(event):
934        PickName = self.PatternTree.GetItemText(self.PickId)
935        if PickName not in ['Image Controls','Masks']:
936            return
937        pixelSize = Data['pixelSize']
938        scalex = 1000./pixelSize[0]
939        scaley = 1000./pixelSize[1]
940        pixLimit = Data['pixLimit']
941        if self.itemPicked is None and PickName == 'Image Controls':
942#            sizexy = Data['size']
943            Xpos = event.xdata
944            if not (Xpos and self.ifGetRing):                   #got point out of frame
945                return
946            Ypos = event.ydata
947            if Ypos and not Page.toolbar._active:         #make sure zoom/pan not selected
948                if event.button == 1:
949                    Xpix = Xpos*scalex
950                    Ypix = Ypos*scaley
951                    xpos,ypos,I,J = G2img.ImageLocalMax(self.ImageZ,pixLimit,Xpix,Ypix)
952                    if I and J:
953                        xpos += .5                              #shift to pixel center
954                        ypos += .5
955                        xpos /= scalex                          #convert to mm
956                        ypos /= scaley
957                        Data['ring'].append([xpos,ypos])
958                PlotImage(self)
959            return
960        else:
961            xpos = event.xdata
962            if xpos:                                        #avoid out of frame mouse position
963                ypos = event.ydata
964                if self.ifGetRing:
965                    xypos = [xpos,ypos]
966                    rings = Data['ring']
967                    for ring in rings:
968                        if np.allclose(ring,xypos,.01,0):
969                            rings.remove(ring)                                                                       
970                else:
971                    tth,azm,dsp = G2img.GetTthAzmDsp(xpos,ypos,Data)
972                    itemPicked = str(self.itemPicked)
973                    if 'Line2D' in itemPicked and PickName == 'Image Controls':
974                        if 'line1' in itemPicked:
975                            Data['IOtth'][0] = tth
976                        elif 'line2' in itemPicked:
977                            Data['IOtth'][1] = tth
978                        elif 'line3' in itemPicked and not Data['fullIntegrate']:
979                            Data['LRazimuth'][0] = int(azm)
980                        elif 'line4' in itemPicked and not Data['fullIntegrate']:
981                            Data['LRazimuth'][1] = int(azm)
982                           
983                        if Data['LRazimuth'][1] < Data['LRazimuth'][0]:
984                            Data['LRazimuth'][1] += 360
985                        if  Data['IOtth'][0] > Data['IOtth'][1]:
986                            Data['IOtth'][0],Data['IOtth'][1] = Data['IOtth'][1],Data['IOtth'][0]
987                           
988                        self.InnerTth.SetValue("%8.2f" % (Data['IOtth'][0]))
989                        self.OuterTth.SetValue("%8.2f" % (Data['IOtth'][1]))
990                        self.Lazim.SetValue("%6d" % (Data['LRazimuth'][0]))
991                        self.Razim.SetValue("%6d" % (Data['LRazimuth'][1]))
992                    elif 'Circle' in itemPicked and PickName == 'Masks':
993                        spots = Masks['Points']
994                        newPos = itemPicked.split(')')[0].split('(')[2].split(',')
995                        newPos = np.array([float(newPos[0]),float(newPos[1])])
996                        for spot in spots:
997                            if np.allclose(np.array([spot[:2]]),newPos):
998                                spot[:2] = xpos,ypos
999                        G2imG.UpdateMasks(self,Masks)
1000                    elif 'Line2D' in itemPicked and PickName == 'Masks':
1001                        Obj = self.itemPicked.findobj()
1002                        rings = Masks['Rings']
1003                        arcs = Masks['Arcs']
1004                        polygons = Masks['Polygons']
1005                        for ring in self.ringList:
1006                            if Obj == ring[0]:
1007                                rN = ring[1]
1008                                if ring[2] == 'o':
1009                                    rings[rN][0] = G2img.GetTth(xpos,ypos,Data)-rings[rN][1]/2.
1010                                else:
1011                                    rings[rN][0] = G2img.GetTth(xpos,ypos,Data)+rings[rN][1]/2.
1012                        for arc in self.arcList:
1013                            if Obj == arc[0]:
1014                                aN = arc[1]
1015                                if arc[2] == 'o':
1016                                    arcs[aN][0] = G2img.GetTth(xpos,ypos,Data)-arcs[aN][2]/2
1017                                elif arc[2] == 'i':
1018                                    arcs[aN][0] = G2img.GetTth(xpos,ypos,Data)+arcs[aN][2]/2
1019                                elif arc[2] == 'l':
1020                                    arcs[aN][1][0] = int(G2img.GetAzm(xpos,ypos,Data))
1021                                else:
1022                                    arcs[aN][1][1] = int(G2img.GetAzm(xpos,ypos,Data))
1023                        for poly in self.polyList:
1024                            if Obj == poly[0]:
1025                                ind = self.itemPicked.contains(self.mousePicked)[1]['ind'][0]
1026                                oldPos = np.array([self.mousePicked.xdata,self.mousePicked.ydata])
1027                                pN = poly[1]
1028                                for i,xy in enumerate(polygons[pN]):
1029                                    if np.allclose(np.array([xy]),oldPos,atol=1.0):
1030                                        polygons[pN][i] = xpos,ypos
1031                        G2imG.UpdateMasks(self,Masks)
1032#                    else:                  #keep for future debugging
1033#                        print str(self.itemPicked),event.xdata,event.ydata,event.button
1034                PlotImage(self)
1035            self.itemPicked = None
1036           
1037    try:
1038        plotNum = self.G2plotNB.plotList.index('2D Powder Image')
1039        Page = self.G2plotNB.nb.GetPage(plotNum)
1040        if not newPlot:
1041            Plot = Page.figure.gca()          #get previous powder plot & get limits
1042            xylim = Plot.get_xlim(),Plot.get_ylim()
1043        Page.figure.clf()
1044        Plot = Page.figure.gca()          #get a fresh plot after clf()
1045       
1046    except ValueError,error:
1047        Plot = self.G2plotNB.addMpl('2D Powder Image').gca()
1048        plotNum = self.G2plotNB.plotList.index('2D Powder Image')
1049        Page = self.G2plotNB.nb.GetPage(plotNum)
1050        Page.canvas.mpl_connect('key_press_event', OnImPlotKeyPress)
1051        Page.canvas.mpl_connect('motion_notify_event', OnImMotion)
1052        Page.canvas.mpl_connect('pick_event', OnImPick)
1053        Page.canvas.mpl_connect('button_release_event', OnImRelease)
1054        xylim = []
1055    if not event:                       #event from GUI TextCtrl - don't want focus to change to plot!!!
1056        Page.SetFocus()
1057    Plot.set_title(self.PatternTree.GetItemText(self.Image)[4:])
1058    size,imagefile = self.PatternTree.GetItemPyData(self.Image)
1059    if imagefile != self.oldImagefile:
1060        imagefile = G2IO.CheckImageFile(self,imagefile)
1061        if not imagefile:
1062            self.G2plotNB.Delete('2D Powder Image')
1063            return
1064        self.PatternTree.SetItemPyData(self.Image,[size,imagefile])
1065        self.ImageZ = G2IO.GetImageData(self,imagefile,imageOnly=True)
1066#        print self.ImageZ.shape,self.ImageZ.size,Data['size'] #might be useful debugging line
1067        self.oldImagefile = imagefile
1068
1069    imScale = 1
1070    if len(self.ImageZ) > 1024:
1071        imScale = len(self.ImageZ)/1024
1072    sizexy = Data['size']
1073    pixelSize = Data['pixelSize']
1074    scalex = 1000./pixelSize[0]
1075    scaley = 1000./pixelSize[1]
1076    Xmax = sizexy[0]*pixelSize[0]/1000.
1077    Ymax = sizexy[1]*pixelSize[1]/1000.
1078    xlim = (0,Xmax)
1079    ylim = (Ymax,0)
1080    Imin,Imax = Data['range'][1]
1081    acolor = mpl.cm.get_cmap(Data['color'])
1082    xcent,ycent = Data['center']
1083    Plot.set_xlabel('Image x-axis, mm',fontsize=12)
1084    Plot.set_ylabel('Image y-axis, mm',fontsize=12)
1085    #do threshold mask - "real" mask - others are just bondaries
1086    Zlim = Masks['Thresholds'][1]
1087    wx.BeginBusyCursor()
1088    try:
1089        MA = ma.masked_greater(ma.masked_less(self.ImageZ,Zlim[0]),Zlim[1])
1090        MaskA = ma.getmaskarray(MA)
1091        A = G2img.ImageCompress(MA,imScale)
1092        AM = G2img.ImageCompress(MaskA,imScale)
1093       
1094        ImgM = Plot.imshow(AM,aspect='equal',cmap='Reds',
1095            interpolation='nearest',vmin=0,vmax=2,extent=[0,Xmax,Xmax,0])
1096        Img = Plot.imshow(A,aspect='equal',cmap=acolor,
1097            interpolation='nearest',vmin=Imin,vmax=Imax,extent=[0,Xmax,Ymax,0])
1098        if self.setPoly:
1099            Img.set_picker(True)
1100   
1101        Plot.plot(xcent,ycent,'x')
1102        if Data['showLines']:
1103            LRAzim = Data['LRazimuth']                  #NB: integers
1104            IOtth = Data['IOtth']
1105            wave = Data['wavelength']
1106            dspI = wave/(2.0*sind(IOtth[0]/2.0))
1107            ellI = G2img.GetEllipse(dspI,Data)           #=False if dsp didn't yield an ellipse (ugh! a parabola or a hyperbola)
1108            dspO = wave/(2.0*sind(IOtth[1]/2.0))
1109            ellO = G2img.GetEllipse(dspO,Data)           #Ditto & more likely for outer ellipse
1110            if Data['fullIntegrate']:
1111                Azm = np.array(range(0,361))
1112            else:
1113                Azm = np.array(range(LRAzim[0],LRAzim[1]+1))
1114            if ellI:
1115                xyI = []
1116                for azm in Azm:
1117                    xyI.append(G2img.GetDetectorXY(dspI,azm,Data))
1118                xyI = np.array(xyI)
1119                arcxI,arcyI = xyI.T
1120                Plot.plot(arcxI,arcyI,picker=3)
1121            if ellO:
1122                xyO = []
1123                for azm in Azm:
1124                    xyO.append(G2img.GetDetectorXY(dspO,azm,Data))
1125                xyO = np.array(xyO)
1126                arcxO,arcyO = xyO.T
1127                Plot.plot(arcxO,arcyO,picker=3)
1128            if ellO and ellI and not Data['fullIntegrate']:
1129                Plot.plot([arcxI[0],arcxO[0]],[arcyI[0],arcyO[0]],picker=3)
1130                Plot.plot([arcxI[-1],arcxO[-1]],[arcyI[-1],arcyO[-1]],picker=3)
1131        for xring,yring in Data['ring']:
1132            Plot.plot(xring,yring,'r+',picker=3)
1133        if Data['setRings']:
1134            rings = np.concatenate((Data['rings']),axis=0)
1135            for xring,yring,dsp in rings:
1136                Plot.plot(xring,yring,'r+')           
1137        for ellipse in Data['ellipses']:
1138            cent,phi,[width,height],col = ellipse
1139            Plot.add_artist(Ellipse([cent[0],cent[1]],2*width,2*height,phi,ec=col,fc='none'))
1140            Plot.text(cent[0],cent[1],'+',color=col,ha='center',va='center')
1141        #masks - mask lines numbered after integration limit lines
1142        spots = Masks['Points']
1143        rings = Masks['Rings']
1144        arcs = Masks['Arcs']
1145        polygons = Masks['Polygons']
1146        for x,y,d in spots:
1147            Plot.add_artist(Circle((x,y),radius=d/2,fc='none',ec='r',picker=3))
1148        self.ringList = []
1149        for iring,(tth,thick) in enumerate(rings):
1150            wave = Data['wavelength']
1151            x1,y1 = np.hsplit(np.array(G2img.makeIdealRing(G2img.GetEllipse(Dsp(tth+thick/2.,wave),Data))),2)
1152            x2,y2 = np.hsplit(np.array(G2img.makeIdealRing(G2img.GetEllipse(Dsp(tth-thick/2.,wave),Data))),2)
1153            self.ringList.append([Plot.plot(x1,y1,'r',picker=3),iring,'o'])           
1154            self.ringList.append([Plot.plot(x2,y2,'r',picker=3),iring,'i'])
1155        self.arcList = []
1156        for iarc,(tth,azm,thick) in enumerate(arcs):           
1157            wave = Data['wavelength']
1158            x1,y1 = np.hsplit(np.array(G2img.makeIdealRing(G2img.GetEllipse(Dsp(tth+thick/2.,wave),Data),azm)),2)
1159            x2,y2 = np.hsplit(np.array(G2img.makeIdealRing(G2img.GetEllipse(Dsp(max(0.01,tth-thick/2.),wave),Data),azm)),2)
1160            self.arcList.append([Plot.plot(x1,y1,'r',picker=3),iarc,'o'])           
1161            self.arcList.append([Plot.plot(x2,y2,'r',picker=3),iarc,'i'])
1162            self.arcList.append([Plot.plot([x1[0],x2[0]],[y1[0],y2[0]],'r',picker=3),iarc,'l'])
1163            self.arcList.append([Plot.plot([x1[-1],x2[-1]],[y1[-1],y2[-1]],'r',picker=3),iarc,'u'])
1164        self.polyList = []
1165        for ipoly,polygon in enumerate(polygons):
1166            x,y = np.hsplit(np.array(polygon),2)
1167            self.polyList.append([Plot.plot(x,y,'r+',picker=10),ipoly])
1168            Plot.plot(x,y,'r')           
1169        colorBar = Page.figure.colorbar(Img)
1170        Plot.set_xlim(xlim)
1171        Plot.set_ylim(ylim)
1172        if not newPlot and xylim:
1173            Page.toolbar.push_current()
1174            Plot.set_xlim(xylim[0])
1175            Plot.set_ylim(xylim[1])
1176            xylim = []
1177            Page.toolbar.push_current()
1178            Page.toolbar.draw()
1179        else:
1180            Page.canvas.draw()
1181    finally:
1182        wx.EndBusyCursor()
1183       
1184def PlotIntegration(self,newPlot=False,event=None):
1185           
1186    def OnMotion(event):
1187        Page.canvas.SetToolTipString('')
1188        Page.canvas.SetCursor(wx.CROSS_CURSOR)
1189        azm = event.ydata
1190        tth = event.xdata
1191        if azm and tth:
1192            self.G2plotNB.status.SetFields(\
1193                ['','Detector 2-th =%9.3fdeg, azm = %7.2fdeg'%(tth,azm)])
1194                               
1195    try:
1196        plotNum = self.G2plotNB.plotList.index('2D Integration')
1197        Page = self.G2plotNB.nb.GetPage(plotNum)
1198        if not newPlot:
1199            Plot = Page.figure.gca()          #get previous plot & get limits
1200            xylim = Plot.get_xlim(),Plot.get_ylim()
1201        Page.figure.clf()
1202        Plot = Page.figure.gca()          #get a fresh plot after clf()
1203       
1204    except ValueError,error:
1205        Plot = self.G2plotNB.addMpl('2D Integration').gca()
1206        plotNum = self.G2plotNB.plotList.index('2D Integration')
1207        Page = self.G2plotNB.nb.GetPage(plotNum)
1208        Page.canvas.mpl_connect('motion_notify_event', OnMotion)
1209        Page.views = False
1210        view = False
1211    if not event:
1212        Page.SetFocus()
1213       
1214    Data = self.PatternTree.GetItemPyData(
1215        G2gd.GetPatternTreeItemId(self,self.Image, 'Image Controls'))
1216    image = self.Integrate[0]
1217    xsc = self.Integrate[1]
1218    ysc = self.Integrate[2]
1219    Imin,Imax = Data['range'][1]
1220    acolor = mpl.cm.get_cmap(Data['color'])
1221    Plot.set_title(self.PatternTree.GetItemText(self.Image)[4:])
1222    Plot.set_ylabel('azimuth',fontsize=12)
1223    Plot.set_xlabel('2-theta',fontsize=12)
1224    Img = Plot.imshow(image,cmap=acolor,vmin=Imin,vmax=Imax,interpolation='nearest', \
1225        extent=[ysc[0],ysc[-1],xsc[-1],xsc[0]],aspect='auto')
1226    colorBar = Page.figure.colorbar(Img)
1227    if Data['setRings']:
1228        rings = np.concatenate((Data['rings']),axis=0)
1229        for xring,yring,dsp in rings:
1230            x,y = G2img.GetTthAzm(xring,yring,Data)
1231            Plot.plot(x,y,'r+')
1232    if Data['ellipses']:           
1233        for ellipse in Data['ellipses']:
1234            ring = np.array(G2img.makeIdealRing(ellipse[:3])) #skip color
1235            x,y = np.hsplit(ring,2)
1236            tth,azm = G2img.GetTthAzm(x,y,Data)
1237            Plot.plot(tth,azm,'b,')
1238    if not newPlot:
1239        Page.toolbar.push_current()
1240        Plot.set_xlim(xylim[0])
1241        Plot.set_ylim(xylim[1])
1242        xylim = []
1243        Page.toolbar.push_current()
1244        Page.toolbar.draw()
1245    else:
1246        Page.canvas.draw()
1247               
1248def PlotTRImage(self,tax,tay,taz,newPlot=False):
1249    #a test plot routine - not normally used
1250           
1251    def OnMotion(event):
1252        Page.canvas.SetToolTipString('')
1253        Page.canvas.SetCursor(wx.CROSS_CURSOR)
1254        azm = event.xdata
1255        tth = event.ydata
1256        if azm and tth:
1257            self.G2plotNB.status.SetFields(\
1258                ['','Detector 2-th =%9.3fdeg, azm = %7.2fdeg'%(tth,azm)])
1259                               
1260    try:
1261        plotNum = self.G2plotNB.plotList.index('2D Transformed Powder Image')
1262        Page = self.G2plotNB.nb.GetPage(plotNum)
1263        if not newPlot:
1264            Plot = Page.figure.gca()          #get previous plot & get limits
1265            xylim = Plot.get_xlim(),Plot.get_ylim()
1266        Page.figure.clf()
1267        Plot = Page.figure.gca()          #get a fresh plot after clf()
1268       
1269    except ValueError,error:
1270        Plot = self.G2plotNB.addMpl('2D Transformed Powder Image').gca()
1271        plotNum = self.G2plotNB.plotList.index('2D Transformed Powder Image')
1272        Page = self.G2plotNB.nb.GetPage(plotNum)
1273        Page.canvas.mpl_connect('motion_notify_event', OnMotion)
1274        Page.views = False
1275        view = False
1276    Page.SetFocus()
1277       
1278    Data = self.PatternTree.GetItemPyData(
1279        G2gd.GetPatternTreeItemId(self,self.Image, 'Image Controls'))
1280    Imin,Imax = Data['range'][1]
1281    step = (Imax-Imin)/5.
1282    V = np.arange(Imin,Imax,step)
1283    acolor = mpl.cm.get_cmap(Data['color'])
1284    Plot.set_title(self.PatternTree.GetItemText(self.Image)[4:])
1285    Plot.set_xlabel('azimuth',fontsize=12)
1286    Plot.set_ylabel('2-theta',fontsize=12)
1287    Plot.contour(tax,tay,taz,V,cmap=acolor)
1288    if Data['showLines']:
1289        IOtth = Data['IOtth']
1290        if Data['fullIntegrate']:
1291            LRAzim = [-180,180]
1292        else:
1293            LRAzim = Data['LRazimuth']                  #NB: integers
1294        Plot.plot([LRAzim[0],LRAzim[1]],[IOtth[0],IOtth[0]],picker=True)
1295        Plot.plot([LRAzim[0],LRAzim[1]],[IOtth[1],IOtth[1]],picker=True)
1296        if not Data['fullIntegrate']:
1297            Plot.plot([LRAzim[0],LRAzim[0]],[IOtth[0],IOtth[1]],picker=True)
1298            Plot.plot([LRAzim[1],LRAzim[1]],[IOtth[0],IOtth[1]],picker=True)
1299    if Data['setRings']:
1300        rings = np.concatenate((Data['rings']),axis=0)
1301        for xring,yring,dsp in rings:
1302            x,y = G2img.GetTthAzm(xring,yring,Data)
1303            Plot.plot(y,x,'r+')           
1304    if Data['ellipses']:           
1305        for ellipse in Data['ellipses']:
1306            ring = np.array(G2img.makeIdealRing(ellipse[:3])) #skip color
1307            x,y = np.hsplit(ring,2)
1308            tth,azm = G2img.GetTthAzm(x,y,Data)
1309            Plot.plot(azm,tth,'b,')
1310    if not newPlot:
1311        Page.toolbar.push_current()
1312        Plot.set_xlim(xylim[0])
1313        Plot.set_ylim(xylim[1])
1314        xylim = []
1315        Page.toolbar.push_current()
1316        Page.toolbar.draw()
1317    else:
1318        Page.canvas.draw()
1319       
1320def PlotStructure(self,data):
1321    generalData = data['General']
1322    cell = generalData['Cell'][1:7]
1323    Amat,Bmat = G2lat.cell2AB(cell)         #Amat - crystal to cartesian, Bmat - inverse
1324    A4mat = np.concatenate((np.concatenate((Amat,[[0],[0],[0]]),axis=1),[[0,0,0,1],]),axis=0)
1325    B4mat = np.concatenate((np.concatenate((Bmat,[[0],[0],[0]]),axis=1),[[0,0,0,1],]),axis=0)
1326    Mydir = generalData['Mydir']
1327    atomData = data['Atoms']
1328    drawingData = data['Drawing']
1329    drawAtoms = drawingData['Atoms']
1330    cx,ct,cs = drawingData['atomPtrs']
1331    Wt = [255,255,255]
1332    Rd = [255,0,0]
1333    Gr = [0,255,0]
1334    Bl = [0,0,255]
1335    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]])
1336    uEdges = np.array([
1337        [uBox[0],uBox[1]],[uBox[0],uBox[3]],[uBox[0],uBox[4]],[uBox[1],uBox[2]], 
1338        [uBox[2],uBox[3]],[uBox[1],uBox[5]],[uBox[2],uBox[6]],[uBox[3],uBox[7]], 
1339        [uBox[4],uBox[5]],[uBox[5],uBox[6]],[uBox[6],uBox[7]],[uBox[7],uBox[4]]])
1340    uColors = [Rd,Gr,Bl,Wt, Wt,Wt,Wt,Wt, Wt,Wt,Wt,Wt]
1341    altDown = False
1342    shiftDown = False
1343    ctrlDown = False
1344   
1345    def OnKeyBox(event):
1346        import Image
1347        Draw()                          #make sure plot is fresh!!
1348        mode = cb.GetValue()
1349        Fname = Mydir+'\\'+generalData['Name']+'.'+mode
1350        size = Page.canvas.GetSize()
1351        glPixelStorei(GL_UNPACK_ALIGNMENT, 1)
1352        if mode in ['jpeg',]:
1353            Pix = glReadPixels(0,0,size[0],size[1],GL_RGBA, GL_UNSIGNED_BYTE)
1354            im = Image.new("RGBA", (size[0],size[1]))
1355        else:
1356            Pix = glReadPixels(0,0,size[0],size[1],GL_RGB, GL_UNSIGNED_BYTE)
1357            im = Image.new("RGB", (size[0],size[1]))
1358        im.fromstring(Pix)
1359        im.save(Fname,mode)
1360        cb.SetValue(' Save as:')
1361        self.G2plotNB.status.SetStatusText('Drawing saved to: '+Fname,1)
1362   
1363    def GetTruePosition(xy):
1364        View = glGetIntegerv(GL_VIEWPORT)
1365        Proj = glGetDoublev(GL_PROJECTION_MATRIX)
1366        Model = glGetDoublev(GL_MODELVIEW_MATRIX)
1367        Zmax = 1.
1368        for i,atom in enumerate(drawAtoms):
1369            x,y,z = atom[cx:cx+3]
1370            X,Y,Z = gluProject(x,y,z,Model,Proj,View)
1371            XY = [int(X),int(View[3]-Y)]
1372            if np.allclose(xy,XY,atol=10) and Z < Zmax:
1373                Zmax = Z
1374                SetSelectedAtoms(i)
1375                   
1376    def OnMouseDown(event):
1377        xy = event.GetPosition()
1378        if event.ShiftDown():
1379            GetTruePosition(xy)
1380        else:
1381            drawingData['Rotation'][3] = xy
1382        Draw()
1383       
1384    def OnMouseMove(event):
1385        newxy = event.GetPosition()
1386        if event.ControlDown() and drawingData['showABC']:
1387            if event.LeftIsDown():
1388                ctrlDown = True
1389                SetTestRot(newxy)
1390            elif event.RightIsDown():
1391                SetTestPos(newxy)
1392            elif event.MiddleIsDown():
1393                SetTestRotZ(newxy)
1394               
1395               
1396        if event.Dragging() and not event.ControlDown():
1397            if event.LeftIsDown():
1398                SetRotation(newxy)
1399            elif event.RightIsDown():
1400                SetTranslation(newxy)
1401            elif event.MiddleIsDown():
1402                SetRotationZ(newxy)
1403        Draw()
1404       
1405    def OnMouseWheel(event):
1406        drawingData['cameraPos'] += event.GetWheelRotation()/24
1407        drawingData['cameraPos'] = max(10,min(500,drawingData['cameraPos']))
1408        page = self.dataDisplay.GetSelection()
1409        if page:
1410            if self.dataDisplay.GetPageText(page) == 'Draw Options':
1411                panel = self.dataDisplay.GetPage(page).GetChildren()[0].GetChildren()
1412                names = [child.GetName() for child in panel]
1413                panel[names.index('cameraPos')].SetLabel('Camera Position: '+'%.2f'%(drawingData['cameraPos']))
1414                panel[names.index('cameraSlider')].SetValue(drawingData['cameraPos'])
1415        Draw()
1416           
1417    def SetViewPointText(VP):
1418        page = self.dataDisplay.GetSelection()
1419        if page:
1420            if self.dataDisplay.GetPageText(page) == 'Draw Options':
1421                panel = self.dataDisplay.GetPage(page).GetChildren()[0].GetChildren()
1422                names = [child.GetName() for child in panel]
1423                panel[names.index('viewPoint')].SetValue('%.3f, %.3f, %.3f'%(VP[0],VP[1],VP[2]))
1424           
1425    def ClearSelectedAtoms():
1426        page = self.dataDisplay.GetSelection()
1427        if page:
1428            if self.dataDisplay.GetPageText(page) == 'Draw Atoms':
1429                self.dataDisplay.GetPage(page).ClearSelection()      #this is the Atoms grid in Draw Atoms
1430            elif self.dataDisplay.GetPageText(page) == 'Atoms':
1431                self.dataDisplay.GetPage(page).ClearSelection()      #this is the Atoms grid in Atoms
1432                   
1433    def SetSelectedAtoms(ind):
1434        page = self.dataDisplay.GetSelection()
1435        if page:
1436            if self.dataDisplay.GetPageText(page) == 'Draw Atoms':
1437                self.dataDisplay.GetPage(page).SelectRow(ind)      #this is the Atoms grid in Draw Atoms
1438            elif self.dataDisplay.GetPageText(page) == 'Atoms':
1439                Id = drawAtoms[ind][-2]
1440                for i,atom in enumerate(atomData):
1441                    if atom[-1] == Id:
1442                        self.dataDisplay.GetPage(page).SelectRow(i)      #this is the Atoms grid in Atoms
1443                 
1444    def GetSelectedAtoms():
1445        page = self.dataDisplay.GetSelection()
1446        Ind = []
1447        if page:
1448            if self.dataDisplay.GetPageText(page) == 'Draw Atoms':
1449                Ind = self.dataDisplay.GetPage(page).GetSelectedRows()      #this is the Atoms grid in Draw Atoms
1450            elif self.dataDisplay.GetPageText(page) == 'Atoms':
1451                Ind = self.dataDisplay.GetPage(page).GetSelectedRows()      #this is the Atoms grid in Atoms
1452        return Ind
1453                                       
1454    def OnKey(event):           #on key UP!!
1455        keyCode = event.GetKeyCode()
1456        if keyCode > 255:
1457            keyCode = 0
1458        key,xyz = chr(keyCode),event.GetPosition()
1459        indx = drawingData['selectedAtoms']
1460        if key in ['c','C']:
1461            drawingData['viewPoint'] = [[.5,.5,.5],[0,0]]
1462            drawingData['testPos'] = [[-.1,-.1,-.1],[0.0,0.0,0.0],[0,0]]
1463            drawingData['Rotation'] = [0.0,0.0,0.0,[]]
1464            SetViewPointText(drawingData['viewPoint'][0])
1465        elif key in ['n','N']:
1466            drawAtoms = drawingData['Atoms']
1467            pI = drawingData['viewPoint'][1]
1468            if indx:
1469                pI[0] = indx[pI[1]]
1470                Tx,Ty,Tz = drawAtoms[pI[0]][cx:cx+3]
1471                pI[1] += 1
1472                if pI[1] >= len(indx):
1473                    pI[1] = 0
1474            else:
1475                Tx,Ty,Tz = drawAtoms[pI[0]][cx:cx+3]               
1476                pI[0] += 1
1477                if pI[0] >= len(drawAtoms):
1478                    pI[0] = 0
1479            drawingData['viewPoint'] = [[Tx,Ty,Tz],pI]
1480            SetViewPointText(drawingData['viewPoint'][0])
1481               
1482        elif key in ['p','P']:
1483            drawAtoms = drawingData['Atoms']
1484            pI = drawingData['viewPoint'][1]
1485            if indx:
1486                pI[0] = indx[pI[1]]
1487                Tx,Ty,Tz = drawAtoms[pI[0]][cx:cx+3]
1488                pI[1] -= 1
1489                if pI[1] < 0:
1490                    pI[1] = len(indx)-1
1491            else:
1492                Tx,Ty,Tz = drawAtoms[pI[0]][cx:cx+3]               
1493                pI[0] -= 1
1494                if pI[0] < 0:
1495                    pI[0] = len(drawAtoms)-1
1496            drawingData['viewPoint'] = [[Tx,Ty,Tz],pI]
1497            SetViewPointText(drawingData['viewPoint'][0])           
1498        Draw()
1499           
1500    def SetBackground():
1501        R,G,B,A = Page.camera['backColor']
1502        glClearColor(R,G,B,A)
1503        glClear(GL_COLOR_BUFFER_BIT | GL_DEPTH_BUFFER_BIT)
1504       
1505    def SetLights():
1506        glEnable(GL_DEPTH_TEST)
1507        glShadeModel(GL_SMOOTH)
1508        glEnable(GL_LIGHTING)
1509        glEnable(GL_LIGHT0)
1510        glLightModeli(GL_LIGHT_MODEL_TWO_SIDE,0)
1511        glLightfv(GL_LIGHT0,GL_AMBIENT,[1,1,1,.8])
1512        glLightfv(GL_LIGHT0,GL_DIFFUSE,[1,1,1,1])
1513       
1514    def SetTranslation(newxy):
1515        Tx,Ty,Tz = drawingData['viewPoint'][0]
1516        anglex,angley,anglez,oldxy = drawingData['Rotation']
1517        Rx = G2lat.rotdMat(anglex,0)
1518        Ry = G2lat.rotdMat(angley,1)
1519        Rz = G2lat.rotdMat(anglez,2)
1520        dxy = list(newxy-oldxy)+[0,]
1521        dxy = np.inner(Bmat,np.inner(Rz,np.inner(Ry,np.inner(Rx,dxy))))
1522        Tx -= dxy[0]*0.01
1523        Ty += dxy[1]*0.01
1524        Tz -= dxy[2]*0.01
1525        drawingData['Rotation'][3] = newxy
1526        drawingData['viewPoint'][0] =  Tx,Ty,Tz
1527        SetViewPointText([Tx,Ty,Tz])
1528       
1529    def SetTestPos(newxy):
1530        Tx,Ty,Tz = drawingData['testPos'][0]
1531        anglex,angley,anglez,oldxy = drawingData['Rotation']
1532        Rx = G2lat.rotdMat(anglex,0)
1533        Ry = G2lat.rotdMat(angley,1)
1534        Rz = G2lat.rotdMat(anglez,2)
1535        dxy = list(newxy-oldxy)+[0,]
1536        dxy = np.inner(Rz,np.inner(Ry,np.inner(Rx,dxy)))
1537        Tx += dxy[0]*0.001
1538        Ty -= dxy[1]*0.001
1539        Tz += dxy[2]*0.001
1540        drawingData['Rotation'][3] = newxy
1541        drawingData['testPos'][0] =  Tx,Ty,Tz
1542       
1543    def SetTestRot(newxy):
1544        Txyz = np.array(drawingData['testPos'][0])
1545        oldxy = drawingData['testPos'][2]
1546        Ax,Ay,Az = drawingData['testPos'][1]
1547        Vxyz = np.array(drawingData['viewPoint'][0])
1548        Dxyz = np.inner(Amat,Txyz-Vxyz)
1549        dxy = list(newxy-oldxy)+[0,]
1550        Ax += dxy[1]*0.01
1551        Ay += dxy[0]*0.01
1552        Rx = G2lat.rotdMat(Ax,0)
1553        Ry = G2lat.rotdMat(Ay,1)
1554        Dxyz = np.inner(Ry,np.inner(Rx,Dxyz))       
1555        Dxyz = np.inner(Bmat,Dxyz)+Vxyz
1556        drawingData['testPos'][1] = [Ax,Ay,Az]
1557        drawingData['testPos'][2] = newxy
1558        drawingData['testPos'][0] = Dxyz
1559       
1560    def SetTestRotZ(newxy):
1561        Txyz = np.array(drawingData['testPos'][0])
1562        oldxy = drawingData['testPos'][2]
1563        Ax,Ay,Az = drawingData['testPos'][1]
1564        Vxyz = np.array(drawingData['viewPoint'][0])
1565        Dxyz = np.inner(Amat,Txyz-Vxyz)       
1566        dxy = list(newxy-oldxy)+[0,]
1567        Az += (dxy[0]+dxy[1])*.01
1568        Rz = G2lat.rotdMat(Az,2)
1569        Dxyz = np.inner(Rz,Dxyz)       
1570        Dxyz = np.inner(Bmat,Dxyz)+Vxyz
1571        drawingData['testPos'][1] = [Ax,Ay,Az]
1572        drawingData['testPos'][2] = newxy
1573        drawingData['testPos'][0] = Dxyz
1574                             
1575    def SetRotation(newxy):       
1576        anglex,angley,anglez,oldxy = drawingData['Rotation']
1577        dxy = newxy-oldxy
1578        anglex += dxy[1]*.25
1579        angley += dxy[0]*.25
1580        oldxy = newxy
1581        drawingData['Rotation'] = [anglex,angley,anglez,oldxy]
1582       
1583    def SetRotationZ(newxy):                       
1584        anglex,angley,anglez,oldxy = drawingData['Rotation']
1585        dxy = newxy-oldxy
1586        anglez += (dxy[0]+dxy[1])*.25
1587        oldxy = newxy
1588        drawingData['Rotation'] = [anglex,angley,anglez,oldxy]
1589       
1590    def RenderBox():
1591        glEnable(GL_COLOR_MATERIAL)
1592        glLineWidth(2)
1593        glEnable(GL_BLEND)
1594        glBlendFunc(GL_SRC_ALPHA,GL_ONE_MINUS_SRC_ALPHA)
1595        glEnable(GL_LINE_SMOOTH)
1596        glBegin(GL_LINES)
1597        for line,color in zip(uEdges,uColors):
1598            glColor3ubv(color)
1599            glVertex3fv(line[0])
1600            glVertex3fv(line[1])
1601        glEnd()
1602        glColor4ubv([0,0,0,0])
1603        glDisable(GL_LINE_SMOOTH)
1604        glDisable(GL_BLEND)
1605        glDisable(GL_COLOR_MATERIAL)
1606       
1607    def RenderUnitVectors(x,y,z):
1608        xyz = np.array([x,y,z])
1609        glEnable(GL_COLOR_MATERIAL)
1610        glLineWidth(1)
1611        glPushMatrix()
1612        glTranslate(x,y,z)
1613        glScalef(1/cell[0],1/cell[1],1/cell[2])
1614        glBegin(GL_LINES)
1615        for line,color in zip(uEdges,uColors)[:3]:
1616            glColor3ubv(color)
1617            glVertex3fv(line[0])
1618            glVertex3fv(line[1])
1619        glEnd()
1620        glPopMatrix()
1621        glColor4ubv([0,0,0,0])
1622        glDisable(GL_COLOR_MATERIAL)
1623               
1624    def RenderSphere(x,y,z,radius,color):
1625        glMaterialfv(GL_FRONT_AND_BACK,GL_DIFFUSE,color)
1626        glPushMatrix()
1627        glTranslate(x,y,z)
1628        glMultMatrixf(B4mat.T)
1629        q = gluNewQuadric()
1630        gluSphere(q,radius,20,10)
1631        glPopMatrix()
1632       
1633    def RenderEllipsoid(x,y,z,ellipseProb,E,R4,color):
1634        s1,s2,s3 = E
1635        glMaterialfv(GL_FRONT_AND_BACK,GL_DIFFUSE,color)
1636        glPushMatrix()
1637        glTranslate(x,y,z)
1638        glMultMatrixf(B4mat.T)
1639        glMultMatrixf(R4.T)
1640        glEnable(GL_NORMALIZE)
1641        glScale(s1,s2,s3)
1642        q = gluNewQuadric()
1643        gluSphere(q,ellipseProb,20,10)
1644        glDisable(GL_NORMALIZE)
1645        glPopMatrix()
1646       
1647    def RenderBonds(x,y,z,Bonds,radius,color,slice=20):
1648        glMaterialfv(GL_FRONT_AND_BACK,GL_DIFFUSE,color)
1649        glPushMatrix()
1650        glTranslate(x,y,z)
1651        glMultMatrixf(B4mat.T)
1652        for bond in Bonds:
1653            glPushMatrix()
1654            Dx = np.inner(Amat,bond)
1655            Z = np.sqrt(np.sum(Dx**2))
1656            azm = atan2d(-Dx[1],-Dx[0])
1657            phi = acosd(Dx[2]/Z)
1658            glRotate(-azm,0,0,1)
1659            glRotate(phi,1,0,0)
1660            q = gluNewQuadric()
1661            gluCylinder(q,radius,radius,Z,slice,2)
1662            glPopMatrix()           
1663        glPopMatrix()
1664               
1665    def RenderLines(x,y,z,Bonds,color):
1666        xyz = np.array([x,y,z])
1667        glEnable(GL_COLOR_MATERIAL)
1668        glLineWidth(1)
1669        glColor3fv(color)
1670        glPushMatrix()
1671        glBegin(GL_LINES)
1672        for bond in Bonds:
1673            glVertex3fv(xyz)
1674            glVertex3fv(xyz+bond)
1675        glEnd()
1676        glColor4ubv([0,0,0,0])
1677        glPopMatrix()
1678        glDisable(GL_COLOR_MATERIAL)
1679       
1680    def RenderPolyhedra(x,y,z,Faces,color):
1681        glPushMatrix()
1682        glTranslate(x,y,z)
1683        glMaterialfv(GL_FRONT_AND_BACK,GL_DIFFUSE,color)
1684        glShadeModel(GL_SMOOTH)
1685        glMultMatrixf(B4mat.T)
1686        for face,norm in Faces:
1687            glPolygonMode(GL_FRONT_AND_BACK,GL_FILL)
1688            glFrontFace(GL_CW)
1689            glNormal3fv(norm)
1690            glBegin(GL_TRIANGLES)
1691            for vert in face:
1692                glVertex3fv(vert)
1693            glEnd()
1694        glPopMatrix()
1695       
1696    def RenderBackbone(Backbone,BackboneColor,radius):
1697        glPushMatrix()
1698        glMultMatrixf(B4mat.T)
1699        glEnable(GL_COLOR_MATERIAL)
1700        glShadeModel(GL_SMOOTH)
1701        gleSetJoinStyle(TUBE_NORM_EDGE | TUBE_JN_ANGLE | TUBE_JN_CAP)
1702        glePolyCylinder(Backbone,BackboneColor,radius)
1703        glPopMatrix()       
1704        glDisable(GL_COLOR_MATERIAL)
1705       
1706    def RenderLabel(x,y,z,label,r):       
1707        glPushMatrix()
1708        glTranslate(x,y,z)
1709        glMultMatrixf(B4mat.T)
1710        glDisable(GL_LIGHTING)
1711        glColor3f(0,1.,0)
1712        glRasterPos3f(r,r,r)
1713        for c in list(label):
1714            glutBitmapCharacter(GLUT_BITMAP_8_BY_13,ord(c))
1715        glEnable(GL_LIGHTING)
1716        glPopMatrix()
1717                           
1718    def Draw():
1719        import numpy.linalg as nl
1720        Ind = GetSelectedAtoms()
1721        VS = np.array(Page.canvas.GetSize())
1722        aspect = float(VS[0])/float(VS[1])
1723        cPos = drawingData['cameraPos']
1724        Zclip = drawingData['Zclip']*cPos/200.
1725        anglex,angley,anglez = drawingData['Rotation'][:3]
1726        Tx,Ty,Tz = drawingData['viewPoint'][0]
1727        cx,ct,cs = drawingData['atomPtrs']
1728        bondR = drawingData['bondRadius']
1729        G,g = G2lat.cell2Gmat(cell)
1730        GS = G
1731        GS[0][1] = GS[1][0] = math.sqrt(GS[0][0]*GS[1][1])
1732        GS[0][2] = GS[2][0] = math.sqrt(GS[0][0]*GS[2][2])
1733        GS[1][2] = GS[2][1] = math.sqrt(GS[1][1]*GS[2][2])
1734        ellipseProb = G2lat.criticalEllipse(drawingData['ellipseProb']/100.)
1735       
1736        SetBackground()
1737        glInitNames()
1738        glPushName(0)
1739       
1740        glMatrixMode(GL_PROJECTION)
1741        glLoadIdentity()
1742        glViewport(0,0,VS[0],VS[1])
1743        gluPerspective(20.,aspect,cPos-Zclip,cPos+Zclip)
1744        gluLookAt(0,0,cPos,0,0,0,0,1,0)
1745        SetLights()           
1746           
1747        glMatrixMode(GL_MODELVIEW)
1748        glLoadIdentity()
1749        glRotate(anglez,0,0,1)
1750        glRotate(anglex,cosd(anglez),-sind(anglez),0)
1751        glRotate(angley,sind(anglez),cosd(anglez),0)
1752        glMultMatrixf(A4mat.T)
1753        glTranslate(-Tx,-Ty,-Tz)
1754        if drawingData['unitCellBox']:
1755            RenderBox()
1756        if drawingData['showABC']:
1757#            try:            #temporary fix - not needed further?
1758#                x,y,z = drawingData['testPos'][0]
1759#            except TypeError:
1760#                x,y,z = drawingData['testPos']
1761            x,y,z = drawingData['testPos'][0]
1762            if altDown:
1763                self.G2plotNB.status.SetStatusText('moving test point %.4f,%.4f,%.4f'%(x,y,z),1)
1764            else:
1765                self.G2plotNB.status.SetStatusText('test point %.4f,%.4f,%.4f'%(x,y,z),1)           
1766            RenderUnitVectors(x,y,z)
1767        Backbone = []
1768        BackboneColor = []
1769        time0 = time.time()
1770        for iat,atom in enumerate(drawingData['Atoms']):
1771            x,y,z = atom[cx:cx+3]
1772            Bonds = atom[-2]
1773            Faces = atom[-1]
1774            try:
1775                atNum = generalData['AtomTypes'].index(atom[ct])
1776            except ValueError:
1777                atNum = -1
1778            CL = atom[cs+2]
1779            color = np.array(CL)/255.
1780            if iat in Ind:
1781                color = np.array(Gr)/255.
1782            radius = 0.5
1783            if atom[cs] != '':
1784                glLoadName(atom[-3])                   
1785            if 'balls' in atom[cs]:
1786                vdwScale = drawingData['vdwScale']
1787                ballScale = drawingData['ballScale']
1788                if atNum < 0:
1789                    radius = 0.2
1790                elif 'H' == atom[ct]:
1791                    if drawingData['showHydrogen']:
1792                        if 'vdW' in atom[cs] and atNum >= 0:
1793                            radius = vdwScale*generalData['vdWRadii'][atNum]
1794                        else:
1795                            radius = ballScale*drawingData['sizeH']
1796                    else:
1797                        radius = 0.0
1798                else:
1799                    if 'vdW' in atom[cs]:
1800                        radius = vdwScale*generalData['vdWRadii'][atNum]
1801                    else:
1802                        radius = ballScale*generalData['BondRadii'][atNum]
1803                RenderSphere(x,y,z,radius,color)
1804                if 'sticks' in atom[cs]:
1805                    RenderBonds(x,y,z,Bonds,bondR,color)
1806            elif 'ellipsoids' in atom[cs]:
1807                RenderBonds(x,y,z,Bonds,bondR,color)
1808                if atom[cs+3] == 'A':                   
1809                    Uij = atom[cs+5:cs+11]
1810                    U = np.multiply(G2spc.Uij2U(Uij),GS)
1811                    U = np.inner(Amat,np.inner(U,Amat).T)
1812                    E,R = nl.eigh(U)
1813                    R4 = np.concatenate((np.concatenate((R,[[0],[0],[0]]),axis=1),[[0,0,0,1],]),axis=0)
1814                    E = np.sqrt(E)
1815                    if atom[ct] == 'H' and not drawingData['showHydrogen']:
1816                        pass
1817                    else:
1818                        RenderEllipsoid(x,y,z,ellipseProb,E,R4,color)                   
1819                else:
1820                    if atom[ct] == 'H' and not drawingData['showHydrogen']:
1821                        pass
1822                    else:
1823                        radius = ellipseProb*math.sqrt(abs(atom[cs+4]))
1824                        RenderSphere(x,y,z,radius,color)
1825            elif 'lines' in atom[cs]:
1826                radius = 0.1
1827                RenderLines(x,y,z,Bonds,color)
1828#                RenderBonds(x,y,z,Bonds,0.05,color,6)
1829            elif atom[cs] == 'sticks':
1830                radius = 0.1
1831                RenderBonds(x,y,z,Bonds,bondR,color)
1832            elif atom[cs] == 'polyhedra':
1833                RenderPolyhedra(x,y,z,Faces,color)
1834            elif atom[cs] == 'backbone':
1835                if atom[ct-1].split()[0] in ['C','N']:
1836                    Backbone.append(list(np.inner(Amat,np.array([x,y,z]))))
1837                    BackboneColor.append(list(color))
1838                   
1839            if atom[cs+1] == 'type':
1840                RenderLabel(x,y,z,atom[ct],radius)
1841            elif atom[cs+1] == 'name':
1842                RenderLabel(x,y,z,atom[ct-1],radius)
1843            elif atom[cs+1] == 'number':
1844                RenderLabel(x,y,z,str(iat+1),radius)
1845            elif atom[cs+1] == 'residue' and atom[ct-1] == 'CA':
1846                RenderLabel(x,y,z,atom[ct-4],radius)
1847            elif atom[cs+1] == '1-letter' and atom[ct-1] == 'CA':
1848                RenderLabel(x,y,z,atom[ct-3],radius)
1849            elif atom[cs+1] == 'chain' and atom[ct-1] == 'CA':
1850                RenderLabel(x,y,z,atom[ct-2],radius)
1851        if Backbone:
1852            RenderBackbone(Backbone,BackboneColor,bondR)
1853#        print time.time()-time0
1854        Page.canvas.SwapBuffers()
1855       
1856    def OnSize(event):
1857        Draw()
1858       
1859    def OnFocus(event):
1860        Draw()
1861       
1862    try:
1863        plotNum = self.G2plotNB.plotList.index(generalData['Name'])
1864        Page = self.G2plotNB.nb.GetPage(plotNum)       
1865    except (ValueError,error):
1866        Plot = self.G2plotNB.addOgl(generalData['Name'])
1867        plotNum = self.G2plotNB.plotList.index(generalData['Name'])
1868        Page = self.G2plotNB.nb.GetPage(plotNum)
1869        Page.views = False
1870        view = False
1871        altDown = False
1872    Page.SetFocus()
1873    cb = wx.ComboBox(self.G2plotNB.status,style=wx.CB_DROPDOWN|wx.CB_READONLY,
1874        choices=(' save as:','jpeg','tiff','bmp'))
1875    cb.Bind(wx.EVT_COMBOBOX, OnKeyBox)
1876    cb.SetValue(' save as:')
1877    Page.canvas.Bind(wx.EVT_MOUSEWHEEL, OnMouseWheel)
1878    Page.canvas.Bind(wx.EVT_LEFT_DOWN, OnMouseDown)
1879    Page.canvas.Bind(wx.EVT_RIGHT_DOWN, OnMouseDown)
1880    Page.canvas.Bind(wx.EVT_MIDDLE_DOWN, OnMouseDown)
1881    Page.canvas.Bind(wx.EVT_KEY_UP, OnKey)
1882    Page.canvas.Bind(wx.EVT_MOTION, OnMouseMove)
1883    Page.canvas.Bind(wx.EVT_SIZE, OnSize)
1884    Page.canvas.Bind(wx.EVT_SET_FOCUS, OnFocus)
1885    Page.camera['position'] = drawingData['cameraPos']
1886    Page.camera['viewPoint'] = np.inner(Amat,drawingData['viewPoint'][0])
1887    Page.camera['backColor'] = np.array(list(drawingData['backColor'])+[0,])/255.
1888    Page.canvas.SetCurrent()
1889    Draw()
1890       
Note: See TracBrowser for help on using the repository browser.