在 R 中屏蔽特定区域的 NetCDF 文件时出错?

Getting error while masking the NetCDF file for a particular region in R?

我编写了以下脚本来根据我的 Shapefile 剪辑 netCDF 文件。 NetCDF 文件包含土壤信息。当我 masking 使用 shapefileNetCDF 文件时,出现错误。任何有关如何解决该问题的建议将不胜感激。我的 netCDF 文件很重,这就是为什么我没有在这里附加它所以我将 link 共享到 NetCDF 文件和 Shapefile GoogleDriveLink.

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

NC_File <- "ADD_PROP1_NCRB.nc"
print(nc_open(NC_File))
 

NetCDF 文件的描述

1 variables (excluding dimension variables):
        byte ADD_PROP[lon,lat]   
            long_name: additional property 
            _FillValue: -100
            missing_value: -100

     2 dimensions:
        lon  Size:3377
            standard_name: longitude
            long_name: longitude
            units: degrees_east
            axis: X
        lat  Size:1725
            standard_name: latitude
            long_name: latitude
            units: degrees_north
            axis: Y

使用 brick 函数进行栅格转换

b <- brick(NC_File)
b
> b 
class      : RasterBrick 
dimensions : 1725, 3377, 5825325, 1  (nrow, ncol, ncell, nlayers)
resolution : 0.00833333, 0.00833333  (x, y)
extent     : -117.6917, -89.55003, 45.31663, 59.69162  (xmin, xmax, ymin, ymax)
crs        : +proj=longlat +datum=WGS84 +no_defs 
source     : G:/Nelson_MIP/ForcingFilesFromOneDrive/Soil/GSDE_NCRB/ADD_PROP1_NCRB.nc 
names      : layer 
varname    : ADD_PROP

正在读取我的 shapefile 以屏蔽 NetCDF

SHP <- readOGR("G:/Nelson_MIP/WatershedFile/ForSoilNetCDFprocessing.shp")
SHP
> SHP
class       : SpatialPolygonsDataFrame 
features    : 1 
extent      : -103.7103, -101.7075, 50.83711, 52.35122  (xmin, xmax, ymin, ymax)
crs         : +proj=longlat +datum=NAD83 +no_defs 
variables   : 1
names       : Strahler_O 
value       :          5

将 shapefile 的投影与 NetCDF 文件的投影相匹配

SHP <- spTransform(SHP, crs(NC_File))

使用我的 shapefile 剪辑 NetCDF 文件

Masking <- mask(b, SHP)
> Masking <- mask(b, SHP)
Error in ncvar_get_inner(ncid2use, varid2use, nc$var[[li]]$missval, addOffset,  : 
  Error: variable has 2 dims, but start has 3 entries.  They must match!

这是我在 rasterdevelopment version 中修复的错误。 使用当前版本,您可以使用 raster 而不是 brick 来解决它,因为这是一个单层数据集

library(raster)
b <- raster("ADD_PROP1_NCRB.nc")
s <- shapefile("ForSoilNetCDFprocessing.shp")
s <- spTransform(s, crs(b))
m <- mask(b, s)