Source code for MDAnalysis.topology.TOPParser

# -*- Mode: python; tab-width: 4; indent-tabs-mode:nil; coding:utf-8 -*-
# vim: tabstop=4 expandtab shiftwidth=4 softtabstop=4 
#
# MDAnalysis --- http://www.MDAnalysis.org
# Copyright (c) 2006-2015 Naveen Michaud-Agrawal, Elizabeth J. Denning, Oliver Beckstein
# and contributors (see AUTHORS for the full list)
#
# Released under the GNU Public Licence, v2 or any higher version
#
# Please cite your use of MDAnalysis in published work:
# 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
#

"""
AMBER PRMTOP topology parser
============================

Reads a  AMBER top file to build the system. It uses atom types,
partial charges and masses from the PRMTOP file.

The format is defined in `PARM parameter/topology file specification`_.
The reader tries to detect if it is a newer (AMBER 12?) file format
by looking for the flag "ATOMIC_NUMBER".

.. Note::

   The Amber charge is converted to electron charges as used in
   MDAnalysis and other packages. To get back Amber charges, multiply
   by 18.2223.

.. _`PARM parameter/topology file specification`:
   http://ambermd.org/formats.html#topology

Classes
-------

.. autoclass:: TOPParser
   :members:
   :inherited-members:

"""
from __future__ import absolute_import
from six.moves import range

import functools
from math import ceil

from ..core.AtomGroup import Atom
from ..units import convert
from ..lib.util import openany, FORTRANReader
from ..core import flags
from .base import TopologyReader


