[发明专利]基于光流法的双目图像最大视差获取方法

摘要

基于光流法的双目图像最大视差获取方法，通过对左右摄像机采集到的图像进行一系列处理，得到图像上各点的运动矢量，对运动矢量进行分析判断最后得到最大视差，以便于进行下一步处理。由于光流法对弱纹理及重复纹理，噪声干扰等的影响较不敏感，对光流算法中数据项的梯度算子以及窗口权重函数的改进。考虑到最大视差唯一性的特点，分简单和复杂两次计算光流，有效减小计算量，加快了计算速度。在CLG‑TV金字塔光流法的基础上，对光流算法中平滑项引入全局扩散系数，对光流算法中数据项的梯度算子、窗口权重函数进行改进，使其更加适用于计算双目图像只有水平位移的情况，对金字塔计算方法，分简单和复杂两次计算光流，有效减小计算量，加快了计算速度。

主权项

基于光流法的双目图像最大视差获取方法，其特征在于：本方法采用了光流法来获取最大视差；通过对左右摄像机采集到的图像进行一系列处理，得到图像上各点的运动矢量，对运动矢量进行分析判断最后得到最大视差，以便于进行下一步处理；由于光流法对弱纹理及重复纹理，噪声干扰的影响较不敏感，得到较准确的最大视差值；具体而言，该方法的具体实施过程如下，1)由粗到精的多尺度计算：按照一定的采样标准，把原始图像分为多层，形成金字塔形式的图像结构，原始图像位于金字塔底部，而金字塔顶部则是图像分层后最低分辨率的近似；这样的多分辨率思想在使用光流分析大位移运动时特别适用，同时使用较小的取样因子使计算量减少，加快了计算速度；2)CLG‑TV光流法求解光流场；使用CLG‑TV光流法处理左右图像获取光流运动矢量；3)对光流算法中平滑项引入全局扩散系数，有效地在不连续处进行精确计算；4)对光流算法中数据项的梯度算子以及窗口权重函数的改进；使用一种更加适用于计算最大视差的梯度算子的计算方法，并通过改进的双边滤波器来计算窗口权重函数；5)考虑到最大视差唯一性的特点，分简单和复杂两次计算光流，有效减小计算量，加快了计算速度；实施流程如下，步骤1、使用左右摄像机采集到的双目图像序列输入计算机，首先对输入图像做预处理，如果源图像的分辨率较大先做下采样使图像缩小到VGA级别以减少计算量，由于本方法不使用图像的彩色信息，首先将彩色图转为灰度图，公式如下：I＝0.299R+0.587G+0.114B然后对图像进行高斯平滑操作来消除摄像头采集过程中产生的随机噪声；步骤2、使用CLG‑TV光流法求解光流场，得到初步估计；光流用来描述物体在观测到的图像中的二维运动，光流法运动估计的计算结果为图像中像素的运动矢量；理想情况下，对双目图像计算光流场将得到图像中像素点的水平运动矢量，这些水平运动矢量就等于立体匹配算法得到的各像素点的视差，光流法计算的运动矢量包含有视差信息，应用于双目图像最大视差的获取；光流法的基本假设将亮度约束和平滑约束相结合，既保证了运动矢量计算的准确性，又使得各区域内运动矢量能够平滑过渡，是一种高精度的运动估计方法；实际图像存在复杂多变的情况，同一景物在左右视点下的成像会有很大差异，会对获得精确的视差值造成困难；利用光流法求解最大视差的优势就在于光流法的精度能够达到亚像素级，完全可以胜任于最大视差的计算；光流法对传统方法中比较难以处理的深度不连续、弱纹理、重复纹理和噪声干扰问题的鲁棒性强；CLG‑TV光流法是结合全局与局部的光流法，基于全局的光流法对小运动精确度高，基于局部的光流法对大运动有更好的精度，将二者结合起来增强鲁棒性，CLG‑TV光流法使用的误差计算式将LK方法和HS方法的误差式结合起来；另外，HS和LK的误差计算式都采用了L2型范数；L2型范数有一些缺点：对噪声敏感度高，过度扩张，选择性差，以及不保留运动边界，为克服L2模型的缺点提出了采用L1型范数的方法将误差式改写，其形式为：ECLG-TV=∫Ω[(ΣregionW(x)(Ixu+Iyv+It)2)+λ(||▿u||+||▿v||)]---(1)]]>其中Ω表示一幅图中的所有像素点，region表示以当前像素点为中心的正方形窗口区域，W(x)是上述窗口区域内各点的权重，λ是权重系数，用来调节数据项和平滑项的比例；设I(x,y,t)是图像点(x,y)在时刻t的亮度，u(x,y)和v(x,y)是该点光流的x和y分量，点在t+δt时运动到(x+δx,y+δy)时，亮度保持不变，即I(x+uδt,y+vδt,t+δt)＝I(x,y,t)亮度随x、y、t光滑变化，将上式的左边用泰勒级数一阶展开：I(x,y,t)+δx∂I∂x+δy∂I∂y+δt∂I∂t+e=I(x,y,t)]]>其中e是关于δx、δy、δt的二阶和二阶以上的项；上式两边的I(x,y,t)相互抵消，两边除以δt，并取极限δt→0，得到这就是亮度不变性约束；(1)式中(∑regionW(x)(Ixu+Iyv+It)2)是LK方法的数据项，其本质是一个区域上的所有像素点满足亮度不变性约束Ixu+Iyv+It＝0；其中亮度对空间的梯度Ix和Iy的物理意义是亮度在水平、竖直方向的变化率，亮度对时间的梯度It的物理意义是亮度随时间的变化率，u和v代表像素点分别在水平，竖直方向的变化速度；(1)式中的是L1型的平滑项，由于实际场景中的物体都是连续，其运动矢量理应是大面积相同或相近的，因此设计平滑项这一项来满足这一特点；表示水平运动速度的空间梯度，反映了水平运动速度u在空间上的连续程度，L1型范数是指模的大小，其值越小连续程度越高；竖直方向同理；整个平滑项的物理意义是指求得的整个运动场应满足所有像素点的之和最小，即整个运动场的大部分点是足够连续，足够平滑的；这样就克服了深度不连续、弱纹理及重复纹理和噪声干扰造成的影响；对于亮度对空间的梯度Ix和Iy，采用了Barron算子与horn算子结合的方式；Barron算子是一维的，其形式为1/12(1 8 0 ‑8 ‑1)；horn算子是三维的，其计算方法如下式所示：Ix=14(Ii+1,j,k+Ii+1,j+1,k+Ii+1,j,k+1+Ii+1,j+1,k+1)-14(Ii,j,k+Ii,j+1,k+Ii,j,k+1+Ii,j+1,k+1)]]>Iy=14(Ii,j+1,k+Ii+1,j+1,k+Ii,j+1,k+1+Ii+1,j+1,k+1)-14(Ii,j,k+Ii+1,Lk+Ii,j,k+1+Ii+1,j,k+1)]]>It=14(Ii,j,k+1+Ii+1,j,k+1+Ii,j+1,k+1+Ii+1,j+1,k+1)-14(Ii,j,k+Ii+1,j,k+Ii,j+1,k+Ii+1,j+1,k)]]>Barron算子适合于计算空间梯度，horn算子适合于计算时间梯度；对于亮度对空间的梯度Ix和Iy，将Barron算子与horn算子计算出的梯度进行加权作为最后梯度，亮度对时间的梯度It，将horn算子算出的梯度作为最后的梯度；通过改进的双边滤波器来计算窗口权重函数，增强数据项中对边缘的保留，在双边滤波器的基础上增加阻塞度一项，使得权重正比于空间距离，亮度差别，阻塞度差别：BFW(x)=exp{-λs|x-x0|2σ12-λi|Ix-Ix0|2σ22}λoo(x)o(x0)---(2)]]>其中，式中x是指窗口内x点的空间坐标矢量，x0是指窗口内x0点的空间坐标矢量，Ix是指在x点处灰度值，是指在x0点处灰度值，λs和σ1用于调节空间距离权重，λi和σ2用于调节灰度值差距权重，λo用于调节阻塞度权重，x点处的阻塞度o(x)定义为o(x)=exp{-d(x)2σd2-It(x)2σe2}---(3)]]>其中，It(x)为x点处亮度对时间的梯度值，σd和σe用于调节It(x)与d(x)二者权重，d(x)为运动矢量的单边散度，ux和vx分别表示x点处的水平和垂直运动矢量；d(x)=div(ux,vx),div(ux,vx)<00,otherwise---(4)]]>双边滤波器将高斯空间滤波与距离滤波相结合，起到既去除噪声又保护边缘的作用，将双边滤波结果作为数据项权重表示当一个区域内的像素点离中心点空间距离越远，那么权重越小；当一个区域内的像素点灰度值与中心点灰度值相差越大，那么权重越小；在此基础上，引入阻塞度这一概念，阻塞度表示了当前像素点的运动矢量的杂乱程度，当一个区域内的像素点阻塞度与中心点阻塞度比值越小，那么权重越小；总之加强了对图像细节的精确计算；对光流算法中平滑项引入全局扩散系数，增强平滑项中对边缘的保留；全局扩散系数定义为：D(||▿I||)=e-α||▿I||β---(5)]]>式中α和β用于调节边缘保留强度；利用图像中的亮度信息来获取一扩散系数，扩散系数表明了图像中物体的连续程度，用它来对平滑项作加权处理使得图像中无纹理的连续区域平滑程度高，而边角这种非连续区域平滑程度低；总之进一步加强了弱纹理处的平滑程度；这样，对两项系数改进后得到的误差式为：ECLG-TV=∫Ω[(ΣregionBFW(x)(Ixu+Iyv+It)2)+λ(||D·▿u||+||D·▿v||)]---(6)]]>该误差式的平滑项是TV‑L1模型，其求解属于全变分正则化问题，直接求解是有困难的，现将其分为三部分：ECLG-TV-1=∫Ω[(ΣregionBFW(x)(Ixu^+Iyv^+It)2)+12·θ(u-u^)2+12·θ(v-v^)2]---(7)]]>ETV-u=∫Ω[12·θ(u-u^)2+||D·▿u||]---(8)]]>ETV-v=∫Ω[12·θ(v-v^)2+||D·▿v||]---(9)]]>式中θ用于控制分解比例，使用交替求解的方法来最小化上面三个式子：(7)式的求解是凸函数最优化问题，直接求其极值点作为最优解，作为未知量，将第一个式子分别对求偏导，并令其为零，整理化简得到方程组：1+2θΣBFW·Ix22θΣBFW·IxIy2θΣBFW·IxIy1+2θΣBFW·Iy2.u^v^=u-2θΣBFW·IxItv-2θΣBFW·IyIt---(10)]]>通过求解以上方程组得到第一个式子的解对于TV‑u和TV‑v式，对于这种形式的积分式构造欧拉‑‑拉格朗日方程然后使用原始对偶算法来求解，TV‑u的欧拉‑‑拉格朗日方程为：-div(D·▿u||▿u||)+1λ·(u-u^)=0---(11)]]>使用原始对偶算法，定义对偶向量为pu→=▿u||▿u||⇒pu→·||▿u||-▿u=0,||pu→||≤1---(12)]]>将欧拉—拉格朗日方程改写为-div(D·pu→)+1λ·(u-u^)=0⇒u=λ·div(D·pu→)+u^---(13)]]>再将(13)式代入(12)式可得：pu→·||▿(λ·div(D·pu→)+u^)||-▿(λ·div(D·pu→)+u^)=0⇒pu→·||▿(div(D·pu→)+u^λ)||-▿(div(D·pu→)+u^λ)=0---(14)]]>为求上式中结合方程形式可以使用固定点迭代的方法求解：此外，求解光流场还需要用到两个策略：1)由粗到精的多尺度计算，降低图像的分辨率，得到金字塔图像；由于当像素点位移大于滤波窗口的尺寸时，使用线性化的数据项直接计算光流会显著增加光流的误差，为了适用于计算较大的光流，一种解决方法是使用由粗到精的多尺度计算，较大的光流在低分辨率尺度变得很小，足以满足光流约束方程；本方法使用自动决定层数的CLG‑TV金字塔光流法，原始图像是最底层图像，高一层的图像由比它低一层的图像经过特定的高斯滤波然后进行下采样得到，此处选取每层的图像的长和宽都是比它低一层的图像的四分之一；由粗到精的多尺度计算很好的适应了计算较大的光流的要求，另外，将金字塔底层图像的分辨率从VGA级别降低至较低的分辨率，减少了计算量；2)翘曲迭代方法实际计算中，运动矢量的求解是经过多次翘曲得到的，将每次的计算结果作为右眼图像的自变量数组进行双线性差值翘曲，然后将左眼图像和翘曲后的右眼图像的运算结果更新上次的结果；整个CLG‑TV金字塔光流法的计算流程；最后，将程序的步骤以伪代码的形式列在下面：1)建立图像金字塔I1，I2，并计算各层梯度；2)对各层I1计算全局扩散系数D；3)每层进行如下步骤：a)如果在最高层，将u,v,初始化为0；b)如果不在最高层，将上层的u,v,上采样，作为本层的初始值；c)for i＝1to warpsi)对u,v中值滤波，减少异常值；ii)用双线性插值对I1,I2,进行warp操作iii)for k＝1to iterationsA.使用(10)式计算B.使用(13)式更新u,v；C.使用(15)式更新本步骤在CLG‑TV光流法的基础上，对光流算法中平滑项进行了改进，使用了各向异性的平滑权重，改善了平滑效果；对光流算法中数据项的梯度算子、窗口权重函数进行改进，加强了对图像细节的精确计算；使用由粗到精的多尺度计算，翘曲迭代方法，加强了对大位移运动的鲁棒性；求解光流场仅计算到第3层为止，目的是对最大视差的像素所在位置作一个初步的估计，通过很小的计算量准确地得到最大视差在整幅图像中的大致位置，以待之后进行精确计算；步骤3、对步骤2计算得到的光流场进行处理，找出其中的最大值umax以及最大值所在的位置(iumax,jumax)，在I1中的该位置截取原图像的一部分，并找到I2中经最大值移动后的该位置(iumax+umax·scale,jumax)，其中scale表示真实图像和当前计算的第3层图像的尺度比，在此处截取原图像的一部分，如图3中所示；光流法的精度保证了截取的小图的计算准确度，同时，截取一小幅图片大大减少了数据量使得计算速度提高；步骤4、使用CLG‑TV光流法精确求解光流场，本步骤与步骤2的不同之处在于：本步骤使用的输入图像I1，I2是步骤3输出的截取图像，而不是原图像；为保证精度，本步骤使用的金字塔采样系数为0.75，并计算到第1层为止；本步骤的LK窗口宽度为11，而步骤2中为5；本层的输出为截取图像的光流；步骤5、将步骤4的输出换算为第0层的值后取其水平运动矢量，作出直方图，从中计算最大值时，对直方图进行一阶差分，然后对一阶差分序列从最小值开始考察，取阈值为100，当一阶差分大于该阈值时，将低于该阈值的运动矢量取出，计算其方差D(uψ)，并根据方差做出判断：umax=min(uψ),D(uψ)≤αmax(uψ),otherwise]]>式中uψ表示低于阈值的运动矢量，将此时得到的最大值与步骤3中得到的步骤2计算结果最大值作算术加，得到的结果为整个系统的输出即最大视差。