{
"cells": [
{
"cell_type": "markdown",
"id": "e616f730",
"metadata": {},
"source": [
"# myfit\n",
"Get in memory an histogram from a root file and fit a user defined function.\n",
"Note that a user defined function must always be defined\n",
"as in this example:\n",
" - first parameter: array of variables (in this example only 1-dimension)\n",
" - second parameter: array of parameters\n",
"Note also that in case of user defined functions, one must set\n",
"an initial value for each parameter.\n",
"\n",
"\n",
"\n",
"\n",
"**Author:** Rene Brun \n",
"This notebook tutorial was automatically generated with ROOTBOOK-izer from the macro found in the ROOT repository on Tuesday, May 19, 2026 at 08:25 PM."
]
},
{
"cell_type": "markdown",
"id": "afa2bceb",
"metadata": {},
"source": [
" Definition of a helper function: "
]
},
{
"cell_type": "code",
"execution_count": 1,
"id": "91d26b76",
"metadata": {
"collapsed": false,
"execution": {
"iopub.execute_input": "2026-05-19T20:25:19.472532Z",
"iopub.status.busy": "2026-05-19T20:25:19.472419Z",
"iopub.status.idle": "2026-05-19T20:25:19.477352Z",
"shell.execute_reply": "2026-05-19T20:25:19.476797Z"
}
},
"outputs": [],
"source": [
"%%cpp -d\n",
"\n",
"#include \n",
"#include \n",
"#include \n",
"#include \n",
"#include \n",
"#include \n",
"\n",
"#include \n",
"\n",
"double fitf(double *x, double *par)\n",
"{\n",
" double arg = 0;\n",
" if (par[2] != 0) arg = (x[0] - par[1])/par[2];\n",
"\n",
" double fitval = par[0]*std::exp(-0.5*arg*arg);\n",
" return fitval;\n",
"}"
]
},
{
"cell_type": "code",
"execution_count": 2,
"id": "65f73eda",
"metadata": {
"collapsed": false,
"execution": {
"iopub.execute_input": "2026-05-19T20:25:19.478761Z",
"iopub.status.busy": "2026-05-19T20:25:19.478645Z",
"iopub.status.idle": "2026-05-19T20:25:20.031082Z",
"shell.execute_reply": "2026-05-19T20:25:20.030643Z"
}
},
"outputs": [],
"source": [
"TString dir = gROOT->GetTutorialDir();\n",
"dir.Append(\"/hsimple.C\");\n",
"dir.ReplaceAll(\"/./\",\"/\");\n",
"if (!gInterpreter->IsLoaded(dir.Data())) gInterpreter->LoadMacro(dir.Data());\n",
"TFile *hsimpleFile = (TFile*)gROOT->ProcessLineFast(\"hsimple(1)\");\n",
"if (!hsimpleFile) return;\n",
"\n",
"TCanvas *c1 = new TCanvas(\"c1\",\"the fit canvas\",500,400);\n",
"\n",
"TH1F *hpx = (TH1F*)hsimpleFile->Get(\"hpx\");"
]
},
{
"cell_type": "markdown",
"id": "b11be27d",
"metadata": {},
"source": [
"Creates a Root function based on function fitf above"
]
},
{
"cell_type": "code",
"execution_count": 3,
"id": "22ca5831",
"metadata": {
"collapsed": false,
"execution": {
"iopub.execute_input": "2026-05-19T20:25:20.036303Z",
"iopub.status.busy": "2026-05-19T20:25:20.036173Z",
"iopub.status.idle": "2026-05-19T20:25:20.240613Z",
"shell.execute_reply": "2026-05-19T20:25:20.240207Z"
}
},
"outputs": [],
"source": [
"TF1 *func = new TF1(\"fitf\",fitf,-2,2,3);"
]
},
{
"cell_type": "markdown",
"id": "f387c814",
"metadata": {},
"source": [
"Sets initial values and parameter names"
]
},
{
"cell_type": "code",
"execution_count": 4,
"id": "f2e85877",
"metadata": {
"collapsed": false,
"execution": {
"iopub.execute_input": "2026-05-19T20:25:20.251127Z",
"iopub.status.busy": "2026-05-19T20:25:20.250999Z",
"iopub.status.idle": "2026-05-19T20:25:20.455668Z",
"shell.execute_reply": "2026-05-19T20:25:20.455038Z"
}
},
"outputs": [],
"source": [
"func->SetParameters(100,0,1);\n",
"func->SetParNames(\"Constant\",\"Mean_value\",\"Sigma\");"
]
},
{
"cell_type": "markdown",
"id": "bf80b0fb",
"metadata": {},
"source": [
"Fit histogram in range defined by function"
]
},
{
"cell_type": "code",
"execution_count": 5,
"id": "b520f739",
"metadata": {
"collapsed": false,
"execution": {
"iopub.execute_input": "2026-05-19T20:25:20.457334Z",
"iopub.status.busy": "2026-05-19T20:25:20.457199Z",
"iopub.status.idle": "2026-05-19T20:25:20.662237Z",
"shell.execute_reply": "2026-05-19T20:25:20.661709Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"****************************************\n",
"Minimizer is Minuit2 / Migrad\n",
"Chi2 = 36.7428\n",
"NDf = 47\n",
"Edm = 2.03167e-06\n",
"NCalls = 101\n",
"Constant = 797.969 +/- 6.79742 \n",
"Mean_value = -7.42918e-05 +/- 0.00734861 \n",
"Sigma = 0.998754 +/- 0.0071337 \n"
]
}
],
"source": [
"hpx->Fit(func,\"r\");"
]
},
{
"cell_type": "markdown",
"id": "77ed2b7b",
"metadata": {},
"source": [
"Gets integral of function between fit limits"
]
},
{
"cell_type": "code",
"execution_count": 6,
"id": "36981b68",
"metadata": {
"collapsed": false,
"execution": {
"iopub.execute_input": "2026-05-19T20:25:20.676222Z",
"iopub.status.busy": "2026-05-19T20:25:20.676082Z",
"iopub.status.idle": "2026-05-19T20:25:20.880900Z",
"shell.execute_reply": "2026-05-19T20:25:20.880342Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Integral of function = 1907.36\n"
]
}
],
"source": [
"printf(\"Integral of function = %g\\n\",func->Integral(-2,2));"
]
},
{
"cell_type": "markdown",
"id": "a58a667b",
"metadata": {},
"source": [
"Draw all canvases "
]
},
{
"cell_type": "code",
"execution_count": 7,
"id": "c21db95a",
"metadata": {
"collapsed": false,
"execution": {
"iopub.execute_input": "2026-05-19T20:25:20.882219Z",
"iopub.status.busy": "2026-05-19T20:25:20.882103Z",
"iopub.status.idle": "2026-05-19T20:25:21.108829Z",
"shell.execute_reply": "2026-05-19T20:25:21.108191Z"
}
},
"outputs": [
{
"data": {
"text/html": [
"\n",
"\n",
"\n",
"
\n",
"\n",
"\n",
"\n"
],
"text/plain": [
""
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"%jsroot on\n",
"gROOT->GetListOfCanvases()->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": 5
}