Skip to content

Commit

Permalink
Add kwarg to index TPRParser from 0 or 1 (#3152)
Browse files Browse the repository at this point in the history
* add kwarg

* change default to True for 2.0

* fix tests to index resids from 0 again

* Update package/CHANGELOG

Co-authored-by: Oliver Beckstein <orbeckst@gmail.com>

* change default arg for TPR utilities and remove warnings

* updated to 1.1.0

Co-authored-by: Oliver Beckstein <orbeckst@gmail.com>
  • Loading branch information
lilyminium and orbeckst authored Mar 15, 2021
1 parent 347a331 commit 6282385
Show file tree
Hide file tree
Showing 6 changed files with 40 additions and 6 deletions.
2 changes: 2 additions & 0 deletions package/CHANGELOG
Original file line number Diff line number Diff line change
Expand Up @@ -133,6 +133,8 @@ Enhancements
explicit in the topology (Issue #2468, PR #2775)

Changes
* TPRParser now loads TPR files with `tpr_resid_from_one=True` by default,
which starts TPR resid indexing from 1 (instead of 0 as in 1.x) (Issue #2364, PR #3152)
* Introduces encore specific C compiler arguments to allow for lowering of
optimisations on non-x86 platforms (Issue #1389, PR #3149)
* Continuous integration uses mamba rather than conda to install the
Expand Down
20 changes: 18 additions & 2 deletions package/MDAnalysis/topology/TPRParser.py
Original file line number Diff line number Diff line change
Expand Up @@ -176,12 +176,27 @@ class TPRParser(TopologyReaderBase):
"""
format = 'TPR'

def parse(self, **kwargs):
def parse(self, tpr_resid_from_one=True, **kwargs):
"""Parse a Gromacs TPR file into a MDAnalysis internal topology structure.
Parameters
----------
tpr_resid_from_one: bool (optional)
Toggle whether to index resids from 1 or 0 from TPR files.
TPR files index resids from 0 by default, even though GRO and ITP
files index from 1.
Returns
-------
structure : dict
.. versionchanged:: 1.1.0
Added the ``tpr_resid_from_one`` keyword to control if
resids are indexed from 0 or 1. Default ``False``.
.. versionchanged:: 2.0.0
Changed to ``tpr_resid_from_one=True`` by default.
"""
with openany(self.filename, mode='rb') as infile:
tprf = infile.read()
Expand Down Expand Up @@ -219,7 +234,8 @@ def parse(self, **kwargs):
tpr_utils.ndo_real(data, state_ngtc) # relevant to Berendsen tcoupl_lambda

if th.bTop:
tpr_top = tpr_utils.do_mtop(data, th.fver)
tpr_top = tpr_utils.do_mtop(data, th.fver,
tpr_resid_from_one=tpr_resid_from_one)
else:
msg = f"{self.filename}: No topology found in tpr file"
logger.critical(msg)
Expand Down
6 changes: 5 additions & 1 deletion package/MDAnalysis/topology/tpr/utils.py
Original file line number Diff line number Diff line change
Expand Up @@ -46,6 +46,7 @@
The module also contains the :func:`do_inputrec` to read the TPR header with.
"""

import numpy as np
import xdrlib
import struct
Expand Down Expand Up @@ -284,7 +285,7 @@ def extract_box_info(data, fver):
return obj.Box(box, box_rel, box_v)


def do_mtop(data, fver):
def do_mtop(data, fver, tpr_resid_from_one=False):
# mtop: the topology of the whole system
symtab = do_symtab(data)
do_symstr(data, symtab) # system_name
Expand Down Expand Up @@ -371,6 +372,9 @@ def do_mtop(data, fver):
molnums = np.array(molnums, dtype=np.int32)
segids = np.array(segids, dtype=object)
resids = np.array(resids, dtype=np.int32)
if tpr_resid_from_one:
resids += 1

resnames = np.array(resnames, dtype=object)
(residx, new_resids,
(new_resnames,
Expand Down
2 changes: 1 addition & 1 deletion testsuite/MDAnalysisTests/analysis/test_density.py
Original file line number Diff line number Diff line change
Expand Up @@ -148,7 +148,7 @@ class DensityParameters(object):

@pytest.fixture()
def universe(self):
return mda.Universe(self.topology, self.trajectory)
return mda.Universe(self.topology, self.trajectory, tpr_resid_from_one=False)

class TestDensityAnalysis(DensityParameters):
def check_DensityAnalysis(self, ag, ref_meandensity,
Expand Down
4 changes: 2 additions & 2 deletions testsuite/MDAnalysisTests/core/test_atomselections.py
Original file line number Diff line number Diff line change
Expand Up @@ -485,7 +485,7 @@ class TestSelectionsTPR(object):
@staticmethod
@pytest.fixture(scope='class')
def universe():
return MDAnalysis.Universe(TPR,XTC)
return MDAnalysis.Universe(TPR, XTC, tpr_resid_from_one=False)

@pytest.mark.parametrize('selstr', [
'same fragment as bynum 1',
Expand Down Expand Up @@ -1299,7 +1299,7 @@ def test_mass_sel_warning(u_fake_masses):


@pytest.mark.parametrize("selstr,n_res", [
("resnum -10 to 3", 14),
("resnum -10 to 3", 13),
("resnum -5--3", 3), # select -5 to -3
("resnum -3 : -5", 0), # wrong way around
])
Expand Down
12 changes: 12 additions & 0 deletions testsuite/MDAnalysisTests/topology/test_tprparser.py
Original file line number Diff line number Diff line change
Expand Up @@ -39,6 +39,7 @@
)
from MDAnalysisTests.topology.base import ParserBase
import MDAnalysis.topology.TPRParser
import MDAnalysis as mda

BONDED_TPRS = (
TPR510_bonded,
Expand Down Expand Up @@ -287,3 +288,14 @@ def test_elements():
'H', '', 'Na', 'Na', 'Na', 'Na',
], dtype=object)
assert_equal(topology.elements.values[-20:], reference)


@pytest.mark.parametrize("resid_from_one,resid_addition", [
(False, 0),
(True, 1), # status quo for 2.x
])
def test_resids(resid_from_one, resid_addition):
u = mda.Universe(TPR, tpr_resid_from_one=resid_from_one)
resids = np.arange(len(u.residues)) + resid_addition
assert_equal(u.residues.resids, resids,
err_msg="tpr_resid_from_one kwarg not switching resids")

0 comments on commit 6282385

Please sign in to comment.