マップタイル作成マニュアル

第1.0版

はじめに
1. ラスタデータを扱う際の事前準備
2. ラスタ（RGB値）データのマップタイル化
3. ラスタ（ラスタ値）データのマップタイル化
4. ベクタデータのマップタイル化
1. (1) 事前準備
2. (2) マップタイルの生成
5. ベクトルタイルのスタイルjsonの作成
1. (1) ベクトルタイルの内容を確認する
2. (2) jsonデータを作成する
脚注

はじめに

本マニュアルは、航空レーザ解析の各業務で整備されたデータを県域などの広域で統合し、マップタイルとして配信するためのデータ変換作業の手順をまとめたものです。出典表示『室木直樹「マップタイル作成マニュアル（第1.0版）」』のみで自由にご利用いただけます（CC BY）。

(1) PC環境

使用したPCの諸元は、OS：Windows11、CPU：intel core i7 12th、RAM：64GB（DDR4-3200、32GB×２）、ROM：M.2接続のSSD 1TB×２である。

使用したソフトウェアは、QGIS（ver.3.34）及びそれと併せてインストールされるOSGeo4W Shell (OSGeo for Windows command shell)である。なお、OSGeo4Wは、GDALやPythonを一括で扱えるコマンドプロンプトである。

(2) 使用データ

基本的には、森林GISフォーラムが運用する「森林資源データ解析・管理標準仕様書」に沿って整備されたデータを、林野庁「森林情報に関するオープンデータ 標準仕様書（案）【航空レーザ森林資源解析データ編】」に沿ってマップタイルを作成する手順を解説している。

ただし、両標準仕様書において規定がないが、マップタイルを作成する上で追加的に規定されていることが望ましいと考えられた事項については、作業手順に追加している。追加した作業手順の概要は次の通りである。

航空レーザ解析業務において、業務範囲外のラスタ値にNodataを割り当てることが多い¹が、そのNodataをRGB値の場合は、白（255,255,255）に、ラスタ値の場合は-9999に統一する手順を追加した。Nodataの値を統一する理由は、マップタイルに変換する際に、透過処理や無効化処理を的確に行うためである。例えば、Nodataを意味する白を透過設定する際に、一部範囲でNodataが黒になっていると、その黒は透過されず残ってしまうことを防止するためである。
航空レーザ解析業務において、業務範囲外のラスタ値について透過設定をし、それをアルファバンド（RGBに次ぐ４バンド目）に指定している場合がある。GDALで変換作業を行う上で、３バンドまでしかないラスタデータと４バンドがあるラスタデータが混在するとエラーが発生するため、アルファバンドを指定する作業を追加した。

(3) 利用上の注意

筆者自身が非エンジニアであり、必ずしも、本稿で紹介する手順がもっとも効率的かつ合理的なものとは限らないことをご了承いただきたい。

1. ラスタデータを扱う際の事前準備

(1) 貸与データの欠損・異常を確認する

航空レーザ解析業務においては、ラスタデータにあっては国土基本図図郭（2,500やその４分の１）、ベクタデータにあっては、市町村等の行政区域を単位として作成されることが多い。この結果、特にラスタデータにあっては、大量のファイルを扱うことになるが、貸与されたデータからいきなりマップタイルを作成すると、貸与データに不足があること、Nodataが適切に設定されていない等の異常に事前に気づくことができず、出戻り作業が発生することも少なくない。

他方で、QGISでは、経験上、データファイル数が100を超えてくると動作が重くなることから、ある程度取り回しの良いサイズにラスタデータをマージし、マージされたラスタでエラーがないか確認するほうが効率がよい。このため、以下の手順により、貸与データを中間加工することを推奨する。

(2) tifをgeotiffに変換する

航空レーザ解析業務で作成されたラスタデータは、基本的には、tifとtfwの２つで構成されている。しかしながら、GDAL等では位置情報付きのgeotiffでの作業が前提となることから、gdalwarp²のループ処理により、一括変換を行う。

コマンド例

for %f in (C:\original\*.tif) do gdalwarp -t_srs EPSG:6678 -r bilinear “%f” “C:\Converted\%~nf.tif”

このコマンドは、OSGeo4W Shellで実行するバッチ処理である。C:\original フォルダ内の全.tifファイルを「バイリニア補間」で補正しながら、EPSG:6678の座標系に再投影して、C:\ Convertedフォルダに変換後のデータを保存するもの。

▼解説

for %f in (C:\original\*.tif)
C:\original フォルダ内のすべての.tifファイルを対象にループ処理を行うための前置き。*.tifとすることで、originalフォルダ内の全てのtifを変換の対象とすることができる。フォルダやファイル名に全角文字があっても支障³はなく、任意のディレクトリを使うことができる。

gdalwarp
GDALの投影変換などに使う機能。

-t_srs EPSG:6678
座標系をEPSG:6678に変換することを指示するもの。EPSGコードは任意の数値を与えることができる。

-r bilinear
投影変換の際のリサンプリング（補間）にバイリニア補間を使用することを指示するもの。（なめらかに変換）。bilinearは、周囲のピクセル値を補間に使い、なめらかにするものであり、連続値データに向き、処理速度と品質のバランスに優れる。このほかに、カテゴリカルデータに向く、最近傍補間（nearest）、連続値に向き、さらに品質のよいバイキュービック（cubic）などもある。

“%f”
変換前のファイルパス。今回はループ処理により一括変換を行うこととしており、前置きした%fがファイルパスとして“ ”で指定されている。

“C:\Converted\%~nf.tif”
変換後のファイルパス。%~nfは、変換前のファイル名である%fから拡張子.tifを除いたものを変換後のファイル名とすることを指示するもの。これにより、Convertedフォルダに同一名称のファイルを出力することができる。

▼参考

このようにコマンドプロンプトが流れていけば成功である。

(3) 国土基本図図郭×業務範囲のベクタデータを作成する

