Merge pull request #932 from heinsimon/issue_917_fixing_points_to_sph…

…eres_3D

Fix `points_to_spheres` bug in 3D #bug

src/porespy/tools/_sphere_insertions.py

@@ -36,7 +36,7 @@ def points_to_spheres(im):
     from scipy.spatial import distance_matrix
     if im.ndim == 3:
         x, y, z = np.where(im > 0)
-        coords = np.vstack((x, y, z)).T
+        coords = np.vstack((x, y, z))
     else:
         x, y = np.where(im > 0)
         coords = np.vstack((x, y))
@@ -46,7 +46,10 @@ def points_to_spheres(im):
         dmap[mask] = np.inf
         r = np.around(dmap.min(axis=0)/2, decimals=0).astype(int)
     else:
-        r = im[x, y].flatten()
+        if im.ndim == 3:
+            r = im[x, y, z].flatten()
+        else:
+            r = im[x, y].flatten()
     im_spheres = np.zeros_like(im, dtype=bool)
     im_spheres = _insert_disks_at_points_parallel(
         im_spheres,

test/unit/test_tools.py

@@ -428,6 +428,18 @@ def test_tic_toc(self):
         t = toc(quiet=True)
         assert t > 1
 
+    def test_points_to_spheres_3D(self):
+        im=np.full((41,41,41),0)
+        im[20,20,20]=10
+        res_ps_3D=ps.tools.points_to_spheres(im=im)
+        assert np.sum(res_ps_3D) == 4169
+
+    def test_points_to_spheres_2D(self):
+        im=np.full((41,41),0)
+        im[20,20]=10
+        res_ps_2D=ps.tools.points_to_spheres(im=im)
+        assert np.sum(res_ps_2D) == 317
+
 
 if __name__ == '__main__':
     t = ToolsTest()