Horn-Schunck Optical Flow (Averaging Velocity)

Question

I'm currently trying to implement the HS-method for optical flow but my u and v always seem to have only zeros in them. I can't seem to figure out my error in here:

vid=VideoReader('outback.AVI');
vid.CurrentTime = 1.5;

alpha=1;
iterations=10;    

frame_one = readFrame(vid);
vid.CurrentTime = 1.6;
frame_two = readFrame(vid);

% convert to grayscale

fone_gr = rgb2gray(frame_one);
ftwo_gr = rgb2gray(frame_two);

% construct  for each image
sobelx=[-1 -2 -1; 0 0 0; 1 2 1];
sobely=sobelx';
time=[-1 1];

fx_fone=imfilter(fone_gr, sobelx);
fy_fone=imfilter(fone_gr, sobely);
ft_fone=imfilter(fone_gr, time);

fx_ftwo=imfilter(ftwo_gr, sobelx);
fy_ftwo=imfilter(ftwo_gr, sobely);
ft_ftwo=imfilter(ftwo_gr, time);

Ix=double(fx_fone+fx_ftwo);
Iy=double(fy_fone+fy_ftwo);
It=double(ft_fone+ft_ftwo);

% flow-variables (velocity = 0 assumption)
velocity_kernel=[0 1 0; 1 0 1; 0 1 0];

u = double(0);
v = double(0);

% iteratively solve for u and v
for i=1:iterations
    neighborhood_average_u=conv2(u, velocity_kernel, 'same');
    neighborhood_average_v=conv2(v, velocity_kernel, 'same');

    data_term = (Ix .* neighborhood_average_u + Iy .* neighborhood_average_v + It);
    smoothness_term = alpha^2 + (Ix).^2 + (Iy).^2;

    numerator_u = Ix .* data_term;
    numerator_v = Iy .* data_term;

    u = neighborhood_average_u - ( numerator_u ./ smoothness_term );
    v = neighborhood_average_v - ( numerator_v ./ smoothness_term );
end

u(isnan(u))=0;
v(isnan(v))=0;

figure
imshow(frame_one); hold on;
quiver(u, v, 5, 'color', 'b', 'linewidth', 2);
set(gca, 'YDir', 'reverse');

The only thing I'm not really confident about is the computation of the neighborhood average:

velocity_kernel=[0 1 0; 1 0 1; 0 1 0];

u = double(0);
v = double(0);

[..]
neighborhood_average_u=conv2(u, velocity_kernel, 'same');
neighborhood_average_v=conv2(v, velocity_kernel, 'same');

Wouldn't that always result in a convolution matrix with only zeros?

I thought about changing it to the following, since I need to compute the average velocity using the velocity kernel on each pixel of my images:

velocity_kernel=[0 1 0; 1 0 1; 0 1 0];
u = double(0);
v = double(0);

[..]
neighborhood_average_u=conv2(Ix, velocity_kernel, 'same');
neighborhood_average_v=conv2(Iy, velocity_kernel, 'same');

But I still don't know if that would be the correct way. I followed the instructions on the bottom of this MATLAB page: http://de.mathworks.com/help/vision/ref/opticalflowhs-class.html

Knipser Knipser · Accepted Answer · 2016-06-20T06:39:28

I found this paper with some further explanations and also some matlab code.

They compute the average of u and v as follows:

% initial values 
u = 0; v = 0;

% weighted average kernel 
kernel = [1/12 1/6 1/12; 1/6 0 1/6; 1/12 1/6 1/12];

for i = 1:iterations
  uAvg = conv2( u, kernel 'same' );
  vAvg = conv2( v, kernel 'same' );
  ...
end

Horn-Schunck Optical Flow (Averaging Velocity)

1 Answers