Tutoriel avancé QGIS et Python

Détails
75602

Index de l'article

Qgis3PythonMots-clefs : interopérabilité Python/QGIS, variable, boucle, fonction, condition, liste, argument, chaîne formatée, fonctions géométriques, création de fichiers, génération de cartes, fonctions cartographiques, fonctions standalone, traitement de données cartographiques, workflow, milieu alpin

Enseignant responsable : Georges Hinot

Pré-requis : une bonne connaissance d'un logiciel SIG supportant Python (QGIS ou ArcGIS par exemple) et des problématiques/besoins/enjeux de la géomatique.

Compétences visées : mettre en place un protocole de développement autour d'une problématique géomatique, de la « re-contextualisation du besoin » jusqu'à l'écriture d'un script Python répondant à ce besoin.

Enjeux du cours : rendre les étudiants autonomes dans la création et l'entretien de scripts « automatisants ». Dans ce type de tâche, l’autonomie consiste à savoir interpréter un besoin humain, ré-interpréter ses propres problématiques de programmation et à écrire/chercher/modifier les blocs de code nécessaires.

Programme pédagogique : initiation aux fondamentaux du langage puis mise en situation (générer des cartes à partir de données géographiques). En plus de la simple automatisation des cartographies, le cours nous donnera l'occasion, toujours en Python, d'utiliser des API, de traiter des données et de mettre en place des modes d'affichages customisés sur QGIS. Selon l'avancée du cours, nous pourrons aussi aborder le mode « standalone » permis par le duo QGIS-Python, la création de plugins QGIS et les possibilités de webmapping/analyse/ingénierie de données permises par Python.

Sources et liens divers :

Tutoriel écrit sur QGIS 3.14 'Py' et Python 3.8, mais j'ai pu constater qu'il fonctionne de façon assez similaire sur des versions précédentes, ainsi que sur QGIS 3.16.

Sur des versions plus récentes comme QGIS 3.2 il semble que certaines API n'aient pas encore été portées (wikipedia et wikidata), mais cela ne posera pas de problème pour 99% du tutoriel.

 


Découverte de la console Python dans QGIS

Ce 1er chapitre est inspiré de l'excellent site qgistutorials.

Ouvrez les couches peaks et eau (jointes en bas de cet article) dans un projet QGIS, ainsi que la console Python du logiciel (un bouton à cet effet est disponible dans le menu QGIS).

La couche de points en premier plan, le linéaire en second, avec une symbologie évidente. Ce sont les principaux sommets des Écrins et le réseau hydrographique français*.

* Les sommets sont issus d'une rapide extraction depuis OSM, et le réseau hydrographique est une couche IGN peu à jour. Les sommets sont encodés en UTF-8, pour bien lire les accents, précisez-le éventuellement dans QGIS. Pour augmenter les difficultés, ces 2 couches sont dans des projections différentes.

Sélectionnez la couche de points dans le panneau calque de QGIS, puis dans la console Python entrez :

layer = iface.activeLayer()

Vous venez de créer une variable arbitrairement nommée layer et contenant la couche active (celle que vous avez sélectionnée dans QGIS).

iface est une variable globale, qui n'est autre que le canvas courant de votre projet QGIS, les données visibles dans votre projet. Vous n'avez pas créé cette variable vous-même, mais elle s'est automatiquement créée à l'ouverture de QGIS.

activeLayer est une méthode de Python, une action à réaliser, qui ici va chercher la couche active de votre projet. Nous reviendrons bientôt sur le point (.) une fois que les conditions de notre compréhension seront plus explicites.

Votre 1ère variable existe donc déjà, et vous pouvez en connaître le contenu simplement puisque :

layer

Vous renverra quelques informations sur le contenu de cette variable :

<QgsMapLayer: 'peaks' (ogr)>

