天天看點

CGAL - user manual - 2D triangulation

CGAL - user manual - 2D triangulation

manual 連結:https://doc.cgal.org/latest/Triangulation_2/index.html。本文對manual中不解地方進行了補上了自己不恰當的注解,僅供參考。

目錄

    • 0 準備知識
    • 1 定義
    • 2 表示
      • 2.1 The Set of Faces
      • 2.2 A Representation Based on Faces and Vertices
    • 3 Software Design
    • 4 Basic Triangulations
      • 4.1 Description
      • 4.2 Implementation
      • 4.3 Geometric Traits
      • 4.4 Example of a Basic Triangulation
    • 5 Delaunay Triangulations
      • 5.1 Description
      • 5.2 Geometric Traits
      • 5.3 Implementation
      • 5.4 Example: a Delaunay Terrain
      • 5.6 Example: Print Voronoi Diagram Edges Restricted to a Rectangle
    • 6 Regular Triangulations
    • 7 Constrained triangulations
      • 7.1 The Geometry Traits
      • 7.2 The Face Base Class of a Constrained Triangulation
    • 8 Constrained Delaunay Triangulations
      • 8.1 Example: a Constrained Delaunay Triangulation
    • 9 Constrained Triangulations with a Bidirectional Mapping between Constraints and Subconstraints
    • 10 The Triangulation Hierarchy
    • 11 Flexibility

單純形:單純形是代數拓撲中最基本的概念。單純形是三角形和四面體的一種泛化,一個 k 維單純形是指包含 (k+1)個節點的凸多面體。1維單純形就是線段;2維單純形就是三角形;三維單純形就是四面體;

Voronoi圖:又叫泰森多邊形或Dirichlet圖,它是由一組由連接配接兩鄰點直線的垂直平分線組成的連續多邊形組成;

二維三角剖分可以大緻描述為三角形小平面的集合T,使得:

  • 兩個三角面片,要麼不相連,要麼共享一個較低次元的面(邊或頂點);
  • 集合T中的面片相連,構成鄰接關系;
  • 集合T中面片union的域\(U_T\),沒有奇點??(剖分後,面片連成一片??)

單純形的描述參見原文。

可以将CGAL的2D三角剖分視為平面分區,其邊界面為三角形,并覆寫該組頂點的凸包。 該分區的單個無邊界面是凸包的互補面。 在許多應用程式中,例如柯克帕特裡克(Kirkpatrick)的層次結構或增量Delaunay構造,僅處理三角形面很友善。 是以,将一個稱為無限頂點的虛拟頂點以及無限邊緣*和入射到該三角形的無限面添加到了三角剖分中。 每個無限邊均入射到無限頂點和凸包的頂點。 每個無限面均入射到無限頂點和凸包邊。

這段話的了解:在平面外引入一個頂點,将該頂點和平面上凸包的頂點相連,構成邊和面。示意如下圖中的右圖所示(下圖中的左圖認為虛拟頂點在平面的無限遠處??):
CGAL - user manual - 2D triangulation

于是,三角化後的每個邊都和兩個三角形相鄰。

由于三角化的結果是固定大小複雜度的三角面片的集合。CGAL使用了基于面片和頂點而不是基于邊的三角化表述方式。這樣的表述方式可以減少存儲空間,以及加快算法實作速度。

每個三角面片可以通路自己的頂點和三個相鄰的面。每個頂點可以通路相連的一個面,通過這個面,通路到該頂點所有相鄰的面。

三角面的三個頂點按照逆時針方向,index依次為0,1,2.面的neighbors也一樣依次标0,1,2. 關系如下圖所示:

CGAL - user manual - 2D triangulation

ccw(i)

cw(i)

分别為

(i+1)mod(3)

(i-1)mod(3)

CGAL - user manual - 2D triangulation

CGAL的triangulations類提供了高抽象級别的集合操作,如triangulation中點的位置,點的插入,删除,替換。他們是建構在triangulation資料結構之上。這個資料結構可以看成是面和頂點的容器。這個資料結構同時關注到此資料結構還照顧到三角剖分的所有組合方面。

幾何和連接配接關系的分離展現了有兩個模闆參數:

  • 第一個用來表示geometric traits,提供幾何圖元(點,線段,三角形)。和這些對象的元素操作(預測,建構);
  • 第二個參數用來表示triangulation data structure類,這個類用來表示三角化的面和頂點,以及handles,iterators,circulators,用于通路頂點和面;CGAL提供了

    Triangulation_data_structure_2<Vb,Fb>

    作為預設的triangulation data structure model。這個類有兩個模闆參數,用來表示頂點類和面類。

