Successive FFT2 transforms cause an increase in amplitude even though each individual transform is normalised.
1 次查看(过去 30 天)
显示 更早的评论
I will be very grateful if anyone can help here. I have some code to simulate Fourier Optics and need to achieve this with a series of Fourier transforms over a short distance where I will eventually take account of turbulence at each step.
Unfortunately I have noticed that when I perform repeated FFT2 transforms on my gaussian beam profile, the resulting beam's amplitude if greatly increased even tho I am normalising the transformed function each step of the way.
My code is as follows:
%%Fourier Optics Code: ===================================================
clear; % clear workspace
steps = 2^8; % Resolution (256)
lim = 1; % Field of View (m)
scale = 2/steps; % Normalisation factor for fft & ifft
ii = sqrt(-1); % Neat Imaginary unit
% Set Up Co-ordinate Grids: -----------------------------------------------
% Direct Space Window 2D:
x = linspace(-lim, lim, steps); y = x; [X,Y]=meshgrid(x,y);
Xs = mean(diff(x)); % Sampling Interval (delta x)
Ks = 1/Xs; % Sampling Frequency (per distance)
Kn = Ks/2; % Nyquist Frequency
% K-Space Window 2D:
kx = linspace(-1, 1, steps)*Kn; ky = kx; [KX,KY] = meshgrid(kx,ky);
% Beam Properties: --------------------------------------------------------
lambda = 900E-06; % Wavelenght (meters)
w0 = 0.01; % Beam Wasit (meters)
k0 = 2*pi/lambda; % Wavenumber (meters^-1)
z0 = (k0*w0^2)/2; % Rayleigh Range
% Propogation Distance: ---------------------------------------------------
dz = 5; % meters
zsteps = 100; % number of propagation steps taken
% Neat Veriable for Cylindrical Polars: -----------------------------------
R = sqrt(X.^2+Y.^2); % Radial Magnitude
K = sqrt(KX.^2+KY.^2); % "Radial" K-vector
dGouy = -atan(dz/z0); % Gouy Phase progression due to dz propagation
% Making Beam: ------------------------------------------------------------
func = exp((-1)*(R./w0).^2); % Gaissian Beam (Waist = w0 in plane z = 0)
% Propogate Beam ---------------------------------------------------------
next = func; peaks = zeros(length(zsteps),1);
for k = 1:zsteps
if k == 1; % Ensures Beam in centre of Field of View
Func = fftshift(fft2(next)*scale^2);
else
Func = fft2(next)*scale^2;
end
% Fourier Optics:
H = exp(ii*k0*dz) .* exp(-ii*(dz.*K.^2)./(2*k0)); % Propagator
HF = H .* Func;
next = (ifft2(HF)/scale^2)*dGouy;
peaks(k) = (abs(max(max(next)))); % Record the Max of the output
end
% Plotting =================================================================
fig = figure(1);
subplot(1,3,1);
surf(X,Y,abs(func),'EdgeColor','none');
title('Original Beam');
xlabel('Direct Space FOV (m)');
ylabel('Direct Space FOV (m)');
axis('tight');
view([0,0,1]);
subplot(1,3,2);
surf(X,Y,abs(next),'EdgeColor','none');
title('Propagated Beam');
xlabel('Direct Space FOV (m)');
ylabel('Direct Space FOV (m)');
axis('tight');
subplot(1,3,3);
plot(diff(peaks))
title('Increase in Max Amplitude');
xlabel('Number of Steps in Z');
ylabel('Change in Amplitude');
axis('tight');
which outputs something like this:
As you can see the amplitude of the beam is increasing exponentially with propagation which isn't right! The amplitude should decrease with distance, and the increase shouldn't be from the scaling of the FFT2 and IFFT2 functions as I have tested my scaling method extensively in other programs. If any one is able to shed any light on this issue or notice something I may be doing wrong, I would be extraordinarily greatful!
0 个评论
回答(0 个)
另请参阅
类别
在 Help Center 和 File Exchange 中查找有关 Fourier Analysis and Filtering 的更多信息
产品
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!