(1)に記載したとおり、国土基本図図郭（2,500やその４分の１）のままでは、その後の作業の取り回しが悪いことから、より大きい図郭単位にマージするほうが望ましい。また、航空レーザ解析業務は、１業務で１都道府県全域のデータが整備されていることは稀であり、市町村単位などにより、複数にデータファイルが分かれていることが多い。そのため、それぞれお互いに業務区域からはみ出たデータが残されていると、データの重複が発生し、どちらか一方のデータが隠れてしまうこともあり得る。

そこで、次の手順により、国土基本図図郭（2,500やその４分の１）単位のデータをマージし、マップタイル変換を行うための中間データを作成する。なお、以下では、国土基本図図郭50,000の４分の１単位を例としているが、これは、前掲の林野庁「森林情報に関するオープンデータ 標準仕様書（案）【航空レーザ森林資源解析データ編】」に沿うデータを併せて作成できるためである。この単位によれば、個々のデータファイルは概ね１～４GBの容量になり、QGISの動作やG空間情報センターへのデータ掲載においても支障がないため、推奨できる。

①国土基本図図郭50,000ファイルを作成
QGISにおいて、MIERUNE社が作成したJapanese Grid Meshプラグインを使うなどにより、国土基本図図郭50,000のファイルを作成する。

②国土基本図図郭50,000の４分の１ファイルを作成
QGISにおいて、ベクタ＞調査ツール＞グリッドを作成により、①で作成した図郭範囲に合致される形で水平20km、垂直15kmのポリゴンを作成し、①のコードを基に②のコード名称を作成する。

③業務範囲ポリゴンでクリップ
QGISにおいて、ベクタ＞空間演算ツール＞切り抜く（Clip）により、②で作成したポリゴンを業務範囲ポリゴンでクリップする。

④クリップポリゴンを図郭ごとの分割ポリゴンとして出力
QGISにおいて、ベクタ＞データ管理ツール＞属性でレイヤを分割により、③で作成したクリップポリゴンを図郭単位に切り分けたジオパッケージで保存する。（分割の基準に図郭のコードを指定する。）

▼参考

(4) 図郭ごとのgeotiffを結合し、分割ポリゴンでクリップする

(2)で変換した国土基本図図郭単位のgeotiffをマージし、(3)で作成した分割ポリゴンでクリップする方法は次のとおりである。

① 仮想ラスタの作成

(3)④の分割ポリゴンでクリップするのに十分な範囲のラスタを特定し、あらかじめマージしたgeotiffを作成しておくことは、必要な図郭の選定に手間がかかる上、最終的に使いもしない大容量のgeotiffを作り出すだけ（ストレージの容量をひっ迫するだけ）の無駄な行為である。そのため、在り様のgeotiff全てを使った仮想ラスタ（.vrt）を作成し、最終的にクリップに用いるGeotiffを仮想ラスタとして用意しておく。GDALのgdalbuildvrtを用いる。

コマンド例

gdalbuildvrt -vrtnodata “0 0 0” “C:\Clipped\merge.vrt” “C:\Converted\*.tif”

▼解説

gdalbuildvrt
複数のラスタを仮想的に一つにまとめた.vrtファイルを生成するGDALの機能。

-vrtnodata “0 0 0”
入力ラスタのうちRGB(0, 0, 0)＝黒をNoData（無効値）として扱うための指示。あらかじめQGISで画像を展開し、Nodataとして扱うべきRGB値を把握し、必要なものを入力する。なお、１の処理につき、１のRGB値しか指定できない。（同時に２色をNodataとすることはできない。）

“C:\Clipped\merge.vrt”
vrtファイルを保存するファイルを指定する。フォルダ名もファイル名（merge）も任意。

“C:\Converted\*.tif”
vrtファイルを生成する基となるgeotiffファイルの保存先を指定する。*.tifとすることで、Convertedフォルダ内の全てのtifがマージされる。

▼参考

このようにコマンドプロンプトが流れていけば成功である。

② 仮想ラスタの分割ポリゴンによる分割

上記で作成した仮想ラスタを(3)④で作成した分割ポリゴンでクリップしたgeotiffを作成する。GDALのgdalwarpを用いる。

コマンド例

gdalwarp -cutline “C:\code_10LD2.gpkg” -crop_to_cutline -co COMPRESS=DEFLATE -co ZLEVEL=6 -co BIGTIFF=YES -t_srs EPSG:6678 “C:\Clipped\merge.vrt” “C:\clipped\10LD2.tif”

▼解説

gdalwarp
GDALの機能である。１(2)のgeotiff変換に用いたものと同じである。

-cutline “C:\code_10LD2.gpkg”
仮想ラスタをgeotiffとして切り出す境界を表わすポリゴンを指定する。

-crop_to_cutline
ポリゴンの外側を除去し、ラスタをトリミングする。

-co COMPRESS=DEFLATE
出力ラスタを可逆圧縮し、データ容量を小さくする。なお、DEFLATEのほか、LZW圧縮（COMPRESS=LZW）も選択肢となる。

-co ZLEVEL=6
圧縮率を0〜9の数値で指定する。6は圧縮率と処理速度のバランスがとれた選択肢。数値を大きくするほど圧縮率は高くなるが、処理速度が遅くなる。作業ストレージに余裕があれば非圧縮（=NONEとする）でも構わない。

-co BIGTIFF=YES
出力ラスタのデータ容量が4GBを超える可能性がある場合には指定しておく必要がある。国土基本図図郭50,000の４分の１を基本単位とする場合は、この指定が必要（たまに４GBを超える場合がある）。

-t_srs EPSG:6678
出力ラスタの座標系をEPSG:6678に指定するもの。任意の座標系を選択できる。

“C:\Clipped\merge.vrt”
gdalbuildvrtで作成した仮想ラスタ（VRT）を指定する。

“C:\clipped\10LD2.tif”
クリップした出力ラスタのフォルダ名、ファイル名を指定する。

