MEPP2 Project
curvature.hpp
Go to the documentation of this file.
1 // Copyright (c) 2012-2019 University of Lyon and CNRS (France).
2 // All rights reserved.
3 //
4 // This file is part of MEPP2; you can redistribute it and/or modify
5 // it under the terms of the GNU Lesser General Public License as
6 // published by the Free Software Foundation; either version 3 of
7 // the License, or (at your option) any later version.
8 //
9 // This file is provided AS IS with NO WARRANTY OF ANY KIND, INCLUDING THE
10 // WARRANTY OF DESIGN, MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE.
11 #pragma once
12 
13 #include <set>
14 #include <stack>
15 #include <cmath> // fabs
16 
17 #include <Eigen/Dense> // for Matrix::Identity()
18 
19 #include <boost/graph/graph_traits.hpp>
20 #include <boost/graph/properties.hpp>
21 
22 #include <CGAL/boost/graph/iterator.h> // for circulators
23 #include <CGAL/boost/graph/helpers.h> // for is_border_edge()
24 
25 #include "eigen_val_vect.hpp"
26 
28 
32 
34 
35 namespace FEVV {
36 namespace Filters {
37 
38 template< typename HalfedgeGraph,
39  typename GeometryTraits = FEVV::Geometry_traits< HalfedgeGraph > >
40 struct v_Curv
41 {
42  double KminCurv;
43  double KmaxCurv;
44 
45  typename GeometryTraits::Vector
47  typename GeometryTraits::Vector
49 };
50 
51 template< typename HalfedgeGraph,
52  typename GeometryTraits = FEVV::Geometry_traits< HalfedgeGraph > >
53 double
55 {
56  return std::max(::fabs(input.KminCurv), ::fabs(input.KmaxCurv));
57 }
58 
59 template< typename HalfedgeGraph,
60  typename PointMap,
61  typename GeometryTraits = FEVV::Geometry_traits< HalfedgeGraph > >
62 double
63 triangle_area(const HalfedgeGraph &g,
64  const PointMap &pm,
65  typename boost::graph_traits< HalfedgeGraph >::face_descriptor fd,
66  const GeometryTraits &gt)
67 {
68  typedef typename boost::property_traits< PointMap >::value_type Point3d;
69  //using Vector = typename GeometryTraits::Vector;
70 
71  CGAL::Halfedge_around_face_circulator< HalfedgeGraph > p_halfedge(
72  halfedge(fd, g), g);
73  Point3d p = get(pm, target(*p_halfedge, g));
74  Point3d q = get(pm, target(next(*p_halfedge, g), g));
75  Point3d r = get(pm, target(next(next(*p_halfedge, g), g), g));
76 
77  return FEVV::Operators::Geometry::triangle_area< GeometryTraits >(
78  p, q, r, gt);
79 }
80 
81 //********************************************
82 // 1-ring neighborhood curvature for a vertex
83 //********************************************
84 template< typename HalfedgeGraph,
85  typename PointMap,
86  typename FaceNormalMap,
87  typename GeometryTraits = FEVV::Geometry_traits< HalfedgeGraph >,
88  typename GraphTraits = boost::graph_traits< HalfedgeGraph >,
89  typename VertexDescriptor = typename GraphTraits::vertex_descriptor >
90 void
91 principal_curvature_per_vert(const HalfedgeGraph &g,
92  const PointMap &pm,
93  const FaceNormalMap &f_nm,
94  VertexDescriptor vd,
95  double pp_matrix_sum[3][3],
96  const GeometryTraits &gt)
97 {
98  typedef typename GeometryTraits::Point Point3d;
99  using Vector = typename GeometryTraits::Vector;
100  typedef typename GraphTraits::face_descriptor face_descriptor;
101 
102  double area = 0;
103 
104  // iterate over all edges
105  CGAL::Halfedge_around_target_circulator< HalfedgeGraph > p_halfedge(vd, g),
106  done(p_halfedge);
107  do
108  {
109  // build edge vector and comput its norm
110  Point3d p1 = get(pm, target(*p_halfedge, g));
111  Point3d p2 = get(pm, source(*p_halfedge, g));
112 #if 1 // WAITING for AIF stuff - TODO : REMOVE THIS BLOC LATER...
113  Vector edge = gt.sub_p(p1, p2);
114 #else
115  Vector edge(p1[0] - p2[0], p1[1] - p2[1], p1[2] - p2[2]);
116 #endif
117  double len_edge = gt.length(edge);
118  if(len_edge == 0) // avoid divide by zero
119  continue;
120 
121  // compute (signed) angle between two incident faces, if exists
122  face_descriptor p_facet1 = face(*p_halfedge, g);
123  face_descriptor p_facet2 = face(opposite(*p_halfedge, g), g);
124  assert(p_facet1 != p_facet2);
125  if(p_facet1 == GraphTraits::null_face() ||
126  p_facet2 == GraphTraits::null_face())
127  continue; // border edge
128 
129  area += triangle_area(g, pm, p_facet1, gt);
130 
131  Vector normal1 = f_nm[p_facet1];
132  Vector normal2 = f_nm[p_facet2];
133 
134  // ---
135  Vector cp(normal1[1] * normal2[2] - normal1[2] * normal2[1],
136  normal1[2] * normal2[0] - normal1[0] * normal2[2],
137  normal1[0] * normal2[1] - normal1[1] * normal2[0]);
138  double ps_cp_xedge = cp[0] * edge[0] + cp[1] * edge[1] + cp[2] * edge[2];
139  // ---
140  double sine =
141  ps_cp_xedge /
142  len_edge; // = (CGAL::cross_product(normal1,normal2)*edge)/len_edge; //
143  // TODO BETTER
144  double beta = FEVV::Operators::Geometry::asin(sine);
145 
146  // compute edge * edge^t * coeff, and add it to current matrix
147  double p_vector_edge[3] = {edge[0], edge[1], edge[2]};
148  double pp_matrix[3][3] = {{0, 0, 0}, {0, 0, 0}, {0, 0, 0}};
150  p_vector_edge, pp_matrix, beta / len_edge);
151  FEVV::Math::Matrix::Square::add(pp_matrix, pp_matrix_sum);
152  } while(++p_halfedge != done);
153 
154  for(int i = 0; i < 3; i++)
155  for(int j = 0; j < 3; j++)
156  pp_matrix_sum[i][j] /= area;
157 }
158 
159 //**********************************************
160 // geodesic neighborhood curvature for a vertex
161 //**********************************************
162 template< typename HalfedgeGraph,
163  typename PointMap,
164  typename FaceNormalMap,
165  typename GeometryTraits = FEVV::Geometry_traits< HalfedgeGraph >,
166  typename GraphTraits = boost::graph_traits< HalfedgeGraph >,
167  typename VertexDescriptor = typename GraphTraits::vertex_descriptor >
168 void
169 geodes_principal_curvature_per_vert(const HalfedgeGraph &g,
170  const PointMap &pm,
171  const FaceNormalMap &f_nm,
172  VertexDescriptor vd,
173  double pp_matrix_sum[3][3],
174  double radius,
175  const GeometryTraits &gt)
176 {
177  typedef typename GeometryTraits::Point Point3d;
178  using Vector = typename GeometryTraits::Vector;
179  typedef typename GraphTraits::face_descriptor face_descriptor;
180 
181  std::set< VertexDescriptor > vertices;
182  Point3d o = get(pm, vd);
183  std::stack< VertexDescriptor > s;
184  s.push(vd);
185  vertices.insert(vd);
186  int iter = 0;
187  while(!s.empty())
188  {
189  VertexDescriptor v = s.top();
190  s.pop();
191  Point3d p = get(pm, v);
192  CGAL::Halfedge_around_target_circulator< HalfedgeGraph > h(v, g), done(h);
193  do
194  {
195  Point3d p1 = get(pm, target(*h, g));
196  Point3d p2 = get(pm, source(*h, g));
197 
198  // ---
199 #if 1 // WAITING for AIF stuff - TODO : REMOVE THIS BLOC LATER...
200  Vector vec = gt.sub_p(p2, p1);
201  Vector pm_o = gt.sub_p(p, o);
202 #else
203  Vector vec(p2[0] - p1[0], p2[1] - p1[1], p2[2] - p1[2]);
204  Vector PmO(P[0] - O[0], P[1] - O[1], P[2] - O[2]);
205 #endif
206  double ps_vx_pm_o =
207  vec[0] * pm_o[0] + vec[1] * pm_o[1] + vec[2] * pm_o[2];
208 
209  // ---
210  if(v == vd ||
211  ps_vx_pm_o >
212  0.0) // if (v==pVertex || vec * (P - O) > 0.0) // TODO BETTER
213  {
214  bool isect =
215  FEVV::Operators::Geometry::sphere_clip_vector(o, radius, p, vec, gt);
216 
217  if(!CGAL::is_border_edge(*h, g))
218  {
219  double len_edge = std::sqrt(
220  vec[0] * vec[0] + vec[1] * vec[1] +
221  vec[2] * vec[2]); // = std::sqrt(vec*vec); // TODO BETTER
222  // compute (signed) angle between two incident faces, if exists
223  face_descriptor p_facet1 = face(*h, g);
224  face_descriptor p_facet2 = face(opposite(*h, g), g);
225  assert(p_facet1 != p_facet2);
226  if(p_facet1 == GraphTraits::null_face() ||
227  p_facet2 == GraphTraits::null_face())
228  continue; // border edge
229 
230  Vector normal1 = f_nm[p_facet1];
231  Vector normal2 = f_nm[p_facet2];
232 
233  // ---
234  Vector cp(normal1[1] * normal2[2] - normal1[2] * normal2[1],
235  normal1[2] * normal2[0] - normal1[0] * normal2[2],
236  normal1[0] * normal2[1] - normal1[1] * normal2[0]);
237  double ps_cpx_v = cp[0] * vec[0] + cp[1] * vec[1] + cp[2] * vec[2];
238  // ---
239  double sine =
240  ps_cpx_v / len_edge; // = (CGAL::cross_product(normal1,
241  // normal2)*vec) / len_edge; // TODO BETTER
242  double beta = FEVV::Operators::Geometry::asin(sine);
243 
244  // compute edge * edge^t * coeff, and add it to current matrix
245  double p_vector_edge[3] = {vec[0], vec[1], vec[2]};
246  double pp_matrix[3][3] = {{0, 0, 0}, {0, 0, 0}, {0, 0, 0}};
248  p_vector_edge, pp_matrix, beta / len_edge);
249  FEVV::Math::Matrix::Square::add(pp_matrix, pp_matrix_sum);
250  }
251 
252  if(!isect)
253  {
254  VertexDescriptor w = source(
255  *h, g); // previously : iterator here, NO !!! -> error from GL ???
256  // // Vertex_iterator w=h->opposite()->vertex();
257 
258  assert(w != GraphTraits::null_vertex());
259 
260  if(vertices.find(w) == vertices.end())
261  {
262  vertices.insert(w);
263  s.push(w);
264  }
265  }
266  }
267  } while(++h != done);
268 
269  iter++;
270  }
271 
272  double area = M_PI * radius * radius;
273 
274  for(int i = 0; i < 3; i++)
275  for(int j = 0; j < 3; j++)
276  pp_matrix_sum[i][j] /= area;
277 }
278 
279 
280 //#define DBG_calculate_curvature
281 
282 
306 template< typename HalfedgeGraph,
307  typename VertexCurvatureMap,
308  typename PointMap,
309  typename FaceNormalMap,
310  typename GeometryTraits = FEVV::Geometry_traits< HalfedgeGraph > >
311 void
312 calculate_curvature(const HalfedgeGraph &g,
313  VertexCurvatureMap &v_cm,
314  const PointMap &pm,
315  const FaceNormalMap &f_nm,
316  bool is_geod,
317  double radius,
318  double &min_nrm_min_curvature,
319  double &max_nrm_min_curvature,
320  double &min_nrm_max_curvature,
321  double &max_nrm_max_curvature,
322  const GeometryTraits &gt)
323 {
324  typedef boost::graph_traits< HalfedgeGraph > GraphTraits;
325  typedef typename GraphTraits::vertex_iterator vertex_iterator;
326 
327  using Vector = typename GeometryTraits::Vector;
328 
329  min_nrm_min_curvature = 100000;
330  max_nrm_min_curvature = -100000;
331 
332  min_nrm_max_curvature = 100000;
333  max_nrm_max_curvature = -100000;
334 
335  // ---
336 
337  std::cout << "Asking to calculate curvature" << std::endl;
338 #ifdef DBG_calculate_curvature
339  std::cout << "Real radius = " << radius << std::endl;
340 #endif
341 
342  unsigned int skipped = 0; // number of skipped vertices
343 
344  const auto iterator_pair =
345  vertices(g); // vertices() returns a vertex_iterator pair
346  vertex_iterator vi = iterator_pair.first;
347  vertex_iterator vi_end = iterator_pair.second;
348  for(; vi != vi_end; ++vi)
349  {
350  // skip isolated vertex
351  if(halfedge(*vi, g) == GraphTraits::null_halfedge())
352  {
353  v_cm[*vi].VKmaxCurv = gt.NULL_VECTOR;
354  v_cm[*vi].VKminCurv = gt.NULL_VECTOR;
355  v_cm[*vi].KmaxCurv = 0;
356  v_cm[*vi].KminCurv = 0;
357  skipped++;
358  continue;
359  }
360 
361  bool no_val_pro = false;
362 
363  double pp_matrix_sum[3][3] = {{0, 0, 0}, {0, 0, 0}, {0, 0, 0}};
364 
365  if(is_geod == true) // geodesic neighborhood
367  HalfedgeGraph,
368  PointMap,
369  FaceNormalMap,
370  GeometryTraits,
371  GraphTraits,
372  typename GraphTraits::vertex_descriptor >( g,
373  pm,
374  f_nm,
375  *vi,
376  pp_matrix_sum,
377  radius,
378  gt);
379  else // 1-ring neighborhood
381  HalfedgeGraph,
382  PointMap,
383  FaceNormalMap,
384  GeometryTraits,
385  GraphTraits,
386  typename GraphTraits::vertex_descriptor >( g,
387  pm,
388  f_nm,
389  *vi,
390  pp_matrix_sum,
391  gt);
392 
393  // Eigen values/vectors
394  EigenMatrix cov_mat;
395  EigenvectorsType vect_pro;
396  EigenvalueType valpro;
397 
398  cov_mat(0, 0) = pp_matrix_sum[0][0];
399  cov_mat(0, 1) = pp_matrix_sum[0][1];
400  cov_mat(0, 2) = pp_matrix_sum[0][2];
401  cov_mat(1, 0) = pp_matrix_sum[1][0];
402  cov_mat(1, 1) = pp_matrix_sum[1][1];
403  cov_mat(1, 2) = pp_matrix_sum[1][2];
404  cov_mat(2, 0) = pp_matrix_sum[2][0];
405  cov_mat(2, 1) = pp_matrix_sum[2][1];
406  cov_mat(2, 2) = pp_matrix_sum[2][2];
407 
408 #ifdef DBG_calculate_curvature
409  std::cout << "\n\nCovMat\n\n" << CovMat << std::endl;
410 #endif
411 
412  if(pp_matrix_sum[0][1] == 0 && pp_matrix_sum[0][2] == 0 &&
413  pp_matrix_sum[1][0] == 0 && pp_matrix_sum[1][2] == 0 &&
414  pp_matrix_sum[2][1] == 0 && pp_matrix_sum[2][0] == 0)
415  {
416  for(int i = 0; i < 3; i++)
417  valpro(i) = cov_mat(i, i);
418  vect_pro = EigenMatrix::Identity();
419  }
420  else
421  {
422  if(eigen_val_vect_compute(cov_mat, vect_pro, valpro) == -1)
423  {
424  no_val_pro = true;
425  }
426  }
427 
428 #ifdef DBG_calculate_curvature
429  std::cout << "\n\nValpro\n\n" << Valpro << std::endl;
430  std::cout << "\n\nVectPro\n\n" << VectPro << std::endl;
431 #endif
432 
433  if(!no_val_pro)
434  {
435  eigen_val_vect_sort(valpro, vect_pro);
436  Vector v_kmax_curv(vect_pro(0, 1), vect_pro(1, 1), vect_pro(2, 1));
437  Vector v_kmin_curv(vect_pro(0, 0), vect_pro(1, 0), vect_pro(2, 0));
438 
439 #ifdef DBG_calculate_curvature
440  std::cout << "\n\nValpro sorted\n\n" << Valpro << std::endl;
441  std::cout << "\n\nVectPro sorted\n\n" << VectPro << std::endl;
442  std::cout << "\nVKmaxCurv = " << VKmaxCurv[0] << " " << VKmaxCurv[1]
443  << " " << VKmaxCurv[2] << std::endl;
444  std::cout << "\nVKminCurv = " << VKminCurv[0] << " " << VKminCurv[1]
445  << " " << VKminCurv[2] << std::endl;
446 #endif
447 
448  v_cm[*vi].VKmaxCurv = v_kmax_curv;
449  v_cm[*vi].VKminCurv = v_kmin_curv;
450 
451  if(valpro(0) > 0)
452  {
453  v_cm[*vi].KmaxCurv = valpro(0);
454  v_cm[*vi].KminCurv = valpro(1);
455  }
456  else
457  {
458  v_cm[*vi].KmaxCurv = -valpro(0);
459  v_cm[*vi].KminCurv = -valpro(1);
460  }
461  }
462  else
463  {
464  v_cm[*vi].VKmaxCurv = gt.NULL_VECTOR;
465  v_cm[*vi].VKminCurv = gt.NULL_VECTOR;
466  v_cm[*vi].KmaxCurv = 0;
467  v_cm[*vi].KminCurv = 0;
468  }
469 
470 #ifdef DBG_calculate_curvature
471  std::cout << "\nv_cm[*vi].KmaxCurv = " << v_cm[*vi].KmaxCurv << std::endl;
472  std::cout << "\nv_cm[*vi].KminCurv = " << v_cm[*vi].KminCurv << std::endl;
473 #endif
474 
475  min_nrm_min_curvature =
476  std::min< double >(min_nrm_min_curvature, v_cm[*vi].KminCurv);
477  max_nrm_min_curvature =
478  std::max< double >(max_nrm_min_curvature, v_cm[*vi].KminCurv);
479 
480  min_nrm_max_curvature =
481  std::min< double >(min_nrm_max_curvature, v_cm[*vi].KmaxCurv);
482  max_nrm_max_curvature =
483  std::max< double >(max_nrm_max_curvature, v_cm[*vi].KmaxCurv);
484  }
485 
486  if(skipped > 0)
487  {
488  std::cout << "Curvature warning: "
489  << skipped << " isolated vertices were skipped."
490  << std::endl;
491  }
492 
493  // ---
494 
495  /*
496  int cpt = 0;
497  auto iterator_pair_log = vertices(g); // vertices() returns a vertex_iterator
498  pair vertex_iterator vi_log = iterator_pair_log.first; vertex_iterator
499  vend_log = iterator_pair_log.second; for (; vi_log != vend_log; ++vi_log)
500  {
501  printf("cpt: %4d - Kmin: %.3lf - Kmax: %.3lf | VKmin: (%.3lf, %.3lf,
502  %.3lf) | VKmax: (%.3lf, %.3lf, %.3lf)\n", cpt, v_cm[*vi_log].KminCurv,
503  v_cm[*vi_log].KmaxCurv, v_cm[*vi_log].VKminCurv[0],
504  v_cm[*vi_log].VKminCurv[1], v_cm[*vi_log].VKminCurv[2],
505  v_cm[*vi_log].VKmaxCurv[0], v_cm[*vi_log].VKmaxCurv[1],
506  v_cm[*vi_log].VKmaxCurv[2]); cpt++;
507  }
508  */
509 }
510 
532 template< typename HalfedgeGraph,
533  typename VertexCurvatureMap,
534  typename PointMap,
535  typename FaceNormalMap,
536  typename GeometryTraits = FEVV::Geometry_traits< HalfedgeGraph > >
537 void
538 calculate_curvature(const HalfedgeGraph &g,
539  VertexCurvatureMap &v_cm,
540  const PointMap &pm,
541  const FaceNormalMap &f_nm,
542  bool is_geod,
543  double radius,
544  double &min_nrm_min_curvature,
545  double &max_nrm_min_curvature,
546  double &min_nrm_max_curvature,
547  double &max_nrm_max_curvature)
548 {
549  GeometryTraits gt(g);
550  calculate_curvature< HalfedgeGraph,
551  VertexCurvatureMap,
552  PointMap,
553  FaceNormalMap,
554  GeometryTraits >(g,
555  v_cm,
556  pm,
557  f_nm,
558  is_geod,
559  radius,
560  min_nrm_min_curvature,
561  max_nrm_min_curvature,
562  min_nrm_max_curvature,
563  max_nrm_max_curvature,
564  gt);
565 }
566 
567 } // namespace Filters
568 } // namespace FEVV
FEVV::DataStructures::AIF::vertices
std::pair< typename boost::graph_traits< FEVV::DataStructures::AIF::AIFMesh >::vertex_iterator, typename boost::graph_traits< FEVV::DataStructures::AIF::AIFMesh >::vertex_iterator > vertices(const FEVV::DataStructures::AIF::AIFMesh &sm)
Returns the iterator range of the vertices of the mesh.
Definition: Graph_traits_aif.h:172
EigenMatrix
Eigen::Matrix3d EigenMatrix
Definition: eigen_val_vect.hpp:19
EigenvectorsType
Eigen::Matrix3d EigenvectorsType
Definition: eigen_val_vect.hpp:20
EigenvalueType
Eigen::Vector3d EigenvalueType
Definition: eigen_val_vect.hpp:21
eigen_val_vect_sort
void eigen_val_vect_sort(EigenvalueType &ad, EigenvectorsType &u)
Definition: eigen_val_vect.hpp:62
FEVV::DataStructures::AIF::next
boost::graph_traits< FEVV::DataStructures::AIF::AIFMesh >::halfedge_descriptor next(typename boost::graph_traits< FEVV::DataStructures::AIF::AIFMesh >::halfedge_descriptor h, const FEVV::DataStructures::AIF::AIFMesh &sm)
Returns the next halfedge around its face.
Definition: Graph_traits_aif.h:599
FEVV::Filters::v_Curv::KmaxCurv
double KmaxCurv
The minimum curvature.
Definition: curvature.hpp:43
Vector
AIFMesh::Vector Vector
Definition: Graph_properties_aif.h:22
FEVV::Math::Matrix::Square::vector_times_transpose_mult
static void vector_times_transpose_mult(const ElementType *p_vector, ElementType pp_matrix[][DIM], ElementType coeff)
Compute coef * V x V^t (V is a vertical DIM-D vector)
Definition: MatrixOperations.hpp:976
FEVV::Filters::v_Curv
Definition: curvature.hpp:41
eigen_val_vect.hpp
FEVV::DataStructures::AIF::opposite
boost::graph_traits< FEVV::DataStructures::AIF::AIFMesh >::halfedge_descriptor opposite(typename boost::graph_traits< FEVV::DataStructures::AIF::AIFMesh >::halfedge_descriptor h, const FEVV::DataStructures::AIF::AIFMesh &sm)
Returns the halfedge with source and target swapped.
Definition: Graph_traits_aif.h:625
FEVV::Geometry_traits< HalfedgeGraph >
Point
AIFMesh::Point Point
Definition: Graph_properties_aif.h:21
FEVV::Filters::v_Curv::KminCurv
double KminCurv
The minimum curvature.
Definition: curvature.hpp:42
FEVV::DataStructures::AIF::source
boost::graph_traits< FEVV::DataStructures::AIF::AIFMesh >::vertex_descriptor source(typename boost::graph_traits< FEVV::DataStructures::AIF::AIFMesh >::edge_descriptor e, const FEVV::DataStructures::AIF::AIFMesh &)
Returns the source vertex of e.
Definition: Graph_traits_aif.h:387
FEVV::DataStructures::AIF::edge
boost::graph_traits< FEVV::DataStructures::AIF::AIFMesh >::edge_descriptor edge(typename boost::graph_traits< FEVV::DataStructures::AIF::AIFMesh >::halfedge_descriptor h, const FEVV::DataStructures::AIF::AIFMesh &sm)
Returns the edge corresponding to h and opposite(h).
Definition: Graph_traits_aif.h:345
FEVV::Filters::calculate_curvature
void calculate_curvature(const HalfedgeGraph &g, VertexCurvatureMap &v_cm, const PointMap &pm, const FaceNormalMap &f_nm, bool is_geod, double radius, double &min_nrm_min_curvature, double &max_nrm_min_curvature, double &min_nrm_max_curvature, double &max_nrm_max_curvature, const GeometryTraits &gt)
Calculate the curvature for a mesh.
Definition: curvature.hpp:312
FEVV::Operators::Geometry::sphere_clip_vector
static bool sphere_clip_vector(const typename GeometryTraits::Point &center, double r, const typename GeometryTraits::Point &p, typename GeometryTraits::Vector &v, const GeometryTraits &gt)
Compute the intersection of a sphere (center + radius) with a ray/line (starting point + direction ve...
Definition: ClippingAndIntersection.hpp:39
FEVV::get
FEVV::PCLPointCloudPointMap::value_type get(const FEVV::PCLPointCloudPointMap &pm, FEVV::PCLPointCloudPointMap::key_type key)
Specialization of get(point_map, key) for PCLPointCloud.
Definition: Graph_properties_pcl_point_cloud.h:117
FEVV
Interfaces for plugins These interfaces will be used for different plugins.
Definition: Assert.h:16
AngleOperations.hpp
FEVV::DataStructures::AIF::halfedge
boost::graph_traits< FEVV::DataStructures::AIF::AIFMesh >::halfedge_descriptor halfedge(typename boost::graph_traits< FEVV::DataStructures::AIF::AIFMesh >::vertex_descriptor v, const FEVV::DataStructures::AIF::AIFMesh &sm)
Returns a halfedge with target v.
Definition: Graph_traits_aif.h:296
FEVV::DataStructures::AIF::target
boost::graph_traits< FEVV::DataStructures::AIF::AIFMesh >::vertex_descriptor target(typename boost::graph_traits< FEVV::DataStructures::AIF::AIFMesh >::edge_descriptor e, const FEVV::DataStructures::AIF::AIFMesh &)
Returns the target vertex of e.
Definition: Graph_traits_aif.h:400
Geometry_traits.h
triangles.hpp
FEVV::DataStructures::AIF::AIFVector
Definition: AIFProperties.h:173
ClippingAndIntersection.hpp
FEVV::Filters::v_Curv::VKminCurv
GeometryTraits::Vector VKminCurv
The minimum curvature direction.
Definition: curvature.hpp:46
FEVV::Filters::triangle_area
double triangle_area(const HalfedgeGraph &g, const PointMap &pm, typename boost::graph_traits< HalfedgeGraph >::face_descriptor fd, const GeometryTraits &gt)
Definition: curvature.hpp:63
msdm2::vertex_descriptor
boost::graph_traits< MeshT >::vertex_descriptor vertex_descriptor
Definition: msdm2_surfacemesh.h:33
FEVV::Filters::fabs
double fabs(const v_Curv< HalfedgeGraph > &input)
Definition: curvature.hpp:54
FEVV::DataStructures::AIF::face
boost::graph_traits< FEVV::DataStructures::AIF::AIFMesh >::face_descriptor face(typename boost::graph_traits< FEVV::DataStructures::AIF::AIFMesh >::halfedge_descriptor h, const FEVV::DataStructures::AIF::AIFMesh &)
Returns the face incident to halfedge h.
Definition: Graph_traits_aif.h:664
FEVV::Filters::geodes_principal_curvature_per_vert
void geodes_principal_curvature_per_vert(const HalfedgeGraph &g, const PointMap &pm, const FaceNormalMap &f_nm, VertexDescriptor vd, double pp_matrix_sum[3][3], double radius, const GeometryTraits &gt)
Definition: curvature.hpp:169
Point3d
EnrichedPolyhedron::Point_3 Point3d
Definition: boolops_cpolyhedron_builder.hpp:27
FEVV::Filters::v_Curv::VKmaxCurv
GeometryTraits::Vector VKmaxCurv
The minimum curvature direction.
Definition: curvature.hpp:48
eigen_val_vect_compute
int eigen_val_vect_compute(const EigenMatrix &a, EigenvectorsType &u, EigenvalueType &ad)
Definition: eigen_val_vect.hpp:33
FEVV::Operators::Geometry::asin
static T asin(T sine)
Safe call to the std::asin function.
Definition: AngleOperations.hpp:31
MatrixOperations.hpp
FEVV::Filters::principal_curvature_per_vert
void principal_curvature_per_vert(const HalfedgeGraph &g, const PointMap &pm, const FaceNormalMap &f_nm, VertexDescriptor vd, double pp_matrix_sum[3][3], const GeometryTraits &gt)
Definition: curvature.hpp:91
FEVV::Math::Matrix::Square::add
static void add(const ElementType p_matrix[][DIM], ElementType p_matrix_sum[][DIM])
add two matrices
Definition: MatrixOperations.hpp:1038