In each recursion step we filter an input vector by
multiplication with a quadratic matrix *C*:

(1) |

The matrix is very similar to a convolution matrix. A wraparound boundary condition is implemented in the lower right corner.

Let us imagine the result *m*_{t} of the multiplication of each odd
row with the vector *B* as a weighted local average of adjacent
vector coefficients *B*_{t}.
The multiplication of each even row with the vector *B* results in
a number *w*_{t} which measures the variation of these adjacent vector
coefficients.
Each pair *m*_{t}, *w*_{t} samples the same location in
vector *B*. I refer to *m*_{t}, *w*_{t} as average or lowpass and
detail or highpass coefficients, respectively.
The distance between locations sampled by adjacent
output pairs is two *B*_{t} samples;
for example the pair *m _{1}*,

Obviously, not every set of *c*_{k} filter coefficients has such a
highpass/lowpass characteristic. The coefficients are
correspondingly chosen (see the section ``The lowpass filter coefficients'').
In order for the DWT to be easily invertible, the matrix
*C* must be orthogonal. The inverse matrix is then the transpose of *C*:

(2) |

This orthogonality condition further constrains the choice of the wavelet
filter coefficients *c*_{k}, as I discuss in the section ``The low
pass filter coefficients.''
The subroutine `dwt.rt` executes the recursion step.

# #% subroutine dwt(conj, b,n, c,nop) integer conj, n, nop, i,j, no2 real b(n),c(nop) temporary real tmp(n)no2 = int(n/2) if(conj != 0) Do j = 1, no2{ tmp(2*j-1) = b(j ) tmp(2*j ) = b(j+no2) } call conjnull(conj, 0, b,n, tmp,n)

Do j = 1, n - 1, 2 Do i = 1, nop if(conj == 0){ tmp(j ) = tmp(j ) + c( i ) * b(mod(j+i-2,n) + 1) tmp(j+1) = tmp(j+1) - (-1)**i * c(nop - i + 1) * b(mod(j+i-2,n) + 1) } else{ b(mod(j+i-2,n) + 1) = b(mod(j+i-2,n) + 1) + c( i ) * tmp(j ) b(mod(j+i-2,n) + 1) = b(mod(j+i-2,n) + 1) -(-1)**i * c(nop-i+1) * tmp(j+1) } if(conj == 0) Do j = 1, no2{ b(j ) = tmp(2*j-1) b(j+no2) = tmp(2*j ) } return; end

11/18/1997