压缩感知点云数据压缩

张习民, 余小清， 万旺根， 张 娟

1.上海大学通信与信息工程学院，上海200444

2.上海大学智慧城市研究院，上海200444

3.河南教育学院物理与电子工程学院，郑州450046

三维激光扫描技术以其非接触性、高精度等特点广泛应用于空间信息的获取，但扫描得到的海量点云数据不仅占用了大量的存储空间，而且不利于传输及后期运算处理.因此，在保证物体的形状精度和数据处理效率的前提下有必要减少点云的数据量.目前，专家学者对点云精简及压缩方面的研究主要集中于网格压缩法和基于空间点模型的压缩法.网格压缩法是一种针对网格拓扑连接信息和几何信息的编码方法，可分为单一位率压缩法（single-rate compression）[1]和渐进压缩法（progressive compression）[2].在单一位率压缩中，以网格的拓扑连接信息和几何数据作为一个整体进行编码压缩、传输和解码；渐进式压缩允许多层次细节（level of details,LODs）的传输和重建.网格编码压缩从早期以拓扑连接信息为主的连接信息驱动压缩算法到后来的几何信息驱动型的渐进网格编码算法[3-4]，需要额外对扫描点进行网格化和拓扑连接信息的处理，运算效率不高.基于点模型的压缩方法节省了网格化的开销，处理和渲染更加方便，更易于体现复杂的拓扑结构，加之3D激光扫描设备的大量应用，引起了广泛关注.其中，文献[5]针对体和面的点云数据，利用八叉树进行空间分割，并结合行程编码（run length coding,RLC）和霍夫曼编码（Huffman coding,HC）对点云模型进行压缩编码.文献[6]给出了一种单一位率的压缩编码方法，通过建立预测树对点云模型进行预测和熵编码.文献[7]在点云数据八叉树结构基础上引入双缓存处理技术，实现对点云的差值编码，获得了较高的压缩率，但基于空间点的压缩方法若要保留其显著的特征点，仍需要对点的法矢量进行隐式求解，算法依然复杂.Daribo等[8-11]提出了一种曲线基的点云压缩方法，将静态网格投影到物体的表面，通过相机捕获投影网格点信息，将空间点以曲线基的方式进行表达，应用竞争的预测编码机制对单位率的点云进行压缩.该方法提高了压缩中的率失真性能，但通用性不大.

本文根据激光扫描的特点，借鉴Daribo曲线基的思路将扫描点换为二维矩阵，创新性地利用压缩感知原理直接对点云的坐标位置信息进行压缩,既不必考虑点云的拓扑关系，又不必重建隐式特征信息，从而实现了对点云的保真压缩.

1 基于压缩感知（CS）的点云数据压缩

1.1 压缩感知的原理介绍

压缩感知理论认为：信号只要在某一正交空间具有稀疏性，就能以较低的频率采样，而以较高概率重构该信号[12-14].运用压缩感知原理对信号进行压缩采样的前提是信号必须稀疏，但并非所有的信号都满足该条件，不过大多信号经过适当的变换是可以满足条件的.因此，压缩感知的第一步工作就是对原始信号进行某种正交稀疏变换.首先，考虑一维空间RN中一个有限长度数字信号x，可看成一个N×1的列向量，于是用一组正交基进行稀疏变换

式中，ψ为一个正交基矩阵[ψ1,ψ2,···,ψN].x为空间域，s表示在基ψ下的变换域表示，是N×1的列向量，显然s和x都能代表信号.xψ如果变换后的s含有k个较大的非零值，则可以说信号x是k稀疏的.

其次，用一个观测矩阵φ对信号x进行观测，得到观测值y

式中，φ为一个M×N矩阵.变换后的信号由RN空间压缩到RM空间（M＜N），实现了对原始信号的压缩处理.需要强调的是：采样和压缩在这个过程中是同时进行的.

观测矩阵的选择也非常重要，应满足约束等容条件（restricted isometric property,RIP）

约束等容的等价描述如下：测量矩阵φ和稀疏正交化矩阵ψ的基不相关，能够满足约束等容性条件的测量矩阵可以是独立同分布的高斯矩阵，也可以是伯努利分布的矩阵，且满足

式中，c是一个较小的常量[15].

1.2 信号的重构

若已知测量矢量y、观测矩阵φ、稀疏化正交矩阵ψ，希望重构出k稀疏信号s，即求解如下欠定方程组：

该方程可转化为l0范数问题

式（6）的求解是l0范数下的优化问题，是一个NP问题，可以转换为l1范数的求解问题.

该优化算法是一个凸优化的过程，可以转化为线性规划的问题，可以通过基追踪（basic pursuit,BP）算法、匹配追踪（match pursuit,MP）等算法来完成，本文采用正交匹配追踪（orthogonal matching pursuit,OMP）算法.

2 激光扫描点云的降维处理

图1 点云降维框图Figure 1 Dimensions reduction principle of the point cloud data

