DGtal  1.5.beta
DGtal::ATSolver2D< TKSpace, TLinearAlgebra > Class Template Reference

Aim: This class solves Ambrosio-Tortorelli functional on a two-dimensional digital space (a 2D grid or 2D digital surface) for a piecewise smooth scalar/vector function u represented as one/several 2-form(s) and a discontinuity function v represented as a 0-form. The 2-form(s) u is a regularized approximation of an input vector data g, while v represents the set of discontinuities of u. The norm chosen for u is the $$l_2$$-norm. More...

#include <DGtal/dec/ATSolver2D.h>

## Public Types

enum  CellOutputPolicy { Average , Minimum , Maximum }

typedef TKSpace KSpace

typedef TLinearAlgebra LinearAlgebra

typedef ATSolver2D< KSpace, LinearAlgebraSelf

typedef KSpace::Space Space

typedef Space::RealVector RealVector

typedef RealVector::Component Scalar

typedef KSpace::SCell SCell

typedef KSpace::Cell Cell

typedef KSpace::Surfel Surfel

typedef HyperRectDomain< SpaceDomain

typedef DiscreteExteriorCalculus< 2, dimension, LinearAlgebraCalculus

typedef KSpace::template SurfelMap< double >::Type SmallestEpsilonMap

typedef Calculus::Index Index

typedef Calculus::PrimalForm0 PrimalForm0

typedef Calculus::PrimalForm1 PrimalForm1

typedef Calculus::PrimalForm2 PrimalForm2

typedef Calculus::PrimalIdentity0 PrimalIdentity0

typedef Calculus::PrimalIdentity1 PrimalIdentity1

typedef Calculus::PrimalIdentity2 PrimalIdentity2

typedef Calculus::PrimalDerivative0 PrimalDerivative0

typedef Calculus::PrimalDerivative1 PrimalDerivative1

typedef Calculus::PrimalAntiderivative1 PrimalAntiderivative1

typedef Calculus::PrimalAntiderivative2 PrimalAntiderivative2

typedef Calculus::PrimalHodge0 PrimalHodge0

typedef Calculus::PrimalHodge1 PrimalHodge1

typedef Calculus::PrimalHodge2 PrimalHodge2

typedef KSpace::template SurfelMap< Index >::Type Surfel2IndexMap

typedef EigenLinearAlgebraBackend::SolverSimplicialLDLT LinearAlgebraSolver

typedef DiscreteExteriorCalculusSolver< Calculus, LinearAlgebraSolver, 2, PRIMAL, 2, PRIMALSolverU2

typedef DiscreteExteriorCalculusSolver< Calculus, LinearAlgebraSolver, 0, PRIMAL, 0, PRIMALSolverV0

## Public Member Functions

Standard services
ATSolver2D (ConstAlias< Calculus > aCalculus, int aVerbose=0)

ATSolver2D ()=delete

~ATSolver2D ()=default

ATSolver2D (const ATSolver2D &other)=default

ATSolver2D (ATSolver2D &&other)=default

ATSolver2Doperator= (const ATSolver2D &other)=default

ATSolver2Doperator= (ATSolver2D &&other)=default

Index size (const int order) const

Initialization services
template<typename VectorFieldInput , typename SurfelRangeConstIterator >
void initInputVectorFieldU2 (const VectorFieldInput &input, SurfelRangeConstIterator itB, SurfelRangeConstIterator itE, bool normalize=false)

template<typename ScalarFieldInput , typename SurfelRangeConstIterator >
void initInputScalarFieldU2 (const ScalarFieldInput &input, SurfelRangeConstIterator itB, SurfelRangeConstIterator itE)

template<typename ScalarVector >
Index initInputVectorFieldU2 (const std::map< Surfel, ScalarVector > &input, bool normalize=false)

template<typename Scalar >
Index initInputScalarFieldU2 (const std::map< Surfel, Scalar > &input)

void setUp (double a, double l)

void setUp (double a, double l, const std::map< Surfel, Scalar > &weights)

template<typename AlphaWeights , typename SurfelRangeConstIterator >
void setUp (double a, double l, const AlphaWeights &weights, SurfelRangeConstIterator itB, SurfelRangeConstIterator itE)

void setEpsilon (double e)

Optimization services
bool solveOneAlternateStep ()

bool solveForEpsilon (double eps, double n_oo_max=1e-4, unsigned int iter_max=10)

bool solveGammaConvergence (double eps1=2.0, double eps2=0.25, double epsr=2.0, bool compute_smallest_epsilon_map=false, double n_oo_max=1e-4, unsigned int iter_max=10)

void normalizeU2 ()

std::tuple< double, double, double > diffV0 () const

Access services
std::tuple< double, double, double > checkV0 () const

PrimalForm0 getV0 () const

PrimalForm1 getV1 () const

PrimalForm2 getV2 () const

PrimalForm2 getU2 (Dimension k) const

template<typename VectorFieldOutput , typename SurfelRangeConstIterator >
void getOutputVectorFieldU2 (VectorFieldOutput &output, SurfelRangeConstIterator itB, SurfelRangeConstIterator itE)

template<typename ScalarFieldOutput , typename SurfelRangeConstIterator >
void getOutputScalarFieldU2 (ScalarFieldOutput &output, SurfelRangeConstIterator itB, SurfelRangeConstIterator itE)

template<typename ScalarFieldOutput , typename CellRangeConstIterator >
void getOutputScalarFieldV0 (ScalarFieldOutput &output, CellRangeConstIterator itB, CellRangeConstIterator itE, CellOutputPolicy policy=CellOutputPolicy::Average)

void updateSmallestEpsilonMap (const double threshold=.5)

Interface services
void selfDisplay (std::ostream &out) const

bool isValid () const

## Data Fields

Surfel2IndexMap surfel2idx

SmallestEpsilonMap smallest_epsilon_map

double alpha

double lambda

double epsilon

bool normalize_u2
Indicates whether to normalize U (unit norm) at each iteration or not. More...

int verbose
Tells the verbose level. More...

## Static Public Attributes

static const Dimension dimension = KSpace::dimension

## Protected Member Functions

Hidden services
void initOperators ()
Initializes the operators. More...

## Protected Attributes

CountedConstPtrOrConstPtr< CalculusptrCalculus
A smart (or not) pointer to a calculus object. More...

PrimalDerivative0 primal_D0
the derivative operator for primal 0-forms More...

PrimalDerivative1 primal_D1
the derivative operator for primal 1-forms More...

PrimalDerivative0 M01
The primal vertex to edge average operator. More...

PrimalDerivative1 M12
The primal edge to face average operator. More...

The antiderivative of primal 2-forms. More...

PrimalIdentity2 alpha_Id2
The alpha-weighted identity operator for primal 2-forms (stored for performance) More...

PrimalIdentity0 l_1_over_4e_Id0
The 1/(4epsilon)-weighted identity operator for primal 0-forms (stored for performance) More...

std::vector< PrimalForm2g2
The N-array of input primal 2-forms g. More...

std::vector< PrimalForm2alpha_g2
The alpha-weighted N-array of input primal 2-forms g. More...

std::vector< PrimalForm2u2
The N-array of regularized primal 2-forms u. More...

PrimalForm0 v0

PrimalForm0 former_v0
The primal 0-form v at the previous iteration. More...

PrimalForm0 l_1_over_4e
The primal 0-form lambda/(4epsilon) (stored for performance) More...

## Detailed Description

