数学建模社区-数学中国

标题: FM调制解调系统 [打印本页]
作者: 李秋锐    时间: 2015-2-1 19:22
标题: FM调制解调系统
%FM调制解调系统.m
7 W; h0 P- [8 O# u2 A%频率调制与解调的Matlab演示源程序, n# Q( ^5 {# F0 u8 B) n
%可以任意改原调制信号函数m(t)
% Q# [9 w* J0 w+ H%电子信息工程  杜怀超 091308305  
8 I* I2 s5 v$ R* {: H%·*·*·*·*·*·*·*·*·*·*·*·*·*·*·*·- E" s4 }6 P" O# |( N  v0 V
%*****************初始化******************
. e2 I' U% Z4 z$ ]5 n7 i/ Z! Aecho   off       - o8 ~2 ]% A/ H, o/ q
close  all' P$ V+ Q, ^9 M; i* u4 H
clear  all
& ^4 r. D) n% s( S  Z. tclc/ P, }5 c1 K( c, O3 U  x( o) z
%*****************************************0 M* m2 k8 i9 v5 j
%·*·*·*·*·*·*·*·*·*·*·*·*·*·*·*·9 o5 @8 K, @( J
%****************FM调制*******************. I( v. v& m6 b" Y9 |3 Q( D; {  J4 M* @
dt=0.001;                                   %设定时间步长9 a3 \% r2 K' P* R
t=0:dt:1.5;                                       %产生时间向量: w2 a- L- \- ]# Z. r% E. E; v
am=5;                                      %设定调制信号幅度+ V, q* D# g8 M# n% w5 l
fm=5;                                      %设定调制信号频率
4 j' Y, O/ B$ f% x; t9 \4 w( Qmt=am*cos(2*pi*fm*t);                        %生成调制信号
# q- Q1 n* R" k2 ?# g/ P/ L3 kfc=50;                                      %设定载波频率: w1 X' W* ], P% W$ Z8 V4 W
ct=cos(2*pi*fc*t);                             %生成载波% q; Y9 C/ |. c; c' H. O) s
kf=10;                                      %设定调频指数7 O6 q% h1 f' \, ?
int_mt(1)=0;0 a2 g! b; N* B8 D' w6 o
for i=1:length(t)-1  
: d5 a4 I  [0 {" S0 e' t    int_mt(i+1)=int_mt(i)+mt(i)*dt;               %求信号m(t)的积分
8 a1 T- g5 F, K+ H0 z$ Jend                                          %调制，产生已调信号$ C: d& j5 h" E& p
sfm=am*cos(2*pi*fc*t+2*pi*kf*int_mt);            %调制信号# W* w0 ?# h" B3 g8 H) z
%*****************************************
  i7 O& S9 n! R. W! b%·*·*·*·*·*·*·*·*·*·*·*·*·*·*·*·. w# T3 |5 w/ @7 q4 X/ h
%*************添加高斯白噪声**************
# s1 s* Y/ H' S* ^" D7 Nsn1=10;                                     %设定信躁比(小信噪比)0 f& @& v3 Z; Y9 ]
sn2=30;                                     %设定信躁比(大信噪比)
! e2 e  E4 n$ p! P/ `0 x: Gsn=0;                                       %设定信躁比(无信噪比)
$ A% g4 D. X1 z! F7 N% Mdb=am^2/(2*(10^(sn/10)));                     %计算对应的高斯白躁声的方差
, l; r# u6 T6 t0 S% L1 N! Q' j7 Y' ~n=sqrt(db)*randn(size(t));                      %生成高斯白躁声
- }* h* U# z, {nsfm=n+sfm;                                %生成含高斯白躁声的已调信号（信号通
: A+ G7 j7 i4 p& g" P%过信道传输）
- u' ?1 C! d7 q%*****************************************, a3 j# s' d+ x! Z* [
%·*·*·*·*·*·*·*·*·*·*·*·*·*·*·*·8 }8 T& V. c/ I- {( l
%****************FM解调*******************
) u4 o% k& }/ B3 ifor i=1:length(t)-1                             %接受信号通过微分器处理/ G, j( q$ l* r6 \+ \
    diff_nsfm(i)=(nsfm(i+1)-nsfm(i))./dt;6 T( v1 a* B7 |( [/ ~" A5 p5 j
end; a6 g$ ~4 N- s, R! W- E
diff_nsfmn = abs(hilbert(diff_nsfm));     %hilbert变换，求绝对值得到瞬时幅度（包络检波）
* L6 p+ H% W5 y' ]- g0 }) I2 mzero=(max(diff_nsfmn)-min(diff_nsfmn))/2;
# |9 `8 F8 p0 \& `  ydiff_nsfmn1=diff_nsfmn-zero;
( r1 p. w2 T0 y%*****************************************2 w9 i/ P1 E! l' @+ s9 }
%·*·*·*·*·*·*·*·*·*·*·*·*·*·*·*·
7 I1 T7 J: x: D( X/ m%**************时域到频域转换**************" I) V, G; S- ^  N# d  t' ~4 G
ts=0.001;                                    %抽样间隔. I+ h# w$ ^+ ^
fs=1/ts;                                      %抽样频率1 H/ M7 G" [7 |9 t! }4 R/ P4 p( S
df=0.25;                                     %所需的频率分辨率，用在求傅里叶变换
9 |! A, F- ^) ]1 y) X1 b%　　　　　　　　　　　　　　　　　　　　　　时，它表示FFT的最小频率间隔. e: S  o1 u  y. c3 R
%*****对调制信号m(t)求傅里叶变换*****. r* G9 s  d' |
m=am*cos(2*pi*fm*t);                         %原调信号2 v/ D- S* q' E' b( ~7 G! Q
fs=1/ts;4 }  v( `' v# `& I6 `
if nargin==2
/ H5 ]" c# S5 T, ~" I3 b! \3 S9 d    n1=0;: _. ?, M. i1 [  i7 A1 q8 H
else
5 S& }/ e  F+ j: Z2 B    n1=fs/df;
5 q3 h& M: C1 F8 Q8 Bend8 H: @5 z- M2 Y& r, d
n2=length(m);
0 F: P' u$ @8 y  R$ r3 @1 tn=2^(max(nextpow2(n1),nextpow2(n2)));
- T3 o6 G8 f: Z4 U9 }$ x/ @8 VM=fft(m,n);
- h0 D) }6 C2 R% @  p& d1 im=[m,zeros(1,n-n2)];4 o5 U, a* H5 C+ X! n/ s
df1=fs/n;                                   %以上程序是对调制后的信号u求傅里变换. ~8 k3 O/ g5 b2 M; B& c
M=M/fs;                                   %缩放，便于在频铺图上整体观察
* W4 C6 ], Q; c6 e+ R7 a/ k9 lf=[0:df1:df1*(length(m)-1)]-fs/2;                %时间向量对应的频率向量
. H6 l5 M: l# [+ g! y3 Z$ I7 }' {2 l* Z6 E* q$ l8 o6 L5 f# c
%************对已调信号u求傅里变换**********
* q0 V) N2 r- Cfs=1/ts;
9 m2 ~1 [* a! U" `$ p& yif nargin==2; \" h; G& [4 T) ?" m( j
    n1=0;. U& [; P" a% l# q. Y  p
else! I# y  `/ n7 V& @) e
    n1=fs/df;
# E' u! c* f, G, ]; ~: p7 D- {- o* Oend4 ~, k5 x% M: q% u9 {
n2=length(sfm);' F/ F- |: E. [, e' O
n=2^(max(nextpow2(n1),nextpow2(n2)));2 A' G/ [" b  e- Q4 T
U=fft(sfm,n);
  a, i' W& e7 T" O  H, \% `u=[sfm,zeros(1,n-n2)];
. R# e8 q5 |1 g; \df1=fs/n;                                   %以上是对已调信号u求傅里变换
9 h6 |# Y6 F/ s+ iU=U/fs;                                    %缩放; B3 g3 p3 ]0 X/ [
%******************************************) z$ y! h3 {6 {# y0 Z$ g
%*****************************************" i' \$ Z' D4 q- e  _
%·*·*·*·*·*·*·*·*·*·*·*·*·*·*·*·, z: w; ^3 T! u4 W
%***************显示程序******************% F6 a9 ?' l- o) l* c+ N
disp('按任意键可以看到原调制信号、载波信号和已调信号的曲线')' z! {$ l7 ~" D
pause
, X& R4 w8 t: `! m: s* C%**************figure(1)******************2 F! Y& Z2 H. X: a. U% E
figure(1)
# }0 c6 r- u# E* \subplot(3,1,1);plot(t,mt);                                 %绘制调制信号的时域图
. A4 ~! I, c. ?" _, y2 j' q7 K3 X* |' wxlabel('时间t');1 D" u2 u/ i, o3 a& Y0 e5 J
title('调制信号的时域图');
) D- x1 G5 P8 G* Q5 Y5 K+ Zsubplot(3,1,2);plot(t,ct);                               %绘制载波的时域图
  Q# b+ Q/ C5 P+ Sxlabel('时间t');
2 j9 e1 n2 F2 e; N5 T0 y" s  Utitle('载波的时域图');
- [7 U$ j8 d1 j$ ]subplot(3,1,3);8 w2 x% o9 f1 a- U: z
plot(t,sfm);                                                  %绘制已调信号的时域图8 C9 D& a( F6 U
xlabel('时间t');7 N- y: [% f4 A( d
title('已调信号的时域图');
: ^' ?$ s# g4 j: D: R%******************************************
- P1 J3 @9 r5 P' c4 S# R; Udisp('按任意键可以看到原调制信号和已调信号在频域内的图形')
! x8 z9 F) [' C7 I, U- @2 bpause
3 _8 Q, B( N+ \* n6 _  b4 y5 V%************figure(2)*********************
, W( e( O8 `3 \, H% Z" ?% X" Efigure(2)
3 w* W' I5 B+ z+ ?subplot(2,1,1)
6 l. C5 N, |& S  w% J" Z( K4 x) Cplot(f,abs(fftshift(M)))                        %fftshift:将FFT中的DC分量移到频谱中心, J1 [3 _, J) {% ]+ Z( q2 c4 y
xlabel('频率f')! _( u  q: f. Q+ i# U
title('原调制信号的频谱图')
- Z' n) F& @9 @/ |subplot(2,1,2)
& k% ?5 Q  j2 }% u7 M; `plot(f,abs(fftshift(U)))
/ o9 A: k# p& f  k& W; uxlabel('频率f')3 V7 W+ t- {3 q) ~' l  Y
title('已调信号的频谱图')
# T6 F: U' R7 i; c%******************************************
) B  K6 {( P4 k, ~7 x; {( rdisp('按任意键可以看到原调制信号、无噪声条件下已调信号和解调信号的曲线')) w- a5 ]7 s1 \4 y% B. B
pause; ?+ G# \$ t  N
%**************figure(3)******************& G3 S8 p) Y1 N& @3 D* Y
figure(3)5 Z( j% P' T( Q$ B' O  [
subplot(3,1,1);plot(t,mt);                                  %绘制调制信号的时域图  ~7 s# a# v: q7 Z0 {* r) l
xlabel('时间t');
2 g8 b1 p% T6 W3 a- v% K: jtitle('调制信号的时域图');
' \+ {: e) |- x; ]! g9 Ysubplot(3,1,2);plot(t,sfm);                              %绘制已调信号的时域图( r: a, N& G3 |$ y0 }& n+ }
xlabel('时间t');) g0 o8 h* s: Z( d6 _
title('无噪声条件下已调信号的时域图');+ u" g# f, g; v% `! s  \$ s
nsfm=sfm;                              
) [- w" ]8 v+ Y$ t  u% `for i=1:length(t)-1                           %接受信号通过微分器处理
4 C- s8 D6 ^' ]4 u3 D    diff_nsfm(i)=(nsfm(i+1)-nsfm(i))./dt;
1 f# Z! B& w  T( y9 L( Qend
5 _3 i6 M& \6 X! A$ ddiff_nsfmn = abs(hilbert(diff_nsfm));     %hilbert变换，求绝对值得到瞬时幅度（包络检波）
3 _  w. s) p& v) Czero=(max(diff_nsfmn)-min(diff_nsfmn))/2;
: g8 h7 H! n( s  p: x/ F3 F, ]6 U# j+ pdiff_nsfmn1=diff_nsfmn-zero;+ y! W, M) q* X2 i  M' J
subplot(3,1,3);                             %绘制无噪声条件下解调信号的时域图
! p' Z% u# W  D. o# yplot((1:length(diff_nsfmn1))./1000,diff_nsfmn1./400,'r');
% R3 L) h0 @( q# ~: Fxlabel('时间t');                          9 Y9 z2 Z. u. }0 t. f
title('无噪声条件下解调信号的时域图');! u+ l3 m& K; v9 V" P2 A/ {7 m8 L
%*****************************************
8 i2 c+ z) k) b* Idisp('按任意键可以看到原调制信号、小信噪比高斯白噪声条件下已调信号和解调信号已调信号的曲线')
' R, L9 }6 h# [" Apause3 P# d- z# N  ~, |; e
%**************figure(4)******************
; y7 N% S6 O- w* Q7 Vfigure(4)% P) @0 s$ s. B0 W  e7 `% P
subplot(3,1,1);plot(t,mt);                                  %绘制调制信号的时域图
; W! Q6 W2 C' j' Cxlabel('时间t');
/ [! P/ g( O4 {6 U% J% w9 l5 n! {title('调制信号的时域图');; `5 q4 O+ N9 u
db1=am^2/(2*(10^(sn1/10)));                   %计算对应的小信噪比高斯白躁声的方差
' ~4 B+ N1 G. a2 h. B8 r& Pn1=sqrt(db1)*randn(size(t));                    %生成高斯白躁声. t; N" C3 I/ n7 N/ w/ z
nsfm1=n1+sfm;                              %生成含高斯白躁声的已调信号（信号通
8 b  v* i: G- K! f7 ?: D9 c( I%　　　　　　　　　　　　　　　　　　　　　　%过信道传输）3 d# M0 n" ^. I4 N1 x" i, _
for i=1:length(t)-1                            %接受信号通过微分器处理
2 t+ W2 K) z, h6 O6 [    diff_nsfm1(i)=(nsfm1(i+1)-nsfm1(i))./dt;: h5 I9 l8 _$ C% F3 g0 f
end0 Q: v) }1 D% J! ]( u/ p1 a, Z
diff_nsfmn1 = abs(hilbert(diff_nsfm1));   %hilbert变换，求绝对值得到瞬时幅度（包络检波）
7 L. x- N1 Z5 U2 r5 C$ k; Czero=(max(diff_nsfmn)-min(diff_nsfmn))/2;$ d& i( X/ u! c
diff_nsfmn1=diff_nsfmn1-zero;4 l, F$ g9 V, J# K! x9 u/ {
subplot(3,1,2);- \# h( u& {  w1 f( g5 |6 D  v4 _
plot(1:length(diff_nsfm),diff_nsfm);      %绘制含小信噪比高斯白噪声已调信号的时域图
0 h, f* V# E* y) lxlabel('时间t');# x9 v1 w3 n7 C1 ]
title('含小信噪比高斯白噪声已调信号的时域图');
4 G, N  M$ S3 H! n) k, m  j4 c9 M: Rsubplot(3,1,3);                        %绘制含小信噪比高斯白噪声解调信号的时域图
, L  d# k/ D$ ?' oplot((1:length(diff_nsfmn1))./1000,diff_nsfmn1./400,'r');6 z: @% q- m3 |: q( B5 y; a; @
xlabel('时间t');                          
0 [$ _3 ]2 [7 I  f! _0 k" A0 Qtitle('含小信噪比高斯白噪声解调信号的时域图');
5 k; w! o4 w; N1 p8 U%*****************************************6 ]+ ~9 v5 U: \' ^+ |2 O
disp('按任意键可以看到原调制信号、大信噪比高斯白噪声条件下已调信号和解调信号已调信号的曲线')* n, b( F. T) H3 Q4 T
pause
0 g" m" a7 b% o) t% A3 X%**************figure(5)******************
* n% j& T2 z* yfigure(5)
# ]2 L4 B4 N$ X2 W& B6 ]subplot(3,1,1);plot(t,mt);                                   %绘制调制信号的时域图
1 a' s! k( W+ I. ?  n. Xxlabel('时间t');
4 c8 _; j% u  V: q. ititle('调制信号的时域图');/ ~" w; B' G! M  J4 l
db1=am^2/(2*(10^(sn2/10)));                   %计算对应的大信噪比高斯白躁声的方差; A- z* y6 H; Z# r# n
n1=sqrt(db1)*randn(size(t));                    %生成高斯白躁声1 t2 ~3 L; P8 @" Z1 L- O- C+ @" U
nsfm1=n1+sfm;                              %生成含高斯白躁声的已调信号（信号通过信道传输）7 [5 A3 k$ }# G% ]) O$ |8 b5 i& N
for i=1:length(t)-1                             %接受信号通过微分器处理( Z" r& P2 p  r; g4 B
    diff_nsfm1(i)=(nsfm1(i+1)-nsfm1(i))./dt;# b  a* {: {5 Y& D0 Z7 x
end
3 g( y9 m4 }* B9 Idiff_nsfmn1 = abs(hilbert(diff_nsfm1));            %hilbert变换，求绝对值得到瞬时幅度（包
( t3 }7 w# m. h4 M5 L3 p. v( `%　　　　　　　　　　　　　　　　　　　　　　%络检波）4 U7 ?2 ~) j' T- o3 \
zero=(max(diff_nsfmn)-min(diff_nsfmn))/2;' Z: X' }8 y  H! x" Q
diff_nsfmn1=diff_nsfmn1-zero;. i3 m5 x# }5 j6 e( X! [
subplot(3,1,2);- R# W) n, k; j! B
plot(1:length(diff_nsfm1),diff_nsfm1);             %绘制含大信噪比高斯白噪声已调信号/ G: q  Z, ^8 F$ X1 w9 {7 b
%　　　　　　　　　　　　　　　　　　　　　　%的时域图
" I0 z; ~/ q9 j$ m" T: S# kxlabel('时间t');
3 F- B7 v: K4 x( Z/ V; @7 htitle('含大信噪比高斯白噪声已调信号的时域图');
8 w! R0 i8 g+ n( [& Jsubplot(3,1,3);                                %绘制含大信噪比高斯白噪声解调信号9 z/ z. p$ s6 v: l
%　　　　　　　　　　　　　　　　　　　　　　%的时域图
3 Z6 u/ K2 g/ B+ B2 V. Z% ?plot((1:length(diff_nsfmn1))./1000,diff_nsfmn1./400,'r');
9 A5 H2 y  |7 S- U4 Z; ]5 D0 G2 Uxlabel('时间t');                            ! P0 G2 v9 j- A7 r$ V6 a
title('含大信噪比高斯白噪声解调信号的时域图');
7 A  j# ]; W+ U0 |- ~%*****************************************
- M0 E/ P/ _: u%******************结 束*******************
$ L7 x- g- B  ]! W& g5 M! Z. X6 z0 e

! ~$ z( H+ Q  q9 T
作者: 李秋锐    时间: 2015-2-1 19:23

: z1 q2 `2 z4 l3 b7 k: j, k8 x! Z7 Z) \7 Q3 ?7 x

7 b/ m& d/ M4 Z% z, S) i/ Z- g  Q。。。。。。。。。。。。。。。
9 a* j: K5 [$ W" m6 |
$ b: O6 p7 u& y7 z; ]" \



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