hirax.net::Keywords::「鳥瞰図」のブログ

■ボケたエアーブラシで細かな字がかけるか?　

画像復元を勉強してみたい その2

「宇宙人はどこにいる? - 画像復元を勉強してみたい その1-」ではボケた画像からオリジナルのシャープな画像を復元してみた。前回の話を例えて言うと、

1. 太郎君が細かい字をエアーブラシで書いた。
2. ボケボケのエアーブラシを使ったから、ボケボケの画になった。
3. そのボケボケの画から、太郎君が何を画こうとしたか、考える。

ということであった。
　今回、やってみたいのは以下のようなことである。

1. 太郎君は太いエアーブラシで字を書きたい。
2. しかも細かな字を書きたい。
3. そんなことができるか?

　直感的には、ボケボケのエアーブラシで細かい字など書けないように思う。その直感が正しいか調べてみたい。考え方は前回と同じく、

　出力画像から、ボケ分布でデコンボリューション処理により、オリジナルの画像を計算する。

というやり方である。前回と違うのは出力画像がシャープな画像(先の例で言うと、細かな字)である、という所である。道具は今回もMathematicaを使う。

　出力したい画像ファイルを読み込む。
<< Utilities`BinaryFiles`
StreamFile = OpenReadBinary["E:\jun\private\dekirukana\ufo\ufo.raw"]
ImageData = Table[ ReadBinary[ StreamFile , Byte] ,{x,64},{y,64}];
ListDensityPlot[ImageData,Mesh->False,PlotRange->{0,255}]

これが得たい出力画像である

画像:1

　この細かな字を太いボケボケなエアーブラシで字を書けるか考える。

　まずは、エアーブラシのボケボケ度をつくる。

(*正規分布=ガウス分布によるぼけパラメータを作成する*)
δ=10;
μ=32;

ListPlot3D[NormalBoke,ColorFunction ->Hue,Mesh->False,PlotRange->All]

ガウス分布のボケ(例で言うと、中央にエアーブラシを吹いた場合に相当する)

画像:2

　デコンボリューション用にガウス分布の場所をずらす。
NormalBoke = RotateRight[NormalBoke,32];
NormalBoke = Transpose[ RotateRight[Transpose[NormalBoke],32] ]; (*上へShift*)
ListPlot3D[NormalBoke,ColorFunction ->Hue,Mesh->False,PlotRange->All]

場所をずらしたボケ

画像:3

　出力画像をエアーブラシのボケボケ度でデコンボリューションする。そうすれば、太郎君がどのように画を画けば良いかがわかる。はたして答えはでるのだろうか?

　計算してみると答えが出てしまう。
SharpImage = Re[InverseFourier[ Fourier[ImageData] / Fourier[NormalBoke]] ];
ListDensityPlot[SharpImage/4,Mesh->False,PlotRange->All]

これが計算されたオリジナル画像

画像:4

出力画像を確認したもの

画像:5

　これには実はタネがある。画像:4を鳥瞰図でみると判るが、画像4は正負の値が高周波で並んでいる。
ListPlot3D[SharpImage/4,ColorFunction ->Hue,Mesh->False,PlotRange->All]

計算されたオリジナル画像(画像:4の鳥瞰図表示)

画像:6

　前回のような光学系の例でも、これが何に対応しているかはすぐわかるが、一番分かりやすいのは電荷と電位の例だと思う。
　電荷が周囲につくる電位分布はボケボケの分布である。ところが、金属などを適当に配置して、その金属に電位を印加してやると、鋭い電位分布をつくることができる。つまり、ボケボケの分布から鋭い電位分布を作成してやることができる。こちらなら直感的にもすぐ納得できるだろう。その際には、金属表面に電荷が鋭く集中するのも、よく知っている話だ。

　実感用に電場計算を行った例を以下に示しておく。使った道具はCUPSの電場計算プログラムである。CUPSは教育用のプログラム集である。
　一応、2次元膜の例で、金属を配置し、適当に電位を印加し、電位・電荷量計算を行ってみる。

電位分布を表示したもの