分析FARO激光扫描仪获取点云数据发现：除了在扫描跳变点外，连续相邻数据间的差异很小，有很强的相关性，而这种相关性正是进行数据压缩的前提.如果能根据数据特点对点云进行重新组合，相关性就会更强.由此，可对激光扫描仪所获取的点云进行降维处理，处理的过程见图1.令原始扫描点云的点集为P={P1,P2,···,Pi,Pi+1,···,Pn}，其中Pi(x,y,z,r,g,b)为其中一个扫描点，包含浮点型的位置坐标分量（x,y,z）和整型的色彩信息分量（r,g,b）.为统一数据类型，可对点集中的数据进行截短处理，截掉色彩分量，只保留浮点型的位置坐标信息分量（x,y,z）.这样，截短后一个点可表示为Pi(x,y,z).需要注意的是在截短的数据中，同一点的3个不同分量相关性较小.考虑相邻点的相同坐标分量有更强相关性，为获得更好的压缩效果，可对截短的点云数据重排序，抽取点云中相同坐标分量x（y或z）排为一列，然后截取连续2n个坐标分量数据作为一个列矢量，最后形成大小2nN的矩阵X来存贮新生成的点云数据，本文取2n=128.

3 激光扫描点云的CS压缩

降维处理后点云矩阵X（128×N）的相关性加强，但进行CS压缩还需稀疏化处理.本文选用简单且行之有效的Harr小波基，并用MATLAB语言编写生成三级变换矩阵T（128×128），实现对矩阵X列向量的稀疏处理，变换后的矩阵可表示为

变换后矩阵列向量行号较小的部分相当于低频分量，行号较大的部分相当于高频分量.

接着进行采样，矩阵X经过三级Harr变换后较大系数值（即k值）约在16～20之间，根据式（4）将M取3～4倍的k值（本文中取M为80），生成该矩阵的MATLAB语句为：“φ=rand(M,N)”，M和N即为采样矩阵的行和列.该压缩采样可用公式表示为

式中，Y为压缩后的数据.

4 压缩点云数据的重建与统计滤波

4.1 压缩点云数据的重建

若已知采样值Y和测量矩阵φ，数据恢复采用OMP算法.恢复X1中第n列的步骤如下：

步骤1 令残差向量r=Yn，Yn代表测量值中的第n列；

步骤2 用测量矩阵的各列与r进行内积运算〈r,φi〉，i=1,···,128，找出内积最大的一列φl,l为列号；

步骤3 用最小二乘法获得X1中第n列第l位置上元素值：

步骤4 修正残差值

步骤5 在观测矩阵φ中剔除掉φl列，返回步骤2，确定X1第n列下一个元素值及其位置，直到待定出全部的k个元素，该列的其他位置置0.

这样就完成了对X1矩阵一列数据的恢复，同理可恢复其他各列.整个X1矩阵恢复后，可以通过Harr逆变换得到点云矩阵X.最后通过矩阵的逆排序进行原始点云的恢复.

4.2 点云重建统计滤波

在CS进行稀疏化处理过程中，如果相邻点位置坐标跳变，或在不同扫描线过渡中存在较大坐标差异，都会引起重建后点云数据的稀疏外点噪声干扰，因此有必要对重建后的点云数据进行统计滤波.统计滤波的原理如下：首先计算点云中每个点的k邻域均值距离µ和标准偏差σ，并认为其全局统计结果服从Gauss分布.滤波时将那些落入µ±λσ之外的点作为噪声外点予以滤除，其中λ为倍乘因子.该统计滤波用c语言编程实现.

5 实验结果及分析

实验中数据源为Faro激光扫描仪所获得的点云，采用MATLAB和VC混合编程的方法，PC配置为Intel®CoreTMi5-2400 CPU＠3.10 GHz,4 G内存，Win7操作系统.

首先将一个模型猴子选作压缩处理对象，预处理后获得的原始点云数据见图2（a），包含700 672个点，压缩复原后的点云见图2（b），点数为437 920.在压缩算法中，高斯矩阵采样的比率为80：128，采样压缩比为1.6.为获得压缩的客观评价，本文求出压缩前后各对应点的距离差值平方，求和后再求平均.计算压缩前后的误差公式为

式中，Xi为原始点云数据中一个点的列向量，Ri为重构后点云数据中一个点的列向量，N为点云数据点的总数量.运算后得到的误差为3.4398×10-7m2，可见误差较小.

图2 点云压缩比较Figure 2 Comparison of point cloud data before and after compression

为了验证算法的通用性，换用大场景数据和不同类型的点云进行测试，压缩所对应的数据见表1，实验结果见图3～5.

表1 不同点云压缩对比Table 1 Compression comparison of different point cloud data

从表1中可以看出，本算法压缩比都在2倍以上，压缩后较好地保持了原始点云数据的特征信息，且整体的误差较小，如图2、图4、图5所示.从图2和5中还可以得到：本算法不仅能较好地实现对较平坦的点云的压缩（见图2），而且能够实现对曲率变化较大的点云数据的压缩（见图5），说明算法具有一定的通用性.不足之处是：虽然对密度较大（对应的文件也较大）的点云的压缩效果较好（见图4和5）；密度较小点云则存在部分细节的丢失，见图3复原图方框部分，且误差也较大（0.039 6 m2），这主要是由于激光扫描获得的原始点云密度较小的缘故.

