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