]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGJE/EMCALJetTasks/UserTasks/AliAnalysisTaskFullpAJets.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
move setter to public
[u/mrichter/AliRoot.git] / PWGJE / EMCALJetTasks / UserTasks / AliAnalysisTaskFullpAJets.cxx
CommitLineData
ac6a3f1e 1#include "AliAnalysisTaskFullpAJets.h"
2
3#include <Riostream.h>
4#include <ctime>
5#include <TString.h>
6#include <TChain.h>
7#include <TTree.h>
1cd88f63
ML		 8#include <TH1F.h>
9#include <TH2F.h>
10#include <TH3F.h>
2d2100d5 11#include <THnSparse.h>
ac6a3f1e 12#include <TCanvas.h>
13#include <TList.h>
14#include <TLorentzVector.h>
15#include <TProfile.h>
16#include <TProfile2D.h>
17#include <TProfile3D.h>
18#include <TRandom.h>
19#include <TRandom3.h>
c6202663 20#include <TClonesArray.h>
21#include <TObjArray.h>
ac6a3f1e 22
23#include "AliAnalysisTaskSE.h"
24#include "AliAnalysisManager.h"
25#include "AliStack.h"
26#include "AliESDtrackCuts.h"
27#include "AliESDEvent.h"
28#include "AliESDInputHandler.h"
29#include "AliAODEvent.h"
08b981da 30#include "AliVEvent.h"
ac6a3f1e 31#include "AliMCEvent.h"
c6202663 32#include "AliVTrack.h"
33#include "AliVCluster.h"
ac6a3f1e 34#include "AliEmcalJet.h"
35#include "AliEMCALGeometry.h"
08b981da 36#include "AliEMCALRecoUtils.h"
37#include "AliVCaloCells.h"
24a61909 38#include "AliPicoTrack.h"
c6202663 39#include "Rtypes.h"
ac6a3f1e 40
41ClassImp(AliAnalysisTaskFullpAJets)
42
43//________________________________________________________________________
44AliAnalysisTaskFullpAJets::AliAnalysisTaskFullpAJets() :
45 AliAnalysisTaskSE(),
46
47 fOutput(0),
62620fff 48 flTrack(0),
49 flCluster(0),
ac6a3f1e 50 fhTrackPt(0),
51 fhTrackEta(0),
52 fhTrackPhi(0),
24a61909 53 fhGlobalTrackPt(0),
54 fhGlobalTrackEta(0),
55 fhGlobalTrackPhi(0),
56 fhComplementaryTrackPt(0),
57 fhComplementaryTrackEta(0),
58 fhComplementaryTrackPhi(0),
ac6a3f1e 59 fhClusterPt(0),
60 fhClusterEta(0),
61 fhClusterPhi(0),
c54b626a 62 fhCentrality(0),
ac6a3f1e 63 fhEMCalCellCounts(0),
78246241 64
ac6a3f1e 65 fhTrackEtaPhi(0),
24a61909 66 fhTrackPhiPt(0),
67 fhTrackEtaPt(0),
68 fhGlobalTrackEtaPhi(0),
69 fhGlobalTrackPhiPt(0),
70 fhGlobalTrackEtaPt(0),
71 fhComplementaryTrackEtaPhi(0),
72 fhComplementaryTrackPhiPt(0),
73 fhComplementaryTrackEtaPt(0),
ac6a3f1e 74 fhClusterEtaPhi(0),
24a61909 75 fhClusterPhiPt(0),
76 fhClusterEtaPt(0),
3e43a01f 77 fhEMCalEventMult(0),
78 fhTPCEventMult(0),
79 fhEMCalTrackEventMult(0),
c6202663 80
81 fpEMCalEventMult(0),
82 fpTPCEventMult(0),
78246241 83
ac6a3f1e 84 fpTrackPtProfile(0),
c54b626a 85 fpClusterPtProfile(0),
86
3e43a01f 87 fTPCRawJets(0),
88 fEMCalRawJets(0),
c6202663 89 fRhoChargedCMSScale(0),
62620fff 90 fRhoChargedScale(0),
91 fRhoFull0(0),
92 fRhoFull1(0),
93 fRhoFull2(0),
94 fRhoFullN(0),
95 fRhoFullDijet(0),
96 fRhoFullkT(0),
97 fRhoFullCMS(0),
98 fRhoCharged0(0),
99 fRhoCharged1(0),
100 fRhoCharged2(0),
101 fRhoChargedN(0),
102 fRhoChargedkT(0),
103 fRhoChargedkTScale(0),
104 fRhoChargedCMS(0),
c6202663 105
106 fTPCJet(0),
107 fTPCFullJet(0),
108 fTPCOnlyJet(0),
8daeee93 109 fTPCJetUnbiased(0),
c6202663 110 fTPCkTFullJet(0),
111 fEMCalJet(0),
112 fEMCalFullJet(0),
113 fEMCalPartJet(0),
8daeee93 114 fEMCalPartJetUnbiased(0),
c6202663 115 fEMCalkTFullJet(0),
116
c54b626a 117 fIsInitialized(0),
c6202663 118 fRJET(4),
c54b626a 119 fnEvents(0),
120 fnEventsCharged(0),
121 fnDiJetEvents(0),
08b981da 122 fEvent(0),
123 fRecoUtil(0),
124 fEMCALGeometry(0),
125 fCells(0),
87a5edfe 126 fDoNEF(0),
d2af75be 127 fDoNEFSignalOnly(0),
3e43a01f 128 fSignalTrackBias(0),
62620fff 129 fTrackQA(0),
130 fClusterQA(0),
131 fCalculateRhoJet(0),
c6202663 132 fEMCalPhiMin(1.39626),
133 fEMCalPhiMax(3.26377),
134 fEMCalPhiTotal(1.86750),
135 fEMCalEtaMin(-0.7),
136 fEMCalEtaMax(0.7),
137 fEMCalEtaTotal(1.4),
138 fEMCalArea(2.61450),
c54b626a 139 fTPCPhiMin(0),
c6202663 140 fTPCPhiMax(6.28319),
141 fTPCPhiTotal(6.28319),
142 fTPCEtaMin(-0.9),
143 fTPCEtaMax(0.9),
144 fTPCEtaTotal(1.8),
145 fTPCArea(11.30973),
24a61909 146 fParticlePtLow(0.0),
147 fParticlePtUp(200.0),
20f2d13a 148 fParticlePtBins(200),
c6202663 149 fJetR(0.4),
150 fJetRForRho(0.5),
151 fJetAreaCutFrac(0.6),
152 fJetAreaThreshold(0.30159),
153 fnEMCalCells(12288),
154 fScaleFactor(1.50),
7acc3e04 155 fNColl(7),
c6202663 156 fTrackMinPt(0.15),
157 fClusterMinPt(0.3),
91d0893e 158 fNEFSignalJetCut(0.9),
3aae667c 159 fCentralityTag("V0A"),
c6202663 160 fCentralityBins(10),
c54b626a 161 fCentralityLow(0),
c6202663 162 fCentralityUp(100),
c54b626a 163 fEventCentrality(0),
c6202663 164 fRhoFull(0),
c54b626a 165 fRhoCharged(0),
c54b626a 166 fnTracks(0),
167 fnClusters(0),
08b981da 168 fnCaloClusters(0),
c54b626a 169 fnAKTFullJets(0),
170 fnAKTChargedJets(0),
171 fnKTFullJets(0),
78246241 172 fnKTChargedJets(0),
c54b626a 173 fnBckgClusters(0),
c6202663 174 fTPCJetThreshold(5),
175 fEMCalJetThreshold(5),
176 fVertexWindow(10),
177 fVertexMaxR(1),
d812e269 178 fTrackName(0),
179 fClusName(0),
180 fkTChargedName(0),
181 fAkTChargedName(0),
182 fkTFullName(0),
183 fAkTFullName(0),
c6202663 184 fOrgTracks(0),
185 fOrgClusters(0),
c54b626a 186 fmyAKTFullJets(0),
187 fmyAKTChargedJets(0),
188 fmyKTFullJets(0),
78246241 189 fmyKTChargedJets(0),
c6202663 190 fmyTracks(0),
191 fmyClusters(0),
192 fEMCalRCBckgFluc(0),
193 fTPCRCBckgFluc(0),
194 fEMCalRCBckgFlucSignal(0),
7acc3e04 195 fTPCRCBckgFlucSignal(0),
196 fEMCalRCBckgFlucNColl(0),
197 fTPCRCBckgFlucNColl(0)
ac6a3f1e 198{
199 // Dummy constructor ALWAYS needed for I/O.
8daeee93 200 fVertex[0]=0.0,fVertex[1]=0.0,fVertex[2]=0.0;
ac6a3f1e 201}
202
203//________________________________________________________________________
204AliAnalysisTaskFullpAJets::AliAnalysisTaskFullpAJets(const char *name) :
205 AliAnalysisTaskSE(name),
206
207 fOutput(0),
62620fff 208 flTrack(0),
209 flCluster(0),
ac6a3f1e 210 fhTrackPt(0),
211 fhTrackEta(0),
212 fhTrackPhi(0),
24a61909 213 fhGlobalTrackPt(0),
214 fhGlobalTrackEta(0),
215 fhGlobalTrackPhi(0),
216 fhComplementaryTrackPt(0),
217 fhComplementaryTrackEta(0),
218 fhComplementaryTrackPhi(0),
ac6a3f1e 219 fhClusterPt(0),
220 fhClusterEta(0),
221 fhClusterPhi(0),
c54b626a 222 fhCentrality(0),
ac6a3f1e 223 fhEMCalCellCounts(0),
78246241 224
ac6a3f1e 225 fhTrackEtaPhi(0),
24a61909 226 fhTrackPhiPt(0),
227 fhTrackEtaPt(0),
228 fhGlobalTrackEtaPhi(0),
229 fhGlobalTrackPhiPt(0),
230 fhGlobalTrackEtaPt(0),
231 fhComplementaryTrackEtaPhi(0),
232 fhComplementaryTrackPhiPt(0),
233 fhComplementaryTrackEtaPt(0),
ac6a3f1e 234 fhClusterEtaPhi(0),
24a61909 235 fhClusterPhiPt(0),
236 fhClusterEtaPt(0),
3e43a01f 237 fhEMCalEventMult(0),
238 fhTPCEventMult(0),
239 fhEMCalTrackEventMult(0),
c6202663 240
241 fpEMCalEventMult(0),
242 fpTPCEventMult(0),
78246241 243
ac6a3f1e 244 fpTrackPtProfile(0),
c54b626a 245 fpClusterPtProfile(0),
246
c6202663 247 fTPCRawJets(0),
248 fEMCalRawJets(0),
c6202663 249 fRhoChargedCMSScale(0),
62620fff 250 fRhoChargedScale(0),
251 fRhoFull0(0),
252 fRhoFull1(0),
253 fRhoFull2(0),
254 fRhoFullN(0),
255 fRhoFullDijet(0),
256 fRhoFullkT(0),
257 fRhoFullCMS(0),
258 fRhoCharged0(0),
259 fRhoCharged1(0),
260 fRhoCharged2(0),
261 fRhoChargedN(0),
262 fRhoChargedkT(0),
263 fRhoChargedkTScale(0),
264 fRhoChargedCMS(0),
c6202663 265
266 fTPCJet(0),
267 fTPCFullJet(0),
268 fTPCOnlyJet(0),
8daeee93 269 fTPCJetUnbiased(0),
c6202663 270 fTPCkTFullJet(0),
271 fEMCalJet(0),
272 fEMCalFullJet(0),
273 fEMCalPartJet(0),
8daeee93 274 fEMCalPartJetUnbiased(0),
c6202663 275 fEMCalkTFullJet(0),
276
c54b626a 277 fIsInitialized(0),
c6202663 278 fRJET(4),
c54b626a 279 fnEvents(0),
280 fnEventsCharged(0),
281 fnDiJetEvents(0),
08b981da 282 fEvent(0),
283 fRecoUtil(0),
284 fEMCALGeometry(0),
285 fCells(0),
87a5edfe 286 fDoNEF(0),
d2af75be 287 fDoNEFSignalOnly(0),
3e43a01f 288 fSignalTrackBias(0),
62620fff 289 fTrackQA(0),
290 fClusterQA(0),
291 fCalculateRhoJet(0),
c6202663 292 fEMCalPhiMin(1.39626),
293 fEMCalPhiMax(3.26377),
294 fEMCalPhiTotal(1.86750),
295 fEMCalEtaMin(-0.7),
296 fEMCalEtaMax(0.7),
297 fEMCalEtaTotal(1.4),
298 fEMCalArea(2.61450),
c54b626a 299 fTPCPhiMin(0),
c6202663 300 fTPCPhiMax(6.28319),
301 fTPCPhiTotal(6.28319),
302 fTPCEtaMin(-0.9),
303 fTPCEtaMax(0.9),
304 fTPCEtaTotal(1.8),
305 fTPCArea(11.30973),
24a61909 306 fParticlePtLow(0.0),
307 fParticlePtUp(200.0),
308 fParticlePtBins(2000),
c6202663 309 fJetR(0.4),
310 fJetRForRho(0.5),
311 fJetAreaCutFrac(0.6),
312 fJetAreaThreshold(0.30159),
313 fnEMCalCells(12288),
314 fScaleFactor(1.50),
7acc3e04 315 fNColl(7),
c6202663 316 fTrackMinPt(0.15),
317 fClusterMinPt(0.3),
91d0893e 318 fNEFSignalJetCut(0.9),
3aae667c 319 fCentralityTag("V0A"),
c6202663 320 fCentralityBins(10),
c54b626a 321 fCentralityLow(0),
c6202663 322 fCentralityUp(100),
c54b626a 323 fEventCentrality(0),
c6202663 324 fRhoFull(0),
c54b626a 325 fRhoCharged(0),
c54b626a 326 fnTracks(0),
327 fnClusters(0),
08b981da 328 fnCaloClusters(0),
c54b626a 329 fnAKTFullJets(0),
330 fnAKTChargedJets(0),
331 fnKTFullJets(0),
78246241 332 fnKTChargedJets(0),
c54b626a 333 fnBckgClusters(0),
c6202663 334 fTPCJetThreshold(5),
335 fEMCalJetThreshold(5),
336 fVertexWindow(10),
337 fVertexMaxR(1),
d812e269 338 fTrackName(0),
339 fClusName(0),
340 fkTChargedName(0),
341 fAkTChargedName(0),
342 fkTFullName(0),
343 fAkTFullName(0),
c6202663 344 fOrgTracks(0),
345 fOrgClusters(0),
c54b626a 346 fmyAKTFullJets(0),
347 fmyAKTChargedJets(0),
348 fmyKTFullJets(0),
78246241 349 fmyKTChargedJets(0),
c6202663 350 fmyTracks(0),
351 fmyClusters(0),
352 fEMCalRCBckgFluc(0),
353 fTPCRCBckgFluc(0),
354 fEMCalRCBckgFlucSignal(0),
7acc3e04 355 fTPCRCBckgFlucSignal(0),
356 fEMCalRCBckgFlucNColl(0),
357 fTPCRCBckgFlucNColl(0)
ac6a3f1e 358{
359 // Constructor
360 // Define input and output slots here (never in the dummy constructor)
361 // Input slot #0 works with a TChain - it is connected to the default input container
362 // Output slot #1 writes into a TH1 container
8daeee93 363 fVertex[0]=0.0,fVertex[1]=0.0,fVertex[2]=0.0;
ac6a3f1e 364
c54b626a 365 DefineOutput(1,TList::Class()); // for output list
ac6a3f1e 366}
367
368//________________________________________________________________________
369AliAnalysisTaskFullpAJets::~AliAnalysisTaskFullpAJets()
370{
371 // Destructor. Clean-up the output list, but not the histograms that are put inside
372 // (the list is owner and will clean-up these histograms). Protect in PROOF case.
373 if (fOutput && !AliAnalysisManager::GetAnalysisManager()->IsProofMode())
374 {
375 delete fOutput;
376 }
ac6a3f1e 377}
378
379//________________________________________________________________________
380void AliAnalysisTaskFullpAJets::UserCreateOutputObjects()
381{
382 // Create histograms
383 // Called once (on the worker node)
384 fIsInitialized=kFALSE;
385 fOutput = new TList();
386 fOutput->SetOwner(); // IMPORTANT!
387
ac6a3f1e 388 // Initialize Global Variables
389 fnEvents=0;
390 fnEventsCharged=0;
391 fnDiJetEvents=0;
ac6a3f1e 392
c54b626a 393 // fRJET=4 -> fJetR=0.4 && fRJET=25 -> fJetR=0.25, but for writing files, should be 4 and 25 respectively
394 if (fRJET>10)
395 {
396 fJetR=(Double_t)fRJET/100.0;
397 }
398 else
399 {
400 fJetR=(Double_t)fRJET/10.0;
401 }
c6202663 402 fJetRForRho=0.5;
403
ac6a3f1e 404 fEMCalPhiMin=(80/(double)360)*2*TMath::Pi();
405 fEMCalPhiMax=(187/(double)360)*2*TMath::Pi();
406 fEMCalPhiTotal= fEMCalPhiMax-fEMCalPhiMin;
407 fEMCalEtaMin=-0.7;
408 fEMCalEtaMax=0.7;
409 fEMCalEtaTotal=fEMCalEtaMax-fEMCalEtaMin;
410 fEMCalArea=fEMCalPhiTotal*fEMCalEtaTotal;
411
412 fTPCPhiMin=(0/(double)360)*2*TMath::Pi();
413 fTPCPhiMax=(360/(double)360)*2*TMath::Pi();
414 fTPCPhiTotal= fTPCPhiMax-fTPCPhiMin;
415 fTPCEtaMin=-0.9;
416 fTPCEtaMax=0.9;
417 fTPCEtaTotal=fTPCEtaMax-fTPCEtaMin;
418 fTPCArea=fTPCPhiTotal*fTPCEtaTotal;
24a61909 419
420 fParticlePtLow=0.0;
421 fParticlePtUp=200.0;
422 fParticlePtBins=Int_t(fParticlePtUp-fParticlePtLow);
423
ac6a3f1e 424 fCentralityBins=10;
425 fCentralityLow=0.0;
426 fCentralityUp=100.0;
c54b626a 427 Int_t CentralityBinMult=10;
428
c6202663 429 fJetAreaCutFrac =0.6; // Fudge factor for selecting on jets with threshold Area or higher
ac6a3f1e 430 fJetAreaThreshold=fJetAreaCutFrac*TMath::Pi()*fJetR*fJetR;
c6202663 431 fTPCJetThreshold=5.0; // Threshold required for an Anti-kT Charged jet to be considered a "true" jet in TPC
432 fEMCalJetThreshold=5.0; // Threshold required for an Anti-kT Charged+Neutral jet to be considered a "true" jet in EMCal
ac6a3f1e 433 fVertexWindow=10.0;
434 fVertexMaxR=1.0;
435
c6202663 436 fnBckgClusters=1;
437 fEMCalRCBckgFluc = new Double_t[fnBckgClusters];
438 fTPCRCBckgFluc = new Double_t[fnBckgClusters];
439 fEMCalRCBckgFlucSignal = new Double_t[fnBckgClusters];
440 fTPCRCBckgFlucSignal = new Double_t[fnBckgClusters];
7acc3e04 441 fEMCalRCBckgFlucNColl = new Double_t[fnBckgClusters];
442 fTPCRCBckgFlucNColl = new Double_t[fnBckgClusters];
ac6a3f1e 443 for (Int_t i=0;i<fnBckgClusters;i++)
444 {
c6202663 445 fEMCalRCBckgFluc[i]=0.0;
446 fTPCRCBckgFluc[i]=0.0;
447 fEMCalRCBckgFlucSignal[i]=0.0;
448 fTPCRCBckgFlucSignal[i]=0.0;
7acc3e04 449 fEMCalRCBckgFlucNColl[i]=0.0;
450 fTPCRCBckgFlucNColl[i]=0.0;
ac6a3f1e 451 }
452
ac6a3f1e 453 fnEMCalCells=12288; // sMods 1-10 have 24x48 cells, sMods 11&12 have 8x48 cells...
454
ac6a3f1e 455 Int_t TCBins=100;
62620fff 456 Int_t multBins=200;
457 Double_t multLow=0;
458 Double_t multUp=200;
24a61909 459
1cd88f63 460 fhCentrality = new TH1F("fhCentrality","Event Centrality Distribution",fCentralityBins*CentralityBinMult,fCentralityLow,fCentralityUp);
c6202663 461 fhCentrality->GetXaxis()->SetTitle(fCentralityTag);
c54b626a 462 fhCentrality->GetYaxis()->SetTitle("1/N_{Events}");
463 fhCentrality->Sumw2();
464
62620fff 465 // Track QA Plots
466 if (fTrackQA==kTRUE)
467 {
468 flTrack = new TList();
469 flTrack->SetName("TrackQA");
470
471 // Hybrid Tracks
1cd88f63 472 fhTrackPt = new TH1F("fhTrackPt","p_{T} distribution of tracks in event",10*fParticlePtBins,fParticlePtLow,fParticlePtUp);
62620fff 473 fhTrackPt->GetXaxis()->SetTitle("p_{T} (GeV/c)");
474 fhTrackPt->GetYaxis()->SetTitle("1/N_{Events} dN/dp_{T}");
475 fhTrackPt->Sumw2();
476
1cd88f63 477 fhTrackPhi = new TH1F("fhTrackPhi","#varphi distribution of tracks in event",TCBins,fTPCPhiMin,fTPCPhiMax);
62620fff 478 fhTrackPhi->GetXaxis()->SetTitle("#varphi");
479 fhTrackPhi->GetYaxis()->SetTitle("1/N_{Events} dN/d#varphi");
480 fhTrackPhi->Sumw2();
481
1cd88f63 482 fhTrackEta = new TH1F("fhTrackEta","#eta distribution of tracks in event",TCBins,fTPCEtaMin,fTPCEtaMax);
62620fff 483 fhTrackEta->GetXaxis()->SetTitle("#eta");
484 fhTrackEta->GetYaxis()->SetTitle("1/N_{Events} dN/d#eta");
485 fhTrackEta->Sumw2();
486
1cd88f63 487 fhTrackEtaPhi = new TH2F("fhTrackEtaPhi","#eta-#varphi distribution of tracks in event",TCBins,fTPCEtaMin,fTPCEtaMax,TCBins,fTPCPhiMin,fTPCPhiMax);
62620fff 488 fhTrackEtaPhi->GetXaxis()->SetTitle("#eta");
489 fhTrackEtaPhi->GetYaxis()->SetTitle("#varphi");
490 fhTrackEtaPhi->GetZaxis()->SetTitle("1/N_{Events} dN/d#etad#varphi");
491 fhTrackEtaPhi->Sumw2();
492
1cd88f63 493 fhTrackPhiPt = new TH2F("fhTrackPhiPt","#varphi-p_{T} distribution of tracks in event",TCBins,fTPCPhiMin,fTPCPhiMax,fParticlePtBins,fParticlePtLow,fParticlePtUp);
62620fff 494 fhTrackPhiPt->GetXaxis()->SetTitle("#varphi");
495 fhTrackPhiPt->GetYaxis()->SetTitle("p_{T} (GeV/c)");
496 fhTrackPhiPt->GetZaxis()->SetTitle("1/N_{Events} dN/d#varphidp_{T}");
497 fhTrackPhiPt->Sumw2();
498
1cd88f63 499 fhTrackEtaPt = new TH2F("fhTrackEtaPt","#eta-p_{T} distribution of tracks in event",TCBins,fTPCEtaMin,fTPCEtaMax,fParticlePtBins,fParticlePtLow,fParticlePtUp);
62620fff 500 fhTrackEtaPt->GetXaxis()->SetTitle("#varphi");
501 fhTrackEtaPt->GetYaxis()->SetTitle("p_{T} (GeV/c)");
502 fhTrackEtaPt->GetZaxis()->SetTitle("1/N_{Events} dN/d#etadp_{T}");
503 fhTrackEtaPt->Sumw2();
1cd88f63 504
62620fff 505 // Global Tracks
1cd88f63 506 fhGlobalTrackPt = new TH1F("fhGlobalTrackPt","Global p_{T} distribution of tracks in event",10*fParticlePtBins,fParticlePtLow,fParticlePtUp);
62620fff 507 fhGlobalTrackPt->GetXaxis()->SetTitle("p_{T} (GeV/c)");
508 fhGlobalTrackPt->GetYaxis()->SetTitle("1/N_{Events} dN/dp_{T}");
509 fhGlobalTrackPt->Sumw2();
510
1cd88f63 511 fhGlobalTrackPhi = new TH1F("fhGlobalTrackPhi","Global #varphi distribution of tracks in event",TCBins,fTPCPhiMin,fTPCPhiMax);
62620fff 512 fhGlobalTrackPhi->GetXaxis()->SetTitle("#varphi");
513 fhGlobalTrackPhi->GetYaxis()->SetTitle("1/N_{Events} dN/d#varphi");
514 fhGlobalTrackPhi->Sumw2();
515
1cd88f63 516 fhGlobalTrackEta = new TH1F("fhGlobalTrackEta","Global #eta distribution of tracks in event",TCBins,fTPCEtaMin,fTPCEtaMax);
62620fff 517 fhGlobalTrackEta->GetXaxis()->SetTitle("#eta");
518 fhGlobalTrackEta->GetYaxis()->SetTitle("1/N_{Events} dN/d#eta");
519 fhGlobalTrackEta->Sumw2();
520
1cd88f63 521 fhGlobalTrackEtaPhi = new TH2F("fhGlobalTrackEtaPhi","Global #eta-#varphi distribution of tracks in event",TCBins,fTPCEtaMin,fTPCEtaMax,TCBins,fTPCPhiMin,fTPCPhiMax);
62620fff 522 fhGlobalTrackEtaPhi->GetXaxis()->SetTitle("#eta");
523 fhGlobalTrackEtaPhi->GetYaxis()->SetTitle("#varphi");
524 fhGlobalTrackEtaPhi->GetZaxis()->SetTitle("1/N_{Events} dN/d#etad#varphi");
525 fhGlobalTrackEtaPhi->Sumw2();
526
1cd88f63 527 fhGlobalTrackPhiPt = new TH2F("fhGlobalTrackPhiPt","Global #varphi-p_{T} distribution of tracks in event",TCBins,fTPCPhiMin,fTPCPhiMax,fParticlePtBins,fParticlePtLow,fParticlePtUp);
62620fff 528 fhGlobalTrackPhiPt->GetXaxis()->SetTitle("#varphi");
529 fhGlobalTrackPhiPt->GetYaxis()->SetTitle("p_{T} (GeV/c)");
530 fhGlobalTrackPhiPt->GetZaxis()->SetTitle("1/N_{Events} dN/d#varphidp_{T}");
531 fhGlobalTrackPhiPt->Sumw2();
532
1cd88f63 533 fhGlobalTrackEtaPt = new TH2F("fhGlobalTrackEtaPt","Global #eta-p_{T} distribution of tracks in event",TCBins,fTPCEtaMin,fTPCEtaMax,fParticlePtBins,fParticlePtLow,fParticlePtUp);
62620fff 534 fhGlobalTrackEtaPt->GetXaxis()->SetTitle("#varphi");
535 fhGlobalTrackEtaPt->GetYaxis()->SetTitle("p_{T} (GeV/c)");
536 fhGlobalTrackEtaPt->GetZaxis()->SetTitle("1/N_{Events} dN/d#etadp_{T}");
537 fhGlobalTrackEtaPt->Sumw2();
538
62620fff 539 // Complementary Tracks
1cd88f63 540 fhComplementaryTrackPt = new TH1F("fhComplementaryTrackPt","Complementary p_{T} distribution of tracks in event",10*fParticlePtBins,fParticlePtLow,fParticlePtUp);
62620fff 541 fhComplementaryTrackPt->GetXaxis()->SetTitle("p_{T} (GeV/c)");
542 fhComplementaryTrackPt->GetYaxis()->SetTitle("1/N_{Events} dN/dp_{T}");
543 fhComplementaryTrackPt->Sumw2();
544
1cd88f63 545 fhComplementaryTrackPhi = new TH1F("fhComplementaryTrackPhi","Complementary #varphi distribution of tracks in event",TCBins,fTPCPhiMin,fTPCPhiMax);
62620fff 546 fhComplementaryTrackPhi->GetXaxis()->SetTitle("#varphi");
547 fhComplementaryTrackPhi->GetYaxis()->SetTitle("1/N_{Events} dN/d#varphi");
548 fhComplementaryTrackPhi->Sumw2();
549
1cd88f63 550 fhComplementaryTrackEta = new TH1F("fhComplementaryTrackEta","Complementary #eta distribution of tracks in event",TCBins,fTPCEtaMin,fTPCEtaMax);
62620fff 551 fhComplementaryTrackEta->GetXaxis()->SetTitle("#eta");
552 fhComplementaryTrackEta->GetYaxis()->SetTitle("1/N_{Events} dN/d#eta");
553 fhComplementaryTrackEta->Sumw2();
554
1cd88f63 555 fhComplementaryTrackEtaPhi = new TH2F("fhComplementaryTrackEtaPhi","Complementary #eta-#varphi distribution of tracks in event",TCBins,fTPCEtaMin,fTPCEtaMax,TCBins,fTPCPhiMin,fTPCPhiMax);
62620fff 556 fhComplementaryTrackEtaPhi->GetXaxis()->SetTitle("#eta");
557 fhComplementaryTrackEtaPhi->GetYaxis()->SetTitle("#varphi");
558 fhComplementaryTrackEtaPhi->GetZaxis()->SetTitle("1/N_{Events} dN/d#etad#varphi");
559 fhComplementaryTrackEtaPhi->Sumw2();
560
1cd88f63 561 fhComplementaryTrackPhiPt = new TH2F("fhComplementaryTrackPhiPt","Complementary #varphi-p_{T} distribution of tracks in event",TCBins,fTPCPhiMin,fTPCPhiMax,fParticlePtBins,fParticlePtLow,fParticlePtUp);
62620fff 562 fhComplementaryTrackPhiPt->GetXaxis()->SetTitle("#varphi");
563 fhComplementaryTrackPhiPt->GetYaxis()->SetTitle("p_{T} (GeV/c)");
564 fhComplementaryTrackPhiPt->GetZaxis()->SetTitle("1/N_{Events} dN/d#varphidp_{T}");
565 fhComplementaryTrackPhiPt->Sumw2();
566
1cd88f63 567 fhComplementaryTrackEtaPt = new TH2F("fhComplementaryTrackEtaPt","Complementary #eta-p_{T} distribution of tracks in event",TCBins,fTPCEtaMin,fTPCEtaMax,fParticlePtBins,fParticlePtLow,fParticlePtUp);
62620fff 568 fhComplementaryTrackEtaPt->GetXaxis()->SetTitle("#varphi");
569 fhComplementaryTrackEtaPt->GetYaxis()->SetTitle("p_{T} (GeV/c)");
570 fhComplementaryTrackEtaPt->GetZaxis()->SetTitle("1/N_{Events} dN/d#etadp_{T}");
571 fhComplementaryTrackEtaPt->Sumw2();
572
1cd88f63 573 fhTPCEventMult = new TH2F("fhTPCEventMult","TPC Event Multiplcity vs Centrality",CentralityBinMult*fCentralityBins,fCentralityLow,fCentralityUp,multBins,multLow,multUp);
62620fff 574 fhTPCEventMult->GetXaxis()->SetTitle(fCentralityTag);
575 fhTPCEventMult->GetYaxis()->SetTitle("Multiplicity");
576 fhTPCEventMult->GetZaxis()->SetTitle("1/N_{Events} dN/dCentdN_{Charged}");
577 fhTPCEventMult->Sumw2();
578
1cd88f63 579 fhEMCalTrackEventMult = new TH2F("fhEMCalTrackEventMult","EMCal Track Event Multiplcity vs Centrality",CentralityBinMult*fCentralityBins,fCentralityLow,fCentralityUp,multBins,multLow,multUp);
62620fff 580 fhEMCalTrackEventMult->GetXaxis()->SetTitle(fCentralityTag);
581 fhEMCalTrackEventMult->GetYaxis()->SetTitle("Multiplicity");
582 fhEMCalTrackEventMult->GetZaxis()->SetTitle("1/N_{Events} dN/dCentdN_{Neutral}");
583 fhEMCalTrackEventMult->Sumw2();
584
585 fpTPCEventMult = new TProfile("fpTPCEventMult","TPC Event Multiplcity vs Centrality",CentralityBinMult*fCentralityBins,fCentralityLow,fCentralityUp);
586 fpTPCEventMult->GetXaxis()->SetTitle(fCentralityTag);
587 fpTPCEventMult->GetYaxis()->SetTitle("Multiplicity");
588
589 // QA::2D Energy Density Profiles for Tracks and Clusters
590 fpTrackPtProfile = new TProfile2D("fpTrackPtProfile","2D Profile of track pT density throughout the TPC",TCBins,fTPCEtaMin,fTPCEtaMax,TCBins,fTPCPhiMin,fTPCPhiMax);
591 fpTrackPtProfile->GetXaxis()->SetTitle("#eta");
592 fpTrackPtProfile->GetYaxis()->SetTitle("#varphi");
593 fpTrackPtProfile->GetZaxis()->SetTitle("p_{T} density (GeV/Area)");
594
595 flTrack->Add(fhTrackPt);
596 flTrack->Add(fhTrackEta);
597 flTrack->Add(fhTrackPhi);
598 flTrack->Add(fhTrackEtaPhi);
599 flTrack->Add(fhTrackPhiPt);
600 flTrack->Add(fhTrackEtaPt);
62620fff 601 flTrack->Add(fhGlobalTrackPt);
602 flTrack->Add(fhGlobalTrackEta);
603 flTrack->Add(fhGlobalTrackPhi);
604 flTrack->Add(fhGlobalTrackEtaPhi);
605 flTrack->Add(fhGlobalTrackPhiPt);
606 flTrack->Add(fhGlobalTrackEtaPt);
62620fff 607 flTrack->Add(fhComplementaryTrackPt);
608 flTrack->Add(fhComplementaryTrackEta);
609 flTrack->Add(fhComplementaryTrackPhi);
610 flTrack->Add(fhComplementaryTrackEtaPhi);
611 flTrack->Add(fhComplementaryTrackPhiPt);
612 flTrack->Add(fhComplementaryTrackEtaPt);
62620fff 613 flTrack->Add(fhTPCEventMult);
614 flTrack->Add(fhEMCalTrackEventMult);
615 flTrack->Add(fpTPCEventMult);
616 flTrack->Add(fpTrackPtProfile);
617 fOutput->Add(flTrack);
618 }
c54b626a 619
62620fff 620 if (fClusterQA==kTRUE)
621 {
622 flCluster = new TList();
623 flCluster->SetName("ClusterQA");
624
1cd88f63 625 fhClusterPt = new TH1F("fhClusterPt","p_{T} distribution of clusters in event",10*fParticlePtBins,fParticlePtLow,fParticlePtUp);
62620fff 626 fhClusterPt->GetXaxis()->SetTitle("p_{T} (GeV/c)");
627 fhClusterPt->GetYaxis()->SetTitle("1/N_{Events} dN/dp_{T}");
628 fhClusterPt->Sumw2();
629
1cd88f63 630 fhClusterPhi = new TH1F("fhClusterPhi","#varphi distribution of clusters in event",TCBins,fTPCPhiMin,fTPCPhiMax);
62620fff 631 fhClusterPhi->GetXaxis()->SetTitle("#varphi");
632 fhClusterPhi->GetYaxis()->SetTitle("1/N_{Events} dN/d#varphi");
633 fhClusterPhi->Sumw2();
634
1cd88f63 635 fhClusterEta = new TH1F("fhClusterEta","#eta distribution of clusters in event",TCBins,fTPCEtaMin,fTPCEtaMax);
62620fff 636 fhClusterEta->GetXaxis()->SetTitle("#eta");
637 fhClusterEta->GetYaxis()->SetTitle("1/N_{Events} dN/d#eta");
638 fhClusterEta->Sumw2();
639
1cd88f63 640 fhClusterEtaPhi = new TH2F("fhClusterEtaPhi","#eta-#varphi distribution of clusters in event",TCBins,fEMCalEtaMin,fEMCalEtaMax,TCBins,fEMCalPhiMin,fEMCalPhiMax);
62620fff 641 fhClusterEtaPhi->GetXaxis()->SetTitle("#eta");
642 fhClusterEtaPhi->GetYaxis()->SetTitle("#varphi");
643 fhClusterEtaPhi->GetZaxis()->SetTitle("1/N_{Events} dN/d#etad#varphi");
644 fhClusterEtaPhi->Sumw2();
645
1cd88f63 646 fhClusterPhiPt = new TH2F("fhClusterPhiPt","#varphi-p_{T} distribution of clusters in event",TCBins,fEMCalPhiMin,fEMCalPhiMax,fParticlePtBins,fParticlePtLow,fParticlePtUp);
62620fff 647 fhClusterPhiPt->GetXaxis()->SetTitle("#varphi");
648 fhClusterPhiPt->GetYaxis()->SetTitle("p_{T} (GeV/c)");
649 fhClusterPhiPt->GetZaxis()->SetTitle("1/N_{Events} dN/d#varphidp_{T}");
650 fhClusterPhiPt->Sumw2();
651
1cd88f63 652 fhClusterEtaPt = new TH2F("fhClusterEtaPt","#eta-p_{T} distribution of clusters in event",TCBins,fEMCalEtaMin,fEMCalEtaMax,fParticlePtBins,fParticlePtLow,fParticlePtUp);
62620fff 653 fhClusterEtaPt->GetXaxis()->SetTitle("#varphi");
654 fhClusterEtaPt->GetYaxis()->SetTitle("p_{T} (GeV/c)");
655 fhClusterEtaPt->GetZaxis()->SetTitle("1/N_{Events} dN/d#etadp_{T}");
656 fhClusterEtaPt->Sumw2();
1cd88f63
ML		 657
658 fhEMCalEventMult = new TH2F("fhEMCalEventMult","EMCal Event Multiplcity vs Centrality",CentralityBinMult*fCentralityBins,fCentralityLow,fCentralityUp,multBins,multLow,multUp);
62620fff 659 fhEMCalEventMult->GetXaxis()->SetTitle(fCentralityTag);
660 fhEMCalEventMult->GetYaxis()->SetTitle("Multiplicity");
661 fhEMCalEventMult->GetZaxis()->SetTitle("1/N_{Events} dN/dCentdN_{Neutral}");
662 fhEMCalEventMult->Sumw2();
663
664 fpEMCalEventMult = new TProfile("fpEMCalEventMult","EMCal Event Multiplcity vs Centrality",CentralityBinMult*fCentralityBins,fCentralityLow,fCentralityUp);
665 fpEMCalEventMult->GetXaxis()->SetTitle(fCentralityTag);
666 fpEMCalEventMult->GetYaxis()->SetTitle("Multiplicity");
667
668 fpClusterPtProfile = new TProfile2D("fpClusterPtProfile","2D Profile of cluster pT density throughout the EMCal",TCBins,fEMCalEtaMin,fEMCalEtaMax,TCBins,fEMCalPhiMin,fEMCalPhiMax);
669 fpClusterPtProfile->GetXaxis()->SetTitle("#eta");
670 fpClusterPtProfile->GetYaxis()->SetTitle("#varphi");
671 fpClusterPtProfile->GetZaxis()->SetTitle("p_{T} density (GeV/Area)");
672
1cd88f63 673 fhEMCalCellCounts = new TH1F("fhEMCalCellCounts","Distribtuion of cluster counts across the EMCal",fnEMCalCells,1,fnEMCalCells);
62620fff 674 fhEMCalCellCounts->GetXaxis()->SetTitle("Absoulute Cell Id");
675 fhEMCalCellCounts->GetYaxis()->SetTitle("Counts per Event");
676 fhEMCalCellCounts->Sumw2();
677
678 flCluster->Add(fhClusterPt);
679 flCluster->Add(fhClusterEta);
680 flCluster->Add(fhClusterPhi);
681 flCluster->Add(fhClusterEtaPhi);
682 flCluster->Add(fhClusterPhiPt);
683 flCluster->Add(fhClusterEtaPt);
62620fff 684 flCluster->Add(fhEMCalEventMult);
685 flCluster->Add(fpEMCalEventMult);
686 flCluster->Add(fpClusterPtProfile);
687 flCluster->Add(fhEMCalCellCounts);
688 fOutput->Add(flCluster);
689 }
c6202663 690
62620fff 691 if (fCalculateRhoJet>=0) // Default Rho & Raw Jet Spectra
692 {
693 fEMCalRawJets = new AlipAJetHistos("EMCalRawJets",fCentralityTag);
694
695 fRhoChargedCMSScale = new AlipAJetHistos("RhoChargedCMSScale",fCentralityTag,fDoNEF);
696 fRhoChargedCMSScale->SetSignalTrackPtBias(fSignalTrackBias);
d2af75be 697 fRhoChargedCMSScale->DoNEFSignalOnly(fDoNEFSignalOnly);
62620fff 698 fOutput->Add(fEMCalRawJets->GetOutputHistos());
699 fOutput->Add(fRhoChargedCMSScale->GetOutputHistos());
700
701 }
702 if (fCalculateRhoJet>=1) // Basic Rho & Raw Jet Spectra
703 {
704 fRhoChargedScale = new AlipAJetHistos("RhoChargedScale",fCentralityTag);
705 fRhoChargedScale->SetSignalTrackPtBias(fSignalTrackBias);
706
707 fOutput->Add(fRhoChargedScale->GetOutputHistos());
708 }
709 if (fCalculateRhoJet>=2) // Basic Rho & Raw Jet Spectra
710 {
711 fTPCRawJets = new AlipAJetHistos("TPCRawJets",fCentralityTag);
712 fRhoFull0 = new AlipAJetHistos("RhoFull0",fCentralityTag);
713 fRhoFull0->SetSignalTrackPtBias(fSignalTrackBias);
714 fRhoFull1 = new AlipAJetHistos("RhoFull1",fCentralityTag);
715 fRhoFull1->SetSignalTrackPtBias(fSignalTrackBias);
716 fRhoFull2 = new AlipAJetHistos("RhoFull2",fCentralityTag);
717 fRhoFull2->SetSignalTrackPtBias(fSignalTrackBias);
718 fRhoFullN = new AlipAJetHistos("RhoFullN",fCentralityTag);
719 fRhoFullN->SetSignalTrackPtBias(fSignalTrackBias);
720 fRhoFullDijet = new AlipAJetHistos("RhoFullDijet",fCentralityTag);
721 fRhoFullDijet->SetSignalTrackPtBias(fSignalTrackBias);
722 fRhoFullkT = new AlipAJetHistos("RhoFullkT",fCentralityTag);
723 fRhoFullkT->SetSignalTrackPtBias(fSignalTrackBias);
724 fRhoFullCMS = new AlipAJetHistos("RhoFullCMS",fCentralityTag);
725 fRhoFullCMS->SetSignalTrackPtBias(fSignalTrackBias);
726 fRhoCharged0 = new AlipAJetHistos("RhoCharged0",fCentralityTag);
727 fRhoCharged0->SetSignalTrackPtBias(fSignalTrackBias);
728 fRhoCharged1 = new AlipAJetHistos("RhoCharged1",fCentralityTag);
729 fRhoCharged1->SetSignalTrackPtBias(fSignalTrackBias);
730 fRhoCharged2 = new AlipAJetHistos("RhoCharged2",fCentralityTag);
731 fRhoCharged2->SetSignalTrackPtBias(fSignalTrackBias);
732 fRhoChargedN = new AlipAJetHistos("RhoChargedN",fCentralityTag);
733 fRhoChargedN->SetSignalTrackPtBias(fSignalTrackBias);
734 fRhoChargedkT = new AlipAJetHistos("RhoChargedkT",fCentralityTag);
735 fRhoChargedkT->SetSignalTrackPtBias(fSignalTrackBias);
736 fRhoChargedkTScale = new AlipAJetHistos("RhoChargedkTScale",fCentralityTag);
737 fRhoChargedkTScale->SetSignalTrackPtBias(fSignalTrackBias);
738 fRhoChargedCMS = new AlipAJetHistos("RhoChargedCMS",fCentralityTag);
739 fRhoChargedCMS->SetSignalTrackPtBias(fSignalTrackBias);
740
741 fOutput->Add(fTPCRawJets->GetOutputHistos());
742 fOutput->Add(fRhoFull0->GetOutputHistos());
743 fOutput->Add(fRhoFull1->GetOutputHistos());
744 fOutput->Add(fRhoFull2->GetOutputHistos());
745 fOutput->Add(fRhoFullN->GetOutputHistos());
746 fOutput->Add(fRhoFullDijet->GetOutputHistos());
747 fOutput->Add(fRhoFullkT->GetOutputHistos());
748 fOutput->Add(fRhoFullCMS->GetOutputHistos());
749 fOutput->Add(fRhoCharged0->GetOutputHistos());
750 fOutput->Add(fRhoCharged1->GetOutputHistos());
751 fOutput->Add(fRhoCharged2->GetOutputHistos());
752 fOutput->Add(fRhoChargedN->GetOutputHistos());
753 fOutput->Add(fRhoChargedkT->GetOutputHistos());
754 fOutput->Add(fRhoChargedkTScale->GetOutputHistos());
755 fOutput->Add(fRhoChargedCMS->GetOutputHistos());
756 }
3e43a01f 757
c54b626a 758 fOutput->Add(fhCentrality);
62620fff 759
ac6a3f1e 760 // Post data for ALL output slots >0 here,
761 // To get at least an empty histogram
762 // 1 is the outputnumber of a certain weg of task 1
763 PostData(1, fOutput);
764}
765
766void AliAnalysisTaskFullpAJets::UserExecOnce()
767{
d812e269 768 // Get the event tracks
769 fOrgTracks = dynamic_cast <TClonesArray*>(InputEvent()->FindListObject(fTrackName));
ac6a3f1e 770
d812e269 771 // Get the event caloclusters
772 fOrgClusters = dynamic_cast <TClonesArray*>(InputEvent()->FindListObject(fClusName));
ac6a3f1e 773
774 // Get charged jets
d812e269 775 fmyKTChargedJets = dynamic_cast <TClonesArray*>(InputEvent()->FindListObject(fkTChargedName));
776 fmyAKTChargedJets = dynamic_cast <TClonesArray*>(InputEvent()->FindListObject(fAkTChargedName));
78246241 777
ac6a3f1e 778 // Get the full jets
d812e269 779 fmyKTFullJets = dynamic_cast <TClonesArray*>(InputEvent()->FindListObject(fkTFullName));
780 fmyAKTFullJets = dynamic_cast <TClonesArray*>(InputEvent()->FindListObject(fAkTFullName));
781
ac6a3f1e 782 fIsInitialized=kTRUE;
783}
784//________________________________________________________________________
785void AliAnalysisTaskFullpAJets::UserExec(Option_t *)
786{
787 if (fIsInitialized==kFALSE)
788 {
789 UserExecOnce();
790 }
791
792 // Get pointer to reconstructed event
08b981da 793 fEvent = InputEvent();
794 if (!fEvent)
ac6a3f1e 795 {
796 AliError("Pointer == 0, this can not happen!");
797 return;
798 }
799
08b981da 800 AliESDEvent* esd = dynamic_cast<AliESDEvent*>(fEvent);
801 AliAODEvent* aod = dynamic_cast<AliAODEvent*>(fEvent);
802
803 fRecoUtil = new AliEMCALRecoUtils();
804 fEMCALGeometry = AliEMCALGeometry::GetInstance();
ac6a3f1e 805
806 if (esd)
807 {
c6202663 808 fEventCentrality=esd->GetCentrality()->GetCentralityPercentile(fCentralityTag);
ac6a3f1e 809
810 if (esd->GetPrimaryVertex()->GetNContributors()<1 || (TMath::Abs(esd->GetPrimaryVertex()->GetZ())>fVertexWindow))
811 {
812 return;
813 }
814 if (TMath::Sqrt(TMath::Power(esd->GetPrimaryVertex()->GetX(),2)+TMath::Power(esd->GetPrimaryVertex()->GetY(),2))>fVertexMaxR)
815 {
816 return;
817 }
818
8daeee93 819 esd->GetPrimaryVertex()->GetXYZ(fVertex);
08b981da 820 fCells = (AliVCaloCells*) esd->GetEMCALCells();
821 fnCaloClusters = esd->GetNumberOfCaloClusters();
ac6a3f1e 822 }
823 else if (aod)
824 {
c6202663 825 fEventCentrality=aod->GetCentrality()->GetCentralityPercentile(fCentralityTag);
ac6a3f1e 826
827 if (aod->GetPrimaryVertex()->GetNContributors()<1 || (TMath::Abs(aod->GetPrimaryVertex()->GetZ())>fVertexWindow))
828 {
829 return;
830 }
831 if (TMath::Sqrt(TMath::Power(aod->GetPrimaryVertex()->GetX(),2)+TMath::Power(aod->GetPrimaryVertex()->GetY(),2))>fVertexMaxR)
832 {
833 return;
834 }
835
8daeee93 836 aod->GetPrimaryVertex()->GetXYZ(fVertex);
08b981da 837 fCells = (AliVCaloCells*) aod->GetEMCALCells();
838 fnCaloClusters = aod->GetNumberOfCaloClusters();
ac6a3f1e 839 }
840 else
841 {
842 AliError("Cannot get AOD/ESD event. Rejecting Event");
843 return;
844 }
845
846 // Make sure Centrality isn't exactly 100% (to avoid bin filling errors in profile plots. Make it 99.99
847 if (fEventCentrality>99.99)
848 {
849 fEventCentrality=99.99;
850 }
c54b626a 851 fhCentrality->Fill(fEventCentrality);
ac6a3f1e 852
c6202663 853 TrackCuts();
ac6a3f1e 854 // Reject any event that doesn't have any tracks, i.e. TPC is off
855 if (fnTracks<1)
856 {
25322962
ML		 857 if (fTrackQA==kTRUE)
858 {
859 fhTPCEventMult->Fill(fEventCentrality,0.0);
860 fpTPCEventMult->Fill(fEventCentrality,0.0);
861 fhEMCalTrackEventMult->Fill(fEventCentrality,0.0);
862 }
ac6a3f1e 863 AliWarning("No PicoTracks, Rejecting Event");
864 return;
865 }
d812e269 866
c6202663 867 ClusterCuts();
ac6a3f1e 868 if (fnClusters<1)
869 {
870 AliInfo("No Corrected CaloClusters, using only charged jets");
871
62620fff 872 if (fTrackQA==kTRUE)
873 {
874 TrackHisto();
875 }
ac6a3f1e 876 InitChargedJets();
c6202663 877 GenerateTPCRandomConesPt();
cf074128 878
25322962
ML		 879 if (fClusterQA==kTRUE)
880 {
881 fhEMCalEventMult->Fill(fEventCentrality,0.0);
882 fpEMCalEventMult->Fill(fEventCentrality,0.0);
883 }
884
c6202663 885 // Rho's
62620fff 886 if (fCalculateRhoJet>=2)
887 {
888 EstimateChargedRho0();
889 EstimateChargedRho1();
890 EstimateChargedRho2();
891 EstimateChargedRhoN();
892 EstimateChargedRhokT();
893 EstimateChargedRhoCMS();
894 }
c6202663 895
c6202663 896 DeleteJetData(kFALSE);
ac6a3f1e 897
898 fnEventsCharged++;
c6202663 899
900 PostData(1, fOutput);
ac6a3f1e 901 return;
902 }
d812e269 903
62620fff 904 if (fTrackQA==kTRUE)
905 {
906 TrackHisto();
907 }
908 if (fClusterQA==kTRUE)
909 {
910 ClusterHisto();
911 }
ac6a3f1e 912
913 // Prep the jets
914 InitChargedJets();
915 InitFullJets();
c6202663 916 GenerateTPCRandomConesPt();
917 GenerateEMCalRandomConesPt();
918
62620fff 919 if (fCalculateRhoJet>=0)
ac6a3f1e 920 {
62620fff 921 EstimateChargedRhoCMSScale();
7acc3e04 922 }
62620fff 923 if (fCalculateRhoJet>=1)
924 {
925 EstimateChargedRhoScale();
926 }
927 if (fCalculateRhoJet>=2)
928 {
929 EstimateChargedRho0();
930 EstimateChargedRho1();
931 EstimateChargedRho2();
932 EstimateChargedRhoN();
933 EstimateChargedRhokT();
934 EstimateChargedRhoCMS();
935
936 EstimateFullRho0();
937 EstimateFullRho1();
938 EstimateFullRho2();
939 EstimateFullRhoN();
940 EstimateFullRhokT();
941 EstimateFullRhoCMS();
942
943 EstimateChargedRhokTScale();
944
945 // Dijet
946 if (IsDiJetEvent()==kTRUE)
947 {
948 EstimateFullRhoDijet();
949 }
950 }
951
ac6a3f1e 952 // Delete Dynamic Arrays
c6202663 953 DeleteJetData(kTRUE);
08b981da 954 delete fRecoUtil;
ac6a3f1e 955 fnEvents++;
956
957 PostData(1, fOutput);
958}
959
960//________________________________________________________________________
961void AliAnalysisTaskFullpAJets::Terminate(Option_t *) //specify what you want to have done
962{
963 // Called once at the end of the query. Done nothing here.
ac6a3f1e 964}
965
966/////////////////////////////////////////////////////////////////////////////////////////
967///////////////// User Defined Sub_Routines ///////////////////////////////////////
968/////////////////////////////////////////////////////////////////////////////////////////
969
c6202663 970void AliAnalysisTaskFullpAJets::TrackCuts()
971{
972 // Fill a TObjArray with the tracks from a TClonesArray which grabs the picotracks.
973 Int_t i;
974
975 fmyTracks = new TObjArray();
976 for (i=0;i<fOrgTracks->GetEntries();i++)
977 {
978 AliVTrack* vtrack = (AliVTrack*) fOrgTracks->At(i);
979 if (vtrack->Pt()>=fTrackMinPt)
980 {
981 fmyTracks->Add(vtrack);
982 }
983 }
984 fnTracks = fmyTracks->GetEntries();
985}
986
987void AliAnalysisTaskFullpAJets::ClusterCuts()
988{
989 // Fill a TObjArray with the clusters from a TClonesArray which grabs the caloclusterscorr.
990 Int_t i;
991
992 fmyClusters = new TObjArray();
62620fff 993 TLorentzVector *cluster_vec = new TLorentzVector;
994
20f2d13a 995 if(fOrgClusters)
c6202663 996 {
20f2d13a 997 for (i=0;i<fOrgClusters->GetEntries();i++)
c6202663 998 {
20f2d13a 999 AliVCluster* vcluster = (AliVCluster*) fOrgClusters->At(i);
8daeee93 1000 vcluster->GetMomentum(*cluster_vec,fVertex);
20f2d13a 1001
36cb7ae2 1002 if (cluster_vec->Pt()>=fClusterMinPt && vcluster->IsEMCAL()==kTRUE)
20f2d13a 1003 {
1004 fmyClusters->Add(vcluster);
1005 }
c6202663 1006 }
c6202663 1007 }
1008 fnClusters = fmyClusters->GetEntries();
62620fff 1009 delete cluster_vec;
c6202663 1010}
1011
ac6a3f1e 1012void AliAnalysisTaskFullpAJets::TrackHisto()
1013{
1014 // Fill track histograms: Phi,Eta,Pt
1015 Int_t i,j;
1016 Int_t TCBins=100;
1cd88f63 1017 TH2F *hdummypT= new TH2F("hdummypT","",TCBins,fTPCEtaMin,fTPCEtaMax,TCBins,fTPCPhiMin,fTPCPhiMax); //!
ac6a3f1e 1018
1019 for (i=0;i<fnTracks;i++)
1020 {
62620fff 1021 AliPicoTrack* vtrack = (AliPicoTrack*) fmyTracks->At(i);
ac6a3f1e 1022 fhTrackPt->Fill(vtrack->Pt());
1023 fhTrackEta->Fill(vtrack->Eta());
1024 fhTrackPhi->Fill(vtrack->Phi());
1025 fhTrackEtaPhi->Fill(vtrack->Eta(),vtrack->Phi());
24a61909 1026 fhTrackPhiPt->Fill(vtrack->Phi(),vtrack->Pt());
1027 fhTrackEtaPt->Fill(vtrack->Eta(),vtrack->Pt());
24a61909 1028
1029 // Fill Associated Track Distributions for AOD QA Productions
1030 // Global Tracks
1031 if (vtrack->GetTrackType()==0)
1032 {
1033 fhGlobalTrackPt->Fill(vtrack->Pt());
1034 fhGlobalTrackEta->Fill(vtrack->Eta());
1035 fhGlobalTrackPhi->Fill(vtrack->Phi());
1036 fhGlobalTrackEtaPhi->Fill(vtrack->Eta(),vtrack->Phi());
1037 fhGlobalTrackPhiPt->Fill(vtrack->Phi(),vtrack->Pt());
1038 fhGlobalTrackEtaPt->Fill(vtrack->Eta(),vtrack->Pt());
24a61909 1039 }
1040 // Complementary Tracks
1041 else if (vtrack->GetTrackType()==1)
1042 {
1043 fhComplementaryTrackPt->Fill(vtrack->Pt());
1044 fhComplementaryTrackEta->Fill(vtrack->Eta());
1045 fhComplementaryTrackPhi->Fill(vtrack->Phi());
1046 fhComplementaryTrackEtaPhi->Fill(vtrack->Eta(),vtrack->Phi());
1047 fhComplementaryTrackPhiPt->Fill(vtrack->Phi(),vtrack->Pt());
1048 fhComplementaryTrackEtaPt->Fill(vtrack->Eta(),vtrack->Pt());
24a61909 1049 }
ac6a3f1e 1050 hdummypT->Fill(vtrack->Eta(),vtrack->Phi(),vtrack->Pt());
1051 }
1052 for (i=1;i<=TCBins;i++)
1053 {
1054 for (j=1;j<=TCBins;j++)
1055 {
1056 fpTrackPtProfile->Fill(hdummypT->GetXaxis()->GetBinCenter(i),hdummypT->GetYaxis()->GetBinCenter(j),fTPCArea*TMath::Power(TCBins,-2)*hdummypT->GetBinContent(i,j));
1057 }
1058 }
1059 delete hdummypT;
1060}
1061
1062void AliAnalysisTaskFullpAJets::ClusterHisto()
1063{
1064 // Fill cluster histograms: Phi,Eta,Pt
1065 Int_t i,j;
1066 Int_t TCBins=100;
1cd88f63 1067 TH2F *hdummypT= new TH2F("hdummypT","",TCBins,fEMCalEtaMin,fEMCalEtaMax,TCBins,fEMCalPhiMin,fEMCalPhiMax); //!
ac6a3f1e 1068 Int_t myCellID=-2;
62620fff 1069 TLorentzVector *cluster_vec = new TLorentzVector;
ac6a3f1e 1070
1071 for (i=0;i<fnClusters;i++)
1072 {
1073 AliVCluster* vcluster = (AliVCluster*) fmyClusters->At(i);
8daeee93 1074 vcluster->GetMomentum(*cluster_vec,fVertex);
ac6a3f1e 1075
1076 fhClusterPt->Fill(cluster_vec->Pt());
1077 fhClusterEta->Fill(cluster_vec->Eta());
1078 fhClusterPhi->Fill(cluster_vec->Phi());
1079 fhClusterEtaPhi->Fill(cluster_vec->Eta(),cluster_vec->Phi());
24a61909 1080 fhClusterPhiPt->Fill(cluster_vec->Phi(),cluster_vec->Pt());
1081 fhClusterEtaPt->Fill(cluster_vec->Eta(),cluster_vec->Pt());
ac6a3f1e 1082 hdummypT->Fill(cluster_vec->Eta(),cluster_vec->Phi(),cluster_vec->Pt());
08b981da 1083 fEMCALGeometry->GetAbsCellIdFromEtaPhi(cluster_vec->Eta(),cluster_vec->Phi(),myCellID);
ac6a3f1e 1084 fhEMCalCellCounts->Fill(myCellID);
ac6a3f1e 1085 }
1086 for (i=1;i<=TCBins;i++)
1087 {
1088 for (j=1;j<=TCBins;j++)
1089 {
1090 fpClusterPtProfile->Fill(hdummypT->GetXaxis()->GetBinCenter(i),hdummypT->GetYaxis()->GetBinCenter(j),fEMCalArea*TMath::Power(TCBins,-2)*hdummypT->GetBinContent(i,j));
1091 }
1092 }
ac6a3f1e 1093 delete hdummypT;
62620fff 1094 delete cluster_vec;
ac6a3f1e 1095}
1096
c6202663 1097void AliAnalysisTaskFullpAJets::InitChargedJets()
1098{
1099 // Preliminary Jet Placement and Selection Cuts
1100 Int_t i;
1101
1102 fnAKTChargedJets = fmyAKTChargedJets->GetEntries();
1103 fnKTChargedJets = fmyKTChargedJets->GetEntries();
1104
1105 fTPCJet = new AlipAJetData("fTPCJet",kFALSE,fnAKTChargedJets);
1106 fTPCFullJet = new AlipAJetData("fTPCFullJet",kFALSE,fnAKTChargedJets);
1107 fTPCOnlyJet = new AlipAJetData("fTPCOnlyJet",kFALSE,fnAKTChargedJets);
8daeee93 1108 fTPCJetUnbiased = new AlipAJetData("fTPCJetUnbiased",kFALSE,fnAKTChargedJets);
c6202663 1109
1110 fTPCJet->SetSignalCut(fTPCJetThreshold);
1111 fTPCJet->SetAreaCutFraction(fJetAreaCutFrac);
1112 fTPCJet->SetJetR(fJetR);
3e43a01f 1113 fTPCJet->SetSignalTrackPtBias(fSignalTrackBias);
c6202663 1114 fTPCFullJet->SetSignalCut(fTPCJetThreshold);
1115 fTPCFullJet->SetAreaCutFraction(fJetAreaCutFrac);
1116 fTPCFullJet->SetJetR(fJetR);
3e43a01f 1117 fTPCFullJet->SetSignalTrackPtBias(fSignalTrackBias);
c6202663 1118 fTPCOnlyJet->SetSignalCut(fTPCJetThreshold);
1119 fTPCOnlyJet->SetAreaCutFraction(fJetAreaCutFrac);
1120 fTPCOnlyJet->SetJetR(fJetR);
3e43a01f 1121 fTPCOnlyJet->SetSignalTrackPtBias(fSignalTrackBias);
8daeee93 1122 fTPCJetUnbiased->SetSignalCut(fTPCJetThreshold);
1123 fTPCJetUnbiased->SetAreaCutFraction(fJetAreaCutFrac);
1124 fTPCJetUnbiased->SetJetR(fJetR);
1125 fTPCJetUnbiased->SetSignalTrackPtBias(kFALSE);
c6202663 1126
1127 // Initialize Jet Data
1128 for (i=0;i<fnAKTChargedJets;i++)
1129 {
1130 AliEmcalJet *myJet =(AliEmcalJet*) fmyAKTChargedJets->At(i);
1131
1132 fTPCJet->SetIsJetInArray(IsInTPC(fJetR,myJet->Phi(),myJet->Eta(),kFALSE),i);
1133 fTPCFullJet->SetIsJetInArray(IsInTPC(fJetR,myJet->Phi(),myJet->Eta(),kTRUE),i);
1134 fTPCOnlyJet->SetIsJetInArray(IsInTPCFull(fJetR,myJet->Phi(),myJet->Eta()),i);
8daeee93 1135 fTPCJetUnbiased->SetIsJetInArray(IsInTPC(fJetR,myJet->Phi(),myJet->Eta(),kFALSE),i);
c6202663 1136 }
1137 fTPCJet->InitializeJetData(fmyAKTChargedJets,fnAKTChargedJets);
1138 fTPCFullJet->InitializeJetData(fmyAKTChargedJets,fnAKTChargedJets);
1139 fTPCOnlyJet->InitializeJetData(fmyAKTChargedJets,fnAKTChargedJets);
8daeee93 1140 fTPCJetUnbiased->InitializeJetData(fmyAKTChargedJets,fnAKTChargedJets);
c6202663 1141
1142 // kT Jets
1143 fTPCkTFullJet = new AlipAJetData("fTPCkTFullJet",kFALSE,fnKTChargedJets);
1144 fTPCkTFullJet->SetSignalCut(fTPCJetThreshold);
1145 fTPCkTFullJet->SetAreaCutFraction(0.25*fJetAreaCutFrac);
1146 fTPCkTFullJet->SetJetR(fJetR);
1147
1148 for (i=0;i<fnKTChargedJets;i++)
1149 {
1150 AliEmcalJet *myJet =(AliEmcalJet*) fmyKTChargedJets->At(i);
1151 fTPCkTFullJet->SetIsJetInArray(IsInTPC(fJetR,myJet->Phi(),myJet->Eta(),kTRUE),i);
1152 }
1153 fTPCkTFullJet->InitializeJetData(fmyKTChargedJets,fnKTChargedJets);
1154
1155 // Raw Charged Jet Spectra
62620fff 1156 if (fCalculateRhoJet>=2)
1157 {
1158 fTPCRawJets->FillBSJS(fEventCentrality,0.0,fTPCJetThreshold,fmyAKTChargedJets,fTPCFullJet->GetJets(),fTPCFullJet->GetTotalJets());
1159 }
c6202663 1160}
1161
1162void AliAnalysisTaskFullpAJets::InitFullJets()
1163{
1164 // Preliminary Jet Placement and Selection Cuts
1165 Int_t i;
1166
1167 fnAKTFullJets = fmyAKTFullJets->GetEntries();
1168 fnKTFullJets = fmyKTFullJets->GetEntries();
1169
1170 fEMCalJet = new AlipAJetData("fEMCalJet",kTRUE,fnAKTFullJets);
1171 fEMCalFullJet = new AlipAJetData("fEMCalFullJet",kTRUE,fnAKTFullJets);
1172 fEMCalPartJet = new AlipAJetData("fEMCalPartJet",kTRUE,fnAKTFullJets);
8daeee93 1173 fEMCalPartJetUnbiased = new AlipAJetData("fEMCalPartJetUnbiased",kTRUE,fnAKTFullJets);
c6202663 1174
1175 fEMCalJet->SetSignalCut(fEMCalJetThreshold);
1176 fEMCalJet->SetAreaCutFraction(fJetAreaCutFrac);
1177 fEMCalJet->SetJetR(fJetR);
91d0893e 1178 fEMCalJet->SetNEF(fNEFSignalJetCut);
3e43a01f 1179 fEMCalJet->SetSignalTrackPtBias(fSignalTrackBias);
c6202663 1180 fEMCalFullJet->SetSignalCut(fEMCalJetThreshold);
1181 fEMCalFullJet->SetAreaCutFraction(fJetAreaCutFrac);
1182 fEMCalFullJet->SetJetR(fJetR);
91d0893e 1183 fEMCalFullJet->SetNEF(fNEFSignalJetCut);
3e43a01f 1184 fEMCalFullJet->SetSignalTrackPtBias(fSignalTrackBias);
c6202663 1185 fEMCalPartJet->SetSignalCut(fEMCalJetThreshold);
1186 fEMCalPartJet->SetAreaCutFraction(fJetAreaCutFrac);
1187 fEMCalPartJet->SetJetR(fJetR);
91d0893e 1188 fEMCalPartJet->SetNEF(fNEFSignalJetCut);
3e43a01f 1189 fEMCalPartJet->SetSignalTrackPtBias(fSignalTrackBias);
8daeee93 1190 fEMCalPartJetUnbiased->SetSignalCut(fEMCalJetThreshold);
1191 fEMCalPartJetUnbiased->SetAreaCutFraction(fJetAreaCutFrac);
1192 fEMCalPartJetUnbiased->SetJetR(fJetR);
1193 fEMCalPartJetUnbiased->SetNEF(1.0);
1194 fEMCalPartJetUnbiased->SetSignalTrackPtBias(kFALSE);
c6202663 1195
1196 // Initialize Jet Data
1197 for (i=0;i<fnAKTFullJets;i++)
1198 {
1199 AliEmcalJet *myJet =(AliEmcalJet*) fmyAKTFullJets->At(i);
1200
1201 fEMCalJet->SetIsJetInArray(IsInEMCal(myJet->Phi(),myJet->Eta()),i);
1202 fEMCalFullJet->SetIsJetInArray(IsInEMCalFull(fJetR,myJet->Phi(),myJet->Eta()),i);
1203 fEMCalPartJet->SetIsJetInArray(IsInEMCalPart(fJetR,myJet->Phi(),myJet->Eta()),i);
8daeee93 1204 fEMCalPartJetUnbiased->SetIsJetInArray(IsInEMCalPart(fJetR,myJet->Phi(),myJet->Eta()),i);
c6202663 1205 }
1206 fEMCalJet->InitializeJetData(fmyAKTFullJets,fnAKTFullJets);
1207 fEMCalFullJet->InitializeJetData(fmyAKTFullJets,fnAKTFullJets);
1208 fEMCalPartJet->InitializeJetData(fmyAKTFullJets,fnAKTFullJets);
8daeee93 1209 fEMCalPartJetUnbiased->InitializeJetData(fmyAKTFullJets,fnAKTFullJets);
c6202663 1210
1211 // kT Jets
1212 fEMCalkTFullJet = new AlipAJetData("fEMCalkTFullJet",kTRUE,fnKTFullJets);
1213 fEMCalkTFullJet->SetSignalCut(fEMCalJetThreshold);
1214 fEMCalkTFullJet->SetAreaCutFraction(0.25*fJetAreaCutFrac);
1215 fEMCalkTFullJet->SetJetR(fJetR);
91d0893e 1216 fEMCalkTFullJet->SetNEF(fNEFSignalJetCut);
3e43a01f 1217 fEMCalkTFullJet->SetSignalTrackPtBias(fSignalTrackBias);
c6202663 1218
1219 for (i=0;i<fnKTFullJets;i++)
1220 {
1221 AliEmcalJet *myJet =(AliEmcalJet*) fmyKTFullJets->At(i);
1222 fEMCalkTFullJet->SetIsJetInArray(IsInEMCalFull(fJetR,myJet->Phi(),myJet->Eta()),i);
1223 }
1224 fEMCalkTFullJet->InitializeJetData(fmyKTFullJets,fnKTFullJets);
1225
1226 // Raw Full Jet Spectra
1227 fEMCalRawJets->FillBSJS(fEventCentrality,0.0,fEMCalJetThreshold,fmyAKTFullJets,fEMCalFullJet->GetJets(),fEMCalFullJet->GetTotalJets());
1228}
1229
c6202663 1230void AliAnalysisTaskFullpAJets::GenerateTPCRandomConesPt()
ac6a3f1e 1231{
1232 Int_t i,j;
1233 Double_t E_tracks_total=0.;
ac6a3f1e 1234 TRandom3 u(time(NULL));
1235
c6202663 1236 Double_t Eta_Center=0.5*(fTPCEtaMin+fTPCEtaMax);
1237 Double_t Phi_Center=0.5*(fTPCPhiMin+fTPCPhiMax);
ac6a3f1e 1238 Int_t event_mult=0;
3e43a01f 1239 Int_t event_track_mult=0;
ac6a3f1e 1240 Int_t clus_mult=0;
1241
c6202663 1242 for (i=0;i<fnBckgClusters;i++)
1243 {
1244 fTPCRCBckgFluc[i]=0.0;
1245 fTPCRCBckgFlucSignal[i]=0.0;
7acc3e04 1246 fTPCRCBckgFlucNColl[i]=0.0;
c6202663 1247 }
1248
ac6a3f1e 1249 TLorentzVector *dummy= new TLorentzVector;
c6202663 1250 TLorentzVector *temp_jet= new TLorentzVector;
62620fff 1251 TLorentzVector *track_vec = new TLorentzVector;
1252
c6202663 1253 // First, consider the RC with no spatial restrictions
ac6a3f1e 1254 for (j=0;j<fnBckgClusters;j++)
1255 {
ac6a3f1e 1256 E_tracks_total=0.;
ac6a3f1e 1257
c6202663 1258 dummy->SetPtEtaPhiE(1,u.Uniform(Eta_Center-fJetR,Eta_Center+fJetR),u.Uniform(Phi_Center-fJetR,Phi_Center+fJetR),0);
ac6a3f1e 1259 // Loop over all tracks
1260 for (i=0;i<fnTracks;i++)
1261 {
1262 AliVTrack* vtrack =(AliVTrack*) fmyTracks->At(i);
c6202663 1263 if (IsInTPC(fJetR,vtrack->Phi(),vtrack->Eta(),kFALSE)==kTRUE)
ac6a3f1e 1264 {
ac6a3f1e 1265 track_vec->SetPtEtaPhiE(vtrack->Pt(),vtrack->Eta(),vtrack->Phi(),vtrack->E());
1266 if (dummy->DeltaR(*track_vec)<fJetR)
1267 {
ac6a3f1e 1268 E_tracks_total+=vtrack->Pt();
1269 }
ac6a3f1e 1270 }
1271 }
c6202663 1272 fTPCRCBckgFlucSignal[j]=E_tracks_total;
1273 }
1274
1275 // Now, consider the RC where the vertex of RC is at least 2R away from the leading signal
1276 E_tracks_total=0.0;
8daeee93 1277 if (fTPCJetUnbiased->GetLeadingPt()<0.0)
c6202663 1278 {
1279 temp_jet->SetPtEtaPhiE(1,u.Uniform(Eta_Center-fJetR,Eta_Center+fJetR),u.Uniform(Phi_Center-fJetR,Phi_Center+fJetR),0);
1280 }
1281 else
1282 {
8daeee93 1283 AliEmcalJet *myJet = (AliEmcalJet*) fmyAKTChargedJets->At(fTPCJetUnbiased->GetLeadingIndex());
c6202663 1284 myJet->GetMom(*temp_jet);
1285 }
1286
1287 for (j=0;j<fnBckgClusters;j++)
1288 {
1289 event_mult=0;
3e43a01f 1290 event_track_mult=0;
c6202663 1291 clus_mult=0;
1292 E_tracks_total=0.;
ac6a3f1e 1293
c6202663 1294 dummy->SetPtEtaPhiE(1,u.Uniform(Eta_Center-fJetR,Eta_Center+fJetR),u.Uniform(Phi_Center-fJetR,Phi_Center+fJetR),0);
1295 while (temp_jet->DeltaR(*dummy)<fJetR)
ac6a3f1e 1296 {
c6202663 1297 dummy->SetPtEtaPhiE(1,u.Uniform(Eta_Center-fJetR,Eta_Center+fJetR),u.Uniform(Phi_Center-fJetR,Phi_Center+fJetR),0);
ac6a3f1e 1298 }
c6202663 1299 // Loop over all tracks
1300 for (i=0;i<fnTracks;i++)
ac6a3f1e 1301 {
c6202663 1302 AliVTrack* vtrack =(AliVTrack*) fmyTracks->At(i);
1303 if (IsInTPC(fJetR,vtrack->Phi(),vtrack->Eta(),kFALSE)==kTRUE)
1304 {
1305 event_mult++;
3e43a01f 1306 if (IsInEMCal(vtrack->Phi(),vtrack->Eta()==kTRUE))
1307 {
1308 event_track_mult++;
1309 }
c6202663 1310 track_vec->SetPtEtaPhiE(vtrack->Pt(),vtrack->Eta(),vtrack->Phi(),vtrack->E());
1311 if (dummy->DeltaR(*track_vec)<fJetR)
1312 {
1313 clus_mult++;
1314 E_tracks_total+=vtrack->Pt();
1315 }
c6202663 1316 }
ac6a3f1e 1317 }
c6202663 1318 fTPCRCBckgFluc[j]=E_tracks_total;
ac6a3f1e 1319 }
62620fff 1320 if (fTrackQA==kTRUE)
1321 {
1322 fhTPCEventMult->Fill(fEventCentrality,event_mult);
1323 fpTPCEventMult->Fill(fEventCentrality,event_mult);
1324 fhEMCalTrackEventMult->Fill(fEventCentrality,event_track_mult);
1325 }
1326 if (fCalculateRhoJet>=2)
1327 {
1328 fTPCRawJets->FillDeltaPt(fEventCentrality,0.0,fJetR,fTPCRCBckgFluc,1);
1329 }
ac6a3f1e 1330
7acc3e04 1331 // For the case of partial exclusion, merely allow a superposition of full and no exclusion with probability p=1/Ncoll
1332 Double_t exclusion_prob;
1333 for (j=0;j<fnBckgClusters;j++)
1334 {
1335 exclusion_prob = u.Uniform(0,1);
1336 if (exclusion_prob<(1/fNColl))
1337 {
1338 fTPCRCBckgFlucNColl[j]=fTPCRCBckgFlucSignal[j];
1339 }
1340 else
1341 {
1342 fTPCRCBckgFlucNColl[j]=fTPCRCBckgFluc[j];
1343 }
1344 }
1345
ac6a3f1e 1346 delete dummy;
c6202663 1347 delete temp_jet;
62620fff 1348 delete track_vec;
ac6a3f1e 1349}
1350
c6202663 1351void AliAnalysisTaskFullpAJets::GenerateEMCalRandomConesPt()
ac6a3f1e 1352{
78246241 1353 Int_t i,j;
c6202663 1354 Double_t E_tracks_total=0.;
1355 Double_t E_caloclusters_total=0.;
1356 TRandom3 u(time(NULL));
ac6a3f1e 1357
c6202663 1358 Double_t Eta_Center=0.5*(fEMCalEtaMin+fEMCalEtaMax);
1359 Double_t Phi_Center=0.5*(fEMCalPhiMin+fEMCalPhiMax);
1360 Int_t event_mult=0;
1361 Int_t clus_mult=0;
ac6a3f1e 1362
c6202663 1363 for (i=0;i<fnBckgClusters;i++)
ac6a3f1e 1364 {
c6202663 1365 fEMCalRCBckgFluc[i]=0.0;
1366 fEMCalRCBckgFlucSignal[i]=0.0;
7acc3e04 1367 fEMCalRCBckgFlucNColl[i]=0.0;
ac6a3f1e 1368 }
78246241 1369
c6202663 1370 TLorentzVector *dummy= new TLorentzVector;
1371 TLorentzVector *temp_jet= new TLorentzVector;
62620fff 1372 TLorentzVector *track_vec = new TLorentzVector;
1373 TLorentzVector *cluster_vec = new TLorentzVector;
1374
c6202663 1375 // First, consider the RC with no spatial restrictions
1376 for (j=0;j<fnBckgClusters;j++)
78246241 1377 {
c6202663 1378 E_tracks_total=0.;
1379 E_caloclusters_total=0.;
1380
1381 dummy->SetPtEtaPhiE(1,u.Uniform(Eta_Center-fJetR,Eta_Center+fJetR),u.Uniform(Phi_Center-fJetR,Phi_Center+fJetR),0);
1382 // Loop over all tracks
1383 for (i=0;i<fnTracks;i++)
78246241 1384 {
c6202663 1385 AliVTrack* vtrack =(AliVTrack*) fmyTracks->At(i);
1386 if (IsInEMCal(vtrack->Phi(),vtrack->Eta())==kTRUE)
78246241 1387 {
c6202663 1388 track_vec->SetPtEtaPhiE(vtrack->Pt(),vtrack->Eta(),vtrack->Phi(),vtrack->E());
1389 if (dummy->DeltaR(*track_vec)<fJetR)
78246241 1390 {
c6202663 1391 E_tracks_total+=vtrack->Pt();
78246241 1392 }
78246241 1393 }
1394 }
78246241 1395
c6202663 1396 // Loop over all caloclusters
1397 for (i=0;i<fnClusters;i++)
78246241 1398 {
c6202663 1399 AliVCluster* vcluster = (AliVCluster*) fmyClusters->At(i);
8daeee93 1400 vcluster->GetMomentum(*cluster_vec,fVertex);
c6202663 1401 if (dummy->DeltaR(*cluster_vec)<fJetR)
1402 {
1403 clus_mult++;
1404 E_caloclusters_total+=vcluster->E();
1405 }
78246241 1406 }
c6202663 1407 fEMCalRCBckgFlucSignal[j]=E_tracks_total+E_caloclusters_total;
78246241 1408 }
c6202663 1409
1410 // Now, consider the RC where the vertex of RC is at least 2R away from the leading signal
1411 E_tracks_total=0.;
1412 E_caloclusters_total=0.;
8daeee93 1413 if (fEMCalPartJetUnbiased->GetLeadingPt()<0.0)
78246241 1414 {
c6202663 1415 temp_jet->SetPtEtaPhiE(1,u.Uniform(Eta_Center-fJetR,Eta_Center+fJetR),u.Uniform(Phi_Center-fJetR,Phi_Center+fJetR),0);
78246241 1416 }
c6202663 1417 else
78246241 1418 {
8daeee93 1419 AliEmcalJet *myJet = (AliEmcalJet*) fmyAKTFullJets->At(fEMCalPartJetUnbiased->GetLeadingIndex());
c6202663 1420 myJet->GetMom(*temp_jet);
78246241 1421 }
1422
c6202663 1423 for (j=0;j<fnBckgClusters;j++)
78246241 1424 {
c6202663 1425 event_mult=0;
1426 clus_mult=0;
1427 E_tracks_total=0.;
1428 E_caloclusters_total=0.;
1429
1430 dummy->SetPtEtaPhiE(1,u.Uniform(Eta_Center-fJetR,Eta_Center+fJetR),u.Uniform(Phi_Center-fJetR,Phi_Center+fJetR),0);
1431 while (temp_jet->DeltaR(*dummy)<fJetR)
1432 {
1433 dummy->SetPtEtaPhiE(1,u.Uniform(Eta_Center-fJetR,Eta_Center+fJetR),u.Uniform(Phi_Center-fJetR,Phi_Center+fJetR),0);
1434 }
1435 // Loop over all tracks
1436 for (i=0;i<fnTracks;i++)
1437 {
1438 AliVTrack* vtrack =(AliVTrack*) fmyTracks->At(i);
1439 if (IsInEMCal(vtrack->Phi(),vtrack->Eta())==kTRUE)
1440 {
1441 event_mult++;
c6202663 1442 track_vec->SetPtEtaPhiE(vtrack->Pt(),vtrack->Eta(),vtrack->Phi(),vtrack->E());
1443 if (dummy->DeltaR(*track_vec)<fJetR)
1444 {
1445 clus_mult++;
1446 E_tracks_total+=vtrack->Pt();
1447 }
c6202663 1448 }
1449 }
1450
1451 // Loop over all caloclusters
1452 for (i=0;i<fnClusters;i++)
1453 {
1454 AliVCluster* vcluster = (AliVCluster*) fmyClusters->At(i);
8daeee93 1455 vcluster->GetMomentum(*cluster_vec,fVertex);
c6202663 1456 event_mult++;
1457 if (dummy->DeltaR(*cluster_vec)<fJetR)
1458 {
1459 clus_mult++;
1460 E_caloclusters_total+=vcluster->E();
1461 }
c6202663 1462 }
1463 fEMCalRCBckgFluc[j]=E_tracks_total+E_caloclusters_total;
78246241 1464 }
62620fff 1465 if (fClusterQA==kTRUE)
1466 {
1467 fhEMCalEventMult->Fill(fEventCentrality,event_mult);
1468 fpEMCalEventMult->Fill(fEventCentrality,event_mult);
1469 }
c6202663 1470 fEMCalRawJets->FillDeltaPt(fEventCentrality,0.0,fJetR,fEMCalRCBckgFluc,1);
1471
7acc3e04 1472 // For the case of partial exclusion, merely allow a superposition of full and no exclusion with probability p=1/Ncoll
1473 Double_t exclusion_prob;
1474 for (j=0;j<fnBckgClusters;j++)
1475 {
1476 exclusion_prob = u.Uniform(0,1);
1477 if (exclusion_prob<(1/fNColl))
1478 {
1479 fEMCalRCBckgFlucNColl[j]=fEMCalRCBckgFlucSignal[j];
1480 }
1481 else
1482 {
1483 fEMCalRCBckgFlucNColl[j]=fEMCalRCBckgFluc[j];
1484 }
1485 }
1486
c6202663 1487 delete dummy;
1488 delete temp_jet;
62620fff 1489 delete track_vec;
1490 delete cluster_vec;
ac6a3f1e 1491}
1492
c6202663 1493// Charged Rho's
1494void AliAnalysisTaskFullpAJets::EstimateChargedRho0()
ac6a3f1e 1495{
ac6a3f1e 1496 Int_t i;
c6202663 1497 Double_t E_tracks_total=0.0;
1498 Double_t TPC_rho=0.;
ac6a3f1e 1499
c6202663 1500 // Loop over all tracks
1501 for (i=0;i<fnTracks;i++)
1502 {
1503 AliVTrack* vtrack =(AliVTrack*) fmyTracks->At(i);
1504 if (IsInTPC(fJetR,vtrack->Phi(),vtrack->Eta(),kFALSE)==kTRUE)
1505 {
1506 E_tracks_total+=vtrack->Pt();
1507 }
1508 }
ac6a3f1e 1509
c6202663 1510 // Calculate the mean Background density
1511 TPC_rho=E_tracks_total/fTPCArea;
1512 fRhoCharged=TPC_rho;
ac6a3f1e 1513
c6202663 1514 // Fill Histograms
1515 fRhoCharged0->FillRho(fEventCentrality,TPC_rho);
1516 fRhoCharged0->FillBSJS(fEventCentrality,TPC_rho,fTPCJetThreshold,fmyAKTChargedJets,fTPCJet->GetJets(),fTPCJet->GetTotalJets());
1517 fRhoCharged0->FillDeltaPt(fEventCentrality,TPC_rho,fJetR,fTPCRCBckgFluc,1);
1518 fRhoCharged0->FillDeltaPtSignal(fEventCentrality,TPC_rho,fJetR,fTPCRCBckgFlucSignal,1);
7acc3e04 1519 fRhoCharged0->FillDeltaPtNColl(fEventCentrality,TPC_rho,fJetR,fTPCRCBckgFlucNColl,1);
c6202663 1520 fRhoCharged0->FillBackgroundFluctuations(fEventCentrality,TPC_rho,fJetR);
1521 fRhoCharged0->FillLeadingJetPtRho(fTPCJet->GetLeadingPt(),TPC_rho);
ac6a3f1e 1522
c6202663 1523}
ac6a3f1e 1524
c6202663 1525void AliAnalysisTaskFullpAJets::EstimateChargedRho1()
1526{
1527 Int_t i;
1528 Double_t E_tracks_total=0.0;
1529 Double_t TPC_rho=0.;
1530
62620fff 1531 TLorentzVector *temp_jet= new TLorentzVector;
1532 TLorentzVector *track_vec = new TLorentzVector;
1533
8daeee93 1534 if (fTPCJetUnbiased->GetLeadingPt()>=fTPCJetThreshold)
ac6a3f1e 1535 {
c6202663 1536 AliEmcalJet *myJet =(AliEmcalJet*) fmyAKTChargedJets->At(fTPCJet->GetLeadingIndex());
c6202663 1537 myJet->GetMom(*temp_jet);
ac6a3f1e 1538
c6202663 1539 // Loop over all tracks
1540 for (i=0;i<fnTracks;i++)
1541 {
1542 AliVTrack* vtrack =(AliVTrack*) fmyTracks->At(i);
1543 if (IsInTPC(fJetR,vtrack->Phi(),vtrack->Eta(),kFALSE)==kTRUE)
1544 {
c6202663 1545 track_vec->SetPtEtaPhiE(vtrack->Pt(),vtrack->Eta(),vtrack->Phi(),vtrack->E());
1546 if (temp_jet->DeltaR(*track_vec)>fJetRForRho)
1547 {
1548 E_tracks_total+=vtrack->Pt();
1549 }
c6202663 1550 }
1551 }
ac6a3f1e 1552
c6202663 1553 // Calculate the mean Background density
1554 TPC_rho=E_tracks_total/(fTPCArea-AreaWithinTPC(fJetR,myJet->Eta()));
1555 }
1556 else // i.e. No signal jets -> same as total background density
1557 {
1558 // Loop over all tracks
1559 for (i=0;i<fnTracks;i++)
ac6a3f1e 1560 {
c6202663 1561 AliVTrack* vtrack =(AliVTrack*) fmyTracks->At(i);
1562 if (IsInTPC(fJetR,vtrack->Phi(),vtrack->Eta(),kFALSE)==kTRUE)
ac6a3f1e 1563 {
c6202663 1564 E_tracks_total+=vtrack->Pt();
ac6a3f1e 1565 }
c6202663 1566 }
1567 // Calculate the mean Background density
1568 TPC_rho=E_tracks_total/fTPCArea;
1569 }
62620fff 1570 delete track_vec;
1571 delete temp_jet;
1572
c6202663 1573 // Fill histograms
1574 fRhoCharged1->FillRho(fEventCentrality,TPC_rho);
1575 fRhoCharged1->FillBSJS(fEventCentrality,TPC_rho,fTPCJetThreshold,fmyAKTChargedJets,fTPCFullJet->GetJets(),fTPCFullJet->GetTotalJets());
1576 fRhoCharged1->FillDeltaPt(fEventCentrality,TPC_rho,fJetR,fTPCRCBckgFluc,1);
1577 fRhoCharged1->FillDeltaPtSignal(fEventCentrality,TPC_rho,fJetR,fTPCRCBckgFlucSignal,1);
7acc3e04 1578 fRhoCharged1->FillDeltaPtNColl(fEventCentrality,TPC_rho,fJetR,fTPCRCBckgFlucNColl,1);
c6202663 1579 fRhoCharged1->FillBackgroundFluctuations(fEventCentrality,TPC_rho,fJetR);
1580 fRhoCharged1->FillLeadingJetPtRho(fTPCFullJet->GetLeadingPt(),TPC_rho);
1581}
1582
1583void AliAnalysisTaskFullpAJets::EstimateChargedRho2()
1584{
1585 Int_t i;
1586 Double_t E_tracks_total=0.0;
1587 Double_t TPC_rho=0.;
62620fff 1588
1589 TLorentzVector *temp_jet1= new TLorentzVector;
1590 TLorentzVector *temp_jet2= new TLorentzVector;
1591 TLorentzVector *track_vec = new TLorentzVector;
1592
8daeee93 1593 if ((fTPCJetUnbiased->GetLeadingPt()>=fTPCJetThreshold) && (fTPCJetUnbiased->GetSubLeadingPt()>=fTPCJetThreshold))
c6202663 1594 {
1595 AliEmcalJet *myhJet =(AliEmcalJet*) fmyAKTChargedJets->At(fTPCJet->GetLeadingIndex());
c6202663 1596 myhJet->GetMom(*temp_jet1);
1597
1598 AliEmcalJet *myJet =(AliEmcalJet*) fmyAKTChargedJets->At(fTPCJet->GetSubLeadingIndex());
c6202663 1599 myJet->GetMom(*temp_jet2);
1600
1601 // Loop over all tracks
1602 for (i=0;i<fnTracks;i++)
1603 {
1604 AliVTrack* vtrack =(AliVTrack*) fmyTracks->At(i);
1605 if (IsInTPC(fJetR,vtrack->Phi(),vtrack->Eta(),kFALSE)==kTRUE)
ac6a3f1e 1606 {
c6202663 1607 track_vec->SetPtEtaPhiE(vtrack->Pt(),vtrack->Eta(),vtrack->Phi(),vtrack->E());
1608 if ((temp_jet1->DeltaR(*track_vec)>fJetRForRho) && (temp_jet2->DeltaR(*track_vec)>fJetRForRho))
c54b626a 1609 {
c6202663 1610 E_tracks_total+=vtrack->Pt();
c54b626a 1611 }
c6202663 1612 }
1613 }
c6202663 1614
1615 // Calculate the mean Background density
1616 TPC_rho=E_tracks_total/(fTPCArea-AreaWithinTPC(fJetR,myhJet->Eta())-AreaWithinTPC(fJetR,myJet->Eta()));
1617 }
8daeee93 1618 else if (fTPCJetUnbiased->GetLeadingPt()>=fTPCJetThreshold)
c6202663 1619 {
1620 AliEmcalJet *myJet =(AliEmcalJet*) fmyAKTChargedJets->At(fTPCJet->GetLeadingIndex());
62620fff 1621 myJet->GetMom(*temp_jet1);
c6202663 1622
1623 // Loop over all tracks
1624 for (i=0;i<fnTracks;i++)
1625 {
1626 AliVTrack* vtrack =(AliVTrack*) fmyTracks->At(i);
1627 if (IsInTPC(fJetR,vtrack->Phi(),vtrack->Eta(),kFALSE)==kTRUE)
1628 {
c6202663 1629 track_vec->SetPtEtaPhiE(vtrack->Pt(),vtrack->Eta(),vtrack->Phi(),vtrack->E());
62620fff 1630 if (temp_jet1->DeltaR(*track_vec)>fJetRForRho)
ac6a3f1e 1631 {
c6202663 1632 E_tracks_total+=vtrack->Pt();
ac6a3f1e 1633 }
c6202663 1634 }
1635 }
c6202663 1636
1637 // Calculate the mean Background density
1638 TPC_rho=E_tracks_total/(fTPCArea-AreaWithinTPC(fJetR,myJet->Eta()));
1639 }
1640 else // i.e. No signal jets -> same as total background density
1641 {
1642 // Loop over all tracks
1643 for (i=0;i<fnTracks;i++)
1644 {
1645 AliVTrack* vtrack =(AliVTrack*) fmyTracks->At(i);
1646 if (IsInTPC(fJetR,vtrack->Phi(),vtrack->Eta(),kFALSE)==kTRUE)
1647 {
1648 E_tracks_total+=vtrack->Pt();
1649 }
1650 }
1651
1652 // Calculate the mean Background density
1653 TPC_rho=E_tracks_total/fTPCArea;
1654 }
62620fff 1655 delete temp_jet1;
1656 delete temp_jet2;
1657 delete track_vec;
1658
c6202663 1659 // Fill histograms
1660 fRhoCharged2->FillRho(fEventCentrality,TPC_rho);
1661 fRhoCharged2->FillBSJS(fEventCentrality,TPC_rho,fTPCJetThreshold,fmyAKTChargedJets,fTPCFullJet->GetJets(),fTPCFullJet->GetTotalJets());
1662 fRhoCharged2->FillDeltaPt(fEventCentrality,TPC_rho,fJetR,fTPCRCBckgFluc,1);
1663 fRhoCharged2->FillDeltaPtSignal(fEventCentrality,TPC_rho,fJetR,fTPCRCBckgFlucSignal,1);
7acc3e04 1664 fRhoCharged2->FillDeltaPtNColl(fEventCentrality,TPC_rho,fJetR,fTPCRCBckgFlucNColl,1);
c6202663 1665 fRhoCharged2->FillBackgroundFluctuations(fEventCentrality,TPC_rho,fJetR);
1666 fRhoCharged2->FillLeadingJetPtRho(fTPCFullJet->GetLeadingPt(),TPC_rho);
1667}
e864d416 1668
c6202663 1669void AliAnalysisTaskFullpAJets::EstimateChargedRhoN()
1670{
1671 Int_t i,j;
1672 Bool_t track_away_from_jet;
1673 Double_t E_tracks_total=0.0;
1674 Double_t TPC_rho=0.0;
1675 Double_t jet_area_total=0.0;
1676
62620fff 1677 TLorentzVector *jet_vec= new TLorentzVector;
1678 TLorentzVector *track_vec = new TLorentzVector;
1679
c6202663 1680 // First, sum all tracks within the EMCal that are away from jet(s) above Pt Threshold
1681 for (i=0;i<fnTracks;i++)
1682 {
1683 // First, check if track is in the EMCal!!
1684 AliVTrack* vtrack = (AliVTrack*) fmyTracks->At(i);
1685 if (IsInTPC(fJetR,vtrack->Phi(),vtrack->Eta(),kFALSE)==kTRUE)
1686 {
8daeee93 1687 if (fTPCJetUnbiased->GetTotalSignalJets()<1)
c6202663 1688 {
1689 E_tracks_total+=vtrack->Pt();
1690 }
1691 else
1692 {
1693 track_away_from_jet=kTRUE;
1694 j=0;
c6202663 1695 track_vec->SetPtEtaPhiE(vtrack->Pt(),vtrack->Eta(),vtrack->Phi(),vtrack->E());
8daeee93 1696 while (track_away_from_jet==kTRUE && j<fTPCJetUnbiased->GetTotalSignalJets())
ac6a3f1e 1697 {
8daeee93 1698 AliEmcalJet *myJet =(AliEmcalJet*) fmyAKTChargedJets->At(fTPCJetUnbiased->GetSignalJetIndex(j));
c6202663 1699 myJet->GetMom(*jet_vec);
1700 if (track_vec->DeltaR(*jet_vec)<=fJetRForRho)
1701 {
1702 track_away_from_jet=kFALSE;
1703 }
c6202663 1704 j++;
ac6a3f1e 1705 }
c6202663 1706 if (track_away_from_jet==kTRUE)
ac6a3f1e 1707 {
c6202663 1708 E_tracks_total+=vtrack->Pt();
ac6a3f1e 1709 }
1710 }
1711 }
c6202663 1712 }
1713
1714 // Determine area of all Jets that are within the EMCal
8daeee93 1715 if (fTPCJetUnbiased->GetTotalSignalJets()==0)
c6202663 1716 {
1717 jet_area_total=0.0;
1718 }
1719 else
1720 {
8daeee93 1721 for (i=0;i<fTPCJetUnbiased->GetTotalSignalJets();i++)
c6202663 1722 {
8daeee93 1723 AliEmcalJet *myJet = (AliEmcalJet*) fmyAKTChargedJets->At(fTPCJetUnbiased->GetSignalJetIndex(i));
c6202663 1724 jet_area_total+=AreaWithinTPC(fJetR,myJet->Eta());
1725 }
1726 }
62620fff 1727 delete jet_vec;
1728 delete track_vec;
1729
c6202663 1730 // Calculate Rho
1731 TPC_rho = E_tracks_total/(fTPCArea-jet_area_total);
1732
1733 // Fill Histogram
1734 fRhoChargedN->FillRho(fEventCentrality,TPC_rho);
1735 fRhoChargedN->FillBSJS(fEventCentrality,TPC_rho,fTPCJetThreshold,fmyAKTChargedJets,fTPCFullJet->GetJets(),fTPCFullJet->GetTotalJets());
1736 fRhoChargedN->FillDeltaPt(fEventCentrality,TPC_rho,fJetR,fTPCRCBckgFluc,1);
1737 fRhoChargedN->FillDeltaPtSignal(fEventCentrality,TPC_rho,fJetR,fTPCRCBckgFlucSignal,1);
7acc3e04 1738 fRhoChargedN->FillDeltaPtNColl(fEventCentrality,TPC_rho,fJetR,fTPCRCBckgFlucNColl,1);
c6202663 1739 fRhoChargedN->FillBackgroundFluctuations(fEventCentrality,TPC_rho,fJetR);
1740 fRhoChargedN->FillLeadingJetPtRho(fTPCFullJet->GetLeadingPt(),TPC_rho);
e864d416 1741
c6202663 1742}
e864d416 1743
c6202663 1744void AliAnalysisTaskFullpAJets::EstimateChargedRhoScale()
1745{
1746 Int_t i,j;
1747 Bool_t track_away_from_jet;
1748 Double_t E_tracks_total=0.0;
1749 Double_t TPC_rho=0.0;
1750 Double_t jet_area_total=0.0;
1751
62620fff 1752 TLorentzVector *jet_vec= new TLorentzVector;
1753 TLorentzVector *track_vec = new TLorentzVector;
1754
c6202663 1755 // First, sum all tracks within the EMCal that are away from jet(s) above Pt Threshold
1756 for (i=0;i<fnTracks;i++)
1757 {
1758 // First, check if track is in the EMCal!!
1759 AliVTrack* vtrack = (AliVTrack*) fmyTracks->At(i);
1760 if (IsInTPC(fJetR,vtrack->Phi(),vtrack->Eta(),kFALSE)==kTRUE)
ac6a3f1e 1761 {
8daeee93 1762 if (fTPCJetUnbiased->GetTotalSignalJets()<1)
ac6a3f1e 1763 {
c6202663 1764 E_tracks_total+=vtrack->Pt();
ac6a3f1e 1765 }
c6202663 1766 else
ac6a3f1e 1767 {
c6202663 1768 track_away_from_jet=kTRUE;
1769 j=0;
c6202663 1770 track_vec->SetPtEtaPhiE(vtrack->Pt(),vtrack->Eta(),vtrack->Phi(),vtrack->E());
8daeee93 1771 while (track_away_from_jet==kTRUE && j<fTPCJetUnbiased->GetTotalSignalJets())
ac6a3f1e 1772 {
8daeee93 1773 AliEmcalJet *myJet =(AliEmcalJet*) fmyAKTChargedJets->At(fTPCJetUnbiased->GetSignalJetIndex(j));
c6202663 1774 myJet->GetMom(*jet_vec);
1775 if (track_vec->DeltaR(*jet_vec)<=fJetRForRho)
1776 {
1777 track_away_from_jet=kFALSE;
1778 }
c6202663 1779 j++;
ac6a3f1e 1780 }
c6202663 1781 if (track_away_from_jet==kTRUE)
ac6a3f1e 1782 {
c6202663 1783 E_tracks_total+=vtrack->Pt();
ac6a3f1e 1784 }
1785 }
c6202663 1786 }
1787 }
1788
7acc3e04 1789 // Determine area of all Jets that are within the TPC
8daeee93 1790 if (fTPCJetUnbiased->GetTotalSignalJets()==0)
c6202663 1791 {
1792 jet_area_total=0.0;
1793 }
1794 else
1795 {
8daeee93 1796 for (i=0;i<fTPCJetUnbiased->GetTotalSignalJets();i++)
c6202663 1797 {
8daeee93 1798 AliEmcalJet *myJet = (AliEmcalJet*) fmyAKTChargedJets->At(fTPCJetUnbiased->GetSignalJetIndex(i));
c6202663 1799 jet_area_total+=AreaWithinTPC(fJetR,myJet->Eta());
1800 }
1801 }
62620fff 1802 delete jet_vec;
1803 delete track_vec;
1804
c6202663 1805 // Calculate Rho
1806 TPC_rho = E_tracks_total/(fTPCArea-jet_area_total);
1807 TPC_rho*=fScaleFactor;
1808
1809 // Fill Histogram
1810 fRhoChargedScale->FillRho(fEventCentrality,TPC_rho);
1811 fRhoChargedScale->FillBSJS(fEventCentrality,TPC_rho,fEMCalJetThreshold,fmyAKTFullJets,fEMCalFullJet->GetJets(),fEMCalFullJet->GetTotalJets());
1812 fRhoChargedScale->FillDeltaPt(fEventCentrality,TPC_rho,fJetR,fEMCalRCBckgFluc,1);
1813 fRhoChargedScale->FillDeltaPtSignal(fEventCentrality,TPC_rho,fJetR,fEMCalRCBckgFlucSignal,1);
7acc3e04 1814 fRhoChargedScale->FillDeltaPtNColl(fEventCentrality,TPC_rho,fJetR,fEMCalRCBckgFlucNColl,1);
c6202663 1815 fRhoChargedScale->FillBackgroundFluctuations(fEventCentrality,TPC_rho,fJetR);
1816 fRhoChargedScale->FillLeadingJetPtRho(fEMCalFullJet->GetLeadingPt(),TPC_rho);
e864d416 1817 fRhoChargedScale->FillMiscJetStats(fmyAKTFullJets,fEMCalFullJet->GetJets(),fEMCalFullJet->GetTotalJets(),fOrgTracks,fOrgClusters);
c6202663 1818
e864d416 1819}
62620fff 1820
c6202663 1821void AliAnalysisTaskFullpAJets::EstimateChargedRhokT()
1822{
1823 Int_t i;
1824 Double_t kTRho = 0.0;
1825 Double_t *pTArray = new Double_t[fTPCkTFullJet->GetTotalJets()];
1826 Double_t *RhoArray = new Double_t[fTPCkTFullJet->GetTotalJets()];
1827
1828 for (i=0;i<fTPCkTFullJet->GetTotalJets();i++)
1829 {
1830 AliEmcalJet *myJet =(AliEmcalJet*) fmyKTChargedJets->At(fTPCkTFullJet->GetJetIndex(i));
1831 pTArray[i]=myJet->Pt();
1832 RhoArray[i]=myJet->Pt()/myJet->Area();
1833 }
1834
1835 if (fTPCkTFullJet->GetTotalJets()>=2)
1836 {
1837 kTRho=MedianRhokT(pTArray,RhoArray,fTPCkTFullJet->GetTotalJets());
1838
1839 fRhoChargedkT->FillRho(fEventCentrality,kTRho);
1840 fRhoChargedkT->FillBSJS(fEventCentrality,kTRho,fTPCJetThreshold,fmyAKTChargedJets,fTPCFullJet->GetJets(),fTPCFullJet->GetTotalJets());
1841 fRhoChargedkT->FillDeltaPt(fEventCentrality,kTRho,fJetR,fTPCRCBckgFluc,1);
1842 fRhoChargedkT->FillDeltaPtSignal(fEventCentrality,kTRho,fJetR,fTPCRCBckgFlucSignal,1);
7acc3e04 1843 fRhoChargedkT->FillDeltaPtNColl(fEventCentrality,kTRho,fJetR,fTPCRCBckgFlucNColl,1);
c6202663 1844 fRhoChargedkT->FillBackgroundFluctuations(fEventCentrality,kTRho,fJetR);
1845 fRhoChargedkT->FillLeadingJetPtRho(fTPCFullJet->GetLeadingPt(),kTRho);
1846 }
1847 delete [] RhoArray;
1848 delete [] pTArray;
1849}
1850
1851void AliAnalysisTaskFullpAJets::EstimateChargedRhokTScale()
1852{
1853 Int_t i;
1854 Double_t kTRho = 0.0;
1855 Double_t *pTArray = new Double_t[fTPCkTFullJet->GetTotalJets()];
1856 Double_t *RhoArray = new Double_t[fTPCkTFullJet->GetTotalJets()];
1857
1858 for (i=0;i<fTPCkTFullJet->GetTotalJets();i++)
1859 {
1860 AliEmcalJet *myJet =(AliEmcalJet*) fmyKTChargedJets->At(fTPCkTFullJet->GetJetIndex(i));
1861 pTArray[i]=myJet->Pt();
1862 RhoArray[i]=myJet->Pt()/myJet->Area();
1863 }
1864
1865 if (fTPCkTFullJet->GetTotalJets()>=2)
1866 {
1867 kTRho=MedianRhokT(pTArray,RhoArray,fTPCkTFullJet->GetTotalJets());
1868 kTRho*=fScaleFactor;
1869
1870 fRhoChargedkTScale->FillRho(fEventCentrality,kTRho);
1871 fRhoChargedkTScale->FillBSJS(fEventCentrality,kTRho,fEMCalJetThreshold,fmyAKTFullJets,fEMCalFullJet->GetJets(),fEMCalFullJet->GetTotalJets());
1872 fRhoChargedkTScale->FillDeltaPt(fEventCentrality,kTRho,fJetR,fEMCalRCBckgFluc,1);
1873 fRhoChargedkTScale->FillDeltaPtSignal(fEventCentrality,kTRho,fJetR,fEMCalRCBckgFlucSignal,1);
7acc3e04 1874 fRhoChargedkTScale->FillDeltaPtNColl(fEventCentrality,kTRho,fJetR,fEMCalRCBckgFlucNColl,1);
c6202663 1875 fRhoChargedkTScale->FillBackgroundFluctuations(fEventCentrality,kTRho,fJetR);
1876 fRhoChargedkTScale->FillLeadingJetPtRho(fEMCalFullJet->GetLeadingPt(),kTRho);
1877 }
1878 delete [] RhoArray;
1879 delete [] pTArray;
1880}
1881
1882void AliAnalysisTaskFullpAJets::EstimateChargedRhoCMS()
1883{
1884 Int_t i,k;
1885 Double_t kTRho = 0.0;
1886 Double_t CMSTotalkTArea = 0.0;
1887 Double_t CMSTrackArea = 0.0;
1888 Double_t CMSCorrectionFactor = 1.0;
1889 Double_t *pTArray = new Double_t[fTPCkTFullJet->GetTotalJets()];
1890 Double_t *RhoArray = new Double_t[fTPCkTFullJet->GetTotalJets()];
1891
1892 k=0;
1893 if ((fTPCJet->GetLeadingPt()>=fTPCJetThreshold) && (fTPCJet->GetSubLeadingPt()>=fTPCJetThreshold))
1894 {
1895 AliEmcalJet *myJet1 =(AliEmcalJet*) fmyAKTChargedJets->At(fTPCJet->GetLeadingIndex());
1896 AliEmcalJet *myJet2 =(AliEmcalJet*) fmyAKTChargedJets->At(fTPCJet->GetSubLeadingIndex());
1897
1898 for (i=0;i<fTPCkTFullJet->GetTotalJets();i++)
1899 {
1900 AliEmcalJet *myJet =(AliEmcalJet*) fmyKTChargedJets->At(fTPCkTFullJet->GetJetIndex(i));
ac6a3f1e 1901
c6202663 1902 CMSTotalkTArea+=myJet->Area();
1903 if (myJet->GetNumberOfTracks()>0)
1904 {
1905 CMSTrackArea+=myJet->Area();
1906 }
1907 if (IsJetOverlap(myJet,myJet1,kFALSE)==kFALSE && IsJetOverlap(myJet,myJet2,kFALSE)==kFALSE)
ac6a3f1e 1908 {
c6202663 1909 pTArray[k]=myJet->Pt();
1910 RhoArray[k]=myJet->Pt()/myJet->Area();
1911 k++;
ac6a3f1e 1912 }
1913 }
c6202663 1914 if (k>0)
ac6a3f1e 1915 {
c6202663 1916 kTRho=MedianRhokT(pTArray,RhoArray,k);
ac6a3f1e 1917 }
c6202663 1918 else
ac6a3f1e 1919 {
c6202663 1920 kTRho=0.0;
ac6a3f1e 1921 }
1922 }
c6202663 1923 else if (fTPCJet->GetLeadingPt()>=fTPCJetThreshold)
ac6a3f1e 1924 {
c6202663 1925 AliEmcalJet *myJet1 =(AliEmcalJet*) fmyAKTChargedJets->At(fTPCJet->GetLeadingIndex());
ac6a3f1e 1926
c6202663 1927 for (i=0;i<fTPCkTFullJet->GetTotalJets();i++)
ac6a3f1e 1928 {
c6202663 1929 AliEmcalJet *myJet =(AliEmcalJet*) fmyKTChargedJets->At(fTPCkTFullJet->GetJetIndex(i));
1930
1931 CMSTotalkTArea+=myJet->Area();
1932 if (myJet->GetNumberOfTracks()>0)
1933 {
1934 CMSTrackArea+=myJet->Area();
1935 }
1936 if (IsJetOverlap(myJet,myJet1,kFALSE)==kFALSE)
1937 {
1938 pTArray[k]=myJet->Pt();
1939 RhoArray[k]=myJet->Pt()/myJet->Area();
1940 k++;
1941 }
1942 }
1943 if (k>0)
1944 {
1945 kTRho=MedianRhokT(pTArray,RhoArray,k);
1946 }
1947 else
1948 {
1949 kTRho=0.0;
ac6a3f1e 1950 }
1951 }
c6202663 1952 else
78246241 1953 {
c6202663 1954 for (i=0;i<fTPCkTFullJet->GetTotalJets();i++)
1955 {
1956 AliEmcalJet *myJet =(AliEmcalJet*) fmyKTChargedJets->At(fTPCkTFullJet->GetJetIndex(i));
1957
1958 CMSTotalkTArea+=myJet->Area();
1959 if (myJet->GetNumberOfTracks()>0)
1960 {
1961 CMSTrackArea+=myJet->Area();
1962 }
1963 pTArray[k]=myJet->Pt();
1964 RhoArray[k]=myJet->Pt()/myJet->Area();
1965 k++;
1966 }
1967 if (k>0)
1968 {
1969 kTRho=MedianRhokT(pTArray,RhoArray,k);
1970 }
1971 else
1972 {
1973 kTRho=0.0;
1974 }
78246241 1975 }
c6202663 1976 // Scale CMS Rho by Correction factor
1977 if (CMSTotalkTArea==0.0)
ac6a3f1e 1978 {
c6202663 1979 CMSCorrectionFactor = 1.0;
ac6a3f1e 1980 }
c6202663 1981 else
ac6a3f1e 1982 {
c6202663 1983 //CMSCorrectionFactor = CMSTrackArea/CMSTotalkTArea;
1984 CMSCorrectionFactor = CMSTrackArea/(fTPCPhiTotal*(fTPCEtaTotal-2*fJetR)); // The total physical area should be reduced by the eta cut due to looping over only fully contained kT jets within the TPC
1985 }
1986 kTRho*=CMSCorrectionFactor;
1987 fRhoChargedCMS->FillRho(fEventCentrality,kTRho);
1988 fRhoChargedCMS->FillBSJS(fEventCentrality,kTRho,fTPCJetThreshold,fmyAKTChargedJets,fTPCFullJet->GetJets(),fTPCFullJet->GetTotalJets());
1989 fRhoChargedCMS->FillDeltaPt(fEventCentrality,kTRho,fJetR,fTPCRCBckgFluc,1);
1990 fRhoChargedCMS->FillDeltaPtSignal(fEventCentrality,kTRho,fJetR,fTPCRCBckgFlucSignal,1);
7acc3e04 1991 fRhoChargedCMS->FillDeltaPtNColl(fEventCentrality,kTRho,fJetR,fTPCRCBckgFlucNColl,1);
c6202663 1992 fRhoChargedCMS->FillBackgroundFluctuations(fEventCentrality,kTRho,fJetR);
1993 fRhoChargedCMS->FillLeadingJetPtRho(fTPCFullJet->GetLeadingPt(),kTRho);
1994 delete [] RhoArray;
1995 delete [] pTArray;
1996}
62620fff 1997
c6202663 1998void AliAnalysisTaskFullpAJets::EstimateChargedRhoCMSScale()
1999{
2000 Int_t i,k;
2001 Double_t kTRho = 0.0;
2002 Double_t CMSTotalkTArea = 0.0;
2003 Double_t CMSTrackArea = 0.0;
2004 Double_t CMSCorrectionFactor = 1.0;
2005 Double_t *pTArray = new Double_t[fTPCkTFullJet->GetTotalJets()];
2006 Double_t *RhoArray = new Double_t[fTPCkTFullJet->GetTotalJets()];
2007
2008 k=0;
2009 if ((fTPCJet->GetLeadingPt()>=fTPCJetThreshold) && (fTPCJet->GetSubLeadingPt()>=fTPCJetThreshold))
2010 {
2011 AliEmcalJet *myJet1 =(AliEmcalJet*) fmyAKTChargedJets->At(fTPCJet->GetLeadingIndex());
2012 AliEmcalJet *myJet2 =(AliEmcalJet*) fmyAKTChargedJets->At(fTPCJet->GetSubLeadingIndex());
2013
2014 for (i=0;i<fTPCkTFullJet->GetTotalJets();i++)
2015 {
2016 AliEmcalJet *myJet =(AliEmcalJet*) fmyKTChargedJets->At(fTPCkTFullJet->GetJetIndex(i));
2017
2018 CMSTotalkTArea+=myJet->Area();
2019 if (myJet->GetNumberOfTracks()>0)
2020 {
2021 CMSTrackArea+=myJet->Area();
2022 }
2023 if (IsJetOverlap(myJet,myJet1,kFALSE)==kFALSE && IsJetOverlap(myJet,myJet2,kFALSE)==kFALSE)
2024 {
2025 pTArray[k]=myJet->Pt();
2026 RhoArray[k]=myJet->Pt()/myJet->Area();
2027 k++;
2028 }
2029 }
2030 if (k>0)
2031 {
2032 kTRho=MedianRhokT(pTArray,RhoArray,k);
2033 }
2034 else
2035 {
2036 kTRho=0.0;
2037 }
2038 }
2039 else if (fTPCJet->GetLeadingPt()>=fTPCJetThreshold)
2040 {
2041 AliEmcalJet *myJet1 =(AliEmcalJet*) fmyAKTChargedJets->At(fTPCJet->GetLeadingIndex());
2042
2043 for (i=0;i<fTPCkTFullJet->GetTotalJets();i++)
2044 {
2045 AliEmcalJet *myJet =(AliEmcalJet*) fmyKTChargedJets->At(fTPCkTFullJet->GetJetIndex(i));
2046
2047 CMSTotalkTArea+=myJet->Area();
2048 if (myJet->GetNumberOfTracks()>0)
2049 {
2050 CMSTrackArea+=myJet->Area();
2051 }
2052 if (IsJetOverlap(myJet,myJet1,kFALSE)==kFALSE)
2053 {
2054 pTArray[k]=myJet->Pt();
2055 RhoArray[k]=myJet->Pt()/myJet->Area();
2056 k++;
2057 }
2058 }
2059 if (k>0)
2060 {
2061 kTRho=MedianRhokT(pTArray,RhoArray,k);
2062 }
2063 else
2064 {
2065 kTRho=0.0;
2066 }
2067 }
2068 else
ac6a3f1e 2069 {
c6202663 2070 for (i=0;i<fTPCkTFullJet->GetTotalJets();i++)
2071 {
2072 AliEmcalJet *myJet =(AliEmcalJet*) fmyKTChargedJets->At(fTPCkTFullJet->GetJetIndex(i));
2073
2074 CMSTotalkTArea+=myJet->Area();
2075 if (myJet->GetNumberOfTracks()>0)
2076 {
2077 CMSTrackArea+=myJet->Area();
2078 }
2079 pTArray[k]=myJet->Pt();
2080 RhoArray[k]=myJet->Pt()/myJet->Area();
2081 k++;
2082 }
2083 if (k>0)
78246241 2084 {
c6202663 2085 kTRho=MedianRhokT(pTArray,RhoArray,k);
78246241 2086 }
c6202663 2087 else
c54b626a 2088 {
c6202663 2089 kTRho=0.0;
c54b626a 2090 }
ac6a3f1e 2091 }
c6202663 2092 kTRho*=fScaleFactor;
2093 // Scale CMS Rho by Correction factor
2094 if (CMSTotalkTArea==0.0)
2095 {
2096 CMSCorrectionFactor = 1.0;
2097 }
2098 else
2099 {
c6202663 2100 CMSCorrectionFactor = CMSTrackArea/(fTPCPhiTotal*(fTPCEtaTotal-2*fJetR)); // The total physical area should be reduced by the eta cut due to looping over only fully contained kT jets within the TPC
2101 }
2102 kTRho*=CMSCorrectionFactor;
2103
2104 fRhoChargedCMSScale->FillRho(fEventCentrality,kTRho);
2105 fRhoChargedCMSScale->FillBSJS(fEventCentrality,kTRho,fEMCalJetThreshold,fmyAKTFullJets,fEMCalFullJet->GetJets(),fEMCalFullJet->GetTotalJets());
2106 fRhoChargedCMSScale->FillDeltaPt(fEventCentrality,kTRho,fJetR,fEMCalRCBckgFluc,1);
2107 fRhoChargedCMSScale->FillDeltaPtSignal(fEventCentrality,kTRho,fJetR,fEMCalRCBckgFlucSignal,1);
7acc3e04 2108 fRhoChargedCMSScale->FillDeltaPtNColl(fEventCentrality,kTRho,fJetR,fEMCalRCBckgFlucNColl,1);
c6202663 2109 fRhoChargedCMSScale->FillBackgroundFluctuations(fEventCentrality,kTRho,fJetR);
2110 fRhoChargedCMSScale->FillLeadingJetPtRho(fEMCalFullJet->GetLeadingPt(),kTRho);
08b981da 2111 fRhoChargedCMSScale->DoNEFAnalysis(fNEFSignalJetCut,fEMCalJetThreshold,fmyAKTFullJets,fEMCalFullJet->GetJets(),fEMCalFullJet->GetTotalJets(),fmyClusters,fOrgClusters,fEvent,fEMCALGeometry,fRecoUtil,fCells);
e864d416 2112 fRhoChargedCMSScale->FillMiscJetStats(fmyAKTFullJets,fEMCalFullJet->GetJets(),fEMCalFullJet->GetTotalJets(),fOrgTracks,fOrgClusters);
2113
c6202663 2114 delete [] RhoArray;
2115 delete [] pTArray;
ac6a3f1e 2116}
62620fff 2117
c6202663 2118// Full Rho's
2119void AliAnalysisTaskFullpAJets::EstimateFullRho0()
ac6a3f1e 2120{
2121 Int_t i;
c6202663 2122 Double_t E_tracks_total=0.0;
2123 Double_t E_caloclusters_total=0.0;
2124 Double_t EMCal_rho=0.0;
ac6a3f1e 2125
62620fff 2126 TLorentzVector *cluster_vec = new TLorentzVector;
2127
ac6a3f1e 2128 // Loop over all tracks
2129 for (i=0;i<fnTracks;i++)
2130 {
2131 AliVTrack* vtrack =(AliVTrack*) fmyTracks->At(i);
2132 if (IsInEMCal(vtrack->Phi(),vtrack->Eta())==kTRUE)
2133 {
2134 E_tracks_total+=vtrack->Pt();
2135 }
2136 }
2137
2138 // Loop over all caloclusters
2139 for (i=0;i<fnClusters;i++)
2140 {
2141 AliVCluster* vcluster = (AliVCluster*) fmyClusters->At(i);
8daeee93 2142 vcluster->GetMomentum(*cluster_vec,fVertex);
c6202663 2143 E_caloclusters_total+=cluster_vec->Pt();
ac6a3f1e 2144 }
62620fff 2145 delete cluster_vec;
2146
ac6a3f1e 2147 // Calculate the mean Background density
2148 EMCal_rho=(E_tracks_total+E_caloclusters_total)/fEMCalArea;
c6202663 2149 fRhoFull=EMCal_rho;
ac6a3f1e 2150
c6202663 2151 // Fill Histograms
c6202663 2152 fRhoFull0->FillRho(fEventCentrality,EMCal_rho);
2153 fRhoFull0->FillBSJS(fEventCentrality,EMCal_rho,fEMCalJetThreshold,fmyAKTFullJets,fEMCalFullJet->GetJets(),fEMCalFullJet->GetTotalJets());
2154 fRhoFull0->FillDeltaPt(fEventCentrality,EMCal_rho,fJetR,fEMCalRCBckgFluc,1);
2155 fRhoFull0->FillDeltaPtSignal(fEventCentrality,EMCal_rho,fJetR,fEMCalRCBckgFlucSignal,1);
7acc3e04 2156 fRhoFull0->FillDeltaPtNColl(fEventCentrality,EMCal_rho,fJetR,fEMCalRCBckgFlucNColl,1);
c6202663 2157 fRhoFull0->FillBackgroundFluctuations(fEventCentrality,EMCal_rho,fJetR);
2158 fRhoFull0->FillLeadingJetPtRho(fEMCalFullJet->GetLeadingPt(),EMCal_rho);
ac6a3f1e 2159}
2160
c6202663 2161void AliAnalysisTaskFullpAJets::EstimateFullRho1()
ac6a3f1e 2162{
2163 Int_t i;
c6202663 2164 Double_t E_tracks_total=0.0;
2165 Double_t E_caloclusters_total=0.0;
2166 Double_t EMCal_rho=0.0;
ac6a3f1e 2167
62620fff 2168 TLorentzVector *temp_jet= new TLorentzVector;
2169 TLorentzVector *track_vec = new TLorentzVector;
2170 TLorentzVector *cluster_vec = new TLorentzVector;
2171
8daeee93 2172 if (fEMCalPartJetUnbiased->GetLeadingPt()>=fEMCalJetThreshold)
ac6a3f1e 2173 {
8daeee93 2174 AliEmcalJet *myJet =(AliEmcalJet*) fmyAKTFullJets->At(fEMCalPartJetUnbiased->GetLeadingIndex());
ac6a3f1e 2175 myJet->GetMom(*temp_jet);
2176
2177 // Loop over all tracks
2178 for (i=0;i<fnTracks;i++)
2179 {
2180 AliVTrack* vtrack =(AliVTrack*) fmyTracks->At(i);
2181 if (IsInEMCal(vtrack->Phi(),vtrack->Eta())==kTRUE)
2182 {
ac6a3f1e 2183 track_vec->SetPtEtaPhiE(vtrack->Pt(),vtrack->Eta(),vtrack->Phi(),vtrack->E());
c6202663 2184 if (temp_jet->DeltaR(*track_vec)>fJetRForRho)
ac6a3f1e 2185 {
2186 E_tracks_total+=vtrack->Pt();
2187 }
ac6a3f1e 2188 }
2189 }
2190
2191 // Loop over all caloclusters
2192 for (i=0;i<fnClusters;i++)
2193 {
2194 AliVCluster* vcluster = (AliVCluster*) fmyClusters->At(i);
8daeee93 2195 vcluster->GetMomentum(*cluster_vec,fVertex);
c6202663 2196 if (temp_jet->DeltaR(*cluster_vec)>fJetRForRho)
ac6a3f1e 2197 {
2198 E_caloclusters_total+=vcluster->E();
2199 }
ac6a3f1e 2200 }
ac6a3f1e 2201 // Calculate the mean Background density
c6202663 2202 EMCal_rho=(E_tracks_total+E_caloclusters_total)/(fEMCalArea-AreaWithinEMCal(fJetR,myJet->Phi(),myJet->Eta()));
ac6a3f1e 2203 }
c6202663 2204 else // i.e. No signal jets -> same as total background density
ac6a3f1e 2205 {
2206 // Loop over all tracks
2207 for (i=0;i<fnTracks;i++)
2208 {
2209 AliVTrack* vtrack =(AliVTrack*) fmyTracks->At(i);
2210 if (IsInEMCal(vtrack->Phi(),vtrack->Eta())==kTRUE)
2211 {
2212 E_tracks_total+=vtrack->Pt();
2213 }
2214 }
2215
2216 // Loop over all caloclusters
2217 for (i=0;i<fnClusters;i++)
2218 {
2219 AliVCluster* vcluster = (AliVCluster*) fmyClusters->At(i);
2220 E_caloclusters_total+=vcluster->E();
2221 }
2222 // Calculate the mean Background density
2223 EMCal_rho=(E_tracks_total+E_caloclusters_total)/fEMCalArea;
2224 }
62620fff 2225 delete temp_jet;
2226 delete track_vec;
2227 delete cluster_vec;
2228
c6202663 2229 // Fill histograms
2230 fRhoFull1->FillRho(fEventCentrality,EMCal_rho);
2231 fRhoFull1->FillBSJS(fEventCentrality,EMCal_rho,fEMCalJetThreshold,fmyAKTFullJets,fEMCalFullJet->GetJets(),fEMCalFullJet->GetTotalJets());
2232 fRhoFull1->FillDeltaPt(fEventCentrality,EMCal_rho,fJetR,fEMCalRCBckgFluc,1);
2233 fRhoFull1->FillDeltaPtSignal(fEventCentrality,EMCal_rho,fJetR,fEMCalRCBckgFlucSignal,1);
7acc3e04 2234 fRhoFull1->FillDeltaPtNColl(fEventCentrality,EMCal_rho,fJetR,fEMCalRCBckgFlucNColl,1);
c6202663 2235 fRhoFull1->FillBackgroundFluctuations(fEventCentrality,EMCal_rho,fJetR);
2236 fRhoFull1->FillLeadingJetPtRho(fEMCalFullJet->GetLeadingPt(),EMCal_rho);
ac6a3f1e 2237}
2238
c6202663 2239void AliAnalysisTaskFullpAJets::EstimateFullRho2()
ac6a3f1e 2240{
2241 Int_t i;
c6202663 2242 Double_t E_tracks_total=0.0;
2243 Double_t E_caloclusters_total=0.0;
2244 Double_t EMCal_rho=0.0;
2245
62620fff 2246 TLorentzVector *temp_jet1 = new TLorentzVector;
2247 TLorentzVector *temp_jet2 = new TLorentzVector;
2248 TLorentzVector *track_vec = new TLorentzVector;
2249 TLorentzVector *cluster_vec = new TLorentzVector;
2250
8daeee93 2251 if ((fEMCalPartJetUnbiased->GetLeadingPt()>=fEMCalJetThreshold) && (fEMCalPartJetUnbiased->GetSubLeadingPt()>=fEMCalJetThreshold))
ac6a3f1e 2252 {
8daeee93 2253 AliEmcalJet *myhJet =(AliEmcalJet*) fmyAKTFullJets->At(fEMCalPartJetUnbiased->GetLeadingIndex());
c6202663 2254 myhJet->GetMom(*temp_jet1);
2255
8daeee93 2256 AliEmcalJet *myJet =(AliEmcalJet*) fmyAKTFullJets->At(fEMCalPartJetUnbiased->GetSubLeadingIndex());
c6202663 2257 myJet->GetMom(*temp_jet2);
2258
2259 // Loop over all tracks
2260 for (i=0;i<fnTracks;i++)
2261 {
2262 AliVTrack* vtrack =(AliVTrack*) fmyTracks->At(i);
2263 if (IsInEMCal(vtrack->Phi(),vtrack->Eta())==kTRUE)
2264 {
c6202663 2265 track_vec->SetPtEtaPhiE(vtrack->Pt(),vtrack->Eta(),vtrack->Phi(),vtrack->E());
2266 if ((temp_jet1->DeltaR(*track_vec)>fJetRForRho) && (temp_jet2->DeltaR(*track_vec)>fJetRForRho))
2267 {
2268 E_tracks_total+=vtrack->Pt();
2269 }
c6202663 2270 }
2271 }
2272
2273 // Loop over all caloclusters
2274 for (i=0;i<fnClusters;i++)
2275 {
2276 AliVCluster* vcluster = (AliVCluster*) fmyClusters->At(i);
8daeee93 2277 vcluster->GetMomentum(*cluster_vec,fVertex);
c6202663 2278 if ((temp_jet1->DeltaR(*cluster_vec)>fJetRForRho) && (temp_jet2->DeltaR(*cluster_vec)>fJetRForRho))
2279 {
2280 E_caloclusters_total+=vcluster->E();
2281 }
c6202663 2282 }
62620fff 2283
c6202663 2284 // Calculate the mean Background density
2285 EMCal_rho=(E_tracks_total+E_caloclusters_total)/(fEMCalArea-AreaWithinEMCal(fJetR,myhJet->Phi(),myhJet->Eta())-AreaWithinEMCal(fJetR,myJet->Phi(),myJet->Eta()));
2286 }
8daeee93 2287 else if (fEMCalPartJetUnbiased->GetLeadingPt()>=fEMCalJetThreshold)
c6202663 2288 {
8daeee93 2289 AliEmcalJet *myJet =(AliEmcalJet*) fmyAKTFullJets->At(fEMCalPartJetUnbiased->GetLeadingIndex());
62620fff 2290 myJet->GetMom(*temp_jet1);
ac6a3f1e 2291
c6202663 2292 // Loop over all tracks
2293 for (i=0;i<fnTracks;i++)
ac6a3f1e 2294 {
c6202663 2295 AliVTrack* vtrack =(AliVTrack*) fmyTracks->At(i);
2296 if (IsInEMCal(vtrack->Phi(),vtrack->Eta())==kTRUE)
ac6a3f1e 2297 {
c6202663 2298 track_vec->SetPtEtaPhiE(vtrack->Pt(),vtrack->Eta(),vtrack->Phi(),vtrack->E());
62620fff 2299 if (temp_jet1->DeltaR(*track_vec)>fJetRForRho)
ac6a3f1e 2300 {
c6202663 2301 E_tracks_total+=vtrack->Pt();
ac6a3f1e 2302 }
c6202663 2303 }
2304 }
2305
2306 // Loop over all caloclusters
2307 for (i=0;i<fnClusters;i++)
2308 {
2309 AliVCluster* vcluster = (AliVCluster*) fmyClusters->At(i);
8daeee93 2310 vcluster->GetMomentum(*cluster_vec,fVertex);
62620fff 2311 if (temp_jet1->DeltaR(*cluster_vec)>fJetRForRho)
c6202663 2312 {
2313 E_caloclusters_total+=vcluster->E();
ac6a3f1e 2314 }
2315 }
c6202663 2316 // Calculate the mean Background density
2317 EMCal_rho=(E_tracks_total+E_caloclusters_total)/(fEMCalArea-AreaWithinEMCal(fJetR,myJet->Phi(),myJet->Eta()));
2318 }
2319 else // i.e. No signal jets -> same as total background density
2320 {
2321 // Loop over all tracks
2322 for (i=0;i<fnTracks;i++)
2323 {
2324 AliVTrack* vtrack =(AliVTrack*) fmyTracks->At(i);
2325 if (IsInEMCal(vtrack->Phi(),vtrack->Eta())==kTRUE)
2326 {
2327 E_tracks_total+=vtrack->Pt();
2328 }
2329 }
2330
2331 // Loop over all caloclusters
2332 for (i=0;i<fnClusters;i++)
2333 {
2334 AliVCluster* vcluster = (AliVCluster*) fmyClusters->At(i);
2335 E_caloclusters_total+=vcluster->E();
2336 }
2337 // Calculate the mean Background density
2338 EMCal_rho=(E_tracks_total+E_caloclusters_total)/fEMCalArea;
ac6a3f1e 2339 }
62620fff 2340 delete temp_jet1;
2341 delete temp_jet2;
2342 delete track_vec;
2343 delete cluster_vec;
2344
c6202663 2345 // Fill histograms
2346 fRhoFull2->FillRho(fEventCentrality,EMCal_rho);
2347 fRhoFull2->FillBSJS(fEventCentrality,EMCal_rho,fEMCalJetThreshold,fmyAKTFullJets,fEMCalFullJet->GetJets(),fEMCalFullJet->GetTotalJets());
2348 fRhoFull2->FillDeltaPt(fEventCentrality,EMCal_rho,fJetR,fEMCalRCBckgFluc,1);
2349 fRhoFull2->FillDeltaPtSignal(fEventCentrality,EMCal_rho,fJetR,fEMCalRCBckgFlucSignal,1);
7acc3e04 2350 fRhoFull2->FillDeltaPtNColl(fEventCentrality,EMCal_rho,fJetR,fEMCalRCBckgFlucNColl,1);
c6202663 2351 fRhoFull2->FillBackgroundFluctuations(fEventCentrality,EMCal_rho,fJetR);
2352 fRhoFull2->FillLeadingJetPtRho(fEMCalFullJet->GetLeadingPt(),EMCal_rho);
ac6a3f1e 2353}
2354
c6202663 2355void AliAnalysisTaskFullpAJets::EstimateFullRhoN()
ac6a3f1e 2356{
2357 Int_t i,j;
2358 Bool_t track_away_from_jet;
c6202663 2359 Bool_t cluster_away_from_jet;
2360 Double_t E_tracks_total=0.0;
2361 Double_t E_caloclusters_total=0.0;
2362 Double_t EMCal_rho=0.0;
2363 Double_t jet_area_total=0.0;
ac6a3f1e 2364
62620fff 2365 TLorentzVector *jet_vec= new TLorentzVector;
2366 TLorentzVector *track_vec = new TLorentzVector;
2367 TLorentzVector *cluster_vec = new TLorentzVector;
2368
ac6a3f1e 2369 // First, sum all tracks within the EMCal that are away from jet(s) above Pt Threshold
ac6a3f1e 2370 for (i=0;i<fnTracks;i++)
2371 {
2372 // First, check if track is in the EMCal!!
c6202663 2373 AliVTrack* vtrack = (AliVTrack*) fmyTracks->At(i);
ac6a3f1e 2374 if (IsInEMCal(vtrack->Phi(),vtrack->Eta())==kTRUE)
2375 {
8daeee93 2376 if (fEMCalPartJetUnbiased->GetTotalSignalJets()<1)
ac6a3f1e 2377 {
2378 E_tracks_total+=vtrack->Pt();
2379 }
2380 else
2381 {
2382 track_away_from_jet=kTRUE;
2383 j=0;
ac6a3f1e 2384 track_vec->SetPtEtaPhiE(vtrack->Pt(),vtrack->Eta(),vtrack->Phi(),vtrack->E());
8daeee93 2385 while (track_away_from_jet==kTRUE && j<fEMCalPartJetUnbiased->GetTotalSignalJets())
ac6a3f1e 2386 {
8daeee93 2387 AliEmcalJet *myJet =(AliEmcalJet*) fmyAKTFullJets->At(fEMCalPartJetUnbiased->GetSignalJetIndex(j));
ac6a3f1e 2388 myJet->GetMom(*jet_vec);
c6202663 2389 if (track_vec->DeltaR(*jet_vec)<=fJetRForRho)
ac6a3f1e 2390 {
2391 track_away_from_jet=kFALSE;
2392 }
ac6a3f1e 2393 j++;
2394 }
2395 if (track_away_from_jet==kTRUE)
2396 {
2397 E_tracks_total+=vtrack->Pt();
2398 }
ac6a3f1e 2399 }
2400 }
2401 }
2402
2403 // Next, sum all CaloClusters within the EMCal (obviously all clusters must be within EMCal!!) that are away from jet(s) above Pt Threshold
ac6a3f1e 2404 for (i=0;i<fnClusters;i++)
2405 {
c6202663 2406 AliVCluster* vcluster = (AliVCluster*) fmyClusters->At(i);
2407 if (fEMCalPartJet->GetTotalSignalJets()<1)
ac6a3f1e 2408 {
2409 E_caloclusters_total+=vcluster->E();
2410 }
2411 else
2412 {
2413 cluster_away_from_jet=kTRUE;
2414 j=0;
2415
8daeee93 2416 vcluster->GetMomentum(*cluster_vec,fVertex);
2417 while (cluster_away_from_jet==kTRUE && j<fEMCalPartJetUnbiased->GetTotalSignalJets())
ac6a3f1e 2418 {
8daeee93 2419 AliEmcalJet *myJet =(AliEmcalJet*) fmyAKTFullJets->At(fEMCalPartJetUnbiased->GetSignalJetIndex(j));
ac6a3f1e 2420 myJet->GetMom(*jet_vec);
c6202663 2421 if (cluster_vec->DeltaR(*jet_vec)<=fJetRForRho)
ac6a3f1e 2422 {
2423 cluster_away_from_jet=kFALSE;
2424 }
ac6a3f1e 2425 j++;
2426 }
2427 if (cluster_away_from_jet==kTRUE)
2428 {
2429 E_caloclusters_total+=vcluster->E();
2430 }
ac6a3f1e 2431 }
2432 }
2433
2434 // Determine area of all Jets that are within the EMCal
c6202663 2435 if (fEMCalPartJet->GetTotalSignalJets()==0)
ac6a3f1e 2436 {
c6202663 2437 jet_area_total=0.0;
ac6a3f1e 2438 }
2439 else
2440 {
8daeee93 2441 for (i=0;i<fEMCalPartJetUnbiased->GetTotalSignalJets();i++)
ac6a3f1e 2442 {
8daeee93 2443 AliEmcalJet *myJet = (AliEmcalJet*) fmyAKTFullJets->At(fEMCalPartJetUnbiased->GetSignalJetIndex(i));
ac6a3f1e 2444 jet_area_total+=AreaWithinEMCal(fJetR,myJet->Phi(),myJet->Eta());
2445 }
2446 }
62620fff 2447 delete jet_vec;
2448 delete track_vec;
2449 delete cluster_vec;
2450
ac6a3f1e 2451 // Calculate Rho
2452 EMCal_rho=(E_tracks_total+E_caloclusters_total)/(fEMCalArea-jet_area_total);
ac6a3f1e 2453
c6202663 2454 // Fill Histogram
2455 fRhoFullN->FillRho(fEventCentrality,EMCal_rho);
2456 fRhoFullN->FillBSJS(fEventCentrality,EMCal_rho,fEMCalJetThreshold,fmyAKTFullJets,fEMCalFullJet->GetJets(),fEMCalFullJet->GetTotalJets());
2457 fRhoFullN->FillDeltaPt(fEventCentrality,EMCal_rho,fJetR,fEMCalRCBckgFluc,1);
2458 fRhoFullN->FillDeltaPtSignal(fEventCentrality,EMCal_rho,fJetR,fEMCalRCBckgFlucSignal,1);
7acc3e04 2459 fRhoFullN->FillDeltaPtNColl(fEventCentrality,EMCal_rho,fJetR,fEMCalRCBckgFlucNColl,1);
c6202663 2460 fRhoFullN->FillBackgroundFluctuations(fEventCentrality,EMCal_rho,fJetR);
2461 fRhoFullN->FillLeadingJetPtRho(fEMCalFullJet->GetLeadingPt(),EMCal_rho);
ac6a3f1e 2462}
2463
c6202663 2464void AliAnalysisTaskFullpAJets::EstimateFullRhoDijet()
ac6a3f1e 2465{
78246241 2466 Int_t i;
78246241 2467 Double_t E_tracks_total=0.0;
2468 Double_t E_caloclusters_total=0.0;
2469 Double_t EMCal_rho=0.0;
78246241 2470
2471 // Loop over all tracks
ac6a3f1e 2472 for (i=0;i<fnTracks;i++)
2473 {
78246241 2474 AliVTrack* vtrack =(AliVTrack*) fmyTracks->At(i);
ac6a3f1e 2475 if (IsInEMCal(vtrack->Phi(),vtrack->Eta())==kTRUE)
2476 {
78246241 2477 E_tracks_total+=vtrack->Pt();
ac6a3f1e 2478 }
2479 }
2480
78246241 2481 // Loop over all caloclusters
ac6a3f1e 2482 for (i=0;i<fnClusters;i++)
2483 {
78246241 2484 AliVCluster* vcluster = (AliVCluster*) fmyClusters->At(i);
2485 E_caloclusters_total+=vcluster->E();
ac6a3f1e 2486 }
2487
78246241 2488 // Calculate the mean Background density
2489 EMCal_rho=(E_tracks_total+E_caloclusters_total)/fEMCalArea;
ac6a3f1e 2490
c6202663 2491 // Fill Histograms
2492 fRhoFullDijet->FillRho(fEventCentrality,EMCal_rho);
2493 fRhoFullDijet->FillBSJS(fEventCentrality,EMCal_rho,fEMCalJetThreshold,fmyAKTFullJets,fEMCalFullJet->GetJets(),fEMCalFullJet->GetTotalJets());
2494 fRhoFullDijet->FillDeltaPt(fEventCentrality,EMCal_rho,fJetR,fEMCalRCBckgFluc,1);
2495 fRhoFullDijet->FillDeltaPtSignal(fEventCentrality,EMCal_rho,fJetR,fEMCalRCBckgFlucSignal,1);
7acc3e04 2496 fRhoFullDijet->FillDeltaPtNColl(fEventCentrality,EMCal_rho,fJetR,fEMCalRCBckgFlucNColl,1);
c6202663 2497 fRhoFullDijet->FillBackgroundFluctuations(fEventCentrality,EMCal_rho,fJetR);
2498 fRhoFullDijet->FillLeadingJetPtRho(fEMCalFullJet->GetLeadingPt(),EMCal_rho);
2499}
2500
2501void AliAnalysisTaskFullpAJets::EstimateFullRhokT()
2502{
2503 Int_t i;
2504 Double_t kTRho = 0.0;
2505 Double_t *pTArray = new Double_t[fEMCalkTFullJet->GetTotalJets()];
2506 Double_t *RhoArray = new Double_t[fEMCalkTFullJet->GetTotalJets()];
2507
2508 for (i=0;i<fEMCalkTFullJet->GetTotalJets();i++)
2509 {
2510 AliEmcalJet *myJet =(AliEmcalJet*) fmyKTFullJets->At(fEMCalkTFullJet->GetJetIndex(i));
2511 pTArray[i]=myJet->Pt();
2512 RhoArray[i]=myJet->Pt()/myJet->Area();
2513 }
ac6a3f1e 2514
c6202663 2515 if (fEMCalkTFullJet->GetTotalJets()>0)
ac6a3f1e 2516 {
c6202663 2517 kTRho=MedianRhokT(pTArray,RhoArray,fEMCalkTFullJet->GetTotalJets());
ac6a3f1e 2518 }
c6202663 2519 else
2520 {
2521 kTRho=0.0;
2522 }
2523 fRhoFullkT->FillRho(fEventCentrality,kTRho);
2524 fRhoFullkT->FillBSJS(fEventCentrality,kTRho,fEMCalJetThreshold,fmyAKTFullJets,fEMCalFullJet->GetJets(),fEMCalFullJet->GetTotalJets());
2525 fRhoFullkT->FillDeltaPt(fEventCentrality,kTRho,fJetR,fEMCalRCBckgFluc,1);
2526 fRhoFullkT->FillDeltaPtSignal(fEventCentrality,kTRho,fJetR,fEMCalRCBckgFlucSignal,1);
7acc3e04 2527 fRhoFullkT->FillDeltaPtNColl(fEventCentrality,kTRho,fJetR,fEMCalRCBckgFlucNColl,1);
c6202663 2528 fRhoFullkT->FillBackgroundFluctuations(fEventCentrality,kTRho,fJetR);
2529 fRhoFullkT->FillLeadingJetPtRho(fEMCalFullJet->GetLeadingPt(),kTRho);
2530 delete [] RhoArray;
2531 delete [] pTArray;
ac6a3f1e 2532}
2533
c6202663 2534void AliAnalysisTaskFullpAJets::EstimateFullRhoCMS()
ac6a3f1e 2535{
c6202663 2536 Int_t i,k;
2537 Double_t kTRho = 0.0;
2538 Double_t CMSTotalkTArea = 0.0;
2539 Double_t CMSParticleArea = 0.0;
2540 Double_t CMSCorrectionFactor = 1.0;
2541 Double_t *pTArray = new Double_t[fEMCalkTFullJet->GetTotalJets()];
2542 Double_t *RhoArray = new Double_t[fEMCalkTFullJet->GetTotalJets()];
2543
2544 k=0;
2545 if ((fEMCalPartJet->GetLeadingPt()>=fEMCalJetThreshold) && (fEMCalPartJet->GetSubLeadingPt()>=fEMCalJetThreshold))
2546 {
2547 AliEmcalJet *myJet1 =(AliEmcalJet*) fmyAKTFullJets->At(fEMCalPartJet->GetLeadingIndex());
2548 AliEmcalJet *myJet2 =(AliEmcalJet*) fmyAKTFullJets->At(fEMCalPartJet->GetSubLeadingIndex());
78246241 2549
c6202663 2550 for (i=0;i<fEMCalkTFullJet->GetTotalJets();i++)
78246241 2551 {
c6202663 2552 AliEmcalJet *myJet =(AliEmcalJet*) fmyKTFullJets->At(fEMCalkTFullJet->GetJetIndex(i));
2553
2554 CMSTotalkTArea+=myJet->Area();
2555 if (myJet->GetNumberOfTracks()>0 || myJet->GetNumberOfClusters()>0)
2556 {
2557 CMSParticleArea+=myJet->Area();
2558 }
2559 if (IsJetOverlap(myJet,myJet1,kTRUE)==kFALSE && IsJetOverlap(myJet,myJet2,kFALSE)==kTRUE)
2560 {
2561 pTArray[k]=myJet->Pt();
2562 RhoArray[k]=myJet->Pt()/myJet->Area();
2563 k++;
2564 }
78246241 2565 }
c6202663 2566 if (k>0)
ac6a3f1e 2567 {
c6202663 2568 kTRho=MedianRhokT(pTArray,RhoArray,k);
ac6a3f1e 2569 }
2570 else
2571 {
c6202663 2572 kTRho=0.0;
2573 }
2574 }
2575 else if (fEMCalJet->GetLeadingPt()>=fEMCalJetThreshold)
2576 {
2577 AliEmcalJet *myJet1 =(AliEmcalJet*) fmyAKTFullJets->At(fEMCalJet->GetLeadingIndex());
2578
2579 for (i=0;i<fEMCalkTFullJet->GetTotalJets();i++)
2580 {
2581 AliEmcalJet *myJet =(AliEmcalJet*) fmyKTFullJets->At(fEMCalkTFullJet->GetJetIndex(i));
2582
2583 CMSTotalkTArea+=myJet->Area();
2584 if (myJet->GetNumberOfTracks()>0 || myJet->GetNumberOfClusters()>0)
ac6a3f1e 2585 {
c6202663 2586 CMSParticleArea+=myJet->Area();
ac6a3f1e 2587 }
c6202663 2588 if (IsJetOverlap(myJet,myJet1,kTRUE)==kFALSE)
ac6a3f1e 2589 {
c6202663 2590 pTArray[k]=myJet->Pt();
2591 RhoArray[k]=myJet->Pt()/myJet->Area();
2592 k++;
ac6a3f1e 2593 }
ac6a3f1e 2594 }
c6202663 2595 if (k>0)
ac6a3f1e 2596 {
c6202663 2597 kTRho=MedianRhokT(pTArray,RhoArray,k);
2598 }
2599 else
2600 {
2601 kTRho=0.0;
ac6a3f1e 2602 }
2603 }
c6202663 2604 else
ac6a3f1e 2605 {
c6202663 2606 for (i=0;i<fEMCalkTFullJet->GetTotalJets();i++)
ac6a3f1e 2607 {
c6202663 2608 AliEmcalJet *myJet =(AliEmcalJet*) fmyKTFullJets->At(fEMCalkTFullJet->GetJetIndex(i));
2609
2610 CMSTotalkTArea+=myJet->Area();
2611 if (myJet->GetNumberOfTracks()>0 || myJet->GetNumberOfClusters()>0)
ac6a3f1e 2612 {
c6202663 2613 CMSParticleArea+=myJet->Area();
ac6a3f1e 2614 }
c6202663 2615 pTArray[k]=myJet->Pt();
2616 RhoArray[k]=myJet->Pt()/myJet->Area();
2617 k++;
2618 }
2619 if (k>0)
2620 {
2621 kTRho=MedianRhokT(pTArray,RhoArray,k);
2622 }
2623 else
2624 {
2625 kTRho=0.0;
ac6a3f1e 2626 }
2627 }
c6202663 2628 // Scale CMS Rho by Correction factor
2629 if (CMSTotalkTArea==0.0)
78246241 2630 {
c6202663 2631 CMSCorrectionFactor = 1.0;
78246241 2632 }
c6202663 2633 else
78246241 2634 {
c6202663 2635 //CMSCorrectionFactor = CMSTrackArea/CMSTotalkTArea;
2636 CMSCorrectionFactor = CMSParticleArea/((fEMCalPhiTotal-2*fJetR)*(fEMCalEtaTotal-2*fJetR)); // The total physical area should be reduced by the eta & phi cuts due to looping over only fully contained kT jets within the EMCal
78246241 2637 }
c6202663 2638 kTRho*=CMSCorrectionFactor;
78246241 2639
c6202663 2640 fRhoFullCMS->FillRho(fEventCentrality,kTRho);
2641 fRhoFullCMS->FillBSJS(fEventCentrality,kTRho,fEMCalJetThreshold,fmyAKTFullJets,fEMCalFullJet->GetJets(),fEMCalFullJet->GetTotalJets());
2642 fRhoFullCMS->FillDeltaPt(fEventCentrality,kTRho,fJetR,fEMCalRCBckgFluc,1);
2643 fRhoFullCMS->FillDeltaPtSignal(fEventCentrality,kTRho,fJetR,fEMCalRCBckgFlucSignal,1);
7acc3e04 2644 fRhoFullCMS->FillDeltaPtNColl(fEventCentrality,kTRho,fJetR,fEMCalRCBckgFlucNColl,1);
c6202663 2645 fRhoFullCMS->FillBackgroundFluctuations(fEventCentrality,kTRho,fJetR);
2646 fRhoFullCMS->FillLeadingJetPtRho(fEMCalFullJet->GetLeadingPt(),kTRho);
2647 delete [] RhoArray;
2648 delete [] pTArray;
78246241 2649}
62620fff 2650
c6202663 2651void AliAnalysisTaskFullpAJets::DeleteJetData(Bool_t EMCalOn)
ac6a3f1e 2652{
c6202663 2653 delete fmyTracks;
2654 delete fTPCJet;
2655 delete fTPCFullJet;
2656 delete fTPCOnlyJet;
2657 delete fTPCkTFullJet;
ac6a3f1e 2658 if (EMCalOn==kTRUE)
2659 {
c6202663 2660 delete fmyClusters;
2661 delete fEMCalJet;
2662 delete fEMCalFullJet;
2663 delete fEMCalPartJet;
2664 delete fEMCalkTFullJet;
ac6a3f1e 2665 }
2666}
2667
2668/////////////////////////////////////////////////////////////////////////////////////////
2669///////////////// User Defined Functions ///////////////////////////////////////
2670/////////////////////////////////////////////////////////////////////////////////////////
2671
2672Bool_t AliAnalysisTaskFullpAJets::IsDiJetEvent()
2673{
78246241 2674 // Determine if event contains a di-jet within the detector. Uses charged jets.
2675 // Requires the delta phi of the jets to be 180 +/- 15 degrees.
2676 // Requires both jets to be outside of the EMCal
2677 // Requires both jets to be signal jets
2678
ac6a3f1e 2679 const Double_t dijet_delta_phi=(180/360.)*2*TMath::Pi();
2680 const Double_t dijet_phi_acceptance=0.5*(30/360.)*2*TMath::Pi(); //Input the total acceptance within the paraenthesis to be +/- dijet_phi_acceptance
78246241 2681 Double_t dummy_phi=0.0;
c6202663 2682
2683 if (fTPCOnlyJet->GetTotalSignalJets()>1)
2684 {
2685 AliEmcalJet *myhJet =(AliEmcalJet*) fmyAKTChargedJets->At(fTPCOnlyJet->GetLeadingIndex());
2686 AliEmcalJet *myJet =(AliEmcalJet*) fmyAKTChargedJets->At(fTPCOnlyJet->GetSubLeadingIndex());
2687 dummy_phi=TMath::Min(TMath::Abs(myhJet->Phi()-myJet->Phi()),2*TMath::Pi()-TMath::Abs(myhJet->Phi()-myJet->Phi()));
2688 if (dummy_phi>(dijet_delta_phi-dijet_phi_acceptance))
2689 {
2690 fnDiJetEvents++;
2691 return kTRUE;
ac6a3f1e 2692 }
2693 }
2694 return kFALSE;
2695}
2696
2697Bool_t AliAnalysisTaskFullpAJets::InsideRect(Double_t phi,Double_t phi_min,Double_t phi_max,Double_t eta,Double_t eta_min,Double_t eta_max)
2698{
2699 if (phi>phi_min && phi<phi_max)
2700 {
2701 if (eta>eta_min && eta<eta_max)
2702 {
2703 return kTRUE;
2704 }
2705 }
2706 return kFALSE;
2707}
2708
2709Bool_t AliAnalysisTaskFullpAJets::IsInEMCal(Double_t phi,Double_t eta)
2710{
2711 return InsideRect(phi,fEMCalPhiMin,fEMCalPhiMax,eta,fEMCalEtaMin,fEMCalEtaMax);
2712}
2713
2714Bool_t AliAnalysisTaskFullpAJets::IsInEMCalFull(Double_t r,Double_t phi,Double_t eta)
2715{
2716 return InsideRect(phi,fEMCalPhiMin+r,fEMCalPhiMax-r,eta,fEMCalEtaMin+r,fEMCalEtaMax-r);
2717}
2718
2719Bool_t AliAnalysisTaskFullpAJets::IsInEMCalPart(Double_t r,Double_t phi,Double_t eta)
2720{
2721 return InsideRect(phi,fEMCalPhiMin-r,fEMCalPhiMax+r,eta,fEMCalEtaMin-r,fEMCalEtaMax+r);
2722}
2723
2724Bool_t AliAnalysisTaskFullpAJets::IsInTPCFull(Double_t r,Double_t phi,Double_t eta)
2725{
2726 Bool_t in_EMCal= InsideRect(phi,fEMCalPhiMin-r,fEMCalPhiMax+r,eta,fEMCalEtaMin-r,fEMCalEtaMax+r);
2727 Bool_t in_TPC= InsideRect(phi,fTPCPhiMin,fTPCPhiMax,eta,fTPCEtaMin+r,fTPCEtaMax-r);
2728
2729 if (in_EMCal==kFALSE && in_TPC==kTRUE)
2730 {
2731 return kTRUE;
2732 }
2733 return kFALSE;
2734}
2735
2736Bool_t AliAnalysisTaskFullpAJets::IsInTPC(Double_t r,Double_t phi,Double_t eta,Bool_t Complete)
2737{
2738 if (Complete==kTRUE)
2739 {
2740 return InsideRect(phi,fTPCPhiMin,fTPCPhiMax,eta,fTPCEtaMin+r,fTPCEtaMax-r);
2741 }
2742 return InsideRect(phi,fTPCPhiMin,fTPCPhiMax,eta,fTPCEtaMin,fTPCEtaMax);
2743}
2744
c6202663 2745Bool_t AliAnalysisTaskFullpAJets::IsJetOverlap(AliEmcalJet *jet1,AliEmcalJet *jet2,Bool_t EMCalOn)
2746{
2747 Int_t i,j;
2748 Int_t jetTrack1=0;
2749 Int_t jetTrack2=0;
2750 Int_t jetCluster1=0;
2751 Int_t jetCluster2=0;
2752
2753 for (i=0;i<jet1->GetNumberOfTracks();i++)
2754 {
2755 jetTrack1=jet1->TrackAt(i);
2756 for (j=0;j<jet2->GetNumberOfTracks();j++)
2757 {
2758 jetTrack2=jet2->TrackAt(j);
2759 if (jetTrack1 == jetTrack2)
2760 {
2761 return kTRUE;
2762 }
2763 }
2764 }
2765 if (EMCalOn == kTRUE)
2766 {
2767 for (i=0;i<jet1->GetNumberOfClusters();i++)
2768 {
2769 jetCluster1=jet1->ClusterAt(i);
2770 for (j=0;j<jet2->GetNumberOfClusters();j++)
2771 {
2772 jetCluster2=jet2->ClusterAt(j);
2773 if (jetCluster1 == jetCluster2)
2774 {
2775 return kTRUE;
2776 }
2777 }
2778 }
2779 }
2780 return kFALSE;
2781}
2782
ac6a3f1e 2783Double_t AliAnalysisTaskFullpAJets::AreaWithinTPC(Double_t r,Double_t eta)
2784{
2785 Double_t z;
2786 if (eta<(fTPCEtaMin+r))
2787 {
2788 z=eta-fTPCEtaMin;
2789 }
2790 else if(eta>(fTPCEtaMax-r))
2791 {
2792 z=fTPCEtaMax-eta;
2793 }
2794 else
2795 {
2796 z=r;
2797 }
2798 return r*r*TMath::Pi()-AreaEdge(r,z);
2799}
2800
2801Double_t AliAnalysisTaskFullpAJets::AreaWithinEMCal(Double_t r,Double_t phi,Double_t eta)
2802{
2803 Double_t x,y;
2804
2805 if (phi<(fEMCalPhiMin-r) || phi>(fEMCalPhiMax+r))
2806 {
2807 x=-r;
2808 }
2809 else if (phi<(fEMCalPhiMin+r))
2810 {
2811 x=phi-fEMCalPhiMin;
2812 }
2813 else if (phi>(fEMCalPhiMin+r) && phi<(fEMCalPhiMax-r))
2814 {
2815 x=r;
2816 }
2817 else
2818 {
2819 x=fEMCalPhiMax-phi;
2820 }
2821
2822 if (eta<(fEMCalEtaMin-r) || eta>(fEMCalEtaMax+r))
2823 {
2824 y=-r;
2825 }
2826 else if (eta<(fEMCalEtaMin+r))
2827 {
2828 y=eta-fEMCalEtaMin;
2829 }
2830 else if (eta>(fEMCalEtaMin+r) && eta<(fEMCalEtaMax-r))
2831 {
2832 y=r;
2833 }
2834 else
2835 {
2836 y=fEMCalEtaMax-eta;
2837 }
2838
2839 if (x>=0 && y>=0)
2840 {
2841 if (TMath::Sqrt(x*x+y*y)>=r)
2842 {
2843 return r*r*TMath::Pi()-AreaEdge(r,x)-AreaEdge(r,y);
2844 }
2845 return r*r*TMath::Pi()-AreaEdge(r,x)-AreaEdge(r,y)+AreaOverlap(r,x,y);
2846 }
2847 else if ((x>=r && y<0) || (y>=r && x<0))
2848 {
2849 return r*r*TMath::Pi()-AreaEdge(r,x)-AreaEdge(r,y);
2850 }
2851 else if (x>0 && x<r && y<0)
2852 {
2853 Double_t a=TMath::Sqrt(r*r-x*x);
2854 Double_t b=TMath::Sqrt(r*r-y*y);
2855 if ((x-b)>0)
2856 {
2857 return r*r*TMath::ASin(b/r)+y*b;
2858 }
2859 else
2860 {
2861 return 0.5*x*a+0.5*r*r*TMath::ASin(x/r)+0.5*y*b+x*y+0.5*r*r*TMath::ASin(b/r);
2862 }
2863 }
2864 else if (y>0 && y<r && x<0)
2865 {
2866 Double_t a=TMath::Sqrt(r*r-x*x);
2867 Double_t b=TMath::Sqrt(r*r-y*y);
2868 if ((y-a)>0)
2869 {
2870 return r*r*TMath::ASin(a/r)+x*a;
2871 }
2872 else
2873 {
2874 return 0.5*y*b+0.5*r*r*TMath::ASin(y/r)+0.5*x*a+x*y+0.5*r*r*TMath::ASin(a/r);
2875 }
2876 }
2877 else
2878 {
2879 Double_t a=TMath::Sqrt(r*r-x*x);
2880 Double_t b=TMath::Sqrt(r*r-y*y);
2881 if ((x+b)<0)
2882 {
2883 return 0;
2884 }
2885 else
2886 {
2887 return 0.5*x*a+0.5*r*r*TMath::ASin(x/r)+0.5*y*b+x*y+0.5*r*r*TMath::ASin(b/r);
2888 }
2889 }
2890}
2891
2892Double_t AliAnalysisTaskFullpAJets::AreaEdge(Double_t r,Double_t z)
2893{
2894 Double_t a=TMath::Sqrt(r*r-z*z);
2895 return r*r*TMath::ASin(a/r)-a*z;
2896}
2897
2898Double_t AliAnalysisTaskFullpAJets::AreaOverlap(Double_t r,Double_t x,Double_t y)
2899{
2900 Double_t a=TMath::Sqrt(r*r-x*x);
2901 Double_t b=TMath::Sqrt(r*r-y*y);
2902 return x*y-0.5*(x*a+y*b)+0.5*r*r*(TMath::ASin(b/r)-TMath::ASin(x/r));
2903}
2904
2905Double_t AliAnalysisTaskFullpAJets::TransverseArea(Double_t r,Double_t psi0,Double_t phi,Double_t eta)
2906{
c6202663 2907 Double_t area_left=0;
2908 Double_t area_right=0;
2909 Double_t eta_a=0;
2910 Double_t eta_b=0;
2911 Double_t eta_up=0;
2912 Double_t eta_down=0;
ac6a3f1e 2913
2914 Double_t u=eta-fEMCalEtaMin;
2915 Double_t v=fEMCalEtaMax-eta;
2916
2917 Double_t phi1=phi+u*TMath::Tan(psi0);
2918 Double_t phi2=phi-u*TMath::Tan(psi0);
2919 Double_t phi3=phi+v*TMath::Tan(psi0);
2920 Double_t phi4=phi-v*TMath::Tan(psi0);
2921
2922 //Calculate the Left side area
2923 if (phi1>=fEMCalPhiMax)
2924 {
2925 eta_a=eta-u*((fEMCalPhiMax-phi)/(phi1-phi));
2926 }
2927 if (phi2<=fEMCalPhiMin)
2928 {
2929 eta_b=eta-u*((phi-fEMCalPhiMin)/(phi-phi2));
2930 }
2931
2932 if ((phi1>=fEMCalPhiMax) && (phi2<=fEMCalPhiMin))
2933 {
2934 eta_up=TMath::Max(eta_a,eta_b);
2935 eta_down=TMath::Min(eta_a,eta_b);
2936
2937 area_left=(eta_down-fEMCalEtaMin)*fEMCalPhiTotal + 0.5*(fEMCalPhiTotal+2*(eta-eta_up)*TMath::Tan(psi0))*(eta_up-eta_down) + (eta-eta_up+r)*TMath::Tan(psi0)*(eta-eta_up-r);
2938 }
2939 else if (phi1>=fEMCalPhiMax)
2940 {
2941 area_left=0.5*(fEMCalPhiMax-phi2+2*(eta-eta_a)*TMath::Tan(psi0))*(eta_a-fEMCalEtaMin) + (eta-eta_a+r)*TMath::Tan(psi0)*(eta-eta_a-r);
2942 }
2943 else if (phi2<=fEMCalPhiMin)
2944 {
2945 area_left=0.5*(phi1-fEMCalPhiMin+2*(eta-eta_b)*TMath::Tan(psi0))*(eta_b-fEMCalEtaMin) + (eta-eta_b+r)*TMath::Tan(psi0)*(eta-eta_b-r);
2946 }
2947 else
2948 {
2949 area_left=0.5*(phi1-phi2+2*r*TMath::Tan(psi0))*(u-r);
2950 }
2951
c6202663 2952 // Calculate the Right side area
ac6a3f1e 2953 if (phi3>=fEMCalPhiMax)
2954 {
2955 eta_a=eta+v*((fEMCalPhiMax-phi)/(phi3-phi));
2956 }
2957 if (phi4<=fEMCalPhiMin)
2958 {
2959 eta_b=eta+v*((phi-fEMCalPhiMin)/(phi-phi4));
2960 }
2961
2962 if ((phi3>=fEMCalPhiMax) && (phi4<=fEMCalPhiMin))
2963 {
2964 eta_up=TMath::Max(eta_a,eta_b);
2965 eta_down=TMath::Min(eta_a,eta_b);
2966
2967 area_right=(fEMCalEtaMax-eta_up)*fEMCalPhiTotal + 0.5*(fEMCalPhiTotal+2*(eta_down-eta)*TMath::Tan(psi0))*(eta_down-eta_up) + (eta_down-eta+r)*TMath::Tan(psi0)*(eta_up-eta-r);
2968 }
2969 else if (phi3>=fEMCalPhiMax)
2970 {
2971 area_right=0.5*(fEMCalPhiMax-phi4+2*(eta_a-eta)*TMath::Tan(psi0))*(fEMCalEtaMax-eta_a) + (eta_a-eta+r)*TMath::Tan(psi0)*(eta_a-eta-r);
2972 }
2973 else if (phi4<=fEMCalPhiMin)
2974 {
2975 area_right=0.5*(phi3-fEMCalPhiMin+2*(eta_b-eta)*TMath::Tan(psi0))*(fEMCalEtaMax-eta_b) + (eta_b-eta+r)*TMath::Tan(psi0)*(eta_b-eta-r);
2976 }
2977 else
2978 {
2979 area_right=0.5*(phi3-phi4+2*r*TMath::Tan(psi0))*(v-r);
2980 }
2981 return area_left+area_right;
2982}
c6202663 2983
2984Double_t AliAnalysisTaskFullpAJets::MedianRhokT(Double_t *pTkTEntries, Double_t *RhokTEntries, Int_t nEntries)
2985{
2986 // This function is used to calculate the median Rho kT value. The procedure is:
2987 // - Order the kT cluster array from highest rho value to lowest
2988 // - Exclude highest rho kT cluster
2989 // - Return the median rho value of the remaining subset
2990
2991 // Sort Array
2992 const Double_t rho_min=-9.9999E+99;
2993 Int_t j,k;
2994 Double_t w[nEntries]; // Used for sorting
2995 Double_t smax=rho_min;
2996 Int_t sindex=-1;
2997
2998 Double_t pTmax=0.0;
2999 Int_t pTmaxID=-1;
3000
3001 for (j=0;j<nEntries;j++)
3002 {
3003 w[j]=0.0;
3004 }
3005
3006 for (j=0;j<nEntries;j++)
3007 {
3008 if (pTkTEntries[j]>pTmax)
3009 {
3010 pTmax=pTkTEntries[j];
3011 pTmaxID=j;
3012 }
3013 }
3014
3015 for (j=0;j<nEntries;j++)
3016 {
3017 for (k=0;k<nEntries;k++)
3018 {
3019 if (RhokTEntries[k]>smax)
3020 {
3021 smax=RhokTEntries[k];
3022 sindex=k;
3023 }
3024 }
3025 w[j]=smax;
3026 RhokTEntries[sindex]=rho_min;
3027 smax=rho_min;
3028 sindex=-1;
3029 }
3030 return w[nEntries/2];
3031}
3032
3033
3034// AlipAJetData Class Member Defs
3035// Constructors
d812e269 3036AliAnalysisTaskFullpAJets::AlipAJetData::AlipAJetData() :
3037
3038 fName(0),
3039 fIsJetsFull(0),
3040 fnTotal(0),
3041 fnJets(0),
3042 fnJetsSC(0),
3043 fJetR(0),
3044 fSignalPt(0),
3045 fAreaCutFrac(0.6),
3e43a01f 3046 fNEF(1.0),
3047 fSignalTrackBias(0),
d812e269 3048 fPtMaxIndex(0),
3049 fPtMax(0),
3050 fPtSubLeadingIndex(0),
3051 fPtSubLeading(0),
3052 fJetsIndex(0),
3053 fJetsSCIndex(0),
20f2d13a 3054 fIsJetInArray(0),
3055 fJetMaxChargedPt(0)
c6202663 3056{
3057 fnTotal=0;
3058 // Dummy constructor ALWAYS needed for I/O.
3059}
3060
d812e269 3061AliAnalysisTaskFullpAJets::AlipAJetData::AlipAJetData(const char *name, Bool_t isFull, Int_t nEntries) :
3062
3063 fName(0),
3064 fIsJetsFull(0),
3065 fnTotal(0),
3066 fnJets(0),
3067 fnJetsSC(0),
3068 fJetR(0),
3069 fSignalPt(0),
3070 fAreaCutFrac(0.6),
3e43a01f 3071 fNEF(1.0),
3072 fSignalTrackBias(0),
d812e269 3073 fPtMaxIndex(0),
3074 fPtMax(0),
3075 fPtSubLeadingIndex(0),
3076 fPtSubLeading(0),
3077 fJetsIndex(0),
3078 fJetsSCIndex(0),
20f2d13a 3079 fIsJetInArray(0),
3080 fJetMaxChargedPt(0)
c6202663 3081{
3082 SetName(name);
3083 SetIsJetsFull(isFull);
3084 SetTotalEntries(nEntries);
3085 SetLeading(0,-9.99E+099);
3086 SetSubLeading(0,-9.99E+099);
3087 SetSignalCut(0);
3088 SetAreaCutFraction(0.6);
d812e269 3089 SetJetR(fJetR);
3e43a01f 3090 SetSignalTrackPtBias(0);
c6202663 3091}
3092
3093// Destructor
3094AliAnalysisTaskFullpAJets::AlipAJetData::~AlipAJetData()
3095{
3096 if (fnTotal!=0)
3097 {
3098 SetName("");
3099 SetIsJetsFull(kFALSE);
3100 SetTotalEntries(0);
3101 SetTotalJets(0);
3102 SetTotalSignalJets(0);
3103 SetLeading(0,0);
3104 SetSubLeading(0,0);
3105 SetSignalCut(0);
3106 SetAreaCutFraction(0);
3107 SetJetR(0);
91d0893e 3108 SetNEF(0);
3e43a01f 3109 SetSignalTrackPtBias(kFALSE);
c6202663 3110
3111 delete [] fJetsIndex;
3112 delete [] fJetsSCIndex;
3113 delete [] fIsJetInArray;
20f2d13a 3114 delete [] fJetMaxChargedPt;
c6202663 3115 }
3116}
3117
3118// User Defined Sub-Routines
3119void AliAnalysisTaskFullpAJets::AlipAJetData::InitializeJetData(TClonesArray *jetList, Int_t nEntries)
3120{
3121 Int_t i=0;
3122 Int_t k=0;
3123 Int_t l=0;
3124 Double_t AreaThreshold = fAreaCutFrac*TMath::Pi()*TMath::Power(fJetR,2);
3125
3126 // Initialize Jet Data
3127 for (i=0;i<nEntries;i++)
3128 {
3129 AliEmcalJet *myJet =(AliEmcalJet*) jetList->At(i);
3130
3131 if (fIsJetInArray[i]==kTRUE && myJet->Area()>AreaThreshold)
3132 {
3133 SetJetIndex(i,k);
3134 if (myJet->Pt()>fPtMax)
3135 {
3136 SetSubLeading(fPtMaxIndex,fPtMax);
3137 SetLeading(i,myJet->Pt());
3138 }
3139 else if (myJet->Pt()>fPtSubLeading)
3140 {
20f2d13a 3141 SetSubLeading(i,myJet->Pt());
c6202663 3142 }
91d0893e 3143 // require leading charged constituent to have a pT greater then the signal threshold & Jet NEF to be less then the Signal Jet NEF cut
20f2d13a 3144 fJetMaxChargedPt[i] = myJet->MaxTrackPt();
8daeee93 3145 if (fSignalTrackBias==kTRUE)
3146 {
3147 if (fJetMaxChargedPt[i]>=fSignalPt && myJet->NEF()<=fNEF)
3148 {
3149 SetSignalJetIndex(i,l);
3150 l++;
3151 }
3152 }
3153 else
c6202663 3154 {
8daeee93 3155 if (myJet->Pt()>=fSignalPt && myJet->NEF()<=fNEF)
3156 {
3157 SetSignalJetIndex(i,l);
3158 l++;
3159 }
c6202663 3160 }
3161 k++;
3162 }
3163 }
3164 SetTotalJets(k);
3165 SetTotalSignalJets(l);
3166}
3167
3168// Setters
3169void AliAnalysisTaskFullpAJets::AlipAJetData::SetName(const char *name)
3170{
3171 fName = name;
3172}
3173
3174void AliAnalysisTaskFullpAJets::AlipAJetData::SetIsJetsFull(Bool_t isFull)
3175{
3176 fIsJetsFull = isFull;
3177}
3178
3179void AliAnalysisTaskFullpAJets::AlipAJetData::SetTotalEntries(Int_t nEntries)
3180{
3181 fnTotal = nEntries;
3182 fJetsIndex = new Int_t[fnTotal];
3183 fJetsSCIndex = new Int_t[fnTotal];
3184 fIsJetInArray = new Bool_t[fnTotal];
20f2d13a 3185 fJetMaxChargedPt = new Double_t[fnTotal];
c6202663 3186}
3187
3188void AliAnalysisTaskFullpAJets::AlipAJetData::SetTotalJets(Int_t nJets)
3189{
3190 fnJets = nJets;
3191}
3192
3193void AliAnalysisTaskFullpAJets::AlipAJetData::SetTotalSignalJets(Int_t nSignalJets)
3194{
3195 fnJetsSC = nSignalJets;
3196}
3197
3198void AliAnalysisTaskFullpAJets::AlipAJetData::SetSignalCut(Double_t Pt)
3199{
3200 fSignalPt = Pt;
3201}
3202
3203void AliAnalysisTaskFullpAJets::AlipAJetData::SetLeading(Int_t index, Double_t Pt)
3204{
3205 fPtMaxIndex = index;
3206 fPtMax = Pt;
3207}
3208
3209void AliAnalysisTaskFullpAJets::AlipAJetData::SetSubLeading(Int_t index, Double_t Pt)
3210{
3211 fPtSubLeadingIndex = index;
3212 fPtSubLeading = Pt;
3213}
3214
3215void AliAnalysisTaskFullpAJets::AlipAJetData::SetJetIndex(Int_t index, Int_t At)
3216{
3217 fJetsIndex[At] = index;
3218}
3219
3220void AliAnalysisTaskFullpAJets::AlipAJetData::SetSignalJetIndex(Int_t index, Int_t At)
3221{
3222 fJetsSCIndex[At] = index;
3223}
3224
3225void AliAnalysisTaskFullpAJets::AlipAJetData::SetIsJetInArray(Bool_t isInArray, Int_t At)
3226{
3227 fIsJetInArray[At] = isInArray;
3228}
3229
3230void AliAnalysisTaskFullpAJets::AlipAJetData::SetAreaCutFraction(Double_t areaFraction)
3231{
3232 fAreaCutFrac = areaFraction;
3233}
3234
3235void AliAnalysisTaskFullpAJets::AlipAJetData::SetJetR(Double_t jetR)
3236{
3237 fJetR = jetR;
3238}
3239
91d0893e 3240void AliAnalysisTaskFullpAJets::AlipAJetData::SetNEF(Double_t nef)
3241{
3242 fNEF = nef;
3243}
3244
8daeee93 3245void AliAnalysisTaskFullpAJets::AlipAJetData::SetSignalTrackPtBias(Bool_t chargedBias)
3246{
3247 fSignalTrackBias = chargedBias;
3248}
3249
c6202663 3250// Getters
3251Int_t AliAnalysisTaskFullpAJets::AlipAJetData::GetTotalEntries()
3252{
3253 return fnTotal;
3254}
3255
3256Int_t AliAnalysisTaskFullpAJets::AlipAJetData::GetTotalJets()
3257{
3258 return fnJets;
3259}
3260
3261Int_t AliAnalysisTaskFullpAJets::AlipAJetData::GetTotalSignalJets()
3262{
3263 return fnJetsSC;
3264}
3265
3266Double_t AliAnalysisTaskFullpAJets::AlipAJetData::GetSignalCut()
3267{
3268 return fSignalPt;
3269}
3270
3271Int_t AliAnalysisTaskFullpAJets::AlipAJetData::GetLeadingIndex()
3272{
3273 return fPtMaxIndex;
3274}
3275
3276Double_t AliAnalysisTaskFullpAJets::AlipAJetData::GetLeadingPt()
3277{
3278 return fPtMax;
3279}
3280
3281Int_t AliAnalysisTaskFullpAJets::AlipAJetData::GetSubLeadingIndex()
3282{
3283 return fPtSubLeadingIndex;
3284}
3285
3286Double_t AliAnalysisTaskFullpAJets::AlipAJetData::GetSubLeadingPt()
3287{
3288 return fPtSubLeading;
3289}
3290
3291Int_t AliAnalysisTaskFullpAJets::AlipAJetData::GetJetIndex(Int_t At)
3292{
3293 return fJetsIndex[At];
3294}
3295
3296Int_t AliAnalysisTaskFullpAJets::AlipAJetData::GetSignalJetIndex(Int_t At)
3297{
3298 return fJetsSCIndex[At];
3299}
3300
3301Bool_t AliAnalysisTaskFullpAJets::AlipAJetData::GetIsJetInArray(Int_t At)
3302{
3303 return fIsJetInArray[At];
3304}
3305
20f2d13a 3306Double_t AliAnalysisTaskFullpAJets::AlipAJetData::GetJetMaxChargedPt(Int_t At)
3307{
3308 return fJetMaxChargedPt[At];
3309}
c6202663 3310
91d0893e 3311Double_t AliAnalysisTaskFullpAJets::AlipAJetData::GetNEF()
3312{
3313 return fNEF;
3314}
3315
c6202663 3316// AlipAJetHistos Class Member Defs
3317// Constructors
d812e269 3318AliAnalysisTaskFullpAJets::AlipAJetHistos::AlipAJetHistos() :
3319
3320 fOutput(0),
3321
3322 fh020Rho(0),
3323 fh80100Rho(0),
3324 fhRho(0),
3325 fhRhoCen(0),
3326 fh020BSPt(0),
3327 fh80100BSPt(0),
3328 fhBSPt(0),
3329 fhBSPtCen(0),
3330 fh020BSPtSignal(0),
3331 fh80100BSPtSignal(0),
3332 fhBSPtSignal(0),
3333 fhBSPtCenSignal(0),
3334 fh020DeltaPt(0),
3335 fh80100DeltaPt(0),
3336 fhDeltaPt(0),
3337 fhDeltaPtCen(0),
3338 fh020DeltaPtSignal(0),
3339 fh80100DeltaPtSignal(0),
3340 fhDeltaPtSignal(0),
3341 fhDeltaPtCenSignal(0),
7acc3e04 3342 fh020DeltaPtNColl(0),
3343 fh80100DeltaPtNColl(0),
3344 fhDeltaPtNColl(0),
3345 fhDeltaPtCenNColl(0),
d812e269 3346 fh020BckgFlucPt(0),
3347 fh80100BckgFlucPt(0),
3348 fhBckgFlucPt(0),
3349 fhBckgFlucPtCen(0),
3350
3351 fpRho(0),
3352 fpLJetRho(0),
1cd88f63 3353
e864d416 3354 fhJetPtArea(0),
2d2100d5
MV		 3355 fhJetConstituentPt(0),
3356 fhJetTracksPt(0),
3357 fhJetClustersPt(0),
3358 fhJetConstituentCounts(0),
3359 fhJetTracksCounts(0),
3360 fhJetClustersCounts(0),
1cd88f63 3361
8daeee93 3362 fNEFOutput(0),
1cd88f63
ML		 3363 fhJetNEFInfo(0),
3364 fhJetNEFSignalInfo(0),
3365 fhClusterNEFInfo(0),
3366 fhClusterNEFSignalInfo(0),
8daeee93 3367 fhClusterShapeAll(0),
36cb7ae2 3368 fhClusterPtCellAll(0),
1cd88f63 3369
d812e269 3370 fName(0),
3371 fCentralityTag(0),
3372 fCentralityBins(0),
3373 fCentralityLow(0),
3374 fCentralityUp(0),
3375 fPtBins(0),
3376 fPtLow(0),
3377 fPtUp(0),
3378 fRhoPtBins(0),
3379 fRhoPtLow(0),
3380 fRhoPtUp(0),
3381 fDeltaPtBins(0),
3382 fDeltaPtLow(0),
3383 fDeltaPtUp(0),
3384 fBckgFlucPtBins(0),
3385 fBckgFlucPtLow(0),
3386 fBckgFlucPtUp(0),
3387 fLJetPtBins(0),
3388 fLJetPtLow(0),
7acc3e04 3389 fLJetPtUp(0),
20f2d13a 3390 fRhoValue(0),
3391 fLChargedTrackPtBins(0),
3392 fLChargedTrackPtLow(0),
91d0893e 3393 fLChargedTrackPtUp(0),
8daeee93 3394 fDoNEFQAPlots(0),
d2af75be 3395 fDoNEFSignalOnly(1),
3e43a01f 3396 fSignalTrackBias(0),
91d0893e 3397 fNEFBins(0),
3398 fNEFLow(0),
8daeee93 3399 fNEFUp(0),
2d2100d5
MV		 3400 fnDimJet(0),
3401 fnDimCluster(0),
8daeee93 3402 fEMCalPhiMin(1.39626),
3403 fEMCalPhiMax(3.26377),
3404 fEMCalEtaMin(-0.7),
3405 fEMCalEtaMax(0.7)
3406
c6202663 3407{
3408 // Dummy constructor ALWAYS needed for I/O.
3409}
3410
d812e269 3411AliAnalysisTaskFullpAJets::AlipAJetHistos::AlipAJetHistos(const char *name) :
3412
3413 fOutput(0),
3414
3415 fh020Rho(0),
3416 fh80100Rho(0),
3417 fhRho(0),
3418 fhRhoCen(0),
3419 fh020BSPt(0),
3420 fh80100BSPt(0),
3421 fhBSPt(0),
3422 fhBSPtCen(0),
3423 fh020BSPtSignal(0),
3424 fh80100BSPtSignal(0),
3425 fhBSPtSignal(0),
3426 fhBSPtCenSignal(0),
3427 fh020DeltaPt(0),
3428 fh80100DeltaPt(0),
3429 fhDeltaPt(0),
3430 fhDeltaPtCen(0),
3431 fh020DeltaPtSignal(0),
3432 fh80100DeltaPtSignal(0),
3433 fhDeltaPtSignal(0),
3434 fhDeltaPtCenSignal(0),
7acc3e04 3435 fh020DeltaPtNColl(0),
3436 fh80100DeltaPtNColl(0),
3437 fhDeltaPtNColl(0),
3438 fhDeltaPtCenNColl(0),
d812e269 3439 fh020BckgFlucPt(0),
3440 fh80100BckgFlucPt(0),
3441 fhBckgFlucPt(0),
3442 fhBckgFlucPtCen(0),
3443
3444 fpRho(0),
3445 fpLJetRho(0),
1cd88f63 3446
e864d416 3447 fhJetPtArea(0),
2d2100d5
MV		 3448 fhJetConstituentPt(0),
3449 fhJetTracksPt(0),
3450 fhJetClustersPt(0),
3451 fhJetConstituentCounts(0),
3452 fhJetTracksCounts(0),
3453 fhJetClustersCounts(0),
1cd88f63 3454
8daeee93 3455 fNEFOutput(0),
1cd88f63
ML		 3456 fhJetNEFInfo(0),
3457 fhJetNEFSignalInfo(0),
3458 fhClusterNEFInfo(0),
3459 fhClusterNEFSignalInfo(0),
8daeee93 3460 fhClusterShapeAll(0),
36cb7ae2 3461 fhClusterPtCellAll(0),
1cd88f63 3462
d812e269 3463 fName(0),
3464 fCentralityTag(0),
3465 fCentralityBins(0),
3466 fCentralityLow(0),
3467 fCentralityUp(0),
3468 fPtBins(0),
3469 fPtLow(0),
3470 fPtUp(0),
3471 fRhoPtBins(0),
3472 fRhoPtLow(0),
3473 fRhoPtUp(0),
3474 fDeltaPtBins(0),
3475 fDeltaPtLow(0),
3476 fDeltaPtUp(0),
3477 fBckgFlucPtBins(0),
3478 fBckgFlucPtLow(0),
3479 fBckgFlucPtUp(0),
3480 fLJetPtBins(0),
3481 fLJetPtLow(0),
7acc3e04 3482 fLJetPtUp(0),
20f2d13a 3483 fRhoValue(0),
3484 fLChargedTrackPtBins(0),
3485 fLChargedTrackPtLow(0),
91d0893e 3486 fLChargedTrackPtUp(0),
8daeee93 3487 fDoNEFQAPlots(0),
d2af75be 3488 fDoNEFSignalOnly(1),
3e43a01f 3489 fSignalTrackBias(0),
91d0893e 3490 fNEFBins(0),
3491 fNEFLow(0),
8daeee93 3492 fNEFUp(0),
2d2100d5
MV		 3493 fnDimJet(0),
3494 fnDimCluster(0),
8daeee93 3495 fEMCalPhiMin(1.39626),
3496 fEMCalPhiMax(3.26377),
3497 fEMCalEtaMin(-0.7),
3498 fEMCalEtaMax(0.7)
91d0893e 3499
c6202663 3500{
3501 SetName(name);
3502 SetCentralityTag("V0A");
3503 SetCentralityRange(100,0,100);
3504 SetPtRange(250,-50,200);
3505 SetRhoPtRange(500,0,50);
3506 SetDeltaPtRange(200,-100,100);
3507 SetBackgroundFluctuationsPtRange(100,0,100);
3508 SetLeadingJetPtRange(200,0,200);
20f2d13a 3509 SetLeadingChargedTrackPtRange(100,0,100);
91d0893e 3510 SetNEFRange(100,0,1);
8daeee93 3511 DoNEFQAPlots(kFALSE);
c6202663 3512
3513 Init();
3514}
3515
8daeee93 3516AliAnalysisTaskFullpAJets::AlipAJetHistos::AlipAJetHistos(const char *name, const char *centag, Bool_t doNEF) :
d812e269 3517
3518 fOutput(0),
3519
3520 fh020Rho(0),
3521 fh80100Rho(0),
3522 fhRho(0),
3523 fhRhoCen(0),
3524 fh020BSPt(0),
3525 fh80100BSPt(0),
3526 fhBSPt(0),
3527 fhBSPtCen(0),
3528 fh020BSPtSignal(0),
3529 fh80100BSPtSignal(0),
3530 fhBSPtSignal(0),
3531 fhBSPtCenSignal(0),
3532 fh020DeltaPt(0),
3533 fh80100DeltaPt(0),
3534 fhDeltaPt(0),
3535 fhDeltaPtCen(0),
3536 fh020DeltaPtSignal(0),
3537 fh80100DeltaPtSignal(0),
3538 fhDeltaPtSignal(0),
3539 fhDeltaPtCenSignal(0),
7acc3e04 3540 fh020DeltaPtNColl(0),
3541 fh80100DeltaPtNColl(0),
3542 fhDeltaPtNColl(0),
3543 fhDeltaPtCenNColl(0),
d812e269 3544 fh020BckgFlucPt(0),
3545 fh80100BckgFlucPt(0),
3546 fhBckgFlucPt(0),
3547 fhBckgFlucPtCen(0),
3548
3549 fpRho(0),
3550 fpLJetRho(0),
1cd88f63 3551
e864d416 3552 fhJetPtArea(0),
2d2100d5
MV		 3553 fhJetConstituentPt(0),
3554 fhJetTracksPt(0),
3555 fhJetClustersPt(0),
3556 fhJetConstituentCounts(0),
3557 fhJetTracksCounts(0),
3558 fhJetClustersCounts(0),
1cd88f63 3559
8daeee93 3560 fNEFOutput(0),
1cd88f63
ML		 3561 fhJetNEFInfo(0),
3562 fhJetNEFSignalInfo(0),
3563 fhClusterNEFInfo(0),
3564 fhClusterNEFSignalInfo(0),
8daeee93 3565 fhClusterShapeAll(0),
36cb7ae2 3566 fhClusterPtCellAll(0),
1cd88f63 3567
d812e269 3568 fName(0),
3569 fCentralityTag(0),
3570 fCentralityBins(0),
3571 fCentralityLow(0),
3572 fCentralityUp(0),
3573 fPtBins(0),
3574 fPtLow(0),
3575 fPtUp(0),
3576 fRhoPtBins(0),
3577 fRhoPtLow(0),
3578 fRhoPtUp(0),
3579 fDeltaPtBins(0),
3580 fDeltaPtLow(0),
3581 fDeltaPtUp(0),
3582 fBckgFlucPtBins(0),
3583 fBckgFlucPtLow(0),
3584 fBckgFlucPtUp(0),
3585 fLJetPtBins(0),
3586 fLJetPtLow(0),
7acc3e04 3587 fLJetPtUp(0),
20f2d13a 3588 fRhoValue(0),
3589 fLChargedTrackPtBins(0),
3590 fLChargedTrackPtLow(0),
91d0893e 3591 fLChargedTrackPtUp(0),
8daeee93 3592 fDoNEFQAPlots(0),
d2af75be 3593 fDoNEFSignalOnly(1),
3e43a01f 3594 fSignalTrackBias(0),
91d0893e 3595 fNEFBins(0),
3596 fNEFLow(0),
8daeee93 3597 fNEFUp(0),
2d2100d5
MV		 3598 fnDimJet(0),
3599 fnDimCluster(0),
8daeee93 3600 fEMCalPhiMin(1.39626),
3601 fEMCalPhiMax(3.26377),
3602 fEMCalEtaMin(-0.7),
3603 fEMCalEtaMax(0.7)
91d0893e 3604
c6202663 3605{
3606 SetName(name);
3607 SetCentralityTag(centag);
3608 SetCentralityRange(100,0,100);
3609 SetPtRange(250,-50,200);
3610 SetRhoPtRange(500,0,50);
3611 SetDeltaPtRange(200,-100,100);
3612 SetBackgroundFluctuationsPtRange(100,0,100);
3613 SetLeadingJetPtRange(200,0,200);
20f2d13a 3614 SetLeadingChargedTrackPtRange(100,0,100);
91d0893e 3615 SetNEFRange(100,0,1);
8daeee93 3616 DoNEFQAPlots(doNEF);
3617
c6202663 3618 Init();
3619}
3620
3621// Destructor
3622AliAnalysisTaskFullpAJets::AlipAJetHistos::~AlipAJetHistos()
3623{
3624 if (fOutput)
3625 {
3626 delete fOutput;
3627 }
3628}
3629
3630void AliAnalysisTaskFullpAJets::AlipAJetHistos::Init()
3631{
1cd88f63
ML		 3632 Int_t i;
3633
8daeee93 3634 // Initialize Private Variables
1cd88f63 3635 Int_t TCBins = 100;
8daeee93 3636 fEMCalPhiMin=(80/(double)360)*2*TMath::Pi();
3637 fEMCalPhiMax=(187/(double)360)*2*TMath::Pi();
3638 fEMCalEtaMin=-0.7;
3639 fEMCalEtaMax=0.7;
3640
c6202663 3641 fOutput = new TList();
3642 fOutput->SetOwner();
3643 fOutput->SetName(fName);
3644
3645 TString RhoString="";
3646 TString PtString="";
3647 TString DeltaPtString="";
3648 TString BckgFlucPtString="";
3649 TString CentralityString;
62620fff 3650 CentralityString = Form("Centrality (%s) ",fCentralityTag);
c6202663 3651
3652 // Rho Spectral Plots
3653 RhoString = Form("%d-%d Centrality, Rho Spectrum",0,20);
1cd88f63 3654 fh020Rho = new TH1F("fh020Rho",RhoString,fRhoPtBins,fRhoPtLow,fRhoPtUp);
c6202663 3655 fh020Rho->GetXaxis()->SetTitle("p_{T}/Area (GeV/c)");
3656 fh020Rho->GetYaxis()->SetTitle("1/N_{Events} dN/d#rho");
3657 fh020Rho->Sumw2();
3658
3659 RhoString = Form("%d-%d Centrality, Rho Spectrum",80,100);
1cd88f63 3660 fh80100Rho = new TH1F("fh80100Rho",RhoString,fRhoPtBins,fRhoPtLow,fRhoPtUp);
c6202663 3661 fh80100Rho->GetXaxis()->SetTitle("p_{T}/Area (GeV/c)");
3662 fh80100Rho->GetYaxis()->SetTitle("1/N_{Events} dN/d#rho");
3663 fh80100Rho->Sumw2();
3664
3665 RhoString = Form("%d-%d Centrality, Rho Spectrum",0,100);
1cd88f63 3666 fhRho = new TH1F("fhRho",RhoString,fRhoPtBins,fRhoPtLow,fRhoPtUp);
c6202663 3667 fhRho->GetXaxis()->SetTitle("p_{T}/Area (GeV/c)");
3668 fhRho->GetYaxis()->SetTitle("1/N_{Events} dN/d#rho");
3669 fhRho->Sumw2();
3670
3671 RhoString = "Rho Spectrum vs Centrality";
1cd88f63 3672 fhRhoCen = new TH2F("fhRhoCen",RhoString,fRhoPtBins,fRhoPtLow,fRhoPtUp,fCentralityBins,fCentralityLow,fCentralityUp);
c6202663 3673 fhRhoCen->GetXaxis()->SetTitle("p_{T}/Area (GeV/c)");
08b981da 3674 fhRhoCen->GetYaxis()->SetTitle(Form("%s",CentralityString.Data()));
c6202663 3675 fhRhoCen->GetZaxis()->SetTitle("1/N_{Events} dN/d#rho");
3676 fhRhoCen->Sumw2();
3677
3678 // Background Subtracted Plots
3679 PtString = Form("%d-%d Centrality, Background Subtracted Jet Spectrum",0,20);
1cd88f63 3680 fh020BSPt = new TH1F("fh020BSPt",PtString,fPtBins,fPtLow,fPtUp);
c6202663 3681 fh020BSPt->GetXaxis()->SetTitle("p_{T} - #rhoA (GeV/c)");
85b11075 3682 fh020BSPt->GetYaxis()->SetTitle("1/N_{Events} dN/dp_{T}d#etad#varphi");
c6202663 3683 fh020BSPt->Sumw2();
3684
3685 PtString = Form("%d-%d Centrality, Background Subtracted Jet Spectrum",80,100);
1cd88f63 3686 fh80100BSPt = new TH1F("fh80100BSPt",PtString,fPtBins,fPtLow,fPtUp);
c6202663 3687 fh80100BSPt->GetXaxis()->SetTitle("p_{T} - #rhoA (GeV/c)");
85b11075 3688 fh80100BSPt->GetYaxis()->SetTitle("1/N_{Events} dN/dp_{T}d#etad#varphi");
c6202663 3689 fh80100BSPt->Sumw2();
3690
3691 PtString = Form("%d-%d Centrality, Background Subtracted Jet Spectrum",0,100);
1cd88f63 3692 fhBSPt = new TH1F("fhBSPt",PtString,fPtBins,fPtLow,fPtUp);
c6202663 3693 fhBSPt->GetXaxis()->SetTitle("p_{T} - #rhoA (GeV/c)");
85b11075 3694 fhBSPt->GetYaxis()->SetTitle("1/N_{Events} dN/dp_{T}d#etad#varphi");
c6202663 3695 fhBSPt->Sumw2();
3696
3697 PtString = "Background Subtracted Jet Spectrum vs Centrality";
1cd88f63 3698 fhBSPtCen = new TH2F("fhBSPtCen",PtString,fPtBins,fPtLow,fPtUp,fCentralityBins,fCentralityLow,fCentralityUp);
c6202663 3699 fhBSPtCen->GetXaxis()->SetTitle("p_{T} - #rhoA (GeV/c)");
08b981da 3700 fhBSPtCen->GetYaxis()->SetTitle(Form("%s",CentralityString.Data()));
85b11075 3701 fhBSPtCen->GetZaxis()->SetTitle("1/N_{Events} dN/dp_{T}d#etad#varphi");
c6202663 3702 fhBSPtCen->Sumw2();
62620fff 3703
c6202663 3704 PtString = Form("%d-%d Centrality, Background Subtracted Signal Jet Spectrum",0,20);
1cd88f63 3705 fh020BSPtSignal = new TH1F("fh020BSPtSignal",PtString,fPtBins,fPtLow,fPtUp);
c6202663 3706 fh020BSPtSignal->GetXaxis()->SetTitle("p_{T} - #rhoA (GeV/c)");
85b11075 3707 fh020BSPtSignal->GetYaxis()->SetTitle("1/N_{Events} dN/dp_{T}d#etad#varphi");
c6202663 3708 fh020BSPtSignal->Sumw2();
3709
3710 PtString = Form("%d-%d Centrality, Background Subtracted Signal Jet Spectrum",80,100);
1cd88f63 3711 fh80100BSPtSignal = new TH1F("fh80100BSPtSignal",PtString,fPtBins,fPtLow,fPtUp);
c6202663 3712 fh80100BSPtSignal->GetXaxis()->SetTitle("p_{T} - #rhoA (GeV/c)");
85b11075 3713 fh80100BSPtSignal->GetYaxis()->SetTitle("1/N_{Events} dN/dp_{T}d#etad#varphi");
c6202663 3714 fh80100BSPtSignal->Sumw2();
3715
3716 PtString = Form("%d-%d Centrality, Background Subtracted Signal Jet Spectrum",0,100);
1cd88f63 3717 fhBSPtSignal = new TH1F("fhBSPtSignal",PtString,fPtBins,fPtLow,fPtUp);
c6202663 3718 fhBSPtSignal->GetXaxis()->SetTitle("p_{T} - #rhoA (GeV/c)");
85b11075 3719 fhBSPtSignal->GetYaxis()->SetTitle("1/N_{Events} dN/dp_{T}d#etad#varphi");
c6202663 3720 fhBSPtSignal->Sumw2();
3721
3722 PtString = "Background Subtracted Signal Jet Spectrum vs Centrality";
1cd88f63 3723 fhBSPtCenSignal = new TH2F("fhBSPtCenSignal",PtString,fPtBins,fPtLow,fPtUp,fCentralityBins,fCentralityLow,fCentralityUp);
c6202663 3724 fhBSPtCenSignal->GetXaxis()->SetTitle("p_{T} - #rhoA (GeV/c)");
08b981da 3725 fhBSPtCenSignal->GetYaxis()->SetTitle(Form("%s",CentralityString.Data()));
85b11075 3726 fhBSPtCenSignal->GetZaxis()->SetTitle("1/N_{Events} dN/dp_{T}d#etad#varphi");
c6202663 3727 fhBSPtCenSignal->Sumw2();
3728
3729 // Delta Pt Plots with RC at least 2R away from Leading Signal
3730 DeltaPtString = Form("%d-%d Centrality, #deltap_{T} Spectrum",0,20);
1cd88f63 3731 fh020DeltaPt = new TH1F("fh020DeltaPt",DeltaPtString,fDeltaPtBins,fDeltaPtLow,fDeltaPtUp);
c6202663 3732 fh020DeltaPt->GetXaxis()->SetTitle("#deltap_{T} (GeV/c)");
3733 fh020DeltaPt->GetYaxis()->SetTitle("Probability Density");
3734 fh020DeltaPt->Sumw2();
3735
3736 DeltaPtString = Form("%d-%d Centrality, #deltap_{T} Spectrum",80,100);
1cd88f63 3737 fh80100DeltaPt = new TH1F("fh80100DeltaPt",DeltaPtString,fDeltaPtBins,fDeltaPtLow,fDeltaPtUp);
c6202663 3738 fh80100DeltaPt->GetXaxis()->SetTitle("#deltap_{T} (GeV/c)");
3739 fh80100DeltaPt->GetYaxis()->SetTitle("Probability Density");
3740 fh80100DeltaPt->Sumw2();
3741
3742 DeltaPtString = Form("%d-%d Centrality, #deltap_{T} Spectrum",0,100);
1cd88f63 3743 fhDeltaPt = new TH1F("fhDeltaPt",DeltaPtString,fDeltaPtBins,fDeltaPtLow,fDeltaPtUp);
c6202663 3744 fhDeltaPt->GetXaxis()->SetTitle("#deltap_{T} (GeV/c)");
3745 fhDeltaPt->GetYaxis()->SetTitle("Probability Density");
3746 fhDeltaPt->Sumw2();
3747
3748 DeltaPtString = "#deltap_{T} Spectrum vs Centrality";
1cd88f63 3749 fhDeltaPtCen = new TH2F("fhDeltaPtCen",DeltaPtString,fDeltaPtBins,fDeltaPtLow,fDeltaPtUp,fCentralityBins,fCentralityLow,fCentralityUp);
c6202663 3750 fhDeltaPtCen->GetXaxis()->SetTitle("#deltap_{T} (GeV/c)");
08b981da 3751 fhDeltaPtCen->GetYaxis()->SetTitle(Form("%s",CentralityString.Data()));
c6202663 3752 fhDeltaPtCen->GetZaxis()->SetTitle("Probability Density");
3753 fhDeltaPtCen->Sumw2();
3754
3755 // Delta Pt Plots with no spatial restrictions on RC
3756 DeltaPtString = Form("%d-%d Centrality, #deltap_{T} Spectrum",0,20);
1cd88f63 3757 fh020DeltaPtSignal = new TH1F("fh020DeltaPtSignal",DeltaPtString,fDeltaPtBins,fDeltaPtLow,fDeltaPtUp);
c6202663 3758 fh020DeltaPtSignal->GetXaxis()->SetTitle("#deltap_{T} (GeV/c)");
3759 fh020DeltaPtSignal->GetYaxis()->SetTitle("Probability Density");
3760 fh020DeltaPtSignal->Sumw2();
3761
3762 DeltaPtString = Form("%d-%d Centrality, #deltap_{T} Spectrum",80,100);
1cd88f63 3763 fh80100DeltaPtSignal = new TH1F("fh80100DeltaPtSignal",DeltaPtString,fDeltaPtBins,fDeltaPtLow,fDeltaPtUp);
c6202663 3764 fh80100DeltaPtSignal->GetXaxis()->SetTitle("#deltap_{T} (GeV/c)");
3765 fh80100DeltaPtSignal->GetYaxis()->SetTitle("Probability Density");
3766 fh80100DeltaPtSignal->Sumw2();
3767
3768 DeltaPtString = Form("%d-%d Centrality, #deltap_{T} Spectrum",0,100);
1cd88f63 3769 fhDeltaPtSignal = new TH1F("fhDeltaPtSignal",DeltaPtString,fDeltaPtBins,fDeltaPtLow,fDeltaPtUp);
c6202663 3770 fhDeltaPtSignal->GetXaxis()->SetTitle("#deltap_{T} (GeV/c)");
3771 fhDeltaPtSignal->GetYaxis()->SetTitle("Probability Density");
3772 fhDeltaPtSignal->Sumw2();
3773
3774 DeltaPtString = "#deltap_{T} Spectrum vs Centrality";
1cd88f63 3775 fhDeltaPtCenSignal = new TH2F("fhDeltaPtCenSignal",DeltaPtString,fDeltaPtBins,fDeltaPtLow,fDeltaPtUp,fCentralityBins,fCentralityLow,fCentralityUp);
c6202663 3776 fhDeltaPtCenSignal->GetXaxis()->SetTitle("#deltap_{T} (GeV/c)");
08b981da 3777 fhDeltaPtCenSignal->GetYaxis()->SetTitle(Form("%s",CentralityString.Data()));
c6202663 3778 fhDeltaPtCenSignal->GetZaxis()->SetTitle("Probability Density");
3779 fhDeltaPtCenSignal->Sumw2();
7acc3e04 3780
3781 // Delta Pt Plots with NColl restrictions on RC
3782 DeltaPtString = Form("%d-%d Centrality, #deltap_{T} Spectrum",0,20);
1cd88f63 3783 fh020DeltaPtNColl = new TH1F("fh020DeltaPtNColl",DeltaPtString,fDeltaPtBins,fDeltaPtLow,fDeltaPtUp);
7acc3e04 3784 fh020DeltaPtNColl->GetXaxis()->SetTitle("#deltap_{T} (GeV/c)");
3785 fh020DeltaPtNColl->GetYaxis()->SetTitle("Probability Density");
3786 fh020DeltaPtNColl->Sumw2();
3787
3788 DeltaPtString = Form("%d-%d Centrality, #deltap_{T} Spectrum",80,100);
1cd88f63 3789 fh80100DeltaPtNColl = new TH1F("fh80100DeltaPtNColl",DeltaPtString,fDeltaPtBins,fDeltaPtLow,fDeltaPtUp);
7acc3e04 3790 fh80100DeltaPtNColl->GetXaxis()->SetTitle("#deltap_{T} (GeV/c)");
3791 fh80100DeltaPtNColl->GetYaxis()->SetTitle("Probability Density");
3792 fh80100DeltaPtNColl->Sumw2();
c6202663 3793
7acc3e04 3794 DeltaPtString = Form("%d-%d Centrality, #deltap_{T} Spectrum",0,100);
1cd88f63 3795 fhDeltaPtNColl = new TH1F("fhDeltaPtNColl",DeltaPtString,fDeltaPtBins,fDeltaPtLow,fDeltaPtUp);
7acc3e04 3796 fhDeltaPtNColl->GetXaxis()->SetTitle("#deltap_{T} (GeV/c)");
3797 fhDeltaPtNColl->GetYaxis()->SetTitle("Probability Density");
3798 fhDeltaPtNColl->Sumw2();
3799
3800 DeltaPtString = "#deltap_{T} Spectrum vs Centrality";
1cd88f63 3801 fhDeltaPtCenNColl = new TH2F("fhDeltaPtCenNColl",DeltaPtString,fDeltaPtBins,fDeltaPtLow,fDeltaPtUp,fCentralityBins,fCentralityLow,fCentralityUp);
7acc3e04 3802 fhDeltaPtCenNColl->GetXaxis()->SetTitle("#deltap_{T} (GeV/c)");
08b981da 3803 fhDeltaPtCenNColl->GetYaxis()->SetTitle(Form("%s",CentralityString.Data()));
7acc3e04 3804 fhDeltaPtCenNColl->GetZaxis()->SetTitle("Probability Density");
3805 fhDeltaPtCenNColl->Sumw2();
3806
c6202663 3807 // Background Fluctuations Pt Plots
3808 BckgFlucPtString = Form("%d-%d Centrality, Background Fluctuation p_{T} Spectrum",0,20);
1cd88f63 3809 fh020BckgFlucPt = new TH1F("fh020BckgFlucPt",PtString,fPtBins,fPtLow,fPtUp);
c6202663 3810 fh020BckgFlucPt->GetXaxis()->SetTitle("p_{T} (GeV/c)");
85b11075 3811 fh020BckgFlucPt->GetYaxis()->SetTitle("1/N_{Events} dN/dp_{T}d#etad#varphi");
c6202663 3812 fh020BckgFlucPt->Sumw2();
3813
3814 BckgFlucPtString = Form("%d-%d Centrality, Background Fluctuation p_{T} Spectrum",80,100);
1cd88f63 3815 fh80100BckgFlucPt = new TH1F("fh80100BckgFlucPt",BckgFlucPtString,fBckgFlucPtBins,fBckgFlucPtLow,fBckgFlucPtUp);
c6202663 3816 fh80100BckgFlucPt->GetXaxis()->SetTitle("p_{T} (GeV/c)");
85b11075 3817 fh80100BckgFlucPt->GetYaxis()->SetTitle("1/N_{Events} dN/dp_{T}d#etad#varphi");
c6202663 3818 fh80100BckgFlucPt->Sumw2();
3819
3820 BckgFlucPtString = Form("%d-%d Centrality, Background Fluctuation p_{T} Spectrum",0,100);
1cd88f63 3821 fhBckgFlucPt = new TH1F("fhBckgFlucPt",BckgFlucPtString,fBckgFlucPtBins,fBckgFlucPtLow,fBckgFlucPtUp);
c6202663 3822 fhBckgFlucPt->GetXaxis()->SetTitle("p_{T} (GeV/c)");
85b11075 3823 fhBckgFlucPt->GetYaxis()->SetTitle("1/N_{Events} dN/dp_{T}d#etad#varphi");
c6202663 3824 fhBckgFlucPt->Sumw2();
3825
3826 BckgFlucPtString = "Background Fluctuation p_{T} Spectrum vs Centrality";
1cd88f63 3827 fhBckgFlucPtCen = new TH2F("fhBckgFlucPtCen",BckgFlucPtString,fBckgFlucPtBins,fBckgFlucPtLow,fBckgFlucPtUp,fCentralityBins,fCentralityLow,fCentralityUp);
c6202663 3828 fhBckgFlucPtCen->GetXaxis()->SetTitle("#p_{T} (GeV/c)");
08b981da 3829 fhBckgFlucPtCen->GetYaxis()->SetTitle(Form("%s",CentralityString.Data()));
85b11075 3830 fhBckgFlucPtCen->GetZaxis()->SetTitle("1/N_{Events} dN/dp_{T}d#etad#varphi");
c6202663 3831 fhBckgFlucPtCen->Sumw2();
3832
3833 // Background Density vs Centrality Profile
3834 RhoString = "Background Density vs Centrality";
3835 fpRho = new TProfile("fpRho",RhoString,fCentralityBins,fCentralityLow,fCentralityUp);
08b981da 3836 fpRho->GetXaxis()->SetTitle(Form("%s",CentralityString.Data()));
c6202663 3837 fpRho->GetYaxis()->SetTitle("p_{T}/Area (GeV/c)");
3838
3839 // Background Density vs Leading Jet Profile
3840 fpLJetRho = new TProfile("fpLJetRho","#rho vs Leading Jet p_{T}",fLJetPtBins,fLJetPtLow,fLJetPtUp);
3841 fpLJetRho->GetXaxis()->SetTitle("Leading Jet p_{T}");
3842 fpLJetRho->GetYaxis()->SetTitle("p_{T}/Area (GeV/c)");
8daeee93 3843
e864d416 3844 // Jet pT vs Area
3845 Int_t JetPtAreaBins=200;
3846 Double_t JetPtAreaLow=0.0;
3847 Double_t JetPtAreaUp=2.0;
3848
1cd88f63 3849 fhJetPtArea = new TH2F("fhJetPtArea","Jet Area Distribution",fPtBins,fPtLow,fPtUp,JetPtAreaBins,JetPtAreaLow,JetPtAreaUp);
e864d416 3850 fhJetPtArea->GetXaxis()->SetTitle("p_{T} (GeV/c)");
3851 fhJetPtArea->GetYaxis()->SetTitle("A_{jet}");
3852 fhJetPtArea->GetZaxis()->SetTitle("1/N_{Events} dN/dA_{jet}dp_{T}");
3853 fhJetPtArea->Sumw2();
3854
2d2100d5 3855 // Jet pT vs Constituent pT
1cd88f63 3856 fhJetConstituentPt = new TH2F("fhJetConstituentPt","Jet constituents p_{T} distribution",fPtBins,fPtLow,fPtUp,10*fPtBins,fPtLow,fPtUp);
2d2100d5
MV		 3857 fhJetConstituentPt->GetXaxis()->SetTitle("Jet p_{T} (GeV/c)");
3858 fhJetConstituentPt->GetYaxis()->SetTitle("Constituent p_{T} (GeV/c)");
1cd88f63 3859 fhJetConstituentPt->GetZaxis()->SetTitle("1/N_{Events} dN/dp_{T,jet}dp_{T,track}");
2d2100d5
MV		 3860 fhJetConstituentPt->Sumw2();
3861
1cd88f63 3862 fhJetTracksPt = new TH2F("fhJetTracksPt","Jet constituents Tracks p_{T} distribution",fPtBins,fPtLow,fPtUp,10*fPtBins,fPtLow,fPtUp);
2d2100d5
MV		 3863 fhJetTracksPt->GetXaxis()->SetTitle("Jet p_{T} (GeV/c)");
3864 fhJetTracksPt->GetYaxis()->SetTitle("Constituent Track p_{T} (GeV/c)");
3865 fhJetTracksPt->GetZaxis()->SetTitle("1/N_{Events} dN/dp_{T,jet}dp_{T,track}");
3866 fhJetTracksPt->Sumw2();
3867
1cd88f63 3868 fhJetClustersPt = new TH2F("fhJetClustersPt","Jet constituents Clusters p_{T} distribution",fPtBins,fPtLow,fPtUp,10*fPtBins,fPtLow,fPtUp);
2d2100d5
MV		 3869 fhJetClustersPt->GetXaxis()->SetTitle("Jet p_{T} (GeV/c)");
3870 fhJetClustersPt->GetYaxis()->SetTitle("Constituent Cluster p_{T} (GeV/c)");
3871 fhJetClustersPt->GetZaxis()->SetTitle("1/N_{Events} dN/dp_{T,jet}dp_{T,cluster}");
3872 fhJetClustersPt->Sumw2();
3873
3874 // Jet pT vs Constituent Counts
1cd88f63 3875 fhJetConstituentCounts = new TH2F("fhJetConstituentCounts","Jet constituents distribution",fPtBins,fPtLow,fPtUp,TCBins,0,(Double_t)TCBins);
2d2100d5
MV		 3876 fhJetConstituentCounts->GetXaxis()->SetTitle("Jet p_{T} (GeV/c)");
3877 fhJetConstituentCounts->GetYaxis()->SetTitle("Constituent Count");
3878 fhJetConstituentCounts->GetZaxis()->SetTitle("1/N_{Events} dN/dp_{T,jet}dN_{constituent}");
3879 fhJetConstituentCounts->Sumw2();
3880
1cd88f63 3881 fhJetTracksCounts = new TH2F("fhJetTracksCounts","Jet constituents Tracks distribution",fPtBins,fPtLow,fPtUp,TCBins,0,(Double_t)TCBins);
2d2100d5
MV		 3882 fhJetTracksCounts->GetXaxis()->SetTitle("Jet p_{T} (GeV/c)");
3883 fhJetTracksCounts->GetYaxis()->SetTitle("Constituent Track Count");
3884 fhJetTracksCounts->GetZaxis()->SetTitle("1/N_{Events} dN/dp_{T,jet}dN_{track}");
3885 fhJetTracksCounts->Sumw2();
3886
1cd88f63 3887 fhJetClustersCounts = new TH2F("fhJetClustersCounts","Jet constituents Clusters distribution",fPtBins,fPtLow,fPtUp,TCBins,0,(Double_t)TCBins);
2d2100d5
MV		 3888 fhJetClustersCounts->GetXaxis()->SetTitle("Jet p_{T} (GeV/c)");
3889 fhJetClustersCounts->GetYaxis()->SetTitle("Constituent Cluster Count");
3890 fhJetClustersCounts->GetZaxis()->SetTitle("1/N_{Events} dN/dp_{T,jet}dN_{cluster}");
3891 fhJetClustersCounts->Sumw2();
3892
8daeee93 3893 // Neutral Energy Fraction Histograms & QA
3894 if (fDoNEFQAPlots==kTRUE)
3895 {
3896 fNEFOutput = new TList();
3897 fNEFOutput->SetOwner();
3898 fNEFOutput->SetName("ListNEFQAPlots");
3899
1cd88f63 3900 SetNEFJetDimensions(8); // Order: nef,Jet Pt,Eta,Phi,Centrality,Constituent mult,Charged mult, Neutral mult
2d2100d5
MV		 3901 SetNEFClusterDimensions(10); // Order: nef,Jet Pt,Eta,Phi,Centrality,F_Cross,z_leading,t_cell,deltat_cell,Cell Count
3902
3903 Int_t dimJetBins[fnDimJet];
3904 Double_t dimJetLow[fnDimJet];
3905 Double_t dimJetUp[fnDimJet];
3906
3907 Int_t dimClusterBins[fnDimCluster];
3908 Double_t dimClusterLow[fnDimCluster];
3909 Double_t dimClusterUp[fnDimCluster];
3910
3911 // Establish dimensinality bin counts and bounds
3912 // NEF
3913 dimJetBins[0] = dimClusterBins[0] = fNEFBins;
3914 dimJetLow[0] = dimClusterLow[0] = fNEFLow;
3915 dimJetUp[0] = dimClusterUp[0] = fNEFUp;
3916
3917 // Jet Pt
3918 dimJetBins[1] = dimClusterBins[1] = fPtBins;
3919 dimJetLow[1] = dimClusterLow[1] = fPtLow;
3920 dimJetUp[1] = dimClusterUp[1] = fPtUp;
3921
3922 // Eta-Phi
3923 dimJetBins[2] = dimJetBins[3] = dimClusterBins[2] = dimClusterBins[3] = TCBins;
3924 dimJetLow[2] = dimClusterLow[2] = fEMCalEtaMin;
3925 dimJetUp[2] = dimClusterUp[2] = fEMCalEtaMax;
3926 dimJetLow[3] = dimClusterLow[3] = fEMCalPhiMin;
3927 dimJetUp[3] = dimClusterUp[3] = fEMCalPhiMax;
1cd88f63
ML		 3928
3929 // Centrality
2d2100d5
MV		 3930 dimJetBins[4] = dimClusterBins[4] = fCentralityBins;
3931 dimJetLow[4] = dimClusterLow[4] = fCentralityLow;
3932 dimJetUp[4] = dimClusterUp[4] = fCentralityUp;
1cd88f63
ML		 3933
3934 // Jets Constituent Multiplicity Info {Total,Charged,Neutral}
3935 for (i=5;i<8;i++)
2d2100d5
MV		 3936 {
3937 dimJetBins[i] = TCBins;
3938 dimJetLow[i] = 0.0;
3939 dimJetUp[i] = (Double_t)TCBins;
3940 }
1cd88f63
ML		 3941
3942 // Cluster F_Cross
2d2100d5
MV		 3943 dimClusterBins[5] = TCBins;
3944 dimClusterLow[5] = 0.0;
3945 dimClusterUp[5] = 1.0;
1cd88f63 3946
2d2100d5
MV		 3947 // Cluster z_leading
3948 dimClusterBins[6] = TCBins;
3949 dimClusterLow[6] = 0.0;
3950 dimClusterUp[6] = 1.0;
3951
3952 // Cluster t_cell
3953 dimClusterBins[7] = 400;
3954 dimClusterLow[7] = -2e-07;
3955 dimClusterUp[7] = 2e-07;
3956
3957 // Cluster delta t_cell
3958 dimClusterBins[8] = 100;
3959 dimClusterLow[8] = 0.0;
3960 dimClusterUp[8] = 1e-07;
1cd88f63 3961
2d2100d5
MV		 3962 // Cluster Cell Count
3963 dimClusterBins[9] = TCBins;
3964 dimClusterLow[9] = 0.0;
3965 dimClusterUp[9] = 100.0;
3966
1cd88f63 3967 fhJetNEFSignalInfo = new THnSparseF("fhJetNEFSignalInfo","Signal Jet NEF Information Histogram",fnDimJet,dimJetBins,dimJetLow,dimJetUp);
2d2100d5
MV		 3968 fhJetNEFSignalInfo->Sumw2();
3969
1cd88f63 3970 fhClusterNEFSignalInfo = new THnSparseF("fhClusterNEFSignalInfo","Signal Jet NEF Cluster Information Histogram",fnDimCluster,dimClusterBins,dimClusterLow,dimClusterUp);
2d2100d5 3971 fhClusterNEFSignalInfo->Sumw2();
1cd88f63 3972
8daeee93 3973 // Cluster Shape QA
1cd88f63 3974 fhClusterShapeAll = new TH1F("fhClusterShapeAll","Cluster Shape of all CaloClustersCorr",10*TCBins,0,10*TCBins);
8daeee93 3975 fhClusterShapeAll->GetXaxis()->SetTitle("Cells");
3976 fhClusterShapeAll->GetYaxis()->SetTitle("1/N_{Events} dN/dCells");
3977 fhClusterShapeAll->Sumw2();
36cb7ae2 3978
1cd88f63 3979 fhClusterPtCellAll = new TH2F("fhClusterPtCellAll","Cluster p_{T} vs Cluster Shape of all CaloClustersCorr",fPtBins,fPtLow,fPtUp,10*TCBins,0,10*TCBins);
36cb7ae2 3980 fhClusterPtCellAll->GetXaxis()->SetTitle("p_{T} (GeV/c)");
3981 fhClusterPtCellAll->GetYaxis()->SetTitle("Cells");
3982 fhClusterPtCellAll->GetZaxis()->SetTitle("1/N_{Events} dN/dp_{T}dCells");
3983 fhClusterPtCellAll->Sumw2();
1cd88f63 3984
d2af75be
CY		 3985 if (fDoNEFSignalOnly == kTRUE)
3986 {
3987 fhJetNEFInfo = new THnSparseF("fhJetNEFInfo","Jet NEF Information Histogram",fnDimJet,dimJetBins,dimJetLow,dimJetUp);
3988 fhJetNEFInfo->Sumw2();
3989
3990 fhClusterNEFInfo = new THnSparseF("fhClusterNEFInfo","Jet NEF Cluster Information Histogram",fnDimCluster,dimClusterBins,dimClusterLow,dimClusterUp);
3991 fhClusterNEFInfo->Sumw2();
3992
3993 fNEFOutput->Add(fhJetNEFInfo);
3994 fNEFOutput->Add(fhClusterNEFInfo);
3995
3996 }
2d2100d5 3997 fNEFOutput->Add(fhJetNEFSignalInfo);
2d2100d5 3998 fNEFOutput->Add(fhClusterNEFSignalInfo);
8daeee93 3999 fNEFOutput->Add(fhClusterShapeAll);
36cb7ae2 4000 fNEFOutput->Add(fhClusterPtCellAll);
1cd88f63 4001 fOutput->Add(fNEFOutput);
8daeee93 4002 }
c6202663 4003
4004 // Add Histos & Profiles to List
4005 fOutput->Add(fh020Rho);
4006 fOutput->Add(fh80100Rho);
4007 fOutput->Add(fhRho);
4008 fOutput->Add(fhRhoCen);
4009 fOutput->Add(fh020BSPt);
4010 fOutput->Add(fh80100BSPt);
4011 fOutput->Add(fhBSPt);
4012 fOutput->Add(fhBSPtCen);
4013 fOutput->Add(fh020BSPtSignal);
4014 fOutput->Add(fh80100BSPtSignal);
4015 fOutput->Add(fhBSPtSignal);
4016 fOutput->Add(fhBSPtCenSignal);
4017 fOutput->Add(fh020DeltaPt);
4018 fOutput->Add(fh80100DeltaPt);
4019 fOutput->Add(fhDeltaPt);
4020 fOutput->Add(fhDeltaPtCen);
4021 fOutput->Add(fh020DeltaPtSignal);
4022 fOutput->Add(fh80100DeltaPtSignal);
4023 fOutput->Add(fhDeltaPtSignal);
4024 fOutput->Add(fhDeltaPtCenSignal);
7acc3e04 4025 fOutput->Add(fh020DeltaPtNColl);
4026 fOutput->Add(fh80100DeltaPtNColl);
4027 fOutput->Add(fhDeltaPtNColl);
4028 fOutput->Add(fhDeltaPtCenNColl);
c6202663 4029 fOutput->Add(fh020BckgFlucPt);
4030 fOutput->Add(fh80100BckgFlucPt);
4031 fOutput->Add(fhBckgFlucPt);
4032 fOutput->Add(fhBckgFlucPtCen);
4033 fOutput->Add(fpRho);
4034 fOutput->Add(fpLJetRho);
e864d416 4035 fOutput->Add(fhJetPtArea);
2d2100d5
MV		 4036 fOutput->Add(fhJetConstituentPt);
4037 fOutput->Add(fhJetTracksPt);
4038 fOutput->Add(fhJetClustersPt);
4039 fOutput->Add(fhJetConstituentCounts);
4040 fOutput->Add(fhJetTracksCounts);
4041 fOutput->Add(fhJetClustersCounts);
c6202663 4042}
4043
4044void AliAnalysisTaskFullpAJets::AlipAJetHistos::SetName(const char *name)
4045{
4046 fName = name;
4047}
4048
4049void AliAnalysisTaskFullpAJets::AlipAJetHistos::SetCentralityTag(const char *name)
4050{
4051 fCentralityTag = name;
4052}
4053
4054void AliAnalysisTaskFullpAJets::AlipAJetHistos::SetCentralityRange(Int_t bins, Double_t low, Double_t up)
4055{
4056 fCentralityBins=bins;
4057 fCentralityLow=low;
4058 fCentralityUp=up;
4059}
4060
4061void AliAnalysisTaskFullpAJets::AlipAJetHistos::SetPtRange(Int_t bins, Double_t low, Double_t up)
4062{
4063 fPtBins=bins;
4064 fPtLow=low;
4065 fPtUp=up;
4066}
4067
4068void AliAnalysisTaskFullpAJets::AlipAJetHistos::SetRhoPtRange(Int_t bins, Double_t low, Double_t up)
4069{
4070 fRhoPtBins=bins;
4071 fRhoPtLow=low;
4072 fRhoPtUp=up;
4073}
4074
4075void AliAnalysisTaskFullpAJets::AlipAJetHistos::SetDeltaPtRange(Int_t bins, Double_t low, Double_t up)
4076{
4077 fDeltaPtBins=bins;
4078 fDeltaPtLow=low;
4079 fDeltaPtUp=up;
4080}
4081
4082void AliAnalysisTaskFullpAJets::AlipAJetHistos::SetBackgroundFluctuationsPtRange(Int_t bins, Double_t low, Double_t up)
4083{
4084 fBckgFlucPtBins=bins;
4085 fBckgFlucPtLow=low;
4086 fBckgFlucPtUp=up;
4087}
4088
4089void AliAnalysisTaskFullpAJets::AlipAJetHistos::SetLeadingJetPtRange(Int_t bins, Double_t low, Double_t up)
4090{
4091 fLJetPtBins=bins;
4092 fLJetPtLow=low;
4093 fLJetPtUp=up;
4094}
4095
20f2d13a 4096void AliAnalysisTaskFullpAJets::AlipAJetHistos::SetLeadingChargedTrackPtRange(Int_t bins, Double_t low, Double_t up)
4097{
4098 fLChargedTrackPtBins=bins;
4099 fLChargedTrackPtLow=low;
4100 fLChargedTrackPtUp=up;
4101}
4102
91d0893e 4103void AliAnalysisTaskFullpAJets::AlipAJetHistos::SetNEFRange(Int_t bins, Double_t low, Double_t up)
4104{
4105 fNEFBins=bins;
4106 fNEFLow=low;
4107 fNEFUp=up;
4108}
4109
3e43a01f 4110void AliAnalysisTaskFullpAJets::AlipAJetHistos::SetSignalTrackPtBias(Bool_t chargedBias)
4111{
4112 fSignalTrackBias = chargedBias;
4113}
4114
2d2100d5
MV		 4115void AliAnalysisTaskFullpAJets::AlipAJetHistos::SetNEFJetDimensions(Int_t n)
4116{
4117 fnDimJet = n;
4118}
4119
4120void AliAnalysisTaskFullpAJets::AlipAJetHistos::SetNEFClusterDimensions(Int_t n)
4121{
4122 fnDimCluster = n;
4123}
4124
c6202663 4125void AliAnalysisTaskFullpAJets::AlipAJetHistos::FillRho(Double_t eventCentrality, Double_t rho)
4126{
7acc3e04 4127 fRhoValue = rho;
4128
c6202663 4129 fhRho->Fill(rho);
4130 fhRhoCen->Fill(rho,eventCentrality);
4131 fpRho->Fill(eventCentrality,rho);
4132
4133 if (eventCentrality<=20)
4134 {
4135 fh020Rho->Fill(rho);
4136 }
4137 else if (eventCentrality>=80)
4138 {
4139 fh80100Rho->Fill(rho);
4140 }
4141}
4142
4143void AliAnalysisTaskFullpAJets::AlipAJetHistos::FillBSJS(Double_t eventCentrality, Double_t rho, Double_t signalCut, TClonesArray *jetList, Int_t *indexJetList, Int_t nIndexJetList)
4144{
4145 Int_t i;
4146 Double_t tempPt=0.0;
20f2d13a 4147 Double_t tempChargedHighPt=0.0;
c6202663 4148
4149 for (i=0;i<nIndexJetList;i++)
4150 {
4151 AliEmcalJet *myJet = (AliEmcalJet*) jetList->At(indexJetList[i]);
4152 tempPt=myJet->Pt()-rho*myJet->Area();
20f2d13a 4153 tempChargedHighPt = myJet->MaxTrackPt();
91d0893e 4154
c6202663 4155 fhBSPt->Fill(tempPt);
4156 fhBSPtCen->Fill(tempPt,eventCentrality);
4157 if (eventCentrality<=20)
4158 {
4159 fh020BSPt->Fill(tempPt);
4160 }
4161 else if (eventCentrality>=80)
4162 {
4163 fh80100BSPt->Fill(tempPt);
4164 }
3e43a01f 4165 if (fSignalTrackBias==kTRUE)
c6202663 4166 {
3e43a01f 4167 if (tempChargedHighPt>=signalCut)
c6202663 4168 {
3e43a01f 4169 fhBSPtSignal->Fill(tempPt);
4170 fhBSPtCenSignal->Fill(tempPt,eventCentrality);
4171 if (eventCentrality<=20)
4172 {
4173 fh020BSPtSignal->Fill(tempPt);
4174 }
4175 else if (eventCentrality>=80)
4176 {
4177 fh80100BSPtSignal->Fill(tempPt);
4178 }
c6202663 4179 }
3e43a01f 4180 }
4181 else
4182 {
4183 if (tempPt>=signalCut)
c6202663 4184 {
3e43a01f 4185 fhBSPtSignal->Fill(tempPt);
4186 fhBSPtCenSignal->Fill(tempPt,eventCentrality);
4187 if (eventCentrality<=20)
4188 {
4189 fh020BSPtSignal->Fill(tempPt);
4190 }
4191 else if (eventCentrality>=80)
4192 {
4193 fh80100BSPtSignal->Fill(tempPt);
4194 }
c6202663 4195 }
4196 }
4197 tempPt=0.0;
20f2d13a 4198 tempChargedHighPt=0.0;
c6202663 4199 }
4200}
4201
4202void AliAnalysisTaskFullpAJets::AlipAJetHistos::FillDeltaPt(Double_t eventCentrality, Double_t rho, Double_t jetRadius, Double_t *RCArray, Int_t nRC)
4203{
4204 Int_t i;
4205 Double_t tempPt=0.0;
4206
4207 for (i=0;i<nRC;i++)
4208 {
4209 tempPt=RCArray[i]-rho*TMath::Power(jetRadius,2);
4210 fhDeltaPt->Fill(tempPt);
4211 fhDeltaPtCen->Fill(tempPt,eventCentrality);
4212 if (eventCentrality<=20)
4213 {
4214 fh020DeltaPt->Fill(tempPt);
4215 }
4216 else if (eventCentrality>=80)
4217 {
4218 fh80100DeltaPt->Fill(tempPt);
4219 }
4220 tempPt=0.0;
4221 }
4222}
4223
4224void AliAnalysisTaskFullpAJets::AlipAJetHistos::FillDeltaPtSignal(Double_t eventCentrality, Double_t rho, Double_t jetRadius, Double_t *RCArray, Int_t nRC)
4225{
4226 Int_t i;
4227 Double_t tempPt=0.0;
4228
4229 for (i=0;i<nRC;i++)
4230 {
4231 tempPt=RCArray[i]-rho*TMath::Power(jetRadius,2);
4232 fhDeltaPtSignal->Fill(tempPt);
4233 fhDeltaPtCenSignal->Fill(tempPt,eventCentrality);
4234 if (eventCentrality<=20)
4235 {
4236 fh020DeltaPtSignal->Fill(tempPt);
4237 }
4238 else if (eventCentrality>=80)
4239 {
4240 fh80100DeltaPtSignal->Fill(tempPt);
4241 }
4242 tempPt=0.0;
4243 }
4244}
4245
7acc3e04 4246void AliAnalysisTaskFullpAJets::AlipAJetHistos::FillDeltaPtNColl(Double_t eventCentrality, Double_t rho, Double_t jetRadius, Double_t *RCArray, Int_t nRC)
4247{
4248 Int_t i;
4249 Double_t tempPt=0.0;
4250
4251 for (i=0;i<nRC;i++)
4252 {
4253 tempPt=RCArray[i]-rho*TMath::Power(jetRadius,2);
4254 fhDeltaPtNColl->Fill(tempPt);
4255 fhDeltaPtCenNColl->Fill(tempPt,eventCentrality);
4256 if (eventCentrality<=20)
4257 {
4258 fh020DeltaPtNColl->Fill(tempPt);
4259 }
4260 else if (eventCentrality>=80)
4261 {
4262 fh80100DeltaPtNColl->Fill(tempPt);
4263 }
4264 tempPt=0.0;
4265 }
4266}
4267
c6202663 4268void AliAnalysisTaskFullpAJets::AlipAJetHistos::FillBackgroundFluctuations(Double_t eventCentrality, Double_t rho, Double_t jetRadius)
4269{
4270 Double_t tempPt=0.0;
4271
4272 tempPt=rho*TMath::Power(jetRadius,2);
4273 fhBckgFlucPt->Fill(tempPt);
4274 fhBckgFlucPtCen->Fill(tempPt,eventCentrality);
4275 if (eventCentrality<=20)
4276 {
4277 fh020BckgFlucPt->Fill(tempPt);
4278 }
4279 else if (eventCentrality>=80)
4280 {
4281 fh80100BckgFlucPt->Fill(tempPt);
4282 }
4283}
4284
4285void AliAnalysisTaskFullpAJets::AlipAJetHistos::FillLeadingJetPtRho(Double_t jetPt, Double_t rho)
4286{
4287 fpLJetRho->Fill(jetPt,rho);
4288}
4289
8daeee93 4290void AliAnalysisTaskFullpAJets::AlipAJetHistos::DoNEFQAPlots(Bool_t doNEFAna)
4291{
4292 fDoNEFQAPlots = doNEFAna;
4293}
4294
d2af75be
CY		 4295void AliAnalysisTaskFullpAJets::AlipAJetHistos::DoNEFSignalOnly(Bool_t doNEFSignalOnly)
4296{
4297 fDoNEFSignalOnly = doNEFSignalOnly;
4298}
4299
08b981da 4300void AliAnalysisTaskFullpAJets::AlipAJetHistos::DoNEFAnalysis(Double_t nefCut, Double_t signalCut, TClonesArray *jetList, Int_t *indexJetList, Int_t nIndexJetList, TObjArray *clusterList, TClonesArray *orgClusterList, AliVEvent *event, AliEMCALGeometry *geometry, AliEMCALRecoUtils *recoUtils, AliVCaloCells *cells)
8daeee93 4301{
4302 if (fDoNEFQAPlots==kFALSE)
4303 {
4304 return;
4305 }
1cd88f63 4306
85b11075 4307 Int_t i,j,k;
2d2100d5
MV		 4308
4309 Double_t valJet[fnDimJet];
4310 Double_t valCluster[fnDimJet];
1cd88f63 4311 for (i=0;i<fnDimJet;i++)
2d2100d5
MV		 4312 {
4313 valJet[i]=0.0;
4314 valCluster[i]=0.0;
4315 }
4316
8daeee93 4317 Double_t nef=0.0;
2d2100d5 4318 Double_t jetPt=0.0;
8daeee93 4319 Double_t tempChargedHighPt=0.0;
4320 Double_t eta=0.0;
4321 Double_t phi=0.0;
4322 Int_t totalMult=0;
2d2100d5 4323 Int_t chargedMult=0;
8daeee93 4324 Int_t neutralMult=0;
08b981da 4325 Int_t iSupMod = -1, absId = -1, ieta = -1, iphi = -1;
4326 Bool_t shared = kFALSE;
36cb7ae2 4327
08b981da 4328 Double_t zLeading=0.0;
4329 Double_t eSeed=0.0;
4330 Double_t tCell=0.0;
4331 Double_t eCross=0.0;
4332 Double_t FCross=0.0;
8daeee93 4333
85b11075 4334 Double_t lowTime=9.99e99;
4335 Double_t upTime=-9.99e99;
4336 Int_t tempCellID=0;
4337 Double_t tempCellTime=0.0;
4338
cf074128 4339 Double_t event_centrality = event->GetCentrality()->GetCentralityPercentile(fCentralityTag);
2d2100d5 4340 valJet[4] = valCluster[4] = event_centrality;
cf074128 4341
8daeee93 4342 // First, do Jet QA
4343 for (i=0;i<nIndexJetList;i++)
4344 {
4345 AliEmcalJet *myJet = (AliEmcalJet*) jetList->At(indexJetList[i]);
4346 tempChargedHighPt = myJet->MaxTrackPt();
4347 nef=myJet->NEF();
1cd88f63 4348 jetPt=myJet->Pt();
8daeee93 4349 eta=myJet->Eta();
4350 phi=myJet->Phi();
4351 totalMult=myJet->GetNumberOfConstituents();
1cd88f63 4352 chargedMult = myJet->GetNumberOfTracks();
8daeee93 4353 neutralMult=myJet->GetNumberOfClusters();
08b981da 4354 zLeading=myJet->MaxClusterPt()/myJet->Pt();
8daeee93 4355
1cd88f63 4356 valJet[0] = valCluster[0] = nef;
2d2100d5
MV		 4357 valJet[1] = valCluster[1] = jetPt;
4358 valJet[2] = valCluster[2] = eta;
4359 valJet[3] = valCluster[3] = phi;
4360
4361 valJet[5] = totalMult;
4362 valJet[6] = chargedMult;
4363 valJet[7] = neutralMult;
1cd88f63 4364
2d2100d5
MV		 4365 valCluster[6] = zLeading;
4366
d2af75be
CY		 4367 // Supress filling of this histogram due to memory size of THnSparse when running over large datasets
4368 if (fDoNEFSignalOnly == kFALSE)
4369 {
4370 fhJetNEFInfo->Fill(valJet);
4371 }
2d2100d5
MV		 4372
4373 if (fSignalTrackBias==kTRUE)
4374 {
4375 if (tempChargedHighPt>=signalCut && nef<=nefCut)
4376 {
4377 fhJetNEFSignalInfo->Fill(valJet);
4378 }
4379 }
4380 else
4381 {
4382 if (jetPt>=signalCut && nef<=nefCut)
4383 {
4384 fhJetNEFSignalInfo->Fill(valJet);
4385 }
4386 }
1cd88f63 4387
08b981da 4388 for (j=0;j<neutralMult;j++)
4389 {
4390 AliVCluster* vcluster = (AliVCluster*) orgClusterList->At(myJet->ClusterAt(j));
4391 recoUtils->GetMaxEnergyCell(geometry,cells,vcluster,absId,iSupMod,ieta,iphi,shared);
36cb7ae2 4392 eSeed = cells->GetCellAmplitude(absId);
08b981da 4393 tCell = cells->GetCellTime(absId);
4394 eCross = recoUtils->GetECross(absId,tCell,cells,event->GetBunchCrossNumber());
4395 FCross = 1 - eCross/eSeed;
1cd88f63 4396
85b11075 4397 // Obtain Delta t of Cluster
4398 lowTime=9.99e99;
4399 upTime=-9.99e99;
4400 for (k=0;k<vcluster->GetNCells();k++)
4401 {
4402 tempCellID=vcluster->GetCellAbsId(k);
4403 tempCellTime=cells->GetCellTime(tempCellID);
4404 if (tempCellTime>upTime)
4405 {
4406 upTime=tempCellTime;
4407 }
4408 if (tempCellTime<lowTime)
4409 {
4410 lowTime=tempCellTime;
4411 }
4412 }
1cd88f63 4413 valCluster[5] = FCross;
2d2100d5
MV		 4414 valCluster[7] = tCell;
4415 valCluster[8] = upTime-lowTime;
4416 valCluster[9] = vcluster->GetNCells();
4417
d2af75be
CY		 4418 if (fDoNEFSignalOnly == kFALSE)
4419 {
4420 fhClusterNEFInfo->Fill(valCluster);
4421 }
4422
2d2100d5 4423 if (fSignalTrackBias==kTRUE)
1cd88f63
ML		 4424 {
4425 if (tempChargedHighPt>=signalCut && nef<=nefCut)
3e43a01f 4426 {
1cd88f63 4427 fhClusterNEFSignalInfo->Fill(valCluster);
3e43a01f 4428 }
4429 }
1cd88f63 4430 else
3e43a01f 4431 {
1cd88f63 4432 if (myJet->Pt()>=signalCut && nef<=nefCut)
3e43a01f 4433 {
1cd88f63 4434 fhClusterNEFSignalInfo->Fill(valCluster);
3e43a01f 4435 }
8daeee93 4436 }
1cd88f63 4437 tempCellID=0;
2d2100d5
MV		 4438 tempCellTime=0.0;
4439 eSeed=0.0;
4440 tCell=0.0;
4441 eCross=0.0;
4442 FCross=0.0;
4443 iSupMod=-1,absId=-1,ieta=-1,iphi=-1;
8daeee93 4444 }
1cd88f63 4445
8daeee93 4446 nef=0.0;
1cd88f63 4447 jetPt=0.0;
8daeee93 4448 tempChargedHighPt=0.0;
4449 eta=0.0;
4450 phi=0.0;
4451 totalMult=0;
1cd88f63 4452 chargedMult=0;
8daeee93 4453 neutralMult=0;
08b981da 4454 zLeading=0.0;
8daeee93 4455 }
4456
4457 // Now do Cluster QA
4458 // Finally, Cluster QA
4459 for (i=0;i<clusterList->GetEntries();i++)
4460 {
4461 AliVCluster *vcluster = (AliVCluster*) clusterList->At(i);
4462 fhClusterShapeAll->Fill(vcluster->GetNCells());
36cb7ae2 4463 fhClusterPtCellAll->Fill(vcluster->E(),vcluster->GetNCells());
8daeee93 4464 }
4465}
4466
e864d416 4467void AliAnalysisTaskFullpAJets::AlipAJetHistos::FillMiscJetStats(TClonesArray *jetList, Int_t *indexJetList, Int_t nIndexJetList, TClonesArray *trackList, TClonesArray *clusterList)
4468{
4469 Int_t i,j;
4470
4471 for (i=0;i<nIndexJetList;i++)
4472 {
4473 AliEmcalJet *myJet = (AliEmcalJet*) jetList->At(indexJetList[i]);
4474
4475 fhJetPtArea->Fill(myJet->Pt(),myJet->Area());
1cd88f63 4476 fhJetConstituentCounts->Fill(myJet->Pt(),myJet->GetNumberOfConstituents());
2d2100d5
MV		 4477 fhJetTracksCounts->Fill(myJet->Pt(),myJet->GetNumberOfTracks());
4478 fhJetClustersCounts->Fill(myJet->Pt(),myJet->GetNumberOfClusters());
e864d416 4479 for (j=0;j<myJet->GetNumberOfTracks();j++)
4480 {
4481 AliVTrack *vtrack = (AliVTrack*) myJet->TrackAt(j,trackList);
4482 fhJetConstituentPt->Fill(myJet->Pt(),vtrack->Pt());
1cd88f63 4483 fhJetTracksPt->Fill(myJet->Pt(),vtrack->Pt());
e864d416 4484 }
4485 for (j=0;j<myJet->GetNumberOfClusters();j++)
4486 {
4487 AliVCluster *vcluster = (AliVCluster*) myJet->ClusterAt(j,clusterList);
4488 fhJetConstituentPt->Fill(myJet->Pt(),vcluster->E());
1cd88f63 4489 fhJetClustersPt->Fill(myJet->Pt(),vcluster->E());
e864d416 4490 }
4491 }
4492}
4493
8daeee93 4494TList* AliAnalysisTaskFullpAJets::AlipAJetHistos::GetOutputHistos()
4495{
4496 return fOutput;
4497}
4498
7acc3e04 4499Double_t AliAnalysisTaskFullpAJets::AlipAJetHistos::GetRho()
4500{
4501 return fRhoValue;
4502}
4503