DJI製ドローン撮影写真の高度の補正

概要

• 以下の「融合高度」の説明の通り、DJI製ドローンの高度情報には問題があります。
• 最近ではネットワークRTKを利用した機体がメインになりつつありますが、旧来と現状では以下の融合高度を考慮しなければなりません。
• 最近の機体でも、RTKによって得られる高度は「地球楕円体高」なので、標高を前提とした測量関連やSfMソフトウェアと相性が悪い。
• このようなわけで、高度の補正ツールを自作して活用しています。

参考サイト

• 「楕円体高」「標高」「ジオイド高」などの用語の意味が分からない人は以下のようなサイトを参照してください。初歩的な説明は省略します。
• 地球中心座標系とは
  https://www.jesea.co.jp/about/coordinate/
• ジオイドを知る
  https://www.gsi.go.jp/buturisokuchi/grageo_geoid.html#:~:text=%E5%9B%BD%E5%9C%9F%E5%9C%B0%E7%90%86%E9%99%A2%E3%81%A7%E3%81%AF%E3%80%81%E9%87%8D%E5%8A%9B,%E6%B1%82%E3%82%81%E3%82%8B%E3%81%93%E3%81%A8%E3%81%8C%E3%81%A7%E3%81%8D%E3%81%BE%E3%81%99%E3%80%82

融合高度についてのDJI社の回答(要約)

• RTK対応のドローン機体の場合、写真のExif情報内の高度は楕円体高。
• RTKを利用しないドローン機体の場合、融合高度を記録。

融合高度についての説明

• カメラが機体内部のフライトコントローラーモジュールの「融合高度」を読み取り、画像のEXIF情報内への書込み。（絶対高度）
  ＊「融合高度」は、「気圧計」及び「GPS」により取得された高度データを一定のアルゴリズムにもとづいて算出された数値。

①

• GPSによる高度データは、GPS受信機が衛星からの信号を受信し、同一時刻の受信機と複数の衛星間の擬似距離を得られた上、融合アルゴリズムにより測定された高度値。
• 精度が限られているので、100m以上の誤差が存在する場合がある。（マイナスの値になる場合も）
• 誤差の構成部分として、GPS衛星による誤差(衛星システムの位置、角度/高度の動態変化)、信号伝送による誤差、受信機による誤差が考えられる。

②

• 気圧センサーは風により、誤差が発生する可能性がある。
  例：横風が吹いている場合、気圧計付近の気圧値が急激に変化することにより気圧計による高度誤差が発生する。
  上記以外気温による影響が発生する可能性もある。
• 結果的に高度測定精度に影響を与える可能性がある要素は非常に多い状況。

説明

• Mavic、Inspire2やPhantom4シリーズを含めた機種で、EXIF情報内に記録される高度は大きな誤差が発生している状況は確認済み。しかし、現段階での改善が非常に難しい状況。
• これらの原因により融合高度は参考値としてのみ利用した方がいい。

より正確な高度情報を求める場合

• 相対高度：画像gps情報内の相対高度(Relative Altitude)を確認。
  ＊相対高度(対地高度)は離陸点との高度差を相対的に正確に測れる。
• 絶対高度：GPS経緯度から撮影地点の高度を調べてた上で、相対高度と両者を加算した数値が実際の絶対高度を得られる。

対処方法の検討

• RTK以外の機体がEXIFに書き込む GPSAltitude(絶対高度) は融合高度なので誤差が100mなど大きく、この値を利用することはあまり良くない。
• RTK機体がEXIFに書き込む GPSAltitude(絶対高度) は地球楕円体高で標高ではないので、標高に変換して利用する必要がある。
• SfMソフトウェアが利用する高度はEXIFの GPSAltitude 。
  ＊私がOpenDroneMapで確認した情報です。
• DJI製機体はEXIF標準の GPSAltitude(絶対高度) とは別に RelativeAltitude(相対高度) を拡張属性として記録する。
  ＊私が調べた限り、Parrot製ドローンはこの拡張属性を使用していません。
• 以下の図から２つの対処方法を検討します。

国土地理院の標高による補正

