Python使用scipy.fft进行大学经典的傅立叶变换
作者:刘润森! 发布时间:2022-09-10 20:26:13
傅里叶变换是在高数是一个很重要的知识点,今天将结合Python代码实现傅立叶变换。
傅立叶变换
我们平时是如何去分解一个复杂的问题呢?一个经典的方法就是把这个复杂的问题分解成为多个简单的可操作的子问题, 傅立叶变换也是基于这个思想。
傅里叶分析是研究如何将数学函数分解为一系列更简单的三角函数的领域。傅里叶变换是该领域的一种工具,用于将函数分解为其分量频率。
在本教程中,傅立叶变换是一种工具,可以获取信号并查看其中每个频率的功率。看一看该傅立叶变换中的重要术语:
信号:信号是随时间变化的信息。例如,音频,视频和电压走线都是信号的示例。
频率:频率是某物重复的速度。例如,时钟以1赫兹(Hz)的频率滴答,或每秒重复1次。
功率:功率表示每个频率的强度。
下图是一些正弦波的频率和功率的直观演示:
第一个是低频正弦波,第二个是高频正弦波,第三个是低频低功率正弦波,因此低功率正弦波比其它两个正弦波的峰较小。
时域与频域
时域与频域是查看信号的两种不同方式,即信号的组成频率或随时间变化的信息。
在时域中,信号是随时间(x轴)幅度(y轴)变化的波。您最有可能在时域中查看图表,例如:
这是一些音频的图像,它是一个时域信号。横轴表示时间,纵轴表示振幅。
在频域中,信号表示为一系列频率(x轴),每个频率都具有关联的功率(y轴)。下图是经过傅立叶变换后的上述音频信号:
代码实现正弦波
音频本质上是正弦波。
下面是产生正弦波的代码:
import numpy as np
from matplotlib import pyplot as plt
SAMPLE_RATE = 44100 # 赫兹
DURATION = 5 # 秒
def generate_sine_wave(freq, sample_rate, duration):
x = np.linspace(0, duration, sample_rate * duration, endpoint=False)
frequencies = x * freq
y = np.sin((2 * np.pi) * frequencies)
return x, y
# 产生持续5秒的2赫兹正弦波
x, y = generate_sine_wave(2, SAMPLE_RATE, DURATION)
plt.plot(x, y)
plt.show()
x轴以秒为单位表示时间,并且由于每秒钟的时间都有两个峰值,因此可以看到正弦波每秒振荡两次。
混合音频
下面将两个正弦波,混合音频信号仅包括两个步骤:
将正弦波加在一起,然后进行归一化的操作。
具体实现的代码如下。
_, nice_tone = generate_sine_wave(400, SAMPLE_RATE, DURATION)
_, noise_tone = generate_sine_wave(4000, SAMPLE_RATE, DURATION)
noise_tone = noise_tone * 0.3
mixed_tone = nice_tone + noise_tone
下一步是归一化,或缩放信号以适合目标格式。由于以后将如何存储音频,目标格式为16位整数,范围为-32768到32767:
normalized_tone = np.int16((mixed_tone / mixed_tone.max()) * 32767)
plt.plot(normalized_tone[:1000])
plt.show()
看到的正弦波是生成的400 Hz音调,将上面的正弦波转化为音频,最简单的方法是使用SciPy
的wavfile.write
方法将其存储在WAV
文件中。16位整数是WAV文件的标准数据类型,因此需要将信号标准化为16位整数:
from scipy.io.wavfile import write
# 记住,采样率=44100赫兹是我们的播放率
write("mysinewave.wav", SAMPLE_RATE, normalized_tone)
这个音频听起来音调很高。
完成此步骤后,就当作音频样本了。下一步是使用傅立叶变换消除高音调!
傅立叶变换
现在对生成的音频上使用FFT了。FFT是一种算法,可实现傅立叶变换并可以在时域中为信号计算频谱。
from scipy.fft import fft, fftfreq
# 标准化音调中的样本数
N = SAMPLE_RATE * DURATION
yf = fft(normalized_tone)
xf = fftfreq(N, 1 / SAMPLE_RATE)
plt.plot(xf, np.abs(yf))
plt.show()
我们可以在正频率中看到两个峰值,正频率峰值位于400 Hz和4000 Hz,与之前生成的音频的频率相对应。
计算傅里叶变换
yf = fft(normalized_tone)
xf = fftfreq(N, 1 / SAMPLE_RATE)
上面代码的功能
fft() 计算转换本身。
fftfreq()计算的输出中每个仓中心的频率fft()。没有这个,就无法在频谱上绘制x轴
fft()输出的频谱围绕y轴反射,因此负半部分是正半部分的镜像,我们一般只需计算一半对称值,即可更快地进行傅立叶变换。scipy.fft以的形式实施此速度骇客rfft()。
from scipy.fft import rfft, rfftfreq
# 注意前面多余的“r”
yf = rfft(normalized_tone)
xf = rfftfreq(N, 1 / SAMPLE_RATE)
plt.plot(xf, np.abs(yf))
plt.show()
过滤信号
傅里叶变换的一大优点是它是可逆的,我们可以利用此优势来过滤音频并摆脱高音调频率。
# 最大频率为采样率的一半
points_per_freq = len(xf) / (SAMPLE_RATE / 2)
# 我们的目标频率是4000赫兹 将44100变成4000
target_idx = int(points_per_freq * 4000)
然后,您可以将其设置yf为0目标频率附近的index来摆脱它:
yf[target_idx - 1 : target_idx + 2] = 0
plt.plot(xf, np.abs(yf))
plt.show()
由于只有一个高峰,下面应用傅立叶逆变换返回时域。
应用逆FFT与应用FFT相似:
from scipy.fft import irfft
new_sig = irfft(yf)
plt.plot(new_sig[:1000])
plt.show()
由于您正在使用rfft(),因此需要使用irfft()来应用反函数。但是,如果您使用过fft(),则反函数将是ifft()。现在,您的绘图应如下所示:
现在有一个以400 Hz振荡的正弦波,并且您已经成功地消除了4000 Hz的噪声。
对信号进行归一化,然后再将其写入文件。
norm_new_sig = np.int16(new_sig * (32767 / new_sig.max()))
write("clean.wav", SAMPLE_RATE, norm_new_sig)
来源:https://blog.csdn.net/weixin_44510615/article/details/117391598
猜你喜欢
- 本文主要分析的是web.py库的application.py这个模块中的代码。总的来说,这个模块主要实现了WSGI兼容的接口,以便应用程序能
- 接触编程的朋友都听过正则表达式,在python中叫re模块,属于文字处理服务里面的一个模块。re里面有一个方法叫match,接下来的文章我来
- Python中对象的行为是由它的类型 (Type) 决定的。所谓类型就是支持某些特定的操作。数字对象在任何编程语言中都是基础元素,支持加、减
- 本文实例讲述了从ThinkPHP3.2.3过渡到ThinkPHP5.0学习笔记。分享给大家供大家参考,具体如下:用tp3.2.3做了不少项目
- 记得以前写过一篇文章 php有效的过滤html标签,js代码,css样式标签: <?php $str = preg_replace(
- 表示文字链接最清楚的方式是“蓝色文字+下划线”,这是在浏览器发展过程中形成的。这个问题大家都说过很多次了,我也曾经说过。然而,这样的规范却总
- '去掉字符串头尾的连续的回车和空格 function trimVBcrlf(str) tr
- 写好脚本,注册好服务之后,经测试,ORACLE可以随RHEL启动而启动,但不能随系统关闭而关闭。在网上找答案,发现几乎所有的设置过程帖子都是
- 前期准备1、机器人框架的下载和配置首先需要一个qq机器人框架,我使用的是基于mirai 以及 MiraiGo 开发的go-cqhttp(里面
- my.ini文件[mysqld]max_allowed_packet = 10M
- 今天无意在坛子里看到这样一个求救帖(这里),看了一下,感觉问题比较好解决。但是问题背后的问题却引起了我的反思。把他的页面整理一下看看(为了便
- 一、基本思想本文思想是基于用asp和DOM来读取和存储XML数据,并利用XML数据来存储留言信息,达到同用数据库存储数据的功能。二、XML留
- '*****************************************************************
- 建立合理的索引提高SQL Server的性能在应用系统中,尤其在联机事务处理系统中,对数据查询及处理速度已成为衡量应用系统成败的标准。而采用
- 巨坑,切忌不要轻易删除Linux系统自带版本的Python1.卸载python(防止未卸载干净)rpm -qa|grep python|xa
- python数据与matlab互通SciPy有时候需要利用python进行科学计算,但需要Matlab进行交互式画图,因此需要掌握pytho
- 一、数据库远程管理技术 对于中小型应用,比如一个网站的建设和维护,这种大型应用平台就显得有些尾大不掉,开销也过于庞大。曾经在互联网技术和Ja
- go-ini的分区go-ini的多个配置项通过分区(section)来划分。有默认(空)分区和命名的分区,没有给分区命名就是默认分区,默认分
- 控制资源访问前文提到threading库在多线程时,对同一资源的访问容易导致破坏与丢失数据。为了保证安全的访问一个资源对象,我们需要创建锁。
- with语句会设置一个临时的上下文,交给上下文管理器对象控制,并且负责清理上下问题。这样做能避免错误并减少样板代码,因此API能更安全,更易