Changeset 1682

Timestamp:

Feb 28, 2015 9:16:41 AM (7 years ago)

Author:

vondreele

Message:

do sums of Bragg intensity, and background terms put results into .lst file

Location:

trunk

Files:

• GSASIIpwd.py (modified) (10 diffs)
• GSASIIstrIO.py (modified) (2 diffs)
• GSASIIstrMath.py (modified) (6 diffs)

trunk/GSASIIpwd.py

@@ -525 +525 @@
     cw = np.diff(xdata)
     cw = np.append(cw,cw[-1])
+    sumBk = [0.,0.,0]
     while True:
         key = pfx+'Back;'+str(nBak)
@@ -531 +532 @@
         else:
             break
+#empirical functions
     if bakType in ['chebyschev','cosine']:
         dt = xdata[-1]-xdata[0]    
@@ -536 +538 @@
             key = pfx+'Back;'+str(iBak)
             if bakType == 'chebyschev':
-                yb += parmDict[key]*(2.*(xdata-xdata[0])/dt-1.)**iBak
+                ybi = parmDict[key]*(2.*(xdata-xdata[0])/dt-1.)**iBak
             elif bakType == 'cosine':
-                yb += parmDict[key]*npcosd(xdata*iBak)
+                ybi = parmDict[key]*npcosd(xdata*iBak)
+            yb += ybi
+            sumBk[0] = np.sum(ybi)
     elif bakType in ['lin interpolate','inv interpolate','log interpolate',]:
         if nBak == 1:
@@ -561 +565 @@
             bakInt = si.interp1d(bakPos,bakVals,'linear')
             yb = bakInt(xdata)
+        sumBk[0] = np.sum(yb)
+#Debye function        
     if pfx+'difC' in parmDict:
         ff = 1.
@@ -578 +584 @@
             dbR = parmDict[pfx+'DebyeR;'+str(iD)]
             dbU = parmDict[pfx+'DebyeU;'+str(iD)]
-            yb += ff*dbA*np.sin(q*dbR)*np.exp(-dbU*q**2)/(q*dbR)
+            ybi = ff*dbA*np.sin(q*dbR)*np.exp(-dbU*q**2)/(q*dbR)
+            yb += ybi
+            sumBk[1] += np.sum(ybi)
             iD += 1       
         except KeyError:
             break
+#peaks
     iD = 0
     while True:
@@ -603 +612 @@
                 iFin = np.searchsorted(xdata,pkP+fmax)
             if 'C' in dataType:
-                yb[iBeg:iFin] += pkI*getFCJVoigt3(pkP,pkS,pkG,0.002,xdata[iBeg:iFin])
+                ybi = pkI*getFCJVoigt3(pkP,pkS,pkG,0.002,xdata[iBeg:iFin])
+                yb[iBeg:iFin] += ybi
             else:   #'T'OF
-                yb[iBeg:iFin] += pkI*getEpsVoigt(pkP,1.,1.,pkS,pkG,xdata[iBeg:iFin])
+                ybi = pkI*getEpsVoigt(pkP,1.,1.,pkS,pkG,xdata[iBeg:iFin])
+                yb[iBeg:iFin] += ybi
+            sumBk[2] += np.sum(ybi)
             iD += 1       
         except KeyError:
@@ -612 +624 @@
             print '**** WARNING - backround peak '+str(iD)+' sigma is negative; fix & try again ****'
             break        
-    return yb
+    return yb,sumBk
     
 def getBackgroundDerv(hfx,parmDict,bakType,dataType,xdata):
@@ -836 +848 @@
     'needs a doc string'
     
-    yb = getBackground('',parmDict,bakType,dataType,xdata)
+    yb = getBackground('',parmDict,bakType,dataType,xdata)[0]
     yc = np.zeros_like(yb)
     cw = np.diff(xdata)
@@ -1574 +1586 @@
         
     sigDict = dict(zip(varyList,sig))
-    yb[xBeg:xFin] = getBackground('',parmDict,bakType,dataType,x[xBeg:xFin])
+    yb[xBeg:xFin] = getBackground('',parmDict,bakType,dataType,x[xBeg:xFin])[0]
     yc[xBeg:xFin] = getPeakProfile(dataType,parmDict,x[xBeg:xFin],varyList,bakType)
     yd[xBeg:xFin] = y[xBeg:xFin]-yc[xBeg:xFin]
@@ -1685 +1697 @@
     msg = 'test '
     gplot = plotter.add('FCJ-Voigt, 11BM').gca()
