3
$\begingroup$

There are already questions on this topic here and here. I am trying to calculate the ifft using a DSP that only has a library function for a forward REAL (assumes conjugate symmetry) fft. For a given a signal of length N, the output of the FFT process is length N interleaved real/complex values, where the real part of the nyquist bin is put into the imaginary part of the DC component (always zero for a real signal - so this keeps a constant memory footprint). I'm having trouble fathoming how to calculate the inverse using such a function.

$\endgroup$
  • 1
    $\begingroup$ FFT and IFFT are essentially the same operation. If you can only do real input → complex output, then maybe you can process real and imaginary parts separately and combine the results. $\endgroup$ – endolith Aug 13 '18 at 18:44
3
$\begingroup$

EDIT: This answer is basically a more detailed version of Hilmar's answer, which I noticed too late. Below there are more details about the derivation and the rearrangement of the outputs of the real FFT.


I assume you use the definitions of the FFT and IFFT as commonly used in DSP:

$$\textrm{FFT: }\quad X[k]=\sum_{n=0}^{N-1}x[n]e^{-j2\pi nk/N}\tag{1}$$

$$\textrm{IFFT: }\quad x[n]=\frac{1}{N}\sum_{k=0}^{N-1}X[k]e^{j2\pi nk/N}\tag{2}$$

For a real FFT we assume that $x[n]$ is real-valued. We also assume that the FFT length $N$ is even (otherwise you don't get an FFT bin at Nyquist).

From $(1)$ and $(2)$ we see that

$$\text{IFFT}\big\{X[k]\big\}=\frac{1}{N}\left(\textrm{FFT}\big\{X^*[k]\big\}\right)^*\tag{3}$$

Since your FFT routine only accepts real-valued inputs, we need to write $X[k]$ in terms of its real and imaginary parts:

$$X[k]=X_R[k]+jX_I[k]\tag{4}$$

Furthermore, if we assume that the inverse FFT of $X[k]$ is real-valued, we can rewrite $(3)$ as

$$\begin{align}\text{IFFT}\big\{X[k]\big\}&=\frac{1}{N}\textrm{Re}\Big(\textrm{FFT}\big\{X_R[k]\big\}-j\,\textrm{FFT}\big\{X_I[k]\big\}\Big)\\&=\frac{1}{N}\Big[\textrm{Re}\big(\textrm{FFT}\big\{X_R[k]\big\}\big)+\textrm{Im}\big(\textrm{FFT}\big\{X_I[k]\big\}\big)\Big]\tag{5}\end{align}$$

A little extra inconvenience lies in the fact that a real FFT routine only computes $N/2+1$ values whereas Eq. $(5)$ assumes FFT results of length $N$. This Matlab/Octave script shows how to correctly arrange the results of the individual FFTs:

N = 64;     % must be even!
N2 = N/2;

x = randn(N,1);

X = realfft(x);
XR = real(X); XR = [XR;XR(N2:-1:2)];
XI = imag(X); XI = [XI;-XI(N2:-1:2)];

x2a = real(realfft(XR)) / N; x2b = imag(realfft(XI)) / N;
x2 = [x2a + x2b; x2a(N2:-1:2) - x2b(N2:-1:2)];

In the above script we assume that there's a function realfft.m defined by

function X = realfft(x)

N = length(x);
if rem(N,2), error('N must be even.'); end

X = fft(x,N);
X = X(1:N/2+1);

end
$\endgroup$
  • $\begingroup$ Great! That was fun in C! The realfft in the example really helped debugging. $\endgroup$ – learnvst Aug 13 '18 at 23:02
2
$\begingroup$

You can basically just to an individual FFT of the real and imaginary parts, post process a bit and sum the results. In essence you calculate the even and odd part of the time signal separately.

The Matlab code illustrates the procedure but doesn't include the packing or unpacking. There are also ways to optimize this by using the symmetry properties which are not shown.

% create a random sequnec
nx = 1024;
x0 = randn(nx,1);
% forward FFT
fx = fft(x0);
% calculate inverse by doing seperate FFTs over real and imaginary part
x1 = real(fft(real(fx)))/nx; % this will be real result: even part
x2 = imag(fft(imag(fx)))/nx; % this will be entirely imaginary,: odd part
% just sum them up
x3 = x2+x1;
% calculate the error
err = x3-x0;
fprintf('Error = %6.2f dB\n',10*log10(sum(err.^2)./sum(x0.^2)));
$\endgroup$
0
$\begingroup$

In Short if $X[k]$ represents the $N$-point FFT of an $N$-point signal $x(n)$ , then $N * x(n)$ (sequence $x(n)$ scaled by $N$) can be obtained by performing FFT on the sequence $X[0], X[N-1], X[N-2], \cdots ,X[1]$.

Explanation:

Consider the equation of IFFT

$$IFFT: x[n] = \frac{1}{N} ( \sum_{k=0}^{N−1}X[k]*e^{j2πnk/N} )$$

It can be split and written in the following way

$$ x[n] = \frac{1}{N} ( X[0] + \sum_{k=1}^{N−1}X[k]*e^{j2πnk/N} )$$

$$ x[n] = \frac{1}{N} ( X[0] + \sum_{k=1}^{N−1}X[N - k]*e^{j2πn(N-k)/N} )$$

$$ x[n] = \frac{1}{N} ( X[0] + \sum_{k=1}^{N−1}X[N - k]*e^{j2πn(N)/N}*e^{-j2πn(k)/N})$$

The term $e^{j2πn(N)/N}$ evaluates to one and hence the above can be re written as

$$ x[n] = \frac{1}{N} ( X[0] + \sum_{k=1}^{N−1}X[N - k]*e^{-j2πn(k)/N})$$

$$ x[n] = \frac{1}{N} ( X[0] * e^{-j2π(0)(k)/N} + \sum_{k=1}^{N−1}X[N - k]*e^{-j2πn(k)/N})$$

Consider a series $Y[k]$ defined as follows

$$Y[0] = X[0]$$ and

$$Y[k] = X[N-k]$$ $for$ $k = 1,2,\cdots,N-1$

Then the last of the above equation can be re written as

$$ x[n] = \frac{1}{N} ( Y[0] * e^{-j2π(0)(k)/N} + \sum_{k=1}^{N−1}Y[k]*e^{-j2πn(k)/N})$$

$$ x[n] = \frac{1}{N} ( \sum_{k=0}^{N−1}Y[k]*e^{-j2πn(k)/N})$$

$$ x[n] = \frac{1}{N} ( FFT(Y[k]))$$

Where $Y[k]$ is the series obtained from rearrangement of $X[k]$

$X[0], X[N-1], X[N-2], \cdots, X[1]$

Note:

$$FFT(x_{re}(n) + j * x_{im}(n)) = FFT(x_{re}(n)) + j * FFT(x_im(n))$$

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.