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