Implementation based on Horn’s optical flow algorithm.

function [Vx,Vy] = OpticalFlow(images,alpha,iterations)
%// Calculating optical flow of a sequence of images.
%// images : 3D array that contains a sequence of images. size of images is (imageHeight, imageWidth, frame number)
%// alpha
%// iterations.
[height,width,frames]=size(images);
%//initialzation of u and v
Vx = zeros(height,width);
Vy = zeros(height,width);
for k = 1:frames-1
% //initialization of Ex Ey and Et
Ex = zeros(height-1,width-1,frames-1);
Ey = zeros(height-1,width-1,frames-1);
Et = zeros(height-1,width-1,frames-1);
%//calculating Ex Ey and Et in frame k.
for x = 2:width-1
for y = 2:height-1
Ex(y,x,k) = (images(y+1,x+1,k)-images(y+1,x,k)+images(y,x+1,k)...
-images(y,x,k)+images(y+1,x+1,k+1)-images(y+1,x,k+1)...
+images(y,x+1,k+1)-images(y,x,k+1))/4;
Ey(y,x,k) = (images(y,x,k)-images(y+1,x,k)+images(y,x+1,k)...
-images(y+1,x+1,k)+images(y,x,k+1)-images(y+1,x,k+1)...
+images(y,x+1,k+1)-images(y+1,x+1,k+1))/4;
Et(y,x,k) = (images(y+1,x,k+1)-images(y+1,x,k)+images(y,x,k+1)...
-images(y,x,k)+images(y+1,x+1,k+1)-images(y+1,x+1,k)...
+images(y,x+1,k+1)-images(y,x+1,k))/4;
end
end
for nn = 1:iterations
for x = 2:width-1
for y = 2:height-1
Vxbar = (Vx(y-1,x)+Vx(y,x+1)+Vx(y+1,x)+Vx(y,x-1))/6+...
(Vx(y-1,x-1)+Vx(y-1,x+1)+Vx(y+1,x+1)+Vx(y+1,x-1))/12;
Vybar = (Vy(y-1,x)+Vy(y,x+1)+Vy(y+1,x)+Vy(y,x-1))/6+...
(Vy(y-1,x-1)+Vy(y-1,x+1)+Vy(y+1,x+1)+Vy(y+1,x-1))/12;
%// chapter 12 of Horn's paper
temp = (Ex(y,x,k)*Vxbar+Ey(y,x,k)*Vybar+Et(y,x,k))/(alpha^2 + Ex(y,x,k)^2 + Ey(y,x,k)^2);
%// update u and v
Vx(y,x) = Vxbar-Ex(y,x,k)*temp;
Vy(y,x) = Vybar-Ey(y,x,k)*temp;
end
end
end
end

Advertisements

Like this:

LikeLoading...

Related

This entry was posted on Monday, June 9th, 2008 at 7:26 pm and is filed under computer vision, matlab. You can follow any responses to this entry through the RSS 2.0 feed.
You can leave a response, or trackback from your own site.

can you please explain what you ment by “images”?
if i have 2 images that i want to check the flow between them how do i create one “images” var from them?

I apply this implementation to the images sequence “office” from http://of-eval.sourceforge.net/.
The results are quite inconsistent. I am wondering if I am using this code correctly.
The code I wrote:
%read the image data
X(:,:,1)=imread(‘anim.20.pgm’); % X data type unit8
X(:,:,2)=imread(‘anim.21.pgm’);
%using your implementation
[vxx,vyy]=OpticalFlow(Y,0.1,5);

Is there any special requirements about the input data X?

Hey Bin Liu, this is a detailed implementation of Horn paper, it’s really nice. Have you been able to fix inconsistent results with the one appearing on http://of-eval.sourceforge.net/?

Jason thanks for the shared knowledge; did you able to get expected results on synthetic data?

Hey,
I’m sorry, i’m pretty new to Matlab and i still don’t understand how to use your algorithm.
can you please explain with the next example? :
lets say i have 2 RGB images that i change to grayscale and get the matrixes ‘E1’ and ‘E2’ with size = [240 352].
how do i create the var ‘images’ from the 2 matrixes i’ve got?
and then, do you think that the parameters: the created ‘images’,’alpha’ = 1, and ‘iterations’ = 16, are reasonable?

