Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
xRooNode.cxx
Go to the documentation of this file.
1/*
2 * Project: xRooFit
3 * Author:
4 * Will Buttinger, RAL 2022
5 *
6 * Copyright (c) 2022, CERN
7 *
8 * Redistribution and use in source and binary forms,
9 * with or without modification, are permitted according to the terms
10 * listed in LICENSE (http://roofit.sourceforge.net/license.txt)
11 */
12
13/** \class ROOT::Experimental::XRooFit::xRooNode
14\ingroup xroofit
15
16The xRooNode class is designed to wrap over a TObject and provide functionality to aid with interacting with that
17object, particularly in the case where the object is a RooFit class instance. It is a smart pointer to the object, so
18you have access to all the methods of the object too.
19
20xRooNode is designed to work in both python and C++, but examples below are given in python because that is imagined
21 be the most common way to use the xRooFit API.
22
23-# [Exploring workspaces](\ref exploring-workspaces)
24
25\anchor exploring-workspaces
26## Exploring workspaces
27
28An existing workspace file (either a ROOT file containing a RooWorkspace, or a json HS3 file) can be opened using
29 xRooNode like this:
30
31\code{.py}
32from ROOT.Experimental import XRooFit
33w = XRooFit.xRooNode("workspace.root") # or can use workspace.json for HS3
34\endcode
35
36 You can explore the content of the workspace somewhat like you would a file system: each node contains sub-nodes,
37 which you can interact with to explore ever deeper. The most relevant methods for navigating the workspace and
38exploring the content are:
39
40
41
42
43 */
44
45#include "RVersion.h"
46
47#if ROOT_VERSION_CODE < ROOT_VERSION(6, 27, 00)
48
49#define protected public
50#include "TRootBrowser.h"
52#define private public
53#include "RooAbsArg.h"
54#include "RooWorkspace.h"
57#include "RooProdPdf.h"
58#include "TGFileBrowser.h"
59#include "RooFitResult.h"
60#include "TPad.h"
61#undef private
62#include "RooAddPdf.h"
63#include "RooRealSumPdf.h"
64#include "RooProduct.h"
65#include "RooHistFunc.h"
66#include "RooConstVar.h"
67#include "RooSimultaneous.h"
68#undef protected
69
70#define GETWS(a) a->_myws
71#define GETWSSETS(w) w->_namedSets
72#define GETWSSNAPSHOTS(w) w->_snapshots
73#define GETACTBROWSER(b) b->fActBrowser
74#define GETROOTDIR(b) b->fRootDir
75#define GETLISTTREE(b) b->fListTree
76#define GETDMP(o, m) o->m
77
78#else
79
80#include "RooAbsArg.h"
81#include "RooWorkspace.h"
82#include "RooFitResult.h"
83#include "RooConstVar.h"
84#include "RooHistFunc.h"
85#include "RooRealSumPdf.h"
86#include "RooSimultaneous.h"
87#include "RooAddPdf.h"
88#include "RooProduct.h"
89#include "TPad.h"
93#include "RooProdPdf.h"
94#include "TRootBrowser.h"
95#include "TGFileBrowser.h"
96#include "TF1.h"
98
100{
101 return a->workspace();
102}
104{
105 return w->sets();
106}
108{
109 return w->getSnapshots();
110}
112{
113 return b->GetActBrowser();
114}
116{
117 return b->GetRootDir();
118}
120{
121 return b->GetListTree();
122}
123#define GETDMP(o, m) \
124 *reinterpret_cast<void **>(reinterpret_cast<unsigned char *>(o) + o->Class()->GetDataMemberOffset(#m))
125
126#endif
127
128#include "RooAddition.h"
129
130#include "RooCategory.h"
131#include "RooRealVar.h"
132#include "RooStringVar.h"
133#include "RooBinning.h"
134#include "RooUniformBinning.h"
135
136#include "RooAbsData.h"
137#include "RooDataHist.h"
138#include "RooDataSet.h"
139
140#include "xRooFit/xRooNode.h"
141#include "xRooFit/xRooFit.h"
142
143#include "TH1.h"
144#include "TBrowser.h"
145#include "TROOT.h"
146#include "TQObject.h"
147#include "TAxis.h"
148#include "TGraphAsymmErrors.h"
149#include "TMath.h"
150#include "TPRegexp.h"
151#include "TRegexp.h"
152#include "TExec.h"
153#include "TPaveText.h"
154
155#include "TGListTree.h"
156#include "TGMsgBox.h"
157#include "TGedEditor.h"
158#include "TGMimeTypes.h"
159#include "TH2.h"
160#include "RooExtendPdf.h"
161#include "RooExtendedBinding.h"
162
164
165#include "coutCapture.h"
166
167// #include "RooFitTrees/RooFitResultTree.h"
168// #include "RooFitTrees/RooDataTree.h"
169#include "TFile.h"
170#include "TSystem.h"
171#include "TKey.h"
172#include "TEnv.h"
173#include "TStyle.h"
174
175#if ROOT_VERSION_CODE >= ROOT_VERSION(6, 26, 00)
177#endif
178
179#if ROOT_VERSION_CODE >= ROOT_VERSION(6, 24, 00)
180#include "RooBinSamplingPdf.h"
181#endif
182
183#include "RooPoisson.h"
184#include "RooGaussian.h"
185#include "RooFormulaVar.h"
186#include "RooGenericPdf.h"
187#include "TVectorD.h"
188#include "TStopwatch.h"
189#include "TTimeStamp.h"
190
191#include <csignal>
192
193#include "TCanvas.h"
194#include "THStack.h"
195
196#include "TLegend.h"
197#include "TLegendEntry.h"
198#include "TGraphErrors.h"
199#include "TMultiGraph.h"
200#include "TFrame.h"
201#include "RooProjectedPdf.h"
202#include "TMemFile.h"
203#include "TGaxis.h"
204#include "TPie.h"
205// #include <thread>
206// #include <future>
207
208#if ROOT_VERSION_CODE < ROOT_VERSION(6, 27, 00)
209#include "RooNaNPacker.h"
210#endif
211
213
214xRooNode::InteractiveObject *xRooNode::gIntObj = nullptr;
215std::map<std::string, std::tuple<std::function<double(double, double, double)>, bool>> xRooNode::auxFunctions;
216void xRooNode::SetAuxFunction(const char *title, const std::function<double(double, double, double)> &func,
217 bool symmetrize)
218{
219 auxFunctions[title] = std::make_tuple(func, symmetrize);
220}
221
222template <typename T>
223const T &_or_func(const T &a, const T &b)
224{
225 if (a)
226 return a;
227 return b;
228}
229
230////////////////////////////////////////////////////////////////////////////////
231/// Create new object of type classname, with given name and title, and own-wrap it
232/// i.e. the xRooNode will delete the object when the node (and any that reference it) is destroyed
233///
234/// \param classname : the type of the object to create
235/// \param name : the name to give the object
236/// \param title : the title to give the object
237
238xRooNode::xRooNode(const char *classname, const char *name, const char *title)
239 : xRooNode(name, std::shared_ptr<TObject>(TClass::GetClass(classname)
240 ? reinterpret_cast<TObject *>(TClass::GetClass(classname)->New())
241 : nullptr,
242 [](TObject *o) {
243 if (auto w = dynamic_cast<RooWorkspace *>(o); w) {
244#if ROOT_VERSION_CODE < ROOT_VERSION(6, 27, 00)
245 w->_embeddedDataList.Delete();
246#endif
247 xRooNode(*w, std::make_shared<xRooNode>()).sterilize();
248 }
249 if (o)
250 delete o;
251 }))
252{
253 if (auto a = get<TNamed>(); a)
254 a->SetName(name);
255 SetTitle(title);
256}
257
258xRooNode::xRooNode(const char *name, const std::shared_ptr<TObject> &comp, const std::shared_ptr<xRooNode> &parent)
259 : TNamed(name, ""), fComp(comp), fParent(parent)
260{
261
262 if (!fComp && !fParent && name && strlen(name) > 0) {
263 char *_path = gSystem->ExpandPathName(name);
264 TString pathName = TString(_path);
265 delete[] _path;
266 if (!gSystem->AccessPathName(pathName)) {
267 // if file is json can try to read
268 if (pathName.EndsWith(".json")) {
269#if ROOT_VERSION_CODE >= ROOT_VERSION(6, 26, 00)
270 fComp = std::make_shared<RooWorkspace>("workspace", name);
271 RooJSONFactoryWSTool tool(*get<RooWorkspace>());
274 if (!tool.importJSON(pathName.Data())) {
275 Error("xRooNode", "Error reading json workspace %s", name);
276 fComp.reset();
277 }
279#else
280 Error("xRooNode", "json format workspaces available only in ROOT 6.26 onwards");
281#endif
282 } else {
283
284 // using acquire in the constructor seems to cause a mem leak according to valgrind ... possibly because
285 // (*this) gets called on it before the node is fully constructed
286 auto _file = std::make_shared<TFile>(
287 pathName); // acquire<TFile>(name); // acquire file to ensure stays open while we have the workspace
288 // actually it appears we don't need to keep the file open once we've loaded the workspace, but should be
289 // no harm doing so
290 // otherwise the workspace doesn't saveas
291 auto keys = _file->GetListOfKeys();
292 if (keys) {
293 for (auto &&k : *keys) {
294 auto cl = TClass::GetClass((static_cast<TKey *>(k))->GetClassName());
295 if (cl == RooWorkspace::Class() || cl->InheritsFrom("RooWorkspace")) {
296 fComp.reset(_file->Get<RooWorkspace>(k->GetName()), [](TObject *ws) {
297 // memory leak in workspace, some RooLinkedLists aren't cleared, fixed in ROOT 6.28
298 if (ws) {
299#if ROOT_VERSION_CODE < ROOT_VERSION(6, 27, 00)
300 dynamic_cast<RooWorkspace *>(ws)->_embeddedDataList.Delete();
301#endif
302 xRooNode(*ws, std::make_shared<xRooNode>()).sterilize();
303 delete ws;
304 }
305 });
306 if (fComp) {
307 TNamed::SetNameTitle(fComp->GetName(), fComp->GetTitle());
308 fParent = std::make_shared<xRooNode>(
309 _file); // keep file alive - seems necessary to save workspace again in some cases
310 break;
311 }
312 }
313 }
314 }
315 }
316 } else if (pathName.EndsWith(".root") || pathName.EndsWith(".json")) {
317 throw std::runtime_error(TString::Format("%s does not exist", name));
318 }
319 }
320
321 if (auto _ws = get<RooWorkspace>(); _ws && (!parent || parent->get<TFile>())) {
324 .removeTopic(RooFit::NumIntegration); // stop info message every time
325
326 // check if any of the open files have version numbers greater than our major version
327 // may not read correctly
328 for (auto f : *gROOT->GetListOfFiles()) {
329 if ((dynamic_cast<TFile *>(f)->GetVersion() / 100) > (gROOT->GetVersionInt() / 100)) {
330 Warning("xRooNode", "There is file open with version %d > current version %d ... results may be wrong",
331 dynamic_cast<TFile *>(f)->GetVersion(), gROOT->GetVersionInt());
332 }
333 }
334
335 // load list of colors if there is one
336 if (auto colors = dynamic_cast<TSeqCollection *>(_ws->obj(gROOT->GetListOfColors()->GetName()))) {
337 gROOT->GetListOfColors()->Clear();
338 for (auto col : *colors) {
339 gROOT->GetListOfColors()->Add(col);
340 }
341 }
342
343 // use the datasets if any to 'mark' observables
344 int checkCount = 0;
345 for (auto &d : _ws->allData()) {
346 for (auto &a : *d->get()) {
347 if (auto v = _ws->var(a->GetName()); v) {
348 v->setAttribute("obs");
349 } else if (auto c = _ws->cat(a->GetName()); c) {
350 c->setAttribute("obs");
351 }
352 }
353 // count how many ds are checked ... if none are checked will check the first
354 checkCount += d->TestBit(1 << 20);
355 }
356
357 if (checkCount == 0 && !_ws->allData().empty())
358 _ws->allData().back()->SetBit(1 << 20, true);
359
360 if (auto _set = dynamic_cast<RooArgSet *>(GETWSSNAPSHOTS(_ws).find("NominalParamValues")); _set) {
361 for (auto s : *_set) {
362 if (auto v = dynamic_cast<RooRealVar *>(s); v) {
363 _ws->var(s->GetName())->setStringAttribute("nominal", TString::Format("%f", v->getVal()));
364 }
365 }
366 }
367
368 // also flag global observables ... relies on ModelConfig existences
369 RooArgSet _allGlobs;
370 for (auto &[k, v] : GETWSSETS(_ws)) {
371 if (k == "globalObservables" || TString(k).EndsWith("_GlobalObservables")) {
372 for (auto &s : v) {
373 _allGlobs.add(*s);
374 s->setAttribute("obs");
375 s->setAttribute("global");
376 }
377 } else if (TString(k).EndsWith("_Observables")) {
378 const_cast<RooArgSet &>(v).setAttribAll("obs");
379 } else if (TString(k).EndsWith("_POI")) {
380 for (auto &s : v) {
381 s->setAttribute("poi");
382 auto _v = dynamic_cast<RooRealVar *>(s);
383 if (!_v)
384 continue;
385 // if (!_v->hasRange("physical")) {
386 // _v->setRange("physical", 0, std::numeric_limits<double>::infinity());
387 // // ensure range of poi is also straddling 0
388 // if (_v->getMin() >= 0)
389 // _v->setMin(-1e-5);
390 // }
391 }
392 } else if (TString(k).EndsWith("_NuisParams")) {
393 const_cast<RooArgSet &>(v).setAttribAll("np");
394 }
395 }
396 if (!_allGlobs.empty() && GETWSSETS(_ws).count("globalObservables") == 0) {
397 _ws->defineSet("globalObservables", _allGlobs);
398 }
399
400 // now check if any pars don't have errors defined (not same as error=0) ... if so, use the first pdf (if there is
401 // one) to try setting values from
402 if (!_ws->allPdfs().empty()) {
403 std::set<RooRealVar *> noErrorPars;
404 std::string parNames;
405 for (auto &p : np()) { // infer errors on all floating non-poi parameters
406 auto v = p->get<RooRealVar>();
407 if (!v)
408 continue;
409 if (!v->hasError()) {
410 noErrorPars.insert(v);
411 if (!parNames.empty())
412 parNames += ",";
413 parNames += v->GetName();
414 }
415 }
416 if (!noErrorPars.empty()) {
417 Warning("xRooNode",
418 "Inferring initial errors of %d parameters (%s%s) (give all nuisance parameters an error to avoid "
419 "this msg)",
420 int(noErrorPars.size()), (*noErrorPars.begin())->GetName(), (noErrorPars.size() > 1) ? ",..." : "");
421 // get the first top-level pdf
422 browse();
423 for (auto &a : *this) {
424 if (noErrorPars.empty()) {
425 break;
426 }
427 if (a->fFolder == "!pdfs") {
428 try {
429 auto fr = a->floats().reduced(parNames).fitResult("prefit");
430 if (auto _fr = fr.get<RooFitResult>(); _fr) {
431 std::set<RooRealVar *> foundPars;
432 for (auto &v : noErrorPars) {
433 if (auto arg = dynamic_cast<RooRealVar *>(_fr->floatParsFinal().find(v->GetName()));
434 arg && arg->hasError()) {
435 v->setError(arg->getError());
436 foundPars.insert(v);
437 }
438 }
439 for (auto &v : foundPars) {
440 noErrorPars.erase(v);
441 }
442 }
443 } catch (...) {
444 }
445 }
446 }
447 }
448 }
449 }
450
451 if (strlen(GetTitle()) == 0) {
452 if (fComp) {
453 TNamed::SetTitle(fComp->GetTitle());
454 } else {
455 TNamed::SetTitle(GetName());
456 }
457 }
458}
459
460xRooNode::xRooNode(const TObject &comp, const std::shared_ptr<xRooNode> &parent)
461 : xRooNode(/*[](const TObject& c) {
462c.InheritsFrom("RooAbsArg");
463if (s) {
464return (s->getStringAttribute("alias")) ? s->getStringAttribute("alias") : c.GetName();
465}
466return c.GetName();
467}(comp)*/
468 (comp.InheritsFrom("RooAbsArg") && dynamic_cast<const RooAbsArg *>(&comp)->getStringAttribute("alias"))
469 ? dynamic_cast<const RooAbsArg *>(&comp)->getStringAttribute("alias")
470 : comp.GetName(),
471 std::shared_ptr<TObject>(const_cast<TObject *>(&comp), [](TObject *) {}), parent)
472{
473}
474
475xRooNode::xRooNode(const std::shared_ptr<TObject> &comp, const std::shared_ptr<xRooNode> &parent)
476 : xRooNode(
477 [&]() {
478 if (auto a = std::dynamic_pointer_cast<RooAbsArg>(comp); a && a->getStringAttribute("alias"))
479 return a->getStringAttribute("alias");
480 if (comp)
481 return comp->GetName();
482 return "";
483 }(),
484 comp, parent)
485{
486}
487
489
490void xRooNode::Checked(TObject *obj, bool val)
491{
492 if (obj != this)
493 return;
494
495 // cycle through states:
496 // unhidden and selected: tick, no uline
497 // hidden and unselected: notick, uline
498 // unhidden and unselected: tick, uline
499 if (auto o = get<RooAbsReal>(); o) {
500 if (o->isSelectedComp() && !val) {
501 // deselecting and hiding
502 o->selectComp(val);
503 o->setAttribute("hidden");
504 } else if (!o->isSelectedComp() && !val) {
505 // selecting
506 o->selectComp(!val);
507 } else if (val) {
508 // unhiding but keeping unselected
509 o->setAttribute("hidden", false);
510 }
511 auto item = GetTreeItem(nullptr);
512 item->CheckItem(!o->getAttribute("hidden"));
513 if (o->isSelectedComp()) {
514 item->ClearColor();
515 } else {
516 item->SetColor(kGray);
517 }
518 return;
519 }
520
521 if (auto o = get(); o) {
522 // if (o->TestBit(1<<20)==val) return; // do nothing
523 o->SetBit(1 << 20, val); // TODO: check is 20th bit ok to play with?
524 if (auto fr = get<RooFitResult>(); fr) {
525 if (auto _ws = ws(); _ws) {
526 if (val) {
527 // ensure fit result is in genericObjects list ... if not, add a copy ...
528 if (!_ws->genobj(fr->GetName())) {
529 _ws->import(*fr);
530 if (auto wfr = dynamic_cast<RooFitResult *>(_ws->genobj(fr->GetName()))) {
531 fr = wfr;
532 }
533 }
534 RooArgSet _allVars = _ws->allVars();
535 _allVars = fr->floatParsFinal();
536 _allVars = fr->constPars();
537 for (auto &i : fr->floatParsInit()) {
538 auto v = dynamic_cast<RooRealVar *>(_allVars.find(i->GetName()));
539 if (v)
540 v->setStringAttribute("initVal", TString::Format("%f", dynamic_cast<RooRealVar *>(i)->getVal()));
541 }
542 // uncheck all other fit results
543 for (auto oo : _ws->allGenericObjects()) {
544 if (auto ffr = dynamic_cast<RooFitResult *>(oo); ffr && ffr != fr) {
545 ffr->ResetBit(1 << 20);
546 }
547 }
548 } else
549 _ws->allVars() = fr->floatParsInit();
550 }
551 if (auto item = GetTreeItem(nullptr); item) {
552 // update check marks on siblings
553 if (auto first = item->GetParent()->GetFirstChild()) {
554 do {
555 if (first->HasCheckBox()) {
556 auto _obj = static_cast<xRooNode *>(first->GetUserData());
557 first->CheckItem(_obj->get() && _obj->get()->TestBit(1 << 20));
558 }
559 } while ((first = first->GetNextSibling()));
560 }
561 }
562 }
563 }
564}
565
567{
568 static bool blockBrowse = false;
569 if (blockBrowse)
570 return;
571 if (b == nullptr) {
572 auto b2 = dynamic_cast<TBrowser *>(gROOT->GetListOfBrowsers()->Last());
573 if (!b2 || !b2->GetBrowserImp()) { // no browser imp if browser was closed
574 blockBrowse = true;
575 gEnv->SetValue("X11.UseXft", "no"); // for faster x11
576 gEnv->SetValue("X11.Sync", "no");
577 gEnv->SetValue("X11.FindBestVisual", "no");
578 gEnv->SetValue("Browser.Name", "TRootBrowser"); // forces classic root browser (in 6.26 onwards)
579 gEnv->SetValue("Canvas.Name", "TRootCanvas");
580 b2 = new TBrowser("nodeBrowser", this, "RooFit Browser");
581 blockBrowse = false;
582 } else if (strcmp(b2->GetName(), "nodeBrowser") == 0) {
583 blockBrowse = true;
584 b2->BrowseObject(this);
585 blockBrowse = false;
586 } else {
587 auto _b = dynamic_cast<TGFileBrowser *>(GETACTBROWSER(dynamic_cast<TRootBrowser *>(b2->GetBrowserImp())));
588 if (_b)
589 _b->AddFSDirectory("Workspaces", nullptr, "SetRootDir");
590 /*auto l = Node2::Class()->GetMenuList();
591 auto o = new CustomClassMenuItem(TClassMenuItem::kPopupUserFunction,Node2::Class(),
592 "blah blah blah","BlahBlah",0,"Option_t*",-1,true);
593 //o->SetCall(o,"BlahBlah","Option_t*",-1);
594 l->AddFirst(o);*/
595 // b->BrowseObject(this);
596 _b->GotoDir(nullptr);
597 _b->Add(this, GetName());
598 // b->Add(this);
599 }
600 return;
601 }
602
603 if (auto item = GetTreeItem(b); item) {
604 if (!item->IsOpen() && IsFolder())
605 return; // no need to rebrowse if closing
606 // update check marks on any child items
607 if (auto first = item->GetFirstChild()) {
608 do {
609 if (first->HasCheckBox()) {
610 auto _obj = static_cast<xRooNode *>(first->GetUserData());
611 first->CheckItem(_obj->get() &&
612 (_obj->get()->TestBit(1 << 20) ||
613 (_obj->get<RooAbsArg>() && !_obj->get<RooAbsArg>()->getAttribute("hidden"))));
614 }
615 } while ((first = first->GetNextSibling()));
616 }
617 }
618
619 browse();
620
621 // for top-level pdfs default to having the .vars browsable too
622 if (get<RooAbsPdf>() && fFolder == "!pdfs" && !_IsShowVars_()) {
623 fBrowsables.push_back(std::make_shared<xRooNode>(vars()));
624 }
625
626 if (auto _fr = get<RooFitResult>(); _fr && fBrowsables.empty()) {
627 // have some common drawing options
628 fBrowsables.push_back(std::make_shared<xRooNode>(".Draw(\"pull\")", nullptr, *this));
629 fBrowsables.push_back(std::make_shared<xRooNode>(".Draw(\"corrcolztext\")", nullptr, *this));
630 if (std::unique_ptr<RooAbsCollection>(_fr->floatParsFinal().selectByAttrib("poi", true))->size() == 1) {
631 fBrowsables.push_back(std::make_shared<xRooNode>(".Draw(\"impact\")", nullptr, *this));
632 }
633 }
634
635 if (empty() && fBrowsables.empty()) {
636 try {
637 if (auto s = get<TStyle>()) {
638 s->SetFillAttributes();
639 if (auto ed = dynamic_cast<TGedEditor *>(TVirtualPadEditor::GetPadEditor())) {
640 ed->SetModel(gPad, s, kButton1Down, true);
641 }
642 } else if (TString(GetName()).BeginsWith(".Draw(\"") && fParent) {
643 fParent->Draw(TString(TString(GetName())(7, strlen(GetName()) - 9)) + b->GetDrawOption());
644 } else {
645 Draw(b->GetDrawOption());
646 }
647 } catch (const std::exception &e) {
648 new TGMsgBox(
649 gClient->GetRoot(),
650 (gROOT->GetListOfBrowsers()->At(0))
651 ? dynamic_cast<TGWindow *>(static_cast<TBrowser *>(gROOT->GetListOfBrowsers()->At(0))->GetBrowserImp())
652 : gClient->GetRoot(),
653 "Exception", e.what(),
654 kMBIconExclamation); // deletes self on dismiss?
655 }
656 }
657
658 bool hasFolders = false;
659 if (strlen(GetName()) > 0 && GetName()[0] != '!') { // folders don't have folders
660 for (auto &c : *this) {
661 if (!c->fFolder.empty()) {
662 hasFolders = true;
663 break;
664 }
665 }
666 }
667 // auto _ws = get<RooWorkspace>();
668 if (/*_ws*/ hasFolders) {
669 // organize in folders
670 auto _folders = find(".folders");
671 if (!_folders) {
672 _folders = emplace_back(std::make_shared<xRooNode>(".folders", nullptr, *this));
673 }
674 // ensure entry in folders for every folder type ...
675 for (auto &v : *this) {
676 if (!v->fFolder.empty() && !_folders->find(v->fFolder, false)) {
677 _folders->emplace_back(std::make_shared<xRooNode>(v->fFolder.c_str(), nullptr, *this));
678 }
679 }
680 // now just add all the folders
681 for (auto &v : *_folders) {
682 TString _name = v->GetName();
683 if (_name.BeginsWith('!'))
684 _name = _name(1, _name.Length()); // strip ! from display
685 b->Add(v.get(), _name);
686 }
687 }
688
689 for (auto &v : *this) {
690 if (hasFolders && !v->fFolder.empty())
691 continue; // in the folders
692 if (strcmp(v->GetName(), ".folders") == 0)
693 continue; // never 'browse' the folders property
694 auto _fr = v->get<RooFitResult>();
695 int _checked = (v->get<RooAbsData>() || _fr) ? v->get()->TestBit(1 << 20) : -1;
696 if (_fr && ((_fr->status() == 0 && _fr->numStatusHistory() == 0) || (_fr->floatParsFinal().empty()))) {
697 // this is a "PARTIAL" fit result ... don't allow it to be selected
698 _checked = -1;
699 }
700 if (v->get<RooAbsPdf>() && get<RooSimultaneous>())
701 _checked = !v->get<RooAbsArg>()->getAttribute("hidden");
702 TString _name = v->GetName();
703 if (v->get() && _name.BeginsWith(TString(v->get()->ClassName()) + "::")) {
704 _name = _name(strlen(v->get()->ClassName()) + 2, _name.Length());
705 }
706 if (_name.BeginsWith(".")) {
707 // property node -- display the name of the contained object
708 if (v->get()) {
709 _name = TString::Format("%s: %s::%s", _name.Data(), v->get()->ClassName(),
710 (v->get<RooAbsArg>() && v->get<RooAbsArg>()->getStringAttribute("alias"))
711 ? v->get<RooAbsArg>()->getStringAttribute("alias")
712 : v->get()->GetName());
713 }
714 } else if (v->get() && !v->get<TFile>() && !TString(v->GetName()).BeginsWith('/'))
715 _name = TString::Format("%s::%s", v->get()->ClassName(), _name.Data());
716 if (auto _type = v->GetNodeType(); strlen(_type)) {
717 // decided not to show const values until figure out how to update if value changes
718 /*if (TString(_type)=="Const") _name += TString::Format(" [%s=%g]",_type,v->get<RooConstVar>()->getVal());
719 else*/
720 _name += TString::Format(" [%s]", _type);
721 }
722 if (auto fv = v->get<RooFormulaVar>()) {
723 TString formu = TString::Format(" [%s]", fv->expression());
724 for (size_t i = 0; i < fv->dependents().size(); i++) {
725 formu.ReplaceAll(TString::Format("x[%zu]", i), fv->dependents()[i].GetName());
726 }
727 _name += formu;
728 } else if (auto gv = v->get<RooGenericPdf>()) {
729 TString formu = TString::Format(" [%s]", gv->expression());
730 for (size_t i = 0; i < gv->dependents().size(); i++) {
731 formu.ReplaceAll(TString::Format("x[%zu]", i), gv->dependents()[i].GetName());
732 }
733 _name += formu;
734 }
735 // tool tip defaults to displaying name and title, so temporarily set name to obj name if has one
736 // and set title to the object type
737 TString nameSave(v->TNamed::GetName());
738 TString titleSave(v->TNamed::GetTitle());
739 if (auto o = v->get(); o)
740 v->TNamed::SetNameTitle(o->GetName(), o->ClassName());
741 b->Add(v.get(), _name, _checked);
742 if (auto o = v->get(); o)
743 v->TNamed::SetNameTitle(nameSave, titleSave);
744 if (_checked != -1) {
745 dynamic_cast<TQObject *>(b->GetBrowserImp())
746 ->Connect("Checked(TObject *, bool)", ClassName(), v.get(), "Checked(TObject *, bool)");
747 }
748 if (_fr) {
749 if (_fr->status() || _fr->covQual() != 3) { // snapshots or bad fits
750 v->GetTreeItem(b)->SetColor((_fr->numStatusHistory() || _fr->floatParsFinal().empty()) ? kRed : kBlue);
751 } else if (_fr->numStatusHistory() == 0) { // partial fit result ..
752 v->GetTreeItem(b)->SetColor(kGray);
753 }
754 }
755 if ((v->fFolder == "!np" || v->fFolder == "!poi")) {
756 if (v->get<RooAbsArg>()->getAttribute("Constant")) {
757 v->GetTreeItem(b)->SetColor(kGray);
758 } else
759 v->GetTreeItem(b)->ClearColor();
760 }
761 if (auto _htr = v->get<RooStats::HypoTestResult>(); _htr) {
762 // check for fit statuses
763 if (auto fits = _htr->GetFitInfo()) {
764 for (int i = 0; i < fits->numEntries(); i++) {
765 // if any fit (other than a genFit) is bad, flag point as bad
766 if (fits->get(i)->getCatIndex("type") != 5 && fits->get(i)->getRealValue("status") != 0) {
767 v->GetTreeItem(b)->SetColor(kRed);
768 break;
769 }
770 }
771 } else {
772 v->GetTreeItem(b)->SetColor(kBlue); // unknown fit status
773 }
774 }
775
776 // v.fBrowsers.insert(b);
777 }
778
779 // for pdfs, check for datasets too and add to list
780 /*if (get<RooAbsPdf>()) {
781 auto dsets = datasets();
782 if (!dsets.empty()) {
783 // check if already have .datasets() in browsables
784 bool found(false);
785 for(auto& p : fBrowsables) {
786 if (TString(p->GetName())==".datasets()") {found=true;
787 // add
788 break;
789 }
790 }
791 if (!found) {
792 fBrowsables.push_back(std::make_shared<xRooNode>(dsets));
793 }
794 }
795 }*/
796 // browse the browsables too
797 for (auto &v : fBrowsables) {
798 TString _name = v->GetName();
799 if (_name == ".memory")
800 continue; // hide the memory from browsing, if put in browsables
801 TString nameSave(v->TNamed::GetName());
802 TString titleSave(v->TNamed::GetTitle());
803 if (auto o = v->get(); o)
804 v->TNamed::SetNameTitle(o->GetName(), o->ClassName());
805 b->Add(v.get(), _name, -1);
806 if (auto o = v->get(); o)
807 v->TNamed::SetNameTitle(nameSave, titleSave);
808 }
809
810 b->SetSelected(this);
811}
812
814{
815 if (!set) {
816 // can't remove as causes a crash, need to remove from the browser first
817 /*for(auto itr = fBrowsables.begin(); itr != fBrowsables.end(); ++itr) {
818 if (strcmp((*itr)->GetName(),".vars")==0) {
819 fBrowsables.erase(itr);
820 }
821 }*/
822 } else {
823 auto v = std::make_shared<xRooNode>(vars());
824 fBrowsables.push_back(v);
825 if (auto l = GetListTree(nullptr)) {
826 l->AddItem(GetTreeItem(nullptr), v->GetName(), v.get());
827 }
828 }
829}
830
832{
833 for (auto &b : fBrowsables) {
834 if (strcmp(b->GetName(), ".vars") == 0)
835 return true;
836 }
837 return false;
838}
839
841{
842 if (strlen(GetName()) > 0 && GetName()[0] == '!')
843 return true;
844 if (strlen(GetName()) > 0 && GetName()[0] == '.' && !(TString(GetName()).BeginsWith(".Draw(\"")))
845 return true;
846 if (empty())
847 const_cast<xRooNode *>(this)->browse();
848 return !empty();
849}
850
851class Axis2 : public TAxis {
852
853public:
854 using TAxis::TAxis;
855 double GetBinWidth(Int_t bin) const override
856 {
857 if (auto v = var(); v)
858 return v->getBinWidth(bin - 1, GetName());
859 return 1;
860 }
861 double GetBinLowEdge(Int_t bin) const override
862 {
863 if (auto v = rvar(); v) {
864 return (bin == v->getBinning(GetName()).numBins() + 1) ? v->getBinning(GetName()).binHigh(bin - 2)
865 : v->getBinning(GetName()).binLow(bin - 1);
866 }
867 return bin - 1;
868 }
869 double GetBinUpEdge(Int_t bin) const override
870 {
871 if (auto v = rvar(); v)
872 return (bin == 0) ? v->getBinning(GetName()).binLow(bin) : v->getBinning(GetName()).binHigh(bin - 1);
873 return bin;
874 }
875
876 const char *GetTitle() const override
877 {
878 return (binning() && strlen(binning()->GetTitle())) ? binning()->GetTitle() : GetParent()->GetTitle();
879 }
880 void SetTitle(const char *title) override
881 {
882 if (binning()) {
883 const_cast<RooAbsBinning *>(binning())->SetTitle(title);
884 } else {
885 dynamic_cast<TNamed *>(GetParent())->SetTitle(title);
886 }
887 }
888
889 void Set(Int_t nbins, const double *xbins) override
890 {
891 if (auto v = dynamic_cast<RooRealVar *>(rvar()))
892 v->setBinning(RooBinning(nbins, xbins), GetName());
893 TAxis::Set(nbins, xbins);
894 }
895 void Set(Int_t nbins, const float *xbins) override
896 {
897 std::vector<double> bins(nbins + 1);
898 for (int i = 0; i <= nbins; i++)
899 bins.at(i) = xbins[i];
900 return Set(nbins, &bins[0]);
901 }
902 void Set(Int_t nbins, double xmin, double xmax) override
903 {
904 if (auto v = dynamic_cast<RooRealVar *>(rvar()))
905 v->setBinning(RooUniformBinning(xmin, xmax, nbins), GetName());
906 TAxis::Set(nbins, xmin, xmax);
907 }
908
909 const RooAbsBinning *binning() const { return var()->getBinningPtr(GetName()); }
910
911 Int_t FindFixBin(const char *label) const override { return TAxis::FindFixBin(label); }
912 Int_t FindFixBin(double x) const override { return (binning()) ? (binning()->binNumber(x) + 1) : x; }
913
914private:
915 RooAbsLValue *var() const { return dynamic_cast<RooAbsLValue *>(GetParent()); }
916 RooAbsRealLValue *rvar() const { return dynamic_cast<RooAbsRealLValue *>(GetParent()); }
917};
918
919std::shared_ptr<TObject> xRooNode::getObject(const std::string &name, const std::string &type) const
920{
921 // if (fParent) return fParent->getObject(name);
922
923 if (auto _owned = find(".memory"); _owned) {
924 for (auto &o : *_owned) {
925 if (name == o->GetName()) {
926 if (type.empty() || o->get()->InheritsFrom(type.c_str()))
927 return o->fComp;
928 }
929 }
930 }
931
932 // see if have a provider
933 auto _provider = fProvider;
934 auto _parent = fParent;
935 while (!_provider && _parent) {
936 _provider = _parent->fProvider;
937 _parent = _parent->fParent;
938 }
939 if (_provider)
940 return _provider->getObject(name, type);
941
942 if (ws()) {
943 std::shared_ptr<TObject> out;
944 if (auto arg = ws()->arg(name.c_str()); arg) {
945 auto _tmp = std::shared_ptr<TObject>(arg, [](TObject *) {});
946 if (!type.empty() && arg->InheritsFrom(type.c_str()))
947 return _tmp;
948 if (!out)
949 out = _tmp;
950 }
951 if (auto arg = ws()->data(name.c_str()); arg) {
952 auto _tmp = std::shared_ptr<TObject>(arg, [](TObject *) {});
953 if (!type.empty() && arg->InheritsFrom(type.c_str()))
954 return _tmp;
955 if (!out)
956 out = _tmp;
957 }
958 if (auto arg = ws()->genobj(name.c_str()); arg) {
959 auto _tmp = std::shared_ptr<TObject>(arg, [](TObject *) {});
960 if (!type.empty() && arg->InheritsFrom(type.c_str()))
961 return _tmp;
962 if (!out)
963 out = _tmp;
964 }
965 if (auto arg = ws()->embeddedData(name.c_str()); arg) {
966 auto _tmp = std::shared_ptr<TObject>(arg, [](TObject *) {});
967 if (!type.empty() && arg->InheritsFrom(type.c_str()))
968 return _tmp;
969 if (!out)
970 out = _tmp;
971 }
972 if (auto arg = GETWSSNAPSHOTS(ws()).find(name.c_str()); arg) {
973 auto _tmp = std::shared_ptr<TObject>(arg, [](TObject *) {});
974 if (!type.empty() && arg->InheritsFrom(type.c_str()))
975 return _tmp;
976 if (!out)
977 out = _tmp;
978 }
979 return out;
980 }
981 if (auto arg = get<RooAbsArg>()) {
982 // can try all nodes
983 RooArgSet nodes;
984 arg->treeNodeServerList(&nodes);
985 if (auto server = nodes.find(name.c_str())) {
986 return std::shared_ptr<TObject>(server, [](TObject *) {});
987 }
988 }
989 return nullptr;
990}
991
993{
994 if (fXAxis) {
995 // check if num bins needs update or not
996 if (auto cat = dynamic_cast<RooAbsCategory *>(fXAxis->GetParent());
997 cat && cat->numTypes() != fXAxis->GetNbins()) {
998 fXAxis.reset();
999 } else {
1000 return fXAxis.get();
1001 }
1002 }
1003 RooAbsLValue *x = nullptr;
1004 if (auto a = get<RooAbsArg>(); a && a->isFundamental())
1005 x = dynamic_cast<RooAbsLValue *>(a); // self-axis
1006
1007 auto _parentX = (!x && fParent && !fParent->get<RooSimultaneous>()) ? fParent->GetXaxis() : nullptr;
1008
1009 auto o = get<RooAbsReal>();
1010 if (!o)
1011 return _parentX;
1012
1013 if (auto xName = o->getStringAttribute("xvar"); xName) {
1014 x = dynamic_cast<RooAbsLValue *>(getObject(xName).get());
1015 }
1016
1017 // if xvar has become set equal to an arg and this is a pdf, we will allow a do-over
1018 if (!x) {
1019 // need to choose from dependent fundamentals, in following order:
1020 // parentX (if not a glob), robs, globs, vars, args
1021
1022 if (_parentX && !dynamic_cast<RooAbsArg *>(_parentX->GetParent())->getAttribute("global") &&
1023 (o->dependsOn(*dynamic_cast<RooAbsArg *>(_parentX->GetParent())) || vars().empty())) {
1024 x = dynamic_cast<RooAbsLValue *>(_parentX->GetParent());
1025 } else if (auto _obs = obs(); !_obs.empty()) {
1026 for (auto &v : _obs) {
1027 if (!v->get<RooAbsArg>()->getAttribute("global")) {
1028 x = v->get<RooAbsLValue>();
1029 if (x)
1030 break;
1031 } else if (!x) {
1032 x = v->get<RooAbsLValue>();
1033 }
1034 }
1035 } else if (auto _pars = pars(); !_pars.empty()) {
1036 for (auto &v : _pars) {
1037 if (!v->get<RooAbsArg>()->getAttribute("Constant")) {
1038 x = v->get<RooAbsLValue>();
1039 if (x)
1040 break;
1041 } else if (!x) {
1042 x = v->get<RooAbsLValue>();
1043 }
1044 }
1045 }
1046
1047 if (!x) {
1048 return nullptr;
1049 }
1050 }
1051
1052 if (o != dynamic_cast<TObject *>(x)) {
1053 o->setStringAttribute("xvar", dynamic_cast<TObject *>(x)->GetName());
1054 }
1055
1056 // decide binning to use
1057 TString binningName = o->getStringAttribute("binning");
1058 auto _bnames = x->getBinningNames();
1059 bool hasBinning = false;
1060 for (auto &b : _bnames) {
1061 if (b == binningName) {
1062 hasBinning = true;
1063 break;
1064 }
1065 }
1066 if (!hasBinning) {
1067 // doesn't have binning, so clear binning attribute
1068 // this can happen after Combine of models because binning don't get combined yet (should fix this)
1069 Warning("GetXaxis", "Binning %s not defined on %s - clearing", binningName.Data(),
1070 dynamic_cast<TObject *>(x)->GetName());
1071 o->setStringAttribute("binning", nullptr);
1072 binningName = "";
1073 }
1074
1075 if (binningName == "" && o != dynamic_cast<TObject *>(x)) {
1076 // has var has a binning matching this nodes name then use that
1077 auto __bnames = x->getBinningNames();
1078 for (auto &b : __bnames) {
1079 if (b == GetName())
1080 binningName = GetName();
1081 if (b == o->GetName()) {
1082 binningName = o->GetName();
1083 break;
1084 } // best match
1085 }
1086 if (binningName == "") {
1087 // if we are binned in this var then will define that as a binning
1088 if (/*o->isBinnedDistribution(*dynamic_cast<RooAbsArg *>(x))*/
1089 auto bins = _or_func(
1090 /*o->plotSamplingHint(*dynamic_cast<RooAbsRealLValue
1091 *>(x),-std::numeric_limits<double>::infinity(),std::numeric_limits<double>::infinity())*/
1092 (std::list<double> *)(nullptr),
1093 o->binBoundaries(*dynamic_cast<RooAbsRealLValue *>(x), -std::numeric_limits<double>::infinity(),
1094 std::numeric_limits<double>::infinity()));
1095 bins) {
1096 std::vector<double> _bins;
1097 for (auto &b : *bins) {
1098 if (_bins.empty() || std::abs(_bins.back() - b) > 1e-5 * _bins.back())
1099 _bins.push_back(b);
1100 }
1101 fXAxis = std::make_shared<Axis2>(_bins.size() - 1, &_bins[0]);
1102 // add this binning to the var to avoid recalling ...
1103 if (auto _v = dynamic_cast<RooRealVar *>(x); _v) {
1104 _v->setBinning(RooBinning(_bins.size() - 1, &_bins[0], o->GetName()), o->GetName());
1105 _v->getBinning(o->GetName())
1106 .SetTitle(""); // indicates to use the current var title when building histograms etc
1107 //_v->getBinning(o->GetName()).SetTitle(strlen(dynamic_cast<TObject*>(x)->GetTitle()) ?
1108 // dynamic_cast<TObject*>(x)->GetTitle() : dynamic_cast<TObject*>(x)->GetName());
1109 }
1110 binningName = o->GetName();
1111 delete bins;
1112 } else if (_parentX) {
1113 // use parent axis binning if defined, otherwise we will default
1114 binningName = _parentX->GetName();
1115 }
1116 }
1117 }
1118
1119 if (!fXAxis) {
1120 if (auto r = dynamic_cast<RooAbsRealLValue *>(x); r) {
1121 if (r->getBinning(binningName).isUniform()) {
1122 fXAxis = std::make_shared<Axis2>(x->numBins(binningName), r->getMin(binningName), r->getMax(binningName));
1123 } else {
1124 fXAxis = std::make_shared<Axis2>(x->numBins(binningName), r->getBinning(binningName).array());
1125 }
1126 } else if (auto cat = dynamic_cast<RooCategory *>(x)) {
1127 std::vector<double> bins = {};
1128 for (int i = 0; i <= x->numBins(binningName); i++)
1129 bins.push_back(i);
1130 fXAxis = std::make_shared<Axis2>(x->numBins(binningName), &bins[0]);
1131 // TODO have to load current state of bin labels if was a category (sadly not a virtual method)
1132 int i = 1;
1133 std::map<int, std::string> cats; // fill into a map to preserve index ordering
1134 for (auto &c : *cat) {
1135 if (cat->isStateInRange(binningName, c.first.c_str())) {
1136 cats[c.second] = c.first;
1137 }
1138 }
1139 for (auto &[_, label] : cats) {
1140 fXAxis->SetBinLabel(i++, label.c_str());
1141 }
1142 }
1143 }
1144
1145 fXAxis->SetName(binningName);
1146 fXAxis->SetParent(dynamic_cast<TObject *>(x));
1147 return fXAxis.get();
1148}
1149
1150const char *xRooNode::GetIconName() const
1151{
1152 if (auto o = get(); o) {
1153 if (o->InheritsFrom("RooWorkspace"))
1154 return "TFile";
1155 if (o->InheritsFrom("RooAbsData"))
1156 return "TProfile";
1157 if (o->InheritsFrom("RooSimultaneous"))
1158 return "TH3D";
1159
1160 if (o->InheritsFrom("RooProdPdf"))
1161 return "a.C"; // or nullptr for folder
1162 if (o->InheritsFrom("RooRealSumPdf") || o->InheritsFrom("RooAddPdf"))
1163 return "TH2D";
1164 // if(o->InheritsFrom("RooProduct")) return "TH1D";
1165 if (o->InheritsFrom("RooFitResult")) {
1166 if (!gClient->GetMimeTypeList()->GetIcon("xRooFitRooFitResult", true)) {
1167 gClient->GetMimeTypeList()->AddType("xRooFitRooFitResult", "xRooFitRooFitResult", "package.xpm",
1168 "package.xpm", "->Browse()");
1169 }
1170 return "xRooFitRooFitResult";
1171 }
1172 if (o->InheritsFrom("RooRealVar") || o->InheritsFrom("RooCategory")) {
1173 if (get<RooAbsArg>()->getAttribute("obs")) {
1174 if (!gClient->GetMimeTypeList()->GetIcon("xRooFitObs", true)) {
1175 gClient->GetMimeTypeList()->AddType("xRooFitObs", "xRooFitObs", "x_pic.xpm", "x_pic.xpm", "->Browse()");
1176 }
1177 if (!gClient->GetMimeTypeList()->GetIcon("xRooFitGlobs", true)) {
1178 gClient->GetMimeTypeList()->AddType("xRooFitGlobs", "xRooFitGlobs", "z_pic.xpm", "z_pic.xpm",
1179 "->Browse()");
1180 }
1181 return (get<RooAbsArg>()->getAttribute("global") ? "xRooFitGlobs" : "xRooFitObs");
1182 }
1183 return "TLeaf";
1184 }
1185 if (o->InheritsFrom("TStyle")) {
1186 if (!gClient->GetMimeTypeList()->GetIcon("xRooFitTStyle", true)) {
1187 gClient->GetMimeTypeList()->AddType("xRooFitTStyle", "xRooFitTStyle", "bld_colorselect.xpm",
1188 "bld_colorselect.xpm", "->Browse()");
1189 }
1190 return "xRooFitTStyle";
1191 }
1192 if (o->InheritsFrom("RooConstVar")) {
1193 /*if (!gClient->GetMimeTypeList()->GetIcon("xRooFitRooConstVar",true)) {
1194 gClient->GetMimeTypeList()->AddType("xRooFitRooConstVar", "xRooFitRooConstVar", "stop_t.xpm", "stop_t.xpm",
1195 "->Browse()");
1196 }
1197 return "xRooFitRooConstVar";*/
1198 return "TMethodBrowsable-leaf";
1199 }
1200 if (o->InheritsFrom("RooStats::HypoTestInverterResult")) {
1201 if (!gClient->GetMimeTypeList()->GetIcon("xRooFitScanStyle", true)) {
1202 gClient->GetMimeTypeList()->AddType("xRooFitScanStyle", "xRooFitScanStyle", "f2_s.xpm", "f2_s.xpm",
1203 "->Browse()");
1204 }
1205 return "xRooFitScanStyle";
1206 }
1207 if (o->InheritsFrom("RooStats::HypoTestResult")) {
1208 if (!gClient->GetMimeTypeList()->GetIcon("xRooFitTestStyle", true)) {
1209 gClient->GetMimeTypeList()->AddType("xRooFitTestStyle", "xRooFitTestStyle", "diamond.xpm", "diamond.xpm",
1210 "->Browse()");
1211 }
1212 return "xRooFitTestStyle";
1213 }
1214 if (o->InheritsFrom("RooStats::HistFactory::FlexibleInterpVar"))
1215 return "TBranchElement-folder";
1216 if (o->InheritsFrom("RooAbsPdf")) {
1217 if (!gClient->GetMimeTypeList()->GetIcon("xRooFitPDFStyle", true)) {
1218 gClient->GetMimeTypeList()->AddType("xRooFitPDFStyle", "xRooFitPDFStyle", "pdf.xpm", "pdf.xpm",
1219 "->Browse()");
1220 }
1221 return "xRooFitPDFStyle";
1222 }
1223 if (o->InheritsFrom("RooStats::ModelConfig")) {
1224 if (!gClient->GetMimeTypeList()->GetIcon("xRooFitMCStyle", true)) {
1225 gClient->GetMimeTypeList()->AddType("xRooFitMCStyle", "xRooFitMCStyle", "app_t.xpm", "app_t.xpm",
1226 "->Browse()");
1227 }
1228 return "xRooFitMCStyle";
1229 }
1230 if (auto a = dynamic_cast<RooAbsReal *>(o); a) {
1231 if (auto _ax = GetXaxis();
1232 _ax && (a->isBinnedDistribution(*dynamic_cast<RooAbsArg *>(_ax->GetParent())) ||
1233 (dynamic_cast<RooAbsRealLValue *>(_ax->GetParent()) &&
1234 std::unique_ptr<std::list<double>>(a->binBoundaries(
1235 *dynamic_cast<RooAbsRealLValue *>(_ax->GetParent()), -std::numeric_limits<double>::infinity(),
1236 std::numeric_limits<double>::infinity()))))) {
1237 return "TH1D";
1238 }
1239 return "TF1";
1240 }
1241 return o->ClassName();
1242 }
1243 if (!IsFolder()) {
1244 return "Unknown";
1245 }
1246 return nullptr;
1247}
1248
1249const char *xRooNode::GetNodeType() const
1250{
1251 if (auto o = get(); o && fParent && (fParent->get<RooProduct>() || fParent->get<RooRealSumPdf>())) {
1252 if (o->InheritsFrom("RooStats::HistFactory::FlexibleInterpVar"))
1253 return "Overall";
1254 if (o->InheritsFrom("PiecewiseInterpolation"))
1255 return (dynamic_cast<RooAbsArg *>(o)->getAttribute("density")) ? "DensityHisto" : "Histo";
1256 if (o->InheritsFrom("RooHistFunc"))
1257 return (dynamic_cast<RooAbsArg *>(o)->getAttribute("density")) ? "ConstDensityHisto" : "ConstHisto";
1258 if (o->InheritsFrom("RooBinWidthFunction"))
1259 return "Density";
1260 if (o->InheritsFrom("ParamHistFunc"))
1261 return "Shape";
1262 if (o->InheritsFrom("RooRealVar"))
1263 return "Norm";
1264 if (o->InheritsFrom("RooConstVar"))
1265 return "Const";
1266 }
1267 return "";
1268}
1269
1270xRooNode xRooNode::coords(bool setVals) const
1271{
1272 xRooNode out(".coords", nullptr, *this);
1273 // go up through parents looking for slice obs
1274 auto _p = std::shared_ptr<xRooNode>(const_cast<xRooNode *>(this), [](xRooNode *) {});
1275 while (_p) {
1276 TString pName(_p->GetName());
1277 // following is commented out while still considering, but idea is to include category in coords
1278 /*if (auto s = _p->get<RooSimultaneous>(); s && s->indexCat().InheritsFrom("RooCategory") &&
1279 !out.find(s->indexCat().GetName())) { auto cat = const_cast<RooCategory*>(dynamic_cast<const
1280 RooCategory*>(&s->indexCat()));
1281 // check if we have a pdf for every category ... if not then add to cut
1282 cat->clearRange("coordRange",true);
1283 bool hasMissing = false;
1284 std::string includedStates;
1285 for (auto state : *cat) {
1286 if (!s->getPdf(state.first.c_str())) {
1287 hasMissing = true;
1288 } else {
1289 if (!includedStates.empty()) {
1290 includedStates += ",";
1291 }
1292 includedStates += state.first;
1293 }
1294 }
1295 if (hasMissing) {
1296 if(includedStates.find(",") != std::string::npos) {
1297 cat->addToRange("coordRange",includedStates.c_str());
1298 } else {
1299 cat->setLabel(includedStates);
1300 }
1301 out.emplace_back(std::make_shared<xRooNode>(cat->GetName(),_p->getObject<RooAbsArg>(cat->GetName()),_p));
1302 }
1303 } else*/
1304 if (auto pos = pName.Index('='); pos != -1) {
1305 if (pos > 0 && pName(pos - 1) == '<') {
1306 // should be a range on a real lvalue, of form low<=name<high
1307 double low = TString(pName(0, pos - 1)).Atof();
1308 pName = pName(pos + 1, pName.Length());
1309 double high = TString(pName(pName.Index('<') + 1, pName.Length())).Atof();
1310 pName = pName(0, pName.Index('<'));
1311 if (auto _obs = _p->getObject<RooAbsRealLValue>(pName.Data()); _obs) {
1312 if (setVals) {
1313 _obs->setVal((high + low) / 2.);
1314 static_cast<RooRealVar *>(_obs.get())->setRange("coordRange", low, high);
1315 _obs->setStringAttribute(
1316 "coordRange", "coordRange"); // will need if we allow multi disconnected regions, need comma list
1317 }
1318 out.emplace_back(std::make_shared<xRooNode>(_obs->GetName(), _obs, _p));
1319 } else {
1320 throw std::runtime_error(TString::Format("Unknown observable: %s", pName.Data()));
1321 }
1322
1323 } else if (auto _obs = _p->getObject<RooAbsArg>(pName(0, pos)); _obs) {
1324 if (setVals) {
1325 if (auto _cat = dynamic_cast<RooAbsCategoryLValue *>(_obs.get()); _cat) {
1326 _cat->setLabel(pName(pos + 1, pName.Length()));
1327 } else if (auto _var = dynamic_cast<RooAbsRealLValue *>(_obs.get()); _var) {
1328 _var->setVal(TString(pName(pos + 1, pName.Length())).Atof());
1329 }
1330 }
1331 out.emplace_back(std::make_shared<xRooNode>(_obs->GetName(), _obs, _p));
1332 } else {
1333 throw std::runtime_error("Unknown observable, could not find");
1334 }
1335 }
1336 _p = _p->fParent;
1337 }
1338 return out;
1339}
1340
1341void xRooNode::_Add_(const char *name, const char *opt)
1342{
1343 try {
1344 Add(name, opt);
1345 } catch (const std::exception &e) {
1346 new TGMsgBox(gClient->GetRoot(), gClient->GetRoot(), "Exception", e.what(),
1347 kMBIconExclamation); // deletes self on dismiss?
1348 }
1349}
1350void xRooNode::_Vary_(const char *what)
1351{
1352 try {
1353 Vary(what);
1354 } catch (const std::exception &e) {
1355 new TGMsgBox(gClient->GetRoot(), gClient->GetRoot(), "Exception", e.what(),
1356 kMBIconExclamation); // deletes self on dismiss?
1357 }
1358}
1359
1361{
1362
1363 if (strcmp(GetName(), ".poi") == 0) {
1364 // demote a parameter from being a poi
1365 auto toRemove =
1366 (child.get<RooAbsArg>() || !find(child.GetName())) ? child : xRooNode(find(child.GetName())->fComp);
1367 if (toRemove) {
1368 if (!toRemove.get<RooAbsArg>()->getAttribute("poi")) {
1369 throw std::runtime_error(TString::Format("%s is not a poi", toRemove.GetName()));
1370 }
1371 toRemove.get<RooAbsArg>()->setAttribute("poi", false);
1372 return toRemove;
1373 }
1374 } else if (strcmp(GetName(), ".factors") == 0 || strcmp(GetName(), ".constraints") == 0 ||
1375 strcmp(GetName(), ".components") == 0) {
1376 auto toRemove =
1377 (child.get<RooAbsArg>() || !find(child.GetName())) ? child : xRooNode(find(child.GetName())->fComp);
1378 if (auto p = fParent->get<RooProdPdf>(); p) {
1379 auto pdf = toRemove.get<RooAbsArg>();
1380 if (!pdf)
1381 pdf = p->pdfList().find(child.GetName());
1382 if (!pdf)
1383 throw std::runtime_error(TString::Format("Cannot find %s in %s", child.GetName(), fParent->GetName()));
1384 auto i = p->pdfList().index(*pdf);
1385 if (i >= 0) {
1386#if ROOT_VERSION_CODE < ROOT_VERSION(6, 27, 00)
1387 const_cast<RooArgList &>(p->pdfList()).remove(*pdf);
1388#if ROOT_VERSION_CODE >= ROOT_VERSION(6, 26, 00)
1389 p->_pdfNSetList.erase(p->_pdfNSetList.begin() + i);
1390#else
1391 auto nset = p->_pdfNSetList.At(i);
1392 p->_pdfNSetList.Remove(nset);
1393 delete nset; // I don't think the RooLinkedList owned it so must delete ourself
1394#endif
1395 if (p->_extendedIndex == i)
1396 p->_extendedIndex = -1;
1397 else if (p->_extendedIndex > i)
1398 p->_extendedIndex--;
1399#else
1400 p->removePdfs(RooArgSet(*pdf));
1401#endif
1402 sterilize();
1403 return xRooNode(*pdf);
1404 } else {
1405 throw std::runtime_error(TString::Format("Cannot find %s in %s", child.GetName(), fParent->GetName()));
1406 }
1407 } else if (auto p2 = fParent->get<RooProduct>(); p2) {
1408 auto arg = toRemove.get<RooAbsArg>();
1409 if (!arg)
1410 arg = p2->components().find(child.GetName());
1411 if (!arg)
1412 throw std::runtime_error(TString::Format("Cannot find %s in %s", child.GetName(), fParent->GetName()));
1413 // remove server ... doesn't seem to trigger removal from proxy
1414#if ROOT_VERSION_CODE < ROOT_VERSION(6, 27, 00)
1415 p2->_compRSet.remove(*arg);
1416#else
1417 const_cast<RooArgList &>(p2->realComponents()).remove(*arg);
1418#endif
1419 p2->removeServer(*arg, true);
1420 sterilize();
1421 return xRooNode(*arg);
1422 } else if (fParent->get<RooSimultaneous>()) {
1423 // remove from all channels
1424 bool removed = false;
1425 for (auto &c : fParent->bins()) {
1426 try {
1427 c->constraints().Remove(toRemove);
1428 removed = true;
1429 } catch (std::runtime_error &) { /* wasn't a constraint in channel */
1430 }
1431 }
1432 sterilize();
1433 if (!removed)
1434 throw std::runtime_error(TString::Format("Cannot find %s in %s", child.GetName(), fParent->GetName()));
1435 return toRemove;
1436 } else if (auto p4 = fParent->get<RooRealSumPdf>(); p4) {
1437 auto arg = toRemove.get<RooAbsArg>();
1438 if (!arg)
1439 arg = p4->funcList().find(child.GetName());
1440 if (!arg)
1441 throw std::runtime_error(TString::Format("Cannot find %s in %s", child.GetName(), fParent->GetName()));
1442 // remove, including coef removal ....
1443 auto idx = p4->funcList().index(arg);
1444
1445 if (idx != -1) {
1446
1447 const_cast<RooArgList &>(p4->funcList()).remove(*arg);
1448 p4->removeServer(*arg, true);
1449 // have to be careful removing coef because if shared will end up removing them all!!
1450 std::vector<RooAbsArg *> _coefs;
1451 for (size_t ii = 0; ii < const_cast<RooArgList &>(p4->coefList()).size(); ii++) {
1452 if (ii != size_t(idx))
1453 _coefs.push_back(const_cast<RooArgList &>(p4->coefList()).at(ii));
1454 }
1455 const_cast<RooArgList &>(p4->coefList()).removeAll();
1456 for (auto &a : _coefs)
1457 const_cast<RooArgList &>(p4->coefList()).add(*a);
1458
1459 sterilize();
1460 } else {
1461 throw std::runtime_error(TString::Format("Cannot find %s in %s", child.GetName(), fParent->GetName()));
1462 }
1463 return xRooNode(*arg);
1464 } else if (auto p5 = fParent->get<RooAddPdf>(); p5) {
1465 auto arg = toRemove.get<RooAbsArg>();
1466 if (!arg)
1467 arg = p5->pdfList().find(child.GetName());
1468 if (!arg)
1469 throw std::runtime_error(TString::Format("Cannot find %s in %s", child.GetName(), fParent->GetName()));
1470 // remove, including coef removal ....
1471 auto idx = p5->pdfList().index(arg);
1472
1473 if (idx != -1) {
1474
1475 const_cast<RooArgList &>(p5->pdfList()).remove(*arg);
1476 p5->removeServer(*arg, true);
1477 // have to be careful removing coef because if shared will end up removing them all!!
1478 std::vector<RooAbsArg *> _coefs;
1479 for (size_t ii = 0; ii < const_cast<RooArgList &>(p5->coefList()).size(); ii++) {
1480 if (ii != size_t(idx))
1481 _coefs.push_back(const_cast<RooArgList &>(p5->coefList()).at(ii));
1482 }
1483 const_cast<RooArgList &>(p5->coefList()).removeAll();
1484 for (auto &a : _coefs)
1485 const_cast<RooArgList &>(p5->coefList()).add(*a);
1486
1487 sterilize();
1488 } else {
1489 throw std::runtime_error(TString::Format("Cannot find %s in %s", child.GetName(), fParent->GetName()));
1490 }
1491 return xRooNode(*arg);
1492 } else if (auto p6 = fParent->get<RooAddition>(); p6) {
1493 auto arg = toRemove.get<RooAbsArg>();
1494 if (!arg)
1495 arg = p6->list().find(child.GetName());
1496 if (!arg)
1497 throw std::runtime_error(TString::Format("Cannot find %s in %s", child.GetName(), fParent->GetName()));
1498 // remove server ... doesn't seem to trigger removal from proxy
1499 const_cast<RooArgList &>(p6->list()).remove(*arg);
1500 p6->removeServer(*arg, true);
1501 sterilize();
1502 return xRooNode(*arg);
1503 }
1504 }
1505
1506 if (auto w = get<RooWorkspace>(); w) {
1507 xRooNode out(child.GetName());
1508 auto arg = w->components().find(child.GetName());
1509 if (!arg)
1510 arg = operator[](child.GetName())->get<RooAbsArg>();
1511 if (!arg) {
1512 throw std::runtime_error(TString::Format("Cannot find %s in workspace %s", child.GetName(), GetName()));
1513 }
1514 // check has no clients ... if so, cannot delete
1515 if (arg->hasClients()) {
1516 throw std::runtime_error(
1517 TString::Format("Cannot remove %s from workspace %s, because it has dependencies - first remove from those",
1518 child.GetName(), GetName()));
1519 }
1520 const_cast<RooArgSet &>(w->components()).remove(*arg); // deletes arg
1521 Info("Remove", "Deleted %s from workspace %s", out.GetName(), GetName());
1522 return out;
1523 } else if (get<RooProduct>() || get<RooProdPdf>()) {
1524 return factors().Remove(child);
1525 } else if (get<RooRealSumPdf>() || get<RooAddPdf>() || get<RooAddition>()) {
1526 return components().Remove(child);
1527 }
1528
1529 throw std::runtime_error("Removal not implemented for object type " +
1530 std::string(get() ? get()->ClassName() : "null"));
1531}
1532
1534{
1535
1536 class AutoUpdater {
1537 public:
1538 AutoUpdater(xRooNode &_n) : n(_n) {}
1539 ~AutoUpdater() { n.browse(); }
1540 xRooNode &n;
1541 };
1542 AutoUpdater xxx(*this);
1543
1544 TString sOpt(opt);
1545 bool considerType(sOpt == "+");
1546
1547 if (strlen(GetName()) > 0 && GetName()[0] == '!' && fParent) {
1548 // folder .. pass onto parent and add folder to child folder list
1549 const_cast<xRooNode &>(child).fFolder += GetName();
1550 return fParent->Add(child, opt);
1551 }
1552 // this is how to get the first real parent ... may be useful at some point?
1553 /*auto realParent = fParent;
1554 while(!realParent->get()) {
1555 realParent = realParent->fParent;
1556 if (!realParent) throw std::runtime_error("No parentage");
1557 }*/
1558
1559 // adding to a collection node will incorporate the child into the parent of the collection
1560 // in the appropriate way
1561 if (strcmp(GetName(), ".factors") == 0) {
1562 // multiply the parent
1563 return fParent->Multiply(child, opt);
1564 } else if (strcmp(GetName(), ".components") == 0) {
1565 // add to the parent
1566 return fParent->Add(child, opt);
1567 } else if (strcmp(GetName(), ".variations") == 0) {
1568 // vary the parent
1569 return fParent->Vary(child);
1570 } else if (strcmp(GetName(), ".constraints") == 0) {
1571 // constrain the parent
1572 return fParent->Constrain(child);
1573 } else if (strcmp(GetName(), ".bins") == 0 && fParent->get<RooSimultaneous>()) {
1574 // adding a channel (should adding a 'bin' be an 'Extend' operation?)
1575 return fParent->Vary(child);
1576 } else if ((strcmp(GetName(), ".globs") == 0)) {
1577 if (child.get<RooAbsArg>() || (!child.fComp && getObject<RooAbsArg>(child.GetName()))) {
1578 auto out = (child.get<RooAbsArg>()) ? child.get<RooAbsArg>() : getObject<RooAbsArg>(child.GetName()).get();
1579 out->setAttribute("obs");
1580 out->setAttribute("global");
1581 return xRooNode(*out, *this);
1582 }
1583 throw std::runtime_error("Failed to add global observable");
1584 } else if ((strcmp(GetName(), ".poi") == 0)) {
1585 if (child.get<RooAbsLValue>() || (!child.fComp && getObject<RooAbsLValue>(child.GetName()))) {
1586 auto out = (child.get<RooAbsArg>()) ? child.get<RooAbsArg>() : getObject<RooAbsArg>(child.GetName()).get();
1587 out->setAttribute("poi");
1588 return xRooNode(*out, *this);
1589 } else if (!child.get() && fParent->get<RooWorkspace>()) {
1590 // may be creating poi at same time as adding, try add to parent
1591 auto res = fParent->Add(child);
1592 if (res.get<RooAbsLValue>())
1593 return Add(res);
1594 }
1595 throw std::runtime_error("Failed to add parameter of interest");
1596 } else if ((strcmp(GetName(), ".pars") == 0 || strcmp(GetName(), ".vars") == 0) && fParent->get<RooWorkspace>()) {
1597 // adding a parameter, interpret as factory string unless no "[" then create RooRealVar
1598 TString fac(child.GetName());
1599 if (!fac.Contains("["))
1600 fac += "[1]";
1601 return xRooNode(*fParent->get<RooWorkspace>()->factory(fac), fParent);
1602 } else if (strcmp(GetName(), ".datasets()") == 0) {
1603 // create a dataset - only allowed for pdfs or workspaces
1604 if (auto _ws = ws(); _ws && fParent) {
1605 sOpt.ToLower();
1606 if (!fParent->get<RooAbsPdf>() && (!fParent->get<RooWorkspace>() || sOpt == "asimov")) {
1607 throw std::runtime_error(
1608 "Datasets can only be created for pdfs or workspaces (except if generated dataset, then must be pdf)");
1609 }
1610
1611 if (sOpt == "asimov" || sOpt == "toy") {
1612 // generate expected dataset - note that globs will be frozen at this time
1613 auto _fr = fParent->fitResult();
1614 if (strlen(_fr->GetName()) == 0) { // ensure fit result has a name so that name is saved inside dataset
1615 _fr.get<RooFitResult>()->SetName(TUUID().AsString());
1616 }
1617 auto ds = fParent->generate(_fr, sOpt == "asimov");
1618 if (strlen(child.GetName())) {
1619 ds.SetName(child.GetName());
1620 ds.get<TNamed>()->SetName(child.GetName());
1621 }
1622 if (auto _ds = ds.get<RooAbsData>()) {
1623 _ws->import(*_ds);
1624 }
1625 if (_fr.get<RooFitResult>()->numStatusHistory() == 0) {
1626 if (!GETWSSNAPSHOTS(_ws).find(_fr->GetName())) {
1627 const_cast<RooLinkedList &>(GETWSSNAPSHOTS(_ws)).Add(_fr->Clone());
1628 }
1629 } else if (!_ws->obj(_fr->GetName())) {
1630 _ws->import((*_fr.get<RooFitResult>()));
1631 } // save fr to workspace, for later retrieval
1632 return xRooNode(*_ws->data(ds.GetName()), fParent);
1633 }
1634
1635 auto parentObs = fParent->obs(); // may own globs so keep alive
1636 auto _obs = parentObs.argList();
1637 // put globs in a snapshot
1638 std::unique_ptr<RooAbsCollection> _globs(_obs.selectByAttrib("global", true));
1639 // RooArgSet _tmp; _tmp.add(*_globs);_ws->saveSnapshot(child.GetName(),_tmp);
1640 _obs.remove(*_globs);
1641
1642 // include any coords
1643 _obs.add(coords(false).argList(), true);
1644 // include axis var too, provided it's an observable
1645 if (auto ax = GetXaxis(); ax && dynamic_cast<RooAbsArg *>(ax->GetParent())->getAttribute("obs")) {
1646 _obs.add(*dynamic_cast<RooAbsArg *>(ax->GetParent()));
1647 }
1648 // check if ws already has a dataset with this name, if it does we may need to extend columns
1649 if (auto _d = _ws->data(child.GetName()); _d) {
1650 // add any missing obs
1651 RooArgSet l(_obs);
1652 l.remove(*_d->get(), true, true);
1653 if (!l.empty()) {
1654 auto _dd = dynamic_cast<RooDataSet *>(_d);
1655 if (!_dd)
1656 throw std::runtime_error("Cannot extend dataset with new columns");
1657 for (auto &x : l) {
1658 _dd->addColumn(*x);
1659 }
1660 }
1661 } else {
1662 RooRealVar w("weightVar", "weightVar", 1);
1663 _obs.add(w);
1664 RooDataSet d(child.GetName(), child.GetTitle(), _obs, RooFit::WeightVar("weightVar"));
1665 _ws->import(d);
1666 // seems have to set bits after importing, not before
1667 if (auto __d = _ws->data(child.GetName()))
1668 __d->SetBit(1 << 20, _ws->allData().size() == 1); // sets as selected if is only ds
1669 }
1670 /*if(!_ws->data(child.GetName())) {
1671 RooRealVar w("weightVar", "weightVar", 1);
1672 RooArgSet _obs; _obs.add(w);
1673 RooDataSet d(child.GetName(), child.GetTitle(), _obs, "weightVar");
1674 _ws->import(d);
1675 }*/
1676 auto out = std::shared_ptr<TObject>(_ws->data(child.GetName()), [](TObject *) {});
1677
1678 if (out) {
1679 xRooNode o(out, fParent);
1680 if (child.get<TH1>())
1681 o = *child.get();
1682 return o;
1683 }
1684 }
1685 throw std::runtime_error("Cannot create dataset");
1686 }
1687
1688 if (!get()) {
1689 if (!fParent)
1690 throw std::runtime_error("Cannot add to null object with no parentage");
1691
1692 auto _ref = emplace_back(std::shared_ptr<xRooNode>(&const_cast<xRooNode &>(child), [](TObject *) {}));
1693 try {
1694 fComp = fParent->Add(*this, "+").fComp;
1695 } catch (...) {
1696 resize(size() - 1);
1697 std::rethrow_exception(std::current_exception());
1698 }
1699 resize(size() - 1); // remove the temporarily added node
1700
1701 if (!fComp) {
1702 throw std::runtime_error("No object");
1703 }
1704 }
1705
1706 if (auto p = get<RooAbsData>(); p) {
1707 if (auto bb = getBrowsable(".sourceds"))
1708 bb->Add(child, opt);
1709 if (auto _data = child.get<RooDataSet>()) {
1710 auto ds = dynamic_cast<RooDataSet *>(p);
1711 if (!ds) {
1712 throw std::runtime_error("Can only add datasets to a dataset");
1713 }
1714
1715 // append any missing globs, and check any existing globs have matching values
1716 RooArgList globsToAdd;
1717 auto _globs = globs();
1718 for (auto &glob : child.globs()) {
1719 if (auto g = _globs.find(glob->GetName()); !g) {
1720 globsToAdd.addClone(*glob->get<RooAbsArg>());
1721 } else if (g->GetContent() != glob->GetContent()) {
1722 Warning("Add", "Global observable %s=%g in dataset %s mismatches %s value %g ... ignoring latter",
1723 g->GetName(), g->GetContent(), GetName(), child.GetName(), glob->GetContent());
1724 }
1725 }
1726 // add any existing globs to list then set the list
1727 if (auto _dglobs = p->getGlobalObservables()) {
1728 globsToAdd.addClone(*_dglobs);
1729 } else {
1730 for (auto g : _globs)
1731 globsToAdd.addClone(*g->get<RooAbsArg>());
1732 }
1733 p->setGlobalObservables(globsToAdd);
1734
1735 // append any missing observables to our dataset, then append the dataset
1736
1737 for (auto col : *_data->get()) {
1738 if (!p->get()->contains(*col)) {
1739 ds->addColumn(*col);
1740 }
1741 }
1742 ds->append(*_data);
1743 ds->SetTitle(TString(ds->GetTitle()) + " + " + _data->GetTitle());
1744 SetTitle(TString(GetTitle()) + " + " + child.GetTitle());
1745 return *this;
1746 }
1747 auto _h = child.get<TH1>();
1748 if (!_h) {
1749 throw std::runtime_error("Can only add histogram or dataset to data");
1750 }
1751 auto _pdf = parentPdf();
1752 if (!_pdf)
1753 throw std::runtime_error("Could not find pdf");
1754 auto _ax = _pdf->GetXaxis();
1755 if (!_ax) {
1756 throw std::runtime_error("Cannot determine binning to add data");
1757 }
1758
1759 RooArgSet obs;
1760 obs.add(*dynamic_cast<RooAbsArg *>(_ax->GetParent()));
1761 obs.add(coords().argList()); // will also move obs to coords
1762
1763 // add any missing obs
1764 RooArgSet l(obs);
1765 l.remove(*p->get(), true, true);
1766 if (!l.empty()) {
1767 auto _d = dynamic_cast<RooDataSet *>(p);
1768 if (!_d)
1769 throw std::runtime_error("Cannot extend dataset with new columns");
1770 for (auto &x : l) {
1771 _d->addColumn(*x);
1772 }
1773 }
1774
1775 // before adding, ensure range is good to cover
1776 for (auto &o : obs) {
1777 if (auto v = dynamic_cast<RooRealVar *>(o); v) {
1778 if (auto dv = dynamic_cast<RooRealVar *>(p->get()->find(v->GetName())); dv) {
1779 if (v->getMin() < dv->getMin())
1780 dv->setMin(v->getMin());
1781 if (v->getMax() > dv->getMax())
1782 dv->setMax(v->getMax());
1783 }
1784 } else if (auto c = dynamic_cast<RooCategory *>(o); c) {
1785 if (auto dc = dynamic_cast<RooCategory *>(p->get()->find(c->GetName())); dc) {
1786 if (!dc->hasLabel(c->getCurrentLabel())) {
1787 dc->defineType(c->getCurrentLabel(), c->getCurrentIndex());
1788 }
1789 }
1790 }
1791 }
1792
1793 for (int i = 1; i <= _h->GetNbinsX(); i++) {
1794 if (auto cat = dynamic_cast<RooAbsCategoryLValue *>(_ax->GetParent())) {
1795 if (!_h->GetXaxis()->GetBinLabel(i)) {
1796 throw std::runtime_error(
1797 TString::Format("Categorical observable %s requires bin labels", _ax->GetParent()->GetName()));
1798 } else if (!cat->hasLabel(_h->GetXaxis()->GetBinLabel(i))) {
1799 throw std::runtime_error(TString::Format("Categorical observable %s does not have label %s",
1800 _ax->GetParent()->GetName(), _h->GetXaxis()->GetBinLabel(i)));
1801 } else {
1802 cat->setLabel(_h->GetXaxis()->GetBinLabel(i));
1803 }
1804 } else {
1805 dynamic_cast<RooAbsRealLValue *>(_ax->GetParent())->setVal(_h->GetBinCenter(i));
1806 }
1807 p->add(obs, _h->GetBinContent(i));
1808 }
1809
1810 return *this;
1811 }
1812
1813 if (auto p = get<RooAddPdf>(); p) {
1814 if ((child.get<RooAbsPdf>() || (!child.fComp && getObject<RooAbsPdf>(child.GetName())))) {
1815 auto out = (child.fComp) ? acquire(child.fComp) : getObject<RooAbsArg>(child.GetName());
1816 // don't add a coef if in 'all-extended' mode and this pdf is extendable
1817 auto _pdf = std::dynamic_pointer_cast<RooAbsPdf>(out);
1818 if (!_pdf) {
1819 throw std::runtime_error("Something went wrong with pdf acquisition");
1820 }
1821
1822 if (auto _ax = GetXaxis(); _ax && dynamic_cast<RooAbsRealLValue *>(_ax->GetParent()) &&
1823 _pdf->dependsOn(*static_cast<RooAbsArg *>(_ax->GetParent()))) {
1824 auto _p = _pdf;
1825
1826 if (auto _boundaries = std::unique_ptr<std::list<double>>(_p->binBoundaries(
1827 *dynamic_cast<RooAbsRealLValue *>(_ax->GetParent()), -std::numeric_limits<double>::infinity(),
1828 std::numeric_limits<double>::infinity()));
1829 !_boundaries && _ax->GetNbins() > 0) {
1830#if ROOT_VERSION_CODE >= ROOT_VERSION(6, 24, 00)
1831 Warning("Add", "Adding unbinned pdf %s to binned %s - will wrap with RooBinSamplingPdf(...)",
1832 _p->GetName(), GetName());
1833 _p = acquireNew<RooBinSamplingPdf>(TString::Format("%s_binned", _p->GetName()), _p->GetTitle(),
1834 *dynamic_cast<RooAbsRealLValue *>(_ax->GetParent()), *_p);
1835 _p->setStringAttribute("alias", std::dynamic_pointer_cast<RooAbsArg>(out)->getStringAttribute("alias"));
1836 if (!_p->getStringAttribute("alias"))
1837 _p->setStringAttribute("alias", out->GetName());
1838#else
1839 throw std::runtime_error(
1840 "unsupported addition of unbinned pdf to binned model - please upgrade to at least ROOT 6.24");
1841#endif
1842 _pdf = _p;
1843 }
1844 }
1845
1846 if (!(_pdf->canBeExtended() && p->coefList().empty())) {
1847 // if extended, use an extended binding as the coef
1848 // otherwise e.g. if adding a RooRealSumPdf the stacked histograms will be above the
1849 // actual pdf histogram because the pdf histogram is just normalized down
1850 if (_pdf->canBeExtended()) {
1851 // FIXME: ExtendedBinding needs the obs list passing to it ... should be fixed in RooFit
1852 // until then, this will return "1" and so the pdf's histograms wont be normalized properly in relation
1853 // to stacks of its comps
1854 const_cast<RooArgList &>(p->coefList())
1855 .add(*acquireNew<RooExtendedBinding>(TString::Format("%s_extBind", _pdf->GetName()),
1856 TString::Format("Expected Events of %s", _pdf->GetTitle()),
1857 *_pdf));
1858 } else {
1859
1860 // need to create a coefficient for each existing pdf first, like above
1861 for (auto i = p->coefList().size(); i < p->pdfList().size(); i++) {
1862 const_cast<RooArgList &>(p->coefList())
1863 .add(*acquireNew<RooExtendedBinding>(
1864 TString::Format("%s_extBind", p->pdfList().at(i)->GetName()),
1865 TString::Format("Expected Events of %s", p->pdfList().at(i)->GetTitle()),
1866 *static_cast<RooAbsPdf *>(p->pdfList().at(i))));
1867 }
1868
1869 const_cast<RooArgList &>(p->coefList()).add(*acquire2<RooAbsArg, RooRealVar>("1", "1", 1));
1870 }
1871 // ensure not in no-coef mode any more
1872 *reinterpret_cast<bool *>(reinterpret_cast<unsigned char *>(p) +
1873 p->Class()->GetDataMemberOffset("_allExtendable")) = false;
1874 *reinterpret_cast<bool *>(reinterpret_cast<unsigned char *>(p) +
1875 p->Class()->GetDataMemberOffset("_haveLastCoef")) = true;
1876 }
1877 const_cast<RooArgList &>(p->pdfList()).add(*_pdf);
1878 sterilize();
1879 return xRooNode(*_pdf, *this);
1880 } else if ((child.get<TH1>() || child.get<RooAbsReal>() ||
1881 (!child.get() && getObject<RooAbsReal>(child.GetName()))) &&
1882 !child.get<RooAbsPdf>()) {
1883 RooRealSumPdf *_pdf = nullptr;
1884 bool tooMany(false);
1885 for (auto &pp : factors()) {
1886 if (auto _p = pp->get<RooRealSumPdf>(); _p) {
1887 if (_pdf) {
1888 _pdf = nullptr;
1889 tooMany = true;
1890 break;
1891 } // more than one!
1892 _pdf = _p;
1893 }
1894 }
1895 if (_pdf) {
1896 return xRooNode(*_pdf, *this).Add(child);
1897 } else if (!tooMany) {
1898 // create a RooRealSumPdf to hold the child
1899 auto _sumpdf = Add(*acquireNew<RooRealSumPdf>(TString::Format("%s_samples", p->GetName()),
1900 TString::Format("%s samples", GetTitle()), RooArgList(),
1901 RooArgList(), true));
1902 _sumpdf.get<RooAbsArg>()->setStringAttribute("alias", "samples");
1903 return _sumpdf.Add(child);
1904 }
1905 }
1906 }
1907
1908 if (auto p = get<RooRealSumPdf>(); p) {
1909 std::shared_ptr<TObject> out;
1910 auto cc = child.fComp;
1911 bool isConverted = (cc != child.convertForAcquisition(*this, sOpt));
1912 if (child.get<RooAbsReal>()) {
1913 out = acquire(child.fComp);
1914 if (std::dynamic_pointer_cast<TH1>(cc) && !TString(cc->GetOption()).Contains("nostyle")) {
1915 xRooNode(out, *this).styles(cc.get()); // transfer style if adding a histogram
1916 }
1917 }
1918 if (!child.fComp && getObject<RooAbsReal>(child.GetName())) {
1919 Info("Add", "Adding existing function %s to %s", child.GetName(), p->GetName());
1920 out = getObject<RooAbsReal>(child.GetName());
1921 }
1922
1923 if (!out && !child.fComp) {
1924 std::shared_ptr<RooAbsArg> _func;
1925 // a null node .. so create either a new RooProduct or RooHistFunc if has observables (or no deps but has
1926 // x-axis)
1927 auto _obs = robs();
1928 if (!_obs.empty() || GetXaxis()) {
1929 if (_obs.empty()) {
1930 // using X axis to construct hist
1931 auto _ax = dynamic_cast<Axis2 *>(GetXaxis());
1932 auto t = TH1::AddDirectoryStatus();
1933 TH1::AddDirectory(false);
1934 auto h =
1935 std::make_unique<TH1D>(child.GetName(), child.GetTitle(), _ax->GetNbins(), _ax->binning()->array());
1937 h->GetXaxis()->SetName(TString::Format("%s;%s", _ax->GetParent()->GetName(), _ax->GetName()));
1938 // technically convertForAcquisition has already acquired so no need to re-acquire but should be harmless
1939 _func = std::dynamic_pointer_cast<RooAbsArg>(acquire(xRooNode(*h).convertForAcquisition(*this)));
1940 } else if (_obs.size() == 1) {
1941 // use the single obs to make a TH1D
1942 auto _x = _obs.at(0)->get<RooAbsLValue>();
1943 auto _bnames = _x->getBinningNames();
1944 TString binningName = p->getStringAttribute("binning");
1945 for (auto &b : _bnames) {
1946 if (b == p->GetName()) {
1947 binningName = p->GetName();
1948 break;
1949 }
1950 }
1951 auto t = TH1::AddDirectoryStatus();
1952 TH1::AddDirectory(false);
1953 auto h = std::make_unique<TH1D>(child.GetName(), child.GetTitle(), _x->numBins(binningName),
1954 _x->getBinningPtr(binningName)->array());
1956 h->GetXaxis()->SetName(
1957 TString::Format("%s;%s", dynamic_cast<TObject *>(_x)->GetName(), binningName.Data()));
1958 // technically convertForAcquisition has already acquired so no need to re-acquire but should be harmless
1959 _func = std::dynamic_pointer_cast<RooAbsArg>(acquire(xRooNode(*h).convertForAcquisition(*this)));
1960 Info("Add", "Created densityhisto factor %s (xaxis=%s) for %s", _func->GetName(), _obs.at(0)->GetName(),
1961 p->GetName());
1962 } else {
1963 throw std::runtime_error("Unsupported creation of new component in SumPdf for this many obs");
1964 }
1965 } else {
1966 _func = acquireNew<RooProduct>(TString::Format("%s_%s", p->GetName(), child.GetName()), child.GetTitle(),
1967 RooArgList());
1968 }
1969 _func->setStringAttribute("alias", child.GetName());
1970 out = _func;
1971 }
1972
1973 if (auto _f = std::dynamic_pointer_cast<RooHistFunc>(
1974 (child.get<RooProduct>()) ? child.factors()[child.GetName()]->fComp : out);
1975 _f) {
1976 // adding a histfunc directly to a sumpdf, should be a density
1977 _f->setAttribute("density");
1978 if (_f->getAttribute("autodensity")) {
1979 // need to divide by bin widths first
1980 for (int i = 0; i < _f->dataHist().numEntries(); i++) {
1981 auto bin_pars = _f->dataHist().get(i);
1982 _f->dataHist().set(*bin_pars, _f->dataHist().weight() / _f->dataHist().binVolume(*bin_pars));
1983 }
1984 _f->setAttribute("autodensity", false);
1985 _f->setValueDirty();
1986 }
1987
1988 // promote the axis vars to observables
1989 // can't use original child as might refer to unacquired deps
1990 for (auto &x : xRooNode("tmp", _f).vars()) {
1991 x->get<RooAbsArg>()->setAttribute("obs");
1992 }
1993 if (isConverted) {
1994 Info("Add", "Created %s factor RooHistFunc::%s for %s",
1995 _f->getAttribute("density") ? "densityhisto" : "histo", _f->GetName(), p->GetName());
1996 }
1997 }
1998
1999 if (auto _p = std::dynamic_pointer_cast<RooAbsPdf>(out); _p) {
2000 // adding a pdf to a RooRealSumPdf will replace it with a RooAddPdf and put the RooRealSumPdf inside that
2001 // if pdf is extended will use in the "no coefficients" state, where the expectedEvents are taking from
2002 // the pdf integrals
2003 TString newName(_p->GetName());
2004 newName.ReplaceAll("_samples", "");
2005 newName += "_components";
2006 Warning("Add", "converting samples to components");
2007
2008 if (auto _ax = GetXaxis(); _ax && dynamic_cast<RooAbsRealLValue *>(_ax->GetParent()) &&
2009 _p->dependsOn(*static_cast<RooAbsArg *>(_ax->GetParent()))) {
2010
2011 if (auto _boundaries = std::unique_ptr<std::list<double>>(_p->binBoundaries(
2012 *dynamic_cast<RooAbsRealLValue *>(_ax->GetParent()), -std::numeric_limits<double>::infinity(),
2013 std::numeric_limits<double>::infinity()));
2014 !_boundaries && _ax->GetNbins() > 0) {
2015#if ROOT_VERSION_CODE >= ROOT_VERSION(6, 24, 00)
2016 Warning("Add", "Adding unbinned pdf %s to binned %s - will wrap with RooBinSamplingPdf(...)",
2017 _p->GetName(), GetName());
2018 _p = acquireNew<RooBinSamplingPdf>(TString::Format("%s_binned", _p->GetName()), _p->GetTitle(),
2019 *dynamic_cast<RooAbsRealLValue *>(_ax->GetParent()), *_p);
2020 _p->setStringAttribute("alias", std::dynamic_pointer_cast<RooAbsArg>(out)->getStringAttribute("alias"));
2021 if (!_p->getStringAttribute("alias"))
2022 _p->setStringAttribute("alias", out->GetName());
2023#else
2024 throw std::runtime_error(
2025 "unsupported addition of unbinned pdf to binned model - please upgrade to at least ROOT 6.24");
2026#endif
2027 }
2028 }
2029
2030 // require to be extended to be in coefficient-free mode ...
2031 // otherwise would lose the integral of the sumPdf (can't think of way to have a coef be the integral)
2032 if (!_p->canBeExtended()) {
2033 _p = acquireNew<RooExtendPdf>(TString::Format("%s_extended", _p->GetName()), _p->GetTitle(), *_p,
2034 *acquire2<RooAbsReal, RooRealVar>("1", "1", 1));
2035 }
2036
2037 return *(Replace(*acquireNew<RooAddPdf>(newName, _p->GetTitle(), RooArgList(*p, *_p)))
2038 .browse()[1]); // returns second node.
2039 }
2040
2041 if (auto _f = std::dynamic_pointer_cast<RooAbsReal>(out); _f) {
2042
2043 // todo: if adding a pdf, should actually replace RooRealSumPdf with a RooAddPdf and put
2044 // the sumPdf and *this* pdf inside that pdf
2045 // only exception is the binSamplingPdf below to integrate unbinned functions across bins
2046
2047 if (auto _ax = GetXaxis(); _ax && dynamic_cast<RooAbsRealLValue *>(_ax->GetParent()) &&
2048 _f->dependsOn(*static_cast<RooAbsArg *>(_ax->GetParent()))) {
2049
2050 if (auto _boundaries = std::unique_ptr<std::list<double>>(_f->binBoundaries(
2051 *dynamic_cast<RooAbsRealLValue *>(_ax->GetParent()), -std::numeric_limits<double>::infinity(),
2052 std::numeric_limits<double>::infinity()));
2053 !_boundaries && _ax->GetNbins() > 0) {
2054#if ROOT_VERSION_CODE >= ROOT_VERSION(6, 24, 00)
2055 Warning(
2056 "Add",
2057 "Adding unbinned function %s to binned %s - will wrap with RooRealSumPdf(RooBinSamplingPdf(...))",
2058 _f->GetName(), GetName());
2059 auto sumPdf = acquireNew<RooRealSumPdf>(TString::Format("%s_pdfWrapper", _f->GetName()), _f->GetTitle(),
2060 *_f, *acquire2<RooAbsArg, RooRealVar>("1", "1", 1), true);
2061 sumPdf->setStringAttribute("alias", _f->getStringAttribute("alias"));
2062 if (!sumPdf->getStringAttribute("alias"))
2063 sumPdf->setStringAttribute("alias", out->GetName());
2064 _f = acquireNew<RooBinSamplingPdf>(TString::Format("%s_binned", _f->GetName()), _f->GetTitle(),
2065 *dynamic_cast<RooAbsRealLValue *>(_ax->GetParent()), *sumPdf);
2066 _f->setStringAttribute("alias", std::dynamic_pointer_cast<RooAbsArg>(out)->getStringAttribute("alias"));
2067 if (!_f->getStringAttribute("alias"))
2068 _f->setStringAttribute("alias", out->GetName());
2069#else
2070 throw std::runtime_error(
2071 "unsupported addition of unbinned function to binned model - please upgrade to at least ROOT 6.24");
2072#endif
2073 }
2074 }
2075
2076 const_cast<RooArgList &>(p->coefList()).add(*acquire2<RooAbsArg, RooRealVar>("1", "1", 1));
2077 const_cast<RooArgList &>(p->funcList()).add(*_f);
2078 // inherit binning if we dont have one yet
2079 if (!p->getStringAttribute("binning"))
2080 p->setStringAttribute("binning", _f->getStringAttribute("binning"));
2081
2082 xRooNode _out(_f, *this);
2083 if (auto gf = p->getStringAttribute("global_factors"); gf) {
2084 TStringToken pattern(gf, ";");
2085 while (pattern.NextToken()) {
2086 auto fac = getObject<RooAbsReal>(pattern.Data());
2087 if (!fac) {
2088 throw std::runtime_error(TString::Format("Could not find global factor %s", pattern.Data()));
2089 }
2090 _out.Multiply(fac);
2091 }
2092 }
2093 sterilize();
2094 // clear children for reload and update shared axis
2095 clear();
2096 fXAxis.reset();
2097 p->setStringAttribute("xvar", nullptr);
2098 browse();
2099 return _out;
2100 }
2101 } else if (auto p2 = get<RooProdPdf>(); p2) {
2102 // can "add" to a RooProdPdf provided trying to add a RooAbsReal not a RooAbsPdf and have a zero or 1
2103 // RooRealSumPdf child.convertForAcquisition(*this); - don't convert here because want generated objects named
2104 // after roorealsumpdf
2105 if (child.get<RooAbsPdf>() || (!child.get() && getObject<RooAbsPdf>(child.GetName()))) {
2106 // can add if 0 or 1 RooAddPdf ....
2107 RooAddPdf *_pdf = nullptr;
2108 bool tooMany(false);
2109 for (auto &pp : factors()) {
2110 if (auto _p = pp->get<RooAddPdf>(); _p) {
2111 if (_pdf) {
2112 _pdf = nullptr;
2113 tooMany = true;
2114 break;
2115 } // more than one!
2116 _pdf = _p;
2117 }
2118 }
2119 if (_pdf) {
2120 return xRooNode(*_pdf, *this).Add(child);
2121 } else if (!tooMany) {
2122 auto out = this->operator[]("components")->Add(child);
2123 return out;
2124 }
2125 } else if ((child.get<TH1>() || child.get<RooAbsReal>() ||
2126 (!child.get() && getObject<RooAbsReal>(child.GetName()))) &&
2127 !child.get<RooAbsPdf>()) {
2128 RooRealSumPdf *_pdf = nullptr;
2129 RooAddPdf *_backup = nullptr;
2130 bool tooMany(false);
2131 for (auto &pp : factors()) {
2132 if (auto _p = pp->get<RooRealSumPdf>(); _p) {
2133 if (_pdf) {
2134 _pdf = nullptr;
2135 tooMany = true;
2136 break;
2137 } // more than one!
2138 _pdf = _p;
2139 } else if (auto _p2 = pp->get<RooAddPdf>(); _p2) {
2140 _backup = _p2;
2141 for (auto &_pdfa : pp->components()) {
2142 if (auto _p3 = _pdfa->get<RooRealSumPdf>(); _p3) {
2143 if (_pdf) {
2144 _pdf = nullptr;
2145 tooMany = true;
2146 break;
2147 } // more than one!
2148 _pdf = _p3;
2149 }
2150 }
2151 }
2152 }
2153 if (_pdf) {
2154 return xRooNode(*_pdf, *this).Add(child);
2155 } else if (_backup) {
2156 // added *INSIDE* the addPdf -- will create a RooRealSumPdf to hold it
2157 return xRooNode(*_backup, *this).Add(child);
2158 } else if (!tooMany) {
2159 auto out = this->operator[]("samples")->Add(child);
2160 // clear our x-axis to re-evaluate
2161 fXAxis.reset();
2162 p2->setStringAttribute("xvar", nullptr);
2163 return out;
2164 }
2165 }
2166 } else if (auto s = get<RooSimultaneous>(); s) {
2167
2168 // adding to a simultaneous means adding a bin
2169 return bins().Add(child);
2170
2171 // if the child is a RooAbsPdf can just add it as a new channel using name of pdf as the channel name
2172 // if child is a histogram, will create a RooProdPdf
2173
2174 } else if (auto w = get<RooWorkspace>(); w) {
2175 child.convertForAcquisition(
2176 *this, child.get() ? "" : "func" /* if child is a string, allow it to be passed to factory */);
2177 if (child.get()) {
2178 if (auto _d = child.get<RooAbsData>()) {
2179 // don't use acquire method to import, because that adds datasets as Embeddded
2180 if (!w->import(*_d)) {
2181 // should upgrade vars with any obs from the dataset
2182 if (_d->get()) {
2183 std::unique_ptr<RooAbsCollection>(w->allVars().selectCommon(*_d->get()))->setAttribAll("obs");
2184 }
2185 if (_d->getGlobalObservables()) {
2186 std::unique_ptr<RooAbsCollection> globs(w->allVars().selectCommon(*_d->get()));
2187 globs->setAttribAll("obs");
2188 globs->setAttribAll("global");
2189 }
2190 return xRooNode(child.GetName(), *w->data(child.GetName()), *this);
2191 } else {
2192 throw std::runtime_error(
2193 TString::Format("Could not import dataset %s into workspace %s", child.GetName(), w->GetName())
2194 .Data());
2195 }
2196 } else {
2197 auto out = acquire(child.fComp);
2198 if (out)
2199 return xRooNode(child.GetName(), out, *this);
2200 }
2201 }
2202
2203 if (!child.empty() || child.fFolder == "!pdfs") {
2204 // create a RooSimultaneous using the children as the channels
2205 // children either have "=" in name if specifying channel cat name or otherwise assume
2206 std::string catName = "channelCat";
2207 if (!child.empty()) {
2208 if (TString ss = child.at(0)->GetName(); ss.Contains("=")) {
2209 catName = ss(0, ss.Index('='));
2210 }
2211 }
2212 auto _cat = acquire<RooCategory>(catName.c_str(), catName.c_str());
2213 _cat->setAttribute("obs");
2214 auto out = acquireNew<RooSimultaneous>(child.GetName(), child.GetTitle(), *_cat);
2215 Info("Add", "Created pdf RooSimultaneous::%s in workspace %s", out->GetName(), w->GetName());
2216 return xRooNode(out, *this);
2217 }
2218 }
2219
2220 if (sOpt == "pdf") {
2221 // can only add a pdf to a workspace
2222 if (get<RooWorkspace>()) {
2223 const_cast<xRooNode &>(child).fFolder = "!pdfs";
2224 return Add(child);
2225 }
2226 } else if (sOpt == "channel") {
2227 // can add to a model or to a workspace (creates a RooProdPdf either way)
2228 if (get<RooSimultaneous>()) {
2229 return Vary(child);
2230 } else if (get<RooWorkspace>()) {
2231 std::shared_ptr<TObject> out;
2232 child.convertForAcquisition(*this);
2233 if (child.get<RooAbsPdf>()) {
2234 out = acquire(child.fComp);
2235 } else if (!child.fComp) {
2236 out = acquireNew<RooProdPdf>(child.GetName(),
2237 (strlen(child.GetTitle())) ? child.GetTitle() : child.GetName(), RooArgList());
2238 Info("Add", "Created channel RooProdPdf::%s in workspace %s", out->GetName(), get()->GetName());
2239 }
2240 return xRooNode(out, *this);
2241 }
2242 } else if (sOpt == "sample" || sOpt == "func") {
2243 if (get<RooProdPdf>()) {
2244 auto _mainChild = mainChild();
2245 if (_mainChild.get<RooRealSumPdf>()) {
2246 return _mainChild.Add(child, sOpt == "func" ? "func" : "");
2247 } else {
2248 return (*this)["samples"]->Add(child, sOpt == "func" ? "func" : "");
2249 }
2250 }
2251 } else if (sOpt == "dataset") {
2252 if (get<RooWorkspace>()) {
2253 // const_cast<xRooNode&>(child).fFolder = "!datasets";return Add(child);
2254 return (*this).datasets().Add(child);
2255 }
2256 }
2257
2258 if (considerType) {
2259
2260 // interpret 'adding' here as dependent on the object type ...
2261 if (get<RooSimultaneous>()) {
2262 return bins().Add(child);
2263 } else if (TString(child.GetName()).Contains('=')) {
2264 return variations().Add(child);
2265 } else if (get<RooProduct>() || get<RooProdPdf>()) {
2266 return factors().Add(child);
2267 }
2268 }
2269
2270 // Nov 2022 - removed ability to add placeholders ... could bring back if rediscover need for them
2271 // if (!child.get() && child.empty() && strlen(child.GetName())) {
2272 // // can add a 'placeholder' node, note it will be deleted at the next browse
2273 // xRooNode out(child.GetName(),nullptr,*this);
2274 // out.SetTitle(child.GetTitle());
2275 // emplace_back(std::make_shared<xRooNode>(out));
2276 // // update the parent in the out node so that it's copy of the parent knows it has itself in it
2277 // // actually maybe not want this :-/
2278 // //out.fParent = std::make_shared<Node2>(*this);
2279 // for(auto o : *gROOT->GetListOfBrowsers()) {
2280 // if(auto b = dynamic_cast<TBrowser*>(o); b && b->GetBrowserImp()){
2281 // if(auto _b = dynamic_cast<TGFileBrowser*>(
2282 // dynamic_cast<TRootBrowser*>(b->GetBrowserImp())->fActBrowser ); _b) {
2283 // auto _root = _b->fRootDir;
2284 // if (!_root) _root = _b->fListTree->GetFirstItem();
2285 // if (auto item = _b->fListTree->FindItemByObj(_root,this); item) {
2286 // _b->fListTree->AddItem(item,back()->GetName(),back().get());
2287 // }
2288 // }
2289 // }
2290 // }
2291 // return out;
2292 // }
2293
2294 throw std::runtime_error(TString::Format("Cannot add %s to %s", child.GetName(), GetName()));
2295}
2296
2297std::string xRooNode::GetPath() const
2298{
2299 if (!fParent)
2300 return GetName();
2301 return fParent->GetPath() + "/" + GetName();
2302}
2303
2305{
2306 // std::cout << "deleting " << GetPath() << std::endl;
2307}
2308
2310{
2311 if (auto a = get<RooAbsArg>()) {
2312 a->setAttribute("hidden", set);
2313 // if(auto item = GetTreeItem(nullptr); item) {
2314 // if(set) item->SetColor(kRed);
2315 // else item->ClearColor();
2316 // }
2317 }
2318}
2320{
2321 auto a = get<RooAbsArg>();
2322 if (a)
2323 return a->getAttribute("hidden");
2324 return false;
2325}
2326
2328{
2329
2330 if (get() == rhs.get()) {
2331 // nothing to do because objects are identical
2332 return *this;
2333 }
2334
2335 // Info("Combine","Combining %s into %s",rhs.GetPath().c_str(),GetPath().c_str());
2336
2337 // combine components, factors, and variations ... when there is a name clash will combine on that object
2338 for (auto &c : rhs.components()) {
2339 if (auto _c = components().find(c->GetName()); _c) {
2340 _c->Combine(*c);
2341 } else {
2342 Add(*c);
2343 }
2344 }
2345
2346 for (auto &f : rhs.factors()) {
2347 if (auto _f = factors().find(f->GetName()); _f) {
2348 _f->Combine(*f);
2349 } else {
2350 Multiply(*f);
2351 }
2352 }
2353
2354 for (auto &v : rhs.variations()) {
2355 if (auto _v = variations().find(v->GetName()); _v) {
2356 _v->Combine(*v);
2357 } else {
2358 Vary(*v);
2359 }
2360 }
2361
2362 // todo: Should also transfer over binnings of observables
2363
2364 return *this;
2365}
2366
2367xRooNode xRooNode::shallowCopy(const std::string &name, std::shared_ptr<xRooNode> parent)
2368{
2369 xRooNode out(name.c_str(), nullptr,
2370 parent /*? parent : fParent -- was passing fParent for getObject benefit before fProvider concept*/);
2371 // if(!parent) out.fAcquirer = true;
2372 if (!parent)
2373 out.fProvider = fParent;
2374
2375 auto o = get();
2376 if (!o) {
2377 return out;
2378 }
2379
2380 if (auto s = get<RooSimultaneous>(); s) {
2381 auto chans = bins();
2382 if (!chans.empty()) {
2383 // create a new RooSimultaneous with shallow copies of each channel
2384
2385 std::shared_ptr<RooSimultaneous> pdf = out.acquire<RooSimultaneous>(
2386 name.c_str(), o->GetTitle(), const_cast<RooAbsCategoryLValue &>(s->indexCat()));
2387
2388 for (auto &c : chans) {
2389 TString cName(c->GetName());
2390 cName = cName(cName.Index('=') + 1, cName.Length());
2391 // by passing out as the parent, will ensure out acquires everything created
2392 auto c_copy =
2393 c->shallowCopy(name + "_" + c->get()->GetName(), std::shared_ptr<xRooNode>(&out, [](xRooNode *) {}));
2394 pdf->addPdf(*dynamic_cast<RooAbsPdf *>(c_copy.get()), cName);
2395 }
2396 out.fComp = pdf;
2397 return out;
2398 }
2399 } else if (auto p = dynamic_cast<RooProdPdf *>(o); p) {
2400 // main pdf will be copied too
2401 std::shared_ptr<RooProdPdf> pdf =
2402 std::dynamic_pointer_cast<RooProdPdf>(out.acquire(std::shared_ptr<TObject>(p->Clone(/*name.c_str()*/)), false,
2403 true)); // use clone to copy all attributes etc too
2404 auto main = mainChild();
2405 if (main) {
2406 auto newMain =
2407 std::dynamic_pointer_cast<RooAbsArg>(out.acquire(std::shared_ptr<TObject>(main->Clone()), false, true));
2408 std::cout << newMain << " " << newMain->GetName() << std::endl;
2409 // pdf->replaceServer(*pdf->pdfList().find(main->GetName()), *newMain, true, true);
2410 // const_cast<RooArgList&>(pdf->pdfList()).replace(*pdf->pdfList().find(main->GetName()), *newMain);
2411 pdf->redirectServers(RooArgList(*newMain));
2412 }
2413 out.fComp = pdf;
2414 out.sterilize();
2415 return out;
2416 }
2417
2418 return out;
2419}
2420
2422{
2423 static std::unique_ptr<cout_redirect> capture;
2424 std::string captureStr;
2425 bool doCapture = false;
2426 if (!capture && gROOT->FromPopUp()) { // FromPopUp means user executed from the context menu
2427 capture = std::make_unique<cout_redirect>(captureStr);
2428 doCapture = true;
2429 }
2430
2431 TString sOpt(opt);
2432 int depth = 0;
2433 if (sOpt.Contains("depth=")) {
2434 depth = TString(sOpt(sOpt.Index("depth=") + 6, sOpt.Length())).Atoi();
2435 sOpt.ReplaceAll(TString::Format("depth=%d", depth), "");
2436 }
2437 int indent = 0;
2438 if (sOpt.Contains("indent=")) {
2439 indent = TString(sOpt(sOpt.Index("indent=") + 7, sOpt.Length())).Atoi();
2440 sOpt.ReplaceAll(TString::Format("indent=%d", indent), "");
2441 }
2442 bool _more = sOpt.Contains("m");
2443 if (_more)
2444 sOpt.Replace(sOpt.Index("m"), 1, "");
2445 if (sOpt != "")
2446 _more = true;
2447
2448 if (indent == 0) { // only print self if not indenting (will already be printed above if tree traverse)
2449 std::cout << GetPath();
2450 if (get() && get() != this) {
2451 std::cout << ": ";
2452 if (_more || (get<RooAbsArg>() && get<RooAbsArg>()->isFundamental()) || get<RooConstVar>() ||
2453 get<RooAbsData>() || get<RooProduct>() || get<RooFitResult>()) {
2454 auto _deps = coords(false).argList(); // want to revert coords after print
2455 auto _snap = std::unique_ptr<RooAbsCollection>(_deps.snapshot());
2456 coords(); // move to coords before printing (in case this matters)
2457 get()->Print(sOpt);
2458 if (auto _fr = get<RooFitResult>(); _fr && dynamic_cast<RooStringVar *>(_fr->constPars().find(".log"))) {
2459 std::cout << "Minimization Logs:" << std::endl;
2460 std::cout << dynamic_cast<RooStringVar *>(_fr->constPars().find(".log"))->getVal() << std::endl;
2461 }
2462 _deps.assignValueOnly(*_snap);
2463 // std::cout << std::endl;
2464 } else {
2465 TString _suffix = "";
2466 if (auto _type = GetNodeType(); strlen(_type)) {
2467 // decided not to show const values until figure out how to update if value changes
2468 /*if (TString(_type)=="Const") _name += TString::Format("
2469 [%s=%g]",_type,v->get<RooConstVar>()->getVal()); else*/
2470 _suffix += TString::Format(" [%s]", _type);
2471 }
2472 if (auto fv = get<RooFormulaVar>()) {
2473 TString formu = TString::Format(" [%s]", fv->expression());
2474 for (size_t i = 0; i < fv->dependents().size(); i++) {
2475 formu.ReplaceAll(TString::Format("x[%zu]", i), fv->dependents()[i].GetName());
2476 }
2477 _suffix += formu;
2478 } else if (auto gv = get<RooGenericPdf>()) {
2479 TString formu = TString::Format(" [%s]", gv->expression());
2480 for (size_t i = 0; i < gv->dependents().size(); i++) {
2481 formu.ReplaceAll(TString::Format("x[%zu]", i), gv->dependents()[i].GetName());
2482 }
2483 _suffix += formu;
2484 }
2485 std::cout << get()->ClassName() << "::" << get()->GetName() << _suffix.Data() << std::endl;
2486 }
2487
2488 } else if (!get()) {
2489 std::cout << std::endl;
2490 }
2491 }
2492 const_cast<xRooNode *>(this)->browse();
2493 std::vector<std::string> folderNames;
2494 for (auto &k : *this) {
2495 if (std::find(folderNames.begin(), folderNames.end(), k->fFolder) == folderNames.end()) {
2496 folderNames.push_back(k->fFolder);
2497 }
2498 }
2499 for (auto &f : folderNames) {
2500 int i = 0;
2501 int iindent = indent;
2502 if (!f.empty()) {
2503 for (int j = 0; j < indent; j++)
2504 std::cout << " ";
2505 std::cout << f << std::endl;
2506 iindent += 1;
2507 }
2508 for (auto &k : *this) {
2509 if (k->fFolder != f) {
2510 i++;
2511 continue;
2512 }
2513 for (int j = 0; j < iindent; j++)
2514 std::cout << " ";
2515 std::cout << i++ << ") " << k->GetName() << " : ";
2516 if (k->get()) {
2517 if (_more || (k->get<RooAbsArg>() && k->get<RooAbsArg>()->isFundamental()) || k->get<RooConstVar>() ||
2518 k->get<RooAbsData>() /*|| k->get<RooProduct>()*/) {
2519 auto _deps = k->coords(false).argList();
2520 auto _snap = std::unique_ptr<RooAbsCollection>(_deps.snapshot());
2521 k->coords(); // move to coords before printing (in case this matters)
2522 k->get()->Print(sOpt); // assumes finishes with an endl
2523 _deps.assignValueOnly(*_snap);
2524 } else {
2525 TString _suffix = "";
2526 if (auto _type = k->GetNodeType(); strlen(_type)) {
2527 // decided not to show const values until figure out how to update if value changes
2528 /*if (TString(_type)=="Const") _name += TString::Format("
2529 [%s=%g]",_type,v->get<RooConstVar>()->getVal()); else*/
2530 _suffix += TString::Format(" [%s]", _type);
2531 }
2532 if (auto fv = k->get<RooFormulaVar>()) {
2533 TString formu = TString::Format(" [%s]", fv->expression());
2534 for (size_t j = 0; j < fv->dependents().size(); j++) {
2535 formu.ReplaceAll(TString::Format("x[%zu]", j), fv->dependents()[j].GetName());
2536 }
2537 _suffix += formu;
2538 } else if (auto gv = k->get<RooGenericPdf>()) {
2539 TString formu = TString::Format(" [%s]", gv->expression());
2540 for (size_t j = 0; j < gv->dependents().size(); j++) {
2541 formu.ReplaceAll(TString::Format("x[%zu]", j), gv->dependents()[j].GetName());
2542 }
2543 _suffix += formu;
2544 }
2545 std::cout << k->get()->ClassName() << "::" << k->get()->GetName() << _suffix.Data() << std::endl;
2546 }
2547 if (depth != 0) {
2548 k->Print(sOpt + TString::Format("depth=%dindent=%d", depth - 1, iindent + 1));
2549 }
2550 } else
2551 std::cout << " NULL " << std::endl;
2552 }
2553 }
2554 if (doCapture) {
2555 capture.reset(); // no captureStr has the string to display
2556 // inject line breaks to avoid msgbox being too wide
2557 size_t lastBreak = 0;
2558 std::string captureStrWithBreaks;
2559 for (size_t i = 0; i < captureStr.size(); i++) {
2560 captureStrWithBreaks += captureStr[i];
2561 if (captureStr[i] == '\n') {
2562 lastBreak = i;
2563 }
2564 if (i - lastBreak > 150) {
2565 captureStrWithBreaks += '\n';
2566 lastBreak = i;
2567 }
2568 }
2569 const TGWindow *w =
2570 (gROOT->GetListOfBrowsers()->At(0))
2571 ? dynamic_cast<TGWindow *>(static_cast<TBrowser *>(gROOT->GetListOfBrowsers()->At(0))->GetBrowserImp())
2572 : gClient->GetRoot();
2573 new TGMsgBox(gClient->GetRoot(), w, GetName(),
2574 captureStrWithBreaks.c_str()); //,nullptr,kMBDismiss,nullptr,kVerticalFrame,kTextLeft|kTextCenterY);
2575 }
2576}
2577
2579{
2580 if (!child.get()) {
2581
2582 if (auto v = get<RooRealVar>(); v) {
2583
2584 TString constrType = child.GetName();
2585 double mean = std::numeric_limits<double>::quiet_NaN();
2586 double sigma = mean;
2587 if (constrType.BeginsWith("gaussian(")) {
2588 // extract the mean and stddev parameters
2589 // if only one given, it is the stddev
2590 if (constrType.Contains(",")) {
2591 mean = TString(constrType(9, constrType.Index(',') - 9)).Atof();
2592 sigma = TString(constrType(constrType.Index(',') + 1, constrType.Index(')') - constrType.Index(',') + 1))
2593 .Atof();
2594 } else {
2595 mean = std::numeric_limits<double>::quiet_NaN(); // will use the var current value below to set mean
2596 sigma = TString(constrType(9, constrType.Index(')') - 9)).Atof();
2597 }
2598 constrType = "normal";
2599 } else if (constrType == "normal") {
2600 mean = 0;
2601 sigma = 1;
2602 } else if (constrType == "gaussian") {
2603 // extract parameters from the variable
2604 // use current value and error on v as constraint
2605 if (!v->hasError())
2606 throw std::runtime_error("No error on parameter for gaussian constraint");
2607 sigma = v->getError();
2608 mean = v->getVal();
2609 constrType = "normal";
2610 } else if (constrType == "poisson") {
2611 if (!v->hasError())
2612 throw std::runtime_error("No error on parameter for poisson constraint");
2613 mean = 1;
2614 sigma = pow(v->getVal() / v->getError(), 2);
2615 }
2616
2617 if (constrType == "poisson") {
2618 // use current value and error on v as constraint
2619 double tau_val = sigma;
2620 auto globs = acquire<RooRealVar>(Form("globs_%s", v->GetName()), Form("globs_%s", v->GetName()),
2621 v->getVal() * tau_val, (v->getVal() - 5 * v->getError()) * tau_val,
2622 (v->getVal() + 5 * v->getError()) * tau_val);
2623 globs->setConstant();
2624 globs->setAttribute("obs");
2625 globs->setAttribute("global");
2626 globs->setStringAttribute("nominal", TString::Format("%f", tau_val));
2627 auto tau = acquireNew<RooConstVar>(TString::Format("tau_%s", v->GetName()), "", tau_val);
2628 auto constr = acquireNew<RooPoisson>(
2629 Form("pois_%s", v->GetName()), TString::Format("Poisson Constraint of %s", v->GetTitle()), *globs,
2630 *acquireNew<RooProduct>(TString::Format("mean_%s", v->GetName()),
2631 TString::Format("Poisson Constraint of %s", globs->GetTitle()),
2632 RooArgList(*v, *tau)),
2633 true /* no rounding */);
2634
2635 auto out = Constrain(xRooNode(Form("pois_%s", GetName()), constr));
2636 if (!v->hasError())
2637 v->setError(mean / sqrt(tau_val)); // if v doesnt have an uncert, will put one on it now
2638 Info("Constrain", "Added poisson constraint pdf RooPoisson::%s (tau=%g) for %s", out->GetName(), tau_val,
2639 GetName());
2640 return out;
2641 } else if (constrType == "normal") {
2642
2643 auto globs = acquire<RooRealVar>(Form("globs_%s", v->GetName()), Form("globs_%s", v->GetName()), mean,
2644 mean - 10 * sigma, mean + 10 * sigma);
2645 globs->setAttribute("obs");
2646 globs->setAttribute("global");
2647 globs->setConstant();
2648
2649 globs->setStringAttribute("nominal", TString::Format("%f", mean));
2650 auto constr = acquireNew<RooGaussian>(
2651 Form("gaus_%s", v->GetName()), TString::Format("Gaussian Constraint of %s", v->GetTitle()), *globs, *v,
2652 *acquireNew<RooConstVar>(TString::Format("sigma_%s", v->GetName()), "", sigma));
2653 auto out = Constrain(xRooNode(Form("gaus_%s", GetName()), constr));
2654 if (!v->hasError())
2655 v->setError(sigma); // if v doesnt have an uncert, will put one on it now
2656 Info("Constrain", "Added gaussian constraint pdf RooGaussian::%s (mean=%g,sigma=%g) for %s", out->GetName(),
2657 mean, sigma, GetName());
2658 return out;
2659 }
2660 }
2661 } else if (auto p = child.get<RooAbsPdf>(); p) {
2662
2663 auto _me = get<RooAbsArg>();
2664 if (!_me) {
2665 throw std::runtime_error("Cannot constrain non arg");
2666 }
2667
2668 if (!p->dependsOn(*_me)) {
2669 throw std::runtime_error("Constraint does not depend on constrainee");
2670 }
2671
2672 // find a parent that can swallow this pdf ... either a RooProdPdf or a RooWorkspace
2673 auto x = fParent;
2674 while (x && !x->get<RooProdPdf>() && !x->get<RooSimultaneous>() && !x->get<RooWorkspace>()) {
2675 x = x->fParent;
2676 }
2677 if (!x) {
2678 throw std::runtime_error("Nowhere to put constraint");
2679 }
2680 // get datasets of the swallower, and add glob to any globs lists
2681 auto childGlobs = child.globs();
2682 if (!childGlobs.empty()) {
2683 for (auto d : x->datasets()) {
2684 if (auto globs = d->get<RooAbsData>()->getGlobalObservables()) {
2685 RooArgSet newGlobs(*globs);
2686 newGlobs.add(*childGlobs.get<RooArgList>());
2687 d->get<RooAbsData>()->setGlobalObservables(newGlobs);
2688 }
2689 }
2690 // also add to the workspaces globalObservables lists
2691 if (x->ws()) {
2692 for (auto &[k, v] : GETWSSETS(x->ws())) {
2693 if (k == "globalObservables" || TString(k).EndsWith("_GlobalObservables")) {
2694 const_cast<RooArgSet &>(v).add(*childGlobs.get<RooArgList>());
2695 }
2696 }
2697 }
2698 }
2699 if (auto s = x->get<RooSimultaneous>(); s) {
2700 // put into every channel that features parameter
2701 x->browse();
2702 for (auto &c : *x) {
2703 if (auto a = c->get<RooAbsArg>(); a->dependsOn(*_me))
2704 c->Multiply(child);
2705 }
2706 return child;
2707 } else if (x->get<RooProdPdf>()) {
2708 return x->Multiply(child);
2709 } else {
2710 return x->Add(child, "+");
2711 }
2712 }
2713
2714 throw std::runtime_error(TString::Format("Cannot constrain %s", GetName()));
2715}
2716
2718{
2719
2720 class AutoUpdater {
2721 public:
2722 AutoUpdater(xRooNode &_n) : n(_n) {}
2723 ~AutoUpdater() { n.browse(); }
2724 xRooNode &n;
2725 };
2726 AutoUpdater xxx(*this);
2727
2728 if (fBinNumber != -1) {
2729 // scaling a bin ...
2730 if (child.get<RooAbsReal>()) { // if not child then let fall through to create a child and call self again below
2731 // doing a bin-multiplication .. the parent should have a ParamHistFunc called binFactors
2732 // if it doesn't then create one
2733 auto o = std::dynamic_pointer_cast<RooAbsReal>(acquire(child.fComp));
2734
2735 // get binFactor unless parent is a ParamHistFunc already ...
2736
2737 auto binFactors = (fParent->get<ParamHistFunc>()) ? fParent : fParent->factors().find("binFactors");
2738
2739 // it can happen in a loop over bins() that another node has moved fParent inside a product
2740 // so check for fParent having a client with the ORIGNAME:<name> attribute
2741 if (!binFactors && fParent->get<RooAbsArg>()) {
2742 for (auto c : fParent->get<RooAbsArg>()->clients()) {
2743 if (c->IsA() == RooProduct::Class() &&
2744 c->getAttribute(TString::Format("ORIGNAME:%s", fParent->get()->GetName()))) {
2745 // try getting binFactors out of this
2746 binFactors = xRooNode(*c).factors().find("binFactors");
2747 break;
2748 }
2749 }
2750 }
2751
2752 if (!binFactors) {
2753 fParent
2754 ->Multiply(TString::Format("%s_binFactors",
2755 (fParent->mainChild().get())
2756 ? fParent->mainChild()->GetName()
2757 : (fParent->get() ? fParent->get()->GetName() : fParent->GetName()))
2758 .Data(),
2759 "blankshape")
2760 .SetName("binFactors"); // creates ParamHistFunc with all pars = 1 (shared const)
2761 binFactors = fParent->factors().find("binFactors");
2762 if (!binFactors) {
2763 throw std::runtime_error(
2764 TString::Format("Could not create binFactors in parent %s", fParent->GetName()));
2765 }
2766 // auto phf = binFactors->get<ParamHistFunc>();
2767
2768 // create RooProducts for all the bins ... so that added factors don't affect selves
2769 int i = 1;
2770 for (auto &b : binFactors->bins()) {
2771 auto p = acquireNew<RooProduct>(TString::Format("%s_bin%d", binFactors->get()->GetName(), i),
2772 TString::Format("binFactors of bin %d", i), RooArgList());
2773 p->setStringAttribute("alias", TString::Format("%s=%g", binFactors->GetXaxis()->GetParent()->GetName(),
2774 binFactors->GetXaxis()->GetBinCenter(i)));
2775 b->Multiply(*p);
2776 i++;
2777 }
2778 }
2779 // then scale the relevant bin ... if the relevant bin is a "1" then just drop in our factor (inside a
2780 // RooProduct though, to avoid it getting modified by subsequent multiplies)
2781 auto _bin = binFactors->bins().at(fBinNumber - 1);
2782 if (auto phf = binFactors->get<ParamHistFunc>(); phf && _bin) {
2783#if ROOT_VERSION_CODE < ROOT_VERSION(6, 27, 00)
2784 RooArgList &pSet = phf->_paramSet;
2785#else
2786 RooArgList &pSet = const_cast<RooArgList &>(phf->paramList());
2787#endif
2788 if (strcmp(_bin->GetName(), "1") == 0) {
2789 RooArgList all;
2790 for (std::size_t i = 0; i < pSet.size(); i++) {
2791 if (int(i) != fBinNumber - 1) {
2792 all.add(*pSet.at(i));
2793 } else {
2794 all.add(*o);
2795 }
2796 }
2797 pSet.removeAll();
2798 pSet.add(all);
2799 } else {
2800 _bin->fBinNumber = -1; // to avoid infinite loop
2801 return _bin->Multiply(child, opt);
2802 }
2803 // } else {else if(_bin->get<RooProduct>()) {
2804 // // multiply the element which will just add it as a factor in the rooproduct
2805 // return _bin->Multiply(child,opt);
2806 // } else {
2807 // // not a rooproduct in this bin yet ... so need to replace with a rooproduct and
2808 // multiply that
2809 // // this avoids the undesired behaviour of shared binFactors getting all impacted by
2810 // mulitplies RooArgList all; auto new_p =
2811 // acquireNew<RooProduct>(TString::Format("%s_bin%d",binFactors->get()->GetName(),fBinNumber),TString::Format("binFactors
2812 // of bin %d",fBinNumber),RooArgList(*_bin->get<RooAbsArg>()));
2813 // new_p->setStringAttribute("alias","")
2814 // for (int i = 0; i < phf->_paramSet.size(); i++) {
2815 // if (i != fBinNumber - 1) all.add(*phf->_paramSet.at(i));
2816 // else all.add(*new_p);
2817 // }
2818 // phf->_paramSet.removeAll();
2819 // phf->_paramSet.add(all);
2820 // // now multiply that bin having converted it to RooProduct
2821 // return binFactors->bins().at(fBinNumber - 1)->Multiply(child,opt);
2822 // }
2823 }
2824 return xRooNode(*o, binFactors);
2825 }
2826 } else if (!get() && fParent) {
2827 // try to 'create' object based on parentage
2828 // add child as a temporary child to help with decision making
2829 auto _ref = emplace_back(std::shared_ptr<xRooNode>(&const_cast<xRooNode &>(child), [](TObject *) {}));
2830 try {
2831 fComp = fParent->Add(*this, "+").fComp;
2832 } catch (...) {
2833 resize(size() - 1);
2834 std::rethrow_exception(std::current_exception());
2835 }
2836 resize(size() - 1); // remove the temporarily added node
2837 }
2838
2839 if (!child.get()) {
2840 TString sOpt(opt);
2841 sOpt.ToLower();
2842 if (auto o = getObject<RooAbsReal>(child.GetName())) {
2843 auto out = Multiply(xRooNode(o, child.fParent));
2844 // have to protect bin case where get() is null (could change but then must change logic above too)
2845 if (get()) {
2846 Info("Multiply", "Scaled %s by existing factor %s::%s",
2847 mainChild().get() ? mainChild().get()->GetName() : get()->GetName(), o->ClassName(), o->GetName());
2848 }
2849 return out;
2850 } else if (sOpt == "norm") {
2851 if (TString(child.GetName()).Contains("[") && ws()) {
2852 // assume factory method wanted
2853 auto arg = ws()->factory(child.GetName());
2854 if (arg) {
2855 auto out = Multiply(*arg);
2856 if (get()) {
2857 Info("Multiply", "Scaled %s by new norm factor %s",
2858 mainChild().get() ? mainChild().get()->GetName() : get()->GetName(), out->GetName());
2859 }
2860 return out;
2861 }
2862 throw std::runtime_error(TString::Format("Failed to create new normFactor %s", child.GetName()));
2863 }
2864 auto out = Multiply(RooRealVar(child.GetName(), child.GetTitle(), 1, -1e-5, 100));
2865 if (get()) {
2866 Info("Multiply", "Scaled %s by new norm factor %s",
2867 mainChild().get() ? mainChild().get()->GetName() : get()->GetName(), out->GetName());
2868 }
2869 return out;
2870 } else if (sOpt == "shape" || sOpt == "histo" || sOpt == "blankshape") {
2871 // needs axis defined
2872 if (auto ax = GetXaxis(); ax) {
2873 auto h = std::shared_ptr<TH1>(BuildHistogram(dynamic_cast<RooAbsLValue *>(ax->GetParent()), true));
2874 h->Reset();
2875 for (int i = 1; i <= h->GetNbinsX(); i++) {
2876 h->SetBinContent(i, 1);
2877 }
2878 h->SetMinimum(0);
2879 h->SetMaximum(100);
2880 h->SetName(TString::Format(";%s", child.GetName())); // ; char indicates don't "rename" this thing
2881 h->SetTitle(child.GetTitle());
2882 if (sOpt.Contains("shape"))
2883 h->SetOption(sOpt);
2884 auto out = Multiply(*h);
2885 if (get()) {
2886 Info("Multiply", "Scaled %s by new %s factor %s",
2887 mainChild().get() ? mainChild().get()->GetName() : get()->GetName(), sOpt.Data(), out->GetName());
2888 }
2889 return out;
2890 }
2891 } else if (sOpt == "overall") {
2892 auto out = Multiply(acquireNew<RooStats::HistFactory::FlexibleInterpVar>(
2893 child.GetName(), child.GetTitle(), RooArgList(), 1, std::vector<double>(), std::vector<double>()));
2894 if (get() /* can happen this is null if on a bin node with no shapeFactors*/) {
2895 Info("Multiply", "Scaled %s by new overall factor %s",
2896 mainChild().get() ? mainChild().get()->GetName() : get()->GetName(), out->GetName());
2897 }
2898 return out;
2899 } else if (sOpt == "func" && ws()) {
2900 // need to get way to get dependencies .. can't pass all as causes circular dependencies issues.
2901 if (auto arg = ws()->factory(TString("expr::") + child.GetName())) {
2902 auto out = Multiply(*arg);
2903 if (get() /* can happen this is null if on a bin node with no shapeFactors*/) {
2904 Info("Multiply", "Scaled %s by new func factor %s",
2905 mainChild().get() ? mainChild().get()->GetName() : get()->GetName(), out->GetName());
2906 }
2907 return out;
2908 }
2909 }
2910 }
2911 if (auto h = child.get<TH1>(); h && strlen(h->GetOption()) == 0 && strlen(opt) > 0) {
2912 // put the option in the hist
2913 h->SetOption(opt);
2914 }
2915 if (auto w = get<RooWorkspace>(); w) {
2916 // just acquire
2917 std::shared_ptr<TObject> out;
2918 child.convertForAcquisition(*this);
2919 if (child.get<RooAbsReal>())
2920 out = acquire(child.fComp);
2921 return out;
2922 }
2923
2924 if (strcmp(GetName(), ".coef") == 0) { // covers both .coef and .coefs
2925 // need to add this into the relevant coef ... if its not a RooProduct, replace it with one first
2926 if (auto p = fParent->fParent->get<RooAddPdf>()) {
2927 // may be in no-coef mode ... in which case must create coefs (use "ExtendedBindings" but note that these need
2928 // obs list passing to them
2929 if (p->coefList().empty() && !p->pdfList().empty()) {
2930 for (auto _pdf : p->pdfList()) {
2931 const_cast<RooArgList &>(p->coefList())
2932 .add(*acquireNew<RooExtendedBinding>(TString::Format("%s_extBind", _pdf->GetName()),
2933 TString::Format("Expected Events of %s", _pdf->GetTitle()),
2934 *static_cast<RooAbsPdf *>(_pdf)));
2935 }
2936 Info("Multiply", "Created RooExtendedBinding coefficients for all pdfs of %s so that can multiply coef",
2937 p->GetName());
2938 *reinterpret_cast<bool *>(reinterpret_cast<unsigned char *>(p) +
2939 p->Class()->GetDataMemberOffset("_allExtendable")) = false;
2940 *reinterpret_cast<bool *>(reinterpret_cast<unsigned char *>(p) +
2941 p->Class()->GetDataMemberOffset("_haveLastCoef")) = true;
2942 }
2943 for (size_t i = 0; i < p->pdfList().size(); i++) {
2944 if (p->pdfList().at(i) == fParent->get<RooAbsArg>()) {
2945 auto coefs = p->coefList().at(i);
2946 if (!coefs->InheritsFrom("RooProduct")) {
2947 RooArgList oldCoef;
2948 if (!(strcmp(coefs->GetName(), "1") == 0 || strcmp(coefs->GetName(), "ONE") == 0))
2949 oldCoef.add(*coefs);
2950 auto newCoefs = fParent->acquireNew<RooProduct>(
2951 TString::Format("coefs_%s", fParent->GetName()),
2952 TString::Format("coefficients for %s", fParent->GetName()), oldCoef);
2953 RooArgList oldCoefs;
2954 for (size_t j = 0; j < p->coefList().size(); j++) {
2955 if (i == j) {
2956 oldCoefs.add(*newCoefs);
2957 } else {
2958 oldCoefs.add(*p->coefList().at(j));
2959 }
2960 }
2961 const_cast<RooArgList &>(p->coefList()).removeAll();
2962 const_cast<RooArgList &>(p->coefList()).add(oldCoefs);
2963 coefs = newCoefs.get();
2964 }
2965 return xRooNode(*coefs, fParent).Multiply(child);
2966 }
2967 }
2968 }
2969 throw std::runtime_error("this coefs case is not supported");
2970 }
2971
2972 if (auto p = get<RooProduct>(); p) {
2973 std::shared_ptr<TObject> out;
2974 auto cc = child.fComp;
2975 bool isConverted = (child.convertForAcquisition(*this) != cc);
2976 if (child.get<RooAbsReal>())
2977 out = acquire(child.fComp);
2978
2979 // child may be a histfunc or a rooproduct of a histfunc and a paramhist if has stat errors
2980 if (auto _f = std::dynamic_pointer_cast<RooHistFunc>(
2981 (child.get<RooProduct>()) ? child.factors()[child.GetName()]->fComp : out);
2982 _f && _f->getAttribute("autodensity")) {
2983 // should we flag this as a density? yes if there's no other term marked as the density
2984 bool hasDensity = false;
2985 for (auto &f : factors()) {
2986 if (f->get<RooAbsArg>()->getAttribute("density")) {
2987 hasDensity = true;
2988 break;
2989 }
2990 }
2991 _f->setAttribute("density", !hasDensity && fParent && fParent->get<RooRealSumPdf>());
2992 if (_f->getAttribute("density")) {
2993
2994 // need to divide by bin widths first
2995 for (int i = 0; i < _f->dataHist().numEntries(); i++) {
2996 auto bin_pars = _f->dataHist().get(i);
2997 _f->dataHist().set(*bin_pars, _f->dataHist().weight() / _f->dataHist().binVolume(*bin_pars));
2998 }
2999 _f->setValueDirty();
3000
3001 // promote the axis vars to observables
3002 for (auto &x : xRooNode("tmp", _f).vars()) {
3003 x->get<RooAbsArg>()->setAttribute("obs");
3004 }
3005 }
3006 _f->setAttribute("autodensity", false);
3007 }
3008
3009 if (isConverted && child.get<RooHistFunc>()) {
3010 Info("Multiply", "Created %s factor %s in %s",
3011 child.get<RooAbsArg>()->getAttribute("density") ? "densityhisto" : "histo", child->GetName(),
3012 p->GetName());
3013 } else if (isConverted && child.get<ParamHistFunc>()) {
3014 Info("Multiply", "Created shape factor %s in %s", child->GetName(), p->GetName());
3015 }
3016
3017 if (auto _f = std::dynamic_pointer_cast<RooAbsReal>(out); _f) {
3018#if ROOT_VERSION_CODE < ROOT_VERSION(6, 27, 00)
3019 p->_compRSet.add(*_f);
3020#else
3021 const_cast<RooArgList &>(p->realComponents()).add(*_f);
3022#endif
3023 p->setValueDirty();
3024
3025 browse();
3026 xRooNode _out(_f, *this);
3027 for (auto &_par : _out.pars()) {
3028 if (auto s = _par->get<RooAbsArg>()->getStringAttribute("boundConstraint"); s) {
3029 bool found = false;
3030 for (auto &_constr : _par->constraints()) {
3031 if (strcmp(s, _constr->get()->GetName()) == 0) {
3032 // constraint is already included
3033 found = true;
3034 break;
3035 }
3036 }
3037 if (!found) {
3038 Info("Multiply", "Pulling in %s boundConstraint: %s", _par->GetName(), s);
3039 auto _pdf = getObject<RooAbsPdf>(s);
3040 if (!_pdf) {
3041 throw std::runtime_error("Couldn't find boundConstraint");
3042 }
3043 _par->Constrain(_pdf);
3044 }
3045 }
3046 }
3047 sterilize();
3048 return _out;
3049 }
3050 } else if (auto p2 = get<RooProdPdf>(); p2) {
3051
3052 std::shared_ptr<TObject> out;
3053 child.convertForAcquisition(*this);
3054 if (child.get<RooAbsPdf>()) {
3055 out = acquire(child.fComp);
3056 } else if (child.get<RooAbsReal>() && mainChild().get<RooRealSumPdf>()) {
3057 // cannot multiply a RooProdPdf by a non pdf
3058 throw std::runtime_error(TString::Format("Cannot multiply %s by non-pdf %s", GetName(), child.GetName()));
3059 // return mainChild().Add(child); - nov 2022 - used to do this but now replaced with exception above
3060 } else if (!child.get() || child.get<RooAbsReal>()) {
3061 // need to create or hide inside a sumpdf or rooadpdf
3062 std::shared_ptr<RooAbsPdf> _pdf;
3063 if (!child.get() && strcmp(child.GetName(), "components") == 0) {
3064 auto _sumpdf = acquireNew<RooAddPdf>(
3065 Form("%s_%s", p2->GetName(), child.GetName()),
3066 (strlen(child.GetTitle()) && strcmp(child.GetTitle(), child.GetName())) ? child.GetTitle()
3067 : p2->GetTitle(),
3068 RooArgList() /*, RooArgList() forces coef-mode if we specify this list */);
3069 _pdf = _sumpdf;
3070 } else {
3071 auto _sumpdf = acquireNew<RooRealSumPdf>(
3072 Form("%s_%s", p2->GetName(), child.GetName()),
3073 (strlen(child.GetTitle()) && strcmp(child.GetTitle(), child.GetName())) ? child.GetTitle()
3074 : p2->GetTitle(),
3075 RooArgList(), RooArgList(), true);
3076 _sumpdf->setFloor(true);
3077 _pdf = _sumpdf;
3078 }
3079 _pdf->setStringAttribute("alias", child.GetName());
3080 // transfer axis attributes if present (TODO: should GetXaxis look beyond the immediate parent?)
3081 _pdf->setStringAttribute("xvar", p2->getStringAttribute("xvar"));
3082 _pdf->setStringAttribute("binning", p2->getStringAttribute("binning"));
3083 out = _pdf;
3084 Info("Multiply", "Created %s::%s in channel %s", _pdf->ClassName(), _pdf->GetName(), p2->GetName());
3085 if (child.get<RooAbsReal>())
3086 xRooNode(*out, *this).Add(child);
3087 }
3088
3089 if (auto _pdf = std::dynamic_pointer_cast<RooAbsPdf>(out); _pdf) {
3090#if ROOT_VERSION_CODE < ROOT_VERSION(6, 27, 00)
3091 const_cast<RooArgList &>(p2->pdfList()).add(*_pdf);
3092#if ROOT_VERSION_CODE >= ROOT_VERSION(6, 26, 00)
3093 p2->_pdfNSetList.emplace_back(std::make_unique<RooArgSet>("nset"));
3094#else
3095 p->_pdfNSetList.Add(new RooArgSet("nset"));
3096#endif
3097 if (!p2->canBeExtended() && _pdf->canBeExtended()) {
3098 p2->_extendedIndex = p2->_pdfList.size() - 1;
3099 }
3100#else
3101 p2->addPdfs(RooArgSet(*_pdf));
3102#endif
3103 sterilize();
3104 browse();
3105 return xRooNode(_pdf, *this);
3106 }
3107 } else if (auto p3 = get<RooRealSumPdf>(); p3) {
3108 // multiplying all current and future components
3109 std::shared_ptr<TObject> out;
3110 child.convertForAcquisition(*this);
3111 if (child.get<RooAbsReal>()) {
3112 out = acquire(child.fComp);
3113 for (auto &c : components()) {
3114 c->Multiply(out);
3115 }
3116 TString s = p3->getStringAttribute("global_factors");
3117 if (s != "")
3118 s += ";";
3119 s += out->GetName();
3120 p3->setStringAttribute("global_factors", s);
3121 Info(
3122 "Multiply",
3123 "Flagged %s as a global factor in channel %s (is applied to all current and future samples in the channel)",
3124 out->GetName(), p3->GetName());
3125 return xRooNode(out, *this);
3126 }
3127
3128 } else if (auto p4 = get<RooAbsPdf>(); p4 && !(fParent && fParent->get<RooRealSumPdf>())) {
3129 // multiply the coefs (if this isn't part of a RooAddPdf or RooRealSumPdf then we will eventually throw exception
3130 return coefs().Multiply(child);
3131 } else if (auto p5 = get<RooAbsReal>(); p5 && (!get<RooAbsPdf>() || (fParent && fParent->get<RooRealSumPdf>()))) {
3132 // replace this obj with a RooProduct to allow for multiplication
3133
3134 // get the list of clients BEFORE creating the new interpolation ... seems list of clients is inaccurate after
3135 std::set<RooAbsArg *> cl;
3136 for (auto &arg : p5->clients()) {
3137 cl.insert(arg);
3138 }
3139
3140 // if multiple clients, see if only one client is in parentage route
3141 // if so, then assume thats the only client we should replace in
3142 if (cl.size() > 1) {
3143 if (cl.count(fParent->get<RooAbsArg>()) > 0) {
3144 cl.clear();
3145 cl.insert(fParent->get<RooAbsArg>());
3146 } else {
3147 Warning("Multiply", "Scaling %s that has multiple clients", p5->GetName());
3148 }
3149 }
3150
3151 auto new_p = acquireNew<RooProduct>(TString::Format("prod_%s", p5->GetName()), p5->GetTitle(), RooArgList(*p5));
3152 // copy attributes over
3153 for (auto &a : p5->attributes())
3154 new_p->setAttribute(a.c_str());
3155 for (auto &a : p5->stringAttributes())
3156 new_p->setStringAttribute(a.first.c_str(), a.second.c_str());
3157 if (!new_p->getStringAttribute("alias"))
3158 new_p->setStringAttribute("alias", p5->GetName());
3159 auto old_p = p5;
3160 new_p->setAttribute(Form("ORIGNAME:%s", old_p->GetName())); // used in redirectServers to say what this replaces
3161 for (auto arg : cl) {
3162 arg->redirectServers(RooArgSet(*new_p), false, true);
3163 }
3164
3165 fComp = new_p;
3166 return Multiply(child);
3167 }
3168
3169 // before giving up here, assume user wanted a norm factor type if child is just a name
3170 if (!child.get() && strlen(opt) == 0)
3171 return Multiply(child, "norm");
3172
3173 throw std::runtime_error(
3174 TString::Format("Cannot multiply %s by %s%s", GetPath().c_str(), child.GetName(),
3175 (!child.get() && strlen(opt) == 0) ? " (forgot to specify factor type?)" : ""));
3176}
3177
3179{
3180
3181 auto p5 = get<RooAbsArg>();
3182 if (!p5) {
3183 throw std::runtime_error("Only replacement of RooAbsArg is supported");
3184 }
3185 node.convertForAcquisition(*this, "func");
3186
3187 auto new_p = node.get<RooAbsArg>();
3188 if (!new_p) {
3189 throw std::runtime_error(TString::Format("Cannot replace with %s", node.GetName()));
3190 }
3191 auto out = acquire(node.fComp);
3192 new_p = std::dynamic_pointer_cast<RooAbsArg>(out).get();
3193
3194 std::set<RooAbsArg *> cl;
3195 for (auto &arg : p5->clients()) {
3196 if (arg == new_p)
3197 continue; // do not replace in self ... although redirectServers will prevent that anyway
3198 cl.insert(arg);
3199 }
3200
3201 // if multiple clients, see if only one client is in parentage route
3202 // if so, then assume thats the only client we should replace in
3203 if (cl.size() > 1) {
3204 if (fParent && fParent->get<RooAbsArg>() && cl.count(fParent->get<RooAbsArg>()) > 0) {
3205 cl.clear();
3206 cl.insert(fParent->get<RooAbsArg>());
3207 } else {
3208 std::stringstream clientList;
3209 for (auto c : cl)
3210 clientList << c->GetName() << ",";
3211 Warning("Replace", "Replacing %s in all clients: %s", p5->GetName(), clientList.str().c_str());
3212 }
3213 }
3214
3215 new_p->setAttribute(Form("ORIGNAME:%s", p5->GetName())); // used in redirectServers to say what this replaces
3216 for (auto arg : cl) {
3217 // if RooFormulaVar need to ensure the internal formula has been "constructed" otherwise will try to construct
3218 // it from the original expression that may have old parameter in it.
3219 if (auto p = dynamic_cast<RooFormulaVar *>(arg))
3220 p->ok(); // triggers creation of RooFormula
3221 arg->redirectServers(RooArgSet(*new_p), false, true);
3222 }
3223 return node;
3224}
3225
3227{
3228
3229 class AutoUpdater {
3230 public:
3231 AutoUpdater(xRooNode &_n) : n(_n) {}
3232 ~AutoUpdater() { n.browse(); }
3233 xRooNode &n;
3234 };
3235 AutoUpdater xxx(*this);
3236
3237 if (!get() && fParent) {
3238 // try to 'create' object based on parentage
3239 // add child as a temporary child to help with decision making
3240 auto _ref = emplace_back(std::shared_ptr<xRooNode>(&const_cast<xRooNode &>(child), [](TObject *) {}));
3241 try {
3242 fComp = fParent->Add(*this, "+").fComp;
3243 } catch (...) {
3244 resize(size() - 1);
3245 std::rethrow_exception(std::current_exception());
3246 }
3247 resize(size() - 1); // remove the temporarily added node
3248 }
3249
3250 if (auto p = mainChild(); p) {
3251 // variations applied to the main child if has one
3252 return p.Vary(child);
3253 }
3254
3255 if (auto s = get<RooSimultaneous>(); s && s->indexCat().IsA() == RooCategory::Class()) {
3256 // name is used as cat label
3257 std::string label = child.GetName();
3258 if (auto pos = label.find('='); pos != std::string::npos)
3259 label = label.substr(pos + 1);
3260 if (!s->indexCat().hasLabel(label)) {
3261 static_cast<RooCategory &>(const_cast<RooAbsCategoryLValue &>(s->indexCat())).defineType(label.c_str());
3262 }
3263 std::shared_ptr<TObject> out;
3264 child.convertForAcquisition(*this);
3265 if (child.get<RooAbsPdf>()) {
3266 out = acquire(child.fComp); // may create a channel from a histogram
3267 } else if (!child.fComp) {
3268 out = acquireNew<RooProdPdf>(TString::Format("%s_%s", s->GetName(), label.c_str()),
3269 (strlen(child.GetTitle())) ? child.GetTitle() : label.c_str(), RooArgList());
3270 Info("Vary", "Created channel RooProdPdf::%s in model %s", out->GetName(), s->GetName());
3271 }
3272
3273 if (auto _pdf = std::dynamic_pointer_cast<RooAbsPdf>(out); _pdf) {
3274 s->addPdf(*_pdf, label.c_str());
3275 sterilize();
3276 // clear children for reload and update shared axis
3277 clear();
3278 fXAxis.reset();
3279 browse();
3280 return xRooNode(TString::Format("%s=%s", s->indexCat().GetName(), label.data()), _pdf, *this);
3281 }
3282
3283 } else if (auto p = get<RooStats::HistFactory::FlexibleInterpVar>(); p) {
3284
3285 // child needs to be a constvar ...
3286 child.convertForAcquisition(*this);
3287 auto _c = child.get<RooConstVar>();
3288 if (!_c && child.get()) {
3289 throw std::runtime_error("Only pure consts can be set as variations of a flexible interpvar");
3290 }
3291#if ROOT_VERSION_CODE < ROOT_VERSION(6, 27, 00)
3292 double value = (_c ? _c->getVal() : p->_nominal);
3293 double nomVal = p->_nominal;
3294#else
3295 double value = (_c ? _c->getVal() : p->nominal());
3296 double nomVal = p->nominal();
3297#endif
3298
3299 TString cName(child.GetName());
3300 if (cName == "nominal") {
3301 p->setNominal(value);
3302 return *(this->variations().at(cName.Data()));
3303 }
3304 if (cName.CountChar('=') != 1) {
3305 throw std::runtime_error("unsupported variation form");
3306 }
3307 std::string parName = cName(0, cName.Index('='));
3308 double parVal = TString(cName(cName.Index('=') + 1, cName.Length())).Atof();
3309 if (parVal != 1 && parVal != -1) {
3310 throw std::runtime_error("unsupported variation magnitude");
3311 }
3312 bool high = parVal > 0;
3313
3314 if (parName.empty()) {
3315 p->setNominal(value);
3316 } else {
3317 auto v = fParent->getObject<RooRealVar>(parName);
3318 if (!v)
3319 v = fParent->acquire<RooRealVar>(parName.c_str(), parName.c_str(), -5, 5);
3320 if (!v->hasError())
3321 v->setError(1);
3322
3323 if (!p->findServer(*v)) {
3324#if ROOT_VERSION_CODE < ROOT_VERSION(6, 27, 00)
3325 p->_paramList.add(*v);
3326 p->_low.push_back(0);
3327 p->_high.push_back(0);
3328 p->_interpCode.push_back(4);
3329#else
3330 const_cast<RooListProxy &>(p->variables()).add(*v);
3331 const_cast<std::vector<double> &>(p->low()).push_back(0);
3332 const_cast<std::vector<double> &>(p->high()).push_back(0);
3333 const_cast<std::vector<int> &>(p->interpolationCodes()).push_back(4);
3334#endif
3335 v->setAttribute(Form("SYMMETRIC%s_%s", high ? "+" : "-", GetName())); // flag for symmetrized
3336 }
3337
3338 if (high) {
3339 p->setHigh(*v, value);
3340 if (v->getAttribute(Form("SYMMETRIC+_%s", GetName()))) {
3341 p->setLow(*v, 2 * nomVal - value);
3342 }
3343 v->setAttribute(Form("SYMMETRIC-_%s", GetName()), false);
3344 } else {
3345 p->setLow(*v, value);
3346 if (v->getAttribute(Form("SYMMETRIC-_%s", GetName()))) {
3347 p->setHigh(*v, 2 * nomVal - value);
3348 }
3349 v->setAttribute(Form("SYMMETRIC+_%s", GetName()), false);
3350 }
3351
3352 /*if (!unconstrained && fParent->pars()[v->GetName()].constraints().empty()) {
3353 fParent->pars()[v->GetName()].constraints().add("normal");
3354 }*/
3355 }
3356 return *(this->variations().at(cName.Data()));
3357 } else if (auto p2 = get<PiecewiseInterpolation>(); p2) {
3358 TString cName(child.GetName());
3359 if (cName.CountChar('=') != 1) {
3360 throw std::runtime_error("unsupported variation form");
3361 }
3362 TString parName = cName(0, cName.Index('='));
3363 double parVal = TString(cName(cName.Index('=') + 1, cName.Length())).Atof();
3364 if (parVal != 1 && parVal != -1) {
3365 throw std::runtime_error("unsupported variation magnitude");
3366 }
3367#if ROOT_VERSION_CODE < ROOT_VERSION(6, 27, 00)
3368 RooHistFunc *f = dynamic_cast<RooHistFunc *>(p2->_nominal.absArg());
3369 if (!f) {
3370 throw std::runtime_error(
3371 TString::Format("Interpolating %s instead of RooHistFunc", p2->_nominal.absArg()->ClassName()));
3372 }
3373#else
3374 RooHistFunc *f = dynamic_cast<RooHistFunc *>(const_cast<RooAbsReal *>(p2->nominalHist()));
3375 if (!f) {
3376 throw std::runtime_error(
3377 TString::Format("Interpolating %s instead of RooHistFunc", p2->nominalHist()->ClassName()));
3378 }
3379#endif
3380 RooHistFunc *nomf = f;
3381 RooHistFunc *otherf = nullptr;
3382 size_t i = 0;
3383 for (auto par : p2->paramList()) {
3384 if (parName == par->GetName()) {
3385 f = dynamic_cast<RooHistFunc *>((parVal > 0 ? p2->highList() : p2->lowList()).at(i));
3386 otherf = dynamic_cast<RooHistFunc *>((parVal > 0 ? p2->lowList() : p2->highList()).at(i));
3387 break;
3388 }
3389 i++;
3390 }
3391 if (i == p2->paramList().size() && !child.get<RooAbsReal>()) {
3392
3393 // need to add the parameter
3394 auto v = acquire<RooRealVar>(parName, parName, -5, 5);
3395 if (!v->hasError())
3396 v->setError(1);
3397
3398 std::shared_ptr<RooHistFunc> up(
3399 static_cast<RooHistFunc *>(f->Clone(Form("%s_%s_up", f->GetName(), parName.Data()))));
3400 std::shared_ptr<RooHistFunc> down(
3401 static_cast<RooHistFunc *>(f->Clone(Form("%s_%s_down", f->GetName(), parName.Data()))));
3402 // RooHistFunc doesn't clone it's data hist ... do it ourself (will be cloned again if imported into a ws)
3403#if ROOT_VERSION_CODE < ROOT_VERSION(6, 27, 00)
3404 std::unique_ptr<RooDataHist> h1(
3405 static_cast<RooDataHist *>(f->dataHist().Clone(Form("hist_%s", up->GetName()))));
3406 std::unique_ptr<RooDataHist> h2(
3407 static_cast<RooDataHist *>(f->dataHist().Clone(Form("hist_%s", down->GetName()))));
3408 up->_dataHist = dynamic_cast<RooDataHist *>(f->dataHist().Clone(Form("hist_%s", up->GetName())));
3409 down->_dataHist = dynamic_cast<RooDataHist *>(f->dataHist().Clone(Form("hist_%s", down->GetName())));
3410#else
3411 up->cloneAndOwnDataHist(TString::Format("hist_%s", up->GetName()));
3412 down->cloneAndOwnDataHist(TString::Format("hist_%s", down->GetName()));
3413#endif
3414 auto ups = std::dynamic_pointer_cast<RooHistFunc>(acquire(up, false, true));
3415 auto downs = std::dynamic_pointer_cast<RooHistFunc>(acquire(down, false, true));
3416#if ROOT_VERSION_CODE < ROOT_VERSION(6, 27, 00)
3417 p2->_highSet.add(*ups.get());
3418 p2->_lowSet.add(*downs.get());
3419 p2->_interpCode.push_back(4);
3420 p2->_paramSet.add(*v);
3421#else
3422 const_cast<RooArgList &>(p2->highList()).add(*ups);
3423 const_cast<RooArgList &>(p2->lowList()).add(*downs);
3424 const_cast<std::vector<int> &>(p2->interpolationCodes()).push_back(4);
3425 const_cast<RooArgList &>(p2->paramList()).add(*v);
3426#endif
3427 p2->setValueDirty();
3428 f = ((parVal > 0) ? ups : downs).get();
3429 otherf = ((parVal > 0) ? downs : ups).get();
3430 // start off with everything being symmetric
3431 f->setStringAttribute("symmetrizes", otherf->GetName());
3432 f->setStringAttribute("symmetrize_nominal", nomf->GetName());
3433 otherf->setStringAttribute("symmetrized_by", f->GetName());
3434
3435 // constrain par if required
3436 /*if (!unconstrained && fParent->pars()[v->GetName()].constraints().empty()) {
3437 fParent->pars()[v->GetName()].constraints().add("normal");
3438 }*/
3439 }
3440
3441 // child.convertForAcquisition(*this);
3442 if (f) {
3443 if (child.get())
3444 xRooNode("tmp", *f, *this) = *child.get();
3445 f->setValueDirty();
3446 xRooNode out(*f, *this);
3447 out.sterilize();
3448 return out;
3449 }
3450
3451 } else if (auto p3 = get<RooConstVar>(); p3) {
3452
3453 // never vary the universal consts ... its too dangerous
3454 if (p3->getAttribute("RooRealConstant_Factory_Object")) {
3455 throw std::runtime_error("Cannot vary pure constants");
3456 }
3457
3458 // inject a FlexibleInterpVar ...
3459
3460 // get the list of clients BEFORE creating the new interpolation ... seems list of clients is inaccurate after
3461 std::set<RooAbsArg *> cl;
3462 for (auto &arg : p3->clients()) {
3463 cl.insert(arg);
3464 }
3465 // if multiple clients, see if only one client is in parentage route
3466 // if so, then assume thats the only client we should replace in
3467 if (cl.size() > 1) {
3468 if (cl.count(fParent->get<RooAbsArg>()) > 0) {
3469 cl.clear();
3470 cl.insert(fParent->get<RooAbsArg>());
3471 } else {
3472 Warning("Vary", "Varying %s that has multiple clients", p3->GetName());
3473 }
3474 }
3475 p3->setStringAttribute("origName", p3->GetName());
3476 TString n = p3->GetName();
3477 p3->SetName(Form("%s_nominal", p3->GetName())); // if problems should perhaps not rename here
3478
3479 auto new_p = acquireNew<RooStats::HistFactory::FlexibleInterpVar>(n, p3->GetTitle(), RooArgList(), p3->getVal(),
3480 std::vector<double>(), std::vector<double>());
3481
3482 // copy attributes over
3483 for (auto &a : p3->attributes())
3484 new_p->setAttribute(a.c_str());
3485 for (auto &a : p3->stringAttributes())
3486 new_p->setStringAttribute(a.first.c_str(), a.second.c_str());
3487 // if (!new_p->getStringAttribute("alias")) new_p->setStringAttribute("alias",p->GetName());
3488 auto old_p = p3;
3489 new_p->setAttribute(Form("ORIGNAME:%s", old_p->GetName())); // used in redirectServers to say what this replaces
3490 for (auto arg : cl) {
3491 arg->redirectServers(RooArgSet(*new_p), false, true);
3492 }
3493
3494 fComp = new_p;
3495 return Vary(child);
3496
3497 } else if (auto p4 = get<RooAbsReal>(); p4) {
3498 // inject an interpolation node
3499
3500 // get the list of clients BEFORE creating the new interpolation ... seems list of clients is inaccurate after
3501 std::set<RooAbsArg *> cl;
3502 for (auto &arg : p4->clients()) {
3503 cl.insert(arg);
3504 }
3505 // if multiple clients, see if only one client is in parentage route
3506 // if so, then assume thats the only client we should replace in
3507 if (cl.size() > 1) {
3508 if (cl.count(fParent->get<RooAbsArg>()) > 0) {
3509 cl.clear();
3510 cl.insert(fParent->get<RooAbsArg>());
3511 } else {
3512 Warning("Vary", "Varying %s that has multiple clients", p4->GetName());
3513 }
3514 }
3515 p4->setStringAttribute("origName", p4->GetName());
3516 TString n = p4->GetName();
3517 p4->SetName(Form("%s_nominal", p4->GetName())); // if problems should perhaps not rename here
3518
3519 auto new_p = acquireNew<PiecewiseInterpolation>(n, p4->GetTitle(), *p4, RooArgList(), RooArgList(), RooArgList());
3520
3521 // copy attributes over
3522 for (auto &a : p4->attributes())
3523 new_p->setAttribute(a.c_str());
3524 for (auto &a : p4->stringAttributes())
3525 new_p->setStringAttribute(a.first.c_str(), a.second.c_str());
3526 // if (!new_p->getStringAttribute("alias")) new_p->setStringAttribute("alias",p->GetName());
3527 auto old_p = p4;
3528 new_p->setAttribute(Form("ORIGNAME:%s", old_p->GetName())); // used in redirectServers to say what this replaces
3529 for (auto arg : cl) {
3530 arg->redirectServers(RooArgSet(*new_p), false, true);
3531 }
3532
3533 fComp = new_p;
3534 return Vary(child);
3535 }
3536
3537 Print();
3538 throw std::runtime_error(TString::Format("Cannot vary %s with %s", GetName(), child.GetName()));
3539}
3540
3542{
3544}
3545
3546bool xRooNode::SetContent(double value, const char *par, double val)
3547{
3548 return SetContents(RooConstVar(GetName(), GetTitle(), value), par, val);
3549}
3550
3553 {
3554 if (x && b)
3555 x->setBinning(*b);
3556 if (b)
3557 delete b;
3558 }
3559 RooRealVar *x = nullptr;
3560 RooAbsBinning *b = nullptr;
3561};
3562
3564{
3565
3566 if (!get()) {
3567 fComp = std::shared_ptr<TObject>(const_cast<TObject *>(&o), [](TObject *) {});
3568 if (fParent && !fParent->find(GetName())) {
3569 // either a temporary or a placeholder so need to try genuinely adding
3570 fComp = fParent->Add(*this, "+").fComp;
3571 if (auto a = get<RooAbsArg>(); a && strcmp(a->GetName(), GetName()) && !a->getStringAttribute("alias")) {
3572 a->setStringAttribute("alias", GetName());
3573 }
3574 if (!fComp)
3575 throw std::runtime_error("Cannot determine type");
3576 return *this;
3577 }
3578 }
3579
3580 if (auto h = dynamic_cast<const TH1 *>(&o); h) {
3581 /*auto f = get<RooHistFunc>();
3582 if (!f) {
3583 // if it's a RooProduct locate child with the same name
3584 if (get<RooProduct>()) {
3585 f = factors()[GetName()]->get<RooHistFunc>();
3586 }
3587
3588
3589
3590 }*/
3591 bool _isData = get<RooAbsData>();
3592 BinningRestorer _b;
3593 if (_isData) {
3594 // need to ensure x-axis matches this h
3595 auto ax = GetXaxis();
3596 if (!ax)
3597 throw std::runtime_error("no xaxis");
3598 auto _v = dynamic_cast<RooRealVar *>(ax->GetParent());
3599 if (_v) {
3600 _b.x = _v;
3601 _b.b = dynamic_cast<RooAbsBinning *>(_v->getBinningPtr(nullptr)->Clone());
3602 if (h->GetXaxis()->IsVariableBinSize()) {
3603 _v->setBinning(RooBinning(h->GetNbinsX(), h->GetXaxis()->GetXbins()->GetArray()));
3604 } else {
3605 _v->setBinning(RooUniformBinning(h->GetXaxis()->GetXmin(), h->GetXaxis()->GetXmax(), h->GetNbinsX()));
3606 }
3607 }
3608 }
3609
3610 if (true) {
3611 for (int bin = 1; bin <= h->GetNbinsX(); bin++) {
3612 SetBinContent(bin, h->GetBinContent(bin));
3613 /*double value = h->GetBinContent(bin);
3614 auto bin_pars = f->dataHist().get(bin - 1);
3615 if (f->getAttribute("density")) {
3616 value /= f->dataHist().binVolume(*bin_pars);
3617 }
3618 f->dataHist().set(*bin_pars, value);*/
3619 if (!_isData && h->GetSumw2N() && !SetBinError(bin, h->GetBinError(bin)))
3620 throw std::runtime_error("Failed setting stat error");
3621 }
3622 return *this;
3623 }
3624 } else if (auto _c = dynamic_cast<const RooConstVar *>(&o); _c) {
3625
3626 if (auto a = get<RooAbsArg>();
3627 (a && a->isFundamental()) || get<RooConstVar>() || get<RooStats::HistFactory::FlexibleInterpVar>()) {
3628 SetBinContent(1, _c->getVal());
3629 return *this;
3630 } else if (get<RooAbsData>()) { // try to do assignment to a dataset (usually setting a bin content)
3631 SetBinContent(0, _c->getVal());
3632 return *this;
3633 }
3634 }
3635
3636 throw std::runtime_error("Assignment failed");
3637
3638 /*
3639
3640 if (fParent && !fParent->mk()) {
3641 throw std::runtime_error("mk failure");
3642 }
3643
3644 if (fComp) return *this;
3645
3646 if (o.InheritsFrom("RooAbsArg")) {
3647 fComp = acquire(std::shared_ptr<TObject>(const_cast<TObject*>(&o),[](TObject* o){}));
3648 std::dynamic_pointer_cast<RooAbsArg>(fComp)->setStringAttribute("alias",GetName());
3649 }
3650
3651 if (fComp && fParent) {
3652 fParent->incorporate(fComp);
3653 }
3654
3655
3656 return *this;
3657 */
3658}
3659
3660void xRooNode::_fit_(const char *constParValues)
3661{
3662 try {
3663 auto _pars = pars();
3664 // std::unique_ptr<RooAbsCollection> snap(_pars.argList().snapshot());
3665 TStringToken pattern(constParValues, ",");
3666 std::map<RooAbsRealLValue *, double> valsToSet;
3667 while (pattern.NextToken()) {
3668 auto idx = pattern.Index('=');
3669 TString pat = (idx == -1) ? TString(pattern) : TString(pattern(0, idx));
3670 double val =
3671 (idx == -1) ? std::numeric_limits<double>::quiet_NaN() : TString(pattern(idx + 1, pattern.Length())).Atof();
3672 for (auto p : _pars.argList()) {
3673 if (TString(p->GetName()).Contains(TRegexp(pat, true))) {
3674 p->setAttribute("Constant", true);
3675 if (!std::isnan(val)) {
3676 valsToSet[dynamic_cast<RooAbsRealLValue *>(p)] = val;
3677 // dynamic_cast<RooAbsRealLValue *>(p)->setVal(val); // don't set yet, to allow for asimov dataset
3678 // creation based on current values
3679 }
3680 }
3681 }
3682 }
3683 // use the first selected dataset
3684 auto _dsets = datasets();
3685 TString dsetName = "";
3686 for (auto &d : _dsets) {
3687 if (d->get()->TestBit(1 << 20)) {
3688 dsetName = d->get()->GetName();
3689 break;
3690 }
3691 }
3692 auto _nll = nll(dsetName.Data());
3693 // can now set the values
3694 for (auto [p, v] : valsToSet) {
3695 p->setVal(v);
3696 }
3697 _nll.fitConfigOptions()->SetValue("LogSize", 65536);
3698 _nll.fitConfig()->MinimizerOptions().SetPrintLevel(0);
3699 auto fr = _nll.minimize();
3700 //_pars.argList() = *snap; // restore values - irrelevant as SetFitResult will restore values
3701 if (!fr.get())
3702 throw std::runtime_error("Fit Failed");
3703 SetFitResult(fr.get());
3704 TString statusCodes;
3705 for (unsigned int i = 0; i < fr->numStatusHistory(); i++) {
3706 statusCodes += TString::Format("\n%s = %d", fr->statusLabelHistory(i), fr->statusCodeHistory(i));
3707 }
3708 const TGWindow *w =
3709 (gROOT->GetListOfBrowsers()->At(0))
3710 ? dynamic_cast<TGWindow *>(static_cast<TBrowser *>(gROOT->GetListOfBrowsers()->At(0))->GetBrowserImp())
3711 : gClient->GetRoot();
3712 if (fr->status() != 0) {
3713 new TGMsgBox(gClient->GetRoot(), w, "Fit Finished with Bad Status Code",
3714 TString::Format("%s\nData = %s\nFit Status Code = %d\nCov Quality = %d\n-------------%s",
3715 fr->GetName(), dsetName.Data(), fr->status(), fr->covQual(), statusCodes.Data()),
3717 } else if (fr->covQual() != 3 && _nll.fitConfig()->ParabErrors()) {
3718 new TGMsgBox(gClient->GetRoot(), w, "Fit Finished with Bad Covariance Quality",
3719 TString::Format("%s\nData = %s\nFit Status Code = %d\nCov Quality = %d\n-------------%s",
3720 fr->GetName(), dsetName.Data(), fr->status(), fr->covQual(), statusCodes.Data()),
3722 } else {
3723 new TGMsgBox(gClient->GetRoot(), w, "Fit Finished Successfully",
3724 TString::Format("%s\nData = %s\nFit Status Code = %d\nCov Quality = %d\n-------------%s",
3725 fr->GetName(), dsetName.Data(), fr->status(), fr->covQual(), statusCodes.Data()));
3726 }
3727 } catch (const std::exception &e) {
3728 new TGMsgBox(
3729 gClient->GetRoot(),
3730 (gROOT->GetListOfBrowsers()->At(0))
3731 ? dynamic_cast<TGWindow *>(static_cast<TBrowser *>(gROOT->GetListOfBrowsers()->At(0))->GetBrowserImp())
3732 : gClient->GetRoot(),
3733 "Exception", e.what(), kMBIconExclamation, kMBOk); // deletes self on dismiss?
3734 }
3735}
3736
3737void xRooNode::_generate_(const char *datasetName, bool expected)
3738{
3739 try {
3740 datasets().Add(datasetName, expected ? "asimov" : "toy");
3741 } catch (const std::exception &e) {
3742 new TGMsgBox(
3743 gClient->GetRoot(),
3744 (gROOT->GetListOfBrowsers()->At(0))
3745 ? dynamic_cast<TGWindow *>(static_cast<TBrowser *>(gROOT->GetListOfBrowsers()->At(0))->GetBrowserImp())
3746 : gClient->GetRoot(),
3747 "Exception", e.what(),
3748 kMBIconExclamation); // deletes self on dismiss?
3749 }
3750}
3751
3752void xRooNode::_scan_(const char *what, double nToys, const char *xvar, int nBinsX, double lowX,
3753 double highX /*, const char*, int, double, double*/, const char *constParValues)
3754{
3755 try {
3756 TString sXvar(xvar);
3757 TString sWhat(what);
3758
3759 // use the first selected dataset
3760 auto _dsets = datasets();
3761 TString dsetName = "";
3762 for (auto &d : _dsets) {
3763 if (d->get()->TestBit(1 << 20)) {
3764 dsetName = d->get()->GetName();
3765 break;
3766 }
3767 }
3768 auto _pars = pars();
3769 std::unique_ptr<RooAbsCollection> snap(_pars.argList().snapshot());
3770 TStringToken pattern(constParValues, ",");
3771 while (pattern.NextToken()) {
3772 auto idx = pattern.Index('=');
3773 TString pat = (idx == -1) ? TString(pattern) : TString(pattern(0, idx));
3774 double val =
3775 (idx == -1) ? std::numeric_limits<double>::quiet_NaN() : TString(pattern(idx + 1, pattern.Length())).Atof();
3776 for (auto par : _pars.argList()) {
3777 if (TString(par->GetName()).Contains(TRegexp(pat, true))) {
3778 par->setAttribute("Constant", true);
3779 if (!std::isnan(val)) {
3780 dynamic_cast<RooAbsRealLValue *>(par)->setVal(val);
3781 }
3782 }
3783 }
3784 }
3785 auto hs = nll(dsetName.Data()).hypoSpace(sXvar);
3786 if (nToys) {
3787 sWhat += " toys";
3788 if (nToys > 0) {
3789 sWhat += TString::Format("=%g", nToys);
3790 }
3791 }
3792 hs.SetTitle(sWhat + " scan" + ((dsetName != "") ? TString::Format(" [data=%s]", dsetName.Data()) : ""));
3793 int scanStatus = hs.scan(sWhat + " visualize", nBinsX, lowX, highX);
3794 if (scanStatus != 0) {
3795 new TGMsgBox(
3796 gClient->GetRoot(),
3797 (gROOT->GetListOfBrowsers()->At(0))
3798 ? dynamic_cast<TGWindow *>(static_cast<TBrowser *>(gROOT->GetListOfBrowsers()->At(0))->GetBrowserImp())
3799 : gClient->GetRoot(),
3800 "Scan Finished with Bad Status Code",
3801 TString::Format("%s\nData = %s\nScan Status Code = %d", hs.GetName(), dsetName.Data(), scanStatus),
3803 }
3804 hs.SetName(TUUID().AsString());
3805 if (ws()) {
3806 if (auto res = hs.result())
3807 ws()->import(*res);
3808 }
3809
3810 _pars.argList() = *snap; // restore pars
3811
3812 } catch (const std::exception &e) {
3813 new TGMsgBox(
3814 gClient->GetRoot(),
3815 (gROOT->GetListOfBrowsers()->At(0))
3816 ? dynamic_cast<TGWindow *>(static_cast<TBrowser *>(gROOT->GetListOfBrowsers()->At(0))->GetBrowserImp())
3817 : gClient->GetRoot(),
3818 "Exception", e.what(), kMBIconExclamation);
3819 }
3820}
3821
3822void xRooNode::_SetBinContent_(int bin, double value, const char *par, double parVal)
3823{
3824 try {
3825 SetBinContent(bin, value, strlen(par) > 0 ? par : nullptr, parVal);
3826 } catch (const std::exception &e) {
3827 new TGMsgBox(gClient->GetRoot(), gClient->GetRoot(), "Exception", e.what(),
3828 kMBIconExclamation); // deletes self on dismiss?
3829 }
3830}
3831
3833{
3834 try {
3835#if ROOT_VERSION_CODE > ROOT_VERSION(6, 29, 00)
3836 // if this is a collection of values, populate a TF1 and display as a dialog
3837 if (!get() && TString(GetName()).BeginsWith("!")) {
3838 browse();
3839 RooArgList args;
3840 for (auto a : *this) {
3841 if (auto arg = a->get<RooRealVar>())
3842 args.add(*arg);
3843 }
3844 TF1 f(GetName(), 0.0, 1.0, std::min(int(args.size()), 10));
3845 int i = 0;
3846 int j = 0;
3847 for (auto c : args) {
3848 j++;
3849 if (j < value) {
3850 continue;
3851 }
3852 auto v = dynamic_cast<RooRealVar *>(c);
3853 f.SetParName(i, c->GetName());
3854 if (v) {
3855 f.SetParLimits(i, v->getMin(), v->getMax());
3856 if (v->isConstant())
3857 f.FixParameter(i, v->getVal());
3858 else {
3859 f.SetParameter(i, v->getVal());
3860 f.SetParError(i, v->getError());
3861 }
3862 }
3863 i++;
3864 if (i == 10) {
3865 break; // max 10 pars shown
3866 }
3867 }
3868 int ret = 0;
3870 gClient->GetDefaultRoot(),
3871 (gROOT->GetListOfBrowsers()->At(0))
3872 ? dynamic_cast<TGWindow *>(static_cast<TBrowser *>(gROOT->GetListOfBrowsers()->At(0))->GetBrowserImp())
3873 : gClient->GetDefaultRoot(),
3874 &f, nullptr, &ret);
3875 if (ret) {
3876 // user has changed parameter values etc, propagate back to parameters
3877 for (i = 0; i < f.GetNpar(); i++) {
3878 auto c = args.find(f.GetParName(i));
3879 auto v = dynamic_cast<RooRealVar *>(c);
3880 if (v) {
3881 v->setVal(f.GetParameter(i));
3882 double low, high;
3883 f.GetParLimits(i, low, high);
3884 if (low == high) {
3885 v->setConstant(low); // if low==high==0 then is not marked constant
3886 } else {
3887 v->setRange(low, high);
3888 }
3889 }
3890 }
3891 }
3892 return;
3893 }
3894#endif
3895
3896 if (!SetContent(value))
3897 throw std::runtime_error("Failed to SetContent");
3898 } catch (const std::exception &e) {
3899 new TGMsgBox(gClient->GetRoot(), gClient->GetRoot(), "Exception", e.what(),
3900 kMBIconExclamation); // deletes self on dismiss?
3901 }
3902}
3903
3904bool xRooNode::SetBinContent(int bin, double value, const char *par, double parVal)
3905{
3906
3907 // create if needed
3908 if (!get()) {
3909 if (fParent && !find(GetName())) {
3910 // if have a binning we create a histogram to match it
3911 if (auto ax = GetXaxis(); ax) {
3912 std::shared_ptr<TH1D> h;
3913 auto _b = dynamic_cast<Axis2 *>(ax)->binning();
3914 auto t = TH1::AddDirectoryStatus();
3915 TH1::AddDirectory(false);
3916 if (_b->isUniform()) {
3917 h.reset(new TH1D(GetName(), GetTitle(), _b->numBins(), _b->lowBound(), _b->highBound()));
3918 } else {
3919 h.reset(new TH1D(GetName(), GetTitle(), _b->numBins(), _b->array()));
3920 }
3921 h->SetOption("nostyle"); // don't transfer style when added
3922 h->SetDirectory(nullptr);
3924 h->GetXaxis()->SetName(TString::Format("%s;%s", ax->GetParent()->GetName(), ax->GetName()));
3925 fComp = h;
3926 }
3927 fComp = fParent->Add(*this, "sample").fComp;
3928 }
3929 }
3930
3931 // if it's a RooProduct locate child with the same name
3932 if (get<RooProduct>()) {
3933 return factors()[GetName()]->SetBinContent(bin, value, par, parVal);
3934 }
3935
3936 if (get<RooAbsData>()) {
3937 if (auto _data = get<RooDataSet>(); _data) {
3938 auto _ax = (bin) ? GetXaxis() : nullptr;
3939 if (!_ax && bin) {
3940 throw std::runtime_error("Cannot determine binning to fill data");
3941 }
3942 if (_ax && _ax->GetNbins() < bin) {
3943 throw std::out_of_range(TString::Format("%s range %s only has %d bins", _ax->GetParent()->GetName(),
3944 _ax->GetName(), _ax->GetNbins()));
3945 }
3946 RooArgSet obs;
3947
3948 TString cut = "";
3949
3950 for (auto _c : coords()) { // coords() moves vars to their respective coordinates too
3951 if (auto _cat = _c->get<RooAbsCategoryLValue>(); _cat) {
3952 if (cut != "")
3953 cut += " && ";
3954 cut += TString::Format("%s==%d", _cat->GetName(), _cat->getCurrentIndex());
3955 obs.add(*_cat); // note: if we ever changed coords to return clones, would need to keep coords alive
3956 } else if (auto _rv = _c->get<RooAbsRealLValue>(); _rv) {
3957 // todo: check coordRange is a single range rather than multirange
3958 if (cut != "")
3959 cut += " && ";
3960 cut +=
3961 TString::Format("%s>=%f&&%s<%f", _rv->GetName(), _rv->getMin(_rv->getStringAttribute("coordRange")),
3962 _rv->GetName(), _rv->getMax(_rv->getStringAttribute("coordRange")));
3963 obs.add(*_rv); // note: if we ever changed coords to return clones, would need to keep coords alive
3964 } else {
3965 throw std::runtime_error("SetBinContent of data: Unsupported coordinate type");
3966 }
3967 }
3968
3969 RooFormulaVar cutFormula("cut1", cut, obs); // doing this to avoid complaints about unused vars
3970 RooFormulaVar icutFormula("icut1", TString::Format("!(%s)", cut.Data()), obs);
3971
3972 TString cut2;
3973 if (_ax) {
3974 cut2 = TString::Format("%s >= %f && %s < %f", _ax->GetParent()->GetName(), _ax->GetBinLowEdge(bin),
3975 _ax->GetParent()->GetName(), _ax->GetBinUpEdge(bin));
3976 obs.add(*dynamic_cast<RooAbsArg *>(_ax->GetParent()));
3977 } else {
3978 cut2 = "1==1";
3979 }
3980 RooFormulaVar cutFormula2("cut2", cut + " && " + cut2, obs);
3981 RooFormulaVar icutFormula2("icut2", TString::Format("!(%s && %s)", cut.Data(), cut2.Data()), obs);
3982
3983 // // go up through parents looking for slice obs
3984 // auto _p = fParent;
3985 // while(_p) {
3986 // TString pName(_p->GetName());
3987 // if (auto pos = pName.Index('='); pos != -1) {
3988 // if(auto _obs = _p->getObject<RooAbsLValue>(pName(0,pos)); _obs) {
3989 // if(auto _cat = dynamic_cast<RooAbsCategoryLValue*>(_obs.get()); _cat) {
3990 // _cat->setLabel(pName(pos+1,pName.Length()));
3991 // cut += TString::Format("%s%s==%d", (cut=="")?"":" && ",_cat->GetName(),
3992 // _cat->getCurrentIndex());
3993 // } else if(auto _var = dynamic_cast<RooAbsRealLValue*>(_obs.get()); _var) {
3994 // _var->setVal(TString(pName(pos+1,pName.Length())).Atof());
3995 // // TODO: Cut for this!!
3996 // }
3997 // obs.add(*dynamic_cast<RooAbsArg*>(_obs.get()));
3998 // } else {
3999 // throw std::runtime_error("Unknown observable, could not find");
4000 // }
4001 // }
4002 // _p = _p->fParent;
4003 // }
4004
4005 // add observables to dataset if necessary
4006 RooArgSet l(obs);
4007 l.remove(*_data->get(), true, true);
4008 if (!l.empty()) {
4009 // addColumns method is buggy: https://github.com/root-project/root/issues/8787
4010 // incredibly though, addColumn works??
4011 for (auto &x : l) {
4012 _data->addColumn(*x);
4013 }
4014 // instead create a copy dataset and merge it into current
4015 // cant use merge because it drops weightVar
4016 /*RooDataSet tmp("tmp","tmp",l);
4017 for(int i=0;i<_data->numEntries();i++) tmp.add(l);
4018 _data->merge(&tmp);*/
4019 // delete _data->addColumns(l);
4020 }
4021 // before adding, ensure range is good to cover
4022 for (auto &o : obs) {
4023 if (auto v = dynamic_cast<RooRealVar *>(o); v) {
4024 if (auto dv = dynamic_cast<RooRealVar *>(_data->get()->find(v->GetName())); dv) {
4025 if (v->getMin() < dv->getMin())
4026 dv->setMin(v->getMin());
4027 if (v->getMax() > dv->getMax())
4028 dv->setMax(v->getMax());
4029 }
4030 } else if (auto c = dynamic_cast<RooCategory *>(o); c) {
4031 if (auto dc = dynamic_cast<RooCategory *>(_data->get()->find(c->GetName())); dc) {
4032 if (!dc->hasLabel(c->getCurrentLabel())) {
4033 dc->defineType(c->getCurrentLabel(), c->getCurrentIndex());
4034 }
4035 }
4036 }
4037 }
4038
4039 // using SetBinContent means dataset must take on a binned form at these coordinates
4040 // if number of entries doesnt match number of bins then will 'bin' the data
4041 if (bin) {
4042 if (auto _nentries = std::unique_ptr<RooAbsData>(_data->reduce(cutFormula))->numEntries();
4043 _nentries != _ax->GetNbins()) {
4044 auto _contents = GetBinContents(1, _ax->GetNbins());
4045
4046 if (_nentries > 0) {
4047 Info("SetBinContent", "Binning %s in channel: %s", GetName(), cut.Data());
4048 auto _reduced = std::unique_ptr<RooAbsData>(_data->reduce(icutFormula));
4049 _data->reset();
4050 for (int j = 0; j < _reduced->numEntries(); j++) {
4051 auto _obs = _reduced->get(j);
4052 _data->add(*_obs, _reduced->weight());
4053 }
4054 }
4055 for (int i = 1; i <= _ax->GetNbins(); i++) {
4056 // can skip over the bin we will be setting to save a reduce step below
4057 if (i == bin)
4058 continue;
4059 dynamic_cast<RooAbsLValue *>(_ax->GetParent())->setBin(i - 1, _ax->GetName());
4060 _data->add(obs, _contents.at(i - 1));
4061 }
4062 }
4063 }
4064 // remove existing entries
4065 if (std::unique_ptr<RooAbsData>(_data->reduce(cutFormula2))->numEntries() > 0) {
4066 auto _reduced = std::unique_ptr<RooAbsData>(_data->reduce(icutFormula2));
4067 _data->reset();
4068 for (int j = 0; j < _reduced->numEntries(); j++) {
4069 auto _obs = _reduced->get(j);
4070 _data->add(*_obs, _reduced->weight());
4071 }
4072 }
4073 if (_ax)
4074 dynamic_cast<RooAbsLValue *>(_ax->GetParent())->setBin(bin - 1, _ax->GetName());
4075 _data->add(obs, value);
4076 if (auto bb = getBrowsable(".sourceds"))
4077 return bb->SetBinContent(bin, value, par, parVal); // apply to source ds if we have one
4078 return true;
4079
4080 } else if (get<RooDataHist>()) {
4081 throw std::runtime_error("RooDataHist not supported yet");
4082 }
4083 }
4084
4085 if (auto _varies = variations(); !_varies.empty() || (par && strlen(par))) {
4086 if (!par || strlen(par) == 0) {
4087 return _varies["nominal"]->SetBinContent(bin, value, par, parVal);
4088 } else if (auto it = _varies.find(Form("%s=%g", par, parVal)); it) {
4089 return it->SetBinContent(bin, value);
4090 } else {
4091 // need to create the variation : note - if no variations existed up to now this will inject a new node
4092 // so we should redirect ourself to the new node
4093 // TODO: Do we need to redirect parents?
4094 TString s = Form("%s=%g", par, parVal);
4095 return Vary(s.Data()).SetBinContent(bin, value);
4096 }
4097 }
4098
4099 auto o = get();
4100 if (auto p = dynamic_cast<RooRealVar *>(o); p) {
4101 if (!par || strlen(par) == 0) {
4102 if (p->getMax() < value)
4103 p->setMax(value);
4104 if (p->getMin() > value)
4105 p->setMin(value);
4106 p->setVal(value);
4107 sterilize();
4108 return true;
4109 }
4110
4111 } else if (auto c = dynamic_cast<RooConstVar *>(o); c) {
4112
4113 // if parent is a FlexibleInterpVar, change the value in that .
4114 if (strcmp(c->GetName(), Form("%g", c->getVal())) == 0) {
4115 c->SetNameTitle(Form("%g", value), Form("%g", value));
4116 }
4117#if ROOT_VERSION_CODE < ROOT_VERSION(6, 24, 00)
4118 c->_value = value; // in future ROOT versions there is a changeVal method!
4119#else
4120 c->changeVal(value);
4121#endif
4122
4124 fParent->Vary(*this);
4125 }
4126
4127 sterilize();
4128 return true;
4129 } else if (auto f = dynamic_cast<RooHistFunc *>(o); f) {
4130 auto bin_pars = f->dataHist().get(bin - 1);
4131 if (f->getAttribute("density")) {
4132 value /= f->dataHist().binVolume(*bin_pars);
4133 }
4134 f->dataHist().set(*bin_pars, value);
4135 f->setValueDirty();
4136
4137 if (auto otherfName = f->getStringAttribute("symmetrized_by"); otherfName) {
4138 // broken symmetry, so update flags ...
4139 f->setStringAttribute("symmetrized_by", nullptr);
4140 if (auto x = getObject<RooAbsArg>(otherfName); x) {
4141 x->setStringAttribute("symmetrizes", nullptr);
4142 x->setStringAttribute("symmetrize_nominal", nullptr);
4143 }
4144 } else if (auto otherfName2 = f->getStringAttribute("symmetrizes"); otherfName2) {
4145 auto nomf = getObject<RooHistFunc>(f->getStringAttribute("symmetrize_nominal"));
4146 auto otherf = getObject<RooHistFunc>(otherfName2);
4147 if (nomf && otherf) {
4148 otherf->dataHist().set(*bin_pars, 2 * nomf->dataHist().weight(bin - 1) - value);
4149 otherf->setValueDirty();
4150 }
4151 }
4152 sterilize();
4153 return true;
4154 } else if (auto f2 = dynamic_cast<RooStats::HistFactory::FlexibleInterpVar *>(o); f2) {
4155 // changing nominal value
4156 f2->setNominal(value);
4157 }
4158 throw std::runtime_error(TString::Format("unable to set bin content of %s", GetPath().c_str()));
4159}
4160
4161bool xRooNode::SetBinData(int bin, double value, const xRooNode &data)
4162{
4163 if (data.get<RooAbsData>()) {
4164 // attach as a child before calling datasets(), so that is included in the list
4165 push_back(std::make_shared<xRooNode>(data));
4166 }
4167 auto node = datasets()[data.GetName()];
4168 if (data.get<RooAbsData>()) {
4169 // remove the child we attached
4170 resize(size() - 1);
4171 }
4172 return node->SetBinContent(bin, value);
4173}
4174
4175bool xRooNode::SetData(const TObject &obj, const xRooNode &data)
4176{
4177 if (data.get<RooAbsData>()) {
4178 // attach as a child before calling datasets(), so that is included in the list
4179 push_back(std::make_shared<xRooNode>(data));
4180 }
4181 auto node = datasets()[data.GetName()];
4182 if (data.get<RooAbsData>()) {
4183 // remove the child we attached
4184 resize(size() - 1);
4185 }
4186 return node->SetContents(obj);
4187}
4188
4189bool xRooNode::SetBinError(int bin, double value)
4190{
4191
4192 // if it's a RooProduct locate child with the same name
4193 if (get<RooProduct>()) {
4194 return factors()[GetName()]->SetBinError(bin, value);
4195 }
4196
4197 if (auto _varies = variations(); !_varies.empty()) {
4198 return _varies["nominal"]->SetBinError(bin, value);
4199 }
4200
4201 auto o = get();
4202
4203 if (auto f = dynamic_cast<RooHistFunc *>(o); f) {
4204
4205 // if (f->getAttribute("density")) { value /= f->dataHist().binVolume(*bin_pars); } - commented out because DON'T
4206 // convert .. sumw and sumw2 attributes will be stored not as densities
4207
4208 // NOTE: Can only do this because factors() makes parents of its children it's own parent (it isn't the parent)
4209 // If ever make factors etc part of the parentage then this would need tweaking to get to the true parent
4210 // find first parent that is a RooProduct, that is where the statFactor would live
4211 // stop as soon as we reach pdf object
4212 auto _prodParent = fParent;
4213 while (_prodParent && !_prodParent->get<RooProduct>() && !_prodParent->get<RooAbsPdf>()) {
4214 if (_prodParent->get<PiecewiseInterpolation>() && strcmp(GetName(), "nominal")) {
4215 _prodParent.reset();
4216 break; // only the 'nominal' variation can look for a statFactor outside the variation container
4217 }
4218 _prodParent = _prodParent->fParent;
4219 }
4220 auto _f_stat =
4221 (_prodParent && !_prodParent->get<RooAbsPdf>()) ? _prodParent->factors().find("statFactor") : nullptr;
4222 auto f_stat = (_f_stat) ? _f_stat->get<ParamHistFunc>() : nullptr;
4223 if (_f_stat && _f_stat->get() && !f_stat) {
4224 throw std::runtime_error("stat factor must be a paramhistfunc");
4225 }
4226
4227 // stat uncertainty lives in the "statFactor" factor, each sample has its own one,
4228 // but they can share parameters
4229 if (!f_stat) {
4230 if (value == 0)
4231 return true;
4232 TString parNames;
4233 for (auto &p : xRooNode("tmp", *f, std::shared_ptr<xRooNode>(nullptr)).vars()) {
4234 if (parNames != "")
4235 parNames += ",";
4236 parNames += p->get()->GetName();
4237 }
4238 auto h = std::unique_ptr<TH1>(f->dataHist().createHistogram(parNames
4239#if ROOT_VERSION_CODE >= ROOT_VERSION(6, 27, 00)
4240 ,
4242#endif
4243 ));
4244 h->Reset();
4245 h->SetName("statFactor");
4246 h->SetTitle(TString::Format("StatFactor of %s", f->GetTitle()));
4247 h->SetOption("blankshape");
4248
4249 // multiply parent if is nominal
4250 auto toMultiply = this;
4251 if (strcmp(GetName(), "nominal") == 0 && fParent && fParent->get<PiecewiseInterpolation>())
4252 toMultiply = fParent.get();
4253
4254 f_stat = dynamic_cast<ParamHistFunc *>(toMultiply->Multiply(*h).get());
4255 if (!f_stat) {
4256 throw std::runtime_error("Failed creating stat shapeFactor");
4257 }
4258 }
4259
4260 auto phf = f_stat;
4261
4262 TString prefix = f->getStringAttribute("statPrefix");
4263 if (value && prefix == "") {
4264 // find the first parent that can hold components (RooAddPdf, RooRealSumPdf, RooAddition, RooWorkspace) ... use
4265 // that name for the stat factor
4266 auto _p = fParent;
4267 while (_p && !(_p->get()->InheritsFrom("RooRealSumPdf") || _p->get()->InheritsFrom("RooAddPdf") ||
4268 _p->get()->InheritsFrom("RooWorkspace") || _p->get()->InheritsFrom("RooAddition"))) {
4269 _p = _p->fParent;
4270 }
4271 prefix = TString::Format("stat_%s", (_p && _p->get<RooAbsReal>()) ? _p->get()->GetName() : f->GetName());
4272 }
4273 auto newVar = (value == 0) ? getObject<RooRealVar>("1")
4274 : acquire<RooRealVar>(Form("%s_bin%d", prefix.Data(), bin),
4275 Form("%s_bin%d", prefix.Data(), bin), 1);
4276#if ROOT_VERSION_CODE < ROOT_VERSION(6, 27, 00)
4277 RooArgList &pSet = phf->_paramSet;
4278#else
4279 RooArgList &pSet = const_cast<RooArgList &>(phf->paramList());
4280#endif
4281 auto var = dynamic_cast<RooRealVar *>(&pSet[bin - 1]);
4282
4283 if (newVar.get() != var) {
4284 // need to swap out var for newVar
4285 // replace ith element in list with new func, or inject into RooProduct
4286 RooArgList all;
4287 for (std::size_t i = 0; i < pSet.size(); i++) {
4288 if (int(i) != bin - 1) {
4289 all.add(*pSet.at(i));
4290 } else {
4291 all.add(*newVar);
4292 }
4293 }
4294 pSet.removeAll();
4295 pSet.add(all);
4296 }
4297
4298 xRooNode v((value == 0) ? *var : *newVar, *this);
4299 auto rrv = dynamic_cast<RooRealVar *>(v.get());
4300 if (strcmp(rrv->GetName(), "1") != 0) {
4301 TString origName = (f->getStringAttribute("origName")) ? f->getStringAttribute("origName") : GetName();
4302 rrv->setStringAttribute(Form("sumw2_%s", origName.Data()), TString::Format("%f", pow(value, 2)));
4303 auto bin_pars = f->dataHist().get(bin - 1);
4304 auto _binContent = f->dataHist().weight();
4305 if (f->getAttribute("density")) {
4306 _binContent *= f->dataHist().binVolume(*bin_pars);
4307 }
4308 rrv->setStringAttribute(Form("sumw_%s", origName.Data()), TString::Format("%f", _binContent));
4309 double sumw2 = 0;
4310 double sumw = 0;
4311 for (auto &[s, sv] : rrv->stringAttributes()) {
4312 if (s.find("sumw_") == 0) {
4313 sumw += TString(sv).Atof();
4314 } else if (s.find("sumw2_") == 0) {
4315 sumw2 += TString(sv).Atof();
4316 }
4317 }
4318 if (sumw2 && sumw2 != std::numeric_limits<double>::infinity()) {
4319 double tau = pow(sumw, 2) / sumw2;
4320 rrv->setError((tau < 1e-15) ? 1e15 : (/*rrv->getVal()*/ 1. / sqrt(tau))); // not sure why was rrv->getVal()?
4321 rrv->setConstant(false);
4322 // parameter must be constrained
4323 auto _constr = v.constraints();
4324 // std::cout << " setting constraint " << v.GetName() << " nomin=" << tau << std::endl;
4325 if (_constr.empty()) {
4326 rrv->setStringAttribute("boundConstraint", _constr.Add("poisson").get()->GetName());
4327 } else {
4328 auto _glob = _constr.at(0)->obs().at(0)->get<RooRealVar>();
4329 // TODO: Update any globs snapshots that are designed to match the nominal
4330 _glob->setStringAttribute("nominal", TString::Format("%f", tau));
4331 double _min = tau * (1. - 5. * sqrt(1. / tau));
4332 double _max = tau * (1. + 5. * sqrt(1. / tau));
4333 _glob->setRange(_min, _max);
4334 _glob->setVal(tau);
4335 _constr.at(0)->pp().at(0)->SetBinContent(0, tau);
4336 rrv->setStringAttribute("boundConstraint", _constr.at(0)->get()->GetName());
4337 }
4338 rrv->setRange(std::max((1. - 5. * sqrt(1. / tau)), 1e-15), 1. + 5. * sqrt(1. / tau));
4339 } else {
4340 // remove constraint
4341 if (auto _constr = v.constraints(); !_constr.empty()) {
4342 v.constraints().Remove(*_constr.at(0));
4343 }
4344 // set const if sumw2 is 0 (i.e. no error)
4345 rrv->setVal(1);
4346 rrv->setError(0);
4347 rrv->setConstant(sumw2 == 0);
4348 }
4349 }
4350
4351 return true;
4352 }
4353
4354 throw std::runtime_error(TString::Format("%s SetBinError failed", GetName()));
4355}
4356
4357std::shared_ptr<xRooNode> xRooNode::at(const std::string &name, bool browseResult) const
4358{
4359 auto res = find(name, browseResult);
4360 if (res == nullptr)
4361 throw std::out_of_range(name + " does not exist");
4362 return res;
4363}
4364
4365////////////////////////////////////////////////////////////////////////////////
4366/// The RooWorkspace this node belong to, if any
4367
4369{
4370 if (auto _w = get<RooWorkspace>(); _w)
4371 return _w;
4372 if (auto a = get<RooAbsArg>(); a && GETWS(a)) {
4373 return GETWS(a);
4374 }
4375 if (fParent)
4376 return fParent->ws();
4377 return nullptr;
4378}
4379
4381{
4382
4383 xRooNode out(".constraints", nullptr, *