FM调制解调系统

李秋锐

%FM调制解调系统.m
& b' f9 ?" `% t%频率调制与解调的Matlab演示源程序
%可以任意改原调制信号函数m(t)
4 j3 r$ C0 t/ B%电子信息工程  杜怀超 091308305  
' p4 ~/ ?, {7 q4 R%·*·*·*·*·*·*·*·*·*·*·*·*·*·*·*·
$ U/ |5 K6 n" D! X/ o# @5 ]%*****************初始化******************
echo   off      
6 O! A& {4 G' E0 i/ v! `close  all
clear  all
clc
%*****************************************
%·*·*·*·*·*·*·*·*·*·*·*·*·*·*·*·
%****************FM调制*******************
& U1 w/ q* _& z& W* y, y6 B- fdt=0.001;                                   %设定时间步长
% d. P5 F1 w) bt=0:dt:1.5;                                       %产生时间向量
0 p+ U/ |' k# }+ ^3 P( f) {am=5;                                      %设定调制信号幅度
fm=5;                                      %设定调制信号频率
mt=am*cos(2*pi*fm*t);                        %生成调制信号
fc=50;                                      %设定载波频率
ct=cos(2*pi*fc*t);                             %生成载波
kf=10;                                      %设定调频指数
/ z" x6 ~  w2 R  t, v  Dint_mt(1)=0;
for i=1:length(t)-1  
4 _6 N, N0 G- K3 i1 Z8 B    int_mt(i+1)=int_mt(i)+mt(i)*dt;               %求信号m(t)的积分
) B) @) v% z5 `( B' yend                                          %调制，产生已调信号
sfm=am*cos(2*pi*fc*t+2*pi*kf*int_mt);            %调制信号
%*****************************************
%·*·*·*·*·*·*·*·*·*·*·*·*·*·*·*·
) A8 @+ O) t: K7 w1 r9 V  ~$ K% C%*************添加高斯白噪声**************
  r6 v2 O" O% D* x' i$ xsn1=10;                                     %设定信躁比(小信噪比)
# J9 X  n, h: B8 \sn2=30;                                     %设定信躁比(大信噪比)
sn=0;                                       %设定信躁比(无信噪比)
db=am^2/(2*(10^(sn/10)));                     %计算对应的高斯白躁声的方差
! N5 _" H, [! k3 D8 ~9 zn=sqrt(db)*randn(size(t));                      %生成高斯白躁声
- _: D' C! I3 ^2 |6 W% Mnsfm=n+sfm;                                %生成含高斯白躁声的已调信号（信号通
%过信道传输）
%*****************************************
%·*·*·*·*·*·*·*·*·*·*·*·*·*·*·*·
%****************FM解调*******************
# T* M2 h! b& h, f% l/ V) K' p2 Sfor i=1:length(t)-1                             %接受信号通过微分器处理
    diff_nsfm(i)=(nsfm(i+1)-nsfm(i))./dt;
end
diff_nsfmn = abs(hilbert(diff_nsfm));     %hilbert变换，求绝对值得到瞬时幅度（包络检波）
3 m. _$ Y+ g! @% hzero=(max(diff_nsfmn)-min(diff_nsfmn))/2;
( X- }7 X$ t/ Tdiff_nsfmn1=diff_nsfmn-zero;
. w( ~: j8 w& P6 a9 f$ B%*****************************************
%·*·*·*·*·*·*·*·*·*·*·*·*·*·*·*·
%**************时域到频域转换**************
ts=0.001;                                    %抽样间隔
fs=1/ts;                                      %抽样频率
df=0.25;                                     %所需的频率分辨率，用在求傅里叶变换
1 j3 ~; i3 T& e! l9 d. P& E%　　　　　　　　　　　　　　　　　　　　　　时，它表示FFT的最小频率间隔
%*****对调制信号m(t)求傅里叶变换*****
: J6 O& i, M* k4 w0 t/ Wm=am*cos(2*pi*fm*t);                         %原调信号
fs=1/ts;
if nargin==2
    n1=0;
9 n! Q8 N( K) u& x$ v& gelse
    n1=fs/df;
9 `% o9 e* [/ u( G: N( W: ^) H2 x& Kend
n2=length(m);
n=2^(max(nextpow2(n1),nextpow2(n2)));
8 r0 ~" Y- T" K/ k& r* K; r6 p; HM=fft(m,n);
m=[m,zeros(1,n-n2)];
2 a: ^$ u9 h* n5 Idf1=fs/n;                                   %以上程序是对调制后的信号u求傅里变换
M=M/fs;                                   %缩放，便于在频铺图上整体观察
7 N' t' x/ I5 s$ O% w6 J& }: M/ Zf=[0:df1:df1*(length(m)-1)]-fs/2;                %时间向量对应的频率向量
2 {8 Q5 P$ A- p
. L5 o2 b; S9 J! u%************对已调信号u求傅里变换**********
fs=1/ts;
if nargin==2
. I/ O: ]% p. ?0 i0 Z    n1=0;
* [. B9 v4 r: G! `, velse
    n1=fs/df;
end
n2=length(sfm);
5 ^% o6 D: }5 |" An=2^(max(nextpow2(n1),nextpow2(n2)));
U=fft(sfm,n);
( N3 K7 o% I0 S7 U1 Y; |8 ?  Pu=[sfm,zeros(1,n-n2)];
df1=fs/n;                                   %以上是对已调信号u求傅里变换
U=U/fs;                                    %缩放
8 w) j6 A6 w2 O# O. {& K: J%******************************************
. b  W: A: j  z%*****************************************
( q6 D: r, u: A4 F%·*·*·*·*·*·*·*·*·*·*·*·*·*·*·*·
%***************显示程序******************
disp('按任意键可以看到原调制信号、载波信号和已调信号的曲线')
+ [4 a+ G7 x9 ^3 s5 fpause
%**************figure(1)******************
figure(1)
9 B) S8 E; z: z6 fsubplot(3,1,1);plot(t,mt);                                 %绘制调制信号的时域图
xlabel('时间t');
$ G, W: V8 B6 X0 U. |. Gtitle('调制信号的时域图');
8 ^: I, l- L2 w% O( e& _, bsubplot(3,1,2);plot(t,ct);                               %绘制载波的时域图
0 C5 S( h- l7 ], M* P9 exlabel('时间t');
title('载波的时域图');
subplot(3,1,3);
  c9 p2 V( D# h  Tplot(t,sfm);                                                  %绘制已调信号的时域图
9 y) @! n2 r- M, s9 y. r2 r+ n( ]xlabel('时间t');
# i* l$ ~( T' i3 v6 V3 Mtitle('已调信号的时域图');
7 v6 q: C/ Y4 D+ b%******************************************
1 c: X! M: A: z6 h( v" ], Kdisp('按任意键可以看到原调制信号和已调信号在频域内的图形')
; N! u4 g/ i7 Y/ qpause
1 t+ K- W) u: W! N4 R" }%************figure(2)*********************
) T/ |# n; P! g4 q3 u2 Y* sfigure(2)
subplot(2,1,1)
: }2 w. I$ g3 O0 O) |plot(f,abs(fftshift(M)))                        %fftshift:将FFT中的DC分量移到频谱中心
! z6 L/ d# i: R9 v4 A3 kxlabel('频率f')
title('原调制信号的频谱图')
1 N/ M# }2 I% Ssubplot(2,1,2)
8 T: c0 ?3 Y# D) y! kplot(f,abs(fftshift(U)))
: ~, A4 h% {0 }6 |xlabel('频率f')
title('已调信号的频谱图')
9 p+ O3 T$ K0 r%******************************************
disp('按任意键可以看到原调制信号、无噪声条件下已调信号和解调信号的曲线')
pause
) B, F! @, Z2 @9 f4 G%**************figure(3)******************
figure(3)
/ v3 h% z8 x" }- x2 a; b& e& L+ [subplot(3,1,1);plot(t,mt);                                  %绘制调制信号的时域图
xlabel('时间t');
; I0 v9 g. z2 o* A( t" Ctitle('调制信号的时域图');
subplot(3,1,2);plot(t,sfm);                              %绘制已调信号的时域图
1 j4 {* V- i7 a- z; q+ x0 G0 rxlabel('时间t');
title('无噪声条件下已调信号的时域图');
6 [- e0 G7 p1 Q1 o1 b  X& bnsfm=sfm;                              
for i=1:length(t)-1                           %接受信号通过微分器处理
! N; L( P6 u( b, V" b7 D& s    diff_nsfm(i)=(nsfm(i+1)-nsfm(i))./dt;
4 s+ h- G8 D: [7 c0 m: Pend
diff_nsfmn = abs(hilbert(diff_nsfm));     %hilbert变换，求绝对值得到瞬时幅度（包络检波）
  a4 ]8 j; h. b+ D% `9 Jzero=(max(diff_nsfmn)-min(diff_nsfmn))/2;
5 p9 M. P3 i+ g# O  x& xdiff_nsfmn1=diff_nsfmn-zero;
subplot(3,1,3);                             %绘制无噪声条件下解调信号的时域图
plot((1:length(diff_nsfmn1))./1000,diff_nsfmn1./400,'r');
/ F" ~5 v5 i7 }0 Hxlabel('时间t');                          
. f1 s$ @/ \3 \& {; I) Wtitle('无噪声条件下解调信号的时域图');
%*****************************************
disp('按任意键可以看到原调制信号、小信噪比高斯白噪声条件下已调信号和解调信号已调信号的曲线')
pause
; Q: }" D7 i0 g) L# _8 H%**************figure(4)******************
figure(4)
2 K6 ^% }! F  Z7 Y& ?; B+ osubplot(3,1,1);plot(t,mt);                                  %绘制调制信号的时域图
xlabel('时间t');
% G4 F: ~* ]8 d  W4 g! D$ }title('调制信号的时域图');
+ C7 r6 g1 X% W4 R: o$ G: w3 M+ R- I* J: Ldb1=am^2/(2*(10^(sn1/10)));                   %计算对应的小信噪比高斯白躁声的方差
n1=sqrt(db1)*randn(size(t));                    %生成高斯白躁声
nsfm1=n1+sfm;                              %生成含高斯白躁声的已调信号（信号通
, d! j; M, ^0 A4 m) }%　　　　　　　　　　　　　　　　　　　　　　%过信道传输）
for i=1:length(t)-1                            %接受信号通过微分器处理
    diff_nsfm1(i)=(nsfm1(i+1)-nsfm1(i))./dt;
