目的：
通过编写和运行相关程序，加深对地震子波类型和特征的理解。
内容及要求：
1、编写程序模拟及显示雷克子波，在同一张图上绘制主频为20Hz、35Hz、50Hz的零相位雷克子波波形。
2、分别绘制上述不同频率子波的振幅谱，对比振幅谱最大值对应频率与子波主频关系。

实验原理

时间域雷克子波函数为： ，其中为主频。对雷克子波频谱分析：对不同频率的雷克子波进行傅里叶变换得，所以不同频率的雷克子波振幅谱极大值对应频率即为主频。
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
代码如下
.
.
.
.
.
.
.
1.模拟及显示雷克子波Python代码

/**仅供学习，如有纰漏，望请指正*/1.import numpy as np2.import matplotlib.pyplot as plt3.import scipy4.import math5.6.f1 = 207.f2 = 358.f3 = 509.10.t = np.linspace(-0.055, 0.055, 500)11.12.rick1 = (1-2*(math.pi*f1*t)**2)*math.e**(-1*(math.pi*f1*t)**2)13.14.15.rick2 = (1-2*(math.pi*f2*t)**2)*math.e**(-1*(math.pi*f2*t)**2)16.17.rick3 = (1-2*(math.pi*f3*t)**2)*math.e**(-1*(math.pi*f3*t)**2)18.19.plt.figure(figsize=(10,6))20.fig,ax = plt.subplots()21.ax.plot(t, rick1, '--', label="20Hz")22.ax.plot(t, rick2, label="35Hz")23.ax.plot(t, rick3, ':', label="50Hz")24.ax.set_xlabel("Time")25.ax.set_ylabel("Amp")26.ax.set_title('Zero phase Rayleigh wavelet waveform')27.ax.legend()28.plt.show()

2.上述不同频率子波的振幅谱Python代码

/**仅供学习，如有纰漏，望请指正 --by wxw*/1.import numpy as np2.import matplotlib.pyplot as plt3.import scipy4.import math5.6.f1 = 207.f2 = 358.f3 = 509.10.f = np.linspace(0, 400, 500)11.12.R1 =(2*f**2)/(f1**2*f1*np.sqrt(math.pi))*math.e**(-1*(f/f1)**2)13.R2 =(2*f**2)/(f2**2*f2*np.sqrt(math.pi))*math.e**(-1*(f/f2)**2)14.R3 =(2*f**2)/(f3**2*f3*np.sqrt(math.pi))*math.e**(-1*(f/f3)**2)15.16.plt.figure(figsize=(10,6))17.fig,ax = plt.subplots()18.ax.plot(f, R1, '--', label="20Hz")19.ax.plot(f, R2, label="35Hz")20.ax.plot(f, R3, ':', label="50Hz")21.ax.set_xlabel("Frequency")22.ax.set_ylabel("Amp")23.ax.set_title('Amplitude spectrum')24.ax.legend()25.plt.show()