Vous venez de lire (d'un point de vue système) le contenu de la variable. Si vous sélectionnez maintenant la couche linéaire, et lancez les mêmes commandes, vous lirez bien-entendu un contenu différent pour la même variable. En effet elle contiendra autre chose, la méthode activeLayer allant chercher la couche sélectionnée dans QGIS. Bienvenue dans Python pour QGIS !

Astuce
Depuis la console, nul besoin de retaper les mêmes commandes déjà utilisées. Utilisez la flèche ↑ (flèche du haut) de votre clavier pour remonter la liste des commandes exécutées.

Vous pouvez aussi interroger les quelques 500 méthodes disponibles pour cet objet :

dir(layer)

Ci-dessus, la variable layer est soumise en tant que paramètre (ce mot n'est pas anodin en programmation) de la fonction dir, entre parenthèse. La fonction sait donc sur quoi s'appliquer.

L'une des méthodes disponibles ets la méthode fields(), nous pouvons consulter toutes ces options sur la documentation officielle pyqgis :

https://qgis.org/pyqgis/3.0/core/other/QgsFields.html

Exemples :

Compter le nombre de champs

layer.fields().count()

Afficher les noms des champs

layer.fields().names()

Les listes

La méthode names() nous renvoie les champs entre crochets, c'est donc une liste (ce mot n'est pas anodin en programmation).

Nous pouvons accéder aux élements indicés des listes en ajoutant [NUMÉRO DE L'INDICE]. Pour afficher le 1er champ par exemple :

layer.fields().names()[0]

Astuce
Les éléments d'une liste sont ordonnés de 0 à X, le 1er élément  ayant toujours le numéro 0. Comme en histoire, quand nous sommes en 50 après JC, nous sommes au 1er siècle après JC. Et quand nous avons 20 ans, nous sommes dans notre 21ème année.

Boucle et itération

L'une des méthodes disponibles pour nos couches s'appelle getFeatures. Nous allons nous en servir pour boucler sur les 341 entités contenues dans la couche de points :

layer = iface.activeLayer()
 
for f in layer.getFeatures():
    print(f['name'], f['osm_id'])

Le mot-clé for signifie que nous démarrons une boucle. Chaque contenu de cette boucle sera stocké dans une variable f (arbitrairement nommée ainsi, mais correspondant à une convention Python, entre autre). Le in regarde où lancer la boucle : sur les entités contenues dans notre variable.

Et qu'est-ce que nous y lancerons ? Des instructions itératives.

Astuce
Si à ce stade la boucle ne s'est pas exécutée, c'est que vous devez sortir de la boucle, en tapant sur la touche ↵ (Entrée).

Cheminement

Dans la 2ème ligne notamment, remarquez l'usage du point (.), qui concatène variables et instructions Python (là ou PHP concatène les chaînes).

Ce n'est donc plus une simple concaténation, mais plutôt le signe d'un cheminement du code : Python voit votre variable layer, mais étant suivie d'un . il sait qu'il y a quelque chose à faire dessus. Quoi donc ? Y lancer la méthode qui suit.

Mais au fait, si le point sert à Python pour concaténer des instructions, comment concatène-t-il les chaînes de texte ? Avec le signe +.

Indentation

Observez l'indentation devant la 3ème ligne. L'indentation est partie intégrante du langage Python, elle permet de structurer son code et plus encore, de signifier des choses dans le code.

Ici l'indentation précédant le print signifie qu'il est physiquement contenu dans la boucle.

L'indentation a donc son rôle à jouer. Avec Python une mauvaise indentation soulèvera une erreur d'exécution (cela contrairement à PHP, où l'indentation n'est là que pour aider à la lecture, au mieux).

Après saisie complète du code ci-dessus dans la console, pour exécuter la boucle, il vous faudra taper une dernière fois sur ↵ (Entrée) , pour sortir de la boucle. Les noms des principaux sommets de la zone de glaciers des Écrins, et leur identifiant OSM, défilent dans la fenêtre de résultats.

Astuce
Une erreur courante est parfois soulevée, suite à des copier-coller ou des modifications du code depuis un autre logiciel :
 

Inconsistent use of tabs and spaces in indentation.

Le log est explicite : il s'agit d'un mélange de tabulations et d'espaces dans votre indentation. Visuellement votre code est propre, mais pas pour la machine, qui voit TOUT. En effet si vous avez commencé par utiliser des espaces ou des tabulations pour indenter, alors vous devez vous y tenir. De même si vous utilisez des indentation à 2 espaces, alors 2 indentations feront 4 espaces ; si vous faîtes des indentations à 4 espaces, alors 2 indentations feront 8 espaces. Souvenez-vous en pour la suite de ce tutoriel !

Corrigez donc la ligne concernée et prenez de bonnes habitudes : Notepad possède un paramètre pour signifier quel caractère systématiquement utiliser dans les indentations de la touche Tab du clavier.

1ère fonction géométrique

De même que nous affichons les noms des sommets, nous pouvons les géocoder, en listant leurs coordonnées géographiques (xy) :

for f in layer.getFeatures():
    geom = f.geometry()
    print(geom.asPoint())

Astuce
Si à ce stade la boucle ne s'est pas exécutée, c'est que vous devez sortir de la boucle, en tapant sur la touche ↵ (Entrée).

Nous avons cette fois créé une variable directement dans notre boucle (f).

La méthode geometry récupère certaines des informations spatiales de nos entités, toutes stockées et bien organisées dans notre variable geom, même si nous ne les voyons pas.

Puis le print interroge notre variable sur l'un de ses contenus, asPoint, méthode qui récupére les latitudes et longitudes.

  • Note personnelle : extraire un réseau hydrographique plus complet depuis OSM, ainsi que les sommets tagués différemment (les pointes par exemple).

Écriture de fichiers

Écriture de fichier Text

Nous allons maintenant exporter les noms des sommets, leurs identifiants OSM ainsi que leur coordonnées géographiques dans un joli fichier texte.

Avant cela, récupérez l'emplacement de votre machine où vous souhaitez écrire votre fichier, car nous mentionnerons son chemin complet dans la 1ère ligne ci-dessous, suivi du nom de fichier voulu.

Vous pouvez aussi écrire directement dans le répertoire courant (si vous le connaissez ! Voir plus bas), dans ce cas le seul nom du fichier suffira :

Astuce
Remarquez que nous utilisons de simples slashs pour définir nos chemins, contrairement à Windows, qui utilise des slahs inversés. SI vous récupérez vos chemins via les Propriétés de vos fichiers Windows, vous devez donc remplacer vos slahs inversés par de simples slahs, ou alors doubler vos slahs inversés (\\).

  1. with open('C:/Users/Georges/Downloads/SommetsEcrins.txt', 'w') as file:
  2.     for f in layer.getFeatures():
  3.         geom = f.geometry()
  4.         line = '{};{};{:.4f};{:.4f}\n'.format(f['name'], f['osm_id'], geom.asPoint().y(), geom.asPoint().x())
  5.         file.write(line)

Astuce
Si à ce stade la boucle ne s'est pas exécutée, c'est parce que ...

Chaînes formatées

Observez la syntaxe de la 4ème ligne ci-dessus, un peu particulière.

Dans la variable line nous indiquons d'abord des zones vides ({}), ce sont les zones dans lesquelles nous souhaitons écrire les valeurs des champs. Nous en profitons parfois pour y indiquer un 1er formatage. La notation :.4f signifie qu'ici les valeurs doivent être arrondies à 4 chiffres décimaux.

Nous utilisons également le caractère spécial \n issu des expressions régulières, pour faire en sorte qu'après chaque écriture d'entité, une ligne soit sautée dans le fichier avant d'écrire la suivante.

Nous y retrouvons donc 4 blocs avant le .format().

Ensuite seulement, après le .format(), nous y indiquons ce que les blocs doivent contenir, en respectant le même ordre (le 1er champ mentionné sera contenu dans la 1er bloc, et ainsi de suite).

Enfin, nous lançons l'écriture de notre variable line.

Il y a plusieurs façons de formater des chaînes avec Python. Ici nous utilisons la méthode {} mais nous en verrons d'autres.

Écriture de CSV

Pour écrire un beau fichier CSV nous allons procéder différemment, sachant par ailleurs qu'il existe plusieurs façons de faire.

  1. out_file = open('C:/Users/Georges/Downloads/SommetsEcrins.csv', 'w', encoding='cp1252')
  2. out_file.write("name" + "," + "osm_id" + "," + "y" + "," + "x" + "\n")
  3.  
  4. for f in layer.getFeatures():
  5.     geom = f.geometry()
  6.     line = '{},{},{:.4f},{:.4f}\n'.format(f['name'], f['osm_id'], geom.asPoint().y(), geom.asPoint().x())
  7.     out_file.write(line)
  8.  
  9. out_file.close()

Astuce
Avec le jeu des copier-coller depuis ce site web ou un éditeur de texte gérant différemment les indentations, séparateur Tab, whitespace et saut de ligne, il peut y avoir des erreurs d'interprêtation de Python.

Dans ce cas utilisez plutôt l'éditeur Python de QGIS (survolez les boutons disponibles en haut de la console Python de QGIS, l'éditeur s'ouvrira à droite de la console) et exécutez-y le code ci-dessus. Une meilleure interprêtation du code est à l'œuvre depuis l'éditeur Python.

La 1ère ligne crée puis ouvre (d'un point de vue système) un fichier CSV vide. Nous ajoutons un encodage cp1252 sans doute nécessaire depuis les machines Windows.

La 2nd ligne ajoute les noms de champs. Cela de façon manuelle, à ce niveau vous pouvez les renommer à votre guise.

Les suite des opérations est plus classique, hormis la dernière ligne qui vient fermer le fichier.

Écriture de Shape

Nous allons maintenant créer un shapefile. Pour que ce soit drôle, et ne pas créer un shape de shape, nous allons plutôt exporter une sélection des aiguilles de notre shape des sommets.

Cette fois inutile de sélectionner la couche dans le panneau des calques, nous utiliserons la méthode mapLayersByName qui ira chercher la couche concernée par son nom (1ère ligne).

  1. layer = QgsProject.instance().mapLayersByName('peaks')[0]
  2. layer.selectByExpression('"NAME" LIKE \'Aiguille%\'')
  3.  
  4. root = r'C:/Users/Georges/Downloads/Aiguilles/'
  5. if not os.path.exists(root):
  6.     os.makedirs(root)
  7.  
  8. file = str(root)+'Aiguilles.shp'
  9.  
  10. writer = QgsVectorFileWriter.writeAsVectorFormat\
  11. (layer, file, 'utf-8', driverName='ESRI Shapefile', onlySelected=True)

La 2nd ligne sélectionne nos aiguilles à la manière du langage SQL.

Des lignes 4 à 6 nous créons le dossier qui va contenir le shape. Avant cela nous vérifions qu'il n'existe pas déjà. De plus nous stockons son chemin complet dans une variable dédiée, car elle nous servira juste après à la génération du shape lui-même. Nous aurions pu procéder plus simplement, mais cette partie du code est optimisée pour éviter la redondance. Nous aurions même pu faire plus court en utilisant un chemin relatif, nous y reviendrons plus tard.

Retenez également le paramètre onlySelected=True dans la dernière ligne, qui vous rappelle sans doute des choses.

Observez le slash inversé (\) qui clôt la 10ème ligne. Puisque les sauts de ligne sont interdits dans le code Python, certaines lignes peuvent être atrocement longue, en particulier les chaînes de texte par exemple. Le slah inversé permet de signifier à Python un saut de ligne sans impact dans le code.

Export d'images

Vous pouvez très rapidement exporter l'image du canvas courant, en format png par exemple :

iface.mapCanvas().saveAsImage("C:/Users/Georges/Downloads/SommetsEcrins.png")

Remarquez que Pyhon-QGIS nous exporte également un fichier World File, correspondant notamment à l'étendue géographique de l'image exportée. Si cela ne vous intéresse pas, supprimez-le à la suite de votre export. Comme ici en format jpg cette fois :

iface.mapCanvas().saveAsImage("C:/Users/Georges/Downloads/SommetsEcrins.jpg")
os.remove("../SommetsEcrins.jgw")

Remarquez l'usage d'un chemin relatif dans la suppression du fichier World File ci-dessus. Mon répertoire courant étant situé dans /Download, pour le remonter une fois j'utilise ../ qui me permet de viser directement le répertoire Download sans mentionner le chemin complet.

Une bonne maîtrise de la syntaxe de vos chemins vous permettra de raccourcir votre code.

Répertoire courant

En cas de doute vous pouvez connaître facilement votre répertoire courant en lançant le code ci-dessous.

path = os.getcwd()
print("Le répertoire courant est : " + path)

Changer de répertoire courant

path = os.getcwd()
print("Le répertoire courant est :\n" + path)
 
os.chdir('C:/Users/Georges/Downloads/')
 
path = os.getcwd()
print("Le répertoire courant a changé et c'est maintenant :\n" + path)

Export de PDF

Pour exporter un PDF il faut d'abord avoir créé un layout, une carte. Créez-en un via le Layouts Manager de QGIS. Dans le code ci-dessous mon layout s'appelle Map1.

Pour connaître les layouts disponible dans un projet QGIS, une petite boucle fera l'affaire :

myLayouts = QgsProject.instance().layoutManager()
for layout in myLayouts.printLayouts():
    print(layout.name())

Pour exporter une carte en PDF :

manager = QgsProject.instance().layoutManager()
layout = manager.layoutByName("Map1")
exporter = QgsLayoutExporter(layout)
exporter.exportToPdf("C:/Users/Georges/Downloads/MaJolieCarte1.pdf",QgsLayoutExporter.PdfExportSettings())

Sélection et fonction

Vous verrez que le mot fonction n'est pas anodin en programmation Python ! En revanche quand nous parlons de sélection, nous parlons d'une simple sélection géographique/attributaire.

Sélection sans fonction

Ouvrez les couches data_BDTOPO_V3_Dep05_adresse et chemin_de_fer (jointes en bas de cet article) dans un projet QGIS, ainsi que la console Python du logiciel. Pour plus de clarté, mettez la couche de points en 1er plan et la couche linéaire en 2nd plan.

Exécutez le code suivant, que nous allons ensuite décrypter :

  1. mylayer = QgsProject.instance().mapLayersByName("ADRESSE_05")[0]
  2. myselect = mylayer.getFeatures( QgsFeatureRequest().setFilterExpression ( u'"CODE_POST" = \'05000\'' ) )
  3. mylayer.selectByIds( [ f.id() for f in myselect ] )
  4. iface.mapCanvas().zoomToSelected(mylayer)

La 1ère ligne crée une variable mylayer, qui regarde le projet QGIS courant, notre instance, et y identifie la couche des points ADRESSE_05.

Le 2nd ligne créé une variable myselect, qui à l'intérieur de la 1ère variable contenant la couche spatiale mylayer, identifie des entités de cette couche en fonction d'une expression (CODE_POST = 05000, pour simplifier. Nous y reviendrons). Ce sont les adresses de la ville de Gap et des communes alentours (code postal 05000).

La 3ème ligne sélectionne les entités concernées avec une boucle for. La 4ème, pour notre confort personnel, zoome sur ces entités.

Nous reviendrons plus tard sur la notion de paramètre, mais remarquez que la méthode zoomToSelected prend déjà un paramètre (mylayer), qui est lui même notre 1ère variable. Situé juste après la méthode et entre parenthèse, il mentionne que la méthode doit s'appliquer sur cette couche spatiale.

Les paramètres des méthodes, ou des fonctions, ne sont pas toujours obligatoires (essayez de le supprimer par exemple, en laissant les parenthèses vides, ré-exécutez le code). Mais dans la mesure où nous avons ici plusieurs couches d'ouvertes dans QGIS, il est nécessaire de préciser sur quelle couche le logiciel doit zoomer. Vous avez ici un 1er exemple du caractère dynamique du langage Python (un autre paramètre, dans cette même fonction, aurait agit ailleurs. Nous reviendrons bientôt sur ce point).

Sélection avec fonction mais sans paramètre

Désélectionnez puis dézoomez. Nous allons maintenant faire la même chose mais en passant par une fonction. Comme son nom l'indique, cela devrait nous permettre d'avoir à notre disposition une vraie fonction, un outil en somme, que nous pourrons appeler et utiliser à notre guise. Exécutez le code suivant :

def SelectAndZoom():
    mylayer = QgsProject.instance().mapLayersByName("ADRESSE_05")[0]
    myselect = mylayer.getFeatures( QgsFeatureRequest().setFilterExpression ( u'"CODE_POST" = \'05000\'' ) )
    mylayer.selectByIds( [ f.id() for f in myselect ] )
    iface.mapCanvas().zoomToSelected(mylayer)

Absolument rien ne s'est produit, c'est tout-à-fait normal.

Observez l'indentation. De toute évidence quelque chose est contenu à l'intérieur de ce bout de code.

Nous avons ici créé une fonction arbitrairement nommée SelectAndZoom (vous commencez à comprendre l'intérêt des fonctions...). Nous savons qu'il s'agit d'une fonction car elle est précédée d'un def, qui définit la présence d'une fonction. Nous l'avons créée, nous l'avons même enregistrée dans le cache de notre instance QGIS, mais nous ne l'avons pas exécutée.

Elle est tout de même là, prête à bondir à chaque instant !

Maintenant appelez-là :

SelectAndZoom()

Astuce
Avec les jeu des copier-coller depuis ce site web ou un éditeur de texte gérant différemment les indentations - séparateur tab, whitespace et saut de ligne - il peut y avoir des erreurs d'interprêtation/compilation de Python à l'appel de vos fonctions.

Dans ce cas utilisez plutôt l'éditeur Python de QGIS (survolez les boutons disponibles en haut de la console Python de QGIS, l'éditeur s'ouvrira à droite de la console).

Hop, nous voici de retour à Gap ! Amusez-vous à dézoomer, déselectionner, vous déplacer dans votre projet QGIS ou même ouvrir des tables attributaires... Puis ré-exécutez la fonction en l'appelant simplement (SelectAndZoom()). Elle est toujours là, et opérationnelle, agissant tant dans votre canvas que dans vos données attributaires.

Sélection avec fonction et un paramètre

Maintenant nous n'allons même plus prendre la peine de mentionner la couche cible dans la fonction.

À la place, nous mentionnerons la nécessité d'un paramètre dans le nom même de la fonction, entre parenthèse et arbitrairement nommé layername, et nous remplaçerons le nom de la couche cible par ce même paramètre :

def SelectAndZoom(layername):
    mylayer = QgsProject.instance().mapLayersByName(layername)[0]
    myselect = mylayer.getFeatures( QgsFeatureRequest().setFilterExpression ( u'"CODE_POST" = \'05000\'' ) )
    mylayer.selectByIds( [ f.id() for f in myselect ] )
    iface.mapCanvas().zoomToSelected(mylayer)

La nouvelle fonction du même nom replace la 1ère, et l'appel se fait maintenant avec un argument, le nom de la couche cible :

SelectAndZoom('ADRESSE_05')

Sélection avec fonction et plusieurs paramètres

Maintenant nous n'allons pas nous contenter d'appeler la fonction et la couche cible, mais carrément la couche cible, le champ et la valeur voulue. Attention nous allons devoir jouer avec la notion de chaîne formatée :

def SelectAndZoom(layername, field, value):
    mylayer = QgsProject.instance().mapLayersByName(layername)[0]
    myselect = mylayer.getFeatures( QgsFeatureRequest().setFilterExpression ( u'"%s" = \'%s\'' % (field, value)) )
    mylayer.selectByIds( [ f.id() for f in myselect ] )
    iface.mapCanvas().zoomToSelected(mylayer)

L'appel requiert maintenant 3 arguments :

SelectAndZoom('ADRESSE_05', 'CODE_POST', '05100')

Chaînes formatées

Les chaînes formatées sont pénibles, hein ? Ainsi que les caractères d'échappement. C'est vrai.

u'"%s" = \'%s\'' % (field, value)

Mais qu'est-ce que c'est que ça ???

Pénible mais très performant, cela nous permet d'inclure du code hors le code. Ici Python attendait une chaîne. Mais nous avons voulu lui fournir une expression à la place, une expression qui aille chercher nos paramètres avant de les remplacer par des arguments (les chaînes finales). Il faut donc lui dire attention, ici ce n'est pas du simple texte, regarde d'abord les paramètres.

Il y a plusieurs façons de formater des chaînes avec Python. Ici nous utilisons la méthode %.

Le u précise que notre chaine sera en caractère Unicode.

Langage dynamique

Bien, vous utilisez maintenant pleinement le côté dynamique de Python : en appelant la même fonction mais avec des argument différents, vous opérerez des actions différentes. Essayez par exemple d'appeler une autre couche, avec un autre champ et une autre valeur :

SelectAndZoom('chemin_de_fer', 'NATURE', 'TGV')

Vous avez compris : un formulaire, autrement dit une interface graphique, pourrait fournir ces arguments suite à une saisie humaine. C'est comme ça que les logiciels et applications web fonctionnent !

  • Note personnelle : ici prévoir exposé sur les différents usages possibles de Python.

Récupération d'informations d'entités déjà sélectionnées

Récupérez, dézippez et ouvrez la couche decathlon_france.

Sélectionnez quelques magasins au hasard, puis exécutez la code suivant.

layer_selected = QgsProject.instance().mapLayersByName("decathlon_france")[0]
selected_features = layer_selected.selectedFeatures()
 
for i in selected_features:
    attrs = i.__geo_interface__
    id_mag = i['id']
    print (attrs)
    print (id_mag)

Hop ! La console nous affiche certaines informations isues de notre sélection.

Remarquez la méthode attrs, qui récupère tous les attributs.

 


Sélection intelligente

Téléchargez les couches decathlon_france et iso_iris. Ce sont des points et leurs zones isochrones constituées de portions d'IRIS à moins 15 minutes en voiture de ces points.

Mettez une transparence complète sur les zones isochrones.

Maintenant que vous savez faire des sélections et enregistrer des fonctions dans la mémoire Python, voyons un exemple concret d'usage de fonction dans un SIG.

Imaginez que vous ayez une couche contenant des entités se superposant. Il est malaisé de les examiner sans faire de sélection attributaire, jouer sur la symbologie, les styles ou l'ordre d'affichage... Créons donc une fonction afin de les faire resortir suite à un événement (ici, un clic sur une entité, ou une simple sélection manuelle de quelques entités).

Mais certaines de ces zones se superposent, et les IRIS, tout ou partie, existent autant de fois qu'il appartiennent à une zone. Le code suivant va sélectionner les IRIS d'une zone isochrone après sélection d'un magasin (et donc faire apparaître la zone isochrone concernée).

layer_selected = QgsProject.instance().mapLayersByName("decathlon_france")[0]
layer_to_select = QgsProject.instance().mapLayersByName("iso_iris")[0]
 
def SelectionAuto():
    selected_features = layer_selected.selectedFeatures()
    for i in selected_features:
        attrs = i.__geo_interface__
        id_mag = i['id']
        #print (id_mag)
 
        myselect = layer_to_select.getFeatures( QgsFeatureRequest().setFilterExpression ( u'"id_mag" = \'%s\'' % id_mag ) )
        layer_to_select.selectByIds( [ f.id() for f in myselect ] )
        #iface.mapCanvas().zoomToSelected(layer_to_select)
 
layer_selected.selectionChanged.connect(SelectionAuto)

Dans la suite de cet article, vous apprendrez à faire des symbologies en Python, vous pourrez donc enrichir ce type de fonction.

  • Note personnelle pour plus tard : enrichir la fonction d'une symbologie et d'un ordre d'affichage des entités sélectionnées. L'objectif étant de ne plus avoir à mettre de la transparence sur la couche pour rendre cette fonction utile.

Action enregistrée

Mais à ce stade votre sélection intelligente ne fonctionnera que tant que votre projet QGIS est ouvert. Fermez-le, et le code Python ne sera pas pris en compte à sa prochaine ouverture.

actionPour enregistrer ce code dans votre projet, il faudra l'enregistrer au niveau des Actions d'une couche. Allez dans les propriétés d'une de vos couches (n'importe laquelle pour ce qui nous concerne, mais celle des magasins est la plus à propos), onglet Actions, choisissez le type Python (😄), nommez votre action, décrivez-la brièvement et collez-y le précédent code. Cochez également l'option Canvas Scope afin de la faire apparaître dans le menu de QGIS.

Votre code est désormais activable dans votre projet QGIS sous l'icône Actions, par un clic sur celle-ci, suivi d'un autre sur la couche stockant l'action (avec la petite croix qui apparaîtra). Sélectionnez maintenant un magasin, toute sa zone isochrone sera aussi sélectionnée.

Les 3 lois d'Asimov

Cependant si les machines n'ont a priori pas de sentiment, elles sont tout-de-même capricieuses, et des événements inopportuns peuvent se produire.

En effet votre fonction n'a de sens qu'après sélection d'un seul magasin, d'ailleurs elle ne sélectionne finalement que les entités liées à un seul magasin. Alors que se passe-t-il en cas de sélection de plusieurs magasins ? De plus si jamais vous sélectionnez tous vos magasins en même temps (Ctrl+A), QGIS ne lancera-t-il pas un important calcul, avec risque de plantage ? Faites quelques tests pour comprendre le problème.

Il serait donc bel et bon d'ajouter un comptage, une éventuelle désélection, ainsi que des conditions (nos 3 lois robotiques ! 😄) pour que la sélection ne se déclenche que si un seul magasin a été sélectionné, et désélectionne d'éventuelles entités précédentes. Ceci afin de ne pas induire l'utilisateur en erreur.

Utilisons les méthodes selectedFeatureCount() et removeSelection(), dans des conditions if.

  1. layer_selected = QgsProject.instance().mapLayersByName("decathlon_france")[0]
  2. layer_to_select = QgsProject.instance().mapLayersByName("iso_iris")[0]
  3.  
  4. def SelectionAuto():
  5.     selected_features = layer_selected.selectedFeatures()
  6.     for i in selected_features:
  7.         if layer_selected.selectedFeatureCount() == 1:
  8.             id_mag = i['id']
  9.             myselect = layer_to_select.getFeatures( QgsFeatureRequest().setFilterExpression ( u'"id_mag" = \'%s\'' % id_mag ) )
  10.             layer_to_select.selectByIds( [ f.id() for f in myselect ] )
  11.         if layer_selected.selectedFeatureCount() > 1:
  12.             layer_to_select.removeSelection()
  13.  
  14. layer_selected.selectionChanged.connect(SelectionAuto)

Dans la ligne 7, une 1ère condition if vérifie que le nombre de magasins sélectionnés est égal à 1, cela grâce à la fonction selectedFeatureCount(), qui renvoie le nombre d'entités sélectionnées. Si c'est bien le cas, le reste du contenu de cette condition s'exécute.

À la ligne 11, un autre if se contente de désélectionner (removeSelection()) la zone éventuellement sélectionnée auparavant, à condition qu'il y ait plus d'1 magasin sélectionné.

  • Note personnelle pour plus tard : faire la même chose mais sans fonction Def, dans une action activable au clic dans la table attributaire ou via l'icône action. Ceci afin de ne pas être confronté à une sélection des zones en permanence (impossibilité de désactiver cette action Def).

Variante pour afficher plusieurs zones isochrones après sélection de plusieurs magasins

  1. layer_selected = QgsProject.instance().mapLayersByName("decathlon_france")[0]
  2. layer_to_select = QgsProject.instance().mapLayersByName("iso_iris")[0]
  3.  
  4. def SelectionAuto():
  5.     myList = []
  6.     selected_features = layer_selected.selectedFeatures()
  7.     for i in selected_features:
  8.         id_mag = i['id']
  9.         myList.append(id_mag)
  10.  
  11.     myList = str(myList).replace("[", "").replace("]", "")
  12.     print(myList)
  13.  
  14.     myselect = layer_to_select.getFeatures( QgsFeatureRequest().setFilterExpression ( u'"id_mag" IN (%s)' % myList ) )
  15.     layer_to_select.selectByIds( [ f.id() for f in myselect ] )
  16.  
  17. layer_selected.selectionChanged.connect(SelectionAuto)

 


Générer une carte

Nous allons maintenant générer une carte PDF à partir du canvas QGIS. À chaque déplacement dans le canvas, une nouvelle exécution du script final mettra à jour la carte. 

Environnement de travail

Ouvez les couches peaks, eau, sol, IRIS et troncons_routes, fournies en bas de cet article, dans un projet QGIS vierge, ainsi qu'un fond de carte OpenStreetMap.

Vérifiez la bonne supersposition des couches (les projections), puis appliquez une symbologie évidente.

Enregistrez votre projet puis décochez les couches eau, sol, IRIS et troncons_routes.

Ouvrez également le Layouts Manager. De même pour les besoins de la démonstration, faîtes-en sorte que vos écrans vous permettent de visualiser en même temps votre projet QGIS, console et éditeur Python ouverts, le Layouts Manager, et encore un peu de place pour un futur layout.

Zoomez sur la couche peak. Nous allons exécuter les bribes de code suivantes de façon itérative, en les décryptant une par une. Au terme de ce chapitre, vous aurez un code complet vous permettant de générer une carte simple. Sauvegardez petit-à-petit l'intégralité du code dans un fichier texte.

Layout

Nous générons un layout vide, que vous allez voir apparaître dans le Layouts Manager.

  1. project = QgsProject.instance()
  2. manager = project.layoutManager()
  3. layoutName = "PrintLayout1"
  4.  
  5. # Vérification de la non-existence d'un layout de même nom
  6. layouts_list = manager.printLayouts()
  7. for layout in layouts_list:
  8.     if layout.name() == layoutName:
  9.         manager.removeLayout(layout)
  10.  
  11. # Génération d'un layout vide
  12. layout = QgsPrintLayout(project)
  13. layout.initializeDefaults()
  14. layout.setName(layoutName)
  15.  
  16. manager.addLayout(layout)

Encart cartographique

Nous tirons maintenant une carte vide dans la layout.

  1. # Charger une carte vide
  2. map = QgsLayoutItemMap(layout)
  3. map.setRect(20, 20, 20, 20)
  4.  
  5. # Mettre un canvas basique
  6. rectangle = QgsRectangle(1355502, -46398, 1734534, 137094)
  7. map.setExtent(rectangle)
  8. layout.addLayoutItem(map)

Canvas courant

Nous plaçons le canvas courant dans la carte.

  1. # Mettre finalement le canvas courant
  2. canvas = iface.mapCanvas()
  3. map.setExtent(canvas.extent())
  4.  
  5. layout.addLayoutItem(map)
  6.  
  7. # Redimensionner la carte
  8. map.attemptMove(QgsLayoutPoint(5, 27, QgsUnitTypes.LayoutMillimeters))
  9. map.attemptResize(QgsLayoutSize(220, 178, QgsUnitTypes.LayoutMillimeters))

Légende dynamique

Ci-après le code d'affichage de la légende complète (une légende fixe, statique). Cette portion du code est commentée (encadré dans une balise """, et donc inactive).

En effet juste en dessous, nous allons préférer ne retenir que les couches apparaissant dans la carte (cochées dans le panneau des couches).

  1. """
  2. # Mettre une légende complète
  3. legend = QgsLayoutItemLegend(layout)
  4. legend.setTitle("Légende")
  5. layout.addLayoutItem(legend)
  6. legend.attemptMove(QgsLayoutPoint(246, 5, QgsUnitTypes.LayoutMillimeters))
  7. """
  8.  
  9. # Légende personnalisée
  10. tree_layers = project.layerTreeRoot().children()
  11. checked_layers = [layer.name() for layer in tree_layers if layer.isVisible()]
  12. print(f"Je vais ajouter uniquement {checked_layers} à la légende." )
  13.  
  14. layers_to_remove = [layer for layer in project.mapLayers().values() if layer.name() not in checked_layers]
  15. print(f"Je vais retirer {layers_to_remove} de la légende." )
  16.  
  17. legend = QgsLayoutItemLegend(layout)
  18. legend.setTitle("Légende")
  19.  
  20. layout.addLayoutItem(legend)
  21.  
  22. legend.attemptMove(QgsLayoutPoint(240, 24, QgsUnitTypes.LayoutMillimeters))
  23.  
  24. # Cette ligne permettra de ne pas sortir les couches inutilisées de votre panneau des calques QGIS, mais uniquement de la légende
  25. legend.setAutoUpdateModel(False)
  26.  
  27. m = legend.model()
  28. g = m.rootGroup()
  29. for l in layers_to_remove:
  30.     g.removeLayer(l)
  31.  
  32. legend.adjustBoxSize()

Jouez avec les print et la console python avec votre layout ouvert pour bien comprendre ce qui se passe.

Algorithmique

Pour la génération de notre légende dynamique ci-dessus, nous avons commencer par lister 3 variables.

  • Une contenant toutes les couches (tree_layers)
  • Une autre contenant uniquement les couches affichées (checked_layers)
  • Enfin une dernière ne contenant que les couches non-cochées (layers_to_remove), grâce à une 1ère soustraction des couches cochées

En algorithmique :

layers_to_remove = tree_layers - checked_layers
  • Ensuite nous affichons la légende complète puis en retirons layers_to_remove dans une 2nd soustraction.

Toujours en algorithmique :

my_legend = legend - layers_to_remove
  • Note personnelle : on pourrait aller plus loin dans cette légende dynamique en masquant de la légende les couches dont aucune entité n'apparaît dans le canvas courant (la couche adresse notamment).

Titre

  1. # Titre
  2. title = QgsLayoutItemLabel(layout)
  3. title.setText("Ma Jolie Carte")
  4. title.setFont(QFont("Verdana", 28))
  5. title.adjustSizeToText()
  6.  
  7. layout.addLayoutItem(title)
  8.  
  9. title.attemptMove(QgsLayoutPoint(5, 4, QgsUnitTypes.LayoutMillimeters))

Tout comme avec la légende, vous connaissez maintenant la méthode attemptMove qui vous permet de positioner les éléments dans votre layout.

Les chiffres en argument sont dans l'ordre : la distance à partir du bord gauche du layout, puis la distance à partir du bord haut.

Sous-titre

  1. # Sous-titre
  2. subtitle = QgsLayoutItemLabel(layout)
  3. subtitle.setText("Est-elle belle ?")
  4. subtitle.setFont(QFont("Verdana", 17))
  5. subtitle.adjustSizeToText()
  6.  
  7. layout.addLayoutItem(subtitle)
  8.  
  9. subtitle.attemptMove(QgsLayoutPoint(5, 19, QgsUnitTypes.LayoutMillimeters))

Échelle

  1. # Échelle
  2. scalebar = QgsLayoutItemScaleBar(layout)
  3. scalebar.setStyle('Single Box')
  4. scalebar.setUnits(QgsUnitTypes.DistanceKilometers)
  5. scalebar.setNumberOfSegments(2)
  6. scalebar.setNumberOfSegmentsLeft(0)
  7. scalebar.setUnitsPerSegment(5)
  8. scalebar.setLinkedMap(map)
  9. scalebar.setUnitLabel('km')
  10. scalebar.setFont(QFont('Verdana', 20))
  11. scalebar.update()
  12.  
  13. layout.addLayoutItem(scalebar)
  14.  
  15. scalebar.attemptMove(QgsLayoutPoint(10, 185, QgsUnitTypes.LayoutMillimeters))

Logo

  1. # Logo
  2. Logo = QgsLayoutItemPicture(layout)
  3. Logo.setPicturePath("https://master-geomatique.org/templates/theme3005/images/logo-ucp-cyu.png")
  4.  
  5. layout.addLayoutItem(Logo)
  6.  
  7. Logo.attemptResize(QgsLayoutSize(40, 15, QgsUnitTypes.LayoutMillimeters))
  8. Logo.attemptMove(QgsLayoutPoint(250,4, QgsUnitTypes.LayoutMillimeters))

Vous utilisez ici la méthode attemptResize pour redimensionner l'élément. Les chiffres en argument sont dans l'ordre la largeur, puis la hauteur.

Petit exercice

Maintenant que vous savez mettre une image, ajoutez donc une flèche nord correctement positionnée.

Vous pouvez faire autrement si vous le souhaitez.

Texte

  1. # Texte
  2. TextCustom = QgsLayoutItemLabel(layout)
  3. TextCustom.setText("Le massif des Écrins est un grand massif montagneux des Alpes françaises situé dans les Hautes-Alpes et en Isère.\
  4. \n\nIl abrite d'importants glaciers, tant en nombre qu'en taille et possède deux sommets de plus de 4 000 mètres.\
  5. \n\nIl était autrefois également nommé massif du Pelvoux.\
  6. \n\nL'Oisans (bassin de la Romanche) au nord-ouest, le Champsaur (haut-bassin du Drac) au sud-ouest, et le Briançonnais (bassin de la Guisane) au nord-est recouvrent une partie du massif.")
  7. TextCustom.setFont(QFont("Verdana", 11))
  8.  
  9. layout.addLayoutItem(TextCustom)
  10.  
  11. TextCustom.attemptResize(QgsLayoutSize(60, 120, QgsUnitTypes.LayoutMillimeters))
  12. TextCustom.attemptMove(QgsLayoutPoint(230, 80, QgsUnitTypes.LayoutMillimeters))

Petit exercice

Maintenant que vous savez mettre un texte (un simple label en réalité, comme les titre et sous-titre), ajoutez la source.

  • Note personnelle : on pourrait aller plus loin en allant chercher les sources directement dans les méta-données de nos couches.

Export

Vous connaissez déjà la manœuvre pour un export. Pensez à adapter le chemin utilisé ci-dessous.

  1. # Export PDF
  2. manager = QgsProject.instance().layoutManager()
  3. layout = manager.layoutByName("PrintLayout1")
  4. exporter = QgsLayoutExporter(layout)
  5. exporter.exportToPdf("C:/Users/Georges/Downloads/Ma Belle Carte.pdf",QgsLayoutExporter.PdfExportSettings())

En ouvrant le PDF dans votre navigateur, en vous déplaçant sur la carte, en cochant-décochant certaines couches puis en re-lançant l'intégralité du code ci-dessus, vous verrez votre PDF exporté être mis à jour en fonction.

  • Note personnelle : on pourrait aller plus loin en ajoutant une vraie flèche nord, dynamique, suivant l'orientation du canvas.

Exercice

Maintenant vous savez :

  • Faire une boucle
  • Faire une sélection
  • Exporter une carte

Et bien que se passerait-il si vous mettiez la sélection et l'export dans une boucle ? La production de cartes en série, pardi.

À vous de jouer !


Générer des cartes

Nous allons maintenant générer 341 cartes en PDF, chacune zoomant sur nos 341 sommets.

Il s'agit en quelque sorte d'une reproduction de la fonction Atlas de QGIS, mais nous allons pouvoir aller beaucoup plus loin : appliquer des zooms différents par cartes, des symbologies conditionnelles, ajouter des contenus web, lancer des géo-traitements...

* Pour les plus pressés (🧐) un exemple de code complet est disponible ici.

Sélection en aveugle

Précédemment nous sélectionnions notre sommet, star de notre carte, sur la base de son identifiant. Mais maintenant que nous souhaitons itérer sur toute la couche peaks, nous ne pouvons plus nommément utiliser les identifiants (en fait si : nous pourrions les écrire un par un dans le code, dans un IN, mais ce n'est pas le but du jeu, une entité ajoutée dans les données ne serait pas prise en compte, et une entité supprimée ferait bugger le code).

Il nous faut donc boucler sur les entités sommitales.

Inspirez vous de la Sélection sans fonction plus haut dans ce court, pour créer une sélection d'un sommet sans connaître son identifiant.

Votre sélection des sommets se trouvera donc dans une boucle. Nous ne pourrons donc plus zoomer ou nous déplacer manuellement sur la couche pour générer les canvas de nos cartes.

Pour plus de clarté commencez par mettre votre expression dans une variable à part (5ème ligne ci-dessous). Utilisez des print pour visualiser progressivement le contenu de vos variables.

mylayer = QgsProject.instance().mapLayersByName("peaks")[0]
 
for feat in mylayer.getFeatures():
    id_peak = feat['OSM_ID']
    expr = u"OSM_ID = '{}'".format(id_peak)
    myselect = mylayer.getFeatures( QgsFeatureRequest().setFilterExpression ( expr ))
 
    mylayer.selectByIds( [ f.id() for f in myselect ] )
 
    iface.mapCanvas().zoomToSelected(mylayer)

Bien, vous êtes capable de lister les différents sommets à l'intérieur d'une boucle. 

Profitons-en pour ajouter une variable contenant leurs noms (champ NAME) et affichons là dans un print.

Contrôlons également le niveau de zoom avec la méthode zoomScale().

  1. for feat in mylayer.getFeatures():
  2.     id_peak = feat['OSM_ID']
  3.     peak_name = feat['NAME']
  4.     expr = u"OSM_ID = '{}'".format(id_peak)
  5.     myselect = mylayer.getFeatures( QgsFeatureRequest().setFilterExpression ( expr ))
  6.     print(f"Voici {peak_name}" )
  7.  
  8.     mylayer.selectByIds( [ f.id() for f in myselect ] )
  9.  
  10.     iface.mapCanvas().zoomToSelected(mylayer)
  11.     iface.mapCanvas().zoomScale(23000)

Notez que le code ci-dessous se termine en sélectionnant la dernière entité listée, mais il a bel et bien listé, sélectionné puis zoomé une par une sur chaque entité. Même si cela ne s'est pas vu !

QGIS en a même gardé une trace dans son cache : si vous jouez avec le bouton Zoom Précédent du menu de QGIS, puis examinez son nom en comparant avec votre liste print, vous vous en apercevrez.

Utiliser les boutons QGIS dans le code Python

Petite parenthèse : notez qu'il existe une façon rigolote de zoomer sur les entités sélectionnées, en utilisant la liste des boutons QGIS :

  1. for feat in mylayer.getFeatures():
  2.     id_peak = feat['OSM_ID']
  3.     peak_name = feat['NAME']
  4.     expr = u"OSM_ID = '{}'".format(id_peak)
  5.     myselect = mylayer.getFeatures( QgsFeatureRequest().setFilterExpression ( expr ))
  6.  
  7.     mylayer.selectByIds( [ f.id() for f in myselect ] )
  8.  
  9.     eMenu = iface.viewMenu()
  10.     eMenu.actions()[12].trigger()
  11.  
  12.     iface.mapCanvas().zoomScale(23000)

Ci-dessus notre variable eMenu contient le menu QGIS, et la méthode actions()[12] exécute le douzième bouton (dans mon cas c'est le bouton Zoom sur la sélection).

Utiliser les données dans la carte

Maintenant vous allez adaptez dans cette boucle le code complet du chapitre Générer une carte, afin de nommer votre layout, votre carte, lui mettre un titre, l'id OSM en guise de sous-titre, et un export correctement nommé également.

Vous allez devoir sélectionner des blocs de code et jouer avec le bouton Tab (indentation) et l'alliance Shift+Tab (retrait d'indentation) afin d'indenter correctement votre code dans la nouvelle boucle.

Testez votre code de façon itérative et procédurale :

  • D'abord la création des layouts, allez ensuite vérifier qu'ils existent bien.
  • Ensuite le zoom de chaque layout sur chaque sommet, vérifiez-en quelques uns.
  • ...

Layouts, canvas et légende

Il vous faudra passer votre variable peak_name en chaîne de texte pour éviter des erreurs d'interprêtation de la variable layoutName. Utilisez la méthode str().

N'oubliez pas non plus de tout mettre dans la boucle ! Soit une indentation à ajouter au code ci-dessous.

...
layoutName = str(peak_name)
 
project = QgsProject.instance()
manager = project.layoutManager()
 
# Vérification de la non-existence d'un layout de même nom
layouts_list = manager.printLayouts()
for layout in layouts_list:
    if layout.name() == layoutName:
        manager.removeLayout(layout)
 
# Génération d'un layout vide
layout = QgsPrintLayout(project)
layout.initializeDefaults()
layout.setName(layoutName)
manager.addLayout(layout)
 
# Charger une carte vide
map = QgsLayoutItemMap(layout)
map.setRect(20, 20, 20, 20)
 
# Mettre un canvas basique
rectangle = QgsRectangle(1355502, -46398, 1734534, 137094)
map.setExtent(rectangle)
layout.addLayoutItem(map)
 
# Mettre finalement le canvas courant
canvas = iface.mapCanvas()
map.setExtent(canvas.extent())
layout.addLayoutItem(map)
 
# Redimensionner la carte
map.attemptMove(QgsLayoutPoint(5, 27, QgsUnitTypes.LayoutMillimeters))
map.attemptResize(QgsLayoutSize(220, 178, QgsUnitTypes.LayoutMillimeters))
 
# Légende personnalisée
tree_layers = project.layerTreeRoot().children()
checked_layers = [layer.name() for layer in tree_layers if layer.isVisible()]
 
# print(f"Je vais ajouter uniquement {checked_layers} à la légende." )
 
layers_to_remove = [layer for layer in project.mapLayers().values() if layer.name() not in checked_layers]
 
# print(f"Je vais retirer {layers_to_remove} de la légende." )
 
legend = QgsLayoutItemLegend(layout)
legend.setTitle("Légende")
layout.addLayoutItem(legend)
legend.attemptMove(QgsLayoutPoint(240, 24, QgsUnitTypes.LayoutMillimeters))
 
# Cette ligne permettra de ne pas sortir les couches inutilisées de votre panneau des calques QGIS, mais uniquement de la légende
legend.setAutoUpdateModel(False)
 
m = legend.model()
g = m.rootGroup()
for l in layers_to_remove:
    g.removeLayer(l)
 
legend.adjustBoxSize()
...

Titres et sous-titres

La variable peak_name contenant le nom des sommets, à utiliser dans le titre de nos cartes, à déjà été passée en texte, avec la méthode str(). Elle peut donc être utilisée telle quelle.

Mais ce n'est pas le cas de l'identifiant OSM, que nous voulons ajouter dans le sous-titre. Nous allons le faire mais en y laissant un peu de texte, grâce à une chaîne formatée.

  1. ...
  2. # Titre
  3. title = QgsLayoutItemLabel(layout)
  4. title.setText(layoutName)
  5. title.setFont(QFont("Verdana", 28))
  6. title.adjustSizeToText()
  7. layout.addLayoutItem(title)
  8. title.attemptMove(QgsLayoutPoint(5, 4, QgsUnitTypes.LayoutMillimeters))
  9.  
  10. # Sous-titre
  11. subtitle = QgsLayoutItemLabel(layout)
  12. subtitle.setText("Identifiant OSM : %s" % (str(id_peak)))
  13. subtitle.setFont(QFont("Verdana", 17))
  14. subtitle.adjustSizeToText()
  15. layout.addLayoutItem(subtitle)
  16. subtitle.attemptMove(QgsLayoutPoint(5, 19, QgsUnitTypes.LayoutMillimeters))
  17. ...

Échelle

L'échelle doit être légèrement adaptée, puisque nous avons modifié le niveau de zoom :

...
scalebar.setUnitsPerSegment(1)
...

Logo, légende et texte

Le logo ne change pas. La legende non plus pour l'instant (plus tard nous appliquerons une symbologie dynamique qui impactera la légende).

Le texte aussi fera l'objet d'un chapite dédié, ne le mettez pas pour l'instant.

Export

Là encore, chaîne formatée pour correctement nommer nos carte grâce à la variable layoutName :

  1. ...
  2. # Export PDF
  3. manager = QgsProject.instance().layoutManager()
  4. layout = manager.layoutByName(layoutName)
  5. exporter = QgsLayoutExporter(layout)
  6. exporter.exportToPdf("C:/Users/Georges/Downloads/qgis/scripts/exports/%s.pdf" % (layoutName),QgsLayoutExporter.PdfExportSettings())

Date du jour

Nous allons ajouter la date du jour dans les noms de nos fichiers.

Ouvrez la console de votre machine (CMD si vous êtes sous Windows) et entrez simplement:

python

La console passe en mode Python et affiche quelques informations de version. Maintenant ce code pour récupérer la date du jour :

  1. import datetime
  2. date = datetime.date.today()
  3. str(date)

Le même code va fonctionner dans la console Python de QGIS. La 1ère ligne est une importation d'un module natif, à inclure en haut de votre fichier.

Le date peut être stockée dans la boucle, nous la stockerons dans une variable today.

import datetime
mylayer = QgsProject.instance().mapLayersByName("peaks")[0]
 
for feat in mylayer.getFeatures():
    id_peak = feat['OSM_ID']
    peak_name = feat['NAME']
 
    date = datetime.date.today()
    today = str(date)
...

Avant de l'inclure dans le nom de nos exports dans la chaîne formatée.

...
exporter.exportToPdf\
("C:/Users/Georges/Downloads/qgis/scripts/exports/%s-%s.pdf"\
% (today, layoutName),QgsLayoutExporter.PdfExportSettings())

Bien, mais nous pouvons encore aller plus loin.


Symbologie dynamique

Symbologie simple sur des points

Pour appliquer une symbologie simple sur nos points, les sommets enregistrés dans la variable my_layer :

  1. my_layer = QgsProject.instance().mapLayersByName("Sommets")[0]
  2.  
  3. symbol_peak = QgsMarkerSymbol.createSimple({'name': 'Triangle', 'color': 'green', 'outline_color': 'black', 'size': '4'})
  4. mylayer.renderer().setSymbol(symbol_peak)
  5. mylayer.triggerRepaint()

Symbologie simple sur des lignes

Cette fois on précise la couleur en HTML.

  1. myroads = QgsProject.instance().mapLayersByName("Routes")[0]
  2. myrivers = QgsProject.instance().mapLayersByName("Fleuves")[0]
  3.  
  4. # Symbologie simple de la couche route
  5. symbol = QgsLineSymbol.createSimple({'line_style': 'dash', 'line_width': '0.5', 'color': 'black'})
  6. myroads.renderer().setSymbol(symbol)
  7. myroads.triggerRepaint()
  8.  
  9. # Symbologie simple de la couche eau
  10. symbol = QgsLineSymbol.createSimple({'line_style': 'solid', 'line_width': '0.75', 'color': '#0088CC'})
  11. myrivers.renderer().setSymbol(symbol)
  12. myrivers.triggerRepaint()

Symbologie simple sur des polygones

Cette fois on précise la couleur en RGB (50,165,77,75), afin d'ajouter un peu de transparence (ici 75, l'indice allant de 0 à 255), et on rafraîchit le panneau des couches afin que notre légende y soit aussi prise en compte,  :

  1. symbol = QgsFillSymbol.createSimple({'line_style': 'solid', 'line_width': '0.2', 'color': '50,165,77,75'})
  2. myground.renderer().setSymbol(symbol)
  3. myground.triggerRepaint()
  4. iface.layerTreeView().refreshLayerSymbology(myground.id())

Il est aussi possible d'examiner les propriétés de symbologie d'un layer :

print(my_layer.renderer().symbol().symbolLayers()[0].properties())

Symbologie catégorisée

Pour une symbologie catégorisée, il va nous falloir définir plusieurs symbologies, puis les attribuer en fonction de valeurs attributaires.

Ici nous mettons tous nos sommets en triangle bleu (▲, puis nous allons chercher un sommet unique en se basant sur son nom (Aiguille de la Gandolière, dans le champ NAME). Celui-ci passera en triangle bleu foncé et de plus grande taille (▲).

  1. symbol_peak = QgsMarkerSymbol.createSimple({'name': 'Triangle', 'color': '#0088CC', 'outline_color': 'black', 'size': '4'})
  2. symbol_peak_selected = QgsMarkerSymbol.createSimple({'name': 'Triangle', 'color': 'blue', 'outline_color': 'black', 'size': '7'})
  3.  
  4. color_peak = QgsRendererCategory(None,symbol_peak,"Sommet",True)
  5. color_peak_selected = QgsRendererCategory("Aiguille de la Gandolière",symbol_peak_selected,"Aiguille de la Gandolière",True)
  6.  
  7. renderer = QgsCategorizedSymbolRenderer("NAME", [color_peak,color_peak_selected])
  8.  
  9. mylayer.setRenderer(renderer)
  10. mylayer.triggerRepaint()

Les 2 premières lignes définissent les symbologies que nous utiliserons, puis dans les lignes 4 et 5, la méthode QgsRendererCategory prend 4 arguments, dans l'ordre :

  • La valeur recherchée dans le champ pour l'attribution de la symbologie
  • La symbologie
  • Le label à afficher dans la légende pour cette catégorie
  • Un booléen, pour l'application ou non de cette symbologie

Puis c'est la ligne 7 qui vient exécuter les conditions, avec 2 arguments :

  • Le nom de champ où rechercher les valeurs
  • Un tableau contenant les conditions à appliquer

Symbologie dynamique

Nous n'avons plus qu'à utiliser une variable pour gérer la symbologie dans la boucle de génération de nos cartes.

La variable contenant le nom du sommet sélectionné existe déjà dans notre code, et a déjà été passé en texte (layoutName). C'est que nous utiliserons, à la fois dans l'argument value (pour filtrer le champ NAME sur cette valeur) et dans l'argument label (pour l'afficher dans la légende).

  1. ...
  2. for feat in mylayer.getFeatures():
  3. ...
  4.     peak_name = feat['NAME']
  5. ...
  6.     layoutName = str(peak_name)
  7. ...
  8.     # Symbologie catégorisée dynamique
  9.     symbol_peak = QgsMarkerSymbol.createSimple({'name': 'Triangle', 'color': '#0088CC', 'outline_color': 'black', 'size': '4'})
  10.     symbol_peak_selected = QgsMarkerSymbol.createSimple({'name': 'Triangle', 'color': 'blue', 'outline_color': 'black', 'size': '7'})
  11.  
  12.     color_peak = QgsRendererCategory(None,symbol_peak,"Autres sommets",True)
  13.     color_peak_selected = QgsRendererCategory(layoutName,symbol_peak_selected,layoutName,True)
  14.  
  15.     renderer = QgsCategorizedSymbolRenderer("NAME", [color_peak,color_peak_selected])
  16.  
  17.     mylayer.setRenderer(renderer)
  18.     mylayer.triggerRepaint()
  19.  
  20.     # Charger une carte vide
  21.     map = QgsLayoutItemMap(layout)
  22.     map.setRect(20, 20, 20, 20)
  23. ...

Logique de programmation

Pour fonctionner, la symbologie dynamique doit donc se trouver :

  • Après la création de la variable peak_name (sinon la symbologie dynamique ne saura pas quel sommet prendre, puisque peak_name devient ensuite layoutName, variable utilisée dans l'une des conditions de notre symbologie)
  • Après la création de la variable layoutName (sinon elle ne saura pas l'utiliser en tant que chaîne de texte)
  • Avant la mise à jour du canvas (ce point n'est pas réellement nécessaire, mais il correspond à une logique de programmation et de lecture/compréhension de notre code)
  • Avant la mise à jour de la légende (puisque la symbologie apparaît dans la legénde. Là encore il n'est pas certain que cela soit nécessaire, mais ce n'est que pure logique itérative, procédurale)

 

  • Note personnelle : ici prévoir exposé sur les grandes logiques de programmation existantes.

Votre code doit pouvoir être aisément lu et compris par vous-même ou quelqu'un d'autre. Un jour vous souhaiterez peut-être le ré-utiliser et l'enrichir. Respecter une logique de prommation et bien commenter votre code n'est pas anodin.

  • Note personnelle : gérer les sauts de ligne et la taille de la légende quand le nom du sommet selectionné est trop long (déborde de la carte).

OK, mais nous pouvons aller encore un peu plus loin dans les contenus de nos cartes.

 


Contenus web

Vous vous rappelez l'encart de texte au sujet du Massif des Écrins, accompagnant la carte que nous avons générée plus haut ? En dur dans le code, nous l'avons volontairement mis de côté au moment de générer nos cartes en série. En effet afficher toujours le même texte dans nos 341 cartes n'avait pas d'intérêt.

Et bien nous allons maintenant enrichir nos cartes avec des contenus dynamiques, en lien avec les données cartographiées, qui iront remplir cet encart de texte. Et nous irons chercher ces contenus sur le web.

Plusieurs méthodes existent, APIs et web-scraping notamment, souvent avec l'utilisation de fichier JSON ou XML.

  • Note personnelle : ici prévoir exposé sur les APIs

Il nous faudra donc faire un lien avec nos données, et ce de façon dynamique bien sûr.

Dans ce chapitre, nous ferons nos tests directement depuis la console Windows, ouvrez donc la console Python dans l'interface de commande Windows (CMD + python).

 


Installer et utiliser PIP

Petite parenthèse : Pip est un installateur d'extensions pour Python, vous allez sans doute devoir vous y intéresser, et lui-même l'installer pour la suite de ce tutoriel. Je vous laisse faire.

Module requests

Si besoin d'installer requests (dans le shell Windows, pas Python) :

pip install requests

 


Exemples de fichiers XML

OSM

Exemple de lien vers le fichier XML de la base de données OSM, suffixé par l'id OSM (que justement, nous possèdons pour tous nos sommets) :

import requests
 
MyOsmId = '582066938'
r = requests.get("https://www.openstreetmap.org/api/0.6/node/"+MyOsmId)
 
print(r.text)

Notez l'importation nécessaire du module en haut de fichier.

Clé-valeur

Nos sommets n'ont pas de donnée attributaire sur leur altitude. Mais il semble que les XML d'OSM la connaissent ! En effet le fichier contient un tableau associatif avec un ensemble de clé-valeur dont ele (elevation en anglais).

Nous allons parser le fichier XML afin de récupérer l'altitude dans une jolie variable.

Parsing

Grâce à une condition if, nous allons récupérer l'altitude, à condition que le fichier XLM contienne le tag ele, sinon nous n'affichons rien.

import requests
import xml.etree.ElementTree
 
MyOsmId = '582066938'
r = requests.get("https://www.openstreetmap.org/api/0.6/node/"+MyOsmId)
root = xml.etree.ElementTree.fromstring(r.content)
 
for child in root.iter('tag'):
    if child.attrib['k'] == 'ele':
        print (child.attrib['v']+" mètres")

Nous sommes donc capable de récupérer l'altitude, bien, maintenant utilisons-là dans le titre de nos cartes générées :

  1. import datetime
  2. import requests
  3. import xml.etree.ElementTree
  4.  
  5. mylayer = QgsProject.instance().mapLayersByName("peaks")[0]
  6. ...
  7. for feat in mylayer.getFeatures():
  8. ...
  9.     # Récupérer altitude sur API OSM
  10.     r = requests.get("https://www.openstreetmap.org/api/0.6/node/"+id_peak)
  11.     root = xml.etree.ElementTree.fromstring(r.content)
  12.  
  13.     for child in root.iter('tag'):
  14.         if child.attrib['k'] == 'ele':
  15.             altitude = ", "+child.attrib['v']+" mètres"
  16.  
  17.     # Titre
  18.     title = QgsLayoutItemLabel(layout)
  19.     title.setText(layoutName+altitude)
  20.     title.setFont(QFont("Verdana", 28))
  21.     title.adjustSizeToText()
  22.     layout.addLayoutItem(title)
  23.     title.attemptMove(QgsLayoutPoint(5, 4, QgsUnitTypes.LayoutMillimeters))
  24. ...

Précédé d'une virgule, notre variable altitude n'a plus qu'à aller se concaténer à la suite de notre variable layoutName dans le titre de nos cartes. Ainsi les enregistrements pour lesquels nous ne connaissons pas l'atitude n'afficheront pas de caractère inopportun à la suite du titre.

  • Note personnelle : on pourrait aller plus loin en ajoutant dans les sources des cartes le nom du contributeur OSM ayant saisi ou modifié l'entité sommitale cible.

Encodage

Il serait bon d'encoder notre caractère spécial è dans le mot mètres, car ici il est affiché en tant que pure chaîne, et dans certains cas spéciaux nous pourrions avoir de mauvaises suprises.

import html
...
            altitude = ", "+child.attrib['v']+html.unescape(" m&amp;egrave;tres")
...
  • Note personnelle : faire aussi l'encodage du fichier lui-même, et parler de l'encodage des fichiers Python.

Wikidata

Certains sommets possèdent un id Wikidata (dans le champ OTHER_TAGS). Exemple de lien vers donnée Wikidata, suffixé par un id Wikidata :

Pour tester, allez chercher le contenu d'une page page HTML d'osm.wikidata.link :

  1. import requests
  2. r = requests.get("https://osm.wikidata.link/Q726652")
  3. print(r.text)

Un peu lent mais bien fourni. Les page wikidata.org se prêtent tout de même mieux au web-scraping (les contenus y sont mieux organisés, à l'intérieur de jolies balises).

Module requests-html

Si besoin, installez requests-html (dans le shell Windows, pas Python) :

pip install requests-html
Astuce
Pour une installation à l'université de Cergy (Python 2) :
py -m pip install requests-html

Récupérer des données Wikidata

Les contenus y sont si bien organisé que nous allons pouvoir récupérer le lien vers la page Wikipédia du sommet, et même cibler celle en français.

Dans la console Python de Windows :

  1. import requests
  2. from requests_html import *
  3.  
  4. session = HTMLSession()
  5.  
  6. MyWikiId = 'Q726652'
  7. r = session.get("https://www.wikidata.org/wiki/"+MyWikiId)
  8.  
  9. my_content = r.html.find('.wikibase-sitelinkview-link-frwiki a', first=True)
  10. #print(my_content)
  11.  
  12. my_link = my_content.xpath('//a/@href')
  13. print(my_link)

 Nous nettoyons ensuite le lien pour le rendre utilisable par le module requests :

  1. import requests
  2. from requests_html import *
  3.  
  4. session = HTMLSession()
  5.  
  6. MyWikiId = 'Q726652'
  7. r = session.get("https://www.wikidata.org/wiki/"+MyWikiId)
  8. my_content = r.html.find('.wikibase-sitelinkview-link-frwiki a', first=True)
  9.  
  10. #print(my_content)
  11. my_link = my_content.xpath('//a/@href')
  12.  
  13. #print(my_link)
  14.  
  15. my_link_string = str(my_link).replace('[', "").replace(']', "").replace("'", "")
  16. print(my_link_string)
  17.  
  18. r2 = session.get(my_link_string)
  19.  
  20. print(r2.text)

Mais c'est finalement les pages Wikipédia dont les contenus ne sont plus suffisament bien organisés... Ce serait sans doute tout de même possible en multipliant les conditions, mais pas très efficace.

L'API Wikipedia

Essayons avec la bibliothèque wikipedia de Python, faîte expressément pour ça !

Si besoin de l'installer (dans le shell Windows, et non Python) :

pip install wikipedia
Astuce
Pour une installation à l'université de Cergy (Python 2) :
py -m pip install wikipedia

Et maintenant dans la console Python de l'invite de commande Windows :

  1. import wikipedia
  2.  
  3. wikipedia.set_lang("fr")
  4.  
  5. summary_fr = wikipedia.summary("Aiguille Dibona", sentences=10)
  6. my_page = wikipedia.page("Aiguille Dibona")
  7.  
  8. print(summary_fr)
  9. print (my_page.url)

Héhé, c'est bien plus simple hein ?

Astuce
Un exemple complet d'usage de l'API Wikipedia est décrit ici : https://hg-map.fr/8-blog/90

L'API Wikidata

https://pypi.org/project/Wikidata/

  1. from wikidata.client import Client
  2.  
  3. client = Client()
  4. entity = client.get('Q1617977', load=True)
  5.  
  6. print(entity)
  7. print(entity.description)
  8. print(entity.data)
  9.  
  10. url_fr = entity.data['sitelinks']['frwiki']['url']
  11.  
  12. print(url_fr)
  13.  
  14. image_prop = client.get('P18')
  15. image = entity[image_prop]
  16.  
  17. print(image.image_resolution)
  18. print(image.image_url)
  19.  
  20. print(entity.__dict__)

Utilisation de modules Python tiers dans QGIS

Mais nous avons un problème. En effet la bibliothèque wikipedia de Python n'est pas forcément accessible depuis QGIS, même si vous l'avez déjà installé sur Windows.

Testez la seule importation du module wikipedia depuis la console Python de QGIS pour en avoir le cœur net :

import wikipedia

C'est parce que la version de Python utilisée par QGIS n'est pas la même que celle de votre machine, puisque sous Windows, QGIS embarque sa propre installation de Python.

 


Installer un module Python sur QGIS

Nous ouvrons donc cette nouvelle parenthèse pour installer un module Python dans QGIS cette fois : il vous faudra passer par le shell d'OSGeo4W.

Suivez ce tutoriel pour installer un module sur la version de Python utilisée par votre QGIS sous Windows, ou plus simplement :

Dans le shell d'OSGeo4W :

py3_env
python -m pip install wikipedia
Astuce
Pour une installation à l'université de Cergy (droits administrateur)
python -m pip install wikipedia --user

Vous aurez peut-être besoin de rallumer QGIS pour prendre en compte votre nouvelle bibliothèque.

Une fois la bibliotheque wikipedia installée, les codes précédents exécuté depuis QGIS fonctionneront également.

Bien, nous n'avons plus qu'à l'adapter à notre code de génération de cartes.

 


Vérification des résultats

Nous ajoutons également une comparaison in afin de vérifier la pertinence du résultat. En effet nous cherchons ici à partir de simples mots-clés dans le moteur de recherche de Wikipédia, et les résulats peuvent parfois être surprenants.

  1. import wikipedia
  2.  
  3. wikipedia.set_lang("fr")
  4. myWikiContent = wikipedia.summary("L'Ourson", sentences=3)
  5.  
  6. # Verifier la pertinence du resultat
  7. if str('Massif des Écrins').lower() in str(myWikiContent).lower():
  8.     pass
  9.  
  10. else:
  11.     myWikiContent = 'Oooouuuuppssss !!!'
  12.  
  13. print(myWikiContent)

Condition

Puis nous mettons le texte dans une condition if else avant d'afficher le texte. Peut-être pas nécessaire mais évite d'encombrer la carte avec des blocs vides.

Le contenu du texte est maintenant une concaténation de nos nouvelles variables :

  1.   ...
  2.   #Texte
  3.     if myWikiContent is None:
  4.         None
  5.     else:
  6.         TextCustom = QgsLayoutItemLabel(layout)
  7.         TextCustom.setText(myWikiContent+"\n\n"+myWikiLink)
  8.         TextCustom.setFont(QFont("Verdana", 11))
  9.         layout.addLayoutItem(TextCustom)
  10.         TextCustom.attemptMove(QgsLayoutPoint(230, 100, QgsUnitTypes.LayoutMillimeters))
  11.         TextCustom.attemptResize(QgsLayoutSize(60, 100, QgsUnitTypes.LayoutMillimeters))
  12.     ...

À partir d'un projet QGIS vide

Code complet pour générer les cartes à partir d'un projet QGIS vide (modifiez les chemins, vérifié sur QGIS 3.14.16 (Pi), préférez l'éditeur pour l'exécution)

  • Note personnelle : mieux gérer la largeur de la légende, ici certaines noms de sommets longs sont coupés.
  • Note personnelle : mieux gérer les liens Wikipédia, pour l'instant il mène à la page d'acceil de Wikipédia quand l'URL complète saute une ligne.

 

Aiguille Dibona

 


Autres exemples d'APIs

Ce chapitre est une parenthèse dans ce tutoriel, pour présenter quelques exemples d'exploitation d'APIs.

Un XML (http://mob.u-strasbg.fr/)

L'université de Strasbourg entretient une API fournissant des POIs, exemple :

La forme du XML généré est plus classique que ceux d'OSM :

<rss version="2.0">
  <poi>
    <id>1</id>
    <name>Central</name>
    <latitude>48.5803337</latitude>
    <longitude>7.7655637</longitude>
    <floor/>
    <type>campus</type>
    <description/>
    <url>http://fr.wikipedia.org/wiki/Campus_central_de_Strasbourg</url>
    <url_image/>
    <tags/>
    <capacity>0</capacity>
  </poi>
 
  <poi>
    <id>2</id>
...

Pour extraire les données de ce type d'XML :

import requests
import xml.etree.ElementTree
 
r = requests.get('http://mob.u-strasbg.fr/cgi-bin/odudsApi.py?method=search&amp;id=allpois')
root = xml.etree.ElementTree.fromstring(r.content)
for my_poi in root.findall('poi'):
    my_name = my_poi.find('name').text
    print(my_name)

Afficher les points à partir de leurs coordonnées

Le code suivant, dans un projet QGIS vide, affiche un fond de carte OSM, crée puis affiche une couche virtuelle de points, avec 2 champs, peuplée des points issus d'un XML. Toujours à partir de l'API de l'université de Strasbourg :

  1. import requests
  2. import xml.etree.ElementTree
  3.  
  4. # VIDER PROJET
  5. QgsProject.instance().removeAllMapLayers()
  6. iface.mapCanvas().refresh()
  7.  
  8. # PARAMETRES OSM
  9. urlWithParams = "type=xyz&amp;url=http://tile.openstreetmap.org/{z}/{x}/{y}.png"
  10. osm = QgsRasterLayer(urlWithParams, "OpenStreetMap", "wms")
  11.  
  12. # CREER UNE COUCHE VIRTUELLE DE POINTS
  13. point_vector = QgsVectorLayer("Point", "my_points", "memory")
  14.  
  15. # AFFICHER LES COUCHES
  16. QgsProject.instance().addMapLayer(osm)
  17. QgsProject.instance().addMapLayer(point_vector)
  18.  
  19. # AJOUTER CHAMPS DANS point_vector
  20. from qgis.PyQt.QtCore import QVariant
  21. pr = point_vector.dataProvider()
  22. pr.addAttributes([QgsField("id", QVariant.Int), QgsField("name", QVariant.String)])
  23. point_vector.updateFields()
  24.  
  25. # PARSER LE FICHIER
  26. r = requests.get('http://mob.u-strasbg.fr/cgi-bin/odudsApi.py?method=search&amp;id=allpois')
  27. root = xml.etree.ElementTree.fromstring(r.content)
  28.  
  29. for my_poi in root.findall('poi'):
  30.     my_id = int(my_poi.find('id').text)
  31.     my_name = str(my_poi.find('name').text)
  32.     my_longitude = float(my_poi.find('longitude').text)
  33.     my_latitude = float(my_poi.find('latitude').text)
  34.     #print(my_id, my_name, my_longitude, my_latitude)
  35.  
  36.     # EDITER COUCHE VIRTUELLE
  37.     f = QgsFeature()
  38.     f.setGeometry(QgsGeometry.fromPointXY(QgsPointXY(my_longitude, my_latitude)))
  39.     f.setAttributes([my_id, my_name])
  40.     pr.addFeature(f)
  41.     point_vector.updateExtents()
  42.     QgsProject.instance().addMapLayer(point_vector)
  43.  
  44. # RAFRAICHIR
  45. iface.mapCanvas().refresh()
  46.  
  47. # CRS
  48. my_crs=QgsCoordinateReferenceSystem(3857)
  49. QgsProject.instance().setCrs(my_crs)

Aller chercher l'image sous condition

...
for my_poi in root.findall('poi'):
    my_id = int(my_poi.find('id').text)
    my_name = str(my_poi.find('name').text)
 
    my_divers = str(my_poi.find('url_image').text)
 
    if my_divers == "None":
        my_divers = ""
    elif my_divers != "None":
        my_divers = "https://mob.u-strasbg.fr/geoloc/url_image/" + my_divers
 
    my_longitude = float(my_poi.find('longitude').text)
...

 

Un JSON (data.culture.gouv.fr, monuments historiques)

Le site data.culture.gouv.fr propose plusieurs APIs, dont les données sont fournis en JSON.

Exemple avec la liste des immeubles parisiens protégés au titre des monuments historiques :

Quand vous êtes confronté à ce type de fichiers, dont la visualisation en ligne n'est pas toujours aisée, il existe ce type d'outil en ligne qui vous permettra d'en comprendre la structure ou de tester votre fichier :

Sur jsonviewer utilisez le bouton Load JSON data et copiez-y l'URL (supprimez le « s » du https de votre URL si besoin). Sur geojson.io copiez vos données GEOJSON dans l'encart prévu, ou ajoutez l'adresse de votre fichier GEOJSON comme dans l'exemple ci-dessus (à la place de [URL]).

Pour parser ce type de JSON, et en extraire les données :

  1. import json
  2. import requests
  3.  
  4. response = requests.get("https://data.culture.gouv.fr/api/records/1.0/search/?dataset=liste-des-immeubles-proteges-au-titre-des-monuments-historiques&amp;q=&amp;facet=reg&amp;facet=dpt_lettre&amp;refine.dpt_lettre=Paris")
  5. my_file = json.loads(response.text)
  6.  
  7. my_json = my_file["records"]
  8.  
  9. for my_block in my_json:
  10.     my_fields = my_block["fields"]
  11.     print(my_fields["wadrs"])

Les lignes 4 et 5 en font un JSON lisible pour Python.

Les données qui nous intéressent sont dans le bloc records, nous l'isolons donc avec la ligne 7.

La ligne 9 boucle sur le bloc records, et pour chaque sous-blocs fields affiche les données du tag wadrs.

Un JSON (data.culture.gouv.fr, événements publics cibul)

https://public.opendatasoft.com/api/records/1.0/search/?dataset=evenements-publics-cibul&q=&rows=1354&facet=tags&facet=placename&facet=department&facet=region&facet=city&facet=date_start&facet=date_end&facet=pricing_info&facet=updated_at&facet=city_district&refine.tags=musique&refine.tags=gratuit&refine.department=Paris

Récupérer les informations :

import json
import requests
 
response = requests.get(
"https://public.opendatasoft.com/api/records/1.0/search/?dataset=evenements-publics-cibul&amp;q=&amp;rows=1354&amp;facet=tags&amp;facet=placename&amp;facet=department&amp;facet=region&amp;facet=city&amp;facet=date_start&amp;facet=date_end&amp;facet=pricing_info&amp;facet=updated_at&amp;facet=city_district&amp;refine.tags=musique&amp;refine.tags=gratuit&amp;refine.department=Paris")
 
my_file = json.loads(response.text)
my_json = my_file["records"]
 
for my_block in my_json:
    my_fields = my_block["fields"]
    my_uid = my_fields["uid"]
    my_title = my_fields["title"]
 
    my_geometry = my_block["geometry"]
    my_coordinates = my_geometry["coordinates"]
 
    print("my_uid :",my_uid, "my_coordinates :", my_coordinates, "my_title :", my_title)

Récupérer proprement les coordonnées

Pour utiliser les points, il faudra d'abord les splitter :

    my_geometry = my_block["geometry"]
    my_coordinates = my_geometry["coordinates"]
 
    my_lat = str(my_coordinates).split(", ")[0].replace("[", "")
    my_lon = str(my_coordinates).split(", ")[1].replace("]", "")

Afficher les points à partir de leurs coordonnées

import requests
import json
 
# VIDER PROJET
QgsProject.instance().removeAllMapLayers()
iface.mapCanvas().refresh()
 
# CREER UNE COUCHE VIRTUELLE DE POINTS
point_vector = QgsVectorLayer("Point", "my_points", "memory")
 
# AFFICHER LES COUCHES
QgsProject.instance().addMapLayer(point_vector)
 
# AJOUTER CHAMPS DANS point_vector
from qgis.PyQt.QtCore import QVariant
pr = point_vector.dataProvider()
pr.addAttributes([QgsField("uid", QVariant.Int), QgsField("title", QVariant.String)])
point_vector.updateFields()
 
# PARSER LE FICHIER
response = requests.get('https://public.opendatasoft.com/api/records/1.0/search/?dataset=evenements-publics-cibul&amp;q=&amp;rows=1354&amp;facet=tags&amp;facet=placename&amp;facet=department&amp;facet=region&amp;facet=city&amp;facet=date_start&amp;facet=date_end&amp;facet=pricing_info&amp;facet=updated_at&amp;facet=city_district&amp;refine.tags=musique&amp;refine.tags=gratuit&amp;refine.department=Paris')
 
my_file = json.loads(response.text)
my_json = my_file["records"]
 
for my_block in my_json:
    my_fields = my_block["fields"]
    my_uid = my_fields["uid"]
    my_title = my_fields["title"]
 
    my_geometry = my_block["geometry"]
    my_coordinates = my_geometry["coordinates"]
    my_lat = str(my_coordinates).split(", ")[0].replace("[", "")
    my_lon = str(my_coordinates).split(", ")[1].replace("]", "")
 
    print("my_uid :", my_uid,"my_coordinates :", my_coordinates,"my_lat :", my_lat,"my_lon :", my_lon,"my_title :", my_title)
 
    # EDITER COUCHE VIRTUELLE
    f = QgsFeature()
    f.setGeometry(QgsGeometry.fromPointXY(QgsPointXY(float(my_lon), float(my_lat))))
    f.setAttributes([my_uid, my_title])
    pr.addFeature(f)
    point_vector.updateExtents()
    QgsProject.instance().addMapLayer(point_vector)
 
# RAFRAICHIR
iface.mapCanvas().refresh()

 

Un GEOJSON

Pour ouvrir un fichier GEOJSON, ce ne sera pas compliqué (source : webgeodatavore) :

iface.addVectorLayer('https://d2ad6b4ur7yvpq.cloudfront.net/naturalearth-3.3.0/ne_50m_populated_places.geojson', 'Populated places', 'ogr')

Mais pour ouvrir des données GEOJSON brutes (construites à la volée par une API par exemple, ou simplement non-accessibles par URL), alors vous devrez d'abord les enregistrer dans un fichier, puis ouvir ce fichier.

Un GEOJSON (API zones isochrones https://openrouteservice.org)

Pensez à vérifier que la structure d'usage des APIs d'openrouteservice n'a pas changé ! Des exemples de codes peuvent être générés sur leur site.

Exemple avec l'excellente API de zones isochrones d'openrouteservice (votre obtiendrez une clé gratuite après une rapide inscription) :

import requests
 
body = {"locations":[[8.681495,49.41461],[8.686507,49.41943]],"range":[300,200]}
headers = {
    'Accept': 'application/json, application/geo+json, application/gpx+xml, img/png; charset=utf-8',
    'Authorization': 'VOTRE CLE API',
    'Content-Type': 'application/geo+json; charset=utf-8'
}
call = requests.post('https://api.openrouteservice.org/v2/isochrones/driving-car', json=body, headers=headers)
 
#print(call.status_code, call.reason)
#print(call.text)
 
with open('C:/Users/Georges/Downloads/myfile.geoson', 'w', encoding='utf-8') as outfile:
    outfile.write(call.text)
 
iface.addVectorLayer('C:/Users/Georges/Downloads/myfile.geoson', 'Isochrones', 'ogr')

Un autre GEOJSON (API directions https://openrouteservice.org)

import requests
 
headers = {
'Accept': 'application/json, application/geo+json, application/gpx+xml, img/png; charset=utf-8',
}
call = requests.get('https://api.openrouteservice.org/v2/directions/driving-car?api_key=VOTRE CLE API&amp;start=8.681495,49.41461&amp;end=8.686507,49.41943', headers=headers)
 
#print(call.status_code, call.reason)
#print(call.text)
 
with open('C:/Users/Georges/Downloads/trajet.geoson', 'w', encoding='utf-8') as outfile:
    outfile.write(call.text)
iface.addVectorLayer('C:/Users/Georges/Downloads/trajet.geoson', 'trajet', 'ogr')

Et pour boucler sur plusieurs mode de déplacement, et créer autant de fichier geojson, et les afficher : 

import requests
 
MesModes = ['foot-hiking', 'driving-car', 'wheelchair']
 
for i in MesModes:
 
    Chemin = 'C:/Users/Georges/Downloads/trajet_' + i + '.geoson'
 
    headers = {
    'Accept': 'application/json, application/geo+json, application/gpx+xml, img/png; charset=utf-8',
    }
    call = requests.get('https://api.openrouteservice.org/v2/directions/' + i + '?api_key=VOTRE CLE API&amp;start=2.011776,49.048384&amp;end=2.077019,49.038914', headers=headers)
 
    with open(Chemin, 'w', encoding='utf-8') as outfile:
        outfile.write(call.text)
    iface.addVectorLayer(Chemin, 'trajet', 'ogr')

 

Obsolète


Exercice

Testez le lien Wikipédia que nous avons mis dans nos cartes. Ils ne fonctionnent pas toujours. en effet nos PDF interprêtent ces liens en tant que chaînes de texte, et les sauts de ligne viennent parfois casser le lien.

Vous avez déjà eu des cours d'HTML l'année dernière, et bien sûr vous vous passionnez pour les solutions web.

Et bien à la place d'afficher le lien complet, qui ne marche même pas, mettez-nous un joli nom de lien fonctionnel !

  • Note personnelle : mettre le code d'affichage propre du lien en correction.

Standalone

Alors ça c'est toute une histoire... Essayez cette page ci-dessous (accès restreint) ou contactez-moi en privé !

https://hg-map.fr/pyqgis-tips

 


Joindre un CSV

layer = QgsProject.instance().mapLayersByName('simple_countries')[0]
 
# IMPORTER CSV
csv_path = 'file:///C:/Users/Georges/Downloads/temp_QGIS/countries_conversion.csv?delimiter=,'
my_csv1 = QgsVectorLayer(csv_path, 'mes_pays', 'delimitedtext')
QgsProject.instance().addMapLayer(my_csv1)
 
shpField='COUNTRY_HB'
csvField='country'
myJoin = QgsVectorLayerJoinInfo()
myJoin.setJoinFieldName(csvField)
myJoin.setTargetFieldName(shpField)
myJoin.setJoinLayerId(my_csv1.id())
myJoin.setUsingMemoryCache(True)
myJoin.setJoinLayer(my_csv1)
layer.addJoin(myJoin)
 
layer.dataProvider().forceReload()

 


Fonctions géométriques

Nous allons lister dans un fichier CSV les adresses des Hautes-Alpes se trouvant à 2kms des limites d'un glacier.

Au terme de cet exercice, il nous faut un code complet permettant de générer ce fichier à partir d'un projet QGIS vide, et d'y inclure le nombre d'adresses concernées.

Buffer

En guise d'entrainement, commencez par créer une zone tampon de 5 mètres autour de nos sommets.

  1. vlayer1 = QgsVectorLayer("C:/Users/Georges/Downloads/qgis/peaks/peaks.shp", "Mes sommets", "ogr")
  2. QgsProject.instance().addMapLayer(vlayer1)
  3.  
  4. bufferPath = "C:/Users/Georges/Downloads/qgis/buffers/buffer1.shp"
  5.  
  6. processing.run('native:buffer', {"INPUT": vlayer1, "DISTANCE": 5, \
  7. "SEGMENTS": 10, "END_CAP_STYLE": 0, "DISSOLVE": False, "OUTPUT": bufferPath})
  8.  
  9. vbuffer1 = QgsVectorLayer(bufferPath, "Buffer 1", "ogr")
  10. QgsProject.instance().addMapLayer(vbuffer1)
  11.  
  12. # Passer Buffer en bas
  13. root = QgsProject.instance().layerTreeRoot()
  14. myAlayer = root.findLayer(vbuffer1.id())
  15. myClone = myAlayer.clone()
  16. parent = myAlayer.parent()
  17. parent.insertChildNode(-1, myClone)
  18. parent.removeChildNode(myAlayer)

Les 2 premières lignes appellent la couche peaks, la 4ème stocke le chemin de notre future zone tampon.

C'est les lignes 6 et 7 qui à elles-seules créent le buffer grâce à l'instruction native/buffer. D'autres paramètres sont disponibles.

Les lignes 9 et 10 appellent ensuite notre couche contenant la zone tampon dans le canvas. Le reste du code est là pour passer les tampons sous les sommets, afin d'y voir clair.

Mais le code ci-dessus fonctionne très mal n'est-ce pas ?

  • Essayez de tricher
  • Essayez de comprendre pourquoi ça ne marche pas
  • Inspectez les unités géométriques utilisées dans votre canvas
  • Contrôlez vos logs dans View / Panels / Log Messages

Utilisez plutôt le code suivant pour créer des buffer de 2 kms autour de nos 108 000 adresses des Hautes-Alpes.

  1. vlayer2 = QgsVectorLayer("C:/Users/Georges/Downloads/qgis/data_BDTOPO_V3_Dep05_adresse/data_BDTOPO_V3_Dep05_adresse/ADRESSE_05.shp", "Addresses", "ogr")
  2. QgsProject.instance().addMapLayer(vlayer2)
  3.  
  4. bufferPath = "C:/Users/Georges/Downloads/qgis/buffers/buffer1.shp"
  5.  
  6. processing.run('native:buffer', {"INPUT": vlayer2, "DISTANCE": 2000, \
  7. "SEGMENTS": 10, "END_CAP_STYLE": 0, "DISSOLVE": False, "OUTPUT": bufferPath})
  8.  
  9. vbuffer1 = QgsVectorLayer(bufferPath, "Buffer 1", "ogr")
  10. QgsProject.instance().addMapLayer(vbuffer1)
  11.  
  12. # Passer Buffer en bas
  13. root = QgsProject.instance().layerTreeRoot()
  14. myAlayer = root.findLayer(vbuffer1.id())
  15. myClone = myAlayer.clone()
  16. parent = myAlayer.parent()
  17. parent.insertChildNode(-1, myClone)
  18. parent.removeChildNode(myAlayer)

Glaciers Alpins

Nous allons compter combien d'adresses se trouvent à 2 kms de distance des limites d'un glacier.

Téléchargez les glaciers contenus dans l'emprise de nos adresses Hautes-Alpines en utilisant l'extension QuickOSM (filtre natural=glaciers). Une couche sur ce modèle est disponible dans les pièces jointes de cet article.

Enregistrez-en un shape dans le même système de projection que vos adresses/tampons.

Pour contrôler la projection d'un layer.

{layer}.crs()

Intersection

L'algorithme de calcul d'intersection se présente de façon assez similaire que celle du buffer. D'autres paramètres sont disponibles.

  1. processing.run('qgis:intersection', {\
  2. "INPUT": {LayerInput},\
  3. "OVERLAY": {LayerOverlay},\
  4. "INPUT_FIELDS": {[list]},\
  5. "OVERLAY_FIELDS": {[list]},\
  6. "OVERLAY_FIELDS_PREFIX": {string},\
  7. "OUTPUT": {File}})

Introduisez l'algorithme ci-dessus dans votre code, exécutez-le puis zoomez, ci-dessous, sur le hameau de la Bérarde, au cœur de la zone de glaciers des Écrins, pour contrôler vos résultats.

  1. # Zoom sur la Berarde
  2. myselect = layerAdresse.getFeatures( QgsFeatureRequest().setFilterExpression ( u'"NOM_1" = \'LA BERARDE\'' ) )
  3. layerAdresse.selectByIds( [ f.id() for f in myselect ] )
  4. iface.mapCanvas().zoomToSelected(layerAdresse)
  5. iface.mapCanvas().zoomScale(100000)

Vos intersections doivent exister pour chaque adresse intersectées (présence de doublons opportuns).

Si votre vérification est bonne (même nombre d'entités intersectées dans un groupe de zones tampons proches que d'adresses au sein de ce groupe), alors il ne vous reste plus qu'à compter le nombre total d'intersections.

{Layer}.featureCount()

Export CSV

Vous savez déjà exporter un fichier CSV (voir chapitre Exporter des fichiers plus haut dans ce court), à vous de jouer !

N'oubliez pas d'inclure le nombre d'adresses concernées dans le nom du fichier.

Nettoyage

La couche de nos zone tampons et celle de nos intersections sont à présent inutiles. Nous avons fait le job, nous pouvons donc supprimez ces couches mais aussi fermer QGIS.

  1. from osgeo import ogr
  2. ...
  3. # Nettoyage
  4. QgsProject.instance().removeMapLayer(layerBuffer)
  5. QgsProject.instance().removeMapLayer(intersection)
  6.  
  7. driver = ogr.GetDriverByName("ESRI Shapefile")
  8.  
  9. driver.DeleteDataSource(bufferPath)
  10. driver.DeleteDataSource(interstectPath)
  11.  
  12. os._exit(0)

Il semble toute fois que la suppression totale des shapes soit impossible à réaliser à la suite du code complet. Sans doute une mémoire cache à vider dans QGIS ou Python.

Code complet de la génération des buffers, des intersections et du fichier CSV (modifiez les chemins, fonctionne dans l'éditeur Python exclusivement)

  • Note personnelle : aller plus loin en se passant de la génération des shapes des buffers (couche virtuelle à la place), des glaciers (API Overpass d'OSM) et des intersections (couche virtuelles). Afin d'obtenir le simple fichier CSV des adresses à 2 kms d'un glacier, et leur nombre.

Classification

#You can use any of these classification method classes:
#QgsClassificationQuantile()
#QgsClassificationEqualInterval()
#QgsClassificationJenks()
#QgsClassificationPrettyBreaks()
#QgsClassificationLogarithmic()
#QgsClassificationStandardDeviation()

Intervalles égaux

  1. # Set layer name and desired paremeters
  2. layer_name = 'simple_countries'
  3. ramp_name = 'Spectral'
  4. value_field = 'population'
  5. num_classes = 5
  6. classification_method = QgsClassificationEqualInterval()
  7.  
  8. layer = QgsProject().instance().mapLayersByName(layer_name)[0]
  9.  
  10. # change format settings as necessary
  11. format = QgsRendererRangeLabelFormat()
  12. format.setFormat("%1 - %2")
  13. format.setPrecision(2)
  14. # format.setTrimTrailingZeroes(True)
  15.  
  16. default_style = QgsStyle().defaultStyle()
  17. color_ramp = default_style.colorRamp(ramp_name)
  18.  
  19. renderer = QgsGraduatedSymbolRenderer()
  20. renderer.setClassAttribute(value_field)
  21. renderer.setClassificationMethod(classification_method)
  22. renderer.setLabelFormat(format)
  23. renderer.updateClasses(layer, num_classes)
  24. renderer.updateColorRamp(color_ramp)
  25.  
  26. layer.setRenderer(renderer)
  27. layer.triggerRepaint()

 

 


Exercice

Maintenant vous savez :

  • Exporter des cartes en série
  • Intersecter des couches géographiques

Et bien ajoutez donc dans le sous-titre de vos cartes la commune où se trouve le sommet sélectionné.

À vous de jouer !

Vous allez devoir utiliser la couche des IRIS (disponible en pièce jointe de cet article) et l'intersecter avec vos sommets :

  1. ...
  2. iris = QgsVectorLayer("C:/Users/Georges/Downloads/qgis/data_IRIS_2019/data_IRIS_2019/iris.shp", "IRIS", "ogr")
  3. ...
  4. QgsProject.instance().addMapLayer(iris)
  5. ...
  6. # Intersecter les sommets et les IRIS
  7. peaks_iris_path = "C:/Users/Georges/Downloads/qgis/exports/peaks_iris.shp"
  8.  
  9. processing.run('qgis:intersection', {\
  10. "INPUT": peaks,\
  11. "OVERLAY": iris,\
  12. "INPUT_FIELDS": ["OSM_ID", "NAME", "OTHER_TAGS"],\
  13. "OVERLAY_FIELDS": ["CODE_IRIS", "NOM_COM"],\
  14. "OVERLAY_FIELDS_PREFIX": "",\
  15. "OUTPUT": peaks_iris_path})
  16.  
  17. peaks_iris = QgsVectorLayer(peaks_iris_path, "Sommets", "ogr")
  18. QgsProject.instance().addMapLayer(peaks_iris)
  19.  
  20. # Retirer ancienne couche des sommets et celle des IRIS
  21. QgsProject.instance().removeMapLayer(peaks)
  22. QgsProject.instance().removeMapLayer(iris)

Puis appeler le nouveau champ disponible (NOM_COM) et modifier la chaîne formatée qui remplit le sous-titre :

...
    commune = feat['NOM_COM']
    ...
    subtitle.setText("Identifiant OSM : {} - Commune : {}".format((str(id_peak)), (commune)))
...

Notez qu'au lieu de l'ancien %, nous utilisons ici une méthode de formatage plus moderne (.format()) 

Code complet de l'export des cartes avec les noms de commune en sous-titre et les liens aux APIs (modifiez les chemins, testez ce code de préférence dans l'éditeur Python, pas la console).


Écriture de points

Dans un shape

Le code ci-dessous est inspiré de l'excellent site opensourceoptions.com.

En partant du principe que vous avez ouvert un shape nommé my_points, en 4326, contenant un champ id et un champ text, le code suivant crée un point sur l'Université de Cergy-Pontoise :

layers = QgsProject.instance().mapLayersByName('my_points')[0]
layer = QgsVectorLayer(layers.dataProvider().dataSourceUri(), '', 'ogr')
caps = layer.dataProvider().capabilities()
 
if caps &amp; QgsVectorDataProvider.AddFeatures:
    feat = QgsFeature(layer.fields())
    feat.setAttributes([1, 'Université de Cergy-Pontoise'])
    feat.setGeometry(QgsGeometry.fromPointXY(QgsPointXY(49.0297879, 2.0616833)))
    res, outFeats = layer.dataProvider().addFeatures([feat])

Dans un shape, éditer seulement un champ

with edit(mylayer):
    for feature in mylayer.getFeatures():
        feature['myfield'] = 'Nouvelle valeur'
        mylayer.updateFeature(feature)

Dans une couche virtuelle ou un GEOJSON

Le code ci-dessous est inspiré de l'excellent site opensourceoptions.com.

Nous créons une couche virtuelle de points, y ajoutons deux champs id et name, avant d'y écrire un 1er point :

point_vector = QgsVectorLayer("Point", "my_points", "memory")
 
QgsProject.instance().addMapLayer(point_vector)
 
from qgis.PyQt.QtCore import QVariant
pr = point_vector.dataProvider()
pr.addAttributes([QgsField("id", QVariant.Int), QgsField("name", QVariant.String)])
point_vector.updateFields()
 
f = QgsFeature()
f.setGeometry(QgsGeometry.fromPointXY(QgsPointXY(49.0297879, 2.0616833)))
f.setAttributes([1,'Université de Cergy-Pontoise'])
pr.addFeature(f)
point_vector.updateExtents()
QgsProject.instance().addMapLayer(point_vector)

 


Fusion de couches

Un cas d'école : un collaborateur récupère un grand nombre de fichiers de données et souhaite les exploiter proprement. Il s'agit de plusieurs centaines de fichiers cadastraux allemands.

Téléchargez, dézippez et ouvrez la couche Koln gml.zip.

Le but est d'automatiser :

  1. L'ouverture d'un grand nombre de fichiers GML
  2. Leur conversion en shapes
  3. Leur édition de structure
  4. Enfin, leur fusion en un seul shape.
  1. import os
  2. import shutil
  3. from os import listdir
  4. from os.path import isfile, join
  5.  
  6. # On supprime les futures variables utilisées, car certaines semblent gêner les exécutions multiples du script
  7. my_var = ['ShapesOuvertsNoms', 'LayersShapesOuverts', 's', 'n', 'mon_shape', 'repertoireShapes', 'pathShape', 'nom', 'writer']
  8. for v in my_var:
  9. if v in globals():
  10. del v
  11.  
  12. monRepertoire = 'C:/Users/Georges/Downloads/Koln/lines/'
  13.  
  14. # On liste les GML en triant d'éventuels fichiers non GML (GFS par exemple)
  15. mesGML_ = [f for f in listdir(monRepertoire) if isfile(join(monRepertoire, f))]
  16. mesGML = [val for val in mesGML_ if val.endswith('.gml')]
  17.  
  18. QgsProject.instance().removeAllMapLayers()
  19. iface.mapCanvas().refresh()
  20.  
  21. NbreGmlValid = 0
  22. NbreGmlNoValid = 0
  23.  
  24. # On va stocker les noms des GML correctement ouverts mais aussi leurs appels en layers
  25. GmlOuvertsNoms = []
  26. LayersGmlOuverts = []
  27.  
  28. # Boucle pour ouvrir tous les GML sur QGIS
  29. for gml in mesGML :
  30.     mon_gml = QgsVectorLayer(monRepertoire+gml, gml, 'ogr')
  31.  
  32.     # On précise la projection
  33.     crs = mon_gml.crs()
  34.     crs.createFromId(25832)
  35.     mon_gml.setCrs(crs)
  36.     QgsProject.instance().addMapLayer(mon_gml)
  37.  
  38.     # On compte et liste les GML correctement ouverts mais aussi leurs appels en layers
  39.     if mon_gml.isValid():
  40.         print(gml, 'est valide ! OKLM.')
  41.         NbreGmlValid += 1
  42.         GmlOuvertsNoms.append(gml)
  43.         LayersGmlOuverts.append(mon_gml)
  44.  
  45.         # On compte les GML en erreur
  46.         if not mon_gml.isValid():
  47.         print(gml, 'est invalide ! Nous le supprimons.')
  48.         QgsProject.instance().removeMapLayer(mon_gml.id())
  49.         NbreGmlNoValid += 1
  50.  
  51. print('Nombre de GML correctement ouverts :',NbreGmlValid)
  52. print('Nombre de GML invalides et supprimés :',NbreGmlNoValid)
  53.  
  54. # EXPORT EN SHAPE
  55. repertoireShapes = 'shapes'
  56. pathShape = os.path.join(monRepertoire, repertoireShapes)
  57.  
  58. # On va stocker les noms des shapes correctement ouverts mais aussi leurs appels en layers
  59. ShapesOuvertsNoms = []
  60. LayersShapesOuverts = []
  61.  
  62. # Création d'un répertoire pour accueillir nos shapes
  63. if not os.path.exists(pathShape):
  64.     os.makedirs(pathShape)
  65. else:
  66.     shutil.rmtree(pathShape)
  67.     os.makedirs(pathShape)
  68.  
  69. # Boucle pour exporter nos GML en shape
  70. for s,n in zip(LayersGmlOuverts, GmlOuvertsNoms):
  71.     nom = n.replace('.gml', '')
  72.     writer = QgsVectorFileWriter.writeAsVectorFormat(s, pathShape+'/' + nom + '.shp', 'utf-8', driverName='ESRI Shapefile', onlySelected=False)
  73.  
  74.     mon_shape = QgsVectorLayer(pathShape+'/' + nom + '.shp', nom, 'ogr')
  75.  
  76.     ShapesOuvertsNoms.append(nom)
  77.     LayersShapesOuverts.append(mon_shape)
  78.  
  79. # OUVRIR LES SHAPES
  80. QgsProject.instance().removeAllMapLayers()
  81. for x in LayersShapesOuverts:
  82.     QgsProject.instance().addMapLayer(x)
  83.  
  84. # On va supprimer tous les champs de tous nos shapes sauf 2 champs
  85. # Lister les champs d'un layer si besoin et tester avec my_vectorlayer = iface.activeLayer()
  86. # D'abord on crée la liste de champs à supprimer
  87. for my_shape in LayersShapesOuverts:
  88.     prov = my_shape.dataProvider()
  89.     field_names = [field.name() for field in prov.fields()]
  90.     listeChampSuppression = []
  91.  
  92.     for a in field_names:
  93.         listeChampSuppression.append(a)
  94.  
  95.     listeChampSuppression.remove('gml_id')
  96.     listeChampSuppression.remove('measuredHe')
  97.  
  98.     # Puis on supprime la liste de champs à supprimer
  99.     for c in listeChampSuppression:
  100.         with edit(my_shape):
  101.             my_field = my_shape.fields().indexFromName(c)
  102.             my_shape.dataProvider().deleteAttributes([my_field])
  103.         my_shape.updateFields()
  104.  
  105. # Fusion des shapes
  106. parameters = {'LAYERS': ShapesOuvertsNoms,'CRS': 'EPSG:25832','OUTPUT': monRepertoire+'results.shp'}
  107. processing.run("qgis:mergevectorlayers", parameters)
  108. maFusion = QgsVectorLayer(monRepertoire+'results.shp', 'Ma jolie fusion', 'ogr')
  109. QgsProject.instance().addMapLayer(maFusion)

Calculatrice de champs pour Python

Ci-dessous un exemple d'utilisation de cette calculatrice Python à la documentation hasardeuse.

Dans Global expression :

def myfunction(myparam):
    if myparam == 'Aiguille Dibona':
        return 2

Puis dans Formula:

Observez que le champ NAME est appelé en le mettant entre chevrons <>.

value = myfunction(<NAME>)

 

 

 

 

 

Liens ou pièces jointes
Accéder à cette adresse URL (https://hg-map.fr/extern/data/shapes/france/chemin_de_fer.zip)chemin_de_fer.zip[ ]0 Ko
Télécharger ce fichier (data_BDTOPO_V3_Dep05_adresse.zip)data_BDTOPO_V3_Dep05_adresse.zip[ ]3889 Ko
Télécharger ce fichier (data_IRIS_2019.zip)data_IRIS_2019.zip[ ]45905 Ko
Télécharger ce fichier (decathlon_france.zip)decathlon_france.zip[308 magasins Décathlon français depuis OSM le 27 décembre 2020]11 Ko
Accéder à cette adresse URL (https://hg-map.fr/extern/data/shapes/france/eau.zip)eau.zip[ ]0 Ko
Télécharger ce fichier (glaciers.zip)glaciers.zip[ ]231 Ko
Télécharger ce fichier (iso_iris.zip)iso_iris.zip[Des zones isochrones à 15 minutes autour de 308 POIs.]12125 Ko
Télécharger ce fichier (Koln GML.zip)Koln gml.zip[ ]2818 Ko
Télécharger ce fichier (peaks.zip)peaks.zip[ ]14 Ko
Télécharger ce fichier (peaks_selection.zip)peaks_selection.zip[ ]1 Ko
Télécharger ce fichier (simple_countries.zip)simple_countries.zip[ ]1880 Ko
Accéder à cette adresse URL (https://hg-map.fr/extern/data/shapes/france/sol.zip)sol.zip[ ]0 Ko
Accéder à cette adresse URL (https://hg-map.fr/extern/data/shapes/france/troncons_routes.zip)troncons_routes.zip[ ]0 Ko