Gerchberg–Saxton Algorithm

Question

Soong Wei Kang 2020-8-3

0
链接

此问题的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/574315-gerchberg-saxton-algorithm

回答：昕跃张 2021-2-8

Hi, i am new to fourier transform and I get this codes from file exhange, but there is 1 line of code that i don't understand, it is "error = [error; sum(sum(abs(1.32*abs(C) - abs(Target))))];"; can anyone explain it to me? Thank you.

These are the codes:

%Author :Musa AYDIN
%e-mail :maydin@fsm.edu.tr
%Fatih Sultan Mehmet Vakif University,
%Department of Computer Engineering,
%Turkey / Istanbul
%  Gerchberg亡axton Algorithm
%  an example an  of the gs algorithm is shown in this example
%--------------------------------------------------------
% pseudo code of GS algorithm
% Gerchberg亡axton Algorithm(Source, Target, Retrieved_Phase)
%  A = IFT(Target)
%  while error criterion is not satisfied
%    B = Amplitude(Source) * exp(i*Phase(A))
%    C = FT(B)
%    D = Amplitude(Target) * exp(i*Phase(C))
%    A = IFT(D)
%  end while
%  Retrieved_Phase = Phase(A)
%---------------------------------------------------------
% %-------------------------------------------------------
% calculate input beam or input intensity,
% in this example a gaussian beam was selected, 
clear all; close all;
tic;
	x = linspace(-10,10,256);
	y = linspace(-10,10,256);
	[X,Y] = meshgrid(x,y);
	x0 = 0;     		% center
	y0 = 0;     		% center
	sigma = 2; 			% beam waist
	A = 1;      		% peak of the beam 
	res = ((X-x0).^2 + (Y-y0).^2)./(2*sigma^2);
	input_intensity = A  * exp(-res);
	surf(input_intensity);
	shading interp 
%---------------------------------------------------------
	Target=rgb2gray(imread('Sample_Image'));
	Target=double(Target);
	A = fftshift(ifft2(fftshift(Target)));
	error = [];
	iteration_num = 100;
	%hologram = |objectWave + referenceWave|.^2
for i=1:iteration_num
  B = abs(input_intensity) .* exp(1i*angle(A));
  C = fftshift(fft2(fftshift(B)));
  D = abs(Target) .* exp(1i*angle(C));
  A = fftshift(ifft2(fftshift(D)));
  error = [error; sum(sum(abs(1.32*abs(C) - abs(Target))))];   
end
	figure
	subplot(2,1,1);
	imshow(Target);
	title('Original image')
	subplot(2,1,2);
	imagesc(abs(C))               %last pattern
	title('reconstructed image');
	figure
	i = 1:1:i;
	plot(i,(error'));
	title('Error');
	figure
	imagesc(abs(C)) %last pattern
	title('reconstructed image');
	
	toc;