1、信號與系統(tǒng)實驗指導(dǎo)手冊山東大學(xué)控制科學(xué)與工程學(xué)院2005-1-26目錄實驗一信號的產(chǎn)生與運算2實驗二離散線性不變系統(tǒng)分析8實驗三周期信號的傅立葉級數(shù)表示11實驗四離散時間傅立葉變換12實驗五連續(xù)時間系統(tǒng)分析16實驗六采樣與重建17實驗一信號的產(chǎn)生與運算-、 信號的產(chǎn)生1.內(nèi)容：產(chǎn)生并畫出以下信號：a>單位沖擊函數(shù)b、單位階躍函數(shù)c、正弦波d、周期二角波和鋸齒波e、周期方波程序:a>subplot (2, 3, 1);%將圖片分為2行3列6部分，在第一部分顯示stem(0, 1);%在零點畫離散信號，幅度為一titlec單位沖激信號')b、subplot (2, 3, 2);

2、x=-10:0.01:10;%x范圍為-10到10,步長為0.01y=(x>=0);plot (x, y);%按上述程序要求畫出連續(xù)信號朋標(biāo)為(x,y)axis(-10, 10,0,2);%顯示坐標(biāo)范圍為橫坐標(biāo)to至10,縱坐標(biāo)0至2titlec單位階躍信號')c> subplot (2, 3, 3); x=-pi:pi/20:pi; plot(x, sin(x); title('正弦波')d、subplot (2, 3, 4); x=-2*pi:0. 01:2*pi;plot (x, sawtooth(x, 0. 5); titlec周期三角波')

3、e、subplot (2, 3, 5); x=-2*pi:0. 01:2*pi; plot (x, saw to oth(x, 1); titlec周期鋸齒波')f、subplot (2, 3, 6);x二-3*pi:0. 01:3*pi; plot (x, square(x); axis(-10, 10,0,2);titlec周期方波')圖像:10.80.60.40.20j單位沖激信號-1 0 1單位階躍信號正眩波10.50-0.5-1周期三角波10.50-0.5-1周期鋸齒波-10 0 10-10 0 10周期方波2. 內(nèi)容：產(chǎn)生并畫出以下離散序列：a> 離散止余弦序

4、列：xn) = 3cos(0.1n + 龍/3) + 2sin(0.57?),0 < n <101/7 = 3b單位釆樣序列= 0/工3f 1 ,2 > 3c、單位階躍序列u(n-3) = '-0,n<3d、實數(shù)指數(shù)序列x(n) = (0.9)0<n<10程序：functionx=impseq(n0, nl, n2)%產(chǎn)生單位采樣序列n=nl:n2;x= (n-no) =0;functionx, n=stepseq (no, nl, n2)%產(chǎn)生單位階躍序列 n=nl:n2;x= (n-no)>=0;a. subplot (2, 2, 1);n

5、=0:10;x=3*cos(0* pi*n+pi/3)+2*sin(0.5*pi*n); stem(n,x)grid ontitlec離散余弦序列')b. subplot (2, 2, 2); x=impseq(3,0j0); stem(n,x); grid on;titlec單位采樣序列')c. subplot (2, 2, 3);x= stepseq(3,0j0); stem(n,x);grid on;titlec單位階躍序列')d. subplot (2, 2, 4);n=0:10;x=(0.9).an;stem(n,x)grid ontitlec實數(shù)指數(shù)序列&#

