function [angle,c] = ctopo2(f)
%
% ctopo2: compute the curl of the 90 degree rotated gradient vectors
% On input:
%   f:  2D image function
% On output:
%   c:  curl of the 90 degree rotated gradient (equivalent to the negative
%       Laplacian
% Call:
%      [angle,c] = ctopo2(f)
%

[fc,fr] = gradient(f);
[fcc,fcr] = gradient(fc);
[frc,frr] = gradient(fr);

c = -fcc-frr;
angle = posori(atan2(fcr,fcc-frr)/2);
%cmask = c>(0.1*max(max(c)));
cmask = c>0;
[m,n] = size(c);
for i1 = 2:m-1
  for j1 = 2:n-1
    if cmask(i1,j1)>0 
      if (locmaxdir(f,i1,j1,angle(i1,j1)) ...
        | locmaxdir(f,i1,j1,rem((angle(i1,j1)+pi/2),2*pi)))
      else
        cmask(i1,j1) = 0;
      end;
    end
  end
end;

c = c.*cmask;
