]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - HLT/TPCLib/AliHLTTPCClusterFinder.cxx
- added TPCLib to doxygen docu, code corrections according to documentation and effC++
[u/mrichter/AliRoot.git] / HLT / TPCLib / AliHLTTPCClusterFinder.cxx
index 1535b1924282c01f0806376b0c2dc94ce1c58e70..162036448a691a5f4018e3e1a5c98911e40c4be5 100644 (file)
@@ -1,14 +1,32 @@
 // @(#) $Id$
-
-// Author: Anders Vestbo <mailto:vestbo@fi.uib.no>, 
-//         Constantin Loizides <mailto:loizides@ikf.uni-frankfurt.de>
-//         Jochen Thaeder <mailto:thaeder@kip.uni-heidelberg.de>
-
-//*-- Copyright &copy ALICE HLT Group
-
+// Original: AliL3ClustFinderNew.cxx,v 1.29 2005/06/14 10:55:21 cvetan Exp 
+
+/**************************************************************************
+ * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
+ *                                                                        *
+ * Authors: Anders Vestbo <mailto:vestbo@fi.uib.no>,                     *
+ *          Constantin Loizides <mailto:loizides@ikf.uni-frankfurt.de>    *
+ *          Jochen Thaeder <mailto:thaeder@kip.uni-heidelberg.de>         *
+ *          for The ALICE Off-line Project.                               *
+ *                                                                        *
+ * Permission to use, copy, modify and distribute this software and its   *
+ * documentation strictly for non-commercial purposes is hereby granted   *
+ * without fee, provided that the above copyright notice appears in all   *
+ * copies and that both the copyright notice and this permission notice   *
+ * appear in the supporting documentation. The authors make no claims     *
+ * about the suitability of this software for any purpose. It is          *
+ * provided "as is" without express or implied warranty.                  *
+ **************************************************************************/
+
+/** @file   AliHLTTPCClusterFinder.cxx
+    @author Anders Vestbo <mailto:vestbo@fi.uib.no>,                
+           Constantin Loizides <mailto:loizides@ikf.uni-frankfurt.de>
+           Jochen Thaeder <mailto:thaeder@kip.uni-heidelberg.de>     
+    @date   
+    @brief  Cluster Finder for the TPC
+*/
 
 #include "AliHLTTPCDigitReader.h"
-#include "AliHLTTPCStandardIncludes.h"
 #include "AliHLTTPCRootTypes.h"
 #include "AliHLTTPCLogging.h"
 #include "AliHLTTPCClusterFinder.h"
@@ -16,6 +34,7 @@
 #include "AliHLTTPCTransform.h"
 #include "AliHLTTPCSpacePointData.h"
 #include "AliHLTTPCMemHandler.h"
+#include "AliHLTTPCPad.h"
 
 #if __GNUC__ >= 3
 using namespace std;
@@ -95,21 +114,61 @@ using namespace std;
 ClassImp(AliHLTTPCClusterFinder)
 
 AliHLTTPCClusterFinder::AliHLTTPCClusterFinder()
+  :
+  fMatch(1),
+  fThreshold(10),
+  fSignalThreshold(-1),
+  fOccupancyLimit(1.0),
+  fXYErr(0.2),
+  fZErr(0.3),
+  fDeconvPad(kTRUE),
+  fDeconvTime(kTRUE),
+  fStdout(kFALSE),
+  fCalcerr(kTRUE),
+  fRawSP(kFALSE),
+  fFirstRow(0),
+  fLastRow(0),
+  fDigitReader(NULL)
 {
   //constructor
-  fMatch = 1;
-  fThreshold = 10;
-  fXYErr = 0.2;
-  fZErr = 0.3;
-  fDeconvPad = kTRUE;
-  fDeconvTime = kTRUE;
-  fStdout = kFALSE;
-  fCalcerr = kTRUE;
-  fRawSP = kFALSE;
-  fFirstRow=0;
-  fLastRow=0;
-  fDigitReader = 0;
+}
+
+AliHLTTPCClusterFinder::AliHLTTPCClusterFinder(const AliHLTTPCClusterFinder& src)
+  :
+  fMatch(src.fMatch),
+  fThreshold(src.fThreshold),
+  fSignalThreshold(src.fSignalThreshold),
+  fOccupancyLimit(src.fOccupancyLimit),
+  fXYErr(src.fXYErr),
+  fZErr(src.fZErr),
+  fDeconvPad(src.fDeconvPad),
+  fDeconvTime(src.fDeconvTime),
+  fStdout(src.fStdout),
+  fCalcerr(src.fCalcerr),
+  fRawSP(src.fRawSP),
+  fFirstRow(src.fFirstRow),
+  fLastRow(src.fLastRow),
+  fDigitReader(src.fDigitReader)
+{
+}
 
