Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
FitHelpers.cxx
Go to the documentation of this file.
1/// \cond ROOFIT_INTERNAL
2
3/*
4 * Project: RooFit
5 * Authors:
6 * Jonas Rembser, CERN 2023
7 *
8 * Copyright (c) 2023, CERN
9 *
10 * Redistribution and use in source and binary forms,
11 * with or without modification, are permitted according to the terms
12 * listed in LICENSE (http://roofit.sourceforge.net/license.txt)
13 */
14
15#include "FitHelpers.h"
16
17#include <RooAbsData.h>
18#include <RooAbsPdf.h>
19#include <RooAbsReal.h>
20#include <RooAddition.h>
21#include <RooBatchCompute.h>
22#include <RooBinSamplingPdf.h>
23#include <RooCategory.h>
24#include <RooCmdConfig.h>
25#include <RooConstraintSum.h>
26#include <RooDataHist.h>
27#include <RooDataSet.h>
28#include <RooDerivative.h>
29#include <RooFit/Evaluator.h>
32#include <RooFitResult.h>
33#include <RooFuncWrapper.h>
34#include <RooLinkedList.h>
35#include <RooMinimizer.h>
36#include <RooRealVar.h>
37#include <RooSimultaneous.h>
38#include <RooFormulaVar.h>
39
40#include <Math/CholeskyDecomp.h>
41
42#include "ConstraintHelpers.h"
43#include "RooEvaluatorWrapper.h"
44#include "RooFitImplHelpers.h"
46
47#ifdef ROOFIT_LEGACY_EVAL_BACKEND
48#include "RooChi2Var.h"
49#include "RooNLLVar.h"
50#include "RooXYChi2Var.h"
51#endif
52
53namespace {
54
55constexpr int extendedFitDefault = 2;
56
57////////////////////////////////////////////////////////////////////////////////
58/// Use the asymptotically correct approach to estimate errors in the presence of weights.
59/// This is slower but more accurate than `SumW2Error`. See also https://arxiv.org/abs/1911.01303).
60/// Applies the calculated covaraince matrix to the RooMinimizer and returns
61/// the quality of the covariance matrix.
62/// See also the documentation of RooAbsPdf::fitTo(), where this function is used.
63/// \param[in] minimizer The RooMinimizer to get the fit result from. The state
64/// of the minimizer will be altered by this function: the covariance
65/// matrix caltulated here will be applied to it via
66/// RooMinimizer::applyCovarianceMatrix().
67/// \param[in] data The dataset that was used for the fit.
69{
70 RooFormulaVar logpdf("logpdf", "log(pdf)", "log(@0)", pdf);
71 RooArgSet obs;
72 logpdf.getObservables(data.get(), obs);
73
74 // Warning if the dataset is binned. TODO: in some cases,
75 // people also use RooDataSet to encode binned data,
76 // e.g. for simultaneous fits. It would be useful to detect
77 // this in this future as well.
78 if (dynamic_cast<RooDataHist const *>(&data)) {
79 oocoutW(&pdf, InputArguments)
80 << "RooAbsPdf::fitTo(" << pdf.GetName()
81 << ") WARNING: Asymptotic error correction is requested for a binned data set. "
82 "This method is not designed to handle binned data. A standard chi2 fit will likely be more suitable.";
83 };
84
85 // Calculated corrected errors for weighted likelihood fits
86 std::unique_ptr<RooFitResult> rw(minimizer.save());
87 // Weighted inverse Hessian matrix
88 const TMatrixDSym &matV = rw->covarianceMatrix();
89 oocoutI(&pdf, Fitting)
90 << "RooAbsPdf::fitTo(" << pdf.GetName()
91 << ") Calculating covariance matrix according to the asymptotically correct approach. If you find this "
92 "method useful please consider citing https://arxiv.org/abs/1911.01303.\n";
93
94 // Initialise matrix containing first derivatives
95 int nFloatPars = rw->floatParsFinal().size();
97 for (int k = 0; k < nFloatPars; k++) {
98 for (int l = 0; l < nFloatPars; l++) {
99 num(k, l) = 0.0;
100 }
101 }
102
103 // Create derivative objects
104 std::vector<std::unique_ptr<RooDerivative>> derivatives;
105 const RooArgList &floated = rw->floatParsFinal();
107 logpdf.getParameters(data.get(), allparams);
108 std::unique_ptr<RooArgSet> floatingparams{allparams.selectByAttrib("Constant", false)};
109
110 const double eps = 1.0e-4;
111
112 // Calculate derivatives of logpdf
113 for (const auto paramresult : floated) {
114 auto paraminternal = static_cast<RooRealVar *>(floatingparams->find(*paramresult));
115 assert(floatingparams->find(*paramresult)->IsA() == RooRealVar::Class());
116 double error = static_cast<RooRealVar *>(paramresult)->getError();
117 derivatives.emplace_back(logpdf.derivative(*paraminternal, obs, 1, eps * error));
118 }
119
120 // Calculate derivatives for number of expected events, needed for extended ML fit
121 RooAbsPdf *extended_pdf = dynamic_cast<RooAbsPdf *>(&pdf);
122 std::vector<double> diffs_expected(floated.size(), 0.0);
123 if (extended_pdf && extended_pdf->expectedEvents(obs) != 0.0) {
124 for (std::size_t k = 0; k < floated.size(); k++) {
125 const auto paramresult = static_cast<RooRealVar *>(floated.at(k));
126 auto paraminternal = static_cast<RooRealVar *>(floatingparams->find(*paramresult));
127
128 *paraminternal = paramresult->getVal();
129 double error = paramresult->getError();
130 paraminternal->setVal(paramresult->getVal() + eps * error);
131 double expected_plus = log(extended_pdf->expectedEvents(obs));
132 paraminternal->setVal(paramresult->getVal() - eps * error);
133 double expected_minus = log(extended_pdf->expectedEvents(obs));
134 *paraminternal = paramresult->getVal();
135 double diff = (expected_plus - expected_minus) / (2.0 * eps * error);
136 diffs_expected[k] = diff;
137 }
138 }
139
140 // Loop over data
141 for (int j = 0; j < data.numEntries(); j++) {
142 // Sets obs to current data point, this is where the pdf will be evaluated
143 obs.assign(*data.get(j));
144 // Determine first derivatives
145 std::vector<double> diffs(floated.size(), 0.0);
146 for (std::size_t k = 0; k < floated.size(); k++) {
147 const auto paramresult = static_cast<RooRealVar *>(floated.at(k));
148 auto paraminternal = static_cast<RooRealVar *>(floatingparams->find(*paramresult));
149 // first derivative to parameter k at best estimate point for this measurement
150 double diff = derivatives[k]->getVal();
151 // need to reset to best fit point after differentiation
152 *paraminternal = paramresult->getVal();
153 diffs[k] = diff;
154 }
155
156 // Fill numerator matrix
157 for (std::size_t k = 0; k < floated.size(); k++) {
158 for (std::size_t l = 0; l < floated.size(); l++) {
159 num(k, l) += data.weightSquared() * (diffs[k] + diffs_expected[k]) * (diffs[l] + diffs_expected[l]);
160 }
161 }
162 }
163 num.Similarity(matV);
164
165 // Propagate corrected errors to parameters objects
166 minimizer.applyCovarianceMatrix(num);
167
168 // The derivatives are found in RooFit and not with the minimizer (e.g.
169 // minuit), so the quality of the corrected covariance matrix corresponds to
170 // the quality of the original covariance matrix
171 return rw->covQual();
172}
173
174////////////////////////////////////////////////////////////////////////////////
175/// Apply correction to errors and covariance matrix. This uses two covariance
176/// matrices, one with the weights, the other with squared weights, to obtain
177/// the correct errors for weighted likelihood fits.
178/// Applies the calculated covaraince matrix to the RooMinimizer and returns
179/// the quality of the covariance matrix.
180/// See also the documentation of RooAbsPdf::fitTo(), where this function is used.
181/// \param[in] minimizer The RooMinimizer to get the fit result from. The state
182/// of the minimizer will be altered by this function: the covariance
183/// matrix caltulated here will be applied to it via
184/// RooMinimizer::applyCovarianceMatrix().
185/// \param[in] nll The NLL object that was used for the fit.
186int calcSumW2CorrectedCovariance(RooAbsReal const &pdf, RooMinimizer &minimizer, RooAbsReal &nll)
187{
188 // Calculated corrected errors for weighted likelihood fits
189 std::unique_ptr<RooFitResult> rw{minimizer.save()};
190 nll.applyWeightSquared(true);
191 oocoutI(&pdf, Fitting) << "RooAbsPdf::fitTo(" << pdf.GetName()
192 << ") Calculating sum-of-weights-squared correction matrix for covariance matrix\n";
193 minimizer.hesse();
194 std::unique_ptr<RooFitResult> rw2{minimizer.save()};
195 nll.applyWeightSquared(false);
196
197 // Apply correction matrix
198 const TMatrixDSym &matV = rw->covarianceMatrix();
199 TMatrixDSym matC = rw2->covarianceMatrix();
201 if (!decomp) {
202 oocoutE(&pdf, Fitting) << "RooAbsPdf::fitTo(" << pdf.GetName()
203 << ") ERROR: Cannot apply sum-of-weights correction to covariance matrix: correction "
204 "matrix calculated with weight-squared is singular\n";
205 return -1;
206 }
207
208 // replace C by its inverse
209 decomp.Invert(matC);
210 // the class lies about the matrix being symmetric, so fill in the
211 // part above the diagonal
212 for (int i = 0; i < matC.GetNrows(); ++i) {
213 for (int j = 0; j < i; ++j) {
214 matC(j, i) = matC(i, j);
215 }
216 }
217 matC.Similarity(matV);
218 // C now contains V C^-1 V
219 // Propagate corrected errors to parameters objects
220 minimizer.applyCovarianceMatrix(matC);
221
222 return std::min(rw->covQual(), rw2->covQual());
223}
224
225/// Configuration struct for RooAbsPdf::minimizeNLL with all the default values
226/// that also should be taked as the default values for RooAbsPdf::fitTo.
227struct MinimizerConfig {
228 double recoverFromNaN = 10.;
229 int optConst = 2;
230 int verbose = 0;
231 int doSave = 0;
232 int doTimer = 0;
233 int printLevel = 1;
234 int strategy = 1;
235 int initHesse = 0;
236 int hesse = 1;
237 int minos = 0;
238 int numee = 10;
239 int doEEWall = 1;
240 int doWarn = 1;
241 int doSumW2 = -1;
242 int doAsymptotic = -1;
243 int maxCalls = -1;
244 int doOffset = -1;
245 int parallelize = 0;
246 bool enableParallelGradient = true;
247 bool enableParallelDescent = false;
248 bool timingAnalysis = false;
249 const RooArgSet *minosSet = nullptr;
250 std::string minType;
251 std::string minAlg = "minuit";
252};
253
255{
256 // Process automatic extended option
259 if (ext) {
260 oocoutI(&pdf, Minimization)
261 << "p.d.f. provides expected number of events, including extended term in likelihood." << std::endl;
262 }
263 return ext;
264 }
265 // If Extended(false) was explicitly set, but the pdf MUST be extended, then
266 // it's time to print an error. This happens when you're fitting a RooAddPdf
267 // with coefficient that represent yields, and without the additional
268 // constraint these coefficients are degenerate because the RooAddPdf
269 // normalizes itself. Nothing correct can come out of this.
270 if (extendedCmdArg == 0) {
272 std::string errMsg = "You used the Extended(false) option on a pdf where the fit MUST be extended! "
273 "The parameters are not well defined and you're getting nonsensical results.";
274 oocoutE(&pdf, InputArguments) << errMsg << std::endl;
275 }
276 }
277 return extendedCmdArg;
278}
279
280/// To set the fitrange attribute of the PDF and custom ranges for the
281/// observables so that RooPlot can automatically plot the fitting range.
282void resetFitrangeAttributes(RooAbsArg &pdf, RooAbsData const &data, std::string const &baseName, const char *rangeName,
283 bool splitRange)
284{
285 // Clear possible range attributes from previous fits.
286 pdf.removeStringAttribute("fitrange");
287
288 // No fitrange was specified, so we do nothing. Or "SplitRange" is used, and
289 // then there are no uniquely defined ranges for the observables (as they
290 // are different in each category).
291 if (!rangeName || splitRange)
292 return;
293
294 RooArgSet observables;
295 pdf.getObservables(data.get(), observables);
296
297 std::string fitrangeValue;
298 auto subranges = ROOT::Split(rangeName, ",");
299 for (auto const &subrange : subranges) {
300 if (subrange.empty())
301 continue;
302 std::string fitrangeValueSubrange = std::string("fit_") + baseName;
303 if (subranges.size() > 1) {
305 }
307 for (RooAbsArg *arg : observables) {
308
309 if (arg->isCategory())
310 continue;
311 auto &observable = static_cast<RooRealVar &>(*arg);
312
313 observable.setRange(fitrangeValueSubrange.c_str(), observable.getMin(subrange.c_str()),
314 observable.getMax(subrange.c_str()));
315 }
316 }
317 pdf.setStringAttribute("fitrange", fitrangeValue.substr(0, fitrangeValue.size() - 1).c_str());
318}
319
320std::unique_ptr<RooAbsArg> createSimultaneousNLL(RooSimultaneous const &simPdf, bool isSimPdfExtended,
321 std::string const &rangeName, RooFit::OffsetMode offset)
322{
323 RooAbsCategoryLValue const &simCat = simPdf.indexCat();
324
325 // Prepare the NLL terms for each component
327 for (auto const &catState : simCat) {
328 std::string const &catName = catState.first;
330
331 // If the channel is not in the selected range of the category variable, we
332 // won't create an for NLL this channel.
333 if (!rangeName.empty()) {
334 // Only the RooCategory supports ranges, not the other
335 // RooAbsCategoryLValue-derived classes.
336 auto simCatAsRooCategory = dynamic_cast<RooCategory const *>(&simCat);
337 if (simCatAsRooCategory && !simCatAsRooCategory->isStateInRange(rangeName.c_str(), catIndex)) {
338 continue;
339 }
340 }
341
342 if (RooAbsPdf *pdf = simPdf.getPdf(catName.c_str())) {
343 auto name = std::string("nll_") + pdf->GetName();
344 std::unique_ptr<RooArgSet> observables{
345 std::unique_ptr<RooArgSet>(pdf->getVariables())->selectByAttrib("__obs__", true)};
346 // In a simultaneous fit, it is allowed that only a subset of the pdfs
347 // are extended. Therefore, we have to make sure that we don't request
348 // extended NLL objects for channels that can't be extended.
349 const bool isPdfExtended = isSimPdfExtended && pdf->extendMode() != RooAbsPdf::CanNotBeExtended;
350 auto nll = std::make_unique<RooFit::Detail::RooNLLVarNew>(name.c_str(), name.c_str(), *pdf, *observables,
352 // Rename the special variables
353 nll->setPrefix(std::string("_") + catName + "_");
354 nllTerms.addOwned(std::move(nll));
355 }
356 }
357
358 for (auto *nll : static_range_cast<RooFit::Detail::RooNLLVarNew *>(nllTerms)) {
359 nll->setSimCount(nllTerms.size());
360 }
361
362 // Time to sum the NLLs
363 auto nll = std::make_unique<RooAddition>("mynll", "mynll", nllTerms);
364 nll->addOwnedComponents(std::move(nllTerms));
365 return nll;
366}
367
368std::unique_ptr<RooAbsReal> createNLLNew(RooAbsPdf &pdf, RooAbsData &data, std::unique_ptr<RooAbsReal> &&constraints,
369 std::string const &rangeName, RooArgSet const &projDeps, bool isExtended,
370 double integrateOverBinsPrecision, RooFit::OffsetMode offset)
371{
372 if (constraints) {
373 // The computation graph for the constraints is very small, no need to do
374 // the tracking of clean and dirty nodes here.
375 constraints->setOperMode(RooAbsArg::ADirty);
376 }
377
378 RooArgSet observables;
379 pdf.getObservables(data.get(), observables);
380 observables.remove(projDeps, true, true);
381
382 oocxcoutI(&pdf, Fitting) << "RooAbsPdf::fitTo(" << pdf.GetName()
383 << ") fixing normalization set for coefficient determination to observables in data"
384 << "\n";
385 pdf.fixAddCoefNormalization(observables, false);
386
387 // Deal with the IntegrateBins argument
389 std::unique_ptr<RooAbsPdf> wrappedPdf = RooBinSamplingPdf::create(pdf, data, integrateOverBinsPrecision);
391 if (wrappedPdf) {
392 binSamplingPdfs.addOwned(std::move(wrappedPdf));
393 }
394 // Done dealing with the IntegrateBins option
395
397
398 auto simPdf = dynamic_cast<RooSimultaneous *>(&finalPdf);
399 if (simPdf) {
400 simPdf->wrapPdfsInBinSamplingPdfs(data, integrateOverBinsPrecision);
401 nllTerms.addOwned(createSimultaneousNLL(*simPdf, isExtended, rangeName, offset));
402 } else {
403 nllTerms.addOwned(std::make_unique<RooFit::Detail::RooNLLVarNew>("RooNLLVarNew", "RooNLLVarNew", finalPdf,
404 observables, isExtended, offset));
405 }
406 if (constraints) {
407 nllTerms.addOwned(std::move(constraints));
408 }
409
410 std::string nllName = std::string("nll_") + pdf.GetName() + "_" + data.GetName();
411 auto nll = std::make_unique<RooAddition>(nllName.c_str(), nllName.c_str(), nllTerms);
412 nll->addOwnedComponents(std::move(binSamplingPdfs));
413 nll->addOwnedComponents(std::move(nllTerms));
414
415 return nll;
416}
417
418} // namespace
419
420namespace RooFit {
421namespace FitHelpers {
422
424{
425 // Default-initialized instance of MinimizerConfig to get the default
426 // minimizer parameter values.
428
429 pc.defineDouble("RecoverFromUndefinedRegions", "RecoverFromUndefinedRegions", 0, minimizerDefaults.recoverFromNaN);
430 pc.defineInt("optConst", "Optimize", 0, minimizerDefaults.optConst);
431 pc.defineInt("verbose", "Verbose", 0, minimizerDefaults.verbose);
432 pc.defineInt("doSave", "Save", 0, minimizerDefaults.doSave);
433 pc.defineInt("doTimer", "Timer", 0, minimizerDefaults.doTimer);
434 pc.defineInt("printLevel", "PrintLevel", 0, minimizerDefaults.printLevel);
435 pc.defineInt("strategy", "Strategy", 0, minimizerDefaults.strategy);
436 pc.defineInt("initHesse", "InitialHesse", 0, minimizerDefaults.initHesse);
437 pc.defineInt("hesse", "Hesse", 0, minimizerDefaults.hesse);
438 pc.defineInt("minos", "Minos", 0, minimizerDefaults.minos);
439 pc.defineInt("numee", "PrintEvalErrors", 0, minimizerDefaults.numee);
440 pc.defineInt("doEEWall", "EvalErrorWall", 0, minimizerDefaults.doEEWall);
441 pc.defineInt("doWarn", "Warnings", 0, minimizerDefaults.doWarn);
442 pc.defineInt("doSumW2", "SumW2Error", 0, minimizerDefaults.doSumW2);
443 pc.defineInt("doAsymptoticError", "AsymptoticError", 0, minimizerDefaults.doAsymptotic);
444 pc.defineInt("maxCalls", "MaxCalls", 0, minimizerDefaults.maxCalls);
445 pc.defineInt("doOffset", "OffsetLikelihood", 0, 0);
446 pc.defineInt("parallelize", "Parallelize", 0, 0); // Three parallelize arguments
447 pc.defineInt("enableParallelGradient", "ParallelGradientOptions", 0, 0);
448 pc.defineInt("enableParallelDescent", "ParallelDescentOptions", 0, 0);
449 pc.defineInt("timingAnalysis", "TimingAnalysis", 0, 0);
450 pc.defineString("mintype", "Minimizer", 0, minimizerDefaults.minType.c_str());
451 pc.defineString("minalg", "Minimizer", 1, minimizerDefaults.minAlg.c_str());
452 pc.defineSet("minosSet", "Minos", 0, minimizerDefaults.minosSet);
453}
454
455////////////////////////////////////////////////////////////////////////////////
456/// Minimizes a given NLL variable by finding the optimal parameters with the
457/// RooMinimzer. The NLL variable can be created with RooAbsPdf::createNLL.
458/// If you are looking for a function that combines likelihood creation with
459/// fitting, see RooAbsPdf::fitTo.
460/// \param[in] nll The negative log-likelihood variable to minimize.
461/// \param[in] data The dataset that was also used for the NLL. It's a necessary
462/// parameter because it is used in the asymptotic error correction.
463/// \param[in] cfg Configuration struct with all the configuration options for
464/// the RooMinimizer. These are a subset of the options that you can
465/// also pass to RooAbsPdf::fitTo via the RooFit command arguments.
466std::unique_ptr<RooFitResult> minimize(RooAbsReal &pdf, RooAbsReal &nll, RooAbsData const &data, RooCmdConfig const &pc)
467{
468 MinimizerConfig cfg;
469 cfg.recoverFromNaN = pc.getDouble("RecoverFromUndefinedRegions");
470 cfg.optConst = pc.getInt("optConst");
471 cfg.verbose = pc.getInt("verbose");
472 cfg.doSave = pc.getInt("doSave");
473 cfg.doTimer = pc.getInt("doTimer");
474 cfg.printLevel = pc.getInt("printLevel");
475 cfg.strategy = pc.getInt("strategy");
476 cfg.initHesse = pc.getInt("initHesse");
477 cfg.hesse = pc.getInt("hesse");
478 cfg.minos = pc.getInt("minos");
479 cfg.numee = pc.getInt("numee");
480 cfg.doEEWall = pc.getInt("doEEWall");
481 cfg.doWarn = pc.getInt("doWarn");
482 cfg.doSumW2 = pc.getInt("doSumW2");
483 cfg.doAsymptotic = pc.getInt("doAsymptoticError");
484 cfg.maxCalls = pc.getInt("maxCalls");
485 cfg.minosSet = pc.getSet("minosSet");
486 cfg.minType = pc.getString("mintype", "");
487 cfg.minAlg = pc.getString("minalg", "minuit");
488 cfg.doOffset = pc.getInt("doOffset");
489 cfg.parallelize = pc.getInt("parallelize");
490 cfg.enableParallelGradient = pc.getInt("enableParallelGradient");
491 cfg.enableParallelDescent = pc.getInt("enableParallelDescent");
492 cfg.timingAnalysis = pc.getInt("timingAnalysis");
493
494 // Determine if the dataset has weights
495 bool weightedData = data.isNonPoissonWeighted();
496
497 std::string msgPrefix = std::string{"RooAbsPdf::fitTo("} + pdf.GetName() + "): ";
498
499 // Warn user that a method to determine parameter uncertainties should be provided if weighted data is offered
500 if (weightedData && cfg.doSumW2 == -1 && cfg.doAsymptotic == -1) {
501 oocoutW(&pdf, InputArguments) << msgPrefix <<
502 R"(WARNING: a likelihood fit is requested of what appears to be weighted data.
503 While the estimated values of the parameters will always be calculated taking the weights into account,
504 there are multiple ways to estimate the errors of the parameters. You are advised to make an
505 explicit choice for the error calculation:
506 - Either provide SumW2Error(true), to calculate a sum-of-weights-corrected HESSE error matrix
507 (error will be proportional to the number of events in MC).
508 - Or provide SumW2Error(false), to return errors from original HESSE error matrix
509 (which will be proportional to the sum of the weights, i.e., a dataset with <sum of weights> events).
510 - Or provide AsymptoticError(true), to use the asymptotically correct expression
511 (for details see https://arxiv.org/abs/1911.01303)."
512)";
513 }
514
515 if (cfg.minos && (cfg.doSumW2 == 1 || cfg.doAsymptotic == 1)) {
516 oocoutE(&pdf, InputArguments)
517 << msgPrefix
518 << " sum-of-weights and asymptotic error correction do not work with MINOS errors. Not fitting.\n";
519 return nullptr;
520 }
521 if (cfg.doAsymptotic == 1 && cfg.minos) {
522 oocoutW(&pdf, InputArguments) << msgPrefix << "WARNING: asymptotic correction does not apply to MINOS errors\n";
523 }
524
525 // avoid setting both SumW2 and Asymptotic for uncertainty correction
526 if (cfg.doSumW2 == 1 && cfg.doAsymptotic == 1) {
527 oocoutE(&pdf, InputArguments) << msgPrefix
528 << "ERROR: Cannot compute both asymptotically correct and SumW2 errors.\n";
529 return nullptr;
530 }
531
532 // Instantiate RooMinimizer
534 minimizerConfig.enableParallelGradient = cfg.enableParallelGradient;
535 minimizerConfig.enableParallelDescent = cfg.enableParallelDescent;
536 minimizerConfig.parallelize = cfg.parallelize;
537 minimizerConfig.timingAnalysis = cfg.timingAnalysis;
538 minimizerConfig.offsetting = cfg.doOffset;
540
541 m.setMinimizerType(cfg.minType);
542 m.setEvalErrorWall(cfg.doEEWall);
543 m.setRecoverFromNaNStrength(cfg.recoverFromNaN);
544 m.setPrintEvalErrors(cfg.numee);
545 if (cfg.maxCalls > 0)
546 m.setMaxFunctionCalls(cfg.maxCalls);
547 if (cfg.printLevel != 1)
548 m.setPrintLevel(cfg.printLevel);
549 if (cfg.optConst)
550 m.optimizeConst(cfg.optConst); // Activate constant term optimization
551 if (cfg.verbose)
552 m.setVerbose(true); // Activate verbose options
553 if (cfg.doTimer)
554 m.setProfile(true); // Activate timer options
555 if (cfg.strategy != 1)
556 m.setStrategy(cfg.strategy); // Modify fit strategy
557 if (cfg.initHesse)
558 m.hesse(); // Initialize errors with hesse
559 m.minimize(cfg.minType.c_str(), cfg.minAlg.c_str()); // Minimize using chosen algorithm
560 if (cfg.hesse)
561 m.hesse(); // Evaluate errors with Hesse
562
563 int corrCovQual = -1;
564
565 if (m.getNPar() > 0) {
566 if (cfg.doAsymptotic == 1)
567 corrCovQual = calcAsymptoticCorrectedCovariance(pdf, m, data); // Asymptotically correct
568 if (cfg.doSumW2 == 1)
570 }
571
572 if (cfg.minos)
573 cfg.minosSet ? m.minos(*cfg.minosSet) : m.minos(); // Evaluate errs with Minos
574
575 // Optionally return fit result
576 std::unique_ptr<RooFitResult> ret;
577 if (cfg.doSave) {
578 auto name = std::string("fitresult_") + pdf.GetName() + "_" + data.GetName();
579 auto title = std::string("Result of fit of p.d.f. ") + pdf.GetName() + " to dataset " + data.GetName();
580 ret = std::unique_ptr<RooFitResult>{m.save(name.c_str(), title.c_str())};
581 if ((cfg.doSumW2 == 1 || cfg.doAsymptotic == 1) && m.getNPar() > 0)
582 ret->setCovQual(corrCovQual);
583 }
584
585 if (cfg.optConst)
586 m.optimizeConst(0);
587 return ret;
588}
589
590std::unique_ptr<RooAbsReal> createNLL(RooAbsPdf &pdf, RooAbsData &data, const RooLinkedList &cmdList)
591{
592 auto baseName = std::string("nll_") + pdf.GetName() + "_" + data.GetName();
593
594 // Select the pdf-specific commands
595 RooCmdConfig pc("RooAbsPdf::createNLL(" + std::string(pdf.GetName()) + ")");
596
597 pc.defineString("rangeName", "RangeWithName", 0, "", true);
598 pc.defineString("addCoefRange", "SumCoefRange", 0, "");
599 pc.defineString("globstag", "GlobalObservablesTag", 0, "");
600 pc.defineString("globssource", "GlobalObservablesSource", 0, "data");
601 pc.defineDouble("rangeLo", "Range", 0, -999.);
602 pc.defineDouble("rangeHi", "Range", 1, -999.);
603 pc.defineInt("splitRange", "SplitRange", 0, 0);
604 pc.defineInt("ext", "Extended", 0, extendedFitDefault);
605 pc.defineInt("numcpu", "NumCPU", 0, 1);
606 pc.defineInt("interleave", "NumCPU", 1, 0);
607 pc.defineInt("verbose", "Verbose", 0, 0);
608 pc.defineInt("optConst", "Optimize", 0, 0);
609 pc.defineInt("cloneData", "CloneData", 0, 2);
610 pc.defineSet("projDepSet", "ProjectedObservables", 0, nullptr);
611 pc.defineSet("cPars", "Constrain", 0, nullptr);
612 pc.defineSet("glObs", "GlobalObservables", 0, nullptr);
613 pc.defineInt("doOffset", "OffsetLikelihood", 0, 0);
614 pc.defineSet("extCons", "ExternalConstraints", 0, nullptr);
615 pc.defineInt("EvalBackend", "EvalBackend", 0, static_cast<int>(RooFit::EvalBackend::defaultValue()));
616 pc.defineDouble("IntegrateBins", "IntegrateBins", 0, -1.);
617 pc.defineMutex("Range", "RangeWithName");
618 pc.defineMutex("GlobalObservables", "GlobalObservablesTag");
619 pc.defineInt("ModularL", "ModularL", 0, 0);
620
621 // New style likelihoods define parallelization through Parallelize(...) on fitTo or attributes on
622 // RooMinimizer::Config.
623 pc.defineMutex("ModularL", "NumCPU");
624
625 // New style likelihoods define offsetting on minimizer, not on likelihood
626 pc.defineMutex("ModularL", "OffsetLikelihood");
627
628 // Process and check varargs
629 pc.process(cmdList);
630 if (!pc.ok(true)) {
631 return nullptr;
632 }
633
634 if (pc.getInt("ModularL")) {
635 int lut[3] = {2, 1, 0};
637 static_cast<RooFit::TestStatistics::RooAbsL::Extended>(lut[pc.getInt("ext")])};
638
642
643 if (auto tmp = pc.getSet("cPars"))
644 cParsSet.add(*tmp);
645
646 if (auto tmp = pc.getSet("extCons"))
647 extConsSet.add(*tmp);
648
649 if (auto tmp = pc.getSet("glObs"))
650 glObsSet.add(*tmp);
651
652 const std::string rangeName = pc.getString("globstag", "", false);
653
655 builder.Extended(ext)
656 .ConstrainedParameters(cParsSet)
657 .ExternalConstraints(extConsSet)
658 .GlobalObservables(glObsSet)
659 .GlobalObservablesTag(rangeName.c_str());
660
661 return std::make_unique<RooFit::TestStatistics::RooRealL>("likelihood", "", builder.build());
662 }
663
664 // Decode command line arguments
665 const char *rangeName = pc.getString("rangeName", nullptr, true);
666 const char *addCoefRangeName = pc.getString("addCoefRange", nullptr, true);
667 const bool ext = interpretExtendedCmdArg(pdf, pc.getInt("ext"));
668
669 int splitRange = pc.getInt("splitRange");
670 int optConst = pc.getInt("optConst");
671 int cloneData = pc.getInt("cloneData");
672 auto offset = static_cast<RooFit::OffsetMode>(pc.getInt("doOffset"));
673
674 // If no explicit cloneData command is specified, cloneData is set to true if optimization is activated
675 if (cloneData == 2) {
677 }
678
679 if (pc.hasProcessed("Range")) {
680 double rangeLo = pc.getDouble("rangeLo");
681 double rangeHi = pc.getDouble("rangeHi");
682
683 // Create range with name 'fit' with above limits on all observables
684 RooArgSet obs;
685 pdf.getObservables(data.get(), obs);
686 for (auto arg : obs) {
687 RooRealVar *rrv = dynamic_cast<RooRealVar *>(arg);
688 if (rrv)
689 rrv->setRange("fit", rangeLo, rangeHi);
690 }
691
692 // Set range name to be fitted to "fit"
693 rangeName = "fit";
694 }
695
696 // Set the fitrange attribute of th PDF, add observables ranges for plotting
698
699 RooArgSet projDeps;
700 auto tmp = pc.getSet("projDepSet");
701 if (tmp) {
702 projDeps.add(*tmp);
703 }
704
705 const std::string globalObservablesSource = pc.getString("globssource", "data", false);
706 if (globalObservablesSource != "data" && globalObservablesSource != "model") {
707 std::string errMsg = "RooAbsPdf::fitTo: GlobalObservablesSource can only be \"data\" or \"model\"!";
708 oocoutE(&pdf, InputArguments) << errMsg << std::endl;
709 throw std::invalid_argument(errMsg);
710 }
711 const bool takeGlobalObservablesFromData = globalObservablesSource == "data";
712
713 // Lambda function to create the correct constraint term for a PDF. In old
714 // RooFit, we use this PDF itself as the argument, for the new BatchMode
715 // we're passing a clone.
716 auto createConstr = [&]() -> std::unique_ptr<RooAbsReal> {
717 return createConstraintTerm(baseName + "_constr", // name
718 pdf, // pdf
719 data, // data
720 pc.getSet("cPars"), // Constrain RooCmdArg
721 pc.getSet("extCons"), // ExternalConstraints RooCmdArg
722 pc.getSet("glObs"), // GlobalObservables RooCmdArg
723 pc.getString("globstag", nullptr, true), // GlobalObservablesTag RooCmdArg
724 takeGlobalObservablesFromData); // From GlobalObservablesSource RooCmdArg
725 };
726
727 auto evalBackend = static_cast<RooFit::EvalBackend::Value>(pc.getInt("EvalBackend"));
728
729 // Construct BatchModeNLL if requested
731
732 // Set the normalization range. We need to do it now, because it will be
733 // considered in `compileForNormSet`.
734 std::string oldNormRange;
735 if (pdf.normRange()) {
736 oldNormRange = pdf.normRange();
737 }
738 pdf.setNormRange(rangeName);
739
741 pdf.getObservables(data.get(), normSet);
742 normSet.remove(projDeps, true, true);
743
744 pdf.setAttribute("SplitRange", splitRange);
745 pdf.setStringAttribute("RangeName", rangeName);
746
748 ctx.setLikelihoodMode(true);
749 std::unique_ptr<RooAbsArg> head = pdf.compileForNormSet(normSet, ctx);
750 std::unique_ptr<RooAbsPdf> pdfClone = std::unique_ptr<RooAbsPdf>{static_cast<RooAbsPdf *>(head.release())};
751
752 // reset attributes
753 pdf.setAttribute("SplitRange", false);
754 pdf.setStringAttribute("RangeName", nullptr);
755
756 // Reset the normalization range
757 pdf.setNormRange(oldNormRange.c_str());
758
759 if (addCoefRangeName) {
760 oocxcoutI(&pdf, Fitting) << "RooAbsPdf::fitTo(" << pdf.GetName()
761 << ") fixing interpretation of coefficients of any component to range "
762 << addCoefRangeName << "\n";
763 pdfClone->fixAddCoefRange(addCoefRangeName, false);
764 }
765
766 std::unique_ptr<RooAbsReal> compiledConstr;
767 if (std::unique_ptr<RooAbsReal> constr = createConstr()) {
769 compiledConstr->addOwnedComponents(std::move(constr));
770 }
771
772 auto nll = createNLLNew(*pdfClone, data, std::move(compiledConstr), rangeName ? rangeName : "", projDeps, ext,
773 pc.getDouble("IntegrateBins"), offset);
774
775 std::unique_ptr<RooAbsReal> nllWrapper;
776
779 bool createGradient = evalBackend == RooFit::EvalBackend::Value::Codegen;
780 auto simPdf = dynamic_cast<RooSimultaneous const *>(pdfClone.get());
781 nllWrapper = std::make_unique<RooFit::Experimental::RooFuncWrapper>("nll_func_wrapper", "nll_func_wrapper",
782 *nll, &data, simPdf, createGradient);
783 if (createGradient)
784 static_cast<Experimental::RooFuncWrapper &>(*nllWrapper).createGradient();
785 } else {
786 nllWrapper = std::make_unique<RooEvaluatorWrapper>(
787 *nll, &data, evalBackend == RooFit::EvalBackend::Value::Cuda, rangeName ? rangeName : "", pdfClone.get(),
788 takeGlobalObservablesFromData);
789 }
790
791 nllWrapper->addOwnedComponents(std::move(nll));
792 nllWrapper->addOwnedComponents(std::move(pdfClone));
793 return nllWrapper;
794 }
795
796 std::unique_ptr<RooAbsReal> nll;
797
798#ifdef ROOFIT_LEGACY_EVAL_BACKEND
799 bool verbose = pc.getInt("verbose");
800
801 int numcpu = pc.getInt("numcpu");
802 int numcpu_strategy = pc.getInt("interleave");
803 // strategy 3 works only for RooSimultaneous.
804 if (numcpu_strategy == 3 && !pdf.InheritsFrom("RooSimultaneous")) {
805 oocoutW(&pdf, Minimization) << "Cannot use a NumCpu Strategy = 3 when the pdf is not a RooSimultaneous, "
806 "falling back to default strategy = 0"
807 << std::endl;
808 numcpu_strategy = 0;
809 }
811
813 RooAbsPdf &actualPdf = binnedLInfo.binnedPdf ? *binnedLInfo.binnedPdf : pdf;
814
815 // Construct NLL
818 cfg.addCoefRangeName = addCoefRangeName ? addCoefRangeName : "";
819 cfg.nCPU = numcpu;
820 cfg.interleave = interl;
821 cfg.verbose = verbose;
822 cfg.splitCutRange = static_cast<bool>(splitRange);
823 cfg.cloneInputData = static_cast<bool>(cloneData);
824 cfg.integrateOverBinsPrecision = pc.getDouble("IntegrateBins");
825 cfg.binnedL = binnedLInfo.isBinnedL;
826 cfg.takeGlobalObservablesFromData = takeGlobalObservablesFromData;
827 cfg.rangeName = rangeName ? rangeName : "";
828 auto nllVar = std::make_unique<RooNLLVar>(baseName.c_str(), "-log(likelihood)", actualPdf, data, projDeps, ext, cfg);
829 nllVar->enableBinOffsetting(offset == RooFit::OffsetMode::Bin);
830 nll = std::move(nllVar);
832
833 // Include constraints, if any, in likelihood
834 if (std::unique_ptr<RooAbsReal> constraintTerm = createConstr()) {
835
836 // Even though it is technically only required when the computation graph
837 // is changed because global observables are taken from data, it is safer
838 // to clone the constraint model in general to reset the normalization
839 // integral caches and avoid ASAN build failures (the PDF of the main
840 // measurement is cloned too anyway, so not much overhead). This can be
841 // reconsidered after the caching of normalization sets by pointer is changed
842 // to a more memory-safe solution.
843 constraintTerm = RooHelpers::cloneTreeWithSameParameters(*constraintTerm, data.get());
844
845 // Redirect the global observables to the ones from the dataset if applicable.
846 constraintTerm->setData(data, false);
847
848 // The computation graph for the constraints is very small, no need to do
849 // the tracking of clean and dirty nodes here.
850 constraintTerm->setOperMode(RooAbsArg::ADirty);
851
852 auto orignll = std::move(nll);
853 nll = std::make_unique<RooAddition>((baseName + "_with_constr").c_str(), "nllWithCons",
854 RooArgSet(*orignll, *constraintTerm));
855 nll->addOwnedComponents(std::move(orignll), std::move(constraintTerm));
856 }
857
858 if (optConst) {
859 nll->constOptimizeTestStatistic(RooAbsArg::Activate, optConst > 1);
860 }
861
863 nll->enableOffsetting(true);
864 }
865#else
866 throw std::runtime_error("RooFit was not built with the legacy evaluation backend");
867#endif
868
869 return nll;
870}
871
872std::unique_ptr<RooAbsReal> createChi2(RooAbsReal &real, RooAbsData &data, const RooLinkedList &cmdList)
873{
874#ifdef ROOFIT_LEGACY_EVAL_BACKEND
875 const bool isDataHist = dynamic_cast<RooDataHist const *>(&data);
876
877 RooCmdConfig pc("createChi2(" + std::string(real.GetName()) + ")");
878
879 pc.defineInt("numcpu", "NumCPU", 0, 1);
880 pc.defineInt("verbose", "Verbose", 0, 0);
881
883
884 if (isDataHist) {
885 // Construct Chi2
887 std::string baseName = "chi2_" + std::string(real.GetName()) + "_" + data.GetName();
888
889 // Clear possible range attributes from previous fits.
890 real.removeStringAttribute("fitrange");
891
892 pc.defineInt("etype", "DataError", 0, (Int_t)RooDataHist::Auto);
893 pc.defineInt("extended", "Extended", 0, extendedFitDefault);
894 pc.defineInt("split_range", "SplitRange", 0, 0);
895 pc.defineDouble("integrate_bins", "IntegrateBins", 0, -1);
896 pc.defineString("addCoefRange", "SumCoefRange", 0, "");
897 pc.allowUndefined();
898
899 pc.process(cmdList);
900 if (!pc.ok(true)) {
901 return nullptr;
902 }
903
904 bool extended = false;
905 if (auto pdf = dynamic_cast<RooAbsPdf const *>(&real)) {
906 extended = interpretExtendedCmdArg(*pdf, pc.getInt("extended"));
907 }
908
909 RooDataHist::ErrorType etype = static_cast<RooDataHist::ErrorType>(pc.getInt("etype"));
910
911 const char *rangeName = pc.getString("rangeName", nullptr, true);
912 const char *addCoefRangeName = pc.getString("addCoefRange", nullptr, true);
913
914 cfg.rangeName = rangeName ? rangeName : "";
915 cfg.nCPU = pc.getInt("numcpu");
917 cfg.verbose = static_cast<bool>(pc.getInt("verbose"));
918 cfg.cloneInputData = false;
919 cfg.integrateOverBinsPrecision = pc.getDouble("integrate_bins");
920 cfg.addCoefRangeName = addCoefRangeName ? addCoefRangeName : "";
921 cfg.splitCutRange = static_cast<bool>(pc.getInt("split_range"));
922 auto chi2 = std::make_unique<RooChi2Var>(baseName.c_str(), baseName.c_str(), real,
923 static_cast<RooDataHist &>(data), extended, etype, cfg);
924
926
927 return chi2;
928 } else {
929 pc.defineInt("integrate", "Integrate", 0, 0);
930 pc.defineObject("yvar", "YVar", 0, nullptr);
931 pc.defineString("rangeName", "RangeWithName", 0, "", true);
932 pc.defineInt("interleave", "NumCPU", 1, 0);
933
934 // Process and check varargs
935 pc.process(cmdList);
936 if (!pc.ok(true)) {
937 return nullptr;
938 }
939
940 // Decode command line arguments
941 bool integrate = pc.getInt("integrate");
942 RooRealVar *yvar = static_cast<RooRealVar *>(pc.getObject("yvar"));
943 const char *rangeName = pc.getString("rangeName", nullptr, true);
944 Int_t numcpu = pc.getInt("numcpu");
945 Int_t numcpu_strategy = pc.getInt("interleave");
946 // strategy 3 works only for RooSimultaneous.
947 if (numcpu_strategy == 3 && !real.InheritsFrom("RooSimultaneous")) {
948 oocoutW(&real, Minimization) << "Cannot use a NumCpu Strategy = 3 when the pdf is not a RooSimultaneous, "
949 "falling back to default strategy = 0"
950 << std::endl;
951 numcpu_strategy = 0;
952 }
954 bool verbose = pc.getInt("verbose");
955
956 cfg.rangeName = rangeName ? rangeName : "";
957 cfg.nCPU = numcpu;
958 cfg.interleave = interl;
959 cfg.verbose = verbose;
960 cfg.verbose = false;
961
962 std::string name = "chi2_" + std::string(real.GetName()) + "_" + data.GetName();
963
964 return std::make_unique<RooXYChi2Var>(name.c_str(), name.c_str(), real, static_cast<RooDataSet &>(data), yvar,
965 integrate, cfg);
966 }
967#else
968 throw std::runtime_error("createChi2() is not supported without the legacy evaluation backend");
969 return nullptr;
970#endif
971}
972
973std::unique_ptr<RooFitResult> fitTo(RooAbsReal &real, RooAbsData &data, const RooLinkedList &cmdList, bool chi2)
974{
975 const bool isDataHist = dynamic_cast<RooDataHist const *>(&data);
976
977 RooCmdConfig pc("fitTo(" + std::string(real.GetName()) + ")");
978
980 std::string nllCmdListString;
981 if (!chi2) {
982 nllCmdListString = "ProjectedObservables,Extended,Range,"
983 "RangeWithName,SumCoefRange,NumCPU,SplitRange,Constrained,Constrain,ExternalConstraints,"
984 "CloneData,GlobalObservables,GlobalObservablesSource,GlobalObservablesTag,"
985 "EvalBackend,IntegrateBins,ModularL";
986
987 if (!cmdList.FindObject("ModularL") || static_cast<RooCmdArg *>(cmdList.FindObject("ModularL"))->getInt(0) == 0) {
988 nllCmdListString += ",OffsetLikelihood";
989 }
990 } else {
991 auto createChi2DataHistCmdArgs = "Range,RangeWithName,NumCPU,Optimize,IntegrateBins,ProjectedObservables,"
992 "AddCoefRange,SplitRange,DataError,Extended";
993 auto createChi2DataSetCmdArgs = "YVar,Integrate,RangeWithName,NumCPU,Verbose";
995 }
996
998
999 pc.defineDouble("prefit", "Prefit", 0, 0);
1001
1002 // Process and check varargs
1003 pc.process(fitCmdList);
1004 if (!pc.ok(true)) {
1005 return nullptr;
1006 }
1007
1008 // TimingAnalysis works only for RooSimultaneous.
1009 if (pc.getInt("timingAnalysis") && !real.InheritsFrom("RooSimultaneous")) {
1010 oocoutW(&real, Minimization) << "The timingAnalysis feature was built for minimization with RooSimultaneous "
1011 "and is not implemented for other PDF's. Please create a RooSimultaneous to "
1012 "enable this feature."
1013 << std::endl;
1014 }
1015
1016 // Decode command line arguments
1017 double prefit = pc.getDouble("prefit");
1018
1019 if (prefit != 0) {
1020 size_t nEvents = static_cast<size_t>(prefit * data.numEntries());
1021 if (prefit > 0.5 || nEvents < 100) {
1022 oocoutW(&real, InputArguments) << "PrefitDataFraction should be in suitable range."
1023 << "With the current PrefitDataFraction=" << prefit
1024 << ", the number of events would be " << nEvents << " out of "
1025 << data.numEntries() << ". Skipping prefit..." << std::endl;
1026 } else {
1027 size_t step = data.numEntries() / nEvents;
1028
1029 RooDataSet tiny("tiny", "tiny", *data.get(), data.isWeighted() ? RooFit::WeightVar() : RooCmdArg());
1030
1031 for (int i = 0; i < data.numEntries(); i += step) {
1032 const RooArgSet *event = data.get(i);
1033 tiny.add(*event, data.weight());
1034 }
1036 pc.filterCmdList(tinyCmdList, "Prefit,Hesse,Minos,Verbose,Save,Timer");
1039
1042
1043 fitTo(real, tiny, tinyCmdList, chi2);
1044 }
1045 }
1046
1048 if (pc.getInt("parallelize") != 0 || pc.getInt("enableParallelGradient") || pc.getInt("enableParallelDescent")) {
1049 // Set to new style likelihood if parallelization is requested
1052 }
1053
1054 std::unique_ptr<RooAbsReal> nll;
1055 if (chi2) {
1056 nll = std::unique_ptr<RooAbsReal>{isDataHist ? real.createChi2(static_cast<RooDataHist &>(data), nllCmdList)
1057 : real.createChi2(static_cast<RooDataSet &>(data), nllCmdList)};
1058 } else {
1059 nll = std::unique_ptr<RooAbsReal>{dynamic_cast<RooAbsPdf &>(real).createNLL(data, nllCmdList)};
1060 }
1061
1062 return RooFit::FitHelpers::minimize(real, *nll, data, pc);
1063}
1064
1065} // namespace FitHelpers
1066} // namespace RooFit
1067
1068/// \endcond
header file containing the templated implementation of matrix inversion routines for use with ROOT's ...
std::unique_ptr< RooAbsReal > createConstraintTerm(std::string const &name, RooAbsPdf const &pdf, RooAbsData const &data, RooArgSet const *constrainedParameters, RooArgSet const *externalConstraints, RooArgSet const *globalObservables, const char *globalObservablesTag, bool takeGlobalObservablesFromData)
Create the parameter constraint sum to add to the negative log-likelihood.
ROOT::RRangeCast< T, false, Range_t > static_range_cast(Range_t &&coll)
#define oocoutW(o, a)
#define oocoutE(o, a)
#define oocoutI(o, a)
#define oocxcoutI(o, a)
int Int_t
Definition RtypesCore.h:45
ROOT::Detail::TRangeCast< T, true > TRangeDynCast
TRangeDynCast is an adapter class that allows the typed iteration through a TCollection.
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void data
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t Float_t Float_t Float_t Int_t Int_t UInt_t UInt_t Rectangle_t Int_t Int_t Window_t TString Int_t GCValues_t GetPrimarySelectionOwner GetDisplay GetScreen GetColormap GetNativeEvent const char const char dpyName wid window const char font_name cursor keysym reg const char only_if_exist regb h Point_t winding char text const char depth char const char Int_t count const char ColorStruct_t color const char Pixmap_t Pixmap_t PictureAttributes_t attr const char char ret_data h unsigned char height h offset
char name[80]
Definition TGX11.cxx:110
class to compute the Cholesky decomposition of a matrix
Common abstract base class for objects that represent a value and a "shape" in RooFit.
Definition RooAbsArg.h:79
void setStringAttribute(const Text_t *key, const Text_t *value)
Associate string 'value' to this object under key 'key'.
RooFit::OwningPtr< RooArgSet > getObservables(const RooArgSet &set, bool valueOnly=true) const
Given a set of possible observables, return the observables that this PDF depends on.
void removeStringAttribute(const Text_t *key)
Delete a string attribute with a given key.
RooFit::OwningPtr< RooArgSet > getVariables(bool stripDisconnected=true) const
Return RooArgSet with all variables (tree leaf nodes of expression tree)
void setAttribute(const Text_t *name, bool value=true)
Set (default) or clear a named boolean attribute of this object.
Abstract base class for objects that represent a discrete value that can be set from the outside,...
virtual bool remove(const RooAbsArg &var, bool silent=false, bool matchByNameOnly=false)
Remove the specified argument from our list.
virtual bool add(const RooAbsArg &var, bool silent=false)
Add the specified argument to list.
void assign(const RooAbsCollection &other) const
Sets the value, cache and constant attribute of any argument in our set that also appears in the othe...
Abstract base class for binned and unbinned datasets.
Definition RooAbsData.h:57
Abstract interface for all probability density functions.
Definition RooAbsPdf.h:40
std::unique_ptr< RooAbsArg > compileForNormSet(RooArgSet const &normSet, RooFit::Detail::CompileContext &ctx) const override
void setNormRange(const char *rangeName)
@ CanBeExtended
Definition RooAbsPdf.h:212
@ MustBeExtended
Definition RooAbsPdf.h:212
@ CanNotBeExtended
Definition RooAbsPdf.h:212
const char * normRange() const
Definition RooAbsPdf.h:250
virtual ExtendMode extendMode() const
Returns ability of PDF to provide extended likelihood terms.
Definition RooAbsPdf.h:216
Abstract base class for objects that represent a real value and implements functionality common to al...
Definition RooAbsReal.h:59
virtual void fixAddCoefNormalization(const RooArgSet &addNormSet=RooArgSet(), bool force=true)
Fix the interpretation of the coefficient of any RooAddPdf component in the expression tree headed by...
static void setEvalErrorLoggingMode(ErrorLoggingMode m)
Set evaluation error logging mode.
RooArgList is a container object that can hold multiple RooAbsArg objects.
Definition RooArgList.h:22
RooArgSet is a container object that can hold multiple RooAbsArg objects.
Definition RooArgSet.h:24
RooArgSet * selectByAttrib(const char *name, bool value) const
Use RooAbsCollection::selectByAttrib(), but return as RooArgSet.
Definition RooArgSet.h:144
static std::unique_ptr< RooAbsPdf > create(RooAbsPdf &pdf, RooAbsData const &data, double precision)
Creates a wrapping RooBinSamplingPdf if appropriate.
Object to represent discrete states.
Definition RooCategory.h:28
Named container for two doubles, two integers two object points and three string pointers that can be...
Definition RooCmdArg.h:26
Int_t getInt(Int_t idx) const
Definition RooCmdArg.h:87
Configurable parser for RooCmdArg named arguments.
void defineMutex(const char *head, Args_t &&... tail)
Define arguments where any pair is mutually exclusive.
bool process(const RooCmdArg &arg)
Process given RooCmdArg.
bool hasProcessed(const char *cmdName) const
Return true if RooCmdArg with name 'cmdName' has been processed.
double getDouble(const char *name, double defaultValue=0.0) const
Return double property registered with name 'name'.
bool defineDouble(const char *name, const char *argName, int doubleNum, double defValue=0.0)
Define double property name 'name' mapped to double in slot 'doubleNum' in RooCmdArg with name argNam...
RooArgSet * getSet(const char *name, RooArgSet *set=nullptr) const
Return RooArgSet property registered with name 'name'.
bool defineSet(const char *name, const char *argName, int setNum, const RooArgSet *set=nullptr)
Define TObject property name 'name' mapped to object in slot 'setNum' in RooCmdArg with name argName ...
bool ok(bool verbose) const
Return true of parsing was successful.
bool defineObject(const char *name, const char *argName, int setNum, const TObject *obj=nullptr, bool isArray=false)
Define TObject property name 'name' mapped to object in slot 'setNum' in RooCmdArg with name argName ...
const char * getString(const char *name, const char *defaultValue="", bool convEmptyToNull=false) const
Return string property registered with name 'name'.
bool defineString(const char *name, const char *argName, int stringNum, const char *defValue="", bool appendMode=false)
Define double property name 'name' mapped to double in slot 'stringNum' in RooCmdArg with name argNam...
bool defineInt(const char *name, const char *argName, int intNum, int defValue=0)
Define integer property name 'name' mapped to integer in slot 'intNum' in RooCmdArg with name argName...
void allowUndefined(bool flag=true)
If flag is true the processing of unrecognized RooCmdArgs is not considered an error.
int getInt(const char *name, int defaultValue=0) const
Return integer property registered with name 'name'.
RooLinkedList filterCmdList(RooLinkedList &cmdInList, const char *cmdNameList, bool removeFromInList=true) const
Utility function to filter commands listed in cmdNameList from cmdInList.
TObject * getObject(const char *name, TObject *obj=nullptr) const
Return TObject property registered with name 'name'.
Container class to hold N-dimensional binned data.
Definition RooDataHist.h:40
Container class to hold unbinned data.
Definition RooDataSet.h:33
static Value & defaultValue()
A RooFormulaVar is a generic implementation of a real-valued object, which takes a RooArgList of serv...
Collection class for internal use, storing a collection of RooAbsArg pointers in a doubly linked list...
Wrapper class around ROOT::Math::Minimizer that provides a seamless interface between the minimizer f...
RooFit::OwningPtr< RooFitResult > save(const char *name=nullptr, const char *title=nullptr)
Save and return a RooFitResult snapshot of current minimizer status.
int hesse()
Execute HESSE.
void applyCovarianceMatrix(TMatrixDSym const &V)
Apply results of given external covariance matrix.
Variable that can be changed from the outside.
Definition RooRealVar.h:37
static TClass * Class()
void setRange(const char *name, double min, double max)
Set a fit or plotting range.
Facilitates simultaneous fitting of multiple PDFs to subsets of a given dataset.
const char * GetName() const override
Returns name of object.
Definition TNamed.h:47
virtual Bool_t InheritsFrom(const char *classname) const
Returns kTRUE if object inherits from class "classname".
Definition TObject.cxx:530
RooCmdArg WeightVar(const char *name="weight", bool reinterpretAsWeight=false)
RooCmdArg Hesse(bool flag=true)
RooCmdArg ModularL(bool flag=false)
RooCmdArg PrintLevel(Int_t code)
std::vector< std::string > Split(std::string_view str, std::string_view delims, bool skipEmpty=false)
Splits a string at each character in delims.
double nll(double pdf, double weight, int binnedL, int doBinOffset)
Definition MathFuncs.h:377
std::unique_ptr< T > compileForNormSet(T const &arg, RooArgSet const &normSet)
The namespace RooFit contains mostly switches that change the behaviour of functions of PDFs (or othe...
Definition JSONIO.h:26
OffsetMode
For setting the offset mode with the Offset() command argument to RooAbsPdf::fitTo()
std::unique_ptr< T > cloneTreeWithSameParameters(T const &arg, RooArgSet const *observables=nullptr)
Clone RooAbsArg object and reattach to original parameters.
BinnedLOutput getBinnedL(RooAbsPdf const &pdf)
std::string rangeName
Stores the configuration parameters for RooAbsTestStatistic.
Config argument to RooMinimizer constructor.
TMarker m
Definition textangle.C:8
TLine l
Definition textangle.C:4