-    gplot.plot(xdata,getBackground('',parmDict0,bakType,'PXC',xdata))   
+    gplot.plot(xdata,getBackground('',parmDict0,bakType,'PXC',xdata)[0])   
     gplot.plot(xdata,getPeakProfile(parmDict0,xdata,varyList,bakType))
     fplot = plotter.add('FCJ-Voigt, Ka1+2').gca()
-    fplot.plot(xdata,getBackground('',parmDict1,bakType,'PXC',xdata))   
+    fplot.plot(xdata,getBackground('',parmDict1,bakType,'PXC',xdata)[0])   
     fplot.plot(xdata,getPeakProfile(parmDict1,xdata,varyList,bakType))
     

trunk/GSASIIstrIO.py

@@ -2475 +2475 @@
                     print >>pFile,' Final refinement RF, RF^2 = %.2f%%, %.2f%% on %d reflections'   \
                         %(Histogram['Residuals'][pfx+'Rf'],Histogram['Residuals'][pfx+'Rf^2'],Histogram['Residuals'][pfx+'Nref'])
+                    print >>pFile,' Bragg intensity sum = %.3g'%(Histogram['Residuals'][pfx+'sumInt'])
                 
                     if pfx+'Scale' in PhFrExtPOSig:
@@ -2851 +2852 @@
                 print >>pFile,ptstr
                 print >>pFile,sigstr
+        sumBk = np.array(Histogram['sumBk'])
+        print >>pFile,' Background sums: empirical %.3g, Debye %.3g, peaks %.3g, Total %.3g'    \
+            %(sumBk[0],sumBk[1],sumBk[2],np.sum(sumBk))
         
     def PrintInstParmsSig(Inst,instSig):

trunk/GSASIIstrMath.py

@@ -1772 +1772 @@
                 sumFosq = 0.0
                 sumdFsq = 0.0
+                sumInt = 0.0
                 for refl in refDict['RefList']:
                     if 'C' in calcControls[hfx+'histType']:
@@ -1785 +1786 @@
                         elif iBeg < iFin:
                             yp[iBeg:iFin] = refl[11+im]*refl[9+im]*G2pwd.getFCJVoigt3(refl[5+im],refl[6+im],refl[7+im],shl,ma.getdata(x[iBeg:iFin]))    #>90% of time spent here
+                            sumInt += refl[11+im]*refl[9+im]
                             if Ka2:
                                 pos2 = refl[5+im]+lamRatio*tand(refl[5+im]/2.0)       # + 360/pi * Dlam/lam * tan(th)
@@ -1791 +1793 @@
                                 iFin2 = min(np.searchsorted(x,pos2+fmax),xF)
                                 yp[iBeg2:iFin2] += refl[11+im]*refl[9+im]*kRatio*G2pwd.getFCJVoigt3(pos2,refl[6+im],refl[7+im],shl,ma.getdata(x[iBeg2:iFin2]))        #and here
+                                sumInt += refl[11+im]*refl[9+im]*kRatio
                             refl[8+im] = np.sum(np.where(ratio[iBeg:iFin2]>0.,yp[iBeg:iFin2]*ratio[iBeg:iFin2]/(refl[11+im]*(1.+kRatio)),0.0))
+                                
                     elif 'T' in calcControls[hfx+'histType']:
                         yp = np.zeros_like(yb)
@@ -1800 +1804 @@
                             yp[iBeg:iFin] = refl[11+im]*refl[9+im]*G2pwd.getEpsVoigt(refl[5+im],refl[12+im],refl[13+im],refl[6+im],refl[7+im],ma.getdata(x[iBeg:iFin]))  #>90% of time spent here
                             refl[8+im] = np.sum(np.where(ratio[iBeg:iFin]>0.,yp[iBeg:iFin]*ratio[iBeg:iFin]/refl[11+im],0.0))
+                            sumInt += refl[11+im]*refl[9+im]
                     Fo = np.sqrt(np.abs(refl[8+im]))
                     Fc = np.sqrt(np.abs(refl[9]+im))
@@ -1812 +1817 @@
                     Histogram['Residuals'][phfx+'Rf'] = 100.
                     Histogram['Residuals'][phfx+'Rf^2'] = 100.
+                Histogram['Residuals'][phfx+'sumInt'] = sumInt
                 Histogram['Residuals'][phfx+'Nref'] = len(refDict['RefList'])
                 Histogram['Residuals']['hId'] = hId
@@ -1852 +1858 @@
     hfx = ':%d:'%(hId)
     bakType = calcControls[hfx+'bakType']
-    yb = G2pwd.getBackground(hfx,parmDict,bakType,calcControls[hfx+'histType'],x)
+    yb,Histogram['sumBk'] = G2pwd.getBackground(hfx,parmDict,bakType,calcControls[hfx+'histType'],x)
     yc = np.zeros_like(yb)
     cw = np.diff(x)