基于法向集差异的点云法向估计算法

Point Cloud Normal Estimation Algorithm Based on Normal Set Difference

摘要: 法向量作为三维点云必不可少的属性，在基于点云的绘制方法和表面重建算法中起到了重要作用。利用法向差异对邻域点进行筛选赋权是法向估计中一种常用的策略，但是在现有的法向估计方法中，一般使用的是单一的、不准确的初始输入法向来对邻域点的局部结构进行刻画，使得筛选的结果并不精准。对此，本文提出了一种利用法向集之间的差异对邻域点进行筛选赋权的法向估计方法。首先通过对单一的、粗糙的初始法向进行扰动筛选，构造一个更为准确的法向集，然后根据法向集之间的差异对邻域点进行筛选赋权，最后利用加权最小二乘拟合得到估计法向。实验结果表明，该方法在不同噪声水平下均能在尖锐特征处得到更好的法向估计效果。

Abstract: The normal vector, as an essential attribute of three-dimensional point cloud, plays a significant role in point cloud-based rendering methods and surface reconstruction algorithms. Utilizing the differences in normal vectors to filter and weight neighboring points is a commonly used strategy in normal estimation. However, in the existing normal estimation methods, a single and inaccurate initial normal vector is generally used to depict the local structure of neighboring points, resulting in imprecise filtering results. To address this issue, this paper proposes a normal estimation method that uses the differences between normal sets to filter and weight neighboring points. Firstly, a more accurate normal set is constructed by perturbing a single and inaccurate initial normal vector. Then, the neighboring points are filtered and weighted based on the differences between the normal sets. Finally, the estimated normal is obtained through weighted least squares fitting. Experimental results show that this method can achieve better normal estimation effects at sharp features under different noise levels.

1. 引言

点云数据以其高精度和丰富的几何信息，为三维建模、物体识别、场景重建等任务提供了坚实的数据基础。随着三维扫描技术的飞速发展，点云数据在众多领域得到了广泛应用，包括工业制造、文化遗产保护、机器人导航、虚拟现实以及医学成像等。然而，在实际应用中，点云数据的处理面临着诸多挑战，其中法向估计是关键且基础的一环。点云法向，即点云中每个点处的曲面法线方向，不仅是点云数据的重要几何属性，更直接影响着后续三维处理任务的效果。但是，由于所采集的点云数据存在不均匀采样、噪声等问题，如何得到准确的法向估计结果成为当前研究中的一个重要课题。

针对上述难题，研究者们提出了多种法向估计的方法，主要分为四类：基于Voronoi/Delaunay的方法[1]-[3]；基于法向修正的方法[4]-[6]；基于深度学习的方法[7]-[9]；基于局部拟合的方法。

基于Voronoi/Delaunay的方法[1]-[3]的核心思想是利用计算几何中的Voronoi图和Delaunay三角剖分技术来捕捉点云中的局部几何和拓扑信息，从而构建法向信息。这种基于Voronoi/Delaunay的方法在处理相对规整的点云数据时表现良好，但当处理大噪声和具有尖锐特征的点云数据时，会出现法向估计不准确的问题[10]。

基于法向修正的方法[4]-[6]的核心思想是通过引入额外的几何信息或优化策略，对初始法向进行修正，从而得到更好的法向估计效果。该方法通过迭代修正得到更加贴合实际几何形态的法向，但该算法的局限性在于依赖初始法向的准确性，尤其在尖锐特征处初始法向计算不准确时，法向估计结果也不够准确，难以还原原始模型尖锐几何特征[11]。

基于深度学习的方法[7]-[9]的核心思想是利用神经网络学习点云局部几何特征。这种方法明显降低了计算复杂度，能够自动学习点云的全局和局部特征，使得模型能够更好地捕捉点云的几何结构，适用于复杂的点云场景。同时，网络通过大量训练数据可自适应识别噪声模式，在非均匀采样或噪声干扰情况下仍能保持鲁棒性。但深度学习模型的性能高度依赖训练数据的质量和数量，如果训练数据不足或存在偏差，模型的泛化能力可能会受到限制，导致在新的点云数据上表现不佳[9]。

