2013年12月19日

特異値分解を用いた連立一次方程式の計算

連立一次方程式を解くとき、逆行列を計算するのではなく、特異値分解をする方が良い、というのをどこかで見た気がしますので、試してみることにします。

今回例題として、高さマップから局所的な2次曲面を計算してみます。
高さ情報をz座標と見立てて、3次元の座標F(x,y,z)=0となる等値面と考えると、
3変数の2次曲面の方程式は、以下のように表すことができます。

svd_007.png

これは2次形式で表現するのが普通なのだと思いますが、とりあえず連立方程式の
形にしてみます。
ある注目画素を中心とするMxMの近傍画素の座標を、 x,y,zそれぞれ-1〜1の範囲に
変換します。そうすると、N=MxM個の連立方程式を作ることができます。

svd_008.png

これを行列の形で表すと、Ax=0の形になります。

svd_009.png

今回は、それぞれの行をノルム1に正規化したものを行列Aとします。
(正規化したほうが近似精度が上がったため)

svd_011.png

この場合、前回のホモグラフィ行列のときにも用いましたが、特異値分解して
特異値が最も小さい特異ベクトルが誤差の二乗が最小になる解になるらしいです。

svd_010.png

Eigenを用いた場合、最小特異値を持つ特異ベクトルは、Vの一番右端に
来ているようですので、そのベクトルがa,b,…,jの値となります。

これで10個のパラメータが求まりましたので、画像上ではx,yを決めてzを求める方法、
3次元空間上では視点とレイを決めて、レイと2次曲面の衝突判定をすることで、
計算した2次曲面で近似した場合の結果を得ることができます。

後、たまに近傍と大きく異なる推定結果の点が出てきましたので、10個のパラメータを
画素ごとに求めて、5x5のガウスフィルタをかけてパラメータの平均化?をしました。

実行結果
各画素における2次曲面を、近傍5x5の高さ情報から計算する、ということをやっています。
画像の拡大とか勾配・法線ベクトルの計算などを試してみないとこの結果だけでは特に意味は無いと思います。
(特異値分解で連立一次方程式を解くテストという意味合いしかありません)

2次曲面の特性上、エッジがはっきりしているようなものは近似しにくいと思いましたので、
入力画像には5x5のガウスフィルタを1回掛けてあります。

1. 円内が1.0、それ以外が0.0の高さを持つマップを近似した結果は以下のようになりました。
近傍5x5の25x10の行列を作って、2次曲面の係数10個を計算しました。
パラメータのガウスフィルタ0回適用です。

svd_002.png svd_014.png
     原画像         近似画像

PSNRが66ぐらいで、エッジ周辺にノイズがみられます。

2. 高さマップを近似した結果は以下のようになりました。
近傍5x5、ガウスフィルタ2回適用です。

svd_004.png svd_013.png
     原画像         近似画像

PSNRが101ぐらいで、かなり再現できています。エッジが少ないからでしょうか。

3. おまけとして、Lenna画像をグレースケール化したものを試してみました。
こちらも近傍5x5で、ガウスフィルタ1回適用です。

svd_006.png svd_012.png
     原画像         近似画像

PSNRが84ぐらいで、まずまずですが、やはりエッジ付近にノイズがみられます。
ベジエ曲面やスプライン曲面ならうまくいくのでしょうか。

次はKinectで取得した点群データをこの2次曲面近似してみようかなと思います。
Kinect持ってないですが。
どこかにKinectで取得した点群データのデータベースとかないんでしょうかね。
posted by シンドラー at 22:44 | Comment(0) | TrackBack(0) | Image Processing | このブログの読者になる | 更新情報をチェックする
この記事へのコメント
コメントを書く
お名前:

メールアドレス:

ホームページアドレス:

コメント:

認証コード: [必須入力]


※画像の中の文字を半角で入力してください。
※ブログオーナーが承認したコメントのみ表示されます。
この記事へのトラックバックURL
http://blog.seesaa.jp/tb/383104769
※ブログオーナーが承認したトラックバックのみ表示されます。

この記事へのトラックバック