这次我们来看一个很有意思的模型。NASA 最新发布的 3D 模型不是某个 AI 绘画工具也不是本地视频生成模型而是一个把地球真实形状“变形放大”之后的三维可视化地球表面并不是完美球体而是一个凹凸不平、看起来像“土豆”的大块头。因为这个特征外媒几乎都用同一个词来形容它——Lumpy Mess。这个模型的本质是“大地水准面”可视化英文叫 Geoid。它是从 GRACE、GOCE 等卫星重力观测数据反演出来的一张全球等势面代表地球重力场在“平均海面”意义下的真实起伏。NASA 把它做成 3D 模型后地球看起来不再是光滑的椭球而是一个有隆起、有凹陷、低处坑洼、高处鼓包的不规则表面。这篇文章会围绕这个模型做一次技术拆解数据从哪来、用什么格式、怎么在本地读取、如何在 Python 里生成 3D 可视化、能做哪些验证和批处理、以及实际工作时最容易踩的坑。如果你平时处理遥感数据、做地理可视化或者想研究卫星重力数据这篇可以直接收藏。1. 核心能力速览先把这个模型的关键信息整理成一张表。需要说明的是NASA 发布页面的具体版本、分辨率会随时间更新下面的描述以常见公开版本为基础具体参数请以数据源页面为准。能力项说明模型类型地球重力场 / 大地水准面三维模型核心数据来源GRACE、GRACE-FO、GOCE 等卫星重力观测结合地面和海洋测量数据主要功能展示大地水准面起伏、全球重力异常、地球三维不规则形状数据格式netCDF、GeoTIFF、txt / csv xyz 网格部分版本提供 OBJ/GLTF 网格常见分辨率1 度、0.5 度、0.25 度、0.1 度视发布版本而定推荐处理工具Python、GMT、QGIS、ParaView、PyVista三维可视化方式将高度或重力异常值做垂直夸张后映射到椭球面上是否支持批量任务支持可批量处理多时相重力场数据并合成动画是否支持接口调用可自行封装 REST API官方主要以数据下载为主适合场景地球科学展示、测绘基准、海洋面研究、卫星轨道计算、地理教学从材料看这个项目不是“双击启动”的本地软件而是一整套公开数据和处理流程。它的价值在于你不需要自己发射卫星也能用公开重力场模型在本地重现 NASA 的可视化效果。2. 适用场景与使用边界这类模型适合谁第一类是做地球物理和大地测量的人。重力场模型是确定全球高程基准、研究地壳结构、反演地下密度分布的基础数据。第二类是搞遥感影像或 GIS 开发的人。Geoid 数据可以用来做高程转换比如 GPS 测到的是椭球高真实海拔需要减去大地水准面差距。第三类是科普和技术展示。NASA 这种“土豆地球”可视化适合做互动界面、教学课件或数据大屏。但它也有非常明确的使用边界。最重要的边界是这个凹凸不平的表面不是真实地形。真实地球表面有山脉、海沟、大陆那是地形学概念。而 Geoid 是一个重力等势面可以理解为“全球静止海平面穿越大陆延伸后形成的形状”。它和参考椭球之间的差异通常在正负一百米左右肉眼根本看不出来。NASA 的 3D 模型为了让观众“看见”这种差异会把高度放大几百到几千倍。所以你看到的是夸张后的重力场示意不是地球真实比例。另一个边界是数据合规。NASA 公开的重力场模型大多是开放共享的但它依赖的局部精细数据、商业数据可能受许可约束。使用时要注意数据源页面给出的引用规范和授权范围。如果涉及本国高精度局部重力数据更要遵守数据管理政策不能随意传播或使用。3. 关键技术原理卫星重力数据是怎么变成“土豆地球”的先把原理讲清楚后面的代码才有意义。地球的重力场可以用一组球谐系数来描述类似傅里叶展开。球谐展开的阶数越高能表达的空间细节越丰富。低阶系数描述地球的整体扁率和大尺度质量分布高阶系数刻画山脉、海沟、密度异常引起的重力变化。GRACE 和 GRACE-FO 卫星主要是通过测量两颗卫星之间的距离变化反演地球重力场的长期变化特别适合监测地表水、冰盖、地下水等质量迁移。GOCE 卫星则通过测量重力梯度得到更高分辨率的静态重力场。把这些数据融合再加上地面重力测量、海洋测高数据就能构建出全球大地水准面模型。模型里面真正展示的表面通常是一个规则经纬度网格。每个网格点上有三个关键数值经度 longitude纬度 latitude大地水准面高度 geoid height单位一般是米Geoid height 表示相对于参考椭球面的起伏。正值表示该处重力等势面高于参考椭球负值表示低于参考椭球。在 3D 可视化时把每个经纬度点放到对应椭球位置再沿径向向外或向内偏移一个高度值同时乘以一个夸张系数就得到了我们看到的高低不平地球。NASA 的新模型之所以叫“Lumpy Mess”就是因为它把这种起伏显示得非常明显。从视觉上看赤道附近会鼓出来南亚、东南亚有较明显的隆起印度洋部分区域有凹陷整体看起来确实不太“圆”。4. 环境准备与数据获取在本地复现这套可视化不需要很夸张的硬件。CPU 就能跑内存和显卡取决于数据分辨率和渲染方式。如果只处理 1 度网格4GB 内存足够如果处理 0.1 度全球网格内存需求会明显上升。建议环境如下操作系统Windows 10/11、Ubuntu 20.04 以上、macOS 均可Python 版本3.9 以上核心库numpy、matplotlib、xarray、netCDF4、pyvista、cartopy可选工具GMT、QGIS、ParaView磁盘空间至少预留 10GB用于存放模型文件和处理结果数据获取推荐几个渠道NASA Earthdata需要注册账号下载 GRACE/GRACE-FO 数据ICGEMInternational Centre for Global Earth Models可直接在线计算 Geoid 格网GFZ、CSR、JPL三家机构发布不同版本的 GRACE 时变重力场如果你只想快速验证可以找一个公开的 Geoid 模型文件比如 EGM2008、GOCO 系列、XGM2019 的网格数据。文件扩展名常见为.nc、.grd、.txt、.xyz。安装 Python 依赖pip install numpy matplotlib xarray netCDF4 pyvista cartopy如果读取 netCDF 时遇到时间相关依赖缺失可以补装pip install cftime5. 本地数据处理与三维可视化下面给出一套通用的处理和可视化流程。实际数据字段名、单位、坐标变量可能因文件而异代码里我会标注需要替换的位置。5.1 读取 Geoid 格网数据假设你拿到的是 netCDF 格式用 xarray 打开import xarray as xr import numpy as np # 打开 Geoid 数据实际路径需要替换 ds xr.open_dataset(geoid_model.nc) print(ds) # 取出经纬度和高度数组 # 常见变量名是 geoid、h、z建议先打印变量列表确认 lat np.array(ds[lat]) lon np.array(ds[lon]) geoid np.array(ds[geoid])如果你的数据是文本网格也可以用 pandas 读取import pandas as pd df pd.read_csv(geoid_xyz.txt, delim_whitespaceTrue, headerNone, names[lon, lat, geoid]) lon df[lon].to_numpy() lat df[lat].to_numpy() geoid df[geoid].to_numpy()这里最容易出错的是经纬度顺序。很多地理数据文件遵循“经度、纬度、值”但也有的文件是“纬度、经度、值”。处理前务必确认。5.2 生成二维平面图先画一张二维伪彩图确认数据范围是否正常import matplotlib.pyplot as plt import cartopy.crs as ccrs fig plt.figure(figsize(12, 6)) ax plt.axes(projectionccrs.PlateCarree()) pcm ax.pcolormesh(lon, lat, geoid, cmapRdBu_r, shadingauto) ax.coastlines(linewidth0.5) ax.set_title(Geoid Height (m)) plt.colorbar(pcm, axax, orientationhorizontal, pad0.05) plt.savefig(geoid_2d.png, dpi150)这一步如果正常你会在海面区域看到明显的重力异常带比如印度洋低值区、西太平洋高值区。5.3 生成“土豆地球”三维模型三维可视化的核心思路先把经纬度网格映射到三维球面坐标然后在径向叠加 Geoid 高度值。为了看到形变夸张系数建议从 1000 倍开始。import numpy as np import pyvista as pv # 将经纬度转成弧度 lon_r np.deg2rad(lon) lat_r np.deg2rad(lat) # 使用 WGS84 椭球参数这里简化为平均半径 R 6371000.0 # 网格化经纬度 lon_grid, lat_grid np.meshgrid(lon_r, lat_r) geoid_grid geoid # 确保 geoid 维度与 lon_grid 一致 # 夸张系数可视情况调整 exaggeration 1000.0 # 径向距离 地球半径 夸张后的 geoid 高度 r R geoid_grid * exaggeration # 转笛卡尔坐标 x r * np.cos(lat_grid) * np.cos(lon_grid) y r * np.cos(lat_grid) * np.sin(lon_grid) z r * np.sin(lat_grid) # 构造 PyVista 网格 points np.column_stack([x.ravel(), y.ravel(), z.ravel()]) cloud pv.PolyData(points) # 通过 Delaunay 三维采样构建表面或用结构化网格更高效 grid pv.StructuredGrid() grid.points points grid.dimensions lon_grid.shape # 添加高度值作为标量用于着色 grid.point_data[geoid] geoid_grid.ravel() # 绘制 plotter pv.Plotter() plotter.add_mesh(grid, scalarsgeoid, cmapRdBu_r, show_edgesFalse) plotter.add_axes() plotter.show()这个代码是通用模板。实际使用中结构化网格要求 lon 和 lat 是严格递增的一维数组并且维度顺序与 reshape 一致。如果遇到索引错位建议先用二维伪彩图比对数据方向再做三维映射。5.4 判断可视化是否成功判断标准有三点三维模型整体保持球形但表面有明显起伏。颜色分布与二维伪彩图一致高值区在图上呈现红色隆起低值区呈蓝色凹陷。旋转模型时能看到赤道附近略微膨胀两极方向偏平且局部区域有明显的鼓包。如果模型出现撕裂、飞点先检查经纬度网格是否单调再检查 geoid 值里是否有 NaN 或异常大值。6. 功能测试与效果验证只做出一张图还不够建议按下面几个维度验证数据和处理流程。6.1 与参考椭球对比生成一组基准球面半径等于地球平均半径颜色设成灰色再叠加 Geoid 模型。这样可以直观看到形变大小sphere pv.Sphere(radiusR, theta_resolution200, phi_resolution200) sphere.point_data[geoid_ref] np.zeros(sphere.n_points) plotter pv.Plotter() plotter.add_mesh(sphere, colorlightgray, opacity0.3) plotter.add_mesh(grid, scalarsgeoid, cmapRdBu_r) plotter.show()如果 Geoid 模型确实是夸张后的你会看到红色部分向外鼓出蓝色部分向里塌陷灰色球体夹在中间。6.2 沿特定纬度做剖面为了验证数据不是“乱画出来的”可以沿某一条经度或纬度画剖面曲线。例如取北纬 30 度上的 Geoid 高度横轴是经度纵轴是高度import matplotlib.pyplot as plt # 假设 lat 和 geoid 都是二维数组 # 找到最接近 30°N 的纬度索引 target_lat 30.0 lat_index np.argmin(np.abs(lat - target_lat)) plt.figure(figsize(10, 4)) plt.plot(lon, geoid[lat_index, :]) plt.xlabel(Longitude) plt.ylabel(Geoid Height (m)) plt.title(fGeoid Profile at Latitude {target_lat}N) plt.grid(True) plt.savefig(geoid_profile.png, dpi150)剖面图应该呈现连续、有物理意义的起伏而不是随机噪声。6.3 批量处理多时相数据GRACE 数据是按月发布的如果你拿到多个时期的 Geoid 或重力异常网格可以批量渲染成动画观察地表质量迁移。批量处理脚本的关键是循环遍历文件并固定投影参数、色标范围否则动画会出现颜色跳动。for f in geoid_*.nc; do python render_frame.py --input $f --output frame_${f%.nc}.png done生成 PNG 帧后用 ffmpeg 合成视频ffmpeg -framerate 10 -i frame_%04d.png -c:v libx264 -pix_fmt yuv420p geoid_animation.mp4批量任务最容易出现的问题是部分数据缺失、NaN 值、坐标范围不一致。建议在脚本里加一个“数据完整性检查”步骤先统计每个文件的 min、max、mean再做渲染。7. 接口 API 与任务服务化如果你的目标不是自己做一张图而是把 Geoid 数据接入业务系统可以考虑把数据处理模块封装成 REST API。下面是一个 FastAPI 示例骨架提供“查询指定窗口内 Geoid 统计值”的接口。实际部署时需要替换数据加载逻辑。from fastapi import FastAPI, Query import xarray as xr import numpy as np app FastAPI() # 全局加载一次避免每次请求都读磁盘 ds xr.open_dataset(geoid_model.nc) app.get(/geoid/stats) def geoid_stats( lon_min: float Query(...), lon_max: float Query(...), lat_min: float Query(...), lat_max: float Query(...) ): # 用 xarray 切片 region ds.sel( lonslice(lon_min, lon_max), latslice(lat_min, lat_max) ) geoid region[geoid].values return { min: float(np.nanmin(geoid)), max: float(np.nanmax(geoid)), mean: float(np.nanmean(geoid)), std: float(np.nanstd(geoid)), count: int(np.count_nonzero(~np.isnan(geoid))) } if __name__ __main__: import uvicorn uvicorn.run(app, host127.0.0.1, port8000)启动接口服务uvicorn main:app --host 127.0.0.1 --port 8000调用接口curl http://127.0.0.1:8000/geoid/stats?lon_min80lon_max120lat_min10lat_max50返回结果示例{ min: -32.5, max: 48.2, mean: 12.6, std: 18.4, count: 120000 }需要注意上面的返回数字只是演示实际数值取决于你加载的模型版本和区域。接口服务化以后可以继续扩展切片下载、PNG 渲染、批量任务队列等功能。批量任务建议加任务 ID、日志文件和失败重试机制不要在主进程里跑长时间渲染。8. 资源占用与性能观察这个模型和常见的 AI 推理不同它主要吃内存、磁盘和一定的 GPU 渲染资源。数据占用的估算公式很简单全球网格点数等于纬度点数乘以经度点数。比如 1 度网格是 180×36064800 个点0.5 度是 360×720259200 个点0.25 度是 720×14401036800 个点。每个点如果是 float64占 8 字节100 万个点大概是 8MB。乍一看不大但如果有多个时间步、多个变量加上经纬度坐标和属性整体数据量会成倍增加。三维渲染时PyVista 会把网格点构造成曲面。顶点数到几十万时一般显卡都能流畅旋转如果到了几百万点建议先做下采样否则旋转会卡顿显存占用也会明显上升。降低资源占用的常用方法下载或加载数据前先裁剪目标区域把 float64 转成 float32三维渲染时降低球体分辨率剖面图用二维画不需要开三维窗口批量渲染时关闭显卡抗锯齿降低截图分辨率观察资源占用Windows 可以用任务管理器Linux 可以用htop和nvidia-smi。如果显存不足优先减网格点数而不是换显卡。9. 常见问题与排查方法下面把最容易遇到的问题列成排查表。问题现象可能原因排查方式解决方案netCDF 文件打开失败缺少 netCDF4 或 cftime打印报错信息安装 netCDF4、cftime经纬度维度顺序不对数据文件是 lat/lon代码按 lon/lat 读取查看变量 shape转置数组或调整索引三维模型出现飞点Geoid 值含 NaN 或极大异常值检查 min、max、NaN 数量清洗数据插值或剔除异常值球体表面出现撕裂经纬度网格不是严格递增打印 lat、lon 数组按升序重排网格颜色标尺不明显Geoid 范围太小夸张系数不足查看数据 min/max增大夸张系数或调整配色范围显存不足三维网格点太多查看显卡占用下采样、裁剪区域、减少球体精度批量任务卡住单个文件读取太慢或内存溢出打印日志定位到具体文件增加超时控制、任务日志、跳过坏文件剖面曲线不连续网格中间有缺测值检查 NaN 分布插值或使用掩膜最值得关注的是“坐标方向错误”。地理数据里lat有时从南到北有时从北到南。如果生成的三维模型高低分布和官方图相反大概率是纬度方向反了。先画二维图对比再进三维能省很多时间。10. 最佳实践与使用建议做完这套可视化流程后面再碰到类似的重力场、地形、海洋测高数据都可以沿用同一套思路。第一数据版本要记录清楚。模型文件名称、来源机构、发布时间、网格分辨率、单位、夸张系数最好写进一个 config 文件或 README。别人拿到你的代码和图片时能知道你到底画的是什么。model: name: XGM2019-like geoid source: ICGEM / NASA Earthdata resolution_deg: 0.5 vertical_exaggeration: 1000 unit: meter output: format: png dpi: 150 colorbar: RdBu_r第二三维可视化和平面的“平面图”要配套使用。平面图适合看全球重力异常分布三维模型适合看整体形变感受。发布或汇报时两者结合更有说服力。第三不要轻视合规问题。全球公开重力模型作为开放数据可以用于研究和教学但商用前要确认许可条款。如果你的工作涉及高精度局部重力数据务必遵守主管部门的使用规定不要随意扩散。第四涉及可视化成果时要在图注里明确“高度经过放大展示”。否则读者可能会误以为地球真实形状就是那么凹凸不平这在科普上是误导。第五工程化批量任务时一定要给每个任务加日志。GRACE 数据一个月一个文件批量处理几十个月的数据很容易中途失败。记录每个文件的处理状态、输出路径、异常原因是稳定跑完批量任务的关键。11. 总结与下一步这个 NASA 新 3D 模型最值得尝试的地方并不是它本身有多炫酷而是它把抽象的重力场变成了肉眼可见的形状。只要学会读取 Geoid 网格、映射到三维球面、加一个垂直夸张系数你完全可以复现类似的“土豆地球”。建议先从 1 度或 0.5 度公开网格开始用 matplotlib 画二维图再用 PyVista 生成三维模型。先验证数据方向和数据范围再逐步加入夸张系数、剖面曲线、批量动画和 API 服务。最容易踩的坑有两个一是把 Geoid 当成真实地形二是在三维映射时忽略经纬度顺序。这两个问题想清楚整套流程基本不会卡住。后续可以继续向两个方向扩展一个是时变方向用 GRACE/GRACE-FO 月重力场做时间序列动画观察地下水、极地冰盖变化另一个是拟合方向把 Geoid 模型和全球地形、卫星测高、GNSS 高程数据做联合分析解决高精度高程转换问题。