----------------------------------------------------------------------
[u/mrichter/AliRoot.git] / PHOS / AliPHOSClusterizerv1.cxx
1 /**************************************************************************
2  * Copyright(c) 1998-1999, 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 //_________________________________________________________________________
17 // A brief description of the class
18 //*-- Author : Yves Schutz  SUBATECH 
19 //////////////////////////////////////////////////////////////////////////////
20
21 // --- ROOT system ---
22
23 #include "TMath.h" 
24
25 // --- Standard library ---
26
27 #include "iostream.h"
28
29 // --- AliRoot header files ---
30
31 #include "AliPHOSClusterizerv1.h"
32 #include "AliPHOSDigit.h"
33 #include "AliPHOSEmcRecPoint.h"
34 #include "AliPHOSPpsdRecPoint.h"
35 #include "AliPHOSv0.h" 
36 #include "AliRun.h" 
37
38 ClassImp(AliPHOSClusterizerv1)
39
40 //____________________________________________________________________________
41 AliPHOSClusterizerv1::AliPHOSClusterizerv1()
42 {
43   fA                       = 0.;
44   fB                       = 0.01 ;
45   fNumberOfEmcClusters     = 0 ; 
46   fNumberOfPpsdClusters    = 0 ; 
47   fEmcClusteringThreshold  = 0.1;   
48   fEmcEnergyThreshold      = 0.01;    
49   fPpsdClusteringThreshold = 0.00000015; 
50   fPpsdEnergyThreshold     = 0.0000001;  
51   fW0                      = 4.5 ;
52   fLocMaxCut               = 0.06 ;
53 }
54
55 //____________________________________________________________________________
56 Int_t AliPHOSClusterizerv1::AreNeighbours(AliPHOSDigit * d1, AliPHOSDigit * d2)
57 {
58   // neigbours are defined as digits having at least common virtix
59   // The order of A and B in AreNeighbours(A,B) is important: first (A) should be digit 
60   // in cluster, which compared with digits, which not clasterized yet  
61   Int_t rv = 0 ; 
62
63   AliPHOSGeometry * geom = AliPHOSGeometry::GetInstance() ;
64
65   Int_t relid1[4] ; 
66   geom->AbsToRelNumbering(d1->GetId(), relid1) ; 
67
68   Int_t relid2[4] ; 
69   geom->AbsToRelNumbering(d2->GetId(), relid2) ; 
70  
71   if ( (relid1[0] == relid2[0]) && (relid1[1]==relid2[1]) ) { // inside the same PHOS module and the same PPSD Module 
72     Int_t RowDiff = TMath::Abs( relid1[2] - relid2[2] ) ;  
73     Int_t ColDiff = TMath::Abs( relid1[3] - relid2[3] ) ;  
74     
75     if (( ColDiff<=1 )  && ( RowDiff <= 1 )){
76       rv = 1 ; 
77     }
78     else {
79       if((relid2[2] > relid1[2]) && (relid2[3] > relid1[3]+1)) 
80         rv = 2; //  Difference in row numbers is too large to look further 
81     }
82
83   } 
84   else {
85     
86     if( (relid1[0] < relid2[0]) || (relid1[1] < relid2[1]) )  
87       rv=2 ;
88
89   }
90   
91   return rv ; 
92 }
93
94 //____________________________________________________________________________
95 void AliPHOSClusterizerv1::FillandSort(const DigitsList * dl, TObjArray * tl) 
96 {
97   // copies the digits with energy above thershold and sorts the list
98   // according to increasing Id number
99
100
101   AliPHOSGeometry * geom = AliPHOSGeometry::GetInstance() ;
102   Int_t relid[4] ;  
103   
104   TIter next(dl) ; 
105   AliPHOSDigit * digit ;
106   
107   while ( (digit = (AliPHOSDigit *)next()) ) { 
108     
109     Int_t id    = digit->GetId() ; 
110     Float_t ene = Calibrate(digit->GetAmp()) ; 
111     geom->AbsToRelNumbering(id, relid) ;
112     
113     if(relid[1]==0){ // EMC
114       if ( ene > fEmcEnergyThreshold ) 
115         tl->Add(digit) ;
116     }
117
118     else { //Ppsd
119       if ( ene > fPpsdEnergyThreshold )
120         tl->Add(digit) ;        
121     }
122
123   }
124   tl->Sort() ; 
125 }
126
127 //____________________________________________________________________________
128 void AliPHOSClusterizerv1:: GetNumberOfClustersFound(Int_t * numb) 
129 {
130  numb[0] = fNumberOfEmcClusters ; 
131  numb[1] = fNumberOfPpsdClusters ; 
132 }
133
134 //____________________________________________________________________________
135 Bool_t AliPHOSClusterizerv1::IsInEmc(AliPHOSDigit * digit) 
136 {
137    Bool_t rv = kFALSE ; 
138
139   AliPHOSGeometry * geom = AliPHOSGeometry::GetInstance() ;  
140
141   Int_t relid[4] ; 
142   geom->AbsToRelNumbering(digit->GetId(), relid) ; 
143
144   if ( relid[1] == 0  )
145     rv = kTRUE; 
146
147   return rv ; 
148 }
149
150 //____________________________________________________________________________
151 void AliPHOSClusterizerv1::MakeClusters(const DigitsList * dl, RecPointsList * emcl, RecPointsList * ppsdl)
152 {
153
154   // Fill and sort the working digits list
155   TObjArray TempoDigitsList( dl->GetEntries() ) ;
156   this->FillandSort(dl, &TempoDigitsList) ; 
157
158
159
160   // Clusterization starts  
161   TIter nextdigit(&TempoDigitsList) ; 
162   AliPHOSDigit * digit ; 
163   Bool_t NotRemoved = kTRUE ;
164
165   while ( (digit = (AliPHOSDigit *)nextdigit()) ) { // scan over the list of digits
166     AliPHOSRecPoint * clu ; 
167
168     int * ClusterDigitsList[dl->GetEntries()] ;   
169     Int_t index ;
170     if (( (this->IsInEmc(digit)) && (Calibrate(digit->GetAmp()) > fEmcClusteringThreshold ) ) || 
171         ( (!this->IsInEmc(digit)) && (Calibrate(digit->GetAmp()) > fPpsdClusteringThreshold ) ) ) {
172   
173       Int_t iDigitInCluster = 0 ; 
174
175       if  (this->IsInEmc(digit) ) {   
176         new ((*emcl)[fNumberOfEmcClusters]) AliPHOSEmcRecPoint(fW0, fLocMaxCut) ;// start a new EMC RecPoint
177
178         clu = (AliPHOSEmcRecPoint *) (*emcl)[fNumberOfEmcClusters++] ; 
179  
180         clu->AddDigit(*digit,Calibrate(digit->GetAmp())) ; 
181
182         ClusterDigitsList[iDigitInCluster++] = (int* ) digit ;  
183         TempoDigitsList.Remove(digit) ; 
184       }
185
186       else { 
187         new ((*ppsdl)[fNumberOfPpsdClusters]) AliPHOSPpsdRecPoint() ;// start a new PPSD cluster
188         clu =  (AliPHOSPpsdRecPoint *) ppsdl->At(fNumberOfPpsdClusters++)  ;  
189         clu->AddDigit(*digit,0.) ;      
190         ClusterDigitsList[iDigitInCluster++] = (int* ) digit ;  
191         TempoDigitsList.Remove(digit) ; 
192         nextdigit.Reset() ;
193         
194         //Here we remove resting EMC digits, which can not make cluster
195         if(NotRemoved){ 
196           
197           while( (digit = (AliPHOSDigit *)nextdigit()) ){
198             
199             if(IsInEmc(digit)) TempoDigitsList.Remove(digit) ;
200             else 
201               break ;
202             
203           }// while digit  
204           
205         } // if NotRemoved 
206         
207       } // else        
208       
209       nextdigit.Reset() ;
210
211       AliPHOSDigit * digitN ; 
212       index = 0 ;
213       while (index < iDigitInCluster){ // scan over digits already in claster 
214         digit = (AliPHOSDigit *) ClusterDigitsList[index++]  ;      
215
216         while ( (digitN = (AliPHOSDigit *)nextdigit()) ) { // scan over the reduced list of digits 
217           Int_t ineb = AreNeighbours(digit, digitN);   //  call (digit,digitN) in THAT oder !!!!!
218           switch (ineb ) {
219           case 0 :   // not a neibors
220             break ;      
221           case 1 :    // Are neibors 
222             clu->AddDigit(*digitN,Calibrate(digitN->GetAmp())) ;
223             ClusterDigitsList[iDigitInCluster++] =(int*) digitN ; 
224             TempoDigitsList.Remove(digitN) ;
225             break ;
226           case 2 :  // to far from each other
227             goto endofloop;   
228           } // switch
229           
230         } // while digitN
231
232       endofloop: ;
233         nextdigit.Reset() ; 
234         
235       } // loop over cluster     
236  
237    }  //below energy theshold  
238   
239   } // while digit
240
241   TempoDigitsList.Clear() ; 
242 }
243
244 //____________________________________________________________________________
245 void AliPHOSClusterizerv1::PrintParameters() 
246 {
247   cout << "PHOS Clusterizer version 1 :" << endl 
248        << "                       EMC  Clustering threshold = " << fEmcClusteringThreshold << endl
249        << "                       EMC  Energy threshold     = " << fEmcEnergyThreshold << endl                  
250        << "                      PPSD  Clustering threshold = " << fPpsdClusteringThreshold << endl
251        << "                      PPSD  Energy threshold     = " << fPpsdEnergyThreshold << endl ;                
252 }