Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
RooDataHist.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 RooDataHist.cxx
19\class RooDataHist
20\ingroup Roofitcore
21
22The RooDataHist is a container class to hold N-dimensional binned data. Each bin's central
23coordinates in N-dimensional space are represented by a RooArgSet containing RooRealVar, RooCategory
24or RooStringVar objects, thus data can be binned in real and/or discrete dimensions.
25
26There is an unbinned equivalent, RooDataSet.
27
28### Inspecting a datahist
29Inspect a datahist using Print() to get the coordinates and `weight()` to get the bin contents:
30```
31datahist->Print("V");
32datahist->get(0)->Print("V"); std::cout << "w=" << datahist->weight(0) << std::endl;
33datahist->get(1)->Print("V"); std::cout << "w=" << datahist->weight(1) << std::endl;
34...
35```
36
37### Plotting data.
38See RooAbsData::plotOn().
39
40### Creating a datahist using RDataFrame
41\see RooAbsDataHelper, rf408_RDataFrameToRooFit.C
42
43**/
44
45#include "RooDataHist.h"
46
47#include "Riostream.h"
48#include "RooMsgService.h"
50#include "RooAbsLValue.h"
51#include "RooArgList.h"
52#include "RooRealVar.h"
53#include "RooMath.h"
54#include "RooBinning.h"
55#include "RooPlot.h"
56#include "RooHistError.h"
57#include "RooCategory.h"
58#include "RooCmdConfig.h"
59#include "RooLinkedListIter.h"
60#include "RooTreeDataStore.h"
61#include "RooVectorDataStore.h"
62#include "RooTrace.h"
63#include "RooFormulaVar.h"
64#include "RooFormula.h"
65#include "RooUniformBinning.h"
66
67#include <ROOT/RSpan.hxx>
68#include <ROOT/StringUtils.hxx>
69
70#include "TAxis.h"
71#include "TH1.h"
72#include "TTree.h"
73#include "TBuffer.h"
74#include "TMath.h"
75#include "Math/Util.h"
76
77using namespace std;
78
80
81
82////////////////////////////////////////////////////////////////////////////////
83/// Default constructor
84
86{
88}
89
90
91std::unique_ptr<RooAbsDataStore>
93{
95 ? static_cast<std::unique_ptr<RooAbsDataStore>>(std::make_unique<RooTreeDataStore>(name, title, vars))
96 : static_cast<std::unique_ptr<RooAbsDataStore>>(std::make_unique<RooVectorDataStore>(name, title, vars));
97}
98
99
100////////////////////////////////////////////////////////////////////////////////
101/// Constructor of an empty data hist from a RooArgSet defining the dimensions
102/// of the data space. The range and number of bins in each dimensions are taken
103/// from getMin()getMax(),getBins() of each RooAbsArg representing that
104/// dimension.
105///
106/// For real dimensions, the fit range and number of bins can be set independently
107/// of the plot range and number of bins, but it is advisable to keep the
108/// ratio of the plot bin width and the fit bin width an integer value.
109/// For category dimensions, the fit ranges always comprises all defined states
110/// and each state is always has its individual bin
111///
112/// To effectively bin real dimensions with variable bin sizes,
113/// construct a RooThresholdCategory of the real dimension to be binned variably.
114/// Set the thresholds at the desired bin boundaries, and construct the
115/// data hist as a function of the threshold category instead of the real variable.
116RooDataHist::RooDataHist(RooStringView name, RooStringView title, const RooArgSet& vars, const char* binningName) :
117 RooAbsData(name,title,vars)
118{
119 // Initialize datastore
121
122 initialize(binningName) ;
123
125
126 appendToDir(this,true) ;
128}
129
130
131
132////////////////////////////////////////////////////////////////////////////////
133/// Constructor of a data hist from an existing data collection (binned or unbinned)
134/// The RooArgSet 'vars' defines the dimensions of the histogram.
135/// The range and number of bins in each dimensions are taken
136/// from getMin(), getMax(), getBins() of each argument passed.
137///
138/// For real dimensions, the fit range and number of bins can be set independently
139/// of the plot range and number of bins, but it is advisable to keep the
140/// ratio of the plot bin width and the fit bin width an integer value.
141/// For category dimensions, the fit ranges always comprises all defined states
142/// and each state is always has its individual bin
143///
144/// To effectively bin real dimensions with variable bin sizes,
145/// construct a RooThresholdCategory of the real dimension to be binned variably.
146/// Set the thresholds at the desired bin boundaries, and construct the
147/// data hist as a function of the threshold category instead of the real variable.
148///
149/// If the constructed data hist has less dimensions that in source data collection,
150/// all missing dimensions will be projected.
151
153 RooDataHist(name,title,vars)
154{
155 add(data,static_cast<const RooFormulaVar*>(nullptr),wgt);
156}
157
158
159
160////////////////////////////////////////////////////////////////////////////////
161/// Constructor of a data hist from a map of TH1,TH2 or TH3 that are collated into a x+1 dimensional
162/// RooDataHist where the added dimension is a category that labels the input source as defined
163/// in the histMap argument. The state names used in histMap must correspond to predefined states
164/// 'indexCat'
165///
166/// The RooArgList 'vars' defines the dimensions of the histogram.
167/// The ranges and number of bins are taken from the input histogram and must be the same in all histograms
168
170 std::map<string,TH1*> histMap, double wgt) :
171 RooAbsData(name,title,RooArgSet(vars,&indexCat))
172{
173 // Initialize datastore
175
176 importTH1Set(vars, indexCat, histMap, wgt, false) ;
177
180}
181
182
183
184////////////////////////////////////////////////////////////////////////////////
185/// Constructor of a data hist from a map of RooDataHists that are collated into a x+1 dimensional
186/// RooDataHist where the added dimension is a category that labels the input source as defined
187/// in the histMap argument. The state names used in histMap must correspond to predefined states
188/// 'indexCat'
189///
190/// The RooArgList 'vars' defines the dimensions of the histogram.
191/// The ranges and number of bins are taken from the input histogram and must be the same in all histograms
192
194 std::map<string,RooDataHist*> dhistMap, double wgt) :
195 RooAbsData(name,title,RooArgSet(vars,&indexCat))
196{
197 // Initialize datastore
199
200 importDHistSet(vars, indexCat, dhistMap, wgt) ;
201
204}
205
206
207
208////////////////////////////////////////////////////////////////////////////////
209/// Constructor of a data hist from an TH1,TH2 or TH3
210/// The RooArgSet 'vars' defines the dimensions of the histogram. The ranges
211/// and number of bins are taken from the input histogram, and the corresponding
212/// values are set accordingly on the arguments in 'vars'
213
214RooDataHist::RooDataHist(RooStringView name, RooStringView title, const RooArgList& vars, const TH1* hist, double wgt) :
215 RooAbsData(name,title,vars)
216{
217 // Initialize datastore
219
220 // Check consistency in number of dimensions
221 if (vars.getSize() != hist->GetDimension()) {
222 std::stringstream errorMsgStream;
223 errorMsgStream << "RooDataHist::ctor(" << GetName() << ") ERROR: dimension of input histogram must match "
224 << "number of dimension variables";
225 const std::string errorMsg = errorMsgStream.str();
226 coutE(InputArguments) << errorMsg << std::endl;
227 throw std::invalid_argument(errorMsg);
228 }
229
230 importTH1(vars,*hist,wgt, false) ;
231
234}
235
236
237
238////////////////////////////////////////////////////////////////////////////////
239/// Constructor of a binned dataset from a RooArgSet defining the dimensions
240/// of the data space. The range and number of bins in each dimensions are taken
241/// from getMin() getMax(),getBins() of each RooAbsArg representing that
242/// dimension.
243///
244/// <table>
245/// <tr><th> Optional Argument <th> Effect
246/// <tr><td> Import(TH1&, bool impDens) <td> Import contents of the given TH1/2/3 into this binned dataset. The
247/// ranges and binning of the binned dataset are automatically adjusted to
248/// match those of the imported histogram.
249///
250/// Please note: for TH1& with unequal binning _only_,
251/// you should decide if you want to import the absolute bin content,
252/// or the bin content expressed as density. The latter is default and will
253/// result in the same histogram as the original TH1. For certain types of
254/// bin contents (containing efficiencies, asymmetries, or ratio is general)
255/// you should import the absolute value and set impDens to false
256///
257///
258/// <tr><td> Weight(double) <td> Apply given weight factor when importing histograms
259///
260/// <tr><td> Index(RooCategory&) <td> Prepare import of multiple TH1/1/2/3 into a N+1 dimensional RooDataHist
261/// where the extra discrete dimension labels the source of the imported histogram
262/// If the index category defines states for which no histogram is be imported
263/// the corresponding bins will be left empty.
264///
265/// <tr><td> Import(const char*, TH1&) <td> Import a THx to be associated with the given state name of the index category
266/// specified in Index(). If the given state name is not yet defined in the index
267/// category it will be added on the fly. The import command can be specified
268/// multiple times.
269/// <tr><td> Import(map<string,TH1*>&) <td> As above, but allows specification of many imports in a single operation
270/// <tr><td> `GlobalObservables(const RooArgSet&)` <td> Define the set of global observables to be stored in this RooDataHist.
271/// A snapshot of the passed RooArgSet is stored, meaning the values wont't change unexpectedly.
272/// </table>
273///
274
275RooDataHist::RooDataHist(RooStringView name, RooStringView title, const RooArgList& vars, const RooCmdArg& arg1, const RooCmdArg& arg2, const RooCmdArg& arg3,
276 const RooCmdArg& arg4,const RooCmdArg& arg5,const RooCmdArg& arg6,const RooCmdArg& arg7,const RooCmdArg& arg8) :
277 RooAbsData(name,title,RooArgSet(vars,(RooAbsArg*)RooCmdConfig::decodeObjOnTheFly("RooDataHist::RooDataHist", "IndexCat",0,nullptr,arg1,arg2,arg3,arg4,arg5,arg6,arg7,arg8)))
278{
279 // Initialize datastore
281
282 // Define configuration for this method
283 RooCmdConfig pc("RooDataHist::ctor(" + std::string(GetName()) + ")");
284 pc.defineObject("impHist","ImportHisto",0) ;
285 pc.defineInt("impDens","ImportHisto",0) ;
286 pc.defineObject("indexCat","IndexCat",0) ;
287 pc.defineObject("impSliceData","ImportDataSlice",0,nullptr,true) ; // array
288 pc.defineString("impSliceState","ImportDataSlice",0,"",true) ; // array
289 pc.defineDouble("weight","Weight",0,1) ;
290 pc.defineObject("dummy1","ImportDataSliceMany",0) ;
291 pc.defineSet("glObs","GlobalObservables",0,nullptr) ;
292 pc.defineMutex("ImportHisto","ImportDataSlice");
293 pc.defineDependency("ImportDataSlice","IndexCat") ;
294
296 l.Add((TObject*)&arg1) ; l.Add((TObject*)&arg2) ;
297 l.Add((TObject*)&arg3) ; l.Add((TObject*)&arg4) ;
298 l.Add((TObject*)&arg5) ; l.Add((TObject*)&arg6) ;
299 l.Add((TObject*)&arg7) ; l.Add((TObject*)&arg8) ;
300
301 // Process & check varargs
302 pc.process(l) ;
303 if (!pc.ok(true)) {
304 throw std::invalid_argument("Invalid command arguments passed to RooDataHist constructor!");
305 }
306
307 if(pc.getSet("glObs")) setGlobalObservables(*pc.getSet("glObs"));
308
309 TH1* impHist = static_cast<TH1*>(pc.getObject("impHist")) ;
310 bool impDens = pc.getInt("impDens") ;
311 double initWgt = pc.getDouble("weight") ;
312 RooCategory* indexCat = static_cast<RooCategory*>(pc.getObject("indexCat")) ;
313 const char* impSliceNames = pc.getString("impSliceState","",true) ;
314 const RooLinkedList& impSliceHistos = pc.getObjectList("impSliceData") ;
315
316
317 if (impHist) {
318
319 // Initialize importing contents from TH1
320 importTH1(vars,*impHist,initWgt, impDens) ;
321
322 } else if (indexCat) {
323
324
325 // Initialize importing mapped set of RooDataHists and TH1s
326 std::map<std::string,RooDataHist*> dmap ;
327 std::map<std::string,TH1*> hmap ;
328 auto hiter = impSliceHistos.begin() ;
329 for (const auto& token : ROOT::Split(impSliceNames, ",", /*skipEmpty=*/true)) {
330 if(auto dHist = dynamic_cast<RooDataHist*>(*hiter)) {
331 dmap[token] = dHist;
332 }
333 if(auto hHist = dynamic_cast<TH1*>(*hiter)) {
334 hmap[token] = hHist;
335 }
336 ++hiter;
337 }
338 if(!dmap.empty() && !hmap.empty()) {
339 std::stringstream errorMsgStream;
340 errorMsgStream << "RooDataHist::ctor(" << GetName() << ") ERROR: you can't import mix of TH1 and RooDataHist";
341 const std::string errorMsg = errorMsgStream.str();
342 coutE(InputArguments) << errorMsg << std::endl;
343 throw std::invalid_argument(errorMsg);
344 }
345 if (!dmap.empty()) {
346 importDHistSet(vars,*indexCat,dmap,initWgt);
347 }
348 if (!hmap.empty()) {
349 importTH1Set(vars,*indexCat,hmap,initWgt,false);
350 }
351
352
353 } else {
354
355 // Initialize empty
356 initialize() ;
357 appendToDir(this,true) ;
358
359 }
360
363
364}
365
366
367
368
369////////////////////////////////////////////////////////////////////////////////
370/// Import data from given TH1/2/3 into this RooDataHist
371
372void RooDataHist::importTH1(const RooArgList& vars, const TH1& histo, double wgt, bool doDensityCorrection)
373{
374 // Adjust binning of internal observables to match that of input THx
375 Int_t offset[3]{0, 0, 0};
376 adjustBinning(vars, histo, offset) ;
377
378 // Initialize internal data structure
379 initialize() ;
380 appendToDir(this,true) ;
381
382 // Define x,y,z as 1st, 2nd and 3rd observable
383 RooRealVar* xvar = (RooRealVar*) _vars.find(vars.at(0)->GetName()) ;
384 RooRealVar* yvar = (RooRealVar*) (vars.at(1) ? _vars.find(vars.at(1)->GetName()) : nullptr ) ;
385 RooRealVar* zvar = (RooRealVar*) (vars.at(2) ? _vars.find(vars.at(2)->GetName()) : nullptr ) ;
386
387 // Transfer contents
388 Int_t xmin(0),ymin(0),zmin(0) ;
389 RooArgSet vset(*xvar) ;
390 xmin = offset[0] ;
391 if (yvar) {
392 vset.add(*yvar) ;
393 ymin = offset[1] ;
394 }
395 if (zvar) {
396 vset.add(*zvar) ;
397 zmin = offset[2] ;
398 }
399
400 Int_t ix(0),iy(0),iz(0) ;
401 for (ix=0 ; ix < xvar->getBins() ; ix++) {
402 xvar->setBin(ix) ;
403 if (yvar) {
404 for (iy=0 ; iy < yvar->getBins() ; iy++) {
405 yvar->setBin(iy) ;
406 if (zvar) {
407 for (iz=0 ; iz < zvar->getBins() ; iz++) {
408 zvar->setBin(iz) ;
409 double bv = doDensityCorrection ? binVolume(vset) : 1;
410 add(vset,bv*histo.GetBinContent(ix+1+xmin,iy+1+ymin,iz+1+zmin)*wgt,bv*TMath::Power(histo.GetBinError(ix+1+xmin,iy+1+ymin,iz+1+zmin)*wgt,2)) ;
411 }
412 } else {
413 double bv = doDensityCorrection ? binVolume(vset) : 1;
414 add(vset,bv*histo.GetBinContent(ix+1+xmin,iy+1+ymin)*wgt,bv*TMath::Power(histo.GetBinError(ix+1+xmin,iy+1+ymin)*wgt,2)) ;
415 }
416 }
417 } else {
418 double bv = doDensityCorrection ? binVolume(vset) : 1 ;
419 add(vset,bv*histo.GetBinContent(ix+1+xmin)*wgt,bv*TMath::Power(histo.GetBinError(ix+1+xmin)*wgt,2)) ;
420 }
421 }
422
423}
424
425namespace {
426bool checkConsistentAxes(const TH1* first, const TH1* second) {
427 return first->GetDimension() == second->GetDimension()
428 && first->GetNbinsX() == second->GetNbinsX()
429 && first->GetNbinsY() == second->GetNbinsY()
430 && first->GetNbinsZ() == second->GetNbinsZ()
431 && first->GetXaxis()->GetXmin() == second->GetXaxis()->GetXmin()
432 && first->GetXaxis()->GetXmax() == second->GetXaxis()->GetXmax()
433 && (first->GetNbinsY() == 1 || (first->GetYaxis()->GetXmin() == second->GetYaxis()->GetXmin()
434 && first->GetYaxis()->GetXmax() == second->GetYaxis()->GetXmax() ) )
435 && (first->GetNbinsZ() == 1 || (first->GetZaxis()->GetXmin() == second->GetZaxis()->GetXmin()
436 && first->GetZaxis()->GetXmax() == second->GetZaxis()->GetXmax() ) );
437}
438
439bool hasConsistentLayoutAndBinning(RooDataHist const& h1, RooDataHist const& h2) {
440
441 // Relative tolerance for bin boundary comparison
442 constexpr double tolerance = 1e-6;
443
444 auto const& vars1 = *h1.get();
445 auto const& vars2 = *h2.get();
446
447 // Check if number of variables and names is consistent
448 if(!vars1.hasSameLayout(vars2)) {
449 return false;
450 }
451
452 for(std::size_t iVar = 0; iVar < vars1.size(); ++iVar) {
453 auto * var1 = dynamic_cast<RooRealVar*>(vars1[iVar]);
454 auto * var2 = dynamic_cast<RooRealVar*>(vars2[iVar]);
455
456 // Check if variables are consistently real-valued
457 if((!var1 && var2) || (var1 && !var2)) return false;
458
459 // Not a real-valued variable
460 if(!var1) continue;
461
462 // Now check the binning
463 auto const& bng1 = var1->getBinning();
464 auto const& bng2 = var2->getBinning();
465
466 // Compare bin numbers
467 if(bng1.numBins() != bng2.numBins()) return false;
468
469 std::size_t nBins = bng1.numBins();
470
471 // Compare bin boundaries
472 for(std::size_t iBin = 0; iBin < nBins; ++iBin) {
473 double v1 = bng1.binLow(iBin);
474 double v2 = bng2.binLow(iBin);
475 if(std::abs((v1 - v2) / v1) > tolerance) return false;
476 }
477 double v1 = bng1.binHigh(nBins - 1);
478 double v2 = bng2.binHigh(nBins - 1);
479 if(std::abs((v1 - v2) / v1) > tolerance) return false;
480 }
481 return true;
482}
483}
484
485
486////////////////////////////////////////////////////////////////////////////////
487/// Import data from given set of TH1/2/3 into this RooDataHist. The category indexCat labels the sources
488/// in the constructed RooDataHist. The stl map provides the mapping between the indexCat state labels
489/// and the import source
490
491void RooDataHist::importTH1Set(const RooArgList& vars, RooCategory& indexCat, std::map<string,TH1*> hmap, double wgt, bool doDensityCorrection)
492{
493 RooCategory* icat = (RooCategory*) _vars.find(indexCat.GetName()) ;
494
495 TH1* histo(nullptr) ;
496 bool init(false) ;
497 for (const auto& hiter : hmap) {
498 // Store pointer to first histogram from which binning specification will be taken
499 if (!histo) {
500 histo = hiter.second;
501 } else {
502 if (!checkConsistentAxes(histo, hiter.second)) {
503 coutE(InputArguments) << "Axes of histogram " << hiter.second->GetName() << " are not consistent with first processed "
504 << "histogram " << histo->GetName() << std::endl;
505 throw std::invalid_argument("Axes of inputs for RooDataHist are inconsistent");
506 }
507 }
508 // Define state labels in index category (both in provided indexCat and in internal copy in dataset)
509 if (!indexCat.hasLabel(hiter.first)) {
510 indexCat.defineType(hiter.first) ;
511 coutI(InputArguments) << "RooDataHist::importTH1Set(" << GetName() << ") defining state \"" << hiter.first << "\" in index category " << indexCat.GetName() << endl ;
512 }
513 if (!icat->hasLabel(hiter.first)) {
514 icat->defineType(hiter.first) ;
515 }
516 }
517
518 // Check consistency in number of dimensions
519 if (histo && (vars.getSize() != histo->GetDimension())) {
520 coutE(InputArguments) << "RooDataHist::importTH1Set(" << GetName() << "): dimension of input histogram must match "
521 << "number of continuous variables" << endl ;
522 throw std::invalid_argument("Inputs histograms for RooDataHist are not compatible with dimensions of variables.");
523 }
524
525 // Copy bins and ranges from THx to dimension observables
526 Int_t offset[3] ;
527 adjustBinning(vars,*histo,offset) ;
528
529 // Initialize internal data structure
530 if (!init) {
531 initialize() ;
532 appendToDir(this,true) ;
533 init = true ;
534 }
535
536 // Define x,y,z as 1st, 2nd and 3rd observable
537 RooRealVar* xvar = (RooRealVar*) _vars.find(vars.at(0)->GetName()) ;
538 RooRealVar* yvar = (RooRealVar*) (vars.at(1) ? _vars.find(vars.at(1)->GetName()) : nullptr ) ;
539 RooRealVar* zvar = (RooRealVar*) (vars.at(2) ? _vars.find(vars.at(2)->GetName()) : nullptr ) ;
540
541 // Transfer contents
542 Int_t xmin(0),ymin(0),zmin(0) ;
543 RooArgSet vset(*xvar) ;
544 double volume = xvar->getMax()-xvar->getMin() ;
545 xmin = offset[0] ;
546 if (yvar) {
547 vset.add(*yvar) ;
548 ymin = offset[1] ;
549 volume *= (yvar->getMax()-yvar->getMin()) ;
550 }
551 if (zvar) {
552 vset.add(*zvar) ;
553 zmin = offset[2] ;
554 volume *= (zvar->getMax()-zvar->getMin()) ;
555 }
556 double avgBV = volume / numEntries() ;
557
558 Int_t ic(0),ix(0),iy(0),iz(0) ;
559 for (ic=0 ; ic < icat->numBins(nullptr) ; ic++) {
560 icat->setBin(ic) ;
561 histo = hmap[icat->getCurrentLabel()] ;
562 for (ix=0 ; ix < xvar->getBins() ; ix++) {
563 xvar->setBin(ix) ;
564 if (yvar) {
565 for (iy=0 ; iy < yvar->getBins() ; iy++) {
566 yvar->setBin(iy) ;
567 if (zvar) {
568 for (iz=0 ; iz < zvar->getBins() ; iz++) {
569 zvar->setBin(iz) ;
570 double bv = doDensityCorrection ? binVolume(vset)/avgBV : 1;
571 add(vset,bv*histo->GetBinContent(ix+1+xmin,iy+1+ymin,iz+1+zmin)*wgt,bv*TMath::Power(histo->GetBinError(ix+1+xmin,iy+1+ymin,iz+1+zmin)*wgt,2)) ;
572 }
573 } else {
574 double bv = doDensityCorrection ? binVolume(vset)/avgBV : 1;
575 add(vset,bv*histo->GetBinContent(ix+1+xmin,iy+1+ymin)*wgt,bv*TMath::Power(histo->GetBinError(ix+1+xmin,iy+1+ymin)*wgt,2)) ;
576 }
577 }
578 } else {
579 double bv = doDensityCorrection ? binVolume(vset)/avgBV : 1;
580 add(vset,bv*histo->GetBinContent(ix+1+xmin)*wgt,bv*TMath::Power(histo->GetBinError(ix+1+xmin)*wgt,2)) ;
581 }
582 }
583 }
584
585}
586
587
588
589////////////////////////////////////////////////////////////////////////////////
590/// Import data from given set of TH1/2/3 into this RooDataHist. The category indexCat labels the sources
591/// in the constructed RooDataHist. The stl map provides the mapping between the indexCat state labels
592/// and the import source
593
594void RooDataHist::importDHistSet(const RooArgList& /*vars*/, RooCategory& indexCat, std::map<std::string,RooDataHist*> dmap, double initWgt)
595{
596 auto* icat = static_cast<RooCategory*>(_vars.find(indexCat.GetName()));
597
598 RooDataHist* dhistForBinning = nullptr;
599
600 for (const auto& diter : dmap) {
601
602 std::string const& label = diter.first;
603 RooDataHist* dhist = diter.second ;
604
605 if(!dhistForBinning) {
606 dhistForBinning = dhist;
607 }
608 else {
609 if(!hasConsistentLayoutAndBinning(*dhistForBinning, *dhist)) {
610 coutE(InputArguments) << "Layout or binning of histogram " << dhist->GetName() << " is not consistent with first processed "
611 << "histogram " << dhistForBinning->GetName() << std::endl;
612 throw std::invalid_argument("Layout or binning of inputs for RooDataHist is inconsistent");
613 }
614 }
615
616 // Define state labels in index category (both in provided indexCat and in internal copy in dataset)
617 if (!indexCat.hasLabel(label)) {
618 indexCat.defineType(label) ;
619 coutI(InputArguments) << "RooDataHist::importDHistSet(" << GetName() << ") defining state \"" << label << "\" in index category " << indexCat.GetName() << endl ;
620 }
621 if (!icat->hasLabel(label)) {
622 icat->defineType(label) ;
623 }
624 }
625
626 // adjust the binning of the created histogram
627 for(auto * theirVar : dynamic_range_cast<RooRealVar*>(dhistForBinning->_vars)) {
628 auto * ourVar = dynamic_cast<RooRealVar*>(_vars.find(theirVar->GetName()));
629 if(!theirVar || !ourVar) continue;
630 ourVar->setBinning(theirVar->getBinning());
631 }
632
633 initialize() ;
634 appendToDir(this,true) ;
635
636
637 for (const auto& diter : dmap) {
638 std::string const& label = diter.first;
639 RooDataHist* dhist = diter.second ;
640
641 icat->setLabel(label.c_str()) ;
642
643 // Transfer contents
644 for (Int_t i=0 ; i<dhist->numEntries() ; i++) {
645 _vars.assign(*dhist->get(i)) ;
646 add(_vars,dhist->weight()*initWgt, pow(dhist->weightError(SumW2),2) ) ;
647 }
648
649 }
650}
651
652////////////////////////////////////////////////////////////////////////////////
653/// Helper doing the actual work of adjustBinning().
654
655void RooDataHist::_adjustBinning(RooRealVar &theirVar, const TAxis &axis,
656 RooRealVar *ourVar, Int_t *offset)
657{
658 const std::string ourVarName(ourVar->GetName() ? ourVar->GetName() : ""), ownName(GetName() ? GetName() : "");
659 // RooRealVar is derived from RooAbsRealLValue which is itself
660 // derived from RooAbsReal and a virtual class RooAbsLValue
661 // supplying setter functions, check if ourVar is indeed derived
662 // as real
663 if (!dynamic_cast<RooAbsReal *>(ourVar)) {
664 coutE(InputArguments) << "RooDataHist::adjustBinning(" << ownName << ") ERROR: dimension " << ourVarName
665 << " must be real\n";
666 throw std::logic_error("Incorrect type object (" + ourVarName +
667 ") passed as argument to RooDataHist::_adjustBinning. Please report this issue.");
668 }
669
670 const double xlo = theirVar.getMin();
671 const double xhi = theirVar.getMax();
672
673 if (axis.GetXbins()->GetArray()) {
674 RooBinning xbins(axis.GetNbins(), axis.GetXbins()->GetArray());
675
676 const double tolerance = 1e-6 * xbins.averageBinWidth();
677
678 // Adjust xlo/xhi to nearest boundary
679 const double xloAdj = xbins.binLow(xbins.binNumber(xlo + tolerance));
680 const double xhiAdj = xbins.binHigh(xbins.binNumber(xhi - tolerance));
681 xbins.setRange(xloAdj, xhiAdj);
682
683 theirVar.setBinning(xbins);
684
685 if (true || std::abs(xloAdj - xlo) > tolerance || std::abs(xhiAdj - xhi) > tolerance) {
686 coutI(DataHandling) << "RooDataHist::adjustBinning(" << ownName << "): fit range of variable " << ourVarName
687 << " expanded to nearest bin boundaries: [" << xlo << "," << xhi << "] --> [" << xloAdj
688 << "," << xhiAdj << "]"
689 << "\n";
690 }
691
692 ourVar->setBinning(xbins);
693
694 if (offset) {
695 *offset = xbins.rawBinNumber(xloAdj + tolerance);
696 }
697 } else {
698 RooBinning xbins(axis.GetXmin(), axis.GetXmax());
699 xbins.addUniform(axis.GetNbins(), axis.GetXmin(), axis.GetXmax());
700
701 const double tolerance = 1e-6 * xbins.averageBinWidth();
702
703 // Adjust xlo/xhi to nearest boundary
704 const double xloAdj = xbins.binLow(xbins.binNumber(xlo + tolerance));
705 const double xhiAdj = xbins.binHigh(xbins.binNumber(xhi - tolerance));
706 xbins.setRange(xloAdj, xhiAdj);
707 theirVar.setRange(xloAdj, xhiAdj);
708
709 if (std::abs(xloAdj - xlo) > tolerance || std::abs(xhiAdj - xhi) > tolerance) {
710 coutI(DataHandling) << "RooDataHist::adjustBinning(" << ownName << "): fit range of variable " << ourVarName
711 << " expanded to nearest bin boundaries: [" << xlo << "," << xhi << "] --> [" << xloAdj
712 << "," << xhiAdj << "]"
713 << "\n";
714 }
715
716 RooUniformBinning xbins2(xloAdj, xhiAdj, xbins.numBins());
717 ourVar->setBinning(xbins2);
718
719 if (offset) {
720 *offset = xbins.rawBinNumber(xloAdj + tolerance);
721 }
722 }
723}
724
725////////////////////////////////////////////////////////////////////////////////
726/// Adjust binning specification on first and optionally second and third
727/// observable to binning in given reference TH1. Used by constructors
728/// that import data from an external TH1.
729/// Both the variables in vars and in this RooDataHist are adjusted.
730/// @param vars List with variables that are supposed to have their binning adjusted.
731/// @param href Reference histogram that dictates the binning
732/// @param offset If not nullptr, a possible bin count offset for the axes x,y,z is saved here as Int_t[3]
733
734void RooDataHist::adjustBinning(const RooArgList& vars, const TH1& href, Int_t* offset)
735{
736 auto xvar = static_cast<RooRealVar*>( _vars.find(*vars.at(0)) );
737 _adjustBinning(*static_cast<RooRealVar*>(vars.at(0)), *href.GetXaxis(), xvar, offset ? &offset[0] : nullptr);
738
739 if (vars.at(1)) {
740 auto yvar = static_cast<RooRealVar*>(_vars.find(*vars.at(1)));
741 if (yvar)
742 _adjustBinning(*static_cast<RooRealVar*>(vars.at(1)), *href.GetYaxis(), yvar, offset ? &offset[1] : nullptr);
743 }
744
745 if (vars.at(2)) {
746 auto zvar = static_cast<RooRealVar*>(_vars.find(*vars.at(2)));
747 if (zvar)
748 _adjustBinning(*static_cast<RooRealVar*>(vars.at(2)), *href.GetZaxis(), zvar, offset ? &offset[2] : nullptr);
749 }
750
751}
752
753
754namespace {
755/// Clone external weight arrays, unless the external array is nullptr.
756void cloneArray(double*& ours, const double* theirs, std::size_t n) {
757 if (ours) delete[] ours;
758 ours = nullptr;
759 if (!theirs) return;
760 ours = new double[n];
761 std::copy(theirs, theirs+n, ours);
762}
763
764/// Allocate and initialise an array with desired size and values.
765void initArray(double*& arr, std::size_t n, double val) {
766 if (arr) delete[] arr;
767 arr = nullptr;
768 if (n == 0) return;
769 arr = new double[n];
770 std::fill(arr, arr+n, val);
771}
772}
773
774
775////////////////////////////////////////////////////////////////////////////////
776/// Initialization procedure: allocate weights array, calculate
777/// multipliers needed for N-space to 1-dim array jump table,
778/// and fill the internal tree with all bin center coordinates
779
780void RooDataHist::initialize(const char* binningName, bool fillTree)
781{
782 _lvvars.clear();
783 _lvbins.clear();
784
785 // Fill array of LValue pointers to variables
786 for (unsigned int i = 0; i < _vars.size(); ++i) {
787 if (binningName) {
788 RooRealVar* rrv = dynamic_cast<RooRealVar*>(_vars[i]);
789 if (rrv) {
790 rrv->setBinning(rrv->getBinning(binningName));
791 }
792 }
793
794 auto lvarg = dynamic_cast<RooAbsLValue*>(_vars[i]);
795 assert(lvarg);
796 _lvvars.push_back(lvarg);
797
798 const RooAbsBinning* binning = lvarg->getBinningPtr(nullptr);
799 _lvbins.emplace_back(binning ? binning->clone() : nullptr);
800 }
801
802
803 // Allocate coefficients array
804 _idxMult.resize(_vars.size()) ;
805
806 _arrSize = 1 ;
807 unsigned int n = 0u;
808 for (const auto var : _vars) {
809 auto arg = dynamic_cast<const RooAbsLValue*>(var);
810 assert(arg);
811
812 // Calculate sub-index multipliers for master index
813 for (unsigned int i = 0u; i<n; i++) {
814 _idxMult[i] *= arg->numBins() ;
815 }
816 _idxMult[n++] = 1 ;
817
818 // Calculate dimension of weight array
819 _arrSize *= arg->numBins() ;
820 }
821
822 // Allocate and initialize weight array if necessary
823 if (!_wgt) {
824 initArray(_wgt, _arrSize, 0.);
825 delete[] _errLo; _errLo = nullptr;
826 delete[] _errHi; _errHi = nullptr;
827 delete[] _sumw2; _sumw2 = nullptr;
828 initArray(_binv, _arrSize, 0.);
829
830 // Refill array pointers in data store when reading
831 // from Streamer
832 if (!fillTree) {
834 }
835 }
836
837 if (!fillTree) return ;
838
839 // Fill TTree with bin center coordinates
840 // Calculate plot bins of components from master index
841
842 for (Int_t ibin=0 ; ibin < _arrSize ; ibin++) {
843 Int_t j(0), idx(0), tmp(ibin) ;
844 double theBinVolume(1) ;
845 for (auto arg2 : _lvvars) {
846 idx = tmp / _idxMult[j] ;
847 tmp -= idx*_idxMult[j++] ;
848 arg2->setBin(idx) ;
849 theBinVolume *= arg2->getBinWidth(idx) ;
850 }
851 _binv[ibin] = theBinVolume ;
852
853 fill() ;
854 }
855
856
857}
858
859
860////////////////////////////////////////////////////////////////////////////////
861
863{
864 if (!_binbounds.empty()) return;
865 for (auto& it : _lvbins) {
866 _binbounds.push_back(std::vector<double>());
867 if (it) {
868 std::vector<double>& bounds = _binbounds.back();
869 bounds.reserve(2 * it->numBins());
870 for (Int_t i = 0; i < it->numBins(); ++i) {
871 bounds.push_back(it->binLow(i));
872 bounds.push_back(it->binHigh(i));
873 }
874 }
875 }
876}
877
878
879////////////////////////////////////////////////////////////////////////////////
880/// Copy constructor
881
882RooDataHist::RooDataHist(const RooDataHist& other, const char* newname) :
883 RooAbsData(other,newname), RooDirItem(), _arrSize(other._arrSize), _idxMult(other._idxMult), _pbinvCache(other._pbinvCache)
884{
885 // Allocate and initialize weight array
886 assert(_arrSize == other._arrSize);
887 cloneArray(_wgt, other._wgt, other._arrSize);
888 cloneArray(_errLo, other._errLo, other._arrSize);
889 cloneArray(_errHi, other._errHi, other._arrSize);
890 cloneArray(_binv, other._binv, other._arrSize);
891 cloneArray(_sumw2, other._sumw2, other._arrSize);
892
893 // Fill array of LValue pointers to variables
894 for (const auto rvarg : _vars) {
895 auto lvarg = dynamic_cast<RooAbsLValue*>(rvarg);
896 assert(lvarg);
897 _lvvars.push_back(lvarg);
898 const RooAbsBinning* binning = lvarg->getBinningPtr(nullptr);
899 _lvbins.emplace_back(binning ? binning->clone() : nullptr) ;
900 }
901
903
904 appendToDir(this,true) ;
905}
906
907
908////////////////////////////////////////////////////////////////////////////////
909/// Implementation of RooAbsData virtual method that drives the RooAbsData::reduce() methods
910
911std::unique_ptr<RooAbsData> RooDataHist::reduceEng(const RooArgSet& varSubset, const RooFormulaVar* cutVar, const char* cutRange,
912 std::size_t nStart, std::size_t nStop)
913{
914 checkInit() ;
915 RooArgSet myVarSubset;
916 _vars.selectCommon(varSubset, myVarSubset);
917 auto rdh = std::make_unique<RooDataHist>(GetName(), GetTitle(), myVarSubset);
918
919 RooFormulaVar* cloneVar = nullptr;
920 std::unique_ptr<RooArgSet> tmp;
921 if (cutVar) {
922 tmp = std::make_unique<RooArgSet>();
923 // Deep clone cutVar and attach clone to this dataset
924 if (RooArgSet(*cutVar).snapshot(*tmp)) {
925 coutE(DataHandling) << "RooDataHist::reduceEng(" << GetName() << ") Couldn't deep-clone cut variable, abort," << endl ;
926 return nullptr;
927 }
928 cloneVar = static_cast<RooFormulaVar*>(tmp->find(*cutVar));
929 cloneVar->attachDataSet(*this) ;
930 }
931
932 double lo,hi ;
933 const std::size_t nevt = nStop < static_cast<std::size_t>(numEntries()) ? nStop : static_cast<std::size_t>(numEntries());
934 for (auto i=nStart; i<nevt ; i++) {
935 const RooArgSet* row = get(i) ;
936
937 bool doSelect(true) ;
938 if (cutRange) {
939 for (const auto arg : *row) {
940 if (!arg->inRange(cutRange)) {
941 doSelect = false ;
942 break ;
943 }
944 }
945 }
946 if (!doSelect) continue ;
947
948 if (!cloneVar || cloneVar->getVal()) {
949 weightError(lo,hi,SumW2) ;
950 rdh->add(*row,weight(),lo*lo) ;
951 }
952 }
953
954 return rdh ;
955}
956
957
958
959////////////////////////////////////////////////////////////////////////////////
960/// Destructor
961
963{
964 delete[] _wgt;
965 delete[] _errLo;
966 delete[] _errHi;
967 delete[] _sumw2;
968 delete[] _binv;
969
970 removeFromDir(this) ;
972}
973
974
975
976
977////////////////////////////////////////////////////////////////////////////////
978/// Calculate bin number of the given coordinates. If only a subset of the internal
979/// coordinates are passed, the missing coordinates are taken at their current value.
980/// \param[in] coord Variables that are representing the coordinates.
981/// \param[in] fast If the variables in `coord` and the ones of the data hist have the
982/// same size and layout, `fast` can be set to skip checking that all variables are
983/// present in `coord`.
984Int_t RooDataHist::getIndex(const RooAbsCollection& coord, bool fast) const {
985 checkInit() ;
986 return calcTreeIndex(coord, fast);
987}
988
990 bool correctForBinSize) const
991{
992 std::string weightName = ctx.makeValidVarName(klass->GetName()) + "_HistWeights";
993
994 std::string arrayDecl = "double " + weightName + "[" + std::to_string(_arrSize) + "] = {";
995 for (Int_t i = 0; i < _arrSize; i++) {
996 arrayDecl += " " + std::to_string(_wgt[i]) + (correctForBinSize ? " / " + std::to_string(_binv[i]) : "") + ",";
997 }
998 arrayDecl.back() = ' ';
999 arrayDecl += "};\n";
1000
1001 ctx.addToCodeBody(arrayDecl, true);
1002
1003 return weightName;
1004}
1005
1007 const RooAbsCollection &coords, bool reverse) const
1008{
1009 assert(coords.size() == _vars.size());
1010
1011 std::string idxName = "idx_" + ctx.getTmpVarName();
1012 ctx.addToCodeBody(klass, "unsigned int " + idxName + " = 0;\n");
1013 int idxMult = 1;
1014
1015 for (std::size_t i = 0; i < _vars.size(); ++i) {
1016
1017 std::size_t iVar = reverse ? _vars.size() - 1 - i : i;
1018 const RooAbsArg *internalVar = _vars[iVar];
1019 const RooAbsArg *theVar = coords[iVar];
1020
1021 const RooAbsBinning *binning = _lvbins[iVar].get();
1022 if (!binning) {
1023 coutE(InputArguments) << "RooHistPdf::weight(" << GetName()
1024 << ") ERROR: Code Squashing currently does not support category values." << std::endl;
1025 return "";
1026 } else if (!dynamic_cast<RooUniformBinning const *>(binning)) {
1027 coutE(InputArguments) << "RooHistPdf::weight(" << GetName()
1028 << ") ERROR: Code Squashing currently only supports uniformly binned cases."
1029 << std::endl;
1030 return "";
1031 }
1032
1033 std::string const &bin = ctx.buildCall("RooFit::Detail::EvaluateFuncs::getUniformBinning", binning->lowBound(),
1034 binning->highBound(), *theVar, binning->numBins());
1035 ctx.addToCodeBody(klass, idxName + " += " + std::to_string(idxMult) + " * " + bin + ";\n");
1036
1037 // Use RooAbsLValue here because it also generalized to categories, which
1038 // is useful in the future. dynamic_cast because it's a cross-cast.
1039 idxMult *= dynamic_cast<RooAbsLValue const *>(internalVar)->numBins();
1040 }
1041
1042 return idxName;
1043}
1044
1045////////////////////////////////////////////////////////////////////////////////
1046/// Calculate the bin index corresponding to the coordinates passed as argument.
1047/// \param[in] coords Coordinates. If `fast == false`, these can be partial.
1048/// \param[in] fast Promise that the coordinates in `coords` have the same order
1049/// as the internal coordinates. In this case, values are looked up only by index.
1050std::size_t RooDataHist::calcTreeIndex(const RooAbsCollection& coords, bool fast) const
1051{
1052 // With fast, caller promises that layout of `coords` is identical to our internal `vars`.
1053 // Previously, this was verified with an assert in debug mode like this:
1054 //
1055 // assert(!fast || coords.hasSameLayout(_vars));
1056 //
1057 // However, there are usecases where the externally provided `coords` have
1058 // different names than the internal variables, even though they correspond
1059 // to each other. For example, if the observables in the computation graph
1060 // are renamed with `redirectServers`. Hence, we can't do a meaningful assert
1061 // here.
1062
1063 if (&_vars == &coords)
1064 fast = true;
1065
1066 std::size_t masterIdx = 0;
1067
1068 for (unsigned int i=0; i < _vars.size(); ++i) {
1069 const RooAbsArg* internalVar = _vars[i];
1070 const RooAbsBinning* binning = _lvbins[i].get();
1071
1072 // Find the variable that we need values from.
1073 // That's either the variable directly from the external coordinates
1074 // or we find the external one that has the same name as "internalVar".
1075 const RooAbsArg* theVar = fast ? coords[i] : coords.find(*internalVar);
1076 if (!theVar) {
1077 // Variable is not in external coordinates. Use current internal value.
1078 theVar = internalVar;
1079 }
1080 // If fast is on, users promise that the sets have the same layout:
1081 //
1082 // assert(!fast || strcmp(internalVar->GetName(), theVar->GetName()) == 0);
1083 //
1084 // This assert is commented out for the same reasons that applied to the
1085 // other assert explained above.
1086
1087 if (binning) {
1088 assert(dynamic_cast<const RooAbsReal*>(theVar));
1089 const double val = static_cast<const RooAbsReal*>(theVar)->getVal();
1090 masterIdx += _idxMult[i] * binning->binNumber(val);
1091 } else {
1092 // We are a category. No binning.
1093 assert(dynamic_cast<const RooAbsCategoryLValue*>(theVar));
1094 auto cat = static_cast<const RooAbsCategoryLValue*>(theVar);
1095 masterIdx += _idxMult[i] * cat->getBin(static_cast<const char*>(nullptr));
1096 }
1097 }
1098
1099 return masterIdx ;
1100}
1101
1102
1103////////////////////////////////////////////////////////////////////////////////
1104/// Back end function to plotting functionality. Plot RooDataHist on given
1105/// frame in mode specified by plot options 'o'. The main purpose of
1106/// this function is to match the specified binning on 'o' to the
1107/// internal binning of the plot observable in this RooDataHist.
1108/// \see RooAbsData::plotOn() for plotting options.
1110{
1111 checkInit() ;
1112 if (o.bins) return RooAbsData::plotOn(frame,o) ;
1113
1114 if(!frame) {
1115 coutE(InputArguments) << ClassName() << "::" << GetName() << ":plotOn: frame is null" << endl;
1116 return nullptr;
1117 }
1118 auto var= static_cast<RooAbsRealLValue*>(frame->getPlotVar());
1119 if(!var) {
1120 coutE(InputArguments) << ClassName() << "::" << GetName()
1121 << ":plotOn: frame does not specify a plot variable" << endl;
1122 return nullptr;
1123 }
1124
1125 auto dataVar = static_cast<RooRealVar*>(_vars.find(*var));
1126 if (!dataVar) {
1127 coutE(InputArguments) << ClassName() << "::" << GetName()
1128 << ":plotOn: dataset doesn't contain plot frame variable" << endl;
1129 return nullptr;
1130 }
1131
1132 o.bins = &dataVar->getBinning() ;
1133 return RooAbsData::plotOn(frame,o) ;
1134}
1135
1136
1137////////////////////////////////////////////////////////////////////////////////
1138/// A vectorized version of interpolateDim for boundary safe quadratic
1139/// interpolation of one dimensional histograms.
1140///
1141/// \param[out] output An array of interpolated weights corresponding to the
1142/// values in xVals.
1143/// \param[in] xVals An array of event coordinates for which the weights should be
1144/// calculated.
1145/// \param[in] correctForBinSize Enable the inverse bin volume correction factor.
1146/// \param[in] cdfBoundaries Enable the special boundary condition for a cdf:
1147/// Underflow bins are assumed to have weight zero and
1148/// overflow bins have weight one. Otherwise, the
1149/// histogram is mirrored at the boundaries for the
1150/// interpolation.
1151
1152void RooDataHist::interpolateQuadratic(double* output, std::span<const double> xVals,
1153 bool correctForBinSize, bool cdfBoundaries)
1154{
1155 const std::size_t nBins = numEntries();
1156 const std::size_t nEvents = xVals.size();
1157
1158 RooAbsBinning const& binning = *_lvbins[0];
1159 // Reuse the output buffer for bin indices and zero-initialize it
1160 auto binIndices = reinterpret_cast<int*>(output + nEvents) - nEvents;
1161 std::fill(binIndices, binIndices + nEvents, 0);
1162 binning.binNumbers(xVals.data(), binIndices, nEvents);
1163
1164 // Extend coordinates and weights with one extra point before the first bin
1165 // and one extra point after the last bin. This means the original histogram
1166 // bins span elements 1 to nBins in coordsExt and weightsExt
1167 std::vector<double> coordsExt(nBins+3);
1168 double* binCoords = coordsExt.data() + 2;
1169 binCoords[0] = binning.lowBound() + 0.5*_binv[0];
1170 for (std::size_t binIdx = 1; binIdx < nBins ; ++binIdx) {
1171 if (binning.isUniform()) {
1172 double binWidth = _binv[0];
1173 binCoords[binIdx] = binIdx*binWidth + binCoords[0];
1174 }
1175 else {
1176 double binCentDiff = 0.5*_binv[binIdx-1] + 0.5*_binv[binIdx];
1177 binCoords[binIdx] = binCoords[binIdx-1] + binCentDiff;
1178 }
1179 }
1180
1181 std::vector<double> weightsExt(nBins+3);
1182 // Fill weights for bins that are inside histogram boundaries
1183 for (std::size_t binIdx = 0; binIdx < nBins; ++binIdx) {
1184 weightsExt[binIdx+2] = correctForBinSize ? _wgt[binIdx] / _binv[binIdx] : _wgt[binIdx];
1185 }
1186
1187 if (cdfBoundaries) {
1188 coordsExt[0] = - 1e-10 + binning.lowBound();
1189 weightsExt[0] = 0.;
1190
1191 coordsExt[1] = binning.lowBound();
1192 weightsExt[1] = 0.;
1193
1194 coordsExt[nBins+2] = binning.highBound();
1195 weightsExt[nBins+2] = 1.;
1196 }
1197 else {
1198 // Mirror first two bins and last bin
1199 coordsExt[0] = binCoords[1] - 2*_binv[0] - _binv[1];
1200 weightsExt[0] = weightsExt[3];
1201
1202 coordsExt[1] = binCoords[0] - _binv[0];
1203 weightsExt[1] = weightsExt[2];
1204
1205 coordsExt[nBins+2] = binCoords[nBins-1] + _binv[nBins-1];
1206 weightsExt[nBins+2] = weightsExt[nBins+1];
1207 }
1208
1209 // We use the current bin center and two bin centers on the left for
1210 // interpolation if xVal is to the left of the current bin center
1211 for (std::size_t i = 0; i < nEvents ; ++i) {
1212 double xVal = xVals[i];
1213 std::size_t binIdx = binIndices[i] + 2;
1214
1215 // If xVal is to the right of the current bin center, shift all bin
1216 // coordinates one step to the right and use that for the interpolation
1217 if (xVal > coordsExt[binIdx]) {
1218 binIdx += 1;
1219 }
1220
1221 double x1 = coordsExt[binIdx-2];
1222 double y1 = weightsExt[binIdx-2];
1223
1224 double x2 = coordsExt[binIdx-1];
1225 double y2 = weightsExt[binIdx-1];
1226
1227 double x3 = coordsExt[binIdx];
1228 double y3 = weightsExt[binIdx];
1229
1230 // Evaluate a few repeated factors
1231 double quotient = (x3-x1) / (x2-x1);
1232 double x1Sqrd = x1*x1;
1233 double x3Sqrd = x3*x3;
1234 // Solve coefficients in system of three quadratic equations!
1235 double secondCoeff = (y3 - y1 - (y2-y1) * quotient) / (x3Sqrd - x1Sqrd - (x2*x2 - x1Sqrd) * quotient);
1236 double firstCoeff = (y3 - y1 - secondCoeff*(x3Sqrd - x1Sqrd)) / (x3-x1);
1237 double zerothCoeff = y1 - secondCoeff*x1Sqrd - firstCoeff*x1;
1238 // Get the interpolated weight using the equation of a second degree polynomial
1239 output[i] = secondCoeff * xVal * xVal + firstCoeff * xVal + zerothCoeff;
1240 }
1241}
1242
1243
1244////////////////////////////////////////////////////////////////////////////////
1245/// A vectorized version of interpolateDim for boundary safe linear
1246/// interpolation of one dimensional histograms.
1247///
1248/// \param[out] output An array of interpolated weights corresponding to the
1249/// values in xVals.
1250/// \param[in] xVals An array of event coordinates for which the weights should be
1251/// calculated.
1252/// \param[in] correctForBinSize Enable the inverse bin volume correction factor.
1253/// \param[in] cdfBoundaries Enable the special boundary condition for a cdf:
1254/// Underflow bins are assumed to have weight zero and
1255/// overflow bins have weight one. Otherwise, the
1256/// histogram is mirrored at the boundaries for the
1257/// interpolation.
1258
1259void RooDataHist::interpolateLinear(double* output, std::span<const double> xVals,
1260 bool correctForBinSize, bool cdfBoundaries)
1261{
1262 const std::size_t nBins = numEntries();
1263 const std::size_t nEvents = xVals.size();
1264
1265 RooAbsBinning const& binning = *_lvbins[0];
1266 // Reuse the output buffer for bin indices and zero-initialize it
1267 auto binIndices = reinterpret_cast<int*>(output + nEvents) - nEvents;
1268 std::fill(binIndices, binIndices + nEvents, 0);
1269 binning.binNumbers(xVals.data(), binIndices, nEvents);
1270
1271 // Extend coordinates and weights with one extra point before the first bin
1272 // and one extra point after the last bin. This means the original histogram
1273 // bins span elements 1 to nBins in coordsExt and weightsExt
1274 std::vector<double> coordsExt(nBins+2);
1275 double* binCoords = coordsExt.data() + 1;
1276 binCoords[0] = binning.lowBound() + 0.5*_binv[0];
1277 for (std::size_t binIdx = 1; binIdx < nBins ; ++binIdx) {
1278 if (binning.isUniform()) {
1279 double binWidth = _binv[0];
1280 binCoords[binIdx] = binIdx*binWidth + binCoords[0];
1281 }
1282 else {
1283 double binCentDiff = 0.5*_binv[binIdx-1] + 0.5*_binv[binIdx];
1284 binCoords[binIdx] = binCoords[binIdx-1] + binCentDiff;
1285 }
1286 }
1287
1288 std::vector<double> weightsExt(nBins+2);
1289 // Fill weights for bins that are inside histogram boundaries
1290 for (std::size_t binIdx = 0; binIdx < nBins; ++binIdx) {
1291 weightsExt[binIdx+1] = correctForBinSize ? _wgt[binIdx] / _binv[binIdx] : _wgt[binIdx];
1292 }
1293
1294 // Fill weights for bins that are outside histogram boundaries
1295 if (cdfBoundaries) {
1296 coordsExt[0] = binning.lowBound();
1297 weightsExt[0] = 0.;
1298 coordsExt[nBins+1] = binning.highBound();
1299 weightsExt[nBins+1] = 1.;
1300 }
1301 else {
1302 // Mirror first and last bins
1303 coordsExt[0] = binCoords[0] - _binv[0];
1304 weightsExt[0] = weightsExt[1];
1305 coordsExt[nBins+1] = binCoords[nBins-1] + _binv[nBins-1];
1306 weightsExt[nBins+1] = weightsExt[nBins];
1307 }
1308
1309 // Interpolate between current bin center and one bin center to the left
1310 // if xVal is to the left of the current bin center
1311 for (std::size_t i = 0; i < nEvents ; ++i) {
1312 double xVal = xVals[i];
1313 std::size_t binIdx = binIndices[i] + 1;
1314
1315 // If xVal is to the right of the current bin center, interpolate between
1316 // current bin center and one bin center to the right instead
1317 if (xVal > coordsExt[binIdx]) { binIdx += 1; }
1318
1319 double x1 = coordsExt[binIdx-1];
1320 double y1 = weightsExt[binIdx-1];
1321 double x2 = coordsExt[binIdx];
1322 double y2 = weightsExt[binIdx];
1323
1324 // Find coefficients by solving a system of two linear equations
1325 double firstCoeff = (y2-y1) / (x2-x1);
1326 double zerothCoeff = y1 - firstCoeff * x1;
1327 // Get the interpolated weight using the equation of a straight line
1328 output[i] = firstCoeff * xVal + zerothCoeff;
1329 }
1330}
1331
1332
1333////////////////////////////////////////////////////////////////////////////////
1334/// A vectorized version of RooDataHist::weight() for one dimensional histograms
1335/// with up to one dimensional interpolation.
1336/// \param[out] output An array of weights corresponding the values in xVals.
1337/// \param[in] xVals An array of coordinates for which the weights should be
1338/// calculated.
1339/// \param[in] intOrder Interpolation order; 0th and 1st order are supported.
1340/// \param[in] correctForBinSize Enable the inverse bin volume correction factor.
1341/// \param[in] cdfBoundaries Enable the special boundary condition for a cdf:
1342/// Underflow bins are assumed to have weight zero and
1343/// overflow bins have weight one. Otherwise, the
1344/// histogram is mirrored at the boundaries for the
1345/// interpolation.
1346
1347void RooDataHist::weights(double* output, std::span<double const> xVals, int intOrder, bool correctForBinSize, bool cdfBoundaries)
1348{
1349 auto const nEvents = xVals.size();
1350
1351 if (intOrder == 0) {
1352 RooAbsBinning const& binning = *_lvbins[0];
1353
1354 // Reuse the output buffer for bin indices and zero-initialize it
1355 auto binIndices = reinterpret_cast<int*>(output + nEvents) - nEvents;
1356 std::fill(binIndices, binIndices + nEvents, 0);
1357 binning.binNumbers(xVals.data(), binIndices, nEvents);
1358
1359 for (std::size_t i=0; i < nEvents; ++i) {
1360 auto binIdx = binIndices[i];
1361 output[i] = correctForBinSize ? _wgt[binIdx] / _binv[binIdx] : _wgt[binIdx];
1362 }
1363 }
1364 else if (intOrder == 1) {
1365 interpolateLinear(output, xVals, correctForBinSize, cdfBoundaries);
1366 }
1367 else if (intOrder == 2) {
1368 interpolateQuadratic(output, xVals, correctForBinSize, cdfBoundaries);
1369 }
1370 else {
1371 // Higher dimensional scenarios not yet implemented
1372 coutE(InputArguments) << "RooDataHist::weights(" << GetName() << ") interpolation in "
1373 << intOrder << " dimensions not yet implemented" << std::endl ;
1374 // Fall back to 1st order interpolation
1375 weights(output, xVals, 1, correctForBinSize, cdfBoundaries);
1376 }
1377}
1378
1379
1380////////////////////////////////////////////////////////////////////////////////
1381/// A faster version of RooDataHist::weight that assumes the passed arguments
1382/// are aligned with the histogram variables.
1383/// \param[in] bin Coordinates for which the weight should be calculated.
1384/// Has to be aligned with the internal histogram variables.
1385/// \param[in] intOrder Interpolation order, i.e. how many neighbouring bins are
1386/// used for the interpolation. If zero, the bare weight for
1387/// the bin enclosing the coordinatesis returned.
1388/// \param[in] correctForBinSize Enable the inverse bin volume correction factor.
1389/// \param[in] cdfBoundaries Enable the special boundary condition for a cdf:
1390/// underflow bins are assumed to have weight zero and
1391/// overflow bins have weight one. Otherwise, the
1392/// histogram is mirrored at the boundaries for the
1393/// interpolation.
1394
1395double RooDataHist::weightFast(const RooArgSet& bin, Int_t intOrder, bool correctForBinSize, bool cdfBoundaries)
1396{
1397 checkInit() ;
1398
1399 // Handle illegal intOrder values
1400 if (intOrder<0) {
1401 coutE(InputArguments) << "RooDataHist::weight(" << GetName() << ") ERROR: interpolation order must be positive" << endl ;
1402 return 0 ;
1403 }
1404
1405 // Handle no-interpolation case
1406 if (intOrder==0) {
1407 const auto idx = calcTreeIndex(bin, true);
1408 return correctForBinSize ? _wgt[idx] / _binv[idx] : _wgt[idx];
1409 }
1410
1411 // Handle all interpolation cases
1412 return weightInterpolated(bin, intOrder, correctForBinSize, cdfBoundaries);
1413}
1414
1415
1416////////////////////////////////////////////////////////////////////////////////
1417/// Return the weight at given coordinates with optional interpolation.
1418/// \param[in] bin Coordinates for which the weight should be calculated.
1419/// \param[in] intOrder Interpolation order, i.e. how many neighbouring bins are
1420/// used for the interpolation. If zero, the bare weight for
1421/// the bin enclosing the coordinatesis returned.
1422/// \param[in] correctForBinSize Enable the inverse bin volume correction factor.
1423/// \param[in] cdfBoundaries Enable the special boundary condition for a cdf:
1424/// underflow bins are assumed to have weight zero and
1425/// overflow bins have weight one. Otherwise, the
1426/// histogram is mirrored at the boundaries for the
1427/// interpolation.
1428/// \param[in] oneSafe Ignored.
1429
1430double RooDataHist::weight(const RooArgSet& bin, Int_t intOrder, bool correctForBinSize, bool cdfBoundaries, bool /*oneSafe*/)
1431{
1432 checkInit() ;
1433
1434 // Handle illegal intOrder values
1435 if (intOrder<0) {
1436 coutE(InputArguments) << "RooDataHist::weight(" << GetName() << ") ERROR: interpolation order must be positive" << endl ;
1437 return 0 ;
1438 }
1439
1440 // Handle no-interpolation case
1441 if (intOrder==0) {
1442 const auto idx = calcTreeIndex(bin, false);
1443 return correctForBinSize ? _wgt[idx] / _binv[idx] : _wgt[idx];
1444 }
1445
1446 // Handle all interpolation cases
1447 _vars.assignValueOnly(bin) ;
1448
1449 return weightInterpolated(_vars, intOrder, correctForBinSize, cdfBoundaries);
1450}
1451
1452
1453////////////////////////////////////////////////////////////////////////////////
1454/// Return the weight at given coordinates with interpolation.
1455/// \param[in] bin Coordinates for which the weight should be calculated.
1456/// Has to be aligned with the internal histogram variables.
1457/// \param[in] intOrder Interpolation order, i.e. how many neighbouring bins are
1458/// used for the interpolation.
1459/// \param[in] correctForBinSize Enable the inverse bin volume correction factor.
1460/// \param[in] cdfBoundaries Enable the special boundary condition for a cdf:
1461/// underflow bins are assumed to have weight zero and
1462/// overflow bins have weight one. Otherwise, the
1463/// histogram is mirrored at the boundaries for the
1464/// interpolation.
1465
1466double RooDataHist::weightInterpolated(const RooArgSet& bin, int intOrder, bool correctForBinSize, bool cdfBoundaries) {
1467 VarInfo const& varInfo = getVarInfo();
1468
1469 const auto centralIdx = calcTreeIndex(bin, true);
1470
1471 double wInt{0} ;
1472 if (varInfo.nRealVars == 1) {
1473
1474 // buffer needs to be 2 x (interpolation order + 1), with the factor 2 for x and y.
1475 _interpolationBuffer.resize(2 * intOrder + 2);
1476
1477 // 1-dimensional interpolation
1478 auto const& realX = static_cast<RooRealVar const&>(*bin[varInfo.realVarIdx1]);
1479 wInt = interpolateDim(varInfo.realVarIdx1, realX.getVal(), centralIdx, intOrder, correctForBinSize, cdfBoundaries) ;
1480
1481 } else if (varInfo.nRealVars == 2) {
1482
1483 // buffer needs to be 2 x 2 x (interpolation order + 1), with one factor 2
1484 // for x and y, and the other for the number of dimensions.
1485 _interpolationBuffer.resize(4 * intOrder + 4);
1486
1487 // 2-dimensional interpolation
1488 auto const& realX = static_cast<RooRealVar const&>(*bin[varInfo.realVarIdx1]);
1489 auto const& realY = static_cast<RooRealVar const&>(*bin[varInfo.realVarIdx2]);
1490 double xval = realX.getVal() ;
1491 double yval = realY.getVal() ;
1492
1493 RooAbsBinning const& binningY = realY.getBinning();
1494
1495 int ybinC = binningY.binNumber(yval) ;
1496 int ybinLo = ybinC-intOrder/2 - ((yval<binningY.binCenter(ybinC))?1:0) ;
1497 int ybinM = binningY.numBins() ;
1498
1499 auto idxMultY = _idxMult[varInfo.realVarIdx2];
1500 auto offsetIdx = centralIdx - idxMultY * ybinC;
1501
1502 // Use a class-member buffer to avoid repeated heap allocations.
1503 double * yarr = _interpolationBuffer.data() + 2 * intOrder + 2; // add offset to skip part reserved for other dim
1504 double * xarr = yarr + intOrder + 1;
1505 for (int i=ybinLo ; i<=intOrder+ybinLo ; i++) {
1506 int ibin ;
1507 if (i>=0 && i<ybinM) {
1508 // In range
1509 ibin = i ;
1510 xarr[i-ybinLo] = binningY.binCenter(ibin) ;
1511 } else if (i>=ybinM) {
1512 // Overflow: mirror
1513 ibin = 2*ybinM-i-1 ;
1514 xarr[i-ybinLo] = 2*binningY.highBound()-binningY.binCenter(ibin) ;
1515 } else {
1516 // Underflow: mirror
1517 ibin = -i -1;
1518 xarr[i-ybinLo] = 2*binningY.lowBound()-binningY.binCenter(ibin) ;
1519 }
1520 auto centralIdxX = offsetIdx + idxMultY * ibin;
1521 yarr[i-ybinLo] = interpolateDim(varInfo.realVarIdx1,xval,centralIdxX,intOrder,correctForBinSize,false) ;
1522 }
1523
1524 if (gDebug>7) {
1525 cout << "RooDataHist interpolating data is" << endl ;
1526 cout << "xarr = " ;
1527 for (int q=0; q<=intOrder ; q++) cout << xarr[q] << " " ;
1528 cout << " yarr = " ;
1529 for (int q=0; q<=intOrder ; q++) cout << yarr[q] << " " ;
1530 cout << endl ;
1531 }
1532 wInt = RooMath::interpolate(xarr,yarr,intOrder+1,yval) ;
1533
1534 } else {
1535
1536 // Higher dimensional scenarios not yet implemented
1537 coutE(InputArguments) << "RooDataHist::weight(" << GetName() << ") interpolation in "
1538 << varInfo.nRealVars << " dimensions not yet implemented" << endl ;
1539 return weightFast(bin,0,correctForBinSize,cdfBoundaries) ;
1540
1541 }
1542
1543 return wInt ;
1544}
1545
1546
1548 if (!_errLo || !_errHi) {
1549 initArray(_errLo, _arrSize, -1.);
1550 initArray(_errHi, _arrSize, -1.);
1552 }
1553}
1554
1555
1556////////////////////////////////////////////////////////////////////////////////
1557/// Return the asymmetric errors on the current weight.
1558/// \see weightError(ErrorType) const for symmetric error.
1559/// \param[out] lo Low error.
1560/// \param[out] hi High error.
1561/// \param[in] etype Type of error to compute. May throw if not supported.
1562/// Supported errors are
1563/// - `Poisson` Default. Asymmetric Poisson errors (68% CL).
1564/// - `SumW2` The square root of the sum of weights. (Symmetric).
1565/// - `None` Return zero.
1566void RooDataHist::weightError(double& lo, double& hi, ErrorType etype) const
1567{
1568 checkInit() ;
1569
1570 switch (etype) {
1571
1572 case Auto:
1573 throw std::invalid_argument("RooDataHist::weightError(" + std::string(GetName()) + ") error type Auto not allowed here");
1574 break ;
1575
1576 case Expected:
1577 throw std::invalid_argument("RooDataHist::weightError(" + std::string(GetName()) + ") error type Expected not allowed here");
1578 break ;
1579
1580 case Poisson:
1581 if (_errLo && _errLo[_curIndex] >= 0.0) {
1582 // Weight is preset or precalculated
1583 lo = _errLo[_curIndex];
1584 hi = _errHi[_curIndex];
1585 return ;
1586 }
1587
1588 // We didn't track asymmetric errors so far, so now we need to allocate
1590
1591 // Calculate poisson errors
1592 double ym,yp ;
1594 _errLo[_curIndex] = weight()-ym;
1595 _errHi[_curIndex] = yp-weight();
1596 lo = _errLo[_curIndex];
1597 hi = _errHi[_curIndex];
1598 return ;
1599
1600 case SumW2:
1601 lo = std::sqrt(weightSquared(_curIndex));
1602 hi = lo;
1603 return ;
1604
1605 case None:
1606 lo = 0 ;
1607 hi = 0 ;
1608 return ;
1609 }
1610}
1611
1612
1613// wve adjust for variable bin sizes
1614
1615////////////////////////////////////////////////////////////////////////////////
1616/// Perform boundary safe 'intOrder'-th interpolation of weights in dimension 'dim'
1617/// at current value 'xval'
1618
1619/// \param[in] iDim Index of the histogram dimension along which to interpolate.
1620/// \param[in] xval Value of histogram variable at dimension `iDim` for which
1621/// we want to interpolate the histogram weight.
1622/// \param[in] centralIdx Index of the bin that the point at which we
1623/// interpolate the histogram weight falls into
1624/// (can be obtained with `RooDataHist::calcTreeIndex`).
1625/// \param[in] intOrder Interpolation order, i.e. how many neighbouring bins are
1626/// used for the interpolation.
1627/// \param[in] correctForBinSize Enable the inverse bin volume correction factor.
1628/// \param[in] cdfBoundaries Enable the special boundary condition for a cdf:
1629/// underflow bins are assumed to have weight zero and
1630/// overflow bins have weight one. Otherwise, the
1631/// histogram is mirrored at the boundaries for the
1632/// interpolation.
1633double RooDataHist::interpolateDim(int iDim, double xval, size_t centralIdx, int intOrder, bool correctForBinSize, bool cdfBoundaries)
1634{
1635 auto const& binning = static_cast<RooRealVar&>(*_vars[iDim]).getBinning();
1636
1637 // Fill workspace arrays spanning interpolation area
1638 int fbinC = binning.binNumber(xval) ;
1639 int fbinLo = fbinC-intOrder/2 - ((xval<binning.binCenter(fbinC))?1:0) ;
1640 int fbinM = binning.numBins() ;
1641
1642 auto idxMult = _idxMult[iDim];
1643 auto offsetIdx = centralIdx - idxMult * fbinC;
1644
1645 // Use a class-member buffer to avoid repeated heap allocations.
1646 double * yarr = _interpolationBuffer.data();
1647 double * xarr = yarr + intOrder + 1;
1648
1649 for (int i=fbinLo ; i<=intOrder+fbinLo ; i++) {
1650 int ibin ;
1651 if (i>=0 && i<fbinM) {
1652 // In range
1653 ibin = i ;
1654 xarr[i-fbinLo] = binning.binCenter(ibin) ;
1655 auto idx = offsetIdx + idxMult * ibin;
1656 yarr[i - fbinLo] = _wgt[idx];
1657 if (correctForBinSize) yarr[i-fbinLo] /= _binv[idx] ;
1658 } else if (i>=fbinM) {
1659 // Overflow: mirror
1660 ibin = 2*fbinM-i-1 ;
1661 if (cdfBoundaries) {
1662 xarr[i-fbinLo] = binning.highBound()+1e-10*(i-fbinM+1) ;
1663 yarr[i-fbinLo] = 1.0 ;
1664 } else {
1665 auto idx = offsetIdx + idxMult * ibin;
1666 xarr[i-fbinLo] = 2*binning.highBound()-binning.binCenter(ibin) ;
1667 yarr[i - fbinLo] = _wgt[idx];
1668 if (correctForBinSize)
1669 yarr[i - fbinLo] /= _binv[idx];
1670 }
1671 } else {
1672 // Underflow: mirror
1673 ibin = -i - 1 ;
1674 if (cdfBoundaries) {
1675 xarr[i-fbinLo] = binning.lowBound()-ibin*(1e-10) ;
1676 yarr[i-fbinLo] = 0.0 ;
1677 } else {
1678 auto idx = offsetIdx + idxMult * ibin;
1679 xarr[i-fbinLo] = 2*binning.lowBound()-binning.binCenter(ibin) ;
1680 yarr[i - fbinLo] = _wgt[idx];
1681 if (correctForBinSize)
1682 yarr[i - fbinLo] /= _binv[idx];
1683 }
1684 }
1685 }
1686 return RooMath::interpolate(xarr,yarr,intOrder+1,xval) ;
1687}
1688
1689
1690
1691
1692////////////////////////////////////////////////////////////////////////////////
1693/// Increment the bin content of the bin enclosing the given coordinates.
1694///
1695/// \param[in] row Coordinates of the bin.
1696/// \param[in] wgt Increment by this weight.
1697/// \param[in] sumw2 Optionally, track the sum of squared weights. If a value > 0 or
1698/// a weight != 1. is passed for the first time, a vector for the squared weights will be allocated.
1699void RooDataHist::add(const RooArgSet& row, double wgt, double sumw2)
1700{
1701 checkInit() ;
1702
1703 if ((sumw2 > 0. || wgt != 1.) && !_sumw2) {
1704 // Receiving a weighted entry. SumW2 != sumw from now on.
1705 _sumw2 = new double[_arrSize];
1706 std::copy(_wgt, _wgt+_arrSize, _sumw2);
1707
1709 }
1710
1711 const auto idx = calcTreeIndex(row, false);
1712
1713 _wgt[idx] += wgt ;
1714 if (_sumw2) _sumw2[idx] += (sumw2 > 0 ? sumw2 : wgt*wgt);
1715
1716 _cache_sum_valid = false;
1717}
1718
1719
1720
1721////////////////////////////////////////////////////////////////////////////////
1722/// Set a bin content.
1723/// \param[in] row Coordinates of the bin to be set.
1724/// \param[in] wgt New bin content.
1725/// \param[in] wgtErrLo Low error of the bin content.
1726/// \param[in] wgtErrHi High error of the bin content.
1727void RooDataHist::set(const RooArgSet& row, double wgt, double wgtErrLo, double wgtErrHi)
1728{
1729 checkInit() ;
1730
1732
1733 const auto idx = calcTreeIndex(row, false);
1734
1735 _wgt[idx] = wgt ;
1736 _errLo[idx] = wgtErrLo ;
1737 _errHi[idx] = wgtErrHi ;
1738
1739 _cache_sum_valid = false;
1740}
1741
1742
1743
1744////////////////////////////////////////////////////////////////////////////////
1745/// Set bin content of bin that was last loaded with get(std::size_t).
1746/// \param[in] binNumber Optional bin number to set. If empty, currently active bin is set.
1747/// \param[in] wgt New bin content.
1748/// \param[in] wgtErr Error of the new bin content. If the weight need not have an error, use 0. or a negative number.
1749void RooDataHist::set(std::size_t binNumber, double wgt, double wgtErr) {
1750 checkInit() ;
1751
1752 if (wgtErr > 0. && !_sumw2) {
1753 // Receiving a weighted entry. Need to track sumw2 from now on:
1754 cloneArray(_sumw2, _wgt, _arrSize);
1755
1757 }
1758
1759 _wgt[binNumber] = wgt ;
1760 if (_errLo) _errLo[binNumber] = wgtErr;
1761 if (_errHi) _errHi[binNumber] = wgtErr;
1762 if (_sumw2) _sumw2[binNumber] = wgtErr*wgtErr;
1763
1764 _cache_sum_valid = false ;
1765}
1766
1767
1768////////////////////////////////////////////////////////////////////////////////
1769/// Set bin content of bin that was last loaded with get(std::size_t).
1770/// \param[in] wgt New bin content.
1771/// \param[in] wgtErr Optional error of the bin content.
1772void RooDataHist::set(double wgt, double wgtErr) {
1773 if (_curIndex == std::numeric_limits<std::size_t>::max()) {
1774 _curIndex = calcTreeIndex(_vars, true) ;
1775 }
1776
1777 set(_curIndex, wgt, wgtErr);
1778}
1779
1780
1781////////////////////////////////////////////////////////////////////////////////
1782/// Set a bin content.
1783/// \param[in] row Coordinates to compute the bin from.
1784/// \param[in] wgt New bin content.
1785/// \param[in] wgtErr Optional error of the bin content.
1786void RooDataHist::set(const RooArgSet& row, double wgt, double wgtErr) {
1787 set(calcTreeIndex(row, false), wgt, wgtErr);
1788}
1789
1790
1791
1792////////////////////////////////////////////////////////////////////////////////
1793/// Add all data points contained in 'dset' to this data set with given weight.
1794/// Optional cut string expression selects the data points to be added and can
1795/// reference any variable contained in this data set
1796
1797void RooDataHist::add(const RooAbsData& dset, const char* cut, double wgt)
1798{
1799 RooFormulaVar cutVar("select",cut,*dset.get()) ;
1800 add(dset,&cutVar,wgt) ;
1801}
1802
1803
1804
1805////////////////////////////////////////////////////////////////////////////////
1806/// Add all data points contained in 'dset' to this data set with given weight.
1807/// Optional RooFormulaVar pointer selects the data points to be added.
1808
1809void RooDataHist::add(const RooAbsData& dset, const RooFormulaVar* cutVar, double wgt)
1810{
1811 checkInit() ;
1812
1813 RooFormulaVar* cloneVar = nullptr;
1814 std::unique_ptr<RooArgSet> tmp;
1815 if (cutVar) {
1816 // Deep clone cutVar and attach clone to this dataset
1817 tmp = std::make_unique<RooArgSet>();
1818 if(RooArgSet(*cutVar).snapshot(*tmp)) {
1819 coutE(DataHandling) << "RooDataHist::add(" << GetName() << ") Couldn't deep-clone cut variable, abort," << endl ;
1820 return ;
1821 }
1822
1823 cloneVar = (RooFormulaVar*) tmp->find(*cutVar) ;
1824 cloneVar->attachDataSet(dset) ;
1825 }
1826
1827
1828 Int_t i ;
1829 for (i=0 ; i<dset.numEntries() ; i++) {
1830 const RooArgSet* row = dset.get(i) ;
1831 if (!cloneVar || cloneVar->getVal()) {
1832 add(*row,wgt*dset.weight(), wgt*wgt*dset.weightSquared()) ;
1833 }
1834 }
1835
1836 _cache_sum_valid = false ;
1837}
1838
1839
1840
1841////////////////////////////////////////////////////////////////////////////////
1842/// Return the sum of the weights of all bins in the histogram.
1843///
1844/// \param[in] correctForBinSize Multiply the sum of weights in each bin
1845/// with the N-dimensional bin volume, making the return value
1846/// the integral over the function represented by this histogram.
1847/// \param[in] inverseBinCor Divide by the N-dimensional bin volume.
1848double RooDataHist::sum(bool correctForBinSize, bool inverseBinCor) const
1849{
1850 checkInit() ;
1851
1852 // Check if result was cached
1853 const CacheSumState_t cache_code = !correctForBinSize ? kNoBinCorrection : (inverseBinCor ? kInverseBinCorr : kCorrectForBinSize);
1854 if (_cache_sum_valid == static_cast<Int_t>(cache_code)) {
1855 return _cache_sum ;
1856 }
1857
1859 for (Int_t i=0; i < _arrSize; i++) {
1860 const double theBinVolume = correctForBinSize ? (inverseBinCor ? 1/_binv[i] : _binv[i]) : 1.0 ;
1861 kahanSum += _wgt[i] * theBinVolume;
1862 }
1863
1864 // Store result in cache
1865 _cache_sum_valid = cache_code;
1866 _cache_sum = kahanSum.Sum();
1867
1868 return kahanSum.Sum();
1869}
1870
1871
1872
1873////////////////////////////////////////////////////////////////////////////////
1874/// Return the sum of the weights of a multi-dimensional slice of the histogram
1875/// by summing only over the dimensions specified in sumSet.
1876///
1877/// The coordinates of all other dimensions are fixed to those given in sliceSet
1878///
1879/// If correctForBinSize is specified, the sum of weights
1880/// is multiplied by the M-dimensional bin volume, (M = N(sumSet)),
1881/// making the return value the integral over the function
1882/// represented by this histogram
1883
1884double RooDataHist::sum(const RooArgSet& sumSet, const RooArgSet& sliceSet, bool correctForBinSize, bool inverseBinCor)
1885{
1886 checkInit() ;
1887
1888 RooArgSet varSave ;
1889 varSave.addClone(_vars) ;
1890
1891 RooArgSet sliceOnlySet(sliceSet);
1892 sliceOnlySet.remove(sumSet,true,true) ;
1893
1894 _vars.assign(sliceOnlySet);
1895 std::vector<double> const * pbinv = nullptr;
1896
1897 if(correctForBinSize && inverseBinCor) {
1898 pbinv = &calculatePartialBinVolume(sliceOnlySet);
1899 } else if(correctForBinSize && !inverseBinCor) {
1900 pbinv = &calculatePartialBinVolume(sumSet);
1901 }
1902
1903 // Calculate mask and reference plot bins for non-iterating variables
1904 std::vector<bool> mask(_vars.size());
1905 std::vector<int> refBin(_vars.size());
1906
1907 for (unsigned int i = 0; i < _vars.size(); ++i) {
1908 const RooAbsArg* arg = _vars[i];
1909 const RooAbsLValue* argLv = _lvvars[i]; // Same as above, but cross-cast
1910
1911 if (sumSet.find(*arg)) {
1912 mask[i] = false ;
1913 } else {
1914 mask[i] = true ;
1915 refBin[i] = argLv->getBin();
1916 }
1917 }
1918
1919 // Loop over entire data set, skipping masked entries
1921 for (Int_t ibin=0; ibin < _arrSize; ++ibin) {
1922
1923 std::size_t tmpibin = ibin;
1924 bool skip(false) ;
1925
1926 // Check if this bin belongs in selected slice
1927 for (unsigned int ivar = 0; !skip && ivar < _vars.size(); ++ivar) {
1928 const Int_t idx = tmpibin / _idxMult[ivar] ;
1929 tmpibin -= idx*_idxMult[ivar] ;
1930 if (mask[ivar] && idx!=refBin[ivar])
1931 skip = true ;
1932 }
1933
1934 if (!skip) {
1935 const double theBinVolume = correctForBinSize ? (inverseBinCor ? 1/(*pbinv)[ibin] : (*pbinv)[ibin] ) : 1.0 ;
1936 total += _wgt[ibin] * theBinVolume;
1937 }
1938 }
1939
1940 _vars.assign(varSave) ;
1941
1942 return total.Sum();
1943}
1944
1945////////////////////////////////////////////////////////////////////////////////
1946/// Return the sum of the weights of a multi-dimensional slice of the histogram
1947/// by summing only over the dimensions specified in sumSet.
1948///
1949/// The coordinates of all other dimensions are fixed to those given in sliceSet
1950///
1951/// If correctForBinSize is specified, the sum of weights
1952/// is multiplied by the M-dimensional bin volume, (M = N(sumSet)),
1953/// or the fraction of it that falls inside the range rangeName,
1954/// making the return value the integral over the function
1955/// represented by this histogram.
1956///
1957/// If correctForBinSize is not specified, the weights are multiplied by the
1958/// fraction of the bin volume that falls inside the range, i.e. a factor of
1959/// binVolumeInRange/totalBinVolume.
1960
1961double RooDataHist::sum(const RooArgSet& sumSet, const RooArgSet& sliceSet,
1962 bool correctForBinSize, bool inverseBinCor,
1963 const std::map<const RooAbsArg*, std::pair<double, double> >& ranges,
1964 std::function<double(int)> getBinScale)
1965{
1966 checkInit();
1968 RooArgSet varSave;
1969 varSave.addClone(_vars);
1970 {
1971 RooArgSet sliceOnlySet(sliceSet);
1972 sliceOnlySet.remove(sumSet, true, true);
1973 _vars.assign(sliceOnlySet);
1974 }
1975
1976 // Calculate mask and reference plot bins for non-iterating variables,
1977 // and get ranges for iterating variables
1978 std::vector<bool> mask(_vars.size());
1979 std::vector<int> refBin(_vars.size());
1980 std::vector<double> rangeLo(_vars.size(), -std::numeric_limits<double>::infinity());
1981 std::vector<double> rangeHi(_vars.size(), +std::numeric_limits<double>::infinity());
1982
1983 for (std::size_t i = 0; i < _vars.size(); ++i) {
1984 const RooAbsArg* arg = _vars[i];
1985 const RooAbsLValue* argLV = _lvvars[i]; // Same object as above, but cross cast
1986
1987 RooAbsArg* sumsetv = sumSet.find(*arg);
1988 RooAbsArg* slicesetv = sliceSet.find(*arg);
1989 mask[i] = !sumsetv;
1990 if (mask[i]) {
1991 assert(argLV);
1992 refBin[i] = argLV->getBin();
1993 }
1994
1995 auto it = ranges.find(sumsetv ? sumsetv : slicesetv);
1996 if (ranges.end() != it) {
1997 rangeLo[i] = it->second.first;
1998 rangeHi[i] = it->second.second;
1999 }
2000 }
2001
2002 // Loop over entire data set, skipping masked entries
2004 for (Int_t ibin = 0; ibin < _arrSize; ++ibin) {
2005 // Check if this bin belongs in selected slice
2006 bool skip{false};
2007 for (int ivar = 0, tmp = ibin; !skip && ivar < int(_vars.size()); ++ivar) {
2008 const Int_t idx = tmp / _idxMult[ivar];
2009 tmp -= idx*_idxMult[ivar];
2010 if (mask[ivar] && idx!=refBin[ivar]) skip = true;
2011 }
2012
2013 if (skip) continue;
2014
2015 // Work out bin volume
2016 // It's not necessary to figure out the bin volume for the slice-only set explicitly here.
2017 // We need to loop over the sumSet anyway to get the partial bin containment correction,
2018 // so we can get the slice-only set volume later by dividing _binv[ibin] / binVolumeSumSetFull.
2019 double binVolumeSumSetFull = 1.;
2020 double binVolumeSumSetInRange = 1.;
2021 for (Int_t ivar = 0, tmp = ibin; ivar < (int)_vars.size(); ++ivar) {
2022 const Int_t idx = tmp / _idxMult[ivar];
2023 tmp -= idx*_idxMult[ivar];
2024
2025 // If the current variable is not in the sumSet, it should not be considered for the bin volume
2026 const auto arg = _vars[ivar];
2027 if (!sumSet.find(*arg)) {
2028 continue;
2029 }
2030
2031 if (_binbounds[ivar].empty()) continue;
2032 const double binLo = _binbounds[ivar][2 * idx];
2033 const double binHi = _binbounds[ivar][2 * idx + 1];
2034 if (binHi < rangeLo[ivar] || binLo > rangeHi[ivar]) {
2035 // bin is outside of allowed range - effective bin volume is zero
2036 binVolumeSumSetInRange = 0.;
2037 break;
2038 }
2039
2040 binVolumeSumSetFull *= binHi - binLo;
2041 binVolumeSumSetInRange *= std::min(rangeHi[ivar], binHi) - std::max(rangeLo[ivar], binLo);
2042 }
2043 const double corrPartial = binVolumeSumSetInRange / binVolumeSumSetFull;
2044 if (0. == corrPartial) continue;
2045 const double corr = correctForBinSize ? (inverseBinCor ? binVolumeSumSetFull / _binv[ibin] : binVolumeSumSetFull ) : 1.0;
2046 total += getBinScale(ibin)*(_wgt[ibin] * corr * corrPartial);
2047 }
2048
2049 _vars.assign(varSave);
2050
2051 return total.Sum();
2052}
2053
2054
2055
2056////////////////////////////////////////////////////////////////////////////////
2057/// Fill the transient cache with partial bin volumes with up-to-date
2058/// values for the partial volume specified by observables 'dimSet'
2059
2060const std::vector<double>& RooDataHist::calculatePartialBinVolume(const RooArgSet& dimSet) const
2061{
2062 // The code bitset has all bits set to one whose position corresponds to arguments in dimSet.
2063 // It is used as the key for the bin volume caching hash map.
2064 int code{0};
2065 {
2066 int i{0} ;
2067 for (auto const& v : _vars) {
2068 code += ((dimSet.find(*v) ? 1 : 0) << i) ;
2069 ++i;
2070 }
2071 }
2072
2073 auto& pbinv = _pbinvCache[code];
2074 if(!pbinv.empty()) {
2075 return pbinv;
2076 }
2077 pbinv.resize(_arrSize);
2078
2079 // Calculate plot bins of components from master index
2080 std::vector<bool> selDim(_vars.size());
2081 for (std::size_t i = 0; i < selDim.size(); ++i) {
2082 selDim[i] = (code >> i) & 1 ;
2083 }
2084
2085 // Recalculate partial bin volume cache
2086 for (Int_t ibin=0; ibin < _arrSize ;ibin++) {
2087 Int_t idx(0), tmp(ibin) ;
2088 double theBinVolume(1) ;
2089 for (unsigned int j=0; j < _lvvars.size(); ++j) {
2090 const RooAbsLValue* arg = _lvvars[j];
2091 assert(arg);
2092
2093 idx = tmp / _idxMult[j];
2094 tmp -= idx*_idxMult[j];
2095 if (selDim[j]) {
2096 theBinVolume *= arg->getBinWidth(idx) ;
2097 }
2098 }
2099 pbinv[ibin] = theBinVolume ;
2100 }
2101
2102 return pbinv;
2103}
2104
2105
2106////////////////////////////////////////////////////////////////////////////////
2107/// Sum the weights of all bins.
2110}
2111
2112
2113
2114////////////////////////////////////////////////////////////////////////////////
2115/// Return the sum of weights in all entries matching cutSpec (if specified)
2116/// and in named range cutRange (if specified)
2117/// Return the
2118
2119double RooDataHist::sumEntries(const char* cutSpec, const char* cutRange) const
2120{
2121 checkInit() ;
2122
2123 if (cutSpec==nullptr && cutRange==nullptr) {
2124 return sumEntries();
2125 } else {
2126
2127 // Setup RooFormulaVar for cutSpec if it is present
2128 std::unique_ptr<RooFormula> select;
2129 if (cutSpec) {
2130 select = std::make_unique<RooFormula>("select",cutSpec,*get());
2131 }
2132
2133 // Otherwise sum the weights in the event
2134 ROOT::Math::KahanSum<> kahanSum;
2135 for (Int_t i=0; i < _arrSize; i++) {
2136 get(i) ;
2137 if ((select && select->eval() == 0.) || (cutRange && !_vars.allInRange(cutRange)))
2138 continue;
2139
2140 kahanSum += weight(i);
2141 }
2142
2143 return kahanSum.Sum();
2144 }
2145}
2146
2147
2148
2149////////////////////////////////////////////////////////////////////////////////
2150/// Reset all bin weights to zero
2151
2153{
2154 // WVE DO NOT CALL RooTreeData::reset() for binned
2155 // datasets as this will delete the bin definitions
2156
2157 std::fill(_wgt, _wgt + _arrSize, 0.);
2158 delete[] _errLo; _errLo = nullptr;
2159 delete[] _errHi; _errHi = nullptr;
2160 delete[] _sumw2; _sumw2 = nullptr;
2161
2163
2164 _cache_sum_valid = false;
2165}
2166
2167
2168
2169////////////////////////////////////////////////////////////////////////////////
2170/// Load bin `binNumber`, and return an argset with the coordinates of the bin centre.
2171/// \note The argset is owned by this data hist, and this function has a side effect, because
2172/// it alters the currently active bin.
2173const RooArgSet* RooDataHist::get(Int_t binNumber) const
2174{
2175 checkInit() ;
2176 _curIndex = binNumber;
2177
2178 return RooAbsData::get(_curIndex);
2179}
2180
2181
2182
2183////////////////////////////////////////////////////////////////////////////////
2184/// Return a RooArgSet with whose coordinates denote the bin centre of the bin
2185/// enclosing the point in `coord`.
2186/// \note The argset is owned by this data hist, and this function has a side effect, because
2187/// it alters the currently active bin.
2188const RooArgSet* RooDataHist::get(const RooArgSet& coord) const {
2189 return get(calcTreeIndex(coord, false));
2190}
2191
2192
2193
2194////////////////////////////////////////////////////////////////////////////////
2195/// Return the volume of the bin enclosing coordinates 'coord'.
2196double RooDataHist::binVolume(const RooArgSet& coord) const {
2197 checkInit() ;
2198 return _binv[calcTreeIndex(coord, false)] ;
2199}
2200
2201
2202////////////////////////////////////////////////////////////////////////////////
2203/// Create an iterator over all bins in a slice defined by the subset of observables
2204/// listed in sliceArg. The position of the slice is given by otherArgs
2205
2207{
2208 // Update to current position
2209 _vars.assign(otherArgs) ;
2210 _curIndex = calcTreeIndex(_vars, true);
2211
2212 RooAbsArg* intArg = _vars.find(sliceArg) ;
2213 if (!intArg) {
2214 coutE(InputArguments) << "RooDataHist::sliceIterator() variable " << sliceArg.GetName() << " is not part of this RooDataHist" << endl ;
2215 return nullptr ;
2216 }
2217 return new RooDataHistSliceIter(*this,*intArg) ;
2218}
2219
2220
2221////////////////////////////////////////////////////////////////////////////////
2222/// Change the name of the RooDataHist
2223
2225{
2226 if (_dir) _dir->GetList()->Remove(this);
2227 // We need to use the function from RooAbsData, because it already overrides TNamed::SetName
2229 if (_dir) _dir->GetList()->Add(this);
2230}
2231
2232
2233////////////////////////////////////////////////////////////////////////////////
2234/// Change the title of this RooDataHist
2235
2236void RooDataHist::SetNameTitle(const char *name, const char* title)
2237{
2238 SetName(name);
2239 SetTitle(title);
2240}
2241
2242
2243////////////////////////////////////////////////////////////////////////////////
2244/// Print value of the dataset, i.e. the sum of weights contained in the dataset
2245
2246void RooDataHist::printValue(ostream& os) const
2247{
2248 os << numEntries() << " bins (" << sumEntries() << " weights)" ;
2249}
2250
2251
2252
2253
2254////////////////////////////////////////////////////////////////////////////////
2255/// Print argument of dataset, i.e. the observable names
2256
2257void RooDataHist::printArgs(ostream& os) const
2258{
2259 os << "[" ;
2260 bool first(true) ;
2261 for (const auto arg : _vars) {
2262 if (first) {
2263 first=false ;
2264 } else {
2265 os << "," ;
2266 }
2267 os << arg->GetName() ;
2268 }
2269 os << "]" ;
2270}
2271
2272
2273
2274////////////////////////////////////////////////////////////////////////////////
2275/// Returns true if dataset contains entries with a non-integer weight.
2276
2278{
2279 for (Int_t i=0; i < _arrSize; ++i) {
2280 const double wgt = _wgt[i];
2281 double intpart;
2282 if (std::abs(std::modf(wgt, &intpart)) > 1.E-10)
2283 return true;
2284 }
2285
2286 return false;
2287}
2288
2289
2290////////////////////////////////////////////////////////////////////////////////
2291/// Print the details on the dataset contents
2292
2293void RooDataHist::printMultiline(ostream& os, Int_t content, bool verbose, TString indent) const
2294{
2295 RooAbsData::printMultiline(os,content,verbose,indent) ;
2296
2297 os << indent << "Binned Dataset " << GetName() << " (" << GetTitle() << ")" << endl ;
2298 os << indent << " Contains " << numEntries() << " bins with a total weight of " << sumEntries() << endl;
2299
2300 if (!verbose) {
2301 os << indent << " Observables " << _vars << endl ;
2302 } else {
2303 os << indent << " Observables: " ;
2305 }
2306
2307 if(verbose) {
2308 if (!_cachedVars.empty()) {
2309 os << indent << " Caches " << _cachedVars << endl ;
2310 }
2311 }
2312}
2313
2314void RooDataHist::printDataHistogram(ostream& os, RooRealVar* obs) const
2315{
2316 for(Int_t i=0; i<obs->getBins(); ++i){
2317 this->get(i);
2318 obs->setBin(i);
2319 os << this->weight() << " +/- " << this->weightSquared() << endl;
2320 }
2321}
2322
2323
2324////////////////////////////////////////////////////////////////////////////////
2325/// Stream an object of class RooDataHist.
2327 if (R__b.IsReading()) {
2328
2329 UInt_t R__s, R__c;
2330 Version_t R__v = R__b.ReadVersion(&R__s, &R__c);
2331
2332 if (R__v > 2) {
2333 R__b.ReadClassBuffer(RooDataHist::Class(),this,R__v,R__s,R__c);
2334 R__b.CheckByteCount(R__s, R__c, RooDataHist::IsA());
2335 initialize(nullptr, false);
2336 } else {
2337
2338 // Legacy dataset conversion happens here. Legacy RooDataHist inherits from RooTreeData
2339 // which in turn inherits from RooAbsData. Manually stream RooTreeData contents on
2340 // file here and convert it into a RooTreeDataStore which is installed in the
2341 // new-style RooAbsData base class
2342
2343 // --- This is the contents of the streamer code of RooTreeData version 2 ---
2344 UInt_t R__s1, R__c1;
2345 Version_t R__v1 = R__b.ReadVersion(&R__s1, &R__c1); if (R__v1) { }
2346
2348 TTree* X_tree(nullptr) ; R__b >> X_tree;
2349 RooArgSet X_truth ; X_truth.Streamer(R__b);
2350 TString X_blindString ; X_blindString.Streamer(R__b);
2351 R__b.CheckByteCount(R__s1, R__c1, TClass::GetClass("RooTreeData"));
2352 // --- End of RooTreeData-v1 streamer
2353
2354 // Construct RooTreeDataStore from X_tree and complete initialization of new-style RooAbsData
2355 _dstore = std::make_unique<RooTreeDataStore>(X_tree,_vars);
2356 _dstore->SetName(GetName()) ;
2357 _dstore->SetTitle(GetTitle()) ;
2358 _dstore->checkInit() ;
2359
2361 R__b >> _arrSize;
2362 delete [] _wgt;
2363 _wgt = new double[_arrSize];
2365 delete [] _errLo;
2366 _errLo = new double[_arrSize];
2368 delete [] _errHi;
2369 _errHi = new double[_arrSize];
2371 delete [] _sumw2;
2372 _sumw2 = new double[_arrSize];
2374 delete [] _binv;
2375 _binv = new double[_arrSize];
2376 RooArgSet tmpSet;
2377 tmpSet.Streamer(R__b);
2378 double tmp;
2379 R__b >> tmp; //_curWeight;
2380 R__b >> tmp; //_curWgtErrLo;
2381 R__b >> tmp; //_curWgtErrHi;
2382 R__b >> tmp; //_curSumW2;
2383 R__b >> tmp; //_curVolume;
2384 R__b >> _curIndex;
2385 R__b.CheckByteCount(R__s, R__c, RooDataHist::IsA());
2386 }
2387
2388 } else {
2389
2391 }
2392}
2393
2394
2395////////////////////////////////////////////////////////////////////////////////
2396/// Return event weights of all events in range [first, first+len).
2397/// If cacheValidEntries() has been called, out-of-range events will have a weight of 0.
2398std::span<const double> RooDataHist::getWeightBatch(std::size_t first, std::size_t len, bool sumW2 /*=false*/) const {
2399 return {(sumW2 && _sumw2 ? _sumw2 : _wgt) + first, len};
2400}
2401
2402
2403////////////////////////////////////////////////////////////////////////////////
2404/// Hand over pointers to our weight arrays to the data store implementation.
2406 _dstore->setExternalWeightArray(_wgt, _errLo, _errHi, _sumw2);
2407}
2408
2409
2410////////////////////////////////////////////////////////////////////////////////
2411/// Return reference to VarInfo struct with cached histogram variable
2412/// information that is frequently used for histogram weights retrieval.
2413///
2414/// If the `_varInfo` struct was not initialized yet, it will be initialized in
2415/// this function.
2417
2418 if(_varInfo.initialized) return _varInfo;
2419
2420 auto& info = _varInfo;
2421
2422 {
2423 // count the number of real vars and get their indices
2424 info.nRealVars = 0;
2425 size_t iVar = 0;
2426 for (const auto real : _vars) {
2427 if (dynamic_cast<RooRealVar*>(real)) {
2428 if(info.nRealVars == 0) info.realVarIdx1 = iVar;
2429 if(info.nRealVars == 1) info.realVarIdx2 = iVar;
2430 ++info.nRealVars;
2431 }
2432 ++iVar;
2433 }
2434 }
2435
2436 {
2437 // assert that the variables are either real values or categories
2438 for (unsigned int i=0; i < _vars.size(); ++i) {
2439 if (_lvbins[i].get()) {
2440 assert(dynamic_cast<const RooAbsReal*>(_vars[i]));
2441 } else {
2442 assert(dynamic_cast<const RooAbsCategoryLValue*>(_vars[i]));
2443 }
2444 }
2445 }
2446
2447 info.initialized = true;
2448
2449 return info;
2450}
#define e(i)
Definition RSha256.hxx:103
#define coutI(a)
#define coutE(a)
#define TRACE_DESTROY
Definition RooTrace.h:24
#define TRACE_CREATE
Definition RooTrace.h:23
int Int_t
Definition RtypesCore.h:45
short Version_t
Definition RtypesCore.h:65
#define ClassImp(name)
Definition Rtypes.h:377
static void indent(ostringstream &buf, int indent_level)
static unsigned int total
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void data
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t Float_t Float_t Float_t Int_t Int_t UInt_t UInt_t Rectangle_t mask
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t Float_t Float_t Float_t Int_t Int_t UInt_t UInt_t Rectangle_t Int_t Int_t Window_t TString Int_t GCValues_t GetPrimarySelectionOwner GetDisplay GetScreen GetColormap GetNativeEvent const char const char dpyName wid window const char font_name cursor keysym reg const char only_if_exist regb h Point_t winding char text const char depth char const char Int_t count const char ColorStruct_t color const char Pixmap_t Pixmap_t PictureAttributes_t attr const char char ret_data h unsigned char height h offset
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t Float_t Float_t Float_t Int_t Int_t UInt_t UInt_t Rectangle_t Int_t Int_t Window_t TString Int_t GCValues_t GetPrimarySelectionOwner GetDisplay GetScreen GetColormap GetNativeEvent const char const char dpyName wid window const char font_name cursor keysym reg const char only_if_exist regb h Point_t winding char text const char depth char const char Int_t count const char ColorStruct_t color const char Pixmap_t Pixmap_t PictureAttributes_t attr const char char ret_data h unsigned char height h Atom_t Int_t ULong_t ULong_t unsigned char prop_list Atom_t Atom_t Atom_t Time_t UChar_t len
Option_t Option_t TPoint TPoint const char x2
Option_t Option_t TPoint TPoint const char x1
Option_t Option_t TPoint TPoint const char y2
Option_t Option_t TPoint TPoint const char y1
char name[80]
Definition TGX11.cxx:110
float xmin
#define hi
float * q
float ymin
Int_t gDebug
Definition TROOT.cxx:595
The Kahan summation is a compensated summation algorithm, which significantly reduces numerical error...
Definition Util.h:122
static KahanSum< T, N > Accumulate(Iterator begin, Iterator end, T initialValue=T{})
Iterate over a range and return an instance of a KahanSum.
Definition Util.h:211
RooAbsArg is the common abstract base class for objects that represent a value and a "shape" in RooFi...
Definition RooAbsArg.h:80
void attachDataSet(const RooAbsData &set)
Replace server nodes with names matching the dataset variable names with those data set variables,...
RooAbsBinning is the abstract base class for RooRealVar binning definitions.
int binNumber(double x) const
Returns the bin number corresponding to the value x.
virtual void binNumbers(double const *x, int *bins, std::size_t n, int coef=1) const =0
Compute the bin indices for multiple values of x.
virtual double binCenter(Int_t bin) const =0
Int_t numBins() const
Return number of bins.
virtual bool isUniform() const
virtual double highBound() const =0
virtual double lowBound() const =0
virtual RooAbsBinning * clone(const char *name=nullptr) const =0
RooAbsCategoryLValue is the common abstract base class for objects that represent a discrete value th...
Int_t numBins(const char *rangeName=nullptr) const override
Return the number of fit bins ( = number of types )
void setBin(Int_t ibin, const char *rangeName=nullptr) override
Set category to i-th fit bin, which is the i-th registered state.
Int_t getBin(const char *=nullptr) const override
Get the index of the plot bin for the current value of this category.
bool hasLabel(const std::string &label) const
Check if a state with name label exists.
virtual const char * getCurrentLabel() const
Return label string of current state.
RooAbsCollection is an abstract container object that can hold multiple RooAbsArg objects.
virtual bool remove(const RooAbsArg &var, bool silent=false, bool matchByNameOnly=false)
Remove the specified argument from our list.
RooAbsCollection & assignValueOnly(const RooAbsCollection &other, bool forceIfSizeOne=false)
Sets the value of any argument in our set that also appears in the other set.
bool allInRange(const char *rangeSpec) const
Return true if all contained object report to have their value inside the specified range.
Int_t getSize() const
Return the number of elements in the collection.
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
virtual RooAbsArg * addClone(const RooAbsArg &var, bool silent=false)
Add a clone of the specified argument to list.
bool selectCommon(const RooAbsCollection &refColl, RooAbsCollection &outColl) const
Create a subset of the current collection, consisting only of those elements that are contained as we...
RooAbsArg * find(const char *name) const
Find object with given name in list.
RooAbsData is the common abstract base class for binned and unbinned datasets.
Definition RooAbsData.h:57
virtual double weight() const =0
virtual const RooArgSet * get() const
Definition RooAbsData.h:101
void printMultiline(std::ostream &os, Int_t contents, bool verbose=false, TString indent="") const override
Interface for detailed printing of object.
void SetName(const char *name) override
Set the name of the TNamed.
void setGlobalObservables(RooArgSet const &globalObservables)
Sets the global observables stored in this data.
void checkInit() const
virtual RooPlot * plotOn(RooPlot *frame, const RooCmdArg &arg1={}, const RooCmdArg &arg2={}, const RooCmdArg &arg3={}, const RooCmdArg &arg4={}, const RooCmdArg &arg5={}, const RooCmdArg &arg6={}, const RooCmdArg &arg7={}, const RooCmdArg &arg8={}) const
static StorageType defaultStorageType
Definition RooAbsData.h:312
std::unique_ptr< RooAbsDataStore > _dstore
Data storage implementation.
Definition RooAbsData.h:351
virtual void fill()
RooArgSet _vars
Dimensions of this data set.
Definition RooAbsData.h:348
RooArgSet _cachedVars
! External variables cached with this data set
Definition RooAbsData.h:349
virtual Int_t numEntries() const
Return number of entries in dataset, i.e., count unweighted entries.
virtual double weightSquared() const =0
void Streamer(TBuffer &) override
Stream an object of class RooAbsData.
Abstract base class for objects that are lvalues, i.e.
virtual Int_t getBin(const char *rangeName=nullptr) const =0
virtual double getBinWidth(Int_t i, const char *rangeName=nullptr) const =0
RooAbsRealLValue is the common abstract base class for objects that represent a real value that may a...
virtual Int_t getBins(const char *name=nullptr) const
Get number of bins of currently defined range.
virtual double getMax(const char *name=nullptr) const
Get maximum of currently defined range.
void setBin(Int_t ibin, const char *rangeName=nullptr) override
Set value to center of bin 'ibin' of binning 'rangeName' (or of default binning if no range is specif...
virtual double getMin(const char *name=nullptr) const
Get minimum of currently defined range.
RooAbsReal is the common abstract base class for objects that represent a real value and implements f...
Definition RooAbsReal.h:59
double getVal(const RooArgSet *normalisationSet=nullptr) const
Evaluate object.
Definition RooAbsReal.h:103
RooArgList is a container object that can hold multiple RooAbsArg objects.
Definition RooArgList.h:22
RooAbsArg * at(Int_t idx) const
Return object at given index, or nullptr if index is out of range.
Definition RooArgList.h:110
RooArgSet is a container object that can hold multiple RooAbsArg objects.
Definition RooArgSet.h:55
RooArgSet * snapshot(bool deepCopy=true) const
Use RooAbsCollection::snapshot(), but return as RooArgSet.
Definition RooArgSet.h:178
void Streamer(TBuffer &) override
Stream an object of class TObject.
Class RooBinning is an implements RooAbsBinning in terms of an array of boundary values,...
Definition RooBinning.h:27
RooCategory is an object to represent discrete states.
Definition RooCategory.h:28
bool defineType(const std::string &label)
Define a state with given name.
RooCmdArg is a named container for two doubles, two integers two object points and three string point...
Definition RooCmdArg.h:26
Class RooCmdConfig is a configurable parser for RooCmdArg named arguments.
void defineMutex(const char *head, Args_t &&... tail)
Define arguments where any pair is mutually exclusive.
bool process(const RooCmdArg &arg)
Process given RooCmdArg.
double getDouble(const char *name, double defaultValue=0.0) const
Return double property registered with name 'name'.
void defineDependency(const char *refArgName, const char *neededArgName)
Define that processing argument name refArgName requires processing of argument named neededArgName t...
bool defineDouble(const char *name, const char *argName, int doubleNum, double defValue=0.0)
Define double property name 'name' mapped to double in slot 'doubleNum' in RooCmdArg with name argNam...
RooArgSet * getSet(const char *name, RooArgSet *set=nullptr) const
Return RooArgSet property registered with name 'name'.
bool defineSet(const char *name, const char *argName, int setNum, const RooArgSet *set=nullptr)
Define TObject property name 'name' mapped to object in slot 'setNum' in RooCmdArg with name argName ...
bool ok(bool verbose) const
Return true of parsing was successful.
bool defineObject(const char *name, const char *argName, int setNum, const TObject *obj=nullptr, bool isArray=false)
Define TObject property name 'name' mapped to object in slot 'setNum' in RooCmdArg with name argName ...
const char * getString(const char *name, const char *defaultValue="", bool convEmptyToNull=false) const
Return string property registered with name 'name'.
bool defineString(const char *name, const char *argName, int stringNum, const char *defValue="", bool appendMode=false)
Define double property name 'name' mapped to double in slot 'stringNum' in RooCmdArg with name argNam...
const RooLinkedList & getObjectList(const char *name) const
Return list of objects registered with name 'name'.
bool defineInt(const char *name, const char *argName, int intNum, int defValue=0)
Define integer property name 'name' mapped to integer in slot 'intNum' in RooCmdArg with name argName...
int getInt(const char *name, int defaultValue=0) const
Return integer property registered with name 'name'.
TObject * getObject(const char *name, TObject *obj=nullptr) const
Return TObject property registered with name 'name'.
The RooDataHist is a container class to hold N-dimensional binned data.
Definition RooDataHist.h:39
std::span< const double > getWeightBatch(std::size_t first, std::size_t len, bool sumW2=false) const override
Return event weights of all events in range [first, first+len).
void interpolateQuadratic(double *output, std::span< const double > xVals, bool correctForBinSize, bool cdfBoundaries)
A vectorized version of interpolateDim for boundary safe quadratic interpolation of one dimensional h...
double sum(bool correctForBinSize, bool inverseCorr=false) const
Return the sum of the weights of all bins in the histogram.
void weights(double *output, std::span< double const > xVals, int intOrder, bool correctForBinSize, bool cdfBoundaries)
A vectorized version of RooDataHist::weight() for one dimensional histograms with up to one dimension...
Int_t _cache_sum_valid
! Is cache sum valid? Needs to be Int_t instead of CacheSumState_t for subclasses.
double interpolateDim(int iDim, double xval, size_t centralIdx, int intOrder, bool correctForBinSize, bool cdfBoundaries)
Perform boundary safe 'intOrder'-th interpolation of weights in dimension 'dim' at current value 'xva...
double weightSquared() const override
Return squared weight of last bin that was requested with get().
friend class RooDataHistSliceIter
void importTH1(const RooArgList &vars, const TH1 &histo, double initWgt, bool doDensityCorrection)
Import data from given TH1/2/3 into this RooDataHist.
void printDataHistogram(std::ostream &os, RooRealVar *obs) const
static TClass * Class()
TClass * IsA() const override
void SetNameTitle(const char *name, const char *title) override
Change the title of this RooDataHist.
double _cache_sum
! Cache for sum of entries ;
void initialize(const char *binningName=nullptr, bool fillTree=true)
Initialization procedure: allocate weights array, calculate multipliers needed for N-space to 1-dim a...
std::string calculateTreeIndexForCodeSquash(RooAbsArg const *klass, RooFit::Detail::CodeSquashContext &ctx, const RooAbsCollection &coords, bool reverse=false) const
VarInfo _varInfo
!
Int_t getIndex(const RooAbsCollection &coord, bool fast=false) const
Calculate bin number of the given coordinates.
void add(const RooArgSet &row, double wgt=1.0) override
Add wgt to the bin content enclosed by the coordinates passed in row.
Definition RooDataHist.h:66
const std::vector< double > & calculatePartialBinVolume(const RooArgSet &dimSet) const
Fill the transient cache with partial bin volumes with up-to-date values for the partial volume speci...
static std::unique_ptr< RooAbsDataStore > makeDefaultDataStore(RooStringView name, RooStringView title, RooArgSet const &vars)
double weightInterpolated(const RooArgSet &bin, int intOrder, bool correctForBinSize, bool cdfBoundaries)
Return the weight at given coordinates with interpolation.
std::unordered_map< int, std::vector< double > > _pbinvCache
! Cache for arrays of partial bin volumes
void checkBinBounds() const
void initializeAsymErrArrays() const
double weight(std::size_t i) const
Return weight of i-th bin.
void set(std::size_t binNumber, double weight, double wgtErr)
Set bin content of bin that was last loaded with get(std::size_t).
void weightError(double &lo, double &hi, ErrorType etype=Poisson) const override
Return the asymmetric errors on the current weight.
double * _errHi
[_arrSize] High-side error on weight array
void importTH1Set(const RooArgList &vars, RooCategory &indexCat, std::map< std::string, TH1 * > hmap, double initWgt, bool doDensityCorrection)
Import data from given set of TH1/2/3 into this RooDataHist.
void adjustBinning(const RooArgList &vars, const TH1 &href, Int_t *offset=nullptr)
Adjust binning specification on first and optionally second and third observable to binning in given ...
double * _binv
[_arrSize] Bin volume array
RooDataHist()
Default constructor.
ULong64_t _curIndex
Current index.
double weightFast(const RooArgSet &bin, int intOrder, bool correctForBinSize, bool cdfBoundaries)
A faster version of RooDataHist::weight that assumes the passed arguments are aligned with the histog...
RooPlot * plotOn(RooPlot *frame, PlotOpt o) const override
Back end function to plotting functionality.
double weight() const override
Return weight of last bin that was requested with get().
std::vector< std::vector< double > > _binbounds
! list of bin bounds per dimension
void printArgs(std::ostream &os) const override
Print argument of dataset, i.e. the observable names.
void importDHistSet(const RooArgList &vars, RooCategory &indexCat, std::map< std::string, RooDataHist * > dmap, double initWgt)
Import data from given set of TH1/2/3 into this RooDataHist.
void _adjustBinning(RooRealVar &theirVar, const TAxis &axis, RooRealVar *ourVar, Int_t *offset)
Helper doing the actual work of adjustBinning().
void printMultiline(std::ostream &os, Int_t content, bool verbose=false, TString indent="") const override
Print the details on the dataset contents.
double * _sumw2
[_arrSize] Sum of weights^2
TIterator * sliceIterator(RooAbsArg &sliceArg, const RooArgSet &otherArgs)
Create an iterator over all bins in a slice defined by the subset of observables listed in sliceArg.
Int_t calcTreeIndex() const
Legacy overload to calculate the tree index from the current value of _vars.
~RooDataHist() override
Destructor.
bool isNonPoissonWeighted() const override
Returns true if dataset contains entries with a non-integer weight.
std::string declWeightArrayForCodeSquash(RooAbsArg const *klass, RooFit::Detail::CodeSquashContext &ctx, bool correctForBinSize) const
std::vector< RooAbsLValue * > _lvvars
! List of observables casted as RooAbsLValue
void SetName(const char *name) override
Change the name of the RooDataHist.
std::vector< std::unique_ptr< const RooAbsBinning > > _lvbins
! List of used binnings associated with lvalues
void Streamer(TBuffer &) override
Stream an object of class RooDataHist.
std::vector< double > _interpolationBuffer
! Buffer to contain values used for weight interpolation
std::vector< Int_t > _idxMult
void registerWeightArraysToDataStore() const
Hand over pointers to our weight arrays to the data store implementation.
void reset() override
Reset all bin weights to zero.
double * _errLo
[_arrSize] Low-side error on weight array
double * _wgt
[_arrSize] Weight array
void printValue(std::ostream &os) const override
Print value of the dataset, i.e. the sum of weights contained in the dataset.
VarInfo const & getVarInfo()
Return reference to VarInfo struct with cached histogram variable information that is frequently used...
std::unique_ptr< RooAbsData > reduceEng(const RooArgSet &varSubset, const RooFormulaVar *cutVar, const char *cutRange=nullptr, std::size_t nStart=0, std::size_t nStop=std::numeric_limits< std::size_t >::max()) override
Implementation of RooAbsData virtual method that drives the RooAbsData::reduce() methods.
const RooArgSet * get() const override
Get bin centre of current bin.
Definition RooDataHist.h:76
void interpolateLinear(double *output, std::span< const double > xVals, bool correctForBinSize, bool cdfBoundaries)
A vectorized version of interpolateDim for boundary safe linear interpolation of one dimensional hist...
double binVolume() const
Return volume of current bin.
double sumEntries() const override
Sum the weights of all bins.
RooDirItem is a utility base class for RooFit objects that are to be attached to ROOT directories.
Definition RooDirItem.h:22
virtual void Streamer(TBuffer &)
void removeFromDir(TObject *obj)
Remove object from directory it was added to.
TDirectory * _dir
! Associated directory
Definition RooDirItem.h:33
void appendToDir(TObject *obj, bool forceMemoryResident=false)
Append object to directory.
A class to maintain the context for squashing of RooFit models into code.
std::string buildCall(std::string const &funcname, Args_t const &...args)
Build the code to call the function with name funcname, passing some arguments.
std::string makeValidVarName(TString in) const
Transform a string into a valid C++ variable name by replacing forbidden.
void addToCodeBody(RooAbsArg const *klass, std::string const &in)
Adds the input string to the squashed code body.
std::string getTmpVarName()
Get a unique variable name to be used in the generated code.
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
Return a confidence interval for the expected number of events given n observed (unweighted) events.
RooLinkedList is an collection class for internal use, storing a collection of RooAbsArg pointers in ...
RooLinkedListIterImpl begin() const
static double interpolate(double yArr[], Int_t nOrder, double x)
Definition RooMath.cxx:78
A RooPlot is a plot frame and a container for graphics objects within that frame.
Definition RooPlot.h:43
RooAbsRealLValue * getPlotVar() const
Definition RooPlot.h:137
virtual void printStream(std::ostream &os, Int_t contents, StyleOption style, TString indent="") const
Print description of object on ostream, printing contents set by contents integer,...
RooRealVar represents a variable that can be changed from the outside.
Definition RooRealVar.h:37
const RooAbsBinning & getBinning(const char *name=nullptr, bool verbose=true, bool createOnTheFly=false) const override
Return binning definition with name.
void setBinning(const RooAbsBinning &binning, const char *name=nullptr)
Add given binning under name 'name' with this variable.
void setRange(const char *name, double min, double max)
Set a fit or plotting range.
The RooStringView is a wrapper around a C-style string that can also be constructed from a std::strin...
RooUniformBinning is an implementation of RooAbsBinning that provides a uniform binning in 'n' bins b...
const Double_t * GetArray() const
Definition TArrayD.h:43
Class to manage histogram axis.
Definition TAxis.h:31
const TArrayD * GetXbins() const
Definition TAxis.h:136
Double_t GetXmax() const
Definition TAxis.h:140
Double_t GetXmin() const
Definition TAxis.h:139
Int_t GetNbins() const
Definition TAxis.h:125
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
virtual Int_t CheckByteCount(UInt_t startpos, UInt_t bcnt, const TClass *clss)=0
virtual void ReadFastArray(Bool_t *b, Int_t n)=0
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
static TClass * GetClass(const char *name, Bool_t load=kTRUE, Bool_t silent=kFALSE)
Static method returning pointer to TClass of the specified class name.
Definition TClass.cxx:2968
virtual TList * GetList() const
Definition TDirectory.h:222
TH1 is the base class of all histogram classes in ROOT.
Definition TH1.h:58
TAxis * GetZaxis()
Definition TH1.h:324
virtual Int_t GetNbinsY() const
Definition TH1.h:296
virtual Double_t GetBinError(Int_t bin) const
Return value of error associated to bin number bin.
Definition TH1.cxx:8932
virtual Int_t GetNbinsZ() const
Definition TH1.h:297
virtual Int_t GetDimension() const
Definition TH1.h:281
TAxis * GetXaxis()
Definition TH1.h:322
virtual Int_t GetNbinsX() const
Definition TH1.h:295
TAxis * GetYaxis()
Definition TH1.h:323
virtual Double_t GetBinContent(Int_t bin) const
Return content of bin number bin.
Definition TH1.cxx:5028
Iterator abstract base class.
Definition TIterator.h:30
void Add(TObject *obj) override
Definition TList.h:81
TObject * Remove(TObject *obj) override
Remove object from the list.
Definition TList.cxx:822
virtual void SetTitle(const char *title="")
Set the title of the TNamed.
Definition TNamed.cxx:164
const char * GetName() const override
Returns name of object.
Definition TNamed.h:47
const char * GetTitle() const override
Returns title of object.
Definition TNamed.h:48
Mother of all ROOT objects.
Definition TObject.h:41
virtual const char * ClassName() const
Returns name of class to which the object belongs.
Definition TObject.cxx:207
Basic string class.
Definition TString.h:139
virtual void Streamer(TBuffer &)
Stream a string object.
Definition TString.cxx:1391
A TTree represents a columnar dataset.
Definition TTree.h:79
const Int_t n
Definition legend1.C:16
TH1F * h1
Definition legend1.C:5
std::vector< std::string > Split(std::string_view str, std::string_view delims, bool skipEmpty=false)
Splits a string at each character in delims.
LongDouble_t Power(LongDouble_t x, LongDouble_t y)
Returns x raised to the power y.
Definition TMath.h:721
Definition first.py:1
RooAbsBinning * bins
Definition RooAbsData.h:181
Structure to cache information on the histogram variable that is frequently used for histogram weight...
TLine l
Definition textangle.C:4
static void output()