### template<typename TKSpace, typename TLinearAlgebra = EigenLinearAlgebraBackend> class DGtal::ATSolver2D< TKSpace, TLinearAlgebra >

Aim: This class solves Ambrosio-Tortorelli functional on a two-dimensional digital space (a 2D grid or 2D digital surface) for a piecewise smooth scalar/vector function u represented as one/several 2-form(s) and a discontinuity function v represented as a 0-form. The 2-form(s) u is a regularized approximation of an input vector data g, while v represents the set of discontinuities of u. The norm chosen for u is the $$l_2$$-norm.

Description of template class 'ATSolver2D'

Template Parameters
 TKSpace any model of CCellularGridSpaceND, e.g KhalimskySpaceND TLinearAlgebra any back-end for performing linear algebra, default is EigenLinearAlgebraBackend.
// Typical use (with appropriate definitions for types and variables).
typedef DiscreteExteriorCalculusFactory<EigenLinearAlgebraBackend> CalculusFactory;
const auto calculus = CalculusFactory::createFromNSCells<2>( surfels.begin(), surfels.end() );
ATSolver2D< KSpace > at_solver(calculus, 1);
at_solver.initInputVectorFieldU2( normals, surfels.cbegin(), surfels.cend() );
at_solver.setUp( alpha_at, lambda_at );
at_solver.solveGammaConvergence( 2.0, 0.5, 2.0 );
at_solver.getOutputVectorFieldU2( normals, surfels.cbegin(), surfels.cend() );
PolyCalculus * calculus
exampleSurfaceATNormals.cpp

Definition at line 90 of file ATSolver2D.h.

## ◆ Calculus

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
 typedef DiscreteExteriorCalculus<2,dimension, LinearAlgebra> DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::Calculus

Definition at line 114 of file ATSolver2D.h.

## ◆ Cell

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
 typedef KSpace::Cell DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::Cell

Definition at line 111 of file ATSolver2D.h.

## ◆ Domain

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
 typedef HyperRectDomain DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::Domain

Definition at line 113 of file ATSolver2D.h.

## ◆ Index

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
 typedef Calculus::Index DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::Index

Definition at line 116 of file ATSolver2D.h.

## ◆ KSpace

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
 typedef TKSpace DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::KSpace

Definition at line 95 of file ATSolver2D.h.

## ◆ LinearAlgebra

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
 typedef TLinearAlgebra DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::LinearAlgebra

Definition at line 96 of file ATSolver2D.h.

## ◆ LinearAlgebraSolver

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
 typedef EigenLinearAlgebraBackend::SolverSimplicialLDLT DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::LinearAlgebraSolver

Definition at line 138 of file ATSolver2D.h.

## ◆ PrimalAntiderivative1

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
 typedef Calculus::PrimalAntiderivative1 DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::PrimalAntiderivative1

Definition at line 125 of file ATSolver2D.h.

## ◆ PrimalAntiderivative2

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
 typedef Calculus::PrimalAntiderivative2 DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::PrimalAntiderivative2

Definition at line 126 of file ATSolver2D.h.

## ◆ PrimalDerivative0

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
 typedef Calculus::PrimalDerivative0 DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::PrimalDerivative0

Definition at line 123 of file ATSolver2D.h.

## ◆ PrimalDerivative1

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
 typedef Calculus::PrimalDerivative1 DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::PrimalDerivative1

Definition at line 124 of file ATSolver2D.h.

## ◆ PrimalForm0

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
 typedef Calculus::PrimalForm0 DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::PrimalForm0

Definition at line 117 of file ATSolver2D.h.

## ◆ PrimalForm1

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
 typedef Calculus::PrimalForm1 DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::PrimalForm1

Definition at line 118 of file ATSolver2D.h.

## ◆ PrimalForm2

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
 typedef Calculus::PrimalForm2 DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::PrimalForm2

Definition at line 119 of file ATSolver2D.h.

## ◆ PrimalHodge0

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
 typedef Calculus::PrimalHodge0 DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::PrimalHodge0

Definition at line 127 of file ATSolver2D.h.

## ◆ PrimalHodge1

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
 typedef Calculus::PrimalHodge1 DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::PrimalHodge1

Definition at line 128 of file ATSolver2D.h.

## ◆ PrimalHodge2

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
 typedef Calculus::PrimalHodge2 DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::PrimalHodge2

Definition at line 129 of file ATSolver2D.h.

## ◆ PrimalIdentity0

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
 typedef Calculus::PrimalIdentity0 DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::PrimalIdentity0

Definition at line 120 of file ATSolver2D.h.

## ◆ PrimalIdentity1

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
 typedef Calculus::PrimalIdentity1 DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::PrimalIdentity1

Definition at line 121 of file ATSolver2D.h.

## ◆ PrimalIdentity2

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
 typedef Calculus::PrimalIdentity2 DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::PrimalIdentity2

Definition at line 122 of file ATSolver2D.h.

## ◆ RealVector

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
 typedef Space::RealVector DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::RealVector

Definition at line 108 of file ATSolver2D.h.

## ◆ Scalar

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
 typedef RealVector::Component DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::Scalar

Definition at line 109 of file ATSolver2D.h.

## ◆ SCell

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
 typedef KSpace::SCell DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::SCell

Definition at line 110 of file ATSolver2D.h.

## ◆ Self

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
 typedef ATSolver2D< KSpace, LinearAlgebra > DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::Self

Definition at line 97 of file ATSolver2D.h.

## ◆ SmallestEpsilonMap

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
 typedef KSpace::template SurfelMap::Type DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::SmallestEpsilonMap

Definition at line 115 of file ATSolver2D.h.

## ◆ SolverU2

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
 typedef DiscreteExteriorCalculusSolver DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::SolverU2

Definition at line 139 of file ATSolver2D.h.

## ◆ SolverV0

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
 typedef DiscreteExteriorCalculusSolver DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::SolverV0

Definition at line 140 of file ATSolver2D.h.

## ◆ Space

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
 typedef KSpace::Space DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::Space

Definition at line 107 of file ATSolver2D.h.

## ◆ Surfel

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
 typedef KSpace::Surfel DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::Surfel

Definition at line 112 of file ATSolver2D.h.

## ◆ Surfel2IndexMap

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
 typedef KSpace::template SurfelMap::Type DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::Surfel2IndexMap

Definition at line 130 of file ATSolver2D.h.

## ◆ CellOutputPolicy

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>

Specifies how to merge the different values of 0-form v at cell vertices when outputing the 0-form v for a range of cells (either pointels, linels, surfels).

Enumerator
Average

compute average values at cell vertices

Minimum

compute minimum value at cell vertices,

Maximum

compute maximum value at cell vertices

Definition at line 103 of file ATSolver2D.h.

103  { Average,
104  Minimum,
105  Maximum,
106  };
@ Maximum
compute maximum value at cell vertices
Definition: ATSolver2D.h:105
@ Average
compute average values at cell vertices
Definition: ATSolver2D.h:103
@ Minimum
compute minimum value at cell vertices,
Definition: ATSolver2D.h:104

## ◆ ATSolver2D() [1/4]

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
 DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::ATSolver2D ( ConstAlias< Calculus > aCalculus, int aVerbose = 0 )
inline

Prepare an AT-solver from a valid calculus.

Parameters
 aCalculus any valid calculus aVerbose tells how the solver displays computing information: 0 none, 1 more, 2 even more...
DiscreteExteriorCalculusFactory for creating calculus objects.

