2014-11-23

Embedded Zerotree Wavelet Encoding

Embedded Zerotree Wavelet Encoding 情報圧縮画像圧縮圧縮 EZW

資料
EZWの名前の由来
- Wavelet transformsを使う(Wavelets のＷ)
- Progressive encodingと呼ばれる。それは、粗い圧縮から細かい圧縮へと連続的に変化させられることを指す。bitの列に変換し、そのbit列の桁数のどこで止めても良い。Embedded encodingとも言う(EmbeddedのＥ)
- Zerotreeは解像度の粗密にまたがった共通性を保管するデータ構造につけられた名前(ZerotreeのＺ)
EZWの特徴
- ２次元画像用に作られたが任意次元に使える
- 自然画像の特徴として、次の特徴がある。Wavelet transformsでは、ハイパス（詳細系）とローパス（近似系）に分けるという処理をするが、自然画像では、そのローパス部分に情報が多い
- Waveletsに付与される係数の大きい方が大事なので、係数の降順でエンコーディングする
- あるWaveletsフィルタに対して、かけるたびに係数を出して、大きい係数が得られれば適切に入れて、小さい係数は後回しにする
- 解像度の粗いところから細かいところに順次進んでいく
- Wavelet 係数が粗い解像度と細かい解像度とで似ていることを利用して、圧縮する。このスケール間の関係を納めるのがZerotree
- ２次元（画像）の場合は、ある一つの要素は、その下部に2x2=4要素を従えている(d次元なら $2^d$ を従えたzerotreeになるはず)
- 圧縮情報の保管にあたっては、配置ルールをあらかじめ定めることで、bit情報が納められた位置が何の情報なのかを別途保管する必要を排している
ZeroTreeの仕組み
- ある閾値を定めて、自身と自身の子孫とがすべての値（の絶対値）が閾値より小さいときは、自身を"zerotreeのルート"と記録して、その子孫に関する情報は記録しないことにする
- もし、自身の値は閾値を超えないが、子孫の中には超えるものがあるようなら、「自身は大したことないけど、子孫には注意」というフラグを立てる
- 自身がしょぼけりゃ、子孫もしょぼいことが多いので、これをやると、圧縮できる、というやりくち
閾値をだんだん下げていくときのしくみ
- 閾値を下げるときに、ある閾値で記録するに足るだけの値を持っていたときは、そこで記録をする。その代り、そのセルの値を減らしてやって、閾値を下げたときの処理において、しょぼいとみなされる確率を上げる。これによって、何度も閾値を変えて記録しながら、記録しすぎを省く
エンコードした後、デコードするのに必要な情報
- (データ次元とエンコード前の全平均値があるとよい)
- Wavelet transformのレベル数(同じWaveletsを使うのであれば、これだけでＯＫ)
アルゴリズム
- 初期閾値
  - 初期位置値の決め方
    - $t_0 = 2^{log_2(MAX(|\gamma(x,y)|))}$ 。 $\gamma(x,y)$ は係数
- メインのループ処理
  - 初期閾値から始めて、２つのパスで処理していく
  - 分けたら、閾値を半分にする
  - 最小閾値までやったら終了
MATLABでの実装例がこちらにある
のアレイを倍々のサイズ拡大していく
- これっていうのは、1からまでの自然数を一辺の長さのd次元アレイとして配置する問題
  - 第一ステップとして、 $2^d$ の配置をd桁の二進法で決め、それを、 $(2^d)-$ 進法のk桁で表すようなもの
  - アレイの座標をハンドリングする

# mortonorder
my.ezw.mortonorder <- function(N,d=2){
	if((log(N,2^d) %% 1) !=0){
		print("Size of X is not appropriate")
		return(0)
	}
	k <- log(N,2^d)
	ret <- matrix(0,N,d)
	init.array <- array(1:(2^d),rep(2,d))
	init.array.add <- which(init.array > 0, arr.ind=TRUE)
	ret[1:(2^d),] <- init.array.add
	cnt <- 2^d+1
	cnt2 <- 2^d+1
	tmpN <- 2^d
	for(i in 2:k){
		tmpK <- 2^(i-1)
		for(j in 2:length(init.array.add[,1])){
			tmp <- init.array.add[j,]
			ret[cnt2:(cnt2+tmpN-1),] <- t(t(ret[1:(cnt-1),]) + (tmp-1) * tmpK)
			cnt2 <- cnt2+tmpN
		}
		cnt <- cnt2
		tmpN <- cnt -1
	}
	ret
}
my.ezw.mortonorder(64^2,2)
my.ezw.mortonorder((2^3)^3,3)

- - 次いで、親子関係を作る

my.ezw.tree <- function(N,d){
	k <- log(N,2^d)
	el <- cbind(rep(1,(2^d-1)),2:(2^d))
	scan <- my.ezw.mortonorder(N,d)
	scan.v <- apply(t(scan-1) * (2^k)^(0:(d-1)),2,sum) + 1
	for(i in 2:length(scan[,1])){
		tmp <- scan[i,]
		tmp.list <- list()
		for(j in 1:length(tmp)){
			tmp.list[[j]] <- c(-1,0) + tmp[j] *2
		}
		tmp.expgrid <- as.matrix(expand.grid(tmp.list))
		tmp.v <- apply(t(tmp.expgrid-1)*(2^k)^(0:(d-1)),2,sum) + 1
		if(max(tmp.expgrid)<=(2^k)){
			tmp.s <- which(scan.v %in% tmp.v)
			el <- rbind(el,cbind(rep(i,length(tmp.s)),tmp.s))
		}
	}
	el
}
d <- 2
k<-4
N <- (2^d)^k
scan <- my.ezw.mortonorder(N,d)
plot(scan,col=1:N)
plot(scan,col=1:N,pch=20)

el <- my.ezw.tree(N,d)
library(igraph)
g <- graph.edgelist(el)
plot(g)
degree(g,mode="out")
sh.p <- shortest.paths(g,mode="out")
sh.p[which(sh.p==Inf)] <- 0
col <- rep(1,N)
col[which(sh.p[8,] > 0)] <- 2
plot(scan,pch=20,col=col)