The conjugate of extracting all the patches is adding them back. Because of overlaps, the conjugate is not the inverse. In many applications, inverse parcelling is required, i.e. patching things back together seamlessly. This can be done with weighting functions. You can have any weighting function you wish and I will provide you a reconstruction like
(1) |
# make wall weights from window weights. # subroutine mkwallwt( k1,k2, windwt, w1,w2, wallwt, n1,n2) integer i1,i2, k1,k2, w1,w2, n1,n2 real windwt( w1,w2), wallwt( n1,n2) call null( wallwt, n1*n2) do i1= 1, k1 { do i2= 1, k2 { call patch( 1, 1, i1,i2, k1,k2, wallwt, n1,n2, windwt,w1,w2 ) }} do i2= 1, n2 { do i1= 1, n1 { if( wallwt(i1,i2) != 0.) wallwt(i1,i2) = 1. / wallwt(i1,i2) }} return; end
No matrices are needed to show this method succeeds because data values are never mixed with one another. An equation for any reconstructed data value as a function of the original value d and the weights w_{i} that hit d is .
To demonstrate the program, I begin by making a random weighting function for each window with subroutine mkrandwt(). The only requirement is that the numbers be positive. Zero valued weights are OK too, but you would lose data there unless it was picked up by an overlapping window. The general strategy allows different weights in different windows. That flexibility adds clutter, however, so here we work with the same weighting function in each window.
# Put positive random numbers in a window. # subroutine mkrandwt( windwt, w1,w2) integer i1,i2, w1,w2, iseed real windwt( w1,w2), rand01 iseed = 1992 do i1= 1, w1 { do i2= 1, w2 { windwt(i1,i2) = rand01( iseed) }} return; end
Typically, we apply some linear operator, say in each patch. Mathematically, this is:
(2) |
# idempatch -- cut up section and paste back together, idempotently. # subroutine idempatch( conj, add, w1,w2, k1,k2, data, n1,n2, fout) integer i1,i2, conj, add, w1,w2, k1,k2, n1,n2 real data( n1,n2),fout( n1,n2) temporary real pdat(w1,w2), copy(n1,n2) temporary real pfou(w1,w2), wallwt(n1,n2), windwt(w1,w2) call conjnull( conj, add, data, n1*n2, fout, n1*n2) call mkrandwt( windwt, w1,w2) call mkwallwt( k1,k2, windwt, w1,w2, wallwt, n1,n2) if( conj == 0 ) { do i2= 1, k2 { do i1= 1, k1 { call patch( 0, 0, i1,i2, k1,k2, data,n1,n2, pdat,w1,w2 ) call ident( 0, 0, 1.0, w1*w2, pdat, pfou ) call diag( 0, 0, windwt, w1*w2, pfou, pfou ) call patch( 1, 1, i1,i2, k1,k2, copy,n1,n2, pfou,w1,w2 ) }} call diag( 0, 1, wallwt, n1*n2, copy, fout ) } else { call diag( 1, 0, wallwt, n1*n2, copy, fout ) do i2= 1, k2 { do i1= 1, k1 { call patch( 0, 0, i1,i2, k1,k2, copy,n1,n2, pfou,w1,w2 ) call diag( 1, 0, windwt, w1*w2, pfou, pfou ) call ident( 1, 0, 1.0, w1*w2, pdat, pfou ) call patch( 1, 1, i1,i2, k1,k2, data,n1,n2, pdat,w1,w2 ) }} } return; end
Figure 2 shows the result when a plane of identical values is passed through the idempatch() subroutine. The result is constant on the 2-axis which confirms that despite the random function as a weight, idempatch() performs as advertised and that there is adequate sampling on the 2-axis. On the 1-axis the output is constant except for being zero in gaps because the windows do not overlap on the 1-axis. Curiously, the conjugate operation produces the same result, i.e. we can do the operations in reverse order.
idempatch
Figure 2 A plane of identical values passed through subroutine idempatch(). Results are shown for the same parameters as Figure 1. |