Integer type:  int32  int64  nag_int  show int32  show int32  show int64  show int64  show nag_int  show nag_int

Chapter Contents
Chapter Introduction
NAG Toolbox

# NAG Toolbox: nag_wav_3d_sngl_fwd (c09fa)

## Purpose

nag_wav_3d_sngl_fwd (c09fa) computes the three-dimensional discrete wavelet transform (DWT) at a single level. The initialization function nag_wav_3d_init (c09ac) must be called first to set up the DWT options.

## Syntax

[c, icomm, ifail] = c09fa(n, fr, a, lenc, icomm, 'm', m)
[c, icomm, ifail] = nag_wav_3d_sngl_fwd(n, fr, a, lenc, icomm, 'm', m)

## Description

nag_wav_3d_sngl_fwd (c09fa) computes the three-dimensional DWT of a given input three-dimensional data array, considered as a number of two-dimensional frames, at a single level. For a chosen wavelet filter pair, the output coefficients are obtained by applying convolution and downsampling by two to the input, A$A$, first over columns, next over rows and finally across frames. The three-dimensional approximation coefficients are produced by the low pass filter over columns, rows and frames. In addition there are 7 sets of three-dimensional detail coefficients, each corresponding to a different order of low pass and high pass filters (see the C09 Chapter Introduction). All coefficients are packed into a single array. To reduce distortion effects at the ends of the data array, several end extension methods are commonly used. Those provided are: periodic or circular convolution end extension, half-point symmetric end extension, whole-point symmetric end extension and zero end extension. The total number, nct${n}_{\mathrm{ct}}$, of coefficients computed is returned by the initialization function nag_wav_3d_init (c09ac).

## References

Daubechies I (1992) Ten Lectures on Wavelets SIAM, Philadelphia

## Parameters

### Compulsory Input Parameters

1:     n – int64int32nag_int scalar
The second dimension of the input data: the number of columns of each two-dimensional frame.
Constraint: this must be the same as the value n passed to the initialization function nag_wav_3d_init (c09ac).
2:     fr – int64int32nag_int scalar
The third dimension of the input data: the number of two-dimensional frames.
Constraint: this must be the same as the value fr passed to the initialization function nag_wav_3d_init (c09ac).
3:     a(lda,sda,fr) – double array
lda, the first dimension of the array, must satisfy the constraint ldam$\mathit{lda}\ge {\mathbf{m}}$.
The m$m$ by n$n$ by fr$\mathit{fr}$ input three-dimensional array A$A$.
4:     lenc – int64int32nag_int scalar
The dimension of the array c as declared in the (sub)program from which nag_wav_3d_sngl_fwd (c09fa) is called.
Constraint: lencnct${\mathbf{lenc}}\ge {n}_{\mathrm{ct}}$, where nct${n}_{\mathrm{ct}}$ is the total number of wavelet coefficients, as returned by nag_wav_3d_init (c09ac).
5:     icomm(260$260$) – int64int32nag_int array
Contains details of the discrete wavelet transform and the problem dimension as setup in the call to the initialization function nag_wav_3d_init (c09ac).

### Optional Input Parameters

1:     m – int64int32nag_int scalar
Default: The first dimension of the array a.
The first dimension of the input data: the number of rows of each two-dimensional frame.
Constraint: this must be the same as the value m passed to the initialization function nag_wav_3d_init (c09ac).

lda sda

### Output Parameters