▼参考

このようにコマンドプロンプトが流れていけば成功である。

2. ラスタ（RGB値）データのマップタイル化

微地形表現図や林相識別図、簡易オルソ画像などのRGB値で構成されるラスタデータをマップタイル（XYZタイル）に変換する手順を記載する。

マップタイル化する前に処理すべきこととして、（１）Nodataの値を統一しておき、タイル生成時に当該Nodataを透過設定できるようにしておくか、（２）透過設定したいRGB値をアルファバンドに設定しておくか、の２パターンの処理がある。（１）と（２）は、必ずしも両方をする必要はなく、貸与データを踏まえながら任意に一方を選択すればよい。例えば、既にアルファバンドが指定されている貸与データがある、Nodataに割り当てられているRGB値が白・黒以外もあるなどパターンが多い、という場合には、全てのデータにアルファバンドを指定するほうが作業量は少ない。他方で、ほとんどの貸与データにおいて、Nodataに白色が充てられており、一部のデータでNodataを白色に置き換えるだけで済む場合は、（１）を選択するほうがよい。一方で、以下のような具体事例を踏まえると、あらかじめ、不用な部分は透明にしてある（２）のほうが後々のエラーを抑制でき、無難であると考えられる。

▼参考｜好ましくない公開データの例

(1) Nodataとする値を統一する手順

① Nodataに指定されている値の確認

GDALのgdalinfoを用い、Nodataが設定されているか、設定されている場合は、どのRGB値がNodataになっているかを確認する。

コマンド例

gdalinfo “C:\Converted\10LD1924.tif”

▼参考

② 特定のRGB値をRGB（255,255,255）＝白に変換し、それをNodataにするpythonコードを作成する

以下のコードをテキストエディタにコピーし、C:\clippedに、RGB値を変換したいgeotiffが保存されているフォルダのパスを、C:\clipped\convertedにRGBを変換した後のgeotiffを保存したいフォルダのパスを入力し、mask = (r == 0) & (g == 0) & (b == 0)に変換したいRGB値を入力し、pythonコードとして保存する⁴。

テキストエディタがない場合は、メモ帳アプリにおいて、ファイル名称replace_rgbとする場合は、replace_rgb.txtとして保存した上で、その拡張子.pyに書き換え、replace_rgb.pyとして保存することにより作成できる。

コード例

import os
from osgeo import gdal
import numpy as np

input_dir = r”C:\clipped“
output_dir = r”C:\clipped\converted“
os.makedirs(output_dir, exist_ok=True)

for filename in os.listdir(input_dir):
    if filename.lower().endswith((“.tif”, “.tiff”)):
        input_path = os.path.join(input_dir, filename)
        output_path = os.path.join(output_dir, filename)
        print(f”処理中: {filename}”)

      ds = gdal.Open(input_path)
        cols = ds.RasterXSize
        rows = ds.RasterYSize
      r = ds.GetRasterBand(1).ReadAsArray()
        g = ds.GetRasterBand(2).ReadAsArray()
        b = ds.GetRasterBand(3).ReadAsArray()

      mask = (r == 0) & (g == 0) & (b == 0)
        r[mask] = 255
        g[mask] = 255
        b[mask] = 255

      driver = gdal.GetDriverByName(“GTiff”)
        out_ds = driver.Create(output_path, cols, rows, 3, gdal.GDT_Byte,                                options=[ “COMPRESS=DEFLATE”, “ZLEVEL=6” , “BIGTIFF=YES”])
        out_ds.SetGeoTransform(ds.GetGeoTransform())
        out_ds.SetProjection(ds.GetProjection())
      out_ds.GetRasterBand(1).WriteArray(r)
        out_ds.GetRasterBand(2).WriteArray(g)
        out_ds.GetRasterBand(3).WriteArray(b)

      for i in range(1, 4):
            out_ds.GetRasterBand(i).SetNoDataValue(255)
      out_ds.FlushCache()
        del out_ds
        del ds

print(“処理完了”)

③ OSGeo4Wで作成したpythonコードを実行する

コマンド例

python “C:\clipped\replace_rgb.py”

▼解説

python
OSGeo4Wにおいて、pythonを使う場合は、最初にpythonと記述する。

“C:\clipped\replace_rgb.py”
②で作成したpythonコードの保存先をフルパスで入力する。ファイルパスは、当該ファイルを右クリックすることで調べることができる（Windows10の場合は、SHIFTを押しながら当該ファイルを右クリックすることで調べることができる。）

▼参考

(2) 透過率をアルファバンドに設定する手順

gdalwarpのルート処理により、フォルダ内のtifを一括して、NodataであるRGB(255,255,255)をアルファバンドに指定し、透過率を記録させる手順は次のとおり。

コマンド例

for %f in (“C:\clipped\*.tif”) do gdalwarp -dstalpha -srcnodata “255 255 255” -co COMPRESS=DEFLATE -co ZLEVEL=6 -co BIGTIFF=YES “%f” “C:\clipped\alpha\%~nxf”

▼解説

for %f in (“C:\clipped\*.tif”)
C:\clippedフォルダ内のすべての.tifファイルに対してループ処理を行う前置き。

gdalwarp
GDALの機能である。１(2)のgeotiff変換に用いたものと同じである。

-dstalpha
アルファバンドに透明度を設定する指示。

-srcnodata “255 255 255”
RGB(255,255,255)＝白を透明として設定する指示。

“%f”
変換前のファイルパス。今回はループ処理により一括変換を行うこととしており、前置きした%fがファイルパスとして“ ”で指定されている。

“C:\clipped\alpha\%~nxf”
変換後のファイルパス。%~nfは、変換前のファイル名である%fから拡張子.tifを除いたものを変換後のファイル名とすることを指示するもの。これにより、clippedフォルダ内のalphaフォルダに同一名称のファイルを出力することができる。

▼参考

(3) マップタイルの作成

① Pythonコード「gdal2tiles.py」の保存場所を確認

