Implement fixed due to Boris Batyunya. Including a fix for a bug with the
[u/mrichter/AliRoot.git] / ITS / AliITSClusterFinderSPD.cxx
index 28be4826125453483b933c2abda0a7051f71ce05..d7cf5fa2c6153fb0a64bd9db985f41a4e3c823c0 100644 (file)
@@ -15,6 +15,8 @@
 
 
 #include "AliITSClusterFinderSPD.h"
+#include "AliITSMapA1.h"
+#include "AliITS.h"
 #include "AliRun.h"
 
 
@@ -44,13 +46,21 @@ AliITSClusterFinderSPD::AliITSClusterFinderSPD()
   fDigits=0;
   fClusters=0;
   fNclusters=0;
+  fMap=0;
   SetDx();
   SetDz();
-  SetMap();
   SetNCells();
   
 }
 
+//_____________________________________________________________________________
+AliITSClusterFinderSPD::~AliITSClusterFinderSPD()
+{
+  // destructor
+  if (fMap) delete fMap;
+
+
+}
 //__________________________________________________________________________
 AliITSClusterFinderSPD::AliITSClusterFinderSPD(const AliITSClusterFinderSPD &source){
   //     Copy Constructor 
@@ -82,8 +92,9 @@ AliITSClusterFinderSPD&
 void AliITSClusterFinderSPD::SetMap()
 {
   // set map
+
   if(!fMap) fMap=new AliITSMapA1(fSegmentation,fDigits);
-  
+
 }
 
 //_____________________________________________________________________________
@@ -101,6 +112,7 @@ void AliITSClusterFinderSPD::Find1DClusters()
   
   // read in digits -> do not apply threshold 
   // signal in fired pixels is always 1
+
   fMap->FillMap();
   
   Int_t nofFoundClusters = 0;
@@ -122,13 +134,13 @@ void AliITSClusterFinderSPD::Find1DClusters()
        id = it+m;
        if(id >= fMaxNofSamples) break;    // ! no possible for the fadc 
        
-       if(fMap->GetHitIndex(k,id)) {   // start of the cluster
+       if(fMap->TestHit(k,id) == kUnused) {   // start of the cluster
          lclx += 1;
          if(lclx == 1) xstart = id;
          
        }
        
-       if(lclx > 0 && !fMap->GetHitIndex(k,id)) {  
+       if(lclx > 0 && fMap->TestHit(k,id) == kEmpty) {  
          // end of cluster if a gap exists
          xstop = id-1;
          ilcl = 1;
@@ -143,14 +155,14 @@ void AliITSClusterFinderSPD::Find1DClusters()
       if(id >= fMaxNofSamples && lclx == 0) break; // the x row finished
       
       if(id < fMaxNofSamples && ilcl == 0 && lclx > 0) {  
-       // cluster end is outside of the window,
+                                  // cluster end is outside of the window,
        mmax += 5;                 // increase mmax and repeat the cluster
-       // finding
+                                  // finding
        it -= 1;
       }
       
       if(id >= fMaxNofSamples && lclx > 0) {  // the x row finished but
-       xstop = fMaxNofSamples - 1;          // the end cluster exists
+       xstop = fMaxNofSamples - 1;           // the end cluster exists
        ilcl = 1;
       } 
       
@@ -161,11 +173,10 @@ void AliITSClusterFinderSPD::Find1DClusters()
        mmax = 10;
            nofFoundClusters++;
            Float_t clusterCharge = 0.;
-            // get this from segmentation when this will be implemented
             Float_t zpitch = fSegmentation->Dpz(k+1); 
            Float_t clusterZ, dummyX; 
             Int_t dummy=0;
-            fSegmentation->GetCellCxz(dummy,k,dummyX,clusterZ);
+            fSegmentation->GetPadCxz(dummy,k,dummyX,clusterZ);
             Float_t zstart = clusterZ - 0.5*zpitch;
             Float_t zstop = clusterZ + 0.5*zpitch;
            Float_t clusterX = 0.;
@@ -177,7 +188,7 @@ void AliITSClusterFinderSPD::Find1DClusters()
             Int_t its;
            for(its=xstart; its<=xstop; its++) {
               Int_t firedpixel=0;
-              if (fMap->GetHitIndex(k,its)) firedpixel=1; 
+              if (fMap->GetHitIndex(k,its)>=0) firedpixel=1; 
              clusterCharge += firedpixel;
               clusterX +=its + 0.5;
            }
@@ -204,6 +215,8 @@ void AliITSClusterFinderSPD::Find1DClusters()
       }    // new cluster (ilcl=1)
     } // X direction loop (it)
   } // Z direction loop (k)
+
+  //fMap->ClearMap();
   return;
   
 }