基于局部拟合的方法的核心思想是利用当前点的邻域点拟合一个几何模型(如平面、二次曲面等)，从而估计该点的法向。Hoppe等[12]是这一类方法的先驱，他们提出的基于主成分分析(Principal Component Analysis, PCA)的局部切平面拟合方法，为后续研究提供了理论基础。该方法将邻域点协方差矩阵的最小特征值所对应的特征向量作为法向，具备计算简单、效率高的优势。但由于PCA方法将邻域内所有点拟合为一个平面，因此难以准确捕捉复杂曲面的真实局部几何结构，使法向估计结果不准确。对此，Guennebaud等[13]和Cazals等[14]采用球面或二次曲面拟合来替代平面近似，更加适配局部非平面的几何结构。除了改进拟合模型外，也有学者采用了对邻域点进行权重分配的方法。基于“距离当前点越近的邻域点，越能反映该点真实的曲面形态”的思想，Pauly等[15]提出了以距离为基础的加权策略，对与当前点距离更近的邻域点赋予更高的权重，有效抑制了距离当前点较远的邻域点的影响。然而这些方法都采用了各向同性的邻域，因此不能从本质上解决尖锐特征处邻域点跨越边界引起的平滑问题[10]。针对这一问题，研究者们引入了鲁棒统计技术以区分不同曲面上的点。基于“拟合残差较大的邻域点更有可能与当前点位于不同曲面上”的思想，Li等[10]和Mederos等[16]利用拟合残差对邻域点进行筛选和加权，其中Li等[10]提出了鲁棒法向估计算法(Robust Normal Estimation, RNE)。上述方法在邻域的筛选与加权上都以距离和残差为基础，没有充分考虑到点在局部结构的相似性，导致估计出的法向在尖锐特征处存在偏差，效果并不理想[17]。

在近些年的邻域筛选与加权工作当中，利用法向差异来刻画邻域点的局部结构受到了诸多关注。Öztireli等[18]将法向差异作为一个独立的权重引入迭代优化过程。在此基础上，Wang等[17]将法向差异作为核心的邻域加权依据，赋予与当前点法向差异较小的邻域点更高的权重值。Zhang等[19]和Liu等[20]分别提出了低秩子空间聚类(Low-Rank Representation, LRR)算法和LRR的加速版本(Low-Rank Representation fast, LRRfast)，该类方法将子空间分割引入到法向估计中，通过法向差异构造关系矩阵，为子空间分割做引导。Zhang等[21]利用点对的初始法向差异来刻画两点在同一平面的可靠性，设计了协同一致性投票机制(Pair Consistency Voting, PCV)。上述方法利用法向差异刻画点的局部结构，进一步抑制了与当前点位于不同曲面的邻域点的作用，有效提高了法向估计结果。但由于尖锐特征处往往包含多个不连续曲面，想要准确估计初始法向本身是困难的，因此对于复杂曲面或几何结构变化剧烈的区域，仅依靠单一的、不准确的初始输入法向来计算法向差异，得到的法向估计结果有待进一步提高。

针对上述问题，本文提出了基于法向集差异的点云法向估计算法。首先，与以往方法不同，本文通过对点云中每一点构建一个法向集来刻画其局部结构。具体而言，就是在计算初始法向后，设计扰动策略生成多个候选法向，对这些候选法向进行筛选得到法向集。这样定义的法向集避免了仅依赖单一的、不准确的初始法向导致的法向估计偏差问题，能够更稳健地刻画点的局部几何特征。其次提出一种新的计算方法来度量当前点与邻域点的法向集差异，用以评估当前点与其邻域点之间局部结构的一致性。最后根据法向集差异为每个邻域点分配权重，并利用加权最小二乘拟合得到所求法向。实验结果表明，与利用单一的、不准确的初始法向计算法向差异对邻域点进行赋权的方法相比，本文方法在尖锐特征处均得到更好的法向估计结果。

图1为本算法进行筛选扰动示意图。

Figure 1. Schematic diagram of the perturbation and screening performed by this algorithm

图1. 本算法进行扰动筛选示意图

2. 算法

给定噪声点云 ($N$ 为采样点个数，$p_i=\left\{x_i,y_i,z_i\right\}$ )，对于点云中任意一点$p_i$ ，它的邻域为$N_{p_i}=\left\{p_i^j|p_i^j是距离p_i最近的K个邻域点\right\}$ 。本文提出的基于法向集差异的点云法向估计算法包括“生成法向集–计算法向集差异–赋权拟合”三个过程。

