]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - TPC/AliTPCclustererMI.cxx
cout removed
[u/mrichter/AliRoot.git] / TPC / AliTPCclustererMI.cxx
index 2859f24511a0675d7e75d709300c063d216e8f8c..8b481f7b859fcd6a1b14c4fe3753599f10582d6d 100644 (file)
@@ -82,6 +82,8 @@
 #include "AliTPCTransform.h"
 #include "AliTPCclustererMI.h"
 
+using std::cerr;
+using std::endl;
 ClassImp(AliTPCclustererMI)
 
 
@@ -121,7 +123,8 @@ AliTPCclustererMI::AliTPCclustererMI(const AliTPCParam* par, const AliTPCRecoPar
   fUseHLTClusters(4),
   fAllBins(NULL),
   fAllSigBins(NULL),
-  fAllNSigBins(NULL)
+  fAllNSigBins(NULL),
+  fHLTClusterAccess(NULL)
 {
   //
   // COSNTRUCTOR
@@ -165,59 +168,7 @@ AliTPCclustererMI::AliTPCclustererMI(const AliTPCParam* par, const AliTPCRecoPar
     fAllNSigBins[iRow]=0;
   }
 }
-//______________________________________________________________
-AliTPCclustererMI::AliTPCclustererMI(const AliTPCclustererMI &param)
-              :TObject(param),
-  fBins(0),
-  fSigBins(0),
-  fNSigBins(0),
-  fLoop(0),
-  fMaxBin(0),
-  fMaxTime(0),
-  fMaxPad(0),
-  fSector(-1),
-  fRow(-1),
-  fSign(0),
-  fRx(0),
-  fPadWidth(0),
-  fPadLength(0),
-  fZWidth(0),
-  fPedSubtraction(kFALSE),
-  fEventHeader(0),
-  fTimeStamp(0),
-  fEventType(0),
-  fInput(0),
-  fOutput(0),
-  fOutputArray(0),
-  fOutputClonesArray(0),
-  fRowCl(0),
-  fRowDig(0),
-  fParam(0),
-  fNcluster(0),
-  fNclusters(0),
-  fDebugStreamer(0),
-  fRecoParam(0),
-  fBDumpSignal(kFALSE),
-  fBClonesArray(kFALSE),
-  fUseHLTClusters(4),
-  fAllBins(NULL),
-  fAllSigBins(NULL),
-  fAllNSigBins(NULL)
-{
-  //
-  // dummy
-  //
-  fMaxBin = param.fMaxBin;
-}
-//______________________________________________________________
-AliTPCclustererMI & AliTPCclustererMI::operator =(const AliTPCclustererMI & param)
-{
-  //
-  // assignment operator - dummy
-  //
-  fMaxBin=param.fMaxBin;
-  return (*this);
-}
+
 //______________________________________________________________
 AliTPCclustererMI::~AliTPCclustererMI(){
   //
@@ -247,6 +198,7 @@ AliTPCclustererMI::~AliTPCclustererMI(){
   delete [] fAllBins;
   delete [] fAllSigBins;
   delete [] fAllNSigBins;
+  if (fHLTClusterAccess) delete fHLTClusterAccess;
 }
 
 void AliTPCclustererMI::SetInput(TTree * tree)
@@ -319,11 +271,25 @@ void AliTPCclustererMI::MakeCluster(Int_t k,Int_t max,Float_t *bins, UInt_t /*m*
 AliTPCclusterMI &c) 
 {
   //
+  //  Make cluster: characterized by position ( mean-  COG) , shape (RMS) a charge, QMax and Q tot
+  //  Additional correction:
+  //  a) To correct for charge below threshold, in the +1 neghborhood to the max charge charge 
+  //       is extrapolated using gaussian approximation assuming given cluster width.. 
+  //       Additional empirical factor is used to account for the charge fluctuation (kVirtualChargeFactor). 
+  //       Actual value of the  kVirtualChargeFactor should obtained minimimizing residuals between the cluster
+  //       and track interpolation.
+  //  b.) For space points with extended shape (in comparison with expected using parameterization) clusters are 
+  //      unfoded     
+  //  
+  //  NOTE. Actual/Empirical  values for correction are hardwired in the code.
+  //
+  // Input paramters for function:
   //  k    - Make cluster at position k  
   //  bins - 2 D array of signals mapped to 1 dimensional array - 
   //  max  - the number of time bins er one dimension
-  //  c    - refernce to cluster to be filled
+  //  c    - reference to cluster to be filled
   //
+  Double_t kVirtualChargeFactor=0.5;
   Int_t i0=k/max;  //central pad
   Int_t j0=k%max;  //central time bin
 
@@ -349,7 +315,7 @@ AliTPCclusterMI &c)
        Float_t ratio = TMath::Exp(-1.2*TMath::Abs(di)/sigmay2)*TMath::Exp(-1.2*TMath::Abs(dj)/sigmaz2);
        amp = ((matrix[2][0]-2)*(matrix[2][0]-2)/(matrix[-di+2][-dj]+2))*ratio;
        if (amp>2) amp = 2;
-       vmatrix[2+di][2+dj]=amp;
+       vmatrix[2+di][2+dj]= kVirtualChargeFactor*amp;
        vmatrix[2+2*di][2+2*dj]=0;
        if ( (di*dj)!=0){       
          //DIAGONAL ELEMENTS
@@ -720,6 +686,7 @@ void AliTPCclustererMI::Digits2Clusters()
     fParam->Dump();
     fRecoParam->Dump();
   }
+  fRowDig = NULL;
 
   //-----------------------------------------------------------------
   // Use HLT clusters
@@ -1624,7 +1591,7 @@ Int_t AliTPCclustererMI::ReadHLTClusters()
   // used in Digits2Clusters
   //
 
-  TObject* pClusterAccess=NULL;
+  if (!fHLTClusterAccess) {
   TClass* pCl=NULL;
   ROOT::NewFunc_t pNewFunc=NULL;
   do {
@@ -1640,34 +1607,37 @@ Int_t AliTPCclustererMI::ReadHLTClusters()
     AliError("unable to create instance of AliHLTTPCClusterAccessHLTOUT");
     return -2;
   }
-  pClusterAccess=reinterpret_cast<TObject*>(p);
-  if (!pClusterAccess) {
-    AliError("instance not of type TObject");
-    return -3 ;
+  fHLTClusterAccess=reinterpret_cast<TObject*>(p);
   }
 
+  TObject* pClusterAccess=fHLTClusterAccess;
+
   const Int_t kNIS = fParam->GetNInnerSector();
   const Int_t kNOS = fParam->GetNOuterSector();
   const Int_t kNS = kNIS + kNOS;
   fNclusters  = 0;
   
+  // make sure that all clusters from the previous event are cleared
+  pClusterAccess->Clear("event");
   for(fSector = 0; fSector < kNS; fSector++) {
 
     Int_t iResult = 1;
     TString param("sector="); param+=fSector;
-    pClusterAccess->Clear();
+    // prepare for next sector
+    pClusterAccess->Clear("sector");
     pClusterAccess->Execute("read", param, &iResult);
     if (iResult < 0) {
       return iResult;
       AliError("HLT Clusters can not be found");
     }
 
-    if (pClusterAccess->FindObject("clusterarray")==NULL) {
+    TObject* pObj=pClusterAccess->FindObject("clusterarray");
+    if (pObj==NULL) {
       AliError("HLT clusters requested, but not cluster array not present");
       return -4;
     }
 
-    TClonesArray* clusterArray=dynamic_cast<TClonesArray*>(pClusterAccess->FindObject("clusterarray"));
+    TObjArray* clusterArray=dynamic_cast<TClonesArray*>(pObj);
     if (!clusterArray) {
       AliError("HLT cluster array is not of class type TClonesArray");
       return -5;
@@ -1716,7 +1686,7 @@ Int_t AliTPCclustererMI::ReadHLTClusters()
     fNclusters+=nClusterSector;
   } // for(fSector = 0; fSector < kNS; fSector++) {
 
-  delete pClusterAccess;
+  pClusterAccess->Clear("event");
 
   Info("Digits2Clusters", "Number of converted HLT clusters : %d", fNclusters);