function resultRGB=marchingsquares(img, K)
    resultRGB = img;

    % convert to grayscale
    grayimage = rgb2gray(img);

    binaryimage = grayimage;

    binaryimage(grayimage>K)=255;
    binaryimage(grayimage<=K)=0;

    imshow(binaryimage);

    x=0;
    y=0;
    
    [width, height] = size(binaryimage);
    
    x=width-1;
    y=height-1;

    for i = 1:width-1
        breakout = false;
        for j=1:height-1
            if(~(binaryimage(i,j,1) == binaryimage(i+1,j,1) && binaryimage(i+1,j,1)==binaryimage(i+1,j+1,1) && binaryimage(i+1,j+1,1)==binaryimage(i,j+1,1)))
                x=i; y=j;
                breakout = true;
                break;
            end
        end
        if(breakout)
            break;
        end
    end

    if(x==width-1 && y==height-1)
        error('Error. \nNo border found.');
    end

    startx=x;
    starty=y;

    finished = false;
    xdirection = 0;
    ydirection = 0;

    prevx = x;
    prevy = y;

    figure, image(img)
    %figure, image(binaryimage)
    axis off
    hold on
    
    while finished==false
        n=getContourNumbering(binaryimage, x,y);
        switch n
			% going RIGHT with these cases:
			% +---+---+		+---+---+   +---+---+   +---+---+
			% | 2 | 3 |		|   | 3 |   | 2 | 3 |   | 2 | 3 |
			% +---+---+		+---+---+   +---+---+   +---+---+
			% | 5 | 7 |		|   |   |   |   |   |   | 5 |   |
			% +---+---+		+---+---+  	+---+---+  	+---+---+
            case {1, 3, 6, 30}
               xdirection = 1;
               ydirection = 0;
            % going LEFT with these cases:
			%	+---+---+   +---+---+   +---+---+
			%	|   |   |   |   |   |   |   | 3 |
			%	+---+---+   +---+---+   +---+---+
			%	| 5 |   |   | 5 | 7 |   | 5 | 7 |
			%	+---+---+  	+---+---+  	+---+---+
            case {35, 105, 5}
               xdirection = -1;
               ydirection = 0;
            % going DOWN with these cases:
			%	+---+---+   +---+---+   +---+---+
			%	|   |   |   |   | 3 |   | 2 | 3 |
			%	+---+---+   +---+---+   +---+---+
			%	|   | 7 |   |   | 7 |   |   | 7 |
			%	+---+---+  	+---+---+  	+---+---+
            case {21, 7, 42}
                xdirection = 0;
                ydirection = 1;
            % going UP with these cases:
            %	+---+---+   +---+---+   +---+---+
            %	| 2 |   |   | 2 |   |   | 2 |   |
            %	+---+---+   +---+---+   +---+---+
            %	|   |   |   | 5 |   |   | 5 | 7 |
            %	+---+---+  	+---+---+  	+---+---+
            case {2, 10, 70}
                xdirection = 0;
                ydirection = -1;
            % special saddle point case 1:
			%	+---+---+ 
			%	|   | 3 | 
			%	+---+---+
			%	| 5 |   |
			%	+---+---+
			% going LEFT if coming from UP
			% else going RIGHT 
            case 15
                if (xdirection==0&&ydirection==-1)
                    xdirection=-1;
                    ydirection=0;
                else
                    xdirection=1;
                    ydirection=0;
                end
            % special saddle point case 2:
			%	+---+---+ 
			%	| 2 |   | 
			%	+---+---+
			%	|   | 7 |
			%	+---+---+
			% going UP if coming from RIGHT
			% else going DOWN
            case {14}
                if (xdirection==1 && ydirection==0)
                    xdirection = 0;
                    ydirection = -1;
                else
                    xdirection=0;
                    ydirection=1;
                end
            otherwise % case n = 210, which is ambiguous
               % dont change direction and just keep going
        end
        x = x+xdirection;
        y = y+ydirection;

        line([prevy,y],[prevx,x],'Color','r','LineWidth',2)

        prevx=x;
        prevy=y;

        if (x==startx && y==starty)
            finished = true;
        end
    end

    hold off
end

function result=getContourNumbering(img, x, y)
    result = 1;
    if(img(x,y,1)==255) result = result*2; end
    if(img(x+1,y,1)==255) result = result*3; end
    if(img(x,y+1,1)==255) result = result*5; end
    if(img(x+1,y+1,1)==255) result = result*7; end
end