Changeset 1597

Timestamp:

Dec 6, 2014 3:12:32 PM (8 years ago)

Author:

vondreele

Message:

Incommensurate Pawley refinement works, including refinement of modulation vector.

Location:

trunk

Files:

• GSASIIlattice.py (modified) (1 diff)
• GSASIIstrIO.py (modified) (5 diffs)
• GSASIIstrMath.py (modified) (10 diffs)

trunk/GSASIIlattice.py

@@ -467 +467 @@
     return np.inner(H,np.inner(G,H))
     
+def calc_rDsqSS(H,A,vec):
+    'needs doc string'
+    rdsq = calc_rDsq(H[:3]+(H[3]*vec).T,A)
+    return rdsq
+       
 def calc_rDsqZ(H,A,Z,tth,lam):
     'needs doc string'

trunk/GSASIIstrIO.py

@@ -562 +562 @@
         print >>pFile,' Process histograms in reverse order: ',Controls['Reverse Seq']
     
-def GetPawleyConstr(SGLaue,PawleyRef,pawleyVary):
+def GetPawleyConstr(SGLaue,PawleyRef,im,pawleyVary):
     'needs a doc string'
 #    if SGLaue in ['-1','2/m','mmm']:
@@ -571 +571 @@
         refI = int(varyI.split(':')[-1])
         ih,ik,il = PawleyRef[refI][:3]
-        dspI = PawleyRef[refI][4]
+        dspI = PawleyRef[refI][4+im]
         for varyJ in pawleyVary[i+1:]:
             refJ = int(varyJ.split(':')[-1])
             jh,jk,jl = PawleyRef[refJ][:3]
-            dspJ = PawleyRef[refJ][4]
+            dspJ = PawleyRef[refJ][4+im]
             if SGLaue in ['4/m','4/mmm']:
                 isum = ih**2+ik**2
@@ -631 +631 @@
 def modVary(pfx,SSGData):
     vary = []
-    for i,item in SSGData['modSymb']:
+    for i,item in enumerate(SSGData['modSymb']):
         if item in ['a','b','g']:
             vary.append(pfx+'mV%d'%(i))
@@ -1148 +1148 @@
                 if refl[5+im]:
                     pawleyVary.append(pfx+'PWLref:'+str(i))
-            GetPawleyConstr(SGData['SGLaue'],PawleyRef,pawleyVary)      #does G2mv.StoreEquivalence
+            GetPawleyConstr(SGData['SGLaue'],PawleyRef,im,pawleyVary)      #does G2mv.StoreEquivalence
             phaseVary += pawleyVary
                 
@@ -1607 +1607 @@
                 for var in ['mV0','mV1','mV2']:
                     namstr += '%12s'%(pfx+var)
-                    ptstr += '%12.4f'%(parmDict[pfx+var])
+                    ptstr += '%12.6f'%(parmDict[pfx+var])
                     if pfx+var in sigDict:
-                        sigstr += '%12.4f'%(sigDict[pfx+var])
+                        sigstr += '%12.6f'%(sigDict[pfx+var])
                     else:
                         sigstr += 12*' '

trunk/GSASIIstrMath.py

@@ -1370 +1370 @@
     return sigDict,gamDict
         
-def GetReflPos(refl,im,wave,A,hfx,calcControls,parmDict):
+def GetReflPos(refl,im,wave,A,pfx,hfx,calcControls,parmDict):
     'Needs a doc string'
-    h,k,l = refl[:3]
-    d = 1./np.sqrt(G2lat.calc_rDsq(np.array([h,k,l]),A))
-
+    if im:
+        h,k,l,m = refl[:4]
+        vec = np.array([parmDict[pfx+'mV0'],parmDict[pfx+'mV1'],parmDict[pfx+'mV2']])
+        d = 1./np.sqrt(G2lat.calc_rDsqSS(np.array([h,k,l,m]),A,vec))
+    else:
+        h,k,l = refl[:3]
+        d = 1./np.sqrt(G2lat.calc_rDsq(np.array([h,k,l]),A))
     refl[4+im] = d
     if 'C' in calcControls[hfx+'histType']:
