Script Python pour parser un fichier laz et afficher son emprise dans un fichier html

Avant de pouvoir jouer le script il faut au préalable disposer des librairies nécessaires.

!pip install laspy lazrs

ensuite on va lancer le script : 

info laz 

Que va faire le script ?

Le script analyse un fichier LiDAR (.laz) et génère une carte interactive affichant l'emprise de la dalle. Voici ce qu'il fait :

1. Lecture du fichier .laz : Le script demande à l'utilisateur de saisir le chemin complet du fichier .laz et ouvre ce fichier pour en extraire des informations.

2. Extraction des informations :

   • Version du fichier
   • Nombre de points
   • ID du format de point
   • Attributs des points
   • Facteurs d'échelle
   • Offsets
   • Limites (X, Y, Z)
   • Informations de projection (tronquées si trop longues)
   • Heure de début et de fin d'acquisition
   • Dimensions de la dalle (largeur, hauteur, plage d'altitude)
   • Surface de la dalle
   • Densité de points
   • Point le plus bas et le plus haut en Z
   • Comptage des points par classe (sol, végétation, bâtiment, etc.)

3. Transformation des coordonnées : Les coordonnées des coins de la dalle sont transformées en latitude et longitude.

4. Création de la carte interactive :

   • Utilisation de Folium pour créer une carte centrée sur l'emprise de la dalle.
   • Ajout de l'emprise de la dalle sous forme de polygone sur la carte.
   • Affichage des informations extraites dans une infobulle lorsque l'emprise est cliquée.

5. Affichage des informations et sauvegarde de la carte :

   • Les informations extraites sont affichées dans la console.
   • La carte interactive est sauvegardée dans le même répertoire que le fichier .laz sous le nom lidar_extent_map.html.

Quand il est lancé alors il faut juste lui glisser le chemin vers le fichier et le script va commencer à regarder la structure.

Cela peut prendre un peu de temps car les fichiers .laz de l'IGN sont lourds. Le laz est la version compressée du las donc quand le script tourne il créée un backend pour ensuite allez lire les datas. Donc sur des fichiers qui pèsent 200 Mo ça nécessite un peu de patience.

Une fois fini il affiche les informations sur le fichier .laz et créée un fichier que j'ai appelé  lidar_extent_map.html.

 

Exemple de sortie sur la dalle : LHD_FXX_0417_6927_PTS_C_LAMB93_IGN69.copc.laz : 

Veuillez entrer le chemin complet et le nom du fichier .laz à analyser : /content/drive/MyDrive/Colab Notebooks/projets/lidar_hd/donnees/LHD_FXX_0417_6927_PTS_C_LAMB93_IGN69.copc.laz

Nom du fichier : LHD_FXX_0417_6927_PTS_C_LAMB93_IGN69.copc.laz

Version : 1.4 Nombre de points : 10981064

ID du format de point : 6

Attributs des points :

- X

- Y

- Z

- intensity

- return_number

- number_of_returns

- synthetic

- key_point

- withheld

- overlap

- scanner_channel

- scan_direction_flag

- edge_of_flight_line

- classification

- user_data

- scan_angle

- point_source_id

- gps_time

Facteurs d'échelle : [0.01 0.01 0.01]

Offsets : [0. 0. 0.]

Limites :

X : [417000.0, 417999.99]

Y : [6926000.0, 6926999.99]

Z : [-1.92, 59.72]

PROJCS["RGF93 v1 / Lambert-93",GEOGCS["RGF93 v1",DATUM["Reseau_Geodesique_Francais_1993_v1",SPHEROID["GRS 1980",6378137,298.257222101,AUTHORITY["EPSG","7019"]],AUTHORITY["EPSG","6171"]],PRIMEM["Greenwich",0,AUTHORITY["EPSG","8901"]],UNIT["degree",0.0174532925199433,AUTHORITY["EPSG","9122"]],AUTHORITY["EPSG","4171"]],PROJECTION["Lambert_Conformal_Conic_2SP"],PARAMETER["latitude_of_origin",46.5],PARAMETER["central_meridian",3],PARAMETER["standard_parallel_1",49],PARAMETER["standard_parallel_2",44],PARAMETER["false_easting",700000],PARAMETER["false_northing",6600000],UNIT["metre",1,AUTHORITY["EPSG","9001"]],AXIS["Easting",EAST],AXIS["Northing",NORTH],AUTHORITY["EPSG","2154"]]

Heure de début d'acquisition : 2022-09-12 04:44:21.803614

Heure de fin d'acquisition : 2022-09-12 13:32:11.685788

Largeur de la dalle (m) : 999.9899999999907

Hauteur de la dalle (m) : 999.9900000002235

Plage d'altitude de la dalle (m) : 61.64

Surface de la dalle (m²) : 999980.0001002142

Densité de points (points/m²) : 10.981283624572011

Point le plus bas (Z) : -1.92

Point le plus haut (Z) : 59.72

Classes et nombre de points associés :

Classe 1 (Non classé): 73894 points

Classe 2 (Sol): 4725803 points

Classe 3 (Végétation basse): 186544 points

Classe 4 (Végétation moyenne): 512541 points

Classe 5 (Végétation haute): 1279158 points

Classe 6 (Bâtiment): 83470 points

Classe 9 (Eau): 4118231 points

Classe 65 (Artefact): 994 points

Classe 67 (User-defined or Reserved): 429 points

Carte interactive créée : /content/drive/MyDrive/Colab Notebooks/projets/lidar_hd/donnees/lidar_extent_map.html

 

Il a créé le fichier html dans le même path.

Attention ! concernant les APRS il s'agit de ceux propres aux spécifications de l'IGN dans ce cas précis. Si vous testez un lidar d'une autre provenance, possible que cela ne fonctionne pas de façon optimum

Quand vous ouvrez le fichier html vous devriez avoir ceci :

et au clic sur la dalle : 

Toutes les informations concernant cette dalle

Un autre exemple avec une autre dalle près de Saint Guillhem le Désert : 

et vérification avec Diffusion LiDAR HD (ign.fr)

Voilà c'est pratique quand même quand on est loin du site et/ou de son SIG comme je le mentionnais au début ça peut permettre de voir s'il n'y a pas de souci dans les coordonnées géographiques par exemple