1. CMP道集加载与道头关键字校验

为什么先做道头校验

SEGY 道集里真正决定“这一道属于哪个 CMP、偏移距多少、坐标在哪”的信息，全部存放在每道开头的道头（trace header）中。速度谱拾取依赖 CDP 把道归入面元、依赖 offset 建立炮检距序列、依赖坐标判断面元几何是否合理。如果道头缺失、全零、越界或同一 CMP 内 CDP 不一致，后续拾取就会把不同面元的道混在一起，速度谱出现虚假能量团。因此第一步不是看波形，而是把道头当作“数据身份证”逐项核对。

关键道头字与校验口径

按 SEG-Y 标准（Rev1）常用字节位置，本例关注以下字段（字节号从 1 开始计，读取时按大端/小端由文件头 Data sample format 与字节序标志确定，这里假设大端）：

校验分三层：存在性（字段是否被写入、是否全零）、范围性（数值是否落在工区合理区间）、一致性（同一 CMP 内 CDP 是否相同、偏移距是否重复、面元与 CDP 映射是否稳定）。三层都通过，才认为道头可用于后续速度谱拾取。

读取与校验的实现

下面用 Python 直接按字节偏移读取道头，不依赖第三方地震库，便于在任何环境复现。核心是把 3200 字节文本头 + 400 字节二进制文件头之后的每道按 240 字节道头 + ns×4 字节样点 切分，只解析道头部分。

import struct
from collections import defaultdict

def read_segy_headers(path, ns=None, endian='>'):
    """只读道头，返回每道的关键字段列表。endian='>' 大端，'<' 小端。"""
    with open(path, 'rb') as f:
        f.seek(3200)
        bin_header = f.read(400)
        # 二进制文件头中 samples per trace 位于 3221-3222 字节（相对文件头 20-21）
        ns_file = struct.unpack(endian + 'H', bin_header[20:22])[0]
        ns = ns or ns_file
        trace_len = 240 + ns * 4
        traces = []
        idx = 0
        while True:
            th = f.read(240)
            if len(th) < 240:
                break
            cdp     = struct.unpack(endian + 'i', th[20:24])[0]
            offset  = struct.unpack(endian + 'i', th[36:40])[0]
            cdp_x   = struct.unpack(endian + 'i', th[180:184])[0]
            cdp_y   = struct.unpack(endian + 'i', th[184:188])[0]
            bin_in  = struct.unpack(endian + 'i', th[188:192])[0]
            bin_xl  = struct.unpack(endian + 'i', th[192:196])[0]
            traces.append(dict(idx=idx, cdp=cdp, offset=offset,
                               x=cdp_x, y=cdp_y, bin_in=bin_in, bin_xl=bin_xl))
            f.seek(trace_len - 240, 1)   # 跳过样点
            idx += 1
    return traces, ns

读取完成后做三层校验。存在性用“全零/缺失”判断；范围性用可配置阈值；一致性按 CMP 分组统计。

def validate_headers(traces, x_range, y_range, off_max):
    issues = []
    by_cmp = defaultdict(list)
    for t in traces:
        # 存在性
        if t['cdp'] == 0:
            issues.append((t['idx'], 'CDP_MISSING', t['cdp']))
        if t['x'] == 0 and t['y'] == 0:
            issues.append((t['idx'], 'COORD_ZERO', (t['x'], t['y'])))
        # 范围性
        if not (x_range[0] <= t['x'] <= x_range[1]):
            issues.append((t['idx'], 'X_OUT_OF_RANGE', t['x']))
        if not (y_range[0] <= t['y'] <= y_range[1]):
            issues.append((t['idx'], 'Y_OUT_OF_RANGE', t['y']))
        if t['offset'] < 0 or t['offset'] > off_max:
            issues.append((t['idx'], 'OFFSET_RANGE', t['offset']))
        by_cmp[t['cdp']].append(t)

    # 一致性：同 CMP 内 CDP 应唯一（按定义），偏移距不应重复
    for cdp, group in by_cmp.items():
        offs = [g['offset'] for g in group]
        if len(offs) != len(set(offs)):
            issues.append((cdp, 'OFFSET_DUP', len(offs) - len(set(offs))))
        bins = {(g['bin_in'], g['bin_xl']) for g in group}
        if len(bins) > 1:
            issues.append((cdp, 'BIN_INCONSISTENT', bins))
    return issues, by_cmp

最后输出每个字段的数值范围与一致性摘要，作为“数据体检报告”：

def summarize(traces, by_cmp):
    def rng(key):
        vals = [t[key] for t in traces]
        return min(vals), max(vals)
    print('道数:', len(traces), ' CMP 数:', len(by_cmp))
    print('CDP 范围:', rng('cdp'))
    print('偏移距范围:', rng('offset'))
    print('X 范围:', rng('x'), ' Y 范围:', rng('y'))
    print('每 CMP 道数分布:', sorted({len(v) for v in by_cmp.values()}))

校验结果如何判读

把上面的输出整理成一张“道头健康度”视图，能直观看出问题集中在哪一层。下图把读取、三层校验与判读结论串成一条数据流，并标出典型异常对应的处理动作。

实操要点与通过条件

• 字节序先确认：SEGY 有大小端之分，读出的 CDP 若出现巨大负数或乱序，先检查 endian 参数，而不是怀疑数据本身。
• ns 以文件头为准：道头读取时用二进制文件头的 samples per trace 计算跳道长度，避免逐道读样点导致速度慢。
• 阈值来自工区：坐标范围、最大偏移距应取自观测系统设计值，不要用固定常数；示例中的 25–2475 m 仅为演示。
• 通过条件：CDP 无零值、坐标无全零、偏移距非负且不重复、同 CMP 内面元唯一、每 CMP 道数接近设计道数。全部满足才进入速度谱拾取。