Changeset 4291 for trunk/GSASIImath.py
 Timestamp:
 Feb 10, 2020 10:06:35 AM (3 years ago)
 File:

 1 edited
Legend:
 Unmodified
 Added
 Removed

trunk/GSASIImath.py
r4289 r4291 639 639 Ids.append(Atoms[j][cia+8]) 640 640 return Neigh,[OId,Ids] 641 642 def FindOctahedron(results): 643 Octahedron = np.array([[1.,0,0],[0,1.,0],[0,0,1.],[1.,0,0],[0,1.,0],[0,0,1.]]) 644 Polygon = np.array([result[3] for result in results]) 645 Dists = np.array([np.sqrt(np.sum(axis**2)) for axis in Polygon]) 646 bond = np.mean(Dists) 647 std = np.std(Dists) 648 Norms = Polygon/Dists[:,nxs] 649 Tilts = acosd(np.dot(Norms,Octahedron[0])) 650 iAng = np.argmin(Tilts) 651 Qavec = np.cross(Norms[iAng],Octahedron[0]) 652 QA = AVdeg2Q(Tilts[iAng],Qavec) 653 newNorms = prodQVQ(QA,Norms) 654 Rots = acosd(np.dot(newNorms,Octahedron[1])) 655 jAng = np.argmin(Rots) 656 Qbvec = np.cross(Norms[jAng],Octahedron[1]) 657 QB = AVdeg2Q(Rots[jAng],Qbvec) 658 QQ = prodQQ(QA,QB) 659 newNorms = prodQVQ(QQ,Norms) 660 dispVecs = np.array([norm[:,nxs]Octahedron.T for norm in newNorms]) 661 disp = np.sqrt(np.sum(dispVecs**2,axis=1)) 662 dispids = np.argmin(disp,axis=1) 663 vecDisp = np.array([dispVecs[i,:,dispid] for i,dispid in enumerate(dispids)]) 664 Disps = np.array([disp[i,dispid] for i,dispid in enumerate(dispids)]) 665 meanDisp = np.mean(Disps) 666 stdDisp = np.std(Disps) 667 A,V = Q2AVdeg(QQ) 668 return bond,std,meanDisp,stdDisp,A,V,vecDisp 669 670 def FindTetrahedron(results): 671 Tetrahedron = np.array([[1.,1,1],[1,1,1],[1,1,1],[1,1,1]])/np.sqrt(3.) 672 Polygon = np.array([result[3] for result in results]) 673 Dists = np.array([np.sqrt(np.sum(axis**2)) for axis in Polygon]) 674 bond = np.mean(Dists) 675 std = np.std(Dists) 676 Norms = Polygon/Dists[:,nxs] 677 Tilts = acosd(np.dot(Norms,Tetrahedron[0])) 678 iAng = np.argmin(Tilts) 679 Qavec = np.cross(Norms[iAng],Tetrahedron[0]) 680 QA = AVdeg2Q(Tilts[iAng],Qavec) 681 newNorms = prodQVQ(QA,Norms) 682 Rots = acosd(np.dot(newNorms,Tetrahedron[1])) 683 jAng = np.argmin(Rots) 684 Qbvec = np.cross(Norms[jAng],Tetrahedron[1]) 685 QB = AVdeg2Q(Rots[jAng],Qbvec) 686 QQ = prodQQ(QA,QB) 687 newNorms = prodQVQ(QQ,Norms) 688 dispVecs = np.array([norm[:,nxs]Tetrahedron.T for norm in newNorms]) 689 disp = np.sqrt(np.sum(dispVecs**2,axis=1)) 690 dispids = np.argmin(disp,axis=1) 691 vecDisp = np.array([dispVecs[i,:,dispid] for i,dispid in enumerate(dispids)]) 692 Disps = np.array([disp[i,dispid] for i,dispid in enumerate(dispids)]) 693 meanDisp = np.mean(Disps) 694 stdDisp = np.std(Disps) 695 A,V = Q2AVdeg(QQ) 696 return bond,std,meanDisp,stdDisp,A,V,vecDisp 641 697 642 698 def FindAllNeighbors(phase,FrstName,AtNames,notName=''): … … 645 701 Atoms = phase['Atoms'] 646 702 atNames = [atom[ct1] for atom in Atoms] 703 atTypes = [atom[ct] for atom in Atoms] 647 704 Cell = General['Cell'][1:7] 648 705 Amat,Bmat = G2lat.cell2AB(Cell) … … 650 707 indices = (1,0,1) 651 708 Units = np.array([[h,k,l] for h in indices for k in indices for l in indices]) 652 atTypes = General['AtomTypes']709 AtTypes = General['AtomTypes'] 653 710 Radii = np.array(General['BondRadii']) 654 711 DisAglCtls = General['DisAglCtls'] 655 712 radiusFactor = DisAglCtls['Factors'][0] 656 AtInfo = dict(zip( atTypes,Radii)) #or General['BondRadii']713 AtInfo = dict(zip(AtTypes,Radii)) #or General['BondRadii'] 657 714 Orig = atNames.index(FrstName) 658 715 OId = Atoms[Orig][cia+8] … … 681 738 else: 682 739 Topstr = ' +(%4d)'%(Top) 683 Neigh.append([AtNames[iA]+Topstr, dist[iU]])740 Neigh.append([AtNames[iA]+Topstr,atTypes[iA],dist[iU],dx[iU]]) 684 741 Ids.append(Atoms[iA][cia+8]) 685 742 return Neigh,[OId,Ids]
Note: See TracChangeset
for help on using the changeset viewer.