source: trunk/GSASIIplot.py @ 1789

Last change on this file since 1789 was 1789, checked in by vondreele, 7 years ago

fix to POData & peak picking

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