Coiflet小波

Coiflet小波是英格丽·多贝西应罗纳德·科夫曼的要求所设计的一种离散小波。Coiflet小波的调整函式（scaling function）及小波函式（wavelet function）能同时拥有高消失动量，且其波形接近对称，常被用于数字信号处理。^[1]

Coiflet with two vanishing moments

性质

• 定义

一个正交小波（英语：Orthogonal_wavelet）（orthogonal wavelet）系统若符合以下条件，则称为Generalized Coiflet Wavelet(GOC)。

${\displaystyle {\begin{array}{lcl}\\{\mathcal {M_{\tilde {\psi }}}}(0,l]=0&{\mbox{for }}l{\mbox{ =0,1,...,L-1}}\\{\mathcal {M_{\tilde {\phi }}}}(t_{0},l]=\delta [l]&{\mbox{for }}l{\mbox{ =0,1,...,L-1}}\\\end{array}}}$

其中${\displaystyle {\psi }}$为小波函式，${\displaystyle {\phi }}$为调整函式

• 滤波器长度与消失动量

若${\displaystyle N}$为滤波器长度，而${\displaystyle L}$为系统消失动量，则最小的${\displaystyle N}$为：

${\displaystyle N={\begin{cases}3L-1,&{\mbox{if }}L{\mbox{ is odd}}\\3L,&{\mbox{if }}L{\mbox{ is even}}\\\end{cases}}}$

• 近线性相位（Near-linear phase）滤波器

当${\displaystyle |w|}$够小时，GOC的低通滤波器的频率响应拥有渐进型态（asymptotic form）：

${\displaystyle \angle {H(w)}=-t_{0}w+C\cdot {w^{L'}}+O(w^{L'+2})}$

