Repository: ComputationalCryoEM/ASPIRE-Python

Issue:
Fix CoordinateSource for 11028
User reported having trouble loading 11028, which is a centers coordinate source based on MRC files.  There are two issues I found so far loading the data.

First, the mrc header is i) incorrectly reporting 3 values ii) reporting 2 distinct values via `mrcfile`, one with a small floating point error.  This can be seen in the mrcfile(...).header.  ASPIRE-Python code would have ignored (i), but (ii) causes ASPIRE's sanity check to fail and raises an error.  If we expect to load files with bad header data, the code will need to be more permissive.  Lets try emitting a log warning, and if that is annoying loading this data, an actual warning. 

In the diff below I changed the value to vx.y because vx.x appeared to have a small error.  Unsure what to do there in general. I suspect it does not matter outside of populating `pixel_size` which could be overrode later...

```
diff --git a/src/aspire/image/image.py b/src/aspire/image/image.py
index 92f43caa..58e3534e 100644
--- a/src/aspire/image/image.py
+++ b/src/aspire/image/image.py
@@ -866,9 +866,7 @@ class Image:
         # Convert from recarray to single values,
         #   checks uniformity.
         if isinstance(vx, np.recarray):
-            if vx.x != vx.y:
-                raise ValueError(f"Voxel sizes are not uniform: {vx}")
-            vx = vx.x
+            vx = vx.y
```

