3.7.3. Radial Distribution Functions — MDAnalysis.analysis.rdf
¶
Tools for calculating pair distribution functions
- # TODO
- Structure factor?
- Coordination number
-
class
MDAnalysis.analysis.rdf.
InterRDF
(g1, g2, nbins=75, range=(0.0, 15.0), exclusion_block=None, start=None, stop=None, step=None)[source]¶ Intermolecular pair distribution function
InterRDF(g1, g2, nbins=75, range=(0.0, 15.0))
Parameters: - g1 – First AtomGroup
- g2 – Second AtomGroup
- Keywords –
- -------- –
- nbins – Number of bins in the histogram [75]
- range – The size of the RDF [0.0, 15.0]
- exclusion_block – A tuple representing the tile to exclude from the distance array. [None]
- start – The frame to start at [0]
- stop – The frame to end at [-1]
- step – The step size through the trajectory in frames [0]
Example
First create the InterRDF object, by supplying two AtomGroups then use the run method
rdf = InterRDF(ag1, ag2) rdf.run()Results are available through the .bins and .rdf attributes
plt.plot(rdf.bins, rdf.rdf)The exclusion_block keyword allows the masking of pairs from within the same molecule. For example, if there are 7 of each atom in each molecule, the exclusion mask (7, 7) can be used.
New in version 0.13.0.