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