Adding possibility to change the tracking window in ctgTheta
[u/mrichter/AliRoot.git] / TPC / AliTPCclustererMI.cxx
index 82cc2d80cf32e9c91cb48cb91bf7566af3f3651d..d546627eff1dc3d2fcd6eecbc43a6461bcfefdcb 100644 (file)
  * provided "as is" without express or implied warranty.                  *
  **************************************************************************/
 
-/*
-$Log$
-Revision 1.4  2003/03/05 11:16:15  kowal2
-Logs added
-
-*/
-
-
-
+/* $Id$ */
 
 //-------------------------------------------------------
 //          Implementation of the TPC clusterer
@@ -29,16 +21,21 @@ Logs added
 //   Origin: Marian Ivanov 
 //-------------------------------------------------------
 
+#include "AliTPCReconstructor.h"
 #include "AliTPCclustererMI.h"
 #include "AliTPCclusterMI.h"
 #include <TObjArray.h>
 #include <TFile.h>
 #include "AliTPCClustersArray.h"
 #include "AliTPCClustersRow.h"
+#include "AliTPCRawStream.h"
 #include "AliDigits.h"
 #include "AliSimDigits.h"
 #include "AliTPCParam.h"
-#include "AliTPCBuffer160.h"
+#include "AliRawReader.h"
+#include "AliTPCRawStream.h"
+#include "AliRunLoader.h"
+#include "AliLoader.h"
 #include "Riostream.h"
 #include <TTree.h>
 
@@ -102,7 +99,7 @@ Float_t  AliTPCclustererMI::GetSigmaZ2(Int_t iz){
   return res;
 }
 
