5 // Author: R.Reed, M.Connors
7 #include "AliAnalysisTaskScale.h"
9 #include <TClonesArray.h>
12 #include <TLorentzVector.h>
15 #include "AliEMCALGeometry.h"
17 #include "AliVCluster.h"
18 #include "AliVParticle.h"
19 #include "AliVTrack.h"
20 #include "AliParticleContainer.h"
21 #include "AliClusterContainer.h"
23 ClassImp(AliAnalysisTaskScale)
25 //________________________________________________________________________
26 AliAnalysisTaskScale::AliAnalysisTaskScale() :
27 AliAnalysisTaskEmcal("AliAnalysisTaskScale", kTRUE),
32 fHistPtEMCALvsCent(0),
35 fHistDeltaScalevsCent(0),
36 fHistScaleEmcalvsCent(0),
37 fHistScale2EmcalvsCent(0),
38 fHistDeltaScale2EmcalvsCent(0),
39 fHistChScalevsCent(0),
40 fHistChScale2EmcalvsCent(0),
41 fHistPtTPCvsNtrack(0),
42 fHistPtEMCALvsNtrack(0),
44 fHistScalevsNtrack(0),
45 fHistDeltaScalevsNtrack(0),
46 fHistScaleEmcalvsNtrack(0),
47 fHistScale2EmcalvsNtrack(0),
48 fHistChScalevsNtrack(0),
49 fHistChScale2EmcalvsNtrack(0),
50 fHistTrackPtvsCent(0),
51 fHistClusterPtvsCent(0),
53 fHistClusterEtaPhi(0),
54 fHistScalevsScale2Emcal(0),
55 fHistScalevsScaleEmcal(0),
56 fHistScaleEmcalvsScale2Emcal(0),
60 // Default constructor.
62 SetMakeGeneralHistograms(kTRUE);
65 //________________________________________________________________________
66 AliAnalysisTaskScale::AliAnalysisTaskScale(const char *name) :
67 AliAnalysisTaskEmcal(name, kTRUE),
72 fHistPtEMCALvsCent(0),
75 fHistDeltaScalevsCent(0),
76 fHistScaleEmcalvsCent(0),
77 fHistScale2EmcalvsCent(0),
78 fHistDeltaScale2EmcalvsCent(0),
79 fHistChScalevsCent(0),
80 fHistChScale2EmcalvsCent(0),
81 fHistPtTPCvsNtrack(0),
82 fHistPtEMCALvsNtrack(0),
84 fHistScalevsNtrack(0),
85 fHistDeltaScalevsNtrack(0),
86 fHistScaleEmcalvsNtrack(0),
87 fHistScale2EmcalvsNtrack(0),
88 fHistChScalevsNtrack(0),
89 fHistChScale2EmcalvsNtrack(0),
90 fHistTrackPtvsCent(0),
91 fHistClusterPtvsCent(0),
93 fHistClusterEtaPhi(0),
94 fHistScalevsScale2Emcal(0),
95 fHistScalevsScaleEmcal(0),
96 fHistScaleEmcalvsScale2Emcal(0),
102 SetMakeGeneralHistograms(kTRUE);
105 //________________________________________________________________________
106 void AliAnalysisTaskScale::UserCreateOutputObjects()
108 // Create my user objects.
110 AliAnalysisTaskEmcal::UserCreateOutputObjects();
112 //Get track and particle container
113 fTracksCont = GetParticleContainer(0);
114 fCaloClustersCont = GetClusterContainer(0);
115 if(fTracksCont) fTracksCont->SetClassName("AliVTrack");
116 if(fCaloClustersCont) fCaloClustersCont->SetClassName("AliVCluster");
119 fHistPtTPCvsCent = new TH2F("fHistPtTPCvsCent", "fHistPtTPCvsCent", 101, -1, 100, 750, 0, 1500);
120 fHistPtTPCvsCent->GetXaxis()->SetTitle("Centrality (%)");
121 fHistPtTPCvsCent->GetYaxis()->SetTitle("#sum p_{T,track}^{TPC} GeV/c");
122 fHistPtTPCvsCent->GetZaxis()->SetTitle("counts");
123 fOutput->Add(fHistPtTPCvsCent);
125 fHistPtEMCALvsCent = new TH2F("fHistPtEMCALvsCent", "fHistPtEMCALvsCent", 101, -1, 100, 250, 0, 500);
126 fHistPtEMCALvsCent->GetXaxis()->SetTitle("Centrality (%)");
127 fHistPtEMCALvsCent->GetYaxis()->SetTitle("#sum p_{T,track}^{EMCal} GeV/c");
128 fHistPtEMCALvsCent->GetZaxis()->SetTitle("counts");
129 fOutput->Add(fHistPtEMCALvsCent);
131 fHistEtvsCent = new TH2F("fHistEtvsCent", "fHistEtvsCent", 101, -1, 100, 250, 0, 500);
132 fHistEtvsCent->GetXaxis()->SetTitle("Centrality (%)");
133 fHistEtvsCent->GetYaxis()->SetTitle("#sum E_{T,cluster} GeV");
134 fHistEtvsCent->GetZaxis()->SetTitle("counts");
135 fOutput->Add(fHistEtvsCent);
137 fHistScalevsCent = new TH2F("fHistScalevsCent", "fHistScalevsCent", 101, -1, 100, 500, 0, 5);
138 fHistScalevsCent->GetXaxis()->SetTitle("Centrality (%)");
139 fHistScalevsCent->GetYaxis()->SetTitle("s_{TPC} = (#sum E_{T,cluster} + #sum p_{T,track}^{TPC}) / #sum p_{T,track}^{TPC}");
140 fHistScalevsCent->GetZaxis()->SetTitle("counts");
141 fOutput->Add(fHistScalevsCent);
143 fHistDeltaScalevsCent = new TH2F("fHistDeltaScalevsCent", "fHistDeltaScalevsCent", 101, -1, 100, 500, -2.5, 2.5);
144 fHistDeltaScalevsCent->GetXaxis()->SetTitle("Centrality (%)");
145 fHistDeltaScalevsCent->GetYaxis()->SetTitle("s_{TPC}-s^{old}");
146 fHistDeltaScalevsCent->GetZaxis()->SetTitle("counts");
147 fOutput->Add(fHistDeltaScalevsCent);
149 fHistScaleEmcalvsCent= new TH2F("fHistScaleEmcalvsCent", "fHistScaleEmcalvsCent", 101, -1, 100, 500, 0, 5);
150 fHistScaleEmcalvsCent->GetXaxis()->SetTitle("Centrality (%)");
151 fHistScaleEmcalvsCent->GetYaxis()->SetTitle("s_{EMC}");
152 fHistScaleEmcalvsCent->GetZaxis()->SetTitle("counts");
153 fOutput->Add(fHistScaleEmcalvsCent);
155 fHistScale2EmcalvsCent = new TH2F("fHistScale2EmcalvsCent", "fHistScale2EmcalvsCent", 101, -1, 100, 500, 0, 5);
156 fHistScale2EmcalvsCent->GetXaxis()->SetTitle("Centrality (%)");
157 fHistScale2EmcalvsCent->GetYaxis()->SetTitle("s_{2 #times EMC}");
158 fHistScale2EmcalvsCent->GetZaxis()->SetTitle("counts");
159 fOutput->Add(fHistScale2EmcalvsCent);
161 fHistDeltaScale2EmcalvsCent = new TH2F("fHistDeltaScale2EmcalvsCent", "fHistDeltaScale2EmcalvsCent", 101, -1, 100, 500, -2.5, 2.5);
162 fHistDeltaScale2EmcalvsCent->GetXaxis()->SetTitle("Centrality (%)");
163 fHistDeltaScale2EmcalvsCent->GetYaxis()->SetTitle("s_{2 #times EMC}-s^{old}");
164 fHistDeltaScale2EmcalvsCent->GetZaxis()->SetTitle("counts");
165 fOutput->Add(fHistDeltaScale2EmcalvsCent);
167 fHistChScalevsCent = new TH2F("fHistChScalevsCent", "fHistChScalevsCent", 101, -1, 100, 500, 0, 5);
168 fHistChScalevsCent->GetXaxis()->SetTitle("Centrality (%)");
169 fHistChScalevsCent->GetYaxis()->SetTitle("s_{TPC}^{ch}");
170 fHistChScalevsCent->GetZaxis()->SetTitle("counts");
171 fOutput->Add(fHistChScalevsCent);
173 fHistChScale2EmcalvsCent = new TH2F("fHistChScale2EmcalvsCent", "fHistChScale2EmcalvsCent", 101, -1, 100, 500, 0, 5);
174 fHistChScale2EmcalvsCent->GetXaxis()->SetTitle("Centrality (%)");
175 fHistChScale2EmcalvsCent->GetYaxis()->SetTitle("s_{2 #times EMC}^{ch}");
176 fHistChScale2EmcalvsCent->GetZaxis()->SetTitle("counts");
177 fOutput->Add(fHistChScale2EmcalvsCent);
179 fHistPtTPCvsNtrack = new TH2F("fHistPtTPCvsNtrack", "fHistPtTPCvsNtrack", 800, 0, 4000, 750, 0, 1500);
180 fHistPtTPCvsNtrack->GetXaxis()->SetTitle("No. of tracks");
181 fHistPtTPCvsNtrack->GetYaxis()->SetTitle("#sum p_{T,track}^{TPC}");
182 fHistPtTPCvsNtrack->GetZaxis()->SetTitle("counts");
183 fOutput->Add(fHistPtTPCvsNtrack);
185 fHistPtEMCALvsNtrack = new TH2F("fHistPtEMCALvsNtrack", "fHistPtEMCALvsNtrack", 800, 0, 4000, 500, 0, 1000);
186 fHistPtEMCALvsNtrack->GetXaxis()->SetTitle("No. of tracks");
187 fHistPtEMCALvsNtrack->GetYaxis()->SetTitle("#sum p_{T,track}^{EMCal}");
188 fHistPtEMCALvsNtrack->GetZaxis()->SetTitle("counts");
189 fOutput->Add(fHistPtEMCALvsNtrack);
191 fHistEtvsNtrack = new TH2F("fHistEtvsNtrack", "fHistEtvsNtrack", 800, 0, 4000, 500, 0, 1000);
192 fHistEtvsNtrack->GetXaxis()->SetTitle("No. of tracks");
193 fHistEtvsNtrack->GetYaxis()->SetTitle("#sum E_{T,cluster}");
194 fHistEtvsNtrack->GetZaxis()->SetTitle("counts");
195 fOutput->Add(fHistEtvsNtrack);
197 fHistScalevsNtrack = new TH2F("fHistScalevsNtrack", "fHistScalevsNtrack", 800, 0, 4000, 500, 0, 5);
198 fHistScalevsNtrack->GetXaxis()->SetTitle("No. of tracks");
199 fHistScalevsNtrack->GetYaxis()->SetTitle("s_{TPC}");
200 fHistScalevsNtrack->GetZaxis()->SetTitle("counts");
201 fOutput->Add(fHistScalevsNtrack);
203 fHistDeltaScalevsNtrack = new TH2F("fHistDeltaScalevsNtrack", "fHistDeltaScalevsNtrack", 800, 0, 4000, 500, -2.5, 2.5);
204 fHistDeltaScalevsNtrack->GetXaxis()->SetTitle("No. of tracks");
205 fHistDeltaScalevsNtrack->GetYaxis()->SetTitle("s_{TPC}-s^{old}");
206 fHistDeltaScalevsNtrack->GetZaxis()->SetTitle("counts");
207 fOutput->Add(fHistDeltaScalevsNtrack);
209 fHistScaleEmcalvsNtrack = new TH2F("fHistScaleEmcalvsNtrack", "fHistScaleEmcalvsNtrack", 800, 0, 4000, 500, 0, 5);
210 fHistScaleEmcalvsNtrack->GetXaxis()->SetTitle("No. of tracks");
211 fHistScaleEmcalvsNtrack->GetYaxis()->SetTitle("s_{EMC}");
212 fHistScaleEmcalvsNtrack->GetZaxis()->SetTitle("counts");
213 fOutput->Add(fHistScaleEmcalvsNtrack);
215 fHistScale2EmcalvsNtrack = new TH2F("fHistScale2EmcalvsNtrack","fHistScale2EmcalvsNtrack", 800, 0, 4000, 500, 0, 5);
216 fHistScale2EmcalvsNtrack->GetXaxis()->SetTitle("No. of tracks");
217 fHistScale2EmcalvsNtrack->GetYaxis()->SetTitle("s_{2 #times EMC}");
218 fHistScale2EmcalvsNtrack->GetZaxis()->SetTitle("counts");
219 fOutput->Add(fHistScale2EmcalvsNtrack);
221 fHistChScalevsNtrack = new TH2F("fHistChScalevsNtrack", "fHistChScalevsNtrack", 800, 0, 4000, 500, 0, 5);
222 fHistChScalevsNtrack->GetXaxis()->SetTitle("No. of tracks");
223 fHistChScalevsNtrack->GetYaxis()->SetTitle("s_{TPC}^{ch}");
224 fHistChScalevsNtrack->GetZaxis()->SetTitle("counts");
225 fOutput->Add(fHistChScalevsNtrack);
227 fHistChScale2EmcalvsNtrack = new TH2F("fHistChScale2EmcalvsNtrack", "fHistChScale2EmcalvsNtrack", 800, 0, 4000, 500, 0, 5);
228 fHistChScale2EmcalvsNtrack->GetXaxis()->SetTitle("No. of tracks");
229 fHistChScale2EmcalvsNtrack->GetYaxis()->SetTitle("s_{2 #times EMC}^{ch}");
230 fHistChScale2EmcalvsNtrack->GetZaxis()->SetTitle("counts");
231 fOutput->Add(fHistChScale2EmcalvsNtrack);
233 fHistTrackPtvsCent = new TH2F("fHistTrackPtvsCent", "fHistTrackPtvsCent", 101, -1, 100, 500, 0, 100);
234 fHistTrackPtvsCent->GetXaxis()->SetTitle("Centrality (%)");
235 fHistTrackPtvsCent->GetYaxis()->SetTitle("p_{T,track} GeV/c");
236 fHistTrackPtvsCent->GetZaxis()->SetTitle("counts");
237 fOutput->Add(fHistTrackPtvsCent);
239 fHistClusterPtvsCent = new TH2F("fHistClusterPtvsCent", "fHistClusterPtvsCent", 101, -1, 100, 500, 0, 100);
240 fHistClusterPtvsCent->GetXaxis()->SetTitle("Centrality (%)");
241 fHistClusterPtvsCent->GetYaxis()->SetTitle("E_{T,cluster} GeV");
242 fHistClusterPtvsCent->GetZaxis()->SetTitle("counts");
243 fOutput->Add(fHistClusterPtvsCent);
245 fHistTrackEtaPhi = new TH2F("fHistTrackEtaPhi", "fHistTrackEtaPhi", 100, -1.0, 1.0, 101, 0, 2.02*TMath::Pi());
246 fHistTrackEtaPhi->GetXaxis()->SetTitle("#eta");
247 fHistTrackEtaPhi->GetYaxis()->SetTitle("#phi");
248 fHistTrackEtaPhi->GetZaxis()->SetTitle("counts");
249 fOutput->Add(fHistTrackEtaPhi);
251 fHistClusterEtaPhi = new TH2F("fHistClusterEtaPhi", "fHistClusterEtaPhi", 100, -1.0, 1.0, 101, 0, 2.02*TMath::Pi());
252 fHistClusterEtaPhi->GetXaxis()->SetTitle("#eta");
253 fHistClusterEtaPhi->GetYaxis()->SetTitle("#phi");
254 fHistClusterEtaPhi->GetZaxis()->SetTitle("counts");
255 fOutput->Add(fHistClusterEtaPhi);
257 fHistScalevsScale2Emcal = new TH2F("fHistScalevsScale2Emcal", "fHistScalevsScale2Emcal",500, 0, 5, 500, 0, 5);
258 fHistScalevsScale2Emcal->GetXaxis()->SetTitle("s_{TPC}");
259 fHistScalevsScale2Emcal->GetYaxis()->SetTitle("s_{2 #times EMC}");
260 fHistScalevsScale2Emcal->GetZaxis()->SetTitle("counts");
261 fOutput->Add(fHistScalevsScale2Emcal);
263 fHistScalevsScaleEmcal = new TH2F("fHistScalevsScaleEmcal", "fHistScalevsScaleEmcal", 500, 0, 5, 500, 0, 5);
264 fHistScalevsScaleEmcal->GetXaxis()->SetTitle("s_{TPC}");
265 fHistScalevsScaleEmcal->GetYaxis()->SetTitle("s_{EMC}");
266 fHistScalevsScaleEmcal->GetZaxis()->SetTitle("counts");
267 fOutput->Add(fHistScalevsScaleEmcal);
269 fHistScaleEmcalvsScale2Emcal = new TH2F("fHistScaleEmcalvsScale2Emcal", "fHistScaleEmcalvsScale2Emcal", 500, 0, 5, 500, 0, 5);
270 fHistScaleEmcalvsScale2Emcal->GetXaxis()->SetTitle("s_{EMC}");
271 fHistScaleEmcalvsScale2Emcal->GetYaxis()->SetTitle("s_{2 #times EMC}");
272 fHistScaleEmcalvsScale2Emcal->GetZaxis()->SetTitle("counts");
273 fOutput->Add(fHistScaleEmcalvsScale2Emcal);
275 PostData(1, fOutput);
278 //________________________________________________________________________
279 Double_t AliAnalysisTaskScale::GetScaleFactor(Double_t cent)
281 // Get scale function.
285 scale = fScaleFunction->Eval(cent);
289 //________________________________________________________________________
290 Bool_t AliAnalysisTaskScale::FillHistograms()
292 // Execute on each event.
294 Double_t EmcalMinPhi = fGeom->GetArm1PhiMin() * TMath::DegToRad();
295 Double_t EmcalMaxPhi = fGeom->GetArm1PhiMax() * TMath::DegToRad();
296 if(InputEvent()->GetRunNumber()>=177295 && InputEvent()->GetRunNumber()<=197470) { //small SM masked in 2012 and 2013
300 const Double_t EmcalWidth = (EmcalMaxPhi-EmcalMinPhi)/2.0;
303 Double_t ptEMCAL = 0;
304 Double_t ptEMCAL2 = 0;
306 const Int_t Ntracks = fTracksCont->GetNAcceptedParticles();
308 fTracksCont->ResetCurrentID();
309 AliVTrack *track = NULL;
310 while((track = static_cast<AliVTrack*>(fTracksCont->GetNextAcceptParticle()))) {
311 fHistTrackPtvsCent->Fill(fCent,track->Pt());
312 fHistTrackEtaPhi->Fill(track->Eta(),track->Phi());
313 ptTPC += track->Pt();
314 if ((track->Phi() > (EmcalMaxPhi+EmcalWidth)) || (track->Phi() < (EmcalMinPhi-EmcalWidth))) continue;
315 ptEMCAL2 += track->Pt();
316 if ((track->Phi() > EmcalMaxPhi) || (track->Phi() < EmcalMinPhi)) continue;
317 ptEMCAL += track->Pt();
325 if (fCaloClustersCont) {
326 fCaloClustersCont->ResetCurrentID();
327 AliVCluster *c = NULL;
328 while((c = fCaloClustersCont->GetNextAcceptCluster())) {
329 TLorentzVector nPart;
330 c->GetMomentum(nPart, fVertex);
332 fHistClusterPtvsCent->Fill(fCent, nPart.Pt());
333 fHistClusterEtaPhi->Fill(nPart.Eta(), nPart.Phi());
339 Double_t scalecalc = -1;
340 if (ptEMCAL > 0 && Et > 0 && ptTPC > 0)
341 scalecalc = ((Et + ptEMCAL) / fEmcalArea) * (fTpcArea / ptTPC);
342 const Double_t scale = GetScaleFactor(fCent);
343 Double_t scalecalcemcal = -1;
345 scalecalcemcal = (Et+ptEMCAL)/ptEMCAL;
346 Double_t scalecalcemcal2 = -1;
347 Double_t Chscalecalcemcal2 = -1;
349 scalecalcemcal2 = 2*(Et+ptEMCAL)/ptEMCAL2;
350 Chscalecalcemcal2 = 2*ptEMCAL/ptEMCAL2;}
351 const Double_t Chscalecalcemcal = ((ptEMCAL) / fEmcalArea) * (fTpcArea / ptTPC);
353 fHistScaleEmcalvsCent->Fill(fCent,scalecalcemcal);
354 fHistScale2EmcalvsCent->Fill(fCent,scalecalcemcal2);
355 fHistChScalevsCent->Fill(fCent,Chscalecalcemcal);
356 fHistChScale2EmcalvsCent->Fill(fCent,Chscalecalcemcal2);
357 fHistScaleEmcalvsNtrack->Fill(Ntracks,scalecalcemcal);
358 fHistScale2EmcalvsNtrack->Fill(Ntracks,scalecalcemcal2);
359 fHistChScalevsNtrack->Fill(Ntracks,Chscalecalcemcal);
360 fHistChScale2EmcalvsNtrack->Fill(Ntracks,Chscalecalcemcal2);
361 fHistPtTPCvsCent->Fill(fCent, ptTPC);
362 fHistPtEMCALvsCent->Fill(fCent, ptEMCAL);
363 fHistEtvsCent->Fill(fCent, Et);
364 fHistScalevsCent->Fill(fCent, scalecalc);
365 fHistDeltaScalevsCent->Fill(fCent, scalecalc - scale);
366 fHistDeltaScale2EmcalvsCent->Fill(fCent, scalecalcemcal2 - scale);
367 fHistPtTPCvsNtrack->Fill(Ntracks, ptTPC);
368 fHistPtEMCALvsNtrack->Fill(Ntracks, ptEMCAL);
369 fHistEtvsNtrack->Fill(Ntracks, Et);
370 fHistScalevsNtrack->Fill(Ntracks, scalecalc);
371 fHistDeltaScalevsNtrack->Fill(Ntracks, scalecalc - scale);
372 fHistScalevsScale2Emcal->Fill(scalecalc,scalecalcemcal2);
373 fHistScalevsScaleEmcal->Fill(scalecalc,scalecalcemcal);
374 fHistScaleEmcalvsScale2Emcal->Fill(scalecalcemcal,scalecalcemcal2);
380 //________________________________________________________________________
381 void AliAnalysisTaskScale::ExecOnce()
383 AliAnalysisTaskEmcal::ExecOnce();
385 const Double_t EmcalMinEta = fGeom->GetArm1EtaMin();
386 const Double_t EmcalMaxEta = fGeom->GetArm1EtaMax();
387 Double_t EmcalMinPhi = fGeom->GetArm1PhiMin() * TMath::DegToRad();
388 Double_t EmcalMaxPhi = fGeom->GetArm1PhiMax() * TMath::DegToRad();
389 if(InputEvent()->GetRunNumber()>=177295 && InputEvent()->GetRunNumber()<=197470) { //small SM masked in 2012 and 2013
394 fEmcalArea = (EmcalMaxPhi - EmcalMinPhi) * (EmcalMinEta - EmcalMaxEta);
396 AliParticleContainer *partCont = GetParticleContainer(0);
398 AliError(Form("%s: No particle container found! Assuming tpc area = 1...",GetName()));
403 Float_t TpcMaxPhi = partCont->GetParticlePhiMax();
404 Float_t TpcMinPhi = partCont->GetParticlePhiMin();
406 if (TpcMaxPhi > TMath::Pi()*2) TpcMaxPhi = TMath::Pi()*2;
407 if (TpcMinPhi < 0) TpcMinPhi = 0;
409 fTpcArea = (TpcMaxPhi - TpcMinPhi) * (EmcalMinEta - EmcalMaxEta);
411 if (fTracksCont && fTracksCont->GetArray() == 0) fTracksCont = 0;
412 if (fCaloClustersCont && fCaloClustersCont->GetArray() == 0) fCaloClustersCont = 0;