快捷導(dǎo)航

Python?scipy利用快速傅里葉變換實(shí)現(xiàn)濾波

更新時(shí)間：2024年01月22日 09:02:25 作者：微小冷

這篇文章主要為大家詳細(xì)介紹了Python?scipy如何利用快速傅里葉變換實(shí)現(xiàn)濾波,文中的示例代碼講解詳細(xì),感興趣的小伙伴可以跟隨小編一起學(xué)習(xí)一下

fft模塊簡(jiǎn)介

scipy官網(wǎng)宣稱，fftpack模塊將不再更新，或許不久之后將被廢棄，也就是說fft將是唯一的傅里葉變換模塊。

Fourier變換極其逆變換在數(shù)學(xué)上的定義如下

下表整理出一部分與Fourier變換相關(guān)的函數(shù)，其中FFT為快速Fourier變換(Fast Fourier Transform)；DCT為離散余弦變換(Discrete Cosine Transform)；DST為離散正弦變換(discrete sine transform)，另外，函數(shù)的前綴和后綴有如下含義

i表示逆變換；
2, n分別表示2維和n維

	正變換	逆變換
通用	fft, fft2, fftn	ifft, ifft2, ifftn
實(shí)數(shù)域	rfft, rfft2, rfftn	irfft, irfft2, irfftn
厄米對(duì)稱	hfft, hfft2, hfftn	ihfft, ihfft2, ihfftn
離散余弦變換	dct, dctn	idct, idctn
離散正弦變換	dst, dstn	idst, idstn
漢克爾變換	fht	ifht
移動(dòng)零頻	fftshift	ifftshift
DFT采樣頻率	fftfreq	ifftfreq

fft函數(shù)示例

在數(shù)值計(jì)算中，一切輸入輸出均為離散值，所以實(shí)際上用到的是離散Fourier變換，即DFT，其功能是將離散的采樣變換為一個(gè)離散的頻譜分布。

下面將手動(dòng)創(chuàng)建一組采樣點(diǎn)，并添加一點(diǎn)噪聲，然后通過FFT獲取其頻域信息。

import numpy as np
from scipy import fft

PI = np.pi*2
fs = 60     #采樣頻率
T = 100     #采樣周期數(shù)
N = fs*T    #采樣點(diǎn)
t = np.linspace(0, T, N)
noise = 2 * np.random.randn(*t.shape)
s = 2 * np.sin(PI * t) + 3 * np.sin(22 * PI * t) + noise
F = fft.fft(s)
f = fft.fftfreq(N, 1.0/fs)

其中，t為時(shí)間序列，s為模擬的采樣點(diǎn)，F是Fourier變換之后的結(jié)果。但由于fft默認(rèn)是在復(fù)數(shù)域上的，故而可以查看其實(shí)部、虛部、模和輻角的值。

下面對(duì)采樣點(diǎn)以及Fourier變換的結(jié)果進(jìn)行繪制

import matplotlib.pyplot as plt
fig = plt.figure()
ax = fig.add_subplot(2,2,1)
ax.plot(t, s)
ax.set_title("t vs s")
f_abs = np.abs(F)
ax = fig.add_subplot(2,2,2)
ax.plot(f, f_abs)
ax.set_title("fs vs |F|")

xlims = [[0,2], [21,23]]
for i, xlim in enumerate(xlims):
    ax = fig.add_subplot(2,2,3+i)
    ax.stem(f, f_abs)
    ax.set_title("fs vs |F|")
    ax.set_xlim(xlim)

plt.show()

結(jié)果為

即f=1和f=22處被篩選了出來。

濾波

有了這個(gè)，就可以在頻域上對(duì)數(shù)據(jù)進(jìn)行濾波，其思路是，對(duì)f_abs中的值進(jìn)行閾值分割，例如，只篩選出低頻部分，然后看一下濾波效果

fig = plt.figure(1)
f_filt = F * (np.abs(f) < 2)
s_filt = fft.ifft(f_filt)
ax = fig.add_subplot()
ax.plot(t, s, lw=0.2)
ax.plot(t, s_filt.real, lw=2)
ax.set_title("threshold=2")
ax.set_xlim([0,10])
plt.show()

效果如下

到此這篇關(guān)于Python scipy利用快速傅里葉變換實(shí)現(xiàn)濾波的文章就介紹到這了,更多相關(guān)Python scipy實(shí)現(xiàn)濾波內(nèi)容請(qǐng)搜索腳本之家以前的文章或繼續(xù)瀏覽下面的相關(guān)文章希望大家以后多多支持腳本之家！

您可能感興趣的文章: