source: trunk/GSASIIplot.py @ 1787

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

remove user reject HKL selection from Reflection List
implement rule bases user reject for HKLF reflections in Controls
user rejection changes mul to -mul; reflection list shows this in red
implement tool tip on 3Dhkl plots showing hkl indices on each point encountered

  • Property svn:eol-style set to native
  • Property svn:keywords set to Date Author Revision URL Id
File size: 224.1 KB
Line 
1# -*- coding: utf-8 -*-
2'''
3*GSASIIplot: plotting routines*
4===============================
5
6'''
7########### SVN repository information ###################
8# $Date: 2015-04-13 17:59:34 +0000 (Mon, 13 Apr 2015) $
9# $Author: vondreele $
10# $Revision: 1787 $
11# $URL: trunk/GSASIIplot.py $
12# $Id: GSASIIplot.py 1787 2015-04-13 17:59:34Z 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: 1787 $")
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            data[0]['delOffset'] = -ypos/Ymax
1244            G2frame.itemPicked = None
1245            G2frame.PatternTree.SetItemPyData(PickId,data)
1246            PlotPatterns(G2frame,plotType=plottype)
1247            return
1248        Parms,Parms2 = G2frame.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(G2frame,G2frame.PatternId, 'Instrument Parameters'))
1249        xpos = event.xdata
1250        if G2frame.PatternTree.GetItemText(PickId) in ['Peak List','Limits'] and xpos:
1251            lines = []
1252            for line in G2frame.Lines: 
1253                lines.append(line.get_xdata()[0])
1254            try:
1255                lineNo = lines.index(G2frame.itemPicked.get_xdata()[0])
1256            except ValueError:
1257                lineNo = -1
1258            if  lineNo in [0,1] or lineNo in exclLines:
1259                LimitId = G2gd.GetPatternTreeItemId(G2frame,G2frame.PatternId, 'Limits')
1260                limits = G2frame.PatternTree.GetItemPyData(LimitId)
1261                id = lineNo/2+1
1262                id2 = lineNo%2
1263                if G2frame.plotStyle['qPlot'] and 'PWDR' in plottype:
1264                    limits[id][id2] = G2lat.Dsp2pos(Parms,2.*np.pi/xpos)
1265                elif G2frame.plotStyle['dPlot'] and 'PWDR' in plottype:
1266                    limits[id][id2] = G2lat.Dsp2pos(Parms,xpos)
1267                else:
1268                    limits[id][id2] = xpos
1269                if id > 1 and limits[id][0] > limits[id][1]:
1270                        limits[id].reverse()
1271                limits[1][0] = min(max(limits[0][0],limits[1][0]),limits[1][1])
1272                limits[1][1] = max(min(limits[0][1],limits[1][1]),limits[1][0])
1273                G2frame.PatternTree.SetItemPyData(LimitId,limits)
1274                if G2frame.PatternTree.GetItemText(G2frame.PickId) == 'Limits':
1275                    G2pdG.UpdateLimitsGrid(G2frame,limits,plottype)
1276            elif lineNo > 1:
1277                PeakId = G2gd.GetPatternTreeItemId(G2frame,G2frame.PatternId, 'Peak List')
1278                peaks = G2frame.PatternTree.GetItemPyData(PeakId)
1279                if event.button == 3:
1280                    del peaks['peaks'][lineNo-2]
1281                else:
1282                    if G2frame.plotStyle['qPlot']:
1283                        peaks['peaks'][lineNo-2][0] = G2lat.Dsp2pos(Parms,2.*np.pi/xpos)
1284                    elif G2frame.plotStyle['dPlot']:
1285                        peaks['peaks'][lineNo-2][0] = G2lat.Dsp2pos(Parms,xpos)
1286                    else:
1287                        peaks['peaks'][lineNo-2][0] = xpos
1288                    peaks['sigDict'] = {}        #no longer valid
1289                G2frame.PatternTree.SetItemPyData(PeakId,peaks)
1290                G2pdG.UpdatePeakGrid(G2frame,peaks)
1291        elif G2frame.PatternTree.GetItemText(PickId) in ['Models',] and xpos:
1292            lines = []
1293            for line in G2frame.Lines: 
1294                lines.append(line.get_xdata()[0])
1295            try:
1296                lineNo = lines.index(G2frame.itemPicked.get_xdata()[0])
1297            except ValueError:
1298                lineNo = -1
1299            if  lineNo in [0,1]:
1300                LimitId = G2gd.GetPatternTreeItemId(G2frame,G2frame.PatternId, 'Limits')
1301                data = G2frame.PatternTree.GetItemPyData(LimitId)
1302                data[1][lineNo] = xpos
1303                data[1][0] = min(max(data[0][0],data[1][0]),data[1][1])
1304                data[1][1] = max(min(data[0][1],data[1][1]),data[1][0])
1305                G2frame.PatternTree.SetItemPyData(LimitId,data)       
1306        elif (G2frame.PatternTree.GetItemText(PickId) == 'Reflection Lists' or \
1307            'PWDR' in G2frame.PatternTree.GetItemText(PickId)) and xpos:
1308            Phases = G2frame.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(G2frame,PatternId,'Reflection Lists'))
1309            pick = str(G2frame.itemPicked).split('(')[1].strip(')')
1310            if 'line' not in pick:       #avoid data points, etc.
1311                data = G2frame.PatternTree.GetItemPyData(PickId)
1312                num = Phases.keys().index(pick)
1313                if num:
1314                    data[0]['refDelt'] = -(event.ydata-Pattern[0]['refOffset'])/(num*Ymax)
1315                else:       #1st row of refl ticks
1316                    data[0]['refOffset'] = event.ydata
1317                G2frame.PatternTree.SetItemPyData(PickId,data)
1318        PlotPatterns(G2frame,plotType=plottype)
1319        G2frame.itemPicked = None   
1320
1321    try:
1322        plotNum = G2frame.G2plotNB.plotList.index('Powder Patterns')
1323        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
1324        if not newPlot:
1325            Plot = Page.figure.gca()          #get previous powder plot & get limits
1326            xylim = Plot.get_xlim(),Plot.get_ylim()
1327        Page.figure.clf()
1328        Plot = Page.figure.gca()          #get a fresh plot after clf()
1329    except ValueError:
1330        if plottype == 'SASD':
1331            G2frame.logPlot = True
1332            G2frame.ErrorBars = True
1333        newPlot = True
1334        G2frame.Cmax = 1.0
1335        Plot = G2frame.G2plotNB.addMpl('Powder Patterns').gca()
1336        plotNum = G2frame.G2plotNB.plotList.index('Powder Patterns')
1337        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
1338        Page.canvas.mpl_connect('key_press_event', OnPlotKeyPress)
1339        Page.canvas.mpl_connect('motion_notify_event', OnMotion)
1340        Page.canvas.mpl_connect('pick_event', OnPick)
1341        Page.canvas.mpl_connect('button_release_event', OnRelease)
1342        Page.canvas.mpl_connect('button_press_event',OnPress)
1343    if plottype == 'PWDR':  # avoids a very nasty clash with KILL_FOCUS in SASD TextCtrl?
1344        Page.SetFocus()
1345    G2frame.G2plotNB.status.DestroyChildren()
1346    if G2frame.Contour:
1347        Page.Choice = (' key press','d: lower contour max','u: raise contour max','o: reset contour max',
1348            'i: interpolation method','s: color scheme','c: contour off')
1349    else:
1350        if G2frame.logPlot:
1351            if 'PWDR' in plottype:
1352                if G2frame.SinglePlot:
1353                    Page.Choice = (' key press','n: log(I) off',
1354                        'c: contour on','q: toggle q plot','t: toggle d-spacing plot',
1355                            'm: toggle multidata plot','w: toggle divide by sig','+: no selection')
1356                else:
1357                    Page.Choice = (' key press','n: log(I) off',
1358                        'd: offset down','l: offset left','r: offset right','u: offset up','o: reset offset',
1359                        'c: contour on','q: toggle q plot','t: toggle d-spacing plot',
1360                        'm: toggle multidata plot','w: toggle divide by sig','+: no selection')
1361            elif 'SASD' in plottype:
1362                if G2frame.SinglePlot:
1363                    Page.Choice = (' key press','b: toggle subtract background file','n: semilog on',
1364                        'q: toggle S(q) plot','m: toggle multidata plot','w: toggle (Io-Ic)/sig plot','+: no selection')
1365                else:
1366                    Page.Choice = (' key press','b: toggle subtract background file','n: semilog on',
1367                        'd: offset down','l: offset left','r: offset right','u: offset up','o: reset offset',
1368                        'q: toggle S(q) plot','m: toggle multidata plot','w: toggle (Io-Ic)/sig plot','+: no selection')
1369        else:
1370            if 'PWDR' in plottype:
1371                if G2frame.SinglePlot:
1372                    Page.Choice = (' key press',
1373                        'b: toggle subtract background','n: log(I) on','s: toggle sqrt plot','c: contour on',
1374                        'q: toggle q plot','t: toggle d-spacing plot','m: toggle multidata plot',
1375                        'w: toggle divide by sig','+: no selection')
1376                else:
1377                    Page.Choice = (' key press','l: offset left','r: offset right','d: offset down',
1378                        'u: offset up','o: reset offset','b: toggle subtract background','n: log(I) on','c: contour on',
1379                        'q: toggle q plot','t: toggle d-spacing plot','m: toggle multidata plot',
1380                        'w: toggle divide by sig','+: no selection')
1381            elif 'SASD' in plottype:
1382                if G2frame.SinglePlot:
1383                    Page.Choice = (' key press','b: toggle subtract background file','n: loglog on','e: toggle error bars',
1384                        'q: toggle S(q) plot','m: toggle multidata plot','w: toggle (Io-Ic)/sig plot','+: no selection')
1385                else:
1386                    Page.Choice = (' key press','b: toggle subtract background file','n: loglog on','e: toggle error bars',
1387                        'd: offset down','l: offset left','r: offset right','u: offset up','o: reset offset',
1388                        'q: toggle S(q) plot','m: toggle multidata plot','w: toggle (Io-Ic)/sig plot','+: no selection')
1389
1390    Page.keyPress = OnPlotKeyPress   
1391    PickId = G2frame.PickId
1392    PatternId = G2frame.PatternId
1393    colors=['b','g','r','c','m','k']
1394    Lines = []
1395    exclLines = []
1396    if G2frame.SinglePlot and PatternId:
1397        Pattern = G2frame.PatternTree.GetItemPyData(PatternId)
1398        Pattern.append(G2frame.PatternTree.GetItemText(PatternId))
1399        PlotList = [Pattern,]
1400        Parms,Parms2 = G2frame.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(G2frame,
1401            G2frame.PatternId, 'Instrument Parameters'))
1402        Sample = G2frame.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(G2frame,G2frame.PatternId, 'Sample Parameters'))
1403        ParmList = [Parms,]
1404        SampleList = [Sample,]
1405        Title = Pattern[-1]
1406    else:       
1407        Title = os.path.split(G2frame.GSASprojectfile)[1]
1408        PlotList = []
1409        ParmList = []
1410        SampleList = []
1411        item, cookie = G2frame.PatternTree.GetFirstChild(G2frame.root)
1412        while item:
1413            if plottype in G2frame.PatternTree.GetItemText(item):
1414                Pattern = G2frame.PatternTree.GetItemPyData(item)
1415                if len(Pattern) < 3:                    # put name on end if needed
1416                    Pattern.append(G2frame.PatternTree.GetItemText(item))
1417                if 'Offset' not in Pattern[0]:     #plot offset data
1418                    Pattern[0].update({'Offset':[0.0,0.0],'delOffset':0.02,'refOffset':-1.0,'refDelt':0.01,})
1419                PlotList.append(Pattern)
1420                ParmList.append(G2frame.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(G2frame,
1421                    item,'Instrument Parameters'))[0])
1422                SampleList.append(G2frame.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(G2frame,
1423                    item, 'Sample Parameters')))
1424            item, cookie = G2frame.PatternTree.GetNextChild(G2frame.root, cookie)               
1425    lenX = 0
1426    if PickId:
1427        if G2frame.PatternTree.GetItemText(PickId) in ['Reflection Lists']:
1428            Phases = G2frame.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(G2frame,PatternId,'Reflection Lists'))
1429            HKL = []
1430            if Phases: 
1431                try:
1432                    for peak in Phases[G2frame.RefList]['RefList']:
1433                        if len(peak) > 15:
1434                            HKL.append(peak[:7])    #SS reflection list - need peak[:7]
1435                        else:
1436                            HKL.append(peak[:6])
1437                except TypeError:   #old style patch
1438                    for peak in Phases[G2frame.RefList]:
1439                        HKL.append(peak[:6])                   
1440                HKL = np.array(HKL)
1441        else:
1442            HKL = np.array(G2frame.HKL)
1443    Ymax = None
1444    for Pattern in PlotList:
1445        xye = Pattern[1]
1446        if xye[1] is None: continue
1447        if Ymax is None: Ymax = max(xye[1])
1448        Ymax = max(Ymax,max(xye[1]))
1449    if Ymax is None: return # nothing to plot
1450    offsetX = Pattern[0]['Offset'][1]
1451    offsetY = Pattern[0]['Offset'][0]
1452    if G2frame.logPlot:
1453        Title = 'log('+Title+')'
1454    Plot.set_title(Title)
1455    if G2frame.plotStyle['qPlot'] or 'SASD' in plottype and not G2frame.Contour:
1456        Plot.set_xlabel(r'$Q, \AA^{-1}$',fontsize=16)
1457    elif G2frame.plotStyle['dPlot'] and 'PWDR' in plottype and not G2frame.Contour:
1458        Plot.set_xlabel(r'$d, \AA$',fontsize=16)
1459    else:
1460        if 'C' in ParmList[0]['Type'][0]:       
1461            Plot.set_xlabel(r'$\mathsf{2\theta}$',fontsize=16)
1462        else:
1463            if G2frame.Contour:
1464                Plot.set_xlabel(r'Channel no.',fontsize=16)           
1465            else:
1466                Plot.set_xlabel(r'$TOF, \mathsf{\mu}$s',fontsize=16)           
1467    if G2frame.Weight:
1468        if 'PWDR' in plottype:
1469            Plot.set_ylabel(r'$\mathsf{I/\sigma(I)}$',fontsize=16)
1470        elif 'SASD' in plottype:
1471            Plot.set_ylabel(r'$\mathsf{\Delta(I)/\sigma(I)}$',fontsize=16)
1472    else:
1473        if 'C' in ParmList[0]['Type'][0]:
1474            if 'PWDR' in plottype:
1475                if G2frame.plotStyle['sqrtPlot']:
1476                    Plot.set_ylabel(r'$\sqrt{Intensity}$',fontsize=16)
1477                else:
1478                    Plot.set_ylabel(r'$Intensity$',fontsize=16)
1479            elif 'SASD' in plottype:
1480                if G2frame.sqPlot:
1481                    Plot.set_ylabel(r'$S(Q)=I*Q^{4}$',fontsize=16)
1482                else:
1483                    Plot.set_ylabel(r'$Intensity, cm^{-1}$',fontsize=16)
1484        else:
1485            if G2frame.plotStyle['sqrtPlot']:
1486                Plot.set_ylabel(r'$\sqrt{Normalized\ intensity}$',fontsize=16)
1487            else:
1488                Plot.set_ylabel(r'$Normalized\ intensity$',fontsize=16)
1489    if G2frame.Contour:
1490        ContourZ = []
1491        ContourY = []
1492        Nseq = 0
1493    for N,Pattern in enumerate(PlotList):
1494        Parms = ParmList[N]
1495        Sample = SampleList[N]
1496        if 'C' in Parms['Type'][0]:
1497            wave = G2mth.getWave(Parms)
1498        else:
1499            difC = Parms['difC'][1]
1500        ifpicked = False
1501        LimitId = 0
1502        if Pattern[1] is None: continue # skip over uncomputed simulations
1503        xye = ma.array(ma.getdata(Pattern[1]))
1504        Zero = Parms.get('Zero',[0.,0.])[1]
1505        if PickId:
1506            ifpicked = Pattern[2] == G2frame.PatternTree.GetItemText(PatternId)
1507            LimitId = G2gd.GetPatternTreeItemId(G2frame,G2frame.PatternId,'Limits')
1508            limits = G2frame.PatternTree.GetItemPyData(LimitId)
1509            excls = limits[2:]
1510            for excl in excls:
1511                xye[0] = ma.masked_inside(xye[0],excl[0],excl[1])
1512        if G2frame.plotStyle['qPlot'] and 'PWDR' in plottype:
1513            Id = G2gd.GetPatternTreeItemId(G2frame,G2frame.root, Pattern[2])
1514            X = 2.*np.pi/G2lat.Pos2dsp(Parms,xye[0])
1515        elif G2frame.plotStyle['dPlot'] and 'PWDR' in plottype:
1516            Id = G2gd.GetPatternTreeItemId(G2frame,G2frame.root, Pattern[2])
1517            X = G2lat.Pos2dsp(Parms,xye[0])
1518        else:
1519            X = xye[0]
1520        if not lenX:
1521            lenX = len(X)
1522        if 'PWDR' in plottype:
1523            if G2frame.plotStyle['sqrtPlot']:
1524                olderr = np.seterr(invalid='ignore') #get around sqrt(-ve) error
1525                Y = np.where(xye[1]>=0.,np.sqrt(xye[1]),-np.sqrt(-xye[1]))
1526                np.seterr(invalid=olderr['invalid'])
1527            else:
1528                Y = xye[1]+offsetY*N*Ymax/100.0
1529        elif 'SASD' in plottype:
1530            B = xye[5]
1531            if G2frame.sqPlot:
1532                Y = xye[1]*Sample['Scale'][0]*(1.05)**(offsetY*N)*X**4
1533            else:
1534                Y = xye[1]*Sample['Scale'][0]*(1.05)**(offsetY*N)
1535        if LimitId and ifpicked:
1536            limits = np.array(G2frame.PatternTree.GetItemPyData(LimitId))
1537            lims = limits[1]
1538            if G2frame.plotStyle['qPlot'] and 'PWDR' in plottype:
1539                lims = 2.*np.pi/G2lat.Pos2dsp(Parms,lims)
1540            elif G2frame.plotStyle['dPlot'] and 'PWDR' in plottype:
1541                lims = G2lat.Pos2dsp(Parms,lims)
1542            Lines.append(Plot.axvline(lims[0],color='g',dashes=(5,5),picker=3.))   
1543            Lines.append(Plot.axvline(lims[1],color='r',dashes=(5,5),picker=3.))
1544            for i,item in enumerate(limits[2:]):
1545                Lines.append(Plot.axvline(item[0],color='m',dashes=(5,5),picker=3.))   
1546                Lines.append(Plot.axvline(item[1],color='m',dashes=(5,5),picker=3.))
1547                exclLines += [2*i+2,2*i+3]
1548        if G2frame.Contour:           
1549            if lenX == len(X):
1550                ContourY.append(N)
1551                ContourZ.append(Y)
1552                if 'C' in ParmList[0]['Type'][0]:       
1553                    ContourX = X
1554                else: #'T'OF
1555                    ContourX = range(lenX)
1556                Nseq += 1
1557                Plot.set_ylabel('Data sequence',fontsize=12)
1558        else:
1559            if 'SASD' in plottype and G2frame.logPlot:
1560                X *= (1.01)**(offsetX*N)
1561            else:
1562                X += offsetX*.005*N
1563            Xum = ma.getdata(X)
1564            DifLine = ['']
1565            if ifpicked:
1566                if G2frame.plotStyle['sqrtPlot']:
1567                    olderr = np.seterr(invalid='ignore') #get around sqrt(-ve) error
1568                    Z = np.where(xye[3]>=0.,np.sqrt(xye[3]),-np.sqrt(-xye[3]))
1569                    np.seterr(invalid=olderr['invalid'])
1570                else:
1571                    Z = xye[3]+offsetY*N*Ymax/100.0
1572                if 'PWDR' in plottype:
1573                    if G2frame.plotStyle['sqrtPlot']:
1574                        olderr = np.seterr(invalid='ignore') #get around sqrt(-ve) error
1575                        W = np.where(xye[4]>=0.,np.sqrt(xye[4]),-np.sqrt(-xye[4]))
1576                        np.seterr(invalid=olderr['invalid'])
1577                        D = np.where(xye[5],(Y-Z),0.)-Ymax*Pattern[0]['delOffset']
1578                    else:
1579                        W = xye[4]+offsetY*N*Ymax/100.0
1580                        D = xye[5]-Ymax*Pattern[0]['delOffset']  #powder background
1581                elif 'SASD' in plottype:
1582                    if G2frame.sqPlot:
1583                        W = xye[4]*X**4
1584                        Z = xye[3]*X**4
1585                        B = B*X**4
1586                    else:
1587                        W = xye[4]
1588                    if G2frame.SubBack:
1589                        YB = Y-B
1590                        ZB = Z
1591                    else:
1592                        YB = Y
1593                        ZB = Z+B
1594                    Plot.set_yscale("log",nonposy='mask')
1595                    if np.any(W>0.):
1596                        Plot.set_ylim(bottom=np.min(np.trim_zeros(W))/2.,top=np.max(Y)*2.)
1597                    else:
1598                        Plot.set_ylim(bottom=np.min(np.trim_zeros(YB))/2.,top=np.max(Y)*2.)
1599                if G2frame.logPlot:
1600                    if 'PWDR' in plottype:
1601                        Plot.set_yscale("log",nonposy='mask')
1602                        Plot.plot(X,Y,colors[N%6]+'+',picker=3.,clip_on=False)
1603                        Plot.plot(X,Z,colors[(N+1)%6],picker=False)
1604                        Plot.plot(X,W,colors[(N+2)%6],picker=False)     #background
1605                    elif 'SASD' in plottype:
1606                        Plot.set_xscale("log",nonposx='mask')
1607                        Ibeg = np.searchsorted(X,limits[1][0])
1608                        Ifin = np.searchsorted(X,limits[1][1])
1609                        if G2frame.Weight:
1610                            Plot.set_yscale("linear")
1611                            DS = (YB-ZB)*np.sqrt(xye[2])
1612                            Plot.plot(X[Ibeg:Ifin],DS[Ibeg:Ifin],colors[(N+3)%6],picker=False)
1613                            Plot.axhline(0.,color=wx.BLACK)
1614                            Plot.set_ylim(bottom=np.min(DS[Ibeg:Ifin])*1.2,top=np.max(DS[Ibeg:Ifin])*1.2)                                                   
1615                        else:
1616                            Plot.set_yscale("log",nonposy='mask')
1617                            if G2frame.ErrorBars:
1618                                if G2frame.sqPlot:
1619                                    Plot.errorbar(X,YB,yerr=X**4*Sample['Scale'][0]*np.sqrt(1./(Pattern[0]['wtFactor']*xye[2])),
1620                                        ecolor=colors[N%6],picker=3.,clip_on=False)
1621                                else:
1622                                    Plot.errorbar(X,YB,yerr=Sample['Scale'][0]*np.sqrt(1./(Pattern[0]['wtFactor']*xye[2])),
1623                                        ecolor=colors[N%6],picker=3.,clip_on=False)
1624                            else:
1625                                Plot.plot(X,YB,colors[N%6]+'+',picker=3.,clip_on=False)
1626                            Plot.plot(X,W,colors[(N+2)%6],picker=False)     #const. background
1627                            Plot.plot(X,ZB,colors[(N+1)%6],picker=False)
1628                elif G2frame.Weight and 'PWDR' in plottype:
1629                    DY = xye[1]*np.sqrt(xye[2])
1630                    Ymax = max(DY)
1631                    DZ = xye[3]*np.sqrt(xye[2])
1632                    DS = xye[5]*np.sqrt(xye[2])-Ymax*Pattern[0]['delOffset']
1633                    ObsLine = Plot.plot(X,DY,colors[N%6]+'+',picker=3.,clip_on=False)         #Io/sig(Io)
1634                    Plot.plot(X,DZ,colors[(N+1)%6],picker=False)                    #Ic/sig(Io)
1635                    DifLine = Plot.plot(X,DS,colors[(N+3)%6],picker=1.)                    #(Io-Ic)/sig(Io)
1636                    Plot.axhline(0.,color=wx.BLACK)
1637                else:
1638                    if G2frame.SubBack:
1639                        if 'PWDR' in plottype:
1640                            Plot.plot(Xum,Y-W,colors[N%6]+'+',picker=False,clip_on=False)  #Io-Ib
1641                            Plot.plot(X,Z-W,colors[(N+1)%6],picker=False)               #Ic-Ib
1642                        else:
1643                            Plot.plot(X,YB,colors[N%6]+'+',picker=3.,clip_on=False)
1644                            Plot.plot(X,ZB,colors[(N+1)%6],picker=False)
1645                    else:
1646                        if 'PWDR' in plottype:
1647                            ObsLine = Plot.plot(Xum,Y,colors[N%6]+'+',picker=3.,clip_on=False)    #Io
1648                            Plot.plot(X,Z,colors[(N+1)%6],picker=False)                 #Ic
1649                        else:
1650                            Plot.plot(X,YB,colors[N%6]+'+',picker=3.,clip_on=False)
1651                            Plot.plot(X,ZB,colors[(N+1)%6],picker=False)
1652                    if 'PWDR' in plottype:
1653                        Plot.plot(X,W,colors[(N+2)%6],picker=False)                 #Ib
1654                        DifLine = Plot.plot(X,D,colors[(N+3)%6],picker=1.)                 #Io-Ic
1655                    Plot.axhline(0.,color=wx.BLACK)
1656                Page.canvas.SetToolTipString('')
1657                if PickId:
1658                    if G2frame.PatternTree.GetItemText(PickId) == 'Peak List':
1659                        tip = 'On data point: Pick peak - L or R MB. On line: L-move, R-delete'
1660                        Page.canvas.SetToolTipString(tip)
1661                        data = G2frame.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(G2frame,PatternId, 'Peak List'))
1662                        for item in data['peaks']:
1663                            if G2frame.plotStyle['qPlot']:
1664                                Lines.append(Plot.axvline(2.*np.pi/G2lat.Pos2dsp(Parms,item[0]),color=colors[N%6],picker=2.))
1665                            elif G2frame.plotStyle['dPlot']:
1666                                Lines.append(Plot.axvline(G2lat.Pos2dsp(Parms,item[0]),color=colors[N%6],picker=2.))
1667                            else:
1668                                Lines.append(Plot.axvline(item[0],color=colors[N%6],picker=2.))
1669                    if G2frame.PatternTree.GetItemText(PickId) == 'Limits':
1670                        tip = 'On data point: Lower limit - L MB; Upper limit - R MB. On limit: MB down to move'
1671                        Page.canvas.SetToolTipString(tip)
1672                        data = G2frame.LimitsTable.GetData()
1673                       
1674            else:   #not picked
1675                if G2frame.logPlot:
1676                    if 'PWDR' in plottype:
1677                        Plot.semilogy(X,Y,colors[N%6],picker=False,nonposy='mask')
1678                    elif 'SASD' in plottype:
1679                        Plot.semilogy(X,Y,colors[N%6],picker=False,nonposy='mask')
1680                else:
1681                    if 'PWDR' in plottype:
1682                        Plot.plot(X,Y,colors[N%6],picker=False)
1683                    elif 'SASD' in plottype:
1684                        Plot.loglog(X,Y,colors[N%6],picker=False,nonposy='mask')
1685                        Plot.set_ylim(bottom=np.min(np.trim_zeros(Y))/2.,top=np.max(Y)*2.)
1686                           
1687                if G2frame.logPlot and 'PWDR' in plottype:
1688                    Plot.set_ylim(bottom=np.min(np.trim_zeros(Y))/2.,top=np.max(Y)*2.)
1689    if PickId and not G2frame.Contour:
1690        Parms,Parms2 = G2frame.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(G2frame,PatternId, 'Instrument Parameters'))
1691        if G2frame.PatternTree.GetItemText(PickId) in ['Index Peak List','Unit Cells List']:
1692            peaks = np.array((G2frame.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(G2frame,PatternId, 'Index Peak List'))))
1693            if not len(peaks): return # are there any peaks?
1694            for peak in peaks[0]:
1695                if peak[2]:
1696                    if G2frame.plotStyle['qPlot']:
1697                        Plot.axvline(2.*np.pi/G2lat.Pos2dsp(Parms,peak[0]),color='b')
1698                    if G2frame.plotStyle['dPlot']:
1699                        Plot.axvline(G2lat.Pos2dsp(Parms,peak[0]),color='b')
1700                    else:
1701                        Plot.axvline(peak[0],color='b')
1702            for hkl in G2frame.HKL:
1703                clr = 'r'
1704                if len(hkl) > 6 and hkl[3]:
1705                    clr = 'g'
1706                if G2frame.plotStyle['qPlot']:
1707                    Plot.axvline(2.*np.pi/G2lat.Pos2dsp(Parms,hkl[-2]),color=clr,dashes=(5,5))
1708                if G2frame.plotStyle['dPlot']:
1709                    Plot.axvline(G2lat.Pos2dsp(Parms,hkl[-2]),color=clr,dashes=(5,5))
1710                else:
1711                    Plot.axvline(hkl[-2],color=clr,dashes=(5,5))
1712        elif G2frame.PatternTree.GetItemText(PickId) in ['Reflection Lists'] or \
1713            'PWDR' in G2frame.PatternTree.GetItemText(PickId):
1714            refColors=['b','r','c','g','m','k']
1715            Phases = G2frame.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(G2frame,PatternId,'Reflection Lists'))
1716            for pId,phase in enumerate(Phases):
1717                try:    #patch for old style reflection lists
1718                    peaks = Phases[phase]['RefList']
1719                except TypeError:
1720                    peaks = Phases[phase]
1721                if not len(peaks):
1722                    continue
1723                if Phases[phase].get('Super',False):
1724                    peak = np.array([[peak[5],peak[6]] for peak in peaks])
1725                else:
1726                    peak = np.array([[peak[4],peak[5]] for peak in peaks])
1727                pos = Pattern[0]['refOffset']-pId*Ymax*Pattern[0]['refDelt']*np.ones_like(peak)
1728                if G2frame.plotStyle['qPlot']:
1729                    Plot.plot(2*np.pi/peak.T[0],pos,refColors[pId%6]+'|',mew=1,ms=8,picker=3.,label=phase)
1730                elif G2frame.plotStyle['dPlot']:
1731                    Plot.plot(peak.T[0],pos,refColors[pId%6]+'|',mew=1,ms=8,picker=3.,label=phase)
1732                else:
1733                    Plot.plot(peak.T[1],pos,refColors[pId%6]+'|',mew=1,ms=8,picker=3.,label=phase)
1734            if len(Phases):
1735                handles,legends = Plot.get_legend_handles_labels()  #got double entries in the legends for some reason
1736                if handles:
1737                    Plot.legend(handles[::2],legends[::2],title='Phases',loc='best')    #skip every other one
1738           
1739    if G2frame.Contour:
1740        acolor = mpl.cm.get_cmap(G2frame.ContourColor)
1741        Img = Plot.imshow(ContourZ,cmap=acolor,vmin=0,vmax=Ymax*G2frame.Cmax,interpolation=G2frame.Interpolate, 
1742            extent=[ContourX[0],ContourX[-1],ContourY[0],ContourY[-1]],aspect='auto',origin='lower')
1743        Page.figure.colorbar(Img)
1744    else:
1745        G2frame.Lines = Lines
1746    if not newPlot:
1747        Page.toolbar.push_current()
1748        Plot.set_xlim(xylim[0])
1749        Plot.set_ylim(xylim[1])
1750#        xylim = []
1751        Page.toolbar.push_current()
1752        Page.toolbar.draw()
1753    else:
1754        Page.canvas.draw()
1755    olderr = np.seterr(invalid='ignore') #ugh - this removes a matplotlib error for mouse clicks in log plots
1756    # and sqrt(-ve) in np.where usage               
1757#    G2frame.Pwdr = True
1758   
1759################################################################################
1760##### PlotDeltSig
1761################################################################################
1762           
1763def PlotDeltSig(G2frame,kind):
1764    'needs a doc string'
1765    try:
1766        plotNum = G2frame.G2plotNB.plotList.index('Error analysis')
1767        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
1768        Page.figure.clf()
1769        Plot = Page.figure.gca()          #get a fresh plot after clf()
1770    except ValueError:
1771        newPlot = True
1772        G2frame.Cmax = 1.0
1773        Plot = G2frame.G2plotNB.addMpl('Error analysis').gca()
1774        plotNum = G2frame.G2plotNB.plotList.index('Error analysis')
1775        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
1776    Page.Choice = None
1777    PatternId = G2frame.PatternId
1778    Pattern = G2frame.PatternTree.GetItemPyData(PatternId)
1779    Pattern.append(G2frame.PatternTree.GetItemText(PatternId))
1780    wtFactor = Pattern[0]['wtFactor']
1781    if kind == 'PWDR':
1782        limits = G2frame.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(G2frame,PatternId, 'Limits'))[1]
1783        xye = np.array(Pattern[1])
1784        xmin = np.searchsorted(xye[0],limits[0])
1785        xmax = np.searchsorted(xye[0],limits[1])
1786        DS = xye[5][xmin:xmax]*np.sqrt(wtFactor*xye[2][xmin:xmax])
1787    elif kind == 'HKLF':
1788        refl = Pattern[1]
1789        DS = []
1790        for ref in refl:
1791            if ref[6] > 0.:
1792                DS.append((ref[5]-ref[7])/ref[6])
1793    Page.SetFocus()
1794    G2frame.G2plotNB.status.DestroyChildren()
1795    DS.sort()
1796    EDS = np.zeros_like(DS)
1797    DX = np.linspace(0.,1.,num=len(DS),endpoint=True)
1798    oldErr = np.seterr(invalid='ignore')    #avoid problem at DX==0
1799    T = np.sqrt(np.log(1.0/DX**2))
1800    top = 2.515517+0.802853*T+0.010328*T**2
1801    bot = 1.0+1.432788*T+0.189269*T**2+0.001308*T**3
1802    EDS = np.where(DX>0,-(T-top/bot),(T-top/bot))
1803    low1 = np.searchsorted(EDS,-1.)
1804    hi1 = np.searchsorted(EDS,1.)
1805    slp,intcp = np.polyfit(EDS[low1:hi1],DS[low1:hi1],deg=1)
1806    frac = 100.*(hi1-low1)/len(DS)
1807    G2frame.G2plotNB.status.SetStatusText(  \
1808        'Over range -1. to 1. :'+' slope = %.3f, intercept = %.3f for %.2f%% of the fitted data'%(slp,intcp,frac),1)
1809    Plot.set_title('Normal probability for '+Pattern[-1])
1810    Plot.set_xlabel(r'expected $\mathsf{\Delta/\sigma}$',fontsize=14)
1811    Plot.set_ylabel(r'observed $\mathsf{\Delta/\sigma}$',fontsize=14)
1812    Plot.plot(EDS,DS,'r+',label='result')
1813    Plot.plot([-2,2],[-2,2],'k',dashes=(5,5),label='ideal')
1814    Plot.legend(loc='upper left')
1815    np.seterr(invalid='warn')
1816    Page.canvas.draw()
1817       
1818################################################################################
1819##### PlotISFG
1820################################################################################
1821           
1822def PlotISFG(G2frame,newPlot=False,type=''):
1823    ''' Plotting package for PDF analysis; displays I(q), S(q), F(q) and G(r) as single
1824    or multiple plots with waterfall and contour plots as options
1825    '''
1826    if not type:
1827        type = G2frame.G2plotNB.plotList[G2frame.G2plotNB.nb.GetSelection()]
1828    if type not in ['I(Q)','S(Q)','F(Q)','G(R)']:
1829        return
1830    superMinusOne = unichr(0xaf)+unichr(0xb9)
1831   
1832    def OnPlotKeyPress(event):
1833        newPlot = False
1834        if event.key == 'u':
1835            if G2frame.Contour:
1836                G2frame.Cmax = min(1.0,G2frame.Cmax*1.2)
1837            elif Pattern[0]['Offset'][0] < 100.:
1838                Pattern[0]['Offset'][0] += 1.
1839        elif event.key == 'd':
1840            if G2frame.Contour:
1841                G2frame.Cmax = max(0.0,G2frame.Cmax*0.8)
1842            elif Pattern[0]['Offset'][0] > 0.:
1843                Pattern[0]['Offset'][0] -= 1.
1844        elif event.key == 'l':
1845            Pattern[0]['Offset'][1] -= 1.
1846        elif event.key == 'r':
1847            Pattern[0]['Offset'][1] += 1.
1848        elif event.key == 'o':
1849            Pattern[0]['Offset'] = [0,0]
1850        elif event.key == 'c':
1851            newPlot = True
1852            G2frame.Contour = not G2frame.Contour
1853            if not G2frame.Contour:
1854                G2frame.SinglePlot = False
1855                Pattern[0]['Offset'] = [0.,0.]
1856        elif event.key == 's':
1857            if G2frame.Contour:
1858                choice = [m for m in mpl.cm.datad.keys() if not m.endswith("_r")]
1859                choice.sort()
1860                dlg = wx.SingleChoiceDialog(G2frame,'Select','Color scheme',choice)
1861                if dlg.ShowModal() == wx.ID_OK:
1862                    sel = dlg.GetSelection()
1863                    G2frame.ContourColor = choice[sel]
1864                else:
1865                    G2frame.ContourColor = 'Paired'
1866                dlg.Destroy()
1867            else:
1868                G2frame.SinglePlot = not G2frame.SinglePlot               
1869        elif event.key == 'i':                  #for smoothing contour plot
1870            choice = ['nearest','bilinear','bicubic','spline16','spline36','hanning',
1871               'hamming','hermite','kaiser','quadric','catrom','gaussian','bessel',
1872               'mitchell','sinc','lanczos']
1873            dlg = wx.SingleChoiceDialog(G2frame,'Select','Interpolation',choice)
1874            if dlg.ShowModal() == wx.ID_OK:
1875                sel = dlg.GetSelection()
1876                G2frame.Interpolate = choice[sel]
1877            else:
1878                G2frame.Interpolate = 'nearest'
1879            dlg.Destroy()
1880        elif event.key == 't' and not G2frame.Contour:
1881            G2frame.Legend = not G2frame.Legend
1882        PlotISFG(G2frame,newPlot=newPlot,type=type)
1883       
1884    def OnKeyBox(event):
1885        if G2frame.G2plotNB.nb.GetSelection() == G2frame.G2plotNB.plotList.index(type):
1886            event.key = cb.GetValue()[0]
1887            cb.SetValue(' key press')
1888            wx.CallAfter(OnPlotKeyPress,event)
1889        Page.canvas.SetFocus() # redirect the Focus from the button back to the plot
1890                       
1891    def OnMotion(event):
1892        xpos = event.xdata
1893        if xpos:                                        #avoid out of frame mouse position
1894            ypos = event.ydata
1895            Page.canvas.SetCursor(wx.CROSS_CURSOR)
1896            try:
1897                if G2frame.Contour:
1898                    G2frame.G2plotNB.status.SetStatusText('R =%.3fA pattern ID =%5d'%(xpos,int(ypos)),1)
1899                else:
1900                    G2frame.G2plotNB.status.SetStatusText('R =%.3fA %s =%.2f'%(xpos,type,ypos),1)                   
1901            except TypeError:
1902                G2frame.G2plotNB.status.SetStatusText('Select '+type+' pattern first',1)
1903   
1904    xylim = []
1905    try:
1906        plotNum = G2frame.G2plotNB.plotList.index(type)
1907        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
1908        if not newPlot:
1909            Plot = Page.figure.gca()          #get previous plot & get limits
1910            xylim = Plot.get_xlim(),Plot.get_ylim()
1911        Page.figure.clf()
1912        Plot = Page.figure.gca()
1913    except ValueError:
1914        newPlot = True
1915        G2frame.Cmax = 1.0
1916        Plot = G2frame.G2plotNB.addMpl(type).gca()
1917        plotNum = G2frame.G2plotNB.plotList.index(type)
1918        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
1919        Page.canvas.mpl_connect('key_press_event', OnPlotKeyPress)
1920        Page.canvas.mpl_connect('motion_notify_event', OnMotion)
1921   
1922    Page.SetFocus()
1923    G2frame.G2plotNB.status.DestroyChildren()
1924    if G2frame.Contour:
1925        Page.Choice = (' key press','d: lower contour max','u: raise contour max',
1926            'i: interpolation method','s: color scheme','c: contour off')
1927    else:
1928        Page.Choice = (' key press','l: offset left','r: offset right','d: offset down','u: offset up',
1929            'o: reset offset','t: toggle legend','c: contour on','s: toggle single plot')
1930    Page.keyPress = OnPlotKeyPress
1931    PatternId = G2frame.PatternId
1932    PickId = G2frame.PickId
1933    Plot.set_title(type)
1934    if type == 'G(R)':
1935        Plot.set_xlabel(r'$R,\AA$',fontsize=14)
1936    else:
1937        Plot.set_xlabel(r'$Q,\AA$'+superMinusOne,fontsize=14)
1938    Plot.set_ylabel(r''+type,fontsize=14)
1939    colors=['b','g','r','c','m','k']
1940    name = G2frame.PatternTree.GetItemText(PatternId)[4:]
1941    Pattern = []   
1942    if G2frame.SinglePlot:
1943        name = G2frame.PatternTree.GetItemText(PatternId)
1944        name = type+name[4:]
1945        Id = G2gd.GetPatternTreeItemId(G2frame,PatternId,name)
1946        Pattern = G2frame.PatternTree.GetItemPyData(Id)
1947        if Pattern:
1948            Pattern.append(name)
1949        PlotList = [Pattern,]
1950    else:
1951        PlotList = []
1952        item, cookie = G2frame.PatternTree.GetFirstChild(G2frame.root)
1953        while item:
1954            if 'PDF' in G2frame.PatternTree.GetItemText(item):
1955                name = type+G2frame.PatternTree.GetItemText(item)[4:]
1956                Id = G2gd.GetPatternTreeItemId(G2frame,item,name)
1957                Pattern = G2frame.PatternTree.GetItemPyData(Id)
1958                if Pattern:
1959                    Pattern.append(name)
1960                    PlotList.append(Pattern)
1961            item, cookie = G2frame.PatternTree.GetNextChild(G2frame.root, cookie)
1962    PDFdata = G2frame.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(G2frame,PatternId, 'PDF Controls'))
1963    numbDen = G2pwd.GetNumDensity(PDFdata['ElList'],PDFdata['Form Vol'])
1964    Xb = [0.,10.]
1965    Yb = [0.,-40.*np.pi*numbDen]
1966    Ymax = 1.0
1967    lenX = 0
1968    for Pattern in PlotList:
1969        try:
1970            xye = Pattern[1]
1971        except IndexError:
1972            return
1973        Ymax = max(Ymax,max(xye[1]))
1974    offset = Pattern[0]['Offset'][0]*Ymax/100.0
1975    if G2frame.Contour:
1976        ContourZ = []
1977        ContourY = []
1978        Nseq = 0
1979    for N,Pattern in enumerate(PlotList):
1980        xye = Pattern[1]
1981        if PickId:
1982            ifpicked = Pattern[2] == G2frame.PatternTree.GetItemText(PatternId)
1983        X = xye[0]
1984        if not lenX:
1985            lenX = len(X)           
1986        Y = xye[1]+offset*N
1987        if G2frame.Contour:
1988            if lenX == len(X):
1989                ContourY.append(N)
1990                ContourZ.append(Y)
1991                ContourX = X
1992                Nseq += 1
1993                Plot.set_ylabel('Data sequence',fontsize=12)
1994        else:
1995            X = xye[0]+Pattern[0]['Offset'][1]*.005*N
1996            if ifpicked:
1997                Plot.plot(X,Y,colors[N%6]+'+',picker=3.,clip_on=False)
1998                Page.canvas.SetToolTipString('')
1999            else:
2000                if G2frame.Legend:
2001                    Plot.plot(X,Y,colors[N%6],picker=False,label='Azm:'+Pattern[2].split('=')[1])
2002                else:
2003                    Plot.plot(X,Y,colors[N%6],picker=False)
2004            if type == 'G(R)':
2005                Plot.plot(Xb,Yb,color='k',dashes=(5,5))
2006            elif type == 'F(Q)':
2007                Plot.axhline(0.,color=wx.BLACK)
2008            elif type == 'S(Q)':
2009                Plot.axhline(1.,color=wx.BLACK)
2010    if G2frame.Contour:
2011        acolor = mpl.cm.get_cmap(G2frame.ContourColor)
2012        Img = Plot.imshow(ContourZ,cmap=acolor,vmin=0,vmax=Ymax*G2frame.Cmax,interpolation=G2frame.Interpolate, 
2013            extent=[ContourX[0],ContourX[-1],ContourY[0],ContourY[-1]],aspect='auto',origin='lower')
2014        Page.figure.colorbar(Img)
2015    elif G2frame.Legend:
2016        Plot.legend(loc='best')
2017    if not newPlot:
2018        Page.toolbar.push_current()
2019        Plot.set_xlim(xylim[0])
2020        Plot.set_ylim(xylim[1])
2021        xylim = []
2022        Page.toolbar.push_current()
2023        Page.toolbar.draw()
2024    else:
2025        Page.canvas.draw()
2026       
2027################################################################################
2028##### PlotCalib
2029################################################################################
2030           
2031def PlotCalib(G2frame,Inst,XY,Sigs,newPlot=False):
2032    '''plot of CW or TOF peak calibration
2033    '''
2034    def OnMotion(event):
2035        xpos = event.xdata
2036        if xpos:                                        #avoid out of frame mouse position
2037            ypos = event.ydata
2038            Page.canvas.SetCursor(wx.CROSS_CURSOR)
2039            try:
2040                G2frame.G2plotNB.status.SetStatusText('X =%9.3f %s =%9.3g'%(xpos,Title,ypos),1)                   
2041            except TypeError:
2042                G2frame.G2plotNB.status.SetStatusText('Select '+Title+' pattern first',1)
2043            found = []
2044            wid = 1
2045            view = Page.toolbar._views.forward()
2046            if view:
2047                view = view[0][:2]
2048                wid = view[1]-view[0]
2049            found = XY[np.where(np.fabs(XY.T[0]-xpos) < 0.005*wid)]
2050            if len(found):
2051                pos = found[0][1]
2052                if 'C' in Inst['Type'][0]: 
2053                    Page.canvas.SetToolTipString('position=%.4f'%(pos))
2054                else:
2055                    Page.canvas.SetToolTipString('position=%.2f'%(pos))
2056            else:
2057                Page.canvas.SetToolTipString('')
2058
2059    Title = 'Position calibration'
2060    try:
2061        plotNum = G2frame.G2plotNB.plotList.index(Title)
2062        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
2063        if not newPlot:
2064            Plot = Page.figure.gca()
2065            xylim = Plot.get_xlim(),Plot.get_ylim()
2066        Page.figure.clf()
2067        Plot = Page.figure.gca()
2068    except ValueError:
2069        newPlot = True
2070        Plot = G2frame.G2plotNB.addMpl(Title).gca()
2071        plotNum = G2frame.G2plotNB.plotList.index(Title)
2072        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
2073        Page.canvas.mpl_connect('motion_notify_event', OnMotion)
2074   
2075    Page.Choice = None
2076    Page.SetFocus()
2077    G2frame.G2plotNB.status.DestroyChildren()
2078    Plot.set_title(Title)
2079    Plot.set_xlabel(r'd-spacing',fontsize=14)
2080    if 'C' in Inst['Type'][0]:
2081        Plot.set_ylabel(r'$\mathsf{\Delta(2\theta)}$',fontsize=14)
2082    else:
2083        Plot.set_ylabel(r'$\mathsf{\Delta}T/T$',fontsize=14)
2084    for ixy,xyw in enumerate(XY):
2085        if len(xyw) > 2:
2086            X,Y,W = xyw
2087        else:
2088            X,Y = xyw
2089            W = 0.
2090        Yc = G2lat.Dsp2pos(Inst,X)
2091        if 'C' in Inst['Type'][0]:
2092            Y = Y-Yc
2093            E = Sigs[ixy]
2094            bin = W/2.
2095        else:
2096            Y = (Y-Yc)/Yc
2097            E = Sigs[ixy]/Yc
2098            bin = W/(2.*Yc)
2099        if E:
2100            Plot.errorbar(X,Y,ecolor='k',yerr=E)
2101        if ixy:
2102            Plot.plot(X,Y,'kx',picker=3)
2103        else:
2104            Plot.plot(X,Y,'kx',label='peak')
2105        if W:
2106            if ixy:
2107                Plot.plot(X,bin,'b+')
2108            else:
2109                Plot.plot(X,bin,'b+',label='bin width')
2110            Plot.plot(X,-bin,'b+')
2111        Plot.axhline(0.,color='r',linestyle='--')
2112    Plot.legend(loc='best')
2113    if not newPlot:
2114        Page.toolbar.push_current()
2115        Plot.set_xlim(xylim[0])
2116        Plot.set_ylim(xylim[1])
2117#        xylim = []
2118        Page.toolbar.push_current()
2119        Page.toolbar.draw()
2120    else:
2121        Page.canvas.draw()
2122
2123################################################################################
2124##### PlotXY
2125################################################################################
2126           
2127def PlotXY(G2frame,XY,XY2=None,labelX=None,labelY=None,newPlot=False,Title=''):
2128    '''simple plot of xy data, used for diagnostic purposes
2129    '''
2130    def OnMotion(event):
2131        xpos = event.xdata
2132        if xpos:                                        #avoid out of frame mouse position
2133            ypos = event.ydata
2134            Page.canvas.SetCursor(wx.CROSS_CURSOR)
2135            try:
2136                G2frame.G2plotNB.status.SetStatusText('X =%9.3f %s =%9.3f'%(xpos,Title,ypos),1)                   
2137            except TypeError:
2138                G2frame.G2plotNB.status.SetStatusText('Select '+Title+' pattern first',1)
2139
2140    try:
2141        plotNum = G2frame.G2plotNB.plotList.index(Title)
2142        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
2143        if not newPlot:
2144            Plot = Page.figure.gca()
2145            xylim = Plot.get_xlim(),Plot.get_ylim()
2146        Page.figure.clf()
2147        Plot = Page.figure.gca()
2148    except ValueError:
2149        newPlot = True
2150        Plot = G2frame.G2plotNB.addMpl(Title).gca()
2151        plotNum = G2frame.G2plotNB.plotList.index(Title)
2152        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
2153        Page.canvas.mpl_connect('motion_notify_event', OnMotion)
2154   
2155    Page.Choice = None
2156    Page.SetFocus()
2157    G2frame.G2plotNB.status.DestroyChildren()
2158    Plot.set_title(Title)
2159    if labelX:
2160        Plot.set_xlabel(r''+labelX,fontsize=14)
2161    else:
2162        Plot.set_xlabel(r'X',fontsize=14)
2163    if labelY:
2164        Plot.set_ylabel(r''+labelY,fontsize=14)
2165    else:
2166        Plot.set_ylabel(r'Y',fontsize=14)
2167    colors=['b','g','r','c','m','k']
2168    for ixy,xy in enumerate(XY):
2169        X,Y = xy
2170        Plot.plot(X,Y,colors[ixy%6]+'+',picker=False)
2171    if len(XY2):
2172        for ixy,xy in enumerate(XY2):
2173            X,Y = xy
2174            Plot.plot(X,Y,colors[ixy%6],picker=False)
2175    if not newPlot:
2176        Page.toolbar.push_current()
2177        Plot.set_xlim(xylim[0])
2178        Plot.set_ylim(xylim[1])
2179        xylim = []
2180        Page.toolbar.push_current()
2181        Page.toolbar.draw()
2182    else:
2183        Page.canvas.draw()
2184
2185################################################################################
2186##### PlotStrain
2187################################################################################
2188           
2189def PlotStrain(G2frame,data,newPlot=False):
2190    '''plot of strain data, used for diagnostic purposes
2191    '''
2192    def OnMotion(event):
2193        xpos = event.xdata
2194        if xpos:                                        #avoid out of frame mouse position
2195            ypos = event.ydata
2196            Page.canvas.SetCursor(wx.CROSS_CURSOR)
2197            try:
2198                G2frame.G2plotNB.status.SetStatusText('d-spacing =%9.5f Azimuth =%9.3f'%(ypos,xpos),1)                   
2199            except TypeError:
2200                G2frame.G2plotNB.status.SetStatusText('Select Strain pattern first',1)
2201
2202    try:
2203        plotNum = G2frame.G2plotNB.plotList.index('Strain')
2204        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
2205        if not newPlot:
2206            Plot = Page.figure.gca()
2207            xylim = Plot.get_xlim(),Plot.get_ylim()
2208        Page.figure.clf()
2209        Plot = Page.figure.gca()
2210    except ValueError:
2211        newPlot = True
2212        Plot = G2frame.G2plotNB.addMpl('Strain').gca()
2213        plotNum = G2frame.G2plotNB.plotList.index('Strain')
2214        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
2215        Page.canvas.mpl_connect('motion_notify_event', OnMotion)
2216   
2217    Page.Choice = None
2218    Page.SetFocus()
2219    G2frame.G2plotNB.status.DestroyChildren()
2220    Plot.set_title('Strain')
2221    Plot.set_ylabel(r'd-spacing',fontsize=14)
2222    Plot.set_xlabel(r'Azimuth',fontsize=14)
2223    colors=['b','g','r','c','m','k']
2224    for N,item in enumerate(data['d-zero']):
2225        Y,X = np.array(item['ImtaObs'])         #plot azimuth as X & d-spacing as Y
2226        Plot.plot(X,Y,colors[N%6]+'+',picker=False)
2227        Y,X = np.array(item['ImtaCalc'])
2228        Plot.plot(X,Y,colors[N%6],picker=False)
2229    if not newPlot:
2230        Page.toolbar.push_current()
2231        Plot.set_xlim(xylim[0])
2232        Plot.set_ylim(xylim[1])
2233        xylim = []
2234        Page.toolbar.push_current()
2235        Page.toolbar.draw()
2236    else:
2237        Page.canvas.draw()
2238       
2239################################################################################
2240##### PlotSASDSizeDist
2241################################################################################
2242           
2243def PlotSASDSizeDist(G2frame):
2244   
2245    def OnPageChanged(event):
2246        PlotText = G2frame.G2plotNB.nb.GetPageText(G2frame.G2plotNB.nb.GetSelection())
2247        if 'Powder' in PlotText:
2248            PlotPatterns(G2frame,plotType='SASD',newPlot=True)
2249        elif 'Size' in PlotText:
2250            PlotSASDSizeDist(G2frame)
2251   
2252    def OnMotion(event):
2253        xpos = event.xdata
2254        if xpos:                                        #avoid out of frame mouse position
2255            ypos = event.ydata
2256            Page.canvas.SetCursor(wx.CROSS_CURSOR)
2257            try:
2258                G2frame.G2plotNB.status.SetStatusText('diameter =%9.3f f(D) =%9.3g'%(xpos,ypos),1)                   
2259            except TypeError:
2260                G2frame.G2plotNB.status.SetStatusText('Select Strain pattern first',1)
2261
2262    try:
2263        plotNum = G2frame.G2plotNB.plotList.index('Size Distribution')
2264        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
2265        Page.figure.clf()
2266        Plot = Page.figure.gca()          #get a fresh plot after clf()
2267    except ValueError:
2268        newPlot = True
2269        Plot = G2frame.G2plotNB.addMpl('Size Distribution').gca()
2270        G2frame.G2plotNB.Bind(wx.aui.EVT_AUINOTEBOOK_PAGE_CHANGED,OnPageChanged)
2271        plotNum = G2frame.G2plotNB.plotList.index('Size Distribution')
2272        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
2273        Page.canvas.mpl_connect('motion_notify_event', OnMotion)
2274    Page.Choice = None
2275    Page.SetFocus()
2276    PatternId = G2frame.PatternId
2277    data = G2frame.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(G2frame,PatternId, 'Models'))
2278    Bins,Dbins,BinMag = data['Size']['Distribution']
2279    Plot.set_title('Size Distribution')
2280    Plot.set_xlabel(r'$D, \AA$',fontsize=14)
2281    Plot.set_ylabel(r'$Volume distribution f(D)$',fontsize=14)
2282    if data['Size']['logBins']:
2283        Plot.set_xscale("log",nonposy='mask')
2284        Plot.set_xlim([np.min(2.*Bins)/2.,np.max(2.*Bins)*2.])
2285    Plot.bar(2.*Bins-Dbins,BinMag,2.*Dbins,facecolor='white')       #plot diameters
2286    if 'Size Calc' in data:
2287        Rbins,Dist = data['Size Calc']
2288        for i in range(len(Rbins)):
2289            if len(Rbins[i]):
2290                Plot.plot(2.*Rbins[i],Dist[i])       #plot diameters
2291    Page.canvas.draw()
2292
2293################################################################################
2294##### PlotPowderLines
2295################################################################################
2296           
2297def PlotPowderLines(G2frame):
2298    ''' plotting of powder lines (i.e. no powder pattern) as sticks
2299    '''
2300    global HKL
2301
2302    def OnMotion(event):
2303        xpos = event.xdata
2304        if xpos:                                        #avoid out of frame mouse position
2305            Page.canvas.SetCursor(wx.CROSS_CURSOR)
2306            G2frame.G2plotNB.status.SetFields(['','2-theta =%9.3f '%(xpos,)])
2307            if G2frame.PickId and G2frame.PatternTree.GetItemText(G2frame.PickId) in ['Index Peak List','Unit Cells List']:
2308                found = []
2309                if len(HKL):
2310                    view = Page.toolbar._views.forward()[0][:2]
2311                    wid = view[1]-view[0]
2312                    found = HKL[np.where(np.fabs(HKL.T[-1]-xpos) < 0.002*wid)]
2313                if len(found):
2314                    h,k,l = found[0][:3] 
2315                    Page.canvas.SetToolTipString('%d,%d,%d'%(int(h),int(k),int(l)))
2316                else:
2317                    Page.canvas.SetToolTipString('')
2318
2319    try:
2320        plotNum = G2frame.G2plotNB.plotList.index('Powder Lines')
2321        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
2322        Page.figure.clf()
2323        Plot = Page.figure.gca()
2324    except ValueError:
2325        Plot = G2frame.G2plotNB.addMpl('Powder Lines').gca()
2326        plotNum = G2frame.G2plotNB.plotList.index('Powder Lines')
2327        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
2328        Page.canvas.mpl_connect('motion_notify_event', OnMotion)
2329       
2330    Page.Choice = None
2331    Page.SetFocus()
2332    Plot.set_title('Powder Pattern Lines')
2333    Plot.set_xlabel(r'$\mathsf{2\theta}$',fontsize=14)
2334    PickId = G2frame.PickId
2335    PatternId = G2frame.PatternId
2336    peaks = G2frame.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(G2frame,PatternId, 'Index Peak List'))[0]
2337    for peak in peaks:
2338        Plot.axvline(peak[0],color='b')
2339    HKL = np.array(G2frame.HKL)
2340    for hkl in G2frame.HKL:
2341        Plot.axvline(hkl[-1],color='r',dashes=(5,5))
2342    xmin = peaks[0][0]
2343    xmax = peaks[-1][0]
2344    delt = xmax-xmin
2345    xlim = [max(0,xmin-delt/20.),min(180.,xmax+delt/20.)]
2346    Plot.set_xlim(xlim)
2347    Page.canvas.draw()
2348    Page.toolbar.push_current()
2349
2350################################################################################
2351##### PlotPeakWidths
2352################################################################################
2353           
2354def PlotPeakWidths(G2frame):
2355    ''' Plotting of instrument broadening terms as function of 2-theta
2356    Seen when "Instrument Parameters" chosen from powder pattern data tree
2357    '''
2358#    sig = lambda Th,U,V,W: 1.17741*math.sqrt(U*tand(Th)**2+V*tand(Th)+W)*math.pi/18000.
2359#    gam = lambda Th,X,Y: (X/cosd(Th)+Y*tand(Th))*math.pi/18000.
2360    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.)
2361#    gamFW2 = lambda s,g: math.sqrt(s**2+(0.4654996*g)**2)+.5345004*g  #Ubaldo Bafile - private communication
2362    PatternId = G2frame.PatternId
2363    limitID = G2gd.GetPatternTreeItemId(G2frame,PatternId, 'Limits')
2364    if limitID:
2365        limits = G2frame.PatternTree.GetItemPyData(limitID)[:2]
2366    else:
2367        return
2368    Parms,Parms2 = G2frame.PatternTree.GetItemPyData( \
2369        G2gd.GetPatternTreeItemId(G2frame,PatternId, 'Instrument Parameters'))
2370    if 'C' in Parms['Type'][0]:
2371        lam = G2mth.getWave(Parms)
2372    else:
2373        difC = Parms['difC'][0]
2374    try:  # PATCH: deal with older peak lists, before changed to dict to implement TOF
2375        peaks = G2frame.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(G2frame,PatternId, 'Peak List'))['peaks']
2376    except TypeError:
2377        print "Your peak list needs reformatting...",
2378        item = G2gd.GetPatternTreeItemId(G2frame,PatternId, 'Peak List')
2379        G2frame.PatternTree.SelectItem(item) 
2380        item = G2gd.GetPatternTreeItemId(G2frame,PatternId, 'Instrument Parameters')
2381        G2frame.PatternTree.SelectItem(item)
2382        print "done"
2383        return
2384    try:
2385        plotNum = G2frame.G2plotNB.plotList.index('Peak Widths')
2386        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
2387        Page.figure.clf()
2388        Plot = Page.figure.gca()
2389    except ValueError:
2390        Plot = G2frame.G2plotNB.addMpl('Peak Widths').gca()
2391        plotNum = G2frame.G2plotNB.plotList.index('Peak Widths')
2392        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
2393    Page.Choice = None
2394    Page.SetFocus()
2395   
2396    Page.canvas.SetToolTipString('')
2397    colors=['b','g','r','c','m','k']
2398    X = []
2399    Y = []
2400    Z = []
2401    W = []
2402    if 'C' in Parms['Type'][0]:
2403        Plot.set_title('Instrument and sample peak widths')
2404        Plot.set_xlabel(r'$Q, \AA^{-1}$',fontsize=14)
2405        Plot.set_ylabel(r'$\Delta Q/Q, \Delta d/d$',fontsize=14)
2406        try:
2407            Xmin,Xmax = limits[1]
2408            X = np.linspace(Xmin,Xmax,num=101,endpoint=True)
2409            Q = 4.*np.pi*npsind(X/2.)/lam
2410            Z = np.ones_like(X)
2411            data = G2mth.setPeakparms(Parms,Parms2,X,Z)
2412            s = 1.17741*np.sqrt(data[4])*np.pi/18000.
2413            g = data[6]*np.pi/18000.
2414            G = G2pwd.getgamFW(g,s)
2415            Y = s/nptand(X/2.)
2416            Z = g/nptand(X/2.)
2417            W = G/nptand(X/2.)
2418            Plot.plot(Q,Y,color='r',label='Gaussian')
2419            Plot.plot(Q,Z,color='g',label='Lorentzian')
2420            Plot.plot(Q,W,color='b',label='G+L')
2421           
2422            fit = G2mth.setPeakparms(Parms,Parms2,X,Z,useFit=True)
2423            sf = 1.17741*np.sqrt(fit[4])*np.pi/18000.
2424            gf = fit[6]*np.pi/18000.
2425            Gf = G2pwd.getgamFW(gf,sf)
2426            Yf = sf/nptand(X/2.)
2427            Zf = gf/nptand(X/2.)
2428            Wf = Gf/nptand(X/2.)
2429            Plot.plot(Q,Yf,color='r',dashes=(5,5),label='Gaussian fit')
2430            Plot.plot(Q,Zf,color='g',dashes=(5,5),label='Lorentzian fit')
2431            Plot.plot(Q,Wf,color='b',dashes=(5,5),label='G+L fit')
2432           
2433            X = []
2434            Y = []
2435            Z = []
2436            W = []
2437            V = []
2438            for peak in peaks:
2439                X.append(4.0*math.pi*sind(peak[0]/2.0)/lam)
2440                try:
2441                    s = 1.17741*math.sqrt(peak[4])*math.pi/18000.
2442                except ValueError:
2443                    s = 0.01
2444                g = peak[6]*math.pi/18000.
2445                G = G2pwd.getgamFW(g,s)
2446                Y.append(s/tand(peak[0]/2.))
2447                Z.append(g/tand(peak[0]/2.))
2448                W.append(G/tand(peak[0]/2.))
2449            if len(peaks):
2450                Plot.plot(X,Y,'+',color='r',label='G peak')
2451                Plot.plot(X,Z,'+',color='g',label='L peak')
2452                Plot.plot(X,W,'+',color='b',label='G+L peak')
2453            Plot.legend(loc='best')
2454            Page.canvas.draw()
2455        except ValueError:
2456            print '**** ERROR - default U,V,W profile coefficients yield sqrt of negative value at 2theta =', \
2457                '%.3f'%(2*theta)
2458            G2frame.G2plotNB.Delete('Peak Widths')
2459    else:   #'T'OF
2460        Plot.set_title('Instrument and sample peak coefficients')
2461        Plot.set_xlabel(r'$Q, \AA^{-1}$',fontsize=14)
2462        Plot.set_ylabel(r'$\alpha, \beta, \Delta Q/Q, \Delta d/d$',fontsize=14)
2463        Xmin,Xmax = limits[1]
2464        T = np.linspace(Xmin,Xmax,num=101,endpoint=True)
2465        Z = np.ones_like(T)
2466        data = G2mth.setPeakparms(Parms,Parms2,T,Z)
2467        ds = T/difC
2468        Q = 2.*np.pi/ds
2469        A = data[4]
2470        B = data[6]
2471        S = 1.17741*np.sqrt(data[8])/T
2472        G = data[10]/T
2473        Plot.plot(Q,A,color='r',label='Alpha')
2474        Plot.plot(Q,B,color='g',label='Beta')
2475        Plot.plot(Q,S,color='b',label='Gaussian')
2476        Plot.plot(Q,G,color='m',label='Lorentzian')
2477
2478        fit = G2mth.setPeakparms(Parms,Parms2,T,Z)
2479        ds = T/difC
2480        Q = 2.*np.pi/ds
2481        Af = fit[4]
2482        Bf = fit[6]
2483        Sf = 1.17741*np.sqrt(fit[8])/T
2484        Gf = fit[10]/T
2485        Plot.plot(Q,Af,color='r',dashes=(5,5),label='Alpha fit')
2486        Plot.plot(Q,Bf,color='g',dashes=(5,5),label='Beta fit')
2487        Plot.plot(Q,Sf,color='b',dashes=(5,5),label='Gaussian fit')
2488        Plot.plot(Q,Gf,color='m',dashes=(5,5),label='Lorentzian fit')
2489       
2490        T = []
2491        A = []
2492        B = []
2493        S = []
2494        G = []
2495        W = []
2496        Q = []
2497        V = []
2498        for peak in peaks:
2499            T.append(peak[0])
2500            A.append(peak[4])
2501            B.append(peak[6])
2502            Q.append(2.*np.pi*difC/peak[0])
2503            S.append(1.17741*np.sqrt(peak[8])/peak[0])
2504            G.append(peak[10]/peak[0])
2505           
2506       
2507        Plot.plot(Q,A,'+',color='r',label='Alpha peak')
2508        Plot.plot(Q,B,'+',color='g',label='Beta peak')
2509        Plot.plot(Q,S,'+',color='b',label='Gaussian peak')
2510        Plot.plot(Q,G,'+',color='m',label='Lorentzian peak')
2511        Plot.legend(loc='best')
2512        Page.canvas.draw()
2513
2514   
2515################################################################################
2516##### PlotSizeStrainPO
2517################################################################################
2518           
2519def PlotSizeStrainPO(G2frame,data,hist='',Start=False):
2520    '''Plot 3D mustrain/size/preferred orientation figure. In this instance data is for a phase
2521    '''
2522   
2523    PatternId = G2frame.PatternId
2524    generalData = data['General']
2525    SGData = generalData['SGData']
2526    SGLaue = SGData['SGLaue']
2527    if Start:                   #initialize the spherical harmonics qlmn arrays
2528        ptx.pyqlmninit()
2529        Start = False
2530#    MuStrKeys = G2spc.MustrainNames(SGData)
2531    cell = generalData['Cell'][1:]
2532    A,B = G2lat.cell2AB(cell[:6])
2533    Vol = cell[6]
2534    useList = data['Histograms']
2535    phase = generalData['Name']
2536    plotType = generalData['Data plot type']
2537    plotDict = {'Mustrain':'Mustrain','Size':'Size','Preferred orientation':'Pref.Ori.'}
2538    for ptype in plotDict:
2539        G2frame.G2plotNB.Delete(ptype)
2540    if plotType in ['None']:
2541        return       
2542    if hist == '':
2543        hist = useList.keys()[0]
2544    numPlots = len(useList)
2545
2546    if plotType in ['Mustrain','Size']:
2547        Plot = mp3d.Axes3D(G2frame.G2plotNB.add3D(plotType))
2548    else:
2549        Plot = G2frame.G2plotNB.addMpl(plotType).gca()       
2550    plotNum = G2frame.G2plotNB.plotList.index(plotType)
2551    Page = G2frame.G2plotNB.nb.GetPage(plotNum)
2552    Page.Choice = None
2553    Page.SetFocus()
2554    G2frame.G2plotNB.status.SetStatusText('',1)
2555    if not Page.IsShown():
2556        Page.Show()
2557   
2558    PHI = np.linspace(0.,360.,30,True)
2559    PSI = np.linspace(0.,180.,30,True)
2560    X = np.outer(npsind(PHI),npsind(PSI))
2561    Y = np.outer(npcosd(PHI),npsind(PSI))
2562    Z = np.outer(np.ones(np.size(PHI)),npcosd(PSI))
2563    try:        #temp patch instead of 'mustrain' for old files with 'microstrain'
2564        coeff = useList[hist][plotDict[plotType]]
2565    except KeyError:
2566        return
2567    if plotType in ['Mustrain','Size']:
2568        if coeff[0] == 'isotropic':
2569            X *= coeff[1][0]
2570            Y *= coeff[1][0]
2571            Z *= coeff[1][0]                               
2572        elif coeff[0] == 'uniaxial':
2573           
2574            def uniaxCalc(xyz,iso,aniso,axes):
2575                Z = np.array(axes)
2576                cp = abs(np.dot(xyz,Z))
2577                sp = np.sqrt(1.-cp**2)
2578                R = iso*aniso/np.sqrt((iso*cp)**2+(aniso*sp)**2)
2579                return R*xyz
2580               
2581            iso,aniso = coeff[1][:2]
2582            axes = np.inner(A,np.array(coeff[3]))
2583            axes /= nl.norm(axes)
2584            Shkl = np.array(coeff[1])
2585            XYZ = np.dstack((X,Y,Z))
2586            XYZ = np.nan_to_num(np.apply_along_axis(uniaxCalc,2,XYZ,iso,aniso,axes))
2587            X,Y,Z = np.dsplit(XYZ,3)
2588            X = X[:,:,0]
2589            Y = Y[:,:,0]
2590            Z = Z[:,:,0]
2591       
2592        elif coeff[0] == 'ellipsoidal':
2593           
2594            def ellipseCalc(xyz,E,R):
2595                XYZ = xyz*E.T
2596                return np.inner(XYZ.T,R)
2597               
2598            S6 = coeff[4]
2599            Sij = G2lat.U6toUij(S6)
2600            E,R = nl.eigh(Sij)
2601            XYZ = np.dstack((X,Y,Z))
2602            XYZ = np.nan_to_num(np.apply_along_axis(ellipseCalc,2,XYZ,E,R))
2603            X,Y,Z = np.dsplit(XYZ,3)
2604            X = X[:,:,0]
2605            Y = Y[:,:,0]
2606            Z = Z[:,:,0]
2607           
2608        elif coeff[0] == 'generalized':
2609           
2610            def genMustrain(xyz,SGData,A,Shkl):
2611                uvw = np.inner(A.T,xyz)
2612                Strm = np.array(G2spc.MustrainCoeff(uvw,SGData))
2613                Sum = np.sum(np.multiply(Shkl,Strm))
2614                Sum = np.where(Sum > 0.01,Sum,0.01)
2615                Sum = np.sqrt(Sum)
2616                return Sum*xyz
2617               
2618            Shkl = np.array(coeff[4])
2619            if np.any(Shkl):
2620                XYZ = np.dstack((X,Y,Z))
2621                XYZ = np.nan_to_num(np.apply_along_axis(genMustrain,2,XYZ,SGData,A,Shkl))
2622                X,Y,Z = np.dsplit(XYZ,3)
2623                X = X[:,:,0]
2624                Y = Y[:,:,0]
2625                Z = Z[:,:,0]
2626                   
2627        if np.any(X) and np.any(Y) and np.any(Z):
2628            errFlags = np.seterr(all='ignore')
2629            Plot.plot_surface(X,Y,Z,rstride=1,cstride=1,color='g',linewidth=1)
2630            np.seterr(all='ignore')
2631            xyzlim = np.array([Plot.get_xlim3d(),Plot.get_ylim3d(),Plot.get_zlim3d()]).T
2632            XYZlim = [min(xyzlim[0]),max(xyzlim[1])]
2633            Plot.set_xlim3d(XYZlim)
2634            Plot.set_ylim3d(XYZlim)
2635            Plot.set_zlim3d(XYZlim)
2636            Plot.set_aspect('equal')
2637        if plotType == 'Size':
2638            Plot.set_title('Crystallite size for '+phase+'\n'+coeff[0]+' model')
2639            Plot.set_xlabel(r'X, $\mu$m')
2640            Plot.set_ylabel(r'Y, $\mu$m')
2641            Plot.set_zlabel(r'Z, $\mu$m')
2642        else:   
2643            Plot.set_title(r'$\mu$strain for '+phase+'\n'+coeff[0]+' model')
2644            Plot.set_xlabel(r'X, $\mu$strain')
2645            Plot.set_ylabel(r'Y, $\mu$strain')
2646            Plot.set_zlabel(r'Z, $\mu$strain')
2647    else:
2648        h,k,l = generalData['POhkl']
2649        if coeff[0] == 'MD':
2650            print 'March-Dollase preferred orientation plot'
2651       
2652        else:
2653            PH = np.array(generalData['POhkl'])
2654            phi,beta = G2lat.CrsAng(PH,cell[:6],SGData)
2655            SHCoef = {}
2656            for item in coeff[5]:
2657                L,N = eval(item.strip('C'))
2658                SHCoef['C%d,0,%d'%(L,N)] = coeff[5][item]                       
2659            ODFln = G2lat.Flnh(Start,SHCoef,phi,beta,SGData)
2660            X = np.linspace(0,90.0,26)
2661            Y = G2lat.polfcal(ODFln,'0',X,0.0)
2662            Plot.plot(X,Y,color='k',label=str(PH))
2663            Plot.legend(loc='best')
2664            Plot.set_title('Axial distribution for HKL='+str(PH)+' in '+phase+'\n'+hist)
2665            Plot.set_xlabel(r'$\psi$',fontsize=16)
2666            Plot.set_ylabel('MRD',fontsize=14)
2667    Page.canvas.draw()
2668   
2669################################################################################
2670##### PlotTexture
2671################################################################################
2672           
2673def PlotTexture(G2frame,data,Start=False):
2674    '''Pole figure, inverse pole figure, 3D pole distribution and 3D inverse pole distribution
2675    plotting.
2676    dict generalData contains all phase info needed which is in data
2677    '''
2678
2679    shModels = ['cylindrical','none','shear - 2/m','rolling - mmm']
2680    SamSym = dict(zip(shModels,['0','-1','2/m','mmm']))
2681    PatternId = G2frame.PatternId
2682    generalData = data['General']
2683    SGData = generalData['SGData']
2684    textureData = generalData['SH Texture']
2685    G2frame.G2plotNB.Delete('Texture')
2686    if not textureData['Order']:
2687        return                  #no plot!!
2688    SHData = generalData['SH Texture']
2689    SHCoef = SHData['SH Coeff'][1]
2690    cell = generalData['Cell'][1:7]
2691    Amat,Bmat = G2lat.cell2AB(cell)
2692    sq2 = 1.0/math.sqrt(2.0)
2693   
2694    def rp2xyz(r,p):
2695        z = npcosd(r)
2696        xy = np.sqrt(1.-z**2)
2697        return xy*npsind(p),xy*npcosd(p),z
2698           
2699    def OnMotion(event):
2700        SHData = data['General']['SH Texture']
2701        if event.xdata and event.ydata:                 #avoid out of frame errors
2702            xpos = event.xdata
2703            ypos = event.ydata
2704            if 'Inverse' in SHData['PlotType']:
2705                r = xpos**2+ypos**2
2706                if r <= 1.0:
2707                    if 'equal' in G2frame.Projection: 
2708                        r,p = 2.*npasind(np.sqrt(r)*sq2),npatan2d(ypos,xpos)
2709                    else:
2710                        r,p = 2.*npatand(np.sqrt(r)),npatan2d(ypos,xpos)
2711                    ipf = G2lat.invpolfcal(IODFln,SGData,np.array([r,]),np.array([p,]))
2712                    xyz = np.inner(Amat.T,np.array([rp2xyz(r,p)]))
2713                    y,x,z = list(xyz/np.max(np.abs(xyz)))
2714                   
2715                    G2frame.G2plotNB.status.SetFields(['',
2716                        'psi =%9.3f, beta =%9.3f, MRD =%9.3f hkl=%5.2f,%5.2f,%5.2f'%(r,p,ipf,x,y,z)])
2717                                   
2718            elif 'Axial' in SHData['PlotType']:
2719                pass
2720               
2721            else:                       #ordinary pole figure
2722                z = xpos**2+ypos**2
2723                if z <= 1.0:
2724                    z = np.sqrt(z)
2725                    if 'equal' in G2frame.Projection: 
2726                        r,p = 2.*npasind(z*sq2),npatan2d(ypos,xpos)
2727                    else:
2728                        r,p = 2.*npatand(z),npatan2d(ypos,xpos)
2729                    pf = G2lat.polfcal(ODFln,SamSym[textureData['Model']],np.array([r,]),np.array([p,]))
2730                    G2frame.G2plotNB.status.SetFields(['','phi =%9.3f, gam =%9.3f, MRD =%9.3f'%(r,p,pf)])
2731   
2732    try:
2733        plotNum = G2frame.G2plotNB.plotList.index('Texture')
2734        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
2735        Page.figure.clf()
2736        Plot = Page.figure.gca()
2737        if not Page.IsShown():
2738            Page.Show()
2739    except ValueError:
2740        Plot = G2frame.G2plotNB.addMpl('Texture').gca()
2741        plotNum = G2frame.G2plotNB.plotList.index('Texture')
2742        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
2743        Page.canvas.mpl_connect('motion_notify_event', OnMotion)
2744
2745    Page.Choice = None
2746    Page.SetFocus()
2747    G2frame.G2plotNB.status.SetFields(['',''])   
2748    PH = np.array(SHData['PFhkl'])
2749    phi,beta = G2lat.CrsAng(PH,cell,SGData)
2750    ODFln = G2lat.Flnh(Start,SHCoef,phi,beta,SGData)
2751    if not np.any(ODFln):
2752        return
2753    PX = np.array(SHData['PFxyz'])
2754    gam = atan2d(PX[0],PX[1])
2755    xy = math.sqrt(PX[0]**2+PX[1]**2)
2756    xyz = math.sqrt(PX[0]**2+PX[1]**2+PX[2]**2)
2757    psi = asind(xy/xyz)
2758    IODFln = G2lat.Glnh(Start,SHCoef,psi,gam,SamSym[textureData['Model']])
2759    if 'Axial' in SHData['PlotType']:
2760        X = np.linspace(0,90.0,26)
2761        Y = G2lat.polfcal(ODFln,SamSym[textureData['Model']],X,0.0)
2762        Plot.plot(X,Y,color='k',label=str(SHData['PFhkl']))
2763        Plot.legend(loc='best')
2764        Plot.set_title('Axial distribution for HKL='+str(SHData['PFhkl']))
2765        Plot.set_xlabel(r'$\psi$',fontsize=16)
2766        Plot.set_ylabel('MRD',fontsize=14)
2767       
2768    else:       
2769        npts = 201
2770        if 'Inverse' in SHData['PlotType']:
2771            X,Y = np.meshgrid(np.linspace(1.,-1.,npts),np.linspace(-1.,1.,npts))
2772            R,P = np.sqrt(X**2+Y**2).flatten(),npatan2d(X,Y).flatten()
2773            if 'equal' in G2frame.Projection:
2774                R = np.where(R <= 1.,2.*npasind(R*sq2),0.0)
2775            else:
2776                R = np.where(R <= 1.,2.*npatand(R),0.0)
2777            Z = np.zeros_like(R)
2778            Z = G2lat.invpolfcal(IODFln,SGData,R,P)
2779            Z = np.reshape(Z,(npts,npts))
2780            try:
2781                CS = Plot.contour(Y,X,Z,aspect='equal')
2782                Plot.clabel(CS,fontsize=9,inline=1)
2783            except ValueError:
2784                pass
2785            Img = Plot.imshow(Z.T,aspect='equal',cmap=G2frame.ContourColor,extent=[-1,1,-1,1])
2786            Page.figure.colorbar(Img)
2787            Plot.set_title('Inverse pole figure for XYZ='+str(SHData['PFxyz']))
2788            Plot.set_xlabel(G2frame.Projection.capitalize()+' projection')
2789                       
2790        else:
2791            X,Y = np.meshgrid(np.linspace(1.,-1.,npts),np.linspace(-1.,1.,npts))
2792            R,P = np.sqrt(X**2+Y**2).flatten(),npatan2d(X,Y).flatten()
2793            if 'equal' in G2frame.Projection:
2794                R = np.where(R <= 1.,2.*npasind(R*sq2),0.0)
2795            else:
2796                R = np.where(R <= 1.,2.*npatand(R),0.0)
2797            Z = np.zeros_like(R)
2798            Z = G2lat.polfcal(ODFln,SamSym[textureData['Model']],R,P)
2799            Z = np.reshape(Z,(npts,npts))
2800            try:
2801                CS = Plot.contour(Y,X,Z,aspect='equal')
2802                Plot.clabel(CS,fontsize=9,inline=1)
2803            except ValueError:
2804                pass
2805            Img = Plot.imshow(Z.T,aspect='equal',cmap=G2frame.ContourColor,extent=[-1,1,-1,1])
2806            Page.figure.colorbar(Img)
2807            Plot.set_title('Pole figure for HKL='+str(SHData['PFhkl']))
2808            Plot.set_xlabel(G2frame.Projection.capitalize()+' projection')
2809    Page.canvas.draw()
2810
2811################################################################################
2812##### Plot Modulation
2813################################################################################
2814
2815def ModulationPlot(G2frame,data,atom,ax,off=0):
2816    global Off,Atom,Ax
2817    Off = off
2818    Atom = atom
2819    Ax = ax
2820    def OnMotion(event):
2821        xpos = event.xdata
2822        if xpos:                                        #avoid out of frame mouse position
2823            ypos = event.ydata
2824            Page.canvas.SetCursor(wx.CROSS_CURSOR)
2825            try:
2826                G2frame.G2plotNB.status.SetStatusText('t =%9.3f %s =%9.3f'%(xpos,GkDelta+Ax,ypos),1)                   
2827            except TypeError:
2828                G2frame.G2plotNB.status.SetStatusText('Select '+Title+' pattern first',1)
2829   
2830    def OnPlotKeyPress(event):
2831        global Off,Atom,Ax
2832        newPlot = False
2833        if event.key == '0':
2834            Off = 0
2835        elif event.key in ['+','=']:
2836            Off += 1
2837        elif event.key == '-':
2838            Off -= 1
2839        elif event.key in ['l','r',] and mapData['Flip']:
2840            roll = 1
2841            if  event.key == 'l':
2842                roll = -1
2843            rho = Map['rho']
2844            Map['rho'] = np.roll(rho,roll,axis=3)
2845        wx.CallAfter(ModulationPlot,G2frame,data,Atom,Ax,Off)
2846
2847    try:
2848        plotNum = G2frame.G2plotNB.plotList.index('Modulation')
2849        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
2850        Page.figure.clf()
2851        Plot = Page.figure.gca()
2852        if not Page.IsShown():
2853            Page.Show()
2854    except ValueError:
2855        Plot = G2frame.G2plotNB.addMpl('Modulation').gca()
2856        plotNum = G2frame.G2plotNB.plotList.index('Modulation')
2857        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
2858        Page.canvas.mpl_connect('motion_notify_event', OnMotion)
2859        Page.canvas.mpl_connect('key_press_event', OnPlotKeyPress)
2860   
2861    Page.SetFocus()
2862    General = data['General']
2863    cx,ct,cs,cia = General['AtomPtrs']
2864    mapData = General['Map']
2865    if mapData['Flip']:
2866        Page.Choice = ['+: shift up','-: shift down','0: reset shift','l: move left','r: move right']
2867    else:
2868        Page.Choice = ['+: shift up','-: shift down','0: reset shift']
2869    Page.keyPress = OnPlotKeyPress
2870    Map = General['4DmapData']
2871    MapType = Map['MapType']
2872    rhoSize = np.array(Map['rho'].shape)
2873    atxyz = np.array(atom[cx:cx+3])
2874    waveType = atom[-1]['SS1']['waveType']
2875    Spos = atom[-1]['SS1']['Spos']
2876    tau = np.linspace(0.,2.,101)
2877    wave = np.zeros((3,101))
2878    if len(Spos):
2879        scof = []
2880        ccof = []
2881        for i,spos in enumerate(Spos):
2882            if waveType in ['Sawtooth','ZigZag'] and not i:
2883                Toff = spos[0][0]
2884                slopes = np.array(spos[0][1:])
2885                if waveType == 'Sawtooth':
2886                    wave = G2mth.posSawtooth(tau,Toff,slopes)
2887                elif waveType == 'ZigZag':
2888                    wave = G2mth.posZigZag(tau,Toff,slopes)
2889            else:
2890                scof.append(spos[0][:3])
2891                ccof.append(spos[0][3:])
2892        wave += G2mth.posFourier(tau,np.array(scof),np.array(ccof))
2893    if mapData['Flip']:
2894        Title = 'Charge flip'
2895    else:
2896        Title = MapType
2897    Title += ' map for atom '+atom[0]+    \
2898        ' at %.4f %.4f %.4f'%(atxyz[0],atxyz[1],atxyz[2])
2899    ix = -np.array(np.rint(rhoSize[:3]*atxyz),dtype='i')
2900    ix += (rhoSize[:3]/2)
2901    ix = ix%rhoSize[:3]
2902    rho = np.roll(np.roll(np.roll(Map['rho'],ix[0],axis=0),ix[1],axis=1),ix[2],axis=2)
2903    ix = rhoSize[:3]/2
2904    ib = 4
2905    if Ax == 'x':
2906        slab = np.sum(np.sum(rho[:,ix[1]-ib:ix[1]+ib,ix[2]-ib:ix[2]+ib,:],axis=2),axis=1)
2907        Plot.plot(tau,wave[0])
2908    elif Ax == 'y':
2909        slab = np.sum(np.sum(rho[ix[0]-ib:ix[0]+ib,:,ix[2]-ib:ix[2]+ib,:],axis=2),axis=0)
2910        Plot.plot(tau,wave[1])
2911    elif Ax == 'z':
2912        slab = np.sum(np.sum(rho[ix[0]-ib:ix[0]+ib,ix[1]-ib:ix[1]+ib,:,:],axis=1),axis=0)
2913        Plot.plot(tau,wave[2])
2914    Plot.set_title(Title)
2915    Plot.set_xlabel('t')
2916    Plot.set_ylabel(r'$\mathsf{\Delta}$%s'%(Ax))
2917    Slab = np.hstack((slab,slab,slab))
2918    Plot.contour(Slab[:,:21],20,extent=(0.,2.,-.5+Off*.005,.5+Off*.005))
2919    Page.canvas.draw()
2920   
2921################################################################################
2922##### PlotCovariance
2923################################################################################
2924           
2925def PlotCovariance(G2frame,Data):
2926    'needs a doc string'
2927    if not Data:
2928        print 'No covariance matrix available'
2929        return
2930    varyList = Data['varyList']
2931    values = Data['variables']
2932    Xmax = len(varyList)
2933    covMatrix = Data['covMatrix']
2934    sig = np.sqrt(np.diag(covMatrix))
2935    xvar = np.outer(sig,np.ones_like(sig))
2936    covArray = np.divide(np.divide(covMatrix,xvar),xvar.T)
2937    title = ' for\n'+Data['title']
2938    newAtomDict = Data.get('newAtomDict',{})
2939    G2frame.G2plotNB.Delete('Covariance')
2940   
2941
2942    def OnPlotKeyPress(event):
2943        newPlot = False
2944        if event.key == 's':
2945            choice = [m for m in mpl.cm.datad.keys() if not m.endswith("_r")]
2946            choice.sort()
2947            dlg = wx.SingleChoiceDialog(G2frame,'Select','Color scheme',choice)
2948            if dlg.ShowModal() == wx.ID_OK:
2949                sel = dlg.GetSelection()
2950                G2frame.VcovColor = choice[sel]
2951            else:
2952                G2frame.VcovColor = 'RdYlGn'
2953            dlg.Destroy()
2954        PlotCovariance(G2frame,Data)
2955
2956    def OnMotion(event):
2957        if event.button:
2958            ytics = imgAx.get_yticks()
2959            ytics = np.where(ytics<len(varyList),ytics,-1)
2960            ylabs = [np.where(0<=i ,varyList[int(i)],' ') for i in ytics]
2961            imgAx.set_yticklabels(ylabs)           
2962        if event.xdata and event.ydata:                 #avoid out of frame errors
2963            xpos = int(event.xdata+.5)
2964            ypos = int(event.ydata+.5)
2965            if -1 < xpos < len(varyList) and -1 < ypos < len(varyList):
2966                if xpos == ypos:
2967                    value = values[xpos]
2968                    name = varyList[xpos]
2969                    if varyList[xpos] in newAtomDict:
2970                        name,value = newAtomDict[name]                       
2971                    msg = '%s value = %.4g, esd = %.4g'%(name,value,sig[xpos])
2972                else:
2973                    msg = '%s - %s: %5.3f'%(varyList[xpos],varyList[ypos],covArray[xpos][ypos])
2974                Page.canvas.SetToolTipString(msg)
2975                G2frame.G2plotNB.status.SetFields(['',msg])
2976               
2977    try:
2978        plotNum = G2frame.G2plotNB.plotList.index('Covariance')
2979        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
2980        Page.figure.clf()
2981        Plot = Page.figure.gca()
2982        if not Page.IsShown():
2983            Page.Show()
2984    except ValueError:
2985        Plot = G2frame.G2plotNB.addMpl('Covariance').gca()
2986        plotNum = G2frame.G2plotNB.plotList.index('Covariance')
2987        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
2988        Page.canvas.mpl_connect('motion_notify_event', OnMotion)
2989        Page.canvas.mpl_connect('key_press_event', OnPlotKeyPress)
2990    Page.Choice = ['s: to change colors']
2991    Page.keyPress = OnPlotKeyPress
2992    Page.SetFocus()
2993    G2frame.G2plotNB.status.SetFields(['',''])   
2994    acolor = mpl.cm.get_cmap(G2frame.VcovColor)
2995    Img = Plot.imshow(covArray,aspect='equal',cmap=acolor,interpolation='nearest',origin='lower',
2996        vmin=-1.,vmax=1.)
2997    imgAx = Img.get_axes()
2998    ytics = imgAx.get_yticks()
2999    ylabs = [varyList[int(i)] for i in ytics[:-1]]
3000    imgAx.set_yticklabels(ylabs)
3001    colorBar = Page.figure.colorbar(Img)
3002    Plot.set_title('V-Cov matrix'+title)
3003    Plot.set_xlabel('Variable number')
3004    Plot.set_ylabel('Variable name')
3005    Page.canvas.draw()
3006   
3007################################################################################
3008##### PlotTorsion
3009################################################################################
3010
3011def PlotTorsion(G2frame,phaseName,Torsion,TorName,Names=[],Angles=[],Coeff=[]):
3012    'needs a doc string'
3013   
3014    global names
3015    names = Names
3016    sum = np.sum(Torsion)
3017    torsion = np.log(2*Torsion+1.)/sum
3018    tMin = np.min(torsion)
3019    tMax = np.max(torsion)
3020    torsion = 3.*(torsion-tMin)/(tMax-tMin)
3021    X = np.linspace(0.,360.,num=45)
3022   
3023    def OnPick(event):
3024        ind = event.ind[0]
3025        msg = 'atoms:'+names[ind]
3026        Page.canvas.SetToolTipString(msg)
3027        try:
3028            page = G2frame.dataDisplay.GetSelection()
3029        except:
3030            return
3031        if G2frame.dataDisplay.GetPageText(page) == 'Torsion restraints':
3032            torGrid = G2frame.dataDisplay.GetPage(page).Torsions
3033            torGrid.ClearSelection()
3034            for row in range(torGrid.GetNumberRows()):
3035                if names[ind] in torGrid.GetCellValue(row,0):
3036                    torGrid.SelectRow(row)
3037            torGrid.ForceRefresh()
3038               
3039    def OnMotion(event):
3040        if event.xdata and event.ydata:                 #avoid out of frame errors
3041            xpos = event.xdata
3042            ypos = event.ydata
3043            msg = 'torsion,energy: %5.3f %5.3f'%(xpos,ypos)
3044            Page.canvas.SetToolTipString(msg)
3045
3046    try:
3047        plotNum = G2frame.G2plotNB.plotList.index('Torsion')
3048        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
3049        Page.figure.clf()
3050        Plot = Page.figure.gca()
3051        if not Page.IsShown():
3052            Page.Show()
3053    except ValueError:
3054        Plot = G2frame.G2plotNB.addMpl('Torsion').gca()
3055        plotNum = G2frame.G2plotNB.plotList.index('Torsion')
3056        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
3057        Page.canvas.mpl_connect('pick_event', OnPick)
3058        Page.canvas.mpl_connect('motion_notify_event', OnMotion)
3059   
3060    Page.SetFocus()
3061    G2frame.G2plotNB.status.SetFields(['','Use mouse LB to identify torsion atoms'])
3062    Plot.plot(X,torsion,'b+')
3063    if len(Coeff):
3064        X2 = np.linspace(0.,360.,45)
3065        Y2 = np.array([-G2mth.calcTorsionEnergy(x,Coeff)[1] for x in X2])
3066        Plot.plot(X2,Y2,'r')
3067    if len(Angles):
3068        Eval = np.array([-G2mth.calcTorsionEnergy(x,Coeff)[1] for x in Angles])
3069        Plot.plot(Angles,Eval,'ro',picker=5)
3070    Plot.set_xlim((0.,360.))
3071    Plot.set_title('Torsion angles for '+TorName+' in '+phaseName)
3072    Plot.set_xlabel('angle',fontsize=16)
3073    Plot.set_ylabel('Energy',fontsize=16)
3074    Page.canvas.draw()
3075   
3076################################################################################
3077##### PlotRama
3078################################################################################
3079
3080def PlotRama(G2frame,phaseName,Rama,RamaName,Names=[],PhiPsi=[],Coeff=[]):
3081    'needs a doc string'
3082
3083    global names
3084    names = Names
3085    rama = np.log(2*Rama+1.)
3086    ramaMax = np.max(rama)
3087    rama = np.reshape(rama,(45,45))
3088    global Phi,Psi
3089    Phi = []
3090    Psi = []
3091
3092    def OnPlotKeyPress(event):
3093        newPlot = False
3094        if event.key == 's':
3095            choice = [m for m in mpl.cm.datad.keys() if not m.endswith("_r")]
3096            choice.sort()
3097            dlg = wx.SingleChoiceDialog(G2frame,'Select','Color scheme',choice)
3098            if dlg.ShowModal() == wx.ID_OK:
3099                sel = dlg.GetSelection()
3100                G2frame.RamaColor = choice[sel]
3101            else:
3102                G2frame.RamaColor = 'RdYlGn'
3103            dlg.Destroy()
3104        PlotRama(G2frame,phaseName,Rama)
3105
3106    def OnPick(event):
3107        ind = event.ind[0]
3108        msg = 'atoms:'+names[ind]
3109        Page.canvas.SetToolTipString(msg)
3110        try:
3111            page = G2frame.dataDisplay.GetSelection()
3112        except:
3113            return
3114        if G2frame.dataDisplay.GetPageText(page) == 'Ramachandran restraints':
3115            ramaGrid = G2frame.dataDisplay.GetPage(page).Ramas
3116            ramaGrid.ClearSelection()
3117            for row in range(ramaGrid.GetNumberRows()):
3118                if names[ind] in ramaGrid.GetCellValue(row,0):
3119                    ramaGrid.SelectRow(row)
3120            ramaGrid.ForceRefresh()
3121
3122    def OnMotion(event):
3123        if event.xdata and event.ydata:                 #avoid out of frame errors
3124            xpos = event.xdata
3125            ypos = event.ydata
3126            msg = 'phi/psi: %5.3f %5.3f'%(xpos,ypos)
3127            Page.canvas.SetToolTipString(msg)
3128           
3129    try:
3130        plotNum = G2frame.G2plotNB.plotList.index('Ramachandran')
3131        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
3132        Page.figure.clf()
3133        Plot = Page.figure.gca()
3134        if not Page.IsShown():
3135            Page.Show()
3136    except ValueError:
3137        Plot = G2frame.G2plotNB.addMpl('Ramachandran').gca()
3138        plotNum = G2frame.G2plotNB.plotList.index('Ramachandran')
3139        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
3140        Page.canvas.mpl_connect('pick_event', OnPick)
3141        Page.canvas.mpl_connect('motion_notify_event', OnMotion)
3142        Page.canvas.mpl_connect('key_press_event', OnPlotKeyPress)
3143
3144    Page.Choice = ['s: to change colors']
3145    Page.keyPress = OnPlotKeyPress
3146    Page.SetFocus()
3147    G2frame.G2plotNB.status.SetFields(['','Use mouse LB to identify phi/psi atoms'])
3148    acolor = mpl.cm.get_cmap(G2frame.RamaColor)
3149    if RamaName == 'All' or '-1' in RamaName:
3150        if len(Coeff): 
3151            X,Y = np.meshgrid(np.linspace(-180.,180.,45),np.linspace(-180.,180.,45))
3152            Z = np.array([-G2mth.calcRamaEnergy(x,y,Coeff)[1] for x,y in zip(X.flatten(),Y.flatten())])
3153            Plot.contour(X,Y,np.reshape(Z,(45,45)))
3154        Img = Plot.imshow(rama,aspect='equal',cmap=acolor,interpolation='nearest',
3155            extent=[-180,180,-180,180],origin='lower')
3156        if len(PhiPsi):
3157            Phi,Psi = PhiPsi.T
3158            Phi = np.where(Phi>180.,Phi-360.,Phi)
3159            Psi = np.where(Psi>180.,Psi-360.,Psi)
3160            Plot.plot(Phi,Psi,'ro',picker=5)
3161        Plot.set_xlim((-180.,180.))
3162        Plot.set_ylim((-180.,180.))
3163    else:
3164        if len(Coeff): 
3165            X,Y = np.meshgrid(np.linspace(0.,360.,45),np.linspace(0.,360.,45))
3166            Z = np.array([-G2mth.calcRamaEnergy(x,y,Coeff)[1] for x,y in zip(X.flatten(),Y.flatten())])
3167            Plot.contour(X,Y,np.reshape(Z,(45,45)))
3168        Img = Plot.imshow(rama,aspect='equal',cmap=acolor,interpolation='nearest',
3169            extent=[0,360,0,360],origin='lower')
3170        if len(PhiPsi):
3171            Phi,Psi = PhiPsi.T
3172            Plot.plot(Phi,Psi,'ro',picker=5)
3173        Plot.set_xlim((0.,360.))
3174        Plot.set_ylim((0.,360.))
3175    Plot.set_title('Ramachandran for '+RamaName+' in '+phaseName)
3176    Plot.set_xlabel(r'$\phi$',fontsize=16)
3177    Plot.set_ylabel(r'$\psi$',fontsize=16)
3178    colorBar = Page.figure.colorbar(Img)
3179    Page.canvas.draw()
3180
3181
3182################################################################################
3183##### PlotSeq
3184################################################################################
3185def PlotSelectedSequence(G2frame,ColumnList,TableGet,SelectX,fitnum=None,fitvals=None):
3186    '''Plot a result from a sequential refinement
3187
3188    :param wx.Frame G2frame: The main GSAS-II tree "window"
3189    :param list ColumnList: list of int values corresponding to columns
3190      selected as y values
3191    :param function TableGet: a function that takes a column number
3192      as argument and returns the column label, the values and there ESDs (or None)
3193    :param function SelectX: a function that returns a selected column
3194      number (or None) as the X-axis selection
3195    '''
3196    global Title,xLabel,yLabel
3197    xLabel = yLabel = Title = ''
3198    def OnMotion(event):
3199        if event.xdata and event.ydata:                 #avoid out of frame errors
3200            xpos = event.xdata
3201            ypos = event.ydata
3202            msg = '%5.3f %.6g'%(xpos,ypos)
3203            Page.canvas.SetToolTipString(msg)
3204
3205    def OnKeyPress(event):
3206        global Title,xLabel,yLabel
3207        if event.key == 's':
3208            G2frame.seqXaxis = G2frame.seqXselect()
3209            Draw()
3210        elif event.key == 't':
3211            dlg = G2gd.MultiStringDialog(G2frame,'Set titles & labels',[' Title ',' x-Label ',' y-Label '],
3212                [Title,xLabel,yLabel])
3213            if dlg.Show():
3214                Title,xLabel,yLabel = dlg.GetValues()
3215            dlg.Destroy()
3216            Draw()
3217           
3218    def Draw():
3219        global Title,xLabel,yLabel
3220        Page.SetFocus()
3221        G2frame.G2plotNB.status.SetStatusText('press s to select X axis, t to change titles',1)
3222        Plot.clear()
3223        if G2frame.seqXaxis is not None:   
3224            xName,X,Xsig = Page.seqTableGet(G2frame.seqXaxis)
3225        else:
3226            X = np.arange(0,G2frame.SeqTable.GetNumberRows(),1)
3227            xName = 'Data sequence number'
3228        for col in Page.seqYaxisList:
3229            name,Y,sig = Page.seqTableGet(col)
3230            # deal with missing (None) values
3231            Xnew = []
3232            Ynew = []
3233            Ysnew = []
3234            for i in range(len(X)):
3235                if X[i] is None or Y[i] is None: continue
3236                Xnew.append(X[i])
3237                Ynew.append(Y[i])
3238                if sig: Ysnew.append(sig[i])
3239            if Ysnew:
3240                if G2frame.seqReverse and not G2frame.seqXaxis:
3241                    Ynew = Ynew[::-1]
3242                    Ysnew = Ysnew[::-1]
3243                Plot.errorbar(Xnew,Ynew,yerr=Ysnew,label=name)
3244            else:
3245                if G2frame.seqReverse and not G2frame.seqXaxis:
3246                    Ynew = Ynew[::-1]
3247                Plot.plot(Xnew,Ynew)
3248                Plot.plot(Xnew,Ynew,'o',label=name)
3249        if Page.fitvals: # TODO: deal with fitting of None values
3250            if G2frame.seqReverse and not G2frame.seqXaxis:
3251                Page.fitvals = Page.fitvals[::-1]
3252            Plot.plot(X,Page.fitvals,label='Fit')
3253           
3254        Plot.legend(loc='best')
3255        print Title,xLabel,yLabel
3256        if Title:
3257            Plot.set_title(Title)
3258        else:
3259            Plot.set_title('')
3260        if xLabel:
3261            Plot.set_xlabel(xLabel)
3262        else:
3263            Plot.set_xlabel(xName)
3264        if yLabel:
3265            Plot.set_ylabel(yLabel)
3266        else:
3267            Plot.set_ylabel('Parameter values')
3268        Page.canvas.draw()           
3269           
3270    G2frame.seqXselect = SelectX
3271    try:
3272        G2frame.seqXaxis
3273    except:
3274        G2frame.seqXaxis = None
3275
3276    if fitnum is None:
3277        label = 'Sequential refinement'
3278    else:
3279        label = 'Parametric fit #'+str(fitnum+1)
3280    try:
3281        plotNum = G2frame.G2plotNB.plotList.index(label)
3282        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
3283        Page.figure.clf()
3284        Plot = Page.figure.gca()
3285        if not Page.IsShown():
3286            Page.Show()
3287    except ValueError:
3288        Plot = G2frame.G2plotNB.addMpl(label).gca()
3289        plotNum = G2frame.G2plotNB.plotList.index(label)
3290        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
3291        Page.canvas.mpl_connect('key_press_event', OnKeyPress)
3292        Page.canvas.mpl_connect('motion_notify_event', OnMotion)
3293    Page.Choice = ['s - select x-axis','t - change titles',]
3294    Page.keyPress = OnKeyPress
3295    Page.seqYaxisList = ColumnList
3296    Page.seqTableGet = TableGet
3297    Page.fitvals = fitvals
3298       
3299    Draw()
3300    G2frame.G2plotNB.nb.SetSelection(plotNum) # raises plot tab
3301               
3302################################################################################
3303##### PlotExposedImage & PlotImage
3304################################################################################
3305           
3306def PlotExposedImage(G2frame,newPlot=False,event=None):
3307    '''General access module for 2D image plotting
3308    '''
3309    plotNo = G2frame.G2plotNB.nb.GetSelection()
3310    if G2frame.G2plotNB.nb.GetPageText(plotNo) == '2D Powder Image':
3311        PlotImage(G2frame,newPlot,event,newImage=True)
3312    elif G2frame.G2plotNB.nb.GetPageText(plotNo) == '2D Integration':
3313        PlotIntegration(G2frame,newPlot,event)
3314
3315def OnStartMask(G2frame):
3316    '''Initiate the start of a Frame or Polygon map
3317
3318    :param wx.Frame G2frame: The main GSAS-II tree "window"
3319    :param str eventkey: a single letter ('f' or 'p') that
3320      determines what type of mask is created.   
3321    '''
3322    Masks = G2frame.PatternTree.GetItemPyData(
3323        G2gd.GetPatternTreeItemId(G2frame,G2frame.Image, 'Masks'))
3324    if G2frame.MaskKey == 'f':
3325        Masks['Frames'] = []
3326    elif G2frame.MaskKey == 'p':
3327        Masks['Polygons'].append([])
3328    elif G2frame.MaskKey == 's':
3329        Masks['Points'].append([])
3330    elif G2frame.MaskKey == 'a':
3331        Masks['Arcs'].append([])
3332    elif G2frame.MaskKey == 'r':
3333        Masks['Rings'].append([])
3334    G2imG.UpdateMasks(G2frame,Masks)
3335    PlotImage(G2frame,newImage=True)
3336   
3337def OnStartNewDzero(G2frame):
3338    '''Initiate the start of adding a new d-zero to a strain data set
3339
3340    :param wx.Frame G2frame: The main GSAS-II tree "window"
3341    :param str eventkey: a single letter ('a') that
3342      triggers the addition of a d-zero.   
3343    '''
3344    G2frame.dataFrame.GetStatusBar().SetStatusText('Add strain ring active - LB pick d-zero value',0)
3345    G2frame.PickId = G2gd.GetPatternTreeItemId(G2frame,G2frame.Image, 'Stress/Strain')
3346    data = G2frame.PatternTree.GetItemPyData(G2frame.PickId)
3347    return data
3348
3349def PlotImage(G2frame,newPlot=False,event=None,newImage=True):
3350    '''Plot of 2D detector images as contoured plot. Also plot calibration ellipses,
3351    masks, etc.
3352    '''
3353    from matplotlib.patches import Ellipse,Arc,Circle,Polygon
3354    import numpy.ma as ma
3355    Dsp = lambda tth,wave: wave/(2.*npsind(tth/2.))
3356    global Data,Masks,StrSta
3357    colors=['b','g','r','c','m','k']
3358    Data = G2frame.PatternTree.GetItemPyData(
3359        G2gd.GetPatternTreeItemId(G2frame,G2frame.Image, 'Image Controls'))
3360# patch
3361    if 'invert_x' not in Data:
3362        Data['invert_x'] = False
3363        Data['invert_y'] = True
3364# end patch
3365    Masks = G2frame.PatternTree.GetItemPyData(
3366        G2gd.GetPatternTreeItemId(G2frame,G2frame.Image, 'Masks'))
3367    try:    #may be absent
3368        StrSta = G2frame.PatternTree.GetItemPyData(
3369            G2gd.GetPatternTreeItemId(G2frame,G2frame.Image, 'Stress/Strain'))
3370    except TypeError:   #is missing
3371        StrSta = {}
3372
3373    def OnImMotion(event):
3374        Page.canvas.SetToolTipString('')
3375        sizexy = Data['size']
3376        FlatBkg = Data.get('Flat Bkg',0.)
3377        if event.xdata and event.ydata and len(G2frame.ImageZ):                 #avoid out of frame errors
3378            Page.canvas.SetToolTipString('%8.2f %8.2fmm'%(event.xdata,event.ydata))
3379            Page.canvas.SetCursor(wx.CROSS_CURSOR)
3380            item = G2frame.itemPicked
3381            pixelSize = Data['pixelSize']
3382            scalex = 1000./pixelSize[0]
3383            scaley = 1000./pixelSize[1]
3384            if item and G2frame.PatternTree.GetItemText(G2frame.PickId) == 'Image Controls':
3385                if 'Text' in str(item):
3386                    Page.canvas.SetToolTipString('%8.3f %8.3fmm'%(event.xdata,event.ydata))
3387                else:
3388                    xcent,ycent = Data['center']
3389                    xpos = event.xdata-xcent
3390                    ypos = event.ydata-ycent
3391                    tth,azm = G2img.GetTthAzm(event.xdata,event.ydata,Data)
3392                    if 'line3' in  str(item) or 'line4' in str(item) and not Data['fullIntegrate']:
3393                        Page.canvas.SetToolTipString('%6d deg'%(azm))
3394                    elif 'line1' in  str(item) or 'line2' in str(item):
3395                        Page.canvas.SetToolTipString('%8.3f deg'%(tth))                           
3396            else:
3397                xpos = event.xdata
3398                ypos = event.ydata
3399                xpix = xpos*scalex
3400                ypix = ypos*scaley
3401                Int = 0
3402                if (0 <= xpix <= sizexy[0]) and (0 <= ypix <= sizexy[1]):
3403                    Int = G2frame.ImageZ[ypix][xpix]-int(FlatBkg)
3404                tth,azm,D,dsp = G2img.GetTthAzmDsp(xpos,ypos,Data)
3405                Q = 2.*math.pi/dsp
3406                fields = ['','Detector 2-th =%9.3fdeg, dsp =%9.3fA, Q = %6.5fA-1, azm = %7.2fdeg, I = %6d'%(tth,dsp,Q,azm,Int)]
3407                if G2frame.MaskKey in ['p','f']:
3408                    fields[1] = 'Polygon/frame mask pick - LB next point, RB close polygon'
3409                elif G2frame.StrainKey:
3410                    fields[0] = 'd-zero pick active'
3411                G2frame.G2plotNB.status.SetFields(fields)
3412
3413    def OnImPlotKeyPress(event):
3414        try:
3415            PickName = G2frame.PatternTree.GetItemText(G2frame.PickId)
3416        except TypeError:
3417            return
3418        if PickName == 'Masks':
3419            if event.key in ['l','p','f','s','a','r']:
3420                G2frame.MaskKey = event.key
3421                OnStartMask(G2frame)
3422                PlotImage(G2frame,newPlot=False)
3423               
3424        elif PickName == 'Stress/Strain':
3425            if event.key in ['a',]:
3426                G2frame.StrainKey = event.key
3427                StrSta = OnStartNewDzero(G2frame)
3428                PlotImage(G2frame,newPlot=False)
3429               
3430        elif PickName == 'Image Controls':
3431            if event.key in ['c',]:
3432                Xpos = event.xdata
3433                if not Xpos:            #got point out of frame
3434                    return
3435                Ypos = event.ydata
3436                dlg = wx.MessageDialog(G2frame,'Are you sure you want to change the center?',
3437                    'Center change',style=wx.OK|wx.CANCEL)
3438                try:
3439                    if dlg.ShowModal() == wx.ID_OK:
3440                        print 'move center to: ',Xpos,Ypos
3441                        Data['center'] = [Xpos,Ypos]
3442                        G2imG.UpdateImageControls(G2frame,Data,Masks)
3443                        PlotImage(G2frame,newPlot=False)
3444                finally:
3445                    dlg.Destroy()
3446                return
3447            elif event.key == 'l':
3448                G2frame.logPlot = not G2frame.logPlot
3449            elif event.key in ['x',]:
3450                Data['invert_x'] = not Data['invert_x']
3451            elif event.key in ['y',]:
3452                Data['invert_y'] = not Data['invert_y']
3453            PlotImage(G2frame,newPlot=True)
3454           
3455    def OnKeyBox(event):
3456        if G2frame.G2plotNB.nb.GetSelection() == G2frame.G2plotNB.plotList.index('2D Powder Image'):
3457            event.key = cb.GetValue()[0]
3458            cb.SetValue(' key press')
3459            if event.key in ['l','s','a','r','p','x','y']:
3460                wx.CallAfter(OnImPlotKeyPress,event)
3461        Page.canvas.SetFocus() # redirect the Focus from the button back to the plot
3462                       
3463    def OnImPick(event):
3464        if G2frame.PatternTree.GetItemText(G2frame.PickId) not in ['Image Controls','Masks']:
3465            return
3466        if G2frame.itemPicked is not None: return
3467        G2frame.itemPicked = event.artist
3468        G2frame.mousePicked = event.mouseevent
3469       
3470    def OnImRelease(event):
3471        try:
3472            PickName = G2frame.PatternTree.GetItemText(G2frame.PickId)
3473        except TypeError:
3474            return
3475        if PickName not in ['Image Controls','Masks','Stress/Strain']:
3476            return
3477        pixelSize = Data['pixelSize']
3478        FlatBkg = Data.get('Flat Bkg',0.)
3479        scalex = 1000./pixelSize[0]
3480        scaley = 1000./pixelSize[1]
3481#        pixLimit = Data['pixLimit']    #can be too tight
3482        pixLimit = 20       #this makes the search box 40x40 pixels
3483        if G2frame.itemPicked is None and PickName == 'Image Controls' and len(G2frame.ImageZ):
3484            Xpos = event.xdata
3485            if not (Xpos and G2frame.ifGetRing):                   #got point out of frame
3486                return
3487            Ypos = event.ydata
3488            if Ypos and not Page.toolbar._active:         #make sure zoom/pan not selected
3489                if event.button == 1:
3490                    Xpix = Xpos*scalex
3491                    Ypix = Ypos*scaley
3492                    xpos,ypos,I,J = G2img.ImageLocalMax(G2frame.ImageZ-FlatBkg,pixLimit,Xpix,Ypix)
3493                    if I and J:
3494                        xpos += .5                              #shift to pixel center
3495                        ypos += .5
3496                        xpos /= scalex                          #convert to mm
3497                        ypos /= scaley
3498                        Data['ring'].append([xpos,ypos])
3499                elif event.button == 3:
3500                    G2frame.dataFrame.GetStatusBar().SetStatusText('Calibrating...',0)
3501                    if G2img.ImageCalibrate(G2frame,Data):
3502                        G2frame.dataFrame.GetStatusBar().SetStatusText('Calibration successful - Show ring picks to check',0)
3503                        print 'Calibration successful'
3504                    else:
3505                        G2frame.dataFrame.GetStatusBar().SetStatusText('Calibration failed - Show ring picks to diagnose',0)
3506                        print 'Calibration failed'
3507                    G2frame.ifGetRing = False
3508                    G2imG.UpdateImageControls(G2frame,Data,Masks)
3509                    return
3510                PlotImage(G2frame,newImage=False)
3511            return
3512        elif G2frame.MaskKey and PickName == 'Masks':
3513            Xpos,Ypos = [event.xdata,event.ydata]
3514            if not Xpos or not Ypos or Page.toolbar._active:  #got point out of frame or zoom/pan selected
3515                return
3516            if G2frame.MaskKey == 's' and event.button == 1:
3517                Masks['Points'][-1] = [Xpos,Ypos,1.]
3518                G2frame.MaskKey = ''               
3519            elif G2frame.MaskKey == 'r' and event.button == 1:
3520                tth = G2img.GetTth(Xpos,Ypos,Data)
3521                Masks['Rings'][-1] = [tth,0.1]
3522                G2frame.MaskKey = ''               
3523            elif G2frame.MaskKey == 'a' and event.button == 1:
3524                tth,azm = G2img.GetTthAzm(Xpos,Ypos,Data)
3525                azm = int(azm)               
3526                Masks['Arcs'][-1] = [tth,[azm-5,azm+5],0.1]
3527                G2frame.MaskKey = ''               
3528            elif G2frame.MaskKey =='p':
3529                polygon = Masks['Polygons'][-1]
3530                if len(polygon) > 2 and event.button == 3:
3531                    x0,y0 = polygon[0]
3532                    polygon.append([x0,y0])
3533                    G2frame.MaskKey = ''
3534                    G2frame.G2plotNB.status.SetFields(['','Polygon closed'])
3535                else:
3536                    G2frame.G2plotNB.status.SetFields(['','New polygon point: %.1f,%.1f'%(Xpos,Ypos)])
3537                    polygon.append([Xpos,Ypos])
3538            elif G2frame.MaskKey =='f':
3539                frame = Masks['Frames']
3540                if len(frame) > 2 and event.button == 3:
3541                    x0,y0 = frame[0]
3542                    frame.append([x0,y0])
3543                    G2frame.MaskKey = ''
3544                    G2frame.G2plotNB.status.SetFields(['','Frame closed'])
3545                else:
3546                    G2frame.G2plotNB.status.SetFields(['','New frame point: %.1f,%.1f'%(Xpos,Ypos)])
3547                    frame.append([Xpos,Ypos])
3548            G2imG.UpdateMasks(G2frame,Masks)
3549            PlotImage(G2frame,newImage=False)
3550        elif PickName == 'Stress/Strain' and G2frame.StrainKey:
3551            Xpos,Ypos = [event.xdata,event.ydata]
3552            if not Xpos or not Ypos or Page.toolbar._active:  #got point out of frame or zoom/pan selected
3553                return
3554            dsp = G2img.GetDsp(Xpos,Ypos,Data)
3555            StrSta['d-zero'].append({'Dset':dsp,'Dcalc':0.0,'pixLimit':10,'cutoff':0.5,
3556                'ImxyObs':[[],[]],'ImtaObs':[[],[]],'ImtaCalc':[[],[]],'Emat':[1.0,1.0,1.0]})
3557            R,r = G2img.MakeStrStaRing(StrSta['d-zero'][-1],G2frame.ImageZ-FlatBkg,Data)
3558            if not len(R):
3559                del StrSta['d-zero'][-1]
3560                G2frame.ErrorDialog('Strain peak selection','WARNING - No points found for this ring selection')
3561            StrSta['d-zero'] = G2mth.sortArray(StrSta['d-zero'],'Dset',reverse=True)
3562            G2frame.StrainKey = ''
3563            G2imG.UpdateStressStrain(G2frame,StrSta)
3564            PlotStrain(G2frame,StrSta)
3565            PlotImage(G2frame,newPlot=False)           
3566        else:
3567            Xpos,Ypos = [event.xdata,event.ydata]
3568            if not Xpos or not Ypos or Page.toolbar._active:  #got point out of frame or zoom/pan selected
3569                return
3570            if G2frame.ifGetRing:                          #delete a calibration ring pick
3571                xypos = [Xpos,Ypos]
3572                rings = Data['ring']
3573                for ring in rings:
3574                    if np.allclose(ring,xypos,.01,0):
3575                        rings.remove(ring)
3576            else:
3577                tth,azm,dsp = G2img.GetTthAzmDsp(Xpos,Ypos,Data)[:3]
3578                itemPicked = str(G2frame.itemPicked)
3579                if 'Line2D' in itemPicked and PickName == 'Image Controls':
3580                    if 'line1' in itemPicked:
3581                        Data['IOtth'][0] = max(tth,0.001)
3582                    elif 'line2' in itemPicked:
3583                        Data['IOtth'][1] = tth
3584                    elif 'line3' in itemPicked:
3585                        Data['LRazimuth'][0] = int(azm)
3586                    elif 'line4' in itemPicked and not Data['fullIntegrate']:
3587                        Data['LRazimuth'][1] = int(azm)
3588                   
3589                    Data['LRazimuth'][0] %= 360
3590                    Data['LRazimuth'][1] %= 360
3591                    if Data['LRazimuth'][0] > Data['LRazimuth'][1]:
3592                        Data['LRazimuth'][1] += 360                       
3593                    if Data['fullIntegrate']:
3594                        Data['LRazimuth'][1] = Data['LRazimuth'][0]+360
3595                       
3596                    if  Data['IOtth'][0] > Data['IOtth'][1]:
3597                        Data['IOtth'][0],Data['IOtth'][1] = Data['IOtth'][1],Data['IOtth'][0]
3598                       
3599                    G2frame.InnerTth.SetValue("%8.2f" % (Data['IOtth'][0]))
3600                    G2frame.OuterTth.SetValue("%8.2f" % (Data['IOtth'][1]))
3601                    G2frame.Lazim.SetValue("%6d" % (Data['LRazimuth'][0]))
3602                    G2frame.Razim.SetValue("%6d" % (Data['LRazimuth'][1]))
3603                elif 'Circle' in itemPicked and PickName == 'Masks':
3604                    spots = Masks['Points']
3605                    newPos = itemPicked.split(')')[0].split('(')[2].split(',')
3606                    newPos = np.array([float(newPos[0]),float(newPos[1])])
3607                    for spot in spots:
3608                        if spot and np.allclose(np.array([spot[:2]]),newPos):
3609                            spot[:2] = Xpos,Ypos
3610                    G2imG.UpdateMasks(G2frame,Masks)
3611                elif 'Line2D' in itemPicked and PickName == 'Masks':
3612                    Obj = G2frame.itemPicked.findobj()
3613                    rings = Masks['Rings']
3614                    arcs = Masks['Arcs']
3615                    polygons = Masks['Polygons']
3616                    frame = Masks['Frames']
3617                    for ring in G2frame.ringList:
3618                        if Obj == ring[0]:
3619                            rN = ring[1]
3620                            if ring[2] == 'o':
3621                                rings[rN][0] = G2img.GetTth(Xpos,Ypos,Data)-rings[rN][1]/2.
3622                            else:
3623                                rings[rN][0] = G2img.GetTth(Xpos,Ypos,Data)+rings[rN][1]/2.
3624                    for arc in G2frame.arcList:
3625                        if Obj == arc[0]:
3626                            aN = arc[1]
3627                            if arc[2] == 'o':
3628                                arcs[aN][0] = G2img.GetTth(Xpos,Ypos,Data)-arcs[aN][2]/2
3629                            elif arc[2] == 'i':
3630                                arcs[aN][0] = G2img.GetTth(Xpos,Ypos,Data)+arcs[aN][2]/2
3631                            elif arc[2] == 'l':
3632                                arcs[aN][1][0] = int(G2img.GetAzm(Xpos,Ypos,Data))
3633                            else:
3634                                arcs[aN][1][1] = int(G2img.GetAzm(Xpos,Ypos,Data))
3635                    for poly in G2frame.polyList:   #merging points problem here?
3636                        if Obj == poly[0]:
3637                            ind = G2frame.itemPicked.contains(G2frame.mousePicked)[1]['ind'][0]
3638                            oldPos = np.array([G2frame.mousePicked.xdata,G2frame.mousePicked.ydata])
3639                            pN = poly[1]
3640                            for i,xy in enumerate(polygons[pN]):
3641                                if np.allclose(np.array([xy]),oldPos,atol=1.0):
3642                                    if event.button == 1:
3643                                        polygons[pN][i] = Xpos,Ypos
3644                                    elif event.button == 3:
3645                                        polygons[pN].insert(i,[Xpos,Ypos])
3646                                        break
3647                    if frame:
3648                        oldPos = np.array([G2frame.mousePicked.xdata,G2frame.mousePicked.ydata])
3649                        for i,xy in enumerate(frame):
3650                            if np.allclose(np.array([xy]),oldPos,atol=1.0):
3651                                if event.button == 1:
3652                                    frame[i] = Xpos,Ypos
3653                                elif event.button == 3:
3654                                    frame.insert(i,[Xpos,Ypos])
3655                                    break
3656                    G2imG.UpdateMasks(G2frame,Masks)
3657#                else:                  #keep for future debugging
3658#                    print str(G2frame.itemPicked),event.xdata,event.ydata,event.button
3659            PlotImage(G2frame,newImage=True)
3660            G2frame.itemPicked = None
3661           
3662    try:
3663        plotNum = G2frame.G2plotNB.plotList.index('2D Powder Image')
3664        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
3665        if not newPlot:
3666            Plot = Page.figure.gca()          #get previous powder plot & get limits
3667            xylim = Plot.get_xlim(),Plot.get_ylim()
3668        if newImage:
3669            Page.figure.clf()
3670            Plot = Page.figure.gca()          #get a fresh plot after clf()
3671    except ValueError:
3672        Plot = G2frame.G2plotNB.addMpl('2D Powder Image').gca()
3673        plotNum = G2frame.G2plotNB.plotList.index('2D Powder Image')
3674        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
3675        Page.canvas.mpl_connect('key_press_event', OnImPlotKeyPress)
3676        Page.canvas.mpl_connect('motion_notify_event', OnImMotion)
3677        Page.canvas.mpl_connect('pick_event', OnImPick)
3678        Page.canvas.mpl_connect('button_release_event', OnImRelease)
3679        xylim = []
3680    Page.Choice = None
3681    if not event:                       #event from GUI TextCtrl - don't want focus to change to plot!!!
3682        Page.SetFocus()
3683    Title = G2frame.PatternTree.GetItemText(G2frame.Image)[4:]
3684    G2frame.G2plotNB.status.DestroyChildren()
3685    if G2frame.logPlot:
3686        Title = 'log('+Title+')'
3687    Plot.set_title(Title)
3688    try:
3689        if G2frame.PatternTree.GetItemText(G2frame.PickId) in ['Image Controls',]:
3690            Page.Choice = (' key press','l: log(I) on','x: flip x','y: flip y',)
3691            if G2frame.logPlot:
3692                Page.Choice[1] = 'l: log(I) off'
3693            Page.keyPress = OnImPlotKeyPress
3694        elif G2frame.PatternTree.GetItemText(G2frame.PickId) in ['Masks',]:
3695            Page.Choice = (' key press','l: log(I) on','s: spot mask','a: arc mask','r: ring mask',
3696                'p: polygon mask','f: frame mask',)
3697            if G2frame.logPlot:
3698                Page.Choice[1] = 'l: log(I) off'
3699            Page.keyPress = OnImPlotKeyPress
3700        elif G2frame.PatternTree.GetItemText(G2frame.PickId) in ['Stress/Strain',]:
3701            Page.Choice = (' key press','a: add new ring',)
3702            Page.keyPress = OnImPlotKeyPress
3703    except TypeError:
3704        pass
3705    size,imagefile = G2frame.PatternTree.GetItemPyData(G2frame.Image)
3706    dark = Data.get('dark image',[0,''])
3707    if dark[0]:
3708        darkfile = G2frame.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(G2frame, 
3709            G2frame.root,dark[0]))[1]
3710    if imagefile != G2frame.oldImagefile:
3711        imagefile = G2IO.CheckImageFile(G2frame,imagefile)
3712        if not imagefile:
3713            G2frame.G2plotNB.Delete('2D Powder Image')
3714            return
3715        G2frame.PatternTree.SetItemPyData(G2frame.Image,[size,imagefile])
3716        G2frame.ImageZ = G2IO.GetImageData(G2frame,imagefile,imageOnly=True)
3717        if dark[0]:
3718            darkImg = G2IO.GetImageData(G2frame,darkfile,imageOnly=True)
3719            G2frame.ImageZ += dark[1]*darkImg
3720        G2frame.oldImagefile = imagefile
3721
3722    imScale = 1
3723    if len(G2frame.ImageZ) > 1024:
3724        imScale = len(G2frame.ImageZ)/1024
3725    sizexy = Data['size']
3726    pixelSize = Data['pixelSize']
3727    scalex = 1000./pixelSize[0]
3728    scaley = 1000./pixelSize[1]
3729    Xmax = sizexy[0]*pixelSize[0]/1000.
3730    Ymax = sizexy[1]*pixelSize[1]/1000.
3731    xlim = (0,Xmax)
3732    ylim = (Ymax,0)
3733    Imin,Imax = Data['range'][1]
3734    acolor = mpl.cm.get_cmap(Data['color'])
3735    xcent,ycent = Data['center']
3736    Plot.set_xlabel('Image x-axis, mm',fontsize=12)
3737    Plot.set_ylabel('Image y-axis, mm',fontsize=12)
3738    #do threshold mask - "real" mask - others are just bondaries
3739    Zlim = Masks['Thresholds'][1]
3740    FlatBkg = Data.get('Flat Bkg',0.0)
3741    wx.BeginBusyCursor()
3742    try:
3743           
3744        if newImage:                   
3745            Imin,Imax = Data['range'][1]
3746            MA = ma.masked_greater(ma.masked_less(G2frame.ImageZ,Zlim[0]+FlatBkg),Zlim[1]+FlatBkg)
3747            MaskA = ma.getmaskarray(MA)
3748            A = G2img.ImageCompress(MA,imScale)-FlatBkg
3749            AM = G2img.ImageCompress(MaskA,imScale)
3750            if G2frame.logPlot:
3751                A = np.where(A>Imin,np.where(A<Imax,A,0),0)
3752                A = np.where(A>0,np.log(A),0)
3753                AM = np.where(AM>0,np.log(AM),0)
3754                Imin,Imax = [np.amin(A),np.amax(A)]
3755            ImgM = Plot.imshow(AM,aspect='equal',cmap='Reds',
3756                interpolation='nearest',vmin=0,vmax=2,extent=[0,Xmax,Ymax,0])
3757            Img = Plot.imshow(A,aspect='equal',cmap=acolor,
3758                interpolation='nearest',vmin=Imin,vmax=Imax,extent=[0,Xmax,Ymax,0])
3759   
3760        Plot.plot(xcent,ycent,'x')
3761        #G2frame.PatternTree.GetItemText(item)
3762        if Data['showLines']:
3763            LRAzim = Data['LRazimuth']                  #NB: integers
3764            Nazm = Data['outAzimuths']
3765            delAzm = float(LRAzim[1]-LRAzim[0])/Nazm
3766            AzmthOff = Data['azmthOff']
3767            IOtth = Data['IOtth']
3768            wave = Data['wavelength']
3769            dspI = wave/(2.0*sind(IOtth[0]/2.0))
3770            ellI = G2img.GetEllipse(dspI,Data)           #=False if dsp didn't yield an ellipse (ugh! a parabola or a hyperbola)
3771            dspO = wave/(2.0*sind(IOtth[1]/2.0))
3772            ellO = G2img.GetEllipse(dspO,Data)           #Ditto & more likely for outer ellipse
3773            Azm = np.array(range(LRAzim[0],LRAzim[1]+1))-AzmthOff
3774            if ellI:
3775                xyI = []
3776                for azm in Azm:
3777                    xy = G2img.GetDetectorXY(dspI,azm,Data)
3778                    if np.any(xy):
3779                        xyI.append(xy)
3780                if len(xyI):
3781                    xyI = np.array(xyI)
3782                    arcxI,arcyI = xyI.T
3783                    Plot.plot(arcxI,arcyI,picker=3)
3784            if ellO:
3785                xyO = []
3786                for azm in Azm:
3787                    xy = G2img.GetDetectorXY(dspO,azm,Data)
3788                    if np.any(xy):
3789                        xyO.append(xy)
3790                if len(xyO):
3791                    xyO = np.array(xyO)
3792                    arcxO,arcyO = xyO.T               
3793                    Plot.plot(arcxO,arcyO,picker=3)
3794            if ellO and ellI:
3795                Plot.plot([arcxI[0],arcxO[0]],[arcyI[0],arcyO[0]],picker=3)
3796                Plot.plot([arcxI[-1],arcxO[-1]],[arcyI[-1],arcyO[-1]],picker=3)
3797            for i in range(Nazm):
3798                cake = LRAzim[0]+i*delAzm-AzmthOff
3799                if Data['centerAzm']:
3800                    cake += delAzm/2.
3801                ind = np.searchsorted(Azm,cake)
3802                Plot.plot([arcxI[ind],arcxO[ind]],[arcyI[ind],arcyO[ind]],color='k',dashes=(5,5))
3803                   
3804        if G2frame.PatternTree.GetItemText(G2frame.PickId) in 'Image Controls':
3805            for xring,yring in Data['ring']:
3806                Plot.plot(xring,yring,'r+',picker=3)
3807            if Data['setRings']:
3808                N = 0
3809                for ring in Data['rings']:
3810                    xring,yring = np.array(ring).T[:2]
3811                    Plot.plot(xring,yring,'.',color=colors[N%6])
3812                    N += 1           
3813            for ellipse in Data['ellipses']:      #what about hyperbola?
3814                cent,phi,[width,height],col = ellipse
3815                if width > 0:       #ellipses
3816                    Plot.add_artist(Ellipse([cent[0],cent[1]],2*width,2*height,phi,ec=col,fc='none'))
3817                    Plot.text(cent[0],cent[1],'+',color=col,ha='center',va='center')
3818        if G2frame.PatternTree.GetItemText(G2frame.PickId) in 'Stress/Strain':
3819            for N,ring in enumerate(StrSta['d-zero']):
3820                xring,yring = ring['ImxyObs']
3821                Plot.plot(xring,yring,colors[N%6]+'.')
3822        #masks - mask lines numbered after integration limit lines
3823        spots = Masks['Points']
3824        rings = Masks['Rings']
3825        arcs = Masks['Arcs']
3826        polygons = Masks['Polygons']
3827        if 'Frames' not in Masks:
3828            Masks['Frames'] = []
3829        frame = Masks['Frames']
3830        for spot in spots:
3831            if spot:
3832                x,y,d = spot
3833                Plot.add_artist(Circle((x,y),radius=d/2,fc='none',ec='r',picker=3))
3834        G2frame.ringList = []
3835        for iring,ring in enumerate(rings):
3836            if ring:
3837                tth,thick = ring
3838                wave = Data['wavelength']
3839                xy1 = []
3840                xy2 = []
3841                Azm = np.linspace(0,362,181)
3842                for azm in Azm:
3843                    xy1.append(G2img.GetDetectorXY(Dsp(tth+thick/2.,wave),azm,Data))      #what about hyperbola
3844                    xy2.append(G2img.GetDetectorXY(Dsp(tth-thick/2.,wave),azm,Data))      #what about hyperbola
3845                x1,y1 = np.array(xy1).T
3846                x2,y2 = np.array(xy2).T
3847                G2frame.ringList.append([Plot.plot(x1,y1,'r',picker=3),iring,'o'])           
3848                G2frame.ringList.append([Plot.plot(x2,y2,'r',picker=3),iring,'i'])
3849        G2frame.arcList = []
3850        for iarc,arc in enumerate(arcs):
3851            if arc:
3852                tth,azm,thick = arc           
3853                wave = Data['wavelength']
3854                xy1 = []
3855                xy2 = []
3856                aR = azm[0],azm[1],azm[1]-azm[0]
3857                if azm[1]-azm[0] > 180:
3858                    aR[2] /= 2
3859                Azm = np.linspace(aR[0],aR[1],aR[2])
3860                for azm in Azm:
3861                    xy1.append(G2img.GetDetectorXY(Dsp(tth+thick/2.,wave),azm,Data))      #what about hyperbola
3862                    xy2.append(G2img.GetDetectorXY(Dsp(tth-thick/2.,wave),azm,Data))      #what about hyperbola
3863                x1,y1 = np.array(xy1).T
3864                x2,y2 = np.array(xy2).T
3865                G2frame.arcList.append([Plot.plot(x1,y1,'r',picker=3),iarc,'o'])           
3866                G2frame.arcList.append([Plot.plot(x2,y2,'r',picker=3),iarc,'i'])
3867                G2frame.arcList.append([Plot.plot([x1[0],x2[0]],[y1[0],y2[0]],'r',picker=3),iarc,'l'])
3868                G2frame.arcList.append([Plot.plot([x1[-1],x2[-1]],[y1[-1],y2[-1]],'r',picker=3),iarc,'u'])
3869        G2frame.polyList = []
3870        for ipoly,polygon in enumerate(polygons):
3871            if polygon:
3872                x,y = np.hsplit(np.array(polygon),2)
3873                G2frame.polyList.append([Plot.plot(x,y,'r+',picker=10),ipoly])
3874                Plot.plot(x,y,'r')           
3875        G2frame.frameList = []
3876        if frame:
3877            x,y = np.hsplit(np.array(frame),2)
3878            G2frame.frameList.append([Plot.plot(x,y,'g+',picker=10),0])
3879            Plot.plot(x,y,'g')           
3880        if newImage:
3881            colorBar = Page.figure.colorbar(Img)
3882        Plot.set_xlim(xlim)
3883        Plot.set_ylim(ylim)
3884        if Data['invert_x']:
3885            Plot.invert_xaxis()
3886        if Data['invert_y']:
3887            Plot.invert_yaxis()
3888        if not newPlot and xylim:
3889            Page.toolbar.push_current()
3890            Plot.set_xlim(xylim[0])
3891            Plot.set_ylim(xylim[1])
3892            xylim = []
3893            Page.toolbar.push_current()
3894            Page.toolbar.draw()
3895            # patch for wx 2.9 on Mac, to force a redraw
3896            i,j= wx.__version__.split('.')[0:2]
3897            if int(i)+int(j)/10. > 2.8 and 'wxOSX' in wx.PlatformInfo:
3898                Page.canvas.draw()
3899        else:
3900            Page.canvas.draw()
3901    finally:
3902        wx.EndBusyCursor()
3903       
3904################################################################################
3905##### PlotIntegration
3906################################################################################
3907           
3908def PlotIntegration(G2frame,newPlot=False,event=None):
3909    '''Plot of 2D image after image integration with 2-theta and azimuth as coordinates
3910    '''
3911           
3912    def OnMotion(event):
3913        Page.canvas.SetToolTipString('')
3914        Page.canvas.SetCursor(wx.CROSS_CURSOR)
3915        azm = event.ydata
3916        tth = event.xdata
3917        if azm and tth:
3918            G2frame.G2plotNB.status.SetFields(\
3919                ['','Detector 2-th =%9.3fdeg, azm = %7.2fdeg'%(tth,azm)])
3920                               
3921    try:
3922        plotNum = G2frame.G2plotNB.plotList.index('2D Integration')
3923        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
3924        if not newPlot:
3925            Plot = Page.figure.gca()          #get previous plot & get limits
3926            xylim = Plot.get_xlim(),Plot.get_ylim()
3927        Page.figure.clf()
3928        Plot = Page.figure.gca()          #get a fresh plot after clf()
3929       
3930    except ValueError:
3931        Plot = G2frame.G2plotNB.addMpl('2D Integration').gca()
3932        plotNum = G2frame.G2plotNB.plotList.index('2D Integration')
3933        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
3934        Page.canvas.mpl_connect('motion_notify_event', OnMotion)
3935        Page.views = False
3936        view = False
3937    Page.Choice = None
3938    if not event:
3939        Page.SetFocus()
3940       
3941    Data = G2frame.PatternTree.GetItemPyData(
3942        G2gd.GetPatternTreeItemId(G2frame,G2frame.Image, 'Image Controls'))
3943    image = G2frame.Integrate[0]
3944    xsc = G2frame.Integrate[1]
3945    ysc = G2frame.Integrate[2]
3946    Imin,Imax = Data['range'][1]
3947    acolor = mpl.cm.get_cmap(Data['color'])
3948    Plot.set_title(G2frame.PatternTree.GetItemText(G2frame.Image)[4:])
3949    Plot.set_ylabel('azimuth',fontsize=12)
3950    Plot.set_xlabel('2-theta',fontsize=12)
3951    Img = Plot.imshow(image,cmap=acolor,vmin=Imin,vmax=Imax,interpolation='nearest', \
3952        extent=[ysc[0],ysc[-1],xsc[-1],xsc[0]],aspect='auto')
3953    colorBar = Page.figure.colorbar(Img)
3954#    if Data['ellipses']:           
3955#        for ellipse in Data['ellipses']:
3956#            x,y = np.array(G2img.makeIdealRing(ellipse[:3])) #skip color
3957#            tth,azm = G2img.GetTthAzm(x,y,Data)
3958##            azm = np.where(azm < 0.,azm+360,azm)
3959#            Plot.plot(tth,azm,'b,')
3960    if not newPlot:
3961        Page.toolbar.push_current()
3962        Plot.set_xlim(xylim[0])
3963        Plot.set_ylim(xylim[1])
3964        xylim = []
3965        Page.toolbar.push_current()
3966        Page.toolbar.draw()
3967    else:
3968        Page.canvas.draw()
3969               
3970################################################################################
3971##### PlotTRImage
3972################################################################################
3973           
3974def PlotTRImage(G2frame,tax,tay,taz,newPlot=False):
3975    '''a test plot routine - not normally used
3976    ''' 
3977           
3978    def OnMotion(event):
3979        Page.canvas.SetToolTipString('')
3980        Page.canvas.SetCursor(wx.CROSS_CURSOR)
3981        azm = event.xdata
3982        tth = event.ydata
3983        if azm and tth:
3984            G2frame.G2plotNB.status.SetFields(\
3985                ['','Detector 2-th =%9.3fdeg, azm = %7.2fdeg'%(tth,azm)])
3986                               
3987    try:
3988        plotNum = G2frame.G2plotNB.plotList.index('2D Transformed Powder Image')
3989        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
3990        if not newPlot:
3991            Plot = Page.figure.gca()          #get previous plot & get limits
3992            xylim = Plot.get_xlim(),Plot.get_ylim()
3993        Page.figure.clf()
3994        Plot = Page.figure.gca()          #get a fresh plot after clf()
3995       
3996    except ValueError:
3997        Plot = G2frame.G2plotNB.addMpl('2D Transformed Powder Image').gca()
3998        plotNum = G2frame.G2plotNB.plotList.index('2D Transformed Powder Image')
3999        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
4000        Page.canvas.mpl_connect('motion_notify_event', OnMotion)
4001        Page.views = False
4002        view = False
4003    Page.Choice = None
4004    Page.SetFocus()
4005       
4006    Data = G2frame.PatternTree.GetItemPyData(
4007        G2gd.GetPatternTreeItemId(G2frame,G2frame.Image, 'Image Controls'))
4008    Imin,Imax = Data['range'][1]
4009    step = (Imax-Imin)/5.
4010    V = np.arange(Imin,Imax,step)
4011    acolor = mpl.cm.get_cmap(Data['color'])
4012    Plot.set_title(G2frame.PatternTree.GetItemText(G2frame.Image)[4:])
4013    Plot.set_xlabel('azimuth',fontsize=12)
4014    Plot.set_ylabel('2-theta',fontsize=12)
4015    Plot.contour(tax,tay,taz,V,cmap=acolor)
4016    if Data['showLines']:
4017        IOtth = Data['IOtth']
4018        if Data['fullIntegrate']:
4019            LRAzim = [-180,180]
4020        else:
4021            LRAzim = Data['LRazimuth']                  #NB: integers
4022        Plot.plot([LRAzim[0],LRAzim[1]],[IOtth[0],IOtth[0]],picker=True)
4023        Plot.plot([LRAzim[0],LRAzim[1]],[IOtth[1],IOtth[1]],picker=True)
4024        if not Data['fullIntegrate']:
4025            Plot.plot([LRAzim[0],LRAzim[0]],[IOtth[0],IOtth[1]],picker=True)
4026            Plot.plot([LRAzim[1],LRAzim[1]],[IOtth[0],IOtth[1]],picker=True)
4027    if Data['setRings']:
4028        rings = np.concatenate((Data['rings']),axis=0)
4029        for xring,yring,dsp in rings:
4030            x,y = G2img.GetTthAzm(xring,yring,Data)
4031            Plot.plot(y,x,'r+')           
4032    if Data['ellipses']:           
4033        for ellipse in Data['ellipses']:
4034            ring = np.array(G2img.makeIdealRing(ellipse[:3])) #skip color
4035            x,y = np.hsplit(ring,2)
4036            tth,azm = G2img.GetTthAzm(x,y,Data)
4037            Plot.plot(azm,tth,'b,')
4038    if not newPlot:
4039        Page.toolbar.push_current()
4040        Plot.set_xlim(xylim[0])
4041        Plot.set_ylim(xylim[1])
4042        xylim = []
4043        Page.toolbar.push_current()
4044        Page.toolbar.draw()
4045    else:
4046        Page.canvas.draw()
4047       
4048################################################################################
4049##### PlotStructure
4050################################################################################
4051           
4052def PlotStructure(G2frame,data,firstCall=False):
4053    '''Crystal structure plotting package. Can show structures as balls, sticks, lines,
4054    thermal motion ellipsoids and polyhedra
4055    '''
4056
4057    def FindPeaksBonds(XYZ):
4058        rFact = data['Drawing'].get('radiusFactor',0.85)    #data['Drawing'] could be empty!
4059        Bonds = [[] for x in XYZ]
4060        for i,xyz in enumerate(XYZ):
4061            Dx = XYZ-xyz
4062            dist = np.sqrt(np.sum(np.inner(Dx,Amat)**2,axis=1))
4063            IndB = ma.nonzero(ma.masked_greater(dist,rFact*2.2))
4064            for j in IndB[0]:
4065                Bonds[i].append(Dx[j]/2.)
4066                Bonds[j].append(-Dx[j]/2.)
4067        return Bonds
4068
4069    # PlotStructure initialization here
4070    ForthirdPI = 4.0*math.pi/3.0
4071    generalData = data['General']
4072    cell = generalData['Cell'][1:7]
4073    Vol = generalData['Cell'][7:8][0]
4074    Amat,Bmat = G2lat.cell2AB(cell)         #Amat - crystal to cartesian, Bmat - inverse
4075    Gmat,gmat = G2lat.cell2Gmat(cell)
4076    A4mat = np.concatenate((np.concatenate((Amat,[[0],[0],[0]]),axis=1),[[0,0,0,1],]),axis=0)
4077    B4mat = np.concatenate((np.concatenate((Bmat,[[0],[0],[0]]),axis=1),[[0,0,0,1],]),axis=0)
4078    SGData = generalData['SGData']
4079    if generalData['Type'] in ['modulated','magnetic']:
4080        SSGData = generalData['SSGData']
4081    Mydir = generalData['Mydir']
4082    atomData = data['Atoms']
4083    mapPeaks = []
4084    drawingData = data['Drawing']
4085    if not drawingData:
4086        return          #nothing setup, nothing to draw   
4087    if 'Map Peaks' in data:
4088        mapPeaks = np.array(data['Map Peaks'])
4089        peakMax = 100.
4090        if len(mapPeaks):
4091            peakMax = np.max(mapPeaks.T[0])
4092    resRBData = data['RBModels'].get('Residue',[])
4093    vecRBData = data['RBModels'].get('Vector',[])
4094    rbAtmDict = {}
4095    for rbObj in resRBData+vecRBData:
4096        exclList = ['X' for i in range(len(rbObj['Ids']))]
4097        rbAtmDict.update(dict(zip(rbObj['Ids'],exclList)))
4098    testRBObj = data.get('testRBObj',{})
4099    rbObj = testRBObj.get('rbObj',{})
4100    MCSA = data.get('MCSA',{})
4101    mcsaModels = MCSA.get('Models',[])
4102    if len(mcsaModels) > 1:
4103        XYZs,Types = G2mth.UpdateMCSAxyz(Bmat,MCSA)
4104        mcsaXYZ = []
4105        mcsaTypes = []
4106        neqv = 0
4107        for xyz,atyp in zip(XYZs,Types):
4108            equiv = G2spc.GenAtom(xyz,SGData,All=True,Move=False)
4109            neqv = max(neqv,len(equiv))
4110            for item in equiv:
4111                mcsaXYZ.append(item[0]) 
4112                mcsaTypes.append(atyp)
4113        mcsaXYZ = np.array(mcsaXYZ)
4114        mcsaTypes = np.array(mcsaTypes)
4115        nuniq = mcsaXYZ.shape[0]/neqv
4116        mcsaXYZ = np.reshape(mcsaXYZ,(nuniq,neqv,3))
4117        mcsaTypes = np.reshape(mcsaTypes,(nuniq,neqv))
4118        cent = np.fix(np.sum(mcsaXYZ+2.,axis=0)/nuniq)-2
4119        cent[0] = [0,0,0]   #make sure 1st one isn't moved
4120        mcsaXYZ = np.swapaxes(mcsaXYZ,0,1)-cent[:,np.newaxis,:]
4121        mcsaTypes = np.swapaxes(mcsaTypes,0,1)
4122        mcsaXYZ = np.reshape(mcsaXYZ,(nuniq*neqv,3))
4123        mcsaTypes = np.reshape(mcsaTypes,(nuniq*neqv))                       
4124        mcsaBonds = FindPeaksBonds(mcsaXYZ)       
4125    drawAtoms = drawingData.get('Atoms',[])
4126    mapData = {}
4127    flipData = {}
4128    rhoXYZ = []
4129    showBonds = False
4130    if 'Map' in generalData:
4131        mapData = generalData['Map']
4132        showBonds = mapData.get('Show bonds',False)
4133    if 'Flip' in generalData:
4134        flipData = generalData['Flip']                       
4135    Wt = np.array([255,255,255])
4136    Rd = np.array([255,0,0])
4137    Gr = np.array([0,255,0])
4138    wxGreen = wx.Colour(0,255,0)
4139    Bl = np.array([0,0,255])
4140    Or = np.array([255,128,0])
4141    wxOrange = wx.Colour(255,128,0)
4142    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]])
4143    uEdges = np.array([
4144        [uBox[0],uBox[1]],[uBox[0],uBox[3]],[uBox[0],uBox[4]],[uBox[1],uBox[2]], 
4145        [uBox[2],uBox[3]],[uBox[1],uBox[5]],[uBox[2],uBox[6]],[uBox[3],uBox[7]], 
4146        [uBox[4],uBox[5]],[uBox[5],uBox[6]],[uBox[6],uBox[7]],[uBox[7],uBox[4]]])
4147    mD = 0.1
4148    mV = np.array([[[-mD,0,0],[mD,0,0]],[[0,-mD,0],[0,mD,0]],[[0,0,-mD],[0,0,mD]]])
4149    mapPeakVecs = np.inner(mV,Bmat)
4150
4151    backColor = np.array(list(drawingData['backColor'])+[0,])
4152    Bc = np.array(list(drawingData['backColor']))
4153    uColors = [Rd,Gr,Bl,Wt-Bc, Wt-Bc,Wt-Bc,Wt-Bc,Wt-Bc, Wt-Bc,Wt-Bc,Wt-Bc,Wt-Bc]
4154    altDown = False
4155    shiftDown = False
4156    ctrlDown = False
4157    G2frame.tau = 0.
4158   
4159    def OnKeyBox(event):
4160        mode = cb.GetValue()
4161        if mode in ['jpeg','bmp','tiff',]:
4162            try:
4163                import Image as Im
4164            except ImportError:
4165                try:
4166                    from PIL import Image as Im
4167                except ImportError:
4168                    print "PIL/pillow Image module not present. Cannot save images without this"
4169                    raise Exception("PIL/pillow Image module not found")
4170            Fname = os.path.join(Mydir,generalData['Name']+'.'+mode)
4171            size = Page.canvas.GetSize()
4172            glPixelStorei(GL_UNPACK_ALIGNMENT, 1)
4173            if mode in ['jpeg',]:
4174                Pix = glReadPixels(0,0,size[0],size[1],GL_RGBA, GL_UNSIGNED_BYTE)
4175                im = Im.new("RGBA", (size[0],size[1]))
4176            else:
4177                Pix = glReadPixels(0,0,size[0],size[1],GL_RGB, GL_UNSIGNED_BYTE)
4178                im = Im.new("RGB", (size[0],size[1]))
4179            im.fromstring(Pix)
4180            im.save(Fname,mode)
4181            cb.SetValue(' save as/key:')
4182            G2frame.G2plotNB.status.SetStatusText('Drawing saved to: '+Fname,1)
4183        else:
4184            event.key = cb.GetValue()[0]
4185            cb.SetValue(' save as/key:')
4186            wx.CallAfter(OnKey,event)
4187        Page.canvas.SetFocus() # redirect the Focus from the button back to the plot
4188
4189    def OnKey(event):           #on key UP!!
4190        try:
4191            keyCode = event.GetKeyCode()
4192            if keyCode > 255:
4193                keyCode = 0
4194            key = chr(keyCode)
4195        except AttributeError:       #if from OnKeyBox above
4196            key = str(event.key).upper()
4197        indx = drawingData['selectedAtoms']
4198        cx,ct = drawingData['atomPtrs'][:2]
4199        if key in ['C']:
4200            drawingData['viewPoint'] = [[.5,.5,.5],[0,0]]
4201            drawingData['viewDir'] = [0,0,1]
4202            drawingData['oldxy'] = []
4203            V0 = np.array([0,0,1])
4204            V = np.inner(Amat,V0)
4205            V /= np.sqrt(np.sum(V**2))
4206            A = np.arccos(np.sum(V*V0))
4207            Q = G2mth.AV2Q(A,[0,1,0])
4208            drawingData['Quaternion'] = Q
4209            SetViewPointText(drawingData['viewPoint'][0])
4210            SetViewDirText(drawingData['viewDir'])
4211            Q = drawingData['Quaternion']
4212            G2frame.G2plotNB.status.SetStatusText('New quaternion: %.2f+, %.2fi+ ,%.2fj+, %.2fk'%(Q[0],Q[1],Q[2],Q[3]),1)
4213        elif key in ['N']:
4214            drawAtoms = drawingData['Atoms']
4215            if not len(drawAtoms):      #no atoms
4216                return
4217            pI = drawingData['viewPoint'][1]
4218            if not len(pI):
4219                pI = [0,0]
4220            if indx:
4221                pI[0] = indx[pI[1]]
4222                Tx,Ty,Tz = drawAtoms[pI[0]][cx:cx+3]
4223                pI[1] += 1
4224                if pI[1] >= len(indx):
4225                    pI[1] = 0
4226            else:
4227                Tx,Ty,Tz = drawAtoms[pI[0]][cx:cx+3]               
4228                pI[0] += 1
4229                if pI[0] >= len(drawAtoms):
4230                    pI[0] = 0
4231            drawingData['viewPoint'] = [[Tx,Ty,Tz],pI]
4232            SetViewPointText(drawingData['viewPoint'][0])
4233            G2frame.G2plotNB.status.SetStatusText('View point at atom '+drawAtoms[pI[0]][ct-1]+str(pI),1)
4234               
4235        elif key in ['P']:
4236            drawAtoms = drawingData['Atoms']
4237            if not len(drawAtoms):      #no atoms
4238                return
4239            pI = drawingData['viewPoint'][1]
4240            if not len(pI):
4241                pI = [0,0]
4242            if indx:
4243                pI[0] = indx[pI[1]]
4244                Tx,Ty,Tz = drawAtoms[pI[0]][cx:cx+3]
4245                pI[1] -= 1
4246                if pI[1] < 0:
4247                    pI[1] = len(indx)-1
4248            else:
4249                Tx,Ty,Tz = drawAtoms[pI[0]][cx:cx+3]               
4250                pI[0] -= 1
4251                if pI[0] < 0:
4252                    pI[0] = len(drawAtoms)-1
4253            drawingData['viewPoint'] = [[Tx,Ty,Tz],pI]
4254            SetViewPointText(drawingData['viewPoint'][0])           
4255            G2frame.G2plotNB.status.SetStatusText('View point at atom '+drawAtoms[pI[0]][ct-1]+str(pI),1)
4256        elif key in ['U','D','L','R'] and mapData['Flip'] == True:
4257            dirDict = {'U':[0,1],'D':[0,-1],'L':[-1,0],'R':[1,0]}
4258            SetMapRoll(dirDict[key])
4259            if 'rho' in generalData.get('4DmapData',{}):
4260                Set4DMapRoll(dirDict[key])
4261            SetPeakRoll(dirDict[key])
4262            SetMapPeaksText(mapPeaks)
4263        elif key in ['+','-','=','0'] and generalData['Type'] in ['modulated','magnetic']:
4264            if key == '0':
4265                G2frame.tau = 0.
4266            elif key in ['+','=']:
4267                G2frame.tau += 0.05
4268            elif key == '-':
4269                G2frame.tau -= 0.05
4270            G2frame.tau %= 1.   #force 0-1 range
4271        Draw('key')
4272           
4273    def GetTruePosition(xy,Add=False):
4274        View = glGetIntegerv(GL_VIEWPORT)
4275        Proj = glGetDoublev(GL_PROJECTION_MATRIX)
4276        Model = glGetDoublev(GL_MODELVIEW_MATRIX)
4277        Zmax = 1.
4278        if Add:
4279            Indx = GetSelectedAtoms()
4280        if G2frame.dataDisplay.GetPageText(getSelection()) == 'Map peaks':
4281            for i,peak in enumerate(mapPeaks):
4282                x,y,z = peak[1:4]
4283                X,Y,Z = gluProject(x,y,z,Model,Proj,View)
4284                XY = [int(X),int(View[3]-Y)]
4285                if np.allclose(xy,XY,atol=10) and Z < Zmax:
4286                    Zmax = Z
4287                    try:
4288                        Indx.remove(i)
4289                        ClearSelectedAtoms()
4290                        for id in Indx:
4291                            SetSelectedAtoms(id,Add)
4292                    except:
4293                        SetSelectedAtoms(i,Add)
4294        else:
4295            cx = drawingData['atomPtrs'][0]
4296            for i,atom in enumerate(drawAtoms):
4297                x,y,z = atom[cx:cx+3]