第一步：生成法向集。首先利用PCA算法得到当前点$p_i$ 的初始法向$n_i^0$ 。其次设计了合适的扰动机制，得到$m$ 个扰动向量$t_i$ ，分别与初始法向$n_i^0$ 相加，生成$m$ 个候选法向，从而得到$m$ 个候选平面。最后通过计算点$p_i$ 的邻域$N_{p_i}$ 到候选平面的距离筛选出$k_{top}$ 个最适合的平面，将这$k_{top}$ 个平面的法向视为当前点$p_i$ 的法向集合，记为$N\left(p_i\right)$ 。

第二步：计算法向集差异。在生成每个点的法向集后，设计了合适的度量来刻画当前点$p_i$ 的法向集$N\left(p_i\right)$ 与其邻域点$p_i^j$ 的法向集$N\left(p_i^j\right)$ 之间的差异性。

第三步：赋权拟合。通过法向集差异对邻域中的每个邻域点赋权，并利用加权最小二乘拟合得到最终的法向$n_i$ 。

图2为算法流程图。

Figure 2. Algorithm flow chart

图2. 算法流程图

2.1. 生成法向集

2.1.1. 生成候选法向集

对于点$p_i$ ，协方差矩阵$C_i$ 的定义如下：

$C_i=\frac{1}{K}\left[\begin{array}{c}{p}_i^1-\overline{p_i}\\ ⋮\\ p_i^K-\overline{p_i}\end{array}\right]\cdot {\left[\begin{array}{c}{p}_i^1-\overline{p_i}\\ ⋮\\ p_i^K-\overline{p_i}\end{array}\right]}^{\text{T}},p_i^j\in N_{p_i},$(1)

其中，$N_{p_i}$ 是点$p_i$ 的邻域，$K$ 是$N_{p_i}$ 的大小，$\overline{p_i}=\frac{1}{K}\left({\displaystyle \sum _{j=1}^K{p}_i^j}\right)$ 是局部邻域$N_{p_i}$ 的平均值。

对协方差矩阵$C_i$ 进行特征分解，得到特征值${\lambda }_1$ ，${\lambda }_2$ ，${\lambda }_3$ 与它们对应的特征向量$v_1$ ，$v_2$ ，$v_3$ ，其中$0\le {\lambda }_1\le {\lambda }_2\le {\lambda }_3$ 。${\lambda }_1$ 对应的特征向量$v_1$ 就是PCA算法计算的初始法向，记为$n_i^0$ 。

由于$n_i^0$ 在尖锐特征处不够准确，所以本文设计了合适的扰动机制得到候选法向集。具体来说，定义扰动幅度$\lambda =0.15+0.05\xi$ ，其中，0.15为基础扰动幅度，$0.05\xi$ 为随机扰动分量，$\xi$ 是一个标准正态分布的随机数。设置扰动系数$\mu =\lambda \cdot \zeta$ ，其中$\zeta$ 为二维的标准正态分布随机向量。将特征向量$\left(v_1,v_2\right)$ 设置为扰动基向量，生成扰动向量$t_i=\mu \cdot \left(v_1,v_2\right)$ 。将点$p_i$ 的初始法向$n_i^0$ 与扰动向量$t_i$ 相加，再进行归一化，得到点$p_i$ 的一个候选法向$n_i^r$ 。重复上述过程生成$m$ 个候选法向，得到点$p_i$ 的候选法向集

$N_{\theta }\left(p_i\right)=\left\{n_i^r|r=1\text{,}\cdots ,m\right\}$ 。

2.1.2. 构造法向集

有了候选法向集$N_{\theta }\left(p_i\right)$ 后，在此基础上构造候选平面集${\Theta }_i=\left\{{\theta }_i^r|r=1\text{,}\cdots ,m\right\}$ ，其中候选平面${\theta }_i^r$ 是经过邻域$N_{p_i}$ 的几何中心点$\overline{p_i}$ 并以$n_i^r$ 为法向的平面。对于候选平面${\theta }_i^r$ ，计算邻域点$p_i^j$ 到该平面的距离$d_i^{\left(j\right)r}$ ，$j=1,\cdots ,K$ 。对距离${\left\{d_i^{\left(j\right)r}\right\}}_{j=1}^K$ 进行升序排序，得到有序序列${d_i^{\left(j\right)r}|}_{\left(1\right)}\le {d_i^{\left(j\right)r}|}_{\left(2\right)}\le \cdots \le {d_i^{\left(j\right)r}|}_{\left(K\right)}$ ，选取序列中10%至90%的距离值并计算其平均值$d_i^r$ ，将$d_i^r$ 记作邻域$N_{p_i}$ 到候选平面${\theta }_i^r$ 的截断距离。

