5.19. Base classes — MDAnalysis.coordinates.base

Derive other Timestep, Reader and Writer classes from the classes in this module. The derived classes must follow the Trajectory API in MDAnalysis.coordinates.__init__.

class MDAnalysis.coordinates.base.Timestep[source]
__init__(n_atoms, **kwargs)[source]

Create a Timestep, representing a frame of a trajectory

Arguments:
n_atoms

The total number of atoms this Timestep describes

Keywords:
positions

Whether this Timestep has position information [True]

velocities

Whether this Timestep has velocity information [False]

forces

Whether this Timestep has force information [False]

dt

The time difference between frames (ps). If time is set, then dt will be ignored.

time_offset

The starting time from which to calculate time (ps)

Changed in version 0.11.0: Added keywords for positions, velocities and forces Can add and remove position/velocity/force information by using the has_* attribute

classmethod from_coordinates(positions, velocities=None, forces=None, **kwargs)[source]

Create an instance of this Timestep, from coordinate data

New in version 0.11.0.

classmethod from_timestep(other, **kwargs)[source]

Create a copy of another Timestep, in the format of this Timestep

New in version 0.11.0.

_init_unitcell()[source]

Create custom datastructure for _unitcell.

n_atoms

A read only view of the number of atoms this Timestep has

Changed in version 0.11.0: Changed to read only property

time

The time in ps of this timestep

This is calculated as:

time = ts.data['time_offset'] + ts.time

Or, if the trajectory doesn’t provide time information:

time = ts.data['time_offset'] + ts.frame * ts.dt

New in version 0.11.0.

dt

The time difference in ps between timesteps

Note

This defaults to 1.0 ps in the absence of time data

New in version 0.11.0.

positions

A record of the positions of all atoms in this Timestep

Setting this attribute will add positions to the Timestep if they weren’t originally present.

Returns:

A numpy.ndarray of shape (n_atoms, 3) of position data for each atom

Raises:
MDAnalysis.exceptions.NoDataError

If the Timestep has no position data

Changed in version 0.11.0: Now can raise NoDataError when no position data present

velocities

A record of the velocities of all atoms in this Timestep

Setting this attribute will add velocities to the Timestep if they weren’t originally present.

Returns:

A numpy.ndarray of shape (n_atoms, 3) of velocity data for each atom

Raises:
MDAnalysis.exceptions.NoDataError

When the Timestep does not contain velocity information

New in version 0.11.0.

forces

A record of the forces of all atoms in this Timestep

Setting this attribute will add forces to the Timestep if they weren’t originally present.

Returns:

A numpy.ndarray of shape (n_atoms, 3) of force data for each atom

Raises:
MDAnalysis.NoDataError

When the Timestep does not contain force information

New in version 0.11.0.

has_positions

A boolean of whether this Timestep has position data

This can be changed to True or False to allocate space for or remove the data.

New in version 0.11.0.

has_velocities

A boolean of whether this Timestep has velocity data

This can be changed to True or False to allocate space for or remove the data.

New in version 0.11.0.

has_forces

A boolean of whether this Timestep has force data

This can be changed to True or False to allocate space for or remove the data.

New in version 0.11.0.

_pos

numpy.ndarray of dtype float32 of shape (n_atoms, 3) and internal FORTRAN order, holding the raw cartesian coordinates (in MDAnalysis units, i.e. Å).

Note

Normally one does not directly access _pos but uses the coordinates() method of an AtomGroup but sometimes it can be faster to directly use the raw coordinates. Any changes to this array are immediately reflected in atom positions. If the frame is written to a new trajectory then the coordinates are changed. If a new trajectory frame is loaded, then all contents of _pos are overwritten.

_velocities

numpy.ndarray of dtype float32. of shape (n_atoms, 3), holding the raw velocities (in MDAnalysis units, i.e. typically Å/ps).

Note

Normally velocities are accessed through the velocities or the velocities() method of an AtomGroup

_velocities only exists if the has_velocities flag is True

New in version 0.7.5.

_forces

numpy.ndarray of dtype float32. of shape (n_atoms, 3), holding the forces

_forces only exists if has_forces is True

New in version 0.11.0: Added as optional to Timestep

dimensions

unitcell dimensions (A, B, C, alpha, beta, gamma)

lengths a, b, c are in the MDAnalysis length unit (Å), and angles are in degrees.

Setting dimensions will populate the underlying native format description of box dimensions

triclinic_dimensions

The unitcell dimensions represented as triclinic vectors

Returns:A (3, 3) numpy.ndarray of unit cell vectors

For example:

>>> ts.dimensions
array([ 13.,  14.,  15.,  90.,  90.,  90.], dtype=float32)
>>> ts.triclinic_dimensions
array([[ 13.,   0.,   0.],
       [  0.,  14.,   0.],
       [  0.,   0.,  15.]], dtype=float32)

Setting the attribute also works:

>>> ts.triclinic_dimensions = [[15, 0, 0], [5, 15, 0], [5, 5, 15]]
>>> ts.dimensions
array([ 15.        ,  15.81138802,  16.58312416,  67.58049774,
        72.45159912,  71.56504822], dtype=float32)

New in version 0.11.0.

volume

volume of the unitcell

__getitem__(atoms)[source]

Get a selection of coordinates

