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"
45#include "RooNLLVarNew.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.
68int calcAsymptoticCorrectedCovariance(RooAbsReal &pdf, RooMinimizer &minimizer, RooAbsData const &data)
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();
96 TMatrixDSym num(nFloatPars);
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();
106 RooArgSet allparams;
107 logpdf.getParameters(data.get(), allparams);
108 std::unique_ptr<RooArgSet> floatingparams{static_cast<RooArgSet *>(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
254bool interpretExtendedCmdArg(RooAbsPdf const &pdf, int extendedCmdArg)
255{
256 // Process automatic extended option
257 if (extendedCmdArg == extendedFitDefault) {
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) {
304 fitrangeValueSubrange += "_" + subrange;
305 }
306 fitrangeValue += fitrangeValueSubrange + ",";
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 isExtended,
321 std::string const &rangeName, RooFit::OffsetMode offset)
322{
323 RooAbsCategoryLValue const &simCat = simPdf.indexCat();
324
325 // Prepare the NLL terms for each component
326 RooArgList nllTerms;
327 for (auto const &catState : simCat) {
328 std::string const &catName = catState.first;
329 RooAbsCategory::value_type catIndex = catState.second;
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 static_cast<RooArgSet *>(std::unique_ptr<RooArgSet>(pdf->getVariables())->selectByAttrib("__obs__", true)));
346 auto nll = std::make_unique<RooNLLVarNew>(name.c_str(), name.c_str(), *pdf, *observables, isExtended, offset);
347 // Rename the special variables
348 nll->setPrefix(std::string("_") + catName + "_");
349 nllTerms.addOwned(std::move(nll));
350 }
351 }
352
353 for (auto *nll : static_range_cast<RooNLLVarNew *>(nllTerms)) {
354 nll->setSimCount(nllTerms.size());
355 }
356
357 // Time to sum the NLLs
358 auto nll = std::make_unique<RooAddition>("mynll", "mynll", nllTerms);
359 nll->addOwnedComponents(std::move(nllTerms));
360 return nll;
361}
362
363std::unique_ptr<RooAbsReal> createNLLNew(RooAbsPdf &pdf, RooAbsData &data, std::unique_ptr<RooAbsReal> &&constraints,
364 std::string const &rangeName, RooArgSet const &projDeps, bool isExtended,
365 double integrateOverBinsPrecision, RooFit::OffsetMode offset)
366{
367 if (constraints) {
368 // The computation graph for the constraints is very small, no need to do
369 // the tracking of clean and dirty nodes here.
370 constraints->setOperMode(RooAbsArg::ADirty);
371 }
372
373 RooArgSet observables;
374 pdf.getObservables(data.get(), observables);
375 observables.remove(projDeps, true, true);
376
377 oocxcoutI(&pdf, Fitting) << "RooAbsPdf::fitTo(" << pdf.GetName()
378 << ") fixing normalization set for coefficient determination to observables in data"
379 << "\n";
380 pdf.fixAddCoefNormalization(observables, false);
381
382 // Deal with the IntegrateBins argument
383 RooArgList binSamplingPdfs;
384 std::unique_ptr<RooAbsPdf> wrappedPdf = RooBinSamplingPdf::create(pdf, data, integrateOverBinsPrecision);
385 RooAbsPdf &finalPdf = wrappedPdf ? *wrappedPdf : pdf;
386 if (wrappedPdf) {
387 binSamplingPdfs.addOwned(std::move(wrappedPdf));
388 }
389 // Done dealing with the IntegrateBins option
390
391 RooArgList nllTerms;
392
393 auto simPdf = dynamic_cast<RooSimultaneous *>(&finalPdf);
394 if (simPdf) {
395 simPdf->wrapPdfsInBinSamplingPdfs(data, integrateOverBinsPrecision);
396 nllTerms.addOwned(createSimultaneousNLL(*simPdf, isExtended, rangeName, offset));
397 } else {
398 nllTerms.addOwned(
399 std::make_unique<RooNLLVarNew>("RooNLLVarNew", "RooNLLVarNew", finalPdf, observables, isExtended, offset));
400 }
401 if (constraints) {
402 nllTerms.addOwned(std::move(constraints));
403 }
404
405 std::string nllName = std::string("nll_") + pdf.GetName() + "_" + data.GetName();
406 auto nll = std::make_unique<RooAddition>(nllName.c_str(), nllName.c_str(), nllTerms);
407 nll->addOwnedComponents(std::move(binSamplingPdfs));
408 nll->addOwnedComponents(std::move(nllTerms));
409
410 return nll;
411}
412
413} // namespace
414
415namespace RooFit {
416namespace FitHelpers {
417
418void defineMinimizationOptions(RooCmdConfig &pc)
419{
420 // Default-initialized instance of MinimizerConfig to get the default
421 // minimizer parameter values.
422 MinimizerConfig minimizerDefaults;
423
424 pc.defineDouble("RecoverFromUndefinedRegions", "RecoverFromUndefinedRegions", 0, minimizerDefaults.recoverFromNaN);
425 pc.defineInt("optConst", "Optimize", 0, minimizerDefaults.optConst);
426 pc.defineInt("verbose", "Verbose", 0, minimizerDefaults.verbose);
427 pc.defineInt("doSave", "Save", 0, minimizerDefaults.doSave);
428 pc.defineInt("doTimer", "Timer", 0, minimizerDefaults.doTimer);
429 pc.defineInt("printLevel", "PrintLevel", 0, minimizerDefaults.printLevel);
430 pc.defineInt("strategy", "Strategy", 0, minimizerDefaults.strategy);
431 pc.defineInt("initHesse", "InitialHesse", 0, minimizerDefaults.initHesse);
432 pc.defineInt("hesse", "Hesse", 0, minimizerDefaults.hesse);
433 pc.defineInt("minos", "Minos", 0, minimizerDefaults.minos);
434 pc.defineInt("numee", "PrintEvalErrors", 0, minimizerDefaults.numee);
435 pc.defineInt("doEEWall", "EvalErrorWall", 0, minimizerDefaults.doEEWall);
436 pc.defineInt("doWarn", "Warnings", 0, minimizerDefaults.doWarn);
437 pc.defineInt("doSumW2", "SumW2Error", 0, minimizerDefaults.doSumW2);
438 pc.defineInt("doAsymptoticError", "AsymptoticError", 0, minimizerDefaults.doAsymptotic);
439 pc.defineInt("maxCalls", "MaxCalls", 0, minimizerDefaults.maxCalls);
440 pc.defineInt("doOffset", "OffsetLikelihood", 0, 0);
441 pc.defineInt("parallelize", "Parallelize", 0, 0); // Three parallelize arguments
442 pc.defineInt("enableParallelGradient", "ParallelGradientOptions", 0, 0);
443 pc.defineInt("enableParallelDescent", "ParallelDescentOptions", 0, 0);
444 pc.defineInt("timingAnalysis", "TimingAnalysis", 0, 0);
445 pc.defineString("mintype", "Minimizer", 0, minimizerDefaults.minType.c_str());
446 pc.defineString("minalg", "Minimizer", 1, minimizerDefaults.minAlg.c_str());
447 pc.defineSet("minosSet", "Minos", 0, minimizerDefaults.minosSet);
448}
449
450////////////////////////////////////////////////////////////////////////////////
451/// Minimizes a given NLL variable by finding the optimal parameters with the
452/// RooMinimzer. The NLL variable can be created with RooAbsPdf::createNLL.
453/// If you are looking for a function that combines likelihood creation with
454/// fitting, see RooAbsPdf::fitTo.
455/// \param[in] nll The negative log-likelihood variable to minimize.
456/// \param[in] data The dataset that was also used for the NLL. It's a necessary
457/// parameter because it is used in the asymptotic error correction.
458/// \param[in] cfg Configuration struct with all the configuration options for
459/// the RooMinimizer. These are a subset of the options that you can
460/// also pass to RooAbsPdf::fitTo via the RooFit command arguments.
461std::unique_ptr<RooFitResult> minimize(RooAbsReal &pdf, RooAbsReal &nll, RooAbsData const &data, RooCmdConfig const &pc)
462{
463 MinimizerConfig cfg;
464 cfg.recoverFromNaN = pc.getDouble("RecoverFromUndefinedRegions");
465 cfg.optConst = pc.getInt("optConst");
466 cfg.verbose = pc.getInt("verbose");
467 cfg.doSave = pc.getInt("doSave");
468 cfg.doTimer = pc.getInt("doTimer");
469 cfg.printLevel = pc.getInt("printLevel");
470 cfg.strategy = pc.getInt("strategy");
471 cfg.initHesse = pc.getInt("initHesse");
472 cfg.hesse = pc.getInt("hesse");
473 cfg.minos = pc.getInt("minos");
474 cfg.numee = pc.getInt("numee");
475 cfg.doEEWall = pc.getInt("doEEWall");
476 cfg.doWarn = pc.getInt("doWarn");
477 cfg.doSumW2 = pc.getInt("doSumW2");
478 cfg.doAsymptotic = pc.getInt("doAsymptoticError");
479 cfg.maxCalls = pc.getInt("maxCalls");
480 cfg.minosSet = pc.getSet("minosSet");
481 cfg.minType = pc.getString("mintype", "");
482 cfg.minAlg = pc.getString("minalg", "minuit");
483 cfg.doOffset = pc.getInt("doOffset");
484 cfg.parallelize = pc.getInt("parallelize");
485 cfg.enableParallelGradient = pc.getInt("enableParallelGradient");
486 cfg.enableParallelDescent = pc.getInt("enableParallelDescent");
487 cfg.timingAnalysis = pc.getInt("timingAnalysis");
488
489 // Determine if the dataset has weights
490 bool weightedData = data.isNonPoissonWeighted();
491
492 std::string msgPrefix = std::string{"RooAbsPdf::fitTo("} + pdf.GetName() + "): ";
493
494 // Warn user that a method to determine parameter uncertainties should be provided if weighted data is offered
495 if (weightedData && cfg.doSumW2 == -1 && cfg.doAsymptotic == -1) {
496 oocoutW(&pdf, InputArguments) << msgPrefix <<
497 R"(WARNING: a likelihood fit is requested of what appears to be weighted data.
498 While the estimated values of the parameters will always be calculated taking the weights into account,
499 there are multiple ways to estimate the errors of the parameters. You are advised to make an
500 explicit choice for the error calculation:
501 - Either provide SumW2Error(true), to calculate a sum-of-weights-corrected HESSE error matrix
502 (error will be proportional to the number of events in MC).
503 - Or provide SumW2Error(false), to return errors from original HESSE error matrix
504 (which will be proportional to the sum of the weights, i.e., a dataset with <sum of weights> events).
505 - Or provide AsymptoticError(true), to use the asymptotically correct expression
506 (for details see https://arxiv.org/abs/1911.01303)."
507)";
508 }
509
510 if (cfg.minos && (cfg.doSumW2 == 1 || cfg.doAsymptotic == 1)) {
511 oocoutE(&pdf, InputArguments)
512 << msgPrefix
513 << " sum-of-weights and asymptotic error correction do not work with MINOS errors. Not fitting.\n";
514 return nullptr;
515 }
516 if (cfg.doAsymptotic == 1 && cfg.minos) {
517 oocoutW(&pdf, InputArguments) << msgPrefix << "WARNING: asymptotic correction does not apply to MINOS errors\n";
518 }
519
520 // avoid setting both SumW2 and Asymptotic for uncertainty correction
521 if (cfg.doSumW2 == 1 && cfg.doAsymptotic == 1) {
522 oocoutE(&pdf, InputArguments) << msgPrefix
523 << "ERROR: Cannot compute both asymptotically correct and SumW2 errors.\n";
524 return nullptr;
525 }
526
527 // Instantiate RooMinimizer
528 RooMinimizer::Config minimizerConfig;
529 minimizerConfig.enableParallelGradient = cfg.enableParallelGradient;
530 minimizerConfig.enableParallelDescent = cfg.enableParallelDescent;
531 minimizerConfig.parallelize = cfg.parallelize;
532 minimizerConfig.timingAnalysis = cfg.timingAnalysis;
533 minimizerConfig.offsetting = cfg.doOffset;
534 RooMinimizer m(nll, minimizerConfig);
535
536 m.setMinimizerType(cfg.minType);
537 m.setEvalErrorWall(cfg.doEEWall);
538 m.setRecoverFromNaNStrength(cfg.recoverFromNaN);
539 m.setPrintEvalErrors(cfg.numee);
540 if (cfg.maxCalls > 0)
541 m.setMaxFunctionCalls(cfg.maxCalls);
542 if (cfg.printLevel != 1)
543 m.setPrintLevel(cfg.printLevel);
544 if (cfg.optConst)
545 m.optimizeConst(cfg.optConst); // Activate constant term optimization
546 if (cfg.verbose)
547 m.setVerbose(true); // Activate verbose options
548 if (cfg.doTimer)
549 m.setProfile(true); // Activate timer options
550 if (cfg.strategy != 1)
551 m.setStrategy(cfg.strategy); // Modify fit strategy
552 if (cfg.initHesse)
553 m.hesse(); // Initialize errors with hesse
554 m.minimize(cfg.minType.c_str(), cfg.minAlg.c_str()); // Minimize using chosen algorithm
555 if (cfg.hesse)
556 m.hesse(); // Evaluate errors with Hesse
557
558 int corrCovQual = -1;
559
560 if (m.getNPar() > 0) {
561 if (cfg.doAsymptotic == 1)
562 corrCovQual = calcAsymptoticCorrectedCovariance(pdf, m, data); // Asymptotically correct
563 if (cfg.doSumW2 == 1)
564 corrCovQual = calcSumW2CorrectedCovariance(pdf, m, nll);
565 }
566
567 if (cfg.minos)
568 cfg.minosSet ? m.minos(*cfg.minosSet) : m.minos(); // Evaluate errs with Minos
569
570 // Optionally return fit result
571 std::unique_ptr<RooFitResult> ret;
572 if (cfg.doSave) {
573 auto name = std::string("fitresult_") + pdf.GetName() + "_" + data.GetName();
574 auto title = std::string("Result of fit of p.d.f. ") + pdf.GetName() + " to dataset " + data.GetName();
575 ret = std::unique_ptr<RooFitResult>{m.save(name.c_str(), title.c_str())};
576 if ((cfg.doSumW2 == 1 || cfg.doAsymptotic == 1) && m.getNPar() > 0)
577 ret->setCovQual(corrCovQual);
578 }
579
580 if (cfg.optConst)
581 m.optimizeConst(0);
582 return ret;
583}
584
585std::unique_ptr<RooAbsReal> createNLL(RooAbsPdf &pdf, RooAbsData &data, const RooLinkedList &cmdList)
586{
587 auto baseName = std::string("nll_") + pdf.GetName() + "_" + data.GetName();
588
589 // Select the pdf-specific commands
590 RooCmdConfig pc("RooAbsPdf::createNLL(" + std::string(pdf.GetName()) + ")");
591
592 pc.defineString("rangeName", "RangeWithName", 0, "", true);
593 pc.defineString("addCoefRange", "SumCoefRange", 0, "");
594 pc.defineString("globstag", "GlobalObservablesTag", 0, "");
595 pc.defineString("globssource", "GlobalObservablesSource", 0, "data");
596 pc.defineDouble("rangeLo", "Range", 0, -999.);
597 pc.defineDouble("rangeHi", "Range", 1, -999.);
598 pc.defineInt("splitRange", "SplitRange", 0, 0);
599 pc.defineInt("ext", "Extended", 0, extendedFitDefault);
600 pc.defineInt("numcpu", "NumCPU", 0, 1);
601 pc.defineInt("interleave", "NumCPU", 1, 0);
602 pc.defineInt("verbose", "Verbose", 0, 0);
603 pc.defineInt("optConst", "Optimize", 0, 0);
604 pc.defineInt("cloneData", "CloneData", 0, 2);
605 pc.defineSet("projDepSet", "ProjectedObservables", 0, nullptr);
606 pc.defineSet("cPars", "Constrain", 0, nullptr);
607 pc.defineSet("glObs", "GlobalObservables", 0, nullptr);
608 pc.defineInt("doOffset", "OffsetLikelihood", 0, 0);
609 pc.defineSet("extCons", "ExternalConstraints", 0, nullptr);
610 pc.defineInt("EvalBackend", "EvalBackend", 0, static_cast<int>(RooFit::EvalBackend::defaultValue()));
611 pc.defineDouble("IntegrateBins", "IntegrateBins", 0, -1.);
612 pc.defineMutex("Range", "RangeWithName");
613 pc.defineMutex("GlobalObservables", "GlobalObservablesTag");
614 pc.defineInt("ModularL", "ModularL", 0, 0);
615
616 // New style likelihoods define parallelization through Parallelize(...) on fitTo or attributes on
617 // RooMinimizer::Config.
618 pc.defineMutex("ModularL", "NumCPU");
619
620 // New style likelihoods define offsetting on minimizer, not on likelihood
621 pc.defineMutex("ModularL", "OffsetLikelihood");
622
623 // Process and check varargs
624 pc.process(cmdList);
625 if (!pc.ok(true)) {
626 return nullptr;
627 }
628
629 if (pc.getInt("ModularL")) {
630 int lut[3] = {2, 1, 0};
632 static_cast<RooFit::TestStatistics::RooAbsL::Extended>(lut[pc.getInt("ext")])};
633
634 RooArgSet cParsSet;
635 RooArgSet extConsSet;
636 RooArgSet glObsSet;
637
638 if (auto tmp = pc.getSet("cPars"))
639 cParsSet.add(*tmp);
640
641 if (auto tmp = pc.getSet("extCons"))
642 extConsSet.add(*tmp);
643
644 if (auto tmp = pc.getSet("glObs"))
645 glObsSet.add(*tmp);
646
647 const std::string rangeName = pc.getString("globstag", "", false);
648
650 builder.Extended(ext)
651 .ConstrainedParameters(cParsSet)
652 .ExternalConstraints(extConsSet)
653 .GlobalObservables(glObsSet)
654 .GlobalObservablesTag(rangeName.c_str());
655
656 return std::make_unique<RooFit::TestStatistics::RooRealL>("likelihood", "", builder.build());
657 }
658
659 // Decode command line arguments
660 const char *rangeName = pc.getString("rangeName", nullptr, true);
661 const char *addCoefRangeName = pc.getString("addCoefRange", nullptr, true);
662 const bool ext = interpretExtendedCmdArg(pdf, pc.getInt("ext"));
663
664 int splitRange = pc.getInt("splitRange");
665 int optConst = pc.getInt("optConst");
666 int cloneData = pc.getInt("cloneData");
667 auto offset = static_cast<RooFit::OffsetMode>(pc.getInt("doOffset"));
668
669 // If no explicit cloneData command is specified, cloneData is set to true if optimization is activated
670 if (cloneData == 2) {
671 cloneData = optConst;
672 }
673
674 if (pc.hasProcessed("Range")) {
675 double rangeLo = pc.getDouble("rangeLo");
676 double rangeHi = pc.getDouble("rangeHi");
677
678 // Create range with name 'fit' with above limits on all observables
679 RooArgSet obs;
680 pdf.getObservables(data.get(), obs);
681 for (auto arg : obs) {
682 RooRealVar *rrv = dynamic_cast<RooRealVar *>(arg);
683 if (rrv)
684 rrv->setRange("fit", rangeLo, rangeHi);
685 }
686
687 // Set range name to be fitted to "fit"
688 rangeName = "fit";
689 }
690
691 // Set the fitrange attribute of th PDF, add observables ranges for plotting
692 resetFitrangeAttributes(pdf, data, baseName, rangeName, splitRange);
693
694 RooArgSet projDeps;
695 auto tmp = pc.getSet("projDepSet");
696 if (tmp) {
697 projDeps.add(*tmp);
698 }
699
700 const std::string globalObservablesSource = pc.getString("globssource", "data", false);
701 if (globalObservablesSource != "data" && globalObservablesSource != "model") {
702 std::string errMsg = "RooAbsPdf::fitTo: GlobalObservablesSource can only be \"data\" or \"model\"!";
703 oocoutE(&pdf, InputArguments) << errMsg << std::endl;
704 throw std::invalid_argument(errMsg);
705 }
706 const bool takeGlobalObservablesFromData = globalObservablesSource == "data";
707
708 // Lambda function to create the correct constraint term for a PDF. In old
709 // RooFit, we use this PDF itself as the argument, for the new BatchMode
710 // we're passing a clone.
711 auto createConstr = [&](bool removeConstraintsFromPdf = false) -> std::unique_ptr<RooAbsReal> {
712 return createConstraintTerm(baseName + "_constr", // name
713 pdf, // pdf
714 data, // data
715 pc.getSet("cPars"), // Constrain RooCmdArg
716 pc.getSet("extCons"), // ExternalConstraints RooCmdArg
717 pc.getSet("glObs"), // GlobalObservables RooCmdArg
718 pc.getString("globstag", nullptr, true), // GlobalObservablesTag RooCmdArg
719 takeGlobalObservablesFromData, // From GlobalObservablesSource RooCmdArg
720 removeConstraintsFromPdf);
721 };
722
723 auto evalBackend = static_cast<RooFit::EvalBackend::Value>(pc.getInt("EvalBackend"));
724
725 // Construct BatchModeNLL if requested
726 if (evalBackend != RooFit::EvalBackend::Value::Legacy) {
727
728 // Set the normalization range. We need to do it now, because it will be
729 // considered in `compileForNormSet`.
730 std::string oldNormRange;
731 if (pdf.normRange()) {
732 oldNormRange = pdf.normRange();
733 }
734 pdf.setNormRange(rangeName);
735
736 RooArgSet normSet;
737 pdf.getObservables(data.get(), normSet);
738 normSet.remove(projDeps, true, true);
739
740 pdf.setAttribute("SplitRange", splitRange);
741 pdf.setStringAttribute("RangeName", rangeName);
742
744 ctx.setLikelihoodMode(true);
745 std::unique_ptr<RooAbsArg> head = pdf.compileForNormSet(normSet, ctx);
746 std::unique_ptr<RooAbsPdf> pdfClone = std::unique_ptr<RooAbsPdf>{static_cast<RooAbsPdf *>(head.release())};
747
748 // reset attributes
749 pdf.setAttribute("SplitRange", false);
750 pdf.setStringAttribute("RangeName", nullptr);
751
752 // Reset the normalization range
753 pdf.setNormRange(oldNormRange.c_str());
754
755 if (addCoefRangeName) {
756 oocxcoutI(&pdf, Fitting) << "RooAbsPdf::fitTo(" << pdf.GetName()
757 << ") fixing interpretation of coefficients of any component to range "
758 << addCoefRangeName << "\n";
759 pdfClone->fixAddCoefRange(addCoefRangeName, false);
760 }
761
762 std::unique_ptr<RooAbsReal> compiledConstr;
763 if (std::unique_ptr<RooAbsReal> constr = createConstr()) {
764 compiledConstr = RooFit::Detail::compileForNormSet(*constr, *data.get());
765 compiledConstr->addOwnedComponents(std::move(constr));
766 }
767
768 auto nll = createNLLNew(*pdfClone, data, std::move(compiledConstr), rangeName ? rangeName : "", projDeps, ext,
769 pc.getDouble("IntegrateBins"), offset);
770
771 std::unique_ptr<RooAbsReal> nllWrapper;
772
773 if (evalBackend == RooFit::EvalBackend::Value::Codegen ||
775 bool createGradient = evalBackend == RooFit::EvalBackend::Value::Codegen;
776 auto simPdf = dynamic_cast<RooSimultaneous const *>(pdfClone.get());
777 nllWrapper = std::make_unique<RooFit::Experimental::RooFuncWrapper>("nll_func_wrapper", "nll_func_wrapper",
778 *nll, &data, simPdf, createGradient);
779 if (createGradient)
780 static_cast<Experimental::RooFuncWrapper &>(*nllWrapper).createGradient();
781 } else {
782 nllWrapper = std::make_unique<RooEvaluatorWrapper>(
783 *nll, &data, evalBackend == RooFit::EvalBackend::Value::Cuda, rangeName ? rangeName : "", pdfClone.get(),
784 takeGlobalObservablesFromData);
785 }
786
787 nllWrapper->addOwnedComponents(std::move(nll));
788 nllWrapper->addOwnedComponents(std::move(pdfClone));
789 return nllWrapper;
790 }
791
792 std::unique_ptr<RooAbsReal> nll;
793
794#ifdef ROOFIT_LEGACY_EVAL_BACKEND
795 bool verbose = pc.getInt("verbose");
796
797 int numcpu = pc.getInt("numcpu");
798 int numcpu_strategy = pc.getInt("interleave");
799 // strategy 3 works only for RooSimultaneous.
800 if (numcpu_strategy == 3 && !pdf.InheritsFrom("RooSimultaneous")) {
801 oocoutW(&pdf, Minimization) << "Cannot use a NumCpu Strategy = 3 when the pdf is not a RooSimultaneous, "
802 "falling back to default strategy = 0"
803 << std::endl;
804 numcpu_strategy = 0;
805 }
806 RooFit::MPSplit interl = (RooFit::MPSplit)numcpu_strategy;
807
808 auto binnedLInfo = RooHelpers::getBinnedL(pdf);
809 RooAbsPdf &actualPdf = binnedLInfo.binnedPdf ? *binnedLInfo.binnedPdf : pdf;
810
811 // Construct NLL
814 cfg.addCoefRangeName = addCoefRangeName ? addCoefRangeName : "";
815 cfg.nCPU = numcpu;
816 cfg.interleave = interl;
817 cfg.verbose = verbose;
818 cfg.splitCutRange = static_cast<bool>(splitRange);
819 cfg.cloneInputData = static_cast<bool>(cloneData);
820 cfg.integrateOverBinsPrecision = pc.getDouble("IntegrateBins");
821 cfg.binnedL = binnedLInfo.isBinnedL;
822 cfg.takeGlobalObservablesFromData = takeGlobalObservablesFromData;
823 cfg.rangeName = rangeName ? rangeName : "";
824 auto nllVar = std::make_unique<RooNLLVar>(baseName.c_str(), "-log(likelihood)", actualPdf, data, projDeps, ext, cfg);
825 nllVar->enableBinOffsetting(offset == RooFit::OffsetMode::Bin);
826 nll = std::move(nllVar);
828
829 // Include constraints, if any, in likelihood
830 if (std::unique_ptr<RooAbsReal> constraintTerm = createConstr()) {
831
832 // Even though it is technically only required when the computation graph
833 // is changed because global observables are taken from data, it is safer
834 // to clone the constraint model in general to reset the normalization
835 // integral caches and avoid ASAN build failures (the PDF of the main
836 // measurement is cloned too anyway, so not much overhead). This can be
837 // reconsidered after the caching of normalization sets by pointer is changed
838 // to a more memory-safe solution.
839 constraintTerm = RooHelpers::cloneTreeWithSameParameters(*constraintTerm, data.get());
840
841 // Redirect the global observables to the ones from the dataset if applicable.
842 constraintTerm->setData(data, false);
843
844 // The computation graph for the constraints is very small, no need to do
845 // the tracking of clean and dirty nodes here.
846 constraintTerm->setOperMode(RooAbsArg::ADirty);
847
848 auto orignll = std::move(nll);
849 nll = std::make_unique<RooAddition>((baseName + "_with_constr").c_str(), "nllWithCons",
850 RooArgSet(*orignll, *constraintTerm));
851 nll->addOwnedComponents(std::move(orignll), std::move(constraintTerm));
852 }
853
854 if (optConst) {
855 nll->constOptimizeTestStatistic(RooAbsArg::Activate, optConst > 1);
856 }
857
859 nll->enableOffsetting(true);
860 }
861#else
862 throw std::runtime_error("RooFit was not built with the legacy evaluation backend");
863#endif
864
865 return nll;
866}
867
868std::unique_ptr<RooAbsReal> createChi2(RooAbsReal &real, RooAbsData &data, const RooLinkedList &cmdList)
869{
870#ifdef ROOFIT_LEGACY_EVAL_BACKEND
871 const bool isDataHist = dynamic_cast<RooDataHist const *>(&data);
872
873 RooCmdConfig pc("createChi2(" + std::string(real.GetName()) + ")");
874
875 pc.defineInt("numcpu", "NumCPU", 0, 1);
876 pc.defineInt("verbose", "Verbose", 0, 0);
877
879
880 if (isDataHist) {
881 // Construct Chi2
883 std::string baseName = "chi2_" + std::string(real.GetName()) + "_" + data.GetName();
884
885 // Clear possible range attributes from previous fits.
886 real.removeStringAttribute("fitrange");
887
888 pc.defineInt("etype", "DataError", 0, (Int_t)RooDataHist::Auto);
889 pc.defineInt("extended", "Extended", 0, extendedFitDefault);
890 pc.defineInt("split_range", "SplitRange", 0, 0);
891 pc.defineDouble("integrate_bins", "IntegrateBins", 0, -1);
892 pc.defineString("addCoefRange", "SumCoefRange", 0, "");
893 pc.allowUndefined();
894
895 pc.process(cmdList);
896 if (!pc.ok(true)) {
897 return nullptr;
898 }
899
900 bool extended = false;
901 if (auto pdf = dynamic_cast<RooAbsPdf const *>(&real)) {
902 extended = interpretExtendedCmdArg(*pdf, pc.getInt("extended"));
903 }
904
905 RooDataHist::ErrorType etype = static_cast<RooDataHist::ErrorType>(pc.getInt("etype"));
906
907 const char *rangeName = pc.getString("rangeName", nullptr, true);
908 const char *addCoefRangeName = pc.getString("addCoefRange", nullptr, true);
909
910 cfg.rangeName = rangeName ? rangeName : "";
911 cfg.nCPU = pc.getInt("numcpu");
913 cfg.verbose = static_cast<bool>(pc.getInt("verbose"));
914 cfg.cloneInputData = false;
915 cfg.integrateOverBinsPrecision = pc.getDouble("integrate_bins");
916 cfg.addCoefRangeName = addCoefRangeName ? addCoefRangeName : "";
917 cfg.splitCutRange = static_cast<bool>(pc.getInt("split_range"));
918 auto chi2 = std::make_unique<RooChi2Var>(baseName.c_str(), baseName.c_str(), real,
919 static_cast<RooDataHist &>(data), extended, etype, cfg);
920
922
923 return chi2;
924 } else {
925 pc.defineInt("integrate", "Integrate", 0, 0);
926 pc.defineObject("yvar", "YVar", 0, nullptr);
927 pc.defineString("rangeName", "RangeWithName", 0, "", true);
928 pc.defineInt("interleave", "NumCPU", 1, 0);
929
930 // Process and check varargs
931 pc.process(cmdList);
932 if (!pc.ok(true)) {
933 return nullptr;
934 }
935
936 // Decode command line arguments
937 bool integrate = pc.getInt("integrate");
938 RooRealVar *yvar = static_cast<RooRealVar *>(pc.getObject("yvar"));
939 const char *rangeName = pc.getString("rangeName", nullptr, true);
940 Int_t numcpu = pc.getInt("numcpu");
941 Int_t numcpu_strategy = pc.getInt("interleave");
942 // strategy 3 works only for RooSimultaneous.
943 if (numcpu_strategy == 3 && !real.InheritsFrom("RooSimultaneous")) {
944 oocoutW(&real, Minimization) << "Cannot use a NumCpu Strategy = 3 when the pdf is not a RooSimultaneous, "
945 "falling back to default strategy = 0"
946 << std::endl;
947 numcpu_strategy = 0;
948 }
949 RooFit::MPSplit interl = (RooFit::MPSplit)numcpu_strategy;
950 bool verbose = pc.getInt("verbose");
951
952 cfg.rangeName = rangeName ? rangeName : "";
953 cfg.nCPU = numcpu;
954 cfg.interleave = interl;
955 cfg.verbose = verbose;
956 cfg.verbose = false;
957
958 std::string name = "chi2_" + std::string(real.GetName()) + "_" + data.GetName();
959
960 return std::make_unique<RooXYChi2Var>(name.c_str(), name.c_str(), real, static_cast<RooDataSet &>(data), yvar,
961 integrate, cfg);
962 }
963#else
964 throw std::runtime_error("createChi2() is not supported without the legacy evaluation backend");
965 return nullptr;
966#endif
967}
968
969std::unique_ptr<RooFitResult> fitTo(RooAbsReal &real, RooAbsData &data, const RooLinkedList &cmdList, bool chi2)
970{
971 const bool isDataHist = dynamic_cast<RooDataHist const *>(&data);
972
973 RooCmdConfig pc("fitTo(" + std::string(real.GetName()) + ")");
974
975 RooLinkedList fitCmdList(cmdList);
976 std::string nllCmdListString;
977 if (!chi2) {
978 nllCmdListString = "ProjectedObservables,Extended,Range,"
979 "RangeWithName,SumCoefRange,NumCPU,SplitRange,Constrained,Constrain,ExternalConstraints,"
980 "CloneData,GlobalObservables,GlobalObservablesSource,GlobalObservablesTag,"
981 "EvalBackend,IntegrateBins,ModularL";
982
983 if (!cmdList.FindObject("ModularL") || static_cast<RooCmdArg *>(cmdList.FindObject("ModularL"))->getInt(0) == 0) {
984 nllCmdListString += ",OffsetLikelihood";
985 }
986 } else {
987 auto createChi2DataHistCmdArgs = "Range,RangeWithName,NumCPU,Optimize,IntegrateBins,ProjectedObservables,"
988 "AddCoefRange,SplitRange,DataError,Extended";
989 auto createChi2DataSetCmdArgs = "YVar,Integrate,RangeWithName,NumCPU,Verbose";
990 nllCmdListString += isDataHist ? createChi2DataHistCmdArgs : createChi2DataSetCmdArgs;
991 }
992
993 RooLinkedList nllCmdList = pc.filterCmdList(fitCmdList, nllCmdListString.c_str());
994
995 pc.defineDouble("prefit", "Prefit", 0, 0);
996 defineMinimizationOptions(pc);
997
998 // Process and check varargs
999 pc.process(fitCmdList);
1000 if (!pc.ok(true)) {
1001 return nullptr;
1002 }
1003
1004 // TimingAnalysis works only for RooSimultaneous.
1005 if (pc.getInt("timingAnalysis") && !real.InheritsFrom("RooSimultaneous")) {
1006 oocoutW(&real, Minimization) << "The timingAnalysis feature was built for minimization with RooSimultaneous "
1007 "and is not implemented for other PDF's. Please create a RooSimultaneous to "
1008 "enable this feature."
1009 << std::endl;
1010 }
1011
1012 // Decode command line arguments
1013 double prefit = pc.getDouble("prefit");
1014
1015 if (prefit != 0) {
1016 size_t nEvents = static_cast<size_t>(prefit * data.numEntries());
1017 if (prefit > 0.5 || nEvents < 100) {
1018 oocoutW(&real, InputArguments) << "PrefitDataFraction should be in suitable range."
1019 << "With the current PrefitDataFraction=" << prefit
1020 << ", the number of events would be " << nEvents << " out of "
1021 << data.numEntries() << ". Skipping prefit..." << std::endl;
1022 } else {
1023 size_t step = data.numEntries() / nEvents;
1024
1025 RooDataSet tiny("tiny", "tiny", *data.get(), data.isWeighted() ? RooFit::WeightVar() : RooCmdArg());
1026
1027 for (int i = 0; i < data.numEntries(); i += step) {
1028 const RooArgSet *event = data.get(i);
1029 tiny.add(*event, data.weight());
1030 }
1031 RooLinkedList tinyCmdList(cmdList);
1032 pc.filterCmdList(tinyCmdList, "Prefit,Hesse,Minos,Verbose,Save,Timer");
1033 RooCmdArg hesse_option = RooFit::Hesse(false);
1034 RooCmdArg print_option = RooFit::PrintLevel(-1);
1035
1036 tinyCmdList.Add(&hesse_option);
1037 tinyCmdList.Add(&print_option);
1038
1039 fitTo(real, tiny, tinyCmdList, chi2);
1040 }
1041 }
1042
1043 RooCmdArg modularL_option;
1044 if (pc.getInt("parallelize") != 0 || pc.getInt("enableParallelGradient") || pc.getInt("enableParallelDescent")) {
1045 // Set to new style likelihood if parallelization is requested
1046 modularL_option = RooFit::ModularL(true);
1047 nllCmdList.Add(&modularL_option);
1048 }
1049
1050 std::unique_ptr<RooAbsReal> nll;
1051 if (chi2) {
1052 nll = std::unique_ptr<RooAbsReal>{isDataHist ? real.createChi2(static_cast<RooDataHist &>(data), nllCmdList)
1053 : real.createChi2(static_cast<RooDataSet &>(data), nllCmdList)};
1054 } else {
1055 nll = std::unique_ptr<RooAbsReal>{dynamic_cast<RooAbsPdf &>(real).createNLL(data, nllCmdList)};
1056 }
1057
1058 return RooFit::FitHelpers::minimize(real, *nll, data, pc);
1059}
1060
1061} // namespace FitHelpers
1062} // namespace RooFit
1063
1064/// \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, bool removeConstraintsFromPdf)
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
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
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void data
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:77
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,...
RooAbsCollection * selectByAttrib(const char *name, bool value) const
Create a subset of the current collection, consisting only of those elements with the specified attri...
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...
Storage_t::size_type size() const
RooAbsArg * first() const
virtual bool addOwned(RooAbsArg &var, bool silent=false)
Add an argument and transfer the ownership to the collection.
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
virtual double expectedEvents(const RooArgSet *nset) const
Return expected number of events to be used in calculation of extended likelihood.
void setNormRange(const char *rangeName)
@ CanBeExtended
Definition RooAbsPdf.h:213
@ MustBeExtended
Definition RooAbsPdf.h:213
const char * normRange() const
Definition RooAbsPdf.h:251
virtual ExtendMode extendMode() const
Returns ability of PDF to provide extended likelihood terms.
Definition RooAbsPdf.h:217
Abstract base class for objects that represent a real value and implements functionality common to al...
Definition RooAbsReal.h:59
double getVal(const RooArgSet *normalisationSet=nullptr) const
Evaluate object.
Definition RooAbsReal.h:103
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.
virtual RooFit::OwningPtr< RooAbsReal > createChi2(RooDataHist &data, const RooLinkedList &cmdList)
RooArgList is a container object that can hold multiple RooAbsArg objects.
Definition RooArgList.h:22
RooAbsArg * at(Int_t idx) const
Return object at given index, or nullptr if index is out of range.
Definition RooArgList.h:110
RooArgSet is a container object that can hold multiple RooAbsArg objects.
Definition RooArgSet.h:55
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:86
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:57
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...
virtual void Add(TObject *arg)
TObject * FindObject(const char *name) const override
Return pointer to object with given name.
Wrapper class around ROOT::Fit:Fitter that provides a seamless interface between the minimizer functi...
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.
RooAbsPdf * getPdf(RooStringView catName) const
Return the p.d.f associated with the given index category name.
void wrapPdfsInBinSamplingPdfs(RooAbsData const &data, double precision)
Wraps the components of this RooSimultaneous in RooBinSamplingPdfs.
const RooAbsCategoryLValue & indexCat() const
Int_t GetNrows() const
TMatrixTSym< Element > & Similarity(const TMatrixT< Element > &n)
Calculate B * (*this) * B^T , final matrix will be (nrowsb x nrowsb) This is a similarity transform w...
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:525
RooCmdArg WeightVar(const char *name="weight", bool reinterpretAsWeight=false)
RooCmdArg Hesse(bool flag=true)
RooCmdArg ModularL(bool flag=false)
RooCmdArg PrintLevel(Int_t code)
RVec< PromoteType< T > > log(const RVec< T > &v)
Definition RVec.hxx:1841
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:353
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