从WMS层提取数据

问题描述 投票:0回答:1

我试图从WMS层提取数据。

作为一个例子,我想分析一下我的区域是否触及了Natura2000区域以及Natura2000区域的具体细节。

Natura2000区域的WMS层位于:https://geodata.nationaalgeoregister.nl/natura2000/ows?service=WMS&request=GetLegendGraphic&format=image/png&width=20&height=20&layer=natura2000

假设我的区域是一个围绕某个x坐标和y坐标的半径为7500米的圆;

我试图用传单包完成这个,但它似乎更像是一个显示信息的工具,而不是分析信息。

x.WGS=6.662226
y.WGS=52.53206

leaflet() %>% setView(x.WGS, y.WGS, zoom = 11) %>%   
  addTiles() %>% 
  addMarkers(lng = x.WGS, lat = y.WGS)%>%
  addWMSTiles(
    "https://geodata.nationaalgeoregister.nl/natura2000/ows?service=WMS&request=GetLegendGraphic&format=image/png&width=20&height=20&layer=natura2000",
  layers = "natura2000",
  options = WMSTileOptions(format = "image/png", transparent = TRUE),
  attribution = "") %>%
  addCircles(lng = x.WGS, lat = y.WGS, weight = 1,
  radius = 7500)

我希望它能归还两件事。我的地区是否触及任何Natura2000区域?换句话说,它们的名称是什么?如果我在Qgis中加载WMS层,那么Natura2000区域的名称应该在naam_n2k下。

在我的例子中,答案应该是我的区域有两个Natura2000区域,这些Natrua2000区域的名称是Vecht-en Beneden-Reggegebied和Engbertsdijksvenen。

r leaflet wms r-leaflet
1个回答
0
投票

传单包是一个用于访问传单功能的R包,它是用于显示交互式地图的最流行的开源JavaScript库之一。更多信息:https://rstudio.github.io/leaflet/

正如@IvanSanchez所说,Web地图服务(WMS)是一个可以在您自己的应用程序中加载的地图。将此视为地理参考图像/图片,没有关于实际多边形,形状等的信息。 如果要对多边形进行某些分析或检索信息,则可能需要Web要素服务(WFS)。

首先,这是正确加载WMS的代码。

# corrected WMS version
library(leaflet)
library(leaflet.extras) # for WMS legend

x.WGS=6.662226
y.WGS=52.53206

leaflet() %>% setView(x.WGS, y.WGS, zoom = 11) %>%
  addTiles() %>%
  addMarkers(lng = x.WGS, lat = y.WGS)%>%
  addWMSTiles(
    baseUrl = "https://geodata.nationaalgeoregister.nl/natura2000/ows?service=WMS",
    layers = "natura2000",
    options = WMSTileOptions(format = "image/png", transparent = TRUE),
    attribution = "") %>%
  addCircles(lng = x.WGS, lat = y.WGS, weight = 1,
             radius = 7500) %>%
  addWMSLegend(uri ="https://geodata.nationaalgeoregister.nl/natura2000/ows?service=WMS&request=GetLegendGraphic&format=image%2Fpng&width=20&height=20&layer=natura2000")

其次,回答你的问题:我们可以将带有几何过滤器的CQL查询添加到WFS请求中。 请注意,这在您的特定情况下是可行的,但这不是所有WFS服务的标准。

# Since you want to do spatial analysis on features, we load the sf package
# more info: https://r-spatial.github.io/sf/
library(sf)

# we convert your point coordinates to the native CRS of the WFS layer:
  # make point and assign WGS84 CRS
  x <- st_point(c(x.WGS,y.WGS))
  x <- st_sfc(x) %>% st_set_crs(NA) %>% st_set_crs(4326)

  # transform to EPSG::28992 // Amersfoort RD is default CRS for the wfs layer
  x_RD <- st_transform(x,28992)

# The WFS allows spatial filtering (this is not always the case)
# so you can use a cql filter with a distance around your point:
# cql_filter=dwithin(natura2000:geom,point(241514 505696.5),7500,meters)
# combining and encoding the cql filter for a valid WFS url:
WFS_url <- paste0("http://geodata.nationaalgeoregister.nl/natura2000/wfs?",
                  "&service=wfs&version=2.0.0&request=GetFeature&",
                  "typeName=natura2000:natura2000&cql_filter=",
                  URLencode(
                    paste0("dwithin(natura2000:geom,",
                            st_as_text(x_RD),
                           ",7500, meters)"),
                    reserved = FALSE),
                  "&outputFormat=application/json"
                  )

# get WFS feature
natura2000 <- st_read(WFS_url)

# Transform to WGS84
natura2000_wgs84 <- st_transform(natura2000,4326)

# load data in leaflet
leaflet() %>%
  addTiles() %>%
  addMarkers(lng = x.WGS, lat = y.WGS) %>%
  addCircles(lng = x.WGS, lat = y.WGS, weight = 1,
             radius = 7500) %>%
  addPolygons(data = natura2000_wgs84,
              label = ~naam_n2k,
              popup = ~naam_n2k)
© www.soinside.com 2019 - 2024. All rights reserved.