Definition at line 202 of file ATSolver2D.h.

203  : ptrCalculus( aCalculus ),
207  g2(), alpha_g2(), u2(), v0( *ptrCalculus ), former_v0( *ptrCalculus ),
208  l_1_over_4e( *ptrCalculus ), verbose( aVerbose )
209  {
210  if ( verbose >= 2 )
211  trace.info() << "[ATSolver::ATSolver] " << *ptrCalculus << std::endl;
212  initOperators();
213  const auto size2 = ptrCalculus->kFormLength( 2, PRIMAL );
214  for ( Index index = 0; index < size2; ++index) {
215  const auto& calculus_cell = ptrCalculus->getSCell( 2, PRIMAL, index );
216  surfel2idx[ calculus_cell ] = index;
217  }
218  }
std::vector< PrimalForm2 > alpha_g2
The alpha-weighted N-array of input primal 2-forms g.
Definition: ATSolver2D.h:162
std::vector< PrimalForm2 > u2
The N-array of regularized primal 2-forms u.
Definition: ATSolver2D.h:164
PrimalIdentity2 alpha_Id2
The alpha-weighted identity operator for primal 2-forms (stored for performance)
Definition: ATSolver2D.h:156
int verbose
Tells the verbose level.
Definition: ATSolver2D.h:191
PrimalDerivative1 M12
The primal edge to face average operator.
Definition: ATSolver2D.h:152
PrimalForm0 l_1_over_4e
The primal 0-form lambda/(4epsilon) (stored for performance)
Definition: ATSolver2D.h:171
The antiderivative of primal 2-forms.
Definition: ATSolver2D.h:154
PrimalForm0 v0
Definition: ATSolver2D.h:167
void initOperators()
Initializes the operators.
Definition: ATSolver2D.h:976
Surfel2IndexMap surfel2idx
Definition: ATSolver2D.h:175
PrimalDerivative1 primal_D1
the derivative operator for primal 1-forms
Definition: ATSolver2D.h:148
PrimalIdentity0 l_1_over_4e_Id0
The 1/(4epsilon)-weighted identity operator for primal 0-forms (stored for performance)
Definition: ATSolver2D.h:158
PrimalForm0 former_v0
The primal 0-form v at the previous iteration.
Definition: ATSolver2D.h:169
std::vector< PrimalForm2 > g2
The N-array of input primal 2-forms g.
Definition: ATSolver2D.h:160
PrimalDerivative0 primal_D0
the derivative operator for primal 0-forms
Definition: ATSolver2D.h:146
CountedConstPtrOrConstPtr< Calculus > ptrCalculus
A smart (or not) pointer to a calculus object.
Definition: ATSolver2D.h:144
PrimalDerivative0 M01
The primal vertex to edge average operator.
Definition: ATSolver2D.h:150
std::ostream & info()
SMesh::Index Index
Trace trace
Definition: Common.h:153
@ PRIMAL
Definition: Duality.h:61
unsigned int index(DGtal::uint32_t n, unsigned int b)
Definition: testBits.cpp:44

## ◆ ATSolver2D() [2/4]

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
 DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::ATSolver2D ( )
delete

Default constructor.

## ◆ ~ATSolver2D()

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
 DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::~ATSolver2D ( )
default

Destructor.

## ◆ ATSolver2D() [3/4]

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
 DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::ATSolver2D ( const ATSolver2D< TKSpace, TLinearAlgebra > & other )
default

Copy constructor.

Parameters
 other the object to clone.

## ◆ ATSolver2D() [4/4]

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
 DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::ATSolver2D ( ATSolver2D< TKSpace, TLinearAlgebra > && other )
default

Move constructor.

Parameters
 other the object to move.

## ◆ checkV0()

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
 std::tuple DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::checkV0 ( ) const
inline

Debug method for checking if v is a scalar field between 0 and 1.

Returns
the tuple (min(v), average(v), max(v))

Definition at line 723 of file ATSolver2D.h.

724  {
725  const double m1 = v0.myContainer.minCoeff();
726  const double m2 = v0.myContainer.maxCoeff();
727  const double ma = v0.myContainer.mean();
728  if ( verbose >= 1 )
729  trace.info() << "0-form v (should be in [0,1]): min=" << m1 << " avg=" << ma << " max=" << m2 << std::endl;
730  return std::make_tuple( m1, m2, ma );
731  }

## ◆ diffV0()

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
 std::tuple DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::diffV0 ( ) const
inline

Computes the norms loo, l2, l1 of (v - former_v), i.e. the evolution of discontinuity function v.

Returns
a tuple (n_infty,n_2,n_1) giving the loo/l2/l1-norm of (v - former_v)

Definition at line 702 of file ATSolver2D.h.

703  {
704  PrimalForm0 delta = v0 - former_v0;
705  delta.myContainer = delta.myContainer.cwiseAbs();
706  const double n_oo = delta.myContainer.maxCoeff();
707  const double n_2 = std::sqrt(delta.myContainer.squaredNorm()/delta.myContainer.size());
708  const double n_1 = delta.myContainer.mean();
709  return std::make_tuple( n_oo, n_2, n_1 );
710  }
Calculus::PrimalForm0 PrimalForm0
Definition: ATSolver2D.h:117
Container myContainer
Definition: KForm.h:131

## ◆ getOutputScalarFieldU2()

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
template<typename ScalarFieldOutput , typename SurfelRangeConstIterator >
 void DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::getOutputScalarFieldU2 ( ScalarFieldOutput & output, SurfelRangeConstIterator itB, SurfelRangeConstIterator itE )
inline

Given a range of surfels [itB,itE), returns in output the regularized scalar field u.

Template Parameters
 ScalarFieldOutput the type of scalar field for output values (RandomAccess container) SurfelRangeConstIterator the type of iterator for traversing a range of surfels
Parameters
 [out] output the vector of output values (a scalar field), which should be of size length(itB,itE) itB the start of the range of surfels. itE past the end of the range of surfels.

Definition at line 798 of file ATSolver2D.h.

800  {
801  ASSERT( u2.size() == 1 && "[ATSolver2D::getOutputScalarFieldU2] "
802  "You try to output a scalar field from a vector field." );
803  Index i = 0;
804  for ( auto it = itB; it != itE; ++it, ++i )
805  {
806  Index idx = surfel2idx[ *it ];
807  output[ i ] = u2[ 0 ].myContainer( idx );
808  }
809  }

## ◆ getOutputScalarFieldV0()

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
template<typename ScalarFieldOutput , typename CellRangeConstIterator >
 void DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::getOutputScalarFieldV0 ( ScalarFieldOutput & output, CellRangeConstIterator itB, CellRangeConstIterator itE, CellOutputPolicy policy = CellOutputPolicy::Average )
inline

Given a range of pointels, linels or 2-cells [itB,itE), returns in output the feature vector v (the average of v for linels/surfels).

Template Parameters
 ScalarFieldOutput the type of scalar field for output values (RandomAccess container) CellRangeConstIterator the type of iterator for traversing a range of cells
Parameters
 [out] output the vector of output scalar values (a scalar field), which should be of size length(itB,itE) [in] itB the start of the range of cells. [in] itE past the end of the range of cells. [in] policy the chosen policy for outputing v values for a given cell.

Definition at line 824 of file ATSolver2D.h.

