用python写一段JONSWAP 谱的代码，输出功率谱密度和波浪角频率

在Python中，我们可以使用matplotlib库来绘制功率谱密度（Power Spectral Density, PSD）图，而JONSWAP谱是一种描述海洋波浪能量分布的常用模型。下面是一个简单的示例，展示如何计算和绘制JONSWAP谱：

import numpy as np
import matplotlib.pyplot as plt

# JONSWAP谱参数（可以根据实际情况调整）
beta = 3.3  # 非线性系数
gamma0 = 0.7  # 最大波高比例
T_p = 8  # 平均周期 (秒)
Hs = 4  # 破碎波高 (米)
g = 9.81  # 重力加速度

# 计算谱宽和谱峰频率
f_peak = 1 / T_p  # 波峰频率
sigma = beta * g * Hs**2 / ((2 * gamma0 * f_peak)**5)  # 谱宽

# 创建波数范围
num_freq = 1000
freq_range = np.linspace(0, 1/f_peak, num_freq)

# 计算JONSWAP谱
SpectralDensity = 1 / (np.sqrt(2 * np.pi) * sigma) * np.exp(-((freq_range - f_peak) / sigma)**2) * (freq_range**5)

# 绘制PSD
plt.figure(figsize=(10, 6))
plt.plot(freq_range, SpectralDensity, label='JONSWAP Power Spectrum')
plt.xlabel('Frequency (Hz)')
plt.ylabel('PSD (m$^2$/s)')
plt.title('JONSWAP Wave Power Spectral Density')
plt.grid(True)
plt.legend()
plt.show()