PLATEAUのLOD2有りのCityGMLがWin版のQGISでは読めるのに、Mac版では読めないので調査！

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/203944/e124ae45-7565-419a-cf91-c4540b6945af.png)

みなさん！PLATEAUで**GIS**やってますかー！！！！

PLATEAUに利用されているCityGMLは地理空間情報（GIS）データなのですが、2Dのデータと3Dのデータを混在させることができるため、今まで接点がなかった分野にも広がっている上、整備済みの都市もどんどん増えています。

また、先人が「[鼻血を出しながら](https://qiita.com/Yfuruchin/items/aa3a847db26c7f378d6e)」調査をしてくれたおかげでLOD1の建物データがQGISでもサクッと見れるようになっていますね！

どの段階から読み込めるようになっていたのかは検証していないですが、win版のQGIS v3.26の時点ではLOD2を含むCityGMLもドラッグ&ドロップでサクッと表示できるようになっています！素晴らしい！！！
（※ただし、言うことを聞いてくれるのは建物ポリゴンのみのようで、他の地物型はすんなり読み込めません）

が、**超絶悲報**です！！！！

Mac版のQGISではv3.28でも読めないし、なんなら最新のナイトリーリリースでも読み込めません！！！

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/203944/117343d4-4035-ae01-3d6e-74b476a162cb.png)
（**何も読み込まれないしこの後、背景地図も表示できなくなってしまう**）

さらに、LOD1しかないCityGMLを開いても、XYが反転され縦長に表示されてしまいます！
なんでやねん！！！つらい！！！

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/203944/350e7757-5381-b086-a762-66f3b43e24ff.png)

## なんでmac版だけうまく開けないのか

### XYが反転する問題

そもそも、QGISでデータを開く際には大部分を[GDAL](https://gdal.org/index.html)という大規模なGISデータ変換ツールに依存しています。

このため、「**データが読み込めない！変な位置に表示される！**」などの理由は大抵GDALに原因があります。

Win版のQGISの概要を見てみると、このように表示されますね。

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/203944/b18b1531-8bfc-85a1-e2e9-d96cb9a89596.png)

QGISv3.26.0のGDLA/OGRバージョンが3.5.0となっています。

が！！！！！！！

Mac版の概要を見てみましょう。

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/203944/a4dfcc60-2728-8a74-0ad4-aff669e6e8c3.png)

**QGISはv3.28.2**ですが**GDAL/OGRは3.3.2**であることがわかります。

原因これやんけ！！！