マップタイルの作成では、pythonコードであるgdal2tiles.pyを使用する。これは、QGISでマップタイルを作成する際にバックグラウンドで使われているコードであるが、その保存場所がQGISのバージョンによって異なる場合があるため、OSGeo4Wで場所を確認する。（のちの作業でフルパスを指定する必要があるため、調べておく必要。）

コマンド例

where gdal2tiles.py

▼解説

where
OSGeo4Wにおいて、最初にwhereと記述し、検索したいキーワードを入れると保存場所を教えてくれる。

▼参考

② Nodataの値が統一され、アルファバンドが設定されていない場合

２(1)の作業を経て、Nodataの値を統一し、アルファバンドが設定されていない場合のデータからマップタイルを作成する手順は、２(1)で作成したgeotiff一式を仮想ラスタで便宜上つなぎ合わせた上で、その仮想ラスタの範囲に対して、マップタイルの作成を指示する流れとなる。

コマンド例

gdalbuildvrt -vrtnodata “255 255 255” “C:\maptiles\merge.vrt” “C:\clipped\*.tif”

▼解説

gdalbuildvrt
複数のラスタを仮想的に一つにまとめた.vrtファイルを生成するGDALの機能。

-vrtnodata “255 255 255”
入力ラスタのうちRGB(255,255,255)＝白をNoData（無効値）として扱うための指示。ここで指定したRGB値が透過設定される。

“C:\maptiles\merge.vrt”
vrtファイルを保存するファイルを指定する。フォルダ名もファイル名（merge）も任意。

“C:\clipped\*.tif”
vrtファイルを生成する基となるgeotiffファイルの保存先を指定する。*.tifとすることで、clippedフォルダ内の全てのtifがマージされる。

次に、以下のコマンドを入力する。

コマンド例

python “C:\Program Files\QGIS 3.22.16\apps\Python39\Scripts\gdal2tiles.py” –tiledriver=WEBP –webp-quality=75 –xyz -r antialias -z 8-18 –s_srs EPSG:6678 –exclude –resume –processes=12 “C:\maptiles\merge.vrt” “C:\maptiles”

▼解説

python “C:\…gdal2tiles.py”
gdal2tiles.py を Python 経由で実行する指示。①でwhereで検索した結果のとおり、フルパスでpythonコードの保存先を指定すること。

–tiledriver=WEBP
タイルの拡張子にwebPを指定するもの。webpは、広く一般に使われているpngと同様に透過設定もできる上、軽量かつ高圧縮が可能なため推奨する。

–webp-quality=75
webpの圧縮品質を指定するもの。0〜100を選択でき、一般に75〜85程度がバランスがよい。なお、pngを選択した場合は、これは記述する必要がない。

–xyz
XYZタイル形式を指定するもの。これを指定しない場合は、WMTS（Web Map Tile Service）形式になってしまう。

-r antialias
リサンプリング方法にアンチエイリアスを指定するもの。微地形表現図や林相識別図などのいわゆる画像系は、アンチエイリアスを選択しておくと、画像の滑らかさを確保できる。

-z 8-18
ズームレベルを指定するもの。ZL-18まであれば、基本的には十分な解像度を得られる。日本全国を対象とするタイルを生成する場合は、ZL-5から作成することも検討。

–s_srs EPSG:6678
タイルの生成に用いるgeotiffの座標系を指定するもの。指定がない場合は、EPSG:3857（WEBメルカトル）とされてしまう。

–exclude
保存先に指定されているフォルダに、既に存在する出力タイルの生成をスキップさせるための指示。

–resume
作業が中断してしまった場合に、前回の中断したところから処理を再開することを指示するもの。

–processes=12
並列処理に使用するCPUプロセス数を指示するもの（マルチコアで高速処理を図ることができる）。数字は個々のPC環境に応じて選択すること。

“C:\maptiles\merge.vrt”
マップタイルの基となる仮想ラスタ（gdalbuildvrtで生成したもの）参照先を指定するもの。

“C:\maptiles”
生成したマップタイルの保存先を指定するもの。

▼参考

③ Nodataの値がアルファバンドで透過設定されている場合

２(2)の作業を経て、アルファバンドが設定されたgeotiffに基づき、マップタイルを作成する場合は、既に透過設定が行われており、仮想ラスタを生成する過程で、特定のRGB値について透過させなさいという指示を与える必要がないため、コマンドの一部を省略することができる。これ以外は、上記と同様である。

アルファバンドの設定がない場合	アルファバンドの設定がある場合
gdalbuildvrt -vrtnodata “255 255 255” “C:\maptiles\merge.vrt” “C:\clipped\*.tif”	gdalbuildvrt “C:\maptiles\merge.vrt” “C:\clipped\*.tif”

3. ラスタ（ラスタ値）データのマップタイル化

DEMやDCHMなどのラスタ値が格納されたラスタデータをデータPNG型のマップタイルに変換する手順を記載する。

マップタイル化する前に処理すべきこととして、第１章で記載した事前準備のほか、Nodataの値を統一しておき、タイル生成時に当該Nodataを無効化処理できるようにしておくことが重要である。第２章においてもNodataの統一について言及しているが、この第３章で紹介する手順（ラスタ値が格納されたラスタデータ）については、Nodataを-9999に統一することが異なる。データPNGを生成する過程において必ずしもNodataを-9999とすることが必須ではないものの、以下に記載する計算式（標高値をRGB値に変換する過程）においてエラーが生じず、通常の標高値との混同が生じないものの代表例として-9999を採用することを推奨するものである。

標高値のRGB値への変換にあっては、いくつかの計算式があるが、本稿においては、Terrain-RGBとPNG標高タイルの２種類の生成について手順を紹介する。

Terrain-RGB
Terrain-RGBは、Mapbox社で考案された以下の変換式が採用されている。QGISやMaplibre GL JSなどで表示でき、FOSS4Gとの相性がよい。記述できる標高値は、RGB（0,0,0）の-10,000mからRGB（255,255,255）の1,667,722mである。

