傅里叶变换的python实现

周期信号的频谱
  为了能既方便又明白地表示一个信号在不同频率下的幅值和相位,可以采用成为频谱图的表示方法。
  在傅里叶分析中,把各个分量的幅度|Fn|或 Cn 随着频率nω1的变化称为信号的幅度谱。
  而把各个分量的相位 φn 随角频率 nω1 变化称为信号的相位谱。
  幅度谱和相位谱通称为信号的频谱。
  三角形式的傅里叶级数频率为非负的,对应的频谱一般称为单边谱;指数形式的傅里叶级数频率为整个实轴,所以称为双边谱。
  下面以周期信号函数作为示范,看看傅里叶级别函数应该怎么画相位谱和幅度谱
  周期函数:
  在这里插入图片描述

函数相应的分量幅度:
在这里插入图片描述

最终傅里叶级数函数的单边图、双边图、相位谱、幅度谱,如下图所示:
这里写图片描述
在Python中,fft(快速傅里叶变换)算法是在NumPy库中实现的。要使用fft函数,需要首先导入NumPy库。以下是示例代码:

import numpy as np

# 创建一个测试信号
t = np.arange(0, 1, 0.001)  # 时间轴
f = 5  # 频率
x = np.sin(2 * np.pi * f * t)

# 进行快速傅里叶变换
X = np.fft.fft(x)

# 获取频谱信息
frequencies = np.fft.fftfreq(len(t), 1/1000)  # 获取频率轴
amplitudes = np.abs(X)  # 获取振幅谱

# 输出结果
for freq, amp in zip(frequencies, amplitudes):
    print(“频率: {:.2f} Hz, 振幅: {:.2f}.format(freq, amp))

在上述代码中,我们首先导入numpy库,并创建一个测试信号(正弦波)。然后,使用np.fft.fft函数对信号进行快速傅里叶变换,得到频谱信息。最后,我们通过np.fft.fftfreq函数获取频率轴,并使用np.abs函数获取振幅谱。在输出结果时,我们遍历频率和振幅,并将其格式化输出。

值得注意的是,如果需要进行逆变换,可以使用np.fft.ifft函数。此外,还可以使用np.fft.fftshift函数将频谱进行频移,以方便在频率轴上显示。
Python的FFT(快速傅立叶变换)函数可以在SciPy包中找到。SciPy是一个基于Python的科学计算库,提供了许多数值计算、优化、插值和统计函数。FFT是其中一个核心功能之一。

下面是关于Python的FFT的一些详细信息:

  1. 导入包:要使用Python中的FFT函数,首先需要导入相关的包。可以使用以下代码导入SciPy包中的FFT函数:
  from scipy.fft import fft

上述代码将从SciPy包中导入FFT函数,该函数是基于FFTPACK库实现的。

  1. FFT函数:scipy.fft.fft函数用于计算一维实输入序列的傅立叶变换。FFT函数采用一个一维数组作为输入,并返回一个包含频率分量的复数数组。以下是FFT函数的基本语法:
   scipy.fft.fft(x, n=None, axis=-1, norm=None)

