5
$\begingroup$

I'm trying to grasp the "ah-ha!" moment for the DFT/FFT. One of the points I struggle with is: if the original time signal $x[n]$ is real, then the values of the DFT $X[k]$ are (correct me if I'm wrong) complex-conjugate between the intervals $[-N/2; 0]$ and $[0; N/2]$, or, equivalently, $[0, N/2]$ and $[N/2; N]$. Therefore, I have read that the mirrored part of the spectrum is "not needed" or "redundant", since it provides no new information.

My question is: if we throw away half of the spectrum, and then reconstruct the signal, what do we obtain? The exact same signal, or do we lose a scaling factor (due to summing only half of the spectrum?)

Trying by myself, in Julia, I see there exist the rfft() function, which specifically assumes a real-valued input signal and returns only the first half of the bins (see picture).enter image description here If I then call ifft() to reconstruct the signal, I get these signals:

  • If I used the usual fft(), I get the original signal, except for minor ($\approx 10^{-16}$) imaginary components due to (I think) rounding errors
  • If I used rfft(), I get a signal with imaginary components, and if I take its $Re()$ part, I get a signal of half the length and period (like it is "contracted" along the time axis, effectively having half of the time points), but same "height". For reference, I leave a picture:comparison of original signal and the one reconstructed from the rfft

So my questions are:

  • is, in fact, the mirrored spectrum needed to reconstruct the signal?
  • why, when using only half of it, I needed to take the real part, and not the absolute value?
  • I think the reason why I get a "contracted" signal is simply because rfft() of course outputs a vector of half the length of the original signal. Is this right? And if yes, does this also explains the halved period?
$\endgroup$
1

4 Answers 4

11
$\begingroup$

is, in fact, the mirrored spectrum needed to reconstruct the signal?

No, and yes.

If I give you a spectrum with points from $n = 0$ to $n = \frac N 2$ and I tell you the original signal was real, then you can either

  • regenerate the negative frequency points and use ifft to regenerate the original waveform, or,
  • Use a "half spectrum" FFT algorithm that "knows" what you're doing -- these are a hair more efficient than the first choice, because they've got the fact that the result should be real baked in.

why, when using only half of it, I needed to take the real part, and not the absolute value?