The second issue is my fault! Shame :D. (https://github.com/ComputationalCryoEM/ASPIRE-Python/commit/fe5a0095e432dec09c71d9a35789a625001a5396).  I previously extended this code to load `TIFF` data, which changed coordinates sources from using `mrcfile` to `aspire.Image`.    All the logic for opening multiple files and dealing with pixel sizes was consolidated to `Image`.  It would seem since that time we have not attempted to load any non-square micrographs.  `Image` raises on non-square data.  We do want to keep the ability to load multiple formats, since that was also an old request, but I don't want opening code duplicated all around.  The lines below are the lines I changed to get current user's data to load.

```
diff --git a/src/aspire/source/coordinates.py b/src/aspire/source/coordinates.py
index 21fec57a..f4f2bf1a 100644
--- a/src/aspire/source/coordinates.py
+++ b/src/aspire/source/coordinates.py
@@ -308,7 +308,7 @@ class CoordinateSource(ImageSource, ABC):
 
         mrc_shapes = np.zeros((self.num_micrographs, 2), dtype=int)
         for i, mrc in enumerate(self.mrc_paths):
-            mrc_shapes[i, :] = Image.load(mrc).resolution
+            mrc_shapes[i, :] = mrcfile.open(mrc)._data.shape
 
         return mrc_shapes
 
@@ -469,7 +469,7 @@ class CoordinateSource(ImageSource, ABC):
         # their origin micrograph
         for mrc_index, coord_list in grouped.items():
             # Load file as 2D numpy array.
-            arr = Image.load(self.mrc_paths[mrc_index]).asnumpy()[0]
+            arr = mrcfile.open(self.mrc_paths[mrc_index])._data
 
             # create iterable of the coordinates in this mrc
             # we don't need to worry about exhausting this iter
```

We're not using any `Image` properties....  so one easy way to do this might be to have an internal `Image._load` that returns the raw image data and pixel size for any of the supported extensions.  Then the `load` method just calls that and returns `Image(im, pixel_size)` line . The two source lines change to call the `_load` method.  This is way less likely to cause trouble than attempting non-square images at this time.  Only a few lines change, no new stuff, just breakup the existing code into two functions.

Of course we should add a test(s) for the rectangular case so that this doesn't regress again.

Discussion:
@j-c-c: Resolved by #1254 

Solution patch:
diff --git a/src/aspire/image/image.py b/src/aspire/image/image.py
index 92f43caa..0395126b 100644
--- a/src/aspire/image/image.py
+++ b/src/aspire/image/image.py
@@ -598,6 +598,10 @@ class Image:
         if self.stack_ndim > 1:
             raise NotImplementedError("`save` is currently limited to 1D image stacks.")
 
+        data = self._data.astype(np.float32)
+        if self.n_images == 1:
+            data = data[0]
+
         if overwrite is None and os.path.exists(mrcs_filepath):
             # If the file exists, append a timestamp to the old file and rename it
             _ = rename_with_timestamp(mrcs_filepath)
@@ -606,13 +610,13 @@ class Image:
 
         with mrcfile.new(mrcs_filepath, overwrite=overwrite) as mrc:
             # original input format (the image index first)
-            mrc.set_data(self._data.astype(np.float32))
+            mrc.set_data(data)
             # Note assigning voxel_size must come after `set_data`
             if self.pixel_size is not None:
                 mrc.voxel_size = self.pixel_size
 
     @staticmethod
-    def load(filepath, dtype=None):
+    def _load_raw(filepath, dtype=None):
         """
         Load raw data from supported files.
 
@@ -621,7 +625,9 @@ class Image:
         :param filepath: File path (string).
         :param dtype: Optionally force cast to `dtype`.
              Default dtype is inferred from the file contents.
-        :return: numpy array of image data.
+        :returns:
+            - numpy array of image data.
+            - pixel size
         """
 
         # Get the file extension
@@ -640,6 +646,23 @@ class Image:
         if dtype is not None:
             im = im.astype(dtype, copy=False)
 
+        return im, pixel_size
+
+    @staticmethod
+    def load(filepath, dtype=None):
+        """
+        Load raw data from supported files.
+
+        Currently MRC and TIFF are supported.
+
+        :param filepath: File path (string).
+        :param dtype: Optionally force cast to `dtype`.
+             Default dtype is inferred from the file contents.
+        :return: Image instance
+        """
+        # Load raw data from filepath with pixel size
+        im, pixel_size = Image._load_raw(filepath, dtype=dtype)
+
         # Return as Image instance
         return Image(im, pixel_size=pixel_size)
 
@@ -867,7 +890,7 @@ class Image:
         #   checks uniformity.
         if isinstance(vx, np.recarray):
             if vx.x != vx.y:
-                raise ValueError(f"Voxel sizes are not uniform: {vx}")
+                logger.warning(f"Voxel sizes are not uniform: {vx}")
             vx = vx.x
 
         # Convert `0` to `None`
diff --git a/src/aspire/source/coordinates.py b/src/aspire/source/coordinates.py
index 21fec57a..495f7a87 100644
--- a/src/aspire/source/coordinates.py
+++ b/src/aspire/source/coordinates.py
@@ -308,7 +308,7 @@ class CoordinateSource(ImageSource, ABC):
 
         mrc_shapes = np.zeros((self.num_micrographs, 2), dtype=int)
         for i, mrc in enumerate(self.mrc_paths):
-            mrc_shapes[i, :] = Image.load(mrc).resolution
+            mrc_shapes[i, :] = Image._load_raw(mrc)[0].shape
 
         return mrc_shapes
 
@@ -469,7 +469,7 @@ class CoordinateSource(ImageSource, ABC):
         # their origin micrograph
         for mrc_index, coord_list in grouped.items():
             # Load file as 2D numpy array.
-            arr = Image.load(self.mrc_paths[mrc_index]).asnumpy()[0]
+            arr = Image._load_raw(self.mrc_paths[mrc_index])[0]
 
             # create iterable of the coordinates in this mrc
             # we don't need to worry about exhausting this iter


Test patch:
diff --git a/tests/test_coordinate_source.py b/tests/test_coordinate_source.py
index c38feff9..5304ef0b 100644
--- a/tests/test_coordinate_source.py
+++ b/tests/test_coordinate_source.py
@@ -1,3 +1,4 @@
+import logging
 import os
 import random
 import shutil
@@ -17,6 +18,8 @@ from aspire.source import BoxesCoordinateSource, CentersCoordinateSource
 from aspire.storage import StarFile
 from aspire.utils import RelionStarFile, importlib_path
 
+logger = logging.getLogger(__name__)
+
 
 class CoordinateSourceTestCase(TestCase):
     def setUp(self):
@@ -696,3 +699,50 @@ class CoordinateSourceTestCase(TestCase):
         self.assertTrue(result_coord.exit_code == 0)
         self.assertTrue(result_star.exit_code == 0)
         self.assertTrue(result_preprocess.exit_code == 0)
+
+
+def create_test_rectangular_micrograph_and_star(tmp_path, voxel_size=(2.0, 2.0, 1.0)):
+    # Create a rectangular micrograph (e.g., 128x256)
+    data = np.random.rand(128, 256).astype(np.float32)
+    mrc_path = tmp_path / "test_micrograph.mrc"
+
+    with mrcfile.new(mrc_path, overwrite=True) as mrc:
+        mrc.set_data(data)
+        mrc.voxel_size = voxel_size
+
+    # Two sample coordinates
+    coordinates = [(50.0, 30.0), (200.0, 100.0)]
+
+    # Write a simple STAR file
+    star_path = tmp_path / "test_coordinates.star"
+    with open(star_path, "w") as f:
+        f.write("data_particles\n\n")
+        f.write("loop_\n")
+        f.write("_rlnCoordinateX #1\n")
+        f.write("_rlnCoordinateY #2\n")
+        for x, y in coordinates:
+            f.write(f"{x:.1f} {y:.1f}\n")
+
+    # Pack files into a list of tuples for consumption by CoordinatSource
+    file_list = [(mrc_path, star_path)]
+
+    return file_list
+
+
+def test_rectangular_coordinate_source(tmp_path):
+    file_list = create_test_rectangular_micrograph_and_star(tmp_path)
+
+    # Check we can instantiate a CoordinateSource with a rectangular micrograph.
+    coord_src = CentersCoordinateSource(file_list, particle_size=32)
+
+    # Check we can access images.
+    _ = coord_src.images[:]
+
+
+def test_coordinate_source_pixel_warning(tmp_path, caplog):
+    # Create micrograph with mismatched pixel dimensions.
+    vx = (2.3, 2.1, 1.0)
+    file_list = create_test_rectangular_micrograph_and_star(tmp_path, voxel_size=vx)
+    with caplog.at_level(logging.WARNING):
+        _ = CentersCoordinateSource(file_list, particle_size=32)
+        assert "Voxel sizes are not uniform" in caplog.text