Search code examples
c++cgal

2D AABBTree in CGAL with custom property map


I am trying to create an AABBTree for checking in which triangle of a 2D mesh some points are. To do that, I am using CGAL, in particular, the Polygon_mesh_processing package and the location functions. However, this feels an impossible task to me.

Since I want to check several points and I care about performance, I want to build an AABB first. To do that, I need the PMP::build_AABB_tree function, which receives a mesh and the AABBTree.

In this particular case, the mesh is of type CGAL::Surface_mesh<K::Point_2>, with a Simple_cartersian kernel. The query points are also 2D. According to documentation, I need to provide a ReadablePropertyMap which converts my 2D points into 3D points.

I do this with the following:

struct Point3VPM {
        using key_type = Mesh::Vertex_index;
        using value_type = VNCS::Sim2D::Kernel::K::Point_3;
        using reference = value_type;
        using category = boost::readable_property_map_tag;

        Point3VPM(const Mesh &m)
            : mesh(std::cref(m))
        {
        }

        friend Point3VPM::value_type get(const Point3VPM &map, Mesh::Vertex_index idx)
        {
            Point p = map.mesh.get().point(idx);
            return {p[0], p[1], 0};
        }

        std::reference_wrapper<const Mesh> mesh;
    };

Here starts my first issue. My AABBTree type is the following:

    using AABBTreeTraits =
        CGAL::AABB_traits<VNCS::Sim2D::Kernel::K, CGAL::AABB_face_graph_triangle_primitive<Mesh, Point3VPM>>;
    using AABBTree = CGAL::AABB_tree<AABBTreeTraits>;

How can I construct my AABBTree? My Point3VPM needs the SurfaceMesh to be able to convert a Vertex_index into a Point_3, but I can't find a way of building the AABBTreeTraits with the SurfaceMesh. CGAL expects a default constructor for the property map but since my property map can't work without th mesh, I need to provide the mesh somehow. How can I fix this?

Any help trying to make this work would be much appreciated

Here you can find the code for a minimal example:

#include <CGAL/Simple_cartesian.h>
#include <CGAL/Surface_mesh.h>
#include <CGAL/AABB_tree.h>
#include <CGAL/AABB_traits.h>
#include <CGAL/AABB_face_graph_triangle_primitive.h>
#include <CGAL/Polygon_mesh_processing/locate.h>

namespace
{
using K = CGAL::Simple_cartesian<double>;
using Point2 = K::Point_2;
using Point3 = K::Point_3;
using Mesh = CGAL::Surface_mesh<Point2>;
namespace PMP = CGAL::Polygon_mesh_processing;

struct Point3VPM {
    using key_type = Mesh::Vertex_index;
    using value_type = Point3;
    using reference = value_type;
    using category = boost::readable_property_map_tag;

    Point3VPM(const Mesh &m)
        : mesh(std::cref(m))
    {
    }

    friend Point3VPM::value_type get(const Point3VPM &map, Mesh::Vertex_index idx)
    {
        Point2 p = map.mesh.get().point(idx);
        return {p[0], p[1], 0};
    }

    std::reference_wrapper<const Mesh> mesh;
};

using AABBTreeTraits = CGAL::AABB_traits<K, CGAL::AABB_face_graph_triangle_primitive<Mesh, Point3VPM>>;
using AABBTree = CGAL::AABB_tree<AABBTreeTraits>;
}  // namespace

int main()
{
    Mesh mesh;

    AABBTree tree;
    PMP::build_AABB_tree(mesh, tree);
    auto location = PMP::locate_with_AABB_tree({0, 0}, tree, mesh);
}

Solution

  • Try the following:

    #include <CGAL/Simple_cartesian.h>
    #include <CGAL/Surface_mesh.h>
    #include <CGAL/AABB_tree.h>
    #include <CGAL/AABB_traits.h>
    #include <CGAL/AABB_face_graph_triangle_primitive.h>
    #include <CGAL/Polygon_mesh_processing/locate.h>
    
    namespace
    {
    using K = CGAL::Simple_cartesian<double>;
    using Point2 = K::Point_2;
    using Point3 = K::Point_3;
    using Mesh = CGAL::Surface_mesh<Point2>;
    namespace PMP = CGAL::Polygon_mesh_processing;
    
    struct Point3VPM {
        using key_type = Mesh::Vertex_index;
        using value_type = Point3;
        using reference = value_type;
        using category = boost::readable_property_map_tag;
    
        Point3VPM()
            : mesh_ptr(nullptr)
        {}
    
        Point3VPM(const Mesh &m)
            : mesh_ptr(&m)
        {}
    
        friend Point3VPM::value_type get(const Point3VPM &map, Mesh::Vertex_index idx)
        {
            Point2 p = map.mesh_ptr->point(idx);
            return {p[0], p[1], 0};
        }
    
        const Mesh* mesh_ptr;
    };
    
    using AABBTreeTraits = CGAL::AABB_traits<K, CGAL::AABB_face_graph_triangle_primitive<Mesh, Point3VPM>>;
    using AABBTree = CGAL::AABB_tree<AABBTreeTraits>;
    }  // namespace
    
    int main()
    {
        Mesh mesh;
        Point3VPM vpm(mesh);
      
        AABBTree tree;
        PMP::build_AABB_tree(mesh, tree, CGAL::parameters::vertex_point_map(vpm));
        auto location = PMP::locate_with_AABB_tree({0, 0}, tree, mesh);
    }
    
    

    Note that you should use CGAL::Exact_predicates_inexact_constructions_kernel instead of Simple_cartesian<double> to get exact predicates.