% Fourier Transform Shifting Operator Demo
close all; clear all;
% Define Initial box
x0 = 97; 
y0 = 97;
xwidth = 64;
ywidth = 64 ;

% Create a white box on a black background image
M= 256; N = 256;
image1 = zeros(M,N);
box = ones(xwidth,ywidth);
%box at centre
image1(x0:x0+xwidth-1,y0:y0+ywidth-1) = box;

figure(1)
imshow(image1)
title('Image 1');

% Define shift
dx= 40; 
dy= 0;


% Shift Image 1
image2 = zeros(M,N);
image2(x0+dx:x0+dx+xwidth-1,y0+dy:y0+dy+ywidth-1) = box;

figure(2)
imshow(image2)
title('Image 2');
% Compute  Fourier Transfoms.

imft1 = fft2(double(image1));
imft2 = fft2(double(image2));


figure(3)
imagesc(abs(fftshift(imft1)));
title('FT of Image 1 (imft1)');

figure(4)
imagesc(abs(fftshift(imft3)));
title('FT of Image 2 (imft2)');


% Define shift in frequency domain
[yF,xF] = meshgrid(-M/2:M/2-1,-N/2:N/2-1);

% Perform the shift
imft3=imft1.*exp(-1i*2*pi.*(xF*dx+yF*dy)/256);

figure(5)
imagesc(abs(fftshift(imft3)));
title('FT of imft3 (imft1 shifted in Fourier Space');


figure(6)
imshow(ifft2(imft3));
title('Inverse FT of imft3');

figure(7)
imagesc(ifft2(imft3)-image2);
title('Residual of IFT of imft3 and image2');


