用R语言从NetCDF文件中提取数据并找到正确的投影

Question

我在这里下载了一个netcf文件，并通过:

library(raster)
b <- brick("precipitationintensity14rcp45modelsensemblemedianannual.nc")

然后我检查了proj4string为了看到投影和如何提供坐标:

> proj4string(b)
[1] "+proj=lcc +lat_0=47.5 +lon_0=13.3299999237061 +lat_1=49 +lat_2=46 +x_0=400000 +y_0=400000 +datum=WGS84 +units=m +no_defs"

为了提取某一点的数据，我使用WGS84格式的长/晚坐标:

coords<-data.frame(lon=47.4, lat=11.5)
coordinates(coords)<-c("lon","lat")

然后我尝试提取数据:

> extract(b$X1981.07.15,coords)
[,1]
[1,]   NA

这让我莫名其妙。因为我期望在这个给定点有一个值。

为了调试我的代码，我尝试比较光栅对象的坐标与我的:

> coordinates(b$X1981.07.15)
x      y
[1,] 111500 571500
[2,] 112500 571500
[3,] 113500 571500

现在我很困惑，因为提供的坐标不匹配proj4string格式(WGES84)。你能帮我一下吗?

1. 为提取数据提供合适的坐标
2. 解释为什么光栅对象中的坐标不像我预期的那样

Answer 1

您需要将点转换为栅格的坐标参考系(而不是反过来!)下面我将演示如何使用"地形"来替换"栅格"。包。

library(terra)
b <- rast("precipitationintensity14rcp45modelsensemblemedianannual.nc")
coords <- data.frame(lon=47.4, lat=11.5)
pts <- vect(coords, crs="+proj=longlat")
#crs(b) = "+proj=lcc +lat_0=47.5 +lon_0=13.3299999237061 ...
lccpts <- project(pts, crs(b))
extract(b$X1981.07.15, lccpts)