JP7440686B1 - 差分処理方法、差分処理装置及びプログラム - Google Patents

Abstract

【課題】データ間の位置ずれの影響を低減してより適切に標高差を算出することができる差分処理方法、差分処理装置及びプログラムを提供する。【解決手段】差分処理方法は、２つの標高値データの一方における単位領域ごとの斜面の方位を特定する方位特定ステップ、２つの標高値データにおける単位領域ごとの標高差を算出する標高差算出ステップ、標高差の代表値を方位ごとに求めて方位別代表値とし、方位と方位別代表値との関係に基づいて、方位に対する標高差の変化を近似した曲線を変化パターンとして算出するパターン算出ステップ、単位領域ごとに、斜面の方位と変化パターンとに応じて単位領域の標高差を補正する補正ステップ、を含む。【選択図】図６

Description

この発明は、差分処理方法、差分処理装置及びプログラムに関する。

航空機からレーザ計測により地表面の高さを計測する航空レーザ測量技術がある。異なる２タイミングにおける地表面の高さを比較することで、斜面の崩壊、侵食や土砂の堆積といった高さの変化を伴う地形変化、特に土砂災害の発生状況に係る情報が得られる。

このときに、２タイミングの地表面の高さに対応付けられる地表の水平位置の情報にずれがあると、正確な地形変化の情報が得られない。特に、傾斜面で水平位置がずれていると、地表の傾斜に応じて元々高さが異なる地点の２タイミングの高さを比較することになり、オフセット値が生じる。特許文献１では、斜面の方位ごとに２タイミングの標高の平均的な差をオフセット値として算出し、このオフセット値を標高差から差し引くことで標高差を補正する技術が開示されている。

特許第７１２９５３６号公報

しかしながら、実際に地形変化が発生して標高の変化が生じている領域の計測結果もオフセット値の算出に利用されると、補正結果が適切な値から外れやすい。これにより、結果的にデータ間の位置ずれの影響が適切に低減されない場合があるという課題がある。

この発明の目的は、データ間の位置ずれの影響を低減してより適切に標高差を算出することができる差分処理方法、差分処理装置及びプログラムを提供することにある。

上記目的を達成するため、本開示は、
２つの標高値データの一方における単位領域ごとの斜面の方位を特定する方位特定ステップ、
前記２つの標高値データにおける前記単位領域ごとの標高差を算出する標高差算出ステップ、
前記標高差の代表値を前記方位ごとに求めて方位別代表値とし、前記方位と前記方位別代表値との関係に基づいて、前記方位に対する標高差の変化を近似した曲線を変化パターンとして算出するパターン算出ステップ、
前記単位領域ごとに、当該単位領域の方位と前記変化パターンとに応じて前記単位領域の標高差を補正する補正ステップ、
を含む差分処理方法である。

本開示によれば、データ間の位置ずれの影響を低減してより適切に標高差を算出することができるという効果がある。

第１実施形態の情報処理装置の機能構成を示すブロック図である。 ２つの標高値データにおける位置ずれについて説明する図である。 斜面の方位角に対する平均標高差の分布の例を示す図である。 斜面の方位の算出について説明する図である。 標高差のばらつきについて説明する図である。 標高差分布出力処理の制御手順を示すフローチャートである。 水平ずれ量と垂直ずれ量との関係を説明する図である。 第２実施形態の標高差分布出力処理の制御手順を示すフローチャートである。

以下、本発明の実施の形態を図面に基づいて説明する。
［第１実施形態］
図１は、第１実施形態の差分処理装置である情報処理装置１の機能構成を示すブロック図である。

情報処理装置１は、通常のコンピュータ（ＰＣ）であってよく、ＣＰＵ１１（Central Processing Unit）と、記憶部１２と、通信部１３と、表示部１４と、操作受付部１５などを備える。

記憶部１２は、揮発性メモリ（ＲＡＭ；Random Access Memory）と不揮発性メモリとを有する。ＲＡＭは、ＣＰＵ１１に作業用のメモリ空間を提供し、一時データを記憶する。一時データには、外部から取得した標高値データ１２２やその処理データなどが含まれる。不揮発性メモリは、プログラム１２１や設定データなどを記憶する。不揮発性メモリは、例えば、フラッシュメモリ及び／又はＨＤＤ（Hard Disk Drive）であってよい。

通信部１３は、インターネット回線やＬＡＮ（Local Area Network）などのネットワークを経由して外部と通信を行う。外部には、標高値データを記憶するデータベース装置などが含まれていてもよい。

表示部１４は、表示画面を有し、ＣＰＵ１１の制御に基づいて表示を行う。表示画面は、特には限られないが、例えば、液晶画面や有機ＥＬ（Electro-Luminescent）画面などである。

操作受付部１５は、ユーザからの操作を受け付けて操作信号をＣＰＵ１１に出力する。操作受付部１５は、特には限られないが、例えば、キーボードやポインティングデバイス（マウスなど）を有する。

次に、本実施形態の差分処理方法について説明する。
情報処理装置１は、２つの標高値データ１２２を取得して記憶部１２に記憶させ、同一の水平位置での標高値Ｚの差分値（標高差ｄＺ）の分布を算出して表示する。２つの標高値データ１２２は、例えば、互いに異なる時期（タイミング）に計測されたものである。これらは、例えば、データベース装置などに記憶されているデータの計測日時を表示部１４により一覧表示させ、当該リスト表示内からユーザが操作受付部１５を介して行う入力操作により選択されることにより取得されるのであってもよい。各標高値データ１２２は、標高差ｄＺの分布を出力する対象となる共通のエリアが含まれていればよく、取得されたデータに対して更に当該共通のエリア（出力対象エリア）の設定が行われる。

