1 /***************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
15 //-----------------------------------------------------//
17 // Source File : PMDDigitizer.cxx, Version 00 //
19 // Date : September 20 2002 //
21 //-----------------------------------------------------//
23 #include <Riostream.h>
27 #include <TGeometry.h>
28 #include <TObjArray.h>
29 #include <TClonesArray.h>
32 #include <TParticle.h>
36 #include "AliPMDhit.h"
38 #include "AliDetector.h"
39 #include "AliRunLoader.h"
40 #include "AliLoader.h"
41 #include "AliConfig.h"
43 #include "AliRunDigitizer.h"
44 #include "AliDigitizer.h"
45 #include "AliHeader.h"
47 #include "AliPMDcell.h"
48 #include "AliPMDsdigit.h"
49 #include "AliPMDdigit.h"
50 #include "AliPMDDigitizer.h"
53 ClassImp(AliPMDDigitizer)
55 AliPMDDigitizer::AliPMDDigitizer() :
66 fSDigits(new TClonesArray("AliPMDsdigit", 1000)),
67 fDigits(new TClonesArray("AliPMDdigit", 1000)),
68 fCell(new TObjArray()),
74 fZPos(361.5)// in units of cm, This is the default position of PMD
76 // Default Constructor
78 for (Int_t i = 0; i < fgkTotUM; i++)
80 for (Int_t j = 0; j < fgkRow; j++)
82 for (Int_t k = 0; k < fgkCol; k++)
86 fPRECounter[i][j][k] = 0;
87 fPRETrackNo[i][j][k] = -1;
88 fCPVTrackNo[i][j][k] = -1;
94 //____________________________________________________________________________
95 AliPMDDigitizer::AliPMDDigitizer(AliRunDigitizer* manager)
96 :AliDigitizer(manager)
98 // ctor which should be used
100 cerr<<"AliPMDDigitizer::AliPMDDigitizer"
101 <<"(AliRunDigitizer* manager) was processed"<<endl;
103 //____________________________________________________________________________
104 AliPMDDigitizer::~AliPMDDigitizer()
106 // Default Destructor
127 //____________________________________________________________________________
128 void AliPMDDigitizer::OpengAliceFile(const char *file, Option_t *option)
130 // Loads galice.root file and corresponding header, kinematics
131 // hits and sdigits or digits depending on the option
133 fRunLoader = AliRunLoader::Open(file,AliConfig::fgkDefaultEventFolderName,
138 Error("Open","Can not open session for file %s.",file);
141 fRunLoader->LoadgAlice();
142 fRunLoader->LoadHeader();
143 fRunLoader->LoadKinematics();
145 gAlice = fRunLoader->GetAliRun();
149 printf("<AliPMDdigitizer::Open> ");
150 printf("AliRun object found on file.\n");
154 printf("<AliPMDdigitizer::Open> ");
155 printf("Could not find AliRun object.\n");
158 fPMD = (AliPMD*)gAlice->GetDetector("PMD");
159 fPMDLoader = fRunLoader->GetLoader("PMDLoader");
160 if (fPMDLoader == 0x0)
162 cerr<<"Hits2Digits : Can not find PMD or PMDLoader\n";
165 const char *cHS = strstr(option,"HS");
166 const char *cHD = strstr(option,"HD");
167 const char *cSD = strstr(option,"SD");
171 fPMDLoader->LoadHits("READ");
172 fPMDLoader->LoadSDigits("recreate");
176 fPMDLoader->LoadHits("READ");
177 fPMDLoader->LoadDigits("recreate");
181 fPMDLoader->LoadSDigits("READ");
182 fPMDLoader->LoadDigits("recreate");
185 //____________________________________________________________________________
186 void AliPMDDigitizer::Hits2SDigits(Int_t ievt)
188 // This reads the PMD Hits tree and assigns the right track number
189 // to a cell and stores in the summable digits tree
191 // cout << " -------- Beginning of Hits2SDigits ----------- " << endl;
193 const Int_t kPi0 = 111;
194 const Int_t kGamma = 22;
202 Float_t xPos, yPos, zPos;
203 Int_t xpad = -1, ypad = -1;
205 Float_t vx = -999.0, vy = -999.0, vz = -999.0;
210 printf("Event Number = %d \n",ievt);
211 Int_t nparticles = fRunLoader->GetHeader()->GetNtrack();
212 printf("Number of Particles = %d \n", nparticles);
213 fRunLoader->GetEvent(ievt);
214 // ------------------------------------------------------- //
215 // Pointer to specific detector hits.
216 // Get pointers to Alice detectors and Hits containers
218 fTreeH = fPMDLoader->TreeH();
220 Int_t ntracks = (Int_t) fTreeH->GetEntries();
221 printf("Number of Tracks in the TreeH = %d \n", ntracks);
223 fTreeS = fPMDLoader->TreeS();
226 fPMDLoader->MakeTree("S");
227 fTreeS = fPMDLoader->TreeS();
229 Int_t bufsize = 16000;
230 fTreeS->Branch("PMDSDigit", &fSDigits, bufsize);
232 if (fPMD) fHits = fPMD->Hits();
234 // Start loop on tracks in the hits containers
236 for (Int_t track=0; track<ntracks;track++)
239 fTreeH->GetEvent(track);
242 npmd = fHits->GetEntriesFast();
243 for (int ipmd = 0; ipmd < npmd; ipmd++)
245 fPMDHit = (AliPMDhit*) fHits->UncheckedAt(ipmd);
246 trackno = fPMDHit->GetTrack();
247 // get kinematics of the particles
249 fParticle = gAlice->GetMCApp()->Particle(trackno);
250 trackpid = fParticle->GetPdgCode();
259 TParticle* mparticle = fParticle;
260 Int_t tracknoOld=0, trackpidOld=0, statusOld = 0;
261 if (mparticle->GetFirstMother() == -1)
263 tracknoOld = trackno;
264 trackpidOld = trackpid;
268 while((imo = mparticle->GetFirstMother()) >= 0)
272 mparticle = gAlice->GetMCApp()->Particle(imo);
273 idmo = mparticle->GetPdgCode();
275 vx = mparticle->Vx();
276 vy = mparticle->Vy();
277 vz = mparticle->Vz();
279 //printf("==> Mother ID %5d %5d %5d Vertex: %13.3f %13.3f %13.3f\n", igen, imo, idmo, vx, vy, vz);
280 //fprintf(fpw1,"==> Mother ID %5d %5d %5d Vertex: %13.3f %13.3f %13.3f\n", igen, imo, idmo, vx, vy, vz);
281 if ((idmo == kGamma || idmo == -11 || idmo == 11) && vx == 0. && vy == 0. && vz == 0.)
289 if (idmo == kPi0 && vx == 0. && vy == 0. && vz == 0.)
298 if (idmo == kPi0 && vx == 0. && vy == 0. && vz == 0.)
310 mtrackno = tracknoOld;
311 mtrackpid = trackpidOld;
317 edep = fPMDHit->GetEnergy();
318 Int_t vol1 = fPMDHit->GetVolume(1); // Column
319 Int_t vol2 = fPMDHit->GetVolume(2); // Row
320 Int_t vol3 = fPMDHit->GetVolume(3); // UnitModule
321 Int_t vol6 = fPMDHit->GetVolume(6); // SuperModule
323 // -----------------------------------------//
324 // For Super Module 1 & 2 //
325 // nrow = 96, ncol = 48 //
326 // For Super Module 3 & 4 //
327 // nrow = 48, ncol = 96 //
328 // -----------------------------------------//
330 smnumber = (vol6-1)*6 + vol3;
332 if (vol6 == 1 || vol6 == 2)
337 else if (vol6 == 3 || vol6 == 4)
343 //cout << "zpos = " << zPos << " edep = " << edep << endl;
345 Float_t zposition = TMath::Abs(zPos);
346 if (zposition < fZPos)
351 else if (zposition > fZPos)
356 Int_t smn = smnumber - 1;
357 Int_t ixx = xpad - 1;
358 Int_t iyy = ypad - 1;
361 fPRE[smn][ixx][iyy] += edep;
362 fPRECounter[smn][ixx][iyy]++;
364 fPMDcell = new AliPMDcell(mtrackno,smn,ixx,iyy,edep);
366 fCell->Add(fPMDcell);
370 fCPV[smn][ixx][iyy] += edep;
371 fCPVTrackNo[smn][ixx][iyy] = mtrackno;
375 } // Track Loop ended
377 TrackAssignment2Cell();
384 for (Int_t idet = 0; idet < 2; idet++)
386 for (Int_t ism = 0; ism < fgkTotUM; ism++)
388 for (Int_t jrow = 0; jrow < fgkRow; jrow++)
390 for (Int_t kcol = 0; kcol < fgkCol; kcol++)
394 deltaE = fPRE[ism][jrow][kcol];
395 trno = fPRETrackNo[ism][jrow][kcol];
400 deltaE = fCPV[ism][jrow][kcol];
401 trno = fCPVTrackNo[ism][jrow][kcol];
406 AddSDigit(trno,detno,ism,jrow,kcol,deltaE);
414 fPMDLoader->WriteSDigits("OVERWRITE");
417 // cout << " -------- End of Hits2SDigit ----------- " << endl;
419 //____________________________________________________________________________
421 void AliPMDDigitizer::Hits2Digits(Int_t ievt)
423 // This reads the PMD Hits tree and assigns the right track number
424 // to a cell and stores in the digits tree
426 const Int_t kPi0 = 111;
427 const Int_t kGamma = 22;
435 Float_t xPos, yPos, zPos;
436 Int_t xpad = -1, ypad = -1;
438 Float_t vx = -999.0, vy = -999.0, vz = -999.0;
442 printf("Event Number = %d \n",ievt);
444 Int_t nparticles = fRunLoader->GetHeader()->GetNtrack();
445 printf("Number of Particles = %d \n", nparticles);
446 fRunLoader->GetEvent(ievt);
447 // ------------------------------------------------------- //
448 // Pointer to specific detector hits.
449 // Get pointers to Alice detectors and Hits containers
451 fPMD = (AliPMD*)gAlice->GetDetector("PMD");
452 fPMDLoader = fRunLoader->GetLoader("PMDLoader");
454 if (fPMDLoader == 0x0)
456 cerr<<"Hits2Digits method : Can not find PMD or PMDLoader\n";
458 fTreeH = fPMDLoader->TreeH();
459 Int_t ntracks = (Int_t) fTreeH->GetEntries();
460 printf("Number of Tracks in the TreeH = %d \n", ntracks);
461 fPMDLoader->LoadDigits("recreate");
462 fTreeD = fPMDLoader->TreeD();
465 fPMDLoader->MakeTree("D");
466 fTreeD = fPMDLoader->TreeD();
468 Int_t bufsize = 16000;
469 fTreeD->Branch("PMDDigit", &fDigits, bufsize);
471 if (fPMD) fHits = fPMD->Hits();
473 // Start loop on tracks in the hits containers
475 for (Int_t track=0; track<ntracks;track++)
478 fTreeH->GetEvent(track);
482 npmd = fHits->GetEntriesFast();
483 for (int ipmd = 0; ipmd < npmd; ipmd++)
485 fPMDHit = (AliPMDhit*) fHits->UncheckedAt(ipmd);
486 trackno = fPMDHit->GetTrack();
488 // get kinematics of the particles
490 fParticle = gAlice->GetMCApp()->Particle(trackno);
491 trackpid = fParticle->GetPdgCode();
500 TParticle* mparticle = fParticle;
501 Int_t tracknoOld=0, trackpidOld=0, statusOld = 0;
502 if (mparticle->GetFirstMother() == -1)
504 tracknoOld = trackno;
505 trackpidOld = trackpid;
510 while((imo = mparticle->GetFirstMother()) >= 0)
514 mparticle = gAlice->GetMCApp()->Particle(imo);
515 idmo = mparticle->GetPdgCode();
517 vx = mparticle->Vx();
518 vy = mparticle->Vy();
519 vz = mparticle->Vz();
521 //printf("==> Mother ID %5d %5d %5d Vertex: %13.3f %13.3f %13.3f\n", igen, imo, idmo, vx, vy, vz);
522 //fprintf(fpw1,"==> Mother ID %5d %5d %5d Vertex: %13.3f %13.3f %13.3f\n", igen, imo, idmo, vx, vy, vz);
523 if ((idmo == kGamma || idmo == -11 || idmo == 11) && vx == 0. && vy == 0. && vz == 0.)
531 if (idmo == kPi0 && vx == 0. && vy == 0. && vz == 0.)
540 if (idmo == kPi0 && vx == 0. && vy == 0. && vz == 0.)
552 mtrackno = tracknoOld;
553 mtrackpid = trackpidOld;
560 edep = fPMDHit->GetEnergy();
561 Int_t vol1 = fPMDHit->GetVolume(1); // Column
562 Int_t vol2 = fPMDHit->GetVolume(2); // Row
563 Int_t vol3 = fPMDHit->GetVolume(3); // UnitModule
564 Int_t vol6 = fPMDHit->GetVolume(6); // SuperModule
566 // -----------------------------------------//
567 // For Super Module 1 & 2 //
568 // nrow = 96, ncol = 48 //
569 // For Super Module 3 & 4 //
570 // nrow = 48, ncol = 96 //
571 // -----------------------------------------//
573 smnumber = (vol6-1)*6 + vol3;
575 if (vol6 == 1 || vol6 == 2)
580 else if (vol6 == 3 || vol6 == 4)
586 //cout << "-zpos = " << -zPos << endl;
588 Float_t zposition = TMath::Abs(zPos);
590 if (zposition < fZPos)
595 else if (zposition > fZPos)
601 Int_t smn = smnumber - 1;
602 Int_t ixx = xpad - 1;
603 Int_t iyy = ypad - 1;
606 fPRE[smn][ixx][iyy] += edep;
607 fPRECounter[smn][ixx][iyy]++;
609 fPMDcell = new AliPMDcell(mtrackno,smn,ixx,iyy,edep);
611 fCell->Add(fPMDcell);
615 fCPV[smn][ixx][iyy] += edep;
616 fCPVTrackNo[smn][ixx][iyy] = mtrackno;
620 } // Track Loop ended
622 TrackAssignment2Cell();
630 for (Int_t idet = 0; idet < 2; idet++)
632 for (Int_t ism = 0; ism < fgkTotUM; ism++)
634 for (Int_t jrow = 0; jrow < fgkRow; jrow++)
636 for (Int_t kcol = 0; kcol < fgkCol; kcol++)
640 deltaE = fPRE[ism][jrow][kcol];
641 trno = fPRETrackNo[ism][jrow][kcol];
646 deltaE = fCPV[ism][jrow][kcol];
647 trno = fCPVTrackNo[ism][jrow][kcol];
653 AddDigit(trno,detno,ism,jrow,kcol,adc);
657 } // supermodule loop
662 fPMDLoader->WriteDigits("OVERWRITE");
665 // cout << " -------- End of Hits2Digit ----------- " << endl;
667 //____________________________________________________________________________
670 void AliPMDDigitizer::SDigits2Digits(Int_t ievt)
672 // This reads the PMD sdigits tree and converts energy deposition
673 // in a cell to ADC and stores in the digits tree
675 // cout << " -------- Beginning of SDigits2Digit ----------- " << endl;
676 fRunLoader->GetEvent(ievt);
678 fTreeS = fPMDLoader->TreeS();
679 AliPMDsdigit *pmdsdigit;
680 TBranch *branch = fTreeS->GetBranch("PMDSDigit");
681 branch->SetAddress(&fSDigits);
683 fTreeD = fPMDLoader->TreeD();
686 fPMDLoader->MakeTree("D");
687 fTreeD = fPMDLoader->TreeD();
689 Int_t bufsize = 16000;
690 fTreeD->Branch("PMDDigit", &fDigits, bufsize);
692 Int_t trno, det, smn;
696 Int_t nmodules = (Int_t) fTreeS->GetEntries();
698 for (Int_t imodule = 0; imodule < nmodules; imodule++)
700 fTreeS->GetEntry(imodule);
701 Int_t nentries = fSDigits->GetLast();
702 //cout << " nentries = " << nentries << endl;
703 for (Int_t ient = 0; ient < nentries+1; ient++)
705 pmdsdigit = (AliPMDsdigit*)fSDigits->UncheckedAt(ient);
706 trno = pmdsdigit->GetTrackNumber();
707 det = pmdsdigit->GetDetector();
708 smn = pmdsdigit->GetSMNumber();
709 irow = pmdsdigit->GetRow();
710 icol = pmdsdigit->GetColumn();
711 edep = pmdsdigit->GetCellEdep();
714 AddDigit(trno,det,smn,irow,icol,adc);
719 fPMDLoader->WriteDigits("OVERWRITE");
720 // cout << " -------- End of SDigits2Digit ----------- " << endl;
722 //____________________________________________________________________________
723 void AliPMDDigitizer::Exec(Option_t *option)
725 // Does the event merging and digitization
728 const char *cdeb = strstr(option,"deb");
731 cout << "**************** PMD Exec *************** " << endl;
734 fDigits = new TClonesArray("AliPMDdigit", 1000);
736 Int_t ninputs = fManager->GetNinputs();
739 cout << " Number of files = " << ninputs << endl;
743 for (Int_t i = 0; i < ninputs; i++)
745 Int_t troffset = fManager->GetMask(i);
746 MergeSDigits(i, troffset);
749 fRunLoader = AliRunLoader::GetRunLoader(fManager->GetOutputFolderName());
750 fPMD = (AliPMD*)gAlice->GetDetector("PMD");
751 fPMDLoader = fRunLoader->GetLoader("PMDLoader");
752 if (fPMDLoader == 0x0)
754 cerr<<"AliPMDDigitizer::Exec : Can not find PMD or PMDLoader\n";
756 fPMDLoader->LoadDigits("recreate");
757 fTreeD = fPMDLoader->TreeD();
760 fPMDLoader->MakeTree("D");
761 fTreeD = fPMDLoader->TreeD();
763 Int_t bufsize = 16000;
764 fTreeD->Branch("PMDDigit", &fDigits, bufsize);
771 for (Int_t idet = 0; idet < 2; idet++)
773 for (Int_t ism = 0; ism < fgkTotUM; ism++)
775 for (Int_t jrow = 0; jrow < fgkRow; jrow++)
777 for (Int_t kcol = 0; kcol < fgkCol; kcol++)
781 deltaE = fPRE[ism][jrow][kcol];
782 trno = fPRETrackNo[ism][jrow][kcol];
787 deltaE = fCPV[ism][jrow][kcol];
788 trno = fCPVTrackNo[ism][jrow][kcol];
794 AddDigit(trno,detno,ism,jrow,kcol,adc);
800 } // supermodule loop
804 fPMDLoader->WriteDigits("OVERWRITE");
806 //____________________________________________________________________________
808 void AliPMDDigitizer::MergeSDigits(Int_t filenumber, Int_t troffset)
811 fRunLoader = AliRunLoader::GetRunLoader(fManager->GetInputFolderName(filenumber));
812 fPMDLoader = fRunLoader->GetLoader("PMDLoader");
813 fPMDLoader->LoadSDigits("read");
814 fTreeS = fPMDLoader->TreeS();
815 AliPMDsdigit *pmdsdigit;
816 TBranch *branch = fTreeS->GetBranch("PMDSDigit");
817 branch->SetAddress(&fSDigits);
819 Int_t itrackno, idet, ism;
823 Int_t nmodules = (Int_t) fTreeS->GetEntries();
826 cout << " nmodules = " << nmodules << endl;
827 cout << " tr offset = " << troffset << endl;
829 for (Int_t imodule = 0; imodule < nmodules; imodule++)
831 fTreeS->GetEntry(imodule);
832 Int_t nentries = fSDigits->GetLast();
835 cout << " nentries = " << nentries << endl;
837 for (Int_t ient = 0; ient < nentries+1; ient++)
839 pmdsdigit = (AliPMDsdigit*)fSDigits->UncheckedAt(ient);
840 itrackno = pmdsdigit->GetTrackNumber();
841 idet = pmdsdigit->GetDetector();
842 ism = pmdsdigit->GetSMNumber();
843 ixp = pmdsdigit->GetRow();
844 iyp = pmdsdigit->GetColumn();
845 edep = pmdsdigit->GetCellEdep();
849 if (fPRE[ism][ixp][iyp] < edep)
851 fPRETrackNo[ism][ixp][iyp] = troffset + itrackno;
853 fPRE[ism][ixp][iyp] += edep;
857 if (fCPV[ism][ixp][iyp] < edep)
859 fCPVTrackNo[ism][ixp][iyp] = troffset + itrackno;
861 fCPV[ism][ixp][iyp] += edep;
867 // ----------------------------------------------------------------------
868 void AliPMDDigitizer::TrackAssignment2Cell()
871 // This block assigns the cell id when there are
872 // multiple tracks in a cell according to the
875 Bool_t jsort = false;
887 pmdTrack = new Int_t ***[fgkTotUM];
888 pmdEdep = new Float_t ***[fgkTotUM];
889 for (i=0; i<fgkTotUM; i++)
891 pmdTrack[i] = new Int_t **[fgkRow];
892 pmdEdep[i] = new Float_t **[fgkRow];
895 for (i = 0; i < fgkTotUM; i++)
897 for (j = 0; j < fgkRow; j++)
899 pmdTrack[i][j] = new Int_t *[fgkCol];
900 pmdEdep[i][j] = new Float_t *[fgkCol];
904 for (i = 0; i < fgkTotUM; i++)
906 for (j = 0; j < fgkRow; j++)
908 for (k = 0; k < fgkCol; k++)
910 Int_t nn = fPRECounter[i][j][k];
913 pmdTrack[i][j][k] = new Int_t[nn];
914 pmdEdep[i][j][k] = new Float_t[nn];
919 pmdTrack[i][j][k] = new Int_t[nn];
920 pmdEdep[i][j][k] = new Float_t[nn];
922 fPRECounter[i][j][k] = 0;
928 Int_t nentries = fCell->GetEntries();
930 Int_t mtrackno, ism, ixp, iyp;
933 for (i = 0; i < nentries; i++)
935 fPMDcell = (AliPMDcell*)fCell->UncheckedAt(i);
937 mtrackno = fPMDcell->GetTrackNumber();
938 ism = fPMDcell->GetSMNumber();
939 ixp = fPMDcell->GetX();
940 iyp = fPMDcell->GetY();
941 edep = fPMDcell->GetEdep();
942 Int_t nn = fPRECounter[ism][ixp][iyp];
943 // cout << " nn = " << nn << endl;
944 pmdTrack[ism][ixp][iyp][nn] = (Int_t) mtrackno;
945 pmdEdep[ism][ixp][iyp][nn] = edep;
946 fPRECounter[ism][ixp][iyp]++;
953 for (im=0; im<fgkTotUM; im++)
955 for (ix=0; ix<fgkRow; ix++)
957 for (iy=0; iy<fgkCol; iy++)
959 nn = fPRECounter[im][ix][iy];
962 // This block handles if a cell is fired
963 // many times by many tracks
964 status1 = new Int_t[nn];
965 status2 = new Int_t[nn];
966 trnarray = new Int_t[nn];
967 for (iz = 0; iz < nn; iz++)
969 status1[iz] = pmdTrack[im][ix][iy][iz];
971 TMath::Sort(nn,status1,status2,jsort);
972 Int_t trackOld = -99999;
973 Int_t track, trCount = 0;
974 for (iz = 0; iz < nn; iz++)
976 track = status1[status2[iz]];
977 if (trackOld != track)
979 trnarray[trCount] = track;
987 trEdp = new Float_t[trCount];
988 fracEdp = new Float_t[trCount];
989 for (il = 0; il < trCount; il++)
992 track = trnarray[il];
993 for (iz = 0; iz < nn; iz++)
995 if (track == pmdTrack[im][ix][iy][iz])
997 trEdp[il] += pmdEdep[im][ix][iy][iz];
1000 totEdp += trEdp[il];
1003 Float_t fracOld = 0.;
1005 for (il = 0; il < trCount; il++)
1007 fracEdp[il] = trEdp[il]/totEdp;
1008 if (fracOld < fracEdp[il])
1010 fracOld = fracEdp[il];
1014 fPRETrackNo[im][ix][iy] = trnarray[ilOld];
1021 // This only handles if a cell is fired
1022 // by only one track
1024 fPRETrackNo[im][ix][iy] = pmdTrack[im][ix][iy][0];
1029 // This is if no cell is fired
1030 fPRETrackNo[im][ix][iy] = -999;
1036 // Delete all the pointers
1038 for (i = 0; i < fgkTotUM; i++)
1040 for (j = 0; j < fgkRow; j++)
1042 for (k = 0; k < fgkCol; k++)
1044 delete [] pmdTrack[i][j][k];
1045 delete [] pmdEdep[i][j][k];
1050 for (i = 0; i < fgkTotUM; i++)
1052 for (j = 0; j < fgkRow; j++)
1054 delete [] pmdTrack[i][j];
1055 delete [] pmdEdep[i][j];
1059 for (i = 0; i < fgkTotUM; i++)
1061 delete [] pmdTrack[i];
1062 delete [] pmdEdep[i];
1067 // End of the cell id assignment
1070 //____________________________________________________________________________
1071 void AliPMDDigitizer::MeV2ADC(Float_t mev, Float_t & adc) const
1073 // This converts the simulated edep to ADC according to the
1079 //____________________________________________________________________________
1080 void AliPMDDigitizer::AddSDigit(Int_t trnumber, Int_t det, Int_t smnumber,
1081 Int_t irow, Int_t icol, Float_t adc)
1085 TClonesArray &lsdigits = *fSDigits;
1086 new(lsdigits[fNsdigit++]) AliPMDsdigit(trnumber,det,smnumber,irow,icol,adc);
1088 //____________________________________________________________________________
1090 void AliPMDDigitizer::AddDigit(Int_t trnumber, Int_t det, Int_t smnumber,
1091 Int_t irow, Int_t icol, Float_t adc)
1095 TClonesArray &ldigits = *fDigits;
1096 new(ldigits[fNdigit++]) AliPMDdigit(trnumber,det,smnumber,irow,icol,adc);
1098 //____________________________________________________________________________
1100 void AliPMDDigitizer::SetZPosition(Float_t zpos)
1104 //____________________________________________________________________________
1105 Float_t AliPMDDigitizer::GetZPosition() const
1109 //____________________________________________________________________________
1111 void AliPMDDigitizer::ResetCell()
1113 // clears the cell array and also the counter
1117 for (Int_t i = 0; i < fgkTotUM; i++)
1119 for (Int_t j = 0; j < fgkRow; j++)
1121 for (Int_t k = 0; k < fgkCol; k++)
1123 fPRECounter[i][j][k] = 0;
1128 //____________________________________________________________________________
1129 void AliPMDDigitizer::ResetSDigit()
1133 if (fSDigits) fSDigits->Clear();
1135 //____________________________________________________________________________
1136 void AliPMDDigitizer::ResetDigit()
1140 if (fDigits) fDigits->Clear();
1142 //____________________________________________________________________________
1144 void AliPMDDigitizer::ResetCellADC()
1146 // Clears individual cells edep and track number
1147 for (Int_t i = 0; i < fgkTotUM; i++)
1149 for (Int_t j = 0; j < fgkRow; j++)
1151 for (Int_t k = 0; k < fgkCol; k++)
1155 fPRETrackNo[i][j][k] = 0;
1156 fCPVTrackNo[i][j][k] = 0;
1161 //____________________________________________________________________________
1163 void AliPMDDigitizer::UnLoad(Option_t *option)
1165 // Unloads all the root files
1167 const char *cS = strstr(option,"S");
1168 const char *cD = strstr(option,"D");
1170 fRunLoader->UnloadgAlice();
1171 fRunLoader->UnloadHeader();
1172 fRunLoader->UnloadKinematics();
1176 fPMDLoader->UnloadHits();
1180 fPMDLoader->UnloadHits();
1181 fPMDLoader->UnloadSDigits();