DGtal 1.3.0
Loading...
Searching...
No Matches
exampleIntegralInvariantCurvature3D.cpp
Go to the documentation of this file.
1
32#include <iostream>
33#include "DGtal/base/Common.h"
34
35// Shape construction
36#include "DGtal/io/readers/VolReader.h"
37#include "DGtal/images/ImageSelector.h"
38#include "DGtal/images/imagesSetsUtils/SetFromImage.h"
39#include "DGtal/images/SimpleThresholdForegroundPredicate.h"
40#include "DGtal/topology/helpers/Surfaces.h"
41#include "DGtal/topology/LightImplicitDigitalSurface.h"
42#include "DGtal/images/ImageHelper.h"
43#include "DGtal/topology/DigitalSurface.h"
44#include "DGtal/graph/DepthFirstVisitor.h"
45#include "DGtal/graph/GraphVisitorRange.h"
46
48#include "DGtal/geometry/surfaces/estimation/IIGeometricFunctors.h"
49#include "DGtal/geometry/surfaces/estimation/IntegralInvariantVolumeEstimator.h"
50
51// Drawing
52#include "DGtal/io/viewers/Viewer3D.h"
53#include "DGtal/io/colormaps/GradientColorMap.h"
54
56
57using namespace DGtal;
58
60
61int main( int argc, char** argv )
62{
63 if ( argc != 4 )
64 {
65 trace.error() << "Usage: " << argv[0]
66 << " <fileName.vol> <threshold> <radius>" << std::endl;
67 trace.error() << "Example : "<< argv[0] << " Al.150.vol 0 7" << std::endl;
68 return 0;
69 }
70
71 trace.beginBlock ( "Example IntegralInvariantCurvature3D" );
72 trace.info() << "Args:";
73 for ( int i = 0; i < argc; ++i )
74 trace.info() << " " << argv[ i ];
75 trace.info() << std::endl;
76
77 double h = 1.0;
78 unsigned int threshold = std::atoi( argv[ 2 ] );
79
83 typedef LightImplicitDigitalSurface< Z3i::KSpace, ImagePredicate > MyLightImplicitDigitalSurface;
85
86 std::string filename = argv[1];
87 Image image = VolReader<Image>::importVol( filename );
88 ImagePredicate predicate = ImagePredicate( image, threshold );
89
90 Z3i::Domain domain = image.domain();
91
92 Z3i::KSpace KSpaceShape;
93
94 bool space_ok = KSpaceShape.init( domain.lowerBound(), domain.upperBound(), true );
95 if (!space_ok)
96 {
97 trace.error() << "Error in the Khalimsky space construction."<<std::endl;
98 return 2;
99 }
100
102 Z3i::KSpace::Surfel bel = Surfaces< Z3i::KSpace >::findABel( KSpaceShape, predicate, 100000 );
103 MyLightImplicitDigitalSurface LightImplDigSurf( KSpaceShape, predicate, SAdj, bel );
104 MyDigitalSurface digSurf( LightImplDigSurf );
105
107 typedef GraphVisitorRange< Visitor > VisitorRange;
108 typedef VisitorRange::ConstIterator SurfelConstIterator;
109
110 VisitorRange range( new Visitor( digSurf, *digSurf.begin() ) );
111 SurfelConstIterator abegin = range.begin();
112 SurfelConstIterator aend = range.end();
113
116 double radius = std::atof(argv[3]);
117
118 typedef functors::IIMeanCurvature3DFunctor<Z3i::Space> MyIICurvatureFunctor;
120
121 // For computing Gaussian curvature instead, for example, change the two typedef above by :
122 // typedef functors::IIGaussianCurvature3DFunctor<Z3i::Space> MyIICurvatureFunctor;
123 // typedef IntegralInvariantCovarianceEstimator< Z3i::KSpace, ImagePredicate, MyIICurvatureFunctor > MyIICurvatureEstimator;
124 // and it's done. The following part is exactly the same.
125
126 typedef MyIICurvatureFunctor::Value Value;
127
128 MyIICurvatureFunctor curvatureFunctor; // Functor used to convert volume -> curvature
129 curvatureFunctor.init( h, radius ); // Initialisation for a grid step and a given Euclidean radius of convolution kernel
130
131 MyIICurvatureEstimator curvatureEstimator( curvatureFunctor );
132 curvatureEstimator.attach( KSpaceShape, predicate ); // Setting a KSpace and a predicate on the object to evaluate
133 curvatureEstimator.setParams( radius / h ); // Setting the digital radius of the convolution kernel
134 curvatureEstimator.init( h, abegin, aend ); // Initialisation for a given h, and a range of surfels
135
136 std::vector< Value > results;
137 std::back_insert_iterator< std::vector< Value > > resultsIt( results ); // output iterator for results of Integral Invariant curvature computation
138 curvatureEstimator.eval( abegin, aend, resultsIt ); // Computation
140
142 Value min = std::numeric_limits < Value >::max();
143 Value max = std::numeric_limits < Value >::min();
144 for ( unsigned int i = 0; i < results.size(); ++i )
145 {
146 if ( results[ i ] < min )
147 {
148 min = results[ i ];
149 }
150 else if ( results[ i ] > max )
151 {
152 max = results[ i ];
153 }
154 }
155
156 QApplication application( argc, argv );
158 Viewer viewer( KSpaceShape );
159 viewer.setWindowTitle("example Integral Invariant 3D");
160 viewer.show();
161
162 typedef GradientColorMap< Value > Gradient;
163 Gradient cmap_grad( min, max );
164 cmap_grad.addColor( Color( 50, 50, 255 ) );
165 cmap_grad.addColor( Color( 255, 0, 0 ) );
166 cmap_grad.addColor( Color( 255, 255, 10 ) );
167
168 VisitorRange range2( new Visitor( digSurf, *digSurf.begin() ) );
169 abegin = range2.begin();
170
171 Z3i::KSpace::Cell dummy_cell;
172 viewer << SetMode3D( dummy_cell.className(), "Basic" );
173
174 for ( unsigned int i = 0; i < results.size(); ++i )
175 {
176 viewer << CustomColors3D( Color::Black, cmap_grad( results[ i ] ))
177 << KSpaceShape.unsigns( *abegin );
178 ++abegin;
179 }
180
181 viewer << Viewer3D<>::updateDisplay;
182
183 trace.endBlock();
184 return application.exec();
185}
186// //
Structure representing an RGB triple with alpha component.
Definition: Color.h:68
static const Color Black
Definition: Color.h:413
Aim: This class is useful to perform a depth-first exploration of a graph given a starting point or s...
Aim: Represents a set of n-1-cells in a nD space, together with adjacency relation between these cell...
ConstIterator begin() const
Aim: This class template may be used to (linearly) convert scalar values in a given range into a colo...
Aim: Transforms a graph visitor into a single pass input range.
const Point & lowerBound() const
const Point & upperBound() const
Aim: implements association bewteen points lying in a digital domain and values.
Definition: Image.h:70
Aim: This class implement an Integral Invariant estimator which computes for each surfel the volume o...
Aim: This class is a model of CCellularGridSpaceND. It represents the cubical grid as a cell complex,...
bool init(const Point &lower, const Point &upper, bool isClosed)
Specifies the upper and lower bounds for the maximal cells in this space.
Cell unsigns(const SCell &p) const
Creates an unsigned cell from a signed one.
Aim: A model of CDigitalSurfaceContainer which defines the digital surface as the boundary of an impl...
static SCell findABel(const KSpace &K, const PointPredicate &pp, unsigned int nbtries=1000)
Aim: Represent adjacencies between surfel elements, telling if it follows an interior to exterior ord...
void beginBlock(const std::string &keyword="")
std::ostream & error()
std::ostream & info()
double endBlock()
virtual void show()
Overload QWidget method in order to add a call to updateList() method (to ensure that the lists are w...
Aim: Define a simple Foreground predicate thresholding image values given a single thresold....
DigitalSurface< MyDigitalSurfaceContainer > MyDigitalSurface
BreadthFirstVisitor< MyDigitalSurface > Visitor
DGtal is the top-level namespace which contains all DGtal functions and types.
Trace trace
Definition: Common.h:154
std::string className() const
Return the style name used for drawing this object.
Modifier class in a Display3D stream. Useful to choose your own mode for a given class....
Represents a signed cell in a cellular grid space by its Khalimsky coordinates and a boolean value.
static ImageContainer importVol(const std::string &filename, const Functor &aFunctor=Functor())
Aim: A functor Real -> Real that returns the 3d mean curvature by transforming the given volume....
int max(int a, int b)
int main()
Definition: testBits.cpp:56
Domain domain