%% Zadanie przykładowe 1.
clear all; close all; clc;
im = imread('tire.tif');
[S L] = wavedec2(im, 2, 'bior3.3');
S1 = S;
L1 = L(:,1).*L(:,2);
% A2 H2 V2 D2 H1 V1 D1
% L1(1) 3xL1(2) 2xL1(3)

S1(1 + L1(1) + 2*L1(2) : L1(1) + 3*L1(2)) = 0;
S1(1 + L1(1) + 3*L1(2) + 2*L1(3) : L1(1) + 3*L1(2) + 3*L1(3)) = 0;

im_rec = waverec2(S1, L, 'bior3.3');
subplot(211); imshow(im); title('Przed');
subplot(212); imshow(im_rec); title('Po');

%% Zadanie przykładowe 2 & zad2 z https://drive.google.com/drive/u/1/folders/0B7iKO8QPDijJd090bkZzczB4SVE
clear all; close all; clc;
im = phantom(256);
im2 = zeros(210, 210);

[t, r] = size(im2);

for k = 1:21
    x = ((k*10)-10)+1;
    y = 1;

    im2(y : 210, x : x + 9) = k / 21;
end

krok = [1 2 3 5 8 10 12, 15];
N = length(krok);

blad = zeros(N, 1);
for k = 1:N
    kat = 0:krok(k):180;
    [R, xp] = radon(im2, kat);
    im_new = iradon(R, kat);

    for kz = 1:t
        for kx = 1:r
            blad(k) = blad(k) + (im(kz, kx) - im_new(kz, kx)).^2;
        end
    end
    blad(k) = sqrt(blad(k));
end

plot(krok, blad, 'r'); xlabel('Krok'); ylabel('Błąd');

%% Zadanie przykładowe 3
clear all; close all; clc;
im = imread('shape2.bmp');
im2 = im(:,:,1) < 255 & im(:,:,2) < 255 & im(:,:,3) < 255;

aseg = bwlabel(im2);
obw = (aseg == 5) | (aseg == 337) | (aseg == 17);
obw = imfill(obw, 'holes');
im3 = im2 | obw;

im3 = imclose(im3, [1 1; 1 1]);
aseg = bwlabel(im3);
N = max(max(aseg));
props = regionprops(im3, 'all');

ob_z_dziura = zeros(size(im3));
prostokat = zeros(size(im3));
gwiazdka = zeros(size(im3));

for k = 1:N
    if props(k).EulerNumber == 0
        ob_z_dziura = ob_z_dziura | aseg == k;
    end

    pole = bwarea(aseg == k);

    if abs((props(k).BoundingBox(3) * props(k).BoundingBox(4)) / pole - 1) < 0.01
        if abs(props(k).BoundingBox(3) - props(k).BoundingBox(4)) > 1
            prostokat = prostokat | aseg == k;
        end
    end

    UpEx = props(k).Extrema(1) - props(k).Extrema(2);
    LeftEx = props(k).Extrema(15) - props(k).Extrema(16);

    if abs(UpEx) <= 1 && abs(LeftEx) <= 3
        gwiazdka = gwiazdka | aseg == k;
    end
end

im3 = 255 * uint8(im3);
im3 = cat(3, im3, im3, im3);

prostokat = 255 * uint8(prostokat);
gwiazdka = 255 * uint8(gwiazdka);
ob_z_dziura = 255 * uint8(ob_z_dziura);

[z, x, dim] = size(im3);

