DGtal  1.5.beta
DECHelpers.h
1 
17 #pragma once
18 
31 #if defined(DECHelpers_RECURSES)
32 #error Recursive header files inclusion detected in DECHelpers.h
33 #else // defined(DECHelpers_RECURSES)
35 #define DECHelpers_RECURSES
36 
37 #if !defined DECHelpers_h
39 #define DECHelpers_h
40 
42 // Inclusions
43 #include <iostream>
44 #include "DGtal/base/Common.h"
45 #include "DGtal/dec/Duality.h"
46 #include "DGtal/dec/DiscreteExteriorCalculus.h"
48 
49 namespace DGtal
50 {
52  namespace dec_helper
53  {
54 
58  template <typename Calculus, DGtal::Dimension dim, DGtal::Duality duality>
61  {
63  typedef typename Calculus::LinearAlgebraBackend::Triplet Triplet;
64  typedef typename Calculus::Index Index;
65  typedef std::vector<Triplet> Triplets;
66 
67  Triplets triplets;
68  for (Index index=0; index<kform.length(); index++)
69  triplets.emplace_back( index, index, kform.myContainer(index) );
70  Operator ope(kform.myCalculus);
71  ope.myContainer.setFromTriplets(triplets.begin(), triplets.end());
72 
73  return ope;
74  }
75 
80  template <typename Calculus>
82  averageOperator01(const Calculus& calculus)
83  {
84  auto M01 = calculus.template derivative<0,DGtal::PRIMAL>();
85  M01.myContainer = 0.5 * M01.myContainer.cwiseAbs();
86  return M01;
87  }
88 
93  template <typename Calculus>
95  averageOperator01(const Calculus& calculus)
96  {
97  auto M12 = calculus.template derivative<1,DGtal::PRIMAL>();
98  M12.myContainer = 0.25 * M12.myContainer.cwiseAbs();
99  return M12;
100  }
101 
106  template <typename Calculus>
108  averageOperator20(const Calculus& calculus)
109  {
110  BOOST_STATIC_ASSERT( Calculus::dimensionEmbedded == 2 );
111  BOOST_STATIC_ASSERT( Calculus::dimensionAmbient == 3 );
112  using DGtal::PRIMAL;
113  using DGtal::DUAL;
114 
116  typedef typename Calculus::LinearAlgebraBackend::Triplet Triplet;
117  typedef typename Calculus::KSpace KSpace;
118  typedef typename Calculus::Index Index;
119  typedef typename Calculus::Cell Cell;
120  typedef typename Calculus::Scalar Scalar;
121  typedef typename Calculus::Point Point;
123  Operator;
124 
125  const KSpace& kspace = calculus.myKSpace;
126 
127  const std::vector<Point> deltas = {
128  Point(0,1,1), Point(0,-1,1), Point(0,-1,-1), Point(0,1,-1),
129  Point(1,0,1), Point(-1,0,1), Point(-1,0,-1), Point(1,0,-1),
130  Point(1,1,0), Point(-1,1,0), Point(-1,-1,0), Point(1,-1,0)
131  };
132 
133  std::vector<Triplet> triplets;
134  for (Index index_point=0; index_point<calculus.kFormLength(0,PRIMAL); index_point++)
135  {
136  const Cell point = kspace.unsigns(calculus.getSCell(0, PRIMAL, index_point));
137  ASSERT( kspace.uDim(point) == 0 );
138 
139  std::vector<Index> indexes_surfel;
140  for (const Point delta : deltas)
141  {
142  const Cell surfel = kspace.uCell((kspace.uKCoords(point)+delta));
143  ASSERT( kspace.uDim(surfel) == 2 );
144  if (calculus.containsCell(surfel))
145  indexes_surfel.push_back(calculus.getCellIndex(surfel));
146  }
147  ASSERT( indexes_surfel.size() > 2 );
148 
149  const double weight = 1/static_cast<Scalar>(indexes_surfel.size());
150  for (const Index index_surfel : indexes_surfel)
151  triplets.emplace_back( index_point, index_surfel, weight );
152  }
153 
154  SparseMatrix matrix( calculus.kFormLength(0, DGtal::PRIMAL),
155  calculus.kFormLength(2, DGtal::PRIMAL) );
156  matrix.setFromTriplets(triplets.begin(), triplets.end());
157 
158  return Operator(calculus, matrix);
159  }
160 
161  } // namespace dec_helper
162 
163 } // namespace DGtal
164 
165 
166 // //
168 
169 #endif // !defined DECHelpers_h
170 
171 #undef DECHelpers_RECURSES
172 #endif // else defined(DECHelpers_RECURSES)
Aim: KForm represents discrete kforms in the dec package.
Definition: KForm.h:66
const Calculus * myCalculus
Definition: KForm.h:136
Container myContainer
Definition: KForm.h:131
Index length() const
Aim: This class is a model of CCellularGridSpaceND. It represents the cubical grid as a cell complex,...
Cell unsigns(const SCell &p) const
Creates an unsigned cell from a signed one.
Cell uCell(const PreCell &c) const
From an unsigned cell, returns an unsigned cell lying into this Khalismky space.
Dimension uDim(const Cell &p) const
Return the dimension of the cell [p].
const Point & uKCoords(const Cell &c) const
Return its Khalimsky coordinates.
Aim: LinearOperator represents discrete linear operator between discrete kforms in the DEC package.
PolyCalculus * calculus
SMesh::Index Index
DGtal::LinearOperator< Calculus, dim, duality, dim, duality > diagonal(const DGtal::KForm< Calculus, dim, duality > &kform)
Definition: DECHelpers.h:60
DGtal::LinearOperator< Calculus, 0, DGtal::PRIMAL, 1, DGtal::PRIMAL > averageOperator01(const Calculus &calculus)
Definition: DECHelpers.h:82
DGtal::LinearOperator< Calculus, 2, DGtal::PRIMAL, 0, DGtal::PRIMAL > averageOperator20(const Calculus &calculus)
Definition: DECHelpers.h:108
DGtal is the top-level namespace which contains all DGtal functions and types.
@ PRIMAL
Definition: Duality.h:61
@ DUAL
Definition: Duality.h:62
MyPointD Point
Definition: testClone2.cpp:383
KSpace::Cell Cell
EigenLinearAlgebraBackend::SparseMatrix SparseMatrix