もちろん、金属内部では均一な電位である。それを条件に解いているのだから当たり前だが。
　その時の電荷分布を下に示す。金属表面に鋭い電荷分布が生じているのがわかるだろう。
　ここでは大雑把な金属の配置にしてしまったが、格子状の金属配置にして、互い違いに違う極性の電位を印加すれば(細かい字に相当する)、正負の極性の電荷分布が鋭く現れるのは当たり前の話だ。

電荷量分布を表示したもの

　電位、電場、電荷量を一緒に示しておく。

電位(左上)、電場の大きさ(左下)、電荷量(右上)、電場(右下)

　今回の話は、単なる計算上の話である。それに、何かどこかで仮定を間違っているような気もするんだよなぁ。信用度アルファ版だからまぁいいか...

■画像ノイズ解析について考える　

考える理由

　その疑問とは次のようなことである。

• 通常の2D-FTでは、入力データ全領域での周波数解析を行う。従って、単発のパルスのようなノイズはバックグラウンドに埋もれてしまい、結果にはなかなか出てこない。
• 同じ理由で、2D-FTでは位置と周波数解析を同時に行うことができない。(もちろん、短時間フーリエ関数を使えば、そのような測定は行うことができる。)
• また、ホワイトノイズのようなフラットな周波数特性を持つノイズもバックグラウンドを押し上げるだけの効果しか持たないため、解析をしづらい。

2D-FTと2D-Waveletの例

2D-FTの例(左から原画像、2D-FT結果、2D-FT結果の鳥瞰図)

　それでは、同じ画像に2D-Waveletをかけてみる。2D-Waveletの結果は位置と周波数強度分布情報(ホントは違うのだが)が両方出てくる。位置情報が2次元で周波数強度分布情報が1次元であるから、合わせて3次元である。そのため、表示に一工夫いる。
　第一段階として高周波成分から調べてみる。すぐにこの結果の意味がわかるだろうか？

2D-Wavelet例(左が原画像、右が一段階Waveletをかけた結果)

　なお、フーリエ変換では基底関数としてSinが用いられるが、Wavelet変換では基底関数としていろいろな関数を使うことができる。今回はDaubechiesの4次のものを用いている。下がその形である。

Daubechiesの4次のフィルター

ドットのノイズを解析してみる

ドット画像(左が原画像、右がノイズを加えた画像)

２次元離散フーリエ変換を行った結果

２次元離散フーリエ変換の結果を鳥瞰図で示したもの

2D-Waveletによる解析結果(左がノイズ無し、右がノイズ有り)

■分数階微分に基づく画像特性を考えてみたい　

同じ年齢でも大違い

前回、分数階微分の謎 - 線形代数、分数階微分、シュレディンガー方程式の三題話- で分数階微分について調べた。例えば、0.7階微分といった、整数階でない微分である。今回はそれを使った応用を考えてみたい。

人間の視覚というものは明るいものは強く感じることができる。これは当たり前である。そして、それだけでなく、強さが変化している所にも(興味を)強く感じ取るようになっている。岡本安春氏の「Delphiでエンジョイプログラミング」によれば、そのような考えはLaming(1986)がdifferential coupling(差動結合)として発表しているらしい。

ということは、人間が画像を感じる特性というものは、画像強度と画像強度変化(画像強度の一階微分)の中間的なものであると言うことができるかもしれない。とすれば、分数階微分を導入すれば面白い表現ができるかもしれない。
今回は、そういう考えのもとに分数階微分を用いて人間の画像特性について考えてみたい。

まずは、元画像を示す。元画像はガウス分布に基づいて作成されたものである。

元画像とその鳥瞰図

まずは、左の元画像を見て欲しい。どこに強い感じを受けるだろうか？白い部分はもちろんであるが、白と黒の境界部にも強い感じを受けるだろう。ギザギザになっているのはデータが少ないからなので、無視して欲しい。というわけで、人間の視覚画像特性は

• 画像強度
• 画像強度変化(画像強度の一階微分)

元画像から元画像の一階微分までの分数階微分画像

元画像
1/2階微分画像
15/20階微分画像
1階微分画像

白地に黒画像バージョンも示しておく。紙の上の画像に慣れた人にはこちらの方が良いだろう。

元画像から元画像の一階微分までの分数階微分画像(白地背景)

元画像

1/2階微分画像

15/20階微分画像

1階微分画像

なお、今回の画像の作成は次のような手順で行っている。

