hi,
I want to implement an spectral derivative using fft() and ifft().
Theoretically it can be done scaling the fft coefficients by the
factor (2 pi i w)^n, where n is the derivate order, and applying the
ifft() to them.
It is, the derivative of a function f_j(x_j) can be obtained by
derivation of it's spectral representation, so if we have
f_j = sum_{w=-M} ^{M} e^{2 pi i w j} ftilde(w)
with j = 0, ..., 2M, then the f'_j can be obtained by
f'_j = ( 2 pi i w )*sum_{w=-M} ^{M} e^{2 pi i w j} ftilde(w)
in terms of the fft() and ifft() functions it should be:
f'_j = ifft( ( 2 pi i w ) * fft( f_j ) )
with w = -M, ..., M
but it doesn't gives the correct result. I suspect that the problem is
in the order of the returned fft() coefficients.
I found that the coef. ftilde(w) are stored in this way:
[ftilde(0) ftilde(1) ... ftilde(M-1/2) ftilde(M+1/2) ftilde(-M+1/2)
... ftilde(-1)] = fft(f_j)
or, defining N = 2M+1 it holds
[ftilde(0) ftilde(1) ... ftilde(N/2-1) ftilde(N/2) ftilde(-N/2+1)
... ftilde(-1)] = fft(f_j)
is that ok?
Here is my test code (it works ok with this simple function y =
cos(2*pi*x), but not for example with y = x + cos(2*pi*x) ) :
N = 2^7;
h = 1/N;
x = [0:N-1]*h;
y = cos(2*pi*x);
figure; plot (real(y));
Y = fft( y , N );
M = (N-1)/2;
W = (2*pi*i)*[0:M-1/2 0 -M+1/2:-1];
Y2 = W.*Y;
y2 = ifft( Y2 );
figure; plot (real(y2));
pause;
i hope someone can help, thanks in advanced
Gastón