基于ggplot2绘制Hatano投影圆形全球栅格图的问题求助

问题根源 #

1. 投影后的栅格未裁剪掉Hatano圆形范围外的无效数据，导致图中出现多余边缘
2. 用base绘图的plot()赋值给myplot，该函数返回值为NULL，因此ggsave生成空白文件

修改后的代码方案 #

方案1：适配原代码逻辑（terra+base绘图） #

library(terra)
library(viridis)
library(graticule)

setwd("C:/Users/usman/Desktop/a")
# 用terra读取栅格，替代raster包更高效
r <- rast("SH.tif")
# 投影到Hatano坐标系
p <- project(r, "+proj=hatano", method="bilinear")

# 创建Hatano投影下的圆形边界（以地球近似半径6371000米为半径）
hatano_crs <- crs(p)
circle <- vect("POLYGON ((0 6371000, 6371000 0, 0 -6371000, -6371000 0, 0 6371000))", crs=hatano_crs)
# 裁剪栅格，仅保留圆形内的数据
p_cropped <- crop(p, circle, mask=TRUE)

# 绘图并直接保存（base绘图需用设备函数）
png("tile_plot2.png", height=4, width=6.5, units="in", res=150)
plot(graticule(60, 45, hatano_crs), background="azure", mar=c(.2,.2,.2,4), lab.cex=0.5, col="light gray")
plot(p_cropped, add=TRUE, axes=FALSE, plg=list(shrink=.8), col=viridis(25))
dev.off()

方案2：ggplot2生态（tidyterra+ggplot2） #

library(terra)
library(ggplot2)
library(viridis)
library(tidyterra)

setwd("C:/Users/usman/Desktop/a")
r <- rast("SH.tif")
p <- project(r, "+proj=hatano", method="bilinear")

# 创建圆形裁剪边界
hatano_crs <- crs(p)
circle <- vect("POLYGON ((0 6371000, 6371000 0, 0 -6371000, -6371000 0, 0 6371000))", crs=hatano_crs)
p_cropped <- crop(p, circle, mask=TRUE)

# 用ggplot绘制并保存
myplot <- ggplot() +
  geom_spatraster(data=p_cropped) +
  scale_fill_viridis(option="viridis", na.value="transparent") +
  geom_graticule(crs=hatano_crs, breaks=c(-60,0,60), lon_breaks=c(-135,-90,-45,0,45,90,135), color="light gray") +
  coord_sf(crs=hatano_crs, expand=FALSE) +
  theme_void() +
  theme(plot.background = element_rect(fill="azure", color=NA))

ggsave("tile_plot2.png", plot=myplot, height=4, width=6.5, dpi=150)

关键说明 #

• 圆形裁剪：Hatano投影的有效区域为圆形，必须用对应投影下的圆形多边形裁剪栅格，清除无效边缘数据
• 保存逻辑：base绘图需用png()+dev.off()的设备调用方式；ggplot绘图可直接用ggsave，因ggplot()返回有效绘图对象

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