Это может быть вопрос для другого форума, если да, пожалуйста, дайте мне знать. Я заметил, что на тег вейвлета подписано всего 14 человек.
У меня есть элегантный способ расширить разложение вейвлета в pywt (пакет pyWavelets) на несколько измерений. Это должно работать из коробки, если установлен pywt. Тест 1 показывает декомпозицию и рекомпозицию трехмерного массива. Все, что нужно сделать, это увеличить количество измерений, и код будет работать при декомпозиции/перекомпоновке с 4, 6 или даже 18 измерениями данных.
Здесь я заменил функции pywt.wavedec и pywt.waverec. Кроме того, в fn_dec я показываю, как новая функция wavedec работает так же, как и старая.
Однако есть одна загвоздка: он представляет вейвлет-коэффициенты в виде массива той же формы, что и данные. Как следствие, с моими ограниченными знаниями о вейвлетах я смог использовать его только для вейвлетов Хаара. Другие, такие как DB4, например, выводят коэффициенты за края этих строгих границ (это не проблема с текущим представлением коэффициентов в виде списка массивов [CA, CD1 ... CDN]. Еще одна загвоздка в том, что я работал только с 2 ^N реберные кубоиды данных.
Теоретически, я думаю, должна быть возможность убедиться, что "кровотечения" не происходит. Алгоритм такого рода вейвлет-разложения и рекомпозиции обсуждается в «Численных рецептах на языке 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