ts[i]

return coordinates for the i’th atom (0-based)

ts[start:stop:skip]

return an array of coordinates, where start, stop and skip correspond to atom indices, MDAnalysis.core.AtomGroup.Atom.index (0-based)
__eq__(other)[source]

Compare with another Timestep

New in version 0.11.0.

__iter__()[source]

Iterate over coordinates

for x in ts

iterate of the coordinates, atom by atom
copy()[source]

Make an independent (“deep”) copy of the whole Timestep.

copy_slice(sel)[source]

Make a new Timestep containing a subset of the original Timestep.

ts.copy_slice(slice(start, stop, skip)) ts.copy_slice([list of indices])

Returns:A Timestep object of the same type containing all header information and all atom information relevant to the selection.

Note

The selection must be a 0 based slice or array of the atom indices in this Timestep

New in version 0.8.

Changed in version 0.11.0: Reworked to follow new Timestep API. Now will strictly only copy official attributes of the Timestep.

class MDAnalysis.coordinates.base.IObase[source]

Base class bundling common functionality for trajectory I/O.

Changed in version 0.8: Added context manager protocol.

close()[source]

Close the trajectory file.

convert_forces_from_native(force, inplace=True)[source]

In-place conversion of forces array force from native units to base units.

By default, the input force is modified in place and also returned.

New in version 0.7.7.

convert_forces_to_native(force, inplace=True)[source]

In-place conversion of force array force from base units to native units.

By default, the input force is modified in place and also returned.

New in version 0.7.7.

convert_pos_from_native(x, inplace=True)[source]

In-place conversion of coordinate array x from native units to base units.

By default, the input x is modified in place and also returned.

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)[source]

Conversion of coordinate array x from base units to native units.

By default, the input x is modified in place and also returned.

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)[source]

Convert time t from native units to base units.

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.)

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)[source]

Convert time t from base units to native units.

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)[source]

In-place conversion of velocities array v from native units to base units.

By default, the input v is modified in place and also returned.

New in version 0.7.5.

convert_velocities_to_native(v, inplace=True)[source]

In-place conversion of coordinate array v from base units to native units.

By default, the input v is modified in place and also returned.

New in version 0.7.5.

format = None

override to define trajectory format of the reader/writer (DCD, XTC, ...)

units = {'velocity': None, 'length': None, 'time': None}

dict with units of of time and length (and velocity, force, ... for formats that support it)

class MDAnalysis.coordinates.base.Reader(filename, convert_units=None, **kwargs)[source]

Base class for trajectory readers that extends ProtoReader with a __del__() method.

New Readers should subclass Reader and properly implement a close() method, to ensure proper release of resources (mainly file handles). Readers that are inherently safe in this regard should subclass ProtoReader instead.

See the Trajectory API definition in MDAnalysis.coordinates.__init__ for the required attributes and methods.

See also

ProtoReader

Changed in version 0.11.0: Most of the base Reader class definitions were offloaded to ProtoReader so as to allow the subclassing of Readers without a __del__() method. Created init method to create common functionality, all Reader subclasses must now super() through this class. Added attribute _ts_kwargs, which is created in init. Provides kwargs to be passed to Timestep

class MDAnalysis.coordinates.base.ChainReader(filenames, **kwargs)[source]

Reader that concatenates multiple trajectories on the fly.

Known issues

  • Trajectory API attributes exist but most of them only reflect the first trajectory in the list; ChainReader.n_frames, ChainReader.n_atoms, and ChainReader.fixed are properly set, though
  • slicing not implemented
  • time will not necessarily return the true time but just number of frames times a provided time between frames (from the keyword delta)

Changed in version 0.11.0: Frames now 0-based instead of 1-based

Set up the chain reader.

Arguments:
filenames

file name or list of file names; the reader will open all file names and provide frames in the order of trajectories from the list. Each trajectory must contain the same number of atoms in the same order (i.e. they all must belong to the same topology). The trajectory format is deduced from the extension of filename.

Extension: filenames are either single filename or list of file names in either plain file names format or (filename,format) tuple combination

skip

skip step (also passed on to the individual trajectory readers); must be same for all trajectories

delta

The time between frames in MDAnalysis time units if no other information is available. If this is not set then any call to time will raise a ValueError.

kwargs

all other keyword arguments are passed on to each trajectory reader unchanged

Changed in version 0.8: The delta keyword was added.

active_reader

Reader instance from which frames are being read.

frame

Cumulative frame number of the current time step.

rewind()[source]

Set current frame to the beginning.

time

Cumulative time of the current frame in MDAnalysis time units (typically ps).

class MDAnalysis.coordinates.base.Writer[source]

Base class for trajectory writers.

See Trajectory API definition in MDAnalysis.coordinates.__init__ for the required attributes and methods.

convert_dimensions_to_unitcell(ts)[source]

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.

has_valid_coordinates(criteria, x)[source]

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 <= max
Arguments:
criteria

dictionary containing the max and min values in native units

x

np.ndarray of (x, y, z) coordinates of atoms selected to be written out.

Returns:

boolean

write(obj)[source]

Write current timestep, using the supplied obj.

The argument should be a AtomGroup or a Universe or a Timestep instance.

Note

The size of the obj must be the same as the number of atom provided when setting up the trajectory.