AliPhysics  608b256 (608b256)
AliAnalysisTaskRho.cxx
Go to the documentation of this file.
1 // $Id$
2 //
3 // Calculation of rho from a collection of jets.
4 // If scale function is given the scaled rho will be exported
5 // with the name as "fOutRhoName".Apppend("_Scaled").
6 //
7 // Authors: R.Reed, S.Aiola
8 
9 #include "AliAnalysisTaskRho.h"
10 
11 #include <TClonesArray.h>
12 #include <TMath.h>
13 
14 #include "AliAnalysisManager.h"
15 #include "AliEmcalJet.h"
16 #include "AliLog.h"
17 #include "AliRhoParameter.h"
18 #include "AliJetContainer.h"
19 #include "AliParticleContainer.h"
20 #include "AliClusterContainer.h"
21 #include "AliVEventHandler.h"
22 #include "AliAnalysisDataContainer.h"
23 
24 ClassImp(AliAnalysisTaskRho)
25 
26 //________________________________________________________________________
29  fNExclLeadJets(0)
30 {
31  // Constructor.
32 }
33 
34 //________________________________________________________________________
36  AliAnalysisTaskRhoBase(name, histo),
37  fNExclLeadJets(0)
38 {
39  // Constructor.
40 }
41 
42 
43 //________________________________________________________________________
45 {
46  // Run the analysis.
47 
48  fOutRho->SetVal(0);
49  if (fOutRhoScaled)
50  fOutRhoScaled->SetVal(0);
51 
52  if (!fJets)
53  return kFALSE;
54 
55  const Int_t Njets = fJets->GetEntries();
56 
57  Int_t maxJetIds[] = {-1, -1};
58  Float_t maxJetPts[] = { 0, 0};
59 
60  if (fNExclLeadJets > 0) {
61  for (Int_t ij = 0; ij < Njets; ++ij) {
62  AliEmcalJet *jet = static_cast<AliEmcalJet*>(fJets->At(ij));
63  if (!jet) {
64  AliError(Form("%s: Could not receive jet %d", GetName(), ij));
65  continue;
66  }
67 
68  if (!AcceptJet(jet))
69  continue;
70 
71  if (jet->Pt() > maxJetPts[0]) {
72  maxJetPts[1] = maxJetPts[0];
73  maxJetIds[1] = maxJetIds[0];
74  maxJetPts[0] = jet->Pt();
75  maxJetIds[0] = ij;
76  } else if (jet->Pt() > maxJetPts[1]) {
77  maxJetPts[1] = jet->Pt();
78  maxJetIds[1] = ij;
79  }
80  }
81  if (fNExclLeadJets < 2) {
82  maxJetIds[1] = -1;
83  maxJetPts[1] = 0;
84  }
85  }
86 
87  static Double_t rhovec[999];
88  Int_t NjetAcc = 0;
89 
90  // push all jets within selected acceptance into stack
91  for (Int_t iJets = 0; iJets < Njets; ++iJets) {
92 
93  // exlcuding lead jets
94  if (iJets == maxJetIds[0] || iJets == maxJetIds[1])
95  continue;
96 
97  AliEmcalJet *jet = static_cast<AliEmcalJet*>(fJets->At(iJets));
98  if (!jet) {
99  AliError(Form("%s: Could not receive jet %d", GetName(), iJets));
100  continue;
101  }
102 
103  if (!AcceptJet(jet))
104  continue;
105 
106  rhovec[NjetAcc] = jet->Pt() / jet->Area();
107  ++NjetAcc;
108  }
109 
110 
111  if (NjetAcc > 0) {
112  //find median value
113  Double_t rho = TMath::Median(NjetAcc, rhovec);
114  fOutRho->SetVal(rho);
115 
116  if (fOutRhoScaled) {
117  Double_t rhoScaled = rho * GetScaleFactor(fCent);
118  fOutRhoScaled->SetVal(rhoScaled);
119  }
120  }
121 
122  return kTRUE;
123 }
124 
125 //________________________________________________________________________
127  const char* nTracks, const char* nClusters, const char* nRho,
128  Double_t jetradius, UInt_t acceptance, AliJetContainer::EJetType_t jetType, const Bool_t histo,
129  AliJetContainer::ERecoScheme_t rscheme, const char* suffix
130 )
131 {
132  // Get the pointer to the existing analysis manager via the static access method.
133  //==============================================================================
134  AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
135  if (!mgr) { ::Error("AddTaskRho", "No analysis manager to connect to."); return NULL; }
136 
137  // Check the analysis type using the event handlers connected to the analysis manager.
138  //==============================================================================
139  AliVEventHandler* handler = mgr->GetInputEventHandler();
140  if (!handler) { ::Error("AddTaskEmcalJetSpectraQA", "This task requires an input event handler"); return NULL; }
141 
142  enum EDataType_t {
143  kUnknown,
144  kESD,
145  kAOD
146  };
147 
148  EDataType_t dataType = kUnknown;
149 
150  if (handler->InheritsFrom("AliESDInputHandler")) {
151  dataType = kESD;
152  }
153  else if (handler->InheritsFrom("AliAODInputHandler")) {
154  dataType = kAOD;
155  }
156 
157  //-------------------------------------------------------
158  // Init the task and do settings
159  //-------------------------------------------------------
160 
161  TString trackName(nTracks);
162  TString clusName(nClusters);
163 
164  if (trackName == "usedefault") {
165  if (dataType == kESD) {
166  trackName = "Tracks";
167  }
168  else if (dataType == kAOD) {
169  trackName = "tracks";
170  }
171  else {
172  trackName = "";
173  }
174  }
175 
176  if (clusName == "usedefault") {
177  if (dataType == kESD) {
178  clusName = "CaloClusters";
179  }
180  else if (dataType == kAOD) {
181  clusName = "caloClusters";
182  }
183  else {
184  clusName = "";
185  }
186  }
187 
188  TString name("AliAnalysisTaskRho");
189  if (strcmp(suffix,"") != 0) {
190  name += "_";
191  name += suffix;
192  }
193 
194  AliAnalysisTaskRho* mgrTask = (AliAnalysisTaskRho*)(mgr->GetTask(name.Data()));
195  if (mgrTask) return mgrTask;
196 
197  AliAnalysisTaskRho *rhotask = new AliAnalysisTaskRho(name, histo);
198  rhotask->SetOutRhoName(nRho);
199 
200  if (trackName == "mcparticles") {
201  AliMCParticleContainer* mcpartCont = rhotask->AddMCParticleContainer(trackName);
202  mcpartCont->SelectPhysicalPrimaries(kTRUE);
203  }
204  else if (trackName == "tracks" || trackName == "Tracks") {
205  AliTrackContainer* trackCont = rhotask->AddTrackContainer(trackName);
206  }
207  else if (!trackName.IsNull()) {
208  rhotask->AddParticleContainer(trackName);
209  }
210  AliParticleContainer* partCont = rhotask->GetParticleContainer(0);
211 
212  AliClusterContainer *clusterCont = rhotask->AddClusterContainer(clusName);
213  if (clusterCont) {
214  clusterCont->SetClusECut(0.);
215  clusterCont->SetClusPtCut(0.);
216  clusterCont->SetClusHadCorrEnergyCut(0.3);
217  clusterCont->SetDefaultClusterEnergy(AliVCluster::kHadCorr);
218  }
219 
220  AliJetContainer *jetCont = rhotask->AddJetContainer(jetType, AliJetContainer::kt_algorithm, rscheme, jetradius, acceptance, partCont, clusterCont);
221  if (jetCont) jetCont->SetJetPtCut(0);
222 
223  //-------------------------------------------------------
224  // Final settings, pass to manager and set the containers
225  //-------------------------------------------------------
226 
227  mgr->AddTask(rhotask);
228 
229  // Create containers for input/output
230  mgr->ConnectInput(rhotask, 0, mgr->GetCommonInputContainer());
231  if (histo) {
232  TString contname(name);
233  contname += "_histos";
234  AliAnalysisDataContainer *coutput1 = mgr->CreateContainer(contname.Data(),
235  TList::Class(), AliAnalysisManager::kOutputContainer,
236  Form("%s", AliAnalysisManager::GetCommonFileName()));
237  mgr->ConnectOutput(rhotask, 1, coutput1);
238  }
239 
240  return rhotask;
241 }
242 
Double_t Area() const
Definition: AliEmcalJet.h:130
double Double_t
Definition: External.C:58
void SetOutRhoName(const char *name)
Container with name, TClonesArray and cuts for particles.
AliJetContainer * AddJetContainer(const char *n, TString defaultCutType, Float_t jetRadius=0.4)
AliClusterContainer * AddClusterContainer(const char *n)
Create new cluster container and attach it to the task.
Container for particles within the EMCAL framework.
static AliAnalysisTaskRho * AddTaskRhoNew(const char *nTracks="usedefault", const char *nClusters="usedefault", const char *nRho="Rho", Double_t jetradius=0.2, UInt_t acceptance=AliEmcalJet::kTPCfid, AliJetContainer::EJetType_t jetType=AliJetContainer::kChargedJet, const Bool_t histo=kFALSE, AliJetContainer::ERecoScheme_t rscheme=AliJetContainer::pt_scheme, const char *suffix="")
AliParticleContainer * GetParticleContainer(Int_t i=0) const
Get particle container attached to this task.
AliRhoParameter * fOutRhoScaled
output rho object
int Int_t
Definition: External.C:63
void SetJetPtCut(Float_t cut)
unsigned int UInt_t
Definition: External.C:33
float Float_t
Definition: External.C:68
AliParticleContainer * AddParticleContainer(const char *n)
Create new particle container and attach it to the task.
Double_t fCent
!event centrality
AliMCParticleContainer * AddMCParticleContainer(const char *n)
Create new container for MC particles and attach it to the task.
TClonesArray * fJets
! jets
Double_t Pt() const
Definition: AliEmcalJet.h:109
virtual Double_t GetScaleFactor(Double_t cent)
AliTrackContainer * AddTrackContainer(const char *n)
Create new track container and attach it to the task.
void SelectPhysicalPrimaries(Bool_t s)
Represent a jet reconstructed using the EMCal jet framework.
Definition: AliEmcalJet.h:51
void SetClusPtCut(Double_t cut)
Bool_t Run()
Run function. This is the core function of the analysis and contains the user code. Therefore users have to implement this function.
virtual Bool_t AcceptJet(AliEmcalJet *jet, Int_t c=0)
bool Bool_t
Definition: External.C:53
EDataType_t
Switch for the data type.
void SetClusECut(Double_t cut)
void SetDefaultClusterEnergy(Int_t d)
Container structure for EMCAL clusters.
Container for MC-true particles within the EMCAL framework.
Container for jet within the EMCAL jet framework.
void SetClusHadCorrEnergyCut(Double_t cut)