• DJIのExif拡張属性 RelativeAltitude + 国土地理院の標高 = GPSAltitudeとして書き込む。
• 海抜高度=AMSL(Above Mean Sea Level)から飛行中ドローンの高度を求める方法です。
• DJI機が記録する相対高度(RelativeAltitude)を信頼するという懸念点(気圧に左右される)がありますが、古いRTKではない機体では少しはまともな値にできると思っています。
• 対象はDJIのRTKではない機体です。Inspire1/2、Mavic、PhantomのRTKなし、MatriceのRTKなし、など。

国土地理院のジオイドによる補正

• Exif属性のGPSAltitude - 国土地理院のジオイド高 = GPSAltitudeとして書き込む。
• GPSの正確な楕円体高から海水面下の高度=ジオイド高を引くことにより、標高=AMSLから飛行中ドローンの高度を求める方法です。
• 最近のネットワークRTK対応機体の場合、こちらの方が正確なので、望ましいと思います。
• 対象はRTK対応の機体。Matrice 210 V2 RTK、Matrice 300 RTK、など。

国土地理院APIのサンプル解説

• あくまで「サンプル」なので、本筋に関係ないコードは全て省略してあります。
• 実際に使用する場合は、よりエラーや例外処理を入れて信頼できるソースコードにしてください。

国土地理院の標高による補正処理(Pythonサンプル)

# 説明:
# DJIのExif拡張属性 RelativeAltitude + 国土地理院の標高 = GPSAltitudeとして書き込む。

from logging import getLogger, StreamHandler, DEBUG
import json
import sys
import glob
import os, shutil
import exiftool
import requests
import time

def main():
    # ログの設定
    # 引数のチェック
    # 元画像ファイルのバックアップ
    # などの処理。でもサンプルなので省略。

    # json出力用の辞書/json
    photoInfoList = {}

    # photoList に個々の写真のパスが入っている。
    # PyExifToolでXMP(拡張属性)の相対高度を抜き出し、Exifの高度に書き込む
    for i in range(len(photoList)):
        with exiftool.ExifTool() as et:
            metadata = et.get_metadata(photoList[i])

            # lat/lonからelevation(標高)を取得する
            lat = metadata['EXIF:GPSLatitude']
            lon = metadata['EXIF:GPSLongitude']
            # 国土地理院に問い合わせる。 https://www.gis-py.com/entry/elevation-api
            url = "http://cyberjapandata2.gsi.go.jp/general/dem/scripts/getelevation.php" "?lon=%s&lat=%s&outtype=%s" %(lon, lat, "JSON")
            resp = ""
            data = ""
            while True:
                resp = requests.get(url, timeout=30)
                data = resp.json()
                if type(data["elevation"]) is list:    # 通常戻り値jsonのelevationはint/float、時々list([])の場合がある。
                    print('data["elevation"] => ', data["elevation"])
                    print('type(data["elevation"]) => ', type(data["elevation"]))
                    time.sleep(1)
                    print('[debug]type(data["elevation"]) is not int or float => Retry.')
                    continue
                break

            if isinstance(metadata['XMP:RelativeAltitude'], float) or isinstance(metadata['XMP:RelativeAltitude'], int):
                # 相対高度と標高を足す
                alt_asml = metadata['XMP:RelativeAltitude'] + data["elevation"]
                alt_string = str(alt_asml)  # 文字列型に変換
            else:
                # 相対高度と標高を足す
                tmp_string = metadata['XMP:RelativeAltitude'][1:]   # 最初の"+"をカット
                alt_asml = float(tmp_string) + data["elevation"] 
                alt_string = str(alt_asml)  # 文字列型に変換
            write_string = "-GPSAltitude=" + alt_string + " -overwrite_original "
            et.execute(bytes(write_string, encoding="shift_jis"), bytes(photoList[i], encoding="shift_jis"))

            # ログ出力用の辞書/jsonを作成
            # 省略

if __name__ == '__main__':
    try:
        main()
    except Exception as ex:
        print('Exception:')
        print(ex)

国土地理院のジオイドによる補正処理(Pythonサンプル)

# 説明:
# Exif属性のGPSAltitude - 国土地理院のジオイド高 = GPSAltitudeとして書き込む。

