简介：面向地震勘探与工程地质评价领域技术人员及地球物理学习者的面波处理与剖面连接软件包，围绕Love波、Rayleigh波数据预处理、相位解缠、频域与时频分析、面波成像，以及速度模型构建、层析反演、多测线剖面连接和模型校验等环节，提供从观测数据到地下速度剖面的完整工具链，适合具备一定地震资料处理基础、希望系统练习浅层结构反演的用户。压缩包共91个文件，约1.8MB，以gif、htm、bmp等图文说明与界面素材为主，辅以exe可执行程序、sys/vxd驱动、doc使用说明及少量h、dll、lib等开发与加密狗组件，兼顾软件运行、操作参考和示例展示。已有923人学习下载，说明其在面波处理与剖面连接方向具有较高参考价值。使用者可借此熟悉CCSWSWIN与CCSWSMAP的模块组织、参数设置和成果表达，理解噪声去除、反演稳定性与模型分辨率等实际问题，并对照图文资料完成速度剖面解释与地质特征判读。

1. 一条测线跑完，活儿才刚开始：面波处理与剖面连接到底卡在哪

野外收工回来，U 盘里躺着几十个道集文件。单条测线的面波处理其实不难：滤波、提频散、反演，一条龙下来半小时能出一条 Vs 剖面。真正花时间的是剖面连接——相邻测线桩号对不上、深度网格不一致、色标各调各的，拼出来一张图自己都不敢往报告里放。这套面波处理及剖面连接软件对准的就是这个环节：把频散提取、横波速度反演、多测线里程对齐、剖面拼接成图塞进一个桌面工作流，让"几十条测线批量出图"变成点几下的事。它适合岩土勘察、地基与路基检测这类要成图交付的团队，也适合想拿现成流程二次开发的人。

2. 技术选型与整体架构：一个桌面处理工具该长什么样

2.1 从道集文件到拼接剖面：五段式处理链路

这套工具的处理链路拆开看是五段，段与段之间靠中间文件解耦，不是一条黑箱流水线。第一段是道集读入，把 SEG-2、SEG-Y、npy 这些来源不同的文件统一成"道集对象"——一个二维数组加采样间隔、道间距、最小偏移距、桩号四个元信息。第二段是预处理，去直流、去线性趋势、端点加窗、带通滤波，目的是让后面的频谱干净。第三段是频散能量图生成，把时距域数据变到频率-相速度域，图上能量脊的位置就是频散曲线。第四段是拾取与反演，沿着能量脊取点得到实测频散曲线，再反演成分层的 Vs-深度模型。第五段是剖面连接，把每条测线的 Vs 剖面按里程摆到同一坐标系里，插值成规则网格，配色出图。

之所以强调"解耦"，是因为这条链上每一步都可能被反复调。频散拾取不满意要回头改带通范围，反演残差大要回头检查能量图，剖面拼不上要回头核对桩号。如果每改一个参数就得把整条链重跑一遍，一条测线十分钟，几十条就是几小时。工程上常见的做法是把每段的结果落盘——能量图存 npz、频散曲线存 csv、Vs 模型存 json，改下游的时候直接读上游文件，不重算。这套工具的项目结构基本就是按这个思路切的，后面讲目录时会看到core/下每个模块都带自己的save/load。

链路里最容易被低估的是第二段。新手往往直接把原始道集丢进 f-k 变换，出来的能量图糊成一团，然后怀疑算法有问题。实际上多数情况是某几道坏道在拖后腿，或者近偏移距道的近场效应把低频段能量带偏了。预处理阶段留一个道编辑入口——手动剔除坏道、按偏移距范围筛选参与计算的道——比事后在能量图上硬描曲线要靠谱得多。

第五段的输入输出格式也值得先说清楚。剖面连接要的不是"图片"，而是带坐标的散点：每个点有里程、深度、Vs 三个值。工具内部通常把单条测线的反演结果存成三列数组，拼接模块只认这个格式。这样一来，反演算法换掉、或者从别的软件导进来一批结果，只要凑成三列数组就能进拼接流程。这个设计看着不起眼，二次开发时价值很大——你可以只换反演内核，拼接和成图原样复用。

2.2 技术栈取舍：Python + Qt 桌面栈为什么更适合这类工具

