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*>(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 
706  // Save real dimensions of dataset separately
707  for (const auto real : _vars) {
708  if (dynamic_cast<RooAbsReal*>(real)) _realVars.add(*real);
709  }
710 
711  _lvvars.clear();
712  _lvbins.clear();
713 
714  // Fill array of LValue pointers to variables
715  for (unsigned int i = 0; i < _vars.size(); ++i) {
716  if (binningName) {
717  RooRealVar* rrv = dynamic_cast<RooRealVar*>(_vars[i]);
718  if (rrv) {
719  rrv->setBinning(rrv->getBinning(binningName));
720  }
721  }
722 
723  auto lvarg = dynamic_cast<RooAbsLValue*>(_vars[i]);
724  assert(lvarg);
725  _lvvars.push_back(lvarg);
726 
727  const RooAbsBinning* binning = lvarg->getBinningPtr(0);
728  _lvbins.emplace_back(binning ? binning->clone() : nullptr);
729  }
730 
731 
732  // Allocate coefficients array
733  _idxMult.resize(_vars.getSize()) ;
734 
735  _arrSize = 1 ;
736  unsigned int n = 0u;
737  for (const auto var : _vars) {
738  auto arg = dynamic_cast<const RooAbsLValue*>(var);
739  assert(arg);
740 
741  // Calculate sub-index multipliers for master index
742  for (unsigned int i = 0u; i<n; i++) {
743  _idxMult[i] *= arg->numBins() ;
744  }
745  _idxMult[n++] = 1 ;
746 
747  // Calculate dimension of weight array
748  _arrSize *= arg->numBins() ;
749  }
750 
751  // Allocate and initialize weight array if necessary
752  if (!_wgt) {
753  initArray(_wgt, _arrSize, 0.);
754  delete[] _errLo; _errLo = nullptr;
755  delete[] _errHi; _errHi = nullptr;
756  delete[] _sumw2; _sumw2 = nullptr;
757  initArray(_binv, _arrSize, 0.);
758 
759  // Refill array pointers in data store when reading
760  // from Streamer
761  if (!fillTree) {
763  }
764  }
765 
766  if (!fillTree) return ;
767 
768  // Fill TTree with bin center coordinates
769  // Calculate plot bins of components from master index
770 
771  for (Int_t ibin=0 ; ibin < _arrSize ; ibin++) {
772  Int_t j(0), idx(0), tmp(ibin) ;
773  Double_t theBinVolume(1) ;
774  for (auto arg2 : _lvvars) {
775  idx = tmp / _idxMult[j] ;
776  tmp -= idx*_idxMult[j++] ;
777  arg2->setBin(idx) ;
778  theBinVolume *= arg2->getBinWidth(idx) ;
779  }
780  _binv[ibin] = theBinVolume ;
781 
782  fill() ;
783  }
784 
785 
786 }
787 
788 
789 ////////////////////////////////////////////////////////////////////////////////
790 
792 {
793  if (!_binbounds.empty()) return;
794  for (auto& it : _lvbins) {
795  _binbounds.push_back(std::vector<Double_t>());
796  if (it) {
797  std::vector<Double_t>& bounds = _binbounds.back();
798  bounds.reserve(2 * it->numBins());
799  for (Int_t i = 0; i < it->numBins(); ++i) {
800  bounds.push_back(it->binLow(i));
801  bounds.push_back(it->binHigh(i));
802  }
803  }
804  }
805 }
806 
807 
808 ////////////////////////////////////////////////////////////////////////////////
809 /// Copy constructor
810 
811 RooDataHist::RooDataHist(const RooDataHist& other, const char* newname) :
812  RooAbsData(other,newname), RooDirItem(), _arrSize(other._arrSize), _idxMult(other._idxMult), _pbinvCache(other._pbinvCache)
813 {
814  // Allocate and initialize weight array
815  assert(_arrSize == other._arrSize);
816  cloneArray(_wgt, other._wgt, other._arrSize);
817  cloneArray(_errLo, other._errLo, other._arrSize);
818  cloneArray(_errHi, other._errHi, other._arrSize);
819  cloneArray(_binv, other._binv, other._arrSize);
820  cloneArray(_sumw2, other._sumw2, other._arrSize);
821 
822  // Save real dimensions of dataset separately
823  for (const auto arg : _vars) {
824  if (dynamic_cast<RooAbsReal*>(arg) != nullptr) _realVars.add(*arg) ;
825  }
826 
827  // Fill array of LValue pointers to variables
828  for (const auto rvarg : _vars) {
829  auto lvarg = dynamic_cast<RooAbsLValue*>(rvarg);
830  assert(lvarg);
831  _lvvars.push_back(lvarg);
832  const RooAbsBinning* binning = lvarg->getBinningPtr(0);
833  _lvbins.emplace_back(binning ? binning->clone() : 0) ;
834  }
835 
837 
838  appendToDir(this,kTRUE) ;
839 }
840 
841 
842 
843 ////////////////////////////////////////////////////////////////////////////////
844 /// Constructor of a data hist from (part of) an existing data hist. The dimensions
845 /// of the data set are defined by the 'vars' RooArgSet, which can be identical
846 /// to 'dset' dimensions, or a subset thereof. Reduced dimensions will be projected
847 /// in the output data hist. The optional 'cutVar' formula variable can used to
848 /// select the subset of bins to be copied.
849 ///
850 /// For most uses the RooAbsData::reduce() wrapper function, which uses this constructor,
851 /// is the most convenient way to create a subset of an existing data
852 
853 RooDataHist::RooDataHist(const char* name, const char* title, RooDataHist* h, const RooArgSet& varSubset,
854  const RooFormulaVar* cutVar, const char* cutRange, Int_t nStart, Int_t nStop, Bool_t copyCache) :
855  RooAbsData(name,title,varSubset)
856 {
857  // Initialize datastore
858  _dstore = new RooTreeDataStore(name,title,*h->_dstore,_vars,cutVar,cutRange,nStart,nStop,copyCache) ;
859 
860  initialize(0,kFALSE) ;
861 
862  // Copy weight array etc
863  assert(_arrSize == h->_arrSize);
864  cloneArray(_wgt, h->_wgt, _arrSize);
865  cloneArray(_errLo, h->_errLo, _arrSize);
866  cloneArray(_errHi, h->_errHi, _arrSize);
867  cloneArray(_binv, h->_binv, _arrSize);
868  cloneArray(_sumw2, h->_sumw2, _arrSize);
869 
871 
872  appendToDir(this,kTRUE) ;
874 }
875 
876 
877 ////////////////////////////////////////////////////////////////////////////////
878 /// Construct a clone of this dataset that contains only the cached variables
879 
880 RooAbsData* RooDataHist::cacheClone(const RooAbsArg* newCacheOwner, const RooArgSet* newCacheVars, const char* newName)
881 {
882  checkInit() ;
883 
884  RooDataHist* dhist = new RooDataHist(newName?newName:GetName(),GetTitle(),this,*get(),0,0,0,2000000000,kTRUE) ;
885 
886  RooArgSet* selCacheVars = (RooArgSet*) newCacheVars->selectCommon(dhist->_cachedVars) ;
887  dhist->attachCache(newCacheOwner, *selCacheVars) ;
888  delete selCacheVars ;
889 
890  return dhist ;
891 }
892 
893 
894 
895 ////////////////////////////////////////////////////////////////////////////////
896 /// Implementation of RooAbsData virtual method that drives the RooAbsData::reduce() methods
897 
898 RooAbsData* RooDataHist::reduceEng(const RooArgSet& varSubset, const RooFormulaVar* cutVar, const char* cutRange,
899  std::size_t nStart, std::size_t nStop, Bool_t /*copyCache*/)
900 {
901  checkInit() ;
902  RooArgSet* myVarSubset = (RooArgSet*) _vars.selectCommon(varSubset) ;
903  RooDataHist *rdh = new RooDataHist(GetName(), GetTitle(), *myVarSubset) ;
904  delete myVarSubset ;
905 
906  RooFormulaVar* cloneVar = 0;
907  RooArgSet* tmp(0) ;
908  if (cutVar) {
909  // Deep clone cutVar and attach clone to this dataset
910  tmp = (RooArgSet*) RooArgSet(*cutVar).snapshot() ;
911  if (!tmp) {
912  coutE(DataHandling) << "RooDataHist::reduceEng(" << GetName() << ") Couldn't deep-clone cut variable, abort," << endl ;
913  return 0 ;
914  }
915  cloneVar = (RooFormulaVar*) tmp->find(*cutVar) ;
916  cloneVar->attachDataSet(*this) ;
917  }
918 
919  Double_t lo,hi ;
920  const std::size_t nevt = nStop < static_cast<std::size_t>(numEntries()) ? nStop : static_cast<std::size_t>(numEntries());
921  for (auto i=nStart; i<nevt ; i++) {
922  const RooArgSet* row = get(i) ;
923 
924  Bool_t doSelect(kTRUE) ;
925  if (cutRange) {
926  for (const auto arg : *row) {
927  if (!arg->inRange(cutRange)) {
928  doSelect = kFALSE ;
929  break ;
930  }
931  }
932  }
933  if (!doSelect) continue ;
934 
935  if (!cloneVar || cloneVar->getVal()) {
936  weightError(lo,hi,SumW2) ;
937  rdh->add(*row,weight(),lo*lo) ;
938  }
939  }
940 
941  if (cloneVar) {
942  delete tmp ;
943  }
944 
945  return rdh ;
946 }
947 
948 
949 
950 ////////////////////////////////////////////////////////////////////////////////
951 /// Destructor
952 
954 {
955  delete[] _wgt;
956  delete[] _errLo;
957  delete[] _errHi;
958  delete[] _sumw2;
959  delete[] _binv;
960 
961  removeFromDir(this) ;
963 }
964 
965 
966 
967 
968 ////////////////////////////////////////////////////////////////////////////////
969 /// Calculate bin number of the given coordinates. If only a subset of the internal
970 /// coordinates are passed, the missing coordinates are taken at their current value.
971 /// \param[in] coord Variables that are representing the coordinates.
972 /// \param[in] fast If the variables in `coord` and the ones of the data hist have the
973 /// same size and layout, `fast` can be set to skip checking that all variables are
974 /// present in `coord`.
976  checkInit() ;
977  return calcTreeIndex(coord, fast);
978 }
979 
980 
981 
982 
983 ////////////////////////////////////////////////////////////////////////////////
984 /// Calculate the bin index corresponding to the coordinates passed as argument.
985 /// \param[in] coords Coordinates. If `fast == false`, these can be partial.
986 /// \param[in] fast Promise that the coordinates in `coords` have the same order
987 /// as the internal coordinates. In this case, values are looked up only by index.
988 std::size_t RooDataHist::calcTreeIndex(const RooAbsCollection& coords, bool fast) const
989 {
990  // With fast, caller promises that layout of "coords" is identical to our internal "vars"
991  assert(!fast || _vars.size() == coords.size());
992 
993  if (&_vars == &coords)
994  fast = true;
995 
996  std::size_t masterIdx = 0;
997 
998  for (unsigned int i=0; i < _vars.size(); ++i) {
999  const RooAbsArg* internalVar = _vars[i];
1000  const RooAbsBinning* binning = _lvbins[i].get();
1001 
1002  // Find the variable that we need values from.
1003  // That's either the variable directly from the external coordinates
1004  // or we find the external one that has the same name as "internalVar".
1005  const RooAbsArg* theVar = fast ? coords[i] : coords.find(*internalVar);
1006  if (!theVar) {
1007  // Variable is not in external coordinates. Use current internal value.
1008  theVar = internalVar;
1009  }
1010  // If fast is on, users promise that the sets have the same layout
1011  assert(!fast || strcmp(internalVar->GetName(), theVar->GetName()) == 0);
1012 
1013  if (binning) {
1014  assert(dynamic_cast<const RooAbsReal*>(theVar));
1015  const double val = static_cast<const RooAbsReal*>(theVar)->getVal();
1016  masterIdx += _idxMult[i] * binning->binNumber(val);
1017  } else {
1018  // We are a category. No binning.
1019  assert(dynamic_cast<const RooAbsCategoryLValue*>(theVar));
1020  auto cat = static_cast<const RooAbsCategoryLValue*>(theVar);
1021  masterIdx += _idxMult[i] * cat->getBin(static_cast<const char*>(nullptr));
1022  }
1023  }
1024 
1025  return masterIdx ;
1026 }
1027 
1028 
1029 
1030 ////////////////////////////////////////////////////////////////////////////////
1031 /// Debug stuff, should go...
1032 
1034 {
1035  cout << "_arrSize = " << _arrSize << endl ;
1036  for (Int_t i=0 ; i < _arrSize ; i++) {
1037  cout << "wgt[" << i << "] = " << _wgt[i]
1038  << "\tsumw2[" << i << "] = " << (_sumw2 ? _sumw2[i] : -1.)
1039  << "\tvol[" << i << "] = " << _binv[i] << endl ;
1040  }
1041 }
1042 
1043 
1044 
1045 ////////////////////////////////////////////////////////////////////////////////
1046 /// Back end function to plotting functionality. Plot RooDataHist on given
1047 /// frame in mode specified by plot options 'o'. The main purpose of
1048 /// this function is to match the specified binning on 'o' to the
1049 /// internal binning of the plot observable in this RooDataHist.
1050 /// \see RooAbsData::plotOn() for plotting options.
1052 {
1053  checkInit() ;
1054  if (o.bins) return RooAbsData::plotOn(frame,o) ;
1055 
1056  if(0 == frame) {
1057  coutE(InputArguments) << ClassName() << "::" << GetName() << ":plotOn: frame is null" << endl;
1058  return 0;
1059  }
1060  RooAbsRealLValue *var= (RooAbsRealLValue*) frame->getPlotVar();
1061  if(0 == var) {
1062  coutE(InputArguments) << ClassName() << "::" << GetName()
1063  << ":plotOn: frame does not specify a plot variable" << endl;
1064  return 0;
1065  }
1066 
1067  RooRealVar* dataVar = (RooRealVar*) _vars.find(*var) ;
1068  if (!dataVar) {
1069  coutE(InputArguments) << ClassName() << "::" << GetName()
1070  << ":plotOn: dataset doesn't contain plot frame variable" << endl;
1071  return 0;
1072  }
1073 
1074  o.bins = &dataVar->getBinning() ;
1076  return RooAbsData::plotOn(frame,o) ;
1077 }
1078 
1079 
1080 ////////////////////////////////////////////////////////////////////////////////
1081 /// Return the weight at given coordinates with optional interpolation.
1082 /// \param[in] bin Coordinates for which the weight should be calculated.
1083 /// \param[in] intOrder If zero, the bare weight for the bin enclosing the coordinatesis returned.
1084 /// For higher values, the result is interpolated in the real dimensions of the dataset.
1085 /// \param[in] correctForBinSize
1086 /// \param[in] cdfBoundaries
1087 /// \param[in] oneSafe Ignored.
1088 
1089 Double_t RooDataHist::weight(const RooArgSet& bin, Int_t intOrder, Bool_t correctForBinSize, Bool_t cdfBoundaries, Bool_t /*oneSafe*/)
1090 {
1091  checkInit() ;
1092 
1093  // Handle illegal intOrder values
1094  if (intOrder<0) {
1095  coutE(InputArguments) << "RooDataHist::weight(" << GetName() << ") ERROR: interpolation order must be positive" << endl ;
1096  return 0 ;
1097  }
1098 
1099  // Handle no-interpolation case
1100  if (intOrder==0) {
1101  const auto idx = calcTreeIndex(bin, false);
1102  if (correctForBinSize) {
1103  return get_wgt(idx) / _binv[idx];
1104  } else {
1105  return get_wgt(idx);
1106  }
1107  }
1108 
1109  // Handle all interpolation cases
1110  _vars.assignValueOnly(bin) ;
1111 
1112  Double_t wInt(0) ;
1113  if (_realVars.getSize()==1) {
1114 
1115  // 1-dimensional interpolation
1116  const auto real = static_cast<RooRealVar*>(_realVars[static_cast<std::size_t>(0)]);
1117  const RooAbsBinning* binning = real->getBinningPtr(0) ;
1118  wInt = interpolateDim(*real,binning,((RooAbsReal*)bin.find(*real))->getVal(), intOrder, correctForBinSize, cdfBoundaries) ;
1119 
1120  } else if (_realVars.getSize()==2) {
1121 
1122  // 2-dimensional interpolation
1123  const auto realX = static_cast<RooRealVar*>(_realVars[static_cast<std::size_t>(0)]);
1124  const auto realY = static_cast<RooRealVar*>(_realVars[static_cast<std::size_t>(1)]);
1125  Double_t xval = ((RooAbsReal*)bin.find(*realX))->getVal() ;
1126  Double_t yval = ((RooAbsReal*)bin.find(*realY))->getVal() ;
1127 
1128  Int_t ybinC = realY->getBin() ;
1129  Int_t ybinLo = ybinC-intOrder/2 - ((yval<realY->getBinning().binCenter(ybinC))?1:0) ;
1130  Int_t ybinM = realY->numBins() ;
1131 
1132  Int_t i ;
1133  Double_t yarr[10] ;
1134  Double_t xarr[10] ;
1135  const RooAbsBinning* binning = realX->getBinningPtr(0) ;
1136  for (i=ybinLo ; i<=intOrder+ybinLo ; i++) {
1137  Int_t ibin ;
1138  if (i>=0 && i<ybinM) {
1139  // In range
1140  ibin = i ;
1141  realY->setBin(ibin) ;
1142  xarr[i-ybinLo] = realY->getVal() ;
1143  } else if (i>=ybinM) {
1144  // Overflow: mirror
1145  ibin = 2*ybinM-i-1 ;
1146  realY->setBin(ibin) ;
1147  xarr[i-ybinLo] = 2*realY->getMax()-realY->getVal() ;
1148  } else {
1149  // Underflow: mirror
1150  ibin = -i -1;
1151  realY->setBin(ibin) ;
1152  xarr[i-ybinLo] = 2*realY->getMin()-realY->getVal() ;
1153  }
1154  yarr[i-ybinLo] = interpolateDim(*realX,binning,xval,intOrder,correctForBinSize,kFALSE) ;
1155  }
1156 
1157  if (gDebug>7) {
1158  cout << "RooDataHist interpolating data is" << endl ;
1159  cout << "xarr = " ;
1160  for (int q=0; q<=intOrder ; q++) cout << xarr[q] << " " ;
1161  cout << " yarr = " ;
1162  for (int q=0; q<=intOrder ; q++) cout << yarr[q] << " " ;
1163  cout << endl ;
1164  }
1165  wInt = RooMath::interpolate(xarr,yarr,intOrder+1,yval) ;
1166 
1167  } else {
1168 
1169  // Higher dimensional scenarios not yet implemented
1170  coutE(InputArguments) << "RooDataHist::weight(" << GetName() << ") interpolation in "
1171  << _realVars.getSize() << " dimensions not yet implemented" << endl ;
1172  return weight(bin,0) ;
1173 
1174  }
1175 
1176  return wInt ;
1177 }
1178 
1179 
1180 
1181 
1182 ////////////////////////////////////////////////////////////////////////////////
1183 /// Return the error of current weight.
1184 /// \param[out] lo Low error.
1185 /// \param[out] hi High error.
1186 /// \param[in] etype Type of error to compute. May throw if not supported.
1187 /// Supported errors are
1188 /// - `Poisson` Default. Asymmetric Poisson errors (68% CL).
1189 /// - `SumW2` The square root of the sum of weights. (Symmetric).
1190 /// - `None` Return zero.
1192 {
1193  checkInit() ;
1194 
1195  switch (etype) {
1196 
1197  case Auto:
1198  throw std::invalid_argument(Form("RooDataHist::weightError(%s) error type Auto not allowed here",GetName())) ;
1199  break ;
1200 
1201  case Expected:
1202  throw std::invalid_argument(Form("RooDataHist::weightError(%s) error type Expected not allowed here",GetName())) ;
1203  break ;
1204 
1205  case Poisson:
1206  if (get_curWgtErrLo() >= 0) {
1207  // Weight is preset or precalculated
1208  lo = get_curWgtErrLo();
1209  hi = get_curWgtErrHi();
1210  return ;
1211  }
1212 
1213  if (!_errLo || !_errHi) {
1214  // We didn't track asymmetric errors so far, so now we need to allocate
1215  initArray(_errLo, _arrSize, -1.);
1216  initArray(_errHi, _arrSize, -1.);
1218  }
1219 
1220  // Calculate poisson errors
1221  Double_t ym,yp ;
1223  _errLo[_curIndex] = weight()-ym;
1224  _errHi[_curIndex] = yp-weight();
1225  lo = _errLo[_curIndex];
1226  hi = _errHi[_curIndex];
1227  return ;
1228 
1229  case SumW2:
1230  lo = sqrt(get_curSumW2());
1231  hi = lo;
1232  return ;
1233 
1234  case None:
1235  lo = 0 ;
1236  hi = 0 ;
1237  return ;
1238  }
1239 }
1240 
1241 
1242 // wve adjust for variable bin sizes
1243 
1244 ////////////////////////////////////////////////////////////////////////////////
1245 /// Perform boundary safe 'intOrder'-th interpolation of weights in dimension 'dim'
1246 /// at current value 'xval'
1247 
1248 Double_t RooDataHist::interpolateDim(RooRealVar& dim, const RooAbsBinning* binning, Double_t xval, Int_t intOrder, Bool_t correctForBinSize, Bool_t cdfBoundaries)
1249 {
1250  // Fill workspace arrays spanning interpolation area
1251  Int_t fbinC = dim.getBin(*binning) ;
1252  Int_t fbinLo = fbinC-intOrder/2 - ((xval<binning->binCenter(fbinC))?1:0) ;
1253  Int_t fbinM = dim.numBins(*binning) ;
1254 
1255 
1256  Int_t i ;
1257  Double_t yarr[10] ;
1258  Double_t xarr[10] ;
1259  for (i=fbinLo ; i<=intOrder+fbinLo ; i++) {
1260  Int_t ibin ;
1261  if (i>=0 && i<fbinM) {
1262  // In range
1263  ibin = i ;
1264  dim.setBinFast(ibin,*binning) ;
1265  //cout << "INRANGE: dim.getVal(ibin=" << ibin << ") = " << dim.getVal() << endl ;
1266  xarr[i-fbinLo] = dim.getVal() ;
1267  Int_t idx = calcTreeIndex(_vars, true);
1268  yarr[i - fbinLo] = get_wgt(idx);
1269  if (correctForBinSize) yarr[i-fbinLo] /= _binv[idx] ;
1270  } else if (i>=fbinM) {
1271  // Overflow: mirror
1272  ibin = 2*fbinM-i-1 ;
1273  dim.setBinFast(ibin,*binning) ;
1274  //cout << "OVERFLOW: dim.getVal(ibin=" << ibin << ") = " << dim.getVal() << endl ;
1275  if (cdfBoundaries) {
1276  xarr[i-fbinLo] = dim.getMax()+1e-10*(i-fbinM+1) ;
1277  yarr[i-fbinLo] = 1.0 ;
1278  } else {
1279  Int_t idx = calcTreeIndex(_vars, true) ;
1280  xarr[i-fbinLo] = 2*dim.getMax()-dim.getVal() ;
1281  yarr[i - fbinLo] = get_wgt(idx);
1282  if (correctForBinSize)
1283  yarr[i - fbinLo] /= _binv[idx];
1284  }
1285  } else {
1286  // Underflow: mirror
1287  ibin = -i - 1 ;
1288  dim.setBinFast(ibin,*binning) ;
1289  //cout << "UNDERFLOW: dim.getVal(ibin=" << ibin << ") = " << dim.getVal() << endl ;
1290  if (cdfBoundaries) {
1291  xarr[i-fbinLo] = dim.getMin()-ibin*(1e-10) ; ;
1292  yarr[i-fbinLo] = 0.0 ;
1293  } else {
1294  Int_t idx = calcTreeIndex(_vars, true) ;
1295  xarr[i-fbinLo] = 2*dim.getMin()-dim.getVal() ;
1296  yarr[i - fbinLo] = get_wgt(idx);
1297  if (correctForBinSize)
1298  yarr[i - fbinLo] /= _binv[idx];
1299  }
1300  }
1301  //cout << "ibin = " << ibin << endl ;
1302  }
1303 // for (int k=0 ; k<=intOrder ; k++) {
1304 // cout << "k=" << k << " x = " << xarr[k] << " y = " << yarr[k] << endl ;
1305 // }
1306  dim.setBinFast(fbinC,*binning) ;
1307  Double_t ret = RooMath::interpolate(xarr,yarr,intOrder+1,xval) ;
1308  return ret ;
1309 }
1310 
1311 
1312 
1313 
1314 ////////////////////////////////////////////////////////////////////////////////
1315 /// Increment the bin content of the bin enclosing the given coordinates.
1316 ///
1317 /// \param[in] row Coordinates of the bin.
1318 /// \param[in] wgt Increment by this weight.
1319 /// \param[in] sumw2 Optionally, track the sum of squared weights. If a value > 0 or
1320 /// a weight != 1. is passed for the first time, a vector for the squared weights will be allocated.
1321 void RooDataHist::add(const RooArgSet& row, Double_t wgt, Double_t sumw2)
1322 {
1323  checkInit() ;
1324 
1325  if ((sumw2 > 0. || wgt != 1.) && !_sumw2) {
1326  // Receiving a weighted entry. SumW2 != sumw from now on.
1327  _sumw2 = new double[_arrSize];
1328  std::copy(_wgt, _wgt+_arrSize, _sumw2);
1329 
1331  }
1332 
1333  const auto idx = calcTreeIndex(row, false);
1334 
1335  _wgt[idx] += wgt ;
1336  if (_sumw2) _sumw2[idx] += (sumw2 > 0 ? sumw2 : wgt*wgt);
1337 
1338  _cache_sum_valid = false;
1339 }
1340 
1341 
1342 
1343 ////////////////////////////////////////////////////////////////////////////////
1344 /// Set a bin content.
1345 /// \param[in] row Coordinates of the bin to be set.
1346 /// \param[in] wgt New bin content.
1347 /// \param[in] wgtErrLo Low error of the bin content.
1348 /// \param[in] wgtErrHi High error of the bin content.
1349 void RooDataHist::set(const RooArgSet& row, Double_t wgt, Double_t wgtErrLo, Double_t wgtErrHi)
1350 {
1351  checkInit() ;
1352 
1353  if (!_errLo || !_errHi) {
1354  initArray(_errLo, _arrSize, -1.);
1355  initArray(_errHi, _arrSize, -1.);
1357  }
1358 
1359  const auto idx = calcTreeIndex(row, false);
1360 
1361  _wgt[idx] = wgt ;
1362  _errLo[idx] = wgtErrLo ;
1363  _errHi[idx] = wgtErrHi ;
1364 
1365  _cache_sum_valid = false;
1366 }
1367 
1368 
1369 
1370 ////////////////////////////////////////////////////////////////////////////////
1371 /// Set bin content of bin that was last loaded with get(std::size_t).
1372 /// \param[in] binNumber Optional bin number to set. If empty, currently active bin is set.
1373 /// \param[in] wgt New bin content.
1374 /// \param[in] wgtErr Error of the new bin content. If the weight need not have an error, use 0. or a negative number.
1375 void RooDataHist::set(std::size_t binNumber, double wgt, double wgtErr) {
1376  checkInit() ;
1377 
1378  if (wgtErr > 0. && !_sumw2) {
1379  // Receiving a weighted entry. Need to track sumw2 from now on:
1380  cloneArray(_sumw2, _wgt, _arrSize);
1381 
1383  }
1384 
1385  _wgt[binNumber] = wgt ;
1386  if (_errLo) _errLo[binNumber] = wgtErr;
1387  if (_errHi) _errHi[binNumber] = wgtErr;
1388  if (_sumw2) _sumw2[binNumber] = wgtErr*wgtErr;
1389 
1391 }
1392 
1393 
1394 ////////////////////////////////////////////////////////////////////////////////
1395 /// Set bin content of bin that was last loaded with get(std::size_t).
1396 /// \deprecated Prefer set(std::size_t, double, double).
1397 /// \param[in] wgt New bin content.
1398 /// \param[in] wgtErr Optional error of the bin content.
1399 void RooDataHist::set(double wgt, double wgtErr) {
1400  if (_curIndex == std::numeric_limits<std::size_t>::max()) {
1401  _curIndex = calcTreeIndex(_vars, true) ;
1402  }
1403 
1404  set(_curIndex, wgt, wgtErr);
1405 }
1406 
1407 
1408 ////////////////////////////////////////////////////////////////////////////////
1409 /// Set a bin content.
1410 /// \param[in] row Coordinates to compute the bin from.
1411 /// \param[in] wgt New bin content.
1412 /// \param[in] wgtErr Optional error of the bin content.
1413 void RooDataHist::set(const RooArgSet& row, Double_t wgt, Double_t wgtErr) {
1414  set(calcTreeIndex(row, false), wgt, wgtErr);
1415 }
1416 
1417 
1418 
1419 ////////////////////////////////////////////////////////////////////////////////
1420 /// Add all data points contained in 'dset' to this data set with given weight.
1421 /// Optional cut string expression selects the data points to be added and can
1422 /// reference any variable contained in this data set
1423 
1424 void RooDataHist::add(const RooAbsData& dset, const char* cut, Double_t wgt)
1425 {
1426  RooFormulaVar cutVar("select",cut,*dset.get()) ;
1427  add(dset,&cutVar,wgt) ;
1428 }
1429 
1430 
1431 
1432 ////////////////////////////////////////////////////////////////////////////////
1433 /// Add all data points contained in 'dset' to this data set with given weight.
1434 /// Optional RooFormulaVar pointer selects the data points to be added.
1435 
1436 void RooDataHist::add(const RooAbsData& dset, const RooFormulaVar* cutVar, Double_t wgt)
1437 {
1438  checkInit() ;
1439 
1440  RooFormulaVar* cloneVar = 0;
1441  RooArgSet* tmp(0) ;
1442  if (cutVar) {
1443  // Deep clone cutVar and attach clone to this dataset
1444  tmp = (RooArgSet*) RooArgSet(*cutVar).snapshot() ;
1445  if (!tmp) {
1446  coutE(DataHandling) << "RooDataHist::add(" << GetName() << ") Couldn't deep-clone cut variable, abort," << endl ;
1447  return ;
1448  }
1449 
1450  cloneVar = (RooFormulaVar*) tmp->find(*cutVar) ;
1451  cloneVar->attachDataSet(dset) ;
1452  }
1453 
1454 
1455  Int_t i ;
1456  for (i=0 ; i<dset.numEntries() ; i++) {
1457  const RooArgSet* row = dset.get(i) ;
1458  if (!cloneVar || cloneVar->getVal()) {
1459  add(*row,wgt*dset.weight(), wgt*wgt*dset.weightSquared()) ;
1460  }
1461  }
1462 
1463  if (cloneVar) {
1464  delete tmp ;
1465  }
1466 
1468 }
1469 
1470 
1471 
1472 ////////////////////////////////////////////////////////////////////////////////
1473 /// Return the sum of the weights of all bins in the histogram.
1474 ///
1475 /// \param[in] correctForBinSize Multiply the sum of weights in each bin
1476 /// with the N-dimensional bin volume, making the return value
1477 /// the integral over the function represented by this histogram.
1478 /// \param[in] inverseBinCor Divide by the N-dimensional bin volume.
1479 Double_t RooDataHist::sum(bool correctForBinSize, bool inverseBinCor) const
1480 {
1481  checkInit() ;
1482 
1483  // Check if result was cached
1484  const CacheSumState_t cache_code = !correctForBinSize ? kNoBinCorrection : (inverseBinCor ? kInverseBinCorr : kCorrectForBinSize);
1485  if (_cache_sum_valid == static_cast<Int_t>(cache_code)) {
1486  return _cache_sum ;
1487  }
1488 
1490  for (Int_t i=0; i < _arrSize; i++) {
1491  const double theBinVolume = correctForBinSize ? (inverseBinCor ? 1/_binv[i] : _binv[i]) : 1.0 ;
1492  kahanSum += get_wgt(i) * theBinVolume;
1493  }
1494 
1495  // Store result in cache
1496  _cache_sum_valid = cache_code;
1497  _cache_sum = kahanSum;
1498 
1499  return kahanSum;
1500 }
1501 
1502 
1503 
1504 ////////////////////////////////////////////////////////////////////////////////
1505 /// Return the sum of the weights of a multi-dimensional slice of the histogram
1506 /// by summing only over the dimensions specified in sumSet.
1507 ///
1508 /// The coordinates of all other dimensions are fixed to those given in sliceSet
1509 ///
1510 /// If correctForBinSize is specified, the sum of weights
1511 /// is multiplied by the M-dimensional bin volume, (M = N(sumSet)),
1512 /// making the return value the integral over the function
1513 /// represented by this histogram
1514 
1515 Double_t RooDataHist::sum(const RooArgSet& sumSet, const RooArgSet& sliceSet, bool correctForBinSize, bool inverseBinCor)
1516 {
1517  checkInit() ;
1518 
1519  RooArgSet varSave ;
1520  varSave.addClone(_vars) ;
1521 
1522  RooArgSet sliceOnlySet(sliceSet);
1523  sliceOnlySet.remove(sumSet,true,true) ;
1524 
1525  _vars = sliceOnlySet;
1526  std::vector<double> const * pbinv = nullptr;
1527 
1528  if(correctForBinSize && inverseBinCor) {
1529  pbinv = &calculatePartialBinVolume(sliceOnlySet);
1530  } else if(correctForBinSize && !inverseBinCor) {
1531  pbinv = &calculatePartialBinVolume(sumSet);
1532  }
1533 
1534  // Calculate mask and refence plot bins for non-iterating variables
1535  std::vector<bool> mask(_vars.getSize());
1536  std::vector<int> refBin(_vars.getSize());
1537 
1538  for (unsigned int i = 0; i < _vars.size(); ++i) {
1539  const RooAbsArg* arg = _vars[i];
1540  const RooAbsLValue* argLv = _lvvars[i]; // Same as above, but cross-cast
1541 
1542  if (sumSet.find(*arg)) {
1543  mask[i] = false ;
1544  } else {
1545  mask[i] = true ;
1546  refBin[i] = argLv->getBin();
1547  }
1548  }
1549 
1550  // Loop over entire data set, skipping masked entries
1552  for (Int_t ibin=0; ibin < _arrSize; ++ibin) {
1553 
1554  std::size_t tmpibin = ibin;
1555  Bool_t skip(false) ;
1556 
1557  // Check if this bin belongs in selected slice
1558  for (unsigned int ivar = 0; !skip && ivar < _vars.size(); ++ivar) {
1559  const Int_t idx = tmpibin / _idxMult[ivar] ;
1560  tmpibin -= idx*_idxMult[ivar] ;
1561  if (mask[ivar] && idx!=refBin[ivar])
1562  skip = true ;
1563  }
1564 
1565  if (!skip) {
1566  const double theBinVolume = correctForBinSize ? (inverseBinCor ? 1/(*pbinv)[ibin] : (*pbinv)[ibin] ) : 1.0 ;
1567  total += get_wgt(ibin) * theBinVolume;
1568  }
1569  }
1570 
1571  _vars = varSave ;
1572 
1573  return total;
1574 }
1575 
1576 ////////////////////////////////////////////////////////////////////////////////
1577 /// Return the sum of the weights of a multi-dimensional slice of the histogram
1578 /// by summing only over the dimensions specified in sumSet.
1579 ///
1580 /// The coordinates of all other dimensions are fixed to those given in sliceSet
1581 ///
1582 /// If correctForBinSize is specified, the sum of weights
1583 /// is multiplied by the M-dimensional bin volume, (M = N(sumSet)),
1584 /// or the fraction of it that falls inside the range rangeName,
1585 /// making the return value the integral over the function
1586 /// represented by this histogram.
1587 ///
1588 /// If correctForBinSize is not specified, the weights are multiplied by the
1589 /// fraction of the bin volume that falls inside the range, i.e. a factor of
1590 /// binVolumeInRange/totalBinVolume.
1591 
1592 Double_t RooDataHist::sum(const RooArgSet& sumSet, const RooArgSet& sliceSet,
1593  bool correctForBinSize, bool inverseBinCor,
1594  const std::map<const RooAbsArg*, std::pair<double, double> >& ranges,
1595  std::function<double(int)> getBinScale)
1596 {
1597  checkInit();
1598  checkBinBounds();
1599  RooArgSet varSave;
1600  varSave.addClone(_vars);
1601  {
1602  RooArgSet sliceOnlySet(sliceSet);
1603  sliceOnlySet.remove(sumSet, true, true);
1604  _vars = sliceOnlySet;
1605  }
1606 
1607  // Calculate mask and reference plot bins for non-iterating variables,
1608  // and get ranges for iterating variables
1609  std::vector<bool> mask(_vars.getSize());
1610  std::vector<int> refBin(_vars.getSize());
1611  std::vector<double> rangeLo(_vars.getSize(), -std::numeric_limits<Double_t>::infinity());
1612  std::vector<double> rangeHi(_vars.getSize(), +std::numeric_limits<Double_t>::infinity());
1613 
1614  for (std::size_t i = 0; i < _vars.size(); ++i) {
1615  const RooAbsArg* arg = _vars[i];
1616  const RooAbsLValue* argLV = _lvvars[i]; // Same object as above, but cross cast
1617 
1618  RooAbsArg* sumsetv = sumSet.find(*arg);
1619  RooAbsArg* slicesetv = sliceSet.find(*arg);
1620  mask[i] = !sumsetv;
1621  if (mask[i]) {
1622  assert(argLV);
1623  refBin[i] = argLV->getBin();
1624  }
1625 
1626  auto it = ranges.find(sumsetv ? sumsetv : slicesetv);
1627  if (ranges.end() != it) {
1628  rangeLo[i] = it->second.first;
1629  rangeHi[i] = it->second.second;
1630  }
1631  }
1632 
1633  // Loop over entire data set, skipping masked entries
1635  for (Int_t ibin = 0; ibin < _arrSize; ++ibin) {
1636  // Check if this bin belongs in selected slice
1637  bool skip{false};
1638  for (int ivar = 0, tmp = ibin; !skip && ivar < int(_vars.size()); ++ivar) {
1639  const Int_t idx = tmp / _idxMult[ivar];
1640  tmp -= idx*_idxMult[ivar];
1641  if (mask[ivar] && idx!=refBin[ivar]) skip = true;
1642  }
1643 
1644  if (skip) continue;
1645 
1646  // Work out bin volume
1647  // It's not necessary to figure out the bin volume for the slice-only set explicitely here.
1648  // We need to loop over the sumSet anyway to get the partial bin containment correction,
1649  // so we can get the slice-only set volume later by dividing _binv[ibin] / binVolumeSumSetFull.
1650  Double_t binVolumeSumSetFull = 1.;
1651  Double_t binVolumeSumSetInRange = 1.;
1652  for (Int_t ivar = 0, tmp = ibin; ivar < (int)_vars.size(); ++ivar) {
1653  const Int_t idx = tmp / _idxMult[ivar];
1654  tmp -= idx*_idxMult[ivar];
1655 
1656  // If the current variable is not in the sumSet, it should not be considered for the bin volume
1657  const auto arg = _vars[ivar];
1658  if (!sumSet.find(*arg)) {
1659  continue;
1660  }
1661 
1662  if (_binbounds[ivar].empty()) continue;
1663  const Double_t binLo = _binbounds[ivar][2 * idx];
1664  const Double_t binHi = _binbounds[ivar][2 * idx + 1];
1665  if (binHi < rangeLo[ivar] || binLo > rangeHi[ivar]) {
1666  // bin is outside of allowed range - effective bin volume is zero
1667  binVolumeSumSetInRange = 0.;
1668  break;
1669  }
1670 
1671  binVolumeSumSetFull *= binHi - binLo;
1672  binVolumeSumSetInRange *= std::min(rangeHi[ivar], binHi) - std::max(rangeLo[ivar], binLo);
1673  }
1674  const Double_t corrPartial = binVolumeSumSetInRange / binVolumeSumSetFull;
1675  if (0. == corrPartial) continue;
1676  const Double_t corr = correctForBinSize ? (inverseBinCor ? binVolumeSumSetFull / _binv[ibin] : binVolumeSumSetFull ) : 1.0;
1677  total += getBinScale(ibin)*(get_wgt(ibin) * corr * corrPartial);
1678  }
1679 
1680  _vars = varSave;
1681 
1682  return total;
1683 }
1684 
1685 
1686 
1687 ////////////////////////////////////////////////////////////////////////////////
1688 /// Fill the transient cache with partial bin volumes with up-to-date
1689 /// values for the partial volume specified by observables 'dimSet'
1690 
1691 const std::vector<double>& RooDataHist::calculatePartialBinVolume(const RooArgSet& dimSet) const
1692 {
1693  // The code bitset has all bits set to one whose position corresponds to arguments in dimSet.
1694  // It is used as the key for the bin volume caching hash map.
1695  int code{0};
1696  {
1697  int i{0} ;
1698  for (auto const& v : _vars) {
1699  code += ((dimSet.find(*v) ? 1 : 0) << i) ;
1700  ++i;
1701  }
1702  }
1703 
1704  auto& pbinv = _pbinvCache[code];
1705  if(!pbinv.empty()) {
1706  return pbinv;
1707  }
1708  pbinv.resize(_arrSize);
1709 
1710  // Calculate plot bins of components from master index
1711  std::vector<bool> selDim(_vars.getSize());
1712  for (std::size_t i = 0; i < selDim.size(); ++i) {
1713  selDim[i] = (code >> i) & 1 ;
1714  }
1715 
1716  // Recalculate partial bin volume cache
1717  for (Int_t ibin=0; ibin < _arrSize ;ibin++) {
1718  Int_t idx(0), tmp(ibin) ;
1719  Double_t theBinVolume(1) ;
1720  for (unsigned int j=0; j < _lvvars.size(); ++j) {
1721  const RooAbsLValue* arg = _lvvars[j];
1722  assert(arg);
1723 
1724  idx = tmp / _idxMult[j];
1725  tmp -= idx*_idxMult[j];
1726  if (selDim[j]) {
1727  theBinVolume *= arg->getBinWidth(idx) ;
1728  }
1729  }
1730  pbinv[ibin] = theBinVolume ;
1731  }
1732 
1733  return pbinv;
1734 }
1735 
1736 
1737 
1738 ////////////////////////////////////////////////////////////////////////////////
1739 /// Return the number of bins
1740 
1742 {
1743  return RooAbsData::numEntries() ;
1744 }
1745 
1746 
1747 
1748 ////////////////////////////////////////////////////////////////////////////////
1749 /// Sum the weights of all bins.
1751 
1752  if (_maskedWeights.empty()) {
1754  } else {
1756  }
1757 }
1758 
1759 
1760 
1761 ////////////////////////////////////////////////////////////////////////////////
1762 /// Return the sum of weights in all entries matching cutSpec (if specified)
1763 /// and in named range cutRange (if specified)
1764 /// Return the
1765 
1766 Double_t RooDataHist::sumEntries(const char* cutSpec, const char* cutRange) const
1767 {
1768  checkInit() ;
1769 
1770  if (cutSpec==0 && cutRange==0) {
1771  return sumEntries();
1772  } else {
1773 
1774  // Setup RooFormulaVar for cutSpec if it is present
1775  RooFormula* select = 0 ;
1776  if (cutSpec) {
1777  select = new RooFormula("select",cutSpec,*get()) ;
1778  }
1779 
1780  // Otherwise sum the weights in the event
1781  ROOT::Math::KahanSum<> kahanSum;
1782  for (Int_t i=0; i < _arrSize; i++) {
1783  get(i) ;
1784  if ((!_maskedWeights.empty() && _maskedWeights[i] == 0.)
1785  || (select && select->eval() == 0.)
1786  || (cutRange && !_vars.allInRange(cutRange)))
1787  continue;
1788 
1789  kahanSum += weight(i);
1790  }
1791 
1792  if (select) delete select ;
1793 
1794  return kahanSum;
1795  }
1796 }
1797 
1798 
1799 
1800 ////////////////////////////////////////////////////////////////////////////////
1801 /// Reset all bin weights to zero
1802 
1804 {
1805  // WVE DO NOT CALL RooTreeData::reset() for binned
1806  // datasets as this will delete the bin definitions
1807 
1808  std::fill(_wgt, _wgt + _arrSize, 0.);
1809  delete[] _errLo; _errLo = nullptr;
1810  delete[] _errHi; _errHi = nullptr;
1811  delete[] _sumw2; _sumw2 = nullptr;
1812 
1814 
1815  _cache_sum_valid = false;
1816 }
1817 
1818 
1819 
1820 ////////////////////////////////////////////////////////////////////////////////
1821 /// Load bin `binNumber`, and return an argset with the coordinates of the bin centre.
1822 /// \note The argset is owned by this data hist, and this function has a side effect, because
1823 /// it alters the currently active bin.
1824 const RooArgSet* RooDataHist::get(Int_t binNumber) const
1825 {
1826  checkInit() ;
1827  _curIndex = binNumber;
1828 
1829  return RooAbsData::get(_curIndex);
1830 }
1831 
1832 
1833 
1834 ////////////////////////////////////////////////////////////////////////////////
1835 /// Return a RooArgSet with whose coordinates denote the bin centre of the bin
1836 /// enclosing the point in `coord`.
1837 /// \note The argset is owned by this data hist, and this function has a side effect, because
1838 /// it alters the currently active bin.
1839 const RooArgSet* RooDataHist::get(const RooArgSet& coord) const {
1840  return get(calcTreeIndex(coord, false));
1841 }
1842 
1843 
1844 
1845 ////////////////////////////////////////////////////////////////////////////////
1846 /// Return the volume of the bin enclosing coordinates 'coord'.
1848  checkInit() ;
1849  return _binv[calcTreeIndex(coord, false)] ;
1850 }
1851 
1852 
1853 ////////////////////////////////////////////////////////////////////////////////
1854 /// Set all the event weight of all bins to the specified value
1855 
1857 {
1858  for (Int_t i=0 ; i<_arrSize ; i++) {
1859  _wgt[i] = value ;
1860  }
1861 
1863 }
1864 
1865 
1866 
1867 ////////////////////////////////////////////////////////////////////////////////
1868 /// Create an iterator over all bins in a slice defined by the subset of observables
1869 /// listed in sliceArg. The position of the slice is given by otherArgs
1870 
1872 {
1873  // Update to current position
1874  _vars = otherArgs ;
1875  _curIndex = calcTreeIndex(_vars, true);
1876 
1877  RooAbsArg* intArg = _vars.find(sliceArg) ;
1878  if (!intArg) {
1879  coutE(InputArguments) << "RooDataHist::sliceIterator() variable " << sliceArg.GetName() << " is not part of this RooDataHist" << endl ;
1880  return 0 ;
1881  }
1882  return new RooDataHistSliceIter(*this,*intArg) ;
1883 }
1884 
1885 
1886 ////////////////////////////////////////////////////////////////////////////////
1887 /// Change the name of the RooDataHist
1888 
1889 void RooDataHist::SetName(const char *name)
1890 {
1891  if (_dir) _dir->GetList()->Remove(this);
1893  if (_dir) _dir->GetList()->Add(this);
1894 }
1895 
1896 
1897 ////////////////////////////////////////////////////////////////////////////////
1898 /// Change the title of this RooDataHist
1899 
1900 void RooDataHist::SetNameTitle(const char *name, const char* title)
1901 {
1902  if (_dir) _dir->GetList()->Remove(this);
1903  TNamed::SetNameTitle(name,title) ;
1904  if (_dir) _dir->GetList()->Add(this);
1905 }
1906 
1907 
1908 ////////////////////////////////////////////////////////////////////////////////
1909 /// Print value of the dataset, i.e. the sum of weights contained in the dataset
1910 
1911 void RooDataHist::printValue(ostream& os) const
1912 {
1913  os << numEntries() << " bins (" << sumEntries() << " weights)" ;
1914 }
1915 
1916 
1917 
1918 
1919 ////////////////////////////////////////////////////////////////////////////////
1920 /// Print argument of dataset, i.e. the observable names
1921 
1922 void RooDataHist::printArgs(ostream& os) const
1923 {
1924  os << "[" ;
1925  Bool_t first(kTRUE) ;
1926  for (const auto arg : _vars) {
1927  if (first) {
1928  first=kFALSE ;
1929  } else {
1930  os << "," ;
1931  }
1932  os << arg->GetName() ;
1933  }
1934  os << "]" ;
1935 }
1936 
1937 
1938 
1939 ////////////////////////////////////////////////////////////////////////////////
1940 /// Compute which bins of the dataset are part of the currently set fit range.
1942 {
1943  checkInit() ;
1944 
1945  _maskedWeights.assign(_wgt, _wgt + _arrSize);
1946 
1947  for (Int_t i=0; i < _arrSize; ++i) {
1948  get(i) ;
1949 
1950  for (const auto arg : _vars) {
1951  if (!arg->inRange(nullptr)) {
1952  _maskedWeights[i] = 0.;
1953  break;
1954  }
1955  }
1956  }
1957 
1958 }
1959 
1960 
1961 
1962 ////////////////////////////////////////////////////////////////////////////////
1963 /// Returns true if dataset contains entries with a non-integer weight.
1964 
1966 {
1967  for (Int_t i=0; i < _arrSize; ++i) {
1968  const double wgt = _wgt[i];
1969  double intpart;
1970  if (fabs(std::modf(wgt, &intpart)) > 1.E-10)
1971  return true;
1972  }
1973 
1974  return false;
1975 }
1976 
1977 
1978 ////////////////////////////////////////////////////////////////////////////////
1979 /// Print the details on the dataset contents
1980 
1981 void RooDataHist::printMultiline(ostream& os, Int_t content, Bool_t verbose, TString indent) const
1982 {
1984 
1985  os << indent << "Binned Dataset " << GetName() << " (" << GetTitle() << ")" << endl ;
1986  os << indent << " Contains " << numEntries() << " bins with a total weight of " << sumEntries() << endl;
1987 
1988  if (!verbose) {
1989  os << indent << " Observables " << _vars << endl ;
1990  } else {
1991  os << indent << " Observables: " ;
1993  }
1994 
1995  if(verbose) {
1996  if (_cachedVars.getSize()>0) {
1997  os << indent << " Caches " << _cachedVars << endl ;
1998  }
1999  }
2000 }
2001 
2002 void RooDataHist::printDataHistogram(ostream& os, RooRealVar* obs) const
2003 {
2004  for(Int_t i=0; i<obs->getBins(); ++i){
2005  this->get(i);
2006  obs->setBin(i);
2007  os << this->weight() << " +/- " << this->weightSquared() << endl;
2008  }
2009 }
2010 
2011 
2012 ////////////////////////////////////////////////////////////////////////////////
2013 /// Stream an object of class RooDataHist.
2014 void RooDataHist::Streamer(TBuffer &R__b) {
2015  if (R__b.IsReading()) {
2016 
2017  UInt_t R__s, R__c;
2018  Version_t R__v = R__b.ReadVersion(&R__s, &R__c);
2019 
2020  if (R__v > 2) {
2021  R__b.ReadClassBuffer(RooDataHist::Class(),this,R__v,R__s,R__c);
2022  R__b.CheckByteCount(R__s, R__c, RooDataHist::IsA());
2023  initialize(0, false);
2024  } else {
2025 
2026  // Legacy dataset conversion happens here. Legacy RooDataHist inherits from RooTreeData
2027  // which in turn inherits from RooAbsData. Manually stream RooTreeData contents on
2028  // file here and convert it into a RooTreeDataStore which is installed in the
2029  // new-style RooAbsData base class
2030 
2031  // --- This is the contents of the streamer code of RooTreeData version 2 ---
2032  UInt_t R__s1, R__c1;
2033  Version_t R__v1 = R__b.ReadVersion(&R__s1, &R__c1); if (R__v1) { }
2034 
2035  RooAbsData::Streamer(R__b);
2036  TTree* X_tree(0) ; R__b >> X_tree;
2037  RooArgSet X_truth ; X_truth.Streamer(R__b);
2038  TString X_blindString ; X_blindString.Streamer(R__b);
2039  R__b.CheckByteCount(R__s1, R__c1, TClass::GetClass("RooTreeData"));
2040  // --- End of RooTreeData-v1 streamer
2041 
2042  // Construct RooTreeDataStore from X_tree and complete initialization of new-style RooAbsData
2043  _dstore = new RooTreeDataStore(X_tree,_vars) ;
2044  _dstore->SetName(GetName()) ;
2045  _dstore->SetTitle(GetTitle()) ;
2046  _dstore->checkInit() ;
2047 
2048  RooDirItem::Streamer(R__b);
2049  R__b >> _arrSize;
2050  delete [] _wgt;
2051  _wgt = new Double_t[_arrSize];
2052  R__b.ReadFastArray(_wgt,_arrSize);
2053  delete [] _errLo;
2054  _errLo = new Double_t[_arrSize];
2056  delete [] _errHi;
2057  _errHi = new Double_t[_arrSize];
2059  delete [] _sumw2;
2060  _sumw2 = new Double_t[_arrSize];
2062  delete [] _binv;
2063  _binv = new Double_t[_arrSize];
2064  _realVars.Streamer(R__b);
2065  double tmp;
2066  R__b >> tmp; //_curWeight;
2067  R__b >> tmp; //_curWgtErrLo;
2068  R__b >> tmp; //_curWgtErrHi;
2069  R__b >> tmp; //_curSumW2;
2070  R__b >> tmp; //_curVolume;
2071  R__b >> _curIndex;
2072  R__b.CheckByteCount(R__s, R__c, RooDataHist::IsA());
2073  }
2074 
2075  } else {
2076 
2077  R__b.WriteClassBuffer(RooDataHist::Class(),this);
2078  }
2079 }
2080 
2081 
2082 ////////////////////////////////////////////////////////////////////////////////
2083 /// Return event weights of all events in range [first, first+len).
2084 /// If no contiguous structure of weights is stored, an empty batch is be returned.
2085 /// In this case, the single-value weight() needs to be used to retrieve it.
2086 RooSpan<const double> RooDataHist::getWeightBatch(std::size_t first, std::size_t len) const {
2087  return _maskedWeights.empty() ?
2088  RooSpan<const double>{_wgt + first, len} :
2090 }
2091 
2092 
2093 ////////////////////////////////////////////////////////////////////////////////
2094 /// Write information to retrieve data columns into `evalData.spans`.
2095 /// All spans belonging to variables of this dataset are overwritten. Spans to other
2096 /// variables remain intact.
2097 /// \param[out] evalData Store references to all data batches in this struct's `spans`.
2098 /// The key to retrieve an item is the pointer of the variable that owns the data.
2099 /// \param first Index of first event that ends up in the batch.
2100 /// \param len Number of events in each batch.
2101 void RooDataHist::getBatches(RooBatchCompute::RunContext& evalData, std::size_t begin, std::size_t len) const {
2102  for (auto&& batch : store()->getBatches(begin, len).spans) {
2103  evalData.spans[batch.first] = std::move(batch.second);
2104  }
2105 }
2106 
2107 ////////////////////////////////////////////////////////////////////////////////
2108 /// Hand over pointers to our weight arrays to the data store implementation.
2111 }
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:261
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
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:748
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
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
Version_t
short Version_t
Definition: RtypesCore.h:65
RooArgSet::addClone
RooAbsArg * addClone(const RooAbsArg &var, Bool_t silent=kFALSE) override
Add clone of specified element to an owning set.
Definition: RooArgSet.cxx:288
RooMsgService.h
RooDataHist::setAllWeights
void setAllWeights(Double_t value)
Set all the event weight of all bins to the specified value.
Definition: RooDataHist.cxx:1856
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:880
RooDataHist::kInverseBinCorr
@ kInverseBinCorr
Definition: RooDataHist.h:263
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:98
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:68
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:802
RooDataHist::_errLo
double * _errLo
Definition: RooDataHist.h:248
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:263
RooAbsBinning::binNumber
virtual Int_t binNumber(Double_t x) const =0
RooLinkedList::GetSize
Int_t GetSize() const
Definition: RooLinkedList.h:60
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:330
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:1750
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:8903
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:584
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:813
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:2086
RooAbsData::store
RooAbsDataStore * store()
Definition: RooAbsData.h:68
coutI
#define coutI(a)
Definition: RooMsgService.h:30
RooDataHist::_binv
double * _binv
Definition: RooDataHist.h:251
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:70
TBuffer
Buffer base class used for serializing objects.
Definition: TBuffer.h:43
RooAbsData::checkInit
void checkInit() const
Definition: RooAbsData.cxx:2331
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:1191
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:1900
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:2101
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:1479
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:265
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:2109
RooDataHist::weightSquared
Double_t weightSquared() const override
Return squared weight of last bin that was requested with get().
Definition: RooDataHist.h:179
RooDataHist::kCorrectForBinSize
@ kCorrectForBinSize
Definition: RooDataHist.h:263
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:238
RooAbsBinning::binCenter
virtual Double_t binCenter(Int_t bin) const =0
RooDataHist::checkBinBounds
void checkBinBounds() const
Definition: RooDataHist.cxx:791
RooArgSet::snapshot
RooArgSet * snapshot(bool deepCopy=true) const
Use RooAbsCollection::snapshot(), but return as RooArgSet.
Definition: RooArgSet.h:118
RooAbsLValue.h
RooDataHist::_idxMult
std::vector< Int_t > _idxMult
Definition: RooDataHist.h:245
total
static unsigned int total
Definition: TGWin32ProxyDefs.h:40
RooDataHist::_arrSize
Int_t _arrSize
Definition: RooDataHist.h:244
RooDataHist::_realVars
RooArgSet _realVars
Definition: RooDataHist.h:253
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:183
RooDataHist::_sumw2
double * _sumw2
Definition: RooDataHist.h:250
RooAbsData::Expected
@ Expected
Definition: RooAbsData.h:99
RooDataHist::interpolateDim
Double_t interpolateDim(RooRealVar &dim, const RooAbsBinning *binning, Double_t xval, Int_t intOrder, Bool_t correctForBinSize, Bool_t cdfBoundaries)
Perform boundary safe 'intOrder'-th interpolation of weights in dimension 'dim' at current value 'xva...
Definition: RooDataHist.cxx:1248
RooFit::DataHandling
@ DataHandling
Definition: RooGlobalFunc.h:69
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:211
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:1033
RooDataHist::get
const RooArgSet * get() const override
Get bin centre of current bin.
Definition: RooDataHist.h:74
RooAbsRealLValue::numBins
virtual Int_t numBins(const char *rangeName=0) const
Definition: RooAbsRealLValue.h:54
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:259
TRACE_DESTROY
#define TRACE_DESTROY
Definition: RooTrace.h:24
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:898
TH1::GetYaxis
TAxis * GetYaxis()
Definition: TH1.h:321
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)
Cache for sum of entries ;.
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
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
RooDataHist::reset
void reset() override
Reset all bin weights to zero.
Definition: RooDataHist.cxx:1803
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:1922
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
RooArgSet::add
Bool_t add(const RooAbsArg &var, Bool_t silent=kFALSE) override
Add element to non-owning set.
Definition: RooArgSet.cxx:261
RooAbsCollection::selectCommon
RooAbsCollection * selectCommon(const RooAbsCollection &refColl) const
Create a subset of the current collection, consisting only of those elements that are contained as we...
Definition: RooAbsCollection.cxx:700
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:1911
RooLinkedList
RooLinkedList is an collection class for internal use, storing a collection of RooAbsArg pointers in ...
Definition: RooLinkedList.h:35
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:975
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:31
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:590
TRACE_CREATE
#define TRACE_CREATE
Definition: RooTrace.h:23
RooDataHist::_lvbins
std::vector< std::unique_ptr< const RooAbsBinning > > _lvbins
List of observables casted as RooAbsLValue.
Definition: RooDataHist.h:260
RooAbsCollection::size
Storage_t::size_type size() const
Definition: RooAbsCollection.h:165
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:2946
RooCategory.h
RooPrintable::kValue
@ kValue
Definition: RooPrintable.h:33
sqrt
double sqrt(double)
RooDataHist::cacheValidEntries
void cacheValidEntries()
Compute which bins of the dataset are part of the currently set fit range.
Definition: RooDataHist.cxx:1941
RooDataHist::_maskedWeights
std::vector< double > _maskedWeights
Definition: RooDataHist.h:254
RooDataHist::isNonPoissonWeighted
Bool_t isNonPoissonWeighted() const override
Returns true if dataset contains entries with a non-integer weight.
Definition: RooDataHist.cxx:1965
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:1517
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:174
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:206
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:256
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:1981
RooAbsData::PlotOpt::bins
RooAbsBinning * bins
Definition: RooAbsData.h:158
unsigned int
ymin
float ymin
Definition: THbookFile.cxx:95
RooDataHist::SetName
void SetName(const char *name) override
Change the name of the RooDataHist.
Definition: RooDataHist.cxx:1889
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:62
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:1691
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
RooAbsRealLValue::getBin
virtual Int_t getBin(const char *rangeName=0) const
Definition: RooAbsRealLValue.h:53
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:1871
RooDataHist::numEntries
Int_t numEntries() const override
Return the number of bins.
Definition: RooDataHist.cxx:1741
RooAbsDataStore::checkInit
virtual void checkInit() const
Definition: RooAbsDataStore.h:116
RooDataHist::_wgt
double * _wgt
Definition: RooDataHist.h:247
RooPrintable::kExtras
@ kExtras
Definition: RooPrintable.h:33
Double_t
double Double_t
Definition: RtypesCore.h:59
RooAbsRealLValue::setBinFast
virtual void setBinFast(Int_t ibin, const RooAbsBinning &binning)
Set value to center of bin 'ibin' of binning 'rangeName' (or of default binning if no range is specif...
Definition: RooAbsRealLValue.cxx:481
TList::Remove
virtual TObject * Remove(TObject *obj)
Remove object from the list.
Definition: TList.cxx:822
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:953
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:231
RooDataHist::_errHi
double * _errHi
Definition: RooDataHist.h:249
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:239
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:1369
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:263
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
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:264
RooMath.h
RooDataHist::printDataHistogram
void printDataHistogram(std::ostream &os, RooRealVar *obs) const
Definition: RooDataHist.cxx:2002
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:1375
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:240
TAxis::GetNbins
Int_t GetNbins() const
Definition: TAxis.h:121
RooAbsCollection::getSize
Int_t getSize() const
Definition: RooAbsCollection.h:182
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:29
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:258
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