PMD/AliPMDClusterFinder.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 //                                                     //
  18 //           Date   : August 05 2003                   //
  19 //  This reads the file PMD.digits.root(TreeD),        //
  20 //  calls the Clustering algorithm and stores the      //
  21 //  clustering output in PMD.RecPoints.root(TreeR)     //
  22 //                                                     //
  23 //-----------------------------------------------------//
  24
  25 #include <Riostream.h>
  26 #include <TTree.h>
  27 #include <TObjArray.h>
  28 #include <TClonesArray.h>
  29
  30 #include "AliLog.h"
  31 #include "AliRunLoader.h"
  32 #include "AliLoader.h"
  33 #include "AliRawReader.h"
  34
  35 #include "AliPMDdigit.h"
  36 #include "AliPMDClusterFinder.h"
  37 #include "AliPMDClustering.h"
  38 #include "AliPMDClusteringV1.h"
  39 #include "AliPMDcluster.h"
  40 #include "AliPMDrecpoint1.h"
  41 #include "AliPMDrechit.h"
  42 #include "AliPMDRawStream.h"
  43 #include "AliPMDCalibData.h"
  44 #include "AliPMDPedestal.h"
  45 #include "AliPMDddldata.h"
  46
  47 #include "AliDAQ.h"
  48 #include "AliCDBManager.h"
  49 #include "AliCDBEntry.h"
  50
  51
  52
  53 ClassImp(AliPMDClusterFinder)
  54
  55 AliPMDClusterFinder::AliPMDClusterFinder():
  56   fRunLoader(0),
  57   fPMDLoader(0),
  58   fCalibGain(GetCalibGain()),
  59   fCalibPed(GetCalibPed()),
  60   fTreeD(0),
  61   fTreeR(0),
  62   fDigits(new TClonesArray("AliPMDdigit", 1000)),
  63   fRecpoints(new TClonesArray("AliPMDrecpoint1", 1000)),
  64   fRechits(new TClonesArray("AliPMDrechit", 1000)),
  65   fNpoint(0),
  66   fNhit(0),
  67   fDetNo(0),
  68   fEcut(0.)
  69 {
  70 //
  71 // Constructor
  72 //
  73 }
  74 // ------------------------------------------------------------------------- //
  75 AliPMDClusterFinder::AliPMDClusterFinder(AliRunLoader* runLoader):
  76   fRunLoader(runLoader),
  77   fPMDLoader(runLoader->GetLoader("PMDLoader")),
  78   fCalibGain(GetCalibGain()),
  79   fCalibPed(GetCalibPed()),
  80   fTreeD(0),
  81   fTreeR(0),
  82   fDigits(new TClonesArray("AliPMDdigit", 1000)),
  83   fRecpoints(new TClonesArray("AliPMDrecpoint1", 1000)),
  84   fRechits(new TClonesArray("AliPMDrechit", 1000)),
  85   fNpoint(0),
  86   fNhit(0),
  87   fDetNo(0),
  88   fEcut(0.)
  89 {
  90 //
  91 // Constructor
  92 //
  93 }
  94 // ------------------------------------------------------------------------- //
  95 AliPMDClusterFinder::AliPMDClusterFinder(const AliPMDClusterFinder & finder):
  96   TObject(finder),
  97   fRunLoader(0),
  98   fPMDLoader(0),
  99   fCalibGain(GetCalibGain()),
 100   fCalibPed(GetCalibPed()),
 101   fTreeD(0),
 102   fTreeR(0),
 103   fDigits(NULL),
 104   fRecpoints(NULL),
 105   fRechits(NULL),
 106   fNpoint(0),
 107   fNhit(0),
 108   fDetNo(0),
 109   fEcut(0.)
 110 {
 111   // copy constructor
 112   AliError("Copy constructor not allowed");
 113 }
 114 // ------------------------------------------------------------------------- //
 115 AliPMDClusterFinder &AliPMDClusterFinder::operator=(const AliPMDClusterFinder & /*finder*/)
 116 {
 117  // assignment op
 118   AliError("Assignment Operator not allowed");
 119   return *this;
 120 }
 121 // ------------------------------------------------------------------------- //
 122 AliPMDClusterFinder::~AliPMDClusterFinder()
 123 {
 124   // Destructor
 125   if (fDigits)
 126     {
 127         fDigits->Clear();
 128     }
 129   if (fRecpoints)
 130     {
 131       fRecpoints->Clear();
 132     }
 133   if (fRechits)
 134     {
 135       fRechits->Clear();
 136     }
 137
 138 }
 139 // ------------------------------------------------------------------------- //
 140
 141 void AliPMDClusterFinder::Digits2RecPoints(Int_t ievt)
 142 {
 143   // Converts digits to recpoints after running clustering
 144   // algorithm on CPV plane and PREshower plane
 145   //
 146
 147   Int_t    det = 0,smn = 0;
 148   Int_t    xpos,ypos;
 149   Float_t  adc;
 150   Int_t    ismn;
 151   Int_t    idet;
 152   Float_t  clusdata[6];
 153
 154   TObjArray *pmdcont = new TObjArray();
 155   AliPMDClustering *pmdclust = new AliPMDClusteringV1();
 156
 157   pmdclust->SetEdepCut(fEcut);
 158
 159   fRunLoader->GetEvent(ievt);
 160
 161
 162   fTreeD = fPMDLoader->TreeD();
 163   if (fTreeD == 0x0)
 164     {
 165       AliFatal("AliPMDClusterFinder: Can not get TreeD");
 166
 167     }
 168   AliPMDdigit  *pmddigit;
 169   TBranch *branch = fTreeD->GetBranch("PMDDigit");
 170   branch->SetAddress(&fDigits);
 171
 172   ResetRecpoint();
 173
 174   fTreeR = fPMDLoader->TreeR();
 175   if (fTreeR == 0x0)
 176     {
 177       fPMDLoader->MakeTree("R");
 178       fTreeR = fPMDLoader->TreeR();
 179     }
 180
 181   Int_t bufsize = 16000;
 182   TBranch * branch1 = fTreeR->Branch("PMDRecpoint", &fRecpoints, bufsize);
 183   TBranch * branch2 = fTreeR->Branch("PMDRechit", &fRechits, bufsize);
 184
 185   Int_t nmodules = (Int_t) fTreeD->GetEntries();
 186
 187   for (Int_t imodule = 0; imodule < nmodules; imodule++)
 188     {
 189       ResetCellADC();
 190       fTreeD->GetEntry(imodule);
 191       Int_t nentries = fDigits->GetLast();
 192       for (Int_t ient = 0; ient < nentries+1; ient++)
 193         {
 194           pmddigit = (AliPMDdigit*)fDigits->UncheckedAt(ient);
 195
 196           det    = pmddigit->GetDetector();
 197           smn    = pmddigit->GetSMNumber();
 198           xpos   = pmddigit->GetRow();
 199           ypos   = pmddigit->GetColumn();
 200           adc    = pmddigit->GetADC();
 201
 202           // CALIBRATION
 203           Float_t gain = fCalibGain->GetGainFact(det,smn,xpos,ypos);
 204           // printf("adc = %d gain = %f\n",adc,gain);
 205
 206           adc = adc*gain;
 207
 208           //Int_t trno   = pmddigit->GetTrackNumber();
 209           fCellADC[xpos][ypos] = (Double_t) adc;
 210         }
 211
 212       idet = det;
 213       ismn = smn;
 214       pmdclust->DoClust(idet,ismn,fCellADC,pmdcont);
 215
 216       Int_t nentries1 = pmdcont->GetEntries();
 217
 218       AliDebug(1,Form("Total number of clusters/module = %d",nentries1));
 219
 220       for (Int_t ient1 = 0; ient1 < nentries1; ient1++)
 221         {
 222           AliPMDcluster *pmdcl = (AliPMDcluster*)pmdcont->UncheckedAt(ient1);
 223           idet        = pmdcl->GetDetector();
 224           ismn        = pmdcl->GetSMN();
 225           clusdata[0] = pmdcl->GetClusX();
 226           clusdata[1] = pmdcl->GetClusY();
 227           clusdata[2] = pmdcl->GetClusADC();
 228           clusdata[3] = pmdcl->GetClusCells();
 229           clusdata[4] = pmdcl->GetClusSigmaX();
 230           clusdata[5] = pmdcl->GetClusSigmaY();
 231
 232           AddRecPoint(idet,ismn,clusdata);
 233
 234           Int_t ncell = (Int_t) clusdata[3];
 235           for(Int_t ihit = 0; ihit < ncell; ihit++)
 236             {
 237               Int_t celldataX = pmdcl->GetClusCellX(ihit);
 238               Int_t celldataY = pmdcl->GetClusCellY(ihit);
 239               AddRecHit(celldataX, celldataY);
 240             }
 241           branch2->Fill();
 242           ResetRechit();
 243         }
 244       pmdcont->Delete();
 245
 246       branch1->Fill();
 247       ResetRecpoint();
 248
 249     } // modules
 250
 251   ResetCellADC();
 252   fPMDLoader = fRunLoader->GetLoader("PMDLoader");
 253   fPMDLoader->WriteRecPoints("OVERWRITE");
 254
 255   //   delete the pointers
 256   delete pmdclust;
 257   delete pmdcont;
 258
 259 }
 260 // ------------------------------------------------------------------------- //
 261
 262 void AliPMDClusterFinder::Digits2RecPoints(TTree *digitsTree,
 263                                            TTree *clustersTree)
 264 {
 265   // Converts digits to recpoints after running clustering
 266   // algorithm on CPV plane and PREshower plane
 267   //
 268
 269   Int_t    det = 0,smn = 0;
 270   Int_t    xpos,ypos;
 271   Float_t  adc;
 272   Int_t    ismn;
 273   Int_t    idet;
 274   Float_t  clusdata[6];
 275
 276   AliPMDcluster *pmdcl = 0x0;
 277
 278   TObjArray *pmdcont = new TObjArray();
 279   AliPMDClustering *pmdclust = new AliPMDClusteringV1();
 280
 281   pmdclust->SetEdepCut(fEcut);
 282
 283   AliPMDdigit  *pmddigit;
 284   TBranch *branch = digitsTree->GetBranch("PMDDigit");
 285   branch->SetAddress(&fDigits);
 286
 287   ResetRecpoint();
 288
 289   Int_t bufsize = 16000;
 290   TBranch * branch1 = clustersTree->Branch("PMDRecpoint", &fRecpoints, bufsize);
 291   TBranch * branch2 = clustersTree->Branch("PMDRechit", &fRechits, bufsize);
 292
 293   Int_t nmodules = (Int_t) digitsTree->GetEntries();
 294
 295   for (Int_t imodule = 0; imodule < nmodules; imodule++)
 296     {
 297
 298       Int_t totADCMod = 0;
 299       ResetCellADC();
 300       digitsTree->GetEntry(imodule);
 301       Int_t nentries = fDigits->GetLast();
 302       for (Int_t ient = 0; ient < nentries+1; ient++)
 303         {
 304           pmddigit = (AliPMDdigit*)fDigits->UncheckedAt(ient);
 305
 306           det    = pmddigit->GetDetector();
 307           smn    = pmddigit->GetSMNumber();
 308           xpos   = pmddigit->GetRow();
 309           ypos   = pmddigit->GetColumn();
 310           adc    = pmddigit->GetADC();
 311
 312           // Pedestal Subtraction
 313           Int_t   pedmeanrms = fCalibPed->GetPedMeanRms(det,smn,xpos,ypos);
 314           Int_t   pedrms1    = (Int_t) pedmeanrms%1000;
 315           Float_t pedrms     = (Float_t)pedrms1/10.;
 316           Float_t pedmean    = (Float_t) (pedmeanrms - pedrms1)/1000.0;
 317           //printf("%f %f\n",pedmean, pedrms);
 318
 319           Float_t adc1 = adc - (pedmean + 3.0*pedrms);
 320
 321           // CALIBRATION
 322           Float_t gain = fCalibGain->GetGainFact(det,smn,xpos,ypos);
 323           // printf("adc = %d gain = %f\n",adc,gain);
 324
 325           adc = adc1*gain;
 326
 327           //Int_t trno   = pmddigit->GetTrackNumber();
 328           fCellADC[xpos][ypos] = (Double_t) adc;
 329
 330           totADCMod += adc;
 331
 332         }
 333
 334       idet = det;
 335       ismn = smn;
 336
 337       if (totADCMod <= 0) continue;
 338
 339       pmdclust->DoClust(idet,ismn,fCellADC,pmdcont);
 340
 341       Int_t nentries1 = pmdcont->GetEntries();
 342
 343       AliDebug(1,Form("Total number of clusters/module = %d",nentries1));
 344
 345       for (Int_t ient1 = 0; ient1 < nentries1; ient1++)
 346         {
 347             pmdcl = (AliPMDcluster*)pmdcont->UncheckedAt(ient1);
 348           idet        = pmdcl->GetDetector();
 349           ismn        = pmdcl->GetSMN();
 350           clusdata[0] = pmdcl->GetClusX();
 351           clusdata[1] = pmdcl->GetClusY();
 352           clusdata[2] = pmdcl->GetClusADC();
 353           clusdata[3] = pmdcl->GetClusCells();
 354           clusdata[4] = pmdcl->GetClusSigmaX();
 355           clusdata[5] = pmdcl->GetClusSigmaY();
 356
 357           AddRecPoint(idet,ismn,clusdata);
 358
 359           Int_t ncell = (Int_t) clusdata[3];
 360           for(Int_t ihit = 0; ihit < ncell; ihit++)
 361             {
 362               Int_t celldataX = pmdcl->GetClusCellX(ihit);
 363               Int_t celldataY = pmdcl->GetClusCellY(ihit);
 364               AddRecHit(celldataX, celldataY);
 365             }
 366           branch2->Fill();
 367           ResetRechit();
 368         }
 369       pmdcont->Delete();
 370
 371       branch1->Fill();
 372       ResetRecpoint();
 373
 374     } // modules
 375
 376   ResetCellADC();
 377
 378   //   delete the pointers
 379   delete pmdclust;
 380   delete pmdcont;
 381
 382 }
 383 // ------------------------------------------------------------------------- //
 384
 385 void AliPMDClusterFinder::Digits2RecPoints(AliRawReader *rawReader,
 386                                            TTree *clustersTree)
 387 {
 388   // Converts RAW data to recpoints after running clustering
 389   // algorithm on CPV and PREshower plane
 390   //
 391   // This method is called at the time of reconstruction
 392
 393
 394     AliPMDddldata *pmdddl = 0x0;
 395     AliPMDcluster *pmdcl  = 0x0;
 396
 397
 398   Float_t  clusdata[6];
 399   TObjArray pmdddlcont;
 400
 401   TObjArray *pmdcont = new TObjArray();
 402
 403   AliPMDClustering *pmdclust = new AliPMDClusteringV1();
 404
 405   pmdclust->SetEdepCut(fEcut);
 406
 407   ResetRecpoint();
 408
 409   Int_t bufsize = 16000;
 410   TBranch *branch1 = clustersTree->Branch("PMDRecpoint", &fRecpoints, bufsize);
 411
 412   TBranch * branch2 = clustersTree->Branch("PMDRechit", &fRechits, bufsize);
 413
 414   const Int_t kRow = 48;
 415   const Int_t kCol = 96;
 416
 417   Int_t idet = 0;
 418   Int_t iSMN = 0;
 419
 420   Int_t indexDDL = -1;
 421   AliPMDRawStream pmdinput(rawReader);
 422
 423   while ((indexDDL = pmdinput.DdlData(&pmdddlcont)) >=0)
 424   {
 425       if (indexDDL < 4)
 426         {
 427           iSMN = 6;
 428         }
 429       else if (indexDDL >= 4)
 430         {
 431           iSMN = 12;
 432         }
 433       Int_t ***precpvADC;
 434       precpvADC = new int **[iSMN];
 435       for (Int_t i=0; i<iSMN; i++) precpvADC[i] = new int *[kRow];
 436       for (Int_t i=0; i<iSMN;i++)
 437         {
 438           for (Int_t j=0; j<kRow; j++) precpvADC[i][j] = new int [kCol];
 439         }
 440       for (Int_t i = 0; i < iSMN; i++)
 441         {
 442           for (Int_t j = 0; j < kRow; j++)
 443             {
 444               for (Int_t k = 0; k < kCol; k++)
 445                 {
 446                   precpvADC[i][j][k] = 0;
 447                 }
 448             }
 449         }
 450       ResetCellADC();
 451
 452       Int_t indexsmn = 0;
 453       Int_t ientries = pmdddlcont.GetEntries();
 454       for (Int_t ient = 0; ient < ientries; ient++)
 455         {
 456             pmdddl = (AliPMDddldata*)pmdddlcont.UncheckedAt(ient);
 457
 458           Int_t det = pmdddl->GetDetector();
 459           Int_t smn = pmdddl->GetSMN();
 460           //Int_t mcm = pmdddl->GetMCM();
 461           //Int_t chno = pmdddl->GetChannel();
 462           Int_t row = pmdddl->GetRow();
 463           Int_t col = pmdddl->GetColumn();
 464           Int_t sig = pmdddl->GetSignal();
 465
 466           // Pedestal Subtraction
 467           Int_t   pedmeanrms = fCalibPed->GetPedMeanRms(det,smn,row,col);
 468           Int_t   pedrms1    = (Int_t) pedmeanrms%1000;
 469           Float_t pedrms     = (Float_t)pedrms1/10.;
 470           Float_t pedmean    = (Float_t) (pedmeanrms - pedrms1)/1000.0;
 471
 472           //printf("%f %f\n",pedmean, pedrms);
 473
 474           // Float_t sig1 = (Float_t) sig;
 475           Float_t sig1 = (Float_t) sig - (pedmean + 3.0*pedrms);
 476
 477           // CALIBRATION
 478           Float_t gain = fCalibGain->GetGainFact(det,smn,row,col);
 479           //printf("sig = %d gain = %f\n",sig,gain);
 480           sig = (Int_t) (sig1*gain);
 481
 482           if (indexDDL < 4)
 483             {
 484               if (det != 0)
 485                 AliError(Form("*DDL %d and Detector NUMBER %d NOT MATCHING *",
 486                               indexDDL, det));
 487               indexsmn = smn - indexDDL * 6;
 488             }
 489           else if (indexDDL == 4)
 490             {
 491               if (det != 1)
 492                 AliError(Form("*DDL %d and Detector NUMBER %d NOT MATCHING *",
 493                               indexDDL, det));
 494               if (smn < 6)
 495                 {
 496                   indexsmn = smn;
 497                 }
 498               else if (smn >= 18 && smn < 24)
 499                 {
 500                   indexsmn = smn - 12;
 501                 }
 502             }
 503           else if (indexDDL == 5)
 504             {
 505               if (det != 1)
 506                 AliError(Form("*DDL %d and Detector NUMBER %d NOT MATCHING *",
 507                               indexDDL, det));
 508               if (smn >= 6 && smn < 12)
 509                 {
 510                   indexsmn = smn - 6;
 511                 }
 512               else if (smn >= 12 && smn < 18)
 513                 {
 514                   indexsmn = smn - 6;
 515                 }
 516             }
 517           precpvADC[indexsmn][row][col] = sig;
 518         }
 519
 520       pmdddlcont.Delete();
 521
 522       Int_t totAdcMod = 0;
 523
 524       Int_t ismn = 0;
 525       for (indexsmn = 0; indexsmn < iSMN; indexsmn++)
 526         {
 527           ResetCellADC();
 528           totAdcMod = 0;
 529           for (Int_t irow = 0; irow < kRow; irow++)
 530             {
 531               for (Int_t icol = 0; icol < kCol; icol++)
 532                 {
 533                   fCellADC[irow][icol] =
 534                     (Double_t) precpvADC[indexsmn][irow][icol];
 535                   totAdcMod += precpvADC[indexsmn][irow][icol];
 536                 } // row
 537             }     // col
 538
 539           if (indexDDL < 4)
 540             {
 541               ismn = indexsmn + indexDDL * 6;
 542               idet = 0;
 543             }
 544           else if (indexDDL == 4)
 545             {
 546               if (indexsmn < 6)
 547                 {
 548                   ismn = indexsmn;
 549                 }
 550               else if (indexsmn >= 6 && indexsmn < 12)
 551                 {
 552                   ismn = indexsmn + 12;
 553                 }
 554               idet = 1;
 555             }
 556           else if (indexDDL == 5)
 557             {
 558               if (indexsmn < 6)
 559                 {
 560                   ismn = indexsmn + 6;
 561                 }
 562               else if (indexsmn >= 6 && indexsmn < 12)
 563                 {
 564                   ismn = indexsmn + 6;
 565                 }
 566               idet = 1;
 567             }
 568
 569           if (totAdcMod <= 0) continue;
 570           pmdclust->DoClust(idet,ismn,fCellADC,pmdcont);
 571           Int_t nentries1 = pmdcont->GetEntries();
 572
 573           AliDebug(1,Form("Total number of clusters/module = %d",nentries1));
 574
 575           for (Int_t ient1 = 0; ient1 < nentries1; ient1++)
 576             {
 577                 pmdcl = (AliPMDcluster*)pmdcont->UncheckedAt(ient1);
 578               idet        = pmdcl->GetDetector();
 579               ismn        = pmdcl->GetSMN();
 580               clusdata[0] = pmdcl->GetClusX();
 581               clusdata[1] = pmdcl->GetClusY();
 582               clusdata[2] = pmdcl->GetClusADC();
 583               clusdata[3] = pmdcl->GetClusCells();
 584               clusdata[4] = pmdcl->GetClusSigmaX();
 585               clusdata[5] = pmdcl->GetClusSigmaY();
 586
 587               AddRecPoint(idet,ismn,clusdata);
 588
 589               Int_t ncell = (Int_t) clusdata[3];
 590               for(Int_t ihit = 0; ihit < ncell; ihit++)
 591                 {
 592                   Int_t celldataX = pmdcl->GetClusCellX(ihit);
 593                   Int_t celldataY = pmdcl->GetClusCellY(ihit);
 594                   AddRecHit(celldataX, celldataY);
 595                 }
 596               branch2->Fill();
 597               ResetRechit();
 598
 599             }
 600           pmdcont->Delete();
 601
 602           branch1->Fill();
 603           ResetRecpoint();
 604
 605
 606         } // smn
 607
 608       for (Int_t i=0; i<iSMN; i++)
 609         {
 610           for (Int_t j=0; j<kRow; j++) delete [] precpvADC[i][j];
 611         }
 612       for (Int_t i=0; i<iSMN; i++) delete [] precpvADC[i];
 613       delete precpvADC;
 614
 615     } // DDL Loop
 616
 617   ResetCellADC();
 618
 619   //   delete the pointers
 620   delete pmdclust;
 621   delete pmdcont;
 622
 623 }
 624 // ------------------------------------------------------------------------- //
 625
 626 void AliPMDClusterFinder::Digits2RecPoints(Int_t ievt, AliRawReader *rawReader)
 627 {
 628   // Converts RAW data to recpoints after running clustering
 629   // algorithm on CPV and PREshower plane
 630   //
 631
 632   Float_t  clusdata[6];
 633   TObjArray pmdddlcont;
 634   TObjArray *pmdcont = new TObjArray();
 635
 636   AliPMDClustering *pmdclust = new AliPMDClusteringV1();
 637
 638   pmdclust->SetEdepCut(fEcut);
 639
 640   fRunLoader->GetEvent(ievt);
 641
 642   ResetRecpoint();
 643
 644   fTreeR = fPMDLoader->TreeR();
 645   if (fTreeR == 0x0)
 646     {
 647       fPMDLoader->MakeTree("R");
 648       fTreeR = fPMDLoader->TreeR();
 649     }
 650   Int_t bufsize = 16000;
 651   TBranch *branch1 = fTreeR->Branch("PMDRecpoint", &fRecpoints, bufsize);
 652   TBranch *branch2 = fTreeR->Branch("PMDRechit", &fRechits, bufsize);
 653
 654   const Int_t kRow = 48;
 655   const Int_t kCol = 96;
 656
 657   Int_t idet = 0;
 658   Int_t iSMN = 0;
 659
 660   AliPMDRawStream pmdinput(rawReader);
 661   Int_t indexDDL = -1;
 662
 663   while ((indexDDL = pmdinput.DdlData(&pmdddlcont)) >=0) {
 664
 665       if (indexDDL < 4)
 666         {
 667           iSMN = 6;
 668         }
 669       else if (indexDDL >= 4)
 670         {
 671           iSMN = 12;
 672         }
 673       Int_t ***precpvADC;
 674       precpvADC = new int **[iSMN];
 675       for (Int_t i=0; i<iSMN; i++) precpvADC[i] = new int *[kRow];
 676       for (Int_t i=0; i<iSMN;i++)
 677         {
 678           for (Int_t j=0; j<kRow; j++) precpvADC[i][j] = new int [kCol];
 679         }
 680       for (Int_t i = 0; i < iSMN; i++)
 681         {
 682           for (Int_t j = 0; j < kRow; j++)
 683             {
 684               for (Int_t k = 0; k < kCol; k++)
 685                 {
 686                   precpvADC[i][j][k] = 0;
 687                 }
 688             }
 689         }
 690       ResetCellADC();
 691
 692
 693       Int_t indexsmn = 0;
 694       Int_t ientries = pmdddlcont.GetEntries();
 695       for (Int_t ient = 0; ient < ientries; ient++)
 696         {
 697           AliPMDddldata *pmdddl = (AliPMDddldata*)pmdddlcont.UncheckedAt(ient);
 698
 699           Int_t det = pmdddl->GetDetector();
 700           Int_t smn = pmdddl->GetSMN();
 701           //Int_t mcm = pmdddl->GetMCM();
 702           //Int_t chno = pmdddl->GetChannel();
 703           Int_t row = pmdddl->GetRow();
 704           Int_t col = pmdddl->GetColumn();
 705           Int_t sig = pmdddl->GetSignal();
 706
 707           // Pedestal Subtraction
 708           Int_t   pedmeanrms = fCalibPed->GetPedMeanRms(det,smn,row,col);
 709           Int_t   pedrms1    = (Int_t) pedmeanrms%1000;
 710           Float_t pedrms     = (Float_t)pedrms1/10.;
 711           Float_t pedmean    = (Float_t) (pedmeanrms - pedrms1)/1000.0;
 712
 713           //printf("%f %f\n",pedmean, pedrms);
 714
 715           //Float_t sig1 = (Float_t) sig;
 716           Float_t sig1 = (Float_t) sig - (pedmean + 3.0*pedrms);
 717           // CALIBRATION
 718           Float_t gain = fCalibGain->GetGainFact(det,smn,row,col);
 719
 720           //printf("sig = %d gain = %f\n",sig,gain);
 721           sig = (Int_t) (sig1*gain);
 722
 723
 724           if (indexDDL < 4)
 725             {
 726               if (det != 0)
 727                 AliError(Form("*DDL %d and Detector NUMBER %d NOT MATCHING *",
 728                               indexDDL, det));
 729               indexsmn = smn - indexDDL * 6;
 730             }
 731           else if (indexDDL == 4)
 732             {
 733               if (det != 1)
 734                 AliError(Form("*DDL %d and Detector NUMBER %d NOT MATCHING *",
 735                               indexDDL, det));
 736               if (smn < 6)
 737                 {
 738                   indexsmn = smn;
 739                 }
 740               else if (smn >= 18 && smn < 24)
 741                 {
 742                   indexsmn = smn - 12;
 743                 }
 744             }
 745           else if (indexDDL == 5)
 746             {
 747               if (det != 1)
 748                 AliError(Form("*DDL %d and Detector NUMBER %d NOT MATCHING *",
 749                               indexDDL, det));
 750               if (smn >= 6 && smn < 12)
 751                 {
 752                   indexsmn = smn - 6;
 753                 }
 754               else if (smn >= 12 && smn < 18)
 755                 {
 756                   indexsmn = smn - 6;
 757                 }
 758             }
 759           precpvADC[indexsmn][row][col] = sig;
 760
 761         }
 762
 763       pmdddlcont.Delete();
 764
 765       Int_t ismn = 0;
 766       for (indexsmn = 0; indexsmn < iSMN; indexsmn++)
 767         {
 768           ResetCellADC();
 769           for (Int_t irow = 0; irow < kRow; irow++)
 770             {
 771               for (Int_t icol = 0; icol < kCol; icol++)
 772                 {
 773                   fCellADC[irow][icol] =
 774                     (Double_t) precpvADC[indexsmn][irow][icol];
 775                 } // row
 776             }     // col
 777
 778
 779           if (indexDDL < 4)
 780             {
 781               ismn = indexsmn + indexDDL * 6;
 782               idet = 0;
 783             }
 784           else if (indexDDL == 4)
 785             {
 786               if (indexsmn < 6)
 787                 {
 788                   ismn = indexsmn;
 789                 }
 790               else if (indexsmn >= 6 && indexsmn < 12)
 791                 {
 792                   ismn = indexsmn + 12;
 793                 }
 794               idet = 1;
 795             }
 796           else if (indexDDL == 5)
 797             {
 798               if (indexsmn < 6)
 799                 {
 800                   ismn = indexsmn + 6;
 801                 }
 802               else if (indexsmn >= 6 && indexsmn < 12)
 803                 {
 804                   ismn = indexsmn + 6;
 805                 }
 806               idet = 1;
 807             }
 808
 809           pmdclust->DoClust(idet,ismn,fCellADC,pmdcont);
 810           Int_t nentries1 = pmdcont->GetEntries();
 811
 812           AliDebug(1,Form("Total number of clusters/module = %d",nentries1));
 813
 814           for (Int_t ient1 = 0; ient1 < nentries1; ient1++)
 815             {
 816               AliPMDcluster *pmdcl =
 817                 (AliPMDcluster*)pmdcont->UncheckedAt(ient1);
 818               idet        = pmdcl->GetDetector();
 819               ismn        = pmdcl->GetSMN();
 820               clusdata[0] = pmdcl->GetClusX();
 821               clusdata[1] = pmdcl->GetClusY();
 822               clusdata[2] = pmdcl->GetClusADC();
 823               clusdata[3] = pmdcl->GetClusCells();
 824               clusdata[4] = pmdcl->GetClusSigmaX();
 825               clusdata[5] = pmdcl->GetClusSigmaY();
 826
 827               AddRecPoint(idet,ismn,clusdata);
 828
 829               Int_t ncell = (Int_t) clusdata[3];
 830               for(Int_t ihit = 0; ihit < ncell; ihit++)
 831                 {
 832                   Int_t celldataX = pmdcl->GetClusCellX(ihit);
 833                   Int_t celldataY = pmdcl->GetClusCellY(ihit);
 834                   AddRecHit(celldataX, celldataY);
 835                 }
 836               branch2->Fill();
 837               ResetRechit();
 838
 839             }
 840           pmdcont->Delete();
 841
 842           branch1->Fill();
 843           ResetRecpoint();
 844
 845
 846         } // smn
 847
 848       for (Int_t i=0; i<iSMN; i++)
 849         {
 850           for (Int_t j=0; j<kRow; j++) delete [] precpvADC[i][j];
 851         }
 852       for (Int_t i=0; i<iSMN; i++) delete [] precpvADC[i];
 853       delete precpvADC;
 854     } // DDL Loop
 855
 856   ResetCellADC();
 857
 858   fPMDLoader = fRunLoader->GetLoader("PMDLoader");
 859   fPMDLoader->WriteRecPoints("OVERWRITE");
 860
 861   //   delete the pointers
 862   delete pmdclust;
 863   delete pmdcont;
 864
 865 }
 866 // ------------------------------------------------------------------------- //
 867 void AliPMDClusterFinder::SetCellEdepCut(Float_t ecut)
 868 {
 869   fEcut = ecut;
 870 }
 871 // ------------------------------------------------------------------------- //
 872 void AliPMDClusterFinder::AddRecPoint(Int_t idet,Int_t ismn,Float_t *clusdata)
 873 {
 874   // Add Reconstructed points
 875   //
 876   TClonesArray &lrecpoints = *fRecpoints;
 877   AliPMDrecpoint1 *newrecpoint;
 878   newrecpoint = new AliPMDrecpoint1(idet, ismn, clusdata);
 879   new(lrecpoints[fNpoint++]) AliPMDrecpoint1(newrecpoint);
 880   delete newrecpoint;
 881 }
 882 // ------------------------------------------------------------------------- //
 883 void AliPMDClusterFinder::AddRecHit(Int_t celldataX,Int_t celldataY)
 884 {
 885   // Add associated cell hits to the Reconstructed points
 886   //
 887   TClonesArray &lrechits = *fRechits;
 888   AliPMDrechit *newrechit;
 889   newrechit = new AliPMDrechit(celldataX, celldataY);
 890   new(lrechits[fNhit++]) AliPMDrechit(newrechit);
 891   delete newrechit;
 892 }
 893 // ------------------------------------------------------------------------- //
 894 void AliPMDClusterFinder::ResetCellADC()
 895 {
 896   // Reset the individual cell ADC value to zero
 897   //
 898   for(Int_t irow = 0; irow < fgkRow; irow++)
 899     {
 900       for(Int_t icol = 0; icol < fgkCol; icol++)
 901         {
 902           fCellADC[irow][icol] = 0.;
 903         }
 904     }
 905 }
 906 // ------------------------------------------------------------------------- //
 907
 908 void AliPMDClusterFinder::ResetRecpoint()
 909 {
 910   // Clear the list of reconstructed points
 911   fNpoint = 0;
 912   if (fRecpoints) fRecpoints->Clear();
 913 }
 914 // ------------------------------------------------------------------------- //
 915 void AliPMDClusterFinder::ResetRechit()
 916 {
 917   // Clear the list of reconstructed points
 918   fNhit = 0;
 919   if (fRechits) fRechits->Clear();
 920 }
 921 // ------------------------------------------------------------------------- //
 922 void AliPMDClusterFinder::Load()
 923 {
 924   // Load all the *.root files
 925   //
 926   fPMDLoader->LoadDigits("READ");
 927   fPMDLoader->LoadRecPoints("recreate");
 928 }
 929 // ------------------------------------------------------------------------- //
 930 void AliPMDClusterFinder::LoadClusters()
 931 {
 932   // Load all the *.root files
 933   //
 934   fPMDLoader->LoadRecPoints("recreate");
 935 }
 936 // ------------------------------------------------------------------------- //
 937 void AliPMDClusterFinder::UnLoad()
 938 {
 939   // Unload all the *.root files
 940   //
 941   fPMDLoader->UnloadDigits();
 942   fPMDLoader->UnloadRecPoints();
 943 }
 944 // ------------------------------------------------------------------------- //
 945 void AliPMDClusterFinder::UnLoadClusters()
 946 {
 947   // Unload all the *.root files
 948   //
 949   fPMDLoader->UnloadRecPoints();
 950 }
 951 // ------------------------------------------------------------------------- //
 952
 953 AliPMDCalibData* AliPMDClusterFinder::GetCalibGain() const
 954 {
 955   // The run number will be centralized in AliCDBManager,
 956   // you don't need to set it here!
 957   // Added by ZA
 958   AliCDBEntry  *entry = AliCDBManager::Instance()->Get("PMD/Calib/Gain");
 959
 960   if(!entry)  AliFatal("Calibration object retrieval failed! ");
 961
 962   AliPMDCalibData *calibdata=0;
 963   if (entry) calibdata = (AliPMDCalibData*) entry->GetObject();
 964
 965   if (!calibdata)  AliFatal("No calibration data from calibration database !");
 966
 967   return calibdata;
 968 }
 969
 970 // ------------------------------------------------------------------------- //
 971
 972 AliPMDPedestal* AliPMDClusterFinder::GetCalibPed() const
 973 {
 974   // The run number will be centralized in AliCDBManager,
 975   // you don't need to set it here!
 976   AliCDBEntry  *entry = AliCDBManager::Instance()->Get("PMD/Calib/Ped");
 977
 978   if(!entry) AliFatal("Pedestal object retrieval failed!");
 979
 980   AliPMDPedestal *pedestal = 0;
 981   if (entry) pedestal = (AliPMDPedestal*) entry->GetObject();
 982
 983   if (!pedestal)  AliFatal("No pedestal data from pedestal database !");
 984
 985   return pedestal;
 986 }