R语言点数据绘图叠加AnnualTemp RasterLayer图层代码求助

你当前代码的核心问题是点数据用的是WGS84经纬度坐标系，你创建的栅格是EPSG:3460投影坐标系，二者坐标系不匹配无法直接叠加，需先统一坐标系、给空栅格插值温度值后再叠加，完整补充后代码如下：

library(raster)
library(gstat)
library(sp)

MinTemp_plots <- data.frame(
  Station = c( "Labasa", "Laucala", "Lautoka", "Levuka", "Matei", "Matuku",
               "Nabouwalu", "Nacocolevu",  "Nadi", "Nausori",  "Ono-I-Lau",  "Penang",
               "Savusavu",   "UduPoint",   "Viwa", "Yasawa"),
  Latitude = c(-16.43, -18.15, -17.6, -17.68, -16.69, -19.15, -16.98,
               -18.1,  -17.75, -18.03, -20.65, -17.37, -16.78, -16.11,
               -17.14, -16.78 ),
  Longitude = c(179.36, 178.45, 177.45, 178.83, 180, 179.76, 178.7, 177.55,
                177.45, 178.56, 178.7, 178.15, 179.34, 180, 176.93,
                177.5),
  AnnualTempMin = c(1.722, 1.711, 0.042, 0.135, 0.264, 0.276, 0.625, 1.215,  
                    1.522, 0.917,  0.617, 0.072, 0.509, 1.057, 1.201, 0.123))

# 1. 将站点数据转为空间点对象，设置原始经纬度坐标系为WGS84(EPSG:4326)
coordinates(MinTemp_plots) <- ~Longitude+Latitude
proj4string(MinTemp_plots) <- CRS("+init=epsg:4326")

# 2. 将站点投影为和栅格一致的EPSG:3460坐标系
MinTemp_proj <- spTransform(MinTemp_plots, CRS("+init=epsg:3460"))

# 3. 创建空栅格
r <- raster(xmn=1790828.61, xmx=2337149.40, ymn=3577110.39, ymx=4504717.19, res=100000)
crs(r) <- crs("+init=epsg:3460")

# 4. 用反距离加权法插值温度到栅格，得到AnnualTemp对应的栅格层
idw_model <- gstat(id = "AnnualTempMin", formula = AnnualTempMin~1, data = MinTemp_proj)
temp_raster <- interpolate(r, idw_model)
# 去掉插值结果冗余层，只保留温度值
temp_raster <- temp_raster[[1]]

# 5. 叠加绘图：先画栅格再加点
plot(temp_raster, main = "年最低温度分布与站点位置")
points(MinTemp_proj, pch = 16, col = "red")

# 如果你要先画点再加栅格，用add参数即可
# plot(MinTemp_proj, pch = 16, col = "red")
# plot(temp_raster, add = TRUE, alpha = 0.7) # alpha设置透明度避免盖住点

补充说明 #

• 如果你已经有现成的AnnualTemp栅格文件，不需要插值的话，直接用temp_raster <- raster("你的栅格文件路径")读取后，用projectRaster转成EPSG:3460坐标系，再按上面的方法叠加即可
• 插值方法可以根据你的需求替换为克里金、样条函数等其他方案

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