@@ -1389 +1393 @@
     return pos
 
-def GetReflPosDerv(refl,im,wave,A,hfx,calcControls,parmDict):
+def GetReflPosDerv(refl,im,wave,A,pfx,hfx,calcControls,parmDict):
     'Needs a doc string'
     dpr = 180./np.pi
-    h,k,l = refl[:3]
-    dstsq = G2lat.calc_rDsq(np.array([h,k,l]),A)
+    if im:
+        h,k,l,m = refl[:4]
+        vec = np.array([parmDict[pfx+'mV0'],parmDict[pfx+'mV1'],parmDict[pfx+'mV2']])
+        dstsq = G2lat.calc_rDsqSS(np.array([h,k,l,m]),A,vec)
+    else:
+        m = 0
+        h,k,l = refl[:3]        
+        dstsq = G2lat.calc_rDsq(np.array([h,k,l]),A)
     dst = np.sqrt(dstsq)
     dsp = 1./dst
@@ -1402 +1412 @@
         dpdA = np.array([h**2,k**2,l**2,h*k,h*l,k*l])*const*wave/(2.0*dst)
         dpdZ = 1.0
-        const = 9.e-2/(np.pi*parmDict[hfx+'Gonio. radius'])                  #shifts in microns
+        dpdV = np.array([2.*h*A[0]+k*A[3]+l*A[4],2*k*A[1]+h*A[3]+l*A[5],
+            2*l*A[2]+h*A[4]+k*A[5]])*m*const*wave/(2.0*dst)
+        shft = 9.e-2/(np.pi*parmDict[hfx+'Gonio. radius'])                  #shifts in microns
         if 'Bragg' in calcControls[hfx+'instType']:
-            dpdSh = -4.*const*cosd(pos/2.0)
-            dpdTr = -const*sind(pos)*100.0
-            return dpdA,dpdw,dpdZ,dpdSh,dpdTr,0.,0.
+            dpdSh = -4.*shft*cosd(pos/2.0)
+            dpdTr = -shft*sind(pos)*100.0
+            return dpdA,dpdw,dpdZ,dpdSh,dpdTr,0.,0.,dpdV
         else:               #Debye-Scherrer - simple but maybe not right
-            dpdXd = -const*cosd(pos)
-            dpdYd = -const*sind(pos)
-            return dpdA,dpdw,dpdZ,0.,0.,dpdXd,dpdYd
+            dpdXd = -shft*cosd(pos)
+            dpdYd = -shft*sind(pos)
+            return dpdA,dpdw,dpdZ,0.,0.,dpdXd,dpdYd,dpdV
     elif 'T' in calcControls[hfx+'histType']:
         dpdA = -np.array([h**2,k**2,l**2,h*k,h*l,k*l])*parmDict[hfx+'difC']*dsp**3/2.
@@ -1417 +1429 @@
         dpdDA = dsp**2
         dpdDB = 1./dsp
-        return dpdA,dpdZ,dpdDC,dpdDA,dpdDB
+        dpdV = np.array([2.*h*A[0]+k*A[3]+l*A[4],2*k*A[1]+h*A[3]+l*A[5],
+            2*l*A[2]+h*A[4]+k*A[5]])*m**parmDict[hfx+'difC']*dsp**3/2.
+        return dpdA,dpdZ,dpdDC,dpdDA,dpdDB,dpdV
             
 def GetHStrainShift(refl,im,SGData,phfx,hfx,calcControls,parmDict):
@@ -1653 +1667 @@
                     h,k,l = refl[:3]
                 Uniq = np.inner(refl[:3],SGMT)