標高値データ１２２は、適宜なサイズで二次元マトリクス状に区切られた単位領域のそれぞれについて三次元座標が定められたものである。このような標高値データ１２２は、例えば、ＤＥＭデータ（Digital Elevation Model；数値標高モデルなどと呼ばれる）として知られているものであってもよい。このような単位領域は、グリッド、セル又はメッシュなどと称されるが、以下ではメッシュと称する。地表の標高値Ｚは、上空からのレーザ計測、すなわち、航空機などから照射したレーザ光が反射して戻ってくるまでの時間によって得られる。レーザ計測では、レーザ光の反射点（計測点）の位置が三次元座標で得られる。しかしながら、樹木や人工物などによりレーザ光が地表面に届かずに反射したり、吸収、減衰されて反射波が戻ってこなかったりする点もあることから、地表面の計測は、その間隔や密度が不均一となり得る。ＤＥＭデータは、各メッシュ内の計測点の計測データが集約、平均化されて、等間隔かつ等密度なデータとなっている。メッシュの２辺は、緯度経度に沿っていてよく、例えば、各辺が１メートルであってもよい。あるいは、メッシュの各辺は、緯度方向についての長さと経度方向についての長さとが異なっていてもよい。出力対象エリア内の各メッシュは、適宜な座標又は識別番号などにより識別可能とされればよい。

これら２つの標高値データ１２２の水平位置が一致するメッシュの標高の差分により標高差ｄＺを得ることで、原理上は２回の計測日時の間での標高の変化量が得られる。新しい方の計測データから古い方の計測データを差し引くことで、標高差ｄＺの正負と変化方向の符号とが一致する、すなわち、正の値が堆積や流入を表し、負の値が侵食や流出を表す。反対に古い方の計測データから新しい方の計測データが差し引かれてもよい。この場合には、結果出力の前に最終的に符号が反転されてもよい。

しかしながら、計測される位置は、２回の計測の間で必ずしも正確に一致しない。すなわち、２つの標高値データ１２２において水平位置（Ｘ，Ｙ）が同一とされているメッシュであっても実際には水平方向にずれている場合がある。また、標高値Ｚが垂直方向にずれている場合がある。

図２は、２つの標高値データ１２２における位置ずれについて説明する図である。ずれには、図２（ａ）の破線と実線とにより示すような垂直方向へのずれと、図２（ｂ）の破線と実線とで示すような水平方向へのずれとが含まれ得る。これらのずれは、例えば、計測時の位置の特定誤差によるものと考えられる。また例えば、地殻変動によるずれも含まれる場合がある。これらのずれの量は微小であるが、出力対象範囲の標高差全体にオフセット値を生じさせるため地形変化を読み取りづらくさせる。そのため、標高差からずれによるオフセット値を差し引く補正を行うことが望ましい。

局所的な地形の変化、例えば、土砂崩れなどを検出したい場合に、全体に上記のようなずれが生じていると、空間全体に非ゼロの標高差ｄＺの領域が広がることになる。その結果、本来の検出対象である地形変化を一見して見分けるのが難しくなる。したがって、地形変化によらないこのような標高差ｄＺを低減させる必要がある。

上記のうち、図２（ａ）に示したように、垂直ずれ量ｄＨは、出力対象エリア内で一律な値である。一方、水平方向へのずれがあっても、傾斜面の方位角φが水平ずれ量ｄＬのずれ方向に垂直な方向の場合、すなわち、図２（ｂ）の表示面に垂直な方向（前後方向）への位置ずれでは、標高差が生じない。また、方位角φが水平ずれ量ｄＬのずれ方向に平行な方向の場合、すなわち、図２（ｂ）の表示面に平行な方向（左右方向）への位置ずれでは、最も大きな標高差が生じる。このように水平方向のずれに応じた標高差ｄＺは、斜面の方位角φと水平のずれ方向との関係に依存する。

このように、２つの標高値データ１２２の間の標高差には、水平方向へのずれに起因するオフセット値であり方位角φに依存した成分と、垂直方向へのずれに起因する一律なオフセット値とが加算された値となっている。変化パターンは、簡易的には、方位角φごとに標高差ｄＺの平均値を平均標高差ｄＺａ（φ）として算出することで得られる。

平均値をとる範囲となる上記出力対象エリアにおいて、実際の地形変化を生じている範囲が占める割合が大きいと、求める平均値において当該地形変化による標高差が支配的となってしまい、平均値を差し引く補正によって当該地形変化による標高差が相殺されてしまう。したがって、本来の検出対象の地形変化が相殺されないように、出力対象エリアは、地形変化を生じているメッシュの数よりも数桁（３－４桁）以上多いメッシュを含むように定められるのがよい。一方で、大地震による地殻変動などで広範囲な変化が非一様に含まれるほど広大な範囲が定められない方がよい。したがって、出力対象エリアは、例えば、数十～百キロメートル四方よりは小さく定められるのがよい。例えば、法律などで定められた図郭の範囲、例えば、１－２ｋｍ四方の矩形範囲が出力対象エリアとして定められてもよい。また、出力対象エリアは、地形の変化以外の要因で標高が変化するエリア、例えば、川や湖の水面などを除外して定められてもよい。

図３は、斜面の方位角φに対する平均標高差ｄＺａ（φ）の分布の例を示す図である。
実線で示すように、平均標高差ｄＺａ（φ）は、概ね正弦波形状を有し、３６０度で一周期の変化パターンを示している。しかしながら、部分的に正弦波形状から外れている部分（スパイク状の部分）が存在している。スパイク状の部分には、実際の地形変化による標高差が含まれていると考えられる。

