Skip to contents
library(happign)
library(sf)
library(tmap); tmap_mode("view") # Set map to interactive
library(dplyr)
library(ggplot2);sf_use_s2(FALSE) # Avoid problem with spherical geometry
library(purrr)
library(stars)
library(terra)

First choose a zone of interest

For the example we will work with the Camors forest. First of all we need the commune border that can be obtained from the insee code. Fortunately, happign provides a table containing all insee codes (data("cog_2023")). Then, get_apicarto_commune is used to download shape from Apicarto commune.

data("cog_2023")

insee_code <- cog_2023[grepl("^Camors", cog_2023$LIBELLE),1]

borders <- get_apicarto_cadastre(insee_code, type = "commune")

tm_shape(borders)+
   tm_borders()+
   tm_text("nom_com")

An other way of getting borders without apicarto is to use ECQL language to directly query IGN WFS geoservers.

borders2 <- get_wfs(layer = "LIMITES_ADMINISTRATIVES_EXPRESS.LATEST:commune",
                    spatial_filter = "intersects",
                    ecql_filter = "nom_m LIKE 'CAMORS'")
#> Features downloaded : 1

Download cadastral parcel

Cadastral parcels are essential for any forest manager Here how to download it with get_wfs.

layers <- get_layers_metadata("wfs", "parcellaire")
parcellaire_layer <- layers[15,1] # "CADASTRALPARCELS.PARCELLAIRE_EXPRESS:parcelle"

parcellaire <- get_wfs(x = borders,
                       spatial_filter = "intersects",
                       layer = parcellaire_layer)
#> Features downloaded : 1000...2000...3000...3490

tm_shape(borders)+
   tm_borders(col = "red", lwd = 2)+
tm_shape(parcellaire)+
   tm_polygons(alpha = 0)