[u/mrichter/AliRoot.git] / PHOS / AliPHOSJetFinder.cxx.old

/**************************************************************************
 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
 *                                                                        *
 * Author: The ALICE Off-line Project.                                    *
 * Contributors are mentioned in the code where appropriate.              *
 *                                                                        *
 * Permission to use, copy, modify and distribute this software and its   *
 * documentation strictly for non-commercial purposes is hereby granted   *
 * without fee, provided that the above copyright notice appears in all   *
 * copies and that both the copyright notice and this permission notice   *
 * appear in the supporting documentation. The authors make no claims     *
 * about the suitability of this software for any purpose. It is          *
 * provided "as is" without express or implied warranty.                  *
 **************************************************************************/

//_________________________________________________________________________
// C++ version of UA2 and/or Lund jet finding algorithm
//  UA1 jet algorithm from LUND JETSET (LUCELL)
//
//*-- Author : D.Peressounko after UA1 coll. etc
//////////////////////////////////////////////////////////////////////////////

// --- ROOT system ---
#include "TClonesArray.h"
//      #include "TIter.h"
#include "TParticle.h"
#include "TTask.h"
// --- Standard library ---

// --- AliRoot header files ---
#include "AliJet.h"
#include "AliPHOSDigit.h"
#include "AliPHOSGetter.h"
#include "AliPHOSJetFinder.h"
#include "AliPHOSDigitizer.h"

ClassImp(AliPHOSJetFinder)


//____________________________________________________________________________ 
  AliPHOSJetFinder::AliPHOSJetFinder():TNamed("AliPHOSJetFinder","") 
{
  fNJets = 0 ;
  fMode  = 0 ;   //No iterations 
  fStatusCode = -999 ; //no selection

  fConeRad = 1.;   //Radius of jet value?????????
  fMaxConeMove = 0.15 ; //value???????
  fMinConeMove = 0.05 ; //value???????
  fEtSeed = 4. ; //Energy of seed particle value??????????
  fEtMin = 5.;   //minimal energy of jet   value??????????
  fPrecBg = 0.00035 ; //value????????
  fSimGain = 0.;
  fSimPedestal = 0.;

  fParticles = 0;
  fJets = 0 ;
  
}

//____________________________________________________________________________ 
  AliPHOSJetFinder::~AliPHOSJetFinder()
{
  if(fParticles){
    delete fParticles ;
    fParticles = 0 ;
  }
  if(fJets){
    delete fJets ;
    fJets = 0 ;
  } 
}

