Changeset 1926 for trunk/GSASIImath.py
 Timestamp:
 Jul 13, 2015 4:10:42 PM (8 years ago)
 File:

 1 edited
Legend:
 Unmodified
 Added
 Removed

trunk/GSASIImath.py
r1925 r1926 385 385 atTypes = General['AtomTypes'] 386 386 Radii = np.array(General['BondRadii']) 387 DisAglCtls = General['DisAglCtls'] 388 radiusFactor = DisAglCtls['Factors'][0] 387 389 AtInfo = dict(zip(atTypes,Radii)) #or General['BondRadii'] 388 390 Orig = atNames.index(FrstName) 391 OId = Atoms[Orig][cia+8] 389 392 OType = Atoms[Orig][ct] 390 393 XYZ = getAtomXYZ(Atoms,cx) 391 394 Neigh = [] 395 Ids = [] 392 396 Dx = np.inner(Amat,XYZXYZ[Orig]).T 393 397 dist = np.sqrt(np.sum(Dx**2,axis=1)) 394 398 sumR = np.array([AtInfo[OType]+AtInfo[atom[ct]] for atom in Atoms]) 395 IndB = ma.nonzero(ma.masked_greater(dist 0.85*sumR,0.))399 IndB = ma.nonzero(ma.masked_greater(distradiusFactor*sumR,0.)) 396 400 for j in IndB[0]: 397 401 if j != Orig: 398 402 Neigh.append([AtNames[j],dist[j]]) 399 return Neigh 403 Ids.append(Atoms[j][cia+8]) 404 return Neigh,[OId,Ids] 405 406 def AddHydrogens(AtLookUp,General,Atoms,AddHydId): 407 cx,ct,cs,cia = General['AtomPtrs'] 408 Cell = General['Cell'][1:7] 409 Amat,Bmat = G2lat.cell2AB(Cell) 410 nBonds = AddHydId[2]+len(AddHydId[1]) 411 Oatom = GetAtomsById(Atoms,AtLookUp,[AddHydId[0],])[0] 412 OXYZ = np.array(Oatom[cx:cx+3]) 413 Tatoms = GetAtomsById(Atoms,AtLookUp,AddHydId[1]) 414 TXYZ = np.array([tatom[cx:cx+3] for tatom in Tatoms]) #3 x xyz 415 DX = np.inner(Amat,TXYZOXYZ).T 416 if nBonds == 4: 417 if AddHydId[2] == 1: 418 Vec = TXYZOXYZ 419 Len = np.sqrt(np.sum(np.inner(Amat,Vec).T**2,axis=0)) 420 Vec = np.sum(Vec/Len,axis=0) 421 Len = np.sqrt(np.sum(Vec**2)) 422 Hpos = OXYZ0.98*np.inner(Bmat,Vec).T/Len 423 return [Hpos,] 424 elif AddHydId[2] == 2: 425 print 'add 2 H' 426 return [] 427 else: 428 print 'add 3 H' 429 return [] 430 elif nBonds == 3: 431 if AddHydId[2] == 1: 432 Vec = np.sum(TXYZOXYZ,axis=0) 433 Len = np.sqrt(np.sum(np.inner(Amat,Vec).T**2)) 434 Vec = 0.93*Vec/Len 435 Hpos = OXYZ+Vec 436 return [Hpos,] 437 elif AddHydId[2] == 2: 438 print 'add 2 H' 439 return [] 440 else: #2 bonds 441 if 'C' in Oatom[ct]: 442 print 'sp atom',Oatom[ct1] 443 return [] 444 elif 'O' in Oatom[ct]: 445 #idea  construct ring at 0.82 from O atom & find high spot? 446 print 'sp3 atom ',Oatom[ct1] 447 print 'add 1 H' 448 return [] 449 return [] 450 451 452 453 400 454 401 455 def AtomUij2TLS(atomData,atPtrs,Amat,Bmat,rbObj): #unfinished & not used
Note: See TracChangeset
for help on using the changeset viewer.