rdf
#
- biotite.structure.rdf(center, atoms, selection=None, interval=(0, 10), bins=100, box=None, periodic=False)[source]#
Compute the radial distribution function g(r) (RDF) for one or multiple given central positions based on a given system of particles.
- Parameters:
- centerAtom or AtomArray or AtomArrayStack or ndarray, dtype=float
Coordinates or atoms(s) to use as origin(s) for RDF calculation.
If a single
Atom
or andarray
with shape (3,) is given, the RDF is only calculated for this position.If an
AtomArray
or andarray
with shape (n,3) is given, the calculated RDF histogram is an average over n postions.If an
AtomArrayStack
or andarray
with shape (m,n,3) is given, different centers are used for each model m. The calculated RDF histogram is an average over m models and n positions. This requires atoms to be anAtomArrayStack
.
- atomsAtomArray or AtomArrayStack
The distribution is calculated based on these atoms. When an an
AtomArrayStack
is provided, the RDF histogram is averaged over all models. Please not that atoms must have an associated box, unless box is set.- selectionndarray, dtype=bool, shape=(n,), optional
Boolean mask for atoms to limit the RDF calculation to specific atoms.
- intervaltuple, optional
The range in which the RDF is calculated.
- binsint or sequence of scalars or str, optional
Bins for the RDF.
If bins is an int, it defines the number of bins for the given interval.
If bins is a sequence, it defines the bin edges, ignoring the interval parameter. The output bins has the length of this input parameter reduced by one.
If bins is a string, it defines the function used to calculate the bins.
See numpy.histogram() for further details.
- boxndarray, shape=(3,3) or shape=(m,3,3), optional
If this parameter is set, the given box is used instead of the box attribute of atoms. Must have shape (3,3) if atoms is an
AtomArray
or (m,3,3) if atoms is anAtomArrayStack
, respectively.- periodicbool, optional
Defines if periodic boundary conditions are taken into account.
- Returns:
- binsndarray, dtype=float, shape=n
The centers of the histogram bins. The length of the array is given by the bins input parameter.
- rdfndarry, dtype=float, shape=n
RDF values for every bin.
Notes
Since the RDF depends on the average particle density of the system, this function strictly requires an box.
Examples
Calculate the oxygen-oxygen radial distribution function of water. The range of the histogram starts at 0.2 Å, in order to ignore the counts for the density for each oxygen to itself.
>>> from os.path import join >>> waterbox = load_structure(join(path_to_structures, "waterbox.gro")) >>> oxygens = waterbox[:, waterbox.atom_name == 'OW'] >>> bins, g_r = rdf(oxygens, oxygens, bins=49, interval=(0.2, 10), periodic=True)
Print the RDF depending on the radius. Bins are in Å.
>>> for x, y in zip(bins, g_r): ... print(f"{x:.2f} {y:.2f}") 0.30 0.00 0.50 0.00 0.70 0.04 0.90 0.02 1.10 0.03 1.30 0.06 1.50 0.03 1.70 0.04 1.90 0.04 2.10 0.04 2.30 0.04 2.50 0.16 2.70 1.99 2.90 2.22 3.10 1.34 3.30 1.04 3.50 0.97 3.70 0.94 3.90 0.97 4.10 0.94 4.30 0.98 4.50 0.97 4.70 0.96 4.90 0.99 5.10 0.99 5.30 1.02 5.50 1.02 5.70 0.99 5.90 0.98 6.10 0.98 6.30 0.99 6.50 1.02 6.70 1.02 6.90 1.00 7.10 1.01 7.30 1.01 7.50 1.00 7.70 1.01 7.90 0.99 8.10 0.99 8.30 0.99 8.50 0.99 8.70 0.99 8.90 1.00 9.10 1.01 9.30 1.01 9.50 1.00 9.70 1.00 9.90 0.99
Find the radius for the first solvation shell. In this simple case, the density peak is identified by finding the maximum of the function.
>>> peak_position = np.argmax(g_r) >>> print(f"{bins[peak_position]/10:.2f} nm") 0.29 nm