ttn.cpp

三维真实地面模型中利用三角网TIN的编程

// TTN.cpp : Defines the entry point for the console application.
//

#include "stdafx.h"
#include "stdio.h"
#include "stdlib.h"
#include "math.h"
#include "iostream.h"

#define n 6							//要建立三角模型的点数目
static int w[n+1][n+1];				//邻接矩阵

int s[n+1];							//纪录已经被选的点
int x[n+1]={0,1,4,6,9,10,7};		//所有点的x坐标分量
int y[n+1]={0,4,8,3,6,2,10};		//所有点的y坐标分量
int top=0;							//已经被选的点数
int rear=0;							//被选数组的下标
int BX[n];							//每选定一条边其所有的备选点放入
double minlength=10000;				//假设一个最短距离
int temp;							//用来检测当边的扩展次数小于2的时候纪录一个待扩展点
int find3;							//找第一个三角形的第三个点
int findd;							//找到的点号
int flag=0;							//标志是否有备选点0 表示没有找到备选点 1表示找到 

int find(int t,int a[])				//检测是否点重复选择
{	
	int i;
	for(i=1;i<=top;i++)
    {
		if(t==a[i])
			return 1;
	}
	return 0;
}

int detect( int ii,int jj)//检测是否已经被扩展两次
{
	int m=0;
	int k;
	for(k=1;k<=top;k++)
		if(w[s[k]][s[ii]]==1&&w[s[k]][s[jj]]==1)
		{
			m++;
			temp=k;
		}
	return m;
}


int main(int argc, char* argv[])
{
/********   首先找出第一个三角形  ********/
	int record;							//纪录第二个点
	int count;							//用来返回边扩展的次数
	double d;
	double  a1,b1,c1,angle;
	double F1,F2;						//代表两条直线方程
	double maxangle=0;					//假设一个最大的角度值
	int k;
	s[1]=1;								//选第一个点为开始点
	for(int i=2;i<=n;i++)				//找第二个点（利用距离公式）
	{ 
		d=sqrt((x[1]-x[i])*(x[1]-x[i])+(y[1]-y[i])*(y[1]-y[i]));
		//cout<<d<<endl;
		if(d<minlength)
		{
			record=i;
			minlength=d;
		}
	}

	w[1][record]=w[record][1]=1;		//邻接矩阵赋值
	s[2]=record;						//把第二个点的点号放入纪录数组
	top=2;

	for( k=1;k<=n;k++)//找第一个三角形的第三个点
	{
		if(!find(k,s))
		{
			a1=sqrt((x[k]-x[1])*(x[k]-x[1])+(y[k]-y[1])*(y[k]-y[1]));
			b1=sqrt((x[k]-x[record])*(x[k]-x[record])+(y[k]-y[record])*(y[k]-y[record]));
			c1=sqrt((x[record]-x[1])*(x[record]-x[1])+(y[record]-y[1])*(y[record]-y[1]));
			angle=acos((a1*a1+b1*b1-c1*c1)/(2*a1*b1));
			if(angle>maxangle)				//!!!!!!!!!!!!!!!符号变了!!!!!!!!!!!
			{
				maxangle=angle;
				find3=k;
			}
		}
	}
	//cout<<"find3="<<find3<<endl;
	//cout<<maxangle<<endl;
	s[3]=find3;									//把第三个点放入纪录数组
	w[1][find3]=w[find3][1]=1;					//邻接矩阵赋值
	w[record][find3]=w[find3][record]=1;		//邻接矩阵赋值
	top=3;										//点号加1
	cout<<"第一个点"<<1<<"第二个点"<<record<<"第三个点"<<find3<<endl;

 
// 到此已经做好所有准备工作并且已经选了一个三角形接着扩展
	for(i=1;i<=top&&top<n;i++)					//当点没有扩展完继续扩展
		for(int j=1;j<=top&&top<n;j++)
		{ 
			rear=0;
			maxangle=0;
			flag=0;
			if(w[s[i]][s[j]]==0)				//如果找到的两个点之间没有边继续下一次循环
				continue;
			if(w[s[i]][s[j]]==1)				//如果有边
			{
				count=detect(i,j);				//返回该边的扩展次数
				if(count==2)					//扩展两次了不能再扩展 继续下一次
				{
					continue;
				} 
				if(count<2)						//可以扩展执行下列代码
				{
					F1=(y[s[j]]-y[s[i]])*(x[s[temp]]-x[s[i]])-(x[s[j]]-x[s[i]])*(y[s[temp]]-y[s[i]]);
					for(int p=1;p<=n;p++)
					{
						if(!find(p,s))
						{
							F2=(y[s[j]]-y[s[i]])*(x[p]-x[s[i]])-(x[s[j]]-x[s[i]])*(y[p]-y[s[i]]);
							if(F1*F2<0)
							{
								flag=1;
								w[s[i]][p]=w[p][s[i]]=1;
								w[s[j]][p]=w[p][s[j]]=1;
								cout<<"第一个点"<<s[i]<<"第二个点"<<s[j]<<"第三个点"<<p<<endl;
								BX[rear]=p;
							}
						}
					}
					for( k=1;k<=rear;k++)				//在备选数组中找到一个最大角度的点
					{
						if(!find(k,BX))
						{
							a1=sqrt((x[BX[k]]-x[s[i]])*(x[BX[k]]-x[s[i]])+(y[BX[k]]-y[s[i]])*(y[BX[k]]-y[s[i]]));
							b1=sqrt((x[BX[k]]-x[s[j]])*(x[BX[k]]-x[s[j]])+(y[BX[k]]-y[s[j]])*(y[BX[k]]-y[s[j]]));
							c1=sqrt((x[s[i]]-x[s[j]])*(x[s[i]]-x[s[j]])+(y[s[i]]-y[s[j]])*(y[s[i]]-y[s[j]]));
							//cout<<a1<<"  "<<b1<<"   "<<c1<<endl;
							angle=acos((a1*a1+b1*b1-c1*c1)/(2*a1*b1));
							if(angle>maxangle)
							{
								maxangle=angle;
								findd=BX[k];
							}
						}
					} 
					if(flag==1)							//在前提有被选点的情况下，求出的点号存放到纪录数组中
					{
						s[++top]=findd;
						w[s[i]][findd]=w[findd][s[i]]=1;
						w[s[j]][findd]=w[findd][s[j]]=1;
					}
				}
			}
}
	
	return 0;
}