我正在尝试使用 Python 将具有时变截止频率的带通滤波器应用于信号。我目前使用的例程将我的信号划分为等长的时间段,然后对于每个段,我应用具有特定时间参数的过滤器,然后将信号重新合并在一起。这些参数基于预先存在的估计。
我似乎遇到的问题是在应用过滤器后出现的每个时间段的边缘都有“涟漪”。这会导致我的信号不连续,从而干扰我的滤波后数据分析。
我希望有人能告诉我是否有任何现有的例程可以在 Python 中实现具有时变参数的过滤器?或者,我将不胜感激有关如何解决此问题的建议。
编辑——下面添加了我想要做的示例。
假设我有一个信号 x(t)。我想用参数为 (100,200) Hz 的带通滤波器过滤前半部分。后半部分我想用参数 (140, 240) Hz 进行过滤。我遍历 x(t),将过滤器应用于每一半,然后重新组合结果。一些示例代码可能如下所示:
outputArray = np.empty(len(x))
segmentSize = len(x) / 2
filtParams = [(100, 200), (140, 240)]
for i in range(2):
tempData = x[i*segmentSize:(i+1)*segmentSize]
tempFiltered = bandPassFilter(tempData, parameters=filtParams[i])
outputArray[i*segmentSize:(i+1)*segmentSize] = tempFiltered
(为了节省空间,假设我有一个执行带通滤波的函数)。
如您所见,数据段不重叠,而是简单地“粘贴”回新数组中。
编辑 2——我的问题的一些示例代码@HD
首先,感谢您迄今为止的重要投入。audiolazy 包看起来是个很棒的工具。
我认为如果我更详细地描述我的目标会更有用。正如我在其他地方发布的那样,我正在尝试使用希尔伯特变换来提取信号的瞬时频率(IF)。我的数据包含大量噪声,但我对 IF 信号所在的带宽有一个很好的估计。然而,我遇到的一个问题是 IF 通常是非平稳的。因此,使用“静态”滤波器方法,我经常需要使用宽带通区域,以确保捕获所有频率。
以下代码演示了增加滤波器带宽对 IF 信号的影响。它包括一个信号生成功能,一个使用 scipy.signal 包的带通滤波器的实现,以及一种提取过滤后信号中频的方法。
from audiolazy import *
import scipy.signal as sig
import numpy as np
from pylab import *
def sineGenerator( ts, f, rate, noiseLevel=None ):
"""generate a sine tone with time, frequency, sample rate and noise specified"""
fs = np.ones(len(ts)) * f
y = np.sin(2*np.pi*fs*ts)
if noiseLevel: y = y + np.random.randn(len(y))/float(noiseLevel)
return y
def bandPassFilter( y, passFreqs, rate, order ):
"""STATIC bandpass filter using scipy.signal Butterworth filter"""
nyquist = rate / 2.0
Wn = np.array([passFreqs[0]/nyquist, passFreqs[1]/nyquist])
z, p, k = sig.butter(order, Wn, btype='bandpass', output='zpk')
b, a = sig.zpk2tf(z, p, k)
return sig.lfilter(b, a, y)
if __name__ == '__main__':
rate = 1e4
ts = np.arange(0, 10, 1/rate)
# CHANGING THE FILTER AFFECTS THE LEVEL OF NOISE
ys = sineGenerator(ts, 600.0, 1e4, noiseLevel=1.0) # a 600Hz signal with noise
filts = [[500, 700], [550, 650], [580, 620]]
for f in filts:
tempFilt = bandPassFilter( ys, f, rate, order=2 )
tempFreq = instantaneousFrequency( tempFilt, rate )
plot( ts[1:], tempFreq, alpha=.7, label=str(f).strip('[]') )
ylim( 500, 750 )
xlabel( 'time' )
ylabel( 'instantaneous frequency (Hz)' )
legend(frameon=False)
title('changing filter passband and instantaneous frequency')
savefig('changingPassBand.png')
信号中有一个单一的频率分量(600Hz)。图例显示了每种情况下使用的过滤器参数。使用更窄的“静态”过滤器可以提供“更干净”的输出。但是我的滤波器可以有多窄受频率的限制。例如,考虑以下具有两个频率分量的信号(一个在 600Hz,另一个在 650Hz)。
在此示例中,我被迫使用更宽的带通滤波器,这导致额外的噪声潜入 IF 数据。
我的想法是,通过使用时变滤波器,我可以在特定时间增量为我的信号“优化”滤波器。因此,对于上述信号,我可能希望在前 5 秒内过滤 580-620Hz 左右,然后在接下来的 5 秒内过滤 630-670Hz。本质上,我想最大限度地减少最终 IF 信号中的噪声。
根据您发送的示例代码,我编写了一个函数,该函数使用 audiolazy 在信号上实现静态巴特沃斯滤波器。
def audioLazyFilter( y, rate, Wp, Ws ):
"""implement a Butterworth filter using audiolazy"""
s, Hz = sHz(rate)
Wp = Wp * Hz # Bandpass range in rad/sample
Ws = Ws * Hz # Bandstop range in rad/sample
order, new_wp_divpi = sig.buttord(Wp/np.pi, Ws/np.pi, gpass=dB10(.6), gstop=dB10(.1))
ssfilt = sig.butter(order, new_wp_divpi, btype='bandpass')
filt_butter = ZFilter(ssfilt[0].tolist(), ssfilt[1].tolist())
return list(filt_butter(y))
使用此过滤器与希尔伯特变换例程结合获得的 IF 数据与使用 scipy.signal 获得的数据相比非常好:
AL_filtered = audioLazyFilter( ys, rate, np.array([580, 620]), np.array([570, 630]) )
SP_filtered = bandPassFilter( ys, [580, 620], rate, order=2 )
plot(ts[1:], instantaneousFrequency( SP_filtered, 1e4 ), alpha=.75, label='scipy.signal Butterworth filter')
plot(ts[1:], instantaneousFrequency( AL_filtered, 1e4 ), 'r', alpha=.75, label='audiolazy Butterworth filter')
我现在的问题是,我可以将 audiolazy Butterworth 例程与您在原始帖子中描述的随时间变化的属性结合起来吗?