clear;clc;
prompt = "Proszę wprowadzić nazwę pliku (bez rozszerzenia)";
prompt = prompt + "\nPlik (typu .gif) pobrany zostanie z folderu Datasets\n\n";
txt = input(prompt, 's');
name = strcat('..\Datasets\', txt);
A=double(imread(name,'gif'));
n=size(A,1);

%
% Dwa progi e_high i e_low. Jeden służy do wstępnej klasyfikacji krawędzi,
% drugi do wtórnej, wśród sąsiadów już zakwalifikowanej krawędzi
%

e_high=70;
e_low=40;

%
% imgaussfilt to filtr Gaussa, zmiękczający. Drugi parametr reguluje
% głębokość zmiękczania, najczęściej stosuje się 1.5 ... 3. Odpowiada to
% wielokrotnym iteracjom splotu z maską Gaussa 5x5
%

B=imgaussfilt(A,1.5);

%
% Maski Prewitta do obliczania pochodnych cząstkowych
%

Dx=[-1 0 1;-1 0 1;-1 0 1];
Dy=Dx';

Gx=conv2(B,Dx,'same');
Gy=conv2(B,Dy,'same');

%
% G to moduł gradientu, T to zkwantyzowany do 4 wartości kierunek gradientu
%

G = abs(complex(Gx,Gy));
colormap(gray);
name = strcat('Images\', txt, '-gradient.gif');
imwrite(uint8(255 - G), name, 'gif');
T = angle(complex(Gx,Gy));
for i=1:n
  for j=1:n
    if T(i,j) < 0
        T(i,j )= T(i,j) + pi;
    end
    if T(i,j) <= pi/8 || T(i,j) >= 7 * pi / 8
        T(i,j) = 1;         %   kierunek poziomy
    elseif T(i,j) >= 5 * pi / 8
        T(i,j) = 2;         %   kierunek NE-SW
    elseif T(i,j) >= 3 * pi / 8
        T(i,j) = 3;         %   kierunek pionowy
    else
        T(i,j) = 4;         %   kierunek NW-SE
    end
  end
end

%
% E to będą krawędzie. Najpierw krawędziami będą lokalne maksima modułu
% gradientu, w kierunku gradientu. Daje to wstepne krawędzie, bardzo
% cienkie
%

E=zeros(n);

for i=2:n-1
  for j=2:n-1
    if T(i,j) == 1
      u = 0;
      v = 1;
    elseif T(i,j) == 4
      u = 1;
      v = 1;
    elseif T(i,j) == 3
      u = 1;
      v = 0;
    else
      u = -1;
      v = 1;
    end
    if (G(i,j) >= G(i+u,j+v)) && (G(i,j) >= G (i-u,j-v))
      E(i,j) = G(i,j);
    end
  end
end

%
% Zachowujemy maksima do wtórnego progowania
%

F = E;

%
% Wstępne progowanie krawędzi, wyższym progiem
%

for i=2:n-1
  for j=2:n-1
    if E(i,j) >= e_high
        E(i,j) = 255;
    else
        E(i,j) = 0;
    end
  end
end
figure(1);
imshow(255 - E); title('Gołe progowanie');
colormap(gray);
name = strcat('Images\', txt, '-edges-initial.gif');
imwrite(uint8(255 - E), name, 'gif');

%
% Ostatni etap to rekurencyjne dodawanie krawędzi, poprzez łagodniejsze
% progowanie sąsiadów już ustalonych krawędzi.
%

for k = 1 : 3
    E_temp = E;
    for i=2:n-1
        for j=2:n-1
            if E(i,j) == 255
                for u = -1:1
                    for v = -1:1
                        if F(i+u,j+v) >= e_low
                            E_temp(i+u,j+v)=255;
                        end
                    end
                end    
            end
        end
    end
    E = E_temp;
end 
figure(2);
imshow(255 - E); title('Po poprawce');
colormap(gray);
name = strcat('Images\', txt, '-edges-final.gif');
imwrite(uint8(255 - E), name, 'gif');
