function [img] = imFilterFFT_lowpass_RESENI( image, sigma, show )
% =========================================================================
% [img] = imFilterFFT_lowpass( image, sigma, show )
% 
% Funkce provede dolnopropustni filteraci s pouzitim Fourier transformace
% 
% PARAMETRY:
% * image - vstupni sedotonovy obrazek o velikosti [n m], kde n, m je sirka a
%           vyska
% * sigma - rozptyl Gaussova rozlozeni ve frekvencni oblasti
% 
% VYSTUP:
% * img - filtrovany obrazek o stejne velikosti jako image
% 
% PRIKLAD:
%     im = imread('ct_brain.jpg');
%     im = rgb2gray(im);
%     img = imFilterFFT_lowpass( im, 20 );
%     figure;
%     subplot(121); imshow(im, []);
%     subplot(122); imshow(img, []);
% -------------------------------------------------------------------------

if nargin<2,  error('too few parameters');  end;
if nargin<3,  show=0;                       end;

    % preved obrazek do double
    im = double(image);
    % vytvor Fourierův obraz vstupniho obrazku
    image_fft = fft2(image);

    % vytvorte matici G (Gauss filter) o stejné velkosti jako FT obrazku
    filter_fft = zeros(size(image_fft));

    % vytvorte pomocnou omatici, která bude obsahovat vzdaalenosti od
    % stredu matice G vytvorenou vyse pouzijte funkci ndgrid
    [sz1, sz2] = size(filter_fft);
    fft_range1 = (-sz1/2):(sz1-sz1/2-1);
    fft_range2 = (-sz2/2):(sz2-sz2/2-1);
    [ids1,ids2] = ndgrid(fft_range1,fft_range2);
    dist = sqrt(ids1.^2 + ids2.^2);

    % spoctete Gaussovo rozdeleni pomoci predpocitanych vzdalenosti a danym
    % paramtrem sigma -> e^(-2*pi*x^2/sigma^2)
    filter_fft(:,:) = (sigma*sqrt(2*pi))^(-1) * exp(-(dist).^2/(2*sigma^2));
%     filter_fft(:,:) = exp(-2*pi*(dist).^2/sigma^2);

    % zkonvertujeme spektrum z centrovaneho formátu nazpet
    filter_fft = ifftshift(filter_fft);

    % provedeme filtraci nasobenim
    image_filt_fft = image_fft .* filter_fft;

    % rekonstrujeme vysledny obraz po filtraci pomocí IFFT
    img = ifft2(image_filt_fft);
    img = real(img);
    
    if show,
        figure;
        subplot(221), imshow(im, []);
        subplot(223), colormap(getFFTmag(image_fft)); colorbar;
        subplot(222), imshow(img, []);
        subplot(224), colormap(getFFTmag(image_filt_fft)); colorbar;
    end;
    
    function [gray_log] = getFFTmag(image_fft)

        % draw centered Fourier spectrum (log of power)
        base_log = 1e-2;
        image_fftlog = log(base_log + abs(image_fft).^2)-log(base_log);
        fftim_sz1 = size(image_fft,1);
        fftim_sz2 = size(image_fft,2);
        fft_range1 = (-fftim_sz1/2):(fftim_sz1-fftim_sz1/2);
        fft_range2 = (-fftim_sz2/2):(fftim_sz2-fftim_sz2/2);
        imagesc(fft_range1,fft_range2, fftshift(image_fftlog));

        % create colormap
        gray_log = (log(base_log + gray.^2) - log(base_log));
        gray_log = gray_log ./ (log(base_log + 1) - log(base_log));
