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