LidarTo3DPrint

Procédure de conversion d'un fichier de points LIDAR en une impression 3D

Contexte

Tenter de pouvoir observer et caractériser le relief d'un territoire au moyen des relevés 3D réalisés par l'état. Le territoire français est en cours de cartographie au Lidar avec une résolution de 10 points par mètre carré : https://geoservices.ign.fr/lidarhd

Ces données seront utilisées et transformées dans ce projet en une surface imprimable en 3D.

Réalisation

Les grandes étapes de la réalisation sont :

• Télécharger les données et si besoin, extraire une zone ou merger plusieurs fichiers ensembles
• Traiter le fichier et combler les "trous" de points et centrer les coordonnées
• Rassembler les points pour moyenner les mesures
• Créer la surface
• Découper la surface pour un faire un volume fermé
• Ouvrir dans l'outil d'impression 3D.

Télécharger les données

Se rendre sur la page web de téléchargement des données brutes : https://cartes.gouv.fr/telechargement/IGNF_NUAGES-DE-POINTS-LIDAR-HD

Sélectionner une zone et télécharger les données, ici nous sélectionnons une zone couvrant le centre ville de Plomelin

Rassembler les points

Afin de traiter les données LIDAR, nous allons utiliser l'outil QGIS et LAStools. Ces outils vont permettre :

• D'ouvrir le nuage de points
• De retirer les catégories difficiles à transformer en surfaces (les arbres et la végétation haute) pour ne garder que le sol, les voitures et les bâtiments.
• D'assembler plusieurs fichiers, ou de découper un sous-ensemble de points

Ouvrir QGIS et installer le plugin LAStools depuis QGIS, cela va intégrer les outils dans l'IHM de QGIS mais les executables sont à télécharger à part. Pour installer les outils dans QGIS commencer par dérouler le menu Plugins puis cliquer sur Manager et Installer des Plugins

Trouver le paquet LAStools et cliquer sur installer, cela va intégrer les outils dans QGIS.

Se rendre sur la page de LASTools et télécharger les outils compilés pour votre OS.

Décompresser le zip et copier le fichier LASTools dans le disque C:\LAStools puis redémarrer QGIS

Maintenant que les outils de traitement LIDAR sont installés on peut Importer un fichier LIDAR (ou plusieurs). Les fichiers sont importés comme des "Couches" dans le projet.

Dans QGIS, il est possible d’avoir un aperçu des fichiers en ouvrant une nouvelle vue 3D avec le menu View => 3D Map View => New 3D Map View

Sur ce zoom on peut apercevoir un château d’eau ainsi que des maisons et de la végétation. Les couleurs affichées dépendent des catégories. En effet lors du post-traitement effectué sur ces données LIDAR, un algorithme de catégorisation va identifier parmi les points : le sol, les bâtiments, la végétation, les lignes électriques... ou alors si rien n'est identifié placer les points dans Uncategorized.

Pour ne pas imprimer tout le village, nous allons sélectionner une plus petite zone. Pour cela commencer par créer une nouvelle couche temporaire, cette couche va servir à dessiner un polygone, qui sera utilisé pour découper le fichier de points.

