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