冒頭の鼻血記事の中で紹介されていた[QGISでProject PLATEAUのCityGMLの緯度経度が入れ替わる問題を修正した件](https://zenn.dev/mugwort_rc/articles/34a427a9394580)での問題の修正はGDAL**v3.3.3**で取り込まれましたが、Mac版の方は**3.3.2**が利用されているため、ギリッギリ取り込まれていませんね。

つまり、2022/12/23現在では、Mac版ではCityGMLがQGISで表示できないので、「**大人しくWin版を使え**」ということになりますね。

つらE。

まぁとはいえひとまず「XYが反転する問題」については原因が発覚し、待てばいづれは解決する問題だということがわかりました。

が、LOD2が全く表示されない原因はまだわかりません。

### LOD2が全く表示されない

Win版では表示されるようになっているので、おそらくGDALv3.3.2~v3.5.0の間に何かしら修正がなされたのは事実でしょう。

諸々割愛しますが、調査のためにGDALのv3.3.2をインストールしました。

```bash
$ ogr2ogr --version
GDAL 3.3.2, released 2021/09/01
```

GDALをインストールすると、ベクトル形式のGISデータの変換ツールである「ogr2ogr」も一緒にインストールされます。

これを利用してGIS界隈で一般的に利用されるGeoJSON形式への変換を試みると…

```python
❯ ogr2ogr -f "GeoJSON" ./53394567_bldg_6697_2_op.geojson ./53394567_bldg_6697_2_op.gml
Warning 1: Unrecognized geometry type : 1015
ERROR 6: OGR geometry type unsupported as a GeoJSON geometry detected. Feature gets NULL geometry assigned.
ERROR 6: OGR geometry type unsupported as a GeoJSON geometry detected. Feature gets NULL geometry assigned.
ERROR 6: OGR geometry type unsupported as a GeoJSON geometry detected. Feature gets NULL geometry assigned.
ERROR 6: OGR geometry type unsupported as a GeoJSON geometry detected. Feature gets NULL geometry assigned.
ERROR 6: OGR geometry type unsupported as a GeoJSON geometry detected. Feature gets NULL geometry assigned.
ERROR 6: OGR geometry type unsupported as a GeoJSON geometry detected. Feature gets NULL geometry assigned.
ERROR 6: OGR geometry type unsupported as a GeoJSON geometry detected. Feature gets NULL geometry assigned.
...
```

このように全ての地物が読み取れず座標がNULLになってしまいました。

気になる点は1行目の`Warning 1: Unrecognized geometry type : 1015`の部分ですね。

GDALのバージョンを3.5.0に上げてみましょう。

```bash
❯ ogr2ogr --version                                                                   
GDAL 3.5.0, released 2022/05/10
```

先ほどと同じコマンドを実行してみます。

```bash
❯ ogr2ogr -f "GeoJSON" ./53394567_bldg_6697_2_op.geojson ./53394567_bldg_6697_2_op.gml
ERROR 6: OGR geometry type unsupported as a GeoJSON geometry detected. Feature gets NULL geometry assigned.
ERROR 6: OGR geometry type unsupported as a GeoJSON geometry detected. Feature gets NULL geometry assigned.
ERROR 6: OGR geometry type unsupported as a GeoJSON geometry detected. Feature gets NULL geometry assigned.
ERROR 6: OGR geometry type unsupported as a GeoJSON geometry detected. Feature gets NULL geometry assigned.
ERROR 6: OGR geometry type unsupported as a GeoJSON geometry detected. Feature gets NULL geometry assigned.
ERROR 6: OGR geometry type unsupported as a GeoJSON geometry detected. Feature gets NULL geometry assigned.
ERROR 6: OGR geometry type unsupported as a GeoJSON geometry detected. Feature gets NULL geometry assigned.
ERROR 6: OGR geometry type unsupported as a GeoJSON geometry detected. Feature gets NULL geometry assigned.
ERROR 6: OGR geometry type unsupported as a GeoJSON geometry detected. Feature gets NULL geometry assigned.
ERROR 6: OGR geometry type unsupported as a GeoJSON geometry detected. Feature gets NULL geometry assigned.
ERROR 6: OGR geometry type unsupported as a GeoJSON geometry detected. Feature gets NULL geometry assigned.
ERROR 6: OGR geometry type unsupported as a GeoJSON geometry detected. Feature gets NULL geometry assigned.
...
```

同じエラーのように見えますが、`Warning 1: Unrecognized geometry type : 1015`が消えており、ogr2ogrは読み取れているが、GeoJSONとして割り当てることができていないように見えます。

もう一つ有名なGISデータのフォーマットとして、GeoPackageというフォーマットも存在します。
こちらはバイナリのため、GeoJSONよりも軽量かつ高速です。

GeoJSONのジオメトリとしてサポートされていない、というエラーが出ているので、GeoPackageへの変換を試みてみましょう。

```bash
❯ ogr2ogr -f "GPKG" ./53394567_bldg_6697_2_op.gpkg ./53394567_bldg_6697_2_op.gml      
Warning 1: The output driver does not natively support IntegerList type for field key. Misconversion can happen. -mapFieldType can be used to control field type conversion.
Warning 1: The output driver does not natively support IntegerList type for field codeValue. Misconversion can happen. -mapFieldType can be used to control field type conversion.
Warning 1: Registering non-standard gpkg_geom_POLYHEDRALSURFACE extension
Warning 1: The output driver does not natively support StringList type for field imageURI. Misconversion can happen. -mapFieldType can be used to control field type conversion.
Warning 1: The output driver does not natively support StringList type for field mimeType. Misconversion can happen. -mapFieldType can be used to control field type conversion.
Warning 1: The output driver does not natively support StringList type for field textureCoordinates. Misconversion can happen. -mapFieldType can be used to control field type conversion.
ERROR 1: failed to execute insert : NOT NULL constraint failed: Appearance.gml_id
ERROR 1: Unable to write feature 0 from layer Appearance.
ERROR 1: Terminating translation prematurely after failed
translation of layer Appearance (use -skipfailures to skip errors)
ERROR 1: sqlite3_exec(CREATE TRIGGER "trigger_insert_feature_count_Building" AFTER INSERT ON "Building" BEGIN UPDATE gpkg_ogr_contents SET feature_count = feature_count + 1 WHERE lower(table_name) = lower('Building'); END;) failed: no such table: main.Building
ERROR 1: sqlite3_exec(CREATE TRIGGER "trigger_delete_feature_count_Building" AFTER DELETE ON "Building" BEGIN UPDATE gpkg_ogr_contents SET feature_count = feature_count - 1 WHERE lower(table_name) = lower('Building'); END;) failed: no such table: main.Building
ERROR 1: sqlite3_exec(CREATE TRIGGER "trigger_insert_feature_count_Appearance" AFTER INSERT ON "Appearance" BEGIN UPDATE gpkg_ogr_contents SET feature_count = feature_count + 1 WHERE lower(table_name) = lower('Appearance'); END;) failed: no such table: main.Appearance
ERROR 1: sqlite3_exec(CREATE TRIGGER "trigger_delete_feature_count_Appearance" AFTER DELETE ON "Appearance" BEGIN UPDATE gpkg_ogr_contents SET feature_count = feature_count - 1 WHERE lower(table_name) = lower('Appearance'); END;) failed: no such table: main.Appearance
```

今度は全く別のエラーが出てきてしまいましたね…

ただ、諦めるのはまだ早いです。

エラーで停止している箇所をよく見るとLOD2特有のAppearance（テクスチャ）が読み込めずに怒られていますね。
さらに`translation of layer Appearance (use -skipfailures to skip errors)
`とも書かれています。

GISでデータを扱う以上はテクスチャの情報は不要ですので、skipしてしまいましょう。

すると、Appearanceに関するエラーは発生していますが、目的の`53394567_bldg_6697_2_op.gpkg`が出力されていますね！

```bash
❯ ogr2ogr -f "GPKG" ./53394567_bldg_6697_2_op.gpkg ./53394567_bldg_6697_2_op.gml -skipfailures      
Warning 1: The output driver does not natively support IntegerList type for field key. Misconversion can happen. -mapFieldType can be used to control field type conversion.
Warning 1: The output driver does not natively support IntegerList type for field codeValue. Misconversion can happen. -mapFieldType can be used to control field type conversion.
Warning 1: Registering non-standard gpkg_geom_POLYHEDRALSURFACE extension
Warning 1: The output driver does not natively support StringList type for field imageURI. Misconversion can happen. -mapFieldType can be used to control field type conversion.
Warning 1: The output driver does not natively support StringList type for field mimeType. Misconversion can happen. -mapFieldType can be used to control field type conversion.
Warning 1: The output driver does not natively support StringList type for field textureCoordinates. Misconversion can happen. -mapFieldType can be used to control field type conversion.
ERROR 1: failed to execute insert : NOT NULL constraint failed: Appearance.gml_id
```

mac版のQGISで開いても地物がちゃんと表示できているのがわかります。

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/203944/bb98c430-399b-91c1-477d-18d296bb211f.png)

