目标跟踪的点跟踪技术(1)

来源：互联网发布：python 打开文件夹编辑：程序博客网时间：2024/04/30 12:36

目标跟踪技术按实际跟踪对象可以分为点跟踪和块跟踪。

所谓点跟踪就是在初始图像帧的目标上找一些具有跟踪价值的点，用点周围的一小块区域的特征对其进行描述，在后续的图像帧中根据特征描述寻找这些点移动到的新位置。

这里需要解决三个问题。特征点选择、特征点描述和特征点匹配。

首先谈一下特征点选择。用于跟踪的特征点周围的纹理应该是复杂的。如果在一片白墙上选择一个特征点，那么显然不可能知道下一帧图像中这个点到底在什么地方。图像周围的纹理越丰富，越容易和其他点区别开。

（该图片来自Richard Szeliski著作Computer Vision: Algorithms and Applications，这里移动的不是景物而是摄像机，原理是一样的）

另外，某点周围如果是比较单纯的边缘划分，也不容易进行跟踪用途，因为仅凭一条边缘的变化无法了解物体的实际移动方向，这一问题称为孔径问题(aperture problem)，详情可以参考以下链接。

http://web.mit.edu/persci/demos/Motion&Form/demos/one-square/one-square.html

因此，比较理想的特征点应该在所谓“角点”的位置，即以该点为中心至少有两个方向有明显的色彩或亮度变化。

这一特性被Harris描述为公式形式。

对于微小的移动(x, y)，周围区域的所有点移动前后的亮度差之和为

$S(x,y) = \sum_u \sum_v w(u,v) \, \left( I(u+x,v+y) - I(u,v)\right)^2$

其中

$w_{u,v}= \exp{-(u^2 + v^2)/2\sigma^2}$

利用泰勒公式，近似有

$I(u+x,v+y) \approx I(u,v) + I_x(u,v)x+I_y(u,v)y$

$I_x$ ， $I_y$ 为对应的偏导。

应用以上近似关系，S可以写成

$S(x,y) \approx \begin{pmatrix} x & y \end{pmatrix} A \begin{pmatrix} x \\ y \end{pmatrix},$

$A = \sum_u \sum_v w(u,v) \begin{bmatrix}I_x^2 & I_x I_y \\I_x I_y & I_y^2 \end{bmatrix}=\begin{bmatrix}\langle I_x^2 \rangle & \langle I_x I_y \rangle\\\langle I_x I_y \rangle & \langle I_y^2 \rangle\end{bmatrix}$

根据之前分析，S应该受到x，y两个方向的变化敏感，因此A的两个特征值都应较大。

可以设置一个阈值，特征值大于这个阈值的点就是用于跟踪的特征点。

因为计算特征值运算较复杂，实际计算以下值

$M_c = \lambda_1 \lambda_2 - \kappa \, (\lambda_1 + \lambda_2)^2= \operatorname{det}(A) - \kappa \, \operatorname{trace}^2(A)$

另外使用局部极大抑制来减少候选点的数量。

以下是程序实现：

function pts = harris( I, k, thres, gsize, sigma )if nargin < 2    k = 0.06;endif nargin < 3    thres = 0.1;endif nargin < 4    gsize = 7;endif nargin < 5    sigma = 2;endif ndims(I) > 2    I = rgb2gray(I);endI = double(I);[height, width] = size(I);hx = [-1 0 1];Ix = imfilter(I, hx, 'replicate');hy = [-1; 0; 1];Iy = imfilter(I, hy, 'replicate');Ix2 = Ix .^2;Iy2 = Iy .^2;Ixy = Ix .* Iy;h = fspecial('gaussian', [gsize, gsize], sigma);Ix2 = filter2(h, Ix2);Iy2 = filter2(h, Iy2);Ixy = filter2(h, Ixy);R = zeros(height, width);Rmax = 0;for i = 1 : height    for j = 1 : width        M = [Ix2(i, j) Ixy(i, j); Ixy(i, j) Iy2(i, j)];                     R(i,j) = det(M) - k * (trace(M)) ^ 2;                        if R(i,j) > Rmax            Rmax = R(i, j);        end;    end;end;thres = Rmax * thres;%non-maximum supressionlocalMax = zeros(height, width);supRadius = 1;for y = supRadius + 1 : height - supRadius    for x = supRadius + 1 : width - supRadius        localMax(y, x) = isMax(R, supRadius, x, y, thres);    endend[ptsc, ptsr] = find(localMax == 1);pts = [ptsr, ptsc];endfunction b = isMax(R, radius, x, y, thres)    if R(y, x) < thres        b = 0;        return;    end     for i = y - radius : y + radius        for j = x - radius : x + radius            if R(y, x) < R(i, j)                b = 0;                return;            end        end    end    b = 1;end

效果如图

0 0