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