AliPhysics  master (3d17d9d)
AliAnalysisTaskRhoBaseDev.cxx
Go to the documentation of this file.
1 /************************************************************************************
2  * Copyright (C) 2017, Copyright Holders of the ALICE Collaboration *
3  * All rights reserved. *
4  * *
5  * Redistribution and use in source and binary forms, with or without *
6  * modification, are permitted provided that the following conditions are met: *
7  * * Redistributions of source code must retain the above copyright *
8  * notice, this list of conditions and the following disclaimer. *
9  * * Redistributions in binary form must reproduce the above copyright *
10  * notice, this list of conditions and the following disclaimer in the *
11  * documentation and/or other materials provided with the distribution. *
12  * * Neither the name of the <organization> nor the *
13  * names of its contributors may be used to endorse or promote products *
14  * derived from this software without specific prior written permission. *
15  * *
16  * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND *
17  * ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED *
18  * WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE *
19  * DISCLAIMED. IN NO EVENT SHALL ALICE COLLABORATION BE LIABLE FOR ANY *
20  * DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES *
21  * (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; *
22  * LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND *
23  * ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT *
24  * (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS *
25  * SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. *
26  ************************************************************************************/
27 #include <TFile.h>
28 #include <TF1.h>
29 #include <TH1F.h>
30 #include <TH2F.h>
31 #include <TClonesArray.h>
32 #include <TGrid.h>
33 
34 #include <AliLog.h>
35 #include <AliVEventHandler.h>
36 #include <AliAnalysisManager.h>
37 
38 #include "AliRhoParameter.h"
39 #include "AliEmcalJet.h"
40 #include "AliParticleContainer.h"
41 #include "AliClusterContainer.h"
42 #include "AliJetContainer.h"
43 
45 
47 
50  fOutRhoName(),
51  fOutRhoScaledName(),
52  fRhoFunction(0),
53  fScaleFunction(0),
54  fAttachToEvent(kTRUE),
55  fTaskConfigured(kFALSE),
56  fOutRho(0),
57  fOutRhoScaled(0),
58  fHistRhoVsCent(0),
59  fHistRhoVsLeadJetPt(),
60  fHistLeadJetPtVsCent(),
61  fHistLeadJetPtDensityVsCent(),
62  fHistTotJetAreaVsCent(),
63  fHistLeadJetNconstVsCent(),
64  fHistLeadJetNconstVsPt(),
65  fHistNjetVsCent(),
66  fHistNjetVsNtrack(),
67  fHistRhoVsLeadTrackPt(0),
68  fHistRhoVsNtrack(0),
69  fHistLeadTrackPtVsCent(0),
70  fHistNtrackVsCent(0),
71  fHistRhoVsLeadClusterE(0),
72  fHistRhoVsNcluster(0),
73  fHistLeadClusterEVsCent(0),
74  fHistNclusterVsCent(0),
75  fHistRhoScaledVsCent(0),
76  fHistRhoScaledVsNtrack(0),
77  fHistRhoScaledVsNcluster(0)
78 {
79 }
80 
82  AliAnalysisTaskJetUE(name, histo),
83  fOutRhoName(),
85  fRhoFunction(0),
86  fScaleFunction(0),
87  fAttachToEvent(kTRUE),
88  fTaskConfigured(kFALSE),
89  fOutRho(0),
90  fOutRhoScaled(0),
91  fHistRhoVsCent(0),
101  fHistRhoVsNtrack(0),
111 {
113 }
114 
116 {
117  if (!fCreateHisto) return;
118 
120 
121  ::Info("UserCreateOutputObjects", "CreateOutputObjects of task %s", GetName());
122 
123  TString name;
124 
125  Int_t maxTracks = 6000;
126  Double_t maxRho = 500;
127  Int_t nRhoBins = 500;
128 
129  if (fForceBeamType == kpp) {
130  maxRho = 50;
131  maxTracks = 200;
132  }
133  else if (fForceBeamType == kpA) {
134  maxRho = 200;
135  maxTracks = 500;
136  }
137 
138  Int_t nPtBins = TMath::CeilNint(fMaxPt / fPtBinWidth);
139 
140  fHistRhoVsCent = new TH2F("fHistRhoVsCent", "fHistRhoVsCent", 100, 0, 100, nRhoBins, 0, maxRho);
141  fHistRhoVsCent->GetXaxis()->SetTitle("Centrality (%)");
142  fHistRhoVsCent->GetYaxis()->SetTitle("#rho (GeV/#it{c} #times rad^{-1})");
143  fOutput->Add(fHistRhoVsCent);
144 
145  if (fParticleCollArray.size() > 0) {
146  fHistRhoVsNtrack = new TH2F("fHistRhoVsNtrack", "fHistRhoVsNtrack", 200, 0, maxTracks, nRhoBins, 0, maxRho);
147  fHistRhoVsNtrack->GetXaxis()->SetTitle("No. of tracks");
148  fHistRhoVsNtrack->GetYaxis()->SetTitle("#rho (GeV/#it{c} #times rad^{-1})");
150 
151  fHistNtrackVsCent = new TH2F("fHistNtrackVsCent", "fHistNtrackVsCent", 100, 0, 100, 200, 0, maxTracks);
152  fHistNtrackVsCent->GetXaxis()->SetTitle("Centrality (%)");
153  fHistNtrackVsCent->GetYaxis()->SetTitle("No. of tracks");
155 
156  fHistRhoVsLeadTrackPt = new TH2F("fHistRhoVsLeadTrackPt", "fHistRhoVsLeadTrackPt", nPtBins, 0, fMaxPt, nRhoBins, 0, maxRho);
157  fHistRhoVsLeadTrackPt->GetXaxis()->SetTitle("#it{p}_{T,track} (GeV/c)");
158  fHistRhoVsLeadTrackPt->GetYaxis()->SetTitle("#rho (GeV/#it{c} #times rad^{-1})");
160 
161  fHistLeadTrackPtVsCent = new TH2F("fHistLeadTrackPtVsCent", "fHistLeadTrackPtVsCent", 100, 0, 100, nPtBins, 0, fMaxPt);
162  fHistLeadTrackPtVsCent->GetXaxis()->SetTitle("Centrality (%)");
163  fHistLeadTrackPtVsCent->GetYaxis()->SetTitle("#it{p}_{T,track} (GeV/c)");
165  }
166 
167  if (fClusterCollArray.size()>0) {
168  fHistRhoVsNcluster = new TH2F("fHistRhoVsNcluster", "fHistRhoVsNcluster", 50, 0, maxTracks / 4, nRhoBins, 0, maxRho);
169  fHistRhoVsNcluster->GetXaxis()->SetTitle("No. of clusters");
170  fHistRhoVsNcluster->GetYaxis()->SetTitle("#rho (GeV/#it{c} #times rad^{-1})");
172 
173  fHistNclusterVsCent = new TH2F("fHistNclusterVsCent", "fHistNclusterVsCent", 100, 0, 100, 50, 0, maxTracks / 4);
174  fHistNclusterVsCent->GetXaxis()->SetTitle("Centrality (%)");
175  fHistNclusterVsCent->GetYaxis()->SetTitle("No. of clusters");
177 
178  fHistRhoVsLeadClusterE = new TH2F("fHistRhoVsLeadClusterE", "fHistRhoVsLeadClusterE", nPtBins, 0, fMaxPt, nRhoBins, 0, maxRho);
179  fHistRhoVsLeadClusterE->GetXaxis()->SetTitle("#it{p}_{T,track} (GeV/c)");
180  fHistRhoVsLeadClusterE->GetYaxis()->SetTitle("#rho (GeV/#it{c} #times rad^{-1})");
182 
183  fHistLeadClusterEVsCent = new TH2F("fHistLeadClusterEVsCent", "fHistLeadClusterEVsCent", 100, 0, 100, nPtBins, 0, fMaxPt);
184  fHistLeadClusterEVsCent->GetXaxis()->SetTitle("Centrality (%)");
185  fHistLeadClusterEVsCent->GetYaxis()->SetTitle("#it{p}_{T,track} (GeV/c)");
187  }
188 
189  for (auto jetCont : fJetCollArray) {
190  name = TString::Format("%s_fHistRhoVsLeadJetPt", jetCont.first.c_str());
191  fHistRhoVsLeadJetPt[jetCont.first] = new TH2F(name, name, nPtBins, 0, fMaxPt, nRhoBins, 0, maxRho);
192  fHistRhoVsLeadJetPt[jetCont.first]->GetXaxis()->SetTitle("#it{p}_{T,jet} (GeV/c)");
193  fHistRhoVsLeadJetPt[jetCont.first]->GetYaxis()->SetTitle("#rho (GeV/#it{c} #times rad^{-1})");
194  fOutput->Add(fHistRhoVsLeadJetPt[jetCont.first]);
195 
196  name = TString::Format("%s_fHistLeadJetPtVsCent", jetCont.first.c_str());
197  fHistLeadJetPtVsCent[jetCont.first] = new TH2F(name, name, 100, 0, 100, nPtBins, 0, fMaxPt);
198  fHistLeadJetPtVsCent[jetCont.first]->GetXaxis()->SetTitle("Centrality (%)");
199  fHistLeadJetPtVsCent[jetCont.first]->GetYaxis()->SetTitle("#it{p}_{T,jet} (GeV/c)");
200  fOutput->Add(fHistLeadJetPtVsCent[jetCont.first]);
201 
202  name = TString::Format("%s_fHistLeadJetPtDensityVsCent", jetCont.first.c_str());
203  fHistLeadJetPtDensityVsCent[jetCont.first] = new TH2F(name, name, 100, 0, 100, nPtBins, 0, fMaxPt*2);
204  fHistLeadJetPtDensityVsCent[jetCont.first]->GetXaxis()->SetTitle("Centrality (%)");
205  fHistLeadJetPtDensityVsCent[jetCont.first]->GetYaxis()->SetTitle("#it{p}_{T,jet} / #it{A}_{jet} (GeV/#it{c})");
206  fOutput->Add(fHistLeadJetPtDensityVsCent[jetCont.first]);
207 
208  name = TString::Format("%s_fHistLeadJetNconstVsCent", jetCont.first.c_str());
209  fHistLeadJetNconstVsCent[jetCont.first] = new TH2F(name, name, 100, 0, 100, 150, -0.5, 149.5);
210  fHistLeadJetNconstVsCent[jetCont.first]->GetXaxis()->SetTitle("Centrality (%)");
211  fHistLeadJetNconstVsCent[jetCont.first]->GetYaxis()->SetTitle("No. of constituents");
212  fOutput->Add(fHistLeadJetNconstVsCent[jetCont.first]);
213 
214  if (fCentBins.size() > 1) {
215  fHistLeadJetNconstVsPt[jetCont.first] = new TH2*[fCentBins.size()-1];
216  for (Int_t i = 0; i < fCentBins.size()-1; i++) {
217  name = TString::Format("%s_fHistJetNconstVsPt_Cent%d_%d", jetCont.first.c_str(), TMath::FloorNint(fCentBins[i]), TMath::FloorNint(fCentBins[i+1]));
218  fHistLeadJetNconstVsPt[jetCont.first][i] = new TH2F(name, name, nPtBins, 0, fMaxPt, 150, -0.5, 149.5);
219  fHistLeadJetNconstVsPt[jetCont.first][i]->GetXaxis()->SetTitle("#it{p}_{T,jet} (GeV/#it{c})");
220  fHistLeadJetNconstVsPt[jetCont.first][i]->GetYaxis()->SetTitle("No. of constituents");
221  fOutput->Add(fHistLeadJetNconstVsPt[jetCont.first][i]);
222  }
223  }
224 
225  name = TString::Format("%s_fHistTotJetAreaVsCent", jetCont.first.c_str());
226  fHistTotJetAreaVsCent[jetCont.first] = new TH2F(name, name, 100, 0, 100, 500, 0, 15);
227  fHistTotJetAreaVsCent[jetCont.first]->GetXaxis()->SetTitle("Centrality (%)");
228  fHistTotJetAreaVsCent[jetCont.first]->GetYaxis()->SetTitle("Jet area");
229  fOutput->Add(fHistTotJetAreaVsCent[jetCont.first]);
230 
231  name = TString::Format("%s_fHistNjetVsCent", jetCont.first.c_str());
232  fHistNjetVsCent[jetCont.first] = new TH2F(name, name, 100, 0, 100, 150, -0.5, 149.5);
233  fHistNjetVsCent[jetCont.first]->GetXaxis()->SetTitle("Centrality (%)");
234  fHistNjetVsCent[jetCont.first]->GetYaxis()->SetTitle("No. of jets");
235  fOutput->Add(fHistNjetVsCent[jetCont.first]);
236 
237  if (fParticleCollArray.size() > 0) {
238  name = TString::Format("%s_fHistNjetVsNtrack", jetCont.first.c_str());
239  fHistNjetVsNtrack[jetCont.first] = new TH2F(name, name, 200, 0, maxTracks, 150, -0.5, 149.5);
240  fHistNjetVsNtrack[jetCont.first]->GetXaxis()->SetTitle("No. of tracks");
241  fHistNjetVsNtrack[jetCont.first]->GetYaxis()->SetTitle("No. of jets");
242  fOutput->Add(fHistNjetVsNtrack[jetCont.first]);
243  }
244  }
245 
246  if (fScaleFunction) {
247  fHistRhoScaledVsCent = new TH2F("fHistRhoScaledVsCent", "fHistRhoScaledVsCent", 100, 0, 100, nRhoBins, 0, maxRho);
248  fHistRhoScaledVsCent->GetXaxis()->SetTitle("Centrality (%)");
249  fHistRhoScaledVsCent->GetYaxis()->SetTitle("#rho_{scaled} (GeV/#it{c} #times rad^{-1})");
251 
252  if (fParticleCollArray.size() > 0) {
253  fHistRhoScaledVsNtrack = new TH2F("fHistRhoScaledVsNtrack", "fHistRhoScaledVsNtrack", 200, 0, maxTracks, nRhoBins, 0, maxRho);
254  fHistRhoScaledVsNtrack->GetXaxis()->SetTitle("No. of tracks");
255  fHistRhoScaledVsNtrack->GetYaxis()->SetTitle("#rho (GeV/#it{c} #times rad^{-1})");
257  }
258 
259  if (fClusterCollArray.size() > 0) {
260  fHistRhoScaledVsNcluster = new TH2F("fHistRhoScaledVsNcluster", "fHistRhoScaledVsNcluster", 50, 0, maxTracks / 4, nRhoBins, 0, maxRho);
261  fHistRhoScaledVsNcluster->GetXaxis()->SetTitle("No. of clusters");
262  fHistRhoScaledVsNcluster->GetYaxis()->SetTitle("#rho_{scaled} (GeV/#it{c} #times rad^{-1})");
264  }
265  }
266 }
267 
269 {
270  Double_t rho = GetRhoFactor(fCent);
271  fOutRho->SetVal(rho);
272 }
273 
275 {
276  fOutRho->SetVal(0);
277  if (fOutRhoScaled) fOutRhoScaled->SetVal(0);
278 
280 
281  CalculateRho();
282 
283  if (fScaleFunction) {
284  Double_t rhoScaled = fOutRho->GetVal() * GetScaleFactor(fCent);
285  fOutRhoScaled->SetVal(rhoScaled);
286  }
287 
288  return kTRUE;
289 }
290 
292 {
293  fHistRhoVsCent->Fill(fCent, fOutRho->GetVal());
294 
295  if (fLeadingParticle) {
297  fHistRhoVsLeadTrackPt->Fill(fLeadingParticle->Pt(), fOutRho->GetVal());
298  }
299 
300  if (fLeadingCluster) {
302  fHistRhoVsLeadClusterE->Fill(fLeadingCluster->E(), fOutRho->GetVal());
303  }
304 
307 
308 
309  if (fHistRhoVsNtrack) fHistRhoVsNtrack->Fill(fNtracks, fOutRho->GetVal());
311 
312  if (fOutRhoScaled) {
313  fHistRhoScaledVsCent->Fill(fCent, fOutRhoScaled->GetVal());
316  }
317 
318  for (auto jetCont : fJetCollArray) {
319  fHistTotJetAreaVsCent[jetCont.first]->Fill(fCent, fTotJetArea[jetCont.first]);
320  if (fLeadingJet[jetCont.first]) {
321  fHistLeadJetPtVsCent[jetCont.first]->Fill(fCent, fLeadingJet[jetCont.first]->Pt());
322  fHistLeadJetPtDensityVsCent[jetCont.first]->Fill(fCent, fLeadingJet[jetCont.first]->Pt() / fLeadingJet[jetCont.first]->Area());
323  fHistLeadJetNconstVsCent[jetCont.first]->Fill(fCent, fLeadingJet[jetCont.first]->GetNumberOfConstituents());
324  fHistRhoVsLeadJetPt[jetCont.first]->Fill(fLeadingJet[jetCont.first]->Pt(), fOutRho->GetVal());
325  if (fCentBin >=0 && fCentBin < fCentBins.size()-1 && fHistLeadJetNconstVsPt[jetCont.first]) fHistLeadJetNconstVsPt[jetCont.first][fCentBin]->Fill(fLeadingJet[jetCont.first]->Pt(), fLeadingJet[jetCont.first]->GetNumberOfConstituents());
326  }
327  if (fHistNjetVsCent[jetCont.first]) fHistNjetVsCent[jetCont.first]->Fill(fCent, fNjets[jetCont.first]);
328 
329  if (fHistNjetVsNtrack[jetCont.first]) fHistNjetVsNtrack[jetCont.first]->Fill(fNtracks, fNjets[jetCont.first]);
330  }
331 
332  return kTRUE;
333 }
334 
336 {
337  if (!fOutRho) {
339 
340  if (fAttachToEvent) {
341  if (!(InputEvent()->FindListObject(fOutRhoName))) {
342  InputEvent()->AddObject(fOutRho);
343  } else {
344  AliFatal(Form("%s: Container with same name %s already present. Aborting", GetName(), fOutRhoName.Data()));
345  return;
346  }
347  }
348  }
349 
350  if (fScaleFunction && !fOutRhoScaled) {
352 
353  if (fAttachToEvent) {
354  if (!(InputEvent()->FindListObject(fOutRhoScaledName))) {
355  InputEvent()->AddObject(fOutRhoScaled);
356  } else {
357  AliFatal(Form("%s: Container with same name %s already present. Aborting", GetName(), fOutRhoScaledName.Data()));
358  return;
359  }
360  }
361  }
362 
364 
366 }
367 
369 {
370  Double_t rho = 0;
371  if (fRhoFunction) rho = fRhoFunction->Eval(cent);
372  return rho;
373 }
374 
376 {
377  Double_t scale = 1;
378  if (fScaleFunction) scale = fScaleFunction->Eval(cent);
379  return scale;
380 }
381 
382 TF1* AliAnalysisTaskRhoBaseDev::LoadRhoFunction(const char* path, const char* name)
383 {
384  TString fname(path);
385  if (fname.BeginsWith("alien://")) {
386  TGrid::Connect("alien://");
387  }
388 
389  TFile* file = TFile::Open(path);
390 
391  if (!file || file->IsZombie()) {
392  ::Error("AddTaskRho", "Could not open scale function file");
393  return 0;
394  }
395 
396  TF1* sfunc = dynamic_cast<TF1*>(file->Get(name));
397 
398  if (sfunc) {
399  ::Info("AliAnalysisTaskRhoBaseDev::LoadRhoFunction", "Scale function %s loaded from file %s.", name, path);
400  }
401  else {
402  ::Error("AliAnalysisTaskRhoBaseDev::LoadRhoFunction", "Scale function %s not found in file %s.", name, path);
403  return 0;
404  }
405 
406  fScaleFunction = static_cast<TF1*>(sfunc->Clone());
407 
408  file->Close();
409  delete file;
410 
411  return fScaleFunction;
412 }
413 
415 {
416  // Get the pointer to the existing analysis manager via the static access method.
417  AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
418  if (!mgr) {
419  ::Error("AliAnalysisTaskRhoBaseDev::AddTaskRhoBaseDev", "No analysis manager to connect to.");
420  return nullptr;
421  }
422 
423  // Check the analysis type using the event handlers connected to the analysis manager.
424  AliVEventHandler* handler = mgr->GetInputEventHandler();
425  if (!handler) {
426  ::Error("AliAnalysisTaskRhoBaseDev::AddTaskRhoBaseDev", "This task requires an input event handler");
427  return nullptr;
428  }
429 
431 
432  if (handler->InheritsFrom("AliESDInputHandler")) {
433  dataType = kESD;
434  }
435  else if (handler->InheritsFrom("AliAODInputHandler")) {
436  dataType = kAOD;
437  }
438 
439  // Init the task and do settings
440  if (trackName == "usedefault") {
441  if (dataType == kESD) {
442  trackName = "Tracks";
443  }
444  else if (dataType == kAOD) {
445  trackName = "tracks";
446  }
447  else {
448  trackName = "";
449  }
450  }
451 
452  if (clusName == "usedefault") {
453  if (dataType == kESD) {
454  clusName = "CaloClusters";
455  }
456  else if (dataType == kAOD) {
457  clusName = "caloClusters";
458  }
459  else {
460  clusName = "";
461  }
462  }
463 
464  TString name("AliAnalysisTaskRhoBaseDev");
465  if (!suffix.IsNull()) {
466  name += "_";
467  name += suffix;
468  }
469 
470  AliAnalysisTaskRhoBaseDev* mgrTask = dynamic_cast<AliAnalysisTaskRhoBaseDev*>(mgr->GetTask(name.Data()));
471  if (mgrTask) {
472  ::Warning("AliAnalysisTaskRhoBaseDev::AddTaskRhoBaseDev", "Not adding the task again, since a task with the same name '%s' already exists", name.Data());
473  return mgrTask;
474  }
475 
476  AliAnalysisTaskRhoBaseDev* rhotask = new AliAnalysisTaskRhoBaseDev(name, histo);
477  rhotask->SetOutRhoName(nRho);
478 
479  AliParticleContainer* partCont = rhotask->AddParticleContainer(trackName.Data());
480  AliClusterContainer *clusterCont = rhotask->AddClusterContainer(clusName.Data());
481  if (clusterCont) {
482  clusterCont->SetClusECut(0.);
483  clusterCont->SetClusPtCut(0.);
484  clusterCont->SetClusHadCorrEnergyCut(0.3);
485  clusterCont->SetDefaultClusterEnergy(AliVCluster::kHadCorr);
486  }
487 
488  AliJetContainer *jetCont = new AliJetContainer(jetType, AliJetContainer::kt_algorithm, rscheme, jetradius, partCont, clusterCont);
489  if (jetCont) {
490  jetCont->SetJetPtCut(0);
491  jetCont->SetJetAcceptanceType(acceptance);
492  jetCont->SetName("Background");
493  rhotask->AdoptJetContainer(jetCont);
494  }
495 
496  // Final settings, pass to manager and set the containers
497  mgr->AddTask(rhotask);
498 
499  // Create containers for input/output
500  mgr->ConnectInput(rhotask, 0, mgr->GetCommonInputContainer());
501  if (histo) {
502  TString contname(name);
503  contname += "_histos";
504  AliAnalysisDataContainer *coutput1 = mgr->CreateContainer(contname.Data(),
505  TList::Class(), AliAnalysisManager::kOutputContainer,
506  Form("%s", AliAnalysisManager::GetCommonFileName()));
507  mgr->ConnectOutput(rhotask, 1, coutput1);
508  }
509 
510  return rhotask;
511 }
TString fOutRhoName
name of output rho object
std::map< std::string, Int_t > fNjets
!number of jets
double Double_t
Definition: External.C:58
EDataType_t
Switch for the data type.
Definition: External.C:236
TF1 * fScaleFunction
pre-computed scale factor as a function of centrality
void SetName(const char *n)
Set the name of the class of the objets inside the underlying array.
AliVCluster * fLeadingCluster
!leading cluster
TString fOutRhoScaledName
name of output scaled rho object
virtual Double_t GetRhoFactor(Double_t cent)
Access rho per centrality.
TH2 * fHistRhoVsCent
!rho vs. centrality
Int_t fNtracks
!number of tracks
Base class for a task that calculates the UE.
Bool_t fAttachToEvent
whether or not attach rho to the event objects list
AliClusterContainer * AddClusterContainer(EMCAL_STRINGVIEW branchName, EMCAL_STRINGVIEW contName="")
std::map< std::string, TH2 * > fHistRhoVsLeadJetPt
!rho vs. leading jet pt
Bool_t fTaskConfigured
!kTRUE if the task is properly configured
Bool_t Run()
Run the analysis.
std::map< std::string, TH2 * > fHistLeadJetPtVsCent
!leading jet pt vs. centrality
std::vector< double > fCentBins
how many centrality bins
std::map< std::string, TH2 * > fHistTotJetAreaVsCent
!total area covered by jets vs. centrality
TH2 * fHistNtrackVsCent
!no. of tracks vs. centrality
TH2 * fHistRhoVsLeadClusterE
!rho vs. leading cluster energy
Container for particles within the EMCAL framework.
AliParticleContainer * AddParticleContainer(EMCAL_STRINGVIEW branchName, EMCAL_STRINGVIEW contName="")
EBeamType_t fForceBeamType
forced beam type
std::map< std::string, TH2 * > fHistLeadJetPtDensityVsCent
!leading jet area vs. centrality
std::map< std::string, TH2 ** > fHistLeadJetNconstVsPt
!leading jet constituents vs. pt
Bool_t fCreateHisto
whether or not create histograms
AliAnalysisTaskRhoBaseDev()
Default constructor. Needed by ROOT I/O.
int Int_t
Definition: External.C:63
void SetJetPtCut(Float_t cut)
void SetOutRhoName(const char *name)
unsigned int UInt_t
Definition: External.C:33
Float_t fMaxPt
Histogram pt limit.
TF1 * LoadRhoFunction(const char *path, const char *name)
Load the scale function from a file.
TH2 * fHistRhoScaledVsCent
!rhoscaled vs. centrality
TH2 * fHistLeadClusterEVsCent
!leading cluster energy vs. centrality
Int_t fCentBin
!event centrality bin
TH2 * fHistRhoScaledVsNcluster
!rhoscaled vs. no. of clusters
std::map< std::string, TH2 * > fHistLeadJetNconstVsCent
!leading jet constituents vs. cent
virtual Double_t GetScaleFactor(Double_t cent)
Get scale factor.
static AliAnalysisTaskRhoBaseDev * AddTaskRhoBaseDev(TString nTracks="usedefault", TString nClusters="usedefault", TString nRho="Rho", Double_t jetradius=0.2, UInt_t acceptance=AliEmcalJet::kTPCfid, AliJetContainer::EJetType_t jetType=AliJetContainer::kChargedJet, AliJetContainer::ERecoScheme_t rscheme=AliJetContainer::pt_scheme, Bool_t histo=kTRUE, TString suffix="")
Create an instance of this class and add it to the analysis manager.
TH2 * fHistRhoVsNcluster
!rho vs. no. of clusters
std::map< std::string, TH2 * > fHistNjetVsNtrack
!no. of jets vs. no. of tracks
TH2 * fHistRhoVsLeadTrackPt
!rho vs. leading track pt
Base class for a task that studies the UE.
Float_t fPtBinWidth
Histogram pt bin width.
Int_t fNclusters
!number of clusters
AliVParticle * fLeadingParticle
!leading particle
Definition: External.C:220
TH2 * fHistNclusterVsCent
!no. of cluster vs. centrality
TH2 * fHistRhoVsNtrack
!rho vs. no. of tracks
std::map< std::string, AliParticleContainer * > fParticleCollArray
particle/track collection array
TFile * file
TList with histograms for a given trigger.
TF1 * fRhoFunction
pre-computed rho as a function of centrality
virtual void CalculateEventProperties()
TH2 * fHistRhoScaledVsNtrack
!rhoscaled vs. no. of tracks
void SetClusPtCut(Double_t cut)
std::map< std::string, AliEmcalJet * > fLeadingJet
!leading jet
void AdoptJetContainer(AliJetContainer *cont)
std::map< std::string, AliJetContainer * > fJetCollArray
jet collection array
virtual void CalculateRho()
Calculates the average background using a given parametrization as a function of centrality.
void SetJetAcceptanceType(UInt_t type)
bool Bool_t
Definition: External.C:53
AliRhoParameter * fOutRhoScaled
!output scaled rho object
void SetClusECut(Double_t cut)
void SetDefaultClusterEnergy(Int_t d)
std::map< std::string, AliClusterContainer * > fClusterCollArray
cluster collection array
Container structure for EMCAL clusters.
Bool_t FillHistograms()
Fill histograms.
std::map< std::string, Double_t > fTotJetArea
!total area covered by jets
AliRhoParameter * fOutRho
!output rho object
Double_t fCent
!event centrality
Container for jet within the EMCAL jet framework.
std::map< std::string, TH2 * > fHistNjetVsCent
!no. of jets vs. centrality
TH2 * fHistLeadTrackPtVsCent
!leading track pt vs. centrality
void UserCreateOutputObjects()
Creating user output.
void SetClusHadCorrEnergyCut(Double_t cut)