企业官网建设流程全解析

简介：一套PyGIMLi示例脚本，聚焦地质建模与反演，面向地球物理、地质工程领域的初学者和研究者。压缩包内共11个Python文件，大小约19KB，全部为可直接运行的.py脚本，并按基础、建模、反演三大主题分类组织。内容涵盖网格创建与插值、网格质量检查、有限元正演（热传导问题）以及直流电阻率一维反演等典型场景；每个脚本对应一个独立任务，如指数拟合、多项式拟合、一维块状反演等，便于边运行边理解PyGIMLi的调用流程和参数设置。示例还展示了软件在几何建模、多种数值方法（有限元、有限差分）与反演框架上的核心用法，适合作为系统学习的起始模板。已有1069人学习，尤其适合希望通过实例快速掌握几何建模、数值求解、反演算法与可视化技巧的读者，也可作为实际项目中的参考模板。

1. 遇到 pyGIMLi 例子：从一头雾水到能跑通第一行反演

如果你刚接触 geophysical 反演，大概率会在某个深夜搜到 pyGIMLi 这个词，然后被一堆 example 脚本和术语劝退。pyGIMLi 是一个用于二维/三维电阻率、地震波、电磁法等地球物理正演与反演的开源库，它的例子不是在教你怎么调用某个函数，而是在展示一条完整的处理链路：网格生成、正演计算、合成数据加噪、反演参数设置、结果可视化，最后才是你真正要的那张地电断面图。看懂这些例子，等于拿到了从数据到模型的钥匙。

这篇文章把 pyGIMLi 的例子拆开讲，按「例子在解决什么问题 → 怎么把例子跑起来 → 参数怎么调 → 反演结果怎么判断 → 坑在哪 → 还能怎么玩」的顺序走一遍，全程用可复现的命令和代码，不绕弯子。适合刚入行想快速上手的学生，也适合接到勘探任务需要短期内用 pyGIMLi 出图的一线从业者。读完你能独立改例子、换数据、调反演参数，并且知道哪些地方容易翻车。

2. pyGIMLi 例子到底在讲什么：核心对象和一条完整链路

2.1 例子背后的三个核心对象：网格、数据和反演框架

pyGIMLi 的例子里几乎都会同时出现三个东西：网格（Mesh）、数据容器，以及反演算法入口。理解这三个对象，比背函数名重要得多。网格负责描述模型空间，数据容器装着你在野外测到的视电阻率或地震走时，反演框架则负责把数据映射到模型上。

以最常见的电阻率例子为例，你会看到类似pg.meshtools.createGrid创建网格，然后用pg.DataContainer承载电极编号和测量值，再用pg.RobustInversion或pg.Inversion做反演。这三个对象不是孤立存在的，它们通过一个「正演算子」串起来，pyGIMLi 里叫pg.ERTModelling，它负责根据当前模型算出理论观测数据，也就是正演，反演过程就是在不断调整网格单元的电导率，让正演结果逼近实测数据。

我个人的经验是，先别急着跑完整例子，把这三个对象在例子里各自干了什么标记出来，然后你会发现 pyGIMLi 所有例子的骨架几乎都一样。这个共性就是它最值钱的地方：学会一个，等于学会了这一类。

2.2 例子的标准链路：从网格生成到反演出图

pyGIMLi 官方仓库里的例子，无论领域怎么变，链路都是固定的五步：

第一步，生成网格。大多数例子用pg.meshtools.createWorld或者直接读入一个已有的网格文件。第二步，加载或仿真数据。电气法例子里常用pg.ERTManager里的simulate方法，地震走时例子则用pg.SeismicManager。第三步，设置反演参数，包括正则化权重、Z 权重、迭代次数等。第四步，调用inv.run()启动反演。第五步，用pg.show可视化模型。

这里要提醒的是，很多例子里的「数据」是为了演示生成的合成数据，并不是你野外的实测数据。合成数据的好处是真实模型已知，可以用来验证反演代码对不对，但代价是它太干净了，野外数据里的电极接触电阻、地形起伏、坏点这些问题它都没有。所以当你把例子跑通以后，一定要意识到：例子只是模板，替换成自己的数据才是开始。

2.3 最小可运行例子：一个二维电阻率反演的完整代码

这里直接给一个我在学习和给新人讲解时最常用的最小例子，它基于 pyGIMLi 的 ERT（电阻率层析成像）模块，代码很短，但足够说明整条链路怎么闭合。

