Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Skip observables contained in particle groups #4615

Merged
Merged
Show file tree
Hide file tree
Changes from 15 commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
1 change: 1 addition & 0 deletions package/AUTHORS
Original file line number Diff line number Diff line change
Expand Up @@ -240,6 +240,7 @@ Chronological list of authors
- Sampurna Mukherjee
- Leon Wehrhan
- Valerij Talagayev
- Fabian Zills



Expand Down
4 changes: 3 additions & 1 deletion package/CHANGELOG
Original file line number Diff line number Diff line change
Expand Up @@ -17,11 +17,13 @@ The rules for this file:
??/??/?? IAlibay, HeetVekariya, marinegor, lilyminium, RMeli,
ljwoods2, aditya292002, pstaerk, PicoCentauri, BFedder,
tyler.je.reddy, SampurnaM, leonwehrhan, kainszs, orionarcher,
yuxuanzhuang
yuxuanzhuang, PythonFZ

* 2.8.0

Fixes
* Do not raise an Error reading H5MD files with datasets like
`observables/<particle>/<property>` (part of Issue #4598, PR #4615)
* Fix PSFParser error when encoutering string-like resids
* (Issue #2053, Issue #4189 PR #4582)
* Fix `MDAnalysis.analysis.align.AlignTraj` not accepting writer kwargs
Expand Down
24 changes: 20 additions & 4 deletions package/MDAnalysis/coordinates/H5MD.py
Original file line number Diff line number Diff line change
Expand Up @@ -207,6 +207,7 @@
:members:

"""
import warnings

import numpy as np
import MDAnalysis as mda
Expand Down Expand Up @@ -673,10 +674,25 @@ def _read_frame(self, frame):
def _copy_to_data(self):
"""assigns values to keys in data dictionary"""

if 'observables' in self._file:
for key in self._file['observables'].keys():
self.ts.data[key] = self._file['observables'][key][
'value'][self._frame]
if "observables" in self._file:
orbeckst marked this conversation as resolved.
Show resolved Hide resolved
for key in self._file["observables"].keys():
try:
# if has value as subkey read directly into data
if "value" in self._file["observables"][key]:
self.ts.data[key] = self._file["observables"][key][
"value"
][self._frame]
# if value is not a subkey, read dict of subkeys
else:
for subkey in self._file["observables"][key].keys():
self.ts.data[key + "/" + subkey] = self._file["observables"][key][
subkey
]["value"][self._frame]
except KeyError:
warnings.warn(
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

If it's not legal, shouldn't we be failing? Or is the reasoning that we don't overly care about observables as opposed to x/v/f?

Perhaps add a comment explaining your rationale.

Also, given that this is not tested, add a #pragma: nocover (or whatever you have to do to show that this code block isn't covered). It will tell anyone else reading the code that it's un-tested.

Tested would be better, of course, but if we don't care about observables then I don't insist on a test.

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I'd be in favor of raising an Error here as well. Shall I create an additional h5md file which contains wrong data to test against this?

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

That would be ideal!

f"Could not read {key} from observables group. "
"Possibly not a legal H5MD observable specification."
)

# pulls 'time' and 'step' out of first available parent group
for name, value in self._has.items():
Expand Down
23 changes: 22 additions & 1 deletion testsuite/MDAnalysisTests/coordinates/test_h5md.py
Original file line number Diff line number Diff line change
Expand Up @@ -7,11 +7,12 @@
from MDAnalysis.coordinates.H5MD import HAS_H5PY
if HAS_H5PY:
import h5py
from MDAnalysis.coordinates.H5MD import H5MDReader
from MDAnalysis.exceptions import NoDataError
from MDAnalysisTests import make_Universe
from MDAnalysisTests.datafiles import (H5MD_xvf, TPR_xvf, TRR_xvf,
COORDINATES_TOPOLOGY,
COORDINATES_H5MD)
COORDINATES_H5MD, H5MD_energy)
from MDAnalysisTests.coordinates.base import (MultiframeReaderTest,
BaseReference, BaseWriterTest,
assert_timestep_almost_equal)
Expand Down Expand Up @@ -894,3 +895,23 @@ def test_writer_no_h5py(self, Writer, outfile):
u.atoms.n_atoms) as W:
for ts in u.trajectory:
W.write(universe)

@pytest.mark.skipif(not HAS_H5PY, reason="h5py not installed")
class TestH5MDReaderWithObservables(object):
"""Read H5MD file with 'observables/atoms/energy'."""

prec = 3
ext = 'h5md'


def test_read_h5md_issue4598(self):
"""Read a H5MD file with observables.

The reader will ignore the 'observables/atoms/energy'.
"""

u = mda.Universe.empty(n_atoms=108, trajectory=True)
reader = H5MDReader(H5MD_energy)
u.trajectory = reader
for ts in u.trajectory:
assert "atoms/energy" in ts.data
Binary file added testsuite/MDAnalysisTests/data/cu.h5
Binary file not shown.
2 changes: 2 additions & 0 deletions testsuite/MDAnalysisTests/datafiles.py
Original file line number Diff line number Diff line change
Expand Up @@ -77,6 +77,7 @@
"GRO_sameresid_diffresname", # Case where two residues share the same resid
"PDB_xvf", "TPR_xvf", "TRR_xvf", # Gromacs coords/veloc/forces (cobrotoxin, OPLS-AA, Gromacs 4.5.5 tpr)
"H5MD_xvf", # TPR_xvf + TRR_xvf converted to h5md format
"H5MD_energy", # H5MD trajectory with observables/atoms/energy
"XVG_BZ2", # Compressed xvg file about cobrotoxin
"PDB_xlserial",
"TPR400", "TPR402", "TPR403", "TPR404", "TPR405", "TPR406", "TPR407",
Expand Down Expand Up @@ -372,6 +373,7 @@
TPR_xvf = (_data_ref / 'cobrotoxin.tpr').as_posix()
TRR_xvf = (_data_ref / 'cobrotoxin.trr').as_posix()
H5MD_xvf = (_data_ref / 'cobrotoxin.h5md').as_posix()
H5MD_energy = (_data_ref / 'cu.h5').as_posix()
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

is it the convention for other h5md files out there to use the *.h5 extension rather than *.h5md?

Copy link
Contributor Author

@PythonFZ PythonFZ Jun 19, 2024

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I can rename the file to *.h5md. I've seen both and personally prefer *.h5 because I was lazy and did not configure the VS Code H5WEB extension to also read *.h5md. Afaik H5MD is the only format specification for H5 files, so there is no risk of confusion.

Edit: I haven't seen a recommendation for the *.h5md suffix in the H5MD definition so I would actually prefer keeping the more general *.h5.

XVG_BZ2 = (_data_ref / 'cobrotoxin_protein_forces.xvg.bz2').as_posix()

XPDB_small = (_data_ref / '5digitResid.pdb').as_posix()
Expand Down
Loading