Ghosh and Jain (1) showed how to charactize a simple, closed, piece-wise linear curve as a Fourier series in a single parametric variable. Such a formulation allows (see Tsow and Land below):

• interpolation of curves by sampling the parameter more densely.
• morphing of one curve into another by locating matching parametric points on the two curves and smoothly interpolating the Fourier components.
• simplification of curves by removing Fourier components. This maintains the overall shape, but eliminates sharp corners.
• location of the curve in the plane as the DC components of the series.
• ability to add band-limited noise to a curve.
• ability to connect serveral curves (possibly representing serial sections) to make a 3D object.

We intend to use the Fourier fit scheme as a way to simplify and characterize the outlines of moving animals. The animals will be videotaped. The videotape will then be analysed to extract the outline of the animal, Fourier fit the outline, then use the low bandwidth components to characterize the motion. For example, the DC components give the position of the animal. The fundamental frequency coefficients of the fit give the ellipse (2) which best fits the animal outline.

Code

A matlab program was written to animate a test object, extract its outline, then do the Fourier fit.

This program Fourier fits a square:

%Fourier fit of motion

clear all
figure(1)
clf; clc
set(gcf,'position',[50,600,400,200])
axis([0 400 0 200])
% axis tight

set(gcf,'units','pixels');
set(gca,'units','pixels');
w_pos = get(gcf, 'position');
set(gca, 'position', [0 0 w_pos(3) w_pos(4)]);

%number of Fourier components to use
%vary this to make a better fit
N = 11;

%make an rectangle

figure(1)
cla;axis off

rectangle('Position', [30 30 100 100]);
im = getframe(gcf);
ed = im.cdata;
% figure(2)
% imshow(ed)

%order the edge pixels
%find them in the array
[xpos,ypos] = find(ed(:,:,1)==0);
npts = length(xpos);
%keep track of visit order
ordering=[1:npts] * 0;
visited=[1:npts] * 0;
%pick the first point on the edge and calculate the normal
%move the next pixel which is:
%--on the edge
%--closeest
%now just move to the nearest pixel which has not been visited
nextindex=1; ordering(1)=1; visited(1)=1;
count = 1;
while (count<npts)
currentindex=nextindex;
currentmin=inf;
for i=1:npts
d = (xpos(currentindex)-xpos(i)).^2 + (ypos(currentindex)-ypos(i)).^2 ;
if (d<currentmin) & (visited(i)==0)
currentmin = d;
nextindex = i;
end
end
%at this point, we have the index to the nearest point
visited(nextindex)=1;
count = count + 1;
ordering(count)=nextindex;
end

%define a piece-wise linear, closed curve
clear p
p(:,1) = ypos(ordering);
p(:,2) = xpos(ordering);
K = length(p);
%close the loop
p(K+1,:) = p(1,:);

%get fourier transform of the closed curve
%Ghosh and Jain, IEEE conputer graphics and appl, 1994
dx = diff(p(:,1));
dy = diff(p(:,2));
dL = sqrt( dx.^2 + dy.^2 );
L = sum(dL);
Lcum = cumsum(dL);

%get a0, the x-offset
temp = p(1,1)*Lcum(1)+ dx(1)*Lcum(1)^2/dL(1)*0.5; %for i=1
for i=2:K
temp = temp + ...
(p(i,1)-dx(i)*Lcum(i-1)/dL(i))*(Lcum(i)-Lcum(i-1)) + ...
0.5*dx(i)/dL(i)*(Lcum(i)^2 - Lcum(i-1)^2);
end
a0=temp/L;

%get c0, the y-offset
temp = p(1,2)*Lcum(1)+ dy(1)*Lcum(1)^2/dL(1)*0.5; %for i=1
for i=2:K
temp = temp + ...
(p(i,2)-dy(i)*Lcum(i-1)/dL(i))*(Lcum(i)-Lcum(i-1)) + ...
0.5*dy(i)/dL(i)*(Lcum(i)^2 - Lcum(i-1)^2);
end
c0=temp/L;

for n=1:N
tempa = dx(1)/dL(1)*(cos(2*pi*n*Lcum(1)/L)-1);
tempb = dx(1)/dL(1)*sin(2*pi*n*Lcum(1)/L);
tempc = dy(1)/dL(1)*(cos(2*pi*n*Lcum(1)/L)-1);
tempd = dy(1)/dL(1)*sin(2*pi*n*Lcum(1)/L);
for i=2:K
tempa = tempa + ...
dx(i)/dL(i)*(cos(2*pi*n*Lcum(i)/L)-cos(2*pi*n*Lcum(i-1)/L));
tempb = tempb + ...
dx(i)/dL(i)*(sin(2*pi*n*Lcum(i)/L)-sin(2*pi*n*Lcum(i-1)/L));
tempc = tempc + ...
dy(i)/dL(i)*(cos(2*pi*n*Lcum(i)/L)-cos(2*pi*n*Lcum(i-1)/L));
tempd = tempd + ...
dy(i)/dL(i)*(sin(2*pi*n*Lcum(i)/L)-sin(2*pi*n*Lcum(i-1)/L));
end
a(n) = tempa * L /(2*pi^2*n^2);
b(n) = tempb * L /(2*pi^2*n^2);...