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