source: trunk/GSASIIplot.py @ 1763

Last change on this file since 1763 was 1763, checked in by vondreele, 8 years ago

new item in imgGUI - Flat Bkg; subtracted from image on display & integration

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