情報処理装置１は、平均標高差ｄＺａ（φ）のφに応じた分布から、スパイク状の部分を除いた標高差の変化を近似した曲線を、変化パターンとして求める。変化パターンは、例えば、上記のように周期２π（３６０度）の正弦関数などで表され得る。あるいは、変化パターンは、３次～５次程度の多項式近似で表されてもよい。多項式近似のパラメータは、例えば最小二乗法などのフィッティングにより、すなわち回帰的に求められてもよいし、元の差分値の分布に対してローパスフィルタを適用することで求められてもよい。

図４は、斜面の方位の算出について説明する図である。
以下では、メッシュの二次元配列方向に沿った水平面上の２軸方向をｘ軸、ｙ軸とし、垂直方向（高さ方向）をｚ軸とする。
斜面の方位は、斜面の法線方向の水平成分の向きにより表される。すなわち、斜面の方位は、当該斜面のｘ軸方向についての傾きとｙ軸方向についての傾きの大きさ（正負を含む）の比により特定される。
メッシュＡ０（Ａ０は、出力対象エリア内の任意のメッシュである）の傾きは、当該メッシュＡ０を中心とする３×３のメッシュＡ０～Ａ８の標高値Ｚ（Ａｋ）（ｋ＝０～８）に基づいて算出される。例えば、ｘ方向（例えば経度方向）についての標高値の変化は、数式１～６に基づいて数式７のように求められる。
Ｚ（Ｐ０）＝（Ｚ（Ａ１）＋Ｚ（Ａ２））／２ … （数式１）
Ｚ（Ｐ１）＝（Ｚ（Ａ２）＋Ｚ（Ａ３））／２ … （数式２）
Ｚ（Ａ０ｘ＋）＝（Ｚ（Ｐ０）＋Ｚ（Ｐ１））／２
＝（Ｚ（Ａ１）＋２・Ｚ（Ａ２）＋Ｚ（Ａ３））／４ … （数式３）
Ｚ（Ｐ２）＝（Ｚ（Ａ５）＋Ｚ（Ａ６））／２ … （数式４）
Ｚ（Ｐ３）＝（Ｚ（Ａ６）＋Ｚ（Ａ７））／２ … （数式５）
Ｚ（Ａ０ｘ－）＝（Ｚ（Ｐ２）＋Ｚ（Ｐ３））／２
＝（Ｚ（Ａ５）＋２・Ｚ（Ａ６）＋Ｚ（Ａ７））／４ … （数式６）
ｄＺ（Ａ０）／ｄｘ＝（Ｚ（Ａ０ｘ＋）＋Ｚ（Ａ０））／２－（Ｚ（Ａ０）＋Ｚ（Ａ０ｘ－））／２
＝（Ｚ（Ａ０ｘ＋）－Ｚ（Ａ０ｘ－））／２ … （数式７）
また、ｙ方向（例えば緯度方向）についての標高値の変化は、例えば、数式８～１３に基づいて数式１４のように求められる。
Ｚ（Ｐ４）＝（Ｚ（Ａ１）＋Ｚ（Ａ８））／２ … （数式８）
Ｚ（Ｐ５）＝（Ｚ（Ａ８）＋Ｚ（Ａ７））／２ … （数式９）
Ｚ（Ａ０ｙ＋）＝（Ｚ（Ｐ４）＋Ｚ（Ｐ５））／２
＝（Ｚ（Ａ１）＋２・Ｚ（Ａ８）＋Ｚ（Ａ７））／４ … （数式１０）
Ｚ（Ｐ６）＝（Ｚ（Ａ４）＋Ｚ（Ａ５））／２ … （数式１１）
Ｚ（Ｐ７）＝（Ｚ（Ａ３）＋Ｚ（Ａ４））／２ … （数式１２）
Ｚ（Ａ０ｙ－）＝（Ｚ（Ｐ６）＋Ｚ（Ｐ７））／２
＝（Ｚ（Ａ３）＋２・Ｚ（Ａ４）＋Ｚ（Ａ５））／４ … （数式１３）
ｄＺ（Ａ０）／ｄｙ＝（（Ｚ（Ａ０ｙ＋）＋Ｚ（Ａ０））／２－（Ｚ（Ａ０）＋Ｚ（Ａ０ｙ－））／２
＝（Ｚ（Ａ０ｙ＋）－Ｚ（Ａ０ｙ－））／２ … （数式１４）

数式７及び数式１４に基づいて、各メッシュの斜面の方位角φ（傾斜面の法線方向の方位）は、以下の数式１５で表される。
φ＝１８０／π・ｔａｎ^－１（（ｄＺ（Ａ０）／ｄｘ）／（ｄＺ（Ａ０）／ｄｙ）） … （数式１５）
逆正接により求められる方位角φは、±９０の範囲であるので、ｄＺ（Ａ０）／ｄｘの符号に応じて±１８０の範囲に拡張されればよい。また、マイナス側の数値に３６０を加算して一周期ずらし、０～３６０の範囲としてもよい。このようにして斜面の方位（方位角）が表される。ここでは、方位角φの基準方向（φ＝０）は真北方向であるが、これに限られない。また、斜面の方位角φは、その他の周知の方法で求められてもよい。

また、数式７及び数式１４に基づいて、各メッシュの斜面の傾斜角度θが以下の数式１６で表される。
θ＝１８０／π・ｔａｎ^－１（（ｄＺ（Ａ０）／ｄｘ）^２＋（ｄＺ（Ａ０）／ｄｙ）^２）^１／２） … （数式１６）

なお、（（ｄＺ（Ａ０）／ｄｘ）^２＋（ｄＺ（Ａ０）／ｄｙ）^２）^１／２がほぼゼロである（基準値未満の）メッシュは、平地であるとして、上記方位角及び傾斜角度の算出対象、並びに以下の処理から除外されてよい。

