点击上方“阿昆的科研日常”，关注我的公众号。
设为“星标”，原创技术文章第一时间推送。
在之前的文章中，分享过Matlab基于KD树的邻域点搜索方法：
在此基础上，进一步分享一下基于KD树的离散点密度特征提取方法。
先来看一下成品效果：
点云密度特征一般用单位面积/体积内的离散点数量表示。
其中，二维(2-Dimension, 2D)平面密度可以通过柱状邻域提取；三维(3-Dimension, 3D)体密度可以通过球邻域提取。
function density_2D = density2D_KD(data,radius)% 功能：利用KD树提取离散点2D密度特征% 输入：data - 原始数据(m*3)% 输出：planes - 拟合所得平面参数M = size(data,1);density_2D = zeros(M,1);idx = rangesearch(data(:,1:2),data(:,1:2),radius,'Distance','euclidean','NSMethod','kdtree');for i = 1:Mdensity_2D(i,1) = length(idx{i})/(pi*radius^2);endend
function density_3D = density3D_KD(data,radius)% 功能：利用KD树提取离散点3D密度特征% 输入：data - 原始数据(m*3)% 输出：planes - 拟合所得平面参数M = size(data,1);density_3D = zeros(M,1);idx = rangesearch(data(:,1:3),data(:,1:3),radius,'Distance','euclidean','NSMethod','kdtree');for i = 1:Mdensity_3D(i,1) = length(idx{i})/(4/3*pi*radius^3);endend
为了检测密度特征提取的效果，采用某地机载LiDAR实测点云数据进行验证：
%% 数据准备% 读取数据data = load('data1.txt');%% 点云密度特征提取% 定义半径radius = 1;% 2D平面密度density_2D = density2D_KD(data(:,1:2),radius);% 3D体密度density_3D = density3D_KD(data(:,1:3),radius);%% 原始点云可视化% 窗口尺寸设置（单位：厘米）figureUnits = 'centimeters';figureWidth = 2*10;figureHeight = 1*10;figureHandle1 = figure;set(gcf, 'Units', figureUnits, 'Position', [0 0 figureWidth figureHeight]);t1 = tiledlayout(1,2,'TileSpacing','compact','Padding','compact');nexttilescatter(data(:,1),data(:,2),7,data(:,3),'filled')hTitle1 = title('原始点云(俯视)');set(gca,'xtick',[])set(gca,'ytick',[])set(gca,'ztick',[])axis equal off tightcolormap(map)colorbarnexttilescatter3(data(:,1),data(:,2),data(:,3),5,data(:,3),'filled')hTitle2 = title('原始点云(侧视)');set(gca,'xtick',[])set(gca,'ytick',[])set(gca,'ztick',[])axis equal off tightview(-1.157471311476228e+02,53.408837245676686)colorbar% 细节优化set([hTitle1,hTitle2], 'FontName', '微软雅黑', 'FontSize', 12, 'FontWeight' , 'bold')% 背景颜色set(gcf,'Color',[1 1 1])% 图片输出figW = figureWidth;figH = figureHeight;set(figureHandle1,'PaperUnits',figureUnits);set(figureHandle1,'PaperPosition',[0 0 figW figH]);fileout = 'test0';print(figureHandle1,[fileout,'.png'],'-r300','-dpng');%% 平面、体密度提取结果可视化验证% 窗口尺寸设置（单位：厘米）figureUnits = 'centimeters';figureWidth = 2.35*10;figureHeight = 1*10;figureHandle2 = figure;set(gcf, 'Units', figureUnits, 'Position', [0 0 figureWidth figureHeight]);t2 = tiledlayout(1,2,'TileSpacing','compact','Padding','compact');nexttilescatter(data(:,1),data(:,2),15,density_2D(:,1),'filled')hTitle1 = title('基于KD树的2D平面密度特征可视化');set(gca,'xtick',[])set(gca,'ytick',[])set(gca,'ztick',[])axis equal off tightcolormap(map)colorbarnexttilescatter(data(:,1),data(:,2),15,density_3D(:,1),'filled')hTitle2 = title('基于KD树的3D体密度特征可视化');set(gca,'xtick',[])set(gca,'ytick',[])set(gca,'ztick',[])axis equal off tightcolorbar% 细节优化set([hTitle1,hTitle2], 'FontName', '微软雅黑', 'FontSize', 12, 'FontWeight' , 'bold')% 背景颜色set(gcf,'Color',[1 1 1])% 图片输出figW = figureWidth;figH = figureHeight;set(figureHandle2,'PaperUnits',figureUnits);set(figureHandle2,'PaperPosition',[0 0 figW figH]);fileout = 'test';print(figureHandle2,[fileout,'.png'],'-r300','-dpng');
其中，为了区分不同对象，从Matlab配色神器TheColor的SCI颜色库中选择渐变色：
% 颜色定义map = TheColor('sci',2068);map = flipud(map);
(点击上图查看TheColor具体功能)
获取方式：公众号（阿昆的科研日常）后台回复 TC
最终结果如下：
以上。
如果你觉得我的分享对你有帮助的话，欢迎大家在这里点赞、在看、分享。当然，也欢迎大家在这里打赏。互动越多，更新越快哦~
声明：本公众号的所有原创内容，在未经允许的情况下，不得用于任何商业用途，违者必究。
