Matlab三维离散点云法向量与特征值的PCA计算方法

具体代码如下：

function [normal_vector,EVs] = EV(knbpts)
% knbpts：邻近点
% normal_vector:单位法向量
% EVs ：单位特征值
P = knbpts(:,1:3);
[m,~] = size(P);
% 计算协方差矩阵
P = P-ones(m,1)*(sum(P,1)/m);
C = P.'*P./(m-1);  
% 计算特征值与特征向量
[V, D] = eig(C);
% 最小特征值对应的特征向量为法向量
s1 = D(1,1);
s2 = D(2,2);
s3 = D(3,3);
if (s1 <= s2 && s1 <= s3)
    normal_vector(1,:) = V(:,1)/norm(V(:,1));
elseif (s2 <= s1 && s2 <= s3)
    normal_vector(1,:) = V(:,2)/norm(V(:,2));
elseif (s3 <= s1 && s3 <= s2)
    normal_vector(1,:) = V(:,3)/norm(V(:,3));
end 
% 单位特征值
epsilon_to_add = 1e-8;
% EVs(1) > EVs(2) > EVs(3)
EVs(1,:) = [D(3,3) D(2,2) D(1,1)];
if EVs(1,3) <= 0
    EVs(1,3) = epsilon_to_add;
    if EVs(1,2) <= 0
        EVs(1,2) = epsilon_to_add;
        if EVs(1,1) <= 0
            EVs(1,1) = epsilon_to_add;
        end
    end
end
sum_EVs = EVs(1) + EVs(2) + EVs(3);    
EVs = EVs/sum_EVs;
end

为了评估法向量提取效果，利用一ISPRS数据进行验证。

结果如下：

懒得动手的朋友，测试数据与示例代码见：

Matlab三维点云法向量与特征值的简易提取方法

法向量可视化代码见：

matlab normal plot.rar

邻近点搜索代码见：

Matlab建立KD树搜索三维点云的邻近点（柱状邻域、球状邻域与KNN）