Mac版のQGISではGDALのバージョンが3.3.2であるため、`Unrecognized geometry type : 1015`というのが理由（つまりGDALの地物として読み込めていない）で表示することができず、逆にWin版ではGDALのバージョンが3.5.0であるため地物型としてはサポートされているため表示可能、であるように見えます。

具体的にGDALの実装を確認して直接的な原因を解明〜というわけにはいかなかったのですが、少なくとも**GDALのバージョンが3.5.0以上であればLOD2の地物が入っていてもちゃんと読み込める**ことがわかりました。

（※ただし、読み込めたとしても、QGISで表示しているのはあくまでLOD0、つまり平面形状のみです。LOD2がCityGML上に格納されていたとしても、それらを無視してLOD0のデータを読み取って、表示できるようになった、というのが今回の調査結果になります。）

ということで「LOD2が全く表示されない」についても原因が発覚し、待てばいづれは解決する問題だということがわかりました。

## Pythonでも触ってみた

ちなみに、僕はよくPythonでGISデータを解析するんですが、PLATEAUのCityGMLを読み取ることはできるんでしょうか？

環境構築などは割愛させていただきますが、ちょっと触ってみました。

Pythonで地理的ベクトルデータを操作する際にはほぼ必須、とも言えるGeoPandasを利用してデータを読み込んでみましょう。

