使用PyGMT与Pandas绘制月球南极高程图的问题排查

问题背景 #

• 数据：日本SELENE（昵称Kaguya）月球探测器南极高程数据，覆盖经度0-360°、纬度约-79°至-90°，已完成清理并导入Pandas DataFrame
• 异常现象：使用PyGMT绘制的等距方位投影（E0/-90/20c）高程图，中心区域及边缘出现空白灰色区域

原始绘图代码 #

def plotGrid(df):
    startTimeOA = time.time()
    grid = pygmt.xyz2grd(
        data=df,
        projection="E0/-90/20c",
        region='0/360/-90/-80',
        spacing=[0.03125, 0.0078125]
    )
    print(f'Took {time.time() - startTimeOA:06.3f} seconds to create grid from relevant Pandas dataFrame.')
    startTime = time.time()
    fig = pygmt.Figure()
    fig.grdimage(grid=grid)
    fig.show()
    print(f'Took {time.time() - startTimeOA:06.3f} seconds overall to plot Lunar South Pole Map '
          f'from relevant LALT data.')
    print('Plotted Long Lat Alt elevation data in pyGMT for relevant LALT data.\n')
    return fig

数据样例（经度/纬度/高程，单位：度/度/千米） #

0.015625 -79.00390625 2.114
0.046875 -79.00390625 2.123
0.078125 -79.00390625 2.133
0.109375 -79.00390625 2.144
0.140625 -79.00390625 2.160
0.171875 -79.00390625 2.177
0.203125 -79.00390625 2.195
0.234375 -79.00390625 2.213
0.265625 -79.00390625 2.230

问题原因及解决方案 #

1. 核心问题分析 #

• 区域与数据不匹配：region='0/360/-90/-80'包含了无数据的南极点（-90°）区域，且0-360°经度在极地投影下会出现衔接空白
• 网格间距方向错误：原始数据纬度从-79°向-90°递减，但spacing中纬度增量设为正值，导致网格生成时无法匹配数据点，出现径向空白条带
• 无插值填充：未设置插值方法，数据间隙直接显示为空白

2. 修复后的代码 #

import pygmt
import time

def plotGrid(df):
    startTimeOA = time.time()
    # 将0-360°经度转换为-180-180°，适配极地投影的经度衔接
    df["longitude"] = df["longitude"].apply(lambda x: x - 360 if x > 180 else x)
    grid = pygmt.xyz2grd(
        data=df,
        region='-180/180/-90/-79',  # 贴合实际数据覆盖范围
        spacing=[0.03125, -0.0078125],  # 纬度增量设为负值，匹配数据递减方向
        method="linear"  # 线性插值填充数据间隙
    )
    print(f'Took {time.time() - startTimeOA:06.3f} seconds to create grid from relevant Pandas dataFrame.')
    startTime = time.time()
    fig = pygmt.Figure()
    fig.grdimage(grid=grid, projection="E0/-90/20c")
    # 可选：添加颜色条，增强可读性
    fig.colorbar(position="JMR+o0.5c/0c+w10c/0.5c", frame=["x+lElevation", "y+lkm"])
    fig.show()
    print(f'Took {time.time() - startTimeOA:06.3f} seconds overall to plot Lunar South Pole Map '
          f'from relevant LALT data.')
    print('Plotted Long Lat Alt elevation data in pyGMT for relevant LALT data.\n')
    return fig

3. 关键修复说明 #

• 经度转换：解决0°/360°衔接处的空白问题，让投影下的经度数据连续
• 区域调整：贴合实际数据的纬度范围，减少无数据区域带来的空白
• 间距方向修正：纬度增量设为负值，确保网格点与数据点的纬度方向一致，消除径向空白条带
• 插值填充：method="linear"让PyGMT自动填充数据间隙，消除局部空白

内容的提问来源于stack exchange，提问作者Python Programmer