数学建模社区-数学中国

标题: FM调制解调系统 [打印本页]
作者: 李秋锐    时间: 2015-2-1 19:22
标题: FM调制解调系统
%FM调制解调系统.m
% r+ n* x7 L& s/ I%频率调制与解调的Matlab演示源程序+ K, M4 p3 Y! t9 l% W( t
%可以任意改原调制信号函数m(t)
/ D$ H4 r! H0 C% n5 q! e4 }%电子信息工程  杜怀超 091308305  ' C- d- P# D# j8 {! f; U
%·*·*·*·*·*·*·*·*·*·*·*·*·*·*·*·
/ S+ c0 ~0 u" U% h  i) G  s$ S%*****************初始化******************
2 s& v. N- n  ]7 v* `6 H* _echo   off      
/ A2 y. N" B) r2 p' K! Dclose  all
1 ?" c. o  F* X2 E7 eclear  all( r5 B& k  t$ M9 z, l/ x9 a5 H
clc
. H0 K6 y5 y0 W2 X* W" t%*****************************************
% ^' U2 p$ [! P%·*·*·*·*·*·*·*·*·*·*·*·*·*·*·*·
3 x/ u2 G9 F) x  R0 q7 ^! D%****************FM调制*******************
& X6 k7 e* W* p1 O2 w/ S4 zdt=0.001;                                   %设定时间步长/ [; L* B  z$ v
t=0:dt:1.5;                                       %产生时间向量6 a$ t8 t2 m$ ]1 U" Q, D9 Z# p; U
am=5;                                      %设定调制信号幅度* s3 [0 z. C) u5 a/ \, N2 l" ~
fm=5;                                      %设定调制信号频率) ~/ d# F. ]* s0 i8 X, V$ j
mt=am*cos(2*pi*fm*t);                        %生成调制信号
8 v! @+ x; s4 Efc=50;                                      %设定载波频率4 ^' X) Y9 B6 l; w3 B
ct=cos(2*pi*fc*t);                             %生成载波1 h. j9 I% i: u5 p
kf=10;                                      %设定调频指数
5 N3 g) D0 I7 \  Z, aint_mt(1)=0;* T5 n8 \1 ^9 x* g/ G
for i=1:length(t)-1  - l" {& G1 b" ?1 s
    int_mt(i+1)=int_mt(i)+mt(i)*dt;               %求信号m(t)的积分
. S' T& [; r8 @+ @* zend                                          %调制，产生已调信号
$ v3 N; Y/ C5 z- L3 Osfm=am*cos(2*pi*fc*t+2*pi*kf*int_mt);            %调制信号
& d$ I2 Q* T- b5 `, y: [%*****************************************
& P* n! ?' w4 S2 r4 J& c%·*·*·*·*·*·*·*·*·*·*·*·*·*·*·*·2 c% w5 q" K# H9 \) n& h6 W8 K; d
%*************添加高斯白噪声**************
4 b* _7 q; d$ @3 ~# ]sn1=10;                                     %设定信躁比(小信噪比): x1 F5 d+ S& C/ u3 _  }
sn2=30;                                     %设定信躁比(大信噪比)
5 j1 f4 r4 x( o! P  ~% I. P7 x8 L0 Csn=0;                                       %设定信躁比(无信噪比)
+ }" g$ c, t4 o; x+ }2 w- I# udb=am^2/(2*(10^(sn/10)));                     %计算对应的高斯白躁声的方差9 k% t* q$ p. l7 y- H5 c
n=sqrt(db)*randn(size(t));                      %生成高斯白躁声
4 e% A; d! p5 ?  X* nnsfm=n+sfm;                                %生成含高斯白躁声的已调信号（信号通% N; s9 a: m! w; e# d$ T# _- G
%过信道传输）
4 U0 d, z+ W0 S* k: R0 _' E3 M%*****************************************
) y  t1 ]9 i* Y. d%·*·*·*·*·*·*·*·*·*·*·*·*·*·*·*·
) U, z$ q7 Q3 r- ^  [3 c; F6 t%****************FM解调*******************, o0 _* P; g* h, F: j, P
for i=1:length(t)-1                             %接受信号通过微分器处理
# d: ]! h, ~  d6 O, F. |* y7 G, l( h    diff_nsfm(i)=(nsfm(i+1)-nsfm(i))./dt;$ j/ }+ r  }& T
end9 r& h: `. O6 \+ Q
diff_nsfmn = abs(hilbert(diff_nsfm));     %hilbert变换，求绝对值得到瞬时幅度（包络检波）
4 |: @! Y7 q6 t" ~0 t) v: V- D. {zero=(max(diff_nsfmn)-min(diff_nsfmn))/2;
% \7 m0 E! x  q1 q( Z7 F8 N4 jdiff_nsfmn1=diff_nsfmn-zero;
1 @1 |: [8 O4 G5 j%*****************************************) y& B- s! I# g
%·*·*·*·*·*·*·*·*·*·*·*·*·*·*·*·
9 g: n) j6 r+ Z# S%**************时域到频域转换**************
1 w5 J0 \3 F9 B% @" ^$ [& ^ts=0.001;                                    %抽样间隔
9 j0 r2 {- h& i) q8 C8 X/ y6 v  Afs=1/ts;                                      %抽样频率
6 |: p$ F6 x7 {6 g$ C) A7 T; kdf=0.25;                                     %所需的频率分辨率，用在求傅里叶变换
9 V0 W, W, d. }6 Y%　　　　　　　　　　　　　　　　　　　　　　时，它表示FFT的最小频率间隔& |% b3 G" C( m( p* C
%*****对调制信号m(t)求傅里叶变换*****
$ ~; s7 g1 J. F  [" _; Vm=am*cos(2*pi*fm*t);                         %原调信号$ c  n4 [0 m/ \( }8 X
fs=1/ts;
: s; ~# c' Y7 |2 u# n" Zif nargin==2- p$ n( m( h  m9 O# M! m$ R* o2 b
    n1=0;
# |* i( L; C' f. |* b4 Ielse
4 |, J2 P8 a( P  |+ X7 _    n1=fs/df;. f) w5 T9 y' f# Z, W
end
/ F7 z6 U) C( W# D4 d4 ~1 \; jn2=length(m);
+ f( d+ d* O  H% \n=2^(max(nextpow2(n1),nextpow2(n2)));8 Y# R# o4 N9 F: u
M=fft(m,n);- I) e1 _/ v; i5 L
m=[m,zeros(1,n-n2)];! _/ X5 T8 H- {4 L2 t# s0 m& K' L
df1=fs/n;                                   %以上程序是对调制后的信号u求傅里变换
' l: z8 \! p: `$ iM=M/fs;                                   %缩放，便于在频铺图上整体观察$ z8 X  V& Y# t) {  Q
f=[0:df1:df1*(length(m)-1)]-fs/2;                %时间向量对应的频率向量( I' p$ D0 v8 E, y& H
  @" n( @% L: s7 `. E
%************对已调信号u求傅里变换**********
$ c# U8 x. t% c; }! [5 `" Ofs=1/ts;
9 Q- S3 p+ f. W: B' Vif nargin==2
9 F2 H0 A$ ~) k& ^) _    n1=0;& `" |8 w2 w, T" M2 Q/ A- R
else
) `$ F! d/ g: J8 s    n1=fs/df;
; g4 E- J$ J4 Q7 Tend' O. r: N5 ^. Q, E0 N" O
n2=length(sfm);5 X+ y) B* a- K8 [; `
n=2^(max(nextpow2(n1),nextpow2(n2)));
5 n: A. n, X3 L$ hU=fft(sfm,n);4 k) e! g: U1 J- d/ }
u=[sfm,zeros(1,n-n2)];
2 B' g1 j' |# D; e% \: vdf1=fs/n;                                   %以上是对已调信号u求傅里变换
7 V' x6 E7 f( i% S' B/ ~, F, {U=U/fs;                                    %缩放
4 K) |! _2 v3 x( b2 N%******************************************( C9 M% `8 @/ f4 n
%*****************************************. d' _& _0 @# ?: I
%·*·*·*·*·*·*·*·*·*·*·*·*·*·*·*·; u9 j& `  a9 O. G% i
%***************显示程序******************
7 q& b( p% M6 i8 e9 j. P8 X$ `4 jdisp('按任意键可以看到原调制信号、载波信号和已调信号的曲线')
9 H3 v% ]. Q9 v; Ppause2 h( Y  g1 b, u0 h8 Y, i, o# c# H- u
%**************figure(1)******************
/ o1 A; D, X2 }$ ffigure(1)
, }6 }6 I: {4 H# B1 t! o3 H3 _subplot(3,1,1);plot(t,mt);                                 %绘制调制信号的时域图
/ o3 l6 S6 Z+ K# E0 b# jxlabel('时间t');1 L2 \+ R8 G+ R2 h
title('调制信号的时域图');5 _5 o5 M2 D. a! j. l
subplot(3,1,2);plot(t,ct);                               %绘制载波的时域图
* d! |* I& i9 M' y$ fxlabel('时间t');
3 |" [" T# w" E) Z$ Btitle('载波的时域图');6 g$ Q2 {* A( N; s- ?
subplot(3,1,3);
) n% d1 [) t8 u8 E* H' Iplot(t,sfm);                                                  %绘制已调信号的时域图
9 K- {, r7 d1 Y0 Vxlabel('时间t');
1 |* H$ g! A" O  i2 j- f" mtitle('已调信号的时域图');
7 y1 K0 e( K3 B( }1 T( \7 ]%******************************************
0 x+ ]: S& J/ G6 p1 s( O5 C3 idisp('按任意键可以看到原调制信号和已调信号在频域内的图形')
2 `- i5 U# ^4 fpause' o2 y; n8 q! E2 U; _- ~
%************figure(2)*********************$ X$ `5 m- {% z+ ?0 \6 {
figure(2)8 n: ]6 ~) i7 p
subplot(2,1,1)( h0 e( z; v' ?3 t1 x# T) F
plot(f,abs(fftshift(M)))                        %fftshift:将FFT中的DC分量移到频谱中心
6 f8 G+ ^* U! K, Z- k& S2 e& qxlabel('频率f')( l& Q" h/ g- k
title('原调制信号的频谱图')
" v% ^' V  B. L- Z8 Z4 psubplot(2,1,2)
: M9 T/ T" `( C( v6 E$ mplot(f,abs(fftshift(U)))/ G/ j/ l* j9 C9 [) C- G
xlabel('频率f')
+ \1 v$ K! S! E4 L# dtitle('已调信号的频谱图')4 K9 G2 v+ q! G3 m7 K
%******************************************
8 j/ v* D6 Y, {+ y9 Q1 s/ odisp('按任意键可以看到原调制信号、无噪声条件下已调信号和解调信号的曲线')
3 ]: s- f* g: t& k% hpause
$ ^7 ]* A: K$ ?3 p! @%**************figure(3)******************
0 Z+ s) {2 E1 @" Wfigure(3)9 R0 i; H% N" o
subplot(3,1,1);plot(t,mt);                                  %绘制调制信号的时域图
; _% n: s: x4 O' Hxlabel('时间t');6 |$ S8 l* u' T+ q. G$ u; x
title('调制信号的时域图');
( u5 `! v0 f- u  osubplot(3,1,2);plot(t,sfm);                              %绘制已调信号的时域图
3 }. L' t% L# e( H. X! W) z( uxlabel('时间t');
1 |1 a  [* |1 g' ^( a7 \# _title('无噪声条件下已调信号的时域图');. |; b% Y2 W: a' [( n
nsfm=sfm;                              
7 N) @8 W; w# G! n. _for i=1:length(t)-1                           %接受信号通过微分器处理
7 \( R1 \9 u: P! z; m9 V! D3 n% E    diff_nsfm(i)=(nsfm(i+1)-nsfm(i))./dt;
, y1 M. ^& v$ P; R8 u9 \- @end
# v+ l/ }: ^# `& U+ K  j$ w% |diff_nsfmn = abs(hilbert(diff_nsfm));     %hilbert变换，求绝对值得到瞬时幅度（包络检波）
( U- E; `' g5 @$ o% v/ }+ D& Azero=(max(diff_nsfmn)-min(diff_nsfmn))/2;1 O2 t( _' a4 W9 a; A
diff_nsfmn1=diff_nsfmn-zero;
+ C; y; L5 R- S% csubplot(3,1,3);                             %绘制无噪声条件下解调信号的时域图
' x. g0 g# S1 E) F; y3 `+ @9 c. rplot((1:length(diff_nsfmn1))./1000,diff_nsfmn1./400,'r');; `" x0 ^/ ^8 j. @8 L& p
xlabel('时间t');                          & Z3 r! s7 z* e, }- \$ g. N
title('无噪声条件下解调信号的时域图');+ c" i8 X# D9 F2 c1 f
%*****************************************
4 w7 q$ u3 z0 Y" ?' C7 wdisp('按任意键可以看到原调制信号、小信噪比高斯白噪声条件下已调信号和解调信号已调信号的曲线')3 S. x- o6 ^/ i9 L9 d
pause  t, _; _, E, W( x; g! k. j
%**************figure(4)******************0 X1 Y3 \4 k( }8 i
figure(4)
. r- c5 `( e, ~. hsubplot(3,1,1);plot(t,mt);                                  %绘制调制信号的时域图+ s8 m4 w& n/ ^0 L$ t1 G: Q
xlabel('时间t');
# n4 r4 O8 N# M$ ]/ X- Qtitle('调制信号的时域图');
8 t, B' i, x) K9 i( h- Cdb1=am^2/(2*(10^(sn1/10)));                   %计算对应的小信噪比高斯白躁声的方差6 Y) A' L; G/ w0 w
n1=sqrt(db1)*randn(size(t));                    %生成高斯白躁声
0 P& |$ H+ g" Tnsfm1=n1+sfm;                              %生成含高斯白躁声的已调信号（信号通
0 o) W+ f. [" C0 S%　　　　　　　　　　　　　　　　　　　　　　%过信道传输）/ A$ {- l7 X8 K
for i=1:length(t)-1                            %接受信号通过微分器处理" M7 R& ?- Y+ U; ]
    diff_nsfm1(i)=(nsfm1(i+1)-nsfm1(i))./dt;  n0 A  r( s7 U! t
end
  h1 }6 b/ E+ y7 Q" h1 }4 [2 C! j& Pdiff_nsfmn1 = abs(hilbert(diff_nsfm1));   %hilbert变换，求绝对值得到瞬时幅度（包络检波）. u9 p1 Q# D8 p% v4 e+ \  ~
zero=(max(diff_nsfmn)-min(diff_nsfmn))/2;
$ A) M0 ]3 J$ N9 V# }diff_nsfmn1=diff_nsfmn1-zero;
( K. z! J6 Q/ m3 a4 Gsubplot(3,1,2);
6 ]2 X( A! \! h/ B- K  S4 A9 J+ Tplot(1:length(diff_nsfm),diff_nsfm);      %绘制含小信噪比高斯白噪声已调信号的时域图( }2 E8 G# @- N6 {' [) D. t
xlabel('时间t');
" y& B. D3 t/ E9 Utitle('含小信噪比高斯白噪声已调信号的时域图');
' f) T4 N, @, }8 Xsubplot(3,1,3);                        %绘制含小信噪比高斯白噪声解调信号的时域图' k& B3 }( O  u  u/ g$ a8 [+ X
plot((1:length(diff_nsfmn1))./1000,diff_nsfmn1./400,'r');# i3 Z( X* C2 J/ M1 P
xlabel('时间t');                          2 i2 I8 ?# l. \+ n. M
title('含小信噪比高斯白噪声解调信号的时域图');
3 W: y& W) L7 ^" s- h8 C# J) z%*****************************************
) Y/ h% m& `, [! Ldisp('按任意键可以看到原调制信号、大信噪比高斯白噪声条件下已调信号和解调信号已调信号的曲线')# y8 J/ w0 p' ~/ ]7 t
pause
! L* u) b6 X- H. z%**************figure(5)******************
% v6 `! P  [) Z& _( xfigure(5)
0 F( ]9 z& k4 Tsubplot(3,1,1);plot(t,mt);                                   %绘制调制信号的时域图
- ]; {1 x: K) P, j. bxlabel('时间t');
: x( M8 X& v* y6 H, Ytitle('调制信号的时域图');; ^( {5 K0 J- b1 |& ]7 X5 R
db1=am^2/(2*(10^(sn2/10)));                   %计算对应的大信噪比高斯白躁声的方差' W' {( B* T8 ?  c, u
n1=sqrt(db1)*randn(size(t));                    %生成高斯白躁声$ S, q! y0 J/ R' n) P5 q& Q/ U5 }
nsfm1=n1+sfm;                              %生成含高斯白躁声的已调信号（信号通过信道传输）) r1 P* ?0 {+ G: w4 [6 y
for i=1:length(t)-1                             %接受信号通过微分器处理
% l, [7 A! p' z. [3 D    diff_nsfm1(i)=(nsfm1(i+1)-nsfm1(i))./dt;
4 N- s0 X3 y0 ^! v0 R0 L; ^, uend7 G8 r2 f& ?+ D/ R( d
diff_nsfmn1 = abs(hilbert(diff_nsfm1));            %hilbert变换，求绝对值得到瞬时幅度（包
- J! S6 B; N" {. c' |( o%　　　　　　　　　　　　　　　　　　　　　　%络检波）3 P+ O; b2 b! e0 B( ~
zero=(max(diff_nsfmn)-min(diff_nsfmn))/2;2 x  N! `! A& G3 S: u! B: d/ i
diff_nsfmn1=diff_nsfmn1-zero;
1 u' N+ Q1 l  q) C+ o$ [2 ?subplot(3,1,2);
. w3 `1 _6 f+ _7 H7 \plot(1:length(diff_nsfm1),diff_nsfm1);             %绘制含大信噪比高斯白噪声已调信号: w. ^& j6 u) s; ~8 w
%　　　　　　　　　　　　　　　　　　　　　　%的时域图' d, E6 D0 _7 D. T& U
xlabel('时间t');4 e' F3 O+ |) Y& K3 o( Z
title('含大信噪比高斯白噪声已调信号的时域图');& e9 y5 k; I% p' Q
subplot(3,1,3);                                %绘制含大信噪比高斯白噪声解调信号- Y: j9 `0 F  D' X$ l0 b
%　　　　　　　　　　　　　　　　　　　　　　%的时域图
! R: f0 t7 R5 K1 J! V# {* N, c. @. Zplot((1:length(diff_nsfmn1))./1000,diff_nsfmn1./400,'r');$ D9 i  }% O( v, D, N( \# e
xlabel('时间t');                            ' h& Z! K' E5 `0 H: O* ~
title('含大信噪比高斯白噪声解调信号的时域图');
! C5 H$ ?- {/ G7 A%*****************************************
, l' m1 n3 i4 O* o2 c%******************结 束*******************& t( e! ~5 I+ z% c6 E5 s

% }% o& V$ `8 c% }/ A, z; u
! f1 `- q4 }; v3 ]
作者: 李秋锐    时间: 2015-2-1 19:23
  N  P1 `  p$ K

- ~/ d2 j3 e# D$ n( y2 ~$ f. S$ s9 l% z/ Q2 y
。。。。。。。。。。。。。。。
! [6 b$ S: j% J1 R
" c+ E! h9 N, U* H" u9 h8 y) v



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