Coding convention obeyed
[u/mrichter/AliRoot.git] / PHOS / AliPHOSJet.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 /* $Id$ */
17
18 //_________________________________________________________________________
19 // Class to calculate jet chararacteristics
20 //
21 // This class implements for PHOS a jet finder for PHOS. It depends on a 
22 // energy seed
23 // minimum energy, cone radius and movement of the cone.
24 //*-- Author :  D.Peressounko
25 //////////////////////////////////////////////////////////////////////////////
26
27 // --- ROOT system ---
28 #include "TParticle.h"
29
30 // --- Standard library ---
31
32 // --- AliRoot header files ---
33 #include "AliPHOSJet.h"
34 //  #include "AliPHOSGetter.h"
35
36 ClassImp(AliPHOSJet)
37   
38 //____________________________________________________________________________ 
39 AliPHOSJet::AliPHOSJet():TObject() {
40   //Inicilize members
41   fNpart = 0 ;
42   fList = 0 ;
43   //  fMode = 0 ;
44   fConeRad = 0 ;
45   fMaxConeMove = 0 ;
46   fMinConeMove = 0 ;
47   
48   fSumEnergy = 0 ;
49   fSumEta = 0 ;
50   fSumPhi = 0 ; 
51   fEnergy = 0 ;
52   fEta = 0 ;
53   fPhi = 0 ; 
54   fLEnergy = 0 ;
55   fLEta = 0 ;
56   fLPhi = 0 ; 
57 }
58
59 //____________________________________________________________________________ 
60 AliPHOSJet::~AliPHOSJet(){
61   //dtor
62   if(fList){
63     delete fList ;
64     fList = 0 ;
65   }
66   
67 }
68 //____________________________________________________________________________ 
69 void AliPHOSJet::AddParticle(const TParticle * p,const Int_t index){
70   //adds particle to jet. Calculates change in jet direction, 
71   //due to addition of this particle and if it is smaller, than fMaxDev, 
72   //add particle, axcept new direction and return true.
73   
74   fSumEnergy+=p->Energy() ;
75   if(p->Pt()/p->Energy() > 0.001)
76     fSumEta+=p->Eta()*p->Energy() ;
77   else
78     fSumEta+=100.*p->Pz() ; 
79   fSumPhi+=p->Phi()*p->Energy() ;    
80   
81   //check if this a leading particle?
82   if(fLEnergy < p->Energy()){
83     fLEnergy = p->Energy() ;
84     if(p->Pt()/p->Energy() > 0.001)
85       fLEta = p->Eta() ;
86     else 
87       fLEta = 100.*p->Pz()/p->Energy() ;
88     fLPhi = p->Phi() ;
89   }
90
91   if(index >=0){ //add index to list of indexes
92     if(!fList){
93       fList = new TArrayI(100) ;
94     }
95     if(fList->GetSize()<=fNpart+1){
96       TArrayI * tmp = new TArrayI(fNpart*2) ;
97       tmp->Adopt(fList->GetSize(),fList->GetArray()) ; //note, we should not delete old array!
98       fList=tmp ;
99     }
100     fList->AddAt(index,fNpart++) ;
101   }
102 }
103 //____________________________________________________________________________ 
104 void AliPHOSJet::AddDigit(const Double_t e,const Double_t eta,const Double_t phi, const Int_t index){
105   //adds particle to jet. Calculates change in jet direction, 
106   //due to addition of this particle and if it is smaller, than fMaxDev, 
107   //add particle, axcept new direction and return true.
108   
109   fSumEnergy+=e ;
110   fSumEta+=eta*e ;
111   fSumPhi+=phi*e ;    
112   
113   //check if this a leading particle?
114   if(fLEnergy < e){
115     fLEnergy = e;
116     fLEta = eta ;
117     fLPhi = phi ;
118   }
119
120   if(index >=0){ //add index to list of indexes
121     if(!fList){
122       fList = new TArrayI(100) ;
123     }
124     if(fList->GetSize()<=fNpart+1){
125       TArrayI * tmp = new TArrayI(fNpart*2) ;
126       tmp->Adopt(fList->GetSize(),fList->GetArray()) ; //note, we should not delete old array!
127       fList=tmp ;
128     }
129     fList->AddAt(index,fNpart++) ;
130   }
131 }
132 // //____________________________________________________________________________ 
133 // Bool_t AliPHOSJet::IsInJet(TParticle * p,Int_t mode)
134 // {
135 //   Double_t dEta ;
136 //   Double_t dPhi ;
137 //   Double_t energy ;
138 //   if(!fEnergy){ //Final values not calculated yet, use intermediate
139 //     if(fSumEnergy==0)
140 //       return kTRUE ; //First particle
141 // note p->Eta() causes fpe! 
142 //     dEta=(p->Eta() - fSumEta/fSumEnergy) ;
143 //     dPhi=(p->Phi() - fSumPhi/fSumEnergy) ;
144 //     energy = fSumEnergy ;
145 //   }
146 //   else{ //calculated with respect to final 
147 //     dEta=(p->Eta() - fEta) ;
148 //     dPhi=(p->Phi() - fPhi) ;
149 //     energy = fEnergy ;    
150 //   }
151   
152 //   switch (mode){
153 //   case 0: 
154 //     return IsInCone(eta,phi) ; //pure geometrical comparison
155 //   case 1: 
156 //     return AcceptConeDeviation(dEta,dPhi,p->Energy() );
157 //   default:
158 //     Error("IsInCone","Unknown mode of cone calculation %d \n",mode );
159 //   }
160 //   return kFALSE ;
161 //}
162 //____________________________________________________________________________ 
163 Bool_t AliPHOSJet::AcceptConeDeviation(const TParticle * p)const
164 { //Calculate cone deviation in case of inclusion of the given
165   //particle to jet. 
166
167   Double_t tmpEnergy = fSumEnergy + p->Energy() ;  
168   Double_t tmpEta = fSumEta ;
169   if(p->Pt()/p->Energy() >0.001)
170     tmpEta+= p->Eta()*p->Energy() ;
171   else
172     tmpEta+= 100.*p->Pz() ;
173   tmpEta = tmpEta/tmpEnergy / (fSumEta/fSumEnergy) - 1 ;
174   Double_t tmpPhi = fSumPhi + p->Phi()*p->Energy() ;
175   tmpPhi = tmpPhi/tmpEnergy /(fSumPhi/fSumEnergy) - 1 ;
176   Double_t dev = TMath::Sqrt(tmpEta*tmpEta + tmpPhi*tmpPhi) ;
177   if(dev < fMaxConeMove)
178     return kTRUE ;
179   else
180     return kFALSE ;
181 }
182 //____________________________________________________________________________ 
183 Bool_t AliPHOSJet::AcceptConeDeviation(const Double_t e,const Double_t eta,const Double_t phi)const
184 { //Calculate cone deviation in case of inclusion of the given
185   //particle to jet. 
186
187   Double_t tmpEnergy = fSumEnergy + e ;
188   Double_t tmpEta = fSumEta + eta*e ;
189   tmpEta = tmpEta/tmpEnergy / (fSumEta/fSumEnergy) ;
190   Double_t tmpPhi = fSumPhi + phi*e ;
191   tmpPhi = tmpPhi/tmpEnergy /(fSumPhi/fSumEnergy) ;
192   Double_t dev = TMath::Sqrt(tmpEta*tmpEta + tmpPhi*tmpPhi) ;
193   if(dev<fMaxConeMove && dev > fMinConeMove)
194     return kTRUE ;
195   else
196     return kFALSE ;
197 }
198 //____________________________________________________________________________ 
199 Bool_t AliPHOSJet::IsInCone(const TParticle * p)const
200 {
201   //Say if  particle is inside the defined cone
202   Double_t dEta ;
203   Double_t dPhi ;
204   if(!fEnergy){ //Final values not calculated yet, use intermediate
205     if(fSumEnergy==0)
206       return kTRUE ; //First particle    
207     if(p->Pt()/p->Energy() > 0.001)
208       dEta=(p->Eta() - fSumEta/fSumEnergy) ;
209     else
210       dEta=(100.*p->Pz()/p->Energy() - fSumEta/fSumEnergy) ;
211     dPhi=(p->Phi() - fSumPhi/fSumEnergy) ;
212   }
213   else{ //calculated with respect to final 
214     if(p->Pt()/p->Energy() > 0.001)
215       dEta=(p->Eta() - fEta) ;
216     else
217       dEta=(100.*p->Pz()/p->Energy() - fEta) ;
218     dPhi=(p->Phi() - fPhi) ;
219   }
220   if(TMath::Sqrt(dEta*dEta + dPhi*dPhi) < fConeRad)
221     return kTRUE ;
222   else
223     return kFALSE ;
224 }
225 //____________________________________________________________________________ 
226 Bool_t AliPHOSJet::IsInCone(const Double_t e,const Double_t eta,const Double_t phi)const
227 {
228   //Says if particle is inside the defined cone
229   Double_t dEta ;
230   Double_t dPhi ;
231   if(!fEnergy){ //Final values not calculated yet, use intermediate
232     if(fSumEnergy==0)
233       return kTRUE ; //First particle    
234     dEta=(eta - fSumEta/fSumEnergy) ;
235     dPhi=(phi - fSumPhi/fSumEnergy) ;
236   }
237   else{ //calculated with respect to final 
238     dEta=(eta - fEta) ;
239     dPhi=(phi - fPhi) ;
240   }
241   if(TMath::Sqrt(dEta*dEta + dPhi*dPhi) < fConeRad)
242     return kTRUE ;
243   else
244     return kFALSE ;
245 }
246 //____________________________________________________________________________ 
247 Double_t AliPHOSJet::DistanceToJet(const TParticle *p)const{
248   //Calculate radius 
249   Double_t dEta ;
250   Double_t dPhi ;
251   if(!fEnergy){ //Final values not calculated yet, use intermediate
252     if(fSumEnergy==0)
253       return kTRUE ; //First particle    
254     if(p->Pt()/p->Energy() > 0.001)
255       dEta=(p->Eta() - fSumEta/fSumEnergy) ;
256     else
257       dEta=(100.*p->Pz()/p->Energy() - fSumEta/fSumEnergy) ;      
258     dPhi=(p->Phi() - fSumPhi/fSumEnergy) ;
259   }
260   else{ //calculated with respect to final 
261     if(p->Pt()/p->Energy() > 0.001)
262       dEta=(p->Eta() - fEta) ;
263     else
264       dEta=(100*p->Pz()/p->Energy() - fEta) ;    
265     dPhi=(p->Phi() - fPhi) ;
266   }
267   return TMath::Sqrt(dEta*dEta + dPhi*dPhi) ;
268
269 }
270 //____________________________________________________________________________ 
271 void AliPHOSJet::CalculateAll(void){
272   //Calculate all jet parameters
273   if(fSumEnergy==0)
274     return  ; //Nothing to calculate    
275   
276   fEta = fSumEta/fSumEnergy ;
277   fPhi = fSumPhi/fSumEnergy ;
278   fEnergy = fSumEnergy ;
279   
280   fSumEnergy = 0. ;
281   fSumEta = 0. ;
282   fSumPhi = 0. ;
283 }
284 //____________________________________________________________________________ 
285 void AliPHOSJet::Print(Option_t * option){
286   //Print jet parameters
287   printf("-------------- AliPHOSJet ------------\n") ;
288   printf(" Energy............. %f \n",fEnergy) ;
289   printf(" Eta................ %f \n",fEta ) ;
290   printf(" Phi................ %f \n",fPhi ) ;
291   printf(" Leading Energy..... %f \n",fLEnergy) ;
292   printf(" Leading Eta........ %f \n",fLEta) ;
293   printf(" Leading Phi........ %f \n",fLPhi) ;
294   printf(" N particles in jet  %d \n",fNpart) ;
295   printf("----------------------------------\n") ;
296 }
297
298
299