1:     c(lenc) – double array
The coefficients of the discrete wavelet transform. The 8$8$ sets of coefficients are stored in the following order: approximation coefficients (LLL) first, followed by 7$7$ sets of detail coefficients: LLH, LHL, LHH, HLL, HLH, HHL, HHH, where L indicates the low pass filter, and H the high pass filter being applied to, respectively, the columns of length m, the rows of length n and then the frames of length fr. Note that for computational efficiency reasons each set of coefficients is stored in the order ncfr × ncm × ncn${n}_{\mathrm{cfr}}×{n}_{\mathrm{cm}}×{n}_{\mathrm{cn}}$ (see output parameters nwcfr, nwct and nwcn in nag_wav_3d_init (c09ac)). See Section [Example] for details of how to access each set of coefficients in order to perform extraction from c following a call to this function, or insertion into c before a call to the three-dimensional inverse routine nag_wav_3d_sngl_inv (c09fb).
2:     icomm(260$260$) – int64int32nag_int array
Contains additional information on the computed transform.
3:     ifail – int64int32nag_int scalar
${\mathrm{ifail}}={\mathbf{0}}$ unless the function detects an error (see [Error Indicators and Warnings]).

## Error Indicators and Warnings

Errors or warnings detected by the function:
ifail = 1${\mathbf{ifail}}=1$
Constraint: fr = fr${\mathbf{fr}}=\mathrm{fr}$, the value of fr on initialization (see nag_wav_3d_init (c09ac)).
Constraint: m = m${\mathbf{m}}=m$, the value of m on initialization (see nag_wav_3d_init (c09ac)).
Constraint: n = n${\mathbf{n}}=n$, the value of n on initialization (see nag_wav_3d_init (c09ac)).
ifail = 2${\mathbf{ifail}}=2$
Constraint: ldam$\mathit{lda}\ge {\mathbf{m}}$.
Constraint: sdan$\mathit{sda}\ge {\mathbf{n}}$.
ifail = 3${\mathbf{ifail}}=3$
Constraint: lencnct${\mathbf{lenc}}\ge {n}_{\mathrm{ct}}$, where nct${n}_{\mathrm{ct}}$ is the number of DWT coefficients returned by nag_wav_3d_init (c09ac) in parameter nwct.
ifail = 6${\mathbf{ifail}}=6$
Either the initialization function nag_wav_3d_init (c09ac) has not been called first or the communication array icomm has been corrupted.
The initialization function was called with wtrans = 'M'${\mathbf{wtrans}}=\text{'M'}$.
ifail = 999${\mathbf{ifail}}=-999$
Dynamic memory allocation failed.

## Accuracy

The accuracy of the wavelet transform depends only on the floating point operations used in the convolution and downsampling and should thus be close to machine precision.

None.

## Example

function nag_wav_3d_sngl_fwd_example
m  = int64(5);
n  = int64(4);
fr = int64(3);
wavnam = 'Haar';
mode = 'half';
wtrans = 'Single Level';
a = zeros(m, n, fr);
a(:, :, 1) = [3, 2, 2, 2;
2, 9, 1, 2;
2, 5, 1, 2;
1, 6, 2, 2;
5, 3, 2, 2];
a(:, :, 2) = [2, 1, 5, 1;
2, 9, 5, 2;
2, 3, 2, 7;
2, 1, 1, 2;
2, 1, 2, 8];
a(:, :, 3) = [3, 1, 4, 1;
1, 1, 2, 1;
4, 1, 7, 2;
3, 2, 1, 5;
1, 1, 2, 2];

% Query wavelet filter dimensions
[lmax, nf, nwct, nwcn, nwcfr, icomm, ifail] = ...
nag_wav_3d_init(wavnam, wtrans, mode, m, n, fr);

nwcm = nwct/(8*nwcn*nwcfr);

% 3D DWT decomposition
[c, icomm, ifail] = nag_wav_3d_sngl_fwd(n, fr, a, nwct, icomm);

d = zeros(nwcm, nwcn, nwcfr);

for cindex = 0:7