计算候选平面集中每个平面的截断距离，得到截断距离集合$D_i=\left\{d_i^r|r=1,\cdots ,m\right\}$ 。选择截断距离最小的$k_{top}$ 个候选平面，这$k_{top}$ 个候选平面的法向构成了当前点$p_i$ 的法向集$N\left(p_i\right)=\left\{n_i^s|s=1,\cdots ,k_{top}\right\}$ 。

2.2. 计算法向集差异

利用2.1的过程为点云中每一点生成法向集。当前点$p_i$ 法向集为$N\left(p_i\right)=\left\{n_i^s|s=1,\cdots ,k_{top}\right\}$ ，其邻域点$p_i^j$ 的法向集为$N\left(p_i^j\right)=\left\{n_i^j{}^t|t=1,\cdots ,k_{top}\right\}$ ，这两个集合中的任意两个法向都可以计算差异：

${\gamma }_{i,j}^{s,t}=\min \left\{{\Vert n_i^s-n_i^j{}^t\Vert }_F^2,{\Vert n_i^s+n_i^j{}^t\Vert }_F^2\right\},n_i^s\in N\left(p_i\right),n_i^j{}^t\in N\left(p_i^j\right),$ (2)

这样得到当前点$p_i$ 与其邻域点$p_i^j$ 法向集的差异集合${\Gamma }_{i,j}=\left\{{\gamma }_{i,j}^{s,t}|s=1,\cdots ,k_{top};t=1,\cdots ,k_{top}\right\}$ ，取该集合中的最小值为当前点$p_i$ 与其邻域点$p_i^j$ 法向集差异，记为${\eta }_i^j$ 。

2.3. 赋权拟合

根据法向集差异${\eta }_i^j$ 将邻域点$p_i^j$ 赋权为：

$w_i^j={\text{e}}^{-\frac{{\left({\eta }_i^j\right)}^2}{{\sigma }^2}},$ (3)

其中，$\sigma =\frac{1}{2}\times \frac{1}{K\cdot k_{top}^2}{\displaystyle \sum _{j=1}^K{\displaystyle \sum _{{\gamma }_{i,j}^{s,t}\in {\Gamma }_{i,j}}{\gamma }_{i,j}^{s,t}}}$ 。

计算当前点$p_i$ 的每个邻域点$p_i^j$ 的权重$w_i^j$ 后，再利用加权最小二乘拟合得到最终法向$n_i$ 。首先，计算加权的邻域点几何中心$\overline{p_i^w}=\frac{{\displaystyle \sum _{j=1}^K{w}_i^j\cdot p_i^j}}{{\displaystyle \sum _{j=1}^K{w}_i^j}}$ ，然后，构建加权协方差矩阵$C_i^w$ ：

$C_i^w=\frac{{\displaystyle \sum _{j=1}^K{w}_i^j\left(p_i^j-\overline{p_i^w}\right)}\cdot {\left(p_i^j-\overline{p_i^w}\right)}^{\text{T}}}{{\displaystyle \sum _{j=1}^K{w}_i^j}},$(4)

最后，加权协方差矩阵$C_i^w$ 最小特征值对应的特征向量为最终所求的法向$n_i$ 。

3. 实验结果及分析

为了评估算法的性能，本文算法与一些具有代表性的算法进行了比较，包括PCA，RNE，LRRfast，PCV。

同时，为了更好地评估运用法向集来刻画局部邻域结构的优势，本文还构造了一个基准算法(法向加权拟合算法)，该算法运用单一的、不准确的初始法向来度量法向差异。首先，利用PCA算法得到点云$P$

