Delaunay三角剖分

Delaunay三角网是俄国数学家B.Delaunay于1934年发现的。关于Delaunay三角网构建的研究有许多，但由于本课题具有数据量大的特征，不宜直接沿用已有构建方法，笔者针对本课题数据特征，研究获得了适应本课题，速度较快的构建方法。Delaunay三角网有一个特性，每个三角网形成的外接圆都不包含其他参考点。利用这一个性质，我们可以直接构成Delaunay三角网：
一、建立第一个三角形
1、判断用来建立TIN的总脚点数，小于3则报错退出。
2、第一点的选择：
链表的第一节点，命名为Pt1;
3、第二点的选择：
A.非Pt1点； B.距Pt1最近
命名为Pt2
4、第三点的选择
A.非Pt1，Pt2点；
B.与Pt1，Pt2点组成的三角形的外接圆内无其他节点；
C.与Pt1，Pt2组成的三角形中的角Pt1Pt3Pt2最大。
命名为Pt3
5、生成三边，加入边表
6、生成第一个三角形，组建三角形表
二、扩展TIN
1、从边表头取一边，要求：该边flag标志为假（只在一个三角形中）
2、从点链表中搜索一点，要求：
A、与边中的Pixel3在边的异侧；
B、该点与边组成的三角形的外接圆内无其他点
C、满足上面两条件的点中角Pt1Pt3Pt2最大的点为Pt3。
3、判断新生成的边，如果边表中没有，则加入边表尾，设定标志flag为假，如果有，则设定该边flag为真。
4、将生成的三角形假如三角形表
5、设定选中的边的标志flag为真
6、转至1，直至边表边的标志flag全部为真。


数据结构：
struct Pixel //脚点数据
{
double x,y,z,g;
bool flag;
};
struct List //数据链表
{
Pixel *pixel;
List *next;
};
struct Line //三角形边
{
Pixel *pixel1; //三角形边一端点
Pixel *pixel2; //三角形边另一端点
Pixel *pixel3; //三角形边所对顶点
bool flag;
};
struct Linelist //三角形边表
{
Line *line;
Linelist *next;
};
struct Triangle //三角形表
{
Line *line1;
Line *line2;
Line *line3;
Triangle *next;
};



以下是我程序中关于建网的部分：
// DelaunayTIN.cpp: implementation of the CDelaunayTIN class.
//
//

#include "stdafx.h"
#include "Reconstruction.h"
#include "DelaunayTIN.h"

#include "MyMath.h"
#include <math.h>
#include "ListControl.h"

#ifdef _DEBUG
#undef THIS_FILE
static char THIS_FILE[]=__FILE__;
#define new DEBUG_NEW
#endif

//
// Construction/Destruction
//

CDelaunayTIN::CDelaunayTIN()
{

}

CDelaunayTIN::~CDelaunayTIN()
{

}