Terrain-RGBの変換式

標高値 = -10000+(R×256×256+G×256+B)×0.1

PNG標高タイル
PNG標高タイルは、国土地理院の標高タイルで採用されている変換式を用いるものである。日本国内のデフォルトとして採用されている変換式のため、他機関データと組み合わせて使用する場合に利点がある。記述できる標高値は、RGB（128,0,0）の-83,886mからRGB（127,255,255）の83,886mである。

PNG標高タイルの変換式

X = 2¹⁶×R+2⁸×G+Bとして
X < 2²³ の場合：標高値 = X×0.01
X > 2²³ の場合：標高値 =(X–2²⁴)×0.01

(1) Nodataとする値を統一する手順

① Nodataに指定されている値を確認する

GDALのgdalinfoを用い、Nodataが設定されているか、設定されている場合は、どのラスタ値がNodataになっているかを確認する。

コマンド例

gdalinfo “C:\geotiff\584146.tiff”

▼参考

② Pythonコード「gdal_calc.py」の保存場所を確認

ラスタ値の変換には、Pythonコードのgdal_calc.pyを用いる。保存場所がQGISのバージョンによって異なる場合があるため、OSGeo4Wで場所を確認する。

コマンド例

where gdal_calc.py

▼解説

where
OSGeo4Wにおいて、最初にwhereと記述し、検索したいキーワードを入れると保存場所を教えてくれる。

▼参考

③ Nodataが指定されていないgeotiffに新たにNodataを指定する

gdalinfoで確認した結果、Nodataが指定されていないgeotiffが処理データに含まれていることが確認された場合は、QGISで展開の上、Nodataとするべき値を特定したうえで、gdal_translateのループ処理により、一括して特定の値をNodataに位置付ける。これは、あくまで既存のラスタ値の一つをNodataに位置付けるものであり、例えば、航空レーザ測量業務の対象範囲外について、あらかじめラスタ値0が与えられていた一方、その数値がNodataに位置付けられていない（ただの0mとして扱われている）場合などに用いるものである。

一方、既に別の値（例えば0）がNodataに位置付けられていて、その数値を-9999に置き換え、-9999をNodataとして用いていくこととする場合は、以下の④で解説するので用いる方法を間違えないよう留意すること。

コマンド例

for %f in (C:\ original\*.tif) do gdal_translate -a_nodata 0 “%f” “C:\Converted\%~nf.tif”

▼解説

for %f in (C:\original\*.tif)
C:\original フォルダ内のすべての.tifファイルを対象にループ処理を行うための前置き。*.tifとすることで、originalフォルダ内の全てのtifを変換の対象とすることができる。

gdal_translate
投影変換を伴わないラスタの変換などに使うGDALの機能。

-a_nodata 0
ラスタ値＝0をNodataに指定する指示。

“%f”
変換前のファイルパス。今回はループ処理により一括変換を行うこととしており、前置きした%fがファイルパスとして“ ”で指定されている。

“C:\Converted\%~nf.tif”
変換後のファイルパス。%~nfは、変換前のファイル名である%fから拡張子.tifを除いたものを変換後のファイル名とすることを指示するもの。これにより、Convertedフォルダに同一名称のファイルを出力することができる。

▼参考

④ Nodataに指定されている値を別の値に置き換えてNodataとして扱う

gdalinfoで確認した結果、別の値（例えば0）が既にNodataに位置付けられている場合において、そのNodataとされているラスタ値（例えば0）を-9999に置き換え、-9999をNodataとして用いていくこととする場合は、Pythonコードのgdal_calc.pyを用いる。③の手順で使ったgdal_translateは、あくまで-9999を新たにNodataに設定することはできたとしても、すでにNodataとされている0は０のままであり、置き換えができるわけではないことに留意する必要がある。

以下のコードをテキストエディタにコピーし、“C:\original”にラスタ値を変換したいgeotiffが保存されているフォルダのパスを、“C:\Converted”にラスタ値の一部を-9999に置き換えたgeotiffを保存したいフォルダのパスを入力し、“C:\Program Files\QGIS 3.34.13\bin\python.exe”と“C:\Program Files\QGIS 3.34.13\apps\Python312\Scripts\gdal_calc.py”にwhereで調べたpython.exe とgdal_calc.pyの保存先を、0に変換したいラスタ値を入力し、pythonコードとして保存する。

テキストエディタがない場合は、メモ帳アプリにおいて、replace_to_nodata.txtとして保存した上で、その拡張子を.pyに置き換え、replace_to_nodata.pyとして保存することにより作成できる。なお、gdal_calc_path = rの後ろは、シングルクォート‘　’で括られていることに注意すること。

コード例

import os
import glob
import subprocess

input_dir = r“C:\original”
output_dir = r“C:\Converted”
gdal_calc_path = r’“C:\Program Files\QGIS 3.34.13\bin\python.exe” “C:\Program Files\QGIS 3.34.13\apps\Python312\Scripts\gdal_calc.py”‘

os.makedirs(output_dir, exist_ok=True)
tiff_files = glob.glob(os.path.join(input_dir, “*.tif”))
for tiff_path in tiff_files:
    filename = os.path.basename(tiff_path)
    output_path = os.path.join(output_dir, filename)

  cmd = f'{gdal_calc_path} -A “{tiff_path}” –outfile=”{output_path}” –calc=”A*(A!=0)+(-9999)*(A==0)” –NoDataValue=-9999 –overwrite’

    print(“Running:”, cmd)
    subprocess.run(cmd, shell=True, check=True)

print(“処理完了”)

⑤ OSGeo4Wで作成したpythonコードを実行する

コマンド例

python “C:\replace_to_nodata.py”

▼解説

python
OSGeo4Wにおいて、pythonを使う場合は、最初にpythonと記述する。