Thank you so much!
and i apologize for driving you nuts here 🙂
Matt

Hi, I believe in order for me to visualize the results I have two options:

1-Apply a velocity threshold and convert Vx and Vy into binary values and then use regionprops to apply a bounding box around the objects that are being tracked.
2-I can plot the optical flow lines on the (I dont know ho exactly) original image.

Has anyone tried this algorithm and visualized the results? I am not quite sure where to feed the Vx and Vy to be able to detect moving objects.

Thank you very much for your code. I have a single motion blurred image. Is it possible to use this code to know the optical flow of a single motion blurred image? How can I do that based on your code since your code is intended for sequence of images?

OK, I have tried this algortihm on a 31-frame intensity video (height*width*31).
One interesting thing is that the output Vx and Vy matrices are of form (height*width), whereas I expected them to be (height*width*30), i.e. giving a different set of Vx-Vy vector for each pixel, comparing each consecutive frame.
You have to put the code in a for loop if you want to achieve a frame by frame optical flow output.
Still, the results (at least for the intentions of object tracking) are not very impressive I must say…

I noticed the same thing, expected a 30 frame output, but it’s just 1 frame. From the code, it seems it’s only returning the optic flow for the last 2 frames, even though it’s calctulating all of it. So I added a ‘k’ in the Vx, Vy assignment at the end of the for loop. What do you think?

The algorithm is to calculate the overall optical flow of all the frames. To get the result you expected, you will have to use two consecutive frames each time and combine the results. Or you can try Diana’s method.

hai guys u r doing a great work i’m not new to matlab but i’m unable to run the code which you have provided.
i need to know how to run it at same time is any image input need to be given and can any 1 specify the software that are used to test optical flow and on wat platform, requirements etc.

Hello! I used this implementation in an application regarding obstacle avoidance. I want to ask you what is the measurement unit for the Vx and Vy components. I need this information to obtain the unit for the time-to-contact. Thanks!

hi all ,I need the code for calculating average angular error for optical flow measurement .Please help me and it’s very urgent .I’m looking forward for it.

those referring to the sourceforge website and finding inconsistencies is because you need to smooth the image first. I used smoothImg() function from Horn Schunk code at mathworks for smoothening before applying this code and I got the results shown at the site. The direction of optical flow is opposite to what is shown though.

what do you mean by term images? if we have a video and we have to compute optical flow between the consecutive frames of that video then how should we input it and then what is the meaning of the argument images which you hav specified in your code.

end
function[Vx,Vy] = OpticalFlow(images,alpha,iterations)
%// Calculating optical flow of a sequence of images.
%// images : 3D array that contains a sequence of images. size of images is (imageHeight, imageWidth, frame number)
%// alpha
%// iterations.
[height,width,frames]=size(images);

%//initialzation of u and v
Vx = zeros(height,width);
Vy = zeros(height,width);

for k = 1:frames-1

% //initialization of Ex Ey and Et
Ex = zeros(height-1,width-1,frames-1);
Ey = zeros(height-1,width-1,frames-1);
Et = zeros(height-1,width-1,frames-1);

%//calculating Ex Ey and Et in frame k.
for x = 2:width-1
for y = 2:height-1
Ex(y,x,k) = (images(y+1,x+1,k)-images(y+1,x,k)+images(y,x+1,k)…
-images(y,x,k)+images(y+1,x+1,k+1)-images(y+1,x,k+1)…
+images(y,x+1,k+1)-images(y,x,k+1))/4;

Et(y,x,k) = (images(y+1,x,k+1)-images(y+1,x,k)+images(y,x,k+1)…
-images(y,x,k)+images(y+1,x+1,k+1)-images(y+1,x+1,k)…
+images(y,x+1,k+1)-images(y,x+1,k))/4;
end
end

for nn = 1:iterations
for x = 2:width-1
for y = 2:height-1

