python基于傅里叶变换的频率滤波－音频降噪

1. 傅里叶简介

法国科学家傅里叶提出，任何一条周期曲线，无论多么跳跃或不规则，都能表示成一组光滑正弦曲线叠加之和。

傅里叶变换即基于傅里叶定理，对一条不规则的曲线进行拆解，从而得到一组光滑正弦曲线函数的过程。

例如：弹钢琴
假设有一时间域函数：y = f(x)，根据傅里叶理论它可以被分解为一系列正弦函数的叠加，这些正弦函数具有不同的振幅A，频率ω或初相位φ。：
$$y=A_{1}sin({\omega }_{1}x+{\varphi }_1)+A_{2}sin({\omega }_{2}x+{\varphi }_2)+A_{2}sin({\omega }_{2}x+{\varphi }_2)+...$$
在信息处理过程中，通常处理步骤是:

• 1.通过傅里叶变换是将时域（即时间域）上的信号转变为频域（即频率域）上的信号，进行相应的滤波去噪、增强、统计分析等处理。
• 2.处理后的信号可以通过傅里叶逆变换在转换成时域信号，从而达到相应的信号处理效果。

2. 示例说明

该示例主要实现通过傅里叶变换实现一段音频数据的降噪处理过程。主要过程如下:

part1. 利用　scipy.io.wavfile.read() 读取音频文件，该接口会返回该音频数据的采用率sample_reate和每个采样点的位移值sigs，根据这组采样点可以绘制音频的时域图，这里为了显示截取了前178个采样点。

part2. 利用numpy.fft.fft()接口完成时域信号到频域信号的转换,并绘制音频的频域图。

part3. 该部分示意一个简单的去噪过程，音频中的高频信号为需要保留的信息，其他都为噪声信号,这里简单的将噪声信号统一设置为０仅保留了高频信号;

part4. 将去噪完的信息通过numpy.fft.ifft()接口逆傅里叶变换成时域信号。

• 测试数据, [nosied.wav] https://pan.baidu.com/s/19vriy2f5XlJid41-4LthbA 提取码: tbii

'''
基于傅里叶变换的频域滤波。
'''
import numpy as np
import numpy.fft as nf
import matplotlib.pyplot as plt
import scipy.io.wavfile as wf

#读取音频文件,将其按照采样率离散化，返回采样率和信号
#sample_reate:采样率(每秒采样个数),　sigs:每个采样位移值。
#================1.原始音频信号,时域信息=================================
sameple_rate,sigs = wf.read('../da_data/noised.wav')
print('采样率:{}'.format(sameple_rate))
print('信号点数量:{}'.format(sigs.size))
sigs = sigs/(2**15)

times = np.arange(len(sigs))/sameple_rate
plt.figure('Filter',facecolor='lightgray')
plt.subplot(221)
plt.title('Time Domain',fontsize=16)
plt.ylabel('Signal',fontsize=12)
plt.grid(linestyle=':')
plt.plot(times[:178],sigs[:178],color='dodgerblue',label='Noised Signal')
plt.legend()

#==================2.转换为频率域信号===================================
#基于傅里叶变换，获取音频频域信息
#绘制音频频域的: 频域/能量图像
freqs = nf.fftfreq(sigs.size, 1/sameple_rate)
complex_arry = nf.fft(sigs)
pows = np.abs(complex_arry)
plt.subplot(222)
plt.title('Frequence Domain',fontsize=16)
plt.ylabel('power',fontsize=12)
plt.grid(linestyle=':')
plt.semilogy(freqs[freqs>0],pows[freqs>0],color='dodgerblue',label='Noised Freq')
plt.legend()

#==============第3步=================================================
#将低能噪声去除后绘制音频频域的: 频率/能量图书
fun_freq = freqs[pows.argmax()] #获取频率域中能量最高的
noised_idx = np.where(freqs != fun_freq)[0] #获取所有噪声的下标
ca = complex_arry[:]
ca[noised_idx] = 0 #高通滤波
filter_pows = np.abs(complex_arry)

plt.subplot(224)
plt.ylabel('power',fontsize=12)
plt.grid(linestyle=':')
plt.plot(freqs[freqs>0],filter_pows[freqs>0],color='dodgerblue',label='Filter Freq')
plt.legend()
#================第4步==============================================
filter_sigs = nf.ifft(ca)
plt.subplot(223)
plt.title('Time Domain',fontsize=16)
plt.ylabel('Signal',fontsize=12)
plt.grid(linestyle=':')
plt.plot(times[:178],filter_sigs[:178],color='dodgerblue',label='Filter Signal')
plt.legend()

#重新生成音频文件
filter_sigs = (filter_sigs*(2**15)).astype('i2')
wf.write('../da_data/filter.wav',sameple_rate,filter_sigs)

plt.tight_layout()
plt.show()