Back to the previous version. The warning was erronously generated
[u/mrichter/AliRoot.git] / TOF / AliTOFanalyzeHits.C
CommitLineData
32bd3eb5 1Int_t AliTOFanalyzeHits(Int_t nevents = 1, Bool_t drawing = kFALSE)
42036c57 2{
3 //
4 // Analyzes the hits and fills QA-histograms
32bd3eb5 5 // Use case:
6 // start aliroot
7 // root [0] .L AliTOFanalyzeHits.C
8 // root [1] AliTOFanalyzeHits()
42036c57 9 //
32bd3eb5 10 // Updated to the new I/O by: A. De Caro, C. Zampolli
42036c57 11
12 Int_t rc = 0;
32bd3eb5 13
42036c57 14 // Define the histograms
15 // x,y,z, rho, tof, padx, padz, sector, plate, strip, (x vs y)
16
17 // hit-map in a plane
18 TH2F *h2hitMap = new TH2F("h2hitMap","Hit Map (projection on the plane)",2500,-12500.,12500.,800,-400.,400.);
19 // time of flight distribution for primaries and secondaries
20 TH1F *htofp = new TH1F("htofp","Time of Flight (primaries)",800,0.,80.);
21 TH1F *htofs = new TH1F("htofs","Time of Flight (secondaries)",800,0.,80.);
22 // momentum when striking the TOF for primaries and secondaries
23 TH1F *htofmomp = new TH1F("htofmomp","Momentum at TOF (primaries)",100,0.,10.);
24 TH1F *htofmoms = new TH1F("htofmoms","Momentum at TOF (secondaries)",100,0.,10.);
25 // TOF hit volumes
26 TH1F *hsector = new TH1F("hsector","Sector",20,0.,20.);
27 TH1F *hplate = new TH1F("hplate","Plate ", 6,0., 6.);
28 TH1F *hstrip = new TH1F("hstrip","Strip ",25,0.,25.);
29 TH1F *hpadz = new TH1F("hpadz","Pad along z ",3,0.,3.);
30 TH1F *hpadx = new TH1F("hpadx","Pad along x",50,0.,50.);
31 // track length when striking the TOF (used by AliTOFT0)
32 TH1F *htrackLenp= new TH1F("htrackLenp","Track Length on TOF for Primaries",800,0.,800.);
33
32bd3eb5 34 // Histograms added to control the right TOF element numbering:
35 // it should be increasing with the azimuthal and polar angles
36
37 TH2F *hmoduleVStheta = new TH2F("hmoduleVStheta", "hmoduleVStheta", 180,0.,180.,6,0,6);
38 TH2F *hsectorVSphi = new TH2F("hsectorVSphi", "hsectorVSphi", 360,0.,360.,19,0,19);
39 TH2F *hstripVStheta = new TH2F("hstripVStheta", "hstripVStheta", 180,0.,180.,25,0,25);
40 TH2F *hpadzVStheta = new TH2F("hpadzVStheta", "hpadzVStheta", 180,0.,180.,3,0,3);
41 TH2F *hpadxVSphi = new TH2F("hpadxVSphi", "hpadxVSphi", 360,0.,360.,49,0,49);
42
43 // Dynamically link some shared libs
44 if (gClassTable->GetID("AliRun") < 0) {
45 gROOT->LoadMacro("loadlibs.C");
46 loadlibs();
47 }
48
49 if (gAlice)
50 {
51 delete gAlice->GetRunLoader();
52 delete gAlice;
53 gAlice = 0x0;
54 }
55
56 AliRunLoader *rl = AliRunLoader::Open("galice.root",AliConfig::fgkDefaultEventFolderName,"read");
57 if (!rl)
58 {
59 cerr<<"Can't load RunLoader from file"<<"!\n";
60 rc = 1;
61 return rc;
62 }
63
64 rl->LoadgAlice();
65 gAlice = rl->GetAliRun();
66
67 if (!gAlice)
68 {
69 cerr << "<AliTOFanalyzeHits> AliRun object not found on file\n ";
70 rc = 2;
71 return rc;
72 }
73
74 // Get the pointer to the TOF detector
75 AliTOF *tof = (AliTOF *) gAlice->GetDetector("TOF");
76 AliLoader *tofl = rl->GetLoader("TOFLoader");
77 if (tof == 0x0 || tofl == 0x0) {
78 cerr << "<AliTOFanalyzeHits> Can not find TOF or TOFLoader\n";
79 rc = 3;
42036c57 80 return rc;
81 }
82
83 Int_t countHits = 0;
32bd3eb5 84
85 rl->LoadHeader();
86
87 for (Int_t ievent=0; ievent<nevents; ievent++) {
88 printf ("Processing event %d \n", ievent);
89 rl->GetEvent(ievent);
90
91 // Get the pointer Hit tree
92 tofl->LoadHits();
93 TTree *hitTree = tofl->TreeH();
94 tof->SetTreeAddress();
95 if (!hitTree) {
96 cout << "<AliTOFanalyzeHits> No TreeH found" << endl;
97 rc = 4;
98 return rc;
99 }
100
101 rl->LoadKinematics();
102 //AliStack* stack = rl->Stack(); // it is not necessary to use the stack!
103
104 // Get the number of entries in the hit tree
105 // (Number of primary particles creating a hit somewhere)
106 Int_t nTrack = (Int_t) hitTree->GetEntries();
107 cout << "<AliTOFanalyzeHits> Found " << nTrack
108 << " primary particles with hits" << endl;
109
110 Int_t nPrimaryOnTof = 0;
111 Int_t nSecondaryOnTof = 0;
112 Int_t nelectron = 0;
113 Int_t npion = 0;
114 Int_t nkaon = 0;
115 Int_t nproton = 0;
116 Int_t nmuon = 0;
42036c57 117
32bd3eb5 118 // Loop through all entries in the tree
119 for (Int_t iTrack = 0; iTrack < nTrack; iTrack++) {
120
121 tof->ResetHits();
122 hitTree->GetEvent(iTrack);
123
124 // Loop through the TOF hits
125 AliTOFhitT0 *hit = (AliTOFhitT0 *) tof->FirstHit(-1);
126 while (hit) {
127
128 countHits++;
129
130 Float_t x = hit->X();
131 Float_t y = hit->Y();
132 Float_t z = hit->Z();
133 Float_t phiAngle=TMath::Pi() + TMath::ATan2(-y,-x);
134 Float_t rhoRadius=TMath::Sqrt(x*x+y*y);
135 Float_t thetaAngle=TMath::Pi() + TMath::ATan2(-rhoRadius,-z);
136 Float_t dummy=rhoRadius*phiAngle;
137 h2hitMap->Fill(dummy,z);
138
139 phiAngle*=180./TMath::Pi();
140 thetaAngle*=180./TMath::Pi();
141
142 Float_t flightTime = hit->GetTof(); // [s]
143 flightTime *= 1.e+09; // convert in [ns]
144 Float_t angle = hit->GetIncA();
145 Float_t tofmom = hit->GetMom(); // [GeV/c]
146 Float_t trackLen = hit->GetLen(); // [cm]
147
148 // TOF hit volumes
149 Int_t sector = hit->GetSector(); // range [1-18]
150 Int_t plate = hit->GetPlate(); // range [1- 5]
151 Int_t strip = hit->GetStrip(); // range [1-20]
152 Int_t padz = hit->GetPadz(); // range [1- 2]
153 Int_t padx = hit->GetPadx(); // range [1-48]
154 // it is QA, then I perform QA!
155 Bool_t isHitBad = (sector<1 || sector>18 || plate<1 || plate >5 || padz<1 || padz>2 || padx<1 || padx>48);
156
157 if (isHitBad) {
158 cout << "<AliTOFanalyzeHits> strange hit found" << endl;
159 rc = 5;
160 return rc;
161 }
162
163 hmoduleVStheta->Fill(thetaAngle,plate);
164 hstripVStheta->Fill(thetaAngle,strip);
165 hsectorVSphi->Fill(phiAngle,sector);
166 hpadxVSphi->Fill(phiAngle,padx);
167 hpadzVStheta->Fill(thetaAngle,padz);
168
169 // filling hit volume histos
170 hsector->Fill(sector);
171 hplate->Fill(plate);
172 hstrip->Fill(strip);
173 hpadx->Fill(padx);
174 hpadz->Fill(padz);
175
176 Int_t track = hit->Track();
177 TParticle *part = gAlice->Particle(track);
178
179 // getting MC info for the current track
180 if (part->GetFirstMother()<0){
181 Int_t icod = TMath::Abs(part->GetPdgCode());
182 switch (icod) {
183 case 211:
184 npion++;
185 break ;
186 case 321:
187 nkaon++;
188 break ;
189 case 2212:
190 nproton++;
191 break ;
192 case 11:
193 nelectron++;
194 break ;
195 case 13:
196 nmuon++;
197 break ;
198 }
199 htofp->Fill(flightTime);
200 htofmomp->Fill(tofmom);
201 htrackLenp->Fill(trackLen);
202 } else {
203 htofs->Fill(flightTime);
204 htofmoms->Fill(tofmom);
205 }
206
207 // go to next hit
208 hit = (AliTOFhitT0 *) tof->NextHit();
42036c57 209
32bd3eb5 210 }
42036c57 211
212 }
213
32bd3eb5 214 tofl->UnloadHits();
215 rl->UnloadKinematics();
42036c57 216
32bd3eb5 217 cout << "<AliTOFanalyzeHits> Found " << countHits << " hits in total" << endl;
218 cout << npion << " primary pions reached the TOF detector" << endl;
219 cout << nkaon << " primary kaons reached the TOF detector" << endl;
220 cout << nproton << " primary protons reached the TOF detector" << endl;
221 cout << nelectron << " primary electrons reached the TOF detector" << endl;
222 cout << nmuon << " primary muons reached the TOF detector" << endl;
42036c57 223
32bd3eb5 224 }
42036c57 225
32bd3eb5 226 rl->UnloadHeader();
227 rl->UnloadgAlice();
228
229 if (drawing) {
230 TCanvas *cHits = new TCanvas("cHits","AliTOFanalyzeHits hit volumes",50,50,900,900);
231 cHits->Divide(3,2);
232 cHits->cd(1);
233 hsector->Draw();
234 cHits->cd(2);
235 hplate->Draw();
236 cHits->cd(3);
237 hstrip->Draw();
238 cHits->cd(4);
239 hpadz->Draw();
240 cHits->cd(5);
241 hpadx->Draw();
242
243 TCanvas *chitmap = new TCanvas("chitmap","AliTOFanalyzeHits Hit Map",50,50,600,600);
244 chitmap->cd();
245 h2hitMap->Draw();
246
247 TCanvas *ctrackLen = new TCanvas("ctrackLen","AliTOFanalyzeHits Track Length for primaries on TOF",50,50,400,400);
248 ctrackLen->cd();
249 htrackLenp->Draw();
250
251 TCanvas *ctofmom = new TCanvas("ctofmom","AliTOFanalyzeHits flight times",50,50,700,700);
252 ctofmom->Divide(2,2);
253 ctofmom->cd(1);
254 gPad->SetLogy();
255 htofp->Draw();
256 ctofmom->cd(2);
257 gPad->SetLogy();
258 htofs->Draw();
259 ctofmom->cd(3);
260 gPad->SetLogy();
261 htofmomp->Draw();
262 ctofmom->cd(4);
263 gPad->SetLogy();
264 htofmoms->Draw();
265 }
266
42036c57 267 // save histos into file TOF_hitsQA.root
268 TFile *fout = new TFile("TOF_hitsQA.root","RECREATE");
269 h2hitMap->Write();
270 htofp->Write();
271 htofs->Write();
272 htofmomp->Write();
273 htofmoms->Write();
274 hsector->Write();
275 hplate->Write();
276 hstrip->Write();
277 hpadz->Write();
278 hpadx->Write();
279 htrackLenp->Write();
32bd3eb5 280
281 hmoduleVStheta->Write();
282 hsectorVSphi->Write();
283 hstripVStheta->Write();
284 hpadzVStheta->Write();
285 hpadxVSphi->Write();
286
42036c57 287 fout->Close();
288
32bd3eb5 289 cout << " Finished AliTOFanalizeHits \n";
290
291 if (gAlice)
292 {
293 delete gAlice->GetRunLoader();
294 delete gAlice;
295 gAlice = 0x0;
296 }
297
42036c57 298 return rc;
299
300}