参数x是输入数据序列,n是输出FFT序列的长度(默认为x的长度),axis是作用的轴(默认为最后一个维度),norm是归一化方式(默认为不归一化)。

  1. FFT输出:FFT函数将计算出的FFT序列作为输出,其中包含了输入序列在不同频率下的振幅和相位信息。输出的数据类型是复数数组,由于FFT对称性,输出数组中的前半部分表示正频率分量,后半部分表示负频率分量。

  2. 频谱分析:使用FFT函数可以进行频谱分析,从而提取信号的频率信息。频率分量对应于FFT输出数组的索引位置,可以通过计算频率分辨率(采样频率除以FFT长度)确定频率轴上的单位。可以使用如下代码使用FFT进行简单的频谱分析:

  import numpy as np
   import matplotlib.pyplot as plt

   # 生成信号
   fs = 1000  # 采样频率
   t = np.arange(0, 1, 1/fs)  # 时间轴
   f = 10  # 信号频率
   x = np.sin(2*np.pi*f*t)  # 输入信号

   # 计算FFT
   X = fft(x)

   # 绘制频谱图
   freq = np.fft.fftfreq(len(x), 1/fs)  # 计算频率轴
   plt.plot(freq, np.abs(X))  # 绘制频率谱
   plt.xlabel(‘Frequency [Hz])
   plt.ylabel(‘Amplitude’)
   plt.show()
  1. 高级FFT功能:除了基本的FFT函数外,SciPy还提供了其他一些与FFT相关的函数,包括rfft函数(用于计算实输入序列的傅立叶变换)、fftshift函数(用于将FFT输出移动到中心)、ifft函数(用于计算逆傅立叶变换)等。这些函数可根据具体需求进行使用。

在实际应用中,FFT的应用非常广泛,包括音频信号处理、图像处理、信号滤波、频谱分析、信号合成等领域。Python的FFT函数提供了简单而高效的FFT实现,可以方便地进行频率域分析和信号处理。

在Python中,FFT(快速傅里叶变换)算法通常位于科学计算库的NumPy包或者SciPy包中。

NumPy是一个基于Python的科学计算库,提供了强大的多维数组对象和广播功能,以及用于快速操作数组的工具。在NumPy中,FFT算法被实现在numpy.fft模块中。

SciPy是一个基于NumPy的库,提供了许多科学计算的函数和工具。在SciPy中,FFT算法也被实现在scipy.fftpack子模块中。

在使用FFT算法之前,需要先安装相应的库。可以通过以下命令在命令行中安装NumPy和SciPy:

pip install numpy
pip install scipy

下面,我们来看一下在实际代码中如何使用FFT算法。

首先,导入所需的库:

import numpy as np
from scipy.fft import fft

接下来,准备输入数据。FFT算法要求输入的数据是一个一维数组,通常表示为复数形式。假设我们有一个实数序列x,我们可以使用NumPy的array函数将其转换为合适的形式:

x = np.array([1.0, 2.0, 1.0, -1.0, 3.0, 0.5])
x_complex = x.astype(complex)

然后,我们可以使用FFT算法对输入数据进行变换:

X = fft(x_complex)

得到的结果X也是一个一维数组,表示了输入数据的频域表示。我们可以使用NumPy的abs函数来计算幅度谱(频域中对应频率的振幅):

X_mag = np.abs(X)

如果我们只对输入数据的一部分进行变换,可以使用NumPy的切片操作:

N = len(x_complex)  # 输入数据的长度
n = 3  # 变换数据的数量
X_partial = fft(x_complex[:n])

除了正向变换之外,我们还可以进行逆变换以恢复原始数据。逆变换的结果将是复数形式的,但实部可以近似等于原始数据。我们可以使用NumPy的real函数来提取实部:

x_recovered = np.real(ifft(X))

以上就是在Python中使用FFT算法的一般步骤和操作流程。根据实际需求,可能还需要进行一些额外的处理和操作,如频域滤波、频谱绘制等。不过,掌握了以上基本方法,就能够在Python中使用FFT算法进行快速傅里叶变换了。

import numpy as np
from scipy.fftpack import fft, fftshift
 
# 创建一个简单的1秒钟的正弦波信号,采样率为44.1kHz
sample_rate = 44100
t = np.linspace(0, 1, sample_rate)
frequency = 440  # 音调的频率为440Hz
 
# 生成正弦波信号
signal = np.sin(2 * np.pi * frequency * t)
 
# 执行傅里叶变换
fft_signal = fft(signal)
 
# 将DC分量移动到频谱的中心
fft_shifted = fftshift(fft_signal)
 
# 计算幅度谱和相位谱
magnitude = np.abs(fft_shifted)
phase = np.angle(fft_shifted)
 
