企业官网建设流程全解析

简介：这份资源面向地震学研究者、地震台网运维人员及相关专业学生，聚焦利用EMR（经验震级关系）方法估算地震台网的最小完整性震级Mc，为台网监测能力评估与布局优化提供可复用的计算工具。压缩包共13个文件，全部为m脚本，整体约17KB，涵盖核心Mc计算、KS检验与bootstrapping误差评估、震级-频度分布拟合、正态分布与泊松概率辅助函数，以及示例运行脚本，构成一套从数据统计到模型验证的完整流程。已有196人学习下载，说明其在同类工具中具有一定参考价值。读者可借助这些脚本输入自有地震目录，完成震级与观测参数经验关系的建立、Mc阈值推算及不确定性评估，并据此判断台网对低震级事件的探测能力，适合作为地震监测能力分析与台站优化研究的基础代码框架。

1. 从一张台网震级底图说起：EMR 方法到底在算什么

如果你手头有一份区域地震台网的观测报告，想回答“这套台网在本地到底能监测到多小的地震”，最直接的办法不是翻规范，而是把 EMR 方法跑一遍。EMR 是“震级—距离”经验关系（Empirical Magnitude–Range）的缩写，核心思路很朴素：用台网历史记录中每个台站实际能检测到的最小震级，按震中距分档统计，再拟合出一条随距离衰减的监测能力曲线。它不依赖理论噪声模型，也不要求你知道每个台基的场地响应，只吃“已经记录到的目录”这一份数据。对做地震台网运维、台址勘选、监测能力评估的人来说，这套方法最大的价值是：结论直接对应“这套台网现在能报出几级地震”，而不是停留在信噪比公式里。标题里的“EMR-example”就是一份可复现的最小算例，把台网目录、震级、震中距三样东西串起来，输出一张监测能力底图。适合谁？适合手里有台网目录、想快速评估监测下限、又不想从零推导检测概率的从业者。

2. EMR 方法的输入准备：目录、震级与震中距怎么对齐

2.1 为什么 EMR 只认“已记录到”的震相

EMR 的统计对象是台网目录里每一条“某台站对某次地震有震相拾取”的记录。它不关心这次地震有没有被定位，只关心这个台站有没有在这个震中距上留下可用的震级。常见做法是：从台网观测报告中导出“台站—地震”配对表，字段至少包含台站代码、震中距、震级、震相类型。这里有个容易翻车的点：很多人直接拿地震目录的震级去配对，但目录震级是台网平均震级，不是单台震级。EMR 要的是单台震级，否则距离分档后震级会被平均掉，曲线整体偏低。我一般会从震相报告里取单台震级，如果只有目录震级，就退而求其次用“台站参与定位”的记录，但要在文档里注明这个近似。

2.2 震中距分档：别用等间距，用等对数

震中距从几公里到几百公里跨度很大，等间距分档会让近距档样本爆炸、远距档样本稀疏。常见做法是按对数等间距分档，比如 0–10 km、10–20 km、20–50 km、50–100 km、100–200 km、200–500 km。每档内统计该档所有记录的最小震级，或者取第 5 百分位震级作为“可检测下限”。为什么不用最小值？因为单条异常记录会把下限拉低，第 5 百分位更稳。下面这段 Python 就是做分档和下限统计的最小骨架。

import numpy as np import pandas as pd # 读取台站-地震配对表，字段：station, dist_km, mag, phase df = pd.read_csv("emr_pairs.csv") # 对数等间距分档边界 bins = [0, 10, 20, 50, 100, 200, 500] labels = ["0-10", "10-20", "20-50", "50-100", "100-200", "200-500"] df["dist_bin"] = pd.cut(df["dist_km"], bins=bins, labels=labels, right=False) # 每档取第5百分位震级作为检测下限 result = df.groupby("dist_bin")["mag"].quantile(0.05).reset_index() result.columns = ["dist_bin", "mag_lower"] print(result)

