C++ 几何计算库

16 篇文章 0 订阅
9 篇文章 0 订阅
代码

#include <iostream>
#include <list>
#include <CGAL/Simple_cartesian.h>
#include <CGAL/AABB_tree.h>
#include <CGAL/AABB_traits.h>
#include <CGAL/AABB_segment_primitive.h>
#include <CGAL/Polygon_2.h>


typedef CGAL::Simple_cartesian<double> K;


// custom point type
struct My_point {
    double m_x;
    double m_y;
    double m_z;

    My_point(const double x,
        const double y,
        const double z)
        : m_x(x), m_y(y), m_z(z) {}
};

// custom triangle type with
// three pointers to points
struct My_triangle {
    My_point* m_pa;
    My_point* m_pb;
    My_point* m_pc;

    My_triangle(My_point* pa,
        My_point* pb,
        My_point* pc)
        : m_pa(pa), m_pb(pb), m_pc(pc) {}
};

// the custom triangles are stored into a vector
typedef std::vector<My_triangle>::const_iterator Iterator;

// The following primitive provides the conversion facilities between
// the custom triangle and point types and the CGAL ones
struct My_triangle_primitive {
public:

    // this is the type of data that the queries returns. For this example
    // we imagine that, for some reasons, we do not want to store the iterators
    // of the vector, but raw pointers. This is to show that the Id type
    // does not have to be the same as the one of the input parameter of the
    // constructor.
    typedef const My_triangle* Id;

    // CGAL types returned
    typedef K::Point_3    Point; // CGAL 3D point type
    typedef K::Triangle_3 Datum; // CGAL 3D triangle type

private:
    Id m_pt; // this is what the AABB tree stores internally

public:
    My_triangle_primitive() {} // default constructor needed

    // the following constructor is the one that receives the iterators from the
    // iterator range given as input to the AABB_tree
    My_triangle_primitive(Iterator it)
        : m_pt(&(*it)) {}

    const Id& id() const { return m_pt; }

    // utility function to convert a custom
    // point type to CGAL point type.
    Point convert(const My_point* p) const
    {
        return Point(p->m_x, p->m_y, p->m_z);
    }

    // on the fly conversion from the internal data to the CGAL types
    Datum datum() const
    {
        return Datum(convert(m_pt->m_pa),
            convert(m_pt->m_pb),
            convert(m_pt->m_pc));
    }

    // returns a reference point which must be on the primitive
    Point reference_point() const
    {
        return convert(m_pt->m_pa);
    }
};

/*
自定义KDTree测试,点相交
*/
int testKDTree()
{
    typedef CGAL::AABB_traits<K, My_triangle_primitive> My_AABB_traits;
    typedef CGAL::AABB_tree<My_AABB_traits> MyTree;
    typedef K::FT FT;

    My_point a(1.0, 0.0, 0.0);
    My_point b(0.0, 1.0, 0.0);
    My_point c(0.0, 0.0, 1.0);
    My_point d(0.0, 0.0, 0.0);

    std::vector<My_triangle> triangles;
    triangles.push_back(My_triangle(&a, &b, &c));
    triangles.push_back(My_triangle(&a, &b, &d));
    triangles.push_back(My_triangle(&a, &d, &c));

    // constructs AABB tree
    MyTree tree(triangles.begin(), triangles.end());

    // counts #intersections
    K::Ray_3 ray_query(K::Point_3(1.0, 0.0, 0.0), K::Point_3(0.0, 1.0, 0.0));
    std::cout << tree.number_of_intersected_primitives(ray_query)
        << " intersections(s) with ray query" << std::endl;

    // computes closest point
    K::Point_3 point_query(2.0, 2.0, 2.0);
    K::Point_3 closest_point = tree.closest_point(point_query);
    std::cerr << "closest point is: " << closest_point << std::endl;

    FT sqd = tree.squared_distance(point_query);
    std::cout << "squared distance: " << sqd << std::endl;
    return EXIT_SUCCESS;

/* 输出
3 intersections(s) with ray query
closest point is: 0.333333 0.333333 0.333333
*/
}

