Changeset 679


Ignore:
Timestamp:
Jul 5, 2012 11:46:52 AM (9 years ago)
Author:
vondreele
Message:

add Debye-Scherrer absorption to refinement

Location:
trunk
Files:
2 edited

Legend:

Unmodified
Added
Removed
  • trunk/GSASIIpwd.py

    r668 r679  
    8181        return 0.0
    8282
    83 def Absorb(Geometry,Abs,Diam,Tth,Phi=0,Psi=0):
     83def Absorb(Geometry,MuR,Tth,Phi=0,Psi=0):
    8484#Calculate sample absorption
    8585#   Geometry: one of 'Cylinder','Bragg-Brentano','Tilting Flat Plate in transmission','Fixed flat plate'
    86 #   Abs: absorption coeff in cm-1
    87 #   Diam: sample thickness/diameter in mm
     86#   MuR: absorption coeff * sample thickness/2 or radius
    8887#   Tth: 2-theta scattering angle - can be numpy array
    8988#   Phi: flat plate tilt angle - future
     
    9291    Cth2 = 1.-Sth2
    9392    if 'Cylinder' in Geometry:      #Lobanov & Alte da Veiga for 2-theta = 0; beam fully illuminates sample
    94         MuR = Abs*Diam/20.0
    9593        if MuR < 3.0:
    9694            T0 = 16.0/(3*np.pi)
     
    116114    elif 'Fixed' in Geometry: #assumes sample plane is perpendicular to incident beam
    117115        # and only defined for 2theta < 90
    118         MuR = Abs*Diam/10.0
    119         T1 = np.exp(-MuR)
    120         T2 = np.exp(-MuR/npcosd(Tth))
    121         Tb = MuR-MuR/npcosd(Tth)
     116        MuT = 2.*MuR
     117        T1 = np.exp(-MuT)
     118        T2 = np.exp(-MuT/npcosd(Tth))
     119        Tb = MuT-MuT/npcosd(Tth)
    122120        return (T2-T1)/Tb
    123121    elif 'Tilting' in Geometry: #assumes symmetric tilt so sample plane is parallel to diffraction vector
    124         MuR = Abs*Diam/10.0
     122        MuT = 2.*MuR
    125123        cth = npcosd(Tth/2.0)
    126         return np.exp(-MuR/cth)/cth
     124        return np.exp(-MuT/cth)/cth
     125       
     126def AbsorbDerv(Geometry,MuR,Tth,Phi=0,Psi=0):
     127    dA = 0.001
     128    AbsP = Absorb(Geometry,MuR+dA,Tth,Phi,Psi)
     129    if MuR:
     130        AbsM = Absorb(Geometry,MuR-dA,Tth,Phi,Psi)
     131        return (AbsP-AbsM)/(2.0*dA)
     132    else:
     133        return (AbsP-1.)/dA
    127134       
    128135def Polarization(Pola,Tth,Azm=0.0):
     
    392399    Tth = xydata['Sample'][1][0]
    393400    dt = (Tth[1]-Tth[0])
    394     xydata['IofQ'][1][1] /= Absorb(data['Geometry'],Abs,data['Diam'],Tth)
     401    MuR = Abs*data['Diam']/20.0
     402    xydata['IofQ'][1][1] /= Absorb(data['Geometry'],MuR,Tth)
    395403    xydata['IofQ'][1][1] /= Polarization(inst['Polariz.'],Tth,Azm=inst['Azimuth'])[0]
    396404    if data['DetType'] == 'Image plate':
  • trunk/GSASIIstruct.py

    r660 r679  
    21232123   
    21242124def SHTXcal(refl,g,pfx,hfx,SGData,calcControls,parmDict):
     2125    #Spherical harmonics texture
    21252126    IFCoup = 'Bragg' in calcControls[hfx+'instType']
    21262127    odfCor = 1.0
     
    21412142   
    21422143def SHTXcalDerv(refl,g,pfx,hfx,SGData,calcControls,parmDict):
     2144    #Spherical harmonics texture derivatives
    21432145    FORPI = 12.5663706143592
    21442146    IFCoup = 'Bragg' in calcControls[hfx+'instType']
     
    21652167   
    21662168def SHPOcal(refl,g,phfx,hfx,SGData,calcControls,parmDict):
     2169    #sphericaql harmonics preferred orientation (cylindrical symmetry only)
    21672170    odfCor = 1.0
    21682171    H = refl[:3]
     
    21852188   
    21862189def SHPOcalDerv(refl,g,phfx,hfx,SGData,calcControls,parmDict):
     2190    #sphericaql harmonics preferred orientation derivatives (cylindrical symmetry only)
    21872191    FORPI = 12.5663706143592
    21882192    odfCor = 1.0
     
    22442248    return POcorr,POderv
    22452249   
     2250def GetAbsorb(refl,hfx,calcControls,parmDict):
     2251    if 'Debye' in calcControls[hfx+'instType']:
     2252        return G2pwd.Absorb('Cylinder',parmDict[hfx+'Absorption'],refl[5],0,0)
     2253    else:
     2254        return 1.0
     2255   
     2256def GetAbsorbDerv(refl,hfx,calcControls,parmDict):
     2257    if 'Debye' in calcControls[hfx+'instType']:
     2258        return G2pwd.AbsorbDerv('Cylinder',parmDict[hfx+'Absorption'],refl[5],0,0)
     2259    else:
     2260        return 0.0
     2261   
    22462262def GetIntensityCorr(refl,G,g,pfx,phfx,hfx,SGData,calcControls,parmDict):
    22472263    Icorr = parmDict[phfx+'Scale']*parmDict[hfx+'Scale']*refl[3]               #scale*multiplicity
     
    22512267    if pfx+'SHorder' in parmDict:
    22522268        Icorr *= SHTXcal(refl,g,pfx,hfx,SGData,calcControls,parmDict)
     2269    Icorr *= GetAbsorb(refl,hfx,calcControls,parmDict)
    22532270    refl[13] = Icorr       
    22542271   
     
    22722289        for i in range(3):
    22732290            dFdSA[i] /= odfCor
    2274     return dIdsh,dIdsp,dIdPola,dIdPO,dFdODF,dFdSA
     2291    dFdAb = GetAbsorbDerv(refl,hfx,calcControls,parmDict)
     2292    return dIdsh,dIdsp,dIdPola,dIdPO,dFdODF,dFdSA,dFdAb
    22752293       
    22762294def GetSampleSigGam(refl,wave,G,GB,phfx,calcControls,parmDict):
     
    27372755            if 'C' in calcControls[hfx+'histType']:        #CW powder
    27382756                h,k,l = refl[:3]
    2739                 dIdsh,dIdsp,dIdpola,dIdPO,dFdODF,dFdSA = GetIntensityDerv(refl,G,g,pfx,phfx,hfx,SGData,calcControls,parmDict)
     2757                dIdsh,dIdsp,dIdpola,dIdPO,dFdODF,dFdSA,dFdAb = GetIntensityDerv(refl,G,g,pfx,phfx,hfx,SGData,calcControls,parmDict)
    27402758                Wd,fmin,fmax = G2pwd.getWidths(refl[5],refl[6],refl[7],shl)
    27412759                iBeg = np.searchsorted(x,refl[5]-fmin)
     
    27902808                    hfx+'I(L2)/I(L1)':[1.0,'L1/L2'],hfx+'Zero':[dpdZ,'pos'],hfx+'Lam':[dpdw,'pos'],
    27912809                    hfx+'Shift':[dpdSh,'pos'],hfx+'Transparency':[dpdTr,'pos'],hfx+'DisplaceX':[dpdX,'pos'],
    2792                     hfx+'DisplaceY':[dpdY,'pos'],}
     2810                    hfx+'DisplaceY':[dpdY,'pos'],hfx+'Absorption':[dFdAb,'int'],}
    27932811                for name in names:
    27942812                    item = names[name]
Note: See TracChangeset for help on using the changeset viewer.