for kz = 1:z
    for kx = 1:x
        if im3(kz, kx, 1) == 255;
            im3(kz, kx, 1) = 0; im3(kz, kx, 2) = 0; im3(kz, kx, 3) = 255;
        end

        if im3(kz, kx, 3) == 0;
            im3(kz, kx, 1) = 255; im3(kz, kx, 2) = 255; im3(kz, kx, 3) = 255;
        end

        if ob_z_dziura(kz, kx, 1) == 255;
            im3(kz, kx, 1) = 0; im3(kz, kx, 2) = 255; im3(kz, kx, 3) = 0;
        end

        if prostokat(kz, kx, 1) == 255;
            im3(kz, kx, 1) = 255; im3(kz, kx, 2) = 0; im3(kz, kx, 3) = 0;
        end

        if gwiazdka(kz, kx, 1) == 255;
            im3(kz, kx, 1) = 255; im3(kz, kx, 2) = 255; im3(kz, kx, 3) = 0;
        end
    end
end

subplot(121); imshow(im);
subplot(122); imshow(im3);

%% Zadanie przykładowe 4
clear all; close all; clc;
vid = zeros(256, 256, 3, 20, 'uint8');
% Środek koła głównego to 128, 128, r = 50
% Małego koła środek to godz 12. -> 178,128
% Małego koła środek końcowy to godz 3. -> 128, 178
% Kąt 270
krok = 255/19;
vid(178, 128, 3, 1) = 1;
for k = 1:20
    kat = -(90 + (k - 1) * 270/19);
    x = round(128 + 50*cosd(kat)); % Środek + promień(kąt), gdzie kąt określa, gdzie w odległości 50px od środka zaczynamy. 0 stopni jest na dodatniej półosi x.
    z = round(128 + 50*sind(kat));
    vid(z, x, 2, k) = (255 - (k - 1) * krok);
    vid(z, x, 1, k) = (0 + (k - 1) * krok);
    vid(:, :, 1:2, k) = imdilate(vid(:, :, 1:2, k), strel('disk', 5 + 2 * (k-1)));
end

implay(vid);

%% Zadanie przykładowe 5
clear all; close all; clc;
im = imread('rice.png');
im2 = imresize(im, [256 * 1.5, 256]); % Ta funkcja używa bikubicznej defaultowo.
im3 = imclose(im2, [1 1 1 1 1; 1 1 1 1 1; 1 1 1 1 1; 1 1 1 1 1; 1 1 1 1 1]);
im4 = imadjust(im3, [0 1], [0.15 0.9]); % Nie wiem, czy to dobrze.

[col, row] = size(im4);
biale = 0;
czarne = 0;

im5 = im4 > 128; % Automatycznie zbinaryzuj?
for kz = 1:col
    for kx = 1:row
        if im5(kz, kx) == 0
            czarne = czarne + 1;
        else
            biale = biale + 1;
        end
    end
end

ratio = biale/czarne; % Pikseli białych jest 2x mniej niż czarnych, ale nie wiem, czy o to chodziło. Raczej niezbyt.

subplot(221); imshow(im);
subplot(222); imshow(im2);
subplot(223); imshow(im3);
subplot(224); imshow(im4);
imshow(im3);
%% Zadanie przykładowe 6
clear all; close all; clc;
im = imread('peppers.png');
im = double(im)/255;

maska = [3 5];
red = im(:,:,1);
green = im(:,:,2);
blue = im(:,:,3);

redf = medfilt2(red, maska);
greenf = medfilt2(green, maska);
bluef = medfilt2(blue, maska);

[col, row, dim] = size(im);
skladowe = zeros(col, 2);
% 1PKT: 1,1
% OstPKT: 384 512
% PolPKT: 192 256
bin = false(col, row);

for kz = 1:col
    for kx = 1:row
        if (kz > kx)
            skladowe(kz, 1) = im(kz, kx, 3);
            skladowe(kz, 2) = im(kz, kx, 2);
        end

        if ((im(kz, kx, 3) > im(kz, kx, 1 )) || (im(kz, kx, 3) > im(kz, kx, 1)))
            bin(kz, kx) = 1;
        end
    end
end

er = imerode(bin, strel('disk', 5));

% Momenty bezwładności:
% I MOMENT
aseg = bwlabel(er);
N = max(max(aseg));

