Fri Feb 1 08:50:59 CST 2013
On Fri, Feb 1, 2013 at 9:50 AM, <firstname.lastname@example.org> wrote:
> On Fri, Feb 1, 2013 at 9:23 AM, François Boulogne <email@example.com> wrote:
>> I need to deconvolve a signal with a filter. I had a look in the
>> documentation. The function exists but the docstring is missing and I'm
>> not satisfied of the result I got from a "simple" example.
>> filter = np.array([0,1,1,1,1,0])
>> step = np.array([0, 0, 0, 0, 0, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
>> 1, 1, 1, 1])
>> # I convolve both
>> res = convolve(step, filter, 'valid')
>> # and it returns a slope as expected
>> array([0, 0, 1, 2, 3, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4])
>> Now, I want to deconvolve.
>> deconvolve(res, filter)
>> # oops, it raises an exception
>> ValueError: BUG: filter coefficient a == 0 not supported yet
>> # well, let's try this
>> deconvolve(res, filter+1e-9)
>> (array([ 0.00000000e+00, 0.00000000e+00, 1.00000000e+09,
>> -9.99999999e+17, 9.99999999e+26, -9.99999999e+35,
>> 9.99999999e+44, -9.99999999e+53, 9.99999999e+62,
>> -9.99999999e+71, 9.99999999e+80, -9.99999999e+89]),
>> array([ 0.00000000e+00, 0.00000000e+00, 1.11022302e-16,
>> -8.27130862e-08, -4.42500000e+01, 5.46901335e+10,
>> 8.27266814e+19, 7.56858250e+28, -8.74285726e+37,
>> 9.99419626e+46, 8.27205507e+55, -8.26933326e+64,
>> 9.99999999e+89, 9.99999999e+89, 9.99999999e+89,
>> 1.00000000e+90, 9.99999999e+80]))
>> It's better but I do not recognize my signal :)
>> 1/ Am I misunderstanding or missing something?
>> 2/ How can I do it correctly?
> not supported, maybe using the numpy polynomial might work for the deconvolution
> from the docstring of lfilter, which is used by deconvolve:
> a : array_like
> The denominator coefficient vector in a 1-D sequence. If ``a``
> is not 1, then both `a` and `b` are normalized by ``a``.
> with the normalization your 1e-9 blows up the calculations.
> (it's been a long time since I tried to figure out deconvolve, and I
> always had 1 in the first position)
>> I also noted that no test exists for deconvolve() :(
>> SciPy-User mailing list
More information about the SciPy-User