このようにして得られた斜面の方位角φごと（ここでは、小数点以下を四捨五入した１度単位）にメッシュがグループ分けされる。方位角φごとに、メッシュのグループにおける標高差ｄＺの平均値である平均標高差ｄＺａ（φ）が求められる。平均標高差ｄＺａ（φ）の分布が曲線近似されて、上記のように正弦波形などで近似した曲線ｆ（φ）が得られる。この曲線ｆ（φ）の値が最大の値となる方位角φが、２つの標高値データの間での水平ずれ方向φ０と推定される。なお、曲線ｆ（φ）の変化の中心がゼロからずれている場合（ｆ（φ）＝ｆ０（φ）＋ｄＨ。ｆ０（φ）は、０を中心とした変化）には、このオフセット値により垂直ずれ量ｄＨが推定される。

各メッシュで得られている標高差ｄＺから、当該メッシュの方位角φに応じた曲線ｆ（φ）の値を差し引くことで、数式１７のように標高差ｄＺを補正した補正標高差ｄＺｍが得られる。
ｄＺｍ＝ｄＺ－ｆ（φ） … （数式１７）

ここで、曲線ｆ（φ）で表される値は、各方位角φの地形変化がないメッシュにおける標高差ｄＺの平均的な値であり、標高差ｄＺには、メッシュごとにさらにばらつきがある。

図５は、標高差のばらつきについて説明する図である。
図５（ａ）に示すように、地形断面の傾斜角度θが一様ではない場合には、ある水平ずれ量ｄＬに対して生じる高さのずれ量が傾斜角度θに依存して変化する。すなわち、傾斜角度θが大きいほどずれ量が大きくなる。したがって、傾斜角度θによらず一律に曲線ｆ（φ）の値を差し引くと、傾斜角度θが小さい場合には、もともと現れていなかったずれ量に対する補正値が過大となって、むしろ不正確な非ゼロの標高差が生じることになる。

図５（ｂ）には、上記の傾斜面の向きが水平ずれ方向φ０と等しい各メッシュの標高差ｄＺを更に傾斜角度θごとに分類して、それぞれ平均標高差ｄＺａ（φ０、θ）（第１仮代表値）を求めた分布が示されている。なお、垂直ずれ量ｄＨがゼロではない場合には、平均標高差ｄＺａ（φ０、θ）から垂直ずれ量ｄＨの影響を差し引いておく。傾斜角度θの分類は、例えば、１度単位であってもよい。上記の通り、白丸で示した平均標高差ｄＺａ（φ０、θ）は、傾斜角度θが大きいほど大きくなり、傾斜角度θが０に近づくとゼロに漸近していく。

この平均標高差ｄＺａ（φ０、θ）から曲線ｆ（φ）のφ＝φ０における値を差し引いた仮補正標高差ｄＺａｍ（φ０、θ）（第２仮代表値）が黒丸で示されている。仮補正標高差ｄＺａｍ（φ０、θ）は、傾斜角度θが３８度以下で正負が逆転している。また、傾斜角度θが２１度以下（２２度未満）では、元の平均標高差ｄＺａ（φ０、θ）よりも仮補正標高差ｄＺａｍ（φ０、θ）の絶対値が大きくなっている。このように補正の結果、むしろずれ量が大きくなる場合には、補正を行う意味がない。情報処理装置１は、曲線ｆ（φ）による補正によりずれ量が小さくなる傾斜角度θの範囲でのみ上記補正値を適用する。

図６は、本実施形態の情報処理装置１で実行される標高差分布出力処理のＣＰＵ１１による制御手順を示すフローチャートである。本実施形態の差分処理方法を含むこの標高差分布出力処理は、例えば、操作受付部１５への所定の入力操作などに応じてプログラム１２１が実行されることで行われる。

ＣＰＵ１１は、比較対象の標高値データを２つ取得する（ステップＳ１０１）。ＣＰＵ１１は、当該２つの標高値データに対し、標高差分布を出力する共通の出力対象エリアを設定する（ステップＳ１０２）。

ＣＰＵ１１は、方位特定手段として、一方の標高値データ、例えば、古い方の標高値データにおける各メッシュ（単位領域ごと）の傾斜面の方位角φを算出し、傾斜角算出手段として、傾斜角度θを算出する（ステップＳ１０３；方位特定ステップ、傾斜角算出ステップ）。ＣＰＵ１１は、標高差算出手段として、２つの標高値データにおける同一メッシュの標高値の差分により標高差ｄＺを算出する（ステップＳ１０４；標高差算出ステップ）。

ＣＰＵ１１は、出力対象エリアの各メッシュを方位ごとにグループ分けして、方位ごとに、グループ分けされたメッシュで算出されている標高差の平均標高差ｄＺａ（φ）を方位別代表値として算出する（ステップＳ１０５）。ＣＰＵ１１は、パターン算出手段として、得られた平均標高差ｄＺａ（φ）の方位角φに応じた変化パターンを曲線ｆ（φ）で近似して求める（ステップＳ１０６；パターン算出ステップ）。近似される曲線ｆ（φ）は、上記のように、正弦関数若しくは余弦関数（オフセット値である垂直ずれ量ｄＨが定数項として含まれていてもよい）、又は３～５次の多項式などである。

ＣＰＵ１１は、曲線ｆ（φ）の値が極値、ここでは絶対値が最大値を取る方位角φを水平ずれ方向φ０（第１方位）として特定する（ステップＳ１０７）。ＣＰＵ１１は、傾斜面の方位が水平ずれ方向φ０と等しいメッシュの標高差ｄＺを傾斜角度θごとに平均した平均標高差ｄＺａ（φ０、θ）をそれぞれ第１仮代表値として算出する（ステップＳ１０８）。