import pygimli as pg from pygimli.physics import ert # 1. 创建一个简单的两层地电模型网格 grid = pg.meshtools.createWorld(start=[-10, 0], end=[10, -10], layers=[-2], world=[False]) # 2. 在网格上定义真实电阻率值，上层 100 欧姆米，下层 10 欧姆米 model_true = pg.solver.parseMapToCellArray([[1, 100.0], [2, 10.0]], grid) grid["res"] = model_true # 3. 生成一个施伦贝格装置的数据容器 data = ert.createData(elecs=pg.utils.grange(-10, 10, dx=1.0), schemeName='schlumberger') # 4. 用正演模拟生成合成视电阻率数据，并加 3% 高斯噪声 data = ert.simulate(mesh=grid, res=model_true, data=data, noiseLevel=0.03, noiseAbs=1e-2) # 5. 创建反演管理器并运行反演 mgr = ert.ERTManager() model = mgr.invert(data, verbose=False) # 6. 输出反演结果并绘图 print(model) pg.show(mesh=model.parent(), data=model, showMesh=True)

这段代码的逻辑说明：createWorld生成一个带有地下一层的矩形网格，parseMapToCellArray给不同区域赋电阻率值，createData用线性排列的 21 个电极生成标准施伦贝格测量方案，simulate调用有限元正演得到视电阻率并添加噪声，ERTManager.invert启动反演优化，最后model对象里就是反演出的每个网格单元的电阻率。参数说明里最重要的是noiseLevel，它直接影响反演难度，建议从 0.01 开始逐步加大，太大反演会失真，太小结果会过拟合。

这段代码我建议你逐行读过一遍再跑，因为它是后续所有电阻率例子的骨架。如果你跑的是地震或直流电法例子，把ert换成对应模块，流程是一样的，差异只在于正演算子和数据容器。

3. 选对例子并改造成自己的数据：从读例子到替换实测数据

3.1 官方例子的目录结构：怎么判断哪个例子适合自己

pyGIMLi 的 examples 目录一般按物探方法分类，常见的有ert（电阻率）、traveltime（地震走时）、tdem（时间域电磁）、inversion（纯反演框架）等。第一次接触的人容易直接挑一个名字看起来顺眼的文件，然后发现里面用到了自己没装好的第三方库，或者数据文件路径对不上。

我的建议是，先看文件名里有没有simulate或synthetic，如果有，说明它自带数据生成逻辑，不依赖外部数据文件，最适合先跑通。再挑一个自己不熟悉的装置类型，比如从wenner换到schlumberger，对比同一条反演链路的差异，会比重复跑同类型例子收获大得多。官方例子的作用是给你提供「最小可用框架」，千万别把它当成标准答案。

3.2 把合成数据换成实测数据的关键改动点

当你拿到自己的野外数据，格式往往是protocol文件或者csv，里面有电极编号、电极坐标和测量的视电阻率值。pyGIMLi 提供了ert.load方法，可以直接读入多种格式，但实测数据里最常见的坑是坐标单位不一致、电极编号从 1 开始而你导入后续索引从 0 开始，还有缺失的电极位置信息。

我一般在替换数据时会做三件事：第一，把数据文件转成 pyGIMLi 的DataContainer，如果原数据是txt的「ABMN R」格式，我写一个很小的解析脚本；第二，检查数据范围，删除视电阻率为负值或明显偏离中位数几个数量级的点；第三，给数据容器附带正确的电极坐标，这一步决定了网格能不能正确覆盖测线。

import pygimli as pg from pygimli.physics import ert # 假设你的实测数据文件是 csv，列依次为：A, B, M, N, rho # 先读成 numpy 数组，再构建 pyGIMLi 数据容器 import numpy as np raw = np.loadtxt('field_data.csv', delimiter=',', skiprows=1) # 提取电极坐标，这里假设电极沿 x 方向从 0 到 40 均匀分布 elecs = pg.utils.grange(0, 40, n=len(set(raw[:, :4].ravel()))) # 创建容器并设置四极装置 data = ert.createData(elecs=elecs, schemeName='custom') # 把实测数据填入容器，注意 pyGIMLi 的 data['rhoa'] 对应视电阻率 data['rhoa'] = raw[:, 4] data['valid'] = np.ones(len(data), dtype=bool) # 过滤掉明显异常值（视电阻率小于 1 或大于 10000 欧姆米） bad = (data['rhoa'] < 1) | (data['rhoa'] > 10000) data.markInvalid(bad) data.removeInvalid() # 打印数据概况检查 print(data)

