数学建模社区-数学中国

标题: FM调制解调系统 [打印本页]
作者: 李秋锐    时间: 2015-2-1 19:22
标题: FM调制解调系统
%FM调制解调系统.m
% Y  }4 X* p. y! I) S3 X%频率调制与解调的Matlab演示源程序
! G* w9 f- t& F3 E6 j- I%可以任意改原调制信号函数m(t)
5 j8 j4 Q  |, @9 {" e0 _+ C+ O1 i%电子信息工程  杜怀超 091308305  ! j3 |7 h( Z% U
%·*·*·*·*·*·*·*·*·*·*·*·*·*·*·*·6 x# z. }3 c, n) b+ A: p3 S8 I9 g
%*****************初始化******************
/ J1 w$ X* a, f% Y) Oecho   off      
5 o, f4 e5 Q9 g) qclose  all
2 n# x/ p( d; {) l- U7 pclear  all
/ i/ ]& l0 k# T9 \# @* qclc
" g7 A3 _: I9 J  X5 X4 b%*****************************************
$ u6 Y9 O( a. H) L%·*·*·*·*·*·*·*·*·*·*·*·*·*·*·*·! S; j! b1 Z' s, w2 z0 M) e7 y
%****************FM调制*******************# ]3 Y2 B, _  ~8 i- j0 R6 M
dt=0.001;                                   %设定时间步长: m% v" _; x9 {0 @; W  A
t=0:dt:1.5;                                       %产生时间向量" R3 P0 N: n) N; S: U8 o
am=5;                                      %设定调制信号幅度. Q4 ^+ m( m  c! [
fm=5;                                      %设定调制信号频率
) o' a3 J( {: qmt=am*cos(2*pi*fm*t);                        %生成调制信号
& o" d' [8 u0 C( zfc=50;                                      %设定载波频率% F6 }, f( C" q7 Z, c% M
ct=cos(2*pi*fc*t);                             %生成载波
: I2 @! K9 U2 Q1 x9 a, \' Y7 xkf=10;                                      %设定调频指数
% j7 X$ O5 C) q2 e* R8 N, Kint_mt(1)=0;0 G3 W) V( z& K8 q
for i=1:length(t)-1  ( b+ F' q0 `) w* g& k
    int_mt(i+1)=int_mt(i)+mt(i)*dt;               %求信号m(t)的积分
* u% }3 t2 E8 y' k* w% K0 }6 jend                                          %调制，产生已调信号
# C$ `0 N5 e3 P2 n8 gsfm=am*cos(2*pi*fc*t+2*pi*kf*int_mt);            %调制信号1 W& I5 O8 ?: c, ?
%*****************************************8 |" U' X0 w1 p* u( I1 p9 j
%·*·*·*·*·*·*·*·*·*·*·*·*·*·*·*·, N$ I  m' e& M5 W8 p
%*************添加高斯白噪声**************
3 G2 N9 v: \" B& O2 Ssn1=10;                                     %设定信躁比(小信噪比)/ G" f8 _' q2 i% [8 C
sn2=30;                                     %设定信躁比(大信噪比)
" t- E7 q' f! k9 v% e1 Fsn=0;                                       %设定信躁比(无信噪比)' z& ^5 o# `2 Q2 k: L9 y* A0 l3 z
db=am^2/(2*(10^(sn/10)));                     %计算对应的高斯白躁声的方差
3 g* }2 y* U6 |* nn=sqrt(db)*randn(size(t));                      %生成高斯白躁声
2 l  O( g( ~2 y/ B& vnsfm=n+sfm;                                %生成含高斯白躁声的已调信号（信号通
8 s5 I7 `9 @( U: w6 J" S4 ?: v%过信道传输）
9 d/ O% V, d0 l( y7 L%*****************************************
' ~( E6 P, i# C1 h! l# z6 x" o8 o%·*·*·*·*·*·*·*·*·*·*·*·*·*·*·*·1 O5 l/ m5 Q' [  `1 z8 k
%****************FM解调*******************- L* s6 ?8 i0 {- b! ]# N
for i=1:length(t)-1                             %接受信号通过微分器处理
& V( n6 ~2 g3 K9 B- n    diff_nsfm(i)=(nsfm(i+1)-nsfm(i))./dt;
( S/ D0 C- K. Jend
  m1 F4 V6 z: ^* j" Hdiff_nsfmn = abs(hilbert(diff_nsfm));     %hilbert变换，求绝对值得到瞬时幅度（包络检波）$ q8 J# ?. e3 J+ z" K0 ~9 d
zero=(max(diff_nsfmn)-min(diff_nsfmn))/2;
) h) ~( T2 }" L) v9 w) w( c  Ddiff_nsfmn1=diff_nsfmn-zero;- G0 p5 y4 t7 y0 ~/ W# F2 Z
%*****************************************. U! n, w' h- l
%·*·*·*·*·*·*·*·*·*·*·*·*·*·*·*·5 n8 E8 z; K- U. U7 h
%**************时域到频域转换**************
0 ^& J' E4 X& y( ^4 @1 V9 O5 Qts=0.001;                                    %抽样间隔( `3 h9 a. J! g0 V/ K, a" g4 L# w9 q
fs=1/ts;                                      %抽样频率% j  Y# U4 T9 Y: D; T( V8 w% w
df=0.25;                                     %所需的频率分辨率，用在求傅里叶变换
* m" t+ O! ^1 w%　　　　　　　　　　　　　　　　　　　　　　时，它表示FFT的最小频率间隔' k6 s: H5 U6 P( t6 ~
%*****对调制信号m(t)求傅里叶变换*****5 N+ a0 f1 J" b: w1 r; s* p
m=am*cos(2*pi*fm*t);                         %原调信号
1 r0 {8 @0 Y2 b; O4 c0 g4 [9 Hfs=1/ts;/ L: |8 k7 g  e% p' `1 r
if nargin==23 G1 A# L4 B0 }0 i: H
    n1=0;
' N5 z5 [' ~! G' z1 T2 oelse
4 \! H* a; c; j8 y8 V  O7 T2 H    n1=fs/df;. J* ^% L/ ^: m! e+ x2 J; L
end
0 l* N  [% C) Nn2=length(m);
7 m8 ]5 S4 A/ o  m" ]( [n=2^(max(nextpow2(n1),nextpow2(n2)));
8 w8 U$ B3 l! p  A" UM=fft(m,n);- ^* f) H4 @% w
m=[m,zeros(1,n-n2)];
- d6 d6 J; Z0 sdf1=fs/n;                                   %以上程序是对调制后的信号u求傅里变换) w9 s) C1 K1 W& r
M=M/fs;                                   %缩放，便于在频铺图上整体观察
! B- e# v2 M9 p; n1 pf=[0:df1:df1*(length(m)-1)]-fs/2;                %时间向量对应的频率向量
* u3 h  t  k$ o% Z. i; D2 i, c. {! u6 z3 {5 ?- C$ _" D
%************对已调信号u求傅里变换**********
6 {. C2 D2 @" g- Gfs=1/ts;# W: L. h; R  B7 G2 N+ g9 p7 W
if nargin==2
9 [. w1 b0 \5 M1 t/ |) f+ S    n1=0;6 S% E* \2 X/ |, b7 f: K/ E; R3 t  O
else
" R6 Y+ a$ i  M3 @8 y    n1=fs/df;* D" g6 u3 M& }* Z1 |' l
end& P/ F* `5 o) k; r+ H3 E! R. X# J
n2=length(sfm);
0 `6 ?  i* E% `% in=2^(max(nextpow2(n1),nextpow2(n2)));
; [' ?7 i/ H3 b& WU=fft(sfm,n);
- m4 x8 N! Z+ }2 t% ?2 |! mu=[sfm,zeros(1,n-n2)];
4 b, O1 O, V0 ?4 hdf1=fs/n;                                   %以上是对已调信号u求傅里变换/ k1 a! L& @5 Q* Z+ w
U=U/fs;                                    %缩放
% U# r; @& c. z6 e2 j%******************************************1 c; I  z1 Z! z' l7 g
%*****************************************
! ~" U/ D$ ~5 o/ `%·*·*·*·*·*·*·*·*·*·*·*·*·*·*·*·6 c0 Z/ g3 _- X  n$ Q
%***************显示程序******************9 C+ P( u& S3 _
disp('按任意键可以看到原调制信号、载波信号和已调信号的曲线')
! D  l1 b1 Q5 _pause' r! E. I8 ^1 h% J
%**************figure(1)******************; _  G% l, X7 [0 O5 P4 h
figure(1)( j$ [7 T7 P9 J' ?- S5 w
subplot(3,1,1);plot(t,mt);                                 %绘制调制信号的时域图
. J; q1 k/ M" xxlabel('时间t');
  d+ T1 p- l3 M* k2 Ptitle('调制信号的时域图');
! l- i- m% J. ~1 P% o6 M2 Ssubplot(3,1,2);plot(t,ct);                               %绘制载波的时域图) t1 Q# K0 U8 h* |: \
xlabel('时间t');
7 h9 V# F6 f" R8 C2 ~7 N! c& xtitle('载波的时域图');
3 I! t: O3 d, N# t; z" U8 w, gsubplot(3,1,3);
7 M5 {7 h" o( b- R% Y  Lplot(t,sfm);                                                  %绘制已调信号的时域图
: ?* {5 f! A, q* x- c- G1 _xlabel('时间t');
2 V% L7 E" h1 T5 R3 gtitle('已调信号的时域图');
9 x0 y3 d) \3 O4 r- T0 X8 n; a%******************************************
$ s6 m2 x. R: q) k0 f7 V* n) A6 Odisp('按任意键可以看到原调制信号和已调信号在频域内的图形')
! w  s0 G  o! |; i. ~! rpause/ ?# ?( x1 R8 |7 C2 C5 e. {8 E
%************figure(2)*********************
! {- ]! K. H$ K0 q+ A, `figure(2)9 z. m0 R  r, \4 X  m
subplot(2,1,1). m( I* T, k1 r! ~: w
plot(f,abs(fftshift(M)))                        %fftshift:将FFT中的DC分量移到频谱中心' {( n5 C. w$ y1 Z  I0 v
xlabel('频率f')+ z1 s1 X; @3 U1 e* K8 |
title('原调制信号的频谱图')- D" p+ h1 [: C9 g$ r8 m
subplot(2,1,2)
- ~& m8 G9 W+ s1 l$ |plot(f,abs(fftshift(U)))! Q' {7 O! q( d2 z7 v" j. V- p
xlabel('频率f')
0 I; w4 y' A! l: o6 a; K8 o* M; t  Gtitle('已调信号的频谱图')% g% O( m' G4 M0 T$ w; S( O' |
%******************************************# `6 ?' v( O1 D/ y9 p) g7 r
disp('按任意键可以看到原调制信号、无噪声条件下已调信号和解调信号的曲线')* x2 j: O! j" f" p! ]6 O* ?3 G) C$ m
pause
" V- N2 B: p; }; H, s5 @9 i0 q%**************figure(3)******************
$ X+ j- b/ G3 Afigure(3): _! U; j/ C" [
subplot(3,1,1);plot(t,mt);                                  %绘制调制信号的时域图( U9 P8 k0 b! z) X
xlabel('时间t');" Y, ^6 e; w, m6 l
title('调制信号的时域图');6 A( X0 U9 p* Y5 v; U9 |) Y! d& p
subplot(3,1,2);plot(t,sfm);                              %绘制已调信号的时域图& G) J' e7 Z7 l, j
xlabel('时间t');
4 N  A4 g0 H) F* ^title('无噪声条件下已调信号的时域图');
: i) ^/ n+ i; p' l) H0 X' t9 Ynsfm=sfm;                              
8 W0 p* Q% |4 f* f' mfor i=1:length(t)-1                           %接受信号通过微分器处理
$ L2 s8 Z  m! e4 X' g5 S# \! ^    diff_nsfm(i)=(nsfm(i+1)-nsfm(i))./dt;. I# Y( s9 ~( m. X3 C6 n
end
4 V7 h) ]  J! |$ C% hdiff_nsfmn = abs(hilbert(diff_nsfm));     %hilbert变换，求绝对值得到瞬时幅度（包络检波）/ w1 `3 \9 H4 X: F0 q- I+ x7 B
zero=(max(diff_nsfmn)-min(diff_nsfmn))/2;6 G1 {- ^$ f8 P# w0 o9 D2 W1 Z4 P2 _  o. H
diff_nsfmn1=diff_nsfmn-zero;
% P! z! K4 l* t8 t/ D; esubplot(3,1,3);                             %绘制无噪声条件下解调信号的时域图' o2 G- P, E. Q% h
plot((1:length(diff_nsfmn1))./1000,diff_nsfmn1./400,'r');
. Y1 x: _$ Y3 _xlabel('时间t');                          9 u# n  G5 J( l* N
title('无噪声条件下解调信号的时域图');
8 z! Y) [( l$ \! V9 ]% W%*****************************************
4 V& ^3 }1 L  s9 jdisp('按任意键可以看到原调制信号、小信噪比高斯白噪声条件下已调信号和解调信号已调信号的曲线')- r# Q5 F9 [) m8 V5 z# w
pause3 }1 R% Z$ e) B& l  j7 [( _
%**************figure(4)******************
* V9 t+ _' ?( afigure(4); i! U  g5 m* ]+ C6 Z
subplot(3,1,1);plot(t,mt);                                  %绘制调制信号的时域图: {4 v& p) ~0 Y- H) U$ r
xlabel('时间t');
: a+ x& Z" j6 C! M' rtitle('调制信号的时域图');2 h6 g3 A0 C: J) v/ I6 ^
db1=am^2/(2*(10^(sn1/10)));                   %计算对应的小信噪比高斯白躁声的方差: I5 j; ]8 s$ b/ l: V
n1=sqrt(db1)*randn(size(t));                    %生成高斯白躁声
7 u5 @! P* ]) r; B5 L7 gnsfm1=n1+sfm;                              %生成含高斯白躁声的已调信号（信号通* E! C8 V8 A' ~& n
%　　　　　　　　　　　　　　　　　　　　　　%过信道传输）7 `+ G3 u! Q( t! w& U
for i=1:length(t)-1                            %接受信号通过微分器处理' T" O' L3 t5 @
    diff_nsfm1(i)=(nsfm1(i+1)-nsfm1(i))./dt;- Q7 H4 y6 O" |2 d3 K4 ^: a$ a
end
9 c5 {0 Y3 g  O! u1 ediff_nsfmn1 = abs(hilbert(diff_nsfm1));   %hilbert变换，求绝对值得到瞬时幅度（包络检波）
( ?  k6 e* e  ~2 a  S" wzero=(max(diff_nsfmn)-min(diff_nsfmn))/2;
0 W5 V- m2 x4 o/ D* Z) W5 m$ x4 ldiff_nsfmn1=diff_nsfmn1-zero;# ?. ]" O( x1 p% j7 Q" z
subplot(3,1,2);* k- m- c" I6 L$ Z
plot(1:length(diff_nsfm),diff_nsfm);      %绘制含小信噪比高斯白噪声已调信号的时域图
- e% z+ r" e5 U/ ?, r+ Txlabel('时间t');7 U+ b, E7 [) v' _% K
title('含小信噪比高斯白噪声已调信号的时域图');
. @3 u; X/ A8 {& isubplot(3,1,3);                        %绘制含小信噪比高斯白噪声解调信号的时域图
( P/ D$ b( m* q5 Nplot((1:length(diff_nsfmn1))./1000,diff_nsfmn1./400,'r');
4 A7 t  w+ ]' ?1 H& t2 U5 lxlabel('时间t');                          - @2 |9 h, l2 M# M" u
title('含小信噪比高斯白噪声解调信号的时域图');! f( P3 h& f4 C
%*****************************************" d; \8 x7 b* R, U8 Y( z- L
disp('按任意键可以看到原调制信号、大信噪比高斯白噪声条件下已调信号和解调信号已调信号的曲线')
  ^$ l6 O9 H( g5 w( l, q& C0 _pause& G% e' O  A- e' N0 V+ k8 W
%**************figure(5)******************
4 K5 R7 }$ L% M% nfigure(5)
7 E- i* X6 U/ a7 U& s' msubplot(3,1,1);plot(t,mt);                                   %绘制调制信号的时域图
, G& i: u6 h. F4 |" P/ O1 j3 \4 `xlabel('时间t');: J1 |/ m* T+ E8 d: p
title('调制信号的时域图');
; T5 g: Z9 x, I# Fdb1=am^2/(2*(10^(sn2/10)));                   %计算对应的大信噪比高斯白躁声的方差! T% Z: I8 h) i, V# Y# o
n1=sqrt(db1)*randn(size(t));                    %生成高斯白躁声
# m, W, z" x' k. lnsfm1=n1+sfm;                              %生成含高斯白躁声的已调信号（信号通过信道传输）
; |. {+ L5 M5 d: F1 Jfor i=1:length(t)-1                             %接受信号通过微分器处理1 E; i3 |0 |* ~, [0 \$ ?, n
    diff_nsfm1(i)=(nsfm1(i+1)-nsfm1(i))./dt;
: v1 t- F0 J8 b# E. Rend3 }( [; ^0 F' G$ ?
diff_nsfmn1 = abs(hilbert(diff_nsfm1));            %hilbert变换，求绝对值得到瞬时幅度（包
' c+ ^) g" n% r! D0 G/ K%　　　　　　　　　　　　　　　　　　　　　　%络检波）
; F2 o- @2 N/ n) L  Wzero=(max(diff_nsfmn)-min(diff_nsfmn))/2;
8 @+ h3 i5 A. s- c( y1 M! C- Bdiff_nsfmn1=diff_nsfmn1-zero;: M! j8 L  T8 I/ h( x
subplot(3,1,2);
% C0 c( J! T+ m' g. r5 T" J, gplot(1:length(diff_nsfm1),diff_nsfm1);             %绘制含大信噪比高斯白噪声已调信号
9 R8 P$ q; P9 L& @+ F%　　　　　　　　　　　　　　　　　　　　　　%的时域图
  Q2 r3 r: {5 s( K8 o& R8 H1 cxlabel('时间t');% G* v3 b/ I8 j7 H8 ^3 E
title('含大信噪比高斯白噪声已调信号的时域图');
' C# c6 e" g  U" ~9 I( T4 j( {! \subplot(3,1,3);                                %绘制含大信噪比高斯白噪声解调信号' C8 j: D  l' |& A/ v) J4 I
%　　　　　　　　　　　　　　　　　　　　　　%的时域图
! \1 K) a) y  P# A& ~; M+ Nplot((1:length(diff_nsfmn1))./1000,diff_nsfmn1./400,'r');2 b4 u) n- u- I/ D) ~
xlabel('时间t');                            , h( H  Z2 b: o2 E7 \/ p
title('含大信噪比高斯白噪声解调信号的时域图');
8 A, q9 Z" w# @%*****************************************
3 {) |' Q) {' R%******************结 束*******************) s+ A7 P: y$ }. o: o! t& V
, [) A7 v/ ^8 x# S) Y, `
& i( }: i, t7 B, v) A; z+ ]
作者: 李秋锐    时间: 2015-2-1 19:23
+ ^4 T7 m, i- v$ z

5 r. m* t( `4 F$ E' z
; \# m: |7 m1 O6 e, |。。。。。。。。。。。。。。。) T6 e9 y; ~2 D1 R, l7 y; Y& ]
' p' a8 c2 y! j( J7 [* [8 Z




欢迎光临 数学建模社区-数学中国 (http://www.madio.net/) Powered by Discuz! X2.5