DGtal  1.4.beta
DGtal::CorrectedNormalCurrentFormula< TRealPoint, TRealVector > Struct Template Reference

Aim: A helper class that provides static methods to compute corrected normal current formulas of curvatures. More...

#include <DGtal/geometry/meshes/CorrectedNormalCurrentFormula.h>

Public Types

typedef TRealPoint RealPoint
 
typedef TRealVector RealVector
 
typedef RealVector::Component Scalar
 
typedef std::vector< RealPointRealPoints
 
typedef std::vector< RealVectorRealVectors
 
typedef SimpleMatrix< Scalar, 3, 3 > RealTensor
 
typedef std::size_t Index
 

Static Public Member Functions

Formulas for mu0 measure
static Scalar mu0ConstantU (const RealPoint &a, const RealPoint &b, const RealPoint &c, const RealVector &u)
 
static Scalar mu0InterpolatedU (const RealPoint &a, const RealPoint &b, const RealPoint &c, const RealVector &ua, const RealVector &ub, const RealVector &uc, bool unit_u=false)
 
static Scalar mu0ConstantU (const RealPoints &pts, const RealVector &u)
 
static Scalar mu0InterpolatedU (const RealPoints &pts, const RealVectors &u, bool unit_u=false)
 
Formulas for mu1 measure
static Scalar mu1ConstantUAtEdge (const RealPoint &a, const RealPoint &b, const RealVector &ur, const RealVector &ul)
 
static Scalar mu1ConstantU (const RealPoint &a, const RealPoint &b, const RealPoint &c, const RealVector &u)
 
static Scalar mu1InterpolatedU (const RealPoint &a, const RealPoint &b, const RealPoint &c, const RealVector &ua, const RealVector &ub, const RealVector &uc, bool unit_u=false)
 
static Scalar mu1ConstantU (const RealPoints &pts, const RealVector &u)
 
static Scalar mu1InterpolatedU (const RealPoints &pts, const RealVectors &u, bool unit_u=false)
 
Formulas for mu2 measure
static Scalar mu2ConstantU (const RealPoint &a, const RealPoint &b, const RealPoint &c, const RealVector &u)
 
static Scalar mu2ConstantUAtVertex (const RealPoint &a, const RealVectors &vu)
 
static Scalar mu2InterpolatedU (const RealPoint &a, const RealPoint &b, const RealPoint &c, const RealVector &ua, const RealVector &ub, const RealVector &uc, bool unit_u=false)
 
static Scalar mu2ConstantU (const RealPoints &pts, const RealVector &u)
 
static Scalar mu2InterpolatedU (const RealPoints &pts, const RealVectors &u, bool unit_u=false)
 
Formulas for muXY measure
static RealTensor muXYConstantU (const RealPoint &a, const RealPoint &b, const RealPoint &c, const RealVector &u)
 
static RealTensor muXYConstantUAtEdge (const RealPoint &a, const RealPoint &b, const RealVector &ur, const RealVector &ul)
 
static RealTensor muXYInterpolatedU (const RealPoint &a, const RealPoint &b, const RealPoint &c, const RealVector &ua, const RealVector &ub, const RealVector &uc, bool unit_u=false)
 
static RealTensor muXYConstantU (const RealPoints &pts, const RealVector &u)
 
static RealTensor muXYInterpolatedU (const RealPoints &pts, const RealVectors &u, bool unit_u=false)
 
Other geometric services
static RealPoint barycenter (const RealPoints &pts)
 
static RealVector averageUnitVector (const RealVectors &vecs)
 

Static Public Attributes

static const Dimension dimension = RealPoint::dimension
 

Detailed Description

template<typename TRealPoint, typename TRealVector>
struct DGtal::CorrectedNormalCurrentFormula< TRealPoint, TRealVector >

Aim: A helper class that provides static methods to compute corrected normal current formulas of curvatures.

Description of class 'CorrectedNormalCurrentFormula'

Template Parameters
TRealPointany model of 3D RealPoint.
TRealVectorany model of 3D RealVector.

Formula for interpolated measures on a triangle with vertices A, B, C, and normal vectors uA, uB, uC:

MU0=-1/6*((uAz + uBz + uCz)*Bx - (uAz + uBz + uCz)*Cx)*Ay + 1/6*((uAz + uBz + uCz)*Ax - (uAz + uBz + uCz)*Cx)*By - 1/6*((uAz + uBz + uCz)*Ax - (uAz + uBz + uCz)*Bx)*Cy + 1/6*((uAy + uBy + uCy)*Bx - (uAy + uBy + uCy)*Cx - (uAx + uBx + uCx)*By + (uAx + uBx + uCx)*Cy)*Az - 1/6*((uAy + uBy + uCy)*Ax - (uAy + uBy + uCy)*Cx - (uAx + uBx + uCx)*Ay + (uAx + uBx + uCx)*Cy)*Bz + 1/6*((uAy + uBy + uCy)*Ax - (uAy + uBy + uCy)*Bx - (uAx + uBx + uCx)*Ay + (uAx + uBx + uCx)*By)*Cz

Let UM=uA+uB+uC.

