# Changeset 477

Ignore:
Timestamp:
Feb 9, 2012 1:23:01 PM (12 years ago)
Message:

final fixes to dist-angle-torsion calc.

Location:
trunk
Files:
5 edited

Unmodified
Removed
• ## trunk/GSASIImath.py

 r470 def calcDist(Ox,Tx,U,inv,C,M,T,Amat): TxT = inv*(np.inner(M,Tx)+T)+C TxT = G2spc.MoveToUnitCell(TxT)+U TxT = inv*(np.inner(M,Tx)+T)+C+U return np.sqrt(np.sum(np.inner(Amat,(TxT-Ox))**2)) return calcAngle(OxA,TxA,TxB,unitA,unitB,invA,CA,MA,TA,invB,CB,MB,TB,Amat),0.0 def GetTorsionSig(Atoms,Amat,SGData,covData={}): XYZ = [] for atom in Atoms: XYZ.append(np.array(atom[3:6])) XYZ = np.inner(Amat,XYZ).T V1 = XYZ[1]-XYZ[0] V2 = XYZ[2]-XYZ[1] V3 = XYZ[3]-XYZ[2] V1 /= np.sqrt(np.sum(V1**2)) V2 /= np.sqrt(np.sum(V2**2)) V3 /= np.sqrt(np.sum(V3**2)) M = np.array([V1,V2,V3]) D = nl.det(M) Ang = 1.0 P12 = np.dot(V1,V2) P13 = np.dot(V1,V3) P23 = np.dot(V2,V3) Tors = acosd((P12*P23-P13)/(np.sqrt(1.-P12**2)*np.sqrt(1.-P23**2)))*D/abs(D) sig = 0.0 def GetTorsionSig(Oatoms,Atoms,Amat,SGData,covData={}): def calcTorsion(Atoms,SyOps,Amat): XYZ = [] for i,atom in enumerate(Atoms): Inv,M,T,C,U = SyOps[i] XYZ.append(np.array(atom[1:4])) XYZ[-1] = Inv*(np.inner(M,np.array(XYZ[-1]))+T)+C+U XYZ[-1] = np.inner(Amat,XYZ[-1]).T V1 = XYZ[1]-XYZ[0] V2 = XYZ[2]-XYZ[1] V3 = XYZ[3]-XYZ[2] V1 /= np.sqrt(np.sum(V1**2)) V2 /= np.sqrt(np.sum(V2**2)) V3 /= np.sqrt(np.sum(V3**2)) M = np.array([V1,V2,V3]) D = nl.det(M) Ang = 1.0 P12 = np.dot(V1,V2) P13 = np.dot(V1,V3) P23 = np.dot(V2,V3) Tors = acosd((P12*P23-P13)/(np.sqrt(1.-P12**2)*np.sqrt(1.-P23**2)))*D/abs(D) return Tors Inv = [] SyOps = [] names = [] for i,atom in enumerate(Oatoms): names += atom[-1] Op,unit = Atoms[i][-1] inv = Op/abs(Op) m,t = SGData['SGOps'][abs(Op)%100-1] c = SGData['SGCen'][abs(Op)/100] SyOps.append([inv,m,t,c,unit]) Tors = calcTorsion(Oatoms,SyOps,Amat) sig = -0.01 if 'covMatrix' in covData: for atom in Atoms: xyz = atom[3:6] Op,Unit = atom[-1] inv = Op/abs(Op) M,T = SGData['SGOps'][abs(Op)%100-1] C = SGData['SGCen'][abs(Op)/100] #reverse inv*(np.inner(M,Tx)+T)+C XYZ = np.inner(nl.inv(M),((xyz-C)*inv-T)) print Op,Unit,xyz,XYZ parmNames = [] dx = .00001 dadx = np.zeros(12) ang = calcTorsion(Oatoms,SyOps,Amat) for i in range(12): ia = i/3 ix = i%3 Oatoms[ia][ix+1] += dx a0 = calcTorsion(Oatoms,SyOps,Amat) Oatoms[ia][ix+1] -= 2*dx dadx[i] = (calcTorsion(Oatoms,SyOps,Amat)-a0)/(2.*dx) covMatrix = covData['covMatrix'] varyList = covData['varyList'] TorVcov = getVCov(names,varyList,covMatrix) sig = np.sqrt(np.inner(dadx,np.inner(TorVcov,dadx))) if sig < 0.01: sig = -0.01 return Tors,sig return str(fmt%(value,int(round(esd*10**(mdec(esd)))))).strip('"') elif esd < 0: return str(round(value,mdec(esd))) return str(round(value,mdec(esd)-1)) else: text = str("%f"%(value))
