在为R中的特定区域屏蔽NetCDF文件时遇到错误？

如何解决在为R中的特定区域屏蔽NetCDF文件时遇到错误？

我有以下编写的脚本，用于根据我的netCDF来裁剪Shapefile文件。 NetCDF文件包含土壤信息。当我使用masking来NetCDF shapefile个文件时，出现错误。任何有关如何解决该问题的建议将不胜感激。我的netCDF文件很重，这就是为什么我没有在这里附加它，所以我共享了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

使用积木函数进行栅格转换

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,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!

解决方法

这是我在raster的{{3}}中修复的错误。在当前版本中，您可以使用raster而不是brick来解决它，因为这是一层数据集

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