[PYTHON] Grundlegende Karteninformationen - Verwenden Sie ein numerisches Höhenmodell

Einführung

Beschreibt die Verwendung des vom National Land Research Institute veröffentlichten numerischen Höhenmodells (DEM). Beschreibt, wie DEM-Daten in 2D-Höhendaten konvertiert werden.

Dateiname

Daten der Basiskarte (numerisches Höhenmodell) 5 m Maschenweite (Höhe) Der Dateiname von hat das folgende Format.

FG-GML-pppp-qq-rr-DEM5X-yyyymmdd.xml

pppp ist die primäre Maschenzahl qq ist die sekundäre Maschenzahl rr ist die tertiäre Maschennummer X ist eine Vermessungsmethode, A: Laservermessung, B: Fotovermessung. yyyymmdd ist das Datum.

--Beispiel

FG-GML-5338-06-89-DEM5A-20130702.xml

Maschennummer

Das primäre Netz ist ein Quadrat mit einer Breite von 1 Grad Länge und einer Höhe von 40 Minuten Breite. Die primäre Maschenzahl wird durch den Breiten- und Längengrad des südwestlichen Endes dieses Quadrats dargestellt. Die primäre Maschennummer wird durch 4 Ziffern dargestellt. Die ersten beiden Ziffern sind das 1,5-fache des Breitengrads. Die letzten beiden Ziffern sind die letzten beiden Ziffern der Länge.

Das sekundäre Netz ist das primäre Netz, das vertikal und horizontal in acht gleiche Teile unterteilt ist. Zahl (Nord-Süd, Ost-West) mit dem südwestlichen Ende als (0,0). Der Nummernkreis ist 0-7. Die sekundäre Maschennummer ist eine zweistellige Nummer, die in der Reihenfolge Nord-Süd und Ost-West angeordnet ist.

Das Tertiärnetz ist das Sekundärnetz, das vertikal und horizontal in 10 gleiche Teile unterteilt ist. Die Nummerierung entspricht der Idee des sekundären Netzes.

image

Dateiinhalt

Es werden nur die Hauptinformationen extrahiert

Es ist ein Array von primären, sekundären und tertiären Maschennummern.

  <mesh>53380689</mesh>

Beschreibt die Anzahl der Gitterzellen, aus denen das Netz besteht. Das niedrige Attribut ist eine Zahl, die das nordwestliche Ende angibt. Es ist immer (0,0). Das hohe Attribut ist eine Zahl, die das südöstliche Ende angibt. Das 5-m-Netz ist auf (224,149) und das 10-m-Netz auf (1124.749) festgelegt.

  <gml:limits>
  	<gml:GridEnvelope>
  		<gml:low>0 0</gml:low>
  		<gml:high>224 149</gml:high>
  	</gml:GridEnvelope>
  </gml:limits>

Die Höhendaten jeder Gitterzelle werden aufgelistet. Die im tupleList-Tag enthaltenen Daten sind die Höhendaten. Wenn keine Daten zur Wasseroberfläche oder Höhe vorhanden sind, beträgt der Höhenwert "-9999".

  <gml:DataBlock>
  	<gml:rangeParameters>

<gml: QuantityList uom = "DEM-Bestandteile" /> </gml:rangeParameters> gml:tupleList (Höhendaten) Beispiel) Bodenoberfläche, 1055.33 </gml:tupleList> </gml:DataBlock>

  			<gml:sequenceRule  order="+x-y">Linear</gml:sequenceRule>

Geben Sie die Sequenz der Daten mit sequenceRule an. "+ X" in der Reihenfolge zeigt an, dass x in der positiven Richtung ist (West → Ost), und "-y" zeigt an, dass y in der negativen Richtung ist (Nord → Süd). Vom nordwestlichen Ende bis zum südöstlichen Ende werden die Daten, die wie unten gezeigt in einer Reihe angeordnet sind, im tupleList-Tag platziert. 1.JPG

  <gml:startPoint>173 0</gml:startPoint>

Geben Sie den Startpunkt der Rasterzelle mit startPoint an. ** Die im DataBlock platzierten Daten sind die Rasterzellennummern (173,0) bis (224,149). ** ** ** ** Bitte beachten Sie, dass die tatsächliche Anzahl der Daten nicht unbedingt (gml: niedrig + 1) * (gml: hoch + 1) ist. ** ** **

Datenabruf

Nachdem wir den Inhalt des XML kennen, werden wir die DEM-Daten als zweidimensionale Höhendaten extrahieren. Ich werde Python verwenden, um es abzurufen. Es ist ein 2D-Ndarray-Typ von Numpy. Die Arraygröße ist auf die Anzahl der Spalten (gml: niedrig + 1) und die Anzahl der Zeilen (gml: hoch + 1) festgelegt. Alle Daten vor startPoint sind auf "-9999" eingestellt.

hyoukou.py


import sys
import re
import numpy as np
import matplotlib.pyplot as plt

with open(sys.argv[1], "r") as f:
    # mesh
    r = re.compile("<mesh>(.+)</mesh>")
    for ln in f:
        m = r.search(ln)
        if m != None:
            mesh = m.group(1)
            break
    # grid len
    r = re.compile("<gml:high>(.+) (.+)</gml:high>")
    for ln in f:
        m = r.search(ln)
        if m != None:
            xlen = int(m.group(1)) + 1
            ylen = int(m.group(2)) + 1
            break
    # start
    r = re.compile("<gml:tupleList>")
    for ln in f:
        m = r.search(ln)
        if m != None:
            break
    # data
    data = []
    r  = re.compile("</gml:tupleList>")
    r2 = re.compile(",(.+)")
    for ln in f:
        m = r.search(ln)
        if m != None:
            break
        else:
            m = r2.search(ln)
            if m != None:
                data.append(float(m.group(1)))
    # start point
    startx = starty = 0
    r = re.compile("<gml:startPoint>(.+) (.+)</gml:startPoint>")
    for ln in f:
        m = r.search(ln)
        if m != None:
            startx = int(m.group(1))
            starty = int(m.group(2))
            break

data2 = np.empty(xlen*ylen)
start_pos = starty*xlen + startx

for i in range(xlen*ylen):
    if i < start_pos:
        data2[i] = -9999.
    else:
        data2[i] = data[i-start_pos]

data = data2.reshape(ylen, xlen)

#print(data)
plt.imshow(data)
plt.colorbar()
plt.show()

So rufen Sie hyoukou.py auf Rufen Sie mit dem DEM-Dateinamen als Argument auf. Zeigen Sie 2D-Höhendaten mit matplotlib als Bild an. Wenn kein Datenverlust vorliegt (-9999.), Kann er normal angezeigt werden.

>python hyoukou.py FG-GML-5338-06-89-DEM5B-20130702.xml

Recommended Posts

Grundlegende Karteninformationen - Verwenden Sie ein numerisches Höhenmodell
Basiskarteninformationen mithilfe der Python-Geotiff-Konvertierung numerischer Höhendaten