本次我们分享使用Matlab进行点云的中值滤波。中值滤波是一种非线性数字滤波技术，最早用于图像处理领域，后被成功应用于点云数据处理。在Matlab环境中，中值滤波通过选择每个点的邻域点集，对这些点的坐标值进行排序，然后用中值替代原始点的坐标，从而有效去除噪声点。
与上一节的均值滤波相比，中值滤波在去除脉冲噪声的同时能够更好地保留点云的几何特征和边缘信息，特别适用于去除孤立的异常点（俗称"飞点"）。
一、Matlab实现的主要流程
1. 数据准备与导入
% 读取点云数据 ptCloud = pcread('input_pointcloud.pcd'); % 或从矩阵创建点云 % xyzPoints = [x, y, z]; % ptCloud = pointCloud(xyzPoints);2. 邻域搜索与构建
Matlab提供了多种邻域搜索方法：
- K近邻搜索：找到距离目标点最近的K个点
- 半径搜索：找到指定半径范围内的所有点
- 球体查询：定义球形邻域进行搜索% K近邻搜索示例 k = 20; % 设置邻域大小 [indices, distances] = findNearestNeighbors(ptCloud, ptCloud.Location, k);3. 中值滤波核心算法
function filteredCloud = medianFilterPointCloud(ptCloud, k) % 获取点的数量 numPoints = ptCloud.Count; % 初始化输出点云 filteredLocation = zeros(numPoints, 3); % 对每个点进行中值滤波 for i = 1:numPoints % 找到当前点的k个最近邻 [indices, ~] = findNearestNeighbors(ptCloud, ptCloud.Location(i,:), k); % 获取邻域点的坐标 neighborPoints = ptCloud.Location(indices, :); % 对每个坐标维度分别计算中值 filteredLocation(i, 1) = median(neighborPoints(:, 1)); filteredLocation(i, 2) = median(neighborPoints(:, 2)); filteredLocation(i, 3) = median(neighborPoints(:, 3)); end % 创建滤波后的点云 filteredCloud = pointCloud(filteredLocation); end4. 自适应中值滤波改进
为提高滤波效果，可采用自适应策略：function filteredCloud = adaptiveMedianFilter(ptCloud, minK, maxK, threshold) numPoints = ptCloud.Count; filteredLocation = ptCloud.Location; for i = 1:numPoints currentK = minK; while currentK <= maxK % 获取当前邻域 [indices, distances] = findNearestNeighbors(ptCloud, ptCloud.Location(i,:), currentK); % 计算邻域内点的坐标方差 neighborPoints = ptCloud.Location(indices, :); variance = var(neighborPoints); % 如果方差小于阈值，应用滤波 if max(variance) < threshold filteredLocation(i, :) = median(neighborPoints, 1); break; end currentK = currentK + 5; % 增加邻域大小 end end filteredCloud = pointCloud(filteredLocation); end5. 性能优化策略
- 并行计算：利用Matlab的并行计算工具箱
- 空间索引：使用KD树或八叉树加速邻域搜索
- GPU加速：对于大规模点云，考虑GPU并行处理% 并行化处理示例 parfor i = 1:numPoints % 并行执行滤波操作 end三、主要应用领域
1. 三维重建预处理
在从多视角图像或激光扫描进行三维重建时，中值滤波能够有效去除由于传感器噪声、匹配错误等产生的异常点，提高重建模型的质量。2. 自动驾驶与环境感知
自动驾驶车辆获取的激光雷达点云数据常含有噪声，中值滤波可以：
- 去除雨雪等天气造成的虚假点
- 清理车辆运动产生的拖尾噪声
- 提高障碍物检测的准确性3. 工业检测与质量控制
在工业产品表面检测中，中值滤波能够：
- 去除扫描设备产生的系统噪声
- 保留产品表面的真实几何特征
- 为后续缺陷检测提供干净的输入数据4. 建筑与文物保护
对于建筑物或文物的三维扫描数据：
- 去除由于反射、遮挡等造成的异常点
- 保持建筑细节和文物特征
- 为数字化保护提供高质量数据5. 医疗影像处理
在医学CT或MRI数据的三维重建中：
- 去除图像重建过程中的伪影
- 保持器官和组织的边界清晰度
- 提高后续诊断和分析的准确性四、优缺点分析
1. 优点：
- 保持边缘：相比线性滤波，中值滤波能更好地保持点云的几何特征
- 去除脉冲噪声：对孤立的异常点特别有效
- 计算简单：算法实现相对简单，易于优化
- 适应性强：可以根据点云特性调整参数2. 缺点：
- 计算量大：对每个点都需要进行邻域搜索和排序
- 参数敏感：邻域大小选择不当可能导致过度平滑或滤波不足
- 对高密度噪声效果有限：当噪声点密度较高时，效果会下降
- 可能改变点密度：滤波后点云的分布特性可能发生变化五、最佳实践建议
1. 参数调优：根据点云密度和噪声特性选择合适的邻域大小
2. 预处理配合：可与其他滤波方法（如统计滤波、半径滤波）结合使用
3. 分块处理：对于大规模点云，采用分块处理提高效率
4. 质量评估：使用点云质量评估指标验证滤波效果
5. 保持备份：滤波前保存原始数据，便于对比和回退通过合理应用Matlab中的中值滤波技术，可以显著提升点云数据质量，为后续的三维建模、分析和应用奠定良好基础。
本次使用的数据是——————窗帘！
一、点云中值滤波的程序
1、最简版
%% 0. 清空环境
clear; clc; close all;
%% 1. 读入点云
[file, path] = uigetfile({'*.ply;*.pcd;*.xyz', '点云文件 (*.ply,*.pcd,*.xyz)'}, ...
'请选择点云');
if file == 0; return; end
fname = fullfile(path, file);
ptCloud = pcread(fname);
N = ptCloud.Count;
fprintf('原始点云有 %d 个点\n', N);
%% 2. 中值滤波（半径 50 mm）
radius = 50;
cloudMedian = pcdMedianFilter(ptCloud, radius);
fprintf('中值滤波(半径=%d mm)后 %d 个点\n', radius, cloudMedian.Count);
%% 3. 可视化
figure('Name', '原始点云', 'NumberTitle', 'off');
pcshow(ptCloud); axis on; view(3);
figure('Name', '中值滤波', 'NumberTitle', 'off');
pcshow(cloudMedian); axis on; view(3);
function out = pcdMedianFilter(pc, radiusMm)
% 中值滤波（模仿 block-wise 风格）
% pc : pointCloud 对象
% radiusMm : 邻域半径，单位 mm
xyz = pc.Location;
radius = radiusMm;
n = size(xyz, 1);
newXYZ = zeros(n, 3, 'single');
kd = createns(xyz, 'NSMethod', 'kdtree');
block = 5000; % 按内存可调
for i = 1:block:n
ir = min(i + block - 1, n);
idxCell = rangesearch(kd, xyz(i:ir, :), radius);
for j = i:ir
idx = idxCell{j-i+1};
newXYZ(j, :) = median(xyz(idx, :), 1); % 中值坐标
end
end
% 重建 pointCloud，保留颜色等信息
out = pointCloud(newXYZ, 'Color', pc.Color, ...
'Normal', pc.Normal);
end
2、GUI版本
function medianFilterGUI
% 中值滤波 GUI —— 2020a 兼容
% 1. 浏览选点云 2. 邻域半径/mm 滑块 3. 实时中值滤波 4. 保存结果
fig = figure('Name','中值滤波工具','NumberTitle','off',...
'MenuBar','none','ToolBar','none','Position',[100 100 1280 720]);
%% ---------- 左侧图像区（78 %） ----------
imgWidth = 0.78;
panelW = imgWidth/2 - 0.01;
pnlOrig = uipanel('Parent',fig,'Units','normalized',...
'FontSize',16,'Position',[0.02 0.02 panelW 0.96],'Title','原始点云');
pnlFilt = uipanel('Parent',fig,'Units','normalized',...
'FontSize',16,'Position',[0.02+panelW+0.01 0.02 panelW 0.96],'Title','中值滤波');
axOrig = axes('Parent',pnlOrig,'Units','normalized','Position',[0.05 0.05 0.90 0.90]);
axFilt = axes('Parent',pnlFilt,'Units','normalized','Position',[0.05 0.05 0.90 0.90]);
%% ---------- 右侧控制区（22 %） ----------
pnlCtrl = uipanel('Parent',fig,'Units','normalized',...
'FontSize',16,'Position',[0.78 0 0.22 1],'Title','控制');
txtH = 0.04;
btnH = 0.06;
gap = 0.02;
yTop = 0.94;
% 1. 浏览
uicontrol('Parent',pnlCtrl,'Style','pushbutton','String','浏览…',...
'FontSize',16,'Units','normalized','Position',[0.05 yTop-btnH 0.90 btnH],...
'Callback',@loadCloud);
yTop = yTop - btnH - gap;
lblInfo = uicontrol('Parent',pnlCtrl,'Style','text','String','未加载点云',...
'FontSize',10,'Units','normalized','Position',[0.05 yTop-txtH 0.90 txtH],...
'HorizontalAlignment','left');
yTop = yTop - txtH - gap;
% 2. radius/mm
uicontrol('Parent',pnlCtrl,'Style','text','String','邻域半径 / mm',...
'FontSize',12,'FontWeight','bold','Units','normalized','Position',[0.05 yTop-txtH 0.90 txtH],...
'HorizontalAlignment','left');
yTop = yTop - txtH - gap;
sliderR = uicontrol('Parent',pnlCtrl,'Style','slider','Min',1,'Max',100,'Value',50,...
'FontSize',16,'Units','normalized','Position',[0.05 yTop-btnH 0.65 btnH],...
'Callback',@refreshFilter);
txtR = uicontrol('Parent',pnlCtrl,'Style','edit','String','50',...
'FontSize',16,'Units','normalized','Position',[0.75 yTop-btnH 0.20 btnH],...
'Callback',@editRCB);
yTop = yTop - btnH - gap;
% 3. 保存
uicontrol('Parent',pnlCtrl,'Style','pushbutton','String','保存滤波结果',...
'FontSize',16,'Units','normalized','Position',[0.05 yTop-btnH 0.90 btnH],...
'Callback',@(s,e)saveCloud(ptCloudFilt));
%% ---------- 数据 ----------
ptCloudOrig = pointCloud.empty;
ptCloudFilt = pointCloud.empty;
%% ---------- 回调 ----------
function loadCloud(~,~)
[file,path] = uigetfile({'*.pcd;*.ply;*.xyz','点云文件'},'选择点云');
if isequal(file,0), return; end
try
ptCloudOrig = pcread(fullfile(path,file));
catch ME
errordlg(ME.message,'读取失败'); return;
end
showPointCloud(axOrig,ptCloudOrig);
N = ptCloudOrig.Count;
set(lblInfo,'String',sprintf('已加载：%s (%d 点)',file,N));
refreshFilter();
end
function refreshFilter(~,~)
if isempty(ptCloudOrig), return; end
r = get(sliderR,'Value');
ptCloudFilt = pcdMedianFilter(ptCloudOrig,r);
showPointCloud(axFilt,ptCloudFilt);
set(txtR,'String',num2str(r));
end
function editRCB(src,~)
v = str2double(get(src,'String'));
if isnan(v), v = 50; end
v = max(1,min(100,v));
set(sliderR,'Value',v);
refreshFilter();
end
function saveCloud(cloud)
if isempty(cloud)
errordlg('请先完成中值滤波','提示'); return;
end
[file,path] = uiputfile({'*.pcd','PCD';'*.ply','PLY';'*.xyz','XYZ'},'保存滤波点云');
if isequal(file,0), return; end
try
pcwrite(cloud,fullfile(path,file),'Precision','double');
msgbox('保存成功！','提示');
catch ME
errordlg(ME.message,'保存失败');
end
end
function showPointCloud(ax,pc)
cla(ax); set(ax,'Color','w');
pcshow(pointCloud(nan(0,3)),'Parent',ax); % 2020a 暖启动
pcshow(pc,'Parent',ax,'MarkerSize',35);
axis(ax,'tight'); grid(ax,'on'); view(ax,3);
end
%% ---------- 中值滤波核心 ----------
function out = pcdMedianFilter(pc,radiusMm)
xyz = pc.Location;
n = size(xyz,1);
radius = radiusMm;
kd = createns(xyz,'NSMethod','kdtree');
block = 10000; % 批查询
xyzMed = zeros(n,3,'like',xyz);
for i = 1:block:n
ir = min(i+block-1,n);
idxCell = rangesearch(kd, xyz(i:ir,:), radius);
for k = i:ir
idx = idxCell{k-i+1};
if numel(idx)>0
xyzMed(k,:) = median(xyz(idx,:),1); % 逐维中值
else
xyzMed(k,:) = xyz(k,:); % 无邻域保持原值
end
end
end
out = pointCloud(xyzMed, ...
'Color', pc.Color, ...
'Normal', pc.Normal);
end
end
二、点云中值滤波的结果
可以看到，随着滑块的变化，窗帘被动态滤波了。感兴趣的童鞋可以试着调一调别的点云呀，绝知此事要躬行！
就酱，下次见^-^
扫一扫
352
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