中每一点的初始法向。然后，计算当前点$p_i$ 与其邻域点$p_i^j$ 的法向差异${\tilde{\gamma }}_i^j=\min \left\{{\Vert {\tilde{n}}_i-{\tilde{n}}_i^j\Vert }_F^2,{\Vert {\tilde{n}}_i+{\tilde{n}}_i^j\Vert }_F^2\right\}$ ，其中，${\tilde{n}}_i$ 为点$p_i$ 的初始法向，${\tilde{n}}_i^j$ 为其邻域点$p_i^j$ 的初始法向。最后利用法向差异对$p_i^j$ 进行赋权：${\tilde{w}}_i^j={\text{e}}^{-\frac{{\left({\tilde{\gamma }}_i^j\right)}^2}{{\sigma }^2}}$，其中，$\sigma =\frac{1}{2}\times \frac{1}{K}{\displaystyle \sum _{j=1}^K{\tilde{\gamma }}_i^j}$ ，并利用加权最小二乘拟合得到所求法向$n_i$ 。

本算法在anchor (40 K)、fandisk (77 K)、mechpart (92 K)、octahedron (50 K)四个模型上进行了测试。所有模型都添加了均值为0，标准差为点云之间平均距离的30%，40%，50%，60%的高斯噪声。

图3展示了在anchor (40 K)模型中选取的同一点，即红色点的近邻点的权重大小，左侧为本算法，右侧为法向加权拟合算法。

Figure 3. The weight magnitudes of the neighboring points of the same point

图3. 同一点的近邻点的权重大小

3.1. 参数的确定

本文算法主要的参数如下：

$K$ ：当前点的邻域大小。

$m$ ：候选法向的数量。

$k_{top}$ ：法向集中的法向数量。

在实验中，$K$ 分别设置为50、100、200、300，$m$ 设置为50，$k_{top}$ 设置为10。其余的算法邻域值$K$ 分别设置为50、100、200、300，其余相关参数采用默认值。

3.2. 评价函数

为了对计算结果进行量化评估，本文使用带有阈值的均方根误差($RMS\_\tau$ )，定义如下：

$RMS\_\tau =\sqrt{\frac{1}{\left|C\right|}{\displaystyle \sum _{P\in C}{v}_P^2}},$ (5)

