function d = nmsdir(f)
%
% nmsdir: directional non-maxima suppression. Retains only those pixels
%         that are a max in 1 of the 4 major directions
% On input:
%     f:  2D image function
% On output:
%     d:  value of f at each local directional max
%

%k2 = (k-1)/2;
[m,n] = size(f);
d = zeros(m,n);
f1 = f.*(f>0);
%for i1 = 1+k2:m-k2
%  for j1 = 1+k2:n-k2
for i1 = 2:m-1
  for j1 = 2:n-1
    val = f1(i1,j1);
	if (val>f1(i1,j+1)) & (val>f1(i1,j1-1))
	  d(i1,j1) = val;
	elseif (val>f1(i1-1,j+1)) & (val>f1(i1+1,j1-1))
	  d(i1,j1) = val;
	elseif (val>f1(i1-1,j1)) & (val>f1(i1+1,j1))
	  d(i1,j1) = val;
	elseif (val>f1(i1-1,j1-1)) & (val>f1(i1+1,j1+1))
	  d(i1,j1) = val;
	end
  end
end