MU0=-1/6*(uMz*Bx - uMz*Cx)*Ay + 1/6*(uMz*Ax - uMz*Cx)*By - 1/6*(uMz*Ax - uMz*Bx)*Cy + 1/6*(uMy*Bx - uMy*Cx - uMx*By + uMx*Cy)*Az - 1/6*(uMy*Ax - uMy*Cx - uMx*Ay + uMx*Cy)*Bz + 1/6*(uMy*Ax - uMy*Bx - uMx*Ay + uMx*By)*Cz

We see by simple computations that MU0 can be written as (uM = UM/3)

MU0=1/2*det( uM, B-A, C-A )

MU1=1/6*((uBy - uCy)*uAz - (uAy + 2*uCy)*uBz + (uAy + 2*uBy)*uCz)*Ax + 1/6*((uBy + 2*uCy)*uAz - (uAy - uCy)*uBz - (2*uAy + uBy)*uCz)*Bx - 1/6*((2*uBy + uCy)*uAz - (2*uAy + uCy)*uBz - (uAy - uBy)*uCz)*Cx - 1/6*((uBx - uCx)*uAz - (uAx + 2*uCx)*uBz + (uAx + 2*uBx)*uCz)*Ay - 1/6*((uBx + 2*uCx)*uAz - (uAx - uCx)*uBz - (2*uAx + uBx)*uCz)*By + 1/6*((2*uBx + uCx)*uAz - (2*uAx + uCx)*uBz - (uAx - uBx)*uCz)*Cy + 1/6*((uBx - uCx)*uAy - (uAx + 2*uCx)*uBy + (uAx + 2*uBx)*uCy)*Az + 1/6*((uBx + 2*uCx)*uAy - (uAx - uCx)*uBy - (2*uAx + uBx)*uCy)*Bz - 1/6*((2*uBx + uCx)*uAy - (2*uAx + uCx)*uBy - (uAx - uBx)*uCy)*Cz

This formula can also be written in a clearer form with determinants:

6*MU1 = det(u_A+u_B+u_C,u_C-u_B,A) + det(u_A+u_B+u_C,u_A-u_C,B) + det(u_A+u_B+u_C,u_B-u_A,C)

