From be8aa9e5e1129179ede7572ab528774273bceb26 Mon Sep 17 00:00:00 2001 From: Hadrien Oliveri Date: Tue, 22 Nov 2016 14:31:42 +0100 Subject: [PATCH 1/2] Added in compute_topomesh_vertex_property_from_faces: possibility to interpolate with angular sectors around a vertex. --- .../property_topomesh_analysis.py | 36 ++++++++++++++----- 1 file changed, 28 insertions(+), 8 deletions(-) diff --git a/src/openalea/cellcomplex/property_topomesh/property_topomesh_analysis.py b/src/openalea/cellcomplex/property_topomesh/property_topomesh_analysis.py index ac62577..d0a51de 100644 --- a/src/openalea/cellcomplex/property_topomesh/property_topomesh_analysis.py +++ b/src/openalea/cellcomplex/property_topomesh/property_topomesh_analysis.py @@ -1392,24 +1392,25 @@ def compute_topomesh_vertex_property_from_faces(topomesh,property_name,weighting * *area*: the weight on the faces is equal to their area * *angle*: the weight of the faces is equal to the incidence angle at the considered vertex (neighborhood=1) * *cotangent*: the weight of the faces is equal to the sum of cotangent of opposite angles (neighborhood=1) - neighborhood (int): + * *angular sector*: the weight of the faces is computed as the angular sector intersecting face barycenter (neighborhood=1) + neighborhood (int): The ring-distance of faces considered around each vertex (1-ring is immediate neighbor faces). - adjacency_sigma (float): + adjacency_sigma (float): The standard deviation of the gaussian weighting using the ring-distance. Returns: None - + Note: The PropertyTopomesh passed as argument is updated. - Example: + Example: >>> from openalea.cellcomplex.property_topomesh.example_topomesh import square_topomesh >>> from openalea.cellcomplex.property_topomesh.property_topomesh_analysis import compute_topomesh_property >>> from openalea.cellcomplex.property_topomesh.property_topomesh_analysis import compute_topomesh_vertex_property_from_faces >>> topomesh = square_topomesh(side_length=1) >>> compute_topomesh_property(topomesh,'normal',2,normal_method='orientation') - >>> compute_topomesh_vertex_property_from_faces(topomesh,'normal',weighting='area',neighborhood=3,adjacency_sigma=1.2) + >>> compute_topomesh_vertex_property_from_faces(topomesh,'normal',weighting='area',neighborhood=3,adjacency_sigma=1.2) """ @@ -1417,7 +1418,7 @@ def compute_topomesh_vertex_property_from_faces(topomesh,property_name,weighting print "--> Computing vertex property from faces" assert topomesh.has_wisp_property(property_name,2,is_computed=True) - assert weighting in ['uniform','area','angle','cotangent'] + assert weighting in ['uniform','area','angle','cotangent','angular sector'] if neighborhood > 1: try: assert weighting in ['uniform','area'] @@ -1436,7 +1437,7 @@ def compute_topomesh_vertex_property_from_faces(topomesh,property_name,weighting vertex_faces = np.concatenate([list(topomesh.regions(0,v,2)) for v in topomesh.wisps(0)]).astype(np.uint32) vertex_face_vertices = np.concatenate([v*np.ones_like(list(topomesh.regions(0,v,2))) for v in topomesh.wisps(0)]).astype(np.uint32) - + compute_topomesh_property(topomesh,'border_neighbors',2) vertex_face_adjacency_matrix[(vertex_face_vertices,vertex_faces)] = 0 @@ -1469,9 +1470,10 @@ def compute_topomesh_vertex_property_from_faces(topomesh,property_name,weighting compute_topomesh_property(topomesh,'vertices',2) compute_topomesh_property(topomesh,'angles',2) vertex_face_face_vertices = topomesh.wisp_property('vertices',2).values(vertex_faces) - + vertex_face_face_angles = topomesh.wisp_property('angles',2).values(vertex_faces) vertex_face_vertex_angles = np.array([angles[vertices == v] for v,vertices,angles in zip(vertex_face_vertices,vertex_face_face_vertices,vertex_face_face_angles)]) + print vertex_face_vertex_angles vertex_face_weight = vertex_face_vertex_angles[:,0] elif weighting == 'cotangent': compute_topomesh_property(topomesh,'vertices',2) @@ -1482,6 +1484,24 @@ def compute_topomesh_vertex_property_from_faces(topomesh,property_name,weighting vertex_face_opposite_angles = np.array([angles[vertices != v] for v,vertices,angles in zip(vertex_face_vertices,vertex_face_face_vertices,vertex_face_face_angles)]) vertex_face_cotangents = 1./np.tan(vertex_face_opposite_angles) vertex_face_weight = vertex_face_cotangents.sum(axis=1) + elif weighting == 'angular sector': + compute_topomesh_property(topomesh,'vertices',2) + compute_topomesh_property(topomesh,'angles',2) + compute_topomesh_property(topomesh,'barycenter',2) + compute_topomesh_property(topomesh,'barycenter',1) + + vertex_face_face_vertices = topomesh.wisp_property('vertices',2).values(vertex_faces) + vertex_face_face_angles = topomesh.wisp_property('angles',2).values(vertex_faces) + vertex_face_barycenter = topomesh.wisp_property('barycenter',2).values(vertex_faces) + vertex_positions = topomesh.wisp_property('barycenter',1).values(vertex_faces) + + face_vertex_positions = np.array([vertex_positions[v] for v in vertex_face_face_vertices]) + face_vertex_barycenter_radii = np.linalg.norm(np.array([[vertex_face_barycenter[f] - v for v in face_vertex_positions[f]] for f in vertex_faces]), axis=1) + + vertex_face_weight = np.array( + [r[vertices == v] * angles[vertices == v] for v, vertices, r, angles + in zip(vertex_face_vertices, vertex_face_face_vertices, + face_vertex_barycenter_radii, vertex_face_face_angles)])[:,0] if neighborhood > 1: vertex_face_weight = vertex_face_weight*vertex_adjacency_gaussian_weight From 6db7f33770443a3da80c9d044a8f7684e39d4fdf Mon Sep 17 00:00:00 2001 From: Guillaume Cerutti Date: Tue, 22 Nov 2016 15:11:08 +0100 Subject: [PATCH 2/2] debugged vertex property computation using angular sector weighting + added test for compute_vertex_property_from_faces --- .../property_topomesh/property_topomesh_analysis.py | 8 ++++---- test/test_topomesh_geometrical_properties.py | 11 +++++++++++ 2 files changed, 15 insertions(+), 4 deletions(-) diff --git a/src/openalea/cellcomplex/property_topomesh/property_topomesh_analysis.py b/src/openalea/cellcomplex/property_topomesh/property_topomesh_analysis.py index d0a51de..0e74b80 100644 --- a/src/openalea/cellcomplex/property_topomesh/property_topomesh_analysis.py +++ b/src/openalea/cellcomplex/property_topomesh/property_topomesh_analysis.py @@ -1488,15 +1488,15 @@ def compute_topomesh_vertex_property_from_faces(topomesh,property_name,weighting compute_topomesh_property(topomesh,'vertices',2) compute_topomesh_property(topomesh,'angles',2) compute_topomesh_property(topomesh,'barycenter',2) - compute_topomesh_property(topomesh,'barycenter',1) vertex_face_face_vertices = topomesh.wisp_property('vertices',2).values(vertex_faces) vertex_face_face_angles = topomesh.wisp_property('angles',2).values(vertex_faces) vertex_face_barycenter = topomesh.wisp_property('barycenter',2).values(vertex_faces) - vertex_positions = topomesh.wisp_property('barycenter',1).values(vertex_faces) + #vertex_positions = topomesh.wisp_property('barycenter',0).values(vertex_face_vertices) - face_vertex_positions = np.array([vertex_positions[v] for v in vertex_face_face_vertices]) - face_vertex_barycenter_radii = np.linalg.norm(np.array([[vertex_face_barycenter[f] - v for v in face_vertex_positions[f]] for f in vertex_faces]), axis=1) + face_vertex_positions = topomesh.wisp_property('barycenter',0).values(vertex_face_face_vertices) + # face_vertex_barycenter_radii = np.linalg.norm(np.array([[vertex_face_barycenter[f] - v for v in face_vertex_positions[f]] for f in vertex_faces]), axis=1) + face_vertex_barycenter_radii = np.linalg.norm(face_vertex_positions - vertex_face_barycenter[:,np.newaxis],axis=1) vertex_face_weight = np.array( [r[vertices == v] * angles[vertices == v] for v, vertices, r, angles diff --git a/test/test_topomesh_geometrical_properties.py b/test/test_topomesh_geometrical_properties.py index 1c438c0..5b8020b 100644 --- a/test/test_topomesh_geometrical_properties.py +++ b/test/test_topomesh_geometrical_properties.py @@ -127,3 +127,14 @@ def test_curvature_property(): assert np.all( np.isclose(topomesh.wisp_property('mean_curvature', 2).values(), 1. / radius, 1e-1)) + +def test_normal_property(): + radius = 10. + + topomesh = sphere_topomesh(radius) + compute_topomesh_property(topomesh, 'normal', 2, + normal_method='orientation') + + for weighting in ['uniform','area','angle','cotangent','angular sector']: + compute_topomesh_vertex_property_from_faces(topomesh, 'normal', weighting=weighting) +