Point Cloud Library (PCL)  1.8.1-dev
poisson.h
1 /*
2  * Software License Agreement (BSD License)
3  *
4  * Copyright (c) 2010, Willow Garage, Inc.
5  * All rights reserved.
6  *
7  * Redistribution and use in source and binary forms, with or without
8  * modification, are permitted provided that the following conditions
9  * are met:
10  *
11  * * Redistributions of source code must retain the above copyright
12  * notice, this list of conditions and the following disclaimer.
13  * * Redistributions in binary form must reproduce the above
14  * copyright notice, this list of conditions and the following
15  * disclaimer in the documentation and/or other materials provided
16  * with the distribution.
17  * * Neither the name of Willow Garage, Inc. nor the names of its
18  * contributors may be used to endorse or promote products derived
19  * from this software without specific prior written permission.
20  *
21  * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
22  * "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
23  * LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS
24  * FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE
25  * COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT,
26  * INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING,
27  * BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES;
28  * LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER
29  * CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
30  * LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN
31  * ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
32  * POSSIBILITY OF SUCH DAMAGE.
33  *
34  * $Id$
35  *
36  */
37 
38 #ifndef PCL_SURFACE_POISSON_H_
39 #define PCL_SURFACE_POISSON_H_
40 
41 #include <pcl/surface/reconstruction.h>
42 
43 namespace pcl
44 {
45  namespace poisson
46  {
47  class CoredVectorMeshData;
48  template <class Real> struct Point3D;
49  }
50 
51  /** \brief The Poisson surface reconstruction algorithm.
52  * \note Code adapted from Misha Kazhdan: http://www.cs.jhu.edu/~misha/Code/PoissonRecon/
53  * \note Based on the paper:
54  * * Michael Kazhdan, Matthew Bolitho, Hugues Hoppe, "Poisson surface reconstruction",
55  * SGP '06 Proceedings of the fourth Eurographics symposium on Geometry processing
56  * \author Alexandru-Eugen Ichim
57  * \ingroup surface
58  */
59  template<typename PointNT>
60  class Poisson : public SurfaceReconstruction<PointNT>
61  {
62  public:
63  typedef boost::shared_ptr<Poisson<PointNT> > Ptr;
64  typedef boost::shared_ptr<const Poisson<PointNT> > ConstPtr;
65 
68 
70 
71  typedef typename pcl::KdTree<PointNT> KdTree;
73 
74  /** \brief Constructor that sets all the parameters to working default values. */
75  Poisson ();
76 
77  /** \brief Destructor. */
78  ~Poisson ();
79 
80  /** \brief Create the surface.
81  * \param[out] output the resultant polygonal mesh
82  */
83  void
85 
86  /** \brief Create the surface.
87  * \param[out] points the vertex positions of the resulting mesh
88  * \param[out] polygons the connectivity of the resulting mesh
89  */
90  void
92  std::vector<pcl::Vertices> &polygons);
93 
94  /** \brief Set the maximum depth of the tree that will be used for surface reconstruction.
95  * \note Running at depth d corresponds to solving on a voxel grid whose resolution is no larger than
96  * 2^d x 2^d x 2^d. Note that since the reconstructor adapts the octree to the sampling density, the specified
97  * reconstruction depth is only an upper bound.
98  * \param[in] depth the depth parameter
99  */
100  inline void
101  setDepth (int depth) { depth_ = depth; }
102 
103  /** \brief Get the depth parameter */
104  inline int
105  getDepth () { return depth_; }
106 
107  inline void
108  setMinDepth (int min_depth) { min_depth_ = min_depth; }
109 
110  inline int
111  getMinDepth () { return min_depth_; }
112 
113  inline void
114  setPointWeight (float point_weight) { point_weight_ = point_weight; }
115 
116  inline float
117  getPointWeight () { return point_weight_; }
118 
119  /** \brief Set the ratio between the diameter of the cube used for reconstruction and the diameter of the
120  * samples' bounding cube.
121  * \param[in] scale the given parameter value
122  */
123  inline void
124  setScale (float scale) { scale_ = scale; }
125 
126  /** Get the ratio between the diameter of the cube used for reconstruction and the diameter of the
127  * samples' bounding cube.
128  */
129  inline float
130  getScale () { return scale_; }
131 
132  /** \brief Set the the depth at which a block Gauss-Seidel solver is used to solve the Laplacian equation
133  * \note Using this parameter helps reduce the memory overhead at the cost of a small increase in
134  * reconstruction time. (In practice, we have found that for reconstructions of depth 9 or higher a subdivide
135  * depth of 7 or 8 can greatly reduce the memory usage.)
136  * \param[in] solver_divide the given parameter value
137  */
138  inline void
139  setSolverDivide (int solver_divide) { solver_divide_ = solver_divide; }
140 
141  /** \brief Get the the depth at which a block Gauss-Seidel solver is used to solve the Laplacian equation */
142  inline int
143  getSolverDivide () { return solver_divide_; }
144 
145  /** \brief Set the depth at which a block iso-surface extractor should be used to extract the iso-surface
146  * \note Using this parameter helps reduce the memory overhead at the cost of a small increase in extraction
147  * time. (In practice, we have found that for reconstructions of depth 9 or higher a subdivide depth of 7 or 8
148  * can greatly reduce the memory usage.)
149  * \param[in] iso_divide the given parameter value
150  */
151  inline void
152  setIsoDivide (int iso_divide) { iso_divide_ = iso_divide; }
153 
154  /** \brief Get the depth at which a block iso-surface extractor should be used to extract the iso-surface */
155  inline int
156  getIsoDivide () { return iso_divide_; }
157 
158  /** \brief Set the minimum number of sample points that should fall within an octree node as the octree
159  * construction is adapted to sampling density
160  * \note For noise-free samples, small values in the range [1.0 - 5.0] can be used. For more noisy samples,
161  * larger values in the range [15.0 - 20.0] may be needed to provide a smoother, noise-reduced, reconstruction.
162  * \param[in] samples_per_node the given parameter value
163  */
164  inline void
165  setSamplesPerNode (float samples_per_node) { samples_per_node_ = samples_per_node; }
166 
167  /** \brief Get the minimum number of sample points that should fall within an octree node as the octree
168  * construction is adapted to sampling density
169  */
170  inline float
171  getSamplesPerNode () { return samples_per_node_; }
172 
173  /** \brief Set the confidence flag
174  * \note Enabling this flag tells the reconstructor to use the size of the normals as confidence information.
175  * When the flag is not enabled, all normals are normalized to have unit-length prior to reconstruction.
176  * \param[in] confidence the given flag
177  */
178  inline void
179  setConfidence (bool confidence) { confidence_ = confidence; }
180 
181  /** \brief Get the confidence flag */
182  inline bool
183  getConfidence () { return confidence_; }
184 
185  /** \brief Enabling this flag tells the reconstructor to output a polygon mesh (rather than triangulating the
186  * results of Marching Cubes).
187  * \param[in] output_polygons the given flag
188  */
189  inline void
190  setOutputPolygons (bool output_polygons) { output_polygons_ = output_polygons; }
191 
192  /** \brief Get whether the algorithm outputs a polygon mesh or a triangle mesh */
193  inline bool
194  getOutputPolygons () { return output_polygons_; }
195 
196  /** \brief Set the degree parameter
197  * \param[in] degree the given degree
198  */
199  inline void
200  setDegree (int degree) { degree_ = degree; }
201 
202  /** \brief Get the degree parameter */
203  inline int
204  getDegree () { return degree_; }
205 
206  /** \brief Set the manifold flag.
207  * \note Enabling this flag tells the reconstructor to add the polygon barycenter when triangulating polygons
208  * with more than three vertices.
209  * \param[in] manifold the given flag
210  */
211  inline void
212  setManifold (bool manifold) { manifold_ = manifold; }
213 
214  /** \brief Get the manifold flag */
215  inline bool
216  getManifold () { return manifold_; }
217 
218  protected:
219  /** \brief Class get name method. */
220  std::string
221  getClassName () const { return ("Poisson"); }
222 
223  private:
224  int depth_;
225  int min_depth_;
226  float point_weight_;
227  float scale_;
228  int solver_divide_;
229  int iso_divide_;
230  float samples_per_node_;
231  bool confidence_;
232  bool output_polygons_;
233 
234  bool no_reset_samples_;
235  bool no_clip_tree_;
236  bool manifold_;
237 
238  int refine_;
239  int kernel_depth_;
240  int degree_;
241  bool non_adaptive_weights_;
242  bool show_residual_;
243  int min_iterations_;
244  float solver_accuracy_;
245 
246  template<int Degree> void
247  execute (poisson::CoredVectorMeshData &mesh,
248  poisson::Point3D<float> &translate,
249  float &scale);
250 
251  public:
252  EIGEN_MAKE_ALIGNED_OPERATOR_NEW
253  };
254 }
255 
256 #ifdef PCL_NO_PRECOMPILE
257 #include <pcl/surface/impl/poisson.hpp>
258 #endif
259 
260 #endif // PCL_SURFACE_POISSON_H_
float getSamplesPerNode()
Get the minimum number of sample points that should fall within an octree node as the octree construc...
Definition: poisson.h:171
bool getManifold()
Get the manifold flag.
Definition: poisson.h:216
int getDegree()
Get the degree parameter.
Definition: poisson.h:204
void setScale(float scale)
Set the ratio between the diameter of the cube used for reconstruction and the diameter of the sample...
Definition: poisson.h:124
void setConfidence(bool confidence)
Set the confidence flag.
Definition: poisson.h:179
void performReconstruction(pcl::PolygonMesh &output)
Create the surface.
Definition: poisson.hpp:152
boost::shared_ptr< const Poisson< PointNT > > ConstPtr
Definition: poisson.h:64
void setSamplesPerNode(float samples_per_node)
Set the minimum number of sample points that should fall within an octree node as the octree construc...
Definition: poisson.h:165
boost::shared_ptr< PointCloud< PointT > > Ptr
Definition: point_cloud.h:428
std::string getClassName() const
Class get name method.
Definition: poisson.h:221
SurfaceReconstruction represents a base surface reconstruction class.
boost::shared_ptr< KdTree< PointT > > Ptr
Definition: kdtree.h:71
bool getConfidence()
Get the confidence flag.
Definition: poisson.h:183
pcl::PointCloud< PointNT >::Ptr PointCloudPtr
Definition: poisson.h:69
void setIsoDivide(int iso_divide)
Set the depth at which a block iso-surface extractor should be used to extract the iso-surface...
Definition: poisson.h:152
int getMinDepth()
Definition: poisson.h:111
void setMinDepth(int min_depth)
Definition: poisson.h:108
float getScale()
Get the ratio between the diameter of the cube used for reconstruction and the diameter of the sample...
Definition: poisson.h:130
The Poisson surface reconstruction algorithm.
Definition: poisson.h:60
void setDegree(int degree)
Set the degree parameter.
Definition: poisson.h:200
float getPointWeight()
Definition: poisson.h:117
void setPointWeight(float point_weight)
Definition: poisson.h:114
void setManifold(bool manifold)
Set the manifold flag.
Definition: poisson.h:212
bool getOutputPolygons()
Get whether the algorithm outputs a polygon mesh or a triangle mesh.
Definition: poisson.h:194
int getIsoDivide()
Get the depth at which a block iso-surface extractor should be used to extract the iso-surface...
Definition: poisson.h:156
Poisson()
Constructor that sets all the parameters to working default values.
Definition: poisson.hpp:64
int getDepth()
Get the depth parameter.
Definition: poisson.h:105
int getSolverDivide()
Get the the depth at which a block Gauss-Seidel solver is used to solve the Laplacian equation...
Definition: poisson.h:143
boost::shared_ptr< Poisson< PointNT > > Ptr
Definition: poisson.h:63
pcl::KdTree< PointNT >::Ptr KdTreePtr
Definition: poisson.h:72
void setSolverDivide(int solver_divide)
Set the the depth at which a block Gauss-Seidel solver is used to solve the Laplacian equation...
Definition: poisson.h:139
KdTree represents the base spatial locator class for kd-tree implementations.
Definition: kdtree.h:56
void setDepth(int depth)
Set the maximum depth of the tree that will be used for surface reconstruction.
Definition: poisson.h:101
~Poisson()
Destructor.
Definition: poisson.hpp:89
void setOutputPolygons(bool output_polygons)
Enabling this flag tells the reconstructor to output a polygon mesh (rather than triangulating the re...
Definition: poisson.h:190
pcl::KdTree< PointNT > KdTree
Definition: poisson.h:71