source: trunk/GSASIIplot.py @ 1958

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

add r-values by SS index
add a fade for frac modulation of atoms in drawing
work on modulation structure factor

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