FIR:有限脉冲响应滤波器。有限说明其脉冲响应是有限的。与IIR相比,它具有线性相位、容易设计的优点。这也就说明,IIR滤波器具有相位不线性,不容易设计的缺点。而另一方面,IIR却拥有FIR所不具有的缺点,那就是设计同样参数的滤波器,FIR比IIR需要更多的参数。这也就说明,要增加DSP的计算量。DSP需要更多的计算时间,对DSP的实时性有影响。以下都是低通滤波器的设计。FIR的设计:FIR滤波器的设计比较简单,就是要设计一个数字滤波器去逼近一个理想的低通滤波器。通常这个理想的低通滤波器在频域上是一个矩形窗。根据傅里叶变换我们可以知道,此函数在时域上是一个采样函数。通常此函数的表达式为:
sa(n)=sin(n∩)/n∏,但是这个采样序列是无限的,计算机是无法对它进行计算的。故我们需要对此采样函数进行截断处理。也就是加一个窗函数。就是传说中的加窗。也就是把这个时域采样序列去乘一个窗函数,就把这个无限的时域采样序列截成了有限个序列值。但是加窗后对此采样序列的频域也产生了影响:此时的频域便不在是一个理想的矩形窗,而是成了一个有过渡带,阻带有波动的低通滤波器。通常根据所加的窗函数的不同,对采样信号加窗后,在频域所得的低通滤波器的阻带衰减也不同。通常我们就是根据此阻带衰减去选择一个合适的窗函数。如矩形窗、汉宁窗、汉明窗、BLACKMAN窗、凯撒窗等。选择一个具体的窗函数之后,根据所设计滤波器的参数来计算所需的阶数、此窗函数的表达式。然后用这个窗函数去和采样序列相乘,就可以得到实际滤波器的脉冲响应。IIR的设计(双线性变换法):IIR的设计理念是这样的:根据所要设计滤波器的参数去确定一个模拟滤波器的传输函数,然后再根据这个传输函数,通过双线性变换、或脉冲响应不变法来进行数字滤波器的设计。它的设计比较复杂,复杂在于它的模拟滤波器传输函数H(s)的确定。这一点我们可以让软件来实现。然后,我们说一下它的具体实现步骤:首先你要先确定你需要一个什么样的滤波器,巴特沃斯型,切比雪夫型,还是其它什么型的滤波器。当你选定一个型号后,你就可以根据设计参数和这个滤波器的计算公式来确定其阶数、传输函数的表达式。通常这个过程中还存在预扭曲的问题(这只是双线性变换法所需要注意的问题,脉冲响应不变法不存在这种问题)。确定H(S)后,就可以通过双线性变换得到其数字域的差分方程。
把原观测平面的磁异常换算到某一高度(上下延拓),把实测的磁异常ΔT换算为Xa,Ya,Za三个分量,或求其水平方向和垂直方向的导数,把斜磁化情况下的磁异常化到地磁极等,这些磁异常转换的方法在解释中得到广泛的应用。对实测的磁异常进行转换处理在频率域实现最为方便、快捷。所谓频率域位场转换是把空间域实测的磁异常通过傅里叶变换得到频谱,再乘以换算因子,反复换算来实现。人们习惯把空间域的磁异常通过傅里叶变换来实现各种换算称为频率域或波数域换算。
在频率域进行磁异常的转换,其最大优点是空间域的褶积关系变为频率域的乘积关系;同时还可以把各种换算统一到一个通用表达式中,从而使磁异常的换算变得简单。另一个优点则是可以从频谱特性出发,形象地讨论各种换算的滤波作用。
下面着重讨论三度异常在频率域中的各种换算,实测磁异常为水平观测面。
(一)延拓
延拓是把原观测面的磁异常通过一定的数学方法换算到高于或低于原观测面的平面上,分为向上延拓与向下延拓。向上延拓是一种常用的处理方法,它的主要用途是削弱局部干扰异常,反映深部异常。我们知道,重磁场随距离的衰减速度与具剩余密度和磁性的地质体体积有关。体积大,重磁场衰减慢;体积小,重磁场衰减快。对于同样大小的地质体,重磁场随距离衰减的速度与地质体埋深有关。埋深大,重磁场衰减慢;埋深小,重磁场衰减快。因此小而浅的地质体重磁场比大而深的地质体重磁场随距离衰减要快得多。这样就可以通过向上延拓来压制局部异常的干扰,反映出深部大的地质体。图3-7-1是内蒙古某地用磁法勘探普查超基性岩的实例。该地区浅部盖有一层不厚的玄武岩,使磁场表现为强烈的跳动。为压制玄武岩的干扰,将磁场向上延拓了500m。由图可知,向上延拓的磁场压制了玄武岩的干扰,同时右侧部分反映了深部的超基性岩磁场。图3-7-2用向上延拓压制了浅部矿体的异常,突出了深部盲矿体产生的低缓异常。
通过向上延拓来研究深部磁性基底构造也是其应用的一个重要方面。如把一个地区航磁资料先化极与向上延拓,消除了浅部磁性体影响后再作磁场分区。
与向上延拓相反,向下延拓时随着延拓深度的加大,一些浅的局部干扰或误差也迅速增大,使延拓曲线发生剧烈跳动,甚至出现振荡而无法利用。为了克服这种影响,往往将圆滑和延拓配合使用,在向下延拓之前要对异常进行圆滑。
利用向下延拓可以分离水平叠加异常。我们知道,磁性体埋深越大,异常显得越宽缓。剖面越接近磁性体,磁异常的范围越接近磁性体边界。例如,对两个相邻的板状体而言,当它接近地表时,实测磁异常可能明显地显示两个峰值。但当埋深大于两个板体的距离时,则其叠加异常将显示为两个宽而平的异常。因此,将叠加的磁异常向下延拓到接近磁性体界面时就可能把各个磁性体的异常分离开来,增强分辨能力(见图3-7-3)。
图3-7-1 用向上延拓压制浅部玄武岩异常的影响
1—玄武岩;2—沉积岩
图3-7-2 用向上延拓压制浅部矿体的异常
利用向下延拓还可以评价低缓异常,低缓异常是指强度和梯度都比较小的异常,显然这是磁性体埋藏较深的标志。低缓异常的某些异常特征是不明显的,用它来进行解释推断有一定困难。解决这一困难的办法就是向下延拓。向下延拓一方面可以突出叠加在区域背景上的局部异常,使之尽量少受区域场的影响;另一方面可以“放大”某些在低缓异常中不够明显的异常特征(如拐点、极值点、零值点等),有利于进一步解释推断。
下面我们来推导频率域延拓公式。
1向上延拓
设场源位于z=H平面以下(H>0),则磁场在z=H平面以上是对x、y、z的连续函数,具有一阶和二阶连续可微的导数。若z=0观测平面上的磁场T (x,y,0)为已知,可以得到向上延拓公式为
地球物理勘探概论
由褶积积分公式可知,上式为T(x,y,0)与 关于变量(x,y)的二维褶积。空间域的褶积与频率域的乘积相对应。下面分别求T(x,y,0)及 的傅里叶变换,设T(x,y,z)对于变量(x,y)的傅里叶变换为ST(u,v,z),有:
地球物理勘探概论
图3-7-3 向下延拓分离水平叠加异常
则
地球物理勘探概论
利用上式可以由已知的T(x,y,0)求出其频谱ST(u,v,0)。进一步求 的傅里叶变换,应用Erdelyi(1954)给出的积分变换表可以得到:
地球物理勘探概论
当z<0时上式成立。由式(1-2-4)、式(1-2-5)并用褶积定理有:
地球物理勘探概论
上式对于z≤0成立。
T(x,y,z)是ST(u,v,z)的反傅里叶变换,即
地球物理勘探概论
式(3-7-6)即为向上延拓的频谱表达式。
2向下延拓
我们还可把这一表达式进一步推广到0<z<H的情况。假如我们已知z=h(0<h<H)平面上的场值,则根据向上延拓公式,可求出向上延拓距离为h的平面(即z=0平面)上的场值:
地球物理勘探概论
根据上面同样的推导方法可以求得:
地球物理勘探概论
所以
地球物理勘探概论
上式对于0<h<H成立。这样就把表示式(3-7-1)扩展到场源以上的-∞<z<H。
式(3-7-6)表明,由z=0平面上的磁场值,求出它的傅里叶变换ST(u,v,0),由它乘以延拓因子 ,z>0时向下延拓,z<0时向上延拓),然后通过反傅里叶变换,即可求出z<H空间磁场的表示式。
(二)导数换算
重磁异常的导数可以突出浅而小的地质体的异常特征而压制区域性深部地质因素的影响,在一定程度上可以划分不同深度和大小异常源产生的叠加异常,且导数的次数越高,这种分辨能力就越强(图2-8-8)。
重磁高阶导数可以将几个互相靠近、埋深相差不大的相邻地质因素引起的叠加异常划分开来,如图3-7-4所示。
图3-7-4 两个相邻球体异常的叠加
这些功能主要是因为导数阶次越高,则异常随中心埋深加大而衰减越快。从水平方向来看,基于同样道理,阶次越高的异常范围越小,因而无论从垂向看或从水平方向看,高阶导数异常的分辨能力提高了。
下面我们来推导重磁异常导数换算的公式。如果令Szx(x,y,z)、Szy(x,y,z)、Szz(x,y,z)及Szxx(x,y,z)、Szyy(x,y,z)、Szzz(x,y,z)分别为Za(x,y,z)对x、y、z的一阶导数及二阶导数的频谱,则由微分定理易于得到:
地球物理勘探概论
地球物理勘探概论
地球物理勘探概论
同理,可以写出:
地球物理勘探概论
由此可知,求磁场的n阶垂向导数的频谱,应乘上的导数因子为[2π(u2+v2)1/2]n;而求磁场沿x方向或y方向的n阶水平导数的频谱,应乘上的导数因子为(2πiu)n或(2πiv)n。
如果求磁场的m阶垂向导数、n阶沿x方向水平导数、l阶沿y方向的导数的频谱(即求 的频谱),应乘上的导数因子为
地球物理勘探概论
进一步设l是实测平面上任一方向,它与x轴的夹角为α,则有:
地球物理勘探概论
两边作傅氏变换并应用微分定理,得知:
STl(u,v,z)=i(2πucosα+2πvsinα)ST(u,v,z)
利用上式即可实现磁场的频率域方向导数计算,以此突出某一方向的异常特征。
(三)ΔT换算Hax、Hay、Za各分量
由磁场与磁位的关系可以得到以下磁场各分量之间的关系式:
地球物理勘探概论
上列式中t0为地磁场方向的单位矢量。
设Sx(u,v,z)、Sy(u,v,z)、Sz(u,v,z)以及ST(u,v,z)分别为Hax(x,y,z)、Hay(x,y,z)、Za(x,y,z)及其ΔT(x,y,z)的频谱。
利用频谱微分定理可得到上列场各分量导数在频率域内相应的换算关系式:
地球物理勘探概论
式中 ,而L0、M0、N0为地磁场单位矢量t0的方向余弦。
(四)化到地磁极
把ΔT化到地磁极的过程包含了ΔT化Za的分量换算和斜磁化Za化垂直磁化Za⊥的磁化方向换算。式(3-7-9)中的第三个式子已经实现了ΔT(x,y,z)与Za(x,y,z)频谱之间的换算,下面我们来进一步推导斜磁化Za化到垂直磁化Za⊥的公式。
磁化方向换算的方法是由斜磁化的磁场Za求垂直磁化方向的磁位U⊥,再由垂直磁化磁位U⊥求垂直磁化的磁场Za⊥。
(1)垂直磁方向的磁位U⊥:
地球物理勘探概论
Za1是原斜磁化方向的磁场。上式表明垂直磁化磁位U⊥是Za1沿着原磁化方向t1反方向的曲线积分。
由傅里叶变换可以写出Za1的频谱表达式:
地球物理勘探概论
则
地球物理勘探概论
(2)由垂直磁化磁位U⊥求垂直磁化磁场Za⊥:
地球物理勘探概论
由(3-7-9)式可以看出,由斜磁化Za化为垂直磁化的转换因子是:
地球物理勘探概论
q1是原磁化方向的方向余弦。
因此,由ΔT化到地磁极的转换因子为
地球物理勘探概论
若不考虑剩磁,即地磁场方向的方向余弦q0与磁化方向的方向余弦一致,则上式可进一步化简为
地球物理勘探概论
由于这种转换相当于把ΔT换算到地磁极的地磁场状态,故称为化到地磁极。
(五)几种换算的滤波作用
1磁异常的转换与电滤波的相似性
我们知道,对无线电技术中的滤波来讲,滤波器的输入和输出可以是随时间变化的电压,线性滤波器的输入εi(t)和输出ε0(t)的关系可用一个褶积型的积分方程式表示:
地球物理勘探概论
函数φ(τ)称为权函数,反映了滤波器的特性,也称滤波器的脉冲响应函数。对褶积积分方程进行傅里叶变换就可得到频率域内输入-输出方程:
地球物理勘探概论
式中:E0(ω)、Ei(ω)和Φ(ω)分别为ε0(t)、εi(t)和Φ(t)的傅里叶变换;ω=2πf称为角频率;Φ(ω)为滤波权函数的频谱,也称为滤波器的频率响应函数。
对于磁异常转换,若设转换前后的二度磁异常分别为Z(x)及Zb(x),则可对空间磁异常的转换写出其一般形式:
地球物理勘探概论
这也是一个褶积方程,φ(x)为参加磁异常换算的函数。
如向上延拓,可表示为
地球物理勘探概论
同样,也可写出频率域内磁异常转换的关系式:
地球物理勘探概论
式中:Sb(ω)、S(ω)、Φ(ω)分别为Zb(x)、Z(x)及φ(x)的频谱。
若将磁异常变换在空间域、频率域的表达式与滤波器在时间域、频率域的输入-输出方程式做一对比,可以发现二者是相似的。转换前后的异常对应于输入、输出电压,转换时参加运算的函数φ(x)则对应于滤波器的权函数,这里所不同的仅仅是电压是时间的函数,而磁异常是空间坐标的函数。因此,我们可以把磁异常的转换看成是对异常的滤波,有关滤波的一些基本理论都能用于磁异常的转换。
2几个常用换算的滤波作用
磁异常的换算相当于一种滤波,滤波器的权函数就反映了换算的特性。在频率域中通过研究权函数的频谱(即频率响应函数)来了解磁场换算的滤波作用,下面以二度异常的几种换算的频率响应为例来说明。
图3-7-5画出了向上延拓、向下延拓、垂向一次导数、垂向二次导数的频率响应函数Φ(ω)随角频率ω变化的图形。
图3-7-5 二度异常几种换算的频率响应曲线
由图可见,向上延拓的频率响应:Φ上延(ω)=e-h|ω|。
当ω=0时,Φ(ω)=0;当ω趋于无穷时, 趋于零。它起到了压制高频让低频通过的作用,相当于一个低通滤波器。
向下延拓、垂向一次导数、垂向二次导数的频率响应分别为eh|ω|、ω、ω2,由图形及分析频率响应的表达式容易知道三者都相当于一个高频放大器。比较二者,可以发现向下延拓的频率响应随ω增大,呈指数规律增加。当ω趋于无穷,Φ趋于无穷,它对高频放大作用特别明显,且是不收敛的,故在计算时要特别注意。在运算时常常加一个圆滑因子,以对高频放大作用进行一些抑制。其次,垂向二次导数的频率响应比垂向一次导数在高频放大作用上更强一些,因为它们随ω增大的规律与导数方次有关。导数次数愈高,愈要注意设法压制高频干扰。
MATLAB 信号处理常用函数
一、 波形产生
函数名 功能
sawtooth 产生锯齿波或三角波Sinc 产生sinc或函数sin(pit)/(pit)
Square 产生方波
Diric 产生Dirichlet或周期sinc函数
二、 滤波器分析和实现
函数名 功能
Abs 求绝对值(幅值)Freqs 模拟滤波器频率响应
Angle 求相角
Freqspace 频率响应中的频率间隔
Conv 求卷积
Freqz 数字滤波器频率响应
Fftfilt 重叠相加法FFT滤波器实现
Grpdelay 平均滤波器延迟(群延迟)
Filter 直接滤波器实现
Impz 数字滤波器的冲激响应
Filtfilt 零相位数字滤波
Zplane 离散系统零极点图
Filtie Filter 函数初始条件选择
三、 线性系统变换
函数名 功能
Convmtx 卷积矩阵Ss2tf 变系统状态空间形式为传递函数形式
Ploy2rc 从多项式系数中计算反射系数
Ss2zp 变系统状态空间形式为零极点增益形式
Rc2ploy 从反射系数中计算多项式系数
Tf2ss 变系统传递函数形式为状态空间形式
Residuez Z变换部分分式展开或留数计算
Tf2zp 变系统传递函数形式为零极点增益形式
Sos2ss 变系统二阶分割形式为状态空间形式
Zp2sos 变系统零极点形式为二阶分割形式
Sos2zp 变系统二阶分割形式为零极点增益形式
Zp2tf 变系统零极点增益形式为传递函数形式
Ss2sos 变系统状态空间形式为二阶分割形式
四、 IIR滤波器设计
Besself Bessel(贝塞尔)模拟滤波器设计Cheby2 Chebyshev(切比雪夫)II型模拟滤波器设计
Butter Butterworth(巴特沃思)模拟滤波器设计
Ellip 椭圆模拟滤波器设计
Cheby1 Chebyshev(切比雪夫)I 型模拟滤波器设计
Yulewalk 递归数字滤波器设计
五、 IIR滤波器阶选择
Buttord Butterworth(巴特沃思)滤波器阶的选择Cheb2ord Chebyshev(切比雪夫)II型滤波器阶的选择
Ehebord Chebyshev(切比雪夫)I 型滤波器阶的选择
Clipord 椭圆滤波器设计阶的选择 模拟原型滤波器设计
Besselap Bessel模拟低通滤波器原型
Cheb2ap Chebyshev(切比雪夫)II型低通滤波器原型
Buttap Butterworth(巴特沃思)模拟低通滤波器原型
Ellipap 椭圆模拟低通滤波器原型
Cheb1ap Chebyshev(切比雪夫)I 型低通滤波器原型
六、 频率变换
Lp2bp 低通到带通模拟滤波器转换Lp2bs 低通到带阻模拟滤波器变换
Lp2hp 低通到高通模拟滤波器变换
Lp2lp 低通到低通模拟滤波器转换
七、 滤波器离散化
Blinear 双线性变换Impinvar 冲激响应不变法
八、 FIR滤波器设计
Fir1 基于窗函数的 FIR 滤波器设计—标准响应Intfilt 内插FIR滤波器设计
Fir2 基于窗函数的 FIR 滤波器设计—任意响应
Remez Firls 最小二乘FIR滤波器设计
Remezord Parks-McCellan 最优 FIR 滤波器 j阶估计
九、 窗函数
Boxcar 矩形窗Hanning Hanning(汉宁)窗
Triang 三角窗
Blackman Blackman(布莱克曼)窗
Bartlett Bartlett(巴特得特)窗
Chebwin Chebyshev(切比雪夫)窗
Hamming Hamming(汉明)窗
Kaiser Kaiser(凯泽)窗
十、 变换
Ctz 线性调频Z变换Fft 一维快速傅里叶变换
Dct 离散余弦变换
Ifft 一维快速傅里叶逆变换
Idct 逆离散余弦变换
Fftshift 重新排列 fft的输出
Dftmtx 离散傅里叶变换矩阵
Hilbert Hilbert(希尔伯特)变换
十一、 统计信号处理
Cov 协方差矩阵Psd 信号功率谱密度(PSD)估计
Xcov 互协方差函数估计
Tfe 从输入输出中估计传递函数
Corrcoef 相关系数矩阵
Periodogram 采用周期图法估计功率谱密度
Xcoor 互相关系数估计
Pwelch 采用 Welch方法估计功率谱密度
Cohere 相关函数平方幅值估计
Rand 生成均匀分布的随机数
Csd 互谱密度估计
Randn 生成正态分布的随机数
十二、 自适应滤波器部分
Adaptfiltlms 最小均方(LMS)自适应算法Adaptfiltrls 递推最小二乘(RLS)自适应算法
Adaptfiltnlms 归一化最小均方(NLMS)自适应算法
十三、 时频分析与小波变换部分
Spectrogram 短时傅里叶变换Idwt 单级离散一维小波逆变换
Waveinfo 介绍小波工具箱中所有小波的信息
Wavedec 多级离散一维小波分解
Cwt 连续一维小波变换
Appcoef 一维小波变换近似系数
Dwt 单级离散一维小波变换
Detcoef 一维小波变换细节系数
十四、 二维信号处理
Conv2 二维卷积Xcorr2 二维互相关参数
Fft2 二维快读傅里叶变换
Dwt2 单级离散二维小波变换
Ifft2 二维逆快速傅里叶变换
Idwt2 单级离散二维小波逆变换
Filter2 二维数字滤波器
Waverec2 多级离散二维小波分解
匹配滤波器是基于信号与噪声比最大化准则设计的,常用于解决信号检测(信号发现)问题,即确定信号是否存在,一般来说,它是由于信号形状严重失真所引起的。此时,假设信号的形状是已知的,并且在滤波过程中的失真不具有实质意义。地球物理数据处理中,信号(异常)的形状可以通过正演问题来确定,也可以通过对邻近具有类似地质条件的异常体上观测值的分析直接确定。在处理地震记录时,信号形状的估计是由原始地震道的自相关或互相关来确定,此外,也存在这样一类震源,像可控源,它们在d性介质中所激发的信号波形是已知的。所以,已知信号形状(即确定异常)的滤波问题,在地球物理勘探中具有实际价值,它的意义与日俱增,是与寻找深部矿体相联系,由于深部矿体异常的确定受到自然界各种各样的噪声干扰,而且噪声的强度往往要大于异常的强度,即通常所说的弱异常提取。
为了求出匹配滤波器的系统函数(频率响应),我们从加性场模型x(n)=s(n)+v(n)入手,其中s(n)信号(异常),v(n)是平稳随机过程。滤波器的输出是非周期信号。假设信号的形状s(n)已知,因此它的傅氏变换S(ejω)(频谱)也是已知的,即
地球物理信息处理基础
滤波器(具有频率响应H(ejω))的输出端信号为
地球物理信息处理基础
滤波器输入端噪声的能量通过功率谱密度Pvv(ejω)来表示
地球物理信息处理基础
考虑到滤波器是线性的,那么在平稳随机过程的变换中,平稳随机过程的功率谱乘以频率响应模的平方(|H(ejω)|2),便得到滤波器的输出端噪声能量的表达式,即
地球物理信息处理基础
滤波器输出端在时刻t0的信噪比为
地球物理信息处理基础
将式(3-17)和式(3-18)代入式(3-19)得
地球物理信息处理基础
在滤波器的输出端y(t0)=so(t0)+no(t0),此时,要求μo最大化。换句话说,要求so(t0)≫no(t0)。利用布尼科夫斯基-舒瓦茨不等式,有
地球物理信息处理基础
那么,式(3-20)为
地球物理信息处理基础
显然,上式是μo的上限,因此,若取
地球物理信息处理基础
则式(3-21)的左端为
地球物理信息处理基础
则式(3-21)的右端为
地球物理信息处理基础
地球物理信息处理基础
故式(3-21)的左右两端相等,由式(3-20)得到最大的信噪比表达式为
地球物理信息处理基础
这样,式(3-24)就是所求的匹配滤波器的频率响应。
不等式(3-22)的物理意义就在于频率间隔(ω,ω+dω)内,有用信号的频谱幅度越大,且噪声的功率谱越小,那么该频带内有用信号通过的程度就越高,此时,根据式(3-24)可推出,信号与噪声的功率谱差别越大,在滤波器的输出端信噪比就越大。
一般可设t0=0,那么根据式(3-23),得
地球物理信息处理基础
如果噪声是不相关的白噪,那么便有Pvv(ejω)=常数,有
H(ejω)=cS(ejω) (3-26)
由式(3-26)可得,滤波器的频率响应完全由信号的频谱确定,换句话说,滤波器的系统函数是与信号的频谱相匹配的,这也就是为什么称它为匹配滤波器。若将式(3-26)变换到时域(空间域,s(n)为实数),那么有
h(t)=cs(-t);h(n)=cs(-n) (3-27)
由式(3-27)同样可得,若噪声不相关,滤波器的单位冲激响应完全由信号的形状确定,即单位冲激响应的形状与信号的形状相匹配。不难看出,式(3-27)中的常数就是噪声的方差倒数,即c=1/σ2。由式(3-25)得
Pvv(ejω)H(ejω)=S(ejω)
根据卷积质,便有
地球物理信息处理基础
若噪声的自相关函数已知,且信号的形状已知,那么式(3-28)为一匹配滤波器的线性方程组,可以写成矩阵的形式
地球物理信息处理基础
在噪声不相关的情况下,有rvv(0)=σ2,rvv(m)=0(m≠0),由此可见,h(n)=s(-n)/σ2。
此外,由μo对h求导,并令其等于零,可以得到最佳匹配滤波器的方程(滤波器的最佳性是相对其它给出较低信噪比的滤波器而言的),根据式(3-9),有
地球物理信息处理基础
地球物理信息处理基础
显然,此时的分母大于零,所以只要取分子为零即可
地球物理信息处理基础
并将其表示成矩阵形式
(sTh)·s·(hTRvvh)-(sTh)2·(Rvvh)=0⇒(hTRvvh)·s-(sTh)·(Rvvh)=0
可以看出,(hTRvvh)和(sTh)为标量,s和Rvvh为列矢量,由此可得
地球物理信息处理基础
这样,又一次得到了式(3-28),不过应注意,信号矢量为s=[s(n),s(n-1),s(n-2),s(n-3),…,s(n-M)]T,下面举一例来说明。
[例3-2]信号与噪声与[例3-1]相同,即[s(0),s(1)]=(3,1),[v(0),v(1)]=(1,0),设计一个最佳匹配滤波器。
解:噪声的自相关矩阵和信号的矢量分别为
地球物理信息处理基础
根据式(3-29)得
地球物理信息处理基础
即h(0)=1,h(1)=3,若考虑归一化,则
地球物理信息处理基础
这时,最大信噪比为
地球物理信息处理基础
最大功率信噪比为
地球物理信息处理基础
[例3-1]中的最佳滤波器(维纳滤波器)的最大信噪比和最大功率信噪比分别为9363、10399,而对于最佳匹配滤波器分别为10、118。由此可见,最佳匹配滤波器同时能保证信噪比和功率信噪比为最大。
图2-6为匹配滤波器应用的另一个例子,由图可见,在相关噪声的干扰下,滤波器的单位冲激响应与信号的形状不再是匹配关系。因此,信号的形状与滤波器的单位冲激响应仅仅在噪声不相关的情况下相互匹配。在第83个观测点处,过滤后的信号最大值与有用信号的中心点相对应,这就是解决在信号形状失真情况下来发现信号的问题。
在地震数据处理中,匹配滤波器也得到了应用,例如用于发现瑞利面波的存在。在各种校正滤波器的设计方法中也广泛应用匹配滤波器,这种滤波器的系统函数是由反滤波器与匹配滤波器系统函数的乘积构成。
图3-2 与噪声相关的信号发现
a—信号s和相关噪声v(n);b—噪声的相关系数;c—滤波器的单位冲激响应;d—信号s与相关噪声v(n)的和x(n);e—匹配滤波器处理结果∑x(i)s(n-i)/σ2;f—信号存在的后验概率
欢迎分享,转载请注明来源:内存溢出
评论列表(0条)