最上層負責根據不同的三角化實作(basic,delaunay,regular,constrained or constrained delaunay)以不同的方式實作三角化的幾何建構。每一種triangulation和不同的類相關聯。下圖給出了不同類之間的關系。

CGAL - user manual - 2D triangulation

Triangulation_2<Traits,Tds>

作為其他2D triangulation類的基類,提供了進行triangulation需要的使用者接口。

該類提供了頂點,面通路的handles,iterators,以及circulators。提供了測試任何要素的無限性的方法,還提供了在給定其頂點的情況下測試特定要素(邊緣或面)的三角剖分中是否存在的方法。提供了定位一個點和三角剖分的位置關系,如,一個頂點是否和三角形的頂點共點,是否位于邊上,位于面上,凸包的外部。

點周遊的示例代碼如下:

#include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
#include <CGAL/Triangulation_2.h>
#include <iostream>
#include <vector>

namespace 
{
    typedef CGAL::Exact_predicates_inexact_constructions_kernel K;
    typedef CGAL::Triangulation_2<K>      Triangulation;
    typedef Triangulation::Vertex_handle  Vertex_handle;
    typedef Triangulation::Point TPoint;
    typedef Triangulation::Finite_vertices_iterator Finite_vertices_iterator;

    void test_triangulation_loop()
    {
        std::vector<TPoint> points = { TPoint(0,0), TPoint(1,0), TPoint(0,1) };
        Triangulation T;
        T.insert(points.begin(), points.end());
        std::cout << "Triangulation_2::Finite_vertices_iterator is like a  Triangulation_2::Vertex_handle\n";

        for (Finite_vertices_iterator it = T.finite_vertices_begin();
            it != T.finite_vertices_end();
            ++it) {
            std::cout << it->point() << std::endl;
        }
    }
}

int main()
{
    test_triangulation_loop();
}

           

上面的insert函數内部主要調用了:

Vertex_handle CGAL::Triangulation_2< Traits, Tds >::insert 	( 	const Point &  	p,
		Face_handle  	f = Face_handle() 
	) 	
           

該函數的解釋如下:

将點p插入到triangulation中,傳回相關的vertex。

如果點p和目前已經存在的重複,那麼直接傳回這個vertex,并且triangulation不會有任何變化。

如果p位于邊上,那麼相關的face,會被分割成為兩份;

如果p嚴格的位于face的内部,那麼這個face會被劃分為三份;

如果p嚴格的位于凸包外部,p将和所有可見的頂點相連,形成新的triangulation。

最後,如果p在仿射外殼之外(在退化的一維或0維三角剖分的情況下),則p被連結到所有其他頂點,進而形成三角剖分,其維數增加一。 最後一個參數f訓示從何處開始的底層定位算法。???這個不了解

定位基于stochastic walk(Olivier Devillers, Sylvain Pion, and Monique Teillaud. Walking in a triangulation. Internat. J. Found. Comput. Sci., 13:181–199, 2002.)實作的。從可選參數給出的面的頂點開始周遊,或者在沒有給出可選參數的情況下,從三角剖分的任意頂點開始。對于Delaunay Triangulations而言,最壞的情況,耗時\(O(n)\),當頂點分布比較均勻的時候,平均耗時隻需要\(O(\sqrt{n})\)。在

Triangulation_hierarchy_2<Traits,Tds>

中實作了更加高效的定位算法。

插入一個點,需要找到這個點對應的face,然後将face進行劃分成新的faces。如果該點落在凸包之外,則通過翻轉來恢複三角剖分。不了解??。和定位不同,插入隻需要花費\(O(1)\)的時間。This bound is only an amortized bound for points located outside the convex hull.不了解???

頂點删除是删除所有相鄰的三角形,并對hole進行重新三角化。删除操作的時間最多和\(d^2\)成正比,其中d是被删除頂點的度。

頂點移位:首先,判斷移位之後,三角形之間的嵌入關系是否保持(這個應該是說移位後,三角形之間不存在交疊區域,對于三角形内點的移位,應該是不能移出該三角形的區域);如果保持,那麼直接将vertex替換成新的位置;否則,在新的位置執行插入操作,并移除舊的位置。