%// chapter 12 of Horn’s paper
temp = (Ex(y,x,k)*Vxbar+Ey(y,x,k)*Vybar+Et(y,x,k))/(alpha^2 + Ex(y,x,k)^2 + Ey(y,x,k)^2);
%// update u and v
Vx(y,x) = Vxbar-Ex(y,x,k)*temp;
Vy(y,x) = Vybar-Ey(y,x,k)*temp;
end
end
end

end

Hi I am not able to run this your code .. there have some error in the function .Can u plz help me out

Try this code it works…..
kindly zoom the figure window to observe the optical flow….

a1=imread(‘TW1.jpg’);
a2=imread(‘TW2.jpg’);

b1 = rgb2gray(a1);
b2 = rgb2gray(a2);

images(:,:,1)=b1;
images(:,:,2)=b2;

[height,width,frames]=size(images);

%//initialzation of u and v
Vx = zeros(height,width);
Vy = zeros(height,width);

for k = 1:frames-1

% //initialization of Ex Ey and Et
Ex = zeros(height-1,width-1,frames-1);
Ey = zeros(height-1,width-1,frames-1);
Et = zeros(height-1,width-1,frames-1);

%//calculating Ex Ey and Et in frame k.
for x = 2:width-1
for y = 2:height-1
Ex(y,x,k) = (images(y+1,x+1,k)-images(y+1,x,k)+images(y,x+1,k)-images(y,x,k)+images(y+1,x+1,k+1)-images(y+1,x,k+1)+images(y,x+1,k+1)-images(y,x,k+1))/4;

Et(y,x,k) = (images(y+1,x,k+1)-images(y+1,x,k)+images(y,x,k+1)-images(y,x,k)+images(y+1,x+1,k+1)-images(y+1,x+1,k)+images(y,x+1,k+1)-images(y,x+1,k))/4;
end
end

%// chapter 12 of Horn’s paper
temp = (Ex(y,x,k)*Vxbar+Ey(y,x,k)*Vybar+Et(y,x,k))/(0.5^2 + Ex(y,x,k)^2 + Ey(y,x,k)^2);
%update u and v
Vx(y,x) = Vxbar-Ex(y,x,k)*temp;
Vy(y,x) = Vybar-Ey(y,x,k)*temp;
end
end
end

With havin so much content and articles do you ever run into
any problems of plagorism or copyright infringement? My site
has a lot of completely unique content I’ve either authored myself
or outsourced but it looks like a lot of it is popping it up all over the web without
my permission. Do you know any methods to help stop content from being ripped off?
I’d really appreciate it.

Usually I do not read article on blogs, but I would like to say that this write-up very pressured me to take a look at and do it! Your writing taste has been surprised me. Thanks, quite nice post.|

Thanks for the code!

can you please explain what you ment by “images”?

if i have 2 images that i want to check the flow between them how do i create one “images” var from them?

cheers!

Matt

i want the program for my research

Matt,

Yes, create an 3-D array in the format images(imageHeight, imageWidth, numberOfFrames)

I apply this implementation to the images sequence “office” from http://of-eval.sourceforge.net/.

The results are quite inconsistent. I am wondering if I am using this code correctly.

The code I wrote:

%read the image data

X(:,:,1)=imread(‘anim.20.pgm’); % X data type unit8

X(:,:,2)=imread(‘anim.21.pgm’);

%using your implementation

[vxx,vyy]=OpticalFlow(Y,0.1,5);

Is there any special requirements about the input data X?

Hey Bin Liu, this is a detailed implementation of Horn paper, it’s really nice. Have you been able to fix inconsistent results with the one appearing on http://of-eval.sourceforge.net/?

Jason thanks for the shared knowledge; did you able to get expected results on synthetic data?

Hey,

I’m sorry, i’m pretty new to Matlab and i still don’t understand how to use your algorithm.

can you please explain with the next example? :

lets say i have 2 RGB images that i change to grayscale and get the matrixes ‘E1’ and ‘E2’ with size = [240 352].

how do i create the var ‘images’ from the 2 matrixes i’ve got?

and then, do you think that the parameters: the created ‘images’,’alpha’ = 1, and ‘iterations’ = 16, are reasonable?

Thank you so much!

and i apologize for driving you nuts here 🙂

