R - 将经度和纬度坐标放入 2x2 网格中

问题描述

我有一个包含经度和纬度的数据框，我想创建一个 0.5x0.5 度的网格，显示哪个纬度、经度落在其中。到目前为止，我已经尝试了几种解决方案，包括在 stackoverflow 上找到的一些解决方案，它们使用 cut 和 expand.grid 以及使用包“sp”的代码，但没有一个对我有用（也许我根本无法实现它们）。

关于如何将数据分组为 0.5x0.5 度网格的任何建议？

纬度	经度
31.602	-39.848
31.675	-39.467
31.747	-39.083
32.152	-36.795
32.218	-36.408
32.285	-36.022
32.348	-35.635
32.412	-35.247
32.475	-34.858
32.535	-34.47
32.595	-34.082
32.677	-33.707
32.763	-33.323

感谢大家付出的时间和努力。

编辑：我最大的努力就是这个片段

library(tidyverse)

pos <- dassem %>% 
  dplyr::select(Latitude,Longitude)

gridx <- seq(from = min(dassem$Longitude),to = max(dassem$Longitude),by = 2)
gridy <- seq(from = min(dassem$Latitude),to = max(dassem$Latitude),by = 2)

xcell <- unlist(lapply(pos$Longitude,function(x) min(which(gridx>x))))
ycell <- unlist(lapply(pos$Latitude,function(y) min(which(gridy>y))))

pos$cell <- (length(gridx) - 1) * ycell + xcell

pos

# A tibble: 45,647 x 3
   Latitude Longitude  cell
      <dbl>     <dbl> <dbl>
 1     51.7     -54.9   638
 2     51.9     -54.5   638
 3     52.1     -54.1   638
 4     52.3     -53.7   639
 5     52.5     -53.2   639
 6     52.7     -52.8   639
 7     52.9     -52.4   639
 8     53.2     -52.0   639

如您所见，它不会返回 2x2 度的网格（我将其设置为 2x2，而不是 0.5x0.5）。

解决方法

library(tidyverse)
library(sf)

df_sf <- df %>%
  st_as_sf(coords = c("lon","lat"),crs = 4326)

grid <- df_sf %>% 
  st_bbox() %>% 
  st_as_sfc() %>% 
  st_make_grid(cellsize = 0.5)

df %>%
  mutate(polygon_id = st_intersects(df_sf,grid) %>% map_int(1))

目前还不清楚您的预期输出是什么，但我会提供几个选项：

dassem <- structure(list(Latitude = c(31.602,31.675,31.747,32.152,32.218,32.285,32.348,32.412,32.475,32.535,32.595,32.677,32.763),Longitude = c(-39.848,-39.467,-39.083,-36.795,-36.408,-36.022,-35.635,-35.247,-34.858,-34.47,-34.082,-33.707,-33.323)),class = "data.frame",row.names = c(NA,-13L))

我将从 1 度网格（任意）开始：

degx <- degy <- 1 # grid sizes
gridx <- seq(min(dassem$Longitude),max(dassem$Longitude) + degx,by = degx)
gridy <- seq(min(dassem$Latitude),max(dassem$Latitude) + degy,by = degy)

dassem %>%
  mutate(
    x = findInterval(Longitude,gridx),y = findInterval(Latitude,gridy)
  )
#    Latitude Longitude x y
# 1    31.602   -39.848 1 1
# 2    31.675   -39.467 1 1
# 3    31.747   -39.083 1 1
# 4    32.152   -36.795 4 1
# 5    32.218   -36.408 4 1
# 6    32.285   -36.022 4 1
# 7    32.348   -35.635 5 1
# 8    32.412   -35.247 5 1
# 9    32.475   -34.858 5 1
# 10   32.535   -34.470 6 1
# 11   32.595   -34.082 6 1
# 12   32.677   -33.707 7 2
# 13   32.763   -33.323 7 2

这提供了每个点的网格索引（如果有兴趣的话）。

如果你不需要单独的索引，你可以组合它们（就像我对那个图像所做的那样）：

dassem %>%
  mutate(
    cell = paste(findInterval(Longitude,findInterval(Latitude,gridy),sep = ",")
  )
#    Latitude Longitude cell
# 1    31.602   -39.848  1,1
# 2    31.675   -39.467  1,1
# 3    31.747   -39.083  1,1
# 4    32.152   -36.795  4,1
# 5    32.218   -36.408  4,1
# 6    32.285   -36.022  4,1
# 7    32.348   -35.635  5,1
# 8    32.412   -35.247  5,1
# 9    32.475   -34.858  5,1
# 10   32.535   -34.470  6,1
# 11   32.595   -34.082  6,1
# 12   32.677   -33.707  7,2
# 13   32.763   -33.323  7,2

geospatial r r