MRC: force read as volume

Author: PlethoraChutney

gridData/core.py

@@ -196,7 +196,7 @@ class Grid(object):
 
     def __init__(self, grid=None, edges=None, origin=None, delta=None,
                  metadata=None, interpolation_spline_order=3,
-                 file_format=None):
+                 file_format=None, **kwargs):
         # file formats are guessed from extension == lower case key
         self._exporters = {
             'DX': self._export_dx,
@@ -238,7 +238,7 @@ def __init__(self, grid=None, edges=None, origin=None, delta=None,
                     filename = str(grid)
 
             if filename is not None:
-                self.load(filename, file_format=file_format)
+                self.load(filename, file_format=file_format, **kwargs)
             else:
                 self._load(grid, edges, metadata, origin, delta)
 
@@ -532,7 +532,7 @@ def _load(
                 "grid={0} edges={1} origin={2} delta={3}".format(
                     grid, edges, origin, delta))
 
-    def load(self, filename, file_format=None):
+    def load(self, filename, file_format=None, **kwargs):
         """Load saved grid and edges from `filename`
 
         The :meth:`load` method calls the class's constructor method and
@@ -545,32 +545,32 @@ def load(self, filename, file_format=None):
             # are not really a file
             raise IOError(errno.ENOENT, "file not found", filename)
         loader = self._get_loader(filename, file_format=file_format)
-        loader(filename)
+        loader(filename, **kwargs)
 
-    def _load_python(self, filename):
+    def _load_python(self, filename, **kwargs):
         with open(filename, 'rb') as f:
             saved = pickle.load(f)
         self._load(grid=saved['grid'],
                    edges=saved['edges'],
                    metadata=saved['metadata'])
 
-    def _load_mrc(self, filename):
+    def _load_mrc(self, filename, **kwargs):
         """Initializes Grid from a MRC/CCP4 file."""
-        mrcfile = mrc.MRC(filename)
+        mrcfile = mrc.MRC(filename, **kwargs)
         grid, edges = mrcfile.histogramdd()
         self._load(grid=grid, edges=edges, metadata=self.metadata)
         # Store header for access from Grid object (undocumented)
         # https://github.com/MDAnalysis/GridDataFormats/pull/100#discussion_r782604833
         self._mrc_header = mrcfile.header.copy()
 
-    def _load_dx(self, filename):
+    def _load_dx(self, filename, **kwargs):
         """Initializes Grid from a OpenDX file."""
         dx = OpenDX.field(0)
         dx.read(filename)
         grid, edges = dx.histogramdd()
         self._load(grid=grid, edges=edges, metadata=self.metadata)
 
-    def _load_plt(self, filename):
+    def _load_plt(self, filename, **kwargs):
         """Initialize Grid from gOpenMol plt file."""
         g = gOpenMol.Plt()
         g.read(filename)

gridData/mrc.py

@@ -43,6 +43,12 @@ class MRC(object):
     ----------
     filename : str (optional)
        input file (or stream), can be compressed
+    is_volume : bool (optional)
+       Force input file to be interpreted as a volume. If None (default),
+       use the MRC file header. If True, force the file to be read as a
+       volume. If False, the file is intepreted as a 2D image stack and
+       raise a ValueError.
+        
 
     Raises
     ------
@@ -86,17 +92,22 @@ class MRC(object):
 
     """
 
-    def __init__(self, filename=None):
+    def __init__(self, filename=None, is_volume=None):
         self.filename = filename
         if filename is not None:
-            self.read(filename)
+            self.read(filename, is_volume)
 
-    def read(self, filename):
+    def read(self, filename, is_volume=None):
         """Populate the instance from the MRC/CCP4 file *filename*."""
         if filename is not None:
             self.filename = filename
         with mrcfile.open(filename) as mrc:
-            if mrc.data is None or len(mrc.data.shape) != 3:                           #pragma: no cover
+            if is_volume is None:
+                is_volume = mrc.is_volume()
+            elif is_volume:
+                # non 3D volumes should always fail, regardless of is_volume value
+                is_volume = mrc.data is not None and len(mrc.data.shape) == 3
+            if not is_volume:                           #pragma: no cover
                 raise ValueError(
                     "MRC file {} is not a volumetric density.".format(filename))
             self.header = h = mrc.header.copy()

gridData/tests/datafiles/__init__.py

@@ -7,6 +7,7 @@
 DX = importlib_resources.files(__name__) / 'test.dx'
 DXGZ = importlib_resources.files(__name__) / 'test.dx.gz'
 CCP4 = importlib_resources.files(__name__) / 'test.ccp4'
+ISPG_0 = importlib_resources.files(__name__) / "ispg_0.mrc"
 # from http://www.ebi.ac.uk/pdbe/coordinates/files/1jzv.ccp4
 # (see issue #57)
 CCP4_1JZV = importlib_resources.files(__name__) / '1jzv.ccp4'

gridData/tests/datafiles/ispg_0.mrc

@@ -113,6 +113,17 @@ def test_triclinic_ValueError():
                        "supported, not"):
         Grid(datafiles.MRC_EMD3001, file_format="MRC")
 
+def test_mrcfile_volume_check():
+    with pytest.raises(ValueError, match="is not a volumetric density"):
+        Grid(datafiles.ISPG_0)
+
+def test_mrcfile_volume_force():
+    Grid(datafiles.ISPG_0, is_volume=True)
+
+def test_mrcfile_force_fail():
+    with pytest.raises(ValueError, match="is not a volumetric density"):
+        Grid(datafiles.CCP4_1JZV, is_volume=False)
+
 class TestGridMRC():
     @pytest.fixture(scope="class")
     def grid(self):