逻辑说明：pd.cut按给定边界分档，right=False表示左闭右开，避免边界值重复计入。quantile(0.05)取第 5 百分位，比最小值抗异常。参数说明：bins要根据你台网的实际震中距分布调整，如果台网最远记录只有 300 km，最后一档写到 500 没意义，改成 300。phase字段可以用来过滤，只保留 Pg/Sg 这类定位震相，避免面波震级混进来。

2.3 单台震级缺失时怎么补

如果震相报告里只有振幅和周期，没有直接给单台震级，就得自己算。常见做法是用地方震震级公式：M = log(A) + 1.11*log(R) + 0.00189*R - 2.09，其中 A 是最大振幅（微米），R 是震中距（公里）。这个公式是华南地区的经验式，换区域要换系数。补完之后一定要做一致性检查：把算出来的单台震级和目录震级做散点，如果系统偏差超过 0.3 级，说明公式系数不匹配，得换。这一步没有捷径，血泪经验是：宁可少用几条记录，也不要用错公式把整条曲线带偏。

3. 用 EMR-example 跑通监测能力曲线：拟合、绘图与参数

3.1 最小可复现流程：从配对表到能力曲线

拿到分档下限后，下一步是拟合一条“震级随距离变化”的曲线。EMR 常用线性拟合：M(R) = a * log10(R) + b，其中 R 是震中距，a 和 b 是拟合系数。为什么用 log10(R)？因为震级本身是对数量纲，距离取对数后线性关系更稳。下面这段代码完成拟合并输出系数。

import numpy as np from scipy.optimize import curve_fit # 取分档中心距离 bin_centers = np.array([5, 15, 35, 75, 150, 350]) mag_lower = result["mag_lower"].values def emr_func(r, a, b): return a * np.log10(r) + b popt, pcov = curve_fit(emr_func, bin_centers, mag_lower) a, b = popt print(f"拟合系数: a={a:.3f}, b={b:.3f}") # 预测 50 km 处的监测下限 print(f"50 km 监测下限: {emr_func(50, a, b):.2f}")

逻辑说明：curve_fit用最小二乘拟合，popt是系数，pcov是协方差。参数说明：bin_centers要和你实际分档一致，如果分档是 0–10、10–20，中心取 5 和 15 没问题；但如果分档是对数等间距，中心应该取几何平均而不是算术平均。拟合完看pcov对角线，如果某个系数方差特别大，说明该档样本太少，考虑合并相邻档。

3.2 绘图时最容易忽略的三个参数

画监测能力图时，横轴是震中距（对数），纵轴是震级。三个参数决定图能不能用：第一，横轴范围要覆盖台网实际监测范围，别默认 0–500 km，如果台网最远 200 km，画到 500 就是空白；第二，纵轴下限要低于最小震级 0.5 级，否则曲线贴底看不清；第三，分档点要画出来，不能只画拟合线，否则审图的人不知道你用了多少数据。下面这段是绘图骨架。

import matplotlib.pyplot as plt plt.figure(figsize=(8, 5)) plt.scatter(bin_centers, mag_lower, color="red", label="分档下限") r_smooth = np.logspace(np.log10(5), np.log10(350), 100) plt.plot(r_smooth, emr_func(r_smooth, a, b), "b-", label="EMR拟合") plt.xscale("log") plt.xlabel("震中距 (km)") plt.ylabel("震级") plt.ylim(mag_lower.min() - 0.5, mag_lower.max() + 0.5) plt.legend() plt.grid(True, which="both", ls="--") plt.savefig("emr_capability.png", dpi=300)

逻辑说明：logspace生成对数等间距的平滑横轴，避免折线感。which="both"让主次网格都显示，方便读对数坐标。参数说明：dpi=300是投稿或报告用的分辨率，如果只是内部看，150 就够。ylim下限减 0.5 是为了让最低点不贴边。

3.3 拟合系数的物理含义与边界

