Pyrocko + PSGRN/PSCMP小问题

1.先看看你的脚本然后诊断 config 文件的问题。问题很明确——YAML 解析config文件时在earthmodel_1d块标量那里报错。大概率是|2缩进指示符和实际内容缩进不匹配。让我先下载脚本看看然后直接诊断fomosto 不在当前环境但不影响修脚本。问题就是 YAMLconfig文件里的|2块标量缩进指示符。分析earthmodel_1d: |2 # ← |2 表示内容缩进为 2 格 0. 3.5 2.0 ... # ← 实际只缩进 1 格矛盾|2告诉 YAML 解析器内容行缩进了 2 个空格但实际内容只缩进 1 个空格。解析器找不到预期的缩进层级后面的tabulated_phases就被当成了映射的 key导致did not find expected key。修复方案把|2改成|自动检测缩进同时确保内容行至少缩进 1 格。修正后的脚本修正完成。改动很小但关键改动对比- earthmodel_1d: |2 # ← |2: 声明内容缩进2格但实际只有1格 earthmodel_1d: | # ← |: 自动检测缩进兼容性最好 **原理** - |2 内容缩进了 2 空格请剥掉 2 空格 → 实际内容只缩了 1 空格 → 解析器读到 0. 变成 .3.5... → YAML 炸了 - | 自动检测最小缩进 → 正确剥掉 1 空格 → 内容完整保留 在你的服务器上运行修正版 bash cd ~/202402351 bash run_experiment_fixed.sh另外如果你跑完fomosto ttt之后fomosto build报错说找不到 psgrn/pscmp 二进制确认一下which psgrn 2/dev/null; which pscmp 2/dev/null 这俩需要从 fomosto-psgrn-pscmp 编译安装后才能用。2.好YAML 解析通过了现在是新问题根因begin阶段定义里有Pv_(cmb)p需要找 CMB核幔边界~2891km 深但你的地球模型只到 200km。cake模块找不到cmb不连续面就炸了。对于 PSGRN/PSCMP近场静态位移/震后变形根本不需要引用 CMB 的远震震相。直接简化 phase 定义改动说明- definition: p,P,p\,P\,Pv_(cmb)p # 引用 CMB (2891km)模型只到 200km → 爆炸 definition: 0.0 # 用固定时间值不依赖深层不连续面begin/end的作用它们只定义格林函数时间窗口的边界最早到达、最晚到达对近场 PSGRN/PSCMP 来说用固定秒数完全够用。在你服务器上重新跑cd ~/2024023551 bash run_experiment_fixed.sh3.如果后面P/p/S/s也报类似错它们依赖特定射线路径也可以简化掉只留begin和end就够了。PSGRN/PSCMP 计算的是近场静态/准静态位移不是动态波形震相表没那么关键。definition: 0.0被 Pyrocko 当成阶段名查找而非常量时间还是触发射线追踪。换个策略——不再整个覆盖 config只修补earthmodel_1d那一行保留fomosto init自动生成的默认 phase 定义它自带适合 psgrn_pscmp 的配置等等Pyrocko 的 YAML 有自定义 tag!pyrocko.gf.meta.ConfigTypeAyaml.safe_load解析不了。去掉验证同时保留你需要的自定义参数深度范围、距离范围等同时把extra/psgrn_pscmp的替换加回来这个文件没有 YAML 问题确认最终脚本策略变了不再整个覆盖 config之前手动写 phase 定义总出问题改为1.fomosto init生成默认配置自带正确的 phase 定义2.用 Python 精准替换earthmodel_1d修掉|23.regex 更新深度/距离参数4.单独写extra/psgrn_pscmp在服务器上跑cd ~/202402351 bash run_experiment_fixed2.sh如果fomosto tt还是报错把完整输出贴过来。另外跑之前确认一下which psgrn which pscmp确保 PSGRN/PSCMP 的二进制已经编译安装好了不然fomosto build也会挂。4.两个独立的问题问题 1: PSGRN —wrong no of distance samples!距离采样参数与 PSGRN 内部期望不匹配。你用 regex 覆盖了fomosto init的默认 distance 值但默认值已经是 PSGRN 兼容的。问题 2: PSCMP —NSCMAX defined too small!PSCMP Fortran 二进制编译时定义了最大快照数上限通常 100-200。你配置了tmax_days: 3650→ 3650 个每日快照远超上限。修复策略不要覆盖 distance 参数保留默认同时减少震后时间窗口和快照数。改动总结# extra/psgrn_pscmp - tmax_days: 3650.0 # 10年 → 3650个快照 → 超过 NSCMAX tmax_days: 100.0 # 100天 → ~100个快照 → 安全 # configPython部分 - 同时覆盖 distance_min/max/delta # PSGRN 不兼容 只覆盖 source_depth receiver_depth # 保留默认 distance在服务器上跑cd ~/202402351 bash run_experiment_fixed2.sh如果后续你想覆盖 10 年的震后变形有两个方案1.增大sampling_interval比如改成 30.0每月一个快照2.重新编译 PSGRN/PSCMP调大NSCMAX参数在 Fortran 源码的pscglob.h里5.问题很清楚。从日志看有两个级联错误根因Step 1 PSGRN 失败STOP Error: wrong no of distance samples!PSGRN 的格林函数生成直接挂了所以 Step 2 PSCMP 找不到格林函数文件uz.ep不存在也跟着崩。问题出在哪你的脚本把distance_min改成了1000.0单位是米但PSGRN 要求距离网格从 0 开始。distance_min应该设为0.0让第一行零偏移距的格林函数被计算。同时注意单位——fomosto config 里距离和深度的单位通常是米m不是千米。修复脚本中这部分 sed 替换config re.sub(rdistance_min:\s*[\d.e-], distance_min: 1000.0, config)改成config re.sub(rdistance_min:\s*[\d.e-], distance_min: 0.0, config)另外也检查一下distance_delta: 5000.0是否和 PSGRN 兼容——某些版本要求 delta 必须是 2 的幂次或特定值。如果还是报wrong no of distance samples试试把distance_delta改成比如1000.0或2000.0同时保持distance_max / distance_delta为整数。要我直接改好脚本给你吗6.Now I see the full picture. Let me get one more piece — the default StoreConfig distance values for psgrn_pscmp.2008a:OK now I have the full picture. The default StoreConfig for psgrn_pscmp.2008a usesdistance_max50km,distance_delta1km. The real problems are:Let me write the definitive fix. Ive traced through the full pyrocko source — the real issue is that the regex-basedearthmodel_1dreplacement corrupts the YAML config, causing fomosto to fall back to bad internal defaults which trigger the wrong no of distance samples error. The fix is touse Pyrockos own APIinstead of regex surgery.Let me also verify the script is correct by checking indentation:找到根因了。问题出在extra/psgrn_pscmp里的gf_distance_spacing: -1.0。看 pyrocko 源码if cg.gf_distance_spacing is None: # 如果是 None用 StoreConfig 的值 gf_distance_spacing fc.distance_delta else: # 否则用它 * km gf_distance_spacing cg.gf_distance_spacing * km # -1.0 * 1000 -100 n_steps_distance int((distance_max - distance_min) / gf_distance_spacing) 1 # int((5000 - 0) / -1000) 1 -49 ← 负数YAML 里写的-1.0意思是用默认但 Pyrocko 不这么理解——它把-1.0当作实际值乘以km1000得到-1000算出来负的步数PSGRN 当然报wrong no of distance samples。修复extra/psgrn_pscmp里删掉gf_distance_spacing和gf_depth_spacing让它们保持Noneoptional 字段不写就是 None。脚本写好了run_experiment_fixed3.sh核心修复找到真正的 bug问题在extra/psgrn_pscmp里的这两行gf_depth_spacing: -1.0 gf_distance_spacing: -1.0 原意是 -1.0 用默认值但 Pyrocko 不这么理解。看源码 python if cg.gf_distance_spacing is None: # None → 用 StoreConfig 的 distance_delta gf_distance_spacing fc.distance_delta else: # 非 None → 直接当值用 gf_distance_spacing cg.gf_distance_spacing * km # -1.0 × 100 -1000 n_steps int((50000 - 0) / -1000) 1 # -49 ← 负数PSGRN 收到负的步数当然报wrong no of distance samples!。改动汇总| 改动 | 原因 | |------| | 删fomosto ttt| psgrn_pscmp 不需要走时表 | | 删gf_distance_spacing/gf_depth_spacing|根因修复— optional 字段不写 None | | 删observation/rectangular_fault_size_factor等冗余字段 | 都有默认值不需要显式写 | |source_depth_min改 1000m (10km) | 避免浅层边界插值问题 | | 保留 distance 默认值 | 默认 0-50km, delta1km 已验证兼容 | | 构建前打印关键参数 | 出问题方便排查 |7.问题找到了。StaticTarget继承自MultiLocation它接受的是lats/lons数组不是标量lat/lon# 你的代码错 gf.StaticTarget(lat31.06, lon103.36, ...) # ← lat, lon 不是合法参数 # 正确写法 gf.StaticTarget(latsnp.array([...]), lonsnp.array([...]), ...)修复的核心改动改动原因StaticTarget(latslats_flat, lonslons_flat, ...)API 要求 1D numpy 数组一个StaticTarget包含全部网格点不再创建 14641 个单独 target结果用static_result.result[displacement.n]等StaticResult.result是 dictkey 是displacement.n/displacement.e/displacement.d结果.reshape(n_lats, n_lons)把 1D 结果重排回 2D 网格脚本已保存到compute_deformation_fixed.py。在你那边复制到/home/lx/2024023551/目录下跑就行cd /home/lx/2024023551 python3 compute_deformation_fixed.py8.找到了错误信息非常明确ValueError: A value (7341.013570101181) in x_new is below the interpolation ranges minimum value (10000.0).Pyrocko 把断层离散化成子补丁时有些子补丁的中心深度到了7341m7.3 km但 store 的source_depth_min10000m查不到格林函数 → NaN。根因source_depth_min设太高了。改成 0 即可。问题确认了Pyrocko 离散化断层时子补丁中心深度到了7341m但 store 只有source_depth_min10000m的格林函数 → 查不到 → NaN。修复source_depth_min从1000改成0。cd psgrn_pscmp_experiment bash run_experiment_v6.sh python3 compute_deformation_fixed.py改动对照| 参数 | v5 | v6 | |------|----| |source_depth_min| 1000.0 (10 km) |0.0 (地表)| | 其他 | 不变 | 不变 |v6 的source_depth_delta5000意味着从 0 到 50 km 每 5 km 一层共 1 层比 v5 的 9 层只多 2 层建库时间增加不多。补充GF store格林函数库就是一个预计算好的地震响应查询表。通俗理解想象你要查地震源在某个深度、距离你多少公里时地面会怎么动。没有 GF store 的话每次都要从头解弹性波方程算一次要几小时甚至几天。有 GF store提前把各种源深度 × 距离组合的响应全算好存成文件。之后查表就行几秒出结果。具体结构my_psgrn_pscmp_gfs_myanmar/ ├── config ← 配置深度范围、距离范围、地球模型 ├── psgrn_functions/ ← 核心PSGRN 算出的格林函数文件 │ ├── gd.cl ← 不同分量 (clcosine-longitudinal, ...) │ ├── gd.ds │ ├── gr.ss │ └── ... ├── index ← 索引快速定位数据 └── extra/psgrn_pscmp ← PSGRN/PSCMP 特有配置config 里定义了网格source_depth: 0, 5, 10, 15, 20, ... 50 km (11层) distance: 0, 3, 6, 9, ... 500 km (167层)每一组 (深度, 距离) 都存了一组格林函数。计算流程断层参数 (strike, dip, rake, slip, ...) ↓ Pyrocko 把断层离散成小补丁 ↓ 对每个补丁 → 查 GF store → 得到该补丁对每个网格点的响应 ↓ 所有补丁的响应叠加 最终位移场为什么之前一直 NaN问题原因distance_max50kmstore 只存了 50km 内的数据远处查不到source_depth_min10kmPyrocko 离散化后有子补丁到 7.3km查不到就像字典里没收录的字查不到就返回 NaN。类比概念类比GF store乘法口诀表建库 (fomosto build)背口诀表费时间查表计算 (engine.process)遇到 7×8 直接查表秒出 56distance_max 太小口诀表只背到 5问 7×8 就不会了2.INDEX_OUT_OF_BOUNDS— 300 km 长的断层最远端的子补丁到接收器的距离超过了 500 km store 上限。比如断层一端离震中 150 km网格角点离震中 ~450 km加起来 600 km 500 km。两个方案选一个方案 A推荐不重建 store把断层拆成多段把 300 km 的断层拆成 3 段 × 100 km每段独立计算再叠加。300 km 拆成 3×100 km 段每段最远子补丁离中心 ~50 km。震中最远的段偏移 ±100 km加上网格对角 ~450 km最大源-接收器距离 ≈ 10050450 600 km还是有点紧。稳妥起见20 km 断层 现有 50km store 刚好够用。改成 2 段分析一下距离2×100 km 断层中心在震中向南北各延伸 100 km网格角点离震中最远 ~454 km角点到最近断层子补丁网格南角到断层南端 ≈ 454 - 10 354 km✅ 500 km之前用的是离震中的距离454100554但 Pyrocko 实际查的是接收器到最近子补丁的距离。2 段 × 100 km 应该够了。直接跑python3 compute_deformation_myanmar.py3.Pyrocko 用的是源中心到接收器的距离不是最近子补丁。段中心偏移 -100km 后网格角点到段中心的距离就 ~490 km逼近 500 km 上限离散化后就超了。最简单的修复缩小网格到 ±2.5°对角 ~380 km不用重建 store。你的文件名是1.2.py我改的是工作区版本。在你的文件里直接改这两行 GRID_LON_RANGE 2.5 GRID_LAT_RANGE 2.5然后跑 bash python3 compute_deformation_myanmar1.2.py或者如果你想保留 ±3° 的大范围重建 store把distance_max改大到7000选一个| 方案 | 操作 | 网格范围 | |------|----------| |A快| 改你文件里两行GRID_LON_RANGE 2.5/GRID_LAT_RANGE 2.5| ±2.5° ≈ 280km | |B完整| 重建 store把 build 脚本里distance_max从50000改成7000再跑bash build_gf_myanmar.sh| ±3° ≈ 330km |推荐先试A马上能出结果。4.问题清楚了。来看两个文件的关键问题1.build_gf_myanmar.sh:distance_max500000(500km) — 但断层本身长 500km 网格延伸远端距离会超出2.compute_deformation_myanmar.py: 只有 2 段 × 100km 200km但你说实际断层 ~500km直接改验证一下修改结果build_gf_myanmar.sh里的实际参数还是 5000需要改现在确认最终修改都正确改好了。总结一下修改内容修改汇总build_gf_myanmar.sh| 参数 | 改前 | 改后 | 原因 | |------|------| |distance_max| 500 (500km) |70000 (700km)| 500km 断层 网格余量最远段到网格边缘可达 ~440km |compute_deformation_myanmar.py| 参数 | 改前 | 改后 | 原因 | |------|------| |N_SEGMENTS| 2 |5| 匹配 ~500km 断层总长 | |GRID_LON_RANGE| 3.0° |4.0°| 500km 断层需要更大网格 | |GRID_LAT_RANGE| 3.0° |4.0°| 同上 | |GRID_SPACING| 0.05° |0.025°| 大断层需要更高分辨率 |使用步骤# 1. 把修改后的文件复制到你的项目录 cp build_gf_myanmar.sh /home/lx/202402351/psgrn_pscmp_experiment/ cp compute_deformation_myanmar.py /home/lx/202402351/psgrn/ # 2. 重新构建 GF store必须重新构建distance_max 从 50→70km cd /home/lx/2024023551/psgrn_pscmp_experiment bash build_gf_myanmar.sh # 3. 运行变形计算 cd /home/lx/202402351/psgrn python3 compute_deformation_myanmar.py⚠️必须重新构建 GF store因为distance_max从 50km 改到了 70km旧的 store 数据不够。构建时间会比之前长一些网格点更多。5.还是同样的错误。问题不是简单的距离超限让我深入分析 pyrocko 的源码逻辑。核心在于pyrocko 的RectangularSource会自动将断层离散化成许多子断层sub-patches子断层点的坐标会超出名义中心位置。再加上multilinear插值需要查询相邻网格点所以即使中心距离在范围内子断层边缘也可能越界。让我做一个更彻底的修复找到根本原因了让我重新计算距离第1段中心 (20.21°N) 到网格西北角 (91.94°E, 26.01°N): dlat 5.8° → 64 km dlon 3.97° (at avg lat ~23.1°) → ~409 km 实际距离 ≈ sqrt(644² 409²) ≈ 763 km !! ← 超过 700km 了之前只算了正北方向644km没算斜角距离_max100km 才够。已经在build_gf_myanmar.sh里改好了。现在给 compute 脚本加上诊断和安全处理Now add a distance diagnostic before eachengine.processcall:验证最终文件好所有修改完成。总结这次修复根本原因之前只算了「正北方向」的距离 644km但实际最远距离是到网格斜角texttext第1段(20.21°N) → 网格西北角(91.94°E, 26.01°N): √(64² 409²) ≈ 763 km 700km ✗最终修改| 文件 | 参数 | 改前 | 改后 | |------|------| |build_gf_myanmar.sh|distance_max| 500km |100km| |build_gf_myanmar.sh|distance_delta| 3km |5km(大范围用粗网格减少存储量) | |build_gf_myanmar.sh|source_depth_delta| 5km |2.5km(更精细) | |compute_deformation_myanmar.py| 诊断 | 无 |每段打印实际最大距离|使用bashbash# 1. 复制到你的目录 # 2. 必须重新构建 GF store cd /home/lx/2024023551/psgrn_pscmp_experiment bash build_gf_myanmar.sh # 3. 运行 cd /home/lx/2024023551/psgrn python3 compute_deformation_myanmar2.py这次构建会更久distance_max 从 50→100km点数翻倍但跑起来后每段会先打印实际距离方便确认是否在 store 范围内。