H-S光流法
一、实验目的
计算出图像序列的光流场,并用有向箭头在图中标出。
二、算法实现
1、构建图像金字塔
金字塔的每一层包括:图像的灰度图、x和y方向的梯度图、相邻两帧灰度图的像素灰度差(
It
)。
a、计算图像x、y、t方向的梯度。 Ix、Iy 计算方法:每个像素沿x、y方向的两个相邻像素的灰度差的一半。 It 计算方法:当前帧与上一帧对应位置的像素灰度差。代码如下:
//x方向梯度
void GradImg::calculateDerivativeX()
{
gradx = cv::Mat::zeros(im.size(),CV_64F);
for(int y = 0; y < im.rows; ++y)
{
for(int x = 0; x < im.cols; ++x)
{
int prev = std::max(x - 1, 0);
int next = std::min(x + 1, im.cols - 1);
gradx.at<double>(y, x) = (double) (im.at<double>(y, next) - im.at<double>(y, prev)) * 0.5f;
}
}
}
//y方向梯度
void GradImg::calculateDerivativeY()
{
grady = cv::Mat::zeros(im.size(),CV_64F);
for(int y = 0; y < im.rows; ++y)
{
for(int x = 0; x < im.cols; ++x)
{
int prev = std::max(y - 1, 0);
int next = std::min(y + 1, im.rows - 1);
grady.at<double>(y, x) = (double) (im.at<double>(next, x) - im.at<double>(prev, x)) * 0.5f;
}
}
}
//t方向梯度
void GradImg::calculateDerivativeT(GradImg& ref)
{
gradt = cv::Mat::zeros(im.size(),CV_64F);
for(int y = 0; y < gradt.rows; ++y)
{
for(int x = 0; x < gradt.cols; ++x)
{
gradt.at<double>(y, x) = (double) (im.at<double>(y, x) - ref.getImg().at<double>(y, x));
}
}
}
b、计算图像金字塔,下采样方法为均匀采样,代码如下:
//均匀采样
cv::Mat HSOpticalFlow::pyrDownMeanSmooth(const cv::Mat& in)
{
cv::Mat out;
out = cv::Mat::zeros(cv::Size(in.size().width / 2, in.size().height / 2),in.type());
for(int y = 0; y < out.rows; ++y)
{
for(int x = 0; x < out.cols; ++x)
{
int x0 = x * 2;
int x1 = x0 + 1;
int y0 = y * 2;
int y1 = y0 + 1;
out.at<double>(y, x) = (double) ( (in.at<double>(y0, x0) + in.at<double>(y0, x1) + in.at<double>(y1, x0) + in.at<double>(y1, x1)) / 4.0f );
}
}
return out.clone();
}
//建图像金字塔
void HSOpticalFlow::buildpyr()
{
//清空存放金字塔的容器
curlevels.clear();
reflevels.clear();
reflevels.push_back(ref);
curlevels.push_back(cur);
cv::Mat refcopy;
cv::Mat curcopy;
refcopy = ref.getImg();
curcopy = cur.getImg();
GradImg gradref;
GradImg gradcur;
for (int i = 0 ; i < py_layers ; i++)
{
cv::Mat refdown = pyrDownMeanSmooth(refcopy);
gradref.setImg(refdown);
reflevels.push_back(gradref);
refcopy = refdown.clone();
cv::Mat curdown = pyrDownMeanSmooth(curcopy);
gradcur.setImg(curdown);
gradcur.calculateDerivativeX();
gradcur.calculateDerivativeY();
gradcur.calculateDerivativeT(gradref);
curlevels.push_back(gradcur);
curcopy = curdown.clone();
}
}
2、计算光流
沿着金字塔,从上往下,依次计算光流。上一层的光流通过插值法传到下一层,作为下一层光流的初值。
这样做的好处:由于光流法基于灰度一致性假设,只能处理帧间运动较小的情况,使用金字塔,可以由粗到精计算大运动情况下的光流。
每个像素的光流计算方法如下:
其中
a
为常数,反映了对图像数据及平滑约束的可信度,当图像数据本身含有较大噪声时,此时需要加大
通过在每一层使用上式迭代计算每个像素的光流,再将该层光流传到下一层,作为下一层光流初值,计算下一层的光流。代码如下:
while(iter < maxiter)
{
iter++;
for(int y = 0; y < im.rows; ++y)
{
for(int x = 0; x < im.cols; ++x)
{
double mm = optow_ub.at<double>(y,x)*dx.at<double>(y,x) +optow_vb.at<double> (y,x)*dy.at<double>(y,x) + dt.at<double>(y,x);
double nn = thea + dx.at<double>(y,x)*dx.at<double>(y,x) + dy.at<double>(y,x)*dy.at<double>(y,x);
optow_u.at<double>(y,x) = optow_ub.at<double>(y,x) - dx.at<double>(y,x)*(mm/nn);
optow_v.at<double>(y,x) = optow_vb.at<double>(y,x) - dy.at<double>(y,x)*(mm/nn);
}}
if ( iter < maxiter )
{
optow_ub = smooth<double>(optow_u);
optow_vb = smooth<double>(optow_v);
}}
3、光流可视化
用箭头表示图像像素的光流,箭头方向表示光流方向,长度表示光流大小。设定一个阈值区间,只显示在该区间内的光流。代码如下:
//在图片的相应像素上画出光流箭头
for (int y = 0 ; y < u.rows ; y++)
for (int x = 0 ; x < u.cols ; x++)
{
if ((abs(u.at<double>(y,x))>th1 && abs(u.at<double>(y,x))<th2) &&
(abs(v.at<double>(y,x))>th1 && abs(v.at<double>(y,x))<th2))
{
cv::Point pStart = Point(x, y);
cv::Point pEnd = Point(x+u.at<double>(y,x) *1.8f, y+v.at<double>(y,x) *1.8f);
drawArrow(rgbcur, pStart, pEnd);
}
}
cv::imshow ( rgb_file, rgbcur.clone() );
三、实验结果
1、Traffic数据集
在图像中显示的光流的上下阈值分别设为:2.6、12.8。
2、outdoor数据集
参考资料
https://www.cnblogs.com/rocbomb/p/3795719.html?utm_source=tuicool&utm_medium=referral
http://blog.csdn.net/u014568921/article/details/46638557