2 /**************************************************************************
3 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
5 * Author: The ALICE Off-line Project. *
6 * Contributors are mentioned in the code where appropriate. *
8 * Permission to use, copy, modify and distribute this software and its *
9 * documentation strictly for non-commercial purposes is hereby granted *
10 * without fee, provided that the above copyright notice appears in all *
11 * copies and that both the copyright notice and this permission notice *
12 * appear in the supporting documentation. The authors make no claims *
13 * about the suitability of this software for any purpose. It is *
14 * provided "as is" without express or implied warranty. *
15 **************************************************************************/
17 //-------------------------------------------------------------
18 // Fill Unit Array with EMCal information
19 // Called by ESD reader for jet analysis
20 // Author: Magali Estienne (magali.estienne@ires.in2p3.fr)
21 //-------------------------------------------------------------
23 // --- Standard library ---
24 #include <Riostream.h>
26 // --- ROOT system ---
28 #include <TLorentzVector.h>
29 #include <TRefArray.h>
31 #include <TGeoManager.h>
33 #include <TClonesArray.h>
36 // --- AliRoot header files ---
37 #include "AliJetFinder.h"
38 #include "AliJetReaderHeader.h"
39 #include "AliJetReader.h"
40 #include "AliJetESDReader.h"
41 #include "AliJetESDReaderHeader.h"
43 #include "AliESDEvent.h"
44 #include "AliJetDummyGeo.h"
45 #include "AliESDCaloCluster.h"
46 #include "AliJetUnitArray.h"
47 #include "AliJetFillUnitArrayEMCalDigits.h"
49 ClassImp(AliJetFillUnitArrayEMCalDigits)
51 //_____________________________________________________________________________
52 AliJetFillUnitArrayEMCalDigits::AliJetFillUnitArrayEMCalDigits()
53 : TTask("AliJetFillUnitArrayEMCalDigits","Fill Unit Array with tpc/its and emcal information"),
80 //_____________________________________________________________________________
81 AliJetFillUnitArrayEMCalDigits::AliJetFillUnitArrayEMCalDigits(AliESDEvent */*esd*/)
82 : TTask("AliJetFillUnitArrayEMCalDigits","Fill Unit Array with tpc/its and emcal information"),
110 //____________________________________________________________________________
111 void AliJetFillUnitArrayEMCalDigits::InitParameters()
113 fNumUnits = fGeom->GetNCells(); // Number of towers in EMCAL
115 fEtaMinCal = fGeom->GetArm1EtaMin();
116 fEtaMaxCal = fGeom->GetArm1EtaMax();
117 fPhiMinCal = fGeom->GetArm1PhiMin();
118 fPhiMaxCal = fGeom->GetArm1PhiMax();
121 if(fDebug>1) printf("\n EMCAL parameters initiated ! \n");
125 //_____________________________________________________________________________
126 AliJetFillUnitArrayEMCalDigits::~AliJetFillUnitArrayEMCalDigits()
131 //_____________________________________________________________________________
132 void AliJetFillUnitArrayEMCalDigits::Exec(Option_t* /*option*/)
138 fDebug = fReaderHeader->GetDebug();
139 fOpt = fReaderHeader->GetDetector();
144 // Get number of clusters from EMCAL
150 Float_t ptMin = fReaderHeader->GetPtCut();
152 // Loop over calo clusters
153 //------------------------------------------------------------------
157 // Total number of EMCAL cluster
158 end = fESD->GetNumberOfCaloClusters();
160 for(Int_t j = beg; j < end; j++) {
161 fClus = fESD->GetCaloCluster(j);
162 if(!fClus->IsEMCAL()) continue;
164 type = fClus->GetClusterType();
165 index = fClus->GetID();
166 nDigitTot = fClus->GetNumberOfDigits();
168 // Keep clusters or pseudo clusters
169 if (type != AliESDCaloCluster::kEMCALClusterv1) continue;
170 // if (type != AliESDCaloCluster::kPseudoCluster) continue;
172 // Get the digit index and the digit information
173 //============================================================
175 // Get number of digits in a cluster
176 Int_t nD = fClus->GetNumberOfDigits();
178 TArrayS *digID = fClus->GetDigitIndex();
179 TArrayS *digEnergy = fClus->GetDigitAmplitude();
180 Float_t *digitEnergy = new Float_t[nD];
181 // Float_t digitEn = 0.;
184 for(Int_t k=0; k<nD; k++) {
186 // Convert energy in GeV
187 Int_t idF = (Int_t)digID->At(k);
188 // Calibration for an energy in GeV
189 digitEnergy[k] = (Float_t)digEnergy->At(k)/500.;
191 // Second method to extract eta, phi positions of a digit
192 //=================================================================
194 Float_t etaD=-10., phiD=-10.;
195 fGeom->EtaPhiFromIndex(idF,etaD,phiD);
196 // fEMCalGrid->GetEtaPhiFromIndex2(idF,phiD,etaD);
197 phiD = ((phiD < 0) ? phiD + 2.* TMath::Pi() : phiD);
199 Float_t etDigit = digitEnergy[k]*TMath::Abs(TMath::Sin(EtaToTheta(etaD)));
201 AliJetUnitArray *uArray = (AliJetUnitArray*)fUnitArray->At(idF);
202 if(uArray->GetUnitEnergy() == 0.) goodDigit++;
204 Float_t unitEnergy = 0.;
206 unitEnergy = uArray->GetUnitEnergy();
208 fRefArray->Add(uArray);
212 uArray->SetUnitEnergy(unitEnergy+etDigit);
214 if(uArray->GetUnitEnergy()<ptMin)
215 uArray->SetUnitCutFlag(kPtSmaller);
217 uArray->SetUnitCutFlag(kPtHigher);
218 if(ok) fNDigitEmcalCut++;
222 uArray->SetUnitDetectorFlag(kAll);
223 else uArray->SetUnitDetectorFlag(kEmcal);
225 // This is for jet multiplicity
226 uArray->SetUnitClusterID(index);
228 if(fDebug > 12) printf("goodDigit : %d\n", goodDigit);
230 } // End loop over digits
232 } // End loop over clusters
238 //_____________________________________________________________________________
239 Float_t AliJetFillUnitArrayEMCalDigits::EtaToTheta(Float_t arg)
241 // return (180./TMath::Pi())*2.*atan(exp(-arg));
242 return 2.*atan(exp(-arg));