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 "AliDigitizationInput.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()),
77 fZPos(361.5) // in units of cm, default position of PMD
79 // Default Constructor
81 for (Int_t i = 0; i < fgkTotUM; i++)
83 for (Int_t j = 0; j < fgkRow; j++)
85 for (Int_t k = 0; k < fgkCol; k++)
89 fCPVCounter[i][j][k] = 0;
90 fPRECounter[i][j][k] = 0;
91 fCPVTrackNo[i][j][k] = -1;
92 fPRETrackNo[i][j][k] = -1;
93 fCPVTrackPid[i][j][k] = -1;
94 fPRETrackPid[i][j][k] = -1;
100 //____________________________________________________________________________
101 AliPMDDigitizer::AliPMDDigitizer(const AliPMDDigitizer& digitizer):
102 AliDigitizer(digitizer),
107 fCalibGain(GetCalibGain()),
108 fCalibPed(GetCalibPed()),
116 fZPos(361.5) // in units of cm, default position of PMD
119 AliError("Copy constructor not allowed ");
122 //____________________________________________________________________________
123 AliPMDDigitizer & AliPMDDigitizer::operator=(const AliPMDDigitizer& /*digitizer*/)
125 // Assignment operator
126 AliError("Assignement operator not allowed ");
130 //____________________________________________________________________________
131 AliPMDDigitizer::AliPMDDigitizer(AliDigitizationInput* digInput):
132 AliDigitizer(digInput),
137 fCalibGain(GetCalibGain()),
138 fCalibPed(GetCalibPed()),
139 fSDigits(new TClonesArray("AliPMDsdigit", 1000)),
140 fDigits(new TClonesArray("AliPMDdigit", 1000)),
146 fZPos(361.5)// in units of cm, This is the default position of PMD
148 // ctor which should be used
151 for (Int_t i = 0; i < fgkTotUM; i++)
153 for (Int_t j = 0; j < fgkRow; j++)
155 for (Int_t k = 0; k < fgkCol; k++)
159 fCPVCounter[i][j][k] = 0;
160 fPRECounter[i][j][k] = 0;
161 fCPVTrackNo[i][j][k] = -1;
162 fPRETrackNo[i][j][k] = -1;
163 fCPVTrackPid[i][j][k] = -1;
164 fPRETrackPid[i][j][k] = -1;
170 //____________________________________________________________________________
171 AliPMDDigitizer::~AliPMDDigitizer()
173 // Default Destructor
191 //____________________________________________________________________________
192 void AliPMDDigitizer::OpengAliceFile(const char *file, Option_t *option)
194 // Loads galice.root file and corresponding header, kinematics
195 // hits and sdigits or digits depending on the option
198 TString evfoldname = AliConfig::GetDefaultEventFolderName();
199 fRunLoader = AliRunLoader::GetRunLoader(evfoldname);
201 fRunLoader = AliRunLoader::Open(file,AliConfig::GetDefaultEventFolderName(), "UPDATE");
205 AliError(Form("Can not open session for file %s.",file));
208 const char *cHS = strstr(option,"HS");
209 const char *cHD = strstr(option,"HD");
210 const char *cSD = strstr(option,"SD");
214 if (!fRunLoader->GetAliRun()) fRunLoader->LoadgAlice();
215 if (!fRunLoader->TreeE()) fRunLoader->LoadHeader();
216 if (!fRunLoader->TreeK()) fRunLoader->LoadKinematics();
218 gAlice = fRunLoader->GetAliRun();
222 AliDebug(1,"Alirun object found");
226 AliError("Could not found Alirun object");
229 fPMD = (AliPMD*)gAlice->GetDetector("PMD");
232 fPMDLoader = fRunLoader->GetLoader("PMDLoader");
233 if (fPMDLoader == 0x0)
235 AliError("Can not find PMDLoader");
241 fPMDLoader->LoadHits("READ");
242 fPMDLoader->LoadSDigits("recreate");
246 fPMDLoader->LoadHits("READ");
247 fPMDLoader->LoadDigits("recreate");
251 fPMDLoader->LoadSDigits("READ");
252 fPMDLoader->LoadDigits("recreate");
255 //____________________________________________________________________________
256 void AliPMDDigitizer::Hits2SDigits(Int_t ievt)
258 // This reads the PMD Hits tree and assigns the right track number
259 // to a cell and stores in the summable digits tree
262 const Int_t kPi0 = 111;
263 const Int_t kGamma = 22;
271 Float_t xPos = 0., yPos = 0., zPos = 0.;
272 Int_t xpad = -1, ypad = -1;
274 Float_t vx = -999.0, vy = -999.0, vz = -999.0;
277 if (!fSDigits) fSDigits = new TClonesArray("AliPMDsdigit", 1000);
280 AliDebug(1,Form("Event Number = %d",ievt));
281 Int_t nparticles = fRunLoader->GetHeader()->GetNtrack();
282 AliDebug(1,Form("Number of Particles = %d",nparticles));
286 fRunLoader->GetEvent(ievt);
287 // ------------------------------------------------------- //
288 // Pointer to specific detector hits.
289 // Get pointers to Alice detectors and Hits containers
291 TTree* treeH = fPMDLoader->TreeH();
293 Int_t ntracks = (Int_t) treeH->GetEntries();
294 AliDebug(1,Form("Number of Tracks in the TreeH = %d", ntracks));
295 TTree* treeS = fPMDLoader->TreeS();
298 fPMDLoader->MakeTree("S");
299 treeS = fPMDLoader->TreeS();
301 Int_t bufsize = 16000;
302 treeS->Branch("PMDSDigit", &fSDigits, bufsize);
304 TClonesArray* hits = 0;
305 if (fPMD) hits = fPMD->Hits();
307 // Start loop on tracks in the hits containers
309 for (Int_t track=0; track<ntracks;track++)
311 gAlice->GetMCApp()->ResetHits();
312 treeH->GetEvent(track);
315 npmd = hits->GetEntriesFast();
316 for (Int_t ipmd = 0; ipmd < npmd; ipmd++)
318 fPMDHit = (AliPMDhit*) hits->UncheckedAt(ipmd);
319 trackno = fPMDHit->GetTrack();
320 // get kinematics of the particles
322 TParticle* mparticle = gAlice->GetMCApp()->Particle(trackno);
323 trackpid = mparticle->GetPdgCode();
324 Int_t ks = mparticle->GetStatusCode();
326 Int_t tracknoOld=0, trackpidOld=0, statusOld = 0;
328 if (mparticle->GetFirstMother() == -1)
330 tracknoOld = trackno;
331 trackpidOld = trackpid;
336 Int_t trnotemp = trackno; // Modified on 25th Nov 2009
337 if(ks==1||(imo = mparticle->GetFirstMother())<0 ){
338 vx = mparticle->Vx();
339 vy = mparticle->Vy();
340 vz = mparticle->Vz();
342 if(trackpid==kGamma||trackpid==11||trackpid==-11||
343 trackpid==kPi0)igstatus=1;
347 while(((imo = mparticle->GetFirstMother()) >= 0)&&
348 (ks = mparticle->GetStatusCode() <1) )
350 mparticle = gAlice->GetMCApp()->Particle(imo);
351 trackpid = mparticle->GetPdgCode();
352 ks = mparticle->GetStatusCode();
353 vx = mparticle->Vx();
354 vy = mparticle->Vy();
355 vz = mparticle->Vz();
357 // Modified on 25th Nov 2009
368 // end of modification on 25th Nov 2009
371 if(trackpid==kGamma||trackpid==11||trackpid==-11||
372 trackpid==kPi0)igstatus=1;
375 trackpid=trackpidOld;
378 //-----------------end of modification----------------
379 Float_t ptime = fPMDHit->GetTime()*1e6; // time in microsec
380 if (ptime < 0. || ptime > 1.2) continue;
386 edep = fPMDHit->GetEnergy();
387 Int_t vol1 = fPMDHit->GetVolume(1); // Column
388 Int_t vol2 = fPMDHit->GetVolume(2); // Row
389 Int_t vol7 = fPMDHit->GetVolume(4); // Serial Module No
392 // -----------------------------------------//
393 // In new geometry after adding electronics //
394 // For Super Module 1 & 2 //
395 // nrow = 48, ncol = 96 //
396 // For Super Module 3 & 4 //
397 // nrow = 96, ncol = 48 //
398 // -----------------------------------------//
406 smnumber = vol7 - 24;
408 Int_t vol8 = smnumber/6 + 1; // fake supermodule
410 if (vol8 == 1 || vol8 == 2)
415 else if (vol8 == 3 || vol8 == 4)
421 AliDebug(2,Form("Zposition = %f Edeposition = %f",zPos,edep));
434 Int_t smn = smnumber;
435 Int_t ixx = xpad - 1;
436 Int_t iyy = ypad - 1;
439 fPRE[smn][ixx][iyy] += edep;
440 fPRECounter[smn][ixx][iyy]++;
442 AliPMDcell* cell = new AliPMDcell(mtrackno,smn,ixx,iyy,edep);
447 fCPV[smn][ixx][iyy] += edep;
448 fCPVCounter[smn][ixx][iyy]++;
449 AliPMDcell* cpvcell = new AliPMDcell(mtrackno,smn,ixx,iyy,edep);
450 fCPVCell.Add(cpvcell);
454 } // Track Loop ended
455 TrackAssignment2CPVCell();
456 TrackAssignment2Cell();
464 for (Int_t idet = 0; idet < 2; idet++)
466 for (Int_t ism = 0; ism < fgkTotUM; ism++)
468 for (Int_t jrow = 0; jrow < fgkRow; jrow++)
470 for (Int_t kcol = 0; kcol < fgkCol; kcol++)
474 deltaE = fPRE[ism][jrow][kcol];
475 trno = fPRETrackNo[ism][jrow][kcol];
480 deltaE = fCPV[ism][jrow][kcol];
481 trno = fCPVTrackNo[ism][jrow][kcol];
487 TParticle *mparticle = gAlice->GetMCApp()->Particle(trno);
488 trpid = mparticle->GetPdgCode();
489 AddSDigit(trno,trpid,detno,ism,jrow,kcol,deltaE);
497 fPMDLoader->WriteSDigits("OVERWRITE");
500 //____________________________________________________________________________
502 void AliPMDDigitizer::Hits2Digits(Int_t ievt)
504 // This reads the PMD Hits tree and assigns the right track number
505 // to a cell and stores in the digits tree
507 const Int_t kPi0 = 111;
508 const Int_t kGamma = 22;
516 Float_t xPos = 0., yPos = 0., zPos = 0.;
517 Int_t xpad = -1, ypad = -1;
519 Float_t vx = -999.0, vy = -999.0, vz = -999.0;
521 if (!fDigits) fDigits = new TClonesArray("AliPMDdigit", 1000);
524 AliDebug(1,Form("Event Number = %d",ievt));
525 Int_t nparticles = fRunLoader->GetHeader()->GetNtrack();
526 AliDebug(1,Form("Number of Particles = %d", nparticles));
528 fRunLoader->GetEvent(ievt);
529 // ------------------------------------------------------- //
530 // Pointer to specific detector hits.
531 // Get pointers to Alice detectors and Hits containers
533 fPMD = (AliPMD*)gAlice->GetDetector("PMD");
534 fPMDLoader = fRunLoader->GetLoader("PMDLoader");
536 if (fPMDLoader == 0x0)
538 AliError("Can not find PMD or PMDLoader");
540 TTree* treeH = fPMDLoader->TreeH();
541 Int_t ntracks = (Int_t) treeH->GetEntries();
542 AliDebug(1,Form("Number of Tracks in the TreeH = %d", ntracks));
543 fPMDLoader->LoadDigits("recreate");
544 TTree* treeD = fPMDLoader->TreeD();
547 fPMDLoader->MakeTree("D");
548 treeD = fPMDLoader->TreeD();
550 Int_t bufsize = 16000;
551 treeD->Branch("PMDDigit", &fDigits, bufsize);
553 TClonesArray* hits = 0;
554 if (fPMD) hits = fPMD->Hits();
556 // Start loop on tracks in the hits containers
558 for (Int_t track=0; track<ntracks;track++)
560 gAlice->GetMCApp()->ResetHits();
561 treeH->GetEvent(track);
565 npmd = hits->GetEntriesFast();
566 for (Int_t ipmd = 0; ipmd < npmd; ipmd++)
568 fPMDHit = (AliPMDhit*) hits->UncheckedAt(ipmd);
569 trackno = fPMDHit->GetTrack();
571 // get kinematics of the particles
573 TParticle* mparticle = gAlice->GetMCApp()->Particle(trackno);
574 trackpid = mparticle->GetPdgCode();
575 Int_t ks = mparticle->GetStatusCode();
577 Int_t tracknoOld=0, trackpidOld=0, statusOld = 0;
578 if (mparticle->GetFirstMother() == -1)
580 tracknoOld = trackno;
581 trackpidOld = trackpid;
587 Int_t trnotemp = trackno; // modified on 25th Nov 2009
588 if(ks==1||(imo = mparticle->GetFirstMother())<0 ){
589 vx = mparticle->Vx();
590 vy = mparticle->Vy();
591 vz = mparticle->Vz();
593 if(trackpid==kGamma||trackpid==11||trackpid==-11||trackpid==kPi0)
598 while(((imo = mparticle->GetFirstMother()) >= 0)&&
599 (ks = mparticle->GetStatusCode() <1) )
601 mparticle = gAlice->GetMCApp()->Particle(imo);
602 trackpid = mparticle->GetPdgCode();
603 ks = mparticle->GetStatusCode();
604 vx = mparticle->Vx();
605 vy = mparticle->Vy();
606 vz = mparticle->Vz();
608 // Modified on 25th Nov 2009
621 if(trackpid==kGamma||trackpid==11||trackpid==-11||trackpid==kPi0)
625 trackpid=trackpidOld;
628 Float_t ptime = fPMDHit->GetTime()*1e6;
629 if (ptime < 0. || ptime > 1.2) continue;
634 edep = fPMDHit->GetEnergy();
635 Int_t vol1 = fPMDHit->GetVolume(1); // Column
636 Int_t vol2 = fPMDHit->GetVolume(2); // Row
637 Int_t vol7 = fPMDHit->GetVolume(4); // Serial Module No
639 // -----------------------------------------//
640 // In new geometry after adding electronics //
641 // For Super Module 1 & 2 //
642 // nrow = 48, ncol = 96 //
643 // For Super Module 3 & 4 //
644 // nrow = 96, ncol = 48 //
645 // -----------------------------------------//
653 smnumber = vol7 - 24;
655 Int_t vol8 = smnumber/6 + 1; // fake supermodule
657 if (vol8 == 1 || vol8 == 2)
662 else if (vol8 == 3 || vol8 == 4)
668 AliDebug(2,Form("ZPosition = %f Edeposition = %f",zPos,edep));
680 Int_t smn = smnumber;
681 Int_t ixx = xpad - 1;
682 Int_t iyy = ypad - 1;
685 fPRE[smn][ixx][iyy] += edep;
686 fPRECounter[smn][ixx][iyy]++;
688 AliPMDcell* cell = new AliPMDcell(mtrackno,smn,ixx,iyy,edep);
694 fCPV[smn][ixx][iyy] += edep;
695 fCPVCounter[smn][ixx][iyy]++;
696 AliPMDcell* cpvcell = new AliPMDcell(mtrackno,smn,ixx,iyy,edep);
697 fCPVCell.Add(cpvcell);
701 } // Track Loop ended
702 TrackAssignment2CPVCell();
703 TrackAssignment2Cell();
713 for (Int_t idet = 0; idet < 2; idet++)
715 for (Int_t ism = 0; ism < fgkTotUM; ism++)
717 for (Int_t jrow = 0; jrow < fgkRow; jrow++)
719 for (Int_t kcol = 0; kcol < fgkCol; kcol++)
723 deltaE = fPRE[ism][jrow][kcol];
724 trno = fPRETrackNo[ism][jrow][kcol];
729 deltaE = fCPV[ism][jrow][kcol];
730 trno = fCPVTrackNo[ism][jrow][kcol];
737 // To decalibrate the adc values
739 gain1 = Gain(idet,ism,jrow,kcol);
742 Int_t adcDecalib = (Int_t)(adc/gain1);
743 adc = (Float_t) adcDecalib;
750 // Pedestal Decalibration
752 fCalibPed->GetPedMeanRms(idet,ism,jrow,kcol);
753 Int_t pedrms1 = (Int_t) pedmeanrms%100;
754 Float_t pedrms = (Float_t)pedrms1/10.;
756 (Float_t) (pedmeanrms - pedrms1)/1000.0;
759 adc += (pedmean + 3.0*pedrms);
761 = gAlice->GetMCApp()->Particle(trno);
762 trpid = mparticle->GetPdgCode();
764 AddDigit(trno,trpid,detno,ism,jrow,kcol,adc);
771 } // supermodule loop
774 fPMDLoader->WriteDigits("OVERWRITE");
778 //____________________________________________________________________________
781 void AliPMDDigitizer::SDigits2Digits(Int_t ievt)
783 // This reads the PMD sdigits tree and converts energy deposition
784 // in a cell to ADC and stores in the digits tree
787 fRunLoader->GetEvent(ievt);
789 TTree* treeS = fPMDLoader->TreeS();
790 AliPMDsdigit *pmdsdigit;
791 TBranch *branch = treeS->GetBranch("PMDSDigit");
794 AliError("PMD Sdigit branch does not exist");
797 if (!fSDigits) fSDigits = new TClonesArray("AliPMDsdigit", 1000);
798 branch->SetAddress(&fSDigits);
800 TTree* treeD = fPMDLoader->TreeD();
803 fPMDLoader->MakeTree("D");
804 treeD = fPMDLoader->TreeD();
806 Int_t bufsize = 16000;
807 if (!fDigits) fDigits = new TClonesArray("AliPMDdigit", 1000);
808 treeD->Branch("PMDDigit", &fDigits, bufsize);
810 Int_t trno = 1, trpid = 0, det = 0, smn = 0;
811 Int_t irow = 0, icol = 0;
812 Float_t edep = 0., adc = 0.;
814 Int_t nmodules = (Int_t) treeS->GetEntries();
815 AliDebug(1,Form("Number of modules = %d",nmodules));
817 for (Int_t imodule = 0; imodule < nmodules; imodule++)
819 treeS->GetEntry(imodule);
820 Int_t nentries = fSDigits->GetLast();
821 AliDebug(2,Form("Number of entries per module = %d",nentries+1));
822 for (Int_t ient = 0; ient < nentries+1; ient++)
824 pmdsdigit = (AliPMDsdigit*)fSDigits->UncheckedAt(ient);
825 trno = pmdsdigit->GetTrackNumber();
826 trpid = pmdsdigit->GetTrackPid();
827 det = pmdsdigit->GetDetector();
828 smn = pmdsdigit->GetSMNumber();
829 irow = pmdsdigit->GetRow();
830 icol = pmdsdigit->GetColumn();
831 edep = pmdsdigit->GetCellEdep();
835 // To decalibrte the adc values
837 Float_t gain1 = Gain(det,smn,irow,icol);
840 Int_t adcDecalib = (Int_t)(adc/gain1);
841 adc = (Float_t) adcDecalib;
847 // Pedestal Decalibration
848 Int_t pedmeanrms = fCalibPed->GetPedMeanRms(det,smn,irow,icol);
849 Int_t pedrms1 = (Int_t) pedmeanrms%100;
850 Float_t pedrms = (Float_t)pedrms1/10.;
851 Float_t pedmean = (Float_t) (pedmeanrms - pedrms1)/1000.0;
854 adc += (pedmean + 3.0*pedrms);
855 AddDigit(trno,trpid,det,smn,irow,icol,adc);
862 fPMDLoader->WriteDigits("OVERWRITE");
865 //____________________________________________________________________________
866 void AliPMDDigitizer::Digitize(Option_t *option)
868 // Does the event merging and digitization
869 const char *cdeb = strstr(option,"deb");
872 AliDebug(100," *** PMD Exec is called ***");
875 Int_t ninputs = fDigInput->GetNinputs();
876 AliDebug(1,Form("Number of files to be processed = %d",ninputs));
879 for (Int_t i = 0; i < ninputs; i++)
881 Int_t troffset = fDigInput->GetMask(i);
882 MergeSDigits(i, troffset);
885 fRunLoader = AliRunLoader::GetRunLoader(fDigInput->GetOutputFolderName());
886 fPMD = (AliPMD*)gAlice->GetDetector("PMD");
887 fPMDLoader = fRunLoader->GetLoader("PMDLoader");
888 if (fPMDLoader == 0x0)
890 AliError("Can not find PMD or PMDLoader");
892 fPMDLoader->LoadDigits("update");
893 TTree* treeD = fPMDLoader->TreeD();
896 fPMDLoader->MakeTree("D");
897 treeD = fPMDLoader->TreeD();
899 Int_t bufsize = 16000;
900 if (!fDigits) fDigits = new TClonesArray("AliPMDdigit", 1000);
901 treeD->Branch("PMDDigit", &fDigits, bufsize);
909 for (Int_t idet = 0; idet < 2; idet++)
911 for (Int_t ism = 0; ism < fgkTotUM; ism++)
913 for (Int_t jrow = 0; jrow < fgkRow; jrow++)
915 for (Int_t kcol = 0; kcol < fgkCol; kcol++)
919 deltaE = fPRE[ism][jrow][kcol];
920 trno = fPRETrackNo[ism][jrow][kcol];
921 trpid = fPRETrackPid[ism][jrow][kcol];
926 deltaE = fCPV[ism][jrow][kcol];
927 trno = fCPVTrackNo[ism][jrow][kcol];
928 trpid = fCPVTrackPid[ism][jrow][kcol];
936 // Gain decalibration
938 Float_t gain1 = Gain(idet,ism,jrow,kcol);
942 Int_t adcDecalib = (Int_t)(adc/gain1);
943 adc = (Float_t) adcDecalib;
949 // Pedestal Decalibration
951 fCalibPed->GetPedMeanRms(idet,ism,jrow,kcol);
952 Int_t pedrms1 = (Int_t) pedmeanrms%100;
953 Float_t pedrms = (Float_t)pedrms1/10.;
955 (Float_t) (pedmeanrms - pedrms1)/1000.0;
958 adc += (pedmean + 3.0*pedrms);
959 AddDigit(trno,trpid,detno,ism,jrow,kcol,adc);
967 } // supermodule loop
969 fPMDLoader->WriteDigits("OVERWRITE");
970 fPMDLoader->UnloadDigits();
973 //____________________________________________________________________________
974 void AliPMDDigitizer::TrackAssignment2CPVCell()
976 // This block assigns the cell id when there are
977 // multiple tracks in a cell according to the
979 // This method added by Ajay
980 Bool_t jsort = false;
982 Int_t i = 0, j = 0, k = 0;
993 cpvTrack = new Int_t ***[fgkTotUM];
994 cpvEdep = new Float_t ***[fgkTotUM];
995 for (i=0; i<fgkTotUM; i++)
997 cpvTrack[i] = new Int_t **[fgkRow];
998 cpvEdep[i] = new Float_t **[fgkRow];
1001 for (i = 0; i < fgkTotUM; i++)
1003 for (j = 0; j < fgkRow; j++)
1005 cpvTrack[i][j] = new Int_t *[fgkCol];
1006 cpvEdep[i][j] = new Float_t *[fgkCol];
1009 for (i = 0; i < fgkTotUM; i++)
1011 for (j = 0; j < fgkRow; j++)
1013 for (k = 0; k < fgkCol; k++)
1015 Int_t nn = fCPVCounter[i][j][k];
1018 cpvTrack[i][j][k] = new Int_t[nn];
1019 cpvEdep[i][j][k] = new Float_t[nn];
1024 cpvTrack[i][j][k] = new Int_t[nn];
1025 cpvEdep[i][j][k] = new Float_t[nn];
1027 fCPVCounter[i][j][k] = 0;
1033 Int_t nentries = fCPVCell.GetEntries();
1035 Int_t mtrackno = 0, ism = 0, ixp = 0, iyp = 0;
1037 for (i = 0; i < nentries; i++)
1039 AliPMDcell* cpvcell = (AliPMDcell*)fCPVCell.UncheckedAt(i);
1041 mtrackno = cpvcell->GetTrackNumber();
1042 ism = cpvcell->GetSMNumber();
1043 ixp = cpvcell->GetX();
1044 iyp = cpvcell->GetY();
1045 edep = cpvcell->GetEdep();
1046 Int_t nn = fCPVCounter[ism][ixp][iyp];
1047 cpvTrack[ism][ixp][iyp][nn] = (Int_t) mtrackno;
1048 cpvEdep[ism][ixp][iyp][nn] = edep;
1049 fCPVCounter[ism][ixp][iyp]++;
1052 Int_t iz = 0, il = 0;
1053 Int_t im = 0, ix = 0, iy = 0;
1055 for (im=0; im<fgkTotUM; im++)
1057 for (ix=0; ix<fgkRow; ix++)
1059 for (iy=0; iy<fgkCol; iy++)
1061 nn = fCPVCounter[im][ix][iy];
1064 // This block handles if a cell is fired
1065 // many times by many tracks
1066 status1 = new Int_t[nn];
1067 status2 = new Int_t[2*nn];
1068 trnarray = new Int_t[nn];
1069 for (iz = 0; iz < nn; iz++)
1071 status1[iz] = cpvTrack[im][ix][iy][iz];
1073 TMath::Sort(nn,status1,status2,jsort);
1074 Int_t trackOld = -99999;
1075 Int_t track, trCount = 0;
1076 for (iz = 0; iz < nn; iz++)
1078 track = status1[status2[iz]];
1079 if (trackOld != track)
1081 trnarray[trCount] = track;
1088 Float_t totEdp = 0.;
1089 trEdp = new Float_t[trCount];
1090 fracEdp = new Float_t[trCount];
1091 for (il = 0; il < trCount; il++)
1094 track = trnarray[il];
1095 for (iz = 0; iz < nn; iz++)
1097 if (track == cpvTrack[im][ix][iy][iz])
1099 trEdp[il] += cpvEdep[im][ix][iy][iz];
1102 totEdp += trEdp[il];
1105 Float_t fracOld = 0.;
1107 for (il = 0; il < trCount; il++)
1109 fracEdp[il] = trEdp[il]/totEdp;
1110 if (fracOld < fracEdp[il])
1112 fracOld = fracEdp[il];
1116 fCPVTrackNo[im][ix][iy] = trnarray[ilOld];
1123 // This only handles if a cell is fired
1124 // by only one track
1126 fCPVTrackNo[im][ix][iy] = cpvTrack[im][ix][iy][0];
1131 // This is if no cell is fired
1132 fCPVTrackNo[im][ix][iy] = -999;
1138 // Delete all the pointers
1140 for (i = 0; i < fgkTotUM; i++)
1142 for (j = 0; j < fgkRow; j++)
1144 for (k = 0; k < fgkCol; k++)
1146 delete []cpvTrack[i][j][k];
1147 delete []cpvEdep[i][j][k];
1152 for (i = 0; i < fgkTotUM; i++)
1154 for (j = 0; j < fgkRow; j++)
1156 delete [] cpvTrack[i][j];
1157 delete [] cpvEdep[i][j];
1161 for (i = 0; i < fgkTotUM; i++)
1163 delete [] cpvTrack[i];
1164 delete [] cpvEdep[i];
1170 // End of the cell id assignment
1173 //____________________________________________________________________________
1175 void AliPMDDigitizer::MergeSDigits(Int_t filenumber, Int_t troffset)
1178 fRunLoader = AliRunLoader::GetRunLoader(fDigInput->GetInputFolderName(filenumber));
1179 fPMDLoader = fRunLoader->GetLoader("PMDLoader");
1180 fPMDLoader->LoadSDigits("read");
1181 TTree* treeS = fPMDLoader->TreeS();
1182 AliPMDsdigit *pmdsdigit;
1183 TBranch *branch = treeS->GetBranch("PMDSDigit");
1184 if (!fSDigits) fSDigits = new TClonesArray("AliPMDsdigit", 1000);
1185 branch->SetAddress(&fSDigits);
1187 Int_t itrackno = 1, itrackpid = 0, idet = 0, ism = 0;
1188 Int_t ixp = 0, iyp = 0;
1190 Int_t nmodules = (Int_t) treeS->GetEntries();
1191 AliDebug(1,Form("Number of Modules in the treeS = %d",nmodules));
1192 AliDebug(1,Form("Track Offset = %d",troffset));
1193 for (Int_t imodule = 0; imodule < nmodules; imodule++)
1195 treeS->GetEntry(imodule);
1196 Int_t nentries = fSDigits->GetLast();
1197 AliDebug(2,Form("Number of Entries per Module = %d",nentries));
1198 for (Int_t ient = 0; ient < nentries+1; ient++)
1200 pmdsdigit = (AliPMDsdigit*)fSDigits->UncheckedAt(ient);
1201 itrackno = pmdsdigit->GetTrackNumber();
1202 itrackpid = pmdsdigit->GetTrackPid();
1203 idet = pmdsdigit->GetDetector();
1204 ism = pmdsdigit->GetSMNumber();
1205 ixp = pmdsdigit->GetRow();
1206 iyp = pmdsdigit->GetColumn();
1207 edep = pmdsdigit->GetCellEdep();
1210 if (fPRE[ism][ixp][iyp] < edep)
1212 fPRETrackNo[ism][ixp][iyp] = troffset + itrackno;
1213 fPRETrackPid[ism][ixp][iyp] = itrackpid;
1215 fPRE[ism][ixp][iyp] += edep;
1219 if (fCPV[ism][ixp][iyp] < edep)
1221 fCPVTrackNo[ism][ixp][iyp] = troffset + itrackno;
1222 fCPVTrackPid[ism][ixp][iyp] = itrackpid;
1224 fCPV[ism][ixp][iyp] += edep;
1230 // ----------------------------------------------------------------------
1231 void AliPMDDigitizer::TrackAssignment2Cell()
1234 // This block assigns the cell id when there are
1235 // multiple tracks in a cell according to the
1236 // energy deposition
1238 Bool_t jsort = false;
1240 Int_t i = 0, j = 0, k = 0;
1249 Float_t ****pmdEdep;
1251 pmdTrack = new Int_t ***[fgkTotUM];
1252 pmdEdep = new Float_t ***[fgkTotUM];
1253 for (i=0; i<fgkTotUM; i++)
1255 pmdTrack[i] = new Int_t **[fgkRow];
1256 pmdEdep[i] = new Float_t **[fgkRow];
1259 for (i = 0; i < fgkTotUM; i++)
1261 for (j = 0; j < fgkRow; j++)
1263 pmdTrack[i][j] = new Int_t *[fgkCol];
1264 pmdEdep[i][j] = new Float_t *[fgkCol];
1268 for (i = 0; i < fgkTotUM; i++)
1270 for (j = 0; j < fgkRow; j++)
1272 for (k = 0; k < fgkCol; k++)
1274 Int_t nn = fPRECounter[i][j][k];
1277 pmdTrack[i][j][k] = new Int_t[nn];
1278 pmdEdep[i][j][k] = new Float_t[nn];
1283 pmdTrack[i][j][k] = new Int_t[nn];
1284 pmdEdep[i][j][k] = new Float_t[nn];
1286 fPRECounter[i][j][k] = 0;
1292 Int_t nentries = fCell.GetEntries();
1294 Int_t mtrackno, ism, ixp, iyp;
1297 for (i = 0; i < nentries; i++)
1299 AliPMDcell* cell = (AliPMDcell*)fCell.UncheckedAt(i);
1301 mtrackno = cell->GetTrackNumber();
1302 ism = cell->GetSMNumber();
1305 edep = cell->GetEdep();
1306 Int_t nn = fPRECounter[ism][ixp][iyp];
1307 pmdTrack[ism][ixp][iyp][nn] = (Int_t) mtrackno;
1308 pmdEdep[ism][ixp][iyp][nn] = edep;
1309 fPRECounter[ism][ixp][iyp]++;
1312 Int_t iz = 0, il = 0;
1313 Int_t im = 0, ix = 0, iy = 0;
1316 for (im=0; im<fgkTotUM; im++)
1318 for (ix=0; ix<fgkRow; ix++)
1320 for (iy=0; iy<fgkCol; iy++)
1322 nn = fPRECounter[im][ix][iy];
1325 // This block handles if a cell is fired
1326 // many times by many tracks
1327 status1 = new Int_t[nn];
1328 status2 = new Int_t[2*nn];
1329 trnarray = new Int_t[nn];
1330 for (iz = 0; iz < nn; iz++)
1332 status1[iz] = pmdTrack[im][ix][iy][iz];
1334 TMath::Sort(nn,status1,status2,jsort);
1335 Int_t trackOld = -99999;
1336 Int_t track, trCount = 0;
1337 for (iz = 0; iz < nn; iz++)
1339 track = status1[status2[iz]];
1340 if (trackOld != track)
1342 trnarray[trCount] = track;
1349 Float_t totEdp = 0.;
1350 trEdp = new Float_t[trCount];
1351 fracEdp = new Float_t[trCount];
1352 for (il = 0; il < trCount; il++)
1355 track = trnarray[il];
1356 for (iz = 0; iz < nn; iz++)
1358 if (track == pmdTrack[im][ix][iy][iz])
1360 trEdp[il] += pmdEdep[im][ix][iy][iz];
1363 totEdp += trEdp[il];
1366 Float_t fracOld = 0.;
1368 for (il = 0; il < trCount; il++)
1370 fracEdp[il] = trEdp[il]/totEdp;
1371 if (fracOld < fracEdp[il])
1373 fracOld = fracEdp[il];
1377 fPRETrackNo[im][ix][iy] = trnarray[ilOld];
1384 // This only handles if a cell is fired
1385 // by only one track
1387 fPRETrackNo[im][ix][iy] = pmdTrack[im][ix][iy][0];
1392 // This is if no cell is fired
1393 fPRETrackNo[im][ix][iy] = -999;
1399 // Delete all the pointers
1401 for (i = 0; i < fgkTotUM; i++)
1403 for (j = 0; j < fgkRow; j++)
1405 for (k = 0; k < fgkCol; k++)
1407 delete [] pmdTrack[i][j][k];
1408 delete [] pmdEdep[i][j][k];
1413 for (i = 0; i < fgkTotUM; i++)
1415 for (j = 0; j < fgkRow; j++)
1417 delete [] pmdTrack[i][j];
1418 delete [] pmdEdep[i][j];
1422 for (i = 0; i < fgkTotUM; i++)
1424 delete [] pmdTrack[i];
1425 delete [] pmdEdep[i];
1430 // End of the cell id assignment
1433 //____________________________________________________________________________
1434 void AliPMDDigitizer::MeV2ADC(Float_t mev, Float_t & adc) const
1436 // This converts the simulated edep to ADC according to the
1438 // PS Test in June 2010, Voltage @ 1300 V
1439 // KeV - ADC conversion for 12bit ADC
1440 // MPV data used for the fit and taken here
1442 // constants are from Test Beam 2010
1444 const Float_t kConstant = 0.612796;
1445 const Float_t kSlope = 130.158;
1447 Float_t adc12bit = kSlope*mev*0.001 + kConstant;
1448 if (adc12bit < 0.) adc12bit = 0.;
1450 //Introducing Readout Resolution for ALICE-PMD
1452 Float_t sigrr = 0.605016 - 0.000273*adc12bit + 6.54e-8*adc12bit*adc12bit;
1453 Float_t adcwithrr = gRandom->Gaus(adc12bit,sigrr);
1459 else if(adcwithrr >= 0. && adcwithrr < 1600.0)
1463 else if (adcwithrr >= 1600.0)
1469 //____________________________________________________________________________
1470 void AliPMDDigitizer::AddSDigit(Int_t trnumber, Int_t trpid, Int_t det,
1471 Int_t smnumber, Int_t irow, Int_t icol,
1476 if (!fSDigits) fSDigits = new TClonesArray("AliPMDsdigit", 1000);
1477 TClonesArray &lsdigits = *fSDigits;
1478 new(lsdigits[fNsdigit++]) AliPMDsdigit(trnumber,trpid,det,smnumber,irow,icol,adc);
1480 //____________________________________________________________________________
1482 void AliPMDDigitizer::AddDigit(Int_t trnumber, Int_t trpid, Int_t det,
1483 Int_t smnumber, Int_t irow, Int_t icol,
1488 if (!fDigits) fDigits = new TClonesArray("AliPMDdigit", 1000);
1489 TClonesArray &ldigits = *fDigits;
1490 new(ldigits[fNdigit++]) AliPMDdigit(trnumber,trpid, det,smnumber,irow,icol,adc);
1492 //____________________________________________________________________________
1494 void AliPMDDigitizer::SetZPosition(Float_t zpos)
1498 //____________________________________________________________________________
1499 Float_t AliPMDDigitizer::GetZPosition() const
1503 //____________________________________________________________________________
1505 void AliPMDDigitizer::ResetCell()
1507 // clears the cell array and also the counter
1512 for (Int_t i = 0; i < fgkTotUM; i++)
1514 for (Int_t j = 0; j < fgkRow; j++)
1516 for (Int_t k = 0; k < fgkCol; k++)
1518 fCPVCounter[i][j][k] = 0;
1519 fPRECounter[i][j][k] = 0;
1524 //____________________________________________________________________________
1525 void AliPMDDigitizer::ResetSDigit()
1529 if (fSDigits) fSDigits->Delete();
1531 //____________________________________________________________________________
1532 void AliPMDDigitizer::ResetDigit()
1536 if (fDigits) fDigits->Delete();
1538 //____________________________________________________________________________
1540 void AliPMDDigitizer::ResetCellADC()
1542 // Clears individual cells edep and track number
1543 for (Int_t i = 0; i < fgkTotUM; i++)
1545 for (Int_t j = 0; j < fgkRow; j++)
1547 for (Int_t k = 0; k < fgkCol; k++)
1551 fCPVTrackNo[i][j][k] = 0;
1552 fPRETrackNo[i][j][k] = 0;
1553 fCPVTrackPid[i][j][k] = -1;
1554 fPRETrackPid[i][j][k] = -1;
1559 //____________________________________________________________________________
1561 void AliPMDDigitizer::UnLoad(Option_t *option)
1563 // Unloads all the root files
1565 const char *cS = strstr(option,"S");
1566 const char *cD = strstr(option,"D");
1568 fRunLoader->UnloadgAlice();
1569 fRunLoader->UnloadHeader();
1570 fRunLoader->UnloadKinematics();
1574 fPMDLoader->UnloadHits();
1578 fPMDLoader->UnloadHits();
1579 fPMDLoader->UnloadSDigits();
1583 //----------------------------------------------------------------------
1584 Float_t AliPMDDigitizer::Gain(Int_t det, Int_t smn, Int_t row, Int_t col) const
1586 // returns of the gain of the cell
1587 // Added this method by ZA
1589 //cout<<" I am here in gain "<<fCalibData<< "smn,row, col "<<smn
1590 //<<" "<<row<<" "<<col<<endl;
1593 AliError("No calibration data loaded from CDB!!!");
1598 GainFact = fCalibGain->GetGainFact(det,smn,row,col);
1601 //----------------------------------------------------------------------
1602 AliPMDCalibData* AliPMDDigitizer::GetCalibGain() const
1604 // The run number will be centralized in AliCDBManager,
1605 // you don't need to set it here!
1606 // Added this method by ZA
1607 // Cleaned up by Alberto
1608 AliCDBEntry *entry = AliCDBManager::Instance()->Get("PMD/Calib/Gain");
1610 if(!entry) AliFatal("Calibration object retrieval failed!");
1612 AliPMDCalibData *calibdata=0;
1613 if (entry) calibdata = (AliPMDCalibData*) entry->GetObject();
1615 if (!calibdata) AliFatal("No calibration data from calibration database !");
1619 //----------------------------------------------------------------------
1620 AliPMDPedestal* AliPMDDigitizer::GetCalibPed() const
1622 // The run number will be centralized in AliCDBManager,
1623 // you don't need to set it here!
1625 AliCDBEntry *entry = AliCDBManager::Instance()->Get("PMD/Calib/Ped");
1627 if(!entry) AliFatal("Pedestal object retrieval failed!");
1629 AliPMDPedestal *pedestal=0;
1630 if (entry) pedestal = (AliPMDPedestal*) entry->GetObject();
1632 if (!pedestal) AliFatal("No pedestal data from calibration database !");