Changeset 1935 for trunk/GSASIImath.py

Timestamp:

Jul 20, 2015 9:25:46 AM (8 years ago)

Author:

vondreele

Message:

apply Uiso (equiv) for H-atom generation/update
revise H-atom addition to allow consecutive addition
deleted H atoms treated correctly

File:

• trunk/GSASIImath.py (modified) (8 diffs)

trunk/GSASIImath.py

@@ -421 +421 @@
     Oatom = GetAtomsById(Atoms,AtLookUp,[AddHydId[0],])[0]
     OXYZ = np.array(Oatom[cx:cx+3])
+    if 'I' in Oatom[cia]:
+        Uiso = Oatom[cia+1]
+    else:
+        Uiso = (Oatom[cia+2]+Oatom[cia+3]+Oatom[cia+4])/3.0       #simple average
+    Uiso = max(Uiso,0.005)                      #set floor!
     Tatoms = GetAtomsById(Atoms,AtLookUp,AddHydId[1])
     TXYZ = np.array([tatom[cx:cx+3] for tatom in Tatoms]) #3 x xyz
@@ -431 +436 @@
             Len = np.sqrt(np.sum(Vec**2))
             Hpos = OXYZ-0.98*np.inner(Bmat,Vec).T/Len
-            return [Hpos,]
+            HU = 1.1*Uiso
+            return [Hpos,],[HU,]
         elif AddHydId[-1] == 2:
             Vec = np.inner(Amat,TXYZ-OXYZ).T
@@ -442 +448 @@
             Hpos = np.array([[-0.97*cosd(54.75),0.97*sind(54.75),0.],
                 [-0.97*cosd(54.75),-0.97*sind(54.75),0.]])
+            HU = 1.2*Uiso*np.ones(2)
             Hpos = np.inner(Bmat,np.inner(iMat,Hpos).T).T+OXYZ
-            return Hpos
+            return Hpos,HU
         else:
             Ratom = GetAtomsById(Atoms,AtLookUp,[AddHydId[2],])[0]
@@ -452 +459 @@
             Hpos = np.array([[a,0.,-b],[a,-b*cosd(30.),0.5*b],[a,b*cosd(30.),0.5*b]])
             Hpos = np.inner(Bmat,np.inner(iMat,Hpos).T).T+OXYZ
-            return Hpos            
+            HU = 1.5*Uiso*np.ones(3)
+            return Hpos,HU          
     elif nBonds == 3:
         if AddHydId[-1] == 1:
@@ -459 +467 @@
             Vec = -0.93*Vec/Len
             Hpos = OXYZ+Vec
-            return [Hpos,]
+            HU = 1.1*Uiso
+            return [Hpos,],[HU,]
         elif AddHydId[-1] == 2:
             Ratom = GetAtomsById(Atoms,AtLookUp,[AddHydId[2],])[0]
@@ -468 +477 @@
             Hpos = [[a,b,0],[a,-b,0]]
             Hpos = np.inner(Bmat,np.inner(iMat,Hpos).T).T+OXYZ
-            return Hpos
+            HU = 1.2*Uiso*np.ones(2)
+            return Hpos,HU
     else:   #2 bonds
         if 'C' in Oatom[ct]:
@@ -475 +485 @@
             Vec = -0.93*Vec/Len
             Hpos = OXYZ+Vec
-            return [Hpos,]
+            HU = 1.1*Uiso
+            return [Hpos,],[HU,]
         elif 'O' in Oatom[ct]:
             mapData = General['Map']
@@ -488 +499 @@
             Rhos = np.array([getRho(pos,mapData) for pos in Hpos])
             imax = np.argmax(Rhos)
-            return [Hpos[imax],]
-    return []
+            HU = 1.5*Uiso
+            return [Hpos[imax],],[HU,]
+    return [],[]
         
 def AtomUij2TLS(atomData,atPtrs,Amat,Bmat,rbObj):   #unfinished & not used