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