%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%% Part 1. Read and plot the original image %%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Read the original image
org_pf = imread("platform.png");
% Convert the image to grayscale and normalize it to the range [0, 1]
org_gray_pf = mean(org_pf, 3) / 255;
[h, w] = size(org_gray_pf);
figure;subplot(2, 3, 4);imshow(org_gray_pf);title("Original Image");
 
org_freq = fftshift(fft2(org_gray_pf)); % The same as circshift(fft2(org_gray_pf), [h/2, w/2]), Shifting the low frequency components to the center. Our definition is equal to (1/(h^2)) * fft2 

% Normalize the frequency component for visualization.
log_org_freq = log(1+abs(org_freq));
max_log_freq = max(max(log_org_freq));
subplot(2, 3, 1);imshow(log_org_freq / max_log_freq);title("Frequency of original image");

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%% Part 2. Read and plot the corrupted image %%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Read the image
pf = imread("noisy_platform.png");
% Convert the image to grayscale and normalize it to the range [0, 1]
gray_pf = mean(pf, 3) / 255;
subplot(2, 3, 5);imshow(gray_pf);title("Corrupted Image");

freq = fftshift(fft2(gray_pf));

% Normalize the frequency component for visualization.
log_freq = log(1+abs(freq));
subplot(2, 3, 2);imshow(log_freq / max_log_freq);title("Frequency of corrupted image");


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%% Part 3. Recover the corrupted image %%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Compute the distance.
[X, Y] = meshgrid([-w/2:w/2-1], [-h/2:h/2-1]); 
dist = X.^2 + Y.^2;
% Define the Gaussian Low Pass Filter.
sigma = 100;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Complete the following lines to define the gaussian low pass filter.
filter = %?????;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Multiply the frequency components by the Gaussian Low Pass Filter.
freq = freq .* filter;

% Visualize the resulting frequency domain.
subplot(2, 3, 3);imshow(log(1+abs(freq)) / max_log_freq);title("Frequency of recovered image");

% Transform back to the spatial domain.
recon = ifft2(ifftshift(freq));

% Plot the recovered image.
subplot(2, 3, 6);imshow(real(recon));title("Recovered Image");


