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