Changeset 4479

Timestamp:

Jun 10, 2020 9:10:31 PM (3 years ago)

Author:

toby

Message:

changes to unique map peaks & plotting after selection

Location:

trunk

Files:

• GSASIImath.py (modified) (2 diffs)
• GSASIIphsGUI.py (modified) (12 diffs)
• GSASIIplot.py (modified) (2 diffs)

trunk/GSASIImath.py

@@ -3934 +3934 @@
                 
 def PeaksUnique(data,Ind):
-    '''Finds the symmetry unique set of peaks from those selected. Works on the 
-    contents of data['Map Peaks']. 
+    '''Finds the symmetry unique set of peaks from those selected. Selects
+    the one closest to the center of the unit cell. 
+    Works on the contents of data['Map Peaks']. Called from OnPeaksUnique in 
+    GSASIIphsGUI.py,
 
     :param data: the phase data structure
     :param list Ind: list of selected peak indices
+
     :returns: the list of symmetry unique peaks from among those given in Ind
-
     '''        
 #    XYZE = np.array([[equiv[0] for equiv in G2spc.GenAtom(xyz[1:4],SGData,Move=True)] for xyz in mapPeaks]) #keep this!!
@@ -3953 +3955 @@
     SGData = generalData['SGData']
     mapPeaks = data['Map Peaks']
-    Indx = {}
-    XYZ = {}
-    for ind in Ind:
-        XYZ[ind] = np.array(mapPeaks[ind][1:4])
-        Indx[ind] = True
+    XYZ = {ind:np.array(mapPeaks[ind][1:4]) for ind in Ind}
+    Indx = [True for ind in Ind]
+    Unique = []
+    # scan through peaks, finding all peaks equivalent to peak ind
     for ind in Ind:
         if Indx[ind]:
             xyz = XYZ[ind]
+            dups = []
             for jnd in Ind:
-                if ind != jnd and Indx[jnd]:                        
+                # only consider peaks we have not looked at before
+                # and were not already found to be equivalent
+                if jnd > ind and Indx[jnd]:                        
                     Equiv = G2spc.GenAtom(XYZ[jnd],SGData,Move=True)
                     xyzs = np.array([equiv[0] for equiv in Equiv])
-                    Indx[jnd] = noDuplicate(xyz,xyzs,Amat)
-    Ind = []
-    for ind in Indx:
-        if Indx[ind]:
-            Ind.append(ind)
-    return Ind
+                    if not noDuplicate(xyz,xyzs,Amat):
+                        Indx[jnd] = False
+                        dups.append(jnd)
+            # select the unique peak closest to cell center
+            cntr = mapPeaks[ind][-1]
+            icntr = ind
+            for jnd in dups:
+                if mapPeaks[jnd][-1] < cntr:
+                    cntr = mapPeaks[jnd][-1]
+                    icntr = jnd
+            Unique.append(icntr)
+    return Unique
 
 ################################################################################

trunk/GSASIIphsGUI.py

@@ -1403 +1403 @@
             'vdwScale':0.67,'ellipseProb':50,'sizeH':0.50,'unitCellBox':True,
             'showABC':True,'selectedAtoms':[],'Atoms':[],'oldxy':[],'magMult':1.0,
-            'bondList':{},'viewDir':[1,0,0],'Plane':[[0,0,1],False,False,0.0,[255,255,0]]}
+            'bondList':{},'viewDir':[1,0,0],'Plane':[[0,0,1],False,False,0.0,[255,255,0]],
+            'peakMoveView':True,'PeakDistRadius':0.0,
+            'atomsExpandRadius':5.,'atomsDistRadius':2.5,
+        }
     for key in defaultDrawing:
         if key not in drawingData: drawingData[key] = defaultDrawing[key]
@@ -8727 +8730 @@
             drawOptions.GetSizer().Clear(True)
         mainSizer = wx.BoxSizer(wx.VERTICAL)
-        mainSizer.Add((5,5),0)
-        mainSizer.Add(wx.StaticText(drawOptions,-1,' Drawing controls:'),0,WACV)
-        mainSizer.Add((5,5),0)        
+        #mainSizer.Add(wx.StaticText(drawOptions,-1,' Drawing controls'),0,
+        #                  wx.ALL|wx.CENTER)
+        #G2G.HorizontalLine(mainSizer,drawOptions)
         mainSizer.Add(SlopSizer(),0)
-        mainSizer.Add((5,5),0)
+        G2G.HorizontalLine(mainSizer,drawOptions)
         mainSizer.Add(ShowSizer(),0,)
-        mainSizer.Add((5,5),0)
+        G2G.HorizontalLine(mainSizer,drawOptions)
         mainSizer.Add(RadSizer(),0,)
-        mainSizer.Add((5,5),0)
+        G2G.HorizontalLine(mainSizer,drawOptions)
         mainSizer.Add(PlaneSizer(),0,)
+        G2G.HorizontalLine(mainSizer,drawOptions)
+        mainSizer.Add(wx.StaticText(drawOptions,wx.ID_ANY,
+                                        'On map peak selection:'),0,WACV)
+        mainSizer.Add(G2G.G2CheckBox(drawOptions,'Move view point',
+                                         drawingData,'peakMoveView'))
+        mapSizer = wx.FlexGridSizer(0,3,5,5)
+        mapSizer.Add(wx.StaticText(drawOptions,wx.ID_ANY,'Show Map points within:'),0,WACV)
+        mapSizer.Add(G2G.ValidatedTxtCtrl(drawOptions,drawingData,'PeakDistRadius',
+                            min=0.0,max=5.0,nDig=(10,1),size=(50,-1)))
+        mapSizer.Add(wx.StaticText(drawOptions,wx.ID_ANY,u"\u212B"),0,WACV)
+        mapSizer.Add(wx.StaticText(drawOptions,wx.ID_ANY,'Show atoms within:'),0,WACV)
+        mapSizer.Add(G2G.ValidatedTxtCtrl(drawOptions,drawingData,'atomsExpandRadius',
+                            min=0.0,max=15.0,nDig=(10,1),size=(50,-1)))
+        mapSizer.Add(wx.StaticText(drawOptions,wx.ID_ANY,u"\u212B"),0,WACV)
+        mapSizer.Add(wx.StaticText(drawOptions,wx.ID_ANY,'Label distance to atoms within:'),0,WACV)
+        mapSizer.Add(G2G.ValidatedTxtCtrl(drawOptions,drawingData,'atomsDistRadius',
+                            min=0.0,max=15.0,nDig=(10,1),size=(50,-1)))
+        mapSizer.Add(wx.StaticText(drawOptions,wx.ID_ANY,u"\u212B"),0,WACV)
+        mainSizer.Add(mapSizer)
         SetPhaseWindow(drawOptions,mainSizer)
 
@@ -11470 +11492 @@
             elif c < 0:                   #only row clicks
                 if event.ControlDown():                    
-                    if r in MapPeaks.GetSelectedRows():
+                    if r in getAtomSelections(MapPeaks):
                         MapPeaks.DeselectRow(r)
                     else:
                         MapPeaks.SelectRow(r,True)
                 elif event.ShiftDown():
-                    indxs = MapPeaks.GetSelectedRows()
+                    indxs = getAtomSelections(MapPeaks)
                     MapPeaks.ClearSelection()
                     ibeg = 0
@@ -11498 +11520 @@
             G2plt.PlotStructure(G2frame,data)                    
             
+        # beginning of FillMapPeaksGrid()
+        #import imp   # debug code
+        #imp.reload(G2plt) # debug code
         G2frame.GetStatusBar().SetStatusText('',1)
         if 'Map Peaks' in data:
@@ -11520 +11545 @@
             mapPeaks = np.array(data['Map Peaks'])
             peakMax = np.amax(mapPeaks.T[0])
-            Ind = MapPeaks.GetSelectedRows()
+            Ind = getAtomSelections(MapPeaks)
             for ind in Ind:
                 mag,x,y,z = mapPeaks[ind][:4]
@@ -11556 +11581 @@
         if 'Map Peaks' in data:
             mapPeaks = data['Map Peaks']
-            Ind = MapPeaks.GetSelectedRows()
+            Ind = getAtomSelections(MapPeaks)
             Ind.sort()
             Ind.reverse()
@@ -11581 +11606 @@
     def OnPeaksEquiv(event):
         if 'Map Peaks' in data:
-            Ind = MapPeaks.GetSelectedRows()
+            Ind = getAtomSelections(MapPeaks)
             if Ind:
                 wx.BeginBusyCursor()
@@ -11609 +11634 @@
         if 'Map Peaks' in data:
             mapPeaks = data['Map Peaks']
-            Ind = MapPeaks.GetSelectedRows()
+            Ind = getAtomSelections(MapPeaks)
             if Ind:
                 wx.BeginBusyCursor()
@@ -11615 +11640 @@
                     Ind = G2mth.PeaksUnique(data,Ind)
                     print (' No. unique peaks: %d Unique peak fraction: %.3f'%(len(Ind),float(len(Ind))/len(mapPeaks)))
+                    tbl = MapPeaks.GetTable().data
+                    tbl[:] = [t for i,t in enumerate(tbl) if i in Ind] + [
+                        t for i,t in enumerate(tbl) if i not in Ind]         
                     for r in range(MapPeaks.GetNumberRows()):
-                        if r in Ind:
+                        if r < len(Ind):
                             MapPeaks.SelectRow(r,addToSelected=True)
                         else:
@@ -11622 +11650 @@
                 finally:
                     wx.EndBusyCursor()
+                MapPeaks.ForceRefresh()
                 G2plt.PlotStructure(G2frame,data)
                 
     def OnPeaksViewPoint(event):
         # set view point
-        indx = MapPeaks.GetSelectedRows()
+        indx = getAtomSelections(MapPeaks)
         if not indx:
             G2frame.ErrorDialog('Set viewpoint','No peaks selected')
@@ -11637 +11666 @@
     def OnPeaksDistVP(event):
         # distance to view point
-        indx = MapPeaks.GetSelectedRows()
+        indx = getAtomSelections(MapPeaks)
         if not indx:
             G2frame.ErrorDialog('Peak distance','No peaks selected')
@@ -11658 +11687 @@
     def OnPeaksDA(event):
         #distance, angle 
-        indx = MapPeaks.GetSelectedRows()
+        indx = getAtomSelections(MapPeaks)
         if len(indx) not in [2,3]:
             G2frame.ErrorDialog('Peak distance/angle','Wrong number of atoms for distance or angle calculation')

trunk/GSASIIplot.py

@@ -8928 +8928 @@
         Zclip = drawingData['Zclip']*cPos/200.
         Q = drawingData['Quaternion']
+        # Move view point on selection of a single peak, if enabled
+        if pageName == 'Map peaks' and len(Ind) == 1 and drawingData['peakMoveView']:
+            ind = Ind[0]
+            mag,x,y,z = mapPeaks[:,:4][ind]
+            drawingData['viewPoint'][0][:] = [x,y,z]
         Tx,Ty,Tz = drawingData['viewPoint'][0]
         cx,ct,cs,ci = drawingData['atomPtrs']
@@ -9078 +9083 @@
         if not FourD and len(rhoXYZ):       #no green dot map for 4D - it's wrong!
             RenderMap(rho,rhoXYZ,indx,Rok)
-        if len(mapPeaks):
+        if len(Ind) == 1 and pageName == 'Map peaks':
+            # one peak has been selected, show as selected by draw options
+            PeakDistRadius = drawingData['PeakDistRadius']
+            atomsExpandRadius = drawingData['atomsExpandRadius']
+            atomsdistRadius = drawingData['atomsDistRadius']
+            ind = Ind[0]
+            mag,x,y,z = mapPeaks[:,:4][ind]
+            RenderMapPeak(x,y,z,Gr,1.0)
+            # distances to other peaks within PeakDistRadius A
+            if PeakDistRadius > 0:
+                XYZ = mapPeaks.T[1:4].T
+                Dx = XYZ-(x,y,z)
+                dist = np.sqrt(np.sum(np.inner(Dx,Amat)**2,axis=1))
+                for i in ma.nonzero(ma.masked_greater(dist,PeakDistRadius))[0]:
+                    lbl = '{:.2f}'.format(dist[i])
+                    RenderLines(x,y,z,[Dx[i]],Gr)
+                    lx,ly,lz = (x,y,z)+(Dx[i]/2)
+                    RenderLabel(lx,ly,lz,lbl,radius,wxGreen,matRot)
+                    [mag,x1,y1,z1] = mapPeaks[:,:4][i]
+                    if mag > 0.:
+                        RenderMapPeak(x1,y1,z1,Wt,mag/peakMax)
+                    else:
+                        RenderMapPeak(x1,y1,z1,Or,-2*mag/peakMax)
+            # find and display atoms within atomsExpandRadius A
+            if atomsExpandRadius > 0:
+                xyzA = np.array((x,y,z))
+                cx,ct,cs,ci = G2phG.getAtomPtrs(data)      
+                cellArray = G2lat.CellBlock(1)
+                radDict = dict(zip(*G2phG.getAtomRadii(data)))
+                Names = []
+                Atoms = []
+                Radii = []
+                Color = []
+                for atomB in data['Atoms']:
+                    xyzB = np.array(atomB[cx:cx+3])
+                    color = np.array(generalData['Color'][generalData['AtomTypes'].index(atomB[ct])])/255.
+                    result = G2spc.GenAtom(xyzB,SGData,False,6*[0.],True)
+                    for item in result:
+                        for xyz in cellArray+np.array(item[0]):
+                            dist = np.sqrt(np.sum(np.inner(Amat,xyz-xyzA)**2))
+                            if 0 < dist <= max(atomsdistRadius,atomsExpandRadius):
+                                ax,ay,az = xyz
+                                RenderSphere(ax,ay,az,radius,color)
+                                Atoms.append(xyz)
+                                Radii.append(radDict[atomB[ct]])
+                                Names.append(atomB[0])
+                                Color.append(color)
+                            if dist < atomsdistRadius:
+                                RenderLines(x,y,z,[xyz-xyzA],Gr)
+                                lx,ly,lz = (xyzA+xyz)/2
+                                lbl = '{:.2f}'.format(dist)
+                                RenderLabel(lx,ly,lz,lbl,radius,wxGreen,matRot)
+                                ax,ay,az = xyz
+                                RenderLabel(ax,ay,az,'  '+atomB[0],radius,wxGreen,matRot)
+                # find bonds
+                bondData = [[] for i in range(len(Atoms))]
+                Atoms = np.array(Atoms)
+                Radii = np.array(Radii)
+                for i,atom in enumerate(Atoms):
+                    Dx = Atoms-atom
+                    sumR = Radii + Radii[i]
+                    dist = ma.masked_less(np.sqrt(np.sum(np.inner(Amat,Dx)**2,axis=0)),0.5)
+                    IndB = ma.nonzero(ma.masked_greater(dist-data['Drawing']['radiusFactor']*sumR,0.))                 #get indices of bonded atoms
+                    for j in IndB[0]:
+                        bondData[i].append(Dx[j]*Radii[i]/sumR[j])
+                        bondData[j].append(-Dx[j]*Radii[j]/sumR[j])
+                # draw bonds
+                bondR = drawingData['bondRadius']
+                for i,xyz in enumerate(Atoms):
+                    ax,ay,az = xyz
+                    RenderBonds(ax,ay,az,bondData[i],bondR,Color[i])
+        elif len(mapPeaks):
             XYZ = mapPeaks.T[1:4].T
             mapBonds = FindPeaksBonds(XYZ)