Interface for reading and writing to a MDTraj HDF5 molecular dynamics trajectory file, whose format is described here.
This is a file-like object, that both reading or writing depending on the mode flag. It implements the context manager protocol, so you can also use it with the python ‘with’ statement.
The format is extremely flexible and high performance. It can hold a wide variety of information about a trajectory, including fields like the temperature and energies. Because it’s built on the fantastic HDF5 library, it’s easily extensible too.
Parameters: | filename : str
mode : {‘r, ‘w’}
force_overwrite : bool
compression : {‘zlib’, None}
|
---|
See also
Attributes
root | Direct access to the root group of the underlying Tables HDF5 file handle. |
title | User-defined title for the data represented in the file |
application | Suite of programs that created the file |
topology | Get the topology out from the file |
randomState | State of the creators internal random number generator at the start of the simulation |
forcefield | Description of the hamiltonian used. |
reference | A published reference that documents the program or parameters used to generate the data |
constraints | Constraints applied to the bond lengths |
Methods
__init__(filename[, mode, force_overwrite, ...]) | |
close() | Close the HDF5 file handle |
flush() | Write all buffered data in the to the disk file. |
read([n_frames, stride, atom_indices]) | Read one or more frames of data from the file |
read_as_traj([n_frames, stride, atom_indices]) | Read a trajectory from the HDF5 file |
seek(offset[, whence]) | Move to a new file position |
tell() | Current file position |
write(coordinates[, time, cell_lengths, ...]) | Write one or more frames of data to the file |
Attributes
application | Suite of programs that created the file |
constraints | Constraints applied to the bond lengths |
distance_unit | |
forcefield | Description of the hamiltonian used. |
randomState | State of the creators internal random number generator at the start of the simulation |
reference | A published reference that documents the program or parameters used to generate the data |
root | Direct access to the root group of the underlying Tables HDF5 file handle. |
title | User-defined title for the data represented in the file |
topology | Get the topology out from the file |
Direct access to the root group of the underlying Tables HDF5 file handle.
This can be used for random or specific access to the underlying arrays on disk
User-defined title for the data represented in the file
Suite of programs that created the file
Get the topology out from the file
Returns: | topology : mdtraj.Topology
|
---|
State of the creators internal random number generator at the start of the simulation
Description of the hamiltonian used. A short, human readable string, like AMBER99sbildn.
A published reference that documents the program or parameters used to generate the data
Constraints applied to the bond lengths
Returns: | constraints : {None, np.array, dtype=[(‘atom1’, ‘<i4’), (‘atom2’, ‘<i4’), (‘distance’, ‘<f4’)])}
|
---|
Read a trajectory from the HDF5 file
Parameters: | n_frames : {int, None}
stride : {int, None}
atom_indices : {int, None}
|
---|---|
Returns: | trajectory : Trajectory
|
Read one or more frames of data from the file
Parameters: | n_frames : {int, None}
stride : {int, None}
atom_indices : {int, None}
|
---|---|
Returns: | frames : namedtuple
|
Notes
If you’d like more flexible access to the data, that is available by using the pytables group directly, which is accessible via the root property on this class.
Write one or more frames of data to the file
This method saves data that is associated with one or more simulation frames. Note that all of the arguments can either be raw numpy arrays or unitted arrays (with simtk.unit.Quantity). If the arrays are unittted, a unit conversion will be automatically done from the supplied units into the proper units for saving on disk. You won’t have to worry about it.
Furthermore, if you wish to save a single frame of simulation data, you can do so naturally, for instance by supplying a 2d array for the coordinates and a single float for the time. This “shape deficiency” will be recognized, and handled appropriately.
Parameters: | coordinates : np.ndarray, shape=(n_frames, n_atoms, 3)
time : np.ndarray, shape=(n_frames,), optional
cell_lengths : np.ndarray, shape=(n_frames, 3), dtype=float32, optional
cell_angles : np.ndarray, shape=(n_frames, 3), dtype=float32, optional
velocities : np.ndarray, shape=(n_frames, n_atoms, 3), optional
kineticEnergy : np.ndarray, shape=(n_frames,), optional
potentialEnergy : np.ndarray, shape=(n_frames,), optional
temperature : np.ndarray, shape=(n_frames,), optional
alchemicalLambda : np.ndarray, shape=(n_frames,), optional
|
---|
Move to a new file position
Parameters: | offset : int
whence : {0, 1, 2}
|
---|
Current file position
Returns: | offset : int
|
---|
Close the HDF5 file handle
Write all buffered data in the to the disk file.