数字信号处理 实习大报告 7.docx
- 文档编号:28637792
- 上传时间:2023-07-19
- 格式:DOCX
- 页数:25
- 大小:524.46KB
数字信号处理 实习大报告 7.docx
《数字信号处理 实习大报告 7.docx》由会员分享,可在线阅读,更多相关《数字信号处理 实习大报告 7.docx(25页珍藏版)》请在冰豆网上搜索。
数字信号处理实习大报告7
一、从给定的程序(文件包Friday.rar)中,选择一个源程序做详细标注。
(目的:
熟悉Matlab程序)
见m文件。
二、能够利用Matlab熟悉地画图,内容包括:
X、Y坐标轴上的label,每幅图上的title,绘画多条曲线时的legend,对图形进行适当的标注等。
(1)在一副图上画出多幅小图;
(2)画出一组二维图形;
(3)画出一组三维图形;(4)画出复数的实部与虚部。
%%%%%%%%%%%%%%
clear
closeall
figure
(1)
x=linspace(0,6);
y1=sin(2*x);
y2=sin(x.^2);
y3=(sin(x)).^2;
subplot(2,2,1);
plot(x,y1);
subplot(2,2,2);
plot(x,y2);
subplot(2,2,3);
plot(x,y3);
subplot(2,2,4);
y4=sin(x);
plot(x,y4);
title('2*2×éͼ')
print-djpeg-r01.jpeg
%%%%%%%%%%%%%%
figure
(2)
t=linspace(-6,6);
f1=t.^2;
f2=-t.^2;
plot(t,f1,'-k');
holdon;
plot(t,f2,':
r');
text(-2,20,'f1=t^2','FontSize',12,'FontWeight','bold');
text(-2,-20,'f2=-t^2','FontSize',12,'FontWeight','bold');
title('¶þάͼ')
print-djpeg-r02.jpeg
%%%%%%%%%%%%%%
figure(3)
[X,Y]=meshgrid(-3:
0.2:
3);
Z=81.4881+1.2877*X+2.9766*Y;
mesh(X,Y,Z)
holdon
plot3(X,Y,Z,'x','MarkerSize',3);
title('Èýάͼ')
print-djpeg-r03.jpeg
%%%%%%%%%%%%%%
figure(4)
b=3+4i;
compass(b);
title('¸´Êýͼ')
print-djpeg-r04.jpeg
三、计算普通褶积与循环褶积,分别使用时间域与频率域两种方法进行正、反演计算,指出循环褶积计算时所存在的边界效应现象;编写一个做相关分析的源程序。
clearall;
closeall;
%普通褶积通过矩阵完成
aa=[2234]';
bb=[13553]';
%%%%%%%%%%%%%%%%
%利用conv函数直接计算,对计算结果进行验证
%%%%%%%%%%%%%%%%
aa_bb=conv(aa,bb);
bb_aa=conv(bb,aa);
%%%%%%%%%%%%
N1=length(aa);
N2=length(bb);
N=N1+N2-1;
A=zeros(N,N2);
B=zeros(N,N1);
forj=1:
N2;
A(j:
j+N1-1,j)=aa;
end
forj=1:
N1;
B(j:
j+N2-1,j)=bb;
end
%%%%%%%%%%%%%
%计算结果
AA=A*bb;
BB=B*aa;
figure
(1);
subplot(3,1,1);stem(aa,'.');title('x(n)');%画图
subplot(3,1,2);stem(bb,'.');title('y(n)');%画图
subplot(3,1,3);stem(AA'.');title('褶积');%画图
print-djpeg-r0secondwork1.jpeg%输出
%%%%%%%%%%%%%%%
%普通褶积通过循环褶积的利用MATLAB中的IFFT和FFT实现
aa1=zeros(N,1);
bb1=zeros(N,1);
aa1(1:
N1)=aa;
bb1(1:
N2)=bb;
AA1=fft(aa1);
BB1=fft(bb1);
CC=AA1.*BB1;
CC1=BB1.*AA1;
cc1=ifft(CC1);
cc=ifft(CC);
%循环褶积通过快速傅立叶变换
w=exp(-sqrt(-1)*2*pi/N);
aa2=zeros(N,1);
bb2=zeros(N,1);
aa2(1:
N1)=aa;
bb2(1:
N2)=bb;
F=zeros(N);
fori=1:
N
forj=1:
N
F(i,j)=power(w,(i-1)*(j-1));
end
end
X1=F*aa2;
Y1=F*bb2;
DD=1/N*conj(F)*(X1.*Y1);
DD1=1/N*conj(F)*(Y1.*X1);
%通过程序写出循环矩阵的系数矩阵,求取循环褶积DD=DD1说明褶积具有可交换性
N=5;
%%%%%%%可以改变N的值来验证N=N1+N2时,循环褶积和线性褶积相等;
C1=zeros(N);
C2=zeros(N);
aa3=zeros(N,1);
bb3=zeros(N,1);
aa3(1:
N1)=aa;
bb3(1:
N2)=bb;%%%%%%%%%%%%%%%%%%%双重循环
forj=1:
N
fori=1:
N
k=i+j-1;
ifk>N
k=k-N;
end
C1(k,j)=aa3(i);
C2(k,j)=bb3(i);
end
end
EE=C1*bb3;
EE1=C2*aa3;
figure
(2);
subplot(3,1,1);stem(aa,'.');title('x(n)');%画图
subplot(3,1,2);stem(bb,'.');title('y(n)');%画图
subplot(3,1,3);stem(EE,'.');title('褶积');%画图
print-djpeg-r0secondwork2.jpeg
边界效应:
两个离散的序列离散x(n)和y(n),他们的长度分别为N1和N2,如果循环褶积的长度N>N1+N2-1,则循环褶积和线性褶积的值相等。
否则,循环褶积和线性褶积的值不相等。
而在循环褶积的系数矩阵中,把第一列离散的序列依次折叠刀第二列,第三列等,这个时候,当循环褶积的N=N1+N2-1时,循环褶积的系数矩阵和线性褶积的系数矩阵不相同,循环褶积系数矩阵边沿是不是零,而线性褶积系数矩阵边沿是零,但是他们运算结果一致,这就是循环褶积在计算过程中的边界效应。
线性相关和循环相关源程序:
clearall;%%%%%%%%%%%%%%%%%%相关分析%%%%%%%%%%%
closeall;
aa=[1,2,3,4]';
bb=[1,2,5,3]';
%%%%%%%%%%%%%%线形相关分析
N1=length(aa);
N2=length(bb);
N=N1+N2-1;
N=4;
A=zeros(N);
aa3=zeros(N,1);
bb3=zeros(N,1);
aa3(1:
N1)=aa;
bb3(1:
N2)=bb;
forj=1:
N;
A(j:
N,N-j+1)=aa3(1:
N-j+1);
end
CC1=A*bb3;
figure
(1)
subplot(3,1,1);stem(aa,'.');title('x(n)');%画图
subplot(3,1,2);stem(bb,'.');title('y(n)');%画图
subplot(3,1,3);stem(CC1,'.');title('线性相关');%画图
print-djpeg-r0secondwork3.jpeg;
%%%%%%%%%%%%%%%循环相关分析
N1=length(aa);
N2=length(bb);
N=N1+N2-1;
N=4;
A1=zeros(N);
aa4=zeros(N,1);
bb4=zeros(N,1);
aa4(1:
N1)=aa;
bb4(1:
N2)=bb;
forj=1:
N;
A1(1:
N-j+1,j)=aa4(j:
N);
A1(N-j+2:
N,j)=aa(1:
j-1);
end
CC2=A1*bb3;
figure
(2)
subplot(3,1,1);stem(aa,'.');title('x(n)');%画图
subplot(3,1,2);stem(bb,'.');title('y(n)');%画图
subplot(3,1,3);stem(CC2,'.');title('循环相关');%画图
print-djpeg-r0secondwork4.jpeg;
四、设计一个病态(矩阵)系统,分析其病态程度;找出对应的解决方法(提示:
添加白噪因子)。
A=[1221]';%定义AB两个列向量
B=[1234]';
N=length(A);%N等于A的长度
F=zeros(N);%定义F、AA、DD都为的全零方阵
AA=zeros(N);
DD=zeros(N);
w=exp(-sqrt(-1)*2*pi/N);%定义w因子
fori=1:
N
forj=1:
N
AA(i,j)=A(mod(i-j,N)+1);%用循环的方式列写A
F(i,j)=power(w,(i-1)*(j-1));
end
end
FA=F*A;%FA为A的频谱
FB=F*B;%FB为B的频谱
C0=AA*B;%C0A与B的褶积
C1=real(conj(F)/N*(FA.*FB));%用F因子的方法算A与B的褶积
C0-C1
FC=F*C0;%FC为C0的频谱
C_FA=conj(FA);%求取FA的共轭
FC2=C_FA.*FC;
FA2=C_FA.*FA+0.001;%添加一个很小的数以保证其频谱大于0
FX=FC2./FA2;
X=conj(F)*FX/N;
BX=inv(AA+0.001*eye(N))*C0;%添加白噪因子
BY=inv(AA)*C0;%无白噪因子
结果如下
BX=%添加白噪因子的结果
1.4996
1.4986
3.4996
3.5006
BY=%没有添加白噪因子,计算不出结果
Inf
Inf
Inf
Inf
五、设计一个一维滤波处理程序(1、分别做低通、高通、带通、带阻等理想滤波器进行处理;2、窗函数)。
%一维理想滤波器设计
closeall
clearall
f1=30;f2=45;%设定滤波器的两个频率,f1为低通滤波器的高通频率、带通滤波器的低通频率,f3为高通滤波器的高截频率、带通滤波器的高通频率
d=0.01;%设定采样间隔
t=0:
d:
1;%建立时间向量
x=sin(2*pi*10*t)+2*sin(2*pi*25*t)+8*cos(2*pi*45*t);%由时间向量构成的含有高、低、中频的原始信号
lx=length(x);%时间向量长度
fc=1/(2*d);%截止频率
f=0:
2*fc/lx:
2*fc;%求出频谱对应的频率序列
FX=fft(x);%求原始信号的频谱
%理想低通滤波器设计
H1=zeros(1,lx);
fori=1:
lx
iff(i)<=f1|f(i)>=(2*fc-f1)
H1(i)=1;
elseH1(i)=0;
end
end
DFX1=FX.*H1;
DX1=ifft(DFX1);
figure
(1);
subplot(2,2,1),plot(H1);axis([0100-0.51.5]);title('低通H1(f)');
subplot(2,2,2),plot(t,x);title('原始信号');
subplot(2,2,3),stem(abs(FX),'.');;title('原始信号频谱');
subplot(2,2,4),stem(abs(DFX1),'.');;title('经低通滤波后的低频信号的频谱');
print-djpeg-r0低通滤波器.jpeg
%理想高通滤波器设计
H2=zeros(1,lx);
fori=1:
lx
iff(i)<=2*f1&&f(i)>=f1
H2(i)=1;
elseH2(i)=0;
end
end
DFX2=FX.*H2;
DX2=ifft(DFX2);
figure
(2);
subplot(2,2,1),plot(H2);axis([0100-0.51.5]);title('高通H2(f)');
subplot(2,2,2),plot(t,x);;title('含有高低中频的原始信号');
subplot(2,2,3),stem(abs(FX),'.');title('原始信号频谱');
subplot(2,2,4),stem(abs(DFX2),'.');title('经高通滤波后的高频信号的频谱');
print-djpeg-r0高通滤波器.jpeg
%理想带通滤波器设计
H3=zeros(1,lx);
fori=1:
lx
iff(i)<=f2&&f(i)>=f1||(f(i)<=(fc+f2-f1)&&(f(i))>=(fc));
H3(i)=1;
elseH3(i)=0;
end
end
DFX3=FX.*H3;
DX3=ifft(DFX3);
figure(3);
subplot(2,2,1),plot(H3);axis([0100-0.51.5]);title('带通H3(f)');
subplot(2,2,2),plot(t,x);title('含有高低中频的原始信号');
subplot(2,2,3),stem(abs(FX),'.');title('原始信号频谱');
subplot(2,2,4),stem(abs(DFX3),'.');title('经带通滤波后的中频信号的频谱');
print-djpeg-r0带通滤波器.jpeg
%理想带阻滤波器设计
f1=20;f2=40;
H4=zeros(1,lx);
fori=1:
lx
if(f(i))<=f2&&(f(i))>=f1||((f(i))<=(f2+fc)&&(f(i))>=(f1+fc))
H4(i)=0;
elseH4(i)=1;
end
end
DFX4=FX.*H4;
DX4=ifft(DFX4);
figure(4);
subplot(2,2,1),plot(H4);axis([0100-0.51.5]);title('带阻H4(f)');
subplot(2,2,2),plot(t,x);title('含有高低中频的原始信号');
subplot(2,2,3),stem(abs(FX),'.');title('原始信号频谱');
subplot(2,2,4),stem(abs(DFX4),'.');title('经带阻滤波后的含有高低频信号的频谱');
print-djpeg-r0带阻滤波器.jpeg
六、设计一个二维滤波处理程序(分别做低通、高通等处理)。
%二维高通滤波#######################
clear;
clc
closeall;
n=200;
e=rand(n,n);
T=peaks(n)+e;%引入的peaks的函数
figure
(1)
size(T);
subplot(221);
mesh(T)
subplot(222)
M=fft(T);
mesh(real(M))
m=1:
70;%设定虑频的宽度
n=110:
200;
M(m,:
)=0;
M(n,:
)=0;
M(:
m)=0;
M(:
n)=0;
subplot(223)
mesh(real(M))
title('二维高通滤波');
subplot(224)
P=ifft(M);
mesh(real(P));
print-djpeg-r0二维高通滤波.jpeg
%二维低通滤波
n=200;
e=rand(n,n);
T=peaks(n)+e;
figure
(2)
subplot(221);
mesh(T)
subplot(222)
M=fft(T);
mesh(real(M))
M(70:
110,:
)=0;
M(:
70:
110)=0;
subplot(223)
mesh(real(M))
title('二维低通滤波');
subplot(224)
H=ifft(M);
mesh(real(H));
print-djpeg-r0二维低通滤波.jpeg
七、验证时间域的循环褶积对应的是频率域的乘积;线性褶积则不然。
closeall;
clearall;
aa=[2234]';
bb=[1355]';
N1=length(aa);
N2=length(bb);
N=N1+N2-1;
%%%%%%%%%%%系数矩阵计算循环褶积%%%%%
aa1=zeros(N,1);
bb1=zeros(N,1);
aa1(1:
N1)=aa;
bb1(1:
N2)=bb;
C1=zeros(N);
C2=zeros(N);
forj=1:
N
fori=1:
N
k=i+j-1;
ifk>N
k=k-N;
end
C1(k,j)=aa1(i);
C2(k,j)=bb1(i);
end
end
EE=C1*bb1
%%%%%%%%利用FFt计算褶积%%%%
aa3=fft(aa1);
bb3=fft(bb1);
M=aa3.*bb3;
EE1=ifft(M)
max=abs(EE1-EE)
%%%%%max=0说明循环褶积时间域和频率域对应%%%
figure
(1);
subplot(4,1,1);stem(aa,'.');title('x(n)');%画图
subplot(4,1,2);stem(bb,'.');title('y(n)');%画图
subplot(4,1,3);stem(EE,'.');title('时间域计算循环褶积');%画图
subplot(4,1,4);stem(EE1,'.');title('频率域计算循环褶积');%画图
print-djpeg-r0sixwork1.jpeg
%%%%%%计算线形褶积%%%%
A=zeros(N,N2);
B=zeros(N,N1);
forj=1:
N2;
A(j:
j+N1-1,j)=aa;
end
forj=1:
N1;
B(j:
j+N2-1,j)=bb;
end
%%%%%%%%%%%%%%%%%%
%计算结果
AA=A*bb;
AA1=fft(aa);
BB1=fft(bb);
CC=AA1.*BB1
AA2=ifft(CC)
figure
(2);
subplot(4,1,1);stem(aa,'.');title('x(n)');%画图
subplot(4,1,2);stem(bb,'.');title('y(n)');%画图
subplot(4,1,3);stem(AA,'.');title('时间域计算普通褶积');%画图
subplot(4,1,4);stem(AA2,'.');title('频率域计算普通褶积');%画图
print-djpeg-r07.jpeg
八、请用通俗、易懂的语言说明数字信号处理中的一种性质、一条定理或一个算例(顺便利用Matlab对其进行实现)。
%信号相加:
信号x(n)=x1(n)+x2(n);值得注意的是,当序列x1(n)
和x2(n)的长度不等或位置不对应时,首先应使两者的位置对齐,然后通过zeros函数左右补0使其长度相等后再相加。
程序如下:
clearall;
n1=0:
3;%设置向量n1
l1=length(n1);
x1=[2,0.50.91];%以x1为纵坐标设置n1对应的值
figure
(1)
stem(n1,x1)
axis([-1802.1])
title('n1')
print-djpeg-r0picture_1.jpeg;
n2=0:
7;%设置向量n2
x2=[00.10.20.30.40.50.60.7];%以x2为纵坐标设置n2对应的值
figure
(2)
stem(n2,x2)
axis([-1800.8])
title('n2')
print-djpeg-r0picture_2jpeg;
n=0:
7;
x1=[x1zeros(1,8-length(n1))];
x2=[zeros(1,8-length(n2)),x2];
x=x1+x2;
figure(3)
stem(n,x)
axis([-1802.1])
title('n3')
print-djpeg-r0picture_3jpeg;
- 配套讲稿:
如PPT文件的首页显示word图标,表示该PPT已包含配套word讲稿。双击word图标可打开word文档。
- 特殊限制:
部分文档作品中含有的国旗、国徽等图片,仅作为作品整体效果示例展示,禁止商用。设计者仅对作品中独创性部分享有著作权。
- 关 键 词:
- 数字信号处理 实习大报告 数字信号 处理 实习 大报