% dm11201
% 信號采樣與重構
% k = 1  臨界采樣  0<k<1  過采樣  k>1  欠采樣
display('Please input the value of k');
k = input('k = ');
wm = 1;                                   %信號帶寬
Ts = k*pi/wm;                             %采樣間隔
ws = 2*pi/Ts;                             %采樣角頻率
wc = ws/2;                                %為簡便，濾波器截止頻率wc=ws/2
n = -10:10;                               %定義采樣點的數量
m = fix(length(n)/2);                     %單邊采樣周期數,fix函數的功能為取整
nTs = n*Ts;                               %計算每個采樣點
dt = 0.05;
t = -m*Ts:dt:m*Ts;
%Sa(t)信號采樣
f = sinc(nTs/pi);                        %計算Sa(t)在采樣點處的函數值，得到采樣信號
%Sa(t)信號重構
fa = f*Ts*wc/pi*sinc((wc/pi)*(ones(length(nTs),1)*t-nTs'*ones(1,length(t)))); 
%Sa(t)包絡線
f2=sinc(t/pi);
subplot(311);
stem(nTs,f);                              %繪制采樣信號
hold on
plot(t,f2,'r:');                          %繪制包絡線
xlabel('t');
title('對sa(t)信號進行采樣');
axis([-m*Ts m*Ts -0.5 1.2]);
hold off
subplot(312);
h1=plot(t,fa);                            %繪制重構信號
hold on
for i=-m*Ts:Ts:m*Ts
    ft = sinc(i/pi)*sinc(wc/pi*(t-i));
    h2=plot(t,ft,'m:');                   %繪制重構信號的各個分量
    hold on
    yy = sinc(i/pi);
    yy1 = yy:yy;
    plot(i,yy1,'o');
end
xlabel('t');
title('重構信號');
axis([-m*Ts m*Ts -0.5 1.2]);
legend([h1,h2],'重構信號','分信號');
hold off
subplot(313);
%計算原信號和重構信號之間的誤差
error=abs(fa-f2);
plot(t,error);
axis([-m*Ts m*Ts min(error) max(error)+0.1*max(error)]);
xlabel('t');
ylabel('error(t)');  
%End