{ "cells": [ { "cell_type": "markdown", "id": "1c795c5f", "metadata": {}, "source": [ "# rf616_morphing\n", "Use Morphing in RooFit.\n", "\n", "This tutorial shows how to use template morphing inside RooFit. As input we have several\n", "Gaussian distributions. The output is one gaussian, with a specific mean value.\n", "Since likelihoods are often used within the framework of morphing, we provide a\n", "way to estimate the negative log likelihood (nll).\n", "\n", "Based on example of Kyle Cranmer https://gist.github.com/cranmer/b67830e46d53d5f7cf2d.\n", "\n", "\n", "\n", "\n", "**Author:** Robin Syring \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:33 PM." ] }, { "cell_type": "code", "execution_count": 1, "id": "cd1ccde4", "metadata": { "collapsed": false, "execution": { "iopub.execute_input": "2026-05-19T20:33:56.154125Z", "iopub.status.busy": "2026-05-19T20:33:56.153995Z", "iopub.status.idle": "2026-05-19T20:33:56.174569Z", "shell.execute_reply": "2026-05-19T20:33:56.174050Z" } }, "outputs": [ { "name": "stderr", "output_type": "stream", "text": [ "[runStaticInitializersOnce]: Failed to materialize symbols: { (main, { __stmts__1, $.cling-module-256.__inits.0, __orc_init_func.cling-module-256, _GLOBAL__sub_I_cling_module_256, _ZN13RooMsgService18setGlobalKillBelowEN6RooFit8MsgLevelE }) }\n" ] } ], "source": [ "%%cpp -d\n", "#include \"RooRealVar.h\"\n", "#include \"RooRealVar.h\"\n", "#include \"RooWorkspace.h\"\n", "#include \"RooGaussian.h\"\n", "#include \"RooUniform.h\"\n", "#include \"RooDataSet.h\"\n", "#include \"RooPlot.h\"\n", "#include \"RooMomentMorphFuncND.h\"\n", "#include \"RooAbsPdf.h\"\n", "\n", "using namespace RooFit;\n", "\n", "const int n_samples = 1000;\n", "\n", "RooMsgService::instance().setGlobalKillBelow(RooFit::WARNING);" ] }, { "cell_type": "markdown", "id": "ffc4a50f", "metadata": {}, "source": [ " Number of samples to fill the histograms\n", "Kills warning massages\n", "Define the morphing routine\n", " " ] }, { "cell_type": "code", "execution_count": 2, "id": "b153315f", "metadata": { "collapsed": false, "execution": { "iopub.execute_input": "2026-05-19T20:33:56.176089Z", "iopub.status.busy": "2026-05-19T20:33:56.175967Z", "iopub.status.idle": "2026-05-19T20:33:56.251438Z", "shell.execute_reply": "2026-05-19T20:33:56.250715Z" } }, "outputs": [ { "name": "stderr", "output_type": "stream", "text": [ "[runStaticInitializersOnce]: Failed to materialize symbols: { (main, { __orc_init_func.cling-module-256 }) }\n" ] } ], "source": [ "%%cpp -d\n", "RooPlot *perform_morphing(RooWorkspace &ws, RooMomentMorphFuncND::Setting setting, double sigma)\n", "{\n", " // Get Variables from the workspace\n", " RooRealVar *x_var = ws.var(\"x\");\n", " RooRealVar *mu_var = ws.var(\"mu\");\n", " RooAbsPdf *gauss = ws.pdf(\"gauss\");\n", "\n", " // Initialize a plot\n", " RooPlot *frame1 = x_var->frame();\n", "\n", " // Define binning for morphing\n", " RooMomentMorphFuncND::Grid grid(RooBinning(4, 0.0, 4.0));\n", "\n", " // Set binning of histograms, has to be customized for optimal results\n", " x_var->setBins(50);\n", "\n", " std::vector parampoints = {0, 1, 2, 3, 4};\n", "\n", " for (auto i : parampoints) {\n", " // Define the sampled gaussians\n", " RooRealVar mu_help(Form(\"mu%d\", i), Form(\"mu%d\", i), i);\n", " // Use * because RooGaussian expects objects no pointers\n", " RooGaussian help(Form(\"g%d\", i), Form(\"g%d\", i), *x_var, mu_help, sigma);\n", " ws.import(help, Silence(true));\n", "\n", " // Fill the histograms use a unique pointer to prevent memory leaks\n", " std::unique_ptr hist1{\n", " dynamic_cast(ws.pdf(Form(\"g%d\", i))->generateBinned(*x_var, 100 * n_samples))};\n", "\n", " // Add the value 1 to each bin\n", " for (int i_bin = 0; i_bin < hist1->numEntries(); ++i_bin) {\n", " const RooArgSet *binContent = hist1->get(i_bin);\n", " hist1->add(*binContent, 1.0);\n", " }\n", "\n", " // Add the pdf to the workspace, the inOrder of 1 is necessary for calculation of the nll\n", " // Adjust it to 0 to see binning\n", " ws.import(RooHistPdf(Form(\"histpdf%d\", i), Form(\"histpdf%d\", i), *x_var, *hist1, 1), Silence(true));\n", "\n", " // Plot and add the pdf to the grid\n", " RooAbsPdf *pdf = ws.pdf(Form(\"histpdf%d\", i));\n", " pdf->plotOn(frame1);\n", " grid.addPdf(*pdf, i);\n", " }\n", "\n", " // Create the morphing\n", " RooMomentMorphFuncND morph_func(\"morpf_func\", \"morph_func\", RooArgList(*mu_var), RooArgList(*x_var), grid, setting);\n", "\n", " // Normalizing the morphed object to be a pdf, set it false to prevent warning messages and gain computational speed\n", " // up\n", " morph_func.setPdfMode();\n", "\n", " // Creating the morphed pdf\n", " RooWrapperPdf morph(\"morph\", \"morph\", morph_func, true);\n", " ws.import(morph, Silence(true));\n", " RooAbsPdf *morph_ = ws.pdf(\"morph\");\n", " morph_->plotOn(frame1, LineColor(kRed));\n", "\n", " return frame1;\n", "}" ] }, { "cell_type": "markdown", "id": "96da959a", "metadata": {}, "source": [ " Define the workspace\n", " " ] }, { "cell_type": "code", "execution_count": 3, "id": "4269d451", "metadata": { "collapsed": false, "execution": { "iopub.execute_input": "2026-05-19T20:33:56.253835Z", "iopub.status.busy": "2026-05-19T20:33:56.253702Z", "iopub.status.idle": "2026-05-19T20:33:56.272807Z", "shell.execute_reply": "2026-05-19T20:33:56.272234Z" } }, "outputs": [], "source": [ "%%cpp -d\n", "std::unique_ptr build_ws(double mu_observed, double sigma)\n", "{\n", " auto ws = std::make_unique();\n", " ws->factory(Form(\"Gaussian::gauss(x[-5,15],mu[%f,0,4], %f)\", mu_observed, sigma));\n", " return ws;\n", "}" ] }, { "cell_type": "markdown", "id": "c9f5ebcd", "metadata": {}, "source": [ "Define the 'observed' mu" ] }, { "cell_type": "code", "execution_count": 4, "id": "01efd0f4", "metadata": { "collapsed": false, "execution": { "iopub.execute_input": "2026-05-19T20:33:56.274236Z", "iopub.status.busy": "2026-05-19T20:33:56.274112Z", "iopub.status.idle": "2026-05-19T20:33:56.594787Z", "shell.execute_reply": "2026-05-19T20:33:56.594066Z" } }, "outputs": [], "source": [ "double mu_observed = 2.5;\n", "double sigma = 1.5;" ] }, { "cell_type": "markdown", "id": "4f52add6", "metadata": {}, "source": [ "Import variables from workspace" ] }, { "cell_type": "code", "execution_count": 5, "id": "510f6300", "metadata": { "collapsed": false, "execution": { "iopub.execute_input": "2026-05-19T20:33:56.596602Z", "iopub.status.busy": "2026-05-19T20:33:56.596482Z", "iopub.status.idle": "2026-05-19T20:33:57.045032Z", "shell.execute_reply": "2026-05-19T20:33:57.044489Z" } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[#1] INFO:ObjectHandling -- RooWorkspace::import() importing dataset genData\n", "[#1] INFO:ObjectHandling -- RooWorkspace::import() importing dataset genData\n", "[#1] INFO:ObjectHandling -- RooWorkSpace::import() changing name of dataset from genData to genData_histpdf1\n", "[#1] INFO:ObjectHandling -- RooWorkspace::import() importing dataset genData\n", "[#1] INFO:ObjectHandling -- RooWorkSpace::import() changing name of dataset from genData to genData_histpdf2\n", "[#1] INFO:ObjectHandling -- RooWorkspace::import() importing dataset genData\n", "[#1] INFO:ObjectHandling -- RooWorkSpace::import() changing name of dataset from genData to genData_histpdf3\n", "[#1] INFO:ObjectHandling -- RooWorkspace::import() importing dataset genData\n", "[#1] INFO:ObjectHandling -- RooWorkSpace::import() changing name of dataset from genData to genData_histpdf4\n", "[#1] INFO:NumericIntegration -- RooRealIntegral::init(histpdf0_MOMENT_1_x_product_Int[x]) using numeric integrator RooIntegrator1D to calculate Int(x)\n", "[#1] INFO:NumericIntegration -- RooRealIntegral::init(histpdf1_MOMENT_1_x_product_Int[x]) using numeric integrator RooIntegrator1D to calculate Int(x)\n", "[#1] INFO:NumericIntegration -- RooRealIntegral::init(histpdf2_MOMENT_1_x_product_Int[x]) using numeric integrator RooIntegrator1D to calculate Int(x)\n", "[#1] INFO:NumericIntegration -- RooRealIntegral::init(histpdf3_MOMENT_1_x_product_Int[x]) using numeric integrator RooIntegrator1D to calculate Int(x)\n", "[#1] INFO:NumericIntegration -- RooRealIntegral::init(histpdf4_MOMENT_1_x_product_Int[x]) using numeric integrator RooIntegrator1D to calculate Int(x)\n", "[#1] INFO:NumericIntegration -- RooRealIntegral::init(histpdf0_MOMENT_2C_x_product_Int[x]) using numeric integrator RooIntegrator1D to calculate Int(x)\n", "[#1] INFO:NumericIntegration -- RooRealIntegral::init(histpdf1_MOMENT_2C_x_product_Int[x]) using numeric integrator RooIntegrator1D to calculate Int(x)\n", "[#1] INFO:NumericIntegration -- RooRealIntegral::init(histpdf2_MOMENT_2C_x_product_Int[x]) using numeric integrator RooIntegrator1D to calculate Int(x)\n", "[#1] INFO:NumericIntegration -- RooRealIntegral::init(histpdf3_MOMENT_2C_x_product_Int[x]) using numeric integrator RooIntegrator1D to calculate Int(x)\n", "[#1] INFO:NumericIntegration -- RooRealIntegral::init(histpdf4_MOMENT_2C_x_product_Int[x]) using numeric integrator RooIntegrator1D to calculate Int(x)\n" ] }, { "name": "stderr", "output_type": "stream", "text": [ "input_line_52:2:2: warning: 'ws' shadows a declaration with the same name in the 'std' namespace; use '::ws' to reference this declaration\n", " std::unique_ptr ws = build_ws(mu_observed, sigma);\n", " ^\n" ] } ], "source": [ "std::unique_ptr ws = build_ws(mu_observed, sigma);\n", "\n", "RooPlot *frame1 = perform_morphing(*ws, RooMomentMorphFuncND::Linear, sigma);\n", "\n", "RooRealVar *x_var = ws->var(\"x\");\n", "RooRealVar *mu_var = ws->var(\"mu\");\n", "RooAbsPdf *gauss = ws->pdf(\"gauss\");\n", "RooDataSet *obs_data = gauss->generate(*x_var, n_samples);" ] }, { "cell_type": "markdown", "id": "3b605699", "metadata": {}, "source": [ "Create the exact negative log likelihood function for Gaussian model" ] }, { "cell_type": "code", "execution_count": 6, "id": "2239147b", "metadata": { "collapsed": false, "execution": { "iopub.execute_input": "2026-05-19T20:33:57.046856Z", "iopub.status.busy": "2026-05-19T20:33:57.046732Z", "iopub.status.idle": "2026-05-19T20:33:57.255933Z", "shell.execute_reply": "2026-05-19T20:33:57.255260Z" } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[#1] INFO:Fitting -- RooAbsPdf::fitTo(gauss_over_gauss_Int[x]) fixing normalization set for coefficient determination to observables in data\n", "[#1] INFO:Fitting -- using generic CPU library compiled with no vectorizations\n", "[#1] INFO:Fitting -- Creation of NLL object took 792.924 μs\n" ] } ], "source": [ "RooAbsReal *nll_gauss = gauss->createNLL(*obs_data);" ] }, { "cell_type": "markdown", "id": "fcafc06f", "metadata": {}, "source": [ "Create the morphed negative log likelihood function" ] }, { "cell_type": "code", "execution_count": 7, "id": "7e8e033a", "metadata": { "collapsed": false, "execution": { "iopub.execute_input": "2026-05-19T20:33:57.257422Z", "iopub.status.busy": "2026-05-19T20:33:57.257303Z", "iopub.status.idle": "2026-05-19T20:33:57.463508Z", "shell.execute_reply": "2026-05-19T20:33:57.463108Z" } }, "outputs": [ { "name": "stderr", "output_type": "stream", "text": [ "input_line_97:2:26: error: reference to 'ws' is ambiguous\n", " RooAbsReal *nll_morph = ws->pdf(\"morph\")->createNLL(*obs_data);\n", " ^\n", "input_line_52:2:32: note: candidate found by name lookup is 'ws'\n", " std::unique_ptr ws = build_ws(mu_observed, sigma);\n", " ^\n", "/usr/lib/gcc/x86_64-redhat-linux/14/../../../../include/c++/14/bits/istream.tcc:1073:5: note: candidate found by name lookup is 'std::ws'\n", " ws(basic_istream<_CharT, _Traits>& __in)\n", " ^\n" ] } ], "source": [ "RooAbsReal *nll_morph = ws->pdf(\"morph\")->createNLL(*obs_data);" ] }, { "cell_type": "markdown", "id": "1f5bc5dd", "metadata": {}, "source": [ "Plot the negative logarithmic summed likelihood" ] }, { "cell_type": "code", "execution_count": 8, "id": "0b8ff86f", "metadata": { "collapsed": false, "execution": { "iopub.execute_input": "2026-05-19T20:33:57.465242Z", "iopub.status.busy": "2026-05-19T20:33:57.465121Z", "iopub.status.idle": "2026-05-19T20:33:57.673938Z", "shell.execute_reply": "2026-05-19T20:33:57.673513Z" } }, "outputs": [ { "name": "stderr", "output_type": "stream", "text": [ "input_line_98:4:37: error: cannot take the address of an rvalue of type 'EColor'\n", "nll_morph->plotOn(frame2, LineColor(kRed), ShiftToZero(), Name(\"morph\"));\n", " ^~~~\n", "Error while creating dynamic expression for:\n", " nll_morph->plotOn(frame2, LineColor(kRed), ShiftToZero(), Name(\"morph\"))\n" ] } ], "source": [ "RooPlot *frame2 = mu_var->frame(Title(\"Negative log Likelihood\"));\n", "nll_gauss->plotOn(frame2, LineColor(kBlue), ShiftToZero(), Name(\"gauss\"));\n", "nll_morph->plotOn(frame2, LineColor(kRed), ShiftToZero(), Name(\"morph\"));\n", "\n", "TCanvas *c = new TCanvas(\"rf616_morphing\", \"rf616_morphing\", 800, 400);\n", "c->Divide(2);\n", "c->cd(1);\n", "gPad->SetLeftMargin(0.15);\n", "frame1->GetYaxis()->SetTitleOffset(1.4);\n", "frame1->Draw();\n", "c->cd(2);\n", "gPad->SetLeftMargin(0.15);\n", "frame2->GetYaxis()->SetTitleOffset(1.8);\n", "frame2->Draw();" ] }, { "cell_type": "markdown", "id": "5ad7f0ca", "metadata": {}, "source": [ "Compute the minimum of the nll via minuit" ] }, { "cell_type": "code", "execution_count": 9, "id": "1b45a045", "metadata": { "collapsed": false, "execution": { "iopub.execute_input": "2026-05-19T20:33:57.675361Z", "iopub.status.busy": "2026-05-19T20:33:57.675238Z", "iopub.status.idle": "2026-05-19T20:33:58.059128Z", "shell.execute_reply": "2026-05-19T20:33:58.058668Z" } }, "outputs": [ { "name": "stderr", "output_type": "stream", "text": [ "[runStaticInitializersOnce]: Failed to materialize symbols: { (main, { _ZN9__gnu_cxxneIPP10RooAbsRealSt6vectorIS2_SaIS2_EEEEbRKNS_17__normal_iteratorIT_T0_EESC_, _ZN9__gnu_cxx17__normal_iteratorIPP10RooAbsRealSt6vectorIS2_SaIS2_EEEC2ERKS3_, _ZN12__cling_N52416__cling_Un1Qu324EPv, _ZN9__gnu_cxx17__normal_iteratorIPP10RooAbsRealSt6vectorIS2_SaIS2_EEEppEv, _ZNK9__gnu_cxx17__normal_iteratorIPP10RooAbsRealSt6vectorIS2_SaIS2_EEE4baseEv, $.cling-module-329.__inits.0, _ZN5cling7runtime8internal15DynamicExprInfoC2EPKcPPvb, __orc_init_func.cling-module-329, _ZN12RooMinimizer6ConfigC1Ev, _ZNSt6vectorIP10RooAbsRealSaIS1_EE3endEv, _ZNK9__gnu_cxx17__normal_iteratorIPP10RooAbsRealSt6vectorIS2_SaIS2_EEEdeEv, _ZN5cling7runtime8internal15DynamicExprInfoC1EPKcPPvb, _ZN12RooMinimizer6ConfigC2Ev, _ZN9__gnu_cxx17__normal_iteratorIPP10RooAbsRealSt6vectorIS2_SaIS2_EEEC1ERKS3_, _GLOBAL__sub_I_cling_module_329, _ZN12RooMinimizer6ConfigD1Ev, _ZNSt6vectorIP10RooAbsRealSaIS1_EE5beginEv, _ZN12__cling_N52424__dynamic__cling_Un1Qu30E, _ZN12__cling_N5244nllsE, _ZNK5cling7runtime8internal15LifetimeHandler9getMemoryEv, cling_module_329_, cling_module_329_.7, _ZN12RooMinimizer6ConfigD2Ev }) }\n", "IncrementalExecutor::executeFunction: symbol '_ZN5cling7runtime8internal15LifetimeHandlerD1Ev' unresolved while linking [cling interface function]!\n", "You are probably missing the definition of cling::runtime::internal::LifetimeHandler::~LifetimeHandler()\n", "Maybe you need to load the corresponding shared library?\n", "IncrementalExecutor::executeFunction: symbol '_ZN5cling7runtime8internal15LifetimeHandlerC1EPNS1_15DynamicExprInfoEPN5clang11DeclContextEPKcPNS_11InterpreterE' unresolved while linking [cling interface function]!\n", "You are probably missing the definition of cling::runtime::internal::LifetimeHandler::LifetimeHandler(cling::runtime::internal::DynamicExprInfo*, clang::DeclContext*, char const*, cling::Interpreter*)\n", "Maybe you need to load the corresponding shared library?\n" ] } ], "source": [ "std::vector nlls = {nll_gauss, nll_morph};\n", "for (auto nll : nlls) {\n", " RooMinimizer minimizer(*nll);\n", " minimizer.setPrintLevel(-1);\n", " minimizer.minimize(\"Minuit2\");\n", " RooFitResult *result = minimizer.save();\n", " result->Print();\n", "}" ] }, { "cell_type": "markdown", "id": "f034de46", "metadata": {}, "source": [ "Draw all canvases " ] }, { "cell_type": "code", "execution_count": 10, "id": "b7f34991", "metadata": { "collapsed": false, "execution": { "iopub.execute_input": "2026-05-19T20:33:58.060658Z", "iopub.status.busy": "2026-05-19T20:33:58.060522Z", "iopub.status.idle": "2026-05-19T20:33:58.291163Z", "shell.execute_reply": "2026-05-19T20:33:58.290574Z" } }, "outputs": [], "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 }