flow1

Download MATLAB source: flow.zip.

% calculate the scene - a field of points
clear
close all
xr = -20:20;
yr = -20:60;
x = ones(length(yr),1)*xr;
y = yr'*ones(1,length(xr));
z = zeros(size(x));
disp('matrix size');
disp(size(x));

matrix size
    81    41

% draw the scene
figure(1);
plot3(x,y,z,'ko');
axis equal
axis off
camproj('perspective');
campos([0 -20 5]);
camtarget([0 -10 0]);
camva(30);
pos = campos();
tgt = camtarget();

% move the camera
delta = [0 0.2 0];
aviobj = avifile('flow2.avi','fps',20,'compression','Cinepak');
for k=1:101
    campos(pos+(k-1)*delta);
    camtarget(tgt+(k-1)*delta);
    drawnow;
    F = getframe(gcf);
    aviobj = addframe(aviobj,F);
end
aviobj = close(aviobj);

% do flow without MATLAB 3D graphpics


figure(2)
up = [0 0 1];
vn = pos - tgt;
vz = vn/norm(vn);
vx = cross(up,vz);
vx = vx/norm(vx);
vy = cross(vz,vx);
R = [vx' vy' -vz'];
T2 = [R  -R*pos'];

% test
disp('transform target point:');
disp(T2*[tgt 1]')

xp = x(:);
yp = y(:);
zp = z(:);

xform = T2*[xp'; yp'; zp'; ones(size(zp'))];
va = 30;
scl = 1/tan(va*pi/180);
idx = find(xform(3,:)>0); % select points in front of camera
xpts = xform(1,idx)./xform(3,idx)*scl;
ypts = xform(2,idx)./xform(3,idx)*scl;

plot(xpts,ypts,'ko');
axis([-0.7 0.7 -0.5 0.5]);

transform target point:
         0
         0
   11.1803

% traditional optical flow diagram
yp = y(:) - 0.2;
xform = T2*[xp'; yp'; zp'; ones(size(zp'))];
xnew = xform(1,idx)./xform(3,idx)*scl;
ynew = xform(2,idx)./xform(3,idx)*scl;

s = 80;
quiver(xpts,ypts,xnew-xpts,ynew-ypts,5,'k');
axis([-0.7 0.7 -0.5 0.5]);