diff --git a/package/CHANGELOG b/package/CHANGELOG index 278e6df0db3..28cbc78ae1b 100644 --- a/package/CHANGELOG +++ b/package/CHANGELOG @@ -79,6 +79,7 @@ Changes * changed the water bridge analysis output format (PR #2087) Fixes + * fixed ChainReader setting format with format keyword (Issue #2334) * fixed lack of check for scaling of NCDFReader velocities (Issue #2323) * fixed PDBReader and PDBWriter newlines for PDB header (Issue #2324) * fixes ProgressMeter issues with older Jupyter Lab versions (Issue #2078) diff --git a/package/MDAnalysis/coordinates/core.py b/package/MDAnalysis/coordinates/core.py index d36a3539879..8f3c0b8aea4 100644 --- a/package/MDAnalysis/coordinates/core.py +++ b/package/MDAnalysis/coordinates/core.py @@ -46,7 +46,7 @@ from ..core._get_readers import get_reader_for, get_writer_for -def reader(filename, **kwargs): +def reader(filename, format=None, **kwargs): """Provide a trajectory reader instance for *filename*. This function guesses the file format from the extension of *filename* and @@ -78,10 +78,14 @@ def reader(filename, **kwargs): if isinstance(filename, tuple): Reader = get_reader_for(filename[0], format=filename[1]) - return Reader(filename[0], **kwargs) + filename = filename[0] else: - Reader = get_reader_for(filename) + Reader = get_reader_for(filename, format=format) + try: return Reader(filename, **kwargs) + except ValueError: + raise TypeError('Unable to read {fn} with {r}.'.format(fn=filename, + r=Reader)) def writer(filename, n_atoms=None, **kwargs): diff --git a/package/MDAnalysis/core/universe.py b/package/MDAnalysis/core/universe.py index b02f462e69d..5594ce4a209 100644 --- a/package/MDAnalysis/core/universe.py +++ b/package/MDAnalysis/core/universe.py @@ -601,7 +601,7 @@ def load_new(self, filename, format=None, in_memory=False, **kwargs): # supply number of atoms for readers that cannot do it for themselves kwargs['n_atoms'] = self.atoms.n_atoms - self.trajectory = reader(filename, **kwargs) + self.trajectory = reader(filename, format=format, **kwargs) if self.trajectory.n_atoms != len(self.atoms): raise ValueError("The topology and {form} trajectory files don't" " have the same number of atoms!\n" diff --git a/testsuite/MDAnalysisTests/coordinates/test_chainreader.py b/testsuite/MDAnalysisTests/coordinates/test_chainreader.py index 0b41aed0eb0..05e78d5b3f1 100644 --- a/testsuite/MDAnalysisTests/coordinates/test_chainreader.py +++ b/testsuite/MDAnalysisTests/coordinates/test_chainreader.py @@ -200,13 +200,29 @@ def test_set_all_format_tuples(self): assert universe.trajectory.n_frames == 21 assert_equal(universe.trajectory.filenames, [PDB, XTC, TRR]) + def test_set_format_tuples_and_format(self): + universe = mda.Universe(GRO, [(PDB, 'pdb'), GRO, GRO, (XTC, 'xtc'), + (TRR, 'trr')], format='gro') + assert universe.trajectory.n_frames == 23 + assert_equal(universe.trajectory.filenames, [PDB, GRO, GRO, XTC, TRR]) + + with pytest.raises(TypeError) as errinfo: + mda.Universe(GRO, [(PDB, 'pdb'), GRO, GRO, (XTC, 'xtc'), + (TRR, 'trr')], format='pdb') + assert 'Unable to read' in str(errinfo.value) + + def test_set_one_format_tuple(self): universe = mda.Universe(PSF, [(PDB_small, 'pdb'), DCD]) assert universe.trajectory.n_frames == 99 def test_set_all_formats(self): - universe = mda.Universe(PSF, [PDB_small, PDB_closed], format='pdb') - assert universe.trajectory.n_frames == 2 + with pytest.raises(TypeError) as errinfo: + mda.Universe(PDB, [PDB, GRO], format='gro') + assert 'Unable to read' in str(errinfo.value) + + universe = mda.Universe(GRO, [PDB, PDB, PDB], format='pdb') + assert_equal(universe.trajectory.filenames, [PDB, PDB, PDB]) def build_trajectories(folder, sequences, fmt='xtc'):