这类地球物理处理工具有两条常见技术路线。一条是 C# WinForms 或 WPF 加商业图表控件，界面做得漂亮，部署也简单，缺点是算法生态薄，做矩阵运算、优化反演得自己啃或者买库。另一条是 Python 加 Qt（PySide6 或 PyQt），界面稍显朴素，但 NumPy、SciPy、Matplotlib 这套生态把数值计算和成图的活全接了，反演直接调优化器，出图直接上 colormap。

我一般会选后者，理由很实际：面波处理的核心计算量在频散能量图（每个频率点沿相速度扫描一遍）和反演（优化器要调几百上千次正演），这两块用 NumPy 向量化和 SciPy 的least_squares、differential_evolution写起来最省事。界面部分反而是次要的，Qt 的QMainWindow加几个pyqtgraph绘图控件足够，能量图需要交互拾取，pyqtgraph的鼠标事件和缩放性能比 Matplotlib 内嵌好一截。

依赖清单通常长这样，版本按你环境里能装上的最新稳定版走，不用死磕：

# requirements.txt 节选 numpy>=1.24 scipy>=1.10 matplotlib>=3.7 pyqtgraph>=0.13 PySide6>=6.5 disba>=0.6 # Rayleigh/Love 波频散正演 h5py>=3.9 # 大数据量道集缓存 pandas>=2.0 # 频散曲线与剖面结果的表格读写

安装和启动就是常规三步，虚拟环境里跑：

python -m venv .venv source .venv/bin/activate # Windows 下是 .venv\Scripts\activate pip install -r requirements.txt python main.py

disba这个库值得单独说一句，它实现了快速标量刚度矩阵法计算层状模型的 Rayleigh 波频散，反演时做正演比从零写快一个量级，而且支持基阶和高阶模式。如果你的项目需要 Love 波或者要算敏感度核，它也能接。不装它的话反演部分就得自己实现正演，工作量翻几倍，不划算。

界面层的组织建议不要把所有窗口堆在一个文件里。主窗口只负责菜单、工具栏和布局，能量图、剖面图、模型表各做成独立 widget，通过信号槽和主窗口通信。Qt 的信号槽机制在这种"一处改动、多处联动"的场景里很好用——拾取点变了发个信号，模型表和图同时刷新，不用手动调刷新函数。

2.3 目录结构：数据层、算法层、界面层怎么切

工程目录按职责分层，算法层不 import 任何 Qt 的东西，这样算法可以脱离界面单独跑批处理脚本：

masw_tool/ ├── main.py # 入口，创建 QApplication 和主窗口 ├── configs/ │ └── default.yaml # 滤波范围、速度扫描区间等默认参数 ├── core/ # 纯算法层，无 UI 依赖 │ ├── gather.py # 道集对象与格式分发读入 │ ├── preprocess.py # 去趋势、加窗、带通 │ ├── dispersion.py # 相移法/f-k 能量图 │ ├── picking.py # 频散曲线拾取与平滑 │ ├── invert.py # 反演封装 │ └── merge.py # 剖面连接与插值 ├── io_utils/ │ ├── seg2.py # SEG-2 解析 │ └── project.py # 工程文件（json）读写 ├── ui/ │ ├── main_window.py │ └── widgets/ │ ├── gather_view.py # 波形显示 │ ├── disp_view.py # 能量图 + 拾取交互 │ └── section_view.py # 拼接剖面显示 └── batch/ └── run_all.py # 无界面批处理，读工程文件批量出图

core/和ui/之间只通过数据对象和信号通信，这条边界守住，后面加批处理、加 Web 服务、换界面框架都不受影响。batch/run_all.py是我特别建议保留的——现场数据动辄几十上百条测线，一条条点界面不现实，写个脚本遍历工程目录、按默认参数跑完、把结果和缩略图落到指定文件夹，第二天上班看结果就行。

configs/default.yaml把参数外置，也是个血泪经验。参数写死在代码里，换个工区就得改代码重新打包；放配置文件里，现场调完存一份新配置，回来直接复用。反演的分层数、速度上下界、滤波频带这些参数，每个地区都不太一样，配置文件是刚需。

3. 核心模块拆解：频散曲线提取与反演的实现链路

3.1 道集读入与预处理：先统一成一种对象

