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