ITS/AliITSClusterFinderV2.cxx

   1 /**************************************************************************
   2  * Copyright(c) 1998-2003, ALICE Experiment at CERN, All rights reserved. *
   3  *                                                                        *
   4  * Author: The ALICE Off-line Project.                                    *
   5  * Contributors are mentioned in the code where appropriate.              *
   6  *                                                                        *
   7  * Permission to use, copy, modify and distribute this software and its   *
   8  * documentation strictly for non-commercial purposes is hereby granted   *
   9  * without fee, provided that the above copyright notice appears in all   *
  10  * copies and that both the copyright notice and this permission notice   *
  11  * appear in the supporting documentation. The authors make no claims     *
  12  * about the suitability of this software for any purpose. It is          *
  13  * provided "as is" without express or implied warranty.                  *
  14  **************************************************************************/
  15 ////////////////////////////////////////////////////////////////////////////
  16 //            Implementation of the ITS clusterer V2 class                //
  17 //                                                                        //
  18 //          Origin: Iouri Belikov, CERN, Jouri.Belikov@cern.ch            //
  19 //                                                                        //
  20 ///////////////////////////////////////////////////////////////////////////
  21
  22
  23 #include "AliRun.h"
  24 #include "AliITSClusterFinderV2.h"
  25 #include "AliITSRecPoint.h"
  26 #include "AliITSgeomTGeo.h"
  27 #include "AliITSDetTypeRec.h"
  28 #include <TParticle.h>
  29 #include "AliMC.h"
  30
  31 ClassImp(AliITSClusterFinderV2)
  32
  33 extern AliRun *gAlice;
  34
  35 AliITSClusterFinderV2::AliITSClusterFinderV2(AliITSDetTypeRec* dettyp):AliITSClusterFinder(dettyp),
  36 fNModules(AliITSgeomTGeo::GetNModules()),
  37 fEvent(0){
  38   //Default constructor
  39 }
  40 /*
  41 //______________________________________________________________________
  42 AliITSClusterFinderV2::AliITSClusterFinderV2(const AliITSClusterFinderV2 &source) : AliITSClusterFinder(source) {
  43   // Copy constructor
  44   // Copies are not allowed. The method is protected to avoid misuse.
  45   Fatal("AliITSClusterFinderV2","Copy constructor not allowed\n");
  46 }
  47 */
  48 //______________________________________________________________________
  49 //AliITSClusterFinderV2& AliITSClusterFinderV2::operator=(const AliITSClusterFinderV2& /* source */){
  50   // Assignment operator
  51   // Assignment is not allowed. The method is protected to avoid misuse.
  52   //Fatal("= operator","Assignment operator not allowed\n");
  53   //return *this;
  54 //}
  55
  56
  57 //______________________________________________________________________
  58 void AliITSClusterFinderV2::CheckLabels2(Int_t lab[10]) {
  59   //------------------------------------------------------------
  60   // Tries to find mother's labels
  61   //------------------------------------------------------------
  62   Int_t nlabels =0;
  63   for (Int_t i=0;i<10;i++) if (lab[i]>=0) nlabels++;
  64   if(nlabels == 0) return; // In case of no labels just exit
  65
  66
  67   Int_t ntracks = gAlice->GetMCApp()->GetNtrack();
  68
  69   for (Int_t i=0;i<10;i++){
  70     Int_t label = lab[i];
  71     if (label>=0 && label<ntracks) {
  72       TParticle *part=(TParticle*)gAlice->GetMCApp()->Particle(label);
  73
  74       if (part->P() < 0.02) {
  75           Int_t m=part->GetFirstMother();
  76           if (m<0) {
  77             continue;
  78           }
  79           if (part->GetStatusCode()>0) {
  80             continue;
  81           }
  82           lab[i]=m;
  83       }
  84       else
  85         if (part->P() < 0.12 && nlabels>3) {
  86           lab[i]=-2;
  87           nlabels--;
  88         }
  89     }
  90     else{
  91       if ( (label>ntracks||label <0) && nlabels>3) {
  92         lab[i]=-2;
  93         nlabels--;
  94       }
  95     }
  96   }
  97   if (nlabels>3){
  98     for (Int_t i=0;i<10;i++){
  99       if (nlabels>3){
 100         Int_t label = lab[i];
 101         if (label>=0 && label<ntracks) {
 102           TParticle *part=(TParticle*)gAlice->GetMCApp()->Particle(label);
 103           if (part->P() < 0.1) {
 104             lab[i]=-2;
 105             nlabels--;
 106           }
 107         }
 108       }
 109     }
 110   }
 111
 112   //compress labels -- if multi-times the same
 113   Int_t lab2[10];
 114   for (Int_t i=0;i<10;i++) lab2[i]=-2;
 115   for (Int_t i=0;i<10  ;i++){
 116     if (lab[i]<0) continue;
 117     for (Int_t j=0;j<10 &&lab2[j]!=lab[i];j++){
 118       if (lab2[j]<0) {
 119         lab2[j]= lab[i];
 120         break;
 121       }
 122     }
 123   }
 124   for (Int_t j=0;j<10;j++) lab[j]=lab2[j];
 125
 126 }
 127
 128 //______________________________________________________________________
 129 void AliITSClusterFinderV2::AddLabel(Int_t lab[10], Int_t label) {
 130
 131   //add label to the cluster
 132
 133   if(label<0) return; // In case of no label just exit
 134
 135   Int_t ntracks = gAlice->GetMCApp()->GetNtrack();
 136   if (label>ntracks) return;
 137   for (Int_t i=0;i<10;i++){
 138     //    if (label<0) break;
 139     if (lab[i]==label) break;
 140     if (lab[i]<0) {
 141       lab[i]= label;
 142       break;
 143     }
 144   }
 145 }
 146
 147
 148 //______________________________________________________________________
 149 void AliITSClusterFinderV2::
 150 FindCluster(Int_t k,Int_t maxz,AliBin *bins,Int_t &n,Int_t *idx) {
 151   //------------------------------------------------------------
 152   // returns an array of indices of digits belonging to the cluster
 153   // (needed when the segmentation is not regular)
 154   //------------------------------------------------------------
 155   if (n<200) idx[n++]=bins[k].GetIndex();
 156   bins[k].Use();
 157
 158   if (bins[k-maxz].IsNotUsed()) FindCluster(k-maxz,maxz,bins,n,idx);
 159   if (bins[k-1   ].IsNotUsed()) FindCluster(k-1   ,maxz,bins,n,idx);
 160   if (bins[k+maxz].IsNotUsed()) FindCluster(k+maxz,maxz,bins,n,idx);
 161   if (bins[k+1   ].IsNotUsed()) FindCluster(k+1   ,maxz,bins,n,idx);
 162   /*
 163   if (bins[k-maxz-1].IsNotUsed()) FindCluster(k-maxz-1,maxz,bins,n,idx);
 164   if (bins[k-maxz+1].IsNotUsed()) FindCluster(k-maxz+1,maxz,bins,n,idx);
 165   if (bins[k+maxz-1].IsNotUsed()) FindCluster(k+maxz-1,maxz,bins,n,idx);
 166   if (bins[k+maxz+1].IsNotUsed()) FindCluster(k+maxz+1,maxz,bins,n,idx);
 167   */
 168 }
 169
 170 //______________________________________________________________________
 171 Bool_t AliITSClusterFinderV2::IsMaximum(Int_t k,Int_t max,const AliBin *bins) {
 172   //------------------------------------------------------------
 173   //is this a local maximum ?
 174   //------------------------------------------------------------
 175   UShort_t q=bins[k].GetQ();
 176   if (q==1023) return kFALSE;
 177   if (bins[k-max].GetQ() > q) return kFALSE;
 178   if (bins[k-1  ].GetQ() > q) return kFALSE;
 179   if (bins[k+max].GetQ() > q) return kFALSE;
 180   if (bins[k+1  ].GetQ() > q) return kFALSE;
 181   if (bins[k-max-1].GetQ() > q) return kFALSE;
 182   if (bins[k+max-1].GetQ() > q) return kFALSE;
 183   if (bins[k+max+1].GetQ() > q) return kFALSE;
 184   if (bins[k-max+1].GetQ() > q) return kFALSE;
 185   return kTRUE;
 186 }
 187
 188 //______________________________________________________________________
 189 void AliITSClusterFinderV2::
 190 FindPeaks(Int_t k,Int_t max,AliBin *b,Int_t *idx,UInt_t *msk,Int_t& n) {
 191   //------------------------------------------------------------
 192   //find local maxima
 193   //------------------------------------------------------------
 194   if (n<31)
 195   if (IsMaximum(k,max,b)) {
 196     idx[n]=k; msk[n]=(2<<n);
 197     n++;
 198   }
 199   b[k].SetMask(0);
 200   if (b[k-max].GetMask()&1) FindPeaks(k-max,max,b,idx,msk,n);
 201   if (b[k-1  ].GetMask()&1) FindPeaks(k-1  ,max,b,idx,msk,n);
 202   if (b[k+max].GetMask()&1) FindPeaks(k+max,max,b,idx,msk,n);
 203   if (b[k+1  ].GetMask()&1) FindPeaks(k+1  ,max,b,idx,msk,n);
 204 }
 205
 206 //______________________________________________________________________
 207 void AliITSClusterFinderV2::
 208 MarkPeak(Int_t k, Int_t max, AliBin *bins, UInt_t m) {
 209   //------------------------------------------------------------
 210   //mark this peak
 211   //------------------------------------------------------------
 212   UShort_t q=bins[k].GetQ();
 213
 214   bins[k].SetMask(bins[k].GetMask()|m);
 215
 216   if (bins[k-max].GetQ() <= q)
 217      if ((bins[k-max].GetMask()&m) == 0) MarkPeak(k-max,max,bins,m);
 218   if (bins[k-1  ].GetQ() <= q)
 219      if ((bins[k-1  ].GetMask()&m) == 0) MarkPeak(k-1  ,max,bins,m);
 220   if (bins[k+max].GetQ() <= q)
 221      if ((bins[k+max].GetMask()&m) == 0) MarkPeak(k+max,max,bins,m);
 222   if (bins[k+1  ].GetQ() <= q)
 223      if ((bins[k+1  ].GetMask()&m) == 0) MarkPeak(k+1  ,max,bins,m);
 224 }
 225
 226 //______________________________________________________________________
 227 void AliITSClusterFinderV2::
 228 MakeCluster(Int_t k,Int_t max,AliBin *bins,UInt_t m,AliITSRecPoint &c) {
 229   //------------------------------------------------------------
 230   //make cluster using digits of this peak
 231   //------------------------------------------------------------
 232   Float_t q=(Float_t)bins[k].GetQ();
 233   Int_t i=k/max, j=k-i*max;
 234
 235   c.SetQ(c.GetQ()+q);
 236   c.SetY(c.GetY()+i*q);
 237   c.SetZ(c.GetZ()+j*q);
 238   c.SetSigmaY2(c.GetSigmaY2()+i*i*q);
 239   c.SetSigmaZ2(c.GetSigmaZ2()+j*j*q);
 240
 241   bins[k].SetMask(0xFFFFFFFE);
 242
 243   if (bins[k-max].GetMask() == m) MakeCluster(k-max,max,bins,m,c);
 244   if (bins[k-1  ].GetMask() == m) MakeCluster(k-1  ,max,bins,m,c);
 245   if (bins[k+max].GetMask() == m) MakeCluster(k+max,max,bins,m,c);
 246   if (bins[k+1  ].GetMask() == m) MakeCluster(k+1  ,max,bins,m,c);
 247 }