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