827  {
828  const KSpace& K = ptrCalculus->myKSpace;
829  const Dimension k = K.uDim( *itB );
830  ASSERT( k <= 2 );
831  Index i = 0;
832  if ( k == 0 )
833  {
834  for ( auto it = itB; it != itE; ++it, ++i )
835  {
836  const Cell pointel = *it;
837  const Index idx = ptrCalculus->getCellIndex( pointel );
838  output[ i ] = v0.myContainer( idx );
839  }
840  }
841  else if ( k == 1 )
842  {
843  for ( auto it = itB; it != itE; ++it, ++i )
844  {
845  const Cell linel = *it;
846  const Dimension d = * K.uDirs( linel );
847  const Cell p0 = K.uIncident( linel, d, false );
848  const Cell p1 = K.uIncident( linel, d, true );
849  const Index idx0 = ptrCalculus->getCellIndex( p0 );
850  const Index idx1 = ptrCalculus->getCellIndex( p1 );
851  switch (policy) {
852  case CellOutputPolicy::Average: output[ i ] = 0.5 * ( v0.myContainer( idx0 ) + v0.myContainer( idx1 ) );
853  break;
854  case CellOutputPolicy::Minimum: output[ i ] = std::min( v0.myContainer( idx0 ), v0.myContainer( idx1 ) );
855  break;
856  case CellOutputPolicy::Maximum: output[ i ] = std::max( v0.myContainer( idx0 ), v0.myContainer( idx1 ) );
857  break;
858  }
859  }
860  }
861  else if ( k == 2 )
862  {
863  for ( auto it = itB; it != itE; ++it, ++i )
864  {
865  const Cell face = *it;
866  const Dimension d = * K.uDirs( face );
867  const Cell l0 = K.uIncident( face, d, false );
868  const Cell l1 = K.uIncident( face, d, true );
869  const Dimension j = * K.uDirs( l0 );
870  const Cell p00 = K.uIncident( l0, j, false );
871  const Cell p01 = K.uIncident( l0, j, true );
872  const Cell p10 = K.uIncident( l1, j, false );
873  const Cell p11 = K.uIncident( l1, j, true );
874  const Index idx00 = ptrCalculus->getCellIndex( p00 );
875  const Index idx01 = ptrCalculus->getCellIndex( p01 );
876  const Index idx10 = ptrCalculus->getCellIndex( p10 );
877  const Index idx11 = ptrCalculus->getCellIndex( p11 );
878  switch (policy) {
879  case CellOutputPolicy::Average:
880  output[ i ] = 0.25 * ( v0.myContainer( idx00 ) + v0.myContainer( idx01 )
881  + v0.myContainer( idx10 ) + v0.myContainer( idx11 ) );
882  break;
883  case CellOutputPolicy::Minimum:
884  output[ i ] = std::min( std::min( v0.myContainer( idx00 ), v0.myContainer( idx01 ) ),
885  std::min( v0.myContainer( idx10 ), v0.myContainer( idx11 ) ) );
886  break;
887  case CellOutputPolicy::Maximum:
888  output[ i ] = std::max( std::max( v0.myContainer( idx00 ), v0.myContainer( idx01 ) ),
889  std::max( v0.myContainer( idx10 ), v0.myContainer( idx11 ) ) );
890  break;
891  }
892  }
893  }
894  }
DGtal::uint32_t Dimension
Definition: Common.h:136
int max(int a, int b)
KSpace K
KSpace::Cell Cell

References K, max(), and DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::ptrCalculus.

## ◆ getOutputVectorFieldU2()

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
template<typename VectorFieldOutput , typename SurfelRangeConstIterator >
 void DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::getOutputVectorFieldU2 ( VectorFieldOutput & output, SurfelRangeConstIterator itB, SurfelRangeConstIterator itE )
inline

Given a range of surfels [itB,itE), returns in output the regularized vector field u.

Template Parameters
 VectorFieldOutput the type of vector field for output values (RandomAccess container) SurfelRangeConstIterator the type of iterator for traversing a range of surfels
Parameters
 [out] output the vector of output values (a scalar or vector field depending on input). itB the start of the range of surfels. itE past the end of the range of surfels.

Definition at line 772 of file ATSolver2D.h.

774  {
775  const Dimension N = u2.size();
776  Index i = 0;
777  for ( auto it = itB; it != itE; ++it, ++i )
778  {
779  Index idx = surfel2idx[ *it ];
780  ASSERT( output[ i ].size() >= N );
781  for ( Dimension k = 0; k < N; ++k )
782  output[ i ][ k ] = u2[ k ].myContainer( idx );
783  }
784  }
Index size(const int order) const
Definition: ATSolver2D.h:258

Referenced by DGtal::ShortcutsGeometry< TKSpace >::getATVectorFieldApproximation(), and main().

## ◆ getU2()

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
 PrimalForm2 DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::getU2 ( Dimension k ) const
inline
Parameters
 k an integer such that 0 <= k < u2.size()
Returns
the k-th piecewise smooth function u as a primal 2-form.

Definition at line 754 of file ATSolver2D.h.

755  {
756  return u2[ k ];
757  }

## ◆ getV0()

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
 PrimalForm0 DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::getV0 ( ) const
inline
Returns
the discontinuity function v as a primal 0-form.

Definition at line 735 of file ATSolver2D.h.

736  {
737  return v0;
738  }

## ◆ getV1()

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
 PrimalForm1 DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::getV1 ( ) const
inline
Returns
the discontinuity function v as a primal 1-form.

Definition at line 741 of file ATSolver2D.h.

742  {
743  return M01*v0;
744  }

## ◆ getV2()

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
 PrimalForm2 DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::getV2 ( ) const
inline
Returns
the discontinuity function u as a primal 2-form.

Definition at line 747 of file ATSolver2D.h.

748  {
749  return M12*M01*v0;
750  }

## ◆ initInputScalarFieldU2() [1/2]

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
template<typename ScalarFieldInput , typename SurfelRangeConstIterator >
 void DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::initInputScalarFieldU2 ( const ScalarFieldInput & input, SurfelRangeConstIterator itB, SurfelRangeConstIterator itE )
inline

Given a range of surfels [itB,itE) and an input scalar field, initializes the AT 2-forms u and g. The 0-form v is itself initialized to 1 everywhere.

Template Parameters
 ScalarFieldInput the type of scalar field for input values (RandomAccess container) SurfelRangeConstIterator the type of iterator for traversing a range of surfels
Parameters
 [in] input the input scalar field (a vector of scalar values) itB the start of the range of surfels. itE past the end of the range of surfels.

Definition at line 330 of file ATSolver2D.h.

332  {
333  if ( verbose >= 1 )
334  trace.beginBlock( "[ATSolver2D::initInputVectorFieldU2] Initializing input data" );
335  ASSERT( ! input.empty() );
336  if ( verbose >= 2 ) trace.info() << "input g as one 2-form." << std::endl;
337  g2 = std::vector<PrimalForm2>( 1, PrimalForm2( *ptrCalculus ) );
338  alpha_g2 = std::vector<PrimalForm2>( 1, PrimalForm2( *ptrCalculus ) );
339  Index i = 0;
340  for ( auto it = itB; it != itE; ++it, ++i )
341  {
342  Index idx = surfel2idx[ *it ];
343  g2[ 0 ].myContainer( idx ) = input[ i ];
344  }
345  // u = g at the beginning
346  if ( verbose >= 2 )
347  trace.info() << "Unknown u[:] = g[:] at beginning." << std::endl;
348  u2 = g2;
349  // v = 1 at the beginning
350  if ( verbose >= 2 ) trace.info() << "Unknown v = 1" << std::endl;
352  if ( verbose >= 1 ) trace.endBlock();
353  normalize_u2 = false;
354  }
Calculus::PrimalForm2 PrimalForm2
Definition: ATSolver2D.h:119
bool normalize_u2
Indicates whether to normalize U (unit norm) at each iteration or not.
Definition: ATSolver2D.h:189
static KForm< TCalculus, order, duality > ones(ConstAlias< Calculus > calculus)
void beginBlock(const std::string &keyword="")
double endBlock()

