Essayez d'imaginer les données d'élévation du National Land Research Institute avec Python
Un programme pour convertir le modèle numérique d'élévation du National Land Research Institute en PNG avec Python créé par le passe-temps de mon supérieur () lors de mon travail précédent.
C'est du gâchis de l'endormir, alors je laisse ça au prochain apprenant
Inutile de dire que le code ci-dessous n'est pas le "meilleur code"!
De plus, je ne suis pas un spécialiste de Python, donc j'écris en pensant: "Je serais heureux si quelqu'un le refactorisait."
* Supplément (2020/02/14) Ce code est mauvais car il a été écrit par moi la première année. S'il vous plaît soyez prudente.
Ce que vous voulez faire
Je souhaite obtenir les données numériques d'élévation du National Land Research Institute et en faire une image en échelle de gris à l'aide de Python
Nous allons générer un fichier PNG pour chaque fichier de données
Veuillez vous référer à l'article Qiita suivant pour l'explication du fichier xml qui peut être supprimé car il est très facile à comprendre.
https://qiita.com/tobira-code/items/43a23362f356198adce2
Ou plutôt, avec cet article, je n'ai pas à écrire un nouvel article ...
Préparer les données
Immédiatement, sélectionnez n'importe quel maillage sur la page de téléchargement du National Land Research Institute et obtenez-le.
https://fgd.gsi.go.jp/download/menu.php
Cliquez sur le bouton "Aller à la sélection de fichier" dans "Informations cartographiques de base Modèle d'élévation numérique" pour passer à l'écran de téléchargement.
L'inscription de membre est requise pour télécharger les données!
Cette fois, nous utiliserons les données de "Ligne de contour de carte de terrain" de "Maillage 10m", donc sélectionnez le bouton radio et cochez le bouton sur le côté gauche en tant que tel, et cliquez sur n'importe quel maillage.
À partir du fichier Zip déposé, décompressez le fichier "FG-GML-0000-0-dem10b-yyyymmdd.xml" et placez-le dans un répertoire approprié.
Ceci termine la préparation des données.
Le code est ci-dessous.
import os
import re
import sys
import numpy as np
from PIL import Image
import xml.etree.ElementTree as ET
#entrée de fichier xml
DATA = sys.argv[1]
tree = ET.parse(DATA)
root = tree.getroot()
namespace = {
'xml': 'http://fgd.gsi.go.jp/spec/2008/FGD_GMLSchema',
'gml': 'http://www.opengis.net/gml/3.2'
}
dem = root.find('xml:DEM', namespace)
#Numéro de maillage
mesh = dem.find('xml:mesh', namespace).text
#Nombre de cellules disposées(La valeur réelle est ajoutée par 1)
high = dem.find('./xml:coverage/gml:gridDomain/gml:Grid/gml:limits/gml:GridEnvelope/gml:high', namespace).text.split(' ')
highX = int(high[0]) + 1
highY = int(high[1]) + 1
#Réglage de la taille de l'image(Le nombre de données==Nombre de pixels)
imgSize = highX * highY
#Tableau de données d'altitude
dem_text = dem.find('./xml:coverage/gml:rangeSet/gml:DataBlock/gml:tupleList', namespace).text
data = re.findall(',(.*)\n', dem_text)
dataNp = np.empty(imgSize)
for i in range(len(data)):
if(data[i] == "-9999.00"):
dataNp[i] = 0
else:
dataNp[i] = float(data[i])
#Coordonnées de début des données
startPoint = dem.find('./xml:coverage/gml:coverageFunction/gml:GridFunction/gml:startPoint', namespace).text.split(' ')
startPointX = int(startPoint[0])
startPointY = int(startPoint[1])
startPosition = startPointY * highX + startPointX
##Lancer la génération d'image
#Lorsque le nombre de données est insuffisant(S'il y a des espaces au-dessus et en dessous de l'image)
if(len(dataNp) < imgSize):
add = []
#Lorsque les données en bas de l'image sont insuffisantes
if(startPosition == 0):
for i in range(imgSize - len(dataNp)):
add.append(0)
dataNp.extend(add)
#Lorsque les données en haut et en bas de l'image sont insuffisantes
else:
for i in range(startPosition):
add.append(0)
dataNp[0:0] = add
add = []
for i in range(imgSize - len(dataNp) - len(add)):
add.append(0)
dataNp.extend(add)
#Conversion d'entiers 8 bits de données
dataNp = (dataNp / 15).astype(np.uint8) #Divisez dataNp par 15 pour correspondre au point le plus élevé du mont Fuji à 255
data = dataNp.reshape(highY, highX)
#Imagerie des données numériques d'élévation
pilImg = Image.fromarray(np.uint8(data))
pilImg = pilImg.resize((int(highX), int(highY)), Image.LANCZOS) # NEAREST
canvas = Image.new('RGB', (highX, highY), (0, 0, 0))
canvas.paste(pilImg, (0, 0))
canvas.save('dem.png', 'PNG', quality=100)
L'exécution est OK avec python .py file .xml file .
Supplément
Rechercher diverses données
Vous pouvez utiliser la méthode ci-dessus pour obtenir le premier numéro de maillage, le nombre de tableaux de cellules, les données d'altitude, etc. Vous pouvez le rechercher avec une expression régulière en tournant la boucle comme .
Pour moi, je préfère chercher en tournant la boucle, donc je le recommande car c'est plus facile, mais ensuite je couvrirai cet article ...
"DataNp / 15" qui apparaît dans la seconde moitié
Je suis désolé pour le nombre magique de gorigori, mais c'est un peu ennuyeux, alors je le laisse tel quel.
Dans le modèle d'élévation numérique, l'élévation de chaque point à l'exception de la surface de l'eau est stockée sous forme de valeur numérique (bien que je ne sache pas si c'est le cas pour tous les MNT). Le point au sommet du mont Fuji doit être stocké dans un fichier xml tel que 3776.0.
Si cette valeur est directement convertie en image PNG en échelle de gris, l'image en noir et blanc sera représentée par une valeur comprise entre 0 et 255, de sorte que tous les terrains avec une altitude de 255 m ou plus seront complètement blancs. Afin d'éviter cela, le code ci-dessus est ajusté en divisant par 15 afin de le rendre 255 ou moins basé sur l'altitude de 3770m au sommet du mont Fuji.
Donc, si vous voulez exprimer le mont Hakodate, qui a une altitude de 333m, de 0 à 255 au lieu du sommet du mont Fuji, vous devez le diviser par 1,3.
Ici, il semble facile de vérifier la valeur maximale du fichier xml et de définir la valeur à diviser automatiquement.
Enfin
Comme ça, quand j'ai converti Hakodate en PNG, c'est devenu comme ça.
L'altitude de la ville est presque 0, et j'ai souvent l'impression que l'altitude du mont Hakodate monte soudainement.
Cependant, n'est-ce pas gênant?
Donc, la prochaine fois, j'écrirai un article qui lit plusieurs fichiers, les trie par numéro de maillage et fait une vue d'ensemble.
Cela vous donnera une image assez intéressante.
Je pense que cet article devrait être revu dans un proche avenir et réécrit pour le rendre plus facile à comprendre. .. .. (Peut être pas)
Recommended Posts