创建具有四个角坐标的矩形多边形 shapefile

Create rectangular polygon shapefiles with four corner coordinates

我有很多观察,每个观察都有四个角坐标来确定一个矩形区域,例如x_min、x_max、y_min和y_max。我想使用四个角坐标创建矩形多边形,然后将它们导出为 shapefile 以便在 ArcGIS 中使用。

以下是一些观察结果。

ID     xmax        xmin         ymax       ymin
1   -6987035.01 -8831307.63  6160489.27  4818866.55
2   -7457581.36 -7918649.51  5705536.08  5370130.4
3   -7527291.93 -7988360.08  5623289.83  5287884.15
4   -7632927.9  -7863461.98  5533399.89  5365697.05

根据坐标,矩形的左下点应该是(xmin,ymin),矩形的左上点应该是(xmin,ymax),矩形的右下点应该是应为 (xmax,ymin),矩形的右上角应为 (xmax, ymax)。

如何根据这四个坐标制作多边形文件并将其导出为可在 ArcGIS 中使用的 shapefile?

我不确定这是否适合你

plot(range(unlist(df[startsWith(names(df), "x")])),
  range(unlist(df[startsWith(names(df), "y")])),
  type = "n", xlab = "", ylab = ""
)
with(df, rect(xmin, ymin, xmax, ymax))

给予

试试下面的方法。

我将您的数据重新格式化为坐标矩阵,假设每行代表一个多边形的边界。

bounds <- matrix(
  c(
    -6987035.01, -8831307.63,  6160489.27,  4818866.55,
    -7457581.36, -7918649.51,  5705536.08,  5370130.4,
    -7527291.93, -7988360.08,  5623289.83,  5287884.15,
    -7632927.9,  -7863461.98,  5533399.89,  5365697.05
  ), nrow = 4, byrow = TRUE, dimnames = list(1:4, c('xmax','xmin','ymax','ymin'))
)

library(rgdal)
library(raster)

polygons_list = apply(bounds, 1, function(x){
  ## Loop throught rows to create extent objects
  out = extent(x[c(2, 1, 4, 3)]) ## x[c(2, 1, 4, 3)] to reoder the row to match the order required by raster::extent function. The order should xmin, xmax, ymin, ymax
  ## convert the extent objects to spatial polygons
  out = as(out, 'SpatialPolygons')
  ## convert the spatial polygons to spatial polygons dataframe
  out = as(out, 'SpatialPolygonsDataFrame')
  ## set CRS 
  crs(out) <- '+proj=merc +a=6378137 +b=6378137 +lat_ts=0.0 +lon_0=0.0 +x_0=0.0 +y_0=0 +k=1.0 +units=m +nadgrids=@null +no_defs'
  out
})

## conbining the list of spatial polygons dataframes into one spatial polygons dataframe
polygons_list = do.call(rbind, polygons_list)
plot(polygons_list)

## export as shapefiles for future use in ArcGIS
writeOGR(polygons_list, getwd(), "polygons_list", driver="ESRI Shapefile", overwrite_layer=FALSE)