## ◆ initInputScalarFieldU2() [2/2]

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
template<typename Scalar >
 Index DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::initInputScalarFieldU2 ( const std::map< Surfel, Scalar > & input )
inline

Given a map Surfel -> Scalar, initializes forms g, u and v of the AT solver. Note that there is only one 2-form u/g.

Parameters
 input any map Surfel -> Scalar
Returns
the number of cells that were initialized.
Template Parameters
 Scalar any type representing a scalar (float, double)

Definition at line 411 of file ATSolver2D.h.

412  {
413  if ( verbose >= 1 ) trace.beginBlock( "[ATSolver2D::initScalarInput] Initializing input data" );
414  if ( verbose >= 2 ) trace.info() << "discontinuity 0-form v = 1." << std::endl;
416  if ( verbose >= 2 ) trace.info() << "input g as one 2-form." << std::endl;
417  g2 = std::vector<PrimalForm2>( 1, PrimalForm2( *ptrCalculus ) );
418  alpha_g2 = std::vector<PrimalForm2>( 1, PrimalForm2( *ptrCalculus ) );
419  const Scalar zero;
420  Index nbok = 0;
421  for ( Index index = 0; index < size(2); index++)
422  {
423  const SCell& cell = g2[ 0 ].getSCell( index );
424  const auto it = input.find( cell );
425  const auto n = ( it != input.end() ) ? *it : zero;
426  nbok += ( it != input.end() ) ? 1 : 0;
427  g2[ 0 ].myContainer( index ) = n;
428  }
429  // u = g at the beginning
430  if ( verbose >= 2 ) trace.info() << "Unknown u[:] = g[:] at beginning." << std::endl;
431  u2 = g2;
432  // v = 1 at the beginning
433  if ( verbose >= 2 ) trace.info() << "Unknown v = 1" << std::endl;
435  if ( verbose >= 1 ) trace.endBlock();
436  normalize_u2 = false;
437  return nbok;
438  }
Z3i::SCell SCell

## ◆ initInputVectorFieldU2() [1/2]

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
template<typename ScalarVector >
 Index DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::initInputVectorFieldU2 ( const std::map< Surfel, ScalarVector > & input, bool normalize = false )
inline

Given a map Surfel -> ScalarVector, initializes forms g, u and v of the AT solver. Note that there are as many 2-forms u/g as the number of dimensions of ScalarVector.

Parameters
 input any map Surfel -> ScalarVector normalize when 'true', the input is supposed to be a unit vector field and the solver will output a unit regularized vector field at the end of each minimization step.
Returns
the number of cells that were initialized.
Template Parameters
 ScalarVector any type representing a vector/array of scalars (float, double)

Definition at line 371 of file ATSolver2D.h.

373  {
374  if ( verbose >= 1 ) trace.beginBlock( "[ATSolver2D::initVectorInput] Initializing input data" );
375  if ( verbose >= 2 ) trace.info() << "discontinuity 0-form v = 1." << std::endl;
376  const Dimension N = ScalarVector().size();
377  if ( verbose >= 2 ) trace.info() << "input g as " << N << " 2-forms." << std::endl;
378  g2 = std::vector<PrimalForm2>( N, PrimalForm2( *ptrCalculus ) );
379  alpha_g2 = std::vector<PrimalForm2>( N, PrimalForm2( *ptrCalculus ) );
380  const ScalarVector zero;
381  Index nbok = 0;
382  for ( Index index = 0; index < size(2); index++)
383  {
384  const SCell& cell = g2[ 0 ].getSCell( index );
385  const auto it = input.find( cell );
386  const auto n = ( it != input.end() ) ? it->second : zero;
387  nbok += ( it != input.end() ) ? 1 : 0;
388  for ( Dimension k = 0; k < N; ++k )
389  g2[ k ].myContainer( index ) = n[ k ];
390  }
391  // u = g at the beginning
392  if ( verbose >= 2 )
393  trace.info() << "Unknown u[:] = g[:] at beginning." << std::endl;
394  u2 = g2;
395  // v = 1 at the beginning
396  if ( verbose >= 2 ) trace.info() << "Unknown v = 1" << std::endl;
398  if ( verbose >= 1 ) trace.endBlock();
399  normalize_u2 = normalize;
400  return nbok;
401  }

## ◆ initInputVectorFieldU2() [2/2]

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
template<typename VectorFieldInput , typename SurfelRangeConstIterator >
 void DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::initInputVectorFieldU2 ( const VectorFieldInput & input, SurfelRangeConstIterator itB, SurfelRangeConstIterator itE, bool normalize = false )
inline

Given a range of surfels [itB,itE) and an input vector field, initializes the AT 2-forms u and g. The 0-form v is itself initialized to 1 everywhere.

Template Parameters
 VectorFieldInput the type of vector field for input values (RandomAccess container) SurfelRangeConstIterator the type of iterator for traversing a range of surfels
Parameters
 [in] input the input vector field (a vector of vector values) [in] itB the start of the range of surfels. [in] itE past the end of the range of surfels. [in] normalize when 'true', the input is supposed to be a unit vector field and the solver will output a unit regularized vector field at the end of each minimization step.

Definition at line 288 of file ATSolver2D.h.

291  {
292  if ( verbose >= 1 )
293  trace.beginBlock( "[ATSolver2D::initInputVectorFieldU2] Initializing input data" );
294  ASSERT( ! input.empty() );
295  const Dimension N = input[ 0 ].size();
296  if ( verbose >= 2 ) trace.info() << "input g as " << N << " 2-forms." << std::endl;
297  g2 = std::vector<PrimalForm2>( N, PrimalForm2( *ptrCalculus ) );
298  alpha_g2 = std::vector<PrimalForm2>( N, PrimalForm2( *ptrCalculus ) );
299  Index i = 0;
300  for ( auto it = itB; it != itE; ++it, ++i )
301  {
302  Index idx = surfel2idx[ *it ];
303  for ( Dimension k = 0; k < N; ++k )
304  g2[ k ].myContainer( idx ) = input[ i ][ k ];
305  }
306  // u = g at the beginning
307  if ( verbose >= 2 )
308  trace.info() << "Unknown u[:] = g[:] at beginning." << std::endl;
309  u2 = g2;
310  // v = 1 at the beginning
311  if ( verbose >= 2 ) trace.info() << "Unknown v = 1" << std::endl;
313  if ( verbose >= 1 ) trace.endBlock();
314  normalize_u2 = normalize;
315  }

Referenced by DGtal::ShortcutsGeometry< TKSpace >::getATVectorFieldApproximation(), and main().

## ◆ initOperators()

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
 void DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::initOperators ( )
inlineprotected

Initializes the operators.

Definition at line 976 of file ATSolver2D.h.

