ROOT logo
#ifndef ALIITSCLUSTERFINDER_H
#define ALIITSCLUSTERFINDER_H
 
/* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
 * See cxx source for full Copyright notice                               */
 
/* $Id$ */

////////////////////////////////////////////////
//  ITS Cluster Finder Class                  //
//                                            //
//                                            //
////////////////////////////////////////////////

#include <TObject.h>
#include <TClonesArray.h>
#include "AliLog.h"

class AliITSMap;
class AliITSresponse;
class AliITSsegmentation;
class AliITSdigit;
class AliITSRecPoint;
class AliITSDetTypeRec;
class AliRawReader;
class TArrayI;

using std::istream;

//----------------------------------------------------------------------
class AliITSClusterFinder :public TObject{
  public:
    AliITSClusterFinder(); // Default constructor
    // Standard Constructor
    AliITSClusterFinder(AliITSDetTypeRec* dettyp);
    AliITSClusterFinder(AliITSDetTypeRec* dettyp,TClonesArray *digits);// Standard+ Constructor
    virtual ~AliITSClusterFinder(); // Destructor
    //
    // Do the Reconstruction.
    virtual void FindRawClusters(Int_t /*mod*/)=0; // Finds cluster of digits.
    virtual void RawdataToClusters(AliRawReader* /*rawReader*/) {
      AliError("Method not implemented in this class ");
    }

    // Digit
    virtual void SetDigits(TClonesArray *itsDigits) {// set digits
        fDigits=itsDigits;fNdigits = fDigits->GetEntriesFast();}
    virtual AliITSdigit* GetDigit(Int_t i){ // Returns ith digit
        return (AliITSdigit*) fDigits->UncheckedAt(i);}
    virtual TClonesArray* Digits(){return fDigits;}// Gets fDigits
    virtual Int_t   NDigits() const {return fNdigits;}// Get Number of Digits

    // Set fClusters up
    virtual void SetClusters(TClonesArray *itsClusters){// set clusters
        fClusters = itsClusters;}
    // Get fCluters
    virtual TClonesArray* Clusters(){return fClusters;}
    // Returns the present number of enteries
    virtual Int_t NClusters()const {return fClusters->GetEntriesFast();}

    virtual void SetModule(Int_t module){fModule = module;}// Set module number
    virtual Int_t GetModule()const {return fModule;}// Returns module number

    void SetEvent(Int_t event) { fEvent=event; }

    // Others
    virtual void  SetMap(AliITSMap *m) {fMap=m;}// map
    AliITSMap* Map(){return fMap;}// map
    virtual Int_t GetNPeaks() const {return fNPeaks;}// returns fNPeaks
    //
    virtual Bool_t IsNeighbor(TObjArray *digs,Int_t i,Int_t j[]) const;
        // Set max. cluster size ; bigger clusters will be rejected

    // IO functions
    void Print(ostream *os) const; // Class ascii print function
    void Read(istream *os);  // Class ascii read function
    virtual void Print(Option_t *option="") const {TObject::Print(option);}
    virtual Int_t Read(const char *name) {return TObject::Read(name);}

    virtual void SetDetTypeRec(AliITSDetTypeRec* dtr) {fDetTypeRec=dtr;}
    AliITSDetTypeRec* GetDetTypeRec() const {return fDetTypeRec;}

    void InitGeometry(); 
    //
    Int_t    GetNClusters()                               const {return fNClusters;}
    void     SetRawID2ClusID(TArrayI *arr)                      {fRawID2ClusID = arr;} 
    TArrayI* GetRawID2ClusID()                            const {return fRawID2ClusID;} 
    // 
  protected:
  class Ali1Dcluster {
  public:
    void SetY(Float_t y) {fY=y;}
    void SetQ(Float_t q) {fQ=q;}
    void SetNd(Int_t n)  {fNd=n;}
    void SetLabels(Int_t *lab) {fLab[0]=lab[0];fLab[1]=lab[1];fLab[2]=lab[2];}
    Float_t GetY() const {return fY;}
    Float_t GetQ() const {return fQ;}
    Int_t GetNd()const {return fNd;}
    Int_t GetLabel(Int_t lab) const { return fLab[lab]; }
  protected:
    Float_t fY; //cluster position
    Float_t fQ; //cluster charge
    Int_t fNd;  //number of digits
    Int_t fLab[3]; //track label
  };
  class AliBin {
  public:
  AliBin():fIndex(0),fMask(0xFFFFFFFE),fRawID(-1),fQ(0){}
    void SetIndex(UInt_t idx) {fIndex=idx;}
    void SetQ(UShort_t q)  {fQ=q;}
    void SetMask(UInt_t m) {fMask=m;}
    void SetRawID(Int_t id) {fRawID=id;}
    void Reset() {fIndex=0; fMask=0xFFFFFFFE; fQ=0; fRawID=-1;}

