clear; close all; clc;

img = imread('pavian.png'); % Choose any picture you like here
img = rgb2gray(img);

%% Display image
figure; subplot(131);
imshow(img);
title('Ausgangsbild');

%% FFT
fft2d = fft2(img);

% Display the FFT image.
subplot(132);
imshow(log(fftshift(abs(fft2d))), []);
title('2D-FFT');

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

% Approach: We want to remove the low frequencies (high pass filter) to
% keep the areas of the image where the lightness changes relative to
% nearby pixels

% Shift FFT so low frequencies are in the center
F = fftshift(fft2d);

% Get image size
[M, N] = size(F);

% Radius of how much to remove
r = 25;

% We cut out the center of the fft low frequencies
H = ones(M, N);
H(M/2-r:M/2+r, N/2-r:N/2+r) = 0;

% Apply filter
F_filtered = F .* H;

% Shift back
fft2d = ifftshift(F_filtered);

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

%% IFFT to obtain filtered image
img_filtered = ifft2(fft2d);
subplot(133);
imshow(real(img_filtered), []);
title('Gefiltertes Ausgangssignal');