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>
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"
46 #include "AliCDBManager.h"
47 #include "AliCDBStorage.h"
48 #include "AliCDBEntry.h"
52 #include "AliPMDhit.h"
53 #include "AliPMDcell.h"
54 #include "AliPMDsdigit.h"
55 #include "AliPMDdigit.h"
56 #include "AliPMDCalibData.h"
57 #include "AliPMDPedestal.h"
58 #include "AliPMDDigitizer.h"
61 ClassImp(AliPMDDigitizer)
63 AliPMDDigitizer::AliPMDDigitizer() :
68 fCalibGain(GetCalibGain()),
69 fCalibPed(GetCalibPed()),
76 fZPos(361.5) // in units of cm, default position of PMD
78 // Default Constructor
80 for (Int_t i = 0; i < fgkTotUM; i++)
82 for (Int_t j = 0; j < fgkRow; j++)
84 for (Int_t k = 0; k < fgkCol; k++)
88 fPRECounter[i][j][k] = 0;
89 fPRETrackNo[i][j][k] = -1;
90 fCPVTrackNo[i][j][k] = -1;
96 //____________________________________________________________________________
97 AliPMDDigitizer::AliPMDDigitizer(const AliPMDDigitizer& digitizer):
98 AliDigitizer(digitizer),
103 fCalibGain(GetCalibGain()),
104 fCalibPed(GetCalibPed()),
111 fZPos(361.5) // in units of cm, default position of PMD
114 AliError("Copy constructor not allowed ");
117 //____________________________________________________________________________
118 AliPMDDigitizer & AliPMDDigitizer::operator=(const AliPMDDigitizer& /*digitizer*/)
120 // Assignment operator
121 AliError("Assignement operator not allowed ");
125 //____________________________________________________________________________
126 AliPMDDigitizer::AliPMDDigitizer(AliRunDigitizer* manager):
127 AliDigitizer(manager),
132 fCalibGain(GetCalibGain()),
133 fCalibPed(GetCalibPed()),
134 fSDigits(new TClonesArray("AliPMDsdigit", 1000)),
135 fDigits(new TClonesArray("AliPMDdigit", 1000)),
140 fZPos(361.5)// in units of cm, This is the default position of PMD
142 // ctor which should be used
145 for (Int_t i = 0; i < fgkTotUM; i++)
147 for (Int_t j = 0; j < fgkRow; j++)
149 for (Int_t k = 0; k < fgkCol; k++)
153 fPRECounter[i][j][k] = 0;
154 fPRETrackNo[i][j][k] = -1;
155 fCPVTrackNo[i][j][k] = -1;
161 //____________________________________________________________________________
162 AliPMDDigitizer::~AliPMDDigitizer()
164 // Default Destructor
181 //____________________________________________________________________________
182 void AliPMDDigitizer::OpengAliceFile(const char *file, Option_t *option)
184 // Loads galice.root file and corresponding header, kinematics
185 // hits and sdigits or digits depending on the option
188 TString evfoldname = AliConfig::GetDefaultEventFolderName();
189 fRunLoader = AliRunLoader::GetRunLoader(evfoldname);
191 fRunLoader = AliRunLoader::Open(file,AliConfig::GetDefaultEventFolderName(), "UPDATE");
195 AliError(Form("Can not open session for file %s.",file));
198 const char *cHS = strstr(option,"HS");
199 const char *cHD = strstr(option,"HD");
200 const char *cSD = strstr(option,"SD");
204 if (!fRunLoader->GetAliRun()) fRunLoader->LoadgAlice();
205 if (!fRunLoader->TreeE()) fRunLoader->LoadHeader();
206 if (!fRunLoader->TreeK()) fRunLoader->LoadKinematics();
208 gAlice = fRunLoader->GetAliRun();
212 AliDebug(1,"Alirun object found");
216 AliError("Could not found Alirun object");
219 fPMD = (AliPMD*)gAlice->GetDetector("PMD");
222 fPMDLoader = fRunLoader->GetLoader("PMDLoader");
223 if (fPMDLoader == 0x0)
225 AliError("Can not find PMDLoader");
231 fPMDLoader->LoadHits("READ");
232 fPMDLoader->LoadSDigits("recreate");
236 fPMDLoader->LoadHits("READ");
237 fPMDLoader->LoadDigits("recreate");
241 fPMDLoader->LoadSDigits("READ");
242 fPMDLoader->LoadDigits("recreate");
245 //____________________________________________________________________________
246 void AliPMDDigitizer::Hits2SDigits(Int_t ievt)
248 // This reads the PMD Hits tree and assigns the right track number
249 // to a cell and stores in the summable digits tree
252 const Int_t kPi0 = 111;
253 const Int_t kGamma = 22;
261 Float_t xPos, yPos, zPos;
262 Int_t xpad = -1, ypad = -1;
264 Float_t vx = -999.0, vy = -999.0, vz = -999.0;
267 if (!fSDigits) fSDigits = new TClonesArray("AliPMDsdigit", 1000);
270 AliDebug(1,Form("Event Number = %d",ievt));
271 Int_t nparticles = fRunLoader->GetHeader()->GetNtrack();
272 AliDebug(1,Form("Number of Particles = %d",nparticles));
273 fRunLoader->GetEvent(ievt);
274 // ------------------------------------------------------- //
275 // Pointer to specific detector hits.
276 // Get pointers to Alice detectors and Hits containers
278 TTree* treeH = fPMDLoader->TreeH();
280 Int_t ntracks = (Int_t) treeH->GetEntries();
281 AliDebug(1,Form("Number of Tracks in the TreeH = %d", ntracks));
282 TTree* treeS = fPMDLoader->TreeS();
285 fPMDLoader->MakeTree("S");
286 treeS = fPMDLoader->TreeS();
288 Int_t bufsize = 16000;
289 treeS->Branch("PMDSDigit", &fSDigits, bufsize);
291 TClonesArray* hits = 0;
292 if (fPMD) hits = fPMD->Hits();
294 // Start loop on tracks in the hits containers
296 for (Int_t track=0; track<ntracks;track++)
299 treeH->GetEvent(track);
302 npmd = hits->GetEntriesFast();
303 for (int ipmd = 0; ipmd < npmd; ipmd++)
305 fPMDHit = (AliPMDhit*) hits->UncheckedAt(ipmd);
306 trackno = fPMDHit->GetTrack();
307 // get kinematics of the particles
309 TParticle* mparticle = gAlice->GetMCApp()->Particle(trackno);
310 trackpid = mparticle->GetPdgCode();
319 Int_t tracknoOld=0, trackpidOld=0, statusOld = 0;
320 if (mparticle->GetFirstMother() == -1)
322 tracknoOld = trackno;
323 trackpidOld = trackpid;
327 while((imo = mparticle->GetFirstMother()) >= 0)
331 mparticle = gAlice->GetMCApp()->Particle(imo);
332 idmo = mparticle->GetPdgCode();
334 vx = mparticle->Vx();
335 vy = mparticle->Vy();
336 vz = mparticle->Vz();
338 //printf("==> Mother ID %5d %5d %5d Vertex: %13.3f %13.3f %13.3f\n", igen, imo, idmo, vx, vy, vz);
339 //fprintf(fpw1,"==> Mother ID %5d %5d %5d Vertex: %13.3f %13.3f %13.3f\n", igen, imo, idmo, vx, vy, vz);
340 if ((idmo == kGamma || idmo == -11 || idmo == 11) && vx == 0. && vy == 0. && vz == 0.)
348 if (idmo == kPi0 && vx == 0. && vy == 0. && vz == 0.)
357 if (idmo == kPi0 && vx == 0. && vy == 0. && vz == 0.)
369 mtrackno = tracknoOld;
370 mtrackpid = trackpidOld;
376 edep = fPMDHit->GetEnergy();
377 Int_t vol1 = fPMDHit->GetVolume(1); // Column
378 Int_t vol2 = fPMDHit->GetVolume(2); // Row
379 Int_t vol7 = fPMDHit->GetVolume(7); // UnitModule
380 Int_t vol8 = fPMDHit->GetVolume(8); // SuperModule
383 // -----------------------------------------//
384 // In new geometry after adding electronics //
385 // For Super Module 1 & 2 //
386 // nrow = 48, ncol = 96 //
387 // For Super Module 3 & 4 //
388 // nrow = 96, ncol = 48 //
389 // -----------------------------------------//
393 smnumber = (vol8-1)*6 + vol7;
395 if (vol8 == 1 || vol8 == 2)
400 else if (vol8 == 3 || vol8 == 4)
406 AliDebug(2,Form("Zposition = %f Edeposition = %f",zPos,edep));
407 //Float_t zposition = TMath::Abs(zPos);
413 else if (zPos > fZPos)
418 Int_t smn = smnumber - 1;
419 Int_t ixx = xpad - 1;
420 Int_t iyy = ypad - 1;
423 fPRE[smn][ixx][iyy] += edep;
424 fPRECounter[smn][ixx][iyy]++;
426 AliPMDcell* cell = new AliPMDcell(mtrackno,smn,ixx,iyy,edep);
431 fCPV[smn][ixx][iyy] += edep;
432 fCPVTrackNo[smn][ixx][iyy] = mtrackno;
436 } // Track Loop ended
438 TrackAssignment2Cell();
445 for (Int_t idet = 0; idet < 2; idet++)
447 for (Int_t ism = 0; ism < fgkTotUM; ism++)
449 for (Int_t jrow = 0; jrow < fgkRow; jrow++)
451 for (Int_t kcol = 0; kcol < fgkCol; kcol++)
455 deltaE = fPRE[ism][jrow][kcol];
456 trno = fPRETrackNo[ism][jrow][kcol];
461 deltaE = fCPV[ism][jrow][kcol];
462 trno = fCPVTrackNo[ism][jrow][kcol];
467 AddSDigit(trno,detno,ism,jrow,kcol,deltaE);
475 fPMDLoader->WriteSDigits("OVERWRITE");
478 //____________________________________________________________________________
480 void AliPMDDigitizer::Hits2Digits(Int_t ievt)
482 // This reads the PMD Hits tree and assigns the right track number
483 // to a cell and stores in the digits tree
485 const Int_t kPi0 = 111;
486 const Int_t kGamma = 22;
494 Float_t xPos, yPos, zPos;
495 Int_t xpad = -1, ypad = -1;
497 Float_t vx = -999.0, vy = -999.0, vz = -999.0;
499 if (!fDigits) fDigits = new TClonesArray("AliPMDdigit", 1000);
502 AliDebug(1,Form("Event Number = %d",ievt));
503 Int_t nparticles = fRunLoader->GetHeader()->GetNtrack();
504 AliDebug(1,Form("Number of Particles = %d", nparticles));
505 fRunLoader->GetEvent(ievt);
506 // ------------------------------------------------------- //
507 // Pointer to specific detector hits.
508 // Get pointers to Alice detectors and Hits containers
510 fPMD = (AliPMD*)gAlice->GetDetector("PMD");
511 fPMDLoader = fRunLoader->GetLoader("PMDLoader");
513 if (fPMDLoader == 0x0)
515 AliError("Can not find PMD or PMDLoader");
517 TTree* treeH = fPMDLoader->TreeH();
518 Int_t ntracks = (Int_t) treeH->GetEntries();
519 AliDebug(1,Form("Number of Tracks in the TreeH = %d", ntracks));
520 fPMDLoader->LoadDigits("recreate");
521 TTree* treeD = fPMDLoader->TreeD();
524 fPMDLoader->MakeTree("D");
525 treeD = fPMDLoader->TreeD();
527 Int_t bufsize = 16000;
528 treeD->Branch("PMDDigit", &fDigits, bufsize);
530 TClonesArray* hits = 0;
531 if (fPMD) hits = fPMD->Hits();
533 // Start loop on tracks in the hits containers
535 for (Int_t track=0; track<ntracks;track++)
538 treeH->GetEvent(track);
542 npmd = hits->GetEntriesFast();
543 for (int ipmd = 0; ipmd < npmd; ipmd++)
545 fPMDHit = (AliPMDhit*) hits->UncheckedAt(ipmd);
546 trackno = fPMDHit->GetTrack();
548 // get kinematics of the particles
550 TParticle* mparticle = gAlice->GetMCApp()->Particle(trackno);
551 trackpid = mparticle->GetPdgCode();
560 Int_t tracknoOld=0, trackpidOld=0, statusOld = 0;
561 if (mparticle->GetFirstMother() == -1)
563 tracknoOld = trackno;
564 trackpidOld = trackpid;
569 while((imo = mparticle->GetFirstMother()) >= 0)
573 mparticle = gAlice->GetMCApp()->Particle(imo);
574 idmo = mparticle->GetPdgCode();
576 vx = mparticle->Vx();
577 vy = mparticle->Vy();
578 vz = mparticle->Vz();
580 //printf("==> Mother ID %5d %5d %5d Vertex: %13.3f %13.3f %13.3f\n", igen, imo, idmo, vx, vy, vz);
581 //fprintf(fpw1,"==> Mother ID %5d %5d %5d Vertex: %13.3f %13.3f %13.3f\n", igen, imo, idmo, vx, vy, vz);
582 if ((idmo == kGamma || idmo == -11 || idmo == 11) && vx == 0. && vy == 0. && vz == 0.)
590 if (idmo == kPi0 && vx == 0. && vy == 0. && vz == 0.)
599 if (idmo == kPi0 && vx == 0. && vy == 0. && vz == 0.)
611 mtrackno = tracknoOld;
612 mtrackpid = trackpidOld;
618 edep = fPMDHit->GetEnergy();
619 Int_t vol1 = fPMDHit->GetVolume(1); // Column
620 Int_t vol2 = fPMDHit->GetVolume(2); // Row
621 Int_t vol7 = fPMDHit->GetVolume(7); // UnitModule
622 Int_t vol8 = fPMDHit->GetVolume(8); // SuperModule
625 // -----------------------------------------//
626 // In new geometry after adding electronics //
627 // For Super Module 1 & 2 //
628 // nrow = 48, ncol = 96 //
629 // For Super Module 3 & 4 //
630 // nrow = 96, ncol = 48 //
631 // -----------------------------------------//
633 smnumber = (vol8-1)*6 + vol7;
635 if (vol8 == 1 || vol8 == 2)
640 else if (vol8 == 3 || vol8 == 4)
646 AliDebug(2,Form("ZPosition = %f Edeposition = %d",zPos,edep));
647 //Float_t zposition = TMath::Abs(zPos);
654 else if (zPos > fZPos)
660 Int_t smn = smnumber - 1;
661 Int_t ixx = xpad - 1;
662 Int_t iyy = ypad - 1;
665 fPRE[smn][ixx][iyy] += edep;
666 fPRECounter[smn][ixx][iyy]++;
668 AliPMDcell* cell = new AliPMDcell(mtrackno,smn,ixx,iyy,edep);
674 fCPV[smn][ixx][iyy] += edep;
675 fCPVTrackNo[smn][ixx][iyy] = mtrackno;
679 } // Track Loop ended
681 TrackAssignment2Cell();
689 for (Int_t idet = 0; idet < 2; idet++)
691 for (Int_t ism = 0; ism < fgkTotUM; ism++)
693 for (Int_t jrow = 0; jrow < fgkRow; jrow++)
695 for (Int_t kcol = 0; kcol < fgkCol; kcol++)
699 deltaE = fPRE[ism][jrow][kcol];
700 trno = fPRETrackNo[ism][jrow][kcol];
705 deltaE = fCPV[ism][jrow][kcol];
706 trno = fCPVTrackNo[ism][jrow][kcol];
713 // To decalibrate the adc values
715 gain1 = Gain(idet,ism,jrow,kcol);
718 Int_t adcDecalib = (Int_t)(adc/gain1);
719 adc = (Float_t) adcDecalib;
726 // Pedestal Decalibration
728 fCalibPed->GetPedMeanRms(idet,ism,jrow,kcol);
729 Int_t pedrms1 = (Int_t) pedmeanrms%1000;
730 Float_t pedrms = (Float_t)pedrms1/10.;
732 (Float_t) (pedmeanrms - pedrms1)/1000.0;
733 //printf("%f %f\n",pedmean, pedrms);
736 adc += (pedmean + 3.0*pedrms);
737 AddDigit(trno,detno,ism,jrow,kcol,adc);
744 } // supermodule loop
747 fPMDLoader->WriteDigits("OVERWRITE");
751 //____________________________________________________________________________
754 void AliPMDDigitizer::SDigits2Digits(Int_t ievt)
756 // This reads the PMD sdigits tree and converts energy deposition
757 // in a cell to ADC and stores in the digits tree
760 fRunLoader->GetEvent(ievt);
762 TTree* treeS = fPMDLoader->TreeS();
763 AliPMDsdigit *pmdsdigit;
764 TBranch *branch = treeS->GetBranch("PMDSDigit");
767 AliError("PMD Sdigit branch does not exist");
770 if (!fSDigits) fSDigits = new TClonesArray("AliPMDsdigit", 1000);
771 branch->SetAddress(&fSDigits);
773 TTree* treeD = fPMDLoader->TreeD();
776 fPMDLoader->MakeTree("D");
777 treeD = fPMDLoader->TreeD();
779 Int_t bufsize = 16000;
780 if (!fDigits) fDigits = new TClonesArray("AliPMDdigit", 1000);
781 treeD->Branch("PMDDigit", &fDigits, bufsize);
783 Int_t trno, det, smn;
787 Int_t nmodules = (Int_t) treeS->GetEntries();
788 AliDebug(1,Form("Number of modules = %d",nmodules));
790 for (Int_t imodule = 0; imodule < nmodules; imodule++)
792 treeS->GetEntry(imodule);
793 Int_t nentries = fSDigits->GetLast();
794 AliDebug(2,Form("Number of entries per module = %d",nentries+1));
795 for (Int_t ient = 0; ient < nentries+1; ient++)
797 pmdsdigit = (AliPMDsdigit*)fSDigits->UncheckedAt(ient);
798 trno = pmdsdigit->GetTrackNumber();
799 det = pmdsdigit->GetDetector();
800 smn = pmdsdigit->GetSMNumber();
801 irow = pmdsdigit->GetRow();
802 icol = pmdsdigit->GetColumn();
803 edep = pmdsdigit->GetCellEdep();
808 // To decalibrte the adc values
810 Float_t gain1 = Gain(det,smn,irow,icol);
813 Int_t adcDecalib = (Int_t)(adc/gain1);
814 adc = (Float_t) adcDecalib;
820 // Pedestal Decalibration
821 Int_t pedmeanrms = fCalibPed->GetPedMeanRms(det,smn,irow,icol);
822 Int_t pedrms1 = (Int_t) pedmeanrms%1000;
823 Float_t pedrms = (Float_t)pedrms1/10.;
824 Float_t pedmean = (Float_t) (pedmeanrms - pedrms1)/1000.0;
825 //printf("%f %f\n",pedmean, pedrms);
828 adc += (pedmean + 3.0*pedrms);
829 AddDigit(trno,det,smn,irow,icol,adc);
836 fPMDLoader->WriteDigits("OVERWRITE");
839 //____________________________________________________________________________
840 void AliPMDDigitizer::Exec(Option_t *option)
842 // Does the event merging and digitization
843 const char *cdeb = strstr(option,"deb");
846 AliDebug(100," *** PMD Exec is called ***");
849 Int_t ninputs = fManager->GetNinputs();
850 AliDebug(1,Form("Number of files to be processed = %d",ninputs));
853 for (Int_t i = 0; i < ninputs; i++)
855 Int_t troffset = fManager->GetMask(i);
856 MergeSDigits(i, troffset);
859 fRunLoader = AliRunLoader::GetRunLoader(fManager->GetOutputFolderName());
860 fPMD = (AliPMD*)gAlice->GetDetector("PMD");
861 fPMDLoader = fRunLoader->GetLoader("PMDLoader");
862 if (fPMDLoader == 0x0)
864 AliError("Can not find PMD or PMDLoader");
866 fPMDLoader->LoadDigits("update");
867 TTree* treeD = fPMDLoader->TreeD();
870 fPMDLoader->MakeTree("D");
871 treeD = fPMDLoader->TreeD();
873 Int_t bufsize = 16000;
874 if (!fDigits) fDigits = new TClonesArray("AliPMDdigit", 1000);
875 treeD->Branch("PMDDigit", &fDigits, bufsize);
882 for (Int_t idet = 0; idet < 2; idet++)
884 for (Int_t ism = 0; ism < fgkTotUM; ism++)
886 for (Int_t jrow = 0; jrow < fgkRow; jrow++)
888 for (Int_t kcol = 0; kcol < fgkCol; kcol++)
892 deltaE = fPRE[ism][jrow][kcol];
893 trno = fPRETrackNo[ism][jrow][kcol];
898 deltaE = fCPV[ism][jrow][kcol];
899 trno = fCPVTrackNo[ism][jrow][kcol];
907 // Gain decalibration
909 Float_t gain1 = Gain(idet,ism,jrow,kcol);
913 Int_t adcDecalib = (Int_t)(adc/gain1);
914 adc = (Float_t) adcDecalib;
920 // Pedestal Decalibration
922 fCalibPed->GetPedMeanRms(idet,ism,jrow,kcol);
923 Int_t pedrms1 = (Int_t) pedmeanrms%1000;
924 Float_t pedrms = (Float_t)pedrms1/10.;
926 (Float_t) (pedmeanrms - pedrms1)/1000.0;
927 //printf("%f %f\n",pedmean, pedrms);
930 adc += (pedmean + 3.0*pedrms);
931 AddDigit(trno,detno,ism,jrow,kcol,adc);
939 } // supermodule loop
941 fPMDLoader->WriteDigits("OVERWRITE");
942 fPMDLoader->UnloadDigits();
945 //____________________________________________________________________________
947 void AliPMDDigitizer::MergeSDigits(Int_t filenumber, Int_t troffset)
950 fRunLoader = AliRunLoader::GetRunLoader(fManager->GetInputFolderName(filenumber));
951 fPMDLoader = fRunLoader->GetLoader("PMDLoader");
952 fPMDLoader->LoadSDigits("read");
953 TTree* treeS = fPMDLoader->TreeS();
954 AliPMDsdigit *pmdsdigit;
955 TBranch *branch = treeS->GetBranch("PMDSDigit");
956 if (!fSDigits) fSDigits = new TClonesArray("AliPMDsdigit", 1000);
957 branch->SetAddress(&fSDigits);
959 Int_t itrackno, idet, ism;
962 Int_t nmodules = (Int_t) treeS->GetEntries();
963 AliDebug(1,Form("Number of Modules in the treeS = %d",nmodules));
964 AliDebug(1,Form("Track Offset = %d",troffset));
965 for (Int_t imodule = 0; imodule < nmodules; imodule++)
967 treeS->GetEntry(imodule);
968 Int_t nentries = fSDigits->GetLast();
969 AliDebug(2,Form("Number of Entries per Module = %d",nentries));
970 for (Int_t ient = 0; ient < nentries+1; ient++)
972 pmdsdigit = (AliPMDsdigit*)fSDigits->UncheckedAt(ient);
973 itrackno = pmdsdigit->GetTrackNumber();
974 idet = pmdsdigit->GetDetector();
975 ism = pmdsdigit->GetSMNumber();
976 ixp = pmdsdigit->GetRow();
977 iyp = pmdsdigit->GetColumn();
978 edep = pmdsdigit->GetCellEdep();
981 if (fPRE[ism][ixp][iyp] < edep)
983 fPRETrackNo[ism][ixp][iyp] = troffset + itrackno;
985 fPRE[ism][ixp][iyp] += edep;
989 if (fCPV[ism][ixp][iyp] < edep)
991 fCPVTrackNo[ism][ixp][iyp] = troffset + itrackno;
993 fCPV[ism][ixp][iyp] += edep;
999 // ----------------------------------------------------------------------
1000 void AliPMDDigitizer::TrackAssignment2Cell()
1003 // This block assigns the cell id when there are
1004 // multiple tracks in a cell according to the
1005 // energy deposition
1007 Bool_t jsort = false;
1017 Float_t ****pmdEdep;
1019 pmdTrack = new Int_t ***[fgkTotUM];
1020 pmdEdep = new Float_t ***[fgkTotUM];
1021 for (i=0; i<fgkTotUM; i++)
1023 pmdTrack[i] = new Int_t **[fgkRow];
1024 pmdEdep[i] = new Float_t **[fgkRow];
1027 for (i = 0; i < fgkTotUM; i++)
1029 for (j = 0; j < fgkRow; j++)
1031 pmdTrack[i][j] = new Int_t *[fgkCol];
1032 pmdEdep[i][j] = new Float_t *[fgkCol];
1036 for (i = 0; i < fgkTotUM; i++)
1038 for (j = 0; j < fgkRow; j++)
1040 for (k = 0; k < fgkCol; k++)
1042 Int_t nn = fPRECounter[i][j][k];
1045 pmdTrack[i][j][k] = new Int_t[nn];
1046 pmdEdep[i][j][k] = new Float_t[nn];
1051 pmdTrack[i][j][k] = new Int_t[nn];
1052 pmdEdep[i][j][k] = new Float_t[nn];
1054 fPRECounter[i][j][k] = 0;
1060 Int_t nentries = fCell.GetEntries();
1062 Int_t mtrackno, ism, ixp, iyp;
1065 for (i = 0; i < nentries; i++)
1067 AliPMDcell* cell = (AliPMDcell*)fCell.UncheckedAt(i);
1069 mtrackno = cell->GetTrackNumber();
1070 ism = cell->GetSMNumber();
1073 edep = cell->GetEdep();
1074 Int_t nn = fPRECounter[ism][ixp][iyp];
1075 pmdTrack[ism][ixp][iyp][nn] = (Int_t) mtrackno;
1076 pmdEdep[ism][ixp][iyp][nn] = edep;
1077 fPRECounter[ism][ixp][iyp]++;
1084 for (im=0; im<fgkTotUM; im++)
1086 for (ix=0; ix<fgkRow; ix++)
1088 for (iy=0; iy<fgkCol; iy++)
1090 nn = fPRECounter[im][ix][iy];
1093 // This block handles if a cell is fired
1094 // many times by many tracks
1095 status1 = new Int_t[nn];
1096 status2 = new Int_t[nn];
1097 trnarray = new Int_t[nn];
1098 for (iz = 0; iz < nn; iz++)
1100 status1[iz] = pmdTrack[im][ix][iy][iz];
1102 TMath::Sort(nn,status1,status2,jsort);
1103 Int_t trackOld = -99999;
1104 Int_t track, trCount = 0;
1105 for (iz = 0; iz < nn; iz++)
1107 track = status1[status2[iz]];
1108 if (trackOld != track)
1110 trnarray[trCount] = track;
1117 Float_t totEdp = 0.;
1118 trEdp = new Float_t[trCount];
1119 fracEdp = new Float_t[trCount];
1120 for (il = 0; il < trCount; il++)
1123 track = trnarray[il];
1124 for (iz = 0; iz < nn; iz++)
1126 if (track == pmdTrack[im][ix][iy][iz])
1128 trEdp[il] += pmdEdep[im][ix][iy][iz];
1131 totEdp += trEdp[il];
1134 Float_t fracOld = 0.;
1136 for (il = 0; il < trCount; il++)
1138 fracEdp[il] = trEdp[il]/totEdp;
1139 if (fracOld < fracEdp[il])
1141 fracOld = fracEdp[il];
1145 fPRETrackNo[im][ix][iy] = trnarray[ilOld];
1152 // This only handles if a cell is fired
1153 // by only one track
1155 fPRETrackNo[im][ix][iy] = pmdTrack[im][ix][iy][0];
1160 // This is if no cell is fired
1161 fPRETrackNo[im][ix][iy] = -999;
1167 // Delete all the pointers
1169 for (i = 0; i < fgkTotUM; i++)
1171 for (j = 0; j < fgkRow; j++)
1173 for (k = 0; k < fgkCol; k++)
1175 delete [] pmdTrack[i][j][k];
1176 delete [] pmdEdep[i][j][k];
1181 for (i = 0; i < fgkTotUM; i++)
1183 for (j = 0; j < fgkRow; j++)
1185 delete [] pmdTrack[i][j];
1186 delete [] pmdEdep[i][j];
1190 for (i = 0; i < fgkTotUM; i++)
1192 delete [] pmdTrack[i];
1193 delete [] pmdEdep[i];
1198 // End of the cell id assignment
1201 //____________________________________________________________________________
1202 void AliPMDDigitizer::MeV2ADC(Float_t mev, Float_t & adc) const
1204 // This converts the simulated edep to ADC according to the
1209 // PS Test in September 2003
1210 // MeV - ADC conversion for 10bit ADC
1212 const Float_t kConstant = 7.181;
1213 const Float_t kErConstant = 0.6899;
1214 const Float_t kSlope = 35.93;
1215 const Float_t kErSlope = 0.306;
1217 //gRandom->SetSeed();
1219 Float_t cons = gRandom->Gaus(kConstant,kErConstant);
1220 Float_t slop = gRandom->Gaus(kSlope,kErSlope);
1222 Float_t adc10bit = slop*mev*0.001 + cons;
1226 Int_t adc12bit = (Int_t) (4.0*adc10bit);
1230 adc = (Float_t) adc12bit;
1232 else if (adc12bit >= 3000)
1238 //____________________________________________________________________________
1239 void AliPMDDigitizer::AddSDigit(Int_t trnumber, Int_t det, Int_t smnumber,
1240 Int_t irow, Int_t icol, Float_t adc)
1244 if (!fSDigits) fSDigits = new TClonesArray("AliPMDsdigit", 1000);
1245 TClonesArray &lsdigits = *fSDigits;
1246 new(lsdigits[fNsdigit++]) AliPMDsdigit(trnumber,det,smnumber,irow,icol,adc);
1248 //____________________________________________________________________________
1250 void AliPMDDigitizer::AddDigit(Int_t trnumber, Int_t det, Int_t smnumber,
1251 Int_t irow, Int_t icol, Float_t adc)
1255 if (!fDigits) fDigits = new TClonesArray("AliPMDdigit", 1000);
1256 TClonesArray &ldigits = *fDigits;
1257 new(ldigits[fNdigit++]) AliPMDdigit(trnumber,det,smnumber,irow,icol,adc);
1259 //____________________________________________________________________________
1261 void AliPMDDigitizer::SetZPosition(Float_t zpos)
1265 //____________________________________________________________________________
1266 Float_t AliPMDDigitizer::GetZPosition() const
1270 //____________________________________________________________________________
1272 void AliPMDDigitizer::ResetCell()
1274 // clears the cell array and also the counter
1278 for (Int_t i = 0; i < fgkTotUM; i++)
1280 for (Int_t j = 0; j < fgkRow; j++)
1282 for (Int_t k = 0; k < fgkCol; k++)
1284 fPRECounter[i][j][k] = 0;
1289 //____________________________________________________________________________
1290 void AliPMDDigitizer::ResetSDigit()
1294 if (fSDigits) fSDigits->Delete();
1296 //____________________________________________________________________________
1297 void AliPMDDigitizer::ResetDigit()
1301 if (fDigits) fDigits->Delete();
1303 //____________________________________________________________________________
1305 void AliPMDDigitizer::ResetCellADC()
1307 // Clears individual cells edep and track number
1308 for (Int_t i = 0; i < fgkTotUM; i++)
1310 for (Int_t j = 0; j < fgkRow; j++)
1312 for (Int_t k = 0; k < fgkCol; k++)
1316 fPRETrackNo[i][j][k] = 0;
1317 fCPVTrackNo[i][j][k] = 0;
1322 //____________________________________________________________________________
1324 void AliPMDDigitizer::UnLoad(Option_t *option)
1326 // Unloads all the root files
1328 const char *cS = strstr(option,"S");
1329 const char *cD = strstr(option,"D");
1331 fRunLoader->UnloadgAlice();
1332 fRunLoader->UnloadHeader();
1333 fRunLoader->UnloadKinematics();
1337 fPMDLoader->UnloadHits();
1341 fPMDLoader->UnloadHits();
1342 fPMDLoader->UnloadSDigits();
1346 //----------------------------------------------------------------------
1347 Float_t AliPMDDigitizer::Gain(Int_t det, Int_t smn, Int_t row, Int_t col) const
1349 // returns of the gain of the cell
1350 // Added this method by ZA
1352 //cout<<" I am here in gain "<<fCalibData<< "smn,row, col "<<smn
1353 //<<" "<<row<<" "<<col<<endl;
1356 AliError("No calibration data loaded from CDB!!!");
1361 GainFact = fCalibGain->GetGainFact(det,smn,row,col);
1364 //----------------------------------------------------------------------
1365 AliPMDCalibData* AliPMDDigitizer::GetCalibGain() const
1367 // The run number will be centralized in AliCDBManager,
1368 // you don't need to set it here!
1369 // Added this method by ZA
1370 // Cleaned up by Alberto
1371 AliCDBEntry *entry = AliCDBManager::Instance()->Get("PMD/Calib/Gain");
1373 if(!entry) AliFatal("Calibration object retrieval failed!");
1375 AliPMDCalibData *calibdata=0;
1376 if (entry) calibdata = (AliPMDCalibData*) entry->GetObject();
1378 if (!calibdata) AliFatal("No calibration data from calibration database !");
1382 //----------------------------------------------------------------------
1383 AliPMDPedestal* AliPMDDigitizer::GetCalibPed() const
1385 // The run number will be centralized in AliCDBManager,
1386 // you don't need to set it here!
1388 AliCDBEntry *entry = AliCDBManager::Instance()->Get("PMD/Calib/Ped");
1390 if(!entry) AliFatal("Pedestal object retrieval failed!");
1392 AliPMDPedestal *pedestal=0;
1393 if (entry) pedestal = (AliPMDPedestal*) entry->GetObject();
1395 if (!pedestal) AliFatal("No pedestal data from calibration database !");