    void Use() {fMask&=0xFFFFFFFE;}
    Bool_t IsNotUsed() const {return (fMask&1);}
    Bool_t IsUsed() const {return !(IsNotUsed());}

    UInt_t   GetIndex() const {return fIndex;}
    UShort_t GetQ()     const {return fQ;}
    UInt_t   GetMask()  const {return fMask;}
    Int_t    GetRawID() const {return fRawID;}
  protected:
    UInt_t fIndex; //digit index
    UInt_t fMask; //peak mask
    Int_t  fRawID; // ID of raw word (used for embedding)
    UShort_t fQ;  //signal
  };
  void MakeCluster(Int_t k,Int_t max,AliBin *bins,UInt_t m,AliITSRecPoint &c);
  static Bool_t IsMaximum(Int_t k, Int_t max, const AliBin *bins);
  static void FindPeaks(Int_t k,Int_t m,AliBin*b,Int_t*idx,UInt_t*msk,Int_t&n);
  static void MarkPeak(Int_t k, Int_t max, AliBin *bins, UInt_t m);
  static void FindCluster(Int_t k,Int_t maxz,AliBin *bins,Int_t &n,Int_t *idx);

  static void CheckLabels2(Int_t lab[10]);
  static void AddLabel(Int_t lab[10], Int_t label);      

  // data members      

  Int_t              fModule;        //! Module number to be reconstuctted
  TClonesArray       *fDigits;       //! digits 
  Int_t              fNdigits;       //! num of digits 
  
  AliITSDetTypeRec* fDetTypeRec; //ITS object for reconstruction
  TClonesArray       *fClusters;     //! Array of clusters
  AliITSMap          *fMap;          //! map
  Int_t              fNPeaks;        //! NPeaks  
  // Data members needed to fill AliCluster objects
  Int_t fNdet[2200];           // detector index  
  Int_t fNlayer[2200];         // detector layer
  
  Int_t fNModules;             // total number of modules    
  Int_t fEvent;                //event number
  Int_t fZmin;   // minimum channel in Zloc
  Int_t fZmax;   // maximum channel in Zloc
  Int_t fXmin;   // minimum channel in Xloc
  Int_t fXmax;   // maximum channel in Xloc 
  //
  UInt_t fNClusters; // total number of clusters found
  //
  TArrayI* fRawID2ClusID;        //! optional array to store raw word ID -> ClusID for embedding (not owned)
  
  AliITSClusterFinder(const AliITSClusterFinder &source); // copy constructor
  // assignment operator
  AliITSClusterFinder& operator=(const AliITSClusterFinder &source);
  