Matt

Bin Liu answered my question.

i will check if i have the same problem as he stated above and let you know.

thanks 🙂

Matt

Hi

m doing a project on Video Segmentation [ Intelligent Video Monitoring System ] n currently working on optical flow.

Plz tell how to get SEQUENCE OF IMAGES..??

Ben, m new to Matlab.. can u Plz tell what is Y in ur code snippet.?

For the image issue asked before: use the follwing snippet.

cat(3,Image0,Image1);

Hi, I believe in order for me to visualize the results I have two options:

1-Apply a velocity threshold and convert Vx and Vy into binary values and then use regionprops to apply a bounding box around the objects that are being tracked.

2-I can plot the optical flow lines on the (I dont know ho exactly) original image.

Has anyone tried this algorithm and visualized the results? I am not quite sure where to feed the Vx and Vy to be able to detect moving objects.

thanks,

berkan

Dear Sir,

Thank you very much for your code. I have a single motion blurred image. Is it possible to use this code to know the optical flow of a single motion blurred image? How can I do that based on your code since your code is intended for sequence of images?

I would like to thank for your kind help.

Berkan,

You can use matlab function “quiver” to plot vectors.

Samira,

My code is intended for sequence of images. By the way, how do you get the motion vectors from a single image?

Hello Samira,

What you are looking for is not optical flow. In order to deal with blurring (deblur), you can try “deconvblind” function at the VIP toolbox.

hth.

OK, I have tried this algortihm on a 31-frame intensity video (height*width*31).

One interesting thing is that the output Vx and Vy matrices are of form (height*width), whereas I expected them to be (height*width*30), i.e. giving a different set of Vx-Vy vector for each pixel, comparing each consecutive frame.

You have to put the code in a for loop if you want to achieve a frame by frame optical flow output.

Still, the results (at least for the intentions of object tracking) are not very impressive I must say…

Berkan,

I noticed the same thing, expected a 30 frame output, but it’s just 1 frame. From the code, it seems it’s only returning the optic flow for the last 2 frames, even though it’s calctulating all of it. So I added a ‘k’ in the Vx, Vy assignment at the end of the for loop. What do you think?

Vx(y,x,k) = Vxbar-Ex(y,x,k)*temp;

Vy(y,x,k) = Vybar-Ey(y,x,k)*temp;

D.

The algorithm is to calculate the overall optical flow of all the frames. To get the result you expected, you will have to use two consecutive frames each time and combine the results. Or you can try Diana’s method.

Hope this helps.

hey,i did’t understand how can i get input arguement images ….m new to matlab

how to plot final result & how to choose alpha & iterations also use this [height,width,frames]=size(images); statement

how can i view the optical flow .I try :

quiver(Vx,Vy) but no resultat

please help me .Thank u a lot.

%% to run the code above

a1=imread(‘ur image 1’);

a2=imread(‘ur image 2’);

ai1 = rgb2gray(a1);

ai2 = rgb2gray(a2);

X(:,:,1)=ai1;

X(:,:,2)=ai2;

[Vx,Vy] = OpticalFlow1(X,0.1,5);

%% to view result

quiver(Vx,Vy);

%% for Avi to frames or seq …

a=aviread(‘ur video.avi’);

s=size(a);

for u=s(1):s(2)

f=a(u).cdata;

ai = rgb2gray(f);

X(:,:,u)=ai;

end

hope it will be helpful !!!

hello!I run the program.But I use ‘quiver’ to see image as arrows!How can do it?

Thanks

Hi again!I know my message made u confused.Sorry.I meant I used ‘quiver’ but I couldn’t see image with arrows.How I can see it?

Thanks

how did u define the function OpticalFlow1

hai guys u r doing a great work i’m not new to matlab but i’m unable to run the code which you have provided.

i need to know how to run it at same time is any image input need to be given and can any 1 specify the software that are used to test optical flow and on wat platform, requirements etc.

dont consider it as spam plzzzzzzz

i desperately need it

Hello! I used this implementation in an application regarding obstacle avoidance. I want to ask you what is the measurement unit for the Vx and Vy components. I need this information to obtain the unit for the time-to-contact. Thanks!

