1. No category

画像処理論
第9回
佐藤 洋一
生産技術研究所
情報理工学系研究科電子情報学専攻／情報学環
http://www.hci.iis.u-tokyo.ac.jp/~ysato/
前回講義内容の復習

2値画像の処理


モルフォロジーによる各種操作
クラスタリング

階層的クラスタリング

分割最適化クラスタリング

スペクトラルクラスタリング
適用例
input
thinning
differencing, thinning end-point thinning,
hit-or-miss
end point thinning
dilation
dilation
Dilationオペレータ（膨張）
 dilation

of A by B
B：構造要素
A ⊕ B = {x ( Bˆ )  A ≠ φ}
x
A
ˆ =B
B
A⊕ B
Erosionオペレータ（収縮）
 erosion
of A by B
A − B = {x ( B) x ⊆ A}
Openingオペレータ
 opening

of A by B
細い所を削除する性質
A  B = ( A − B) ⊕ B
Closing
 closing

of A by B
細い溝を埋める性質
A • B = ( A ⊕ B) − B
Hit-or-Miss Transform

AとACのErosion

Erosionの結果を統合→Xの検出
A ∗ B = ( A − X )  [ AC − (W − X )]
その他の操作

Region filling
X k = ( X k −1 ⊕ B)  AC

Thinning
A ⊗ {B} = (( (( A ⊗ B1 ) ⊗ B 2 ) ) ⊗ B n )
A ⊗ B = A − ( A ∗ B) 特徴空間におけるクラスタリング

分類対象の特徴量を表現する空間


特徴空間（R-G平面）
色で分類⇒(r,g,b)
色と場所で分類⇒（x,y,r,g,b)
特徴空間への射影

クラスタリングとは、特徴空間内の分布から群（クラスタ）を
見つける処理
階層的クラスタリング

最初にN個のクラスタが与えられたとし、クラスタの数がCに
減少するまで最も類似した２つを融合

Cを事前に決定

クラスタ間の類似度の選択
－重心間距離法
－群平均法
－最短距離法
－最長距離法
－ウォード法
分割最適化クラスタリング： K-平均法
スペクトラルクラスタリングの概要
Data
Similarities
Block-Detection
ノード⇒クラスタの結合度と要素間類似度


k個のノードのあるクラスタへの結合度をk次元ベクトルwで
表現
ノード⇒クラスタ結合度と要素間類似度の合計
T
w Aw
∑{[ノードiの結合度]×∑{[ノードiとノードjの類似度] ×[ノードjの結合度]}}
類似度行列A
重みベクトルw
ノード1
ノードk
結合度と類似度の最大化

ノード→クラスタの結合度、要素間の類似度の合計の最大化
T
w Aw の最大化

重みベクトルの正規化を条件として
wT Aw + λ ( wT w − 1) の最大化

wに関する微分＝0より固有値問題に帰着
Aw = λw
Spectral Space

Can put items into blocks by eigenvectors:
1
1
.2
0
.71
0
1
1
0
-.2
.69
-.14
.2
0
1
1
.14
.69
0
-.2
1
1
0
.71
e1
e2

Clusters clear regardless of row ordering:
1
.2
1
0
.71
0
.2
1
0
1
.14
.69
1
0
1
-.2
.69
-.14
0
1
-.2
1
0
.71
e1
e2
e1
e2
e1
e2
スケールの影響
データ分布
固有値列
小←大
スペクトラルクラスタリングの欠点

類似度行列の複数の固有値が等しい場合
固有値列
データ分布
４つの固有ベクトル
Normalized Cuts


Current criterion evaluates within cluster similarity, but not
across cluster difference
Instead, we’d like to maximize the within cluster similarity
compared to the across cluster difference

Write graph as V, one cluster as A and the other as B

Maximize
N assoc =
V
assoc( A, A) assoc(B, B )
+
assoc( A,V ) assoc( B,V )
assoc(A,A) : the sum of weights of all edges within A
assoc(A,V): the sum of weights of all edges that have one end in A
J. Shi and J. Malik, “Normalized Cuts and Image Segmentation, IEEE Trans. PAMI, 22(8), 2000.
Nassoc and Ncut

Consider Ncut instead of Nassoc
N assoc
assoc( A, A) assoc(B, B )
=
+
assoc( A,V ) assoc( B,V )
N cut
V
cut ( A, B )
cut ( A, B )
=
+
assoc( A,V ) assoc( B,V )
Cut(A,B) : the sum of weights of all edges that have one end in A and the other in B

