]> git.uio.no Git - u/mrichter/AliRoot.git/blob - PMD/AliPMDDigitizer.cxx
Bug fix (Y. Belikov)
[u/mrichter/AliRoot.git] / PMD / AliPMDDigitizer.cxx
1 /***************************************************************************
2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3  *                                                                        *
4  * Author: The ALICE Off-line Project.                                    *
5  * Contributors are mentioned in the code where appropriate.              *
6  *                                                                        *
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 //-----------------------------------------------------//
16 //                                                     //
17 //  Source File : PMDDigitizer.cxx, Version 00         //
18 //                                                     //
19 //  Date   : September 20 2002                         //
20 //                                                     //
21 //-----------------------------------------------------//
22
23 #include <Riostream.h>
24 #include <TBRIK.h>
25 #include <TNode.h>
26 #include <TTree.h>
27 #include <TGeometry.h>
28 #include <TObjArray.h>
29 #include <TClonesArray.h>
30 #include <TFile.h>
31 #include <TNtuple.h>
32 #include <TParticle.h>
33 #include <TRandom.h>
34
35 #include "AliLog.h"
36 #include "AliRun.h"
37 #include "AliHit.h"
38 #include "AliDetector.h"
39 #include "AliRunLoader.h"
40 #include "AliLoader.h"
41 #include "AliConfig.h"
42 #include "AliMagF.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"
49 #include "AliMC.h"
50
51 #include "AliPMD.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"
59
60
61 ClassImp(AliPMDDigitizer)
62
63 AliPMDDigitizer::AliPMDDigitizer() :
64   fRunLoader(0),
65   fPMDHit(0),
66   fPMD(0),
67   fPMDLoader(0),
68   fCalibGain(GetCalibGain()),
69   fCalibPed(GetCalibPed()),
70   fSDigits(0),
71   fDigits(0),
72   fCPVCell(0),
73   fCell(0),
74   fNsdigit(0),
75   fNdigit(0),
76   fDetNo(0),
77   fZPos(361.5)   // in units of cm, default position of PMD
78 {
79   // Default Constructor
80   //
81   for (Int_t i = 0; i < fgkTotUM; i++)
82     {
83       for (Int_t j = 0; j < fgkRow; j++)
84         {
85           for (Int_t k = 0; k < fgkCol; k++)
86             {
87               fCPV[i][j][k] = 0.;
88               fPRE[i][j][k] = 0.;
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             }
94         }
95     }
96
97 }
98 //____________________________________________________________________________
99 AliPMDDigitizer::AliPMDDigitizer(const AliPMDDigitizer& digitizer):
100   AliDigitizer(digitizer),
101   fRunLoader(0),
102   fPMDHit(0),
103   fPMD(0),
104   fPMDLoader(0),
105   fCalibGain(GetCalibGain()),
106   fCalibPed(GetCalibPed()),
107   fSDigits(0),
108   fDigits(0),
109   fCPVCell(0),
110   fCell(0),
111   fNsdigit(0),
112   fNdigit(0),
113   fDetNo(0),
114   fZPos(361.5)   // in units of cm, default position of PMD
115 {
116   // copy constructor
117   AliError("Copy constructor not allowed ");
118   
119 }
120 //____________________________________________________________________________
121 AliPMDDigitizer & AliPMDDigitizer::operator=(const AliPMDDigitizer& /*digitizer*/)
122 {
123   // Assignment operator
124   AliError("Assignement operator not allowed ");
125
126   return *this;
127 }
128 //____________________________________________________________________________
129 AliPMDDigitizer::AliPMDDigitizer(AliRunDigitizer* manager):
130   AliDigitizer(manager),
131   fRunLoader(0),
132   fPMDHit(0),
133   fPMD(0),
134   fPMDLoader(0),
135   fCalibGain(GetCalibGain()),
136   fCalibPed(GetCalibPed()),
137   fSDigits(new TClonesArray("AliPMDsdigit", 1000)),
138   fDigits(new TClonesArray("AliPMDdigit", 1000)),
139   fCPVCell(0),
140   fCell(0),
141   fNsdigit(0),
142   fNdigit(0),
143   fDetNo(0),
144   fZPos(361.5)// in units of cm, This is the default position of PMD
145 {
146   // ctor which should be used
147
148
149   for (Int_t i = 0; i < fgkTotUM; i++)
150     {
151       for (Int_t j = 0; j < fgkRow; j++)
152         {
153           for (Int_t k = 0; k < fgkCol; k++)
154             {
155               fCPV[i][j][k] = 0.;
156               fPRE[i][j][k] = 0.;
157               fCPVCounter[i][j][k] =  0; 
158               fPRECounter[i][j][k] =  0;
159               fCPVTrackNo[i][j][k] = -1;
160               fPRETrackNo[i][j][k] = -1;
161             }
162         }
163     }
164 }
165
166 //____________________________________________________________________________
167 AliPMDDigitizer::~AliPMDDigitizer()
168 {
169   // Default Destructor
170   //
171   if (fSDigits) {
172     fSDigits->Delete();
173     delete fSDigits;
174     fSDigits=0;
175   }
176   if (fDigits) {
177     fDigits->Delete();
178     delete fDigits;
179     fDigits=0;
180   }
181   fCPVCell.Delete();
182   fCell.Delete();
183 }
184 //
185 // Member functions
186 //
187 //____________________________________________________________________________
188 void AliPMDDigitizer::OpengAliceFile(const char *file, Option_t *option)
189 {
190   // Loads galice.root file and corresponding header, kinematics
191   // hits and sdigits or digits depending on the option
192   //
193
194   TString evfoldname = AliConfig::GetDefaultEventFolderName();
195   fRunLoader = AliRunLoader::GetRunLoader(evfoldname);
196   if (!fRunLoader)
197       fRunLoader = AliRunLoader::Open(file,AliConfig::GetDefaultEventFolderName(), "UPDATE");
198   
199   if (!fRunLoader)
200    {
201      AliError(Form("Can not open session for file %s.",file));
202    }
203
204   const char *cHS = strstr(option,"HS");
205   const char *cHD = strstr(option,"HD");
206   const char *cSD = strstr(option,"SD");
207   
208   if(cHS || cHD)
209     {
210       if (!fRunLoader->GetAliRun()) fRunLoader->LoadgAlice();
211       if (!fRunLoader->TreeE()) fRunLoader->LoadHeader();
212       if (!fRunLoader->TreeK()) fRunLoader->LoadKinematics();
213   
214       gAlice = fRunLoader->GetAliRun();
215   
216       if (gAlice)
217         {
218           AliDebug(1,"Alirun object found");
219         }
220       else
221         {
222           AliError("Could not found Alirun object");
223         }
224   
225       fPMD  = (AliPMD*)gAlice->GetDetector("PMD");
226     }
227
228   fPMDLoader = fRunLoader->GetLoader("PMDLoader");
229   if (fPMDLoader == 0x0)
230     {
231       AliError("Can not find PMDLoader");
232     }
233
234
235   if (cHS)
236     {
237       fPMDLoader->LoadHits("READ");
238       fPMDLoader->LoadSDigits("recreate");
239     }
240   else if (cHD)
241     {
242       fPMDLoader->LoadHits("READ");
243       fPMDLoader->LoadDigits("recreate");
244     }
245   else if (cSD)
246     {
247       fPMDLoader->LoadSDigits("READ");
248       fPMDLoader->LoadDigits("recreate");
249     }
250 }
251 //____________________________________________________________________________
252 void AliPMDDigitizer::Hits2SDigits(Int_t ievt)
253 {
254   // This reads the PMD Hits tree and assigns the right track number
255   // to a cell and stores in the summable digits tree
256   //
257
258   const Int_t kPi0 = 111;
259   const Int_t kGamma = 22;
260   Int_t npmd;
261   Int_t trackno;
262   Int_t smnumber;
263   Int_t trackpid;
264   Int_t mtrackno;
265   Int_t mtrackpid;
266
267   Float_t xPos, yPos, zPos;
268   Int_t xpad = -1, ypad = -1;
269   Float_t edep;
270   Float_t vx = -999.0, vy = -999.0, vz = -999.0;
271
272
273   if (!fSDigits) fSDigits = new TClonesArray("AliPMDsdigit", 1000);
274   ResetSDigit();
275
276   AliDebug(1,Form("Event Number = %d",ievt));
277   Int_t nparticles = fRunLoader->GetHeader()->GetNtrack();
278   AliDebug(1,Form("Number of Particles = %d",nparticles));
279   fRunLoader->GetEvent(ievt);
280   // ------------------------------------------------------- //
281   // Pointer to specific detector hits.
282   // Get pointers to Alice detectors and Hits containers
283
284   TTree* treeH = fPMDLoader->TreeH();
285   
286   Int_t ntracks    = (Int_t) treeH->GetEntries();
287   AliDebug(1,Form("Number of Tracks in the TreeH = %d", ntracks));
288   TTree* treeS = fPMDLoader->TreeS();
289   if (treeS == 0x0)
290     {
291       fPMDLoader->MakeTree("S");
292       treeS = fPMDLoader->TreeS();
293     }
294   Int_t bufsize = 16000;
295   treeS->Branch("PMDSDigit", &fSDigits, bufsize);
296   
297   TClonesArray* hits = 0;
298   if (fPMD) hits = fPMD->Hits();
299
300   // Start loop on tracks in the hits containers
301
302   for (Int_t track=0; track<ntracks;track++)
303     {
304       gAlice->ResetHits();
305       treeH->GetEvent(track);
306       if (fPMD)
307         {
308           npmd = hits->GetEntriesFast();
309           for (Int_t ipmd = 0; ipmd < npmd; ipmd++)
310             {
311               fPMDHit = (AliPMDhit*) hits->UncheckedAt(ipmd);
312               trackno = fPMDHit->GetTrack();
313               //  get kinematics of the particles
314
315               TParticle* mparticle = gAlice->GetMCApp()->Particle(trackno);
316               trackpid  = mparticle->GetPdgCode();
317               Int_t  ks = mparticle->GetStatusCode();
318               Int_t imo;
319               Int_t tracknoOld=0, trackpidOld=0, statusOld = 0;
320               
321               if (mparticle->GetFirstMother() == -1)
322                 {
323                   tracknoOld  = trackno;
324                   trackpidOld = trackpid;
325                   statusOld   = -1;
326                 }
327               Int_t igstatus = 0;
328               //------------------modified by Mriganka ----------------------
329               if(ks==1||(imo = mparticle->GetFirstMother())<0 ){
330                 vx = mparticle->Vx();
331                 vy = mparticle->Vy();
332                 vz = mparticle->Vz();
333                 
334                 if(trackpid==kGamma||trackpid==11||trackpid==-11||
335                    trackpid==kPi0)igstatus=1;
336               }
337               
338               
339               while(((imo = mparticle->GetFirstMother()) >= 0)&& 
340                     (ks = mparticle->GetStatusCode() <1) )
341                 {
342                   mparticle =  gAlice->GetMCApp()->Particle(imo);
343                   trackpid = mparticle->GetPdgCode();
344                   ks = mparticle->GetStatusCode();
345                   vx = mparticle->Vx();
346                   vy = mparticle->Vy();
347                   vz = mparticle->Vz();
348                   
349                   trackno=imo;
350                                 
351                     }
352          
353               if(trackpid==kGamma||trackpid==11||trackpid==-11||
354                  trackpid==kPi0)igstatus=1;
355               mtrackpid=trackpid;
356               mtrackno=trackno;
357               trackpid=trackpidOld;
358               trackno=tracknoOld;
359               
360               //-----------------end of modification----------------
361               xPos = fPMDHit->X();
362               yPos = fPMDHit->Y();
363               zPos = fPMDHit->Z();
364               
365               edep       = fPMDHit->GetEnergy();
366               Int_t vol1 = fPMDHit->GetVolume(1); // Column
367               Int_t vol2 = fPMDHit->GetVolume(2); // Row
368               Int_t vol7 = fPMDHit->GetVolume(7); // UnitModule
369               Int_t vol8 = fPMDHit->GetVolume(8); // SuperModule
370
371
372               // -----------------------------------------//
373               // In new geometry after adding electronics //
374               // For Super Module 1 & 2                   //
375               //  nrow = 48, ncol = 96                    //
376               // For Super Module 3 & 4                   //
377               //  nrow = 96, ncol = 48                    //
378               // -----------------------------------------//
379
380
381               
382               smnumber = (vol8-1)*6 + vol7;
383
384               if (vol8 == 1 || vol8 == 2)
385                 {
386                   xpad = vol2;
387                   ypad = vol1;
388                 }
389               else if (vol8 == 3 || vol8 == 4)
390                 {
391                   xpad = vol1;
392                   ypad = vol2;
393                 }
394
395               AliDebug(2,Form("Zposition = %f Edeposition = %f",zPos,edep));
396               //Float_t zposition = TMath::Abs(zPos);
397               if (zPos < fZPos)
398                 {
399                   // CPV
400                   fDetNo = 1;
401                 }
402               else if (zPos > fZPos)
403                 {
404                   // PMD
405                   fDetNo = 0;
406                 }
407               Int_t smn = smnumber - 1;
408               Int_t ixx = xpad     - 1;
409               Int_t iyy = ypad     - 1;
410               if (fDetNo == 0)
411                 {
412                   fPRE[smn][ixx][iyy] += edep;
413                   fPRECounter[smn][ixx][iyy]++;
414
415                   AliPMDcell* cell = new AliPMDcell(mtrackno,smn,ixx,iyy,edep);
416                   fCell.Add(cell);
417                 }
418               else if(fDetNo == 1)
419                 {
420                   fCPV[smn][ixx][iyy] += edep;
421                   fCPVCounter[smn][ixx][iyy]++;
422                   AliPMDcell* cpvcell = new AliPMDcell(mtrackno,smn,ixx,iyy,edep); 
423                   fCPVCell.Add(cpvcell);
424                 }
425             }
426         }
427     } // Track Loop ended
428   TrackAssignment2CPVCell();
429   TrackAssignment2Cell();
430   ResetCell();
431
432   Float_t deltaE      = 0.;
433   Int_t   detno       = 0;
434   Int_t   trno        = -1;
435
436   for (Int_t idet = 0; idet < 2; idet++)
437     {
438       for (Int_t ism = 0; ism < fgkTotUM; ism++)
439         {
440           for (Int_t jrow = 0; jrow < fgkRow; jrow++)
441             {
442               for (Int_t kcol = 0; kcol < fgkCol; kcol++)
443                 {
444                   if (idet == 0)
445                     {
446                       deltaE = fPRE[ism][jrow][kcol];
447                       trno   = fPRETrackNo[ism][jrow][kcol];
448                       detno = 0;
449                     }
450                   else if (idet == 1)
451                     {
452                       deltaE = fCPV[ism][jrow][kcol];
453                       trno   = fCPVTrackNo[ism][jrow][kcol];
454                       detno = 1;
455                     }
456                   if (deltaE > 0.)
457                     {
458                       AddSDigit(trno,detno,ism,jrow,kcol,deltaE);
459                     }
460                 }
461             }
462           treeS->Fill();
463           ResetSDigit();
464         }
465     }
466   fPMDLoader->WriteSDigits("OVERWRITE");
467   ResetCellADC();
468 }
469 //____________________________________________________________________________
470
471 void AliPMDDigitizer::Hits2Digits(Int_t ievt)
472 {
473   // This reads the PMD Hits tree and assigns the right track number
474   // to a cell and stores in the digits tree
475   //
476   const Int_t kPi0 = 111;
477   const Int_t kGamma = 22;
478   Int_t npmd;
479   Int_t trackno;
480   Int_t smnumber;
481   Int_t trackpid;
482   Int_t mtrackno;
483   Int_t mtrackpid;
484
485   Float_t xPos, yPos, zPos;
486   Int_t xpad = -1, ypad = -1;
487   Float_t edep;
488   Float_t vx = -999.0, vy = -999.0, vz = -999.0;
489
490   if (!fDigits) fDigits = new TClonesArray("AliPMDdigit", 1000);
491   ResetDigit();
492
493   AliDebug(1,Form("Event Number =  %d",ievt));
494   Int_t nparticles = fRunLoader->GetHeader()->GetNtrack();
495   AliDebug(1,Form("Number of Particles = %d", nparticles));
496   fRunLoader->GetEvent(ievt);
497   // ------------------------------------------------------- //
498   // Pointer to specific detector hits.
499   // Get pointers to Alice detectors and Hits containers
500
501   fPMD  = (AliPMD*)gAlice->GetDetector("PMD");
502   fPMDLoader = fRunLoader->GetLoader("PMDLoader");
503
504   if (fPMDLoader == 0x0)
505     {
506       AliError("Can not find PMD or PMDLoader");
507     }
508   TTree* treeH = fPMDLoader->TreeH();
509   Int_t ntracks    = (Int_t) treeH->GetEntries();
510   AliDebug(1,Form("Number of Tracks in the TreeH = %d", ntracks));
511   fPMDLoader->LoadDigits("recreate");
512   TTree* treeD = fPMDLoader->TreeD();
513   if (treeD == 0x0)
514     {
515       fPMDLoader->MakeTree("D");
516       treeD = fPMDLoader->TreeD();
517     }
518   Int_t bufsize = 16000;
519   treeD->Branch("PMDDigit", &fDigits, bufsize);
520   
521   TClonesArray* hits = 0;
522   if (fPMD) hits = fPMD->Hits();
523
524   // Start loop on tracks in the hits containers
525
526   for (Int_t track=0; track<ntracks;track++)
527     {
528       gAlice->ResetHits();
529       treeH->GetEvent(track);
530       
531       if (fPMD)
532         {
533           npmd = hits->GetEntriesFast();
534           for (Int_t ipmd = 0; ipmd < npmd; ipmd++)
535             {
536               fPMDHit = (AliPMDhit*) hits->UncheckedAt(ipmd);
537               trackno = fPMDHit->GetTrack();
538               
539               //  get kinematics of the particles
540               
541               TParticle* mparticle = gAlice->GetMCApp()->Particle(trackno);
542               trackpid  = mparticle->GetPdgCode();
543               Int_t  ks = mparticle->GetStatusCode();
544               Int_t imo;
545               Int_t tracknoOld=0, trackpidOld=0, statusOld = 0;
546               if (mparticle->GetFirstMother() == -1)
547                 {
548                   tracknoOld  = trackno;
549                   trackpidOld = trackpid;
550                   statusOld   = -1;
551                 }
552
553               Int_t igstatus = 0;
554               //-----------------------modified by Mriganka ------------------
555               if(ks==1||(imo = mparticle->GetFirstMother())<0 ){
556                 vx = mparticle->Vx();
557                 vy = mparticle->Vy();
558                 vz = mparticle->Vz();
559                 
560                 if(trackpid==kGamma||trackpid==11||trackpid==-11||trackpid==kPi0)
561                   igstatus=1;
562               }
563               
564               
565               while(((imo = mparticle->GetFirstMother()) >= 0)&& 
566                     (ks = mparticle->GetStatusCode() <1) )
567                 {
568                   mparticle =  gAlice->GetMCApp()->Particle(imo);
569                   trackpid = mparticle->GetPdgCode();
570                   ks = mparticle->GetStatusCode();
571                   vx = mparticle->Vx();
572                   vy = mparticle->Vy();
573                   vz = mparticle->Vz();
574                   
575                   trackno=imo;
576                                 
577                     }
578          
579               if(trackpid==kGamma||trackpid==11||trackpid==-11||trackpid==kPi0)
580                 igstatus=1;
581               mtrackpid=trackpid;
582               mtrackno=trackno;
583               trackpid=trackpidOld;
584               trackno=tracknoOld;
585               
586               //-----------------end of modification----------------
587               xPos = fPMDHit->X();
588               yPos = fPMDHit->Y();
589               zPos = fPMDHit->Z();
590               edep       = fPMDHit->GetEnergy();
591               Int_t vol1 = fPMDHit->GetVolume(1); // Column
592               Int_t vol2 = fPMDHit->GetVolume(2); // Row
593               Int_t vol7 = fPMDHit->GetVolume(7); // UnitModule
594               Int_t vol8 = fPMDHit->GetVolume(8); // SuperModule
595
596
597               // -----------------------------------------//
598               // In new geometry after adding electronics //
599               // For Super Module 1 & 2                   //
600               //  nrow = 48, ncol = 96                    //
601               // For Super Module 3 & 4                   //
602               //  nrow = 96, ncol = 48                    //
603               // -----------------------------------------//
604               
605               smnumber = (vol8-1)*6 + vol7;
606
607               if (vol8 == 1 || vol8 == 2)
608                 {
609                   xpad = vol2;
610                   ypad = vol1;
611                 }
612               else if (vol8 == 3 || vol8 == 4)
613                 {
614                   xpad = vol1;
615                   ypad = vol2;
616                 }
617
618               AliDebug(2,Form("ZPosition = %f Edeposition = %d",zPos,edep));
619               //Float_t zposition = TMath::Abs(zPos);
620
621               if (zPos < fZPos)
622                 {
623                   // CPV
624                   fDetNo = 1;
625                 }
626               else if (zPos > fZPos)
627                 {
628                   // PMD
629                   fDetNo = 0;
630                 }
631
632               Int_t smn = smnumber - 1;
633               Int_t ixx = xpad     - 1;
634               Int_t iyy = ypad     - 1;
635               if (fDetNo == 0)
636                 {
637                   fPRE[smn][ixx][iyy] += edep;
638                   fPRECounter[smn][ixx][iyy]++;
639
640                   AliPMDcell* cell = new AliPMDcell(mtrackno,smn,ixx,iyy,edep);
641
642                   fCell.Add(cell);
643                 }
644               else if(fDetNo == 1)
645                 {
646                   fCPV[smn][ixx][iyy] += edep;
647                   fCPVCounter[smn][ixx][iyy]++;
648                   AliPMDcell* cpvcell = new AliPMDcell(mtrackno,smn,ixx,iyy,edep); 
649                   fCPVCell.Add(cpvcell);
650                 }
651             }
652         }
653     } // Track Loop ended
654   TrackAssignment2CPVCell();
655   TrackAssignment2Cell();
656   ResetCell();
657
658   Float_t gain1;
659   Float_t adc;
660   Float_t deltaE = 0.;
661   Int_t detno = 0;
662   Int_t trno = 1;
663   for (Int_t idet = 0; idet < 2; idet++)
664   {
665       for (Int_t ism = 0; ism < fgkTotUM; ism++)
666       {
667           for (Int_t jrow = 0; jrow < fgkRow; jrow++)
668           {
669               for (Int_t kcol = 0; kcol < fgkCol; kcol++)
670               {
671                   if (idet == 0)
672                   {
673                       deltaE = fPRE[ism][jrow][kcol];
674                       trno   = fPRETrackNo[ism][jrow][kcol];
675                       detno = 0;
676                   }
677                   else if (idet == 1)
678                   {
679                       deltaE = fCPV[ism][jrow][kcol];
680                       trno   = fCPVTrackNo[ism][jrow][kcol];
681                       detno = 1;
682                   }
683                   if (deltaE > 0.)
684                   {
685                       MeV2ADC(deltaE,adc);
686
687                       // To decalibrate the adc values
688                       //
689                       gain1 = Gain(idet,ism,jrow,kcol);
690                       if (gain1 != 0.)
691                       {
692                           Int_t adcDecalib = (Int_t)(adc/gain1);
693                           adc = (Float_t) adcDecalib;
694                       }
695                       else if(gain1 == 0.)
696                       {
697                           adc = 0.;
698                       }
699
700                       // Pedestal Decalibration
701                       Int_t   pedmeanrms = 
702                           fCalibPed->GetPedMeanRms(idet,ism,jrow,kcol);
703                       Int_t   pedrms1    = (Int_t) pedmeanrms%100;
704                       Float_t pedrms     = (Float_t)pedrms1/10.;
705                       Float_t pedmean    = 
706                           (Float_t) (pedmeanrms - pedrms1)/1000.0;
707                       if (adc > 0.)
708                       {
709                           adc += (pedmean + 3.0*pedrms);
710                           AddDigit(trno,detno,ism,jrow,kcol,adc);
711                       }
712                   }
713               } // column loop
714           } // row    loop
715           treeD->Fill();
716           ResetDigit();
717       } // supermodule loop
718   } // detector loop
719   
720   fPMDLoader->WriteDigits("OVERWRITE");
721   ResetCellADC();
722   
723 }
724 //____________________________________________________________________________
725
726
727 void AliPMDDigitizer::SDigits2Digits(Int_t ievt)
728 {
729   // This reads the PMD sdigits tree and converts energy deposition
730   // in a cell to ADC and stores in the digits tree
731   //
732
733   fRunLoader->GetEvent(ievt);
734
735   TTree* treeS = fPMDLoader->TreeS();
736   AliPMDsdigit  *pmdsdigit;
737   TBranch *branch = treeS->GetBranch("PMDSDigit");
738   if(!branch)
739     {
740       AliError("PMD Sdigit branch does not exist");
741       return;
742     }
743   if (!fSDigits) fSDigits = new TClonesArray("AliPMDsdigit", 1000);
744   branch->SetAddress(&fSDigits);
745
746   TTree* treeD = fPMDLoader->TreeD();
747   if (treeD == 0x0)
748     {
749       fPMDLoader->MakeTree("D");
750       treeD = fPMDLoader->TreeD();
751     }
752   Int_t bufsize = 16000;
753   if (!fDigits) fDigits = new TClonesArray("AliPMDdigit", 1000);
754   treeD->Branch("PMDDigit", &fDigits, bufsize);
755
756   Int_t   trno, det, smn;
757   Int_t   irow, icol;
758   Float_t edep, adc;
759
760   Int_t nmodules = (Int_t) treeS->GetEntries();
761   AliDebug(1,Form("Number of modules = %d",nmodules));
762
763   for (Int_t imodule = 0; imodule < nmodules; imodule++)
764     {
765       treeS->GetEntry(imodule);
766       Int_t nentries = fSDigits->GetLast();
767       AliDebug(2,Form("Number of entries per module = %d",nentries+1));
768       for (Int_t ient = 0; ient < nentries+1; ient++)
769         {
770           pmdsdigit = (AliPMDsdigit*)fSDigits->UncheckedAt(ient);
771           trno   = pmdsdigit->GetTrackNumber();
772           det    = pmdsdigit->GetDetector();
773           smn    = pmdsdigit->GetSMNumber();
774           irow   = pmdsdigit->GetRow();
775           icol   = pmdsdigit->GetColumn();
776           edep   = pmdsdigit->GetCellEdep();
777
778           MeV2ADC(edep,adc);
779
780           // To decalibrte the adc values
781           //
782           Float_t gain1 = Gain(det,smn,irow,icol);
783           if (gain1 != 0.)
784           {
785             Int_t adcDecalib = (Int_t)(adc/gain1);
786             adc = (Float_t) adcDecalib;
787           }
788           else if(gain1 == 0.)
789           {
790               adc = 0.;
791           }
792           // Pedestal Decalibration
793           Int_t   pedmeanrms = fCalibPed->GetPedMeanRms(det,smn,irow,icol);
794           Int_t   pedrms1    = (Int_t) pedmeanrms%100;
795           Float_t pedrms     = (Float_t)pedrms1/10.;
796           Float_t pedmean    = (Float_t) (pedmeanrms - pedrms1)/1000.0;
797           if(adc > 0.)
798           {
799               adc += (pedmean + 3.0*pedrms);
800               AddDigit(trno,det,smn,irow,icol,adc);
801           }
802
803         }
804       treeD->Fill();
805       ResetDigit();
806     }
807   fPMDLoader->WriteDigits("OVERWRITE");
808
809 }
810 //____________________________________________________________________________
811 void AliPMDDigitizer::Exec(Option_t *option)
812 {
813   // Does the event merging and digitization
814   const char *cdeb = strstr(option,"deb");
815   if(cdeb)
816     {
817       AliDebug(100," *** PMD Exec is called ***");
818     }
819
820   Int_t ninputs = fManager->GetNinputs();
821   AliDebug(1,Form("Number of files to be processed = %d",ninputs));
822   ResetCellADC();
823
824   for (Int_t i = 0; i < ninputs; i++)
825     {
826       Int_t troffset = fManager->GetMask(i);
827       MergeSDigits(i, troffset);
828     }
829
830   fRunLoader = AliRunLoader::GetRunLoader(fManager->GetOutputFolderName());
831   fPMD  = (AliPMD*)gAlice->GetDetector("PMD");
832   fPMDLoader = fRunLoader->GetLoader("PMDLoader");
833   if (fPMDLoader == 0x0)
834     {
835       AliError("Can not find PMD or PMDLoader");
836     }
837   fPMDLoader->LoadDigits("update");
838   TTree* treeD = fPMDLoader->TreeD();
839   if (treeD == 0x0)
840     {
841       fPMDLoader->MakeTree("D");
842       treeD = fPMDLoader->TreeD();
843     }
844   Int_t bufsize = 16000;
845   if (!fDigits) fDigits = new TClonesArray("AliPMDdigit", 1000);
846   treeD->Branch("PMDDigit", &fDigits, bufsize);
847
848   Float_t adc;
849   Float_t deltaE = 0.;
850   Int_t detno = 0;
851   Int_t trno = 1;
852
853   for (Int_t idet = 0; idet < 2; idet++)
854     {
855       for (Int_t ism = 0; ism < fgkTotUM; ism++)
856         {
857           for (Int_t jrow = 0; jrow < fgkRow; jrow++)
858             {
859               for (Int_t kcol = 0; kcol < fgkCol; kcol++)
860                 {
861                   if (idet == 0)
862                     {
863                       deltaE = fPRE[ism][jrow][kcol];
864                       trno   = fPRETrackNo[ism][jrow][kcol];
865                       detno = 0;
866                     }
867                   else if (idet == 1)
868                     {
869                       deltaE = fCPV[ism][jrow][kcol];
870                       trno   = fCPVTrackNo[ism][jrow][kcol];
871                       detno = 1;
872                     }
873                   if (deltaE > 0.)
874                     {
875                       MeV2ADC(deltaE,adc);
876
877                       //
878                       // Gain decalibration
879                       //
880                       Float_t gain1 = Gain(idet,ism,jrow,kcol);
881
882                       if (gain1 != 0.)
883                       {
884                           Int_t adcDecalib = (Int_t)(adc/gain1);
885                           adc = (Float_t) adcDecalib;
886                       }
887                       else if(gain1 == 0.)
888                       {
889                           adc = 0.;
890                       }
891                       // Pedestal Decalibration
892                       Int_t   pedmeanrms = 
893                           fCalibPed->GetPedMeanRms(idet,ism,jrow,kcol);
894                       Int_t   pedrms1    = (Int_t) pedmeanrms%100;
895                       Float_t pedrms     = (Float_t)pedrms1/10.;
896                       Float_t pedmean    = 
897                           (Float_t) (pedmeanrms - pedrms1)/1000.0;
898                       if (adc > 0.)
899                       {
900                           adc += (pedmean + 3.0*pedrms);
901                           AddDigit(trno,detno,ism,jrow,kcol,adc);
902                       }
903
904                     }
905                 } // column loop
906             } // row    loop
907           treeD->Fill();
908           ResetDigit();
909         } // supermodule loop
910     } // detector loop
911   fPMDLoader->WriteDigits("OVERWRITE");
912   fPMDLoader->UnloadDigits();
913   ResetCellADC();
914 }
915 //____________________________________________________________________________
916 void AliPMDDigitizer::TrackAssignment2CPVCell()
917 {
918   // This block assigns the cell id when there are
919   // multiple tracks in a cell according to the
920   // energy deposition
921   // This method added by Ajay
922   Bool_t jsort = false;
923
924   Int_t i, j, k;
925
926   Int_t   *status1;
927   Int_t   *status2;
928   Int_t   *trnarray;  
929   Float_t *fracEdp;
930   Float_t *trEdp;
931   
932   Int_t   ****cpvTrack;
933   Float_t ****cpvEdep;
934
935   cpvTrack = new Int_t ***[fgkTotUM];
936   cpvEdep  = new Float_t ***[fgkTotUM];
937   for (i=0; i<fgkTotUM; i++)
938     {
939       cpvTrack[i] = new Int_t **[fgkRow];
940       cpvEdep[i]  = new Float_t **[fgkRow];
941     }
942
943   for (i = 0; i < fgkTotUM; i++)
944     {
945       for (j = 0; j < fgkRow; j++)
946         {
947           cpvTrack[i][j] = new Int_t *[fgkCol];
948           cpvEdep[i][j]  = new Float_t *[fgkCol];
949         }
950     }
951   for (i = 0; i < fgkTotUM; i++)
952     {
953       for (j = 0; j < fgkRow; j++)
954         {
955           for (k = 0; k < fgkCol; k++)
956             {
957               Int_t nn = fCPVCounter[i][j][k];
958               if(nn > 0)
959                 {
960                   cpvTrack[i][j][k] = new Int_t[nn];
961                   cpvEdep[i][j][k] = new Float_t[nn];
962                 }
963               else
964                 {
965                   nn = 1;
966                   cpvTrack[i][j][k] = new Int_t[nn];
967                   cpvEdep[i][j][k] = new Float_t[nn];
968                 }                     
969               fCPVCounter[i][j][k] = 0;
970             }
971         }
972     }
973
974
975   Int_t nentries = fCPVCell.GetEntries();
976  
977   Int_t   mtrackno, ism, ixp, iyp;
978   Float_t edep;
979   for (i = 0; i < nentries; i++)
980     {
981       AliPMDcell* cpvcell = (AliPMDcell*)fCPVCell.UncheckedAt(i);
982       
983       mtrackno = cpvcell->GetTrackNumber();
984       ism      = cpvcell->GetSMNumber();
985       ixp      = cpvcell->GetX();
986       iyp      = cpvcell->GetY();
987       edep     = cpvcell->GetEdep();
988       Int_t nn = fCPVCounter[ism][ixp][iyp];
989       cpvTrack[ism][ixp][iyp][nn] = (Int_t) mtrackno;
990       cpvEdep[ism][ixp][iyp][nn] = edep;
991       fCPVCounter[ism][ixp][iyp]++;
992     }
993   
994   Int_t iz, il;
995   Int_t im, ix, iy;
996   Int_t nn;
997   for (im=0; im<fgkTotUM; im++)
998     {
999       for (ix=0; ix<fgkRow; ix++)
1000         {
1001           for (iy=0; iy<fgkCol; iy++)
1002             {
1003               nn = fCPVCounter[im][ix][iy];
1004               if (nn > 1)
1005                 {
1006                   // This block handles if a cell is fired
1007                   // many times by many tracks
1008                   status1  = new Int_t[nn];
1009                   status2  = new Int_t[nn];
1010                   trnarray = new Int_t[nn];
1011                   for (iz = 0; iz < nn; iz++)
1012                     {
1013                       status1[iz] = cpvTrack[im][ix][iy][iz];
1014                     }
1015                   TMath::Sort(nn,status1,status2,jsort);
1016                   Int_t trackOld = -99999;
1017                   Int_t track, trCount = 0;
1018                   for (iz = 0; iz < nn; iz++)
1019                     {
1020                       track = status1[status2[iz]];
1021                       if (trackOld != track)
1022                         {
1023                           trnarray[trCount] = track;
1024                           trCount++;
1025                         }                             
1026                       trackOld = track;
1027                     }
1028                   delete [] status1;
1029                   delete [] status2;
1030                   Float_t totEdp = 0.;
1031                   trEdp = new Float_t[trCount];
1032                   fracEdp = new Float_t[trCount];
1033                   for (il = 0; il < trCount; il++)
1034                     {
1035                       trEdp[il] = 0.;
1036                       track = trnarray[il];
1037                       for (iz = 0; iz < nn; iz++)
1038                         {
1039                           if (track == cpvTrack[im][ix][iy][iz])
1040                             {
1041                               trEdp[il] += cpvEdep[im][ix][iy][iz];
1042                             }
1043                         }
1044                       totEdp += trEdp[il];
1045                     }
1046                   Int_t ilOld = 0;
1047                   Float_t fracOld = 0.;
1048                   
1049                   for (il = 0; il < trCount; il++)
1050                     {
1051                       fracEdp[il] = trEdp[il]/totEdp;
1052                       if (fracOld < fracEdp[il])
1053                         {
1054                           fracOld = fracEdp[il];
1055                           ilOld = il;
1056                         }
1057                     }
1058                   fCPVTrackNo[im][ix][iy] = trnarray[ilOld];
1059                   delete [] fracEdp;
1060                   delete [] trEdp;
1061                   delete [] trnarray;
1062                 }
1063               else if (nn == 1)
1064                 {
1065                   // This only handles if a cell is fired
1066                   // by only one track
1067                   
1068                   fCPVTrackNo[im][ix][iy] = cpvTrack[im][ix][iy][0];
1069                   
1070                 }
1071               else if (nn ==0)
1072                 {
1073                   // This is if no cell is fired
1074                   fCPVTrackNo[im][ix][iy] = -999;
1075                 }
1076             } // end of iy
1077         } // end of ix
1078     } // end of im
1079   
1080   // Delete all the pointers
1081   
1082  for (i = 0; i < fgkTotUM; i++)
1083     {
1084       for (j = 0; j < fgkRow; j++)
1085         {
1086           for (k = 0; k < fgkCol; k++)
1087             {
1088               delete []cpvTrack[i][j][k];
1089               delete []cpvEdep[i][j][k];
1090             }
1091         }
1092     }
1093  
1094   for (i = 0; i < fgkTotUM; i++)
1095     {
1096       for (j = 0; j < fgkRow; j++)
1097         {
1098           delete [] cpvTrack[i][j];
1099           delete [] cpvEdep[i][j];
1100         }
1101     }
1102   
1103   for (i = 0; i < fgkTotUM; i++)
1104     {
1105       delete [] cpvTrack[i];
1106       delete [] cpvEdep[i];
1107     }
1108   delete [] cpvTrack;
1109   delete [] cpvEdep;
1110   
1111   // 
1112   // End of the cell id assignment
1113   //
1114 }
1115 //____________________________________________________________________________
1116
1117 void AliPMDDigitizer::MergeSDigits(Int_t filenumber, Int_t troffset)
1118 {
1119   // merging sdigits
1120   fRunLoader = AliRunLoader::GetRunLoader(fManager->GetInputFolderName(filenumber));
1121   fPMDLoader = fRunLoader->GetLoader("PMDLoader");
1122   fPMDLoader->LoadSDigits("read");
1123   TTree* treeS = fPMDLoader->TreeS();
1124   AliPMDsdigit  *pmdsdigit;
1125   TBranch *branch = treeS->GetBranch("PMDSDigit");
1126   if (!fSDigits) fSDigits = new TClonesArray("AliPMDsdigit", 1000);
1127   branch->SetAddress(&fSDigits);
1128
1129   Int_t   itrackno, idet, ism;
1130   Int_t   ixp, iyp;
1131   Float_t edep;
1132   Int_t nmodules = (Int_t) treeS->GetEntries();
1133   AliDebug(1,Form("Number of Modules in the treeS = %d",nmodules));
1134   AliDebug(1,Form("Track Offset = %d",troffset));
1135   for (Int_t imodule = 0; imodule < nmodules; imodule++)
1136     {
1137       treeS->GetEntry(imodule);
1138       Int_t nentries = fSDigits->GetLast();
1139       AliDebug(2,Form("Number of Entries per Module = %d",nentries));
1140       for (Int_t ient = 0; ient < nentries+1; ient++)
1141         {
1142           pmdsdigit = (AliPMDsdigit*)fSDigits->UncheckedAt(ient);
1143           itrackno  = pmdsdigit->GetTrackNumber();
1144           idet      = pmdsdigit->GetDetector();
1145           ism       = pmdsdigit->GetSMNumber();
1146           ixp       = pmdsdigit->GetRow();
1147           iyp       = pmdsdigit->GetColumn();
1148           edep      = pmdsdigit->GetCellEdep();
1149           if (idet == 0)
1150             {
1151               if (fPRE[ism][ixp][iyp] < edep)
1152                 {
1153                   fPRETrackNo[ism][ixp][iyp] = troffset + itrackno;
1154                 }
1155               fPRE[ism][ixp][iyp] += edep;
1156             }
1157           else if (idet == 1)
1158             {
1159               if (fCPV[ism][ixp][iyp] < edep)
1160                 {
1161                   fCPVTrackNo[ism][ixp][iyp] = troffset + itrackno;
1162                 }
1163               fCPV[ism][ixp][iyp] += edep;
1164             }
1165         }
1166     }
1167
1168 }
1169 // ----------------------------------------------------------------------
1170 void AliPMDDigitizer::TrackAssignment2Cell()
1171 {
1172   // 
1173   // This block assigns the cell id when there are
1174   // multiple tracks in a cell according to the
1175   // energy deposition
1176   //
1177   Bool_t jsort = false;
1178
1179   Int_t i, j, k;
1180
1181   Int_t   *status1;
1182   Int_t   *status2;
1183   Int_t   *trnarray;
1184   Float_t *fracEdp;
1185   Float_t *trEdp;
1186   
1187   Int_t   ****pmdTrack;
1188   Float_t ****pmdEdep;
1189
1190   pmdTrack = new Int_t ***[fgkTotUM];
1191   pmdEdep  = new Float_t ***[fgkTotUM];
1192   for (i=0; i<fgkTotUM; i++)
1193     {
1194       pmdTrack[i] = new Int_t **[fgkRow];
1195       pmdEdep[i]  = new Float_t **[fgkRow];
1196     }
1197
1198   for (i = 0; i < fgkTotUM; i++)
1199     {
1200       for (j = 0; j < fgkRow; j++)
1201         {
1202           pmdTrack[i][j] = new Int_t *[fgkCol];
1203           pmdEdep[i][j]  = new Float_t *[fgkCol];
1204         }
1205     }
1206   
1207   for (i = 0; i < fgkTotUM; i++)
1208     {
1209       for (j = 0; j < fgkRow; j++)
1210         {
1211           for (k = 0; k < fgkCol; k++)
1212             {
1213               Int_t nn = fPRECounter[i][j][k];
1214               if(nn > 0)
1215                 {
1216                   pmdTrack[i][j][k] = new Int_t[nn];
1217                   pmdEdep[i][j][k] = new Float_t[nn];
1218                 }
1219               else
1220                 {
1221                   nn = 1;
1222                   pmdTrack[i][j][k] = new Int_t[nn];
1223                   pmdEdep[i][j][k] = new Float_t[nn];
1224                 }
1225               fPRECounter[i][j][k] = 0;
1226             }
1227         }
1228     }
1229
1230
1231   Int_t nentries = fCell.GetEntries();
1232
1233   Int_t   mtrackno, ism, ixp, iyp;
1234   Float_t edep;
1235
1236   for (i = 0; i < nentries; i++)
1237     {
1238       AliPMDcell* cell = (AliPMDcell*)fCell.UncheckedAt(i);
1239       
1240       mtrackno = cell->GetTrackNumber();
1241       ism      = cell->GetSMNumber();
1242       ixp      = cell->GetX();
1243       iyp      = cell->GetY();
1244       edep     = cell->GetEdep();
1245       Int_t nn = fPRECounter[ism][ixp][iyp];
1246       pmdTrack[ism][ixp][iyp][nn] = (Int_t) mtrackno;
1247       pmdEdep[ism][ixp][iyp][nn] = edep;
1248       fPRECounter[ism][ixp][iyp]++;
1249     }
1250   
1251   Int_t iz, il;
1252   Int_t im, ix, iy;
1253   Int_t nn;
1254   
1255   for (im=0; im<fgkTotUM; im++)
1256     {
1257       for (ix=0; ix<fgkRow; ix++)
1258         {
1259           for (iy=0; iy<fgkCol; iy++)
1260             {
1261               nn = fPRECounter[im][ix][iy];
1262               if (nn > 1)
1263                 {
1264                   // This block handles if a cell is fired
1265                   // many times by many tracks
1266                   status1  = new Int_t[nn];
1267                   status2  = new Int_t[nn];
1268                   trnarray = new Int_t[nn];
1269                   for (iz = 0; iz < nn; iz++)
1270                     {
1271                       status1[iz] = pmdTrack[im][ix][iy][iz];
1272                     }
1273                   TMath::Sort(nn,status1,status2,jsort);
1274                   Int_t trackOld = -99999;
1275                   Int_t track, trCount = 0;
1276                   for (iz = 0; iz < nn; iz++)
1277                     {
1278                       track = status1[status2[iz]];
1279                       if (trackOld != track)
1280                         {
1281                           trnarray[trCount] = track;
1282                           trCount++;
1283                         }
1284                       trackOld = track;
1285                     }
1286                   delete [] status1;
1287                   delete [] status2;
1288                   Float_t totEdp = 0.;
1289                   trEdp = new Float_t[trCount];
1290                   fracEdp = new Float_t[trCount];
1291                   for (il = 0; il < trCount; il++)
1292                     {
1293                       trEdp[il] = 0.;
1294                       track = trnarray[il];
1295                       for (iz = 0; iz < nn; iz++)
1296                         {
1297                           if (track == pmdTrack[im][ix][iy][iz])
1298                             {
1299                               trEdp[il] += pmdEdep[im][ix][iy][iz];
1300                             }
1301                         }
1302                       totEdp += trEdp[il];
1303                     }
1304                   Int_t ilOld = 0;
1305                   Float_t fracOld = 0.;
1306                   
1307                   for (il = 0; il < trCount; il++)
1308                     {
1309                       fracEdp[il] = trEdp[il]/totEdp;
1310                       if (fracOld < fracEdp[il])
1311                         {
1312                           fracOld = fracEdp[il];
1313                           ilOld = il;
1314                         }
1315                     }
1316                   fPRETrackNo[im][ix][iy] = trnarray[ilOld];
1317                   delete [] fracEdp;
1318                   delete [] trEdp;
1319                   delete [] trnarray;
1320                 }
1321               else if (nn == 1)
1322                 {
1323                   // This only handles if a cell is fired
1324                   // by only one track
1325                   
1326                   fPRETrackNo[im][ix][iy] = pmdTrack[im][ix][iy][0];
1327                   
1328                 }
1329               else if (nn ==0)
1330                 {
1331                   // This is if no cell is fired
1332                   fPRETrackNo[im][ix][iy] = -999;
1333                 }
1334             } // end of iy
1335         } // end of ix
1336     } // end of im
1337   
1338   // Delete all the pointers
1339   
1340   for (i = 0; i < fgkTotUM; i++)
1341     {
1342       for (j = 0; j < fgkRow; j++)
1343         {
1344           for (k = 0; k < fgkCol; k++)
1345             {
1346               delete [] pmdTrack[i][j][k];
1347               delete [] pmdEdep[i][j][k];
1348             }
1349         }
1350     }
1351   
1352   for (i = 0; i < fgkTotUM; i++)
1353     {
1354       for (j = 0; j < fgkRow; j++)
1355         {
1356           delete [] pmdTrack[i][j];
1357           delete [] pmdEdep[i][j];
1358         }
1359     }
1360   
1361   for (i = 0; i < fgkTotUM; i++)
1362     {
1363       delete [] pmdTrack[i];
1364       delete [] pmdEdep[i];
1365     }
1366   delete [] pmdTrack;
1367   delete [] pmdEdep;
1368   // 
1369   // End of the cell id assignment
1370   //
1371 }
1372 //____________________________________________________________________________
1373 void AliPMDDigitizer::MeV2ADC(Float_t mev, Float_t & adc) const
1374 {
1375   // This converts the simulated edep to ADC according to the
1376   // Test Beam Data
1377   //PS Test in September 2003 and 2006
1378   // KeV - ADC conversion for 12bit ADC
1379   // Modified by Ajay
1380  
1381   const Float_t kConstant   =  0.07;
1382   const Float_t kErConstant =  0.1;
1383   const Float_t kSlope      = 76.0;
1384   const Float_t kErSlope    =  5.0;
1385   
1386   Float_t cons   = gRandom->Gaus(kConstant,kErConstant);
1387   Float_t slop   = gRandom->Gaus(kSlope,kErSlope);
1388   
1389   Float_t adc12bit = slop*mev*0.001 + cons;
1390   
1391                            
1392   if(adc12bit < 1600.0)
1393     {
1394       adc = (Float_t) adc12bit;
1395     }
1396   else if (adc12bit >= 1600.0)
1397     {
1398       adc = 1600.0;
1399     }
1400 }
1401 //____________________________________________________________________________
1402 void AliPMDDigitizer::AddSDigit(Int_t trnumber, Int_t det, Int_t smnumber,
1403                                 Int_t irow, Int_t icol, Float_t adc)
1404 {
1405   // Add SDigit
1406   //
1407   if (!fSDigits) fSDigits = new TClonesArray("AliPMDsdigit", 1000);
1408   TClonesArray &lsdigits = *fSDigits;
1409   new(lsdigits[fNsdigit++])  AliPMDsdigit(trnumber,det,smnumber,irow,icol,adc);
1410 }
1411 //____________________________________________________________________________
1412
1413 void AliPMDDigitizer::AddDigit(Int_t trnumber, Int_t det, Int_t smnumber,
1414                                Int_t irow, Int_t icol, Float_t adc)
1415 {
1416   // Add Digit
1417   //
1418   if (!fDigits) fDigits = new TClonesArray("AliPMDdigit", 1000);
1419   TClonesArray &ldigits = *fDigits;
1420   new(ldigits[fNdigit++]) AliPMDdigit(trnumber,det,smnumber,irow,icol,adc);
1421 }
1422 //____________________________________________________________________________
1423
1424 void AliPMDDigitizer::SetZPosition(Float_t zpos)
1425 {
1426   fZPos = zpos;
1427 }
1428 //____________________________________________________________________________
1429 Float_t AliPMDDigitizer::GetZPosition() const
1430 {
1431   return fZPos;
1432 }
1433 //____________________________________________________________________________
1434
1435 void AliPMDDigitizer::ResetCell()
1436 {
1437   // clears the cell array and also the counter
1438   //  for each cell
1439   //
1440   fCPVCell.Delete();
1441   fCell.Delete();
1442   for (Int_t i = 0; i < fgkTotUM; i++)
1443     {
1444       for (Int_t j = 0; j < fgkRow; j++)
1445         {
1446           for (Int_t k = 0; k < fgkCol; k++)
1447             {
1448               fCPVCounter[i][j][k] = 0; 
1449               fPRECounter[i][j][k] = 0;
1450             }
1451         }
1452     }
1453 }
1454 //____________________________________________________________________________
1455 void AliPMDDigitizer::ResetSDigit()
1456 {
1457   // Clears SDigits
1458   fNsdigit = 0;
1459   if (fSDigits) fSDigits->Delete();
1460 }
1461 //____________________________________________________________________________
1462 void AliPMDDigitizer::ResetDigit()
1463 {
1464   // Clears Digits
1465   fNdigit = 0;
1466   if (fDigits) fDigits->Delete();
1467 }
1468 //____________________________________________________________________________
1469
1470 void AliPMDDigitizer::ResetCellADC()
1471 {
1472   // Clears individual cells edep and track number
1473   for (Int_t i = 0; i < fgkTotUM; i++)
1474     {
1475       for (Int_t j = 0; j < fgkRow; j++)
1476         {
1477           for (Int_t k = 0; k < fgkCol; k++)
1478             {
1479               fCPV[i][j][k] = 0.;
1480               fPRE[i][j][k] = 0.;
1481               fCPVTrackNo[i][j][k] = 0;
1482               fPRETrackNo[i][j][k] = 0;
1483             }
1484         }
1485     }
1486 }
1487 //____________________________________________________________________________
1488
1489 void AliPMDDigitizer::UnLoad(Option_t *option)
1490 {
1491   // Unloads all the root files
1492   //
1493   const char *cS = strstr(option,"S");
1494   const char *cD = strstr(option,"D");
1495
1496   fRunLoader->UnloadgAlice();
1497   fRunLoader->UnloadHeader();
1498   fRunLoader->UnloadKinematics();
1499
1500   if (cS)
1501     {
1502       fPMDLoader->UnloadHits();
1503     }
1504   if (cD)
1505     {
1506       fPMDLoader->UnloadHits();
1507       fPMDLoader->UnloadSDigits();
1508     }
1509 }
1510
1511 //----------------------------------------------------------------------
1512 Float_t AliPMDDigitizer::Gain(Int_t det, Int_t smn, Int_t row, Int_t col) const
1513 {
1514   // returns of the gain of the cell
1515   // Added this method by ZA
1516
1517   //cout<<" I am here in gain "<<fCalibData<< "smn,row, col "<<smn
1518   //<<" "<<row<<" "<<col<<endl;
1519
1520   if(!fCalibGain) {
1521     AliError("No calibration data loaded from CDB!!!");
1522     return 1;
1523   }
1524
1525   Float_t GainFact;
1526   GainFact = fCalibGain->GetGainFact(det,smn,row,col);
1527   return GainFact;
1528 }
1529 //----------------------------------------------------------------------
1530 AliPMDCalibData* AliPMDDigitizer::GetCalibGain() const
1531 {
1532   // The run number will be centralized in AliCDBManager,
1533   // you don't need to set it here!
1534   // Added this method by ZA
1535   // Cleaned up by Alberto
1536   AliCDBEntry  *entry = AliCDBManager::Instance()->Get("PMD/Calib/Gain");
1537   
1538   if(!entry) AliFatal("Calibration object retrieval failed!");
1539   
1540   AliPMDCalibData *calibdata=0;
1541   if (entry) calibdata = (AliPMDCalibData*) entry->GetObject();
1542   
1543   if (!calibdata)  AliFatal("No calibration data from calibration database !");
1544   
1545   return calibdata;
1546 }
1547 //----------------------------------------------------------------------
1548 AliPMDPedestal* AliPMDDigitizer::GetCalibPed() const
1549 {
1550   // The run number will be centralized in AliCDBManager,
1551   // you don't need to set it here!
1552
1553   AliCDBEntry  *entry = AliCDBManager::Instance()->Get("PMD/Calib/Ped");
1554   
1555   if(!entry) AliFatal("Pedestal object retrieval failed!");
1556   
1557   AliPMDPedestal *pedestal=0;
1558   if (entry) pedestal = (AliPMDPedestal*) entry->GetObject();
1559   
1560   if (!pedestal)  AliFatal("No pedestal data from calibration database !");
1561   
1562   return pedestal;
1563 }