由不规则离散点集生成TIN与GRID DEM.doc

由不规则离散点集生成TIN与GRID DEM现阶段生成数字高程模型（DEM）的方法较多，如以摄影测量得到的像对为数据源跟踪生成等高线及DEM，由机载激光测距仪记录规则点集后生产数据，也可采用传统的地形图扫描后跟踪等高线，记录一连串离散点集，接着运用各类算法进行处理，最后生成不规则三角网（TIN）与规则格网（GRID）DEM的方法。本文主要介绍的就是以等高线（参考图一）和离散点集为数据源，产生TIN与GRID DEM的技术路线。具体步骤如下：1）跟踪等高线生成离散点集，记录在文本文件中。参考图二和图三。2）读取文本文件中的数据，进行预处理。主要工作是找到XY轴方向上最小最大数值，压缩数据范围，避免数据范围跨度太大或太小，即出现数据分布稠密或稀疏的情况。 while (!_demfile.eof() _demfile point3dXYZi0 point3dXYZi1 point3dXYZi2; point3dXYZi2 = point3dXYZi2 / 2; /因为XY轴在随后调整，因此相应调Z轴数值 if(xMinpoint3dXYZi0) xMin = point3dXYZi0;/得到整个范围的最大与最小数值 if(xMaxpoint3dXYZi1) yMin = point3dXYZi1; if(yMaxpoint3dXYZi2) zMin = point3dXYZi2; if(zMax= 0 & iNeighbor = 0 & iTriangle m_Extent.Get_YRange() ? m_Extent.Get_XRange() : m_Extent.Get_YRange(); /取范围中XY两轴的较大数值xmid= m_Extent.Get_XCenter(); /中心点 ymid= m_Extent.Get_YCenter();PointsnPoints + 0-m_Point0 = xmid - 20 * dmax;/第一个顶点的X数值PointsnPoints + 0-m_Point1 = ymid - dmax;/第一个顶点的Y数值，此数值取得相当小PointsnPoints + 1-m_Point0 = xmid; /第二个顶点的X数值PointsnPoints + 1-m_Point1 = ymid + 20 * dmax;/第二个顶点的Y数值，居中，但偏高PointsnPoints + 2-m_Point0 = xmid + 20 * dmax;/第三个顶点的X数值PointsnPoints + 2-m_Point1 = ymid - dmax;/第三个顶点的Y数值，第三个顶点与第一个顶点高度相同，在X轴上对称Triangles0.p1= nPoints + 0; /三角形的顶点记录的是顶点数组中的索引Triangles0.p2= nPoints + 1;Triangles0.p3= nPoints + 2; complete 0= false; nTriangles = 1; /初始化三角形个数为1 /- /每次增加一个顶点到已存在的TIN格网中 for(i=0; im_Point0; /使用浮点数进行记录 yp= Pointsi-m_Point1; /使用浮点数进行记录 nedge = 0; /对于每个顶点，初始化边的个数为0 /如果顶点位于三角形的外接圆内，则原三角形的三条边被加到一个缓存中，并移走此三角形 for(j=0; jm_Point0;/三角形的三个顶点 y1 = PointsTrianglesj.p1-m_Point1; x2 = PointsTrianglesj.p2-m_Point0; y2 = PointsTrianglesj.p2-m_Point1; x3 = PointsTrianglesj.p3-m_Point0; y3 = PointsTrianglesj.p3-m_Point1; inside = _CircumCircle(xp, yp, x1, y1, x2, y2, x3, y3, &xc, &yc, &r);/判断顶点与外接圆的关系 if( xc + r = emax ) emax+= 100; /增加边的个数 if( (edges = (TTIN_Edge *)realloc(edges, emax * sizeof(TTIN_Edge) = NULL ) /增加新的内存分配 status = 3; if( edges ) free(edges); if( complete ) free(complete); edgesnedge + 0.p1= Trianglesj.p1; /存储边的信息 edgesnedge + 0.p2= Trianglesj.p2; /顶点的索引：1,2; 2,3; 3,1 edgesnedge + 1.p1= Trianglesj.p2; edgesnedge + 1.p2= Trianglesj.p3; edgesnedge + 2.p1= Trianglesj.p3; edgesnedge + 2.p2= Trianglesj.p1; nedge += 3; /对于最外层循环的顶点来说，相关边的个数加上3 Trianglesj= TrianglesnTriangles - 1;/将链表中最后一个三角形赋值过来 complete j= complete nTriangles - 1; nTriangles-; /三角形个数减去 j-; /对应：if( inside ) /结束对所有三角形的循环 /标记被重复使用的边 /若三角形以逆时针方向存储，则内部被多次使用的边存储的顶点顺序是相反的 for(j=0; jnedge-1; j+) for(k=j+1; knedge; k+) if( (edgesj.p1 = edgesk.p2) & (edgesj.p2 = edgesk.p1) ) /同样的一条边存储的方向刚刚相反 edgesj.p1 = -1; /做标记 edgesj.p2 = -1; edgesk.p1 = -1; edgesk.p2 = -1; if( (edgesj.p1 = edgesk.p1) & (edgesj.p2 = edgesk.p2) ) /同样的一条边 edgesj.p1 = -1; edgesj.p2 = -1; edgesk.p1 = -1; edgesk.p2 = -1; /组成新三角形 /忽略所有已经打过标志的边 /所有的边以顺时针方向安排 for(j=0; jnedge; j+) if( edgesj.p1 0 | edgesj.p2 = trimax ) status = 4; if( edges ) free(edges); if( complete ) free(complete); TrianglesnTriangles.p1= edgesj.p1; /记录关联的边的一个顶点索引 TrianglesnTriangles.p2= edgesj.p2;/记记录关联的边的另一个顶点索引 TrianglesnTriangles.p3= i;/对于最外层循环的顶点来说，记录此顶点的索引 complete nTriangles = false; nTriangles+; /对所有的顶点循环完毕 /对应着for(i=0;i;i+) /- /Remove triangles with supertriangle vertices /移走初始定义的大外接三角形 /These are triangles which have a vertex number greater than nPoints /这些三角形的顶点索引超出了顶点个数 for(i=0; i= nPoints |Trianglesi.p2 = nPoints |Trianglesi.p3 = nPoints ) Triang