comp.lang.idl-pvwave archive
Messages from Usenet group comp.lang.idl-pvwave, compiled by Paulo Penteado

Home » Public Forums » archive » More problems with Curvefit
Show: Today's Messages :: Show Polls :: Message Navigator
E-mail to friend 
Switch to threaded view of this topic Create a new topic Submit Reply
More problems with Curvefit [message #35619] Mon, 30 June 2003 07:35 Go to next message
heather.williams is currently offline  heather.williams
Messages: 3
Registered: July 2002
Junior Member
Good afternoon, everyone. I'm having some problems using Curvefit (in
IDL 5.4) to fit my data. I've reviewed the messages which have already
been posted on this subject, and had a good look at the MPCURVEFIT
substitute, but am no wiser.

My code looks like this:

PRO data_fit

;Define the vectors of tabulated:

F18_x = FLTARR(12) & F18_x(0) = 9.055227833
F18_x(1) = 9.908886278 & F18_x(2) = 11.86860889
F18_x(3) = 13.281685 & F18_x(4) = 16.69834393
F18_x(5) = 19.52864256 & F18_x(6) = 23.17273836
F18_x(7) = 28.51793219 & F18_x(8) = 31.23624055
F18_x(9) = 33.53401408 & F18_x(10) = 38.12262897
F18_x(11) = 39.15701348
F18_y = FLTARR(12) & F18_y(0) = 0.108598707
F18_y(1) = 0.329883541 & F18_y(2) = 0.504690343
F18_y(3) = 0.685805013 & F18_y(4) = 0.780161321
F18_y(5) = 0.87284238 & F18_y(6) = 0.890067419
F18_y(7) = 0.907523914 & F18_y(8) = 0.98011631
F18_y(9) = 0.943832957 & F18_y(10) = 0.966238284
F18_y(11) = 1

X = FLTARR(12) & X(*) = F18_x(*) - F18_x(0)
Y = FLTARR(12) & Y(*) = F18_y(*) - F18_x(0)

;Define a vector of weights:
W = 1.0
;Provide an initial guess of the function's parameters:
A = [1.0, 1.0]
;Compute the parameters a0 and a1:
yfit = CURVEFIT(X, Y, W, A, SIGMA_A, FUNCTION_NAME = 'fit_funct')
;Print the parameters, which are returned in A:
PRINT, A

END

PRO fit_funct, X, A, F, PDER

F = (1.0 - EXP(-A[0] * X)) + (1.0 - EXP(-A[1]*X))
; PDER's column dimension is equal to the number of
; elements in xi and its row dimension is equal to
; the number of parameters in the function F:
pder = FLTARR(N_ELEMENTS(X), 2)
; Compute the partial derivatives with respect to
; a0 and place in the first row of PDER:
pder[*, 0] = A[0] * EXP(-A[0] * X)
pder[*, 1] = A[1] * EXP(-A[1] * X)

END


Which looks alright, if not particularly elegant, to me. However, when
I run it, I get this error message (which relates to the line
beginning y_fit =) :

% Operands of matrix multiply have incompatible dimensions: <FLOAT
Array[2, 12]>, <FLOAT Array[1, 2]>.
% Error occurred at: CURVEFIT 269
O:\Rsi\Idl54\lib\curvefit.pro
% DATA_FIT 21 H:\PhD
IDL\Progs\data_fit.pro
% $MAIN$

How do I avoid this error and get the fit to work?

Thanks for your help,
Heather Williams

