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