PLATEAUのCityGMLを鼻血出してQGISで表示する

【2021/9/21 追記】
この記事は2021/3下旬での検証結果で、2021/9/20に @mugwort_rc さんがGDAL本家(!)にFIXをPUSHしてくださり、時期GDALバージョン(3.3.3あたり？）で改善され、それを取り込むちょっと先のQGISでもなおる見込みです（そのころにはこの記事はふつうのXY交換のやりかた紹介記事に転生します😂 ）。 @mugwort_rc さん本当にありがとうございます！！！

すばらしい修正記事はこちら↓

はじめに

「こんなに生やしてどうするのォォ！」でおなじみの、日本全国の3D都市モデルの整備・オープンデータ化プロジェクトPLATEAUが話題になっています。
この記事を書いている前日()にG空間情報センターからCityGML形式データなどがCCBY4.0準拠のオープンライセンスでダウンロードできるようになりました。わたしはCode for Japanのオンラインハッカソンイベント、Social Hack Dayのチーム4「PLATEAU/3D都市モデルを触ってみよう」でご一緒した@tossetoさん@hal_skさんたちとともにこのPLATEAUのCityGMLデータに半日しがみついておりました(死闘記録はこちら)。その結果、チームの皆さんの助けもあって無事QGIS上に表示することができましたので、その結果を紹介します。慣れたらLOD1（豆腐建物）なら５分くらいでQGISに載せれますよ。
もちろんもっとスマートでスッキリできる処理方法もあるかもしれませんが、まずはできましたということで♨

こんなかたに

• CityGMLを地図上でさくっと確認したいQGISな方
• 国交省PJに参加したものの鬼データもらって頭を抱える役所の方
• PLATEAUのスペル覚えるのがむつかしいひと

やりかた

注：けっこう負荷の高い作業をさせますので、なるべくスペックの高いPCで試してくださいね

データのダウンロードと解凍

G空間情報センターダウンロードページの3D都市モデル（Project PLATEAU）東京都23区（CityGML 2020年度）から、閲覧したいエリアの地域基準メッシュ番号に対応したデータを選択します。データに対応する位置は同じページ内の「東京都23区構築範囲図」で閲覧できますし、よくわかんなかったらこのあたりの索引地図から引いても良いです。
なお、範囲図のピンクのエリアであるLOD2レベルの入っている地域はデータが数Gバイトになるのでご注意下さい。また、本項ではLOD2レベルのデータ確認検証は行っていません（たぶんドライバの関係でむづかしいかも）、ごめんなさい。

では、目黒区あたりの533935をDLします。
解凍するとbldg/tran/dem/luse…などのzipがでてきます。

これらは、東京23区における3D都市モデルのための拡張製品仕様によると、下記のような分類になっています。

bldg 建築物、建築物部分、建築物付属物、及びこれらの境界面 
tran 道路、通路
luse 土地利用
urf 都市計画区域、区域区分、地域地区
dem 地形(起伏)
fld 洪水浸水想定区域
tnm 津波浸水想定 
lsld 土砂災害警戒区域 
brid 橋梁
frn 設置物

全部試して「俺のALL目黒区inQGIS」を作ってみるのもなかなかチャレンジングではありますが、まずは建物である「bldg.zip」を更に解凍してみると、大量のgmlファイルが生成されます。
この一つ一つには、更に詳細な地域基準メッシュ番号ごとの建物CityGMLファイルとなっています。

また、appearanceフォルダにはこんな建物のおもしろテクスチャ（表面画像）が入っていますが今回は参考までにくらいでみておきます。

どうやって見るんですか

ここで大きな壁にぶつかります。
メンバーたちはどうすればQGISでCityGMLを開けるか数時間悩んでいました。
CityGMLからCityJSONに変換してもQGISに出ない。有料変換サービスに$7払ったらレシートとエラーログしか返ってこなくてダメ、サーバー立てて3DTiles形式を読む？　Blender経由でFBX形式からなんかに変換する？…
散々悩んだ末、ある瞬間「あれ、これ、QGISドラッグ＆ドロップで読めるやん」ということに気がついたのです。なんという灯台モトクラシー！こんなワクワクあるから検証はやめられない。なお、QGISのバージョンは3.16以上じゃないと読めませんのでご注意を（わたし古いバージョンで試してた…）

これは瀬戸さんが調べてくださったのですが、QGISのバージョン3.16の更新履歴にList available GDAL vector driversとありました。なお、バージョン情報をみると、QGISに搭載されているGDALが3.0.4から3.1.4にバージョンアップされています。GDALは3.1.1から単純なGMLにも対応とあるため、その恩恵がQGISにもやってきたというわけです。

なんか出たけど変

では、とりあえず目黒区の53393585あたりをドラッグ＆ドロップでQGISバージョン3.16以上で開いてみます。なお、appearanceが一緒だとちょっと読み込みに数分かかるのでご注意を。
しばらくするとデータどっちを開くか聞いてくるので、とりあえずBulidingを選択してOKを押します。

投影法の設定などについて聞いてくる場合がありますが基本EPSG:6697でも大丈夫そうなので、OKを押します。
その後ポリゴンレイヤにズームさせるとなんか建物矩形ぽいものが！

しかしここでまた別の問題が。
この位置がどこなのかを確認するためにOpenStreetMapや地理院地図様などの正しい位置が入ったタイル地図を読み込んでこようとするとハングアップしてしまいます。また、正しい市町村ポリゴンなどとも合致しません。ここでまた1時間ほど悩んでいましたがあることに気がついたのです。
XY座標がQGISの表示と逆になっているのです。
原因は元データの不備よりかは、GIS/CADのXYに対する考え方の違いとか仕様の差によるものかもしれませんが、とにかく位置が合った絵を作りたいので（）がっちり探っていません。

読み方や座標変換のPROJ4をなんかすればいいのかもですが今、今回はデータをFIXする方法をとってみました

矯正してみる

この点の"XY-YX問題"は過去何回か食らっていまして、アジャスト的なベクトル変換などをかければ対応できます。
その前に生のGMLだといろいろエラーが出るので、一旦 .gpkg（GeoPackage形式）などに変換しておきます。
SHPに変換もできなくもないのですが、今回の国交省さんCityGMLのワイルドな長くて日本語満載の属性フィールド名ですとSHPがお腹を下すので、*.gpkgとかのほうがおすすめです。

読んできたCityGMLのレイヤを選択後右クリックで地物の保存→任意のフォルダにGeoPackage形式で保存します。

さて、材料が揃ったところで、XY-YX変換を行います。以前のQGISのVer２時代は専用のプラグインがあったようですが、いまでは基本機能のプロセッシングツールボックス内に変換ツールが入っています、開発者のみなさんありがとうございます。
プロセッシングメニュー→プロセッシングツールボックスを表示させてswapと🔍検索窓に打ち込むと座標のXY入れ替えという機能が出てくるので、起動します。

右のヘルプには

このアルゴリズムは、入力ジオメトリの頂点のXとYを入れ替えます。誤って緯度と経度を逆に入力したジオメトリを修復するために使用できます

とありますが、あちらさんもご事情があるでしょうし😅温かい心持ちで、入力レイヤに先程のgpkgを、出力レイヤはとりあえず[一時レイヤを作成]でもよいのでOKを押します（慣れたら「...」から*.gpkgに保存してもよいかと）

変換されたレイヤにズーム後、OpenStreetMapや地理院地図様の地図タイルを敷いてみると、おお！ハマりました。
サンプルのgeojsonはここに置いときます

投影法をEPSG：3857などにして、QGIS2Threejsプラグインなどを使って、建物高さを"measuredHeight"*1.5とかにするといい感じで３Dにしてくれます。今回はLOD1レベルのデータなので、「ちょっとかっこいい豆腐」ではありますが、これだけでもさまざまな活用が見込めそうですね。

現地で測量してくださった方、様々な点群データからCityGMLを生成してくれた方、かっちょいいプラットフォームこさえてくれた方、オープンデータのために各方面への調整に奔走してくれた方、その他たくさんの方の苦労がこのデータに詰まっているようにも思います、お疲れさまでした、ありがとうございます。うやうやしく活用させていただいて、これらの成果がデータを軸とした社会貢献の下地になることを願ってやみません。

おわりに

なお、QGISおじさんたちがヒーヒー言ってGMLともみくちゃ泥まみれになっている間にメンバーのmamixさんはFBX形式から10分くらいでさっさとUnityに取り込んで猫オブジェクト立てたりして遊んでました…orz