//____________________________________________________________________________ 
void  AliPHOSJetFinder::FindJetsFromParticles(const TClonesArray * plist,TObjArray * jetslist) 
{

  TIter next(plist) ;

  TIter nextJet(jetslist) ;

  fNJets = 0 ;
  TParticle * p ;
  AliJet * jet ;
  //In this cicle we find number of jets and define approx. their directions
  //note, that we do not really add particles to jet (index =-1)
  while((p=static_cast<TParticle*>(next()))){
    if(fStatusCode==-999 || p->GetStatusCode()==fStatusCode){
    if(p->Energy() >= fEtSeed){ //Energetic enough
      //cout << "p " << p->Energy() << endl ;
      //cout << "Status "<<fStatusCode<<" "<<p->GetName()<<" " << p->Energy() << " "<<p->Eta()<< " "<<p->Phi()<<endl ;
      Bool_t startnew = kTRUE ;
      //Do not start new jet if part of older jet
      nextJet.Reset() ;
      while((jet=static_cast<AliJet*>(nextJet()))){
	//jet->Print() ;
	if(jet->AcceptConeDeviation(p)){
	  startnew = kFALSE ;
	  //cout << "false" << endl ;
	  break ;
	}
      }
      if(startnew){
	//cout << "new " << endl ;
	jet = new AliJet() ;
	jetslist->Add(jet) ;
	//	jet = static_cast<AliJet*>(jetslist->Last()) ;
	jet->SetConeRadius(fConeRad) ;
	jet->SetMaxConeMove(fMaxConeMove) ;
	//	jet->SetMinConeMove(fMinConeMove) ;
	jet->AddParticle(p,-1) ;
	fNJets++;
      }
    }
    while((jet=static_cast<AliJet*>(nextJet()))){
      if(jet->AcceptConeDeviation(p))
	jet->AddParticle(p,-1) ; //Just recalculate direction of jet
    }
    }
  }
    
  //now calculate directions of jets using collected information
  nextJet.Reset() ;
  while((jet=static_cast<AliJet*>(nextJet()))){
    jet->CalculateAll() ;
    if(jet->Energy() < fEtMin){
      jetslist->Remove(jet) ;
      delete jet ;
    }
  }

  jetslist->Compress() ;
  //And finally, really add particles to jets
  for(Int_t iPart=0; iPart<plist->GetEntries();iPart++){
    p=static_cast<TParticle*>(plist->At(iPart)) ;
    if(fStatusCode == -999 || p->GetStatusCode()==fStatusCode){
    Double_t dist = 999999. ; //big distance
    Int_t iJet = -1 ;
    for(Int_t i=0; i<jetslist->GetEntriesFast();i++){
      jet=static_cast<AliJet*>(jetslist->At(i)) ;
      if(jet->IsInCone(p)){
	Double_t cdist = jet->DistanceToJet(p); 
	if(cdist < dist){
	  dist = cdist ;
	  iJet = i ;
	}
      }
    }
    if(iJet>-1)
      (static_cast<AliJet*>(jetslist->At(iJet)))->AddParticle(p,iPart); //assign particle to closest jet
  }
  }
  
  //Calculate jet parameters 
  nextJet.Reset() ;
  while((jet=static_cast<AliJet*>(nextJet()))){
    jet->CalculateAll() ;
  }

}
//____________________________________________________________________________ 
void AliPHOSJetFinder::FindJetsFromDigits(const TClonesArray * digits, TObjArray * jets){

  if(digits->GetEntries()==0){
    Error("JetsFromDigits","No entries in digits list \n") ;
    return ;
  }


  TClonesArray * copyDigits = new TClonesArray(*digits) ;

  //Remove CPV digits if any
  AliPHOSGeometry * geom = AliPHOSGeometry::GetInstance("GPS2","") ;
  Int_t iDigit ;
  AliPHOSDigit * digit ;
  for(iDigit=copyDigits->GetEntries()-1;iDigit>=0;iDigit--){
    digit=static_cast<AliPHOSDigit *>(copyDigits->At(iDigit)) ;
    if(!geom->IsInEMC(digit->GetId()))
      copyDigits->RemoveAt(iDigit) ;
    else
      break ;
  }
  copyDigits->Compress() ;

  Double_t totalEnergy = 0 ;
  Float_t energy[copyDigits->GetEntries()] ;
  //calculate average energy of digits
  //fill array of energies
  for(iDigit=0;iDigit<copyDigits->GetEntries();iDigit++){
    digit=static_cast<AliPHOSDigit *>(copyDigits->At(iDigit)) ;
    energy[iDigit] = Calibrate(digit) ;
    totalEnergy+=energy[iDigit] ;
  }
  
  //Sort digits in decreasing energy.
  Int_t index[copyDigits->GetEntries()] ;
  TMath::Sort(copyDigits->GetEntries(),energy,index) ;
  
  Double_t eAverage = totalEnergy/copyDigits->GetEntries()  ;
  //remove digits below average energy
  for(iDigit=copyDigits->GetEntries()-1;iDigit>=0;iDigit--){
    digit=static_cast<AliPHOSDigit *>(copyDigits->At(index[iDigit])) ;
    if(energy[index[iDigit]] < eAverage)
      copyDigits->RemoveAt(iDigit) ;
    else
      break ;
  }
  
  
  AliJet * jet ;
  Int_t iIter = 0 ;
  while(iIter < 10){//less than 10 iterations
    
    //while digits above seed
    for(Int_t ind=0;ind<copyDigits->GetEntriesFast();ind++){
      digit=static_cast<AliPHOSDigit*>(copyDigits->At(index[ind])) ;
      if(energy[index[ind]] > fEtSeed && digit){ //start new jet      
 	jet = new AliJet() ;
	Double_t e,eta,phi ;
	CalculateEEtaPhi(digit,e,eta,phi) ;
	jet->AddDigit(e,eta,phi,-1) ;  
	//loop over left digits
	for(Int_t iDigit = 0 ; iDigit < copyDigits->GetEntries() ; iDigit++){
	  if(iDigit!= ind){ //first digit already in jet
	    digit = static_cast<AliPHOSDigit *>(copyDigits->At(iDigit));
	    CalculateEEtaPhi(digit,e,eta,phi) ;	    
	    if(jet->IsInCone(e,eta,phi) && //is cell in cone
	       jet->AcceptConeDeviation(e,eta,phi)){//if cone does not move too much	      
	      jet->AddDigit(e,eta,phi,-1) ;  //accept new direction
	    }
	  }
	}//end of loop over cells
	
	//accept  all anused cells incide cone
	//note, that digits might be returned as anused later
	for(Int_t icell = 0 ; icell < copyDigits->GetEntries() ; icell++){
	  digit = static_cast<AliPHOSDigit *>(copyDigits->At(icell));
	  if(jet->IsInCone(e,eta,phi)){ //is cell in cone
	    CalculateEEtaPhi(digit,e,eta,phi) ;
	    jet->AddDigit(e,eta,phi,digit->GetIndexInList()) ; 
	  }
	}
	
	//Accept Jet with Et > Et_min and remove all belonging digits
	if(jet->Energy()/TMath::CosH(jet->Eta()) > fEtMin){
	  Int_t nIndxs ;
	  const Int_t * indxs = jet->Indexs(nIndxs) ;
	  for(Int_t i=0;i<nIndxs;i++){
	    copyDigits->RemoveAt(indxs[i]) ;
	  }
	  jet->CalculateAll() ;
	  jets->AddAt(jet,fNJets++);
	}
	else{ //remove jet and do not touch digits
	  delete jet ;
	}
      }
      else{ 
	if(energy[index[ind]] < fEtSeed){ // no more digits above threshold left, return from loop
	  break ;
	}
      }
      
      iIter++ ;
      //calculate new energy of backrgound
      Double_t oldTotalEnergy = totalEnergy ;
      totalEnergy = 0 ;
      for(Int_t i=0 ; i<copyDigits->GetEntriesFast() ; i++){
	digit=static_cast<AliPHOSDigit*>(copyDigits->At(index[ind])) ;
	if(digit)
	  totalEnergy+=energy[i] ;
      }
      if(!fMode || (oldTotalEnergy != 0) && 
	 (TMath::Abs(oldTotalEnergy - totalEnergy)/oldTotalEnergy < fPrecBg))
	break ;	
    }
  }

  copyDigits->Delete() ;
  
}
//____________________________________________________________________________ 
Double_t AliPHOSJetFinder::Calibrate(const AliPHOSDigit * digit){ 
//   if(fPedestals || fGains ){  //use calibration data
//     if(!fPedestals || !fGains ){
//       Error("Calibrate","Either Pedestals of Gains not set!") ;
//       return 0 ;
//     }
//     Float_t en=(digit->GetAmp() - fPedestals->Data(digit->GetId)()))*fGains->Data(digit->GetId()) ;
//     if(en>0)
//       return en ;
//     else
//       return 0. ;
//   }
//   else{ //simulation
  if(fSimGain==0){ //read simulation parameters
    AliPHOSGetter * gime = AliPHOSGetter::GetInstance() ;
    if(!gime){
      Error("Calibrate","Can not read Calibration parameters") ;
      return 0 ;
    }
    const TTask * task = gime->Digitizer() ;
    if(strcmp(task->IsA()->GetName(),"AliPHOSDigitizer")==0){
      const AliPHOSDigitizer * dig = static_cast<const AliPHOSDigitizer *>(task) ;
      fSimGain     = dig->GetEMCchannel() ;
      fSimPedestal = dig->GetEMCpedestal();
    }
  }
  
  return fSimPedestal + digit->GetAmp()*fSimGain ;        
  
  //  }  
}
//____________________________________________________________________________ 
void AliPHOSJetFinder::CalculateEEtaPhi(const AliPHOSDigit * d,Double_t &e, Double_t &eta, Double_t &phi){
  e=Calibrate(d) ;
  AliPHOSGeometry * geom = AliPHOSGeometry::GetInstance("GPS2","") ;
  TVector3 pos ;
  geom->RelPosInAlice(d->GetId(), pos) ;
  eta = pos.Eta() ;
  phi = pos.Phi() ;
}
//____________________________________________________________________________ 
void AliPHOSJetFinder::Print(Option_t * option){	
  printf("\n --------------- AliPHOSJetFinder --------------- \n") ;
  printf(" Jets found .........%d \n",fNJets) ;
  printf(" Seed energy cut ....%f \n",fEtSeed) ;
  printf(" Cone radius ........%f \n",fConeRad) ;
  printf(" Minimal cone move ..%f \n",fMinConeMove) ;
  printf(" Maximal cone move ..%f \n",fMaxConeMove) ;
  printf("------------------------------------------------- \n") ;
}
Commit	Line	Data
5240525e	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	// C++ version of UA2 and/or Lund jet finding algorithm
	18	// UA1 jet algorithm from LUND JETSET (LUCELL)
	19	//
	20	//*-- Author : D.Peressounko after UA1 coll. etc
	21	//////////////////////////////////////////////////////////////////////////////
	22
	23	// --- ROOT system ---
	24	#include "TClonesArray.h"
	25	// #include "TIter.h"
	26	#include "TParticle.h"
	27	#include "TTask.h"
	28	// --- Standard library ---
	29
	30	// --- AliRoot header files ---
	31	#include "AliJet.h"
	32	#include "AliPHOSDigit.h"
	33	#include "AliPHOSGetter.h"
	34	#include "AliPHOSJetFinder.h"
	35	#include "AliPHOSDigitizer.h"
	36
	37	ClassImp(AliPHOSJetFinder)
	38
	39
	40	//____________________________________________________________________________
	41	AliPHOSJetFinder::AliPHOSJetFinder():TNamed("AliPHOSJetFinder","")
	42	{
	43	fNJets = 0 ;
	44	fMode = 0 ; //No iterations
	45	fStatusCode = -999 ; //no selection
	46
	47	fConeRad = 1.; //Radius of jet value?????????
	48	fMaxConeMove = 0.15 ; //value???????
	49	fMinConeMove = 0.05 ; //value???????
	50	fEtSeed = 4. ; //Energy of seed particle value??????????
	51	fEtMin = 5.; //minimal energy of jet value??????????
	52	fPrecBg = 0.00035 ; //value????????
	53	fSimGain = 0.;
	54	fSimPedestal = 0.;
	55
	56	fParticles = 0;
	57	fJets = 0 ;
	58
	59	}
	60
	61	//____________________________________________________________________________
	62	AliPHOSJetFinder::~AliPHOSJetFinder()
	63	{
	64	if(fParticles){
65	delete fParticles ;
66	fParticles = 0 ;
67	}
68	if(fJets){
69	delete fJets ;
70	fJets = 0 ;
71	}
72	}
73
74	//____________________________________________________________________________
75	void AliPHOSJetFinder::FindJetsFromParticles(const TClonesArray * plist,TObjArray * jetslist)
76	{
77
78	TIter next(plist) ;
79
80	TIter nextJet(jetslist) ;
81
82	fNJets = 0 ;
83	TParticle * p ;
84	AliJet * jet ;
85	//In this cicle we find number of jets and define approx. their directions
86	//note, that we do not really add particles to jet (index =-1)
87	while((p=static_cast<TParticle*>(next()))){
88	if(fStatusCode==-999 \|\| p->GetStatusCode()==fStatusCode){
89	if(p->Energy() >= fEtSeed){ //Energetic enough
90	//cout << "p " << p->Energy() << endl ;
91	//cout << "Status "<<fStatusCode<<" "<<p->GetName()<<" " << p->Energy() << " "<<p->Eta()<< " "<<p->Phi()<<endl ;
92	Bool_t startnew = kTRUE ;
93	//Do not start new jet if part of older jet
94	nextJet.Reset() ;
95	while((jet=static_cast<AliJet*>(nextJet()))){
96	//jet->Print() ;
97	if(jet->AcceptConeDeviation(p)){
98	startnew = kFALSE ;
99	//cout << "false" << endl ;
100	break ;
101	}
102	}
103	if(startnew){
104	//cout << "new " << endl ;
105	jet = new AliJet() ;
106	jetslist->Add(jet) ;
107	// jet = static_cast<AliJet*>(jetslist->Last()) ;
108	jet->SetConeRadius(fConeRad) ;
109	jet->SetMaxConeMove(fMaxConeMove) ;
110	// jet->SetMinConeMove(fMinConeMove) ;
111	jet->AddParticle(p,-1) ;
112	fNJets++;
113	}
114	}
115	while((jet=static_cast<AliJet*>(nextJet()))){
116	if(jet->AcceptConeDeviation(p))
117	jet->AddParticle(p,-1) ; //Just recalculate direction of jet
118	}
119	}
120	}
121
122	//now calculate directions of jets using collected information
123	nextJet.Reset() ;
124	while((jet=static_cast<AliJet*>(nextJet()))){
125	jet->CalculateAll() ;
126	if(jet->Energy() < fEtMin){
127	jetslist->Remove(jet) ;
128	delete jet ;
129	}
130	}
131
132	jetslist->Compress() ;
133	//And finally, really add particles to jets
134	for(Int_t iPart=0; iPart<plist->GetEntries();iPart++){
135	p=static_cast<TParticle*>(plist->At(iPart)) ;
136	if(fStatusCode == -999 \|\| p->GetStatusCode()==fStatusCode){
137	Double_t dist = 999999. ; //big distance
138	Int_t iJet = -1 ;
139	for(Int_t i=0; i<jetslist->GetEntriesFast();i++){
140	jet=static_cast<AliJet*>(jetslist->At(i)) ;
141	if(jet->IsInCone(p)){
142	Double_t cdist = jet->DistanceToJet(p);
143	if(cdist < dist){
144	dist = cdist ;
145	iJet = i ;
146	}
147	}
148	}
149	if(iJet>-1)
150	(static_cast<AliJet*>(jetslist->At(iJet)))->AddParticle(p,iPart); //assign particle to closest jet
151	}
152	}
153
154	//Calculate jet parameters
155	nextJet.Reset() ;
156	while((jet=static_cast<AliJet*>(nextJet()))){
157	jet->CalculateAll() ;
158	}
159
160	}
161	//____________________________________________________________________________
162	void AliPHOSJetFinder::FindJetsFromDigits(const TClonesArray * digits, TObjArray * jets){
163
164	if(digits->GetEntries()==0){
165	Error("JetsFromDigits","No entries in digits list \n") ;
166	return ;
167	}
168
169
170	TClonesArray * copyDigits = new TClonesArray(*digits) ;
171
172	//Remove CPV digits if any
173	AliPHOSGeometry * geom = AliPHOSGeometry::GetInstance("GPS2","") ;
174	Int_t iDigit ;
175	AliPHOSDigit * digit ;
176	for(iDigit=copyDigits->GetEntries()-1;iDigit>=0;iDigit--){
177	digit=static_cast<AliPHOSDigit *>(copyDigits->At(iDigit)) ;
178	if(!geom->IsInEMC(digit->GetId()))
179	copyDigits->RemoveAt(iDigit) ;
180	else
181	break ;
182	}
183	copyDigits->Compress() ;
184
185	Double_t totalEnergy = 0 ;
186	Float_t energy[copyDigits->GetEntries()] ;
187	//calculate average energy of digits
188	//fill array of energies
189	for(iDigit=0;iDigit<copyDigits->GetEntries();iDigit++){
190	digit=static_cast<AliPHOSDigit *>(copyDigits->At(iDigit)) ;
191	energy[iDigit] = Calibrate(digit) ;
192	totalEnergy+=energy[iDigit] ;
193	}
194
195	//Sort digits in decreasing energy.
196	Int_t index[copyDigits->GetEntries()] ;
197	TMath::Sort(copyDigits->GetEntries(),energy,index) ;
198
199	Double_t eAverage = totalEnergy/copyDigits->GetEntries() ;
200	//remove digits below average energy
201	for(iDigit=copyDigits->GetEntries()-1;iDigit>=0;iDigit--){
202	digit=static_cast<AliPHOSDigit *>(copyDigits->At(index[iDigit])) ;
203	if(energy[index[iDigit]] < eAverage)
204	copyDigits->RemoveAt(iDigit) ;
205	else
206	break ;
207	}
208
209
210	AliJet * jet ;
211	Int_t iIter = 0 ;
212	while(iIter < 10){//less than 10 iterations
213
214	//while digits above seed
215	for(Int_t ind=0;ind<copyDigits->GetEntriesFast();ind++){
216	digit=static_cast<AliPHOSDigit*>(copyDigits->At(index[ind])) ;
217	if(energy[index[ind]] > fEtSeed && digit){ //start new jet
218	jet = new AliJet() ;
219	Double_t e,eta,phi ;
220	CalculateEEtaPhi(digit,e,eta,phi) ;
221	jet->AddDigit(e,eta,phi,-1) ;
222	//loop over left digits
223	for(Int_t iDigit = 0 ; iDigit < copyDigits->GetEntries() ; iDigit++){
224	if(iDigit!= ind){ //first digit already in jet
225	digit = static_cast<AliPHOSDigit *>(copyDigits->At(iDigit));
226	CalculateEEtaPhi(digit,e,eta,phi) ;
227	if(jet->IsInCone(e,eta,phi) && //is cell in cone
228	jet->AcceptConeDeviation(e,eta,phi)){//if cone does not move too much
229	jet->AddDigit(e,eta,phi,-1) ; //accept new direction
230	}
231	}
232	}//end of loop over cells
233
234	//accept all anused cells incide cone
235	//note, that digits might be returned as anused later
236	for(Int_t icell = 0 ; icell < copyDigits->GetEntries() ; icell++){
237	digit = static_cast<AliPHOSDigit *>(copyDigits->At(icell));
238	if(jet->IsInCone(e,eta,phi)){ //is cell in cone
239	CalculateEEtaPhi(digit,e,eta,phi) ;
240	jet->AddDigit(e,eta,phi,digit->GetIndexInList()) ;
241	}
242	}
243
244	//Accept Jet with Et > Et_min and remove all belonging digits
245	if(jet->Energy()/TMath::CosH(jet->Eta()) > fEtMin){
246	Int_t nIndxs ;
247	const Int_t * indxs = jet->Indexs(nIndxs) ;
248	for(Int_t i=0;i<nIndxs;i++){
249	copyDigits->RemoveAt(indxs[i]) ;
250	}
251	jet->CalculateAll() ;
252	jets->AddAt(jet,fNJets++);
253	}
254	else{ //remove jet and do not touch digits
255	delete jet ;
256	}
257	}
258	else{
259	if(energy[index[ind]] < fEtSeed){ // no more digits above threshold left, return from loop
260	break ;
261	}
262	}
263
264	iIter++ ;
265	//calculate new energy of backrgound
266	Double_t oldTotalEnergy = totalEnergy ;
267	totalEnergy = 0 ;
268	for(Int_t i=0 ; i<copyDigits->GetEntriesFast() ; i++){
269	digit=static_cast<AliPHOSDigit*>(copyDigits->At(index[ind])) ;
270	if(digit)
271	totalEnergy+=energy[i] ;
272	}
273	if(!fMode \|\| (oldTotalEnergy != 0) &&
274	(TMath::Abs(oldTotalEnergy - totalEnergy)/oldTotalEnergy < fPrecBg))
275	break ;
276	}
277	}
278
279	copyDigits->Delete() ;
280
281	}
282	//____________________________________________________________________________
283	Double_t AliPHOSJetFinder::Calibrate(const AliPHOSDigit * digit){
284	// if(fPedestals \|\| fGains ){ //use calibration data
285	// if(!fPedestals \|\| !fGains ){
286	// Error("Calibrate","Either Pedestals of Gains not set!") ;
287	// return 0 ;
288	// }
289	// Float_t en=(digit->GetAmp() - fPedestals->Data(digit->GetId)()))*fGains->Data(digit->GetId()) ;
290	// if(en>0)
291	// return en ;
292	// else
293	// return 0. ;
294	// }
295	// else{ //simulation
296	if(fSimGain==0){ //read simulation parameters
297	AliPHOSGetter * gime = AliPHOSGetter::GetInstance() ;
298	if(!gime){
299	Error("Calibrate","Can not read Calibration parameters") ;
300	return 0 ;
301	}
302	const TTask * task = gime->Digitizer() ;
303	if(strcmp(task->IsA()->GetName(),"AliPHOSDigitizer")==0){
304	const AliPHOSDigitizer * dig = static_cast<const AliPHOSDigitizer *>(task) ;
305	fSimGain = dig->GetEMCchannel() ;
306	fSimPedestal = dig->GetEMCpedestal();
307	}
308	}
309
310	return fSimPedestal + digit->GetAmp()*fSimGain ;
311
312	// }
313	}
314	//____________________________________________________________________________
315	void AliPHOSJetFinder::CalculateEEtaPhi(const AliPHOSDigit * d,Double_t &e, Double_t &eta, Double_t &phi){
316	e=Calibrate(d) ;
317	AliPHOSGeometry * geom = AliPHOSGeometry::GetInstance("GPS2","") ;
318	TVector3 pos ;
319	geom->RelPosInAlice(d->GetId(), pos) ;
320	eta = pos.Eta() ;
321	phi = pos.Phi() ;
322	}
323	//____________________________________________________________________________
324	void AliPHOSJetFinder::Print(Option_t * option){
325	printf("\n --------------- AliPHOSJetFinder --------------- \n") ;
326	printf(" Jets found .........%d \n",fNJets) ;
327	printf(" Seed energy cut ....%f \n",fEtSeed) ;
328	printf(" Cone radius ........%f \n",fConeRad) ;
329	printf(" Minimal cone move ..%f \n",fMinConeMove) ;
330	printf(" Maximal cone move ..%f \n",fMaxConeMove) ;
331	printf("------------------------------------------------- \n") ;
332	}