import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
import pywt
from scipy.signal import find_peaks
# 读取数据
file_path = r'D:\wave\precipitation_A.xlsx'
df = pd.read_excel(file_path)
# 提取年份和降水数据
year = df['year'].values
precipitation = df['P'].values
# Morlet小波变换
wavelet = 'cmor'
scales = np.arange(1, 128) # 定义尺度范围
# 进行小波变换
coefficients, frequencies = pywt.cwt(precipitation, scales, wavelet)
# 计算能量
power = (np.abs(coefficients))**2
# 计算平均能量
avg_power = np.mean(power, axis=0)
# 找到频谱中的峰值
peaks, _ = find_peaks(avg_power)
# 绘制Morlet小波分析结果
plt.figure(figsize=(12, 6))
plt.subplot(2, 1, 1)
plt.plot(year, precipitation, label='降水数据')
plt.title('降水数据及Morlet小波变换')
plt.legend()
plt.subplot(2, 1, 2)
plt.plot(frequencies, avg_power, label='平均能量')
plt.plot(frequencies[peaks], avg_power[peaks], 'ro', label='峰值')
plt.title('Morlet小波变换能量谱')
plt.legend()
plt.show()