• ## trunk/GSASIIphsGUI.py

 r476 if atom[cuia] == 'A': Uij = atom[cuij:cuij+6] result = G2spc.GenAtom(XYZ,SGData,False,Uij) result = G2spc.GenAtom(XYZ,SGData,False,Uij,False) for item in result: atom = copy.copy(atomData[ind]) atomData.append(atom[:]) else: result = G2spc.GenAtom(XYZ,SGData,False) result = G2spc.GenAtom(XYZ,SGData,False,Move=False) for item in result: atom = copy.copy(atomData[ind]) event.StopPropagation() def FindAtomIndexByIDs(atomData,IDs): def FindAtomIndexByIDs(atomData,IDs,Draw=True): indx = [] for i,atom in enumerate(atomData): if atom[-2] in IDs: if Draw and atom[-3] in IDs: indx.append(i) elif atom[-1] in IDs: indx.append(i) return indx return TorsionData = {} ocx,oct,ocs,cia = data['General']['AtomPtrs'] drawingData = data['Drawing'] atomData = data['Atoms'] #        for atom in atomData: print atom atomDData = drawingData['Atoms'] #        for atom in atomDData: print atom colLabels = [drawAtoms.GetColLabelValue(c) for c in range(drawAtoms.GetNumberCols())] cx = colLabels.index('x') cn = colLabels.index('Name') cid = colLabels.index('I/A')+8 xyz = [] Oxyz = [] for i,atom in enumerate(atomDData): if i in indx: xyz.append([i,]+atom[cn:cn+2]+atom[cx:cx+4]) #also gets Sym Op id = FindAtomIndexByIDs(atomData,[atom[cid],],False)[0] Oxyz.append([id,]+atomData[id][cx+1:cx+4]) TorsionData['Datoms'] = xyz TorsionData['Oatoms'] = Oxyz generalData = data['General'] TorsionData['Name'] = generalData['Name']
• ## trunk/GSASIIplot.py

 r468 return G2frame.dataDisplay.GetSelection() except AttributeError: print G2frame.dataDisplay.GetLabel() G2frame.G2plotNB.status.SetStatusText('Select this from Phase data window!') return 0 def Draw(): import numpy.linalg as nl Ind = GetSelectedAtoms() VS = np.array(Page.canvas.GetSize())
• ## trunk/GSASIIspc.py

 r459 return new def GenAtom(XYZ,SGData,All=False,Uij=[]): def GenAtom(XYZ,SGData,All=False,Uij=[],Move=True): ''' Generates the equivalent positions for a specified coordinate and space group Cell = [] X = np.array(XYZ) X = MoveToUnitCell(X) if Move: X = MoveToUnitCell(X) for ic,cen in enumerate(SGData['SGCen']): C = np.array(cen)
• ## trunk/GSASIIstruct.py

 r474 TxyzNames = [pfx+'dAx:%d'%(Tatom[0]),pfx+'dAy:%d'%(Tatom[0]),pfx+'dAz:%d'%(Tatom[0])] Xvcov = G2mth.getVCov(OxyzNames+TxyzNames,varyList,covMatrix) result = G2spc.GenAtom(Tatom[3:6],SGData,False) result = G2spc.GenAtom(Tatom[3:6],SGData,False,Move=False) BsumR = (Radii[Oatom[2]][0]+Radii[Tatom[2]][0])*Factor[0] AsumR = (Radii[Oatom[2]][1]+Radii[Tatom[2]][1])*Factor[1] print 80*'*' print '   Torsion angle for phase '+name print 80*'*','\n' print 80*'*' ShowBanner(TorsionData['Name']) SGData = TorsionData['SGData'] SGtext = G2spc.SGPrint(SGData) for line in SGtext: print line Cell = TorsionData['Cell'] Amat,Bmat = G2lat.cell2AB(Cell[:6]) covData = {} pfx = '' if 'covData' in TorsionData: covData = TorsionData['covData'] varyList = covData['varyList'] pfx = str(TorsionData['pId'])+'::' A = G2lat.cell2A(Cell[:6]) cellSig = getCellEsd(pfx,SGData,A,covData) names = [' a = ',' b = ',' c = ',' alpha = ',' beta = ',' gamma = ',' Volume = '] valEsd = [G2mth.ValEsd(Cell[i],cellSig[i],True) for i in range(7)] line = '\n Unit cell:' for name,vals in zip(names,valEsd): line += name+vals print line else: print '\n Unit cell: a = ','%.5f'%(Cell[0]),' b = ','%.5f'%(Cell[1]),' c = ','%.5f'%(Cell[2]), \ ' alpha = ','%.3f'%(Cell[3]),' beta = ','%.3f'%(Cell[4]),' gamma = ', \ '%.3f'%(Cell[5]),' volume = ','%.3f'%(Cell[6]) #find one end of 4 atom string - involved in longest distance dist = {} dist = {} X1 = TorsionData['Datoms'][end] for i,X2 in enumerate(TorsionData['Datoms']): for i,X2 in enumerate(TorsionData['Datoms']): dist[np.sqrt(np.sum(np.inner(Amat,np.array(X2[3:6])-np.array(X1[3:6]))**2))] = i sortdist = dist.keys() sortdist.sort() Datoms = [] Oatoms = [] for d in sortdist: atom = TorsionData['Datoms'][dist[d]] symop[1] = eval(symop[1]) atom[-1] = symop print atom Datoms.append(atom) Tors,sig = G2mth.GetTorsionSig(Datoms,Amat,SGData,covData={}) print ' Torsion: ',G2mth.ValEsd(Tors,sig) oatom = TorsionData['Oatoms'][dist[d]] names = ['','',''] if pfx: names = [pfx+'dAx:'+str(oatom[0]),pfx+'dAy:'+str(oatom[0]),pfx+'dAz:'+str(oatom[0])] oatom += [names,] Oatoms.append(oatom) Tors,sig = G2mth.GetTorsionSig(Oatoms,Datoms,Amat,SGData,covData) print ' Torsion angle for atom sequence: ',[Datoms[i][1] for i in range(4)],'=',G2mth.ValEsd(Tors,sig) print ' NB: Atom sequence determined by interatomic distances' def BestPlane(PlaneData):
Note: See TracChangeset for help on using the changeset viewer.