libxdrfile2, a derivative of the Gromacs libxdrfile library, provides an interface to some high-level functions for XTC/TRR trajectory handling. Only functions required for reading and processing whole trajectories are exposed at the moment; low-level routines to read individual numbers are not provided. In addition, libxdrfile2 exposes functions to allow fast frame indexing and XDR file seeking.
The functions querying the numbers of atoms in a trajectory frame (read_xtc_natoms() and read_trr_natoms()) open a file themselves and only require the file name.
All other functions operate on a XDRFILE object, which is a special file handle for xdr files. Any xdr-based trajectory file (XTC or TRR format) always has to be opened with xdrfile_open(). When done, close the trajectory with xdrfile_close().
The functions fill or read existing arrays of coordinates; they never allocate these arrays themselves. Hence they need to be setup outside libxdrfile2 as numpy arrays. The exception to these are the indexing ones functions that take care of array allocation and transference to a garbage-collectable memory object.
Changed in version 0.8.0: libxdrfile2 is now used instead of libxdrfile. libxdrfile2 is based on libxdrfile but has xdr seeking and indexing capabilities. Unlike libxdrfile before it, libxdrfile2 is distributed under the GNU GENERAL PUBLIC LICENSE, version 2 (or higher).
In the example we read coordinate frames from an existing XTC trajectory:
import numpy as np
from libxdrfile2 import xdrfile_open, xdrfile_close, read_xtc_natoms, read_xtc, DIM, exdrOK
xtc = 'md.xtc'
# get number of atoms
natoms = read_xtc_natoms(xtc)
# allocate coordinate array of the right size and type
# (the type float32 is crucial to match the underlying C-code!!)
x = np.zeros((natoms, DIM), dtype=np.float32)
# allocate unit cell box
box = np.zeros((DIM, DIM), dtype=np.float32)
# open file
XTC = xdrfile_open(xtc, 'r')
# loop through file until return status signifies end or a problem
# (it should become exdrENDOFFILE on the last iteration)
status = exdrOK
while status == exdrOK:
status,step,time,prec = read_xtc(XTC, box, x)
# do something with x
centre = x.mean(axis=0)
print 'Centre of geometry at %(time)g ps: %(centre)r' % vars()
# finally close file
xdrfile_close(XTC)
Note that only the contents of the coordinate and unitcell arrays x and box change.
The module defines a number of constants such as DIM or the Status symbols.
The number of cartesian dimensions for which the underlying C-code was compiled; this is most certainly 3.
A number of symbols are exported; they all start with the letters exdr. Important ones are:
Success of xdr file read/write operation.
xdr file is closed
end of file was reached (response of read_xtc() and read_trr() after the last read frame)
xdrfile_open() cannot find the requested file
Other symbols that are used internally are:
header
string
double precision floating point number
integer
floating point number
unsigned integer
compressed 3D coordinates
magic number
not enough memory to allocate space for a XDR data structure.
Two low-level functions are used to obtain a XDRFILE object (a file handle) to access xdr files such as XTC or TRR trajectories.
The XTC trajectory format is a lossy compression format that only stores coordinates. Compression level is determined by the precision argument to the write_xtc() function. Coordinates (Gromacs uses nm natively) are multiplied by precision and truncated to the integer part. A typical value is 1000.0, which gives an accuracy of 1/100 of an Angstroem.
The advantage of XTC over TRR is its significantly reduced size.
Read the number of atoms natoms from a xtc file fn.
Arguments: |
|
---|---|
Raises: | IOError if the supplied filed is not a XTC or if it is not readable. |
Read through the whole trajectory headers to obtain the total number of frames. The process is speeded up by reading frame headers for the amount of data in the frame, and then skipping directly to the next header. An array of frame offsets is also returned, which can later be used to seek direcly to arbitrary frames in the trajectory.
Arguments: |
|
---|---|
Returns: |
|
Raises: | IOError if the supplied filed is not a XTC or if it is not readable. |
Read the next frame from the opened xtc trajectory into x.
Arguments: |
|
---|---|
Returns: |
|
Write the next frame x to the opened xtc trajectory.
Arguments: |
|
---|---|
Returns: | status, integer status (0 = OK), see the libxdrfile2.exdr* constants under Status symbols for other values) |
TRR is the Gromacs native full-feature trajectory storage format. It can contain position coordinates, velocities and forces, and the lambda value for free energy perturbation calculations. Velocities and forces are optional in the sense that they can be all zero.
Read the number of atoms natoms from a trr file fn.
Arguments: |
|
---|---|
Raises: | IOError if the supplied filed is not a TRR or if it is not readable. |
Read through the whole trajectory headers to obtain the total number of frames. The process is speeded up by reading frame headers for the amount of data in the frame, and then skipping directly to the next header. An array of frame offsets is also returned, which can later be used to seek direcly to arbitrary frames in the trajectory.
Arguments: |
|
---|---|
Returns: |
|
Raises: | IOError if the supplied filed is not a TRR or if it is not readable. |
Read the next frame from the opened trr trajectory into x, v, and f.
Arguments: |
|
---|---|
Returns: |
|
Write the next frame to the opened trr trajectory.
Arguments: |
|
---|
Changed in version 0.8.0: either one of x, v, or f can now be set as a natom,0-DIM numpy array((natom, 0),dtype=nump.float32). This will cause the corresponding property to be skipped when writing to file.
Returns: | status, integer status (0 = OK), see the libxdrfile2.exdr* constants under Status symbols for other values) |
---|