/*
光线追踪,面相交
*/
void testGlyph() {
    typedef K::FT FT;
    typedef K::Segment_3 Segment;
    typedef K::Point_3 Point;

    typedef std::list<Segment> SegmentRange;
    typedef SegmentRange::const_iterator Iterator;
    typedef CGAL::AABB_segment_primitive<K, Iterator> Primitive;
    typedef CGAL::AABB_traits<K, Primitive> Traits;
    typedef CGAL::AABB_tree<Traits> Tree;
    typedef Tree::Point_and_primitive_id Point_and_primitive_id;

    Point a(0.0, 0.0, 1);
    Point b(2.0, 1.0, 1);
    Point c(3.0, 4.0, 1);
    Point d(1.0, 6.0, 1);
    Point e(-1.0, 3.0, 1);

    std::list<Segment> seg;
    seg.push_back(Segment(a, b));
    seg.push_back(Segment(b, c));
    seg.push_back(Segment(c, d));
    seg.push_back(Segment(d, e));
    seg.push_back(Segment(e, a));

    // constructs the AABB tree and the internal search tree for
    // efficient distance computations.
    Tree tree(seg.begin(), seg.end());
    tree.build();

    tree.accelerate_distance_queries();

    // counts #intersections with a segment query
    Segment segment_query(Point(1.0, 0.0, 1), Point(0.0, 7.0, 1));
    std::cout << tree.number_of_intersected_primitives(segment_query)
        << " intersections(s) with segment" << std::endl;

    // computes the closest point from a point query
    Point point_query(1.5, 3.0, 1);
    Point closest = tree.closest_point(point_query);
    std::cerr << "closest point is: " << closest << std::endl;

    Point_and_primitive_id id = tree.closest_point_and_primitive(point_query);
    std::cout << id.second->source() << " " << id.second->target() << std::endl;
}

/*
测试布尔运算
*/

#include <CGAL/Exact_predicates_exact_constructions_kernel.h>
#include <CGAL/Polygon_2.h>
#include <CGAL/Polygon_with_holes_2.h>
#include <CGAL/Polygon_set_2.h>
#include <CGAL/Boolean_set_operations_2.h>
#include <CGAL/Intersections_3/Line_3_Line_3.h>
//-----------------------------------------------------------------------------
// Pretty-print a CGAL polygon.
//
template<class Kernel, class Container>
void print_polygon(const CGAL::Polygon_2<Kernel, Container>& P)
{
    typename CGAL::Polygon_2<Kernel, Container>::Vertex_const_iterator  vit;

    std::cout << "[ " << P.size() << " vertices:";
    for (vit = P.vertices_begin(); vit != P.vertices_end(); ++vit)
        std::cout << " (" << *vit << ')';
    std::cout << " ]" << std::endl;

    return;
}

//-----------------------------------------------------------------------------
// Pretty-print a polygon with holes.
//
template<class Kernel, class Container>
void print_polygon_with_holes
(const CGAL::Polygon_with_holes_2<Kernel, Container>& pwh)
{
    if (!pwh.is_unbounded())
    {
        std::cout << "{ Outer boundary = ";
        print_polygon(pwh.outer_boundary());
    }
    else
        std::cout << "{ Unbounded polygon." << std::endl;

    typename CGAL::Polygon_with_holes_2<Kernel, Container>::
        Hole_const_iterator  hit;
    unsigned int                                                     k = 1;

    std::cout << "  " << pwh.number_of_holes() << " holes:" << std::endl;
    for (hit = pwh.holes_begin(); hit != pwh.holes_end(); ++hit, ++k)
    {
        std::cout << "    Hole #" << k << " = ";
        print_polygon(*hit);
    }
    std::cout << " }" << std::endl;

    return;
}

