python绘制地震动傅里叶谱

最新推荐文章于 2024-12-08 00:13:32 发布

原创

最新推荐文章于 2024-12-08 00:13:32 发布 · 1.7k 阅读

1 ·

CC 4.0 BY-SA版权

文章标签：

#python

本文详细介绍了地震动的傅里叶变换与逆变换，阐述了快速傅里叶变换（FFT）在Python中的应用，并提供了相关代码。傅里叶谱振幅与快速傅里叶变换结果之间的关系被解析，强调了频率、相位的计算以及傅里叶积分在处理非周期现象中的重要性。

先说结论

快速傅里叶变换结果为 $y_k$ ，我们需要的傅里叶谱振幅为 $F_k$ ，二者的关系如下
$F_k =\Delta t\cdot|y_k|$
需要注意： $F_0 = \Delta t/2\cdot|y_0|$ ， $F_{N/2} = \Delta t/2\cdot|y_{N/2}|$

python代码

import numpy as np
from scipy.fft import fft


def fourier(ag, dt):
    n = len(ag) - 1
    fw = fft(ag)
    i_half = n // 2
    fw_half = fw[:i_half + 1]  # 取单边频谱
    ampl = np.abs(fw_half) * dt  # 幅值
    ampl[0] = 0.  # 直流分量置零
    ampl[-1] = ampl[-1] / 2
    freq = np.linspace(0, 1 / (2 * dt), i_half + 1)
    return ampl, freq

原因

1 地震动的傅里叶变换与逆变换

地震动 $x_m$ 是离散值，傅里叶逆变换为
$x_m=\frac{A_{0}}{2}+\sum_{k=1}^{N / 2-1}\left[A_{k} \cos \frac{2 \pi k m}{N }+B_{k} \sin \frac{2 \pi k m}{N }\right]+\frac{A_{N / 2}}{2} \cos \frac{2 \pi(N / 2) m}{N }\tag{1}$
傅里叶变换为
$A_{k}=\frac{2}{N} \sum_{m=0}^{N-1} x_{m} \cos \frac{2 \pi k m}{N} \quad k=0,1,2, \cdots, \frac{N}{2}-1, \frac{N}{2}\tag{2}$