其中${\displaystyle L'=2\left\lfloor L/2\right\rfloor +1}$，${\displaystyle C}$则为常数项。因此${\displaystyle h}$具有渐近线性相位（asymptotic linear phase）的特性：

${\displaystyle \lim _{L\to \infty }\angle {H(w)}=-t_{0}w}$

应用

Coiflet小波除了被使用在影像压缩外，目前也有被使用在电力系统讯号监测上^[2]

系数

下表为coiflet小波调整函式C6~30的系数，小波函式的系数可以借由每两个系数变换一次符号来推导(例如C6 wavelet = {−0.022140543057, 0.102859456942, 0.544281086116, −1.205718913884, 0.477859456942, 0.102859456942})

Coiflets系数

k	C6	C12	C18	C24	C30
-10					-0.0002999290456692
-9					0.0005071055047161
-8				0.0012619224228619	0.0030805734519904
-7				-0.0023044502875399	-0.0058821563280714
-6			-0.0053648373418441	-0.0103890503269406	-0.0143282246988201
-5			0.0110062534156628	0.0227249229665297	0.0331043666129858
-4		0.0231751934774337	0.0331671209583407	0.0377344771391261	0.0398380343959686
-3		-0.0586402759669371	-0.0930155289574539	-0.1149284838038540	-0.1299967565094460
-2	-0.1028594569415370	-0.0952791806220162	-0.0864415271204239	-0.0793053059248983	-0.0736051069489375
-1	0.4778594569415370	0.5460420930695330	0.5730066705472950	0.5873348100322010	0.5961918029174380
0	1.2057189138830700	1.1493647877137300	1.1225705137406600	1.1062529100791000	1.0950165427080700
1	0.5442810861169260	0.5897343873912380	0.6059671435456480	0.6143146193357710	0.6194005181568410
2	-0.1028594569415370	-0.1081712141834230	-0.1015402815097780	-0.0942254750477914	-0.0877346296564723
3	-0.0221405430584631	-0.0840529609215432	-0.1163925015231710	-0.1360762293560410	-0.1492888402656790
4		0.0334888203265590	0.0488681886423339	0.0556272739169390	0.0583893855505615
5		0.0079357672259240	0.0224584819240757	0.0354716628454062	0.0462091445541337
6		-0.0025784067122813	-0.0127392020220977	-0.0215126323101745	-0.0279425853727641
7		-0.0010190107982153	-0.0036409178311325	-0.0080020216899011	-0.0129534995030117
8			0.0015804102019152	0.0053053298270610	0.0095622335982613
9			0.0006593303475864	0.0017911878553906	0.0034387669687710
10			-0.0001003855491065	-0.0008330003901883	-0.0023498958688271
11			-0.0000489314685106	-0.0003676592334273	-0.0009016444801393
12				0.0000881604532320	0.0004268915950172
13				0.0000441656938246	0.0001984938227975
14				-0.0000046098383254	-0.0000582936877724
15				-0.0000025243583600	-0.0000300806359640
16					0.0000052336193200
17					0.0000029150058427
18					-0.0000002296399300
19					-0.0000001358212135

滤波器

这边列出6, 12 ,18点的filter消失点(vanish moment)也就是k/6。

这边的g是低频的滤波器，h是高频的滤波器。

可以看出${\displaystyle h[n]=(-1)^{n}g[-n]}$，这边采用python的语法，h会是g的反序，且n为奇数时要乘上-1

Daubechies filter ${\displaystyle g_{k},h_{k}}$	0	1	2	3	4	5	6	7	8	9	10	11	12	13	14	15	16	17
${\displaystyle h_{6}}$	-0.0727	-0.3378	0.8525	-0.38486	-0.0727	0.0156
${\displaystyle g_{6}}$	-0.0156	-0.0727	0.38486	0.8525	0.3378	-0.0727
${\displaystyle h_{12}}$	0.0010	-0.0026	0.0079	-0.0335	0.0841	-0.1082	0.5897	-1.1494	0.5460	-0.0953	-0.0586	-0.0232
${\displaystyle g_{12}}$	0.0232	0.0586	0.0953	0.5460	1.1494	0.5897	0.1082	0.0841	0.0335	0.0079	0.0026	0.0010
${\displaystyle h_{18}}$	-0.003793	-0.007782	0.023452	0.065771	-0.06112	-0.40517	0.79377	-0.42848	-0.071799	0.082301	0.034555	-0.015880	-0.0090079	0.002574	-0.001117	0.000466	-0.00007098	-0.0000345
${\displaystyle g_{18}}$	-0.0000345	-0.00007098	0.000466	0.001117	-0.002574	-0.0090079	0.015880	0.034555	-0.082301	-0.071799	0.42848	0.79377	0.40517	-0.06112	-0.065771	0.023452	0.007782	-0.003793

Matlab程式

F = coifwavf(W)会回传N=str2num(W)的coiflet小波的调整函式，其中N只能为1~5的整数。^[3]

Python程式

这边根据上面所提供的高频和低频的滤波器去实作出离散小波转换并且应用到影像压缩上面，

透过小波转换良好的性质，尝试实作出连续的影像压缩。

这边单纯是在示范小波转换进行多次压缩和反压缩的结果，

如果不想转换这么多次修改Jmax即可。

import numpy as np
def subsampling(x, d):
    if d == 1:
        y = x[::2, :]
    elif d == 2:
        y = x[:, ::2]
    return y

def upsampling(x, d):
    s = x.shape
    if d == 1:
        y = np.zeros((p * s[0], s[1]))
        y[::2, :] = x
    elif d == 2:
        y = np.zeros((s[0], p * s[1]))
        y[:, ::2] = x
    return y
    
def cconv(x, h, d):
    if d == 2:
        return np.transpose(cconv(np.transpose(x), h, 1))
    y = np.zeros(x.shape)
    p = len(h)
    pc = int(round( float((p - 1) / 2 )))
    for i in range(0, p):
        y = y + h[i] * np.roll(x, i - pc, axis=0)
    return y
    
def DWT_compression(image, h, g): # discrete wavelet transformation
    Jmin = 0
    Jmax = log2(n)-1
    fW = image.copy()
    j = int(np.log2(image.shape[0])-1)
    A = fW[:2**(j+1):,:2**(j+1):]
    for j in arange(Jmax,Jmin-1,-1):
        A = fW[:2**(j+1):,:2**(j+1):]
        for d in arange(1,3):
            Coarse = subsampling(cconv(A,h,d),d)
            Detail = subsampling(cconv(A,g,d),d)
            A = concatenate( (Coarse, Detail), axis=d-1 )
        fW[:2**(j+1):,:2**(j+1):] = A
        j1 = Jmax-j
    return fW
    
def iDWT_decompression(image, fW, h, g): #image is original, fW is after DWT of that.
    Jmin = 0
    Jmax = int(np.log2(image.shape[0])-1)
    f1 = fW.copy()
    j = int(np.log2(image.shape[0])-1)
    A = f1[:2**(j+1):,:2**(j+1):]
    h1 = h[::-1]
    g1 = g[::-1]
    for j in arange(Jmin,Jmax+1): 
        A = f1[:2**(j+1):,:2**(j+1):]
        for d in arange(1,3):
            if d==1:
                Coarse = A[:2**j:,:]
                Detail = A[2**j:2**(j+1):,:]
            else:
                Coarse = A[:,:2**j:]
                Detail = A[:,2**j:2**(j+1):]
            Coarse = cconv(upsampling(Coarse,d),h1,d)
            Detail = cconv(upsampling(Detail,d),g1,d)
            A = Coarse + Detail
            j1 = Jmax-j
            if j1>0 and j1<5:
                imageplot(A, 'j=' + str(int(j)), [2,2,j1])
    f1[:2**(j+1):,:2**(j+1):] = A
    return f1

双正交Coiflet小波

• 定义

一个order为${\displaystyle (L,{\tilde {L}})}$的双正交Coiflet小波需符合以下条件：

${\displaystyle {\begin{array}{lcl}\\{\mathcal {M_{\tilde {\phi }}}}(0,l]=\delta [l]&{\mbox{for }}l{\mbox{ =0,1,...,L-1}}\\{\mathcal {M_{\tilde {\psi }}}}(0,l]=0&{\mbox{for }}l{\mbox{ =0,1,...,L-1}}\\{\mathcal {M_{\phi }}}(0,l]=\delta [l]&{\mbox{for }}l{\mbox{ =0,1,...,L-1}}\\{\mathcal {M_{\psi }}}(0,l]=0&{\mbox{for }}l{\mbox{ =0,1,...,}}{\tilde {L}}{\mbox{-1}}\\\end{array}}}$

• 性质

当${\displaystyle L}$为偶数时，${\displaystyle h}$会对称于原点：${\displaystyle h[n]=h[-n]}$。这让双正交coiflet在图片压缩方面能有较好的峰值信噪比（PSNR）。