为什么不同的栅格索引函数会给我不同的答案?

Why are different raster indexing functions giving me different answers?

我正在尝试在空间 raster 中查找特定位置的索引。但是根据我计算索引的方式,我得到了不同的答案,我无法弄清楚为什么它们不同或哪个是正确的。

请在下面找到一个最小的工作示例:

> library(raster)
> library(tidyverse) 
> 
> ### define raster (grid) & point of interest (poi) 
> lon <- seq(-124.6875, -114.0625, by = 0.625)
> lat <- seq(32.5, 49.0, by = 0.5)
> grid <- expand.grid(lon,lat) %>% 
+     mutate(values = rnorm(nrow(.))) %>% 
+     rasterFromXYZ(crs = "+proj=longlat +datum=WGS84 +no_defs")
> poi <- c(-122.8125, 38.5)
> 
> ### option 1: find row & column indices manually ###
> 
> # row-col index:
> c <- which(lon==poi[1]) #longitude = x = matrix columns
> r <- which(lat==poi[2]) #latitude = y = matrix rows
> c(row=r, col=c)
row col 
 13   4 
> 
> # cell index: 
> cellFromRowCol(grid, row=r, col=c)
[1] 220
> 
> # cell value:
> grid[r,c]
[1] -0.6369638
> 
> ### option 2: use raster functions ###
> 
> # cell index:
> cell <- cellFromXY(grid, xy = poi) 
> cell 
[1] 382
> 
> # row-col index:
> rowColFromCell(grid, cell)
     row col
[1,]  22   4
> 
> # cell value: 
> grid[cell]
[1] -0.8445352
> 
> ## none of these match the results from the manual method! 

就目前而言,Option 1 不正确,Option 2 正确。

Option 1 中的错误来自这样一个事实,即在 raster 中,单元格的编号顺序是从左上角的单元格到右上角的单元格,然后在下一个单元格的左侧继续行,依此类推,直到栅格 lower-right 侧的最后一个像元。因此,如果您想使用 Option 1 并获得正确的结果(即与 Option 2 相同),您需要反转 lat 向量的值(即按降序设置值).

请在下面找到一个代表。

Reprex

library(raster)
library(tidyverse)


### define raster (grid) & point of interest (poi) 
lon <- seq(-124.6875, -114.0625, by = 0.625)
lat <- rev(seq(32.5, 49.0, by = 0.5))        # !!!! Reverse the order of the 'lat' vector
set.seed(1234)                       # !!!!Set the seed for the reprex to be reproducible
grid <- expand.grid(lon,lat) %>% 
  mutate(values = rnorm(nrow(.))) %>% 
  rasterFromXYZ(crs = "+proj=longlat +datum=WGS84 +no_defs")

poi <- c(-122.8125, 38.5)

### option 1: find row & column indices manually ###
 
# row-col index:
c <- which(lon==poi[1]) #longitude = x = matrix columns
r <- which(lat==poi[2]) #latitude = y = matrix rows
c(row=r, col=c)
#> row col 
#>  22   4
 
# cell index: 
cellFromRowCol(grid, row=r, col=c)
#> [1] 382

# cell value:
grid[r,c]
#> [1] -2.651741


### option 2: use raster functions ###
 
# cell index:
cell <- cellFromXY(grid, xy = poi) 
cell 
#> [1] 382

# cell value
grid[cell]
#> [1] -2.651741

reprex package (v2.0.1)

于 2022-03-24 创建