pole_pow = 0; % Nie mam pojęcia, czy to ma być dla wszystkich figur czy całego obrazka.

for k = 1:N
    pole_pow = pole_pow + bwarea(aseg == k);
end

x = 0;
y = 0;
for kz = 1:col
    for kx = 1:row
        if er(kz, kx) == 1
            x = x + kx;
            y = y + kz;
        end
    end
end

M1x = (1/pole_pow)*x;
M1y = (1/pole_pow)*y;

% II MOMENT
x = 0;
y = 0;
for kz = 1:col
    for kx = 1:row
        if er(kz, kx) == 1
            x = x + (kx - M1x)^2;
            y = y + (kz - M1y)^2;
        end
    end
end

M2x = 1/pole_pow*x;
M2y = 1/pole_pow*y;

%% Zadanie przykładowe 7
clear all; close all; clc;
im = imread('pout.tif');
im = double(im)/255;
im = imadjust(im);

IM = fftshift(fft2(im));

[col, row] = size(im);
maska = zeros(size(im));
maska(290/2, 240/2) = 1;
maska = bwdist(maska, 'euclidean');

LP = 1.0 ./ (1 + (maska / 20).^10);
IM2 = IM .* LP;
im_new = ifft2(ifftshift(IM2));

imshow(im_new);

bin = false(col, row);
licznik = 0;
for kz = 1:col
    for kx = 1:row
        if im_new(kz, kx) >= 0.2 && im_new(kz, kx) <= 0.25
            licznik = licznik + 1;
            bin(kz, kx) = 1;
        end
    end
end

aseg = bwlabel(bin);
N = max(max(aseg(:)));

%%
clear all; close all; clc;
im = imread('new-shape.bmp');

bin = (im > 50);
aseg = bwlabel(bin);
N = max(aseg(:));
props = regionprops(bin, 'all');

x = zeros(1, N);
y = zeros(1, N);

for k = 1:N
    fig = (aseg == k);
    pole = bwarea(fig);
    obwod = bwarea(bwperim(fig));
    bwk = (4*pi*pole)/(obwod*obwod);
    if abs(bwk - 1) < 0.05
        x(k) = props(k).BoundingBox(1) - 0.5;
        y(k) = props(k).BoundingBox(2) - 0.5;
    end
end

bin2 = ~bin;
aseg2 = bwlabel(bin2);
N2 = max(aseg2(:));
dziura_kolo = false(size(bin2));

for k = 1:N2
    fig = aseg2 == k;

    for kx = 1:7
        r = x(1, kx);
        c = y(1, kx);

        if c == 0 && r == 0
            continue;
        end

        if fig(c, r) == 1
            dziura_kolo = dziura_kolo | fig;
        end
    end
end

props2 = regionprops(bin2, 'all');
szesciokaty = zeros(size(bin2));

for k = 1:N2
    fig = aseg2 == k;

    temp = props2(k).Extrema(1) - props2(k).Extrema(2);
    if (abs(temp) <= 1)
        temp2 = props2(k).Extrema(5) - props2(k).Extrema(6);
        if (abs(temp2) <= 3)
            temp3 = props2(k).Extrema(11) - props2(k).Extrema(12);
            if (abs(temp3) > 90)
                szesciokaty = szesciokaty | fig;
            end
        end
    end
end

subplot(211); imshow(szesciokaty);
subplot(212); imshow(dziura_kolo);

%% 4 1b.jpg
clear all; close all; clc;
im = imread('pout.tif');
im = imadjust(im);
im = double(im)/255;

im_n = imnoise(im, 'salt & pepper', 0.02);
im_n = imnoise(im_n, 'gaussian', 0.0, 0.001);

K = 2;
M = 3;

im_new = wiener2(im_n, [K, M]);
im_new = medfilt2(im_new, [K, M], 'symmetric');

imshow(im_new);

% Obliczanie błędu:
oc = sqrt(sum(im(:) - im_new(:)).^2);