The cost Ncut is related to Nassoc as Ncut=2- Nassoc .
Finding Minimum Normalized-Cut

It can be shown that min N cut
y T (D − W )y
= min y
y T Dy
such that y (i ) ∈ {1,−b}, 0 < b ≤ 1, and y T D1 = 0


If y(i)=1, then i-th node belongs to one cluster, otherwise the other
cluster.
D is N x N diagonal matrix where
D(i, i ) = ∑ W (i, j )
j

If y is allowed to take real values then the minimization can
be done by solving the generalized eigenvalue system
(D − W )y = λDy
Normalized Cut: Overview

Compute matrices W & D

Solve

for eigen vectors
Use the eigen vector with second smallest eigen value to
bipartition the graph. (The smallest one is always zero.)


(D − W )y = λDy
Threshold for y(i) can be chosen s.t.
y T (D − W )y
y T Dy
is minimized.
Recursively partition the segmented parts if necessary.
Normalized Cutによる分割の例
本日の講義の内容

画像復元

画像復元とは？

ノイズ低減のためのフィルタ

劣化関数の推定に基づいた画像復元
画像復元
画像復元(image restoration)とは？

劣化プロセスに関する知識をもとに、劣化画像から元画像を
推定する処理
劣化画像（ぶれとノイズ）

復元画像
c.f. 画像強調(image enhancement)ではad-hocな処理が多い
劣化と復元プロセスのモデル

劣化関数がlinear & position-invariantとすると
実空間領域 : g ( x, y ) = h( x, y ) * f ( x, y ) + η ( x, y )
周波数領域 : G (u , v) = H (u , v) F (u, v) + N (u , v)
導出の詳細は補足資料参照（GW5.5）
ノイズ低減のための処理

ノイズのみの劣化の場合
g ( x, y ) = f ( x, y ) + η ( x, y )


空間領域のフィルタ

線形フィルタ（移動平均フィルタ、加重平均フィルタ）

非線形フィルタ（メディアンフィルタ、バイラテラルフィルタ）
周波数領域のフィルタ G (u , v) R (u , v) ⇒ Fˆ (u , v)

Low-pass, high-pass filters

Bandreject filters

Notch filters

Optimal notch filters
ノイズで劣化した画像の例
Bandreject Filters

ある大きさの周波数成分のみを削除

理想的なbandreject filterの例

1

R (u , v) = 0

1

if D (u , v) < D0 −
W
2
W
W
≤ D (u , v) < D0 +
2
2
W
if D (u , v) > D0 +
2
if D0 −
D0
W
D(u , v)
Butterworth Bandreject Filters

滑らかに変化するbandreject filterとしての性質

オーダーnのButterworth bandreject filter
R(u , v) =
1
 D(u , v)W 
1+ 
2
2
D
u
v
D
(
,
)
−
0 

2n
n =1
Bandreject filterによる画像復元の例
Notch Filters

Bandreject filterと異なり、ある特定の周波数成分のみを
削除するフィルタ

理想的なNotch filterの例

(u,v)平面上で原点対称の2点のみを削除
0 if D1 (u , v) ≤ D0 or D2 (u , v) ≤ D0
R(u , v) = 
1 otherwise
D1 (u , v)
D2 (u , v)
Butterworth Notch Filters

滑らかに変化するNotch filter
R (u , v) =
1


D
1+ 

D
u
v
D
u
v
(
,
)
(
,
)
2
 1

2
0
n
n=2
Notch filterによる画像復元の例
より複雑なノイズの問題

Bandreject filterやNotch filterでは上手く処理できない場合


ノイズ成分を削除してしまうと元画像自体が大きく劣化する可能
性が大
ノイズ信号を推定し、重みをコントロールしながらノイズを除去
Optimal Notch Filtering

まず最初に、劣化画像からノイズ信号を抽出


劣化画像G(u,v)のスペクトラルを見ながらユーザがNotch filter
R(u,v)を作成
ノイズ信号を復元
spikeを抽出
η ( x, y ) = ℑ−1{R(u , v)G (u, v)}
Optimal Notch Filtering

次に、劣化画像からノイズ成分を除去

単純に除去する代わりに、重みw(x,y)を調整
fˆ ( x, y ) = g ( x, y ) − w( x, y )η ( x, y )

復元画像の近傍領域内の画素値の分散が最小となる重み

導出はGW5.4.4参照
g ( x, y )η ( x, y ) − g ( x, y )η ( x, y )
w( x, y ) =
η 2 ( x, y ) − η 2 ( x, y )
Optimal Notch Filterによる処理の例
劣化関数の推定に基づいた画像復元

