Segmenter une couche d'éléments linéaires de type polygone avec ArcMap

	Niveau	Débutant
	Logiciels utilisés	ArcMap 10.X
	Plateforme	Windows

Quand on fait de l'analyse spatiale dans un SIG, il arrive souvent que l'on ait besoin de découper une couche en fragments homogènes afin d'appliquer géotraitrements ou analyses à chacun de ces fragments comparables entre eux. Dans de nombreux cas on peut utiliser la création d'un maillage pour les couches de type polygones ou d'un fractionnement pour les couches de type polyligne mais comment faire quand on veut segmenter une couche d'éléments linéaires de type polygone (par exemple, la surface des cours d'eau d'un bassin versant, la surface d'un réseau de voirie, etc.) tous les x mètres ?

Voici une méthode permettant cette segmentation en utilisant les outils disponibles dans ArcMap avec l'aide d'extentions disponibles gratuitement.

Les couches en entrée

Pour réaliser cette méthode, il vous faudra une couche de type polyligne correspondant au "squelette" du réseau à segmenter et d'une couche de type polygone correspondant à la surface de ce réseau.

Sachez que si vous ne disposez pas de la couche de type polyligne correspondant au "squelette" du linéaire, il existe à ma connaissance une méthode de squelettisation basée sur les diagrammes de Voronoï détaillée par Jérôme Rolland sur son blog  et un plugin pour le logiciel SIG libre OpenJump : skeletonizer 1.0 

ETAPE 1 : Créer une couche de points suivant le linéaire tous les "x" mètres avec pour attribut l’azimut de la ligne sur laquelle ils sont

1 – Fusionner toutes les entités du linéraire "squelette" de type polyligne pour n’avoir plus qu’une entité (Toolbox gestion de données -> généralisation -> fusionner). Faire ensuite une copie de la couche obtenue.

2 – Diviser le linéaire en tronçon de "x" mètres : menu éditeur -> fractionner : choisir « en parties égales et diviser la longueur total de la ligne par "x" pour avoir des tronçons de "x" mètres (dans notre exemple nous choisiront de faire une segmentation tous les 10 mètres.

3 – Créer des points à chaque tronçon de polylignes (extension ET Geo Wizards de Ian-ko Convert –> Polyline to Point, utiliser l’option « nodes »).

4 – Repartir de la copie du squelette fusionné, créer des points à partir chaque sommet (utiliser le même outil que pour l’étape 1.3 mais avec l’option Vertices).

5 – Découper le linéaire fusionné par les points créée à l’étape 1.4 avec l'outil "Fractionner des lignes par des points" (Toolbox gestion de données -> Entités -> Fractionner des lignes par des points)

6 – Chaque polyligne est en fait devenue une ligne avec une seule orientation (un seul azimut). Nous allons donc, pour chaque ligne, renseigner son azimut dans la table attributaire. L’add-in ET EasyCalculate 10 de Ian-ko nous permet de faire cela (Category : Geometry / Expression : polyline_Get_Azimuth / Point Along = 0)

7 – Faire ensuite une jointure spatiale (clic droit sur la couche => jointures et relations => joindre... => joindre les données d'une autre couche selon l'emplacement)  pour que chaque point de la couche créée à l’étape 1.3 (un point tous les 10 mètres) ait dans sa table attributaire l’azimut de la ligne sur laquelle il se trouve. 

ETAPE 2 : Calculer les coordonnées des points x1 y1 x2 y2 de lignes passant par les points situés tous les 10 mètres et étant perpendiculaire au linéaire voirie

1 – Créer et calculer les champs X et Y de chaque point issu de la couche de l’étape 1.7 (dans la table attributaire clic-droit sur le champ puis choisir "calculer la géométrie").

2 – Créer les champs dx et dy dans lesquels nous allons calculer les différences entre les coordonnées X et Y des l’extrémités des lignes que l’on veut créer et les coordonnées X,Y des points d’origine.
Le calcul à réaliser est issu des formules de trigonométrie pour les triangles-rectangle :
 

 

Ici nous cherchons [AC] (dx) et [BC] (dy) en connaissant  (azimut) et [AB] (la distance d’une demi-ligne, ici je veux tracer des lignes de 20 mètres alors j’utiliserai c = 10)

[AC] = cos (Â) * [AB]
dx = cos(azimut) * 10

[BC] = sin(Â) * [AB]
dy = sin(azimut) * 10

Voici les formules à entrer dans la calculatrice de valeurs de champs avec pour Analyseur Python :

Dx=math.cos((!azimut!*math.pi)/180)*10

Dy=math.sin((!azimut!*math.pi)/180)*10

!azimut! est le nom de champs où sont stockées les azimuts
*math.pi/180 sert à transformer les angles degrée en radian

3 – On crée maintenant les champs x1, y1, x2, y2 pour accueillir les coordonnées des extrémités de nos futures lignes.

Voici les calcul à effectuer avec la calculatrice de valeurs de champ :

x1=x+dx
y1=y-dy
x2=x-dx
y2=y+dy

Nous avons maintenant les coordonnées des lignes que l’on va créer pour couper notre couche linéaire de type polygone correspondant à la surface de notre réseau.

ETAPE 3 : Tracer les lignes et découper le linéaire de type polygone

1 – pour tracer nos lignes nous allons utiliser l'outil "XY vers lignes" dans l'ArcToolbox (Toolbox gestion de données -> Entités -> XY vers lignes). Nous avons donc une nouvelle couche contenant les lignes que nous allons utiliser pour découper notre linéaire de type polygone.

2 – Pour découper il y a trois solutions :

• Soit on utilise l’outil entité vers polygone (toolbox Outils de gestion de données-> entités) mais cet outil nécessite une licence ArcInfo.
• Soit on utilise l’outil Split Polygons with Polylines (ET Geo Wizards dans le menu Polygon) mais cet outil est limité à 100 entités en version gratuite.
• Soit quand comme moi, on a qu'une licence ArcView, on fusionne toutes les entités de la couche de lignes de découpe puis on se met en édition, on sélectionne le polygone de notre couche linéaire de type polygone, puis avec l’outil « découpe de polygone » on se place sur n'importe quelle « ligne de découpe », on clique droit et on choisi « remplacer la construction ».

Conclusion

Nous obtenons le résultat qui nécessite quelques retouches manuelles à certains endroits mais qui, dans l’ensemble, est assez propre.

A noter qu'il existe d'autres moyens qu'ArcMap pour réaliser cette méthode, voir par exemple ce sujet sur gis.stackexchange.com qui utilise le module Python Shapely, signalé par Gene sur le forumSIG.

Site officiel : ArcMap
Autres Liens : Jérôme Rolland - Du diagramme de Voronoï à la Squelettisation de polygone
Autres Liens : Plugin Skeletonizer 1.0 pour Open Jump
Autres Liens : Generates multiple profile lines perpendicular to an input line #