11import numpy as np
22from compas .plugins import plugin
3- import sys
4- import platform
53
64try :
75 from compas_libigl import _meshing
6+
87 HAS_MESHING_MODULE = True
98except ImportError :
109 HAS_MESHING_MODULE = False
@@ -103,29 +102,29 @@ def trimesh_remesh_along_isolines(M, scalars, isovalues):
103102
104103def sequential_isoline_remeshing (V , F , S , isovalues ):
105104 """Pure Python implementation of sequential isoline remeshing.
106-
105+
107106 This is used as a fallback when the C++ implementation is not available.
108107 """
109108 current_V = V .copy ()
110109 current_F = F .copy ()
111110 current_S = S .copy ()
112-
111+
113112 # Sort isovalues to ensure consistent results
114113 sorted_isovalues = sorted (isovalues )
115-
114+
116115 for isovalue in sorted_isovalues :
117116 # Remesh along current isoline
118117 if HAS_MESHING_MODULE :
119118 result = _meshing .trimesh_remesh_along_isoline (current_V , current_F , current_S , isovalue )
120119 else :
121120 result = _pure_python_trimesh_remesh_along_isoline (current_V , current_F , current_S , isovalue )
122-
121+
123122 # Update the current mesh for the next iteration
124123 current_V = result [0 ]
125124 current_F = result [1 ]
126125 # We need to interpolate scalar values for new vertices
127126 current_S = _interpolate_vertex_values (current_V , result [0 ], current_S )
128-
127+
129128 return current_V , current_F , current_S
130129
131130
@@ -134,154 +133,138 @@ def _interpolate_vertex_values(old_vertices, new_vertices, old_values):
134133 # Simple implementation for fallback - in practice you'd want to use proper interpolation
135134 # based on the barycentric coordinates, but this is sufficient for a fallback
136135 new_values = np .zeros (len (new_vertices ))
137-
136+
138137 # For each new vertex, find the closest old vertex and use its value
139138 for i , new_v in enumerate (new_vertices ):
140- min_dist = float (' inf' )
139+ min_dist = float (" inf" )
141140 min_idx = 0
142141 for j , old_v in enumerate (old_vertices ):
143142 dist = np .sum ((new_v - old_v ) ** 2 )
144143 if dist < min_dist :
145144 min_dist = dist
146145 min_idx = j
147146 new_values [i ] = old_values [min_idx ]
148-
147+
149148 return new_values
150149
151150
152151def _pure_python_trimesh_remesh_along_isoline (V , F , S , isovalue ):
153152 """Pure Python implementation of trimesh_remesh_along_isoline.
154-
153+
155154 This is a simplified version that can be used when the C++ implementation is not available.
156155 It does not handle all edge cases correctly but provides basic functionality.
157156 """
158- from compas .geometry import intersect_line_plane
159157 from compas .datastructures import Mesh
160-
158+
161159 # Convert numpy arrays to lists
162160 vertices = V .tolist ()
163161 faces = F .tolist ()
164162 scalars = S .tolist ()
165-
163+
166164 # Create a mesh
167165 mesh = Mesh .from_vertices_and_faces (vertices , faces )
168-
166+
169167 # Assign scalar values to vertices
170168 for i , value in enumerate (scalars ):
171- mesh .vertex_attribute (i , ' scalar' , value )
172-
169+ mesh .vertex_attribute (i , " scalar" , value )
170+
173171 # New geometry to collect
174172 new_vertices = []
175173 new_faces = []
176174 labels = []
177-
178- # Map original vertices to new indices
179- vertex_map = {}
180-
181- # Create a plane from the isovalue
182- # Assuming isovalue is in the z-direction for simplicity
183- plane = ([0 , 0 , isovalue ], [0 , 0 , 1 ])
184-
175+
185176 # Process each original vertex
186177 for i in range (len (vertices )):
187- vertex = vertices [i ]
188178 scalar = scalars [i ]
189-
190- # Add to new vertices and update map
191- new_idx = len (new_vertices )
192- new_vertices .append (vertex )
193- vertex_map [i ] = new_idx
194-
179+
180+ # Add to new vertices
181+ new_vertices .append (vertices [i ])
182+
195183 # Determine label based on scalar value
196184 label = 0 if scalar < isovalue else 1
197185 labels .append (label )
198-
186+
199187 # Process each face to add intersection points
200188 edge_intersections = {}
201189 for face_idx , face in enumerate (faces ):
202- face_vertices = [vertices [i ] for i in face ]
203190 face_scalars = [scalars [i ] for i in face ]
204-
191+
205192 # Check if face crosses the isovalue
206193 crosses = False
207194 for i in range (len (face )):
208- if (face_scalars [i ] < isovalue ) != (face_scalars [(i + 1 ) % len (face )] < isovalue ):
195+ if (face_scalars [i ] < isovalue ) != (face_scalars [(i + 1 ) % len (face )] < isovalue ):
209196 crosses = True
210197 break
211-
198+
212199 if crosses :
213200 # For each edge, compute intersection if it crosses the isovalue
214201 for i in range (len (face )):
215202 v1_idx = face [i ]
216- v2_idx = face [(i + 1 ) % len (face )]
217-
203+ v2_idx = face [(i + 1 ) % len (face )]
204+
218205 # Ensure consistent ordering of edge
219206 edge = tuple (sorted ([v1_idx , v2_idx ]))
220-
207+
221208 # Skip if we've already processed this edge
222209 if edge in edge_intersections :
223210 continue
224-
211+
225212 s1 = scalars [v1_idx ]
226213 s2 = scalars [v2_idx ]
227-
214+
228215 # Check if edge crosses isovalue
229216 if (s1 < isovalue ) != (s2 < isovalue ):
230217 # Calculate intersection point
231218 t = (isovalue - s1 ) / (s2 - s1 )
232219 v1 = vertices [v1_idx ]
233220 v2 = vertices [v2_idx ]
234- intersect_point = [
235- v1 [0 ] + t * (v2 [0 ] - v1 [0 ]),
236- v1 [1 ] + t * (v2 [1 ] - v1 [1 ]),
237- v1 [2 ] + t * (v2 [2 ] - v1 [2 ])
238- ]
239-
221+ intersect_point = [v1 [0 ] + t * (v2 [0 ] - v1 [0 ]), v1 [1 ] + t * (v2 [1 ] - v1 [1 ]), v1 [2 ] + t * (v2 [2 ] - v1 [2 ])]
222+
240223 # Add intersection point to new vertices
241224 intersect_idx = len (new_vertices )
242225 new_vertices .append (intersect_point )
243226 edge_intersections [edge ] = intersect_idx
244-
227+
245228 # Assign label based on position (on the isovalue)
246229 labels .append (0.5 ) # 0.5 indicates on the boundary
247-
230+
248231 # Create new faces
249232 for face_idx , face in enumerate (faces ):
250233 # Extract edges and check if any cross the isovalue
251- edges = [(face [i ], face [(i + 1 ) % len (face )]) for i in range (len (face ))]
234+ edges = [(face [i ], face [(i + 1 ) % len (face )]) for i in range (len (face ))]
252235 sorted_edges = [tuple (sorted (edge )) for edge in edges ]
253236 crossing_edges = [edge for edge in sorted_edges if edge in edge_intersections ]
254-
237+
255238 if len (crossing_edges ) == 0 :
256239 # Face doesn't cross isovalue, add it directly
257- new_face = [vertex_map [ v ] for v in face ]
240+ new_face = [i for i in face ]
258241 new_faces .append (new_face )
259242 elif len (crossing_edges ) == 2 :
260243 # Face crosses isovalue, split into two new faces
261244 # Get intersection points
262245 p1_idx = edge_intersections [crossing_edges [0 ]]
263246 p2_idx = edge_intersections [crossing_edges [1 ]]
264-
247+
265248 # Determine which vertices are on which side
266- vertices_below = [vertex_map [ v ] for v in face if scalars [v ] < isovalue ]
267- vertices_above = [vertex_map [ v ] for v in face if scalars [v ] >= isovalue ]
268-
249+ vertices_below = [i for i in face if scalars [i ] < isovalue ]
250+ vertices_above = [i for i in face if scalars [i ] >= isovalue ]
251+
269252 # Create new faces
270253 if vertices_below :
271254 for v in vertices_below :
272255 new_faces .append ([v , p1_idx , p2_idx ])
273-
256+
274257 if vertices_above :
275258 for v in vertices_above :
276259 new_faces .append ([v , p1_idx , p2_idx ])
277-
260+
278261 return np .array (new_vertices ), np .array (new_faces ), np .array (labels )
279262
280263
281264# Function to split a mesh by a plane
282265def split_mesh_by_plane (mesh , plane_point , plane_normal ):
283266 """Split a mesh into two parts based on a plane.
284-
267+
285268 Parameters
286269 ----------
287270 mesh : compas.datastructures.Mesh
@@ -290,33 +273,33 @@ def split_mesh_by_plane(mesh, plane_point, plane_normal):
290273 A point on the splitting plane.
291274 plane_normal : list
292275 The normal vector of the splitting plane.
293-
276+
294277 Returns
295278 -------
296279 tuple
297280 A tuple containing two meshes (above_mesh, below_mesh).
298281 """
299282 # Normalize the plane normal
300- import numpy as np
301283 from compas .datastructures import Mesh
302- from compas .geometry import normalize_vector , dot_vectors , distance_point_plane
303-
284+ from compas .geometry import distance_point_plane
285+ from compas .geometry import normalize_vector
286+
304287 plane_normal = normalize_vector (plane_normal )
305288 plane = (plane_point , plane_normal )
306-
289+
307290 # Generate scalar field based on signed distance to plane
308- vertices = mesh .vertices_attributes (' xyz' )
291+ vertices = mesh .vertices_attributes (" xyz" )
309292 scalars = [distance_point_plane ((x , y , z ), plane ) for x , y , z in vertices ]
310-
293+
311294 # Use isoline remeshing at the zero distance (plane intersection)
312295 M = mesh .to_vertices_and_faces ()
313296 V2 , F2 , L = trimesh_remesh_along_isoline (M , scalars , 0.0 )
314-
297+
315298 # Create meshes for positive and negative regions
316299 pos_faces = [f for i , f in enumerate (F2 ) if all (L [v ] >= 0 for v in f )]
317300 neg_faces = [f for i , f in enumerate (F2 ) if all (L [v ] <= 0 for v in f )]
318-
301+
319302 above_mesh = Mesh .from_vertices_and_faces (V2 , pos_faces )
320303 below_mesh = Mesh .from_vertices_and_faces (V2 , neg_faces )
321-
304+
322305 return above_mesh , below_mesh
0 commit comments