信号分析与处理实验报告Word文件下载.docx
《信号分析与处理实验报告Word文件下载.docx》由会员分享,可在线阅读,更多相关《信号分析与处理实验报告Word文件下载.docx(20页珍藏版)》请在冰豆网上搜索。
%%
clearall;
closeall;
clc;
%输入数据%
A=input('
输入x(n)序列'
'
s'
);
A=str2num(A);
%A=[1,2,-1,4];
%测试数据%
%校验序列,%
n=length(A);
m=log2(n);
if(fix(m)~=m)
disp('
输入序列长度错误,请重新输入!
'
A=input('
A=str2num(A);
else
输入正确,请运行下一步'
)
end
%码位倒置%
fork=0:
n-1
forj=1:
m%取M位的二进制数%
x1(j)=bitget(k,j);
%倒取出二进制数%
end
x1=num2str(x1);
%将数字序列转化为字符串%
y(k+1)=bin2dec(x1);
%二进制序列转化为十进制数%
clearx1
fork=1:
n
B(k)=A(y(k)+1);
%时间抽取序列%
clearA
%计算%
forL=1:
m%分解为M级进行运算%
LE=2^L;
%第L级群间隔为2^L%
LE1=2^(L-1);
%第L级中共有2^(L-1)个Wn乘数,进行运算蝶运算的两数序号相隔LE1%
W=1;
W1=exp(-1i*pi/LE1);
forR=1:
LE1%针对第R个Wn系数进行一轮蝶运算,共进行LE1次%
forP=R:
LE:
n%每个蝶的大小为LE%
Q=P+LE1;
T=B(Q)*W;
B(Q)=B(P)-T;
B(P)=B(P)+T;
W=W*W1;
B%输出X(k)%
验证结果:
例6-4
b.IFFT
输入X(k)序列'
%A=[6,2+2i,-6,2-2i];
B=conj(B);
%取共轭%
B=B/n%输出x(n)%
六、实验心得与结论
本次实验借助于Matlab软件,我避开了用C平台进行复杂的复数运算,在一定程度上简化了程序,并添加了简单的检错代码,码位倒置我通过查阅资料,使用了一些函数,涉及到十-二进制转换,数字-文本转换,二-文本转换,相对较复杂,蝶运算我参考了书上了流程图,做些许改动就能直接实现。
通过本次实验,我对FFT的过程更加熟悉了,尤其是WN系数的特性,可以通过累乘来实现W系数的更新。
IIR数字滤波器的设计实验
一、实验目的
1.掌握实现模拟信号数字化的编程方法。
2.掌握实现时间抽取快速傅里叶变换(FFT)的编程方法。
3.掌握设计IIR数字滤波器的双线性变换法。
4.加深对信号分析与处理的理解。
1.生成三个不同频率的正弦信号
、
,将它们叠加成一个信号
。
用满足采样定理的采样频率对模拟信号
进行采样,得到离散时间序列
2.编制子程序实现FFT运算,并对
进行FFT运算,观察它的频谱。
3.设计满足已知的设计指标的IIR数字滤波器,并用它对
进行滤波处理。
对滤波后得到的波形再做FFT运算,观察滤波后的信号频谱。
三、实验结果
1.生成的信号
图1原始信号及合成信号的时域波形
图2原始信号的幅度谱
2.滤波器的性能检验
图3模拟滤波器的频率特性
图4数字滤波器的频率特性
3.滤波后的信号
图5滤波前后信号的时域波形对比
图6滤波前后信号的幅度谱对比
附录
附录一:
FFT子程序
function[A]=fftlzr(A)
N=length(A);
M=log2(N);
j=complex(0,1);
iffix(M)~=M
A=[A(:
zeros(2^(fix(M)+1)-N,1)];
I=0;
J=0;
L=1;
whileI~=N-1
ifI<
J
t=A(J+1);
A(J+1)=A(I+1);
A(I+1)=t;
K=N/2;
whileK<
=J
J=J-K;
K=K/2;
J=J+K;
I=I+1;
whileL<
=M
R=1;
LE1=LE/2;
W1=exp(-j*pi/LE1);
whileR<
=LE1
P=R;
whileP<
=N
T=A(Q)*W;
A(Q)=A(P)-T;
A(P)=A(P)+T;
P=P+LE;
R=R+1;
L=L+1;
附录二主程序
clearall
s_f=15;
s_xy=15;
s_t=20;
s_mar=10;
w_f=2.5;
w_axi=2;
%信号产生
T0=8e-2;
Ts=1/800;
N=T0/Ts;
t=0:
Ts:
(N-1)*Ts;
f1=50;
f2=100;
f3=200;
fp=80;
fs=120;
Sig=4*sin(f1*2*pi*t)+4*sin(f2*2*pi*t)+4*sin(f3*2*pi*t);
s_need=4*sin(f1*2*pi*t);
s_1=4*sin(f2*2*pi*t);
s_2=4*sin(f3*2*pi*t);
P_qian=fftlzr(Sig);
figure
(1)
h1=stem(2*pi/(N*Ts).*(0:
N-1),abs(P_qian));
gridon
hx1=xlabel('
模拟角频率/rad/s'
ht1=title('
滤波前幅度频谱'
ha1=gca;
axistight
%滤波器设计
Fs=1/Ts;
wp=2/Ts*tan(fp*2*pi*Ts/2);
ws=2/Ts*tan(fs*2*pi*Ts/2);
rp=1;
rs=15;
[Nl,Wn]=buttord(wp,ws,rp,rs,'
%求阶数及3dB截止频率
[z,p,k]=buttap(Nl);
%求极零点增益
[b,a]=zp2tf(z,p,k);
%系统传递函数
[b2,a2]=lp2lp(b,a,Wn);
%转换为所需低通滤波器
[HS,WS]=freqs(b2,a2,linspace(0,2000,512));
%求0到2000rad/s的频率响应
[hs,wss]=freqs(b2,a2,[wp,ws,[f1,f2,f3].*2*pi]);
%求通带阻带边缘角频率及原始信号频率响应
%画模拟滤波器频谱并标注
figure(6)
subplot(211)
hsf=plot(WS,20*log10(abs(HS)));
holdon
stem(wss,20*log10(abs(hs)),'
filled'
r:
linewidth'
w_f)
text(wss
(1)-100,20*log10(abs(hs
(1)))-10,...
{'
\Omega_p=160\pirad/s'
;
[num2str(20*log10(abs(hs
(1)))),'
dB'
]});
text(wss
(2)-65,20*log10(abs(hs
(2)))-12,...
\Omega_s=240\pirad/s'
[num2str(20*log10(abs(hs
(2)))),'
text(wss(3)-100,20*log10(abs(hs(3)))-10,...
\Omega_1=100\pirad/s'
[num2str(20*log10(abs(hs(3)))),'
text(wss(4)-65,20*log10(abs(hs(4)))-12,...
\Omega_2=200\pirad/s'
[num2str(20*log10(abs(hs(4)))),'
text(wss(5)-70,20*log10(abs(hs(5)))-10,...
\Omega_3=400\pirad/s'
[num2str(20*log10(abs(hs(5)))),'
hxsf=xlabel('
hysf=ylabel('
幅值/dB'
htsf=title('
幅频特性'
hasf=gca;
subplot(212)
hsx=plot(WS,angle(HS));
stem(wss,angle(hs),'
text(wss
(1),angle(hs
(1))+0.7,...
[num2str(angle(hs
(1))),'
rad'
text(wss
(2),angle(hs
(2))+0.7,...
[num2str(angle(hs
(2))),'
text(wss(3)-150,angle(hs(3))-0.8,...
[num2str(angle(hs(3))),'
text(wss(4)-20,angle(hs(4))+0.8,...
[num2str(angle(hs(4))),'
text(wss(5)-50,angle(hs(5))-1.3,...
[num2str(angle(hs(5))),'
hxsx=xlabel('
hysx=ylabel('
相值/rad'
htsx=title('
相频特性'
hasx=gca;
[b3,a3]=bilinear(b2,a2,Fs);
%双线性变换到数字滤波器
[H,W]=freqz(b3,a3,linspace(0,0.5*pi,512));
%求数字滤波器频率响应
[hd,wd]=freqz(b3,a3,[wp,ws,[f1,f2,f3].*2*pi]*Ts);
%求通带阻带边缘数字角频率及原始信号数字频率响应
wd=wd/pi;
%画数字滤波器频谱并标注
figure
(2)
h2=plot(W/pi,20*log10(abs(H)));
stem(wd,20*log10(abs(hd)),'
text(wd
(1)-0.0035,20*log10(abs(hd
(1)))-12,...
{['
\omega_p='
num2str(wd
(1)),'
\pirad'
];
[num2str(20*log10(abs(hd
(1)))),'
text(wd
(2)-0.003,20*log10(abs(hd
(2)))-12,...
\omega_s='
num2str(wd
(2)),'
[num2str(20*log10(abs(hd
(2)))),'
text(wd(3)-0.0035,20*log10(abs(hd(3)))-12,...
\omega_1='
num2str(wd(3)),'
[num2str(20*log10(abs(hd(3)))),'
text(wd(4)-0.003,20*log10(abs(hd(4)))-12,...
\omega_2='
num2str(wd(4)),'
[num2str(20*log10(abs(hd(4)))),'
text(wd(5)-0.003,20*log10(abs(hd(5)))-12,...
\omega_3='
num2str(wd(5)),'
[num2str(20*log10(abs(hd(5)))),'
hx2=xlabel('
归一化数字角频率/\pirad'
hy1=ylabel('
ht2=title('
ha2=gca;
h3=plot(W/pi,angle(H));
grid;
title('
stem(wd,angle(hd),'
text(wd
(1),angle(hd
(1))+1,...
[num2str(angle(hd
(1))),'
text(wd
(2)-0.002,angle(hd
(2))+1,...
[num2str(angle(hd
(2))),'
text(wd(3)-0.003,angle(hd(3))-1,...
[num2str(angle(hd(3))),'
text(wd(4)-0.003,angle(hd(4))-2,...
[num2str(angle(hd(4))),'
text(wd(5)-0.003,angle(hd(5))-1.3,...
[num2str(angle(hd(5))),'
hx3=xlabel('
hy2=ylabel('
相值/\pirad'
ht3=title('
ha3=gca;
%滤波
data=filter(b3,a3,Sig);
figure(3)
h4=plot(t,data,'
k'
t,Sig,'
b'
t,s_need,'
r'
ht4=title('
滤波效果'
xlabel('
t/s'
fontsize'
s_xy);
ylabel('
幅值'
legend('
滤波后波形'
滤波前波形'
50Hz信号'
ha4=gca;
P_hou=fftlzr(data);
figure(4)
h5=stem(2*pi/(N*Ts).*(0:
hx4=xlabel('
ht5=title('
axis([0,1.5e3,0,2e3])
text(f1*2*pi-20,abs(P_qian(f1*N*Ts+1))+200,[num2str(f1)'
Hz'
]);
text(f2*2*pi-20,abs(P_qian(f2*N*Ts+1))+200,[num2str(f2)'
text(f3*2*pi-20,abs(P_qian(f3*N*Ts+1))+200,[num2str(f3)'
ha5=gca;
h6=stem(2*pi/(N*Ts).*(0:
N-1),abs(P_hou));
hx5=xlabel('
ht6=title('
滤波后幅度频谱'
text(f1*2*pi-20,abs(P_hou(f1*N*Ts+1))+200,[num2str(f1)'
text(f2*2*pi-20,abs(P_hou(f2*N*Ts+1))+200,[num2str(f2)'
text(f3*2*pi-20,abs(P_hou(f3*N*Ts+1))+200,[num2str(f3)'
ha6=gca;
figure(5)
h7=plot(t,Sig,t,s_need,'
t,s_1,'
g'
t,s_2,'
c'
hx6=xlabel('
hy3=xlabel('
ht7=title('
生成信号'
set(h7,'
2)
合成信号'
50Hz'
100Hz'
200Hz'
ha7=gca;
%图像处理
h=[h1,h2,h3,h4'
h5,h6,hsf,hsx];
ht=[ht1,ht2,ht3,ht4,ht5,ht6,ht7,htsf,htsx];
hx=[hx1,hx2,hx3,hx4,hx5,hx6,hxsf,hxsx];
hy=[hy1,hy2,hy3,hysf,hysx];
ha=[ha1,ha2,ha3,ha4,ha5,ha6,ha7,hasf,hasx];
set(h([2:
6,9:
10]),'
Linewidth'
set(ha,'
w_axi,'
Fontsize'
s_f);
set([hxhyht],'
s_xy,'
FontWeight'
bold'
set((1:
6),'
outerposition'
get(0,'
screensize'
),'
color'
w'