void testBoolOpeation() {
    typedef CGAL::Exact_predicates_exact_constructions_kernel Kernel;
    typedef Kernel::Point_2                                   Point_2;
    typedef CGAL::Polygon_2<Kernel>                           Polygon_2;
    typedef CGAL::Polygon_with_holes_2<Kernel>                Polygon_with_holes_2;
    typedef CGAL::Polygon_set_2<Kernel>                       Polygon_set_2;
    typedef std::list<Polygon_with_holes_2>                   Pwh_list_2;
    typedef CGAL::Line_3<K>           Line;
    typedef CGAL::Point_3<K>          Point;

    // Construct the two initial polygons and the clipping rectangle.
    Polygon_2 P;
    P.push_back(Point_2(0, 1));
    P.push_back(Point_2(2, 0));
    P.push_back(Point_2(1, 1));
    P.push_back(Point_2(2, 2));

    Polygon_2 Q;
    Q.push_back(Point_2(3, 1));
    Q.push_back(Point_2(1, 2));
    Q.push_back(Point_2(2, 1));
    Q.push_back(Point_2(1, 0));

    Polygon_2 rect;
    rect.push_back(Point_2(0, 0));
    rect.push_back(Point_2(3, 0));
    rect.push_back(Point_2(3, 2));
    rect.push_back(Point_2(0, 2));

    // Perform a sequence of operations.
    Polygon_set_2 S;
    S.insert(P);
    S.join(Q);                   // Compute the union of P and Q.
    S.complement();               // Compute the complement.
    S.intersection(rect);        // Intersect with the clipping rectangle.

    // Print the result.
    std::list<Polygon_with_holes_2> res;
    std::list<Polygon_with_holes_2>::const_iterator it;

    std::cout << "The result contains " << S.number_of_polygons_with_holes()
        << " components:" << std::endl;

    S.polygons_with_holes(std::back_inserter(res));
    for (it = res.begin(); it != res.end(); ++it) {
        std::cout << "--> ";
        print_polygon_with_holes(*it);
    }

    // 交集
    if ((CGAL::do_intersect(P, Q)))
        std::cout << "The two polygons intersect in their interior." << std::endl;
    else
        std::cout << "The two polygons do not intersect." << std::endl;

    // union
      // Compute the union of P and Q.
    Polygon_with_holes_2 unionR;
    if (CGAL::join(P, Q, unionR)) {
        std::cout << "The union: ";
        print_polygon_with_holes(unionR);
    }
    else
        std::cout << "P and Q are disjoint and their union is trivial."
        << std::endl;

    // Compute the intersection of P and Q.
    CGAL::intersection(P, Q, std::back_inserter(res));

    // Compute the symmetric difference of P and Q.
    CGAL::symmetric_difference(P, Q, std::back_inserter(res));

    // 直线相交
    const Line l = Line(Point(0, 0, 0), Point(1, 2, 3));
    const Line l_1 = Line(Point(0, 0, 0), Point(-3, -2, -1));
    const CGAL::Object obj1 = CGAL::intersection(l, l_1);
}

/*
凸包检测
*/
#include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
#include <CGAL/point_generators_3.h>
#include <CGAL/Delaunay_triangulation_3.h>
#include <CGAL/Surface_mesh.h>
#include <CGAL/algorithm.h>
#include <CGAL/convex_hull_3_to_face_graph.h>
void testConvexHull() {
    typedef CGAL::Exact_predicates_inexact_constructions_kernel     K;
    typedef K::Point_3                                              Point_3;
    typedef CGAL::Delaunay_triangulation_3<K>                       Delaunay;
    typedef Delaunay::Vertex_handle                                 Vertex_handle;
    typedef CGAL::Surface_mesh<Point_3>                             Surface_mesh;
    CGAL::Random_points_in_sphere_3<Point_3> gen(100.0);
    std::list<Point_3> points;

    // generate 250 points randomly in a sphere of radius 100.0
    // and insert them into the triangulation
    std::copy_n(gen, 250, std::back_inserter(points));
    Delaunay T;
    T.insert(points.begin(), points.end());

    std::list<Vertex_handle>  vertices;
    T.incident_vertices(T.infinite_vertex(), std::back_inserter(vertices));
    std::cout << "This convex hull of the 250 points has "
        << vertices.size() << " points on it." << std::endl;

    // remove 25 of the input points
    std::list<Vertex_handle>::iterator v_set_it = vertices.begin();
    for (int i = 0; i < 25; i++)
    {
        T.remove(*v_set_it);
        v_set_it++;
    }

    //copy the convex hull of points into a polyhedron and use it
    //to get the number of points on the convex hull
    Surface_mesh chull;
    CGAL::convex_hull_3_to_face_graph(T, chull);

    std::cout << "After removal of 25 points, there are "
        << num_vertices(chull) << " points on the convex hull." << std::endl;
}

