source: trunk/GSASIIplot.py @ 796

Last change on this file since 796 was 796, checked in by vondreele, 10 years ago

revise instrument parameters (again); now two dictionaries
revise key stroke options on plots

  • Property svn:eol-style set to native
  • Property svn:keywords set to Date Author Revision URL Id
File size: 133.0 KB
Line 
1# -*- coding: utf-8 -*-
2#GSASII plotting routines
3########### SVN repository information ###################
4# $Date: 2012-11-07 20:54:48 +0000 (Wed, 07 Nov 2012) $
5# $Author: vondreele $
6# $Revision: 796 $
7# $URL: trunk/GSASIIplot.py $
8# $Id: GSASIIplot.py 796 2012-11-07 20:54:48Z vondreele $
9########### SVN repository information ###################
10import math
11import time
12import copy
13import os.path
14import numpy as np
15import numpy.ma as ma
16import numpy.linalg as nl
17import wx
18import wx.aui
19import wx.glcanvas
20import matplotlib as mpl
21import mpl_toolkits.mplot3d.axes3d as mp3d
22import GSASIIpath
23GSASIIpath.SetVersionNumber("$Revision: 796 $")
24import GSASIIgrid as G2gd
25import GSASIIimage as G2img
26import GSASIIpwd as G2pwd
27import GSASIIIO as G2IO
28import GSASIIpwdGUI as G2pdG
29import GSASIIimgGUI as G2imG
30import GSASIIphsGUI as G2phG
31import GSASIIlattice as G2lat
32import GSASIIspc as G2spc
33import GSASIImath as G2mth
34import pytexture as ptx
35from  OpenGL.GL import *
36from OpenGL.GLU import *
37from OpenGL.GLUT import *
38from OpenGL.GLE import *
39from matplotlib.backends.backend_wx import _load_bitmap
40from matplotlib.backends.backend_wxagg import FigureCanvasWxAgg as Canvas
41from matplotlib.backends.backend_wxagg import NavigationToolbar2Wx as Toolbar
42
43# useful degree trig functions
44sind = lambda x: math.sin(x*math.pi/180.)
45cosd = lambda x: math.cos(x*math.pi/180.)
46tand = lambda x: math.tan(x*math.pi/180.)
47asind = lambda x: 180.*math.asin(x)/math.pi
48acosd = lambda x: 180.*math.acos(x)/math.pi
49atan2d = lambda x,y: 180.*math.atan2(y,x)/math.pi
50atand = lambda x: 180.*math.atan(x)/math.pi
51# numpy versions
52npsind = lambda x: np.sin(x*np.pi/180.)
53npcosd = lambda x: np.cos(x*np.pi/180.)
54npacosd = lambda x: 180.*np.arccos(x)/np.pi
55npasind = lambda x: 180.*np.arcsin(x)/np.pi
56npatand = lambda x: 180.*np.arctan(x)/np.pi
57npatan2d = lambda x,y: 180.*np.arctan2(x,y)/np.pi
58   
59class G2PlotMpl(wx.Panel):   
60    def __init__(self,parent,id=-1,dpi=None,**kwargs):
61        wx.Panel.__init__(self,parent,id=id,**kwargs)
62        mpl.rcParams['legend.fontsize'] = 10
63        self.figure = mpl.figure.Figure(dpi=dpi,figsize=(5,6))
64        self.canvas = Canvas(self,-1,self.figure)
65        self.toolbar = GSASIItoolbar(self.canvas)
66
67        self.toolbar.Realize()
68       
69        sizer=wx.BoxSizer(wx.VERTICAL)
70        sizer.Add(self.canvas,1,wx.EXPAND)
71        sizer.Add(self.toolbar,0,wx.LEFT|wx.EXPAND)
72        self.SetSizer(sizer)
73       
74class G2PlotOgl(wx.Panel):
75    def __init__(self,parent,id=-1,dpi=None,**kwargs):
76        self.figure = wx.Panel.__init__(self,parent,id=id,**kwargs)
77        self.canvas = wx.glcanvas.GLCanvas(self,-1,**kwargs)
78        self.camera = {}
79        sizer=wx.BoxSizer(wx.VERTICAL)
80        sizer.Add(self.canvas,1,wx.EXPAND)
81        self.SetSizer(sizer)
82       
83class G2Plot3D(wx.Panel):
84    def __init__(self,parent,id=-1,dpi=None,**kwargs):
85        wx.Panel.__init__(self,parent,id=id,**kwargs)
86        self.figure = mpl.figure.Figure(dpi=dpi,figsize=(6,6))
87        self.canvas = Canvas(self,-1,self.figure)
88        self.toolbar = GSASIItoolbar(self.canvas)
89
90        self.toolbar.Realize()
91       
92        sizer=wx.BoxSizer(wx.VERTICAL)
93        sizer.Add(self.canvas,1,wx.EXPAND)
94        sizer.Add(self.toolbar,0,wx.LEFT|wx.EXPAND)
95        self.SetSizer(sizer)
96                             
97class G2PlotNoteBook(wx.Panel):
98    def __init__(self,parent,id=-1):
99        wx.Panel.__init__(self,parent,id=id)
100        #so one can't delete a plot page!!
101        self.nb = wx.aui.AuiNotebook(self, \
102            style=wx.aui.AUI_NB_DEFAULT_STYLE ^ wx.aui.AUI_NB_CLOSE_ON_ACTIVE_TAB)
103        sizer = wx.BoxSizer()
104        sizer.Add(self.nb,1,wx.EXPAND)
105        self.SetSizer(sizer)
106        self.status = parent.CreateStatusBar()
107        self.status.SetFieldsCount(2)
108        self.status.SetStatusWidths([150,-1])
109        self.Bind(wx.aui.EVT_AUINOTEBOOK_PAGE_CHANGED, self.OnPageChanged)
110       
111        self.plotList = []
112           
113    def addMpl(self,name=""):
114        page = G2PlotMpl(self.nb)
115        self.nb.AddPage(page,name)
116       
117        self.plotList.append(name)
118       
119        return page.figure
120       
121    def add3D(self,name=""):
122        page = G2Plot3D(self.nb)
123        self.nb.AddPage(page,name)
124       
125        self.plotList.append(name)
126       
127        return page.figure
128       
129    def addOgl(self,name=""):
130        page = G2PlotOgl(self.nb)
131        self.nb.AddPage(page,name)
132       
133        self.plotList.append(name)
134       
135        return page.figure
136       
137    def Delete(self,name):
138        try:
139            item = self.plotList.index(name)
140            del self.plotList[item]
141            self.nb.DeletePage(item)
142        except ValueError:          #no plot of this name - do nothing
143            return     
144               
145    def clear(self):
146        while self.nb.GetPageCount():
147            self.nb.DeletePage(0)
148        self.plotList = []
149        self.status.DestroyChildren()
150       
151    def Rename(self,oldName,newName):
152        try:
153            item = self.plotList.index(oldName)
154            self.plotList[item] = newName
155            self.nb.SetPageText(item,newName)
156        except ValueError:          #no plot of this name - do nothing
157            return     
158       
159    def OnPageChanged(self,event):       
160        if self.plotList:
161            self.status.SetStatusText('Better to select this from GSAS-II data tree',1)
162        self.status.DestroyChildren()                           #get rid of special stuff on status bar
163       
164class GSASIItoolbar(Toolbar):
165    ON_MPL_HELP = wx.NewId()
166    ON_MPL_KEY = wx.NewId()
167    def __init__(self,plotCanvas):
168        Toolbar.__init__(self,plotCanvas)
169        POSITION_OF_CONFIGURE_SUBPLOTS_BTN = 6
170        self.DeleteToolByPos(POSITION_OF_CONFIGURE_SUBPLOTS_BTN)
171        parent = self.GetParent()
172        key = os.path.join(os.path.split(__file__)[0],'key.ico')
173        self.AddSimpleTool(self.ON_MPL_KEY,_load_bitmap(key),'Key press','Select key press')
174        wx.EVT_TOOL(self,self.ON_MPL_KEY,self.OnKey)
175        help = os.path.join(os.path.split(__file__)[0],'help.ico')
176        self.AddSimpleTool(self.ON_MPL_HELP,_load_bitmap(help),'Help on','Show help on')
177        wx.EVT_TOOL(self,self.ON_MPL_HELP,self.OnHelp)
178    def OnHelp(self,event):
179        Page = self.GetParent().GetParent()
180        pageNo = Page.GetSelection()
181        bookmark = Page.GetPageText(pageNo)
182        bookmark = bookmark.strip(')').replace('(','_')
183        G2gd.ShowHelp(bookmark,self.TopLevelParent)
184    def OnKey(self,event):
185        parent = self.GetParent()
186        if parent.Choice:
187            dlg = wx.SingleChoiceDialog(parent,'Select','Key press',list(parent.Choice))
188            if dlg.ShowModal() == wx.ID_OK:
189                sel = dlg.GetSelection()
190                event.key = parent.Choice[sel][0]
191                parent.keyPress(event)
192            dlg.Destroy()
193   
194
195################################################################################
196##### PlotSngl
197################################################################################
198           
199def PlotSngl(self,newPlot=False):
200    '''Single crystal structure factor plotting package - displays zone of reflections as rings proportional
201        to F, F**2, etc. as requested
202    '''
203    from matplotlib.patches import Circle
204    global HKL,HKLF
205
206    def OnSCMotion(event):
207        xpos = event.xdata
208        if xpos:
209            xpos = round(xpos)                                        #avoid out of frame mouse position
210            ypos = round(event.ydata)
211            zpos = Data['Layer']
212            if '100' in Data['Zone']:
213                HKLtxt = '(%3d,%3d,%3d)'%(zpos,xpos,ypos)
214            elif '010' in Data['Zone']:
215                HKLtxt = '(%3d,%3d,%3d)'%(xpos,zpos,ypos)
216            elif '001' in Data['Zone']:
217                HKLtxt = '(%3d,%3d,%3d)'%(xpos,ypos,zpos)
218            Page.canvas.SetToolTipString(HKLtxt)
219            self.G2plotNB.status.SetFields(['HKL = '+HKLtxt,''])
220               
221    def OnSCPick(event):
222        zpos = Data['Layer']
223        pos = event.artist.center
224        if '100' in Data['Zone']:
225            Page.canvas.SetToolTipString('(picked:(%3d,%3d,%3d))'%(zpos,pos[0],pos[1]))
226            hkl = np.array([zpos,pos[0],pos[1]])
227        elif '010' in Data['Zone']:
228            Page.canvas.SetToolTipString('(picked:(%3d,%3d,%3d))'%(pos[0],zpos,pos[1]))
229            hkl = np.array([pos[0],zpos,pos[1]])
230        elif '001' in Data['Zone']:
231            Page.canvas.SetToolTipString('(picked:(%3d,%3d,%3d))'%(pos[0],pos[1],zpos))
232            hkl = np.array([pos[0],pos[1],zpos])
233        h,k,l = hkl
234        hklf = HKLF[np.where(np.all(HKL-hkl == [0,0,0],axis=1))]
235        if len(hklf):
236            Fosq,sig,Fcsq = hklf[0]
237            HKLtxt = '(%3d,%3d,%3d %.2f %.3f %.2f %.2f)'%(h,k,l,Fosq,sig,Fcsq,(Fosq-Fcsq)/(scale*sig))
238            self.G2plotNB.status.SetFields(['','HKL, Fosq, sig, Fcsq, delFsq/sig = '+HKLtxt])
239                                 
240    def OnSCKeyPress(event):
241        print event.key
242
243    try:
244        plotNum = self.G2plotNB.plotList.index('Structure Factors')
245        Page = self.G2plotNB.nb.GetPage(plotNum)
246        if not newPlot:
247            Plot = Page.figure.gca()          #get previous powder plot & get limits
248            xylim = Plot.get_xlim(),Plot.get_ylim()
249        Page.figure.clf()
250        Page.Choice = None
251        Plot = Page.figure.gca()          #get a fresh plot after clf()
252    except ValueError:
253        Plot = self.G2plotNB.addMpl('Structure Factors').gca()
254        plotNum = self.G2plotNB.plotList.index('Structure Factors')
255        Page = self.G2plotNB.nb.GetPage(plotNum)
256#        Page.canvas.mpl_connect('key_press_event', OnSCKeyPress)
257        Page.canvas.mpl_connect('pick_event', OnSCPick)
258        Page.canvas.mpl_connect('motion_notify_event', OnSCMotion)
259    Page.Choice = None
260    Page.SetFocus()
261   
262    Plot.set_aspect(aspect='equal')
263    HKLref = self.PatternTree.GetItemPyData(self.Sngl)[1]
264    Data = self.PatternTree.GetItemPyData( 
265        G2gd.GetPatternTreeItemId(self,self.Sngl, 'HKL Plot Controls'))
266    Type = Data['Type']           
267    scale = Data['Scale']
268    HKLmax = Data['HKLmax']
269    HKLmin = Data['HKLmin']
270    FosqMax = Data['FoMax']
271    FoMax = math.sqrt(FosqMax)
272    ifFc = Data['ifFc']
273    xlabel = ['k, h=','h, k=','h, l=']
274    ylabel = ['l','l','k']
275    zones = ['100','010','001']
276    pzone = [[1,2],[0,2],[0,1]]
277    izone = zones.index(Data['Zone'])
278    Plot.set_title(self.PatternTree.GetItemText(self.Sngl)[5:])
279    HKL = []
280    HKLF = []
281    for refl in HKLref:
282        H = np.array(refl[:3])
283        sig,Fosq,Fcsq = refl[7:10]
284        HKL.append(H)
285        HKLF.append([Fosq,sig,Fcsq])
286        if H[izone] == Data['Layer']:
287            B = 0
288            if Type == 'Fosq':
289                A = scale*Fosq/FosqMax
290                B = scale*Fcsq/FosqMax
291                C = abs(A-B)
292            elif Type == 'Fo':
293                A = scale*math.sqrt(max(0,Fosq))/FoMax
294                B = scale*math.sqrt(max(0,Fcsq))/FoMax
295                C = abs(A-B)
296            elif Type == '|DFsq|/sig':
297                A = abs(Fosq-Fcsq)/(scale*sig)
298            elif Type == '|DFsq|>sig':
299                A = abs(Fosq-Fcsq)/(scale*sig)
300                if A < 1.0: A = 0                   
301            elif Type == '|DFsq|>3sig':
302                A = abs(Fosq-Fcsq)/(scale*sig)
303                if A < 3.0: A = 0                   
304            xy = (H[pzone[izone][0]],H[pzone[izone][1]])
305            if A > 0.0:
306                Plot.add_artist(Circle(xy,radius=A,ec='g',fc='w',picker=3))
307            if B:
308                Plot.add_artist(Circle(xy,radius=B,ec='b',fc='w'))
309                radius = C
310                if radius > 0:
311                    if A > B:
312                        Plot.add_artist(Circle(xy,radius=radius,ec='g',fc='g'))
313                    else:                   
314                        Plot.add_artist(Circle(xy,radius=radius,ec='r',fc='r'))
315    HKL = np.array(HKL,dtype=np.int)
316    HKLF = np.array(HKLF)
317    Plot.set_xlabel(xlabel[izone]+str(Data['Layer']),fontsize=12)
318    Plot.set_ylabel(ylabel[izone],fontsize=12)
319    Plot.set_xlim((HKLmin[pzone[izone][0]],HKLmax[pzone[izone][0]]))
320    Plot.set_ylim((HKLmin[pzone[izone][1]],HKLmax[pzone[izone][1]]))
321    if not newPlot:
322        Page.toolbar.push_current()
323        Plot.set_xlim(xylim[0])
324        Plot.set_ylim(xylim[1])
325        xylim = []
326        Page.toolbar.push_current()
327        Page.toolbar.draw()
328    else:
329        Page.canvas.draw()
330       
331################################################################################
332##### PlotPatterns
333################################################################################
334           
335def PlotPatterns(G2frame,newPlot=False):
336    '''Powder pattern plotting package - displays single or multiple powder patterns as intensity vs
337    2-theta, q or TOF. Can display multiple patterns as "waterfall plots" or contour plots. Log I
338    plotting available.
339    '''
340    global HKL
341
342#    def OnKeyBox(event):
343#        if G2frame.G2plotNB.nb.GetSelection() == G2frame.G2plotNB.plotList.index('Powder Patterns'):
344#            event.key = cb.GetValue()[0]
345#            cb.SetValue(' key press')
346#            wx.CallAfter(OnPlotKeyPress,event)
347#                       
348    def OnPlotKeyPress(event):
349        newPlot = False
350        if event.key == 'w':
351            if G2frame.Weight:
352                G2frame.Weight = False
353            else:
354                G2frame.Weight = True
355                G2frame.SinglePlot = True
356            newPlot = True
357        elif event.key == 'b':
358            if G2frame.SubBack:
359                G2frame.SubBack = False
360            else:
361                G2frame.SubBack = True
362                G2frame.SinglePlot = True               
363        elif event.key == 'n':
364            if G2frame.Contour:
365                pass
366            else:
367                if G2frame.logPlot:
368                    G2frame.logPlot = False
369                else:
370                    G2frame.Offset[0] = 0
371                    G2frame.logPlot = True
372                newPlot = True
373        elif event.key == 'u':
374            if G2frame.Contour:
375                G2frame.Cmax = min(1.0,G2frame.Cmax*1.2)
376            elif G2frame.logPlot:
377                pass
378            elif G2frame.Offset[0] < 100.:
379                G2frame.Offset[0] += 1.
380        elif event.key == 'd':
381            if G2frame.Contour:
382                G2frame.Cmax = max(0.0,G2frame.Cmax*0.8)
383            elif G2frame.logPlot:
384                pass
385            elif G2frame.Offset[0] > 0.:
386                G2frame.Offset[0] -= 1.
387        elif event.key == 'l':
388            G2frame.Offset[1] -= 1.
389        elif event.key == 'r':
390            G2frame.Offset[1] += 1.
391        elif event.key == 'o':
392            G2frame.Offset = [0,0]
393        elif event.key == 'c':
394            newPlot = True
395            if G2frame.Contour:
396                G2frame.Contour = False
397            else:
398                G2frame.Contour = True
399                G2frame.SinglePlot = False
400                G2frame.Offset = [0.,0.]
401        elif event.key == 'q':
402            newPlot = True
403            if G2frame.qPlot:
404                G2frame.qPlot = False
405            else:
406                G2frame.qPlot = True
407        elif event.key == 's':
408            if G2frame.Contour:
409                choice = [m for m in mpl.cm.datad.keys() if not m.endswith("_r")]
410                choice.sort()
411                dlg = wx.SingleChoiceDialog(G2frame,'Select','Color scheme',choice)
412                if dlg.ShowModal() == wx.ID_OK:
413                    sel = dlg.GetSelection()
414                    G2frame.ContourColor = choice[sel]
415                else:
416                    G2frame.ContourColor = 'Paired'
417                dlg.Destroy()
418            else:               
419                if G2frame.SinglePlot:
420                    G2frame.SinglePlot = False
421                else:
422                    G2frame.SinglePlot = True
423            newPlot = True
424        elif event.key == '+':
425            if G2frame.PickId:
426                G2frame.PickId = False
427        elif event.key == 'i':                  #for smoothing contour plot
428            choice = ['nearest','bilinear','bicubic','spline16','spline36','hanning',
429               'hamming','hermite','kaiser','quadric','catrom','gaussian','bessel',
430               'mitchell','sinc','lanczos']
431            dlg = wx.SingleChoiceDialog(G2frame,'Select','Interpolation',choice)
432            if dlg.ShowModal() == wx.ID_OK:
433                sel = dlg.GetSelection()
434                G2frame.Interpolate = choice[sel]
435            else:
436                G2frame.Interpolate = 'nearest'
437            dlg.Destroy()
438           
439        PlotPatterns(G2frame,newPlot=newPlot)
440       
441    def OnMotion(event):
442        xpos = event.xdata
443        if xpos:                                        #avoid out of frame mouse position
444            ypos = event.ydata
445            Page.canvas.SetCursor(wx.CROSS_CURSOR)
446            try:
447                Parms,Parms2 = G2frame.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(G2frame,G2frame.PatternId, 'Instrument Parameters'))
448                if 'C' in Parms['Type'][0]:
449                    wave = G2mth.getWave(Parms)
450                    if G2frame.qPlot:
451                        try:
452                            xpos = 2.0*asind(xpos*wave/(4*math.pi))
453                        except ValueError:      #avoid bad value in asin beyond upper limit
454                            pass
455                    dsp = 0.0
456                    if abs(xpos) > 0.:                  #avoid possible singularity at beam center
457                        dsp = wave/(2.*sind(abs(xpos)/2.0))
458                    if G2frame.Contour:
459                        G2frame.G2plotNB.status.SetStatusText('2-theta =%9.3f d =%9.5f pattern ID =%5d'%(xpos,dsp,int(ypos)),1)
460                    else:
461                        G2frame.G2plotNB.status.SetStatusText('2-theta =%9.3f d =%9.5f Intensity =%9.1f'%(xpos,dsp,ypos),1)
462                else:       #TOF neutrons
463                    dsp = 0.0
464                    difC = Parms['difC'][1]
465                    dsp = xpos/difC             #rough approx.!
466                    if G2frame.Contour:
467                        G2frame.G2plotNB.status.SetStatusText('TOF =%9.3f d =%9.5f pattern ID =%5d'%(xpos,dsp,int(ypos)),1)
468                    else:
469                        G2frame.G2plotNB.status.SetStatusText('TOF =%9.3f d =%9.5f Intensity =%9.1f'%(xpos,dsp,ypos),1)
470                if G2frame.itemPicked:
471                    Page.canvas.SetToolTipString('%9.3f'%(xpos))
472                if G2frame.PickId:
473                    found = []
474                    if G2frame.PatternTree.GetItemText(G2frame.PickId) in ['Index Peak List','Unit Cells List','Reflection Lists'] or \
475                        'PWDR' in G2frame.PatternTree.GetItemText(PickId):
476                        if len(HKL):
477                            view = Page.toolbar._views.forward()[0][:2]
478                            wid = view[1]-view[0]
479                            found = HKL[np.where(np.fabs(HKL.T[5]-xpos) < 0.002*wid)]
480                        if len(found):
481                            h,k,l = found[0][:3] 
482                            Page.canvas.SetToolTipString('%d,%d,%d'%(int(h),int(k),int(l)))
483                        else:
484                            Page.canvas.SetToolTipString('')
485
486            except TypeError:
487                G2frame.G2plotNB.status.SetStatusText('Select PWDR powder pattern first',1)
488                                                   
489    def OnPick(event):
490        if G2frame.itemPicked is not None: return
491        PatternId = G2frame.PatternId
492        try:
493            Parms,Parms2 = G2frame.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(G2frame,G2frame.PatternId, 'Instrument Parameters'))
494        except TypeError:
495            return
496        if 'C' in Parms['Type'][0]:
497            wave = G2mth.getWave(Parms)
498        else:
499            difC = Parms['difC'][1]
500        PickId = G2frame.PickId
501        pick = event.artist
502        mouse = event.mouseevent       
503        xpos = pick.get_xdata()
504        ypos = pick.get_ydata()
505        ind = event.ind
506        xy = list(zip(np.take(xpos,ind),np.take(ypos,ind))[0])
507        if G2frame.PatternTree.GetItemText(PickId) == 'Peak List':
508            if ind.all() != [0]:                                    #picked a data point
509                data = G2frame.PatternTree.GetItemPyData(G2frame.PickId)
510                XY = G2mth.setPeakparms(Parms,Parms2,xy[0],xy[1])
511                data.append(XY)
512                G2pdG.UpdatePeakGrid(G2frame,data)
513                PlotPatterns(G2frame)
514            else:                                                   #picked a peak list line
515                G2frame.itemPicked = pick
516        elif G2frame.PatternTree.GetItemText(PickId) == 'Limits':
517            if ind.all() != [0]:                                    #picked a data point
518                LimitId = G2gd.GetPatternTreeItemId(G2frame,PatternId, 'Limits')
519                data = G2frame.PatternTree.GetItemPyData(LimitId)
520                if 'C' in Parms['Type'][0]:                            #CW data - TOF later in an elif
521                    if G2frame.qPlot:                              #qplot - convert back to 2-theta
522                        xy[0] = 2.0*asind(xy[0]*wave/(4*math.pi))
523                if mouse.button==1:
524                    data[1][0] = min(xy[0],data[1][1])
525                if mouse.button==3:
526                    data[1][1] = max(xy[0],data[1][0])
527                G2frame.PatternTree.SetItemPyData(LimitId,data)
528                G2pdG.UpdateLimitsGrid(G2frame,data)
529                PlotPatterns(G2frame)
530            else:                                                   #picked a limit line
531                G2frame.itemPicked = pick
532        elif G2frame.PatternTree.GetItemText(PickId) == 'Reflection Lists' or \
533            'PWDR' in G2frame.PatternTree.GetItemText(PickId):
534            G2frame.itemPicked = pick
535            pick = str(pick)
536       
537    def OnRelease(event):
538        if G2frame.itemPicked is None: return
539        Parms,Parms2 = G2frame.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(G2frame,G2frame.PatternId, 'Instrument Parameters'))
540        if 'C' in Parms['Type'][0]:
541            wave = G2mth.getWave(Parms)
542        else:
543            difC = Parms['difC'][1]
544        xpos = event.xdata
545        PickId = G2frame.PickId
546        if G2frame.PatternTree.GetItemText(PickId) in ['Peak List','Limits'] and xpos:
547            lines = []
548            for line in G2frame.Lines: 
549                lines.append(line.get_xdata()[0])
550#            print G2frame.itemPicked.get_xdata()
551            lineNo = lines.index(G2frame.itemPicked.get_xdata()[0])
552            if  lineNo in [0,1]:
553                LimitId = G2gd.GetPatternTreeItemId(G2frame,G2frame.PatternId, 'Limits')
554                data = G2frame.PatternTree.GetItemPyData(LimitId)
555                if G2frame.qPlot:
556                    if 'C' in Parms['Type'][0]:
557                        data[1][lineNo] = 2.0*asind(wave*xpos/(4*math.pi))
558                    else:
559                        data[1][lineNo] = 2*math.pi*Parms['difC'][1]/xpos
560                else:
561                    data[1][lineNo] = xpos
562                G2frame.PatternTree.SetItemPyData(LimitId,data)
563                if G2frame.PatternTree.GetItemText(G2frame.PickId) == 'Limits':
564                    G2pdG.UpdateLimitsGrid(G2frame,data)
565            else:
566                PeakId = G2gd.GetPatternTreeItemId(G2frame,G2frame.PatternId, 'Peak List')
567                data = G2frame.PatternTree.GetItemPyData(PeakId)
568#                print 'peaks',xpos
569                if event.button == 3:
570                    del data[lineNo-2]
571                else:
572                    if G2frame.qPlot:
573                        data[lineNo-2][0] = 2.0*asind(wave*xpos/(4*math.pi))
574                    else:
575                        data[lineNo-2][0] = xpos
576                G2frame.PatternTree.SetItemPyData(PeakId,data)
577                G2pdG.UpdatePeakGrid(G2frame,data)
578        elif (G2frame.PatternTree.GetItemText(PickId) == 'Reflection Lists' or \
579            'PWDR' in G2frame.PatternTree.GetItemText(PickId)) and xpos:
580            Phases = G2frame.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(G2frame,PatternId,'Reflection Lists'))
581            pick = str(G2frame.itemPicked).split('(')[1].strip(')')
582            if 'line' not in pick:       #avoid data points, etc.
583                num = Phases.keys().index(pick)
584                if num:
585                    G2frame.refDelt = -(event.ydata-G2frame.refOffset)/(num*Ymax)
586                else:       #1st row of refl ticks
587                    G2frame.refOffset = event.ydata
588        PlotPatterns(G2frame)
589        G2frame.itemPicked = None   
590
591    xylim = []
592    try:
593        plotNum = G2frame.G2plotNB.plotList.index('Powder Patterns')
594        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
595        if not newPlot:
596            Plot = Page.figure.gca()          #get previous powder plot & get limits
597            xylim = Plot.get_xlim(),Plot.get_ylim()
598        Page.figure.clf()
599        Plot = Page.figure.gca()          #get a fresh plot after clf()
600    except ValueError:
601        newPlot = True
602        G2frame.Cmax = 1.0
603        Plot = G2frame.G2plotNB.addMpl('Powder Patterns').gca()
604        plotNum = G2frame.G2plotNB.plotList.index('Powder Patterns')
605        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
606        Page.canvas.mpl_connect('key_press_event', OnPlotKeyPress)
607        Page.canvas.mpl_connect('motion_notify_event', OnMotion)
608        Page.canvas.mpl_connect('pick_event', OnPick)
609        Page.canvas.mpl_connect('button_release_event', OnRelease)
610    Page.SetFocus()
611    G2frame.G2plotNB.status.DestroyChildren()
612    if G2frame.Contour:
613        Page.Choice = (' key press','d: lower contour max','u: raise contour max',
614            'i: interpolation method','s: color scheme','c: contour off')
615    else:
616        if G2frame.logPlot:
617            Page.Choice = (' key press','n: log(I) off','l: offset left','r: offset right',
618                'c: contour on','q: toggle q plot','s: toggle single plot','+: no selection')
619        else:
620            Page.Choice = (' key press','l: offset left','r: offset right','d: offset down',
621                'u: offset up','o: reset offset','b: toggle subtract background','n: log(I) on','c: contour on',
622                'q: toggle q plot','s: toggle single plot','w: toggle divide by sig','+: no selection')
623    Page.keyPress = OnPlotKeyPress
624   
625    PickId = G2frame.PickId
626    PatternId = G2frame.PatternId
627    colors=['b','g','r','c','m','k']
628    Lines = []
629    if G2frame.SinglePlot:
630        Pattern = G2frame.PatternTree.GetItemPyData(PatternId)
631        Pattern.append(G2frame.PatternTree.GetItemText(PatternId))
632        PlotList = [Pattern,]
633        Parms,Parms2 = G2frame.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(G2frame,
634            G2frame.PatternId, 'Instrument Parameters'))
635        ParmList = [Parms,]
636    else:       
637        PlotList = []
638        ParmList = []
639        item, cookie = G2frame.PatternTree.GetFirstChild(G2frame.root)
640        while item:
641            if 'PWDR' in G2frame.PatternTree.GetItemText(item):
642                Pattern = G2frame.PatternTree.GetItemPyData(item)
643                if len(Pattern) < 3:                    # put name on end if needed
644                    Pattern.append(G2frame.PatternTree.GetItemText(item))
645                PlotList.append(Pattern)
646                ParmList.append(G2frame.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(G2frame,
647                    item,'Instrument Parameters'))[0])
648            item, cookie = G2frame.PatternTree.GetNextChild(G2frame.root, cookie)               
649    Ymax = 1.0
650    lenX = 0
651    if PickId:
652        if G2frame.PatternTree.GetItemText(PickId) in ['Reflection Lists']:
653            Phases = G2frame.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(G2frame,PatternId,'Reflection Lists'))
654            HKL = []
655            if Phases:
656                for peak in Phases[G2frame.RefList]:
657                    HKL.append(peak[:6])
658                HKL = np.array(HKL)
659        else:
660            HKL = np.array(G2frame.HKL)
661    for Pattern in PlotList:
662        xye = Pattern[1]
663        Ymax = max(Ymax,max(xye[1]))
664    offset = G2frame.Offset[0]*Ymax/100.0
665    Title = 'Powder Patterns: '+os.path.split(G2frame.GSASprojectfile)[1]
666    if G2frame.logPlot:
667        Title = 'log('+Title+')'
668    Plot.set_title(Title)
669    if G2frame.qPlot:
670        Plot.set_xlabel(r'$q, \AA^{-1}$',fontsize=14)
671    else:
672        if 'C' in ParmList[0]['Type'][0]:       
673            Plot.set_xlabel(r'$\mathsf{2\theta}$',fontsize=14)
674        else:
675            Plot.set_xlabel(r'TOF, $\mathsf{\mu}$s',fontsize=14)           
676    if G2frame.Weight:
677        Plot.set_ylabel(r'$\mathsf{I/\sigma(I)}$',fontsize=14)
678    else:
679        if 'C' in ParmList[0]['Type'][0]:
680            Plot.set_ylabel('Intensity',fontsize=14)
681        else:
682            Plot.set_ylabel('Normalized intensity',fontsize=14)
683    if G2frame.Contour:
684        ContourZ = []
685        ContourY = []
686        Nseq = 0
687    if len(PlotList) < 2:
688        G2frame.Contour = False
689    for N,Pattern in enumerate(PlotList):
690        Parms = ParmList[N]
691        if 'C' in Parms['Type'][0]:
692            wave = G2mth.getWave(Parms)
693        else:
694            difC = Parms['difC'][1]
695        ifpicked = False
696        LimitId = 0
697        xye = np.array(Pattern[1])
698        if PickId:
699            ifpicked = Pattern[2] == G2frame.PatternTree.GetItemText(PatternId)
700            LimitId = G2gd.GetPatternTreeItemId(G2frame,PatternId, 'Limits')
701        if G2frame.qPlot:
702            Id = G2gd.GetPatternTreeItemId(G2frame,G2frame.root, Pattern[2])
703            if 'C' in Parms['Type'][0]:
704                X = 4*np.pi*npsind(xye[0]/2.0)/wave
705            else:
706                X = 2*np.pi*Parms['difC'][1]/xye[0]
707        else:
708            X = xye[0]
709        if not lenX:
710            lenX = len(X)           
711        Y = xye[1]+offset*N
712        if LimitId:
713            limits = np.array(G2frame.PatternTree.GetItemPyData(LimitId))
714            if G2frame.qPlot:
715                if 'C' in Parms['Type'][0]:
716                    limits = 4*np.pi*npsind(limits/2.0)/wave
717                else:
718                    limits = 2*np.pi*difC/limits
719            Lines.append(Plot.axvline(limits[1][0],color='g',dashes=(5,5),picker=3.))   
720            Lines.append(Plot.axvline(limits[1][1],color='r',dashes=(5,5),picker=3.))                   
721        if G2frame.Contour:
722           
723            if lenX == len(X):
724                ContourY.append(N)
725                ContourZ.append(Y)
726                ContourX = X
727                Nseq += 1
728                Plot.set_ylabel('Data sequence',fontsize=12)
729        else:
730            X += G2frame.Offset[1]*.005*N
731            if ifpicked:
732                Z = xye[3]+offset*N
733                W = xye[4]+offset*N
734                D = xye[5]-Ymax*G2frame.delOffset
735                if G2frame.logPlot:
736                    Plot.semilogy(X,Y,colors[N%6]+'+',picker=3.,clip_on=False,nonposy='mask')
737                    Plot.semilogy(X,Z,colors[(N+1)%6],picker=False,nonposy='mask')
738                    Plot.semilogy(X,W,colors[(N+2)%6],picker=False,nonposy='mask')
739                elif G2frame.Weight:
740                    DY = xye[1]*np.sqrt(xye[2])
741                    DYmax = max(DY)
742                    DZ = xye[3]*np.sqrt(xye[2])
743                    DS = xye[5]*np.sqrt(xye[2])-DYmax*G2frame.delOffset
744                    Plot.plot(X,DY,colors[N%6]+'+',picker=3.,clip_on=False)
745                    Plot.plot(X,DZ,colors[(N+1)%6],picker=False)
746                    Plot.plot(X,DS,colors[(N+3)%6],picker=False)
747                    Plot.axhline(0.,color=wx.BLACK)
748                else:
749                    if G2frame.SubBack:
750                        Plot.plot(X,Y-W,colors[N%6]+'+',picker=3.,clip_on=False)
751                        Plot.plot(X,Z-W,colors[(N+1)%6],picker=False)
752                    else:
753                        Plot.plot(X,Y,colors[N%6]+'+',picker=3.,clip_on=False)
754                        Plot.plot(X,Z,colors[(N+1)%6],picker=False)
755                    Plot.plot(X,W,colors[(N+2)%6],picker=False)
756                    Plot.plot(X,D,colors[(N+3)%6],picker=False)
757                    Plot.axhline(0.,color=wx.BLACK)
758                Page.canvas.SetToolTipString('')
759                if G2frame.PatternTree.GetItemText(PickId) == 'Peak List':
760                    tip = 'On data point: Pick peak - L or R MB. On line: L-move, R-delete'
761                    Page.canvas.SetToolTipString(tip)
762                    data = G2frame.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(G2frame,PatternId, 'Peak List'))
763                    for item in data:
764                        if G2frame.qPlot:
765                            if 'C' in Parms['Type'][0]:
766                                Lines.append(Plot.axvline(4*math.pi*sind(item[0]/2.)/wave,color=colors[N%6],picker=2.))
767                            else:
768                                Lines.append(Plot.axvline(2*math.pi*difC/item[0],color=colors[N%6],picker=2.))                               
769                        else:
770                            Lines.append(Plot.axvline(item[0],color=colors[N%6],picker=2.))
771                if G2frame.PatternTree.GetItemText(PickId) == 'Limits':
772                    tip = 'On data point: Lower limit - L MB; Upper limit - R MB. On limit: MB down to move'
773                    Page.canvas.SetToolTipString(tip)
774                    data = G2frame.LimitsTable.GetData()
775            else:
776                if G2frame.logPlot:
777                    Plot.semilogy(X,Y,colors[N%6],picker=False,nonposy='clip')
778                else:
779                    Plot.plot(X,Y,colors[N%6],picker=False)
780    if PickId and not G2frame.Contour:
781        Parms,Parms2 = G2frame.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(G2frame,PatternId, 'Instrument Parameters'))
782        if 'C' in Parms['Type'][0]:
783            wave = G2mth.getWave(Parms)
784        else:
785            difC = Parms['difC'][1]
786        if G2frame.PatternTree.GetItemText(PickId) in ['Index Peak List','Unit Cells List']:
787            peaks = np.array((G2frame.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(G2frame,PatternId, 'Index Peak List'))))
788            for peak in peaks:
789                if G2frame.qPlot:
790                    Plot.axvline(4*np.pi*sind(peak[0]/2.0)/wave,color='b')
791                else:
792                    Plot.axvline(peak[0],color='b')
793            for hkl in G2frame.HKL:
794                if G2frame.qPlot:
795                    Plot.axvline(4*np.pi*sind(hkl[5]/2.0)/wave,color='r',dashes=(5,5))
796                else:
797                    Plot.axvline(hkl[5],color='r',dashes=(5,5))
798        elif G2frame.PatternTree.GetItemText(PickId) in ['Reflection Lists'] or \
799            'PWDR' in G2frame.PatternTree.GetItemText(PickId):
800            Phases = G2frame.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(G2frame,PatternId,'Reflection Lists'))
801            for pId,phase in enumerate(Phases):
802                peaks = Phases[phase]
803                if not peaks:
804                    continue
805                peak = np.array([[peak[4],peak[5]] for peak in peaks])
806                pos = G2frame.refOffset-pId*Ymax*G2frame.refDelt*np.ones_like(peak)
807                if G2frame.qPlot:
808                    Plot.plot(2*np.pi/peak.T[0],pos,colors[pId%6]+'|',mew=1,ms=8,picker=3.,label=phase)
809                else:
810                    Plot.plot(peak.T[1],pos,colors[pId%6]+'|',mew=1,ms=8,picker=3.,label=phase)
811            if len(Phases):
812                handles,legends = Plot.get_legend_handles_labels()  #got double entries in the legends for some reason
813                if handles:
814                    Plot.legend(handles[::2],legends[::2],title='Phases',loc='best')    #skip every other one
815           
816    if G2frame.Contour:
817        acolor = mpl.cm.get_cmap(G2frame.ContourColor)
818        Img = Plot.imshow(ContourZ,cmap=acolor,vmin=0,vmax=Ymax*G2frame.Cmax,interpolation=G2frame.Interpolate, 
819            extent=[ContourX[0],ContourX[-1],ContourY[0],ContourY[-1]],aspect='auto',origin='lower')
820        Page.figure.colorbar(Img)
821    else:
822        G2frame.Lines = Lines
823    if not newPlot:
824        Page.toolbar.push_current()
825        Plot.set_xlim(xylim[0])
826        Plot.set_ylim(xylim[1])
827        xylim = []
828        Page.toolbar.push_current()
829        Page.toolbar.draw()
830    else:
831        Page.canvas.draw()
832#    G2frame.Pwdr = True
833   
834################################################################################
835##### PlotDeltSig
836################################################################################
837           
838def PlotDeltSig(G2frame,kind):
839    try:
840        plotNum = G2frame.G2plotNB.plotList.index('Error analysis')
841        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
842        Page.figure.clf()
843        Plot = Page.figure.gca()          #get a fresh plot after clf()
844    except ValueError:
845        newPlot = True
846        G2frame.Cmax = 1.0
847        Plot = G2frame.G2plotNB.addMpl('Error analysis').gca()
848        plotNum = G2frame.G2plotNB.plotList.index('Error analysis')
849        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
850    Page.Choice = None
851    PatternId = G2frame.PatternId
852    Pattern = G2frame.PatternTree.GetItemPyData(PatternId)
853    Pattern.append(G2frame.PatternTree.GetItemText(PatternId))
854    wtFactor = Pattern[0]['wtFactor']
855    if kind == 'PWDR':
856        limits = G2frame.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(G2frame,PatternId, 'Limits'))[1]
857        xye = np.array(Pattern[1])
858        xmin = np.searchsorted(xye[0],limits[0])
859        xmax = np.searchsorted(xye[0],limits[1])
860        DS = xye[5][xmin:xmax]*np.sqrt(wtFactor*xye[2][xmin:xmax])
861    elif kind == 'HKLF':
862        refl = Pattern[1]
863        DS = []
864        for ref in refl:
865            if ref[6] > 0.:
866                DS.append((ref[5]-ref[7])/ref[6])
867    Page.SetFocus()
868    G2frame.G2plotNB.status.DestroyChildren()
869    DS.sort()
870    EDS = np.zeros_like(DS)
871    DX = np.linspace(0.,1.,num=len(DS),endpoint=True)
872    oldErr = np.seterr(invalid='ignore')    #avoid problem at DX==0
873    T = np.sqrt(np.log(1.0/DX**2))
874    top = 2.515517+0.802853*T+0.010328*T**2
875    bot = 1.0+1.432788*T+0.189269*T**2+0.001308*T**3
876    EDS = np.where(DX>0,-(T-top/bot),(T-top/bot))
877    low1 = np.searchsorted(EDS,-1.)
878    hi1 = np.searchsorted(EDS,1.)
879    slp,intcp = np.polyfit(EDS[low1:hi1],DS[low1:hi1],deg=1)
880    frac = 100.*(hi1-low1)/len(DS)
881    G2frame.G2plotNB.status.SetStatusText(  \
882        'Over range -1. to 1. :'+' slope = %.3f, intercept = %.3f for %.2f%% of the fitted data'%(slp,intcp,frac),1)
883    Plot.set_title('Normal probability for '+Pattern[-1])
884    Plot.set_xlabel(r'expected $\mathsf{\Delta/\sigma}$',fontsize=14)
885    Plot.set_ylabel(r'observed $\mathsf{\Delta/\sigma}$',fontsize=14)
886    Plot.plot(EDS,DS,'r+',label='result')
887    Plot.plot([-2,2],[-2,2],'k',dashes=(5,5),label='ideal')
888    Plot.legend(loc='upper left')
889    np.seterr(invalid='warn')
890    Page.canvas.draw()
891       
892################################################################################
893##### PlotISFG
894################################################################################
895           
896def PlotISFG(G2frame,newPlot=False,type=''):
897    ''' PLotting package for PDF analysis; displays I(q), S(q), F(q) and G(r) as single
898    or multiple plots with waterfall and contour plots as options
899    '''
900    if not type:
901        type = G2frame.G2plotNB.plotList[G2frame.G2plotNB.nb.GetSelection()]
902    if type not in ['I(Q)','S(Q)','F(Q)','G(R)']:
903        return
904    superMinusOne = unichr(0xaf)+unichr(0xb9)
905   
906    def OnPlotKeyPress(event):
907        newPlot = False
908        if event.key == 'u':
909            if G2frame.Contour:
910                G2frame.Cmax = min(1.0,G2frame.Cmax*1.2)
911            elif G2frame.Offset[0] < 100.:
912                G2frame.Offset[0] += 1.
913        elif event.key == 'd':
914            if G2frame.Contour:
915                G2frame.Cmax = max(0.0,G2frame.Cmax*0.8)
916            elif G2frame.Offset[0] > 0.:
917                G2frame.Offset[0] -= 1.
918        elif event.key == 'l':
919            G2frame.Offset[1] -= 1.
920        elif event.key == 'r':
921            G2frame.Offset[1] += 1.
922        elif event.key == 'o':
923            G2frame.Offset = [0,0]
924        elif event.key == 'c':
925            newPlot = True
926            if G2frame.Contour:
927                G2frame.Contour = False
928            else:
929                G2frame.Contour = True
930                G2frame.SinglePlot = False
931                G2frame.Offset = [0.,0.]
932        elif event.key == 's':
933            if G2frame.Contour:
934                choice = [m for m in mpl.cm.datad.keys() if not m.endswith("_r")]
935                choice.sort()
936                dlg = wx.SingleChoiceDialog(G2frame,'Select','Color scheme',choice)
937                if dlg.ShowModal() == wx.ID_OK:
938                    sel = dlg.GetSelection()
939                    G2frame.ContourColor = choice[sel]
940                else:
941                    G2frame.ContourColor = 'Paired'
942                dlg.Destroy()
943            else:               
944                if G2frame.SinglePlot:
945                    G2frame.SinglePlot = False
946                else:
947                    G2frame.SinglePlot = True
948        elif event.key == 'i':                  #for smoothing contour plot
949            choice = ['nearest','bilinear','bicubic','spline16','spline36','hanning',
950               'hamming','hermite','kaiser','quadric','catrom','gaussian','bessel',
951               'mitchell','sinc','lanczos']
952            dlg = wx.SingleChoiceDialog(G2frame,'Select','Interpolation',choice)
953            if dlg.ShowModal() == wx.ID_OK:
954                sel = dlg.GetSelection()
955                G2frame.Interpolate = choice[sel]
956            else:
957                G2frame.Interpolate = 'nearest'
958            dlg.Destroy()
959        elif event.key == 't' and not G2frame.Contour:
960            if G2frame.Legend:
961                G2frame.Legend = False
962            else:
963                G2frame.Legend = True
964           
965           
966        PlotISFG(G2frame,newPlot=newPlot,type=type)
967       
968    def OnKeyBox(event):
969        if G2frame.G2plotNB.nb.GetSelection() == G2frame.G2plotNB.plotList.index(type):
970            event.key = cb.GetValue()[0]
971            cb.SetValue(' key press')
972            wx.CallAfter(OnPlotKeyPress,event)
973                       
974    def OnMotion(event):
975        xpos = event.xdata
976        if xpos:                                        #avoid out of frame mouse position
977            ypos = event.ydata
978            Page.canvas.SetCursor(wx.CROSS_CURSOR)
979            try:
980                if G2frame.Contour:
981                    G2frame.G2plotNB.status.SetStatusText('R =%.3fA pattern ID =%5d'%(xpos,int(ypos)),1)
982                else:
983                    G2frame.G2plotNB.status.SetStatusText('R =%.3fA %s =%.2f'%(xpos,type,ypos),1)                   
984            except TypeError:
985                G2frame.G2plotNB.status.SetStatusText('Select '+type+' pattern first',1)
986   
987    xylim = []
988    try:
989        plotNum = G2frame.G2plotNB.plotList.index(type)
990        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
991        if not newPlot:
992            Plot = Page.figure.gca()          #get previous plot & get limits
993            xylim = Plot.get_xlim(),Plot.get_ylim()
994        Page.figure.clf()
995        Plot = Page.figure.gca()
996    except ValueError:
997        newPlot = True
998        G2frame.Cmax = 1.0
999        Plot = G2frame.G2plotNB.addMpl(type).gca()
1000        plotNum = G2frame.G2plotNB.plotList.index(type)
1001        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
1002        Page.canvas.mpl_connect('key_press_event', OnPlotKeyPress)
1003        Page.canvas.mpl_connect('motion_notify_event', OnMotion)
1004   
1005    Page.SetFocus()
1006    G2frame.G2plotNB.status.DestroyChildren()
1007    if G2frame.Contour:
1008        Page.Choice = (' key press','d: lower contour max','u: raise contour max',
1009            'i: interpolation method','s: color scheme','c: contour off')
1010    else:
1011        Page.Choice = (' key press','l: offset left','r: offset right','d: offset down','u: offset up',
1012            'o: reset offset','t: toggle legend','c: contour on','s: toggle single plot')
1013    Page.keyPress = OnPlotKeyPress
1014    PatternId = G2frame.PatternId
1015    PickId = G2frame.PickId
1016    Plot.set_title(type)
1017    if type == 'G(R)':
1018        Plot.set_xlabel(r'$R,\AA$',fontsize=14)
1019    else:
1020        Plot.set_xlabel(r'$Q,\AA$'+superMinusOne,fontsize=14)
1021    Plot.set_ylabel(r''+type,fontsize=14)
1022    colors=['b','g','r','c','m','k']
1023    name = G2frame.PatternTree.GetItemText(PatternId)[4:]
1024    Pattern = []   
1025    if G2frame.SinglePlot:
1026        name = G2frame.PatternTree.GetItemText(PatternId)
1027        name = type+name[4:]
1028        Id = G2gd.GetPatternTreeItemId(G2frame,PatternId,name)
1029        Pattern = G2frame.PatternTree.GetItemPyData(Id)
1030        if Pattern:
1031            Pattern.append(name)
1032        PlotList = [Pattern,]
1033    else:
1034        PlotList = []
1035        item, cookie = G2frame.PatternTree.GetFirstChild(G2frame.root)
1036        while item:
1037            if 'PDF' in G2frame.PatternTree.GetItemText(item):
1038                name = type+G2frame.PatternTree.GetItemText(item)[4:]
1039                Id = G2gd.GetPatternTreeItemId(G2frame,item,name)
1040                Pattern = G2frame.PatternTree.GetItemPyData(Id)
1041                if Pattern:
1042                    Pattern.append(name)
1043                    PlotList.append(Pattern)
1044            item, cookie = G2frame.PatternTree.GetNextChild(G2frame.root, cookie)
1045    PDFdata = G2frame.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(G2frame,PatternId, 'PDF Controls'))
1046    numbDen = G2pwd.GetNumDensity(PDFdata['ElList'],PDFdata['Form Vol'])
1047    Xb = [0.,10.]
1048    Yb = [0.,-40.*np.pi*numbDen]
1049    Ymax = 1.0
1050    lenX = 0
1051    for Pattern in PlotList:
1052        xye = Pattern[1]
1053        Ymax = max(Ymax,max(xye[1]))
1054    offset = G2frame.Offset[0]*Ymax/100.0
1055    if G2frame.Contour:
1056        ContourZ = []
1057        ContourY = []
1058        Nseq = 0
1059    for N,Pattern in enumerate(PlotList):
1060        xye = Pattern[1]
1061        if PickId:
1062            ifpicked = Pattern[2] == G2frame.PatternTree.GetItemText(PatternId)
1063        X = xye[0]
1064        if not lenX:
1065            lenX = len(X)           
1066        Y = xye[1]+offset*N
1067        if G2frame.Contour:
1068            if lenX == len(X):
1069                ContourY.append(N)
1070                ContourZ.append(Y)
1071                ContourX = X
1072                Nseq += 1
1073                Plot.set_ylabel('Data sequence',fontsize=12)
1074        else:
1075            X = xye[0]+G2frame.Offset[1]*.005*N
1076            if ifpicked:
1077                Plot.plot(X,Y,colors[N%6]+'+',picker=3.,clip_on=False)
1078                Page.canvas.SetToolTipString('')
1079            else:
1080                if G2frame.Legend:
1081                    Plot.plot(X,Y,colors[N%6],picker=False,label='Azm:'+Pattern[2].split('=')[1])
1082                else:
1083                    Plot.plot(X,Y,colors[N%6],picker=False)
1084            if type == 'G(R)':
1085                Plot.plot(Xb,Yb,color='k',dashes=(5,5))
1086            elif type == 'F(Q)':
1087                Plot.axhline(0.,color=wx.BLACK)
1088            elif type == 'S(Q)':
1089                Plot.axhline(1.,color=wx.BLACK)
1090    if G2frame.Contour:
1091        acolor = mpl.cm.get_cmap(G2frame.ContourColor)
1092        Img = Plot.imshow(ContourZ,cmap=acolor,vmin=0,vmax=Ymax*G2frame.Cmax,interpolation=G2frame.Interpolate, 
1093            extent=[ContourX[0],ContourX[-1],ContourY[0],ContourY[-1]],aspect='auto',origin='lower')
1094        Page.figure.colorbar(Img)
1095    elif G2frame.Legend:
1096        Plot.legend(loc='best')
1097    if not newPlot:
1098        Page.toolbar.push_current()
1099        Plot.set_xlim(xylim[0])
1100        Plot.set_ylim(xylim[1])
1101        xylim = []
1102        Page.toolbar.push_current()
1103        Page.toolbar.draw()
1104    else:
1105        Page.canvas.draw()
1106       
1107################################################################################
1108##### PlotXY
1109################################################################################
1110           
1111def PlotXY(G2frame,XY,newPlot=False,type=''):
1112    '''simple plot of xy data, used for diagnostic purposes
1113    '''
1114    def OnMotion(event):
1115        xpos = event.xdata
1116        if xpos:                                        #avoid out of frame mouse position
1117            ypos = event.ydata
1118            Page.canvas.SetCursor(wx.CROSS_CURSOR)
1119            try:
1120                G2frame.G2plotNB.status.SetStatusText('X =%9.3f %s =%9.3f'%(xpos,type,ypos),1)                   
1121            except TypeError:
1122                G2frame.G2plotNB.status.SetStatusText('Select '+type+' pattern first',1)
1123
1124    try:
1125        plotNum = G2frame.G2plotNB.plotList.index(type)
1126        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
1127        if not newPlot:
1128            Plot = Page.figure.gca()
1129            xylim = Plot.get_xlim(),Plot.get_ylim()
1130        Page.figure.clf()
1131        Plot = Page.figure.gca()
1132    except ValueError:
1133        newPlot = True
1134        Plot = G2frame.G2plotNB.addMpl(type).gca()
1135        plotNum = G2frame.G2plotNB.plotList.index(type)
1136        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
1137        Page.canvas.mpl_connect('motion_notify_event', OnMotion)
1138   
1139    Page.Choice = None
1140    Page.SetFocus()
1141    G2frame.G2plotNB.status.DestroyChildren()
1142    Plot.set_title(type)
1143    Plot.set_xlabel(r'X',fontsize=14)
1144    Plot.set_ylabel(r''+type,fontsize=14)
1145    colors=['b','g','r','c','m','k']
1146    Ymax = 1.0
1147    lenX = 0
1148    X,Y = XY[:2]
1149    Ymax = max(Ymax,max(Y))
1150    Plot.plot(X,Y,'k',picker=False)
1151    if not newPlot:
1152        Page.toolbar.push_current()
1153        Plot.set_xlim(xylim[0])
1154        Plot.set_ylim(xylim[1])
1155        xylim = []
1156        Page.toolbar.push_current()
1157        Page.toolbar.draw()
1158    else:
1159        Page.canvas.draw()
1160
1161################################################################################
1162##### PlotPowderLines
1163################################################################################
1164           
1165def PlotPowderLines(G2frame):
1166    ''' plotting of powder lines (i.e. no powder pattern) as sticks
1167    '''
1168    global HKL
1169
1170    def OnMotion(event):
1171        xpos = event.xdata
1172        if xpos:                                        #avoid out of frame mouse position
1173            Page.canvas.SetCursor(wx.CROSS_CURSOR)
1174            G2frame.G2plotNB.status.SetFields(['','2-theta =%9.3f '%(xpos,)])
1175            if G2frame.PickId and G2frame.PatternTree.GetItemText(G2frame.PickId) in ['Index Peak List','Unit Cells List']:
1176                found = []
1177                if len(HKL):
1178                    view = Page.toolbar._views.forward()[0][:2]
1179                    wid = view[1]-view[0]
1180                    found = HKL[np.where(np.fabs(HKL.T[5]-xpos) < 0.002*wid)]
1181                if len(found):
1182                    h,k,l = found[0][:3] 
1183                    Page.canvas.SetToolTipString('%d,%d,%d'%(int(h),int(k),int(l)))
1184                else:
1185                    Page.canvas.SetToolTipString('')
1186
1187    try:
1188        plotNum = G2frame.G2plotNB.plotList.index('Powder Lines')
1189        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
1190        Page.figure.clf()
1191        Plot = Page.figure.gca()
1192    except ValueError:
1193        Plot = G2frame.G2plotNB.addMpl('Powder Lines').gca()
1194        plotNum = G2frame.G2plotNB.plotList.index('Powder Lines')
1195        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
1196        Page.canvas.mpl_connect('motion_notify_event', OnMotion)
1197       
1198    Page.Choice = None
1199    Page.SetFocus()
1200    Plot.set_title('Powder Pattern Lines')
1201    Plot.set_xlabel(r'$\mathsf{2\theta}$',fontsize=14)
1202    PickId = G2frame.PickId
1203    PatternId = G2frame.PatternId
1204    peaks = G2frame.PatternTree.GetItemPyData(G2gd.GetPatternTreeItemId(G2frame,PatternId, 'Index Peak List'))
1205    for peak in peaks:
1206        Plot.axvline(peak[0],color='b')
1207    HKL = np.array(G2frame.HKL)
1208    for hkl in G2frame.HKL:
1209        Plot.axvline(hkl[5],color='r',dashes=(5,5))
1210    xmin = peaks[0][0]
1211    xmax = peaks[-1][0]
1212    delt = xmax-xmin
1213    xlim = [max(0,xmin-delt/20.),min(180.,xmax+delt/20.)]
1214    Plot.set_xlim(xlim)
1215    Page.canvas.draw()
1216    Page.toolbar.push_current()
1217
1218################################################################################
1219##### PlotPeakWidths
1220################################################################################
1221           
1222def PlotPeakWidths(G2frame):
1223    ''' Plotting of instrument broadening terms as function of 2-theta (future TOF)
1224    Seen when "Instrument Parameters" chosen from powder pattern data tree
1225    '''
1226    PatternId = G2frame.PatternId
1227    limitID = G2gd.GetPatternTreeItemId(G2frame,PatternId, 'Limits')
1228    if limitID:
1229        limits = G2frame.PatternTree.GetItemPyData(limitID)
1230    else:
1231        return
1232    Parms,Parms2 = G2frame.PatternTree.GetItemPyData( \
1233        G2gd.GetPatternTreeItemId(G2frame,PatternId, 'Instrument Parameters'))
1234    if 'C' in Parms['Type'][0]:
1235        try:           
1236            lam = Parms['Lam'][1]
1237        except KeyError:
1238            lam = Parms['Lam1'][1]
1239        GU = Parms['U'][0]
1240        GV = Parms['V'][0]
1241        GW = Parms['W'][0]
1242        LX = Parms['X'][0]
1243        LY = Parms['Y'][0]
1244    else:
1245        difC = Parms['difC'][0]
1246        alp = Parms['alpha'][0]
1247        bet0 = Parms['beta-0'][0]
1248        bet1 = Parms['beta-1'][0]
1249        sig = Parms['var-inst'][0]
1250        LX = Parms['X'][0]
1251        LY = Parms['Y'][0]
1252    peakID = G2gd.GetPatternTreeItemId(G2frame,PatternId, 'Peak List')
1253    if peakID:
1254        peaks = G2frame.PatternTree.GetItemPyData(peakID)
1255    else:
1256        peaks = []
1257   
1258    try:
1259        plotNum = G2frame.G2plotNB.plotList.index('Peak Widths')
1260        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
1261        Page.figure.clf()
1262        Plot = Page.figure.gca()
1263    except ValueError:
1264        Plot = G2frame.G2plotNB.addMpl('Peak Widths').gca()
1265        plotNum = G2frame.G2plotNB.plotList.index('Peak Widths')
1266        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
1267    Page.Choice = None
1268    Page.SetFocus()
1269   
1270    Page.canvas.SetToolTipString('')
1271    colors=['b','g','r','c','m','k']
1272    Xmin,Xmax = limits[1]
1273    Xmin = min(0.5,max(Xmin,1))
1274    Xmin /= 2
1275    Xmax /= 2
1276    nPts = 100
1277    delt = (Xmax-Xmin)/nPts
1278    thetas = []
1279    for i in range(nPts):
1280        thetas.append(Xmin+i*delt)
1281    X = []
1282    Y = []
1283    Z = []
1284    W = []
1285    sig = lambda Th,U,V,W: 1.17741*math.sqrt(U*tand(Th)**2+V*tand(Th)+W)*math.pi/18000.
1286    gam = lambda Th,X,Y: (X/cosd(Th)+Y*tand(Th))*math.pi/18000.
1287    gamFW = lambda s,g: math.exp(math.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.)
1288#    gamFW2 = lambda s,g: math.sqrt(s**2+(0.4654996*g)**2)+.5345004*g  #Ubaldo Bafile - private communication
1289    Plot.set_title('Instrument and sample peak widths')
1290    Plot.set_ylabel(r'$\Delta q/q, \Delta d/d$',fontsize=14)
1291    Plot.set_xlabel(r'$q, \AA^{-1}$',fontsize=14)
1292    if 'C' in Parms['Type'][0]:
1293        try:
1294            for theta in thetas:
1295                X.append(4.0*math.pi*sind(theta)/lam)              #q
1296                s = sig(theta,GU,GV,GW)
1297                g = gam(theta,LX,LY)
1298                G = gamFW(g,s)
1299                Y.append(s/tand(theta))
1300                Z.append(g/tand(theta))
1301                W.append(G/tand(theta))
1302            Plot.plot(X,Y,color='r',label='Gaussian')
1303            Plot.plot(X,Z,color='g',label='Lorentzian')
1304            Plot.plot(X,W,color='b',label='G+L')
1305            X = []
1306            Y = []
1307            Z = []
1308            W = []
1309            V = []
1310            for peak in peaks:
1311                X.append(4.0*math.pi*sind(peak[0]/2.0)/lam)
1312                try:
1313                    s = 1.17741*math.sqrt(peak[4])*math.pi/18000.
1314                except ValueError:
1315                    s = 0.01
1316                g = peak[6]*math.pi/18000.
1317                G = gamFW(g,s)
1318                Y.append(s/tand(peak[0]/2.))
1319                Z.append(g/tand(peak[0]/2.))
1320                W.append(G/tand(peak[0]/2.))
1321            Plot.plot(X,Y,'+',color='r',label='G peak')
1322            Plot.plot(X,Z,'+',color='g',label='L peak')
1323            Plot.plot(X,W,'+',color='b',label='G+L peak')
1324            Plot.legend(loc='best')
1325            Page.canvas.draw()
1326        except ValueError:
1327            print '**** ERROR - default U,V,W profile coefficients yield sqrt of negative value at 2theta =', \
1328                '%.3f'%(2*theta)
1329            G2frame.G2plotNB.Delete('Peak Widths')
1330    else:
1331        pass    #for TOF peak parms
1332
1333   
1334################################################################################
1335##### PlotSizeStrainPO
1336################################################################################
1337           
1338def PlotSizeStrainPO(G2frame,data,Start=False):
1339    '''Plot 3D mustrain/size/preferred orientation figure. In this instance data is for a phase
1340    '''
1341   
1342    PatternId = G2frame.PatternId
1343    generalData = data['General']
1344    SGData = generalData['SGData']
1345    SGLaue = SGData['SGLaue']
1346    if Start:                   #initialize the spherical harmonics qlmn arrays
1347        ptx.pyqlmninit()
1348        Start = False
1349#    MuStrKeys = G2spc.MustrainNames(SGData)
1350    cell = generalData['Cell'][1:]
1351    A,B = G2lat.cell2AB(cell[:6])
1352    Vol = cell[6]
1353    useList = data['Histograms']
1354    phase = generalData['Name']
1355    plotType = generalData['Data plot type']
1356    plotDict = {'Mustrain':'Mustrain','Size':'Size','Preferred orientation':'Pref.Ori.'}
1357    for ptype in plotDict:
1358        G2frame.G2plotNB.Delete(ptype)
1359    if plotType in ['None']:
1360        return       
1361
1362    for item in useList:
1363        if useList[item]['Show']:
1364            break
1365    else:
1366        return            #nothing to show!!
1367   
1368    numPlots = len(useList)
1369
1370    if plotType in ['Mustrain','Size']:
1371        Plot = mp3d.Axes3D(G2frame.G2plotNB.add3D(plotType))
1372    else:
1373        Plot = G2frame.G2plotNB.addMpl(plotType).gca()       
1374    plotNum = G2frame.G2plotNB.plotList.index(plotType)
1375    Page = G2frame.G2plotNB.nb.GetPage(plotNum)
1376    Page.Choice = None
1377    Page.SetFocus()
1378    G2frame.G2plotNB.status.SetStatusText('',1)
1379    if not Page.IsShown():
1380        Page.Show()
1381   
1382    for item in useList:
1383        if useList[item]['Show']:
1384            PHI = np.linspace(0.,360.,30,True)
1385            PSI = np.linspace(0.,180.,30,True)
1386            X = np.outer(npsind(PHI),npsind(PSI))
1387            Y = np.outer(npcosd(PHI),npsind(PSI))
1388            Z = np.outer(np.ones(np.size(PHI)),npcosd(PSI))
1389            try:        #temp patch instead of 'mustrain' for old files with 'microstrain'
1390                coeff = useList[item][plotDict[plotType]]
1391            except KeyError:
1392                break
1393            if plotType in ['Mustrain','Size']:
1394                if coeff[0] == 'isotropic':
1395                    X *= coeff[1][0]
1396                    Y *= coeff[1][0]
1397                    Z *= coeff[1][0]                               
1398                elif coeff[0] == 'uniaxial':
1399                   
1400                    def uniaxCalc(xyz,iso,aniso,axes):
1401                        Z = np.array(axes)
1402                        cp = abs(np.dot(xyz,Z))
1403                        sp = np.sqrt(1.-cp**2)
1404                        R = iso*aniso/np.sqrt((iso*cp)**2+(aniso*sp)**2)
1405                        return R*xyz
1406                       
1407                    iso,aniso = coeff[1][:2]
1408                    axes = np.inner(A,np.array(coeff[3]))
1409                    axes /= nl.norm(axes)
1410                    Shkl = np.array(coeff[1])
1411                    XYZ = np.dstack((X,Y,Z))
1412                    XYZ = np.nan_to_num(np.apply_along_axis(uniaxCalc,2,XYZ,iso,aniso,axes))
1413                    X,Y,Z = np.dsplit(XYZ,3)
1414                    X = X[:,:,0]
1415                    Y = Y[:,:,0]
1416                    Z = Z[:,:,0]
1417               
1418                elif coeff[0] == 'ellipsoidal':
1419                   
1420                    def ellipseCalc(xyz,E,R):
1421                        XYZ = xyz*E.T
1422                        return np.inner(XYZ.T,R)
1423                       
1424                    S6 = coeff[4]
1425                    Sij = G2lat.U6toUij(S6)
1426                    E,R = nl.eigh(Sij)
1427                    XYZ = np.dstack((X,Y,Z))
1428                    XYZ = np.nan_to_num(np.apply_along_axis(ellipseCalc,2,XYZ,E,R))
1429                    X,Y,Z = np.dsplit(XYZ,3)
1430                    X = X[:,:,0]
1431                    Y = Y[:,:,0]
1432                    Z = Z[:,:,0]
1433                   
1434                elif coeff[0] == 'generalized':
1435                   
1436                    def genMustrain(xyz,SGData,A,Shkl):
1437                        uvw = np.inner(A.T,xyz)
1438                        Strm = np.array(G2spc.MustrainCoeff(uvw,SGData))
1439                        sum = np.sum(np.multiply(Shkl,Strm))
1440                        sum = np.where(sum > 0.01,sum,0.01)
1441                        sum = np.sqrt(sum)*math.pi/0.018      #centidegrees to radians!
1442                        return sum*xyz
1443                       
1444                    Shkl = np.array(coeff[4])
1445                    if np.any(Shkl):
1446                        XYZ = np.dstack((X,Y,Z))
1447                        XYZ = np.nan_to_num(np.apply_along_axis(genMustrain,2,XYZ,SGData,A,Shkl))
1448                        X,Y,Z = np.dsplit(XYZ,3)
1449                        X = X[:,:,0]
1450                        Y = Y[:,:,0]
1451                        Z = Z[:,:,0]
1452                           
1453                if np.any(X) and np.any(Y) and np.any(Z):
1454                    errFlags = np.seterr(all='ignore')
1455                    Plot.plot_surface(X,Y,Z,rstride=1,cstride=1,color='g',linewidth=1)
1456                    np.seterr(all='ignore')
1457                    xyzlim = np.array([Plot.get_xlim3d(),Plot.get_ylim3d(),Plot.get_zlim3d()]).T
1458                    XYZlim = [min(xyzlim[0]),max(xyzlim[1])]
1459                    Plot.set_xlim3d(XYZlim)
1460                    Plot.set_ylim3d(XYZlim)
1461                    Plot.set_zlim3d(XYZlim)
1462                    Plot.set_aspect('equal')
1463                if plotType == 'Size':
1464                    Plot.set_title('Crystallite size for '+phase+'\n'+coeff[0]+' model')
1465                    Plot.set_xlabel(r'X, $\mu$m')
1466                    Plot.set_ylabel(r'Y, $\mu$m')
1467                    Plot.set_zlabel(r'Z, $\mu$m')
1468                else:   
1469                    Plot.set_title(r'$\mu$strain for '+phase+'\n'+coeff[0]+' model')
1470                    Plot.set_xlabel(r'X, $\mu$strain')
1471                    Plot.set_ylabel(r'Y, $\mu$strain')
1472                    Plot.set_zlabel(r'Z, $\mu$strain')
1473            else:
1474                h,k,l = generalData['POhkl']
1475                if coeff[0] == 'MD':
1476                    print 'March-Dollase preferred orientation plot'
1477               
1478                else:
1479                    PH = np.array(generalData['POhkl'])
1480                    phi,beta = G2lat.CrsAng(PH,cell[:6],SGData)
1481                    SHCoef = {}
1482                    for item in coeff[5]:
1483                        L,N = eval(item.strip('C'))
1484                        SHCoef['C%d,0,%d'%(L,N)] = coeff[5][item]                       
1485                    ODFln = G2lat.Flnh(Start,SHCoef,phi,beta,SGData)
1486                    X = np.linspace(0,90.0,26)
1487                    Y = G2lat.polfcal(ODFln,'0',X,0.0)
1488                    Plot.plot(X,Y,color='k',label=str(PH))
1489                    Plot.legend(loc='best')
1490                    Plot.set_title('Axial distribution for HKL='+str(PH)+' in '+phase)
1491                    Plot.set_xlabel(r'$\psi$',fontsize=16)
1492                    Plot.set_ylabel('MRD',fontsize=14)
1493    Page.canvas.draw()
1494   
1495################################################################################
1496##### PlotTexture
1497################################################################################
1498           
1499def PlotTexture(G2frame,data,Start=False):
1500    '''Pole figure, inverse pole figure, 3D pole distribution and 3D inverse pole distribution
1501    plotting.
1502    dict generalData contains all phase info needed which is in data
1503    '''
1504
1505    shModels = ['cylindrical','none','shear - 2/m','rolling - mmm']
1506    SamSym = dict(zip(shModels,['0','-1','2/m','mmm']))
1507    PatternId = G2frame.PatternId
1508    generalData = data['General']
1509    SGData = generalData['SGData']
1510    textureData = generalData['SH Texture']
1511    G2frame.G2plotNB.Delete('Texture')
1512    if not textureData['Order']:
1513        return                  #no plot!!
1514    SHData = generalData['SH Texture']
1515    SHCoef = SHData['SH Coeff'][1]
1516    cell = generalData['Cell'][1:7]
1517    Amat,Bmat = G2lat.cell2AB(cell)
1518    sq2 = 1.0/math.sqrt(2.0)
1519   
1520    def rp2xyz(r,p):
1521        z = npcosd(r)
1522        xy = np.sqrt(1.-z**2)
1523        return xy*npsind(p),xy*npcosd(p),z
1524           
1525    def OnMotion(event):
1526        SHData = data['General']['SH Texture']
1527        if event.xdata and event.ydata:                 #avoid out of frame errors
1528            xpos = event.xdata
1529            ypos = event.ydata
1530            if 'Inverse' in SHData['PlotType']:
1531                r = xpos**2+ypos**2
1532                if r <= 1.0:
1533                    if 'equal' in G2frame.Projection: 
1534                        r,p = 2.*npasind(np.sqrt(r)*sq2),npatan2d(ypos,xpos)
1535                    else:
1536                        r,p = 2.*npatand(np.sqrt(r)),npatan2d(ypos,xpos)
1537                    ipf = G2lat.invpolfcal(IODFln,SGData,np.array([r,]),np.array([p,]))
1538                    xyz = np.inner(Amat.T,np.array([rp2xyz(r,p)]))
1539                    y,x,z = list(xyz/np.max(np.abs(xyz)))
1540                   
1541                    G2frame.G2plotNB.status.SetFields(['',
1542                        'psi =%9.3f, beta =%9.3f, MRD =%9.3f hkl=%5.2f,%5.2f,%5.2f'%(r,p,ipf,x,y,z)])
1543                                   
1544            elif 'Axial' in SHData['PlotType']:
1545                pass
1546               
1547            else:                       #ordinary pole figure
1548                z = xpos**2+ypos**2
1549                if z <= 1.0:
1550                    z = np.sqrt(z)
1551                    if 'equal' in G2frame.Projection: 
1552                        r,p = 2.*npasind(z*sq2),npatan2d(ypos,xpos)
1553                    else:
1554                        r,p = 2.*npatand(z),npatan2d(ypos,xpos)
1555                    pf = G2lat.polfcal(ODFln,SamSym[textureData['Model']],np.array([r,]),np.array([p,]))
1556                    G2frame.G2plotNB.status.SetFields(['','phi =%9.3f, gam =%9.3f, MRD =%9.3f'%(r,p,pf)])
1557   
1558    try:
1559        plotNum = G2frame.G2plotNB.plotList.index('Texture')
1560        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
1561        Page.figure.clf()
1562        Plot = Page.figure.gca()
1563        if not Page.IsShown():
1564            Page.Show()
1565    except ValueError:
1566        Plot = G2frame.G2plotNB.addMpl('Texture').gca()
1567        plotNum = G2frame.G2plotNB.plotList.index('Texture')
1568        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
1569        Page.canvas.mpl_connect('motion_notify_event', OnMotion)
1570
1571    Page.Choice = None
1572    Page.SetFocus()
1573    G2frame.G2plotNB.status.SetFields(['',''])   
1574    PH = np.array(SHData['PFhkl'])
1575    phi,beta = G2lat.CrsAng(PH,cell,SGData)
1576    ODFln = G2lat.Flnh(Start,SHCoef,phi,beta,SGData)
1577    PX = np.array(SHData['PFxyz'])
1578    gam = atan2d(PX[0],PX[1])
1579    xy = math.sqrt(PX[0]**2+PX[1]**2)
1580    xyz = math.sqrt(PX[0]**2+PX[1]**2+PX[2]**2)
1581    psi = asind(xy/xyz)
1582    IODFln = G2lat.Glnh(Start,SHCoef,psi,gam,SamSym[textureData['Model']])
1583    if 'Axial' in SHData['PlotType']:
1584        X = np.linspace(0,90.0,26)
1585        Y = G2lat.polfcal(ODFln,SamSym[textureData['Model']],X,0.0)
1586        Plot.plot(X,Y,color='k',label=str(SHData['PFhkl']))
1587        Plot.legend(loc='best')
1588        Plot.set_title('Axial distribution for HKL='+str(SHData['PFhkl']))
1589        Plot.set_xlabel(r'$\psi$',fontsize=16)
1590        Plot.set_ylabel('MRD',fontsize=14)
1591       
1592    else:       
1593        npts = 201
1594        if 'Inverse' in SHData['PlotType']:
1595            X,Y = np.meshgrid(np.linspace(1.,-1.,npts),np.linspace(-1.,1.,npts))
1596            R,P = np.sqrt(X**2+Y**2).flatten(),npatan2d(X,Y).flatten()
1597            if 'equal' in G2frame.Projection:
1598                R = np.where(R <= 1.,2.*npasind(R*sq2),0.0)
1599            else:
1600                R = np.where(R <= 1.,2.*npatand(R),0.0)
1601            Z = np.zeros_like(R)
1602            Z = G2lat.invpolfcal(IODFln,SGData,R,P)
1603            Z = np.reshape(Z,(npts,npts))
1604            CS = Plot.contour(Y,X,Z,aspect='equal')
1605            Plot.clabel(CS,fontsize=9,inline=1)
1606            try:
1607                Img = Plot.imshow(Z.T,aspect='equal',cmap=G2frame.ContourColor,extent=[-1,1,-1,1])
1608            except ValueError:
1609                pass
1610            Page.figure.colorbar(Img)
1611            Plot.set_title('Inverse pole figure for XYZ='+str(SHData['PFxyz']))
1612            Plot.set_xlabel(G2frame.Projection.capitalize()+' projection')
1613                       
1614        else:
1615            X,Y = np.meshgrid(np.linspace(1.,-1.,npts),np.linspace(-1.,1.,npts))
1616            R,P = np.sqrt(X**2+Y**2).flatten(),npatan2d(X,Y).flatten()
1617            if 'equal' in G2frame.Projection:
1618                R = np.where(R <= 1.,2.*npasind(R*sq2),0.0)
1619            else:
1620                R = np.where(R <= 1.,2.*npatand(R),0.0)
1621            Z = np.zeros_like(R)
1622            Z = G2lat.polfcal(ODFln,SamSym[textureData['Model']],R,P)
1623            Z = np.reshape(Z,(npts,npts))
1624            try:
1625                CS = Plot.contour(Y,X,Z,aspect='equal')
1626                Plot.clabel(CS,fontsize=9,inline=1)
1627            except ValueError:
1628                pass
1629            Img = Plot.imshow(Z.T,aspect='equal',cmap=G2frame.ContourColor,extent=[-1,1,-1,1])
1630            Page.figure.colorbar(Img)
1631            Plot.set_title('Pole figure for HKL='+str(SHData['PFhkl']))
1632            Plot.set_xlabel(G2frame.Projection.capitalize()+' projection')
1633    Page.canvas.draw()
1634
1635################################################################################
1636##### PlotCovariance
1637################################################################################
1638           
1639def PlotCovariance(G2frame,Data={}):
1640    if not Data:
1641        Data = G2frame.PatternTree.GetItemPyData(
1642            G2gd.GetPatternTreeItemId(G2frame,G2frame.root, 'Covariance'))
1643    if not Data:
1644        print 'No covariance matrix available'
1645        return
1646    varyList = Data['varyList']
1647    values = Data['variables']
1648    Xmax = len(varyList)
1649    covMatrix = Data['covMatrix']
1650    sig = np.sqrt(np.diag(covMatrix))
1651    xvar = np.outer(sig,np.ones_like(sig))
1652    covArray = np.divide(np.divide(covMatrix,xvar),xvar.T)
1653    title = ' for\n'+Data['title']
1654    newAtomDict = Data['newAtomDict']
1655
1656    def OnPlotKeyPress(event):
1657        newPlot = False
1658        if event.key == 's':
1659            choice = [m for m in mpl.cm.datad.keys() if not m.endswith("_r")]
1660            choice.sort()
1661            dlg = wx.SingleChoiceDialog(G2frame,'Select','Color scheme',choice)
1662            if dlg.ShowModal() == wx.ID_OK:
1663                sel = dlg.GetSelection()
1664                G2frame.VcovColor = choice[sel]
1665            else:
1666                G2frame.VcovColor = 'RdYlGn'
1667            dlg.Destroy()
1668        PlotCovariance(G2frame)
1669
1670    def OnMotion(event):
1671        if event.button:
1672            ytics = imgAx.get_yticks()
1673            ytics = np.where(ytics<len(varyList),ytics,-1)
1674            ylabs = [np.where(0<=i ,varyList[int(i)],' ') for i in ytics]
1675            imgAx.set_yticklabels(ylabs)
1676           
1677        if event.xdata and event.ydata:                 #avoid out of frame errors
1678            xpos = int(event.xdata+.5)
1679            ypos = int(event.ydata+.5)
1680            if -1 < xpos < len(varyList) and -1 < ypos < len(varyList):
1681                if xpos == ypos:
1682                    value = values[xpos]
1683                    name = varyList[xpos]
1684                    if varyList[xpos] in newAtomDict:
1685                        name,value = newAtomDict[name]                       
1686                    msg = '%s value = %.4g, esd = %.4g'%(name,value,sig[xpos])
1687                else:
1688                    msg = '%s - %s: %5.3f'%(varyList[xpos],varyList[ypos],covArray[xpos][ypos])
1689                Page.canvas.SetToolTipString(msg)
1690                G2frame.G2plotNB.status.SetFields(['',msg])
1691    try:
1692        plotNum = G2frame.G2plotNB.plotList.index('Covariance')
1693        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
1694        Page.figure.clf()
1695        Plot = Page.figure.gca()
1696        if not Page.IsShown():
1697            Page.Show()
1698    except ValueError:
1699        Plot = G2frame.G2plotNB.addMpl('Covariance').gca()
1700        plotNum = G2frame.G2plotNB.plotList.index('Covariance')
1701        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
1702        Page.canvas.mpl_connect('motion_notify_event', OnMotion)
1703        Page.canvas.mpl_connect('key_press_event', OnPlotKeyPress)
1704
1705    Page.Choice = ['s: to change colors']
1706    Page.keyPress = OnPlotKeyPress
1707    Page.SetFocus()
1708    G2frame.G2plotNB.status.SetFields(['',''])   
1709    acolor = mpl.cm.get_cmap(G2frame.VcovColor)
1710    Img = Plot.imshow(covArray,aspect='equal',cmap=acolor,interpolation='nearest',origin='lower')
1711    imgAx = Img.get_axes()
1712    ytics = imgAx.get_yticks()
1713    ylabs = [varyList[int(i)] for i in ytics[:-1]]
1714    imgAx.set_yticklabels(ylabs)
1715    colorBar = Page.figure.colorbar(Img)
1716    Plot.set_title('V-Cov matrix'+title)
1717    Plot.set_xlabel('Variable number')
1718    Plot.set_ylabel('Variable name')
1719    Page.canvas.draw()
1720   
1721################################################################################
1722##### PlotSeq
1723################################################################################
1724           
1725def PlotSeq(G2frame,SeqData,SeqSig,SeqNames,sampleParm):
1726   
1727    def OnKeyPress(event):
1728        if event.key == 's' and sampleParm:
1729            if G2frame.xAxis:
1730                G2frame.xAxis = False
1731            else:
1732                G2frame.xAxis = True
1733            Draw(False)
1734    try:
1735        plotNum = G2frame.G2plotNB.plotList.index('Sequential refinement')
1736        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
1737        Page.figure.clf()
1738        Plot = Page.figure.gca()
1739        if not Page.IsShown():
1740            Page.Show()
1741    except ValueError:
1742        Plot = G2frame.G2plotNB.addMpl('Sequential refinement').gca()
1743        plotNum = G2frame.G2plotNB.plotList.index('Sequential refinement')
1744        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
1745        Page.canvas.mpl_connect('key_press_event', OnKeyPress)
1746        G2frame.xAxis = False
1747    Page.Choice = ['s to toggle x-axis = sample environment parameter']
1748    Page.keyPress = OnKeyPress
1749       
1750    def Draw(newPlot):
1751        Page.SetFocus()
1752        G2frame.G2plotNB.status.SetFields(['','press '])
1753        if len(SeqData):
1754            Plot.clear()
1755            if G2frame.xAxis:   
1756                xName = sampleParm.keys()[0]
1757                X = sampleParm[xName]
1758            else:
1759                X = np.arange(0,len(SeqData[0]),1)
1760                xName = 'Data sequence number'
1761            for Y,sig,name in zip(SeqData,SeqSig,SeqNames):
1762                Plot.errorbar(X,Y,yerr=sig,label=name)       
1763            Plot.legend(loc='best')
1764            Plot.set_ylabel('Parameter values')
1765            Plot.set_xlabel(xName)
1766            Page.canvas.draw()           
1767    Draw(True)
1768           
1769################################################################################
1770##### PlotExposedImage & PlotImage
1771################################################################################
1772           
1773def PlotExposedImage(G2frame,newPlot=False,event=None):
1774    '''General access module for 2D image plotting
1775    '''
1776    plotNo = G2frame.G2plotNB.nb.GetSelection()
1777    if G2frame.G2plotNB.nb.GetPageText(plotNo) == '2D Powder Image':
1778        PlotImage(G2frame,newPlot,event,newImage=True)
1779    elif G2frame.G2plotNB.nb.GetPageText(plotNo) == '2D Integration':
1780        PlotIntegration(G2frame,newPlot,event)
1781
1782def PlotImage(G2frame,newPlot=False,event=None,newImage=True):
1783    '''Plot of 2D detector images as contoured plot. Also plot calibration ellipses,
1784    masks, etc.
1785    '''
1786    from matplotlib.patches import Ellipse,Arc,Circle,Polygon
1787    import numpy.ma as ma
1788    Dsp = lambda tth,wave: wave/(2.*sind(tth/2.))
1789    global Data,Masks
1790    colors=['b','g','r','c','m','k']
1791    Data = G2frame.PatternTree.GetItemPyData(
1792        G2gd.GetPatternTreeItemId(G2frame,G2frame.Image, 'Image Controls'))
1793    Masks = G2frame.PatternTree.GetItemPyData(
1794        G2gd.GetPatternTreeItemId(G2frame,G2frame.Image, 'Masks'))
1795    try:    #may be absent
1796        StrSta = G2frame.PatternTree.GetItemPyData(
1797            G2gd.GetPatternTreeItemId(G2frame,G2frame.Image, 'Stress/Strain'))
1798    except TypeError:   #is missing
1799        StrSta = {}
1800
1801    def OnImMotion(event):
1802        Page.canvas.SetToolTipString('')
1803        sizexy = Data['size']
1804        if event.xdata and event.ydata and len(G2frame.ImageZ):                 #avoid out of frame errors
1805            Page.canvas.SetCursor(wx.CROSS_CURSOR)
1806            item = G2frame.itemPicked
1807            pixelSize = Data['pixelSize']
1808            scalex = 1000./pixelSize[0]
1809            scaley = 1000./pixelSize[1]
1810            if item and G2frame.PatternTree.GetItemText(G2frame.PickId) == 'Image Controls':
1811                if 'Text' in str(item):
1812                    Page.canvas.SetToolTipString('%8.3f %8.3fmm'%(event.xdata,event.ydata))
1813                else:
1814                    xcent,ycent = Data['center']
1815                    xpos = event.xdata-xcent
1816                    ypos = event.ydata-ycent
1817                    tth,azm = G2img.GetTthAzm(event.xdata,event.ydata,Data)
1818                    if 'line3' in  str(item) or 'line4' in str(item) and not Data['fullIntegrate']:
1819                        Page.canvas.SetToolTipString('%6d deg'%(azm))
1820                    elif 'line1' in  str(item) or 'line2' in str(item):
1821                        Page.canvas.SetToolTipString('%8.3fdeg'%(tth))                           
1822            else:
1823                xpos = event.xdata
1824                ypos = event.ydata
1825                xpix = xpos*scalex
1826                ypix = ypos*scaley
1827                Int = 0
1828                if (0 <= xpix <= sizexy[0]) and (0 <= ypix <= sizexy[1]):
1829                    Int = G2frame.ImageZ[ypix][xpix]
1830                tth,azm,dsp = G2img.GetTthAzmDsp(xpos,ypos,Data)
1831                Q = 2.*math.pi/dsp
1832                if G2frame.setPoly:
1833                    G2frame.G2plotNB.status.SetFields(['','Polygon mask pick - LB next point, RB close polygon'])
1834                else:
1835                    G2frame.G2plotNB.status.SetFields(\
1836                        ['','Detector 2-th =%9.3fdeg, dsp =%9.3fA, Q = %6.5fA-1, azm = %7.2fdeg, I = %6d'%(tth,dsp,Q,azm,Int)])
1837
1838    def OnImPlotKeyPress(event):
1839        try:
1840            PickName = G2frame.PatternTree.GetItemText(G2frame.PickId)
1841        except TypeError:
1842            return
1843        if PickName == 'Masks':
1844            Xpos = event.xdata
1845            if not Xpos:            #got point out of frame
1846                return
1847            Ypos = event.ydata
1848            if event.key == 's':
1849                Masks['Points'].append([Xpos,Ypos,1.])
1850            elif event.key == 'r':
1851                tth = G2img.GetTth(Xpos,Ypos,Data)
1852                Masks['Rings'].append([tth,0.1])
1853            elif event.key == 'a':
1854                tth,azm = G2img.GetTthAzm(Xpos,Ypos,Data)
1855                azm = int(azm)               
1856                Masks['Arcs'].append([tth,[azm-5,azm+5],0.1])
1857            elif event.key == 'p':
1858                G2frame.setPoly = True
1859                Masks['Polygons'].append([])
1860                G2frame.G2plotNB.status.SetFields(['','Polygon mask active - LB pick next point, RB close polygon'])
1861            G2imG.UpdateMasks(G2frame,Masks)
1862        elif PickName == 'Image Controls':
1863            if event.key == 'c':
1864                Xpos = event.xdata
1865                if not Xpos:            #got point out of frame
1866                    return
1867                Ypos = event.ydata
1868                dlg = wx.MessageDialog(G2frame,'Are you sure you want to change the center?',
1869                    'Center change',style=wx.OK|wx.CANCEL)
1870                try:
1871                    if dlg.ShowModal() == wx.ID_OK:
1872                        print 'move center to: ',Xpos,Ypos
1873                        Data['center'] = [Xpos,Ypos]
1874                        G2imG.UpdateImageControls(G2frame,Data,Masks)
1875                finally:
1876                    dlg.Destroy()
1877            elif event.key == 'l':
1878                if G2frame.logPlot:
1879                    G2frame.logPlot = False
1880                else:
1881                    G2frame.logPlot = True
1882        PlotImage(G2frame,newImage=True)
1883           
1884    def OnKeyBox(event):
1885        if G2frame.G2plotNB.nb.GetSelection() == G2frame.G2plotNB.plotList.index('2D Powder Image'):
1886            event.key = cb.GetValue()[0]
1887            cb.SetValue(' key press')
1888            if event.key in 'l':
1889                wx.CallAfter(OnImPlotKeyPress,event)
1890                       
1891    def OnImPick(event):
1892        if G2frame.PatternTree.GetItemText(G2frame.PickId) not in ['Image Controls','Masks']:
1893            return
1894        if G2frame.setPoly:
1895            polygon = Masks['Polygons'][-1]
1896            xpos,ypos = event.mouseevent.xdata,event.mouseevent.ydata
1897            if xpos and ypos:                       #point inside image
1898                if len(polygon) > 2 and event.mouseevent.button == 3:
1899                    x0,y0 = polygon[0]
1900                    polygon.append([x0,y0])
1901                    G2frame.setPoly = False
1902                    G2frame.G2plotNB.status.SetFields(['','Polygon closed - RB drag a vertex to change shape'])
1903                else:
1904                    G2frame.G2plotNB.status.SetFields(['','New polygon point: %.1f,%.1f'%(xpos,ypos)])
1905                    polygon.append([xpos,ypos])
1906                G2imG.UpdateMasks(G2frame,Masks)
1907        else:
1908            if G2frame.itemPicked is not None: return
1909            G2frame.itemPicked = event.artist
1910            G2frame.mousePicked = event.mouseevent
1911       
1912    def OnImRelease(event):
1913        try:
1914            PickName = G2frame.PatternTree.GetItemText(G2frame.PickId)
1915        except TypeError:
1916            return
1917        if PickName not in ['Image Controls','Masks']:
1918            return
1919        pixelSize = Data['pixelSize']
1920        scalex = 1000./pixelSize[0]
1921        scaley = 1000./pixelSize[1]
1922        pixLimit = Data['pixLimit']
1923        if G2frame.itemPicked is None and PickName == 'Image Controls' and len(G2frame.ImageZ):
1924#            sizexy = Data['size']
1925            Xpos = event.xdata
1926            if not (Xpos and G2frame.ifGetRing):                   #got point out of frame
1927                return
1928            Ypos = event.ydata
1929            if Ypos and not Page.toolbar._active:         #make sure zoom/pan not selected
1930                if event.button == 1:
1931                    Xpix = Xpos*scalex
1932                    Ypix = Ypos*scaley
1933                    xpos,ypos,I,J = G2img.ImageLocalMax(G2frame.ImageZ,pixLimit,Xpix,Ypix)
1934                    if I and J:
1935                        xpos += .5                              #shift to pixel center
1936                        ypos += .5
1937                        xpos /= scalex                          #convert to mm
1938                        ypos /= scaley
1939                        Data['ring'].append([xpos,ypos])
1940                elif event.button == 3:
1941                    G2frame.dataFrame.GetStatusBar().SetStatusText('Calibrating...')
1942                    if G2img.ImageCalibrate(G2frame,Data):
1943                        G2frame.dataFrame.GetStatusBar().SetStatusText('Calibration successful - Show ring picks to check')
1944                        print 'Calibration successful'
1945                    else:
1946                        G2frame.dataFrame.GetStatusBar().SetStatusText('Calibration failed - Show ring picks to diagnose')
1947                        print 'Calibration failed'
1948                    G2frame.ifGetRing = False
1949                    G2imG.UpdateImageControls(G2frame,Data,Masks)
1950                    return
1951                PlotImage(G2frame,newImage=False)
1952            return
1953        else:
1954            xpos = event.xdata
1955            if xpos:                                        #avoid out of frame mouse position
1956                ypos = event.ydata
1957                if G2frame.ifGetRing:                          #delete a calibration ring pick
1958                    xypos = [xpos,ypos]
1959                    rings = Data['ring']
1960                    for ring in rings:
1961                        if np.allclose(ring,xypos,.01,0):
1962                            rings.remove(ring)                                                                       
1963                else:
1964                    tth,azm,dsp = G2img.GetTthAzmDsp(xpos,ypos,Data)
1965                    itemPicked = str(G2frame.itemPicked)
1966                    if 'Line2D' in itemPicked and PickName == 'Image Controls':
1967                        if 'line1' in itemPicked:
1968                            Data['IOtth'][0] = max(tth,0.001)
1969                        elif 'line2' in itemPicked:
1970                            Data['IOtth'][1] = tth
1971                        elif 'line3' in itemPicked:
1972                            Data['LRazimuth'][0] = int(azm)
1973                            if Data['fullIntegrate']:
1974                                Data['LRazimuth'][1] = Data['LRazimuth'][0]+360
1975                        elif 'line4' in itemPicked and not Data['fullIntegrate']:
1976                            Data['LRazimuth'][1] = int(azm)
1977                           
1978                        if Data['LRazimuth'][0] > Data['LRazimuth'][1]:
1979                            Data['LRazimuth'][0] -= 360
1980                           
1981                        azLim = np.array(Data['LRazimuth'])   
1982                        if np.any(azLim>360):
1983                            azLim -= 360
1984                            Data['LRazimuth'] = list(azLim)
1985                           
1986                        if  Data['IOtth'][0] > Data['IOtth'][1]:
1987                            Data['IOtth'][0],Data['IOtth'][1] = Data['IOtth'][1],Data['IOtth'][0]
1988                           
1989                        G2frame.InnerTth.SetValue("%8.2f" % (Data['IOtth'][0]))
1990                        G2frame.OuterTth.SetValue("%8.2f" % (Data['IOtth'][1]))
1991                        G2frame.Lazim.SetValue("%6d" % (Data['LRazimuth'][0]))
1992                        G2frame.Razim.SetValue("%6d" % (Data['LRazimuth'][1]))
1993                    elif 'Circle' in itemPicked and PickName == 'Masks':
1994                        spots = Masks['Points']
1995                        newPos = itemPicked.split(')')[0].split('(')[2].split(',')
1996                        newPos = np.array([float(newPos[0]),float(newPos[1])])
1997                        for spot in spots:
1998                            if np.allclose(np.array([spot[:2]]),newPos):
1999                                spot[:2] = xpos,ypos
2000                        G2imG.UpdateMasks(G2frame,Masks)
2001                    elif 'Line2D' in itemPicked and PickName == 'Masks':
2002                        Obj = G2frame.itemPicked.findobj()
2003                        rings = Masks['Rings']
2004                        arcs = Masks['Arcs']
2005                        polygons = Masks['Polygons']
2006                        for ring in G2frame.ringList:
2007                            if Obj == ring[0]:
2008                                rN = ring[1]
2009                                if ring[2] == 'o':
2010                                    rings[rN][0] = G2img.GetTth(xpos,ypos,Data)-rings[rN][1]/2.
2011                                else:
2012                                    rings[rN][0] = G2img.GetTth(xpos,ypos,Data)+rings[rN][1]/2.
2013                        for arc in G2frame.arcList:
2014                            if Obj == arc[0]:
2015                                aN = arc[1]
2016                                if arc[2] == 'o':
2017                                    arcs[aN][0] = G2img.GetTth(xpos,ypos,Data)-arcs[aN][2]/2
2018                                elif arc[2] == 'i':
2019                                    arcs[aN][0] = G2img.GetTth(xpos,ypos,Data)+arcs[aN][2]/2
2020                                elif arc[2] == 'l':
2021                                    arcs[aN][1][0] = int(G2img.GetAzm(xpos,ypos,Data))
2022                                else:
2023                                    arcs[aN][1][1] = int(G2img.GetAzm(xpos,ypos,Data))
2024                        for poly in G2frame.polyList:
2025                            if Obj == poly[0]:
2026                                ind = G2frame.itemPicked.contains(G2frame.mousePicked)[1]['ind'][0]
2027                                oldPos = np.array([G2frame.mousePicked.xdata,G2frame.mousePicked.ydata])
2028                                pN = poly[1]
2029                                for i,xy in enumerate(polygons[pN]):
2030                                    if np.allclose(np.array([xy]),oldPos,atol=1.0):
2031                                        polygons[pN][i] = xpos,ypos
2032                        G2imG.UpdateMasks(G2frame,Masks)
2033#                    else:                  #keep for future debugging
2034#                        print str(G2frame.itemPicked),event.xdata,event.ydata,event.button
2035                PlotImage(G2frame,newImage=True)
2036            G2frame.itemPicked = None
2037           
2038    try:
2039        plotNum = G2frame.G2plotNB.plotList.index('2D Powder Image')
2040        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
2041        if not newPlot:
2042            Plot = Page.figure.gca()          #get previous powder plot & get limits
2043            xylim = Plot.get_xlim(),Plot.get_ylim()
2044        if newImage:
2045            Page.figure.clf()
2046            Plot = Page.figure.gca()          #get a fresh plot after clf()
2047       
2048    except ValueError:
2049        Plot = G2frame.G2plotNB.addMpl('2D Powder Image').gca()
2050        plotNum = G2frame.G2plotNB.plotList.index('2D Powder Image')
2051        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
2052        Page.canvas.mpl_connect('key_press_event', OnImPlotKeyPress)
2053        Page.canvas.mpl_connect('motion_notify_event', OnImMotion)
2054        Page.canvas.mpl_connect('pick_event', OnImPick)
2055        Page.canvas.mpl_connect('button_release_event', OnImRelease)
2056        xylim = []
2057    Page.Choice = None
2058    if not event:                       #event from GUI TextCtrl - don't want focus to change to plot!!!
2059        Page.SetFocus()
2060    Title = G2frame.PatternTree.GetItemText(G2frame.Image)[4:]
2061    G2frame.G2plotNB.status.DestroyChildren()
2062    if G2frame.logPlot:
2063        Title = 'log('+Title+')'
2064    Plot.set_title(Title)
2065    try:
2066        if G2frame.PatternTree.GetItemText(G2frame.PickId) in ['Image Controls',]:
2067            Page.Choice = (' key press','l: log(I) on',)
2068            if G2frame.logPlot:
2069                Page.Choice[1] = 'l: log(I) off'
2070            Page.keyPress = OnImPlotKeyPress
2071    except TypeError:
2072        pass
2073    size,imagefile = G2frame.PatternTree.GetItemPyData(G2frame.Image)
2074    if imagefile != G2frame.oldImagefile:
2075        imagefile = G2IO.CheckImageFile(G2frame,imagefile)
2076        if not imagefile:
2077            G2frame.G2plotNB.Delete('2D Powder Image')
2078            return
2079        G2frame.PatternTree.SetItemPyData(G2frame.Image,[size,imagefile])
2080        G2frame.ImageZ = G2IO.GetImageData(G2frame,imagefile,imageOnly=True)
2081        G2frame.oldImagefile = imagefile
2082
2083    imScale = 1
2084    if len(G2frame.ImageZ) > 1024:
2085        imScale = len(G2frame.ImageZ)/1024
2086    sizexy = Data['size']
2087    pixelSize = Data['pixelSize']
2088    scalex = 1000./pixelSize[0]
2089    scaley = 1000./pixelSize[1]
2090    Xmax = sizexy[0]*pixelSize[0]/1000.
2091    Ymax = sizexy[1]*pixelSize[1]/1000.
2092    xlim = (0,Xmax)
2093    ylim = (Ymax,0)
2094    Imin,Imax = Data['range'][1]
2095    acolor = mpl.cm.get_cmap(Data['color'])
2096    xcent,ycent = Data['center']
2097    Plot.set_xlabel('Image x-axis, mm',fontsize=12)
2098    Plot.set_ylabel('Image y-axis, mm',fontsize=12)
2099    #do threshold mask - "real" mask - others are just bondaries
2100    Zlim = Masks['Thresholds'][1]
2101    wx.BeginBusyCursor()
2102    try:
2103           
2104        if newImage:                   
2105            MA = ma.masked_greater(ma.masked_less(G2frame.ImageZ,Zlim[0]),Zlim[1])
2106            MaskA = ma.getmaskarray(MA)
2107            A = G2img.ImageCompress(MA,imScale)
2108            AM = G2img.ImageCompress(MaskA,imScale)
2109            if G2frame.logPlot:
2110                A = np.where(A>0,np.log(A),0)
2111                AM = np.where(AM>0,np.log(AM),0)
2112                Imin,Imax = [np.amin(A),np.amax(A)]
2113            ImgM = Plot.imshow(AM,aspect='equal',cmap='Reds',
2114                interpolation='nearest',vmin=0,vmax=2,extent=[0,Xmax,Ymax,0])
2115            Img = Plot.imshow(A,aspect='equal',cmap=acolor,
2116                interpolation='nearest',vmin=Imin,vmax=Imax,extent=[0,Xmax,Ymax,0])
2117            if G2frame.setPoly:
2118                Img.set_picker(True)
2119   
2120        Plot.plot(xcent,ycent,'x')
2121        #G2frame.PatternTree.GetItemText(item)
2122        if Data['showLines']:
2123            LRAzim = Data['LRazimuth']                  #NB: integers
2124            Nazm = Data['outAzimuths']
2125            delAzm = float(LRAzim[1]-LRAzim[0])/Nazm
2126            AzmthOff = Data['azmthOff']
2127            IOtth = Data['IOtth']
2128            wave = Data['wavelength']
2129            dspI = wave/(2.0*sind(IOtth[0]/2.0))
2130            ellI = G2img.GetEllipse(dspI,Data)           #=False if dsp didn't yield an ellipse (ugh! a parabola or a hyperbola)
2131            dspO = wave/(2.0*sind(IOtth[1]/2.0))
2132            ellO = G2img.GetEllipse(dspO,Data)           #Ditto & more likely for outer ellipse
2133            Azm = np.array(range(LRAzim[0],LRAzim[1]+1))-AzmthOff
2134            if ellI:
2135                xyI = []
2136                for azm in Azm:
2137                    xyI.append(G2img.GetDetectorXY(dspI,azm-90.,Data))
2138                xyI = np.array(xyI)
2139                arcxI,arcyI = xyI.T
2140                Plot.plot(arcxI,arcyI,picker=3)
2141            if ellO:
2142                xyO = []
2143                for azm in Azm:
2144                    xyO.append(G2img.GetDetectorXY(dspO,azm-90.,Data))
2145                xyO = np.array(xyO)
2146                arcxO,arcyO = xyO.T
2147                Plot.plot(arcxO,arcyO,picker=3)
2148            if ellO and ellI:
2149                Plot.plot([arcxI[0],arcxO[0]],[arcyI[0],arcyO[0]],picker=3)
2150                Plot.plot([arcxI[-1],arcxO[-1]],[arcyI[-1],arcyO[-1]],picker=3)
2151            for i in range(Nazm):
2152                cake = LRAzim[0]+i*delAzm-AzmthOff
2153                if Data['centerAzm']:
2154                    cake += delAzm/2.
2155                ind = np.searchsorted(Azm,cake)
2156                Plot.plot([arcxI[ind],arcxO[ind]],[arcyI[ind],arcyO[ind]],color='k',dashes=(5,5))
2157                   
2158        if G2frame.PatternTree.GetItemText(G2frame.PickId) in 'Image Controls':
2159            for xring,yring in Data['ring']:
2160                Plot.plot(xring,yring,'r+',picker=3)
2161            if Data['setRings']:
2162    #            rings = np.concatenate((Data['rings']),axis=0)
2163                N = 0
2164                for ring in Data['rings']:
2165                    xring,yring = np.array(ring).T[:2]
2166                    Plot.plot(xring,yring,'+',color=colors[N%6])
2167                    N += 1           
2168            for ellipse in Data['ellipses']:
2169                cent,phi,[width,height],col = ellipse
2170                Plot.add_artist(Ellipse([cent[0],cent[1]],2*width,2*height,phi,ec=col,fc='none'))
2171                Plot.text(cent[0],cent[1],'+',color=col,ha='center',va='center')
2172        if G2frame.PatternTree.GetItemText(G2frame.PickId) in 'Stress/Strain':
2173            print 'plot stress/strain stuff'
2174            for ring in StrSta['d-zero']:
2175                xring,yring = ring['ImxyObs']
2176                Plot.plot(xring,yring,'r+')
2177#                for xring,yring in ring['ImxyCalc'].T:
2178#                    Plot.add_artist(Polygon(ring['ImxyCalc'].T,ec='b',fc='none'))
2179#                    Plot.plot(xring,yring)
2180        #masks - mask lines numbered after integration limit lines
2181        spots = Masks['Points']
2182        rings = Masks['Rings']
2183        arcs = Masks['Arcs']
2184        polygons = Masks['Polygons']
2185        for x,y,d in spots:
2186            Plot.add_artist(Circle((x,y),radius=d/2,fc='none',ec='r',picker=3))
2187        G2frame.ringList = []
2188        for iring,(tth,thick) in enumerate(rings):
2189            wave = Data['wavelength']
2190            x1,y1 = np.hsplit(np.array(G2img.makeIdealRing(G2img.GetEllipse(Dsp(tth+thick/2.,wave),Data))),2)
2191            x2,y2 = np.hsplit(np.array(G2img.makeIdealRing(G2img.GetEllipse(Dsp(tth-thick/2.,wave),Data))),2)
2192            G2frame.ringList.append([Plot.plot(x1,y1,'r',picker=3),iring,'o'])           
2193            G2frame.ringList.append([Plot.plot(x2,y2,'r',picker=3),iring,'i'])
2194        G2frame.arcList = []
2195        for iarc,(tth,azm,thick) in enumerate(arcs):           
2196            wave = Data['wavelength']
2197            x1,y1 = np.hsplit(np.array(G2img.makeIdealRing(G2img.GetEllipse(Dsp(tth+thick/2.,wave),Data),azm)),2)
2198            x2,y2 = np.hsplit(np.array(G2img.makeIdealRing(G2img.GetEllipse(Dsp(max(0.01,tth-thick/2.),wave),Data),azm)),2)
2199            G2frame.arcList.append([Plot.plot(x1,y1,'r',picker=3),iarc,'o'])           
2200            G2frame.arcList.append([Plot.plot(x2,y2,'r',picker=3),iarc,'i'])
2201            G2frame.arcList.append([Plot.plot([x1[0],x2[0]],[y1[0],y2[0]],'r',picker=3),iarc,'l'])
2202            G2frame.arcList.append([Plot.plot([x1[-1],x2[-1]],[y1[-1],y2[-1]],'r',picker=3),iarc,'u'])
2203        G2frame.polyList = []
2204        for ipoly,polygon in enumerate(polygons):
2205            x,y = np.hsplit(np.array(polygon),2)
2206            G2frame.polyList.append([Plot.plot(x,y,'r+',picker=10),ipoly])
2207            Plot.plot(x,y,'r')           
2208        if newImage:
2209            colorBar = Page.figure.colorbar(Img)
2210        Plot.set_xlim(xlim)
2211        Plot.set_ylim(ylim)
2212        Plot.invert_yaxis()
2213        if not newPlot and xylim:
2214            Page.toolbar.push_current()
2215            Plot.set_xlim(xylim[0])
2216            Plot.set_ylim(xylim[1])
2217            xylim = []
2218            Page.toolbar.push_current()
2219            Page.toolbar.draw()
2220        else:
2221            Page.canvas.draw()
2222    finally:
2223        wx.EndBusyCursor()
2224       
2225################################################################################
2226##### PlotIntegration
2227################################################################################
2228           
2229def PlotIntegration(G2frame,newPlot=False,event=None):
2230    '''Plot of 2D image after image integration with 2-theta and azimuth as coordinates
2231    '''
2232           
2233    def OnMotion(event):
2234        Page.canvas.SetToolTipString('')
2235        Page.canvas.SetCursor(wx.CROSS_CURSOR)
2236        azm = event.ydata
2237        tth = event.xdata
2238        if azm and tth:
2239            G2frame.G2plotNB.status.SetFields(\
2240                ['','Detector 2-th =%9.3fdeg, azm = %7.2fdeg'%(tth,azm)])
2241                               
2242    try:
2243        plotNum = G2frame.G2plotNB.plotList.index('2D Integration')
2244        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
2245        if not newPlot:
2246            Plot = Page.figure.gca()          #get previous plot & get limits
2247            xylim = Plot.get_xlim(),Plot.get_ylim()
2248        Page.figure.clf()
2249        Plot = Page.figure.gca()          #get a fresh plot after clf()
2250       
2251    except ValueError:
2252        Plot = G2frame.G2plotNB.addMpl('2D Integration').gca()
2253        plotNum = G2frame.G2plotNB.plotList.index('2D Integration')
2254        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
2255        Page.canvas.mpl_connect('motion_notify_event', OnMotion)
2256        Page.views = False
2257        view = False
2258    Page.Choice = None
2259    if not event:
2260        Page.SetFocus()
2261       
2262    Data = G2frame.PatternTree.GetItemPyData(
2263        G2gd.GetPatternTreeItemId(G2frame,G2frame.Image, 'Image Controls'))
2264    image = G2frame.Integrate[0]
2265    xsc = G2frame.Integrate[1]
2266    ysc = G2frame.Integrate[2]
2267    Imin,Imax = Data['range'][1]
2268    acolor = mpl.cm.get_cmap(Data['color'])
2269    Plot.set_title(G2frame.PatternTree.GetItemText(G2frame.Image)[4:])
2270    Plot.set_ylabel('azimuth',fontsize=12)
2271    Plot.set_xlabel('2-theta',fontsize=12)
2272    Img = Plot.imshow(image,cmap=acolor,vmin=Imin,vmax=Imax,interpolation='nearest', \
2273        extent=[ysc[0],ysc[-1],xsc[-1],xsc[0]],aspect='auto')
2274    colorBar = Page.figure.colorbar(Img)
2275    if Data['ellipses']:           
2276        for ellipse in Data['ellipses']:
2277            ring = np.array(G2img.makeIdealRing(ellipse[:3])) #skip color
2278            x,y = np.hsplit(ring,2)
2279            tth,azm = G2img.GetTthAzm(x,y,Data)
2280#            azm = np.where(azm < 0.,azm+360,azm)
2281            Plot.plot(tth,azm,'b,')
2282    if not newPlot:
2283        Page.toolbar.push_current()
2284        Plot.set_xlim(xylim[0])
2285        Plot.set_ylim(xylim[1])
2286        xylim = []
2287        Page.toolbar.push_current()
2288        Page.toolbar.draw()
2289    else:
2290        Page.canvas.draw()
2291               
2292################################################################################
2293##### PlotTRImage
2294################################################################################
2295           
2296def PlotTRImage(G2frame,tax,tay,taz,newPlot=False):
2297    '''a test plot routine - not normally used
2298    ''' 
2299           
2300    def OnMotion(event):
2301        Page.canvas.SetToolTipString('')
2302        Page.canvas.SetCursor(wx.CROSS_CURSOR)
2303        azm = event.xdata
2304        tth = event.ydata
2305        if azm and tth:
2306            G2frame.G2plotNB.status.SetFields(\
2307                ['','Detector 2-th =%9.3fdeg, azm = %7.2fdeg'%(tth,azm)])
2308                               
2309    try:
2310        plotNum = G2frame.G2plotNB.plotList.index('2D Transformed Powder Image')
2311        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
2312        if not newPlot:
2313            Plot = Page.figure.gca()          #get previous plot & get limits
2314            xylim = Plot.get_xlim(),Plot.get_ylim()
2315        Page.figure.clf()
2316        Plot = Page.figure.gca()          #get a fresh plot after clf()
2317       
2318    except ValueError:
2319        Plot = G2frame.G2plotNB.addMpl('2D Transformed Powder Image').gca()
2320        plotNum = G2frame.G2plotNB.plotList.index('2D Transformed Powder Image')
2321        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
2322        Page.canvas.mpl_connect('motion_notify_event', OnMotion)
2323        Page.views = False
2324        view = False
2325    Page.Choice = None
2326    Page.SetFocus()
2327       
2328    Data = G2frame.PatternTree.GetItemPyData(
2329        G2gd.GetPatternTreeItemId(G2frame,G2frame.Image, 'Image Controls'))
2330    Imin,Imax = Data['range'][1]
2331    step = (Imax-Imin)/5.
2332    V = np.arange(Imin,Imax,step)
2333    acolor = mpl.cm.get_cmap(Data['color'])
2334    Plot.set_title(G2frame.PatternTree.GetItemText(G2frame.Image)[4:])
2335    Plot.set_xlabel('azimuth',fontsize=12)
2336    Plot.set_ylabel('2-theta',fontsize=12)
2337    Plot.contour(tax,tay,taz,V,cmap=acolor)
2338    if Data['showLines']:
2339        IOtth = Data['IOtth']
2340        if Data['fullIntegrate']:
2341            LRAzim = [-180,180]
2342        else:
2343            LRAzim = Data['LRazimuth']                  #NB: integers
2344        Plot.plot([LRAzim[0],LRAzim[1]],[IOtth[0],IOtth[0]],picker=True)
2345        Plot.plot([LRAzim[0],LRAzim[1]],[IOtth[1],IOtth[1]],picker=True)
2346        if not Data['fullIntegrate']:
2347            Plot.plot([LRAzim[0],LRAzim[0]],[IOtth[0],IOtth[1]],picker=True)
2348            Plot.plot([LRAzim[1],LRAzim[1]],[IOtth[0],IOtth[1]],picker=True)
2349    if Data['setRings']:
2350        rings = np.concatenate((Data['rings']),axis=0)
2351        for xring,yring,dsp in rings:
2352            x,y = G2img.GetTthAzm(xring,yring,Data)
2353            Plot.plot(y,x,'r+')           
2354    if Data['ellipses']:           
2355        for ellipse in Data['ellipses']:
2356            ring = np.array(G2img.makeIdealRing(ellipse[:3])) #skip color
2357            x,y = np.hsplit(ring,2)
2358            tth,azm = G2img.GetTthAzm(x,y,Data)
2359            Plot.plot(azm,tth,'b,')
2360    if not newPlot:
2361        Page.toolbar.push_current()
2362        Plot.set_xlim(xylim[0])
2363        Plot.set_ylim(xylim[1])
2364        xylim = []
2365        Page.toolbar.push_current()
2366        Page.toolbar.draw()
2367    else:
2368        Page.canvas.draw()
2369       
2370################################################################################
2371##### PlotStructure
2372################################################################################
2373           
2374def PlotStructure(G2frame,data):
2375    '''Crystal structure plotting package. Can show structures as balls, sticks, lines,
2376    thermal motion ellipsoids and polyhedra
2377    '''
2378    ForthirdPI = 4.0*math.pi/3.0
2379    generalData = data['General']
2380    cell = generalData['Cell'][1:7]
2381    Vol = generalData['Cell'][7:8][0]
2382    Amat,Bmat = G2lat.cell2AB(cell)         #Amat - crystal to cartesian, Bmat - inverse
2383    Gmat,gmat = G2lat.cell2Gmat(cell)
2384    A4mat = np.concatenate((np.concatenate((Amat,[[0],[0],[0]]),axis=1),[[0,0,0,1],]),axis=0)
2385    B4mat = np.concatenate((np.concatenate((Bmat,[[0],[0],[0]]),axis=1),[[0,0,0,1],]),axis=0)
2386    Mydir = generalData['Mydir']
2387    atomData = data['Atoms']
2388    mapPeaks = []
2389    if 'Map Peaks' in data:
2390        mapPeaks = np.array(data['Map Peaks'])
2391        peakMax = 100.
2392        if len(mapPeaks):
2393            peakMax = np.max(mapPeaks.T[0])
2394    drawingData = data['Drawing']
2395    try:
2396        drawAtoms = drawingData['Atoms']
2397    except KeyError:
2398        drawAtoms = []
2399    mapData = {}
2400    flipData = {}
2401    rhoXYZ = []
2402    if 'Map' in generalData:
2403        mapData = generalData['Map']
2404    if 'Flip' in generalData:
2405        flipData = generalData['Flip']                       
2406        flipData['mapRoll'] = [0,0,0]
2407    cx,ct,cs,ci = drawingData['atomPtrs']
2408    Wt = np.array([255,255,255])
2409    Rd = np.array([255,0,0])
2410    Gr = np.array([0,255,0])
2411    Bl = np.array([0,0,255])
2412    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]])
2413    uEdges = np.array([
2414        [uBox[0],uBox[1]],[uBox[0],uBox[3]],[uBox[0],uBox[4]],[uBox[1],uBox[2]], 
2415        [uBox[2],uBox[3]],[uBox[1],uBox[5]],[uBox[2],uBox[6]],[uBox[3],uBox[7]], 
2416        [uBox[4],uBox[5]],[uBox[5],uBox[6]],[uBox[6],uBox[7]],[uBox[7],uBox[4]]])
2417    mD = 0.1
2418    mV = np.array([[[-mD,0,0],[mD,0,0]],[[0,-mD,0],[0,mD,0]],[[0,0,-mD],[0,0,mD]]])
2419    mapPeakVecs = np.inner(mV,Bmat)
2420
2421    uColors = [Rd,Gr,Bl,Wt, Wt,Wt,Wt,Wt, Wt,Wt,Wt,Wt]
2422    altDown = False
2423    shiftDown = False
2424    ctrlDown = False
2425   
2426    def OnKeyBox(event):
2427        import Image
2428#        Draw()                          #make sure plot is fresh!!
2429        mode = cb.GetValue()
2430        if mode in ['jpeg','bmp','tiff',]:
2431            Fname = os.path.joint(Mydir,generalData['Name']+'.'+mode)
2432            size = Page.canvas.GetSize()
2433            glPixelStorei(GL_UNPACK_ALIGNMENT, 1)
2434            if mode in ['jpeg',]:
2435                Pix = glReadPixels(0,0,size[0],size[1],GL_RGBA, GL_UNSIGNED_BYTE)
2436                im = Image.new("RGBA", (size[0],size[1]))
2437            else:
2438                Pix = glReadPixels(0,0,size[0],size[1],GL_RGB, GL_UNSIGNED_BYTE)
2439                im = Image.new("RGB", (size[0],size[1]))
2440            im.fromstring(Pix)
2441            im.save(Fname,mode)
2442            cb.SetValue(' save as/key:')
2443            G2frame.G2plotNB.status.SetStatusText('Drawing saved to: '+Fname,1)
2444        else:
2445            event.key = cb.GetValue()[0]
2446            cb.SetValue(' save as/key:')
2447            wx.CallAfter(OnKey,event)
2448
2449    def OnKey(event):           #on key UP!!
2450#        Draw()                          #make sure plot is fresh!!
2451        try:
2452            keyCode = event.GetKeyCode()
2453            if keyCode > 255:
2454                keyCode = 0
2455            key = chr(keyCode)
2456        except AttributeError:       #if from OnKeyBox above
2457            key = str(event.key).upper()
2458        indx = drawingData['selectedAtoms']
2459        if key in ['C']:
2460            drawingData['viewPoint'] = [[.5,.5,.5],[0,0]]
2461            drawingData['viewDir'] = [0,0,1]
2462            drawingData['oldxy'] = []
2463            V0 = np.array([0,0,1])
2464            V = np.inner(Amat,V0)
2465            V /= np.sqrt(np.sum(V**2))
2466            A = np.arccos(np.sum(V*V0))
2467            Q = G2mth.AV2Q(A,[0,1,0])
2468            drawingData['Quaternion'] = Q
2469            SetViewPointText(drawingData['viewPoint'][0])
2470            SetViewDirText(drawingData['viewDir'])
2471            Q = drawingData['Quaternion']
2472            G2frame.G2plotNB.status.SetStatusText('New quaternion: %.2f+, %.2fi+ ,%.2fj+, %.2fk'%(Q[0],Q[1],Q[2],Q[3]),1)
2473        elif key in ['N']:
2474            drawAtoms = drawingData['Atoms']
2475            pI = drawingData['viewPoint'][1]
2476            if indx:
2477                pI[0] = indx[pI[1]]
2478                Tx,Ty,Tz = drawAtoms[pI[0]][cx:cx+3]
2479                pI[1] += 1
2480                if pI[1] >= len(indx):
2481                    pI[1] = 0
2482            else:
2483                Tx,Ty,Tz = drawAtoms[pI[0]][cx:cx+3]               
2484                pI[0] += 1
2485                if pI[0] >= len(drawAtoms):
2486                    pI[0] = 0
2487            drawingData['viewPoint'] = [[Tx,Ty,Tz],pI]
2488            SetViewPointText(drawingData['viewPoint'][0])
2489            G2frame.G2plotNB.status.SetStatusText('View point at atom '+drawAtoms[pI[0]][ct-1]+str(pI),1)
2490               
2491        elif key in ['P']:
2492            drawAtoms = drawingData['Atoms']
2493            pI = drawingData['viewPoint'][1]
2494            if indx:
2495                pI[0] = indx[pI[1]]
2496                Tx,Ty,Tz = drawAtoms[pI[0]][cx:cx+3]
2497                pI[1] -= 1
2498                if pI[1] < 0:
2499                    pI[1] = len(indx)-1
2500            else:
2501                Tx,Ty,Tz = drawAtoms[pI[0]][cx:cx+3]               
2502                pI[0] -= 1
2503                if pI[0] < 0:
2504                    pI[0] = len(drawAtoms)-1
2505            drawingData['viewPoint'] = [[Tx,Ty,Tz],pI]
2506            SetViewPointText(drawingData['viewPoint'][0])           
2507            G2frame.G2plotNB.status.SetStatusText('View point at atom '+drawAtoms[pI[0]][ct-1]+str(pI),1)
2508        elif key in ['U','D','L','R'] and mapData['Flip'] == True:
2509            dirDict = {'U':[0,1],'D':[0,-1],'L':[-1,0],'R':[1,0]}
2510            SetMapRoll(dirDict[key])
2511            SetPeakRoll(dirDict[key])
2512            SetMapPeaksText(mapPeaks)
2513        Draw('key')
2514           
2515    def GetTruePosition(xy,Add=False):
2516        View = glGetIntegerv(GL_VIEWPORT)
2517        Proj = glGetDoublev(GL_PROJECTION_MATRIX)
2518        Model = glGetDoublev(GL_MODELVIEW_MATRIX)
2519        Zmax = 1.
2520        if Add:
2521            Indx = GetSelectedAtoms()
2522        if G2frame.dataDisplay.GetPageText(getSelection()) == 'Map peaks':
2523            for i,peak in enumerate(mapPeaks):
2524                x,y,z = peak[1:4]
2525                X,Y,Z = gluProject(x,y,z,Model,Proj,View)
2526                XY = [int(X),int(View[3]-Y)]
2527                if np.allclose(xy,XY,atol=10) and Z < Zmax:
2528                    Zmax = Z
2529                    try:
2530                        Indx.remove(i)
2531                        ClearSelectedAtoms()
2532                        for id in Indx:
2533                            SetSelectedAtoms(id,Add)
2534                    except:
2535                        SetSelectedAtoms(i,Add)
2536        else:
2537            for i,atom in enumerate(drawAtoms):
2538                x,y,z = atom[cx:cx+3]
2539                X,Y,Z = gluProject(x,y,z,Model,Proj,View)
2540                XY = [int(X),int(View[3]-Y)]
2541                if np.allclose(xy,XY,atol=10) and Z < Zmax:
2542                    Zmax = Z
2543                    try:
2544                        Indx.remove(i)
2545                        ClearSelectedAtoms()
2546                        for id in Indx:
2547                            SetSelectedAtoms(id,Add)
2548                    except:
2549                        SetSelectedAtoms(i,Add)
2550                                       
2551    def OnMouseDown(event):
2552        xy = event.GetPosition()
2553        if event.ShiftDown():
2554            if event.LeftIsDown():
2555                GetTruePosition(xy)
2556            elif event.RightIsDown():
2557                GetTruePosition(xy,True)
2558        else:
2559            drawingData['oldxy'] = list(xy)
2560#        Draw()
2561       
2562    def OnMouseMove(event):
2563        if event.ShiftDown():           #don't want any inadvertant moves when picking
2564            return
2565        newxy = event.GetPosition()
2566        page = getSelection()
2567                               
2568        if event.Dragging() and not event.ControlDown():
2569            if event.LeftIsDown():
2570                SetRotation(newxy)
2571                Q = drawingData['Quaternion']
2572                G2frame.G2plotNB.status.SetStatusText('New quaternion: %.2f+, %.2fi+ ,%.2fj+, %.2fk'%(Q[0],Q[1],Q[2],Q[3]),1)
2573            elif event.RightIsDown():
2574                SetTranslation(newxy)
2575                Tx,Ty,Tz = drawingData['viewPoint'][0]
2576                G2frame.G2plotNB.status.SetStatusText('New view point: %.4f, %.4f, %.4f'%(Tx,Ty,Tz),1)
2577            elif event.MiddleIsDown():
2578                SetRotationZ(newxy)
2579                Q = drawingData['Quaternion']
2580                G2frame.G2plotNB.status.SetStatusText('New quaternion: %.2f+, %.2fi+ ,%.2fj+, %.2fk'%(Q[0],Q[1],Q[2],Q[3]),1)
2581            Draw('move')
2582       
2583    def OnMouseWheel(event):
2584        if event.ShiftDown():
2585            return
2586        drawingData['cameraPos'] += event.GetWheelRotation()/24
2587        drawingData['cameraPos'] = max(10,min(500,drawingData['cameraPos']))
2588        G2frame.G2plotNB.status.SetStatusText('New camera distance: %.2f'%(drawingData['cameraPos']),1)
2589        page = getSelection()
2590        if page:
2591            if G2frame.dataDisplay.GetPageText(page) == 'Draw Options':
2592                panel = G2frame.dataDisplay.GetPage(page).GetChildren()[0].GetChildren()
2593                names = [child.GetName() for child in panel]
2594                panel[names.index('cameraPos')].SetLabel('Camera Position: '+'%.2f'%(drawingData['cameraPos']))
2595                panel[names.index('cameraSlider')].SetValue(drawingData['cameraPos'])
2596            Draw('wheel')
2597       
2598    def getSelection():
2599        try:
2600            return G2frame.dataDisplay.GetSelection()
2601        except AttributeError:
2602            G2frame.G2plotNB.status.SetStatusText('Select this from Phase data window!')
2603            return 0
2604           
2605    def SetViewPointText(VP):
2606        page = getSelection()
2607        if page:
2608            if G2frame.dataDisplay.GetPageText(page) == 'Draw Options':
2609                panel = G2frame.dataDisplay.GetPage(page).GetChildren()[0].GetChildren()
2610                names = [child.GetName() for child in panel]
2611                panel[names.index('viewPoint')].SetValue('%.3f %.3f %.3f'%(VP[0],VP[1],VP[2]))
2612               
2613    def SetViewDirText(VD):
2614        page = getSelection()
2615        if page:
2616            if G2frame.dataDisplay.GetPageText(page) == 'Draw Options':
2617                panel = G2frame.dataDisplay.GetPage(page).GetChildren()[0].GetChildren()
2618                names = [child.GetName() for child in panel]
2619                panel[names.index('viewDir')].SetValue('%.3f %.3f %.3f'%(VD[0],VD[1],VD[2]))
2620               
2621    def SetMapPeaksText(mapPeaks):
2622        page = getSelection()
2623        if page:
2624            if G2frame.dataDisplay.GetPageText(page) == 'Map peaks':
2625                G2frame.MapPeaksTable.SetData(mapPeaks)
2626                panel = G2frame.dataDisplay.GetPage(page).GetChildren()
2627                names = [child.GetName() for child in panel]
2628                panel[names.index('grid window')].Refresh()
2629           
2630    def ClearSelectedAtoms():
2631        page = getSelection()
2632        if page:
2633            if G2frame.dataDisplay.GetPageText(page) == 'Draw Atoms':
2634                G2frame.dataDisplay.GetPage(page).ClearSelection()      #this is the Atoms grid in Draw Atoms
2635            elif G2frame.dataDisplay.GetPageText(page) == 'Map peaks':
2636                G2frame.dataDisplay.GetPage(page).ClearSelection()      #this is the Atoms grid in Atoms
2637            elif G2frame.dataDisplay.GetPageText(page) == 'Atoms':
2638                G2frame.dataDisplay.GetPage(page).ClearSelection()      #this is the Atoms grid in Atoms
2639               
2640                   
2641    def SetSelectedAtoms(ind,Add=False):
2642        page = getSelection()
2643        if page:
2644            if G2frame.dataDisplay.GetPageText(page) == 'Draw Atoms':
2645                G2frame.dataDisplay.GetPage(page).SelectRow(ind,Add)      #this is the Atoms grid in Draw Atoms
2646            elif G2frame.dataDisplay.GetPageText(page) == 'Map peaks':
2647                G2frame.dataDisplay.GetPage(page).SelectRow(ind,Add)                 
2648            elif G2frame.dataDisplay.GetPageText(page) == 'Atoms':
2649                Id = drawAtoms[ind][-2]
2650                for i,atom in enumerate(atomData):
2651                    if atom[-1] == Id:
2652                        G2frame.dataDisplay.GetPage(page).SelectRow(i)      #this is the Atoms grid in Atoms
2653                 
2654    def GetSelectedAtoms():
2655        page = getSelection()
2656        Ind = []
2657        if page:
2658            if G2frame.dataDisplay.GetPageText(page) == 'Draw Atoms':
2659                Ind = G2frame.dataDisplay.GetPage(page).GetSelectedRows()      #this is the Atoms grid in Draw Atoms
2660            elif G2frame.dataDisplay.GetPageText(page) == 'Map peaks':
2661                Ind = G2frame.dataDisplay.GetPage(page).GetSelectedRows()
2662            elif G2frame.dataDisplay.GetPageText(page) == 'Atoms':
2663                Ind = G2frame.dataDisplay.GetPage(page).GetSelectedRows()      #this is the Atoms grid in Atoms
2664        return Ind
2665                                       
2666    def SetBackground():
2667        R,G,B,A = Page.camera['backColor']
2668        glClearColor(R,G,B,A)
2669        glClear(GL_COLOR_BUFFER_BIT | GL_DEPTH_BUFFER_BIT)
2670       
2671    def SetLights():
2672        glEnable(GL_DEPTH_TEST)
2673        glShadeModel(GL_SMOOTH)
2674        glEnable(GL_LIGHTING)
2675        glEnable(GL_LIGHT0)
2676        glLightModeli(GL_LIGHT_MODEL_TWO_SIDE,0)
2677        glLightfv(GL_LIGHT0,GL_AMBIENT,[1,1,1,.8])
2678        glLightfv(GL_LIGHT0,GL_DIFFUSE,[1,1,1,1])
2679       
2680    def GetRoll(newxy,rho):
2681        Q = drawingData['Quaternion']
2682        dxy = G2mth.prodQVQ(G2mth.invQ(Q),np.inner(Bmat,newxy+[0,]))
2683        dxy = np.array(dxy*rho.shape)       
2684        roll = np.where(dxy>0.5,1,np.where(dxy<-.5,-1,0))
2685        return roll
2686               
2687    def SetMapRoll(newxy):
2688        rho = mapData['rho']
2689        roll = GetRoll(newxy,rho)
2690        mapData['rho'] = np.roll(np.roll(np.roll(rho,roll[0],axis=0),roll[1],axis=1),roll[2],axis=2)
2691        drawingData['oldxy'] = list(newxy)
2692       
2693    def SetPeakRoll(newxy):
2694        rho = mapData['rho']
2695        roll = GetRoll(newxy,rho)
2696        steps = 1./np.array(rho.shape)
2697        dxy = roll*steps
2698        for peak in mapPeaks:
2699            peak[1:4] += dxy
2700            peak[1:4] %= 1.
2701            peak[4] = np.sqrt(np.sum(np.inner(Amat,peak[1:4])**2))
2702               
2703    def SetTranslation(newxy):
2704#first get translation vector in screen coords.       
2705        oldxy = drawingData['oldxy']
2706        if not len(oldxy): oldxy = list(newxy)
2707        dxy = newxy-oldxy
2708        drawingData['oldxy'] = list(newxy)
2709        V = np.array([-dxy[0],dxy[1],0.])
2710#then transform to rotated crystal coordinates & apply to view point       
2711        Q = drawingData['Quaternion']
2712        V = np.inner(Bmat,G2mth.prodQVQ(G2mth.invQ(Q),V))
2713        Tx,Ty,Tz = drawingData['viewPoint'][0]
2714        Tx += V[0]*0.01
2715        Ty += V[1]*0.01
2716        Tz += V[2]*0.01
2717        drawingData['oldxy'] = list(newxy)
2718        drawingData['viewPoint'][0] =  Tx,Ty,Tz
2719        SetViewPointText([Tx,Ty,Tz])
2720       
2721    def SetRotation(newxy):
2722#first get rotation vector in screen coords. & angle increment       
2723        oldxy = drawingData['oldxy']
2724        if not len(oldxy): oldxy = list(newxy)
2725        dxy = newxy-oldxy
2726        drawingData['oldxy'] = list(newxy)
2727        V = np.array([dxy[1],dxy[0],0.])
2728        A = 0.25*np.sqrt(dxy[0]**2+dxy[1]**2)
2729# next transform vector back to xtal coordinates via inverse quaternion
2730# & make new quaternion
2731        Q = drawingData['Quaternion']
2732        V = G2mth.prodQVQ(G2mth.invQ(Q),np.inner(Bmat,V))
2733        DQ = G2mth.AVdeg2Q(A,V)
2734        Q = G2mth.prodQQ(Q,DQ)
2735        drawingData['Quaternion'] = Q
2736# finally get new view vector - last row of rotation matrix
2737        VD = np.inner(Bmat,G2mth.Q2Mat(Q)[2])
2738        VD /= np.sqrt(np.sum(VD**2))
2739        drawingData['viewDir'] = VD
2740        SetViewDirText(VD)
2741       
2742    def SetRotationZ(newxy):                       
2743#first get rotation vector (= view vector) in screen coords. & angle increment       
2744        View = glGetIntegerv(GL_VIEWPORT)
2745        cent = [View[2]/2,View[3]/2]
2746        oldxy = drawingData['oldxy']
2747        if not len(oldxy): oldxy = list(newxy)
2748        dxy = newxy-oldxy
2749        drawingData['oldxy'] = list(newxy)
2750        V = drawingData['viewDir']
2751        A = [0,0]
2752        A[0] = dxy[1]*.25
2753        A[1] = dxy[0]*.25
2754        if newxy[0] > cent[0]:
2755            A[0] *= -1
2756        if newxy[1] < cent[1]:
2757            A[1] *= -1       
2758# next transform vector back to xtal coordinates & make new quaternion
2759        Q = drawingData['Quaternion']
2760        V = np.inner(Amat,V)
2761        Qx = G2mth.AVdeg2Q(A[0],V)
2762        Qy = G2mth.AVdeg2Q(A[1],V)
2763        Q = G2mth.prodQQ(Q,Qx)
2764        Q = G2mth.prodQQ(Q,Qy)
2765        drawingData['Quaternion'] = Q
2766
2767    def RenderBox():
2768        glEnable(GL_COLOR_MATERIAL)
2769        glLineWidth(2)
2770        glEnable(GL_BLEND)
2771        glBlendFunc(GL_SRC_ALPHA,GL_ONE_MINUS_SRC_ALPHA)
2772        glEnable(GL_LINE_SMOOTH)
2773        glBegin(GL_LINES)
2774        for line,color in zip(uEdges,uColors):
2775            glColor3ubv(color)
2776            glVertex3fv(line[0])
2777            glVertex3fv(line[1])
2778        glEnd()
2779        glColor4ubv([0,0,0,0])
2780        glDisable(GL_LINE_SMOOTH)
2781        glDisable(GL_BLEND)
2782        glDisable(GL_COLOR_MATERIAL)
2783       
2784    def RenderUnitVectors(x,y,z):
2785        xyz = np.array([x,y,z])
2786        glEnable(GL_COLOR_MATERIAL)
2787        glLineWidth(1)
2788        glPushMatrix()
2789        glTranslate(x,y,z)
2790        glScalef(1/cell[0],1/cell[1],1/cell[2])
2791        glBegin(GL_LINES)
2792        for line,color in zip(uEdges,uColors)[:3]:
2793            glColor3ubv(color)
2794            glVertex3fv(-line[1]/2.)
2795            glVertex3fv(line[1]/2.)
2796        glEnd()
2797        glPopMatrix()
2798        glColor4ubv([0,0,0,0])
2799        glDisable(GL_COLOR_MATERIAL)
2800               
2801    def RenderSphere(x,y,z,radius,color):
2802        glMaterialfv(GL_FRONT_AND_BACK,GL_DIFFUSE,color)
2803        glPushMatrix()
2804        glTranslate(x,y,z)
2805        glMultMatrixf(B4mat.T)
2806        q = gluNewQuadric()
2807        gluSphere(q,radius,20,10)
2808        glPopMatrix()
2809       
2810    def RenderDots(XYZ,RC):
2811        glEnable(GL_COLOR_MATERIAL)
2812        XYZ = np.array(XYZ)
2813        glPushMatrix()
2814        for xyz,rc in zip(XYZ,RC):
2815            x,y,z = xyz
2816            r,c = rc
2817            glColor3ubv(c)
2818            glPointSize(r*50)
2819            glBegin(GL_POINTS)
2820            glVertex3fv(xyz)
2821            glEnd()
2822        glPopMatrix()
2823        glColor4ubv([0,0,0,0])
2824        glDisable(GL_COLOR_MATERIAL)
2825       
2826    def RenderSmallSphere(x,y,z,radius,color):
2827        glMaterialfv(GL_FRONT_AND_BACK,GL_DIFFUSE,color)
2828        glPushMatrix()
2829        glTranslate(x,y,z)
2830        glMultMatrixf(B4mat.T)
2831        q = gluNewQuadric()
2832        gluSphere(q,radius,4,2)
2833        glPopMatrix()
2834               
2835    def RenderEllipsoid(x,y,z,ellipseProb,E,R4,color):
2836        s1,s2,s3 = E
2837        glMaterialfv(GL_FRONT_AND_BACK,GL_DIFFUSE,color)
2838        glPushMatrix()
2839        glTranslate(x,y,z)
2840        glMultMatrixf(B4mat.T)
2841        glMultMatrixf(R4.T)
2842        glEnable(GL_NORMALIZE)
2843        glScale(s1,s2,s3)
2844        q = gluNewQuadric()
2845        gluSphere(q,ellipseProb,20,10)
2846        glDisable(GL_NORMALIZE)
2847        glPopMatrix()
2848       
2849    def RenderBonds(x,y,z,Bonds,radius,color,slice=20):
2850        glMaterialfv(GL_FRONT_AND_BACK,GL_DIFFUSE,color)
2851        glPushMatrix()
2852        glTranslate(x,y,z)
2853        glMultMatrixf(B4mat.T)
2854        for bond in Bonds:
2855            glPushMatrix()
2856            Dx = np.inner(Amat,bond)
2857            Z = np.sqrt(np.sum(Dx**2))
2858            if Z:
2859                azm = atan2d(-Dx[1],-Dx[0])
2860                phi = acosd(Dx[2]/Z)
2861                glRotate(-azm,0,0,1)
2862                glRotate(phi,1,0,0)
2863                q = gluNewQuadric()
2864                gluCylinder(q,radius,radius,Z,slice,2)
2865            glPopMatrix()           
2866        glPopMatrix()
2867               
2868    def RenderLines(x,y,z,Bonds,color):
2869        xyz = np.array([x,y,z])
2870        glEnable(GL_COLOR_MATERIAL)
2871        glLineWidth(1)
2872        glColor3fv(color)
2873        glPushMatrix()
2874        glBegin(GL_LINES)
2875        for bond in Bonds:
2876            glVertex3fv(xyz)
2877            glVertex3fv(xyz+bond)
2878        glEnd()
2879        glColor4ubv([0,0,0,0])
2880        glPopMatrix()
2881        glDisable(GL_COLOR_MATERIAL)
2882       
2883    def RenderPolyhedra(x,y,z,Faces,color):
2884        glPushMatrix()
2885        glTranslate(x,y,z)
2886        glMaterialfv(GL_FRONT_AND_BACK,GL_DIFFUSE,color)
2887        glShadeModel(GL_SMOOTH)
2888        glMultMatrixf(B4mat.T)
2889        for face,norm in Faces:
2890            glPolygonMode(GL_FRONT_AND_BACK,GL_FILL)
2891            glFrontFace(GL_CW)
2892            glNormal3fv(norm)
2893            glBegin(GL_TRIANGLES)
2894            for vert in face:
2895                glVertex3fv(vert)
2896            glEnd()
2897        glPopMatrix()
2898
2899    def RenderMapPeak(x,y,z,color,den):
2900        xyz = np.array([x,y,z])
2901        glEnable(GL_COLOR_MATERIAL)
2902        glLineWidth(3)
2903        glColor3fv(color*den/255)
2904        glPushMatrix()
2905        glBegin(GL_LINES)
2906        for vec in mapPeakVecs:
2907            glVertex3fv(vec[0]+xyz)
2908            glVertex3fv(vec[1]+xyz)
2909        glEnd()
2910        glColor4ubv([0,0,0,0])
2911        glPopMatrix()
2912        glDisable(GL_COLOR_MATERIAL)
2913       
2914    def RenderBackbone(Backbone,BackboneColor,radius):
2915        glPushMatrix()
2916        glMultMatrixf(B4mat.T)
2917        glEnable(GL_COLOR_MATERIAL)
2918        glShadeModel(GL_SMOOTH)
2919        gleSetJoinStyle(TUBE_NORM_EDGE | TUBE_JN_ANGLE | TUBE_JN_CAP)
2920        glePolyCylinder(Backbone,BackboneColor,radius)
2921        glPopMatrix()       
2922        glDisable(GL_COLOR_MATERIAL)
2923       
2924    def RenderLabel(x,y,z,label,r):       
2925        glPushMatrix()
2926        glTranslate(x,y,z)
2927        glMultMatrixf(B4mat.T)
2928        glDisable(GL_LIGHTING)
2929        glColor3f(0,1.,0)
2930        glRasterPos3f(r,r,r)
2931        for c in list(label):
2932            glutBitmapCharacter(GLUT_BITMAP_8_BY_13,ord(c))
2933        glEnable(GL_LIGHTING)
2934        glPopMatrix()
2935       
2936    def RenderMap(rho,rhoXYZ,indx,Rok):
2937        cLevel = drawingData['contourLevel']
2938        XYZ = []
2939        RC = []
2940        for i,xyz in enumerate(rhoXYZ):
2941            if not Rok[i]:
2942                x,y,z = xyz
2943                I,J,K = indx[i]
2944                alpha = 1.0
2945                if cLevel < 1.:
2946                    alpha = (abs(rho[I,J,K])/mapData['rhoMax']-cLevel)/(1.-cLevel)
2947                if rho[I,J,K] < 0.:
2948                    XYZ.append(xyz)
2949                    RC.append([0.1*alpha,Rd])
2950                else:
2951                    XYZ.append(xyz)
2952                    RC.append([0.1*alpha,Gr])
2953        RenderDots(XYZ,RC)
2954                           
2955    def Draw(caller=''):
2956#useful debug?       
2957#        if caller:
2958#            print caller
2959# end of useful debug
2960        mapData = generalData['Map']
2961        pageName = ''
2962        page = getSelection()
2963        if page:
2964            pageName = G2frame.dataDisplay.GetPageText(page)
2965        rhoXYZ = []
2966        if len(mapData['rho']):
2967            VP = np.array(drawingData['viewPoint'][0])-np.array([.5,.5,.5])
2968            contLevel = drawingData['contourLevel']*mapData['rhoMax']
2969            if 'delt-F' in mapData['MapType']:
2970                rho = ma.array(mapData['rho'],mask=(np.abs(mapData['rho'])<contLevel))
2971            else:
2972                rho = ma.array(mapData['rho'],mask=(mapData['rho']<contLevel))
2973            steps = 1./np.array(rho.shape)
2974            incre = np.where(VP>=0,VP%steps,VP%steps-steps)
2975            Vsteps = -np.array(VP/steps,dtype='i')
2976            rho = np.roll(np.roll(np.roll(rho,Vsteps[0],axis=0),Vsteps[1],axis=1),Vsteps[2],axis=2)
2977            indx = np.array(ma.nonzero(rho)).T
2978            rhoXYZ = indx*steps+VP-incre
2979            Nc = len(rhoXYZ)
2980            rcube = 2000.*Vol/(ForthirdPI*Nc)
2981            rmax = math.exp(math.log(rcube)/3.)**2
2982            radius = min(drawingData['mapSize']**2,rmax)
2983            view = np.array(drawingData['viewPoint'][0])
2984            Rok = np.sum(np.inner(Amat,rhoXYZ-view).T**2,axis=1)>radius
2985        Ind = GetSelectedAtoms()
2986        VS = np.array(Page.canvas.GetSize())
2987        aspect = float(VS[0])/float(VS[1])
2988        cPos = drawingData['cameraPos']
2989        Zclip = drawingData['Zclip']*cPos/200.
2990        Q = drawingData['Quaternion']
2991        Tx,Ty,Tz = drawingData['viewPoint'][0]
2992        cx,ct,cs,ci = drawingData['atomPtrs']
2993        bondR = drawingData['bondRadius']
2994        G,g = G2lat.cell2Gmat(cell)
2995        GS = G
2996        GS[0][1] = GS[1][0] = math.sqrt(GS[0][0]*GS[1][1])
2997        GS[0][2] = GS[2][0] = math.sqrt(GS[0][0]*GS[2][2])
2998        GS[1][2] = GS[2][1] = math.sqrt(GS[1][1]*GS[2][2])
2999        ellipseProb = G2lat.criticalEllipse(drawingData['ellipseProb']/100.)
3000       
3001        SetBackground()
3002        glInitNames()
3003        glPushName(0)
3004       
3005        glMatrixMode(GL_PROJECTION)
3006        glLoadIdentity()
3007        glViewport(0,0,VS[0],VS[1])
3008        gluPerspective(20.,aspect,cPos-Zclip,cPos+Zclip)
3009        gluLookAt(0,0,cPos,0,0,0,0,1,0)
3010        SetLights()           
3011           
3012        glMatrixMode(GL_MODELVIEW)
3013        glLoadIdentity()
3014        matRot = G2mth.Q2Mat(Q)
3015        matRot = np.concatenate((np.concatenate((matRot,[[0],[0],[0]]),axis=1),[[0,0,0,1],]),axis=0)
3016        glMultMatrixf(matRot.T)
3017        glMultMatrixf(A4mat.T)
3018        glTranslate(-Tx,-Ty,-Tz)
3019        if drawingData['unitCellBox']:
3020            RenderBox()
3021        if drawingData['showABC']:
3022            x,y,z = drawingData['viewPoint'][0]
3023            RenderUnitVectors(x,y,z)
3024        Backbones = {}
3025        BackboneColor = []
3026        time0 = time.time()
3027#        glEnable(GL_BLEND)
3028#        glBlendFunc(GL_SRC_ALPHA,GL_ONE_MINUS_SRC_ALPHA)
3029        for iat,atom in enumerate(drawingData['Atoms']):
3030            x,y,z = atom[cx:cx+3]
3031            Bonds = atom[-2]
3032            Faces = atom[-1]
3033            try:
3034                atNum = generalData['AtomTypes'].index(atom[ct])
3035            except ValueError:
3036                atNum = -1
3037            CL = atom[cs+2]
3038            color = np.array(CL)/255.
3039            if iat in Ind and G2frame.dataDisplay.GetPageText(getSelection()) != 'Map peaks':
3040                color = np.array(Gr)/255.
3041#            color += [.25,]
3042            radius = 0.5
3043            if atom[cs] != '':
3044                try:
3045                    glLoadName(atom[-3])
3046                except: #problem with old files - missing code
3047                    pass                   
3048            if 'balls' in atom[cs]:
3049                vdwScale = drawingData['vdwScale']
3050                ballScale = drawingData['ballScale']
3051                if atNum < 0:
3052                    radius = 0.2
3053                elif 'H' == atom[ct]:
3054                    if drawingData['showHydrogen']:
3055                        if 'vdW' in atom[cs] and atNum >= 0:
3056                            radius = vdwScale*generalData['vdWRadii'][atNum]
3057                        else:
3058                            radius = ballScale*drawingData['sizeH']
3059                    else:
3060                        radius = 0.0
3061                else:
3062                    if 'vdW' in atom[cs]:
3063                        radius = vdwScale*generalData['vdWRadii'][atNum]
3064                    else:
3065                        radius = ballScale*generalData['BondRadii'][atNum]
3066                RenderSphere(x,y,z,radius,color)
3067                if 'sticks' in atom[cs]:
3068                    RenderBonds(x,y,z,Bonds,bondR,color)
3069            elif 'ellipsoids' in atom[cs]:
3070                RenderBonds(x,y,z,Bonds,bondR,color)
3071                if atom[cs+3] == 'A':                   
3072                    Uij = atom[cs+5:cs+11]
3073                    U = np.multiply(G2spc.Uij2U(Uij),GS)
3074                    U = np.inner(Amat,np.inner(U,Amat).T)
3075                    E,R = nl.eigh(U)
3076                    R4 = np.concatenate((np.concatenate((R,[[0],[0],[0]]),axis=1),[[0,0,0,1],]),axis=0)
3077                    E = np.sqrt(E)
3078                    if atom[ct] == 'H' and not drawingData['showHydrogen']:
3079                        pass
3080                    else:
3081                        RenderEllipsoid(x,y,z,ellipseProb,E,R4,color)                   
3082                else:
3083                    if atom[ct] == 'H' and not drawingData['showHydrogen']:
3084                        pass
3085                    else:
3086                        radius = ellipseProb*math.sqrt(abs(atom[cs+4]))
3087                        RenderSphere(x,y,z,radius,color)
3088            elif 'lines' in atom[cs]:
3089                radius = 0.1
3090                RenderLines(x,y,z,Bonds,color)
3091#                RenderBonds(x,y,z,Bonds,0.05,color,6)
3092            elif atom[cs] == 'sticks':
3093                radius = 0.1
3094                RenderBonds(x,y,z,Bonds,bondR,color)
3095            elif atom[cs] == 'polyhedra':
3096                RenderPolyhedra(x,y,z,Faces,color)
3097            elif atom[cs] == 'backbone':
3098                if atom[ct-1].split()[0] in ['C','N']:
3099                    if atom[2] not in Backbones:
3100                        Backbones[atom[2]] = []
3101                    Backbones[atom[2]].append(list(np.inner(Amat,np.array([x,y,z]))))
3102                    BackboneColor.append(list(color))
3103                   
3104            if atom[cs+1] == 'type':
3105                RenderLabel(x,y,z,atom[ct],radius)
3106            elif atom[cs+1] == 'name':
3107                RenderLabel(x,y,z,atom[ct-1],radius)
3108            elif atom[cs+1] == 'number':
3109                RenderLabel(x,y,z,str(iat),radius)
3110            elif atom[cs+1] == 'residue' and atom[ct-1] == 'CA':
3111                RenderLabel(x,y,z,atom[ct-4],radius)
3112            elif atom[cs+1] == '1-letter' and atom[ct-1] == 'CA':
3113                RenderLabel(x,y,z,atom[ct-3],radius)
3114            elif atom[cs+1] == 'chain' and atom[ct-1] == 'CA':
3115                RenderLabel(x,y,z,atom[ct-2],radius)
3116#        glDisable(GL_BLEND)
3117        if len(rhoXYZ):
3118            RenderMap(rho,rhoXYZ,indx,Rok)
3119        if len(mapPeaks):
3120            for ind,[mag,x,y,z,d] in enumerate(mapPeaks):
3121                if ind in Ind and pageName == 'Map peaks':
3122                    RenderMapPeak(x,y,z,Gr,1.0)
3123                else:
3124                    RenderMapPeak(x,y,z,Wt,mag/peakMax)
3125        if Backbones:
3126            for chain in Backbones:
3127                Backbone = Backbones[chain]
3128                RenderBackbone(Backbone,BackboneColor,bondR)
3129#        print time.time()-time0
3130        Page.canvas.SwapBuffers()
3131       
3132    def OnSize(event):
3133        Draw('size')
3134       
3135    def OnFocus(event):         #not needed?? Bind commented out below
3136        Draw('focus')
3137       
3138    try:
3139        plotNum = G2frame.G2plotNB.plotList.index(generalData['Name'])
3140        Page = G2frame.G2plotNB.nb.GetPage(plotNum)       
3141    except ValueError:
3142        Plot = G2frame.G2plotNB.addOgl(generalData['Name'])
3143        plotNum = G2frame.G2plotNB.plotList.index(generalData['Name'])
3144        Page = G2frame.G2plotNB.nb.GetPage(plotNum)
3145        Page.views = False
3146        view = False
3147        altDown = False
3148    Page.SetFocus()
3149    Page.Choice = None
3150    if mapData['Flip']:
3151        choice = [' save as/key:','jpeg','tiff','bmp','u: roll up','d: roll down','l: roll left','r: roll right']
3152    else:
3153        choice = [' save as/key:','jpeg','tiff','bmp','c: center on 1/2,1/2,1/2','n: next','p: previous']
3154    cb = wx.ComboBox(G2frame.G2plotNB.status,style=wx.CB_DROPDOWN|wx.CB_READONLY,choices=choice)
3155    cb.Bind(wx.EVT_COMBOBOX, OnKeyBox)
3156    cb.SetValue(' save as/key:')
3157    Page.canvas.Bind(wx.EVT_MOUSEWHEEL, OnMouseWheel)
3158    Page.canvas.Bind(wx.EVT_LEFT_DOWN, OnMouseDown)
3159    Page.canvas.Bind(wx.EVT_RIGHT_DOWN, OnMouseDown)
3160    Page.canvas.Bind(wx.EVT_MIDDLE_DOWN, OnMouseDown)
3161    Page.canvas.Bind(wx.EVT_KEY_UP, OnKey)
3162    Page.canvas.Bind(wx.EVT_MOTION, OnMouseMove)
3163    Page.canvas.Bind(wx.EVT_SIZE, OnSize)
3164#    Page.canvas.Bind(wx.EVT_SET_FOCUS, OnFocus)
3165    Page.camera['position'] = drawingData['cameraPos']
3166    Page.camera['viewPoint'] = np.inner(Amat,drawingData['viewPoint'][0])
3167    Page.camera['backColor'] = np.array(list(drawingData['backColor'])+[0,])/255.
3168    Page.canvas.SetCurrent()
3169    Draw('main')
3170       
Note: See TracBrowser for help on using the repository browser.