“C:\replace_to_nodata.py”
作成したpythonコードの保存先をフルパスで入力する。ファイルパスは、当該ファイルを右クリックすることで調べることができる（Windows10の場合は、SHIFTを押しながら当該ファイルを右クリックすることで調べることができる。）

▼参考

(2) Terrain-RGBを生成する

① Nodataを無効化処理する（Nodataとして扱っていた-9999をラスタ値-9999として扱う）

Terrain-RGBは、ラスタ値をRGB値に変換したgeotiffを生成し、当該geotiffにおいて、Nodataとすべき値（RGB値）を透過設定しつつ、マップタイルを生成する、という二段階を経る。ただし、ラスタ値をRGB値に変換する上では、Nodataの設定があるgeotiffを扱うことができず、一旦、Nodataとされている値を再びラスタ値として認識できるようにする（Nodataを無効化処理する）必要がある⁵。

Nodataの無効化処理には、GDALのgdal_translateを用い、ループ処理によりフォルダ内のtifを一括して変換する。

コマンド例

for %f in (“C:\converted\*.tif”) do gdal_translate -a_nodata none “%~f” “C:\Disabled\%~nxf”

▼解説

for %f in (“C:\converted\*.tif”)
C:\convertedフォルダ内のすべての.tifファイルを対象にループ処理を行うための前置き。*.tifとすることで、convertedフォルダ内の全てのtifを変換の対象とすることができる。

gdal_translate
投影変換を伴わないラスタの変換などに使うGDALの機能。

-a_nodata none
Nodataを無効化処理する指示。

“%~f”
変換前のファイルパス（フルパス）。今回はループ処理により一括変換を行うこととしており、前置きした%fがファイルパスとして“ ”で指定されている。

“C:\Disabled\%~nxf”
変換後のファイルパス。%~nxfは、変換前のファイル名を変換後のファイル名とすることを指示するもの。これにより、Convertedフォルダに同一名称のファイルを出力することができる。

▼参考

② Pythonパッケージツールの「rio-rgbify」をインストールする

ラスタ値なgeotiffをTerrain-RGBの変換式に基づくRGB値なgeotiffへの変換は、pythonパッケージツールのrio-rgbifyを用いる。rio-rgbifyは、QGISをインストールしただけの初期環境ではインストールされていないため、OSGeo4Wにおいて以下のコマンドにより、インストールを実行する。

コマンド例

pip install rio-rgbify

次に、pythonパッケージツールが保存されている場所を確認するため、OSGeo4Wにおいて以下のコマンドにより、検索を行う。

コマンド例

where pip

▼参考

筆者のPC環境では、「C:\Program Files\QGIS 3.34.13\apps\Python312\Scripts」にあることが確認できる。この検索結果は、以下③において使用するのでメモしておくこと。

③ ラスタ値をRGB値に変換したgeotiffを生成する

OSGeo4Wのコマンドでpythonパッケージツールが起動するよう、以下のコマンドにより環境設定を行う⁶。C:\Users…\Scriptsは、上記で調べた結果を指定する（ただし、\pip.batは含めず、\Script止めとすること）。

コマンド例

set PATH=%PATH%;C:\Users\naoki\AppData\Roaming\Python\Python312\Scripts

▼参考

以下のとおり何も表示されず、次のコマンド入力が行えるようになっていれば、正常に進んでいる。

次に、rio-rgbifyのループ処理により、フォルダ内のtifを一括してTerrain-RGBの変換式によるRGB値化されたgeotiffを生成する。

コマンド例

for %f in (C:\Disabled\*.tif) do rio rgbify -b -10000 -i 0.1 –co COMPRESS=DEFLATE –co ZLEVEL=6 “%f” “C:\rgb\%~nxf”

▼解説

for %f in (C:\Disabled\*.tif)
C:\Disabledフォルダ内のすべての.tifファイルを対象にループ処理を行うための前置き。*.tifとすることで、Disabledフォルダ内の全てのtifを変換の対象とすることができる。

rio rgbify
次の記述の計算式に基づき、ラスタ値をRGB値に変換するpythonパッケージツール。

-b -10000
RGB（0,0,0）の時に-10,000mと規定するベース値を表す。Terrain-RGBの変換式に基づき、-10,000としている。

-i 0.1
RGB値が１変わるたびにラスタ値をどのように動かすかを規定するインターバル値。Terrain-RGBの変換式に基づき、0.1としている。

–co COMPRESS=DEFLATE
RGB値化したgeotiffをDeflate圧縮するように指示するもの。（任意）

–co ZLEVEL=6
Deflate圧縮のレベル（1-9）のうち、レベル6（比較的高圧縮かつ処理が早め）を設定するもの。（任意）

“%f”
変換前のファイルパス（フルパス）。今回はループ処理により一括変換を行うこととしており、前置きした%fがファイルパスとして“ ”で指定されている。

“C:\rgb\%~nxf”
変換後のファイルパス。%~nxfは、変換前のファイル名を変換後のファイル名とすることを指示するもの。これにより、rgbフォルダに同一名称のファイルを出力することができる。

▼参考

④ マップタイルの作成

微地形表現図などのマップタイル化に用いたものと同様、pythonコードであるgdal2tiles.pyを使用する。ただし、以下のようにコマンド入力で数点異なるところがあるため、留意すること。

コマンド例

gdalbuildvrt -vrtnodata “0 0 10” “C:\maptiles\merge.vrt” “C:\rgb\*.tif”

▼解説

gdalbuildvrt
複数のラスタを仮想的に一つにまとめた.vrtファイルを生成するGDALの機能。

-vrtnodata “0 0 10”
入力ラスタのうちRGB(0,0,10)＝NoData（無効値）として扱うための指示。ここで指定したRGB値が透過設定される。RGB（0,0,10）は、Terrain-RGBの変換式によればラスタ値＝-9999であり、もともとNodataであるべきラスタ値を透過させることに相当する。