PhD Student, Manchester PET Centre
Manchester, UK
Re: More problems with Curvefit [message #35689 is a reply to message #35619] Tue, 01 July 2003 10:12 Go to previous message
Heinz Stege is currently offline  Heinz Stege
Messages: 189
Registered: January 2003
Senior Member
On 30 Jun 2003 07:35:26 -0700, heather.williams@physics.cr.man.ac.uk
(Heather Williams) wrote:

> Y = FLTARR(12) & Y(*) = F18_y(*) - F18_x(0)

Hmm, I didn't go into the details. But this line looks a little bit
strange. Do you eventually want the following operation:

Y = FLTARR(12) & Y(*) = F18_y(*) - F18_y(0)

Heinz
Re: More problems with Curvefit [message #35706 is a reply to message #35619] Tue, 01 July 2003 04:54 Go to previous message
heather.williams is currently offline  heather.williams
Messages: 3
Registered: July 2002
Junior Member
Thanks to everyone for your help - it's now working happily, but even
with the recommended changes, this function doesn't seem to match the
data! I shall go in search of a more appropriate model, but any
suggestions from the previous respondents would be welcome.

- Heather
Re: More problems with Curvefit [message #35711 is a reply to message #35619] Mon, 30 June 2003 11:43 Go to previous message
Craig Markwardt is currently offline  Craig Markwardt
Messages: 1869
Registered: November 1996
Senior Member
Paul van Delst <paul.vandelst@noaa.gov> writes:
> Try doing
> w = make_array(12, value = 1.0 )
> or
> w = replicate( 1.0, 12 )
> to actually give you a vector for w. The IDL documentation is a wee bit misleading here as
> the actual words say "For no weighting, set Weightsi = 1.0" where the "i" suffix is an
> indicator that Weights is an array but it's not entirely clear (to me at least.)

As a side note, MPCURVEFIT, and the other "MP" fitting programs can
accept either a scalar or a vector for the uncertainties or weights.
Either way it does the "right thing." And of course the MP fitting
programs do not necessarily require the user to compute function
derivatives.

To Heather Williams, the original poster, you should be aware that
your fitting function is degenerate, since it contains a linear
combination of the *same* basis function. This leads to a singular
normal matrix, and takes a canned routine like CURVEFIT on a trip to
la-la land. This effect is further exacerbated by the choice of
initial conditions (both coefficients equal). MPCURVEFIT is more
robust in this sense, since it will do a singular value decomposition
in the face of a singular matrix.


Happy fitting!
Craig

--
------------------------------------------------------------ --------------
Craig B. Markwardt, Ph.D. EMAIL: craigmnet@cow.physics.wisc.edu
Astrophysics, IDL, Finance, Derivatives | Remove "net" for better response
------------------------------------------------------------ --------------
Re: More problems with Curvefit [message #35715 is a reply to message #35619] Mon, 30 June 2003 10:17 Go to previous message
the_cacc is currently offline  the_cacc
Messages: 104
Registered: October 2001
Senior Member
1) You should put 'fit_funct' at the top of the file so it compiles
first.
2) Change "W = 1.0" to "W = replicate(1.0,12)"
3) pder should be the partial deriv wrt the parameters A - ie.
pder[*,0] = X * EXP(-A[0]*X) and pder[*,1] = X * EXP(-A[1]*X).
4) Looking at your data (PLOT,x,y), it seems that the choice of
fitting function is unlikely to fit well... ever. You need to allow
for a constant offset - the value 9 looks "right". If you can be sure
it is 9, then simply change y to y+9 in the call to curvefit.
5) To see the fit add this line at the end of your code: PLOT,x,y &
fit_funct,x,a,f & OPLOT,x,f-9.

Hope this gives you somewhere to start. I recommend changing your fit
function to deal with the offset issue. More generally, I would also
recommend looking at AMOEBA for fitting a small number of parameters -
very stable, doesn't need derivatives and not fussy about
discontinuities etc.