from logging import getLogger, StreamHandler, DEBUG
import json
import sys
import glob
import os, shutil
import exiftool
import requests
import time

def main():
    # ログの設定
    # 引数のチェック
    # 元画像ファイルのバックアップ
    # などの処理。でもサンプルなので省略。

    # json出力用の辞書/json
    photoInfoList = {}

    # photoList に個々の写真のパスが入っている。
    # PyExifToolでXMPの相対高度を抜き出し、Exifの高度に書き込む
    for i in range(len(photoList)):
        with exiftool.ExifTool() as et:
            metadata = et.get_metadata(photoList[i])

            # lat/lonからgeoid高を取得する
            lat = metadata['EXIF:GPSLatitude']
            lon = metadata['EXIF:GPSLongitude']
            # 国土地理院に問い合わせる。
            url = "http://vldb.gsi.go.jp/sokuchi/surveycalc/geoid/calcgh/cgi/geoidcalc.pl" "?outputType=%s&latitude=%s&longitude=%s" %("json", lat, lon)
            resp = ""
            data = ""
            while True:
                resp = requests.get(url, timeout=30)
                if "50" in str(resp):
                    print("Response Error => Retry.")
                    time.sleep(1)
                    continue
                data = resp.json()
                if type(data["OutputData"]["geoidHeight"]) is list:    # 通常戻り値jsonのgeoidHeightはint/float、時々list([])の場合がある。
                    print('data["OutputData"]["geoidHeight"] => ', data["OutputData"]["geoidHeight"])
                    print('type(data["OutputData"]["geoidHeight"]) => ', type(data["OutputData"]["geoidHeight"]))
                    time.sleep(1)
                    print('[debug]type(data["OutputData"]["geoidHeight"]) is not int or float => Retry.')
                    continue
                break

            alt_asml = float(metadata['EXIF:GPSAltitude']) - float(data["OutputData"]["geoidHeight"]) # GPSAltitudeから国土地理院のジオイド高を引く。
            alt_string = str(alt_asml)  # 文字列型に変換。
            write_string = "-GPSAltitude=" + alt_string + " -overwrite_original "
            et.execute(bytes(write_string, encoding="shift_jis"), bytes(photoList[i], encoding="shift_jis"))


            # ログ出力用の辞書/jsonを作成
            # 省略

    logger.debug(" ")

if __name__ == '__main__':
    try:
        main()
    except Exception as ex:
        print('Exception:')
        print(ex)

国土地理院APIサーバの注意

• GPS座標から標高を求めるAPIサーバ。
  参考サイト：Pythonで国土地理院の標高API を使ってみよう！
  https://www.gis-py.com/entry/elevation-api
• GPS座標からジオイド高を求めるAPIサーバ。
  標高APIと基本的に同じです。GPS座標を渡すと、その地点のジオイド高が返ってきます。
• でも、利用するにあたっては注意があります！

国土地理院APIの動作

• とにかく HTTP Status Code 200 OK を返す。
• でもレスポンスのjsonに標高やジオイド高の数値が入ってない(空)の場合がある！
• おそらく裏方でデーらベースにクエリを投げて値を取得していると思うのですが、なんとクエリが間に合わなくて値が空でもとにかくjsonを返すという厄介な動作です。
• (心の声：私はずるい実装だと思います。クエリが間に合わないで値が無いなら「Value cannot get」とかエラーにすべきです。)
• 気付くのに結構時間かかりました。通信ダンプしたりレスポンスの中身を全部精査して調べてやっと分かりました。
• このようなわけで、HTTP 200 OKでもレスポンスの中身を見てjsonの値が入ってなければ再び同じGPS座標で問い合わせを繰り返しています。
• 詳細はサンプルソースコードを見れば分かると思います。

現状とこれから

• 実際に使用する場合は、これよりももっと機能追加して、DJIの機体種別をEXIFから判別して処理するなど改善したいと思っています。
• ネットワークRTKの機体も増えてきているので、対応したいです。
• 最近のドローン業界はDJI一辺倒ではなくて、他の種類のPixhawk系や国産ドローンなども増えてきているので、更に調査して対応していきたいと思っています。
• それでは。興味のある方は是非！