“C:\maptiles\merge.vrt”
vrtファイルを保存するファイルを指定する。フォルダ名もファイル名（merge）も任意。

“C:\rgb\*.tif”
vrtファイルを生成する基となるgeotiffファイルの保存先を指定する。*.tifとすることで、rgbフォルダ内の全てのtifがマージされる。

コマンド例

python “C:\Program Files\QGIS 3.34.13\apps\Python312\Scripts\gdal2tiles.py” –tiledriver=PNG -r near –xyz -z 8-18 –s_srs EPSG:6668 –exclude –resume –processes=12 “C:\maptiles\merge.vrt” “C:\maptiles”

▼解説

python “C:\…gdal2tiles.py”
gdal2tiles.py を Python 経由で実行する指示。フルパスでpythonコードの保存先を指定すること。

–tiledriver=PNG
タイルの拡張子にpngを指定するもの。

-r near
リサンプリング方法に最近傍を指定するもの。アンチエイリアスなどの計算過程が含まれるもを選択することは、RGB値上で計算が行われ、現実にそぐわない標高値が返されるおそれがあり、お勧めしない。

–xyz
XYZタイル形式を指定するもの。これを指定しない場合は、WMTS（Web Map Tile Service）形式になってしまう。

-z 8-18
ズームレベルを指定するもの。ZL-18まであれば、基本的には十分な解像度を得られる。日本全国を対象とするタイルを生成する場合は、ZL-5から作成することも検討。

–s_srs EPSG:6668
タイルの生成に用いるgeotiffの座標系を指定するもの。指定がない場合は、EPSG:3857（WEBメルカトル）とされてしまう。

–exclude
保存先に指定されているフォルダに、既に存在する出力タイルの生成をスキップさせるための指示。

–resume
作業が中断してしまった場合に、前回の中断したところから処理を再開することを指示するもの。

–processes=12
並列処理に使用するCPUプロセス数を指示するもの（マルチコアで高速処理を図ることができる）。数字は個々のPC環境に応じて選択すること。

“C:\maptiles\merge.vrt”
マップタイルの基となる仮想ラスタ（gdalbuildvrtで生成したもの）参照先を指定するもの。

“C:\maptiles”
生成したマップタイルの保存先を指定するもの。

▼参考

(3) PNG標高タイルを生成する

PNG標高タイルは、全国Q地図の運営者さんが公開しているpythonコードgdal2NPtiles.pyを使い生成する。（１）の手順に基づき、Nodataを-9999に統一したラスタ値なgeotiffが用意できていれば、そのまま、マップタイルの生成に移行できる。

① gdal2NPtiles.pyのダウンロード

GitHubからZIPをダウンロードし、その中にあるgdal2NPtiles.pyを任意の場所に保存する。

▼参考

② gdal2NPtiles.pyの実行

コマンド例

gdalbuildvrt -vrtnodata -9999 “C:\maptiles\merge.vrt” “C:\converted\*.tif”

▼解説

gdalbuildvrt
複数のラスタを仮想的に一つにまとめた.vrtファイルを生成するGDALの機能。

-vrtnodata -9999
入力ラスタのうちラスタ値-9999をNoData（無効値）として扱うための指示。ここで指定したRGB値が透過設定される。

“C:\maptiles\merge.vrt”
vrtファイルを保存するファイルを指定する。フォルダ名もファイル名（merge）も任意。

“C:\ converted \*.tif”
vrtファイルを生成する基となるgeotiffファイルの保存先を指定する。*.tifとすることで、rgbフォルダ内の全てのtifがマージされる。

コマンド例

python “C:\gdal2NPtiles.py” –numerical -z8-18 –xyz –srcnodata=-9999 –processes=12 “C:\maptiles\merge.vrt” “C:\maptiles”

▼解説

python “C:\gdal2NPtiles.py”
gdal2NPtiles.py を Python 経由で実行する指示。フルパスでpythonコードの保存先を指定すること。

–numerical
数値PNGタイル生成モードを有効にする指示。

-z8-18
ズームレベルを指定するもの。ZL-18まであれば、基本的には十分な解像度を得られる。日本全国を対象とするタイルを生成する場合は、ZL-5から作成することも検討。

–xyz
XYZタイル形式を指定するもの。これを指定しない場合は、WMTS（Web Map Tile Service）形式になってしまう。

–srcnodata=-9999
タイルの生成に用いるgeotiffのNodataを-9999と認識させる指示。

–processes=12
並列処理に使用するCPUプロセス数を指示するもの（マルチコアで高速処理を図ることができる）。数字は個々のPC環境に応じて選択すること。

“C:\maptiles\merge.vrt”
マップタイルの基となる仮想ラスタ（gdalbuildvrtで生成したもの）参照先を指定するもの。

“C:\maptiles”
生成したマップタイルの保存先を指定するもの。

▼参考

4. ベクタデータのマップタイル化

森林計画図や樹種ポリゴン等のベクタデータから、ベクトルタイルの作成手順を紹介する。

(1) 事前準備

ラスタタイルの場合は、第１章として事前準備を設けたが、ベクトルタイルの場合は、以下の点をGIS等において、処理しておくだけである。具体の作業手順については、他誌が詳しいので省略する。

１つのベクトルタイルとする範囲について、ベクタデータをクリップ、マージ⁷しておく。
属性フィールドの記載がそのまま使われることから、略称や記号表記について、視認性のよい表現に変換しておく。

(2) マップタイルの生成

QGISのプロセッシングツールのうち、「ベクタタイルを書き出し（XYZ）」を用いる。パラメータの入力は次のとおり入力する。

ファイルのテンプレート
初期入力値　{z}/{x}/{y}.pbf　のままでよい。

入力レイヤ
１つのベクトルタイルとしたいレイヤを全て選択する。（基本的には、上記のとおり、あらかじめ１つのレイヤにマージしておき、それのみを選択することを推奨。）