Ciao.
Re: More problems with Curvefit [message #35716 is a reply to message #35619] Mon, 30 June 2003 09:19 Go to previous message
Paul Van Delst[1] is currently offline  Paul Van Delst[1]
Messages: 1157
Registered: April 2002
Senior Member
Heather Williams wrote:
>
> Good afternoon, everyone. I'm having some problems using Curvefit (in
> IDL 5.4) to fit my data. I've reviewed the messages which have already
> been posted on this subject, and had a good look at the MPCURVEFIT
> substitute, but am no wiser.
>
> My code looks like this:
>
> PRO data_fit
>
> ;Define the vectors of tabulated:
>
> F18_x = FLTARR(12) & F18_x(0) = 9.055227833
> F18_x(1) = 9.908886278 & F18_x(2) = 11.86860889
> F18_x(3) = 13.281685 & F18_x(4) = 16.69834393
> F18_x(5) = 19.52864256 & F18_x(6) = 23.17273836
> F18_x(7) = 28.51793219 & F18_x(8) = 31.23624055
> F18_x(9) = 33.53401408 & F18_x(10) = 38.12262897
> F18_x(11) = 39.15701348
> F18_y = FLTARR(12) & F18_y(0) = 0.108598707
> F18_y(1) = 0.329883541 & F18_y(2) = 0.504690343
> F18_y(3) = 0.685805013 & F18_y(4) = 0.780161321
> F18_y(5) = 0.87284238 & F18_y(6) = 0.890067419
> F18_y(7) = 0.907523914 & F18_y(8) = 0.98011631
> F18_y(9) = 0.943832957 & F18_y(10) = 0.966238284
> F18_y(11) = 1
>
> X = FLTARR(12) & X(*) = F18_x(*) - F18_x(0)
> Y = FLTARR(12) & Y(*) = F18_y(*) - F18_x(0)
>
> ;Define a vector of weights:
> W = 1.0

Try doing

w = make_array(12, value = 1.0 )

or

w = replicate( 1.0, 12 )

to actually give you a vector for w. The IDL documentation is a wee bit misleading here as
the actual words say "For no weighting, set Weightsi = 1.0" where the "i" suffix is an
indicator that Weights is an array but it's not entirely clear (to me at least.)


This made your code work for me.....but I got the following:

% CURVEFIT: Failed to converge- CHISQ increasing without bound.
1.00000 1.00000
% Program caused arithmetic error: Floating overflow

Changing the line

yfit = CURVEFIT(X, Y, W, A, SIGMA_A, FUNCTION_NAME = 'fit_funct')

to use Craig's MPCURVEFIT

yfit = MPCURVEFIT(X, Y, W, A, SIGMA_A, FUNCTION_NAME = 'fit_funct')

gave me the following:

IDL> data_fit
% Compiled module: MPCURVEFIT.
% Compiled module: MPFIT.
Iter 1 CHI-SQUARE = 1215.4340 DOF = 10
P(0) = 1.00000
P(1) = 1.00000
Iter 2 CHI-SQUARE = 1212.4130 DOF = 10
P(0) = 0.887302
P(1) = 0.887302
Iter 3 CHI-SQUARE = 1208.7012 DOF = 10
P(0) = 0.776462
P(1) = 0.776462
Iter 4 CHI-SQUARE = 1203.9254 DOF = 10
P(0) = 0.665368
P(1) = 0.665368
Iter 5 CHI-SQUARE = 1189.1635 DOF = 10
P(0) = 0.448375
P(1) = 0.448375
Iter 6 CHI-SQUARE = 945.99725 DOF = 10
P(0) = 0.0271755
P(1) = 0.0271756
Iter 7 CHI-SQUARE = 440.82562 DOF = 10
P(0) = -0.0570307
P(1) = -0.0570309
-0.0570307 -0.0570309
% Program caused arithmetic error: Floating underflow
% Program caused arithmetic error: Floating overflow
% Program caused arithmetic error: Floating illegal operand

This output just *has* to be more meaningful wrt diagnostics than the curvefit output.

paulv

--
Paul van Delst
CIMSS @ NOAA/NCEP/EMC
Ph: (301)763-8000 x7748
Fax:(301)763-8545
  Switch to threaded view of this topic Create a new topic Submit Reply
Previous Topic: Re: Format Code Question
Next Topic: IDL from crontab going defunct?

-=] Back to Top [=-
[ Syndicate this forum (XML) ] [ RSS ] [ PDF ]

Current Time: Wed Oct 08 13:04:57 PDT 2025

Total time taken to generate the page: 0.00457 seconds