AliPhysics  fe039ad (fe039ad)
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
AliAnalysisTaskEmcalSubjet.cxx
Go to the documentation of this file.
1 #include <AliVParticle.h>
2 
3 #include "THistManager.h"
4 
5 #include "AliAnalysisManager.h"
6 #include "AliVEventHandler.h"
7 
8 #include "AliJetContainer.h"
9 #include "AliEmcalJetFinder.h"
10 
12 //=============================================================================
13 
15 
16 //_____________________________________________________________________________
19 fSubjetRadius(0.),
20 fSubjetAlgorithm(0),
21 fSubjetFinder(nullptr),
22 fHistMgr(nullptr)
23 {
24 //
25 // AliAnalysisTaskEmcalSubjet::AliAnalysisTaskEmcalSubjet() :
26 //
27 }
28 
29 //_____________________________________________________________________________
31 AliAnalysisTaskEmcalJet(name,bHistos),
32 fSubjetRadius(0.1),
33 fSubjetAlgorithm(1),
34 fSubjetFinder(0),
35 fHistMgr(nullptr)
36 {
37 //
38 // AliAnalysisTaskEmcalSubjet::AliAnalysisTaskEmcalSubjet(const char *name, const Bool_t bHistos) :
39 //
40 
42 }
43 
44 //_____________________________________________________________________________
46 {
47 //
48 // AliAnalysisTaskEmcalSubjet::~AliAnalysisTaskEmcalSubjet
49 //
50 
51  if (fHistMgr) { delete fHistMgr; fHistMgr = nullptr; }
52  if (fSubjetFinder) { delete fSubjetFinder; fSubjetFinder = nullptr; }
53 }
54 
55 //_____________________________________________________________________________
57 {
58 //
59 // AliAnalysisTaskEmcalSubjet::UserCreateOutputObjects
60 //
61 
62  fCreateHisto = kTRUE;
64 //=============================================================================
65 
66  if (fSubjetFinder) { delete fSubjetFinder; fSubjetFinder = nullptr; }
67 
71 //=============================================================================
72 
73  if (fHistMgr) { delete fHistMgr; fHistMgr = nullptr; }
74  fHistMgr = new THistManager(GetName());
75 
79 //=============================================================================
80 
81  TObject *p(nullptr);
82  TIter next(fHistMgr->GetListOfHistograms());
83  while ((p = next())) fOutput->Add(p);
84 //=============================================================================
85 
86  return;
87 }
88 
89 //_____________________________________________________________________________
91 {
92 //
93 // Bool_t AliAnalysisTaskEmcalSubjet::Run()
94 //
95 
96  if (!AliAnalysisTaskEmcalJet::Run()) return kFALSE;
97 //=============================================================================
98 
99  TIter next(&fJetCollArray);
100  AliJetContainer *pc(nullptr);
101  while ((pc = static_cast<AliJetContainer*>(next()))) LoopJets(pc);
102 //=============================================================================
103 
104  return kTRUE;
105 }
106 
107 //_____________________________________________________________________________
109 {
110 //
111 // void AliAnalysisTaskEmcalSubjet::LoopJets(AliJetContainer const *pc)
112 //
113  if (!pc) return;
114 //=============================================================================
115 
116  const TString s(pc->GetName());
117  const auto d(pc->GetRhoName().IsNull() ? -1. : pc->GetRhoVal());
118 //=============================================================================
119 
120  UInt_t w(0);
121  for (auto pj : pc->accepted()) if (!(pj->IsGhost())) {
122 
123  w += 1;
124  const auto da(pj->Area());
125  const auto dj(d>0. ? pj->Pt() - d*da : pj->Pt());
126  fHistMgr->FillTH1(Form("%s/hJetPt_%d",s.Data(),fCentBin), dj);
127  fHistMgr->FillTH2(Form("%s/hArea_%d", s.Data(),fCentBin), dj, da);
128 //=============================================================================
129 
130  LoopJetConstis(pj, pc);
131  LoopSubjets(pj, pc);
132  }
133 //=============================================================================
134 
135  fHistMgr->FillTH1(Form("%s/hNjets_%d",s.Data(),fCentBin), w);
136 //=============================================================================
137 
138  return;
139 }
140 
141 //_____________________________________________________________________________
143 {
144 //
145 // void AliAnalysisTaskEmcalSubjet::LoopSubjets(AliEmcalJet const *pj, AliJetContainer const *pc)
146 //
147 
148  if (!(pj && pc)) return;
149  if (!(fSubjetFinder->Filter(const_cast<AliEmcalJet*>(pj),const_cast<AliJetContainer*>(pc),fVertex))) return;
150 //=============================================================================
151 
152  const TString s(Form("%s_sj",pc->GetName()));
153  const auto d(pc->GetRhoName().IsNull() ? -1. : pc->GetRhoVal());
154  const auto dj((d>0.) ? (pj->Pt() - d*(pj->Area())) : pj->Pt());
155 //=============================================================================
156 
157  UInt_t w(0);
158  auto d1(-9999.), d2(-9999.);
159  for (auto i=0; i<fSubjetFinder->GetNumberOfJets(); ++i) {
160  auto p(fSubjetFinder->GetJet(i)); if (!p) continue;
161 // if (p->IsGhost()) continue;
162 
163  w += 1;
164  const auto da(p->Area());
165  const auto dd(d>0. ? (p->Pt() - d*da) : p->Pt());
166  if (dd>d1) { d2 = d1; d1 = dd; } else if (dd>d2) { d2 = dd; }
167 //=============================================================================
168 
169  fHistMgr->FillTH2(Form("%s/hSubjetPt_%d", s.Data(),fCentBin), dd, dj);
170  fHistMgr->FillTH2(Form("%s/hSubjetArea_%d",s.Data(),fCentBin), dd, da);
171  }
172 //=============================================================================
173 
174  fHistMgr->FillTH1(Form("%s/hNsubjets_%d",s.Data(),fCentBin), w);
175  fHistMgr->FillTH2(Form("%s/hL1stSjPt_%d",s.Data(),fCentBin), d1, dj);
176  fHistMgr->FillTH2(Form("%s/hL2ndSjPt_%d",s.Data(),fCentBin), d2, dj);
177  fHistMgr->FillTH2(Form("%s/hDelta_Sj_%d",s.Data(),fCentBin), d1 - d2, dj);
178 //=============================================================================
179 
180  return;
181 }
182 
183 //_____________________________________________________________________________
185 {
186 //
187 // void AliAnalysisTaskEmcalSubjet::LoopJetConstis(AliEmcalJet const *pj, AliJetContainer const *pc)
188 //
189 
190  if (!(pj && pc)) return;
191 //=============================================================================
192 
193  const TString s(Form("%s_consti",pc->GetName()));
194  const auto d(pc->GetRhoName().IsNull() ? -1. : pc->GetRhoVal());
195  const auto dj((d>0.) ? (pj->Pt() - d*(pj->Area())) : pj->Pt());
196 //=============================================================================
197 
198  UInt_t w(0);
199  auto d1(-9999.), d2(-9999.);
200  for (auto i=0; i<pj->GetNumberOfTracks(); ++i) {
201  const auto p(pj->Track(i)); if (!p) continue;
202 
203  w += 1;
204  const auto dd(p->Pt());
205  if (dd>d1) { d2 = d1; d1 = dd; } else if (dd>d2) { d2 = dd; }
206  }
207 //=============================================================================
208 
209 /*TLorentzVector v;
210  for (auto i=0; i<pj->GetNumberOfClusters(); ++i) {
211  auto p(pj->Cluster(i)); if (!p) continue;
212  p->GetMomentum(v,fVertex);
213 
214  w += 1;
215  auto dd(v.Pt());
216  if (dd>d1) { d2 = d1; d1 = dd; } else if (dd>d2) { d2 = dd; }
217  }*/
218 //=============================================================================
219 
220  fHistMgr->FillTH1(Form("%s/hNtracks_%d",s.Data(),fCentBin), w);
221  fHistMgr->FillTH2(Form("%s/hL1stTrk_%d",s.Data(),fCentBin), d1, dj);
222  fHistMgr->FillTH2(Form("%s/hL2ndTrk_%d",s.Data(),fCentBin), d2, dj);
223  fHistMgr->FillTH2(Form("%s/hDeltaTk_%d",s.Data(),fCentBin), d1-d2, dj);
224 //=============================================================================
225 
226  return;
227 }
228 
229 //_____________________________________________________________________________
231 {
232 //
233 // void AliAnalysisTaskEmcalSubjet::CreateHistoJets()
234 //
235 
236  TIter next(&fJetCollArray);
237  AliJetContainer *p(nullptr);
238  while ((p = static_cast<AliJetContainer*>(next()))) {
239  const TString s(p->GetName());
240 
241  if (fHistMgr->FindObject(s)) {
242  AliWarning(Form("Found group name %s in hist manager", s.Data()));
243  continue;
244  }
245 //=============================================================================
246 
248  const auto bRho(!(p->GetRhoName().IsNull()));
249  const auto dMin(bRho ? -fMaxBinPt/2. : fMinBinPt);
250  const auto dMax(bRho ? fMaxBinPt/2. : fMaxBinPt);
251 //=============================================================================
252 
253  const auto nbs(2*fNbins);
254  for (auto i=0; i<fNcentBins; ++i) {
255  fHistMgr->CreateTH1(Form("%s/hNjets_%d",s.Data(),i), "", 500, -0.5, 499.5, "s");
256  fHistMgr->CreateTH1(Form("%s/hJetPt_%d",s.Data(),i), "", fNbins, dMin, dMax, "s");
257  fHistMgr->CreateTH2(Form("%s/hArea_%d", s.Data(),i), "", fNbins, dMin, dMax, nbs, 0., 3., "s");
258  }
259  }
260 //=============================================================================
261 
262  return;
263 }
264 
265 //_____________________________________________________________________________
267 {
268 //
269 // void AliAnalysisTaskEmcalSubjet::CreateHistoSubjets()
270 //
271 
272  TIter next(&fJetCollArray);
273  AliJetContainer *p(nullptr);
274  while ((p = static_cast<AliJetContainer*>(next()))) {
275  const TString s(Form("%s_sj",p->GetName()));
276 
277  if (fHistMgr->FindObject(s)) {
278  AliWarning(Form("Found group name %s in hist manager", s.Data()));
279  continue;
280  }
281 //=============================================================================
282 
284  const auto bRho(!(p->GetRhoName().IsNull()));
285  const auto dMin(bRho ? -fMaxBinPt/2. : fMinBinPt);
286  const auto dMax(bRho ? fMaxBinPt/2. : fMaxBinPt);
287 //=============================================================================
288 
289  const auto nbs(2*fNbins);
290  for (auto i=0; i<fNcentBins; ++i) {
291  fHistMgr->CreateTH1(Form("%s/hNsubjets_%d", s.Data(),i), "", 500, -0.5, 499.5, "s");
292  fHistMgr->CreateTH2(Form("%s/hSubjetPt_%d", s.Data(),i), "", nbs, dMin, dMax, fNbins, dMin, dMax, "s");
293  fHistMgr->CreateTH2(Form("%s/hL1stSjPt_%d", s.Data(),i), "", nbs, dMin, dMax, fNbins, dMin, dMax, "s");
294  fHistMgr->CreateTH2(Form("%s/hL2ndSjPt_%d", s.Data(),i), "", nbs, dMin, dMax, fNbins, dMin, dMax, "s");
295  fHistMgr->CreateTH2(Form("%s/hDelta_Sj_%d", s.Data(),i), "", nbs, dMin, dMax, fNbins, dMin, dMax, "s");
296  fHistMgr->CreateTH2(Form("%s/hSubjetArea_%d",s.Data(),i), "", nbs, dMin, dMax, nbs, 0., 3., "s");
297  }
298  }
299 //=============================================================================
300 
301  return;
302 }
303 
304 //_____________________________________________________________________________
306 {
307 //
308 // void AliAnalysisTaskEmcalSubjet::CreateHistoJetConstis()
309 //
310 
311  TIter next(&fJetCollArray);
312  AliJetContainer *p(nullptr);
313  while ((p = static_cast<AliJetContainer*>(next()))) {
314  const TString s(Form("%s_consti",p->GetName()));
315 
316  if (fHistMgr->FindObject(s)) {
317  AliWarning(Form("Found group name %s in hist manager", s.Data()));
318  continue;
319  }
320 //=============================================================================
321 
323  const auto bRho(!(p->GetRhoName().IsNull()));
324  const auto dMin(bRho ? -fMaxBinPt/2. : fMinBinPt);
325  const auto dMax(bRho ? fMaxBinPt/2. : fMaxBinPt);
326 //=============================================================================
327 
328  const auto nbs(2*fNbins);
329  for (auto i=0; i<fNcentBins; ++i) {
330  fHistMgr->CreateTH1(Form("%s/hNtracks_%d",s.Data(),i), "", 500, -0.5, 499.5, "s");
331  fHistMgr->CreateTH2(Form("%s/hL1stTrk_%d",s.Data(),i), "", nbs, dMin, dMax, fNbins, dMin, dMax, "s");
332  fHistMgr->CreateTH2(Form("%s/hL2ndTrk_%d",s.Data(),i), "", nbs, dMin, dMax, fNbins, dMin, dMax, "s");
333  fHistMgr->CreateTH2(Form("%s/hDeltaTk_%d",s.Data(),i), "", nbs, dMin, dMax, fNbins, dMin, dMax, "s");
334  }
335  }
336 //=============================================================================
337 
338  return;
339 }
340 
341 //_____________________________________________________________________________
343  const TString sClus,
344  const TString sCells)
345 {
346 //
347 // AliAnalysisTaskEmcalSubjet *AliAnalysisTaskEmcalSubjet::AddTask()
348 //
349 
350  auto mgr(AliAnalysisManager::GetAnalysisManager());
351 
352  if (!mgr) {
353  ::Error(Form("AliAnalysisTaskEmcalSubjet::%s",__func__), "No analysis manager to connect to");
354  return nullptr;
355  }
356 //=============================================================================
357 
358  AliVEventHandler *pH(mgr->GetInputEventHandler());
359 
360  if (!pH) {
361  ::Error(Form("AliAnalysisTaskEmcalSubjet::%s",__func__), "This task requires an input event handler");
362  return nullptr;
363  }
364 //=============================================================================
365 
366  enum EDataType_t { kUnknown, kESD, kAOD };
367 
368  EDataType_t wType(kUnknown);
369  if (pH->InheritsFrom("AliESDInputHandler")) wType = kESD;
370  if (pH->InheritsFrom("AliAODInputHandler")) wType = kAOD;
371 
372  if (wType==kUnknown) {
373  ::Error(Form("AliAnalysisTaskEmcalSubjet::%s",__func__), "Unkown data input");
374  return nullptr;
375  }
376 //=============================================================================
377 
378  auto task(new AliAnalysisTaskEmcalSubjet("AliAnalysisTaskEmcalSubjet",kTRUE));
379 //=============================================================================
380 
381  TString sTrkName(sTrks);
382  if (sTrks=="usedefault") {
383  if (wType==kESD) sTrkName = "Tracks";
384  if (wType==kAOD) sTrkName = "tracks";
385  }
386 
387  if (sTrkName=="mcparticles") {
388  task->AddMCParticleContainer(sTrkName);
389  } else if ((sTrkName=="tracks") || (sTrkName=="Tracks")) {
390  task->AddTrackContainer(sTrkName);
391  } else if (!sTrkName.IsNull()) {
392  task->AddParticleContainer(sTrkName);
393  }
394 //=============================================================================
395 
396  TString sClsName(sClus);
397  if (sClus=="usedefault") {
398  if (wType==kESD) sClsName = "CaloClusters";
399  if (wType==kAOD) sClsName = "caloClusters";
400  }
401 
402  task->AddClusterContainer(sClsName);
403 //=============================================================================
404 
405  TString sCellName(sCells);
406  if (sCells=="usedefault") {
407  if (wType==kESD) sCellName = "EMCALCells";
408  if (wType==kAOD) sCellName = "emcalCells";
409  }
410 
411  task->SetCaloCellsName(sCellName);
412 //=============================================================================
413 
414  mgr->AddTask(task);
415  mgr->ConnectInput( task, 0, mgr->GetCommonInputContainer());
416  mgr->ConnectOutput(task, 1, mgr->CreateContainer("listSubjet",
417  AliEmcalList::Class(),
418  AliAnalysisManager::kOutputContainer,
419  AliAnalysisManager::GetCommonFileName()));
420 //=============================================================================
421 
422  return task;
423 }
void SetRadius(Double_t val)
virtual 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.
Double_t Area() const
Definition: AliEmcalJet.h:130
THashList * CreateHistoGroup(const char *groupname)
Create a new group of histograms within a parent group.
Double_t GetRhoVal() const
const TString & GetRhoName() const
void LoopJetConstis(AliEmcalJet const *pj, AliJetContainer const *pc)
Bool_t Filter(AliEmcalJet *pJet, AliJetContainer *pContJets, Double_t dVtx[3])
void FillTH2(const char *hname, double x, double y, double weight=1., Option_t *opt="")
Fill a 2D histogram within the container.
Double_t fMinBinPt
min pt in histograms
Bool_t fGeneralHistograms
whether or not it should fill some general histograms
AliVParticle * Track(Int_t idx) const
Int_t fCentBin
!event centrality bin
UShort_t GetNumberOfTracks() const
Definition: AliEmcalJet.h:139
void LoopSubjets(AliEmcalJet const *pj, AliJetContainer const *pc)
TH2 * CreateTH2(const char *name, const char *title, int nbinsx, double xmin, double xmax, int nbinsy, double ymin, double ymax, Option_t *opt="")
Create a new TH2 within the container.
void SetJetAlgorithm(Int_t val)
TObject * FindObject(const char *name) const
Find an object inside the container.
unsigned int UInt_t
Definition: External.C:33
THashList * GetListOfHistograms() const
Get the list of histograms.
Definition: THistManager.h:671
void LoopJets(AliJetContainer const *pc)
TH1 * CreateTH1(const char *name, const char *title, int nbins, double xmin, double xmax, Option_t *opt="")
Create a new TH1 within the container.
Int_t fNcentBins
how many centrality bins
void FillTH1(const char *hname, double x, double weight=1., Option_t *opt="")
Fill a 1D histogram within the container.
TObjArray fJetCollArray
jet collection array
AliEmcalJet * GetJet(Int_t index)
Double_t Pt() const
Definition: AliEmcalJet.h:109
AliEmcalList * fOutput
!output list
Double_t fMaxBinPt
max pt in histograms
Double_t fVertex[3]
!event vertex
static AliAnalysisTaskEmcalSubjet * AddTask(const TString sTrks="usedefault", const TString sClus="usedefault", const TString sCells="usedefault")
Bool_t fCreateHisto
whether or not create histograms
virtual 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.
Base task in the EMCAL jet framework.
Represent a jet reconstructed using the EMCal jet framework.
Definition: AliEmcalJet.h:51
Container class for histograms.
Definition: THistManager.h:99
void UserCreateOutputObjects()
Main initialization function on the worker.
const AliJetIterableContainer accepted() const
bool Bool_t
Definition: External.C:53
EDataType_t
Switch for the data type.
Container for jet within the EMCAL jet framework.
Int_t fNbins
no. of pt bins