@@ -217,12 +230,11 @@ void  AliITSClusterFinderSPD::GroupClusters()
   // get number of clusters for this module
   Int_t nofClusters = fClusters->GetEntriesFast();
   nofClusters -= fNclusters;
-  //printf("Group: fNclusters nofClusters %d %d\n",fNclusters, nofClusters);
   
   AliITSRawClusterSPD *clusterI;
   AliITSRawClusterSPD *clusterJ;
   
-  Int_t *label = new Int_t [nofClusters];  // activate this for DEC machines
+  Int_t *label=new Int_t[nofClusters];  
   Int_t i,j;
   for(i=0; i<nofClusters; i++) label[i] = 0;
   for(i=0; i<nofClusters; i++) {
@@ -237,8 +249,8 @@ void  AliITSClusterFinderSPD::GroupClusters()
        //    if((clusterI->XStop() == clusterJ->XStart()-1)||(clusterI->XStart()==clusterJ->XStop()+1)) cout<<"!! Diagonal cluster"<<endl;
        /*    
              cout << "clusters " << i << "," << j << " before grouping" << endl;
-             clusterI->Print();
-             clusterJ->Print();
+             clusterI->PrintInfo();
+             clusterJ->PrintInfo();
        */    
        clusterI->Add(clusterJ);
        //        cout << "remove cluster " << j << endl;
@@ -246,25 +258,114 @@ void  AliITSClusterFinderSPD::GroupClusters()
        fClusters->RemoveAt(j);
        /*
          cout << "cluster  " << i << " after grouping" << endl;
-         clusterI->Print();
+         clusterI->PrintInfo();
        */
       }  // pair
     } // J clusters  
     label[i] = 1;
   } // I clusters
   fClusters->Compress();
-  //Int_t totalNofClusters = fClusters->GetEntriesFast();
-  //cout << " Nomber of clusters at the group end ="<< totalNofClusters<<endl;
-  
-  delete [] nofClusters;
-  delete nofClusters;
+  //  Int_t totalNofClusters = fClusters->GetEntriesFast();
+  //  cout << " Nomber of clusters at the group end ="<< totalNofClusters<<endl;
   
+  delete [] label;
+
   return;
   
   
 }
 //_____________________________________________________________________________
 