整個三角化的過程是一個動态的過程,随着點的變換,動态調整。

Geometric Traits類用來提供幾何對象(點,線段,三角形)和這些對象相關的預測,用于建構triangulation。相關的預測有:

  • 兩個點的x或y坐标的比較;
  • 定向測試,計算三個給定點的順序類型;

概念

TriangulationTraits_2

描述了triangulation對geometric traits的需求。

以下程式使用預設的kernel Exact_predicate_inexact_constructions_kernel建立2D點的三角剖分。

#include <fstream>
#include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
#include <CGAL/Triangulation_2.h>
namespace
{
    typedef ::CGAL::Triangulation_2<Kernel>         Triangulation;
    typedef Triangulation::Vertex_circulator Vertex_circulator;
    typedef Triangulation::Point             Point;
    int test_basic_triangulation() {
        std::ifstream in(kExamplePath + "Triangulation_2/data/triangulation_prog1.cin");
        std::istream_iterator<Point> begin(in);
        std::istream_iterator<Point> end;
        Triangulation t;
        t.insert(begin, end);
        Vertex_circulator vc = t.incident_vertices(t.infinite_vertex()),
            done(vc);
        if (vc != 0) {
            do {
                std::cout << vc->point() << std::endl;
            } while (++vc != done);
        }

        return 0;
    }
}
           

t.incident_vertices(t.infinite_vertex()),

用于擷取最外圍的凸多邊形的邊界點。如果需要擷取所有的三角形,如下:

#include <fstream>
#include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
#include <CGAL/Triangulation_2.h>
namespace
{
    typedef ::CGAL::Triangulation_2<Kernel>         Triangulation;
    typedef Triangulation::Vertex_circulator Vertex_circulator;
    typedef Triangulation::Point             Point;
    int test_basic_triangulation() {
        std::ifstream in("D:/Triangulation_2/data/triangulation_prog1.cin");
        std::istream_iterator<Point> begin(in);
        std::istream_iterator<Point> end;
        Triangulation t;
        t.insert(begin, end);

        for (auto faceHandle : t.finite_face_handles())
        {
            std::cout << faceHandle->vertex(0)->point() << ", "
                << faceHandle->vertex(1)->point() << ", "
                << faceHandle->vertex(2)->point() << std::endl;
        }

        return 0;
    }
}
           

輸入為:

1 0
3 2
4 5 
9 8 
7 4 
5 2 
6 3
10 1
           

輸出結果為:

1 0, 3 2, 4 5
7 4, 10 1, 9 8
3 2, 9 8, 4 5
3 2, 7 4, 9 8
3 2, 1 0, 7 4
7 4, 6 3, 10 1
6 3, 1 0, 5 2
10 1, 6 3, 5 2
7 4, 1 0, 6 3
5 2, 1 0, 10 1
           

形成的三角形為:

CGAL - user manual - 2D triangulation

基本原理了解:https://www.cnblogs.com/RenLiQQ/archive/2008/02/06/1065399.html

https://loopvoid.github.io/2019/11/28/Delaunay三角剖分/

Delaunay_triangulation_2<Traits,Tds>

類用于平面上點集的Delaunay三角化實作。Delaunay三角剖分實作以下空圓屬性(也稱為Delaunay屬性):Delaunay三角剖分中任意三角形的外接圓内不包括其他結點。對于沒有四個共圓點子集的點集,Delaunay三角剖分是唯一的,它是點集的Voronoi圖的對偶。該類繼承自

Triangulation_2<Traits,Tds>

。他對應的traits還提供了額外的類型,用于表示voronoi圖。

該類重寫了點insert,move,remove函數。它還具有成員函數(Vertex_handle Delaunay_triangulation_2 :: nearest_vertex(const Point&p))來回答最近的鄰居查詢,并具有成員函數來構造對偶Voronoi圖的元素(頂點和邊)。

DelaunayTriangulationTraits_2

概念對應的模型。特殊的是,這個概念,提供了side_of_oriented_circle用來描述點s,相對于經過p,q,r點的外接圓的位置關系。

首先使用基本三角剖分的插入成員函數,然後執行一系列翻轉操作以恢複Delaunay屬性,進而在Delaunay三角剖分中插入新點。flip的數量為\(O(d)\),其中d為新插入的點的度。由于點随機均勻分布,一旦插入的點的位置确定後,插入的平均時間為\(O(1)\)。

