This question was originally posted in StackOverflow.
I am trying to use French spatial data distributed using WFS. I want to use these data in a Shiny app that allows the user to update the data, so I would really prefer to use WFS and not only to embed the data from a shapefile in the app.
I can get the data but there is something off with it in the sense that the coordinates seem to have been reversed (X instead of Y and Y instead of X) resulting in a strange picture of France. I am not very familiar with the details of spatial analyses and data and I am certainly missing something.
I imported the data using the function st_read which I fed with the link to the WFS service. Then I get the upside-down map of French hydrographic network. When I connect QGIS using the same WFS link, however, the map is correct so I guess that sf and QGIS do not exactly handle things in the same way, but investigating why is beyond my knowledge.
When looking at the bounding box we see that the values returned for the x correspond to what we expect for the y and vice versa.
I had the problem with two different data sets distributed using WFS by the same provider (http://www.sandre.eaufrance.fr/). The first one was displaying points and the data contained XY coordinates so I removed the geometry and transformed the data frame using st_as_sf with the names of the X and Y fields in the coords parameter. This is a dirty workaround but it does not work with my second data set that has a MULTILINESTRING geometry. So I would prefer a clean way of importing/transforming the data in the correct format.
Here is a minimal code example to reproduce the unexpected behavior.
library(dplyr) library(sf) library(ggplot2) data <- st_read("http://services.sandre.eaufrance.fr/geo/zon_FXX?Request=GetCapabilities&SERVICE=WFS&VERSION=1.1.0", layer = "SegClassContinuiteEco_Liste2", stringsAsFactors = FALSE) st_bbox(data) ggplot() + geom_sf(data = data)
The code above gives this map: wrong map
The correct map should look like that (I rotated and flipped the original): correct map
R version 3.5.3 (2019-03-11) Platform: x86_64-w64-mingw32/x64 (64-bit) Running under: Windows >= 8 x64 (build 9200) Matrix products: default locale:  LC_COLLATE=French_France.1252 LC_CTYPE=French_France.1252  LC_MONETARY=French_France.1252 LC_NUMERIC=C  LC_TIME=French_France.1252 attached base packages:  stats graphics grDevices utils datasets methods base other attached packages:  sf_0.7-3 forcats_0.4.0 stringr_1.4.0 dplyr_0.8.0.1 purrr_0.3.2  readr_1.3.1 tidyr_0.8.3 tibble_2.1.1 ggplot2_3.1.1 tidyverse_1.2.1 loaded via a namespace (and not attached):  Rcpp_1.0.1 cellranger_1.1.0 pillar_1.3.1 compiler_3.5.3  plyr_1.8.4 class_7.3-15 tools_3.5.3 jsonlite_1.6  lubridate_1.7.4 gtable_0.3.0 nlme_3.1-137 lattice_0.20-38  pkgconfig_2.0.2 rlang_0.3.4 DBI_1.0.0 cli_1.1.0  rstudioapi_0.10 yaml_2.2.0 haven_2.1.0 e1071_1.7-1  withr_2.1.2 xml2_1.2.0 httr_1.4.0 fs_1.2.7  generics_0.0.2 hms_0.4.2.9001 classInt_0.3-1 grid_3.5.3  tidyselect_0.2.5 glue_1.3.1 R6_2.4.0 fansi_0.4.0  readxl_1.3.1 modelr_0.1.4 magrittr_1.5 usethis_1.5.0  units_0.6-2 backports_1.1.4 scales_1.0.0 rvest_0.3.2  assertthat_0.2.1 colorspace_1.4-1 utf8_1.1.4 stringi_1.4.3  lazyeval_0.2.2 munsell_0.5.0 broom_0.5.2 crayon_1.3.4