source: trunk/Mars/mhbase/MH3.cc@ 10009

Last change on this file since 10009 was 9858, checked in by tbretz, 14 years ago
Implemented TProfile3D in MH3 and MHn and the possibility to skip the reset of the histograms in the SetupFill.
File size: 41.2 KB
Line 
1/* ======================================================================== *\
2!
3! *
4! * This file is part of MARS, the MAGIC Analysis and Reconstruction
5! * Software. It is distributed to you in the hope that it can be a useful
6! * and timesaving tool in analysing Data of imaging Cerenkov telescopes.
7! * It is distributed WITHOUT ANY WARRANTY.
8! *
9! * Permission to use, copy, modify and distribute this software and its
10! * documentation for any purpose is hereby granted without fee,
11! * provided that the above copyright notice appear in all copies and
12! * that both that copyright notice and this permission notice appear
13! * in supporting documentation. It is provided "as is" without express
14! * or implied warranty.
15! *
16!
17!
18! Author(s): Thomas Bretz 2002 <mailto:tbretz@astro.uni-wuerzburg.de>
19!
20! Copyright: MAGIC Software Development, 2000-2010
21!
22!
23\* ======================================================================== */
24
25/////////////////////////////////////////////////////////////////////////////
26//
27// MH3
28//
29// With this histogram you can fill a histogram with up to three
30// variables from Mars parameter containers in an eventloop.
31//
32// In the constructor you can give up to three variables which should be
33// filled in the histogram. Dependend on the number of given variables
34// (data members) a TH1D, TH2D or TH3D is created.
35// Specify the data mamber like the following:
36// "MHillas.fLength"
37// Where MHillas is the name of the parameter container in the parameter
38// list and fLength is the name of the data member which should be filled
39// in the histogram. Assuming that your MHillas container has a different
40// name (MyHillas) the name to give would be:
41// "MyHillas.fLength"
42//
43// If you want to use a different unit for histogramming use SetScaleX,
44// SetScaleY and SetScaleZ.
45//
46//
47// Binning/Binning name
48// =====================
49//
50// The axis binning is retrieved from the parameter list, too. Create a
51// MBinning with the name "Binning" plus the name of your MH3 container
52// plus the axis name ("X", "Y" or "Z") and add it to the parameter list.
53//
54// If the binning should have a different name than the histogram name
55// the binning name can be added to the name, eg.:
56// SetName("MyHistName;MyXBins;MyYBins")
57// Instead of BinningMyHistName[XYZ] the parameter list will be searched
58// for BinningMyXBinning, BinningMyYBins and BinningMyHistNameZ
59//
60// If you don't want to use a MBinning object from the parameter list
61// you can also set one directly, for example
62// MBinning bins(10, 0, 1);
63// mh3.SetBinningX(&bins);
64// You must not delete the MBinning object before the class has been
65// PreProcessed.
66//
67//
68// Axis titles
69// ===========
70//
71// 1) If no other title is given the rule for this axis is used.
72// 2) If the MBinning used for this axis has a non-default Title
73// (MBinning::HasTitle) this title is used for the corresponding axis
74// 3) If the MH3 has a non-default title (MH3::SetTitle called)
75// this title is set as the histogram title. It can be used to overwrite
76// the axis titles. For more information see TH1::SetTitle, eg.
77// SetTitle("MyHist;x[mm];y[cm];Counts");
78//
79//
80// Labels
81// ======
82//
83// To use labels at an axis you have to initialize this for the axis
84// by either calling InitLabels(Labels_t) or setiting a DefaultLabel.
85// For the axis for which the labels have been initialized the
86// number returned by the given corresponding phrase is converted
87// to int with TMath::Nint and used as a label. If you want to replace
88// this id by a named label you can call DefineLabel to do that.
89// Several ids can be replaced by the same label. If you define
90// named labels for every label which was not defined the default
91// is used, if any, otherwise an unnamed label is created.
92//
93// In the case of an axis with labels the axis-title cannot be
94// set via a MBinning, because the MBinning is not evaluated.
95//
96// Please note that for some reason not all combinations of
97// labels, dimensions and weights are available in the root-
98// histogram classes. Please check the MH3::Fill function to see
99// whether your combination is supported.
100//
101//
102// Examples:
103// =========
104//
105// 1) MH3 myhist("MHillas.fLength");
106// myhist.SetName("MyHist");
107// myhist.SetScaleX(geomcam.GetConvMm2Deg()); //convert length to degree
108// MBinning bins("BinningMyHistX", "Title for my x-axis [Hz]");
109// bins.SetEdges(10, 0, 150);
110// plist.AddToList(&bins);
111//
112// 2) MH3 myhist("MHillas.fLength");
113// myhist.SetName("MyHist;MyX");
114// myhist.SetTitle("Histogram Title;X-Title [mm];Counts");
115// MBinning bins("BinningMyX");
116// bins.SetEdges(10, 0, 150);
117// plist.AddToList(&bins);
118//
119// 3) MH3 myhist("MTriggerPatter.GetUnprescaled");
120// myhist.SetWeight("1./MRawRunHeader.GetRunLength");
121// myhist.SetTitle("Rate of the trigger pattern [Hz];Run Number;Trigger Pattern;Rate [Hz]");
122// myhist.InitLabels(MH3::kLabelsXY);
123// myhist.DefaultLabelY("UNKNOWN"); // Lvl1
124// myhist.DefineLabelY( 1, "Trig"); // Lvl1
125// myhist.DefineLabelY( 2, "Cal"); // Cal
126// myhist.DefineLabelY( 4, "Trig"); // Lvl2
127// myhist.DefineLabelY( 8, "Ped"); // Ped
128//
129//
130// Class Version 1:
131// ----------------
132// - MData *fData[3];
133// + MData *fData[4];
134//
135// Class Version 2:
136// ----------------
137// - MDataChain *fData[3]; // Object from which the data is filled
138// + MData *fData[3]; // Object from which the data is filled
139//
140// Class Version 3:
141// ----------------
142// - Byte_t fStyleBits
143// + MBinning fBins[3]
144//
145// Class Version 5:
146// ----------------
147// + TFormula *fConversion
148//
149// Class Version 5:
150// ----------------
151// + MData *fWeight;
152//
153/////////////////////////////////////////////////////////////////////////////
154#include "MH3.h"
155
156#include <ctype.h> // tolower
157#include <stdlib.h> // atoi (Ubuntu 8.10)
158#include <fstream>
159
160#include <TMath.h>
161#include <TFormula.h>
162
163#include <THashList.h>
164#include <TObjString.h>
165
166//#include <TPad.h>
167#include <TStyle.h>
168#include <TCanvas.h>
169
170#include <TH2.h>
171#include <TH3.h>
172#include <TProfile.h>
173#include <TProfile2D.h>
174#include <TProfile3D.h>
175
176#include "MLog.h"
177#include "MLogManip.h"
178
179#include "MString.h"
180
181#include "MParList.h"
182#include "MBinning.h"
183#include "MDataPhrase.h"
184
185ClassImp(MH3);
186
187using namespace std;
188
189const TString MH3::gsDefName = "MH3";
190const TString MH3::gsDefTitle = "Container for a n-D Mars Histogram";
191
192// --------------------------------------------------------------------------
193//
194// Set fStyleBits to 0, Reset fBins to NULL, fScale to 1, set name and title
195// to gsDefName and gsDefTitle and if fHist!=NULL UseCurrentStyle and
196// SetDirectory(0)
197//
198void MH3::Init()
199{
200 fStyleBits = 0;
201
202 fWeight = NULL;
203
204 fData[0] = NULL;
205 fData[1] = NULL;
206 fData[2] = NULL;
207 fData[3] = NULL;
208
209 fBins[0] = NULL;
210 fBins[1] = NULL;
211 fBins[2] = NULL;
212
213 fScale[0] = 1;
214 fScale[1] = 1;
215 fScale[2] = 1;
216
217 fConversion = NULL;
218
219 fName = gsDefName;
220 fTitle = gsDefTitle;
221
222 if (!fHist)
223 return;
224
225 fHist->UseCurrentStyle();
226 fHist->SetDirectory(NULL);
227}
228
229// --------------------------------------------------------------------------
230//
231// Default constructor.
232//
233MH3::MH3(const Int_t dim, Type_t type) : fDimension(dim), fHist(NULL)
234{
235 // FIXME?
236 switch (type)
237 {
238 case kHistogram:
239 if (fDimension>3)
240 fDimension=3;
241 break;
242 case kProfile:
243 case kProfileSpread:
244 fDimension = -TMath::Abs(fDimension);
245 if (fDimension<-2)
246 fDimension = -2;
247 break;
248 }
249
250 switch (fDimension)
251 {
252 case 1:
253 fHist = new TH1D;
254 fHist->SetYTitle("Counts");
255 break;
256 case -1:
257 fHist = new TProfile;
258 fHist->SetYTitle("Average");
259 if (type==kProfileSpread)
260 static_cast<TProfile*>(fHist)->SetErrorOption("s");
261 break;
262 case 2:
263 fHist = new TH2D;
264 fHist->SetZTitle("Counts");
265 break;
266 case -2:
267 fHist = new TProfile2D;
268 fHist->SetZTitle("Average");
269 static_cast<TProfile2D*>(fHist)->BuildOptions(0, 0, type==kProfileSpread?"s":"");
270 break;
271 case 3:
272 fHist = new TH3D;
273 break;
274 case -3:
275 fHist = new TProfile2D;
276 fHist->SetZTitle("Average");
277 static_cast<TProfile2D*>(fHist)->BuildOptions(0, 0, type==kProfileSpread?"s":"");
278 break;
279 }
280
281 Init();
282}
283
284// --------------------------------------------------------------------------
285//
286// Creates an TH1D. memberx is filled into the X-bins. For a more detailed
287// description see the class description above.
288//
289MH3::MH3(const char *memberx, Type_t type) : fDimension(1)
290{
291 fHist = new TH1D;
292 fHist->SetYTitle("Counts");
293
294 Init();
295
296 fData[0] = new MDataPhrase(memberx);
297}
298
299// --------------------------------------------------------------------------
300//
301// Adapt a given histogram
302//
303MH3::MH3(const TH1 &h1) : fDimension(1)
304{
305 if (h1.InheritsFrom(TH3::Class()))
306 fDimension = 3;
307 if (h1.InheritsFrom(TH2::Class()))
308 fDimension = 2;
309
310 if (h1.InheritsFrom(TProfile3D::Class()))
311 {
312 fDimension = -3;
313 *fLog << warn << "WARNING - MH3::MH3(TH1&) does not support TProfile3D." << endl;
314 }
315 if (h1.InheritsFrom(TProfile2D::Class()))
316 fDimension = -2;
317 if (h1.InheritsFrom(TProfile::Class()))
318 fDimension = -1;
319
320 fHist = (TH1*)h1.Clone();
321
322 Init(); // Before without SeUseCurrentStyle!
323
324 switch (fDimension)
325 {
326 case 3:
327 case -2:
328 fData[2] = new MDataPhrase(h1.GetZaxis()->GetTitle());
329 case 2:
330 case -1:
331 fData[1] = new MDataPhrase(h1.GetYaxis()->GetTitle());
332 case 1:
333 fData[0] = new MDataPhrase(h1.GetXaxis()->GetTitle());
334 }
335}
336
337// --------------------------------------------------------------------------
338//
339// Creates an TH2D. memberx is filled into the X-bins. membery is filled
340// into the Y-bins. For a more detailed description see the class
341// description above.
342//
343MH3::MH3(const char *memberx, const char *membery, Type_t type)
344 : fDimension(type==kHistogram?2:-1)
345{
346
347 switch (fDimension)
348 {
349 case 2:
350 fHist = static_cast<TH1*>(new TH2D);
351 break;
352 case -1:
353 fHist = static_cast<TH1*>(new TProfile);
354 if (type==kProfileSpread)
355 static_cast<TProfile*>(fHist)->SetErrorOption("s");
356
357 break;
358 }
359
360 fHist->SetZTitle(fDimension>0?"Counts":"Average");
361
362 Init();
363
364 fData[0] = new MDataPhrase(memberx);
365 fData[1] = new MDataPhrase(membery);
366}
367
368// --------------------------------------------------------------------------
369//
370// Creates an TH3D. memberx is filled into the X-bins. membery is filled
371// into the Y-bins. membery is filled into the Z-bins. For a more detailed
372// description see the class description above.
373//
374MH3::MH3(const char *memberx, const char *membery, const char *memberz, Type_t type)
375 : fDimension(type==kHistogram?3:-2)
376{
377 switch (fDimension)
378 {
379 case 3:
380 fHist = static_cast<TH1*>(new TH3D);
381 break;
382 case -2:
383 fHist = static_cast<TH1*>(new TProfile2D);
384 static_cast<TProfile2D*>(fHist)->BuildOptions(0, 0, type==kProfileSpread?"s":"");
385 }
386
387 Init();
388
389 fData[0] = new MDataPhrase(memberx);
390 fData[1] = new MDataPhrase(membery);
391 fData[2] = new MDataPhrase(memberz);
392}
393
394// --------------------------------------------------------------------------
395//
396// Creates an TH3D. memberx is filled into the X-bins. membery is filled
397// into the Y-bins. membery is filled into the Z-bins. Weight is used as a
398// weight for the profile histogram. For a more detailed description see the
399// class description above.
400//
401MH3::MH3(const char *memberx, const char *membery, const char *memberz, const char *weight, Type_t type)
402 : fDimension(-3)
403{
404 fHist = static_cast<TH1*>(new TProfile3D);
405 static_cast<TProfile3D*>(fHist)->BuildOptions(0, 0, type==kProfileSpread?"s":"");
406
407 Init();
408
409 fData[0] = new MDataPhrase(memberx);
410 fData[1] = new MDataPhrase(membery);
411 fData[2] = new MDataPhrase(memberz);
412 fData[3] = new MDataPhrase(weight);
413}
414
415// --------------------------------------------------------------------------
416//
417// Deletes the histogram
418//
419MH3::~MH3()
420{
421 if (fHist)
422 delete fHist;
423
424 if (fConversion)
425 delete fConversion;
426
427 if (fWeight)
428 delete fWeight;
429
430 for (int i=0; i<4; i++)
431 if (fData[i])
432 delete fData[i];
433
434 for (int i=0; i<3; i++)
435 if (fLabels[i].GetDefault())
436 delete fLabels[i].GetDefault();
437}
438
439// --------------------------------------------------------------------------
440//
441// You can set a weight as a phrase additionally to the one given
442// as an argument to Fill (most likely from MFillH). The two weights
443// are multiplied together.
444//
445void MH3::SetWeight(const char *phrase)
446{
447 if (fWeight)
448 delete fWeight;
449 fWeight = new MDataPhrase(phrase);
450}
451
452// --------------------------------------------------------------------------
453//
454// Set a function which is applied to the histogram before it is displayed.
455// Note, that it only effects the displayed histogram.
456//
457// e.g. SetConversion("sqrt(x)");
458//
459Bool_t MH3::SetConversion(const char *func)
460{
461 if (TString(func).IsNull())
462 {
463 delete fConversion;
464 fConversion = 0;
465 return kTRUE;
466 }
467
468 fConversion = new TFormula;
469
470 // Must have a name otherwise all axis labels disappear like a miracle
471 fConversion->SetName("ConversionFunction");
472 if (fConversion->Compile(func))
473 {
474 *fLog << err << dbginf << "Syntax Error: TFormula::Compile failed for " << func << endl;
475 delete fConversion;
476 fConversion = 0;
477 return kFALSE;
478 }
479
480 gROOT->GetListOfFunctions()->Remove(fConversion);
481
482 return kTRUE;
483}
484
485// --------------------------------------------------------------------------
486//
487// The axis label is centered and the labeling of the axis is initialized.
488//
489// This function must not be called after any label has been created!
490//
491void MH3::InitLabels(TAxis &x) const
492{
493 x.CenterTitle();
494 x.SetBinLabel(1, "");
495 x.LabelsOption("h"); // FIXME: Is "a" thread safe? (Paint and Fill?)
496 x.GetLabels()->Delete();
497}
498
499// --------------------------------------------------------------------------
500//
501// Depending on the bits set the InitLabels(TAxis&) function for
502// the corresponding axes are called. In any case the kCanRebin bit
503// is set.
504//
505// This function must not be called after any label has been created!
506//
507void MH3::InitLabels(Labels_t type) const
508{
509 if (!fHist)
510 return;
511
512 if (type&kLabelsX && fHist->GetXaxis())
513 InitLabels(*fHist->GetXaxis());
514
515 if (type&kLabelsY && fHist->GetYaxis())
516 InitLabels(*fHist->GetYaxis());
517
518 if (type&kLabelsZ && fHist->GetZaxis())
519 InitLabels(*fHist->GetZaxis());
520
521 if (type&kLabelsXYZ)
522 fHist->SetBit(TH1::kCanRebin);
523}
524
525// --------------------------------------------------------------------------
526//
527// Return the corresponding Labels_t describing for which axis
528// axis-labels are switched on.
529//
530MH3::Labels_t MH3::GetLabels() const
531{
532 UInt_t type = kNoLabels;
533 if (fHist->GetXaxis() && fHist->GetXaxis()->GetLabels())
534 type |= kLabelsX;
535 if (fHist->GetYaxis() && fHist->GetYaxis()->GetLabels())
536 type |= kLabelsY;
537 if (fHist->GetZaxis() && fHist->GetZaxis()->GetLabels())
538 type |= kLabelsZ;
539 return (Labels_t)type;
540}
541
542// --------------------------------------------------------------------------
543//
544// Calls the LabelsDeflate from the histogram for all axes.
545// LabelsDeflate will just do nothing if the axis has no labels
546// initialized.
547//
548void MH3::DeflateLabels() const
549{
550 fHist->LabelsDeflate("X");
551 fHist->LabelsDeflate("Y");
552 fHist->LabelsDeflate("Z");
553}
554
555// --------------------------------------------------------------------------
556//
557// Returns the named label corresponding to the given value
558// and the given axis. The names are defined with the
559// DefineLabel-functions. if no name is defined the value
560// is converted to a string with %d and TMath::Nint.
561// If names are defined, but not for the given value, the default
562// label is returned instead. If no default is defined the
563// %d-converted string is returned.
564//
565const char *MH3::GetLabel(Int_t axe, Double_t val) const
566{
567 const Int_t v = TMath::Nint(val);
568
569 if (fLabels[axe].GetSize())
570 {
571 const char *l = fLabels[axe].GetObjName(v);
572 if (l)
573 return l;
574 }
575
576 return Form("%d", v);
577}
578
579// --------------------------------------------------------------------------
580//
581// Return the data members used by the data chain to be used in
582// MTask::AddBranchToList
583//
584TString MH3::GetDataMember() const
585{
586 TString str=fData[0]->GetDataMember();
587
588 for (int i=1; i<4; i++)
589 if (fData[i])
590 {
591 str += ";";
592 str += fData[i]->GetDataMember();
593 }
594
595 return str;
596}
597
598// --------------------------------------------------------------------------
599//
600// Setup the Binning for the histograms automatically if the correct
601// instances of MBinning are found in the parameter list
602// For a more detailed description see class description above.
603//
604Bool_t MH3::SetupFill(const MParList *plist)
605{
606 // reset histogram (necessary if the same eventloop is run more than once)
607 if (!TestBit(kDoNotReset))
608 fHist->Reset();
609
610 // Tokenize name into name and binnings names
611 TObjArray *tok = fName.Tokenize(";");
612
613 const TString name = (*tok)[0] ? (*tok)[0]->GetName() : fName.Data();
614
615 TString bx = (*tok)[1] ? (*tok)[1]->GetName() : Form("%sX", name.Data());
616 TString by = (*tok)[2] ? (*tok)[2]->GetName() : Form("%sY", name.Data());
617 TString bz = (*tok)[3] ? (*tok)[3]->GetName() : Form("%sZ", name.Data());
618
619 bx.Prepend("Binning");
620 by.Prepend("Binning");
621 bz.Prepend("Binning");
622
623 delete tok;
624
625 MBinning *binsx = NULL;
626 MBinning *binsy = NULL;
627 MBinning *binsz = NULL;
628
629 const Labels_t labels = GetLabels();
630
631 switch (TMath::Abs(fDimension))
632 {
633 case 3:
634 if (fData[2])
635 fHist->SetZTitle(fData[2]->GetTitle());
636 if (!(labels&kLabelsZ))
637 {
638 binsz = fBins[2] ? fBins[2] : (MBinning*)plist->FindObject(bz, "MBinning");
639 if (!binsz)
640 {
641 *fLog << err << dbginf << "MBinning '" << bz << "' not found... aborting." << endl;
642 return kFALSE;
643 }
644 if (binsz->HasTitle())
645 fHist->SetZTitle(binsz->GetTitle());
646 if (binsz->IsLogarithmic())
647 fHist->SetBit(kIsLogz);
648 }
649 case 2:
650 if (fData[1])
651 fHist->SetYTitle(fData[1]->GetTitle());
652 if (!(labels&kLabelsY))
653 {
654 binsy = fBins[1] ? fBins[1] : (MBinning*)plist->FindObject(by, "MBinning");
655 if (!binsy)
656 {
657 *fLog << err << dbginf << "MBinning '" << by << "' not found... aborting." << endl;
658 return kFALSE;
659 }
660 if (binsy->HasTitle())
661 fHist->SetYTitle(binsy->GetTitle());
662 if (binsy->IsLogarithmic())
663 fHist->SetBit(kIsLogy);
664 }
665 case 1:
666 if (fData[0]!=NULL)
667 fHist->SetXTitle(fData[0]->GetTitle());
668 if (!(labels&kLabelsX))
669 {
670 binsx = fBins[0] ? fBins[0] : (MBinning*)plist->FindObject(bx, "MBinning");
671 if (!binsx)
672 {
673 if (fDimension==1)
674 binsx = (MBinning*)plist->FindObject("Binning"+fName, "MBinning");
675
676 if (!binsx)
677 {
678 *fLog << err << dbginf << "Neither '" << bx << "' nor 'Binning" << fName << "' found... aborting." << endl;
679 return kFALSE;
680 }
681 }
682 if (binsx->HasTitle())
683 fHist->SetXTitle(binsx->GetTitle());
684 if (binsx->IsLogarithmic())
685 fHist->SetBit(kIsLogx);
686 }
687 }
688
689 // PreProcess existing fData members
690 for (int i=0; i<4; i++)
691 if (fData[i] && !fData[i]->PreProcess(plist))
692 return kFALSE;
693
694 if (fWeight && !fWeight->PreProcess(plist))
695 return kFALSE;
696
697 TString title(fDimension>0?"Histogram":"Profile");
698 title += " for ";
699 title += name;
700 title += Form(" (%dD)", TMath::Abs(fDimension));
701
702 fHist->SetName(name);
703 fHist->SetTitle(fTitle==gsDefTitle ? title : fTitle);
704 fHist->SetDirectory(0);
705
706 // This is for the case we have set lables
707 const MBinning def(1, 0, 1);
708 if (!binsx)
709 binsx = const_cast<MBinning*>(&def);
710 if (!binsy)
711 binsy = const_cast<MBinning*>(&def);
712 if (!binsz)
713 binsz = const_cast<MBinning*>(&def);
714
715 // set binning
716 switch (TMath::Abs(fDimension))
717 {
718 case 1:
719 SetBinning(*fHist, *binsx);
720 return kTRUE;
721 case 2:
722 SetBinning(static_cast<TH2&>(*fHist), *binsx, *binsy);
723 return kTRUE;
724 case 3:
725 SetBinning(static_cast<TH3&>(*fHist), *binsx, *binsy, *binsz);
726 return kTRUE;
727 }
728
729 *fLog << err << "ERROR - MH3 has " << TMath::Abs(fDimension) << " dimensions!" << endl;
730 return kFALSE;
731}
732
733// --------------------------------------------------------------------------
734//
735// Set the name of the histogram ant the MH3 container
736//
737void MH3::SetName(const char *name)
738{
739 if (fHist)
740 {
741 if (gPad)
742 {
743 const TString pfx(MString::Format("%sProfX", fHist->GetName()));
744 const TString pfy(MString::Format("%sProfY", fHist->GetName()));
745
746 TProfile *p = 0;
747 if ((p=dynamic_cast<TProfile*>(gPad->FindObject(pfx))))
748 p->SetName(MString::Format("%sProfX", name));
749 if ((p=dynamic_cast<TProfile*>(gPad->FindObject(pfy))))
750 p->SetName(MString::Format("%sProfY", name));
751 }
752
753 fHist->SetName(name);
754 fHist->SetDirectory(0);
755
756 }
757 MParContainer::SetName(name);
758}
759
760// --------------------------------------------------------------------------
761//
762// Set the title of the histogram ant the MH3 container
763//
764void MH3::SetTitle(const char *title)
765{
766 if (fHist)
767 fHist->SetTitle(title);
768 MParContainer::SetTitle(title);
769}
770
771// --------------------------------------------------------------------------
772//
773// Fills the one, two or three data members into our histogram
774//
775Int_t MH3::Fill(const MParContainer *par, const Stat_t ww)
776{
777 // Get Information about labels (UInt_t, to supress warning about
778 // unhandeled cases in switch)
779 const UInt_t type = GetLabels();
780
781 // Get values for axis
782 Double_t x=0;
783 Double_t y=0;
784 Double_t z=0;
785 Double_t t=0;
786 Double_t w=ww;
787
788 switch (fDimension)
789 {
790 case -3:
791 t = fData[3]->GetValue();
792 case -2:
793 case 3:
794 z = fData[2]->GetValue()*fScale[2];
795 case -1:
796 case 2:
797 y = fData[1]->GetValue()*fScale[1];
798 case 1:
799 x = fData[0]->GetValue()*fScale[0];
800 }
801
802 if (fData[4])
803 w *= fData[4]->GetValue();
804
805 // If label option is set, convert value to label
806 TString labelx, labely, labelz;
807 if (type&kLabelsX)
808 labelx = GetLabel(0, x);
809 if (type&kLabelsY)
810 labely = GetLabel(1, y);
811 if (type&kLabelsZ)
812 labelz = GetLabel(2, z);
813
814 // Fill histogram
815 switch (fDimension)
816 {
817 case 3:
818 switch (type)
819 {
820 case kNoLabels:
821 static_cast<TH3*>(fHist)->Fill(x, y, z, w);
822 return kTRUE;
823 case kLabelsX:
824 static_cast<TH3*>(fHist)->Fill(labelx, y, z);
825 return kTRUE;
826 case kLabelsY:
827 static_cast<TH3*>(fHist)->Fill(x, labely, z, w);
828 return kTRUE;
829 case kLabelsZ:
830 static_cast<TH3*>(fHist)->Fill(x, y, labelz, w);
831 return kTRUE;
832 case kLabelsXY:
833 static_cast<TH3*>(fHist)->Fill(labelx, labely, z, w);
834 return kTRUE;
835 case kLabelsXZ:
836 static_cast<TH3*>(fHist)->Fill(labelx, y, labelz, w);
837 return kTRUE;
838 case kLabelsYZ:
839 static_cast<TH3*>(fHist)->Fill(x, labely, labelz, w);
840 return kTRUE;
841 case kLabelsXYZ:
842 static_cast<TH3*>(fHist)->Fill(labelx, labely, labelz, w);
843 return kTRUE;
844 }
845 break;
846 case 2:
847 switch (type)
848 {
849 case kNoLabels:
850 static_cast<TH2*>(fHist)->Fill(x, y, w);
851 return kTRUE;
852 case kLabelsX:
853 static_cast<TH2*>(fHist)->Fill(x, labely, w);
854 return kTRUE;
855 case kLabelsY:
856 static_cast<TH2*>(fHist)->Fill(labelx, y, w);
857 return kTRUE;
858 case kLabelsXY:
859 static_cast<TH2*>(fHist)->Fill(labelx, labely, w);
860 return kTRUE;
861 }
862 break;
863 case 1:
864 switch (type)
865 {
866 case kNoLabels:
867 fHist->Fill(x, w);
868 return kTRUE;
869 case kLabelsX:
870 fHist->Fill(labelx, w);
871 return kTRUE;
872 }
873 break;
874 case -1:
875 switch (type)
876 {
877 case kNoLabels:
878 static_cast<TProfile*>(fHist)->Fill(x, y, w);
879 return kTRUE;
880 case kLabelsX:
881 static_cast<TProfile*>(fHist)->Fill(labelx, y, w);
882 return kTRUE;
883 }
884 break;
885 case -2:
886 switch (type)
887 {
888 case kNoLabels:
889 static_cast<TProfile2D*>(fHist)->Fill(x, y, z, w);
890 return kTRUE;
891 case kLabelsX:
892 static_cast<TProfile2D*>(fHist)->Fill(labelx, y, z);
893 return kTRUE;
894 case kLabelsY:
895 static_cast<TProfile2D*>(fHist)->Fill(x, labely, z);
896 return kTRUE;
897 case kLabelsXY:
898 static_cast<TProfile2D*>(fHist)->Fill(labelx, labely, z);
899 return kTRUE;
900 }
901 break;
902 case -3:
903 switch (type)
904 {
905 case kNoLabels:
906 static_cast<TProfile3D*>(fHist)->Fill(x, y, z, t, w);
907 return kTRUE;
908 default:
909 *fLog << err << "ERROR - Labels not supported in TProfile3D." << endl;
910 }
911 break;
912 }
913
914 *fLog << err << "MH3::Fill: ERROR - A fatal error occured." << endl;
915 return kERROR;
916}
917
918// --------------------------------------------------------------------------
919//
920// If an auto range bit is set the histogram range of the corresponding
921// axis is set to show only the non-empty bins (with a single empty bin
922// on both sides)
923//
924Bool_t MH3::Finalize()
925{
926 DeflateLabels();
927
928 Bool_t autorangex=TESTBIT(fStyleBits, 0);
929 Bool_t autorangey=TESTBIT(fStyleBits, 1);
930 //Bool_t autorangez=TESTBIT(fStyleBits, 2);
931
932 Int_t lo, hi;
933 if (autorangex)
934 {
935 GetRangeX(*fHist, lo, hi);
936 fHist->GetXaxis()->SetRange(lo-2, hi+1);
937 }
938 if (autorangey)
939 {
940 GetRangeY(*fHist, lo, hi);
941 fHist->GetYaxis()->SetRange(lo-2, hi+1);
942 }
943 /*
944 if (autorangez)
945 {
946 GetRangeZ(*fHist, lo, hi);
947 fHist->GetZaxis()->SetRange(lo-2, hi+1);
948 }
949 */
950 return kTRUE;
951}
952
953// --------------------------------------------------------------------------
954//
955// Apply the conversion function to the contents stored in fHist and
956// store the result in h.
957//
958// In the case of a TProfile we keep it a TProfile to keep the mean and
959// rms in y displayed. To get the error correctly we have to reverse
960// the calculation done in TProfile::GetBinError of course.
961//
962void MH3::Convert(TH1 &h) const
963{
964 const Bool_t prof = h.InheritsFrom(TProfile::Class()) || h.InheritsFrom(TProfile2D::Class());
965
966 for (Int_t z=0; z<=h.GetNbinsZ()+1; z++)
967 for (Int_t y=0; y<=h.GetNbinsY()+1; y++)
968 for (Int_t x=0; x<=h.GetNbinsX()+1; x++)
969 {
970 h.SetBinContent(x, y, z, fConversion->Eval(fHist->GetBinContent(x, y, z)));
971
972 if (prof)
973 h.SetBinError(x, y, z, TMath::Hypot(fConversion->Eval(fHist->GetBinContent(x, y, z)),
974 fConversion->Eval(fHist->GetBinError( x, y, z))));
975 else
976 h.SetBinError(x, y, z, fConversion->Eval(fHist->GetBinError(x, y, z)));
977 }
978
979 TProfile *p1 = dynamic_cast<TProfile*>(fHist);
980 if (p1)
981 for (Int_t i=0; i<p1->GetSize(); i++)
982 static_cast<TProfile&>(h).SetBinEntries(i, p1->GetBinEntries(i)>0 ? 1 : 0);
983
984 TProfile *p2 = dynamic_cast<TProfile*>(fHist);
985 if (p2)
986 for (Int_t i=0; i<p2->GetSize(); i++)
987 static_cast<TProfile2D&>(h).SetBinEntries(i, p2->GetBinEntries(i)>0 ? 1 : 0);
988}
989
990// --------------------------------------------------------------------------
991//
992// FIXME
993//
994void MH3::Paint(Option_t *o)
995{
996 TProfile *p=0;
997
998 if (TMath::Abs(fDimension)==2)
999 MH::SetPalette("pretty");
1000
1001 if (fConversion)
1002 {
1003 TH1 *h = 0;
1004 if ((h=dynamic_cast<TH1*>(gPad->FindObject(fHist->GetName()))))
1005 Convert(*h);
1006 }
1007
1008 const TString pfx(MString::Format("%sProfX", fHist->GetName()));
1009 if ((p=dynamic_cast<TProfile*>(gPad->FindObject(pfx))))
1010 {
1011 Int_t col = p->GetLineColor();
1012 p = ((TH2*)fHist)->ProfileX(pfx, -1, -1, "s");
1013 p->SetLineColor(col);
1014 }
1015
1016 const TString pfy(MString::Format("%sProfY", fHist->GetName()));
1017 if ((p=dynamic_cast<TProfile*>(gPad->FindObject(pfy))))
1018 {
1019 Int_t col = p->GetLineColor();
1020 p = ((TH2*)fHist)->ProfileY(pfy, -1, -1, "s");
1021 p->SetLineColor(col);
1022 }
1023/*
1024 if (fHist->TestBit(kIsLogx) && fHist->GetEntries()>0) gPad->SetLogx();
1025 if (fHist->TestBit(kIsLogy) && fHist->GetEntries()>0) gPad->SetLogy();
1026 if (fHist->TestBit(kIsLogz) && fHist->GetEntries()>0) gPad->SetLogz();
1027 */
1028}
1029
1030// --------------------------------------------------------------------------
1031//
1032// If Xmax is < 3000*Xmin SetMoreLogLabels is called. If Xmax<5000
1033// the exponent is switched off (SetNoExponent)
1034//
1035void MH3::HandleLogAxis(TAxis &axe) const
1036{
1037 if (axe.GetXmax()>3000*axe.GetXmin())
1038 return;
1039
1040 axe.SetMoreLogLabels();
1041 if (axe.GetXmax()<5000)
1042 axe.SetNoExponent();
1043}
1044
1045// --------------------------------------------------------------------------
1046//
1047// Creates a new canvas and draws the histogram into it.
1048//
1049// Possible options are:
1050// PROFX: Draw a x-profile into the histogram (for 2D histograms only)
1051// PROFY: Draw a y-profile into the histogram (for 2D histograms only)
1052// ONLY: Draw the profile histogram only (for 2D histograms only)
1053// BLUE: Draw the profile in blue color instead of the histograms
1054// line color
1055//
1056// If the kIsLog?-Bit is set the axis is displayed lkogarithmically.
1057// eg this is set when applying a logarithmic MBinning
1058//
1059// Be careful: The histogram belongs to this object and won't get deleted
1060// together with the canvas.
1061//
1062void MH3::Draw(Option_t *opt)
1063{
1064 TVirtualPad *pad = gPad ? gPad : MakeDefCanvas(fHist);
1065 pad->SetBorderMode(0);
1066 pad->SetGridx();
1067 pad->SetGridy();
1068
1069 if (fHist->TestBit(kIsLogx))
1070 {
1071 pad->SetLogx();
1072 HandleLogAxis(*fHist->GetXaxis());
1073 }
1074 if (fHist->TestBit(kIsLogy))
1075 {
1076 pad->SetLogy();
1077 HandleLogAxis(*fHist->GetYaxis());
1078 }
1079 if (fHist->TestBit(kIsLogz))
1080 {
1081 pad->SetLogz();
1082 HandleLogAxis(*fHist->GetZaxis());
1083 }
1084
1085 fHist->SetFillStyle(4000);
1086
1087 TString str(opt);
1088 str.ToLower();
1089
1090 if (GetLabels())
1091 {
1092 if (str.IsNull() && fDimension==2)
1093 str = "colz";
1094
1095 if (str.Contains("box", TString::kIgnoreCase) && fDimension==2)
1096 fHist->SetLineColor(kBlue);
1097 }
1098
1099 const Bool_t only = str.Contains("only") && TMath::Abs(fDimension)==2;
1100 const Bool_t same = str.Contains("same") && TMath::Abs(fDimension)<3;
1101 const Bool_t blue = str.Contains("blue") && TMath::Abs(fDimension)==2;
1102 const Bool_t profx = str.Contains("profx") && TMath::Abs(fDimension)==2;
1103 const Bool_t profy = str.Contains("profy") && TMath::Abs(fDimension)==2;
1104
1105 str.ReplaceAll("only", "");
1106 str.ReplaceAll("blue", "");
1107 str.ReplaceAll("profx", "");
1108 str.ReplaceAll("profy", "");
1109 str.ReplaceAll(" ", "");
1110
1111 if (same && TMath::Abs(fDimension)==1)
1112 {
1113 fHist->SetLineColor(kBlue);
1114 fHist->SetMarkerColor(kBlue);
1115 }
1116
1117
1118 TH1 *h = fHist;
1119
1120 if (fConversion)
1121 {
1122 h = static_cast<TH1*>(fHist->Clone());
1123
1124 h->SetDirectory(0);
1125 h->SetBit(kCanDelete);
1126
1127 Convert(*h);
1128 }
1129
1130 // FIXME: We may have to remove all our own options from str!
1131 if (!only)
1132 h->Draw(str);
1133
1134 AppendPad();
1135
1136 TProfile *p=0;
1137 if (profx)
1138 {
1139 const TString pfx(MString::Format("%sProfX", h->GetName()));
1140
1141 if (same && (p=dynamic_cast<TProfile*>(gPad->FindObject(pfx))))
1142 *fLog << warn << "TProfile " << pfx << " already in pad." << endl;
1143
1144 p = ((TH2*)h)->ProfileX(pfx, -1, -1, "s");
1145 p->UseCurrentStyle();
1146 p->SetLineColor(blue ? kBlue : h->GetLineColor());
1147 p->SetBit(kCanDelete);
1148 p->SetDirectory(NULL);
1149 p->SetXTitle(h->GetXaxis()->GetTitle());
1150 p->SetYTitle(h->GetYaxis()->GetTitle());
1151 p->Draw(only&&!same?"":"same");
1152 }
1153 if (profy)
1154 {
1155 const TString pfy(MString::Format("%sProfY", h->GetName()));
1156
1157 if (same && (p=dynamic_cast<TProfile*>(gPad->FindObject(pfy))))
1158 *fLog << warn << "TProfile " << pfy << " already in pad." << endl;
1159
1160 p = ((TH2*)h)->ProfileY(pfy, -1, -1, "s");
1161 p->UseCurrentStyle();
1162 p->SetLineColor(blue ? kBlue : h->GetLineColor());
1163 p->SetBit(kCanDelete);
1164 p->SetDirectory(NULL);
1165 p->SetYTitle(h->GetXaxis()->GetTitle());
1166 p->SetXTitle(h->GetYaxis()->GetTitle());
1167 p->Draw(only&&!same?"":"same");
1168 }
1169
1170 //AppendPad("log");
1171}
1172
1173// --------------------------------------------------------------------------
1174//
1175// Implementation of SavePrimitive. Used to write the call to a constructor
1176// to a macro. In the original root implementation it is used to write
1177// gui elements to a macro-file.
1178//
1179void MH3::StreamPrimitive(ostream &out) const
1180{
1181 TString name = GetUniqueName();
1182
1183 out << " MH3 " << name << "(\"";
1184 out << fData[0]->GetRule() << "\"";
1185 if (fDimension>1 || fDimension<0)
1186 out << ", \"" << fData[1]->GetRule() << "\"";
1187 if (fDimension>2 || fDimension<-1)
1188 out << ", \"" << fData[2]->GetRule() << "\"";
1189
1190 out << ");" << endl;
1191
1192 if (fName!=gsDefName)
1193 out << " " << name << ".SetName(\"" << fName << "\");" << endl;
1194
1195 if (fTitle!=gsDefTitle)
1196 out << " " << name << ".SetTitle(\"" << fTitle << "\");" << endl;
1197
1198 if (fWeight)
1199 out << " " << name << ".SetWeight(\"" << fWeight->GetRule() << "\");" << endl;
1200
1201 switch (fDimension)
1202 {
1203 case -3:
1204 if (fScale[3]!=1)
1205 out << " " << name << ".SetScaleT(" << fScale[3] << ");" << endl;
1206 case -2:
1207 case 3:
1208 if (fScale[2]!=1)
1209 out << " " << name << ".SetScaleZ(" << fScale[2] << ");" << endl;
1210 case -1:
1211 case 2:
1212 if (fScale[1]!=1)
1213 out << " " << name << ".SetScaleY(" << fScale[1] << ");" << endl;
1214 case 1:
1215 if (fScale[0]!=1)
1216 out << " " << name << ".SetScaleX(" << fScale[0] << ");" << endl;
1217 }
1218}
1219
1220MH3::Type_t MH3::GetType() const
1221{
1222 switch (fDimension)
1223 {
1224 case -1:
1225 return TString(static_cast<TProfile*>(fHist)->GetErrorOption())=="s" ? kProfileSpread : kProfile;
1226 case -2:
1227 return TString(static_cast<TProfile2D*>(fHist)->GetErrorOption())=="s" ? kProfileSpread : kProfile;
1228 case -3:
1229 return TString(static_cast<TProfile3D*>(fHist)->GetErrorOption())=="s" ? kProfileSpread : kProfile;
1230 }
1231 return kHistogram;
1232}
1233
1234// --------------------------------------------------------------------------
1235//
1236// Used to rebuild a MH3 object of the same type (data members,
1237// dimension, ...)
1238//
1239MParContainer *MH3::New() const
1240{
1241 // FIXME: TREAT THE NEW OPTIONS CORRECTLY (PROFILE, LABELS)
1242
1243 MH3 *h = NULL;
1244
1245 if (fData[0] == NULL)
1246 h=new MH3(fDimension);
1247 else
1248 switch (fDimension)
1249 {
1250 case 1:
1251 h=new MH3(fData[0]->GetRule());
1252 break;
1253 case 2:
1254 case -1:
1255 h=new MH3(fData[0]->GetRule(), fData[1]->GetRule(), GetType());
1256 break;
1257 case 3:
1258 case -2:
1259 h=new MH3(fData[0]->GetRule(), fData[1]->GetRule(), fData[2]->GetRule(), GetType());
1260 break;
1261 case -3:
1262 h=new MH3(fData[0]->GetRule(), fData[1]->GetRule(), fData[2]->GetRule(), fData[3]->GetRule());
1263 break;
1264 }
1265
1266 switch (fDimension)
1267 {
1268 case -2:
1269 case 3:
1270 h->SetScaleZ(fScale[2]);
1271 case 2:
1272 case -1:
1273 h->SetScaleY(fScale[1]);
1274 case 1:
1275 h->SetScaleX(fScale[0]);
1276 }
1277
1278 if (fData[4])
1279 h->SetWeight(fData[4]->GetRule());
1280
1281 return h;
1282}
1283
1284// --------------------------------------------------------------------------
1285//
1286// FIXME
1287//
1288TString MH3::GetRule(const Char_t axis) const
1289{
1290 switch (tolower(axis))
1291 {
1292 case 'x':
1293 case 'X':
1294 return fData[0] ? fData[0]->GetRule() : TString("");
1295 case 'y':
1296 case 'Y':
1297 return fData[1] ? fData[1]->GetRule() : TString("");
1298 case 'z':
1299 case 'Z':
1300 return fData[2] ? fData[2]->GetRule() : TString("");
1301 case 't':
1302 case 'T':
1303 return fData[3] ? fData[3]->GetRule() : TString("");
1304 case 'w':
1305 case 'W':
1306 return fWeight ? fWeight->GetRule() : TString("");
1307 default:
1308 return "<n/a>";
1309 }
1310}
1311
1312// --------------------------------------------------------------------------
1313//
1314// Returns the total number of bins in a histogram (excluding under- and
1315// overflow bins)
1316//
1317Int_t MH3::GetNbins() const
1318{
1319 Int_t num = 1;
1320
1321 switch (TMath::Abs(fDimension))
1322 {
1323 case 3:
1324 num *= fHist->GetNbinsZ()+2;
1325 case 2:
1326 num *= fHist->GetNbinsY()+2;
1327 case 1:
1328 num *= fHist->GetNbinsX()+2;
1329 }
1330
1331 return num;
1332}
1333
1334// --------------------------------------------------------------------------
1335//
1336// Returns the total number of bins in a histogram (excluding under- and
1337// overflow bins) Return -1 if bin is underflow or overflow
1338//
1339Int_t MH3::FindFixBin(Double_t x, Double_t y, Double_t z) const
1340{
1341 const TAxis &axex = *fHist->GetXaxis();
1342 const TAxis &axey = *fHist->GetYaxis();
1343 const TAxis &axez = *fHist->GetZaxis();
1344
1345 Int_t binz = 0;
1346 Int_t biny = 0;
1347 Int_t binx = 0;
1348
1349 switch (fDimension)
1350 {
1351 case 3:
1352 case -2:
1353 binz = axez.FindFixBin(z);
1354 if (binz>axez.GetLast() || binz<axez.GetFirst())
1355 return -1;
1356 case 2:
1357 case -1:
1358 biny = axey.FindFixBin(y);
1359 if (biny>axey.GetLast() || biny<axey.GetFirst())
1360 return -1;
1361 case 1:
1362 binx = axex.FindFixBin(x);
1363 if (binx<axex.GetFirst() || binx>axex.GetLast())
1364 return -1;
1365 }
1366
1367 const Int_t nx = fHist->GetNbinsX()+2;
1368 const Int_t ny = fHist->GetNbinsY()+2;
1369
1370 return binx + nx*(biny +ny*binz);
1371}
1372
1373// --------------------------------------------------------------------------
1374//
1375// Return the MObjLookup corresponding to the axis/character.
1376// Note that only lower-case charecters (x, y, z) are supported.
1377// If for the axis no labels were set, the corresponding
1378// InitLabels is called.
1379//
1380MObjLookup *MH3::GetLabels(char axe)
1381{
1382 if (!fHist)
1383 return 0;
1384
1385 TAxis *x = 0;
1386
1387 switch (axe)
1388 {
1389 case 'x':
1390 x = fHist->GetXaxis();
1391 break;
1392 case 'y':
1393 x = fHist->GetYaxis();
1394 break;
1395 case 'z':
1396 x = fHist->GetZaxis();
1397 break;
1398 }
1399
1400 if (!x)
1401 return 0;
1402
1403 const Int_t idx = axe-'x';
1404
1405 if (!x->GetLabels())
1406 switch (idx)
1407 {
1408 case 0:
1409 InitLabels(kLabelsX);
1410 break;
1411 case 1:
1412 InitLabels(kLabelsY);
1413 break;
1414 case 2:
1415 InitLabels(kLabelsZ);
1416 break;
1417 }
1418
1419 return &fLabels[idx];
1420}
1421
1422// --------------------------------------------------------------------------
1423//
1424// Set a default label which is used if no other is found in the list
1425// of labels. if a default was set already it is overwritten. If the
1426// axis has not yet been initialized to use labels it it now.
1427//
1428void MH3::DefaultLabel(char axe, const char *name)
1429{
1430 MObjLookup *arr = GetLabels(axe);
1431 if (!arr)
1432 return;
1433
1434 if (arr->GetDefault())
1435 {
1436 delete arr->GetDefault();
1437 arr->SetDefault(0);
1438 }
1439
1440 if (name)
1441 arr->SetDefault(new TObjString(name));
1442}
1443
1444// --------------------------------------------------------------------------
1445//
1446// Define a name for a label. More than one label can have the same
1447// name. If the axis has not yet been initialized to use labels
1448// it it now.
1449//
1450void MH3::DefineLabel(char axe, Int_t label, const char *name)
1451{
1452 MObjLookup *arr = GetLabels(axe);
1453
1454 if (!arr || !name)
1455 return;
1456
1457 if (arr->GetObj(label)!=arr->GetDefault())
1458 return;
1459
1460 arr->Add(label, name);
1461}
1462
1463// --------------------------------------------------------------------------
1464//
1465// Define names for labels, like
1466// 1=Trig;2=Cal;4=Ped;8=Lvl2
1467// More than one label can have the same name. If the axis has not
1468// yet been initialized to use labels it it now.
1469//
1470// A default cannot be set here. Use DefaultLabel instead.
1471//
1472void MH3::DefineLabels(char axe, const TString &labels)
1473{
1474 TObjArray *arr = labels.Tokenize(';');
1475
1476 for (int i=0; i<arr->GetEntries(); i++)
1477 {
1478 const char *s = (*arr)[0]->GetName();
1479 const char *v = strchr(s, '=');
1480
1481 if (v)
1482 DefineLabel(axe, atoi(s), v+1);
1483 }
1484
1485 delete arr;
1486}
1487
1488void MH3::RecursiveRemove(TObject *obj)
1489{
1490 if (obj==fHist)
1491 fHist = 0;
1492
1493 MH::RecursiveRemove(obj);
1494}
Note: See TracBrowser for help on using the repository browser.