13. GSASIIpwd: Powder calculations module¶
-
GSASIIpwd.
Absorb
(Geometry, MuR, Tth, Phi=0, Psi=0)[source]¶ Calculate sample absorption :param str Geometry: one of ‘Cylinder’,’Bragg-Brentano’,’Tilting Flat Plate in transmission’,’Fixed flat plate’ :param float MuR: absorption coeff * sample thickness/2 or radius :param Tth: 2-theta scattering angle - can be numpy array :param float Phi: flat plate tilt angle - future :param float Psi: flat plate tilt axis - future
-
GSASIIpwd.
CalcPDF
(data, inst, limits, xydata)[source]¶ Computes I(Q), S(Q) & G(r) from Sample, Bkg, etc. diffraction patterns loaded into dict xydata; results are placed in xydata. Calculation parameters are found in dicts data and inst and list limits. The return value is at present an empty list.
-
GSASIIpwd.
Dict2Values
(parmdict, varylist)[source]¶ Use before call to leastsq to setup list of values for the parameters in parmdict, as selected by key in varylist
-
GSASIIpwd.
DoPeakFit
(FitPgm, Peaks, Background, Limits, Inst, Inst2, data, fixback=None, prevVaryList=[], oneCycle=False, controls=None, wtFactor=1.0, dlg=None, noFit=False)[source]¶ Called to perform a peak fit, refining the selected items in the peak table as well as selected items in the background.
Parameters: - FitPgm (str) – type of fit to perform. At present this is ignored.
- Peaks (list) – a list of peaks. Each peak entry is a list with 8 values: four values followed by a refine flag where the values are: position, intensity, sigma (Gaussian width) and gamma (Lorentzian width). From the Histogram/”Peak List” tree entry, dict item “peaks”
- Background (list) – describes the background. List with two items. Item 0 specifies a background model and coefficients. Item 1 is a dict. From the Histogram/Background tree entry.
- Limits (list) – min and max x-value to use
- Inst (dict) – Instrument parameters
- Inst2 (dict) – more Instrument parameters
- data (numpy.array) – a 5xn array. data[0] is the x-values, data[1] is the y-values, data[2] are weight values, data[3], [4] and [5] are calc, background and difference intensities, respectively.
- fixback (array) – fixed background array; same size as data[0-5]
- prevVaryList (list) – Used in sequential refinements to override the variable list. Defaults as an empty list.
- oneCycle (bool) – True if only one cycle of fitting should be performed
- controls (dict) – a dict specifying two values, Ftol = controls[‘min dM/M’] and derivType = controls[‘deriv type’]. If None default values are used.
- wtFactor (float) – weight multiplier; = 1.0 by default
- dlg (wx.Dialog) – A dialog box that is updated with progress from the fit. Defaults to None, which means no updates are done.
- noFit (bool) – When noFit is True, a refinement is not performed. Default is False.
-
GSASIIpwd.
GetAsfMean
(ElList, Sthl2)[source]¶ Calculate various scattering factor terms for PDF calcs
Parameters: - ElList (dict) – element dictionary contains scattering factor coefficients, etc.
- Sthl2 (np.array) – numpy array of sin theta/lambda squared values
Returns: mean(f^2), mean(f)^2, mean(compton)
-
GSASIIpwd.
GetPDFfitAtomVar
(Phase, RMCPdict)[source]¶ Find dict of independent “@n” variables for PDFfit in atom constraints
-
GSASIIpwd.
MEMupdateReflData
(prfName, data, reflData)[source]¶ Update reflection data with new Fosq, phase result from Dysnomia
Parameters: - prfName (str) – phase.mem file name
- reflData (list) – GSAS-II reflection data
-
GSASIIpwd.
MakefullrmcRun
(pName, Phase, RMCPdict)[source]¶ Creates a script to run fullrmc. Returns the name of the file that was created.
-
GSASIIpwd.
PhaseWtSum
(G2frame, histo)[source]¶ Calculate sum of phase mass*phase fraction for PWDR data (exclude magnetic phases)
Parameters: - G2frame – GSASII main frame structure
- histo (str) – histogram name
Returns: sum(scale*mass) for phases in histo
-
GSASIIpwd.
Polarization
(Pola, Tth, Azm=0.0)[source]¶ Calculate angle dependent x-ray polarization correction (not scaled correctly!)
Parameters: - Pola – polarization coefficient e.g 1.0 fully polarized, 0.5 unpolarized
- Azm – azimuthal angle e.g. 0.0 in plane of polarization - can be numpy array
- Tth – 2-theta scattering angle - can be numpy array which (if either) of these is “right”?
Returns: (pola, dpdPola) - both 2-d arrays * pola = ((1-Pola)*npcosd(Azm)**2+Pola*npsind(Azm)**2)*npcosd(Tth)**2+ (1-Pola)*npsind(Azm)**2+Pola*npcosd(Azm)**2 * dpdPola: derivative needed for least squares
-
GSASIIpwd.
SetBackgroundParms
(Background)[source]¶ Loads background parameters into dicts/lists to create varylist & parmdict
-
GSASIIpwd.
StackSim
(Layers, ctrls, scale=0.0, background={}, limits=[], inst={}, profile=[])[source]¶ Simulate powder or selected area diffraction pattern from stacking faults using DIFFaX
Parameters: - Layers (dict) –
dict with following items
{'Laue':'-1','Cell':[False,1.,1.,1.,90.,90.,90,1.], 'Width':[[10.,10.],[False,False]],'Toler':0.01,'AtInfo':{}, 'Layers':[],'Stacking':[],'Transitions':[]}
- ctrls (str) – controls string to be written on DIFFaX controls.dif file
- scale (float) – scale factor
- background (dict) – background parameters
- limits (list) – min/max 2-theta to be calculated
- inst (dict) – instrument parameters dictionary
- profile (list) – powder pattern data
Note that parameters all updated in place
- Layers (dict) –
-
GSASIIpwd.
SurfaceRough
(SRA, SRB, Tth)[source]¶ Suortti (J. Appl. Cryst, 5,325-331, 1972) surface roughness correction :param float SRA: Suortti surface roughness parameter :param float SRB: Suortti surface roughness parameter :param float Tth: 2-theta(deg) - can be numpy array
-
GSASIIpwd.
SurfaceRoughDerv
(SRA, SRB, Tth)[source]¶ Suortti surface roughness correction derivatives :param float SRA: Suortti surface roughness parameter (dimensionless) :param float SRB: Suortti surface roughness parameter (dimensionless) :param float Tth: 2-theta(deg) - can be numpy array :return list: [dydSRA,dydSRB] derivatives to be used for intensity derivative
-
GSASIIpwd.
Transmission
(Geometry, Abs, Diam)[source]¶ Calculate sample transmission
Parameters: - Geometry (str) – one of ‘Cylinder’,’Bragg-Brentano’,’Tilting flat plate in transmission’,’Fixed flat plate’
- Abs (float) – absorption coeff in cm-1
- Diam (float) – sample thickness/diameter in mm
-
GSASIIpwd.
Values2Dict
(parmdict, varylist, values)[source]¶ Use after call to leastsq to update the parameter dictionary with values corresponding to keys in varylist
-
GSASIIpwd.
abeles
(kz, depth, rho, irho=0, sigma=0)[source]¶ Optical matrix form of the reflectivity calculation. O.S. Heavens, Optical Properties of Thin Solid Films
Reflectometry as a function of kz for a set of slabs.
Parameters: - kz – float[n] (1/Ang). Scattering vector, \(2\pi\sin(\theta)/\lambda\). This is \(\tfrac12 Q_z\).
- depth – float[m] (Ang). thickness of each layer. The thickness of the incident medium and substrate are ignored.
- rho – float[n,k] (1e-6/Ang^2) Real scattering length density for each layer for each kz
- irho – float[n,k] (1e-6/Ang^2) Imaginary scattering length density for each layer for each kz Note: absorption cross section mu = 2 irho/lambda for neutrons
- sigma – float[m-1] (Ang) interfacial roughness. This is the roughness between a layer and the previous layer. The sigma array should have m-1 entries.
Slabs are ordered with the surface SLD at index 0 and substrate at index -1, or reversed if kz < 0.
-
GSASIIpwd.
ellipseSize
(H, Sij, GB)[source]¶ Implements r=1/sqrt(sum((1/S)*(q.v)^2) per note from Alexander Brady
-
GSASIIpwd.
ellipseSizeDerv
(H, Sij, GB)[source]¶ Implements r=1/sqrt(sum((1/S)*(q.v)^2) derivative per note from Alexander Brady
-
GSASIIpwd.
factorize
(num)[source]¶ Provide prime number factors for integer num :returns: dictionary of prime factors (keys) & power for each (data)
-
class
GSASIIpwd.
fcjde_gen
(*args, **kwargs)[source]¶ Finger-Cox-Jephcoat D(2phi,2th) function for S/L = H/L Ref: J. Appl. Cryst. (1994) 27, 892-900.
Parameters: - x – array -1 to 1
- t – 2-theta position of peak
- s – sum(S/L,H/L); S: sample height, H: detector opening, L: sample to detector opening distance
- dx – 2-theta step size in deg
Returns: for fcj.pdf
T = x*dx+t
s = S/L+H/L
if x < 0:
fcj.pdf = [1/sqrt({cos(T)**2/cos(t)**2}-1) - 1/s]/|cos(T)|
if x >= 0: fcj.pdf = 0
-
GSASIIpwd.
findPDFfit
()[source]¶ Find if PDFfit2 is installed (may be local to GSAS-II). Does the following: :returns: the full path to a python executable or None, if it was not found.
-
GSASIIpwd.
findfullrmc
()[source]¶ Find where fullrmc is installed. Tries the following:
- Returns the Config var ‘fullrmc_exec’, if defined. No check is done that the interpreter has fullrmc
- The current Python interpreter if fullrmc can be imported and fullrmc is version 5+
- The path is checked for a fullrmc image as named by Bachir
Returns: the full path to a python executable that is assumed to have fullrmc installed or None, if it was not found.
-
GSASIIpwd.
getBackground
(pfx, parmDict, bakType, dataType, xdata, fixback=None)[source]¶ Computes the background from vars pulled from gpx file or tree.
-
GSASIIpwd.
getBackgroundDerv
(hfx, parmDict, bakType, dataType, xdata, fixback=None)[source]¶ needs a doc string
-
GSASIIpwd.
getEpsVoigt
(pos, alp, bet, sig, gam, xdata)[source]¶ Compute the double exponential Pseudo-Voigt convolution function for a neutron TOF & CW pink peak in external Fortran routine
-
GSASIIpwd.
getFCJVoigt
(pos, intens, sig, gam, shl, xdata)[source]¶ Compute the Finger-Cox-Jepcoat modified Voigt function for a CW powder peak by direct convolution. This version is not used.
-
GSASIIpwd.
getFCJVoigt3
(pos, sig, gam, shl, xdata)[source]¶ Compute the Finger-Cox-Jepcoat modified Pseudo-Voigt function for a CW powder peak in external Fortran routine
param pos: peak position in degrees param sig: Gaussian variance in centideg^2 param gam: Lorentzian width in centideg param shl: axial divergence parameter (S+H)/L param xdata: array; profile points for peak to be calculated; bounded by 20FWHM to 50FWHM (or vv)
returns: array: calculated peak function at each xdata returns: integral of peak; nominally = 1.0
-
GSASIIpwd.
getFWHM
(pos, Inst)[source]¶ Compute total FWHM from Thompson, Cox & Hastings (1987) , J. Appl. Cryst. 20, 79-83 via getgamFW(g,s).
Parameters: - pos – float peak position in deg 2-theta or tof in musec
- Inst – dict instrument parameters
Returns float: total FWHM of pseudoVoigt in deg or musec
-
GSASIIpwd.
getHKLpeak
(dmin, SGData, A, Inst=None, nodup=False)[source]¶ Generates allowed by symmetry reflections with d >= dmin NB: GenHKLf & checkMagextc return True for extinct reflections
Parameters: - dmin – minimum d-spacing
- SGData – space group data obtained from SpcGroup
- A – lattice parameter terms A1-A6
- Inst – instrument parameter info
Returns: HKLs: np.array hkl, etc for allowed reflections
-
GSASIIpwd.
getPeakProfile
(dataType, parmDict, xdata, fixback, varyList, bakType)[source]¶ Computes the profiles from multiple peaks for individual peak fitting for powder patterns. NB: not used for Rietveld refinement
-
GSASIIpwd.
getPeakProfileDerv
(dataType, parmDict, xdata, fixback, varyList, bakType)[source]¶ Computes the profile derivatives for a powder pattern for single peak fitting
return: np.array([dMdx1,dMdx2,…]) in same order as varylist = backVary,insVary,peakVary order
NB: not used for Rietveld refinement
-
GSASIIpwd.
getPsVoigt
(pos, sig, gam, xdata)[source]¶ Compute the simple Pseudo-Voigt function for a small angle Bragg peak in external Fortran routine
param pos: peak position in degrees param sig: Gaussian variance in centideg^2 param gam: Lorentzian width in centideg
returns: array: calculated peak function at each xdata returns: integral of peak; nominally = 1.0
-
GSASIIpwd.
getWidthsCW
(pos, sig, gam, shl)[source]¶ Compute the peak widths used for computing the range of a peak for constant wavelength data. On low-angle side, 50 FWHM are used, on high-angle side 75 are used, high angle side extended for axial divergence (for peaks above 90 deg, these are reversed.)
Parameters: - pos – peak position; 2-theta in degrees
- sig – Gaussian peak variance in centideg^2
- gam – Lorentzian peak width in centidegrees
- shl – axial divergence parameter (S+H)/L
Returns: widths; [Gaussian sigma, Lorentzian gamma] in degrees, and low angle, high angle ends of peak; 20 FWHM & 50 FWHM from position reversed for 2-theta > 90 deg.
-
GSASIIpwd.
getWidthsED
(pos, sig)[source]¶ Compute the peak widths used for computing the range of a peak for energy dispersive data. On low-energy side, 20 FWHM are used, on high-energy side 20 are used
Parameters: - pos – peak position; energy in keV (not used)
- sig – Gaussian peak variance in keV^2
Returns: widths; [Gaussian sigma] in keV, and low angle, high angle ends of peak; 20 FWHM & 50 FWHM from position
-
GSASIIpwd.
getWidthsTOF
(pos, alp, bet, sig, gam)[source]¶ Compute the peak widths used for computing the range of a peak for constant wavelength data. 50 FWHM are used on both sides each extended by exponential coeff.
param pos: peak position; TOF in musec (not used) param alp,bet: TOF peak exponential rise & decay parameters param sig: Gaussian peak variance in musec^2 param gam: Lorentzian peak width in musec
returns: widths; [Gaussian sigma, Lornetzian gamma] in musec returns: low TOF, high TOF ends of peak; 50FWHM from position
-
GSASIIpwd.
getdEpsVoigt
(pos, alp, bet, sig, gam, xdata)[source]¶ Compute the double exponential Pseudo-Voigt convolution function derivatives for a neutron TOF & CW pink peak in external Fortran routine
-
GSASIIpwd.
getdFCJVoigt3
(pos, sig, gam, shl, xdata)[source]¶ Compute analytic derivatives the Finger-Cox-Jepcoat modified Pseudo-Voigt function for a CW powder peak
param pos: peak position in degrees param sig: Gaussian variance in centideg^2 param gam: Lorentzian width in centideg param shl: axial divergence parameter (S+H)/L param xdata: array; profile points for peak to be calculated; bounded by 20FWHM to 50FWHM (or vv)
returns: arrays: function and derivatives of pos, sig, gam, & shl
-
GSASIIpwd.
getdPsVoigt
(pos, sig, gam, xdata)[source]¶ Compute the simple Pseudo-Voigt function derivatives for a small angle Bragg peak peak in external Fortran routine
param pos: peak position in degrees param sig: Gaussian variance in centideg^2 param gam: Lorentzian width in centideg
returns: arrays: function and derivatives of pos, sig & gam NB: the pos derivative has the opposite sign compared to that in other profile functions
-
GSASIIpwd.
getgamFW
(g, s)[source]¶ Compute total FWHM from Thompson, Cox & Hastings (1987), J. Appl. Cryst. 20, 79-83 lambda fxn needs FWHM for both Gaussian & Lorentzian components
Parameters: - g – float Lorentzian gamma = FWHM(L)
- s – float Gaussian sig
Returns float: total FWHM of pseudoVoigt
-
GSASIIpwd.
makeFFTsizeList
(nmin=1, nmax=1023, thresh=15)[source]¶ Provide list of optimal data sizes for FFT calculations
Parameters: - nmin (int) – minimum data size >= 1
- nmax (int) – maximum data size > nmin
- thresh (int) – maximum prime factor allowed
Returns: list of data sizes where the maximum prime factor is < thresh
-
GSASIIpwd.
makeMEMfile
(data, reflData, MEMtype, DYSNOMIA)[source]¶ make Dysnomia .mem file of reflection data, etc.
Parameters: - data (dict) – GSAS-II phase data
- reflData (list) – GSAS-II reflection data
- MEMtype (int) – 1 for neutron data with negative scattering lengths 0 otherwise
- DYSNOMIA (str) – path to dysnomia.exe
-
GSASIIpwd.
makePRFfile
(data, MEMtype)[source]¶ makes Dysnomia .prf control file from Dysnomia GUI controls
Parameters: - data (dict) – GSAS-II phase data
- MEMtype (int) – 1 for neutron data with negative scattering lengths 0 otherwise
Returns str: name of Dysnomia control file
-
GSASIIpwd.
peakInstPrmMode
= True¶ Determines the mode used for peak fitting. When peakInstPrmMode=True peak width parameters are computed from the instrument parameters (UVW,… or alpha,… etc) unless the individual parameter is refined. This allows the instrument parameters to be refined. When peakInstPrmMode=False, the instrument parameters are not used and cannot be refined. The default is peakFitMode=True. This is changed only in
setPeakInstPrmMode()
, which is called only fromGSASIIscriptable
.
-
GSASIIpwd.
setPeakInstPrmMode
(normal=True)[source]¶ Determines the mode used for peak fitting. If normal=True (default) peak width parameters are computed from the instrument parameters unless the individual parameter is refined. If normal=False, peak widths are used as supplied for each peak. At present this is called only in G2scriptable.
Note that normal=True unless this routine is called. Also, instrument parameters can only be refined with normal=True.
Parameters: normal (bool) – setting to apply to global variable peakInstPrmMode