source: trunk/GSASIIplot.py @ 1799

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

trap missing data error
numpy array sph.harm. texture fit function & derivative routines - speedier
change pole figure titles

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