% Decide which coefficient type we are considering and advance the
% pointer locc to the first element of that 3D array in C.
switch (cindex)
case {0}
fprintf('Approximation coefficients (LLL)\n');
locc = 1;
case {1}
fprintf('Detail coefficients (LLH)\n');
% Advance pointer past approximation coefficients
locc = nwcm*nwcn*nwcfr + 1;
case {2}
fprintf('Detail coefficients (LHL)\n');
% Advance pointer past approximation coefficients and 1 set of
% detail coefficients
locc = 2*nwcm*nwcn*nwcfr + 1;
case {3}
fprintf('Detail coefficients (LHH)\n');
% Advance pointer past approximation coefficients and 2 sets of
% detail coefficients
locc = 3*nwcm*nwcn*nwcfr + 1;
case {4}
fprintf('Detail coefficients (HLL)\n');
% Advance pointer past approximation coefficients and 3 sets of
% detail coefficients
locc = 4*nwcm*nwcn*nwcfr + 1;
case {5}
fprintf('Detail coefficients (HLH)\n');
% Advance pointer past approximation coefficients and 4 sets of
% detail coefficients
locc = 5*nwcm*nwcn*nwcfr + 1;
case {6}
fprintf('Detail coefficients (HHL)\n');
% Advance pointer past approximation coefficients and 5 sets of
% detail coefficients
locc = 6*nwcm*nwcn*nwcfr + 1;
case {7}
fprintf('Detail coefficients (HHH)\n');
% Advance pointer past approximation coefficients and 6 sets of
% detail coefficients
locc = 7*nwcm*nwcn*nwcfr + 1;
end

for k = 1:nwcfr
for j = 1:nwcn
for i = 1:nwcm
i1 = locc - 1 + (j-1)*nwcfr*nwcm + (i-1)*nwcfr + k;
d(i,j,k) = c(i1);
end
end
end

for j = 1:nwcfr
if (j==1)
fprintf('Coefficients        Frame 1');
else
fprintf('          Frame %d', j);
end
end
fprintf('\n');
d2 = reshape(d, nwcm, nwcn*nwcfr);
for i = 1:nwcm
if i == 1
fprintf('    %d         ', cindex);
else
fprintf('              ');
end
for j=1:nwcn*nwcfr
fprintf('%8.4f ', d2(i,j));
end
fprintf('\n');
end
end

% 3D DWT reconstruction
[b, ifail] = nag_wav_3d_sngl_inv(m, n, fr, c, icomm);

fprintf('\nOutput Data          b : \n');
% Convert to int16 to get more compact output
for i=1:nwcm
fprintf('Frame %d :\n', i);
disp(b(:, :, i));
end

Approximation coefficients (LLL)
Coefficients        Frame 1          Frame 2
0          10.6066   7.0711   4.2426   5.6569
7.7782   6.7175   7.0711  10.6066
7.7782   9.8995   2.8284   5.6569
Detail coefficients (LLH)
Coefficients        Frame 1          Frame 2
1           0.7071  -2.1213   0.0000   0.0000
2.1213  -1.7678   0.0000   0.0000
3.5355  -4.2426   0.0000   0.0000
Detail coefficients (LHL)
Coefficients        Frame 1          Frame 2
2          -4.2426   2.1213   1.4142   2.8284
-2.8284  -2.4749   2.8284   0.7071
2.1213  -4.2426   0.0000   0.0000
Detail coefficients (LHH)
Coefficients        Frame 1          Frame 2
3           0.0000  -2.8284   0.0000   0.0000
-2.8284   1.7678   0.0000   0.0000
0.7071   4.2426   0.0000   0.0000
Detail coefficients (HLL)
Coefficients        Frame 1          Frame 2
4          -4.9497   0.0000   1.4142   1.4142
0.7071   1.7678  -0.0000   2.1213
0.0000   0.0000   0.0000   0.0000
Detail coefficients (HLH)
Coefficients        Frame 1          Frame 2
5           0.7071   0.7071   0.0000   0.0000
-0.7071  -2.4749   0.0000   0.0000
0.0000   0.0000   0.0000   0.0000
Detail coefficients (HHL)
Coefficients        Frame 1          Frame 2
6           5.6569   0.7071   1.4142   1.4142
0.0000  -1.7678   1.4142   6.3640
0.0000   0.0000   0.0000   0.0000
Detail coefficients (HHH)
Coefficients        Frame 1          Frame 2
7           0.0000   0.0000   0.0000   0.0000
1.4142   1.0607   0.0000   0.0000
0.0000   0.0000   0.0000   0.0000

