MATLAB Answers

Fourier Transform: IFFT Not Giving Original Image -- All Values Off By Same Constant

10 views (last 30 days)
Corey Hoydic
Corey Hoydic on 4 Feb 2020
Commented: Hiro Yoshino on 5 Feb 2020
I am attempting to write a 2D FFT "round trip" algorithm for a real-valued, symmetric Covariance function (hence, why there is no imaginary component in the transform). I am taking the FFT of an image and immediately taking the IFFT thereafter. Although the shape of the original image is preserved, every value is off by some constant (the constant is 1.0785 in the example case below, but changes for different images). Is there some scaling factor that I am not considering? Note that omitting the value (1/(2*pi))^2 in the forward transform makes the discrepancy even larger. I have included two images and my code below. I have been stuck on this problem longer than I would like to admit -- I dearly appreciate any help. Thank you!
and here is the function:
function [CovIFFT,fomega] = FFT2DDebug(h1,h2,Cov) % h1 is x lags, h2 is y lags (enter -h:h)
%% Calculate omega based on lag and total lags
for i=1:numel(h1)
omega1(i,1)=2*pi*h1(i)/numel(h1);
end
for i=1:numel(h2)
omega2(i,1)=2*pi*h2(i)/numel(h2);
end
%convert covariance to column
covcount=1;
for i=1:numel(h1)
for j=1:numel(h2)
Covcol(covcount,1)=Cov(j,i);
covcount=covcount+1;
end
end
%FFT
fomega=NaN(numel(omega2),numel(omega1));
for ox=1:numel(omega1)
for oy=1:numel(omega2)
sum1=0;
covcount=1;
for hx=1:numel(h1)
for hy=1:numel(h2)
sum1=((1/(2*pi))^2*cos(omega1(ox)*h1(hx)+omega2(oy)*h2(hy)))*Covcol(covcount,1)+sum1;
covcount=covcount+1;
end
end
fomega(oy,ox)=sum1;
end
end
%convert density spectrum to column for smoothing
omegacount=1;
for i=1:numel(omega1)
for j=1:numel(omega2)
fomegacol(omegacount,1)=fomega(j,i);
omegacount=omegacount+1;
end
end
%IFFT
CovIFFT=NaN(numel(omega2),numel(omega1));
for hx=1:numel(h1)
for hy=1:numel(h2)
sum1=0;
omegacount=1;
for ox=1:numel(omega1)
for oy=1:numel(omega2)
sum1=((1/(2*pi))^2*cos(omega1(ox)*h1(hx)+omega2(oy)*h2(hy)))*fomegacol(omegacount,1)+sum1;
omegacount=omegacount+1;
end
end
CovIFFT(hy,hx)=sum1;
end
end
end

  0 Comments

Sign in to comment.

Answers (1)

Hiro Yoshino
Hiro Yoshino on 4 Feb 2020
Why not using fft2?:
I would say it is educational to write it from scratch but no one would do it nowadays to just use it for something else.
If you understand the principle well, you could easily fix it.
From my experience, 99% do not really understand the story behind since 99% of the application are non complex signals.
As you point out, the constant does matter. It works as an adjustor between Time and Frequency domains.

  2 Comments

Corey Hoydic
Corey Hoydic on 4 Feb 2020
I am not using fft2 because the indexing appears to be off -- every real-valued signal I input results in a complex-valued Fourier transform, which should not be the case. Perhaps this is because I am working with a spatial signal that ranges from x,y coordinates of -20 to 20, whereas fft treats those as 1 to 41.
My work aims to employ Bochner's theorem:
which states that a real-valued function in the frequency domain is also real in the spatial domain. The work I am attempting to do requires some processing of the signal in the frequency domain, so I require this signal to be real.
You say that "the constant does matter. It works as an adjustor between Time and Frequency domains." Are you saying that this constant is an error on my part or a necessary result of performing the forward and back transforms?
Thank you for your help.
Hiro Yoshino
Hiro Yoshino on 5 Feb 2020
I see your point.
As for the indexing, you can shift the coordinates by using fftshift.
Let me correct the comments I made, yes the constant does matter but not sure if it is the case with Bochner's theorem.

Sign in to comment.

Sign in to answer this question.


Translated by