filters - Upsample data using FFTs. How is this exactly done?

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.

Answer

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.

---

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.')