使用 R 中的 Raster 包查找单元格中存在的附近行和列的值

Question

我想查找单元格中存在的行和列的值。我已经编写了一个代码来获取值，如下所示：

x<-raster()
    values(x)<-1:ncell(x)
class       : RasterLayer 
dimensions  : 180, 360, 64800  (nrow, ncol, ncell)
resolution  : 1, 1  (x, y)
extent      : -180, 180, -90, 90  (xmin, xmax, ymin, ymax)
coord. ref. : +proj=longlat +datum=WGS84 +ellps=WGS84 +towgs84=0,0,0 
data source : in memory
names       : layer 
values      : 1, 64800  (min, max)

下面给出的另一个文件包含经纬度信息。因此，我们必须在栅格对象中找到这些点，然后以最后给出的格式存储值。

long<-c(106.61291,-81.97224,-84.4277,-97.66631,-72.68604,-73.1247,-119.05662,
   -86.75413,-108.5377,-86.67737,-116.22231,-71.00956,-91.15146,-73.1516)
lat<-c(35.04333,33.37378,33.64073,30.19743,41.93887,41.1642,35.43291,33.56243,
       45.8037,36.12632,43.56582,42.36561,30.53236,44.4707)

xy<-data.frame(long,lat)
coordinates(xy)=c('long','lat')
proj4string(xy)=proj4string(x)
points<-cellFromXY(x,xy)

so , points 包含 xy 在栅格对象 x 中存在的所有像元编号。

As given below:
[1] 19727 20259 20256 21323 17388 17387 19501 20254 15912 19174 16624 17029 21329 16307

我想找到以点为单位的每个单元格的行和列的值，如下所示： 例如从“点”中获取单元格编号 19727

rowColFromCell(x,19727)
  row col
  [1,]  55 287

  getValues(x,row=55)[287]
  value:16109

我想要 9*9 网格的值，其中给定的行和列即 55,287 位于中心。

我想对点中给出的所有单元格做同样的事情。

最终输出应如下所示：

 long        lat    cell_1  cell_2  cell_3 cell _4 ......... cell_9
106.61291 35.04333     16109  ...     ....   ....              ....

提前致谢！

Answer 1

我们可以定义一个函数，根据栅格层的分辨率来获取每个点的 9*9 网格的坐标。之后，我们可以使用 raster 包中的 extract 函数来提取值。在以下示例中，xy2 是最终输出。列5 是来自中心单元格的值。所有其他列显示来自相邻 9*9 网格的值。

library(raster)

# Create the example raster
r <- raster()
values(r) <- 1:ncell(r)

# Create a data frame with long-lat coordinates
long<-c(106.61291,-81.97224,-84.4277,-97.66631,-72.68604,-73.1247,-119.05662,
        -86.75413,-108.5377,-86.67737,-116.22231,-71.00956,-91.15146,-73.1516)
lat<-c(35.04333,33.37378,33.64073,30.19743,41.93887,41.1642,35.43291,33.56243,
       45.8037,36.12632,43.56582,42.36561,30.53236,44.4707)

xy <- data.frame(long, lat)

# Design a function to calculate the coordinates of the surrounding eight cells
coord_nine <- function(long, lat, res){
  cell1 <- c(long - res[1], lat + res[2])
  cell2 <- c(long, lat + res[2])
  cell3 <- c(long + res[1], lat + res[2])
  cell4 <- c(long - res[1], lat)
  cell5 <- c(long, lat)
  cell6 <- c(long + res[1], lat)
  cell7 <- c(long - res[1], lat - res[2])
  cell8 <- c(long, lat - res[2])
  cell9 <- c(long + res[1], lat - res[2])

  coord_list <- list(cell1, cell2, cell3, cell4, cell5, 
                     cell6, cell7, cell8, cell9)

  coord <- do.call(rbind, coord_list)
  return(coord)
}

# Use lapply to loop through all points to get the 9-cell coordinates for each of them
coord_list <- lapply(1:nrow(xy), function(i){
  coord_nine(xy[i, 1], xy[i, 2], res = res(r))
})

# Extract the values based on the coordinates in coord_list
ex_values <- sapply(coord_list, function(coords){
  extract(r, coords)
})

# Combine the results to the original xy data frame
xy2 <- cbind(xy, t(ex_values))

# View the results
xy2
         long      lat     1     2     3     4     5     6     7     8     9
1   106.61291 35.04333 19366 19367 19368 19726 19727 19728 20086 20087 20088
2   -81.97224 33.37378 19898 19899 19900 20258 20259 20260 20618 20619 20620
3   -84.42770 33.64073 19895 19896 19897 20255 20256 20257 20615 20616 20617
4   -97.66631 30.19743 20962 20963 20964 21322 21323 21324 21682 21683 21684
5   -72.68604 41.93887 17027 17028 17029 17387 17388 17389 17747 17748 17749
6   -73.12470 41.16420 17026 17027 17028 17386 17387 17388 17746 17747 17748
7  -119.05662 35.43291 19140 19141 19142 19500 19501 19502 19860 19861 19862
8   -86.75413 33.56243 19893 19894 19895 20253 20254 20255 20613 20614 20615
9  -108.53770 45.80370 15551 15552 15553 15911 15912 15913 16271 16272 16273
10  -86.67737 36.12632 18813 18814 18815 19173 19174 19175 19533 19534 19535
11 -116.22231 43.56582 16263 16264 16265 16623 16624 16625 16983 16984 16985
12  -71.00956 42.36561 16668 16669 16670 17028 17029 17030 17388 17389 17390
13  -91.15146 30.53236 20968 20969 20970 21328 21329 21330 21688 21689 21690
14  -73.15160 44.47070 15946 15947 15948 16306 16307 16308 16666 16667 16668