ＣＰＵ１１は、得られた平均標高差ｄＺａ（φ０、θ）の各々から曲線ｆ（φ０）の値を差し引いて仮補正標高差ｄＺａｍ（φ０、θ）（第２仮代表値）を算出し、当該仮補正標高差ｄＺａｍ（φ０、θ）と、元の平均標高差ｄＺａ（φ０、θ）とを比較する。ＣＰＵ１１は、比較の結果、補正により標高差が大きくならない、すなわち、元の平均標高差の絶対値｜ｄＺａ（φ０、θ）｜＜仮補正標高差の絶対値｜ｄＺａｍ（φ０、θ）｜である傾斜角度θの範囲を不適切補正範囲として定める（ステップＳ１０９；範囲特定ステップ）。ＣＰＵ１１は、補正手段として、各メッシュの標高差ｄＺのうち傾斜角度θが不適切補正範囲内にないものについて、当該メッシュの斜面の方位角φに応じた曲線ｆ（φ）の値を差し引くことで、それぞれ補正標高差ｄＺｍを算出する（ステップＳ１１０；補正ステップ）。ＣＰＵ１１は、不適切補正範囲内にあるメッシュについては、元の標高差ｄＺをそのまま補正標高差ｄＺｍに代入してもよい。

ＣＰＵ１１は、得られた補正標高差ｄＺｍをその位置データとともに出力する（ステップＳ１１１）。出力は、画像データであってもよい。例えば、予め、補正標高差ｄＺｍの値と色や濃淡との対応関係が設定されて、ＣＰＵ１１は、各メッシュの補正標高差ｄＺｍの値を色や濃度に変換して出力してもよい。変換された色や濃度を各メッシュの位置に対応する画像上の位置に割り当てることで、補正標高差ｄＺｍの分布に応じたカラー画像や濃淡画像が得られる。また、これらの画像は、同一領域の航空写真画像と重ねられて透過表示されてもよい。この出力において、ＣＰＵ１１は、単に補正標高差ｄＺｍの分布を出力するだけではなくてもよい。例えば、補正標高差ｄＺｍが基準を超えている範囲の輪郭を検出して、当該輪郭の位置情報などが併せて出力されてもよい。また、ＣＰＵ１１は、この輪郭内での補正標高差ｄＺｍを積分して、合計の土砂の移動量を概算出力してもよい。ＣＰＵ１１は、標高差分布出力処理を終了する。

［第２実施形態］
次に、第２実施形態の情報処理装置１について説明する。
第２実施形態の情報処理装置１の機能構成は、第１実施形態の情報処理装置１と同一であるので、同一の符号を用いることとして説明を省略する。

第２実施形態の情報処理装置１における差分値補正の処理内容について説明する。
上記第１実施形態で示したように、標高差ｄＺの方位角φごとの平均標高差ｄＺａ（φ）から得られた曲線ｆ（φ）が極大となる方向により、水平ずれ方向φ０が特定される。

図７は、水平ずれ方向φ０を含む垂直面における破線斜面から実線斜面への水平ずれ量ｄＬと垂直ずれ量ｄＨとの関係を説明する図である。
斜面の方位角φが水平ずれ方向φ０と等しいメッシュの平均傾斜角度θ０（φ０）を傾斜角度の代表値として、この平均傾斜角度θ０（φ０）に基づいて、水平ずれ量ｄＬは、以下の数式１８のように推定される。
ｄＬ＝（ｄＺａ（φ０）－ｄＨ０）／ｔａｎ（θ０（φ０）） … （数式１８）

斜面の方位角φが水平ずれ方向φ０と等しいメッシュにおいて、任意の傾斜角度θの斜面では、垂直方向の位置ずれ量ｄｄＺは、数式１９で表される。
ｄｄＺ（φ０、θ）＝ｄＬ・ｔａｎ（θ）＋ｄＨ０ … （数式１９）

任意の斜面の方位角φでは、上記水平方向ずれに対する位置ずれ量であるｄＬ・ｔａｎ（θ）のφ＝φ０の面へ投影された分の位置ずれが生じる。したがって、任意のθ、φの組み合わせのメッシュにおける垂直方向の位置ずれ量ｄｄＺ
は、数式２０で表される。
ｄｄＺ（φ、θ）＝ｄＬ・ｃоｓ（φ－φ０）・ｔａｎ（θ）＋ｄＨ０ … （数式２０）

もとの標高差ｄＺから上記の標高差ｄｄＺ（φ、θ）を差し引いて補正された補正標高差ｄＺｍが再算出されて得られる。

図８は、第２実施形態の情報処理装置１が実行する標高差分布出力処理の制御手順を示すフローチャートである。
この標高差分布出力処理は、第１実施形態の標高差分布出力処理におけるステップＳ１０８～Ｓ１１０の処理の代わりにステップＳ１２１～Ｓ１２３の処理が含まれている。その他の処理は同一であり、同一の処理内容は同一の符号で示されることとして説明は省略される。

ステップＳ１０７の処理の後、ＣＰＵ１１は、斜面方位が水平ずれ方向φ０であるメッシュの平均傾斜角度θ０（φ０）を算出する（ステップＳ１２１；基準値特定ステップ）。ＣＰＵ１１は、上述の数式１８により、水平ずれ量ｄＬを算出する（ステップＳ１２２；位置ずれ算出ステップ）。ＣＰＵ１１は、各メッシュの標高差ｄＺから数式２０で得られるｄＺ（φ、θ）を差し引いて、補正標高差ｄＺｍを算出する（ステップＳ１２３）。ＣＰＵ１１の処理は、ステップＳ１１１へ移行する。

このように、位置ずれの補正量を斜面の傾斜に応じて調整することで、一律な変化の場合と比較して、より柔軟な補正量が得られる。

