如何构建两口钻井地层连接的完整地质剖面？

我有一个存储两口钻井数据的列表，数据格式为[深度起始值, 深度终止值, 地层描述]，具体数据如下：

wells = [
    [[0, -4, 'soil'], [-4, -8, 'clay'], [-8, -12, 'gravel'], [-12, -20, 'sand'], [-20, -24, 'basalts']],
    [[0, -4, 'soil'], [-4, -16, 'sand'], [-16, -20, 'galka'], [-20, -32, 'pebble'], [-32, -36, 'basalts']]
]

需求说明 #

需要绘制所有地层间的连线，生成包含全部地层的完整地质剖面：

• 正确连接两口井中的相同地层
• 准确识别仅存在于单井的地层，并完成符合地质逻辑的过渡连接

问题现状 #

现有基于matplotlib的代码仅能连接描述完全相同的地层区间，无法实现正确的地层剖面连接效果。

改进思路与实现代码 #

核心问题在于不能仅通过地层名称匹配连接，需要结合地层的序列位置（上下层关系）处理：同名地层作为连接锚点，夹在锚点之间的独有地层，需连接到对应锚点的边界；首尾独有地层同理处理。

以下是改进后的完整代码：

import matplotlib.pyplot as plt

def plot_geological_profile(wells):
    # 定义两口井的x坐标
    well_x = [1, 3]
    fig, ax = plt.subplots(figsize=(8, 10))
    
    # 1. 绘制两口井的地层柱状图
    for idx, well in enumerate(wells):
        x = well_x[idx]
        for layer in well:
            top, bottom, name = layer
            # 绘制地层填充矩形
            ax.fill_between([x-0.2, x+0.2], top, bottom, alpha=0.7, label=name if idx==0 else "")
            # 添加地层名称标注
            ax.text(x, (top+bottom)/2, name.strip(), ha='center', va='center')
    
    # 2. 提取结构化地层数据并匹配同名锚点
    well1_layers = [(layer[0], layer[1], layer[2].strip()) for layer in wells[0]]
    well2_layers = [(layer[0], layer[1], layer[2].strip()) for layer in wells[1]]
    
    anchor_pairs = []
    # 匹配同名地层作为锚点，避免重复匹配
    for w1_layer in well1_layers.copy():
        w1_top, w1_bottom, w1_name = w1_layer
        for w2_layer in well2_layers.copy():
            w2_top, w2_bottom, w2_name = w2_layer
            if w1_name == w2_name:
                anchor_pairs.append((w1_layer, w2_layer))
                well1_layers.remove(w1_layer)
                well2_layers.remove(w2_layer)
                break
    
    # 3. 连接同名锚点的顶、底边界
    for w1_layer, w2_layer in anchor_pairs:
        w1_top, w1_bottom, _ = w1_layer
        w2_top, w2_bottom, _ = w2_layer
        ax.plot([well_x[0], well_x[1]], [w1_top, w2_top], 'k--', linewidth=1)
        ax.plot([well_x[0], well_x[1]], [w1_bottom, w2_bottom], 'k--', linewidth=1)
    
    # 4. 处理单井独有地层的过渡连接
    # 处理第一口井剩余的独有地层
    for layer in well1_layers:
        top, bottom, _ = layer
        # 匹配上层锚点（当前地层顶对应锚点的底）
        upper_anchor = next((p for p in anchor_pairs if p[0][1] == top), None)
        # 匹配下层锚点（当前地层底对应锚点的顶）
        lower_anchor = next((p for p in anchor_pairs if p[0][0] == bottom), None)
        if upper_anchor and lower_anchor:
            ax.plot([well_x[0], well_x[1]], [top, upper_anchor[1][1]], 'k--', linewidth=1)
            ax.plot([well_x[0], well_x[1]], [bottom, lower_anchor[1][0]], 'k--', linewidth=1)
    
    # 处理第二口井剩余的独有地层
    for layer in well2_layers:
        top, bottom, _ = layer
        upper_anchor = next((p for p in anchor_pairs if p[1][1] == top), None)
        lower_anchor = next((p for p in anchor_pairs if p[1][0] == bottom), None)
        if upper_anchor and lower_anchor:
            ax.plot([well_x[0], well_x[1]], [upper_anchor[0][1], top], 'k--', linewidth=1)
            ax.plot([well_x[0], well_x[1]], [lower_anchor[0][0], bottom], 'k--', linewidth=1)
    
    # 可视化优化
    ax.set_ylim(-38, 2)
    ax.set_xlim(0, 4)
    ax.set_xticks(well_x)
    ax.set_xticklabels(['Well 1', 'Well 2'])
    ax.set_ylabel('Depth (m)')
    ax.set_title('Geological Profile')
    # 去重图例
    handles, labels = plt.gca().get_legend_handles_labels()
    by_label = dict(zip(labels, handles))
    ax.legend(by_label.values(), by_label.keys(), loc='upper right')
    plt.grid(axis='y', linestyle='--')
    plt.show()

# 调用生成剖面
plot_geological_profile(wells)

代码说明 #

1. 柱状图绘制：先完成两口井的地层柱状展示，清晰呈现单井地层结构。
2. 锚点匹配：优先匹配同名地层作为剖面连接的基准，确保核心地层的对应关系正确。
3. 独有地层处理：针对单井独有的地层，通过上下相邻的锚点边界完成过渡连接，符合地质剖面的渐变逻辑。
4. 可视化优化：添加图例、网格、坐标标注，提升剖面的可读性和专业性。

内容的提问来源于stack exchange，提问作者Edward Zotov