source: trunk/GSASIIplot.py @ 766

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

map peaks now show with gray scale proportional to density

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