977  {
978  if ( verbose >= 1 ) trace.beginBlock( "[ATSolver2D::initOperators] Solver initialization" );
979  if ( verbose >= 2 ) trace.info() << "derivative of primal 0-forms: primal_D0" << std::endl;
980  primal_D0 = ptrCalculus->template derivative<0,PRIMAL>();
981  if ( verbose >= 2 ) trace.info() << "derivative of primal 1-forms: primal_D1" << std::endl;
982  primal_D1 = ptrCalculus->template derivative<1,PRIMAL>();
983  if ( verbose >= 2 ) trace.info() << "antiderivative of primal 2-forms: primal_AD2" << std::endl;
984  primal_AD2 = ptrCalculus->template antiderivative<2,PRIMAL>();
985  if ( verbose >= 2 ) trace.info() << "vertex to edge average operator: M01" << std::endl;
986  M01 = primal_D0;
987  M01.myContainer = .5 * M01.myContainer.cwiseAbs();
988  if ( verbose >= 2 ) trace.info() << "edge to face average operator: M12" << std::endl;
989  M12 = primal_D1;
990  M12.myContainer = .25 * M12.myContainer.cwiseAbs();
991  if ( verbose >= 1 ) trace.endBlock();
992  }

## ◆ isValid()

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
 bool DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::isValid ( ) const
inline

Checks the validity/consistency of the object.

Returns
'true' if the object is valid, 'false' otherwise.

Definition at line 962 of file ATSolver2D.h.

963  {
964  return true;
965  }

## ◆ normalizeU2()

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
 void DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::normalizeU2 ( )
inline

Forces the normalization of the vector u, meaning for all index i, $$\sum_{k=0}^{K-1} u[k][i]^2 = 1$$. Can be useful in some applications where you are looking for unitary vector field.

Definition at line 684 of file ATSolver2D.h.

685  {
686  for ( Index index = 0; index < size( 2 ); index++)
687  {
688  double n2 = 0.0;
689  for ( unsigned int d = 0; d < u2.size(); ++d )
690  n2 += u2[ d ].myContainer( index ) * u2[ d ].myContainer( index );
691  double norm = sqrt( n2 );
692  if (norm == 0.0) continue;
693  for ( unsigned int d = 0; d < u2.size(); ++d )
694  u2[ d ].myContainer( index ) /= norm;
695  }
696  }

## ◆ operator=() [1/2]

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
 ATSolver2D& DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::operator= ( ATSolver2D< TKSpace, TLinearAlgebra > && other )
default

Move assignment operator.

Parameters
 other the object to move.
Returns
a reference on 'this'.

## ◆ operator=() [2/2]

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
 ATSolver2D& DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::operator= ( const ATSolver2D< TKSpace, TLinearAlgebra > & other )
default

Copy assignment operator.

Parameters
 other the object to copy.
Returns
a reference on 'this'.

## ◆ selfDisplay()

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
 void DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::selfDisplay ( std::ostream & out ) const
inline

Writes/Displays the object on an output stream.

Parameters
 out the output stream where the object is written.

Definition at line 949 of file ATSolver2D.h.

950  {
951  auto cv = checkV0();
952  out << "[ATSolver2D] v is between min/avg/max:"
953  << std::get<0>(cv) << "/"
954  << std::get<1>(cv) << "/"
955  << std::get<2>(cv) << std::endl;
956  }
std::tuple< double, double, double > checkV0() const
Definition: ATSolver2D.h:723

## ◆ setEpsilon()

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
 void DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::setEpsilon ( double e )
inline

Initializes the epsilon parameter of AT and precomputes the assaociated forms and operators.

Parameters
 e the epsilon parameter in AT

Definition at line 525 of file ATSolver2D.h.

526  {
527  epsilon = e;
529  l_1_over_4e_Id0 = (lambda/4./epsilon)*ptrCalculus->template identity<0, PRIMAL>();
530  }

## ◆ setUp() [1/3]

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
 void DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::setUp ( double a, double l )
inline

Initializes the alpha and lambda parameters of AT.

Parameters
 a the global alpha parameter l the global lambda parameter
Note
all 2-cells have the same weight for the data term.

Definition at line 444 of file ATSolver2D.h.

445  {
446  const Dimension N = (Dimension)g2.size();
447  alpha = a;
448  lambda = l;
449  alpha_Id2 = alpha * ptrCalculus->template identity<2, PRIMAL>();
450  for ( Dimension k = 0; k < N; ++k )
451  alpha_g2[ k ] = alpha * g2[ k ];
452  }

## ◆ setUp() [2/3]

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
template<typename AlphaWeights , typename SurfelRangeConstIterator >
 void DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::setUp ( double a, double l, const AlphaWeights & weights, SurfelRangeConstIterator itB, SurfelRangeConstIterator itE )
inline

Initializes the alpha and lambda parameters of AT, with weights on the 2-cells for the data terms.

Template Parameters
 AlphaWeights the type of RandomAccess container for alpha weight values SurfelRangeConstIterator the type of iterator for traversing a range of surfels
Parameters
 [in] a the global alpha parameter [in] l the global lambda parameter [in] weights the vector of alpha weights for each surfel of the range [itB,itE) [in] itB the start of the range of surfels. [in] itE past the end of the range of surfels.
Note
Useful for inpainting applications or for adaptive piecewise smooth reconstruction.

Definition at line 499 of file ATSolver2D.h.

502  {
503  const Dimension N = g2.size();
504  alpha = a;
505  lambda = l;
506  PrimalForm2 w_form( *ptrCalculus );
507  if ( verbose >= 2 )
508  trace.info() << "Using variable weights for fitting (alpha term)" << std::endl;
509  for ( Dimension k = 0; k < N; ++k )
510  alpha_g2[ k ] = alpha * g2[ k ];
511  Index i = 0;
512  for ( auto it = itB; it != itE; ++it, ++i )
513  {
514  const Index idx = surfel2idx[ *it ];
515  const Scalar w = weights[ i ];
516  w_form.myContainer( idx ) = w;
517  for ( Dimension k = 0; k < N; ++k )
518  alpha_g2[ k ].myContainer( idx ) *= w;
519  }
520  alpha_Id2 = alpha * diagonal( w_form );
521  }
DGtal::LinearOperator< Calculus, dim, duality, dim, duality > diagonal(const DGtal::KForm< Calculus, dim, duality > &kform)

## ◆ setUp() [3/3]

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
 void DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::setUp ( double a, double l, const std::map< Surfel, Scalar > & weights )
inline

Initializes the alpha and lambda parameters of AT, with weights on the 2-cells for the data terms.

Parameters
 a the global alpha parameter l the global lambda parameter weights the map Surfel -> Scalar that gives the weight of each 2-cell in the data terms.
Note
Useful for inpainting applications or for adaptive piecewise smooth reconstruction.

Definition at line 463 of file ATSolver2D.h.

464  {
465  const Dimension N = g2.size();
466  alpha = a;
467  lambda = l;
468  PrimalForm2 w_form( *ptrCalculus );
469  if ( verbose >= 2 )
470  trace.info() << "Using variable weights for fitting (alpha term)" << std::endl;
471  for ( Dimension k = 0; k < N; ++k )
472  alpha_g2[ k ] = alpha * g2[ k ];
473  for ( Index index = 0; index < size( 2 ); index++)
474  {
475  const SCell& cell = g2[ 0 ].getSCell( index );
476  const Scalar& w = weights[ cell ];
477  w_form.myContainer( index ) = w;
478  for ( Dimension k = 0; k < N; ++k )
479  alpha_g2[ k ].myContainer( index ) *= w;
480  }
481  alpha_Id2 = alpha * diagonal( w_form );
482  }

