# -*- Mode: python; tab-width: 4; indent-tabs-mode:nil; coding:utf-8 -*-
# vim: tabstop=4 expandtab shiftwidth=4 softtabstop=4
#
# MDAnalysis --- https://www.mdanalysis.org
# Copyright (c) 2006-2017 The MDAnalysis Development Team and contributors
# (see the file AUTHORS for the full list of names)
#
# Released under the GNU Public Licence, v2 or any higher version
#
# Please cite your use of MDAnalysis in published work:
#
# R. J. Gowers, M. Linke, J. Barnoud, T. J. E. Reddy, M. N. Melo, S. L. Seyler,
# D. L. Dotson, J. Domanski, S. Buchoux, I. M. Kenney, and O. Beckstein.
# MDAnalysis: A Python package for the rapid analysis of molecular dynamics
# simulations. In S. Benthall and S. Rostrup editors, Proceedings of the 15th
# Python in Science Conference, pages 102-109, Austin, TX, 2016. SciPy.
#
# N. Michaud-Agrawal, E. J. Denning, T. B. Woolf, and O. Beckstein.
# MDAnalysis: A Toolkit for the Analysis of Molecular Dynamics Simulations.
# J. Comput. Chem. 32 (2011), 2319--2327, doi:10.1002/jcc.21787
#
# TPR parser and tpr support module
# Copyright (c) 2011 Zhuyi Xue
# Released under the GNU Public Licence, v2
"""Gromacs portable run input TPR format parser
============================================
Classes
-------
.. autoclass:: TPRParser
:members:
See Also
--------
:mod:`MDAnalysis.topology.tpr`
Development notes
-----------------
The TPR reader is a pure-python implementation of a basic TPR
parser. Currently the following sections of the topology are parsed:
* Atoms: number, name, type, resname, resid, segid, mass, charge,
[residue, segment, radius, bfactor, resnum, moltype]
* Bonds
* Angels
* Dihedrals
* Impropers
This tpr parser is written according to the following files
- :file:`{gromacs_dir}/src/kernel/gmxdump.c`
- :file:`{gromacs_dir}/src/gmxlib/tpxio.c` (the most important one)
- :file:`{gromacs_dir}/src/gmxlib/gmxfio_rw.c`
- :file:`{gromacs_dir}/src/gmxlib/gmxfio_xdr.c`
- :file:`{gromacs_dir}/include/gmxfiofio.h`
or their equivalent in more recent versions of Gromacs.
The function :func:`read_tpxheader` is based on the
`TPRReaderDevelopment`_ notes. Functions with names starting with
``read_`` or ``do_`` are trying to be similar to those in
:file:`gmxdump.c` or :file:`tpxio.c`, those with ``extract_`` are new.
Wherever ``fver_err(fver)`` is used, it means the tpx version problem
has not been solved. Versions prior to Gromacs 4.0.x are not supported.
.. Links
.. _Gromacs: http://www.gromacs.org
.. _`Gromacs manual`: http://manual.gromacs.org/documentation/5.1/manual-5.1.pdf
.. _TPR file: http://manual.gromacs.org/current/online/tpr.html
.. _`Issue Tracker`: https://github.com/MDAnalysis/mdanalysis/issues
.. _`Issue 2`: https://github.com/MDAnalysis/mdanalysis/issues/2
.. _`Issue 463`: https://github.com/MDAnalysis/mdanalysis/pull/463
.. _TPRReaderDevelopment: https://github.com/MDAnalysis/mdanalysis/wiki/TPRReaderDevelopment
"""
from __future__ import absolute_import
__author__ = "Zhuyi Xue"
__copyright__ = "GNU Public Licence, v2"
import xdrlib
from . import guessers
from ..lib.util import openany
from .tpr import utils as tpr_utils
from .base import TopologyReaderBase
from ..core.topologyattrs import Resnums
import logging
logger = logging.getLogger("MDAnalysis.topology.TPRparser")
[docs]class TPRParser(TopologyReaderBase):
"""Read topology information from a Gromacs_ `TPR file`_.
.. _Gromacs: http://www.gromacs.org
.. _TPR file: http://manual.gromacs.org/current/online/tpr.html
"""
format = 'TPR'
[docs] def parse(self, **kwargs):
"""Parse a Gromacs TPR file into a MDAnalysis internal topology structure.
Returns
-------
structure : dict
"""
with openany(self.filename, mode='rb') as infile:
tprf = infile.read()
data = xdrlib.Unpacker(tprf)
try:
th = tpr_utils.read_tpxheader(data) # tpxheader
except EOFError:
msg = "{0}: Invalid tpr file or cannot be recognized".format(self.filename)
logger.critical(msg)
raise IOError(msg)
self._log_header(th)
V = th.fver # since it's used very often
state_ngtc = th.ngtc # done init_state() in src/gmxlib/tpxio.c
if th.bBox:
tpr_utils.extract_box_info(data, V)
if state_ngtc > 0 and V >= 28:
if V < 69: # redundancy due to different versions
tpr_utils.ndo_real(data, state_ngtc)
tpr_utils.ndo_real(data, state_ngtc) # relevant to Berendsen tcoupl_lambda
if V < 26:
tpr_utils.fileVersion_err(V)
if th.bTop:
tpr_top = tpr_utils.do_mtop(data, V)
else:
msg = "{0}: No topology found in tpr file".format(self.filename)
logger.critical(msg)
raise IOError(msg)
tpr_top.add_TopologyAttr(Resnums(tpr_top.resids.values.copy()))
return tpr_top
# THE FOLLOWING CODE IS WORKING FOR TPX VERSION 58, BUT SINCE THESE INFO IS
# NOT INTERESTED, SO IT'S NOT COVERED IN ALL VERSIONS. PARSING STOPS HERE.
# if th.bX:
# ndo_rvec(data, th.natoms)
# if th.bV:
# ndo_rvec(data, th.natoms)
# if th.bF:
# ndo_rvec(data, th.natoms)
# not useful at the moment
# ePBC = -1;
# bPeriodicMols = False
# if th.bIr:
# # update
# data.unpack_int() # ePBC
# data.unpack_bool() # bPeriodicMols
# # 17 < 23. and ir (ir is from the c code, seems not apply here
# if th.fgen < setting.tpx_generation:
# # a crazily long (670 lines) function in c, slightly better here
# # (240 lines), so put it in setting.py
# utils.do_inputrec(data)
def _log_header(self, th):
logger.info("Gromacs version : {0}".format(th.ver_str))
logger.info("tpx version : {0}".format(th.fver))
logger.info("tpx generation : {0}".format(th.fgen))
logger.info("tpx precision : {0}".format(th.precision))
logger.info("tpx file_tag : {0}".format(th.file_tag))
logger.info("tpx natoms : {0}".format(th.natoms))
logger.info("tpx ngtc : {0}".format(th.ngtc))
logger.info("tpx fep_state : {0}".format(th.fep_state))
logger.info("tpx lambda : {0}".format(th.lamb))
logger.debug("tpx bIr (input record): {0}".format(th.bIr))
logger.debug("tpx bTop : {0}".format(th.bTop))
logger.debug("tpx bX : {0}".format(th.bX))
logger.debug("tpx bV : {0}".format(th.bV))
logger.debug("tpx bF : {0}".format(th.bF))
logger.debug("tpx bBox : {0}".format(th.bBox))