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>
37 #include "AliPMDhit.h"
39 #include "AliDetector.h"
40 #include "AliRunLoader.h"
41 #include "AliLoader.h"
42 #include "AliConfig.h"
44 #include "AliRunDigitizer.h"
45 #include "AliDigitizer.h"
46 #include "AliHeader.h"
48 #include "AliPMDcell.h"
49 #include "AliPMDsdigit.h"
50 #include "AliPMDdigit.h"
51 #include "AliPMDDigitizer.h"
54 ClassImp(AliPMDDigitizer)
56 AliPMDDigitizer::AliPMDDigitizer() :
67 fZPos(361.5)// in units of cm, This is the default position of PMD
69 // Default Constructor
71 for (Int_t i = 0; i < fgkTotUM; i++)
73 for (Int_t j = 0; j < fgkRow; j++)
75 for (Int_t k = 0; k < fgkCol; k++)
79 fPRECounter[i][j][k] = 0;
80 fPRETrackNo[i][j][k] = -1;
81 fCPVTrackNo[i][j][k] = -1;
87 //____________________________________________________________________________
88 AliPMDDigitizer::AliPMDDigitizer(AliRunDigitizer* manager)
89 :AliDigitizer(manager),
94 fSDigits(new TClonesArray("AliPMDsdigit", 1000)),
95 fDigits(new TClonesArray("AliPMDdigit", 1000)),
100 fZPos(361.5)// in units of cm, This is the default position of PMD
102 // ctor which should be used
104 for (Int_t i = 0; i < fgkTotUM; i++)
106 for (Int_t j = 0; j < fgkRow; j++)
108 for (Int_t k = 0; k < fgkCol; k++)
112 fPRECounter[i][j][k] = 0;
113 fPRETrackNo[i][j][k] = -1;
114 fCPVTrackNo[i][j][k] = -1;
119 //____________________________________________________________________________
120 AliPMDDigitizer::~AliPMDDigitizer()
122 // Default Destructor
139 //____________________________________________________________________________
140 void AliPMDDigitizer::OpengAliceFile(const char *file, Option_t *option)
142 // Loads galice.root file and corresponding header, kinematics
143 // hits and sdigits or digits depending on the option
146 TString evfoldname = AliConfig::GetDefaultEventFolderName();
147 fRunLoader = AliRunLoader::GetRunLoader(evfoldname);
149 fRunLoader = AliRunLoader::Open(file,AliConfig::GetDefaultEventFolderName(),
154 AliError(Form("Can not open session for file %s.",file));
157 if (!fRunLoader->GetAliRun()) fRunLoader->LoadgAlice();
158 if (!fRunLoader->TreeE()) fRunLoader->LoadHeader();
159 if (!fRunLoader->TreeK()) fRunLoader->LoadKinematics();
161 gAlice = fRunLoader->GetAliRun();
165 AliDebug(1,"Alirun object found");
169 AliError("Could not found Alirun object");
172 fPMD = (AliPMD*)gAlice->GetDetector("PMD");
173 fPMDLoader = fRunLoader->GetLoader("PMDLoader");
174 if (fPMDLoader == 0x0)
176 AliError("Can not find PMDLoader");
179 const char *cHS = strstr(option,"HS");
180 const char *cHD = strstr(option,"HD");
181 const char *cSD = strstr(option,"SD");
185 fPMDLoader->LoadHits("READ");
186 fPMDLoader->LoadSDigits("recreate");
190 fPMDLoader->LoadHits("READ");
191 fPMDLoader->LoadDigits("recreate");
195 fPMDLoader->LoadSDigits("READ");
196 fPMDLoader->LoadDigits("recreate");
199 //____________________________________________________________________________
200 void AliPMDDigitizer::Hits2SDigits(Int_t ievt)
202 // This reads the PMD Hits tree and assigns the right track number
203 // to a cell and stores in the summable digits tree
206 const Int_t kPi0 = 111;
207 const Int_t kGamma = 22;
215 Float_t xPos, yPos, zPos;
216 Int_t xpad = -1, ypad = -1;
218 Float_t vx = -999.0, vy = -999.0, vz = -999.0;
221 if (!fSDigits) fSDigits = new TClonesArray("AliPMDsdigit", 1000);
224 AliDebug(1,Form("Event Number = %d",ievt));
225 Int_t nparticles = fRunLoader->GetHeader()->GetNtrack();
226 AliDebug(1,Form("Number of Particles = %d",nparticles));
227 fRunLoader->GetEvent(ievt);
228 // ------------------------------------------------------- //
229 // Pointer to specific detector hits.
230 // Get pointers to Alice detectors and Hits containers
232 TTree* treeH = fPMDLoader->TreeH();
234 Int_t ntracks = (Int_t) treeH->GetEntries();
235 AliDebug(1,Form("Number of Tracks in the TreeH = %d", ntracks));
236 TTree* treeS = fPMDLoader->TreeS();
239 fPMDLoader->MakeTree("S");
240 treeS = fPMDLoader->TreeS();
242 Int_t bufsize = 16000;
243 treeS->Branch("PMDSDigit", &fSDigits, bufsize);
245 TClonesArray* hits = 0;
246 if (fPMD) hits = fPMD->Hits();
248 // Start loop on tracks in the hits containers
250 for (Int_t track=0; track<ntracks;track++)
253 treeH->GetEvent(track);
256 npmd = hits->GetEntriesFast();
257 for (int ipmd = 0; ipmd < npmd; ipmd++)
259 fPMDHit = (AliPMDhit*) hits->UncheckedAt(ipmd);
260 trackno = fPMDHit->GetTrack();
261 // get kinematics of the particles
263 TParticle* mparticle = gAlice->GetMCApp()->Particle(trackno);
264 trackpid = mparticle->GetPdgCode();
273 Int_t tracknoOld=0, trackpidOld=0, statusOld = 0;
274 if (mparticle->GetFirstMother() == -1)
276 tracknoOld = trackno;
277 trackpidOld = trackpid;
281 while((imo = mparticle->GetFirstMother()) >= 0)
285 mparticle = gAlice->GetMCApp()->Particle(imo);
286 idmo = mparticle->GetPdgCode();
288 vx = mparticle->Vx();
289 vy = mparticle->Vy();
290 vz = mparticle->Vz();
292 //printf("==> Mother ID %5d %5d %5d Vertex: %13.3f %13.3f %13.3f\n", igen, imo, idmo, vx, vy, vz);
293 //fprintf(fpw1,"==> Mother ID %5d %5d %5d Vertex: %13.3f %13.3f %13.3f\n", igen, imo, idmo, vx, vy, vz);
294 if ((idmo == kGamma || idmo == -11 || idmo == 11) && vx == 0. && vy == 0. && vz == 0.)
302 if (idmo == kPi0 && vx == 0. && vy == 0. && vz == 0.)
311 if (idmo == kPi0 && vx == 0. && vy == 0. && vz == 0.)
323 mtrackno = tracknoOld;
324 mtrackpid = trackpidOld;
330 edep = fPMDHit->GetEnergy();
331 Int_t vol1 = fPMDHit->GetVolume(1); // Column
332 Int_t vol2 = fPMDHit->GetVolume(2); // Row
333 Int_t vol3 = fPMDHit->GetVolume(3); // UnitModule
334 Int_t vol6 = fPMDHit->GetVolume(6); // SuperModule
336 // -----------------------------------------//
337 // For Super Module 1 & 2 //
338 // nrow = 96, ncol = 48 //
339 // For Super Module 3 & 4 //
340 // nrow = 48, ncol = 96 //
341 // -----------------------------------------//
343 smnumber = (vol6-1)*6 + vol3;
345 if (vol6 == 1 || vol6 == 2)
350 else if (vol6 == 3 || vol6 == 4)
356 AliDebug(2,Form("Zposition = %f Edeposition = %f",zPos,edep));
357 Float_t zposition = TMath::Abs(zPos);
358 if (zposition < fZPos)
363 else if (zposition > fZPos)
368 Int_t smn = smnumber - 1;
369 Int_t ixx = xpad - 1;
370 Int_t iyy = ypad - 1;
373 fPRE[smn][ixx][iyy] += edep;
374 fPRECounter[smn][ixx][iyy]++;
376 AliPMDcell* cell = new AliPMDcell(mtrackno,smn,ixx,iyy,edep);
382 fCPV[smn][ixx][iyy] += edep;
383 fCPVTrackNo[smn][ixx][iyy] = mtrackno;
387 } // Track Loop ended
389 TrackAssignment2Cell();
396 for (Int_t idet = 0; idet < 2; idet++)
398 for (Int_t ism = 0; ism < fgkTotUM; ism++)
400 for (Int_t jrow = 0; jrow < fgkRow; jrow++)
402 for (Int_t kcol = 0; kcol < fgkCol; kcol++)
406 deltaE = fPRE[ism][jrow][kcol];
407 trno = fPRETrackNo[ism][jrow][kcol];
412 deltaE = fCPV[ism][jrow][kcol];
413 trno = fCPVTrackNo[ism][jrow][kcol];
418 AddSDigit(trno,detno,ism,jrow,kcol,deltaE);
426 fPMDLoader->WriteSDigits("OVERWRITE");
430 //____________________________________________________________________________
432 void AliPMDDigitizer::Hits2Digits(Int_t ievt)
434 // This reads the PMD Hits tree and assigns the right track number
435 // to a cell and stores in the digits tree
437 const Int_t kPi0 = 111;
438 const Int_t kGamma = 22;
446 Float_t xPos, yPos, zPos;
447 Int_t xpad = -1, ypad = -1;
449 Float_t vx = -999.0, vy = -999.0, vz = -999.0;
451 if (!fDigits) fDigits = new TClonesArray("AliPMDdigit", 1000);
454 AliDebug(1,Form("Event Number = %d",ievt));
455 Int_t nparticles = fRunLoader->GetHeader()->GetNtrack();
456 AliDebug(1,Form("Number of Particles = %d", nparticles));
457 fRunLoader->GetEvent(ievt);
458 // ------------------------------------------------------- //
459 // Pointer to specific detector hits.
460 // Get pointers to Alice detectors and Hits containers
462 fPMD = (AliPMD*)gAlice->GetDetector("PMD");
463 fPMDLoader = fRunLoader->GetLoader("PMDLoader");
465 if (fPMDLoader == 0x0)
467 AliError("Can not find PMD or PMDLoader");
469 TTree* treeH = fPMDLoader->TreeH();
470 Int_t ntracks = (Int_t) treeH->GetEntries();
471 AliDebug(1,Form("Number of Tracks in the TreeH = %d", ntracks));
472 fPMDLoader->LoadDigits("recreate");
473 TTree* treeD = fPMDLoader->TreeD();
476 fPMDLoader->MakeTree("D");
477 treeD = fPMDLoader->TreeD();
479 Int_t bufsize = 16000;
480 treeD->Branch("PMDDigit", &fDigits, bufsize);
482 TClonesArray* hits = 0;
483 if (fPMD) hits = fPMD->Hits();
485 // Start loop on tracks in the hits containers
487 for (Int_t track=0; track<ntracks;track++)
490 treeH->GetEvent(track);
494 npmd = hits->GetEntriesFast();
495 for (int ipmd = 0; ipmd < npmd; ipmd++)
497 fPMDHit = (AliPMDhit*) hits->UncheckedAt(ipmd);
498 trackno = fPMDHit->GetTrack();
500 // get kinematics of the particles
502 TParticle* mparticle = gAlice->GetMCApp()->Particle(trackno);
503 trackpid = mparticle->GetPdgCode();
512 Int_t tracknoOld=0, trackpidOld=0, statusOld = 0;
513 if (mparticle->GetFirstMother() == -1)
515 tracknoOld = trackno;
516 trackpidOld = trackpid;
521 while((imo = mparticle->GetFirstMother()) >= 0)
525 mparticle = gAlice->GetMCApp()->Particle(imo);
526 idmo = mparticle->GetPdgCode();
528 vx = mparticle->Vx();
529 vy = mparticle->Vy();
530 vz = mparticle->Vz();
532 //printf("==> Mother ID %5d %5d %5d Vertex: %13.3f %13.3f %13.3f\n", igen, imo, idmo, vx, vy, vz);
533 //fprintf(fpw1,"==> Mother ID %5d %5d %5d Vertex: %13.3f %13.3f %13.3f\n", igen, imo, idmo, vx, vy, vz);
534 if ((idmo == kGamma || idmo == -11 || idmo == 11) && vx == 0. && vy == 0. && vz == 0.)
542 if (idmo == kPi0 && vx == 0. && vy == 0. && vz == 0.)
551 if (idmo == kPi0 && vx == 0. && vy == 0. && vz == 0.)
563 mtrackno = tracknoOld;
564 mtrackpid = trackpidOld;
571 edep = fPMDHit->GetEnergy();
572 Int_t vol1 = fPMDHit->GetVolume(1); // Column
573 Int_t vol2 = fPMDHit->GetVolume(2); // Row
574 Int_t vol3 = fPMDHit->GetVolume(3); // UnitModule
575 Int_t vol6 = fPMDHit->GetVolume(6); // SuperModule
577 // -----------------------------------------//
578 // For Super Module 1 & 2 //
579 // nrow = 96, ncol = 48 //
580 // For Super Module 3 & 4 //
581 // nrow = 48, ncol = 96 //
582 // -----------------------------------------//
584 smnumber = (vol6-1)*6 + vol3;
586 if (vol6 == 1 || vol6 == 2)
591 else if (vol6 == 3 || vol6 == 4)
597 AliDebug(2,Form("ZPosition = %f Edeposition = %d",zPos,edep));
598 Float_t zposition = TMath::Abs(zPos);
600 if (zposition < fZPos)
605 else if (zposition > fZPos)
611 Int_t smn = smnumber - 1;
612 Int_t ixx = xpad - 1;
613 Int_t iyy = ypad - 1;
616 fPRE[smn][ixx][iyy] += edep;
617 fPRECounter[smn][ixx][iyy]++;
619 AliPMDcell* cell = new AliPMDcell(mtrackno,smn,ixx,iyy,edep);
625 fCPV[smn][ixx][iyy] += edep;
626 fCPVTrackNo[smn][ixx][iyy] = mtrackno;
630 } // Track Loop ended
632 TrackAssignment2Cell();
640 for (Int_t idet = 0; idet < 2; idet++)
642 for (Int_t ism = 0; ism < fgkTotUM; ism++)
644 for (Int_t jrow = 0; jrow < fgkRow; jrow++)
646 for (Int_t kcol = 0; kcol < fgkCol; kcol++)
650 deltaE = fPRE[ism][jrow][kcol];
651 trno = fPRETrackNo[ism][jrow][kcol];
656 deltaE = fCPV[ism][jrow][kcol];
657 trno = fCPVTrackNo[ism][jrow][kcol];
663 AddDigit(trno,detno,ism,jrow,kcol,adc);
669 } // supermodule loop
672 fPMDLoader->WriteDigits("OVERWRITE");
676 //____________________________________________________________________________
679 void AliPMDDigitizer::SDigits2Digits(Int_t ievt)
681 // This reads the PMD sdigits tree and converts energy deposition
682 // in a cell to ADC and stores in the digits tree
685 fRunLoader->GetEvent(ievt);
687 TTree* treeS = fPMDLoader->TreeS();
688 AliPMDsdigit *pmdsdigit;
689 TBranch *branch = treeS->GetBranch("PMDSDigit");
692 AliError("PMD Sdigit branch does not exist");
695 if (!fSDigits) fSDigits = new TClonesArray("AliPMDsdigit", 1000);
696 branch->SetAddress(&fSDigits);
698 TTree* treeD = fPMDLoader->TreeD();
701 fPMDLoader->MakeTree("D");
702 treeD = fPMDLoader->TreeD();
704 Int_t bufsize = 16000;
705 if (!fDigits) fDigits = new TClonesArray("AliPMDdigit", 1000);
706 treeD->Branch("PMDDigit", &fDigits, bufsize);
708 Int_t trno, det, smn;
712 Int_t nmodules = (Int_t) treeS->GetEntries();
713 AliDebug(1,Form("Number of modules = %d",nmodules));
715 for (Int_t imodule = 0; imodule < nmodules; imodule++)
717 treeS->GetEntry(imodule);
718 Int_t nentries = fSDigits->GetLast();
719 AliDebug(2,Form("Number of entries per module = %d",nentries+1));
720 for (Int_t ient = 0; ient < nentries+1; ient++)
722 pmdsdigit = (AliPMDsdigit*)fSDigits->UncheckedAt(ient);
723 trno = pmdsdigit->GetTrackNumber();
724 det = pmdsdigit->GetDetector();
725 smn = pmdsdigit->GetSMNumber();
726 irow = pmdsdigit->GetRow();
727 icol = pmdsdigit->GetColumn();
728 edep = pmdsdigit->GetCellEdep();
731 AddDigit(trno,det,smn,irow,icol,adc);
736 fPMDLoader->WriteDigits("OVERWRITE");
739 //____________________________________________________________________________
740 void AliPMDDigitizer::Exec(Option_t *option)
742 // Does the event merging and digitization
743 const char *cdeb = strstr(option,"deb");
746 AliDebug(100," *** PMD Exec is called ***");
749 Int_t ninputs = fManager->GetNinputs();
750 AliDebug(1,Form("Number of files to be processed = %d",ninputs));
753 for (Int_t i = 0; i < ninputs; i++)
755 Int_t troffset = fManager->GetMask(i);
756 MergeSDigits(i, troffset);
759 fRunLoader = AliRunLoader::GetRunLoader(fManager->GetOutputFolderName());
760 fPMD = (AliPMD*)gAlice->GetDetector("PMD");
761 fPMDLoader = fRunLoader->GetLoader("PMDLoader");
762 if (fPMDLoader == 0x0)
764 AliError("Can not find PMD or PMDLoader");
766 fPMDLoader->LoadDigits("update");
767 TTree* treeD = fPMDLoader->TreeD();
770 fPMDLoader->MakeTree("D");
771 treeD = fPMDLoader->TreeD();
773 Int_t bufsize = 16000;
774 if (!fDigits) fDigits = new TClonesArray("AliPMDdigit", 1000);
775 treeD->Branch("PMDDigit", &fDigits, bufsize);
782 for (Int_t idet = 0; idet < 2; idet++)
784 for (Int_t ism = 0; ism < fgkTotUM; ism++)
786 for (Int_t jrow = 0; jrow < fgkRow; jrow++)
788 for (Int_t kcol = 0; kcol < fgkCol; kcol++)
792 deltaE = fPRE[ism][jrow][kcol];
793 trno = fPRETrackNo[ism][jrow][kcol];
798 deltaE = fCPV[ism][jrow][kcol];
799 trno = fCPVTrackNo[ism][jrow][kcol];
805 AddDigit(trno,detno,ism,jrow,kcol,adc);
811 } // supermodule loop
813 fPMDLoader->WriteDigits("OVERWRITE");
814 fPMDLoader->UnloadDigits();
817 //____________________________________________________________________________
819 void AliPMDDigitizer::MergeSDigits(Int_t filenumber, Int_t troffset)
822 fRunLoader = AliRunLoader::GetRunLoader(fManager->GetInputFolderName(filenumber));
823 fPMDLoader = fRunLoader->GetLoader("PMDLoader");
824 fPMDLoader->LoadSDigits("read");
825 TTree* treeS = fPMDLoader->TreeS();
826 AliPMDsdigit *pmdsdigit;
827 TBranch *branch = treeS->GetBranch("PMDSDigit");
828 if (!fSDigits) fSDigits = new TClonesArray("AliPMDsdigit", 1000);
829 branch->SetAddress(&fSDigits);
831 Int_t itrackno, idet, ism;
835 Int_t nmodules = (Int_t) treeS->GetEntries();
836 AliDebug(1,Form("Number of Modules in the treeS = %d",nmodules));
837 AliDebug(1,Form("Track Offset = %d",troffset));
838 for (Int_t imodule = 0; imodule < nmodules; imodule++)
840 treeS->GetEntry(imodule);
841 Int_t nentries = fSDigits->GetLast();
842 AliDebug(2,Form("Number of Entries per Module = %d",nentries));
843 for (Int_t ient = 0; ient < nentries+1; ient++)
845 pmdsdigit = (AliPMDsdigit*)fSDigits->UncheckedAt(ient);
846 itrackno = pmdsdigit->GetTrackNumber();
847 idet = pmdsdigit->GetDetector();
848 ism = pmdsdigit->GetSMNumber();
849 ixp = pmdsdigit->GetRow();
850 iyp = pmdsdigit->GetColumn();
851 edep = pmdsdigit->GetCellEdep();
855 if (fPRE[ism][ixp][iyp] < edep)
857 fPRETrackNo[ism][ixp][iyp] = troffset + itrackno;
859 fPRE[ism][ixp][iyp] += edep;
863 if (fCPV[ism][ixp][iyp] < edep)
865 fCPVTrackNo[ism][ixp][iyp] = troffset + itrackno;
867 fCPV[ism][ixp][iyp] += edep;
873 // ----------------------------------------------------------------------
874 void AliPMDDigitizer::TrackAssignment2Cell()
877 // This block assigns the cell id when there are
878 // multiple tracks in a cell according to the
881 Bool_t jsort = false;
893 pmdTrack = new Int_t ***[fgkTotUM];
894 pmdEdep = new Float_t ***[fgkTotUM];
895 for (i=0; i<fgkTotUM; i++)
897 pmdTrack[i] = new Int_t **[fgkRow];
898 pmdEdep[i] = new Float_t **[fgkRow];
901 for (i = 0; i < fgkTotUM; i++)
903 for (j = 0; j < fgkRow; j++)
905 pmdTrack[i][j] = new Int_t *[fgkCol];
906 pmdEdep[i][j] = new Float_t *[fgkCol];
910 for (i = 0; i < fgkTotUM; i++)
912 for (j = 0; j < fgkRow; j++)
914 for (k = 0; k < fgkCol; k++)
916 Int_t nn = fPRECounter[i][j][k];
919 pmdTrack[i][j][k] = new Int_t[nn];
920 pmdEdep[i][j][k] = new Float_t[nn];
925 pmdTrack[i][j][k] = new Int_t[nn];
926 pmdEdep[i][j][k] = new Float_t[nn];
928 fPRECounter[i][j][k] = 0;
934 Int_t nentries = fCell.GetEntries();
936 Int_t mtrackno, ism, ixp, iyp;
939 for (i = 0; i < nentries; i++)
941 AliPMDcell* cell = (AliPMDcell*)fCell.UncheckedAt(i);
943 mtrackno = cell->GetTrackNumber();
944 ism = cell->GetSMNumber();
947 edep = cell->GetEdep();
948 Int_t nn = fPRECounter[ism][ixp][iyp];
949 pmdTrack[ism][ixp][iyp][nn] = (Int_t) mtrackno;
950 pmdEdep[ism][ixp][iyp][nn] = edep;
951 fPRECounter[ism][ixp][iyp]++;
958 for (im=0; im<fgkTotUM; im++)
960 for (ix=0; ix<fgkRow; ix++)
962 for (iy=0; iy<fgkCol; iy++)
964 nn = fPRECounter[im][ix][iy];
967 // This block handles if a cell is fired
968 // many times by many tracks
969 status1 = new Int_t[nn];
970 status2 = new Int_t[nn];
971 trnarray = new Int_t[nn];
972 for (iz = 0; iz < nn; iz++)
974 status1[iz] = pmdTrack[im][ix][iy][iz];
976 TMath::Sort(nn,status1,status2,jsort);
977 Int_t trackOld = -99999;
978 Int_t track, trCount = 0;
979 for (iz = 0; iz < nn; iz++)
981 track = status1[status2[iz]];
982 if (trackOld != track)
984 trnarray[trCount] = track;
992 trEdp = new Float_t[trCount];
993 fracEdp = new Float_t[trCount];
994 for (il = 0; il < trCount; il++)
997 track = trnarray[il];
998 for (iz = 0; iz < nn; iz++)
1000 if (track == pmdTrack[im][ix][iy][iz])
1002 trEdp[il] += pmdEdep[im][ix][iy][iz];
1005 totEdp += trEdp[il];
1008 Float_t fracOld = 0.;
1010 for (il = 0; il < trCount; il++)
1012 fracEdp[il] = trEdp[il]/totEdp;
1013 if (fracOld < fracEdp[il])
1015 fracOld = fracEdp[il];
1019 fPRETrackNo[im][ix][iy] = trnarray[ilOld];
1026 // This only handles if a cell is fired
1027 // by only one track
1029 fPRETrackNo[im][ix][iy] = pmdTrack[im][ix][iy][0];
1034 // This is if no cell is fired
1035 fPRETrackNo[im][ix][iy] = -999;
1041 // Delete all the pointers
1043 for (i = 0; i < fgkTotUM; i++)
1045 for (j = 0; j < fgkRow; j++)
1047 for (k = 0; k < fgkCol; k++)
1049 delete [] pmdTrack[i][j][k];
1050 delete [] pmdEdep[i][j][k];
1055 for (i = 0; i < fgkTotUM; i++)
1057 for (j = 0; j < fgkRow; j++)
1059 delete [] pmdTrack[i][j];
1060 delete [] pmdEdep[i][j];
1064 for (i = 0; i < fgkTotUM; i++)
1066 delete [] pmdTrack[i];
1067 delete [] pmdEdep[i];
1072 // End of the cell id assignment
1075 //____________________________________________________________________________
1076 void AliPMDDigitizer::MeV2ADC(Float_t mev, Float_t & adc) const
1078 // This converts the simulated edep to ADC according to the
1083 // PS Test in September 2003
1084 // MeV - ADC conversion for 10bit ADC
1086 const Float_t kConstant = 7.181;
1087 const Float_t kErConstant = 0.6899;
1088 const Float_t kSlope = 35.93;
1089 const Float_t kErSlope = 0.306;
1091 //gRandom->SetSeed();
1093 Float_t cons = gRandom->Gaus(kConstant,kErConstant);
1094 Float_t slop = gRandom->Gaus(kSlope,kErSlope);
1096 Float_t adc10bit = slop*mev*0.001 + cons;
1100 Int_t adc12bit = (Int_t) (4.0*adc10bit);
1104 adc = (Float_t) adc12bit;
1106 else if (adc12bit >= 3000)
1112 //____________________________________________________________________________
1113 void AliPMDDigitizer::AddSDigit(Int_t trnumber, Int_t det, Int_t smnumber,
1114 Int_t irow, Int_t icol, Float_t adc)
1118 if (!fSDigits) fSDigits = new TClonesArray("AliPMDsdigit", 1000);
1119 TClonesArray &lsdigits = *fSDigits;
1120 new(lsdigits[fNsdigit++]) AliPMDsdigit(trnumber,det,smnumber,irow,icol,adc);
1122 //____________________________________________________________________________
1124 void AliPMDDigitizer::AddDigit(Int_t trnumber, Int_t det, Int_t smnumber,
1125 Int_t irow, Int_t icol, Float_t adc)
1129 if (!fDigits) fDigits = new TClonesArray("AliPMDdigit", 1000);
1130 TClonesArray &ldigits = *fDigits;
1131 new(ldigits[fNdigit++]) AliPMDdigit(trnumber,det,smnumber,irow,icol,adc);
1133 //____________________________________________________________________________
1135 void AliPMDDigitizer::SetZPosition(Float_t zpos)
1139 //____________________________________________________________________________
1140 Float_t AliPMDDigitizer::GetZPosition() const
1144 //____________________________________________________________________________
1146 void AliPMDDigitizer::ResetCell()
1148 // clears the cell array and also the counter
1152 for (Int_t i = 0; i < fgkTotUM; i++)
1154 for (Int_t j = 0; j < fgkRow; j++)
1156 for (Int_t k = 0; k < fgkCol; k++)
1158 fPRECounter[i][j][k] = 0;
1163 //____________________________________________________________________________
1164 void AliPMDDigitizer::ResetSDigit()
1168 if (fSDigits) fSDigits->Delete();
1170 //____________________________________________________________________________
1171 void AliPMDDigitizer::ResetDigit()
1175 if (fDigits) fDigits->Delete();
1177 //____________________________________________________________________________
1179 void AliPMDDigitizer::ResetCellADC()
1181 // Clears individual cells edep and track number
1182 for (Int_t i = 0; i < fgkTotUM; i++)
1184 for (Int_t j = 0; j < fgkRow; j++)
1186 for (Int_t k = 0; k < fgkCol; k++)
1190 fPRETrackNo[i][j][k] = 0;
1191 fCPVTrackNo[i][j][k] = 0;
1196 //____________________________________________________________________________
1198 void AliPMDDigitizer::UnLoad(Option_t *option)
1200 // Unloads all the root files
1202 const char *cS = strstr(option,"S");
1203 const char *cD = strstr(option,"D");
1205 fRunLoader->UnloadgAlice();
1206 fRunLoader->UnloadHeader();
1207 fRunLoader->UnloadKinematics();
1211 fPMDLoader->UnloadHits();
1215 fPMDLoader->UnloadHits();
1216 fPMDLoader->UnloadSDigits();