Re: Seeking guidance on implementing a user defined equation in ROOT

From: Rene Brun <Rene.Brun_at_cern.ch>
Date: Thu, 19 Feb 2009 09:16:06 +0100


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/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

Received on Thu Feb 19 2009 - 09:16:32 CET

This archive was generated by hypermail 2.2.0 : Tue Feb 24 2009 - 11:50:02 CET