hi all ,I need the code for calculating average angular error for optical flow measurement .Please help me and it’s very urgent .I’m looking forward for it.

Hi I’m new to matlab and I’m trying calculate the velocity of a falling object.

Why do I get this error when I use optical flow?

??? Undefined function or method ‘OpticalFlow’ for input arguments of

type ‘uint8’.

Please help me

i’m sorry, but this Optical flow is not robust…

those referring to the sourceforge website and finding inconsistencies is because you need to smooth the image first. I used smoothImg() function from Horn Schunk code at mathworks for smoothening before applying this code and I got the results shown at the site. The direction of optical flow is opposite to what is shown though.

Ok i will try it again, Thanks alot

what do you mean by term images? if we have a video and we have to compute optical flow between the consecutive frames of that video then how should we input it and then what is the meaning of the argument images which you hav specified in your code.

im1 = double(imread(‘C:\Users\User\Desktop\Measure_Optical Flow\1.bmp’))/256;

im2 = double(imread(‘C:\Users\User\Desktop\Measure_Optical Flow\2.bmp’))/256;

figure; imshow(im1);

figure; imshow(im2);

a1=imread(‘C:\Users\User\Desktop\Measure_Optical Flow\1.bmp’);

a2=imread(‘C:\Users\User\Desktop\Measure_Optical Flow\2.bmp’);

ai1 = rgb2gray(a1);

ai2 = rgb2gray(a2);

X(:,:,1)=ai1;

X(:,:,2)=ai2;

[Vx,Vy] = OpticalFlow1(X,0.1,5);

%% to view result

quiver(Vx,Vy);

%% for Avi to frames or seq …

%a=aviread();

s=size(a);

for u=s(1):s(2)

f=a(u).cdata;

ai = rgb2gray(f);

X(:,:,u)=ai;

end

function[Vx,Vy] = OpticalFlow(images,alpha,iterations)

%// Calculating optical flow of a sequence of images.

%// images : 3D array that contains a sequence of images. size of images is (imageHeight, imageWidth, frame number)

%// alpha

%// iterations.

[height,width,frames]=size(images);

%//initialzation of u and v

Vx = zeros(height,width);

Vy = zeros(height,width);

for k = 1:frames-1

% //initialization of Ex Ey and Et

Ex = zeros(height-1,width-1,frames-1);

Ey = zeros(height-1,width-1,frames-1);

Et = zeros(height-1,width-1,frames-1);

%//calculating Ex Ey and Et in frame k.

for x = 2:width-1

for y = 2:height-1

Ex(y,x,k) = (images(y+1,x+1,k)-images(y+1,x,k)+images(y,x+1,k)…

-images(y,x,k)+images(y+1,x+1,k+1)-images(y+1,x,k+1)…

+images(y,x+1,k+1)-images(y,x,k+1))/4;

Ey(y,x,k) = (images(y,x,k)-images(y+1,x,k)+images(y,x+1,k)…

-images(y+1,x+1,k)+images(y,x,k+1)-images(y+1,x,k+1)…

+images(y,x+1,k+1)-images(y+1,x+1,k+1))/4;

Et(y,x,k) = (images(y+1,x,k+1)-images(y+1,x,k)+images(y,x,k+1)…

-images(y,x,k)+images(y+1,x+1,k+1)-images(y+1,x+1,k)…

+images(y,x+1,k+1)-images(y,x+1,k))/4;

end

end

for nn = 1:iterations

for x = 2:width-1

for y = 2:height-1

Vxbar = (Vx(y-1,x)+Vx(y,x+1)+Vx(y+1,x)+Vx(y,x-1))/6+…

(Vx(y-1,x-1)+Vx(y-1,x+1)+Vx(y+1,x+1)+Vx(y+1,x-1))/12;

Vybar = (Vy(y-1,x)+Vy(y,x+1)+Vy(y+1,x)+Vy(y,x-1))/6+…

(Vy(y-1,x-1)+Vy(y-1,x+1)+Vy(y+1,x+1)+Vy(y+1,x-1))/12;

%// chapter 12 of Horn’s paper