なお、本発明は、上記実施の形態に限られるものではなく、様々な変更が可能である。
例えば、上記実施の形態では、方位角φや傾斜角度θを１度単位で分割するものとして説明したが、これに限られない。他の角度単位であってもよい。例えば、１度単位の分割では各角度のデータ数が足りない場合には、分割する角度幅が広く３－５度程度に定められてもよい。

また、上記実施の形態では、正弦関数又は３～５次の多項式による近似により曲線ｆ（φ）を求めるとしたが、平均標高差ｄＺａ（φ）の細かい変動を反映し過ぎない範囲で適切に近似可能な他の関数であってもよい。多項式である場合には、曲線ｆ（φ）は、３次以上の関数であることが好ましい。

また、斜面の方位角φ及び傾斜角度θの求め方は、上記説明したものに限られない。他の公知の計算方法で求められてもよい。

また、上記実施の形態では、各方位角φごと、各傾斜角度θごとのメッシュの標高差の代表値として平均値を用いるものとして説明したが、これに限られない。代表値は、例えば、中央値や最頻値などであってもよい。ただし、最頻値が用いられる場合には、傾斜角度が基準値未満、すなわち平坦又は平坦に近いメッシュの頻度は除外されてもよい。

また、上記実施の形態では、水平ずれ方向φ０に斜面が向いているメッシュの傾斜角度θに応じた平均標高差ｄＺａ（φ０、θ）に基づいて、補正の有無が定められた。しかしながら、他の特定の方位角φを向いた斜面のメッシュにより補正の有無が定められてもよい。方位が水平ずれ方向φ０から垂直な方向に近くなると、上記の通り、標高差ｄＺに見かけ上のずれが生じなくなっていくので、水平ずれ方向φ０に近い向きであることが好ましい。また、水平ずれ方向φ０と反対向きである場合には、正負が反転するので、これに合わせて基準が定められる。

また、曲線ｆ（φ）が１８０度ごとに反対称の場合には、不適切補正範囲の設定時などには、１８０度反対の２方向のデータが符号を反転させて足し合わせて利用されてもよい。すなわち、可能な場合には、第１方位は、一方向に限られなくてもよい。

また、曲線ｆ（φ）が多項式近似されたものであって、必ずしも１８０度ごとに反対称ではない場合には、極値を取る２つの方向が特定されてもよい。特定された２つの極値を取る方位角φが１８０度差ではない場合には、一方の方位角と他方の方位角の１８０度反対側の方位角との平均方向が水平ずれ方向φ０として定められてもよい。

また、上記実施の形態では、補正により標高差が大きくなる傾斜角度では補正を行わない例を説明したが、これに限られない。情報処理装置１は、平均標高差ｄＺａ（φ０）が小さい場合、例えば、ＤＥＭの計測誤差などと比較して小さい場合などには、そもそも補正を行わないように判定してもよい。

また、不適切補正範囲の設定は、曲線ｆ（φ）を求めずに、直接ある方位φ０及び傾斜角度θの各々における平均標高差ｄＺａ（φ０、θ）（第１仮代表値）と、補正標高差ｄＺａｍ２（φ０）＝ｄＺａ（φ０、θ）－ｄＺａ（φ０）（第２仮代表値）との各絶対値の大小に応じて定められてもよい。

また、上記の標高差分布出力処理は、単一の情報処理装置による処理に限られない。複数の情報処理装置により分散処理されてもよい。

また、以上の説明では、本発明の差分値処理の制御に係るプログラム１２１を記憶するコンピュータ読み取り可能な媒体としてＨＤＤ、フラッシュメモリなどの不揮発性メモリなどからなる記憶部１２を例に挙げて説明したが、これらに限定されない。その他のコンピュータ読み取り可能な媒体として、ＭＲＡＭなどの他の不揮発性メモリや、ＣＤ－ＲＯＭ、ＤＶＤディスクなどの可搬型記録媒体を適用することが可能である。また、本発明に係るプログラムのデータを、通信回線を介して提供する媒体として、キャリアウェーブ（搬送波）も本発明に適用される。
その他、上記実施の形態で示した具体的な構成、処理動作の内容及び手順などは、本発明の趣旨を逸脱しない範囲において適宜変更可能である。本発明の範囲は、特許請求の範囲に記載した発明の範囲とその均等の範囲を含む。

以上のように、本実施形態の差分処理方法は、次のステップを含む。２つの標高値データ１２２の一方におけるメッシュごとの斜面の方位角φを特定する方位特定ステップ。２つの標高値データ１２２におけるメッシュごとの標高差ｄＺを算出する標高差算出ステップ。標高差ｄＺの代表値としての平均値を方位角φごとに求めて方位別代表値（平均標高差ｄＺａ（φ））とし、方位角φと平均標高差ｄＺａ（φ）との関係に基づいて、方位角φに対する標高差ｄＺの変化を近似した曲線ｆ（φ）を変化パターンとして算出するパターン算出ステップ。メッシュごとに、当該メッシュの方位角φと曲線ｆ（φ）とに応じて各メッシュの標高差ｄＺを補正して補正標高差ｄＺｍを得る補正ステップ。
このように、２つの標高値データ１２２の位置ずれに起因する見かけ上の標高差を差し引く補正をする際に、差分処理方法では、水平位置ずれに係る斜面の方位角φに対する依存性を考慮する。この方位依存性において、水平位置ずれに対して想定される略正弦波形状の標高差（ずれ）の傾向から外れる標高差は、見かけ上のものではなく、実際に生じているものである可能性がある。このような本来検出対象の標高差を排除するのは、本末転倒である。本実施形態の差分処理方法では、標高差の各方位についての変化を３６０度の範囲全体で近似的に表した曲線ｆ（φ）を変化パターンとして求める。これにより、実際に生じている標高差が差し引かれて正しくユーザに認識されないおそれを低減させることができる。よって、この差分処理方法は、データ間の位置ずれの影響を低減して、より適切に標高差を算出することができる。このようにより適切に算出された標高差は、例えば、土砂崩れ、侵食、堆積、堰堤への土砂の流入など、自然災害の把握や防災事業の計画などに好適に利用され得る。

