Transient convolution

Matrix multiplication and transpose multiplication still fit easily in the same computational framework when the matrix has a special form, such as  

$$\bold y \eq \left[ \begin{array} {c} y_1 \\ y_2 \\ y... ..._4 \\ x_5 \\ x_6 \end{array} \right] \eq \bold B \bold x$$ (4)

Notice that columns of equation (4) all contain the same signal, but with different shifts. This signal is called the filter's impulse response.

Equation (4) could be rewritten as  

$$\bold y \eq \left[ \begin{array} {c} y_1 \\ y_2 \\ y... ..._1 \\ b_2 \\ b_3 \end{array} \right] \eq \bold X \bold b$$ (5)

In applications we can choose between $\bold y = \bold X \bold b$and $\bold y=\bold B\bold x$.In one case the output $\bold y$is dual to the filter $\bold b$,and in the other case the output $\bold y$is dual to the input $\bold x$.Sometimes we must solve for $\bold b$ and sometimes for $\bold x$;so sometimes we use equation (5) and sometimes (4). Such solutions begin from the adjoints. The adjoint of (4) is  

$$\left[ \begin{array} {c} \hat x_1 \\ \hat x_2 \\ \hat x... ...y_4 \\ y_5 \\ y_6 \\ y_7 \\ y_8 \end{array} \right]$$ (6)

The adjoint crosscorrelates with the filter instead of convolving with it (because the filter is backwards). Notice that each row in equation (6) contains all the filter coefficients and there are no rows where the filter somehow uses zero values off the ends of the data as we saw earlier. In some applications it is important not to assume zero values beyond the interval where where is given.

The adjoint of (5) crosscorrelates a fixed portion of filter input across a variable portion of filter output.  

$$\left[ \begin{array} {c} \hat b_1 \\ \hat b_2 \\ \hat b... ..._4 \\ y_5 \\ y_6 \\ y_7 \\ y_8 \end{array} \right]$$ (7)

Module tcai1 is used for $\bold y=\bold B\bold x$and module tcaf1 is used for $\bold y = \bold X \bold b$. 

module tcai1 {            # Transient Convolution Adjoint Input 1-D. yy(m1+n1)
real, dimension (:), pointer :: bb
integer                      :: nx
#%  _init( bb, nx)
#%  _lop ( xx, yy)
integer b, x, y
if( size(yy) < nx + size(bb) - 1 ) call erexit('tcai')
do b= 1, size(bb) {
do x= 1, nx       {                 y = x + b - 1
        if( adj)
                        xx(x) = xx(x) + yy(y) * bb(b)
	else	
                        yy(y) = yy(y) + xx(x) * bb(b)
        }}
}

 

module tcaf1 {              # Transient Convolution, Adjoint is the Filter, 1-D
real, dimension (:), pointer :: xx
#%  _init( xx)
#%  _lop ( bb,  yy)
integer        x,   b,   y
if( size(yy) < size(xx) + size(bb) - 1 )   call erexit('tcaf')
do b= 1, size(bb) {
do x= 1, size(xx) {                      y = x + b - 1
        if( adj)
                        bb(b) = bb(b) + yy(y) * xx(x)
	else	
                        yy(y) = yy(y) + bb(b) * xx(x)
        } }
}