道格拉斯—普克法(Douglas—Peucker)简化线算法报告

线简化算法程序设计报告矢量数据是GIS中，使用非常普遍的一种数据类型。在使用中，有时需要对矢量数据进行压缩，矢量数据压缩的目的是删除冗余数据，减少数据的存贮量，节省存贮空间，加快后继处理的速度。矢量数据的压缩方法常用的有道格拉斯—普克法、垂距法、光栏法。本文主要讨论道格拉斯—普克法，运用该算法的思想，用C语言于TC20中编写一个小程序，实现对既有线的简化。首先，简要介绍下道格拉斯—普克法(Douglas—Peucker)的核心思想：基本思路(如图1)：对每一条曲线的首末端点连一条线，求所有点到该直线的距离，并找出最大距离值dmax，用dmax与限差D相比：（1）若dmax＜D，这条曲线上的中间点全部舍去；（2）若dmax≥D，保留dmax对应的坐标点，并以该点为界，把曲线分为两部分，对这两部分重复使用该方法。图1由该算法的基本思路可知，该算法是递归的，而且算法的核心是求得点到直线的距离。根据以上分析，结合老师给出的数据，编写出了程序代码。经过调试，能实现线简化的功能。使用TC20做为程序的开发环境，主程序流程图如下所示：该算法的主要功能函数是Simp()，该函数的功能是根据输入的限差LimitDis确定线上的点是保留还是去除。Simp()函数的实现思想是，将曲线上的各点与曲线的端点连接，组成一个三角形，如下：先根据点到点的距离公式，求出a,b,c的值。然后再应用余弦公式cosA=(b*b+c*c-a*a)/(2*b*c);cosB=(a*a+c*c-b*b)/(2*a*c);sinA=(float)sqrt(1-cosA*cosA);再根据cosA，cosB，sinA的值，算出距离。对每个点进行上述的计算，求出每个点到AB的距离，取最大值和LimitDis比较，如果大于限差，那么从C点将曲线分为二段，重复第一步。如果小于，则去除AB间的所有点。该函数的流程图如下：用Savetostruct()函数将坐标数据读入结构数组以文字形式显示读入的坐标数据调用DrawLine()函数以图形形式显示坐标数据根据输入的LimitDis值，调用Simp()函数简化线根据输入的LimitDis值，调用Simp()函数简化线分别以文字和图形的形式显示简化后的结果CABbca算法源代码：#includestdio.h#includemath.h#includegraphics.h计算AB的距离:c计算i点到A,B的距离，得到b,a;根据a,b,c求sinA,cosA,cosB判断cosA，cosB是否为零curdis=b*sinAYNcurdis=b或者a判断curdis=MaxDisCurdis赋值给MaxDis,i赋值给MaxNOY判断MaxDis=LimitDisSimp(A点,MaxNO点);Simp(MaxNo点,B点);Simp(maxNO,p2);Y删除AB间的所有点N/*---------------------------------------------------------*//*----------------definethepointsstruct-----------------*//*---------------------------------------------------------*/typedefstructPoint{floatx,y;charRes;};FILE*fp;intiNum;floatLimitDis;structPointPt[30];/*-------------------------------------------------------*//*----------savethepointstothestructtype------------*//*-------------------------------------------------------*/voidSavetostruct(){floatpx,py;intj=0;if((fp=fopen(C:\\data.txt,r))==NULL){printf(Can'tOpentheFile!!!\n);exit(0);}while(!feof(fp)){fscanf(fp,%f%f,&px,&py);Pt[iNum].x=px;Pt[iNum].y=py;Pt[iNum].Res='T';iNum++;}printf(ThesearetheOriginalPointsbeforeSimplized:);for(j=0;j=iNum-1;j++){printf(\nOriginalPt-%dx=%fy=%f\n,j,Pt[j].x,Pt[j].y);}printf(\nToatalpoints:%d\n,j);};/*-------------------------------------------------------*//*--------------------Simplizetheline------------------*//*-------------------------------------------------------*/voidSimp(p1,p2){floata,b,c,cosA,cosB,sinA,maxdis,curdis;inti=0,maxNO=0;floatp2pdis();if((p2-p1)=2){maxdis=0.00;c=p2pdis(p1,p2);i=p1+1;while(ip2){curdis=0.00;b=p2pdis(p1,i);a=p2pdis(p2,i);cosA=(b*b+c*c-a*a)/(2*b*c);cosB=(a*a+c*c-b*b)/(2*a*c);sinA=(float)sqrt(1-cosA*cosA);if((cosA==0)||(cosB==0)){if(cosA==0){curdis=b;}else{curdis=a;}}else{curdis=b*sinA;}if(maxdis=curdis){maxdis=curdis;maxNO=i;}i++;}/*endwhile*/if(maxdis=LimitDis){Simp(p1,maxNO);Simp(maxNO,p2);}else{Delpt(p1,p2);}}/*endif*/}/*endSimp()*//*-------------------------------------------------------*//*------------------Distanceof2pts--------------------*//*-------------------------------------------------------*/floatp2pdis(pa,pb){floatd;d=(float)sqrt((Pt[pa].x-Pt[pb].x)*(Pt[pa].x-Pt[pb].x)+(Pt[pa].y-Pt[pb].y)*(Pt[pa].y-Pt[pb].y));returnd;}/*-------------------------------------------------------*//*-----------------Delpoint-----------------------------*//*-------------------------------------------------------*/Delpt(a,b){intc=a+1;while(cb){Pt[c].Res='F';c++;}}/*-------------------------------------------------------*//*-----------------DrawLine-----------------------------*//*-------------------------------------------------------*/DrawLine(){intdriver,mode,j,iSum=0,lastx,lasty;driver=DETECT;mode=0;initgraph(&driver,&mode,);setcolor(15);for(j=0;j=iNum;j++){if(Pt[j].Res=='T'){iSum++;if(iSum2){lastx=Pt[j].x;lasty=Pt[j].y;}else{line(lastx,lasty,Pt[j].x,Pt[j].y);lastx=Pt[j].x;lasty=Pt[j].y;}putpixel(Pt[j].x,Pt[j].y,RED);}}getch();restorecrtmode();}/*-------------------------------------------------------*//*---------------------MainFunction---------------------*//*-------------------------------------------------------*/main(){intk;clrscr();Savetostruct();printf(\nPressanykeytoseetheOrignalLineingraphicsmode:\n);getch();DrawLine();printf(Plsinputanumberasthelimitdistance:\n);scanf(%f,&LimitDis);Simp(0,(iNum-1));printf(ThesearethePointsafterSimplized:\n);for(k=0;k=iNum-1;k++){if(Pt[k].Res=='T')printf(ResultPt-%dx=%fy=%f\n,k,Pt[k].x,Pt[k].y);}printf(\nPressanykeytoseethesimplizedresultingraphicsmode:\n);getch();DrawLine();