Sélectionner une géométrie de type polygon (pour faire la sélection) et choisir la même géodésique que celle du nuage de points (vérifier qu'il n'y ait pas de logo de type point d'exclamation à côté du menu déroulant).

Une nouvelle couche apparaît maintenant en bas a gauche, nous allons y dessiner un polygone qui contiendra uniquement les points à garder.

Cliquer sur "Add Polygon" pour ajouter un polygone et sélectionner une région a garder pour impression 3D.

Utiliser l’outil d’ajout de polygone (click gauche pour les points du polygone, click droit pour terminer)

Chercher l’outil de processing « Point Cloud Clip »

Sélectionner la couche de points source ainsi que le polygone

Cela va sélectionner les points contenus dans le polygone et créer une nouvelle couche. Ne cocher que la couche « Clipped » pour retenir les points sélectionnés.

Si au contraire on souhaite assembler plusieurs fichiers de points pour représenter une grande surface, cela reste possible avec la fonction Point cloud data management ==> Merge

Désélectionner les catégories non souhaitées

Il est possible maintenant d’utiliser les classes de ces données lidar pour retirer ce que l’on ne va pas imprimer. Ici je ne garde que le sol et les bâtiments

En cochant / décochant des cases, il n’est possible que d’afficher ou de masquer des points. Cependant pour les étapes suivantes nous souhaitons les retirer du jeu de données avant de transférer le fichier dans un autre logiciel

Faire un click droit sur la couche et cliquer sur Filter… créer un filtre en ne gardant que les classifications 2 et 6

Faire un click droit sur la couche, Export => Save As et enregistrer dans un fichier LAS

Vous n’aurez que les points associés aux catégories filtrées dans ce fichier

Recommencer l’opération et cette fois ci enregistrer les données au format « Comma Separated Values »

Combler les trous

Afin de combler les trous un utilisateur à créer un programme en C# permettant de traiter un nuage de points LIDAR. Son procédé date de 2018 ainsi que les outils qu'il utilisait, il décrit cela sur son blog, et propose une archive du code qu'il à utilisé sur un dépot Github. Son projet à été cloné ici au cas où il disparaîtrais.

L'idée est de l'utiliser afin de combler les vides, par exemple ici à gauche, le nuage de points originel, et à droite le nuage "corrigé"

On peut voir que sur le nuage originel, pris par avion, certaines facades de bâtiments sont vides de points. Cela rend la reconstruction compliquée...

L'outil créer par marian42 permet de combler ces vides en appellant les lignes de commande bash

:: This script will search through pointcloud data in the data folder and generate a 3D printable model.
:: Please update these parameters before using it:

:: The center of the square to extract, using the same coordinate system as the XYZ data supplied.
SET x=166578.91
SET y=6783615.12

:: The size of the square to extract, in meters
SET size=100

:: Add some space between the lowest part of the surface and the bottom of the 3D mesh. In meters.
SET verticaloffset=4

:: The directory where the program will look for .xyz files
SET datadirectory=data

:: Make sure Meshlab and Blender are installed and put their executable locations here.
SET meshlab="C:\Program Files\VCG\MeshLab\meshlabserver.exe"
SET blender="C:\Program Files\Blender Foundation\Blender 2.82\blender.exe"

:: Nothing to configure below this.

pointcloudtool.exe extract %datadirectory% pointcloud.xyz %x% %y% %size% || goto :error

pointcloudtool.exe fix pointcloud.xyz heightmap.xyz|| goto :error

Créer la surface

Maintenant que le nuage de points est corrigé, on cherche à en créer la surface, pour cela ouvrir « Meshlab » et importer le fichier traité par le programme de marian42. Le fichier xyz comporte les latitudes, longitudes (centrées autour du point défini dans le script batch) et l'élévation séparées par une virgule.

En regardant en haut a droite on peut observer que nous avons beaucoup de points mais aucune face. Il faut reconstruire les faces d’après tous les points en utilisant le bon filtre et potentiellement des corrections manuelles.

Ces points ont été relevés par avion et ont un certain bruit de mesure. Partir de ce nuage pour reconstruire le modèle va créer des incohérences étranges dues a une mauvaise interprétation des normales (qu'est ce qui est "dedans" et ce qui est "dehors"). Sur l'image ci-dessous une partie du toit du building de gauche à été inversée comme "dehors"

Pour réduire ces effets, on peut moyenner les points. Cela va induire une perte de "résolution" et de détails sur le modèle, mais permet d'éviter les formes incohérentes.

Pour cela on va simplifier le nuage de points en regroupant les points, je trouve que de garder 0,1% sur ces relevés donne de bons résultats. Avec 0,09% j'ai toujours des problèmes donc 0,1 est peut être le minimum pour ne pas avoir de défauts sans tout autant perdre trop de résolution.

Il faut ensuite lancer une estimation des normales. Cet algorithmes va tenter d'identifier ce qui est dedans et ce qui est dehors, ainsi que le vecteur orthogonal aux points. Je trouve qu'une valeur qui fonctionne bien est d'utiliser une quantité de 40 points par normale.

Finalement on peut utiliser la reconstruction de Poisson pour identifier la surface qui correspond le mieux a tous ces points, toutes ces normales et leur orientation. Ici je trouve qu'une valeu qui fonctionne bien est d'utiliser les paramètres suivants:

Enregistrer maintenant la surface au format STL en utilisant le menu File => Export As et l’appeler "mesh.stl"

Découper la surface

Le reste de l’algorithme de marian42 va maintenant créer un cube de côté égale au paramètre donné dans le script, découper la surface suivant ce cube, et aussi créer "un fond" pour fermer le modèle et en faire un solide. Pour cela appeler ce code batch

Ce bout de code nécessite d'avoir installé Blender sur son ordinateur

:: This script will search through pointcloud data in the data folder and generate a 3D printable model.
:: Please update these parameters before using it:

:: The center of the square to extract, using the same coordinate system as the XYZ data supplied.
SET x=166578.91
SET y=6783615.12

:: The size of the square to extract, in meters
SET size=100

:: Add some space between the lowest part of the surface and the bottom of the 3D mesh. In meters.
SET verticaloffset=4

SET blender="C:\Program Files\Blender Foundation\Blender 2.82\blender.exe"

pointcloudtool.exe makeSolid mesh.stl mesh-solid.stl cube-solid.stl %size% %verticaloffset% || goto :error

%blender% -b -P intersect.py -- mesh.stl cube.stl mesh.stl || goto :error

Cela va créer un modèle STL importable dans Bambu Studio

Ouvrir l'objet pour impression 3D

Ouvrir Bambu Studio et importer le STL, le redimensionner et lancer l'impression.