]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - PMD/AliPMDClusterFinder.cxx
Added files for the reference version based on an old version of Anders'
[u/mrichter/AliRoot.git] / PMD / AliPMDClusterFinder.cxx
index 5c9f3709da6a4fb69940364f9de820227cd1b4ac..beac2228f7d303e056b492f8f8b51d2ea39e8912 100644 (file)
 //-----------------------------------------------------//
 
 #include <Riostream.h>
-#include <TMath.h>
-#include <TBRIK.h>
-#include <TNode.h>
 #include <TTree.h>
-#include <TGeometry.h>
 #include <TObjArray.h>
 #include <TClonesArray.h>
-#include <TFile.h>
-#include <TNtuple.h>
-#include <TParticle.h>
 
-#include "AliRun.h"
-#include "AliPMD.h"
-#include "AliDetector.h"
+#include "AliLog.h"
 #include "AliRunLoader.h"
 #include "AliLoader.h"
-#include "AliHeader.h"
 #include "AliRawReader.h"
 
 #include "AliPMDdigit.h"
 #include "AliPMDrecpoint1.h"
 #include "AliPMDRawStream.h"
 
+
+
 ClassImp(AliPMDClusterFinder)
 
+AliPMDClusterFinder::AliPMDClusterFinder():
+  fRunLoader(0),
+  fPMDLoader(0),
+  fTreeD(0),
+  fTreeR(0),
+  fDigits(new TClonesArray("AliPMDdigit", 1000)),
+  fRecpoints(new TClonesArray("AliPMDrecpoint1", 1000)),
+  fNpoint(0),
+  fEcut(0.)
+{
+//
+// Constructor
+//
+}
+// ------------------------------------------------------------------------- //
 AliPMDClusterFinder::AliPMDClusterFinder(AliRunLoader* runLoader):
   fRunLoader(runLoader),
   fPMDLoader(runLoader->GetLoader("PMDLoader")),
@@ -59,7 +66,6 @@ AliPMDClusterFinder::AliPMDClusterFinder(AliRunLoader* runLoader):
   fDigits(new TClonesArray("AliPMDdigit", 1000)),
   fRecpoints(new TClonesArray("AliPMDrecpoint1", 1000)),
   fNpoint(0),
