Shannon vs dirichlet kernel interpolation method for signal reconstruction

Shannon vs dirichlet kernel interpolation method for signal reconstruction

I am currently studying fourier transform, and especially the way that the
signal could be reconstructed from its spectrum.
In many lectures, I have seen the shannon interpolation method to
reconstruct the signal: $ f(t) = \sum_{n=-N}^{N} f(nT_{e})
sinc\left(\frac{t-nT_{e}}{T_{e}}\right) $
with $ T_{e} $ being the sampling period
I did not tried to do the demonstration, I read that it could be achieved
using poisson summation formula but I did not investigated much more on
it.
But in other lectures, I have found a demonstration, that uses the
dirichlet kernel: $ f(t) = \frac{1}{2N+1} \sum_{n=0}^{2N} f[nT_{e}]
D\left( 2\pi F_{e} \left( t - \frac{nT_{e}}{2N+1} \right) \right) $
Where D which is defined as follow:
$ D(t) = \frac{sin((N+\frac{1}{2})t)}{sin(\frac{t}{2})} $
I tried to re-do the demonstration, and I think I achieved to do it using
fourier coefficient expression
So I decided to try these two versions on matlab to assess the accuracy of
each of them. But the dirichlet version give me strange results
(discontinuity).
I would like to know if comparing these two methods is the right things to
do,or if there is a fundamental difference between the that I have not
seen.
Here is my matlab code, I would be glad to have a feedback on it :
clc;
clear all;
close all;
Fe = 3000;
Te = 1/Fe;
Nech = 100;
F1 = 500;
F2 = 1000;
FMax = 1500;
time = [0:Te:(Nech-1)*Te];
timeDiscrete = [1:1:Nech];
frequency = (timeDiscrete/Nech)*Fe;
signal = cos(2*pi*F1*(time))+cos(2*pi*F2*(time))+cos(2*pi*FMax*(time));
%Compute the FFT
spectrum=zeros(Nech);
for k = timeDiscrete
for l = timeDiscrete
spectrum(k) = spectrum(k) + signal(l)*exp(-2*pi*j*l*k/Nech);
end
end
%Compute de inverse FFT
reconstruction=zeros(Nech);
for k = timeDiscrete
for l = timeDiscrete
reconstruction(k) = reconstruction(k) +
spectrum(l)*exp(2*pi*j*l*k/Nech);
end
end
reconstruction=reconstruction/Nech;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%% Now interpolation will take place %%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
Finterp = 18000;
Tinterp = 1/Finterp;
TimeInterp = [0:Tinterp:(Nech-1)*Te];
[m,n] = size(TimeInterp);
NechInterp = n;
TimeInterpDiscrete = [1:1:NechInterp];
%Compute original signal value without any interpolation
signalResampled =
cos(2*pi*F1*(TimeInterp))+cos(2*pi*F2*(TimeInterp))+cos(2*pi*FMax*(TimeInterp));
%Compute original signal interpolation through shannon interpolation method
interp_shannon=zeros(NechInterp);
for k = TimeInterpDiscrete
for l = timeDiscrete
interp_shannon(k) = interp_shannon(k) +
reconstruction(l)*sinc(Fe*(TimeInterp(k)-time(l)));
end
end
%Compute original signal interpolation through dirichlet kernel
interpolation method
interp_dirichlet=zeros(NechInterp);
for k = TimeInterpDiscrete
for l = timeDiscrete
x = 2*pi*Fe*(TimeInterp(k)-time(l))/NechInterp;
interp_dirichlet(k) = interp_dirichlet(k) +
reconstruction(l)*(sin((NechInterp/2 + 0.5)*x)/sin(0.5*x));
end
end
interp_dirichlet = interp_dirichlet/NechInterp;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%% Let's print out the result %%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% figure(1);
% plot(time,signal);
% hold on;
% plot(time,real(reconstruction),'r');
%
% figure(2);
% plot(frequency(1:Nech/2),abs(spectrum(1:Nech/2)));
figure(3);
% Ground truth : deterministic siganl is recomputed
plot(TimeInterp,signalResampled,'g');
hold on;
% linear interpolation between subsampled points (matlab tracing tool)
plot(time,real(reconstruction),'r');
hold on;
% Shannon interpolation method
plot(TimeInterp,real(interp_shannon),'b');
hold on;
% Dirichlet kernel interpolation method
plot(TimeInterp,real(interp_dirichlet),'k');
Thank you for your help