Statistics: assert that images have the same size

makeStatistics() assumes that image[0,0] and mask[0,0] are the same position.
If mask is bigger than image, but all the unmasked area falls outside the
image given this alignment, then attempting to calculate statistics will
produce a NaN. We avoid this possibility by refusing to proceed if the image
and mask do not have the same dimensions.

include/lsst/afw/math/Statistics.h

@@ -416,6 +416,7 @@ class ImageImposter {
     x_iterator row_end(int) const { return _v.end(); }
     int getWidth() const { return _v.size(); }
     int getHeight() const { return 1; }
+    afw::geom::Extent2I getDimensions() const { return afw::geom::Extent2I(getWidth(), getHeight()); }
 
     bool empty() const { return _v.empty(); }
 private:

src/math/Statistics.cc

@@ -763,6 +763,24 @@ namespace {
 
     template<typename T>
     bool isEmpty(afwImage::Image<T> const& im) { return (im.getWidth() == 0 && im.getHeight() == 0); }
+
+    // Asserts that image dimensions are equal
+    template<typename ImageT1, typename ImageT2>
+    void checkDimensions(ImageT1 const& image1, ImageT2 const& image2)
+    {
+        if (image1.getDimensions() != image2.getDimensions()) {
+            throw LSST_EXCEPT(pexExceptions::InvalidParameterError,
+                              (boost::format("Image sizes don't match: %s vs %s") %
+                               image1.getDimensions() % image2.getDimensions()).str());
+        }
+    }
+
+    // Overloads for MaskImposter (which doesn't have a size)
+    template<typename ImageT, typename PixelT>
+    void checkDimensions(ImageT const& image1, afwMath::MaskImposter<PixelT> const& image2) {}
+    template<typename ImageT, typename PixelT>
+    void checkDimensions(afwMath::MaskImposter<PixelT> const& image1, ImageT const& image2) {}
+
 }
 
 template<typename ImageT, typename MaskT, typename VarianceT, typename WeightT>
@@ -803,6 +821,11 @@ void afwMath::Statistics::doStatistics(
     if (_n == 0) {
         throw LSST_EXCEPT(pexExceptions::InvalidParameterError, "Image contains no pixels");
     }
+    checkDimensions(img, msk);
+    checkDimensions(img, var);
+    if (sctrl.getWeighted()) {
+        checkDimensions(img, weights);
+    }
 
     // Check that an int's large enough to hold the number of pixels
     assert(img.getWidth()*static_cast<double>(img.getHeight()) < std::numeric_limits<int>::max());

tests/statistics.py

@@ -53,7 +53,7 @@
 
 #-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-
 
-class StatisticsTestCase(unittest.TestCase):
+class StatisticsTestCase(utilsTests.TestCase):
     """A test case for Statistics"""
     def setUp(self):
         self.val = 10
@@ -458,6 +458,23 @@ def testMeanClipSingleValue(self):
         stats = afwMath.makeStatistics(img, afwMath.MEANCLIP)
         self.assertEqual(stats.getValue(), 0)
 
+    def testMismatch(self):
+        """Test that we get an exception when there's a size mismatch"""
+        scale = 5
+        dims = self.image.getDimensions()
+        mask = afwImage.MaskU(dims*scale)
+        mask.set(0xFF)
+        ctrl = afwMath.StatisticsControl()
+        ctrl.setAndMask(0xFF)
+        # If it didn't raise, this would result in a NaN (the image data is completely masked).
+        self.assertRaises(lsst.pex.exceptions.InvalidParameterError, afwMath.makeStatistics,
+                          self.image, mask, afwMath.MEDIAN, ctrl)
+        subMask = afwImage.MaskU(mask, afwGeom.Box2I(afwGeom.Point2I(dims*(scale - 1)), dims))
+        subMask.set(0)
+        # Using subMask is successful.
+        self.assertEqual(afwMath.makeStatistics(self.image, subMask, afwMath.MEDIAN, ctrl).getValue(),
+                         self.val)
+
 #-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-
 
 def suite():