在 R 中使用 mapSpain 将加那利群岛 sf 数据重新投影到西班牙边界地图的相同投影
Reproject Canary Islands sf data to the same projection of Spain boundary map using mapSpain in R
我有一个包含西班牙一些干旱数据的 sf 对象，我想将这些数据投影到传单地图中。数据来源为：https://www.miteco.gob.es/es/biodiversidad/temas/desertificacion-restauracion/lucha-contra-la-desertificacion/lch_pand_descargas.aspx
我正在使用库 mapSpain 中的函数 esp_get_country() 来获取国家/地区的边界。
问题在于，在 sf 对象中，加那利群岛远离边界地图。如何从 sf 投影加那利群岛多边形，因为它们在边界图中？
我尝试了函数 st_crs() 从边界图中获取 crs 和函数 st_transform() 重新投影 sf 对象，但它没有做任何事情。
我使用的代码是这样的：
# get the boundary map from package 'mapSpain'
spain_coun_sf <- esp_get_country()
# import the drought data from Canary Islands (extraction of the whole drought data)
drought_data_canarias <- readOGR(
dsn = path,
layer = "pand_c",
verbose = FALSE
)
drought_data_canarias <- st_as_sf(drought_data_canarias)
drought_data_canarias <- st_transform(drought_data_canarias , crs = st_crs(spain_coun_sf))
leaflet(spain_coun_sf) %>% addPolygons() %>% addPolygons(data = deser_sf_canarias)
回复你后面的评论：
("Canary map far away from boundary map")
esp_get_country() 图书馆 mapSpain returns 一张扭曲的地图：加那利群岛和西班牙大陆之间的距离被故意缩短了。因此，如果您从具有真实坐标的单独空间特征绘制金丝雀，它们将不会与扭曲的地图重合。
如果您想使用 esp_get_country 返回的功能，您可以执行以下操作：
- 将spain_coun_sf拆分为单独的多边形：
spain_fragments <-
spain_boundary %>%
st_cast('POLYGON') %>% ## split MULTIPOLYGON into POLYGONS
## create polygon IDs for later identification
mutate(polygon_id = as.factor(row_number())
- 使用交互式
ggplotly 绘图识别加那利多边形：
library(plotly)
ggplot() +
geom_sf(data = spain_fragments , aes(label = polygon_id))
ggplotly()
## I got polygon IDs c(17:18, 20:23, 26:28, 33)
## but better check yourself
- 通过选择 spain_fragments:
的相关行来创建仅包含 Canaries 的新 sf
the_canaries <-
spain_fragments %>%
filter(polygon_id %in% c(17:18, 20:23, 26:28, 33)) %>%
st_union
## check:
ggplot() +
geom_sf(data = spain_fragments) +
geom_sf(data = the_canaries, col = 'red')
...现在您有了每个加那利群岛的边界以加入干旱数据。
地图的创建者说西班牙。
您可以使用选项 moveCAN = FALSE 避免加那利群岛的位移。这在包的所有功能中都可用，并且所有文档包都有关于此的特定部分：
https://ropenspain.github.io/mapSpain/reference/esp_get_country.html#displacing-the-canary-islands
Displacing the Canary Islands
While moveCAN is useful for visualization, it would alter the actual geographic position of the Canary Islands. When using the output for spatial analysis or using tiles (e.g. with esp_getTiles() or addProviderEspTiles()) this option should be set to FALSE in order to get the actual coordinates, instead of the modified ones.
在此处查看完整的代表。由于您仅使用 Canarias 图层，我建议您使用 esp_get_ccaa("Canarias", moveCAN = FALSE) 来仅关注该 CCAA：
library(mapSpain)
library(ggplot2)
spain_coun_sf_nomoved <- esp_get_country()
ggplot(spain_coun_sf_nomoved) +
geom_sf()
# Deactivate
spain_coun_sf <- esp_get_country(moveCAN = FALSE)
ggplot(spain_coun_sf) +
geom_sf()
# Your reprex, use CCAA instead to focus in Canarias
canarias <- esp_get_ccaa("Canarias", moveCAN = FALSE)
# Download file
tempfile <- tempfile(fileext = ".zip")
download.file("https://www.miteco.gob.es/es/biodiversidad/servicios/banco-datos-naturaleza/0904712280184ea3_tcm30-177178.zip", tempfile)
unzip(tempfile, junkpaths = TRUE, exdir = tempdir())
library(sf)
#> Linking to GEOS 3.9.1, GDAL 3.2.1, PROJ 7.2.1; sf_use_s2() is TRUE
file_path <- list.files(tempdir(), pattern = "shp$", full.names = TRUE)
drought_data_canarias <- st_read(file_path)
#> Reading layer `pand_c' from data source
#> `C:\Users\XXXX\RtmpCOGdcp\pand_c.shp' using driver `ESRI Shapefile'
#> Simple feature collection with 9642 features and 3 fields
#> Geometry type: POLYGON
#> Dimension: XY
#> Bounding box: xmin: 188300.1 ymin: 3060674 xmax: 653718.9 ymax: 3235418
#> Projected CRS: WGS 84 / UTM zone 28N
# Project everything to 4326 for leaflet
canarias_4326 <- st_transform(canarias, 4326)
drought_data_canarias_4326 <- st_transform(drought_data_canarias, 4326)
ggplot(canarias_4326) +
geom_sf() +
geom_sf(data=drought_data_canarias_4326, aes(fill=DESER_CLA), color=NA) +
scale_fill_viridis_c()
library(leaflet)
leaflet(canarias_4326) %>%
addTiles() %>%
addPolygons() %>%
addPolygons(data = drought_data_canarias_4326)
由 reprex package (v2.0.1)
于 2022-05-10 创建
我有一个包含西班牙一些干旱数据的 sf 对象，我想将这些数据投影到传单地图中。数据来源为：https://www.miteco.gob.es/es/biodiversidad/temas/desertificacion-restauracion/lucha-contra-la-desertificacion/lch_pand_descargas.aspx
我正在使用库 mapSpain 中的函数 esp_get_country() 来获取国家/地区的边界。
问题在于，在 sf 对象中，加那利群岛远离边界地图。如何从 sf 投影加那利群岛多边形，因为它们在边界图中？
我尝试了函数 st_crs() 从边界图中获取 crs 和函数 st_transform() 重新投影 sf 对象，但它没有做任何事情。
我使用的代码是这样的：
# get the boundary map from package 'mapSpain'
spain_coun_sf <- esp_get_country()
# import the drought data from Canary Islands (extraction of the whole drought data)
drought_data_canarias <- readOGR(
dsn = path,
layer = "pand_c",
verbose = FALSE
)
drought_data_canarias <- st_as_sf(drought_data_canarias)
drought_data_canarias <- st_transform(drought_data_canarias , crs = st_crs(spain_coun_sf))
leaflet(spain_coun_sf) %>% addPolygons() %>% addPolygons(data = deser_sf_canarias)
回复你后面的评论：
("Canary map far away from boundary map")
esp_get_country() 图书馆 mapSpain returns 一张扭曲的地图：加那利群岛和西班牙大陆之间的距离被故意缩短了。因此，如果您从具有真实坐标的单独空间特征绘制金丝雀，它们将不会与扭曲的地图重合。
如果您想使用 esp_get_country 返回的功能，您可以执行以下操作：
- 将spain_coun_sf拆分为单独的多边形：
spain_fragments <-
spain_boundary %>%
st_cast('POLYGON') %>% ## split MULTIPOLYGON into POLYGONS
## create polygon IDs for later identification
mutate(polygon_id = as.factor(row_number())
- 使用交互式
ggplotly绘图识别加那利多边形：
library(plotly)
ggplot() +
geom_sf(data = spain_fragments , aes(label = polygon_id))
ggplotly()
## I got polygon IDs c(17:18, 20:23, 26:28, 33)
## but better check yourself
- 通过选择 spain_fragments: 的相关行来创建仅包含 Canaries 的新
sf
the_canaries <-
spain_fragments %>%
filter(polygon_id %in% c(17:18, 20:23, 26:28, 33)) %>%
st_union
## check:
ggplot() +
geom_sf(data = spain_fragments) +
geom_sf(data = the_canaries, col = 'red')
...现在您有了每个加那利群岛的边界以加入干旱数据。
地图的创建者说西班牙。
您可以使用选项 moveCAN = FALSE 避免加那利群岛的位移。这在包的所有功能中都可用，并且所有文档包都有关于此的特定部分：
https://ropenspain.github.io/mapSpain/reference/esp_get_country.html#displacing-the-canary-islands
Displacing the Canary Islands
While
moveCANis useful for visualization, it would alter the actual geographic position of the Canary Islands. When using the output for spatial analysis or using tiles (e.g. withesp_getTiles()oraddProviderEspTiles()) this option should be set toFALSEin order to get the actual coordinates, instead of the modified ones.
在此处查看完整的代表。由于您仅使用 Canarias 图层，我建议您使用 esp_get_ccaa("Canarias", moveCAN = FALSE) 来仅关注该 CCAA：
library(mapSpain)
library(ggplot2)
spain_coun_sf_nomoved <- esp_get_country()
ggplot(spain_coun_sf_nomoved) +
geom_sf()
# Deactivate
spain_coun_sf <- esp_get_country(moveCAN = FALSE)
ggplot(spain_coun_sf) +
geom_sf()
# Your reprex, use CCAA instead to focus in Canarias
canarias <- esp_get_ccaa("Canarias", moveCAN = FALSE)
# Download file
tempfile <- tempfile(fileext = ".zip")
download.file("https://www.miteco.gob.es/es/biodiversidad/servicios/banco-datos-naturaleza/0904712280184ea3_tcm30-177178.zip", tempfile)
unzip(tempfile, junkpaths = TRUE, exdir = tempdir())
library(sf)
#> Linking to GEOS 3.9.1, GDAL 3.2.1, PROJ 7.2.1; sf_use_s2() is TRUE
file_path <- list.files(tempdir(), pattern = "shp$", full.names = TRUE)
drought_data_canarias <- st_read(file_path)
#> Reading layer `pand_c' from data source
#> `C:\Users\XXXX\RtmpCOGdcp\pand_c.shp' using driver `ESRI Shapefile'
#> Simple feature collection with 9642 features and 3 fields
#> Geometry type: POLYGON
#> Dimension: XY
#> Bounding box: xmin: 188300.1 ymin: 3060674 xmax: 653718.9 ymax: 3235418
#> Projected CRS: WGS 84 / UTM zone 28N
# Project everything to 4326 for leaflet
canarias_4326 <- st_transform(canarias, 4326)
drought_data_canarias_4326 <- st_transform(drought_data_canarias, 4326)
ggplot(canarias_4326) +
geom_sf() +
geom_sf(data=drought_data_canarias_4326, aes(fill=DESER_CLA), color=NA) +
scale_fill_viridis_c()
library(leaflet)
leaflet(canarias_4326) %>%
addTiles() %>%
addPolygons() %>%
addPolygons(data = drought_data_canarias_4326)
由 reprex package (v2.0.1)
于 2022-05-10 创建