## ◆ size()

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
 Index DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::size ( const int order ) const
inline
Parameters
 order the dimension of cells (0,1,2)
Returns
the number of cells with dimension order

Definition at line 258 of file ATSolver2D.h.

259  {
260  return ptrCalculus->kFormLength(order, PRIMAL);
261  }

## ◆ solveForEpsilon()

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
 bool DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::solveForEpsilon ( double eps, double n_oo_max = 1e-4, unsigned int iter_max = 10 )
inline

Solves the alternate minimization of AT for a given eps. Solves for u then for v till convergence.

Parameters
 eps the epsilon parameter at which AT is solved. n_oo_max the alternate minimization will stop when the loo-norm of $$v^{k+1} - v^k$$ is below this bound. iter_max the alternate minimization will stop when the number of minimization steps exceeds iter_max.
Returns
true if everything went fine, false if there was a problem in the optimization.
Note
Use diffV0 to check if you are close to a critical point of AT.

Definition at line 605 of file ATSolver2D.h.

608  {
609  (void)n_oo_max;//paramerter not used
610
611  bool ok = true;
612  if ( verbose >= 1 ) {
613  std::ostringstream sstr;
614  sstr << "******* Solving AT for epsilon = " << eps << " **********";
615  trace.beginBlock( sstr.str() );
616  }
617  setEpsilon( eps );
618  for ( unsigned int i = 0; i < iter_max; ++i )
619  {
620  if ( verbose >= 1 )
621  trace.info() << "---------- Iteration "
622  << i << "/" << iter_max << " ---------------" << std::endl;
624  auto diffs_v = diffV0();
625  if ( verbose >= 1 ) {
626  trace.info() << "Variation |v^k+1 - v^k|_oo = " << std::get<0>( diffs_v )
627  << std::endl;
628  if ( verbose >= 2 ) {
629  trace.info() << "Variation |v^k+1 - v^k|_2 = " << std::get<1>( diffs_v )
630  << std::endl;
631  trace.info() << "Variation |v^k+1 - v^k|_1 = " << std::get<2>( diffs_v )
632  << std::endl;
633  }
634  }
635  if ( std::get<0>( diffs_v ) < 1e-4 ) break;
636  }
637  if ( verbose >= 1 ) trace.endBlock();
638  return ok;
639  }
bool solveOneAlternateStep()
Definition: ATSolver2D.h:546
void setEpsilon(double e)
Definition: ATSolver2D.h:525
std::tuple< double, double, double > diffV0() const
Definition: ATSolver2D.h:702

## ◆ solveGammaConvergence()

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
 bool DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::solveGammaConvergence ( double eps1 = 2.0, double eps2 = 0.25, double epsr = 2.0, bool compute_smallest_epsilon_map = false, double n_oo_max = 1e-4, unsigned int iter_max = 10 )
inline

Solves AT by progressively decreasing epsilon from eps1 to eps2. AT is solved with solveForEpsilon at each epsilon.

Parameters
 eps1 the first epsilon parameter at which AT is solved. eps2 the last epsilon parameter at which AT is solved. epsr the ratio (>1) used to decrease progressively epsilon. compute_smallest_epsilon_map when 'true' determines for each surfel the smallest epsilon for which it is a discontinuity. n_oo_max the alternate minimization will stop when the loo-norm of $$v^{k+1} - v^k$$ is below this bound. iter_max the alternate minimization will stop when the number of minimization steps exceeds iter_max.
Returns
true if everything went fine, false if there was a problem in the optimization.

Definition at line 657 of file ATSolver2D.h.

663  {
664  bool ok = true;
665  if ( epsr <= 1.0 ) epsr = 2.0;
666  if ( verbose >= 1 )
667  trace.beginBlock( "#### Solve AT by Gamma-convergence ##########" );
668  if ( compute_smallest_epsilon_map ) smallest_epsilon_map.clear();
669  for ( double eps = eps1; eps >= eps2; eps /= epsr )
670  {
671  solveForEpsilon( eps, n_oo_max, iter_max );
672  if ( compute_smallest_epsilon_map )
674  }
675  if ( verbose >= 1 )
676  trace.endBlock();
677  return ok;
678  }
bool solveForEpsilon(double eps, double n_oo_max=1e-4, unsigned int iter_max=10)
Definition: ATSolver2D.h:605
void updateSmallestEpsilonMap(const double threshold=.5)
Definition: ATSolver2D.h:903
SmallestEpsilonMap smallest_epsilon_map
Definition: ATSolver2D.h:178

## ◆ solveOneAlternateStep()

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
 bool DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::solveOneAlternateStep ( )
inline

Solves one step of the alternate minimization of AT. Solves for u then for v.

Returns
true if everything went fine, false if there was a problem in the optimization.
Note
Use diffV0 to check if you are close to a critical point of AT.

Definition at line 546 of file ATSolver2D.h.

547  {
548  bool solve_ok = true;
549  if ( verbose >= 1 ) trace.beginBlock("Solving for u as a 2-form");
550  PrimalForm1 v1_squared = M01*v0;
551  v1_squared.myContainer.array() = v1_squared.myContainer.array().square();
552  const PrimalIdentity2 ope_u2 = alpha_Id2
553  + primal_AD2.transpose() * dec_helper::diagonal( v1_squared ) * primal_AD2;
554
555  if ( verbose >= 2 ) trace.info() << "Prefactoring matrix U associated to u" << std::endl;
556  SolverU2 solver_u2;
557  solver_u2.compute( ope_u2 );
558  for ( Dimension d = 0; d < u2.size(); ++d )
559  {
560  if ( verbose >= 2 ) trace.info() << "Solving U u[" << d << "] = a g[" << d << "]" << std::endl;
561  u2[ d ] = solver_u2.solve( alpha_g2[ d ] );
562  if ( verbose >= 2 ) trace.info() << " => " << ( solver_u2.isValid() ? "OK" : "ERROR" )
563  << " " << solver_u2.myLinearAlgebraSolver.info() << std::endl;
564  solve_ok = solve_ok && solver_u2.isValid();
565  }
566  if ( normalize_u2 ) normalizeU2();
567  if ( verbose >= 1 ) trace.endBlock();
568  if ( verbose >= 1 ) trace.beginBlock("Solving for v");
569  former_v0 = v0;
570  PrimalForm1 squared_norm_d_u2 = PrimalForm1::zeros(*ptrCalculus);
571  for ( Dimension d = 0; d < u2.size(); ++d )
572  squared_norm_d_u2.myContainer.array() += (primal_AD2 * u2[ d ] ).myContainer.array().square();
573  trace.info() << "build metric u2" << std::endl;
574  const PrimalIdentity0 ope_v0 = l_1_over_4e_Id0
576  + M01.transpose() * dec_helper::diagonal( squared_norm_d_u2 ) * M01;
577
578  if ( verbose >= 2 ) trace.info() << "Prefactoring matrix V associated to v" << std::endl;
579  SolverV0 solver_v0;
580  solver_v0.compute( ope_v0 );
581  if ( verbose >= 2 ) trace.info() << "Solving V v = l/4e * 1" << std::endl;
582  v0 = solver_v0.solve( l_1_over_4e );
583  if ( verbose >= 2 ) trace.info() << " => " << ( solver_v0.isValid() ? "OK" : "ERROR" )
584  << " " << solver_v0.myLinearAlgebraSolver.info() << std::endl;
585  solve_ok = solve_ok && solver_v0.isValid();
586  if ( verbose >= 1 ) trace.endBlock();
587  return solve_ok;
588  }
Calculus::PrimalForm1 PrimalForm1
Definition: ATSolver2D.h:118
Calculus::PrimalIdentity2 PrimalIdentity2
Definition: ATSolver2D.h:122
Calculus::PrimalIdentity0 PrimalIdentity0
Definition: ATSolver2D.h:120
DiscreteExteriorCalculusSolver< Calculus, LinearAlgebraSolver, 0, PRIMAL, 0, PRIMAL > SolverV0
Definition: ATSolver2D.h:140
DiscreteExteriorCalculusSolver< Calculus, LinearAlgebraSolver, 2, PRIMAL, 2, PRIMAL > SolverU2
Definition: ATSolver2D.h:139
DiscreteExteriorCalculusSolver & compute(const Operator &linear_operator)
static KForm< TCalculus, order, duality > zeros(ConstAlias< Calculus > calculus)
TransposedLinearOperator transpose() const
DGtal::LinearOperator< Calculus, dim, duality, dim, duality > diagonal(const DGtal::KForm< Calculus, dim, duality > &kform)
Definition: DECHelpers.h:60

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
 void DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::updateSmallestEpsilonMap ( const double threshold = .5 )