+void AliITSClusterFinderSPD::TracksInCluster()
+{
+  
+  // Find tracks creating one cluster
+
+  // get number of clusters for this module
+  Int_t nofClusters = fClusters->GetEntriesFast();
+  nofClusters -= fNclusters;
+
+  Int_t i, ix, iz, jx, jz, xstart, xstop, zstart, zstop, nclx, nclz;
+  //  Int_t signal, track0, track1, track2;
+  Int_t trmax = 100;
+  Int_t cltracks[trmax], itr, tracki, ii, is, js, ie, ntr, tr0, tr1, tr2;
+
+  for(i=0; i<nofClusters; i++) { 
+    ii = 0;
+    memset(cltracks,-1,sizeof(int)*trmax);
+    tr0=tr1=tr2=-1;
+
+    AliITSRawClusterSPD *clusterI = (AliITSRawClusterSPD*) fClusters->At(i);
+
+    nclx = clusterI->NclX();
+    nclz = clusterI->NclZ();
+    xstart = clusterI->XStartf();
+    xstop = clusterI->XStopf();
+    zstart = clusterI->Zend()-nclz+1;
+    zstop = clusterI->Zend();
+
+    Int_t ind; 
+
+     for(iz=0; iz<nclz; iz++) { 
+         jz = zstart + iz;
+       for(ix=0; ix<nclx; ix++) { 
+        jx = xstart + ix;
+        ind = fMap->GetHitIndex(jz,jx);
+        if(ind == 0 && iz >= 0 && ix > 0) {
+          continue;
+         }
+        if(ind == 0 && iz > 0 && ix >= 0) {
+          continue;
+         }
+        if(ind == 0 && iz == 0 && ix == 0 && i > 0) {
+          continue;
+         }
+
+        AliITSdigitSPD *dig = (AliITSdigitSPD*)fMap->GetHit(jz,jx);
+       /*
+         signal=dig->fSignal;
+         track0=dig->fTracks[0];
+         track1=dig->fTracks[1];
+         track2=dig->fTracks[2];
+       */
+          for(itr=0; itr<3; itr++) { 
+           tracki = dig->fTracks[itr];
+            if(tracki >= 0) {
+             ii += 1;
+             cltracks[ii-1] = tracki;
+            }
+         }
+       } // ix pixel
+     }  // iz pixel
+     for(is=0; is<trmax; is++) { 
+         if(cltracks[is]<0) continue;
+       for(js=is+1; js<trmax; js++) { 
+         if(cltracks[js]<0) continue;
+         if(cltracks[js]==cltracks[is]) cltracks[js]=-5;
+       }
+     }
+
+     ntr = 0;
+     for(ie=0; ie<trmax; ie++) { 
+       if(cltracks[ie] >= 0) {
+        ntr=ntr+1;
+        if(ntr==1) tr0=cltracks[ie];
+        if(ntr==2) tr1=cltracks[ie];
+        if(ntr==3) tr2=cltracks[ie];
+       }
+     }
+     // if delta ray only
+     if(ntr == 0) ntr = 1;
+
+     clusterI->SetNTracks(ntr);
+     clusterI->SetTracks(tr0,tr1,tr2);
+
+  } // I cluster
+
+}
+//_____________________________________________________________________________
+
 void AliITSClusterFinderSPD::GetRecPoints()
 {
   // get rec points
@@ -273,25 +374,35 @@ void AliITSClusterFinderSPD::GetRecPoints()
   // get number of clusters for this module
   Int_t nofClusters = fClusters->GetEntriesFast();
   nofClusters -= fNclusters;
-  //printf("GetRecP: fNclusters nofClusters %d %d\n",fNclusters, nofClusters);
-  
   const Float_t kconv = 1.0e-4;
   const Float_t kRMSx = 12.0*kconv; // microns -> cm ITS TDR Table 1.3
   const Float_t kRMSz = 70.0*kconv; // microns -> cm ITS TDR Table 1.3
 
+  Float_t spdLength = fSegmentation->Dz();
+  Float_t spdWidth = fSegmentation->Dx();
+
   Int_t i;
+  Int_t track0, track1, track2;
+
   for(i=0; i<nofClusters; i++) { 
+
     AliITSRawClusterSPD *clusterI = (AliITSRawClusterSPD*) fClusters->At(i);
+    clusterI->GetTracks(track0, track1, track2); 
     AliITSRecPoint rnew;
-    rnew.SetX(clusterI->X()*kconv);
-    rnew.SetZ(clusterI->Z()*kconv);
+
+    rnew.SetX((clusterI->X() - spdWidth/2)*kconv);
+    rnew.SetZ((clusterI->Z() - spdLength/2)*kconv);
     rnew.SetQ(1.);
     rnew.SetdEdX(0.);
     rnew.SetSigmaX2(kRMSx*kRMSx);
     rnew.SetSigmaZ2(kRMSz*kRMSz);
-    rnew.SetProbability(1.);
+    rnew.fTracks[0]=track0;
+    rnew.fTracks[1]=track1;
+    rnew.fTracks[2]=track2;
     iTS->AddRecPoint(rnew);
   } // I clusters
+
+  fMap->ClearMap();
   
 }
 //_____________________________________________________________________________
@@ -301,6 +412,8 @@ void AliITSClusterFinderSPD::FindRawClusters()
   // find raw clusters
   Find1DClusters();
   GroupClusters();
+  TracksInCluster();
   GetRecPoints();
+
 }