-  fDebug(0),
   fEcut(0.)
 {
 //
@@ -99,21 +105,24 @@ void AliPMDClusterFinder::Digits2RecPoints(Int_t ievt)
 
   TObjArray *pmdcont = new TObjArray();
   AliPMDClustering *pmdclust = new AliPMDClustering();
-  pmdclust->SetDebug(fDebug);
+
   pmdclust->SetEdepCut(fEcut);
 
   fRunLoader->GetEvent(ievt);
-  //cout << " ***** Beginning::Digits2RecPoints *****" << endl;
+
+
   fTreeD = fPMDLoader->TreeD();
   if (fTreeD == 0x0)
     {
-      cout << " Can not get TreeD" << endl;
+      AliFatal("AliPMDClusterFinder: Can not get TreeD");
+
     }
   AliPMDdigit  *pmddigit;
   TBranch *branch = fTreeD->GetBranch("PMDDigit");
   branch->SetAddress(&fDigits);
 
   ResetRecpoint();
+
   fTreeR = fPMDLoader->TreeR();
   if (fTreeR == 0x0)
     {
@@ -125,7 +134,7 @@ void AliPMDClusterFinder::Digits2RecPoints(Int_t ievt)
   fTreeR->Branch("PMDRecpoint", &fRecpoints, bufsize); 
 
   Int_t nmodules = (Int_t) fTreeD->GetEntries();
-  
+
   for (Int_t imodule = 0; imodule < nmodules; imodule++)
     {
       ResetCellADC();
@@ -149,7 +158,9 @@ void AliPMDClusterFinder::Digits2RecPoints(Int_t ievt)
       pmdclust->DoClust(idet,ismn,fCellADC,pmdcont);
       
       Int_t nentries1 = pmdcont->GetEntries();
-//      cout << " nentries1 = " << nentries1 << endl;
+
+      AliDebug(1,Form("Total number of clusters/module = %d",nentries1));
+
       for (Int_t ient1 = 0; ient1 < nentries1; ient1++)
        {
          AliPMDcluster *pmdcl = (AliPMDcluster*)pmdcont->UncheckedAt(ient1);
@@ -160,7 +171,7 @@ void AliPMDClusterFinder::Digits2RecPoints(Int_t ievt)
          clusdata[2] = pmdcl->GetClusADC();
          clusdata[3] = pmdcl->GetClusCells();
          clusdata[4] = pmdcl->GetClusRadius();
-         
+
          AddRecPoint(idet,ismn,clusdata);
        }
       pmdcont->Clear();
@@ -178,107 +189,366 @@ void AliPMDClusterFinder::Digits2RecPoints(Int_t ievt)
   delete pmdclust;
   delete pmdcont;
     
-  //  cout << " ***** End::Digits2RecPoints *****" << endl;
 }
 // ------------------------------------------------------------------------- //
 
-void AliPMDClusterFinder::Digits2RecPoints(Int_t ievt, AliRawReader *rawReader)
+void AliPMDClusterFinder::Digits2RecPoints(AliRawReader *rawReader,
+                                          TTree *clustersTree)
 {
-  // Converts digits to recpoints after running clustering
-  // algorithm on CPV plane and PREshower plane
+  // Converts RAW data to recpoints after running clustering
+  // algorithm on CPV and PREshower plane
   //
 
   Float_t  clusdata[5];
 
+  TObjArray *pmdcont = new TObjArray();
+  AliPMDClustering *pmdclust = new AliPMDClustering();
+
+  pmdclust->SetEdepCut(fEcut);
+
+  ResetRecpoint();
 
-  rawReader->Reset();
-  AliPMDRawStream pmdinput(rawReader);
+  Int_t bufsize = 16000;
+  clustersTree->Branch("PMDRecpoint", &fRecpoints, bufsize); 
+
+  const Int_t kDDL = 6;
+  const Int_t kRow = 48;
+  const Int_t kCol = 96;
+
+  Int_t idet = 0;
+  Int_t iSMN = 0;
+  
+  for (Int_t indexDDL = 0; indexDDL < kDDL; indexDDL++)
+    {
+      if (indexDDL < 4)
+       {
+         iSMN = 6;
+       }
+      else if (indexDDL >= 4)
+       {
+         iSMN = 12;
+       }
+      Int_t ***precpvADC;
+      precpvADC = new int **[iSMN];
+      for (Int_t i=0; i<iSMN; i++) precpvADC[i] = new int *[kRow];
+      for (Int_t i=0; i<iSMN;i++)
+       {
+         for (Int_t j=0; j<kRow; j++) precpvADC[i][j] = new int [kCol];
+       }
+      for (Int_t i = 0; i < iSMN; i++)
+       {
+         for (Int_t j = 0; j < kRow; j++)
+           {
+             for (Int_t k = 0; k < kCol; k++)
+               {
+                 precpvADC[i][j][k] = 0;
+               }
+           }
+       }
+      ResetCellADC();
+      rawReader->Reset();
+      AliPMDRawStream pmdinput(rawReader);
+      rawReader->Select(12, indexDDL, indexDDL);
+      while(pmdinput.Next())
+       {
+         Int_t det = pmdinput.GetDetector();
+         Int_t smn = pmdinput.GetSMN();
+         //Int_t mcm = pmdinput.GetMCM();
+         //Int_t chno = pmdinput.GetChannel();
+         Int_t row = pmdinput.GetRow();
+         Int_t col = pmdinput.GetColumn();
+         Int_t sig = pmdinput.GetSignal();
+         
+         Int_t indexsmn = 0;
+
+         if (indexDDL < 4)
+           {
+             if (det != 0)
+               AliError(Form("*DDL %d and Detector NUMBER %d NOT MATCHING *",
+                             indexDDL, det));
+             indexsmn = smn - indexDDL * 6;
+           }
+         else if (indexDDL == 4)
+           {
+             if (det != 1)
+               AliError(Form("*DDL %d and Detector NUMBER %d NOT MATCHING *",
+                             indexDDL, det));
+             if (smn < 6)
+               {
+                 indexsmn = smn;
+               }
+             else if (smn >= 12 && smn < 18)
+               {
+                 indexsmn = smn - 6;
+               }
+           }
+         else if (indexDDL == 5)
+           {
+             if (det != 1)
+               AliError(Form("*DDL %d and Detector NUMBER %d NOT MATCHING *",
+                             indexDDL, det));
+             if (smn >= 6 && smn < 12)
+               {
+                 indexsmn = smn - 6;
+               }
+             else if (smn >= 18 && smn < 24)
+               {
+                 indexsmn = smn - 12;
+               }
+           }         
+         precpvADC[indexsmn][row][col] = sig;
+       } // while loop
+
+      Int_t ismn = 0;
+      for (Int_t indexsmn = 0; indexsmn < iSMN; indexsmn++)
+       {
+         ResetCellADC();
+         for (Int_t irow = 0; irow < kRow; irow++)
+           {
+             for (Int_t icol = 0; icol < kCol; icol++)
+               {
+                 fCellADC[irow][icol] = 
+                   (Double_t) precpvADC[indexsmn][irow][icol];
+               } // row
+           }     // col
+         if (indexDDL < 4)
+           {
+             ismn = indexsmn + indexDDL * 6;
+             idet = 0;
+           }
+         else if (indexDDL == 4)
+           {
+             if (indexsmn < 6)
+               {
+                 ismn = indexsmn;
+               }
+             else if (indexsmn >= 6 && indexsmn < 12)
+               {
+                 ismn = indexsmn + 6;
+               }
+             idet = 1;
+           }
+         else if (indexDDL == 5)
+           {
+             if (indexsmn < 6)
+               {
+                 ismn = indexsmn + 6;
+               }
+             else if (indexsmn >= 6 && indexsmn < 12)
+               {
+                 ismn = indexsmn + 12;
+               }
+             idet = 1;
+           }
+
+
+         pmdclust->DoClust(idet,ismn,fCellADC,pmdcont);
+         Int_t nentries1 = pmdcont->GetEntries();
+
+         AliDebug(1,Form("Total number of clusters/module = %d",nentries1));
+
+         for (Int_t ient1 = 0; ient1 < nentries1; ient1++)
+           {
+             AliPMDcluster *pmdcl = 
+               (AliPMDcluster*)pmdcont->UncheckedAt(ient1);
+             idet        = pmdcl->GetDetector();
+             ismn        = pmdcl->GetSMN();
+             clusdata[0] = pmdcl->GetClusX();
+             clusdata[1] = pmdcl->GetClusY();
+             clusdata[2] = pmdcl->GetClusADC();
+             clusdata[3] = pmdcl->GetClusCells();
+             clusdata[4] = pmdcl->GetClusRadius();
+
+             AddRecPoint(idet,ismn,clusdata);
+           }
+         pmdcont->Clear();
+         
+         //fTreeR->Fill();
+         clustersTree->Fill();
+         ResetRecpoint();
+
+
+       } // smn
+
+      for (Int_t i=0; i<iSMN; i++)
+       {
+         for (Int_t j=0; j<kRow; j++) delete [] precpvADC[i][j];
+       }
+      for (Int_t i=0; i<iSMN; i++) delete [] precpvADC[i];
+      delete precpvADC;
+    } // DDL Loop
+  
+  ResetCellADC();
+  
+  //   delete the pointers
+  delete pmdclust;
+  delete pmdcont;
+
+}
+// ------------------------------------------------------------------------- //
+
+void AliPMDClusterFinder::Digits2RecPoints(Int_t ievt, AliRawReader *rawReader)
+{
+  // Converts RAW data to recpoints after running clustering
+  // algorithm on CPV and PREshower plane
+  //
+
+  Float_t  clusdata[5];
 
   TObjArray *pmdcont = new TObjArray();
   AliPMDClustering *pmdclust = new AliPMDClustering();
-  pmdclust->SetDebug(fDebug);
+
   pmdclust->SetEdepCut(fEcut);
 
   fRunLoader->GetEvent(ievt);
 
   ResetRecpoint();
+
   fTreeR = fPMDLoader->TreeR();
   if (fTreeR == 0x0)
     {
       fPMDLoader->MakeTree("R");
       fTreeR = fPMDLoader->TreeR();
     }
-
   Int_t bufsize = 16000;
   fTreeR->Branch("PMDRecpoint", &fRecpoints, bufsize); 
-  const Int_t kDet = 2;
-  const Int_t kSMN = 24;
+
+  const Int_t kDDL = 6;
   const Int_t kRow = 48;
   const Int_t kCol = 96;
-  Int_t preADC[kSMN][kRow][kCol];
-  Int_t cpvADC[kSMN][kRow][kCol];
 
-  for (Int_t i = 0; i < kSMN; i++)
+  Int_t idet = 0;
+  Int_t iSMN = 0;
+  
+  for (Int_t indexDDL = 0; indexDDL < kDDL; indexDDL++)
     {
-      for (Int_t j = 0; j < kRow; j++)
+      if (indexDDL < 4)
        {
-         for (Int_t k = 0; k < kCol; k++)
-           {
-             preADC[i][j][k] = 0;
-             cpvADC[i][j][k] = 0;
-           }
+         iSMN = 6;
        }
-    }
-  ResetCellADC();
-
-  while(pmdinput.Next())
-    {
-      Int_t det = pmdinput.GetDetector();
-      Int_t smn = pmdinput.GetSMN();
-      //Int_t mcm = pmdinput.GetMCM();
-      //Int_t chno = pmdinput.GetChannel();
-      Int_t row = pmdinput.GetRow();
-      Int_t col = pmdinput.GetColumn();
-      Int_t sig = pmdinput.GetSignal();
-
-      //cout << " det = " << det << " smn = " << smn 
-      //   << " row = " << row << " col = " << col
-      //   << " sig = " << sig << endl;
-
-      if (det == 0)
+      else if (indexDDL >= 4)
        {
-         preADC[smn][row][col] = sig;
+         iSMN = 12;
        }
-      else if (det == 1)
+      Int_t ***precpvADC;
+      precpvADC = new int **[iSMN];
+      for (Int_t i=0; i<iSMN; i++) precpvADC[i] = new int *[kRow];
+      for (Int_t i=0; i<iSMN;i++)
        {
-         cpvADC[smn][row][col] = sig;
+         for (Int_t j=0; j<kRow; j++) precpvADC[i][j] = new int [kCol];
        }
+      for (Int_t i = 0; i < iSMN; i++)
+       {
+         for (Int_t j = 0; j < kRow; j++)
+           {
+             for (Int_t k = 0; k < kCol; k++)
+               {
+                 precpvADC[i][j][k] = 0;
+               }
+           }
+       }
+      ResetCellADC();
+      rawReader->Reset();
+      AliPMDRawStream pmdinput(rawReader);
+      rawReader->Select(12, indexDDL, indexDDL);
+      while(pmdinput.Next())
+       {
+         Int_t det = pmdinput.GetDetector();
+         Int_t smn = pmdinput.GetSMN();
+         //Int_t mcm = pmdinput.GetMCM();
+         //Int_t chno = pmdinput.GetChannel();
+         Int_t row = pmdinput.GetRow();
+         Int_t col = pmdinput.GetColumn();
+         Int_t sig = pmdinput.GetSignal();
+         
+         Int_t indexsmn = 0;
 
-    } // while loop
-
-  for (Int_t idet = 0; idet < kDet; idet++)
-    {
-      for (Int_t ismn = 0; ismn < kSMN; ismn++)
+         if (indexDDL < 4)
+           {
+             if (det != 0)
+               AliError(Form("*DDL %d and Detector NUMBER %d NOT MATCHING *",
+                             indexDDL, det));
+             indexsmn = smn - indexDDL * 6;
+           }
+         else if (indexDDL == 4)
+           {
+             if (det != 1)
+               AliError(Form("*DDL %d and Detector NUMBER %d NOT MATCHING *",
+                             indexDDL, det));
+             if (smn < 6)
+               {
+                 indexsmn = smn;
+               }
+             else if (smn >= 12 && smn < 18)
+               {
+                 indexsmn = smn - 6;
+               }
+           }
+         else if (indexDDL == 5)
+           {
+             if (det != 1)
+               AliError(Form("*DDL %d and Detector NUMBER %d NOT MATCHING *",
+                             indexDDL, det));
+             if (smn >= 6 && smn < 12)
+               {
+                 indexsmn = smn - 6;
+               }
+             else if (smn >= 18 && smn < 24)
+               {
+                 indexsmn = smn - 12;
+               }
+           }         
+         precpvADC[indexsmn][row][col] = sig;
+       } // while loop
+
+      Int_t ismn = 0;
+      for (Int_t indexsmn = 0; indexsmn < iSMN; indexsmn++)
        {
+         ResetCellADC();
          for (Int_t irow = 0; irow < kRow; irow++)
            {
              for (Int_t icol = 0; icol < kCol; icol++)
                {
-                 if (idet == 0)
-                   {
-                     fCellADC[irow][icol] = 
-                       (Double_t) preADC[ismn][irow][icol];
-                   }
-                 else if (idet == 1)
-                   {
-                     fCellADC[irow][icol] = 
-                       (Double_t) cpvADC[ismn][irow][icol];
-                   }
+                 fCellADC[irow][icol] = 
+                   (Double_t) precpvADC[indexsmn][irow][icol];
                } // row
            }     // col
+         if (indexDDL < 4)
+           {
+             ismn = indexsmn + indexDDL * 6;
+             idet = 0;
+           }
+         else if (indexDDL == 4)
+           {
+             if (indexsmn < 6)
+               {
+                 ismn = indexsmn;
+               }
+             else if (indexsmn >= 6 && indexsmn < 12)
+               {
+                 ismn = indexsmn + 6;
+               }
+             idet = 1;
+           }
+         else if (indexDDL == 5)
+           {
+             if (indexsmn < 6)
+               {
+                 ismn = indexsmn + 6;
+               }
+             else if (indexsmn >= 6 && indexsmn < 12)
+               {
+                 ismn = indexsmn + 12;
+               }
+             idet = 1;
+           }
+
+
          pmdclust->DoClust(idet,ismn,fCellADC,pmdcont);
          Int_t nentries1 = pmdcont->GetEntries();
-         //      cout << " nentries1 = " << nentries1 << endl;
+
+         AliDebug(1,Form("Total number of clusters/module = %d",nentries1));
+
          for (Int_t ient1 = 0; ient1 < nentries1; ient1++)
            {
              AliPMDcluster *pmdcl = 
@@ -290,18 +560,27 @@ void AliPMDClusterFinder::Digits2RecPoints(Int_t ievt, AliRawReader *rawReader)
              clusdata[2] = pmdcl->GetClusADC();
              clusdata[3] = pmdcl->GetClusCells();
              clusdata[4] = pmdcl->GetClusRadius();
-             
+
              AddRecPoint(idet,ismn,clusdata);
            }
          pmdcont->Clear();
          
          fTreeR->Fill();
          ResetRecpoint();
-         
+
+
        } // smn
-    }     // det (pre, cpv)
 
+      for (Int_t i=0; i<iSMN; i++)
+       {
+         for (Int_t j=0; j<kRow; j++) delete [] precpvADC[i][j];
+       }
+      for (Int_t i=0; i<iSMN; i++) delete [] precpvADC[i];
+      delete precpvADC;
+    } // DDL Loop
+  
   ResetCellADC();
+  
   fPMDLoader = fRunLoader->GetLoader("PMDLoader");  
   fPMDLoader->WriteRecPoints("OVERWRITE");
 
@@ -309,7 +588,6 @@ void AliPMDClusterFinder::Digits2RecPoints(Int_t ievt, AliRawReader *rawReader)
   delete pmdclust;
   delete pmdcont;
 
-  //  cout << " ***** End::Digits2RecPoints :: Raw *****" << endl;
 }
 // ------------------------------------------------------------------------- //
 void AliPMDClusterFinder::SetCellEdepCut(Float_t ecut)
@@ -317,11 +595,6 @@ void AliPMDClusterFinder::SetCellEdepCut(Float_t ecut)
   fEcut = ecut;
 }
 // ------------------------------------------------------------------------- //
-void AliPMDClusterFinder::SetDebug(Int_t idebug)
-{
-  fDebug = idebug;
-}
-// ------------------------------------------------------------------------- //
 void AliPMDClusterFinder::AddRecPoint(Int_t idet,Int_t ismn,Float_t *clusdata)
 {
   // Add Reconstructed points
@@ -362,6 +635,13 @@ void AliPMDClusterFinder::Load()
   fPMDLoader->LoadRecPoints("recreate");
 }
 // ------------------------------------------------------------------------- //
+void AliPMDClusterFinder::LoadClusters()
+{
+  // Load all the *.root files
+  //
+  fPMDLoader->LoadRecPoints("recreate");
+}
+// ------------------------------------------------------------------------- //
 void AliPMDClusterFinder::UnLoad()
 {
   // Unload all the *.root files
@@ -370,3 +650,10 @@ void AliPMDClusterFinder::UnLoad()
   fPMDLoader->UnloadRecPoints();
 }
 // ------------------------------------------------------------------------- //
+void AliPMDClusterFinder::UnLoadClusters()
+{
+  // Unload all the *.root files
+  //
+  fPMDLoader->UnloadRecPoints();
+}
+// ------------------------------------------------------------------------- //