On Thu, Mar 4, 2010 at 7:29 PM, Raul Miller <[email protected]> wrote: > Can you show me how to deal with the 1 dimensional case using > this style of approach, so I can focus on the differences between > the two implementations?
Let's take these numeric examples from http://rosettacode.org/wiki/Deconvolution/1D a=: _8 _9 _3 _1 _6 7 b=: _3 _6 _1 8 _6 3 _1 _9 _9 3 _2 5 2 _2 _7 _1 +//. a */ b 24 75 71 _34 3 22 _45 23 245 25 52 25 _67 _96 96 31 55 36 29 _43 _7 We can make a matrix bm from b such that the product of bm with a (thought of as a row vector) is the convolution c of a and b. This matrix bm does not depend on a except that it depends on its length. convmat=: 4 :'|:(-i.x)|."0 _(<:x+#y){.y' ]bm=: (#a) convmat b _3 0 0 0 0 0 _6 _3 0 0 0 0 _1 _6 _3 0 0 0 8 _1 _6 _3 0 0 _6 8 _1 _6 _3 0 3 _6 8 _1 _6 _3 _1 3 _6 8 _1 _6 _9 _1 3 _6 8 _1 _9 _9 _1 3 _6 8 3 _9 _9 _1 3 _6 _2 3 _9 _9 _1 3 5 _2 3 _9 _9 _1 2 5 _2 3 _9 _9 _2 2 5 _2 3 _9 _7 _2 2 5 _2 3 _1 _7 _2 2 5 _2 0 _1 _7 _2 2 5 0 0 _1 _7 _2 2 0 0 0 _1 _7 _2 0 0 0 0 _1 _7 0 0 0 0 0 _1 ]c=: +/"1 a *"1 bm 24 75 71 _34 3 22 _45 23 245 25 52 25 _67 _96 96 31 55 36 29 _43 _7 Now it's clear that we can get a back from b and c: we just divide c with the matrix bm from the right. bm-: (>:c-&#b) convmat b 1 ]av=: c %. bm _8 _9 _3 _1 _6 7 a-:av 1 Deconvolution in 1D done. The only problem is that the matrix is not square, so the equation system has too many equations. This is no wonder, for you don't get all vectors as a convolution of b. The two dimensional case works the same: we just construct a matrix bm so that the row vector (,a) multiplied by bm is the vector (,c) where c is the 2D convolution of a and b. One more note. The earlier deconvolution solution I gave fails if b starts with a zero, whereas this solution works in that degenerate case too. deconv=: 4 :('r=.$0'; 'while. x<:&#y'; 'do. r=.r,t=.y%&{.x'; 'y=.}.y-/@:,:t*x'; 'end.'; 'r') 0 1 3 deconv 0 1 13 30 0 _ 0 1 13 30 %. 2 convmat 0 1 3 1 10 +//. 1 10 */ 0 1 3 0 1 13 30 Similarly, the earlier solution I gave for 2D deconvolution fails in some degenerate special cases, but this matrix division one should work. Ambrus ---------------------------------------------------------------------- For information about J forums see http://www.jsoftware.com/forums.htm
