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