【Python】GPXファイルからXYZファイルへの変換

はじめに

 山形大学大学院修士1年のHagianです。この記事は、Python Advent Calendar 2022 の15日目となります。

やりたかったこと

1. GPXファイル（GPSやGNSSの情報をXML形式で記録したもの）から緯度・経度・高さを読み込む
2. UTM（ユニバーサル横メルカトル）座標に変換する
3. 変換した座標データをXYZファイル（X, Y, Z方向の基準点からの距離をスペース区切りで記録したもの）で書き出す

GPXとXYZについて

①GPX

 GPXは英語で、GPS Exchange Formatといい、アプリケーション間でGPSやGNSSの情報をやり取りするために用いられるフォーマットです。初耳だ、何やら難しそう… と思われる方もいらっしゃるかもしれませんが、一言でまとめるとXML形式 です。GPXでは軌跡やルート、ポイントを記録することができます。

実際に使われている例として登山者向けアプリの「YAMAP」、「ヤマレコ」や、フィットネスアプリとして有名な「Strava」等でアクティビティの軌跡を出力する際に、この形式を選択することができます。

以下にGPXファイルの中身の例を示します。富士登山を想定して作ったデータであり、架空のものですので実在しません（笑）

Subashiri.gpx

<?xml version='1.0' encoding='UTF-8'?>
<gpx creator="YAMAP - https://yamap.com" xsi:schemaLocation="http://www.topografix.com/GPX/1/1 http://www.topografix.com/GPX/1/1/gpx.xsd" version="1.1" xmlns="http://www.topografix.com/GPX/1/1" xmlns:xsi="http://www.w3.org/2001/XMLSchema-instance">
<metadata>
    <bounds maxlat="35.369134" maxlon="138.778500" minlat="35.366435" minlon="138.739171"/>
</metadata>
<trk>
    <name>須走口5合目〜8合目</name>
    <number>1</number>
    <trkseg>
        <trkpt lat="35.366435" lon="138.778500"><ele>1966.2</ele><time>2022-09-10T03:48:17Z</time></trkpt>
        <trkpt lat="35.366872" lon="138.778276"><ele>1969.2</ele><time>2022-09-10T03:50:50Z</time></trkpt>
        <trkpt lat="35.367560" lon="138.777829"><ele>1975.7</ele><time>2022-09-10T03:55:16Z</time></trkpt>
        .
        .
        .
        <trkpt lat="35.369134" lon="138.739171"><ele>3362.7</ele><time>2022-09-10T07:45:39Z</time></trkpt>
    </trkseg>
</trk>
</gpx>

②XYZ

 XYZはアルファベット最後の3文字…ではなくて、XYZの座標を記録したファイルのことです。単位はメートルで、ある位置（基準点）からの距離を表します。XYZファイルはあまり一般的ではない形式ですが、中身はいたって簡単でCSVの区切り文字がカンマから半角スペースになっているだけです。主にGISソフトで用いられるほか、シミュレーションなどをする際の3次元プロットの元データとしても用いられます。
以下に、XYZファイルの一例（先述したGPXファイルをXYZに変換したもの）を示します。

298175.47862468986 3915945.848919105 1966.2
298156.21326954296 3915994.7826505005 1969.2
298117.31068591564 3916072.0151253897 1975.7
294608.64727569465 3916326.204215371 3362.7

左からX、Y、Zの順に並んでいます。

実装方法

開発環境

• Python 3.10.8
• Visual Studio Code 1.74.0

使用ライブラリ

• gpxpy

 簡単に使うことのできるGPXファイルパーサーです。

• pyproj

 地理的情報を扱う際に、座標系・測地系の変換や方位角・距離の計算を簡単に行えるライブラリです。

緯度・経度・高さをUTM座標に変換

 まず先述したgpxpyを使って、GPXファイルから緯度・経度・高さのみを読み込みます。下記リンクを参考にさせていただきました。

