Replacing Header with Id
[u/mrichter/AliRoot.git] / STEER / AliMCQA.cxx
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 /* $Id$ */
17
18 ///////////////////////////////////////////////////////////////////////////////
19 //                                                                           //
20 // Quality Assurance class                                                   //
21 // Provides standard histograms for every element of the detector            //
22 //                                                                           //
23 ///////////////////////////////////////////////////////////////////////////////
24
25
26 #include <strings.h>
27
28 #include <TBrowser.h>
29 #include <TCanvas.h>
30 #include <TDatabasePDG.h>
31 #include <TExec.h>
32 #include <TH1.h>
33 #include <TList.h>
34 #include <TLorentzVector.h>
35 #include <TMath.h>
36 #include <TMath.h>
37 #include <TObjArray.h>
38 #include <TPad.h>
39 #include <TPaveLabel.h>
40 #include <TROOT.h>
41
42 #include "AliMCQA.h"
43 #include "AliModule.h"
44 #include "AliRun.h"
45
46 ClassImp(AliMCQA)
47
48 //_______________________________________________________________________
49 AliMCQA::AliMCQA():
50   fNdets(0),
51   fNvolumes(0),
52   fQAList(0),
53   fOldId(0),
54   fDetDone(0),
55   fQAHist(0),
56   fVolNames(0),
57   fModNames(0),
58   fMPaveLabel(0),
59   fVPaveLabel(0)
60 {
61   //
62   // Default constructor
63   //
64 }
65
66 //_______________________________________________________________________
67 AliMCQA::AliMCQA(const AliMCQA &qa):
68   TObject(qa),
69   fNdets(0),
70   fNvolumes(0),
71   fQAList(0),
72   fOldId(0),
73   fDetDone(0),
74   fQAHist(0),
75   fVolNames(0),
76   fModNames(0),
77   fMPaveLabel(0),
78   fVPaveLabel(0)
79 {
80   //
81   // Copy constructor
82   //
83   qa.Copy(*this);
84 }
85
86 //_______________________________________________________________________
87 AliMCQA::AliMCQA(Int_t ndets):
88   fNdets(ndets),
89   fNvolumes(gMC->NofVolumes()),
90   fQAList(new TObjArray(ndets)),
91   fOldId(0),
92   fDetDone(new Int_t[ndets]),
93   fQAHist(new TObjArray(2)),
94   fVolNames(new TObjArray(fNvolumes)),
95   fModNames(new TObjArray(fNdets)),
96   fMPaveLabel(0),
97   fVPaveLabel(0)
98 {
99   //
100   // Constructor, creates the list of lists of histograms
101   //
102   TList *list;
103   TH1F* h;
104   Int_t i;
105   
106   TObjArray &hist = *fQAList;
107
108   char title[100];
109   //
110   TObjArray &mods = *(gAlice->Modules());
111   TList *dir = gDirectory->GetList();
112   for (i=0; i<ndets; i++) {
113     hist[i] = list = new TList();
114     AliModule *mod = dynamic_cast<AliModule*>(mods[i]);
115
116     // Energy Spectrum
117     sprintf(title,"Spectrum entering: %s ",mod->GetName());
118     list->Add(new TH1F("hEnIn",strcat(title,mod->GetTitle()),100,-4,2));
119     dir->Remove(dir->FindObject("hEnIn"));
120
121     sprintf(title,"Spectrum exiting: %s ",mod->GetName());
122     list->Add(new TH1F("hEnOut",strcat(title,mod->GetTitle()),100,-4,2));
123     dir->Remove(dir->FindObject("hEnOut"));
124
125     // Z position
126     sprintf(title,"Z coordinate entering: %s ",mod->GetName());
127     list->Add(new TH1F("hZIn",strcat(title,mod->GetTitle()),100,
128                        mod->ZMin(),mod->ZMax()));
129     dir->Remove(dir->FindObject("hZIn"));
130
131     sprintf(title,"Z coordinate exiting: %s ",mod->GetName());
132     list->Add(new TH1F("hZOut",strcat(title,mod->GetTitle()),100,
133                        mod->ZMin(),mod->ZMax()));
134     dir->Remove(dir->FindObject("hZOut"));
135   }
136   //
137   gROOT->GetListOfBrowsables()->Add(this,"AliMCQA");
138
139   //
140   // Global QA histograms
141   //
142   
143   fQAHist->Add(new TH1F("hMCVcalls","Monte Carlo calls per volume",
144                         fNvolumes, 0.5, fNvolumes+0.5));
145   h = dynamic_cast<TH1F*>(dir->FindObject("hMCVcalls"));
146   h->GetListOfFunctions()->Add(new TExec("ex","gAlice->GetMCQA()->AddVolumeName()"));
147   dir->Remove(dir->FindObject("hMCVcalls"));
148   //
149   // Build list of volume names
150   //
151   for(i=0;i<fNvolumes;++i) {
152     AliModule *mod = dynamic_cast<AliModule*>
153       ((*gAlice->Modules())[gAlice->DetFromMate(gMC->VolId2Mate(i+1))]);
154     (*fVolNames)[i]=new TNamed(gMC->VolName(i+1),mod->GetName());
155   }
156
157   fQAHist->Add(new TH1F("hMCMcalls","Monte Carlo calls per module",
158                         fNdets, -0.5, fNdets-0.5));
159   h = dynamic_cast<TH1F*>(dir->FindObject("hMCMcalls"));
160    h->GetListOfFunctions()->Add(new TExec("ex","gAlice->GetMCQA()->AddModuleName()"));
161
162   dir->Remove(dir->FindObject("hMCMcalls"));
163   //
164   // Build list of module names
165   //
166   for(i=0;i<fNdets;++i) 
167     (*fModNames)[i]=
168       new TNamed((dynamic_cast<AliModule*>((*gAlice->Modules())[i]))->GetName(),"");
169 }
170
171 //_______________________________________________________________________
172 void AliMCQA::Copy(AliMCQA &) const
173 {
174   Fatal("Copy ctor","Not implemented!\n");
175 }
176
177 //_______________________________________________________________________
178 AliMCQA::~AliMCQA() 
179 {
180   //
181   // Destructor
182   //
183   gROOT->GetListOfBrowsables()->Remove(this);
184   if (fQAList) {
185     fQAList->Delete();
186     delete fQAList;
187     fQAList=0;
188   }
189   if (fQAHist) {
190     fQAHist->Delete();
191     delete fQAHist;
192     fQAHist=0;
193   }
194   if (fVolNames) {
195     fVolNames->Delete();
196     delete fVolNames;
197     fVolNames=0;
198   }
199   if (fModNames) {
200     fModNames->Delete();
201     delete fModNames;
202     fModNames=0;
203   }
204   delete [] fDetDone;
205   delete fMPaveLabel;
206   delete fVPaveLabel;
207 }
208
209 //_______________________________________________________________________
210 void AliMCQA::Browse(TBrowser *b)
211 {
212   //
213   // Called when the item "Run" is clicked on the left pane
214   // of the Root browser.
215   // It displays the Root Trees and all detectors.
216   //
217   TH1 *hist;
218   //
219   // Global histos first
220   //
221   TIter global(fQAHist);
222   while((hist = dynamic_cast<TH1*>(global())))
223     b->Add(hist,hist->GetTitle());
224   //
225   // Module histograms now
226   //
227   TIter next(fQAList);
228   TList *histos;
229   while((histos = dynamic_cast<TList*>(next()))) {
230     TIter next1(histos);
231     while((hist = dynamic_cast<TH1*>(next1())))
232       b->Add(hist,hist->GetTitle());
233   }
234 }
235
236 //_______________________________________________________________________
237 void AliMCQA::PreTrack()
238 {
239   //
240   // Called before each track
241   //
242   fOldId=-1;
243   for(Int_t i=0; i<fNdets; i++) fDetDone[i]=0;
244 }
245
246 //_______________________________________________________________________
247 void AliMCQA::StepManager(Int_t id)
248 {
249   //
250   // Called at each step
251   //
252   TH1F *hist;
253   Int_t copy;
254   //
255   // Fill Global histograms first
256   //
257
258
259   static TH1F* mcvcalls = dynamic_cast<TH1F*>(fQAHist->FindObject("hMCVcalls"));
260   mcvcalls->Fill(gMC->CurrentVolID(copy));
261   static TH1F* mcmcalls = dynamic_cast<TH1F*>(fQAHist->FindObject("hMCMcalls"));
262   mcmcalls->Fill(id);
263
264   //
265   // Now the step manager histograms
266   //
267   if(fOldId != id) {
268     static  Double_t mpi0=0;
269     static  Double_t mpip=0;
270     static  Double_t mpim=0;
271     static  Double_t mep=0;
272     static  Double_t mem=0;
273     Double_t mass=0;
274     Int_t num = gMC->TrackPid();
275
276     switch (num) {
277     case 111: 
278       if (mpi0==0) mpi0=gAlice->PDGDB()->GetParticle(num)->Mass(); 
279       mass=mpi0;
280       break;
281     case 211:
282       if (mpip==0) mpip=gAlice->PDGDB()->GetParticle(num)->Mass(); 
283       mass=mpip;
284       break;
285     case -211:
286       if (mpim==0) mpim=gAlice->PDGDB()->GetParticle(num)->Mass(); 
287       mass=mpim;
288       break;
289     case 11:
290       if (mep==0) mep=gAlice->PDGDB()->GetParticle(num)->Mass(); 
291       mass=mep;
292       break;
293     case -11:
294       if (mem==0) mem=gAlice->PDGDB()->GetParticle(num)->Mass(); 
295       mass=mem;
296       break;
297     default:
298       mass =gAlice->PDGDB()->GetParticle(num)->Mass();
299       break; 
300     }
301
302     static TLorentzVector p, x;
303     gMC->TrackMomentum(p);
304     gMC->TrackPosition(x);
305     Double_t energy = TMath::Max(p[3]-mass,1.e-12);
306     if(fOldId > -1) {
307       if(!fDetDone[fOldId] && !gMC->IsNewTrack()) {
308         TList *histold = dynamic_cast<TList*>((*fQAList)[fOldId]);
309         hist = dynamic_cast<TH1F*>(histold->FindObject("hEnOut"));
310         hist->Fill(TMath::Log10(energy));
311         hist = dynamic_cast<TH1F*>(histold->FindObject("hZOut"));
312         hist->Fill(x[2]);
313         fDetDone[fOldId]=1;
314       }
315     }
316     if(!fDetDone[id] && !gMC->IsNewTrack()) {
317       TList *histnew = dynamic_cast<TList*>((*fQAList)[id]);
318       hist = dynamic_cast<TH1F*>(histnew->FindObject("hEnIn"));
319       hist->Fill(TMath::Log10(energy));
320       hist = dynamic_cast<TH1F*>(histnew->FindObject("hZIn"));
321       hist->Fill(x[2]);
322     }
323     fOldId=id;
324   }
325 }
326
327 //_______________________________________________________________________
328 void AliMCQA::AddModuleName()
329 {
330   //
331   // Add function DrawModuleName to the module frequency histogram pad
332   //
333   static TVirtualPad* oldpad=0;
334   if(oldpad!=gPad) {
335     gPad->GetCanvas()->FeedbackMode(kTRUE);
336     if(gPad) gPad->AddExec("ex","gAlice->GetMCQA()->DrawModuleName()");
337     DrawPaveLabel(fMPaveLabel);
338     oldpad=gPad;
339   }
340 }
341
342 //_______________________________________________________________________
343 void AliMCQA::AddVolumeName()
344 {
345   //
346   // Add function DrawVolumeName to the volume frequency histogram pad
347   //
348   static TVirtualPad* oldpad=0;
349   if(oldpad!=gPad) {
350     gPad->GetCanvas()->FeedbackMode(kTRUE);
351     if(gPad) gPad->AddExec("ex","gAlice->GetMCQA()->DrawVolumeName()");
352     DrawPaveLabel(fVPaveLabel);
353     oldpad=gPad;
354   }
355 }
356
357 //_______________________________________________________________________
358 void AliMCQA::DrawPaveLabel(TPaveLabel *&pv)
359 {
360   //
361   // Draws the PaveLabel with the meaning of the bin
362   //
363   float uxmin = gPad->GetUxmin();
364   float uxmax = gPad->GetUxmax();
365   float uymin = gPad->GetUymin();
366   float uymax = gPad->GetUymax();
367   float lx = uxmax-uxmin;
368   float ly = uymax-uymin;
369
370   if(pv) delete pv;
371   pv 
372     = new TPaveLabel(uxmin+0.05*lx,uymax-(0.05+0.1)*ly,
373                      uxmin+(0.05+0.2)*lx,uymax-0.05*ly,
374                      "");
375   pv->Draw();
376 }
377
378 //_______________________________________________________________________
379 Int_t AliMCQA::GetHBin(const char* hname)
380 {
381   //
382   // Get the bin where the cursor is
383   //
384   TList *dir = gDirectory->GetList();
385   TH1 *h=dynamic_cast<TH1*>(dir->FindObject(hname));
386   
387
388   int px = gPad->GetEventX();
389   Float_t upx = gPad->AbsPixeltoX(px);
390   Float_t x = gPad->PadtoX(upx);
391     
392   return h->GetXaxis()->FindBin(x);
393 }
394
395 //_______________________________________________________________________
396 void AliMCQA::DrawModuleName()
397 {
398   //
399   // Writes the name of the module of the bin where the cursor is
400   //
401   TObject *select = gPad->GetSelected();
402   if(!select) return;
403   
404   Int_t binx = GetHBin("hMCMcalls");
405   if(0<binx && binx<=fNdets) {
406     char lab[15];
407     strcpy(lab,dynamic_cast<TNamed*>((*fModNames)[binx-1])->GetName());
408     fMPaveLabel->SetLabel(lab);
409   
410     gPad->Modified();
411     gPad->Update();
412   }
413 }
414
415 //_______________________________________________________________________
416 void AliMCQA::DrawVolumeName()
417 {
418   //
419   // Writes the name of the volume:module of the bin where the cursor is
420   //
421   TObject *select = gPad->GetSelected();
422   if(!select) return;
423
424   Int_t binx = GetHBin("hMCVcalls");
425   if(0<binx && binx<=fNvolumes) {
426     char lab[20];
427     sprintf(lab,"%s: %s",dynamic_cast<TNamed*>((*fVolNames)[binx-1])->GetName(),
428             dynamic_cast<TNamed*>((*fVolNames)[binx-1])->GetTitle());
429     fVPaveLabel->SetLabel(lab);
430     
431     gPad->Modified();
432     gPad->Update();
433   }
434 }
435
436