source: branches/Mars_McMismatchStudy/mhist/MHCompProb.cc@ 19062

Last change on this file since 19062 was 8076, checked in by tbretz, 18 years ago
*** empty log message ***
File size: 5.9 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): Abelardo Moralejo <mailto:moralejo@pd.infn.it>
19! Author(s): Thomas Bretz, 5/2002 <mailto:moralejo@pd.infn.it>
20!
21! Copyright: MAGIC Software Development, 2000-2002
22!
23!
24\* ======================================================================== */
25
26///////////////////////////////////////////////////////////////////////
27//
28// MHCompProb
29//
30// This class contains different histograms of the Hillas parameters
31// and composite probabilities based on them.
32//
33///////////////////////////////////////////////////////////////////////
34
35#include "MHCompProb.h"
36
37#include <TH2.h>
38#include <TPad.h>
39#include <TText.h>
40#include <TStyle.h>
41#include <TCanvas.h>
42#include <TProfile.h>
43
44#include "MLog.h"
45#include "MLogManip.h"
46
47#include "MParList.h"
48#include "MBinning.h"
49#include "MDataPhrase.h"
50
51#include "MMcEvt.hxx"
52
53ClassImp(MHCompProb);
54
55using namespace std;
56
57// --------------------------------------------------------------------------
58//
59// Setup histograms
60//
61MHCompProb::MHCompProb(Int_t nbins, const char *name, const char *title)
62 : fNumLoop(0)
63{
64 //
65 // set the name and title of this object
66 //
67 fName = name ? name : "MHCompProb";
68 fTitle = title ? title : "Gamma/Hadron Separation Quality Histograms";
69
70 fData = new TList;
71 fRules = new TList;
72 fHists = new TList;
73 fHistVar = new TList;
74
75 fData->SetOwner();
76 fRules->SetOwner();
77 fHists->SetOwner();
78 fHistVar->SetOwner();
79}
80
81// --------------------------------------------------------------------------
82//
83// Delete the histograms
84//
85MHCompProb::~MHCompProb()
86{
87 delete fData;
88 delete fRules;
89 delete fHists;
90 delete fHistVar;
91}
92
93// --------------------------------------------------------------------------
94//
95//
96//
97void MHCompProb::Add(const char *rule, Int_t n, Float_t min, Float_t max)
98{
99 MDataPhrase &chain = *new MDataPhrase(rule);
100 fData->Add(&chain);
101
102 TNamed &name = *new TNamed(rule, "");
103 fRules->Add(&name);
104
105 TH1D &hist = *new TH1D(Form("Hist_%s", rule), rule, n, min, max);
106 hist.SetXTitle(rule);
107 hist.SetYTitle("Counts");
108 hist.SetDirectory(NULL);
109 fHists->Add(&hist);
110
111 TH1D &varhist = *new TH1D;
112 varhist.SetName(Form("Var_%s", rule));
113 varhist.SetTitle(rule);
114 varhist.SetXTitle(rule);
115 varhist.SetYTitle("Counts");
116 varhist.SetDirectory(NULL);
117 fHistVar->Add(&varhist);
118}
119
120// --------------------------------------------------------------------------
121//
122//
123//
124Bool_t MHCompProb::SetupFill(const MParList *plist)
125{
126 if (fData->GetSize()==0)
127 {
128 *fLog << err << "No data members spcified for usage... aborting." << endl;
129 return kFALSE;
130 }
131
132 TIter Next(fData);
133 MData *data=NULL;
134 while ((data=(MData*)Next()))
135 if (!data->PreProcess(plist))
136 return kFALSE;
137
138 return kTRUE;
139}
140
141// --------------------------------------------------------------------------
142//
143//
144//
145void MHCompProb::Fill(TList &list)
146{
147 MData *data = NULL;
148
149 TIter NextD(fData);
150 TIter NextH(&list);
151
152 while ((data=(MData*)NextD()))
153 {
154 TH1D *hist = (TH1D*)NextH();
155 hist->Fill(data->GetValue());
156 }
157}
158
159// --------------------------------------------------------------------------
160//
161//
162//
163Bool_t MHCompProb::Fill(const MParContainer *par, const Stat_t w)
164{
165 const MMcEvt &mcevt = *(MMcEvt*)par;
166
167 switch (fNumLoop)
168 {
169 case 0: // First loop : fill the fixed-bin histograms with gammas.
170 if (mcevt.GetPartId() == MMcEvt::kGAMMA)
171 Fill(*fHists);
172 return kTRUE;
173
174 case 1: // Second Loop: fill the variable-bin histograms with protons.
175 if (mcevt.GetPartId() != MMcEvt::kGAMMA)
176 Fill(*fHistVar);
177 return kTRUE;
178 default:
179 *fLog << err << "Error - Invalid Loop Number... aborting." << endl;
180 return kFALSE;
181 }
182}
183
184// --------------------------------------------------------------------------
185//
186//
187//
188Bool_t MHCompProb::Finalize()
189{
190 switch (fNumLoop++)
191 {
192 case 0:
193 *fLog << inf << "Finished filling fixed bin size histograms with gamma-data." << endl;
194 SetBinningHistVar();
195 fHists->Delete();
196 return kTRUE;
197 case 1:
198 *fLog << inf << "Finished filling variable bin size histogram with proton data." << endl;
199 return kTRUE;
200 default:
201 *fLog << err << "Error - Invalid Loop Number... aborting." << endl;
202 return kFALSE;
203 }
204}
205
206// --------------------------------------------------------------------------
207//
208//
209//
210void MHCompProb::SetBinningHistVar()
211{
212 Int_t nedges = 51; // Number of bins in variable-bin histograms.
213
214 TIter NextH(fHists);
215 TIter NextV(fHistVar);
216 TH1D *hist = NULL;
217 while ((hist=(TH1D*)NextH()))
218 {
219 Int_t n = hist->GetNbinsX();
220
221 TArrayD edges(nedges);
222
223 edges[0] = hist->GetBinLowEdge(1);
224 edges[nedges-1] = hist->GetBinLowEdge(n+1);
225
226 Float_t newwidth = hist->Integral(1, n)/nedges;
227
228 Int_t jbin = 1;
229 for (Int_t j=1; j<n && jbin<nedges-1; j++)
230 {
231 if (hist->Integral(1, j) <= jbin*newwidth)
232 continue;
233
234 edges[jbin++] = hist->GetBinLowEdge(j+1);
235 }
236
237 MBinning bins;
238 bins.SetEdges(edges);
239
240 SetBinning((TH1D*)NextV(), &bins);
241 }
242}
243
Note: See TracBrowser for help on using the repository browser.