这段代码里有个关键的参数是schemeName='custom'，当你没有标准装置方案时用它最省事，pyGIMLi 不会强加电极配置约束。markInvalid和removeInvalid是处理坏点的标准动作，实测数据里经常有接地电阻不好导致的飞点，不处理会让反演结果出现大片假异常。注意我这里的电极坐标用grange均匀生成，如果你的实测电极间距不等，这个假设就错了，需要从数据文件里真正读取每个电极的坐标。

3.3 为什么要复用反演管理器而不是每次重建

很多例子最后都会出现mgr = ert.ERTManager()然后mgr.invert(data)。有的新手会问，为什么不直接用pg.Inversion类，非要套一个管理器。管理器的价值在于：它封装了正演算子、数据变换、误差加权、默认边界条件和可视化等一堆反演外围逻辑，你只改数据就可以迭代反演，不需要重复配置。

实际使用中我更推荐这种模式：数据预处理一次，生成多个测线数据文件，然后循环喂给同一个ERTManager，每次invert会复用上一次的网格初始猜测，对多测线联合成像特别有用。但注意，复用管理器有一个副作用：上一次反演得到的模型会成为下一次的初始模型。如果不是连续测线，记得调用mgr.reset()或者在创建管理器时指定mgr = ert.ERTManager(verbose=True)，避免上一测线的残差影响本次结果。

4. 反演参数逐个拆解：从默认值开始调，而不是冲动改

4.1 正则化参数 lambda 和 Z 权重：控制平滑度的两个旋钮

pyGIMLi 的大部分例子默认使用二阶平滑约束正则化，通过pg.meshtools生成的网格配合盒状约束，使反演结果在地层边界处仍能保持清晰。反演控制器里常用的参数是lam（正则化权重）和zWeight（垂直与水平平滑比例）。如果你直接用mgr.invert(data)，库会启动一个自动选择 lambda 的机制，从较大值开始逐渐降低，直到拟合残差达到目标。

这看起来省心，但在实际数据上，自动选项往往偏保守，结果过于平滑，分辨率不够。我的做法是，先把自动反演跑一遍，记下它最终选择的 lambda 值，然后把它乘以 0.3 到 0.5 作为手动初始值，再看模型细节是否增多且不发散。zWeight 的作用是调整垂直和水平方向平滑的权重，默认值在 pyGIMLi 里通常是 1，意味着各向同性平滑。对测线方向拉长的二维结构，我会把 zWeight 设为 0.5 到 0.8，让垂直方向有更多变化空间。

inv = mgr.inversion # 获取反演对象 inv.setLambda(20.0) # 手动指定正则化权重，而不是自动选择 inv.setZWeight(0.8) # 垂直方向平滑比水平方向略弱 # 设置反演停止条件：相对残差变化小于 0.01 或达到 15 次迭代 inv.setMaxIter(15) inv.setChiLimit(1.01) mgr.invert(data, startModel=pg.meshtools.createGrid(...))

参数说明：setLambda(20.0)在我的同类型数据集上是从自动值 54 调下来的，模型细节明显增加但仍在稳定范围；setChiLimit控制卡方目标，1.01 对含 3% 噪声的数据比较合理，如果噪声更大可以放宽到 1.05。这里要注意startModel必须是一个网格或者与待求模型长度一致的数组，否则反演会直接报错。

4.2 数据误差模型：比正则化更值得先调的东西

pyGIMLi 默认认为数据是高斯独立同分布的，误差由绝对误差和相对误差组成。在 ERT 模块里，simulate函数生成合成数据时加的noiseLevel=0.03和noiseAbs=1e-2，到了反演阶段，pyGIMLi 会默认用数据本身的中位数来估计误差比例。这个假设在大部分野外数据上是脆弱的。

如果你用的是实测数据，强烈建议在反演前显式设置数据误差。pyGIMLi 中data['err']是一个可选字段，如果没有，反演内部会用经验公式生成。我踩过最直观的坑是，不同电极排列的测量数据误差差异很大，比如偶极-偶极排列在长极距时信噪比低，如果统一用一个固定相对误差，长极距的噪点会被放大成浅部假异常。

