如何用Python对CSV水位时间序列做FFT并绘制cycles per day频率图

前置依赖安装 #

首先确认你已经安装了所需的第三方库，未安装可执行以下命令安装：

pip install pandas numpy matplotlib scipy

步骤1：导入依赖库 #

import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
from scipy.signal import detrend

# 可选配置：解决matplotlib中文显示异常问题
plt.rcParams['font.sans-serif'] = ['SimHei']
plt.rcParams['axes.unicode_minus'] = False

步骤2：读取并预处理CSV数据 #

# 读取CSV文件，大文件可通过dtype参数指定字段类型降低内存占用
df = pd.read_csv('你的水位数据文件路径.csv', parse_dates=['datetime'])

# 时间序列一致性校验：统一为15分钟固定间隔，缺失值用线性插值填充
df = df.set_index('datetime').resample('15T').interpolate(method='linear').reset_index()

# 提取水位序列，去除线性趋势避免低频干扰
wl_series = detrend(df['WL (mm)'].values)

步骤3：FFT计算 #

# 采样间隔转换为天单位：15分钟 = 15/(24*60) 天
dt = 15 / (24 * 60)
# 采样频率，单位为 cycles per day（每天采样次数）
Fs = 1 / dt
# 序列长度
N = len(wl_series)

# 执行FFT变换
yf = np.fft.fft(wl_series)
# 计算频率轴，单位 cycles per day
xf = np.fft.fftfreq(N, d=dt)

# 提取正频率段（实信号FFT结果对称，仅保留正频率部分即可）
positive_mask = xf > 0
freq_cpd = xf[positive_mask]
amplitude = 2.0 / N * np.abs(yf[positive_mask])

步骤4：绘制频率分布图 #

plt.figure(figsize=(12, 6))
plt.plot(freq_cpd, amplitude, linewidth=1)

plt.xlabel('频率 (cycles per day)', fontsize=12)
plt.ylabel('振幅 (mm)', fontsize=12)
plt.title('水位时间序列FFT频率分布图', fontsize=14)

plt.grid(alpha=0.3)
# 可选：添加常用潮汐频率参考线
# plt.axvline(x=1, color='orange', linestyle='--', label='全日潮 (1 cpd)')
# plt.axvline(x=2, color='red', linestyle='--', label='半日潮 (2 cpd)')
# plt.legend()

# x轴最大限制到奈奎斯特频率（采样频率的一半，为有效最高频率）
plt.xlim(0, Fs/2)
plt.show()

注意事项 #

• 原始数据噪声较大时可先对水位序列做平滑处理后再执行FFT，得到的频谱图会更干净
• numpy的FFT实现对大数据量优化较好，普通PC处理百万级序列无明显卡顿

内容的提问来源于stack exchange，提问作者Engr M Faysal