利用双目立体视觉技术测量降雨诱发滑坡三维变形的方法

1.本发明属于计算机视觉和摄影测量领域，涉及一种滑坡体变形测量方法，具体涉及一种利用双目立体视觉技术测量降雨诱发滑坡三维变形的方法。


背景技术：

2.降雨是滑坡发生最重要的触发因素和动力来源，降雨对边坡岩土体起加载作用，增加土体自重并产生动静水压力，同时降雨会侵蚀坡脚，改变边坡结构，在我国很多滑坡具有大雨大滑、小雨小滑、无雨不滑的特点。在降雨和渗透的作用下，斜坡平衡遭受破坏而产生滑动的现象，包括缓慢长期的斜坡变形和突然急剧的爆发过程。相较于一般滑坡，降雨诱发的滑坡，分布最广，发生频率最高(占滑坡总数的90％)，集中发生在雨季及暴雨之后，长时间持续降雨还有可能会引发二次滑坡，测量难度更高，因此对降雨诱发滑坡区的无接触变形测量工作就显得极为重要。
3.传统滑坡体变形测量方法主要是大地测量法。大地测量法以人工监测为主要手段，通过全站仪导线和电磁波测距三角高程方法进行变形监测。大地测量法虽然积累了大量的经验，但其监测频率较低，监测空档期长，需要人员到现场观测，工作量大，在降雨条件下，测量人员安全也难以保证。近年来，随着科技的发展，一些新的方法不断用于滑坡检测技术，如gnss技术、三维激光扫描、无人机数字摄影测量、insar技术和botdr技术，gnss技术通过接收空中卫星信号测距进行定位，优点是不受通视条件制约精度高、不受天气干扰、可以连续长时间进行实时监测，但其设备造价高；三维激光扫描技术是继gnss技术之后发展起来的一门新兴测绘技术，该技术利用激光测距原理，通过记录滑坡表面点的三维坐标实现滑坡体变形测量，优点是无需事先埋设监测点，实现了无接触测量，精度高，可以穿透植被，但其设备造价较高；数字摄影测量技术通过对被监测区域进行摄影，用摄影测量软件系统获取被监测区域的三维表面点云模型，来实现滑坡体的变形测量，该技术数据获取效率高、运行成本低、响应速度快，但精度较低；insar技术是以合成孔径雷达的相位信息为依据获取地表三维信息的一项技术，通过单轨模式或重复轨道模式获取地面上同一地区的复图像对，由复杂图像对上的干涉条纹图获得精确位置信息，实现滑坡体变形测量，具冇成本低、无需建立监测网、覆盖范围大、时间短、分辨率高等优势，但受环境影响大；botdr技术是一种新型的分布式光纤传感技术，通过光纤的微弯损耗实现滑坡体的变形测量，具有高灵敏度、适用范围广的优点，但光纤布设在滑坡体表面，易断裂。立体视觉技术在滑坡体变形测量中也有一定应用，但其测量精度较低，而且受环境影响大，难以保证在暴雨天气下工作，不适用于降雨诱发滑坡体的变形测量。


技术实现要素：

