Matplotlib Contour/Contourf 在北极/日期变更线的 Cartopy 奇点

Question

我正在尝试生成显示 RF link 到北极上方卫星的大气衰减值的热图，但我对 Matplotlib contour/contourf 函数完成的插值有疑问。

contourf 函数完成的线性插值在 N.Pole 周围效果不佳，因为我怀疑它不知道在从（-180 度到 +180 度）的值之间进行插值 - 即越过日期变更线，或越过极点。

有什么关于生成热图的不同方法的建议，以避免中心出现这个可怕的洞？！

下面的代码生成绘图。

import cartopy.crs as ccrs
import cartopy.feature


plt.figure(figsize=(10,10))

# Initialise Cartopy Axes.
proj=ccrs.LambertAzimuthalEqualArea(central_longitude=0, central_latitude=90)

ax = plt.axes(projection = proj)

ax.set_extent([-180,180,45,90], ccrs.PlateCarree())

ax.add_feature(cartopy.feature.LAND)
ax.add_feature(cartopy.feature.OCEAN)
ax.add_feature(cartopy.feature.COASTLINE)
ax.add_feature(cartopy.feature.BORDERS, linestyle=':')
ax.gridlines(ls=":",color="grey",lw=0.5)


x0,x1 = attenuation_df.lon.min(), attenuation_df.lon.max()
y0,y1 = attenuation_df.lat.min(), attenuation_df.lat.max()

x,y = np.linspace(x0,x1,1000), np.linspace(y0,y1,1000)
X,Y = np.meshgrid(x,y)

Z = scipy.interpolate.griddata(
    attenuation_df[["lon","lat"]],
    attenuation_df["attenuation"],
    (X,Y),
    method="linear",
)

plt.contourf(X,Y,Z,transform=ccrs.PlateCarree(),alpha=0.5)
plt.colorbar(shrink=0.5)

plt.title("Attenuation")

plt.show()

Attenuation_df 是一个 Pandas 数据框，其中包含大约 3500 个采样点的衰减值，这些采样点在全球范围内等距分布。这是样本点的位置：

这里是 attenuation_df 的 header:

	lon	lat	attenuation
0	-30.8538	48.8813	0.860307
1	-29.0448	49.5026	0.783662
2	-27.2358	50.1317	0.720165
3	-32.6628	48.2676	0.947662
4	37.4226	46.0322	0.27495

attenuation_df 的 csv link 在这里：https://pastebin.com/NYA1jFgt

Answer 1

一个解决方案是将您的数据重新投影到不同的坐标系，我的建议是使用 Polar Stereographic system。但是，以北极为中心的大“洞”不是来自正在使用的坐标系，而是数据集中存在一些 nan，因此您首先必须删除这些值。

这是一个可行的解决方案：

from pyproj import Proj

# Define a pyproj function to reproject data
def coordinate_conv(x, y, inverse = True):
    p = Proj('+proj=stere +lat_0=90 +lat_ts=70 +lon_0=-45 +k=1 +x_0=0 +y_0=0 +a=6378273 +b=6356889.449 +units=m +no_defs')
    return p(x, y, inverse = inverse)

# Drop null values
attenuation_df.dropna(how = 'any', inplace = True)

# Reproject data
rpjx, rpjy = coordinate_conv(attenuation_df.lon, attenuation_df.lat, False)
rpj_cord = pd.DataFrame({'x': rpjx, 'y': rpjy})

# Interpoolate data
x,y = np.linspace(rpjx.min(),rpjx.max(),1000), np.linspace(rpjy.min(),rpjy.max(),1000)
X,Y = np.meshgrid(x,y)

Z = interpolate.griddata(
    rpj_cord,
    attenuation_df["attenuation"],
    (X,Y),
    method="linear",
)

# Figure
plt.figure(figsize=(10,10))

# Initialise Cartopy Axes.
proj=ccrs.LambertAzimuthalEqualArea(central_longitude=0, central_latitude=90)

ax = plt.axes(projection = proj)

ax.set_extent([-180,180,45,90], ccrs.PlateCarree())

ax.add_feature(cartopy.feature.LAND)
ax.add_feature(cartopy.feature.OCEAN)
ax.add_feature(cartopy.feature.COASTLINE)
ax.add_feature(cartopy.feature.BORDERS, linestyle=':')
ax.gridlines(ls=":",color="grey",lw=0.5)

kw = dict(central_latitude=90, central_longitude=-45, true_scale_latitude=70)
plt.contourf(X,Y,Z, transform=ccrs.Stereographic(**kw),alpha=0.5)
plt.colorbar(shrink=0.5)

plt.title("Attenuation")

这是输出图：

Matplotlib Contour/Contourf 在北极/日期变更线的 Cartopy 奇点

Matplotlib Contour/Contourf in Cartopy singularity at North Pole/ dateline

python

matplotlib

geospatial

geo

cartopy