Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
THnBase.cxx
Go to the documentation of this file.
1// @(#)root/hist:$Id$
2// Author: Axel Naumann (2011-12-20)
3
4/*************************************************************************
5 * Copyright (C) 1995-2012, Rene Brun and Fons Rademakers. *
6 * All rights reserved. *
7 * *
8 * For the licensing terms see $ROOTSYS/LICENSE. *
9 * For the list of contributors see $ROOTSYS/README/CREDITS. *
10 *************************************************************************/
11
12#include "THnBase.h"
13
14#include "TAxis.h"
15#include "TBrowser.h"
16#include "TError.h"
17#include "TClass.h"
18#include "TF1.h"
19#include "TH1D.h"
20#include "TH2D.h"
21#include "TH3D.h"
22#include "THn.h"
23#include "THnSparse.h"
24#include "TMath.h"
25#include "TRandom.h"
26#include "TVirtualPad.h"
27#include "THashList.h"
28#include "TObjString.h"
29
30#include "HFitInterface.h"
31#include "Fit/DataRange.h"
32#include "Fit/SparseData.h"
35
36
37/** \class THnBase
38 \ingroup Hist
39Multidimensional histogram base.
40Defines common functionality and interfaces for THn, THnSparse.
41*/
42
44
45////////////////////////////////////////////////////////////////////////////////
46/// Construct a THnBase with "dim" dimensions,
47/// "nbins" holds the number of bins for each dimension;
48/// "xmin" and "xmax" the minimal and maximal value for each dimension.
49/// The arrays "xmin" and "xmax" can be NULL; in that case SetBinEdges()
50/// must be called for each dimension.
51
52THnBase::THnBase(const char* name, const char* title, Int_t dim,
53 const Int_t* nbins, const Double_t* xmin, const Double_t* xmax):
54TNamed(name, title), fNdimensions(dim), fAxes(dim), fBrowsables(dim),
55fEntries(0), fTsumw(0), fTsumw2(-1.), fTsumwx(dim), fTsumwx2(dim),
56fIntegral(0), fIntegralStatus(kNoInt)
57{
58 for (Int_t i = 0; i < fNdimensions; ++i) {
59 TAxis* axis = new TAxis(nbins[i], xmin ? xmin[i] : 0., xmax ? xmax[i] : 1.);
60 axis->SetName(TString::Format("axis%d", i));
61 fAxes.AddAtAndExpand(axis, i);
62 }
63 SetTitle(title);
65}
66
67THnBase::THnBase(const char *name, const char *title, Int_t dim, const Int_t *nbins,
68 const std::vector<std::vector<double>> &xbins)
69 : TNamed(name, title), fNdimensions(dim), fAxes(dim), fBrowsables(dim), fEntries(0), fTsumw(0), fTsumw2(-1.),
70 fTsumwx(dim), fTsumwx2(dim), fIntegral(0), fIntegralStatus(kNoInt)
71{
72 if (Int_t(xbins.size()) != fNdimensions) {
73 Error("THnBase", "Mismatched number of dimensions %d with number of bin edge vectors %zu", fNdimensions,
74 xbins.size());
75 }
76 for (Int_t i = 0; i < fNdimensions; ++i) {
77 if (Int_t(xbins[i].size()) != (nbins[i] + 1)) {
78 Error("THnBase", "Mismatched number of bins %d with number of bin edges %zu", nbins[i], xbins[i].size());
79 }
80 TAxis *axis = new TAxis(nbins[i], xbins[i].data());
81 axis->SetName(TString::Format("axis%d", i));
82 fAxes.AddAtAndExpand(axis, i);
83 }
84 SetTitle(title);
86}
87
89 : TNamed(other), fNdimensions(other.fNdimensions), fAxes(fNdimensions), fBrowsables(fNdimensions),
90 fEntries(other.fEntries), fTsumw(other.fTsumw), fTsumw2(other.fTsumw2), fTsumwx(other.fTsumwx),
91 fTsumwx2(other.fTsumwx2), fIntegral(other.fIntegral), fIntegralStatus(other.fIntegralStatus)
92{
93
94 for (Int_t i = 0; i < fNdimensions; ++i) {
95 TAxis *axis = new TAxis(*static_cast<TAxis *>(other.fAxes[i]));
96 fAxes.AddAtAndExpand(axis, i);
97 }
99}
100
102{
103
104 if (this == &other)
105 return *this;
106
107 TNamed::operator=(other);
111 fEntries = other.fEntries;
112 fTsumw = other.fTsumw;
113 fTsumw2 = other.fTsumw2;
114 fTsumwx = other.fTsumwx;
115 fTsumwx2 = other.fTsumwx2;
116 fIntegral = other.fIntegral;
118
119 for (Int_t i = 0; i < fNdimensions; ++i) {
120 TAxis *axis = new TAxis(*static_cast<TAxis *>(other.fAxes[i]));
121 fAxes.AddAtAndExpand(axis, i);
122 }
123 fAxes.SetOwner();
124
125 return *this;
126}
127
129 : TNamed(other), fNdimensions(other.fNdimensions), fAxes(other.fAxes), fBrowsables(fNdimensions),
130 fEntries(other.fEntries), fTsumw(other.fTsumw), fTsumw2(other.fTsumw2), fTsumwx(other.fTsumwx),
131 fTsumwx2(other.fTsumwx2), fIntegral(std::move(other.fIntegral)), fIntegralStatus(other.fIntegralStatus)
132{
133
134 other.fAxes.SetOwner(false);
135 other.fAxes.Clear();
136 fAxes.SetOwner();
137}
138
140{
141
142 if (this == &other)
143 return *this;
144
145 TNamed::operator=(other);
146 fNdimensions = other.fNdimensions;
147 fAxes = other.fAxes;
149 fEntries = other.fEntries;
150 fTsumw = other.fTsumw;
151 fTsumw2 = other.fTsumw2;
152 fTsumwx = other.fTsumwx;
153 fTsumwx2 = other.fTsumwx2;
154 fIntegral = std::move(other.fIntegral);
155 fIntegralStatus = other.fIntegralStatus;
156
157 other.fAxes.SetOwner(false);
158 other.fAxes.Clear();
159 fAxes.SetOwner();
160
161 return *this;
162}
163
164////////////////////////////////////////////////////////////////////////////////
165/// Destruct a THnBase
166
168 if (fIntegralStatus != kNoInt)
169 fIntegral.clear();
170}
171
172
173////////////////////////////////////////////////////////////////////////////////
174/// Create a new THnBase object that is of the same type as *this,
175/// but with dimensions and bins given by axes.
176/// If keepTargetAxis is true, the axes will keep their original xmin / xmax,
177/// else they will be restricted to the range selected (first / last).
178
179THnBase* THnBase::CloneEmpty(const char* name, const char* title,
180 const TObjArray* axes, Bool_t keepTargetAxis) const
181{
182 THnBase* ret = (THnBase*)IsA()->New();
183 Int_t chunkSize = 1024 * 16;
185 chunkSize = ((const THnSparse*)this)->GetChunkSize();
186 }
187 ret->Init(name, title, axes, keepTargetAxis, chunkSize);
188 return ret;
189}
190
191
192////////////////////////////////////////////////////////////////////////////////
193/// Initialize axes and name.
194
195void THnBase::Init(const char* name, const char* title,
196 const TObjArray* axes, Bool_t keepTargetAxis,
197 Int_t chunkSize /*= 1024 * 16*/)
198{
199 SetNameTitle(name, title);
200 if (!axes) {
201 ::Error("THnBase::Init", "Input parameter `axes` is null, no axes were provided at initialization");
202 return;
203 }
204 TIter iAxis(axes);
205 const TAxis* axis = nullptr;
206 Int_t pos = 0;
207 Int_t *nbins = new Int_t[axes->GetEntriesFast()];
208 while ((axis = (TAxis*)iAxis())) {
209 if (!axis) {
210 ::Error("THnBase::Init", "Input parameter `axes` has a null element in the array, cannot create new axis");
211 continue;
212 }
213 TAxis* reqaxis = new TAxis(*axis);
214 if (!keepTargetAxis && axis->TestBit(TAxis::kAxisRange)) {
215 Int_t binFirst = axis->GetFirst();
216 // The lowest egde of the underflow is meaningless.
217 if (binFirst == 0)
218 binFirst = 1;
219 Int_t binLast = axis->GetLast();
220 // The overflow edge is implicit.
221 if (binLast > axis->GetNbins())
222 binLast = axis->GetNbins();
223 Int_t nBins = binLast - binFirst + 1;
224 if (axis->GetXbins()->GetSize()) {
225 // non-uniform bins:
226 reqaxis->Set(nBins, axis->GetXbins()->GetArray() + binFirst - 1);
227 } else {
228 // uniform bins:
229 reqaxis->Set(nBins, axis->GetBinLowEdge(binFirst), axis->GetBinUpEdge(binLast));
230 }
231 reqaxis->ResetBit(TAxis::kAxisRange);
232 }
233
234 nbins[pos] = reqaxis->GetNbins();
235 fAxes.AddAtAndExpand(reqaxis, pos++);
236 }
237 fAxes.SetOwner();
238
240 InitStorage(nbins, chunkSize);
243 delete [] nbins;
244}
245
246
247////////////////////////////////////////////////////////////////////////////////
248/// Create an empty histogram with name and title with a given
249/// set of axes. Create a TH1D/TH2D/TH3D, depending on the number
250/// of elements in axes.
251
252TH1* THnBase::CreateHist(const char* name, const char* title,
253 const TObjArray* axes,
254 Bool_t keepTargetAxis ) const {
255 if (!axes) {
256 ::Error("THnBase::CreateHist", "Input parameter `axes` is null, no axes were provided at creation");
257 return nullptr;
258 }
259 const int ndim = axes->GetSize();
260
261 TH1* hist = nullptr;
262 // create hist with dummy axes, we fix them later.
263 if (ndim == 1)
264 hist = new TH1D(name, title, 1, 0., 1.);
265 else if (ndim == 2)
266 hist = new TH2D(name, title, 1, 0., 1., 1, 0., 1.);
267 else if (ndim == 3)
268 hist = new TH3D(name, title, 1, 0., 1., 1, 0., 1., 1, 0., 1.);
269 else {
270 Error("CreateHist", "Cannot create histogram %s with %d dimensions!", name, ndim);
271 return nullptr;
272 }
273
274 TAxis* hax[3] = {hist->GetXaxis(), hist->GetYaxis(), hist->GetZaxis()};
275 for (Int_t d = 0; d < ndim; ++d) {
276 TAxis* reqaxis = (TAxis*)(*axes)[d];
277 if (!reqaxis) {
278 ::Error("THnBase::CreateHist", "Input parameter `axes` has a null element in the position %d of the array, cannot create new axis", d);
279 continue;
280 }
281 hax[d]->SetTitle(reqaxis->GetTitle());
282 if (!keepTargetAxis && reqaxis->TestBit(TAxis::kAxisRange)) {
283 // axis cannot extend to underflow/overflows (fix ROOT-8781)
284 Int_t binFirst = std::max(reqaxis->GetFirst(),1);
285 Int_t binLast = std::min(reqaxis->GetLast(), reqaxis->GetNbins() );
286 Int_t nBins = binLast - binFirst + 1;
287 if (reqaxis->GetXbins()->GetSize()) {
288 // non-uniform bins:
289 hax[d]->Set(nBins, reqaxis->GetXbins()->GetArray() + binFirst - 1);
290 } else {
291 // uniform bins:
292 hax[d]->Set(nBins, reqaxis->GetBinLowEdge(binFirst), reqaxis->GetBinUpEdge(binLast));
293 }
294 } else {
295 if (reqaxis->GetXbins()->GetSize()) {
296 // non-uniform bins:
297 hax[d]->Set(reqaxis->GetNbins(), reqaxis->GetXbins()->GetArray());
298 } else {
299 // uniform bins:
300 hax[d]->Set(reqaxis->GetNbins(), reqaxis->GetXmin(), reqaxis->GetXmax());
301 }
302 // Copy the axis labels if needed.
303 THashList* labels = reqaxis->GetLabels();
304 if (labels) {
305 TIter iL(labels);
306 Int_t i = 1;
307 while (auto lb = static_cast<TObjString *>(iL())) {
308 hax[d]->SetBinLabel(i,lb->String().Data());
309 i++;
310 }
311 }
312 }
313 }
314
315 hist->Rebuild();
316
317 return hist;
318}
319
320////////////////////////////////////////////////////////////////////////////////
321/// Create a THn / THnSparse object from a histogram deriving from TH1.
322
323THnBase* THnBase::CreateHnAny(const char* name, const char* title,
324 const TH1* h, Bool_t sparse, Int_t chunkSize)
325{
326 if (!h) {
327 ::Error("THnBase::CreateHnAny", "Input parameter `h` is null, no histogram was provided upon creation");
328 return nullptr;
329 }
330 // Get the dimension of the TH1
331 int ndim = h->GetDimension();
332
333 if (ndim > 3) {
334 ::Error("THnBase::CreateHnAny", "Can only handle 1..3 dimensional histograms!");
335 return nullptr;
336 }
337
338 // Axis properties
339 int nbins[3] = {0,0,0};
340 double minRange[3] = {0.,0.,0.};
341 double maxRange[3] = {0.,0.,0.};
342 const TAxis* axis[3] = { h->GetXaxis(), h->GetYaxis(), h->GetZaxis() };
343 for (int i = 0; i < ndim; ++i) {
344 nbins[i] = axis[i]->GetNbins();
345 minRange[i] = axis[i]->GetXmin();
346 maxRange[i] = axis[i]->GetXmax();
347 }
348
349 // Create the corresponding THnSparse, depending on the storage
350 // type of the TH1. The class name will be "TH??\0" where the first
351 // ? is 1,2 or 3 and the second ? indicates the storage as C, S,
352 // I, L, F or D.
353 THnBase* s = nullptr;
354 const char* cname( h->ClassName() );
355 if (cname[0] == 'T' && cname[1] == 'H'
356 && cname[2] >= '1' && cname[2] <= '3' && cname[4] == 0) {
357
358#define R__THNBCASE(TAG) \
359if (sparse) { \
360s = new _NAME2_(THnSparse,TAG)(name, title, ndim, nbins, \
361minRange, maxRange, chunkSize); \
362} else { \
363s = new _NAME2_(THn,TAG)(name, title, ndim, nbins, \
364minRange, maxRange); \
365} \
366break;
367
368 switch (cname[3]) {
369 case 'F': R__THNBCASE(F);
370 case 'D': R__THNBCASE(D);
371 case 'I': R__THNBCASE(I);
372 case 'L': R__THNBCASE(L);
373 case 'S': R__THNBCASE(S);
374 case 'C': R__THNBCASE(C);
375 }
376#undef R__THNBCASE
377 }
378 if (!s) {
379 ::Warning("THnSparse::CreateHnAny", "Unknown Type of Histogram");
380 return nullptr;
381 }
382
383 for (int i = 0; i < ndim; ++i) {
384 s->GetAxis(i)->SetTitle(axis[i]->GetTitle());
385 }
386
387 // Get the array to know the number of entries of the TH1
388 const TArray *array = dynamic_cast<const TArray*>(h);
389 if (!array) {
390 ::Warning("THnSparse::CreateHnAny", "Unknown Type of Histogram");
391 return nullptr;
392 }
393
394 s->Add(h);
395 return s;
396}
397
398
399////////////////////////////////////////////////////////////////////////////////
400/// Create a THnSparse (if "sparse") or THn from "hn", possibly
401/// converting THn <-> THnSparse.
402
403THnBase* THnBase::CreateHnAny(const char* name, const char* title,
404 const THnBase* hn, Bool_t sparse,
405 Int_t chunkSize /*= 1024 * 16*/)
406{
407 if (!hn) {
408 ::Error("THnBase::CreateHnAny", "Input parameter `hn` is null, no histogram was provided upon creation");
409 return nullptr;
410 }
411 TClass* type = nullptr;
412 if (hn->InheritsFrom(THnSparse::Class())) {
413 if (sparse) type = hn->IsA();
414 else {
415 char bintype;
416 if (hn->InheritsFrom(THnSparseD::Class())) bintype = 'D';
417 else if (hn->InheritsFrom(THnSparseF::Class())) bintype = 'F';
418 else if (hn->InheritsFrom(THnSparseL::Class())) bintype = 'L';
419 else if (hn->InheritsFrom(THnSparseI::Class())) bintype = 'I';
420 else if (hn->InheritsFrom(THnSparseS::Class())) bintype = 'S';
421 else if (hn->InheritsFrom(THnSparseC::Class())) bintype = 'C';
422 else {
423 hn->Error("CreateHnAny", "Type %s not implemented; please inform the ROOT team!",
424 hn->IsA()->GetName());
425 return nullptr;
426 }
427 type = TClass::GetClass(TString::Format("THn%c", bintype));
428 }
429 } else if (hn->InheritsFrom(THn::Class())) {
430 if (!sparse) type = hn->IsA();
431 else {
432 char bintype = 0;
433 if (hn->InheritsFrom(THnD::Class())) bintype = 'D';
434 else if (hn->InheritsFrom(THnF::Class())) bintype = 'F';
435 else if (hn->InheritsFrom(THnC::Class())) bintype = 'C';
436 else if (hn->InheritsFrom(THnS::Class())) bintype = 'S';
437 else if (hn->InheritsFrom(THnI::Class())) bintype = 'I';
438 else if (hn->InheritsFrom(THnL::Class()) || hn->InheritsFrom(THnL64::Class())) bintype = 'L';
439 if (bintype) {
440 type = TClass::GetClass(TString::Format("THnSparse%c", bintype));
441 }
442 }
443 } else {
444 hn->Error("CreateHnAny", "Unhandled type %s, not deriving from THn nor THnSparse!",
445 hn->IsA()->GetName());
446 return nullptr;
447 }
448 if (!type) {
449 hn->Error("CreateHnAny", "Unhandled type %s, please inform the ROOT team!",
450 hn->IsA()->GetName());
451 return nullptr;
452 }
453
454 THnBase* ret = (THnBase*)type->New();
455 ret->Init(name, title, hn->GetListOfAxes(),
456 kFALSE /*keepTargetAxes*/, chunkSize);
457
458 ret->Add(hn);
459 return ret;
460}
461
462
463////////////////////////////////////////////////////////////////////////////////
464/// Fill the THnBase with the bins of hist that have content
465/// or error != 0.
466
467void THnBase::Add(const TH1* hist, Double_t c /*=1.*/)
468{
469 if (!hist) {
470 ::Error("THnBase::Add", "Input parameter `hist` is null, no histogram was provided");
471 return;
472 }
473 Long64_t nbins = hist->GetNcells();
474 int x[3] = {0,0,0};
475 for (int i = 0; i < nbins; ++i) {
476 double value = hist->GetBinContent(i);
477 double error = hist->GetBinError(i);
478 if (!value && !error) continue;
479 hist->GetBinXYZ(i, x[0], x[1], x[2]);
481 SetBinError(x, error * c);
482 }
483}
484
485////////////////////////////////////////////////////////////////////////////////
486/// Fit a THnSparse with function f
487///
488/// since the data is sparse by default a likelihood fit is performed
489/// merging all the regions with empty bins for better performance efficiency
490///
491/// Since the THnSparse is not drawn no graphics options are passed
492/// Here is the list of possible options
493///
494/// = "I" Use integral of function in bin instead of value at bin center
495/// = "X" Use chi2 method (default is log-likelihood method)
496/// = "U" Use a User specified fitting algorithm (via SetFCN)
497/// = "Q" Quiet mode (minimum printing)
498/// = "V" Verbose mode (default is between Q and V)
499/// = "E" Perform better Errors estimation using Minos technique
500/// = "B" Use this option when you want to fix one or more parameters
501/// and the fitting function is like "gaus", "expo", "poln", "landau".
502/// = "M" More. Improve fit results
503/// = "R" Use the Range specified in the function range
504
506{
507
508 Foption_t fitOption;
509
510 if (!TH1::FitOptionsMake(option,fitOption)) return 0;
511
512 // The function used to fit cannot be stored in a THnSparse. It
513 // cannot be drawn either. Perhaps in the future.
514 fitOption.Nostore = true;
515 // Use likelihood fit if not specified
516 if (!fitOption.Chi2) fitOption.Like = true;
517 // create range and minimizer options with default values
519 for ( int i = 0; i < GetNdimensions(); ++i ) {
520 TAxis *axis = GetAxis(i);
521 range.AddRange(i, axis->GetXmin(), axis->GetXmax());
522 }
524
525 return ROOT::Fit::FitObject(this, f , fitOption , minOption, goption, range);
526}
527
528////////////////////////////////////////////////////////////////////////////////
529/// \brief THnBase::GetBinCenter
530/// \param idx an array of bin index in each dimension.
531/// \return vector of bin centers in each dimension; empty in case of error.
532/// \note Throws error if size is different from nDimensions.
533/// \sa GetAxis(dim)::GetBinCenter(idx) as an alternative
534std::vector<Double_t> THnBase::GetBinCenter(const std::vector<Int_t> &idx) const
535{
536 if (idx.size() != static_cast<decltype(idx.size())>(fNdimensions)) {
537 Error("THnBase::GetBinCenter",
538 "Mismatched number of dimensions %d with bin index vector size %zu, returning an empty vector.",
539 fNdimensions, idx.size());
540 return {};
541 }
542 std::vector<Double_t> centers(fNdimensions);
543 std::generate(centers.begin(), centers.end(), [i = 0, &idx, this]() mutable {
544 auto bincenter = GetAxis(i)->GetBinCenter(idx[i]);
545 i++;
546 return bincenter;
547 });
548 return centers;
549}
550
551////////////////////////////////////////////////////////////////////////////////
552/// Generate an n-dimensional random tuple based on the histogrammed
553/// distribution. If subBinRandom, the returned tuple will be additionally
554/// randomly distributed within the randomized bin, using a flat
555/// distribution.
556
557void THnBase::GetRandom(Double_t *rand, Bool_t subBinRandom /* = kTRUE */)
558{
559 // check whether the integral array is valid
562
563 // generate a random bin
564 Double_t p = gRandom->Rndm();
565 Long64_t idx = TMath::BinarySearch(GetNbins() + 1, fIntegral.data(), p);
566 const Int_t nStaticBins = 40;
567 Int_t bin[nStaticBins];
568 Int_t* pBin = bin;
569 if (GetNdimensions() > nStaticBins) {
570 pBin = new Int_t[GetNdimensions()];
571 }
572 GetBinContent(idx, pBin);
573
574 // convert bin coordinates to real values
575 for (Int_t i = 0; i < fNdimensions; i++) {
576 rand[i] = GetAxis(i)->GetBinCenter(pBin[i]);
577
578 // randomize the vector within a bin
579 if (subBinRandom)
580 rand[i] += (gRandom->Rndm() - 0.5) * GetAxis(i)->GetBinWidth(pBin[i]);
581 }
582 if (pBin != bin) {
583 delete [] pBin;
584 }
585
586 return;
587}
588
589////////////////////////////////////////////////////////////////////////////////
590/// Check whether bin coord is in range, as defined by TAxis::SetRange().
591
593{
594 Int_t min = 0;
595 Int_t max = 0;
596 for (Int_t i = 0; i < fNdimensions; ++i) {
597 TAxis *axis = GetAxis(i);
598 if (!axis->TestBit(TAxis::kAxisRange)) continue;
599 min = axis->GetFirst();
600 max = axis->GetLast();
601 if (coord[i] < min || coord[i] > max)
602 return kFALSE;
603 }
604 return kTRUE;
605}
606
607////////////////////////////////////////////////////////////////////////////////
608/// Project all bins into a ndim-dimensional THn / THnSparse (whatever
609/// *this is) or if (ndim < 4 and !wantNDim) a TH1/2/3 histogram,
610/// keeping only axes in dim (specifying ndim dimensions).
611/// If "option" contains:
612/// - "E" errors will be calculated.
613/// - "A" ranges of the target axes will be ignored.
614/// - "O" original axis range of the target axes will be
615/// kept, but only bins inside the selected range
616/// will be filled.
617
619 Bool_t wantNDim,
620 Option_t* option /*= ""*/) const
621{
623 name +="_proj";
624
625 for (Int_t d = 0; d < ndim; ++d) {
626 name += "_";
627 name += dim[d];
628 }
629
630 TString title(GetTitle());
631 Ssiz_t posInsert = title.First(';');
632 if (posInsert == kNPOS) {
633 title += " projection ";
634 for (Int_t d = 0; d < ndim; ++d)
635 title += GetAxis(dim[d])->GetTitle();
636 } else {
637 for (Int_t d = ndim - 1; d >= 0; --d) {
638 title.Insert(posInsert, GetAxis(d)->GetTitle());
639 if (d)
640 title.Insert(posInsert, ", ");
641 }
642 title.Insert(posInsert, " projection ");
643 }
644
645 TObjArray newaxes(ndim);
646 for (Int_t d = 0; d < ndim; ++d) {
647 newaxes.AddAt(GetAxis(dim[d]),d);
648 }
649
650 THnBase* hn = nullptr;
651 TH1* hist = nullptr;
652 TObject* ret = nullptr;
653
654 Bool_t* hadRange = nullptr;
655 Bool_t ignoreTargetRange = (option && (strchr(option, 'A') || strchr(option, 'a')));
656 Bool_t keepTargetAxis = ignoreTargetRange || (option && (strchr(option, 'O') || strchr(option, 'o')));
657 if (ignoreTargetRange) {
658 hadRange = new Bool_t[ndim];
659 for (Int_t d = 0; d < ndim; ++d){
660 TAxis *axis = GetAxis(dim[d]);
661 hadRange[d] = axis->TestBit(TAxis::kAxisRange);
663 }
664 }
665
666 if (wantNDim)
667 ret = hn = CloneEmpty(name, title, &newaxes, keepTargetAxis);
668 else
669 ret = hist = CreateHist(name, title, &newaxes, keepTargetAxis);
670
671 if (keepTargetAxis) {
672 // make the whole axes visible, i.e. unset the range
673 if (wantNDim) {
674 for (Int_t d = 0; d < ndim; ++d) {
675 hn->GetAxis(d)->SetRange(0, 0);
676 }
677 } else {
678 hist->GetXaxis()->SetRange(0, 0);
679 hist->GetYaxis()->SetRange(0, 0);
680 hist->GetZaxis()->SetRange(0, 0);
681 }
682 }
683
684 Bool_t haveErrors = GetCalculateErrors();
685 Bool_t wantErrors = haveErrors || (option && (strchr(option, 'E') || strchr(option, 'e')));
686
687 Int_t* bins = new Int_t[ndim];
688 Long64_t myLinBin = 0;
689
690 std::unique_ptr<ROOT::Internal::THnBaseBinIter> iter{CreateIter(true /*use axis range*/)};
691
692 while ((myLinBin = iter->Next()) >= 0) {
693 Double_t v = GetBinContent(myLinBin);
694
695 for (Int_t d = 0; d < ndim; ++d) {
696 bins[d] = iter->GetCoord(dim[d]);
697 if (!keepTargetAxis && GetAxis(dim[d])->TestBit(TAxis::kAxisRange)) {
698 Int_t binOffset = GetAxis(dim[d])->GetFirst();
699 // Don't subtract even more if underflow is alreday included:
700 if (binOffset > 0) --binOffset;
701 bins[d] -= binOffset;
702 }
703 }
704
705 Long64_t targetLinBin = -1;
706 if (!wantNDim) {
707 if (ndim == 1) targetLinBin = bins[0];
708 else if (ndim == 2) targetLinBin = hist->GetBin(bins[0], bins[1]);
709 else if (ndim == 3) targetLinBin = hist->GetBin(bins[0], bins[1], bins[2]);
710 } else {
711 targetLinBin = hn->GetBin(bins, kTRUE /*allocate*/);
712 }
713
714 if (wantErrors) {
715 Double_t err2 = 0.;
716 if (haveErrors) {
717 err2 = GetBinError2(myLinBin);
718 } else {
719 err2 = v;
720 }
721 if (wantNDim) {
722 hn->AddBinError2(targetLinBin, err2);
723 } else {
724 Double_t preverr = hist->GetBinError(targetLinBin);
725 hist->SetBinError(targetLinBin, TMath::Sqrt(preverr * preverr + err2));
726 }
727 }
728
729 // only _after_ error calculation, or sqrt(v) is taken into account!
730 if (wantNDim)
731 hn->AddBinContent(targetLinBin, v);
732 else
733 hist->AddBinContent(targetLinBin, v);
734 }
735
736 delete [] bins;
737
738 if (wantNDim) {
739 hn->SetEntries(fEntries);
740 } else {
741 if (!iter->HaveSkippedBin()) {
742 hist->SetEntries(fEntries);
743 } else {
744 // re-compute the entries
745 // in case of error calculation (i.e. when Sumw2() is set)
746 // use the effective entries for the entries
747 // since this is the only way to estimate them
748 hist->ResetStats();
749 Double_t entries = hist->GetEffectiveEntries();
750 if (!wantErrors) {
751 // to avoid numerical rounding
752 entries = TMath::Floor(entries + 0.5);
753 }
754 hist->SetEntries(entries);
755 }
756 }
757
758 if (hadRange) {
759 // reset kAxisRange bit:
760 for (Int_t d = 0; d < ndim; ++d)
761 GetAxis(dim[d])->SetBit(TAxis::kAxisRange, hadRange[d]);
762
763 delete [] hadRange;
764 }
765
766 return ret;
767}
768
769////////////////////////////////////////////////////////////////////////////////
770/// Scale contents and errors of this histogram by c:
771/// this = this * c
772/// It does not modify the histogram's number of entries.
773
775{
776
777 Double_t nEntries = GetEntries();
778 // Scale the contents & errors
779 Bool_t haveErrors = GetCalculateErrors();
780 Long64_t i = 0;
781 std::unique_ptr<ROOT::Internal::THnBaseBinIter> iter{CreateIter(false)};
782 while ((i = iter->Next()) >= 0) {
783 // Get the content of the bin from the current histogram
785 SetBinContent(i, c * v);
786 if (haveErrors) {
787 Double_t err2 = GetBinError2(i);
788 SetBinError2(i, c * c * err2);
789 }
790 }
791 SetEntries(nEntries);
792}
793
794////////////////////////////////////////////////////////////////////////////////
795/// Add() implementation for both rebinned histograms and those with identical
796/// binning. See THnBase::Add().
797
799{
800 if (fNdimensions != h->GetNdimensions()) {
801 Warning("RebinnedAdd", "Different number of dimensions, cannot carry out operation on the histograms");
802 return;
803 }
804
805 // Trigger error calculation if h has it
806 if (!GetCalculateErrors() && h->GetCalculateErrors())
807 Sumw2();
808 Bool_t haveErrors = GetCalculateErrors();
809
810 Double_t* x = nullptr;
811 if (rebinned) {
812 x = new Double_t[fNdimensions];
813 }
814 Int_t* coord = new Int_t[fNdimensions];
815
816 // Expand the exmap if needed, to reduce collisions
817 Long64_t numTargetBins = GetNbins() + h->GetNbins();
818 Reserve(numTargetBins);
819
820 Long64_t i = 0;
821 std::unique_ptr<ROOT::Internal::THnBaseBinIter> iter{h->CreateIter(false)};
822 // Add to this whatever is found inside the other histogram
823 while ((i = iter->Next(coord)) >= 0) {
824 // Get the content of the bin from the second histogram
825 Double_t v = h->GetBinContent(i);
826
827 Long64_t mybinidx = -1;
828 if (rebinned) {
829 // Get the bin center given a coord
830 for (Int_t j = 0; j < fNdimensions; ++j)
831 x[j] = h->GetAxis(j)->GetBinCenter(coord[j]);
832
833 mybinidx = GetBin(x, kTRUE /* allocate*/);
834 } else {
835 mybinidx = GetBin(coord, kTRUE /*allocate*/);
836 }
837
838 if (haveErrors) {
839 Double_t err2 = h->GetBinError2(i) * c * c;
840 AddBinError2(mybinidx, err2);
841 }
842 // only _after_ error calculation, or sqrt(v) is taken into account!
843 AddBinContent(mybinidx, c * v);
844 }
845
846 delete [] coord;
847 delete [] x;
848
849 // add also the statistics
850 fTsumw += c * h->fTsumw;
851 if (haveErrors) {
852 fTsumw2 += c * c * h->fTsumw2;
853 if (h->fTsumwx.fN == fNdimensions && h->fTsumwx2.fN == fNdimensions) {
854 for (Int_t d = 0; d < fNdimensions; ++d) {
855 fTsumwx[d] += c * h->fTsumwx[d];
856 fTsumwx2[d] += c * h->fTsumwx2[d];
857 }
858 }
859 }
860
861 Double_t nEntries = GetEntries() + c * h->GetEntries();
862 SetEntries(nEntries);
863}
864
865////////////////////////////////////////////////////////////////////////////////
866/// Add contents of h scaled by c to this histogram:
867/// this = this + c * h
868/// Note that if h has Sumw2 set, Sumw2 is automatically called for this
869/// if not already set.
870
872{
873 // Check consistency of the input
874 if (!CheckConsistency(h, "Add")) return;
875
877}
878
879////////////////////////////////////////////////////////////////////////////////
880/// Add contents of h scaled by c to this histogram:
881/// this = this + c * h
882/// Note that if h has Sumw2 set, Sumw2 is automatically called for this
883/// if not already set.
884/// In contrast to Add(), RebinnedAdd() does not require consistent binning of
885/// this and h; instead, each bin's center is used to determine the target bin.
886
888{
889 AddInternal(h, c, kTRUE);
890}
891
892
893////////////////////////////////////////////////////////////////////////////////
894/// Merge this with a list of THnBase's. All THnBase's provided
895/// in the list must have the same bin layout!
896
898{
899 if (!list) return 0;
900 if (list->IsEmpty()) return (Long64_t)GetEntries();
901
902 Long64_t sumNbins = GetNbins();
903 TIter iter(list);
904 const TObject* addMeObj = nullptr;
905 while ((addMeObj = iter())) {
906 const THnBase* addMe = dynamic_cast<const THnBase*>(addMeObj);
907 if (addMe) {
908 sumNbins += addMe->GetNbins();
909 }
910 }
911 Reserve(sumNbins);
912
913 iter.Reset();
914 while ((addMeObj = iter())) {
915 const THnBase* addMe = dynamic_cast<const THnBase*>(addMeObj);
916 if (!addMe)
917 Error("Merge", "Object named %s is not THnBase! Skipping it.",
918 addMeObj->GetName());
919 else
920 Add(addMe);
921 }
922 return (Long64_t)GetEntries();
923}
924
925
926////////////////////////////////////////////////////////////////////////////////
927/// Multiply this histogram by histogram h
928/// this = this * h
929/// Note that if h has Sumw2 set, Sumw2 is automatically called for this
930/// if not already set.
931
933{
934 // Check consistency of the input
935 if(!CheckConsistency(h, "Multiply"))return;
936
937 // Trigger error calculation if h has it
938 Bool_t wantErrors = kFALSE;
939 if (GetCalculateErrors() || h->GetCalculateErrors())
940 wantErrors = kTRUE;
941
942 if (wantErrors) Sumw2();
943
944 Double_t nEntries = GetEntries();
945 // Now multiply the contents: in this case we have the intersection of the sets of bins
946 Int_t* coord = new Int_t[fNdimensions];
947 Long64_t i = 0;
948 std::unique_ptr<ROOT::Internal::THnBaseBinIter> iter{CreateIter(false)};
949 // Add to this whatever is found inside the other histogram
950 while ((i = iter->Next(coord)) >= 0) {
951 // Get the content of the bin from the current histogram
953 // Now look at the bin with the same coordinates in h
954 Long64_t idxh = h->GetBin(coord);
955 Double_t v2 = 0.;
956 if (idxh >= 0) v2 = h->GetBinContent(idxh);
957 SetBinContent(i, v1 * v2);
958 if (wantErrors) {
959 Double_t err1 = GetBinError(i) * v2;
960 Double_t err2 = 0.;
961 if (idxh >= 0) err2 = h->GetBinError(idxh) * v1;
962 SetBinError(i, TMath::Sqrt((err2 * err2 + err1 * err1)));
963 }
964 }
965 SetEntries(nEntries);
966
967 delete [] coord;
968}
969
970////////////////////////////////////////////////////////////////////////////////
971/// Performs the operation: this = this*c*f1
972/// if errors are defined, errors are also recalculated.
973///
974/// Only bins inside the function range are recomputed.
975/// IMPORTANT NOTE: If you intend to use the errors of this histogram later
976/// you should call Sumw2 before making this operation.
977/// This is particularly important if you fit the histogram after
978/// calling Multiply()
979
981{
982 Int_t* coord = new Int_t[fNdimensions];
984
985 Bool_t wantErrors = GetCalculateErrors();
986 if (wantErrors) Sumw2();
987
988 Long64_t i = 0;
989 std::unique_ptr<ROOT::Internal::THnBaseBinIter> iter{CreateIter(false)};
990 // Add to this whatever is found inside the other histogram
991 while ((i = iter->Next(coord)) >= 0) {
993
994 // Get the bin coordinates given an index array
995 for (Int_t j = 0; j < fNdimensions; ++j)
996 x[j] = GetAxis(j)->GetBinCenter(coord[j]);
997
998 if (!f->IsInside(x))
999 continue;
1001
1002 // Evaluate function at points
1003 Double_t fvalue = f->EvalPar(x, nullptr);
1004
1005 SetBinContent(i, c * fvalue * value);
1006 if (wantErrors) {
1007 Double_t error = GetBinError(i);
1008 SetBinError(i, c * fvalue * error);
1009 }
1010 }
1011
1012 delete [] x;
1013 delete [] coord;
1014}
1015
1016////////////////////////////////////////////////////////////////////////////////
1017/// Divide this histogram by h
1018/// this = this/(h)
1019/// Note that if h has Sumw2 set, Sumw2 is automatically called for
1020/// this if not already set.
1021/// The resulting errors are calculated assuming uncorrelated content.
1022
1024{
1025 // Check consistency of the input
1026 if (!CheckConsistency(h, "Divide"))return;
1027
1028 // Trigger error calculation if h has it
1029 Bool_t wantErrors=GetCalculateErrors();
1030 if (!GetCalculateErrors() && h->GetCalculateErrors())
1031 wantErrors=kTRUE;
1032
1033 // Remember original histogram statistics
1034 Double_t nEntries = fEntries;
1035
1036 if (wantErrors) Sumw2();
1037 Bool_t didWarn = kFALSE;
1038
1039 // Now divide the contents: also in this case we have the intersection of the sets of bins
1040 Int_t* coord = new Int_t[fNdimensions];
1041 Long64_t i = 0;
1042 std::unique_ptr<ROOT::Internal::THnBaseBinIter> iter{CreateIter(false)};
1043 // Add to this whatever is found inside the other histogram
1044 while ((i = iter->Next(coord)) >= 0) {
1045 // Get the content of the bin from the first histogram
1047 // Now look at the bin with the same coordinates in h
1048 Long64_t hbin = h->GetBin(coord);
1049 Double_t v2 = h->GetBinContent(hbin);
1050 if (!v2) {
1051 v1 = 0.;
1052 v2 = 1.;
1053 if (!didWarn) {
1054 Warning("Divide(h)", "Histogram h has empty bins - division by zero! Setting bin to 0.");
1055 didWarn = kTRUE;
1056 }
1057 }
1058 SetBinContent(i, v1 / v2);
1059 if (wantErrors) {
1060 Double_t err1 = GetBinError(i) * v2;
1061 Double_t err2 = h->GetBinError(hbin) * v1;
1062 Double_t b22 = v2 * v2;
1063 Double_t err = (err1 * err1 + err2 * err2) / (b22 * b22);
1064 SetBinError2(i, err);
1065 }
1066 }
1067 delete [] coord;
1068 SetEntries(nEntries);
1069}
1070
1071////////////////////////////////////////////////////////////////////////////////
1072/// Replace contents of this histogram by multiplication of h1 by h2
1073/// this = (c1*h1)/(c2*h2)
1074/// Note that if h1 or h2 have Sumw2 set, Sumw2 is automatically called for
1075/// this if not already set.
1076/// The resulting errors are calculated assuming uncorrelated content.
1077/// However, if option ="B" is specified, Binomial errors are computed.
1078/// In this case c1 and c2 do not make real sense and they are ignored.
1079
1081{
1082
1083 TString opt = option;
1084 opt.ToLower();
1085 Bool_t binomial = kFALSE;
1086 if (opt.Contains("b")) binomial = kTRUE;
1087
1088 // Check consistency of the input
1089 if (!CheckConsistency(h1, "Divide") || !CheckConsistency(h2, "Divide"))return;
1090 if (!c2) {
1091 Error("Divide","Coefficient of dividing histogram cannot be zero");
1092 return;
1093 }
1094
1095 Reset();
1096
1097 // Trigger error calculation if h1 or h2 have it
1098 if (!GetCalculateErrors() && (h1->GetCalculateErrors()|| h2->GetCalculateErrors() != 0))
1099 Sumw2();
1100
1101 // Count filled bins
1102 Long64_t nFilledBins=0;
1103
1104 // Now divide the contents: we have the intersection of the sets of bins
1105
1106 Int_t* coord = new Int_t[fNdimensions];
1107 memset(coord, 0, sizeof(Int_t) * fNdimensions);
1108 Bool_t didWarn = kFALSE;
1109
1110 Long64_t i = 0;
1111 std::unique_ptr<ROOT::Internal::THnBaseBinIter> iter{h1->CreateIter(false)};
1112 // Add to this whatever is found inside the other histogram
1113 while ((i = iter->Next(coord)) >= 0) {
1114 // Get the content of the bin from the first histogram
1116 // Now look at the bin with the same coordinates in h2
1117 Long64_t h2bin = h2->GetBin(coord);
1118 Double_t v2 = h2->GetBinContent(h2bin);
1119 if (!v2) {
1120 v1 = 0.;
1121 v2 = 1.;
1122 if (!didWarn) {
1123 Warning("Divide(h1, h2)", "Histogram h2 has empty bins - division by zero! Setting bin to 0.");
1124 didWarn = kTRUE;
1125 }
1126 }
1127 nFilledBins++;
1128 Long64_t myBin = GetBin(coord);
1129 SetBinContent(myBin, c1 * v1 / c2 / v2);
1130 if(GetCalculateErrors()){
1131 Double_t err1 = h1->GetBinError(i);
1132 Double_t err2 = h2->GetBinError(h2bin);
1133 Double_t errSq = 0.;
1134 if (binomial) {
1135 if (v1 != v2) {
1136 Double_t w = v1 / v2;
1137 err2 *= w;
1138 errSq = TMath::Abs( ( (1. - 2.*w) * err1 * err1 + err2 * err2 ) / (v2 * v2) );
1139 }
1140 } else {
1141 c1 *= c1;
1142 c2 *= c2;
1143 Double_t b22 = v2 * v2 * c2;
1144 err1 *= v2;
1145 err2 *= v1;
1146 errSq = c1 * c2 * (err1 * err1 + err2 * err2) / (b22 * b22);
1147 }
1148 SetBinError2(myBin, errSq);
1149 }
1150 }
1151
1152 delete [] coord;
1153 SetFilledBins(nFilledBins);
1154
1155 // Set as entries in the result histogram the entries in the numerator
1157}
1158
1159////////////////////////////////////////////////////////////////////////////////
1160/// Consistency check on (some of) the parameters of two histograms (for operations).
1161
1162Bool_t THnBase::CheckConsistency(const THnBase *h, const char *tag) const
1163{
1164 if (fNdimensions != h->GetNdimensions()) {
1165 Warning(tag, "Different number of dimensions, cannot carry out operation on the histograms");
1166 return kFALSE;
1167 }
1168 for (Int_t dim = 0; dim < fNdimensions; dim++){
1169 if (GetAxis(dim)->GetNbins() != h->GetAxis(dim)->GetNbins()) {
1170 Warning(tag, "Different number of bins on axis %i, cannot carry out operation on the histograms", dim);
1171 return kFALSE;
1172 }
1173 }
1174 return kTRUE;
1175}
1176
1177////////////////////////////////////////////////////////////////////////////////
1178/// Set the axis # of bins and bin limits on dimension idim
1179
1180void THnBase::SetBinEdges(Int_t idim, const Double_t* bins)
1181{
1182 TAxis* axis = (TAxis*) fAxes[idim];
1183 axis->Set(axis->GetNbins(), bins);
1184}
1185
1186////////////////////////////////////////////////////////////////////////////////
1187/// Change (i.e. set) the title.
1188///
1189/// If title is in the form "stringt;string0;string1;string2 ..."
1190/// the histogram title is set to stringt, the title of axis0 to string0,
1191/// of axis1 to string1, of axis2 to string2, etc, just like it is done
1192/// for TH1/TH2/TH3.
1193/// To insert the character ";" in one of the titles, one should use "#;"
1194/// or "#semicolon".
1195
1196void THnBase::SetTitle(const char *title)
1197{
1198 fTitle = title;
1199 fTitle.ReplaceAll("#;",2,"#semicolon",10);
1200
1201 Int_t endHistTitle = fTitle.First(';');
1202 if (endHistTitle >= 0) {
1203 // title contains a ';' so parse the axis titles
1204 Int_t posTitle = endHistTitle + 1;
1205 Int_t lenTitle = fTitle.Length();
1206 Int_t dim = 0;
1207 while (posTitle > 0 && posTitle < lenTitle && dim < fNdimensions){
1208 Int_t endTitle = fTitle.Index(";", posTitle);
1209 TString axisTitle = fTitle(posTitle, endTitle - posTitle);
1210 axisTitle.ReplaceAll("#semicolon", 10, ";", 1);
1211 GetAxis(dim)->SetTitle(axisTitle);
1212 dim++;
1213 if (endTitle > 0)
1214 posTitle = endTitle + 1;
1215 else
1216 posTitle = -1;
1217 }
1218 // Remove axis titles from histogram title
1219 fTitle.Remove(endHistTitle, lenTitle - endHistTitle);
1220 }
1221
1222 fTitle.ReplaceAll("#semicolon", 10, ";", 1);
1223
1224}
1225
1226////////////////////////////////////////////////////////////////////////////////
1227/// Combine the content of "group" neighboring bins into
1228/// a new bin and return the resulting THnBase.
1229/// For group=2 and a 3 dimensional histogram, all "blocks"
1230/// of 2*2*2 bins will be put into a bin.
1231
1233{
1234 Int_t* ngroup = new Int_t[GetNdimensions()];
1235 for (Int_t d = 0; d < GetNdimensions(); ++d)
1236 ngroup[d] = group;
1237 THnBase* ret = RebinBase(ngroup);
1238 delete [] ngroup;
1239 return ret;
1240}
1241
1242////////////////////////////////////////////////////////////////////////////////
1243/// Combine the content of "group" neighboring bins for each dimension
1244/// into a new bin and return the resulting THnBase.
1245/// For group={2,1,1} and a 3 dimensional histogram, pairs of x-bins
1246/// will be grouped.
1247
1249{
1250 Int_t ndim = GetNdimensions();
1251 TString name(GetName());
1252 for (Int_t d = 0; d < ndim; ++d)
1253 name += TString::Format("_%d", group[d]);
1254
1255
1256 TString title(GetTitle());
1257 Ssiz_t posInsert = title.First(';');
1258 if (posInsert == kNPOS) {
1259 title += " rebin ";
1260 for (Int_t d = 0; d < ndim; ++d)
1261 title += TString::Format("{%d}", group[d]);
1262 } else {
1263 for (Int_t d = ndim - 1; d >= 0; --d)
1264 title.Insert(posInsert, TString::Format("{%d}", group[d]));
1265 title.Insert(posInsert, " rebin ");
1266 }
1267
1268 TObjArray newaxes(ndim);
1269 newaxes.SetOwner();
1270 for (Int_t d = 0; d < ndim; ++d) {
1271 newaxes.AddAt(new TAxis(*GetAxis(d) ),d);
1272 if (group[d] > 1) {
1273 TAxis* newaxis = (TAxis*) newaxes.At(d);
1274 Int_t newbins = (newaxis->GetNbins() + group[d] - 1) / group[d];
1275 if (newaxis->GetXbins() && newaxis->GetXbins()->GetSize()) {
1276 // variable bins
1277 Double_t *edges = new Double_t[newbins + 1];
1278 for (Int_t i = 0; i < newbins + 1; ++i)
1279 if (group[d] * i <= newaxis->GetNbins())
1280 edges[i] = newaxis->GetXbins()->At(group[d] * i);
1281 else edges[i] = newaxis->GetXmax();
1282 newaxis->Set(newbins, edges);
1283 delete [] edges;
1284 } else {
1285 newaxis->Set(newbins, newaxis->GetXmin(), newaxis->GetXmax());
1286 }
1287 }
1288 }
1289
1290 THnBase* h = CloneEmpty(name.Data(), title.Data(), &newaxes, kTRUE);
1291 Bool_t haveErrors = GetCalculateErrors();
1292 Bool_t wantErrors = haveErrors;
1293
1294 Int_t* bins = new Int_t[ndim];
1295 Int_t* coord = new Int_t[fNdimensions];
1296
1297 Long64_t i = 0;
1298 std::unique_ptr<ROOT::Internal::THnBaseBinIter> iter{CreateIter(false)};
1299 while ((i = iter->Next(coord)) >= 0) {
1301 for (Int_t d = 0; d < ndim; ++d) {
1302 bins[d] = TMath::CeilNint( (double) coord[d]/group[d] );
1303 }
1304 Long64_t idxh = h->GetBin(bins, kTRUE /*allocate*/);
1305
1306 if (wantErrors) {
1307 // wantErrors == haveErrors, thus:
1308 h->AddBinError2(idxh, GetBinError2(i));
1309 }
1310
1311 // only _after_ error calculation, or sqrt(v) is taken into account!
1312 h->AddBinContent(idxh, v);
1313 }
1314
1315 delete [] bins;
1316 delete [] coord;
1317
1318 h->SetEntries(fEntries);
1319
1320 return h;
1321
1322}
1323
1324////////////////////////////////////////////////////////////////////////////////
1325/// Clear the histogram
1326
1327void THnBase::ResetBase(Option_t * /*option = ""*/)
1328{
1329 fEntries = 0.;
1330 fTsumw = 0.;
1331 fTsumw2 = -1.;
1332 if (fIntegralStatus != kNoInt) {
1333 fIntegral.clear();
1335 }
1336}
1337
1338////////////////////////////////////////////////////////////////////////////////
1339/// \brief Compute integral (sum of counts) of histogram in all dimensions
1340/// \param respectAxisRange if false, count all bins including under/overflows,
1341/// if true, restrict sum to the user-set axis range
1342/// \sa Projection(0)::Integral() as alternative
1343/// \note this function is different from ComputeIntegral, that is a normalized
1344/// cumulative sum
1345Double_t THnBase::Integral(Bool_t respectAxisRange) const
1346{
1347 Long64_t myLinBin = 0;
1348 std::unique_ptr<ROOT::Internal::THnBaseBinIter> iter{CreateIter(respectAxisRange)};
1349 Double_t sum = 0.;
1350 while ((myLinBin = iter->Next()) >= 0) {
1351 sum += GetBinContent(myLinBin);
1352 }
1353 return sum;
1354}
1355
1356////////////////////////////////////////////////////////////////////////////////
1357/// Compute integral (normalized cumulative sum of bins) w/o under/overflows
1358/// The result is stored in fIntegral and used by the GetRandom functions.
1359/// This function is automatically called by GetRandom when the fIntegral
1360/// array does not exist or when the number of entries in the histogram
1361/// has changed since the previous call to GetRandom.
1362/// The resulting integral is normalized to 1.
1363/// \return 1 if success, 0 if integral is zero
1364
1366{
1367 // delete old integral
1368 if (fIntegralStatus != kNoInt) {
1369 fIntegral.clear();
1371 }
1372
1373 // check number of bins
1374 if (GetNbins() == 0) {
1375 Error("ComputeIntegral", "The histogram must have at least one bin.");
1376 return 0.;
1377 }
1378
1379 // allocate integral array
1380 fIntegral.resize(GetNbins() + 1);
1381 fIntegral[0] = 0.;
1382
1383 // fill integral array with contents of regular bins (non over/underflow)
1384 Int_t* coord = new Int_t[fNdimensions];
1385 Long64_t i = 0;
1386 std::unique_ptr<ROOT::Internal::THnBaseBinIter> iter{CreateIter(false)};
1387 while ((i = iter->Next(coord)) >= 0) {
1389
1390 // check whether the bin is regular
1391 bool regularBin = true;
1392 for (Int_t dim = 0; dim < fNdimensions; dim++) {
1393 if (coord[dim] < 1 || coord[dim] > GetAxis(dim)->GetNbins()) {
1394 regularBin = false;
1395 break;
1396 }
1397 }
1398
1399 // if outlier, count it with zero weight
1400 if (!regularBin) v = 0.;
1401
1402 fIntegral[i + 1] = fIntegral[i] + v;
1403 }
1404 delete [] coord;
1405
1406 // check sum of weights
1407 if (fIntegral[GetNbins()] == 0.) {
1408 Error("ComputeIntegral", "Integral = 0, no hits in histogram bins (excluding over/underflow).");
1409 fIntegral.clear();
1410 return 0.;
1411 }
1412
1413 // normalize the integral array
1414 for (Long64_t j = 0; j <= GetNbins(); ++j)
1415 fIntegral[j] /= fIntegral[GetNbins()];
1416
1417 // set status to valid
1419 return fIntegral[GetNbins()];
1420}
1421
1422////////////////////////////////////////////////////////////////////////////////
1423/// Print bin with linex index "idx".
1424/// For valid options see PrintBin(Long64_t idx, Int_t* bin, Option_t* options).
1425
1426void THnBase::PrintBin(Long64_t idx, Option_t* options) const
1427{
1428 Int_t* coord = new Int_t[fNdimensions];
1429 PrintBin(idx, coord, options);
1430 delete [] coord;
1431}
1432
1433////////////////////////////////////////////////////////////////////////////////
1434/// Print one bin. If "idx" is != -1 use that to determine the bin,
1435/// otherwise (if "idx" == -1) use the coordinate in "bin".
1436/// If "options" contains:
1437/// - '0': only print bins with an error or content != 0
1438/// Return whether the bin was printed (depends on options)
1439
1441{
1442 Double_t v = -42;
1443 if (idx == -1) {
1444 idx = GetBin(bin);
1445 v = GetBinContent(idx);
1446 } else {
1447 v = GetBinContent(idx, bin);
1448 }
1449
1450 Double_t err = 0.;
1451 if (GetCalculateErrors()) {
1452 if (idx != -1) {
1453 err = GetBinError(idx);
1454 }
1455 }
1456
1457 if (v == 0. && err == 0. && options && strchr(options, '0')) {
1458 // suppress zeros, and we have one.
1459 return kFALSE;
1460 }
1461
1462 TString coord;
1463 for (Int_t dim = 0; dim < fNdimensions; ++dim) {
1464 coord += bin[dim];
1465 coord += ',';
1466 }
1467 coord.Remove(coord.Length() - 1);
1468
1469 if (GetCalculateErrors()) {
1470 Printf("Bin at (%s) = %g (+/- %g)", coord.Data(), v, err);
1471 } else {
1472 Printf("Bin at (%s) = %g", coord.Data(), v);
1473 }
1474
1475 return kTRUE;
1476}
1477
1478////////////////////////////////////////////////////////////////////////////////
1479/// Print "howmany" entries starting at "from". If "howmany" is -1, print all.
1480/// If "options" contains:
1481/// - 'x': print in the order of axis bins, i.e. (0,0,...,0), (0,0,...,1),...
1482/// - '0': only print bins with content != 0
1483
1484void THnBase::PrintEntries(Long64_t from /*=0*/, Long64_t howmany /*=-1*/,
1485 Option_t* options /*=0*/) const
1486{
1487 if (from < 0) from = 0;
1488 if (howmany == -1) howmany = GetNbins();
1489
1490 Int_t* bin = new Int_t[fNdimensions];
1491
1492 if (options && (strchr(options, 'x') || strchr(options, 'X'))) {
1493 Int_t* nbins = new Int_t[fNdimensions];
1494 for (Int_t dim = fNdimensions - 1; dim >= 0; --dim) {
1495 nbins[dim] = GetAxis(dim)->GetNbins();
1496 bin[dim] = from % nbins[dim];
1497 from /= nbins[dim];
1498 }
1499
1500 for (Long64_t i = 0; i < howmany; ++i) {
1501 if (!PrintBin(-1, bin, options) || !strchr(options, '0'))
1502 ++howmany;
1503 // Advance to next bin in last dimension:
1504 ++bin[fNdimensions - 1];
1505 // If overflow, fix it by adding to previous dimension
1506 for (Int_t dim = fNdimensions - 1; dim >= 0; --dim) {
1507 if (bin[dim] >= nbins[dim] + 2) { // include under/overflow bins
1508 bin[dim] = 0;
1509 if (dim > 0) {
1510 ++bin[dim - 1];
1511 } else {
1512 howmany = -1; // aka "global break"
1513 }
1514 }
1515 }
1516 }
1517 delete [] nbins;
1518 } else {
1519 for (Long64_t i = from; i < from + howmany; ++i) {
1520 if (!PrintBin(i, bin, options))
1521 ++howmany;
1522 }
1523 }
1524 delete [] bin;
1525}
1526
1527////////////////////////////////////////////////////////////////////////////////
1528/// Print a THnBase. If "option" contains:
1529/// - 'a': print axis details
1530/// - 'm': print memory usage
1531/// - 's': print statistics
1532/// - 'c': print its content, too (this can generate a LOT of output!)
1533/// Other options are forwarded to PrintEntries().
1534
1535void THnBase::Print(Option_t* options) const
1536{
1537 Bool_t optAxis = options && (strchr(options, 'A') || (strchr(options, 'a')));
1538 Bool_t optMem = options && (strchr(options, 'M') || (strchr(options, 'm')));
1539 Bool_t optStat = options && (strchr(options, 'S') || (strchr(options, 's')));
1540 Bool_t optContent = options && (strchr(options, 'C') || (strchr(options, 'c')));
1541
1542 Printf("%s (*0x%zx): \"%s\" \"%s\"", IsA()->GetName(), (size_t)this, GetName(), GetTitle());
1543 Printf(" %d dimensions, %g entries in %lld filled bins", GetNdimensions(), GetEntries(), GetNbins());
1544
1545 if (optAxis) {
1546 for (Int_t dim = 0; dim < fNdimensions; ++dim) {
1547 TAxis* axis = GetAxis(dim);
1548 Printf(" axis %d \"%s\": %d bins (%g..%g), %s bin sizes", dim,
1549 axis->GetTitle(), axis->GetNbins(), axis->GetXmin(), axis->GetXmax(),
1550 (axis->GetXbins() ? "variable" : "fixed"));
1551 }
1552 }
1553
1554 if (optStat) {
1555 Printf(" %s error calculation", (GetCalculateErrors() ? "with" : "without"));
1556 if (GetCalculateErrors()) {
1557 Printf(" Sum(w)=%g, Sum(w^2)=%g", GetSumw(), GetSumw2());
1558 for (Int_t dim = 0; dim < fNdimensions; ++dim) {
1559 Printf(" axis %d: Sum(w*x)=%g, Sum(w*x^2)=%g", dim, GetSumwx(dim), GetSumwx2(dim));
1560 }
1561 }
1562 }
1563
1564 if (optMem && InheritsFrom(THnSparse::Class())) {
1565 const THnSparse* hsparse = dynamic_cast<const THnSparse*>(this);
1566 Printf(" coordinates stored in %d chunks of %d entries\n %g of bins filled using %g of memory compared to an array",
1567 hsparse->GetNChunks(), hsparse->GetChunkSize(),
1568 hsparse->GetSparseFractionBins(), hsparse->GetSparseFractionMem());
1569 }
1570
1571 if (optContent) {
1572 Printf(" BIN CONTENT:");
1573 PrintEntries(0, -1, options);
1574 }
1575}
1576
1577
1578////////////////////////////////////////////////////////////////////////////////
1579/// Browse a THnSparse: create an entry (ROOT::THnSparseBrowsable) for each
1580/// dimension.
1581
1583{
1584 if (fBrowsables.IsEmpty()) {
1585 for (Int_t dim = 0; dim < fNdimensions; ++dim) {
1587 }
1589 }
1590
1591 for (Int_t dim = 0; dim < fNdimensions; ++dim) {
1592 b->Add(fBrowsables[dim]);
1593 }
1594}
1595
1596
1597
1598
1599/** \class ROOT::Internal::THnBaseBinIter
1600 Iterator over THnBase bins (internal implementation).
1601*/
1602
1603/// Destruct a bin iterator.
1604
1606 // Not much to do, but pin vtable
1607}
1608
1609
1610
1611/** \class THnIter
1612 Iterator over THnBase bins
1613*/
1614
1616
1618 // Destruct a bin iterator.
1619 delete fIter;
1620}
1621
1622
1623/** \class ROOT::Internal::THnBaseBrowsable
1624 TBrowser helper for THnBase.
1625*/
1626
1627
1629
1630////////////////////////////////////////////////////////////////////////////////
1631/// Construct a THnBaseBrowsable.
1632
1634fHist(hist), fAxis(axis), fProj(nullptr)
1635{
1636 TString axisName = hist->GetAxis(axis)->GetName();
1637 if (axisName.IsNull()) {
1638 axisName = TString::Format("axis%d", axis);
1639 }
1640
1641 SetNameTitle(axisName,
1642 TString::Format("Projection on %s of %s", axisName.Data(),
1643 hist->IsA()->GetName()).Data());
1644}
1645
1646////////////////////////////////////////////////////////////////////////////////
1647/// Destruct a THnBaseBrowsable.
1648
1650{
1651 delete fProj;
1652}
1653
1654////////////////////////////////////////////////////////////////////////////////
1655/// Browse an axis of a THnBase, i.e. draw its projection.
1656
1658{
1659 if (!fProj) {
1660 fProj = fHist->Projection(fAxis);
1661 }
1662 fProj->Draw(b ? b->GetDrawOption() : "");
1663 gPad->Update();
1664}
1665
#define d(i)
Definition RSha256.hxx:102
#define b(i)
Definition RSha256.hxx:100
#define f(i)
Definition RSha256.hxx:104
#define c(i)
Definition RSha256.hxx:101
#define h(i)
Definition RSha256.hxx:106
size_t size(const MatrixT &matrix)
retrieve the size of a square matrix
int Int_t
Definition RtypesCore.h:45
constexpr Bool_t kFALSE
Definition RtypesCore.h:94
constexpr Ssiz_t kNPOS
Definition RtypesCore.h:117
long long Long64_t
Definition RtypesCore.h:69
constexpr Bool_t kTRUE
Definition RtypesCore.h:93
const char Option_t
Definition RtypesCore.h:66
#define ClassImp(name)
Definition Rtypes.h:382
winID h TVirtualViewer3D TVirtualGLPainter p
Option_t Option_t option
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void data
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t Float_t Float_t Float_t Int_t Int_t UInt_t UInt_t Rectangle_t Int_t Int_t Window_t TString Int_t GCValues_t GetPrimarySelectionOwner GetDisplay GetScreen GetColormap GetNativeEvent const char const char dpyName wid window const char font_name cursor keysym reg const char only_if_exist regb h Point_t winding char text const char depth char const char Int_t count const char cname
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void value
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t Float_t Float_t Float_t Int_t Int_t UInt_t UInt_t Rectangle_t Int_t Int_t Window_t TString Int_t GCValues_t GetPrimarySelectionOwner GetDisplay GetScreen GetColormap GetNativeEvent const char const char dpyName wid window const char font_name cursor keysym reg const char only_if_exist regb h Point_t winding char text const char depth char const char Int_t count const char ColorStruct_t color const char Pixmap_t Pixmap_t PictureAttributes_t attr const char char ret_data h unsigned char height h Atom_t Int_t ULong_t ULong_t unsigned char prop_list Atom_t Atom_t Atom_t Time_t type
char name[80]
Definition TGX11.cxx:110
float xmin
float xmax
#define R__THNBCASE(TAG)
R__EXTERN TRandom * gRandom
Definition TRandom.h:62
void Printf(const char *fmt,...)
Formats a string in a circular formatting buffer and prints the string.
Definition TString.cxx:2503
#define gPad
class describing the range in the coordinates it supports multiple range in a coordinate.
Definition DataRange.h:35
void AddRange(unsigned int icoord, double xmin, double xmax)
add a range [xmin,xmax] for the new coordinate icoord Adding a range does not delete existing one,...
Definition DataRange.cxx:94
virtual ~THnBaseBinIter()
Destruct a bin iterator.
Definition THnBase.cxx:1605
TBrowser helper for THnBase.
Definition THnBase.h:300
~THnBaseBrowsable() override
Destruct a THnBaseBrowsable.
Definition THnBase.cxx:1649
THnBaseBrowsable(THnBase *hist, Int_t axis)
Construct a THnBaseBrowsable.
Definition THnBase.cxx:1633
void Browse(TBrowser *b) override
Browse an axis of a THnBase, i.e. draw its projection.
Definition THnBase.cxx:1657
Double_t At(Int_t i) const
Definition TArrayD.h:79
void Set(Int_t n) override
Set size of this array to n doubles.
Definition TArrayD.cxx:106
const Double_t * GetArray() const
Definition TArrayD.h:43
Abstract array base class.
Definition TArray.h:31
Int_t GetSize() const
Definition TArray.h:47
Class to manage histogram axis.
Definition TAxis.h:31
virtual void SetBinLabel(Int_t bin, const char *label)
Set label for bin.
Definition TAxis.cxx:886
const char * GetTitle() const override
Returns title of object.
Definition TAxis.h:135
virtual Double_t GetBinCenter(Int_t bin) const
Return center of bin.
Definition TAxis.cxx:478
const TArrayD * GetXbins() const
Definition TAxis.h:136
Double_t GetXmax() const
Definition TAxis.h:140
@ kAxisRange
Definition TAxis.h:65
virtual Double_t GetBinLowEdge(Int_t bin) const
Return low edge of bin.
Definition TAxis.cxx:518
virtual void Set(Int_t nbins, Double_t xmin, Double_t xmax)
Initialize axis with fix bins.
Definition TAxis.cxx:794
Int_t GetLast() const
Return last bin on the axis i.e.
Definition TAxis.cxx:469
Double_t GetXmin() const
Definition TAxis.h:139
Int_t GetNbins() const
Definition TAxis.h:125
virtual void SetRange(Int_t first=0, Int_t last=0)
Set the viewing range for the axis using bin numbers.
Definition TAxis.cxx:1052
virtual Double_t GetBinWidth(Int_t bin) const
Return bin width.
Definition TAxis.cxx:540
virtual Double_t GetBinUpEdge(Int_t bin) const
Return up edge of bin.
Definition TAxis.cxx:528
Int_t GetFirst() const
Return first bin on the axis i.e.
Definition TAxis.cxx:458
THashList * GetLabels() const
Definition TAxis.h:121
Using a TBrowser one can browse all ROOT objects.
Definition TBrowser.h:37
TClass instances represent classes, structs and namespaces in the ROOT type system.
Definition TClass.h:81
void * New(ENewType defConstructor=kClassNew, Bool_t quiet=kFALSE) const
Return a pointer to a newly allocated object of this class.
Definition TClass.cxx:5047
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:3037
Collection abstract base class.
Definition TCollection.h:65
virtual void SetOwner(Bool_t enable=kTRUE)
Set whether this collection is the owner (enable==true) of its content.
virtual Bool_t IsEmpty() const
virtual Int_t GetSize() const
Return the capacity of the collection, i.e.
1-Dim function class
Definition TF1.h:233
static void RejectPoint(Bool_t reject=kTRUE)
Static function to set the global flag to reject points the fgRejectPoint global flag is tested by al...
Definition TF1.cxx:3683
Provides an indirection to the TFitResult class and with a semantics identical to a TFitResult pointe...
1-D histogram with a double per channel (see TH1 documentation)
Definition TH1.h:670
TH1 is the base class of all histogram classes in ROOT.
Definition TH1.h:59
virtual Double_t GetEffectiveEntries() const
Number of effective entries of the histogram.
Definition TH1.cxx:4448
virtual void Rebuild(Option_t *option="")
Using the current bin info, recompute the arrays for contents and errors.
Definition TH1.cxx:7108
TAxis * GetZaxis()
Definition TH1.h:326
virtual void AddBinContent(Int_t bin)
Increment bin content by 1.
Definition TH1.cxx:1268
virtual Double_t GetBinError(Int_t bin) const
Return value of error associated to bin number bin.
Definition TH1.cxx:9084
static Int_t FitOptionsMake(Option_t *option, Foption_t &Foption)
Decode string choptin and fill fitOption structure.
Definition TH1.cxx:4673
TAxis * GetXaxis()
Definition TH1.h:324
virtual void GetBinXYZ(Int_t binglobal, Int_t &binx, Int_t &biny, Int_t &binz) const
Return binx, biny, binz corresponding to the global bin number globalbin see TH1::GetBin function abo...
Definition TH1.cxx:4995
virtual Int_t GetNcells() const
Definition TH1.h:300
virtual Int_t GetBin(Int_t binx, Int_t biny=0, Int_t binz=0) const
Return Global bin number corresponding to binx,y,z.
Definition TH1.cxx:4982
virtual void SetBinError(Int_t bin, Double_t error)
Set the bin Error Note that this resets the bin eror option to be of Normal Type and for the non-empt...
Definition TH1.cxx:9227
TAxis * GetYaxis()
Definition TH1.h:325
virtual Double_t GetEntries() const
Return the current number of entries.
Definition TH1.cxx:4423
virtual void ResetStats()
Reset the statistics including the number of entries and replace with values calculated from bin cont...
Definition TH1.cxx:7923
virtual Double_t GetBinContent(Int_t bin) const
Return content of bin number bin.
Definition TH1.cxx:5082
virtual void SetEntries(Double_t n)
Definition TH1.h:391
2-D histogram with a double per channel (see TH1 documentation)
Definition TH2.h:357
3-D histogram with a double per channel (see TH1 documentation)
Definition TH3.h:363
THashList implements a hybrid collection class consisting of a hash table and a list to store TObject...
Definition THashList.h:34
Multidimensional histogram base.
Definition THnBase.h:43
virtual void Sumw2()=0
virtual ROOT::Internal::THnBaseBinIter * CreateIter(Bool_t respectAxisRange) const =0
void SetEntries(Double_t entries)
Definition THnBase.h:192
virtual void SetFilledBins(Long64_t)
Definition THnBase.h:105
void Browse(TBrowser *b) override
Browse a THnSparse: create an entry (ROOT::THnSparseBrowsable) for each dimension.
Definition THnBase.cxx:1582
virtual void InitStorage(Int_t *nbins, Int_t chunkSize)=0
Double_t GetBinError(const Int_t *idx) const
Definition THnBase.h:187
void SetBinError(const Int_t *idx, Double_t e)
Definition THnBase.h:189
Double_t fEntries
Number of entries, spread over chunks.
Definition THnBase.h:48
Double_t GetSumw2() const
Definition THnBase.h:214
Bool_t IsInRange(Int_t *coord) const
Check whether bin coord is in range, as defined by TAxis::SetRange().
Definition THnBase.cxx:592
TFitResultPtr Fit(TF1 *f1, Option_t *option="", Option_t *goption="")
Fit a THnSparse with function f.
Definition THnBase.cxx:505
void AddBinContent(const Int_t *x, Double_t v=1.)
Definition THnBase.h:191
void Scale(Double_t c)
Scale contents and errors of this histogram by c: this = this * c It does not modify the histogram's ...
Definition THnBase.cxx:774
virtual void SetBinError2(Long64_t bin, Double_t e2)=0
TObjArray * GetListOfAxes()
Definition THnBase.h:128
TH1 * CreateHist(const char *name, const char *title, const TObjArray *axes, Bool_t keepTargetAxis) const
Create an empty histogram with name and title with a given set of axes.
Definition THnBase.cxx:252
Bool_t PrintBin(Long64_t idx, Int_t *coord, Option_t *options) const
Print one bin.
Definition THnBase.cxx:1440
Double_t GetSumwx(Int_t dim) const
Definition THnBase.h:215
@ kNoInt
Definition THnBase.h:55
@ kValidInt
Definition THnBase.h:56
Bool_t CheckConsistency(const THnBase *h, const char *tag) const
Consistency check on (some of) the parameters of two histograms (for operations).
Definition THnBase.cxx:1162
Int_t GetNdimensions() const
Definition THnBase.h:140
~THnBase() override
Destruct a THnBase.
Definition THnBase.cxx:167
static THnBase * CreateHnAny(const char *name, const char *title, const TH1 *h1, Bool_t sparse, Int_t chunkSize=1024 *16)
Create a THn / THnSparse object from a histogram deriving from TH1.
Definition THnBase.cxx:323
TArrayD fTsumwx2
Total sum of weight*X*X for each dimension.
Definition THnBase.h:52
void ResetBase(Option_t *option="")
Clear the histogram.
Definition THnBase.cxx:1327
virtual Long64_t GetNbins() const =0
THnBase * RebinBase(Int_t group) const
Combine the content of "group" neighboring bins into a new bin and return the resulting THnBase.
Definition THnBase.cxx:1232
TClass * IsA() const override
Definition THnBase.h:292
TObjArray fAxes
Axes of the histogram.
Definition THnBase.h:46
void SetBinEdges(Int_t idim, const Double_t *bins)
Set the axis # of bins and bin limits on dimension idim.
Definition THnBase.cxx:1180
Bool_t GetCalculateErrors() const
Definition THnBase.h:141
void AddInternal(const THnBase *h, Double_t c, Bool_t rebinned)
Add() implementation for both rebinned histograms and those with identical binning.
Definition THnBase.cxx:798
void SetTitle(const char *title) override
Change (i.e.
Definition THnBase.cxx:1196
void PrintEntries(Long64_t from=0, Long64_t howmany=-1, Option_t *options=nullptr) const
Print "howmany" entries starting at "from".
Definition THnBase.cxx:1484
enum THnBase::@83 fIntegralStatus
! status of integral
Double_t GetBinContent(const Int_t *idx) const
Definition THnBase.h:197
void Add(const THnBase *h, Double_t c=1.)
Add contents of h scaled by c to this histogram: this = this + c * h Note that if h has Sumw2 set,...
Definition THnBase.cxx:871
Double_t fTsumw2
Total sum of weights squared; -1 if no errors are calculated.
Definition THnBase.h:50
virtual void Reserve(Long64_t)
Definition THnBase.h:104
THnBase()
Definition THnBase.h:61
Double_t GetEntries() const
Definition THnBase.h:138
std::vector< Double_t > GetBinCenter(const std::vector< Int_t > &idx) const
THnBase::GetBinCenter.
Definition THnBase.cxx:534
virtual Double_t GetBinError2(Long64_t linidx) const =0
virtual void Reset(Option_t *option="")=0
void SetBinContent(const Int_t *idx, Double_t v)
Definition THnBase.h:207
void Multiply(const THnBase *h)
Multiply this histogram by histogram h this = this * h Note that if h has Sumw2 set,...
Definition THnBase.cxx:932
Double_t Integral(Bool_t respectAxisRange) const
Compute integral (sum of counts) of histogram in all dimensions.
Definition THnBase.cxx:1345
void Divide(const THnBase *h)
Divide this histogram by h this = this/(h) Note that if h has Sumw2 set, Sumw2 is automatically calle...
Definition THnBase.cxx:1023
void Print(Option_t *option="") const override
Print a THnBase.
Definition THnBase.cxx:1535
Double_t GetSumwx2(Int_t dim) const
Definition THnBase.h:216
TObjArray fBrowsables
! Browser-helpers for each axis
Definition THnBase.h:47
Double_t GetSumw() const
Definition THnBase.h:213
TObject * ProjectionAny(Int_t ndim, const Int_t *dim, Bool_t wantNDim, Option_t *option="") const
Project all bins into a ndim-dimensional THn / THnSparse (whatever *this is) or if (ndim < 4 and !...
Definition THnBase.cxx:618
void GetRandom(Double_t *rand, Bool_t subBinRandom=kTRUE)
Generate an n-dimensional random tuple based on the histogrammed distribution.
Definition THnBase.cxx:557
TAxis * GetAxis(Int_t dim) const
Definition THnBase.h:130
THnBase * CloneEmpty(const char *name, const char *title, const TObjArray *axes, Bool_t keepTargetAxis) const
Create a new THnBase object that is of the same type as *this, but with dimensions and bins given by ...
Definition THnBase.cxx:179
Long64_t Merge(TCollection *list)
Merge this with a list of THnBase's.
Definition THnBase.cxx:897
Double_t fTsumw
Total sum of weights.
Definition THnBase.h:49
virtual Long64_t GetBin(const Int_t *idx) const =0
Double_t ComputeIntegral()
Compute integral (normalized cumulative sum of bins) w/o under/overflows The result is stored in fInt...
Definition THnBase.cxx:1365
virtual void AddBinError2(Long64_t bin, Double_t e2)=0
void RebinnedAdd(const THnBase *h, Double_t c=1.)
Add contents of h scaled by c to this histogram: this = this + c * h Note that if h has Sumw2 set,...
Definition THnBase.cxx:887
Int_t fNdimensions
Number of dimensions.
Definition THnBase.h:45
void Init(const char *name, const char *title, const TObjArray *axes, Bool_t keepTargetAxis, Int_t chunkSize=1024 *16)
Initialize axes and name.
Definition THnBase.cxx:195
THnBase & operator=(const THnBase &other)
Definition THnBase.cxx:101
std::vector< Double_t > fIntegral
! vector with bin weight sums
Definition THnBase.h:53
TArrayD fTsumwx
Total sum of weight*X for each dimension.
Definition THnBase.h:51
Iterator over THnBase bins.
Definition THnBase.h:333
~THnIter() override
Definition THnBase.cxx:1617
static TClass * Class()
Efficient multidimensional histogram.
Definition THnSparse.h:37
Double_t GetSparseFractionBins() const
Return the amount of filled bins over all bins.
Double_t GetSparseFractionMem() const
Return the amount of used memory over memory that would be used by a non-sparse n-dimensional histogr...
Int_t GetNChunks() const
Definition THnSparse.h:89
Int_t GetChunkSize() const
Definition THnSparse.h:88
static TClass * Class()
static TClass * Class()
static TClass * Class()
void Reset()
The TNamed class is the base class for all named ROOT classes.
Definition TNamed.h:29
virtual void SetTitle(const char *title="")
Set the title of the TNamed.
Definition TNamed.cxx:164
const char * GetName() const override
Returns name of object.
Definition TNamed.h:47
const char * GetTitle() const override
Returns title of object.
Definition TNamed.h:48
TString fTitle
Definition TNamed.h:33
void Clear(Option_t *option="") override
Set name and title to empty strings ("").
Definition TNamed.cxx:64
virtual void SetName(const char *name)
Set the name of the TNamed.
Definition TNamed.cxx:140
TNamed & operator=(const TNamed &rhs)
TNamed assignment operator.
Definition TNamed.cxx:51
virtual void SetNameTitle(const char *name, const char *title)
Set all the TNamed parameters (name and title).
Definition TNamed.cxx:154
An array of TObjects.
Definition TObjArray.h:31
Int_t GetEntriesFast() const
Definition TObjArray.h:58
void AddAt(TObject *obj, Int_t idx) override
Add object at position ids.
virtual void AddAtAndExpand(TObject *obj, Int_t idx)
Add object at position idx.
TObject * At(Int_t idx) const override
Definition TObjArray.h:164
Bool_t IsEmpty() const override
Definition TObjArray.h:65
Collectable string class.
Definition TObjString.h:28
Mother of all ROOT objects.
Definition TObject.h:41
virtual const char * GetName() const
Returns name of object.
Definition TObject.cxx:444
R__ALWAYS_INLINE Bool_t TestBit(UInt_t f) const
Definition TObject.h:199
virtual void Warning(const char *method, const char *msgfmt,...) const
Issue warning message.
Definition TObject.cxx:979
void SetBit(UInt_t f, Bool_t set)
Set or unset the user status bits as specified in f.
Definition TObject.cxx:786
virtual Bool_t InheritsFrom(const char *classname) const
Returns kTRUE if object inherits from class "classname".
Definition TObject.cxx:530
virtual void Error(const char *method, const char *msgfmt,...) const
Issue error message.
Definition TObject.cxx:993
void ResetBit(UInt_t f)
Definition TObject.h:198
Double_t Rndm() override
Machine independent random number generator.
Definition TRandom.cxx:559
Basic string class.
Definition TString.h:139
Ssiz_t Length() const
Definition TString.h:417
void ToLower()
Change string to lower-case.
Definition TString.cxx:1182
TString & Insert(Ssiz_t pos, const char *s)
Definition TString.h:661
Ssiz_t First(char c) const
Find first occurrence of a character c.
Definition TString.cxx:538
const char * Data() const
Definition TString.h:376
TString & ReplaceAll(const TString &s1, const TString &s2)
Definition TString.h:704
Bool_t IsNull() const
Definition TString.h:414
TString & Remove(Ssiz_t pos)
Definition TString.h:685
static TString Format(const char *fmt,...)
Static method which formats a string using a printf style format descriptor and return a TString.
Definition TString.cxx:2378
Bool_t Contains(const char *pat, ECaseCompare cmp=kExact) const
Definition TString.h:632
Ssiz_t Index(const char *pat, Ssiz_t i=0, ECaseCompare cmp=kExact) const
Definition TString.h:651
return c1
Definition legend1.C:41
Double_t x[n]
Definition legend1.C:17
TH1F * h1
Definition legend1.C:5
return c2
Definition legend2.C:14
#define F(x, y, z)
#define I(x, y, z)
TFitResultPtr FitObject(TH1 *h1, TF1 *f1, Foption_t &option, const ROOT::Math::MinimizerOptions &moption, const char *goption, ROOT::Fit::DataRange &range)
fitting function for a TH1 (called from TH1::Fit)
Definition HFitImpl.cxx:972
Double_t Floor(Double_t x)
Rounds x downward, returning the largest integral value that is not greater than x.
Definition TMath.h:680
Double_t Sqrt(Double_t x)
Returns the square root of x.
Definition TMath.h:662
Int_t CeilNint(Double_t x)
Returns the nearest integer of TMath::Ceil(x).
Definition TMath.h:674
Long64_t BinarySearch(Long64_t n, const T *array, T value)
Binary search in an array of n values to locate value.
Definition TMathBase.h:347
Short_t Abs(Short_t d)
Returns the absolute value of parameter Short_t d.
Definition TMathBase.h:123
int Like
Definition Foption.h:34
int Nostore
Definition Foption.h:41
int Chi2
Definition Foption.h:32
static uint64_t sum(uint64_t i)
Definition Factory.cxx:2345