temp = (Ex(y,x,k)*Vxbar+Ey(y,x,k)*Vybar+Et(y,x,k))/(alpha^2 + Ex(y,x,k)^2 + Ey(y,x,k)^2);

%// update u and v

Vx(y,x) = Vxbar-Ex(y,x,k)*temp;

Vy(y,x) = Vybar-Ey(y,x,k)*temp;

end

end

end

end

Hi I am not able to run this your code .. there have some error in the function .Can u plz help me out

HI,i did’t understand how can i get input argument images …. I am new to mat lab… Plz help me…

I am getting NaN values for Vx and Vy. What should i do,i want to identify stationary and moving pixels in the images. Please help….

Try this code it works…..

kindly zoom the figure window to observe the optical flow….

a1=imread(‘TW1.jpg’);

a2=imread(‘TW2.jpg’);

b1 = rgb2gray(a1);

b2 = rgb2gray(a2);

images(:,:,1)=b1;

images(:,:,2)=b2;

[height,width,frames]=size(images);

%//initialzation of u and v

Vx = zeros(height,width);

Vy = zeros(height,width);

for k = 1:frames-1

% //initialization of Ex Ey and Et

Ex = zeros(height-1,width-1,frames-1);

Ey = zeros(height-1,width-1,frames-1);

Et = zeros(height-1,width-1,frames-1);

%//calculating Ex Ey and Et in frame k.

for x = 2:width-1

for y = 2:height-1

Ex(y,x,k) = (images(y+1,x+1,k)-images(y+1,x,k)+images(y,x+1,k)-images(y,x,k)+images(y+1,x+1,k+1)-images(y+1,x,k+1)+images(y,x+1,k+1)-images(y,x,k+1))/4;

Ey(y,x,k) = (images(y,x,k)-images(y+1,x,k)+images(y,x+1,k)-images(y+1,x+1,k)+images(y,x,k+1)-images(y+1,x,k+1)+images(y,x+1,k+1)-images(y+1,x+1,k+1))/4;

Et(y,x,k) = (images(y+1,x,k+1)-images(y+1,x,k)+images(y,x,k+1)-images(y,x,k)+images(y+1,x+1,k+1)-images(y+1,x+1,k)+images(y,x+1,k+1)-images(y,x+1,k))/4;

end

end

for nn = 1:5

for x = 2:width-1

for y = 2:height-1

Vxbar = (Vx(y-1,x)+Vx(y,x+1)+Vx(y+1,x)+Vx(y,x-1))/6+(Vx(y-1,x-1)+Vx(y-1,x+1)+Vx(y+1,x+1)+Vx(y+1,x-1))/12;

Vybar = (Vy(y-1,x)+Vy(y,x+1)+Vy(y+1,x)+Vy(y,x-1))/6+(Vy(y-1,x-1)+Vy(y-1,x+1)+Vy(y+1,x+1)+Vy(y+1,x-1))/12;

%// chapter 12 of Horn’s paper

temp = (Ex(y,x,k)*Vxbar+Ey(y,x,k)*Vybar+Et(y,x,k))/(0.5^2 + Ex(y,x,k)^2 + Ey(y,x,k)^2);

%update u and v

Vx(y,x) = Vxbar-Ex(y,x,k)*temp;

Vy(y,x) = Vybar-Ey(y,x,k)*temp;

end

end

end

end

quiver(Vx,Vy);

[…] MATLAB examples https://chi3x10.wordpress.com/2008/06/09/optical-flow-in-matlab/ […]

heloo guys..i am doing some work to detect motion of a vehicle in opposite direction using matlab..can any body help me regarding this!

Thanks extremely helpful. Will certainly share site with my pals

With havin so much content and articles do you ever run into

any problems of plagorism or copyright infringement? My site

has a lot of completely unique content I’ve either authored myself

or outsourced but it looks like a lot of it is popping it up all over the web without

my permission. Do you know any methods to help stop content from being ripped off?

I’d really appreciate it.

Usually I do not read article on blogs, but I would like to say that this write-up very pressured me to take a look at and do it! Your writing taste has been surprised me. Thanks, quite nice post.|