Rapid fitting of exponentials
meverest
// The following function finds A, Tau, and B from an exponential decay by // (1) taking the FFT, and // (2) using the complicated function for alpha found in Rev. Sci. Instrum. 79, 023108 (2008) Function FFTFit(RDTWave, A, Tau, B) wave RDTWave Variable &A, &Tau, &B FFT/OUT=1/DEST=RDTWave_FFT RDTWave variable tm = pnt2x(RDTWave, (numpnts(RDTWave)-1)) - pnt2x(RDTWave, 0) variable omega = pnt2x(RDTWave_FFT,1)*2*pi variable alpha = ( 1/ deltax(RDTWave) )* ln(cos(2*pi/numpnts(RDTWave))+(real(RDTWave_FFT[1])/imag(RDTWave_FFT[1]))*sin(2*pi/numpnts(RDTWave) ) ) Tau = 1/alpha A = (alpha^2 + omega^2)*(Real(RDTWave_FFT[1])*deltax(RDTWave))/( alpha*(1 - exp(-alpha*tm)) ) B = (RDTWave_FFT[0]/numpnts(RDTWave) ) - ( A*(1 - exp(-alpha*tm) )/(alpha*tm) ) end
// The following function finds A, Tau, and B from an exponential decay by // the method of corrected successive integration // Rev. Sci. Instrum. 79, 023108 (2008) // Function DSIFit(RDTWave, A, Tau, B) wave RDTWave Variable &A, &Tau, &B variable Status = 0 make /O/D/N=3 CoefMat=nan make /O/D/N=3 DataMat=nan make/O/D/N=(3,3) ArrayMat=nan make /O/D/N=(numpnts(RDTWave)) twave CopyScales RDTWave twave twave = x // S_I Integrate RDTWave /D=RDT_Int variable /D S_I = sum(RDT_Int) // S_II variable /D S_II = MatrixDot(RDT_Int, RDT_Int) // S_f variable /D S_f =sum(RDTWave) // S_tI variable /D S_tI = MatrixDot(twave, RDT_Int) // S_fI variable /D S_fI = MatrixDot(RDTWave, RDT_Int) // S_ft variable /D S_ft = MatrixDot(RDTWave, twave) // S_t variable /D S_t = deltax(RDTWave)*(numpnts(RDTWave)-1)*(numpnts(RDTWave))/2 // S_tt variable /D S_tt = (deltax(RDTWave)^2)*(numpnts(RDTWave)-1)*(numpnts(RDTWave))*(2*(numpnts(RDTWave)-1)+1)/6 ArrayMat = {{numpnts(RDTWave), S_I, S_t},{S_I, S_II, S_tI},{S_t, S_tI, S_tt}} DataMat={S_f, S_fI, S_ft} MatrixOp /O CoefMat =Inv(ArrayMat) x DataMat Tau = Deltax(RDTWave)/ln(1 - coefmat[1]*Deltax(RDTWave)) // Rectangular Integration B = -CoefMat[2]/CoefMat[1] A = exp(-Deltax(RDTWave)/Tau)*CoefMat[0] - B Return Status end
Forum
Support
Gallery
Igor Pro 9
Learn More
Igor XOP Toolkit
Learn More
Igor NIDAQ Tools MX
Learn More
What are A and B? Amplitude and vertical offset? Oh- a little experimentation shows that's correct.
I tried it on this:
Make/D junk
setscale/I x 0,1,junk
junk += gnoise(.1) // pretty noisy
•printFFTFit(junk)
A= 1.37807 Tau= 0.0749016 B= -0.0283905
•printDSIFit(junk)
A= 0.984629 Tau= 0.102131 B= -0.0311092
Using Igor's built-in fit I got
y0 =-0.030051 ± 0.0128
A =0.99762 ± 0.0575
tau =0.10196 ± 0.00998
The FFT fit seems to be insensitive to where the exponential starts (that is, setscale/I x 1,2,junk gives the same answer), but DSIFit gives quite different answers.
John Weeks
WaveMetrics, Inc.
support@wavemetrics.com
April 3, 2008 at 12:30 pm - Permalink