]> git.uio.no Git - u/mrichter/AliRoot.git/blame - STEER/AliMCQA.cxx
Corrections to destructors
[u/mrichter/AliRoot.git] / STEER / AliMCQA.cxx
CommitLineData
65fb704d 1/**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3 * *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
6 * *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
15
16/*
17$Log$
e460afec 18Revision 1.3 2000/12/18 14:16:31 alibrary
19HP compatibility fix
20
af2e514f 21Revision 1.2 2000/12/18 11:33:48 alibrary
22New call frequence histograms per module and volume
23
8409d8c9 24Revision 1.1 2000/11/30 07:12:48 alibrary
25Introducing new Rndm and QA classes
26
65fb704d 27*/
28
29///////////////////////////////////////////////////////////////////////////////
30// //
31// //
32///////////////////////////////////////////////////////////////////////////////
33
8409d8c9 34
65fb704d 35#include <strings.h>
36
37#include "TObjArray.h"
38#include "TH1.h"
39#include "TList.h"
40#include "TROOT.h"
41#include "TBrowser.h"
42#include "TMath.h"
43#include "TLorentzVector.h"
44#include "TDatabasePDG.h"
45#include "TMath.h"
8409d8c9 46#include "TPad.h"
47#include "TExec.h"
48#include "TPaveLabel.h"
49#include "TCanvas.h"
65fb704d 50
51#include "AliMCQA.h"
52#include "AliRun.h"
53#include "AliModule.h"
54#include "AliMC.h"
55
56ClassImp(AliMCQA)
57
58
59//_____________________________________________________________________________
8409d8c9 60AliMCQA::AliMCQA() : fQAList(0), fDetDone(0), fQAHist(0), fVolNames(0),
61 fModNames(0),fMPaveLabel(0),fVPaveLabel(0)
62{
63 //
64 // Default constructor
65 //
66}
67
68//_____________________________________________________________________________
69AliMCQA::AliMCQA(Int_t ndets) : fMPaveLabel(0),fVPaveLabel(0)
65fb704d 70{
71 //
72 // Constructor, creates the list of lists of histograms
73 //
74 TList *list;
8409d8c9 75 TH1F* h;
af2e514f 76 Int_t i;
8409d8c9 77
78 fNdets=ndets;
65fb704d 79
80 fQAList = new TObjArray(ndets);
81 TObjArray &hist = *fQAList;
82
83 char title[100];
84 //
85 TObjArray &mods = *(gAlice->Modules());
86 AliModule *mod;
87 TList *dir = gDirectory->GetList();
af2e514f 88 for (i=0; i<ndets; i++) {
65fb704d 89 hist[i] = list = new TList();
90 mod = (AliModule *) mods[i];
91
92 // Energy Spectrum
93 sprintf(title,"Spectrum entering: %s ",mod->GetName());
94 list->Add(new TH1F("hEnIn",strcat(title,mod->GetTitle()),100,-4,2));
95 dir->Remove(dir->FindObject("hEnIn"));
96
97 sprintf(title,"Spectrum exiting: %s ",mod->GetName());
98 list->Add(new TH1F("hEnOut",strcat(title,mod->GetTitle()),100,-4,2));
99 dir->Remove(dir->FindObject("hEnOut"));
100
101 // Z position
102 sprintf(title,"Z coordinate entering: %s ",mod->GetName());
103 list->Add(new TH1F("hZIn",strcat(title,mod->GetTitle()),100,
104 mod->ZMin(),mod->ZMax()));
105 dir->Remove(dir->FindObject("hZIn"));
106
107 sprintf(title,"Z coordinate exiting: %s ",mod->GetName());
108 list->Add(new TH1F("hZOut",strcat(title,mod->GetTitle()),100,
109 mod->ZMin(),mod->ZMax()));
110 dir->Remove(dir->FindObject("hZOut"));
111 }
112 //
113 gROOT->GetListOfBrowsables()->Add(this,"AliMCQA");
114
115 fDetDone = new Int_t[fNdets];
8409d8c9 116
117 //
118 // Global QA histograms
119 //
120 fQAHist = new TObjArray(2);
121 fNvolumes=gMC->NofVolumes();
122
123 fQAHist->Add(new TH1F("hMCVcalls","Monte Carlo calls per volume",
124 fNvolumes, 0.5, fNvolumes+0.5));
125 h = (TH1F*) dir->FindObject("hMCVcalls");
126 h->GetListOfFunctions()->Add(new TExec("ex","gAlice->GetMCQA()->AddVolumeName()"));
127 dir->Remove(dir->FindObject("hMCVcalls"));
128 //
129 // Build list of volume names
130 //
8409d8c9 131 fVolNames=new TObjArray(fNvolumes);
132 for(i=0;i<fNvolumes;++i) {
133 AliModule *mod = (AliModule*)
134 (*gAlice->Modules())[gAlice->DetFromMate(gMC->VolId2Mate(i+1))];
135 (*fVolNames)[i]=new TNamed(gMC->VolName(i+1),mod->GetName());
136 }
137
138 fQAHist->Add(new TH1F("hMCMcalls","Monte Carlo calls per module",
139 fNdets, -0.5, fNdets-0.5));
140 h = (TH1F*) dir->FindObject("hMCMcalls");
141 h->GetListOfFunctions()->Add(new TExec("ex","gAlice->GetMCQA()->AddModuleName()"));
142
143 dir->Remove(dir->FindObject("hMCMcalls"));
144 //
145 // Build list of module names
146 //
147 fModNames=new TObjArray(fNdets);
148 for(i=0;i<fNdets;++i)
149 (*fModNames)[i]=
150 new TNamed(((AliModule *)(*gAlice->Modules())[i])->GetName(),"");
65fb704d 151}
152
e460afec 153//_____________________________________________________________________________
154
155AliMCQA::~AliMCQA() {
156 if (fQAList) {
157 fQAList->Delete();
158 delete fQAList;
159 fQAList=0;
160 }
161 if (fQAHist) {
162 fQAHist->Delete();
163 delete fQAHist;
164 fQAHist=0;
165 }
166 if (fVolNames) {
167 fVolNames->Delete();
168 delete fVolNames;
169 fVolNames=0;
170 }
171 if (fModNames) {
172 fModNames->Delete();
173 delete fModNames;
174 fModNames=0;
175 }
176}
177
65fb704d 178//_____________________________________________________________________________
179void AliMCQA::Browse(TBrowser *b)
180{
181 //
182 // Called when the item "Run" is clicked on the left pane
183 // of the Root browser.
184 // It displays the Root Trees and all detectors.
185 //
8409d8c9 186 TH1 *hist;
187 //
188 // Global histos first
189 //
190 TIter global(fQAHist);
191 while((hist = (TH1*)global()))
192 b->Add(hist,hist->GetTitle());
193 //
194 // Module histograms now
195 //
65fb704d 196 TIter next(fQAList);
197 TList *histos;
65fb704d 198 while((histos = (TList*)next())) {
199 TIter next1(histos);
8409d8c9 200 while((hist = (TH1*)next1()))
65fb704d 201 b->Add(hist,hist->GetTitle());
65fb704d 202 }
203}
204
205//_____________________________________________________________________________
206void AliMCQA::PreTrack()
207{
8409d8c9 208 //
209 // Called before each track
210 //
65fb704d 211 fOldId=-1;
212 for(Int_t i=0; i<fNdets; i++) fDetDone[i]=0;
213}
214
215//_____________________________________________________________________________
216void AliMCQA::StepManager(Int_t id)
217{
8409d8c9 218 //
219 // Called at each step
220 //
221 TH1F *hist;
222 Int_t copy;
223 //
224 // Fill Global histograms first
225 //
226 hist = (TH1F*) fQAHist->FindObject("hMCVcalls");
227 hist->Fill(gMC->CurrentVolID(copy));
228 hist = (TH1F*) fQAHist->FindObject("hMCMcalls");
229 hist->Fill(id);
230 //
231 // Now the step manager histograms
232 //
65fb704d 233 if(fOldId != id) {
65fb704d 234 TLorentzVector p, x;
235 gMC->TrackMomentum(p);
236 gMC->TrackPosition(x);
237 Double_t energy = TMath::Max(
238 p[3]-gAlice->PDGDB()->GetParticle(gMC->TrackPid())->Mass(),1.e-12);
239 if(fOldId > -1) {
240 if(!fDetDone[fOldId] && !gMC->IsNewTrack()) {
241 TList *histold = (TList*) (*fQAList)[fOldId];
242 hist = (TH1F*) histold->FindObject("hEnOut");
243 hist->Fill(TMath::Log10(energy));
244 hist = (TH1F*) histold->FindObject("hZOut");
245 hist->Fill(x[2]);
246 fDetDone[fOldId]=1;
247 }
248 }
249 if(!fDetDone[id] && !gMC->IsNewTrack()) {
250 TList *histnew = (TList*) (*fQAList)[id];
251 hist = (TH1F*) histnew->FindObject("hEnIn");
252 hist->Fill(TMath::Log10(energy));
253 hist = (TH1F*) histnew->FindObject("hZIn");
254 hist->Fill(x[2]);
255 }
256 fOldId=id;
257 }
258}
8409d8c9 259
260//_____________________________________________________________________________
261void AliMCQA::AddModuleName()
262{
263 //
264 // Add function DrawModuleName to the module frequency histogram pad
265 //
266 static TVirtualPad* oldpad=0;
267 if(oldpad!=gPad) {
268 gPad->GetCanvas()->FeedbackMode(kTRUE);
269 if(gPad) gPad->AddExec("ex","gAlice->GetMCQA()->DrawModuleName()");
270 DrawPaveLabel(fMPaveLabel);
271 oldpad=gPad;
272 }
273}
274
275//_____________________________________________________________________________
276void AliMCQA::AddVolumeName()
277{
278 //
279 // Add function DrawVolumeName to the volume frequency histogram pad
280 //
281 static TVirtualPad* oldpad=0;
282 if(oldpad!=gPad) {
283 gPad->GetCanvas()->FeedbackMode(kTRUE);
284 if(gPad) gPad->AddExec("ex","gAlice->GetMCQA()->DrawVolumeName()");
285 DrawPaveLabel(fVPaveLabel);
286 oldpad=gPad;
287 }
288}
289
290//_____________________________________________________________________________
291void AliMCQA::DrawPaveLabel(TPaveLabel *&pv)
292{
293 //
294 // Draws the PaveLabel with the meaning of the bin
295 //
296 float uxmin = gPad->GetUxmin();
297 float uxmax = gPad->GetUxmax();
298 float uymin = gPad->GetUymin();
299 float uymax = gPad->GetUymax();
300 float lx = uxmax-uxmin;
301 float ly = uymax-uymin;
302
303 if(pv) delete pv;
304 pv
305 = new TPaveLabel(uxmin+0.05*lx,uymax-(0.05+0.1)*ly,
306 uxmin+(0.05+0.2)*lx,uymax-0.05*ly,
307 "");
308 pv->Draw();
309}
310
311//_____________________________________________________________________________
312Int_t AliMCQA::GetHBin(const char* hname)
313{
314 //
315 // Get the bin where the cursor is
316 //
317 TList *dir = gDirectory->GetList();
318 TH1 *h=(TH1*)dir->FindObject(hname);
319
320
321 int px = gPad->GetEventX();
322 Float_t upx = gPad->AbsPixeltoX(px);
323 Float_t x = gPad->PadtoX(upx);
324
325 return h->GetXaxis()->FindBin(x);
326}
327
328//_____________________________________________________________________________
329void AliMCQA::DrawModuleName()
330{
331 //
332 // Writes the name of the module of the bin where the cursor is
333 //
334 TObject *select = gPad->GetSelected();
335 if(!select) return;
336
337 Int_t binx = GetHBin("hMCMcalls");
338 if(0<binx && binx<=fNdets) {
339 char lab[15];
340 strcpy(lab,((TNamed*)(*fModNames)[binx-1])->GetName());
341 fMPaveLabel->SetLabel(lab);
342
343 gPad->Modified();
344 gPad->Update();
345 }
346}
347
348//_____________________________________________________________________________
349void AliMCQA::DrawVolumeName()
350{
351 //
352 // Writes the name of the volume:module of the bin where the cursor is
353 //
354 TObject *select = gPad->GetSelected();
355 if(!select) return;
356
357 Int_t binx = GetHBin("hMCVcalls");
358 if(0<binx && binx<=fNvolumes) {
359 char lab[20];
360 sprintf(lab,"%s: %s",((TNamed*)(*fVolNames)[binx-1])->GetName(),
361 ((TNamed*)(*fVolNames)[binx-1])->GetTitle());
362 fVPaveLabel->SetLabel(lab);
363
364 gPad->Modified();
365 gPad->Update();
366 }
367}
368
369