また、差分処理方法は、さらに、以下のステップを含む。２つの標高値データ１２２の一方におけるメッシュごとに、斜面の傾斜角度θを算出する傾斜角算出ステップ。水平ずれ方向φ０といったある第１方位における標高差ｄＺの代表値として平均標高差ｄＺａ（φ０、θ）を傾斜角度θごとに求め、第１方位（水平ずれ方向φ０）における傾斜角度θごとに、変化パターン（曲線ｆ（φ））に基づいて仮補正標高差ｄＺａｍ（φ０、θ）を算出し、平均標高差ｄＺａ（φ０、θ）と仮補正標高差ｄＺａｍ（φ０、θ）との絶対値大小に応じて補正が適切ではない傾斜角度θの範囲を特定する範囲特定ステップ。上記補正ステップでは、傾斜角度θが上記範囲にあるメッシュの標高差ｄＺを補正しない。
このように、補正のしすぎにより、特に傾斜角度θが小さい領域でむしろ非ゼロの見かけ上の標高差が現れるのを抑制することができる。つまり、この差分処理方法によれば、平坦地などで広く見かけ上の標高差が生じて、実際の標高差の発生状況が見にくくなるのを抑制することができる。よって、本実施形態の差分処理方法は、より適切に見かけ上の標高差を低減させることができる。

また、第１方位は、変化パターンを示す曲線ｆ（φ）の値が極値をとる方位であってもよい。この方位は、垂直ずれ量ｄＨを差し引いた曲線ｆ０（φ）であっても同じである。すなわち、水平ずれの影響が大きい方位を基準として基準方向が定められる。また、位置ずれの傾斜角度θに対する依存性を見る場合にも、方位角φに依存した変化パターンの振幅が大きい方が精度よく求めやすくなる。

一方で、差分処理方法は、以下のステップを含むものであってもよい。２つの標高値データ１２２の一方におけるメッシュごとの斜面の方位を特定する方位特定ステップ。２つの標高値データ１２２におけるメッシュごとの標高差ｄＺを算出する標高差算出ステップ。方位角φごとの標高差ｄＺの代表値（平均標高差ｄＺａ（φ））を変化パターンとして算出するパターン算出ステップ。メッシュごとに、当該メッシュの方位角φと変化パターンｄＺａ（φ）とに応じてこのメッシュの標高差ｄＺを補正して補正標高差ｄＺｍを得る補正ステップ。２つの標高値データ１２２の一方におけるメッシュごとに、斜面の傾斜角度θを算出する傾斜角算出ステップ。ある第１方位（水平ずれ方向φ０）における標高差ｄＺの代表値として平均標高差ｄＺａ（φ０、θ）を傾斜角度θごとに第１仮代表値として求め、水平ずれ方向φ０における傾斜角度θごとに、変化パターンに基づいて補正された標高差の第２仮代表値（補正標高差ｄＺａｍ２（φ０、θ）＝ｄＺａ（φ０、θ）－ｄＺａ（φ０））を算出し、平均標高差ｄＺａ（φ０、θ）の絶対値と補正標高差ｄＺａｍ２（φ０、θ）の各絶対値との大小に応じて補正が適切ではない傾斜角度θの範囲を特定する範囲特定ステップ。補正ステップでは、傾斜角度θが、補正が適切ではない範囲にあるメッシュの標高差ｄＺを補正しない。
このように、補正しない傾斜角度θの範囲は、変化パターンとして近似曲線ｆ（φ）を算出しない場合でも定められ得る。差分処理方法は、データ間の位置ずれの影響を低減して、より適切に標高差を算出することができる。

また、差分処理方法は、さらに、以下のステップを含み得る。平均標高差ｄＺａ（φ）のうち絶対値が最大の方位を水平ずれ方向φ０として特定し、この水平ずれ方向φ０に斜面が向いたメッシュにおける平均傾斜角度θ０（φ０）を求める基準値特定ステップ。変化パターンを示す曲線ｆ（φ）の振幅（ｆ（φ０）－ｄＨ）と、平均傾斜角度θ０（φ０）とに基づいて、２つの標高値データ１２２間の水平ずれ量ｄＬを算出する位置ずれ算出ステップ。補正ステップでは、メッシュごとに、当該メッシュの方位角φ、傾斜角度θ及び水平ずれ量に基づいて、２つの標高値データ１２２の間での標高差を再算出することで補正標高差ｄＺｍを得る。
このように、方位角φに応じて変化する位置ずれ量を考慮する際に、更に斜面の傾斜角度θの影響も考慮する。したがって、この差分処理方法では、傾斜角度θに応じた補正のし過ぎや補正不足を低減させ、より確実に見かけ上の標高差を低減することができる。したがって、この差分処理方法は、２つの標高値データ間の地形変化がユーザに分かりやすいデータを得ることができる。

また、本実施形態の差分処理装置である情報処理装置１は、ＣＰＵ１１を備える。ＣＰＵ１１は、方位特定手段として、２つの標高値データ１２２の一方におけるメッシュごとの斜面の方位角φを特定する。ＣＰＵ１１は、標高差算出手段として、２つの標高値データ１２２におけるメッシュごとの標高差ｄＺを算出する。ＣＰＵ１１は、パターン算出手段として、標高差ｄＺの方位角φごとの平均値である平均標高差ｄＺａ（φ）を求め、方位角φと平均標高差ｄＺａ（φ）との関係に基づいて、方位角φに対する標高差ｄＺの変化を近似した曲線ｆ（φ）を変化パターンをとして算出する。ＣＰＵ１１は、補正手段として、メッシュごとに、当該メッシュの方位角φと曲線ｆ（φ）とに応じて標高差ｄＺを補正して補正標高差ｄＺｍを得る。
このような情報処理装置１によれば、標高値データ１２２間の位置ずれの影響を低減して、より適切に標高差を求めることができる。