-                refl[5+im] = GetReflPos(refl,im,wave,A,hfx,calcControls,parmDict)         #corrected reflection position
+                refl[5+im] = GetReflPos(refl,im,wave,A,pfx,hfx,calcControls,parmDict)         #corrected reflection position
                 Lorenz = 1./(2.*sind(refl[5+im]/2.)**2*cosd(refl[5+im]/2.))           #Lorentz correction
 #                refl[5+im] += GetHStrainShift(refl,im,SGData,phfx,hfx,calcControls,parmDict)               #apply hydrostatic strain shift
@@ -1696 +1710 @@
                 h,k,l = refl[:3]
                 Uniq = np.inner(refl[:3],SGMT)
-                refl[5+im] = GetReflPos(refl,im,0.0,A,hfx,calcControls,parmDict)         #corrected reflection position
+                refl[5+im] = GetReflPos(refl,im,0.0,A,pfx,hfx,calcControls,parmDict)         #corrected reflection position
                 Lorenz = sind(parmDict[hfx+'2-theta']/2)*refl[4+im]**4                                                #TOF Lorentz correction
 #                refl[5+im] += GetHStrainShift(refl,im,SGData,phfx,hfx,calcControls,parmDict)               #apply hydrostatic strain shift
@@ -1821 +1835 @@
         time0 = time.time()
         for iref,refl in enumerate(refDict['RefList']):
-            h,k,l = refl[:3]
+            if im:
+                h,k,l,m = refl[:4]
+            else:
+                h,k,l = refl[:3]
             Uniq = np.inner(refl[:3],SGMT)
             if 'T' in calcControls[hfx+'histType']:
@@ -1882 +1899 @@
                     pass
             if 'C' in calcControls[hfx+'histType']:
-                dpdA,dpdw,dpdZ,dpdSh,dpdTr,dpdX,dpdY = GetReflPosDerv(refl,im,wave,A,hfx,calcControls,parmDict)
+                dpdA,dpdw,dpdZ,dpdSh,dpdTr,dpdX,dpdY,dpdV = GetReflPosDerv(refl,im,wave,A,pfx,hfx,calcControls,parmDict)
                 names = {hfx+'Scale':[dIdsh,'int'],hfx+'Polariz.':[dIdpola,'int'],phfx+'Scale':[dIdsp,'int'],
                     hfx+'U':[tanth**2,'sig'],hfx+'V':[tanth,'sig'],hfx+'W':[1.0,'sig'],
@@ -1895 +1912 @@
                     names.update({hfx+'Absorption':[dFdAb,'int'],})
             else:   #'T'OF
-                dpdA,dpdZ,dpdDC,dpdDA,dpdDB = GetReflPosDerv(refl,im,0.0,A,hfx,calcControls,parmDict)
+                dpdA,dpdZ,dpdDC,dpdDA,dpdDB,dpdV = GetReflPosDerv(refl,im,0.0,A,pfx,hfx,calcControls,parmDict)
                 names = {hfx+'Scale':[dIdsh,'int'],phfx+'Scale':[dIdsp,'int'],
                     hfx+'difC':[dpdDC,'pos'],hfx+'difA':[dpdDA,'pos'],hfx+'difB':[dpdDB,'pos'],
@@ -1961 +1978 @@
                     if Ka2:
                         depDerivDict[name][iBeg2:iFin2] += dDijDict[name]*dervDict2['pos']
+            for i,name in enumerate([pfx+'mV0',pfx+'mV1',pfx+'mV2']):
+                if name in varylist:
+                    dMdv[varylist.index(name)][iBeg:iFin] += dpdV[i]*dervDict['pos']
+                    if Ka2:
+                        dMdv[varylist.index(name)][iBeg2:iFin2] += dpdV[i]*dervDict2['pos']
+                elif name in dependentVars:
+                    depDerivDict[name][iBeg:iFin] += dpdV[i]*dervDict['pos']
+                    if Ka2:
+                        depDerivDict[name][iBeg2:iFin2] += dpdV[i]*dervDict2['pos']
             if 'C' in calcControls[hfx+'histType']:
                 sigDict,gamDict = GetSampleSigGamDerv(refl,im,wave,G,GB,SGData,hfx,phfx,calcControls,parmDict)