source: trunk/GSASIIplot.py @ 1743

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

fix to data selection issue
formatting error for histogram scale
change ':' to ';' for SASD parameter names

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