]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - ITS/AliITSclustererV2.cxx
Modifications to reduce the size in memory (E.Crescio)
[u/mrichter/AliRoot.git] / ITS / AliITSclustererV2.cxx
index 9b10f5017a0b5923c275bc8bb3d5f58e619dae3d..767099ceb1de0c383fd35509776f8d4b71f0ff1f 100644 (file)
@@ -10,6 +10,7 @@
 
 #include "AliITSclustererV2.h"
 #include "AliITSclusterV2.h"
+#include "AliITSDetTypeRec.h"
 #include "AliRawReader.h"
 #include "AliITSRawStreamSPD.h"
 #include "AliITSRawStreamSDD.h"
@@ -91,6 +92,7 @@ AliITSclustererV2::AliITSclustererV2(const AliITSgeom *geom) {
   fHlSSD=2.00;
   fTanP=0.0275;
   fTanN=0.0075;
+
 }
 
 
@@ -118,6 +120,10 @@ Int_t AliITSclustererV2::Digits2Clusters(TTree *dTree, TTree *cTree) {
   else branch->SetAddress(&clusters);
 
   Int_t mmax=(Int_t)dTree->GetEntries();
+  if (mmax!=fNModules) {
+    Error("Digits2Clusters","Number of entries != number of modules !");
+    return 1;
+  }
 
   for (fI=0; fI<mmax; fI++) {
     dTree->GetEvent(fI);
@@ -219,13 +225,13 @@ void AliITSclustererV2::Digits2Clusters(AliRawReader* rawReader) {
 //**** Fast clusters *******************************
 #include "TParticle.h"
 
-#include "AliITS.h"
+//#include "AliITS.h"
 #include "AliITSmodule.h"
 #include "AliITSRecPoint.h"
 #include "AliITSsimulationFastPoints.h"
 #include "AliITSRecPoint.h"
 
-/*
+
 static void CheckLabels(Int_t lab[3]) {
   //------------------------------------------------------------
   // Tries to find mother's labels
@@ -252,12 +258,15 @@ static void CheckLabels(Int_t lab[3]) {
         else ;//cerr<<"CheckLabels : No empty labels !\n";
     }
 }
-*/
+
+/*
 static void CheckLabels(Int_t lab[3]) {
   //------------------------------------------------------------
   // Tries to find mother's labels
   //------------------------------------------------------------
 
+  if(lab[0]<0 && lab[1]<0 && lab[2]<0) return; // In case of no labels just exit
+
   Int_t ntracks = gAlice->GetMCApp()->GetNtrack();
   for (Int_t i=0;i<3;i++){
     Int_t label = lab[i];
@@ -277,15 +286,16 @@ static void CheckLabels(Int_t lab[3]) {
   }
   
 }
-
+*/
 static void CheckLabels2(Int_t lab[10]) {
   //------------------------------------------------------------
   // Tries to find mother's labels
   //------------------------------------------------------------
-
-  Int_t ntracks = gAlice->GetMCApp()->GetNtrack();
   Int_t nlabels =0; 
   for (Int_t i=0;i<10;i++) if (lab[i]>=0) nlabels++;
+  if(nlabels == 0) return; // In case of no labels just exit
+
+  Int_t ntracks = gAlice->GetMCApp()->GetNtrack();
 
   for (Int_t i=0;i<10;i++){
     Int_t label = lab[i];
@@ -346,10 +356,13 @@ static void CheckLabels2(Int_t lab[10]) {
 }
 
 static void AddLabel(Int_t lab[10], Int_t label) {
+
+  if(label<0) return; // In case of no label just exit
+
   Int_t ntracks = gAlice->GetMCApp()->GetNtrack();
   if (label>ntracks) return;
   for (Int_t i=0;i<10;i++){
-    if (label<0) break;
+    //    if (label<0) break;
     if (lab[i]==label) break;
     if (lab[i]<0) {
       lab[i]= label;
@@ -383,55 +396,6 @@ void AliITSclustererV2::RecPoints2Clusters
   }  
 } 
 
-Int_t AliITSclustererV2::Hits2Clusters(TTree *hTree, TTree *cTree) {
-  //------------------------------------------------------------
-  // This function creates ITS clusters
-  //------------------------------------------------------------
-  if (!gAlice) {
-     Error("Hits2Clusters","gAlice==0 !");
-     return 1;
-  }
-
-  AliITS *its  = (AliITS*)gAlice->GetModule("ITS");
-  if (!its) { 
-     Error("Hits2Clusters","Can't find the ITS !"); 
-     return 2; 
-  }
-  AliITSgeom *geom=its->GetITSgeom();
-  Int_t mmax=geom->GetIndexMax();
-
-  its->InitModules(-1,mmax);
-  its->FillModules(hTree,0);
-
-  TClonesArray *clusters=new TClonesArray("AliITSclusterV2",1000);
-  TBranch *branch=cTree->GetBranch("Clusters");
-  if (!branch) cTree->Branch("Clusters",&clusters);
-  else branch->SetAddress(&clusters);
-
-  static TClonesArray *points=its->RecPoints();
-  AliITSsimulationFastPoints sim;
-  Int_t ncl=0;
-  for (Int_t m=0; m<mmax; m++) {
-    AliITSmodule *mod=its->GetModule(m);      
-    sim.CreateFastRecPoints(mod,m,gRandom);      
-
-    RecPoints2Clusters(points, m, clusters);
-    its->ResetRecPoints();
-
-    ncl+=clusters->GetEntriesFast();
-    cTree->Fill();
-    clusters->Clear();
-  }
-
-  Info("Hits2Clusters","Number of found fast clusters : %d",ncl);
-
-  //cTree->Write();
-
-  delete clusters;
-
-  return 0;
-}
-
 //***********************************
 
 #ifndef V1
@@ -524,7 +488,7 @@ FindClustersSPD(const TClonesArray *digits, TClonesArray *clusters) {
        for (Int_t iy=ymin; iy<=ymax;iy+=idy){
         //
         Int_t ndigits =0;
-        Double_t y=0.,z=0.,q=0.;        
+        Float_t y=0.,z=0.,q=0.;         
         for (Int_t l=0; l<ni; l++) {
           d=(AliITSdigitSPD*)digits->UncheckedAt(idx[l]);
           if (zmax-zmin>=idz || ymax-ymin>=idy){
@@ -532,7 +496,7 @@ FindClustersSPD(const TClonesArray *digits, TClonesArray *clusters) {
             if (TMath::Abs( d->GetCoord1()-iz)>0.75*idz) continue;
           }
           ndigits++;
-          Double_t qq=d->GetSignal();
+          Float_t qq=d->GetSignal();
           y+=qq*fYSPD[d->GetCoord2()]; z+=qq*fZSPD[d->GetCoord1()]; q+=qq;   
          
         }     
@@ -576,7 +540,9 @@ void AliITSclustererV2::FindClustersSPD(AliITSRawStream* input,
   Int_t kNzBins = fNzSPD + 2;
   Int_t kNyBins = fNySPD + 2;
   Int_t kMaxBin = kNzBins * kNyBins;
-  AliBin* bins = NULL;
+  AliBin *binsSPD = new AliBin[kMaxBin];
+  AliBin *binsSPDInit = new AliBin[kMaxBin];
+  AliBin *bins = NULL;
 
   // read raw data input stream
   while (kTRUE) {
@@ -610,9 +576,6 @@ void AliITSclustererV2::FindClustersSPD(AliITSRawStream* input,
          Int_t ymax = ymin;
          Int_t zmin = (idxBins[0] % kNzBins) - 1;
          Int_t zmax = zmin;
-         Float_t y = 0.;
-         Float_t z = 0.;
-         Float_t q = 0.;
          for (Int_t idx = 0; idx < nBins; idx++) {
            Int_t iy = (idxBins[idx] / kNzBins) - 1;
            Int_t iz = (idxBins[idx] % kNzBins) - 1;
@@ -620,39 +583,59 @@ void AliITSclustererV2::FindClustersSPD(AliITSRawStream* input,
            if (ymax < iy) ymax = iy;
            if (zmin > iz) zmin = iz;
            if (zmax < iz) zmax = iz;
-
-           Float_t qBin = bins[idxBins[idx]].GetQ();
-           y += qBin * fYSPD[iy]; 
-           z += qBin * fZSPD[iz]; 
-           q += qBin;   
          }
-         y /= q; 
-         z /= q;
-         y -= fHwSPD; 
-         z -= fHlSPD;
-
-         Float_t hit[5];  // y, z, sigma(y)^2, sigma(z)^2, charge
-         hit[0] = -y-fYshift[iModule]; 
-         if (iModule <= fLastSPD1) hit[0] = -hit[0];
-         hit[1] = z+fZshift[iModule];
-         hit[2] = fYpitchSPD*fYpitchSPD/12.;
-         hit[3] = fZ1pitchSPD*fZ1pitchSPD/12.;
-//       hit[4] = q;
-         hit[4] = (zmax-zmin+1)*100 + (ymax-ymin+1);
-
-         CheckLabels(label);
-         Int_t info[3]={0,0,0};
-         new (clusters[iModule]->AddrAt(nClusters)) 
-           AliITSclusterV2(label, hit,info); 
-         nClusters++;
+
+         Int_t idy = 0; // max 2 clusters
+         if ((iModule <= fLastSPD1) &&idy<3) idy=3;
+         if ((iModule > fLastSPD1) &&idy<4) idy=4; 
+         Int_t idz =3;
+         for (Int_t iiz=zmin; iiz<=zmax;iiz+=idz)
+           for (Int_t iiy=ymin; iiy<=ymax;iiy+=idy){
+             //
+             Int_t ndigits =0;
+             Float_t y=0.,z=0.,q=0.;    
+             for (Int_t idx = 0; idx < nBins; idx++) {
+               Int_t iy = (idxBins[idx] / kNzBins) - 1;
+               Int_t iz = (idxBins[idx] % kNzBins) - 1;
+               if (zmax-zmin>=idz || ymax-ymin>=idy){
+                 if (TMath::Abs(iy-iiy)>0.75*idy) continue;
+                 if (TMath::Abs(iz-iiz)>0.75*idz) continue;
+               }
+               ndigits++;
+               Float_t qBin = bins[idxBins[idx]].GetQ();
+               y += qBin * fYSPD[iy]; 
+               z += qBin * fZSPD[iz]; 
+               q += qBin;   
+             }
+             if (ndigits==0) continue;
+             y /= q; 
+             z /= q;
+             y -= fHwSPD; 
+             z -= fHlSPD;
+
+             Float_t hit[5];  // y, z, sigma(y)^2, sigma(z)^2, charge
+             hit[0] = -(-y+fYshift[iModule]); 
+             if (iModule <= fLastSPD1) hit[0] = -hit[0];
+             hit[1] = -z+fZshift[iModule];
+             hit[2] = fYpitchSPD*fYpitchSPD/12.;
+             hit[3] = fZ1pitchSPD*fZ1pitchSPD/12.;
+             //          hit[4] = q;
+             hit[4] = (zmax-zmin+1)*100 + (ymax-ymin+1);
+             //          CheckLabels(label);
+             Int_t info[3]={ymax-ymin+1,zmax-zmin+1,fNlayer[iModule]};
+             new (clusters[iModule]->AddrAt(nClusters)) 
+               AliITSclusterV2(label, hit,info); 
+             nClusters++;
+           }
        }
 
        nClustersSPD += nClusters;
-       delete bins;
+       bins = NULL;
       }
 
       if (!next) break;
-      bins = new AliBin[kMaxBin];
+      bins = binsSPD;
+      memcpy(binsSPD,binsSPDInit,sizeof(AliBin)*kMaxBin);
     }
 
     // fill the current digit into the bins array
@@ -662,6 +645,9 @@ void AliITSclustererV2::FindClustersSPD(AliITSRawStream* input,
     bins[index].SetQ(1);
   }
 
+  delete [] binsSPDInit;
+  delete [] binsSPD;
+
   Info("FindClustersSPD", "found clusters in ITS SPD: %d", nClustersSPD);
 }
 
@@ -807,12 +793,14 @@ FindClustersSDD(AliBin* bins[2], Int_t nMaxBin, Int_t nzBins,
               if (dj>maxj) maxj=dj;
               if (dj<minj) minj=dj;
               //
-              if (TMath::Abs(di)<2&&TMath::Abs(dj)<2){
-                AliITSdigitSDD* d=(AliITSdigitSDD*)digits->UncheckedAt(b->GetIndex());
-                for (Int_t itrack=0;itrack<10;itrack++){
-                  Int_t track = (d->GetTracks())[itrack];
-                  if (track>=0) {
-                    AddLabel(milab, track); 
+              if(digits) {
+                if (TMath::Abs(di)<2&&TMath::Abs(dj)<2){
+                  AliITSdigitSDD* d=(AliITSdigitSDD*)digits->UncheckedAt(b->GetIndex());
+                  for (Int_t itrack=0;itrack<10;itrack++){
+                    Int_t track = (d->GetTracks())[itrack];
+                    if (track>=0) {
+                      AddLabel(milab, track); 
+                    }
                   }
                 }
               }
@@ -926,6 +914,7 @@ FindClustersSDD(const TClonesArray *digits, TClonesArray *clusters) {
 
   delete[] bins[0];
   delete[] bins[1];
+
 }
 
 void AliITSclustererV2::FindClustersSDD(AliITSRawStream* input, 
@@ -937,7 +926,10 @@ void AliITSclustererV2::FindClustersSDD(AliITSRawStream* input,
   Int_t nClustersSDD = 0;
   Int_t kNzBins = fNzSDD + 2;
   Int_t kMaxBin = kNzBins * (fNySDD+2);
-  AliBin* bins[2] = {NULL, NULL};
+  AliBin *binsSDDInit = new AliBin[kMaxBin];
+  AliBin *binsSDD1 = new AliBin[kMaxBin];
+  AliBin *binsSDD2 = new AliBin[kMaxBin];
+  AliBin *bins[2] = {NULL, NULL};
 
   // read raw data input stream
   while (kTRUE) {
@@ -952,25 +944,32 @@ void AliITSclustererV2::FindClustersSDD(AliITSRawStream* input,
        FindClustersSDD(bins, kMaxBin, kNzBins, NULL, clusters[iModule]);
        Int_t nClusters = clusters[iModule]->GetEntriesFast();
        nClustersSDD += nClusters;
-       delete[] bins[0];
-       delete[] bins[1];
+       bins[0] = bins[1] = NULL;
       }
 
       if (!next) break;
-      bins[0] = new AliBin[kMaxBin];
-      bins[1] = new AliBin[kMaxBin];
+      bins[0]=binsSDD1;
+      bins[1]=binsSDD2;
+      memcpy(binsSDD1,binsSDDInit,sizeof(AliBin)*kMaxBin);
+      memcpy(binsSDD2,binsSDDInit,sizeof(AliBin)*kMaxBin);
     }
 
     // fill the current digit into the bins array
-    Int_t iz = input->GetCoord1()+1;
-    Int_t side = ((iz <= fNzSDD) ? 0 : 1);
-    iz -= side*fNzSDD;
-    Int_t index = (input->GetCoord2()+1) * kNzBins + iz;
-    bins[side][index].SetQ(input->GetSignal());
-    bins[side][index].SetMask(1);
-    bins[side][index].SetIndex(index);
+    if(input->GetSignal()>=3) {
+      Int_t iz = input->GetCoord1()+1;
+      Int_t side = ((iz <= fNzSDD) ? 0 : 1);
+      iz -= side*fNzSDD;
+      Int_t index = (input->GetCoord2()+1) * kNzBins + iz;
+      bins[side][index].SetQ(input->GetSignal());
+      bins[side][index].SetMask(1);
+      bins[side][index].SetIndex(index);
+    }
   }
 
+  delete[] binsSDD1;
+  delete[] binsSDD2;
+  delete[] binsSDDInit;
+
   Info("FindClustersSDD", "found clusters in ITS SDD: %d", nClustersSDD);
 }
 
@@ -998,6 +997,7 @@ FindClustersSSD(Ali1Dcluster* neg, Int_t nn,
   Int_t cused2[3000];
   for (Int_t i=0;i<3000;i++) {cnegative[i]=0; cused1[i]=0;}
   for (Int_t i=0;i<3000;i++) {cpositive[i]=0; cused2[i]=0;}
+  for (Int_t i=0;i<30000;i++){negativepair[i]=positivepair[i]=0;}
   static Short_t pairs[1000][1000];
   memset(pairs,0,sizeof(Short_t)*1000000);
 //   Short_t ** pairs = new Short_t*[1000];
@@ -1546,15 +1546,15 @@ FindClustersSSD(const TClonesArray *alldigits, TClonesArray *clusters) {
   c[*n].SetNd(nd);
   c[*n].SetLabels(lab);
   //Split suspiciously big cluster
-  if (nd>3) {
-     c[*n].SetY(y/q-0.5*nd);
+  if (nd>4 && nd<25) {
+     c[*n].SetY(y/q-0.25*nd);
      c[*n].SetQ(0.5*q);
      (*n)++;
      if (*n==MAX) {
         Error("FindClustersSSD","Too many 1D clusters !");
         return;
      }
-     c[*n].SetY(y/q+0.5*nd);
+     c[*n].SetY(y/q+0.25*nd);
      c[*n].SetQ(0.5*q);
      c[*n].SetNd(nd);
      c[*n].SetLabels(lab);
@@ -1578,20 +1578,23 @@ void AliITSclustererV2::FindClustersSSD(AliITSRawStream* input,
   const Int_t MAX = 1000;
   Ali1Dcluster clusters1D[2][MAX];
   Int_t nClusters[2] = {0, 0};
+  Int_t lab[3]={-2,-2,-2};
   Float_t q = 0.;
   Float_t y = 0.;
   Int_t nDigits = 0;
   Int_t prevStrip = -1;
   Int_t prevFlag = -1;
+  Int_t prevModule = -1;
 
   // read raw data input stream
   while (kTRUE) {
     Bool_t next = input->Next();
 
+    if(input->GetSignal()<3 && next) continue;
     // check if a new cluster starts
     Int_t strip = input->GetCoord2();
     Int_t flag = input->GetCoord1();
-    if ((!next || input->IsNewModule() ||
+    if ((!next || (input->GetModuleID() != prevModule)||
         (strip-prevStrip > 1) || (flag != prevFlag)) &&
        (nDigits > 0)) {
       if (nClusters[prevFlag] == MAX) {
@@ -1602,26 +1605,28 @@ void AliITSclustererV2::FindClustersSSD(AliITSRawStream* input,
       cluster.SetY(y/q);
       cluster.SetQ(q);
       cluster.SetNd(nDigits);
+      cluster.SetLabels(lab);
 
       //Split suspiciously big cluster
-      if (nDigits > 3) {
-       cluster.SetY(y/q - 0.5*nDigits);
+      if (nDigits > 4&&nDigits < 25) {
+       cluster.SetY(y/q - 0.25*nDigits);
         cluster.SetQ(0.5*q);
        if (nClusters[prevFlag] == MAX) {
          Error("FindClustersSSD", "Too many 1D clusters !");
          return;
        }
        Ali1Dcluster& cluster2 = clusters1D[prevFlag][nClusters[prevFlag]++];
-       cluster2.SetY(y/q + 0.5*nDigits);
+       cluster2.SetY(y/q + 0.25*nDigits);
        cluster2.SetQ(0.5*q);
        cluster2.SetNd(nDigits);
+       cluster2.SetLabels(lab);
       }
       y = q = 0.;
       nDigits = 0;
     }
 
-    if (!next || input->IsNewModule()) {
-      Int_t iModule = input->GetPrevModuleID();
+    if (!next || (input->GetModuleID() != prevModule)) {
+      Int_t iModule = prevModule;
 
       // when all data from a module was read, search for clusters
       if (prevFlag >= 0) {
@@ -1645,6 +1650,8 @@ void AliITSclustererV2::FindClustersSSD(AliITSRawStream* input,
     nDigits++;
     prevStrip = strip;
     prevFlag = flag;
+    prevModule = input->GetModuleID();
+
   }
 
   Info("FindClustersSSD", "found clusters in ITS SSD: %d", nClustersSSD);
@@ -1653,7 +1660,7 @@ void AliITSclustererV2::FindClustersSSD(AliITSRawStream* input,
 #else   //V1
 
 #include "AliITSDetType.h"
-
+#include "AliITS.h"
 #include "AliITSsegmentationSPD.h"
 #include "AliITSClusterFinderSPD.h"