MODWT IMODWT code without fft

3 Ansichten (letzte 30 Tage)
Emiliano Rosso
Emiliano Rosso am 21 Mär. 2023
Bearbeitet: Emiliano Rosso am 12 Apr. 2023
Hi.
I studied the maximal overlap wavelet transform and its properties on "Wavelet Methods for Time Series Analysis by Donald B. Percival, Andrew T. Walden " and I saw that the application of the fft is performed only after the development of the algorithm for speed up the code. To better understand how it works I need to develop this algorithm without using the fft. Is there any library that does this?
Thank you!
Now i found the modwt code :
function [coeffs] = NOFFTmodwt(x, filters, level)
% filters: 'db1', 'db2', 'haar', ...
% return: see matlab
% filter
N=length(x);
[Lo_D,Hi_D] = wfilters(filters);
g_t = Lo_D / sqrt(2);
h_t = Hi_D / sqrt(2);
wavecoeff= cell(1,5);
v_j_1 = x;
for j = 1:level
w = circular_convolve_d(h_t, v_j_1, j);
v_j_1 = circular_convolve_d(g_t, v_j_1, j);
wavecoeff{j} = w;
end
wavecoeff{level+1} = v_j_1;
coeffs = vertcat(wavecoeff{:});
coeffs(1:level,1:N)=-coeffs(1:level,1:N);
function w_j = circular_convolve_d(h_t, v_j_1, j)
% jth level decomposition
% h_t: \tilde{h} = h / sqrt(2)
% v_j_1: v_{j-1}, the (j-1)th scale coefficients
% return: w_j (or v_j)
N = length(v_j_1);
L = length(h_t);
w_j = zeros(1,N);
l = 0:L-1;
for t = 1:N
index = mod(t -1 - 2^(j-1)*l, N) + 1;
v_p = v_j_1(index);
w_j(t) = dot(h_t, v_p);
end
It's freely translated from the Python function modwtpy which you can find here :
It's equal in results to the Matlab modwt function , border sx oriented for the 'haar' , shit invariante , same coeffs, ecc..
Now I must find the same for the imodwt. this is what i managed to do but it doesn't work!
function [x] = NOFFTimodwt(coeffs, filters)
% inverse modwt
% filter
[Lo_R,Hi_R] = wfilters(filters);
g_t = Lo_R / sqrt(2);
h_t = Hi_R / sqrt(2);
level = size(coeffs, 1) - 1;
v_j = coeffs(end, :);
for j = level:-1:1
%w_j = circshift(coeffs(j,:),[0,2^(j-1)]);
%v_j = circular_convolve_s(h_t, g_t, w_j, v_j, j);
v_j = circular_convolve_s(h_t, g_t, coeffs(j,:), v_j, j);
end
x = v_j;
function v_j_1 = circular_convolve_s(h_t, g_t, w_j, v_j, j)
% (j-1)th level synthesis from w_j, w_j
% see function circular_convolve_d
N = length(v_j);
L = length(h_t);
v_j_1 = zeros(1,N);
l = 0:L-1;
for t = 1:N
index = mod(t -1 + 2^(j-1)*l, N) + 1;
%index = mod(t + 2^(j-1)*l - 1, N) + 1;
w_p = w_j(index);
v_p = v_j(index);
v_j_1(t) = dot(h_t, w_p) + dot(g_t, v_p);
end
Any suggestions?
Thanks!

Antworten (0)

Kategorien

Mehr zu Filter Banks finden Sie in Help Center und File Exchange

Tags

Produkte


Version

R2020b

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by