]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGJE/EMCALJetTasks/UserTasks/AliAnalysisTaskEmcalJetTriggerQA.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
timing calibration for 2013 period added
[u/mrichter/AliRoot.git] / PWGJE / EMCALJetTasks / UserTasks / AliAnalysisTaskEmcalJetTriggerQA.cxx
CommitLineData
0531c107 1//
2// Jet trigger QA analysis task.
3//
4// Author: M.Verweij
5
6#include <TClonesArray.h>
7#include <TH1F.h>
8#include <TH2F.h>
9#include <TH3F.h>
10#include <THnSparse.h>
11#include <TList.h>
12#include <TLorentzVector.h>
13
14#include "AliVCluster.h"
15#include "AliVTrack.h"
16#include "AliEmcalJet.h"
17#include "AliRhoParameter.h"
18#include "AliLog.h"
19#include "AliAnalysisUtils.h"
20#include "AliEmcalParticle.h"
21#include "AliAODCaloTrigger.h"
22#include "AliEMCALGeometry.h"
23#include "AliVCaloCells.h"
9e5eee5d 24#include "AliJetContainer.h"
25#include "AliClusterContainer.h"
0531c107 26
27#include "AliAnalysisTaskEmcalJetTriggerQA.h"
28
29ClassImp(AliAnalysisTaskEmcalJetTriggerQA)
30
31//________________________________________________________________________
32AliAnalysisTaskEmcalJetTriggerQA::AliAnalysisTaskEmcalJetTriggerQA() :
9e5eee5d 33 AliAnalysisTaskEmcalJetDev("AliAnalysisTaskEmcalJetTriggerQA", kTRUE),
0531c107 34 fDebug(kFALSE),
35 fUseAnaUtils(kTRUE),
36 fAnalysisUtils(0),
0531c107 37 fTriggerClass(""),
38 fBitJ1((1<<8)),
39 fBitJ2((1<<11)),
9e5eee5d 40 fContainerFull(0),
41 fContainerCharged(1),
0531c107 42 fMaxPatchEnergy(0),
43 fTriggerType(-1),
44 fNFastOR(16),
45 fhNEvents(0),
46 fh3PtEtaPhiJetFull(0),
47 fh3PtEtaPhiJetCharged(0),
48 fh2NJetsPtFull(0),
49 fh2NJetsPtCharged(0),
50 fh3PtEtaAreaJetFull(0),
51 fh3PtEtaAreaJetCharged(0),
52 fh2PtNConstituentsCharged(0),
53 fh2PtNConstituents(0),
54 fh2PtMeanPtConstituentsCharged(0),
55 fh2PtMeanPtConstituentsNeutral(0),
56 fh2PtNEF(0),
3c709670 57 fh3NEFEtaPhi(0),
58 fh2NEFNConstituentsCharged(0),
59 fh2NEFNConstituentsNeutral(0),
0531c107 60 fh2Ptz(0),
3c709670 61 fh2PtzCharged(0),
0531c107 62 fh2PtLeadJet1VsLeadJet2(0),
63 fh3EEtaPhiCluster(0),
64 fh3PtLeadJet1VsPatchEnergy(0),
65 fh3PtLeadJet2VsPatchEnergy(0),
992189ac 66 fh3PatchEnergyEtaPhiCenter(0),
67 fh2CellEnergyVsTime(0)
0531c107 68{
69 // Default constructor.
70
71
72 SetMakeGeneralHistograms(kTRUE);
73}
74
75//________________________________________________________________________
76AliAnalysisTaskEmcalJetTriggerQA::AliAnalysisTaskEmcalJetTriggerQA(const char *name) :
9e5eee5d 77 AliAnalysisTaskEmcalJetDev(name, kTRUE),
0531c107 78 fDebug(kFALSE),
79 fUseAnaUtils(kTRUE),
80 fAnalysisUtils(0),
0531c107 81 fTriggerClass(""),
82 fBitJ1((1<<8)),
83 fBitJ2((1<<11)),
9e5eee5d 84 fContainerFull(0),
85 fContainerCharged(1),
0531c107 86 fMaxPatchEnergy(0),
87 fTriggerType(-1),
88 fNFastOR(16),
89 fhNEvents(0),
90 fh3PtEtaPhiJetFull(0),
91 fh3PtEtaPhiJetCharged(0),
92 fh2NJetsPtFull(0),
93 fh2NJetsPtCharged(0),
94 fh3PtEtaAreaJetFull(0),
95 fh3PtEtaAreaJetCharged(0),
96 fh2PtNConstituentsCharged(0),
97 fh2PtNConstituents(0),
98 fh2PtMeanPtConstituentsCharged(0),
99 fh2PtMeanPtConstituentsNeutral(0),
100 fh2PtNEF(0),
3c709670 101 fh3NEFEtaPhi(0),
102 fh2NEFNConstituentsCharged(0),
103 fh2NEFNConstituentsNeutral(0),
0531c107 104 fh2Ptz(0),
3c709670 105 fh2PtzCharged(0),
0531c107 106 fh2PtLeadJet1VsLeadJet2(0),
107 fh3EEtaPhiCluster(0),
108 fh3PtLeadJet1VsPatchEnergy(0),
109 fh3PtLeadJet2VsPatchEnergy(0),
992189ac 110 fh3PatchEnergyEtaPhiCenter(0),
111 fh2CellEnergyVsTime(0)
0531c107 112{
113 // Standard constructor.
114
115 SetMakeGeneralHistograms(kTRUE);
116}
117
118//________________________________________________________________________
119AliAnalysisTaskEmcalJetTriggerQA::~AliAnalysisTaskEmcalJetTriggerQA()
120{
121 // Destructor.
122}
123
124//----------------------------------------------------------------------------------------------
125void AliAnalysisTaskEmcalJetTriggerQA::InitOnce() {
126 //
127 // only initialize once
128 //
129
130 // Initialize analysis util class for vertex selection
131 if(fUseAnaUtils) {
132 fAnalysisUtils = new AliAnalysisUtils();
133 fAnalysisUtils->SetMinVtxContr(2);
134 fAnalysisUtils->SetMaxVtxZ(10.);
135 }
136}
137
138//________________________________________________________________________
139Bool_t AliAnalysisTaskEmcalJetTriggerQA::SelectEvent() {
140 //
141 // Decide if event should be selected for analysis
142 //
143
144 fhNEvents->Fill(0.5);
145
146 if(fAnalysisUtils) {
147 if(!fAnalysisUtils->IsVertexSelected2013pA(InputEvent()))
148 return kFALSE;
149
150 fhNEvents->Fill(2.5);
151
152 if(fAnalysisUtils->IsPileUpEvent(InputEvent()))
153 return kFALSE;
154 }
155 else{
156 if(fUseAnaUtils)
157 AliError("fAnalysisUtils not initialized. Call AliAnalysisTaskEmcalJetTriggerQA::InitOnce()");
158 }
159
160 fhNEvents->Fill(3.5);
161
162 if(!fTriggerClass.IsNull()) {
163 //Check if requested trigger was fired
164 TString firedTrigClass = InputEvent()->GetFiredTriggerClasses();
165 if(!firedTrigClass.Contains(fTriggerClass))
166 return kFALSE;
167 }
168
169 fhNEvents->Fill(1.5);
170
171 return kTRUE;
172
173}
174
175//________________________________________________________________________
176void AliAnalysisTaskEmcalJetTriggerQA::FindTriggerPatch() {
177
178 //Code to get position of trigger
179
180 if(!fGeom)
181 fGeom = AliEMCALGeometry::GetInstance();
182
183
184 TString firedTrigClass = InputEvent()->GetFiredTriggerClasses();
185
186 AliAODCaloTrigger *trg = dynamic_cast<AliAODCaloTrigger*>(InputEvent()->GetCaloTrigger("EMCAL"));
187 trg->Reset();
188
189 int col, row; //FASTOR position
190
191 Int_t nPatchNotEmpty = 0; //counter number of patches which are not empty
192 fMaxPatchEnergy = 0.;
193 while (trg->Next())
194 {
195
196 trg->GetPosition(col, row); //col (0 to 63), row (0 to 47)
197
198 if (col > -1 && row > -1)
199 {
200 Int_t id = -1; //FASTOR index
201 Int_t cellIndex[4] = {-1};
202 fGeom->GetAbsFastORIndexFromPositionInEMCAL(col,row,id); //phi is column, eta is row
203 if(id<0) continue;
204
205 fGeom->GetCellIndexFromFastORIndex(id,cellIndex);
206
207 Int_t ts = 0;
208 trg->GetL1TimeSum(ts);
209
210 Float_t ampTrg = 0.;
211 trg->GetAmplitude(ampTrg);
212
213 //L1 analysis
214 Int_t bit = 0;
215 trg->GetTriggerBits(bit);
216
217 Bool_t bTrigJ = TestFilterBit(bit,fBitJ1|fBitJ2);
218 if(!bTrigJ)
219 continue;
220
221 Bool_t bTrigJ1 = TestFilterBit(bit,fBitJ1);
222 Bool_t bTrigJ2 = TestFilterBit(bit,fBitJ2);
223
224 if(bTrigJ1) fTriggerType = 0;
225 if(bTrigJ2) fTriggerType = 1;
226 if(!bTrigJ1 && !bTrigJ2) fTriggerType = -1;
227
228 Double_t minPhiPatch = 10.;
229 Double_t maxPhiPatch = -10.;
230 Double_t minEtaPatch = 10.;
231 Double_t maxEtaPatch = -10.;
232
233 //Get energy in trigger patch 8x8 FASTOR
234 Double_t patchEnergy = 0.;
235 Double_t sumAmp = 0.;
236 // const Int_t nFastOR = 8;//16;
237 for(Int_t fastrow = 0; fastrow<fNFastOR; fastrow++) {
238 for(Int_t fastcol = 0; fastcol<fNFastOR; fastcol++) {
239 Int_t nrow = row+fastrow;
240 Int_t ncol = col+fastcol;
241 fGeom->GetAbsFastORIndexFromPositionInEMCAL(ncol,nrow,id);
242
243 if(id<0) {
244 AliWarning(Form("%s: id smaller than 0 %d",GetName(),id));
245 continue;
246 }
247
248 fGeom->GetCellIndexFromFastORIndex(id,cellIndex);
249 for(Int_t icell=0; icell<4; icell++) {
250
251 if(!fGeom->CheckAbsCellId(cellIndex[icell])) continue;
252
253 Double_t amp =0., time = 0., efrac = 0;
254 Int_t mclabel = -1;
255 Short_t absId = -1;
256 Int_t nSupMod = -1, nModule = -1, nIphi = -1, nIeta = -1;
257
258 fGeom->GetCellIndex(cellIndex[icell], nSupMod, nModule, nIphi, nIeta );
259
260 fCaloCells->GetCell(cellIndex[icell], absId, amp, time,mclabel,efrac);
261
262 Double_t eta, phi;
263 fGeom->EtaPhiFromIndex(cellIndex[icell], eta, phi);
264
265 if(phi<minPhiPatch) minPhiPatch = phi;
266 if(phi>maxPhiPatch) maxPhiPatch = phi;
267 if(eta<minEtaPatch) minEtaPatch = eta;
268 if(eta>maxEtaPatch) maxEtaPatch = eta;
269
270 sumAmp+=fCaloCells->GetAmplitude(cellIndex[icell]);
271 patchEnergy+=amp;
272
273 }//cells in fastor loop
274 }//fastor col loop
275 }//fastor row loop
276
277 Double_t etaCent = minEtaPatch + 0.5*(maxEtaPatch-minEtaPatch);
278 Double_t phiCent = minPhiPatch + 0.5*(maxPhiPatch-minPhiPatch);
279 fh3PatchEnergyEtaPhiCenter->Fill(patchEnergy,etaCent,phiCent);
280
281 if(patchEnergy>0)
282 if(fDebug>2) AliInfo(Form("%s: patch edges eta: %f - %f phi: %f - %f ampTrg: %f patchEnergy: %f sumAmp: %f",GetName(),minEtaPatch,maxEtaPatch,minPhiPatch,maxPhiPatch,ampTrg,patchEnergy,sumAmp));
283
284 if(patchEnergy>0) nPatchNotEmpty++;
285
286 if(patchEnergy>fMaxPatchEnergy) fMaxPatchEnergy=patchEnergy;
287
288 }
289 }
290
291}
292
293//________________________________________________________________________
294void AliAnalysisTaskEmcalJetTriggerQA::LoadExtraBranches() {
295 //
9e5eee5d 296 // load extra brances
0531c107 297 //
298
0531c107 299}
300
301//________________________________________________________________________
302void AliAnalysisTaskEmcalJetTriggerQA::UserCreateOutputObjects()
303{
304 // Create user output.
305
306 InitOnce();
307
9e5eee5d 308 AliAnalysisTaskEmcalJetDev::UserCreateOutputObjects();
0531c107 309
310 Bool_t oldStatus = TH1::AddDirectoryStatus();
311 TH1::AddDirectory(kFALSE);
312
313 fhNEvents = new TH1F("fhNEvents","fhNEvents;selection;N_{evt}",5,0,5);
314 fOutput->Add(fhNEvents);
315
316 Int_t nBinsPt = 120;
317 Double_t minPt = -20.;
318 Double_t maxPt = 100.;
8daeee93 319 Int_t nBinsEta = 100;
0531c107 320 Double_t minEta = -1.;
321 Double_t maxEta = 1.;
8daeee93 322 Int_t nBinsPhi = 18*8;
0531c107 323 Double_t minPhi = 0.;
324 Double_t maxPhi = TMath::TwoPi();
325 Int_t nBinsArea = 100;
326 Double_t minArea = 0.;
327 Double_t maxArea = 1.;
328
329 fh3PtEtaPhiJetFull = new TH3F("fh3PtEtaPhiJetFull","fh3PtEtaPhiJetFull;#it{p}_{T}^{jet};#eta;#varphi",nBinsPt,minPt,maxPt,nBinsEta,minEta,maxEta,nBinsPhi,minPhi,maxPhi);
330 fOutput->Add(fh3PtEtaPhiJetFull);
331
332 fh3PtEtaPhiJetCharged = new TH3F("fh3PtEtaPhiJetCharged","fh3PtEtaPhiJetCharged;#it{p}_{T}^{jet};#eta;#varphi",nBinsPt,minPt,maxPt,nBinsEta,minEta,maxEta,nBinsPhi,minPhi,maxPhi);
333 fOutput->Add(fh3PtEtaPhiJetCharged);
334
335 fh2NJetsPtFull = new TH2F("fh2NJetsPtFull","fh2NJetsPtFull;N_{jets};#it{p}_{T}^{jet}",20,-0.5,19.5,nBinsPt,minPt,maxPt);
336 fOutput->Add(fh2NJetsPtFull);
337
338 fh2NJetsPtCharged = new TH2F("fh2NJetsPtCharged","fh2NJetsPtCharged;N_{jets};#it{p}_{T}^{jet}",20,-0.5,19.5,nBinsPt,minPt,maxPt);
339 fOutput->Add(fh2NJetsPtCharged);
340
341 fh3PtEtaAreaJetFull = new TH3F("fh3PtEtaAreaJetFull","fh3PtEtaAreaJetFull;#it{p}_{T}^{jet};#eta;A",nBinsPt,minPt,maxPt,nBinsEta,minEta,maxEta,nBinsArea,minArea,maxArea);
342 fOutput->Add(fh3PtEtaAreaJetFull);
343
344 fh3PtEtaAreaJetCharged = new TH3F("fh3PtEtaAreaJetCharged","fh3PtEtaAreaJetCharged;#it{p}_{T}^{jet};#eta;A",nBinsPt,minPt,maxPt,nBinsEta,minEta,maxEta,nBinsArea,minArea,maxArea);
345 fOutput->Add(fh3PtEtaAreaJetCharged);
346
347 Int_t nBinsConst =100;
348 Double_t minConst = 0.;
349 Double_t maxConst = 100.;
350
351 Int_t nBinsMeanPt = 200;
352 Double_t minMeanPt = 0.;
353 Double_t maxMeanPt = 20.;
354
355 Int_t nBinsNEF = 101;
356 Double_t minNEF = 0.;
357 Double_t maxNEF = 1.01;
358
359 Int_t nBinsz = 101;
360 Double_t minz = 0.;
361 Double_t maxz = 1.01;
362
363 Int_t nBinsECluster =100;
364 Double_t minECluster = 0.;
365 Double_t maxECluster = 100.;
366
367
368 fh2PtNConstituentsCharged = new TH2F("fh2PtNConstituentsCharged","fh2PtNConstituentsCharged;#it{p}_{T}^{jet};N_{charged constituents}",nBinsPt,minPt,maxPt,nBinsConst,minConst,maxConst);
369 fOutput->Add(fh2PtNConstituentsCharged);
370
371 fh2PtNConstituents = new TH2F("fh2PtNConstituents","fh2PtNConstituents;#it{p}_{T}^{jet};N_{constituents}",nBinsPt,minPt,maxPt,nBinsConst,minConst,maxConst);
372 fOutput->Add(fh2PtNConstituents);
373
374 fh2PtMeanPtConstituentsCharged = new TH2F("fh2PtMeanPtConstituentsCharged","fh2PtMeanPtConstituentsCharged;#it{p}_{T}^{jet};charged #langle #it{p}_{T} #rangle",nBinsPt,minPt,maxPt,nBinsMeanPt,minMeanPt,maxMeanPt);
375 fOutput->Add(fh2PtMeanPtConstituentsCharged);
376
377 fh2PtMeanPtConstituentsNeutral = new TH2F("fh2PtMeanPtConstituentsNeutral","fh2PtMeanPtConstituentsNeutral;#it{p}_{T}^{jet};neutral langle #it{p}_{T} #rangle",nBinsPt,minPt,maxPt,nBinsMeanPt,minMeanPt,maxMeanPt);
378 fOutput->Add(fh2PtMeanPtConstituentsNeutral);
379
380 fh2PtNEF = new TH2F("fh2PtNEF","fh2PtNEF;#it{p}_{T}^{jet};NEF",nBinsPt,minPt,maxPt,nBinsNEF,minNEF,maxNEF);
381 fOutput->Add(fh2PtNEF);
382
3c709670 383 fh3NEFEtaPhi = new TH3F("fh3NEFEtaPhi","fh3NEFEtaPhi;NEF;#eta;#varphi",nBinsNEF,minNEF,maxNEF,nBinsEta,minEta,maxEta,nBinsPhi,minPhi,maxPhi);
384 fOutput->Add(fh3NEFEtaPhi);
385
386 fh2NEFNConstituentsCharged = new TH2F("fh2NEFNConstituentsCharged","fh2NEFNConstituentsCharged;NEF;N_{charged constituents}",nBinsNEF,minNEF,maxNEF,nBinsConst,minConst,maxConst);
387 fOutput->Add(fh2NEFNConstituentsCharged);
388
389 fh2NEFNConstituentsNeutral = new TH2F("fh2NEFNConstituentsNeutral","fh2NEFNConstituentsNeutral;NEF;N_{clusters}",nBinsNEF,minNEF,maxNEF,nBinsConst,minConst,maxConst);
390 fOutput->Add(fh2NEFNConstituentsNeutral);
391
0531c107 392 fh2Ptz = new TH2F("fh2Ptz","fh2Ptz;#it{p}_{T}^{jet};z=p_{t,trk}^{proj}/p_{jet}",nBinsPt,minPt,maxPt,nBinsz,minz,maxz);
393 fOutput->Add(fh2Ptz);
394
3c709670 395 fh2PtzCharged = new TH2F("fh2PtzCharged","fh2Ptz;#it{p}_{T}^{ch jet};z=p_{t,trk}^{proj}/p_{ch jet}",nBinsPt,minPt,maxPt,nBinsz,minz,maxz);
396 fOutput->Add(fh2PtzCharged);
397
0531c107 398 fh2PtLeadJet1VsLeadJet2 = new TH2F("fh2PtLeadJet1VsLeadJet2","fh2PtLeadJet1VsLeadJet2;#it{p}_{T}^{jet 1};#it{p}_{T}^{jet 2}",nBinsPt,minPt,maxPt,nBinsPt,minPt,maxPt);
399 fOutput->Add(fh2PtLeadJet1VsLeadJet2);
400
992189ac 401 fh3EEtaPhiCluster = new TH3F("fh3EEtaPhiCluster","fh3EEtaPhiCluster;E_{clus};#eta;#phi",nBinsECluster,minECluster,maxECluster,nBinsEta,minEta,maxEta,nBinsPhi,minPhi,maxPhi);
0531c107 402 fOutput->Add(fh3EEtaPhiCluster);
403
404 fh3PtLeadJet1VsPatchEnergy = new TH3F("fh3PtLeadJet1VsPatchEnergy","fh3PtLeadJet1VsPatchEnergy;#it{p}_{T}^{jet 1};Amplitude_{patch};trig type",nBinsPt,minPt,maxPt,nBinsPt,minPt,maxPt,2,-0.5,1.5);
405 fOutput->Add(fh3PtLeadJet1VsPatchEnergy);
406 fh3PtLeadJet2VsPatchEnergy = new TH3F("fh3PtLeadJet2VsPatchEnergy","fh3PtLeadJet2VsPatchEnergy;#it{p}_{T}^{jet 1};Amplitude_{patch};trig type",nBinsPt,minPt,maxPt,nBinsPt,minPt,maxPt,2,-0.5,1.5);
407 fOutput->Add(fh3PtLeadJet2VsPatchEnergy);
408
409 fh3PatchEnergyEtaPhiCenter = new TH3F("fh3PatchEnergyEtaPhiCenter","fh3PatchEnergyEtaPhiCenter;E_{patch};#eta;#phi",nBinsPt,minPt,maxPt,nBinsEta,minEta,maxEta,nBinsPhi,minPhi,maxPhi);
410 fOutput->Add(fh3PatchEnergyEtaPhiCenter);
411
992189ac 412 fh2CellEnergyVsTime = new TH2F("fh2CellEnergyVsTime","fh2CellEnergyVsTime;E_{cell};time",100,0.,100.,700,-400,1000);
413 fOutput->Add(fh2CellEnergyVsTime);
414
0531c107 415
416 // =========== Switch on Sumw2 for all histos ===========
417 for (Int_t i=0; i<fOutput->GetEntries(); ++i) {
418 TH1 *h1 = dynamic_cast<TH1*>(fOutput->At(i));
419 if (h1){
420 h1->Sumw2();
421 continue;
422 }
423 TH2 *h2 = dynamic_cast<TH2*>(fOutput->At(i));
424 if (h2){
425 h2->Sumw2();
426 continue;
427 }
428 TH3 *h3 = dynamic_cast<TH3*>(fOutput->At(i));
429 if (h3){
430 h3->Sumw2();
431 continue;
432 }
433 THnSparse *hn = dynamic_cast<THnSparse*>(fOutput->At(i));
434 if(hn)hn->Sumw2();
435 }
436
437 TH1::AddDirectory(oldStatus);
438
439 PostData(1, fOutput); // Post data for ALL output slots > 0 here.
440}
441
442//________________________________________________________________________
443Bool_t AliAnalysisTaskEmcalJetTriggerQA::FillHistograms()
444{
445 // Fill histograms.
0531c107 446
9e5eee5d 447 AliClusterContainer *clusCont = GetClusterContainer(0);
448
449 if (clusCont) {
450 Int_t nclusters = clusCont->GetNClusters();
0531c107 451 for (Int_t ic = 0; ic < nclusters; ic++) {
9e5eee5d 452 AliVCluster *cluster = static_cast<AliVCluster*>(clusCont->GetCluster(ic));
0531c107 453 if (!cluster) {
454 AliError(Form("Could not receive cluster %d", ic));
455 continue;
456 }
9e5eee5d 457 if (!cluster->IsEMCAL()) {
458 AliDebug(11,Form("%s: Cluster is not emcal",GetName()));
459 continue;
460 }
0531c107 461
462 TLorentzVector lp;
463 cluster->GetMomentum(lp, const_cast<Double_t*>(fVertex));
464
465 //Fill eta,phi,E of clusters here
466 fh3EEtaPhiCluster->Fill(lp.E(),lp.Eta(),lp.Phi());
467 }
468 }
469
992189ac 470 if(fCaloCells) {
471 const Short_t nCells = fCaloCells->GetNumberOfCells();
472
473 for(Int_t iCell=0; iCell<nCells; ++iCell) {
474 Short_t cellId = fCaloCells->GetCellNumber(iCell);
475 Double_t cellE = fCaloCells->GetCellAmplitude(cellId);
476 Double_t cellT = fCaloCells->GetCellTime(cellId);
477
478 AliDebug(2,Form("cell energy = %f time = %f",cellE,cellT*1e9));
479 fh2CellEnergyVsTime->Fill(cellE,cellT*1e9);
480
481 }
482
483 }
484
0531c107 485 Double_t ptLeadJet1 = 0.;
486 Double_t ptLeadJet2 = 0.;
487
488 TArrayI *nJetsArr = new TArrayI(fh2NJetsPtFull->GetNbinsY()+1);
489 nJetsArr->Reset(0);
490 nJetsArr->Set(fh2NJetsPtFull->GetNbinsY()+1);
491
9e5eee5d 492 if (GetJetContainer(fContainerFull)) {
493 const Int_t njets = GetNJets(fContainerFull);
0531c107 494 for (Int_t ij = 0; ij < njets; ij++) {
495
9e5eee5d 496 AliEmcalJet* jet = GetAcceptJetFromArray(ij,fContainerFull);
497 if (!jet)
498 continue; //jet not selected
499
0531c107 500 Double_t jetPt = jet->Pt();
501 if(jetPt>ptLeadJet1) ptLeadJet1=jetPt;
502 fh3PtEtaPhiJetFull->Fill(jetPt,jet->Eta(),jet->Phi());
503 fh3PtEtaAreaJetFull->Fill(jetPt,jet->Eta(),jet->Area());
504
505 //count jets above certain pT threshold
506 Int_t ptbin = fh2NJetsPtFull->GetYaxis()->FindBin(jetPt);
507 for(Int_t iptbin = ptbin; iptbin<=fh2NJetsPtFull->GetNbinsY(); iptbin++)
508 nJetsArr->AddAt(nJetsArr->At(iptbin)+1,iptbin);
509
510 fh2PtNConstituentsCharged->Fill(jetPt,jet->GetNumberOfTracks());
511 fh2PtNConstituents->Fill(jetPt,jet->GetNumberOfConstituents());
512
513 fh2PtNEF->Fill(jetPt,jet->NEF());
3c709670 514 fh3NEFEtaPhi->Fill(jet->NEF(),jet->Eta(),jet->Phi());
515 fh2NEFNConstituentsCharged->Fill(jet->NEF(),jet->GetNumberOfTracks());
516 fh2NEFNConstituentsNeutral->Fill(jet->NEF(),jet->GetNumberOfClusters());
0531c107 517
518 AliVParticle *vp;
519 Double_t sumPtCh = 0.;
520 for(Int_t icc=0; icc<jet->GetNumberOfTracks(); icc++) {
521 vp = static_cast<AliVParticle*>(jet->TrackAt(icc, fTracks));
522 if(!vp) continue;
523 fh2Ptz->Fill(jetPt,GetZ(vp,jet));
0531c107 524 sumPtCh+=vp->Pt();
0531c107 525 }
526
527 if(jet->GetNumberOfTracks()>0)
528 fh2PtMeanPtConstituentsCharged->Fill(jetPt,sumPtCh/(double)(jet->GetNumberOfTracks()) );
529
530
531 AliVCluster *vc = 0x0;
532 Double_t sumPtNe = 0.;
9e5eee5d 533 if (clusCont) {
534 for(Int_t icc=0; icc<jet->GetNumberOfClusters(); icc++) {
535 vc = static_cast<AliVCluster*>(clusCont->GetCluster(icc));
536 if(!vc) continue;
537
538 TLorentzVector lp;
539 vc->GetMomentum(lp, const_cast<Double_t*>(fVertex));
540 sumPtNe+=lp.Pt();
541
542 }
0531c107 543
9e5eee5d 544 if(jet->GetNumberOfClusters()>0)
545 fh2PtMeanPtConstituentsNeutral->Fill(jetPt,sumPtNe/(double)(jet->GetNumberOfClusters()) );
0531c107 546 }
0531c107 547 } //full jet loop
548
549 for(Int_t i=1; i<=fh2NJetsPtFull->GetNbinsY(); i++) {
550 Int_t nJetsInEvent = nJetsArr->At(i);
551 fh2NJetsPtFull->Fill(nJetsInEvent,fh2NJetsPtFull->GetYaxis()->GetBinCenter(i));
552 }
553
554 }
555
556 //Reset array to zero to also count charged jets
557 nJetsArr->Reset(0);
558
559 //Loop over charged jets
9e5eee5d 560 if (GetJetContainer(fContainerCharged)) {
561 const Int_t njets = GetNJets(fContainerCharged);
562 for (Int_t ij = 0; ij < njets; ij++) {
0531c107 563
9e5eee5d 564 AliEmcalJet* jet = GetAcceptJetFromArray(ij,fContainerCharged);
565 if (!jet)
566 continue; //jet not selected
0531c107 567
0531c107 568 Double_t jetPt = jet->Pt();
569 if(jetPt>ptLeadJet2) ptLeadJet2=jetPt;
570 fh3PtEtaPhiJetCharged->Fill(jetPt,jet->Eta(),jet->Phi());
571 fh3PtEtaAreaJetCharged->Fill(jetPt,jet->Eta(),jet->Area());
3c709670 572
573 AliVParticle *vp;
574 for(Int_t icc=0; icc<jet->GetNumberOfTracks(); icc++) {
575 vp = static_cast<AliVParticle*>(jet->TrackAt(icc, fTracks));
576 if(!vp) continue;
577 fh2PtzCharged->Fill(jetPt,GetZ(vp,jet));
578 }
0531c107 579
580 //count jets above certain pT threshold
581 Int_t ptbin = fh2NJetsPtCharged->GetYaxis()->FindBin(jetPt);
582 for(Int_t iptbin = ptbin; iptbin<=fh2NJetsPtCharged->GetNbinsY(); iptbin++)
583 nJetsArr->AddAt(nJetsArr->At(iptbin)+1,iptbin);
584
585 }//ch jet loop
586 for(Int_t i=1; i<=fh2NJetsPtCharged->GetNbinsY(); i++) {
587 Int_t nJetsInEvent = nJetsArr->At(i);
588 fh2NJetsPtCharged->Fill(nJetsInEvent,fh2NJetsPtCharged->GetYaxis()->GetBinCenter(i));
589 }
590 }
591
9e5eee5d 592 if(GetJetContainer(fContainerFull) && GetJetContainer(fContainerCharged)) {
0531c107 593 fh2PtLeadJet1VsLeadJet2->Fill(ptLeadJet1,ptLeadJet2);
594 }
595
596 fh3PtLeadJet1VsPatchEnergy->Fill(ptLeadJet1,fMaxPatchEnergy,fTriggerType);
597 fh3PtLeadJet2VsPatchEnergy->Fill(ptLeadJet2,fMaxPatchEnergy,fTriggerType);
598
599 if(nJetsArr) delete nJetsArr;
600
601 return kTRUE;
602}
603
604//________________________________________________________________________
605Bool_t AliAnalysisTaskEmcalJetTriggerQA::Run()
606{
607 // Run analysis code here, if needed. It will be executed before FillHistograms().
608
609 //Check if event is selected (vertex & pile-up)
610 if(!SelectEvent())
611 return kFALSE;
612
613 LoadExtraBranches();
614
609c7a0d 615 if(!fTriggerClass.IsNull())
616 FindTriggerPatch();
0531c107 617
618 return kTRUE; // If return kFALSE FillHistogram() will NOT be executed.
619}
620
621//_______________________________________________________________________
622void AliAnalysisTaskEmcalJetTriggerQA::Terminate(Option_t *)
623{
624 // Called once at the end of the analysis.
625}
626//________________________________________________________________________
627Double_t AliAnalysisTaskEmcalJetTriggerQA::GetZ(const AliVParticle *trk, const AliEmcalJet *jet) const
628{
629 // Get Z of constituent trk
630
631 return GetZ(trk->Px(),trk->Py(),trk->Pz(),jet->Px(),jet->Py(),jet->Pz());
632}
633
634//________________________________________________________________________
635Double_t AliAnalysisTaskEmcalJetTriggerQA::GetZ(const Double_t trkPx, const Double_t trkPy, const Double_t trkPz, const Double_t jetPx, const Double_t jetPy, const Double_t jetPz) const
636{
637 //
638 // Get the z of a constituent inside of a jet
639 //
640
641 Double_t pJetSq = jetPx*jetPx+jetPy*jetPy+jetPz*jetPz;
642
643 if(pJetSq>0.)
644 return (trkPx*jetPx+trkPy*jetPy+trkPz*jetPz)/pJetSq;
645 else {
646 AliWarning(Form("%s: strange, pjet*pjet seems to be zero pJetSq: %f",GetName(), pJetSq));
647 return 0;
648 }
649}