最近空間解析を行う機会が多くなり、ジオプロセシング ツールをあれこれ使ってなんだかんだ処理をしています。
地理座標系の標高データ
傾斜角ツール（ArcGIS Pro Help より引用）
これは標高データから傾斜角を求めるツールですが、地理座標系が設定されている標高データをそのまま使用すると正常な処理結果になりません。なぜなら、地理座標系は XY座標の単位が「度」だからです。Z座標（高さ）はメートル（距離）なのに対して XY 座標が度（角度）だと正しい解析結果とはなりません。
処理を実行する前に XY座標と Z座標の軸を同じ距離単位に指定する必要があります。日本では通常高さはメートルで示しますので、XY座標の単位がメートルであり、且つ距離が正しい投影座標系を使用するべきです。メルカトル図法のような距離が正確ではない投影座標系も正しい結果とはなりません。全世界で任意の二地点間の距離を正しく表現できる投影法は存在しないので、局所であれば UTM座標系や平面直角座標系（日本国内の場合）を使用し、日本全国程度の範囲であればランベルト正角円錐図法やボンヌ図法などを使用します。
しかし、都度標高データを投影座標系に変換するのも手間です。全世界で様々な場所を処理する場合はなおさらです。
Z 係数の設定
傾斜角や陰影起伏などのツールには「Z 係数 (Z factor) 」というパラメーターがあり、XY 値に対する Z 値の倍率を設定することができます。
［傾斜角］ツールのパラメーター
水平方向も垂直方向も同じ距離で計算するならば 1 を設定すれば良いのですが、たとえば投影座標系の距離単位がメートルなのに対して高さ方向の単位がフィートの場合は、「1メートル = 0.3048フィート」としてZ 係数を 0.3048 に設定します。航空関係では距離はマイル、高さはフィートが使用されます。投影座標系の単位がマイルなら「1マイル= 5280フィート」として Z 係数を 5280 に設定します。
このように一律で係数をかければ計算できる場合は問題ないのですが、場所によって 1度の距離は異なるため 1度相当の距離を一律に求めることはできません。以前の記事を読んでいただくとよく分かりますが、経度は南北極点に近づくほど 0に近くなりますし、緯度は南北極点が最大長で赤道に近づくほど若干短くなります。
Z 係数換算表
先日空間解析のスペシャリストから Esri の方が書かれた ArcGIS Blog の記事を紹介してもらいました。
「1m＝XX度」を緯度別に求めたものです。
この記事によると、以下のような緯度の場合 Z 係数にそれぞれの値を入れて計算すれば正しい傾斜角や陰影起伏が求められるというものです。表は北半球に住んでいる日本人向けに 0度を下にしました。
| 緯度（南北） | Z 係数（メートル） | Z 係数（フィート） |
|---|---|---|
| 80° | 0.00005156 | 0.00001571 |
| 70° | 0.00002619 | 0.00000798 |
| 60° | 0.00001792 | 0.00000546 |
| 50° | 0.00001395 | 0.00000425 |
| 40° | 0.00001171 | 0.00000357 |
| 30° | 0.00001036 | 0.00000316 |
| 20° | 0.00000956 | 0.00000291 |
| 10° | 0.00000912 | 0.00000278 |
| 0° | 0.00000898 | 0.00000273 |
※参考記事は Longitude（経度）と書かれていますが、上の表で Latitude（緯度）における東西方向の 1m （もしくは 1フィート）が何度かを求めたのものです。同じ緯度なら経度間の距離は同じです。
つまり東京付近である北緯35度ならば 0.00001036 ～ 0.00001171 を設定すれば良いということになります。地理座標系の標高データを空間解析で処理する際は、必要に応じてこの Z 係数換算表を使用してください。厳密には緯度が 1度違うと Z 係数も変わってきます。参考記事には計算式が書かれていなかったので、中間の 35度の場合はどんな値になるか正確な値は分かりませんでした。これは後日実測してみます。
2018年10月21日追記
1度単位で Z係数を知りたい場合は、こちらの記事にある表の「経度1m角 (度) 」列をご利用ください。
新バージョンでの回避方法
ArcGIS Desktop 10.5 から、傾斜角ツールでは地理座標系の場合による Z 係数の計算が不要になりました。ArcGIS Pro はどのバージョンからか不明ですが、少なくともバージョン 2.2 では測地線に対応しています。各バージョンのヘルプをご覧ください。
- 傾斜角 (ArcMap 10.4)
- 傾斜角 (ArcMap 10.6 current)
- 傾斜角 (ArcGIS Pro 2.2 current)
※ArcGIS Pro には過去バージョンのヘルプは存在しません
ArcGIS Desktop 10.5 から、傾斜角ツールのパラメーターに［メソッド］が追加され、平面 (PLANAR) と測地線 (GEODESIC) が選択できるようになりました。メソッドで測地線を選択すれば、回転楕円体に基づいた傾斜角を求めるので、標高データの座標系が地理座標系や距離が正確でないメルカトル図法であっても正しい傾斜角を求められます。ただし、現時点（ArcGIS Desktop 10.6、ArcGIS Pro 2.2）で Z 係数パラメーターを持つすべてのツールに［メソッド］はが追加されている訳ではありません。たとえば陰影起伏ツールはまだ平面しか考慮されませんので、上記の計算を行うか、距離が正確な投影座標系を使用する必要があります。
その他の活用方法
「1度＝XXm」が判れば、ラスターデータのエクスポートにも活用できます。地理座標系のラスターレイヤーを右クリックし、データのエクスポートを実行すると、セルサイズの最適値に毎回悩みます。ArcGIS が適当なデフォルト値を設定してくれますがそれが妥当なものとは限りません。
換算表を使うと、たとえば 10m 解像度のメッシュで緯度 30度なら 0.0001036 を指定すれば良いことが判ります。
まとめ
地理座標系やメルカトルのような距離が正確でない投影座標系の高さを処理するツールを実行する際は、Z 係数に注意しましょう。地理座標系の場合は上記の表に基づき Z 係数を設定するか、ツールが対応している場合はメソッドで［測地線］を指定してください。
ArcScene でも同様の考慮が必要で XY 座標と Z座標の単位はそろえなければ正しい表示になりません。学生の頃である ArcGIS 8.x の時代だと 3D 表示は ArcScene でしか実現できませんでした。ArcScene で地理座標系のデータをうっかりそのまま立ち上げようとすると、Z値が度相当（100m の建物が 100度）となってしまい、途方もない縦長の表示にしてしまう失敗をよくやっていました。我々はこれを「鉛筆シーン」と命名していましたが、ArcGlobe や ArcGIS Pro のシーンは地理座標系のデータであっても高さを正しく表示してくれるので最近「鉛筆シーン」を見かけることはなくなりました。
