auto r = ROOT::R::TRInterface::Instance();
A simple example using random points with Gaussian distribution. We use the interactive ROOT visualisation, JSROOT.
%jsroot on
auto c = new TCanvas("c");
TH1D h1("h1", "h1", 256, -4, 4);
h1.FillRandom("gaus");
ROOT::Math::MinimizerOptions::SetDefaultMinimizer("RMinimizer","BFGS");
h1.Fit("gaus");
c->Draw();
Calling R with command result <- optim( initialparams, minfunction,method='BFGS',control = list(ndeps=stepsizes,maxit=1000,trace=0,abstol=1.000000e-02),hessian=TRUE) Value at minimum =247.533 **************************************** Minimizer is RMinimizer / BFGS Chi2 = 247.533 NDf = 197 NCalls = 110 Constant = 68.7203 +/- 4.08237 Mean = 0.00966072 +/- 0.0118658 Sigma = 0.943347 +/- 0.03857 (limited)
The next example creates an exponential fit. The idea is to create a set of numbers x,y plus noise using ROOT, pass them to R, fit the data using nls(Nonlinear Least Squares), obtain fitted coefficient and then plot data, known function and fitted function with ROOT
// draw a frame to define the range
TMultiGraph mg;
const Int_t n = 24;
Double_t x[n] ;
Double_t y[n] ;
Generate random points from $X^3$ with noise
TRandom rg;
rg.SetSeed(520);
for (Int_t i = 0; i < n; i++) {
x[i] = rg.Uniform(0, 1);
y[i] = TMath::Power(x[i], 3) + rg.Gaus() * 0.06;
}
%%cpp -d
void SetPointStyle(TGraph& g, Int_t color)
{
g.SetMarkerColor(color);
g.SetFillColor(kWhite);
g.SetMarkerStyle(8);
g.SetMarkerSize(1);
}
TGraph gr1(n,x,y);
SetPointStyle(gr1, kBlue);
mg.Add(&gr1);
TF1 f_known("f_known","pow(x,3)",0,1);
TGraph gr2(&f_known);
SetPointStyle(gr2, kRed);
mg.Add(&gr2);
passing data to R for fitting
TVectorD vx(n, x);
TVectorD vy(n, y);
r["x"]=vx;
r["y"]=vy;
r<<"ds<-data.frame(x=x,y=y)";
r << "m <- nls(y ~ I(x^power),data = ds, start = list(power = 1),trace = T)";
1.506598 : 1 0.3263132 : 1.885851 0.09498773 : 2.719307 0.0764783 : 3.086762 0.07622443 : 3.136312 0.07622369 : 3.139027 0.07622369 : 3.139151
Double_t power;
r["summary(m)$coefficients[1]"] >> power;
Plot fitted function
TCanvas c1("c1", "Curve Fitting", 700, 500);
c1.SetGrid();
TF1 f_fitted("f_fitted","pow(x,[0])",0,1);
f_fitted.SetParameter(0,power);
TGraph gr3(&f_fitted);
SetPointStyle(gr3, kGreen);
mg.Add(&gr3);
mg.Draw("ap");
TLegend leg(0.1,0.6,0.5,0.9,"brNDC");
leg.SetHeader("Fitting x^power");
leg.AddEntry(&gr1, "Points with gaussian noise to be fitted", "P");
leg.AddEntry(&gr2, "Known function x^3", "P");
TString fmsg;
fmsg.Form(" \"Green\" Fitted function with power=%.4lf",power);
leg.AddEntry(&gr3, fmsg, "P");
leg.Draw();
c1.Draw();