DGtal 1.3.0
Loading...
Searching...
No Matches
VoronoiCovarianceMeasureOnDigitalSurface.ih
1/**
2 * This program is free software: you can redistribute it and/or modify
3 * it under the terms of the GNU Lesser General Public License as
4 * published by the Free Software Foundation, either version 3 of the
5 * License, or (at your option) any later version.
6 *
7 * This program is distributed in the hope that it will be useful,
8 * but WITHOUT ANY WARRANTY; without even the implied warranty of
9 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
10 * GNU General Public License for more details.
11 *
12 * You should have received a copy of the GNU General Public License
13 * along with this program. If not, see <http://www.gnu.org/licenses/>.
14 *
15 **/
16
17/**
18 * @file VoronoiCovarianceMeasureOnDigitalSurface.ih
19 * @author Jacques-Olivier Lachaud (\c jacques-olivier.lachaud@univ-savoie.fr )
20 * Laboratory of Mathematics (CNRS, UMR 5127), University of Savoie, France
21 *
22 * @date 2014/02/13
23 *
24 * Implementation of inline methods defined in VoronoiCovarianceMeasureOnDigitalSurface.h
25 *
26 * This file is part of the DGtal library.
27 */
28
29
30//////////////////////////////////////////////////////////////////////////////
31#include <cstdlib>
32#include "DGtal/topology/CanonicSCellEmbedder.h"
33#include "DGtal/math/ScalarFunctors.h"
34#include "DGtal/geometry/surfaces/estimation/LocalEstimatorFromSurfelFunctorAdapter.h"
35#include "DGtal/geometry/surfaces/estimation/estimationFunctors/ElementaryConvolutionNormalVectorEstimator.h"
36#include "DGtal/geometry/volumes/distance/LpMetric.h"
37
38//////////////////////////////////////////////////////////////////////////////
39
40///////////////////////////////////////////////////////////////////////////////
41// IMPLEMENTATION of inline methods.
42///////////////////////////////////////////////////////////////////////////////
43
44///////////////////////////////////////////////////////////////////////////////
45// ----------------------- Standard services ------------------------------
46
47//-----------------------------------------------------------------------------
48template <typename TDigitalSurfaceContainer, typename TSeparableMetric, typename TKernelFunction>
49inline
50DGtal::VoronoiCovarianceMeasureOnDigitalSurface<TDigitalSurfaceContainer, TSeparableMetric, TKernelFunction>::
51~VoronoiCovarianceMeasureOnDigitalSurface()
52{
53}
54//-----------------------------------------------------------------------------
55template <typename TDigitalSurfaceContainer, typename TSeparableMetric, typename TKernelFunction>
56inline
57DGtal::VoronoiCovarianceMeasureOnDigitalSurface<TDigitalSurfaceContainer, TSeparableMetric, TKernelFunction>::
58VoronoiCovarianceMeasureOnDigitalSurface( ConstAlias< Surface > _surface,
59 Surfel2PointEmbedding _surfelEmbedding,
60 Scalar _R, Scalar _r,
61 KernelFunction chi_r,
62 Scalar t, Metric aMetric, bool verbose )
63 : mySurface( _surface ), mySurfelEmbedding( _surfelEmbedding ), myChi( chi_r ),
64 myVCM( _R, _r, aMetric, verbose ), myRadiusTrivial( t )
65{
66 if ( verbose ) trace.beginBlock( "Computing VCM on digital surface." );
67 const KSpace & ks = mySurface->container().space();
68 std::vector<Point> vectPoints;
69
70 // Get points.
71 if ( verbose ) trace.beginBlock( "Getting points." );
72 std::set<Point> pointSet;
73 for ( ConstIterator it = mySurface->begin(), itE = mySurface->end(); it != itE; ++it )
74 getPoints( std::inserter( pointSet, pointSet.begin() ), *it );
75 vectPoints.resize( pointSet.size() );
76 std::copy( pointSet.begin(), pointSet.end(), vectPoints.begin() );
77 pointSet.clear();
78 if ( verbose ) trace.endBlock();
79
80 // Compute Voronoi Covariance Matrix for all points.
81 myVCM.init( vectPoints.begin(), vectPoints.end() );
82
83 // Compute VCM( chi_r ) for each point.
84 if ( verbose ) trace.beginBlock ( "Integrating VCM( chi_r(p) ) for each point." );
85 int i = 0;
86 // HatPointFunction< Point, Scalar > chi_r( 1.0, r );
87 for ( typename std::vector<Point>::const_iterator it = vectPoints.begin(), itE = vectPoints.end();
88 it != itE; ++it )
89 {
90 if ( verbose ) trace.progressBar( ++i, vectPoints.size() );
91 Point p = *it;
92 MatrixNN measure = myVCM.measure( myChi, p );
93 // On diagonalise le résultat.
94 EigenStructure & evcm = myPt2EigenStructure[ p ];
95 LinearAlgebraTool::getEigenDecomposition( measure, evcm.vectors, evcm.values );
96 }
97 myVCM.clean(); // free some memory.
98 if ( verbose ) trace.endBlock();
99
100 if ( verbose ) trace.beginBlock ( "Computing average orientation for each surfel." );
101 typedef functors::HatFunction<Scalar> Functor;
102 Functor fct( 1.0, myRadiusTrivial );
103 LpMetric<Space> l2(2.0); //L2 metric in R^3 for surface propagation.
104 typedef functors::ElementaryConvolutionNormalVectorEstimator< Surfel, CanonicSCellEmbedder<KSpace> >
105 SurfelFunctor;
106 typedef LocalEstimatorFromSurfelFunctorAdapter< DigitalSurfaceContainer, LpMetric<Space>, SurfelFunctor, Functor>
107 NormalEstimator;
108
109 CanonicSCellEmbedder<KSpace> canonic_embedder( ks );
110 SurfelFunctor surfelFct( canonic_embedder, 1.0 );
111 NormalEstimator estimator;
112 estimator.attach( *mySurface);
113 estimator.setParams( l2, surfelFct, fct , myRadiusTrivial);
114 estimator.init( 1.0, mySurface->begin(), mySurface->end());
115 i = 0;
116 std::vector<Point> pts;
117 int surf_size = mySurface->size();
118 for ( ConstIterator it = mySurface->begin(), itE = mySurface->end(); it != itE; ++it )
119 {
120 if ( verbose ) trace.progressBar(++i, surf_size );
121 Surfel s = *it;
122 Normals & normals = mySurfel2Normals[ s ];
123 // get rough estimation of normal
124 normals.trivialNormal = estimator.eval( it );
125 // get points associated with surfel s
126 getPoints( std::back_inserter( pts ), s );
127 for ( typename std::vector<Point>::const_iterator itPts = pts.begin(), itPtsE = pts.end();
128 itPts != itPtsE; ++itPts )
129 {
130 Point p = *itPts;
131 const EigenStructure& evcm = myPt2EigenStructure[ p ];
132 VectorN n = evcm.vectors.column( Space::dimension-1 );
133 if ( n.dot( normals.trivialNormal ) < 0 ) normals.vcmNormal -= n;
134 else normals.vcmNormal += n;
135 }
136 if ( pts.size() > 1 ) normals.vcmNormal /= pts.size();
137 pts.clear();
138 }
139 if ( verbose ) trace.endBlock();
140
141 if ( verbose ) trace.endBlock();
142}
143
144//-----------------------------------------------------------------------------
145template <typename TDigitalSurfaceContainer, typename TSeparableMetric, typename TKernelFunction>
146inline
147DGtal::CountedConstPtrOrConstPtr< typename DGtal::VoronoiCovarianceMeasureOnDigitalSurface<TDigitalSurfaceContainer, TSeparableMetric, TKernelFunction>::Surface >
148DGtal::VoronoiCovarianceMeasureOnDigitalSurface<TDigitalSurfaceContainer, TSeparableMetric, TKernelFunction>::
149surface() const
150{
151 return mySurface;
152}
153//-----------------------------------------------------------------------------
154template <typename TDigitalSurfaceContainer, typename TSeparableMetric, typename TKernelFunction>
155inline
156DGtal::Surfel2PointEmbedding
157DGtal::VoronoiCovarianceMeasureOnDigitalSurface<TDigitalSurfaceContainer, TSeparableMetric, TKernelFunction>::
158surfelEmbedding() const
159{
160 return mySurfelEmbedding;
161}
162//-----------------------------------------------------------------------------
163template <typename TDigitalSurfaceContainer, typename TSeparableMetric, typename TKernelFunction>
164inline
165typename DGtal::VoronoiCovarianceMeasureOnDigitalSurface<TDigitalSurfaceContainer, TSeparableMetric, TKernelFunction>::Scalar
166DGtal::VoronoiCovarianceMeasureOnDigitalSurface<TDigitalSurfaceContainer, TSeparableMetric, TKernelFunction>::
167R() const
168{
169 return myVCM.R();
170}
171//-----------------------------------------------------------------------------
172template <typename TDigitalSurfaceContainer, typename TSeparableMetric, typename TKernelFunction>
173inline
174typename DGtal::VoronoiCovarianceMeasureOnDigitalSurface<TDigitalSurfaceContainer, TSeparableMetric, TKernelFunction>::Scalar
175DGtal::VoronoiCovarianceMeasureOnDigitalSurface<TDigitalSurfaceContainer, TSeparableMetric, TKernelFunction>::
176r() const
177{
178 return myVCM.r();
179}
180//-----------------------------------------------------------------------------
181template <typename TDigitalSurfaceContainer, typename TSeparableMetric, typename TKernelFunction>
182inline
183typename DGtal::VoronoiCovarianceMeasureOnDigitalSurface<TDigitalSurfaceContainer, TSeparableMetric, TKernelFunction>::Scalar
184DGtal::VoronoiCovarianceMeasureOnDigitalSurface<TDigitalSurfaceContainer, TSeparableMetric, TKernelFunction>::
185radiusTrivial() const
186{
187 return myRadiusTrivial;
188}
189//-----------------------------------------------------------------------------
190template <typename TDigitalSurfaceContainer, typename TSeparableMetric, typename TKernelFunction>
191inline
192const typename DGtal::VoronoiCovarianceMeasureOnDigitalSurface<TDigitalSurfaceContainer, TSeparableMetric, TKernelFunction>::Surfel2Normals&
193DGtal::VoronoiCovarianceMeasureOnDigitalSurface<TDigitalSurfaceContainer, TSeparableMetric, TKernelFunction>::
194mapSurfel2Normals() const
195{
196 return mySurfel2Normals;
197}
198//-----------------------------------------------------------------------------
199template <typename TDigitalSurfaceContainer, typename TSeparableMetric, typename TKernelFunction>
200inline
201const typename DGtal::VoronoiCovarianceMeasureOnDigitalSurface<TDigitalSurfaceContainer, TSeparableMetric, TKernelFunction>::Point2EigenStructure&
202DGtal::VoronoiCovarianceMeasureOnDigitalSurface<TDigitalSurfaceContainer, TSeparableMetric, TKernelFunction>::
203mapPoint2ChiVCM() const
204{
205 return myPt2EigenStructure;
206}
207
208//-----------------------------------------------------------------------------
209template <typename TDigitalSurfaceContainer, typename TSeparableMetric, typename TKernelFunction>
210inline
211bool
212DGtal::VoronoiCovarianceMeasureOnDigitalSurface<TDigitalSurfaceContainer, TSeparableMetric, TKernelFunction>::
213getChiVCMEigenvalues( VectorN& values, Surfel s ) const
214{
215 std::vector<Point> pts;
216 getPoints( std::back_inserter( pts ), s );
217 bool ok = true;
218 int i = 0;
219 values = VectorN(); // Setting values to 0 before averaging.
220 for ( typename std::vector<Point>::const_iterator itPts = pts.begin(), itPtsE = pts.end();
221 itPts != itPtsE; ++itPts, ++i )
222 {
223 Point p = *itPts;
224 typename Point2EigenStructure::const_iterator itEigen = myPt2EigenStructure.find( p );
225 if ( itEigen == myPt2EigenStructure.end() )
226 {
227 ok = false;
228 break;
229 }
230 const EigenStructure& evcm = itEigen->second;
231 values += evcm.values;
232 }
233 if ( i > 1 ) values /= i;
234 return ok;
235}
236
237//-----------------------------------------------------------------------------
238template <typename TDigitalSurfaceContainer, typename TSeparableMetric, typename TKernelFunction>
239inline
240bool
241DGtal::VoronoiCovarianceMeasureOnDigitalSurface<TDigitalSurfaceContainer, TSeparableMetric, TKernelFunction>::
242getChiVCMEigenStructure( VectorN& values, MatrixNN& vectors, Surfel s ) const
243{
244 std::vector<Point> pts;
245 getPoints( std::back_inserter( pts ), s );
246 bool ok = true;
247 int i = 0;
248 values = VectorN(); // Setting values to 0 before averaging.
249 vectors = MatrixNN(); // Setting values to 0 before averaging.
250 for ( typename std::vector<Point>::const_iterator itPts = pts.begin(), itPtsE = pts.end();
251 itPts != itPtsE; ++itPts, ++i )
252 {
253 Point p = *itPts;
254 typename Point2EigenStructure::const_iterator itEigen = myPt2EigenStructure.find( p );
255 if ( itEigen == myPt2EigenStructure.end() )
256 {
257 ok = false;
258 break;
259 }
260 const EigenStructure& evcm = itEigen->second;
261 values += evcm.values;
262 vectors += evcm.vectors;
263 }
264 if ( i > 1 ) {
265 values /= i;
266 vectors /= i;
267 }
268 return ok;
269}
270
271//-----------------------------------------------------------------------------
272template <typename TDigitalSurfaceContainer, typename TSeparableMetric, typename TKernelFunction>
273template <typename PointOutputIterator>
274inline
275PointOutputIterator
276DGtal::VoronoiCovarianceMeasureOnDigitalSurface<TDigitalSurfaceContainer, TSeparableMetric, TKernelFunction>::
277getPoints( PointOutputIterator outIt, Surfel s ) const
278{
279 BOOST_CONCEPT_ASSERT(( boost::OutputIterator< PointOutputIterator, Point > ));
280 const KSpace & ks = mySurface->container().space();
281 Dimension k = ks.sOrthDir( s );
282 switch ( mySurfelEmbedding ) {
283 case Pointels:
284 {
285 typename KSpace::Cells faces = ks.uFaces( ks.unsigns( s ) );
286 for ( typename KSpace::Cells::const_iterator it = faces.begin(), itE = faces.end();
287 it != itE; ++it )
288 {
289 if ( ks.uDim( *it ) == 0 ) // get only pointels (cell of dim 0)
290 *outIt++ = ks.uCoords( *it );
291 }
292 // Dimension i = (k+1)%3;
293 // Dimension j = (i+1)%3;
294 // SCell l1 = ks.sIncident( s, i, true );
295 // SCell l2 = ks.sIncident( s, i, false );
296 // *outIt++ = ks.sCoords( ks.sIncident( l1, j, true ) );
297 // *outIt++ = ks.sCoords( ks.sIncident( l1, j, false ) );
298 // *outIt++ = ks.sCoords( ks.sIncident( l2, j, true ) );
299 // *outIt++ = ks.sCoords( ks.sIncident( l2, j, false ) );
300 }
301 break;
302 case InnerSpel:
303 *outIt++ = ks.sCoords( ks.sDirectIncident( s, k ) );
304 break;
305 case OuterSpel:
306 *outIt++ = ks.sCoords( ks.sIndirectIncident( s, k ) );
307 break;
308 }
309 return outIt;
310}
311
312///////////////////////////////////////////////////////////////////////////////
313// Interface - public :
314
315/**
316 * Writes/Displays the object on an output stream.
317 * @param out the output stream where the object is written.
318 */
319template <typename TDigitalSurfaceContainer, typename TSeparableMetric, typename TKernelFunction>
320inline
321void
322DGtal::VoronoiCovarianceMeasureOnDigitalSurface<TDigitalSurfaceContainer, TSeparableMetric, TKernelFunction>::
323selfDisplay ( std::ostream & out ) const
324{
325 out << "[VoronoiCovarianceMeasureOnDigitalSurface"
326 << " #pts=" << myPt2EigenStructure.size()
327 << " #surf=" << mySurfel2Normals.size()
328 << "]";
329}
330
331/**
332 * Checks the validity/consistency of the object.
333 * @return 'true' if the object is valid, 'false' otherwise.
334 */
335template <typename TDigitalSurfaceContainer, typename TSeparableMetric, typename TKernelFunction>
336inline
337bool
338DGtal::VoronoiCovarianceMeasureOnDigitalSurface<TDigitalSurfaceContainer, TSeparableMetric, TKernelFunction>::
339isValid() const
340{
341 return true;
342}
343
344
345
346///////////////////////////////////////////////////////////////////////////////
347// Implementation of inline functions //
348
349template <typename TDigitalSurfaceContainer, typename TSeparableMetric, typename TKernelFunction>
350inline
351std::ostream&
352DGtal::operator<< ( std::ostream & out,
353 const VoronoiCovarianceMeasureOnDigitalSurface<TDigitalSurfaceContainer, TSeparableMetric, TKernelFunction> & object )
354{
355 object.selfDisplay( out );
356 return out;
357}
358
359// //
360///////////////////////////////////////////////////////////////////////////////
361
362
Generated on Sun Nov 27 2022 15:15:39 for DGtal by doxygen 1.9.5