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 <TMath.h>
  27 #include <TBRIK.h>
  28 #include <TNode.h>
  29 #include <TTree.h>
  30 #include <TGeometry.h>
  31 #include <TObjArray.h>
  32 #include <TClonesArray.h>
  33 #include <TFile.h>
  34 #include <TNtuple.h>
  35 #include <TParticle.h>
  36
  37 #include "AliRun.h"
  38 #include "AliPMD.h"
  39 #include "AliDetector.h"
  40 #include "AliRunLoader.h"
  41 #include "AliLoader.h"
  42 #include "AliHeader.h"
  43 #include "AliRawReader.h"
  44
  45 #include "AliPMDdigit.h"
  46 #include "AliPMDClusterFinder.h"
  47 #include "AliPMDClustering.h"
  48 #include "AliPMDcluster.h"
  49 #include "AliPMDrecpoint1.h"
  50 #include "AliPMDRawStream.h"
  51
  52 ClassImp(AliPMDClusterFinder)
  53
  54 AliPMDClusterFinder::AliPMDClusterFinder(AliRunLoader* runLoader):
  55   fRunLoader(runLoader),
  56   fPMDLoader(runLoader->GetLoader("PMDLoader")),
  57   fTreeD(0),
  58   fTreeR(0),
  59   fDigits(new TClonesArray("AliPMDdigit", 1000)),
  60   fRecpoints(new TClonesArray("AliPMDrecpoint1", 1000)),
  61   fNpoint(0),
  62   fDebug(0),
  63   fEcut(0.)
  64 {
  65 //
  66 // Constructor
  67 //
  68 }
  69 // ------------------------------------------------------------------------- //
  70 AliPMDClusterFinder::~AliPMDClusterFinder()
  71 {
  72   // Destructor
  73   if (fDigits)
  74     {
  75       fDigits->Delete();
  76       delete fDigits;
  77       fDigits=0;
  78     }
  79   if (fRecpoints)
  80     {
  81       fRecpoints->Delete();
  82       delete fRecpoints;
  83       fRecpoints=0;
  84     }
  85 }
  86 // ------------------------------------------------------------------------- //
  87
  88 void AliPMDClusterFinder::Digits2RecPoints(Int_t ievt)
  89 {
  90   // Converts digits to recpoints after running clustering
  91   // algorithm on CPV plane and PREshower plane
  92   //
  93   Int_t    det = 0,smn = 0;
  94   Int_t    xpos,ypos;
  95   Float_t  adc;
  96   Int_t    ismn;
  97   Int_t    idet;
  98   Float_t  clusdata[5];
  99
 100   TObjArray *pmdcont = new TObjArray();
 101   AliPMDClustering *pmdclust = new AliPMDClustering();
 102   pmdclust->SetDebug(fDebug);
 103   pmdclust->SetEdepCut(fEcut);
 104
 105   fRunLoader->GetEvent(ievt);
 106   //cout << " ***** Beginning::Digits2RecPoints *****" << endl;
 107
 108   fTreeD = fPMDLoader->TreeD();
 109   if (fTreeD == 0x0)
 110     {
 111       cout << " Can not get TreeD" << endl;
 112     }
 113   AliPMDdigit  *pmddigit;
 114   TBranch *branch = fTreeD->GetBranch("PMDDigit");
 115   branch->SetAddress(&fDigits);
 116
 117   ResetRecpoint();
 118
 119   fTreeR = fPMDLoader->TreeR();
 120   if (fTreeR == 0x0)
 121     {
 122       fPMDLoader->MakeTree("R");
 123       fTreeR = fPMDLoader->TreeR();
 124     }
 125
 126   Int_t bufsize = 16000;
 127   fTreeR->Branch("PMDRecpoint", &fRecpoints, bufsize);
 128
 129   Int_t nmodules = (Int_t) fTreeD->GetEntries();
 130
 131   for (Int_t imodule = 0; imodule < nmodules; imodule++)
 132     {
 133       ResetCellADC();
 134       fTreeD->GetEntry(imodule);
 135       Int_t nentries = fDigits->GetLast();
 136       for (Int_t ient = 0; ient < nentries+1; ient++)
 137         {
 138           pmddigit = (AliPMDdigit*)fDigits->UncheckedAt(ient);
 139
 140           det    = pmddigit->GetDetector();
 141           smn    = pmddigit->GetSMNumber();
 142           xpos   = pmddigit->GetRow();
 143           ypos   = pmddigit->GetColumn();
 144           adc    = pmddigit->GetADC();
 145           //Int_t trno   = pmddigit->GetTrackNumber();
 146           fCellADC[xpos][ypos] = (Double_t) adc;
 147         }
 148
 149       idet = det;
 150       ismn = smn;
 151       pmdclust->DoClust(idet,ismn,fCellADC,pmdcont);
 152
 153       Int_t nentries1 = pmdcont->GetEntries();
 154 //      cout << " nentries1 = " << nentries1 << endl;
 155       for (Int_t ient1 = 0; ient1 < nentries1; ient1++)
 156         {
 157           AliPMDcluster *pmdcl = (AliPMDcluster*)pmdcont->UncheckedAt(ient1);
 158           idet        = pmdcl->GetDetector();
 159           ismn        = pmdcl->GetSMN();
 160           clusdata[0] = pmdcl->GetClusX();
 161           clusdata[1] = pmdcl->GetClusY();
 162           clusdata[2] = pmdcl->GetClusADC();
 163           clusdata[3] = pmdcl->GetClusCells();
 164           clusdata[4] = pmdcl->GetClusRadius();
 165
 166           AddRecPoint(idet,ismn,clusdata);
 167         }
 168       pmdcont->Clear();
 169
 170       fTreeR->Fill();
 171       ResetRecpoint();
 172
 173     } // modules
 174
 175   ResetCellADC();
 176   fPMDLoader = fRunLoader->GetLoader("PMDLoader");
 177   fPMDLoader->WriteRecPoints("OVERWRITE");
 178
 179   //   delete the pointers
 180   delete pmdclust;
 181   delete pmdcont;
 182
 183   //  cout << " ***** End::Digits2RecPoints *****" << endl;
 184 }
 185 // ------------------------------------------------------------------------- //
 186
 187 void AliPMDClusterFinder::Digits2RecPoints(Int_t ievt, AliRawReader *rawReader)
 188 {
 189   // Converts RAW data to recpoints after running clustering
 190   // algorithm on CPV and PREshower plane
 191   //
 192
 193   Float_t  clusdata[5];
 194
 195   TObjArray *pmdcont = new TObjArray();
 196   AliPMDClustering *pmdclust = new AliPMDClustering();
 197   pmdclust->SetDebug(fDebug);
 198   pmdclust->SetEdepCut(fEcut);
 199
 200   fRunLoader->GetEvent(ievt);
 201
 202   ResetRecpoint();
 203
 204   fTreeR = fPMDLoader->TreeR();
 205   if (fTreeR == 0x0)
 206     {
 207       fPMDLoader->MakeTree("R");
 208       fTreeR = fPMDLoader->TreeR();
 209     }
 210   Int_t bufsize = 16000;
 211   fTreeR->Branch("PMDRecpoint", &fRecpoints, bufsize);
 212
 213   const Int_t kDDL = 6;
 214   const Int_t kRow = 48;
 215   const Int_t kCol = 96;
 216
 217   Int_t idet = 0;
 218   Int_t iSMN = 0;
 219
 220   for (Int_t indexDDL = 0; indexDDL < kDDL; indexDDL++)
 221     {
 222       if (indexDDL < 4)
 223         {
 224           iSMN = 6;
 225         }
 226       else if (indexDDL >= 4)
 227         {
 228           iSMN = 12;
 229         }
 230       Int_t ***precpvADC;
 231       precpvADC = new int **[iSMN];
 232       for (Int_t i=0; i<iSMN; i++) precpvADC[i] = new int *[kRow];
 233       for (Int_t i=0; i<iSMN;i++)
 234         {
 235           for (Int_t j=0; j<kRow; j++) precpvADC[i][j] = new int [kCol];
 236         }
 237       for (Int_t i = 0; i < iSMN; i++)
 238         {
 239           for (Int_t j = 0; j < kRow; j++)
 240             {
 241               for (Int_t k = 0; k < kCol; k++)
 242                 {
 243                   precpvADC[i][j][k] = 0;
 244                 }
 245             }
 246         }
 247       ResetCellADC();
 248       rawReader->Reset();
 249       AliPMDRawStream pmdinput(rawReader);
 250       rawReader->Select(12, indexDDL, indexDDL);
 251       while(pmdinput.Next())
 252         {
 253           Int_t det = pmdinput.GetDetector();
 254           Int_t smn = pmdinput.GetSMN();
 255           //Int_t mcm = pmdinput.GetMCM();
 256           //Int_t chno = pmdinput.GetChannel();
 257           Int_t row = pmdinput.GetRow();
 258           Int_t col = pmdinput.GetColumn();
 259           Int_t sig = pmdinput.GetSignal();
 260
 261           Int_t indexsmn = 0;
 262
 263           if (indexDDL < 4)
 264             {
 265               if (det != 0)
 266                 printf(" *** DDL %d and Detector NUMBER %d NOT MATCHING *** ",
 267                        indexDDL, det);
 268               indexsmn = smn - indexDDL * 6;
 269             }
 270           else if (indexDDL == 4)
 271             {
 272               if (det != 1)
 273                 printf(" *** DDL %d and Detector NUMBER %d NOT MATCHING *** ",
 274                        indexDDL, det);
 275               if (smn < 6)
 276                 {
 277                   indexsmn = smn;
 278                 }
 279               else if (smn >= 12 && smn < 18)
 280                 {
 281                   indexsmn = smn - 6;
 282                 }
 283             }
 284           else if (indexDDL == 5)
 285             {
 286               if (det != 1)
 287                 printf(" *** DDL %d and Detector NUMBER %d NOT MATCHING *** ",
 288                        indexDDL, det);
 289               if (smn >= 6 && smn < 12)
 290                 {
 291                   indexsmn = smn - 6;
 292                 }
 293               else if (smn >= 18 && smn < 24)
 294                 {
 295                   indexsmn = smn - 12;
 296                 }
 297             }
 298           precpvADC[indexsmn][row][col] = sig;
 299         } // while loop
 300
 301       Int_t ismn = 0;
 302       for (Int_t indexsmn = 0; indexsmn < iSMN; indexsmn++)
 303         {
 304           ResetCellADC();
 305           for (Int_t irow = 0; irow < kRow; irow++)
 306             {
 307               for (Int_t icol = 0; icol < kCol; icol++)
 308                 {
 309                   fCellADC[irow][icol] =
 310                     (Double_t) precpvADC[indexsmn][irow][icol];
 311                 } // row
 312             }     // col
 313           if (indexDDL < 4)
 314             {
 315               ismn = indexsmn + indexDDL * 6;
 316               idet = 0;
 317             }
 318           else if (indexDDL == 4)
 319             {
 320               if (indexsmn < 6)
 321                 {
 322                   ismn = indexsmn;
 323                 }
 324               else if (indexsmn >= 6 && indexsmn < 12)
 325                 {
 326                   ismn = indexsmn + 6;
 327                 }
 328               idet = 1;
 329             }
 330           else if (indexDDL == 5)
 331             {
 332               if (indexsmn < 6)
 333                 {
 334                   ismn = indexsmn + 6;
 335                 }
 336               else if (indexsmn >= 6 && indexsmn < 12)
 337                 {
 338                   ismn = indexsmn + 12;
 339                 }
 340               idet = 1;
 341             }
 342
 343
 344           pmdclust->DoClust(idet,ismn,fCellADC,pmdcont);
 345           Int_t nentries1 = pmdcont->GetEntries();
 346           //      cout << " nentries1 = " << nentries1 << endl;
 347           for (Int_t ient1 = 0; ient1 < nentries1; ient1++)
 348             {
 349               AliPMDcluster *pmdcl =
 350                 (AliPMDcluster*)pmdcont->UncheckedAt(ient1);
 351               idet        = pmdcl->GetDetector();
 352               ismn        = pmdcl->GetSMN();
 353               clusdata[0] = pmdcl->GetClusX();
 354               clusdata[1] = pmdcl->GetClusY();
 355               clusdata[2] = pmdcl->GetClusADC();
 356               clusdata[3] = pmdcl->GetClusCells();
 357               clusdata[4] = pmdcl->GetClusRadius();
 358
 359               AddRecPoint(idet,ismn,clusdata);
 360             }
 361           pmdcont->Clear();
 362
 363           fTreeR->Fill();
 364           ResetRecpoint();
 365
 366
 367         } // smn
 368
 369       for (Int_t i=0; i<iSMN; i++)
 370         {
 371           for (Int_t j=0; j<kRow; j++) delete [] precpvADC[i][j];
 372         }
 373       for (Int_t i=0; i<iSMN; i++) delete [] precpvADC[i];
 374       delete precpvADC;
 375     } // DDL Loop
 376
 377   ResetCellADC();
 378
 379   fPMDLoader = fRunLoader->GetLoader("PMDLoader");
 380   fPMDLoader->WriteRecPoints("OVERWRITE");
 381
 382   //   delete the pointers
 383   delete pmdclust;
 384   delete pmdcont;
 385
 386   //  cout << " ***** End::Digits2RecPoints :: Raw *****" << endl;
 387 }
 388 // ------------------------------------------------------------------------- //
 389 void AliPMDClusterFinder::SetCellEdepCut(Float_t ecut)
 390 {
 391   fEcut = ecut;
 392 }
 393 // ------------------------------------------------------------------------- //
 394 void AliPMDClusterFinder::SetDebug(Int_t idebug)
 395 {
 396   fDebug = idebug;
 397 }
 398 // ------------------------------------------------------------------------- //
 399 void AliPMDClusterFinder::AddRecPoint(Int_t idet,Int_t ismn,Float_t *clusdata)
 400 {
 401   // Add Reconstructed points
 402   //
 403   TClonesArray &lrecpoints = *fRecpoints;
 404   AliPMDrecpoint1 *newrecpoint;
 405   newrecpoint = new AliPMDrecpoint1(idet, ismn, clusdata);
 406   new(lrecpoints[fNpoint++]) AliPMDrecpoint1(newrecpoint);
 407   delete newrecpoint;
 408 }
 409 // ------------------------------------------------------------------------- //
 410 void AliPMDClusterFinder::ResetCellADC()
 411 {
 412   // Reset the individual cell ADC value to zero
 413   //
 414   for(Int_t irow = 0; irow < fgkRow; irow++)
 415     {
 416       for(Int_t icol = 0; icol < fgkCol; icol++)
 417         {
 418           fCellADC[irow][icol] = 0.;
 419         }
 420     }
 421 }
 422 // ------------------------------------------------------------------------- //
 423
 424 void AliPMDClusterFinder::ResetRecpoint()
 425 {
 426   // Clear the list of reconstructed points
 427   fNpoint = 0;
 428   if (fRecpoints) fRecpoints->Clear();
 429 }
 430 // ------------------------------------------------------------------------- //
 431 void AliPMDClusterFinder::Load()
 432 {
 433   // Load all the *.root files
 434   //
 435   fPMDLoader->LoadDigits("READ");
 436   fPMDLoader->LoadRecPoints("recreate");
 437 }
 438 // ------------------------------------------------------------------------- //
 439 void AliPMDClusterFinder::LoadClusters()
 440 {
 441   // Load all the *.root files
 442   //
 443   fPMDLoader->LoadRecPoints("recreate");
 444 }
 445 // ------------------------------------------------------------------------- //
 446 void AliPMDClusterFinder::UnLoad()
 447 {
 448   // Unload all the *.root files
 449   //
 450   fPMDLoader->UnloadDigits();
 451   fPMDLoader->UnloadRecPoints();
 452 }
 453 // ------------------------------------------------------------------------- //
 454 void AliPMDClusterFinder::UnLoadClusters()
 455 {
 456   // Unload all the *.root files
 457   //
 458   fPMDLoader->UnloadRecPoints();
 459 }
 460 // ------------------------------------------------------------------------- //