clear all;
close all;
A=1;
fm=1;           
fc=2;           %2Hz;
N_sample=8;
N=100;           %碼元數
Ts=1;             %1Baud/s

dt=Ts/fc/N_sample;%波形采樣間隔
T=Ts-dt;
t=0:dt:N*T;
Lt=length(t);

%產生二進制信源
d=sign(randn(1,N));
dd=sigexpand((d+1)/2,fc*N_sample);
gt=ones(1,fc*N_sample);%NRZ波形

figure(1)
subplot(221);%輸入NRZ信號波形（單極性）
d_NRZ=conv(dd,gt);
plot(t,d_NRZ(1:length(t)));
axis([0 15 0 1.2]);ylabel('輸入信號');   

subplot(222);%輸入NRZ頻譜
[f,d_NRZf]=T2F(t,d_NRZ(1:length(t)));
plot(f,10*log10(abs(d_NRZf).^2/T));  
axis([-2 2 -50 20]);ylabel('輸入信號功率譜密度（dB/Hz）');

%2FSK
%s_2fsk=Acos(2*pi*fc*t+int(2*d_NRZ-1));
sd_2fsk=2*d_NRZ-1;

s_2fsk=A*cos(2*pi*fc*t+2*pi*sd_2fsk(1:length(t)).*t);
subplot(223)
plot(t,s_2fsk);
axis([0 15 -1.2 1.2]);xlabel('t');ylabel('FSK');
subplot(224)
[f,s_2fskf]=T2F(t,s_2fsk);
plot(f,10*log10(abs(s_2fskf).^2/T));
axis([-fc-4 fc+4 -50 20]);xlabel('f');ylabel('2FSK功率譜密度（dB/Hz）');

%N0=0.01;
B=2*fm;
%noise=noise_nb(fc,B,N0,t);
%s_2fsk=s_2fsk+noise
figure(2)
subplot(311)
plot(t,s_2fsk);
title('fsk調制信號');
xlabel('t');

rt=s_2fsk.*(2*pi*fc*t+2*pi*sd_2fsk(1:length(t)).*t);
rt=rt-mean(rt);
[f,rf]=T2F(t,rt);
[t,rt]=lpf(f,rf,2*fm);

subplot(312)
plot(t,rt);
title('t');
subplot(313)
[f,sf]=T2F(t,s_2fsk);
psf=(abs(sf).^2)/T;
plot(f,psf);
axis([-2*fc 2*fc 0 max(psf)]);
title('fsk功率譜');
xlabel('f');