Skip to content

Commit

Permalink
cleanup related to BioPDB removal (#777)
Browse files Browse the repository at this point in the history
- Changed tests to eliminate known failures caused by BioPython
- Updated CHANGELOGs and fixed wrong version number in AtomGroup.py,
- fixed indentation issue in PDB.py
- Fixed doc references to version strings and the permissive flags,
- got rid of extraneous text in PrimitivePDBParser, fixed scope of warnings
- Used boolean property of collections as suggested by QC
  • Loading branch information
jdetle authored and orbeckst committed Apr 22, 2016
1 parent ecc8239 commit 5e13207
Show file tree
Hide file tree
Showing 7 changed files with 54 additions and 44 deletions.
49 changes: 27 additions & 22 deletions package/MDAnalysis/coordinates/PDB.py
Original file line number Diff line number Diff line change
Expand Up @@ -56,15 +56,13 @@
Implementations
---------------
PDB I/O is available in the form of the":ref:`permissive<permissive>`"
Reader/Writers.
PDB I/O is available in the form of the Simple PDB Reader/Writers.
..deprecated:: 0.14.1
Only available in the form of the ":ref:":ref:`permissive<permissive>`"
Readers and Writers, see below
..deprecated:: 0.15.0
Readers amd writers solely available in the form of
Simple Readers and Writers, see below.
.. _permissive:
Simple (permissive) PDB Reader and Writer
Simple PDB Reader and Writer
-----------------------------------------
A pure-Python implementation for PDB files commonly encountered in MD
simulations comes under the names :class:`PDBReader` and
Expand Down Expand Up @@ -132,7 +130,7 @@
:members:
..deprecated:: 0.14.1
..deprecated:: 0.15.0
Setting the flag "permissive_pdb_reader" in :data:`MDAnalysis.core.flags`
(see :ref:`flags-label`) to ``False``::
Expand Down Expand Up @@ -262,7 +260,7 @@ def __init__(self, filename, **kwargs):
with util.openany(filename, 'rt') as pdbfile:
for i, line in enumerate(pdbfile):
line = line.strip() # Remove extra spaces
if len(line) == 0: # Skip line if empty
if not line: # Skip line if empty
continue
record = line[:6].strip()

Expand Down Expand Up @@ -326,7 +324,7 @@ def __init__(self, filename, **kwargs):
self.convert_pos_from_native(self.ts._unitcell[:3]) # in-place ! (only lengths)

# No 'MODEL' entries
if len(frames) == 0:
if not frames:
frames[0] = 0

self.frames = frames
Expand Down Expand Up @@ -736,12 +734,12 @@ def _update_frame(self, obj):
Attributes initialized/updated:
* :attr:`PrimitivePDBWriter.obj` (the entity that provides topology information *and*
* :attr:`PDBWriter.obj` (the entity that provides topology information *and*
coordinates, either a :class:`~MDAnalysis.core.AtomGroup.AtomGroup` or a whole
:class:`~MDAnalysis.core.AtomGroup.Universe`)
* :attr:`PrimitivePDBWriter.trajectory` (the underlying trajectory
* :attr:`PDBWriter.trajectory` (the underlying trajectory
:class:`~MDAnalysis.coordinates.base.Reader`)
* :attr:`PrimitivePDBWriter.timestep` (the underlying trajectory
* :attr:`PDBWriter.timestep` (the underlying trajectory
:class:`~MDAnalysis.coordinates.base.Timestep`)
Before calling :meth:`write_next_timestep` this method **must** be
Expand Down Expand Up @@ -855,14 +853,14 @@ def write_next_timestep(self, ts=None, **kwargs):
:Keywords:
*ts*
:class:`base.Timestep` object containing coordinates to be written to trajectory file;
if ``None`` then :attr:`PrimitivePDBWriter.ts`` is tried.
if ``None`` then :attr:`PDBWriter.ts`` is tried.
*multiframe*
``False``: write a single frame (default); ``True`` behave as a trajectory writer
.. Note::
Before using this method with another :class:`base.Timestep` in the *ts*
argument, :meth:`PrimitivePDBWriter._update_frame` *must* be called
argument, :meth:`PDBWriter._update_frame` *must* be called
with the :class:`~MDAnalysis.core.AtomGroup.AtomGroup.Universe` as
its argument so that topology information can be gathered.
'''
Expand Down Expand Up @@ -1027,7 +1025,7 @@ def END(self):
Only a single END record is written. Calling END multiple times has no
effect. Because :meth:`~PDBWriter.close` also calls this
method right before closing the file it is recommended to *not* call
:meth:`~PrimitivePDBWriter.END` explicitly.
:meth:`~PDBWriter.END` explicitly.
.. _END: http://www.wwpdb.org/documentation/format32/sect11.html#END
Expand Down Expand Up @@ -1055,16 +1053,23 @@ def CONECT(self, conect):
conect = "".join(conect)
self.pdbfile.write(self.fmt['CONECT'].format(conect))

warnings.warn('PrimitivePDBReader is identical to the PDBReader,'
'it is deprecated in favor of the shorter name',
category=DeprecationWarning)

class PrimitivePDBReader(PDBReader):
def __init__(self, filename, *args, **kwargs):
warnings.warn('PrimitivePDBReader is identical to the PDBReader,'
' it is deprecated in favor of the shorter name'
' removal targeted for version 0.16.0',
category=DeprecationWarning)
super(PrimitivePDBReader, self).__init__(filename, *args, **kwargs)
format = 'Permissive_PDB'

warnings.warn('PrimitivePDBWriter is now identical to the PDBWriter,'
'it is deprecated in favor of the shorter name',
category=DeprecationWarning)
class PrimitivePDBWriter(PDBWriter):
def __init__(self, filename, *args, **kwargs):
warnings.warn('PrimitivePDBWriter is identical to the Writer,'
'it is deprecated in favor of the shorter name'
' removal targeted for version 0.16.0',
category=DeprecationWarning)
super(PrimitivePDBWriter, self).__init__(filename, *args, **kwargs)
format = 'Permissive_PDB'

class ExtendedPDBReader(PDBReader):
Expand Down
12 changes: 8 additions & 4 deletions package/MDAnalysis/core/AtomGroup.py
Original file line number Diff line number Diff line change
Expand Up @@ -1429,7 +1429,7 @@ def ids(self):
"""
out = np.array([atom.serial for atom in self._atoms])

if not any(out):
out = np.array([atom.id for atom in self._atoms])

Expand Down Expand Up @@ -3334,7 +3334,7 @@ def split(self, level):
return [AtomGroup([a]) for a in self]

if level in ('resid', 'segid'):
warnings.warn("'resid' or 'segid' are no longer allowed levels "
warnings.warn("'resid' or 'segid' are no longer allowed levels "
"in version 0.16.0; instead give "
"'residue' or 'segment', respectively.",
DeprecationWarning)
Expand Down Expand Up @@ -4196,8 +4196,12 @@ def __init__(self, *args, **kwargs):
*permissive*
currently only relevant to PDBs, Readers and Writers and Parsers
returned by setting ::permisive = False are depecrated in favor
returned by setting ::permisive = False are deprecated in favor
of their ::permissive = True counterparts
..deprecated:: 0.15.0
Readers amd writers solely available in the form of
Simple Readers and Writers, see below. Permissive flag
will soon be removed.
*topology_format*
provide the file format of the topology file; ``None`` guesses it from the file
Expand Down Expand Up @@ -4260,7 +4264,7 @@ def __init__(self, *args, **kwargs):
.. versionchanged:: 0.11.0
Added the *is_anchor* and *anchor_name* keywords for finer behavior
control when unpickling instances of :class:`MDAnalysis.core.AtomGroup.AtomGroup`.
.. versionchanged:: 0.14.1
.. versionchanged:: 0.15.0
*permissive* set to ``True`` and set to ``False`` both yield permissive
PDB readers, writers and parsers
"""
Expand Down
2 changes: 1 addition & 1 deletion package/MDAnalysis/topology/PDBParser.py
Original file line number Diff line number Diff line change
Expand Up @@ -100,7 +100,7 @@ def _parseatoms(self):
resid_prev = 0 # resid looping hack
for i, line in enumerate(f):
line = line.strip() # Remove extra spaces
if len(line) == 0: # Skip line if empty
if not line: # Skip line if empty
continue
record = line[:6].strip()

Expand Down
19 changes: 9 additions & 10 deletions package/MDAnalysis/topology/PrimitivePDBParser.py
Original file line number Diff line number Diff line change
Expand Up @@ -45,9 +45,8 @@
:members:
:inherited-members:
..deprecated:: 0.14.1
PrimitivePDBParser now has the same behavior as PDBParser, and will be
absent from future versions
..deprecated:: 0.15.0
PDBParser has been replaced with PrimitivePDBParser.
"""

from __future__ import absolute_import, print_function
Expand All @@ -62,17 +61,17 @@
from .base import TopologyReader


warnings.warn('PrimitivePDBParser is identical to the PDBParser,'
'it is deprecated in favor of the shorter name',
category=DeprecationWarning)

class PrimitivePDBParser(PDBParser.PDBParser):
def __init__(self, *args, **kwargs):
warnings.warn('PrimitivePDBParser is identical to the PDBParser,'
' it is deprecated in favor of the shorter name',
category=DeprecationWarning)
super(PDBParser.PDBParser, self).__init__(*args, **kwargs)
format = 'Permissive_PDB'

warnings.warn('PrimitivePDBParser is identical to the PDBParser,'
'it is deprecated in favor of the shorter name,'
'please use the parse_conect function of that class',
category=DeprecationWarning)
def _parse_conect(conect):

"""parse a CONECT record from pdbs
Parameters
Expand Down
1 change: 1 addition & 0 deletions testsuite/CHANGELOG
Original file line number Diff line number Diff line change
Expand Up @@ -16,6 +16,7 @@ and https://github.com/MDAnalysis/mdanalysis/wiki/UnitTests
??/??/16 orbeckst, jbarnoud, pedrishi, fiona-naughton, jdetle
* 0.15.0
- Added test for weighted rmsd (issue #814)
- removed biopython PDB parser for coordinates and topology (Issue #777)
- metadata update: link download_url to GitHub releases so that
Depsy recognizes contributors (issue #749) and added
@richardjgowers as maintainer
Expand Down
13 changes: 7 additions & 6 deletions testsuite/MDAnalysisTests/coordinates/test_pdb.py
Original file line number Diff line number Diff line change
Expand Up @@ -30,18 +30,19 @@ def setUp(self):
# http://www.wwpdb.org/documentation/format32/sect9.html#ATOM
self.prec = 3

def test_uses_Biopython(self):

def test_uses_PDBReader(self):
from MDAnalysis.coordinates.PDB import PDBReader

assert_(isinstance(self.universe.trajectory, PDBReader),
"failed to choose Biopython PDBReader")
"failed to choose PDBReader")


@knownfailure("Biopython PDB reader does not parse CRYST1", AssertionError)
def test_dimensions(self):
assert_almost_equal(
self.universe.trajectory.ts.dimensions, RefAdKSmall.ref_unitcell,
self.prec,
"Biopython reader failed to get unitcell dimensions from CRYST1")
"PDBReader failed to get unitcell dimensions from CRYST1")


class _PDBMetadata(TestCase, Ref4e43):
Expand Down Expand Up @@ -712,11 +713,11 @@ def setUp(self):
# http://www.wwpdb.org/documentation/format32/sect9.html#ATOM
self.prec = 3

def test_uses_Biopython(self):
def test_uses_PDBReader(self):
from MDAnalysis.coordinates.PDB import PDBReader

assert_(isinstance(self.universe.trajectory, PDBReader),
"failed to choose Biopython PDBReader")
"failed to choose PDBReader")


class TestPDBWriterOccupancies(object):
Expand Down
2 changes: 1 addition & 1 deletion testsuite/MDAnalysisTests/test_streamio.py
Original file line number Diff line number Diff line change
Expand Up @@ -333,7 +333,7 @@ def test_PrimitivePDBReader(self):
u = MDAnalysis.Universe(streamData.as_NamedStream('PDB'))
assert_equal(u.atoms.n_atoms, self.ref_n_atoms)

@knownfailure()

def test_PDBReader(self):
try:
u = MDAnalysis.Universe(streamData.as_NamedStream('PDB'), permissive=False)
Expand Down

0 comments on commit 5e13207

Please sign in to comment.