Re: Seeking further guidance on implementing a user defined equation in ROOT (More)

From: Rene Brun <Rene.Brun_at_cern.ch>
Date: Mon, 23 Feb 2009 09:17:29 +0100


Your input value for alpha_k is 0.221243. The best fit is obtained with alpha_k= 0.225598
Up to you ::)

Rene Brun

Russell Leslie wrote:
> Dear Rene,
>
> Thank you for the response and the advice that not fixing alpha_k
> would give a better fit.
>
> Unfortunately, in my current approach to this problem alpha_k is the
> defining physical parameter of the situation. It is not clear that
> the model would remain physically realistic if the value of alpha_k
> was to change substantially.
>
> I will explore this issue further.
>
> Regards
>
> Russell Leslie
> Nuclear Physics
> R.S.Phys.S.E.
> ANU
>
>
>
> ----- Original Message -----
> From: Rene Brun <Rene.Brun_at_cern.ch>
> Date: Monday, February 23, 2009 6:59 pm
> Subject: Re: [ROOT] Seeking further guidance on implementing a user
> defined equation in ROOT (More)
> To: Russell Leslie <u4504546_at_anu.edu.au>
> Cc: roottalk_at_lxbuild091.cern.ch
>
> > ROOT considers that the fit fails (despite computing good
> > parameters)
> > because your data points do not have correct errors associated
> > with each
> > bin. By the way, when you do not fix alpha_k you get a better fit.
> >
> > Rene Brun
> >
> > Russell Leslie wrote:
> > > Dear Rene,
> > >
> > > Thank you very much for the scripted fitting function and the
> > advice
> > > to check the equation (in the definition of function "fint"
> > the
> > > parameter "Lambdas" had been replaced with parameter
> > "lambda_k" ).
> > >
> > > Your advice to fix one or more parameters was well
> > taken. The
> > > parameters alpha_k and econst are fixed for any individual
> > case
> > > (varying only between cases). I also set upper and lower
> > bounds for
> > > the non-fixed parameters (as none of the parameters can take
> > negative
> > > values).
> > >
> > > One thing that puzzles me is that the end result of the fit
> > looks
> > > pretty well perfect, the values of the various parameters are
> > > physically plausible and the curves appear to match well, but
> > root
> > > gave the following message:
> > >
> > > Warning in <TH1F::Fit>: Abnormal termination of minimization.
> > > FCN=0.0598415 FROM MIGRAD
> > STATUS=FAILED 458
> > CALLS
> > > 459 TOTAL
> > > EDM=0.000252546 STRATEGY= 1 ERR MATRIX
> > > NOT POS-DEF
> > > EXT
> > PARAMETER APPROXIMATE STEP FIRST
> > > NO. NAME
> > VALUE ERROR SIZE DERIVATIVE
> > > 1
> > alpha_k 2.21243e-
> > 01 fixed
> > > 2
> > gamma 4.86939e-
> > 02 2.20292e-01 -0.00000e+00 -1.76568e-02
> > > 3 lambda_k
> > 2.54645e-01 4.59892e-02 -0.00000e+00 -
> > 2.42907e-01
> > > 4
> > tau_c 5.79966e-
> > 02 3.65613e-01 0.00000e+00 -
> > 8.21423e-02
> > > 5
> > econst 1.00000e+00 fixed
> > > 6
> > Lambdas 7.83294e-02
> > 1.37364e-01 -0.00000e+00 3.41404e-02
> > >
> > > Should I worry that root considers the fit a failure - or
> > should I
> > > just accept the results of the fit.
> > >
> > > Regards
> > >
> > > Russell Leslie
> > > Nuclear Physics
> > > R.S.Phys.S.E.
> > > ANU
> > >
> > >
> > > ----- Original Message -----
> > > From: Rene Brun <Rene.Brun_at_cern.ch>
> > > Date: Saturday, February 21, 2009 12:51 am
> > > Subject: Re: [ROOT] Seeking further uidance on implementing a
> > user
> > > defined equation in ROOT (More)
> > > To: Russell Leslie <u4504546_at_anu.edu.au>
> > > Cc: roottalk_at_lxbuild091.cern.ch
> > >
> > > > In my previouis mail I forgot to mention 2 problems
> > > > -you should provide a good error estimate for
> > the bin
> > > > contents of your
> > > > histogram. The default being to use the sqrt(bin content) does
> > > > not make
> > > > much sense in your case. I used the option "w" that means use
> > > > same
> > > > weight=1 for all bins.
> > > >
> > > > -you should fix a few parameters. In my previous mail, all
> > > > parameters
> > > > were free. in the attachment
> > > > you will find an example where one of your parameters is fixed,
> > > > you can
> > > > add more.
> > > >
> > > > Rene Brun
> > > >
> > > > Russell Leslie wrote:
> > > > > Dear Rene,
> > > > >
> > > > > I seem to be missing some vital step in using a user defined
> > > > equation
> > > > > to fit to a histogram.
> > > > >
> > > > > I am trying to use the function gks from the script gks.C
> > > > (that you
> > > > > produced for me) as a fitting function for the Monte Carlo
> > > > data that
> > > > > we produce from a FORTRAN program. I have attached a
> > > > copy of a root
> > > > > file containing a typical set of such data (file MCGcalc.root).
> > > > >
> > > > > I have compiled your script (using the command ".x gks.C+")
> > > > and
> > > > > produced a shared object gks_C.so.
> > > > >
> > > > > I have loaded gks_C.so (using the command ".L gks_C.so") and
> > > > then
> > > > > tried to fit using various permutations of commands such as
> > > > >
> > > > > G2values.Fit("gks") but everything I have tried so far gives
> > > > rise to
> > > > > messages similar to "Unknown function: gks" (with the
> > > > occasional seg
> > > > > fault).
> > > > >
> > > > > I have been similarly unsuccessful using the FitPanel GUI.
> > > > >
> > > > > Any guidance you can give on what I might be doing wrong would
> > > > be
> > > > > greatly appreciated.
> > > > >
> > > > > Regards
> > > > >
> > > > > Russell Leslie
> > > > > Nuclear Physics
> > > > > R.S.Phys.S.E.
> > > > > ANU
> > > > >
> > > > >
> > > > > ----- Original Message -----
> > > > > From: Russell Leslie <russell.leslie_at_gmail.com>
> > > > > Date: Friday, February 20, 2009 8:07 pm
> > > > > Subject: Fwd: [ROOT] Seeking guidance on implementing a user
> > > > defined
> > > > > equation in ROOT
> > > > > To: u4504546_at_anu.edu.au
> > > > >
> > > > > >
> > > > > >
> > > > > > ---------- Forwarded message ----------
> > > > > > From: *Rene Brun* <Rene.Brun_at_cern.ch>
> > > > > > Date: Thu, Feb 19, 2009 at 7:16 PM
> > > > > > Subject: Re: [ROOT] Seeking guidance on implementing a user
> > > > defined
> > > > > equation in ROOT
> > > > > > To: Russell Leslie <u4504546_at_anu.edu.au
> > > > > <javascript:main.compose('new', 't=u4504546_at_anu.edu.au')>>
> > > > > > Cc: roottalk_at_lxbuild091.cern.ch
> > > > <javascript:main.compose('new',
> > > > > 't=roottalk_at_lxbuild091.cern.ch')>
> > > > > >
> > > > > >
> > > > > > Try the function gks.C in attachment
> > > > > > To run do:
> > > > > > root > .x gks.C+
> > > > > >
> > > > > > Rene Brun
> > > > > >
> > > > > > Russell Leslie wrote:
> > > > >
> > > > > > Dear all,
> > > > > > NOTE - message resent because
> > > > I sent this originally from
> > > > > another account and it did not
> > turn up
> > > > on the roottalk list.
> > > > > > I am still feeling my way
> > > > forward on ROOT, but I have found it
> > > > > to be a very useful tool-set.
> > > > > > One problem I am currently working
> > > > on is to try to implement the
> > > > > following equations in ROOT (LaTex
> > > > form of the equations below).
> > > > > >
> > > > > > $$
> > > > > > G_{k}(t) = e^{-\Lambda^{*}
> > > > t}G_{k}^{(fluct.)}(t)+ \int_{0}^{t}
> > > > > \Lambda^{*} G_{k}^{(fluct.)}(t')\cdotp
> > > > e^{-\Lambda^{*}
> > > > > t'}G_{k}^{(stat.)} (t-t')dt'
> > > > > > $$
> > > > > > Where Gk(stat) is give by the
> > > > following equation
> > > > > > $$
> > > > > > G_{k}^{(stat.)} (t) = \left\langle
> > > > \alpha_{k} \right\rangle + (
> > > > > 1 - \left\langle
> > > > \alpha_{k}\right\rangle )e^{- \Gamma t}
> > > > > > $$
> > > > > > and Gk(fluct) is given by the
> > > > this equation
> > > > > > $$
> > > > > > G_{k}^{(fluct.)}(t) = ( 1 -
> > > > \lambda_{k}\tau_{c} )e^{ - \lambda
> > > > > _{k} t }
> > > > > > $$
> > > > > > The equations are from the
> > > > manual for the Coulex code GOSIA and
> > > > > I want to use them as the fitting
> > > > function for the output of a
> > > > > Monte Carlo simulation that we use.
> > > > > > I have had no problems
> > > > producing the curves for Gk(stat) and
> > > > > Gk(fluct) - but I have not been able
> > > > to implement the integral
> > > > > term in the Gk equation (the first
> > > > equation above). The code I
> > > > > have so far is given below
> > > > > > {
> > > > > > //Calculate Gks
> > > > > > //parameter[0]
> > = alpha_k
> > > > > > //parameter[1]
> > = gamma
> > > > > > TF1 *Gks = new
> > > > TF1("Gks","[0]+(1-[0])*exp(-[1]*x)",0,15);
> > > > > > Gks-
> > > > >SetParameters(0.221243,0.4);>
> > > > > //Calculate Gkf
> > > > > > //parameter[0] =
> > > > lambda_k> >
> > > > //parameter[1] = tau_c
> > > > > > TF1 *Gkf = new
> > > > TF1("Gkf","(1-[0]*[1])*exp(-[0]*x)",0,15);
> > > > > > Gkf-
> > > > >SetParameters(0.172,0.5);>
> > > > > //Calculate expLamba*t
> > > > > > //parameter[0] =
> > > > 1 //parameter[1] = Lambda*
> > > > > > TF1 *eLst = new
> > > > TF1("eLst","[0]*exp(-[1]*x)",0,15);
> > > > > > eLst-
> > > > >SetParameters(1,0.0345);>
> > > > > //Combine Gks and Gkf to for
> > Gkcalc> > > > TF1
> > *Gkcalc = new
> > > > TF1("Gkcalc","Gkf*eLst*[0]*[1]",0,15); //
> > > > > This is the bit I can't get!!!!!!!!
> > > > > >
> > > > //Create a canvas to draw Gkcalc
> > > > > > TCanvas
> > *Gkcalc1 = new
> > > > >
> > > >
> > TCanvas("Gkcalc1","Gkcalc1",10,10,1400,1200);> > //Draw and
> > > format Gkcalc
> > > > > > Gkcalc->Draw();
> > > > > > Gkcalc-
> > > > >SetLineWidth(2.0);> >
> > > > Gkcalc->SetLineColor(kBlue);
> > > > > > Gkcalc-
> > >GetXaxis()-
> > > > >SetTitle("time");> >
> > > > Gkcalc->GetYaxis()->SetTitle("G_{k}(t)");
> > > > > > Gkcalc-
> > > > >SetTitle("G_{k}(t) = e^{-#Lambda^{*}
> > > > > t}G_{k}^{(fluct.)}(t)+ #int_{0}^{t}
> > > > #Lambda^{*}>
> > G_{k}^{(fluct.)}(t').e^{-
> > > > #Lambda^{*} t'}G_{k}^{(stat.)} (t-t')dt' ");
> > > > > > Gkcalc-
> > >GetYaxis()-
> > > > >SetLabelFont(22);> >
> > > > Gkcalc->GetYaxis()->SetTitleFont(22);
> > > > > > Gkcalc-
> > >GetXaxis()-
> > > > >SetLabelFont(22);> >
> > > > Gkcalc->GetXaxis()->SetTitleFont(22);
> > > > > > //Add other lines
> > > > > > Gks->Draw("same");
> > > > > > Gks-
> > > > >SetLineWidth(2.0);> >
> > > > Gks->SetLineColor(kRed);
> > > > > > Gkf->Draw("same");
> > > > > > Gkf-
> > > > >SetLineWidth(2.0);> >
> > > > Gkf->SetLineColor(kGreen);
> > > > > > eLst->Draw("same");
> > > > > > eLst-
> > > > >SetLineWidth(2.0);> >
> > > > eLst->SetLineColor(kMagenta);
> > > > > > TLegend
> > *legend =
> > > > new TLegend(.75,.80,.95,.95);
> > > > > > legend-
> > > > >AddEntry(Gkcalc,"G_{k}(t)");>
> > > > > legend->AddEntry(Gks,"G_{k}^{(stat)}(t)");
> > > > > > legend-
> > > > >AddEntry(Gkf,"G_{k}^{(fluct)}(t)");>
> > > > > legend->AddEntry(Gks,"e^{-#Lambda^{*} t}");
> > > > > > legend->Draw();
> > > > > > }
> > > > > > My root details as as follows ROOT
> > > > 5.18/00b> (branches/v5-18-00-
> > > > patches_at_22563 <mailto:branches
> > > > > <javascript:main.compose('new',
> > > > > 't=branches')>/v5-18-00-
> > > > patches_at_22563>, Oct 19 2008, 22:04:00 on
> > > > > linuxx8664gcc4.3)
> > > > >
> > > > > > Any help or guidance you can
> > > > give on this would be greatly
> > > > > appreciated.
> > > > > > Russell Leslie
> > > > > > Nuclear Physics
> > > > > > R.S.Phys.S.E.
> > > > > > ANU
> > > > >
> > > > > >
> > > > > >
> > > > >
> > > > >
> > > > > Russell Leslie
> > > > > Nuclear Physics
> > > > > R.S.Phys.S.E.
> > > > > ANU
> > > >
> > >
> > > Russell Leslie
> > > Nuclear Physics
> > > R.S.Phys.S.E.
> > > ANU
> >
>
> Russell Leslie
> Nuclear Physics
> R.S.Phys.S.E.
> ANU
Received on Mon Feb 23 2009 - 09:17:53 CET

This archive was generated by hypermail 2.2.0 : Mon Feb 23 2009 - 17:50:01 CET