数字信号处理合成地震记录fft.docx
- 文档编号:2280695
- 上传时间:2022-10-28
- 格式:DOCX
- 页数:10
- 大小:74.60KB
数字信号处理合成地震记录fft.docx
《数字信号处理合成地震记录fft.docx》由会员分享,可在线阅读,更多相关《数字信号处理合成地震记录fft.docx(10页珍藏版)》请在冰豆网上搜索。
数字信号处理合成地震记录fft
数字信号处理实验报告
实验一、地震子波波形显示及一维地震记录合成
一、实验目的
1、认识地震子波(以雷克子波为例),对子波有直观的认识。
2、利用线性褶积公式合成一维地震记录。
二、实验内容
1、雷克子波:
(零相位子波)、
(最小相位子波),
其中代表子波的中心频率,代表子波宽度,随着的增大,子波能量后移,当=7时,最小相位子波可视为混合相位子波,这里取=25Hz,=3;
2、根据公式编程实现零相位子波、最小相位子波的波形显示;
3、设计反射系数(n=500),其中,,,,,其它为0;
4、应用褶积公式合成一维地震记录,并图形显示;
5、根据所学知识对实验结果进行分析。
三、实验结果:
1、零相位子波:
(1)程序源代码:
%编写零相位子波
t=0.002;
r=3;
fm=25;
forn=1:
51
w(n)=exp(-(2*pi*fm/r)^2*(t*n)^2)*cos(2*pi*fm*t*n);
end
plot(w)
xlabel('n')
ylabel('w')
title('零相位子波')
(2)图像:
2、最小相位子波:
(1)程序源代码:
%最小相位子波
t=0.002;
r=3;
fm=25;
forn=1:
51
w(n)=exp(-(2*pi*fm/r)^2*(t*n)^2)*sin(2*pi*fm*t*n);
end
plot(w)
xlabel('n');ylabel('w');
title('最小相位子波')
(2)图像:
3、对比不同时的波形图
(1)程序:
t=0.002;
r=3;
fm=25;
forn=1:
51
w1(n)=exp(-(2*pi*fm/r)^2*(t*n)^2)*cos(2*pi*fm*t*n);
end
r=4;
forn=1:
51
w2(n)=exp(-(2*pi*fm/r)^2*(t*n)^2)*cos(2*pi*fm*t*n);
end
r=5;
forn=1:
51
w3(n)=exp(-(2*pi*fm/r)^2*(t*n)^2)*cos(2*pi*fm*t*n);
end
subplot(1,3,1),plot(w1);
axis([0,55,-0.7,1]);xlabel('n');title('r=3时');
subplot(1,3,2),plot(w2);
axis([0,55,-0.7,1]);xlabel('n');title('r=4时');
subplot(1,3,3),plot(w3);
axis([0,55,-0.7,1]);xlabel('n');title('r=5时');
(2)图像:
(3)分析:
代表子波宽度,随着的增大,子波能量后移。
4、一维地震记录:
(1)零相位子波程序:
t=0.002;
r=3;
fm=25;
forn=1:
51
w(n)=exp(-(2*pi*fm/r)^2*(t*n)^2)*cos(2*pi*fm*t*n);
end
%设置反射系数
r=zeros(500);
r(100)=1.0;
r(200)=-0.7;
r(300)=0.5;
r(400)=0.4;
r(500)=0.6;
%编写褶积公式
f=zeros(1,550);
forn=1:
550
form=1:
500
if(1<=(n-m)&&(n-m)<=51)
f(n)=f(n)+r(m)*w(n-m);
end
end
end
plot(f)
(2)零相位子波图像:
(3)最小相位子波程序:
t=0.002;
r=3;
fm=25;
forn=1:
51
w(n)=exp(-(2*pi*fm/r)^2*(t*n)^2)*sin(2*pi*fm*t*n);
end
%设置反射系数
r=zeros(1,500);
r(100)=1.0;
r(200)=-0.7;
r(300)=0.5;
r(400)=0.4;
r(500)=0.6;
%编写褶积公式
f=zeros(1,550);
forn=1:
550
form=1:
500
if(1<=(n-m)&&(n-m)<=51)
f(n)=f(n)+r(m)*w(n-m);
end
end
end
plot(f)
(4)最小相位子波图像:
(5)对比零相位子波和最小相位子波的一维地震记录:
放大如下图:
局部放大可发现,最小相位子波比零相位子波的地震记录要滞后。
最小相位子波的能量要稍小于零相位子波的能量。
程序:
t=0.002;
r=3;
fm=25;
forn=1:
51
w1(n)=exp(-(2*pi*fm/r)^2*(t*n)^2)*cos(2*pi*fm*t*n);
w2(n)=exp(-(2*pi*fm/r)^2*(t*n)^2)*sin(2*pi*fm*t*n);
end
r=zeros(1,500);
r(100)=1.0;
r(200)=-0.7;
r(300)=0.5;
r(400)=0.4;
r(500)=0.6;
f1=zeros(1,550);
f2=zeros(1,550);
forn=1:
550
form=1:
500
if(0<(n-m)&&(n-m)<=51)
f1(n)=f1(n)+r(m)*w1(n-m);
f2(n)=f2(n)+r(m)*w2(n-m);
end
end
end
plot(f1)
hold
plot(f2,'r')
title('对比最小相位和零相位子波的一维地震记录')
四、结果分析:
1、离散雷克子波,取采样间隔一般是2毫秒或者1毫秒,在这里我取了2毫秒。
采样点数取51个。
2、代表子波宽度,随着的增大,子波能量后移。
3、褶积公式:
中,r的长度为M=500,w的长度为N=51,则f的长度为N+M-1=550。
所以计算过程中应该用二重循环,外层是n从1到550,内层是m从1到500.同时循环过程要满足1<=n-m<=51。
实验二、带通滤波及频谱分析
一、实验内容
对复合频率信号进行频谱分析,并根据其振幅谱设计带陷滤波器,滤掉某些频率成分。
二、实验步骤
1、设计某一信号,包含多种频率成分(可用雷克子波);
2、将离散,并应用fft进行频谱分析,绘出振幅谱;
3、分析振幅谱有什么特点,在频率域设计带陷滤波器(可加斜坡),以消除某频段(大于62.5Hz)的频率成分,并显示滤波后的振幅谱。
(要求绘出滤波器图形)
4、将滤波后的信号反变换回时间域,并绘出信号曲线,观察其与原信号的差别。
5、根据所学知识对实验结果进行分析。
三、实验过程
以主频为25Hz的雷克子波为例。
设置的子波长度为200,滤波器长度为51。
(1)不加斜坡的滤波器
分析:
关于N/2=128对称。
但是吉卜斯现象波动明显。
程序如下:
clf
t=0.002;
r=3;
fm=25;
forn=1:
200
x(n)=exp(-(2*pi*fm/r)^2*(t*n)^2)*cos(2*pi*fm*t*n);
end
figure
(1)
plot(x)
title('主频为25Hz的雷克子波')
forn=1:
256
forn=1:
200
x(n)=exp(-(2*pi*fm/r)^2*(t*n)^2)*cos(2*pi*fm*t*n);
end
forn=201:
256
x(n)=0;
end
end
s=fft(x);
real_s=real(s);
imag_s=imag(s);
p=sqrt(real_s.^2+imag_s.^2);
figure
(2)
plot(p)
title('快速傅立叶变换后的振幅谱')
ylabel('振幅')
fx=62.5;%消除大于62.5Hz的成分
N=256;
ff=1/(N*t);
k1=fx/ff;
k2=N-k1;
forn=1:
256
forn=1:
k1
r(n)=1;
end
forn=k1-1:
k2
r(n)=0;
end
forn=k2+1:
N
r(n)=1;
end
xx(n)=r(n)*p(n);
end
figure(3)
plot(r);
title('滤波器图形')
xxx=abs(xx);
figure(4)
plot(xxx);
title('滤波后的振幅谱')
ylabel('振幅')
s_s=ifft(xxx);
figure(5)
plot(real(s_s))
title('滤波后的信号曲线')
(2)加斜坡的滤波器
程序:
除滤波器处的程序不同外,其余相同,不做赘述。
下面附上加斜坡的滤波器的程序如下:
forn=1:
256
forn=1:
k1-5
r(n)=1;
end
forn=k1-4:
k1-1
r(n)=8-k1/4;
end
forn=k1:
k2
r(n)=0;
end
forn=k2:
k2+4
r(n)=k2/4-56;
end
forn=k2+5:
N
r(n)=1;
end
xx(n)=r(n)*p(n);
end
分析,在0和1过渡时没有直接跳跃,而是设计了一个线性函数,即加了一个斜坡过渡。
吉卜斯效应:
当频率特性曲线是不连续函数而对滤波因子去有限项时,导致对应的频率特征曲线是一波动的曲线,频率域滤波算子反生畸变。
(3)对比两种滤波器
(4)对比两种滤波后的信号曲线
(5)对比不同斜率的斜坡的滤波作用:
- 配套讲稿:
如PPT文件的首页显示word图标,表示该PPT已包含配套word讲稿。双击word图标可打开word文档。
- 特殊限制:
部分文档作品中含有的国旗、国徽等图片,仅作为作品整体效果示例展示,禁止商用。设计者仅对作品中独创性部分享有著作权。
- 关 键 词:
- 数字信号 处理 合成 地震 记录 fft