整套流程的第一道关是把五花八门的文件格式归一。SEG-2 是面波仪最常见的输出，SEG-Y 偶尔遇到，试验阶段可能还有 npy 或 csv。工程上的做法是定义一个Gather数据类，所有读入函数最终都返回它，下游代码只认这个类型。

from dataclasses import dataclass import numpy as np @dataclass class Gather: data: np.ndarray # (nch, ns) float32，按道优先 dt: float # 采样间隔，单位 s dx: float # 道间距，单位 m offset0: float # 最小偏移距，单位 m chainage: float = 0.0 # 该测线起点桩号，单位 m @property def offsets(self): # 各道相对震源的偏移距 return self.offset0 + np.arange(self.data.shape[0]) * self.dx

offset0和dx这两个参数不是装饰品，后面算相移、算理论频散都要用。野外记录里这两个值通常写在道头或班报上，读 SEG-2 的时候顺手取出来；实在取不到就让用户在界面上填，填错了能量图会整体偏移，这在避坑章节里还会碰到。chainage是给剖面连接预留的，单条测线处理时用不上，拼图时必须。

预处理函数把去直流、去线性趋势、端点加窗、带通滤波串起来：

def preprocess(g: Gather, fmin=5.0, fmax=80.0, taper=0.05): x = g.data.astype(np.float64).copy() x -= x.mean(axis=1, keepdims=True) # 去直流分量 # 逐道去线性趋势，压制长周期漂移 t = np.arange(x.shape[1], dtype=np.float64) for i in range(x.shape[0]): p = np.polyfit(t, x[i], 1) x[i] -= np.polyval(p, t) # 两端余弦窗，减小 FFT 泄漏 n = x.shape[1] m = max(1, int(n * taper)) w = np.ones(n) ramp = 0.5 - 0.5 * np.cos(np.pi * np.arange(m) / m) w[:m] = ramp w[-m:] = ramp[::-1] x *= w # 频域带通：直接置零频带外的谱 X = np.fft.rfft(x, axis=1) f = np.fft.rfftfreq(n, d=g.dt) band = (f >= fmin) & (f <= fmax) X[:, ~band] = 0.0 return np.fft.irfft(X, n=n, axis=1)

fmin和fmax的选择直接决定反演深度范围，这是调参的主战场。经验上低频端受检波器主频和测线长度限制，4.5 Hz 检波器配几十米偏移距，可靠下限大概在 5 到 8 Hz；高频端受道间距和空间假频限制，道间距 1 m 时 60 到 80 Hz 就到顶了。taper取 0.05 是通用值，数据两端有明显突变时可以加到 0.1，代价是有效信号也损失一点。

带通这里用的是"硬置零"，频率响应是矩形，时域会有振铃。要求高的场合可以换成 Butterworth 零相位滤波（scipy.signal.filtfilt），代价是要处理相位补偿。我一般先用硬置零快速看效果，定了频带再换正式滤波。

3.2 相移法生成频散能量图：能量脊就是频散曲线

从道集到频散曲线，中间这步是整个工具的技术核心。常用两种做法：二维 FFT 得到 f-k 谱，或者相移法（phase shift）直接扫相速度。f-k 谱的横轴是波数，要换算成相速度还得除以频率，分辨率受测线长度限制；相移法直接把相速度作为扫描变量，出的图横轴就是频率、纵轴就是相速度，拾取更直观。这套工具用的是后者。

原理不复杂：对某个频率 f，如果有一列波以相速度 c 沿测线传播，那么相邻道之间它的相位差是 2πf·dx/c。把每道在该频率的复频谱乘上对应的相位校正因子再叠加，如果 c 猜对了，各道同相叠加、幅值最大；猜错了就互相抵消。对每个频率扫一遍相速度，得到一张能量图。