6、39;)圖像:-%離散余弦序列510p)v9()(<)<)i<>0-2-42單位釆樣序列單位階躍序列實數(shù)指數(shù)序列4內(nèi)容：產(chǎn)生復(fù)數(shù)值信號:x(ri) = e+jq3)n -10 < n < 10衣四個子圖中畫出其幅度、相位、實部和虛部的波形。 程序：>> n=-10:1:10;alpha=-0. 1+0. 3j;» x=exp(alpha*n);» subplot (2, 2,1);stem(n, real (x) ; title(，實部');xlabel c n');» subplot (2, 2,

7、2); stem(n, imag(x) ; title('虛部');xlabel (' n');>> subplot (2, 2, 3); stem(n, abs (x); title c 振幅');xlabel (' n');» subplot(2, 2, 4); stem(n, (180/pi) *angle(x) :titlec 相位');xlabel c n*); 圖像：實部)06(ar> 虛部-1050510n 振幅n相位二、對序列的運算1.相加：在matlab中可川算術(shù)運算符“ + ”實現(xiàn)。

8、然而坷)和x2(n)的長度必須相等。我們必 須首先給殲)和x2(h)以適當(dāng)?shù)膮?shù)使他們有相同的位置向量斤。這需要注意matlab的 下標(biāo)運算。特別是邏輯與運算，關(guān)系運算如“s”和“=”，以及要用find函數(shù)來使xl(n) 和x2(n)具有相同的長度。下而稱為sidadd的函數(shù)演示了這些運算。function(y, n)=sigadd(xl, nl, x2, n2)%實現(xiàn) y(n)=xl(n)+x2(n)%y, n=sigadd(xl, nl, x2, n2)%y=在包括nl和n2的n上求和序列%n= (min (mi n (nl), min (n2) :max (max (nl), max (

9、n2) ;%y (n)的長度 yl=zeros(l, length(n) ;y2=yl;%初始化yl (find(n>=min(nl)&(n<=max(nl)=l) )=xl;%具有 y (n)的長度 xly2 (f ind(n>=min(n2) )&(n<=max (n2)=l) )=x2;%具有 y (n)的長度 x2y二yl+y2;%序列相加2. 相乘：這是對對應(yīng)采樣值z間的相乘，它rh數(shù)組運算符實現(xiàn)。但它受到“ + ”運算符同 樣的限制。因此必須建立一個與sigadd相仿的sigmult函數(shù)：function(y, n)=sigmul t (xl

10、, nl, x2, n2)%實現(xiàn) y(n)=xl(n)*x2(n)%y, n=sigadd(xl, nl, x2, n2)%y=在包括nl和n2的n上求和序列%n=min(min (nl), min(n2): max (max (nl), max (n2) ;%y (n)的長度 yl=zeros(l, length(n) ;y2=yl;%初始化yl (find(n>=min(nl)&(n<=max (nl) =1) =xl;%具有 y(n)的長度 xly2(find(n>=min(n2)&(n<=max(n2)=l)=x2;%具有 y (n)的長度 x2

11、y二yl. *y2;%序列相乘3. 移位：在這個運算屮，兀)的每一個樣木都移動k個周期，移位后的序列y(n) oy(n) = x(n - k)若令m = n-k,因而n = m + k,上述運算化成：y(m + £) = x(m)實現(xiàn)這一運算的函數(shù)如下：function(y, n)二sigshift (x, m, no)%實現(xiàn) y (n) =xl (n-no)%n=m+no;y=x;4. 折疊：在這個運算中，x(n)的每個樣本都對« = 0翻轉(zhuǎn)，得到一個折疊后的序列y(n)。y(z1) = x(-n)在matlab中，這一運算對采樣值由fliplr (x)函數(shù)來實現(xiàn)，而對采

12、樣位置則由一fliplr (n)得到，這可從函數(shù)sigfold中看到。function(y, n)二sigfold(x, n)%實現(xiàn) y (n) =x (-n)%y=fliplr(x);n=-fliplr(n);5信號能壘：序列x(h)的能量山下式給出：0000ex =工卜-00一8有限長序列的能量可用以下的matlab命令求得：>>ex=sum(x. *conj(x) ;%種方法>>ex=sum(abs (x). 2) ;%另一-種方法內(nèi)容：畫出以下各序列在給定區(qū)間的波形圖a. x(n) = 2§n + t)-s(n-4),-5 < n < 5b

13、. x(n) = nu(n) -u(n 10) + 10e3(n-10) un -10) un - 20),0 < n < 20c. x(n) = cos(0.04n) + 0.2vv(/?),0 < n < 50 其中w(n)為具有零均值及單位方差的高斯隨機序列。d. x(n) = ,5,4,3,2,1,5,4,3,2,1,5,4,3,2,1 ,一 10 < n <9 程序：functionx =impseq (no, n 1, n2)%產(chǎn)生單位采樣序列函數(shù)n=nl: n2;x=(n-no)=0;functionx, n二stepseq (no, nl,

14、n2)%產(chǎn)生單位階躍序列函數(shù)n=nl:n2;x=(n-no)>二0;a. n二-5:5;>> x二2*impseq(-2, -5, 5)-impseq(4, -5, 5);>> subplot (2, 2, 1);>> stem(n, x);» titlec (1)的序列')；>> xlabel (' n ) ;ylabel (' x (n)');b. n二0:20;>> xl=n. * (stepseq (0, 0, 20) -stepseq (10, 0, 20);>>

15、x2=10*exp (-0. 3*(nt0). * (stepseq (10, 0, 20)-stepseq (20, 0, 20);>> x=xl+x2;>> subplot(2, 2, 2): stem(n, x); title c (2)的序列)；c. n二0:50;>> x二cos(0. 04*pi*n)+0. 2*randn(size(n);>> subplot(2, 2, 3) ;stem(n, x); titlec* (3)的序列');d. n =-10:9 ;x=5, 4, 3, 2, 1;>> xtilde=

16、x，*ones (1, 4);>> xtilde=(xtilde(:)* ;>> subplot (2, 2, 4): stem(n, xtilde); title c (4)的序列');圖像:(1)的序列(2)的序列(3)的序列(4)的序列實驗二離散線性不變系統(tǒng)分析1.內(nèi)容：設(shè)(斤)=(0.9)n w(n),輸入 x(n) = u(n)-u(n -10),求系統(tǒng)輸出 y(n) = h(n) * x(n)程序：function u, n二impulse(no, nl, n2) ; %從 nl 至 n2,在 no 處產(chǎn)生一個脈沖 n=nl:n2;u二(n-n0)=

17、0;function x, n=step(no, nl, n2) ;%從 nl 至 n2,在 no 處以后產(chǎn)牛階躍信號n=nl:n2;x=(n-n0)>=0;function y, ny=conv_m(x, nx, h, nh) ;%信號處理的改卷積程序nyb=nx(l)+nh(l) ;nye=nx(length(x)+ nh仃ength(h);ny二nyb:nye;y=conv (x, h);n二-5:50;ul二stepseq (0, -5, 50);u2=stepseq(10, -5, 50); x二ul-u2;%表示維數(shù)相同的向量對應(yīng)位置相乘，*指矩陣相乘h=(0. 9). b)

18、. *ul; figure(l);subplot (3, 1, 1);stem(n, x,'.');axis(-5, 50, 0, 2);ylabel (' xn'); subplot (3, 1, 2); stem (n, h,'.'); axis (-5, 50, 0, 2); ylabel (' h n');y, ny=conv_m(x, n, h, n); subplot (3, 1,3); stem (ny, y,'.'); ylabel (，yn');圖像：2.內(nèi)容:設(shè)離散系統(tǒng)對由下列差分方程表

19、示y(n) - y(n -1) + 0.9y(n - 2) = x(n)a、計算n=-20,100上的沖激響應(yīng)b、計算n=-20100上的階越響應(yīng)程序：functionx=impseq (no, nl, n2)%產(chǎn)生單位采樣序列 n二nl :n2;x=(n-no)=0;functionx, n二stepseq(no, nl, n2)%產(chǎn)生單位階躍序列 n=nl:n2;x=(n-no) >=0;a=l, -1, 0. 9;b=l;x=impseq (0, -20, 120);n=-20:120; h=filter(b, a, x);%差分方程響應(yīng)端系數(shù)矩陣為a,輸入端為b,輸入信號為 x,

20、求響應(yīng)subplot (2, 1, 1); stem (n, h, ' ');xl abel (' n'); ylabel (，hnj，);x=stepseq(0, -20, 120); s=f i 1 ter (b, a, x);subplot (2, 1, 2); stem (n, s, t );xlabel (' n');ylabel (' sn');圖象：n實驗三周期信號的傅立葉級數(shù)表示內(nèi)容：川傅立葉級數(shù)合成方波圖形步驟:方波是一個奇諧信號,rh傅里葉級數(shù)可知，它是由無窮個奇次諧波分量合成的，現(xiàn)在我們用圖形的方式來表示它的

21、合成。我們只考察從t二0到t=10s這段時間內(nèi)的信號。 讓我們先畫出方波奇諧信號。t = 0:. 1:10; y = sin(t) ; s二(y>0)*0. 91;subplot (2, 3, 1);plot (t, s), xlabel (* 時間(秒)')ylabel (* 方波波形'),axis(0 10 0 1), grid然后畫出基波分量。subplot (2, 3, 2);plot (t, y), xlabel c 時間(秒)')，ylabel ('波形')，grid現(xiàn)在，將三次諧波加到基波z上，并畫出結(jié)果，如圖所示。subplot (

22、2, 3, 3);y = sin(t) + sin(3*t)/3;plot (t, y), xlabel ('時間(秒)，ylabel (* 波形')，grid再將一次、三次、五次、七次和九次諧波加在一起。subplot (2, 3, 4);y = sin (t) + sin(3*t)/3 + sin(5*t)/5 + sin(7*t)/7 + sin (9*t)/9;plot (t, y), xlabel (* 時間(秒),ylabel (* 波形')，grid最后，我們合并從基波到十九次諧波的各奇次諧波分量。subplot (2, 3, 5);t = 0:. 05:

23、10; y二0;for k=l:2:19y = y + sin(k*t)/k;end plot (t, y), xlabel ('時間(秒)'),ylabel ('波形') grid將上述波形分別畫在一幅圖中，可以看出它們逼近方波的過程。注意“吉布斯現(xiàn)象”。 周期信號傅里葉級數(shù)在信號的連續(xù)點收于該信號，在不連續(xù)點收斂于信號左右極限的平均 值。如果我們用，周期信號傅里葉級數(shù)的部分和來近似周期信號，在不連續(xù)點附近將，會出 現(xiàn)起伏和超量。在實際中，如果應(yīng)用這種近似，就應(yīng)該選擇足夠大的n,以保證這些起伏擁 有的能量可以忽略。subplot (2, 3, 6);x = z

24、eros (size (t);t = 0:. 02:3. 14;y = zeros (10, max (size (t); for k二1:2:19x 二 x + sin(k*t)/k; y(k+l)/2, :) = x;endplot (t, y(l:2:9, :),), title(，方波的組成:吉布斯現(xiàn)象') xlabel c 時間(秒)'),ylabelc 波形')；set (gca,' xliin , 0 3. 142) grid圖像：-1o10.50-0.5510時間（秒）0.5-0.50510時間（秒）0-110-0.50.50510時間（秒）方波

25、的組成:吉布斯現(xiàn)象實驗四離散時間傅立葉變換1內(nèi)容：求以下有限吋寬序列的付里葉變換x(r°)a. 已知x(n) = (0.9小”門)0<x<10b. 已知兀()=2"105x510程序：nl=0:1:10;x=0. 9*exp(j*pi/3). *nl;k二200:200;w= (pi/100)*k;x=x* (exp (-j*pi/100). 7nf *k); magx=abs(x) ;% 幅值angx=angle(x) ;% 相角figure(l);subplot (2, 2, 1);plot (w/pi, magx);xlabel (' f'

26、);ylabel ('i xi |'); subplot (2, 2, 2); plot (w/pi, angx/pi);xlabel (' f'); ylabel ('xi');n2=-10:l:10;y二2. ”n2;y=y*(exp(-j*pi/100). *(n2 *k); magy=abs(y);angy=anglc(y);figure(l);subplot (2, 2, 3);plot (w/pi, magy); xlabel (' f');ylabel (' x2 ');subplot (2, 2,

27、4);plot (w/pi, angy/pi); xlabelc f');ylabelc <x2f);圖象：os 5 1000 < 77 < 1002.內(nèi)容：設(shè)有兩序歹ij x(n) = cos(n/2) y(n) = e/,/r/4x(n) 試比較 x(ej),y(e) 程序：n=0:100;x=cos(pi*n/2);k=-100:100;w二pi/100*k;x=x*(exp(-j*pi/100)(n'*k); magx=abs(x);angx=angle(x); subplot (2, 2, 1);plot (w/pi, magx);xlabol (&

28、#39; frequcnce*pi);ylabel ('|x|'); subplot (2, 2, 2);plot (w/pi, angx);xlabelc frenquence*pi); ylabel ('x');m二0:100;y=exp(j*pi*m/4) *x; l=-100:100;wl二pi/100*l;y=y*(exp(-j*pi/100).magy=abs(y);angy=angle(y); subplot (2, 2, 3);plot (wl/pi, magy);xlabol (' frenqucnce*pi"); ylabe

29、l ('i y|');subplot (2, 2, 4);plot (wl/pi, angy);xlabel (，frenquency*pi'); ylabel c <y");圖像:604020-10°21210-1-2mill imiuhi inni in-2-10 1ox023.內(nèi)容：頭序列x()= sin(加/2)-5 << 10求出x(ej(0)的實部和虛部，同時求出奇、偶部相應(yīng)的xe(ej(0xo(ejm) 程序：n=-5:10;x=sin(pi*n/2);k=-100:100;w二pi/100*k;x=x* (exp (

30、-j*pi/100). * (n，*k);subplot (2, 1, 1);plot (w*pi, x);xlabelc frequency * pi');yl abel (' real & x e');subplot (2, 1, 2);plot (w*pi, -j*x);xlabelc frequency * pi');ylabelc imagry & x o,)圖象:-2-10-8-6-4-20246810frequency * pi-1086420246810frequency * pi5 0 50xo3-1010實驗五連續(xù)時間系統(tǒng)分析

31、1.內(nèi)容：設(shè)有一個穩(wěn)定的lt1系統(tǒng)可以由下面的微分方程描述a. dy(t)/dt + 3y(/) = 3x(z)b. 3d2y(t)/dr + 4dy(t) / dt + y(t) = d2x(t) / dt2 +5r 程序：%afigure(l);a=3;b二3,1;freqs (b, a); %此函數(shù)計算連續(xù)時間系統(tǒng)的頻率響應(yīng)%bfigure (2);c=5, 0, 1:d二1,4, 3; freqs (d, c);圖像： pnlcbe乏io410210°10'210'210-11010101010100010frequency (rad/s)80604020-2

32、10-10 1 210 10 10103104frequency (rad/s)1 010 102 o o o pnlcbe 乏10-210frequency (rad/s)-15050050-100-2o1010 10frequency (rad/s)10實驗六采樣與重建1 內(nèi)容：設(shè)xa二凹"a. 求其付里葉變換xa(ejf0)b. 分別用頻率 fs=5000hz、fs=1000hz 對 xa(f)進(jìn)行取樣，給出 xej(0x2(ei<0) 曲線c. 車建信號xal(txa2(t),并且與xq作比較說明取樣頻率對信號的重建影響。 程序：)t二0. 00005;t=-0. 00

33、5: dt: 0. 0045;xa=exp (t000*abs (t);%設(shè)置連續(xù)時間 fs二5000hz, ts二0. 2msts 二0.0002;n二25:1:25;x=exp(-1000*abs(n*ts);%做fourier變換k=500;k二0:l:k;w二pi*k/k;x二x*exp(-j*n' *w);x二real (x);w=-fliplr(w), w(2:k+l);x=fliplr(x),x(2:k+l);figure(2);subplot (2, 1, 1);plot (t*1000, xa);xlabel (' t inmsec.');ylabel

34、 (' xl (t)');titlec* discrete signal');hold on;stem(n*ts*1000, x);hold off;axis(-5 5 0 1);subplot (2, 1,2);plot (w/pi, x);xlabelfreq in pirads units');ylabel (' xi (w)');title(，discrete-time fourier transform'); %gtext (' ts=0. 2msec');%di serete-time si gnalts=0.

35、001;n=-25:1:25;%fs=1000hz, ts=lms x=exp(-1000*abs(n*ts);%diserete-time fourier transformk二500;k二0:l:k;w二pi*k/k;x二x*exp (-j*n，*w);x=real(x);w=-fliplr(w),w(2:k+l);x=fliplr(x),x(2:k+l);figure (3);subplot (2, 1, 1);plot (t*1000, xa);xlabel (' t inmsec');ylabel (' xl (t)');t itle('di

36、screte signal');hold on;stem(n*ts*l000, x);hold off;axis (-5 5 0 1);subplot (2, 1, 2);plot (w/pi, x);xlabel (* freq in pi rads units');ylabel (' xl (w)');title(，discrete-time fourier transform'); %gtext (' ts=lmsec，);dt=o.00005;t=-0. 005:dt:0. 005;xa=exp(1000*abs(t);wmax=2*pi*2000;k二500;k二0:l:k;w二k*wmax/k;xa=xa*exp(-j*t' *w) *dt;%xa=real(xa);w=-flipl