+AliHLTTPCClusterFinder& AliHLTTPCClusterFinder::operator=(const AliHLTTPCClusterFinder& src)
+{
+  fMatch=src.fMatch;
+  fThreshold=src.fThreshold;
+  fSignalThreshold=src.fSignalThreshold;
+  fOccupancyLimit=src.fOccupancyLimit;
+  fXYErr=src.fXYErr;
+  fZErr=src.fZErr;
+  fDeconvPad=src.fDeconvPad;
+  fDeconvTime=src.fDeconvTime;
+  fStdout=src.fStdout;
+  fCalcerr=src.fCalcerr;
+  fRawSP=src.fRawSP;
+  fFirstRow=src.fFirstRow;
+  fLastRow=src.fLastRow;
+  fDigitReader=src.fDigitReader;
+  return (*this);
 }
 
 AliHLTTPCClusterFinder::~AliHLTTPCClusterFinder()
@@ -156,18 +215,21 @@ void AliHLTTPCClusterFinder::ProcessDigits()
   bool readValue = true;
   Int_t newRow = 0;    
   Int_t rowOffset = 0;
-  UChar_t pad;
-  UShort_t time,newTime=0;
-  UInt_t charge,newPad=0;
+  UShort_t time=0,newTime=0;
+  UInt_t pad=0,newPad=0;
+  AliHLTTPCSignal_t charge=0;
 
   fNClusters = 0;
 
   // initialize block for reading packed data
-  fDigitReader->InitBlock(fPtr,fSize,fFirstRow,fLastRow);
+  fDigitReader->InitBlock(fPtr,fSize,fFirstRow,fLastRow,fCurrentPatch,fCurrentSlice);
   readValue = fDigitReader->Next();
 
-  if (!readValue)return;
+  // Matthias 08.11.2006 the following return would cause termination without writing the
+  // ClusterData and thus would block the component. I just wnt to have the commented line
+  // here for information
+  //if (!readValue)return;
+
   pad = fDigitReader->GetPad();
   time = fDigitReader->GetTime();
   fCurrentRow = fDigitReader->GetRow();
@@ -178,9 +240,11 @@ void AliHLTTPCClusterFinder::ProcessDigits()
   fCurrentRow += rowOffset;
 
   UInt_t lastpad = 123456789;
-  AliClusterData *pad1[5000]; //2 lists for internal memory=2pads
-  AliClusterData *pad2[5000]; //2 lists for internal memory=2pads
-  AliClusterData clusterlist[10000]; //Clusterlist
+  const Int_t kPadArraySize=5000;
+  const Int_t kClusterListSize=10000;
+  AliClusterData *pad1[kPadArraySize]; //2 lists for internal memory=2pads
+  AliClusterData *pad2[kPadArraySize]; //2 lists for internal memory=2pads
+  AliClusterData clusterlist[kClusterListSize]; //Clusterlist
 
   AliClusterData **currentPt;  //List of pointers to the current pad
   AliClusterData **previousPt; //List of pointers to the previous pad
@@ -188,6 +252,22 @@ void AliHLTTPCClusterFinder::ProcessDigits()
   previousPt = pad1;
   UInt_t nprevious=0,ncurrent=0,ntotal=0;
 