最小ズーム値
ラスタタイルと異なり、ベクタデータの性質により１タイルあたりのデータ容量が変動する。特に、ズームレベルが小さい場合、１タイルあたりのデータ容量が大きくなり、WEB利用に向かいない場合もあるので、注意が必要である。

最大ズーム値
ラスタタイルと同様、ズームレベルを大きくすればするほど、タイルのサイズが小さくなる。元のベクタデータのポリゴン等のサイズがそれほど小さくない場合は、ズームレベルに大きな値を設定しなくても十分な場合もある。そのような場合に、ズームレベルを大きくすると、１のポリゴンが不必要に四角く切り分けられて、かえって見づらい場合もある。

領域、出力フォルダ
適宜指定する。

▼参考

▼参考｜林班ポリゴン（１ポリゴンあたり0.3km²/属性5列）の例

5. ベクトルタイルのスタイルjsonの作成

ベクトルタイルは、ラスタタイルと異なり、属性情報に基づきスタイル設定を任意で設定できることがメリットである一方、ベクトルタイルの利用者がQGIS等において、一からスタイル設定を行うのはやや手間である。そこで、基本的なスタイル設定はベクトルタイルの配信と合わせ、style.jsonも配信しておくことが望まれる。style.jsonは、テキストエディタで記述する。

(1) ベクトルタイルの内容を確認する

QGISでベクトルタイルを接続し、地物情報を確認する。選択したポリゴンの名称を確認し、今回の事例の場合は、style.jsonの記述に必ず必要なレイヤ名称（rinpan_tohoku ― 林班）と、ラベル表示の設定に使いたい属性名称（林班主番）をメモしておく。

▼参考

(2) jsonデータを作成する

style.jsonの基本骨格は次のとおりである。

以下の例は、
・レイヤ名（rinpan_tohoku — 林班）とされていたベクトルデータから生成し、
・file:///C:/maptiles/{z}/{x}/{y}.pbfとして保存されているベクトルタイルについて、
・ズームレベル8-18は、透過率50％の深緑でポリゴンを塗りつぶし、
・ズームレベル8-18は、透過率0％、1.0ポイント幅の深緑の実線でポリゴンを括り、
・ズームレベル14-18は、白で縁取りされた16ポイントの黒文字で林班主番をラベル表示する
という記述がされている。このほかにも、属性値によって塗りつぶしの色を変えるなどのアレンジができるが、林野庁が公開している樹種ポリゴンや崩壊箇所等判読結果のstyle.jsonなども参考に学習してほしい。

コード例

{
“version”: 8,
“name”: “林班ポリゴン初期設定”,
“sources”: {
“rinpan”: {
“type”: “vector”,
“tiles”: [
“file:///C:/maptiles/{z}/{x}/{y}.pbf”
],
“minzoom”: 8,
“maxzoom”: 18
}
},
“layers”: [
{
“id”: “rinpan-fill”,
“type”: “fill”,
“source”: “rinpan”,
“source-layer”: “rinpan_tohoku — 林班”,
“minzoom”: 8,
“maxzoom”: 18,
“paint”: {
“fill-color”: “rgba(51, 160, 44, 0.5)”,
“fill-outline-color”: “rgba(51, 160, 44, 1.0)“
}
},
{
“id”: “rinpan-label”,
“type”: “symbol”,
“source”: “rinpan”,
“source-layer”: “rinpan_tohoku — 林班“,
“minzoom”: 14,
“maxzoom”: 18,
“layout”: {
“text-field”: [“get”, “林班主番”],
“text-font”: [“Noto Sans Regular”, “Arial Unicode MS”],
“text-size”: 16,
“text-anchor”: “center”,
“symbol-placement”: “point”,
“text-allow-overlap”: false
},
“paint”: {
“text-color”: “#000000”,
“text-halo-color”: “#ffffff”,
“text-halo-width”: 1
}
}
]
}

▼参考｜上記style.jsonの表示例

脚注

RGB値は、白とすることが多いようであるが、黒（0,0,0）が設定されている場合もある。ラスタ値は、±3.4e¹⁸や、-32,768が設定されている場合もあり、航空レーザ解析業務では統一されていない。 ↩︎
単なるtifのgeotiff変換は、gdal_translateで実行するものであるが、経験上、貸与データには複数の座標系が混在していることが多々あり、（事実上、投影変換が空振りになる可能性はあるものの、）geotiffへの変換とともに、投影変換も可能なgdalwarpを使用することにした。 ↩︎
括弧が禁則である。基本的には全角文字を混ぜないように習慣化してほしい。 ↩︎
Pythonを使わずに、GDALだけで処理することもできるが、手順が増えるため、一括処理できるPythonコードを作成することとした。 ↩︎
そうであるならば、前処理において-9999をNodataとしなければよいとの考え方もあるが、Nodataに何が設定されているかの確認の習慣をつけることと、後段で紹介するPNG標高タイルの生成において重要となるため、このような手順の紹介としている。 ↩︎
PCのシステムプロパティから、この環境変数を経常指定することもできるようであるが、筆者の環境ではそれがうまくいかず、rio-rgbifyを使用するたびに、このコマンド処理を行っている。 ↩︎
必ずしも１つのベクタデータとなっている必要はないが、第５章で紹介するstyle.jsonの記述が複雑となることから、マージしておくことを推奨する。 ↩︎

ズームレベル小（8）の例	ズームレベル大（18）の例

（QGISでの表示）

１タイルあたり1-2MBくらいであれば、通信上の支障も少なく、タイルでポリゴンより十分に大きい。	タイルがポリゴンに比べ小さすぎ、データ容量も小さすぎる。この例では、ZL-18は不要。

変換前	変換後

Nodataの値を統一せず、白を透過設定してマップタイルを作成した（Nodataであるはずの黒が残った）例	Nodataの値は白に統一されていたが、白を透過設定せずにマップタイルを作成した例

変換前	変換後

変換前	変換後

gdalinfo の結果	QGISのレイヤプロパティ