CGAL 5.6 - 2D Generalized Barycentric Coordinates
Barycentric_coordinates_2/affine_coordinates.cpp
#include <Eigen/Core>
#include <Eigen/Dense>
#include <CGAL/barycenter.h>
#include <CGAL/Simple_cartesian.h>
#include <boost/iterator/transform_iterator.hpp>
// Typedefs.
using Point_2 = Kernel::Point_2;
using VectorXd = Eigen::VectorXd;
using MatrixXd = Eigen::MatrixXd;
using Output_iterator =
std::back_insert_iterator< std::vector<double> >;
int main() {
// Create a set of vertices.
const std::vector<Point_2> vertices = {
Point_2(0.0, 0.0), Point_2(0.75, 0.25), Point_2(0.5, 0.5), Point_2(0.4, -0.2) };
// Create a set of query points.
const std::vector<Point_2> queries = {
Point_2(0.2, 0.2), Point_2(0.3, 0.3), Point_2(0.4, 0.4) };
// Create a lambda function with affine coordinates.
// This implementation is based on the following paper:
// S. Waldron. Affine generalized barycentric coordinates.
// Jaen Journal on Approximation, 3(2):209-226, 2011.
// This function is a model of the `AnalyticWeights_2` concept.
const auto affine = [&](
const Point_2& query,
Output_iterator coordinates) {
const std::size_t n = vertices.size();
const auto lambda = [](const Point_2& p){ return std::make_pair(p, 1.0); };
const Point_2 b = CGAL::barycenter(
boost::make_transform_iterator(vertices.begin(), lambda),
boost::make_transform_iterator(vertices.end() , lambda), Kernel());
MatrixXd V(2, n);
for (std::size_t i = 0; i < n; ++i) {
V(0, i) = vertices[i].x() - b.x();
V(1, i) = vertices[i].y() - b.y();
}
const auto A = V.adjoint();
const auto mat = V * A;
const auto inv = mat.inverse();
Point_2 diff; VectorXd vec(2);
for (std::size_t i = 0; i < n; ++i) {
const double x = query.x() - b.x();
const double y = query.y() - b.y();
diff = Point_2(x, y);
vec(0) = V(0, i);
vec(1) = V(1, i);
const auto res = inv * vec;
*(coordinates++) =
diff.x() * res(0) + diff.y() * res(1) + 1.0 / double(n);
}
};
// Compute affine coordinates for all query points.
std::cout << std::endl << "affine coordinates (all queries): " << std::endl << std::endl;
std::vector<double> coordinates;
coordinates.reserve(4);
for (const auto& query : queries) {
coordinates.clear();
affine(query, std::back_inserter(coordinates));
for (std::size_t i = 0; i < coordinates.size() - 1; ++i) {
std::cout << coordinates[i] << ", ";
}
std::cout << coordinates[coordinates.size() - 1] << std::endl;
}
std::cout << std::endl;
return EXIT_SUCCESS;
}