Author: | Oliver Beckstein |
---|---|
Year: | 2012 |
Copyright: | GNU Public License v2 |
New in version 0.7.7.
The module contains code to analyze root mean square quantities such as the RMSD or RMSF (not implemented yet).
This module uses the fast QCP algorithm [Theobald2005] to calculate the root mean square distance (RMSD) between two coordinate sets (as implemented in MDAnalysis.core.qcprot.CalcRMSDRotationalMatrix()).
When using this module in published work please cite [Theobald2005].
See also
Returns RMSD between two coordinate sets a and b.
a and b are arrays of the coordinates of N atoms of shape N*3 as generated by, e.g., MDAnalysis.core.AtomGroup.AtomGroup.coordinates().
The weights can be an array of length N, containing e.g. masses for a weighted RMSD calculation.
The function uses Douglas Theobald’s fast QCP algorithm [Theobald2005] to calculate the RMSD.
>>> u = Universe(PSF,DCD)
>>> bb = u.selectAtoms('backbone')
>>> A = bb.coordinates() # coordinates of first frame
>>> u.trajectory[-1] # forward to last frame
>>> B = bb.coordinates() # coordinates of last frame
>>> rmsd(A,B)
6.8342494129169804
Class to perform RMSD analysis on a trajectory.
Run the analysis with RMSD.run(), which stores the results in the array RMSD.rmsd:
frame time (ps) RMSD (A)
This class uses Douglas Theobald’s fast QCP algorithm [Theobald2005] to calculate the RMSD.
New in version 0.7.7.
Setting up the RMSD analysis.
The RMSD will be computed between select and reference for all frames in the trajectory in universe.
Arguments : |
|
---|
Results are stored in this N×3 numpy.ndarray array, (frame, time (ps), RMSD (Å)).
Perform RMSD analysis on the trajectory.
A number of parameters can be changed from the defaults. The result is stored as the array RMSD.rmsd.
Keywords : |
|
---|
Save RMSD from RMSD.rmsd to text file filename.
If filename is not supplied then the default provided to the constructor is used.
The data are saved with numpy.savetxt().