b2b13ee2 |
1 | //-----------------------------------------------------// |
2 | // // |
3 | // Source File : PMDDigitization.cxx, Version 00 // |
4 | // // |
5 | // Date : September 20 2002 // |
6 | // // |
7 | //-----------------------------------------------------// |
8 | |
9 | #include <TBRIK.h> |
10 | #include <TNode.h> |
11 | #include <TTree.h> |
12 | #include <TGeometry.h> |
13 | #include <TObjArray.h> |
14 | #include <TClonesArray.h> |
15 | #include <TFile.h> |
16 | #include <TNtuple.h> |
17 | #include <TParticle.h> |
18 | |
19 | #include "AliRun.h" |
20 | #include "AliPMD.h" |
21 | #include "AliHit.h" |
22 | #include "AliDetector.h" |
23 | #include "AliRunLoader.h" |
24 | #include "AliLoader.h" |
25 | #include "AliConfig.h" |
26 | #include "AliMagF.h" |
27 | #include "AliRunDigitizer.h" |
28 | #include "AliHeader.h" |
29 | |
30 | #include "AliPMDcell.h" |
31 | #include "AliPMDsdigit.h" |
32 | #include "AliPMDdigit.h" |
33 | #include "AliPMDDigitizer.h" |
34 | #include "AliPMDClustering.h" |
35 | #include "AliPMDContainer.h" |
36 | #include "AliPMDrecpoint.h" |
37 | |
38 | |
39 | ClassImp(AliPMDDigitizer) |
40 | // |
41 | // Constructor |
42 | // |
43 | AliPMDDigitizer::AliPMDDigitizer() |
44 | { |
45 | if (!fSDigits) fSDigits = new TClonesArray("AliPMDsdigit", 1000); |
46 | fNsdigit = 0; |
47 | if (!fDigits) fDigits = new TClonesArray("AliPMDdigit", 1000); |
48 | fNdigit = 0; |
49 | |
50 | for (Int_t i = 0; i < fTotSM; i++) |
51 | { |
52 | for (Int_t j = 0; j < fNCell; j++) |
53 | { |
54 | for (Int_t k = 0; k < fNCell; k++) |
55 | { |
56 | fCPV[i][j][k] = 0.; |
57 | fPMD[i][j][k] = 0.; |
58 | } |
59 | } |
60 | } |
61 | |
62 | if (!fCell) fCell = new TObjArray(); |
63 | |
64 | fZPos = 361.5; // in units of cm, This is the default position of PMD |
65 | |
66 | } |
67 | AliPMDDigitizer::~AliPMDDigitizer() |
68 | { |
69 | delete fSDigits; |
70 | delete fDigits; |
71 | delete fCell; |
72 | } |
73 | // |
74 | // Member functions |
75 | // |
76 | void AliPMDDigitizer::OpengAliceFile(Char_t *file, Option_t *option) |
77 | { |
78 | |
79 | fRunLoader = AliRunLoader::Open(file,AliConfig::fgkDefaultEventFolderName, |
80 | "UPDATE"); |
81 | |
82 | if (!fRunLoader) |
83 | { |
84 | Error("Open","Can not open session for file %s.",file); |
85 | } |
86 | |
87 | fRunLoader->LoadgAlice(); |
88 | fRunLoader->LoadHeader(); |
89 | fRunLoader->LoadKinematics(); |
90 | |
91 | gAlice = fRunLoader->GetAliRun(); |
92 | |
93 | if (gAlice) |
94 | { |
95 | printf("<AliPMDdigitizer::Open> "); |
96 | printf("AliRun object found on file.\n"); |
97 | } |
98 | else |
99 | { |
100 | printf("<AliPMDdigitizer::Open> "); |
101 | printf("Could not find AliRun object.\n"); |
102 | } |
103 | |
104 | PMD = (AliPMD*)gAlice->GetDetector("PMD"); |
105 | pmdloader = fRunLoader->GetLoader("PMDLoader"); |
106 | if (pmdloader == 0x0) |
107 | { |
108 | cerr<<"Hits2Digits : Can not find PMD or PMDLoader\n"; |
109 | } |
110 | |
111 | const char *cHS = strstr(option,"HS"); |
112 | const char *cHD = strstr(option,"HD"); |
113 | const char *cSD = strstr(option,"SD"); |
114 | |
115 | if (cHS) |
116 | { |
117 | pmdloader->LoadHits("READ"); |
118 | pmdloader->LoadSDigits("recreate"); |
119 | } |
120 | else if (cHD) |
121 | { |
122 | pmdloader->LoadHits("READ"); |
123 | pmdloader->LoadDigits("recreate"); |
124 | } |
125 | else if (cSD) |
126 | { |
127 | pmdloader->LoadSDigits("READ"); |
128 | pmdloader->LoadDigits("recreate"); |
129 | } |
130 | |
131 | } |
132 | void AliPMDDigitizer::Hits2SDigits(Int_t ievt) |
133 | { |
134 | cout << " -------- Beginning of Hits2SDigits ----------- " << endl; |
135 | |
136 | Int_t kPi0 = 111; |
137 | Int_t kGamma = 22; |
138 | Int_t npmd; |
139 | Int_t trackno; |
140 | Int_t cellnumber; |
141 | Int_t smnumber1; |
142 | Int_t smnumber; |
143 | Int_t trackpid; |
144 | Int_t mtrackno; |
145 | Int_t mtrackpid; |
146 | |
147 | Float_t xPos, yPos, zPos; |
148 | Int_t xpad, ypad; |
149 | Float_t edep; |
150 | Float_t vx = -999.0, vy = -999.0, vz = -999.0; |
151 | |
152 | |
153 | ResetSDigit(); |
154 | |
155 | printf("Event Number = %d \n",ievt); |
156 | Int_t nparticles = fRunLoader->GetHeader()->GetNtrack(); |
157 | printf("Number of Particles = %d \n", nparticles); |
158 | fRunLoader->GetEvent(ievt); |
159 | Particles = gAlice->Particles(); |
160 | // ------------------------------------------------------- // |
161 | // Pointer to specific detector hits. |
162 | // Get pointers to Alice detectors and Hits containers |
163 | |
164 | treeH = pmdloader->TreeH(); |
165 | |
166 | Int_t ntracks = (Int_t) treeH->GetEntries(); |
167 | printf("Number of Tracks in the TreeH = %d \n", ntracks); |
168 | |
169 | treeS = pmdloader->TreeS(); |
170 | if (treeS == 0x0) |
171 | { |
172 | pmdloader->MakeTree("S"); |
173 | treeS = pmdloader->TreeS(); |
174 | } |
175 | Int_t bufsize = 16000; |
176 | treeS->Branch("PMDSDigit", &fSDigits, bufsize); |
177 | |
178 | if (PMD) PMDhits = PMD->Hits(); |
179 | |
180 | // Start loop on tracks in the hits containers |
181 | |
182 | |
183 | for (Int_t track=0; track<ntracks;track++) |
184 | { |
185 | gAlice->ResetHits(); |
186 | treeH->GetEvent(track); |
187 | |
188 | if (PMD) |
189 | { |
190 | npmd = PMDhits->GetEntriesFast(); |
191 | for (int ipmd = 0; ipmd < npmd; ipmd++) |
192 | { |
193 | pmdHit = (AliPMDhit*) PMDhits->UncheckedAt(ipmd); |
194 | trackno = pmdHit->GetTrack(); |
195 | |
196 | // get kinematics of the particles |
197 | |
198 | particle = gAlice->Particle(trackno); |
199 | trackpid = particle->GetPdgCode(); |
200 | |
201 | Int_t igatr = -999; |
202 | Int_t ichtr = -999; |
203 | |
204 | Int_t imo; |
205 | Int_t igen = 0; |
206 | Int_t id_mo = -999; |
207 | TParticle* mparticle = particle; |
208 | |
209 | while((imo = mparticle->GetFirstMother()) >= 0) |
210 | { |
211 | igen++; |
212 | mparticle = gAlice->Particle(imo); |
213 | id_mo = mparticle->GetPdgCode(); |
214 | |
215 | vx = mparticle->Vx(); |
216 | vy = mparticle->Vy(); |
217 | vz = mparticle->Vz(); |
218 | |
219 | //printf("==> Mother ID %5d %5d %5d Vertex: %13.3f %13.3f %13.3f\n", igen, imo, id_mo, vx, vy, vz); |
220 | //fprintf(ftest1,"==> Mother ID %5d %5d %5d Vertex: %13.3f %13.3f %13.3f\n", igen, imo, id_mo, vx, vy, vz); |
221 | |
222 | if (id_mo == kGamma && vx == 0. && vy == 0. && vz == 0.) |
223 | { |
224 | igatr = imo; |
225 | } |
226 | ichtr = imo; |
227 | } |
228 | |
229 | if (id_mo == kPi0 && vx == 0. && vy == 0. && vz == 0.) |
230 | { |
231 | mtrackno = igatr; |
232 | mtrackpid = kGamma; |
233 | } |
234 | else |
235 | { |
236 | mtrackno = ichtr; |
237 | mtrackpid = id_mo; |
238 | } |
239 | |
240 | xPos = pmdHit->X(); |
241 | yPos = pmdHit->Y(); |
242 | zPos = pmdHit->Z(); |
243 | cellnumber = pmdHit->fVolume[1]; |
244 | smnumber1 = pmdHit->fVolume[4]; |
245 | edep = pmdHit->fEnergy; |
246 | |
247 | if (smnumber1 > 3 && smnumber1 <= 6) |
248 | { |
249 | Int_t ny = (cellnumber-1)/48 + 1; |
250 | Int_t nx = cellnumber - (ny-1)*48; |
251 | Int_t nx1 = nx + 24; |
252 | Int_t ny1 = ny + 48; |
253 | Int_t ncell1 = (ny1 - 1)* 72 + nx1; |
254 | cellnumber = ncell1; |
255 | } |
256 | smnumber = Convert2RealSMNumber(smnumber1); |
257 | ypad = (cellnumber - 1)/fNCell + 1; |
258 | xpad = cellnumber - (ypad-1)*fNCell; |
259 | |
260 | //cout << "-zpos = " << -zPos << endl; |
261 | |
262 | if (-zPos < fZPos) |
263 | { |
264 | // CPV |
265 | fDetNo = 1; |
266 | } |
267 | else if (-zPos > fZPos) |
268 | { |
269 | // PMD |
270 | fDetNo = 0; |
271 | } |
272 | |
273 | if (fDetNo == 0) |
274 | { |
275 | fPMD[smnumber-1][xpad-1][ypad-1] += edep; |
276 | fPMDCounter[smnumber-1][xpad-1][ypad-1]++; |
277 | Int_t smn = smnumber - 1; |
278 | Int_t ixx = xpad - 1; |
279 | Int_t iyy = ypad - 1; |
280 | |
281 | pmdcell = new AliPMDcell(mtrackno,smn,ixx,iyy,edep); |
282 | |
283 | fCell->Add(pmdcell); |
284 | } |
285 | else if(fDetNo == 1) |
286 | { |
287 | fCPV[smnumber-1][xpad-1][ypad-1] += edep; |
288 | fCPVTrackNo[smnumber-1][xpad-1][ypad-1] = mtrackno; |
289 | } |
290 | } |
291 | } |
292 | } // Track Loop ended |
293 | |
294 | |
295 | TrackAssignment2Cell(); |
296 | |
297 | ResetCell(); |
298 | |
299 | Float_t deltaE = 0.; |
300 | Int_t detno = 0; |
301 | Int_t count_digit = 0; |
302 | Int_t trno = -1; |
303 | Int_t cellno = 0; |
304 | |
305 | for (Int_t idet = 0; idet < 2; idet++) |
306 | { |
307 | for (Int_t ism = 0; ism < fTotSM; ism++) |
308 | { |
309 | for (Int_t jrow = 0; jrow < fNCell; jrow++) |
310 | { |
311 | for (Int_t kcol = 0; kcol < fNCell; kcol++) |
312 | { |
313 | cellno = jrow + kcol*fNCell; |
314 | if (idet == 0) |
315 | { |
316 | deltaE = fPMD[ism][jrow][kcol]; |
317 | trno = fPMDTrackNo[ism][jrow][kcol]; |
318 | detno = 0; |
319 | } |
320 | else if (idet == 1) |
321 | { |
322 | deltaE = fCPV[ism][jrow][kcol]; |
323 | trno = fCPVTrackNo[ism][jrow][kcol]; |
324 | detno = 1; |
325 | } |
326 | if (deltaE > 0.) |
327 | { |
328 | count_digit++; |
329 | AddSDigit(trno,detno,ism,cellno,deltaE); |
330 | } |
331 | } |
332 | } |
333 | treeS->Fill(); |
334 | ResetSDigit(); |
335 | } |
336 | } |
337 | pmdloader->WriteSDigits("OVERWRITE"); |
338 | |
339 | ResetCellADC(); |
340 | |
341 | // cout << " -------- End of Hits2SDigit ----------- " << endl; |
342 | } |
343 | |
344 | void AliPMDDigitizer::Hits2Digits(Int_t ievt) |
345 | { |
346 | Int_t kPi0 = 111; |
347 | Int_t kGamma = 22; |
348 | Int_t npmd; |
349 | Int_t trackno; |
350 | Int_t cellnumber; |
351 | Int_t smnumber1; |
352 | Int_t smnumber; |
353 | Int_t trackpid; |
354 | Int_t mtrackno; |
355 | Int_t mtrackpid; |
356 | |
357 | Float_t xPos, yPos, zPos; |
358 | Int_t xpad, ypad; |
359 | Float_t edep; |
360 | Float_t vx = -999.0, vy = -999.0, vz = -999.0; |
361 | |
362 | |
363 | ResetDigit(); |
364 | |
365 | printf("Event Number = %d \n",ievt); |
366 | |
367 | Int_t nparticles = fRunLoader->GetHeader()->GetNtrack(); |
368 | printf("Number of Particles = %d \n", nparticles); |
369 | fRunLoader->GetEvent(ievt); |
370 | Particles = gAlice->Particles(); |
371 | // ------------------------------------------------------- // |
372 | // Pointer to specific detector hits. |
373 | // Get pointers to Alice detectors and Hits containers |
374 | |
375 | PMD = (AliPMD*)gAlice->GetDetector("PMD"); |
376 | pmdloader = fRunLoader->GetLoader("PMDLoader"); |
377 | |
378 | if (pmdloader == 0x0) |
379 | { |
380 | cerr<<"Hits2Digits method : Can not find PMD or PMDLoader\n"; |
381 | } |
382 | treeH = pmdloader->TreeH(); |
383 | Int_t ntracks = (Int_t) treeH->GetEntries(); |
384 | printf("Number of Tracks in the TreeH = %d \n", ntracks); |
385 | pmdloader->LoadDigits("recreate"); |
386 | treeD = pmdloader->TreeD(); |
387 | if (treeD == 0x0) |
388 | { |
389 | pmdloader->MakeTree("D"); |
390 | treeD = pmdloader->TreeD(); |
391 | } |
392 | Int_t bufsize = 16000; |
393 | treeD->Branch("PMDDigit", &fDigits, bufsize); |
394 | |
395 | if (PMD) PMDhits = PMD->Hits(); |
396 | |
397 | // Start loop on tracks in the hits containers |
398 | |
399 | for (Int_t track=0; track<ntracks;track++) |
400 | { |
401 | gAlice->ResetHits(); |
402 | treeH->GetEvent(track); |
403 | |
404 | if (PMD) |
405 | { |
406 | npmd = PMDhits->GetEntriesFast(); |
407 | for (int ipmd = 0; ipmd < npmd; ipmd++) |
408 | { |
409 | pmdHit = (AliPMDhit*) PMDhits->UncheckedAt(ipmd); |
410 | trackno = pmdHit->GetTrack(); |
411 | |
412 | // get kinematics of the particles |
413 | |
414 | particle = gAlice->Particle(trackno); |
415 | trackpid = particle->GetPdgCode(); |
416 | |
417 | Int_t igatr = -999; |
418 | Int_t ichtr = -999; |
419 | |
420 | Int_t imo; |
421 | Int_t igen = 0; |
422 | Int_t id_mo = -999; |
423 | TParticle* mparticle = particle; |
424 | |
425 | while((imo = mparticle->GetFirstMother()) >= 0) |
426 | { |
427 | igen++; |
428 | mparticle = gAlice->Particle(imo); |
429 | id_mo = mparticle->GetPdgCode(); |
430 | |
431 | vx = mparticle->Vx(); |
432 | vy = mparticle->Vy(); |
433 | vz = mparticle->Vz(); |
434 | |
435 | //printf("==> Mother ID %5d %5d %5d Vertex: %13.3f %13.3f %13.3f\n", igen, imo, id_mo, vx, vy, vz); |
436 | //fprintf(ftest1,"==> Mother ID %5d %5d %5d Vertex: %13.3f %13.3f %13.3f\n", igen, imo, id_mo, vx, vy, vz); |
437 | |
438 | if (id_mo == kGamma && vx == 0. && vy == 0. && vz == 0.) |
439 | { |
440 | igatr = imo; |
441 | } |
442 | ichtr = imo; |
443 | } |
444 | |
445 | if (id_mo == kPi0 && vx == 0. && vy == 0. && vz == 0.) |
446 | { |
447 | mtrackno = igatr; |
448 | mtrackpid = kGamma; |
449 | } |
450 | else |
451 | { |
452 | mtrackno = ichtr; |
453 | mtrackpid = id_mo; |
454 | } |
455 | |
456 | xPos = pmdHit->X(); |
457 | yPos = pmdHit->Y(); |
458 | zPos = pmdHit->Z(); |
459 | cellnumber = pmdHit->fVolume[1]; |
460 | smnumber1 = pmdHit->fVolume[4]; |
461 | edep = pmdHit->fEnergy; |
462 | |
463 | if (smnumber1 > 3 && smnumber1 <= 6) |
464 | { |
465 | Int_t ny = (cellnumber-1)/48 + 1; |
466 | Int_t nx = cellnumber - (ny-1)*48; |
467 | Int_t nx1 = nx + 24; |
468 | Int_t ny1 = ny + 48; |
469 | Int_t ncell1 = (ny1 - 1)* 72 + nx1; |
470 | cellnumber = ncell1; |
471 | } |
472 | |
473 | smnumber = Convert2RealSMNumber(smnumber1); |
474 | ypad = (cellnumber - 1)/fNCell + 1; |
475 | xpad = cellnumber - (ypad-1)*fNCell; |
476 | |
477 | //cout << "-zpos = " << -zPos << endl; |
478 | if (-zPos < fZPos) |
479 | { |
480 | // CPV |
481 | fDetNo = 1; |
482 | } |
483 | else if (-zPos > fZPos) |
484 | { |
485 | // PMD |
486 | fDetNo = 0; |
487 | } |
488 | |
489 | if(fDetNo == 1) |
490 | { |
491 | fCPV[smnumber-1][xpad-1][ypad-1] += edep; |
492 | } |
493 | else if (fDetNo == 0) |
494 | { |
495 | fPMD[smnumber-1][xpad-1][ypad-1] += edep; |
496 | fPMDCounter[smnumber-1][xpad-1][ypad-1]++; |
497 | } |
498 | } |
499 | } |
500 | } // Track Loop ended |
501 | |
502 | TrackAssignment2Cell(); |
503 | ResetCell(); |
504 | |
505 | Float_t deltaE = 0.; |
506 | Int_t detno = 0; |
507 | Int_t count_digit = 0; |
508 | Int_t trno = 1; |
509 | Int_t cellno; |
510 | |
511 | for (Int_t idet = 0; idet < 2; idet++) |
512 | { |
513 | for (Int_t ism = 0; ism < fTotSM; ism++) |
514 | { |
515 | for (Int_t jrow = 0; jrow < fNCell; jrow++) |
516 | { |
517 | for (Int_t kcol = 0; kcol < fNCell; kcol++) |
518 | { |
519 | cellno = jrow + kcol*fNCell; |
520 | if (idet == 0) |
521 | { |
522 | deltaE = fPMD[ism][jrow][kcol]; |
523 | detno = 0; |
524 | } |
525 | else if (idet == 1) |
526 | { |
527 | deltaE = fCPV[ism][jrow][kcol]; |
528 | detno = 1; |
529 | } |
530 | if (deltaE > 0.) |
531 | { |
532 | count_digit++; |
533 | AddDigit(trno,detno,ism,cellno,deltaE); |
534 | } |
535 | } // column loop |
536 | } // row loop |
537 | } // supermodule loop |
538 | treeD->Fill(); |
539 | ResetDigit(); |
540 | } // detector loop |
541 | |
542 | pmdloader->WriteDigits("OVERWRITE"); |
543 | |
544 | ResetCellADC(); |
545 | |
546 | // cout << " -------- End of Hits2Digit ----------- " << endl; |
547 | } |
548 | |
549 | |
550 | void AliPMDDigitizer::SDigits2Digits(Int_t ievt) |
551 | { |
552 | // cout << " -------- Beginning of SDigits2Digit ----------- " << endl; |
553 | fRunLoader->GetEvent(ievt); |
554 | |
555 | treeS = pmdloader->TreeS(); |
556 | AliPMDsdigit *pmdsdigit; |
557 | TBranch *branch = treeS->GetBranch("PMDSDigit"); |
558 | branch->SetAddress(&fSDigits); |
559 | |
560 | treeD = pmdloader->TreeD(); |
561 | if (treeD == 0x0) |
562 | { |
563 | pmdloader->MakeTree("D"); |
564 | treeD = pmdloader->TreeD(); |
565 | } |
566 | Int_t bufsize = 16000; |
567 | treeD->Branch("PMDDigit", &fDigits, bufsize); |
568 | |
569 | Int_t trno, det, smn; |
570 | Int_t cellno; |
571 | Float_t edep, adc; |
572 | |
573 | Int_t nmodules = (Int_t) treeS->GetEntries(); |
574 | |
575 | for (Int_t imodule = 0; imodule < nmodules; imodule++) |
576 | { |
577 | treeS->GetEntry(imodule); |
578 | Int_t nentries = fSDigits->GetLast(); |
579 | //cout << " nentries = " << nentries << endl; |
580 | for (Int_t ient = 0; ient < nentries+1; ient++) |
581 | { |
582 | pmdsdigit = (AliPMDsdigit*)fSDigits->UncheckedAt(ient); |
583 | trno = pmdsdigit->GetTrackNumber(); |
584 | det = pmdsdigit->GetDetector(); |
585 | smn = pmdsdigit->GetSMNumber(); |
586 | cellno = pmdsdigit->GetCellNumber(); |
587 | edep = pmdsdigit->GetCellEdep(); |
588 | |
589 | MeV2ADC(edep,adc); |
590 | |
591 | AddDigit(trno,det,smn,cellno,adc); |
592 | } |
593 | treeD->Fill(); |
594 | ResetDigit(); |
595 | } |
596 | pmdloader->WriteDigits("OVERWRITE"); |
597 | // cout << " -------- End of SDigits2Digit ----------- " << endl; |
598 | } |
599 | |
600 | void AliPMDDigitizer::TrackAssignment2Cell() |
601 | { |
602 | // To be checked again |
603 | // |
604 | // This blocks assign the cell id when there is |
605 | // multiple tracks in a cell according to the |
606 | // energy deposition |
607 | // |
608 | |
609 | Int_t i, j, k; |
610 | |
611 | Float_t *frac_edp; |
612 | Float_t *tr_edp; |
613 | Int_t *status; |
614 | Int_t ****PMDTrack; |
615 | Float_t ****PMDEdep; |
616 | |
617 | PMDTrack = new Int_t ***[27]; |
618 | PMDEdep = new Float_t ***[27]; |
619 | for (i=0; i<fTotSM; i++) |
620 | { |
621 | PMDTrack[i] = new Int_t **[72]; |
622 | PMDEdep[i] = new Float_t **[72]; |
623 | } |
624 | |
625 | for (i = 0; i < fTotSM; i++) |
626 | { |
627 | for (j = 0; j < fNCell; j++) |
628 | { |
629 | PMDTrack[i][j] = new Int_t *[72]; |
630 | PMDEdep[i][j] = new Float_t *[72]; |
631 | } |
632 | } |
633 | |
634 | for (i = 0; i < fTotSM; i++) |
635 | { |
636 | for (j = 0; j < fNCell; j++) |
637 | { |
638 | for (k = 0; k < fNCell; k++) |
639 | { |
640 | Int_t nn = fPMDCounter[i][j][k]; |
641 | if(nn > 0) |
642 | { |
643 | PMDTrack[i][j][k] = new Int_t[nn]; |
644 | PMDEdep[i][j][k] = new Float_t[nn]; |
645 | } |
646 | else |
647 | { |
648 | nn = 1; |
649 | PMDTrack[i][j][k] = new Int_t[nn]; |
650 | PMDEdep[i][j][k] = new Float_t[nn]; |
651 | } |
652 | fPMDCounter[i][j][k] = 0; |
653 | } |
654 | } |
655 | } |
656 | |
657 | |
658 | Int_t nentries = fCell->GetEntries(); |
659 | |
660 | Int_t mtrackno, ism, ixp, iyp; |
661 | Float_t edep; |
662 | |
663 | for (i = 0; i < nentries; i++) |
664 | { |
665 | pmdcell = (AliPMDcell*)fCell->UncheckedAt(i); |
666 | |
667 | mtrackno = pmdcell->GetTrackNumber(); |
668 | ism = pmdcell->GetSMNumber(); |
669 | ixp = pmdcell->GetX(); |
670 | iyp = pmdcell->GetY(); |
671 | edep = pmdcell->GetEdep(); |
672 | |
673 | Int_t nn = fPMDCounter[ism][ixp][iyp]; |
674 | |
675 | // cout << " nn = " << nn << endl; |
676 | |
677 | PMDTrack[ism][ixp][iyp][nn] = (Int_t) mtrackno; |
678 | PMDEdep[ism][ixp][iyp][nn] = edep; |
679 | fPMDCounter[ism][ixp][iyp]++; |
680 | } |
681 | |
682 | Int_t iz, il; |
683 | Int_t im, ix, iy; |
684 | Int_t nn; |
685 | |
686 | for (im=0; im<27; im++) |
687 | { |
688 | for (ix=0; ix<72; ix++) |
689 | { |
690 | for (iy=0; iy<72; iy++) |
691 | { |
692 | nn = fPMDCounter[im][ix][iy]; |
693 | if (nn > 1) |
694 | { |
695 | // This block handles if a cell is fired |
696 | // many times by many tracks |
697 | |
698 | status = new Int_t[nn]; |
699 | for (iz = 0; iz < nn; iz++) |
700 | { |
701 | status[iz] = PMDTrack[im][ix][iy][iz]; |
702 | } |
703 | sort(status,status+nn); |
704 | Int_t track_old = -99999; |
705 | Int_t track, tr_count = 0; |
706 | for (iz = 0; iz < nn; iz++) |
707 | { |
708 | track = status[iz]; |
709 | if (track_old != track) |
710 | { |
711 | tr_count++; |
712 | vjunkTRN.push_back(track); |
713 | } |
714 | track_old = track; |
715 | } |
716 | delete status; |
717 | Float_t tot_edp = 0.; |
718 | tr_edp = new Float_t[tr_count]; |
719 | frac_edp = new Float_t[tr_count]; |
720 | for (il = 0; il < tr_count; il++) |
721 | { |
722 | tr_edp[il] = 0.; |
723 | track = vjunkTRN[il]; |
724 | for (iz = 0; iz < nn; iz++) |
725 | { |
726 | if (track == PMDTrack[im][ix][iy][iz]) |
727 | { |
728 | tr_edp[il] += PMDEdep[im][ix][iy][iz]; |
729 | } |
730 | } |
731 | tot_edp += tr_edp[il]; |
732 | } |
733 | vjunkTRN.clear(); |
734 | Int_t il_old = 0; |
735 | Float_t frac_old = 0.; |
736 | |
737 | for (il = 0; il < tr_count; il++) |
738 | { |
739 | frac_edp[il] = tr_edp[il]/tot_edp; |
740 | if (frac_old < frac_edp[il]) |
741 | { |
742 | frac_old = frac_edp[il]; |
743 | il_old = il; |
744 | } |
745 | } |
746 | delete frac_edp; |
747 | delete tr_edp; |
748 | |
749 | fPMDTrackNo[im][ix][iy] = vjunkTRN[il_old]; |
750 | } |
751 | else if (nn == 1) |
752 | { |
753 | // This only handles if a cell is fired |
754 | // by only one track |
755 | |
756 | fPMDTrackNo[im][ix][iy] = PMDTrack[im][ix][iy][0]; |
757 | |
758 | } |
759 | else if (nn ==0) |
760 | { |
761 | // This is if no cell is fired |
762 | fPMDTrackNo[im][ix][iy] = -999; |
763 | } |
764 | } // end of iy |
765 | } // end of ix |
766 | } // end of im |
767 | |
768 | // Delete all the pointers |
769 | |
770 | for (i = 0; i < fTotSM; i++) |
771 | { |
772 | for (j = 0; j < fNCell; j++) |
773 | { |
774 | for (k = 0; k < fNCell; k++) |
775 | { |
776 | delete [] PMDTrack[i][j][k]; |
777 | delete [] PMDEdep[i][j][k]; |
778 | } |
779 | } |
780 | } |
781 | |
782 | for (i = 0; i < fTotSM; i++) |
783 | { |
784 | for (j = 0; j < fNCell; j++) |
785 | { |
786 | delete [] PMDTrack[i][j]; |
787 | delete [] PMDEdep[i][j]; |
788 | } |
789 | } |
790 | |
791 | for (i = 0; i < fTotSM; i++) |
792 | { |
793 | delete [] PMDTrack[i]; |
794 | delete [] PMDEdep[i]; |
795 | } |
796 | delete PMDTrack; |
797 | delete PMDEdep; |
798 | // |
799 | // End of the cell id assignment |
800 | // |
801 | } |
802 | |
803 | |
804 | void AliPMDDigitizer::MeV2ADC(Float_t mev, Float_t & adc) |
805 | { |
806 | // To be done |
807 | |
808 | adc = mev*1.; |
809 | } |
810 | void AliPMDDigitizer::AddSDigit(Int_t trnumber, Int_t det, Int_t smnumber, |
811 | Int_t cellnumber, Float_t adc) |
812 | { |
813 | TClonesArray &lsdigits = *fSDigits; |
814 | AliPMDsdigit *newcell; |
815 | newcell = new AliPMDsdigit(trnumber,det,smnumber,cellnumber,adc); |
816 | new(lsdigits[fNsdigit++]) AliPMDsdigit(newcell); |
817 | delete newcell; |
818 | } |
819 | |
820 | void AliPMDDigitizer::AddDigit(Int_t trnumber, Int_t det, Int_t smnumber, |
821 | Int_t cellnumber, Float_t adc) |
822 | { |
823 | TClonesArray &ldigits = *fDigits; |
824 | AliPMDdigit *newcell; |
825 | newcell = new AliPMDdigit(trnumber,det,smnumber,cellnumber,adc); |
826 | new(ldigits[fNdigit++]) AliPMDdigit(newcell); |
827 | delete newcell; |
828 | } |
829 | |
830 | Int_t AliPMDDigitizer::Convert2RealSMNumber(Int_t smnumber1) |
831 | { |
832 | Int_t smnumber = -999; |
833 | |
834 | if (smnumber1==1) smnumber = 1; |
835 | if (smnumber1==2) smnumber = 10; |
836 | if (smnumber1==3) smnumber = 19; |
837 | if (smnumber1==4) smnumber = 1; |
838 | if (smnumber1==5) smnumber = 10; |
839 | if (smnumber1==6) smnumber = 19; |
840 | if (smnumber1==7) smnumber = 2; |
841 | if (smnumber1==8) smnumber = 3; |
842 | if (smnumber1==9) smnumber = 4; |
843 | if (smnumber1==10) smnumber = 5; |
844 | if (smnumber1==11) smnumber = 6; |
845 | if (smnumber1==12) smnumber = 7; |
846 | if (smnumber1==13) smnumber = 8; |
847 | if (smnumber1==14) smnumber = 9; |
848 | if (smnumber1==15) smnumber = 11; |
849 | if (smnumber1==16) smnumber = 12; |
850 | if (smnumber1==17) smnumber = 13; |
851 | if (smnumber1==18) smnumber = 14; |
852 | if (smnumber1==19) smnumber = 15; |
853 | if (smnumber1==20) smnumber = 16; |
854 | if (smnumber1==21) smnumber = 17; |
855 | if (smnumber1==22) smnumber = 18; |
856 | if (smnumber1==23) smnumber = 20; |
857 | if (smnumber1==24) smnumber = 21; |
858 | if (smnumber1==25) smnumber = 22; |
859 | if (smnumber1==26) smnumber = 23; |
860 | if (smnumber1==27) smnumber = 24; |
861 | if (smnumber1==28) smnumber = 25; |
862 | if (smnumber1==29) smnumber = 26; |
863 | if (smnumber1==30) smnumber = 27; |
864 | |
865 | return smnumber; |
866 | } |
867 | void AliPMDDigitizer::SetZPosition(Float_t zpos) |
868 | { |
869 | fZPos = zpos; |
870 | } |
871 | Float_t AliPMDDigitizer::GetZPosition() const |
872 | { |
873 | return fZPos; |
874 | } |
875 | |
876 | void AliPMDDigitizer::ResetCell() |
877 | { |
878 | fCell->Clear(); |
879 | for (Int_t i = 0; i < fTotSM; i++) |
880 | { |
881 | for (Int_t j = 0; j < fNCell; j++) |
882 | { |
883 | for (Int_t k = 0; k < fNCell; k++) |
884 | { |
885 | fPMDCounter[i][j][k] = 0; |
886 | } |
887 | } |
888 | } |
889 | } |
890 | void AliPMDDigitizer::ResetSDigit() |
891 | { |
892 | fNsdigit = 0; |
893 | if (fSDigits) fSDigits->Clear(); |
894 | } |
895 | void AliPMDDigitizer::ResetDigit() |
896 | { |
897 | fNdigit = 0; |
898 | if (fDigits) fDigits->Clear(); |
899 | } |
900 | |
901 | void AliPMDDigitizer::ResetCellADC() |
902 | { |
903 | for (Int_t i = 0; i < fTotSM; i++) |
904 | { |
905 | for (Int_t j = 0; j < fNCell; j++) |
906 | { |
907 | for (Int_t k = 0; k < fNCell; k++) |
908 | { |
909 | fCPV[i][j][k] = 0.; |
910 | fPMD[i][j][k] = 0.; |
911 | } |
912 | } |
913 | } |
914 | } |
915 | |
916 | void AliPMDDigitizer::UnLoad(Option_t *option) |
917 | { |
918 | const char *cS = strstr(option,"S"); |
919 | const char *cD = strstr(option,"D"); |
920 | |
921 | fRunLoader->UnloadgAlice(); |
922 | fRunLoader->UnloadHeader(); |
923 | fRunLoader->UnloadKinematics(); |
924 | |
925 | if (cS) |
926 | { |
927 | pmdloader->UnloadHits(); |
928 | } |
929 | if (cD) |
930 | { |
931 | pmdloader->UnloadHits(); |
932 | pmdloader->UnloadSDigits(); |
933 | } |
934 | } |