2011年12月7日水曜日

Wavelet その2

メキシカンハット関数を使った、ウェーブレット解析(ただの演算ですが)の続きです。

昨日UPした通り、(1-x^2・・・から始めることにしました。
2次元関数はこのような感じです。


マニュアルではスケール S を1で固定し、DEM のグリッドサイズで波長が変わるといった内容です。λ=4グリッドといったところですね。実務的ですね。
x-aが影響するところまでを考えると、最低でも中心から左右に4グリッドずつ必要です。従って、フィルタは9*9のマトリックスになります。

早速、SAGA で User Defined Filter があるか調べてみました。結果、装備されていました。しかし、3*3までなので適用できません。仕方ないですね。
Surfer では9*9でも可能ですから、こちらでやってみましょう。ただし、Filter による演算結果は、最後にマトリックスの合計値で割る設定になっています。そのため、マニュアルの積分に合わせるには、後で Grid Math で合計値を掛け戻してやる必要があります。また、Filter の edge の処理設定もしないといけません。

演算結果を確認してみますと、あれ?といった感じの仕上がりです。全くマニュの図とは異なっています。係数のオーダーも大きく違います。チェックしましたが計算過程は問題ありません。
ということは、(1-x^2・・・ではなくて、(2-x^2が正しいということでしょうか?

(2-x^2・・・の2次元関数は以下のようになります。振幅の最大値は大きくなりますが、同じような形状です。


Surfer で演算し直し、できた GRD ファイルを SAGA で読み込んだ結果がこちら。


オーダーもあっています。(2-X^2・・・で問題ないようでした。
最終的に、地すべり地形を判読しやすいと実感したら、式の導出まで追いかけてみましょうか。このまま使うのは気持ち悪いですからね。


とりあえず、地上開度、ウェーブレット解析、両方ともクリアーです。
地上開度もよく読んでみますと、1グリッドで計算した場合も書かれていますので、最大値を取るアルゴを組まなくてもよくなっています。非常に実務的と言いますか、汎用ソフトの利用を前提に、適用しやすいよう単純化して書かれています。技術的に良いか悪いかは分かりませんが、裾野を広げることはできるでしょう。


あとは地上開度図とウェーブレット解析図の合成のみです。
こちらはまた後日。

0 件のコメント:

コメントを投稿