source: trunk/GSASIIplot.py @ 1334

Last change on this file since 1334 was 1334, checked in by vondreele, 8 years ago

auto reverse plot if seq refinement done in reverse

  • Property svn:eol-style set to native
  • Property svn:keywords set to Date Author Revision URL Id
File size: 183.3 KB
Line 
1# -*- coding: utf-8 -*-
2'''
3*GSASIIplot: plotting routines*
4===============================
5
6'''
7########### SVN repository information ###################
8# $Date: 2014-05-07 18:47:58 +0000 (Wed, 07 May 2014) $
9# $Author: vondreele $
10# $Revision: 1334 $
11# $URL: trunk/GSASIIplot.py $
12# $Id: GSASIIplot.py 1334 2014-05-07 18:47:58Z vondreele $
13########### SVN repository information ###################
14import math
15import time
16import copy
17import sys
18import os.path
19import numpy as np
20import numpy.ma as ma
21import numpy.linalg as nl
22import wx
23import wx.aui
24import wx.glcanvas
25import matplotlib as mpl
26import mpl_toolkits.mplot3d.axes3d as mp3d
27import GSASIIpath
28GSASIIpath.SetVersionNumber("$Revision: 1334 $")
29import GSASIIgrid as G2gd
30import GSASIIimage as G2img
31import GSASIIpwd as G2pwd
32import GSASIIIO as G2IO
33import GSASIIpwdGUI as G2pdG
34import GSASIIimgGUI as G2imG
35import GSASIIphsGUI as G2phG
36import GSASIIlattice as G2lat
37import GSASIIspc as G2spc
38import GSASIImath as G2mth
39import pytexture as ptx
40from  OpenGL.GL import *
41from OpenGL.GLU import *
42from OpenGL.GLE import *
43import gltext
44from matplotlib.backends.backend_wx import _load_bitmap
45from matplotlib.backends.backend_wxagg import FigureCanvasWxAgg as Canvas
46from matplotlib.backends.backend_wxagg import NavigationToolbar2Wx as Toolbar
47
48# useful degree trig functions
49sind = lambda x: math.sin(x*math.pi/180.)
50cosd = lambda x: math.cos(x*math.pi/180.)
51tand = lambda x: math.tan(x*math.pi/180.)
52asind = lambda x: 180.*math.asin(x)/math.pi
53acosd = lambda x: 180.*math.acos(x)/math.pi
54atan2d = lambda x,y: 180.*math.atan2(y,x)/math.pi
55atand = lambda x: 180.*math.atan(x)/math.pi
56# numpy versions
57npsind = lambda x: np.sin(x*np.pi/180.)
58npcosd = lambda x: np.cos(x*np.pi/180.)
59nptand = lambda x: np.tan(x*np.pi/180.)
60npacosd = lambda x: 180.*np.arccos(x)/np.pi
61npasind = lambda x: 180.*np.arcsin(x)/np.pi
62npatand = lambda x: 180.*np.arctan(x)/np.pi
63npatan2d = lambda x,y: 180.*np.arctan2(x,y)/np.pi
64   
65class G2PlotMpl(wx.Panel):   
66    'needs a doc string'
67    def __init__(self,parent,id=-1,dpi=None,**kwargs):
68        wx.Panel.__init__(self,parent,id=id,**kwargs)
69        mpl.rcParams['legend.fontsize'] = 10
70        self.figure = mpl.figure.Figure(dpi=dpi,figsize=(5,6))
71        self.canvas = Canvas(self,-1,self.figure)
72        self.toolbar = GSASIItoolbar(self.canvas)
73
74        self.toolbar.Realize()
75       
76        sizer=wx.BoxSizer(wx.VERTICAL)
77        sizer.Add(self.canvas,1,wx.EXPAND)
78        sizer.Add(self.toolbar,0,wx.LEFT|wx.EXPAND)
79        self.SetSizer(sizer)
80       
81class G2PlotOgl(wx.Panel):
82    'needs a doc string'
83    def __init__(self,parent,id=-1,dpi=None,**kwargs):
84        self.figure = wx.Panel.__init__(self,parent,id=id,**kwargs)
85        if 'win' in sys.platform:           #Windows already double buffered
86            self.canvas = wx.glcanvas.GLCanvas(self,-1,**kwargs)
87        else:                               #fix from Jim Hester for X systems
88            attribs = (wx.glcanvas.WX_GL_DOUBLEBUFFER,)         
89            self.canvas = wx.glcanvas.GLCanvas(self,-1,attribList=attribs,**kwargs)
90        self.camera = {}
91        sizer=wx.BoxSizer(wx.VERTICAL)
92        sizer.Add(self.canvas,1,wx.EXPAND)
93        self.SetSizer(sizer)
94       
95class G2Plot3D(wx.Panel):
96    'needs a doc string'
97    def __init__(self,parent,id=-1,dpi=None,**kwargs):
98        wx.Panel.__init__(self,parent,id=id,**kwargs)
99        self.figure = mpl.figure.Figure(dpi=dpi,figsize=(6,6))
100        self.canvas = Canvas(self,-1,self.figure)
101        self.toolbar = GSASIItoolbar(self.canvas)
102
103        self.toolbar.Realize()
104       
105        sizer=wx.BoxSizer(wx.VERTICAL)
106        sizer.Add(self.canvas,1,wx.EXPAND)
107        sizer.Add(self.toolbar,0,wx.LEFT|wx.EXPAND)
108        self.SetSizer(sizer)
109                             
110class G2PlotNoteBook(wx.Panel):
111    'create a tabbed window for plotting'
112    def __init__(self,parent,id=-1):
113        wx.Panel.__init__(self,parent,id=id)
114        #so one can't delete a plot page!!
115        self.nb = wx.aui.AuiNotebook(self, \
116            style=wx.aui.AUI_NB_DEFAULT_STYLE ^ wx.aui.AUI_NB_CLOSE_ON_ACTIVE_TAB)
117        sizer = wx.BoxSizer()
118        sizer.Add(self.nb,1,wx.EXPAND)
119        self.SetSizer(sizer)
120        self.status = parent.CreateStatusBar()
121        self.status.SetFieldsCount(2)
122        self.status.SetStatusWidths([150,-1])
123        self.Bind(wx.aui.EVT_AUINOTEBOOK_PAGE_CHANGED, self.OnPageChanged)
124        self.nb.Bind(wx.EVT_KEY_UP,self.OnNotebookKey)
125       
126        self.plotList = []
127           
128    def OnNotebookKey(self,event):
129        '''Called when a keystroke event gets picked up by the notebook window
130        rather the child. This is not expected, but somehow it does sometimes
131        on the Mac and perhaps Linux.
132
133        Assume that the page associated with the currently displayed tab
134        has a child, .canvas; give that child the focus and pass it the event.
135        '''
136        try:
137            Page = self.nb.GetPage(self.nb.GetSelection())
138        except ValueError: # occurs with no plot tabs
139            return
140        try:
141            Page.canvas.SetFocus()
142            wx.PostEvent(Page.canvas,event)
143        except AttributeError:
144            pass
145
146    def addMpl(self,name=""):
147        'Add a tabbed page with a matplotlib plot'
148        page = G2PlotMpl(self.nb)
149        self.nb.AddPage(page,name)
150       
151        self.plotList.append(name)
152       
153        return page.figure
154       
155    def add3D(self,name=""):
156        'Add a tabbed page with a 3D plot'
157        page = G2Plot3D(self.nb)
158        self.nb.AddPage(page,name)
159       
160        self.plotList.append(name)
161       
162        return page.figure
163       
164    def addOgl(self,name=""):
165        'Add a tabbed page with an openGL plot'
166        page = G2PlotOgl(self.nb)
167        self.nb.AddPage(page,name)
168       
169        self.plotList.append(name)
170       
171        return page.figure
172       
173    def Delete(self,name):
174        'delete a tabbed page'
175        try:
176            item = self.plotList.index(name)
177            del self.plotList[item]
178            self.nb.DeletePage(item)
179        except ValueError:          #no plot of this name - do nothing
180            return     
181               
182    def clear(self):
183        'clear all pages from plot window'
184        while self.nb.GetPageCount():
185            self.nb.DeletePage(0)
186        self.plotList = []
187        self.status.DestroyChildren()
188       
189    def Rename(self,oldName,newName):
190        'rename a tab'
191        try:
192            item = self.plotList.index(oldName)
193            self.plotList[item] = newName
194            self.nb.SetPageText(item,newName)
195        except ValueError:          #no plot of this name - do nothing
196            return     
197       
198    def OnPageChanged(self,event):
199        'respond to someone pressing a tab on the plot window'
200        if self.plotList:
201            self.status.SetStatusText('Better to select this from GSAS-II data tree',1)
202        self.status.DestroyChildren()                           #get rid of special stuff on status bar
203       
204class GSASIItoolbar(Toolbar):
205    'needs a doc string'
206    ON_MPL_HELP = wx.NewId()
207    ON_MPL_KEY = wx.NewId()
208    def __init__(self,plotCanvas):
209        Toolbar.__init__(self,plotCanvas)
210        POSITION_OF_CONFIGURE_SUBPLOTS_BTN = 6
211        self.DeleteToolByPos(POSITION_OF_CONFIGURE_SUBPLOTS_BTN)
212        parent = self.GetParent()
213        key = os.path.join(os.path.split(__file__)[0],'key.ico')
214        self.AddSimpleTool(self.ON_MPL_KEY,_load_bitmap(key),'Key press','Select key press')
215        wx.EVT_TOOL(self,self.ON_MPL_KEY,self.OnKey)
216        help = os.path.join(os.path.split(__file__)[0],'help.ico')
217        self.AddSimpleTool(self.ON_MPL_HELP,_load_bitmap(help),'Help on','Show help on')
218        wx.EVT_TOOL(self,self.ON_MPL_HELP,self.OnHelp)
219    def OnHelp(self,event):
220        'needs a doc string'
221        Page = self.GetParent().GetParent()
222        pageNo = Page.GetSelection()
223        bookmark = Page.GetPageText(pageNo)
224        bookmark = bookmark.strip(')').replace('(','_')
225        G2gd.ShowHelp(bookmark,self.TopLevelParent)
226    def OnKey(self,event):
227        'needs a doc string'
228        parent = self.GetParent()
229        if parent.Choice:
230            dlg = wx.SingleChoiceDialog(parent,'Select','Key press',list(parent.Choice))
231            if dlg.ShowModal() == wx.ID_OK:
232                sel = dlg.GetSelection()
233                event.key = parent.Choice[sel][0]
234                parent.keyPress(event)
235            dlg.Destroy()
236           
237################################################################################
238##### PlotSngl
239################################################################################
240           
241def PlotSngl(G2frame,newPlot=False,Data=None,hklRef=None,Title=''):
242    '''Single crystal structure factor plotting package - displays zone of reflections as rings proportional
243        to F, F**2, etc. as requested
244    '''
245    from matplotlib.patches import Circle,CirclePolygon
246    global HKL,HKLF
247   
248    def OnSCKeyPress(event):
249        i = zones.index(Data['Zone'])
250        newPlot = False
251        typeChoice = {'1':'|DFsq|>sig','3':'|DFsq|>3sig','w':'|DFsq|/sig','f':'Fo','s':'Fosq'}
252        if event.key == 'h':
253            Data['Zone'] = '100'
254            newPlot = True
255        elif event.key == 'k':
256            Data['Zone'] = '010'
257            newPlot = True
258        elif event.key == 'l':
259            Data['Zone'] = '001'
260            newPlot = True
261        elif event.key == 'u':
262            Data['Scale'] *= 1.1
263        elif event.key == 'd':
264            Data['Scale'] /= 1.1
265        elif event.key == '+':
266            Data['Layer'] = min(Data['Layer']+1,HKLmax[i])
267        elif event.key == '-':
268            Data['Layer'] = max(Data['Layer']-1,HKLmin[i])
269        elif event.key == '0':
270            Data['Layer'] = 0
271        elif event.key in typeChoice and 'HKLF' in Name:
272            Data['Type'] = typeChoice[event.key]           
273            newPlot = True
274        PlotSngl(G2frame,newPlot,Data,hklRef,Title)
275
276    def OnSCMotion(event):
277        xpos = event.xdata
278        if xpos:
279            xpos = round(xpos)                                        #avoid out of frame mouse position
280            ypos = round(event.ydata)
281            zpos = Data['Layer']
282            if '100' in Data['Zone']:
283                HKLtxt = '(%3d,%3d,%3d)'%(zpos,xpos,ypos)
284            elif '010' in Data['Zone']:
285                HKLtxt = '(%3d,%3d,%3d)'%(xpos,zpos,ypos)
286            elif '001' in Data['Zone']:
287                HKLtxt = '(%3d,%3d,%3d)'%(xpos,ypos,zpos)
288            Page.canvas.SetToolTipString(HKLtxt)
289            G2frame.G2plotNB.status.SetStatusText('HKL = '+HKLtxt,0)
290            G2frame.G2plotNB.status.SetStatusText('Use K-box to set plot controls',1)
291               
292    def OnSCPress(event):
293        zpos = Data['Layer']
294        xpos = event.xdata
295        if xpos:
296            pos = int(round(event.xdata)),int(round(event.ydata))
297            if '100' in Data['Zone']:
298                Page.canvas.SetToolTipString('(picked:(%3d,%3d,%3d))'%(zpos,pos[0],pos[1]))
299                hkl = np.array([zpos,pos[0],pos[1]])
300            elif '010' in Data['Zone']:
301                Page.canvas.SetToolTipString('(picked:(%3d,%3d,%3d))'%(pos[0],zpos,pos[1]))
302                hkl = np.array([pos[0],zpos,pos[1]])
303            elif '001' in Data['Zone']:
304                Page.canvas.SetToolTipString('(picked:(%3d,%3d,%3d))'%(pos[0],pos[1],zpos))
305                hkl = np.array([pos[0],pos[1],zpos])
306            h,k,l = hkl
307            hklf = HKLF[np.where(np.all(HKL-hkl == [0,0,0],axis=1))]
308            if len(hklf):
309                Fosq,sig,Fcsq = hklf[0]
310                HKLtxt = '( %.2f %.3f %.2f %.2f)'%(Fosq,sig,Fcsq,(Fosq-Fcsq)/(scale*sig))
311                G2frame.G2plotNB.status.SetStatusText('Fosq, sig, Fcsq, delFsq/sig = '+HKLtxt,1)
312                                 
313    Name = G2frame.PatternTree.GetItemText(G2frame.PatternId)
314    if not Title:
315        Title = Name
316    try:
317        plotNum = G2frame.G2plotNB.plotList.index('Structure Factors')
318        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
319        if not newPlot:
320            Plot = Page.figure.gca()          #get previous powder plot & get limits
321            xylim = Plot.get_xlim(),Plot.get_ylim()
322        Page.figure.clf()
323        Plot = Page.figure.gca()          #get a fresh plot after clf()
324    except ValueError:
325        Plot = G2frame.G2plotNB.addMpl('Structure Factors').gca()
326        plotNum = G2frame.G2plotNB.plotList.index('Structure Factors')
327        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
328        Page.canvas.mpl_connect('button_press_event', OnSCPress)
329        Page.canvas.mpl_connect('motion_notify_event', OnSCMotion)
330        Page.canvas.mpl_connect('key_press_event', OnSCKeyPress)
331        Page.keyPress = OnSCKeyPress
332        if 'HKLF' in Name:
333            Page.Choice = (' key press','u: increase scale','d: decrease scale',
334                'f: select Fo','s: select Fosq','w: select |DFsq|/sig',
335                '1: select |DFsq|>sig','3 select |DFsq|>3sig',
336                'h: select 100 zone','k: select 010 zone','l: select 001 zone',
337                '+: increase index','-: decrease index','0: zero layer',)
338        else:   
339            Page.Choice = (' key press','u: increase scale','d: decrease scale',
340                'h: select 100 zone','k: select 010 zone','l: select 001 zone',
341                '+: increase index','-: decrease index','0: zero layer',)
342    Page.SetFocus()
343   
344    G2frame.G2plotNB.status.SetStatusText('Use K-box to set plot controls',1)
345    Plot.set_aspect(aspect='equal')
346   
347    Type = Data['Type']           
348    scale = Data['Scale']
349    HKLmax = Data['HKLmax']
350    HKLmin = Data['HKLmin']
351    FosqMax = Data['FoMax']
352    FoMax = math.sqrt(FosqMax)
353    xlabel = ['k, h=','h, k=','h, l=']
354    ylabel = ['l','l','k']
355    zones = ['100','010','001']
356    pzone = [[1,2],[0,2],[0,1]]
357    izone = zones.index(Data['Zone'])
358    Plot.set_title(Data['Type']+' '+Title)
359    HKL = []
360    HKLF = []
361    time0 = time.time()
362    for refl in hklRef:
363        H = np.array(refl[:3])
364        if 'HKLF' in Name:
365            Fosq,sig,Fcsq = refl[5:8]
366        else:
367            Fosq,sig,Fcsq = refl[8],1.0,refl[9]
368        HKL.append(H)
369        HKLF.append([Fosq,sig,Fcsq])
370        if H[izone] == Data['Layer']:
371            A = 0
372            B = 0
373            if Type == 'Fosq':
374                A = scale*Fosq/FosqMax
375                B = scale*Fcsq/FosqMax
376                C = abs(A-B)
377            elif Type == 'Fo':
378                A = scale*math.sqrt(max(0,Fosq))/FoMax
379                B = scale*math.sqrt(max(0,Fcsq))/FoMax
380                C = abs(A-B)
381            elif Type == '|DFsq|/sig':
382                if sig > 0.:
383                    A = (Fosq-Fcsq)/(3*sig)
384                B = 0
385            elif Type == '|DFsq|>sig':
386                if sig > 0.:
387                    A = (Fosq-Fcsq)/(3*sig)
388                if abs(A) < 1.0: A = 0
389                B = 0                   
390            elif Type == '|DFsq|>3sig':
391                if sig > 0.:
392                    A = (Fosq-Fcsq)/(3*sig)
393                if abs(A) < 3.0: A = 0
394                B = 0                   
395            xy = (H[pzone[izone][0]],H[pzone[izone][1]])
396            if Type in ['|DFsq|/sig','|DFsq|>sig','|DFsq|>3sig']:
397                if A > 0.0:
398                    Plot.add_artist(Circle(xy,radius=A,ec='g',fc='w'))
399                else:
400                    Plot.add_artist(Circle(xy,radius=-A,ec='r',fc='w'))
401            else:
402                if A > 0.0 and A > B:
403                    Plot.add_artist(Circle(xy,radius=A,ec='g',fc='w'))
404                if B:
405                    Plot.add_artist(Circle(xy,radius=B,ec='b',fc='w'))
406                    if A < B:
407                        Plot.add_artist(Circle(xy,radius=A,ec='g',fc='w'))
408                    radius = C
409                    if radius > 0:
410                        if A > B:
411                            Plot.add_artist(Circle(xy,radius=radius,ec='g',fc='g'))
412                        else:                   
413                            Plot.add_artist(Circle(xy,radius=radius,ec='r',fc='r'))
414#    print 'plot time: %.3f'%(time.time()-time0)
415    HKL = np.array(HKL,dtype=np.int)
416    HKLF = np.array(HKLF)
417    Plot.set_xlabel(xlabel[izone]+str(Data['Layer']),fontsize=12)
418    Plot.set_ylabel(ylabel[izone],fontsize=12)
419    if not newPlot:
420        Page.toolbar.push_current()
421        Plot.set_xlim(xylim[0])
422        Plot.set_ylim(xylim[1])
423#        xylim = []
424        Page.toolbar.push_current()
425        Page.toolbar.draw()
426    else:
427        Plot.set_xlim((HKLmin[pzone[izone][0]],HKLmax[pzone[izone][0]]))
428        Plot.set_ylim((HKLmin[pzone[izone][1]],HKLmax[pzone[izone][1]]))
429        Page.canvas.draw()
430       
431################################################################################
432##### PlotPatterns
433################################################################################
434           
435def PlotPatterns(G2frame,newPlot=False,plotType='PWDR'):
436    '''Powder pattern plotting package - displays single or multiple powder patterns as intensity vs
437    2-theta, q or TOF. Can display multiple patterns as "waterfall plots" or contour plots. Log I
438    plotting available.
439    '''
440    global HKL
441    global exclLines
442    plottype = plotType
443    def OnPlotKeyPress(event):
444        newPlot = False
445        if event.key == 'w' and 'PWDR' in plottype:
446            G2frame.ErrorBars = not G2frame.ErrorBars
447            if not G2frame.Weight:
448                G2frame.SinglePlot = True
449            newPlot = True
450        elif event.key == 'e' and 'SASD' in plottype:
451            G2frame.ErrorBars = not G2frame.ErrorBars
452        elif event.key == 'b':
453            G2frame.SubBack = not G2frame.SubBack
454            if not G2frame.SubBack:
455                G2frame.SinglePlot = True               
456        elif event.key == 'n':
457            if G2frame.Contour:
458                pass
459            else:
460                G2frame.logPlot = not G2frame.logPlot
461                if not G2frame.logPlot:
462                    G2frame.Offset[0] = 0
463                newPlot = True
464        elif event.key == 'u':
465            if G2frame.Contour:
466                G2frame.Cmax = min(1.0,G2frame.Cmax*1.2)
467            elif G2frame.Offset[0] < 100.:
468                G2frame.Offset[0] += 1.
469        elif event.key == 'd':
470            if G2frame.Contour:
471                G2frame.Cmax = max(0.0,G2frame.Cmax*0.8)
472            elif G2frame.Offset[0] > 0.:
473                G2frame.Offset[0] -= 1.
474        elif event.key == 'l':
475            G2frame.Offset[1] -= 1.
476        elif event.key == 'r':
477            G2frame.Offset[1] += 1.
478        elif event.key == 'o':
479            G2frame.Cmax = 1.0
480            G2frame.Offset = [0,0]
481        elif event.key == 'c' and 'PWDR' in plottype:
482            newPlot = True
483            G2frame.Contour = not G2frame.Contour
484            if not G2frame.Contour:
485                G2frame.SinglePlot = False
486                G2frame.Offset = [0.,0.]
487        elif event.key == 'q': 
488            if 'PWDR' in plottype:
489                newPlot = True
490                G2frame.qPlot = not G2frame.qPlot
491            elif 'SASD' in plottype:
492                newPlot = True
493                G2frame.sqPlot = not G2frame.sqPlot       
494        elif event.key == 's':
495            if G2frame.Contour:
496                choice = [m for m in mpl.cm.datad.keys() if not m.endswith("_r")]
497                choice.sort()
498                dlg = wx.SingleChoiceDialog(G2frame,'Select','Color scheme',choice)
499                if dlg.ShowModal() == wx.ID_OK:
500                    sel = dlg.GetSelection()
501                    G2frame.ContourColor = choice[sel]
502                else:
503                    G2frame.ContourColor = 'Paired'
504                dlg.Destroy()
505            else:
506                G2frame.SinglePlot = not G2frame.SinglePlot               
507            newPlot = True
508        elif event.key == '+':
509            if G2frame.PickId:
510                G2frame.PickId = False
511        elif event.key == 'i' and G2frame.Contour:                  #for smoothing contour plot
512            choice = ['nearest','bilinear','bicubic','spline16','spline36','hanning',
513               'hamming','hermite','kaiser','quadric','catrom','gaussian','bessel',
514               'mitchell','sinc','lanczos']
515            dlg = wx.SingleChoiceDialog(G2frame,'Select','Interpolation',choice)
516            if dlg.ShowModal() == wx.ID_OK:
517                sel = dlg.GetSelection()
518                G2frame.Interpolate = choice[sel]
519            else:
520                G2frame.Interpolate = 'nearest'
521            dlg.Destroy()
522           
523        wx.CallAfter(PlotPatterns,G2frame,newPlot=newPlot,plotType=plottype)
524       
525    def OnMotion(event):
526        xpos = event.xdata
527        if xpos:                                        #avoid out of frame mouse position
528            ypos = event.ydata
529            Page.canvas.SetCursor(wx.CROSS_CURSOR)
530            try:
531                Parms,Parms2 = G2frame.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(G2frame,G2frame.PatternId, 'Instrument Parameters'))
532                if 'C' in Parms['Type'][0]:
533                    wave = G2mth.getWave(Parms)
534                    if G2frame.qPlot and 'PWDR' in plottype:
535                        try:
536                            xpos = 2.0*asind(xpos*wave/(4*math.pi))
537                        except ValueError:      #avoid bad value in asin beyond upper limit
538                            pass
539                    dsp = 0.0
540                    if abs(xpos) > 0.:
541                        if 'PWDR' in plottype:                  #avoid possible singularity at beam center
542                            dsp = wave/(2.*sind(abs(xpos)/2.0))
543                            q = 2.*np.pi/dsp
544                        elif 'SASD' in plottype:
545                            dsp = 2*np.pi/xpos
546                            q = xpos
547                    if G2frame.Contour:
548                        if 'PWDR' in plottype:
549                            G2frame.G2plotNB.status.SetStatusText('2-theta =%9.3f d =%9.5f q = %9.5f pattern ID =%5d'%(xpos,dsp,q,int(ypos)),1)
550                        elif 'SASD' in plottype:
551                            G2frame.G2plotNB.status.SetStatusText('d =%9.5f q = %9.5f pattern ID =%5d'%(dsp,q,int(ypos)),1)
552                    else:
553                        if 'PWDR' in plottype:
554                            G2frame.G2plotNB.status.SetStatusText('2-theta =%9.3f d =%9.5f q = %9.5f Intensity =%9.2f'%(xpos,dsp,q,ypos),1)
555                        elif 'SASD' in plottype:
556                            G2frame.G2plotNB.status.SetStatusText('d =%9.5f q = %9.5f Intensity =%12.5g'%(dsp,q,ypos),1)
557                else:       #TOF neutrons
558                    dsp = 0.0
559                    difC = Parms['difC'][1]
560                    dsp = xpos/difC             #rough approx.!
561                    q = 2.*np.pi/dsp
562                    if G2frame.Contour:
563                        G2frame.G2plotNB.status.SetStatusText('TOF =%9.3f d =%9.5f q = %9.5f pattern ID =%5d'%(xpos,dsp,q,int(ypos)),1)
564                    else:
565                        G2frame.G2plotNB.status.SetStatusText('TOF =%9.3f d =%9.5f q = %9.5f Intensity =%9.2f'%(xpos,dsp,q,ypos),1)
566                if G2frame.itemPicked:
567                    Page.canvas.SetToolTipString('%9.5f'%(xpos))
568                if G2frame.PickId:
569                    found = []
570                    if G2frame.PatternTree.GetItemText(G2frame.PickId) in ['Index Peak List','Unit Cells List','Reflection Lists'] or \
571                        'PWDR' in G2frame.PatternTree.GetItemText(PickId):
572                        if len(HKL):
573                            view = Page.toolbar._views.forward()[0][:2]
574                            wid = view[1]-view[0]
575                            found = HKL[np.where(np.fabs(HKL.T[5]-xpos) < 0.002*wid)]
576                        if len(found):
577                            h,k,l = found[0][:3] 
578                            Page.canvas.SetToolTipString('%d,%d,%d'%(int(h),int(k),int(l)))
579                        else:
580                            Page.canvas.SetToolTipString('')
581
582            except TypeError:
583                G2frame.G2plotNB.status.SetStatusText('Select '+plottype+' powder pattern first',1)
584               
585    def OnPress(event): #ugh - this removes a matplotlib error for mouse clicks in log plots                 
586        olderr = np.seterr(invalid='ignore')
587                                                   
588    def OnPick(event):
589        if G2frame.itemPicked is not None: return
590        PatternId = G2frame.PatternId
591        try:
592            Parms,Parms2 = G2frame.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(G2frame,G2frame.PatternId, 'Instrument Parameters'))
593        except TypeError:
594            return
595        if 'C' in Parms['Type'][0]:
596            wave = G2mth.getWave(Parms)
597        else:
598            difC = Parms['difC'][1]
599        PickId = G2frame.PickId
600        pick = event.artist
601        mouse = event.mouseevent
602        xpos = pick.get_xdata()
603        ypos = pick.get_ydata()
604        ind = event.ind
605        xy = list(zip(np.take(xpos,ind),np.take(ypos,ind))[0])
606        if G2frame.PatternTree.GetItemText(PickId) == 'Peak List':
607            if ind.all() != [0]:                                    #picked a data point
608                data = G2frame.PatternTree.GetItemPyData(G2frame.PickId)
609                if 'C' in Parms['Type'][0]:                            #CW data - TOF later in an elif
610                    if G2frame.qPlot:                              #qplot - convert back to 2-theta
611                        xy[0] = 2.0*asind(xy[0]*wave/(4*math.pi))
612                XY = G2mth.setPeakparms(Parms,Parms2,xy[0],xy[1])
613                data.append(XY)
614                G2pdG.UpdatePeakGrid(G2frame,data)
615                PlotPatterns(G2frame,plotType=plottype)
616            else:                                                   #picked a peak list line
617                G2frame.itemPicked = pick
618        elif G2frame.PatternTree.GetItemText(PickId) == 'Limits':
619            if ind.all() != [0]:                                    #picked a data point
620                LimitId = G2gd.GetPatternTreeItemId(G2frame,PatternId, 'Limits')
621                data = G2frame.PatternTree.GetItemPyData(LimitId)
622                if 'C' in Parms['Type'][0]:                            #CW data - TOF later in an elif
623                    if G2frame.qPlot and 'PWDR' in plottype:                              #qplot - convert back to 2-theta
624                        xy[0] = 2.0*asind(xy[0]*wave/(4*math.pi))
625                if G2frame.ifGetExclude:
626                    excl = [0,0]
627                    excl[0] = max(data[1][0],min(xy[0],data[1][1]))
628                    excl[1] = excl[0]+0.1
629                    data.append(excl)
630                    G2frame.ifGetExclude = False
631                else:
632                    if mouse.button==1:
633                        data[1][0] = min(xy[0],data[1][1])
634                    if mouse.button==3:
635                        data[1][1] = max(xy[0],data[1][0])
636                G2frame.PatternTree.SetItemPyData(LimitId,data)
637                G2pdG.UpdateLimitsGrid(G2frame,data,plottype)
638                wx.CallAfter(PlotPatterns,G2frame,plotType=plottype)
639            else:                                                   #picked a limit line
640                G2frame.itemPicked = pick
641        elif G2frame.PatternTree.GetItemText(PickId) == 'Models':
642            if ind.all() != [0]:                                    #picked a data point
643                LimitId = G2gd.GetPatternTreeItemId(G2frame,PatternId, 'Limits')
644                data = G2frame.PatternTree.GetItemPyData(LimitId)
645                if mouse.button==1:
646                    data[1][0] = min(xy[0],data[1][1])
647                if mouse.button==3:
648                    data[1][1] = max(xy[0],data[1][0])
649                G2frame.PatternTree.SetItemPyData(LimitId,data)
650                wx.CallAfter(PlotPatterns,G2frame,plotType=plottype)
651            else:                                                   #picked a limit line
652                G2frame.itemPicked = pick
653        elif G2frame.PatternTree.GetItemText(PickId) == 'Reflection Lists' or \
654            'PWDR' in G2frame.PatternTree.GetItemText(PickId):
655            G2frame.itemPicked = pick
656            pick = str(pick)
657       
658    def OnRelease(event):
659        if G2frame.itemPicked is None: return
660        Parms,Parms2 = G2frame.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(G2frame,G2frame.PatternId, 'Instrument Parameters'))
661        if 'C' in Parms['Type'][0]:
662            wave = G2mth.getWave(Parms)
663        else:
664            difC = Parms['difC'][1]
665        xpos = event.xdata
666        PickId = G2frame.PickId
667        if G2frame.PatternTree.GetItemText(PickId) in ['Peak List','Limits'] and xpos:
668            lines = []
669            for line in G2frame.Lines: 
670                lines.append(line.get_xdata()[0])
671            try:
672                lineNo = lines.index(G2frame.itemPicked.get_xdata()[0])
673            except ValueError:
674                lineNo = -1
675            if  lineNo in [0,1] or lineNo in exclLines:
676                LimitId = G2gd.GetPatternTreeItemId(G2frame,G2frame.PatternId, 'Limits')
677                data = G2frame.PatternTree.GetItemPyData(LimitId)
678                id = lineNo/2+1
679                id2 = lineNo%2
680                if G2frame.qPlot and 'PWDR' in plottype:
681                    if 'C' in Parms['Type'][0]:
682                        data[id][id2] = 2.0*asind(wave*xpos/(4*math.pi))
683                    else:
684                        data[id][id2] = 2*math.pi*Parms['difC'][1]/xpos
685                else:
686                    data[id][id2] = xpos
687                if id > 1 and data[id][0] > data[id][1]:
688                        data[id].reverse()
689                data[1][0] = min(max(data[0][0],data[1][0]),data[1][1])
690                data[1][1] = max(min(data[0][1],data[1][1]),data[1][0])
691                G2frame.PatternTree.SetItemPyData(LimitId,data)
692                if G2frame.PatternTree.GetItemText(G2frame.PickId) == 'Limits':
693                    G2pdG.UpdateLimitsGrid(G2frame,data,plottype)
694            elif lineNo > 1:
695                PeakId = G2gd.GetPatternTreeItemId(G2frame,G2frame.PatternId, 'Peak List')
696                data = G2frame.PatternTree.GetItemPyData(PeakId)
697                if event.button == 3:
698                    del data[lineNo-2]
699                else:
700                    if G2frame.qPlot:
701                        data[lineNo-2][0] = 2.0*asind(wave*xpos/(4*math.pi))
702                    else:
703                        data[lineNo-2][0] = xpos
704                G2frame.PatternTree.SetItemPyData(PeakId,data)
705                G2pdG.UpdatePeakGrid(G2frame,data)
706        elif G2frame.PatternTree.GetItemText(PickId) in ['Models',] and xpos:
707            lines = []
708            for line in G2frame.Lines: 
709                lines.append(line.get_xdata()[0])
710            try:
711                lineNo = lines.index(G2frame.itemPicked.get_xdata()[0])
712            except ValueError:
713                lineNo = -1
714            if  lineNo in [0,1]:
715                LimitId = G2gd.GetPatternTreeItemId(G2frame,G2frame.PatternId, 'Limits')
716                data = G2frame.PatternTree.GetItemPyData(LimitId)
717                data[1][lineNo] = xpos
718                data[1][0] = min(max(data[0][0],data[1][0]),data[1][1])
719                data[1][1] = max(min(data[0][1],data[1][1]),data[1][0])
720                G2frame.PatternTree.SetItemPyData(LimitId,data)       
721        elif (G2frame.PatternTree.GetItemText(PickId) == 'Reflection Lists' or \
722            'PWDR' in G2frame.PatternTree.GetItemText(PickId)) and xpos:
723            Phases = G2frame.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(G2frame,PatternId,'Reflection Lists'))
724            pick = str(G2frame.itemPicked).split('(')[1].strip(')')
725            if 'line' not in pick:       #avoid data points, etc.
726                num = Phases.keys().index(pick)
727                if num:
728                    G2frame.refDelt = -(event.ydata-G2frame.refOffset)/(num*Ymax)
729                else:       #1st row of refl ticks
730                    G2frame.refOffset = event.ydata
731        PlotPatterns(G2frame,plotType=plottype)
732        G2frame.itemPicked = None   
733
734    try:
735        plotNum = G2frame.G2plotNB.plotList.index('Powder Patterns')
736        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
737        if not newPlot:
738            Plot = Page.figure.gca()          #get previous powder plot & get limits
739            xylim = Plot.get_xlim(),Plot.get_ylim()
740        Page.figure.clf()
741        Plot = Page.figure.gca()          #get a fresh plot after clf()
742    except ValueError:
743        if plottype == 'SASD':
744            G2frame.logPlot = True
745            G2frame.ErrorBars = True
746        newPlot = True
747        G2frame.Cmax = 1.0
748        Plot = G2frame.G2plotNB.addMpl('Powder Patterns').gca()
749        plotNum = G2frame.G2plotNB.plotList.index('Powder Patterns')
750        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
751        Page.canvas.mpl_connect('key_press_event', OnPlotKeyPress)
752        Page.canvas.mpl_connect('motion_notify_event', OnMotion)
753        Page.canvas.mpl_connect('pick_event', OnPick)
754        Page.canvas.mpl_connect('button_release_event', OnRelease)
755        Page.canvas.mpl_connect('button_press_event',OnPress)
756    if plottype == 'PWDR':  # avoids a very nasty clash with KILL_FOCUS in SASD TextCtrl?
757        Page.SetFocus()
758    G2frame.G2plotNB.status.DestroyChildren()
759    if G2frame.Contour:
760        Page.Choice = (' key press','d: lower contour max','u: raise contour max','o: reset contour max',
761            'i: interpolation method','s: color scheme','c: contour off')
762    else:
763        if G2frame.logPlot:
764            if 'PWDR' in plottype:
765                Page.Choice = (' key press','n: log(I) off','l: offset left','r: offset right','o: reset offset',
766                    'c: contour on','q: toggle q plot','s: toggle single plot','+: no selection')
767            elif 'SASD' in plottype:
768                Page.Choice = (' key press','b: toggle subtract background','n: semilog on',
769                    'd: offset down','l: offset left','r: offset right','u: offset up','o: reset offset',
770                    'q: toggle S(q) plot','s: toggle single plot','+: no selection')
771        else:
772            if 'PWDR' in plottype:
773                Page.Choice = (' key press','l: offset left','r: offset right','d: offset down',
774                    'u: offset up','o: reset offset','b: toggle subtract background','n: log(I) on','c: contour on',
775                    'q: toggle q plot','s: toggle single plot','w: toggle divide by sig','+: no selection')
776            elif 'SASD' in plottype:
777                Page.Choice = (' key press','b: toggle subtract background','n: loglog on','e: toggle error bars',
778                    'd: offset down','l: offset left','r: offset right','u: offset up','o: reset offset',
779                    'q: toggle S(q) plot','s: toggle single plot','+: no selection')
780
781    Page.keyPress = OnPlotKeyPress   
782    PickId = G2frame.PickId
783    PatternId = G2frame.PatternId
784    colors=['b','g','r','c','m','k']
785    Lines = []
786    exclLines = []
787    if G2frame.SinglePlot:
788        Pattern = G2frame.PatternTree.GetItemPyData(PatternId)
789        Pattern.append(G2frame.PatternTree.GetItemText(PatternId))
790        PlotList = [Pattern,]
791        Parms,Parms2 = G2frame.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(G2frame,
792            G2frame.PatternId, 'Instrument Parameters'))
793        Sample = G2frame.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(G2frame,G2frame.PatternId, 'Sample Parameters'))
794        ParmList = [Parms,]
795        SampleList = [Sample,]
796        Title = Pattern[-1]
797    else:       
798        Title = os.path.split(G2frame.GSASprojectfile)[1]
799        PlotList = []
800        ParmList = []
801        SampleList = []
802        item, cookie = G2frame.PatternTree.GetFirstChild(G2frame.root)
803        while item:
804            if plottype in G2frame.PatternTree.GetItemText(item):
805                Pattern = G2frame.PatternTree.GetItemPyData(item)
806                if len(Pattern) < 3:                    # put name on end if needed
807                    Pattern.append(G2frame.PatternTree.GetItemText(item))
808                PlotList.append(Pattern)
809                ParmList.append(G2frame.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(G2frame,
810                    item,'Instrument Parameters'))[0])
811                SampleList.append(G2frame.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(G2frame,
812                    item, 'Sample Parameters')))
813            item, cookie = G2frame.PatternTree.GetNextChild(G2frame.root, cookie)               
814    lenX = 0
815    if PickId:
816        if G2frame.PatternTree.GetItemText(PickId) in ['Reflection Lists']:
817            Phases = G2frame.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(G2frame,PatternId,'Reflection Lists'))
818            HKL = []
819            if Phases:
820                try:
821                    for peak in Phases[G2frame.RefList]['RefList']:
822                        HKL.append(peak[:6])
823                except TypeError:
824                    for peak in Phases[G2frame.RefList]:
825                        HKL.append(peak[:6])                   
826                HKL = np.array(HKL)
827        else:
828            HKL = np.array(G2frame.HKL)
829    Ymax = None
830    for Pattern in PlotList:
831        xye = Pattern[1]
832        if xye[1] is None: continue
833        if Ymax is None: Ymax = max(xye[1])
834        Ymax = max(Ymax,max(xye[1]))
835    if Ymax is None: return # nothing to plot
836    offset = G2frame.Offset[0]*Ymax/100.0
837    if G2frame.logPlot:
838        Title = 'log('+Title+')'
839    Plot.set_title(Title)
840    if G2frame.qPlot or 'SASD' in plottype:
841        Plot.set_xlabel(r'$Q, \AA^{-1}$',fontsize=16)
842    else:
843        if 'C' in ParmList[0]['Type'][0]:       
844            Plot.set_xlabel(r'$\mathsf{2\theta}$',fontsize=16)
845        else:
846            Plot.set_xlabel(r'$TOF, \mathsf{\mu}$s',fontsize=16)           
847    if G2frame.Weight:
848        Plot.set_ylabel(r'$\mathsf{I/\sigma(I)}$',fontsize=16)
849    else:
850        if 'C' in ParmList[0]['Type'][0]:
851            if 'PWDR' in plottype:
852                Plot.set_ylabel('$Intensity$',fontsize=16)
853            elif 'SASD' in plottype:
854                if G2frame.sqPlot:
855                    Plot.set_ylabel('$S(Q)=I*Q^{4}$',fontsize=16)
856                else:
857                    Plot.set_ylabel('$Intensity, cm^{-1}$',fontsize=16)
858        else:
859            Plot.set_ylabel('Normalized intensity',fontsize=16)
860    if G2frame.Contour:
861        ContourZ = []
862        ContourY = []
863        Nseq = 0
864    if len(PlotList) < 2:
865        G2frame.Contour = False
866    for N,Pattern in enumerate(PlotList):
867        Parms = ParmList[N]
868        Sample = SampleList[N]
869        if 'C' in Parms['Type'][0]:
870            wave = G2mth.getWave(Parms)
871        else:
872            difC = Parms['difC'][1]
873        ifpicked = False
874        LimitId = 0
875        if Pattern[1] is None: continue # skip over uncomputed simulations
876        xye = ma.array(ma.getdata(Pattern[1]))
877        Zero = Parms.get('Zero',[0.,0.])[1]
878        if PickId:
879            ifpicked = Pattern[2] == G2frame.PatternTree.GetItemText(PatternId)
880            LimitId = G2gd.GetPatternTreeItemId(G2frame,PatternId, 'Limits')
881            limits = np.array(G2frame.PatternTree.GetItemPyData(LimitId))
882            excls = limits[2:]
883            for excl in excls:
884                xye[0] = ma.masked_inside(xye[0],excl[0],excl[1])
885        if G2frame.qPlot and 'PWDR' in plottype:
886            Id = G2gd.GetPatternTreeItemId(G2frame,G2frame.root, Pattern[2])
887            if 'C' in Parms['Type'][0]:
888                X = 4*np.pi*npsind((xye[0]-Zero)/2.0)/wave
889            else:
890                X = 2*np.pi*Parms['difC'][1]/(xye[0]-Zero)
891        else:
892            X = xye[0]-Zero
893        if not lenX:
894            lenX = len(X)
895        if 'PWDR' in plottype:
896            Y = xye[1]+offset*N
897        elif 'SASD' in plottype:
898            B = xye[5]
899            if G2frame.sqPlot:
900                Y = xye[1]*Sample['Scale'][0]*(1.005)**(offset*N)*X**4
901            else:
902                Y = xye[1]*Sample['Scale'][0]*(1.005)**(offset*N)
903        if LimitId and ifpicked:
904            limits = np.array(G2frame.PatternTree.GetItemPyData(LimitId))
905            if G2frame.qPlot and 'PWDR' in plottype:
906                if 'C' in Parms['Type'][0]:
907                    limits = 4*np.pi*npsind(limits/2.0)/wave
908                else:
909                    limits = 2*np.pi*difC/limits
910            Lines.append(Plot.axvline(limits[1][0],color='g',dashes=(5,5),picker=3.))   
911            Lines.append(Plot.axvline(limits[1][1],color='r',dashes=(5,5),picker=3.))
912            for i,item in enumerate(limits[2:]):
913                Lines.append(Plot.axvline(item[0],color='m',dashes=(5,5),picker=3.))   
914                Lines.append(Plot.axvline(item[1],color='m',dashes=(5,5),picker=3.))
915                exclLines += [2*i+2,2*i+3]
916        if G2frame.Contour:
917           
918            if lenX == len(X):
919                ContourY.append(N)
920                ContourZ.append(Y)
921                ContourX = X
922                Nseq += 1
923                Plot.set_ylabel('Data sequence',fontsize=12)
924        else:
925            if 'SASD' in plottype and G2frame.logPlot:
926                X *= (1.01)**(G2frame.Offset[1]*N)
927            else:
928                X += G2frame.Offset[1]*.005*N
929            Xum = ma.getdata(X)
930            if ifpicked:
931                Z = xye[3]+offset*N
932                if 'PWDR' in plottype:  #powder background
933                    W = xye[4]+offset*N
934                    D = xye[5]-Ymax*G2frame.delOffset
935                elif 'SASD' in plottype:
936                    if G2frame.sqPlot:
937                        W = xye[4]*X**4
938                        Z = xye[3]*X**4
939                        B = B*X**4
940                    else:
941                        W = xye[4]
942                    if G2frame.SubBack:
943                        YB = Y-B
944                        ZB = Z
945                    else:
946                        YB = Y
947                        ZB = Z+B
948                    Plot.set_yscale("log",nonposy='mask')
949                    Plot.set_ylim(bottom=np.min(np.trim_zeros(Y))/2.,top=np.max(Y)*2.)
950                if G2frame.logPlot:
951                    if 'PWDR' in plottype:
952                        Plot.set_yscale("log",nonposy='mask')
953                        Plot.plot(X,Y,colors[N%6]+'+',picker=3.,clip_on=False)
954                        Plot.plot(X,Z,colors[(N+1)%6],picker=False)
955                        Plot.plot(X,W,colors[(N+2)%6],picker=False)
956                    elif 'SASD' in plottype:
957                        Plot.set_xscale("log",nonposx='mask')
958                        Plot.set_yscale("log",nonposy='mask')
959                        if G2frame.ErrorBars:
960                            if G2frame.sqPlot:
961                                Plot.errorbar(X,YB,yerr=X**4*Sample['Scale'][0]*np.sqrt(1./(Pattern[0]['wtFactor']*xye[2])),
962                                    ecolor=colors[N%6],picker=3.,clip_on=False)
963                            else:
964                                Plot.errorbar(X,YB,yerr=Sample['Scale'][0]*np.sqrt(1./(Pattern[0]['wtFactor']*xye[2])),
965                                    ecolor=colors[N%6],picker=3.,clip_on=False)
966                        else:
967                            Plot.plot(X,YB,colors[N%6]+'+',picker=3.,clip_on=False)
968                        Plot.plot(X,W,colors[(N+2)%6],picker=False)
969                        Plot.plot(X,ZB,colors[(N+1)%6],picker=False)
970                elif G2frame.Weight and 'PWDR' in plottype:
971                    DY = xye[1]*np.sqrt(xye[2])
972                    DYmax = max(DY)
973                    DZ = xye[3]*np.sqrt(xye[2])
974                    DS = xye[5]*np.sqrt(xye[2])-DYmax*G2frame.delOffset
975                    Plot.plot(X,DY,colors[N%6]+'+',picker=3.,clip_on=False)
976                    Plot.plot(X,DZ,colors[(N+1)%6],picker=False)
977                    Plot.plot(X,DS,colors[(N+3)%6],picker=False)
978                    Plot.axhline(0.,color=wx.BLACK)
979                else:
980                    if G2frame.SubBack:
981                        Plot.plot(Xum,Y-W,colors[N%6]+'+',picker=3.,clip_on=False)
982                        Plot.plot(X,Z-W,colors[(N+1)%6],picker=False)
983                    else:
984                        Plot.plot(Xum,Y,colors[N%6]+'+',picker=3.,clip_on=False)
985                        Plot.plot(X,Z,colors[(N+1)%6],picker=False)
986                        if G2frame.logPlot:
987                            Plot.set_yscale("log",nonposy='mask')
988                            Plot.set_ylim(bottom=np.min(np.trim_zeros(Y))/2.,top=np.max(Y)*2.)
989                    if 'PWDR' in plottype:
990                        Plot.plot(X,W,colors[(N+2)%6],picker=False)
991                        Plot.plot(X,D,colors[(N+3)%6],picker=False)
992                    Plot.axhline(0.,color=wx.BLACK)
993                Page.canvas.SetToolTipString('')
994                if G2frame.PatternTree.GetItemText(PickId) == 'Peak List':
995                    tip = 'On data point: Pick peak - L or R MB. On line: L-move, R-delete'
996                    Page.canvas.SetToolTipString(tip)
997                    data = G2frame.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(G2frame,PatternId, 'Peak List'))
998                    for item in data:
999                        if G2frame.qPlot:
1000                            if 'C' in Parms['Type'][0]:
1001                                Lines.append(Plot.axvline(4*math.pi*sind(item[0]/2.)/wave,color=colors[N%6],picker=2.))
1002                            else:
1003                                Lines.append(Plot.axvline(2*math.pi*difC/item[0],color=colors[N%6],picker=2.))                               
1004                        else:
1005                            Lines.append(Plot.axvline(item[0],color=colors[N%6],picker=2.))
1006                if G2frame.PatternTree.GetItemText(PickId) == 'Limits':
1007                    tip = 'On data point: Lower limit - L MB; Upper limit - R MB. On limit: MB down to move'
1008                    Page.canvas.SetToolTipString(tip)
1009                    data = G2frame.LimitsTable.GetData()
1010            else:
1011                if G2frame.logPlot:
1012                    if 'PWDR' in plottype:
1013                        Plot.semilogy(X,Y,colors[N%6],picker=False,nonposy='mask')
1014                    elif 'SASD' in plottype:
1015                        Plot.semilogy(X,Y,colors[N%6],picker=False,nonposy='mask')
1016                else:
1017                    if 'PWDR' in plottype:
1018                        Plot.plot(X,Y,colors[N%6],picker=False)
1019                    elif 'SASD' in plottype:
1020                        Plot.loglog(X,Y,colors[N%6],picker=False,nonposy='mask')
1021                        Plot.set_ylim(bottom=np.min(np.trim_zeros(Y))/2.,top=np.max(Y)*2.)
1022            if G2frame.logPlot:
1023                Plot.set_ylim(bottom=np.min(np.trim_zeros(Y))/2.,top=np.max(Y)*2.)
1024    if PickId and not G2frame.Contour:
1025        Parms,Parms2 = G2frame.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(G2frame,PatternId, 'Instrument Parameters'))
1026        if 'C' in Parms['Type'][0]:
1027            wave = G2mth.getWave(Parms)
1028        else:
1029            difC = Parms['difC'][1]
1030        if G2frame.PatternTree.GetItemText(PickId) in ['Index Peak List','Unit Cells List']:
1031            peaks = np.array((G2frame.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(G2frame,PatternId, 'Index Peak List'))))
1032            for peak in peaks:
1033                if G2frame.qPlot:
1034                    Plot.axvline(4*np.pi*sind(peak[0]/2.0)/wave,color='b')
1035                else:
1036                    Plot.axvline(peak[0],color='b')
1037            for hkl in G2frame.HKL:
1038                if G2frame.qPlot:
1039                    Plot.axvline(4*np.pi*sind(hkl[5]/2.0)/wave,color='r',dashes=(5,5))
1040                else:
1041                    Plot.axvline(hkl[5],color='r',dashes=(5,5))
1042        elif G2frame.PatternTree.GetItemText(PickId) in ['Reflection Lists'] or \
1043            'PWDR' in G2frame.PatternTree.GetItemText(PickId):
1044            refColors=['b','r','c','g','m','k']
1045            Phases = G2frame.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(G2frame,PatternId,'Reflection Lists'))
1046            for pId,phase in enumerate(Phases):
1047                try:    #patch for old style reflection lists
1048                    peaks = Phases[phase]['RefList']
1049                except TypeError:
1050                    peaks = Phases[phase]
1051                if not len(peaks):
1052                    continue
1053                peak = np.array([[peak[4],peak[5]] for peak in peaks])
1054                pos = G2frame.refOffset-pId*Ymax*G2frame.refDelt*np.ones_like(peak)
1055                if G2frame.qPlot:
1056                    Plot.plot(2*np.pi/peak.T[0],pos,refColors[pId%6]+'|',mew=1,ms=8,picker=3.,label=phase)
1057                else:
1058                    Plot.plot(peak.T[1],pos,refColors[pId%6]+'|',mew=1,ms=8,picker=3.,label=phase)
1059            if len(Phases):
1060                handles,legends = Plot.get_legend_handles_labels()  #got double entries in the legends for some reason
1061                if handles:
1062                    Plot.legend(handles[::2],legends[::2],title='Phases',loc='best')    #skip every other one
1063           
1064    if G2frame.Contour:
1065        acolor = mpl.cm.get_cmap(G2frame.ContourColor)
1066        Img = Plot.imshow(ContourZ,cmap=acolor,vmin=0,vmax=Ymax*G2frame.Cmax,interpolation=G2frame.Interpolate, 
1067            extent=[ContourX[0],ContourX[-1],ContourY[0],ContourY[-1]],aspect='auto',origin='lower')
1068        Page.figure.colorbar(Img)
1069    else:
1070        G2frame.Lines = Lines
1071    if not newPlot:
1072        Page.toolbar.push_current()
1073        Plot.set_xlim(xylim[0])
1074        Plot.set_ylim(xylim[1])
1075#        xylim = []
1076        Page.toolbar.push_current()
1077        Page.toolbar.draw()
1078    else:
1079        Page.canvas.draw()
1080    olderr = np.seterr(invalid='warn') #ugh - this removes a matplotlib error for mouse clicks in log plots                 
1081#    G2frame.Pwdr = True
1082   
1083################################################################################
1084##### PlotDeltSig
1085################################################################################
1086           
1087def PlotDeltSig(G2frame,kind):
1088    'needs a doc string'
1089    try:
1090        plotNum = G2frame.G2plotNB.plotList.index('Error analysis')
1091        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
1092        Page.figure.clf()
1093        Plot = Page.figure.gca()          #get a fresh plot after clf()
1094    except ValueError:
1095        newPlot = True
1096        G2frame.Cmax = 1.0
1097        Plot = G2frame.G2plotNB.addMpl('Error analysis').gca()
1098        plotNum = G2frame.G2plotNB.plotList.index('Error analysis')
1099        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
1100    Page.Choice = None
1101    PatternId = G2frame.PatternId
1102    Pattern = G2frame.PatternTree.GetItemPyData(PatternId)
1103    Pattern.append(G2frame.PatternTree.GetItemText(PatternId))
1104    wtFactor = Pattern[0]['wtFactor']
1105    if kind == 'PWDR':
1106        limits = G2frame.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(G2frame,PatternId, 'Limits'))[1]
1107        xye = np.array(Pattern[1])
1108        xmin = np.searchsorted(xye[0],limits[0])
1109        xmax = np.searchsorted(xye[0],limits[1])
1110        DS = xye[5][xmin:xmax]*np.sqrt(wtFactor*xye[2][xmin:xmax])
1111    elif kind == 'HKLF':
1112        refl = Pattern[1]
1113        DS = []
1114        for ref in refl:
1115            if ref[6] > 0.:
1116                DS.append((ref[5]-ref[7])/ref[6])
1117    Page.SetFocus()
1118    G2frame.G2plotNB.status.DestroyChildren()
1119    DS.sort()
1120    EDS = np.zeros_like(DS)
1121    DX = np.linspace(0.,1.,num=len(DS),endpoint=True)
1122    oldErr = np.seterr(invalid='ignore')    #avoid problem at DX==0
1123    T = np.sqrt(np.log(1.0/DX**2))
1124    top = 2.515517+0.802853*T+0.010328*T**2
1125    bot = 1.0+1.432788*T+0.189269*T**2+0.001308*T**3
1126    EDS = np.where(DX>0,-(T-top/bot),(T-top/bot))
1127    low1 = np.searchsorted(EDS,-1.)
1128    hi1 = np.searchsorted(EDS,1.)
1129    slp,intcp = np.polyfit(EDS[low1:hi1],DS[low1:hi1],deg=1)
1130    frac = 100.*(hi1-low1)/len(DS)
1131    G2frame.G2plotNB.status.SetStatusText(  \
1132        'Over range -1. to 1. :'+' slope = %.3f, intercept = %.3f for %.2f%% of the fitted data'%(slp,intcp,frac),1)
1133    Plot.set_title('Normal probability for '+Pattern[-1])
1134    Plot.set_xlabel(r'expected $\mathsf{\Delta/\sigma}$',fontsize=14)
1135    Plot.set_ylabel(r'observed $\mathsf{\Delta/\sigma}$',fontsize=14)
1136    Plot.plot(EDS,DS,'r+',label='result')
1137    Plot.plot([-2,2],[-2,2],'k',dashes=(5,5),label='ideal')
1138    Plot.legend(loc='upper left')
1139    np.seterr(invalid='warn')
1140    Page.canvas.draw()
1141       
1142################################################################################
1143##### PlotISFG
1144################################################################################
1145           
1146def PlotISFG(G2frame,newPlot=False,type=''):
1147    ''' Plotting package for PDF analysis; displays I(q), S(q), F(q) and G(r) as single
1148    or multiple plots with waterfall and contour plots as options
1149    '''
1150    if not type:
1151        type = G2frame.G2plotNB.plotList[G2frame.G2plotNB.nb.GetSelection()]
1152    if type not in ['I(Q)','S(Q)','F(Q)','G(R)']:
1153        return
1154    superMinusOne = unichr(0xaf)+unichr(0xb9)
1155   
1156    def OnPlotKeyPress(event):
1157        newPlot = False
1158        if event.key == 'u':
1159            if G2frame.Contour:
1160                G2frame.Cmax = min(1.0,G2frame.Cmax*1.2)
1161            elif G2frame.Offset[0] < 100.:
1162                G2frame.Offset[0] += 1.
1163        elif event.key == 'd':
1164            if G2frame.Contour:
1165                G2frame.Cmax = max(0.0,G2frame.Cmax*0.8)
1166            elif G2frame.Offset[0] > 0.:
1167                G2frame.Offset[0] -= 1.
1168        elif event.key == 'l':
1169            G2frame.Offset[1] -= 1.
1170        elif event.key == 'r':
1171            G2frame.Offset[1] += 1.
1172        elif event.key == 'o':
1173            G2frame.Offset = [0,0]
1174        elif event.key == 'c':
1175            newPlot = True
1176            G2frame.Contour = not G2frame.Contour
1177            if not G2frame.Contour:
1178                G2frame.SinglePlot = False
1179                G2frame.Offset = [0.,0.]
1180        elif event.key == 's':
1181            if G2frame.Contour:
1182                choice = [m for m in mpl.cm.datad.keys() if not m.endswith("_r")]
1183                choice.sort()
1184                dlg = wx.SingleChoiceDialog(G2frame,'Select','Color scheme',choice)
1185                if dlg.ShowModal() == wx.ID_OK:
1186                    sel = dlg.GetSelection()
1187                    G2frame.ContourColor = choice[sel]
1188                else:
1189                    G2frame.ContourColor = 'Paired'
1190                dlg.Destroy()
1191            else:
1192                G2frame.SinglePlot = not G2frame.SinglePlot               
1193        elif event.key == 'i':                  #for smoothing contour plot
1194            choice = ['nearest','bilinear','bicubic','spline16','spline36','hanning',
1195               'hamming','hermite','kaiser','quadric','catrom','gaussian','bessel',
1196               'mitchell','sinc','lanczos']
1197            dlg = wx.SingleChoiceDialog(G2frame,'Select','Interpolation',choice)
1198            if dlg.ShowModal() == wx.ID_OK:
1199                sel = dlg.GetSelection()
1200                G2frame.Interpolate = choice[sel]
1201            else:
1202                G2frame.Interpolate = 'nearest'
1203            dlg.Destroy()
1204        elif event.key == 't' and not G2frame.Contour:
1205            G2frame.Legend = not G2frame.Legend
1206        PlotISFG(G2frame,newPlot=newPlot,type=type)
1207       
1208    def OnKeyBox(event):
1209        if G2frame.G2plotNB.nb.GetSelection() == G2frame.G2plotNB.plotList.index(type):
1210            event.key = cb.GetValue()[0]
1211            cb.SetValue(' key press')
1212            wx.CallAfter(OnPlotKeyPress,event)
1213        Page.canvas.SetFocus() # redirect the Focus from the button back to the plot
1214                       
1215    def OnMotion(event):
1216        xpos = event.xdata
1217        if xpos:                                        #avoid out of frame mouse position
1218            ypos = event.ydata
1219            Page.canvas.SetCursor(wx.CROSS_CURSOR)
1220            try:
1221                if G2frame.Contour:
1222                    G2frame.G2plotNB.status.SetStatusText('R =%.3fA pattern ID =%5d'%(xpos,int(ypos)),1)
1223                else:
1224                    G2frame.G2plotNB.status.SetStatusText('R =%.3fA %s =%.2f'%(xpos,type,ypos),1)                   
1225            except TypeError:
1226                G2frame.G2plotNB.status.SetStatusText('Select '+type+' pattern first',1)
1227   
1228    xylim = []
1229    try:
1230        plotNum = G2frame.G2plotNB.plotList.index(type)
1231        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
1232        if not newPlot:
1233            Plot = Page.figure.gca()          #get previous plot & get limits
1234            xylim = Plot.get_xlim(),Plot.get_ylim()
1235        Page.figure.clf()
1236        Plot = Page.figure.gca()
1237    except ValueError:
1238        newPlot = True
1239        G2frame.Cmax = 1.0
1240        Plot = G2frame.G2plotNB.addMpl(type).gca()
1241        plotNum = G2frame.G2plotNB.plotList.index(type)
1242        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
1243        Page.canvas.mpl_connect('key_press_event', OnPlotKeyPress)
1244        Page.canvas.mpl_connect('motion_notify_event', OnMotion)
1245   
1246    Page.SetFocus()
1247    G2frame.G2plotNB.status.DestroyChildren()
1248    if G2frame.Contour:
1249        Page.Choice = (' key press','d: lower contour max','u: raise contour max',
1250            'i: interpolation method','s: color scheme','c: contour off')
1251    else:
1252        Page.Choice = (' key press','l: offset left','r: offset right','d: offset down','u: offset up',
1253            'o: reset offset','t: toggle legend','c: contour on','s: toggle single plot')
1254    Page.keyPress = OnPlotKeyPress
1255    PatternId = G2frame.PatternId
1256    PickId = G2frame.PickId
1257    Plot.set_title(type)
1258    if type == 'G(R)':
1259        Plot.set_xlabel(r'$R,\AA$',fontsize=14)
1260    else:
1261        Plot.set_xlabel(r'$Q,\AA$'+superMinusOne,fontsize=14)
1262    Plot.set_ylabel(r''+type,fontsize=14)
1263    colors=['b','g','r','c','m','k']
1264    name = G2frame.PatternTree.GetItemText(PatternId)[4:]
1265    Pattern = []   
1266    if G2frame.SinglePlot:
1267        name = G2frame.PatternTree.GetItemText(PatternId)
1268        name = type+name[4:]
1269        Id = G2gd.GetPatternTreeItemId(G2frame,PatternId,name)
1270        Pattern = G2frame.PatternTree.GetItemPyData(Id)
1271        if Pattern:
1272            Pattern.append(name)
1273        PlotList = [Pattern,]
1274    else:
1275        PlotList = []
1276        item, cookie = G2frame.PatternTree.GetFirstChild(G2frame.root)
1277        while item:
1278            if 'PDF' in G2frame.PatternTree.GetItemText(item):
1279                name = type+G2frame.PatternTree.GetItemText(item)[4:]
1280                Id = G2gd.GetPatternTreeItemId(G2frame,item,name)
1281                Pattern = G2frame.PatternTree.GetItemPyData(Id)
1282                if Pattern:
1283                    Pattern.append(name)
1284                    PlotList.append(Pattern)
1285            item, cookie = G2frame.PatternTree.GetNextChild(G2frame.root, cookie)
1286    PDFdata = G2frame.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(G2frame,PatternId, 'PDF Controls'))
1287    numbDen = G2pwd.GetNumDensity(PDFdata['ElList'],PDFdata['Form Vol'])
1288    Xb = [0.,10.]
1289    Yb = [0.,-40.*np.pi*numbDen]
1290    Ymax = 1.0
1291    lenX = 0
1292    for Pattern in PlotList:
1293        xye = Pattern[1]
1294        Ymax = max(Ymax,max(xye[1]))
1295    offset = G2frame.Offset[0]*Ymax/100.0
1296    if G2frame.Contour:
1297        ContourZ = []
1298        ContourY = []
1299        Nseq = 0
1300    for N,Pattern in enumerate(PlotList):
1301        xye = Pattern[1]
1302        if PickId:
1303            ifpicked = Pattern[2] == G2frame.PatternTree.GetItemText(PatternId)
1304        X = xye[0]
1305        if not lenX:
1306            lenX = len(X)           
1307        Y = xye[1]+offset*N
1308        if G2frame.Contour:
1309            if lenX == len(X):
1310                ContourY.append(N)
1311                ContourZ.append(Y)
1312                ContourX = X
1313                Nseq += 1
1314                Plot.set_ylabel('Data sequence',fontsize=12)
1315        else:
1316            X = xye[0]+G2frame.Offset[1]*.005*N
1317            if ifpicked:
1318                Plot.plot(X,Y,colors[N%6]+'+',picker=3.,clip_on=False)
1319                Page.canvas.SetToolTipString('')
1320            else:
1321                if G2frame.Legend:
1322                    Plot.plot(X,Y,colors[N%6],picker=False,label='Azm:'+Pattern[2].split('=')[1])
1323                else:
1324                    Plot.plot(X,Y,colors[N%6],picker=False)
1325            if type == 'G(R)':
1326                Plot.plot(Xb,Yb,color='k',dashes=(5,5))
1327            elif type == 'F(Q)':
1328                Plot.axhline(0.,color=wx.BLACK)
1329            elif type == 'S(Q)':
1330                Plot.axhline(1.,color=wx.BLACK)
1331    if G2frame.Contour:
1332        acolor = mpl.cm.get_cmap(G2frame.ContourColor)
1333        Img = Plot.imshow(ContourZ,cmap=acolor,vmin=0,vmax=Ymax*G2frame.Cmax,interpolation=G2frame.Interpolate, 
1334            extent=[ContourX[0],ContourX[-1],ContourY[0],ContourY[-1]],aspect='auto',origin='lower')
1335        Page.figure.colorbar(Img)
1336    elif G2frame.Legend:
1337        Plot.legend(loc='best')
1338    if not newPlot:
1339        Page.toolbar.push_current()
1340        Plot.set_xlim(xylim[0])
1341        Plot.set_ylim(xylim[1])
1342        xylim = []
1343        Page.toolbar.push_current()
1344        Page.toolbar.draw()
1345    else:
1346        Page.canvas.draw()
1347       
1348################################################################################
1349##### PlotXY
1350################################################################################
1351           
1352def PlotXY(G2frame,XY,XY2=None,labelX=None,labelY=None,newPlot=False,type=''):
1353    '''simple plot of xy data, used for diagnostic purposes
1354    '''
1355    def OnMotion(event):
1356        xpos = event.xdata
1357        if xpos:                                        #avoid out of frame mouse position
1358            ypos = event.ydata
1359            Page.canvas.SetCursor(wx.CROSS_CURSOR)
1360            try:
1361                G2frame.G2plotNB.status.SetStatusText('X =%9.3f %s =%9.3f'%(xpos,type,ypos),1)                   
1362            except TypeError:
1363                G2frame.G2plotNB.status.SetStatusText('Select '+type+' pattern first',1)
1364
1365    try:
1366        plotNum = G2frame.G2plotNB.plotList.index(type)
1367        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
1368        if not newPlot:
1369            Plot = Page.figure.gca()
1370            xylim = Plot.get_xlim(),Plot.get_ylim()
1371        Page.figure.clf()
1372        Plot = Page.figure.gca()
1373    except ValueError:
1374        newPlot = True
1375        Plot = G2frame.G2plotNB.addMpl(type).gca()
1376        plotNum = G2frame.G2plotNB.plotList.index(type)
1377        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
1378        Page.canvas.mpl_connect('motion_notify_event', OnMotion)
1379   
1380    Page.Choice = None
1381    Page.SetFocus()
1382    G2frame.G2plotNB.status.DestroyChildren()
1383    Plot.set_title(type)
1384    if xLabel:
1385        Plot.set_xlabel(r''+xLabel,fontsize=14)
1386    else:
1387        Plot.set_xlabel(r'X',fontsize=14)
1388    if yLabel:
1389        Plot.set_ylabel(r''+yLabel,fontsize=14)
1390    else:
1391        Plot.set_ylabel(r'Y',fontsize=14)
1392    colors=['b','g','r','c','m','k']
1393    for ixy,xy in enumerate(XY):
1394        X,Y = xy
1395        Plot.plot(X,Y,color(ixy%6)+'+',picker=False)
1396    if len(XY2):
1397        for ixy,xy in enumerate(XY2):
1398            X,Y = xy
1399            Plot.plot(X,Y,color(ixy%6),picker=False)
1400    if not newPlot:
1401        Page.toolbar.push_current()
1402        Plot.set_xlim(xylim[0])
1403        Plot.set_ylim(xylim[1])
1404        xylim = []
1405        Page.toolbar.push_current()
1406        Page.toolbar.draw()
1407    else:
1408        Page.canvas.draw()
1409
1410################################################################################
1411##### PlotStrain
1412################################################################################
1413           
1414def PlotStrain(G2frame,data,newPlot=False):
1415    '''plot of strain data, used for diagnostic purposes
1416    '''
1417    def OnMotion(event):
1418        xpos = event.xdata
1419        if xpos:                                        #avoid out of frame mouse position
1420            ypos = event.ydata
1421            Page.canvas.SetCursor(wx.CROSS_CURSOR)
1422            try:
1423                G2frame.G2plotNB.status.SetStatusText('d-spacing =%9.5f Azimuth =%9.3f'%(ypos,xpos),1)                   
1424            except TypeError:
1425                G2frame.G2plotNB.status.SetStatusText('Select Strain pattern first',1)
1426
1427    try:
1428        plotNum = G2frame.G2plotNB.plotList.index('Strain')
1429        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
1430        if not newPlot:
1431            Plot = Page.figure.gca()
1432            xylim = Plot.get_xlim(),Plot.get_ylim()
1433        Page.figure.clf()
1434        Plot = Page.figure.gca()
1435    except ValueError:
1436        newPlot = True
1437        Plot = G2frame.G2plotNB.addMpl('Strain').gca()
1438        plotNum = G2frame.G2plotNB.plotList.index('Strain')
1439        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
1440        Page.canvas.mpl_connect('motion_notify_event', OnMotion)
1441   
1442    Page.Choice = None
1443    Page.SetFocus()
1444    G2frame.G2plotNB.status.DestroyChildren()
1445    Plot.set_title('Strain')
1446    Plot.set_ylabel(r'd-spacing',fontsize=14)
1447    Plot.set_xlabel(r'Azimuth',fontsize=14)
1448    colors=['b','g','r','c','m','k']
1449    for N,item in enumerate(data['d-zero']):
1450        Y,X = np.array(item['ImtaObs'])         #plot azimuth as X & d-spacing as Y
1451        Plot.plot(X,Y,colors[N%6]+'+',picker=False)
1452        Y,X = np.array(item['ImtaCalc'])
1453        Plot.plot(X,Y,colors[N%6],picker=False)
1454    if not newPlot:
1455        Page.toolbar.push_current()
1456        Plot.set_xlim(xylim[0])
1457        Plot.set_ylim(xylim[1])
1458        xylim = []
1459        Page.toolbar.push_current()
1460        Page.toolbar.draw()
1461    else:
1462        Page.canvas.draw()
1463       
1464################################################################################
1465##### PlotSASDSizeDist
1466################################################################################
1467           
1468def PlotSASDSizeDist(G2frame):
1469   
1470    def OnPageChanged(event):
1471        PlotText = G2frame.G2plotNB.nb.GetPageText(G2frame.G2plotNB.nb.GetSelection())
1472        if 'Powder' in PlotText:
1473            PlotPatterns(G2frame,plotType='SASD',newPlot=True)
1474        elif 'Size' in PlotText:
1475            PlotSASDSizeDist(G2frame)
1476   
1477    def OnMotion(event):
1478        xpos = event.xdata
1479        if xpos:                                        #avoid out of frame mouse position
1480            ypos = event.ydata
1481            Page.canvas.SetCursor(wx.CROSS_CURSOR)
1482            try:
1483                G2frame.G2plotNB.status.SetStatusText('diameter =%9.3f f(D) =%9.3g'%(xpos,ypos),1)                   
1484            except TypeError:
1485                G2frame.G2plotNB.status.SetStatusText('Select Strain pattern first',1)
1486
1487    try:
1488        plotNum = G2frame.G2plotNB.plotList.index('Size Distribution')
1489        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
1490        Page.figure.clf()
1491        Plot = Page.figure.gca()          #get a fresh plot after clf()
1492    except ValueError:
1493        newPlot = True
1494        Plot = G2frame.G2plotNB.addMpl('Size Distribution').gca()
1495        G2frame.G2plotNB.Bind(wx.aui.EVT_AUINOTEBOOK_PAGE_CHANGED,OnPageChanged)
1496        plotNum = G2frame.G2plotNB.plotList.index('Size Distribution')
1497        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
1498        Page.canvas.mpl_connect('motion_notify_event', OnMotion)
1499    Page.Choice = None
1500    Page.SetFocus()
1501    PatternId = G2frame.PatternId
1502    data = G2frame.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(G2frame,PatternId, 'Models'))
1503    Bins,Dbins,BinMag = data['Size']['Distribution']
1504    Plot.set_title('Size Distribution')
1505    Plot.set_xlabel(r'$D, \AA$',fontsize=14)
1506    Plot.set_ylabel(r'$Volume distribution f(D)$',fontsize=14)
1507    if data['Size']['logBins']:
1508        Plot.set_xscale("log",nonposy='mask')
1509        Plot.set_xlim([np.min(2.*Bins)/2.,np.max(2.*Bins)*2.])
1510    Plot.bar(2.*Bins-Dbins,BinMag,2.*Dbins,facecolor='white')       #plot diameters
1511    if 'Size Calc' in data:
1512        Rbins,Dist = data['Size Calc']
1513        for i in range(len(Rbins)):
1514            if len(Rbins[i]):
1515                Plot.plot(2.*Rbins[i],Dist[i])       #plot diameters
1516    Page.canvas.draw()
1517
1518################################################################################
1519##### PlotPowderLines
1520################################################################################
1521           
1522def PlotPowderLines(G2frame):
1523    ''' plotting of powder lines (i.e. no powder pattern) as sticks
1524    '''
1525    global HKL
1526
1527    def OnMotion(event):
1528        xpos = event.xdata
1529        if xpos:                                        #avoid out of frame mouse position
1530            Page.canvas.SetCursor(wx.CROSS_CURSOR)
1531            G2frame.G2plotNB.status.SetFields(['','2-theta =%9.3f '%(xpos,)])
1532            if G2frame.PickId and G2frame.PatternTree.GetItemText(G2frame.PickId) in ['Index Peak List','Unit Cells List']:
1533                found = []
1534                if len(HKL):
1535                    view = Page.toolbar._views.forward()[0][:2]
1536                    wid = view[1]-view[0]
1537                    found = HKL[np.where(np.fabs(HKL.T[5]-xpos) < 0.002*wid)]
1538                if len(found):
1539                    h,k,l = found[0][:3] 
1540                    Page.canvas.SetToolTipString('%d,%d,%d'%(int(h),int(k),int(l)))
1541                else:
1542                    Page.canvas.SetToolTipString('')
1543
1544    try:
1545        plotNum = G2frame.G2plotNB.plotList.index('Powder Lines')
1546        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
1547        Page.figure.clf()
1548        Plot = Page.figure.gca()
1549    except ValueError:
1550        Plot = G2frame.G2plotNB.addMpl('Powder Lines').gca()
1551        plotNum = G2frame.G2plotNB.plotList.index('Powder Lines')
1552        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
1553        Page.canvas.mpl_connect('motion_notify_event', OnMotion)
1554       
1555    Page.Choice = None
1556    Page.SetFocus()
1557    Plot.set_title('Powder Pattern Lines')
1558    Plot.set_xlabel(r'$\mathsf{2\theta}$',fontsize=14)
1559    PickId = G2frame.PickId
1560    PatternId = G2frame.PatternId
1561    peaks = G2frame.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(G2frame,PatternId, 'Index Peak List'))
1562    for peak in peaks:
1563        Plot.axvline(peak[0],color='b')
1564    HKL = np.array(G2frame.HKL)
1565    for hkl in G2frame.HKL:
1566        Plot.axvline(hkl[5],color='r',dashes=(5,5))
1567    xmin = peaks[0][0]
1568    xmax = peaks[-1][0]
1569    delt = xmax-xmin
1570    xlim = [max(0,xmin-delt/20.),min(180.,xmax+delt/20.)]
1571    Plot.set_xlim(xlim)
1572    Page.canvas.draw()
1573    Page.toolbar.push_current()
1574
1575################################################################################
1576##### PlotPeakWidths
1577################################################################################
1578           
1579def PlotPeakWidths(G2frame):
1580    ''' Plotting of instrument broadening terms as function of 2-theta
1581    Seen when "Instrument Parameters" chosen from powder pattern data tree
1582    '''
1583#    sig = lambda Th,U,V,W: 1.17741*math.sqrt(U*tand(Th)**2+V*tand(Th)+W)*math.pi/18000.
1584#    gam = lambda Th,X,Y: (X/cosd(Th)+Y*tand(Th))*math.pi/18000.
1585    gamFW = lambda s,g: np.exp(np.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.)
1586#    gamFW2 = lambda s,g: math.sqrt(s**2+(0.4654996*g)**2)+.5345004*g  #Ubaldo Bafile - private communication
1587    PatternId = G2frame.PatternId
1588    limitID = G2gd.GetPatternTreeItemId(G2frame,PatternId, 'Limits')
1589    if limitID:
1590        limits = G2frame.PatternTree.GetItemPyData(limitID)[:2]
1591    else:
1592        return
1593    Parms,Parms2 = G2frame.PatternTree.GetItemPyData( \
1594        G2gd.GetPatternTreeItemId(G2frame,PatternId, 'Instrument Parameters'))
1595    if 'C' in Parms['Type'][0]:
1596        lam = G2mth.getWave(Parms)
1597    else:
1598        difC = Parms['difC'][0]
1599    peakID = G2gd.GetPatternTreeItemId(G2frame,PatternId, 'Peak List')
1600    if peakID:
1601        peaks = G2frame.PatternTree.GetItemPyData(peakID)
1602    else:
1603        peaks = []
1604   
1605    try:
1606        plotNum = G2frame.G2plotNB.plotList.index('Peak Widths')
1607        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
1608        Page.figure.clf()
1609        Plot = Page.figure.gca()
1610    except ValueError:
1611        Plot = G2frame.G2plotNB.addMpl('Peak Widths').gca()
1612        plotNum = G2frame.G2plotNB.plotList.index('Peak Widths')
1613        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
1614    Page.Choice = None
1615    Page.SetFocus()
1616   
1617    Page.canvas.SetToolTipString('')
1618    colors=['b','g','r','c','m','k']
1619    X = []
1620    Y = []
1621    Z = []
1622    W = []
1623    if 'C' in Parms['Type'][0]:
1624        Plot.set_title('Instrument and sample peak widths')
1625        Plot.set_xlabel(r'$Q, \AA^{-1}$',fontsize=14)
1626        Plot.set_ylabel(r'$\Delta Q/Q, \Delta d/d$',fontsize=14)
1627        try:
1628            Xmin,Xmax = limits[1]
1629            X = np.linspace(Xmin,Xmax,num=101,endpoint=True)
1630            Q = 4.*np.pi*npsind(X/2.)/lam
1631            Z = np.ones_like(X)
1632            data = G2mth.setPeakparms(Parms,Parms2,X,Z)
1633            s = 1.17741*np.sqrt(data[4])*np.pi/18000.
1634            g = data[6]*np.pi/18000.
1635            G = G2pwd.getgamFW(g,s)
1636            Y = s/nptand(X/2.)
1637            Z = g/nptand(X/2.)
1638            W = G/nptand(X/2.)
1639            Plot.plot(Q,Y,color='r',label='Gaussian')
1640            Plot.plot(Q,Z,color='g',label='Lorentzian')
1641            Plot.plot(Q,W,color='b',label='G+L')
1642           
1643            fit = G2mth.setPeakparms(Parms,Parms2,X,Z,useFit=True)
1644            sf = 1.17741*np.sqrt(fit[4])*np.pi/18000.
1645            gf = fit[6]*np.pi/18000.
1646            Gf = G2pwd.getgamFW(gf,sf)
1647            Yf = sf/nptand(X/2.)
1648            Zf = gf/nptand(X/2.)
1649            Wf = Gf/nptand(X/2.)
1650            Plot.plot(Q,Yf,color='r',dashes=(5,5),label='Gaussian fit')
1651            Plot.plot(Q,Zf,color='g',dashes=(5,5),label='Lorentzian fit')
1652            Plot.plot(Q,Wf,color='b',dashes=(5,5),label='G+L fit')
1653           
1654            X = []
1655            Y = []
1656            Z = []
1657            W = []
1658            V = []
1659            for peak in peaks:
1660                X.append(4.0*math.pi*sind(peak[0]/2.0)/lam)
1661                try:
1662                    s = 1.17741*math.sqrt(peak[4])*math.pi/18000.
1663                except ValueError:
1664                    s = 0.01
1665                g = peak[6]*math.pi/18000.
1666                G = G2pwd.getgamFW(g,s)
1667                Y.append(s/tand(peak[0]/2.))
1668                Z.append(g/tand(peak[0]/2.))
1669                W.append(G/tand(peak[0]/2.))
1670            if len(peaks):
1671                Plot.plot(X,Y,'+',color='r',label='G peak')
1672                Plot.plot(X,Z,'+',color='g',label='L peak')
1673                Plot.plot(X,W,'+',color='b',label='G+L peak')
1674            Plot.legend(loc='best')
1675            Page.canvas.draw()
1676        except ValueError:
1677            print '**** ERROR - default U,V,W profile coefficients yield sqrt of negative value at 2theta =', \
1678                '%.3f'%(2*theta)
1679            G2frame.G2plotNB.Delete('Peak Widths')
1680    else:
1681        Plot.set_title('Instrument and sample peak coefficients')
1682        Plot.set_xlabel(r'$Q, \AA^{-1}$',fontsize=14)
1683        Plot.set_ylabel(r'$\alpha, \beta, \Delta Q/Q, \Delta d/d$',fontsize=14)
1684        Xmin,Xmax = limits[1]
1685        T = np.linspace(Xmin,Xmax,num=101,endpoint=True)
1686        Z = np.ones_like(T)
1687        data = G2mth.setPeakparms(Parms,Parms2,T,Z)
1688        ds = T/difC
1689        Q = 2.*np.pi/ds
1690        A = data[4]
1691        B = data[6]
1692        S = 1.17741*np.sqrt(data[8])/T
1693        G = data[10]/T
1694        Plot.plot(Q,A,color='r',label='Alpha')
1695        Plot.plot(Q,B,color='g',label='Beta')
1696        Plot.plot(Q,S,color='b',label='Gaussian')
1697        Plot.plot(Q,G,color='m',label='Lorentzian')
1698
1699        fit = G2mth.setPeakparms(Parms,Parms2,T,Z)
1700        ds = T/difC
1701        Q = 2.*np.pi/ds
1702        Af = fit[4]
1703        Bf = fit[6]
1704        Sf = 1.17741*np.sqrt(fit[8])/T
1705        Gf = fit[10]/T
1706        Plot.plot(Q,Af,color='r',dashes=(5,5),label='Alpha fit')
1707        Plot.plot(Q,Bf,color='g',dashes=(5,5),label='Beta fit')
1708        Plot.plot(Q,Sf,color='b',dashes=(5,5),label='Gaussian fit')
1709        Plot.plot(Q,Gf,color='m',dashes=(5,5),label='Lorentzian fit')
1710       
1711        T = []
1712        A = []
1713        B = []
1714        S = []
1715        G = []
1716        W = []
1717        Q = []
1718        V = []
1719        for peak in peaks:
1720            T.append(peak[0])
1721            A.append(peak[4])
1722            B.append(peak[6])
1723            Q.append(2.*np.pi*difC/peak[0])
1724            S.append(1.17741*np.sqrt(peak[8])/peak[0])
1725            G.append(peak[10]/peak[0])
1726           
1727       
1728        Plot.plot(Q,A,'+',color='r',label='Alpha peak')
1729        Plot.plot(Q,B,'+',color='g',label='Beta peak')
1730        Plot.plot(Q,S,'+',color='b',label='Gaussian peak')
1731        Plot.plot(Q,G,'+',color='m',label='Lorentzian peak')
1732        Plot.legend(loc='best')
1733        Page.canvas.draw()
1734
1735   
1736################################################################################
1737##### PlotSizeStrainPO
1738################################################################################
1739           
1740def PlotSizeStrainPO(G2frame,data,Start=False):
1741    '''Plot 3D mustrain/size/preferred orientation figure. In this instance data is for a phase
1742    '''
1743   
1744    PatternId = G2frame.PatternId
1745    generalData = data['General']
1746    SGData = generalData['SGData']
1747    SGLaue = SGData['SGLaue']
1748    if Start:                   #initialize the spherical harmonics qlmn arrays
1749        ptx.pyqlmninit()
1750        Start = False
1751#    MuStrKeys = G2spc.MustrainNames(SGData)
1752    cell = generalData['Cell'][1:]
1753    A,B = G2lat.cell2AB(cell[:6])
1754    Vol = cell[6]
1755    useList = data['Histograms']
1756    phase = generalData['Name']
1757    plotType = generalData['Data plot type']
1758    plotDict = {'Mustrain':'Mustrain','Size':'Size','Preferred orientation':'Pref.Ori.'}
1759    for ptype in plotDict:
1760        G2frame.G2plotNB.Delete(ptype)
1761    if plotType in ['None']:
1762        return       
1763
1764    for item in useList:
1765        if useList[item]['Show']:
1766            break
1767    else:
1768        return            #nothing to show!!
1769   
1770    numPlots = len(useList)
1771
1772    if plotType in ['Mustrain','Size']:
1773        Plot = mp3d.Axes3D(G2frame.G2plotNB.add3D(plotType))
1774    else:
1775        Plot = G2frame.G2plotNB.addMpl(plotType).gca()       
1776    plotNum = G2frame.G2plotNB.plotList.index(plotType)
1777    Page = G2frame.G2plotNB.nb.GetPage(plotNum)
1778    Page.Choice = None
1779    Page.SetFocus()
1780    G2frame.G2plotNB.status.SetStatusText('',1)
1781    if not Page.IsShown():
1782        Page.Show()
1783   
1784    for item in useList:
1785        if useList[item]['Show']:
1786            PHI = np.linspace(0.,360.,30,True)
1787            PSI = np.linspace(0.,180.,30,True)
1788            X = np.outer(npsind(PHI),npsind(PSI))
1789            Y = np.outer(npcosd(PHI),npsind(PSI))
1790            Z = np.outer(np.ones(np.size(PHI)),npcosd(PSI))
1791            try:        #temp patch instead of 'mustrain' for old files with 'microstrain'
1792                coeff = useList[item][plotDict[plotType]]
1793            except KeyError:
1794                break
1795            if plotType in ['Mustrain','Size']:
1796                if coeff[0] == 'isotropic':
1797                    X *= coeff[1][0]
1798                    Y *= coeff[1][0]
1799                    Z *= coeff[1][0]                               
1800                elif coeff[0] == 'uniaxial':
1801                   
1802                    def uniaxCalc(xyz,iso,aniso,axes):
1803                        Z = np.array(axes)
1804                        cp = abs(np.dot(xyz,Z))
1805                        sp = np.sqrt(1.-cp**2)
1806                        R = iso*aniso/np.sqrt((iso*cp)**2+(aniso*sp)**2)
1807                        return R*xyz
1808                       
1809                    iso,aniso = coeff[1][:2]
1810                    axes = np.inner(A,np.array(coeff[3]))
1811                    axes /= nl.norm(axes)
1812                    Shkl = np.array(coeff[1])
1813                    XYZ = np.dstack((X,Y,Z))
1814                    XYZ = np.nan_to_num(np.apply_along_axis(uniaxCalc,2,XYZ,iso,aniso,axes))
1815                    X,Y,Z = np.dsplit(XYZ,3)
1816                    X = X[:,:,0]
1817                    Y = Y[:,:,0]
1818                    Z = Z[:,:,0]
1819               
1820                elif coeff[0] == 'ellipsoidal':
1821                   
1822                    def ellipseCalc(xyz,E,R):
1823                        XYZ = xyz*E.T
1824                        return np.inner(XYZ.T,R)
1825                       
1826                    S6 = coeff[4]
1827                    Sij = G2lat.U6toUij(S6)
1828                    E,R = nl.eigh(Sij)
1829                    XYZ = np.dstack((X,Y,Z))
1830                    XYZ = np.nan_to_num(np.apply_along_axis(ellipseCalc,2,XYZ,E,R))
1831                    X,Y,Z = np.dsplit(XYZ,3)
1832                    X = X[:,:,0]
1833                    Y = Y[:,:,0]
1834                    Z = Z[:,:,0]
1835                   
1836                elif coeff[0] == 'generalized':
1837                   
1838                    def genMustrain(xyz,SGData,A,Shkl):
1839                        uvw = np.inner(A.T,xyz)
1840                        Strm = np.array(G2spc.MustrainCoeff(uvw,SGData))
1841                        sum = np.sum(np.multiply(Shkl,Strm))
1842                        sum = np.where(sum > 0.01,sum,0.01)
1843                        sum = np.sqrt(sum)*math.pi/.18      #centidegrees to radians!
1844                        return sum*xyz
1845                       
1846                    Shkl = np.array(coeff[4])
1847                    if np.any(Shkl):
1848                        XYZ = np.dstack((X,Y,Z))
1849                        XYZ = np.nan_to_num(np.apply_along_axis(genMustrain,2,XYZ,SGData,A,Shkl))
1850                        X,Y,Z = np.dsplit(XYZ,3)
1851                        X = X[:,:,0]
1852                        Y = Y[:,:,0]
1853                        Z = Z[:,:,0]
1854                           
1855                if np.any(X) and np.any(Y) and np.any(Z):
1856                    errFlags = np.seterr(all='ignore')
1857                    Plot.plot_surface(X,Y,Z,rstride=1,cstride=1,color='g',linewidth=1)
1858                    np.seterr(all='ignore')
1859                    xyzlim = np.array([Plot.get_xlim3d(),Plot.get_ylim3d(),Plot.get_zlim3d()]).T
1860                    XYZlim = [min(xyzlim[0]),max(xyzlim[1])]
1861                    Plot.set_xlim3d(XYZlim)
1862                    Plot.set_ylim3d(XYZlim)
1863                    Plot.set_zlim3d(XYZlim)
1864                    Plot.set_aspect('equal')
1865                if plotType == 'Size':
1866                    Plot.set_title('Crystallite size for '+phase+'\n'+coeff[0]+' model')
1867                    Plot.set_xlabel(r'X, $\mu$m')
1868                    Plot.set_ylabel(r'Y, $\mu$m')
1869                    Plot.set_zlabel(r'Z, $\mu$m')
1870                else:   
1871                    Plot.set_title(r'$\mu$strain for '+phase+'\n'+coeff[0]+' model')
1872                    Plot.set_xlabel(r'X, $\mu$strain')
1873                    Plot.set_ylabel(r'Y, $\mu$strain')
1874                    Plot.set_zlabel(r'Z, $\mu$strain')
1875            else:
1876                h,k,l = generalData['POhkl']
1877                if coeff[0] == 'MD':
1878                    print 'March-Dollase preferred orientation plot'
1879               
1880                else:
1881                    PH = np.array(generalData['POhkl'])
1882                    phi,beta = G2lat.CrsAng(PH,cell[:6],SGData)
1883                    SHCoef = {}
1884                    for item in coeff[5]:
1885                        L,N = eval(item.strip('C'))
1886                        SHCoef['C%d,0,%d'%(L,N)] = coeff[5][item]                       
1887                    ODFln = G2lat.Flnh(Start,SHCoef,phi,beta,SGData)
1888                    X = np.linspace(0,90.0,26)
1889                    Y = G2lat.polfcal(ODFln,'0',X,0.0)
1890                    Plot.plot(X,Y,color='k',label=str(PH))
1891                    Plot.legend(loc='best')
1892                    Plot.set_title('Axial distribution for HKL='+str(PH)+' in '+phase)
1893                    Plot.set_xlabel(r'$\psi$',fontsize=16)
1894                    Plot.set_ylabel('MRD',fontsize=14)
1895    Page.canvas.draw()
1896   
1897################################################################################
1898##### PlotTexture
1899################################################################################
1900           
1901def PlotTexture(G2frame,data,Start=False):
1902    '''Pole figure, inverse pole figure, 3D pole distribution and 3D inverse pole distribution
1903    plotting.
1904    dict generalData contains all phase info needed which is in data
1905    '''
1906
1907    shModels = ['cylindrical','none','shear - 2/m','rolling - mmm']
1908    SamSym = dict(zip(shModels,['0','-1','2/m','mmm']))
1909    PatternId = G2frame.PatternId
1910    generalData = data['General']
1911    SGData = generalData['SGData']
1912    textureData = generalData['SH Texture']
1913    G2frame.G2plotNB.Delete('Texture')
1914    if not textureData['Order']:
1915        return                  #no plot!!
1916    SHData = generalData['SH Texture']
1917    SHCoef = SHData['SH Coeff'][1]
1918    cell = generalData['Cell'][1:7]
1919    Amat,Bmat = G2lat.cell2AB(cell)
1920    sq2 = 1.0/math.sqrt(2.0)
1921   
1922    def rp2xyz(r,p):
1923        z = npcosd(r)
1924        xy = np.sqrt(1.-z**2)
1925        return xy*npsind(p),xy*npcosd(p),z
1926           
1927    def OnMotion(event):
1928        SHData = data['General']['SH Texture']
1929        if event.xdata and event.ydata:                 #avoid out of frame errors
1930            xpos = event.xdata
1931            ypos = event.ydata
1932            if 'Inverse' in SHData['PlotType']:
1933                r = xpos**2+ypos**2
1934                if r <= 1.0:
1935                    if 'equal' in G2frame.Projection: 
1936                        r,p = 2.*npasind(np.sqrt(r)*sq2),npatan2d(ypos,xpos)
1937                    else:
1938                        r,p = 2.*npatand(np.sqrt(r)),npatan2d(ypos,xpos)
1939                    ipf = G2lat.invpolfcal(IODFln,SGData,np.array([r,]),np.array([p,]))
1940                    xyz = np.inner(Amat.T,np.array([rp2xyz(r,p)]))
1941                    y,x,z = list(xyz/np.max(np.abs(xyz)))
1942                   
1943                    G2frame.G2plotNB.status.SetFields(['',
1944                        'psi =%9.3f, beta =%9.3f, MRD =%9.3f hkl=%5.2f,%5.2f,%5.2f'%(r,p,ipf,x,y,z)])
1945                                   
1946            elif 'Axial' in SHData['PlotType']:
1947                pass
1948               
1949            else:                       #ordinary pole figure
1950                z = xpos**2+ypos**2
1951                if z <= 1.0:
1952                    z = np.sqrt(z)
1953                    if 'equal' in G2frame.Projection: 
1954                        r,p = 2.*npasind(z*sq2),npatan2d(ypos,xpos)
1955                    else:
1956                        r,p = 2.*npatand(z),npatan2d(ypos,xpos)
1957                    pf = G2lat.polfcal(ODFln,SamSym[textureData['Model']],np.array([r,]),np.array([p,]))
1958                    G2frame.G2plotNB.status.SetFields(['','phi =%9.3f, gam =%9.3f, MRD =%9.3f'%(r,p,pf)])
1959   
1960    try:
1961        plotNum = G2frame.G2plotNB.plotList.index('Texture')
1962        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
1963        Page.figure.clf()
1964        Plot = Page.figure.gca()
1965        if not Page.IsShown():
1966            Page.Show()
1967    except ValueError:
1968        Plot = G2frame.G2plotNB.addMpl('Texture').gca()
1969        plotNum = G2frame.G2plotNB.plotList.index('Texture')
1970        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
1971        Page.canvas.mpl_connect('motion_notify_event', OnMotion)
1972
1973    Page.Choice = None
1974    Page.SetFocus()
1975    G2frame.G2plotNB.status.SetFields(['',''])   
1976    PH = np.array(SHData['PFhkl'])
1977    phi,beta = G2lat.CrsAng(PH,cell,SGData)
1978    ODFln = G2lat.Flnh(Start,SHCoef,phi,beta,SGData)
1979    PX = np.array(SHData['PFxyz'])
1980    gam = atan2d(PX[0],PX[1])
1981    xy = math.sqrt(PX[0]**2+PX[1]**2)
1982    xyz = math.sqrt(PX[0]**2+PX[1]**2+PX[2]**2)
1983    psi = asind(xy/xyz)
1984    IODFln = G2lat.Glnh(Start,SHCoef,psi,gam,SamSym[textureData['Model']])
1985    if 'Axial' in SHData['PlotType']:
1986        X = np.linspace(0,90.0,26)
1987        Y = G2lat.polfcal(ODFln,SamSym[textureData['Model']],X,0.0)
1988        Plot.plot(X,Y,color='k',label=str(SHData['PFhkl']))
1989        Plot.legend(loc='best')
1990        Plot.set_title('Axial distribution for HKL='+str(SHData['PFhkl']))
1991        Plot.set_xlabel(r'$\psi$',fontsize=16)
1992        Plot.set_ylabel('MRD',fontsize=14)
1993       
1994    else:       
1995        npts = 201
1996        if 'Inverse' in SHData['PlotType']:
1997            X,Y = np.meshgrid(np.linspace(1.,-1.,npts),np.linspace(-1.,1.,npts))
1998            R,P = np.sqrt(X**2+Y**2).flatten(),npatan2d(X,Y).flatten()
1999            if 'equal' in G2frame.Projection:
2000                R = np.where(R <= 1.,2.*npasind(R*sq2),0.0)
2001            else:
2002                R = np.where(R <= 1.,2.*npatand(R),0.0)
2003            Z = np.zeros_like(R)
2004            Z = G2lat.invpolfcal(IODFln,SGData,R,P)
2005            Z = np.reshape(Z,(npts,npts))
2006            CS = Plot.contour(Y,X,Z,aspect='equal')
2007            Plot.clabel(CS,fontsize=9,inline=1)
2008            try:
2009                Img = Plot.imshow(Z.T,aspect='equal',cmap=G2frame.ContourColor,extent=[-1,1,-1,1])
2010            except ValueError:
2011                pass
2012            Page.figure.colorbar(Img)
2013            Plot.set_title('Inverse pole figure for XYZ='+str(SHData['PFxyz']))
2014            Plot.set_xlabel(G2frame.Projection.capitalize()+' projection')
2015                       
2016        else:
2017            X,Y = np.meshgrid(np.linspace(1.,-1.,npts),np.linspace(-1.,1.,npts))
2018            R,P = np.sqrt(X**2+Y**2).flatten(),npatan2d(X,Y).flatten()
2019            if 'equal' in G2frame.Projection:
2020                R = np.where(R <= 1.,2.*npasind(R*sq2),0.0)
2021            else:
2022                R = np.where(R <= 1.,2.*npatand(R),0.0)
2023            Z = np.zeros_like(R)
2024            Z = G2lat.polfcal(ODFln,SamSym[textureData['Model']],R,P)
2025            Z = np.reshape(Z,(npts,npts))
2026            try:
2027                CS = Plot.contour(Y,X,Z,aspect='equal')
2028                Plot.clabel(CS,fontsize=9,inline=1)
2029            except ValueError:
2030                pass
2031            Img = Plot.imshow(Z.T,aspect='equal',cmap=G2frame.ContourColor,extent=[-1,1,-1,1])
2032            Page.figure.colorbar(Img)
2033            Plot.set_title('Pole figure for HKL='+str(SHData['PFhkl']))
2034            Plot.set_xlabel(G2frame.Projection.capitalize()+' projection')
2035    Page.canvas.draw()
2036
2037################################################################################
2038##### PlotCovariance
2039################################################################################
2040           
2041def PlotCovariance(G2frame,Data):
2042    'needs a doc string'
2043    if not Data:
2044        print 'No covariance matrix available'
2045        return
2046    varyList = Data['varyList']
2047    values = Data['variables']
2048    Xmax = len(varyList)
2049    covMatrix = Data['covMatrix']
2050    sig = np.sqrt(np.diag(covMatrix))
2051    xvar = np.outer(sig,np.ones_like(sig))
2052    covArray = np.divide(np.divide(covMatrix,xvar),xvar.T)
2053    title = ' for\n'+Data['title']
2054    newAtomDict = Data.get('newAtomDict',{})
2055   
2056
2057    def OnPlotKeyPress(event):
2058        newPlot = False
2059        if event.key == 's':
2060            choice = [m for m in mpl.cm.datad.keys() if not m.endswith("_r")]
2061            choice.sort()
2062            dlg = wx.SingleChoiceDialog(G2frame,'Select','Color scheme',choice)
2063            if dlg.ShowModal() == wx.ID_OK:
2064                sel = dlg.GetSelection()
2065                G2frame.VcovColor = choice[sel]
2066            else:
2067                G2frame.VcovColor = 'RdYlGn'
2068            dlg.Destroy()
2069        PlotCovariance(G2frame,Data)
2070
2071    def OnMotion(event):
2072        #there is a problem here - reports wrong values
2073        if event.button:
2074            ytics = imgAx.get_yticks()
2075            ytics = np.where(ytics<len(varyList),ytics,-1)
2076            ylabs = [np.where(0<=i ,varyList[int(i)],' ') for i in ytics]
2077            imgAx.set_yticklabels(ylabs)           
2078        if event.xdata and event.ydata:                 #avoid out of frame errors
2079            xpos = int(event.xdata+.5)
2080            ypos = int(event.ydata+.5)
2081            if -1 < xpos < len(varyList) and -1 < ypos < len(varyList):
2082                if xpos == ypos:
2083                    value = values[xpos]
2084                    name = varyList[xpos]
2085                    if varyList[xpos] in newAtomDict:
2086                        name,value = newAtomDict[name]                       
2087                    msg = '%s value = %.4g, esd = %.4g'%(name,value,sig[xpos])
2088                else:
2089                    msg = '%s - %s: %5.3f'%(varyList[xpos],varyList[ypos],covArray[xpos][ypos])
2090                Page.canvas.SetToolTipString(msg)
2091                G2frame.G2plotNB.status.SetFields(['',msg])
2092               
2093    try:
2094        plotNum = G2frame.G2plotNB.plotList.index('Covariance')
2095        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
2096        Page.figure.clf()
2097        Plot = Page.figure.gca()
2098        if not Page.IsShown():
2099            Page.Show()
2100    except ValueError:
2101        Plot = G2frame.G2plotNB.addMpl('Covariance').gca()
2102        plotNum = G2frame.G2plotNB.plotList.index('Covariance')
2103        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
2104        Page.canvas.mpl_connect('motion_notify_event', OnMotion)
2105        Page.canvas.mpl_connect('key_press_event', OnPlotKeyPress)
2106    Page.Choice = ['s: to change colors']
2107    Page.keyPress = OnPlotKeyPress
2108    Page.SetFocus()
2109    G2frame.G2plotNB.status.SetFields(['',''])   
2110    acolor = mpl.cm.get_cmap(G2frame.VcovColor)
2111    Img = Plot.imshow(covArray,aspect='equal',cmap=acolor,interpolation='nearest',origin='lower',
2112        vmin=-1.,vmax=1.)
2113    imgAx = Img.get_axes()
2114    ytics = imgAx.get_yticks()
2115    ylabs = [varyList[int(i)] for i in ytics[:-1]]
2116    imgAx.set_yticklabels(ylabs)
2117    colorBar = Page.figure.colorbar(Img)
2118    Plot.set_title('V-Cov matrix'+title)
2119    Plot.set_xlabel('Variable number')
2120    Plot.set_ylabel('Variable name')
2121    Page.canvas.draw()
2122   
2123################################################################################
2124##### PlotTorsion
2125################################################################################
2126
2127def PlotTorsion(G2frame,phaseName,Torsion,TorName,Names=[],Angles=[],Coeff=[]):
2128    'needs a doc string'
2129   
2130    global names
2131    names = Names
2132    sum = np.sum(Torsion)
2133    torsion = np.log(2*Torsion+1.)/sum
2134    tMin = np.min(torsion)
2135    tMax = np.max(torsion)
2136    torsion = 3.*(torsion-tMin)/(tMax-tMin)
2137    X = np.linspace(0.,360.,num=45)
2138   
2139    def OnPick(event):
2140        ind = event.ind[0]
2141        msg = 'atoms:'+names[ind]
2142        Page.canvas.SetToolTipString(msg)
2143        try:
2144            page = G2frame.dataDisplay.GetSelection()
2145        except:
2146            return
2147        if G2frame.dataDisplay.GetPageText(page) == 'Torsion restraints':
2148            torGrid = G2frame.dataDisplay.GetPage(page).Torsions
2149            torGrid.ClearSelection()
2150            for row in range(torGrid.GetNumberRows()):
2151                if names[ind] in torGrid.GetCellValue(row,0):
2152                    torGrid.SelectRow(row)
2153            torGrid.ForceRefresh()
2154               
2155    def OnMotion(event):
2156        if event.xdata and event.ydata:                 #avoid out of frame errors
2157            xpos = event.xdata
2158            ypos = event.ydata
2159            msg = 'torsion,energy: %5.3f %5.3f'%(xpos,ypos)
2160            Page.canvas.SetToolTipString(msg)
2161
2162    try:
2163        plotNum = G2frame.G2plotNB.plotList.index('Torsion')
2164        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
2165        Page.figure.clf()
2166        Plot = Page.figure.gca()
2167        if not Page.IsShown():
2168            Page.Show()
2169    except ValueError:
2170        Plot = G2frame.G2plotNB.addMpl('Torsion').gca()
2171        plotNum = G2frame.G2plotNB.plotList.index('Torsion')
2172        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
2173        Page.canvas.mpl_connect('pick_event', OnPick)
2174        Page.canvas.mpl_connect('motion_notify_event', OnMotion)
2175   
2176    Page.SetFocus()
2177    G2frame.G2plotNB.status.SetFields(['','Use mouse LB to identify torsion atoms'])
2178    Plot.plot(X,torsion,'b+')
2179    if len(Coeff):
2180        X2 = np.linspace(0.,360.,45)
2181        Y2 = np.array([-G2mth.calcTorsionEnergy(x,Coeff)[1] for x in X2])
2182        Plot.plot(X2,Y2,'r')
2183    if len(Angles):
2184        Eval = np.array([-G2mth.calcTorsionEnergy(x,Coeff)[1] for x in Angles])
2185        Plot.plot(Angles,Eval,'ro',picker=5)
2186    Plot.set_xlim((0.,360.))
2187    Plot.set_title('Torsion angles for '+TorName+' in '+phaseName)
2188    Plot.set_xlabel('angle',fontsize=16)
2189    Plot.set_ylabel('Energy',fontsize=16)
2190    Page.canvas.draw()
2191   
2192################################################################################
2193##### PlotRama
2194################################################################################
2195
2196def PlotRama(G2frame,phaseName,Rama,RamaName,Names=[],PhiPsi=[],Coeff=[]):
2197    'needs a doc string'
2198
2199    global names
2200    names = Names
2201    rama = np.log(2*Rama+1.)
2202    ramaMax = np.max(rama)
2203    rama = np.reshape(rama,(45,45))
2204    global Phi,Psi
2205    Phi = []
2206    Psi = []
2207
2208    def OnPlotKeyPress(event):
2209        newPlot = False
2210        if event.key == 's':
2211            choice = [m for m in mpl.cm.datad.keys() if not m.endswith("_r")]
2212            choice.sort()
2213            dlg = wx.SingleChoiceDialog(G2frame,'Select','Color scheme',choice)
2214            if dlg.ShowModal() == wx.ID_OK:
2215                sel = dlg.GetSelection()
2216                G2frame.RamaColor = choice[sel]
2217            else:
2218                G2frame.RamaColor = 'RdYlGn'
2219            dlg.Destroy()
2220        PlotRama(G2frame,phaseName,Rama)
2221
2222    def OnPick(event):
2223        ind = event.ind[0]
2224        msg = 'atoms:'+names[ind]
2225        Page.canvas.SetToolTipString(msg)
2226        try:
2227            page = G2frame.dataDisplay.GetSelection()
2228        except:
2229            return
2230        if G2frame.dataDisplay.GetPageText(page) == 'Ramachandran restraints':
2231            ramaGrid = G2frame.dataDisplay.GetPage(page).Ramas
2232            ramaGrid.ClearSelection()
2233            for row in range(ramaGrid.GetNumberRows()):
2234                if names[ind] in ramaGrid.GetCellValue(row,0):
2235                    ramaGrid.SelectRow(row)
2236            ramaGrid.ForceRefresh()
2237
2238    def OnMotion(event):
2239        if event.xdata and event.ydata:                 #avoid out of frame errors
2240            xpos = event.xdata
2241            ypos = event.ydata
2242            msg = 'phi/psi: %5.3f %5.3f'%(xpos,ypos)
2243            Page.canvas.SetToolTipString(msg)
2244           
2245    try:
2246        plotNum = G2frame.G2plotNB.plotList.index('Ramachandran')
2247        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
2248        Page.figure.clf()
2249        Plot = Page.figure.gca()
2250        if not Page.IsShown():
2251            Page.Show()
2252    except ValueError:
2253        Plot = G2frame.G2plotNB.addMpl('Ramachandran').gca()
2254        plotNum = G2frame.G2plotNB.plotList.index('Ramachandran')
2255        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
2256        Page.canvas.mpl_connect('pick_event', OnPick)
2257        Page.canvas.mpl_connect('motion_notify_event', OnMotion)
2258        Page.canvas.mpl_connect('key_press_event', OnPlotKeyPress)
2259
2260    Page.Choice = ['s: to change colors']
2261    Page.keyPress = OnPlotKeyPress
2262    Page.SetFocus()
2263    G2frame.G2plotNB.status.SetFields(['','Use mouse LB to identify phi/psi atoms'])
2264    acolor = mpl.cm.get_cmap(G2frame.RamaColor)
2265    if RamaName == 'All' or '-1' in RamaName:
2266        if len(Coeff): 
2267            X,Y = np.meshgrid(np.linspace(-180.,180.,45),np.linspace(-180.,180.,45))
2268            Z = np.array([-G2mth.calcRamaEnergy(x,y,Coeff)[1] for x,y in zip(X.flatten(),Y.flatten())])
2269            Plot.contour(X,Y,np.reshape(Z,(45,45)))
2270        Img = Plot.imshow(rama,aspect='equal',cmap=acolor,interpolation='nearest',
2271            extent=[-180,180,-180,180],origin='lower')
2272        if len(PhiPsi):
2273            Phi,Psi = PhiPsi.T
2274            Phi = np.where(Phi>180.,Phi-360.,Phi)
2275            Psi = np.where(Psi>180.,Psi-360.,Psi)
2276            Plot.plot(Phi,Psi,'ro',picker=5)
2277        Plot.set_xlim((-180.,180.))
2278        Plot.set_ylim((-180.,180.))
2279    else:
2280        if len(Coeff): 
2281            X,Y = np.meshgrid(np.linspace(0.,360.,45),np.linspace(0.,360.,45))
2282            Z = np.array([-G2mth.calcRamaEnergy(x,y,Coeff)[1] for x,y in zip(X.flatten(),Y.flatten())])
2283            Plot.contour(X,Y,np.reshape(Z,(45,45)))
2284        Img = Plot.imshow(rama,aspect='equal',cmap=acolor,interpolation='nearest',
2285            extent=[0,360,0,360],origin='lower')
2286        if len(PhiPsi):
2287            Phi,Psi = PhiPsi.T
2288            Plot.plot(Phi,Psi,'ro',picker=5)
2289        Plot.set_xlim((0.,360.))
2290        Plot.set_ylim((0.,360.))
2291    Plot.set_title('Ramachandran for '+RamaName+' in '+phaseName)
2292    Plot.set_xlabel(r'$\phi$',fontsize=16)
2293    Plot.set_ylabel(r'$\psi$',fontsize=16)
2294    colorBar = Page.figure.colorbar(Img)
2295    Page.canvas.draw()
2296
2297
2298################################################################################
2299##### PlotSeq
2300################################################################################
2301def PlotSelectedSequence(G2frame,ColumnList,TableGet,SelectX,fitnum=None,fitvals=None):
2302    '''Plot a result from a sequential refinement
2303
2304    :param wx.Frame G2frame: The main GSAS-II tree "window"
2305    :param list ColumnList: list of int values corresponding to columns
2306      selected as y values
2307    :param function TableGet: a function that takes a column number
2308      as argument and returns the column label, the values and there ESDs (or None)
2309    :param function SelectX: a function that returns a selected column
2310      number (or None) as the X-axis selection
2311    '''
2312    G2frame.seqXselect = SelectX
2313    try:
2314        G2frame.seqXaxis
2315    except:
2316        G2frame.seqXaxis = None
2317
2318    def OnMotion(event):
2319        if event.xdata and event.ydata:                 #avoid out of frame errors
2320            xpos = event.xdata
2321            ypos = event.ydata
2322            msg = '%5.3f %.6g'%(xpos,ypos)
2323            Page.canvas.SetToolTipString(msg)
2324
2325    def OnKeyPress(event):
2326        if event.key == 's':
2327            G2frame.seqXaxis = G2frame.seqXselect()
2328            Draw()
2329           
2330    if fitnum is None:
2331        label = 'Sequential refinement'
2332    else:
2333        label = 'Parametric fit #'+str(fitnum+1)
2334    try:
2335        plotNum = G2frame.G2plotNB.plotList.index(label)
2336        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
2337        Page.figure.clf()
2338        Plot = Page.figure.gca()
2339        if not Page.IsShown():
2340            Page.Show()
2341    except ValueError:
2342        Plot = G2frame.G2plotNB.addMpl(label).gca()
2343        plotNum = G2frame.G2plotNB.plotList.index(label)
2344        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
2345        Page.canvas.mpl_connect('key_press_event', OnKeyPress)
2346        Page.canvas.mpl_connect('motion_notify_event', OnMotion)
2347    Page.Choice = ['s to select plot x-axis',]
2348    Page.keyPress = OnKeyPress
2349    Page.seqYaxisList = ColumnList
2350    Page.seqTableGet = TableGet
2351    Page.fitvals = fitvals
2352       
2353    def Draw():
2354        Page.SetFocus()
2355        G2frame.G2plotNB.status.SetStatusText('press s to select X axis',1)
2356        Plot.clear()
2357        if G2frame.seqXaxis is not None:   
2358            xName,X,Xsig = Page.seqTableGet(G2frame.seqXaxis)
2359        else:
2360            X = np.arange(0,G2frame.SeqTable.GetNumberRows(),1)
2361            xName = 'Data sequence number'
2362        for col in Page.seqYaxisList:
2363            name,Y,sig = Page.seqTableGet(col)
2364            if sig:
2365                if G2frame.seqReverse and not G2frame.seqXaxis:
2366                    Y = Y[::-1]
2367                    sig = sig[::-1]
2368                Plot.errorbar(X,Y,yerr=sig,label=name)
2369            else:
2370                if G2frame.seqReverse and not G2frame.seqXaxis:
2371                    Y = Y[::-1]
2372                Plot.plot(X,Y)
2373                Plot.plot(X,Y,'o',label=name)
2374        if Page.fitvals:
2375            if G2frame.seqReverse and not G2frame.seqXaxis:
2376                Page.fitvals = Page.fitvals[::-1]
2377            Plot.plot(X,Page.fitvals,label='Fit')
2378           
2379        Plot.legend(loc='best')
2380        Plot.set_ylabel('Parameter values')
2381        Plot.set_xlabel(xName)
2382        Page.canvas.draw()           
2383    Draw()
2384               
2385################################################################################
2386##### PlotExposedImage & PlotImage
2387################################################################################
2388           
2389def PlotExposedImage(G2frame,newPlot=False,event=None):
2390    '''General access module for 2D image plotting
2391    '''
2392    plotNo = G2frame.G2plotNB.nb.GetSelection()
2393    if G2frame.G2plotNB.nb.GetPageText(plotNo) == '2D Powder Image':
2394        PlotImage(G2frame,newPlot,event,newImage=True)
2395    elif G2frame.G2plotNB.nb.GetPageText(plotNo) == '2D Integration':
2396        PlotIntegration(G2frame,newPlot,event)
2397
2398def OnStartMask(G2frame):
2399    '''Initiate the start of a Frame or Polygon map
2400
2401    :param wx.Frame G2frame: The main GSAS-II tree "window"
2402    :param str eventkey: a single letter ('f' or 'p') that
2403      determines what type of mask is created.   
2404    '''
2405    Masks = G2frame.PatternTree.GetItemPyData(
2406        G2gd.GetPatternTreeItemId(G2frame,G2frame.Image, 'Masks'))
2407    if G2frame.MaskKey == 'f':
2408        Masks['Frames'] = []
2409    elif G2frame.MaskKey == 'p':
2410        Masks['Polygons'].append([])
2411    elif G2frame.MaskKey == 's':
2412        Masks['Points'].append([])
2413    elif G2frame.MaskKey == 'a':
2414        Masks['Arcs'].append([])
2415    elif G2frame.MaskKey == 'r':
2416        Masks['Rings'].append([])
2417    G2imG.UpdateMasks(G2frame,Masks)
2418    PlotImage(G2frame,newImage=True)
2419   
2420def OnStartNewDzero(G2frame):
2421    '''Initiate the start of adding a new d-zero to a strain data set
2422
2423    :param wx.Frame G2frame: The main GSAS-II tree "window"
2424    :param str eventkey: a single letter ('a') that
2425      triggers the addition of a d-zero.   
2426    '''
2427    G2frame.dataFrame.GetStatusBar().SetStatusText('Add strain ring active - LB pick d-zero value',1)
2428    StrSta = G2frame.PatternTree.GetItemPyData(
2429        G2gd.GetPatternTreeItemId(G2frame,G2frame.Image, 'Stress/Strain'))
2430
2431def PlotImage(G2frame,newPlot=False,event=None,newImage=True):
2432    '''Plot of 2D detector images as contoured plot. Also plot calibration ellipses,
2433    masks, etc.
2434    '''
2435    from matplotlib.patches import Ellipse,Arc,Circle,Polygon
2436    import numpy.ma as ma
2437    Dsp = lambda tth,wave: wave/(2.*npsind(tth/2.))
2438    global Data,Masks
2439    colors=['b','g','r','c','m','k']
2440    Data = G2frame.PatternTree.GetItemPyData(
2441        G2gd.GetPatternTreeItemId(G2frame,G2frame.Image, 'Image Controls'))
2442# patch
2443    if 'invert_x' not in Data:
2444        Data['invert_x'] = False
2445        Data['invert_y'] = True
2446# end patch
2447    Masks = G2frame.PatternTree.GetItemPyData(
2448        G2gd.GetPatternTreeItemId(G2frame,G2frame.Image, 'Masks'))
2449    try:    #may be absent
2450        StrSta = G2frame.PatternTree.GetItemPyData(
2451            G2gd.GetPatternTreeItemId(G2frame,G2frame.Image, 'Stress/Strain'))
2452    except TypeError:   #is missing
2453        StrSta = {}
2454
2455    def OnImMotion(event):
2456        Page.canvas.SetToolTipString('')
2457        sizexy = Data['size']
2458        if event.xdata and event.ydata and len(G2frame.ImageZ):                 #avoid out of frame errors
2459            Page.canvas.SetToolTipString('%8.2f %8.2fmm'%(event.xdata,event.ydata))
2460            Page.canvas.SetCursor(wx.CROSS_CURSOR)
2461            item = G2frame.itemPicked
2462            pixelSize = Data['pixelSize']
2463            scalex = 1000./pixelSize[0]
2464            scaley = 1000./pixelSize[1]
2465            if item and G2frame.PatternTree.GetItemText(G2frame.PickId) == 'Image Controls':
2466                if 'Text' in str(item):
2467                    Page.canvas.SetToolTipString('%8.3f %8.3fmm'%(event.xdata,event.ydata))
2468                else:
2469                    xcent,ycent = Data['center']
2470                    xpos = event.xdata-xcent
2471                    ypos = event.ydata-ycent
2472                    tth,azm = G2img.GetTthAzm(event.xdata,event.ydata,Data)
2473                    if 'line3' in  str(item) or 'line4' in str(item) and not Data['fullIntegrate']:
2474                        Page.canvas.SetToolTipString('%6d deg'%(azm))
2475                    elif 'line1' in  str(item) or 'line2' in str(item):
2476                        Page.canvas.SetToolTipString('%8.3f deg'%(tth))                           
2477            else:
2478                xpos = event.xdata
2479                ypos = event.ydata
2480                xpix = xpos*scalex
2481                ypix = ypos*scaley
2482                Int = 0
2483                if (0 <= xpix <= sizexy[0]) and (0 <= ypix <= sizexy[1]):
2484                    Int = G2frame.ImageZ[ypix][xpix]
2485                tth,azm,D,dsp = G2img.GetTthAzmDsp(xpos,ypos,Data)
2486                Q = 2.*math.pi/dsp
2487                fields = ['','Detector 2-th =%9.3fdeg, dsp =%9.3fA, Q = %6.5fA-1, azm = %7.2fdeg, I = %6d'%(tth,dsp,Q,azm,Int)]
2488                if G2frame.MaskKey in ['p','f']:
2489                    fields[1] = 'Polygon/frame mask pick - LB next point, RB close polygon'
2490                elif G2frame.StrainKey:
2491                    fields[0] = 'd-zero pick active'
2492                G2frame.G2plotNB.status.SetFields(fields)
2493
2494    def OnImPlotKeyPress(event):
2495        try:
2496            PickName = G2frame.PatternTree.GetItemText(G2frame.PickId)
2497        except TypeError:
2498            return
2499        if PickName == 'Masks':
2500            if event.key in ['l','p','f','s','a','r']:
2501                G2frame.MaskKey = event.key
2502                OnStartMask(G2frame)
2503                PlotImage(G2frame,newPlot=False)
2504               
2505        elif PickName == 'Stress/Strain':
2506            if event.key in ['a',]:
2507                G2frame.StrainKey = event.key
2508                OnStartNewDzero(G2frame)
2509                PlotImage(G2frame,newPlot=False)
2510               
2511        elif PickName == 'Image Controls':
2512            if event.key in ['c',]:
2513                Xpos = event.xdata
2514                if not Xpos:            #got point out of frame
2515                    return
2516                Ypos = event.ydata
2517                dlg = wx.MessageDialog(G2frame,'Are you sure you want to change the center?',
2518                    'Center change',style=wx.OK|wx.CANCEL)
2519                try:
2520                    if dlg.ShowModal() == wx.ID_OK:
2521                        print 'move center to: ',Xpos,Ypos
2522                        Data['center'] = [Xpos,Ypos]
2523                        G2imG.UpdateImageControls(G2frame,Data,Masks)
2524                        PlotImage(G2frame,newPlot=False)
2525                finally:
2526                    dlg.Destroy()
2527                return
2528            elif event.key == 'l':
2529                G2frame.logPlot = not G2frame.logPlot
2530            elif event.key in ['x',]:
2531                Data['invert_x'] = not Data['invert_x']
2532            elif event.key in ['y',]:
2533                Data['invert_y'] = not Data['invert_y']
2534            PlotImage(G2frame,newPlot=True)
2535           
2536    def OnKeyBox(event):
2537        if G2frame.G2plotNB.nb.GetSelection() == G2frame.G2plotNB.plotList.index('2D Powder Image'):
2538            event.key = cb.GetValue()[0]
2539            cb.SetValue(' key press')
2540            if event.key in ['l','s','a','r','p','x','y']:
2541                wx.CallAfter(OnImPlotKeyPress,event)
2542        Page.canvas.SetFocus() # redirect the Focus from the button back to the plot
2543                       
2544    def OnImPick(event):
2545        if G2frame.PatternTree.GetItemText(G2frame.PickId) not in ['Image Controls','Masks']:
2546            return
2547        if G2frame.itemPicked is not None: return
2548        G2frame.itemPicked = event.artist
2549        G2frame.mousePicked = event.mouseevent
2550       
2551    def OnImRelease(event):
2552        try:
2553            PickName = G2frame.PatternTree.GetItemText(G2frame.PickId)
2554        except TypeError:
2555            return
2556        if PickName not in ['Image Controls','Masks','Stress/Strain']:
2557            return
2558        pixelSize = Data['pixelSize']
2559        scalex = 1000./pixelSize[0]
2560        scaley = 1000./pixelSize[1]
2561        pixLimit = Data['pixLimit']
2562        if G2frame.itemPicked is None and PickName == 'Image Controls' and len(G2frame.ImageZ):
2563            Xpos = event.xdata
2564            if not (Xpos and G2frame.ifGetRing):                   #got point out of frame
2565                return
2566            Ypos = event.ydata
2567            if Ypos and not Page.toolbar._active:         #make sure zoom/pan not selected
2568                if event.button == 1:
2569                    Xpix = Xpos*scalex
2570                    Ypix = Ypos*scaley
2571                    xpos,ypos,I,J = G2img.ImageLocalMax(G2frame.ImageZ,pixLimit,Xpix,Ypix)
2572                    if I and J:
2573                        xpos += .5                              #shift to pixel center
2574                        ypos += .5
2575                        xpos /= scalex                          #convert to mm
2576                        ypos /= scaley
2577                        Data['ring'].append([xpos,ypos])
2578                elif event.button == 3:
2579                    G2frame.dataFrame.GetStatusBar().SetStatusText('Calibrating...',0)
2580                    if G2img.ImageCalibrate(G2frame,Data):
2581                        G2frame.dataFrame.GetStatusBar().SetStatusText('Calibration successful - Show ring picks to check',0)
2582                        print 'Calibration successful'
2583                    else:
2584                        G2frame.dataFrame.GetStatusBar().SetStatusText('Calibration failed - Show ring picks to diagnose',0)
2585                        print 'Calibration failed'
2586                    G2frame.ifGetRing = False
2587                    G2imG.UpdateImageControls(G2frame,Data,Masks)
2588                    return
2589                PlotImage(G2frame,newImage=False)
2590            return
2591        elif G2frame.MaskKey and PickName == 'Masks':
2592            Xpos,Ypos = [event.xdata,event.ydata]
2593            if not Xpos or not Ypos or Page.toolbar._active:  #got point out of frame or zoom/pan selected
2594                return
2595            if G2frame.MaskKey == 's' and event.button == 1:
2596                Masks['Points'][-1] = [Xpos,Ypos,1.]
2597                G2frame.MaskKey = ''               
2598            elif G2frame.MaskKey == 'r' and event.button == 1:
2599                tth = G2img.GetTth(Xpos,Ypos,Data)
2600                Masks['Rings'][-1] = [tth,0.1]
2601                G2frame.MaskKey = ''               
2602            elif G2frame.MaskKey == 'a' and event.button == 1:
2603                tth,azm = G2img.GetTthAzm(Xpos,Ypos,Data)
2604                azm = int(azm)               
2605                Masks['Arcs'][-1] = [tth,[azm-5,azm+5],0.1]
2606                G2frame.MaskKey = ''               
2607            elif G2frame.MaskKey =='p':
2608                polygon = Masks['Polygons'][-1]
2609                if len(polygon) > 2 and event.button == 3:
2610                    x0,y0 = polygon[0]
2611                    polygon.append([x0,y0])
2612                    G2frame.MaskKey = ''
2613                    G2frame.G2plotNB.status.SetFields(['','Polygon closed'])
2614                else:
2615                    G2frame.G2plotNB.status.SetFields(['','New polygon point: %.1f,%.1f'%(Xpos,Ypos)])
2616                    polygon.append([Xpos,Ypos])
2617            elif G2frame.MaskKey =='f':
2618                frame = Masks['Frames']
2619                if len(frame) > 2 and event.button == 3:
2620                    x0,y0 = frame[0]
2621                    frame.append([x0,y0])
2622                    G2frame.MaskKey = ''
2623                    G2frame.G2plotNB.status.SetFields(['','Frame closed'])
2624                else:
2625                    G2frame.G2plotNB.status.SetFields(['','New frame point: %.1f,%.1f'%(Xpos,Ypos)])
2626                    frame.append([Xpos,Ypos])
2627            G2imG.UpdateMasks(G2frame,Masks)
2628            PlotImage(G2frame,newImage=False)
2629        elif PickName == 'Stress/Strain' and G2frame.StrainKey:
2630            Xpos,Ypos = [event.xdata,event.ydata]
2631            if not Xpos or not Ypos or Page.toolbar._active:  #got point out of frame or zoom/pan selected
2632                return
2633            dsp = G2img.GetDsp(Xpos,Ypos,Data)
2634            StrSta['d-zero'].append({'Dset':dsp,'Dcalc':0.0,'pixLimit':10,'cutoff':1.0,
2635                'ImxyObs':[[],[]],'ImtaObs':[[],[]],'ImtaCalc':[[],[]],'Emat':[1.0,1.0,1.0]})
2636            R,r = G2img.MakeStrStaRing(StrSta['d-zero'][-1],G2frame.ImageZ,Data)
2637            if not len(R):
2638                del StrSta['d-zero'][-1]
2639                G2frame.ErrorDialog('Strain peak selection','WARNING - No points found for this ring selection')
2640            StrSta['d-zero'] = G2mth.sortArray(StrSta['d-zero'],'Dset',reverse=True)
2641            G2frame.StrainKey = ''
2642            G2imG.UpdateStressStrain(G2frame,StrSta)
2643            PlotStrain(G2frame,StrSta)
2644            PlotImage(G2frame,newPlot=False)           
2645        else:
2646            Xpos,Ypos = [event.xdata,event.ydata]
2647            if not Xpos or not Ypos or Page.toolbar._active:  #got point out of frame or zoom/pan selected
2648                return
2649            if G2frame.ifGetRing:                          #delete a calibration ring pick
2650                xypos = [Xpos,Ypos]
2651                rings = Data['ring']
2652                for ring in rings:
2653                    if np.allclose(ring,xypos,.01,0):
2654                        rings.remove(ring)
2655            else:
2656                tth,azm,dsp = G2img.GetTthAzmDsp(Xpos,Ypos,Data)[:3]
2657                itemPicked = str(G2frame.itemPicked)
2658                if 'Line2D' in itemPicked and PickName == 'Image Controls':
2659                    if 'line1' in itemPicked:
2660                        Data['IOtth'][0] = max(tth,0.001)
2661                    elif 'line2' in itemPicked:
2662                        Data['IOtth'][1] = tth
2663                    elif 'line3' in itemPicked:
2664                        Data['LRazimuth'][0] = int(azm)
2665                    elif 'line4' in itemPicked and not Data['fullIntegrate']:
2666                        Data['LRazimuth'][1] = int(azm)
2667                   
2668                    Data['LRazimuth'][0] %= 360
2669                    Data['LRazimuth'][1] %= 360
2670                    if Data['LRazimuth'][0] > Data['LRazimuth'][1]:
2671                        Data['LRazimuth'][1] += 360                       
2672                    if Data['fullIntegrate']:
2673                        Data['LRazimuth'][1] = Data['LRazimuth'][0]+360
2674                       
2675                    if  Data['IOtth'][0] > Data['IOtth'][1]:
2676                        Data['IOtth'][0],Data['IOtth'][1] = Data['IOtth'][1],Data['IOtth'][0]
2677                       
2678                    G2frame.InnerTth.SetValue("%8.2f" % (Data['IOtth'][0]))
2679                    G2frame.OuterTth.SetValue("%8.2f" % (Data['IOtth'][1]))
2680                    G2frame.Lazim.SetValue("%6d" % (Data['LRazimuth'][0]))
2681                    G2frame.Razim.SetValue("%6d" % (Data['LRazimuth'][1]))
2682                elif 'Circle' in itemPicked and PickName == 'Masks':
2683                    spots = Masks['Points']
2684                    newPos = itemPicked.split(')')[0].split('(')[2].split(',')
2685                    newPos = np.array([float(newPos[0]),float(newPos[1])])
2686                    for spot in spots:
2687                        if np.allclose(np.array([spot[:2]]),newPos):
2688                            spot[:2] = Xpos,Ypos
2689                    G2imG.UpdateMasks(G2frame,Masks)
2690                elif 'Line2D' in itemPicked and PickName == 'Masks':
2691                    Obj = G2frame.itemPicked.findobj()
2692                    rings = Masks['Rings']
2693                    arcs = Masks['Arcs']
2694                    polygons = Masks['Polygons']
2695                    frame = Masks['Frames']
2696                    for ring in G2frame.ringList:
2697                        if Obj == ring[0]:
2698                            rN = ring[1]
2699                            if ring[2] == 'o':
2700                                rings[rN][0] = G2img.GetTth(Xpos,Ypos,Data)-rings[rN][1]/2.
2701                            else:
2702                                rings[rN][0] = G2img.GetTth(Xpos,Ypos,Data)+rings[rN][1]/2.
2703                    for arc in G2frame.arcList:
2704                        if Obj == arc[0]:
2705                            aN = arc[1]
2706                            if arc[2] == 'o':
2707                                arcs[aN][0] = G2img.GetTth(Xpos,Ypos,Data)-arcs[aN][2]/2
2708                            elif arc[2] == 'i':
2709                                arcs[aN][0] = G2img.GetTth(Xpos,Ypos,Data)+arcs[aN][2]/2
2710                            elif arc[2] == 'l':
2711                                arcs[aN][1][0] = int(G2img.GetAzm(Xpos,Ypos,Data))
2712                            else:
2713                                arcs[aN][1][1] = int(G2img.GetAzm(Xpos,Ypos,Data))
2714                    for poly in G2frame.polyList:   #merging points problem here?
2715                        if Obj == poly[0]:
2716                            ind = G2frame.itemPicked.contains(G2frame.mousePicked)[1]['ind'][0]
2717                            oldPos = np.array([G2frame.mousePicked.xdata,G2frame.mousePicked.ydata])
2718                            pN = poly[1]
2719                            for i,xy in enumerate(polygons[pN]):
2720                                if np.allclose(np.array([xy]),oldPos,atol=1.0):
2721                                    if event.button == 1:
2722                                        polygons[pN][i] = Xpos,Ypos
2723                                    elif event.button == 3:
2724                                        polygons[pN].insert(i,[Xpos,Ypos])
2725                                        break
2726                    if frame:
2727                        oldPos = np.array([G2frame.mousePicked.xdata,G2frame.mousePicked.ydata])
2728                        for i,xy in enumerate(frame):
2729                            if np.allclose(np.array([xy]),oldPos,atol=1.0):
2730                                if event.button == 1:
2731                                    frame[i] = Xpos,Ypos
2732                                elif event.button == 3:
2733                                    frame.insert(i,[Xpos,Ypos])
2734                                    break
2735                    G2imG.UpdateMasks(G2frame,Masks)
2736#                else:                  #keep for future debugging
2737#                    print str(G2frame.itemPicked),event.xdata,event.ydata,event.button
2738            PlotImage(G2frame,newImage=True)
2739            G2frame.itemPicked = None
2740           
2741    try:
2742        plotNum = G2frame.G2plotNB.plotList.index('2D Powder Image')
2743        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
2744        if not newPlot:
2745            Plot = Page.figure.gca()          #get previous powder plot & get limits
2746            xylim = Plot.get_xlim(),Plot.get_ylim()
2747        if newImage:
2748            Page.figure.clf()
2749            Plot = Page.figure.gca()          #get a fresh plot after clf()
2750    except ValueError:
2751        Plot = G2frame.G2plotNB.addMpl('2D Powder Image').gca()
2752        plotNum = G2frame.G2plotNB.plotList.index('2D Powder Image')
2753        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
2754        Page.canvas.mpl_connect('key_press_event', OnImPlotKeyPress)
2755        Page.canvas.mpl_connect('motion_notify_event', OnImMotion)
2756        Page.canvas.mpl_connect('pick_event', OnImPick)
2757        Page.canvas.mpl_connect('button_release_event', OnImRelease)
2758        xylim = []
2759    Page.Choice = None
2760    if not event:                       #event from GUI TextCtrl - don't want focus to change to plot!!!
2761        Page.SetFocus()
2762    Title = G2frame.PatternTree.GetItemText(G2frame.Image)[4:]
2763    G2frame.G2plotNB.status.DestroyChildren()
2764    if G2frame.logPlot:
2765        Title = 'log('+Title+')'
2766    Plot.set_title(Title)
2767    try:
2768        if G2frame.PatternTree.GetItemText(G2frame.PickId) in ['Image Controls',]:
2769            Page.Choice = (' key press','l: log(I) on','x: flip x','y: flip y',)
2770            if G2frame.logPlot:
2771                Page.Choice[1] = 'l: log(I) off'
2772            Page.keyPress = OnImPlotKeyPress
2773        elif G2frame.PatternTree.GetItemText(G2frame.PickId) in ['Masks',]:
2774            Page.Choice = (' key press','l: log(I) on','s: spot mask','a: arc mask','r: ring mask',
2775                'p: polygon mask','f: frame mask',)
2776            if G2frame.logPlot:
2777                Page.Choice[1] = 'l: log(I) off'
2778            Page.keyPress = OnImPlotKeyPress
2779        elif G2frame.PatternTree.GetItemText(G2frame.PickId) in ['Stress/Strain',]:
2780            Page.Choice = (' key press','a: add new ring',)
2781            Page.keyPress = OnImPlotKeyPress
2782    except TypeError:
2783        pass
2784    size,imagefile = G2frame.PatternTree.GetItemPyData(G2frame.Image)
2785    if imagefile != G2frame.oldImagefile:
2786        imagefile = G2IO.CheckImageFile(G2frame,imagefile)
2787        if not imagefile:
2788            G2frame.G2plotNB.Delete('2D Powder Image')
2789            return
2790        G2frame.PatternTree.SetItemPyData(G2frame.Image,[size,imagefile])
2791        G2frame.ImageZ = G2IO.GetImageData(G2frame,imagefile,imageOnly=True)
2792        G2frame.oldImagefile = imagefile
2793
2794    imScale = 1
2795    if len(G2frame.ImageZ) > 1024:
2796        imScale = len(G2frame.ImageZ)/1024
2797    sizexy = Data['size']
2798    pixelSize = Data['pixelSize']
2799    scalex = 1000./pixelSize[0]
2800    scaley = 1000./pixelSize[1]
2801    Xmax = sizexy[0]*pixelSize[0]/1000.
2802    Ymax = sizexy[1]*pixelSize[1]/1000.
2803    xlim = (0,Xmax)
2804    ylim = (Ymax,0)
2805    Imin,Imax = Data['range'][1]
2806    acolor = mpl.cm.get_cmap(Data['color'])
2807    xcent,ycent = Data['center']
2808    Plot.set_xlabel('Image x-axis, mm',fontsize=12)
2809    Plot.set_ylabel('Image y-axis, mm',fontsize=12)
2810    #do threshold mask - "real" mask - others are just bondaries
2811    Zlim = Masks['Thresholds'][1]
2812    wx.BeginBusyCursor()
2813    try:
2814           
2815        if newImage:                   
2816            Imin,Imax = Data['range'][1]
2817            MA = ma.masked_greater(ma.masked_less(G2frame.ImageZ,Zlim[0]),Zlim[1])
2818            MaskA = ma.getmaskarray(MA)
2819            A = G2img.ImageCompress(MA,imScale)
2820            AM = G2img.ImageCompress(MaskA,imScale)
2821            if G2frame.logPlot:
2822                A = np.where(A>Imin,np.where(A<Imax,A,0),0)
2823                A = np.where(A>0,np.log(A),0)
2824                AM = np.where(AM>0,np.log(AM),0)
2825                Imin,Imax = [np.amin(A),np.amax(A)]
2826            ImgM = Plot.imshow(AM,aspect='equal',cmap='Reds',
2827                interpolation='nearest',vmin=0,vmax=2,extent=[0,Xmax,Ymax,0])
2828            Img = Plot.imshow(A,aspect='equal',cmap=acolor,
2829                interpolation='nearest',vmin=Imin,vmax=Imax,extent=[0,Xmax,Ymax,0])
2830   
2831        Plot.plot(xcent,ycent,'x')
2832        #G2frame.PatternTree.GetItemText(item)
2833        if Data['showLines']:
2834            LRAzim = Data['LRazimuth']                  #NB: integers
2835            Nazm = Data['outAzimuths']
2836            delAzm = float(LRAzim[1]-LRAzim[0])/Nazm
2837            AzmthOff = Data['azmthOff']
2838            IOtth = Data['IOtth']
2839            wave = Data['wavelength']
2840            dspI = wave/(2.0*sind(IOtth[0]/2.0))
2841            ellI = G2img.GetEllipse(dspI,Data)           #=False if dsp didn't yield an ellipse (ugh! a parabola or a hyperbola)
2842            dspO = wave/(2.0*sind(IOtth[1]/2.0))
2843            ellO = G2img.GetEllipse(dspO,Data)           #Ditto & more likely for outer ellipse
2844            Azm = np.array(range(LRAzim[0],LRAzim[1]+1))-AzmthOff
2845            if ellI:
2846                xyI = []
2847                for azm in Azm:
2848                    xy = G2img.GetDetectorXY(dspI,azm,Data)
2849                    if np.any(xy):
2850                        xyI.append(xy)
2851                if len(xyI):
2852                    xyI = np.array(xyI)
2853                    arcxI,arcyI = xyI.T
2854                    Plot.plot(arcxI,arcyI,picker=3)
2855            if ellO:
2856                xyO = []
2857                for azm in Azm:
2858                    xy = G2img.GetDetectorXY(dspO,azm,Data)
2859                    if np.any(xy):
2860                        xyO.append(xy)
2861                if len(xyO):
2862                    xyO = np.array(xyO)
2863                    arcxO,arcyO = xyO.T               
2864                    Plot.plot(arcxO,arcyO,picker=3)
2865            if ellO and ellI:
2866                Plot.plot([arcxI[0],arcxO[0]],[arcyI[0],arcyO[0]],picker=3)
2867                Plot.plot([arcxI[-1],arcxO[-1]],[arcyI[-1],arcyO[-1]],picker=3)
2868            for i in range(Nazm):
2869                cake = LRAzim[0]+i*delAzm-AzmthOff
2870                if Data['centerAzm']:
2871                    cake += delAzm/2.
2872                ind = np.searchsorted(Azm,cake)
2873                Plot.plot([arcxI[ind],arcxO[ind]],[arcyI[ind],arcyO[ind]],color='k',dashes=(5,5))
2874                   
2875        if G2frame.PatternTree.GetItemText(G2frame.PickId) in 'Image Controls':
2876            for xring,yring in Data['ring']:
2877                Plot.plot(xring,yring,'r+',picker=3)
2878            if Data['setRings']:
2879                N = 0
2880                for ring in Data['rings']:
2881                    xring,yring = np.array(ring).T[:2]
2882                    Plot.plot(xring,yring,'.',color=colors[N%6])
2883                    N += 1           
2884            for ellipse in Data['ellipses']:      #what about hyperbola?
2885                cent,phi,[width,height],col = ellipse
2886                if width > 0:       #ellipses
2887                    Plot.add_artist(Ellipse([cent[0],cent[1]],2*width,2*height,phi,ec=col,fc='none'))
2888                    Plot.text(cent[0],cent[1],'+',color=col,ha='center',va='center')
2889        if G2frame.PatternTree.GetItemText(G2frame.PickId) in 'Stress/Strain':
2890            for N,ring in enumerate(StrSta['d-zero']):
2891                xring,yring = ring['ImxyObs']
2892                Plot.plot(xring,yring,colors[N%6]+'.')
2893        #masks - mask lines numbered after integration limit lines
2894        spots = Masks['Points']
2895        rings = Masks['Rings']
2896        arcs = Masks['Arcs']
2897        polygons = Masks['Polygons']
2898        if 'Frames' not in Masks:
2899            Masks['Frames'] = []
2900        frame = Masks['Frames']
2901        for spot in spots:
2902            if spot:
2903                x,y,d = spot
2904                Plot.add_artist(Circle((x,y),radius=d/2,fc='none',ec='r',picker=3))
2905        G2frame.ringList = []
2906        for iring,ring in enumerate(rings):
2907            if ring:
2908                tth,thick = ring
2909                wave = Data['wavelength']
2910                xy1 = []
2911                xy2 = []
2912                Azm = np.linspace(0,362,181)
2913                for azm in Azm:
2914                    xy1.append(G2img.GetDetectorXY(Dsp(tth+thick/2.,wave),azm,Data))      #what about hyperbola
2915                    xy2.append(G2img.GetDetectorXY(Dsp(tth-thick/2.,wave),azm,Data))      #what about hyperbola
2916                x1,y1 = np.array(xy1).T
2917                x2,y2 = np.array(xy2).T
2918                G2frame.ringList.append([Plot.plot(x1,y1,'r',picker=3),iring,'o'])           
2919                G2frame.ringList.append([Plot.plot(x2,y2,'r',picker=3),iring,'i'])
2920        G2frame.arcList = []
2921        for iarc,arc in enumerate(arcs):
2922            if arc:
2923                tth,azm,thick = arc           
2924                wave = Data['wavelength']
2925                xy1 = []
2926                xy2 = []
2927                aR = azm[0],azm[1],azm[1]-azm[0]
2928                if azm[1]-azm[0] > 180:
2929                    aR[2] /= 2
2930                Azm = np.linspace(aR[0],aR[1],aR[2])
2931                for azm in Azm:
2932                    xy1.append(G2img.GetDetectorXY(Dsp(tth+thick/2.,wave),azm,Data))      #what about hyperbola
2933                    xy2.append(G2img.GetDetectorXY(Dsp(tth-thick/2.,wave),azm,Data))      #what about hyperbola
2934                x1,y1 = np.array(xy1).T
2935                x2,y2 = np.array(xy2).T
2936                G2frame.arcList.append([Plot.plot(x1,y1,'r',picker=3),iarc,'o'])           
2937                G2frame.arcList.append([Plot.plot(x2,y2,'r',picker=3),iarc,'i'])
2938                G2frame.arcList.append([Plot.plot([x1[0],x2[0]],[y1[0],y2[0]],'r',picker=3),iarc,'l'])
2939                G2frame.arcList.append([Plot.plot([x1[-1],x2[-1]],[y1[-1],y2[-1]],'r',picker=3),iarc,'u'])
2940        G2frame.polyList = []
2941        for ipoly,polygon in enumerate(polygons):
2942            if polygon:
2943                x,y = np.hsplit(np.array(polygon),2)
2944                G2frame.polyList.append([Plot.plot(x,y,'r+',picker=10),ipoly])
2945                Plot.plot(x,y,'r')           
2946        G2frame.frameList = []
2947        if frame:
2948            x,y = np.hsplit(np.array(frame),2)
2949            G2frame.frameList.append([Plot.plot(x,y,'g+',picker=10),0])
2950            Plot.plot(x,y,'g')           
2951        if newImage:
2952            colorBar = Page.figure.colorbar(Img)
2953        Plot.set_xlim(xlim)
2954        Plot.set_ylim(ylim)
2955        if Data['invert_x']:
2956            Plot.invert_xaxis()
2957        if Data['invert_y']:
2958            Plot.invert_yaxis()
2959        if not newPlot and xylim:
2960            Page.toolbar.push_current()
2961            Plot.set_xlim(xylim[0])
2962            Plot.set_ylim(xylim[1])
2963            xylim = []
2964            Page.toolbar.push_current()
2965            Page.toolbar.draw()
2966        else:
2967            Page.canvas.draw()
2968    finally:
2969        wx.EndBusyCursor()
2970       
2971################################################################################
2972##### PlotIntegration
2973################################################################################
2974           
2975def PlotIntegration(G2frame,newPlot=False,event=None):
2976    '''Plot of 2D image after image integration with 2-theta and azimuth as coordinates
2977    '''
2978           
2979    def OnMotion(event):
2980        Page.canvas.SetToolTipString('')
2981        Page.canvas.SetCursor(wx.CROSS_CURSOR)
2982        azm = event.ydata
2983        tth = event.xdata
2984        if azm and tth:
2985            G2frame.G2plotNB.status.SetFields(\
2986                ['','Detector 2-th =%9.3fdeg, azm = %7.2fdeg'%(tth,azm)])
2987                               
2988    try:
2989        plotNum = G2frame.G2plotNB.plotList.index('2D Integration')
2990        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
2991        if not newPlot:
2992            Plot = Page.figure.gca()          #get previous plot & get limits
2993            xylim = Plot.get_xlim(),Plot.get_ylim()
2994        Page.figure.clf()
2995        Plot = Page.figure.gca()          #get a fresh plot after clf()
2996       
2997    except ValueError:
2998        Plot = G2frame.G2plotNB.addMpl('2D Integration').gca()
2999        plotNum = G2frame.G2plotNB.plotList.index('2D Integration')
3000        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
3001        Page.canvas.mpl_connect('motion_notify_event', OnMotion)
3002        Page.views = False
3003        view = False
3004    Page.Choice = None
3005    if not event:
3006        Page.SetFocus()
3007       
3008    Data = G2frame.PatternTree.GetItemPyData(
3009        G2gd.GetPatternTreeItemId(G2frame,G2frame.Image, 'Image Controls'))
3010    image = G2frame.Integrate[0]
3011    xsc = G2frame.Integrate[1]
3012    ysc = G2frame.Integrate[2]
3013    Imin,Imax = Data['range'][1]
3014    acolor = mpl.cm.get_cmap(Data['color'])
3015    Plot.set_title(G2frame.PatternTree.GetItemText(G2frame.Image)[4:])
3016    Plot.set_ylabel('azimuth',fontsize=12)
3017    Plot.set_xlabel('2-theta',fontsize=12)
3018    Img = Plot.imshow(image,cmap=acolor,vmin=Imin,vmax=Imax,interpolation='nearest', \
3019        extent=[ysc[0],ysc[-1],xsc[-1],xsc[0]],aspect='auto')
3020    colorBar = Page.figure.colorbar(Img)
3021#    if Data['ellipses']:           
3022#        for ellipse in Data['ellipses']:
3023#            x,y = np.array(G2img.makeIdealRing(ellipse[:3])) #skip color
3024#            tth,azm = G2img.GetTthAzm(x,y,Data)
3025##            azm = np.where(azm < 0.,azm+360,azm)
3026#            Plot.plot(tth,azm,'b,')
3027    if not newPlot:
3028        Page.toolbar.push_current()
3029        Plot.set_xlim(xylim[0])
3030        Plot.set_ylim(xylim[1])
3031        xylim = []
3032        Page.toolbar.push_current()
3033        Page.toolbar.draw()
3034    else:
3035        Page.canvas.draw()
3036               
3037################################################################################
3038##### PlotTRImage
3039################################################################################
3040           
3041def PlotTRImage(G2frame,tax,tay,taz,newPlot=False):
3042    '''a test plot routine - not normally used
3043    ''' 
3044           
3045    def OnMotion(event):
3046        Page.canvas.SetToolTipString('')
3047        Page.canvas.SetCursor(wx.CROSS_CURSOR)
3048        azm = event.xdata
3049        tth = event.ydata
3050        if azm and tth:
3051            G2frame.G2plotNB.status.SetFields(\
3052                ['','Detector 2-th =%9.3fdeg, azm = %7.2fdeg'%(tth,azm)])
3053                               
3054    try:
3055        plotNum = G2frame.G2plotNB.plotList.index('2D Transformed Powder Image')
3056        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
3057        if not newPlot:
3058            Plot = Page.figure.gca()          #get previous plot & get limits
3059            xylim = Plot.get_xlim(),Plot.get_ylim()
3060        Page.figure.clf()
3061        Plot = Page.figure.gca()          #get a fresh plot after clf()
3062       
3063    except ValueError:
3064        Plot = G2frame.G2plotNB.addMpl('2D Transformed Powder Image').gca()
3065        plotNum = G2frame.G2plotNB.plotList.index('2D Transformed Powder Image')
3066        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
3067        Page.canvas.mpl_connect('motion_notify_event', OnMotion)
3068        Page.views = False
3069        view = False
3070    Page.Choice = None
3071    Page.SetFocus()
3072       
3073    Data = G2frame.PatternTree.GetItemPyData(
3074        G2gd.GetPatternTreeItemId(G2frame,G2frame.Image, 'Image Controls'))
3075    Imin,Imax = Data['range'][1]
3076    step = (Imax-Imin)/5.
3077    V = np.arange(Imin,Imax,step)
3078    acolor = mpl.cm.get_cmap(Data['color'])
3079    Plot.set_title(G2frame.PatternTree.GetItemText(G2frame.Image)[4:])
3080    Plot.set_xlabel('azimuth',fontsize=12)
3081    Plot.set_ylabel('2-theta',fontsize=12)
3082    Plot.contour(tax,tay,taz,V,cmap=acolor)
3083    if Data['showLines']:
3084        IOtth = Data['IOtth']
3085        if Data['fullIntegrate']:
3086            LRAzim = [-180,180]
3087        else:
3088            LRAzim = Data['LRazimuth']                  #NB: integers
3089        Plot.plot([LRAzim[0],LRAzim[1]],[IOtth[0],IOtth[0]],picker=True)
3090        Plot.plot([LRAzim[0],LRAzim[1]],[IOtth[1],IOtth[1]],picker=True)
3091        if not Data['fullIntegrate']:
3092            Plot.plot([LRAzim[0],LRAzim[0]],[IOtth[0],IOtth[1]],picker=True)
3093            Plot.plot([LRAzim[1],LRAzim[1]],[IOtth[0],IOtth[1]],picker=True)
3094    if Data['setRings']:
3095        rings = np.concatenate((Data['rings']),axis=0)
3096        for xring,yring,dsp in rings:
3097            x,y = G2img.GetTthAzm(xring,yring,Data)
3098            Plot.plot(y,x,'r+')           
3099    if Data['ellipses']:           
3100        for ellipse in Data['ellipses']:
3101            ring = np.array(G2img.makeIdealRing(ellipse[:3])) #skip color
3102            x,y = np.hsplit(ring,2)
3103            tth,azm = G2img.GetTthAzm(x,y,Data)
3104            Plot.plot(azm,tth,'b,')
3105    if not newPlot:
3106        Page.toolbar.push_current()
3107        Plot.set_xlim(xylim[0])
3108        Plot.set_ylim(xylim[1])
3109        xylim = []
3110        Page.toolbar.push_current()
3111        Page.toolbar.draw()
3112    else:
3113        Page.canvas.draw()
3114       
3115################################################################################
3116##### PlotStructure
3117################################################################################
3118           
3119def PlotStructure(G2frame,data):
3120    '''Crystal structure plotting package. Can show structures as balls, sticks, lines,
3121    thermal motion ellipsoids and polyhedra
3122    '''
3123
3124    def FindPeaksBonds(XYZ):
3125        rFact = data['Drawing']['radiusFactor']
3126        Bonds = [[] for x in XYZ]
3127        for i,xyz in enumerate(XYZ):
3128            Dx = XYZ-xyz
3129            dist = np.sqrt(np.sum(np.inner(Dx,Amat)**2,axis=1))
3130            IndB = ma.nonzero(ma.masked_greater(dist,rFact*2.2))
3131            for j in IndB[0]:
3132                Bonds[i].append(Dx[j]/2.)
3133                Bonds[j].append(-Dx[j]/2.)
3134        return Bonds
3135
3136    ForthirdPI = 4.0*math.pi/3.0
3137    generalData = data['General']
3138    cell = generalData['Cell'][1:7]
3139    Vol = generalData['Cell'][7:8][0]
3140    Amat,Bmat = G2lat.cell2AB(cell)         #Amat - crystal to cartesian, Bmat - inverse
3141    Gmat,gmat = G2lat.cell2Gmat(cell)
3142    A4mat = np.concatenate((np.concatenate((Amat,[[0],[0],[0]]),axis=1),[[0,0,0,1],]),axis=0)
3143    B4mat = np.concatenate((np.concatenate((Bmat,[[0],[0],[0]]),axis=1),[[0,0,0,1],]),axis=0)
3144    SGData = generalData['SGData']
3145    Mydir = generalData['Mydir']
3146    atomData = data['Atoms']
3147    mapPeaks = []
3148    drawingData = data['Drawing']   
3149    if 'Map Peaks' in data:
3150        mapPeaks = np.array(data['Map Peaks'])
3151        peakMax = 100.
3152        if len(mapPeaks):
3153            peakMax = np.max(mapPeaks.T[0])
3154    resRBData = data['RBModels'].get('Residue',[])
3155    vecRBData = data['RBModels'].get('Vector',[])
3156    rbAtmDict = {}
3157    for rbObj in resRBData+vecRBData:
3158        exclList = ['X' for i in range(len(rbObj['Ids']))]
3159        rbAtmDict.update(dict(zip(rbObj['Ids'],exclList)))
3160    testRBObj = data.get('testRBObj',{})
3161    rbObj = testRBObj.get('rbObj',{})
3162    MCSA = data.get('MCSA',{})
3163    mcsaModels = MCSA.get('Models',[])
3164    if mcsaModels:
3165        XYZs,Types = G2mth.UpdateMCSAxyz(Bmat,MCSA)
3166        mcsaXYZ = []
3167        mcsaTypes = []
3168        for xyz,atyp in zip(XYZs,Types):
3169            for item in G2spc.GenAtom(xyz,SGData):
3170                mcsaXYZ.append(item[0]) 
3171                mcsaTypes.append(atyp)
3172        mcsaBonds = FindPeaksBonds(mcsaXYZ)       
3173    drawAtoms = drawingData.get('Atoms',[])
3174    mapData = {}
3175    flipData = {}
3176    rhoXYZ = []
3177    showBonds = False
3178    if 'Map' in generalData:
3179        mapData = generalData['Map']
3180        showBonds = mapData.get('Show bonds',False)
3181    if 'Flip' in generalData:
3182        flipData = generalData['Flip']                       
3183        flipData['mapRoll'] = [0,0,0]
3184    Wt = np.array([255,255,255])
3185    Rd = np.array([255,0,0])
3186    Gr = np.array([0,255,0])
3187    wxGreen = wx.Color(0,255,0)
3188    Bl = np.array([0,0,255])
3189    Or = np.array([255,128,0])
3190    wxOrange = wx.Color(255,128,0)
3191    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]])
3192    uEdges = np.array([
3193        [uBox[0],uBox[1]],[uBox[0],uBox[3]],[uBox[0],uBox[4]],[uBox[1],uBox[2]], 
3194        [uBox[2],uBox[3]],[uBox[1],uBox[5]],[uBox[2],uBox[6]],[uBox[3],uBox[7]], 
3195        [uBox[4],uBox[5]],[uBox[5],uBox[6]],[uBox[6],uBox[7]],[uBox[7],uBox[4]]])
3196    mD = 0.1
3197    mV = np.array([[[-mD,0,0],[mD,0,0]],[[0,-mD,0],[0,mD,0]],[[0,0,-mD],[0,0,mD]]])
3198    mapPeakVecs = np.inner(mV,Bmat)
3199
3200    uColors = [Rd,Gr,Bl,Wt, Wt,Wt,Wt,Wt, Wt,Wt,Wt,Wt]
3201    altDown = False
3202    shiftDown = False
3203    ctrlDown = False
3204   
3205    def OnKeyBox(event):
3206        import Image
3207#        Draw()                          #make sure plot is fresh!!
3208        mode = cb.GetValue()
3209        if mode in ['jpeg','bmp','tiff',]:
3210            Fname = os.path.join(Mydir,generalData['Name']+'.'+mode)
3211            size = Page.canvas.GetSize()
3212            glPixelStorei(GL_UNPACK_ALIGNMENT, 1)
3213            if mode in ['jpeg',]:
3214                Pix = glReadPixels(0,0,size[0],size[1],GL_RGBA, GL_UNSIGNED_BYTE)
3215                im = Image.new("RGBA", (size[0],size[1]))
3216            else:
3217                Pix = glReadPixels(0,0,size[0],size[1],GL_RGB, GL_UNSIGNED_BYTE)
3218                im = Image.new("RGB", (size[0],size[1]))
3219            im.fromstring(Pix)
3220            im.save(Fname,mode)
3221            cb.SetValue(' save as/key:')
3222            G2frame.G2plotNB.status.SetStatusText('Drawing saved to: '+Fname,1)
3223        else:
3224            event.key = cb.GetValue()[0]
3225            cb.SetValue(' save as/key:')
3226            wx.CallAfter(OnKey,event)
3227        Page.canvas.SetFocus() # redirect the Focus from the button back to the plot
3228
3229    def OnKey(event):           #on key UP!!
3230#        Draw()                          #make sure plot is fresh!!
3231        try:
3232            keyCode = event.GetKeyCode()
3233            if keyCode > 255:
3234                keyCode = 0
3235            key = chr(keyCode)
3236        except AttributeError:       #if from OnKeyBox above
3237            key = str(event.key).upper()
3238        indx = drawingData['selectedAtoms']
3239        cx,ct = drawingData['atomPtrs'][:2]
3240        if key in ['C']:
3241            drawingData['viewPoint'] = [[.5,.5,.5],[0,0]]
3242            drawingData['viewDir'] = [0,0,1]
3243            drawingData['oldxy'] = []
3244            V0 = np.array([0,0,1])
3245            V = np.inner(Amat,V0)
3246            V /= np.sqrt(np.sum(V**2))
3247            A = np.arccos(np.sum(V*V0))
3248            Q = G2mth.AV2Q(A,[0,1,0])
3249            drawingData['Quaternion'] = Q
3250            SetViewPointText(drawingData['viewPoint'][0])
3251            SetViewDirText(drawingData['viewDir'])
3252            Q = drawingData['Quaternion']
3253            G2frame.G2plotNB.status.SetStatusText('New quaternion: %.2f+, %.2fi+ ,%.2fj+, %.2fk'%(Q[0],Q[1],Q[2],Q[3]),1)
3254        elif key in ['N']:
3255            drawAtoms = drawingData['Atoms']
3256            if not len(drawAtoms):      #no atoms
3257                return
3258            pI = drawingData['viewPoint'][1]
3259            if not len(pI):
3260                pI = [0,0]
3261            if indx:
3262                pI[0] = indx[pI[1]]
3263                Tx,Ty,Tz = drawAtoms[pI[0]][cx:cx+3]
3264                pI[1] += 1
3265                if pI[1] >= len(indx):
3266                    pI[1] = 0
3267            else:
3268                Tx,Ty,Tz = drawAtoms[pI[0]][cx:cx+3]               
3269                pI[0] += 1
3270                if pI[0] >= len(drawAtoms):
3271                    pI[0] = 0
3272            drawingData['viewPoint'] = [[Tx,Ty,Tz],pI]
3273            SetViewPointText(drawingData['viewPoint'][0])
3274            G2frame.G2plotNB.status.SetStatusText('View point at atom '+drawAtoms[pI[0]][ct-1]+str(pI),1)
3275               
3276        elif key in ['P']:
3277            drawAtoms = drawingData['Atoms']
3278            if not len(drawAtoms):      #no atoms
3279                return
3280            pI = drawingData['viewPoint'][1]
3281            if not len(pI):
3282                pI = [0,0]
3283            if indx:
3284                pI[0] = indx[pI[1]]
3285                Tx,Ty,Tz = drawAtoms[pI[0]][cx:cx+3]
3286                pI[1] -= 1
3287                if pI[1] < 0:
3288                    pI[1] = len(indx)-1
3289            else:
3290                Tx,Ty,Tz = drawAtoms[pI[0]][cx:cx+3]               
3291                pI[0] -= 1
3292                if pI[0] < 0:
3293                    pI[0] = len(drawAtoms)-1
3294            drawingData['viewPoint'] = [[Tx,Ty,Tz],pI]
3295            SetViewPointText(drawingData['viewPoint'][0])           
3296            G2frame.G2plotNB.status.SetStatusText('View point at atom '+drawAtoms[pI[0]][ct-1]+str(pI),1)
3297        elif key in ['U','D','L','R'] and mapData['Flip'] == True:
3298            dirDict = {'U':[0,1],'D':[0,-1],'L':[-1,0],'R':[1,0]}
3299            SetMapRoll(dirDict[key])
3300            SetPeakRoll(dirDict[key])
3301            SetMapPeaksText(mapPeaks)
3302        Draw('key')
3303           
3304    def GetTruePosition(xy,Add=False):
3305        View = glGetIntegerv(GL_VIEWPORT)
3306        Proj = glGetDoublev(GL_PROJECTION_MATRIX)
3307        Model = glGetDoublev(GL_MODELVIEW_MATRIX)
3308        Zmax = 1.
3309        if Add:
3310            Indx = GetSelectedAtoms()
3311        if G2frame.dataDisplay.GetPageText(getSelection()) == 'Map peaks':
3312            for i,peak in enumerate(mapPeaks):
3313                x,y,z = peak[1:4]
3314                X,Y,Z = gluProject(x,y,z,Model,Proj,View)
3315                XY = [int(X),int(View[3]-Y)]
3316                if np.allclose(xy,XY,atol=10) and Z < Zmax:
3317                    Zmax = Z
3318                    try:
3319                        Indx.remove(i)
3320                        ClearSelectedAtoms()
3321                        for id in Indx:
3322                            SetSelectedAtoms(id,Add)
3323                    except:
3324                        SetSelectedAtoms(i,Add)
3325        else:
3326            cx = drawingData['atomPtrs'][0]
3327            for i,atom in enumerate(drawAtoms):
3328                x,y,z = atom[cx:cx+3]
3329                X,Y,Z = gluProject(x,y,z,Model,Proj,View)
3330                XY = [int(X),int(View[3]-Y)]
3331                if np.allclose(xy,XY,atol=10) and Z < Zmax:
3332                    Zmax = Z
3333                    try:
3334                        Indx.remove(i)
3335                        ClearSelectedAtoms()
3336                        for id in Indx:
3337                            SetSelectedAtoms(id,Add)
3338                    except:
3339                        SetSelectedAtoms(i,Add)
3340                                       
3341    def OnMouseDown(event):
3342        xy = event.GetPosition()
3343        if event.ShiftDown():
3344            if event.LeftIsDown():
3345                GetTruePosition(xy)
3346            elif event.RightIsDown():
3347                GetTruePosition(xy,True)
3348        else:
3349            drawingData['oldxy'] = list(xy)
3350       
3351    def OnMouseMove(event):
3352        if event.ShiftDown():           #don't want any inadvertant moves when picking
3353            return
3354        newxy = event.GetPosition()
3355                               
3356        if event.Dragging():
3357            if event.AltDown() and rbObj:
3358                if event.LeftIsDown():
3359                    SetRBRotation(newxy)
3360                    Q = rbObj['Orient'][0]
3361                    G2frame.G2plotNB.status.SetStatusText('New quaternion: %.2f+, %.2fi+ ,%.2fj+, %.2fk'%(Q[0],Q[1],Q[2],Q[3]),1)
3362                elif event.RightIsDown():
3363                    SetRBTranslation(newxy)
3364                    Tx,Ty,Tz = rbObj['Orig'][0]
3365                    G2frame.G2plotNB.status.SetStatusText('New view point: %.4f, %.4f, %.4f'%(Tx,Ty,Tz),1)
3366                elif event.MiddleIsDown():
3367                    SetRBRotationZ(newxy)
3368                    Q = rbObj['Orient'][0]
3369                    G2frame.G2plotNB.status.SetStatusText('New quaternion: %.2f+, %.2fi+ ,%.2fj+, %.2fk'%(Q[0],Q[1],Q[2],Q[3]),1)
3370                Draw('move')
3371            elif not event.ControlDown():
3372                if event.LeftIsDown():
3373                    SetRotation(newxy)
3374                    Q = drawingData['Quaternion']
3375                    G2frame.G2plotNB.status.SetStatusText('New quaternion: %.2f+, %.2fi+ ,%.2fj+, %.2fk'%(Q[0],Q[1],Q[2],Q[3]),1)
3376                elif event.RightIsDown():
3377                    SetTranslation(newxy)
3378                    Tx,Ty,Tz = drawingData['viewPoint'][0]
3379                    G2frame.G2plotNB.status.SetStatusText('New view point: %.4f, %.4f, %.4f'%(Tx,Ty,Tz),1)
3380                elif event.MiddleIsDown():
3381                    SetRotationZ(newxy)
3382                    Q = drawingData['Quaternion']
3383                    G2frame.G2plotNB.status.SetStatusText('New quaternion: %.2f+, %.2fi+ ,%.2fj+, %.2fk'%(Q[0],Q[1],Q[2],Q[3]),1)
3384                Draw('move')
3385           
3386       
3387    def OnMouseWheel(event):
3388        if event.ShiftDown():
3389            return
3390        drawingData['cameraPos'] += event.GetWheelRotation()/24.
3391        drawingData['cameraPos'] = max(10,min(500,drawingData['cameraPos']))
3392        G2frame.G2plotNB.status.SetStatusText('New camera distance: %.2f'%(drawingData['cameraPos']),1)
3393        page = getSelection()
3394        if page:
3395            if G2frame.dataDisplay.GetPageText(page) == 'Draw Options':
3396                G2frame.dataDisplay.cameraPosTxt.SetLabel('Camera Position: '+'%.2f'%(drawingData['cameraPos']))
3397                G2frame.dataDisplay.cameraSlider.SetValue(drawingData['cameraPos'])
3398            Draw('wheel')
3399       
3400    def getSelection():
3401        try:
3402            return G2frame.dataDisplay.GetSelection()
3403        except AttributeError:
3404            G2frame.G2plotNB.status.SetStatusText('Select this from Phase data window!',1)
3405            return 0
3406           
3407    def SetViewPointText(VP):
3408        page = getSelection()
3409        if page:
3410            if G2frame.dataDisplay.GetPageText(page) == 'Draw Options':
3411                G2frame.dataDisplay.viewPoint.SetValue('%.3f %.3f %.3f'%(VP[0],VP[1],VP[2]))
3412               
3413    def SetRBOrigText():
3414        page = getSelection()
3415        if page:
3416            if G2frame.dataDisplay.GetPageText(page) == 'RB Models':
3417                for i,sizer in enumerate(testRBObj['Sizers']['Xsizers']):
3418                    sizer.SetValue('%8.5f'%(testRBObj['rbObj']['Orig'][0][i]))
3419                   
3420    def SetRBOrienText():
3421        page = getSelection()
3422        if page:
3423            if G2frame.dataDisplay.GetPageText(page) == 'RB Models':
3424                for i,sizer in enumerate(testRBObj['Sizers']['Osizers']):
3425                    sizer.SetValue('%8.5f'%(testRBObj['rbObj']['Orient'][0][i]))
3426               
3427    def SetViewDirText(VD):
3428        page = getSelection()
3429        if page:
3430            if G2frame.dataDisplay.GetPageText(page) == 'Draw Options':
3431                G2frame.dataDisplay.viewDir.SetValue('%.3f %.3f %.3f'%(VD[0],VD[1],VD[2]))
3432               
3433    def SetMapPeaksText(mapPeaks):
3434        page = getSelection()
3435        if page:
3436            if G2frame.dataDisplay.GetPageText(page) == 'Map peaks':
3437                G2frame.MapPeaksTable.SetData(mapPeaks)
3438                panel = G2frame.dataDisplay.GetPage(page).GetChildren()
3439                names = [child.GetName() for child in panel]
3440                panel[names.index('grid window')].Refresh()
3441           
3442    def ClearSelectedAtoms():
3443        page = getSelection()
3444        if page:
3445            if G2frame.dataDisplay.GetPageText(page) == 'Draw Atoms':
3446                G2frame.dataDisplay.GetPage(page).ClearSelection()      #this is the Atoms grid in Draw Atoms
3447            elif G2frame.dataDisplay.GetPageText(page) == 'Map peaks':
3448                G2frame.dataDisplay.GetPage(page).ClearSelection()      #this is the Atoms grid in Atoms
3449            elif G2frame.dataDisplay.GetPageText(page) == 'Atoms':
3450                G2frame.dataDisplay.GetPage(page).ClearSelection()      #this is the Atoms grid in Atoms
3451               
3452                   
3453    def SetSelectedAtoms(ind,Add=False):
3454        page = getSelection()
3455        if page:
3456            if G2frame.dataDisplay.GetPageText(page) == 'Draw Atoms':
3457                G2frame.dataDisplay.GetPage(page).SelectRow(ind,Add)      #this is the Atoms grid in Draw Atoms
3458            elif G2frame.dataDisplay.GetPageText(page) == 'Map peaks':
3459                G2frame.dataDisplay.GetPage(page).SelectRow(ind,Add)                 
3460            elif G2frame.dataDisplay.GetPageText(page) == 'Atoms':
3461                Id = drawAtoms[ind][-3]
3462                for i,atom in enumerate(atomData):
3463                    if atom[-1] == Id:
3464                        G2frame.dataDisplay.GetPage(page).SelectRow(i)      #this is the Atoms grid in Atoms
3465                 
3466    def GetSelectedAtoms():
3467        page = getSelection()
3468        Ind = []
3469        if page:
3470            if G2frame.dataDisplay.GetPageText(page) == 'Draw Atoms':
3471                Ind = G2frame.dataDisplay.GetPage(page).GetSelectedRows()      #this is the Atoms grid in Draw Atoms
3472            elif G2frame.dataDisplay.GetPageText(page) == 'Map peaks':
3473                Ind = G2frame.dataDisplay.GetPage(page).GetSelectedRows()
3474            elif G2frame.dataDisplay.GetPageText(page) == 'Atoms':
3475                Ind = G2frame.dataDisplay.GetPage(page).GetSelectedRows()      #this is the Atoms grid in Atoms
3476        return Ind
3477                                       
3478    def SetBackground():
3479        R,G,B,A = Page.camera['backColor']
3480        glClearColor(R,G,B,A)
3481        glClear(GL_COLOR_BUFFER_BIT | GL_DEPTH_BUFFER_BIT)
3482       
3483    def SetLights():
3484        glEnable(GL_DEPTH_TEST)
3485        glShadeModel(GL_SMOOTH)
3486        glEnable(GL_LIGHTING)
3487        glEnable(GL_LIGHT0)
3488        glLightModeli(GL_LIGHT_MODEL_TWO_SIDE,0)
3489        glLightfv(GL_LIGHT0,GL_AMBIENT,[1,1,1,.8])
3490        glLightfv(GL_LIGHT0,GL_DIFFUSE,[1,1,1,1])
3491       
3492    def GetRoll(newxy,rho):
3493        Q = drawingData['Quaternion']
3494        dxy = G2mth.prodQVQ(G2mth.invQ(Q),np.inner(Bmat,newxy+[0,]))
3495        dxy = np.array(dxy*rho.shape)       
3496        roll = np.where(dxy>0.5,1,np.where(dxy<-.5,-1,0))
3497        return roll
3498               
3499    def SetMapRoll(newxy):
3500        rho = mapData['rho']
3501        roll = GetRoll(newxy,rho)
3502        mapData['rho'] = np.roll(np.roll(np.roll(rho,roll[0],axis=0),roll[1],axis=1),roll[2],axis=2)
3503        drawingData['oldxy'] = list(newxy)
3504       
3505    def SetPeakRoll(newxy):
3506        rho = mapData['rho']
3507        roll = GetRoll(newxy,rho)
3508        steps = 1./np.array(rho.shape)
3509        dxy = roll*steps
3510        for peak in mapPeaks:
3511            peak[1:4] += dxy
3512            peak[1:4] %= 1.
3513            peak[4] = np.sqrt(np.sum(np.inner(Amat,peak[1:4])**2))
3514               
3515    def SetTranslation(newxy):
3516#first get translation vector in screen coords.       
3517        oldxy = drawingData['oldxy']
3518        if not len(oldxy): oldxy = list(newxy)
3519        dxy = newxy-oldxy
3520        drawingData['oldxy'] = list(newxy)
3521        V = np.array([-dxy[0],dxy[1],0.])
3522#then transform to rotated crystal coordinates & apply to view point       
3523        Q = drawingData['Quaternion']
3524        V = np.inner(Bmat,G2mth.prodQVQ(G2mth.invQ(Q),V))
3525        Tx,Ty,Tz = drawingData['viewPoint'][0]
3526        Tx += V[0]*0.01
3527        Ty += V[1]*0.01
3528        Tz += V[2]*0.01
3529        drawingData['viewPoint'][0] =  Tx,Ty,Tz
3530        SetViewPointText([Tx,Ty,Tz])
3531       
3532    def SetRBTranslation(newxy):
3533#first get translation vector in screen coords.       
3534        oldxy = drawingData['oldxy']
3535        if not len(oldxy): oldxy = list(newxy)
3536        dxy = newxy-oldxy
3537        drawingData['oldxy'] = list(newxy)
3538        V = np.array([-dxy[0],dxy[1],0.])
3539#then transform to rotated crystal coordinates & apply to RB origin       
3540        Q = drawingData['Quaternion']
3541        V = np.inner(Bmat,G2mth.prodQVQ(G2mth.invQ(Q),V))
3542        Tx,Ty,Tz = rbObj['Orig'][0]
3543        Tx -= V[0]*0.01
3544        Ty -= V[1]*0.01
3545        Tz -= V[2]*0.01
3546        rbObj['Orig'][0] =  Tx,Ty,Tz
3547        SetRBOrigText()
3548       
3549    def SetRotation(newxy):
3550#first get rotation vector in screen coords. & angle increment       
3551        oldxy = drawingData['oldxy']
3552        if not len(oldxy): oldxy = list(newxy)
3553        dxy = newxy-oldxy
3554        drawingData['oldxy'] = list(newxy)
3555        V = np.array([dxy[1],dxy[0],0.])
3556        A = 0.25*np.sqrt(dxy[0]**2+dxy[1]**2)
3557# next transform vector back to xtal coordinates via inverse quaternion
3558# & make new quaternion
3559        Q = drawingData['Quaternion']
3560        V = G2mth.prodQVQ(G2mth.invQ(Q),np.inner(Bmat,V))
3561        DQ = G2mth.AVdeg2Q(A,V)
3562        Q = G2mth.prodQQ(Q,DQ)
3563        drawingData['Quaternion'] = Q
3564# finally get new view vector - last row of rotation matrix
3565        VD = np.inner(Bmat,G2mth.Q2Mat(Q)[2])
3566        VD /= np.sqrt(np.sum(VD**2))
3567        drawingData['viewDir'] = VD
3568        SetViewDirText(VD)
3569       
3570    def SetRotationZ(newxy):                       
3571#first get rotation vector (= view vector) in screen coords. & angle increment       
3572        View = glGetIntegerv(GL_VIEWPORT)
3573        cent = [View[2]/2,View[3]/2]
3574        oldxy = drawingData['oldxy']
3575        if not len(oldxy): oldxy = list(newxy)
3576        dxy = newxy-oldxy
3577        drawingData['oldxy'] = list(newxy)
3578        V = drawingData['viewDir']
3579        A = [0,0]
3580        A[0] = dxy[1]*.25
3581        A[1] = dxy[0]*.25
3582        if newxy[0] > cent[0]:
3583            A[0] *= -1
3584        if newxy[1] < cent[1]:
3585            A[1] *= -1       
3586# next transform vector back to xtal coordinates & make new quaternion
3587        Q = drawingData['Quaternion']
3588        V = np.inner(Amat,V)
3589        Qx = G2mth.AVdeg2Q(A[0],V)
3590        Qy = G2mth.AVdeg2Q(A[1],V)
3591        Q = G2mth.prodQQ(Q,Qx)
3592        Q = G2mth.prodQQ(Q,Qy)
3593        drawingData['Quaternion'] = Q
3594
3595    def SetRBRotation(newxy):
3596#first get rotation vector in screen coords. & angle increment       
3597        oldxy = drawingData['oldxy']
3598        if not len(oldxy): oldxy = list(newxy)
3599        dxy = newxy-oldxy
3600        drawingData['oldxy'] = list(newxy)
3601        V = np.array([dxy[1],dxy[0],0.])
3602        A = 0.25*np.sqrt(dxy[0]**2+dxy[1]**2)
3603# next transform vector back to xtal coordinates via inverse quaternion
3604# & make new quaternion
3605        Q = rbObj['Orient'][0]              #rotate RB to Cart
3606        QC = drawingData['Quaternion']      #rotate Cart to drawing
3607        V = G2mth.prodQVQ(G2mth.invQ(QC),V)
3608        V = G2mth.prodQVQ(G2mth.invQ(Q),V)
3609        DQ = G2mth.AVdeg2Q(A,V)
3610        Q = G2mth.prodQQ(Q,DQ)
3611        rbObj['Orient'][0] = Q
3612        SetRBOrienText()
3613       
3614    def SetRBRotationZ(newxy):                       
3615#first get rotation vector (= view vector) in screen coords. & angle increment       
3616        View = glGetIntegerv(GL_VIEWPORT)
3617        cent = [View[2]/2,View[3]/2]
3618        oldxy = drawingData['oldxy']
3619        if not len(oldxy): oldxy = list(newxy)
3620        dxy = newxy-oldxy
3621        drawingData['oldxy'] = list(newxy)
3622        V = drawingData['viewDir']
3623        A = [0,0]
3624        A[0] = dxy[1]*.25
3625        A[1] = dxy[0]*.25
3626        if newxy[0] < cent[0]:
3627            A[0] *= -1
3628        if newxy[1] > cent[1]:
3629            A[1] *= -1       
3630# next transform vector back to RB coordinates & make new quaternion
3631        Q = rbObj['Orient'][0]              #rotate RB to cart
3632        V = np.inner(Amat,V)
3633        V = -G2mth.prodQVQ(G2mth.invQ(Q),V)
3634        Qx = G2mth.AVdeg2Q(A[0],V)
3635        Qy = G2mth.AVdeg2Q(A[1],V)
3636        Q = G2mth.prodQQ(Q,Qx)
3637        Q = G2mth.prodQQ(Q,Qy)
3638        rbObj['Orient'][0] = Q
3639        SetRBOrienText()
3640
3641    def RenderBox():
3642        glEnable(GL_COLOR_MATERIAL)
3643        glLineWidth(2)
3644        glEnable(GL_BLEND)
3645        glBlendFunc(GL_SRC_ALPHA,GL_ONE_MINUS_SRC_ALPHA)
3646        glEnable(GL_LINE_SMOOTH)
3647        glBegin(GL_LINES)
3648        for line,color in zip(uEdges,uColors):
3649            glColor3ubv(color)
3650            glVertex3fv(line[0])
3651            glVertex3fv(line[1])
3652        glEnd()
3653        glColor4ubv([0,0,0,0])
3654        glDisable(GL_LINE_SMOOTH)
3655        glDisable(GL_BLEND)
3656        glDisable(GL_COLOR_MATERIAL)
3657       
3658    def RenderUnitVectors(x,y,z):
3659        xyz = np.array([x,y,z])
3660        glEnable(GL_COLOR_MATERIAL)
3661        glLineWidth(1)
3662        glPushMatrix()
3663        glTranslate(x,y,z)
3664        glScalef(1/cell[0],1/cell[1],1/cell[2])
3665        glBegin(GL_LINES)
3666        for line,color in zip(uEdges,uColors)[:3]:
3667            glColor3ubv(color)
3668            glVertex3fv(-line[1]/2.)
3669            glVertex3fv(line[1]/2.)
3670        glEnd()
3671        glPopMatrix()
3672        glColor4ubv([0,0,0,0])
3673        glDisable(GL_COLOR_MATERIAL)
3674               
3675    def RenderSphere(x,y,z,radius,color):
3676        glMaterialfv(GL_FRONT_AND_BACK,GL_DIFFUSE,color)
3677        glPushMatrix()
3678        glTranslate(x,y,z)
3679        glMultMatrixf(B4mat.T)
3680        q = gluNewQuadric()
3681        gluSphere(q,radius,20,10)
3682        glPopMatrix()
3683       
3684    def RenderDots(XYZ,RC):
3685        glEnable(GL_COLOR_MATERIAL)
3686        XYZ = np.array(XYZ)
3687        glPushMatrix()
3688        for xyz,rc in zip(XYZ,RC):
3689            x,y,z = xyz
3690            r,c = rc
3691            glColor3ubv(c)
3692            glPointSize(r*50)
3693            glBegin(GL_POINTS)
3694            glVertex3fv(xyz)
3695            glEnd()
3696        glPopMatrix()
3697        glColor4ubv([0,0,0,0])
3698        glDisable(GL_COLOR_MATERIAL)
3699       
3700    def RenderSmallSphere(x,y,z,radius,color):
3701        glMaterialfv(GL_FRONT_AND_BACK,GL_DIFFUSE,color)
3702        glPushMatrix()
3703        glTranslate(x,y,z)
3704        glMultMatrixf(B4mat.T)
3705        q = gluNewQuadric()
3706        gluSphere(q,radius,4,2)
3707        glPopMatrix()
3708               
3709    def RenderEllipsoid(x,y,z,ellipseProb,E,R4,color):
3710        s1,s2,s3 = E
3711        glMaterialfv(GL_FRONT_AND_BACK,GL_DIFFUSE,color)
3712        glPushMatrix()
3713        glTranslate(x,y,z)
3714        glMultMatrixf(B4mat.T)
3715        glMultMatrixf(R4.T)
3716        glEnable(GL_NORMALIZE)
3717        glScale(s1,s2,s3)
3718        q = gluNewQuadric()
3719        gluSphere(q,ellipseProb,20,10)
3720        glDisable(GL_NORMALIZE)
3721        glPopMatrix()
3722       
3723    def RenderBonds(x,y,z,Bonds,radius,color,slice=20):
3724        glMaterialfv(GL_FRONT_AND_BACK,GL_DIFFUSE,color)
3725        glPushMatrix()
3726        glTranslate(x,y,z)
3727        glMultMatrixf(B4mat.T)
3728        for bond in Bonds:
3729            glPushMatrix()
3730            Dx = np.inner(Amat,bond)
3731            Z = np.sqrt(np.sum(Dx**2))
3732            if Z:
3733                azm = atan2d(-Dx[1],-Dx[0])
3734                phi = acosd(Dx[2]/Z)
3735                glRotate(-azm,0,0,1)
3736                glRotate(phi,1,0,0)
3737                q = gluNewQuadric()
3738                gluCylinder(q,radius,radius,Z,slice,2)
3739            glPopMatrix()           
3740        glPopMatrix()
3741               
3742    def RenderLines(x,y,z,Bonds,color):
3743        glShadeModel(GL_FLAT)
3744        xyz = np.array([x,y,z])
3745        glEnable(GL_COLOR_MATERIAL)
3746        glLineWidth(1)
3747        glColor3fv(color)
3748        glPushMatrix()
3749        glBegin(GL_LINES)
3750        for bond in Bonds:
3751            glVertex3fv(xyz)
3752            glVertex3fv(xyz+bond)
3753        glEnd()
3754        glColor4ubv([0,0,0,0])
3755        glPopMatrix()
3756        glDisable(GL_COLOR_MATERIAL)
3757        glShadeModel(GL_SMOOTH)
3758       
3759    def RenderPolyhedra(x,y,z,Faces,color):
3760        glShadeModel(GL_FLAT)
3761        glPushMatrix()
3762        glTranslate(x,y,z)
3763        glMaterialfv(GL_FRONT_AND_BACK,GL_DIFFUSE,color)
3764        glShadeModel(GL_SMOOTH)
3765        glMultMatrixf(B4mat.T)
3766        for face,norm in Faces:
3767            glPolygonMode(GL_FRONT_AND_BACK,GL_FILL)
3768            glFrontFace(GL_CW)
3769            glNormal3fv(norm)
3770            glBegin(GL_TRIANGLES)
3771            for vert in face:
3772                glVertex3fv(vert)
3773            glEnd()
3774        glPopMatrix()
3775        glShadeModel(GL_SMOOTH)
3776
3777    def RenderMapPeak(x,y,z,color,den):
3778        glShadeModel(GL_FLAT)
3779        xyz = np.array([x,y,z])
3780        glEnable(GL_COLOR_MATERIAL)
3781        glLineWidth(3)
3782        glColor3fv(color*den/255)
3783        glPushMatrix()
3784        glBegin(GL_LINES)
3785        for vec in mapPeakVecs:
3786            glVertex3fv(vec[0]+xyz)
3787            glVertex3fv(vec[1]+xyz)
3788        glEnd()
3789        glColor4ubv([0,0,0,0])
3790        glPopMatrix()
3791        glDisable(GL_COLOR_MATERIAL)
3792        glShadeModel(GL_SMOOTH)
3793       
3794    def RenderBackbone(Backbone,BackboneColor,radius):
3795        glPushMatrix()
3796        glMultMatrixf(B4mat.T)
3797        glEnable(GL_COLOR_MATERIAL)
3798        glShadeModel(GL_SMOOTH)
3799        gleSetJoinStyle(TUBE_NORM_EDGE | TUBE_JN_ANGLE | TUBE_JN_CAP)
3800        glePolyCylinder(Backbone,BackboneColor,radius)
3801        glPopMatrix()       
3802        glDisable(GL_COLOR_MATERIAL)
3803       
3804    def RenderLabel(x,y,z,label,r,color,matRot):
3805        '''
3806        color wx.Color object
3807        '''       
3808        glPushMatrix()
3809        glTranslate(x,y,z)
3810        glMultMatrixf(B4mat.T)
3811        glDisable(GL_LIGHTING)
3812        glRasterPos3f(0,0,0)
3813        glMultMatrixf(matRot)
3814        glRotate(180,1,0,0)             #fix to flip about x-axis
3815        text = gltext.Text(text=label,font=Font,foreground=color)
3816        text.draw_text(scale=0.025)
3817        glEnable(GL_LIGHTING)
3818        glPopMatrix()
3819       
3820    def RenderMap(rho,rhoXYZ,indx,Rok):
3821        glShadeModel(GL_FLAT)
3822        cLevel = drawingData['contourLevel']
3823        XYZ = []
3824        RC = []
3825        for i,xyz in enumerate(rhoXYZ):
3826            if not Rok[i]:
3827                x,y,z = xyz
3828                I,J,K = indx[i]
3829                alpha = 1.0
3830                if cLevel < 1.:
3831                    alpha = (abs(rho[I,J,K])/mapData['rhoMax']-cLevel)/(1.-cLevel)
3832                if rho[I,J,K] < 0.:
3833                    XYZ.append(xyz)
3834                    RC.append([0.1*alpha,Rd])
3835                else:
3836                    XYZ.append(xyz)
3837                    RC.append([0.1*alpha,Gr])
3838        RenderDots(XYZ,RC)
3839        glShadeModel(GL_SMOOTH)
3840                           
3841    def Draw(caller=''):
3842#useful debug?       
3843#        if caller:
3844#            print caller
3845# end of useful debug
3846        mapData = generalData['Map']
3847        pageName = ''
3848        page = getSelection()
3849        if page:
3850            pageName = G2frame.dataDisplay.GetPageText(page)
3851        rhoXYZ = []
3852        if len(mapData['rho']):
3853            VP = np.array(drawingData['viewPoint'][0])-np.array([.5,.5,.5])
3854            contLevel = drawingData['contourLevel']*mapData['rhoMax']
3855            if 'delt-F' in mapData['MapType']:
3856                rho = ma.array(mapData['rho'],mask=(np.abs(mapData['rho'])<contLevel))
3857            else:
3858                rho = ma.array(mapData['rho'],mask=(mapData['rho']<contLevel))
3859            steps = 1./np.array(rho.shape)
3860            incre = np.where(VP>=0,VP%steps,VP%steps-steps)
3861            Vsteps = -np.array(VP/steps,dtype='i')
3862            rho = np.roll(np.roll(np.roll(rho,Vsteps[0],axis=0),Vsteps[1],axis=1),Vsteps[2],axis=2)
3863            indx = np.array(ma.nonzero(rho)).T
3864            rhoXYZ = indx*steps+VP-incre
3865            Nc = len(rhoXYZ)
3866            rcube = 2000.*Vol/(ForthirdPI*Nc)
3867            rmax = math.exp(math.log(rcube)/3.)**2
3868            radius = min(drawingData['mapSize']**2,rmax)
3869            view = np.array(drawingData['viewPoint'][0])
3870            Rok = np.sum(np.inner(Amat,rhoXYZ-view).T**2,axis=1)>radius
3871        Ind = GetSelectedAtoms()
3872        VS = np.array(Page.canvas.GetSize())
3873        aspect = float(VS[0])/float(VS[1])
3874        cPos = drawingData['cameraPos']
3875        Zclip = drawingData['Zclip']*cPos/200.
3876        Q = drawingData['Quaternion']
3877        Tx,Ty,Tz = drawingData['viewPoint'][0]
3878        cx,ct,cs,ci = drawingData['atomPtrs']
3879        bondR = drawingData['bondRadius']
3880        G,g = G2lat.cell2Gmat(cell)
3881        GS = G
3882        GS[0][1] = GS[1][0] = math.sqrt(GS[0][0]*GS[1][1])
3883        GS[0][2] = GS[2][0] = math.sqrt(GS[0][0]*GS[2][2])
3884        GS[1][2] = GS[2][1] = math.sqrt(GS[1][1]*GS[2][2])
3885        ellipseProb = G2lat.criticalEllipse(drawingData['ellipseProb']/100.)
3886       
3887        SetBackground()
3888        glInitNames()
3889        glPushName(0)
3890       
3891        glMatrixMode(GL_PROJECTION)
3892        glLoadIdentity()
3893        glViewport(0,0,VS[0],VS[1])
3894        gluPerspective(20.,aspect,cPos-Zclip,cPos+Zclip)
3895        gluLookAt(0,0,cPos,0,0,0,0,1,0)
3896        SetLights()           
3897           
3898        glMatrixMode(GL_MODELVIEW)
3899        glLoadIdentity()
3900        matRot = G2mth.Q2Mat(Q)
3901        matRot = np.concatenate((np.concatenate((matRot,[[0],[0],[0]]),axis=1),[[0,0,0,1],]),axis=0)
3902        glMultMatrixf(matRot.T)
3903        glMultMatrixf(A4mat.T)
3904        glTranslate(-Tx,-Ty,-Tz)
3905        if drawingData['unitCellBox']:
3906            RenderBox()
3907        if drawingData['showABC']:
3908            x,y,z = drawingData['viewPoint'][0]
3909            RenderUnitVectors(x,y,z)
3910        Backbones = {}
3911        BackboneColor = []
3912        time0 = time.time()
3913#        glEnable(GL_BLEND)
3914#        glBlendFunc(GL_SRC_ALPHA,GL_ONE_MINUS_SRC_ALPHA)
3915        for iat,atom in enumerate(drawingData['Atoms']):
3916            x,y,z = atom[cx:cx+3]
3917            Bonds = atom[-2]
3918            Faces = atom[-1]
3919            try:
3920                atNum = generalData['AtomTypes'].index(atom[ct])
3921            except ValueError:
3922                atNum = -1
3923            CL = atom[cs+2]
3924            atColor = np.array(CL)/255.
3925            if drawingData['showRigidBodies'] and atom[ci] in rbAtmDict:
3926                bndColor = Or
3927            else:
3928                bndColor = atColor
3929            if iat in Ind and G2frame.dataDisplay.GetPageText(getSelection()) != 'Map peaks':
3930                atColor = np.array(Gr)/255.
3931#            color += [.25,]
3932            radius = 0.5
3933            if atom[cs] != '':
3934                try:
3935                    glLoadName(atom[-3])
3936                except: #problem with old files - missing code
3937                    pass                   
3938            if 'balls' in atom[cs]:
3939                vdwScale = drawingData['vdwScale']
3940                ballScale = drawingData['ballScale']
3941                if atNum < 0:
3942                    radius = 0.2
3943                elif 'H' == atom[ct]:
3944                    if drawingData['showHydrogen']:
3945                        if 'vdW' in atom[cs] and atNum >= 0:
3946                            radius = vdwScale*generalData['vdWRadii'][atNum]
3947                        else:
3948                            radius = ballScale*drawingData['sizeH']
3949                    else:
3950                        radius = 0.0
3951                else:
3952                    if 'vdW' in atom[cs]:
3953                        radius = vdwScale*generalData['vdWRadii'][atNum]
3954                    else:
3955                        radius = ballScale*generalData['BondRadii'][atNum]
3956                RenderSphere(x,y,z,radius,atColor)
3957                if 'sticks' in atom[cs]:
3958                    RenderBonds(x,y,z,Bonds,bondR,bndColor)
3959            elif 'ellipsoids' in atom[cs]:
3960                RenderBonds(x,y,z,Bonds,bondR,bndColor)
3961                if atom[cs+3] == 'A':                   
3962                    Uij = atom[cs+5:cs+11]
3963                    U = np.multiply(G2spc.Uij2U(Uij),GS)
3964                    U = np.inner(Amat,np.inner(U,Amat).T)
3965                    E,R = nl.eigh(U)
3966                    R4 = np.concatenate((np.concatenate((R,[[0],[0],[0]]),axis=1),[[0,0,0,1],]),axis=0)
3967                    E = np.sqrt(E)
3968                    if atom[ct] == 'H' and not drawingData['showHydrogen']:
3969                        pass
3970                    else:
3971                        RenderEllipsoid(x,y,z,ellipseProb,E,R4,atColor)                   
3972                else:
3973                    if atom[ct] == 'H' and not drawingData['showHydrogen']:
3974                        pass
3975                    else:
3976                        radius = ellipseProb*math.sqrt(abs(atom[cs+4]))
3977                        RenderSphere(x,y,z,radius,atColor)
3978            elif 'lines' in atom[cs]:
3979                radius = 0.1
3980                RenderLines(x,y,z,Bonds,bndColor)
3981#                RenderBonds(x,y,z,Bonds,0.05,color,6)
3982            elif atom[cs] == 'sticks':
3983                radius = 0.1
3984                RenderBonds(x,y,z,Bonds,bondR,bndColor)
3985            elif atom[cs] == 'polyhedra':
3986                RenderPolyhedra(x,y,z,Faces,atColor)
3987            elif atom[cs] == 'backbone':
3988                if atom[ct-1].split()[0] in ['C','N']:
3989                    if atom[2] not in Backbones:
3990                        Backbones[atom[2]] = []
3991                    Backbones[atom[2]].append(list(np.inner(Amat,np.array([x,y,z]))))
3992                    BackboneColor.append(list(atColor))
3993                   
3994            if atom[cs+1] == 'type':
3995                RenderLabel(x,y,z,'  '+atom[ct],radius,wxGreen,matRot)
3996            elif atom[cs+1] == 'name':
3997                RenderLabel(x,y,z,'  '+atom[ct-1],radius,wxGreen,matRot)
3998            elif atom[cs+1] == 'number':
3999                RenderLabel(x,y,z,'  '+str(iat),radius,wxGreen,matRot)
4000            elif atom[cs+1] == 'residue' and atom[ct-1] == 'CA':
4001                RenderLabel(x,y,z,'  '+atom[ct-4],radius,wxGreen,matRot)
4002            elif atom[cs+1] == '1-letter' and atom[ct-1] == 'CA':
4003                RenderLabel(x,y,z,'  '+atom[ct-3],radius,wxGreen,matRot)
4004            elif atom[cs+1] == 'chain' and atom[ct-1] == 'CA':
4005                RenderLabel(x,y,z,'  '+atom[ct-2],radius,wxGreen,matRot)
4006#        glDisable(GL_BLEND)
4007        if len(rhoXYZ):
4008            RenderMap(rho,rhoXYZ,indx,Rok)
4009        if len(mapPeaks):
4010            XYZ = mapPeaks.T[1:4].T
4011            mapBonds = FindPeaksBonds(XYZ)
4012            for ind,[mag,x,y,z,d] in enumerate(mapPeaks):
4013                if ind in Ind and pageName == 'Map peaks':
4014                    RenderMapPeak(x,y,z,Gr,1.0)
4015                else:
4016                    RenderMapPeak(x,y,z,Wt,mag/peakMax)
4017                if showBonds:
4018                    RenderLines(x,y,z,mapBonds[ind],Wt)
4019        if len(testRBObj) and pageName == 'RB Models':
4020            XYZ = G2mth.UpdateRBXYZ(Bmat,testRBObj['rbObj'],testRBObj['rbData'],testRBObj['rbType'])[0]
4021            rbBonds = FindPeaksBonds(XYZ)
4022            for ind,[x,y,z] in enumerate(XYZ):
4023                aType = testRBObj['rbAtTypes'][ind]
4024                name = '  '+aType+str(ind)
4025                color = np.array(testRBObj['AtInfo'][aType][1])
4026                RenderSphere(x,y,z,0.2,color/255.)
4027                RenderBonds(x,y,z,rbBonds[ind],0.03,Gr)
4028                RenderLabel(x,y,z,name,0.2,wxOrange,matRot)
4029        if len(mcsaModels) > 1 and pageName == 'MC/SA':             #skip the default MD entry
4030            for ind,[x,y,z] in enumerate(mcsaXYZ):
4031                aType = mcsaTypes[ind]
4032                name = '  '+aType+str(ind)
4033                color = np.array(MCSA['AtInfo'][aType][1])
4034                RenderSphere(x,y,z,0.2,color/255.)
4035                RenderBonds(x,y,z,mcsaBonds[ind],0.03,Gr)
4036                RenderLabel(x,y,z,name,0.2,wxOrange,matRot)
4037        if Backbones:
4038            for chain in Backbones:
4039                Backbone = Backbones[chain]
4040                RenderBackbone(Backbone,BackboneColor,bondR)
4041#        print time.time()-time0
4042        Page.canvas.SwapBuffers()
4043       
4044    def OnSize(event):
4045        Draw('size')
4046       
4047    def OnFocus(event):         #not needed?? Bind commented out below
4048        Draw('focus')
4049       
4050    try:
4051        plotNum = G2frame.G2plotNB.plotList.index(generalData['Name'])
4052        Page = G2frame.G2plotNB.nb.GetPage(plotNum)       
4053    except ValueError:
4054        Plot = G2frame.G2plotNB.addOgl(generalData['Name'])
4055        plotNum = G2frame.G2plotNB.plotList.index(generalData['Name'])
4056        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
4057        Page.views = False
4058        view = False
4059        altDown = False
4060    Font = Page.GetFont()
4061    Page.SetFocus()
4062    Page.Choice = None
4063    if mapData['Flip']:
4064        choice = [' save as/key:','jpeg','tiff','bmp','c: center on 1/2,1/2,1/2',
4065            'u: roll up','d: roll down','l: roll left','r: roll right']
4066    else:
4067        choice = [' save as/key:','jpeg','tiff','bmp','c: center on 1/2,1/2,1/2','n: next','p: previous']
4068    cb = wx.ComboBox(G2frame.G2plotNB.status,style=wx.CB_DROPDOWN|wx.CB_READONLY,choices=choice)
4069    cb.Bind(wx.EVT_COMBOBOX, OnKeyBox)
4070    cb.SetValue(' save as/key:')
4071    Page.canvas.Bind(wx.EVT_MOUSEWHEEL, OnMouseWheel)
4072    Page.canvas.Bind(wx.EVT_LEFT_DOWN, OnMouseDown)
4073    Page.canvas.Bind(wx.EVT_RIGHT_DOWN, OnMouseDown)
4074    Page.canvas.Bind(wx.EVT_MIDDLE_DOWN, OnMouseDown)
4075    Page.canvas.Bind(wx.EVT_KEY_UP, OnKey)
4076    Page.canvas.Bind(wx.EVT_MOTION, OnMouseMove)
4077    Page.canvas.Bind(wx.EVT_SIZE, OnSize)
4078#    Page.canvas.Bind(wx.EVT_SET_FOCUS, OnFocus)
4079    Page.camera['position'] = drawingData['cameraPos']
4080    Page.camera['viewPoint'] = np.inner(Amat,drawingData['viewPoint'][0])
4081    Page.camera['backColor'] = np.array(list(drawingData['backColor'])+[0,])/255.
4082    try:
4083        Page.canvas.SetCurrent()
4084    except:
4085        pass
4086    Draw('main')
4087       
4088################################################################################
4089#### Plot Rigid Body
4090################################################################################
4091
4092def PlotRigidBody(G2frame,rbType,AtInfo,rbData,defaults):
4093    '''RB plotting package. Can show rigid body structures as balls & sticks
4094    '''
4095
4096    def FindBonds(XYZ):
4097        rbTypes = rbData['rbTypes']
4098        Radii = []
4099        for Atype in rbTypes:
4100            Radii.append(AtInfo[Atype][0])
4101            if Atype == 'H':
4102                Radii[-1] = 0.5
4103        Radii = np.array(Radii)
4104        Bonds = [[] for i in range(len(Radii))]
4105        for i,xyz in enumerate(XYZ):
4106            Dx = XYZ-xyz
4107            dist = np.sqrt(np.sum(Dx**2,axis=1))
4108            sumR = Radii[i]+Radii
4109            IndB = ma.nonzero(ma.masked_greater(dist-0.85*sumR,0.))
4110            for j in IndB[0]:
4111                Bonds[i].append(Dx[j]*Radii[i]/sumR[j])
4112                Bonds[j].append(-Dx[j]*Radii[j]/sumR[j])
4113        return Bonds
4114                       
4115    Wt = np.array([255,255,255])
4116    Rd = np.array([255,0,0])
4117    Gr = np.array([0,255,0])
4118    Bl = np.array([0,0,255])
4119    uBox = np.array([[0,0,0],[1,0,0],[0,1,0],[0,0,1]])
4120    uEdges = np.array([[uBox[0],uBox[1]],[uBox[0],uBox[2]],[uBox[0],uBox[3]]])
4121    uColors = [Rd,Gr,Bl]
4122    if rbType == 'Vector':
4123        atNames = [str(i)+':'+Ty for i,Ty in enumerate(rbData['rbTypes'])]
4124        XYZ = np.array([[0.,0.,0.] for Ty in rbData['rbTypes']])
4125        for imag,mag in enumerate(rbData['VectMag']):
4126            XYZ += mag*rbData['rbVect'][imag]
4127        Bonds = FindBonds(XYZ)
4128    elif rbType == 'Residue':
4129#        atNames = [str(i)+':'+Ty for i,Ty in enumerate(rbData['atNames'])]
4130        atNames = rbData['atNames']
4131        XYZ = np.copy(rbData['rbXYZ'])      #don't mess with original!
4132        Seq = rbData['rbSeq']
4133        for ia,ib,ang,mv in Seq:
4134            va = XYZ[ia]-XYZ[ib]
4135            Q = G2mth.AVdeg2Q(ang,va)
4136            for im in mv:
4137                vb = XYZ[im]-XYZ[ib]
4138                vb = G2mth.prodQVQ(Q,vb)
4139                XYZ[im] = XYZ[ib]+vb
4140        Bonds = FindBonds(XYZ)
4141    elif rbType == 'Z-matrix':
4142        pass
4143
4144#    def SetRBOrigin():
4145#        page = getSelection()
4146#        if page:
4147#            if G2frame.dataDisplay.GetPageText(page) == 'Rigid bodies':
4148#                G2frame.MapPeaksTable.SetData(mapPeaks)
4149#                panel = G2frame.dataDisplay.GetPage(page).GetChildren()
4150#                names = [child.GetName() for child in panel]
4151#                panel[names.index('grid window')].Refresh()
4152           
4153    def OnMouseDown(event):
4154        xy = event.GetPosition()
4155        defaults['oldxy'] = list(xy)
4156
4157    def OnMouseMove(event):
4158        newxy = event.GetPosition()
4159                               
4160        if event.Dragging():
4161            if event.LeftIsDown():
4162                SetRotation(newxy)
4163                Q = defaults['Quaternion']
4164                G2frame.G2plotNB.status.SetStatusText('New quaternion: %.2f+, %.2fi+ ,%.2fj+, %.2fk'%(Q[0],Q[1],Q[2],Q[3]),1)
4165#            elif event.RightIsDown():
4166#                SetRBOrigin(newxy)
4167            elif event.MiddleIsDown():
4168                SetRotationZ(newxy)
4169                Q = defaults['Quaternion']
4170                G2frame.G2plotNB.status.SetStatusText('New quaternion: %.2f+, %.2fi+ ,%.2fj+, %.2fk'%(Q[0],Q[1],Q[2],Q[3]),1)
4171            Draw('move')
4172       
4173    def OnMouseWheel(event):
4174        defaults['cameraPos'] += event.GetWheelRotation()/24
4175        defaults['cameraPos'] = max(10,min(500,defaults['cameraPos']))
4176        G2frame.G2plotNB.status.SetStatusText('New camera distance: %.2f'%(defaults['cameraPos']),1)
4177        Draw('wheel')
4178       
4179    def SetBackground():
4180        R,G,B,A = Page.camera['backColor']
4181        glClearColor(R,G,B,A)
4182        glClear(GL_COLOR_BUFFER_BIT | GL_DEPTH_BUFFER_BIT)
4183       
4184    def SetLights():
4185        glEnable(GL_DEPTH_TEST)
4186        glShadeModel(GL_FLAT)
4187        glEnable(GL_LIGHTING)
4188        glEnable(GL_LIGHT0)
4189        glLightModeli(GL_LIGHT_MODEL_TWO_SIDE,0)
4190        glLightfv(GL_LIGHT0,GL_AMBIENT,[1,1,1,.8])
4191        glLightfv(GL_LIGHT0,GL_DIFFUSE,[1,1,1,1])
4192       
4193#    def SetRBOrigin(newxy):
4194##first get translation vector in screen coords.       
4195#        oldxy = defaults['oldxy']
4196#        if not len(oldxy): oldxy = list(newxy)
4197#        dxy = newxy-oldxy
4198#        defaults['oldxy'] = list(newxy)
4199#        V = np.array([dxy[0],-dxy[1],0.])/100.
4200#        Q = defaults['Quaternion']
4201#        V = G2mth.prodQVQ(G2mth.invQ(Q),V)
4202#        rbData['rbXYZ'] += V
4203#        PlotRigidBody(G2frame,rbType,AtInfo,rbData,defaults)
4204#               
4205    def SetRotation(newxy):
4206#first get rotation vector in screen coords. & angle increment       
4207        oldxy = defaults['oldxy']
4208        if not len(oldxy): oldxy = list(newxy)
4209        dxy = newxy-oldxy
4210        defaults['oldxy'] = list(newxy)
4211        V = np.array([dxy[1],dxy[0],0.])
4212        A = 0.25*np.sqrt(dxy[0]**2+dxy[1]**2)
4213# next transform vector back to xtal coordinates via inverse quaternion
4214# & make new quaternion
4215        Q = defaults['Quaternion']
4216        V = G2mth.prodQVQ(G2mth.invQ(Q),V)
4217        DQ = G2mth.AVdeg2Q(A,V)
4218        Q = G2mth.prodQQ(Q,DQ)
4219        defaults['Quaternion'] = Q
4220# finally get new view vector - last row of rotation matrix
4221        VD = G2mth.Q2Mat(Q)[2]
4222        VD /= np.sqrt(np.sum(VD**2))
4223        defaults['viewDir'] = VD
4224       
4225    def SetRotationZ(newxy):                       
4226#first get rotation vector (= view vector) in screen coords. & angle increment       
4227        View = glGetIntegerv(GL_VIEWPORT)
4228        cent = [View[2]/2,View[3]/2]
4229        oldxy = defaults['oldxy']
4230        if not len(oldxy): oldxy = list(newxy)
4231        dxy = newxy-oldxy
4232        defaults['oldxy'] = list(newxy)
4233        V = defaults['viewDir']
4234        A = [0,0]
4235        A[0] = dxy[1]*.25
4236        A[1] = dxy[0]*.25
4237        if newxy[0] > cent[0]:
4238            A[0] *= -1
4239        if newxy[1] < cent[1]:
4240            A[1] *= -1       
4241# next transform vector back to xtal coordinates & make new quaternion
4242        Q = defaults['Quaternion']
4243        Qx = G2mth.AVdeg2Q(A[0],V)
4244        Qy = G2mth.AVdeg2Q(A[1],V)
4245        Q = G2mth.prodQQ(Q,Qx)
4246        Q = G2mth.prodQQ(Q,Qy)
4247        defaults['Quaternion'] = Q
4248
4249    def RenderUnitVectors(x,y,z):
4250        xyz = np.array([x,y,z])
4251        glEnable(GL_COLOR_MATERIAL)
4252        glLineWidth(1)
4253        glPushMatrix()
4254        glTranslate(x,y,z)
4255        glBegin(GL_LINES)
4256        for line,color in zip(uEdges,uColors):
4257            glColor3ubv(color)
4258            glVertex3fv(-line[1])
4259            glVertex3fv(line[1])
4260        glEnd()
4261        glPopMatrix()
4262        glColor4ubv([0,0,0,0])
4263        glDisable(GL_COLOR_MATERIAL)
4264               
4265    def RenderSphere(x,y,z,radius,color):
4266        glMaterialfv(GL_FRONT_AND_BACK,GL_DIFFUSE,color)
4267        glPushMatrix()
4268        glTranslate(x,y,z)
4269        q = gluNewQuadric()
4270        gluSphere(q,radius,20,10)
4271        glPopMatrix()
4272       
4273    def RenderBonds(x,y,z,Bonds,radius,color,slice=20):
4274        glMaterialfv(GL_FRONT_AND_BACK,GL_DIFFUSE,color)
4275        glPushMatrix()
4276        glTranslate(x,y,z)
4277        for Dx in Bonds:
4278            glPushMatrix()
4279            Z = np.sqrt(np.sum(Dx**2))
4280            if Z:
4281                azm = atan2d(-Dx[1],-Dx[0])
4282                phi = acosd(Dx[2]/Z)
4283                glRotate(-azm,0,0,1)
4284                glRotate(phi,1,0,0)
4285                q = gluNewQuadric()
4286                gluCylinder(q,radius,radius,Z,slice,2)
4287            glPopMatrix()           
4288        glPopMatrix()
4289               
4290    def RenderLabel(x,y,z,label,matRot):       
4291        glPushMatrix()
4292        glTranslate(x,y,z)
4293        glDisable(GL_LIGHTING)
4294        glRasterPos3f(0,0,0)
4295        glMultMatrixf(matRot)
4296        glRotate(180,1,0,0)             #fix to flip about x-axis
4297        text = gltext.TextElement(text=label,font=Font,foreground=wx.WHITE)
4298        text.draw_text(scale=0.025)
4299        glEnable(GL_LIGHTING)
4300        glPopMatrix()
4301       
4302    def Draw(caller=''):
4303#useful debug?       
4304#        if caller:
4305#            print caller
4306# end of useful debug
4307        cPos = defaults['cameraPos']
4308        VS = np.array(Page.canvas.GetSize())
4309        aspect = float(VS[0])/float(VS[1])
4310        Zclip = 500.0
4311        Q = defaults['Quaternion']
4312        SetBackground()
4313        glInitNames()
4314        glPushName(0)
4315       
4316        glMatrixMode(GL_PROJECTION)
4317        glLoadIdentity()
4318        glViewport(0,0,VS[0],VS[1])
4319        gluPerspective(20.,aspect,1.,500.)
4320        gluLookAt(0,0,cPos,0,0,0,0,1,0)
4321        SetLights()           
4322           
4323        glMatrixMode(GL_MODELVIEW)
4324        glLoadIdentity()
4325        matRot = G2mth.Q2Mat(Q)
4326        matRot = np.concatenate((np.concatenate((matRot,[[0],[0],[0]]),axis=1),[[0,0,0,1],]),axis=0)
4327        glMultMatrixf(matRot.T)
4328        RenderUnitVectors(0.,0.,0.)
4329        radius = 0.2
4330        for iat,atom in enumerate(XYZ):
4331            x,y,z = atom
4332            CL = AtInfo[rbData['rbTypes'][iat]][1]
4333            color = np.array(CL)/255.
4334            RenderSphere(x,y,z,radius,color)
4335            RenderBonds(x,y,z,Bonds[iat],0.05,color)
4336            RenderLabel(x,y,z,'  '+atNames[iat],matRot)
4337        Page.canvas.SwapBuffers()
4338
4339    def OnSize(event):
4340        Draw('size')
4341       
4342    try:
4343        plotNum = G2frame.G2plotNB.plotList.index('Rigid body')
4344        Page = G2frame.G2plotNB.nb.GetPage(plotNum)       
4345    except ValueError:
4346        Plot = G2frame.G2plotNB.addOgl('Rigid body')
4347        plotNum = G2frame.G2plotNB.plotList.index('Rigid body')
4348        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
4349        Page.views = False
4350        view = False
4351        altDown = False
4352    Page.SetFocus()
4353    Font = Page.GetFont()
4354    Page.canvas.Bind(wx.EVT_MOUSEWHEEL, OnMouseWheel)
4355    Page.canvas.Bind(wx.EVT_LEFT_DOWN, OnMouseDown)
4356    Page.canvas.Bind(wx.EVT_RIGHT_DOWN, OnMouseDown)
4357    Page.canvas.Bind(wx.EVT_MIDDLE_DOWN, OnMouseDown)
4358    Page.canvas.Bind(wx.EVT_MOTION, OnMouseMove)
4359    Page.canvas.Bind(wx.EVT_SIZE, OnSize)
4360    Page.camera['position'] = defaults['cameraPos']
4361    Page.camera['backColor'] = np.array([0,0,0,0])
4362    Page.canvas.SetCurrent()
4363    Draw('main')
Note: See TracBrowser for help on using the repository browser.