source: trunk/GSASIIplot.py @ 1656

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

show succession of 3D HKL plots when switching histograms with arrow keys

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