#!/usr/bin/env python
##############################################################################
#
# pdffit2 by DANSE Diffraction group
# Simon J. L. Billinge
# (c) 2006 trustees of the Michigan State University.
# All rights reserved.
#
# File coded by: Chris Farros, Pavol Juhas
#
# See AUTHORS.txt for a list of people who contributed.
# See LICENSE.txt for license information.
#
##############################################################################
"""PdfFit class for fitting a structural model to PDF data."""
from __future__ import print_function
import itertools
import numbers
# Load pdffit2 and output modules to the current namespace.
# Note that "import diffpy.pdffit2.output as output" would
# crash with AttributeError when executed during imports of
# parent packages.
from diffpy.pdffit2 import output, pdffit2
# helper routines
def _format_value_std(value, stdev):
"""Convert value to a string with standard deviation in brackets.
value -- a number
stdev -- standard deviation. Ignored when small compared to value.
Return string.
"""
if stdev > abs(value) * 1e-8:
s = "%g (%g)" % (value, stdev)
elif str(stdev) == "nan":
s = "%g (NaN)" % value
else:
s = "%g" % value
return s
def _format_bond_length(dij, ddij, ij1, symij):
"""Return string with formatted bond length info for a pair of atoms.
dij -- distance between atoms i and j
ddij -- standard deviation of dij. Ignored when small relative to dij.
ij1 -- tuple of atom indices starting at 1
symij -- tuple of atom symbols
Return formatted string.
"""
w_smbidx = 10
w_equals = 30
s0 = "%s (#%i)" % (symij[0], ij1[0])
s1 = "%s (#%i)" % (symij[1], ij1[1])
leader0 = " " + s0.ljust(w_smbidx) + " - " + s1 + " "
leader1 = leader0.ljust(w_equals) + "= "
s = leader1 + _format_value_std(dij, ddij) + " A"
return s
def _convertCallable(var):
"""Convert an object to the result of its call when callable.
var -- string or callable object that returns string
Return var or var().
"""
if callable(var):
rv = var()
else:
rv = var
return rv
# constants
__intro_message__ = """
******************************************************************************
* P D F F I T Version %(version)s *
* %(date)s *
* -------------------------------------------------------------------------- *
* (c) 1998-2007 Trustees of the Michigan State University. *
* (c) 2008-%(year)s Trustees of the Columbia University *
* in the city of New York. *
* *
* Authors: *
* Thomas Proffen - Email: tproffen@lanl.gov *
* Jacques Bloch - Email: bloch@pa.msu.edu *
* Christopher Farrow - Email: clf2121@columbia.edu *
* Pavol Juhas - Email: pjuhas@bnl.gov *
* Simon Billinge - Email: sb2896@columbia.edu *
******************************************************************************
"""
##############################################################################
[docs]
class PdfFit(object):
"""Class for handling PdfFit calculations and refinements.
Attributes
----------
stru_files : list
The list to store structure files.
data_files : list
The list to store data files.
"""
# constants and enumerators from pdffit.h:
# selection of all atoms
selalias = {"ALL": -1}
# constraint type identifiers
FCON = {"USER": 0, "IDENT": 1, "FCOMP": 2, "FSQR": 3}
# scattering type identifiers
Sctp = {"X": 0, "N": 1}
def _exportAll(self, namespace):
"""_exportAll(self, namespace) --> Export all 'public' class methods
into namespace.
This function allows for a module-level PdfFit object which
doesn't have to be referenced when calling a method. This
function makes old (python) scripts compatible with this class.
At the top of the script, create a pdffit object, and then call
this method. Usually, namespace = locals().
"""
# string aliases (var = "var")
for a in itertools.chain(self.selalias, self.FCON, self.Sctp):
exec("%s = %r" % (a, a), namespace)
public = [
a
for a in dir(self)
if "__" not in a and a not in ["_handle", "_exportAll", "selalias", "FCON", "Sctp"]
]
for funcname in public:
namespace[funcname] = getattr(self, funcname)
return
[docs]
def intro():
"""Show introductory message."""
import re
from diffpy.pdffit2 import __date__, __version__
date = __date__[:10]
d = {"version": __version__, "date": date, "year": date[:4] or "2019"}
msg = __intro_message__ % d
def filler(mx):
return mx.group(0).rstrip(" *").ljust(77) + "*"
msg_ljust = re.sub("(?m)^(.{1,77}|.{79}.*)$", filler, msg)
print(msg_ljust, file=output.stdout)
return
intro = staticmethod(intro)
[docs]
def add_structure(self, stru):
"""add_structure(stru) --> Add new structure to PdfFit instance.
stru -- instance of Structure class from diffpy.structure.
No return value.
Raises pdffit2.structureError when stru contains unknown
atom species.
"""
s = stru.writeStr("pdffit")
self.read_struct_string(s)
return
[docs]
def read_struct(self, struct):
"""read_struct(struct) --> Read structure from file into memory.
struct -- name of file from which to read structure
Raises:
pdffit2.calculationError when a lattice cannot be created from the
given structure
pdffit2.structureError when a structure file is malformed
IOError when the file cannot be read from the disk
"""
pdffit2.read_struct(self._handle, struct)
self.stru_files.append(struct)
return
[docs]
def read_struct_string(self, struct, name=""):
"""read_struct_string(struct, name = "") --> Read structure from a
string into memory.
struct -- string containing the contents of the structure file
name -- tag with which to label structure
Raises:
pdffit2.calculationError when a lattice cannot be created from the
given structure
pdffit2.structureError when a structure file is malformed
"""
pdffit2.read_struct_string(self._handle, struct)
self.stru_files.append(name)
return
[docs]
def read_data(self, data, stype, qmax, qdamp):
"""read_data(data, stype, qmax, qdamp) --> Read pdf data from file into
memory.
data -- name of file from which to read data
stype -- 'X' (xray) or 'N' (neutron)
qmax -- Q-value cutoff used in PDF calculation.
Use qmax=0 to neglect termination ripples.
qdamp -- instrumental Q-resolution factor
Raises: IOError when the file cannot be read from disk
"""
pdffit2.read_data(self._handle, data, stype.encode(), qmax, qdamp)
self.data_files.append(data)
return
[docs]
def read_data_string(self, data, stype, qmax, qdamp, name=""):
"""read_data_string(data, stype, qmax, qdamp, name = "") --> Read pdf
data from a string into memory.
data -- string containing the contents of the data file
stype -- 'X' (xray) or 'N' (neutron)
qmax -- Q-value cutoff used in PDF calculation.
Use qmax=0 to neglect termination ripples.
qdamp -- instrumental Q-resolution factor
name -- tag with which to label data
"""
pdffit2.read_data_string(self._handle, data, stype.encode(), qmax, qdamp, name)
name = data
self.data_files.append(name)
return
[docs]
def read_data_lists(self, stype, qmax, qdamp, r_data, Gr_data, dGr_data=None, name="list"):
"""read_data_lists(stype, qmax, qdamp, r_data, Gr_data, dGr_data =
None, name = "list") --> Read pdf data into memory from lists.
All lists must be of the same length.
stype -- 'X' (xray) or 'N' (neutron)
qmax -- Q-value cutoff used in PDF calculation.
Use qmax=0 to neglect termination ripples.
qdamp -- instrumental Q-resolution factor
r_data -- list of r-values
Gr_data -- list of G(r) values
dGr_data -- list of G(r) uncertainty values
name -- tag with which to label data
Raises: ValueError when the data lists are of different length
"""
pdffit2.read_data_arrays(self._handle, stype.encode(), qmax, qdamp, r_data, Gr_data, dGr_data, name)
self.data_files.append(name)
return
[docs]
def pdfrange(self, iset, rmin, rmax):
"""Pdfrange(iset, rmin, rmax) --> Set the range of the fit.
iset -- data set to consider
rmin -- minimum r-value of fit
rmax -- maximum r-value of fit
Raises: ValueError for bad input values
"""
pdffit2.pdfrange(self._handle, iset, rmin, rmax)
return
[docs]
def reset(self):
"""Reset() --> Clear all stored fit, structure, and parameter data."""
self.stru_files = []
self.data_files = []
pdffit2.reset(self._handle)
return
[docs]
def alloc(self, stype, qmax, qdamp, rmin, rmax, bin):
"""Alloc(stype, qmax, qdamp, rmin, rmax, bin) --> Allocate space for a
PDF calculation.
The structure from which to calculate the PDF must first be imported
with the read_struct() or read_struct_string() method.
stype -- 'X' (xray) or 'N' (neutron)
qmax -- Q-value cutoff used in PDF calculation.
Use qmax=0 to neglect termination ripples.
qdamp -- instrumental Q-resolution factor
rmin -- minimum r-value of calculation
rmax -- maximum r-value of calculation
bin -- number of data points in calculation
Raises:
ValueError for bad input values
pdffit.unassignedError when no structure has been loaded
"""
pdffit2.alloc(self._handle, stype.encode(), qmax, qdamp, rmin, rmax, bin)
return
[docs]
def calc(self):
"""Calc() --> Calculate the PDF of the imported structure.
Space for the calculation must first be allocated with the alloc()
method.
Raises:
pdffit2.calculationError when allocated space cannot
accommodate calculation
pdffit.unassignedError when space for calculation has not been
allocated
"""
pdffit2.calc(self._handle)
return
[docs]
def refine(self, toler=0.00000001):
"""Refine(toler = 0.00000001) --> Fit the theory to the imported data.
toler -- tolerance of the fit
Raises:
pdffit2.calculationError when the model pdf cannot be calculated
pdffit2.constraintError when refinement fails due to bad
constraint
pdffit2.unassigedError when a constraint used but never initialized
using setpar()
"""
step = 0
finished = 0
while not finished:
finished = pdffit2.refine_step(self._handle, toler)
step += 1
return
[docs]
def refine_step(self, toler=0.00000001):
"""refine_step(toler = 0.00000001) --> Run a single step of the fit.
toler -- tolerance of the fit
Raises:
pdffit2.calculationError when the model pdf cannot be calculated
pdffit2.constraintError when refinement fails due to bad
constraint
pdffit2.unassigedError when a constraint used but never initialized
using setpar()
Returns: 1 (0) if refinement is (is not) finished
"""
self.finished = pdffit2.refine_step(self._handle, toler)
return self.finished
[docs]
def save_pdf(self, iset, fname):
"""save_pdf(iset, fname) --> Save calculated or fitted PDF to file.
iset -- data set to save
Raises:
IOError if file cannot be saved
pdffit2.unassignedError if the data set is undefined
"""
pdffit2.save_pdf(self._handle, iset, fname)
return
[docs]
def save_pdf_string(self, iset):
"""save_pdf_string(iset) --> Save calculated or fitted PDF to string.
iset -- data set to save
Raises:
pdffit2.unassignedError if the data set is undefined
Returns: string containing contents of save file
"""
pdffilestring = pdffit2.save_pdf(self._handle, iset, "")
return pdffilestring
[docs]
def save_dif(self, iset, fname):
"""save_dif(iset, fname) --> Save data and fitted PDF difference to
file.
iset -- data set to save
Raises:
IOError if file cannot be saved
pdffit2.unassignedError if the data set is undefined
"""
pdffit2.save_dif(self._handle, iset, fname)
return
[docs]
def save_dif_string(self, iset):
"""save_dif_string(iset) --> Save data and fitted PDF difference to
string.
iset -- data set to save
Raises:
pdffit2.unassignedError if the data set is undefined
Returns: string containing contents of save file
"""
diffilestring = pdffit2.save_dif(self._handle, iset, "")
return diffilestring
[docs]
def save_res(self, fname):
"""save_res(fname) --> Save fit-specific data to file.
Raises:
IOError if file cannot be saved
pdffit2.unassignedError if there is no refinement data to save
"""
pdffit2.save_res(self._handle, fname)
return
[docs]
def save_res_string(self):
"""save_res_string() --> Save fit-specific data to a string.
Raises:
pdffit2.unassignedError if there is no refinement data to save
Returns: string containing contents of save file
"""
resfilestring = pdffit2.save_res(self._handle, "")
return resfilestring
[docs]
def get_structure(self, ip):
"""get_structure(ip) --> Get a copy of specified phase data.
ip -- index of existing PdfFit phase starting from 1
Return Structure object from diffpy.structure.
Raise pdffit2.unassignedError if phase ip is undefined.
"""
from diffpy.structure import PDFFitStructure
s = self.save_struct_string(ip)
stru = PDFFitStructure()
stru.readStr(s, "pdffit")
return stru
[docs]
def save_struct(self, ip, fname):
"""save_struct(ip, fname) --> Save structure resulting from fit to
file.
ip -- phase to save
Raises:
IOError if file cannot be saved
pdffit2.unassignedError if the data set is undefined
"""
pdffit2.save_struct(self._handle, ip, fname)
return
[docs]
def save_struct_string(self, ip):
"""save_struct(ip) --> Save structure resulting from fit to string.
ip -- phase to save
Raises:
pdffit2.unassignedError if phase ip is undefined.
Returns: string containing contents of save file
"""
structfilestring = pdffit2.save_struct(self._handle, ip, "")
return structfilestring
[docs]
def show_struct(self, ip):
"""show_struct(ip) --> Print structure resulting from fit.
ip -- phase to display
Raises: pdffit2.unassignedError if the phase is undefined
"""
pdffit2.show_struct(self._handle, ip)
return
[docs]
def constrain(self, var, par, fcon=None):
"""Constrain(var, par[, fcon]) --> Constrain a variable to a parameter.
A variable can be constrained to a number or equation string.
var -- variable to constrain, such as x(1)
par -- parameter which to constrain the variable. This can be
an integer or an equation string containing a reference
to another parameter. Equation strings use standard c++
syntax. The value of a constrained parameter is accessed
as @p in an equation string, where p is the parameter.
e.g.
>>> constrain(x(1), 1)
>>> constrain(x(2), "0.5+@1")
fcon -- 'USER', 'IDENT', 'FCOMP', or 'FSQR'
this is an optional parameter, and I don't know how it is
used!
Raises:
pdffit2.constraintError if a constraint is bad
pdffit2.unassignedError if variable does not yet exist
ValueError if variable index does not exist (e.g. lat(7))
"""
var_ref = self.__getRef(var)
varnc = _convertCallable(var)
if fcon:
fc = self.FCON[fcon]
pdffit2.constrain_int(self._handle, var_ref, varnc, par, fc)
elif isinstance(par, str):
pdffit2.constrain_str(self._handle, var_ref, varnc, par)
else:
pdffit2.constrain_int(self._handle, var_ref, varnc, par)
return
[docs]
def setpar(self, par, val):
"""Setpar(par, val) --> Set value of constrained parameter.
val -- Either a numerical value or a reference to a variable
Raises:
pdffit2.unassignedError when variable is yet to be assigned
"""
# people do not use parenthesis, e.g., "setpar(3, qdamp)"
# in such case val is a reference to PdfFit method
val = _convertCallable(val)
try:
val = float(val)
pdffit2.setpar_dbl(self._handle, par, val)
except ValueError:
var_ref = self.__getRef(val)
pdffit2.setpar_RV(self._handle, par, var_ref)
return
[docs]
def setvar(self, var, val):
"""Setvar(var, val) --> Set the value of a variable.
Raises:
pdffit2.unassignedError if variable does not yet exist
ValueError if variable index does not exist (e.g. lat(7))
"""
var_ref = self.__getRef(var)
pdffit2.setvar(self._handle, var_ref, val)
return
[docs]
def getvar(self, var):
"""Getvar(var) --> Get stored value of a variable.
Raises:
pdffit2.unassignedError if variable does not yet exist
ValueError if variable index does not exist (e.g. lat(7))
"""
var_ref = self.__getRef(var)
retval = pdffit2.getvar(self._handle, var_ref)
return retval
[docs]
def getrw(self):
"""Getrw() --> Get normalized total error of the fit rw.
getrw calculates total fit error summed for all datasets in the
fit.
Return float.
"""
rw = pdffit2.getrw(self._handle)
return rw
[docs]
def getcrw(self):
"""Getcrw() --> Get cumulative Rw for the current dataset.
Cumulative Rw is a list of Rw partial sums cost values evaluated against
observed PDF data in the error sums evaluated against
the r-points in the fit.
Raises: pdffit2.unassignedError if no data exists
Returns: List of crw points, equidistant in r or empty list
if the refine function has not been called yet.
"""
crw = pdffit2.getcrw(self._handle)
return crw
[docs]
def getR(self):
"""GetR() --> Get r-points used in the fit.
This function should only be called after data has been loaded or
calculated. Before a refinement, the list of r-points will reflect the
data. Afterwards, they will reflect the fit range.
Raises: pdffit2.unassignedError if no data exists
Returns: List of equidistance r-points used in fit.
"""
R = pdffit2.getR(self._handle)
return R
[docs]
def getpdf_fit(self):
"""getpdf_fit() --> Get fitted PDF.
This function should only be called after a refinement or refinement
step has been done.
Raises: pdffit2.unassignedError if no data exists
Returns: List of fitted points, equidistant in r.
"""
pdfdata = pdffit2.getpdf_fit(self._handle)
return pdfdata
[docs]
def getpdf_obs(self):
"""getpdf_obs() --> Get observed PDF.
This function should only be called after data has been loaded or
calculated. Before a refinement, the list of r-points will reflect the
data. Afterwards, they will reflect the fit range.
Raises: pdffit2.unassignedError if no data exists
Returns: List of data points, equidistant in r.
"""
pdfdata = pdffit2.getpdf_obs(self._handle)
return pdfdata
[docs]
def getpdf_diff(self):
"""Obtain difference between observed and fitted PDF.
This function should only be called after data has been loaded or
calculated. Before a refinement, the list of r-points will reflect the
data. Afterwards, they will reflect the fit range.
Raises: pdffit2.unassignedError if no data exists
Returns: List of data points, equidistant in r.
"""
Gdiff = pdffit2.getpdf_diff(self._handle)
return Gdiff
[docs]
def get_atoms(self, ip=None):
"""get_atoms() --> Get element symbols of all atoms in the structure.
ip -- index of phase to get the elements from (starting from 1)
when ip is not given, use current phase
This function should only be called after a structure has been loaded.
Raises: pdffit2.unassignedError if no structure exists
Returns: List of atom names in structure.
"""
if ip is None:
rv = pdffit2.get_atoms(self._handle)
else:
rv = pdffit2.get_atoms(self._handle, ip)
return rv
[docs]
def get_atom_types(self, ip=None):
"""get_atom_types() --> Ordered unique element symbols in the
structure.
ip -- index of phase to get the elements from (starting from 1)
when ip is not given, use current phase
This function should only be called after a structure has been loaded.
Raises:
pdffit2.unassignedError if no structure exists
Returns: List of unique atom symbols as they occur in structure.
"""
if ip is None:
rv = pdffit2.get_atom_types(self._handle)
else:
rv = pdffit2.get_atom_types(self._handle, ip)
return rv
[docs]
def getpar(self, par):
"""Getpar(par) --> Get value of parameter.
Raises: ValueError if parameter does not exists
"""
return pdffit2.getpar(self._handle, par)
[docs]
def fixpar(self, par):
"""Fixpar(par) --> Fix a parameter.
Fixed parameters are not fitted in a refinement. Passed parameter
can be 'ALL', in which case all parameters are fixed.
Raises: pdffit.unassignedError when parameter has not been assigned
"""
if isinstance(par, str) and par.upper() in self.selalias:
par = self.selalias[par.upper()]
pdffit2.fixpar(self._handle, par)
return
[docs]
def freepar(self, par):
"""Freepar(par) --> Free a parameter.
Freed parameters are fitted in a refinement. Passed parameter
can be 'ALL', in which case all parameters are freed.
Raises: pdffit.unassignedError when parameter has not been assigned
"""
if isinstance(par, str) and par.upper() in self.selalias:
par = self.selalias[par.upper()]
pdffit2.freepar(self._handle, par)
return
[docs]
def setphase(self, ip):
"""Setphase(ip) --> Switch to phase ip.
ip -- index of the phase starting at 1.
All parameters assigned after this method is called refer only to the
current phase.
Raises: pdffit.unassignedError when phase does not exist
"""
pdffit2.setphase(self._handle, ip)
return
[docs]
def setdata(self, iset):
"""Setdata(iset) --> Set the data set in focus.
iset -- integer index of data set starting at 1.
Raises: pdffit.unassignedError when data set does not exist
"""
pdffit2.setdata(self._handle, iset)
return
[docs]
def psel(self, ip):
"""Psel(ip) --> Include phase ip in calculation of total PDF.
psel('ALL') selects all phases for PDF calculation.
Raises: pdffit2.unassignedError if selected phase does not exist
"""
if isinstance(ip, str) and ip.upper() in self.selalias:
ip = self.selalias[ip.upper()]
pdffit2.psel(self._handle, ip)
return
[docs]
def pdesel(self, ip):
"""Pdesel(ip) --> Exclude phase ip from calculation of total PDF.
pdesel('ALL') excludes all phases from PDF calculation.
Raises: pdffit2.unassignedError if selected phase does not exist
"""
if isinstance(ip, str) and ip.upper() in self.selalias:
ip = self.selalias[ip.upper()]
pdffit2.pdesel(self._handle, ip)
return
[docs]
def selectAtomType(self, ip, ijchar, symbol, flag):
"""Configure partial PDF - mark the specified atom type in phase ip
as included or excluded as a first or second in pair for distance
evaluation.
ip -- phase index starting at 1
ijchar -- 'i' or 'j' for first or second in pair
symbol -- element symbol
flag -- bool flag, True for selection, False for exclusion
Raises:
pdffit2.unassignedError if selected phase does not exist
ValueError for invalid value of ijchar
"""
pdffit2.selectAtomType(self._handle, ip, ijchar.encode(), symbol, flag)
return
[docs]
def selectAtomIndex(self, ip, ijchar, aidx, flag):
"""Configure partial PDF - mark the atom of given index in phase ip
as included or excluded as a first or second in pair for distance
evaluation.
ip -- phase index starting at 1
ijchar -- 'i' or 'j' for first or second in pair
aidx -- integer index of atom starting at 1
flag -- bool flag, True for selection, False for exclusion
Raises:
pdffit2.unassignedError if selected phase does not exist
ValueError if atom index or ijchar are invalid
"""
pdffit2.selectAtomIndex(self._handle, ip, ijchar.encode(), aidx, flag)
return
[docs]
def selectAll(self, ip, ijchar):
"""Configure partial PDF - include all atoms of phase ip as first or
second element in pair for distance evaluation.
ip -- phase index starting at 1
ijchar -- 'i' or 'j' for first or second in pair
Raises:
pdffit2.unassignedError if selected phase does not exist
ValueError if ijchar is invalid
"""
pdffit2.selectAll(self._handle, ip, ijchar.encode())
return
[docs]
def selectNone(self, ip, ijchar):
"""Configure partial PDF - exclude all atoms of phase ip from first
or second element of pair distance evaluation.
ip -- phase index starting at 1
ijchar -- 'i' or 'j' for first or second in pair
Raises:
pdffit2.unassignedError if selected phase does not exist
ValueError if ijchar is invalid
"""
pdffit2.selectNone(self._handle, ip, ijchar.encode())
return
[docs]
def bang(self, i, j, k):
"""Bang(i, j, k) --> Show bond angle defined by atoms i, j, k.
No return value. Use bond_angle() to get the result.
Raises: ValueError if selected atom(s) does not exist
pdffit.unassignedError when no structure has been loaded
"""
angle, stdev = pdffit2.bond_angle(self._handle, i, j, k)
# indices should be already checked here by bond_angle
atom_symbols = self.get_atoms()
leader = " %s (#%i) - %s (#%i) - %s (#%i) = " % (
atom_symbols[i - 1],
i,
atom_symbols[j - 1],
j,
atom_symbols[k - 1],
k,
)
s = leader + _format_value_std(angle, stdev) + " degrees"
print(s, file=output.stdout)
return
[docs]
def bond_angle(self, i, j, k):
"""bond_angle(i, j, k) --> bond angle defined by atoms i, j, k. Angle
is calculated using the shortest ji and jk lengths with respect to
periodic boundary conditions.
i, j, k -- atom indices starting at 1
Return a tuple of (angle, angle_error), both values are in degrees.
Raises: ValueError if selected atom(s) does not exist
pdffit.unassignedError when no structure has been loaded
"""
rv = pdffit2.bond_angle(self._handle, i, j, k)
return rv
[docs]
def blen(self, *args):
"""Blen(i, j) --> Show bond length defined by atoms i and j.
i -- index of the first atom starting at 1
j -- index of the second atom starting at 1
No return value. Use bond_length_atoms() to retrieve result.
Second form:
blen(a1, a2, lb, ub) --> Show sorted lengths of all a1-a2 bonds.
a1 -- symbol of the first element in pair or "ALL"
a2 -- symbol of the second element in pair or "ALL"
lb -- lower bond length boundary
ub -- upper bond length boundary
No return value. Use bond_length_types() to retrieve results.
Raises: ValueError if selected atom(s) does not exist
pdffit.unassignedError when no structure has been loaded
"""
# first form
if len(args) == 2:
dij, ddij = self.bond_length_atoms(*args[0:2])
atom_symbols = self.get_atoms()
ij = (args[0], args[1])
# indices were already checked in bond_length_atoms call
assert (0 <= min(ij) - 1) and (max(ij) - 1 < len(atom_symbols))
symij = (atom_symbols[ij[0] - 1].upper(), atom_symbols[ij[1] - 1].upper())
print(_format_bond_length(dij, ddij, ij, symij), file=output.stdout)
# second form
elif len(args) == 4:
a1, a2, lb, ub = args
try:
atom_types = self.get_atom_types()
if isinstance(a1, numbers.Integral):
a1 = atom_types[a1 - 1]
if isinstance(a2, numbers.Integral):
a2 = atom_types[a2 - 1]
except IndexError:
# index of non-existent atom type
return
# arguments are OK here, get bond length dictionary
bld = pdffit2.bond_length_types(self._handle, a1, a2, lb, ub)
s = "(%s,%s) bond lengths in [%gA,%gA] for current phase :" % (a1, a2, lb, ub)
print(s, file=output.stdout)
atom_symbols = self.get_atoms()
npts = len(bld["dij"])
for idx in range(npts):
dij = bld["dij"][idx]
ddij = bld["ddij"][idx]
ij0 = bld["ij0"][idx]
ij1 = bld["ij1"][idx]
symij = (atom_symbols[ij0[0]], atom_symbols[ij0[1]])
s = _format_bond_length(dij, ddij, ij1, symij)
print(s, file=output.stdout)
print(file=output.stdout)
if not bld["dij"]:
print(" *** No pairs found ***", file=output.stdout)
else:
emsg = "blen() takes 2 or 4 arguments (%i given)" % len(args)
raise TypeError(emsg)
# done
return
[docs]
def bond_length_atoms(self, i, j):
"""bond_length_atoms(i, j) --> shortest distance between atoms i, j.
Periodic boundary conditions are applied to find the shortest bond.
i -- index of the first atom starting at 1
j -- index of the second atom starting at 1
Return a tuple of (distance, distance_error).
Raises: ValueError if selected atom(s) does not exist
pdffit.unassignedError when no structure has been loaded.
"""
rv = pdffit2.bond_length_atoms(self._handle, i, j)
return rv
[docs]
def bond_length_types(self, a1, a2, lb, ub):
"""bond_length_types(a1, a2, lb, ub) --> get all a1-a2 distances.
a1 -- symbol of the first element in pair or "ALL"
a2 -- symbol of the second element in pair or "ALL"
lb -- lower bond length boundary
ub -- upper bond length boundary
Return a dictionary of distance data containing
dij : list of bond lengths within given bounds
ddij : list of bond length standard deviations
ij0 : pairs of atom indices starting from 0
ij1 : pairs of atom indices starting from 1
Raises: ValueError if selected atom(s) does not exist
pdffit.unassignedError when no structure has been loaded.
"""
rv = pdffit2.bond_length_types(self._handle, a1, a2, lb, ub)
return rv
[docs]
def show_scat(self, stype):
"""show_scat(stype) --> Print scattering length for all atoms in the
current phase.
stype -- 'X' (xray) or 'N' (neutron).
Raises: pdffit2.unassignedError if no phase exists
"""
print(self.get_scat_string(stype), file=output.stdout)
return
[docs]
def get_scat_string(self, stype):
"""get_scat_string(stype) --> Get string with scattering factors of all
atoms in the current phase.
stype -- 'X' (xray) or 'N' (neutron).
Raises:
pdffit2.unassignedError if no phase exists
Returns: string with all scattering factors.
"""
return pdffit2.get_scat_string(self._handle, stype.encode())
[docs]
def get_scat(self, stype, element):
"""get_scat(stype, element) --> Get active scattering factor for given
element. If scattering factor has been changed using set_scat the
result may depend on the active phase. When no phase has been loaded,
return the standard value.
stype -- 'X' (xray) or 'N' (neutron).
element -- case-insensitive element symbol such as "Na" or "CL"
Return float.
Raises:
ValueError if element is not known.
"""
rv = pdffit2.get_scat(self._handle, stype.encode(), element)
return rv
[docs]
def set_scat(self, stype, element, value):
"""set_scat(stype, element, value) --> Set custom scattering factor for
given element. The new scattering factor applies only for the current
phase, in other phases it keeps its default value.
stype -- 'X' (xray) or 'N' (neutron).
element -- case-insensitive element symbol such as "Na" or "CL"
value -- new value of scattering factor
No return value.
Raises:
pdffit2.unassignedError if no phase exists.
ValueError if element is not known.
See also reset_scat, get_scat.
"""
pdffit2.set_scat(self._handle, stype.encode(), element, value)
return
[docs]
def reset_scat(self, element):
"""reset_scat(stype, element) --> Reset scattering factors for given
element to their standard values. The reset_scat applies only for the
current phase.
element -- case-insensitive element symbol such as "Na" or "CL"
Raises:
pdffit2.unassignedError if no phase exists
ValueError if element is not known.
"""
pdffit2.reset_scat(self._handle, element)
return
[docs]
def num_atoms(self):
"""num_atoms() --> Get number of atoms in current phase.
Raises: pdffit2.unassignedError if no atoms exist
"""
return pdffit2.num_atoms(self._handle)
[docs]
def num_phases(self):
"""num_phases() --> Number of phases loaded in PdfFit instance.
Use setphase to bring a specific phase in focus.
Return integer.
"""
n = pdffit2.num_phases(self._handle)
return n
[docs]
def num_datasets(self):
"""num_datasets() --> Number of datasets loaded in PdfFit instance.
Use setdata to bring a specific dataset in focus.
Return integer.
"""
n = pdffit2.num_datasets(self._handle)
return n
[docs]
def phase_fractions(self):
"""phase_fractions() --> relative phase fractions for current dataset.
Convert phase scale factors to relative phase fractions given the
scattering type of current dataset.
Return a dictionary of phase fractions with following keys:
"atom" -- list of fractions normalized to atom count
"stdatom" -- errors of atom count fractions
"cell" -- list of fractions normalized to unit cell count
"stdcell" -- errors of unit cell count fractions
"mass" -- list of relative weight fractions
"stdmass" -- errors of relative weight fractions
Raises: pdffit2.unassignedError if no dataset exists.
"""
return pdffit2.phase_fractions(self._handle)
# Begin refinable variables.
[docs]
def lat(n):
"""Lat(n) --> Get reference to lattice variable n.
n can be an integer or a string representing the lattice
variable. 1 <==> 'a' 2 <==> 'b' 3 <==> 'c' 4 <==> 'alpha' 5 <==>
'beta' 6 <==> 'gamma'
"""
LatParams = {"a": 1, "b": 2, "c": 3, "alpha": 4, "beta": 5, "gamma": 6}
if isinstance(n, str):
n = LatParams[n]
return "lat(%i)" % n
lat = staticmethod(lat)
[docs]
def x(i):
"""X(i) --> Get reference to x-value of atom i."""
return "x(%i)" % i
x = staticmethod(x)
[docs]
def y(i):
"""Y(i) --> Get reference to y-value of atom i."""
return "y(%i)" % i
y = staticmethod(y)
[docs]
def z(i):
"""Z(i) --> Get reference to z-value of atom i."""
return "z(%i)" % i
z = staticmethod(z)
[docs]
def u11(i):
"""U11(i) --> Get reference to U(1,1) for atom i.
U is the anisotropic thermal factor tensor.
"""
return "u11(%i)" % i
u11 = staticmethod(u11)
[docs]
def u22(i):
"""U22(i) --> Get reference to U(2,2) for atom i.
U is the anisotropic thermal factor tensor.
"""
return "u22(%i)" % i
u22 = staticmethod(u22)
[docs]
def u33(i):
"""U33(i) --> Get reference to U(3,3) for atom i.
U is the anisotropic thermal factor tensor.
"""
return "u33(%i)" % i
u33 = staticmethod(u33)
[docs]
def u12(i):
"""U12(i) --> Get reference to U(1,2) for atom i.
U is the anisotropic thermal factor tensor.
"""
return "u12(%i)" % i
u12 = staticmethod(u12)
[docs]
def u13(i):
"""U13(i) --> Get reference to U(1,3) for atom i.
U is the anisotropic thermal factor tensor.
"""
return "u13(%i)" % i
u13 = staticmethod(u13)
[docs]
def u23(i):
"""U23(i) --> Get reference to U(2,3) for atom i.
U is the anisotropic thermal factor tensor.
"""
return "u23(%i)" % i
u23 = staticmethod(u23)
[docs]
def occ(i):
"""Occ(i) --> Get reference to occupancy of atom i."""
return "occ(%i)" % i
occ = staticmethod(occ)
[docs]
def pscale():
"""Pscale() --> Get reference to pscale.
pscale is the fraction of the total structure that the current
phase represents.
"""
return "pscale"
pscale = staticmethod(pscale)
[docs]
def sratio():
"""Sratio() --> Get reference to sigma ratio.
The sigma ratio determines the reduction in the Debye-Waller
factor for distances below rcut.
"""
return "sratio"
sratio = staticmethod(sratio)
[docs]
def delta1():
"""Delta1() --> Get reference to 1/R peak sharpening factor."""
return "delta1"
delta1 = staticmethod(delta1)
[docs]
def delta2():
"""Delta2() --> Reference to (1/R^2) sharpening factor.
The phenomenological correlation constant in the Debye-Waller
factor. The (1/R^2) peak sharpening factor.
"""
return "delta2"
delta2 = staticmethod(delta2)
[docs]
def dscale():
"""Dscale() --> Get reference to dscale.
The data scale factor.
"""
return "dscale"
dscale = staticmethod(dscale)
[docs]
def qdamp():
"""Qdamp() --> Get reference to qdamp.
Qdamp controls PDF damping due to instrument Q-resolution.
"""
return "qdamp"
qdamp = staticmethod(qdamp)
[docs]
def qbroad():
"""Qbroad() --> Get reference to qbroad.
Quadratic peak broadening factor.
"""
return "qbroad"
qbroad = staticmethod(qbroad)
[docs]
def spdiameter():
"""Spdiameter() --> Get reference to spdiameter (phase property).
Diameter value for the spherical particle PDF correction.
Spherical envelope is not applied when spdiameter equals 0.
"""
return "spdiameter"
spdiameter = staticmethod(spdiameter)
[docs]
def stepcut():
"""Stepcut() --> Get reference to stepcut (phase property).
stepcut is cutoff radius for empirical step-function PDF
envelope. stepcut can be used to approximate loss of pair
correlations in amorphous phase. stepcut cannot be refined.
Step cutoff is not applied when stepcut equals 0.
"""
return "stepcut"
stepcut = staticmethod(stepcut)
[docs]
def rcut():
"""Rcut() --> Get reference to rcut.
rcut is the value of r below which peak sharpening, defined by
the sigma ratio (sratio), applies. rcut cannot be refined.
"""
return "rcut"
rcut = staticmethod(rcut)
# End refinable variables.
def __init__(self, create_intro=True):
"""Initialize the Pdffit class, create a new PdfFit object.
Parameters
----------
create_intro : bool, optional
The flag to control the display of an introduction message.
If True, display an introduction message, else not. Default is True.
"""
self.stru_files = []
self.data_files = []
self._handle = pdffit2.create()
if create_intro:
self.intro()
return
def __getRef(self, var_string):
"""Return the actual reference to the variable in the var_string.
This function must be called before trying to actually reference an
internal variable. See the constrain method for an example.
Raises:
pdffit2.unassignedError if variable is not yet assigned
ValueError if variable index does not exist (e.g. lat(7))
"""
var_string = _convertCallable(var_string)
arg_int = None
try:
method_string, arg_string = var_string.split("(")
method_string = method_string.strip()
arg_int = int(arg_string.strip(")").strip())
except ValueError: # There is no arg_string
method_string = var_string.strip()
f = getattr(pdffit2, method_string)
if arg_int is None:
retval = f(self._handle)
else:
retval = f(self._handle, arg_int)
return retval
# End of class PdfFit
# End of file