移除在三角剖分中調用移除,然後以滿足Delaunay準則的方式重新三角剖分所建立的hole。删除度為d的頂點需要的時間為\(O(d^2)\)。當待删除的點的度小于等于7的時候,可以通過特殊的方式減少時間(Olivier Devillers. Vertex removal in two dimensional Delaunay triangulation: Asymptotic complexity is pointless. Research Report 7104, INRIA, 2009.)。

頂點移位:首先,判斷移位之後,三角形之間的嵌入關系是否保持(這個應該是說移位後,三角形之間不存在交疊區域,對于三角形内點的移位,應該是不能移出該三角形的區域);如果保持,那麼直接将vertex替換成新的位置,恢複Delaunay特性,需要\(O(d)\);否則,在新的位置執行插入操作,并移除舊的位置,此時最壞的複雜度為\(O(n)\),但對于在機關正方形中均勻分布的頂點隻需要\(O(1 + \delta\sqrt{n})\), 其中\(\delta\)為新的位置和舊的位置之間的歐幾裡得距離。

對于點的位置的查找,一個點的最鄰近的查找最壞需要\(O(n)\)的時間,但是随機均勻分布的頂點隻需要\(O(1)\)。

#include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
#include <CGAL/Delaunay_triangulation_2.h>
#include <fstream>

namespace 
{
    typedef ::CGAL::Exact_predicates_inexact_constructions_kernel K;
    typedef ::CGAL::Delaunay_triangulation_2<K>  DTriangulation;
    typedef DTriangulation::Edge_iterator  Edge_iterator;
    typedef DTriangulation::Point          Point;
    int test_voronoi_diagram()
    {
        std::ifstream in(kExamplePath + "Triangulation_2/data/voronoi.cin");
        std::istream_iterator<Point> begin(in);
        std::istream_iterator<Point> end;
        DTriangulation T;
        T.insert(begin, end);
        int ns = 0;
        int nr = 0;

        Edge_iterator eit = T.edges_begin();
        for (; eit != T.edges_end(); ++eit) {
            CGAL::Object o = T.dual(eit);
            if (auto seg = CGAL::object_cast<K::Segment_2>(&o)) {
                std::cout << "seg: " << seg->point(0) << ", " << seg->point(1) << std::endl;
                ++ns; 
            }
            else if (auto ray = CGAL::object_cast<K::Ray_2>(&o)) { 
                std::cout << "ray: " << ray->point(0) << ", " << ray->point(1) << std::endl;
                ++nr; 
            }
        }
        std::cout << "The Voronoi diagram has " << ns << " finite edges "
            << " and " << nr << " rays" << std::endl;
        return 0;
    }
}

int main()
{
    test_voronoi_diagram();
    return 0;
}
           
ray: 0 1, -4 1
seg: 0 1, 0.5 0.5
ray: 0 1, 2 3
ray: 0.5 0.5, 2.5 0.5
ray: 0.5 0.5, 0.5 -1.5
The Voronoi diagram has 1 finite edges  and 4 rays
           

繪制得到的結果如下:

CGAL - user manual - 2D triangulation

其中藍色的點為輸出,紅色的點會輸出,紅色的線為線段,紫色的線為射線。

将Delaunay triangulation得到的結果限制到矩形框内。

該示例詳細見原文。其中主要的一個函數,如下,用來将dual形成的Voronoi圖輸出到stream中,示例中借助stream,将範圍限制在了矩形内:

template < class Stream>
  Stream& draw_dual(Stream & ps)
  {
    Finite_edges_iterator eit= this->finite_edges_begin();
    for (; eit != this->finite_edges_end(); ++eit) {
      Object o = dual(eit);
      typename Geom_traits::Line_2  l;
      typename Geom_traits::Ray_2   r;
      Segment s;
      if(CGAL::assign(s,o)) ps << s;
      if(CGAL::assign(r,o)) ps << r;
      if(CGAL::assign(l,o)) ps << l;
    }
    return ps;
  }
           

TODO:當需要應用到該概念的時候再進行擴充

限制三角剖分是一組點的三角剖分,必須在其邊緣之間包括一組給定的連接配接點的折線。 折線稱為限制。 相應的邊稱為限制邊。

CGAL - user manual - 2D triangulation

通過上圖了解,就是說,已知一些點,和一些相連的邊,然後進行三角化的過程,就是constrained triangulations。

