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