AliPhysics  8d00e07 (8d00e07)
AliAnalysisTaskTrackDensityData.cxx
Go to the documentation of this file.
1 /*
2  * AliAnalysisTaskTrackDensityData.cxx
3  *
4  * Created on: Mar 11, 2016
5  * Author: markus
6  */
7 #include <THistManager.h>
8 #include <TMath.h>
9 
10 #include "AliEmcalJet.h"
11 #include "AliJetContainer.h"
12 #include "AliTrackContainer.h"
15 
17 
18 namespace EMCalTriggerPtAnalysis {
19 
22  fHistos(NULL),
23  fTrackSelection(NULL),
24  fBinHandler(NULL),
25  fNameJetContainer(""),
26  fNameTrackContainer("")
27 {
28 }
29 
31  AliAnalysisTaskEmcalJet(name, true),
32  fHistos(NULL),
33  fTrackSelection(NULL),
37 {
39  double defaultJetRadii[] = {0, 0.01, 0.02, 0.05, 0.1, 0.2, 0.4},
40  defaultPtMinSteps[] = {0.5, 1, 2, 5, 10, 20},
41  defaultJetPtBins[] = {20, 40, 60, 80, 100, 150, 200, 1000};
42 
43  fBinHandler->SetLinearBinning("jetpt", 100, 0, 100.);
44  fBinHandler->SetLinearBinning("contributors", 101, -0.5, 100.5);
45  fBinHandler->SetBinning("jetradii", sizeof(defaultJetRadii)/sizeof(double) - 1, defaultJetRadii);
46  fBinHandler->SetBinning("ptmin", sizeof(defaultPtMinSteps)/sizeof(double) - 1, defaultPtMinSteps);
47  fBinHandler->SetBinning("jetlarge", sizeof(defaultJetPtBins)/sizeof(double) - 1, defaultJetPtBins);
48 }
49 
50 
52  if(fHistos) delete fHistos;
53  if(fBinHandler) delete fBinHandler;
54 }
55 
58  fHistos = new THistManager("histos");
60 
62  bininit.Create(fBinHandler);
63 
64  const TBinning *trackptbinning = fBinHandler->GetBinning("pt"),
65  *jetptbinning = fBinHandler->GetBinning("jetpt"),
66  *contributorbinning = fBinHandler->GetBinning("contributors"),
67  *jetradii = fBinHandler->GetBinning("jetradii"),
68  *ptminsteps = fBinHandler->GetBinning("ptmin"),
69  *jetptlarge = fBinHandler->GetBinning("jetlarge");
70 
71  fHistos->CreateTH1("hTrackPtSel", "Pt spectrum of selected tracks", *trackptbinning);
72  fHistos->CreateTH1("hTrackPtSelEvent", "Pt spectrum of selected tracks (directly from the input event)", *trackptbinning);
73  fHistos->CreateTH2("hJetMultiplicity", "Multiplicity of particles in jets", *jetptbinning, *contributorbinning);
74  fHistos->CreateTH2("hParticlePtJet", "Correlation between track pt and jet pt", *jetptbinning, *trackptbinning);
75 
76  TArrayD jrad, ptms, jptl;
77  jetradii->CreateBinEdges(jrad);
78  ptminsteps->CreateBinEdges(ptms);
79  jetptlarge->CreateBinEdges(jptl);
80  for(int irad = 0; irad < jrad.GetSize()-1; irad++){
81  for(int ptstep = 0; ptstep <ptms.GetSize(); ptstep++){
82  fHistos->CreateTH2(Form("trackDensityJet_r%d_%d_minpt%d",
83  static_cast<int>(jrad[irad] * 100.),
84  static_cast<int>(jrad[irad+1] * 100.),
85  static_cast<int>(ptms[ptstep] * 10.)),
86  Form("Density of tracks with p_{t} > %f GeV/c in r [%.2f, %.2f]; p_{t, jet} (GeV/c); Number of tracks",
87  ptms[ptstep],
88  jrad[irad],
89  jrad[irad+1]),
90  200, 0., 200.,
91  102, -0.5, 100.5);
92  }
93  for(int jetptbin = 0 ; jetptbin < jptl.GetSize()-1; jetptbin++){
94  fHistos->CreateTH2(Form("trackDensityParticle_r%d_%d_jetpt%d_%d",
95  static_cast<int>(jrad[irad] * 100.),
96  static_cast<int>(jrad[irad+1] * 100.),
97  static_cast<int>(jptl[jetptbin]),
98  static_cast<int>(jptl[jetptbin+1])),
99  Form("Density of tracks in jet with p_{t} [%.1f, %.1f] in r[%.2f,%2f]",
100  jptl[jetptbin],
101  jptl[jetptbin+1],
102  jrad[irad],
103  jrad[irad+1]),
104  *trackptbinning, *contributorbinning);
105  }
106  }
107 
108  for(THistManager::iterator it = fHistos->begin(); it != fHistos->end(); it++){
109  fOutput->Add(*it);
110  }
111 
112  /*
113  for(auto it : *fHistos){
114  fOutput->Add(it);
115  }
116  */
117 }
118 
120  // Loop over jets
123  AliEmcalJet *myjet = NULL;
124  AliVParticle *jetparticle = NULL;
125  TArrayD particlePtBinning, jetradii, ptminsteps;
126  fBinHandler->GetBinning("pt")->CreateBinEdges(particlePtBinning);
127  fBinHandler->GetBinning("jetradii")->CreateBinEdges(jetradii);
128  fBinHandler->GetBinning("ptmin")->CreateBinEdges(ptminsteps);
129  for(int ipart = 0; ipart < InputEvent()->GetNumberOfTracks(); ipart++){
130  jetparticle = InputEvent()->GetTrack(ipart);
131  if(TMath::Abs(jetparticle->Eta()) > 0.8) continue;
132  if(!fTrackSelection->IsTrackAccepted(static_cast<AliVTrack *>(jetparticle))) continue;
133  fHistos->FillTH1("hTrackPtSelEvent", TMath::Abs(jetparticle->Pt()));
134  }
135  const AliTrackIterableContainer accepted_tracks = tcont->accepted();
136  for(AliTrackIterableContainer::iterator trackiter = accepted_tracks.begin(); trackiter != accepted_tracks.end(); ++trackiter){
137  jetparticle = *trackiter;
138  if(TMath::Abs(jetparticle->Eta()) > 0.8) continue;
139  if(!fTrackSelection->IsTrackAccepted(static_cast<AliVTrack *>(jetparticle))) continue;
140  fHistos->FillTH1("hTrackPtSel", TMath::Abs(jetparticle->Pt()));
141  }
142  const AliJetIterableContainer accepted_jets = jcont->accepted();
143  for(AliJetIterableContainer::iterator jetiter = accepted_jets.begin(); jetiter != accepted_jets.end(); ++jetiter){
144  myjet = *jetiter;
145  fHistos->FillTH2("hJetMultiplicity", myjet->Pt(), myjet->GetNumberOfConstituents());
146 
147  for(int iconst = 0; iconst < myjet->GetNumberOfTracks(); iconst++){
148  jetparticle = myjet->TrackAt(iconst, tcont->GetArray());
149  if(TMath::Abs(jetparticle->Eta()) > 0.8) continue;
150  if(!fTrackSelection->IsTrackAccepted(static_cast<AliVTrack *>(jetparticle))) continue;
151  fHistos->FillTH2("hParticlePtJet", myjet->Pt(), jetparticle->Pt());
152  }
153 
154  for(int irad = 0 ; irad < jetradii.GetSize()-1; irad++){
155  for(int ptstep = 0; ptstep < ptminsteps.GetSize(); ptstep++){
156  fHistos->FillTH2(Form("trackDensityJet_r%d_%d_minpt%d", static_cast<int>(jetradii[irad] * 100.), static_cast<int>(jetradii[irad+1] * 100.), static_cast<int>(ptminsteps[ptstep] * 10.)),
157  TMath::Abs(myjet->Pt()), GetParticleMultiplicity(*myjet, *tcont, ptminsteps[ptstep], 10000., jetradii[irad], jetradii[irad+1]));
158  }
159  double jetptmin, jetptmax;
160  FindJetPtBin(myjet, jetptmin, jetptmax);
161  if(jetptmin > 0 && jetptmax > 0){
162  for(int ptstep = 0; ptstep < particlePtBinning.GetSize()-1; ptstep++){
163  double mean = (particlePtBinning[ptstep] + particlePtBinning[ptstep+1])/2.;
164  fHistos->FillTH2(Form("trackDensityParticle_r%d_%d_jetpt%d_%d", static_cast<int>(jetradii[irad] * 100.),static_cast<int>(jetradii[irad+1] * 100.),
165  static_cast<int>(jetptmin), static_cast<int>(jetptmax)), mean, GetParticleMultiplicity(*myjet, *tcont, particlePtBinning[ptstep], particlePtBinning[+1], jetradii[irad], jetradii[irad+1]));
166  }
167  }
168  }
169 
170  }
171  return true;
172 }
173 
174 int AliAnalysisTaskTrackDensityData::GetParticleMultiplicity(const AliEmcalJet &jet, const AliParticleContainer &partcont, double ptmin, double ptmax, double rmin, double rmax) const {
175  AliDebug(1, Form("Next jet: %s\n", jet.toString().Data()));
176  TLorentzVector jetaxis(jet.Px(), jet.Py(), jet.Pz(), jet.E());
177  int nselected = 0;
178  AliVParticle *jetparticle(NULL);
179  for(int ipart = 0; ipart < jet.GetNumberOfTracks(); ipart++){
180  jetparticle = static_cast<AliVParticle *>(jet.TrackAt(ipart, partcont.GetArray()));
181  if(TMath::Abs(jetparticle->Eta()) > 0.8) continue;
182  if(!fTrackSelection->IsTrackAccepted(static_cast<AliVTrack *>(jetparticle))) continue;
183  double partpt = TMath::Abs(jetparticle->Pt());
184  if(partpt >= ptmin && partpt < ptmax){
185  TLorentzVector partvector(jetparticle->Px(), jetparticle->Py(), jetparticle->Pz(), jetparticle->E());
186  double r = TMath::Abs(jetaxis.DeltaR(partvector));
187  if(r >= rmin && r < rmax) nselected++;
188  }
189  }
190  return nselected;
191 }
192 
193 void AliAnalysisTaskTrackDensityData::FindJetPtBin(const AliEmcalJet *const jet, double &ptmin, double &ptmax) const {
194  TArrayD jetptlarge;
195  const TBinning *jptlbin = fBinHandler->GetBinning("jetlarge");
196  jptlbin->CreateBinEdges(jetptlarge);
197  ptmin = ptmax = -1;
198  double jetpt = TMath::Abs(jet->Pt());
199  for(int ptstep = 0; ptstep < jetptlarge.GetSize() - 1; ptstep++){
200  if(jetpt >= jetptlarge[ptstep] && jetpt < jetptlarge[ptstep+1]){
201  ptmin = jetptlarge[ptstep];
202  ptmax = jetptlarge[ptstep+1];
203  break;
204  }
205  }
206 }
207 
208 
209 } /* namespace EMCalTriggerPtAnalysis */
AliJetContainer * GetJetContainer(Int_t i=0) const
bidirectional stl iterator over the EMCAL iterable container
Double_t Py() const
Definition: AliEmcalJet.h:107
Container with name, TClonesArray and cuts for particles.
int GetParticleMultiplicity(const AliEmcalJet &jet, const AliParticleContainer &partcont, double ptmin, double ptmax, double rmin, double rmax) const
void FillTH2(const char *hname, double x, double y, double weight=1., Option_t *opt="")
Fill a 2D histogram within the container.
virtual void CreateBinEdges(TArrayD &binedges) const =0
Double_t E() const
Definition: AliEmcalJet.h:119
Interface for binnings used by the histogram handler.
Definition: TBinning.h:21
Container for particles within the EMCAL framework.
virtual bool Run()
Run function. This is the core function of the analysis and contains the user code. Therefore users have to implement this function.
Int_t TrackAt(Int_t idx) const
Definition: AliEmcalJet.h:160
UShort_t GetNumberOfTracks() const
Definition: AliEmcalJet.h:139
Double_t Px() const
Definition: AliEmcalJet.h:106
Global binning handler used by several analysis components.
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.
const Double_t ptmax
TString toString() const
TH1 * CreateTH1(const char *name, const char *title, int nbins, double xmin, double xmax, Option_t *opt="")
Create a new TH1 within the container.
void ReleaseOwner()
Definition: THistManager.h:254
const Double_t ptmin
void SetBinning(const char *dimname, int nbins, const double *binning)
stl-iterator for the histogram manager
Definition: THistManager.h:115
void FillTH1(const char *hname, double x, double weight=1., Option_t *opt="")
Fill a 1D histogram within the container.
Global binning definition for the high- charged particle analysis.
Double_t Pt() const
Definition: AliEmcalJet.h:109
iterator begin() const
Create forward iterator starting at the beginning of the container.
Definition: THistManager.h:730
virtual PWG::EMCAL::AliEmcalTrackSelResultPtr IsTrackAccepted(AliVTrack *const trk)=0
Interface for track selection code.
AliEmcalList * fOutput
!output list
AliEMCalTriggerBinningComponent * fBinHandler
EMCAL track selection.
Analysis of high- tracks in triggered events.
AliTrackContainer * GetTrackContainer(Int_t i=0) const
iterator end() const
Create forward iterator starting at the end of the container.
Definition: THistManager.h:738
void SetMakeGeneralHistograms(Bool_t g)
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
Double_t Pz() const
Definition: AliEmcalJet.h:108
void UserCreateOutputObjects()
Main initialization function on the worker.
const AliJetIterableContainer accepted() const
Container for jet within the EMCAL jet framework.
void SetLinearBinning(const char *dirname, int nbins, double min, double max)
void FindJetPtBin(const AliEmcalJet *const jet, double &ptmin, double &ptmax) const