4.针对传统和现代滑坡体变形测量方法的缺点，本发明给出了一种利用双目立体视觉技术测量降雨诱发滑坡三维变形的方法。该测量方法成本低廉，简单方便，精度高，时间快，可以连续测量，适用于条件恶劣的作业环境。
5.本发明的目的是通过以下技术方案实现的：
6.一种利用双目立体视觉技术测量降雨诱发滑坡三维变形的方法，包括如下步骤：
7.步骤一、构建双目立体视觉系统并标定
8.步骤一一、构建双目立体视觉系统对测量滑坡区域进行信息采集，建立摄像机坐标系，其中：
9.双目立体视觉系统由左ccd摄像机和右ccd摄像机组成；
10.摄像机坐标系为左ccd摄像机坐标系和右ccd摄像机坐标系；
11.步骤一二、对双目立体视觉系统进行标定，得到双目立体视觉系统的标定参数，其中：
12.标定方法采用张正友的棋盘格标定法，分别对左ccd摄像机和右ccd摄像机标定，求取摄像机内参数矩阵、畸变变量以及两台摄像机相对于标定物的旋转矩阵和平移向量畸变变量；根据两台ccd 摄像机之间的位置关系，通过计算得到旋转矩阵和平移矩阵；通过立体校正把实际中非共面行对准的两幅图像，校正成共面行对准(共面行对准：两摄像机图像平面在同一平面上，且同一点投影到两个摄像机图像平面时，应该在两个像素坐标系的同一行)，将实际的双目系统校正为理想的双目系统；
13.标定参数包括左ccd摄像机坐标系和右ccd摄像机坐标系的内参数和外参数、旋转矩阵和平移向量；
14.步骤二、图像采集与图像处理
15.步骤二一、在被测点处进行标志布设，标志点采用同心圆形，通过同心圆提供多几何约束进行圆心定位，其中：
16.标志点的设计宜根据具体的工程情况进行考虑，选择与环境对比度高的回光反射材料制作而成；
17.步骤二二、对标志点布设区域进行图像采集和图像处理，得到清晰的标志图像，其中：
18.图像采集的环境可能为阴天或夜晚，因此在阳光不足时会在图像采集的同时开启闪光灯；
19.降雨诱发滑坡的图像采集环境多为雨天，由于降雨诱发滑坡的发生环境存在雨痕和雨雾等干扰因素，导致成像质量低，影响变形测量准确性，需要对采集到的原始图像进行去雨痕和雨雾处理，得到清晰的标志图像；
20.步骤三、标志点的识别和定位
21.步骤三一、对清晰的标志图像进行滤波、去噪和二值化，得到圆形标识点目标参数；
22.步骤三二、对单个圆形标志点采用周长、面积和椭圆度三个形状准则进行约束，去除干扰因素，进行椭圆拟合，得到单个圆形标志点的圆心精确坐标，最后通过聚类得到圆心像素坐标代表值；
23.步骤三三、对各标志点进行双目重建，得到不同时间段各标志点于世界坐标系的三维坐标，通过数据融合，定量获得滑坡测点在一段时间内的水平位移、垂直位移和三维空间位移，实现降雨诱发滑坡三维变形的测量，进而预测滑坡变形发展趋势，其中：
24.通过canny算法进行边缘点粗定位，使用zernike矩算法进行边缘的亚像素定位，
以提高定位精度；通过最小二乘法拟合椭圆方程得到各同心圆圆心坐标；最后通过中心聚类算法，获得圆心像素坐标的代表值，即可得到标志点在世界坐标系内的三维坐标。
25.相比于现有技术，本发明具有如下优点：
26.本发明属于间接测量和非接触测量，可以确定测点的三维空间数据，并通过长时间连续监测得到测点的水平位移、垂直位移和三维空间位移，能够很好的用于降雨诱发滑坡体变形过程的定量化研究。
附图说明
27.图1为利用双目立体视觉技术测量降雨诱发滑坡三维变形的流程图。
具体实施方式
28.下面结合附图对本发明的技术方案作进一步的说明，但并不局限于此，凡是对本发明技术方案进行修改或者等同替换，而不脱离本发明技术方案的精神和范围，均应涵盖在本发明的保护范围中。
29.本发明提供了一种利用双目立体视觉技术测量降雨诱发滑坡三维变形的系统，所述系统包括标志点、标志点支架、光源、左ccd 摄像机、右ccd摄像机和计算机，其中：
30.所述标志点由反光板和反光膜组成；
31.所述标志点支架采用钢结构，根据现场情况在测点处进行架设，用以支撑标志点；
32.所述左ccd摄像机、右ccd摄像机均与计算机相连，用以对标志点区域进行摄像，并将图像传递给计算机供其处理；
33.所述计算机包括图像处理模块、标志图像识别模块、双目重建模块和数据融合模块；
34.所述图像处理模块由图像去雨雾处理模块和图像去雨痕处理模块组成，降雨诱发滑坡的拍摄环境多为雨天，雨水在下落过程中以雨痕的形式出现在图像中，在落地后以雨雾的形式出现在图像中，影响标志点识别的精度，通过图像去雨雾处理模块和图像去雨痕处理模块对雨环境下采集到的图像进行去雨雾和去雨痕处理，使其清晰；
35.所述标志图像识别模块与图像处理模块相连，对图像处理模块处理完成后的图像进行标志点识别，得到标志点圆心像素坐标信息；
36.所述双目重建项模块基于双目立体视觉原理，根据标志点圆心像素坐标通过计算获得测点的三维坐标；
37.所述数据融合模块与双目重建模块相连，得到测点的水平位移、垂直位移和三维空间位移，并将其以图表的形式显示给测量员。
38.一种利用双目立体视觉技术测量降雨诱发滑坡三维变形的智能方法，包括如下步骤：
39.1、构建双目立体视觉系统，并对双目立体视觉系统进行标定。
40.将左ccd摄像机和右ccd摄像机固定，尽量保持前向平行布置，根据测试要求选择不同焦距范围的光学镜头，确保标志点所在区域在摄像机的视野范围内。采用张正友的棋盘格标定法，分别对左ccd 摄像机和右ccd摄像机标定，求取ccd摄像机内参数矩阵、畸变变量以及两台ccd摄像机相对于标定物的旋转矩阵和平移向量畸变变量；根据两台ccd摄像机
之间的位置关系，得到旋转矩阵和平移矩阵；通过立体校正把实际中非共面行对准的两幅图像，校正成共面行对准(共面行对准：两摄像机图像平面在同一平面上，且同一点投影到两个摄像机图像平面时，应该在两个像素坐标系的同一行)，将实际的双目系统校正为理想的双目系统。
41.2、图像采集与图像处理
42.对标志点所在区域以一定频率进行连续拍摄，若环境亮度不足则需在摄像机拍摄的同时给出闪光，将采集到的图像上传至图像存储模块，降雨型滑坡的拍摄环境多存在雨痕和雨雾等干扰因素，还需通过图像去雾算法和图像去雨算法对图像进行一定的处理，突出图像细节，提升对比度，使之看起来更加清晰。
43.3、标志点的识别和定位
44.对采集到的标志点图像完成滤波、去噪和二值化，得到圆形标识点目标参数；之后通过canny算法进行边缘点粗定位，使用zernike 矩算法进行边缘的亚像素定位，以提高定位精度；通过最小二乘法拟合椭圆方程得到各同心圆圆心坐标；最后通过中心聚类算法，获得圆心坐标的代表值，即可得到标志点在世界坐标系内的三维坐标；记录不同时间段各标志点的三维坐标可以实现对滑坡体变形的测量。
45.各模块工作原理如下：
46.1、标志点
47.选用同心圆形标志图案，透视变形后为椭圆，易于识别；多个同心圆圆心坐标进行中心聚类得到的结果更加准确；易于从背景中分离，根据具体工程选用和环境对比度高的颜色，并覆以反光膜更加易于识别。
48.2、图像采集模块
49.2.1摄像机的成像过程
50.根据针孔相机模型，空间中任意一点p(x,y,z)在世界坐标系中的坐标值pw(xw,yw,zw)可变换为像素平面坐标系中的坐标值p(u,v)。
[0051][0052][0053][0054]
式中，r是从世界坐标系到相机坐标系的旋转矩阵，(r1、r2、r3)、 (r4、r5、r6)和(r7、r8、r9)为相机坐标系的基底向量在世界坐标系中的表示；t是从世界坐标系到相机坐标系的平移向量，t
x
、ty、tz在数值上与世界坐标系原点在相机坐标系中的坐标相等；zc是比例因子；b是照相机的固有参数矩阵。
[0055]
由于镜像畸变，实际图像点与理想图像点之间存在一定的偏差。考虑到径向畸变和切向畸变，本发明采用了包含四个参数的畸变模型，畸变后的归一化坐标值可以用pd(xd,yd)表示为：
[0056][0057]
式中，(xu,yu)是图像点在未畸变条件下的归一化坐标，(xd,yd) 是畸变条件下图像点的归一化坐标，归一化坐标就相当于在z的方向上当z＝1时用一个平面截断，这时光心与p点的连线在该面上的点即为p点的归一化点坐标；k1、k2为径向畸变系数；p1、p2为切向畸变系数。公式(1)中的矩阵b和k1、k2、p1、p2统称为相机的固有参数。
[0058]
像平面坐标系中点p(xd,yd)可通过像素坐标系上的坐标值p(u,v) 获得：
[0059][0060]
式中，dx、dy为每一个像素在x轴和y轴上的物理尺寸(mm)；u0为镜头光轴z与图像平面的交点o1在像素坐标系轴上的坐标值；v0为o1在像素坐标系v轴上的坐标值。
[0061]
2.2摄像机标定过程
[0062]
选择无畸变的针孔照相机模型作为计算模型，重建左右图像上棋盘格角点的图像坐标，采用最大似然估计法求解投影矩阵，利用旋转矩阵的正交性得到左右ccd摄像机的内在参数和外部参数。计算出 ccd摄像机的内参数和外参数初值。然后，考虑畸变的影响，将图像上n个角点的坐标值(ui,vi)利用畸变模型计算这n个角点的实际坐标值(ui,vi)，然后利用公式(6)进行数据拟合和迭代：
[0063][0064]
2.3立体视觉测量
[0065]
由于左右图像上对应点的匹配关系，可以根据两台ccd摄像机的几何约束，得到世界坐标系中p点的三维坐标，方程如下所示：
[0066][0067]
式中，f
l
和fr分别为左ccd摄像机和右ccd摄像机的有效焦距， (x
l
，y
l
)为左ccd摄像机像平面坐标系测点的坐标，(xr，yr)为右ccd摄像机像平面坐标系测点的坐标。
[0068]
3、图像处理模块
[0069]
雨滴颗粒较大，且下落速度快，在下落过程中，通常会以雨痕的形式出现在图像中，可能会遮挡标志图像，影响信息的识别；在雨滴落到地面后，会产生类似气溶胶的水雾，对光线造成折射的效果，同样影响标志信息的识别。
[0070]
本发明通过雨雾去除和雨痕去除两方面完成对于雨环境图像的处理，增强了降雨型滑坡体变形测量的准确性。
[0071]
3.1图像去雨雾算法
[0072]
雨雾环境会使采集到的图像模糊，影响目标定位精度。本发明采用暗通道先验去雾法对雾环境下拍摄的图像进行去雾化处理。
[0073]
图像去雨雾过程就是根据有雾图像i(x)求解无雾图像j(x)的过程。从公式(8)可以看出，基于物理模型的去雨雾算法本质是根据已知的有雨雾图像i(x)求出透射率t(x)和全局大气光成分a(a》0)，进而得到复原图像j(x)(场景光成分)。
[0074]
i(x)＝j(x)t(x)+a(1-t(x))
ꢀꢀꢀꢀ
(8)；
[0075]
式中，a、i、j三个向量是共面的，并且末端点共线，并可以推出t(x)的表达式：
[0076][0077][0078]
式中，ω(x)表示以x为中心的局部块，c为r、g、b三个通道， jc表示图像j的彩色通道图像。即先求出r、g、b通道中像素的最小值，然后将其存进和原始图像尺寸相同的灰度图中，再对该灰度图进行最小值滤波。
[0079]
若j是无雨雾图像，除了天空区域，式(10)的值趋近与0。
[0080]
在暗通道先验去雾中，可以通过于暗通道图像来从有雨雾图像中获取a值，其基本步骤是：首先从暗通道图像中获得前0.1％亮度大小的像素位置。其次在这些像素位置中，在有雨雾图像i(x)中寻找对应的具有最高亮度像素位置的亮度值，即可作为a值，式(8)可以推出：
[0081][0082]
假设局部块ω(x)是连续的，块的传输率为式(11)对每个通道去最小值运算，则有：
[0083][0084]
由于式(10)趋近于0，a》0，所以有：
[0085][0086]
把式(13)带入(12)，得到：
[0087][0088]
由于天空在无限远处且传输率趋于0，所以式(13)此时也是可以满足的。同时考虑到大气透视现象的存在，如果完全移除雨雾，图像看起来会不自然并且可能丢失深度信息，所以需要对远处的对象保留少量的雨雾，于是在式(14)中引入了恒参ω，一般取值为0.95，式(14)变为：
[0089][0090]
采用了软抠图算法来优化透射率，优化后的精确的透射率用t(x) 来表示，估算透射率与精确透射率的差值可表示为：
[0091][0092]
式中，λ时一个正则化参数，l为抠图拉普拉斯矩阵，l矩阵的第 (i,j)项定义为：
[0093][0094]
式中，ii和ij是输入图像在i和j处像素的颜色；δ
ij
为克罗内克函数，μk和∑k表示区域ωk内的颜色的均值矩阵和协方差矩阵；u3是一个3
×
3的矩阵；ε是一个正则化参数；|ωk|表示区域ωk内的像素个数。
[0095][0096]
得到透射率图t后，可以根据雾霾模型公式恢复出场景的亮度，但是，当t(x)接近0的时候，直接衰减项j(x)t(x)也是接近0，直接恢复j会容易受噪声干扰，因此需要对t做一个下界t0的限制，一般取0.1，则无雨雾图像可表示为：
[0097][0098]
3.2图像去雨痕算法
[0099]
本发明采用基于形态成分分析(mca)的图像分解方法来完成雨痕的去除，主要是利用了字典学习来去雨，其主要过程如下：
[0100]
首先对输入取雨雾图像j进行双边滤波，如公式(20)，将图像分为低频图像j
lf
和高频图像j
hf
，低频图像保留了原始图像的大致信息，这部分中不含雨，因此保留不变，而高频图像则包含了原始图像中的雨纹信息以及背景纹理信息；
[0101]
j＝j
lf
+j
hf
ꢀꢀꢀ
(20)
[0102]
对高频图像进行稀疏编码和字典学习，字典学习问题可以表述为公式(21)，这里的字典d
hf
包括了雨的字典d
hf-r
与图像边界的字典d
hf-g
，需要进一步区分；
[0103][0104]
式中，λ表示正则化参数，yh为第h个纹理块的训练样本， h＝1,2,
…
p，θh表示相对于d
hf
稀疏系数的向量，得到j
hf
基于d
hf
的稀疏表示。
[0105]
字典的划分，雨所在的块的梯度方向一般比较一致，图像边界所在的块表现出来的梯度就比较杂乱，对每个字典元素提取其梯度方向直方图作为特征，使用k-means进行聚类，方差比较大的则为图像边界的字典d
hf-g
，另外一类则为雨的字典d
hf-r
。
[0106]
通过公式(22)获得字典的稀疏系数，l表示非零系数θ
hhf
的最大数。将雨的字典所对应的稀疏系数设为0，利用稀疏系数和字典重构回图像，得到高频图像中的图像细节，再与低频图像相加，即为去雨过后的处理图像，如公式(23)。
[0107][0108][0109]
4、标志图像识别模块
[0110]
本发明通过提取圆形标识点圆心坐标实现特征的识别和定位，涉及圆形标识点的边缘检测、特征识别及中心拟合算法。
[0111]
首先对左右ccd摄像机采集的同步图像进行简单处理分析，完成滤波去噪和二值化，得到圆形标识点目标参数，包括形状参数、偏心率、球状度、圆形度以及边缘长度等，以便进行标志图像的边缘检测。
[0112]
利用canny算法进行初步的边缘点粗定位，所得到单一的边缘为下一步精确定位提供了方便；然后，使用zernike矩算法进行边缘的亚像素定位，以提高定位精度。两者的结合充分发挥了两者的优点，可快速精确地定位圆形标识点的边缘。
[0113]
圆形标志在相机成像平面上所呈现的是一个椭圆，本发明选择椭圆拟合法进行定位。利用最小二乘法拟合椭圆方程，根据获得的图像亚像素边缘进行椭圆拟合，从而定位每一个同心圆的中心坐标，为后续中心坐标聚类提供数据基础。
[0114]
椭圆方程为：
[0115]
ax2+bxy+cy2+dx+ey+f＝0
ꢀꢀꢀ
(24)。
[0116]
将椭圆的一般形式转换为矩阵形式：
[0117]
f(u,v)＝u*v＝1
ꢀꢀꢀ
(25)；
[0118]
式中，u＝[a b c de 1]，v＝[x
2 xy y
2 x y 1]t，在上述问题中引入约束条件||u||＝1，可建立目标函数：
[0119][0120]
进行最优化求解可得椭圆中心坐标：
[0121][0122]
根据上述边缘检测与拟合可获得多个同心圆圆心点，记为 x＝{x1,x2,x3,...,xn}，本发明将根据上述圆心利用k-means聚类的方法实现圆心代表值的提取。k-means算法的目标是将n个对象依据对象间的相似性聚集到指定的k个类簇中，每个对象属于且仅属于一个其到类簇中心距离最小的类簇中。对于k-means，首先需要初始化k个聚类中心{c1,c2,c3,...,ck}，1《k≤n，然后计算每一个对象到每一个聚类中心的欧式距离，如下：
[0123][0124]
式中，xi表示第i个对象1≤i≤n，cj表示第j个聚类中心的1≤j≤k， x
it
表示第i个对象的第t个属性，1≤t≤2，c
jt
表示第j个聚类中心的第 t个属性。依次比较每一个对象到每一个聚类中心的距离，将对象分配到距离最近的聚类中心的类簇中，得到k个类簇{s1,s2,s3,...,sk}。
[0125]
k-means算法用中心定义了类簇的原型，类簇中心就是类簇内所有对象在各个维度的均值，其计算公式如下：
[0126][0127]
式中，c
l
表示第l个聚类中心，1≤j≤k，|sl|表示第l个类簇中对象的个数，xi表示第l个类簇中第i个对象，1≤i≤|s
l
|。
[0128]
聚类中心坐标即为标志点圆心像素坐标代表值，通过双目重建可获得测点的三维坐标，记录不同时间段的三维坐标代表制，即可获得测点的水平位移和垂直位移，实现利用双目立体视觉技术测量降雨诱发滑坡三维变形的智能方法。