# Extraction du réseau bocager avec HedgeTools
Note méthodologique pour extraire le bocage des rasters arborés issus de notre travail acharné de digitalisation de vignettes :)
---
## Vue d'ensemble
Cette chaîne de traitement que je vous propose et que je vais appliquer aussi à toute la Normandie pour mon travail repose sur une partie de la méthodologie développé par Marques *et al.* (2024) "HedgeTools : un outil d’analyse spatiale dédié à l’évaluation de la multifonctionnalité des haies."
Les traitements sont plus rapides pour les buffers et j'utilise le serveur de calcul pour la Normandie.
<br>
La différence réside dans le fait que l'on saute l'étape de différenciation de la couverture arborée qu'ils ont développé à partir d'orthoimages infrarouge et du modèle numérique de surface.
On part d'un **raster binaire** (végétation arborée = 1, reste = 0) d'une résolution de 5x5m pour arriver à un **réseau linéaire de haies**. L'ensemble des traitements sont reproductible sur n'importe quel raster arboré en entrée.
A noter que le plugin Hedgetools propose aussi d'autres mesures du paysage bocager à partir des **polygones** de haies.
```
RASTER ARBORÉ → POLYGONES → PRÉTRAITEMENT → CLASSIFICATION → LINÉAIRE → EXTRACTION → FILTRE
.tif .gpkg .gpkg .gpkg .gpkg .gpkg .gpkg
```
----
Les rasters arborés fusionnés de la Normandie sont disponibles sur le serveur de partage du labo IDEES (K:\Labo\111 - ECOTONES).
Si vous n'y avez pas accès, il faut demander à Julien Frérot qui vous fournir le .bat à éxecuter pour y accéder.
### Les 5 étapes
| Étape | Ce qu'on fait | Résultat |
| ------------------ | -------------------------------------------- | ----------------- |
| **Polygonisation** | Convertir les pixels en polygones | Polygones bruts |
| **Prétraitement** | Nettoyer les géométries | Polygones propres |
| **Classification** | Séparer forêt / haie / bosquet / arbre isolé | 4 couches |
| **Linéarisation** | Extraire l'axe médian des haies | Réseau linéaire |
| **Simplification/Extraction** | Différencier les haies des vergers | Réseau linéaire |
---
## Prérequis
- **QGIS** (3.28+) avec le plugin **HedgeTools** installé -> Télécharger le zip du plugin (https://plugins.qgis.org/plugins/hedge_tools/) -> extension qgis -> installer depuis un zip
- Un raster binaire fusionné de toute la Normandie de végétation arborée fourni par T.Gaubert :)
---
## Organisation des fichiers
```
mon_projet/
├── data/
│ ├── raster/
│ │ └── arbore.tif ← Raster d'entrée
│ └── vecteur/
│ ├── buffers.gpkg ← Zones d'étude
│ ├── buffers_polygonises/ ← Étape 1
│ └── buffers_preprocessed/ ← Étape 2
└── results/
└── buffers_ecotone_class/
├── forest/ ← Étape 3
├── hedge/ ← Étape 3
├── grove/ ← Étape 3
├── lonetree/ ← Étape 3
├── lineaire_brut/ ← Étape 4
└── lineaire_filtre/ ← Étape 5
```
---
## 1. Polygonisation
### 1.1 Principe
Pour chaque buffer (zone d'étude) :
1. Découpe le raster sur l'emprise du buffer
2. Masque avec la géométrie du buffer
3. Convertit les pixels = 1 en polygones
#### Pourquoi GDAL et pas R ?
La polygonisation avec GDAL en ligne de commande ou via Qgis est **plus rapide** qu'avec R sur de gros rasters. Les deux méthodes font appel aux mêmes algorithme GDAL, mais R doit passer par des couches *d'abstraction* pour convertir le code C/C++ natif en R et convertir les objets de C vers R. Pour certains traitement SIG massif il est conseillé de passer directement par GDAL.
### 1.2 Execution
1. Ouvrir le terminal OSGeo4W Shell, natif avec l'installation qgis
2. Lancer le script :
```bash
python chemin_vers_le_script/01_polygonise_buffers.py
```
### Paramètres à modifier (en haut du script)
```python
RASTER_INPUT = "chemin/vers/arbore.tif" # Votre raster
BUFFERS_INPUT = "chemin/vers/buffers.gpkg" # Vos zones d'étude
OUTPUT_DIR = "chemin/vers/sortie" # Dossier de sortie
```
---
## 2. Prétraitement
### 2.1 Principe
Les polygones issus d'un raster sont "pixelisés" avec des bords en escalier. Le prétraitement :
- **Supprime les petits polygones** (< 50 m²) → bruit
- **Lisse les contours** (Douglas-Peucker 2m) → formes naturelles, réduit le nombre de vertices
- **Bouche les petits trous** (< 50 m²) → artefacts
### 2.2 Execution
1. Ouvrir QGIS
2. Ouvrir la console Python : `Ctrl+Alt+P`
3. Coller :
```python
exec(open("C:/chemin/vers/03_preprocess_buffers.py").read())
```
### Paramètres
```python
DELETE_AREA = 50 # m² - polygones plus petits = supprimés
SIMPLIFY_TOL = 2.0 # m - tolérance de lissage
HOLE_AREA = 50 # m² - trous plus petits = bouchés
```
---
## 3. Classification
### 3.1 Principe
Pour ces calculs morphométriques, HedgeTools se base sur une approche inspirée de Touya et al. (2010), qui séparent la végétation arborée en classes selon la **taille** et la **forme** en deux étapes :
- Une première étape qui distingue les noyaux forestiers des autres objets selon un processus **d'ouverture morphologique**. L'étape d'érosion supprime par défaut toutes les surfaces de tous les objets dans un buffer de 15m, les haies et bosquets sont ainsi supprimés. L'étape de dilatation rend la surface originelle aux forêts qui ont été conservées.

- Une deuxième étape de classification des surfaces éliminés par la première étape, elle repose sur une combinaison de **mesures géométriques** :
- Taux de convexité
- Élongation du MBR
- Compacité (indice de Miller)
- Surface seuil pour distinguer les bosquets du reste
- la classification s'opère selon ce schéma de séquence (Touya et al, 2010)

---
| Classe | Critère | Description |
| --------------- | -------------------- | --------------------- |
| **Forêt** | Aire > 2000 m² | Grands massifs boisés |
| **Haie** | Forme allongée | Éléments linéaires |
| **Bosquet** | 200 < Aire < 2000 m² | Petits boisements |
| **Arbre isolé** | Aire < 200 m² | Arbres individuels |
### L'algorithme
1. **Érosion-dilatation** (20m) : sépare les éléments connectés
2. **Tri par taille** : forêt vs bosquet vs arbre isolé
3. **Tri par forme** : ce qui reste = haies (forme allongée)
### 3.2 Execution
Dans la console Python QGIS :
```python
exec(open("C:/chemin/vers/05_classify_buffers.py").read())
```
### Paramètres
```python
OPENING = 10 # m - largeur max des haies (érosion-dilatation)
FOREST_AREA = 2000 # m² - seuil forêt
GROVE_AREA = 200 # m² - seuil bosquet
```
> **OPENING = 20m** signifie que deux éléments séparés de moins de 20m seront considérés comme connectés. C'est aussi la largeur max typique d'une haie.
---
## 4. Extraction linéaire
### 4.1 Principe
On extrait l'**axe médian** des polygones de haies pour obtenir un réseau linéaire. Je développe pas plus cette partie là on est dans des calculs de géométrie un peu trop complexes pour moi que je n'ai pas encore creusé
```
Polygone haie Axe médian
┌─────────────┐
│ │ → ────────────
└─────────────┘
```
### 4.2 Execution
Dans la console Python QGIS :
```python
exec(open("C:/chemin/vers/07_lineaire_buffers.py").read())
```
### 4.3 Résumé visuel

### Paramètres
```python
MIN_WIDTH = -1 # -1 = pas de filtre (garder toutes les largeurs)
DANGLE = 30 # m
```
---
## 5. Simplification du linéaire, isolation des vergers
### 5.1 Problématique
les vergers "parasitent" le résultat. Ils sont surtout le signe que l'algorithme fonctionne bien : ce sont des arbres alignées. Mais ils s'organisent de manière plus irrégulière que les haies et sont plus courts.
````
Linéaire haie Linéaire verger
─────────────┐ ─────┐
└────── │
│
│
└────┐
│
│
└───┘
````
### 5.2 : Principe
Pour une même longueur de haies, le verger fait souvent plus de virages (constat visuel). Il s'agit alors d'estimer cette *tortuosité* ou *sinuosité*. On peut sommer ces changements de direction en mesurant leur angle et en le divisant par la longueur du segment via la métrique *Average Cumulative Angle* (Bíl et al., 2018).
>Cette étape est entièrement **expérimental** et **empirique**, pour bien faire il faudra vérifier la sortie de ce traitement avec nos vignettes digitalisées manuellement.
Toutes choses égales par ailleurs, la meilleure solution serait d'éliminer les vergers manuellement post traitement pour obtenir un résultat 100% fiable.
----
Cette étape se déroule en trois temps :
- Dans un premier temps, on cherche à réduire le nombre de vertices en simplifiant la géométrie du vecteur via l'algorithme de Douglas-Peucker (tolérance de 2 m), comme le préconise HedgeTools. Cette simplification permet d'éliminer les micro-variations liées à la vectorisation tout en conservant la forme générale du linéaire.
- Dans un deuxième temps, on calcule l'indice de giration *G*, adapté de l'*Average Cumulative Angle* de Bíl et al. (2018). Pour chaque linéaire, on mesure la somme des valeurs absolues des changements de direction entre segments consécutifs, normalisée par la longueur totale :
$$ACA = \frac{\sum_{i=1}^{n-2} |\alpha_i|}{L}$$
où :
- $\alpha_i$ représente le changement d'angle (en radians) entre le segment reliant les vertices $i-1$ et $i$, et le segment reliant les vertices $i$ et $i+1$
- $n$ correspond au nombre de vertices de la polyligne
- $L$ désigne la longueur totale du linéaire (en mètres)
L'indice *G* est exprimé en rad/m. Les valeurs élevées caractérisent des formes tortueuses avec de nombreux changements de direction (typiques des vergers en alignement), tandis que les valeurs faibles indiquent des structures plus linéaires (haies bocagères).
- Enfin, on mesure la longueur des segments de haie pour filter les petits linéaires qu'on retrouve dans les vergers (Déjà fait avec la BDHaie pour mon travail, après oeil expert de Daniel)
### 5.3 Execution
Toujours dans la console python de Qgis
```python
exec(open("C:/chemin/vers/05_classify_buffers.py").read())
```
**Il s'agit ensuite** de filtrer les haies en fonction de cet indice de Giration. Après vérification visuel avec Aurélia, on a proposé ce filtre là pour les buffers:
- Pour l'indice G : **0.06**
- Pour la longueur des haies à filtrer : **30 ou 50m**
Pour ces deux derniers points je vous laisser voir ça entre vous :)
Voici une commande rapide pour filtrer l'intégralité de vos buffers que vous aurez mis dans Qgis (**ctrl+alt+p** pour ouvrir la commande python)
````{p}
exec("""
filtre = '"giration" < 0.06 AND "longueur" > 50'
for layer in QgsProject.instance().mapLayers().values():
if layer.type() == QgsMapLayer.VectorLayer:
if "giration" in [f.name() for f in layer.fields()]:
layer.setSubsetString(filtre)
print(f"{layer.name()} - filtre applique")
""")
````
Libre à vous de modifier ces paramètres et de contrôler visuellement. Il ne vous reste plus qu'à extraire les entités sélectionnées pour vos buffers.
## Bibliographie
Bíl, M., Andrášik, R., Sedoník, J., & Cícha, V. (2018). ROCA – An ArcGIS toolbox for road alignment identification and horizontal curve radii computation. PLOS ONE, 13(12), e0208407. https://doi.org/10.1371/journal.pone.0208407
Marquès, G., Villierme, L., Boissonnat, J.-B., Guébin, G., Lang, M., Monteil, C., & Sheeren, D. (2024). HedgeTools : un outil d'analyse spatiale dédié à l'évaluation de la multifonctionnalité des haies. Mappemonde, 138. https://doi.org/10.4000/12jjm
Touya, G., Duchêne, C., & Mustière, S. (2010). Généralisation et intégration pour un fond vert commun entre l'IFN et l'IGN. Revue Internationale de Géomatique, 20(1), 65–86.