如何扩展pyWavelets以处理N维数据?
|
这可能是另一个论坛的问题,如果是这样,请告诉我。我注意到只有14个人遵循小波标签。
我在这里有一种优雅的方法,可以将pywt(pyWavelets包)中的小波分解扩展到多个维度。如果安装了pywt,这应该是开箱即用的。测试1显示了3D阵列的分解和重组。所有要做的事情就是增加维数,并且代码将可以对4、6甚至18个维的数据进行分解/重组。
我在这里替换了pywt.wavedec和pywt.waverec函数。另外,在fn_dec中,我展示了新的wavedec函数如何像旧函数一样工作。
但是有一个问题:它把小波系数表示为与数据形状相同的数组。结果,由于我对小波的了解有限,所以我只能将其用于Haar小波。例如DB4之类的其他代码,例如在此严格边界的边缘上渗出系数(当前以数组列表[CA,CD1 ... CDN]表示的系数不是问题。另一个问题是,我仅使用2 ^ N个数据的边缘长方体。
从理论上讲,我认为应该可以确保不会发生“出血”。 William Press,Saul A teukolsky,William T. Vetterling和Brian P. Flannery(第二版)在“ C中的数值接收”中讨论了用于这种小波分解和重组的算法。尽管此算法假定在边缘反射而不是其他形式的边缘扩展(例如zpd),但该方法足够通用,可以用于其他形式的扩展。
关于如何将此工作扩展到其他小波的任何建议?
注意:此查询也发布在http://groups.google.com/group/pywavelets上
谢谢,
阿城
import pywt
import sys
import numpy as np
def waveFn(wavelet):
if not isinstance(wavelet, pywt.Wavelet):
return pywt.Wavelet(wavelet)
else:
return wavelet
# given a single dimensional array ... returns the coefficients.
def wavedec(data, wavelet, mode=\'sym\'):
wavelet = waveFn(wavelet)
dLen = len(data)
coeffs = np.zeros_like(data)
level = pywt.dwt_max_level(dLen, wavelet.dec_len)
a = data
end_idx = dLen
for idx in xrange(level):
a, d = pywt.dwt(a, wavelet, mode)
begin_idx = end_idx/2
coeffs[begin_idx:end_idx] = d
end_idx = begin_idx
coeffs[:end_idx] = a
return coeffs
def waverec(data, wavelet, mode=\'sym\'):
wavelet = waveFn(wavelet)
dLen = len(data)
level = pywt.dwt_max_level(dLen, wavelet.dec_len)
end_idx = 1
a = data[:end_idx] # approximation ... also the original data
d = data[end_idx:end_idx*2]
for idx in xrange(level):
a = pywt.idwt(a, d, wavelet, mode)
end_idx *= 2
d = data[end_idx:end_idx*2]
return a
def fn_dec(arr):
return np.array(map(lambda row: reduce(lambda x,y : np.hstack((x,y)), pywt.wavedec(row, \'haar\', \'zpd\')), arr))
# return np.array(map(lambda row: row*2, arr))
if __name__ == \'__main__\':
test = 1
np.random.seed(10)
wavelet = waveFn(\'haar\')
if test==0:
# SIngle dimensional test.
a = np.random.randn(1,8)
print \"original values A\"
print a
print \"decomposition of A by method in pywt\"
print fn_dec(a)
print \" decomposition of A by my method\"
coeffs = wavedec(a[0], \'haar\', \'zpd\')
print coeffs
print \"recomposition of A by my method\"
print waverec(coeffs, \'haar\', \'zpd\')
sys.exit()
if test==1:
a = np.random.randn(4,4,4)
# 2 D test
print \"original value of A\"
print a
# decompose the signal into wavelet coefficients.
dimensions = a.shape
for dim in dimensions:
a = np.rollaxis(a, 0, a.ndim)
ndim = a.shape
#a = fn_dec(a.reshape(-1, dim))
a = np.array(map(lambda row: wavedec(row, wavelet), a.reshape(-1, dim)))
a = a.reshape(ndim)
print \" decomposition of signal into coefficients\"
print a
# re-composition of the coefficients into original signal
for dim in dimensions:
a = np.rollaxis(a, 0, a.ndim)
ndim = a.shape
a = np.array(map(lambda row: waverec(row, wavelet), a.reshape(-1, dim)))
a = a.reshape(ndim)
print \"recomposition of coefficients to signal\"
print a
没有找到相关结果
已邀请:
1 个回复
末钉蹈泰唬
其中
是对两个维度(LL)都应用了近似变换的系数数组,而
是对第一个维度应用了细节变换的系数阵列,对第二维(HL)应用了细节变换的近似数组(与dwt2输出进行比较)。 基于此,将其扩展到多级案例应该相当容易。 这是我对分解部分的看法:https://gist.github.com/934166。 我也想解决您在问题中提到的一个问题: 有一个问题: 将小波系数表示为 与数据形状相同的数组。 在我看来,将结果表示为与输入数据具有相同形状/大小的数组的方法是有害的。这使整个事情变得不必要地复杂,难以理解和使用,因为无论如何您都必须进行假设或维护带有索引的二级数据结构,以便能够访问输出数组中的系数并执行逆变换(请参见Matlab的文档以获取相关信息)。 wavedec / waverec)。 同样,即使它在纸张上效果很好,但由于您提到的问题,它也不总是适合实际应用:大多数情况下,输入数据的大小不是2 ^ n,并且使用小波滤波器对信号进行卷积的抽取结果较大“存储空间”,进而可能导致数据丢失和重建不完善。 为了避免这些问题,我建议使用更自然的数据结构来表示结果数据层次结构,例如Python的列表,字典和元组(如果可用)。