+  /* quick implementation of baseline calculation and zero suppression
+     open a pad object for each pad and delete it after processing.
+     later a list of pad objects with base line history can be used
+     The whole thing only works if we really get unprocessed raw data, if
+     the data is already zero suppressed, there might be gaps in the time
+     bins.
+   */
+  Int_t gatingGridOffset=50;
+  AliHLTTPCPad baseline(gatingGridOffset, AliHLTTPCTransform::GetNTimeBins());
+  // just to make later conversion to a list of objects easier
+  AliHLTTPCPad* pCurrentPad=NULL;
+  if (fSignalThreshold>=0) {
+    pCurrentPad=&baseline;
+    baseline.SetThreshold(fSignalThreshold);
+  }
+
   while ( readValue ){   // Reads through all digits in block
 
     if(pad != lastpad){
@@ -229,17 +309,63 @@ void AliHLTTPCClusterFinder::ProcessDigits()
       lastwas_falling = 0;
     }
 
+    while(1){ //Loop over time bins of current pad
+      // read all the values for one pad at once to calculate the base line
+      if (pCurrentPad) {
+       if (!pCurrentPad->IsStarted()) {
+         //HLTDebug("reading data for pad %d, padrow %d", fDigitReader->GetPad(), fDigitReader->GetRow()+rowOffset);
+         pCurrentPad->SetID(fDigitReader->GetRow()+rowOffset,fDigitReader->GetPad());
+         if ((pCurrentPad->StartEvent())>=0) {
+           do {
+             if ((fDigitReader->GetRow()+rowOffset)!=pCurrentPad->GetRowNumber()) break;
+             if (fDigitReader->GetPad()!=pCurrentPad->GetPadNumber()) break;
+             pCurrentPad->SetRawData(fDigitReader->GetTime(), fDigitReader->GetSignal());
+             //HLTDebug("set raw data to pad: bin %d charge %d", fDigitReader->GetTime(), fDigitReader->GetSignal());
+           } while ((readValue = fDigitReader->Next())!=0);
+         }
+         pCurrentPad->CalculateBaseLine(AliHLTTPCTransform::GetNTimeBins()/2);
+         if (pCurrentPad->Next(kTRUE/*do zero suppression*/)==0) {
+           //HLTDebug("no data available after zero suppression");
+           pCurrentPad->StopEvent();
+           pCurrentPad->ResetHistory();
+           break;
+         }
+         time=pCurrentPad->GetCurrentPosition();
+         if (time>pCurrentPad->GetSize()) {
+           HLTError("invalid time bin for pad");
+           break;
+         }
+       }
+      }
+
+      if (pCurrentPad) {
+       Float_t occupancy=pCurrentPad->GetOccupancy();
+       //HLTDebug("pad %d occupancy level: %f", pCurrentPad->GetPadNumber(), occupancy);
+       if ( occupancy < fOccupancyLimit ) {
+         charge = pCurrentPad->GetCorrectedData();
+       } else {
+         charge = 0;
+         //HLTDebug("ignoring pad %d with occupancy level %f", pCurrentPad->GetPadNumber(), occupancy);
+       }
+      } else {
+       charge = fDigitReader->GetSignal();
+      }
+      //HLTDebug("get next charge value: position %d charge %d", time, charge);
+
+
+      // CHARGE DEBUG
+      if (fDigitReader->GetRow() == 90){
+/////    LOG(AliHLTTPCLog::kFatal,"AliHLTTPCClusterFinder::Row","row90")  << "PAD=" <<  fDigitReader->GetPad() << "  TIME=" <<  fDigitReader->GetTime() 
+         //                                       << "  SIGNAL=" <<  fDigitReader->GetSignal() << ENDLOG;
 
-    // LOOP OVER CURRENR SEQUENCE
-    while(1){ //Loop over current
-      charge = fDigitReader->GetSignal();
+      }
 
       if(time >= AliHLTTPCTransform::GetNTimeBins()){
-       LOG(AliHLTTPCLog::kFatal,"AliHLTTPCClusterFinder::ProcessRow","Digits")
-         <<"Timebin out of range "<<(Int_t)time<<ENDLOG;
+       HLTWarning("Timebin (%d) out of range (%d)", time, AliHLTTPCTransform::GetNTimeBins());
        break;
       }
-      
+
+
       //Get the current ADC-value
       if(fDeconvTime){
 
@@ -259,14 +385,31 @@ void AliHLTTPCClusterFinder::ProcessDigits()
       seqaverage += time*charge;
       seqerror += time*time*charge;
       
+      if (pCurrentPad) {
+       
+       if((pCurrentPad->Next(kTRUE/*do zero suppression*/))==0) {
+         pCurrentPad->StopEvent();
+         pCurrentPad->ResetHistory();
+         if(readValue) {
+           newPad = fDigitReader->GetPad();
+           newTime = fDigitReader->GetTime();
+           newRow = fDigitReader->GetRow() + rowOffset;
+         }
+         break;
+       }
+
+       newPad=pCurrentPad->GetPadNumber();
+       newTime=pCurrentPad->GetCurrentPosition();
+       newRow=pCurrentPad->GetRowNumber();
+      } else {
       readValue = fDigitReader->Next();
-      
       //Check where to stop:
       if(!readValue) break; //No more value
 
       newPad = fDigitReader->GetPad();
       newTime = fDigitReader->GetTime();
       newRow = fDigitReader->GetRow() + rowOffset;
+      }
 
       if(newPad != pad)break; //new pad
       if(newTime != time+1) break; //end of sequence
@@ -276,6 +419,12 @@ void AliHLTTPCClusterFinder::ProcessDigits()
 
     }//end loop over sequence
 
+    //HLTDebug("ended time bin sequence loop: seqcharge=%d readValue=%d", seqcharge, readValue);
+    //HLTDebug("pad=%d newpad=%d current row=%d newrow=%d", pad, newPad, fCurrentRow, newRow);
+    if (seqcharge<=0) {
+      // with active zero suppression zero values are possible
+      continue;
+    }
 
     //Calculate mean of sequence:
     Int_t seqmean=0;
@@ -294,7 +443,7 @@ void AliHLTTPCClusterFinder::ProcessDigits()
 
 
     //Compare with results on previous pad:
-    for(UInt_t p=0; p<nprevious; p++){
+    for(UInt_t p=0; p<nprevious && p<kPadArraySize && ncurrent<kPadArraySize; p++){
       
       //dont merge sequences on the same pad twice
       if(previousPt[p]->fLastMergedPad==pad) continue;
@@ -336,7 +485,7 @@ void AliHLTTPCClusterFinder::ProcessDigits()
     } //Loop over results on previous pad.
 
 
-    if(newcluster){
+    if(newcluster && ncurrent<kPadArraySize){
       //Start a new cluster. Add it to the clusterlist, and update
       //the list of pointers to clusters in current pad.
       //current pad will be previous pad on next pad.
@@ -368,13 +517,17 @@ void AliHLTTPCClusterFinder::ProcessDigits()
   
     // to prevent endless loop  
     if(time >= AliHLTTPCTransform::GetNTimeBins()){
-      LOG(AliHLTTPCLog::kFatal,"AliHLTTPCClusterFinder::ProcessRow","Digits")
-       <<"Timebin out of range "<<(Int_t)time<<ENDLOG;
+      HLTWarning("Timebin (%d) out of range (%d)", time, AliHLTTPCTransform::GetNTimeBins());
       break;
     }
 
 
     if(!readValue) break; //No more value
+    
+    if (ntotal>=kClusterListSize || ncurrent>=kPadArraySize) {
+      HLTWarning("pad array size exceeded ntotal=%d ncurrent=%d, skip rest of the data", ntotal, ncurrent);
+      break;
+    }
 
     if(fCurrentRow != newRow){
       WriteClusters(ntotal,clusterlist);
@@ -392,14 +545,12 @@ void AliHLTTPCClusterFinder::ProcessDigits()
 
     pad = newPad;
     time = newTime;
-  
+
   } // END while(readValue)
 
   WriteClusters(ntotal,clusterlist);
 
-  LOG(AliHLTTPCLog::kInformational,"AliHLTTPCClusterFinder::ProcessDigits","Space points") 
-    << "ClusterFinder found " << fNClusters << " clusters in slice " << fCurrentSlice << " patch " 
-    << fCurrentPatch << ENDLOG;
+  HLTInfo("ClusterFinder found %d clusters in slice %d patch %d", fNClusters, fCurrentSlice, fCurrentPatch);
 
 } // ENDEND
 
@@ -420,6 +571,11 @@ void AliHLTTPCClusterFinder::WriteClusters(Int_t nclusters,AliClusterData *list)
       Float_t ftime =(Float_t)list[j].fTime / list[j].fTotalCharge;
       Float_t ftime2=fZErr*fZErr;  //fixed given error
 
+
+   
+     
+
+
       if(fCalcerr) { //calc the errors, otherwice take the fixed error 
        Int_t patch = AliHLTTPCTransform::GetPatch(fCurrentRow);
        UInt_t q2=list[j].fTotalCharge*list[j].fTotalCharge;
@@ -487,6 +643,9 @@ void AliHLTTPCClusterFinder::WriteClusters(Int_t nclusters,AliClusterData *list)
       fSpacePointData[counter].fSigmaY2 = fpad2;
       fSpacePointData[counter].fSigmaZ2  = ftime2;
 
+      fSpacePointData[counter].fUsed = kFALSE;         // only used / set in AliHLTTPCDisplay
+      fSpacePointData[counter].fTrackN = -1;           // only used / set in AliHLTTPCDisplay
+
       Int_t patch=fCurrentPatch;
       if(patch==-1) patch=0; //never store negative patch number
       fSpacePointData[counter].fID = counter