6.3. DCD trajectory I/O — MDAnalysis.coordinates.DCD
¶
Classes to read and write DCD binary trajectories, the format used by CHARMM, NAMD, and also LAMMPS. Trajectories can be read regardless of system-endianness as this is auto-detected.
Generally, DCD trajectories produced by any code can be read (with the
DCDReader
) although there can be issues with the unitcell (simulation
box) representation (see DCDReader.dimensions
). DCDs can also be
written but the DCDWriter
follows recent NAMD/VMD convention for the
unitcell but still writes AKMA time. Reading and writing these trajectories
within MDAnalysis will work seamlessly but if you process those trajectories
with other tools you might need to watch out that time and unitcell dimensions
are correctly interpreted.
See also
MDAnalysis.coordinates.LAMMPS
- module provides a more flexible DCD reader/writer.
6.3.1. Classes¶
-
class
MDAnalysis.coordinates.DCD.
DCDReader
(filename, convert_units=True, dt=None, **kwargs)[source]¶ Reader for the DCD format.
DCD is used by NAMD, CHARMM and LAMMPS as the default trajectory format. The DCD file format is not well defined. In particular, NAMD and CHARMM use it differently. Currently, MDAnalysis tries to guess the correct format for the unitcell representation but it can be wrong. Check the unitcell dimensions, especially for triclinic unitcells (see Issue 187). DCD trajectories produced by CHARMM and NAMD( >2.5) record time in AKMA units. If other units have been recorded (e.g., ps) then employ the configurable :class:MDAnalysis.coordinates.LAMMPS.DCDReader and set the time unit as a optional argument. You can find a list of units used in the DCD formats on the MDAnalysis wiki.
MDAnalysis always uses
(*A*, *B*, *C*, *alpha*, *beta*, *gamma*)
to represent the unit cell. Lengths A, B, C are in the MDAnalysis length unit (Å), and angles are in degrees.The ordering of the angles in the unitcell is the same as in recent versions of VMD’s DCDplugin (2013), namely the X-PLOR DCD format: The original unitcell is read as
[A, gamma, B, beta, alpha, C]
from the DCD file. If any of these values are < 0 or if any of the angles are > 180 degrees then it is assumed it is a new-style CHARMM unitcell (at least since c36b2) in which box vectors were recorded.Warning
The DCD format is not well defined. Check your unit cell dimensions carefully, especially when using triclinic boxes. Different software packages implement different conventions and MDAnalysis is currently implementing the newer NAMD/VMD convention and tries to guess the new CHARMM one. Old CHARMM trajectories might give wrong unitcell values. For more details see Issue 187.
Parameters: Changed in version 0.17.0: Changed to use libdcd.pyx library and removed the correl function
-
OtherWriter
(filename, **kwargs)¶ Returns a writer appropriate for filename.
Sets the default keywords start, step and dt (if available). n_atoms is always set from
Reader.n_atoms
.See also
Reader.Writer()
-
add_auxiliary
(auxname, auxdata, format=None, **kwargs)¶ Add auxiliary data to be read alongside trajectory.
Auxiliary data may be any data timeseries from the trajectory additional to that read in by the trajectory reader. auxdata can be an
AuxReader
instance, or the data itself as e.g. a filename; in the latter case an appropriateAuxReader
is guessed from the data/file format. An appropriate format may also be directly provided as a key word argument.On adding, the AuxReader is initially matched to the current timestep of the trajectory, and will be updated when the trajectory timestep changes (through a call to
next()
or jumping timesteps withtrajectory[i]
).The representative value(s) of the auxiliary data for each timestep (as calculated by the
AuxReader
) are stored in the current timestep in thets.aux
namespace under auxname; e.g. to add additional pull force data stored in pull-force.xvg:u = MDAnalysis.Universe(PDB, XTC) u.trajectory.add_auxiliary('pull', 'pull-force.xvg')
The representative value for the current timestep may then be accessed as
u.trajectory.ts.aux.pull
oru.trajectory.ts.aux['pull']
.See also
Note
Auxiliary data is assumed to be time-ordered, with no duplicates. See the Auxiliary API.
-
add_transformations
(*transformations)¶ Add all transformations to be applied to the trajectory.
This function take as list of transformations as an argument. These transformations are functions that will be called by the Reader and given a
Timestep
object as argument, which will be transformed and returned to the Reader. The transformations can be part of thetransformations
module, or created by the user, and are stored as a list transformations. This list can only be modified once, and further calls of this function will raise an exception.u = MDAnalysis.Universe(topology, coordinates) workflow = [some_transform, another_transform, this_transform] u.trajectory.add_transformations(*workflow)
The transformations are applied in the order given in the list transformations, i.e., the first transformation is the first or innermost one to be applied to the
Timestep
. The example above would be equivalent tofor ts in u.trajectory: ts = this_transform(another_transform(some_transform(ts)))
Parameters: transform_list (list) – list of all the transformations that will be applied to the coordinates in the order given in the list See also
-
aux_list
¶ Lists the names of added auxiliary data.
-
check_slice_indices
(start, stop, step)¶ Check frame indices are valid and clip to fit trajectory.
The usage follows standard Python conventions for
range()
but see the warning below.Parameters: - start (int or None) – Starting frame index (inclusive).
None
corresponds to the default of 0, i.e., the initial frame. - stop (int or None) – Last frame index (exclusive).
None
corresponds to the default of n_frames, i.e., it includes the last frame of the trajectory. - step (int or None) – step size of the slice,
None
corresponds to the default of 1, i.e, include every frame in the range start, stop.
Returns: start, stop, step – Integers representing the slice
Return type: Warning
The returned values start, stop and step give the expected result when passed in
range()
but gives unexpected behavior when passed in aslice
whenstop=None
andstep=-1
This can be a problem for downstream processing of the output from this method. For example, slicing of trajectories is implemented by passing the values returned by
check_slice_indices()
torange()
range(start, stop, step)
and using them as the indices to randomly seek to. On the other hand, in
MDAnalysis.analysis.base.AnalysisBase
the values returned bycheck_slice_indices()
are used to splice the trajectory by creating aslice
instanceslice(start, stop, step)
This creates a discrepancy because these two lines are not equivalent:
range(10, -1, -1) # [10, 9, 8, 7, 6, 5, 4, 3, 2, 1, 0] range(10)[slice(10, -1, -1)] # []
- start (int or None) – Starting frame index (inclusive).
-
convert_forces_from_native
(force, inplace=True)¶ Conversion of forces array force from native to base units
Parameters: - force (array_like) – Forces to transform
- inplace (bool (optional)) – Whether to modify the array inplace, overwriting previous data
Note
By default, the input force is modified in place and also returned. In-place operations improve performance because allocating new arrays is avoided.
New in version 0.7.7.
-
convert_forces_to_native
(force, inplace=True)¶ Conversion of force array force from base to native units.
Parameters: - force (array_like) – Forces to transform
- inplace (bool (optional)) – Whether to modify the array inplace, overwriting previous data
Note
By default, the input force is modified in place and also returned. In-place operations improve performance because allocating new arrays is avoided.
New in version 0.7.7.
-
convert_pos_from_native
(x, inplace=True)¶ Conversion of coordinate array x from native units to base units.
Parameters: - x (array_like) – Positions to transform
- inplace (bool (optional)) – Whether to modify the array inplace, overwriting previous data
Note
By default, the input x is modified in place and also returned. In-place operations improve performance because allocating new arrays is avoided.
Changed in version 0.7.5: Keyword inplace can be set to
False
so that a modified copy is returned unless no conversion takes place, in which case the reference to the unmodified x is returned.
-
convert_pos_to_native
(x, inplace=True)¶ Conversion of coordinate array x from base units to native units.
Parameters: - x (array_like) – Positions to transform
- inplace (bool (optional)) – Whether to modify the array inplace, overwriting previous data
Note
By default, the input x is modified in place and also returned. In-place operations improve performance because allocating new arrays is avoided.
Changed in version 0.7.5: Keyword inplace can be set to
False
so that a modified copy is returned unless no conversion takes place, in which case the reference to the unmodified x is returned.
-
convert_time_from_native
(t, inplace=True)¶ Convert time t from native units to base units.
Parameters: - t (array_like) – Time values to transform
- inplace (bool (optional)) – Whether to modify the array inplace, overwriting previous data
Note
By default, the input t is modified in place and also returned (although note that scalar values t are passed by value in Python and hence an in-place modification has no effect on the caller.) In-place operations improve performance because allocating new arrays is avoided.
Changed in version 0.7.5: Keyword inplace can be set to
False
so that a modified copy is returned unless no conversion takes place, in which case the reference to the unmodified x is returned.
-
convert_time_to_native
(t, inplace=True)¶ Convert time t from base units to native units.
Parameters: - t (array_like) – Time values to transform
- inplace (bool, optional) – Whether to modify the array inplace, overwriting previous data
Note
By default, the input t is modified in place and also returned. (Also note that scalar values t are passed by value in Python and hence an in-place modification has no effect on the caller.)
Changed in version 0.7.5: Keyword inplace can be set to
False
so that a modified copy is returned unless no conversion takes place, in which case the reference to the unmodified x is returned.
-
convert_velocities_from_native
(v, inplace=True)¶ Conversion of velocities array v from native to base units
Parameters: - v (array_like) – Velocities to transform
- inplace (bool (optional)) – Whether to modify the array inplace, overwriting previous data
Note
By default, the input v is modified in place and also returned. In-place operations improve performance because allocating new arrays is avoided.
New in version 0.7.5.
-
convert_velocities_to_native
(v, inplace=True)¶ Conversion of coordinate array v from base to native units
Parameters: - v (array_like) – Velocities to transform
- inplace (bool (optional)) – Whether to modify the array inplace, overwriting previous data
Note
By default, the input v is modified in place and also returned. In-place operations improve performance because allocating new arrays is avoided.
New in version 0.7.5.
-
copy
()¶ Return independent copy of this Reader.
New Reader will have its own file handle and can seek/iterate independently of the original.
Will also copy the current state of the Timestep held in the original Reader
-
dimensions
¶ unitcell dimensions (A, B, C, alpha, beta, gamma)
-
dt
¶ timestep between frames
-
frame
¶ Frame number of the current time step.
This is a simple short cut to
Timestep.frame
.
-
get_aux_attribute
(auxname, attrname)¶ Get the value of attrname from the auxiliary auxname
Parameters: See also
-
get_aux_descriptions
(auxnames=None)¶ Get descriptions to allow reloading the specified auxiliaries.
If no auxnames are provided, defaults to the full list of added auxiliaries.
Passing the resultant description to
add_auxiliary()
will allow recreation of the auxiliary. e.g., to duplicate all auxiliaries into a second trajectory:descriptions = trajectory_1.get_aux_descriptions() for aux in descriptions: trajectory_2.add_auxiliary(**aux)
Returns: List of dictionaries of the args/kwargs describing each auxiliary. Return type: list
-
iter_as_aux
(auxname)¶ Iterate through timesteps for which there is at least one assigned step from the auxiliary auxname within the cutoff specified in auxname.
See also
-
iter_auxiliary
(auxname, start=None, stop=None, step=None, selected=None)¶ Iterate through the auxiliary auxname independently of the trajectory.
Will iterate over the specified steps of the auxiliary (defaults to all steps). Allows to access all values in an auxiliary, including those out of the time range of the trajectory, without having to also iterate through the trajectory.
After interation, the auxiliary will be repositioned at the current step.
Parameters: - auxname (str) – Name of the auxiliary to iterate over.
- stop, step) ((start,) – Options for iterating over a slice of the auxiliary.
- selected (lst | ndarray, optional) – List of steps to iterate over.
Yields: AuxStep
objectSee also
-
n_frames
¶ number of frames in trajectory
-
next
()¶ Forward one step to next frame.
-
next_as_aux
(auxname)¶ Move to the next timestep for which there is at least one step from the auxiliary auxname within the cutoff specified in auxname.
This allows progression through the trajectory without encountering
NaN
representative values (unless these are specifically part of the auxiliary data).If the auxiliary cutoff is not set, where auxiliary steps are less frequent (
auxiliary.dt > trajectory.dt
), this allows progression at the auxiliary pace (rounded to nearest timestep); while if the auxiliary steps are more frequent, this will work the same as callingnext()
.See the Auxiliary API.
See also
-
static
parse_n_atoms
(filename, **kwargs)[source]¶ Read the coordinate file and deduce the number of atoms
Returns: n_atoms – the number of atoms in the coordinate file Return type: int Raises: NotImplementedError
– when the number of atoms can’t be deduced
-
rename_aux
(auxname, new)¶ Change the name of the auxiliary auxname to new.
Provided there is not already an auxiliary named new, the auxiliary name will be changed in ts.aux namespace, the trajectory’s list of added auxiliaries, and in the auxiliary reader itself.
Parameters: Raises: ValueError
– If the name new is already in use by an existing auxiliary.
-
rewind
()¶ Position at beginning of trajectory
-
set_aux_attribute
(auxname, attrname, new)¶ Set the value of attrname in the auxiliary auxname.
Parameters: See also
-
time
¶ Time of the current frame in MDAnalysis time units (typically ps).
This is either read straight from the Timestep, or calculated as time =
Timestep.frame
*Timestep.dt
-
timeseries
(asel=None, start=None, stop=None, step=None, order='afc')[source]¶ Return a subset of coordinate data for an AtomGroup
Parameters: - asel (
AtomGroup
) – TheAtomGroup
to read the coordinates from. Defaults to None, in which case the full set of coordinate data is returned. - start (int (optional)) – Begin reading the trajectory at frame index start (where 0 is the
index of the first frame in the trajectory); the default
None
starts at the beginning. - stop (int (optional)) – End reading the trajectory at frame index stop-1, i.e, stop is
excluded. The trajectory is read to the end with the default
None
. - step (int (optional)) – Step size for reading; the default
None
is equivalent to 1 and means to read every frame. - order (str (optional)) – the order/shape of the return data array, corresponding to (a)tom, (f)rame, (c)oordinates all six combinations of ‘a’, ‘f’, ‘c’ are allowed ie “fac” - return array where the shape is (frame, number of atoms, coordinates)
Changed in version 1.0.0: skip and format keywords have been removed.
- asel (
-
totaltime
¶ Total length of the trajectory
The time is calculated as
(n_frames - 1) * dt
, i.e., we assume that the first frame no time as elapsed. Thus, a trajectory with two frames will be considered to have a length of a single time step dt and a “trajectory” with a single frame will be reported as length 0.
-
transformations
¶ Returns the list of transformations
-
-
class
MDAnalysis.coordinates.DCD.
DCDWriter
(filename, n_atoms, convert_units=True, step=1, dt=1, remarks='', nsavc=1, istart=0, **kwargs)[source]¶ DCD Writer class
The writer follows recent NAMD/VMD convention for the unitcell (box lengths in Å and angle-cosines,
[A, cos(gamma), B, cos(beta), cos(alpha), C]
) and writes positions in Å and time in AKMA time units.Note
When writing out timesteps without
dimensions
(i.e. setNone
) theDCDWriter
will write out a zeroed unitcell (i.e.[0, 0, 0, 0, 0, 0]
). As this behaviour is poorly defined, it may not match the expectations of other software.Parameters: - filename (str) – filename of trajectory
- n_atoms (int) – number of atoms to be written
- convert_units (bool (optional)) – convert from MDAnalysis units to format specific units
- step (int (optional)) – number of steps between frames to be written
- dt (float (optional)) – time between two frames. If
None
guess from first writtenTimeStep
- remarks (str (optional)) – remarks to be stored in DCD. Shouldn’t be more then 240 characters
- nsavc (int (optional)) – DCD files can also store the number of integrator time steps that correspond to the interval between two frames as nsavc (i.e., every how many MD steps is a frame saved to the DCD). By default, this number is just set to one and this should be sufficient for almost all cases but if required, nsavc can be changed.
- istart (int (optional)) – starting frame number in integrator timesteps. CHARMM defaults to
nsavc, i.e., start at frame number 1 = istart / nsavc. The value
None
will set istart to nsavc (the CHARMM default). The MDAnalysis default is 0 so that the frame number and time of the first frame is 0. - **kwargs (dict) – General writer arguments
-
convert_dimensions_to_unitcell
(ts, inplace=True)¶ Read dimensions from timestep ts and return appropriate unitcell.
The default is to return
[A,B,C,alpha,beta,gamma]
; if this is not appropriate then this method has to be overriden.
-
convert_forces_from_native
(force, inplace=True)¶ Conversion of forces array force from native to base units
Parameters: - force (array_like) – Forces to transform
- inplace (bool (optional)) – Whether to modify the array inplace, overwriting previous data
Note
By default, the input force is modified in place and also returned. In-place operations improve performance because allocating new arrays is avoided.
New in version 0.7.7.
-
convert_forces_to_native
(force, inplace=True)¶ Conversion of force array force from base to native units.
Parameters: - force (array_like) – Forces to transform
- inplace (bool (optional)) – Whether to modify the array inplace, overwriting previous data
Note
By default, the input force is modified in place and also returned. In-place operations improve performance because allocating new arrays is avoided.
New in version 0.7.7.
-
convert_pos_from_native
(x, inplace=True)¶ Conversion of coordinate array x from native units to base units.
Parameters: - x (array_like) – Positions to transform
- inplace (bool (optional)) – Whether to modify the array inplace, overwriting previous data
Note
By default, the input x is modified in place and also returned. In-place operations improve performance because allocating new arrays is avoided.
Changed in version 0.7.5: Keyword inplace can be set to
False
so that a modified copy is returned unless no conversion takes place, in which case the reference to the unmodified x is returned.
-
convert_pos_to_native
(x, inplace=True)¶ Conversion of coordinate array x from base units to native units.
Parameters: - x (array_like) – Positions to transform
- inplace (bool (optional)) – Whether to modify the array inplace, overwriting previous data
Note
By default, the input x is modified in place and also returned. In-place operations improve performance because allocating new arrays is avoided.
Changed in version 0.7.5: Keyword inplace can be set to
False
so that a modified copy is returned unless no conversion takes place, in which case the reference to the unmodified x is returned.
-
convert_time_from_native
(t, inplace=True)¶ Convert time t from native units to base units.
Parameters: - t (array_like) – Time values to transform
- inplace (bool (optional)) – Whether to modify the array inplace, overwriting previous data
Note
By default, the input t is modified in place and also returned (although note that scalar values t are passed by value in Python and hence an in-place modification has no effect on the caller.) In-place operations improve performance because allocating new arrays is avoided.
Changed in version 0.7.5: Keyword inplace can be set to
False
so that a modified copy is returned unless no conversion takes place, in which case the reference to the unmodified x is returned.
-
convert_time_to_native
(t, inplace=True)¶ Convert time t from base units to native units.
Parameters: - t (array_like) – Time values to transform
- inplace (bool, optional) – Whether to modify the array inplace, overwriting previous data
Note
By default, the input t is modified in place and also returned. (Also note that scalar values t are passed by value in Python and hence an in-place modification has no effect on the caller.)
Changed in version 0.7.5: Keyword inplace can be set to
False
so that a modified copy is returned unless no conversion takes place, in which case the reference to the unmodified x is returned.
-
convert_velocities_from_native
(v, inplace=True)¶ Conversion of velocities array v from native to base units
Parameters: - v (array_like) – Velocities to transform
- inplace (bool (optional)) – Whether to modify the array inplace, overwriting previous data
Note
By default, the input v is modified in place and also returned. In-place operations improve performance because allocating new arrays is avoided.
New in version 0.7.5.
-
convert_velocities_to_native
(v, inplace=True)¶ Conversion of coordinate array v from base to native units
Parameters: - v (array_like) – Velocities to transform
- inplace (bool (optional)) – Whether to modify the array inplace, overwriting previous data
Note
By default, the input v is modified in place and also returned. In-place operations improve performance because allocating new arrays is avoided.
New in version 0.7.5.
-
has_valid_coordinates
(criteria, x)¶ Returns
True
if all values are within limit values of their formats.Due to rounding, the test is asymmetric (and min is supposed to be negative):
min < x <= maxParameters: - criteria (dict) – dictionary containing the max and min values in native units
- x (numpy.ndarray) –
(x, y, z)
coordinates of atoms selected to be written out
Returns: Return type:
-
write
(obj)¶ Write current timestep, using the supplied obj.
Parameters: obj ( AtomGroup
orUniverse
) – write coordinate information associate with objNote
The size of the obj must be the same as the number of atoms provided when setting up the trajectory.
Changed in version 2.0.0: Deprecated support for Timestep argument to write has now been removed. Use AtomGroup or Universe as an input instead.