def dispersion_image(g: Gather, data: np.ndarray, freqs=None, cs=None): """ 相移法频散能量图 g : 道集元信息（提供 dx、dt） data : (nch, ns) 预处理后的波形 freqs : 频率扫描轴，Hz；None 时取 FFT 频率轴 cs : 相速度扫描轴，m/s """ nch, ns = data.shape if freqs is None: freqs = np.fft.rfftfreq(ns, d=g.dt) if cs is None: cs = np.linspace(50.0, 1200.0, 240) spec = np.fft.rfft(data, axis=1) # (nch, nf) 复频谱 # 按道归一化振幅，避免强能量道主导叠加结果 amp = np.abs(spec) amp[amp < 1e-12] = 1e-12 spec_n = spec / amp x = g.offsets[:, None] # (nch, 1) 各道偏移距 energy = np.zeros((len(freqs), len(cs)), dtype=np.float64) for j, c in enumerate(cs): # 沿偏移距的相位校正：猜对 c 时各道同相 phase = np.exp(-1j * 2.0 * np.pi * freqs[None, :] * x / c) stacked = np.sum(spec_n * phase, axis=0) / nch energy[:, j] = np.abs(stacked) return freqs, cs, energy

cs的扫描范围要覆盖场地可能的相速度，土层场地 50 到 1200 m/s 基本够用，岩石场地往上加。扫描点数 240 是分辨率和耗时的折中，点数太少能量脊会断。spec_n那步振幅归一化很关键——野外数据近道能量往往比远道强一个量级，不归一化的话叠加结果被近道绑架，低频段能量脊会歪。

能量图出来后，横轴频率、纵轴相速度，真正的频散曲线就是图上那条最亮的脊。拾取的做法是在每个频率切片上找能量峰值，再加一个人工修正的环节——自动拾取在低频段经常跳到高阶模式或者被噪声带偏，得让用户在图上点几下纠正。拾取结果存成(freq, phase_vel)两列，后面反演用。

还有一个参数值得提：参与计算的道数。近偏移距道受近场效应影响，远偏移距道信噪比低，实际参与叠加的道通常取中间一段。界面上给个道号范围选择，比全局改参数方便。

3.3 频散曲线反演：把相速度变成 Vs 剖面

反演是这套工具里最"玄学"的一步，也是最容易被误用的一步。输入是一条实测频散曲线（频率 vs 相速度），输出是分层模型的横波速度结构。本质是个非线性最小二乘问题：调整各层的 Vs 和厚度，让正演出来的理论频散曲线逼近实测曲线。

正演用disba，反演用 SciPy 的least_squares，参数化方式是把模型写成"各层 Vs + 各层厚度"的向量：

import numpy as np from disba import PhaseDispersion from scipy.optimize import least_squares def build_model(params, n_layers): """params = [vs_1..vs_n, h_1..h_{n-1}]，单位均为 m/s 和 m""" vs = params[:n_layers] h = np.abs(params[n_layers:]) # 厚度取正 vp = vs * 1.85 # 泊松比约 0.29 的经验换算 rho = np.full(n_layers, 2.0) # 密度取常数 g/cm3 return h, vp, vs, rho def residual(params, freqs, c_obs, n_layers): h, vp, vs, rho = build_model(params, n_layers) # disba 的单位是 km 和 km/s，这里统一换算 thickness = np.append(h, 0.0) / 1000.0 model = (vp / 1000.0, vs / 1000.0, rho) try: pd = PhaseDispersion(thickness, model) period = 1.0 / freqs _, c_pred = pd(period, mode=0, wave="rayleigh") except Exception: return np.full_like(freqs, 1e3) # 正演失败给大残差 return (c_pred * 1000.0 - c_obs) / c_obs # 相对残差

build_model里vp = vs * 1.85是泊松比约 0.29 时的经验关系，饱和土层泊松比接近 0.5、比值能到 3 以上，干燥土层比值小一些。这个系数在缺少纵波资料时可以当固定值，有 PS 测井或者折射资料的话直接代入更准。密度对 Rayleigh 波频散的影响相对小，取常数 2.0 g/cm³ 一般够用。

反演调用就是给初值和边界：

def invert(freqs, c_obs, n_layers=5): # 初值：Vs 从 150 到 600 线性分布，厚度从浅到深递增 vs0 = np.linspace(150.0, 600.0, n_layers) h0 = np.linspace(2.0, 8.0, n_layers - 1) x0 = np.concatenate([vs0, h0]) lower = np.concatenate([np.full(n_layers, 80.0), np.full(n_layers - 1, 0.5)]) upper = np.concatenate([np.full(n_layers, 1200.0), np.full(n_layers - 1, 30.0)]) res = least_squares(residual, x0, bounds=(lower, upper), args=(freqs, c_obs, n_layers), max_nfev=200, x_scale="jac") return res.x, res.cost, res.optimality

