数字集成系统ii作业

1、仿真一：证明所有样本函数的变换相等。变换的平方的均值与自相关函数的%验证样本函数 clc;clear; x=zeros(1000,1000); a1=zeros(1,1999); y=;%求 1000 个样本函数for k=1:1000变换平方和均值的频谱等于自相关函数的频谱变换x(k,:)=randn(1, 1000);y(k,:)=fft(x(k,:);end%求 1000 个样本函数for k=1:1000if k=1y1=abs(y(k,:).2;else变换的平方和y1=y1+abs(y(k,:).2;endend y=y1;y=y/1000;%即 1000 个样本函数k=fix(l

2、ength(y)/2); fs=1000;f=fs*(0:k-1)/length(y);变换平方和的均值%画出 1000 个样本函数figure(1) subplot(2,1,1)plot(f,abs(y(1:k)/length(y)变换平方和均值的频谱title(样本函数变换平方和均值的频谱)xlabel(Hz)%求 1000 个样本函数各自的自相关函数FFT for k=1:1000a,b=xcorr(x(k,:),unbiased); a1=a1+a;yr=fft(a,length(y);%对自相关函数做if k=1c=yr;elsec=c+yr;endend yr=c/1000;变换%

3、求 1000 个样本自相关函数FFT 的平均值，作为平稳序列的自相关函数的FFTa1=a1/1000;%求 1000 个样本自相关函数的平均值，作为平稳序列的自相关函数%画出平稳序列自相关函数的频谱 kr=fix(length(yr)/2);f=fs*(0:kr-1)/length(yr); subplot(2,1,2) plot(f,abs(yr(1:kr)title(自相关函数的频谱)xlabel(Hz)%画出自相关函数图形 figure(2)plot(b,a1)title(自相关函数)%验证样本函数变换平方和均值的频谱等于自相关函数的频谱figure(3) k=1000;kr=lengt

4、h(y);error=abs(y(1:k)/length(y)-abs(yr(1:kr); plot(error)title(error)%输出最大误差和平均误差 display(最大绝对误差max_error：) max_error=max(abs(error)display(平均误差mean_error：) mean_error=mean(error)display(误差的方差var_error：)var_error=var(error)输出结果：最大绝对误差max_error： max_error =0.1289平均误差mean_error：mean_error =-0.0015误差的方

5、差var_error：var_error=0.0020仿真二：证明所有样本函数的 z 变换的平方的均值与自相关函数的z 变换相等。%根据z 变换的定义计算 z 变换，即对x(n)z(-n)求和 clc;clear;%初始化z1=zeros(1,1000);%存随机序列z 变换 z2=zeros(1,1000);%存自相关函数z 变换 a1=zeros(1,1999);y=zeros(1,1000);%随机序列z 变换平方的yy=zeros(1,1000);%求 1000 个样本函数Z 变换for k=1:10 x=randn(1,1000);%生成随机序列 n=1;a,b=xcorr(x,un

6、biased);%求自相关函数 a1=a1+a;for z=linspace(1.1,2.1,1000) for i=1:1000z1(n)=z1(n)+x(i)*z(-i);%对随机序列z 变换z2(n)=z2(n)+a(i)*z(-i);%对自相关函数z 变换end n=n+1;end y=y+z1.2; yy=yy+z2;endy=y/10;%求 10 个样本 z 变换平方的平均值figure(1) subplot(2,1,1) plot(y)title(样本z 变换平方的均值)yy=yy/10;%求 10 个样本自相关函数 z 变换的平均值，作为平稳序列自相关函数的z 变换subplo

7、t(2,1,2) plot(abs(yy)title(样本自相关函数的z 变换)figure(2) plot(b,a1/10)title(样本的自相关函数)输出结果：输出结果误差比较大，以下是选取误差较小的图仿真三：证明cross-spectra 等于各自变换的乘积%验证cross-spectra 等于各自clc;clear; x=zeros(10,1000); y=zeros(10,1000); xx=zeros(1,1000); yy=zeros(1,1000); for k=1:10 x(k,:)=randn(1,1000);y(k,:)=randn(1,1000); x+x(k,:);

8、yy=yy+y(k,:);end%x=randn(1,1000);%y=randn(1,1000);/10;y=yy/10; X=fft(x);Y=fft(y); XY=X.*Y;a,b=xcorr(x,y,unbiased);变换的乘积r=fft(a,length(XY);%对自相关函数做变换figure(1) k=fix(length(XY)/2); fs=1000;f=fs*(0:k-1)/length(XY); subplot(2,1,1) plot(f,abs(XY(1:k)/length(XY)title(各自变换的乘积)kr=fix(length(r)/2); f=fs*(0:k

9、r-1)/length(r); subplot(2,1,2) plot(f,abs(r(1:kr) title(cross-spectra)%验证cross-spectra 等于各自 figure(2) k=length(XY);kr=length(r);error=abs(XY(1:k)/length(XY)-abs(r(1:kr); plot(error)title(error)变换的乘积%输出最大误差和平均误差 display(最大绝对误差max_error：) max_error=max(abs(error)display(平均误差mean_error：) mean_error=mea

10、n(error)display(误差的方差var_error：) var_error=var(error)输出结果：最大绝对误差max_error： max_error =0.5811平均误差mean_error： mean_error =-0.1137误差的方差var_error： var_error =0.0143仿真四：白噪声、平滑滤波%白噪声、平滑滤波（中值滤波） clc;clear;N=10000;x=randn(1,N);%产生均值为 0，方差为 1 的figure(1)随机序列，也就是白噪声序列plot(x);%白噪声title(白噪声)n=N/10;a,b=xcorr(x,n,

11、biased); %求自相关函数 figure(2)plot(b,a)%画出自相关函数 title(自相关函数)yr1=fft(a);%求自相关函数 y=abs(yr1);%求功率谱 figure(3)plot(y);title(功率谱密度) axis(0 2*n -5 5)变换%中值滤波 yr2=y;d=6;k=length(yr2)-d; for i=1:kxr=yr2(i:i+d); yr2(i)=median(xr);end figure(4)plot(yr2(1:length(yr2)-d)%画出中值滤波后的功率谱密度axis(0 2*n -5 5)title(中值滤波后的功率谱密度)display(中值滤波前的均值和方