| Re: Gaussian Curve Fitting [message #7619 is a reply to message #7578] |
Mon, 02 December 1996 00:00   |
Peter Mason
Messages: 145
Registered: June 1996
|
Senior Member |
|
|
On Fri, 29 Nov 1996, Professor Michael Alport wrote:
> Can anyone give me their benefit of their experience with Gaussian
> curve fitting and (I suspect the primary culprit) curvefit.pro.
> These routines are just not robust enough. 1D Gaussfit sometimes
> suddenly fits very badly - due to subtle differences in the input data.
> Furthermore, IDL 4.01 help refers to the keyword "NTERMS" but, my
> gaussfit.pro does not accept this. Cutting and running the Gauss2Dfit
> example from the help documentation caused my 100MHz pentium to hang -
> at least after 30 minutes I rebooted.
> The Gauss2Dfit documentation refers to the need to have the data span at
> least 5-8 halfwidths - any other pearls of wisdom?
If you just want to fit a single Gaussian to an array of positive values,
you could do it directly instead of using curvefit, which implements an
iterative non-linear least-squares fitting algorithm.
e.g., Fit a Gaussian to an array y: find H, C and W which minimise
y - ( H * exp( -((x-C)/W)^2 ) ) in a least-squares sense. (X being some
wavelength array or such, with N elements.)
a=fltarr(N,3)+1.0
a(0,1)=x & a(0,2)=x*x ;the Gaussian is just a 2nd order poly in log space
b=alog(y>1.0e-6) ;..this is why Y must be positive
svdc,a,x,u,v,/column
p=svsol(u,x,v,b,/column)
W = abs(1.0/p(2))
C = 0.5*p(1)*W ;centre
H = exp(p(0)+C*C/W) ;height
W = sqrt(W) ;width
If you are fitting multiple Gaussians then you can't get a direct solution, of
course, and have to use something like Curvefit.
Now I've seen this posted in this group before, and sure have found the same:
algorithms like Curvefit's are extremely sensitive to starting values - the
initial guess of the parameters' values. There usually are several local
minima (non-optimal solutions), and the algorithm can easily become trapped in
one and never find the optimal solution if it doesn't start off fairly close
to the optimal solution.
If you use Gaussfit, it may be worth your while to check how it finds the
initial guess, given your data. (You may want to modify the code, or try
smoothing your data before fitting.)
My copy of Gaussfit also doesn't support NTERMS (IDL 4.0.1, DEC ALPHA/OSF).
Now Gauss2dfit calls Gaussfit with NTERMS=4, and naturally this crashes IDL.
(You should have got an error message in your PC's IDL console when you ran
the Gauss2dfit example.)
I'd say a bug report is in order.
Peter Mason
|
|
|
|
comp.lang.idl-pvwave archive
Members
Search
Help
Login
Home
