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%1000;
 704                       Float_t pedrms     = (Float_t)pedrms1/10.;
 705                       Float_t pedmean    =
 706                           (Float_t) (pedmeanrms - pedrms1)/1000.0;
 707                       //printf("%f %f\n",pedmean, pedrms);
 708                       if (adc > 0.)
 709                       {
 710                           adc += (pedmean + 3.0*pedrms);
 711                           AddDigit(trno,detno,ism,jrow,kcol,adc);
 712                       }
 713                   }
 714               } // column loop
 715           } // row    loop
 716           treeD->Fill();
 717           ResetDigit();
 718       } // supermodule loop
 719   } // detector loop
 720
 721   fPMDLoader->WriteDigits("OVERWRITE");
 722   ResetCellADC();
 723
 724 }
 725 //____________________________________________________________________________
 726
 727
 728 void AliPMDDigitizer::SDigits2Digits(Int_t ievt)
 729 {
 730   // This reads the PMD sdigits tree and converts energy deposition
 731   // in a cell to ADC and stores in the digits tree
 732   //
 733
 734   fRunLoader->GetEvent(ievt);
 735
 736   TTree* treeS = fPMDLoader->TreeS();
 737   AliPMDsdigit  *pmdsdigit;
 738   TBranch *branch = treeS->GetBranch("PMDSDigit");
 739   if(!branch)
 740     {
 741       AliError("PMD Sdigit branch does not exist");
 742       return;
 743     }
 744   if (!fSDigits) fSDigits = new TClonesArray("AliPMDsdigit", 1000);
 745   branch->SetAddress(&fSDigits);
 746
 747   TTree* treeD = fPMDLoader->TreeD();
 748   if (treeD == 0x0)
 749     {
 750       fPMDLoader->MakeTree("D");
 751       treeD = fPMDLoader->TreeD();
 752     }
 753   Int_t bufsize = 16000;
 754   if (!fDigits) fDigits = new TClonesArray("AliPMDdigit", 1000);
 755   treeD->Branch("PMDDigit", &fDigits, bufsize);
 756
 757   Int_t   trno, det, smn;
 758   Int_t   irow, icol;
 759   Float_t edep, adc;
 760
 761   Int_t nmodules = (Int_t) treeS->GetEntries();
 762   AliDebug(1,Form("Number of modules = %d",nmodules));
 763
 764   for (Int_t imodule = 0; imodule < nmodules; imodule++)
 765     {
 766       treeS->GetEntry(imodule);
 767       Int_t nentries = fSDigits->GetLast();
 768       AliDebug(2,Form("Number of entries per module = %d",nentries+1));
 769       for (Int_t ient = 0; ient < nentries+1; ient++)
 770         {
 771           pmdsdigit = (AliPMDsdigit*)fSDigits->UncheckedAt(ient);
 772           trno   = pmdsdigit->GetTrackNumber();
 773           det    = pmdsdigit->GetDetector();
 774           smn    = pmdsdigit->GetSMNumber();
 775           irow   = pmdsdigit->GetRow();
 776           icol   = pmdsdigit->GetColumn();
 777           edep   = pmdsdigit->GetCellEdep();
 778
 779           MeV2ADC(edep,adc);
 780
 781
 782           // To decalibrte the adc values
 783           //
 784           Float_t gain1 = Gain(det,smn,irow,icol);
 785           if (gain1 != 0.)
 786           {
 787               Int_t adcDecalib = (Int_t)(adc/gain1);
 788               adc = (Float_t) adcDecalib;
 789           }
 790           else if(gain1 == 0.)
 791           {
 792               adc = 0.;
 793           }
 794           // Pedestal Decalibration
 795           Int_t   pedmeanrms = fCalibPed->GetPedMeanRms(det,smn,irow,icol);
 796           Int_t   pedrms1    = (Int_t) pedmeanrms%1000;
 797           Float_t pedrms     = (Float_t)pedrms1/10.;
 798           Float_t pedmean    = (Float_t) (pedmeanrms - pedrms1)/1000.0;
 799           //printf("%f %f\n",pedmean, pedrms);
 800           if(adc > 0.)
 801           {
 802               adc += (pedmean + 3.0*pedrms);
 803               AddDigit(trno,det,smn,irow,icol,adc);
 804           }
 805
 806         }
 807       treeD->Fill();
 808       ResetDigit();
 809     }
 810   fPMDLoader->WriteDigits("OVERWRITE");
 811
 812 }
 813 //____________________________________________________________________________
 814 void AliPMDDigitizer::Exec(Option_t *option)
 815 {
 816   // Does the event merging and digitization
 817   const char *cdeb = strstr(option,"deb");
 818   if(cdeb)
 819     {
 820       AliDebug(100," *** PMD Exec is called ***");
 821     }
 822
 823   Int_t ninputs = fManager->GetNinputs();
 824   AliDebug(1,Form("Number of files to be processed = %d",ninputs));
 825   ResetCellADC();
 826
 827   for (Int_t i = 0; i < ninputs; i++)
 828     {
 829       Int_t troffset = fManager->GetMask(i);
 830       MergeSDigits(i, troffset);
 831     }
 832
 833   fRunLoader = AliRunLoader::GetRunLoader(fManager->GetOutputFolderName());
 834   fPMD  = (AliPMD*)gAlice->GetDetector("PMD");
 835   fPMDLoader = fRunLoader->GetLoader("PMDLoader");
 836   if (fPMDLoader == 0x0)
 837     {
 838       AliError("Can not find PMD or PMDLoader");
 839     }
 840   fPMDLoader->LoadDigits("update");
 841   TTree* treeD = fPMDLoader->TreeD();
 842   if (treeD == 0x0)
 843     {
 844       fPMDLoader->MakeTree("D");
 845       treeD = fPMDLoader->TreeD();
 846     }
 847   Int_t bufsize = 16000;
 848   if (!fDigits) fDigits = new TClonesArray("AliPMDdigit", 1000);
 849   treeD->Branch("PMDDigit", &fDigits, bufsize);
 850
 851   Float_t adc;
 852   Float_t deltaE = 0.;
 853   Int_t detno = 0;
 854   Int_t trno = 1;
 855
 856   for (Int_t idet = 0; idet < 2; idet++)
 857     {
 858       for (Int_t ism = 0; ism < fgkTotUM; ism++)
 859         {
 860           for (Int_t jrow = 0; jrow < fgkRow; jrow++)
 861             {
 862               for (Int_t kcol = 0; kcol < fgkCol; kcol++)
 863                 {
 864                   if (idet == 0)
 865                     {
 866                       deltaE = fPRE[ism][jrow][kcol];
 867                       trno   = fPRETrackNo[ism][jrow][kcol];
 868                       detno = 0;
 869                     }
 870                   else if (idet == 1)
 871                     {
 872                       deltaE = fCPV[ism][jrow][kcol];
 873                       trno   = fCPVTrackNo[ism][jrow][kcol];
 874                       detno = 1;
 875                     }
 876                   if (deltaE > 0.)
 877                     {
 878                       MeV2ADC(deltaE,adc);
 879
 880                       //
 881                       // Gain decalibration
 882                       //
 883                       Float_t gain1 = Gain(idet,ism,jrow,kcol);
 884
 885                       if (gain1 != 0.)
 886                       {
 887                           Int_t adcDecalib = (Int_t)(adc/gain1);
 888                           adc = (Float_t) adcDecalib;
 889                       }
 890                       else if(gain1 == 0.)
 891                       {
 892                           adc = 0.;
 893                       }
 894                       // Pedestal Decalibration
 895                       Int_t   pedmeanrms =
 896                           fCalibPed->GetPedMeanRms(idet,ism,jrow,kcol);
 897                       Int_t   pedrms1    = (Int_t) pedmeanrms%1000;
 898                       Float_t pedrms     = (Float_t)pedrms1/10.;
 899                       Float_t pedmean    =
 900                           (Float_t) (pedmeanrms - pedrms1)/1000.0;
 901                       //printf("%f %f\n",pedmean, pedrms);
 902                       if (adc > 0.)
 903                       {
 904                           adc += (pedmean + 3.0*pedrms);
 905                           AddDigit(trno,detno,ism,jrow,kcol,adc);
 906                       }
 907
 908                     }
 909                 } // column loop
 910             } // row    loop
 911           treeD->Fill();
 912           ResetDigit();
 913         } // supermodule loop
 914     } // detector loop
 915   fPMDLoader->WriteDigits("OVERWRITE");
 916   fPMDLoader->UnloadDigits();
 917   ResetCellADC();
 918 }
 919 //____________________________________________________________________________
 920 void AliPMDDigitizer::TrackAssignment2CPVCell()
 921 {
 922   // This block assigns the cell id when there are
 923   // multiple tracks in a cell according to the
 924   // energy deposition
 925   // This method added by Ajay
 926   Bool_t jsort = false;
 927
 928   Int_t i, j, k;
 929
 930   Int_t   *status1;
 931   Int_t   *status2;
 932   Int_t   *trnarray;
 933   Float_t *fracEdp;
 934   Float_t *trEdp;
 935
 936   Int_t   ****cpvTrack;
 937   Float_t ****cpvEdep;
 938
 939   cpvTrack = new Int_t ***[fgkTotUM];
 940   cpvEdep  = new Float_t ***[fgkTotUM];
 941   for (i=0; i<fgkTotUM; i++)
 942     {
 943       cpvTrack[i] = new Int_t **[fgkRow];
 944       cpvEdep[i]  = new Float_t **[fgkRow];
 945     }
 946
 947   for (i = 0; i < fgkTotUM; i++)
 948     {
 949       for (j = 0; j < fgkRow; j++)
 950         {
 951           cpvTrack[i][j] = new Int_t *[fgkCol];
 952           cpvEdep[i][j]  = new Float_t *[fgkCol];
 953         }
 954     }
 955   for (i = 0; i < fgkTotUM; i++)
 956     {
 957       for (j = 0; j < fgkRow; j++)
 958         {
 959           for (k = 0; k < fgkCol; k++)
 960             {
 961               Int_t nn = fCPVCounter[i][j][k];
 962               if(nn > 0)
 963                 {
 964                   cpvTrack[i][j][k] = new Int_t[nn];
 965                   cpvEdep[i][j][k] = new Float_t[nn];
 966                 }
 967               else
 968                 {
 969                   nn = 1;
 970                   cpvTrack[i][j][k] = new Int_t[nn];
 971                   cpvEdep[i][j][k] = new Float_t[nn];
 972                 }
 973               fCPVCounter[i][j][k] = 0;
 974             }
 975         }
 976     }
 977
 978
 979   Int_t nentries = fCPVCell.GetEntries();
 980
 981   Int_t   mtrackno, ism, ixp, iyp;
 982   Float_t edep;
 983   for (i = 0; i < nentries; i++)
 984     {
 985       AliPMDcell* cpvcell = (AliPMDcell*)fCPVCell.UncheckedAt(i);
 986
 987       mtrackno = cpvcell->GetTrackNumber();
 988       ism      = cpvcell->GetSMNumber();
 989       ixp      = cpvcell->GetX();
 990       iyp      = cpvcell->GetY();
 991       edep     = cpvcell->GetEdep();
 992       Int_t nn = fCPVCounter[ism][ixp][iyp];
 993       cpvTrack[ism][ixp][iyp][nn] = (Int_t) mtrackno;
 994       cpvEdep[ism][ixp][iyp][nn] = edep;
 995       fCPVCounter[ism][ixp][iyp]++;
 996     }
 997
 998   Int_t iz, il;
 999   Int_t im, ix, iy;
1000   Int_t nn;
1001   for (im=0; im<fgkTotUM; im++)
1002     {
1003       for (ix=0; ix<fgkRow; ix++)
1004         {
1005           for (iy=0; iy<fgkCol; iy++)
1006             {
1007               nn = fCPVCounter[im][ix][iy];
1008               if (nn > 1)
1009                 {
1010                   // This block handles if a cell is fired
1011                   // many times by many tracks
1012                   status1  = new Int_t[nn];
1013                   status2  = new Int_t[nn];
1014                   trnarray = new Int_t[nn];
1015                   for (iz = 0; iz < nn; iz++)
1016                     {
1017                       status1[iz] = cpvTrack[im][ix][iy][iz];
1018                     }
1019                   TMath::Sort(nn,status1,status2,jsort);
1020                   Int_t trackOld = -99999;
1021                   Int_t track, trCount = 0;
1022                   for (iz = 0; iz < nn; iz++)
1023                     {
1024                       track = status1[status2[iz]];
1025                       if (trackOld != track)
1026                         {
1027                           trnarray[trCount] = track;
1028                           trCount++;
1029                         }
1030                       trackOld = track;
1031                     }
1032                   delete [] status1;
1033                   delete [] status2;
1034                   Float_t totEdp = 0.;
1035                   trEdp = new Float_t[trCount];
1036                   fracEdp = new Float_t[trCount];
1037                   for (il = 0; il < trCount; il++)
1038                     {
1039                       trEdp[il] = 0.;
1040                       track = trnarray[il];
1041                       for (iz = 0; iz < nn; iz++)
1042                         {
1043                           if (track == cpvTrack[im][ix][iy][iz])
1044                             {
1045                               trEdp[il] += cpvEdep[im][ix][iy][iz];
1046                             }
1047                         }
1048                       totEdp += trEdp[il];
1049                     }
1050                   Int_t ilOld = 0;
1051                   Float_t fracOld = 0.;
1052
1053                   for (il = 0; il < trCount; il++)
1054                     {
1055                       fracEdp[il] = trEdp[il]/totEdp;
1056                       if (fracOld < fracEdp[il])
1057                         {
1058                           fracOld = fracEdp[il];
1059                           ilOld = il;
1060                         }
1061                     }
1062                   fCPVTrackNo[im][ix][iy] = trnarray[ilOld];
1063                   delete [] fracEdp;
1064                   delete [] trEdp;
1065                   delete [] trnarray;
1066                 }
1067               else if (nn == 1)
1068                 {
1069                   // This only handles if a cell is fired
1070                   // by only one track
1071
1072                   fCPVTrackNo[im][ix][iy] = cpvTrack[im][ix][iy][0];
1073
1074                 }
1075               else if (nn ==0)
1076                 {
1077                   // This is if no cell is fired
1078                   fCPVTrackNo[im][ix][iy] = -999;
1079                 }
1080             } // end of iy
1081         } // end of ix
1082     } // end of im
1083
1084   // Delete all the pointers
1085
1086  for (i = 0; i < fgkTotUM; i++)
1087     {
1088       for (j = 0; j < fgkRow; j++)
1089         {
1090           for (k = 0; k < fgkCol; k++)
1091             {
1092               delete []cpvTrack[i][j][k];
1093               delete []cpvEdep[i][j][k];
1094             }
1095         }
1096     }
1097
1098   for (i = 0; i < fgkTotUM; i++)
1099     {
1100       for (j = 0; j < fgkRow; j++)
1101         {
1102           delete [] cpvTrack[i][j];
1103           delete [] cpvEdep[i][j];
1104         }
1105     }
1106
1107   for (i = 0; i < fgkTotUM; i++)
1108     {
1109       delete [] cpvTrack[i];
1110       delete [] cpvEdep[i];
1111     }
1112   delete [] cpvTrack;
1113   delete [] cpvEdep;
1114
1115   //
1116   // End of the cell id assignment
1117   //
1118 }
1119 //____________________________________________________________________________
1120
1121 void AliPMDDigitizer::MergeSDigits(Int_t filenumber, Int_t troffset)
1122 {
1123   // merging sdigits
1124   fRunLoader = AliRunLoader::GetRunLoader(fManager->GetInputFolderName(filenumber));
1125   fPMDLoader = fRunLoader->GetLoader("PMDLoader");
1126   fPMDLoader->LoadSDigits("read");
1127   TTree* treeS = fPMDLoader->TreeS();
1128   AliPMDsdigit  *pmdsdigit;
1129   TBranch *branch = treeS->GetBranch("PMDSDigit");
1130   if (!fSDigits) fSDigits = new TClonesArray("AliPMDsdigit", 1000);
1131   branch->SetAddress(&fSDigits);
1132
1133   Int_t   itrackno, idet, ism;
1134   Int_t   ixp, iyp;
1135   Float_t edep;
1136   Int_t nmodules = (Int_t) treeS->GetEntries();
1137   AliDebug(1,Form("Number of Modules in the treeS = %d",nmodules));
1138   AliDebug(1,Form("Track Offset = %d",troffset));
1139   for (Int_t imodule = 0; imodule < nmodules; imodule++)
1140     {
1141       treeS->GetEntry(imodule);
1142       Int_t nentries = fSDigits->GetLast();
1143       AliDebug(2,Form("Number of Entries per Module = %d",nentries));
1144       for (Int_t ient = 0; ient < nentries+1; ient++)
1145         {
1146           pmdsdigit = (AliPMDsdigit*)fSDigits->UncheckedAt(ient);
1147           itrackno  = pmdsdigit->GetTrackNumber();
1148           idet      = pmdsdigit->GetDetector();
1149           ism       = pmdsdigit->GetSMNumber();
1150           ixp       = pmdsdigit->GetRow();
1151           iyp       = pmdsdigit->GetColumn();
1152           edep      = pmdsdigit->GetCellEdep();
1153           if (idet == 0)
1154             {
1155               if (fPRE[ism][ixp][iyp] < edep)
1156                 {
1157                   fPRETrackNo[ism][ixp][iyp] = troffset + itrackno;
1158                 }
1159               fPRE[ism][ixp][iyp] += edep;
1160             }
1161           else if (idet == 1)
1162             {
1163               if (fCPV[ism][ixp][iyp] < edep)
1164                 {
1165                   fCPVTrackNo[ism][ixp][iyp] = troffset + itrackno;
1166                 }
1167               fCPV[ism][ixp][iyp] += edep;
1168             }
1169         }
1170     }
1171
1172 }
1173 // ----------------------------------------------------------------------
1174 void AliPMDDigitizer::TrackAssignment2Cell()
1175 {
1176   //
1177   // This block assigns the cell id when there are
1178   // multiple tracks in a cell according to the
1179   // energy deposition
1180   //
1181   Bool_t jsort = false;
1182
1183   Int_t i, j, k;
1184
1185   Int_t   *status1;
1186   Int_t   *status2;
1187   Int_t   *trnarray;
1188   Float_t *fracEdp;
1189   Float_t *trEdp;
1190
1191   Int_t   ****pmdTrack;
1192   Float_t ****pmdEdep;
1193
1194   pmdTrack = new Int_t ***[fgkTotUM];
1195   pmdEdep  = new Float_t ***[fgkTotUM];
1196   for (i=0; i<fgkTotUM; i++)
1197     {
1198       pmdTrack[i] = new Int_t **[fgkRow];
1199       pmdEdep[i]  = new Float_t **[fgkRow];
1200     }
1201
1202   for (i = 0; i < fgkTotUM; i++)
1203     {
1204       for (j = 0; j < fgkRow; j++)
1205         {
1206           pmdTrack[i][j] = new Int_t *[fgkCol];
1207           pmdEdep[i][j]  = new Float_t *[fgkCol];
1208         }
1209     }
1210
1211   for (i = 0; i < fgkTotUM; i++)
1212     {
1213       for (j = 0; j < fgkRow; j++)
1214         {
1215           for (k = 0; k < fgkCol; k++)
1216             {
1217               Int_t nn = fPRECounter[i][j][k];
1218               if(nn > 0)
1219                 {
1220                   pmdTrack[i][j][k] = new Int_t[nn];
1221                   pmdEdep[i][j][k] = new Float_t[nn];
1222                 }
1223               else
1224                 {
1225                   nn = 1;
1226                   pmdTrack[i][j][k] = new Int_t[nn];
1227                   pmdEdep[i][j][k] = new Float_t[nn];
1228                 }
1229               fPRECounter[i][j][k] = 0;
1230             }
1231         }
1232     }
1233
1234
1235   Int_t nentries = fCell.GetEntries();
1236
1237   Int_t   mtrackno, ism, ixp, iyp;
1238   Float_t edep;
1239
1240   for (i = 0; i < nentries; i++)
1241     {
1242       AliPMDcell* cell = (AliPMDcell*)fCell.UncheckedAt(i);
1243
1244       mtrackno = cell->GetTrackNumber();
1245       ism      = cell->GetSMNumber();
1246       ixp      = cell->GetX();
1247       iyp      = cell->GetY();
1248       edep     = cell->GetEdep();
1249       Int_t nn = fPRECounter[ism][ixp][iyp];
1250       pmdTrack[ism][ixp][iyp][nn] = (Int_t) mtrackno;
1251       pmdEdep[ism][ixp][iyp][nn] = edep;
1252       fPRECounter[ism][ixp][iyp]++;
1253     }
1254
1255   Int_t iz, il;
1256   Int_t im, ix, iy;
1257   Int_t nn;
1258
1259   for (im=0; im<fgkTotUM; im++)
1260     {
1261       for (ix=0; ix<fgkRow; ix++)
1262         {
1263           for (iy=0; iy<fgkCol; iy++)
1264             {
1265               nn = fPRECounter[im][ix][iy];
1266               if (nn > 1)
1267                 {
1268                   // This block handles if a cell is fired
1269                   // many times by many tracks
1270                   status1  = new Int_t[nn];
1271                   status2  = new Int_t[nn];
1272                   trnarray = new Int_t[nn];
1273                   for (iz = 0; iz < nn; iz++)
1274                     {
1275                       status1[iz] = pmdTrack[im][ix][iy][iz];
1276                     }
1277                   TMath::Sort(nn,status1,status2,jsort);
1278                   Int_t trackOld = -99999;
1279                   Int_t track, trCount = 0;
1280                   for (iz = 0; iz < nn; iz++)
1281                     {
1282                       track = status1[status2[iz]];
1283                       if (trackOld != track)
1284                         {
1285                           trnarray[trCount] = track;
1286                           trCount++;
1287                         }
1288                       trackOld = track;
1289                     }
1290                   delete [] status1;
1291                   delete [] status2;
1292                   Float_t totEdp = 0.;
1293                   trEdp = new Float_t[trCount];
1294                   fracEdp = new Float_t[trCount];
1295                   for (il = 0; il < trCount; il++)
1296                     {
1297                       trEdp[il] = 0.;
1298                       track = trnarray[il];
1299                       for (iz = 0; iz < nn; iz++)
1300                         {
1301                           if (track == pmdTrack[im][ix][iy][iz])
1302                             {
1303                               trEdp[il] += pmdEdep[im][ix][iy][iz];
1304                             }
1305                         }
1306                       totEdp += trEdp[il];
1307                     }
1308                   Int_t ilOld = 0;
1309                   Float_t fracOld = 0.;
1310
1311                   for (il = 0; il < trCount; il++)
1312                     {
1313                       fracEdp[il] = trEdp[il]/totEdp;
1314                       if (fracOld < fracEdp[il])
1315                         {
1316                           fracOld = fracEdp[il];
1317                           ilOld = il;
1318                         }
1319                     }
1320                   fPRETrackNo[im][ix][iy] = trnarray[ilOld];
1321                   delete [] fracEdp;
1322                   delete [] trEdp;
1323                   delete [] trnarray;
1324                 }
1325               else if (nn == 1)
1326                 {
1327                   // This only handles if a cell is fired
1328                   // by only one track
1329
1330                   fPRETrackNo[im][ix][iy] = pmdTrack[im][ix][iy][0];
1331
1332                 }
1333               else if (nn ==0)
1334                 {
1335                   // This is if no cell is fired
1336                   fPRETrackNo[im][ix][iy] = -999;
1337                 }
1338             } // end of iy
1339         } // end of ix
1340     } // end of im
1341
1342   // Delete all the pointers
1343
1344   for (i = 0; i < fgkTotUM; i++)
1345     {
1346       for (j = 0; j < fgkRow; j++)
1347         {
1348           for (k = 0; k < fgkCol; k++)
1349             {
1350               delete [] pmdTrack[i][j][k];
1351               delete [] pmdEdep[i][j][k];
1352             }
1353         }
1354     }
1355
1356   for (i = 0; i < fgkTotUM; i++)
1357     {
1358       for (j = 0; j < fgkRow; j++)
1359         {
1360           delete [] pmdTrack[i][j];
1361           delete [] pmdEdep[i][j];
1362         }
1363     }
1364
1365   for (i = 0; i < fgkTotUM; i++)
1366     {
1367       delete [] pmdTrack[i];
1368       delete [] pmdEdep[i];
1369     }
1370   delete [] pmdTrack;
1371   delete [] pmdEdep;
1372   //
1373   // End of the cell id assignment
1374   //
1375 }
1376 //____________________________________________________________________________
1377 void AliPMDDigitizer::MeV2ADC(Float_t mev, Float_t & adc) const
1378 {
1379   // This converts the simulated edep to ADC according to the
1380   // Test Beam Data
1381   //PS Test in September 2003 and 2006
1382   // KeV - ADC conversion for 12bit ADC
1383   // Modified by Ajay
1384   const Float_t kConstant   = 9.0809;
1385   const Float_t kErConstant = 1.6763;
1386   const Float_t kSlope      = 128.348;
1387   const Float_t kErSlope    = 0.4703;
1388
1389   Float_t cons   = gRandom->Gaus(kConstant,kErConstant);
1390   Float_t slop   = gRandom->Gaus(kSlope,kErSlope);
1391
1392   Float_t adc12bit = slop*mev*0.001 + cons;
1393
1394   if(adc12bit < 1600.0)
1395     {
1396       adc = (Float_t) adc12bit;
1397     }
1398   else if (adc12bit >= 1600.0)
1399     {
1400       adc = 1600.0;
1401     }
1402 }
1403 //____________________________________________________________________________
1404 void AliPMDDigitizer::AddSDigit(Int_t trnumber, Int_t det, Int_t smnumber,
1405                                 Int_t irow, Int_t icol, Float_t adc)
1406 {
1407   // Add SDigit
1408   //
1409   if (!fSDigits) fSDigits = new TClonesArray("AliPMDsdigit", 1000);
1410   TClonesArray &lsdigits = *fSDigits;
1411   new(lsdigits[fNsdigit++])  AliPMDsdigit(trnumber,det,smnumber,irow,icol,adc);
1412 }
1413 //____________________________________________________________________________
1414
1415 void AliPMDDigitizer::AddDigit(Int_t trnumber, Int_t det, Int_t smnumber,
1416                                Int_t irow, Int_t icol, Float_t adc)
1417 {
1418   // Add Digit
1419   //
1420   if (!fDigits) fDigits = new TClonesArray("AliPMDdigit", 1000);
1421   TClonesArray &ldigits = *fDigits;
1422   new(ldigits[fNdigit++]) AliPMDdigit(trnumber,det,smnumber,irow,icol,adc);
1423 }
1424 //____________________________________________________________________________
1425
1426 void AliPMDDigitizer::SetZPosition(Float_t zpos)
1427 {
1428   fZPos = zpos;
1429 }
1430 //____________________________________________________________________________
1431 Float_t AliPMDDigitizer::GetZPosition() const
1432 {
1433   return fZPos;
1434 }
1435 //____________________________________________________________________________
1436
1437 void AliPMDDigitizer::ResetCell()
1438 {
1439   // clears the cell array and also the counter
1440   //  for each cell
1441   //
1442   fCPVCell.Delete();
1443   fCell.Delete();
1444   for (Int_t i = 0; i < fgkTotUM; i++)
1445     {
1446       for (Int_t j = 0; j < fgkRow; j++)
1447         {
1448           for (Int_t k = 0; k < fgkCol; k++)
1449             {
1450               fCPVCounter[i][j][k] = 0;
1451               fPRECounter[i][j][k] = 0;
1452             }
1453         }
1454     }
1455 }
1456 //____________________________________________________________________________
1457 void AliPMDDigitizer::ResetSDigit()
1458 {
1459   // Clears SDigits
1460   fNsdigit = 0;
1461   if (fSDigits) fSDigits->Delete();
1462 }
1463 //____________________________________________________________________________
1464 void AliPMDDigitizer::ResetDigit()
1465 {
1466   // Clears Digits
1467   fNdigit = 0;
1468   if (fDigits) fDigits->Delete();
1469 }
1470 //____________________________________________________________________________
1471
1472 void AliPMDDigitizer::ResetCellADC()
1473 {
1474   // Clears individual cells edep and track number
1475   for (Int_t i = 0; i < fgkTotUM; i++)
1476     {
1477       for (Int_t j = 0; j < fgkRow; j++)
1478         {
1479           for (Int_t k = 0; k < fgkCol; k++)
1480             {
1481               fCPV[i][j][k] = 0.;
1482               fPRE[i][j][k] = 0.;
1483               fCPVTrackNo[i][j][k] = 0;
1484               fPRETrackNo[i][j][k] = 0;
1485             }
1486         }
1487     }
1488 }
1489 //____________________________________________________________________________
1490
1491 void AliPMDDigitizer::UnLoad(Option_t *option)
1492 {
1493   // Unloads all the root files
1494   //
1495   const char *cS = strstr(option,"S");
1496   const char *cD = strstr(option,"D");
1497
1498   fRunLoader->UnloadgAlice();
1499   fRunLoader->UnloadHeader();
1500   fRunLoader->UnloadKinematics();
1501
1502   if (cS)
1503     {
1504       fPMDLoader->UnloadHits();
1505     }
1506   if (cD)
1507     {
1508       fPMDLoader->UnloadHits();
1509       fPMDLoader->UnloadSDigits();
1510     }
1511 }
1512
1513 //----------------------------------------------------------------------
1514 Float_t AliPMDDigitizer::Gain(Int_t det, Int_t smn, Int_t row, Int_t col) const
1515 {
1516   // returns of the gain of the cell
1517   // Added this method by ZA
1518
1519   //cout<<" I am here in gain "<<fCalibData<< "smn,row, col "<<smn
1520   //<<" "<<row<<" "<<col<<endl;
1521
1522   if(!fCalibGain) {
1523     AliError("No calibration data loaded from CDB!!!");
1524     return 1;
1525   }
1526
1527   Float_t GainFact;
1528   GainFact = fCalibGain->GetGainFact(det,smn,row,col);
1529   return GainFact;
1530 }
1531 //----------------------------------------------------------------------
1532 AliPMDCalibData* AliPMDDigitizer::GetCalibGain() const
1533 {
1534   // The run number will be centralized in AliCDBManager,
1535   // you don't need to set it here!
1536   // Added this method by ZA
1537   // Cleaned up by Alberto
1538   AliCDBEntry  *entry = AliCDBManager::Instance()->Get("PMD/Calib/Gain");
1539
1540   if(!entry) AliFatal("Calibration object retrieval failed!");
1541
1542   AliPMDCalibData *calibdata=0;
1543   if (entry) calibdata = (AliPMDCalibData*) entry->GetObject();
1544
1545   if (!calibdata)  AliFatal("No calibration data from calibration database !");
1546
1547   return calibdata;
1548 }
1549 //----------------------------------------------------------------------
1550 AliPMDPedestal* AliPMDDigitizer::GetCalibPed() const
1551 {
1552   // The run number will be centralized in AliCDBManager,
1553   // you don't need to set it here!
1554
1555   AliCDBEntry  *entry = AliCDBManager::Instance()->Get("PMD/Calib/Ped");
1556
1557   if(!entry) AliFatal("Pedestal object retrieval failed!");
1558
1559   AliPMDPedestal *pedestal=0;
1560   if (entry) pedestal = (AliPMDPedestal*) entry->GetObject();
1561
1562   if (!pedestal)  AliFatal("No pedestal data from calibration database !");
1563
1564   return pedestal;
1565 }