我有一个PCL点云(3D(，我想将其转换为地面三角网(2.5D(，然后对点(2D(进行采样，以在投影到三角网上时找到它们的高程。为了做到这一点，我一直在使用CGAL Delaunay三角测量类，它基本上运行得很好！
我能够使用建立在三角网_vertex_base_with_info_2上的Delaunay_triangulation_2来实现这一点，并创建一个好看的三角网。我还编写了一个函数，使用CGAL locate((函数提取二维空间中任意点的面和顶点，这样，如果点投影到三角网上，我就可以对点的高度进行插值。我需要信息字段来保存一个索引，使我能够将三角测量中的顶点与PCL点云结构中的点相关联。
然而，当使用基本的三角测量类时，locate((函数很慢(从三角测量中的任意顶点开始随机行走(，因为我必须对云中的每个查询点进行插值(以估计投影高度(，所以这是目前我整个管道中最慢的部分。因此，我考虑使用Triangulation Hierarchy类来提高效率。
我不知道如何让Triangulation_hierarchy类使用带信息的顶点库，我认为我只是做错了一些愚蠢的事情。下面是一个模拟示例，显示了我使用简单的三角测量结构(没有层次结构(的慢速解决方案，它确实有效：
#include <chrono>
#include <pcl/point_types.h>
#include <pcl/filters/passthrough.h>
#include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
#include <CGAL/Triangulation_vertex_base_with_info_2.h>
#include <CGAL/Delaunay_triangulation_2.h>
typedef CGAL::Exact_predicates_inexact_constructions_kernel K;
typedef CGAL::Triangulation_vertex_base_with_info_2<unsigned int, K> Vb;
typedef CGAL::Triangulation_data_structure_2<Vb> Tds;
typedef CGAL::Delaunay_triangulation_2<K, Tds> Delaunay;
typedef Delaunay::Point_2 CGALPoint;
typedef Delaunay::Face_handle Face_handle;
// This is templated on a PCL PointCloud Ptr - input cloud is basically a 3D vector of XYZ points
template <typename CloudType>
void delaunayTriangulation(CloudType input_cloud, Delaunay& triangulation)
{
std::cout << "Performing Delaunay triangulation on cloud of size " << input_cloud->points.size() << std::endl;
// Convert ground minima cloud to CGAL vector of points
std::vector< std::pair<CGALPoint, unsigned> > minima_vec;
for(std::size_t i=0; i<input_cloud->points.size(); i++)
{
minima_vec.push_back(std::make_pair(CGALPoint(input_cloud->points[i].x,input_cloud->points[i].y), i));
}
// Generate Delaunay Triangulation for ground minima
triangulation = Delaunay(minima_vec.begin(), minima_vec.end());
std::cout << " Number of vertices in Delaunay: " << triangulation.number_of_vertices() << std::endl;
std::cout << " Number of faces in Delaunay: " << triangulation.number_of_faces() << std::endl;
}
int main()
{
// Generate a starting point cloud with random points
pcl::PointCloud<pcl::PointXYZ>::Ptr cloud(new pcl::PointCloud<pcl::PointXYZ>);
for(int i=0; i<500; i++)
for(int j=0; j<500; j++)
{
// Generate points which are gridded + a bit of noise in XY, and random Z
pcl::PointXYZ point;
point.x = i + (std::rand()%100)/100.0;
point.y = j + (std::rand()%100)/100.0;
point.z = std::rand();
cloud->points.push_back(point);
}
// Get the ground triangulation
Delaunay triangulation;
delaunayTriangulation(cloud, triangulation);
// Locate the containing face for a bunch of random points
std::cout << "Starting to search for faces..." << std::endl;
auto start_time = std::chrono::high_resolution_clock::now();
for(int i=0; i<3000000; i++)
{
// Random point with X and Y between 0 and 500
CGALPoint test_point((std::rand() % 5000)/10.0, (std::rand() % 5000)/10.0);
Face_handle face = triangulation.locate(test_point);
// here we would do some math using the vertices located above
}
auto stop_time = std::chrono::high_resolution_clock::now();
float duration = std::chrono::duration_cast<std::chrono::microseconds>( stop_time - start_time ).count();
std::cout << "Finished searching for faces - took " << duration/1000000 << std::endl;
}
如果我改为使用在Delaunay_Triangulation_2类型之上构建的Triangulation_hierarchy_2对象，它将不允许我插入包含信息字段的点对-只有当我自己使用点向量构建对象时，它才会编译：
#include <chrono>
#include <pcl/point_types.h>
#include <pcl/filters/passthrough.h>
#include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
#include <CGAL/Triangulation_vertex_base_with_info_2.h>
#include <CGAL/Delaunay_triangulation_2.h>
#include <CGAL/Triangulation_hierarchy_2.h>
typedef CGAL::Exact_predicates_inexact_constructions_kernel K;
typedef CGAL::Triangulation_vertex_base_with_info_2<unsigned int, K> Vbb;
typedef CGAL::Triangulation_hierarchy_vertex_base_2<Vbb> Vb;
typedef CGAL::Triangulation_data_structure_2<Vb> Tds;
typedef CGAL::Delaunay_triangulation_2<K, Tds> Delaunay;
typedef Delaunay::Point_2 CGALPoint;
typedef Delaunay::Face_handle Face_handle;
typedef CGAL::Triangulation_hierarchy_2<Delaunay> Delaunay_hierarchy;
// This is templated on a PCL PointCloud Ptr - input cloud is basically a 3D vector of XYZ points
template <typename CloudType>
void delaunayTriangulation(CloudType input_cloud, Delaunay_hierarchy& triangulation)
{
std::cout << "Performing Delaunay triangulation on cloud of size " << input_cloud->points.size() << std::endl;
// Convert ground minima cloud to CGAL vector of points
std::vector<CGALPoint> minima_vec_simple;
for(std::size_t i=0; i<input_cloud->points.size(); i++)
{
minima_vec_simple.push_back(CGALPoint(input_cloud->points[i].x,input_cloud->points[i].y));
}
// Generate Delaunay Triangulation for ground minima
triangulation = Delaunay_hierarchy(minima_vec_simple.begin(), minima_vec_simple.end());
std::cout << " Number of vertices in Delaunay: " << triangulation.number_of_vertices() << std::endl;
std::cout << " Number of faces in Delaunay: " << triangulation.number_of_faces() << std::endl;
}
int main()
{
// Generate a starting point cloud with random points
pcl::PointCloud<pcl::PointXYZ>::Ptr cloud(new pcl::PointCloud<pcl::PointXYZ>);
for(int i=0; i<500; i++)
for(int j=0; j<500; j++)
{
// Generate points which are gridded + a bit of noise in XY, and random Z
pcl::PointXYZ point;
point.x = i + (std::rand()%100)/100.0;
point.y = j + (std::rand()%100)/100.0;
point.z = std::rand();
cloud->points.push_back(point);
}
// Get the ground triangulation
Delaunay_hierarchy triangulation;
delaunayTriangulation(cloud, triangulation);
// Locate the containing face for a bunch of random points
std::cout << "Starting to search for faces..." << std::endl;
auto start_time = std::chrono::high_resolution_clock::now();
for(int i=0; i<3000000; i++)
{
// Random point with X and Y between 0 and 500
CGALPoint test_point((std::rand() % 5000)/10.0, (std::rand() % 5000)/10.0);
Face_handle face = triangulation.locate(test_point);
// here we would do some math using the vertices located above
}
auto stop_time = std::chrono::high_resolution_clock::now();
float duration = std::chrono::duration_cast<std::chrono::microseconds>( stop_time - start_time ).count();
std::cout << "Finished searching for faces - took " << duration/1000000 << std::endl;
}
CGAL中的所有typedef和模板层对我来说都有点不透明——有没有什么方法可以让我设置这些层，让我可以基于Triangulation_hierarchy_vertex_base_with_info_2构建类似Triangulation_hierarchy_2的东西？第二种类型看起来并不存在。
谢谢！
Ok-我最终尝试了一些路由，并对每个路由的执行时间进行了一些基本的基准测试。这些是基于我的问题中代码中的设置：
- 三角网，25000个点，主要分布在XY中0到500的网格上，XY中有一点噪声，z值随机
- 在X/Y中用3000000个随机值介于0和500之间的点测试云
我尝试了什么：
- CGAL locate((使用简单的三角结构，任意起始面->79.3秒
- CGAL locate((使用三角形层次结构，任意起始面->4.74秒
- CGAL locate((使用最近相邻顶点给定的起始面，使用顶点输入云上的K-D树找到-->3.41秒
为此，在创建三角测量后，我在面上迭代，并在输入云顶点索引和三角测量面句柄之间进行映射：
std::vector<Face_handle> face_mapping(cloud->points.size());
std::vector<bool> faces_filled(cloud->points.size(), false);
auto start_time_list = std::chrono::high_resolution_clock::now();
// Iterate over all faces in triangulation
for (Face_handle face : triangulation.finite_face_handles())
// Iterate over 3 vertices for each face
for(int i=0; i<3; i++)
{
int index = uint32_t(face->vertex(i)->info());
if(!faces_filled[index])
face_mapping[index] = face;
}
然后，当我们运行点定位搜索时，我们可以做这样的事情：
std::vector<int> nearest_indices;
std::vector<float> nearest_dists;
pcl::KdTreeFLANN<pcl::Point2DGround> tree;
pcl::PointCloud<pcl::Point2DGround>::Ptr cloud_2d(new pcl::PointCloud<pcl::Point2DGround>);
copyPointCloud3D(cloud, cloud_2d);
tree.setInputCloud(cloud_2d);
for(int i=0; i<3000000; i++)
{
// Random point with X and Y between 0 and 500
CGALPoint test_point((std::rand() % 5000)/10.0, (std::rand() % 5000)/10.0);
// Get closest vertex in triangulation
pcl::Point2DGround pcl_point;
pcl_point.x = test_point.x();
pcl_point.y = test_point.y();
tree.nearestKSearch(pcl_point, 1, nearest_indices, nearest_dists);
Face_handle face = triangulation.locate(test_point, face_mapping[nearest_indices[0]]);
// here we would do some math using the vertices located above
}
其中pcl:：Point2DGround必须是某种自定义点类型，其中表示仅为2D(以使搜索树正常工作(。
我最终没有对要测试的传入点进行任何排序，因为上面的方法有效，而且我要定位的测试点((-ed总是比三角网顶点本身多得多，所以我认为对它们强制排序可能会更昂贵。
我认为这是我的情况下最简单/最快的解决方案，所以我会接受这个！感谢@marcglisse和@andreasfabri的评论。