其中，$v_P=\{\begin{array}{l}\langle {\text{n}}_{P,ref},{\text{n}}_{P,est}\rangle ,若\langle {\text{n}}_{P,ref},{\text{n}}_{P,est}\rangle <\tau \\ \frac{\pi }{2},否则\end{array}$ ，$\langle {\text{n}}_{P,ref},{\text{n}}_{P,est}\rangle$ 是法向${\text{n}}_{P,ref}$ 与${\text{n}}_{P,est}$ 之间的夹角，${\text{n}}_{P,ref}$ 是点$P$ 的真实法向，${\text{n}}_{P,est}$ 是点$P$ 的估计法向。本文中取$\tau ={10}^{\circ }$ 。

3.3. 结果评估

表1展示了四种邻域尺度下的不同算法在四种噪声尺度模型下$RMS\_\tau$ 的均值比较，其中粗体表示最优值。图4进一步展示了各算法在四种模型下的$RMS\_\tau$ 值的比较。图5为本文算法、法向加权拟合算法、PCA算法在噪声水平均为30%的四种模型下的可视化呈现，左侧为本算法，中间为法向加权拟合算法，右侧为PCA算法。为进行可视化比较，本算法与法向加权拟合算法均选取了$RMS\_\tau$ 值最小时的邻域大小，而PCA算法则采用了与法向加权拟合算法相同的邻域大小。下文将基于实验结果，对不同算法的表现进行综合比较与讨论。

如表1与图4所示，在四种不同噪声尺度模型与多种邻域尺度设置下，与法向加权拟合、PCA、RNE算法相比，本算法能够得到更低的$RMS\_\tau$ 值，展现出良好的整体稳定性与对噪声的鲁棒性。在部分模型(如anchor (40 K))、部分邻域尺度(如$K=50$ )情况下，本算法取得了与LRRfast、PCV算法相近甚至略优的结果，说明本算法在特定条件下具备一定的有效性。

Table 1. Comparison of the mean values of $RMS\_\tau$ under four noise scale models for different algorithms at four neighborhood scales (bold font indicates the optimal value)

表1. 四种邻域尺度下的不同算法在四种噪声尺度模型下$RMS\_\tau$ 的均值比较(粗体表示最优值)

模型	邻域大小	本算法	法向加权拟合	PCA	RNE	LRRfast	PCV
anchor (40 K)	50	0.3636	0.5784	0.7132	0.9690	0.5218	0.3855
	100	0.2759	0.3929	0.8583	0.6148	0.3499	0.2977
	200	0.3394	0.4050	1.0222	0.5645	0.4436	0.3925
	300	0.4491	0.5209	1.1249	0.5689	0.5276	0.4427
fandisk (77 K)	50	0.3281	0.6063	0.6324	1.0080	0.5150	0.4228
	100	0.2145	0.3485	0.7364	0.4645	0.2496	0.2206
	200	0.2427	0.3146	0.8668	0.3722	0.2048	0.2013
	300	0.2853	0.3616	0.9505	0.3553	0.2659	0.2384
mechpart (92 K)	50	0.4836	0.7219	0.6658	0.9695	0.6066	0.5559
	100	0.3444	0.4890	0.7758	0.6148	0.3717	0.3572
	200	0.3514	0.4087	0.9143	0.5005	0.3216	0.3049
	300	0.3883	0.4111	1.0083	0.4652	0.3508	0.3085
octahedron (50 K)	50	0.4332	0.7775	0.5546	1.1195	0.6067	0.5691
	100	0.1981	0.3974	0.6432	0.4832	0.2719	0.2497
	200	0.1643	0.3148	0.7536	0.3534	0.1214	0.0962
	300	0.1725	0.3449	0.8261	0.3047	0.1262	0.0698
总平均		0.3146	0.4621	0.8154	0.6080	0.3659	0.3195

Figure 4. Comparison of $RMS\_\tau$ of four models under different algorithms

图4. 不同算法在四种模型下的$RMS\_\tau$ 值的比较

如图5可视化结果所示，在噪声的影响下，PCA算法仍然可以准确地估计出光滑曲面的法向，但是在模型的尖锐特征处估计出的法向准确性较低。这是由于在尖锐特征处，特征点的邻域点可能采样自多个光滑曲面，所以PCA算法选择整个邻域进行平面拟合时，易受采样自其他光滑曲面的点影响，导致估计出的法向过于光滑，误差较大。这样，法向加权拟合算法利用PCA算法估计的单一的、不准确的初始法向来加权拟合，在尖锐特征处得到的法向估计结果不够准确。而本算法在不同噪声环境下均能有效抑制跨曲面干扰，保持了尖锐区域的几何特征，与法向加权拟合算法相比，本算法提高了在尖锐特征处的法向估计的准确性。

Figure 5. The visualizations of this algorithm, the normal-weighted fitting algorithm and the PCA algorithm under different models

图5. 本算法、法向加权拟合算法、PCA算法在不同模型下的可视化

Figure 6. The display of the normal sets and weights of the selected edge point and its two neighboring points in the octahedron (50 K) model

图6. 本算法在octahedron (50 K)模型中选取的边点及其两邻域点的法向集与权重展示

本文方法依靠当前点与其邻域点之间的法向集差异进行加权，以此判断当前点与其邻域点是否属于同一曲面结构。在octahedron (50 K)模型中，棱边两侧平面的夹角较大，几何结构表现为弱特征，导致在局部邻域内难以有效区分尖锐特征。如图6所示，当前点，即红色点与其位于不同平面的邻域点之间的法向集高度接近，使得不同平面上的邻域点被赋予相近的权重，在加权拟合过程中难以有效区分不同的曲面结构，最终使尖锐特征处的估计结果出现过光滑现象。

4. 结论

本文提出了一种基于法向集差异的法向估计算法。本算法通过对单一的、不准确的初始法向设计合适的扰动筛选机制，为每一点生成法向集，并利用法向集差异对邻域点进行赋权筛选。实验结果表明，与法向加权拟合算法相比，本算法在尖锐特征处可以得到较好的法向估计结果，可行性较强。然而，从整体准确性来看，由于法向集构造的局限性，本算法对弱特征的感知不够敏感，在弱特征模型下的鲁棒性不及LRRfast与PCV算法，比如在octahedron (50 K)模型上误差明显偏高。在未来工作中，我们会构造更合理的候选法向集的构造方式，更好地刻画邻域特征，同时可以利用本算法对其他利用单一的、不准确的初始法向度量法向差异的算法进行改进，如LRRfast、PCV算法。

NOTES

^*通讯作者。