目录
2、scipy库介绍+IIR滤波器设计（含GUI）+绘制频谱响应
四、IIR低通滤波器、IIR高通滤波器、IIR带通滤波器、IIR带阻滤波器
一、参考文献
1、scipy接口
2、scipy库介绍+IIR滤波器设计（含GUI）+绘制频谱响应
https://www.bilibili.com/video/BV1D64y1C7fv?from=search&seid=5485332069106154425&spm_id_from=333.337.0.03、读写播放wav音频文件
https://www.cnblogs.com/LXP-Never/p/10078200.html4、利用设计好的滤波器对信号滤波
二、代码
1、音频资源
【python数字信号处理】scipy库设计滤波器（IIR为例）、绘制滤波器频谱响应、IIR滤波器滤波、读写wav音频文件-编解码文档类资源-CSDN下载相关文章：【python数字信号处理】——scipy库设计滤波器（IIR为例）、绘制滤波器频谱响应、更多下载资源、学习资料请访问CSDN下载频道.
https://download.csdn.net/download/qq_45769063/85017358
音频也可以是自己网上找的一段wav文件?
2、滤波器设计代码
#*-----------------设置一个低通IIR滤波器----------------------*#
# 导入模块
import matplotlib.pyplot as plt
from scipy import signal
import numpy as np
import matplotlib
matplotlib.use("QT5Agg")
matplotlib.rcParams["font.family"] = "STSong"
# 设置参数
'''
这里的参数均是归一化后的
通带【0,0.1】，阻带【0,0.2】，因为是低通滤波，因此需要保证通带范围在阻带范围内，反之亦然,如设计一个高通滤波，通带【0.4-1】，则阻带边缘值必须小于0.4;
通带最大损耗dB和阻带最小损耗dB值分别表示在对应的带中滤波时衰减的数值
'''
## 通带边缘频率
wp = 0.1
## 阻带边缘频率
ws = 0.2
## 通带最大损耗dB，表示在通带中的信号值最大的衰减不得高于1dB
gpass = 1
## 阻带最小损耗dB，表示在通带中的信号值最大的衰减不得低于40dB
gstop = 40
# 利用设置好的参数设计一个IIR滤波器
'''
滤波器说白了骑士就是一个传递函数，输出可以选择分子分母的形式来输出吧b,a，该比值即为传递函数
根据传递函数即可计算频率响应，其中w表示频率，h表示对应频率的衰减值
由于计算得到的衰减值有负数，负数的衰减值是没有意义的，因此需要对负数衰减值进行求绝对值
a = [ 0.01173816 -0.03292858 0.0451046 -0.03292858 0.01173816]
b = [ 1. -3.60478527 4.97876142 -3.11559024 0.74467021]
habs = [ 8.91250938e-01+0.00000000e+00j 8.91153448e-01-1.87257159e-02j 。。。]
'''
b,a = signal.iirdesign(wp,ws,gpass,gstop,False,ftype="ellip",output='ba')
print(b,"\n",a)
## 计算滤波器的频率响应
w,h = signal.freqz(b,a)
habs = np.abs(h)
print(w,"\n",h,"\n",habs)
# 绘图查看频率响应图
plt.figure()
plt.title("IIR信号滤波")
plt.subplot(321)
plt.xlabel("频率")
plt.ylabel("衰减")
plt.plot(w,habs)3、音频读取以及快速傅里叶变换FFT代码
# 读取音频文件
import wave
# 打开WAV音频
f = wave.open(r"D:\ParttimeJobs\DigitSignalProcess\resource\sound.wav", "rb")
# 读取格式信息
'''
(声道数、量化位数、采样频率、采样点数、压缩类型、压缩类型的描述)
# (nchannels, sampwidth, framerate, nframes, comptype, compname)
nchannels, sampwidth, framerate, nframes: 1 2 44100 62976
'''
params = f.getparams()
nchannels, sampwidth, framerate, nframes = params[:4]
print("nchannels, sampwidth, framerate, nframes:",nchannels, sampwidth, framerate, nframes)
# 读取nframes个数据，返回字符串格式
str_data = f.readframes(nframes)
f.close()
# 显示原信号波形和频谱
'''
wave_data_shape: (62976,)
'''
## 将字符串转换为数组，得到一维的short类型的数组
wave_data = np.fromstring(str_data, dtype=np.short)
print("wave_data_shape:",wave_data.shape)
plt.subplot(323)
plt.plot(wave_data)
plt.xlabel("采样点")
plt.ylabel("频率")
# 对原信号离散傅里叶变换
wave_data_xf = np.fft.fft(wave_data)
wave_data_xbs = np.abs(wave_data_xf)
plt.subplot(324)
plt.plot(wave_data_xbs[:len(wave_data_xbs)//2]) # 对称的取一半即可
plt.xlabel("频率")
plt.ylabel("幅值")
4、利用IIR滤波器对音频信号滤波显示时域频域波形
# 滤波
'''
filtedData_shape: (62976,)
'''
filtedData = signal.filtfilt(b, a, wave_data)
print("filtedData_shape:",filtedData.shape)
plt.subplot(325)
plt.plot(filtedData)
plt.xlabel("采样点")
plt.ylabel("频率")
# # 对滤波后信号离散傅里叶变换
filtedData_xf = np.fft.fft(filtedData)
filtedData_xbs = np.abs(filtedData_xf)
plt.subplot(326)
plt.plot(filtedData_xbs[:len(filtedData_xbs)//2])
plt.xlabel("频率")
plt.ylabel("幅值")
plt.tight_layout()
plt.show()5、写入滤波结束后的代码
# 将滤波后的信号写入音频
import struct
# wav文件写入
## 获取待写入的信号，数组型，打开wav文件
outData = filtedData
outwave = wave.open(r"D:\ParttimeJobs\DigitSignalProcess\resource\soundfilterlowpass.wav", 'wb')
## 设置音频参数
nchannels = nchannels # 通道数设置为1
sampwidth = sampwidth # 量化位数设置为2
framerate = framerate # 采样频率8000
nframes = len(outData) # 采样点数
comptype = "NONE"
compname = "not compressed"
outwave.setparams((nchannels, sampwidth, framerate, nframes,comptype, compname))
# ## 对每一帧进行打包写入（逐帧写入）
# for i in outData:
# outwave.writeframes(struct.pack('h', int(i * 64000 / 2)))
# # struct.pack(FMT, V1)将V1的值转换为FMT格式字符串
outwave.writeframes(filtedData.tobytes())
outwave.close()6、完整代码
#*-----------------设置一个低通IIR滤波器----------------------*#
# 导入模块
import matplotlib.pyplot as plt
from scipy import signal
import numpy as np
import matplotlib
# matplotlib.use("QT5Agg")
matplotlib.rcParams["font.family"] = "STSong"
# 设置参数
'''
这里的参数均是归一化后的
通带【0,0.1】，阻带【0,0.2】，因为是低通滤波，因此需要保证通带范围在阻带范围内，反之亦然,如设计一个高通滤波，通带【0.4-1】，则阻带边缘值必须小于0.4;
通带最大损耗dB和阻带最小损耗dB值分别表示在对应的带中滤波时衰减的数值
'''
## 通带边缘频率
wp = 0.1
## 阻带边缘频率
ws = 0.2
## 通带最大损耗dB，表示在通带中的信号值最大的衰减不得高于1dB
gpass = 1
## 阻带最小损耗dB，表示在通带中的信号值最大的衰减不得低于40dB
gstop = 40
# 利用设置好的参数设计一个IIR滤波器
'''
滤波器说白了骑士就是一个传递函数，输出可以选择分子分母的形式来输出吧b,a，该比值即为传递函数
根据传递函数即可计算频率响应，其中w表示频率，h表示对应频率的衰减值
由于计算得到的衰减值有负数，负数的衰减值是没有意义的，因此需要对负数衰减值进行求绝对值
a = [ 0.01173816 -0.03292858 0.0451046 -0.03292858 0.01173816]
b = [ 1. -3.60478527 4.97876142 -3.11559024 0.74467021]
habs = [ 8.91250938e-01+0.00000000e+00j 8.91153448e-01-1.87257159e-02j 。。。]
'''
b,a = signal.iirdesign(wp,ws,gpass,gstop,False,ftype="ellip",output='ba')
print(b,"\n",a)
## 计算滤波器的频率响应
w,h = signal.freqz(b,a)
habs = np.abs(h)
print(w,"\n",h,"\n",habs)
# 绘图查看频率响应图
plt.figure()
plt.title("IIR信号滤波")
plt.subplot(321)
plt.xlabel("频率")
plt.ylabel("衰减")
plt.plot(w,habs)
# 读取音频文件
import wave
# 打开WAV音频
f = wave.open(r"D:\ParttimeJobs\DigitSignalProcess\resource\sound.wav", "rb")
# 读取格式信息
'''
(声道数、量化位数、采样频率、采样点数、压缩类型、压缩类型的描述)
# (nchannels, sampwidth, framerate, nframes, comptype, compname)
nchannels, sampwidth, framerate, nframes: 1 2 44100 62976
'''
params = f.getparams()
nchannels, sampwidth, framerate, nframes = params[:4]
print("nchannels, sampwidth, framerate, nframes:",nchannels, sampwidth, framerate, nframes)
# 读取nframes个数据，返回字符串格式
str_data = f.readframes(nframes)
f.close()
# 显示原信号波形和频谱
'''
wave_data_shape: (62976,)
'''
## 将字符串转换为数组，得到一维的short类型的数组
wave_data = np.fromstring(str_data, dtype=np.short)
print("wave_data_shape:",wave_data.shape)
plt.subplot(323)
plt.plot(wave_data)
plt.xlabel("采样点")
plt.ylabel("频率")
# 对原信号离散傅里叶变换
wave_data_xf = np.fft.fft(wave_data)
wave_data_xbs = np.abs(wave_data_xf)
plt.subplot(324)
plt.plot(wave_data_xbs) # 对称的取一半即可
plt.xlabel("频率")
plt.ylabel("幅值")
# 滤波
'''
filtedData_shape: (62976,)
'''
filtedData = signal.filtfilt(b, a, wave_data)
print("filtedData_shape:",filtedData.shape)
plt.subplot(325)
plt.plot(filtedData)
plt.xlabel("采样点")
plt.ylabel("频率")
# # 对滤波后信号离散傅里叶变换
filtedData_xf = np.fft.fft(filtedData)
filtedData_xbs = np.abs(filtedData_xf)
plt.subplot(326)
plt.plot(filtedData_xbs)
plt.xlabel("频率")
plt.ylabel("幅值")
plt.tight_layout()
plt.show()
# 将滤波后的信号写入音频
import struct
# wav文件写入
## 获取待写入的信号，数组型，打开wav文件
outData = filtedData
outwave = wave.open(r"D:\ParttimeJobs\DigitSignalProcess\resource\soundfilterlowpass.wav", 'wb')
## 设置音频参数
nchannels = nchannels # 通道数设置为1
sampwidth = sampwidth # 量化位数设置为2
framerate = framerate # 采样频率8000
nframes = len(outData) # 采样点数
comptype = "NONE"
compname = "not compressed"
outwave.setparams((nchannels, sampwidth, framerate, nframes,comptype, compname))
# ## 对每一帧进行打包写入（逐帧写入）
# for i in outData:
# outwave.writeframes(struct.pack('h', int(i * 64000 / 2)))
# # struct.pack(FMT, V1)将V1的值转换为FMT格式字符串
outwave.writeframes(filtedData.tobytes())
outwave.close()
三、结果
四、IIR低通滤波器、IIR高通滤波器、IIR带通滤波器、IIR带阻滤波器
改变IIR滤波器设计代码中的wp和ws值即可，具体如下：
1、IIR低通滤波器
## 通带边缘频率
wp = 0.1
## 阻带边缘频率
ws = 0.2保证0<=wp<ws<=1
2、IIR高通滤波器
## 通带边缘频率
wp = 0.2
## 阻带边缘频率
ws = 0.1保证0<=ws<wp<=1
?滤波后音频：soundfilterhighpass.wav
?3、IIR带通滤波器
# wp = [wp1,wp2]
# ws = [ws1,ws2]
wp = [0.2,0.3]
ws = [0.1,0.4]?保证：
0 <= ws1 < wp1 < wp2 < ws2 <= 1
?滤波后音频：soundfilterbandpass.wav
?4、IIR带阻滤波器
# wp = [wp1,wp2]
# ws = [ws1,ws2]
wp = [0.1,0.4]
ws = [0.2,0.3]?保证：
0 <= wp1 < ws1< ws2 < wp2 <= 1
?滤波后音频：soundfilterbandstop.wav
?
?