From: kk KKsingh on
I just revised some thing in my code so posted here , It is based on following paper
http://docs.google.com/leaf?id=0B9lyGDKrglBfOGEyYzNiZjAtNWE2Zi00OWU1LThjN2YtNjgyZDFhOWFiOGYz&sort=name&layout=list&num=50

Just see eqn 4,5 and 12 .......Am i doing some thing wrong

%For decimation please type as .20 for 20%;

%This is to test code


clear all;
close all;
N=1024; %Number of Sample points
Fs=1000;%Sampling Frequency
dt=1/Fs;%Sampling interval
fo=50;% Frequency of sin wave
t=0:dt:dt*(N-1);%Sampling Period
omega=2*pi*fo;
S=sin(omega*t);




%Applying FFT
f=fftshift(fft(S));

%Applying the decimation kernal

pd=input('Enter the percent decimation')
if pd == 0
W = ones(size(S));
elseif pd ==1
W = zeros(size(S));
else
W = rand(size(S));
W = W + (1-mean(W) - pd);
mean(W)
W = round(W);
end

S1=W.*S;
%finding zeros from the signal
index=find(S1==0)

%Removing zeros from the signal
S2=S1;
t2=t;
S2(index)=[];
t2(index)=[];


N1=numel(S2);

freq_NDFT=Fs/N1*(-N1/2:N1/2-1);
S2=S2(:);
t2=t2(:);
freq_NDFT=freq_NDFT(:);


%Forward Fourier Transform
deltak=1/max(t2)-min(t2);
DFT=exp(-2*pi*i*deltak*freq_NDFT*t2');
X_DFT=DFT*S2;

%Inverse Discrete Fourier Transform
IDFT=DFT'*(X_DFT.*freq_NDFT);






%Now Deciding weights as peer dujhdham


%------------------------------------------------------
%Schonvillie NDFT is solved using Reinmann sum weights are
%used here
%NDFT Accoring to Schonvillie ยจ
N3=numel(S2);
Period=max(t2)-min(t2);
deltax_n=zeros(N1,1);
deltax_n(1,:)=((t2(2)-(t2(N1)-Period)))/2;
deltax_n(N1,:)=((t2(1)+Period)-(t2(N3-1)))/2;
for k=2:N1-1;
deltax_n(k,:)= ((t2(k+1)-t2(k-1))/2);
end

%Now applying weights to the DFT
X_DFTw=DFT*(S2.*deltax_n)

%Bandlimiting the signal
X_DFTw(1:200,:) = 0;
X_DFTw((end-199):end,:) = 0;

%Making inverse kernal and Applying to wighted DFT
X_DFT_least= inv(DFT*diag(deltax_n)*DFT')*X_DFTw;




figure(1)
subplot(311)
plot(freq_NDFT,abs(X_DFT))
legend('NDFT')
subplot(312)
plot(freq_NDFT,abs(X_DFTw))
legend('NDFT with Weights')
subplot(313)
plot(freq_NDFT,abs(X_DFT_least))
legend('Least Squares NDFT')

IDFT = exp(2*pi*1i*t'*freq_NDFT');
sigout = IDFT*X_DFT_least
figure(2)
plot(t,real(sigout),t,real(S1),t,real(S));

I will really appreciate your help !
From: Jan Simon on
Dear KKsingh,

> Just see eqn 4,5 and 12 .......Am i doing some thing wrong

I don't know and I cannot know. Do *you* have the impression, that something is wrong? Errors? Differences between expectations and results?

Kind regards, Jan