Contents

%function map = find_edges(im,zm)

find edges

m0 = [ 1 2 1; 2 4 2; 1 2 1];
m1 = [ -1 0 1; -2 0 2; -1 0 1];
m2 = [ -1 -2 -1; 0 0 0; 1 2 1];
a0 = imfilter(im,m0);
a1 = imfilter(im,m1);
a2 = imfilter(im,m2);
slope = a1.*a1+a2.*a2;
[r c] = size(im);
slope(1,:)=0;
slope(:,1)=0;
slope(r,:)=0;
slope(:,c)=0;
imshow(slope);
d = 16*zm-a0;

imshow(d);
idx = find(abs(d)<sqrt(slope));
x = ones(r,1)*(1:c);
y = (1:r)'*ones(1,c);
map = zeros(r,c);
map(idx)=1;
u = x(idx) + (a1(idx)./slope(idx)).*d(idx)/2;
v = y(idx) + (a2(idx)./slope(idx)).*d(idx)/2;
x = x(idx);
y = y(idx);
clear m0 m1 m2;
clear r c idx;

edge linking

ymin = min(y);
ymax = max(y);
nseg = 0;
label = zeros(size(x));
next = zeros(size(x));
merge=3;
for yi=ymin:ymax
    ndx = find(y==yi);
for j=1:length(ndx)
    connect=0;
    k=ndx(j);
    ulink = 0;
    vlink = 0;
    umin=1.5;
    vmin=1.5;
    for n=1:nseg
        % just look at open edges
        ka = seg(n).a;
        kb = seg(n).b;

        if (seg(n).length>0 && kb>0)
            dx1 = u(ka)-u(k);
            dy1 = v(ka)-v(k);
            dist = sqrt(dx1*dx1+dy1*dy1);
            if (dist<umin)
                ulink = n;
                umin = dist;
            end
            dx2 = u(kb)-u(k);
            dy2 = v(kb)-v(k);
            dist = sqrt(dx2*dx2+dy2*dy2);
            if (dist<vmin)
                vlink = n;
                vmin = dist;
            end
        end
    end

    % append to existing edge
    if (vlink>0)
        connect=connect+1;
        kb = seg(vlink).b;
        next(kb)=k;
        seg(vlink).b = k;
        label(k)=vlink;
        seg(vlink).length = seg(vlink).length + 1;
    elseif (ulink>0)
        connect=connect+1;
        %fprintf('prepending to segment %d (%g)\n',ulink,umin);
        ka = seg(ulink).a;
        next(k)=ka;
        seg(ulink).a = k;
        label(k)=ulink;
        seg(ulink).length = seg(ulink).length + 1;
    end
    % join two existing edges
    if ((ulink>0) && (vlink>0))
        if (ulink~=vlink)
            %fprintf('line %d ',yi-ymin+1);
            %fprintf('linking segments %d [%d] and %d [%d]\n',ulink,seg(ulink).length,vlink,seg(vlink).length);
            %fprintf('segment %d a %d b %d length %d\n',ulink,seg(ulink).a,seg(ulink).b,seg(ulink).length);
            %fprintf('segment %d a %d b %d length %d\n',vlink,seg(vlink).a,seg(vlink).b,seg(vlink).length);
            %fprintf('link element: %d\n',k);
            kdx = find(label==vlink);
            label(kdx)=ulink;
            seg(ulink).length = seg(ulink).length+seg(vlink).length;
            seg(vlink).length = 0;
            next(k)=seg(ulink).a;
            seg(ulink).a = seg(vlink).a;
        end
    end
    % create a new segment
    if (connect==0)
        nseg = nseg + 1;
        seg(nseg) = struct('a',{k},'b',{k},'length',{1},'color',{'b'});
        label(k)=nseg;
    end
end

eliminate empty segments

if (mod(yi-ymin,4)==0)
    merge_edges;
else
    trim=0;
    trim_edges;
end
segment 14 is self-linking
segment 14 is self-linking
end
clear dx1 dx2 dy1 dy2;
clear connect dist j k len n ka kb kdx ndx ulink vlink;
clear x y yi ymin ymax umin vmin;
merge_edges;
show_segments(seg,nseg);
%imshow(map);
draw_edgels(seg,nseg,u,v,next);
segment 29 is self-linking
   n    a    b  len
   1 1669 1769  221
   2 1933 1937  223
   3 2204 2208  220
   4 2474 2575  221
   5 2737 2843  222
   6 3010 3093  222
   7 3276 3382  221
   8 3545 3550  223
   9 3846 3918  221
  10 4080 4161  221
  11  860  905  157
  12 1214 1229  121
  13  258  275   32
  14 4344 4420  909
  15  178  313   32
  16  171  305   32
  17  942 1002  191
  18  246  334  119
  19 1289 1296   66
  20  254  315    9
  21  954  995    6
  22  712  807   95
  23   75  416  121
  24 1219 1258   42
  25 1169 1221   41
  26  676  793  102
  27  774  774    1
  28   69  407  129

32 allocated segments
4420 total edgels