Coding a function for the Fourier Transform.
显示 更早的评论
I have the following code for a function here,
function F = FourierTransInv(vec,arg)
% This function caluclates the Discrete Fourier Transform or the Inverse
% Fourier Transform.
% The input vec is a vector of length 2^n (n in N with 0) where the entries
% are all real numbers. arg is an input that determines whether or not the
% inverse or regular Fourier Transform is done. An arg of 1 corresponds to
% regular Fourier transform, while an arg of 0 is the inverse.
% The output F is the result of the Fourier or inverse transform.
if ~exist('vec')
error('Please make sure your inputs are valid.');
elseif ~isrow(vec) || ~isnumeric(vec)
error('Please make sure your inputs are valid.');
% Here 2^27 is the maximum because I observed that the difference in time
% between 2^20 and 2^21 is a factor of 2, where 2^20 took 30 seconds. So
% assuming the time increases by a factor of 2 each time, the minimum that
% 2^27 can take is 1 hour which is an extremely long time.
% This can of course be changed easily.
elseif ~ismember(length(vec),2.^(0:27))
error('Please make sure your input is a power of 2.');
elseif ~ismember(arg,[0,1])
arg=1;
warning('Argument is being set to 1 due to invalid input.');
end
N=length(vec);
F=ones(1,N);
switch arg
case 1
w=exp(-2*pi*1i/N);
case 0
w=exp(2*pi*1i/N);
end
switch N
case 1
F(1)=vec(1);
otherwise
Feven=FourierTransInv(vec(2:2:end),arg);
Fodd=FourierTransInv(vec(1:2:end),arg);
for k=1:N/2
F(k)=Fodd(k)+w^(k-1)*Feven(k);
F(k+N/2)=Fodd(k)-w^(k-1)*Feven(k);
end
end
if arg==0
% Division by sqrt(N) because this part always divides twice.
F=F/sqrt(N);
end
end
Now this works for example if I input the vector 
and then run it back through the inversion I get the same vector back. However, I tried some more complicated things such as the one produced with this line of code,
and then run it back through the inversion I get the same vector back. However, I tried some more complicated things such as the one produced with this line of code,n=2^8; L=6*pi; dt=L/n; t=dt*[0:n-1];
y = sin(5*2*pi*t/L)+sin(10*2*pi*t/L);
I've tested many more than this and it all seems to lead to the same conclusion, if I run it through my function, taking the Fourier Transform it agrees with the function fft( . ) everywhere which is great. Then once I run it back through I get back almost the same thing except that every term is multiplied by a factor of 
, why is this? I cannot see anywhere in the code where this would happen, even less so because the simple 
works but not some slightly more complicated vector.
, why is this? I cannot see anywhere in the code where this would happen, even less so because the simple works but not some slightly more complicated vector.采纳的回答
David Goodmanson
2020-11-22
编辑：David Goodmanson
2020-11-22
Hi Anthony,
the only problem here is the normalization. For the inverse transform you are dividing by sqrt(N) every time, (N being the 'local' N). However, the ifft has an overall factor of 1/2^N where 2^N is the length of the initial vector that is being transformed. So you need to divide by 2 at each level of recursion, but not if you have an element of length 2^0 = 1. If you replace
if arg==0
% Division by sqrt(N) because this part always divides twice.
F=F/sqrt(N);
end
with
if arg==0 & N~=1
F=F/2;
end
the results look good.
更多回答（0 个）
类别
在 帮助中心 和 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!
Translated by 
