source: trunk/GSASIIplot.py @ 1630

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

add export routine for tables of x,y0,y1,y2,... for multiple powder data sets - apparently useful for plotting in Excel, etc.
implement SS structure factor calcs for position modulation - Fourier only
working on printing wave results

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