根据另一个 sf 对象计算 sf 对象的面积答案

【问题标题】：Calculate area of sf object based on another sf object根据另一个 sf 对象计算 sf 对象的面积
【发布时间】：2020-11-13 09:05:59
【问题描述】：

我有一个研究区（左），以及解剖它的道路（中）。我想计算两种模式相交所产生的区域（右）——这将产生 5 个子区域，这些子区域的总和为研究区域对象的总面积。

两个对象都是sf 并且具有相同的投影（见下文）。

如何使用sf 库计算子区域？

> str(myarea)
Classes ‘sf’ and 'data.frame':  1 obs. of  4 variables:
 $ Id      : num 0
 $ Habitat : num 1
 $ Area    : num 202
 $ geometry:sfc_POLYGON of length 1; first list element: List of 1
  ..$ : num [1:135, 1:2] 858895 859084 859358 859865 860105 ...
  ..- attr(*, "class")= chr  "XY" "POLYGON" "sfg"
 - attr(*, "sf_column")= chr "geometry"
 - attr(*, "agr")= Factor w/ 3 levels "constant","aggregate",..: NA NA NA
  ..- attr(*, "names")= chr  "Id" "Habitat" "Area"


> str(roads)
Classes ‘sf’ and 'data.frame':  12 obs. of  23 variables:
 $ LINKNO    : Factor w/ 13673 levels "A1-1-1-1","A1-1-1-2",..: 7483 7327 7326 7325 2433 2436 174 2438 2439 2434 ...
 $ ROADNO    : Factor w/ 2981 levels "A1","A104","A104A",..: 1444 1415 1415 1415 246 246 2 247 247 246 ...
 $ START_KM  : num  0 0 0 0 0 0 0 0 0 0 ...
 $ END_KM    : num  14 8 4 1 5 7 16 10 7 6 ...
 $ LENGTHKM  : num  14.07 7.95 4.25 1.32 4.92 ...
 $ STRATDESC : Factor w/ 2908 levels "1005","1699A",..: 61 1778 1112 1112 628 1115 92 1112 1112 628 ...
 $ ENDDESC   : Factor w/ 2990 levels "1","A.P.LINE",..: 58 1667 1668 1668 1025 72 655 56 56 1025 ...
 $ ROADCLASS : Factor w/ 12 levels "A","B","C","D",..: 5 5 5 5 4 4 1 4 4 4 ...
 $ CLASS     : Factor w/ 3 levels "Primary (Trunk)",..: 3 3 3 3 2 2 1 2 2 2 ...
 $ REGION    : Factor w/ 9 levels "0","Central",..: 8 8 8 8 8 8 8 8 8 8 ...
 $ WIDTH     : num  6 5.5 5.5 5.5 5.5 5.5 7 5.5 5.5 6 ...
 $ LANES     : int  2 2 2 2 2 2 2 2 2 2 ...
 $ TYPE      : Factor w/ 10 levels "230","Concrete (Jt-Plain)",..: 7 4 4 4 4 6 7 6 4 10 ...
 $ SURFTYPE  : Factor w/ 6 levels "Gravel","Paved",..: 2 6 6 6 6 6 2 6 6 2 ...
 $ CONDITION : Factor w/ 6 levels "Excellent","Fair",..: 4 6 2 4 2 3 4 4 4 3 ...
 $ AADT      : num  0 0 69 69 50 ...
 $ iso3      : Factor w/ 1 level "KEN": 1 1 1 1 1 1 1 1 1 1 ...
 $ AICD_REG  : int  NA NA NA NA NA NA 1 NA NA NA ...
 $ TAH       : int  NA NA NA NA NA NA 1 NA NA NA ...
 $ WB_Project: int  NA NA NA NA NA NA NA NA NA NA ...
 $ Kamp_Momba: int  NA NA NA NA NA NA 1 NA NA NA ...
 $ Nairo_Addi: int  NA NA NA NA NA NA NA NA NA NA ...
 $ geometry  :sfc_GEOMETRY of length 12; first list element:  'XY' num [1:581, 1:2] 858210 858211 858215 858219 858223 ...
 - attr(*, "sf_column")= chr "geometry"
 - attr(*, "agr")= Factor w/ 3 levels "constant","aggregate",..: NA NA NA NA NA NA NA NA NA NA ...
  ..- attr(*, "names")= chr  "LINKNO" "ROADNO" "START_KM" "END_KM" ...

【问题讨论】：

标签： r sf

【解决方案1】：

您可以使用 lwgeom 库中的 st_split()，然后使用 st_area 计算面积。我只是编了一些数据来演示。请注意，这会将您的初始多边形裁剪为多个多边形。

library(sf)
library(lwgeom)

# Some sample data
poly <- st_bbox(c(xmin = 0, xmax = 15, ymin = 0, ymax=15)) %>% st_as_sfc()
lines <- data.frame(x = seq(-1, 18), y = seq(1,20)) %>%      st_as_sf(coords=c("x","y")) %>% st_union() %>% st_cast('LINESTRING')

初始

# split and add area to ssf
spltply <- lwgeom::st_split(poly, lines) %>% st_collection_extract(c("POLYGON")) %>% st_sf()
spltply['area'] <- st_area(spltply)

拆分后

【讨论】：