[port] tested surf2vol and remeshsurf

Author: fangq

iso2mesh/core.py

@@ -602,8 +602,8 @@ def surf2volz(node, face, xi, yi, zi):
         enum = bcutedges.shape[0]
 
         for j in range(enum):
-            e0 = bcutpos[bcutedges[j, 0], :2]
-            e1 = bcutpos[bcutedges[j, 1], :2]
+            e0 = bcutpos[bcutedges[j, 0] - 1, :2]
+            e1 = bcutpos[bcutedges[j, 1] - 1, :2]
             length = np.ceil(np.sum(np.abs(e1 - e0)) / (np.abs(dx) + np.abs(dy))) + 1
             dd = (e1 - e0) / length
 
@@ -648,7 +648,7 @@ def surf2vol(node, face, xi, yi, zi, **kwargs):
 
     opt = kwargs
     label = opt.get("label", 0)
-    elabel = 1
+    elabel = np.array([1], dtype=np.int32)
 
     if face.shape[1] >= 4:
         elabel = np.unique(face[:, -1])
@@ -671,7 +671,9 @@ def surf2vol(node, face, xi, yi, zi, **kwargs):
 
         im = surf2volz(node[:, :3], fc[:, :3], xi, yi, zi)
         im |= np.moveaxis(surf2volz(node[:, [2, 0, 1]], fc[:, :3], zi, xi, yi), 0, 2)
-        im |= np.moveaxis(surf2volz(node[:, [1, 2, 0]], fc[:, :3], yi, zi, xi), 0, 1)
+        im |= np.moveaxis(
+            surf2volz(node[:, [1, 2, 0]], fc[:, :3], yi, zi, xi), [0, 1], [-2, -1]
+        )
 
         if opt.get("fill", 0) or label:
             im = binary_fill_holes(im)
@@ -680,15 +682,15 @@ def surf2vol(node, face, xi, yi, zi, **kwargs):
 
         img = np.maximum(im.astype(img.dtype), img)
 
-    v2smap = None
     if "v2smap" in kwargs:
         dlen = np.abs([xi[1] - xi[0], yi[1] - yi[0], zi[1] - zi[0]])
         offset = np.min(node, axis=0)
         v2smap = np.eye(4)
         v2smap[:3, :3] = np.diag(np.abs(dlen))
         v2smap[:3, 3] = offset
+        return img.astype(np.uint8), v2smap
 
-    return img, v2smap
+    return img.astype(np.uint8)
 
 
 def binsurface(img, nface=3):
@@ -1205,13 +1207,13 @@ def remeshsurf(node, face, opt):
         dx = opt / 4
 
     x_range = np.arange(p0[0] - dx, p1[0] + dx, dx)
-    y_range = np.arange(p0[1] - dx, p1[1] + dx)
-    z_range = np.arange(p0[2] - dx, p1[2] + dx)
+    y_range = np.arange(p0[1] - dx, p1[1] + dx, dx)
+    z_range = np.arange(p0[2] - dx, p1[2] + dx, dx)
 
     img = surf2vol(node, face, x_range, y_range, z_range)
 
     # Compute surface edges
-    eg = surfedge(face)
+    eg = surfedge(face)[0]
 
     closesize = 0
     if eg.size > 0 and isinstance(opt, dict):

iso2mesh/modify.py

@@ -208,11 +208,12 @@ def qmeshcut(elem, node, value, cutat):
         else:
             dist = eval(cutat[0]) - cutat[1]
         asign = np.where(dist <= 0, -1, 1)
-    elif not isinstance(cutat, (int, float)) and (len(cutat) == 9 or len(cutat) == 4):
-        if len(cutat) == 9:
-            a, b, c, d = getplanefrom3pt(np.array(cutat).reshape(3, 3))
+    elif not isinstance(cutat, (int, float)) and isinstance(cutat, (list, np.ndarray)):
+        cutat = np.array(cutat)
+        if cutat.size == 9:
+            a, b, c, d = getplanefrom3pt(cutat.reshape(3, 3))
         else:
-            a, b, c, d = cutat
+            a, b, c, d = cutat.tolist()
         dist = np.dot(node, np.array([a, b, c])) + d
         asign = np.where(dist >= 0, 1, -1)
     else:
@@ -261,7 +262,7 @@ def qmeshcut(elem, node, value, cutat):
 
     if value.shape[0] == node.shape[0]:
         if isinstance(cutat, (str, list)) or (
-            not isinstance(cutat, (int, float)) and len(cutat) in [4, 9]
+            not isinstance(cutat, (int, float)) and np.array(cutat).size in [4, 9]
         ):
             cutvalue = (
                 value[edges[cutedges, 0] - 1] * cutweight[:, [1]].squeeze()

test/run_test.py

@@ -1217,6 +1217,20 @@ def test_s2m_plc(self):
             sum(elemvolume(no[:, :3], fc[:, :3])), 412.8904758569056, 4
         )
 
+    def test_remeshsurf(self):
+        no, fc, _, _ = v2s(self.dist, 100, 2)
+        no1, fc1 = remeshsurf(no, fc, 1)
+        self.assertAlmostEqual(
+            sum(elemvolume(no[:, :3], fc[:, :3])) * 0.0001, 2.034550093898604, 4
+        )
+
+    # def test_surf2vol(self):
+    #     no, fc, _, _ = v2s(self.dist < 100, 0.5, 2)
+    #     vol = surf2vol(
+    #         no, fc, np.arange(19, 42, 1), np.arange(19, 42, 1), np.arange(19, 42, 1)
+    #     )
+    #     self.assertAlmostEqual(np.sum(vol.astype(np.float32)), 2.034550093898604, 4)
+
 
 @unittest.skipIf(
     (int(mpl_ver[0]), int(mpl_ver[1])) < (3, 6), "Requires Matplotlib 3.6 or higher"