これまではノイズのみを考慮してきたが、劣化関数を取り扱
うにはどうすればよいか？
実空間領域 : g ( x, y ) = h( x, y ) * f ( x, y ) + η ( x, y )
周波数領域 : G (u, v) = H (u, v) F (u, v) + N (u, v)


Estimation of H(u,v) by

Experimentation

Mathematical modeling
劣化関数の推定を伴う画像復元はblind deconvolutionと
呼ばれる
Estimation by Experimentation

劣化画像を撮影した撮影システムに関して、事前に劣化プロ
セスを観察


インパルス応答を直接計測
G (u , v)
H (u , v) =
A
劣化画像のみを与えられた場合には不可能
Estimation by Modeling


劣化プロセスに関する物理モデルにもとづき劣化関数を決
定する方法
例：air turbulence
H (u , v) = exp(−k (u 2 + v 2 ) 5 / 6 )
k: turbulenceの度合い
手振れによる劣化の例

画像の露出時間t (t=0～T)の間に、平行移動x0(t), y0(t)により
ぶれた画像
手振れによる劣化関数
x0 (t ) =
bt
at
, y0 (t ) = の場合
T
T
T
sin[π (ua + vb)]exp(− jπ (ua + vb))
H (u , v) =
π (ua + vb)
（導出はGW5.6.3参照）
Motion Deblurring using Hybrid Imaging

Hybrid Camera（高解像度静止画＋低解像度高フレーム
レート）による手ぶれ画像の補正
Ben-Erza&Nayar, CVPR2003
Overview of Approach
Motion Analysis
y
x
Low-Res. camera
Same time period PSF Estimation
Hi-Res. camera
Deconvolution
Designs for Hybrid Imaging
A rig of two cameras
Using a beam splitter
Using a special chip
Prototype: Rig of Two Cameras
Primary detector
(2048x1536)
Secondary detector
(360x240)
Example - Blurred Night Image
f = 884mm, Exp. Time 4 Sec (> -11 stops)
PSF Estimation from Motion
Low resolution sequence.
0.003
Y (Pixels)
30
10
0.001
10
X (Pixels)
f = 884mm, Exp. Time 4 Sec
60
Deblurred Night Image
f = 884mm, Exp. Time 4 Sec
Comparison
Tripod image (Ground Truth)
Blurred image
Deblurred image
劣化関数を用いた画像復元

推定された劣化関数を用いて画像を復元
G (u , v)
N (u , v)
Fˆ (u , v) =
= F (u , v) +
H (u , v)
H (u , v)


復元画像と元画像とは完全には一致しない
H(u,v)が0に近いと、解が不安定になるという問題

劣化画像に含まれるノイズ成分に影響を受けやすい
画像復元の失敗例
G (u , v)
Fˆ (u , v) =
H (u , v)
元画像
劣化画像
H (u , v) = exp(−k (u 2 + v 2 ) 5 / 6 )
安定化のための方策

画像復元の際に、高周波成分は取り扱わない
G (u , v)
Fˆ (u , v) =
H (u , v)
劣化画像
Wiener Filter

ノイズに関する情報をより積極的に利用したフィルタとして
Wiener Filterが知られている。

元画像と復元画像との最小2乗誤差基準から導出
*


,
)
(
H
u
v
Fˆ (u , v) = 
G (u , v)
2
2
2
 H (u , v) + N (u , v) / F (u , v) 

しかし、ノイズ信号と元画像のパワースペクトルが既知なこと
は稀であるため、定数Kで代用することも多い
 H * (u , v) 
Fˆ (u , v) = 
G (u , v)
2
 H (u , v) + K 
Wiener Filterによる画像復元の例
Air turbulenceによる劣化画像の場合
H (u , v) = exp(−k (u 2 + v 2 ) 5 / 6 )
Wiener Filterによる画像復元の例
Motion blurによる劣化画像の場合
H (u, v) =
劣化画像
inverse filtering
Wiener filtering
T
sin[π (ua + vb)]exp(− jπ (ua + vb))
π (ua + vb)
本日の講義内容のまとめ

画像復元

画像復元とは？

ノイズ低減のためのフィルタ



空間領域

周波数領域
劣化関数の推定に基づいた画像復元

Estimation by experimentation

Estimation by modeling

Inverse filtering

Wiener filtering
参考資料

Forsyth & Ponce, Chapter 14.5

Gonzalez & Woods, Chapter 5.1, 5.4～5.8