3.实例：引力波的绘制

引力波

在物理学中，引力波是因为时空弯曲对外以辐射形式传播的能量。

爱因斯坦基于广义相对论预言了引力波的存在，2016年6月16日，LIGO合作组宣布2015年12月26日03：38：53（UTC），两台不同位置的引力波探测器同时探测到了一个引力波信号。

引力波的绘制

数据源：

• http://python123.io/dv/grawave.html
• http://python123.io/dv/H1_Strain.wav
• http://python123.io/dv/L1_Strain.wav
• http://python123.io/dv/wf_template.txt

import numpy as np
import matplotlib.pyplot as plt
from scipy.io import wavfile

rate_h, hstrain= wavfile.read(r"H1_Strain.wav","rb")  # 读取速率rate和数据strain
rate_l, lstrain= wavfile.read(r"L1_Strain.wav","rb")
reftime, ref_H1 = np.genfromtxt('wf_template.txt').transpose()  # reftime是时间序列，ref_H1为信号数据，genfromtxt执行两个循环，第一个循环是将每一行转换为字符串序列，第二个循环是将每一个字符串序列转换成相应的数据类型，transpose将这个两行矩阵进行转置，在赋予给reftime和ref_H1

# 读取应变数据
htime_interval = 1/rate_h  # 求倒数得到波形的时间间隔
ltime_interval = 1/rate_l
fig = plt.figure(figsize=(12, 6))  # 创建一个大小为12x6的绘图空间

# 丢失信号起始点，绘制H1探测器图像
htime_len = hstrain.shape[0]/rate_h  #mshape读取hstrain第一维度的数据的长度，也就是数据点的个数，在除以rate后，得到函数在坐标上的图长度
htime = np.arange(-htime_len/2, htime_len/2 , htime_interval)  # 创建一个数组代表X轴
plth = fig.add_subplot(221)  # 将图像添加到2x2绘图区域的第一个区域上
plth.plot(htime, hstrain, 'y')  # 画出以时间为x轴，应变数据为Y轴的图像。
plth.set_xlabel('Time (seconds)')
plth.set_ylabel('H1 Strain')
plth.set_title('H1 Strain')

# 绘制L1探测器图像
ltime_len = lstrain.shape[0]/rate_l
ltime = np.arange(-ltime_len/2, ltime_len/2 , ltime_interval)
pltl = fig.add_subplot(222)
pltl.plot(ltime, lstrain, 'g')
pltl.set_xlabel('Time (seconds)')
pltl.set_ylabel('L1 Strain')
pltl.set_title('L1 Strain')

# 绘制templates
pltref = fig.add_subplot(212)
pltref.plot(reftime, ref_H1)
pltref.set_xlabel('Time (seconds)')
pltref.set_ylabel('Template Strain')
pltref.set_title('Template')
fig.tight_layout()  # 自动调整图像的四周边缘

plt.savefig("Gravitational_Waves_Original.png")  # 将图像保存为png
plt.show()
plt.close(fig)

LIGO数据的原始封装形式hdf5，这是高效专用于压缩numpy数组的第三方库，LIGO将数据转换成方便下载的文件