Logo ROOT  
Reference Guide
Loading...
Searching...
No Matches
RooTreeDataStore.cxx
Go to the documentation of this file.
1/*****************************************************************************
2 * Project: RooFit *
3 * Package: RooFitCore *
4 * @(#)root/roofitcore:$Id$
5 * Authors: *
6 * WV, Wouter Verkerke, UC Santa Barbara, verkerke@slac.stanford.edu *
7 * DK, David Kirkby, UC Irvine, dkirkby@uci.edu *
8 * *
9 * Copyright (c) 2000-2005, Regents of the University of California *
10 * and Stanford University. All rights reserved. *
11 * *
12 * Redistribution and use in source and binary forms, *
13 * with or without modification, are permitted according to the terms *
14 * listed in LICENSE (http://roofit.sourceforge.net/license.txt) *
15 *****************************************************************************/
16
17/**
18\file RooTreeDataStore.cxx
19\class RooTreeDataStore
20\ingroup Roofitcore
21
22TTree-backed data storage. When a file is opened before
23creating the data storage, the storage will be file-backed. This reduces memory
24pressure because it allows storing the data in the file and reading it on demand.
25For a completely memory-backed storage, which is faster than the file-backed storage,
26RooVectorDataStore can be used.
27
28With tree-backed storage, the tree can be found in the file with the name
29`RooTreeDataStore_name_title` for a dataset created as
30`RooDataSet("name", "title", ...)`.
31
32\note A file needs to be opened **before** creating the data storage to enable file-backed
33storage.
34```
35TFile outputFile("filename.root", "RECREATE");
36RooAbsData::setDefaultStorageType(RooAbsData::Tree);
37RooDataSet mydata(...);
38```
39
40One can also change between TTree- and std::vector-backed storage using
41RooAbsData::convertToTreeStore() and
42RooAbsData::convertToVectorStore().
43**/
44
45#include "RooTreeDataStore.h"
46
47#include "RooMsgService.h"
48#include "RooFormulaVar.h"
49#include "RooRealVar.h"
50#include "RooHistError.h"
51
52#include "ROOT/StringUtils.hxx"
53
54#include "TTree.h"
55#include "TFile.h"
56#include "TChain.h"
57#include "TDirectory.h"
58#include "TBuffer.h"
59#include "TBranch.h"
60#include "TROOT.h"
61
62#include <iomanip>
63using std::endl, std::list, std::string;
64
65
66
68
69
70
71////////////////////////////////////////////////////////////////////////////////
72
74
75
76
77////////////////////////////////////////////////////////////////////////////////
78/// Constructor to facilitate reading of legacy RooDataSets
79
80RooTreeDataStore::RooTreeDataStore(TTree* t, const RooArgSet& vars, const char* wgtVarName) :
81 RooAbsDataStore("blah","blah",varsNoWeight(vars,wgtVarName)),
82 _tree(t),
83 _defCtor(true),
84 _varsww(vars),
85 _wgtVar(weightVar(vars,wgtVarName))
86{
87}
88
89
90
91
92////////////////////////////////////////////////////////////////////////////////
93
94RooTreeDataStore::RooTreeDataStore(RooStringView name, RooStringView title, const RooArgSet& vars, const char* wgtVarName) :
95 RooAbsDataStore(name,title,varsNoWeight(vars,wgtVarName)),
96 _varsww(vars),
97 _wgtVar(weightVar(vars,wgtVarName))
98{
99 initialize() ;
100}
101
102
103////////////////////////////////////////////////////////////////////////////////
104
105RooTreeDataStore::RooTreeDataStore(RooStringView name, RooStringView title, const RooArgSet& vars, TTree& t, const char* selExpr, const char* wgtVarName) :
106 RooAbsDataStore(name,title,varsNoWeight(vars,wgtVarName)),
107 _varsww(vars),
108 _wgtVar(weightVar(vars,wgtVarName))
109{
110 initialize() ;
111
112 if (selExpr && *selExpr) {
113 // Create a RooFormulaVar cut from given cut expression
114 RooFormulaVar select(selExpr, selExpr, _vars, /*checkVariables=*/false);
115 loadValues(&t,&select);
116 } else {
117 loadValues(&t);
118 }
119}
120
121
122////////////////////////////////////////////////////////////////////////////////
123
124RooTreeDataStore::RooTreeDataStore(RooStringView name, RooStringView title, const RooArgSet& vars, const RooAbsDataStore& ads, const char* selExpr, const char* wgtVarName) :
125 RooAbsDataStore(name,title,varsNoWeight(vars,wgtVarName)),
126 _varsww(vars),
127 _wgtVar(weightVar(vars,wgtVarName))
128{
129 initialize() ;
130
131 if (selExpr && *selExpr) {
132 // Create a RooFormulaVar cut from given cut expression
133 RooFormulaVar select(selExpr, selExpr, _vars, /*checkVariables=*/false);
134 loadValues(&ads,&select);
135 } else {
136 loadValues(&ads);
137 }
138}
139
140
141
142
143////////////////////////////////////////////////////////////////////////////////
144
146 const RooFormulaVar *cutVar, const char *cutRange, Int_t nStart, Int_t nStop,
147 const char *wgtVarName)
148 : RooAbsDataStore(name, title, varsNoWeight(vars, wgtVarName)),
149 _varsww(vars),
150 _wgtVar(weightVar(vars, wgtVarName))
151{
152 // WVE NEED TO ADJUST THIS FOR WEIGHTS
153
154 // Protected constructor for internal use only
155
156 createTree(makeTreeName(), title);
157
158 // Deep clone cutVar and attach clone to this dataset
159 std::unique_ptr<RooFormulaVar> cloneVar;
160 if (cutVar) {
161 cloneVar.reset(static_cast<RooFormulaVar*>(cutVar->cloneTree()));
162 cloneVar->attachDataStore(tds) ;
163 }
164
165 // Constructor from existing data set with list of variables that preserves the cache
166 initialize();
167
168 attachCache(nullptr,(static_cast<RooTreeDataStore&>(tds))._cachedVars) ;
169
170 // WVE copy values of cached variables here!!!
171 _cacheTree->CopyEntries((static_cast<RooTreeDataStore&>(tds))._cacheTree) ;
172 _cacheOwner = nullptr ;
173
174 loadValues(&tds,cloneVar.get(),cutRange,nStart,nStop);
175}
176
177
178std::unique_ptr<RooAbsDataStore> RooTreeDataStore::reduce(RooStringView name, RooStringView title,
179 const RooArgSet& vars, const RooFormulaVar* cutVar, const char* cutRange,
180 std::size_t nStart, std::size_t nStop) {
181 RooArgSet tmp(vars) ;
182 if(_wgtVar && !tmp.contains(*_wgtVar)) {
183 tmp.add(*_wgtVar) ;
184 }
185 const char* wgtVarName = _wgtVar ? _wgtVar->GetName() : nullptr;
186 return std::make_unique<RooTreeDataStore>(name, title, *this, tmp, cutVar, cutRange, nStart, nStop, wgtVarName);
187}
188
189
190////////////////////////////////////////////////////////////////////////////////
191/// Utility function for constructors
192/// Return RooArgSet that is copy of allVars minus variable matching wgtName if specified
193
194RooArgSet RooTreeDataStore::varsNoWeight(const RooArgSet& allVars, const char* wgtName)
195{
196 RooArgSet ret(allVars) ;
197 if(wgtName) {
198 RooAbsArg* wgt = allVars.find(wgtName) ;
199 if (wgt) {
200 ret.remove(*wgt,true,true) ;
201 }
202 }
203 return ret ;
204}
205
206
207
208////////////////////////////////////////////////////////////////////////////////
209/// Utility function for constructors
210/// Return pointer to weight variable if it is defined
211
212RooRealVar* RooTreeDataStore::weightVar(const RooArgSet& allVars, const char* wgtName)
213{
214 if(wgtName) {
215 RooRealVar* wgt = dynamic_cast<RooRealVar*>(allVars.find(wgtName)) ;
216 return wgt ;
217 }
218 return nullptr ;
219}
220
221
222
223
224////////////////////////////////////////////////////////////////////////////////
225/// Initialize cache of dataset: attach variables of cache ArgSet
226/// to the corresponding TTree branches
227
228void RooTreeDataStore::attachCache(const RooAbsArg* newOwner, const RooArgSet& cachedVarsIn)
229{
230 // iterate over the cache variables for this dataset
231 _cachedVars.removeAll() ;
232 for (RooAbsArg * var : cachedVarsIn) {
233 var->attachToTree(*_cacheTree,_defTreeBufSize) ;
234 _cachedVars.add(*var) ;
235 }
236 _cacheOwner = newOwner ;
237
238}
239
240
241
242
243
244
245////////////////////////////////////////////////////////////////////////////////
246
247RooTreeDataStore::RooTreeDataStore(const RooTreeDataStore& other, const char* newname) :
248 RooAbsDataStore(other,newname),
249 _varsww(other._varsww),
250 _wgtVar(other._wgtVar),
255 _curWgt(other._curWgt),
259{
260 initialize() ;
261 loadValues(&other) ;
262}
263
264
265////////////////////////////////////////////////////////////////////////////////
266
267RooTreeDataStore::RooTreeDataStore(const RooTreeDataStore& other, const RooArgSet& vars, const char* newname) :
268 RooAbsDataStore(other,varsNoWeight(vars,other._wgtVar?other._wgtVar->GetName():nullptr),newname),
269 _varsww(vars),
270 _wgtVar(other._wgtVar?weightVar(vars,other._wgtVar->GetName()):nullptr),
275 _curWgt(other._curWgt),
279{
280 initialize() ;
281 loadValues(&other) ;
282}
283
284
285
286
287////////////////////////////////////////////////////////////////////////////////
288/// Destructor
289
291{
292 if (_tree) {
293 delete _tree ;
294 }
295 if (_cacheTree) {
296 delete _cacheTree ;
297 }
298}
299
300
301
302////////////////////////////////////////////////////////////////////////////////
303/// One-time initialization common to all constructor forms. Attach
304/// variables of internal ArgSet to the corresponding TTree branches
305
307{
308 // Recreate (empty) cache tree
310
311 // Attach each variable to the dataset
312 for (auto var : _varsww) {
313 var->attachToTree(*_tree,_defTreeBufSize) ;
314 }
315}
316
317
318
319
320
321////////////////////////////////////////////////////////////////////////////////
322/// Create TTree object that lives in memory, independent of current
323/// location of gDirectory
324
326{
327 if (!_tree) {
328 _tree = new TTree(name.c_str(),title.c_str());
329 _tree->ResetBit(kCanDelete);
330 _tree->ResetBit(kMustCleanup);
331 _tree->SetDirectory(nullptr);
332 }
333
334 TString pwd(gDirectory->GetPath()) ;
335 TString memDir(gROOT->GetName()) ;
336 memDir.Append(":/") ;
337 bool notInMemNow= (pwd!=memDir) ;
338
339 // std::cout << "RooTreeData::createTree pwd=" << pwd << " memDir=" << memDir << " notInMemNow = " << (notInMemNow?"T":"F") << std::endl ;
340
341 if (notInMemNow) {
342 gDirectory->cd(memDir) ;
343 }
344
345 if (!_cacheTree) {
346 _cacheTree = new TTree(TString{name.c_str()} + "_cacheTree", TString{title.c_str()});
347 _cacheTree->SetDirectory(nullptr) ;
348 gDirectory->RecursiveRemove(_cacheTree) ;
349 }
350
351 if (notInMemNow) {
352 gDirectory->cd(pwd) ;
353 }
354
355}
356
357
358
359
360////////////////////////////////////////////////////////////////////////////////
361/// Load values from tree 't' into this data collection, optionally
362/// selecting events using the RooFormulaVar 'select'.
363///
364/// The source tree 't' is cloned to not disturb its branch
365/// structure when retrieving information from it.
366void RooTreeDataStore::loadValues(const TTree *t, const RooFormulaVar* select, const char* /*rangeName*/, Int_t /*nStart*/, Int_t /*nStop*/)
367{
368 // Make our local copy of the tree, so we can safely loop through it.
369 // We need a custom deleter, because if we don't deregister the Tree from the directory
370 // of the original, it tears it down at destruction time!
371 auto deleter = [](TTree* tree){tree->SetDirectory(nullptr); delete tree;};
372 std::unique_ptr<TTree, decltype(deleter)> tClone(static_cast<TTree*>(t->Clone()), deleter);
373 tClone->SetDirectory(t->GetDirectory());
374
375 // Clone list of variables
376 RooArgSet sourceArgSet;
377 _varsww.snapshot(sourceArgSet, false);
378
379 // Check that we have the branches:
380 bool missingBranches = false;
381 for (const auto var : sourceArgSet) {
382 if (!tClone->GetBranch(var->GetName())) {
383 missingBranches = true;
384 coutE(InputArguments) << "Didn't find a branch in Tree '" << tClone->GetName() << "' to read variable '"
385 << var->GetName() << "' from."
386 << "\n\tNote: Name the RooFit variable the same as the branch." << std::endl;
387 }
388 }
389 if (missingBranches) {
390 coutE(InputArguments) << "Cannot import data from TTree '" << tClone->GetName()
391 << "' because some branches are missing !" << std::endl;
392 return;
393 }
394
395 // Attach args in cloned list to cloned source tree
396 for (const auto sourceArg : sourceArgSet) {
397 sourceArg->attachToTree(*tClone,_defTreeBufSize) ;
398 }
399
400 // Redirect formula servers to sourceArgSet
401 std::unique_ptr<RooFormulaVar> selectClone;
402 if (select) {
403 selectClone.reset( static_cast<RooFormulaVar*>(select->cloneTree()) );
404 selectClone->recursiveRedirectServers(sourceArgSet) ;
405 selectClone->setOperMode(RooAbsArg::ADirty,true) ;
406 }
407
408 // Loop over events in source tree
409 Int_t numInvalid(0) ;
410 const Long64_t nevent = tClone->GetEntries();
411 for(Long64_t i=0; i < nevent; ++i) {
412 const auto entryNumber = tClone->GetEntryNumber(i);
413 if (entryNumber<0) break;
414 tClone->GetEntry(entryNumber,1);
415
416 // Copy from source to destination
417 bool allOK(true) ;
418 for (unsigned int j=0; j < sourceArgSet.size(); ++j) {
419 auto destArg = _varsww[j];
420 const auto sourceArg = sourceArgSet[j];
421
422 destArg->copyCache(sourceArg) ;
423 sourceArg->copyCache(destArg) ;
424 if (!destArg->isValid()) {
425 numInvalid++ ;
426 allOK=false ;
427 if (numInvalid < 5) {
428 auto& log = coutI(DataHandling);
429 log << "RooTreeDataStore::loadValues(" << GetName() << ") Skipping event #" << i << " because " << destArg->GetName()
430 << " cannot accommodate the value ";
431 if(sourceArg->isCategory()) {
432 log << static_cast<RooAbsCategory*>(sourceArg)->getCurrentIndex();
433 } else {
434 log << static_cast<RooAbsReal*>(sourceArg)->getVal();
435 }
436 log << std::endl;
437 } else if (numInvalid == 5) {
438 coutI(DataHandling) << "RooTreeDataStore::loadValues(" << GetName() << ") Skipping ..." << std::endl;
439 }
440 break ;
441 }
442 }
443
444 // Does this event pass the cuts?
445 if (!allOK || (selectClone && selectClone->getVal()==0)) {
446 continue ;
447 }
448
449 fill() ;
450 }
451
452 if (numInvalid>0) {
453 coutW(DataHandling) << "RooTreeDataStore::loadValues(" << GetName() << ") Ignored " << numInvalid << " out-of-range events" << std::endl ;
454 }
455
456 SetTitle(t->GetTitle());
457}
458
459
460
461
462
463
464////////////////////////////////////////////////////////////////////////////////
465/// Load values from dataset 't' into this data collection, optionally
466/// selecting events using 'select' RooFormulaVar
467///
468
470 const char* rangeName, std::size_t nStart, std::size_t nStop)
471{
472 // Redirect formula servers to source data row
473 std::unique_ptr<RooFormulaVar> selectClone;
474 if (select) {
475 selectClone.reset( static_cast<RooFormulaVar*>(select->cloneTree()) );
476 selectClone->recursiveRedirectServers(*ads->get()) ;
477 selectClone->setOperMode(RooAbsArg::ADirty,true) ;
478 }
479
480 // Force RDS internal initialization
481 ads->get(0) ;
482
483 // Loop over events in source tree
484 const auto numEntr = static_cast<std::size_t>(ads->numEntries());
485 std::size_t nevent = nStop < numEntr ? nStop : numEntr;
486
487 auto TDS = dynamic_cast<const RooTreeDataStore*>(ads) ;
488 if (TDS) {
489 const_cast<RooTreeDataStore*>(TDS)->resetBuffers();
490 }
491
492 std::vector<std::string> ranges;
493 if (rangeName) {
494 ranges = ROOT::Split(rangeName, ",");
495 }
496
497 for (auto i=nStart; i < nevent ; ++i) {
498 ads->get(i) ;
499
500 // Does this event pass the cuts?
501 if (selectClone && selectClone->getVal()==0) {
502 continue ;
503 }
504
505
506 if (TDS) {
507 _varsww.assignValueOnly(TDS->_varsww) ;
508 } else {
509 _varsww.assignValueOnly(*ads->get()) ;
510 }
511
512 // Check that all copied values are valid
513 bool allValid = true;
514 for (const auto arg : _varsww) {
515 allValid = arg->isValid() && (ranges.empty() || std::any_of(ranges.begin(), ranges.end(),
516 [arg](const std::string& range){return arg->inRange(range.c_str());}) );
517 if (!allValid)
518 break ;
519 }
520
521 if (!allValid) {
522 continue ;
523 }
524
525 _cachedVars.assign(static_cast<RooTreeDataStore const*>(ads)->_cachedVars) ;
526 fill() ;
527 }
528
529 if (TDS) {
530 const_cast<RooTreeDataStore*>(TDS)->restoreAlternateBuffers();
531 }
532
533 SetTitle(ads->GetTitle());
534}
535
536
537////////////////////////////////////////////////////////////////////////////////
538/// Interface function to TTree::Fill
539
541{
542 return _tree->Fill() ;
543}
544
545
546
547////////////////////////////////////////////////////////////////////////////////
548/// Load the n-th data point (n='index') in memory
549/// and return a pointer to the internal RooArgSet
550/// holding its coordinates.
551
553{
554 checkInit() ;
555
556 Int_t ret = _tree->GetEntry(index, 1);
557 _cacheTree->GetEntry(index,1);
558
559 if(!ret) return nullptr;
560
561 if (_doDirtyProp) {
562 // Raise all dirty flags
563 for (auto var : _vars) {
564 var->setValueDirty(); // This triggers recalculation of all clients
565 }
566
567 for (auto var : _cachedVars) {
568 var->setValueDirty(); // This triggers recalculation of all clients, but doesn't recalculate self
569 var->clearValueDirty();
570 }
571 }
572
573 // Update current weight cache
574 if (_extWgtArray) {
575
576 // If external array is specified use that
577 _curWgt = _extWgtArray[index] ;
580 _curWgtErr = sqrt( _extSumW2Array ? _extSumW2Array[index] : _extWgtArray[index] );
581
582 } else if (_wgtVar) {
583
584 // Otherwise look for weight variable
585 _curWgt = _wgtVar->getVal() ;
586 _curWgtErrLo = _wgtVar->getAsymErrorLo() ;
587 _curWgtErrHi = _wgtVar->getAsymErrorHi() ;
588 _curWgtErr = _wgtVar->hasAsymError() ? ((_wgtVar->getAsymErrorHi() - _wgtVar->getAsymErrorLo())/2) : _wgtVar->getError() ;
589
590 } else {
591
592 // Otherwise return 1
593 _curWgt=1.0 ;
594 _curWgtErrLo = 0 ;
595 _curWgtErrHi = 0 ;
596 _curWgtErr = 0 ;
597
598 }
599
600 return &_vars;
601}
602
603
604////////////////////////////////////////////////////////////////////////////////
605/// Return the weight of the n-th data point (n='index') in memory
606
608{
609 return _curWgt ;
610}
611
612
613////////////////////////////////////////////////////////////////////////////////
614
616{
617 if (_extWgtArray) {
618
619 // We have a weight array, use that info
620
621 // Return symmetric error on current bin calculated either from Poisson statistics or from SumOfWeights
622 double lo = 0;
623 double hi = 0;
624 weightError(lo,hi,etype) ;
625 return (lo+hi)/2 ;
626
627 } else if (_wgtVar) {
628
629 // We have a weight variable, use that info
630 if (_wgtVar->hasAsymError()) {
631 return ( _wgtVar->getAsymErrorHi() - _wgtVar->getAsymErrorLo() ) / 2 ;
632 } else {
633 return _wgtVar->getError() ;
634 }
635
636 }
637
638 // We have no weights
639 return 0.0;
640}
641
642
643
644////////////////////////////////////////////////////////////////////////////////
645
646void RooTreeDataStore::weightError(double& lo, double& hi, RooAbsData::ErrorType etype) const
647{
648 if (_extWgtArray) {
649
650 // We have a weight array, use that info
651 switch (etype) {
652
653 case RooAbsData::Auto:
654 throw string(Form("RooDataHist::weightError(%s) error type Auto not allowed here",GetName())) ;
655 break ;
656
658 throw string(Form("RooDataHist::weightError(%s) error type Expected not allowed here",GetName())) ;
659 break ;
660
662 // Weight may be preset or precalculated
663 if (_curWgtErrLo>=0) {
664 lo = _curWgtErrLo ;
665 hi = _curWgtErrHi ;
666 return ;
667 }
668
669 // Otherwise Calculate poisson errors
670 double ym;
671 double yp;
673 lo = weight()-ym ;
674 hi = yp-weight() ;
675 return ;
676
678 lo = _curWgtErr ;
679 hi = _curWgtErr ;
680 return ;
681
682 case RooAbsData::None:
683 lo = 0 ;
684 hi = 0 ;
685 return ;
686 }
687
688 } else if (_wgtVar) {
689
690 // We have a weight variable, use that info
691 if (_wgtVar->hasAsymError()) {
692 hi = _wgtVar->getAsymErrorHi() ;
693 lo = _wgtVar->getAsymErrorLo() ;
694 } else {
695 hi = _wgtVar->getError() ;
696 lo = _wgtVar->getError() ;
697 }
698
699 } else {
700
701 // We are unweighted
702 lo=0 ;
703 hi=0 ;
704
705 }
706}
707
708
709////////////////////////////////////////////////////////////////////////////////
710/// Change name of internal observable named 'from' into 'to'
711
712bool RooTreeDataStore::changeObservableName(const char* from, const char* to)
713{
714 // Find observable to be changed
715 RooAbsArg* var = _vars.find(from) ;
716
717 // Check that we found it
718 if (!var) {
719 coutE(InputArguments) << "RooTreeDataStore::changeObservableName(" << GetName() << " no observable " << from << " in this dataset" << std::endl ;
720 return true ;
721 }
722
723 // Process name change
724 TString oldBranchName = var->cleanBranchName() ;
725 var->SetName(to) ;
726
727 // Change the branch name as well
728 if (_tree->GetBranch(oldBranchName.Data())) {
729
730 // Simple case varName = branchName
731 _tree->GetBranch(oldBranchName.Data())->SetName(var->cleanBranchName().Data()) ;
732
733 // Process any error branch if existing
734 if (_tree->GetBranch(Form("%s_err",oldBranchName.Data()))) {
735 _tree->GetBranch(Form("%s_err",oldBranchName.Data()))->SetName(Form("%s_err",var->cleanBranchName().Data())) ;
736 }
737 if (_tree->GetBranch(Form("%s_aerr_lo",oldBranchName.Data()))) {
738 _tree->GetBranch(Form("%s_aerr_lo",oldBranchName.Data()))->SetName(Form("%s_aerr_lo",var->cleanBranchName().Data())) ;
739 }
740 if (_tree->GetBranch(Form("%s_aerr_hi",oldBranchName.Data()))) {
741 _tree->GetBranch(Form("%s_aerr_hi",oldBranchName.Data()))->SetName(Form("%s_aerr_hi",var->cleanBranchName().Data())) ;
742 }
743
744 } else {
745
746 // Native category case branchNames = varName_idx and varName_lbl
747 if (_tree->GetBranch(Form("%s_idx",oldBranchName.Data()))) {
748 _tree->GetBranch(Form("%s_idx",oldBranchName.Data()))->SetName(Form("%s_idx",var->cleanBranchName().Data())) ;
749 }
750 if (_tree->GetBranch(Form("%s_lbl",oldBranchName.Data()))) {
751 _tree->GetBranch(Form("%s_lbl",oldBranchName.Data()))->SetName(Form("%s_lb",var->cleanBranchName().Data())) ;
752 }
753
754 }
755
756 return false ;
757}
758
759
760
761////////////////////////////////////////////////////////////////////////////////
762/// Add a new column to the data set which holds the pre-calculated values
763/// of 'newVar'. This operation is only meaningful if 'newVar' is a derived
764/// value.
765///
766/// The return value points to the added element holding 'newVar's value
767/// in the data collection. The element is always the corresponding fundamental
768/// type of 'newVar' (e.g. a RooRealVar if 'newVar' is a RooFormulaVar)
769///
770/// Note: This function is explicitly NOT intended as a speed optimization
771/// opportunity for the user. Components of complex PDFs that can be
772/// precalculated with the dataset are automatically identified as such
773/// and will be precalculated when fitting to a dataset
774///
775/// By forcibly precalculating functions with non-trivial Jacobians,
776/// or functions of multiple variables occurring in the data set,
777/// using addColumn(), you may alter the outcome of the fit.
778///
779/// Only in cases where such a modification of fit behaviour is intentional,
780/// this function should be used.
781
783{
784 checkInit() ;
785
786 // Create a fundamental object of the right type to hold newVar values
787 auto valHolder = std::unique_ptr<RooAbsArg>{newVar.createFundamental()}.release();
788 // Sanity check that the holder really is fundamental
789 if(!valHolder->isFundamental()) {
790 coutE(InputArguments) << GetName() << "::addColumn: holder argument is not fundamental: \""
791 << valHolder->GetName() << "\"" << std::endl;
792 return nullptr;
793 }
794
795 // WVE need to reset TTRee buffers to original datamembers here
796 resetBuffers() ;
797
798 // Clone variable and attach to cloned tree
799 std::unique_ptr<RooAbsArg> newVarClone{newVar.cloneTree()};
800 newVarClone->recursiveRedirectServers(_vars,false) ;
801
802 // Attach value place holder to this tree
803 ((RooAbsArg*)valHolder)->attachToTree(*_tree,_defTreeBufSize) ;
804 _vars.add(*valHolder) ;
805 _varsww.add(*valHolder) ;
806
807
808 // Fill values of placeholder
809 for (int i=0 ; i < _tree->GetEntries() ; i++) {
810 get(i) ;
811
812 newVarClone->syncCache(&_vars) ;
813 valHolder->copyCache(newVarClone.get());
814 valHolder->fillTreeBranch(*_tree) ;
815 }
816
817 // WVE need to restore TTRee buffers to previous values here
819
820 if (adjustRange) {
821// // Set range of valHolder to (just) bracket all values stored in the dataset
822// double vlo,vhi ;
823// RooRealVar* rrvVal = dynamic_cast<RooRealVar*>(valHolder) ;
824// if (rrvVal) {
825// getRange(*rrvVal,vlo,vhi,0.05) ;
826// rrvVal->setRange(vlo,vhi) ;
827// }
828 }
829
830 return valHolder ;
831}
832
833
834////////////////////////////////////////////////////////////////////////////////
835/// Merge columns of supplied data set(s) with this data set. All
836/// data sets must have equal number of entries. In case of
837/// duplicate columns the column of the last dataset in the list
838/// prevails
839
840RooAbsDataStore* RooTreeDataStore::merge(const RooArgSet& allVars, list<RooAbsDataStore*> dstoreList)
841{
842 RooTreeDataStore* mergedStore = new RooTreeDataStore("merged","merged",allVars) ;
843
844 Int_t nevt = dstoreList.front()->numEntries() ;
845 for (int i=0 ; i<nevt ; i++) {
846
847 // Cope data from self
848 mergedStore->_vars.assign(*get(i)) ;
849
850 // Copy variables from merge sets
851 for (list<RooAbsDataStore*>::iterator iter = dstoreList.begin() ; iter!=dstoreList.end() ; ++iter) {
852 const RooArgSet* partSet = (*iter)->get(i) ;
853 mergedStore->_vars.assign(*partSet) ;
854 }
855
856 mergedStore->fill() ;
857 }
858 return mergedStore ;
859}
860
861
862
863
864
865////////////////////////////////////////////////////////////////////////////////
866
868{
869 Int_t nevt = other.numEntries() ;
870 for (int i=0 ; i<nevt ; i++) {
871 _vars.assign(*other.get(i)) ;
872 if (_wgtVar) {
873 _wgtVar->setVal(other.weight()) ;
874 }
875
876 fill() ;
877 }
878}
879
880
881////////////////////////////////////////////////////////////////////////////////
882
884{
885 if (_wgtVar) {
886
887 double sum(0);
888 double carry(0);
889 Int_t nevt = numEntries() ;
890 for (int i=0 ; i<nevt ; i++) {
891 get(i) ;
892 // Kahan's algorithm for summing to avoid loss of precision
893 double y = _wgtVar->getVal() - carry;
894 double t = sum + y;
895 carry = (t - sum) - y;
896 sum = t;
897 }
898 return sum ;
899
900 } else if (_extWgtArray) {
901
902 double sum(0);
903 double carry(0);
904 Int_t nevt = numEntries() ;
905 for (int i=0 ; i<nevt ; i++) {
906 // Kahan's algorithm for summing to avoid loss of precision
907 double y = _extWgtArray[i] - carry;
908 double t = sum + y;
909 carry = (t - sum) - y;
910 sum = t;
911 }
912 return sum ;
913
914 } else {
915
916 return numEntries() ;
917
918 }
919}
920
921
922
923
924////////////////////////////////////////////////////////////////////////////////
925
927{
928 return _tree->GetEntries() ;
929}
930
931
932
933////////////////////////////////////////////////////////////////////////////////
934
936{
937 _tree->Reset() ;
938}
939
940
941
942////////////////////////////////////////////////////////////////////////////////
943/// Cache given RooAbsArgs with this tree: The tree is
944/// given direct write access of the args internal cache
945/// the args values is pre-calculated for all data points
946/// in this data collection. Upon a get() call, the
947/// internal cache of 'newVar' will be loaded with the
948/// precalculated value and it's dirty flag will be cleared.
949
950void RooTreeDataStore::cacheArgs(const RooAbsArg* owner, RooArgSet& newVarSet, const RooArgSet* nset, bool /*skipZeroWeights*/)
951{
952 checkInit() ;
953
954 _cacheOwner = owner ;
955
956 std::unique_ptr<RooArgSet> constExprVarSet{newVarSet.selectByAttrib("ConstantExpression", true)};
957
958 bool doTreeFill = (_cachedVars.empty()) ;
959
960 for (RooAbsArg * arg : *constExprVarSet) {
961 // Attach original newVar to this tree
962 arg->attachToTree(*_cacheTree,_defTreeBufSize) ;
963 //arg->recursiveRedirectServers(_vars) ;
964 _cachedVars.add(*arg) ;
965 }
966
967 // WVE need to reset TTRee buffers to original datamembers here
968 //resetBuffers() ;
969
970 // Refill regular and cached variables of current tree from clone
971 for (int i=0 ; i < _tree->GetEntries() ; i++) {
972 get(i) ;
973
974 // Evaluate the cached variables and store the results
975 for (RooAbsArg * arg : *constExprVarSet) {
976 arg->setValueDirty() ;
977 arg->syncCache(nset) ;
978 if (!doTreeFill) {
979 arg->fillTreeBranch(*_cacheTree) ;
980 }
981 }
982
983 if (doTreeFill) {
984 _cacheTree->Fill() ;
985 }
986 }
987
988 // WVE need to restore TTRee buffers to previous values here
989 //restoreAlternateBuffers() ;
990}
991
992
993
994
995////////////////////////////////////////////////////////////////////////////////
996/// Activate or deactivate the branch status of the TTree branch associated
997/// with the given set of dataset observables
998
1000{
1001 for (RooAbsArg * arg : set) {
1002 RooAbsArg* depArg = _vars.find(arg->GetName()) ;
1003 if (!depArg) {
1004 coutE(InputArguments) << "RooTreeDataStore::setArgStatus(" << GetName()
1005 << ") dataset doesn't contain variable " << arg->GetName() << std::endl ;
1006 continue ;
1007 }
1008 depArg->setTreeBranchStatus(*_tree,active) ;
1009 }
1010}
1011
1012
1013
1014////////////////////////////////////////////////////////////////////////////////
1015/// Remove tree with values of cached observables
1016/// and clear list of cached observables
1017
1019{
1020 // Empty list of cached functions
1021 _cachedVars.removeAll() ;
1022
1023 // Delete & recreate cache tree
1024 delete _cacheTree ;
1025 _cacheTree = nullptr ;
1026 createTree(makeTreeName().c_str(), GetTitle());
1027
1028 return ;
1029}
1030
1031
1032
1033
1034////////////////////////////////////////////////////////////////////////////////
1035
1037{
1038 _attachedBuffers.removeAll() ;
1039 for (const auto arg : _varsww) {
1040 RooAbsArg* extArg = extObs.find(arg->GetName()) ;
1041 if (extArg) {
1042 if (arg->getAttribute("StoreError")) {
1043 extArg->setAttribute("StoreError") ;
1044 }
1045 if (arg->getAttribute("StoreAsymError")) {
1046 extArg->setAttribute("StoreAsymError") ;
1047 }
1048 extArg->attachToTree(*_tree) ;
1049 _attachedBuffers.add(*extArg) ;
1050 }
1051 }
1052}
1053
1054
1055
1056////////////////////////////////////////////////////////////////////////////////
1057
1059{
1060 for(RooAbsArg * arg : _varsww) {
1061 arg->attachToTree(*_tree) ;
1062 }
1063}
1064
1065
1066
1067////////////////////////////////////////////////////////////////////////////////
1068
1070{
1071 for(RooAbsArg * arg : _attachedBuffers) {
1072 arg->attachToTree(*_tree) ;
1073 }
1074}
1075
1076
1077
1078////////////////////////////////////////////////////////////////////////////////
1079
1081{
1082 if (_defCtor) {
1083 const_cast<RooTreeDataStore*>(this)->initialize() ;
1084 _defCtor = false ;
1085 }
1086}
1087
1088////////////////////////////////////////////////////////////////////////////////
1089/// Stream an object of class RooTreeDataStore.
1090
1092{
1093 if (R__b.IsReading()) {
1094 UInt_t R__s;
1095 UInt_t R__c;
1096 const Version_t R__v = R__b.ReadVersion(&R__s, &R__c);
1097
1098 R__b.ReadClassBuffer(RooTreeDataStore::Class(), this, R__v, R__s, R__c);
1099
1100 if (!_tree) {
1101 // If the tree has not been deserialised automatically, it is time to load
1102 // it now.
1103 TFile* parent = dynamic_cast<TFile*>(R__b.GetParent());
1104 assert(parent);
1105 parent->GetObject(makeTreeName().c_str(), _tree);
1106 }
1107
1108 initialize();
1109
1110 } else {
1111
1112 TTree* tmpTree = _tree;
1113 auto parent = dynamic_cast<TDirectory*>(R__b.GetParent());
1114 if (_tree && parent) {
1115 // Large trees cannot be written because of the 1Gb I/O limitation.
1116 // Here, we take the tree away from our instance, write it, and continue
1117 // to write the rest of the class normally
1118 auto tmpDir = _tree->GetDirectory();
1119
1120 _tree->SetDirectory(parent);
1121 _tree->FlushBaskets(false);
1122 parent->WriteObject(_tree, makeTreeName().c_str());
1123 _tree->SetDirectory(tmpDir);
1124 _tree = nullptr;
1125 }
1126
1128
1129 _tree = tmpTree;
1130 }
1131}
1132
1133////////////////////////////////////////////////////////////////////////////////
1134/// Generate a name for the storage tree from the name and title of this instance.
1136 std::string title = GetTitle();
1137 std::replace(title.begin(), title.end(), ' ', '_');
1138 std::replace(title.begin(), title.end(), '-', '_');
1139 return std::string("RooTreeDataStore_") + GetName() + "_" + title;
1140}
1141
1142
1143////////////////////////////////////////////////////////////////////////////////
1144/// Get the weights of the events in the range [first, first+len).
1145/// This implementation will fill a vector with every event retrieved one by one
1146/// (even if the weight is constant). Then, it returns a span.
1147std::span<const double> RooTreeDataStore::getWeightBatch(std::size_t first, std::size_t len) const {
1148
1149 if (_extWgtArray) {
1150 return {_extWgtArray + first, len};
1151 }
1152
1153 if (!_weightBuffer) {
1154 _weightBuffer = std::make_unique<std::vector<double>>();
1155 _weightBuffer->reserve(len);
1156
1157 for (int i = 0; i < _tree->GetEntries(); ++i) {
1158 _weightBuffer->push_back(weight(i));
1159 }
1160 }
1161
1162 return {_weightBuffer->data() + first, len};
1163}
#define coutI(a)
#define coutW(a)
#define coutE(a)
char * ret
Definition Rotated.cxx:221
int Int_t
Signed integer 4 bytes (int).
Definition RtypesCore.h:59
short Version_t
Class version identifier (short).
Definition RtypesCore.h:79
unsigned int UInt_t
Unsigned integer 4 bytes (unsigned int).
Definition RtypesCore.h:60
long long Long64_t
Portable signed long integer 8 bytes.
Definition RtypesCore.h:83
return
#define gDirectory
Definition TDirectory.h:385
char name[80]
Definition TGX11.cxx:148
#define hi
#define gROOT
Definition TROOT.h:417
char * Form(const char *fmt,...)
Formats a string in a circular formatting buffer.
Definition TString.cxx:2496
Common abstract base class for objects that represent a value and a "shape" in RooFit.
Definition RooAbsArg.h:76
void SetName(const char *name) override
Set the name of the TNamed.
virtual RooFit::OwningPtr< RooAbsArg > createFundamental(const char *newname=nullptr) const =0
Create a fundamental-type object that stores our type of value.
virtual void attachToTree(TTree &t, Int_t bufSize=32000)=0
Overloadable function for derived classes to implement attachment as branch to a TTree.
virtual void setTreeBranchStatus(TTree &t, bool active)=0
virtual RooAbsArg * cloneTree(const char *newname=nullptr) const
Clone tree expression of objects.
TString cleanBranchName() const
Construct a mangled name from the actual name that is free of any math symbols that might be interpre...
void setAttribute(const Text_t *name, bool value=true)
Set (default) or clear a named boolean attribute of this object.
bool contains(const char *name) const
Check if collection contains an argument with a specific name.
Storage_t const & get() const
Const access to the underlying stl container.
virtual bool add(const RooAbsArg &var, bool silent=false)
Add the specified argument to list.
void assign(const RooAbsCollection &other) const
Sets the value, cache and constant attribute of any argument in our set that also appears in the othe...
Storage_t::size_type size() const
RooAbsArg * find(const char *name) const
Find object with given name in list.
virtual const RooArgSet * get(Int_t index) const =0
bool _doDirtyProp
Switch do (de)activate dirty state propagation when loading a data point.
virtual double weight() const =0
virtual Int_t numEntries() const =0
RooArgSet is a container object that can hold multiple RooAbsArg objects.
Definition RooArgSet.h:24
RooArgSet * selectByAttrib(const char *name, bool value) const
Use RooAbsCollection::selectByAttrib(), but return as RooArgSet.
Definition RooArgSet.h:144
A RooFormulaVar is a generic implementation of a real-valued object, which takes a RooArgList of serv...
static const RooHistError & instance()
Return a reference to a singleton object that is created the first time this method is called.
bool getPoissonInterval(Int_t n, double &mu1, double &mu2, double nSigma=1) const
Calculate a confidence interval for the expected number of events given n observed (unweighted) event...
Variable that can be changed from the outside.
Definition RooRealVar.h:37
The RooStringView is a wrapper around a C-style string that can also be constructed from a std::strin...
void initialize()
One-time initialization common to all constructor forms.
double _curWgtErr
Weight of current event.
double weightError(RooAbsData::ErrorType etype=RooAbsData::Poisson) const override
void resetBuffers() override
double _curWgt
Weight of current event.
void attachCache(const RooAbsArg *newOwner, const RooArgSet &cachedVars) override
Initialize cache of dataset: attach variables of cache ArgSet to the corresponding TTree branches.
double _curWgtErrHi
Weight of current event.
static TClass * Class()
Int_t numEntries() const override
std::string makeTreeName() const
Generate a name for the storage tree from the name and title of this instance.
RooArgSet varsNoWeight(const RooArgSet &allVars, const char *wgtName=nullptr)
Utility function for constructors Return RooArgSet that is copy of allVars minus variable matching wg...
~RooTreeDataStore() override
Destructor.
void createTree(RooStringView name, RooStringView title)
Create TTree object that lives in memory, independent of current location of gDirectory.
const double * _extWgtErrHiArray
! External weight array - high error
void attachBuffers(const RooArgSet &extObs) override
void reset() override
TTree * _cacheTree
! TTree holding the cached function values
RooAbsDataStore * merge(const RooArgSet &allvars, std::list< RooAbsDataStore * > dstoreList) override
Merge columns of supplied data set(s) with this data set.
static Int_t _defTreeBufSize
RooArgSet _attachedBuffers
! Currently attached buffers (if different from _varsww)
Int_t fill() override
Interface function to TTree::Fill.
double sumEntries() const override
std::unique_ptr< RooAbsDataStore > reduce(RooStringView name, RooStringView title, const RooArgSet &vars, const RooFormulaVar *cutVar, const char *cutRange, std::size_t nStart, std::size_t nStop) override
bool _defCtor
! Was object constructed with default ctor?
RooAbsArg * addColumn(RooAbsArg &var, bool adjustRange=true) override
Add a new column to the data set which holds the pre-calculated values of 'newVar'.
double weight() const override
Return the weight of the n-th data point (n='index') in memory.
void loadValues(const TTree *t, const RooFormulaVar *select=nullptr, const char *rangeName=nullptr, Int_t nStart=0, Int_t nStop=2000000000)
Load values from tree 't' into this data collection, optionally selecting events using the RooFormula...
void append(RooAbsDataStore &other) override
std::span< const double > getWeightBatch(std::size_t first, std::size_t len) const override
Get the weights of the events in the range [first, first+len).
const double * _extWgtErrLoArray
! External weight array - low error
void checkInit() const override
std::unique_ptr< std::vector< double > > _weightBuffer
! Buffer for weights in case a batch of values is requested.
const double * _extSumW2Array
! External sum of weights array
void Streamer(TBuffer &) override
Stream an object of class RooTreeDataStore.
RooRealVar * weightVar(const RooArgSet &allVars, const char *wgtName=nullptr)
Utility function for constructors Return pointer to weight variable if it is defined.
const double * _extWgtArray
! External weight array
double _curWgtErrLo
Weight of current event.
bool changeObservableName(const char *from, const char *to) override
Change name of internal observable named 'from' into 'to'.
void resetCache() override
Remove tree with values of cached observables and clear list of cached observables.
void setArgStatus(const RooArgSet &set, bool active) override
Activate or deactivate the branch status of the TTree branch associated with the given set of dataset...
void cacheArgs(const RooAbsArg *owner, RooArgSet &varSet, const RooArgSet *nset=nullptr, bool skipZeroWeights=false) override
Cache given RooAbsArgs with this tree: The tree is given direct write access of the args internal cac...
virtual const RooArgSet * get() const
const RooAbsArg * _cacheOwner
! Object owning cache contents
Buffer base class used for serializing objects.
Definition TBuffer.h:43
virtual Version_t ReadVersion(UInt_t *start=nullptr, UInt_t *bcnt=nullptr, const TClass *cl=nullptr)=0
TObject * GetParent() const
Return pointer to parent of this buffer.
Definition TBuffer.cxx:261
virtual Int_t ReadClassBuffer(const TClass *cl, void *pointer, const TClass *onfile_class=nullptr)=0
Bool_t IsReading() const
Definition TBuffer.h:86
virtual Int_t WriteClassBuffer(const TClass *cl, void *pointer)=0
Describe directory structure in memory.
Definition TDirectory.h:45
void GetObject(const char *namecycle, T *&ptr)
Get an object with proper type checking.
Definition TDirectory.h:213
A file, usually with extension .root, that stores data and code in the form of serialized objects in ...
Definition TFile.h:130
TObject * Clone(const char *newname="") const override
Make a clone of an object using the Streamer facility.
Definition TNamed.cxx:73
virtual void SetTitle(const char *title="")
Set the title of the TNamed.
Definition TNamed.cxx:173
const char * GetName() const override
Returns name of object.
Definition TNamed.h:49
const char * GetTitle() const override
Returns title of object.
Definition TNamed.h:50
@ kCanDelete
if object in a list can be deleted
Definition TObject.h:71
@ kMustCleanup
if object destructor must call RecursiveRemove()
Definition TObject.h:73
Basic string class.
Definition TString.h:138
const char * Data() const
Definition TString.h:384
TString & Append(const char *cs)
Definition TString.h:581
A TTree represents a columnar dataset.
Definition TTree.h:89
TDirectory * GetDirectory() const
Definition TTree.h:509
virtual void SetDirectory(TDirectory *dir)
Change the tree's directory.
Definition TTree.cxx:9220
STL class.
Double_t y[n]
Definition legend1.C:17
std::vector< std::string > Split(std::string_view str, std::string_view delims, bool skipEmpty=false)
Splits a string at each character in delims.
static uint64_t sum(uint64_t i)
Definition Factory.cxx:2338