また、上記差分処理方法に係るプログラム１２１をコンピュータ（情報処理装置１）にインストールして実行することで、容易に汎用のコンピュータで２つの標高値データ１２２間の位置ずれの影響を低減して、より適切な標高差を求めることができる。

１ 情報処理装置
１１ ＣＰＵ
１２ 記憶部
１２１ プログラム
１２２ 標高値データ
１３ 通信部
１４ 表示部
１５ 操作受付部
ｄＨ 垂直ずれ量
ｄＬ 水平ずれ量
ｄＺ 標高差
ｄＺａ 平均標高差
ｄＺａｍ 仮補正標高差
ｄＺｍ 補正標高差
θ 傾斜角度
θ０ 平均傾斜角度
φ 方位角
φ０ 水平ずれ方向

Claims (7)

1. ２つの標高値データの一方における単位領域ごとの斜面の方位を特定する方位特定ステップ、
   前記２つの標高値データにおける前記単位領域ごとの標高差を算出する標高差算出ステップ、
   前記標高差の代表値を前記方位ごとに求めて方位別代表値とし、前記方位と前記方位別代表値との関係に基づいて、前記方位に対する標高差の変化を近似した曲線を変化パターンとして算出するパターン算出ステップ、
   前記単位領域ごとに、当該単位領域の方位と前記変化パターンとに応じて前記単位領域の標高差を補正する補正ステップ、
   を含む差分処理方法。
2. 前記２つの標高値データの一方における単位領域ごとに、前記斜面の傾斜角度を算出する傾斜角算出ステップ、
   ある第１方位における前記標高差の代表値を前記傾斜角度ごとに第１仮代表値として求め、前記第１方位における前記傾斜角度ごとに、前記変化パターンに基づいて補正された前記標高差の第２仮代表値を算出し、前記第１仮代表値と前記第２仮代表値との絶対値の大小に応じて前記補正が適切ではない前記傾斜角度の範囲を特定する範囲特定ステップ、
   を含み、
   前記補正ステップでは、前記傾斜角度が前記範囲にある前記単位領域の前記標高差を補正しない、
   請求項１記載の差分処理方法。
3. ２つの標高値データの一方における単位領域ごとの斜面の方位を特定する方位特定ステップ、
   前記２つの標高値データにおける前記単位領域ごとの標高差を算出する標高差算出ステップ、
   前記方位ごとの前記標高差の代表値の変化パターンを算出するパターン算出ステップ、
   前記単位領域ごとに、当該単位領域の方位と前記変化パターンとに応じて前記単位領域の標高差を補正する補正ステップ、
   前記２つの標高値データの一方における単位領域ごとに、前記斜面の傾斜角度を算出する傾斜角算出ステップ、
   ある第１方位における前記標高差の代表値を前記傾斜角度ごとに第１仮代表値として求め、前記第１方位における前記傾斜角度ごとに、前記変化パターンに基づいて補正された前記標高差の第２仮代表値を算出し、前記第１仮代表値と前記第２仮代表値の大小に応じて前記補正が適切ではない前記傾斜角度の範囲を特定する範囲特定ステップ、
   を含み、
   前記補正ステップでは、前記傾斜角度が前記範囲にある前記単位領域の前記標高差を補正しない、
   差分処理方法。
4. 前記第１方位は、前記変化パターンにおいて極値をとる方位である、請求項２又は３記載の差分処理方法。
5. 前記２つの標高値データの一方における単位領域ごとに、前記斜面の傾斜角度を算出する傾斜角算出ステップ、
   前記方位別代表値のうち絶対値が最大の当該方位別代表値及びその第１方位を特定し、当該第１方位の前記単位領域の傾斜角度の代表値を求める基準値特定ステップ、
   前記変化パターンの振幅と、前記傾斜角度の代表値とに基づいて、前記２つの標高値データ間の水平ずれ量を算出する位置ずれ算出ステップ、
   を含み、
   前記補正ステップでは、前記単位領域ごとに、当該単位領域の前記方位、前記傾斜角度及び前記水平ずれ量に基づいて、前記２つの標高値データの間での標高差を再算出することで前記標高差を補正する、
   請求項１記載の差分処理方法。
6. ２つの標高値データの一方における単位領域ごとの斜面の方位を特定する方位特定手段、
   前記２つの標高値データにおける前記単位領域ごとの標高差を算出する標高差算出手段、
   前記標高差の代表値を前記方位ごとに求めて方位別代表値とし、前記方位と前記方位別代表値との関係に基づいて、前記方位に対する標高差の変化を近似した曲線を変化パターンとして算出するパターン算出手段、
   前記単位領域ごとに、当該単位領域の方位と前記変化パターンとに応じて前記単位領域の標高差を補正する補正手段、
   を備える差分処理装置。
7. コンピュータを、
   ２つの標高値データの一方における単位領域ごとの斜面の方位を特定する方位特定手段、
   前記２つの標高値データにおける前記単位領域ごとの標高差を算出する標高差算出手段、
   前記標高差の代表値を前記方位ごとに求めて方位別代表値とし、前記方位と前記方位別代表値との関係に基づいて、前記方位に対する標高差の変化を近似した曲線を変化パターンとして算出するパターン算出手段、
   前記単位領域ごとに、当該単位領域の方位と前記変化パターンとに応じて前記単位領域の標高差を補正する補正手段、
   として機能させるプログラム。

Images