-void AliTPCclustererMI::MakeCluster(Int_t k,Int_t max,Int_t *bins, UInt_t m,
+void AliTPCclustererMI::MakeCluster(Int_t k,Int_t max,Int_t *bins, UInt_t /*m*/,
 AliTPCclusterMI &c) 
 {
   Int_t i0=k/max;  //central pad
@@ -448,21 +445,20 @@ void AliTPCclustererMI::AddCluster(AliTPCclusterMI &c){
 
 
 //_____________________________________________________________________________
-void AliTPCclustererMI::Digits2Clusters(const AliTPCParam *par, Int_t eventn)
+void AliTPCclustererMI::Digits2Clusters()
 {
   //-----------------------------------------------------------------
   // This is a simple cluster finder.
   //-----------------------------------------------------------------
-  TDirectory *savedir=gDirectory; 
 
if (fInput==0)
-    cerr<<"AliTPC::Digits2Clusters(): input tree not initialised !\n";
 if (!fInput) 
+    Error("Digits2Clusters", "input tree not initialised");
     return;
   }
  
-  if (fOutput==0) {
-     cerr<<"AliTPC::Digits2Clusters(): output tree not initialised !\n";
-     return;
+  if (!fOutput) {
+    Error("Digits2Clusters", "output tree not initialised");
+    return;
   }
 
   AliSimDigits digarr, *dummy=&digarr;
@@ -470,17 +466,14 @@ void AliTPCclustererMI::Digits2Clusters(const AliTPCParam *par, Int_t eventn)
   fInput->GetBranch("Segment")->SetAddress(&dummy);
   Stat_t nentries = fInput->GetEntries();
   
-  fMaxTime=par->GetMaxTBin()+6; // add 3 virtual time bins before and 3   after
+  fMaxTime=fParam->GetMaxTBin()+6; // add 3 virtual time bins before and 3   after
     
-  fParam = par;
-  ((AliTPCParam*)par)->Write(par->GetTitle());
-
   Int_t nclusters  = 0;
   
   for (Int_t n=0; n<nentries; n++) {
     fInput->GetEvent(n);
     Int_t row;
-    if (!par->AdjustSectorRow(digarr.GetID(),fSector,row)) {
+    if (!fParam->AdjustSectorRow(digarr.GetID(),fSector,row)) {
       cerr<<"AliTPC warning: invalid segment ID ! "<<digarr.GetID()<<endl;
       continue;
     }
@@ -492,21 +485,21 @@ void AliTPCclustererMI::Digits2Clusters(const AliTPCParam *par, Int_t eventn)
 
     clrow->SetID(digarr.GetID());
     fOutput->GetBranch("Segment")->SetAddress(&clrow);
-    fRx=par->GetPadRowRadii(fSector,row);
+    fRx=fParam->GetPadRowRadii(fSector,row);
     
     
-    const Int_t kNIS=par->GetNInnerSector(), kNOS=par->GetNOuterSector();
+    const Int_t kNIS=fParam->GetNInnerSector(), kNOS=fParam->GetNOuterSector();
     fZWidth = fParam->GetZWidth();
     if (fSector < kNIS) {
-      fMaxPad = par->GetNPadsLow(row);
+      fMaxPad = fParam->GetNPadsLow(row);
       fSign = (fSector < kNIS/2) ? 1 : -1;
-      fPadLength = par->GetPadPitchLength(fSector,row);
-      fPadWidth = par->GetPadPitchWidth();
+      fPadLength = fParam->GetPadPitchLength(fSector,row);
+      fPadWidth = fParam->GetPadPitchWidth();
     } else {
-      fMaxPad = par->GetNPadsUp(row);
+      fMaxPad = fParam->GetNPadsUp(row);
       fSign = ((fSector-kNIS) < kNOS/2) ? 1 : -1;
-      fPadLength = par->GetPadPitchLength(fSector,row);
-      fPadWidth  = par->GetPadPitchWidth();
+      fPadLength = fParam->GetPadPitchLength(fSector,row);
+      fPadWidth  = fParam->GetPadPitchWidth();
     }
     
     
@@ -518,7 +511,7 @@ void AliTPCclustererMI::Digits2Clusters(const AliTPCParam *par, Int_t eventn)
     if (digarr.First()) //MI change
       do {
        Short_t dig=digarr.CurrentDigit();
-       if (dig<=par->GetZeroSup()) continue;
+       if (dig<=fParam->GetZeroSup()) continue;
        Int_t j=digarr.CurrentRow()+3, i=digarr.CurrentColumn()+3;
        fBins[i*fMaxTime+j]=dig;
       } while (digarr.Next());
@@ -532,142 +525,154 @@ void AliTPCclustererMI::Digits2Clusters(const AliTPCParam *par, Int_t eventn)
     delete[] fBins;      
     delete[] fResBins;  
   }  
-  cerr<<"Number of found clusters : "<<nclusters<<"                        \n";
-  
-  fOutput->Write();
-  savedir->cd();
+
+  Info("Digits2Clusters", "Number of found clusters : %d", nclusters);
 }
 
-void AliTPCclustererMI::Digits2Clusters(const char* fileName)
+void AliTPCclustererMI::Digits2Clusters(AliRawReader* rawReader)
 {
 //-----------------------------------------------------------------
 // This is a cluster finder for raw data.
-//
-// It assumes, that the input file contains uncompressed data 
-// in TPC Altro format without mini headers.
-// It also assumes, that the data is sorted in a such a way that
-// data from different pad rows is not mixed. The order of pads
-// in a pad row can be random.
 //-----------------------------------------------------------------
 
-  TDirectory* savedir = gDirectory; 
-
-  AliTPCBuffer160 buffer(fileName, 0);  // buffer for reading raw data
-  fRowDig = NULL;
   if (!fOutput) {
-    Error("Digits2Clusters", "output tree not initialised !\n");
+    Error("Digits2Clusters", "output tree not initialised");
     return;
   }
-  fOutput->GetCurrentFile()->cd();
-  const_cast<AliTPCParam*>(fParam)->Write();
+
+  rawReader->Reset();
+  AliTPCRawStream input(rawReader);
+
+  fRowDig = NULL;
 
   Int_t nclusters  = 0;
   
   fMaxTime = fParam->GetMaxTBin() + 6; // add 3 virtual time bins before and 3 after
   const Int_t kNIS = fParam->GetNInnerSector();
   const Int_t kNOS = fParam->GetNOuterSector();
+  const Int_t kNS = kNIS + kNOS;
   fZWidth = fParam->GetZWidth();
   Int_t zeroSup = fParam->GetZeroSup();
-  Int_t offset = 1;
 
-  Int_t nWords, iPad, iRow;
-  Int_t prevSector = -1;
-  Int_t prevRow = -1;
   fBins = NULL;
-  while (buffer.ReadTrailerBackward(nWords, iPad, iRow, fSector) == 0) {
-    Int_t nSkip = (4 - (nWords % 4)) % 4;
+  Int_t** splitRows = new Int_t* [kNS*2];
+  Int_t** splitRowsRes = new Int_t* [kNS*2];
+  for (Int_t iSector = 0; iSector < kNS*2; iSector++)
+    splitRows[iSector] = NULL;
+  Int_t iSplitRow = -1;
+
+  Bool_t next = kTRUE;
+  while (next) {
+    next = input.Next();
 
     // when the sector or row number has changed ...
-    if ((fSector != prevSector) || (iRow != prevRow)) {
+    if (input.IsNewRow() || !next) {
 
       // ... find clusters in the previous pad row, and ...
       if (fBins) {
-       FindClusters();
-
-       fOutput->Fill();
-       delete fRowCl;    
-       nclusters += fNcluster;    
-       delete[] fBins;      
-       delete[] fResBins;  
+       if ((iSplitRow < 0) || splitRows[fSector + kNS*iSplitRow]) {
+         fRowCl = new AliTPCClustersRow;
+         fRowCl->SetClass("AliTPCclusterMI");
+         fRowCl->SetArray(1);
+         fRowCl->SetID(fParam->GetIndex(fSector, input.GetPrevRow()));
+         fOutput->GetBranch("Segment")->SetAddress(&fRowCl);
+
+         FindClusters();
+
+         fOutput->Fill();
+         delete fRowCl;    
+         nclusters += fNcluster;    
+         delete[] fBins;
+         delete[] fResBins;
+         if (iSplitRow >= 0) splitRows[fSector + kNS*iSplitRow] = NULL;
+
+       } else if (iSplitRow >= 0) {
+         splitRows[fSector + kNS*iSplitRow] = fBins;
+         splitRowsRes[fSector + kNS*iSplitRow] = fResBins;
+       }
       }
 
-      // ... prepare for the next pad row
-      fRowCl = new AliTPCClustersRow;
-      fRowCl->SetClass("AliTPCclusterMI");
-      fRowCl->SetArray(1);
+      if (!next) break;
 
-      fRowCl->SetID(fParam->GetIndex(fSector, iRow));
-      fOutput->GetBranch("Segment")->SetAddress(&fRowCl);
+      // ... prepare for the next pad row
+      fSector = input.GetSector();
+      Int_t iRow = input.GetRow();
       fRx = fParam->GetPadRowRadii(fSector, iRow);
     
+      iSplitRow = -1;
       if (fSector < kNIS) {
        fMaxPad = fParam->GetNPadsLow(iRow);
        fSign = (fSector < kNIS/2) ? 1 : -1;
+       if (iRow == 30) iSplitRow = 0;
       } else {
        fMaxPad = fParam->GetNPadsUp(iRow);
        fSign = ((fSector-kNIS) < kNOS/2) ? 1 : -1;
+       if (iRow == 27) iSplitRow = 0;
+       else if (iRow == 76) iSplitRow = 1;
       }
       fPadLength = fParam->GetPadPitchLength(fSector, iRow);
       fPadWidth  = fParam->GetPadPitchWidth();
     
       fMaxBin = fMaxTime*(fMaxPad+6);  // add 3 virtual pads  before and 3 after
-      fBins    = new Int_t[fMaxBin];
-      fResBins = new Int_t[fMaxBin];  //fBins with residuals after 1 finder loop 
-      memset(fBins, 0, sizeof(Int_t)*fMaxBin);
-
-      prevSector = fSector;
-      prevRow = iRow;
+      if ((iSplitRow < 0) || !splitRows[fSector + kNS*iSplitRow]) {
+       fBins    = new Int_t[fMaxBin];
+       fResBins = new Int_t[fMaxBin];  //fBins with residuals after 1 finder loop 
+       memset(fBins, 0, sizeof(Int_t)*fMaxBin);
+      } else {
+       fBins    = splitRows[fSector + kNS*iSplitRow];
+       fResBins = splitRowsRes[fSector + kNS*iSplitRow];
+      }
     }
 
     // fill fBins with digits data
-    for (Int_t iWord = 0; iWord < nSkip; iWord++) {
-      if (buffer.GetNextBackWord() != 0x2AA) {
-       Error("Digits2Clusters", "could not read dummy data\n");
-       break;
-      }
-    }
-    while (nWords > 0) {
-      Int_t bunchLength = buffer.GetNextBackWord() - 2;
-      if (bunchLength < 0) {
-       Error("Digits2Clusters", "could not read bunch length\n");
-       break;
-      }
-      Int_t timeBin = buffer.GetNextBackWord();
-      if (timeBin < 0) {
-       Error("Digits2Clusters", "could not read time bin\n");
-       break;
-      }
-      for (Int_t iTimeBin = timeBin; iTimeBin > timeBin-bunchLength; iTimeBin--) {
-       Int_t dig = buffer.GetNextBackWord();
-       if (dig < 0) {
-         Error("Digits2Clusters", "could not read sample amplitude\n");
-         break;
+    if (input.GetSignal() <= zeroSup) continue;
+    Int_t i = input.GetPad() + 3;
+    Int_t j = input.GetTime() + 3;
+    fBins[i*fMaxTime+j] = input.GetSignal();
+  }
+
+  // find clusters in split rows that were skipped until now.
+  // this can happen if the rows were not splitted 
+  for (fSector = 0; fSector < kNS; fSector++)
+    for (Int_t iSplit = 0; iSplit < 2; iSplit++)
+      if (splitRows[fSector + kNS*iSplit]) {
+
+       Int_t iRow = -1;
+       if (fSector < kNIS) {
+         iRow = 30;
+         fMaxPad = fParam->GetNPadsLow(iRow);
+         fSign = (fSector < kNIS/2) ? 1 : -1;
+       } else {
+         if (iSplit == 0) iRow = 27; else iRow = 76;
+         fMaxPad = fParam->GetNPadsUp(iRow);
+         fSign = ((fSector-kNIS) < kNOS/2) ? 1 : -1;
        }
-       dig += offset;
-       if (dig <= zeroSup) continue;
-       Int_t i = iPad + 3;
-       Int_t j = iTimeBin + 3;
-       fBins[i*fMaxTime+j] = dig;
-      }
-      nWords -= 2+bunchLength;
-    }
-  }  
+       fRx = fParam->GetPadRowRadii(fSector, iRow);
+       fPadLength = fParam->GetPadPitchLength(fSector, iRow);
+       fPadWidth  = fParam->GetPadPitchWidth();
 
-  // find clusters in the last pad row
-  if (fBins) {
-    FindClusters();
+       fMaxBin = fMaxTime*(fMaxPad+6);  // add 3 virtual pads  before and 3 after
+       fBins    = splitRows[fSector + kNS*iSplit];
+       fResBins = splitRowsRes[fSector + kNS*iSplit];
 
-    fOutput->Fill();
-    delete fRowCl;    
-    nclusters += fNcluster;    
-    delete[] fBins;      
-    delete[] fResBins;  
-  }
+       fRowCl = new AliTPCClustersRow;
+       fRowCl->SetClass("AliTPCclusterMI");
+       fRowCl->SetArray(1);
+       fRowCl->SetID(fParam->GetIndex(fSector, iRow));
+       fOutput->GetBranch("Segment")->SetAddress(&fRowCl);
 
-  Info("Digits2Clusters", "Number of found clusters : %d\n", nclusters);
+       FindClusters();
 
-  fOutput->Write();
-  savedir->cd();
+       fOutput->Fill();
+       delete fRowCl;    
+       nclusters += fNcluster;    
+       delete[] fBins;
+       delete[] fResBins;
+      }
+
+  delete[] splitRows;
+  delete[] splitRowsRes;
+  Info("Digits2Clusters", "Number of found clusters : %d\n", nclusters);
 }
 
 void AliTPCclustererMI::FindClusters()
@@ -703,7 +708,7 @@ void AliTPCclustererMI::FindClusters()
   //first loop - for "gold cluster" 
   fLoop=1;
   Int_t *b=&fBins[-1]+2*fMaxTime;
-  Int_t crtime = Int_t((fParam->GetZLength()-1.05*fRx)/fZWidth-5);
+  Int_t crtime = Int_t((fParam->GetZLength()-AliTPCReconstructor::GetCtgRange()*fRx)/fZWidth-5);
 
   for (Int_t i=2*fMaxTime; i<fMaxBin-2*fMaxTime; i++) {
     b++;