从 WRF 投影质量网格点并在 R 中添加形状文件

如何解决从 WRF 投影质量网格点并在 R 中添加形状文件

我想将 shapefile 添加到 WRF 地理文件和 WRF 输出文件的绘图中，但我真的很难让投影正确。

我尝试了几种网上找到的方法，例如帖子 here 和 here，以及星星包和 st_transform 以及 eixport 包，但没有任何工作正常。我使用 this code in slighly modified form 作为函数 get_wrf_raster 从 WRF 地理文件中获取投影。输入文件和函数可以在here中找到。

library(ncdf4)
library(rgdal)
library(raster)         

source('WRFnc_to_raster.R') # for function get_wrf_raster

# read shapefile
shp <- readOGR(dsn=".",layer="boundary")
extent <- c(-139.27,-138.97,60.61,60.91)
shp_crop <- crop(shp,extent)

# ORIGINAL ELEVATION DATA (ASTER-DEM)
###############################################
orig <- raster("orig.tif")
orig_crop <- crop(orig,extent)
plot(orig_crop)
plot(shp_crop,add=T)

original_data_in_ll --> Shapefile 匹配高程数据

# PROJECTION TO LCC FROM WRF GEOFILES (LAmberT CONFORMAL CONIC)
###############################################
geo <- get_wrf_raster('geo_em.d04.nc',"HGT_M")
crs_lcc <- crs(geo)
shp_lcc <- spTransform(shp,crs_lcc)
plot(geo)
lines(shp_lcc)

wrf_geofile_in_lcc --> 错

# PROJECTION TO LL FROM WRF GEOFILES (LATITUDE LONGITUDE)
###############################################
crs_ll <- crs(orig)
geo_ll <- projectRaster(geo,crs = crs_ll,method = "ngb")
geo_ll_crop <- crop(geo_ll,extent)
plot(geo_ll_crop)
lines(shp_crop)

wrf_geofile_in_ll --> 错

# PROJECTION FROM WRF OUTFILE
###############################################
wrf <- raster('wrfout.d04.nc',varname="HGT")
# projection from namelist information
crs(wrf) <- "+proj=lcc +lat_1=50 +lat_2=70 +lat_0=60.75889 +lon_0=-139.1246 +datum=wgs84 +ellps=sphere +a=6370000 +b=6370000 +units=m +no_defs"
extent(wrf) <- extent(geo) # ???

1) 你能帮我把投影正确吗？问题是精度、质量网格点（WRF 地理文件中的“HGT_M”）还是我的方法？

2) 如何直接从 wrf 输出文件中获取范围，而不使用地理文件中的角纬度/经度？我是否必须以不同的方式对待 WRF 输出文件（交错）中的“HGT”？

任何帮助将不胜感激。非常感谢！