這兩天一直在研究FFT算法计贰,研究了兩天總算整明白了快速傅里葉變換算法掖棉,今天就寫一篇文章來總結(jié)一下墓律。
FFT和IFFT的Python語言實現(xiàn)源代碼
直接把我用了一個晚上寫好的快速傅里葉變換和快速傅里葉逆變換的Python語言代碼貼出,關鍵部分有注釋幔亥,里面只用到了Python標準庫cmath庫耻讽,因為要計算cos、sin函數(shù)的值帕棉。直接復制到自己的Python程序中就可以直接使用了针肥。
"""
@Author: Sam
@Function: Fast Fourier Transform
@Time: 2020.02.22 16:00
"""
from cmath import sin, cos, pi
class FFT_pack():
def __init__(self, _list=[], N=0): # _list 是傳入的待計算的離散序列,N是序列采樣點數(shù)笤昨,對于本方法祖驱,點數(shù)必須是2^n才可以得到正確結(jié)果
self.list = _list # 初始化數(shù)據(jù)
self.N = N
self.total_m = 0 # 序列的總層數(shù)
self._reverse_list = [] # 位倒序列表
self.output = [] # 計算結(jié)果存儲列表
self._W = [] # 系數(shù)因子列表
for _ in range(len(self.list)):
self._reverse_list.append(self.list[self._reverse_pos(_)])
self.output = self._reverse_list.copy()
for _ in range(self.N):
self._W.append((cos(2 * pi / N) - sin(2 * pi / N) * 1j) ** _) # 提前計算W值,降低算法復雜度
def _reverse_pos(self, num) -> int: # 得到位倒序后的索引
out = 0
bits = 0
_i = self.N
data = num
while (_i != 0):
_i = _i // 2
bits += 1
for i in range(bits - 1):
out = out << 1
out |= (data >> i) & 1
self.total_m = bits - 1
return out
def FFT(self, _list, N, abs=True) -> list: # 計算給定序列的傅里葉變換結(jié)果瞒窒,返回一個列表捺僻,結(jié)果是沒有經(jīng)過歸一化處理的
"""參數(shù)abs=True表示輸出結(jié)果是否取得絕對值"""
self.__init__(_list, N)
for m in range(self.total_m):
_split = self.N // 2 ** (m + 1)
num_each = self.N // _split
for _ in range(_split):
for __ in range(num_each // 2):
temp = self.output[_ * num_each + __]
temp2 = self.output[_ * num_each + __ + num_each // 2] * self._W[__ * 2 ** (self.total_m - m - 1)]
self.output[_ * num_each + __] = (temp + temp2)
self.output[_ * num_each + __ + num_each // 2] = (temp - temp2)
if abs == True:
for _ in range(len(self.output)):
self.output[_] = self.output[_].__abs__()
return self.output
def FFT_normalized(self, _list, N) -> list: # 計算給定序列的傅里葉變換結(jié)果,返回一個列表崇裁,結(jié)果經(jīng)過歸一化處理
self.FFT(_list, N)
max = 0 # 存儲元素最大值
for _ in range(len(self.output)):
if max < self.output[_]:
max = self.output[_]
for _ in range(len(self.output)):
self.output[_] /= max
return self.output
def IFFT(self, _list, N) -> list: # 計算給定序列的傅里葉逆變換結(jié)果匕坯,返回一個列表
self.__init__(_list, N)
for _ in range(self.N):
self._W[_] = (cos(2 * pi / N) - sin(2 * pi / N) * 1j) ** (-_)
for m in range(self.total_m):
_split = self.N // 2 ** (m + 1)
num_each = self.N // _split
for _ in range(_split):
for __ in range(num_each // 2):
temp = self.output[_ * num_each + __]
temp2 = self.output[_ * num_each + __ + num_each // 2] * self._W[__ * 2 ** (self.total_m - m - 1)]
self.output[_ * num_each + __] = (temp + temp2)
self.output[_ * num_each + __ + num_each // 2] = (temp - temp2)
for _ in range(self.N): # 根據(jù)IFFT計算公式對所有計算列表中的元素進行*1/N的操作
self.output[_] /= self.N
self.output[_] = self.output[_].__abs__()
return self.output
def DFT(self, _list, N) -> list: # 計算給定序列的離散傅里葉變換結(jié)果,算法復雜度較大拔稳,返回一個列表葛峻,結(jié)果沒有經(jīng)過歸一化處理
self.__init__(_list, N)
origin = self.list.copy()
for i in range(self.N):
temp = 0
for j in range(self.N):
temp += origin[j] * (((cos(2 * pi / self.N) - sin(2 * pi / self.N) * 1j)) ** (i * j))
self.output[i] = temp.__abs__()
return self.output
if __name__ == '__main__':
list = [1, 2, 3, 4, 5, 6, 7, 8, 0, 0, 0, 0, 0, 0, 0, 0]
a = FFT_pack().FFT(list, 16, False)
print(a)
下面對源代碼的一些關鍵部分進行解釋:
我直接將計算方法的調(diào)用封裝成了一個類,這樣可以很方便進行擴展和調(diào)用巴比。在FFT_pack()
這個類之中术奖,我定義了下面的幾種方法
1.構(gòu)造函數(shù)__init__()
礁遵,目的是初始化需要進行FFT變換的序列,對采樣點數(shù)進行賦值采记、計算分治算法需要進行分組的層數(shù)佣耐,計算旋轉(zhuǎn)因子的系數(shù)列表等。
2._reverse_pos()
方法唧龄,是為了獲得位倒序后的順序兼砖,這個方法是一個不需要外部調(diào)用的方法。
3.FFT()
方法既棺,顧名思義讽挟,是計算給定序列的快速傅里葉變換結(jié)果。里面涉及到四個參數(shù)丸冕,在實際調(diào)用的時候需要傳入三個參數(shù)-list
N
abs
耽梅。_list
參數(shù)是需要進行FFT運算的列表,N
參數(shù)是參加計算的點的個數(shù)晨仑,注意: 這里面N的值必須是2的m次方褐墅,例如N可以是2拆檬、4洪己、8、16竟贯、32答捕、……1024、2048屑那、4096等這樣的數(shù)拱镐,如果填入了別的數(shù),無法得到正確的計算結(jié)果持际。abs
參數(shù)是缺省值為True的參數(shù)沃琅,當abs賦值為True或者使用缺省值的時候,返回的FFT運算結(jié)果是取絕對值以后的序列蜘欲,當abs參數(shù)賦值為False時益眉,返回的FFT運算結(jié)果就是一個復數(shù),含有實部和虛部的值姥份。
4.FFT_normalized()
方法郭脂,用法與FFT()
方法類似,只不過沒有abs參數(shù)澈歉,方法的返回值是經(jīng)過歸一化處理的FFT變換結(jié)果展鸡。
5.IFFT()
方法,返回給定序列的快速傅里葉逆變換的序列埃难。
6.DFT()
方法莹弊,返回給定序列的離散傅里葉變換序列涤久,返回結(jié)果是經(jīng)過取絕對值運算的。這個DFT()方法主要是用來與FFT方法的運算性能進行對比的忍弛,實際使用中還是使用FFT方法拴竹。
下面總結(jié)一下FFT算法
對于每一個鰈形運算最小單元
對于不同旋轉(zhuǎn)因子的種類和選取方法
將原始序列進行位倒序
IFFT算法
根據(jù)IDFT算法的定義式,可以知道逆變換與正變換的差別非常小剧罩,旋轉(zhuǎn)因子表每一項的方次取負值栓拜,運算之后的序列整個乘以1/N即可雏吭。