1. 1次元のガウス分布を作成する。
2. 微分値が正であるような半分の領域を線対称に回転させ、2次元画像を作成する。

今回は

• 画像強度
• 画像強度変化(画像強度の一階微分)

• 電位
• 電界(電位の微分、といっても本来は電位が電界の積分か)

• 人口密度
• 人口密度変化(人口密度の微分)

さて、分数階微分を調べる中で、バナッハ空間についても調べた。調べ始めた時には、聞き覚えもなかったが、調べてみるとヒルベルト空間の導入で登場していた。きれいさっぱり忘れていたようである。
京大数学教室 徳永健一氏のWEB (http://www.kusm.kyoto-u.ac.jp/~kenichi/)
から辿れる「「年齢の本」数学者版」によれば
バナッハがバナッハ空間を提唱したのは30歳の時であるらしい。(http://www.kusm.kyoto-u.ac.jp/~kenichi/age/30.html)
うーん...

■電界計算をしてみたい[有限要素法編その1]　

有限と微小のパン

　今回のサブタイトルは一目瞭然であるが、森博嗣のミステリのタイトルそのままである。

森博嗣 著 「有限と微小のパン」

　何故、「電界計算をしてみたい-有限要素法編その1-」が「有限と微小のパン」に繋がるのか。もちろん、"有限要素法"と"有限と微小のパン"の「有限」をかけた駄洒落ではない。有限要素法を考えるとき、私は森博嗣に足を向けては寝ることができない。それが、なぜかは下の本を見ればわかる。

森博嗣 著 「C言語による有限要素法入門」

　これは、学生時代に有限要素法を勉強するために使った本である。「森 博嗣　著」と書いてあるのがわかるだろうか。いや、まさかこの本の作者がミステリを量産するとは想像もしなかった。ビックリである。講談社ノベルズと森北出版の両方から本を出している人は他にはいそうにない。

　本題と関係のない話はここまでにしておく。今回はMathematicaで有限要素法を用いて静電界計算を行いたい。とりあえず、ソルバーとプリ・プロセッサまでつくる。その応用は続きの回で行いたい。Mathematicaで有限要素法を勉強するには、森北出版の依田　潔　著「Mathematicaによる電磁界シミュレーション入門」を参考にした。任意の電荷配置のPoisson方程式を解くようにしてある。

　今回使用したMathematicaのNotebookをHTMLで出力したものをここに示す。Notebook中でエラーが表示されているところは初期設定の変数をきちんと設定してやれば、エラーは出ないはずである。

次回に詳しく計算モデルの説明を行うので、今回は計算モデルの詳細については記述しない。Notebook内に、モデルの詳細は記述してある。

　このNotebookを使った計算、出力例を以下に示す。

電界計算の例

平行平板電極の間に誘電体層があるモデル	平板電極と三角柱電極の間に誘電体層があるモデル	平板電極と円柱電極の間に誘電体層があるモデル
分割要素	分割要素	分割要素
電位表示(色がきちんとしたhueでないことに注意)	電位表示(色がきちんとしたhueでないことに注意)	電位表示(色がきちんとしたhueでないことに注意)
半分の領域の電位を鳥瞰図にしたもの	半分の領域の電位を鳥瞰図にしたもの	半分の領域の電位を鳥瞰図にしたもの

　Mathematica3.0のHTML出力は大変便利だが、漢字が化けるのが困りものだ。しかも、ちょっと似た漢字に化けてしまうからわかりにくい。今回のNotebook中で化けた漢字を以下に示す。

• 油界　<- 電界
• 堰素　<- 要素
• 誘油　<- 誘電
• 姦み込む　<- 組み込む
• 表傭　<- 表面
• 壓さ　<- 高さ
• 堆心　<- 重心
• 肖似　<- 近似
• 内占　<- 内部
• 傭積　<- 面積
• 回寂　<- 回転
• 進当　<- 適当
• 懷瞰　<- 鳥瞰

中国語みたいな化け方である。しかも、意味としても何か変な化け方である。いつか、この対処方法と理由を考えてみたい。それにしても、週末の遊び道具としてはMathematicaは素晴らしいと思う。

■Photoshop鳥瞰図＆表計算フィルタ　

　ひとまず完成。フィルタダイアログからの直接プリントは今のところできない。（リンク）（リンク）