【问题标题】:Export a polygon shapefile as a list of XY vertex pairs将多边形 shapefile 导出为 XY 顶点对列表
【发布时间】:2019-10-17 18:20:45
【问题描述】:

对于一个连续的多边形 shapefile,我需要一个格式如下的多边形顶点坐标文件作为遗留程序的输入:

polygon1name, AUTO
  -120.408750    34.591250
  -120.398313    34.591250
  -120.396250    34.593313
  -120.396250    34.593354
END
polygon2name, AUTO
  -120.423354    34.641250
  -120.423313    34.641250
  -120.421250    35.643313
  -120.421250    35.647521
END

从一个示例文件看来,旧程序希望这些对按逆时针绘制顺序排列。下面以北卡罗来纳县为例。我希望在如何导出 XY 对和包含 , AUTOEND 件方面获得帮助。

library(tidyverse) #for the %>% pipes and transmute()
library(sf) #for st_read()
library(rmapshaper) #for ms_simplify()

nc <- st_read(system.file("shape/nc.shp", package="sf")) %>%
      transmute(NAME, geometry) %>% #keeps just the county column for simplicity
      ms_simplify(keep = 0.01) #reduces the number of vertices for simplicity
plot(nc)

有什么想法吗?

谢谢。

【问题讨论】:

  • 您可以通过将数据嵌套在NAME 上开始,然后添加类似mutate(coords = map(data, st_coordinates)) 的内容。至于逆时针排序,我不知道,这就是为什么我没有足够的内容来充实一个完整的答案
  • 谢谢@camille

标签: r sf


【解决方案1】:

不确定坐标的逆时针顺序,但这回答了您问题的输出格式部分。

#extract coordinates from sf
coord    <- st_coordinates(nc) %>% 
  as.data.frame() %>%
  group_by( L3 ) %>% 
  mutate(L4 = row_number() )

#extract data from sf
polygons <- st_drop_geometry(nc) %>% 
  mutate( NAME = as.character( NAME ) ) %>%
  rownames_to_column( var = "id" ) %>% 
  mutate( id = as.numeric(id) ) %>%
  #join coordinates
  left_join( coord, by = c("id" = "L3") )

#split polygons-dataframe to list
l <- split( polygons, f = polygons$id )

#extract text needed from each polygon
result <- lapply( l, function(x) {
  paste0 ( paste0( unique( x$NAME ), ", AUTO\n" ),
           paste0( "  ", x$X, "    ", x$Y, collapse = "\n" ),
           "\nEND" )
})

#unlist and write lines
writeLines( unlist(result) )

# Ashe, AUTO
# -81.4727554    36.2343559
# -81.7410736    36.3917847
# -81.6699982    36.5896492
# -81.3452988    36.5728645
# -81.2398911    36.3653641
# -81.4727554    36.2343559
# END
# Alleghany, AUTO
# -81.2398911    36.3653641
# -81.3452988    36.5728645
# -80.9034424    36.5652122
# -80.9563904    36.4037971
# -81.2398911    36.3653641
# END
# Surry, AUTO
# -80.4563446    36.2425575
# -80.874382    36.2338829

更新

对于逆时针部分:查看sf::st_read()check_ring_dir 参数。

check_ring_dir
合乎逻辑的;如果为 TRUE,则检查多边形环方向 必要时进行修正(从上方看:外环 逆时针,孔顺时针)

【讨论】:

  • 更新了关于 st_read() 中顺时针/逆时针参数的信息
  • 非常感谢@Wimpel! st_coordinates() 的新手,我很头疼,为什么 st_coordinates(my_sf_file) 生成 L1 和 L2 而不是 L3 列? my_sf_file &lt;- readOGR(".", "myshapefile) %&gt;% st_as_sf() %&gt;% transmute(Name, geometry) 使我看起来与nc 相同的格式..
  • 更新 - 我想我明白为什么了,my_sf_file 被标记为 POLYGON 而不是 MULTIPOLYGON..
  • 更新——确实——只需要使用st_cast("MULTIPOLYGON")
猜你喜欢
  • 2012-12-05
  • 2013-01-24
  • 1970-01-01
  • 1970-01-01
  • 2013-01-13
  • 1970-01-01
  • 2021-08-15
  • 1970-01-01
  • 2016-05-28
相关资源
最近更新 更多