[SciPy-User] fitting with convolution?

David Baddeley david_baddeley@yahoo.com...
Fri Jul 19 08:57:32 CDT 2013

Hi Petro,

when I run the code I get very similar, although not identical curves, which are roughly what you'd expect. Your exponential kernel that you're convolving with is not band limited, but FFTs (and hence FFT based convolution) are band-limited to whatever your effective nyquist frequency is (strictly speaking out of band frequencies might be aliased into the pass band - but in your case you can just think of the high frequency components as being discarded). This probably applies generally to any discrete way of doing convolution, when compared to an analytical solution. Put simply, you can never sample fine enough  to reproduce all the detail of the initial spikey bit of your exponential (if you increase your sampling frequency you'll get a better match, but you'll never quite get there). If you are referring to the factor of 14.5, it is a result of fft normalisation, and should more accurately be sqrt(N)/pi.

hope this helps,

 From: Petro <x.piter@gmail.com>
To: scipy-user@scipy.org 
Sent: Friday, 12 July 2013 2:13 PM
Subject: [SciPy-User] fitting with convolution?

Hi all,
I try to fir a time-resolved dataset with multiple exponents convoluted 
with a Gaussian instrument response function (IRF).
I had a look how it is done in Origin
There fft_fft_convolution calculates the circular convolution of an 
exponent with IRF.
I have found a similar function for python here:

This convolution also can be calculated analytically as, for example, in 
this package:

def convolutedexp(tau,mu,fwhm,x):
     d = (fwhm/(2*sqrt(2*log(2))))
     return 0.5*exp(-x/tau)*exp((mu+(d**2.)/(2.*tau))/tau)* 

def gaussian(mu,fwhm,x):
    d = (fwhm/(2.*sqrt(2.*log(2.))))
    return exp(-((x-mu)**2.)/(2.*d**2.))

My problem is if I compare analytical and circular convolution they do 
not match:
import numpy
from scipy.special import erf
def cconv(a, b):
     Computes the circular convolution of the (real-valued) vectors a and b.
     return fft.ifft(fft.fft(a) * fft.fft(b)).real

def convolutedexp(tau,mu,fwhm,x):
     d = (fwhm/(2*sqrt(2*log(2))))

def gaussian(mu,fwhm,x):
    d = (fwhm/(2.*sqrt(2.*log(2.))))
    return exp(-((x-mu)**2.)/(2.*d**2.))

t = array(linspace(-10.0,1000.0,2040.0))[:-1]
mu = 0
fwhm = 4.0
tau = 20.0
uf = gaussian(mu,fwhm,t)

vf = exp(-t/tau)
uvf1 = cconv(uf,vf)
uvf2 = convolutedexp(tau,mu,fwhm,t)


My feeling is that I miss something about convolution?
Can anybody give me a hint?

SciPy-User mailing list
-------------- next part --------------
An HTML attachment was scrubbed...
URL: http://mail.scipy.org/pipermail/scipy-user/attachments/20130719/33048e09/attachment.html 

More information about the SciPy-User mailing list