Le portail francophone de la géomatique

import rasterio
import numpy as np
from tqdm import tqdm
import os

# Chemin vers le dossier contenant les fichiers raster
dossier_source = r'chemin/vers/votre/dossier'

# Liste des fichiers raster à traiter
liste = ['chemin/vers/votre/fichier1.tif', 'chemin/vers/votre/fichier2.tif']

for fichier in liste:
    # Chemin vers le fichier de sortie modifié
    output_raster_path = fichier[:-4] + "_modifie.tif"
    
    # Passe si un fichier de sortie de même nom a déjà été généré
    if os.path.exists(output_raster_path):
        print(output_raster_path + " existe déjà !")
        continue
    
    # Lire le raster
    with rasterio.open(fichier) as src:
        # Lire les données raster
        raster_data = src.read(1)

        # Copier les données raster pour effectuer les modifications
        modified_raster_data = raster_data.copy()

        # Tableau pour stocker les groupes
        group_table = np.zeros_like(raster_data, dtype=np.int32)

        # Dictionnaire pour stocker les coordonnées des pixels de contour pour chaque groupe
        contour_pixels = {}

        # Valeur de groupe initial
        current_group = 1

    # Parcourir tous les pixels du raster en utilisant un parcours en profondeur (DFS) itératif
    def dfs(i, j, current_group):
        stack = [(i, j)]
        while stack:
            i, j = stack.pop()
            if (0 <= i < raster_data.shape[0] and
                    0 <= j < raster_data.shape[1] and
                    raster_data[i, j] == 4 and
                    group_table[i, j] == 0):
                group_table[i, j] = current_group
                contour_pixels[current_group].append((i, j))

                # Ajouter les voisins non marqués à la pile
                stack.extend([
                    (i, j-1),  # Gauche
                    (i-1, j-1),  # Haut gauche
                    (i-1, j),  # Haut
                    (i-1, j+1),  # Haut droite
                    (i, j+1),  # Droite
                    (i+1, j+1),  # Bas droite
                    (i+1, j),  # Bas
                    (i+1, j-1)  # Bas gauche
                ])

    # Parcourir tous les pixels du raster
    for i in tqdm(range(raster_data.shape[0]), desc='Processing'):
        for j in range(raster_data.shape[1]):
            if raster_data[i, j] == 4 and group_table[i, j] == 0:
                current_group += 1
                contour_pixels[current_group] = []  # Initialiser la liste des pixels de contour pour ce groupe
                dfs(i, j, current_group)

    # Parcourir les groupes pour vérifier le pourcentage de pixels voisins des pixels de contour
    for group, pixels in tqdm(contour_pixels.items(), desc='Processing contours'):
        contour_array = np.array(pixels)
        contour_i, contour_j = contour_array[:, 0], contour_array[:, 1]

        # Obtenir les indices des voisins (les 8 voisins)
        neighbors_i = np.array([contour_i-1, contour_i-1, contour_i-1, contour_i, contour_i, contour_i+1, contour_i+1, contour_i+1])
        neighbors_j = np.array([contour_j-1, contour_j, contour_j+1, contour_j-1, contour_j+1, contour_j-1, contour_j, contour_j+1])

        # Filtrer les voisins valides
        valid_neighbors = np.logical_and(
            np.logical_and(neighbors_i >= 0, neighbors_i < raster_data.shape[0]),
            np.logical_and(neighbors_j >= 0, neighbors_j < raster_data.shape[1])
        )
        neighbors_i = neighbors_i[valid_neighbors]
        neighbors_j = neighbors_j[valid_neighbors]

        # Obtenir les valeurs des voisins
        neighbor_values = raster_data[neighbors_i, neighbors_j]

        # Compter les voisins ayant la valeur 3
        imper_neighbors = np.count_nonzero(neighbor_values == 3)

        # Compter les voisins ayant une valeur différente de 4 pour obtenir le nombre total de voisins
        total_neighbors = np.count_nonzero(neighbor_values != 4)

        # Calculer le pourcentage
        percentage = (imper_neighbors * 100) / total_neighbors

        if percentage > 60:
            # Changer les valeurs des pixels du groupe de 4 à 3 dans raster_data
            group_indices = np.where(group_table == group)
            raster_data[group_indices] = 3

    # Enregistrer le raster modifié
    with rasterio.open(output_raster_path, 'w', **src.meta) as dst:
        dst.write(raster_data, 1)

print("Programme terminé")