パイクラおじさんの日記

MinecraftでPythonを勉強するおじさんの日記です。

標高データをMatplotlibで表示してみる

# coding: utf-8

import xml.etree.ElementTree as ET
import numpy as np
import matplotlib.pyplot as plot
import sys

GEO_DIR = "FG-GML-{0:04d}-{1:02d}-DEM5A"
GEO_XML = "FG-GML-{0:04d}-{1:02d}-{2:02d}-DEM5A-20161001.xml"

if len(sys.argv) < 4:
    print("xmlplot.py mesh1 mesh2 mesh3")
    exit(-1)

mesh1 = int(sys.argv[1])
mesh2 = int(sys.argv[2])
mesh3 = int(sys.argv[3])

dir = GEO_DIR.format(mesh1, mesh2)
fname = GEO_XML.format(mesh1, mesh2, mesh3)

tree = ET.parse(dir+"/"+fname)
root = tree.getroot()

tl = root.find('./{http://fgd.gsi.go.jp/spec/2008/FGD_GMLSchema}DEM/{http://fgd.gsi.go.jp/spec/2008/FGD_GMLSchema}coverage/{http://www.opengis.net/gml/3.2}rangeSet/{http://www.opengis.net/gml/3.2}DataBlock/{http://www.opengis.net/gml/3.2}tupleList')
if tl is None:
    print("{http://www.opengis.net/gml/3.2}tupleList is not found")
    exit(-1)

lines = tl.text.split()
array = np.zeros(len(lines))
i = 0
for l in lines:
    (t, h) = l.split(",")
    array[i] = float(h)
    i += 1

img_array = array.reshape((150, 225))
plot.imshow(img_array)
plot.colorbar()
plot.show()

注意:reshapeの引数が逆だったので、直しました。

実行すると…

$ python xmlplot.py 5340 22 0

ウインドウが開いて、Fig. 1のようなグラフ(?)が表示されました。

f:id:pycra:20170917214713p:plain
Fig. 1

3次メッシュ番号0はいいけど、1, 2, 3などはFig. 2のように真っ黄色に紫のグラフになりました。

<gml:tupleList>内の「内水面, -9999」の行がいけないらしい。

これは池なんだろう。

f:id:pycra:20170917215902p:plain
Fig. 2

# coding: utf-8

import xml.etree.ElementTree as ET
import numpy as np
import matplotlib.pyplot as plot
import sys

GEO_DIR = "FG-GML-{0:04d}-{1:02d}-DEM5A"
GEO_XML = "FG-GML-{0:04d}-{1:02d}-{2:02d}-DEM5A-20161001.xml"

if len(sys.argv) < 4:
    print("xmlplot.py mesh1 mesh2 mesh3")
    exit(-1)

mesh1 = int(sys.argv[1])
mesh2 = int(sys.argv[2])
mesh3 = int(sys.argv[3])

dir = GEO_DIR.format(mesh1, mesh2)
fname = GEO_XML.format(mesh1, mesh2, mesh3)

tree = ET.parse(dir+"/"+fname)
root = tree.getroot()

tl = root.find('./{http://fgd.gsi.go.jp/spec/2008/FGD_GMLSchema}DEM/{http://fgd.gsi.go.jp/spec/2008/FGD_GMLSchema}coverage/{http://www.opengis.net/gml/3.2}rangeSet/{http://www.opengis.net/gml/3.2}DataBlock/{http://www.opengis.net/gml/3.2}tupleList')
if tl is None:
    print("{http://www.opengis.net/gml/3.2}tupleList is not found")
    exit(-1)

lines = tl.text.split()
array = np.zeros(len(lines))
i = 0
for l in lines:
    (t, h) = l.split(",")
    hval = float(h)
    if hval == -9999:
        array[i] = -9
    else:
        array[i] = hval
    i += 1

img_array = array.reshape((150, 225))
plot.imshow(img_array)
plot.colorbar()
plot.show()

判定を入れて-9999だったら-9にするようにしたらFig. 3のようになった。

f:id:pycra:20170917220737p:plain
Fig. 3

しかし、便利だなぁ〜!すぐに図として確認できる。ただし、みんなレンジが違うので注意。