限制邊的端點當然是三角剖分的頂點。 但是,三角剖分還可以包括其他頂點。 限制三角剖分共有三種形式:

  • 在基本的版本中,Constrained triangulations并不處理相交限制,輸出的限制是不相交的線段(斷點處的相交不算)。受限制的邊緣可能是垂直的或長度為零。
  • 其他兩種情況,是針對限制邊有相交,或交疊,部分交疊的時候。這個時候需要引入額外的點。這個時候有兩種處理方式(TODO:這兩種到底有啥差別):
    • 第一種是,對于條件判斷是準确評估,對于結構計算(如,相交計算)進行近似計算;
    • 第二種,執行準确評估, 和準确計算;

處理的類為:

Constrained_triangulation_2<Traits,Tds,Itag>

,模闆參數

Itag

用來表示如何處理相交限制,可以是:

  • No_intersection_tag

    輸入的限制條件不相交;
  • Exact_predicates_tag

    :準确評估,結構近似計算;
  • Exact_intersections_tag

    :準确評估,結構準确計算;

通常情況下符合概念

TriangulationTraits_2

. 當處理的限制相交的時候,需要滿足概念

ConstrainedTriangulationTraits_2

,該概念增加了兩個線段相交計算。

限制邊存儲在三角化的面中。限制三角剖分的面基礎必須是概念ConstrainedTriangulationFaceBase_2的模型,該模型完善了概念TriangulationFaceBase_2。

關于限制Delaunay三角剖分概念的介紹可以參見:https://blog.csdn.net/alegriabaile/article/details/88816511

受限制的Delaunay三角剖分是在基于限制條件進行三角剖分的時候,盡可能多的保留Delaunay特性。由于受限制的邊緣不一定是Delaunay邊緣,限制Delaunay三角剖分的三角形不一定滿足空圓屬性,但它們滿足較弱的限制空圓屬性。正如上面連結中介紹的,限制的邊,會擋住三角形發現點的視線。那麼一個三角形是限制Delaunay當且僅當,從這個三角形出發往四周看,看不到其他的vertex,如下圖中的箭頭示意尋找點的視線方向,有兩個方向被限制邊擋住了。

CGAL - user manual - 2D triangulation

有三種限制Delaunay三角剖分的方式提供:第一種用于處理最多再斷點處相交的限制邊;另外兩種用來處理限制邊相交的情況。

Constrained_Delaunay_triangulation_2<Traits,Tds,Itag>

用來處理帶限制的Delaunay三角剖分。

#include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
#include <CGAL/Constrained_Delaunay_triangulation_2.h>
#include <cassert>
#include <iostream>

namespace
{
    typedef ::CGAL::Exact_predicates_inexact_constructions_kernel K;
    typedef ::CGAL::Exact_predicates_tag Itag;
    typedef ::CGAL::Constrained_Delaunay_triangulation_2<K, ::CGAL::Default, Itag> CDT;
    typedef CDT::Point Point;
    typedef CDT::Edge Edge;

    void test_cdt()
    {
        CDT cdt;
        int maxNum = 2;
        for (int i = 1; i < maxNum; ++i) {
            cdt.insert_constraint(Point(0, i), Point(maxNum, i));
        }
        for (int j = 1; j < maxNum; ++j) {
            cdt.insert_constraint(Point(j, 0), Point(j, maxNum));
        }

        assert(cdt.is_valid());
        for (auto& e : cdt.finite_edges()) {
            auto seg = cdt.segment(e);
            if (cdt.is_constrained(e)) {
                std::cout << "constrained edge: ";
            }
            else {
                std::cout << "normal edge     : ";
            }
            std::cout << seg.point(0) << ", " << seg.point(1) << std::endl;
        }
    }
}
int main()
{
    test_cdt();
    return 0;
}
           

得到的結果如下:

normal edge     : 2 1, 1 0
normal edge     : 0 1, 1 0
constrained edge: 1 0, 1 1
constrained edge: 1 1, 0 1
constrained edge: 1 1, 1 2
normal edge     : 1 2, 0 1
normal edge     : 1 2, 2 1
constrained edge: 2 1, 1 1
           

圖如下:

CGAL - user manual - 2D triangulation

TODO

版權說明

作者: grassofsky

出處: http://www.cnblogs.com/grass-and-moon

本文版權歸作者,歡迎轉載,但未經作者同意必須保留此段聲明,且在文章頁面明顯位置給出, 原文連結 如有問題, 可郵件([email protected])咨詢.