  ClassDef(AliITSClusterFinder,11) //Class for clustering and reconstruction of space points
};
// Input and output functions for standard C++ input/output.
ostream &operator<<(ostream &os,AliITSClusterFinder &source);
istream &operator>>(istream &os,AliITSClusterFinder &source);
#endif
 AliITSClusterFinder.h:1
 AliITSClusterFinder.h:2
 AliITSClusterFinder.h:3
 AliITSClusterFinder.h:4
 AliITSClusterFinder.h:5
 AliITSClusterFinder.h:6
 AliITSClusterFinder.h:7
 AliITSClusterFinder.h:8
 AliITSClusterFinder.h:9
 AliITSClusterFinder.h:10
 AliITSClusterFinder.h:11
 AliITSClusterFinder.h:12
 AliITSClusterFinder.h:13
 AliITSClusterFinder.h:14
 AliITSClusterFinder.h:15
 AliITSClusterFinder.h:16
 AliITSClusterFinder.h:17
 AliITSClusterFinder.h:18
 AliITSClusterFinder.h:19
 AliITSClusterFinder.h:20
 AliITSClusterFinder.h:21
 AliITSClusterFinder.h:22
 AliITSClusterFinder.h:23
 AliITSClusterFinder.h:24
 AliITSClusterFinder.h:25
 AliITSClusterFinder.h:26
 AliITSClusterFinder.h:27
 AliITSClusterFinder.h:28
 AliITSClusterFinder.h:29
 AliITSClusterFinder.h:30
 AliITSClusterFinder.h:31
 AliITSClusterFinder.h:32
 AliITSClusterFinder.h:33
 AliITSClusterFinder.h:34
 AliITSClusterFinder.h:35
 AliITSClusterFinder.h:36
 AliITSClusterFinder.h:37
 AliITSClusterFinder.h:38
 AliITSClusterFinder.h:39
 AliITSClusterFinder.h:40
 AliITSClusterFinder.h:41
 AliITSClusterFinder.h:42
 AliITSClusterFinder.h:43
 AliITSClusterFinder.h:44
 AliITSClusterFinder.h:45
 AliITSClusterFinder.h:46
 AliITSClusterFinder.h:47
 AliITSClusterFinder.h:48
 AliITSClusterFinder.h:49
 AliITSClusterFinder.h:50
 AliITSClusterFinder.h:51
 AliITSClusterFinder.h:52
 AliITSClusterFinder.h:53
 AliITSClusterFinder.h:54
 AliITSClusterFinder.h:55
 AliITSClusterFinder.h:56
 AliITSClusterFinder.h:57
 AliITSClusterFinder.h:58
 AliITSClusterFinder.h:59
 AliITSClusterFinder.h:60
 AliITSClusterFinder.h:61
 AliITSClusterFinder.h:62
 AliITSClusterFinder.h:63
 AliITSClusterFinder.h:64
 AliITSClusterFinder.h:65
 AliITSClusterFinder.h:66
 AliITSClusterFinder.h:67
 AliITSClusterFinder.h:68
 AliITSClusterFinder.h:69
 AliITSClusterFinder.h:70
 AliITSClusterFinder.h:71
 AliITSClusterFinder.h:72
 AliITSClusterFinder.h:73
 AliITSClusterFinder.h:74
 AliITSClusterFinder.h:75
 AliITSClusterFinder.h:76
 AliITSClusterFinder.h:77
 AliITSClusterFinder.h:78
 AliITSClusterFinder.h:79
 AliITSClusterFinder.h:80
 AliITSClusterFinder.h:81
 AliITSClusterFinder.h:82
 AliITSClusterFinder.h:83
 AliITSClusterFinder.h:84
 AliITSClusterFinder.h:85
 AliITSClusterFinder.h:86
 AliITSClusterFinder.h:87
 AliITSClusterFinder.h:88
 AliITSClusterFinder.h:89
 AliITSClusterFinder.h:90
 AliITSClusterFinder.h:91
 AliITSClusterFinder.h:92
 AliITSClusterFinder.h:93
 AliITSClusterFinder.h:94
 AliITSClusterFinder.h:95
 AliITSClusterFinder.h:96
 AliITSClusterFinder.h:97
 AliITSClusterFinder.h:98
 AliITSClusterFinder.h:99
 AliITSClusterFinder.h:100
 AliITSClusterFinder.h:101
 AliITSClusterFinder.h:102
 AliITSClusterFinder.h:103
 AliITSClusterFinder.h:104
 AliITSClusterFinder.h:105
 AliITSClusterFinder.h:106
 AliITSClusterFinder.h:107
 AliITSClusterFinder.h:108
 AliITSClusterFinder.h:109
 AliITSClusterFinder.h:110
 AliITSClusterFinder.h:111
 AliITSClusterFinder.h:112
 AliITSClusterFinder.h:113
 AliITSClusterFinder.h:114
 AliITSClusterFinder.h:115
 AliITSClusterFinder.h:116
 AliITSClusterFinder.h:117
 AliITSClusterFinder.h:118
 AliITSClusterFinder.h:119
 AliITSClusterFinder.h:120
 AliITSClusterFinder.h:121
 AliITSClusterFinder.h:122
 AliITSClusterFinder.h:123
 AliITSClusterFinder.h:124
 AliITSClusterFinder.h:125
 AliITSClusterFinder.h:126
 AliITSClusterFinder.h:127
 AliITSClusterFinder.h:128
 AliITSClusterFinder.h:129
 AliITSClusterFinder.h:130
 AliITSClusterFinder.h:131
 AliITSClusterFinder.h:132
 AliITSClusterFinder.h:133
 AliITSClusterFinder.h:134
 AliITSClusterFinder.h:135
 AliITSClusterFinder.h:136
 AliITSClusterFinder.h:137
 AliITSClusterFinder.h:138
 AliITSClusterFinder.h:139
 AliITSClusterFinder.h:140
 AliITSClusterFinder.h:141
 AliITSClusterFinder.h:142
 AliITSClusterFinder.h:143
 AliITSClusterFinder.h:144
 AliITSClusterFinder.h:145
 AliITSClusterFinder.h:146
 AliITSClusterFinder.h:147
 AliITSClusterFinder.h:148
 AliITSClusterFinder.h:149
 AliITSClusterFinder.h:150
 AliITSClusterFinder.h:151
 AliITSClusterFinder.h:152
 AliITSClusterFinder.h:153
 AliITSClusterFinder.h:154
 AliITSClusterFinder.h:155
 AliITSClusterFinder.h:156
 AliITSClusterFinder.h:157
 AliITSClusterFinder.h:158
 AliITSClusterFinder.h:159
 AliITSClusterFinder.h:160
 AliITSClusterFinder.h:161
 AliITSClusterFinder.h:162
 AliITSClusterFinder.h:163
 AliITSClusterFinder.h:164
 AliITSClusterFinder.h:165
 AliITSClusterFinder.h:166
 AliITSClusterFinder.h:167
 AliITSClusterFinder.h:168
 AliITSClusterFinder.h:169
 AliITSClusterFinder.h:170
 AliITSClusterFinder.h:171
 AliITSClusterFinder.h:172
 AliITSClusterFinder.h:173