以下の例は読み込んだ緯度・経度・高さをCSVとして書き出すプログラムです。

gpx2utm.py

import csv
import gpxpy

with (open("Subashiri.gpx",mode="r") as gpx_file, open("Subashiri.csv",mode="w") as csv_file):
    coordinate_list = []
    gpx = gpxpy.parse(gpx_file)
    for trk in gpx.tracks:
        for seg in trk.segments:
            for pt in seg.points:
                coordinate_list.append([pt.latitude, pt.longitude, pt.elevation])

 次に、読み込んだ緯度・経度・高さをUTM座標に変換します。UTM（ユニバーサル横メルカトル）については下記リンクを参照ください。

 UTM座標系では経度によって「UTMゾーン」というものが決められています。計算式で求めることもできるのですが、ここではUTMゾーンがわかっている前提で実装します。下記リンクを参考にさせていただきました。

utm2csv.py

import csv
from pyproj import Proj

utm_zone = input("UTMゾーンの値を入力（半角数字）：")
# WGS84楕円体に準拠する場合、ellps="WGS84"とすればよい。
utm_conv = Proj(proj="utm", zone=utm_zone, ellps="GRS80")
out_list = []

with (open("Subashiri.csv", mode="r") as input_file, open("Subashiri_utm.csv", mode="w") as out_file):
    reader = csv.reader(input_file)
    for row in reader:
        coordinate_list = [float(i) for i in row]
        lat = coordinate_list[0]
        lon = coordinate_list[1]
        ele = coordinate_list[2]
        # 経度は東西方向なのでX方向、緯度は南北方向なのでY方向のデータになる
        utmx, utmy = utm_conv(lon, lat)

        # 南半球ならオフセットする
        if lat < 0:
            utmy = utmy + 10000000
        else:
            pass

        out_list.append([utmx, utmy, ele])
    writer = csv.writer(out_file)
    writer.writerows(out_list)

 ここでは準拠する楕円体をGRS80楕円体にしましたが、WGS84楕円体にすることも可能です。

XYZファイルへの書き出し

 XYZファイルは拡張子こそ特殊ですが、区切り文字がスペースになったCSVのようなものです。そこで標準ライブラリのCSVを用いて書き出しを行い、保存するファイルの拡張子を.xyzにする裏技（？）を使います。

実装例

 完成したコードを以下に示します。

gpx2xyz.py

# GPXファイルから緯度経度標高を抽出してXYZファイルに書き出すプログラム
import csv
import gpxpy
from pyproj import Proj

def from_gpx_to_xyz():
    file_name = input("GPXファイルの名称を入力（拡張子除く）：")
    utm_zone = input("UTMゾーンの値を入力（半角数字）：")
    utm_conv = Proj(proj="utm", zone=utm_zone, ellps="GRS80")
    xyz_list = []

    with (open(f"{file_name}.gpx",mode="r") as gpx_file, open(f"{file_name}_from_gpx.xyz",mode="w") as xyz_file):
        gpx = gpxpy.parse(gpx_file)
        for trk in gpx.tracks:
            for seg in trk.segments:
                for pt in seg.points:
        # 経度は東西方向なのでX方向、緯度は南北方向なのでY方向のデータになる
                    utmx, utmy = utm_conv(pt.longitude, pt.latitude)
                    # 南半球ならオフセットする
                    if pt.latitude < 0:
                        utmy = utmy + 10000000
                    else:
                        pass
                    xyz_list.append([utmx, utmy, pt.elevation])
        writer = csv.writer(xyz_file, delimiter=" ")
        writer.writerows(xyz_list)

if __name__ == "__main__":
    from_gpx_to_xyz()

 ファイルの名称とUTMゾーンを指定するだけで、GPXファイルからXYZファイルに一発変換するツールを作ることができました。

おわりに

 実装例のコード、およびテストに用いたGPXファイルとXYZファイルはGithubで公開しています。