正确做法是给每个测量值设定独立误差，一般做法是data.setError(np.abs(data['rhoa']) * 0.02 + 0.005)，或者更稳妥地从重复测量数据里估算标准差。pyGIMLi 反演框架在加权最小二乘里会把误差作为分母，误差越小权重越高，所以低误差数据会对反演产生主导作用。如果你发现自己的反演结果被某个数据点控制，去查一下那个点的误差是否过小。

4.3 网格密度与地形：两个容易被忽视的边界

pyGIMLi 例子里常用createWorld生成规则网格，边界处网格无限向外延伸，避免边界反射。当你把例子改成实测数据时，必须重新检查网格是否覆盖所有电极坐标，并考虑地形。pyGIMLi 支持带地形的网格，常见做法是用pg.meshtools.createSurveyMesh根据电极坐标生成贴合地表的网格。

网格太疏会导致分辨率不够，太密会导致反演矩阵巨大，而且浅层容易出现锯齿。我的经验是，在电极间距的 1/5 到 1/10 之间设置网格细化级别。你不需要手工建网格，用pg.meshtools.createSurveyMesh( data, quality=34)会自动在电极附近加密网格。地形处理方面，如果测线高差超过电极间距的 5%，必须把电极坐标的 z 值传入网格生成函数，否则反演模型会把人造地形起伏当成浅部电性异常。

5. pyGIMLi 例子避坑：我在实际项目里踩过的五个真实问题

5.1 现象：反演迭代发散，每次迭代残差反而增大

原因：最常见的是数据里有负值或零值。pyGIMLi 反演通常采用对数电阻率参数，零或负值取对数直接报错，即使没报错，也会导致正演矩阵奇异。另一种可能是起始模型给了一个极端值，比如startModel给了 1e-6 欧姆米，正演计算浮点溢出。

解决：先做数据清洗，data['rhoa'] <= 0的点一律删掉；再看startModel，用中位数初始化。代码里顺手加一行data.removeInvalid()，能避免 80% 的发散。

5.2 现象：反演结果出现带状条纹，沿测线方向规则排列

原因：电极位置坐标有误，最常见是电极间距不等但用了均匀网格，或者 x 和 z 坐标填反。另一个原因是测量数据存在系统偏差，比如电极接触电阻大导致所有相邻极距的数据整体偏高，反演不得不以条带状模型补偿。

解决：先用pg.show(data)查看数据的伪剖面图，如果伪剖面里有明显的横向突变，优先怀疑数据质量；再检查网格与电极坐标是否严格匹配，使用mgr.model与电极位置比对。

5.3 现象：反演结果与真实地质分层不符，深部电阻率平均化严重

原因：这不是 pyGIMLi 的问题，而是观测装置对深部的灵敏度不足。偶极-偶极装置比温纳装置有更好的横向分辨率，但深部信号弱且误差大。如果沿用例子里的施伦贝格装置，深部约束本来就弱。

解决：换装置类型重新采集数据，或者在反演中加入先验约束。pyGIMLi 支持pg.Region约束某个区域电阻率固定，比如你知道表层是冲积层，就把浅层区域约束在一定范围内，减轻深部平均化。

5.4 现象：运行反演时内存爆炸，电脑直接卡死

原因：网格太密，或者反演求解器默认用了直接法。三维反演时尤其明显，二维问题一般还好。在高分辨率地形起伏大的测线上，网格数量上了百万，直接法求一次正演就吃掉几十 G 内存。

解决：降低网格质量参数，从quality=34调到quality=30，或者限制网格最大面积；如果 pyGIMLi 版本支持，指定求解器为迭代法，常见做法是创建pg.solver.Solver并设置solver='cg'，然后传给反演管理器。不过这会让反演时间变长，需要权衡。

5.5 现象：pg.show显示模型时颜色范围永远被少量极端值霸占，断面图看不出结构

原因：这个是可视化常见误用，pyGIMLi 默认的色标范围是数据的最小最大值，非鲁棒。

解决：手动指定cMin和cMax，一般用数据的 5% 和 95% 分位数，或者直接先对模型取对数再显示。pg.show(mesh, model, cMin=10, cMax=200)就比自动范围实用得多。这个小技巧在实际报告里非常重要，审阅者第一眼看到的就是图，颜色范围不对等于模型白反演了。

6. 把例子进化成自己的研究工具：批量反演与结果质量评估