void test() {
    testKDTree();
    testGlyph();
    testBoolOpeation();
    testConvexHull();
}
输出
3 intersections(s) with ray query
closest point is: 0.333333 0.333333 0.333333
squared distance: 8.33333
2 intersections(s) with segment
closest point is: 2.55 2.65 1
2 1 1 3 4 1
The result contains 2 components:
--> { Outer boundary = [ 10 vertices: (2 2) (1 2) (0 2) (0 1) (0 0) (1 0) (2 0) (3 0) (3 1) (3 2) ]
  1 holes:
    Hole #1 = [ 12 vertices: (3 1) (1.66667 0.333333) (2 0) (1.5 0.25) (1 0) (1.33333 0.333333) (0 1) (1.33333 1.66667) (1 2) (1.5 1.75) (2 2) (1.66667 1.66667) ]
 }
--> { Outer boundary = [ 4 vertices: (1 1) (1.5 0.5) (2 1) (1.5 1.5) ]
  0 holes:
 }
The two polygons intersect in their interior.
The union: { Outer boundary = [ 12 vertices: (1.33333 0.333333) (1 0) (1.5 0.25) (2 0) (1.66667 0.333333) (3 1) (1.66667 1.66667) (2 2) (1.5 1.75) (1 2) (1.33333 1.66667) (0 1) ]
  1 holes:
    Hole #1 = [ 4 vertices: (1.5 1.5) (2 1) (1.5 0.5) (1 1) ]
 }
This convex hull of the 250 points has 88 points on it.
After removal of 25 points, there are 84 points on the convex hull.

GitHub - CGAL/cgal: The public CGAL repository, see the README below


创作不易,小小的支持一下吧!

  • 18
    点赞
  • 2
    收藏
    觉得还不错? 一键收藏
  • 打赏
    打赏
  • 0
    评论
㈠ 点的基本运算 1. 平面上两点之间距离 1 2. 判断两点是否重合 1 3. 矢量叉乘 1 4. 矢量点乘 2 5. 判断点是否在线段上 2 6. 求一点饶某点旋转后的坐标 2 7. 求矢量夹角 2 ㈡ 线段及直线的基本运算 1. 点与线段的关系 3 2. 求点到线段所在直线垂线的垂足 4 3. 点到线段的最近点 4 4. 点到线段所在直线的距离 4 5. 点到折线集的最近距离 4 6. 判断圆是否在多边形内 5 7. 求矢量夹角余弦 5 8. 求线段之间的夹角 5 9. 判断线段是否相交 6 10.判断线段是否相交但不交在端点处 6 11.求线段所在直线的方程 6 12.求直线的斜率 7 13.求直线的倾斜角 7 14.求点关于某直线的对称点 7 15.判断两条直线是否相交及求直线交点 7 16.判断线段是否相交,如果相交返回交点 7 ㈢ 多边形常用算法模块 1. 判断多边形是否简单多边形 8 2. 检查多边形顶点的凸凹性 9 3. 判断多边形是否凸多边形 9 4. 求多边形面积 9 5. 判断多边形顶点的排列方向,方法一 10 6. 判断多边形顶点的排列方向,方法二 10 7. 射线法判断点是否在多边形内 10 8. 判断点是否在凸多边形内 11 9. 寻找点集的graham算法 12 10.寻找点集凸包的卷包裹法 13 11.判断线段是否在多边形内 14 12.求简单多边形的重心 15 13.求凸多边形的重心 17 14.求肯定在给定多边形内的一个点 17 15.求从多边形外一点出发到该多边形的切线 18 16.判断多边形的核是否存在 19 ㈣ 圆的基本运算 1 .点是否在圆内 20 2 .求不共线的三点所确定的圆 21 ㈤ 矩形的基本运算 1.已知矩形三点坐标,求第4点坐标 22 ㈥ 常用算法的描述 22 ㈦ 补充 1.两圆关系: 24 2.判断圆是否在矩形内: 24 3.点到平面的距离: 25 4.点是否在直线同侧: 25 5.镜面反射线: 25 6.矩形包含: 26 7.两圆交点: 27 8.两圆公共面积: 28 9. 圆和直线关系: 29 10. 内切圆: 30 11. 求切点: 31 12. 线段的左右旋: 31
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包

打赏作者

码力码力我爱你

创作不易,小小的支持一下吧!

¥1 ¥2 ¥4 ¥6 ¥10 ¥20
扫码支付:¥1
获取中
扫码支付

您的余额不足,请更换扫码支付或充值

打赏作者

实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值