]>
Commit | Line | Data |
---|---|---|
8e64dd77 | 1 | |
793ff80c | 2 | Int_t AliTRDtest() |
3 | { | |
4 | // | |
5 | // Test macro for the TRD code | |
6 | // | |
7 | ||
8 | Int_t rc = 0; | |
9 | ||
8e64dd77 | 10 | TFile *file; |
11 | ||
793ff80c | 12 | // Initialize the test setup |
5931c105 | 13 | gAlice->Init("$(ALICE_ROOT)/TRD/AliTRDconfig.C"); |
793ff80c | 14 | |
15 | // Run one event and create the hits | |
16 | gAlice->Run(1); | |
17 | ||
1f4ee546 | 18 | //if (gAlice) delete gAlice; |
19 | //file = (TFile *) gROOT->GetListOfFiles()->FindObject("TRD_test.root"); | |
20 | //gAlice = (AliRun *) file->Get("gAlice"); | |
6244debe | 21 | |
793ff80c | 22 | // Analyze the TRD hits |
793ff80c | 23 | if (rc = AliTRDanalyzeHits()) return rc; |
24 | ||
25 | // Run the digitization | |
793ff80c | 26 | if (rc = AliTRDcreateDigits()) return rc; |
27 | ||
1f4ee546 | 28 | // if (gAlice) delete gAlice; |
29 | // file = (TFile *) gROOT->GetListOfFiles()->FindObject("TRD_test.root"); | |
30 | // gAlice = (AliRun *) file->Get("gAlice"); | |
6244debe | 31 | |
a328fff9 | 32 | // Analyze the digits |
33 | if (rc = AliTRDanalyzeDigits()) return rc; | |
793ff80c | 34 | |
1f4ee546 | 35 | // // Create the cluster |
36 | // if (rc = AliTRDcreateCluster()) return rc; | |
4a0fe73c | 37 | |
1f4ee546 | 38 | // if (gAlice) delete gAlice; |
39 | // file = (TFile *) gROOT->GetListOfFiles()->FindObject("TRD_test.root"); | |
40 | // gAlice = (AliRun *) file->Get("gAlice"); | |
4a0fe73c | 41 | |
1f4ee546 | 42 | // // Analyze the cluster |
43 | // if (rc = AliTRDanalyzeCluster()) return rc; | |
4a0fe73c | 44 | |
8e64dd77 | 45 | //file = (TFile *) gROOT->GetListOfFiles()->FindObject("TRD_test.root"); |
46 | //file->Close(); | |
47 | ||
48 | // Find the tracks | |
49 | //if (rc = AliTRDcreateTracks()) return rc; | |
50 | ||
793ff80c | 51 | return rc; |
52 | ||
53 | } | |
fa148e6c | 54 | |
55 | //_____________________________________________________________________________ | |
56 | Int_t AliTRDanalyzeHits() | |
57 | { | |
58 | // | |
59 | // Analyzes the hits and fills QA-histograms | |
60 | // | |
61 | ||
62 | Int_t rc = 0; | |
63 | ||
1f4ee546 | 64 | AliRunLoader *rl = gAlice->GetRunLoader(); |
65 | if (!rl) { | |
66 | cout << "<AliTRDanalyzeHits> No RunLoader found" << endl; | |
fa148e6c | 67 | rc = 1; |
68 | return rc; | |
69 | } | |
1f4ee546 | 70 | |
71 | // Import the Trees for the event nEvent in the file | |
72 | rl->GetEvent(0); | |
73 | ||
74 | AliLoader* loader = rl->GetLoader("TRDLoader"); | |
75 | if (!loader) { | |
76 | cout << "<AliTRDanalyzeHits> No TRDLoader found" << endl; | |
77 | rc = 2; | |
78 | return rc; | |
79 | } | |
fa148e6c | 80 | |
81 | // Get the pointer to the TRD detector | |
82 | AliTRD *trd = (AliTRD *) gAlice->GetDetector("TRD"); | |
83 | if (!trd) { | |
84 | cout << "<AliTRDanalyzeHits> No TRD detector found" << endl; | |
85 | rc = 2; | |
86 | return rc; | |
87 | } | |
88 | ||
89 | // Get the pointer to the geometry object | |
90 | AliTRDgeometryFull *geo; | |
91 | if (trd) { | |
92 | geo = (AliTRDgeometryFull *) trd->GetGeometry(); | |
93 | } | |
94 | else { | |
95 | cout << "<AliTRDanalyzeHits> No TRD geometry found" << endl; | |
96 | rc = 3; | |
97 | return rc; | |
98 | } | |
99 | ||
3551db50 | 100 | AliTRDCommonParam *par = AliTRDCommonParam::Instance(); |
fa148e6c | 101 | |
102 | // Define the histograms | |
103 | TH1F *hQdedx = new TH1F("hQdedx","Charge dedx-hits",100,0.0,1000.0); | |
104 | TH1F *hQtr = new TH1F("hQtr" ,"Charge TR-hits" ,100,0.0,1000.0); | |
105 | ||
106 | Float_t rmin = geo->Rmin(); | |
107 | Float_t rmax = geo->Rmax(); | |
108 | Float_t length = geo->GetChamberLength(0,2); | |
109 | Float_t width = geo->GetChamberWidth(0); | |
110 | Int_t ncol = par->GetColMax(0); | |
111 | Int_t nrow = par->GetRowMax(0,2,13); | |
3551db50 | 112 | Int_t ntime = ((Int_t) (rmax - rmin) / 22.0); |
fa148e6c | 113 | |
114 | TH2F *hZY = new TH2F("hZY" ,"Y vs Z (chamber 0)", nrow,-length/2.,length/2. | |
115 | ,ntime, rmin, rmax); | |
116 | TH2F *hXZ = new TH2F("hXZ" ,"Z vs X (plane 0)" , ncol, -width/2., width/2. | |
117 | , nrow,-length/2.,length/2.); | |
118 | ||
fa148e6c | 119 | // Get the pointer hit tree |
1f4ee546 | 120 | TTree *hitTree = loader->TreeH(); |
fa148e6c | 121 | if (!hitTree) { |
122 | cout << "<AliTRDanalyzeHits> No hit tree found" << endl; | |
123 | rc = 4; | |
124 | return rc; | |
125 | } | |
126 | ||
127 | Int_t countHits = 0; | |
128 | Int_t nBytes = 0; | |
129 | ||
130 | // Get the number of entries in the hit tree | |
131 | // (Number of primary particles creating a hit somewhere) | |
132 | Int_t nTrack = (Int_t) hitTree->GetEntries(); | |
133 | cout << "<AliTRDanalyzeHits> Found " << nTrack | |
134 | << " primary particles with hits" << endl; | |
135 | ||
136 | // Loop through all entries in the tree | |
137 | for (Int_t iTrack = 0; iTrack < nTrack; iTrack++) { | |
138 | ||
139 | gAlice->ResetHits(); | |
140 | nBytes += hitTree->GetEvent(iTrack); | |
141 | ||
fa148e6c | 142 | // Loop through the TRD hits |
143 | Int_t iHit = 0; | |
144 | AliTRDhit *hit = (AliTRDhit *) trd->FirstHit(-1); | |
145 | while (hit) { | |
146 | ||
147 | countHits++; | |
148 | iHit++; | |
149 | ||
150 | Float_t x = hit->X(); | |
151 | Float_t y = hit->Y(); | |
152 | Float_t z = hit->Z(); | |
153 | Float_t q = hit->GetCharge(); | |
154 | Int_t track = hit->Track(); | |
155 | Int_t det = hit->GetDetector(); | |
156 | Int_t plane = geo->GetPlane(det); | |
157 | ||
158 | if (q > 0) { | |
159 | hQdedx->Fill(q); | |
160 | hZY->Fill(z,y); | |
161 | if (plane == 0) { | |
162 | hXZ->Fill(x,z); | |
163 | } | |
164 | } | |
165 | else if (q < 0) { | |
166 | hQtr->Fill(TMath::Abs(q)); | |
167 | } | |
168 | ||
fa148e6c | 169 | hit = (AliTRDhit *) trd->NextHit(); |
170 | ||
171 | } | |
172 | ||
fa148e6c | 173 | } |
174 | ||
175 | cout << "<AliTRDanalyzeHits> Found " << countHits << " hits in total" << endl; | |
176 | ||
177 | TCanvas *cHits = new TCanvas("cHits","AliTRDanalyzeHits",50,50,600,600); | |
178 | cHits->Divide(2,2); | |
179 | cHits->cd(1); | |
180 | hXZ->Draw("COL"); | |
181 | cHits->cd(2); | |
182 | hZY->Draw("COL"); | |
183 | cHits->cd(3); | |
184 | gPad->SetLogy(); | |
185 | hQdedx->Draw(); | |
186 | cHits->cd(4); | |
187 | gPad->SetLogy(); | |
188 | hQtr->Draw(); | |
189 | ||
fa148e6c | 190 | return rc; |
191 | ||
192 | } | |
193 | ||
194 | //_____________________________________________________________________________ | |
195 | Int_t AliTRDcreateDigits() | |
196 | { | |
197 | // | |
198 | // Creates the digits from the hits of the slow simulator | |
199 | // | |
200 | ||
201 | Int_t rc = 0; | |
202 | ||
203 | if (!gAlice) { | |
204 | cout << "<AliTRDcreateDigits> No AliRun object found" << endl; | |
205 | rc = 1; | |
206 | return rc; | |
207 | } | |
208 | gAlice->GetEvent(0); | |
209 | ||
210 | // Create the TRD digitzer | |
211 | AliTRDdigitizer *digitizer = new AliTRDdigitizer("TRDdigitizer" | |
212 | ,"TRD digitizer class"); | |
1f4ee546 | 213 | digitizer->Open("TRD_test.root"); |
fa148e6c | 214 | digitizer->InitDetector(); |
1f4ee546 | 215 | digitizer->InitOutput(0); |
fa148e6c | 216 | |
217 | // Set the parameter | |
218 | digitizer->SetDebug(1); | |
219 | ||
220 | // Define the parameter object | |
221 | // If no external parameter object is defined, | |
222 | // default parameter will be used | |
223 | AliTRDparameter *parameter = new AliTRDparameter("TRDparameter" | |
224 | ,"TRD parameter class"); | |
225 | digitizer->SetParameter(parameter); | |
226 | ||
227 | // Create the digits | |
228 | if (!(digitizer->MakeDigits())) { | |
229 | rc = 2; | |
230 | return rc; | |
231 | } | |
232 | ||
233 | // Write the digits into the input file | |
1f4ee546 | 234 | //if (!(digitizer->MakeBranch())) { |
235 | // rc = 3; | |
236 | // return rc; | |
237 | //} | |
fa148e6c | 238 | |
239 | // Write the digits into the input file | |
240 | if (!(digitizer->WriteDigits())) { | |
241 | rc = 4; | |
242 | return rc; | |
243 | } | |
244 | ||
245 | // Save the parameter object in the AliROOT file | |
246 | if (!(parameter->Write())) { | |
247 | rc = 4; | |
248 | return rc; | |
249 | } | |
250 | ||
251 | return rc; | |
252 | ||
253 | } | |
254 | ||
255 | //_____________________________________________________________________________ | |
256 | Int_t AliTRDanalyzeDigits() | |
257 | { | |
258 | // | |
259 | // Analyzes the digits | |
260 | // | |
261 | ||
262 | Int_t rc = 0; | |
263 | ||
264 | const Int_t kNpla = AliTRDgeometry::Nplan(); | |
265 | ||
266 | if (!gAlice) { | |
267 | cout << "<AliTRDanalyzeDigits> No AliRun object found" << endl; | |
268 | rc = 1; | |
269 | return rc; | |
270 | } | |
fa148e6c | 271 | |
a328fff9 | 272 | AliRunLoader *rl = gAlice->GetRunLoader(); |
273 | if (!rl) { | |
274 | cout << "<AliTRDanalyzeHits> No RunLoader found" << endl; | |
fa148e6c | 275 | rc = 2; |
276 | return rc; | |
277 | } | |
278 | ||
a328fff9 | 279 | // Import the Trees for the event nEvent in the file |
280 | rl->LoadKinematics(); | |
281 | rl->GetEvent(0); | |
282 | rl->GetHeader(); | |
283 | ||
284 | AliLoader* loader = rl->GetLoader("TRDLoader"); | |
285 | if (!loader) { | |
286 | cout << "<AliTRDanalyzeHits> No TRDLoader found" << endl; | |
fa148e6c | 287 | rc = 3; |
288 | return rc; | |
289 | } | |
290 | ||
a328fff9 | 291 | // Get the pointer to the TRD detector |
292 | AliTRD *trd = (AliTRD *) gAlice->GetDetector("TRD"); | |
293 | if (!trd) { | |
294 | cout << "<AliTRDanalyzeDigits> No TRD detector found" << endl; | |
295 | rc = 4; | |
296 | return rc; | |
297 | } | |
298 | ||
299 | // Get the parameter object | |
3551db50 | 300 | AliTRDSimParam *parSim = AliTRDSimParam::Instance(); |
301 | AliTRDCommonParam *parCom = AliTRDCommonParam::Instance(); | |
302 | AliTRDcalibDB *cal = AliTRDcalibDB::Instance(); | |
a328fff9 | 303 | |
fa148e6c | 304 | // Define the histograms |
3551db50 | 305 | Int_t adcRange = ((Int_t) parSim->GetADCoutRange()); |
fa148e6c | 306 | TH1F *hAmpAll = new TH1F("hAmpAll" ,"Amplitude of the digits (all)" |
307 | ,adcRange+1,-0.5,((Float_t) adcRange)+0.5); | |
308 | TH1F *hAmpEl = new TH1F("hAmpEl" ,"Amplitude of the digits (electrons)" | |
309 | ,adcRange+1,-0.5,((Float_t) adcRange)+0.5); | |
310 | TH1F *hAmpPi = new TH1F("hAmpPi" ,"Amplitude of the digits (pions)" | |
311 | ,adcRange+1,-0.5,((Float_t) adcRange)+0.5); | |
312 | TH1F *hAmpNoise = new TH1F("hAmpNoise","Amplitude of the digits (noise)" | |
313 | ,5,-0.5,4.5); | |
314 | ||
315 | // Get the pointer to the geometry object | |
316 | AliTRDgeometry *geo; | |
317 | if (trd) { | |
318 | geo = trd->GetGeometry(); | |
319 | } | |
320 | else { | |
321 | cout << "<AliTRDanalyzeDigits> No TRD geometry found" << endl; | |
a328fff9 | 322 | rc = 5; |
fa148e6c | 323 | return rc; |
324 | } | |
325 | ||
326 | // Create the digits manager | |
327 | AliTRDdigitsManager *digitsManager = new AliTRDdigitsManager(); | |
8e64dd77 | 328 | digitsManager->SetDebug(1); |
fa148e6c | 329 | |
fa148e6c | 330 | // Read the digits from the file |
a328fff9 | 331 | if (!(digitsManager->ReadDigits(loader->TreeD()))) { |
fa148e6c | 332 | cout << "<AliTRDanalyzeDigits> Cannot read the digits" << endl; |
a328fff9 | 333 | rc = 6; |
334 | return rc; | |
335 | } | |
336 | ||
337 | // Get the particle stack | |
338 | AliStack *kineStack = rl->Stack(); | |
339 | if (!kineStack) { | |
340 | cout << "<AliTRDanalyzeDigits> Cannot find the KINE stack" << endl; | |
341 | rc = 7; | |
fa148e6c | 342 | return rc; |
343 | } | |
344 | ||
345 | AliTRDmatrix *matrix; | |
346 | ||
347 | Int_t countDigits = 0; | |
348 | Int_t iSec = trd->GetSensSector(); | |
349 | Int_t iCha = trd->GetSensChamber(); | |
3551db50 | 350 | Int_t timeMax = cal->GetNumberOfTimeBins(); |
fa148e6c | 351 | |
352 | TProfile *hAmpTimeEl = new TProfile("hAmpTimeEl","Amplitude of the digits (electrons)" | |
353 | ,timeMax,-0.5,((Double_t) timeMax)-0.5); | |
354 | TProfile *hAmpTimePi = new TProfile("hAmpTimePi","Amplitude of the digits (pions)" | |
355 | ,timeMax,-0.5,((Double_t) timeMax)-0.5); | |
356 | ||
357 | // Loop over all planes | |
358 | for (Int_t iPla = 0; iPla < kNpla; iPla++) { | |
359 | ||
360 | Int_t iDet = geo->GetDetector(iPla,iCha,iSec); | |
3551db50 | 361 | Int_t rowMax = parCom->GetRowMax(iPla,iCha,iSec); |
362 | Int_t colMax = parCom->GetColMax(iPla); | |
fa148e6c | 363 | |
364 | if (iPla == 0) { | |
365 | matrix = new AliTRDmatrix(rowMax,colMax,timeMax,iSec,iCha,iPla); | |
366 | } | |
367 | ||
368 | // Loop through the detector pixel | |
369 | for (Int_t time = 0; time < timeMax; time++) { | |
370 | for (Int_t col = 0; col < colMax; col++) { | |
371 | for (Int_t row = 0; row < rowMax; row++) { | |
372 | ||
373 | AliTRDdigit *digit = digitsManager->GetDigit(row,col,time,iDet); | |
374 | Int_t amp = digit->GetAmp(); | |
375 | Int_t track0 = digitsManager->GetTrack(0,row,col,time,iDet); | |
376 | Int_t track1 = digitsManager->GetTrack(1,row,col,time,iDet); | |
377 | TParticle *particle = 0; | |
378 | if (track0 > -1) { | |
a328fff9 | 379 | particle = (TParticle *) kineStack->Particle(track0); |
fa148e6c | 380 | } |
381 | ||
382 | if (amp > 0) { | |
383 | countDigits++; | |
384 | if (iPla == 0) { | |
385 | matrix->SetSignal(row,col,time,amp); | |
386 | } | |
387 | } | |
388 | ||
389 | // Total spectrum | |
390 | hAmpAll->Fill(amp); | |
391 | ||
392 | // Noise spectrum | |
393 | if (track0 < 0) { | |
394 | hAmpNoise->Fill(amp); | |
395 | } | |
396 | ||
397 | // Electron digit | |
398 | if ((particle) && (particle->GetPdgCode() == 11) && (track1 < 0)) { | |
399 | hAmpEl->Fill(amp); | |
400 | hAmpTimeEl->Fill(time,amp); | |
401 | } | |
402 | ||
403 | // Pion digit | |
404 | if ((particle) && (particle->GetPdgCode() == -211) && (track1 < 0)) { | |
405 | hAmpPi->Fill(amp); | |
406 | hAmpTimePi->Fill(time,amp); | |
407 | } | |
408 | ||
409 | delete digit; | |
410 | ||
411 | } | |
412 | } | |
413 | } | |
414 | ||
415 | } | |
416 | ||
417 | cout << "<AliTRDanalyzeDigits> Found " << countDigits << " digits in total" << endl; | |
418 | ||
419 | // Display the detector matrix | |
420 | matrix->Draw(); | |
421 | ||
422 | TCanvas *cDigits = new TCanvas("cDigits","AliTRDanalyzeDigits",50,50,600,800); | |
423 | cDigits->Divide(2,3); | |
424 | cDigits->cd(1); | |
425 | gPad->SetLogy(); | |
426 | hAmpAll->SetXTitle("Amplitude (ADC-channels)"); | |
427 | hAmpAll->SetYTitle("Entries"); | |
428 | hAmpAll->Draw(); | |
429 | cDigits->cd(2); | |
430 | gPad->SetLogy(); | |
431 | hAmpNoise->SetXTitle("Amplitude (ADC-channels)"); | |
432 | hAmpNoise->SetYTitle("Entries"); | |
433 | hAmpNoise->Draw(); | |
434 | cDigits->cd(3); | |
435 | gPad->SetLogy(); | |
436 | hAmpEl->SetXTitle("Amplitude (ADC-channels)"); | |
437 | hAmpEl->SetYTitle("Entries"); | |
438 | hAmpEl->Draw(); | |
439 | cDigits->cd(4); | |
440 | gPad->SetLogy(); | |
441 | hAmpPi->SetXTitle("Amplitude (ADC-channels)"); | |
442 | hAmpPi->SetYTitle("Entries"); | |
443 | hAmpPi->Draw(); | |
444 | cDigits->cd(5); | |
445 | hAmpTimeEl->SetXTitle("Timebin number"); | |
446 | hAmpTimeEl->SetYTitle("Mean amplitude"); | |
447 | hAmpTimeEl->Draw("HIST"); | |
448 | cDigits->cd(6); | |
449 | hAmpTimePi->SetXTitle("Timebin number"); | |
450 | hAmpTimePi->SetYTitle("Mean amplitude"); | |
451 | hAmpTimePi->Draw("HIST"); | |
452 | ||
fa148e6c | 453 | return rc; |
454 | ||
455 | } | |
456 | ||
457 | //_____________________________________________________________________________ | |
458 | Int_t AliTRDcreateCluster() | |
459 | { | |
460 | // | |
461 | // Creates the cluster from the digits | |
462 | // | |
463 | ||
464 | Int_t rc = 0; | |
465 | ||
466 | // Create the clusterizer | |
467 | AliTRDclusterizerV1 *clusterizer = new AliTRDclusterizerV1("TRDclusterizer" | |
468 | ,"TRD clusterizer class"); | |
469 | clusterizer->SetVerbose(1); | |
470 | ||
471 | // Define the parameter object | |
472 | // If no external parameter object is defined, | |
473 | // default parameter will be used | |
474 | AliTRDparameter *parameter = new AliTRDparameter("TRDparameter" | |
475 | ,"TRD parameter class"); | |
476 | parameter->SetClusMaxThresh(0); | |
477 | parameter->SetClusSigThresh(0); | |
478 | clusterizer->SetParameter(parameter); | |
479 | ||
480 | // Open the file | |
481 | if (!(clusterizer->Open("TRD_test.root",0))) { | |
482 | rc = 1; | |
483 | return rc; | |
484 | } | |
485 | ||
486 | // Load the digits | |
487 | if (!(clusterizer->ReadDigits())) { | |
488 | rc = 2; | |
489 | return rc; | |
490 | } | |
491 | ||
492 | // Find the cluster | |
493 | if (!(clusterizer->MakeClusters())) { | |
494 | rc = 3; | |
495 | return rc; | |
496 | } | |
497 | ||
498 | // Write the cluster tree into the file | |
499 | if (!(clusterizer->WriteClusters(-1))) { | |
500 | rc = 4; | |
501 | return rc; | |
502 | } | |
503 | ||
504 | // Save the clusterizer class in the file | |
505 | if (!(clusterizer->Write())) { | |
506 | rc = 5; | |
507 | return rc; | |
508 | } | |
509 | ||
510 | return rc; | |
511 | ||
512 | } | |
513 | ||
514 | //_____________________________________________________________________________ | |
515 | Int_t AliTRDanalyzeCluster() | |
516 | { | |
517 | // | |
518 | // Analyzes the cluster | |
519 | // | |
520 | ||
521 | Int_t rc = 0; | |
522 | ||
523 | if (!gAlice) { | |
524 | cout << "<AliTRDanalyzeCluster> No AliRun object found" << endl; | |
525 | rc = 1; | |
526 | return rc; | |
527 | } | |
528 | gAlice->GetEvent(0); | |
529 | ||
530 | // Get the pointer to the TRD detector | |
531 | AliTRD *trd = (AliTRD *) gAlice->GetDetector("TRD"); | |
532 | if (!trd) { | |
533 | cout << "<AliTRDanalyzeCluster> No TRD detector found" << endl; | |
534 | rc = 2; | |
535 | return rc; | |
536 | } | |
537 | ||
538 | // Define the histograms | |
539 | TH1F *hClusAll = new TH1F("hClusAll" ,"Amplitude of the cluster (all)" | |
540 | ,501,-0.5,500.5); | |
541 | TH1F *hClusNoise = new TH1F("hClusNoise","Amplitude of the cluster (noise)" | |
963b25ce | 542 | , 11,-0.5, 10.5); |
fa148e6c | 543 | TH1F *hClusEl = new TH1F("hClusEl" ,"Amplitude of the cluster (electron)" |
544 | ,501,-0.5,500.5); | |
545 | TH1F *hClusPi = new TH1F("hClusPi" ,"Amplitude of the cluster (pion)" | |
546 | ,501,-0.5,500.5); | |
547 | ||
548 | // Get the pointer to the geometry object | |
549 | AliTRDgeometry *geo; | |
550 | if (trd) { | |
551 | geo = trd->GetGeometry(); | |
552 | } | |
553 | else { | |
554 | cout << "<AliTRDanalyzeCluster> No TRD geometry found" << endl; | |
555 | rc = 3; | |
556 | return rc; | |
557 | } | |
558 | ||
559 | // Get the pointer to the hit-tree | |
560 | TFile *file = (TFile *) gROOT->GetListOfFiles()->FindObject("TRD_test.root"); | |
561 | TTree *clusterTree = (TTree *) file->Get("TreeR0_TRD"); | |
562 | if (!(clusterTree)) { | |
563 | cout << "<AliTRDanalyzeCluster> No tree with clusters found" << endl; | |
564 | rc = 4; | |
565 | return rc; | |
566 | } | |
567 | ||
568 | // Get the pointer to the hit container | |
569 | TObjArray *clusterArray = trd->RecPoints(); | |
570 | if (!(clusterArray)) { | |
571 | cout << "<AliTRDanalyzeCluster> No clusterArray found" << endl; | |
572 | rc = 5; | |
573 | return rc; | |
574 | } | |
575 | ||
576 | // Set the branch address | |
577 | clusterTree->GetBranch("TRDcluster")->SetAddress(&clusterArray); | |
578 | Int_t nEntries = clusterTree->GetEntries(); | |
579 | cout << "<AliTRDanalyzeCluster> Number of entries in the cluster tree = " | |
580 | << nEntries | |
581 | << endl; | |
582 | ||
583 | Int_t countCluster = 0; | |
584 | Int_t countOverlap = 0; | |
585 | ||
586 | // Loop through all entries in the tree | |
587 | Int_t nbytes; | |
588 | for (Int_t iEntry = 0; iEntry < nEntries; iEntry++) { | |
589 | ||
590 | // Import the tree | |
591 | nbytes += clusterTree->GetEvent(iEntry); | |
592 | ||
593 | // Get the number of points in the detector | |
594 | Int_t nCluster = clusterArray->GetEntriesFast(); | |
595 | ||
596 | // Loop through all TRD digits | |
597 | for (Int_t iCluster = 0; iCluster < nCluster; iCluster++) { | |
598 | ||
599 | // Get the information for this digit | |
600 | AliTRDcluster *cluster = (AliTRDcluster *) clusterArray->UncheckedAt(iCluster); | |
601 | Int_t detector = cluster->GetDetector(); | |
602 | Int_t sector = geo->GetSector(detector); | |
603 | Int_t plane = geo->GetPlane(detector); | |
604 | Int_t chamber = geo->GetChamber(detector); | |
605 | Float_t energy = cluster->GetQ(); | |
f6ff4cec | 606 | Int_t track0 = cluster->GetLabel(0); |
607 | Int_t track1 = cluster->GetLabel(1); | |
608 | Int_t track2 = cluster->GetLabel(2); | |
fa148e6c | 609 | TParticle *particle = 0; |
610 | if (track0 > -1) { | |
611 | particle = gAlice->Particle(track0); | |
612 | } | |
613 | ||
614 | countCluster++; | |
615 | if (!cluster->Isolated()) countOverlap++; | |
616 | ||
617 | // Total spectrum | |
618 | hClusAll->Fill(energy); | |
619 | ||
620 | if (cluster->Isolated()) { | |
621 | ||
622 | // Noise spectrum | |
623 | if (track0 < 0) { | |
624 | hClusNoise->Fill(energy); | |
625 | } | |
626 | ||
627 | // Electron cluster | |
628 | if ((particle) && (particle->GetPdgCode() == 11) && (track1 < 0)) { | |
629 | hClusEl->Fill(energy); | |
630 | } | |
631 | ||
632 | // Pion cluster | |
633 | if ((particle) && (particle->GetPdgCode() == -211) && (track1 < 0)) { | |
634 | hClusPi->Fill(energy); | |
635 | } | |
636 | ||
637 | } | |
638 | ||
639 | } | |
640 | ||
641 | } | |
642 | ||
643 | cout << "<AliTRDanalyzeCluster> Found " << countCluster << " cluster in total" << endl; | |
644 | cout << "<AliTRDanalyzeCluster> Found " << countOverlap << " overlapping cluster" << endl; | |
645 | cout << endl; | |
646 | ||
647 | TCanvas *cCluster = new TCanvas("cCluster","AliTRDanalyzeCluster",50,50,600,600); | |
648 | cCluster->Divide(2,2); | |
649 | ||
650 | TF1 *fun; | |
651 | cCluster->cd(1); | |
652 | gPad->SetLogy(); | |
653 | hClusAll->Fit("landau","0"); | |
654 | fun = (TF1 *) hClusAll->GetListOfFunctions()->First(); | |
655 | Float_t meanAll = fun->GetParameter(1); | |
656 | hClusAll->Draw(); | |
657 | fun->SetLineColor(2); | |
658 | fun->Draw("SAME"); | |
659 | ||
963b25ce | 660 | cCluster->cd(2); |
661 | gPad->SetLogy(); | |
fa148e6c | 662 | Float_t meanNoise = hClusNoise->GetMean(); |
963b25ce | 663 | hClusNoise->Draw(); |
fa148e6c | 664 | |
665 | cCluster->cd(3); | |
666 | gPad->SetLogy(); | |
667 | hClusEl->Fit("landau","0"); | |
668 | fun = (TF1 *) hClusEl->GetListOfFunctions()->First(); | |
669 | fun->SetLineColor(2); | |
670 | Float_t meanEl = fun->GetParameter(1); | |
671 | hClusEl->Draw(); | |
672 | fun->Draw("SAME"); | |
673 | ||
674 | cCluster->cd(4); | |
675 | gPad->SetLogy(); | |
676 | hClusPi->Fit("landau","0"); | |
677 | fun = (TF1 *) hClusPi->GetListOfFunctions()->First(); | |
678 | fun->SetLineColor(2); | |
679 | Float_t meanPi = fun->GetParameter(1); | |
680 | hClusPi->Draw(); | |
681 | fun->Draw("SAME"); | |
682 | ||
683 | cout << endl; | |
684 | cout << "##################################################################" << endl; | |
685 | cout << " Mean all = " << meanAll << endl; | |
686 | cout << " Mean noise = " << meanNoise << endl; | |
687 | cout << " Mean electrons = " << meanEl << endl; | |
688 | cout << " Mean pions = " << meanPi << endl; | |
689 | cout << "##################################################################" << endl; | |
690 | cout << endl; | |
691 | ||
692 | return rc; | |
693 | ||
694 | } | |
8e64dd77 | 695 | |
696 | //_____________________________________________________________________________ | |
697 | Int_t AliTRDcreateTracks() | |
698 | { | |
699 | // | |
700 | // Creates the tracks | |
701 | // | |
702 | ||
703 | Int_t rc = 0; | |
704 | ||
705 | // Create the tracker | |
706 | AliTRDtracker *tracker = new AliTRDtracker("TRDtracker","TRD tracker"); | |
707 | ||
708 | // Read in the kine tree and the cluster | |
709 | tracker->GetEvent("TRD_test.root","TRD_test.root"); | |
710 | ||
711 | // Find the tracks | |
712 | TH1F *hs = new TH1F("hs","hs",100,0.0,1.0); | |
713 | TH1F *hd = new TH1F("hd","hd",100,0.0,1.0); | |
714 | tracker->Clusters2Tracks(hs,hd); | |
715 | ||
716 | // Store the tracks | |
717 | tracker->WriteTracks("TRD_test.root"); | |
718 | ||
719 | return rc; | |
720 | ||
721 | } |