CDF 9/7 离散小波变换（卷积）

Question

我正在尝试写一个简单的程序，目的是对一个一维列表进行一次离散小波变换，使用的是CDF 9/7小波，然后再把它重建回来。我只是用卷积/滤波器组的方法来理解这个过程。简单来说，就是用一个滤波器对列表进行卷积，得到尺度系数；再用另一个滤波器对列表进行卷积，得到小波系数，但只从每隔一个元素开始。接着进行上采样（也就是在元素之间加零），对小波和尺度系数应用滤波器，把它们加在一起，最终得到原来的列表。

我能让这个方法在Haar小波滤波器上工作，但当我尝试使用CDF 9/7滤波器时，结果却无法还原成原始输入。不过，结果列表和原始列表的和是相同的。

我相信这一定是卷积中的一个很简单的错误，但我就是搞不明白。我尝试了很多卷积的不同组合，比如把滤波器中心放在索引“i”上，而不是从左边开始，但似乎都没有效果……这可能是那种等我搞明白后会让我拍头的错误。

这是我的代码：

import random
import math
length = 128
array = list()
row = list()
scaleCoefficients = list()
waveletCoefficients = list()
reconstruction = list()

def upsample(lst, index):
    if (index % 2 == 0):
        return 0.0
    else:
        return lst[index/2]

for i in range(length):
    array.append(random.random())

## CDF 9/7 Wavelet (doesn't work?)
DWTAnalysisLowpass = [.026749, -.016864, -.078223, .266864, .602949, .266864, -.078223, -.016864, .026749]
for i in range(len(DWTAnalysisLowpass)):
    DWTAnalysisLowpass[i] = math.sqrt(2.0) * DWTAnalysisLowpass[i]
DWTAnalysisHighpass = [0.0, .091272, -.057544, -0.591272, 1.115087, -.591272, -.057544, .091272, 0.0]
for i in range(len(DWTAnalysisHighpass)):
    DWTAnalysisHighpass[i] = 1.0/math.sqrt(2.0) * DWTAnalysisHighpass[i]

DWTSynthesisLowpass = [0.0, -.091272, -.057544, 0.591272, 1.115087, .591272, -.057544, -.091272, 0.0]
for i in range(len(DWTSynthesisLowpass)):
    DWTSynthesisLowpass[i] = 1.0/math.sqrt(2.0) * DWTSynthesisLowpass[i]
DWTSynthesisHighpass = [.026749, .016864, -.078223, -.266864, .602949, -.266864, -.078223, .016864, .026749]
for i in range(len(DWTSynthesisHighpass)):
    DWTSynthesisHighpass[i] = math.sqrt(2.0) * DWTSynthesisHighpass[i]

## Haar Wavelet (Works)
## c = 1.0/math.sqrt(2)
## DWTAnalysisLowpass = [c,c]
## DWTAnalysisHighpass = [c, -c]
## DWTSynthesisLowpass = [c, c]
## DWTSynthesisHighpass = [-c, c]


## Do the forward transform - we only need to do it on half the elements
for i in range(0,length,2):
    newVal = 0.0
    ## Convolve the next j elements
    for j in range(len(DWTAnalysisLowpass)):
        index = i + j
        if(index >= length):
            index = index - length

        newVal = newVal + array[index]*DWTAnalysisLowpass[j]

    scaleCoefficients.append(newVal)

    newVal = 0.0
    for j in range(len(DWTAnalysisHighpass)):
        index = i + j
        if(index >= length):
            index = index - length

        newVal = newVal + array[index]*DWTAnalysisHighpass[j]

    waveletCoefficients.append(newVal)

## Do the inverse transform
for i in range(length):
    newVal = 0.0
    for j in range(len(DWTSynthesisHighpass)):
        index = i + j
        if(index >= length):
            index = index - length

        newVal = newVal + upsample(waveletCoefficients, index)*DWTSynthesisHighpass[j]

    for j in range(len(DWTSynthesisLowpass)):
        index = i + j
        if(index >= length):
            index = index - length

        newVal = newVal + upsample(scaleCoefficients, index)*DWTSynthesisLowpass[j]

    reconstruction.append(newVal)

print sum(reconstruction)
print sum(array)
print reconstruction
print array

顺便提一下，我从这里的附录中获取了滤波器的值：http://www1.cs.columbia.edu/~rso2102/AWR/Files/Overbeck2009AWR.pdf，我也在很多matlab示例代码中见过这些值。

信号处理卷积滤波器离散小波变换 cdf 9/7 上采样尺度系数小波系数

CDF 9/7 离散小波变换（卷积）

1 个回答

撰写回答