diff --git a/pyproject.toml b/pyproject.toml
index 54c9017..29b8bbe 100644
--- a/pyproject.toml
+++ b/pyproject.toml
@@ -6,7 +6,7 @@ build-backend = "setuptools.build_meta"
 
 [project]
 name = "yirgacheffe"
-version = "0.8.6"
+version = "0.8.7"
 description = "Abstraction of gdal datasets for doing basic math operations"
 readme = "README.md"
 authors = [{ name = "Michael Dales", email = "mwd24@cam.ac.uk" }]
diff --git a/tests/test_operators.py b/tests/test_operators.py
index fe73e5c..ab270d6 100644
--- a/tests/test_operators.py
+++ b/tests/test_operators.py
@@ -501,8 +501,57 @@ def test_write_mulitband_raster() -> None:
 
             assert (expected == actual).all()
 
-def test_sum_float32(monkeypatch) -> None:
+def test_save_and_sum_float32(monkeypatch) -> None:
+    random.seed(42)
+    data = []
+    for _ in range(10):
+        row = []
+        for _ in range(10):
+            row.append(random.random())
+        data.append(row)
 
+    data1 = np.array(data, dtype=np.float32)
+    layer1 = RasterLayer(gdal_dataset_with_data((0.0, 0.0), 0.02, data1))
+    assert layer1.datatype == gdal.GDT_Float32
+
+    # Sum forces things to float64
+    expected = np.sum(data1.astype(np.float64))
+
+    with monkeypatch.context() as m:
+        for blocksize in range(1,11):
+            m.setattr(yirgacheffe.constants, "YSTEP", blocksize)
+            with RasterLayer.empty_raster_layer_like(layer1) as store:
+                actual = layer1.save(store, and_sum=True)
+            assert expected == actual
+
+def test_parallel_save_and_sum_float32(monkeypatch) -> None:
+    random.seed(42)
+    data = []
+    for _ in range(10):
+        row = []
+        for _ in range(10):
+            row.append(random.random())
+        data.append(row)
+
+    with tempfile.TemporaryDirectory() as tempdir:
+        path1 = os.path.join(tempdir, "test1.tif")
+        data1 = np.array(data, dtype=np.float32)
+        dataset1 = gdal_dataset_with_data((0.0, 0.0), 0.02, data1, filename=path1)
+        dataset1.Close()
+        layer1 = RasterLayer.layer_from_file(path1)
+        assert layer1.datatype == gdal.GDT_Float32
+
+        # Sum forces things to float64
+        expected = np.sum(data1.astype(np.float64))
+
+        with monkeypatch.context() as m:
+            for blocksize in range(1,11):
+                m.setattr(yirgacheffe.constants, "YSTEP", blocksize)
+                with RasterLayer.empty_raster_layer_like(layer1) as store:
+                    actual = layer1.parallel_save(store, and_sum=True)
+                assert expected == actual
+
+def test_sum_float32(monkeypatch) -> None:
     random.seed(42)
     data = []
     for _ in range(10):
diff --git a/yirgacheffe/operators.py b/yirgacheffe/operators.py
index 351fdef..597f44a 100644
--- a/yirgacheffe/operators.py
+++ b/yirgacheffe/operators.py
@@ -243,7 +243,7 @@ def save(self, destination_layer, and_sum=False, callback=None, band=1):
                 yoffset + destination_window.yoff,
             )
             if and_sum:
-                total += np.sum(chunk)
+                total += np.sum(chunk.astype(np.float64))
         if callback:
             callback(1.0)
 
@@ -316,7 +316,6 @@ def parallel_save(self, destination_layer, and_sum=False, callback=None, paralle
             gdal.GDT_UInt64:  np.dtype('uint64'),
         }[band.DataType]
 
-
         with multiprocessing.Manager() as manager:
             with SharedMemoryManager() as smm:
 
@@ -375,7 +374,7 @@ def parallel_save(self, destination_layer, and_sum=False, callback=None, paralle
                         yoffset + destination_window.yoff,
                     )
                     if and_sum:
-                        total += np.sum(arr[0:step])
+                        total += np.sum(np.array(arr[0:step]).astype(np.float64))
                     sem.release()
                     retired_blocks += 1
                     if callback: