AliPhysics  8bb951a (8bb951a)
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
AliAnalysisTaskEmcalSample.cxx
Go to the documentation of this file.
1 // $Id$
2 //
3 // Emcal sample analysis task.
4 //
5 // Author: S.Aiola, M. Verweij
6 
7 #include <TClonesArray.h>
8 #include <TH1F.h>
9 #include <TH2F.h>
10 #include <TH3F.h>
11 #include <TList.h>
12 #include <TLorentzVector.h>
13 
14 #include "AliVCluster.h"
15 #include "AliAODCaloCluster.h"
16 #include "AliESDCaloCluster.h"
17 #include "AliVTrack.h"
18 #include "AliLog.h"
19 #include "AliParticleContainer.h"
20 #include "AliClusterContainer.h"
21 #include "AliPicoTrack.h"
22 
24 
26 
27 //________________________________________________________________________
30  fHistTracksPt(0),
31  fHistClustersPt(0),
32  fHistPtDEtaDPhiTrackClus(0),
33  fHistPtDEtaDPhiClusTrack(0),
34  fTracksCont(0),
35  fCaloClustersCont(0)
36 {
37  // Default constructor.
38 
39  fHistTracksPt = new TH1*[fNcentBins];
40  fHistClustersPt = new TH1*[fNcentBins];
41 
42  for (Int_t i = 0; i < fNcentBins; i++) {
43  fHistTracksPt[i] = 0;
44  fHistClustersPt[i] = 0;
45  }
46 
47  SetMakeGeneralHistograms(kTRUE);
48 }
49 
50 //________________________________________________________________________
52  AliAnalysisTaskEmcal(name, kTRUE),
53  fHistTracksPt(0),
54  fHistClustersPt(0),
55  fHistPtDEtaDPhiTrackClus(0),
56  fHistPtDEtaDPhiClusTrack(0),
57  fTracksCont(0),
58  fCaloClustersCont(0)
59 {
60  // Standard constructor.
61 
62  fHistTracksPt = new TH1*[fNcentBins];
63  fHistClustersPt = new TH1*[fNcentBins];
64 
65  for (Int_t i = 0; i < fNcentBins; i++) {
66  fHistTracksPt[i] = 0;
67  fHistClustersPt[i] = 0;
68  }
69 
71 }
72 
73 //________________________________________________________________________
75 {
76  // Destructor.
77 }
78 
79 //________________________________________________________________________
81 {
82  // Create user output.
83 
85 
88  fTracksCont->SetClassName("AliVTrack");
89  fCaloClustersCont->SetClassName("AliVCluster");
90 
91  TString histname;
92 
93  for (Int_t i = 0; i < fNcentBins; i++) {
94  if (fParticleCollArray.GetEntriesFast()>0) {
95  histname = "fHistTracksPt_";
96  histname += i;
97  fHistTracksPt[i] = new TH1F(histname.Data(), histname.Data(), fNbins / 2, fMinBinPt, fMaxBinPt / 2);
98  fHistTracksPt[i]->GetXaxis()->SetTitle("p_{T,track} (GeV/c)");
99  fHistTracksPt[i]->GetYaxis()->SetTitle("counts");
100  fOutput->Add(fHistTracksPt[i]);
101  }
102 
103  if (fClusterCollArray.GetEntriesFast()>0) {
104  histname = "fHistClustersPt_";
105  histname += i;
106  fHistClustersPt[i] = new TH1F(histname.Data(), histname.Data(), fNbins / 2, fMinBinPt, fMaxBinPt / 2);
107  fHistClustersPt[i]->GetXaxis()->SetTitle("p_{T,clus} (GeV/c)");
108  fHistClustersPt[i]->GetYaxis()->SetTitle("counts");
109  fOutput->Add(fHistClustersPt[i]);
110  }
111  }
112 
113  histname = "fHistPtDEtaDPhiTrackClus";
114  fHistPtDEtaDPhiTrackClus = new TH3F(histname.Data(),Form("%s;#it{p}_{T}^{track};#Delta#eta;#Delta#varphi",histname.Data()),100,0.,100.,100,-0.1,0.1,100,-0.1,0.1);
116 
117  histname = "fHistPtDEtaDPhiClusTrack";
118  fHistPtDEtaDPhiClusTrack = new TH3F(histname.Data(),Form("%s;#it{p}_{T}^{clus};#Delta#eta;#Delta#varphi",histname.Data()),100,0.,100.,100,-0.1,0.1,100,-0.1,0.1);
120 
121  PostData(1, fOutput); // Post data for ALL output slots > 0 here.
122 }
123 
124 //________________________________________________________________________
126 {
127  // Fill histograms.
128 
129  if (fTracksCont) {
131  AliVTrack *track = static_cast<AliVTrack*>(fTracksCont->GetNextAcceptParticle());
132  while(track) {
133  fHistTracksPt[fCentBin]->Fill(track->Pt());
134  track = static_cast<AliVTrack*>(fTracksCont->GetNextAcceptParticle());
135  }
136  }
137 
138  if (fCaloClustersCont) {
140  AliVCluster *cluster = fCaloClustersCont->GetNextAcceptCluster();
141  while(cluster) {
142  TLorentzVector nPart;
143  cluster->GetMomentum(nPart, fVertex);
144  fHistClustersPt[fCentBin]->Fill(nPart.Pt());
146  }
147  }
148 
150 
151  return kTRUE;
152 }
153 
154 //________________________________________________________________________
156 {
157 
159  return;
160 
161  Double_t deta = 999;
162  Double_t dphi = 999;
163 
164  //Get closest cluster to track
166  AliVTrack *track = static_cast<AliVTrack*>(fTracksCont->GetNextAcceptParticle());
167  while(track) {
168  //Get matched cluster
169  Int_t emc1 = track->GetEMCALcluster();
170  if(fCaloClustersCont && emc1>=0) {
171  AliVCluster *clusMatch = fCaloClustersCont->GetCluster(emc1);
172  if(clusMatch) {
173  AliPicoTrack::GetEtaPhiDiff(track, clusMatch, dphi, deta);
174  fHistPtDEtaDPhiTrackClus->Fill(track->Pt(),deta,dphi);
175  }
176  }
177  track = static_cast<AliVTrack*>(fTracksCont->GetNextAcceptParticle());
178  }
179 
180  //Get closest track to cluster
182  AliVCluster *cluster = fCaloClustersCont->GetNextAcceptCluster();
183  while(cluster) {
184  TLorentzVector nPart;
185  cluster->GetMomentum(nPart, fVertex);
186  fHistClustersPt[fCentBin]->Fill(nPart.Pt());
187 
188  //Get matched track
189  AliVTrack *mt = NULL;
190  AliAODCaloCluster *acl = dynamic_cast<AliAODCaloCluster*>(cluster);
191  if(acl) {
192  if(acl->GetNTracksMatched()>1)
193  mt = static_cast<AliVTrack*>(acl->GetTrackMatched(0));
194  }
195  else {
196  AliESDCaloCluster *ecl = dynamic_cast<AliESDCaloCluster*>(cluster);
197  Int_t im = ecl->GetTrackMatchedIndex();
198  if(fTracksCont && im>=0) {
199  mt = static_cast<AliVTrack*>(fTracksCont->GetParticle(im));
200  }
201  }
202  if(mt) {
203  AliPicoTrack::GetEtaPhiDiff(mt, cluster, dphi, deta);
204  fHistPtDEtaDPhiClusTrack->Fill(nPart.Pt(),deta,dphi);
205 
206  /* //debugging
207  if(mt->IsEMCAL()) {
208  Int_t emc1 = mt->GetEMCALcluster();
209  Printf("current id: %d emc1: %d",fCaloClustersCont->GetCurrentID(),emc1);
210  AliVCluster *clm = fCaloClustersCont->GetCluster(emc1);
211  AliPicoTrack::GetEtaPhiDiff(mt, clm, dphi, deta);
212  Printf("deta: %f dphi: %f",deta,dphi);
213  }
214  */
215  }
217  }
218 }
219 
220 //________________________________________________________________________
222 
224 
225  if (fTracksCont && fTracksCont->GetArray() == 0) fTracksCont = 0;
227 
228 }
229 
230 //________________________________________________________________________
232 {
233  // Run analysis code here, if needed. It will be executed before FillHistograms().
234 
235  return kTRUE; // If return kFALSE FillHistogram() will NOT be executed.
236 }
237 
238 //________________________________________________________________________
240 {
241  // Called once at the end of the analysis.
242 }
TObjArray fClusterCollArray
cluster collection array
virtual AliVParticle * GetNextAcceptParticle()
Base task in the EMCAL framework.
Double_t fMinBinPt
min pt in histograms
Int_t fCentBin
!event centrality bin
AliClusterContainer * fCaloClustersCont
Tracks.
TList * fOutput
!output list
TH3 * fHistPtDEtaDPhiClusTrack
track pt, delta eta, delta phi to matched cluster
TObjArray fParticleCollArray
particle/track collection array
TH1 ** fHistClustersPt
Track pt spectrum.
AliParticleContainer * GetParticleContainer(Int_t i=0) const
TH3 * fHistPtDEtaDPhiTrackClus
Cluster pt spectrum.
virtual AliVParticle * GetParticle(Int_t i=-1) const
Int_t fNcentBins
how many centrality bins
AliClusterContainer * GetClusterContainer(Int_t i=0) const
TClonesArray * GetArray() const
AliVCluster * GetCluster(Int_t i) const
void SetClassName(const char *clname)
ClassImp(AliAnalysisTaskEmcalSample) AliAnalysisTaskEmcalSample
Double_t fMaxBinPt
max pt in histograms
Double_t fVertex[3]
!event vertex
void SetMakeGeneralHistograms(Bool_t g)
AliParticleContainer * fTracksCont
cluster pt, delta eta, delta phi to matched track
static void GetEtaPhiDiff(const AliVTrack *t, const AliVCluster *v, Double_t &phidiff, Double_t &etadiff)
AliVCluster * GetNextAcceptCluster()
void ResetCurrentID(Int_t i=-1)
Int_t fNbins
no. of pt bins