a 值反映监测能力随距离衰减的快慢，a 越大，远距离能监测到的震级越高，说明台网在远距越吃力。b 值是 1 km 处的截距，但实际没有 1 km 的记录，所以 b 只作为拟合参数，不要单独解释。常见范围：区域台网 a 在 1.0–2.0 之间，b 在 -1.0–0.5 之间。如果 a 超过 2.5，检查是不是远距档样本太少导致拟合被拉偏；如果 a 小于 0.5，检查是不是近距档震级普遍偏高，把曲线压平了。这些边界值不是规范，是我跑了几十个台网总结出来的经验区间，超出就要回头查数据。

4. 避坑与排查：EMR 算例里最容易翻车的五件事

4.1 现象：曲线整体比预期低 0.5 级 → 原因：用了目录震级而非单台震级 → 解决：回震相报告取单台震级

这是最常见的翻车。目录震级是多个台站平均后的结果，天然比单台震级平滑且偏低。如果你直接拿目录震级配对，分档下限会被拉低，曲线整体下移。解决方法是回到震相报告，取每个台站自己的震级。如果震相报告里没有单台震级，用振幅周期自己算，但要用对区域公式。

4.2 现象：远距档拟合线翘尾 → 原因：远距档样本太少，被个别大震级记录主导 → 解决：合并档位或设最小样本数

远距档往往只有几条记录，如果这几条恰好都是较大地震，下限就偏高，拟合线在远端上翘。解决方法是设一个最小样本数，比如每档至少 10 条记录，不够就和相邻档合并。合并后档中心要重新算，别直接用原来的。

4.3 现象：近距档下限异常低 → 原因：混入了非定位震相或爆破记录 → 解决：按震相类型和事件类型过滤

近距档如果混入爆破、塌陷或非定位震相，震级会异常低，把下限拉下去。过滤方法是：只保留 Pg/Sg 震相，事件类型只保留天然地震，爆破和塌陷单独处理。如果目录里没有事件类型字段，用震源深度和波形特征辅助判断，但这一步比较费时，建议在数据准备阶段就做好。

4.4 现象：拟合不收敛 → 原因：分档中心距离取了算术平均而非几何平均 → 解决：对数分档用几何中心

对数等间距分档时，档内距离不是均匀分布，算术平均会偏向大值。比如 10–20 km 档，算术平均是 15，但几何平均是 14.1。差别不大，但如果所有档都偏，拟合就会不收敛。解决方法是统一用几何平均：sqrt(下限*上限)。

4.5 现象：换区域后曲线完全不可用 → 原因：震级公式系数没换 → 解决：按区域换地方震震级公式

地方震震级公式是区域性的，华南的系数拿到西北用，系统偏差可能超过 0.5 级。换区域时，先找该区域已发表的地方震震级公式，如果没有，用该区域台网的历史记录重新拟合系数。这一步没有后悔药，只能老老实实做。

5. 进阶技巧：用 EMR 曲线反推台站增补优先级

跑完 EMR 曲线只是第一步，真正有用的是拿它做台站布局决策。具体做法：把现有台网的 EMR 曲线画出来，再模拟“如果某个方位增补一个台站，该方位的监测下限能降多少”。常见做法是按方位角分扇区，每个扇区单独跑 EMR，找出监测能力最弱的扇区，优先增补。下面这个表格是我常用的扇区评估模板，填完就能看出哪个方位最需要补台。

填表时注意：记录数少于 50 的扇区，下限置信度低，先补数据再下结论。达标标准按你的台网任务定，如果任务是监测 1.5 级以上，50 km 下限低于 1.5 才算达标。增补优先级按“下限差值 × 扇区面积”排序，差值越大、面积越大，越优先。最后说一个我自己的习惯：每次跑完 EMR，我都会把分档下限和拟合系数存成 CSV，下次换区域时直接对比，看是数据问题还是方法问题。这个习惯帮我省了很多重复排查的时间。希望帮到你。

本文还有配套的精品资源，点击获取