renaming function to avoid library conflict (discovered in test/generators/TUHKMgen)
[u/mrichter/AliRoot.git] / macros / analHits.C
1 void analHits (const char *filename="galice.root",Int_t evNumber=0, char *opt="Liny"){
2 /////////////////////////////////////////////////////////////////////////
3 //   This macro is a small example of a ROOT macro
4 //   illustrating how to read the output of GALICE
5 //   and fill some histograms.
6 //   
7 //     Root > .L analHits.C               //this loads the macro in memory
8 //     Root > analHits();                 //by default process first event   
9 //     Root > analHits("galice2.root",2); //process third event from 
10 //                                          galice2.root file.
11 //Begin_Html
12 /*
13 <img src="picts/analHits.gif">
14 */
15 //End_Html
16 /////////////////////////////////////////////////////////////////////////
17   if(gAlice){
18     delete gAlice;
19     gAlice=0;
20   }
21   else{
22     // Dynamically link some shared libs
23     if(gClassTable->GetID("AliRun") < 0) {
24       gROOT->LoadMacro("loadlibs.C");
25       loadlibs();
26     } // end if
27   }
28 // Connect the Root Galice file containing Geometry, Kine and Hits
29     TFile *file = (TFile*)gROOT->GetListOfFiles()->FindObject(filename);
30     if(!file) file = new TFile(filename);
31
32 // Get AliRun object from file or create it if not on file
33     if(!gAlice) {
34         gAlice = (AliRun*)file->Get("gAlice");
35         if(gAlice) printf("AliRun object found on file\n");
36         if(!gAlice) gAlice = new AliRun("gAlice","Alice test program");
37     } // end if !gAlice
38       
39 // Set event pointer to this event
40     Int_t nparticles = gAlice->GetEvent(evNumber);
41     if (nparticles <= 0){
42         cout << "No particles found for event " << evNumber;
43         cout << " in file " << filename << endl;
44         return;
45     } // end if nparticles <=0
46
47 // Pointer to specific detector hits.
48     AliFMDhit    *fmdHit;
49     AliITShit    *itsHit;
50     AliMUONHit   *muonHit;
51     AliPHOSHit   *phosHit;
52     AliPMDhit    *pmdHit;
53     AliHMPIDHit   *richHit;
54     AliT0hit  *startHit;
55     AliTOFhit    *tofHit;
56     AliTPChit    *tpcHit;
57     AliTPCTrackHits *tpc2Hit;
58     AliTRDhit    *trdHit;
59     AliCASTORhit *castorHit;
60     AliZDCHit    *zdcHit;
61     AliEMCALHit  *emcalHit;
62
63 // Get pointers to ALL Alice detectors and Hits containers
64     AliFMD    *FMD    = (AliFMD*)    gAlice->GetDetector("FMD");
65     AliITS    *ITS    = (AliITS*)    gAlice->GetDetector("ITS");
66     AliMUON   *MUON   = (AliMUON*)   gAlice->GetDetector("MUON");
67     AliPHOS   *PHOS   = (AliPHOS*)   gAlice->GetDetector("PHOS");
68     AliPMD    *PMD    = (AliPMD*)    gAlice->GetDetector("PMD");
69     AliHMPID   *HMPID   = (AliHMPID*)   gAlice->GetDetector("HMPID");
70     AliT0  *T0  = (AliT0*)  gAlice->GetDetector("T0");
71     AliTOF    *TOF    = (AliTOF*)    gAlice->GetDetector("TOF");
72     AliTPC    *TPC    = (AliTPC*)    gAlice->GetDetector("TPC");
73     AliTPC    *TPC    = (AliTPC*)    gAlice->GetDetector("TPC");
74     AliTRD    *TRD    = (AliTRD*)    gAlice->GetDetector("TRD");
75     AliCASTOR *CASTOR = (AliCASTOR*) gAlice->GetDetector("CASTOR");
76     AliZDC    *ZDC    = (AliZDC*)    gAlice->GetDetector("ZDC");
77     AliEMCAL  *EMCAL  = (AliEMCAL*)  gAlice->GetDetector("EMCAL");
78
79 // Get pointer to the particles
80 //    TClonesArray *Particles = gAlice->Particles();
81 //    TParticle    *part;
82
83     // Create histograms
84     if(FMD)   TH1F *hFMD   = new TH1F("hFMD"   ,"Hit Radius",100,0.,100.);
85     if(ITS)   TH1F *hITS   = new TH1F("hITS"   ,"Ionization",100,0.,3.e-3);
86     if(MUON)  TH1F *hMUON  = new TH1F("hMUON"  ,"Hit Radius",100,0.,500.);
87     if(PHOS)   TH1F *hPHOS  = new TH1F("hPHOS"  ,"Energy Dep.",100,0.,0.5);
88     if(PMD)     TH1F *hPMD   = new TH1F("hPMD"   ,"Energy Dep.",100,0.,1.e+5);
89     if(HMPID)  TH1F *hHMPID  = new TH1F("hHMPID"  ,"Energy loss",100,0.,1.e-5);
90     if(T0) TH1F *hT0 = new TH1F("hT0" ,"Time of Flight",100,0.,10.);
91     if(TOF)   TH1F *hTOF   = new TH1F("hTOF"   ,"Time of Flight",100,0.,1.e-5);
92     if(TPC)   TH1F *hTPC   = new TH1F("hTPC"   ,"Charge",100,0.,70.2);
93     if(TRD)   TH1F *hTRD   = new TH1F("hTRD"   ,"Charge",100,0.,10.);
94     if(CASTOR)TH1F *hCASTOR= new TH1F("hCASTOR","Hit Radius",100,0.,10.);
95     if(ZDC)   TH1F *hZDC   = new TH1F("hZDC"   ,"Energy",100,0.,5.);
96     if(EMCAL) TH1F *hEMCAL = new TH1F("hEMCAL" ,"Energy",100,0.,2.);
97 //    TH1F *hTPAR  = new TH1F("hTPAR" ,"?",6,1,7);  
98     Int_t track,ntracks = gAlice->TreeH()->GetEntries();
99 // Start loop on tracks in the hits containers
100     for(track=0; track<ntracks;track++){
101       //MI change
102       gAlice->ResetHits();
103       gAlice->TreeH()->GetEvent(track);
104       if(FMD){
105             for(fmdHit=(AliFMDhit*)FMD->FirstHit(-1);fmdHit;
106                 fmdHit=(AliFMDhit*)FMD->NextHit()){
107                 hFMD->Fill(TMath::Hypot(fmdHit->X(),fmdHit->Y()));
108             } // end for fmdHit
109         } // end if FMD
110         if(ITS){
111             for(itsHit=(AliITShit*)ITS->FirstHit(-1);itsHit;
112                 itsHit=(AliITShit*)ITS->NextHit()){
113                 if(itsHit->GetIonization()>0.0){//only after a step in the ITS
114                     hITS->Fill(itsHit->GetIonization());
115                 } // end if
116             } // end for itsHit
117         } // end if ITS
118         if(MUON){
119             for(muonHit=(AliMUONHit*)MUON->FirstHit(-1);muonHit;
120                 muonHit=(AliMUONHit*)MUON->NextHit()){
121                 hMUON->Fill(TMath::Hypot(muonHit->X(),muonHit->Y()));
122             } // end for muonHit
123         } // end if MUON
124         if(PHOS){
125             for(phosHit=(AliPHOSHit*)PHOS->FirstHit(-1);phosHit;
126                 phosHit=(AliPHOSHit*)PHOS->NextHit()){
127                 hPHOS->Fill(phosHit->GetEnergy());
128             } // end for phosHit
129         } // end if PHOS
130         if(PMD){
131             for(pmdHit=(AliPMDhit*)PMD->FirstHit(-1);pmdHit;
132                 pmdHit=(AliPMDhit*)PMD->NextHit()){
133                 hPMD->Fill(pmdHit->GetEnergy());
134             } // end for pmdHit
135         } // end if PMD
136         if(HMPID){
137             for(richHit=(AliHMPIDHit*)HMPID->FirstHit(-1);richHit;
138                 richHit=(AliHMPIDHit*)HMPID->NextHit()){
139                 hHMPID->Fill(richHit->fEloss);
140             } // end for richHit
141         } // end if HMPID
142         if(T0){
143             for(startHit=(AliT0hit*)T0->FirstHit(-1);startHit;
144                 startHit=(AliT0hit*)T0->NextHit()){
145                 hT0->Fill(startHit->fTime);
146             } // end for startHit
147         } // end if T0
148         if(TOF){
149             for(tofHit=(AliTOFhit*)TOF->FirstHit(-1);tofHit;
150                 tofHit=(AliTOFhit*)TOF->NextHit()){
151                 hTOF->Fill(tofHit->GetTof());
152             } // end for tofHit
153         } // end if TOF
154
155         if(TRD) {
156             for(trdHit=(AliTRDhit*)TRD->FirstHit(-1);trdHit;
157                 trdHit=(AliTRDhit*)TRD->NextHit()) {
158                 hTRD->Fill((Float_t)(trdHit->GetCharge()));
159             } // end for
160         } // end if TRD
161         if(CASTOR) {
162             for(castorHit=(AliCASTORhit*)CASTOR->FirstHit(-1);castorHit;
163                 castorHit=(AliCASTORhit*)CASTOR->NextHit()) {
164                 hCASTOR->Fill(TMath::Hypot(castorHit->X(),castorHit->Y()));
165             } // end for
166         } // end if CASTOR
167         if(ZDC){
168             for(zdcHit=(AliZDCHit*)ZDC->FirstHit(-1);zdcHit;
169                 zdcHit=(AliZDCHit*)ZDC->NextHit()){
170                 hZDC->Fill(zdcHit->GetEnergy());
171             } // end for zdcdHit
172         } // end if ZDC
173         if(TPC) {        
174           for(tpcHit=(AliTPChit*)TPC->FirstHit(-1);tpcHit;
175               tpcHit=(AliTPChit*)TPC->NextHit()) {
176             hTPC->Fill((Float_t)(tpcHit->fQ));
177           } // end for tpcHit
178         } // end if TPC
179         if(EMCAL) {
180           for(emcalHit=(AliEMCALHit*)EMCAL->FirstHit(-1);emcalHit;
181               emcalHit=(AliEMCALHit*)EMCAL->NextHit()) {
182             hEMCAL->Fill((Float_t)(emcalHit->GetEnergy()));
183           } // end for tpcHit
184         } // end if TPC
185
186     } // end for track
187
188 //Create a canvas, set the view range, show histograms
189     TCanvas *c0 = new TCanvas("c0","Alice Detectors",400,10,600,700);
190     if(opt=="Logy")c0->SetLogy();
191     if(FMD){
192         hFMD->SetFillColor(42);
193         hFMD->Draw();
194         c0->Print("analHitsFMD.ps");
195     } // end if FMD
196     if(ITS){
197         hITS->SetFillColor(42);
198         hITS->Draw();
199         c0->SaveAs("analHitsITS.ps");
200     } // end if ITS
201     if(MUON){
202         hMUON->SetFillColor(42);
203         hMUON->Draw();
204         c0->SaveAs("analHitsMUON.ps");
205     } // end if MUON
206     if(PHOS){
207         hPHOS->SetFillColor(42);
208         hPHOS->Draw();
209         c0->SaveAs("analHitsPHOS.ps");
210     } // end if PHOS
211     if(PMD){
212         hPMD->SetFillColor(42);
213         hPMD->Draw();
214         c0->SaveAs("analHitsPMD.ps");
215     } // end if PMD
216     if(HMPID){
217         hHMPID->SetFillColor(42);
218         hHMPID->Draw();
219         c0->SaveAs("analHitsHMPID.ps");
220     } // end if HMPID
221     if(T0){
222         hT0->SetFillColor(42);
223         hT0->Draw();
224         c0->SaveAs("analHitsT0.ps");
225     } // end if T0
226     if(TOF){
227         hTOF->SetFillColor(42);
228         hTOF->Draw();
229         c0->SaveAs("analHitsTOF.ps");
230     } // end if TOF
231     if(TPC){
232         hTPC->SetFillColor(42);
233         hTPC->Draw();
234         c0->SaveAs("analHitsTPC.ps");
235     } // end if TPC
236     if(TRD){
237         hTRD->SetFillColor(42);
238         hTRD->Draw();
239         c0->SaveAs("analHitsTRD.ps");
240     } // end if TRD
241     if(CASTOR){
242         hCASTOR->SetFillColor(42);
243         hCASTOR->Draw();
244         c0->SaveAs("analHitsCASTOR.ps");
245     } // end if TRD
246     if(ZDC){
247         hZDC->SetFillColor(42);
248         hZDC->Draw();
249         c0->SaveAs("analHitsZDC.ps");
250     } // end if ZDC
251     if(EMCAL){
252         hEMCAL->SetFillColor(42);
253         hEMCAL->Draw();
254         c0->SaveAs("analHitsEMCAL.ps");
255     } // end if ZDC
256
257 // Clean Up
258     /*
259     if(FMD)    delete hFMD;
260     if(ITS)    delete hITS;
261     if(MUON)   delete hMUON;
262     if(PHOS)   delete hPHOS;
263     if(PMD)    delete hPMD;
264     if(HMPID)   delete hHMPID;
265     if(T0)  delete hT0;
266     if(TOF)    delete hTOF;
267     if(TPC)    delete hTPC;
268     if(TRD)    delete hTRD;
269     if(CASTOR) delete hCASTOR;
270     if(ZDC)    delete hZDC;
271     */
272 }