inline

Computes the map that stores for each surfel the smallest epsilon for which the surfel was in the discontinuity zone (more precisely, the surfel has at least two vertices that belongs to the set of discontinuity).

Parameters
 [in] threshold the threshold for discontinuity function v (below u is discontinuous, above u is continuous)

Definition at line 903 of file ATSolver2D.h.

904  {
905  const KSpace& K = ptrCalculus->myKSpace;
906  for ( const SCell surfel : ptrCalculus->template getIndexedSCells<2, PRIMAL>() )
907  {
908  const Cell face = K.unsigns( surfel );
909  const Dimension k1 = * K.uDirs( face );
910  const Cell l0 = K.uIncident( face, k1, false );
911  const Cell l1 = K.uIncident( face, k1, true );
912  const Dimension k2 = * K.uDirs( l0 );
913  const Cell ll0 = K.uIncident( face, k2, false );
914  const Cell ll1 = K.uIncident( face, k2, true );
915  const Cell p00 = K.uIncident( l0, k2, false );
916  const Cell p01 = K.uIncident( l0, k2, true );
917  const Cell p10 = K.uIncident( l1, k2, false );
918  const Cell p11 = K.uIncident( l1, k2, true );
919
920  std::vector<double> features( 4 );
921  features[ 0 ] = v0.myContainer( ptrCalculus->getCellIndex( p00 ) );
922  features[ 1 ] = v0.myContainer( ptrCalculus->getCellIndex( p01 ) );
923  features[ 2 ] = v0.myContainer( ptrCalculus->getCellIndex( p10 ) );
924  features[ 3 ] = v0.myContainer( ptrCalculus->getCellIndex( p11 ) );
925  std::sort( features.begin(), features.end() );
926
927  if ( features[ 1 ] <= threshold )
928  {
929  auto it = smallest_epsilon_map.find( surfel );
930  if ( it != smallest_epsilon_map.end() )
931  it->second = std::min( epsilon, it->second );
932  else smallest_epsilon_map[ surfel ] = epsilon;
933  }
934  }
935  }

## ◆ alpha

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
 double DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::alpha

The global coefficient alpha giving the smoothness of the reconstruction (the smaller, the smoother)

Definition at line 181 of file ATSolver2D.h.

## ◆ alpha_g2

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
 std::vector DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::alpha_g2
protected

The alpha-weighted N-array of input primal 2-forms g.

Definition at line 162 of file ATSolver2D.h.

## ◆ alpha_Id2

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
 PrimalIdentity2 DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::alpha_Id2
protected

The alpha-weighted identity operator for primal 2-forms (stored for performance)

Definition at line 156 of file ATSolver2D.h.

## ◆ dimension

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
 const Dimension DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::dimension = KSpace::dimension
static

Definition at line 99 of file ATSolver2D.h.

## ◆ epsilon

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
 double DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::epsilon

The global coefficient epsilon giving the width of the discontinuities (the smaller, the thinner)

Definition at line 187 of file ATSolver2D.h.

## ◆ former_v0

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
 PrimalForm0 DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::former_v0
protected

The primal 0-form v at the previous iteration.

Definition at line 169 of file ATSolver2D.h.

## ◆ g2

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
 std::vector DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::g2
protected

The N-array of input primal 2-forms g.

Definition at line 160 of file ATSolver2D.h.

## ◆ l_1_over_4e

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
 PrimalForm0 DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::l_1_over_4e
protected

The primal 0-form lambda/(4epsilon) (stored for performance)

Definition at line 171 of file ATSolver2D.h.

## ◆ l_1_over_4e_Id0

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
 PrimalIdentity0 DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::l_1_over_4e_Id0
protected

The 1/(4epsilon)-weighted identity operator for primal 0-forms (stored for performance)

Definition at line 158 of file ATSolver2D.h.

## ◆ lambda

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
 double DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::lambda

The global coefficient lambda giving the length of discontinuities (the smaller, the more discontinuities)

Definition at line 184 of file ATSolver2D.h.

## ◆ M01

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
 PrimalDerivative0 DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::M01
protected

The primal vertex to edge average operator.

Definition at line 150 of file ATSolver2D.h.

## ◆ M12

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
 PrimalDerivative1 DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::M12
protected

The primal edge to face average operator.

Definition at line 152 of file ATSolver2D.h.

## ◆ normalize_u2

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
 bool DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::normalize_u2

Indicates whether to normalize U (unit norm) at each iteration or not.

Definition at line 189 of file ATSolver2D.h.

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
 PrimalAntiderivative2 DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::primal_AD2
protected

The antiderivative of primal 2-forms.

Definition at line 154 of file ATSolver2D.h.

## ◆ primal_D0

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
 PrimalDerivative0 DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::primal_D0
protected

the derivative operator for primal 0-forms

Definition at line 146 of file ATSolver2D.h.

## ◆ primal_D1

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
 PrimalDerivative1 DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::primal_D1
protected

the derivative operator for primal 1-forms

Definition at line 148 of file ATSolver2D.h.

## ◆ ptrCalculus

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
 CountedConstPtrOrConstPtr< Calculus > DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::ptrCalculus
protected

## ◆ smallest_epsilon_map

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
 SmallestEpsilonMap DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::smallest_epsilon_map

The map Surfel -> double telling the smallest epsilon for which the surfel was a discontinuity.

Definition at line 178 of file ATSolver2D.h.

## ◆ surfel2idx

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
 Surfel2IndexMap DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::surfel2idx

## ◆ u2

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
 std::vector DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::u2
protected

The N-array of regularized primal 2-forms u.

Definition at line 164 of file ATSolver2D.h.

## ◆ v0

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
 PrimalForm0 DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::v0
protected

The primal 0-form v representing the set of discontinuities S (more precisely $$1 - \chi_S$$)

Definition at line 167 of file ATSolver2D.h.

## ◆ verbose

template<typename TKSpace , typename TLinearAlgebra = EigenLinearAlgebraBackend>
 int DGtal::ATSolver2D< TKSpace, TLinearAlgebra >::verbose

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