FM调制解调系统

李秋锐

%FM调制解调系统.m
%频率调制与解调的Matlab演示源程序
%可以任意改原调制信号函数m(t)
4 U& h" ?: b9 ]8 L) Z$ {%电子信息工程  杜怀超 091308305  
# {! j4 N8 H, q1 E. s" u4 b%·*·*·*·*·*·*·*·*·*·*·*·*·*·*·*·
2 K2 B1 o2 g9 {6 X! j4 H, K%*****************初始化******************
echo   off      
" a4 x, B0 B! `' W3 F( jclose  all
clear  all
( o" L/ i% ?; `+ M( M, [clc
%*****************************************
%·*·*·*·*·*·*·*·*·*·*·*·*·*·*·*·
%****************FM调制*******************
dt=0.001;                                   %设定时间步长
2 y; A2 Q1 L) d! o1 E8 Lt=0:dt:1.5;                                       %产生时间向量
am=5;                                      %设定调制信号幅度
- `4 G4 M8 p0 j" Nfm=5;                                      %设定调制信号频率
mt=am*cos(2*pi*fm*t);                        %生成调制信号
5 C6 N' e- M* Jfc=50;                                      %设定载波频率
ct=cos(2*pi*fc*t);                             %生成载波
kf=10;                                      %设定调频指数
int_mt(1)=0;
$ w5 Y( |: M+ a5 o2 Y, ?& @: nfor i=1:length(t)-1  
    int_mt(i+1)=int_mt(i)+mt(i)*dt;               %求信号m(t)的积分
% I# B1 K4 a+ E. mend                                          %调制，产生已调信号
) s; Q6 l- o; Y8 D5 s$ o4 ~* Wsfm=am*cos(2*pi*fc*t+2*pi*kf*int_mt);            %调制信号
%*****************************************
%·*·*·*·*·*·*·*·*·*·*·*·*·*·*·*·
%*************添加高斯白噪声**************
( r& m* b- ]& V' dsn1=10;                                     %设定信躁比(小信噪比)
sn2=30;                                     %设定信躁比(大信噪比)
  Q6 s9 ]! e) O, Z4 Esn=0;                                       %设定信躁比(无信噪比)
db=am^2/(2*(10^(sn/10)));                     %计算对应的高斯白躁声的方差
n=sqrt(db)*randn(size(t));                      %生成高斯白躁声
1 d3 P) q+ U2 {. }5 ~$ @nsfm=n+sfm;                                %生成含高斯白躁声的已调信号（信号通
  U; N( z. u8 j9 ~0 p%过信道传输）
%*****************************************
%·*·*·*·*·*·*·*·*·*·*·*·*·*·*·*·
%****************FM解调*******************
for i=1:length(t)-1                             %接受信号通过微分器处理
% q: X7 X4 i3 ], ]    diff_nsfm(i)=(nsfm(i+1)-nsfm(i))./dt;
/ h+ F- C* k! u4 C4 o. Cend
# @7 K0 Z0 B: {4 {9 }# Z! }$ t- c, Udiff_nsfmn = abs(hilbert(diff_nsfm));     %hilbert变换，求绝对值得到瞬时幅度（包络检波）
' p2 F1 O& q' O- q2 b" k  azero=(max(diff_nsfmn)-min(diff_nsfmn))/2;
, i! G1 M- H5 M$ P$ c& L* Qdiff_nsfmn1=diff_nsfmn-zero;
) M+ p; r% j+ U$ `5 _# P%*****************************************
%·*·*·*·*·*·*·*·*·*·*·*·*·*·*·*·
% M8 C0 U; I' R9 `3 `%**************时域到频域转换**************
ts=0.001;                                    %抽样间隔
/ O; D" |8 s+ j+ n2 H, ufs=1/ts;                                      %抽样频率
df=0.25;                                     %所需的频率分辨率，用在求傅里叶变换
%　　　　　　　　　　　　　　　　　　　　　　时，它表示FFT的最小频率间隔
%*****对调制信号m(t)求傅里叶变换*****
m=am*cos(2*pi*fm*t);                         %原调信号
fs=1/ts;
if nargin==2
    n1=0;
5 r, S0 |" `9 I! r7 q5 g( l4 F$ f) Jelse
    n1=fs/df;
4 m7 ]% _' ~) V2 L' n( J( ?end
n2=length(m);
n=2^(max(nextpow2(n1),nextpow2(n2)));
M=fft(m,n);
m=[m,zeros(1,n-n2)];
- Z9 }7 ]9 N& m% D% m1 Cdf1=fs/n;                                   %以上程序是对调制后的信号u求傅里变换
+ `% }+ b- q3 GM=M/fs;                                   %缩放，便于在频铺图上整体观察
5 r+ `1 [0 `0 r  `- Rf=[0:df1:df1*(length(m)-1)]-fs/2;                %时间向量对应的频率向量
: M( Q8 g; S) ~6 L0 g
1 h( B; y4 T* M3 n( s6 {%************对已调信号u求傅里变换**********
fs=1/ts;
if nargin==2
    n1=0;
+ F" S' ^" V! S$ u+ o4 y9 Q1 yelse
    n1=fs/df;
end
n2=length(sfm);
n=2^(max(nextpow2(n1),nextpow2(n2)));
( y- V9 f) o  P& PU=fft(sfm,n);
* s9 a2 o4 _5 a8 k- @! u# Fu=[sfm,zeros(1,n-n2)];
/ g; ^  n. Z8 C* k" fdf1=fs/n;                                   %以上是对已调信号u求傅里变换
U=U/fs;                                    %缩放
* Q& ]4 t3 K% m2 M4 b0 [%******************************************
%*****************************************
%·*·*·*·*·*·*·*·*·*·*·*·*·*·*·*·
%***************显示程序******************
0 u/ l) ^/ V0 U2 F$ q1 Bdisp('按任意键可以看到原调制信号、载波信号和已调信号的曲线')
6 v0 _8 U: _4 d1 x7 x$ ypause
  U' ]8 ]# `) z6 b%**************figure(1)******************
+ D; J3 w+ z, \6 w" Hfigure(1)
subplot(3,1,1);plot(t,mt);                                 %绘制调制信号的时域图
/ s' {/ f& V/ `/ K8 oxlabel('时间t');
title('调制信号的时域图');
subplot(3,1,2);plot(t,ct);                               %绘制载波的时域图
xlabel('时间t');
' t. A' v4 A* x3 G' `" |" Vtitle('载波的时域图');
  w2 N0 X. ?, a" t  S# T1 ^subplot(3,1,3);
plot(t,sfm);                                                  %绘制已调信号的时域图
xlabel('时间t');
* m+ A5 e) ?: `; S4 j1 Otitle('已调信号的时域图');
9 K) T$ E2 b9 @4 R' d2 I%******************************************
disp('按任意键可以看到原调制信号和已调信号在频域内的图形')
pause
%************figure(2)*********************
  p; _3 j6 q- [0 X, H: Zfigure(2)
subplot(2,1,1)
plot(f,abs(fftshift(M)))                        %fftshift:将FFT中的DC分量移到频谱中心
8 W+ N& P( |% L8 y2 W$ }xlabel('频率f')
' M  c1 _& a- b% btitle('原调制信号的频谱图')
+ [% L# i- }' |4 i  c& osubplot(2,1,2)
2 L2 k1 k3 u+ j9 ]0 g$ e# S2 Gplot(f,abs(fftshift(U)))
xlabel('频率f')
2 J7 p8 u8 W6 u  D1 Ititle('已调信号的频谱图')
7 G& {3 Z0 Z$ A# P5 t* F6 [%******************************************
& H/ R. Y2 w" ^* D" q, C( R0 Wdisp('按任意键可以看到原调制信号、无噪声条件下已调信号和解调信号的曲线')
pause
%**************figure(3)******************
& L' Y, `# b) K" {' Ifigure(3)
. a$ H2 x2 D3 z5 l2 n$ ^8 o. Q6 x; Z% qsubplot(3,1,1);plot(t,mt);                                  %绘制调制信号的时域图
+ \- Z7 `  y4 \, L5 H; j2 Sxlabel('时间t');
title('调制信号的时域图');
6 p7 m0 D! ]" G" J) ?subplot(3,1,2);plot(t,sfm);                              %绘制已调信号的时域图
xlabel('时间t');
6 l) ]( O! x; e! j: Mtitle('无噪声条件下已调信号的时域图');
' K; s' q: A3 _8 v4 }: Mnsfm=sfm;                              
  q2 o( _7 d% zfor i=1:length(t)-1                           %接受信号通过微分器处理
    diff_nsfm(i)=(nsfm(i+1)-nsfm(i))./dt;
end
diff_nsfmn = abs(hilbert(diff_nsfm));     %hilbert变换，求绝对值得到瞬时幅度（包络检波）
zero=(max(diff_nsfmn)-min(diff_nsfmn))/2;
9 C( C$ m7 F1 A2 U( X' c  Ediff_nsfmn1=diff_nsfmn-zero;
2 c7 ?( T9 J& G9 v: \0 i% w. ssubplot(3,1,3);                             %绘制无噪声条件下解调信号的时域图
/ n: Q; u: `) V. ~- U0 ?* Bplot((1:length(diff_nsfmn1))./1000,diff_nsfmn1./400,'r');
$ B6 L9 |/ d5 `xlabel('时间t');                          
$ U2 d9 @8 z$ h. S+ D1 }0 I3 W& ctitle('无噪声条件下解调信号的时域图');
& K5 S* E3 b6 y( A%*****************************************
3 U, [# Q- q' r5 ydisp('按任意键可以看到原调制信号、小信噪比高斯白噪声条件下已调信号和解调信号已调信号的曲线')
pause
%**************figure(4)******************
figure(4)
! E) [  e: l/ Q; asubplot(3,1,1);plot(t,mt);                                  %绘制调制信号的时域图
xlabel('时间t');
title('调制信号的时域图');
. i; R; Q6 E* cdb1=am^2/(2*(10^(sn1/10)));                   %计算对应的小信噪比高斯白躁声的方差
n1=sqrt(db1)*randn(size(t));                    %生成高斯白躁声
5 S* y# c$ F0 m) tnsfm1=n1+sfm;                              %生成含高斯白躁声的已调信号（信号通
%　　　　　　　　　　　　　　　　　　　　　　%过信道传输）
for i=1:length(t)-1                            %接受信号通过微分器处理
. t: J( K& D( f) J: p    diff_nsfm1(i)=(nsfm1(i+1)-nsfm1(i))./dt;
end
diff_nsfmn1 = abs(hilbert(diff_nsfm1));   %hilbert变换，求绝对值得到瞬时幅度（包络检波）
zero=(max(diff_nsfmn)-min(diff_nsfmn))/2;
diff_nsfmn1=diff_nsfmn1-zero;
subplot(3,1,2);
5 S6 N! I; @5 g* N( L+ b" J4 c7 Uplot(1:length(diff_nsfm),diff_nsfm);      %绘制含小信噪比高斯白噪声已调信号的时域图
1 K4 p- ~5 ]& Zxlabel('时间t');
" I) }1 h. Y  M1 [title('含小信噪比高斯白噪声已调信号的时域图');
subplot(3,1,3);                        %绘制含小信噪比高斯白噪声解调信号的时域图
plot((1:length(diff_nsfmn1))./1000,diff_nsfmn1./400,'r');
xlabel('时间t');                          
# a+ ?6 [6 a  I3 Ititle('含小信噪比高斯白噪声解调信号的时域图');
  j8 j% A. [- k3 ~2 x0 H) g%*****************************************
disp('按任意键可以看到原调制信号、大信噪比高斯白噪声条件下已调信号和解调信号已调信号的曲线')
pause
%**************figure(5)******************
; \& \: ]' P0 I6 c; y/ p! H. w; ~figure(5)
subplot(3,1,1);plot(t,mt);                                   %绘制调制信号的时域图
7 D8 H2 U+ A+ l5 ]5 exlabel('时间t');
title('调制信号的时域图');
& }# Y, B7 i" d0 }db1=am^2/(2*(10^(sn2/10)));                   %计算对应的大信噪比高斯白躁声的方差
n1=sqrt(db1)*randn(size(t));                    %生成高斯白躁声
nsfm1=n1+sfm;                              %生成含高斯白躁声的已调信号（信号通过信道传输）
for i=1:length(t)-1                             %接受信号通过微分器处理
    diff_nsfm1(i)=(nsfm1(i+1)-nsfm1(i))./dt;
3 o) T" l7 z& S* ^, o3 m9 V5 Tend
diff_nsfmn1 = abs(hilbert(diff_nsfm1));            %hilbert变换，求绝对值得到瞬时幅度（包
%　　　　　　　　　　　　　　　　　　　　　　%络检波）
zero=(max(diff_nsfmn)-min(diff_nsfmn))/2;
diff_nsfmn1=diff_nsfmn1-zero;
5 ^) P4 r) i- }0 j' i* Q. hsubplot(3,1,2);
2 s+ `' g* ^: _1 b; v9 U- uplot(1:length(diff_nsfm1),diff_nsfm1);             %绘制含大信噪比高斯白噪声已调信号
/ S5 g7 f. a' s& D% J  z3 m5 ~%　　　　　　　　　　　　　　　　　　　　　　%的时域图
xlabel('时间t');
3 {& h5 Q7 b) J1 @. Etitle('含大信噪比高斯白噪声已调信号的时域图');
: M. _1 P/ i7 h% d7 M; A8 }* `: }) Z- r8 Hsubplot(3,1,3);                                %绘制含大信噪比高斯白噪声解调信号
% W( q7 e3 p! C3 f4 `! [2 O%　　　　　　　　　　　　　　　　　　　　　　%的时域图
plot((1:length(diff_nsfmn1))./1000,diff_nsfmn1./400,'r');
xlabel('时间t');                           
title('含大信噪比高斯白噪声解调信号的时域图');
%*****************************************
%******************结 束*******************

李秋锐

。。。。。。。。。。。。。。。
5 u% l- `9 N# ~0 M* g/ ^  S% i