2013-12-09

Chapter 4 Conventional DWT step by step　ぱらぱらめくる『Conceptual Wavelets in Digital Signal Procesing』

フィルタを２のべき乗で延ばすのではなく、シグナルの方を縮める方法
２つの注意点
- 同じ値があるとaliasingと呼ばれる変な現象が現れる
- Downsamplingの仕方によって、同一のシグナルから異なる結果が生じる(これを"time-invariant"である、と言う)
  - フィルタをかけた後のベクトルから、半分にダウンサンプリングするときに、偶数番目を残すか奇数番目を残すかで結果が変わる、ということ
  - これを実例で見ておこう

out.S <- my.wkeep.v(conv(c(1,1),S),length(S))
round(out.S)
round(out.S[which(1:length(out.S)%%2==0)])

> out.S <- my.wkeep.v(conv(c(1,1),S),length(S))
> round(out.S)
[1]  80 160 160 160  80   0   0   0
> round(out.S[which(1:length(out.S)%%2==0)])
[1] 160 160   0   0

ダウンサンプリングとアップサンプリングの関数を作っておこう
- アップサンプリングでベクトルの長さが２倍＋１になっていることに注意

my.dyadDWN <- function(x,d=2){
	x[which(1:length(x)%%d==0)]
}
my.dyadUP <- function(x,d=2){
	tmp <- matrix(0,d-1,length(x))
	c(rbind(tmp,x),0)
}

x
my.dyadDWN(x)
my.dyadUP(x)

> x
[1] 80 80 80 80  0  0  0  0
> my.dyadDWN(x)
[1] 80 80  0  0
> my.dyadUP(x)
 [1]  0 80  0 80  0 80  0 80  0  0  0  0  0  0  0  0  0

(1,1)フィルタをかけて、ダウンサンプリングをして、アップサンプリングをしなおして、averagingをすると、定数倍の点を除いて元に戻る

out1 <- conv(c(1,1),S)
out2 <- my.dyadDWN(out1)
out3 <- my.dyadUP(out2)
out4 <- conv(c(1,1),out3)
out5 <- my.wkeep.v(out4,length(S))
S
round(out1)
round(out2)
round(out3)
round(out4)
round(out5)

> S
[1] 80 80 80 80  0  0  0  0
> round(out1)
[1]  80 160 160 160  80   0   0   0   0
> round(out2)
[1] 160 160   0   0
> round(out3)
[1]   0 160   0 160   0   0   0   0   0
> round(out4)
 [1]   0 160 160 160 160   0   0   0   0   0
> round(out5)
[1] 160 160 160 160   0   0   0   0

細かくする方向と粗くする方向の両方で同じようにダウンサンプリング、アップサンプリングをして、最後に合算する

S <- c(rep(80,4),rep(0,4))
out1 <- conv(c(1,1),S)
out2 <- my.dyadDWN(out1)
out3 <- my.dyadUP(out2)
out4 <- conv(c(1,1),out3)
out5 <- my.wkeep.v(out4,length(S))
S
round(out1)
round(out2)
round(out3)
round(out4)
round(out5)

OUT1 <- conv(c(-1,1),S)
OUT2 <- my.dyadDWN(OUT1)
OUT3 <- my.dyadUP(OUT2)
OUT4 <- conv(c(1,-1),OUT3)
OUT5 <- my.wkeep.v(OUT4,length(S))
S
round(OUT1)
round(OUT2)
round(OUT3)
round(OUT4)
round(OUT5)

round(out5 + OUT5)

> S
[1] 80 80 80 80  0  0  0  0
> round(OUT1)
[1] -80   0   0   0  80   0   0   0   0
> round(OUT2)
[1] 0 0 0 0
> round(OUT3)
[1] 0 0 0 0 0 0 0 0 0
> round(OUT4)
 [1] 0 0 0 0 0 0 0 0 0 0
> round(OUT5)
[1] 0 0 0 0 0 0 0 0
> 
> round(out5 + OUT5)
[1] 160 160 160 160   0   0   0   0

シグナルに適当な値を入れても復元できる

S <- rnorm(8)
out1 <- conv(c(1,1),S)
out2 <- my.dyadDWN(out1)
out3 <- my.dyadUP(out2)
out4 <- conv(c(1,1),out3)
out5 <- my.wkeep.v(out4,length(S))
S

OUT1 <- conv(c(-1,1),S)
OUT2 <- my.dyadDWN(OUT1)
OUT3 <- my.dyadUP(OUT2)
OUT4 <- conv(c(1,-1),OUT3)
OUT5 <- my.wkeep.v(OUT4,length(S))
round((out4 + OUT4)/2,5)
S

> round((out4 + OUT4)/2,5)
 [1]  0.00000  0.75621 -1.81127 -1.43747 -1.94896 -0.31869  1.04858  0.62294
 [9] -0.46164  0.00000
> S
[1]  0.7562057 -1.8112699 -1.4374677 -1.9489588 -0.3186942  1.0485760  0.6229445
[8] -0.4616435

この後、省略するけれど、基本的には、アップサンプリング・ダウンサンプリングとフィルタと畳み込みの組み合わせ