1
$\begingroup$

Following this question: Shannon-Nyquist theorem reconstruct 1Hz sine wave from 2 samples

could you explain the algorithm to apply for sinc interpolation to avoid the "sawtooth" effect due to linear interpolation?

(It seems to me that the https://en.wikipedia.org/wiki/Whittaker%E2%80%93Shannon_interpolation_formula Shannon-Whittaker formula would be the one suitable for this?)

Improve this question
$\endgroup$
4
  • $\begingroup$ Have you read this question and answers? dsp.stackexchange.com/questions/59068/… $\endgroup$
    – Cedron Dawg
    Commented Jul 4, 2019 at 15:08
  • $\begingroup$ Well you must admit that this one is a bit discouraging ... $\endgroup$
    – SheppLogan
    Commented Jul 4, 2019 at 15:12
  • $\begingroup$ Discouraging? Your question or the link I gave? Here's another one that may give you some insight: dsp.stackexchange.com/questions/58032/… $\endgroup$
    – Cedron Dawg
    Commented Jul 4, 2019 at 15:20
  • $\begingroup$ I have expanded the answer to add a sinc interpolator simulation to your previous question. $\endgroup$
    – Fat32
    Commented Jul 5, 2019 at 7:18

2 Answers 2

Reset to default
2
$\begingroup$

Here again let me note that exact Nyquist frequency for a pure sine wave should be avoided. Shannon-Nyquist sampling theorem requires that there's no impulse at the Nyquist frequency as a consequence of bandlimitedness, the content at the exact Nyquist frequency is taken to be zero.

Then the following code demonstrates the approximate simulation of an ideal sinc based interpolator applide to (near) critical samples of a pure sine wave. Note that, any finite observation of a signal cannot be bandlimited, so this simulation is not a perfect representation of the true output from an ideal interpolator, nevertheless by chosing signal duration long enough, one can attain an approximately bandlimited signal.

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

With the result of

enter image description here

Improve this answer
$\endgroup$
1
  • 1
    $\begingroup$ Fantastic answer, thanks $\endgroup$
    – SheppLogan
    Commented Jul 5, 2019 at 12:03
1
$\begingroup$

A pure sine wave sampled at exactly the Nyquist frequency is simply an alternating sequence of $C(-1)^n$, where C can be anything from zero to the amplitude of the sine wave depending on where in the cycle the points are sampled. Thus the only plausible reconstruction (which can be done as an infinite series of sincs, see Convergence of periodic sinc interpolation) is a sine wave with the amplitude of C.

The Nyquist is the limit. Anything close is going to give you numerical problems when implemented.

Ced

Your question is actually the same as what I am working on in regards to the resampling of a sampled cycle having an even number of points.

Consider this:

There is a merry-go-round in a playground. There is one child on the merry-go-round and a camera is taking pictures at even intervals. The collection of pictures is examined. It shows the child on opposite sides of the merry-go-round frame by frame. So, does this mean the merry-go-round is going clockwise at one revolution per two pics taken, or counter-clockwise at the same frequency, or is the merry-go-round standing still and the child is running side to side between the pics?

Now, let's introduce aliases.

Couldn't it also be the merry-go-round is going around one and a half times for every picture taken, instead of a half? Or two and a half? In the opposite direction?

Likewise, maybe the child is running back and forth one and a half, or two and a half ....

Who says the assumption of smooth motion is even required? Couldn't the merry-go-round be doing all sort of crazy spinning that just happens to show the child at those positions at those times?

How much faster does the camera have to be to determine any of that? Do you need a second camera?

Improve this answer
$\endgroup$
4
  • $\begingroup$ Interesting analogy, thanks.So the aliases would occur if we took only 2 pictures per second and the merry-go-round would actually be rotating at strictly more than 1 rotation per second, according to the Nyquist-Shannon thm. But in fact there we might not even notice these aliases, i guess if it is say 4 cycles/second, or any integer multiple of 1 ie. any natural number of cycles / second, .. right ? $\endgroup$
    – SheppLogan
    Commented Jul 5, 2019 at 14:51
  • $\begingroup$ Where do "seconds" come from? At Nyquist, you have two samples per cycle, or 1/2 cycle per sample, or $\pi$ radians per sample. If you double the frequency, then you have 1 cycle per sample and the child looks like it isn't moving. If you triple the frequency, you have 1 1/2 cycles per sample which is indistiguishable (an alias of) 1/2 cycle per sample. So, the odd multipliers, not the even ones. Can also be seen as additions of wholes. Nyquist is at "the band limit", so the aliases are not "bandlimited", which is why "bandlimited" is often specified to rule out aliases. $\endgroup$
    – Cedron Dawg
    Commented Jul 5, 2019 at 15:20
  • $\begingroup$ Getting to the heart of my analogy, suppose my samples are in the form of $(-1)^n$, that fits $\cos(n\pi)$, so that is a real valued bandlimited solution. $$ n (samples) \cdot \pi (radians per sample) = n \pi (radians) $$ Suppose though, that my signal is allowed to be complex valued, then $e^{i n\pi}$ and $e^{-i n\pi}$ are also solutions. As is any linear combination of them, including $\cos(n\pi)$. $\endgroup$
    – Cedron Dawg
    Commented Jul 5, 2019 at 15:21
  • $\begingroup$ I can also add any multiple (real or complex) of $\sin(n\pi)$ to your underlying signal and not change the sampled values one bit. The "multiple" can even be a well behaved function, but then I would not longer be "bandlimited". $\endgroup$
    – Cedron Dawg
    Commented Jul 5, 2019 at 15:26

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.