前言
最近闲暇之余，无事喜欢研究一些经典的算法，这两天看了Delauny三角网的一些算法，以其最经典的一种为例进行分享，今天又是反对内卷的一天。
什么是Delauny三角网？
首先说三角网，三角网是由一系列连续三角形构成的网状的平面控制图形，用于描述二、三维物体表面起伏特征的网格状图形，如下图所示：
所谓的Delauny（德劳内）三角网，就是由Delauny三角形所组成的三角网，称为Delauny三角网，关于Delauny三角网的定义，可以参照连接：Delaunay三角网是什么? - 知乎 (zhihu.com)
还是那句话，把定义搞懂了，代码都是小事，不想看定义的可以直接跳到代码部分。
Bowyer-Watson
Delauny三角网和Bowyer-Watson的关系
了解了Delauny三角网是什么，那么什么又是Bowyer-Watson呢？
新手容易搞混哈，要想弄明白也简单，Delauny三角网只是一种理想模型，一种理论目标，并不是一种算法，并不指向具体过程，就比如：领导给你交给你一项任务，最终要形成什么样子的效果，达到什么目的，中间的过程是什么就靠自己去实现了，最终完成目标即可。Delauny三角网就是一个目标，要构建Delauny三角网，算法有很多种，这里笔者用的是Bowyer-Watson算法。
核心过程
假定已生成了连接若干个顶点的 Delaunay 三角网格加入一个新的节点，找出所有外接圆包含新加入节点的三角形，并将这些三角形删除，形成一个空腔空腔的节点与新加入的节点连接，形成新的 Delaunay 三角形网格不断循环直到遍历完所有点、
这里笔者就只写一个核心过程，具体关于Bowyer-Watson构建Delauny三角网的过程可以参照博客：Bowyer-Watson算法-CSDN博客
效果演示
网上关于三角网的算法有很多，现成的库也有，但是笔者为啥非得自己去研究代码呢，这可能是闲得慌吧，哈哈，笔者对于任何东西都不只是想停留在表面，学习原理是进步的一种方式，研究算法的原理过程中所获得的心得是直接看书的两倍甚至以上，每次笔者在研究的过程中都会产生新的灵感与思路，只有学习了原理才能在原来的基础之上加入自己的东西，而不是单靠着一些现成功能的调用就觉得事情很简单，调用谁不会呢。
当然了，我也不提倡什么都研究，如果对自己学习和研究方向有用的东西，可以适当的去研究一下，说不定你们也能获得新的东西。
话说远了哈，就此打住。
代码
主代码：
///
///
/// 德劳内三角剖分（效率较低，待优化）
///
public class DelaunyTriangulations
{
///
/// 外部点集
///
private List
///
/// 超级三角
///
public Triangle SuperTriangle { get; private set; }
public DelaunyTriangulations(List
{
this.Coordinates = Coordinates;
this.SuperTriangle = GetSuperTriangle();
}
///
/// 获取超级三角（这里我生成了一个等边三角）
///
///
private Triangle GetSuperTriangle()
{
//先获取点集的最大外包矩形
IExtent Envelope = new IExtent();//先创建一个空矩形对象
Coordinates = Coordinates.AsEnumerable().OrderBy(pt => pt.Y).ToList();//（Y升序排列）
double Ymin = Coordinates[0].Y, Ymax = Coordinates[Coordinates.Count - 1].Y;//获取最小Y值和最大Y值
Coordinates = Coordinates.AsEnumerable().OrderBy(pt => pt.X).ToList();//(X升序排列）
double Xmin = Coordinates[0].X, Xmax = Coordinates[Coordinates.Count - 1].X;//获取最小X值和最大X值
Envelope.MinPoint = new IPoint(Xmin, Ymin, 0);//赋值最小点（左下角）
Envelope.MaxPoint = new IPoint(Xmax, Ymax, 0);//赋值最大点（右上角）
Envelope.Scale(1.5, 1.5);//防止后续有点集在外部大三角边上的情况
return ISpatialMath.GetEquilateralTriangle(Envelope);//返回矩形的外接大三角（等边三角形）
}
///
/// 构建三角网
///
///
public List
{
if (Coordinates.Count < 3) return new List
List
List
for (int i = 0; i < Coordinates.Count; i++)
{
IPoint curPoint = Coordinates[i];//当前顶点
if (i == 0)
{
ReSplitTriangles.Add(SuperTriangle);//将外部最大三角放入分割列表
ReSplitTriangles = SplitTriangleWith(ReSplitTriangles, curPoint);
DelaunyTriangles.AddRange(ReSplitTriangles);//将分割后的三角形存放至德劳内三角列表
continue;//第一个点先进行分割得到基础的德劳内三角网
}
//判断当前点位于哪些德劳内三角形的外接圆内
//删除所在外接圆的三角形，并用当前点重新分割删除的三角形
//这一步可以借助空间索引筛选
for (int c = 0; c < DelaunyTriangles.Count; c++)
{
Triangle triangle = DelaunyTriangles[c];
if (triangle.Area <= 0)//移除空三角
{
DelaunyTriangles.RemoveAt(c);
c--;
continue;
}
Circle circle = triangle.Tag as Circle;
//当前点到三角形的外接圆圆心距离如果小于等于圆的半径，则认为点在圆内部
if (circle.Center.DistanceTo(curPoint, true) <= circle.Radius)
{
ReSplitTriangles.Add(triangle);//将当前三角纳入重新分割的列表
DelaunyTriangles.RemoveAt(c);//将当前三角从德劳内三角列表移除
c--;
}
}
ReSplitTriangles = SplitTriangleWith(ReSplitTriangles, curPoint);
DelaunyTriangles.AddRange(ReSplitTriangles);//将分割后的三角形存放至德劳内三角列表
if (ReSplitTriangles.Count != 0) ReSplitTriangles.Clear();//清理临时列表
}
FilterDelauny(DelaunyTriangles);//过滤三角网（将包含外部大三角顶点的德劳内三角删除）
return DelaunyTriangles;
}
///
/// 通过点将三角形拆分成多个三角形
///
/// 拆分点
///
private List
{
List
if (Triangles.Count == 0) return triangles;//如果没有待分割的三角形，则返回空值
List
for (int i = 0; i < Triangles.Count; i++)
{
IPoint top = Triangles[i].Top;
IPoint left = Triangles[i].LeftBottom;
IPoint right = Triangles[i].RightBottom;
//如果点集内不存在当前点坐标则加入
var q1 = from pt in Vertex where pt.Equals(top) select pt;
if (q1.Count() == 0)
{
Vertex.Add(top);
}
var q2 = from pt in Vertex where pt.Equals(left) select pt;
if (q2.Count() == 0)
{
Vertex.Add(left);
}
var q3 = from pt in Vertex where pt.Equals(right) select pt;
if (q3.Count() == 0)
{
Vertex.Add(right);
}
}
//将点集按照与拆分点连线的水平方向夹角进行排序（逆时针）
Vertex = Vertex.AsEnumerable().OrderBy(pt => Point.VectorTo(pt).Degree).ToList();
for (int i = 0; i < Vertex.Count; i++)
{
int next = i == Vertex.Count - 1 ? 0 : i + 1;
IPoint left = Vertex[i];
IPoint right = Vertex[next];
Triangle tri = new Triangle(Point, left, right);//创建三角形
CircumCircle circumCircle = new CircumCircle();//调用外接圆工具
tri.Tag = circumCircle.BuildFromTriangle(tri);//绑定外接圆形至三角
triangles.Add(tri);//将分割后的三角形添加到结果
}
Vertex.Clear();//清理临时列表
return triangles;
}
///
/// 过滤德劳内三角（将包含外部大三角顶点的德劳内三角删除）
///
///
private void FilterDelauny(List
{
List
outPoints.Add(SuperTriangle.Top);
outPoints.Add(SuperTriangle.LeftBottom);
outPoints.Add(SuperTriangle.RightBottom);
for (int i = 0; i < Delauny.Count; i++)
{
IPoint top = Delauny[i].Top;
IPoint left = Delauny[i].LeftBottom;
IPoint right = Delauny[i].RightBottom;
var q = from pt in outPoints where pt.Equals(top) || pt.Equals(left) || pt.Equals(right) select pt;
if (q.Count() > 0)
{
Delauny.RemoveAt(i);
i--;
}
}
outPoints.Clear();
}
}
辅助类
Vector
///
/// 向量
///
public class Vector
{
///
/// 向量X方向增量
///
public double X { get; set; }
///
/// 向量Y方向的增量
///
public double Y { get; set; }
///
/// 向量Z方向的增量
///
public double Z { get; set; }
///
/// 向量的模长
///
public double Norm { get { return GetNormOfVector(); } }
public Vector(double X, double Y, double Z = 0)
{
this.X = X;
this.Y = Y;
this.Z = Z;
}
///
/// 计算与其他向量的夹角
///
///
///
public double DegreeTo(Vector OtherVector)
{
// 计算向量长度
double vector1Length = Math.Sqrt(this.X * this.X + this.Y * this.Y + this.Z * this.Z);
double vector2Length = Math.Sqrt(OtherVector.X * OtherVector.X + OtherVector.Y * OtherVector.Y + OtherVector.Z * OtherVector.Z);
// 计算向量点积
double dotProduct = this.X * OtherVector.X + this.Y * OtherVector.Y + this.Z * OtherVector.Z;
// 计算夹角（单位为弧度）
double angle = Math.Acos(dotProduct / (vector1Length * vector2Length));
return angle.RadiansToDegree();
}
///
/// 求向量的模长
///
///
private double GetNormOfVector()
{
return Math.Sqrt(this.X * this.X + this.Y * this.Y + this.Z * this.Z);
}
}
IMath
public static class IMath
{
///
/// 弧度转角度
///
/// 弧度值
///
public static double RadiansToDegree(this double Radians)
{
return (180.00 / Math.PI) * Radians;
}
///
/// 角度转弧度
///
/// 角度值
///
public static double DegreeToRadians(this double Degree)
{
return Degree * Math.PI / 180.00;
}
}
Triangle
public class Triangle
{
///
/// 顶点
///
public IPoint Top { get { return _Top; } set { _Top = value; } }
private IPoint _Top;
///
/// 左底点
///
public IPoint LeftBottom { get { return _LeftBottom; } set { _LeftBottom = value; } }
private IPoint _LeftBottom;
///
/// 右底点
///
public IPoint RightBottom { get { return _RightBottom; } set { _RightBottom = value; } }
private IPoint _RightBottom;
///
/// 左边
///
public ILine LeftSide { get { return new ILine(Top, LeftBottom); } }
///
/// 右边
///
public ILine RightSide { get { return new ILine(Top, RightBottom); } }
///
/// 底边
///
public ILine Bottom { get { return new ILine(LeftBottom, RightBottom); ; } }
///
/// 周长
///
public double Length { get { return LeftSide.Length + RightSide.Length + Bottom.Length; } }
///
/// 面积
///
public double Area { get { return Math.Sqrt((this.Length / 2) * ((this.Length / 2 - RightSide.Length) * (this.Length / 2 - LeftSide.Length) * (this.Length / 2 - Bottom.Length))); } }
///
/// 判断是否有效三角
///
public bool IsValid { get { return !(this.Area == double.NaN || this.Area <= 0); } }
///
/// 外部挂接信息（可自定义对象）
///
public object Tag { get; set; }
///
/// 实例化
///
///
///
///
public Triangle(IPoint Top,IPoint LeftBottom,IPoint RightBottom)
{
this.Top = Top;
this.RightBottom = RightBottom;
this.LeftBottom = LeftBottom;
}
}
IPoint
///
/// 点对象
///
public class IPoint
{
public double X
{
get
{
//无精度限定
if (LeaderConfig.GeoDights==-1)
{
return _x;
}
else
{
return Math.Round(_x, LeaderConfig.GeoDights);
}
}
set
{
_x = value;
}
}
private double _x;
public double Y
{
get
{
if (LeaderConfig.GeoDights==-1)
{
return _y;
}
else
{
return Math.Round(_y, LeaderConfig.GeoDights);
}
}
set
{
_y = value;
}
}
private double _y;
public double Z
{
get
{
if (LeaderConfig.GeoDights==-1)
{
return _z;
}
else
{
return Math.Round(_z, LeaderConfig.GeoDights);
}
}
set
{
_z = value;
}
}
private double _z;
///
/// 实例化
///
///
///
///
public IPoint(double X, double Y, double Z)
{
this.X = X;
this.Y = Y;
this.Z = Z;
}
///
/// 两点距离
///
/// 第二点
/// 是否忽略Z值
///
public double DistanceTo(IPoint SecondPoint,bool IgnoreZvalue=false)
{
if (IgnoreZvalue)
{
return Math.Sqrt(Math.Pow(this.X - SecondPoint.X, 2) + Math.Pow(this.Y - SecondPoint.Y, 2));
}
else
{
return this.VectorTo(SecondPoint).Norm;
}
}
///
/// 获取到另一点的向量
///
///
public Vector VectorTo(IPoint SecondPoint)
{
return new Vector(SecondPoint.X - this.X, SecondPoint.Y - this.Y, SecondPoint.Z - this.Z);
}
///
/// 判断两点是否相等
///
///
///
public bool Equals(IPoint iPoint)
{
return this.X == iPoint.X && this.Y == iPoint.Y && this.Z == iPoint.Z;
}
}
ILine
///
/// 线段
///
public class ILine
{
///
/// 起点
///
public IPoint StartPoint { get { return _startPoint; } set { _startPoint = value; } }
private IPoint _startPoint;
///
/// 终点
///
public IPoint EndPoint { get { return _endPoint; } set { _endPoint = value; } }
private IPoint _endPoint;
public double Length { get { return StartPoint.DistanceTo(EndPoint); } }
///
/// 实例化
///
///
///
public ILine(IPoint start, IPoint end)
{
_startPoint = start;
_endPoint = end;
}
}
CicumCircle
public class CircumCircle
{
///
/// 生成三角形的外接圆
///
///
///
public Circle BuildFromTriangle(Triangle triangle)
{
if (!triangle.IsValid) return null;
IPoint Center = GetCircleCenterWith(triangle.Top, triangle.LeftBottom, triangle.RightBottom);
return new Circle(Center, triangle.Top.DistanceTo(Center,true));
}
///
/// 已知圆上3点获取圆心
///
///
///
///
///
private IPoint GetCircleCenterWith(IPoint P1, IPoint P2, IPoint P3)
{
double a13 = P1.X - P3.X;
double a13_ = P1.X + P3.X;
double b13 = P1.Y - P3.Y;
double b13_ = P1.Y + P3.Y;
double a12 = P1.X - P2.X;
double a12_ = P1.X + P2.X;
double b12 = P1.Y - P2.Y;
double b12_ = P1.Y + P2.Y;
double a12b12_2 = a12 * a12_ + b12 * b12_;
double a13b13_2 = a13 * a13_ + b13 * b13_;
double a13b12 = 2 * a13 * b12;
double a12b13 = 2 * a12 * b13;
if (a12b13 - a13b12 == 0) return new IPoint((P2.X + P1.X) / 2.0, (P2.Y + P1.Y) / 2.0, 0);
double af = a12b13 - a13b12;
double bf = a13b12 - a12b13;
double az = b13 * a12b12_2 - b12 * a13b13_2;
double bz = a13 * a12b12_2 - a12 * a13b13_2;
double a = az / af;
double b = bz / bf;
return new IPoint(a, b, 0);
}
}
Circle
public class Circle:IEntity
{
///
/// 圆心
///
public IPoint Center { get; set; }
///
/// 半径
///
public double Radius { get; set; }
///
/// 面积
///
public double Area { get { return Math.PI * Math.Pow(Radius, 2); } }
///
/// 周长
///
public double Length { get { return Math.PI * Radius * 2; } }
public bool IsValid { get { return Radius > 0; } }
///
/// 点径构圆
///
///
///
public Circle(IPoint Center,double Radius)
{
this.Center = Center;
this.Radius = Radius;
}
///
/// 两点构圆
///
///
///
public Circle(IPoint Center, IPoint PointOnCircle)
{
this.Center = Center;
this.Radius = Center.DistanceTo(PointOnCircle);
}
}
结语
如果要自己改代码，一定得先看懂BowyerWatson算法的过程，重要的事说三遍，不懂的可以私聊我。
参考链接
您阅读本篇文章共花了： 大家都在找：
c#：c#菜鸟教程
算法：算法工程师学什么专业
windows：windows官网
发表评论