n_layers不是越多越好。层数多了反演非唯一性放大，同一套数据能拟出好几种速度结构，残差都很小但剖面形态差很远。4 到 6 层是常见区间，具体看场地地层复杂程度。max_nfev是正演次数上限，200 次对 5 层模型通常够收敛，层数加到 8 以上适当调高。

res.cost和res.optimality两个返回值一定要看，不能只看曲线拟合得好不好。cost是残差平方和的一半，量级应该在几十以内；optimality反映收敛程度，数值大说明还在爬坡。拟合图上也别只看曲线重合度，反演出的模型如果出现 Vs 随深度剧烈跳变、或者层厚趋近边界值，多数是层数给多了或者频带不够宽导致的，得回去改模型参数而不是继续调优化器。

反演的深度分辨率有个物理上限：最大探测深度大致是最大波长的一半到三分之一，而最大波长受最低可靠频率限制。5 Hz 的 Rayleigh 波在 300 m/s 的介质里波长 60 m，探测深度也就 20 到 30 m。指望 5 Hz 的下限去反演 50 米深的基岩面，出来的结果不可信。这个换算关系建议直接做进界面提示，用户填了频带就自动算个预估深度显示出来，比事后解释省事。

4. 剖面连接模块：多测线拼图的坐标系与插值

4.1 里程-深度坐标系：拼图前先把尺子统一

剖面连接听上去就是"把几张图拼一起"，实际上难点全在坐标系。每条测线有自己的起点桩号、自己的方向、自己的测线长度，反演出来的 Vs 点云也是按各自测线的相对位置排的。不把里程统一，拼出来必然是错位的。

工程上的做法是维护一个测线注册表，每条测线登记三项：起点绝对里程、终点绝对里程、方向（正向或反向）。单条测线反演产出的三列数组(沿测线距离, 深度, Vs)里，第一列是相对距离，连接前先换算成绝对里程：

def to_absolute_chainage(dist, start_chainage, direction=1): """dist: 沿测线相对距离（m），direction=1 正向，-1 反向""" if direction == 1: return start_chainage + dist else: # 反向测线：起点桩号对应最大里程 return start_chainage - dist

方向判断不能靠猜。野外布线经常来回折返，同一条路两侧各布一条线、方向相反很常见。注册表里显式记方向，比事后看图发现镜像了再返工强。start_chainage用绝对桩号（比如道路 K3+250 记成 3250.0），同一工区所有测线共用一个里程基准，这是拼图能对上的前提。

里程统一之后还要处理重叠。相邻测线在接头处往往有几十米重叠段，两套反演结果在同一里程上 Vs 值不一样，直接拼会有台阶。常见做法是重叠段按距离加权平均——靠近本测线中心权重高，靠近端点权重低：

def blend_overlap(seg_a, seg_b, overlap_start, overlap_end): """seg_a/seg_b: (chainage, depth, vs) 三元组，在重叠段线性过渡""" ... w = (ch - overlap_start) / (overlap_end - overlap_start) vs_blend = (1 - w) * vs_a + w * vs_b

权重函数用线性过渡最省事，重叠段长的时候可以换余弦窗，过渡更平滑。这一步不做，拼接剖面上会出现明显的竖向条带，报告里很难解释。

4.2 网格插值与空白区补齐

散点拼成图要落到规则网格上，两个方向的网格间距得统一。沿测线方向按里程布点，深度方向按固定步长布点，然后对每个网格节点做插值。SciPy 的griddata能直接干这活：

import numpy as np from scipy.interpolate import griddata def build_section_grid(sections, dx=1.0, dz=0.5): """ sections: list of dict，每项含 'chainage'、'depth'、'vs' 三个等长数组 dx, dz : 沿测线和深度方向的网格步长，单位 m """ pts = np.concatenate([np.column_stack([s["chainage"], s["depth"]]) for s in sections]) vals = np.concatenate([s["vs"] for s in sections]) xg = np.arange(pts[:, 0].min(), pts[:, 0].max() + dx, dx) zg = np.arange(pts[:, 1].min(), pts[:, 1].max() + dz, dz) X, Z = np.meshgrid(xg, zg) # 先用线性插值保住细节，再用最近邻补边角空白 V = griddata(pts, vals, (X, Z), method="linear") nan_mask = np.isnan(V) if nan_mask.any(): V[nan_mask] = griddata(pts, vals, (X[nan_mask], Z[nan_mask]), method="nearest") return X, Z, V