6 ~# x! b) \" k. I) w, @, ^8 |7 Y% Send
! n* K. e7 v: ?  P% `- a& Ddiff_nsfmn1 = abs(hilbert(diff_nsfm1));   %hilbert变换，求绝对值得到瞬时幅度（包络检波）
zero=(max(diff_nsfmn)-min(diff_nsfmn))/2;
diff_nsfmn1=diff_nsfmn1-zero;
subplot(3,1,2);
& V* W- q$ H* i- nplot(1:length(diff_nsfm),diff_nsfm);      %绘制含小信噪比高斯白噪声已调信号的时域图
7 f. I# J/ @2 B! J2 R0 ~' K5 [/ a  ?xlabel('时间t');
title('含小信噪比高斯白噪声已调信号的时域图');
subplot(3,1,3);                        %绘制含小信噪比高斯白噪声解调信号的时域图
plot((1:length(diff_nsfmn1))./1000,diff_nsfmn1./400,'r');
, G$ R2 ]& m" F6 Y* C: Kxlabel('时间t');                          
) s0 H! @+ T2 S0 xtitle('含小信噪比高斯白噪声解调信号的时域图');
+ R- Z9 j7 `4 r%*****************************************
disp('按任意键可以看到原调制信号、大信噪比高斯白噪声条件下已调信号和解调信号已调信号的曲线')
) g5 y* O( L5 J( xpause
%**************figure(5)******************
figure(5)
subplot(3,1,1);plot(t,mt);                                   %绘制调制信号的时域图
; o+ {) r1 A7 G" U8 {% pxlabel('时间t');
title('调制信号的时域图');
9 x1 \  G/ t* s3 d$ bdb1=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;
end
+ `! K' R& D! ?/ J8 t6 cdiff_nsfmn1 = abs(hilbert(diff_nsfm1));            %hilbert变换，求绝对值得到瞬时幅度（包
%　　　　　　　　　　　　　　　　　　　　　　%络检波）
& E2 ]) v9 }3 B$ \1 ?: Lzero=(max(diff_nsfmn)-min(diff_nsfmn))/2;
diff_nsfmn1=diff_nsfmn1-zero;
) |$ T  y/ j* A6 b# ?4 Dsubplot(3,1,2);
plot(1:length(diff_nsfm1),diff_nsfm1);             %绘制含大信噪比高斯白噪声已调信号
%　　　　　　　　　　　　　　　　　　　　　　%的时域图
) ]. R- L8 k; Q9 o1 p) y1 f# sxlabel('时间t');
  C; P! ?/ `8 Z3 i  ]( Ztitle('含大信噪比高斯白噪声已调信号的时域图');
subplot(3,1,3);                                %绘制含大信噪比高斯白噪声解调信号
%　　　　　　　　　　　　　　　　　　　　　　%的时域图
7 v5 T- ?9 o* R0 X: }  ?4 L7 P$ Dplot((1:length(diff_nsfmn1))./1000,diff_nsfmn1./400,'r');
6 j, o1 M6 }% _xlabel('时间t');                           
' {8 g) S* @9 I" W# w8 x/ i( V  xtitle('含大信噪比高斯白噪声解调信号的时域图');
  q: M" F! p/ ^) A  E6 w%*****************************************
$ z$ {) \8 A4 F$ l1 X/ O%******************结 束*******************
. V9 y+ e1 M! G! E- {

李秋锐

# l( v# E" `& D+ r% C
$ w+ c+ ]% _9 ?" |+ b$ X, P* E。。。。。。。。。。。。。。。

" ]0 S' I# }) T8 \