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