Loading [MathJax]/extensions/tex2jax.js
Logo ROOT  
Reference Guide
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Modules Pages
SimpleFitting.C File Reference

Detailed Description

View in nbviewer Open in SWAN Create an exponential fitting The idea is to create a set of numbers x,y with the function x^3 and some noise from ROOT, fit the function to get the exponent (which must be near 3) and plot the points with noise, the known function and the fitted function

#include<TRInterface.h>
#include<TRandom.h>
TCanvas *SimpleFitting(){
TCanvas *c1 = new TCanvas("c1","Curve Fitting",700,500);
c1->SetGrid();
// draw a frame to define the range
// create the first plot (points with gaussian noise)
const Int_t n = 24;
Double_t y1[n] ;
//Generate the points along a X^3 with noise
TRandom rg;
rg.SetSeed(520);
for (Int_t i = 0; i < n; i++) {
x1[i] = rg.Uniform(0, 1);
y1[i] = TMath::Power(x1[i], 3) + rg.Gaus() * 0.06;
}
TGraph *gr1 = new TGraph(n,x1,y1);
gr1->SetMarkerStyle(8);
gr1->SetMarkerSize(1);
mg->Add(gr1);
// create the second plot
TF1 *f_known=new TF1("f_known","pow(x,3)",0,1);
TGraph *gr2 = new TGraph(f_known);
gr2->SetMarkerStyle(8);
gr2->SetMarkerSize(1);
mg->Add(gr2);
//passing data to Rfot fitting
r["x"]<<TVectorD(n, x1);
r["y"]<<TVectorD(n, y1);
//creating a R data frame
r<<"ds<-data.frame(x=x,y=y)";
//fitting x and y to X^power using Nonlinear Least Squares
r<<"m <- nls(y ~ I(x^power),data = ds, start = list(power = 1),trace = T)";
//getting the exponent
Double_t power;
r["summary(m)$coefficients[1]"]>>power;
TF1 *f_fitted=new TF1("f_fitted","pow(x,[0])",0,1);
f_fitted->SetParameter(0,power);
//plotting the fitted function
TGraph *gr3 = new TGraph(f_fitted);
gr3->SetMarkerStyle(8);
gr3->SetMarkerSize(1);
mg->Add(gr3);
mg->Draw("ap");
//displaying basic results
TPaveText *pt = new TPaveText(0.1,0.6,0.5,0.9,"brNDC");
pt->AddText("Fitting x^power ");
pt->AddText(" \"Blue\" Points with gaussian noise to be fitted");
pt->AddText(" \"Red\" Known function x^3");
TString fmsg;
fmsg.Form(" \"Green\" Fitted function with power=%.4lf",power);
pt->AddText(fmsg);
pt->Draw();
c1->Update();
return c1;
}
ROOT::R::TRInterface & r
Definition: Object.C:4
static const double x1[5]
int Int_t
Definition: RtypesCore.h:43
double Double_t
Definition: RtypesCore.h:57
@ kRed
Definition: Rtypes.h:64
@ kGreen
Definition: Rtypes.h:64
@ kBlue
Definition: Rtypes.h:64
TVectorT< Double_t > TVectorD
Definition: TVectorDfwd.h:22
ROOT R was implemented using the R Project library and the modules Rcpp and RInside
Definition: TRInterface.h:136
static TRInterface & Instance()
static method to get an TRInterface instance reference
virtual void SetFillColor(Color_t fcolor)
Set the fill area color.
Definition: TAttFill.h:37
virtual void SetMarkerColor(Color_t mcolor=1)
Set the marker color.
Definition: TAttMarker.h:38
virtual void SetMarkerStyle(Style_t mstyle=1)
Set the marker style.
Definition: TAttMarker.h:40
virtual void SetMarkerSize(Size_t msize=1)
Set the marker size.
Definition: TAttMarker.h:41
virtual void SetTextAlign(Short_t align=11)
Set the text alignment.
Definition: TAttText.h:41
The Canvas class.
Definition: TCanvas.h:27
1-Dim function class
Definition: TF1.h:210
virtual void SetParameter(Int_t param, Double_t value)
Definition: TF1.h:628
A TGraph is an object made of two arrays X and Y with npoints each.
Definition: TGraph.h:41
A TMultiGraph is a collection of TGraph (or derived) objects.
Definition: TMultiGraph.h:36
A Pave (see TPave) with text, lines or/and boxes inside.
Definition: TPaveText.h:21
virtual TText * AddText(Double_t x1, Double_t y1, const char *label)
Add a new Text line to this pavetext at given coordinates.
Definition: TPaveText.cxx:182
virtual void Draw(Option_t *option="")
Draw this pavetext with its current attributes.
Definition: TPaveText.cxx:233
This is the base class for the ROOT Random number generators.
Definition: TRandom.h:27
virtual Double_t Gaus(Double_t mean=0, Double_t sigma=1)
Samples a random number from the standard Normal (Gaussian) Distribution with the given mean and sigm...
Definition: TRandom.cxx:263
virtual void SetSeed(ULong_t seed=0)
Set the random generator seed.
Definition: TRandom.cxx:597
virtual Double_t Uniform(Double_t x1=1)
Returns a uniform deviate on the interval (0, x1).
Definition: TRandom.cxx:635
Basic string class.
Definition: TString.h:131
void Form(const char *fmt,...)
Formats a string using a printf style format descriptor.
Definition: TString.cxx:2289
TPaveText * pt
return c1
Definition: legend1.C:41
const Int_t n
Definition: legend1.C:16
static constexpr double mg
LongDouble_t Power(LongDouble_t x, LongDouble_t y)
Definition: TMath.h:725
Author
Omar Zapata

Definition in file SimpleFitting.C.