It follows that MU1=1/2( det(uM, u_C-u_B, A) + det(uM, u_A-u_C, B) + det(uM, u_B-u_A, C)

Gaussian curvature measure is

MU2=-1/2*uCx*uBy*uAz + 1/2*uBx*uCy*uAz + 1/2*uCx*uAy*uBz - 1/2*uAx*uCy*uBz - 1/2*uBx*uAy*uCz + 1/2*uAx*uBy*uCz

which is simply MU2=1/2*det( uA, uB, uC )

Anisotropic curvature measure is written as (for X, Y arbitrary vectors, <.|.> the standard Euclidean scalar product).

MUXY = 1/2 det( uM, < Y | uC-uA > X, (B-A)) - det( uM, < Y | uB-uA > X, (C-A))

Definition at line 99 of file CorrectedNormalCurrentFormula.h.

Member Typedef Documentation

◆ Index

template<typename TRealPoint , typename TRealVector >
typedef std::size_t DGtal::CorrectedNormalCurrentFormula< TRealPoint, TRealVector >::Index

Definition at line 107 of file CorrectedNormalCurrentFormula.h.

◆ RealPoint

template<typename TRealPoint , typename TRealVector >
typedef TRealPoint DGtal::CorrectedNormalCurrentFormula< TRealPoint, TRealVector >::RealPoint

Definition at line 101 of file CorrectedNormalCurrentFormula.h.

◆ RealPoints

template<typename TRealPoint , typename TRealVector >
typedef std::vector< RealPoint > DGtal::CorrectedNormalCurrentFormula< TRealPoint, TRealVector >::RealPoints

Definition at line 104 of file CorrectedNormalCurrentFormula.h.

◆ RealTensor

template<typename TRealPoint , typename TRealVector >
typedef SimpleMatrix< Scalar, 3, 3 > DGtal::CorrectedNormalCurrentFormula< TRealPoint, TRealVector >::RealTensor

Definition at line 106 of file CorrectedNormalCurrentFormula.h.

◆ RealVector

template<typename TRealPoint , typename TRealVector >
typedef TRealVector DGtal::CorrectedNormalCurrentFormula< TRealPoint, TRealVector >::RealVector

Definition at line 102 of file CorrectedNormalCurrentFormula.h.

◆ RealVectors

template<typename TRealPoint , typename TRealVector >
typedef std::vector< RealVector > DGtal::CorrectedNormalCurrentFormula< TRealPoint, TRealVector >::RealVectors

Definition at line 105 of file CorrectedNormalCurrentFormula.h.

◆ Scalar

template<typename TRealPoint , typename TRealVector >
typedef RealVector::Component DGtal::CorrectedNormalCurrentFormula< TRealPoint, TRealVector >::Scalar

Definition at line 103 of file CorrectedNormalCurrentFormula.h.

Member Function Documentation

◆ averageUnitVector()

template<typename TRealPoint , typename TRealVector >
static RealVector DGtal::CorrectedNormalCurrentFormula< TRealPoint, TRealVector >::averageUnitVector ( const RealVectors vecs)
inlinestatic

Given a vector of unit vectors, returns their average unit vector.

Parameters
vecsany vector of vectors.
Returns
the average unit vector.

Definition at line 633 of file CorrectedNormalCurrentFormula.h.

634  {
635  RealVector avg;
636  for ( auto v : vecs ) avg += v;
637  auto avg_norm = avg.norm();
638  return avg_norm != 0.0 ? avg / avg_norm : avg;
639  }

Referenced by DGtal::CorrectedNormalCurrentFormula< TRealPoint, TRealVector >::mu0InterpolatedU(), DGtal::CorrectedNormalCurrentFormula< TRealPoint, TRealVector >::mu1InterpolatedU(), DGtal::CorrectedNormalCurrentFormula< TRealPoint, TRealVector >::mu2InterpolatedU(), and DGtal::CorrectedNormalCurrentFormula< TRealPoint, TRealVector >::muXYInterpolatedU().

◆ barycenter()

template<typename TRealPoint , typename TRealVector >
static RealPoint DGtal::CorrectedNormalCurrentFormula< TRealPoint, TRealVector >::barycenter ( const RealPoints pts)
inlinestatic

Given a vector of points, returns its barycenter.

Parameters
ptsany vector of points
Returns
the barycenter of these points.

Definition at line 621 of file CorrectedNormalCurrentFormula.h.

622  {
623  RealPoint b;
624  for ( auto p : pts ) b += p;
625  b /= pts.size();
626  return b;
627  }
PointVector< 3, double > RealPoint

Referenced by DGtal::CorrectedNormalCurrentFormula< TRealPoint, TRealVector >::mu0ConstantU(), DGtal::CorrectedNormalCurrentFormula< TRealPoint, TRealVector >::mu0InterpolatedU(), DGtal::CorrectedNormalCurrentFormula< TRealPoint, TRealVector >::mu1InterpolatedU(), DGtal::CorrectedNormalCurrentFormula< TRealPoint, TRealVector >::mu2InterpolatedU(), and DGtal::CorrectedNormalCurrentFormula< TRealPoint, TRealVector >::muXYInterpolatedU().

◆ mu0ConstantU() [1/2]

template<typename TRealPoint , typename TRealVector >
static Scalar DGtal::CorrectedNormalCurrentFormula< TRealPoint, TRealVector >::mu0ConstantU ( const RealPoint a,
const RealPoint b,
const RealPoint c,
const RealVector u 
)
inlinestatic

Computes mu0 measure (area) of triangle abc given a constant corrected normal vector u.

Parameters
aany point
bany point
cany point
uthe constant corrected normal vector to triangle abc
Returns
the mu0-measure of triangle abc, i.e. its area.

Definition at line 123 of file CorrectedNormalCurrentFormula.h.

126  {
127  return 0.5 * ( b - a ).crossProduct( c - a ).dot( u );
128  }
auto crossProduct(PointVector< 3, LeftEuclideanRing, LeftContainer > const &lhs, PointVector< 3, RightEuclideanRing, RightContainer > const &rhs) -> decltype(DGtal::constructFromArithmeticConversion(lhs, rhs))
Cross product of two 3D Points/Vectors.

References DGtal::crossProduct().

Referenced by DGtal::CorrectedNormalCurrentFormula< TRealPoint, TRealVector >::mu0ConstantU().

◆ mu0ConstantU() [2/2]

template<typename TRealPoint , typename TRealVector >
static Scalar DGtal::CorrectedNormalCurrentFormula< TRealPoint, TRealVector >::mu0ConstantU ( const RealPoints pts,
const RealVector u 
)
inlinestatic

Computes mu0 measure (area) of polygonal face pts given a constant corrected normal vector u.

Parameters
ptsthe (ccw ordered) points forming the vertices of a polygonal face.
uthe constant corrected normal vector to this polygonal face.
Returns
the mu0-measure of the given polygonal face, i.e. its area.

Definition at line 165 of file CorrectedNormalCurrentFormula.h.

166  {
167  if ( pts.size() < 3 ) return 0.0;
168  if ( pts.size() == 3 )
169  return mu0ConstantU( pts[ 0 ], pts[ 1 ], pts[ 2 ], u );
170  const RealPoint b = barycenter( pts );
171  Scalar a = 0.0;
172  for ( Index i = 0; i < pts.size(); i++ )
173  a += mu0ConstantU( b, pts[ i ], pts[ (i+1)%pts.size() ], u );
174  return a;
175  }
SMesh::Index Index
static RealPoint barycenter(const RealPoints &pts)
static Scalar mu0ConstantU(const RealPoint &a, const RealPoint &b, const RealPoint &c, const RealVector &u)

References DGtal::CorrectedNormalCurrentFormula< TRealPoint, TRealVector >::barycenter(), and DGtal::CorrectedNormalCurrentFormula< TRealPoint, TRealVector >::mu0ConstantU().

◆ mu0InterpolatedU() [1/2]

template<typename TRealPoint , typename TRealVector >
static Scalar DGtal::CorrectedNormalCurrentFormula< TRealPoint, TRealVector >::mu0InterpolatedU ( const RealPoint a,
const RealPoint b,
const RealPoint c,
const RealVector ua,
const RealVector ub,
const RealVector uc,
bool  unit_u = false 
)
inlinestatic

Computes mu0 measure (area) of triangle abc given an interpolated corrected normal vector ua, ub, uc.

Parameters
aany point
bany point
cany point
uathe corrected normal vector at point a
ubthe corrected normal vector at point b
ucthe corrected normal vector at point c
unit_uwhen 'true' considers that interpolated corrected normals should be made unitary, otherwise interpolated corrected normals may have smaller norms.
Returns
the mu0-measure of triangle abc, i.e. its area.

Definition at line 143 of file CorrectedNormalCurrentFormula.h.

147  {
148  // MU0=1/2*det( uM, B-A, C-A )
149  // = 1/2 < ( (u_A + u_B + u_C)/3.0 ) | (AB x AC ) >
150  RealVector uM = ( ua+ub+uc ) / 3.0;
151  if ( unit_u )
152  {
153  auto uM_norm = uM.norm();
154  uM = uM_norm == 0.0 ? uM : uM / uM_norm;
155  }
156  return 0.5 * ( b - a ).crossProduct( c - a ).dot( uM );
157  }

References DGtal::crossProduct().

Referenced by DGtal::CorrectedNormalCurrentFormula< TRealPoint, TRealVector >::mu0InterpolatedU().

◆ mu0InterpolatedU() [2/2]

template<typename TRealPoint , typename TRealVector >
static Scalar DGtal::CorrectedNormalCurrentFormula< TRealPoint, TRealVector >::mu0InterpolatedU ( const RealPoints pts,
const RealVectors u,
bool  unit_u = false 
)
inlinestatic

Computes area of polygonal face pts given an interpolated corrected normal vector ua, ub, uc.

Parameters
ptsthe (ccw ordered) points forming the vertices of a polygonal face.
uthe (ccw ordered) normal vectors at the corresponding vertices in pts.
unit_uwhen 'true' considers that interpolated corrected normals should be made unitary, otherwise interpolated corrected normals may have smaller norms.
Returns
the mu0-measure of the given polygonal face, i.e. its area.

Definition at line 186 of file CorrectedNormalCurrentFormula.h.

188  {
189  ASSERT( pts.size() == u.size() );
190  if ( pts.size() < 3 ) return 0.0;
191  if ( pts.size() == 3 )
192  return mu0InterpolatedU( pts[ 0 ], pts[ 1 ], pts[ 2 ],
193  u[ 0 ], u[ 1 ], u[ 2 ], unit_u );
194  const RealPoint b = barycenter( pts );
195  const RealVector ub = averageUnitVector( u );
196  Scalar a = 0.0;
197  for ( Index i = 0; i < pts.size(); i++ )
198  a += mu0InterpolatedU( b, pts[ i ], pts[ (i+1)%pts.size() ],
199  ub, u[ i ], u[ (i+1)%pts.size() ], unit_u );
200  return a;
201  }
static RealVector averageUnitVector(const RealVectors &vecs)
static Scalar mu0InterpolatedU(const RealPoint &a, const RealPoint &b, const RealPoint &c, const RealVector &ua, const RealVector &ub, const RealVector &uc, bool unit_u=false)

References DGtal::CorrectedNormalCurrentFormula< TRealPoint, TRealVector >::averageUnitVector(), DGtal::CorrectedNormalCurrentFormula< TRealPoint, TRealVector >::barycenter(), and DGtal::CorrectedNormalCurrentFormula< TRealPoint, TRealVector >::mu0InterpolatedU().

◆ mu1ConstantU() [1/2]

template<typename TRealPoint , typename TRealVector >
static Scalar DGtal::CorrectedNormalCurrentFormula< TRealPoint, TRealVector >::mu1ConstantU ( const RealPoint a,
const RealPoint b,
const RealPoint c,
const RealVector u 
)
inlinestatic

Computes mu1 measure (twice the mean curvature) of triangle abc given a constant corrected normal vector u.

Parameters
aany point
bany point
cany point
uthe constant corrected normal vector to triangle abc
Returns
the mu1-measure of triangle abc, i.e. twice its mean curvature, always 0.0.

Definition at line 248 of file CorrectedNormalCurrentFormula.h.

251  {
252  (void)a; (void)b; (void)c; (void)u;
253  return 0.0;
254  }

◆ mu1ConstantU() [2/2]

template<typename TRealPoint , typename TRealVector >
static Scalar DGtal::CorrectedNormalCurrentFormula< TRealPoint, TRealVector >::mu1ConstantU ( const RealPoints pts,
const RealVector u 
)
inlinestatic

Computes mu1 measure (twice the mean curvature) of polygonal face pts given a constant corrected normal vector u.

Parameters
ptsthe (ccw ordered) points forming the vertices of a polygonal face.
uthe constant corrected normal vector to this polygonal face.
Returns
the mu1-measure of the given polygonal face, i.e. twice its mean curvature, always 0.0.

Definition at line 291 of file CorrectedNormalCurrentFormula.h.

292  {
293  (void)pts; //not used parameter
294  (void)u;
295  return 0.0;
296  }

◆ mu1ConstantUAtEdge()

template<typename TRealPoint , typename TRealVector >
static Scalar DGtal::CorrectedNormalCurrentFormula< TRealPoint, TRealVector >::mu1ConstantUAtEdge ( const RealPoint a,
const RealPoint b,
const RealVector ur,
const RealVector ul 
)
inlinestatic

Computes mu1 measure (twice the mean curvature) at edge ab given a constant corrected normal vector ur to the right of ab, and a constant corrected normal vector ul to the left of ab.

The formula is \( int_{ab} \Psi \langle e | e_1 \rangle dH^1 \), where e is the unit vector parallel to ab, and \( e_1 \) is the unit vector orthogonal to ur and ul, and \( \Psi \) is the angle between these two vectors.

Parameters
aany point
bany point
urthe constant corrected normal vector to the right of oriented edge ab
ulthe constant corrected normal vector to the left of oriented edge ab
Returns
the mu1-measure of edge ab, i.e. twice its mean curvature.

Definition at line 226 of file CorrectedNormalCurrentFormula.h.

229  {
230  RealVector e = b - a;
231  RealVector e1 = ur.crossProduct( ul );
232  const Scalar l1 = std::min( e1.norm(), 1.0 );
233  if ( l1 < 1e-10 ) return 0.0;
234  e1 /= l1;
235  const Scalar Psi = asin( l1 );
236  return -Psi * e.dot( e1 );
237  }

◆ mu1InterpolatedU() [1/2]

template<typename TRealPoint , typename TRealVector >
static Scalar DGtal::CorrectedNormalCurrentFormula< TRealPoint, TRealVector >::mu1InterpolatedU ( const RealPoint a,
const RealPoint b,
const RealPoint c,
const RealVector ua,
const RealVector ub,
const RealVector uc,
bool  unit_u = false 
)
inlinestatic

Computes mu1 measure (twice the mean curvature) of triangle abc given an interpolated corrected normal vector ua, ub, uc.

Parameters
aany point
bany point
cany point
uathe corrected normal vector at point a
ubthe corrected normal vector at point b
ucthe corrected normal vector at point c
unit_uwhen 'true' considers that interpolated corrected normals should be made unitary, otherwise interpolated corrected normals may have smaller norms.
Returns
the mu1-measure of triangle abc, i.e. twice its mean curvature.

Definition at line 271 of file CorrectedNormalCurrentFormula.h.

275  {
276  // MU1=1/2( | uM u_C-u_B A | + | uM u_A-u_C B | + | uM u_B-u_A C |
277  RealVector uM = ( ua+ub+uc ) / 3.0;
278  if ( unit_u ) uM /= uM.norm();
279  return 0.5 * ( uM.crossProduct( uc - ub ).dot( a )
280  + uM.crossProduct( ua - uc ).dot( b )
281  + uM.crossProduct( ub - ua ).dot( c ) );
282  }

Referenced by DGtal::CorrectedNormalCurrentFormula< TRealPoint, TRealVector >::mu1InterpolatedU().

◆ mu1InterpolatedU() [2/2]

template<typename TRealPoint , typename TRealVector >
static Scalar DGtal::CorrectedNormalCurrentFormula< TRealPoint, TRealVector >::mu1InterpolatedU ( const RealPoints pts,
const RealVectors u,
bool  unit_u = false 
)
inlinestatic

Computes mean curvature of polygonal face pts given an interpolated corrected normal vector ua, ub, uc.

Parameters
ptsthe (ccw ordered) points forming the vertices of a polygonal face.
uthe (ccw ordered) normal vectors at the corresponding vertices in pts.
unit_uwhen 'true' considers that interpolated corrected normals should be made unitary, otherwise interpolated corrected normals may have smaller norms.
Returns
the mu1-measure of the given polygonal face, i.e. its mean curvature.

Definition at line 307 of file CorrectedNormalCurrentFormula.h.

309  {
310  ASSERT( pts.size() == u.size() );
311  if ( pts.size() < 3 ) return 0.0;
312  if ( pts.size() == 3 )
313  return mu1InterpolatedU( pts[ 0 ], pts[ 1 ], pts[ 2 ],
314  u[ 0 ], u[ 1 ], u[ 2 ], unit_u );
315  const RealPoint b = barycenter( pts );
316  const RealVector ub = averageUnitVector( u );
317  Scalar a = 0.0;
318  for ( Index i = 0; i < pts.size(); i++ )
319  a += mu1InterpolatedU( b, pts[ i ], pts[ (i+1)%pts.size() ],
320  ub, u[ i ], u[ (i+1)%pts.size() ], unit_u );
321  return a;
322  }
static Scalar mu1InterpolatedU(const RealPoint &a, const RealPoint &b, const RealPoint &c, const RealVector &ua, const RealVector &ub, const RealVector &uc, bool unit_u=false)

References DGtal::CorrectedNormalCurrentFormula< TRealPoint, TRealVector >::averageUnitVector(), DGtal::CorrectedNormalCurrentFormula< TRealPoint, TRealVector >::barycenter(), and DGtal::CorrectedNormalCurrentFormula< TRealPoint, TRealVector >::mu1InterpolatedU().

◆ mu2ConstantU() [1/2]

template<typename TRealPoint , typename TRealVector >
static Scalar DGtal::CorrectedNormalCurrentFormula< TRealPoint, TRealVector >::mu2ConstantU ( const RealPoint a,
const RealPoint b,
const RealPoint c,
const RealVector u 
)
inlinestatic

Computes mu2 measure (Gaussian curvature) of triangle abc given a constant corrected normal vector u.

Parameters
aany point
bany point
cany point
uthe constant corrected normal vector to triangle abc
Returns
the mu2-measure of triangle abc, i.e. its Gaussian curvature, always 0.0.

Definition at line 339 of file CorrectedNormalCurrentFormula.h.

342  {
343  (void)a; (void)b; (void)c; (void)u;
344  return 0.0;
345  }

◆ mu2ConstantU() [2/2]

template<typename TRealPoint , typename TRealVector >
static Scalar DGtal::CorrectedNormalCurrentFormula< TRealPoint, TRealVector >::mu2ConstantU ( const RealPoints pts,
const RealVector u 
)
inlinestatic

Computes mu2 measure (Gaussian curvature) of polygonal face pts given a constant corrected normal vector u.

Parameters
ptsthe (ccw ordered) points forming the vertices of a polygonal face.
uthe constant corrected normal vector to this polygonal face.
Returns
the mu2-measure of the given polygonal face, i.e. its Gaussian curvature, always 0.0.

Definition at line 431 of file CorrectedNormalCurrentFormula.h.

432  {
433  (void) pts; (void) u;
434  return 0.0;
435  }

◆ mu2ConstantUAtVertex()

template<typename TRealPoint , typename TRealVector >
static Scalar DGtal::CorrectedNormalCurrentFormula< TRealPoint, TRealVector >::mu2ConstantUAtVertex ( const RealPoint a,
const RealVectors vu 
)
inlinestatic

Computes mu2 measure (Gaussian curvature) at given vertex a, surrounded by faces with constant corrected normal vectors vu.

Parameters
aany point
vuthe vector of constant corrected normal vectors of the faces incident to a.
Returns
the mu2-measure at this vertex, i.e. its Gaussian curvature.
Precondition
the unit normal vectors should lie in the same hemisphere of the Gauss sphere.

Definition at line 360 of file CorrectedNormalCurrentFormula.h.

362  {
363  typedef SpaceND< dimension > Space;
364  typedef SphericalTriangle<Space> ST;
365  RealVector avg_u;
366  for ( const auto& u : vu ) avg_u += u;
367  const Scalar l = avg_u.norm();
368  if ( l < 1e-10 )
369  {
370  trace.warning() << "[CorrectedNormalCurrentFormula::mu2ConstantUAtVertex]"
371  << " Invalid surrounding normal vectors at vertex "
372  << a << "."
373  << " Unit normal vectors should lie strictly in the same"
374  << " hemisphere." << std::endl;
375  return 0.0;
376  }
377  avg_u /= l;
378  Scalar mu2 = 0.0;
379  const auto n = vu.size();
380  for ( size_t i = 0; i < n; ++i )
381  {
382  ST st( avg_u, vu[ i ], vu[ (i+1) % n ] );
383  mu2 += st.algebraicArea();
384  }
385  return mu2;
386  }
std::ostream & warning()
Trace trace
Definition: Common.h:153

References DGtal::trace, and DGtal::Trace::warning().

◆ mu2InterpolatedU() [1/2]

template<typename TRealPoint , typename TRealVector >
static Scalar DGtal::CorrectedNormalCurrentFormula< TRealPoint, TRealVector >::mu2InterpolatedU ( const RealPoint a,
const RealPoint b,
const RealPoint c,
const RealVector ua,
const RealVector ub,
const RealVector uc,
bool  unit_u = false 
)
inlinestatic

Computes mu2 measure (Gaussian curvature) of triangle abc given an interpolated corrected normal vector ua, ub, uc.

Parameters
aany point
bany point
cany point
uathe corrected normal vector at point a
ubthe corrected normal vector at point b
ucthe corrected normal vector at point c
unit_uwhen 'true' considers that interpolated corrected normals should be made unitary, otherwise interpolated corrected normals may have smaller norms.
Returns
the mu2-measure of triangle abc, i.e. its Gaussian curvature.

Definition at line 403 of file CorrectedNormalCurrentFormula.h.

407  {
408  (void)a; //not used
409  (void)b;
410  (void)c;
411 
412  // Using non unitary interpolated normals give
413  // MU2=1/2*det( uA, uB, uC )
414  // When normals are unitary, it is the area of a spherical triangle.
415  if ( unit_u )
416  {
417  typedef SpaceND< dimension > Space;
418  SphericalTriangle<Space> ST( ua, ub, uc ); // check order.
419  return ST.algebraicArea();
420  }
421  else
422  return 0.5 * ( ua.crossProduct( ub ).dot( uc ) );
423  }

References DGtal::SphericalTriangle< TSpace >::algebraicArea().

Referenced by DGtal::CorrectedNormalCurrentFormula< TRealPoint, TRealVector >::mu2InterpolatedU().

◆ mu2InterpolatedU() [2/2]

template<typename TRealPoint , typename TRealVector >
static Scalar DGtal::CorrectedNormalCurrentFormula< TRealPoint, TRealVector >::mu2InterpolatedU ( const RealPoints pts,
const RealVectors u,
bool  unit_u = false 
)
inlinestatic

Computes Gaussian curvature of polygonal face pts given an interpolated corrected normal vector ua, ub, uc.

Parameters
ptsthe (ccw ordered) points forming the vertices of a polygonal face.
uthe (ccw ordered) normal vectors at the corresponding vertices in pts.
unit_uwhen 'true' considers that interpolated corrected normals should be made unitary, otherwise interpolated corrected normals may have smaller norms.
Returns
the mu2-measure of the given polygonal face, i.e. its Gaussian curvature.

Definition at line 446 of file CorrectedNormalCurrentFormula.h.

448  {
449  ASSERT( pts.size() == u.size() );
450  if ( pts.size() < 3 ) return 0.0;
451  if ( pts.size() == 3 )
452  return mu2InterpolatedU( pts[ 0 ], pts[ 1 ], pts[ 2 ],
453  u[ 0 ], u[ 1 ], u[ 2 ], unit_u );
454  const RealPoint b = barycenter( pts );
455  const RealVector ub = averageUnitVector( u );
456  Scalar a = 0.0;
457  for ( Index i = 0; i < pts.size(); i++ )
458  a += mu2InterpolatedU( b, pts[ i ], pts[ (i+1)%pts.size() ],
459  ub, u[ i ], u[ (i+1)%pts.size() ], unit_u );
460  return a;
461  }
static Scalar mu2InterpolatedU(const RealPoint &a, const RealPoint &b, const RealPoint &c, const RealVector &ua, const RealVector &ub, const RealVector &uc, bool unit_u=false)

References DGtal::CorrectedNormalCurrentFormula< TRealPoint, TRealVector >::averageUnitVector(), DGtal::CorrectedNormalCurrentFormula< TRealPoint, TRealVector >::barycenter(), and DGtal::CorrectedNormalCurrentFormula< TRealPoint, TRealVector >::mu2InterpolatedU().

◆ muXYConstantU() [1/2]

template<typename TRealPoint , typename TRealVector >
static RealTensor DGtal::CorrectedNormalCurrentFormula< TRealPoint, TRealVector >::muXYConstantU ( const RealPoint a,
const RealPoint b,
const RealPoint c,
const RealVector u 
)
inlinestatic

Computes muXY measure (anisotropic curvature) of triangle abc given a constant corrected normal vector u.

Parameters
aany point
bany point
cany point
uthe constant corrected normal vector to triangle abc
Returns
the muXY-measure of triangle abc, i.e. its anisotropic curvature, always 0.0.

Definition at line 478 of file CorrectedNormalCurrentFormula.h.

481  {
482  (void)a; (void)b; (void)c; (void)u;
483  return RealTensor { 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0 };
484  }

◆ muXYConstantU() [2/2]

template<typename TRealPoint , typename TRealVector >
static RealTensor DGtal::CorrectedNormalCurrentFormula< TRealPoint, TRealVector >::muXYConstantU ( const RealPoints pts,
const RealVector u 
)
inlinestatic

Computes muXY measure (anisotropic curvature) of polygonal face pts given a constant corrected normal vector u.

Parameters
ptsthe (ccw ordered) points forming the vertices of a polygonal face.
uthe constant corrected normal vector to this polygonal face.
Returns
the muXY-measure of the given polygonal face, i.e. its anisotropic curvature, always 0.0.

Definition at line 577 of file CorrectedNormalCurrentFormula.h.

578  {
579  (void)pts; (void)u;
580  return RealTensor { 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0 };
581  }

◆ muXYConstantUAtEdge()

template<typename TRealPoint , typename TRealVector >
static RealTensor DGtal::CorrectedNormalCurrentFormula< TRealPoint, TRealVector >::muXYConstantUAtEdge ( const RealPoint a,
const RealPoint b,
const RealVector ur,
const RealVector ul 
)
inlinestatic

Computes muXY measure (anisotropic measure) at edge ab given a constant corrected normal vector ur to the right of ab, and a constant corrected normal vector ul to the left of ab.

The formula is \( int_{ab} \Psi \langle e | e_1 \rangle dH^1 \), where e is the unit vector parallel to ab, and \( e_1 \) is the unit vector orthogonal to ur and ul, and \( \Psi \) is the angle between these two vectors.

Parameters
aany point
bany point
urthe constant corrected normal vector to the right of oriented edge ab
ulthe constant corrected normal vector to the left of oriented edge ab
Returns
the mu1-measure of edge ab, i.e. twice its mean curvature.

Definition at line 502 of file CorrectedNormalCurrentFormula.h.

505  {
506  RealTensor M { 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0 };
507  RealVector e = b - a;
508  RealVector e1 = ur.crossProduct( ul );
509  const Scalar l1 = std::min( e1.norm(), 1.0 );
510  if ( l1 < 1e-10 ) return M;
511  e1 /= l1;
512  const Scalar psi = asin( l1 );
513  const Scalar sin_psi = sin( psi );
514  const Scalar sin_2psi = sin( 2.0 * psi );
515  const RealVector e_x_ur = e. crossProduct( ur );
516  const RealVector e1_x_ur = e1.crossProduct( ur );
517  const RealVector e_x_e1_x_ur = e. crossProduct( e1_x_ur );
518  for ( Dimension i = 0; i < 3; i++ )
519  for ( Dimension j = 0; j < 3; j++ )
520  {
521  Scalar v = (-psi - 0.5*sin_2psi) * e_x_ur[ i ] * e1_x_ur[ j ]
522  + (sin_psi*sin_psi) * ( e_x_ur[ i ] * ur[ j ]
523  - e_x_e1_x_ur[ i ] * e1_x_ur[ j ] )
524  + (psi - 0.5*sin_2psi) * e_x_e1_x_ur[ i ] * ur[ j ];
525  M.setComponent( i, j, -0.5 * v );
526  }
527  return M;
528  }
DGtal::uint32_t Dimension
Definition: Common.h:136

References DGtal::crossProduct().

◆ muXYInterpolatedU() [1/2]

template<typename TRealPoint , typename TRealVector >
static RealTensor DGtal::CorrectedNormalCurrentFormula< TRealPoint, TRealVector >::muXYInterpolatedU ( const RealPoint a,
const RealPoint b,
const RealPoint c,
const RealVector ua,
const RealVector ub,
const RealVector uc,
bool  unit_u = false 
)
inlinestatic

Computes muXY measure (anisotropic curvature) of triangle abc given an interpolated corrected normal vector ua, ub, uc.

Parameters
aany point
bany point
cany point
uathe corrected normal vector at point a
ubthe corrected normal vector at point b
ucthe corrected normal vector at point c
unit_uwhen 'true' considers that interpolated corrected normals should be made unitary, otherwise interpolated corrected normals may have smaller norms.
Returns
the muXY-measure of triangle abc, i.e. its anisotropic curvature.

Definition at line 543 of file CorrectedNormalCurrentFormula.h.

547  {
548  RealTensor T;
549  // MUXY = 1/2 < uM | < Y | uc-ua > X x (b-a) - < Y | ub-ua > X x (c-a) >
550  // MUXY = 1/2 ( < Y | ub-ua > | X uM (c-a) | - < Y | uc-ua > | X uM (b-a) | )
551  RealVector uM = ( ua+ub+uc ) / 3.0;
552  if ( unit_u ) uM /= uM.norm();
553  const RealVector uac = uc - ua;
554  const RealVector uab = ub - ua;
555  const RealVector ab = b - a;
556  const RealVector ac = c - a;
557  for ( Dimension i = 0; i < 3; ++i ) {
558  RealVector X = RealVector::base( i, 1.0 );
559  for ( Dimension j = 0; j < 3; ++j ) {
560  // Since RealVector Y = RealVector::base( j, 1.0 );
561  // < Y | uac > = uac[ j ]
562  const Scalar tij =
563  0.5 * uM.dot( uac[ j ] * X.crossProduct( ab )
564  - uab[ j ] * X.crossProduct( ac ) );
565  T.setComponent( i, j, tij );
566  }
567  }
568  return T;
569  }

References DGtal::SimpleMatrix< TComponent, TM, TN >::setComponent().

Referenced by DGtal::CorrectedNormalCurrentFormula< TRealPoint, TRealVector >::muXYInterpolatedU().

◆ muXYInterpolatedU() [2/2]

template<typename TRealPoint , typename TRealVector >
static RealTensor DGtal::CorrectedNormalCurrentFormula< TRealPoint, TRealVector >::muXYInterpolatedU ( const RealPoints pts,
const RealVectors u,
bool  unit_u = false 
)
inlinestatic

Computes anisotropic curvature of polygonal face pts given an interpolated corrected normal vector ua, ub, uc.

Parameters
ptsthe (ccw ordered) points forming the vertices of a polygonal face.
uthe (ccw ordered) normal vectors at the corresponding vertices in pts.
unit_uwhen 'true' considers that interpolated corrected normals should be made unitary, otherwise interpolated corrected normals may have smaller norms.
Returns
the muXY-measure of the given polygonal face, i.e. its anisotropic curvature.

Definition at line 592 of file CorrectedNormalCurrentFormula.h.

594  {
595  ASSERT( pts.size() == u.size() );
596  if ( pts.size() < 3 ) return RealTensor { 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0 };
597  if ( pts.size() == 3 )
598  return muXYInterpolatedU( pts[ 0 ], pts[ 1 ], pts[ 2 ],
599  u[ 0 ], u[ 1 ], u[ 2 ], unit_u );
600  const RealPoint b = barycenter( pts );
601  const RealVector ub = averageUnitVector( u );
602  RealTensor T = { 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0 };
603  for ( Index i = 0; i < pts.size(); i++ )
604  T += muXYInterpolatedU( b, pts[ i ], pts[ (i+1)%pts.size() ],
605  ub, u[ i ], u[ (i+1)%pts.size() ], unit_u );
606  return T;
607  }
static RealTensor muXYInterpolatedU(const RealPoint &a, const RealPoint &b, const RealPoint &c, const RealVector &ua, const RealVector &ub, const RealVector &uc, bool unit_u=false)

References DGtal::CorrectedNormalCurrentFormula< TRealPoint, TRealVector >::averageUnitVector(), DGtal::CorrectedNormalCurrentFormula< TRealPoint, TRealVector >::barycenter(), and DGtal::CorrectedNormalCurrentFormula< TRealPoint, TRealVector >::muXYInterpolatedU().

Field Documentation

◆ dimension

template<typename TRealPoint , typename TRealVector >
const Dimension DGtal::CorrectedNormalCurrentFormula< TRealPoint, TRealVector >::dimension = RealPoint::dimension
static

Definition at line 108 of file CorrectedNormalCurrentFormula.h.


The documentation for this struct was generated from the following file: