MATLAB 中的高通滤波

Question

有谁知道如何在MATLAB中使用过滤器吗？ 我不是一个爱好者，所以我不关心滚降特性等 - 我有一个以 100 kHz 采样的一维信号向量 x，我想对其执行高通滤波（例如，拒绝任何低于 10Hz 的信号） ）以消除基线漂移。

帮助中描述了巴特沃斯、椭圆和切比雪夫滤波器，但没有简单说明如何实现。

Answer 1

可以使用多种过滤器，过滤器的实际选择取决于您想要实现的目标。既然您提到了巴特沃斯、切比雪夫和椭圆滤波器，我假设您正在寻找一般的 IIR 滤波器。

维基百科是开始阅读不同过滤器及其用途的好地方。例如，巴特沃斯在通带中最大程度平坦，而响应在阻带中滚降。在切比雪夫中，您可以在通带（类型 2）或阻带（类型 1）及更大范围内获得平滑响应，另一个带中存在不规则波纹，最后，在椭圆滤波器中，两个带中都有波纹。下图取自维基百科。

因此，在这三种情况下，您都必须用某些东西来交换其他东西。在巴特沃斯中，不会出现纹波，但频率响应滚降速度较慢。在上图中，需要从 0.4 到大约 0.55 才能达到一半功率。在切比雪夫中，你会得到更陡峭的滚降，但你必须考虑到其中一个频带中存在不规则且更大的波纹，而在椭圆中，你会几乎立即被切断，但在两个频带中都有波纹。

过滤器的选择完全取决于您的应用。您是否想在几乎没有损失的情况下获得干净的信号？然后，您需要能够在通带（巴特沃斯/切比2）中提供平滑响应的东西。您是否试图消除阻带中的频率，并且不介意通带响应中的轻微损失？然后您将需要在阻带中平滑的东西（Cheby1）。您是否需要极其尖锐的截止角，即任何稍微超出通带的东西都会对您的分析产生不利影响？如果是这样，您应该使用椭圆过滤器。

关于 IIR 滤波器需要记住的一点是它们有极点。与 FIR 滤波器不同，您可以增加滤波器的阶数，唯一的影响是滤波器延迟，而增加 IIR 滤波器的阶数将使滤波器不稳定。我所说的不稳定是指极点位于单位圆之外。要了解为什么会这样，您可以阅读有关 IIR 滤波器 的 wiki 文章，尤其是有关稳定性的部分。

为了进一步说明我的观点，请考虑以下带通滤波器。

fpass=[0.05 0.2];%# passband
fstop=[0.045 0.205]; %# frequency where it rolls off to half power
Rpass=1;%# max permissible ripples in stopband (dB)
Astop=40;%# min 40dB attenuation
n=cheb2ord(fpass,fstop,Rpass,Astop);%# calculate minimum filter order to achieve these design requirements

[b,a]=cheby2(n,Astop,fstop);

现在，如果您使用 zplane(b,a) 查看零极点图，您会发现有几个极点 (x) 位于单位圆之外，这使得这种方法不稳定。

从频率响应完全混乱的事实可以明显看出这一点。使用 freqz(b,a) 获取以下内容

获取更多信息如果您想要稳定滤波器满足您的具体设计要求，您需要在 MATLAB 中使用 zpk 方法而不是 ba 来使用二阶滤波器。下面是与上面相同的滤波器的操作方法：

[z,p,k]=cheby2(n,Astop,fstop);
[s,g]=zp2sos(z,p,k);%# create second order sections
Hd=dfilt.df2sos(s,g);%# create a dfilt object.

现在，如果您查看该滤波器的特性，您会发现所有极点都位于单位圆内（因此稳定）并且符合设计要求

方法butter 和 ellip 类似，具有等效的 buttord 和 ellipord。 MATLAB 文档还提供了有关设计滤波器的很好的示例。您可以在这些示例和我的示例的基础上根据您的需要设计过滤器。

要对数据使用过滤器，您可以执行 filter(b,a,data) 或 filter(Hd,data)，具体取决于您最终使用的过滤器。如果您想要零相位失真，请使用filtfilt。但是，这不接受 dfilt 对象。因此，要使用 Hd 进行零相位滤波器，请使用 filtfilthd Mathworks 文件交换网站上提供的文件

编辑

这是对@DarenW 评论的回应。平滑和过滤是两种不同的操作，尽管它们在某些方面相似（移动平均是低通滤波器），但您不能简单地用其中一种代替另一种，除非您可以确定它不会具体应用中要注意。

例如，实现达伦对 0-25kHz 线性调频信号的建议，以 100kHz 采样，这是使用高斯滤波器平滑后的频谱

当然，接近 10Hz 的漂移几乎为零。然而，该操作完全改变了原始信号中频率分量的性质。出现这种差异是因为他们完全忽略了平滑操作的滚降（参见红线），并假设它会平坦为零。如果这是真的，那么减法就会起作用。但遗憾的是，事实并非如此，这就是为什么存在设计滤波器的整个领域的原因。

Answer 2

创建您的过滤器 - 例如使用 [B,A] = butter(N,Wn,'high') 其中 N 是过滤器的顺序 - 如果您不确定这是什么，只需设置它到 10。Wn 是在 0 和 1 之间归一化的截止频率，其中 1 对应于信号采样率的一半。如果您的采样率为 fs，并且您想要 10 Hz 的截止频率，则需要设置 Wn = (10/(fs/2))。

然后，您可以使用 Y = filter(B,A,X) 应用过滤器，其中 X 是您的信号。您还可以查看 filtfilt 函数。

Answer 3

进行这种滤波的一种廉价方法是：

* Make a copy of the signal
* Smooth it.   For a 100KHz signal and wanting to eliminate about 10Hz on down, you'll need to smooth over about 10,000 points.   Use a Gaussian smoother, or a box smoother maybe 1/2 that width twice, or whatever is handy.  (A simple box smoother of total width 10,000 used once may produce unwanted edge effects)
* Subtract the smoothed version from the original.  Baseline drift will be gone.

如果原始信号是尖峰的，您可能需要在大信号之前使用一个短中值滤波器， 而不需要在设计、零点、极点和纹波等方面使脑细胞紧张。更顺畅。

这很容易推广到 2D 图像、3D 体积数据等。