9
$\begingroup$

I am basing my question on this post here, because I would like additional details on it, as I have not had any success in re-creating it.

I would like to simply learn, in no uncertain terms, how to do upsampling via the FFT. There seem to be a lot of details regarding if the signal is even or odd, where exactly to zero-pad, how much exactly, etc, and I am not getting those correct numbers for some reason.

Question: So let us say I have an even or odd length signal, of length $N$. I would like to upsample it to, say, $2N$, using FFTs.

How is this done exactly?

What I tried:

I tried following the post and others on the web but there seems to be very little consistency. It basically amounts to "zero-padding the frequency domain and then IFFTing" but I do not get proper results. In my case, I have a signal of length $100$, and I would like to upsample it to length $200$. So I do the following command:

ifft(fftshift([zeros(1,50) fftshift(fft(signal, 100)) zeros(1,50)]));

This however gives me complex data back, and not a signal that is simply upsampled by a factor of $2$. I am not sure where the problem is.

Thank you.

Improve this question
$\endgroup$
2
  • $\begingroup$ In short, you have been misinterpreting how the spectrum is stored. Thanks for the answer, Peter K. $\endgroup$ Commented Mar 10, 2014 at 19:37
  • $\begingroup$ @user7358 Ok, but why the downvote? $\endgroup$ Commented Mar 10, 2014 at 20:27

2 Answers 2

Reset to default
14
$\begingroup$

One trick, for even-length signals, is what to do with the "middle" sample. Here, I've split it half and half between each side of the FFT.

The other trick is to ensure that you have the right amplitudes in the resampled signal. Here's it's a factor of 2.

Try this in scilab:

x = rand(1,100,'normal'); X = fft(x); XX = 2*[X(1:50) X(51)/2 zeros(1,99) X(51)/2 X(52:100)]; xx = ifft(XX); clf; plot([0:199]/200,xx) plot([0:199]/200,xx,'go') plot([0:99]/100,x,'r.')

which appears to generate the right thing. The red dots are the original samples, and the green circles (and blue connecting lines) are the resampled data.

enter image description here

UPDATE

For the case of an odd number of samples, there are fewer fiddly bits (no splitting the last coefficient):

x2 = rand(1,101,'normal'); X2 = fft(x2); XX2 = 2*[X2(1:51) zeros(1,101) X2(52:101)]; xx2 = ifft(XX2); clf; plot([0:201]/202,xx2) plot([0:201]/202,xx2,'gx') plot([0:100]/101,x2,'r.')
Improve this answer
$\endgroup$
2
  • 1
    $\begingroup$ Thanks Peter. May you also please place an addition to your answer for the odd-cases for completion, based on the same example? Much appreciated! $\endgroup$ Commented Mar 10, 2014 at 14:02
  • 1
    $\begingroup$ @TheGrapeBeyond:Done! $\endgroup$ Commented Mar 10, 2014 at 15:21
0
$\begingroup$

To extend PeterK's answer and refer a derivation, having a Nyquist bin = aliasing! This means no general-purpose recovery (or upsampling) exists, but PeterK's approach for real signals at worst loses an imaginary part of one RFFT bin, which is perhaps the best we can get.

This has an overlookable implication: in some applications, we keep all FFTs as powers of 2 for speed, and do filtering/conv/etc with pow2. The processing pipeline may have subsampling and upsampling steps, or, we may wish to test the quality of our lowpass filtering by comparing recovery against original. What we may miss is, if the downsampled signal has a non-zero Nyquist, we've already aliased. This is the only bin with such property.

Python version, with example in the linked answer:

xf, N = fft(x), len(x) xupf = 2 * hstack([xf[:N//2], xf[N//2]/2, zeros(N - 1), xf[N//2]/2, xf[-(N//2 - 1):]]) xup = ifft(xupf)
Improve this answer
$\endgroup$

Start asking to get answers

Find the answer to your question by asking.

Ask question

Explore related questions

See similar questions with these tags.