Migrate to native Python API for 'splashsurf'.

Author: duburcqa

genesis/engine/entities/particle_entity.py

@@ -247,7 +247,8 @@ def sample(self):
             self._vfaces = np.array([], dtype=gs.np_float)
         elif isinstance(self._morph, gs.options.morphs.MeshSet):
             for i in range(len(self._morph.files)):
-                pos_i, euler_i = map(np.asarray, (self._morph.poss[i], self._morph.eulers[i]))
+                pos_i = np.asarray(self._morph.poss[i], dtype=gs.np_float)
+                euler_i = np.asarray(self._morph.eulers[i], dtype=gs.np_float)
                 quat_i = gs.utils.geom.xyz_to_quat(euler_i, rpy=True, degrees=True)
                 self._vmesh[i].apply_transform(gu.trans_quat_to_T(pos_i, quat_i))
 
@@ -284,7 +285,8 @@ def sample(self):
             origin = np.mean(self._morph.poss, dtype=gs.np_float)
         else:
             # transform vmesh
-            pos, quat = map(np.asarray, (self._morph.pos, self._morph.quat))
+            pos = np.asarray(self._morph.pos, dtype=gs.np_float)
+            quat = np.asarray(self._morph.quat, dtype=gs.np_float)
             self._vmesh.apply_transform(gu.trans_quat_to_T(pos, quat))
             # transform particles
             particles = gu.transform_by_trans_quat(

genesis/engine/materials/MPM/elasto_plastic.py

@@ -12,7 +12,8 @@ class ElastoPlastic(Base):
 
     Note
     ----
-    Default yield ratio comes from the SNOW material in taichi's MPM implementation: https://github.com/taichi-dev/taichi_elements/blob/d19678869a28b09a32ef415b162e35dc929b792d/engine/mpm_solver.py#L434
+    Default yield ratio comes from the SNOW material in taichi's MPM implementation:
+    https://github.com/taichi-dev/taichi_elements/blob/d19678869a28b09a32ef415b162e35dc929b792d/engine/mpm_solver.py#L434
 
     Parameters
     ----------

genesis/options/surfaces.py

@@ -13,12 +13,16 @@
 class Surface(Options):
     """
     Base class for all surfaces types in Genesis.
-    A ``Surface`` object encapsulates all visual information used for rendering an entity or its sub-components (links, geoms, etc.)
-    The surface contains different types textures: diffuse_texture, specular_texture, roughness_texture, metallic_texture, transmission_texture, normal_texture, and emissive_texture. Each one of them is a `gs.textures.Texture` object.
+
+    A ``Surface`` object encapsulates all visual information used for rendering an entity or its sub-components (links,
+    geoms, ...). The surface contains different types textures: diffuse_texture, specular_texture, roughness_texture,
+    metallic_texture, transmission_texture, normal_texture, and emissive_texture. Each one of them is a
+    `gs.textures.Texture` object.
 
     Tip
     ---
-    If any of the textures only has single value (instead of a map), you can use the shortcut attribute (e.g., `color`, `roughness`, `metallic`, `emissive`) instead of creating a texture object.
+    If any of the textures only has single value (instead of a map), you can use the shortcut attribute (e.g., `color`,
+    `roughness`, `metallic`, `emissive`) instead of creating a texture object.
 
     Note
     ----
@@ -49,26 +53,26 @@ class Surface(Options):
     emissive_texture : gs.textures.Texture | None, optional
         Emissive texture of the surface.
     default_roughness : float, optional
-        Default roughness value when `roughness` is not set and the asset does not have a roughness texture. Defaults to 1.0.
+        Default roughness value when `roughness` is not set and the asset does not have a roughness texture. Defaults
+        to 1.0.
     vis_mode : str | None, optional
-        How the entity should be visualized. Possible values are ['visual', 'particle', 'collision', 'sdf', 'recon'].
-
+        How the entity should be visualized, e.g.
         - 'visual': Render the entity's visual geometry.
         - 'collision': Render the entity's collision geometry.
         - 'particle': Render the entity's particle representation (if applicable).
         - 'sdf': Render the reconstructed surface mesh of the entity's sdf.
         - 'recon': Render the reconstructed surface mesh of the entity's particle representation.
-
     smooth : bool, optional
         Whether to smooth face normals by interpolating vertex normals.
     double_sided : bool | None, optional
-        Whether to render both sides of the surface. Useful for non-watertight 2D objects. Defaults to True for Cloth material and False for others.
+        Whether to render both sides of the surface. Useful for non-watertight 2D objects. Defaults to True for Cloth
+        material and False for others.
     beam_angle : float
         The beam angle of emission. Defaults to 180.0.
     normal_diff_clamp : float, optional
         Controls the threshold for computing surface normals by interpolating vertex normals.
     recon_backend : str, optional
-        Backend for surface reconstruction. Possible values are ['splashsurf', 'openvdb'].
+        Backend for surface reconstruction. Possible values are ['splashsurf', 'openvdb']. Defaults to 'splashsurf'.
     generate_foam : bool, optional
         Whether to generate foam particles for visual effects for particle-based entities.
     foam_options : gs.options.FoamOptions, optional
@@ -92,14 +96,12 @@ class Surface(Options):
 
     default_roughness: float = 1.0
 
-    vis_mode: Optional[str] = None  # ['visual', 'particle', 'collision', 'sdf', 'recon']
-    smooth: bool = True  # whether to smooth face normals by interpolating vertex normals
-    double_sided: Optional[bool] = (
-        None  # whether to render both sides of the surface. Defaults to True for Cloth material and False for others.
-    )
+    vis_mode: Optional[str] = None
+    smooth: bool = True
+    double_sided: Optional[bool] = None
     beam_angle: float = 180
     normal_diff_clamp: float = 180
-    recon_backend: str = "splashsurf"  # backend for surface recon
+    recon_backend: str = "splashsurf"
     generate_foam: bool = False
     foam_options: Optional[FoamOptions] = None
 

genesis/utils/particle.py

@@ -6,6 +6,7 @@
 import tempfile
 
 import igl
+import pysplashsurf
 import numpy as np
 import trimesh
 
@@ -347,48 +348,26 @@ def parse_args(backend):
         return trimesh.Trimesh(vertices, faces, process=False)
 
     elif "splashsurf" in backend:
-        if gs.platform != "Linux":
-            gs.raise_exception("Backend 'splashsurf' is only supported on Linux.")
-
-        fd, xyz_path = tempfile.mkstemp(suffix=".xyz")
-        os.close(fd)
-        fd, obj_path = tempfile.mkstemp(suffix=".obj")
-        os.close(fd)
-        positions.astype(np.float32, copy=False).tofile(xyz_path)
-
         # Suggested value is 1.4-1.6, but 1.0 seems more detailed
-        radius_scale = args_dict.get("rscale", 1.0)
-        smooth_iter = args_dict.get("smooth")
-        r = radius * radius_scale
-
-        try:
-            command = ["splashsurf", "reconstruct", xyz_path, f"-r={r}", "-c=0.8", "-l=2.0", "-t=0.6", "-o", obj_path]
-            if smooth_iter is not None:
-                command += [
-                    "--mesh-cleanup=on",
-                    "--mesh-smoothing-weights=on",
-                    f"--mesh-smoothing-iters={int(smooth_iter)}",
-                    "--normals=on",
-                    "--normals-smoothing-iters=10",
-                ]
-
-            result = subprocess.run(map(str, command), capture_output=True, text=True)
-            if result.stdout:
-                gs.logger.debug(result.stdout)
-            if result.stderr:
-                gs.logger.warning(result.stderr)
-            if os.path.getsize(obj_path) == 0:
-                raise OSError("Output OBJ file is empty.")
-
-            # Read the generated OBJ file
-            mesh = trimesh.load_mesh(obj_path)
-            gs.logger.debug(f"[splashsurf]: reconstruct vertices: {mesh.vertices.shape}, {mesh.faces.shape}")
-        except OSError as e:
-            gs.raise_exception_from("Surface reconstruction failed.", e)
-        finally:
-            os.remove(xyz_path)
-            os.remove(obj_path)
-
+        mesh_with_data, _ = pysplashsurf.reconstruction_pipeline(
+            positions,
+            particle_radius=radius * args_dict.get("rscale", 1.0),
+            smoothing_length=2.0,
+            cube_size=0.8,
+            iso_surface_threshold=0.6,
+            mesh_smoothing_weights=True,
+            mesh_smoothing_iters=int(args_dict.get("smooth", 25)),
+            normals_smoothing_iters=10,
+            mesh_cleanup=True,
+            compute_normals=True,
+            enable_multi_threading=True,
+        )
+        mesh = trimesh.Trimesh(
+            vertices=mesh_with_data.mesh.vertices,
+            faces=mesh_with_data.mesh.triangles,
+            face_normals=mesh_with_data.get_point_attribute("normals"),
+        )
+        gs.logger.debug(f"[splashsurf]: reconstruct vertices: {mesh.vertices.shape}, {mesh.faces.shape}")
         return mesh
 
     else:

genesis/utils/path_planing.py

@@ -617,12 +617,11 @@ def plan(
 
 @ti.data_oriented
 class RRTConnect(PathPlanner):
-    def _init_rrt_connect_fields(self, goal_bias=0.1, max_nodes=4000, pos_tol=5e-3, max_step_size=0.05):
+    def _init_rrt_connect_fields(self, goal_bias=0.1, max_nodes=4000, max_step_size=0.05):
         self._is_rrt_connect_init = getattr(self, "_is_rrt_connect_init", False)
         if not self._is_rrt_connect_init:
             self._rrt_goal_bias = goal_bias
             self._rrt_max_nodes = max_nodes
-            self._rrt_pos_tol = pos_tol
             self._rrt_max_step_size = max_step_size
             self._rrt_start_configuration = ti.field(
                 dtype=gs.ti_float, shape=self._solver._batch_shape(self._entity.n_qs)

pyproject.toml

@@ -59,6 +59,8 @@ dependencies = [
     "gs-madrona>=0.0.2; platform_system == 'Linux' and (platform_machine == 'x86_64' or platform_machine == 'AMD64')",
     # Used for mesh simplification
     "fast_simplification>=0.1.12",
+    # Surface reconstruction library for particle data from SPH simulations
+    "pysplashsurf",
 ]
 
 [project.optional-dependencies]