Output Data          b :
Frame 1 :
3.0000    2.0000    2.0000    2.0000
2.0000    9.0000    1.0000    2.0000
2.0000    5.0000    1.0000    2.0000
1.0000    6.0000    2.0000    2.0000
5.0000    3.0000    2.0000    2.0000

Frame 2 :
2.0000    1.0000    5.0000    1.0000
2.0000    9.0000    5.0000    2.0000
2.0000    3.0000    2.0000    7.0000
2.0000    1.0000    1.0000    2.0000
2.0000    1.0000    2.0000    8.0000

Frame 3 :
3.0000    1.0000    4.0000    1.0000
1.0000    1.0000    2.0000    1.0000
4.0000    1.0000    7.0000    2.0000
3.0000    2.0000    1.0000    5.0000
1.0000    1.0000    2.0000    2.0000

function c09fa_example
m  = int64(5);
n  = int64(4);
fr = int64(3);
wavnam = 'Haar';
mode = 'half';
wtrans = 'Single Level';
a = zeros(m, n, fr);
a(:, :, 1) = [3, 2, 2, 2;
2, 9, 1, 2;
2, 5, 1, 2;
1, 6, 2, 2;
5, 3, 2, 2];
a(:, :, 2) = [2, 1, 5, 1;
2, 9, 5, 2;
2, 3, 2, 7;
2, 1, 1, 2;
2, 1, 2, 8];
a(:, :, 3) = [3, 1, 4, 1;
1, 1, 2, 1;
4, 1, 7, 2;
3, 2, 1, 5;
1, 1, 2, 2];

% Query wavelet filter dimensions
[lmax, nf, nwct, nwcn, nwcfr, icomm, ifail] = ...
c09ac(wavnam, wtrans, mode, m, n, fr);

nwcm = nwct/(8*nwcn*nwcfr);

% 3D DWT decomposition
[c, icomm, ifail] = c09fa(n, fr, a, nwct, icomm);

d = zeros(nwcm, nwcn, nwcfr);

for cindex = 0:7

% Decide which coefficient type we are considering and advance the
% pointer locc to the first element of that 3D array in C.
switch (cindex)
case {0}
fprintf('Approximation coefficients (LLL)\n');
locc = 1;
case {1}
fprintf('Detail coefficients (LLH)\n');
% Advance pointer past approximation coefficients
locc = nwcm*nwcn*nwcfr + 1;
case {2}
fprintf('Detail coefficients (LHL)\n');
% Advance pointer past approximation coefficients and 1 set of
% detail coefficients
locc = 2*nwcm*nwcn*nwcfr + 1;
case {3}
fprintf('Detail coefficients (LHH)\n');
% Advance pointer past approximation coefficients and 2 sets of
% detail coefficients
locc = 3*nwcm*nwcn*nwcfr + 1;
case {4}
fprintf('Detail coefficients (HLL)\n');
% Advance pointer past approximation coefficients and 3 sets of
% detail coefficients
locc = 4*nwcm*nwcn*nwcfr + 1;
case {5}
fprintf('Detail coefficients (HLH)\n');
% Advance pointer past approximation coefficients and 4 sets of
% detail coefficients
locc = 5*nwcm*nwcn*nwcfr + 1;
case {6}
fprintf('Detail coefficients (HHL)\n');
% Advance pointer past approximation coefficients and 5 sets of
% detail coefficients
locc = 6*nwcm*nwcn*nwcfr + 1;
case {7}
fprintf('Detail coefficients (HHH)\n');
% Advance pointer past approximation coefficients and 6 sets of
% detail coefficients
locc = 7*nwcm*nwcn*nwcfr + 1;
end