/
//函数名： CreateDelaunayTIN
//编写者： Polaris
//参考资料：
//功能： 用给定的数据链表数据，组建Delaunay不规则三角网
//输入参数：数据链表list;区域范围(XMin,YMin),(XMax,YMax)
//输出参数：不规则三角网首三角形地址
//备注：
/
struct Triangle * CDelaunayTIN::CreateDelaunayTIN(List *list)
{
//组建第一个三角形
CMyMath MyMath;
CListControl ListControl;
struct List *node;
struct Pixel *pt1,*pt2,*pt3;
bool flag;
struct Triangle *TIN;
pt1=list->pixel;
pt2=list->next->pixel;
node=list->next->next;
while(node!=NULL)
{
if(MyMath.Calculate2PtDistanceIn3D
(pt1->x,pt1->y,pt1->z,node->pixel->x,node->pixel->y,node->pixel->z)
<MyMath.Calculate2PtDistanceIn3D
(pt1->x,pt1->y,pt1->z,pt2->x,pt2->y,pt2->z))
{
pt2=node->pixel;
}
node=node->next;
}
node=list->next;
pt3=NULL;
while(node!=NULL)
{
if(node->pixel==pt1 || node->pixel==pt2)
{
node=node->next;
continue;
}
if(pt3==NULL)
{
pt3=node->pixel;
}
else
{
if((pow(MyMath.Calculate2PtDistanceIn2D(pt1->x,pt1->y,node->pixel->x,node->pixel->y),2)+pow(MyMath.Calculate2PtDistanceIn2D(pt2->x,pt2->y,node->pixel->x,node->pixel->y),2)-pow(MyMath.Calculate2PtDistanceIn2D(pt1->x,pt1->y,pt2->x,pt2->y),2))/(2*MyMath.Calculate2PtDistanceIn2D(pt1->x,pt1->y,node->pixel->x,node->pixel->y)*MyMath.Calculate2PtDistanceIn2D(pt2->x,pt2->y,node->pixel->x,node->pixel->y))
<(pow(MyMath.Calculate2PtDistanceIn2D(pt1->x,pt1->y,pt3->x,pt3->y),2)+pow(MyMath.Calculate2PtDistanceIn2D(pt2->x,pt2->y,pt3->x,pt3->y),2)-pow(MyMath.Calculate2PtDistanceIn2D(pt1->x,pt1->y,pt2->x,pt2->y),2))/(2*MyMath.Calculate2PtDistanceIn2D(pt1->x,pt1->y,pt3->x,pt3->y)*MyMath.Calculate2PtDistanceIn2D(pt2->x,pt2->y,pt3->x,pt3->y)))
{
pt3=node->pixel;
}
}
node=node->next;
}
//LineList
Linelist *linehead,*linenode,*linelast;
Line *ln1,*ln2,*ln3;
linenode=new Linelist;
linenode->line=new Line;
linenode->line->pixel1=pt1;
linenode->line->pixel2=pt2;
linenode->line->pixel3=pt3;
linenode->line->flag=false;
linenode->next=NULL;
linehead=linelast=linenode;
ln1=linenode->line;
linenode=new Linelist;
linenode->line=new Line;
linenode->line->pixel1=pt2;
linenode->line->pixel2=pt3;
linenode->line->pixel3=pt1;
linenode->line->flag=false;
linenode->next=NULL;
linelast->next=linenode;
linelast=linenode;
ln2=linenode->line;
linenode=new Linelist;
linenode->line=new Line;
linenode->line->pixel1=pt3;
linenode->line->pixel2=pt1;
linenode->line->pixel3=pt2;
linenode->line->flag=false;
linenode->next=NULL;
linelast->next=linenode;
linelast=linenode;
ln3=linenode->line;
//first Triangle
Triangle *tglhead,*tglnode,*tgllast;
tglnode=new Triangle;
tglnode->line1=ln1;
tglnode->line2=ln2;
tglnode->line3=ln3;
tglnode->next=NULL;
tglhead=tgllast=tglnode;

/*//test
TIN=tglhead;
TINCreateOK=true;
PostMessage(WM_PAINT);
AfxMessageBox("ok");
//test end*/

//expend tin;
Linelist *linetmp,*linetemp;
List *pixeltmp;
double x1,y1,x2,y2,x3,y3;
linetmp=linehead;
while(linetmp!=NULL)
{
if(linetmp->line->flag==true)
{
linetmp=linetmp->next;
continue;
}
ln1=linetmp->line;
pt1=linetmp->line->pixel1;
pt2=linetmp->line->pixel2;
x1=linetmp->line->pixel1->x;
y1=linetmp->line->pixel1->y;
x2=linetmp->line->pixel2->x;
y2=linetmp->line->pixel2->y;
x3=linetmp->line->pixel3->x;
y3=linetmp->line->pixel3->y;

pixeltmp=list;
pt3=NULL;
while(pixeltmp!=NULL)
{
if(pixeltmp->pixel==pt1 || pixeltmp->pixel==pt2)
{
pixeltmp=pixeltmp->next;
continue;
}
if(((y2-y1)*pixeltmp->pixel->x+(x1-x2)*pixeltmp->pixel->y+(x2*y1-x1*y2))
*((y2-y1)*x3+(x1-x2)*y3+(x2*y1-x1*y2))>=0)
{
pixeltmp=pixeltmp->next;
continue;
}
if(pt3==NULL)pt3=pixeltmp->pixel;
else
{
if((pow(MyMath.Calculate2PtDistanceIn2D(pt1->x,pt1->y,pixeltmp->pixel->x,pixeltmp->pixel->y),2)+pow(MyMath.Calculate2PtDistanceIn2D(pt2->x,pt2->y,pixeltmp->pixel->x,pixeltmp->pixel->y),2)-pow(MyMath.Calculate2PtDistanceIn2D(pt1->x,pt1->y,pt2->x,pt2->y),2))/(2*MyMath.Calculate2PtDistanceIn2D(pt1->x,pt1->y,pixeltmp->pixel->x,pixeltmp->pixel->y)*MyMath.Calculate2PtDistanceIn2D(pt2->x,pt2->y,pixeltmp->pixel->x,pixeltmp->pixel->y))
<(pow(MyMath.Calculate2PtDistanceIn2D(pt1->x,pt1->y,pt3->x,pt3->y),2)+pow(MyMath.Calculate2PtDistanceIn2D(pt2->x,pt2->y,pt3->x,pt3->y),2)-pow(MyMath.Calculate2PtDistanceIn2D(pt1->x,pt1->y,pt2->x,pt2->y),2))/(2*MyMath.Calculate2PtDistanceIn2D(pt1->x,pt1->y,pt3->x,pt3->y)*MyMath.Calculate2PtDistanceIn2D(pt2->x,pt2->y,pt3->x,pt3->y)))
{
pt3=pixeltmp->pixel;
}
}
pixeltmp=pixeltmp->next;
}
if(pt3!=NULL)
{
linetemp=linehead;
flag=false;
while(linetemp!=NULL)
{
if((pt1==linetemp->line->pixel1 && pt3==linetemp->line->pixel2)
|| (pt3==linetemp->line->pixel1 && pt1==linetemp->line->pixel2))
{
linetemp->line->flag=true;
flag=true;
break;
}
linetemp=linetemp->next;
}
if(!flag)
{
linenode=new Linelist;
linenode->line=new Line;
linenode->line->pixel1=pt3;
linenode->line->pixel2=pt1;
linenode->line->pixel3=pt2;
linenode->line->flag=false;
linenode->next=NULL;
linelast->next=linenode;
linelast=linenode;
ln2=linenode->line;
}
linetemp=linehead;
flag=false;
while(linetemp!=NULL)
{
if((pt2==linetemp->line->pixel1 && pt3==linetemp->line->pixel2)
|| (pt3==linetemp->line->pixel1 && pt2==linetemp->line->pixel2))
{
linetemp->line->flag=true;
flag=true;
break;
}
linetemp=linetemp->next;
}
if(!flag)
{
linenode=new Linelist;
linenode->line=new Line;
linenode->line->pixel1=pt2;
linenode->line->pixel2=pt3;
linenode->line->pixel3=pt1;
linenode->line->flag=false;
linenode->next=NULL;
linelast->next=linenode;
linelast=linenode;
ln3=linenode->line;
}
tglnode=new Triangle;
tglnode->line1=ln1;
tglnode->line2=ln2;
tglnode->line3=ln3;
tglnode->next=NULL;
tgllast->next=tglnode;
tgllast=tglnode;

/*//test
TIN=tglhead;
TINCreateOK=true;
PostMessage(WM_PAINT);
AfxMessageBox("ok");
//test end*/
}
linetmp->line->flag=true;
linetmp=linetmp->next;
}
TIN=tglhead;
return TIN;