c++ - 给定点,计算法线并在 CGAL 中执行多边形表面重建
问题描述
我有 .xyz 文件,我想在 CGAL 中使用“多边形表面重建”重新创建表面。重建表面需要法线,所以我想首先在 CGAL 中计算这些,然后使用点及其法线作为多边形表面重建的输入。
但是,两者的输出和输入类型不同:
正常估计的输出是向量对的列表,参见https://doc.cgal.org/latest/Point_set_processing_3/index.html#Point_set_processing_3NormalOrientation : 12.1.1 示例;而多边形表面重建的输入是一个point_vector((https://doc.cgal.org/latest/Polygonal_surface_reconstruction/index.html)。
我试图创建一个转换数据格式的“convert_function”,但我不是 C++ 程序员,因此我很难转换。如果有人能帮我解决这个问题,那就太好了,因为我已经坚持了好几个小时了。我的目标是使这个过程尽可能地节省时间,因此非常感谢有关如何解决这个问题的更好的想法。
数据类型:
->正常计算:
typedef std::pair<Point, Vector> PointVectorPair;
->表面重建:
typedef boost::tuple<Point, Vector, int> PNI;
typedef std::vector<PNI> Point_vector;
typedef CGAL::Nth_of_tuple_property_map<0, PNI> Point_map;
typedef CGAL::Nth_of_tuple_property_map<1, PNI> Normal_map;
typedef CGAL::Nth_of_tuple_property_map<2, PNI> Plane_index_map;
我的反转函数不起作用,因为语法不正确:
std::vector<PNI> convert_function(const std::list<PointVectorPair>& list)
{
std::vector<PNI> out;//a temporary object to store the output
for (std::size_t i = 0; i < list.size(); ++i) {
{
out[i].get<0>() = list[i].get<0>();
out[i].get<1>() = list[i].get<1>();
}
return out;
}
再进一步,我使用以下内容:
std::vector<PNI> points2; // store points
points2 = convert_function(points);
现在完整的代码如下:
#include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
#include <CGAL/IO/read_xyz_points.h>
#include <CGAL/IO/Writer_OFF.h>
#include <CGAL/property_map.h>
#include <CGAL/Surface_mesh.h>
#include <CGAL/Shape_detection/Efficient_RANSAC.h>
#include <CGAL/Polygonal_surface_reconstruction.h>
//extras formal calculation
#include <CGAL/compute_average_spacing.h>
#include <CGAL/pca_estimate_normals.h>
#include <CGAL/mst_orient_normals.h>
#include <utility> // defines std::pair
#include <list>
///////////// TRY TO COMBINE NORMAL CALCULATION + RANSAC + POLYFIT
#ifdef CGAL_USE_SCIP // defined (or not) by CMake scripts, do not define by hand
#include <CGAL/SCIP_mixed_integer_program_traits.h>
typedef CGAL::SCIP_mixed_integer_program_traits<double> MIP_Solver;
#elif defined(CGAL_USE_GLPK) // defined (or not) by CMake scripts, do not define by hand
#include <CGAL/GLPK_mixed_integer_program_traits.h>
typedef CGAL::GLPK_mixed_integer_program_traits<double> MIP_Solver;
#endif
#if defined(CGAL_USE_GLPK) || defined(CGAL_USE_SCIP)
#include <CGAL/Timer.h>
#include <fstream>
typedef CGAL::Exact_predicates_inexact_constructions_kernel Kernel;
typedef Kernel::Point_3 Point;
typedef Kernel::Vector_3 Vector;
// Point with normal, and plane index
//
typedef boost::tuple<Point, Vector, int> PNI;
typedef std::vector<PNI> Point_vector;
//bij normals calc: Point with normal vector stored in a std::pair.
//typedef std::pair<Point, Vector> PointVectorPair;
typedef CGAL::Nth_of_tuple_property_map<0, PNI> Point_map;
typedef CGAL::Nth_of_tuple_property_map<1, PNI> Normal_map;
typedef CGAL::Nth_of_tuple_property_map<2, PNI> Plane_index_map;
typedef std::pair<Point, Vector> PointVectorPair;
typedef CGAL::Shape_detection::Efficient_RANSAC_traits<Kernel, Point_vector, Point_map, Normal_map> Traits;
typedef CGAL::Shape_detection::Efficient_RANSAC<Traits> Efficient_ransac;
typedef CGAL::Shape_detection::Plane<Traits> Plane;
typedef CGAL::Shape_detection::Point_to_shape_index_map<Traits> Point_to_shape_index_map;
typedef CGAL::Polygonal_surface_reconstruction<Kernel> Polygonal_surface_reconstruction;
typedef CGAL::Surface_mesh<Point> Surface_mesh;
// Concurrency
typedef CGAL::Parallel_if_available_tag Concurrency_tag;
std::vector<PNI> convert_function(const std::list<PointVectorPair>& list)
{
std::vector<PNI> out;//a temporary object to store the output
for (std::size_t i = 0; i < list.size(); ++i) {
{
out[i].get<0>() = list[i].get<0>();
out[i].get<1>() = list[i].get<1>();
}
return out;
}
/*
* This example first extracts planes from the input point cloud
* (using RANSAC with default parameters) and then reconstructs
* the surface model from the planes.
*/
int main(int argc, char* argv[])
{
const char* fname = (argc > 1) ? argv[1] : "data/47038.xyz";
// Reads a .xyz point set file in points[].
std::list<PointVectorPair> points;
std::ifstream stream(fname);
if (!stream ||
!CGAL::read_xyz_points(stream,
std::back_inserter(points),
CGAL::parameters::point_map(CGAL::First_of_pair_property_map<PointVectorPair>())))
{
std::cerr << "Error: cannot read file " << fname << std::endl;
return EXIT_FAILURE;
}
std::cout << "Calculating normals...";
CGAL::Timer t;
t.start();
// Estimates normals direction.
// Note: pca_estimate_normals() requires a range of points
// as well as property maps to access each point's position and normal.
const int nb_neighbors = 18; // K-nearest neighbors = 3 rings
if (argc > 2 && std::strcmp(argv[2], "-r") == 0) // Use a fixed neighborhood radius
{
// First compute a spacing using the K parameter
double spacing
= CGAL::compute_average_spacing<Concurrency_tag>
(points, nb_neighbors,
CGAL::parameters::point_map(CGAL::First_of_pair_property_map<PointVectorPair>()));
// Then, estimate normals with a fixed radius
CGAL::pca_estimate_normals<Concurrency_tag>
(points,
0, // when using a neighborhood radius, K=0 means no limit on the number of neighbors returns
CGAL::parameters::point_map(CGAL::First_of_pair_property_map<PointVectorPair>()).
normal_map(CGAL::Second_of_pair_property_map<PointVectorPair>()).
neighbor_radius(2. * spacing)); // use 2*spacing as neighborhood radius
}
else // Use a fixed number of neighbors
{
CGAL::pca_estimate_normals<Concurrency_tag>
(points, nb_neighbors,
CGAL::parameters::point_map(CGAL::First_of_pair_property_map<PointVectorPair>()).
normal_map(CGAL::Second_of_pair_property_map<PointVectorPair>()));
}
// Orients normals.
// Note: mst_orient_normals() requires a range of points
// as well as property maps to access each point's position and normal.
std::list<PointVectorPair>::iterator unoriented_points_begin =
CGAL::mst_orient_normals(points, nb_neighbors,
CGAL::parameters::point_map(CGAL::First_of_pair_property_map<PointVectorPair>()).
normal_map(CGAL::Second_of_pair_property_map<PointVectorPair>()));
// Optional: delete points with an unoriented normal
// if you plan to call a reconstruction algorithm that expects oriented normals.
points.erase(unoriented_points_begin, points.end());
/// Convert pointvectorpair to pointvector
std::vector<PNI> points2; // store points
points2 = convert_function(points);
///
/// use points + normals for surface reconstruction
///
// Shape detection
Efficient_ransac ransac;
ransac.set_input(points2);
ransac.add_shape_factory<Plane>();
std::cout << "Extracting planes...";
t.reset();
ransac.detect();
Efficient_ransac::Plane_range planes = ransac.planes();
std::size_t num_planes = planes.size();
std::cout << " Done. " << num_planes << " planes extracted. Time: " << t.time() << " sec." << std::endl;
// Stores the plane index of each point as the third element of the tuple.
Point_to_shape_index_map shape_index_map(points, planes);
for (std::size_t i = 0; i < points.size(); ++i) {
// Uses the get function from the property map that accesses the 3rd element of the tuple.
int plane_index = get(shape_index_map, i);
points2[i].get<2>() = plane_index;
}
//////////////////////////////////////////////////////////////////////////
std::cout << "Generating candidate faces...";
t.reset();
Polygonal_surface_reconstruction algo(
points,
Point_map(),
Normal_map(),
Plane_index_map()
);
std::cout << " Done. Time: " << t.time() << " sec." << std::endl;
//////////////////////////////////////////////////////////////////////////
Surface_mesh model;
std::cout << "Reconstructing...";
t.reset();
if (!algo.reconstruct<MIP_Solver>(model)) {
std::cerr << " Failed: " << algo.error_message() << std::endl;
return EXIT_FAILURE;
}
const std::string& output_file("data/cube_result_47038.off");
std::ofstream output_stream(output_file.c_str());
if (output_stream && CGAL::write_off(output_stream, model)) {
// flush the buffer
output_stream << std::flush;
std::cout << " Done. Saved to " << output_file << ". Time: " << t.time() << " sec." << std::endl;
}
else {
std::cerr << " Failed saving file." << std::endl;
return EXIT_FAILURE;
}
//////////////////////////////////////////////////////////////////////////
// Also stores the candidate faces as a surface mesh to a file
Surface_mesh candidate_faces;
algo.output_candidate_faces(candidate_faces);
const std::string& candidate_faces_file("data/cube_candidate_faces_47038.off");
std::ofstream candidate_stream(candidate_faces_file.c_str());
if (candidate_stream && CGAL::write_off(candidate_stream, candidate_faces)) {
// flush the buffer
output_stream << std::flush;
std::cout << "Candidate faces saved to " << candidate_faces_file << "." << std::endl;
}
return EXIT_SUCCESS;
}
#else
int main(int, char**)
{
std::cerr << "This test requires either GLPK or SCIP.\n";
return EXIT_SUCCESS;
}
#endif // defined(CGAL_USE_GLPK) || defined(CGAL_USE_SCIP)
谢谢,
芭芭拉
解决方案
您尝试在列表中使用 operator[] :
std::vector<PNI> convert_function(const std::list<PointVectorPair>& list)
{
...
for (std::size_t i = 0; i < list.size(); ++i)
{
out[i].get<0>() = list[i].get<0>();
你根本做不到,list[i]
因为std::list
不提供它。
https://en.cppreference.com/w/cpp/container/list
相反,您可以使用:
auto it = list.begin();
和
it != list.end();
和
it++;
结合
(*it).get<0>();
这使您可以从头开始,检测结尾并推进迭代器。结合迭代器的取消引用。
在尝试写入之前,您还需要resize()
您的向量。out
推荐阅读
- python - 通过 array[rows,cols] 规范更新二维数组?
- javascript - 将 .xml 转换为 JSON
- java - 需要修复我的递归代码(get 需要返回 int 或 stackoverflow 错误)
- django - Bootstrap 模态提交按钮在 Django 2.2 中不能与 DeleteView 一起使用
- jenkins-pipeline - 如何在 Jenkinsfile 中将构建部署到不同的 Artifactory?
- javascript - 当我重定向到另一个页面时,本地存储已损坏
- sql-server - 如何显示两个日期之间的日期
- java - 获取亚马逊页面和产品信息的最佳方式
- django - django分析没有模型的文件的heroku限制?
- laravel - 使用外部 Laravel 护照流明 api 进行 Laravel 客户端身份验证