本文中的压缩率可以通过调节采样矩阵的大小进行调节，但一定要满足约束等容条件，否则会使复原后的点云质量下降，文中的采样矩阵为80×128.

图3 上大图书馆点云压缩比较Figure 3 Comparison of SHU library point cloud data before and after compression

图4 书本点云数据压缩比较Figure 4 Comparison of a book point cloud data before and after compression

图5 纸杯点云数据压缩比较Figure 5 Comparison of a paper bottle point cloud data before and after compression

当然，为了去除压缩恢复后点云外部噪声点，可用前面提及的统计滤波对猴子点云数据进行处理，其中滤波均值k=30，标准偏差σ=1，对比结果见图6，左半部分滤波前的点云数据为700 762个，右半部分滤波后的点云数据为619 292，滤除的点数为81 470.由图6中（a）和（b）所标记的绿色方框显然可以看到，滤波前后噪声外点得到明显的控制，从而使恢复后的点云边缘部分更加清晰.

图6 统计滤波前后对比图Figure 6 Comparison of point cloud data after and before statistical f iltering

6 结语

本文结合激光扫描点云数据的特点，创新性地将CS理论用于点云压缩，该算法不需要点云数据的拓扑关系，只需根据扫描数据的先后顺序，直接对数据坐标点进行采样压缩.通过实验对比，可以发现压缩效果良好.但也存在一些需要改进的地方：随着压缩比的增大，噪声点明显增多，于是可以考虑采用更优的采样矩阵及恢复算法加以克服；压缩算法中未考虑彩色信息，这也是下一步需要解决的问题.

[1]ALLIEZ P,DESBRUNM.Valence-driven connectivity encoding for 3D meshes[C]//European Association for Computer Graphics,2001:480-489.

[2]ALLIEZ P,DESBRUN M.Progressive encoding for lossless transmission of triangle meshes[C]//Association for Computing Machinery's Special Interest Group on Computer Graphics and Interactive Techniques,2001:198-205.

[3]PENG J, KUO C C J.Progressive geometry encoder using octree-based space partitioning[C]//Proceedings of the 2004 IEEE International Conference on Multimedia and Exposition,2004:1-4.

[4]PENG J,KUO C C J.Geometry-guided progressive lossless 3D mesh coding with octree(OT)decomposition[C]//Association for Computing Machinery's Special Interest Group on Computer Graphics and Interactive Techniques,2005:609-616.

[5]SIDDIQUI R A,CELASUN I,OCTREE U B.Based compression of volumetric and surface 3D point cloud data[C]//Proceedings of 13th Intl Conference on Virtual Systems and Multimedia,Brisbane,Australia,2007:278-282.

[6]GUMHOLD S,KAMI Z,ISENBURG M,SEIDEL H P.Predictive point-cloud compression[C]//Proceedings of the sixth Israel-Korea Binational Conference,2005:125-129.

[7]KAMMERL J,BLODOw N,RUSU R B,GEDIKLI S,BEETZ M,STEINBACH E.Real-time compression of point cloud streams[C]//Robotics and Automation(ICRA),Saint Pual,2012:778-785.

[8]DARIBO I,FURUKAwA R,SAGAwA R,KAwASAKI H,HIURA S,ASADA N.Point cloud compression for grid-pattern-based 3D Scanning System[C]//Visual Communications and Image Processing(VCIP),2011:1-4.

[9]DARIBO I,FURUKAwA R,SAGAwA R,KAwASAKI H.Dynamic compression of curve-based point cloud[C]//SIVT'11 Proceedings of the 5th Pacif ic Rim conference on Advances in Image and Video Technology-Volume Part II,2011:323-334.

[10]DARIBO I,FURUKAwA R,SAGAwA R,KAwASA K H,HIURAS,ASADAN.Efficient rate-distortion compression of dynamic point cloud for grid-pattern-based 3D scanning systems[C]//3D Research,2012,3(1):1-9.

[11]DARIBOI,FURUKAwA R,SAGAwA R,KAwASAKIH,HIURA S,ASADA N.Curve-based representation of point cloud for efficient compression[C]//The 14th Meeting on Image Recognition and Understanding,2011:1385-1390.

[12]DONOHO D．Compressed sensing[J]．IEEE Transactions on Information Theory,2006,52(4):1289-1306．

[13]CANDE S E,EMAMNUEL J C.Compressive sampling[C]//Proceedings of International Congress of Mathematicians Switzerland:European Mathematical Society Publishing House,2006:1433-1452．

[14]CANDE'S E,ROMBERG J,TAO T．Robust uncertainty principles:exact signal reconstruction from highly incomplete frequency information[J]．IEEE Transactions on Information Theory,2006,52(2):489-509．

[15]BARANIUK R．A lecture on compressive sensing[J].Signal Processing Magazine,2007,24(4):118-121．