close all;
N=128;
%N=256;
f1=200;f0=300;
fc=1000;
fs=8000;

%窗函數(shù)（矩形窗）法設(shè)計(jì)低通濾波器（LPF）
w1=(2*f1)/fs;
w0=(2*f0)/fs;
wL=[0 w1 w0 1];
mL=[1 1 0 0];
h1=fir2(N,wL,mL);
[H1,f]=freqz(h1,1,512,fs);
figure(1);
subplot(2,1,1);
plot(h1);title('LPF單位沖激響應(yīng):h1(n)');
grid on;
subplot(2,1,2);
plot(f,abs(H1));title('LPF頻率特性：H1(f)');
grid on;

%窗函數(shù)（矩形窗）法設(shè)計(jì)帶通濾波器（BPF）
w01=(2*(fc-f0))/fs;
w11=(2*(fc-f1))/fs;
w1r=(2*(fc+f1))/fs;
w0r=(2*(fc+f0))/fs;
fB=[0 w01 w11 w1r w0r 1];
mB=[0 0 1 1 0 0];
h2=fir2(N,fB,mB);
[H2,f]=freqz(h2,1,512,fs);
figure(2);
subplot(2,1,1);
plot(h2);title('BPF單位沖擊響應(yīng):h2(n)');
grid on;
subplot(2,1,2);
plot(f,abs(H2));title('BPF頻率特性：H2(f)');
grid on;

%繪制基帶信號(hào)（0915）
T=ones(1,200);
s=[-T,-T,-T,-T,T,-T,-T,T,-T,-T,-T,T,-T,T,T,-T,];
S_f=abs(fft(s));
figure(3);
subplot(2,1,1);
plot(s);title('基帶信號(hào)：s(n)');
axis([0 3200 -2 2]);
grid on;
subplot(2,1,2);
plot(S_f);title('基帶信號(hào)頻譜：S(f)');
grid on;

%基帶信號(hào)進(jìn)低通濾波
b=conv(s,h1);
B_f=abs(fft(b));
figure(4);
subplot(2,1,1);
plot(b);title('低通濾波后信號(hào)：b(n)=s(n)*h1(n)');
grid on;
subplot(2,1,2);
plot(B_f);title('B(f)');
grid on;

%用濾波后的信號(hào)對(duì)載波進(jìn)行調(diào)制
n=0:3199;
wc=(2*pi*fc)/fs;
y=cos(wc*n);
Y_f=abs(fft(y));
figure(5);
subplot(3,1,1);
plot(n,y);title('載波：y=cos((2*pi*fc/fs)*n)');
subplot(3,1,2);
plot(n,y);
axis([0 100 -1 1]);
grid on;
subplot(3,1,3);
plot(Y_f);title('載波頻譜：Y(f)');
grid on;

b2=b(64:3263);%b2=b(128:3327);
c=b2.*y;
figure(6);
subplot(2,1,1);
plot(c);title('對(duì)b(n)進(jìn)行調(diào)制得到：c(n)=b(n)cos((2*pi*fc/fs)*n)');
grid on;
subplot(2,1,2);
C_f=abs(fft(c));
plot(C_f);title('C(f)');
grid on;

%向已調(diào)信號(hào)中加入白噪聲
Noisy=1*rand(size(c));
d=c+Noisy;
figure(7);
subplot(2,1,1);
plot(d);title('疊加白噪聲：d(n)=c(n)+N(n)');
grid on;
subplot(2,1,2);
D_f=abs(fft(d));
plot(D_f);title('D(f)');
grid on;

%加入白噪聲后的信號(hào)通過(guò)帶通濾波器
e=conv(d,h2);
figure(8);
subplot(2,1,1);
plot(e);title('對(duì)d(n)帶通濾波得到：e(n)=d(n)*h2(n)');
grid on;
subplot(2,1,2);
E_f=abs(fft(e));
plot(E_f);title('E(f)');
axis([0 3300 0 500]);
grid on;

%用相干法對(duì)已調(diào)信號(hào)進(jìn)行解調(diào)（y=cos(wc*n)）
e2=e(64:3263);%e2=e(128:3327);
f=e2.*y;
figure(9);
subplot(2,1,1);
plot(f);title('對(duì)e(n)相干解調(diào)得到：f(n)=e(n)cos(2πfcn/fs)');
grid on;
subplot(2,1,2);
H_f=abs(fft(f));
plot(H_f);title('F(f)');
grid on;

%解調(diào)后信號(hào)通過(guò)低通濾波器
g=conv(f,h1);
g2=g(64:3263);%g2=g(128:3327);
figure(10);
subplot(2,1,1);
plot(g2);title('對(duì)f(n)低通濾波得到：g(n)=f(n)*h1(n)');
grid on;
subplot(2,1,2);
G_f=abs(fft(g2));
plot(G_f);title('G(f)');
grid on;

%判決恢復(fù)后的基帶信號(hào)波形
for p=1:16
    sum=0;
    for q=1:200
        sum=sum+g2((p-1)*200+q);
    end
    if sum>0
        R(p)=1;
           else
       R(p)=-1;
   end
end
for i=1:16
      if R(i)==1
         for j=1:200
              K((i-1)*200+j)=1;
         end
     else
         for j=1:200
              K((i-1)*200+j)=-1;
         end
     end
end
figure(11);
subplot(3,1,1);
plot(R);
axis([0 16 -2 2]);
grid on;
subplot(3,1,2);
plot(K);title('判決后得到的基帶信號(hào)：s1(n)');
axis([0 3200 -2 2]);   
grid on;
subplot(3,1,3);
plot(s);title('基帶信號(hào)：s(n)');
axis([0 3200 -2 2]);   
grid on;