Skip to content

Commit

Permalink
Fix bug in GRO reader (mosdef-hub#775)
Browse files Browse the repository at this point in the history
* parse atom/residue name/id by slicing string

* Update gmso/formats/gro.py

* reorder error checking block, update test

* add warnings if velocity info is present
  • Loading branch information
daico007 authored Jan 15, 2024
1 parent aab43c0 commit c99799a
Show file tree
Hide file tree
Showing 2 changed files with 18 additions and 14 deletions.
30 changes: 17 additions & 13 deletions gmso/formats/gro.py
Original file line number Diff line number Diff line change
Expand Up @@ -42,8 +42,8 @@ def read_gro(filename):
Gro files do not specify connections between atoms, the returned topology
will not have connections between sites either.
Currently this implementation does not support a gro file with more than 1
frame.
Currently this implementation does not support parsing velocities from a gro file or gro file
with more than 1 frame.
All residues and resid information from the gro file are currently lost
when converting to `topology`.
Expand All @@ -57,32 +57,36 @@ def read_gro(filename):
coords = u.nm * np.zeros(shape=(n_atoms, 3))
for row, _ in enumerate(coords):
line = gro_file.readline()
content = line.split()
if not line:
msg = (
"Incorrect number of lines in .gro file. Based on the "
"number in the second line of the file, {} rows of"
"atoms were expected, but at least one fewer was found."
)
raise ValueError(msg.format(n_atoms))
res_id = int(line[:5].strip())
res_name = line[5:10].strip()
atom_name = line[10:15].strip()
atom_id = line[15:20].strip()

res = content[0]
atom_name = content[1]
atom_id = content[2]
positions = line[20:].split()
coords[row] = u.nm * np.array(
[
float(content[3]),
float(content[4]),
float(content[5]),
float(positions[0]),
float(positions[1]),
float(positions[2]),
]
)
site = Atom(name=atom_name, position=coords[row])

r = re.compile("([0-9]+)([a-zA-Z]+)")
m = r.match(res)
site.molecule = (m.group(2), int(m.group(1)))
site.residue = (m.group(2), int(m.group(1)))
site.molecule = (res_name, res_id)
site.residue = (res_name, res_id)
top.add_site(site, update_types=False)

if len(positions) == 6:
warnings.warn(
"Velocity information presents but will not be parsed."
)
top.update_topology()

# Box information
Expand Down
2 changes: 1 addition & 1 deletion gmso/tests/test_gro.py
Original file line number Diff line number Diff line change
Expand Up @@ -39,7 +39,7 @@ def test_read_gro(self):
)

def test_wrong_n_atoms(self):
with pytest.raises(IndexError):
with pytest.raises(ValueError):
Topology.load(get_fn("too_few_atoms.gro"))
with pytest.raises(ValueError):
Topology.load(get_fn("too_many_atoms.gro"))
Expand Down

0 comments on commit c99799a

Please sign in to comment.