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