# 计算频率轴
frequency_axis = np.fft.fftfreq(sample_rate, 1 / sample_rate)
 
# 打印结果
print("频率轴:", frequency_axis[:len(frequency_axis) // 2 + 1])
print("幅度谱:", magnitude[:len(magnitude) // 2 + 1])
print("相位谱:", phase[:len(phase) // 2 + 1])

在这个例子中,我们首先创建了一个简单的正弦波信号,然后使用scipy.fftpack中的fft函数进行傅里叶变换。fftshift用于将频谱的直流分量移动到频谱的中心,以便于更好地观察频率成分。然后,我们计算幅度谱和相位谱,并且生成频率轴。

参考文献:
https://blog.csdn.net/u014792561/article/details/78662381
https://worktile.com/kb/ask/93595.html

附录:

import numpy as np
import matplotlib.pyplot as plt
from scipy.fft import fft
from scipy.integrate import quad

# 示例信号生成
def generate_signal(freq1, freq2, duration, sample_rate):
    t = np.linspace(0, duration, int(sample_rate * duration), endpoint=False)
    signal = (np.sin(2 * np.pi * freq1 * t) + 0.5 * np.sin(2 * np.pi * freq2 * t))
    return t, signal

# 拉普拉斯变换的定义
def laplace_transform(func, t):
    result, _ = quad(lambda s: func(s) * np.exp(-s * t), 0, np.inf)
    return result

# 生成信号
duration = 1.0  # 持续时间1秒
sample_rate = 1000  # 采样率1000Hz
freq1 = 5  # 频率1
freq2 = 50  # 频率2

t, signal = generate_signal(freq1, freq2, duration, sample_rate)

# 傅里叶变换
N = len(signal)
yf = fft(signal)
xf = np.fft.fftfreq(N, 1/sample_rate)[:N//2]
amplitude_spectrum = 2.0/N * np.abs(yf[:N//2])
power_spectrum = amplitude_spectrum**2
phase_spectrum = np.angle(yf[:N//2])

# 拉普拉斯变换示例
t_values = np.linspace(0, duration, 100)
laplace_results = []
for t_val in t_values:
    index = int(t_val * sample_rate)
    if index < N:  # 确保索引不超出范围
        laplace_results.append(laplace_transform(lambda s: signal[index], t_val))
    else:
        laplace_results.append(0)  # 超出范围时设为0或其他处理

image_path = f'filter_image_data2/傅里叶变换+拉普拉斯变换-功率谱.png'
# 可视化
plt.figure(figsize=(15, 15))

# 原始信号
plt.subplot(5, 1, 1)
plt.plot(t, signal)
plt.title('Original Signal')
plt.xlabel('Time [s]')
plt.ylabel('Amplitude')
plt.grid()

# 傅里叶变换结果
plt.subplot(5, 1, 2)
plt.plot(xf, amplitude_spectrum)
plt.title('Magnitude Spectrum')
plt.xlabel('Frequency [Hz]')
plt.ylabel('Magnitude')
plt.grid()

# 功率谱
plt.subplot(5, 1, 3)
plt.plot(xf, power_spectrum)
plt.title('Power Spectrum')
plt.xlabel('Frequency [Hz]')
plt.ylabel('Power')
plt.grid()

# 相位谱
plt.subplot(5, 1, 4)
plt.plot(xf, phase_spectrum)
plt.title('Phase Spectrum')
plt.xlabel('Frequency [Hz]')
plt.ylabel('Phase [radians]')
plt.grid()

# 拉普拉斯变换结果
plt.subplot(5, 1, 5)
plt.plot(t_values, laplace_results)
plt.title('Laplace Transform of the Signal')
plt.xlabel('t')
plt.ylabel('Laplace Transform Value')
plt.grid()

plt.tight_layout()
plt.savefig(image_path)
plt.close()

在这里插入图片描述

评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值