为什么不同的栅格索引函数会给我不同的答案?
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 创建
我正在尝试在空间 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 创建