基于MATLAB的道格拉斯普克算法递归实现
基于MATLAB的道格拉斯普克算法递归实现@TOC
道格拉斯普克算法
（道格拉斯-普克）Douglas-Peukcer算法由D.Douglas和T.Peueker于1973年提出，是线状要素抽稀的经典算法。用它处理大量冗余的几何数据点，既可以达到数据量精简的目的，有可以在很大程度上保留几何形状的骨架。
递归算法
递归是一种代码运行的机制，简单来说就是函数自己调用自己实现有条件的迭代。在道格拉斯普克算法中，处理的原始数据为折线段，首先将折线段有条件的划分为两段（有条件所指的这个条件就是某点距离折线首尾连线距离足够远，也就是说这个点足够突出，抽稀时不可忽视，需要保留！）。以这个需要保留的特征点为分界点，将原始折线分成两段，每一段再调用当前函数继续提取特征点和进一步分段，直到所有特征点提取完毕，停止递归迭代。话不多说，上代码：
递归函数代码（ARecursionFun.m）
%% 设计一个递归函数（递归函数的定义初中学过，注意复习！！！），依次寻找特征点并存储
% 输入两个相邻特征点之间的扫描线pointsTab，特征点表A(A是折线首尾两个端点)
% 输出补充新发现的特征点后的特征点表A
% 函数名称为ARecursionFun（一个递归函数）
function [A] = ARecursionFun(pointsTab,A,Threshold)
[r,~] = size(pointsTab); % 获取扫描线片段上点的个数
if r > 2 % 如果这条扫描线片段上点数大于2则执行操作
Q1 = [pointsTab(1,1);pointsTab(1,2)]; % 起点坐标对的列向量表示（为了便于点到直线距离计算的表示方法）
Q2 = [pointsTab(r,1);pointsTab(r,2)]; % 终点坐标对的列向量表示（作用同上）
% 遍历这个扫描线，依次计算每个点到扫描线起点终点连线的距离==================
for i = 1:1:r
P = [pointsTab(i,1);pointsTab(i,2)]; % 当前点坐标的列向量表示
d(i,1) = abs(det([Q2-Q1,P-Q1]))/norm(Q2-Q1); % 计算点到直线的距离
end
% 计算完毕，每个点到直线的距离存入列向量d中================================
if max(d) > Threshold % 如果距离列向量中最大值大于阈值则进行下述操作
ind = find(all(repmat(max(d),size(d,1),1)==d,2)); % 获取列向量中最大值对应的点的序号
[rA,~] = size(A); % 获取当前特征点表A已存点的个数
A(rA+1,1) = pointsTab(ind,1); % 将这个点作为特征点存储起来（x坐标）
A(rA+1,2) = pointsTab(ind,2); % 将这个点作为特征点存储起来（y坐标）
%{
到上述步骤的时候其实就是一条扫描线上找到了距离扫描线首尾连线最远的那一个点并存储起来
你应该是看得懂的！！！
接下来两句则是该函数，或者说是该类函数（递归函数）的精髓，加几个断点，认真看！！！
%}
pointsTabb = pointsTab(1:ind,:); % 以刚才存储的特征点为界限，起点到该点建立新的片段扫描线
A = ARecursionFun(pointsTabb,A,Threshold); % 函数自身调用进行递归，进一步获取片段内的特征点
pointsTabe = pointsTab(ind:r,:); % 以刚才存储的特征点为界限，该点到终点建立新的片段扫描线
A = ARecursionFun(pointsTabe,A,Threshold); % 函数自身调用进行递归，进一步获取片段内的特征点
end
end
调用函数(main.m)
%% 主函数入口（在该函数界面下点击运行实验）
%{
程序运行结果图中，黑色线为原始扫描线的连线，红色*为获取的特征点
红色线为特征点连线
工作区中，points为扫描线点集，A为特征点集合
%}
clc;clear;close all; %清空所有变量，开始崭新的工作！！
% 模拟一条折线（用的时候替换成你自己的就行）=================================
points(:,1) = 5:5:300; %x值为1到60
points(:,2) = 10 + 3 * rand(60,1); %y为10加一个0到1的随机数
points(25:35,2) = 5 + 3 * rand(11,1); %其中第25到第35个点低一点
% =========================================================================
[r,c] = size(points);
A(1,1) = points(1,1); A(1,2) = points(1,2);
A(2,1) = points(r,1); A(2,2) = points(r,2);
Threshold = 3; %给定阈值
[A] = ARecursionFun(points,A,Threshold); % 递归
A = sortrows(A,1);
figure(1); %创建图层
plot(points(:,1),points(:,2),’-k’); %绘制原始折线
hold on; %保留当前图层的要素
plot(A(:,1),A(:,2),’*-r’); %在原图基础上绘制特征点
title([‘阈值为：’,num2str(Threshold)]);
msgbox ‘运行完毕！’;