dx和dz的选择要看原始数据密度。反演出来的 Vs 点沿深度方向通常每米一个点，dz取 0.5 m 是合理的内插；dx取决于测线间距和沿测线的采样密度，取 1 m 比较稳妥。网格取得比原始数据还粗会丢细节，取得太细只是把插值误差放大，没有实际信息增益。

插值方法这里有个坑：linear在凸包边缘会留 NaN，nearest补出来的边角颜色会跟邻近区域不一致，看着像一块补丁。如果边角区域大，说明测线覆盖本身有缺口，与其硬补，不如把空白区留白并在图上标注"无数据"。这一点在交付报告里反而更专业。

剖面深度方向还要处理一个统一问题。各条测线反演的分层厚度不一样，最大深度也不一样，拼图时深度轴得取所有测线的最大深度，浅的测线在深部留空。这个逻辑放在网格构建阶段做，比在绘图阶段补要干净。

4.3 渲染、色标与导出

拼好的网格用pcolormesh出图，色标选 Vs 的经典彩虹或者分段色阶都行，关键是整个工区所有剖面图用同一套色标范围和分级，否则单看每张都正常，摆一起就没法对比。色标范围建议取全工区 Vs 的 5% 和 95% 分位数，避免个别异常值把色标拉满。

import matplotlib.pyplot as plt from matplotlib.colors import BoundaryNorm def plot_section(X, Z, V, out_png, vmin=100, vmax=800, step=50): levels = np.arange(vmin, vmax + step, step) cmap = plt.get_cmap("jet", len(levels) - 1) norm = BoundaryNorm(levels, cmap.N) fig, ax = plt.subplots(figsize=(12, 4), dpi=200) m = ax.pcolormesh(X, Z, V, cmap=cmap, norm=norm, shading="auto") ax.invert_yaxis() # 深度向下增大 ax.set_xlabel("里程 (m)") ax.set_ylabel("深度 (m)") fig.colorbar(m, ax=ax, label="Vs (m/s)", ticks=levels[::2]) fig.tight_layout() fig.savefig(out_png) plt.close(fig)

shading="auto"让 pcolormesh 按网格中心对齐，剖面不会整体偏移半格。invert_yaxis是 Earth science 的惯例，深度向下，不做这一步图是上下颠倒的，评审一眼就看出来。

导出除了 PNG，工程上通常还要矢量图和原始网格。PDF 或 SVG 给报告排版用，网格数据存成 csv 或 npz 方便后续在其他软件里加工。这几个导出入口都做在界面上，比每次改代码导要实际。

5. 避坑与排查：几个把结果带偏的细节

5.1 数据与格式类问题

现象：频散能量图整体糊成一片，看不到清晰的能量脊。原因通常有两个。一是道间距dx或采样间隔dt填错了，相移计算里的相位校正因子直接跟这两个量挂钩，错一个数量级能量就对不齐；二是原始道集里有几道明显坏道（死道、极性反转、工频干扰严重），叠加时把有效信号压下去了。解决办法是先核对道头或班报里的采集参数，用单道波形图快速扫一眼有没有异常道，有的话在参与计算的道范围里剔掉。极性反转的道肉眼不容易发现，可以算一下各道与参考道的互相关，相关系数为负的翻极性。

现象：反演出的 Vs 值比预期高或低一个数量级，但曲线拟合看着还挺好。这是单位换算的经典翻车。disba内部用 km 和 km/s，输入是 m 和 m/s，差 1000 倍。相对残差的形式让拟合优度对整体缩放不敏感，所以曲线形状对得上，绝对值全错。解决方式是在build_model里统一除以 1000，并且加一个断言：如果反演出的 Vs 均值超过 5000 m/s 或者低于 50 m/s，直接报警而不是静默返回。

