如何制作高通滤波器？

Question

我有一个海平面数据（时间，y，x）的3D数据矩阵，我通过FFT的平方找到了功率谱，但是有一些低频真正占主导地位。我想通过应用高通滤波器来消除这些低频。。。我该怎么做呢？ 数据集和结构/代码示例如下：

这是数据集并创建阵列：

Yearmin = 2018
Yearmax = 2019
year_len = Yearmax - Yearmin + 1.0 # number of years

direcInput = "filepath"
a = s.Dataset(direcInput+"test.nc", mode='r') 

#creating arrays
lat = a.variables["latitude"][:] 
lon = a.variables["longitude"][:] 
time1 = a.variables["time"][:] #DAYS SINCE JAN 1ST 1950
sla = a.variables["sla"][:,:,:] #t, y, x
time = Yearmin + (year_len * (time1 - np.min(time1)) / ( np.max(time1) - np.min(time1))) 


#detrending and normalizing data 
def standardize(y, detrend = True, normalize = True):
    if detrend == True:
        y = signal.detrend(y, axis=0)
    y = (y - np.mean(y, axis=0))
    if normalize == True:
        y = y / np.std(y, axis=0)
    return y

sla_standard = standardize(sla)

print(sla_standard.shape) = (710, 81, 320)


#fft
fft = np.fft.rfft(sla_standard, axis=0)
spec = np.square(abs(fft))

frequencies = (0, nyquist, df)


#PLOTTING THE FREQUENCIES VS SPECTRUM FOR A FEW DIFFERENT SPATIAL LOCATIONS
plt.plot(frequencies, spec[:, 68,85])
plt.plot(frequencies, spec[:, 23,235])
plt.plot(frequencies, spec[:, 39,178])
plt.plot(frequencies, spec[:, 30,149])
plt.xlim(0,.05)
plt.show()

我的目标是对原始时间序列（sla_标准）制作一个高通滤波器，以消除两个非常大的峰值。我应该使用哪种类型的过滤器？谢谢大家!

Answer 1

使用.axes.Axes.set_ylim设置y轴限制

Axes.set_ylim(self, left=None, right=None, emit=True, auto=False, *, ymin=None, ymax=None)

因此，在您的例子中ymin=None，您可以在开始绘图之前将ymax设置为ymax=60000

因此plt.ylim(ymin=None, ymax=60000)

不应该在这里取出数据，因为它“伪造结果”。你真正想要的是放大图表。如果没有事先告知，独立于您阅读图表的人会错误地解释数据。超出图表的峰值是可以接受的，因为每个人都明白这一点

或：

直接替换数组中的某些值（arr）：

arr[arr > ori] = dest

例如，在您的案例ori=60000和dest=1

大于“>；”60k的所有值将替换为1

Answer 2

不同的滤波器：当您声明滤波器作用于信号的频率时。存在不同的过滤器形状，其中一些过滤器具有复杂的表达式，因为它们需要在实时处理（因果）中实现。然而，在您的情况下，您似乎要对数据进行后期处理。您可以使用傅里叶变换，这需要所有数据（非因果）

要选择的滤波器：因此，您可以通过在频率上应用掩码，直接在傅里叶域中执行滤波操作。如果你想移除频率，我建议你使用由0和1组成的二进制掩码。为什么？因为它是你能想到的最简单的过滤器。这是科学相关的国家，你完全删除了一些频率（说出来，并证明它）。然而，更难断言的是，你让一些衰减，让另一些衰减，你任意选择衰减因子

Python实现

signal_fft = np.fft.rfft(sla_standard,axis=0)
mask = np.ones_like(sla_standard)
mask[freq_to_filter,...] = 0.0 # define here the frequencies to filter
filtered_signal = np.fft.irfft(mask*signal_fft,axis=0)