source: trunk/GSASIIplot.py @ 1694

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

some missing stuff in PWDR[0] dictionary
add new fxn to G2img - meanAZm - calc effective azimuth for a azm range for correct polarization calc.
fix some FlexGridSizers? in G2phsGUI
fix '+'/'=' key response for Anaconda vs Enthought versions; one sees '=' the other sees '+'!
fix L/R shifting in waterfall plots.

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