6.1 用脚本批量处理多条测线：从手动跑例子到流程自动化

你手里的野外数据一般不止一条测线。手动改文件名跑例子太慢，我通常会写一个循环，批量读入测线、反演、保存结果图。这一步能最大限度发挥 pyGIMLi 的价值。

import glob from pygimli.physics import ert # 获取所有测线数据文件 files = glob.glob('survey_*.csv') for f in files: data = ert.load(f) # 假定是 pyGIMLi 支持格式 mgr = ert.ERTManager() model = mgr.invert(data, verbose=False) # 保存反演结果模型为 vtk 格式，方便后续在 paraview 里查看 pg.save(f.replace('.csv', '_model.vtk'), model) # 保存反演拟合图 ax = pg.show(mgr.fop.mesh(), model, cMin=10, cMax=200) pg.show(data, ax=ax, cMin=10, cMax=200)

逻辑说明：glob匹配所有测线文件，ert.load自动识别常见格式，然后每条测线单独创建ERTManager避免初始模型串扰。保存为 vtk 是因为后续我想在专业可视化软件里剖切联合显示，pyGIMLi 自带的绘图适合快速验证但不好排版。

运行这段之前，务必在每条测线上跑一次单条数据，确认没有坏点。自动化里最容易出问题的就是某条测线数据有问题导致反演中断，所以在循环外层加try...except记录失败文件，避免整个批处理崩掉。

6.2 用 L-curve 和模型粗糙度判断反演结果可信度

pyGIMLi 反演完成后，除了看图，还需要量化说明结果可信度。mgr.invert的历史记录里保存了每次迭代的残差和模型粗糙度，我常用mgr.inversion.chi2看最终卡方值，如果它远大于 1，说明数据误差估计偏低或者模型不满足数据；如果小于 0.8，可能过拟合了。

一种更直观的验证方法是做「正演验证」：把反演得到的模型作为真实模型，用同样的测量方案重新生成模拟数据，再对比模拟数据和实测数据的拟合程度，也就是交叉验证思路。我经常这么做：用一套网格做反演，用另一套更密或更疏的网格做正演，看结果是否稳健。如果模型对网格变化特别敏感，说明它主要由正则化控制，而不是由数据驱动，解释时要格外小心。

import pygimli as pg from pygimli.physics import ert # 已知反演模型 model_res 和对应的网格 mesh mgr_check = ert.ERTManager() new_data = mgr_check.simulate(mesh=mesh, res=model_res, scheme=data.scheme()) # 对比模拟数据与实测数据的归一化 RMS rms = pg.norm((new_data['rhoa'] - data['rhoa']) / data['rhoa']) / pg.norm(data['rhoa']) print('Normalized RMS:', rms)

这段代码的用途不是重新反演，而是单纯做正演模拟来验证模型。如果模拟数据和实测数据的 RMS 在 0.05 以内，说明这个模型对当前数据是合理的，但如果 RMS 很大，即使反演收敛了，也要怀疑网格或误差设置有问题。

6.3 我最常用的一个技巧：反演前把数据做对数变换

pyGIMLi 内部在反演时默认对电阻率做对数变换，但很多例子在预处理阶段不会主动告诉你。野外视电阻率数据跨度经常从个位数到几千欧姆米，如果不做对数变换，高阻值的数据点会在目标函数中占据主导，低阻地区的细节会被掩盖。

我有个习惯，在把实测数据塞进反演之前，先把rhoa取对数，并对对数数据设置误差，这样反演更稳定。pyGIMLi 里可以通过data.log()复制一份对数数据容器，或者直接在读取后设置data['rhoa'] = np.log(data['rhoa'])，然后反演管理器内部会自动检测。最终显示的时候你再指数还原，注意pg.show的色标标签标成欧姆米即可。

这个技巧对电阻率跨几个数量级的项目帮助很大。我印象最深的是一次在第四系覆盖区做工程勘察，表层干砂电阻率上千，下部含水黏土只有几十，不做对数变换反演，含水层边界基本糊成一片，做了对数变换后分层清晰多了。反演不仅仅是调参数，数据变换本身也是反演的一部分，这个观点在 pyGIMLi 的例子里体现得不明显，需要你自己在实际项目中体会。

希望我讲的这些能让你少走一点弯路，跑通例子只是一小步，用好参数和数据判断才是在这个方向投入时间的真正回报。

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