Because you were misusing the rfft/ifft pair. There should be a matching inverse FFT to rfft somewhere in the library. From a comment: "Looking at the top documentation result for 'Julia FFT' there is an irfft function available which is the inverse of rfft"*. If you use that (and it's implemented the way I expect), then the process will Just Work(TM).

I think the reason why I get a "contracted" signal is simply because rfft() of course outputs a vector of half the length of the original signal. Is this right? And if yes, does this also explains the halved period?

Well, yes, keeping in mind that you misused the rfft output by feeding it to ifft. You either need to regenerate the negative half of the spectrum, or you need to find the correct inverse FFT to go with rfft.

* Thank you, @orlp.

$\endgroup$
5
  • 2
    $\begingroup$ Looking at the top documentation result for "Julia FFT" there is an irfft function available which is the inverse of rfft. $\endgroup$
    – orlp
    Commented Oct 3, 2022 at 10:32
  • 1
    $\begingroup$ Thank you. Your technological distinctiveness has been assimilated into my answer. $\endgroup$
    – TimWescott
    Commented Oct 3, 2022 at 15:53
  • $\begingroup$ :O Copied & pasted this time, as I should have done originally. $\endgroup$
    – TimWescott
    Commented Oct 3, 2022 at 18:20
  • $\begingroup$ Thanks! I think I'm getting it now. But I want to mix your answer with @John Bofarull 's one: If I am using the wrong function, why do I re-obtain the right amplitude of the signal? If I am using half of the frequencies, shouldn't I get a signal halved also in amplitude? (since I superpose only half of the composing sinusoids...?) $\endgroup$ Commented Oct 5, 2022 at 10:34
  • $\begingroup$ I don't know -- I'd have to know how Julia scales their FFT. If the fft routine is scaled so that $\sin \omega t$ has an amplitude of 1 after the FFT, then having a 1 in any bin will give you a sine wave with amplitude 1 after the inverse FFT. $\endgroup$
    – TimWescott
    Commented Oct 5, 2022 at 19:57
4
$\begingroup$

Knowing the input spectrum is real is enough information to know that the -ve frequencies are conjugate duplicates of the +ve frequencies, so can be ignored for analysis, and simply recreated from the +ve frequencies for synthesis. A complex input spectrum yields independent +ve and -ve frequencies.

This is consistent with the information content point of view. The FFT works with complex data. If you have N complex data points input, then you have 2N scalars specified. After a transform, you will still have 2N scalars of information, though it will be frequency rather than time information. 2N scalars needs N complex data, and so all outputs carry information. If instead you have N reals as input, all imaginary parts set to be zero, then you only have N scalars of input information, and the same amount after transformation. As the output of the FFT is complex, you only have N/2 independent complex results. By convention, we take the +ve frequency results, but we could just as easily take the -ve set.

$\endgroup$
2
$\begingroup$

I don’t speak Julia, but in Matlab you might do something like:

x = randn(8,1)
Y = fft(x)
x2=ifft([Y(1:5); conj(Y(4:-1:2))])
max(abs(x-x2))

ans =

1.1102e-16

In plain English: Generate N even samples of a real vector. Do the FFT. Use only the first N/2+1 coefficients. Mirror and conjugate to replace larger coeffs. Inverse FFT and confirm that you lost only rounding errors.

Or replace the last bit with:

x3=ifft([Y(1:5); NaN(3,1)], 'symmetric')
max(abs(x-x3))

ans =

 1.1102e-16

Plain English: with the ‘symmetric’ option, the ifft completely ignored the latter ~half of coeffs and not even NaN there will affect the output.

I do not know if Matlab will call the «real iFFT» function of the FFTW library in this last case, or if it will simply mirror and conjugate behind the scenes. A performance test might reveal that.

See the FFTW documentation here: http://www.fftw.org/fftw2_doc/fftw_2.html

$\endgroup$
0
$\begingroup$

If you throw away half the FFT symmetric sectrum you are halving the power of the signal.

Therefore to preserve signal integrity you have to double Spectrum amplitude to compensate for halving the amount of symmetric data around the Nyquist frequency.

$\endgroup$
3
  • $\begingroup$ That doesn't seem to be true, if you look at the second picture I posted: the amplitude of the signal is preserved. At least I think: did I see something wrong? $\endgroup$ Commented Oct 5, 2022 at 10:27
  • $\begingroup$ I said halving the spectrum halves the signal Power, not the amplitude. Px=sum(x.^2). When I learnt this point I also thought it counterintuitive but let show you otherwise: 1. Given x(t) then X(f)=TF(x(t)) 2. X and x are the same signal x in time, X in frequency 3. Power of X(f) = Power of x(t) otherwise X would be the spectrum of a signal different than x 4. now halve X(f) and bring it back to time X1=X([1:f_half]), x1=invTF(X1) 5. for a starter x1 is no longer real but complex because x real X symmetric also applies f to t, P(x1)=1/2*P(x) 6. now bring back x1 to f and cut X1 again in 2 $\endgroup$ Commented Oct 5, 2022 at 19:02
  • $\begingroup$ 7. keep repeating .. and one ends up with a null power signal. Just try it, measure power of X(f) remove one side and measure power of X again? or bring back X1 (X halved) to time and compare P(x) with P(x1). You will see that the power of x1(t) is 3dB below the power of x(t). I would be happy to write a MATLAB script to illustrate this point in detail, so if you post another question to prove P(X)==P(X1) I will answer it. $\endgroup$ Commented Oct 5, 2022 at 19:06

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

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