function X = pm_haar(x,s)
%
% MATLAB function for pm Haar 1-D DWT algorithm.
% This function receives N (a power of two) real values in x
% and the number of stages (scales) of decomposition required in s,
% computes the Haar DWT, and returns the N DWT coefficients in X.
%
% For the theory behind this algorithm and example input/output,
% please refer the paper:
% Fundamentals of the discrete Haar wavelet transform
% Duraisamy Sundararajan
% dsprelated.com, 2011
% articles/paper section
%
N = length(x); % length of the input vector
b = N;sq2 = sqrt(2);fac = 1;
for ns =1:s % outer loop stepping over stages
p = 1;q = 1;
while(p < b) % inner loop stepping over each butterfly in a stage
r = p + 1;
xmr(q) = x(p) - x(r); xpr(q) = x(p) + x(r);
p = p + 2; q = q + 1;
end
fac = fac / sq2;
x(1:b) = [xpr(1:b/2) xmr(1:b/2) * fac];
b = b / 2;
end
x(1:b) = fac * x(1:b);
X=x;
function x = pm_haar_inv(X,s)
%
% MATLAB function for pm Haar 1-D IDWT algorithm.
% This function receives N (a power of two) Haar DWT coefficients in X
% and the number of stages (scales) decomposed in s,
% computes the Haar IDWT, and returns the N real values in x.
%
% For the theory behind this algorithm and example input/output,
% please refer the paper:
% Fundamentals of the discrete Haar wavelet transform
% Duraisamy Sundararajan
% dsprelated.com, 2011
% articles/paper section
%
N = length(X); % length of the input vector
b = N/2^s;sq2 = sqrt(2);fac=1/(sq2^s);
X(1:b) = X(1:b) * fac;
for ns =1:s % outer loop stepping over stages
p = 1;q = 1;
xpr(1:b) = X(1:b);
xmr(1:b) = X(b+1:2*b) * fac;
while(p <= b) % inner loop stepping over each butterfly in a stage
X(q + 1) = xpr(p) - xmr(p); X(q) = xpr(p) + xmr(p);
p = p + 1; q = q + 2;
end
fac = sq2 * fac;
b = b * 2;
end
x = X;
% Main program for calling functions to compute forward and inverse
% Haar 1-D DWT
% functions pm_haar.m and pm_haar_inv.m are required in the same directory
clear;
% Example input x = [0 1 2 3 4 5 6 7], s = 3
% Example output X = [9.8995, -5.6569, -2, -2,
% -0.7071, -0.7071, -0.7071, -0.7071]
%
x = [0 1 2 3 4 5 6 7];
s = 3; % number of stages of computation required, s = 1,2,...,log2(N),
% where N is the number of elements in the input vector x. N is an
% integral power of two
X = pm_haar(x,s); % Invoke pm_haar.m for
% computing the forward DWT
X % the Haar DWT of input vector x after s stages of computation
x = pm_haar_inv(X,s); % Invoke pm_haar_inv.m
% for computing the IDWT
x % the given input after s stages of inverse computation
% Main program for calling functions to compute forward and inverse
% 2-D Haar DWT
% functions pm_haar.m and pm_haar_inv.m are required in the same directory
clear;
% Example input x2 = [
% 1 2 -1 3
% 2 1 4 3
% 1 1 2 2
% 4 2 1 3
% ]
% scale = 2
% Example output X = [
% 7.7500 -0.7500 0 -1.5000
%-0.2500 -0.7500 1.0000 -1.0000
% 0 -2.5000 -1.0000 -2.5000
%-2.0000 0 -1.0000 1.0000
% ]
x2 = [1 2 -1 3;2 1 4 3;1 1 2 2;4 2 1 3];
N = size(x2,1); % size of each row, N must be an integral power of two
scale = 2; % number of stages (scales) of computation required,
% s = 1,2,...,log2(N)
lim = N; % size of the row or column of the subimage in different stages
% computation of the forward 2-D Haar DWT of x2 by row-column method
for k=1:scale % compute the 2-D DWT scale by scale
%
for j=1:lim
x = x2(j,1:lim); % copy each row of input
X = pm_haar(x,1); % Invoke pm_haar.m to compute each row DWT
x2(j,1:lim) = X;
end
%
for j = 1:lim
x = x2(1:lim,j); % copy each column of the result of row DFTs
X = pm_haar(x',1); % Invoke pm_haar.m to compute each column DWT
x2(1:lim,j) = X';
end
lim = lim / 2; % transform size reduced for each scale
end
x2 % 2-D DWT of x2
% computation of the 2-D Haar IDWT of x2 by row-column method
lim = N/(2^(scale-1)); % size of the row or column of the subimage in
% different stages computation of the 2-D Haar IDWT of x2 by
% row-column method
for k = 1:scale % compute the 2-D IDWT scale by scale
for j = 1:lim
X = x2(j,1:lim); % copy each row of the coefficients
x = pm_haar_inv(X,1); % Invoke pm_haar_inv.m to compute each row IDWTs
x2(j,1:lim) = x;
end
%
for j = 1:lim
X = x2(1:lim,j); % copy each column of the results of row IDWTs
x = pm_haar_inv(X',1); % Invoke pm_haar_inv.m to compute column IDWTs
x2(1:lim,j) = x';
end
lim = lim * 2; % transform size increased for each scale
end
%
x2 % the given input after required stages of IDWT computa