```python
import geopandas as gpd

gdf = gpd.read_file("./53394567_bldg_6697_2_op_no_lod2.gml")
gdf
```

結果は…

```python
UnsupportedGeometryTypeError              Traceback (most recent call last)
Cell In[15], line 1
----> 1 gdf = gpd.read_file("./53394567_bldg_6697_2_op_no_lod2.gml")
      2 gdf

File ~/miniforge3/envs/plateau/lib/python3.9/site-packages/geopandas/io/file.py:259, in _read_file(filename, bbox, mask, rows, engine, **kwargs)
    256     path_or_bytes = filename
    258 if engine == "fiona":
--> 259     return _read_file_fiona(
    260         path_or_bytes, from_bytes, bbox=bbox, mask=mask, rows=rows, **kwargs
    261     )
    262 elif engine == "pyogrio":
    263     return _read_file_pyogrio(
    264         path_or_bytes, bbox=bbox, mask=mask, rows=rows, **kwargs
    265     )

File ~/miniforge3/envs/plateau/lib/python3.9/site-packages/geopandas/io/file.py:351, in _read_file_fiona(path_or_bytes, from_bytes, bbox, mask, rows, where, **kwargs)
    349     f_filt = features
    350 # get list of columns
--> 351 columns = list(features.schema["properties"])
    352 datetime_fields = [
    353     k for (k, v) in features.schema["properties"].items() if v == "datetime"
    354 ]
    355 if kwargs.get("ignore_geometry", False):
...
File fiona/ogrext.pyx:719, in fiona.ogrext.Session.get_schema()

File fiona/_geometry.pyx:81, in fiona._geometry.normalize_geometry_type_code()

UnsupportedGeometryTypeError: 1015
```

あかんやんけ！

と思いましたが、最終行に見慣れた数値が…

`UnsupportedGeometryTypeError: 1015`

先ほどのGDAL v3.3.2で変換を試みた時と同じ数値が書かれているようです。

それも当然の話で、GeoPandasではベクトルデータの変換をfionaというライブラリに依存していますが、fionaはGDALに依存しています。

なので、当然GDALが対応していなければデータの変換はできません。

GeoPandasの最新版（安定版）はv0.12.2で、依存しているfionaのバージョンはv1.8.21になります。

```bash
$ conda list | grep geopandas
geopandas                 0.12.2             pyhd8ed1ab_0    conda-forge
geopandas-base            0.12.2             pyha770c72_0    conda-forge

$ conda list | grep fiona    
fiona                     1.8.21           py39hf45f784_2    conda-forge
```

それよりも新しいv1.8.22がGDAL v3.4.3を利用していると書かれているので、やはりGDALは3.5系以上でなければLOD2を含むCityGMLは取り扱えないようです。
https://github.com/Toblerity/Fiona/releases/tag/1.8.22

なので、Pythonで扱うには今の所すんなりとはいかず、v3.5以上のGDLAをインストールし、ogr2ogrのコマンドを利用して変換してから読み込む必要がありそうです。

ただし、fionaの開発版（v1.9b1）ではGDAL3.5.3が含まれていると書いているので、もう少し辛抱すればPythonでもサクッとCityGMLを読み込んで地理的な解析をすることができそうです。
https://github.com/Toblerity/Fiona/releases/tag/1.9b1

## 終わりに

ということで、Mac版のQGISやPythonでももう少しだけ待てば読み込めるようになるだろうが、現時点では未対応、ということがわかりました。

今回利用したツールは全てOSSですので、皆さんがコミットすることによってリリース時期が早まるかもしれませんし、建物以外の地物型に対応することだってできるかもしれません。

OSSにコミットしてどんどんCityGMLや周辺ツールを良くしていきましょう！