Delaunay三角网,写了用半天,调试BUG用了2天……醉了。
基本思路比较简单,但效率并不是很快。
1. 先生成一个凸包;
2. 只考虑凸包上的点,将凸包环切,生成一个三角网,暂时不考虑Delaunay三角网各种规则。将生成的三角形放进三角形集合 Triangles 中;
3.将其它非凸包上的点全都插入。每插入一个点 ptA,都要判断被插入的点存在于 Triangles 集合的哪个三角形(trianA)之内,并将 ptA 与此三角形的三个点进行连接,删除 trianA,并将新生成的三角形加入到集合 Triangles 中。初始三角网生成结束;
3.1. 若 ptA 处在三角形DEF的DE边上,那么只连接点F与ptA ;
4.遍历三角形集合 Triangles(曾考虑用邻接矩阵,但是使用矩阵的复杂度反而会更高),每遍历到一个三角形DEF,都要遍历DEF的三条边DE,EF, FD,并分别寻找另外一个三角形,此三角形与DEF存在一个公共边;
5.使用LOP(Local Optimization Procedure: 局部优化处理)处理步骤4返回的两个三角形 与DEG。
5.1. 做三角形DEF的外接圆cirDEF,如果点G在圆cirDEF之内,则从Triangles 集合中删除三角形DEF与DEG,并加入两个新三角形:FGD与FGE。实际上,是将两个三角形组成的四边形的对角线对调了。
6.使用新生成的两个三角形FGD与FGE执行步骤5,递归。递归出口为,对三角形进行LOP处理时,没有出现5.1的情况。
Note:
1.暂时存在的问题:
a. 当点数过多时,生成过程中会出现“栈溢出”情况。测试时,2000个点之内是成功的。此情况是由于三角形过多导致的递归过深,需要重新组织算法结构才能避免。
b. 生成时间太长。可采取的办法1是完全采用不同的生成算法,比如分块等等。其次是由于寻找三角形的复杂度过高,达到了O(n),利用空间换时间的办法能接近O(1),只不过多一个成员来存储三角形的拓扑结构,这个比较简单但优化效果有限。
2. 一旦遇到了一些算法上的问题,那么搬出数学知识往往非常有效(尽管我的数学水平一望见底)。使用矩阵的方式能够让算法的实现更加清晰,但复杂度有可能会升高。
*下面使用了https://www.cnblogs.com/rkexy/p/9768475.html其中的代码。*
File: Delaunay.h
1 #pragma once 2 3 #include "ConvexHull.h" 4 #include "Triangle.h" 5 #include "Circle.cpp" 6 7 class TrianIndex final 8 { 9 public: 10 TrianIndex() 11 { 12 _isVaild = false; 13 } 14 15 ~TrianIndex() 16 { 17 _isVaild = false; 18 } 19 20 TrianIndex(unsigned int iA, unsigned int iB, unsigned int iC) 21 { 22 Init(iA, iB, iC); 23 } 24 25 TrianIndex(std::array<unsigned int, 3> pts) : 26 TrianIndex(pts[0], pts[1], pts[2]) 27 { 28 } 29 30 TrianIndex(const TrianIndex& other) 31 { 32 this->_isVaild = other._isVaild; 33 this->_ptIdx = other._ptIdx; 34 } 35 36 unsigned int& Get(int i) 37 { 38 if (i < 0 || i > 2 || !_isVaild) 39 { 40 ErrorThrow("Error Triangle Point Index[0, 2] Or Invaild Triangle: " + std::to_string(i)); 41 return _ptIdx[0]; 42 } 43 44 return _ptIdx[i]; 45 } 46 47 unsigned int& operator[](int i) 48 { 49 return Get(i); 50 } 51 52 const unsigned int Get(int i) const 53 { 54 TrianIndex* pThis = const_cast<TrianIndex*>(this); 55 return pThis->Get(i); 56 } 57 58 const unsigned int operator[](int i) const 59 { 60 TrianIndex* pThis = const_cast<TrianIndex*>(this); 61 return (*pThis)[i]; 62 } 63 64 TrianIndex& operator=(const TrianIndex& other) 65 { 66 this->_isVaild = other._isVaild; 67 this->_ptIdx = other._ptIdx; 68 return *this; 69 } 70 71 bool IsVaild() const 72 { 73 return _isVaild; 74 } 75 76 void SetVaild(bool isVaild) 77 { 78 _isVaild = isVaild; 79 } 80 81 void Init(unsigned int iA, unsigned int iB, unsigned int iC) 82 { 83 _ptIdx[0] = iA; 84 _ptIdx[1] = iB; 85 _ptIdx[2] = iC; 86 _isVaild = true; 87 } 88 89 private: 90 bool _isVaild; 91 std::array<unsigned int, 3> _ptIdx; 92 }; 93 94 //CHECK IT: It Has Been Abandoned. 95 // 96 class TrianMatrix final 97 { 98 public: 99 enum IsConnected : bool 100 { 101 NotConnected = false, 102 Connected = true 103 }; 104 105 TrianMatrix(): 106 TrianMatrix(0) 107 { 108 } 109 110 TrianMatrix(size_t ptCount) 111 { 112 _sign.clear(); 113 _sign.resize(ptCount); 114 for (unsigned int x = 0; x < ptCount; ++x) 115 { 116 _sign[x].resize(ptCount); 117 for (unsigned int y = 0; y < ptCount; ++y) 118 { 119 _sign[x][y] = IsConnected::NotConnected; 120 } 121 } 122 } 123 124 ~TrianMatrix() 125 { 126 _sign.clear(); 127 } 128 129 public: 130 size_t PointCount() const 131 { 132 return _sign.size(); 133 } 134 135 bool IsEmpty() const 136 { 137 return _sign.empty(); 138 } 139 140 // 141 bool GetTrianglesByPoint(unsigned int ptIdx, 142 std: