【发布时间】：2015-11-21 04:03:56
【问题描述】：
我遇到了我需要完成的任务的性能问题。目前的瓶颈之一是从非结构化网格中获取插值字段值。
慢速部分是，给定一个 2D 点和一个非结构化的 2D 网格，找到紧邻该点的网格点。最好能找到它所在的三角形。
现在我正在使用 CGAL，但它太慢了。使用当前的实现，整个任务需要几天才能完成，在高端 CPU 上并行运行。
我认为慢的部分是 CGAL::natural_neighbor_coordinates_2。
#ifndef FIELD_INTERPOLATOR_H
#define FIELD_INTERPOLATOR_H
#include "Vec.h"
#include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
#include <CGAL/Delaunay_triangulation_2.h>
#include <CGAL/Interpolation_traits_2.h>
#include <CGAL/natural_neighbor_coordinates_2.h>
#include <CGAL/interpolation_functions.h>
#include <map>
#include <vector>
typedef CGAL::Exact_predicates_inexact_constructions_kernel Kernel;
typedef CGAL::Delaunay_triangulation_2< Kernel > Delaunay_triangulation;
typedef Kernel::FT FieldType;
typedef Kernel::Point_2 MeshType;
struct FieldInterpolator23 {
Delaunay_triangulation m_triangulation;
std::map< MeshType, FieldType, Kernel::Less_xy_2 > m_vX;
std::map< MeshType, FieldType, Kernel::Less_xy_2 > m_vY;
std::map< MeshType, FieldType, Kernel::Less_xy_2 > m_vZ;
typedef CGAL::Data_access< std::map< MeshType, FieldType, Kernel::Less_xy_2 > > ValueAccess;
FieldInterpolator23() {}
FieldInterpolator23( const std::vector< TN::Vec2 > & mesh, const std::vector< TN::Vec3 > & field )
{
const int N = mesh.size();
for ( int i = 0; i < N; ++i ) {
MeshType p( mesh[i].x(), mesh[i].y() );
m_triangulation.insert( p );
m_vX.insert( std::make_pair( p, field[i].x() ) );
m_vY.insert( std::make_pair( p, field[i].y() ) );
m_vZ.insert( std::make_pair( p, field[i].z() ) );
}
}
void set( const std::vector< TN::Vec2 > & mesh, const std::vector< TN::Vec3 > & field ) {
m_triangulation.clear();
m_vX.clear();
m_vY.clear();
m_vZ.clear();
const int N = mesh.size();
for ( int i = 0; i < N; ++i ) {
MeshType p( mesh[i].x(), mesh[i].y() );
m_triangulation.insert( p );
m_vX.insert( std::make_pair( p, field[i].x() ) );
m_vY.insert( std::make_pair( p, field[i].y() ) );
m_vZ.insert( std::make_pair( p, field[i].z() ) );
}
}
TN::Vec3 operator() ( TN::Vec2 p ) {
MeshType pos( p.x(), p.y() );
std::vector< std::pair< MeshType, FieldType > > coords;
FieldType norm =
CGAL::natural_neighbor_coordinates_2( m_triangulation, pos, std::back_inserter( coords ) ).second;
FieldType resX =
CGAL::linear_interpolation(
coords.begin(),
coords.end(),
norm,
ValueAccess( m_vX )
);
FieldType resY =
CGAL::linear_interpolation(
coords.begin(),
coords.end(),
norm,
ValueAccess( m_vY )
);
FieldType resZ =
CGAL::linear_interpolation(
coords.begin(),
coords.end(),
norm,
ValueAccess( m_vZ )
);
return TN::Vec3( resX, resY, resZ );
}
};
#endif
谁能指出一个可接受的更高性能解决方案的方向，无论是不同的库还是算法？
【问题讨论】：
-
三角剖分的顶点数是多少？您需要在三角剖分中定位多少个查询点？
-
大约有 100 万个网格点，我需要做大约 20 亿次查询。
-
如果你事先知道查询，它有助于对它们进行排序（比如沿着空间填充曲线，CGAL 有相应的功能），然后按此顺序进行点位置查询，传递每个查询以前的位置作为起点。这样，大多数查询只会检查起点是否已经是直角三角形，而其他查询最多需要步行 1 或 2 步才能到达。
-
如果您的查询点在规则网格上，一个快速的方法是三角形填充：创建一个像素对应查询点的图像，并使用三角形 ID 作为填充颜色填充每个三角形.这将使查询成为恒定时间。 （对于如此大量的查询，这在您的情况下可能不切实际 - 除非您可以通过子窗口工作。）
标签： c++ performance computational-geometry cgal