function W = linear_kernel(width);

% to be called in REDUCE

N = (width-1)/2;

% returns a matrix descirbing the optimal 2-D linear kernel
% the width of the filter is width
%not the most efficient algorithm but who cares... it's just a few operations

%first step : compute the 1D optimal linear filter
w = zeros(1,width);
a = -3-sqrt(8);
for i=0:N/2
  w(2*i+N+1) = a^(-i);
  w(-2*i+N+1) = w(2*i+N+1);
end
if (mod(N/2,1)~=0)
  for i=0:N/2
    w(2*i+N+2) = 0.5*(a^(-i)+a^(-i-1));
    w(-2*i+N) = w(2*i+N+2);  
  end
else
  for i=0:N/2-1
    w(2*i+N+2) = 0.5*(a^(-i)+a^(-i-1));
    w(-2*i+N) = w(2*i+N+2);  
  end
end
% normalisation of the filter.
w = w/sqrt(2);
% assume the filter is separable... w(ij)=w(i)w(j)
W = w'*w;
%normalisation so that the energy in the system is conserved...
W = W/sum(sum(W));