现象：相邻测线拼出来的剖面在接头处错开几十米。原因多半是里程基准不统一——有的测线用相对桩号（从 0 开始），有的用绝对桩号（从路段起点算）。核对每条的start_chainage和方向，统一到绝对里程后重拼。反向测线忘了取负号也会造成镜像错位，这类错误在只看单条剖面时完全发现不了，必须拼图才暴露。

5.2 参数与结果类问题

现象：低频段拾取到的频散曲线突然跳到高速，跟中高频段接不上。这是空间假频或高阶模式被误当成基阶。低频波长长，道间距不够密时相邻道相位差超过 π，出现假频，能量脊跳到错误位置。近场效应也会让低频段能量偏离。解决办法是低频段只用远偏移距道参与计算，或者把相速度扫描下限抬高，避开明显不合理的高速分支。拾取时人工检查低频端几个点，别全信自动拾取。

现象：拼接剖面上出现规则竖向条纹或锯齿。两条原因。一是各测线反演输出的深度采样不一致，直接混在一起插值时产生阶梯；二是重叠段没做加权过渡，硬拼接留下接缝。解决方式是在进插值之前把所有测线重采样到统一的深度网格（比如每 0.5 m 一点），重叠段按距离做线性或余弦加权。这两步都不复杂，但不做的话图交出去会被问是不是数据质量问题。

现象：反演残差很小，但 Vs 剖面随深度剧烈跳动，物理上不合理。层数给多了。层数一多，反演的非唯一性被放大，优化器能找到很多残差相近但形态迥异的解。把层数降到 4 到 6 层，加一个平滑约束（相邻层 Vs 差值的惩罚项）通常能压住。另一个可能是频带太窄，可用的频率点太少，约束不住那么多参数，这种情况下减层数比调优化器有效。

6. 二次开发：把它改成你自己项目能用的出图工具

拿到源码后最想改的地方通常是三个：频散提取算法、反演策略、出图样式。三者的改动成本差很多，建议按"先换出图、再换反演、最后动提取"的顺序来，每一步都能独立验证。

出图样式改动最安全，因为不碰数值。plot_section里的色标、分级、坐标轴、标注全在绘图函数里，换成本单位常用的模板（比如加钻孔柱状图叠加、加高程坐标）只需要改这一个函数。验证方式是拿一条已知结果的测线跑一遍，网格数据不变、图变了，就算改对。

反演策略的替换要守住接口。invert()的输入是(freqs, c_obs)，输出是模型参数和残差，内部换成遗传算法、模拟退火、或者神经网络初值加局部优化都行，只要接口不变，上游拾取和下游拼接都不用动。我一般会在换算法时保留原实现，加一个method参数切换，同一批数据两套算法各跑一遍，残差和剖面形态都对比过再决定用哪个。单看残差选算法容易被误导，非唯一性意味着残差小的模型不一定更接近真实地层。

频散提取算法改动风险最高，因为它是整条链的源头，换个变换方法（比如 f-k 换成相移、或者加个多模式联合提取）出来的能量图形态会变，下游拾取习惯也得跟着调。改这步之前先用合成数据验证——按照已知的分层模型正演出理论道集，再走一遍你的提取流程，看能不能还原出理论频散曲线。合成数据能过，再上实测数据。这个习惯救过我很多次，特别是改相位校正公式的时候，实测数据出问题根本分不清是算法错还是数据差。

还有一个容易被忽略的改造点：结果的组织方式。原始工程可能把每条测线的结果存成独立文件，测线一多就散得到处都是。改成"一个工区一个工程文件 + 结果存 HDF5"会好管很多，h5py按/{测线号}/vs组织，读的时候一次遍历全取出来。工程文件里记好参数版本，后面回看某张图是哪套参数跑出来的，不用靠回忆。

批量处理脚本是二次开发里性价比最高的东西。前面提过的batch/run_all.py，模板大概是遍历工程里的测线、读参数、跑预处理到反演、结果落盘、生成缩略图。写一次，之后每次数据更新重跑一遍就行，比在界面上点几十次可靠得多。脚本里加个日志，每条测线的残差、耗时、警告都记下来，出问题能快速定位是哪条测线哪一步。

从那以后我每次改反演相关的东西，都强制先跑一遍合成数据回归，再拿一条实测线对比改动前后的剖面差异，最后才批量出图。这个流程看着慢，但比图交出去被退回来重做快得多。希望帮到你。

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