bash pip install obspy tensorflow

克隆或安装 PhaseNet 的 Python 实现（例如 github.com/wayneweiqiang/PhaseNet）

### 5.2 加载预训练模型

```python
from phasenet.model import PhaseNet

model = PhaseNet()
model.load_weights('model/phasenet.h5')  # 预训练权重
model.compile(loss='categorical_crossentropy', optimizer='adam')

5.3 读取连续波形并进行预处理

from obspy import read
import numpy as np

st = read('example_waveform.mseed')  # 三分量数据
# 假设采样率为 100 Hz，长度至少 30 秒
st.detrend('linear')
st.filter('bandpass', freqmin=1.0, freqmax=20.0)

# 截取 30 秒窗口
start = st[0].stats.starttime
end = start + 30
st.trim(start, end)

# 合并成 (3, 3000) 数组，单位为采样点
data = np.array([tr.data[:3000] for tr in st])
data = data.T  # 转置为 (3000, 3)
data = data[np.newaxis, :, :, np.newaxis]  # 添加 batch 维度和通道维度

5.4 模型推理与震相拾取

prob = model.predict(data)[0, :, :]  # 输出形状 (3000, 3)
p_prob = prob[:, 0]  # P 波概率
s_prob = prob[:, 1]  # S 波概率

# 设定阈值提取到时
def extract_picks(prob_array, threshold=0.5, min_dist=50):
    picks = []
    above = prob_array > threshold
    # 寻找局部极大值...
    # 简单示例：峰值位置
    from scipy.signal import find_peaks
    peaks, _ = find_peaks(prob_array, height=threshold, distance=min_dist)
    return peaks

p_peaks = extract_picks(p_prob, threshold=0.5)
s_peaks = extract_picks(s_prob, threshold=0.5)

print("P波到时样本点:", p_peaks, "S波到时样本点:", s_peaks)
# 转换为时间
p_times = [str(start + t / 100) for t in p_peaks]  # 100 Hz
s_times = [str(start + t / 100) for t in s_peaks]