source: trunk/GSASIIplot.py @ 687

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

histogram wtFactors for PWDR & HKLF (latter not implemented yet)
means that data[0] is now {'wtFactor':1.0} instead of a mostly unused file name
also put in atomIndx for doing restraint calcs.

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