[docs]class TOPParser(TopologyReader): """Reads topology information from an AMBER top file. It uses atom types, partial charges and masses from the PRMTOP file. The format is defined in `PARM parameter/topology file specification`_. The reader tries to detect if it is a newer (AMBER 12?) file format by looking for the flag "ATOMIC_NUMBER". .. _`PARM parameter/topology file specification`: http://ambermd.org/formats.html#topology .. versionchanged:: 0.7.6 parses both amber10 and amber12 formats """ format = ['TOP', 'PRMTOP', 'PARM7']
[docs] def parse(self): """Parse Amber PRMTOP topology file *filename*. :Returns: MDAnalysis internal *structure* dict. """ formatversion = 10 with openany(self.filename, 'rt') as topfile: for line in topfile: if line.startswith("%FLAG ATOMIC_NUMBER"): formatversion = 12 break if formatversion == 12: sections = [ ("ATOM_NAME", 1, 20, self._parseatoms, "_name", 0), ("CHARGE", 1, 5, self._parsesection, "_charge", 0), ("ATOMIC_NUMBER", 1, 10, self._parsesectionint, "_skip", 0), ("MASS", 1, 5, self._parsesection, "_mass", 0), ("ATOM_TYPE_INDEX", 1, 10, self._parsesectionint, "_atom_type", 0), ("NUMBER_EXCLUDED_ATOMS", 1, 10, self._parseskip, "_skip", 8), ("NONBONDED_PARM_INDEX", 1, 10, self._parseskip, "_skip", 8), ("RESIDUE_LABEL", 1, 20, self._parseatoms, "_resname", 11), ("RESIDUE_POINTER", 2, 10, self._parsesectionint, "_respoint", 11), ] #("BOND_FORCE_CONSTANT", 1, 5, self._parseskip,"_skip",8), #("BOND_EQUIL_VALUE", 1, 5, self._parseskip,"_skip",8), #("ANGLE_FORCE_CONSTANT", 1, 5, self._parseskip,"_skip",8), #("ANGLE_EQUIL_VALUE", 1, 5, self._parseskip,"_skip",8), #("DIHEDRAL_FORCE_CONSTANT", 1, 5, self._parseskip,"_skip",8), #("DIHEDRAL_PERIODICITY", 1, 5, self._parseskip,"_skip",8), #("DIHEDRAL_PHASE", 1, 5, self._parseskip,"_skip",8), #("SOLTY", 1, 5, self._parseskip,"_skip",8), #("LENNARD_JONES_ACOEF", 1, 5, self._parseskip,"_skip",8), #("LENNARD_JONES_BCOEF", 1, 5, self._parseskip,"_skip",8), #("BONDS_INC_HYDROGEN", 2, 4, self._parsebond, "_bonds",2), #("ANGLES_INC_HYDROGEN", 3, 3, self._parsesection, "_angles"), #("DIHEDRALS_INC_HYDROGEN", 4, 2, self._parsesection, "_dihe"), #("NIMPHI", 4, 2, self._parsesection, "_impr"), #("NDON", 2, 4, self._parsesection,"_donors"), #("NACC", 2, 4, self._parsesection,"_acceptors"), elif formatversion == 10: sections = [ ("ATOM_NAME", 1, 20, self._parseatoms, "_name", 0), ("CHARGE", 1, 5, self._parsesection, "_charge", 0), ("MASS", 1, 5, self._parsesection, "_mass", 0), ("ATOM_TYPE_INDEX", 1, 10, self._parsesectionint, "_atom_type", 0), ("NUMBER_EXCLUDED_ATOMS", 1, 10, self._parseskip, "_skip", 8), ("NONBONDED_PARM_INDEX", 1, 10, self._parseskip, "_skip", 8), ("RESIDUE_LABEL", 1, 20, self._parseatoms, "_resname", 11), ("RESIDUE_POINTER", 2, 10, self._parsesectionint, "_respoint", 11), ] #("BOND_FORCE_CONSTANT", 1, 5, self._parseskip,"_skip",8), #("BOND_EQUIL_VALUE", 1, 5, self._parseskip,"_skip",8), #("ANGLE_FORCE_CONSTANT", 1, 5, self._parseskip,"_skip",8), #("ANGLE_EQUIL_VALUE", 1, 5, self._parseskip,"_skip",8), #("DIHEDRAL_FORCE_CONSTANT", 1, 5, self._parseskip,"_skip",8), #("DIHEDRAL_PERIODICITY", 1, 5, self._parseskip,"_skip",8), #("DIHEDRAL_PHASE", 1, 5, self._parseskip,"_skip",8), #("SOLTY", 1, 5, self._parseskip,"_skip",8), #("LENNARD_JONES_ACOEF", 1, 5, self._parseskip,"_skip",8), #("LENNARD_JONES_BCOEF", 1, 5, self._parseskip,"_skip",8), #("BONDS_INC_HYDROGEN", 2, 4, self._parsebond, "_bonds",2), #("ANGLES_INC_HYDROGEN", 3, 3, self._parsesection, "_angles"), #("DIHEDRALS_INC_HYDROGEN", 4, 2, self._parsesection, "_dihe")] #("NIMPHI", 4, 2, self._parsesection, "_impr"), #("NDON", 2, 4, self._parsesection,"_donors"), #("NACC", 2, 4, self._parsesection,"_acceptors")] # Open and check top validity # Reading header info POINTERS with openany(self.filename, 'rt') as topfile: next_line = functools.partial(next, topfile) header = next_line() if header[:3] != "%VE": raise ValueError("{0} is not a valid TOP file. %VE Missing in header".format(topfile)) title = next_line().split() if not (title[1] == "TITLE"): raise ValueError("{0} is not a valid TOP file. 'TITLE' missing in header".format(topfile)) while header[:14] != '%FLAG POINTERS': header = next_line() header = next_line() topremarks = [next_line().strip() for i in range(4)] sys_info = [int(k) for i in topremarks for k in i.split()] structure = {} final_structure = {} try: for info in sections: self._parse_sec(sys_info, info, next_line, structure, final_structure) except StopIteration: raise ValueError("The TOP file didn't contain the minimum" " required section of ATOM_NAME") # Completing info respoint to include all atoms in last resid structure["_respoint"].append(sys_info[0]) structure["_respoint"][-1] = structure["_respoint"][-1] + 1 atoms = [None, ]*sys_info[0] j = 0 segid = "SYSTEM" for i in range(sys_info[0]): charge = convert(structure["_charge"][i], 'Amber', flags['charge_unit']) if structure["_respoint"][j] <= i+1 < structure["_respoint"][j+1]: resid = j + 1 resname = structure["_resname"][j] else: j += 1 resid = j + 1 resname = structure["_resname"][j] mass = structure["_mass"][i] atomtype = structure["_atom_type"][i] atomname = structure["_name"][i] #segid = 'SYSTEM' # does not exist in Amber atoms[i] = Atom(i, atomname, atomtype, resname, resid, segid, mass, charge, universe=self._u) final_structure["atoms"] = atoms final_structure["_numatoms"] = sys_info[0] return final_structure
def _parse_sec(self, sys_info, section_info, next_line, structure, final_structure): desc, atoms_per, per_line, parsefunc, data_struc, sect_num = section_info # Get the number num = sys_info[sect_num] if data_struc in ["_resname", "_bond"]: pass else: header = next_line() # Now figure out how many lines to read numlines = int(ceil(float(num)/per_line)) #print data_struc, numlines if parsefunc == self._parsebond: parsefunc(next_line, atoms_per, data_struc, final_structure, numlines) else: parsefunc(next_line, atoms_per, data_struc, structure, numlines) def _parseskip(self, lines, atoms_per, attr, structure, numlines): while (lines()[:5] != "%FLAG"): pass def _parsebond(self, lines, atoms_per, attr, structure, numlines): section = [] # [None,]*numlines for i in range(numlines): l = lines() # Subtract 1 from each number to ensure zero-indexing for the atoms f = map(int, l.split()) fields = [a-1 for a in f] for j in range(0, len(fields), atoms_per): section.append(tuple(fields[j:j+atoms_per])) structure[attr] = section def _parsesectionint(self, lines, atoms_per, attr, structure, numlines): section = [] # [None,]*numlines y = lines().strip("%FORMAT(") y.strip(")") x = FORTRANReader(y) for i in range(numlines): l = lines() # Subtract 1 from each number to ensure zero-indexing for the atoms try: for j in range(len(x.entries)): section.append(int(l[x.entries[j].start:x.entries[j].stop].strip())) except: continue structure[attr] = section def _parsesection(self, lines, atoms_per, attr, structure, numlines): section = [] # [None,]*numlines y = lines().strip("%FORMAT(") y.strip(")") x = FORTRANReader(y) for i in range(numlines): l = lines() # Subtract 1 from each number to ensure zero-indexing for the atoms try: for j in range(0, len(x.entries)): section.append(float(l[x.entries[j].start:x.entries[j].stop].strip())) except: continue structure[attr] = section def _parseatoms(self, lines, atoms_per, attr, structure, numlines): section = [] # [None,]*numlines y = lines().strip("%FORMAT(") y.strip(")") x = FORTRANReader(y) for i in range(numlines): l = lines() # Subtract 1 from each number to ensure zero-indexing for the atoms for j in range(0, len(x.entries)): if l[x.entries[j].start:x.entries[j].stop] != '': #print l[x.entries[j].start:x.entries[j].stop] section.append(l[x.entries[j].start:x.entries[j].stop].strip()) else: continue structure[attr] = section