`rasterize`函数的问题

Question

我需要使用一些 modis 图像来计算 EVI，我已经使用过这段代码并且它有效，但是现在，当我尝试使用

rasterize

函数时，我收到此错误：

Error in (function (cond)  :  error in evaluating the argument 'x' in selecting a method for function 'is.factor': Not compatible with requested type: [type=character; target=double].

我之前在 QGIS 中分离了 EVI 层，因此我正在使用仅包含 EVI 信息的 .tif 文件。

请参阅下面的代码和图像信息。非常欢迎任何帮助。

install.packages("C:/Users/Anto/Downloads/rgdal_1.4-3.tar.gz", repos=NULL, type='source')

library("raster")
library("rgdal")

setwd("C:/Users/Anto/Desktop/Comportamiento_cortisol/GIS_payunia/UTM19S/all_years_spring")

imagen <- raster("20110914_623_UTM19S.tif")   

shape <-readOGR("polig_N_UTM19S.shp")  
zones <- rasterize(shape, imagen, field = "id")  


zones <- rasterize(shape, imagen, field = "id")  

> imagen
class      : RasterLayer 
dimensions : 5120, 4799, 24570880  (nrow, ncol, ncell)
resolution : 393.1105, 224.7721  (x, y)
extent     : -296537.6, 1590000, 5530378, 6681211  (xmin, xmax, ymin, ymax)
crs        : +proj=utm +zone=19 +south +datum=WGS84 +units=m +no_defs 
source     : 20110914_623_UTM19S.tif 
names      : X20110914_623_UTM19S 
values     : -0.1998, 0.9411  (min, max)

> shape
class       : SpatialPolygonsDataFrame 
features    : 2 
extent      : 469983.2, 539187.8, 5984042, 6016778  (xmin, xmax, ymin, ymax)
crs         : +proj=utm +zone=19 +south +datum=WGS84 +units=m +no_defs 
variables   : 14
names       : id,        area, count,              sum,              mean,            median,             stdev,               min,               max, count_1,            sum_1,            mean_1,             min_1,             max_1 
min values  :  1, 202291562.4,  2293, 162.667400129139,  0.07094086355392, 0.068999998271465, 0.011263931890283,   0.0168999992311, 0.114500001072884,    2293, 162.667400129139,  0.07094086355392,   0.0168999992311, 0.114500001072884 
max values  :  2, 513447724.7,  5810,  428.03230074048, 0.073671652451029, 0.073600001633167, 0.012532315292168, 0.032499998807907, 0.116499997675419,    5810,  428.03230074048, 0.073671652451029, 0.032499998807907, 0.116499997675419 
> 

Answer 1

这是一个对我有用的示例，使用 Sentinel-2 带和我数字化的多边形形状文件。 shapefile 有一个整数字段“id”。

library(terra)
# Polygon shapefile
vect_file <- "MyName_shp.shp"
# Sentinel 2 Band
S2_blue  <- "T37NEB_20190416T072621_B02_60m.jp2"
v <- terra::vect(vect_file)
r <- terra::rast(S2_blue)
v
 class       : SpatVector 
 geometry    : polygons 
 dimensions  : 2, 1  (geometries, attributes)
 extent      : 527979.8, 551753.5, 146122.6, 155954.9  (xmin, xmax, ymin, ymax)
 source      : MyName_shp.shp
 coord. ref. : WGS 84 / UTM zone 37N (EPSG:32637) 
 names       :    id
 type        : <int>
 values      :     1
                   2
r
class       : SpatRaster 
dimensions  : 1830, 1830, 1  (nrow, ncol, nlyr)
resolution  : 60, 60  (x, y)
extent      : 499980, 609780, 90240, 200040  (xmin, xmax, ymin, ymax)
coord. ref. : WGS 84 / UTM zone 37N (EPSG:32637) 
source      : T37NEB_20190416T072621_B02_60m.jp2 
name        : T37NEB_20190416T072621_B02_60m 
v_rasterized <- terra::rasterize(v,r, field = "id")
v_rasterized
class       : SpatRaster 
dimensions  : 1830, 1830, 1  (nrow, ncol, nlyr)
resolution  : 60, 60  (x, y)
extent      : 499980, 609780, 90240, 200040  (xmin, xmax, ymin, ymax)
coord. ref. : WGS 84 / UTM zone 37N (EPSG:32637) 
source(s)   : memory
varname     : T37NEB_20190416T072621_B02_60m 
name        : id 
min value   :  1 
max value   :  2 

我想补充一点，这不会为您提供原始栅格的 EVI 值。相反，光栅化矢量将具有每个像素的 id

 字段的值。如果您想将原始栅格裁剪到 shapefile 多边形的范围内，您需要的是 

terra::crop()

 和 

terra::mask()

。

HTH