2
$\begingroup$

Lets say I want to set the minimal sampling rate to reconstruct a 1Hz sine wave, according to the Nyquist-Shannon theorem that states that the maximum recoverable frequency is Fs/2 i.e. we must sample the signal 2 times the maximum frequency.

It seems obvious that the limit is a sampling frequency of 2Hz. Lets even say I interpret it as just more than 2Fs, lets say 3Fs. This would imply that 3 samples is enough.

So, in Matlab i generate:

f=1;%my 1Hz freq
Fs=2*f+1;%=3 =>a bit more than the Nyquist freq
t = linspace(0,1,Fs);
y=sin(2*pi*f*t);
subplot(1,2,1);
plot(t,y);
title('sin(2*pi*1*t)');
subplot(1,2,2);
plot(abs(fft(y)));
title('abs(fft(y))');

The ''sine'' does not even go back up to 0, not even mentioning that even if it did it would be a saw-tooth wave rather than a sine but i guess that is not the problem.

What am I missing, why do I need at least 4 samples rather than 2?

I think this is important to understand the theorem 'in practice'.

enter image description here

Although even with 4 samples the spike in the FT is wrong, it is at 3 rather than 2 (2 is the 1Hz since the first would be the DC freq) enter image description here

enter image description here

In fact I need 5 samples to get finally this spike at 2 in the FFT! Why ?

$\endgroup$
4
$\begingroup$

It's your observation interval which creates the main problem.

Your reasoning based on the Nyquist sampling theorem is ok; of course with a pure sine wave at the exact Nyquist frequency you will have troubles and therefore it's wise to relax the sampling frequency (slightly) above that of Nyquist rate, such 2.2 Hz instead of a strict 2 Hz... So this is one problem you will practically have.

But your main problem that appears on the FFT plot is about the spectral resolution due to short observation interval. Mainly with two samples (or one second observation of a 1 Hz sine wave) your FFT resolution will be limited to even less than one Hz. Please search the site for PSD, spectral resolution, FFT bin frequency to get a better understanding of spectral analysis of practical windowed data.

In order to see sharp frequency peaks (ideally impulses) in the FFT output of your sine wave, you should increase the spectral resolution, which requires you to increase the observation interval.

I have modified and also expanded your code to see the result of ideal sinc based interpolator on the (near) critically sampled data. Note that I've mereley included a digital simulation of ideal sinc based interpolator (not the simulation of analog interpolator) to see that it would indeed reconstruct the pure sinudosidal from its given samples taken at close to Nyquist rate. Note that for the ideal sinc interpolator to work, the original signal should strictly be bandlimited, or at least sufficiently so, which will have a lot of consequences on the success and efficiency of interpolation.

f  = 1;                % 1 Hz. sine wave...
Fs = 4.2*f;            % sampling frequency Fs = 2.2*f  ; a bit more than the Nyquist rate.
Td = 25;               % duration of observation ultimately determines the spectral resolution.
t  = 0:1/Fs:Td;        % observe 25 seconds of this sine wave at Ts = 1/Fs
Td = t(end);           % get the resulting final duration
L  = length(t);        % number of samples in the sequence
M = 2^nextpow2(10*L);  % DFT / FFT length (for smoother spectral display, not better resolution! )

x = sin(2*pi*f*t);     % sinusoidal signal in [0,Td]
%x = x.*hamming(L)';   % hamming window applied for improved spectral display

% Part-II : Approximate a sinc() interpolator :
% ---------------------------------------------
K = 25;                  % expansion factor
xe = zeros(1,K*L);       % expanded signal 
xe(1:K:end) = x;

D = 1024*8;
b = K*fir1(D,1/K);     % ideal lowpass filter for interpolation

y = conv(xe,b);
yi = y(D/2+1:D/2+K*L);

subplot(3,1,1);
plot(t,x);
title(['1 Hz sine wave sampled at Fs = ',num2str(Fs),' Hz, Duration : ', num2str(Td), ' s'])
%xlabel(' time [s]');

subplot(3,1,2);
plot(linspace(-Fs/2,Fs/2-Fs/M,M),fftshift(abs(fft(x,M))));
title(['magnitude of ', num2str(M), '-point DFT / FFT of y[n]']);
%xlabel('Frequency [Hz]');


subplot(3,1,3)
plot(linspace(0,Td,length(yi)),yi);
xlabel('approx simulation of ideal sinc interpolation');

Below is a plot for the result of interpolation from a set of near critical sampling.

enter image description here

And below is the same simulation with a more relaxed sampling, as you can see the interpolator performs much better for this set of improved samples (consequence of better bandlimitedness)

enter image description here

$\endgroup$
  • $\begingroup$ Great answer, thanks! $\endgroup$ – Machupicchu Jul 1 at 5:58
  • $\begingroup$ Maybe just a complementary question about the 'smoothness', why is the 'sawtooth' considered OK? (e.g., for audio, is there some kind if smoothing taking place, e.g. when Fs=44.1KHz or even we could say 8KHz, less, high frequencies will look like these sawtooth (always talking about pure sines), ? $\endgroup$ – Machupicchu Jul 1 at 6:03
  • 1
    $\begingroup$ @Machupicchu there is no sawtooth. It's your mis interpretation of the available sample points, based on an implied simple crude linear interpolation applied to view them on the plot. If you connect the sample points with regard to ideal sinc based interpolation there won't be a sawtooth but a pure sinusoidal appearing. $\endgroup$ – Fat32 Jul 1 at 11:49
  • $\begingroup$ Oh ok thanks, so is it correct to say that I used the "wrong interpolator" to display my wave, and that I should have used a sinc interpolator? $\endgroup$ – Machupicchu Jul 2 at 7:10
  • 1
    $\begingroup$ Thanks for all your answers, great job. $\endgroup$ – Machupicchu Jul 4 at 13:31

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.