for k = 1:nwcfr
for j = 1:nwcn
for i = 1:nwcm
i1 = locc - 1 + (j-1)*nwcfr*nwcm + (i-1)*nwcfr + k;
d(i,j,k) = c(i1);
end
end
end

for j = 1:nwcfr
if (j==1)
fprintf('Coefficients        Frame 1');
else
fprintf('          Frame %d', j);
end
end
fprintf('\n');
d2 = reshape(d, nwcm, nwcn*nwcfr);
for i = 1:nwcm
if i == 1
fprintf('    %d         ', cindex);
else
fprintf('              ');
end
for j=1:nwcn*nwcfr
fprintf('%8.4f ', d2(i,j));
end
fprintf('\n');
end
end

% 3D DWT reconstruction
[b, ifail] = c09fb(m, n, fr, c, icomm);

fprintf('\nOutput Data          b : \n');
% Convert to int16 to get more compact output
for i=1:nwcm
fprintf('Frame %d :\n', i);
disp(b(:, :, i));
end

Approximation coefficients (LLL)
Coefficients        Frame 1          Frame 2
0          10.6066   7.0711   4.2426   5.6569
7.7782   6.7175   7.0711  10.6066
7.7782   9.8995   2.8284   5.6569
Detail coefficients (LLH)
Coefficients        Frame 1          Frame 2
1           0.7071  -2.1213   0.0000   0.0000
2.1213  -1.7678   0.0000   0.0000
3.5355  -4.2426   0.0000   0.0000
Detail coefficients (LHL)
Coefficients        Frame 1          Frame 2
2          -4.2426   2.1213   1.4142   2.8284
-2.8284  -2.4749   2.8284   0.7071
2.1213  -4.2426   0.0000   0.0000
Detail coefficients (LHH)
Coefficients        Frame 1          Frame 2
3           0.0000  -2.8284   0.0000   0.0000
-2.8284   1.7678   0.0000   0.0000
0.7071   4.2426   0.0000   0.0000
Detail coefficients (HLL)
Coefficients        Frame 1          Frame 2
4          -4.9497   0.0000   1.4142   1.4142
0.7071   1.7678  -0.0000   2.1213
0.0000   0.0000   0.0000   0.0000
Detail coefficients (HLH)
Coefficients        Frame 1          Frame 2
5           0.7071   0.7071   0.0000   0.0000
-0.7071  -2.4749   0.0000   0.0000
0.0000   0.0000   0.0000   0.0000
Detail coefficients (HHL)
Coefficients        Frame 1          Frame 2
6           5.6569   0.7071   1.4142   1.4142
0.0000  -1.7678   1.4142   6.3640
0.0000   0.0000   0.0000   0.0000
Detail coefficients (HHH)
Coefficients        Frame 1          Frame 2
7           0.0000   0.0000   0.0000   0.0000
1.4142   1.0607   0.0000   0.0000
0.0000   0.0000   0.0000   0.0000

Output Data          b :
Frame 1 :
3.0000    2.0000    2.0000    2.0000
2.0000    9.0000    1.0000    2.0000
2.0000    5.0000    1.0000    2.0000
1.0000    6.0000    2.0000    2.0000
5.0000    3.0000    2.0000    2.0000

Frame 2 :
2.0000    1.0000    5.0000    1.0000
2.0000    9.0000    5.0000    2.0000
2.0000    3.0000    2.0000    7.0000
2.0000    1.0000    1.0000    2.0000
2.0000    1.0000    2.0000    8.0000

Frame 3 :
3.0000    1.0000    4.0000    1.0000
1.0000    1.0000    2.0000    1.0000
4.0000    1.0000    7.0000    2.0000
3.0000    2.0000    1.0000    5.0000
1.0000    1.0000    2.0000    2.0000