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