关于上期葵花卫星绘图补充

用户11172986

发布于 2024-06-20 16:36:22

670

发布于 2024-06-20 16:36:22

文章被收录于专栏：气python风雨气python风雨

*封面图片由ai生成

昨日有读者提到葵花卫星图是使用imshow绘制的，能不能添加经纬度刻度？小编测试一下是可行的。代码与效果如下

增加经纬度刻度

import cartopy.crs as ccrs
import matplotlib.pyplot as plt

dat = ds2.metpy.parse_cf('channel_0001_scaled_radiance')
geos = dat.metpy.cartopy_crs
x = dat.x
y = dat.y

fig = plt.figure(figsize=(15, 12))
ax = fig.add_subplot(1, 1, 1, projection=geos)
ax.imshow(rgb, origin='upper', extent=(x.min(), x.max(), y.min(), y.max()), transform=geos)

# 添加经纬度网格线
ax.gridlines(draw_labels=True, dms=True, x_inline=False, y_inline=False)

# 添加海岸线和边界
ax.coastlines(resolution='50m', color='white', linewidth=1)
ax.add_feature(ccrs.cartopy.feature.BORDERS, linewidth=1)

plt.title('True Color', loc='left', fontweight='bold', fontsize=15)

plt.show()

只需添加ax.gridlines函数即可。

显示局部地区的经纬度刻度的代码效果如下

#获取投影
dat = ds2.metpy.parse_cf('channel_0001_scaled_radiance')

geos = dat.metpy.cartopy_crs
#转换RGB
vmin = 0
vmax = .5
B1 = b1[1000:4000,1000:3500].clip(vmin, vmax) / vmax * 255
B2 = b2[1000:4000,1000:3500].clip(vmin, vmax) / vmax * 255
B3 = b3[1000:4000,1000:3500].clip(vmin, vmax) / vmax * 255

X = x[1000:4000]
Y = y[1000:3500]
rgb = np.stack((B3, B2, B1), axis=2).astype('uint8')

fig = plt.figure(figsize=(20,12))
ax1 = fig.add_subplot(121, projection=geos)
ax2 = fig.add_subplot(122, projection=geos)
# Plot image
ax1.imshow(rgb,extent=(X[0], X[-1], Y[-1], Y[0]))

ax2.imshow(rgb, origin='upper',extent=(X[0], X[-1], Y[-1], Y[0]), transform=geos)

# Add Coastlines and States
ax2.coastlines(resolution='50m', color='white', linewidth=1)
ax2.add_feature(ccrs.cartopy.feature.STATES, linewidth=1)
gl = ax1.gridlines(draw_labels=True, alpha=0.5, linewidth=1, color='w', linestyle='--')

gl.xlabels_top = False  # 关闭顶端标签
gl.ylabels_right = False  # 关闭右侧标签
gl.xformatter = LONGITUDE_FORMATTER  # x轴设为经度格式
gl.yformatter = LATITUDE_FORMATTER  # y轴设为纬度格式
plt.title('Himawari 8 \n', fontsize=20)

# Adjust tick mark size
plt.tick_params(labelsize=16)