close all;
clear all;

knee = double(imread('knee.png'));


figure()
imshow(knee,[]);

[X,Y] = ginput(1);
X = floor(X);
Y = floor(Y);

stack = zeros(1000,2);
visited = zeros(size(knee));
segmented = zeros(size(knee));
top = 1;

stack(top,1) = X;
stack(top,2) = Y;
visited(X,Y) = 1;
segmented(X,Y) = 1;

prog = 4;

while top ~= 0
    X = stack(top,1);
    Y = stack(top,2);
    top = top - 1;
    neighbours = [X-1 Y+1; X-1 Y; X-1 Y-1; X Y-1; X+1 Y-1; X+1 Y; X+1 Y+1; X Y+1];
    for i = 1:8
        if visited(neighbours(i,1),neighbours(i,2)) == 0
            measure = abs(knee(X,Y)-knee(neighbours(i,1),neighbours(i,2)));
            if measure < prog
             segmented(neighbours(i,1),neighbours(i,2)) = 1;
             top = top + 1;
             stack(top,1) = neighbours(i,1);
             stack(top,2) = neighbours(i,2);
            end
        end
        visited(neighbours(i,1),neighbours(i,2)) = 1;
    end
end

figure();
imshow(segmented,[]);

%%%%%%%%%%%%%%%%%%%

figure();
imshow(knee,[]);

knee = imfilter(knee,fspecial('gaussian',7,7));

figure();
imshow(knee,[]);

[X,Y] = ginput(1);
X = floor(X);
Y = floor(Y);

Icopy = knee;

stack = zeros(1000,2);
visited = zeros(size(knee));
segmented = zeros(size(knee));
top = 1;

stack(top,1) = X;
stack(top,2) = Y;
visited(X,Y) = 1;
segmented(X,Y) = 1;

prog = 15;

n = 1;
mean = knee(X,Y);

while top ~= 0
    X = stack(top,1);
    Y = stack(top,2);
    top = top - 1;
    n = n+1;
    neighbours = [X-1 Y+1; X-1 Y; X-1 Y-1; X Y-1; X+1 Y-1; X+1 Y; X+1 Y+1; X Y+1];
    for i = 1:8
        if visited(neighbours(i,1),neighbours(i,2)) == 0
            mean = (mean*(n-1)+knee(X,Y))/n;
            measure = abs(mean-knee(neighbours(i,1),neighbours(i,2)));
            if measure < prog
             segmented(neighbours(i,1),neighbours(i,2)) = 1;
             top = top + 1;
             stack(top,1) = neighbours(i,1);
             stack(top,2) = neighbours(i,2);
            end
        end
        visited(neighbours(i,1),neighbours(i,2)) = 1;
    end
    if mod(n,800) == 0
        imshow(Icopy.*segmented);
        drawnow;
    end
end

figure();
imshow(segmented,[]);