]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGJE/EMCALJetTasks/UserTasks/AliAnalysisTaskJetShapeDeriv.cxx
bug fix for 2.5% wide centrality bins
[u/mrichter/AliRoot.git] / PWGJE / EMCALJetTasks / UserTasks / AliAnalysisTaskJetShapeDeriv.cxx
CommitLineData
f4b02da3 1//
2// Do subtraction for jet shapes using derivatives arXiv:1211:2811
3//
4// Author: M.Verweij
5
6#include <TClonesArray.h>
7#include <TH1F.h>
8#include <TH2F.h>
9#include <TH3F.h>
10#include <THnSparse.h>
11#include <TF1.h>
12#include <TList.h>
13#include <TLorentzVector.h>
14#include <TProfile.h>
15#include <TChain.h>
16#include <TSystem.h>
17#include <TFile.h>
18#include <TKey.h>
19#include <TTree.h>
20
21#include "AliVCluster.h"
22#include "AliVTrack.h"
23#include "AliEmcalJet.h"
24#include "AliRhoParameter.h"
25#include "AliLog.h"
26#include "AliEmcalParticle.h"
27#include "AliMCEvent.h"
28#include "AliAODEvent.h"
29#include "AliGenPythiaEventHeader.h"
30#include "AliAODMCHeader.h"
31#include "AliAnalysisManager.h"
32#include "AliJetContainer.h"
33#include "AliParticleContainer.h"
34
35#include "AliAnalysisTaskJetShapeDeriv.h"
36
37ClassImp(AliAnalysisTaskJetShapeDeriv)
38
39//________________________________________________________________________
40AliAnalysisTaskJetShapeDeriv::AliAnalysisTaskJetShapeDeriv() :
41 AliAnalysisTaskEmcalJet("AliAnalysisTaskJetShapeDeriv", kTRUE),
42 fContainerBase(0),
43 fMinFractionShared(0),
44 fSingleTrackEmb(kFALSE),
45 fCreateTree(kFALSE),
46 fTreeJetBkg(),
47 fJet1Vec(new TLorentzVector()),
48 fJet2Vec(new TLorentzVector()),
49 fArea(0),
50 fAreaPhi(0),
51 fAreaEta(0),
52 fRho(0),
53 fRhoM(0),
54 fNConst(0),
55 fM1st(0),
56 fM2nd(0),
57 fDeriv1st(0),
58 fDeriv2nd(0),
59 fMatch(0),
60 fh2MSubMatch(0x0),
61 fh2MSubPtRawAll(0x0),
62 fh2MSubPtRawMatch(0x0),
63 fh2MSubPtTrue(0x0),
64 fh2MTruePtTrue(0x0),
65 fh2PtTrueDeltaM(0x0),
0947b103 66 fh2PtTrueDeltaMRel(0x0),
f4b02da3 67 fhnMassResponse(0x0),
68 fh2PtTrueSubFacV1(0x0),
69 fh2PtRawSubFacV1(0x0),
70 fh2PtCorrSubFacV1(0x0),
71 fh2NConstSubFacV1(0x0),
72 fh2PtTrueSubFacV2(0x0),
73 fh2PtRawSubFacV2(0x0),
74 fh2PtCorrSubFacV2(0x0),
75 fh2NConstSubFacV2(0x0)
76{
77 // Default constructor.
78
79 fh2MSubMatch = new TH2F*[fNcentBins];
80 fh2MSubPtRawAll = new TH2F*[fNcentBins];
81 fh2MSubPtRawMatch = new TH2F*[fNcentBins];
82 fh2MSubPtTrue = new TH2F*[fNcentBins];
83 fh2MTruePtTrue = new TH2F*[fNcentBins];
84 fh2PtTrueDeltaM = new TH2F*[fNcentBins];
0947b103 85 fh2PtTrueDeltaMRel = new TH2F*[fNcentBins];
f4b02da3 86 fhnMassResponse = new THnSparse*[fNcentBins];
87 fh2PtTrueSubFacV1 = new TH2F*[fNcentBins];
88 fh2PtRawSubFacV1 = new TH2F*[fNcentBins];
89 fh2PtCorrSubFacV1 = new TH2F*[fNcentBins];
90 fh2NConstSubFacV1 = new TH2F*[fNcentBins];
91 fh2PtTrueSubFacV2 = new TH2F*[fNcentBins];
92 fh2PtRawSubFacV2 = new TH2F*[fNcentBins];
93 fh2PtCorrSubFacV2 = new TH2F*[fNcentBins];
94 fh2NConstSubFacV2 = new TH2F*[fNcentBins];
95
96 for (Int_t i = 0; i < fNcentBins; i++) {
97 fh2MSubMatch[i] = 0;
98 fh2MSubPtRawAll[i] = 0;
99 fh2MSubPtRawMatch[i] = 0;
100 fh2MSubPtTrue[i] = 0;
101 fh2MTruePtTrue[i] = 0;
102 fh2PtTrueDeltaM[i] = 0;
0947b103 103 fh2PtTrueDeltaMRel[i] = 0;
f4b02da3 104 fhnMassResponse[i] = 0;
105 fh2PtTrueSubFacV1[i] = 0;
106 fh2PtRawSubFacV1[i] = 0;
107 fh2PtCorrSubFacV1[i] = 0;
108 fh2NConstSubFacV1[i] = 0;
109 fh2PtTrueSubFacV2[i] = 0;
110 fh2PtRawSubFacV2[i] = 0;
111 fh2PtCorrSubFacV2[i] = 0;
112 fh2NConstSubFacV2[i] = 0;
113 }
114
115 SetMakeGeneralHistograms(kTRUE);
116 DefineOutput(2, TTree::Class());
117}
118
119//________________________________________________________________________
120AliAnalysisTaskJetShapeDeriv::AliAnalysisTaskJetShapeDeriv(const char *name) :
121 AliAnalysisTaskEmcalJet(name, kTRUE),
122 fContainerBase(0),
123 fMinFractionShared(0),
124 fSingleTrackEmb(kFALSE),
125 fCreateTree(kFALSE),
126 fTreeJetBkg(0),
127 fJet1Vec(new TLorentzVector()),
128 fJet2Vec(new TLorentzVector()),
129 fArea(0),
130 fAreaPhi(0),
131 fAreaEta(0),
132 fRho(0),
133 fRhoM(0),
134 fNConst(0),
135 fM1st(0),
136 fM2nd(0),
137 fDeriv1st(0),
138 fDeriv2nd(0),
139 fMatch(0),
140 fh2MSubMatch(0x0),
141 fh2MSubPtRawAll(0x0),
142 fh2MSubPtRawMatch(0x0),
143 fh2MSubPtTrue(0x0),
144 fh2MTruePtTrue(0x0),
145 fh2PtTrueDeltaM(0x0),
0947b103 146 fh2PtTrueDeltaMRel(0x0),
f4b02da3 147 fhnMassResponse(0x0),
148 fh2PtTrueSubFacV1(0x0),
149 fh2PtRawSubFacV1(0x0),
150 fh2PtCorrSubFacV1(0x0),
151 fh2NConstSubFacV1(0x0),
152 fh2PtTrueSubFacV2(0x0),
153 fh2PtRawSubFacV2(0x0),
154 fh2PtCorrSubFacV2(0x0),
155 fh2NConstSubFacV2(0x0)
156{
157 // Standard constructor.
158
159 fh2MSubMatch = new TH2F*[fNcentBins];
160 fh2MSubPtRawAll = new TH2F*[fNcentBins];
161 fh2MSubPtRawMatch = new TH2F*[fNcentBins];
162 fh2MSubPtTrue = new TH2F*[fNcentBins];
163 fh2MTruePtTrue = new TH2F*[fNcentBins];
164 fh2PtTrueDeltaM = new TH2F*[fNcentBins];
0947b103 165 fh2PtTrueDeltaMRel = new TH2F*[fNcentBins];
f4b02da3 166 fhnMassResponse = new THnSparse*[fNcentBins];
167 fh2PtTrueSubFacV1 = new TH2F*[fNcentBins];
168 fh2PtRawSubFacV1 = new TH2F*[fNcentBins];
169 fh2PtCorrSubFacV1 = new TH2F*[fNcentBins];
170 fh2NConstSubFacV1 = new TH2F*[fNcentBins];
171 fh2PtTrueSubFacV2 = new TH2F*[fNcentBins];
172 fh2PtRawSubFacV2 = new TH2F*[fNcentBins];
173 fh2PtCorrSubFacV2 = new TH2F*[fNcentBins];
174 fh2NConstSubFacV2 = new TH2F*[fNcentBins];
175
176 for (Int_t i = 0; i < fNcentBins; i++) {
177 fh2MSubMatch[i] = 0;
178 fh2MSubPtRawAll[i] = 0;
179 fh2MSubPtRawMatch[i] = 0;
180 fh2MSubPtTrue[i] = 0;
181 fh2MTruePtTrue[i] = 0;
182 fh2PtTrueDeltaM[i] = 0;
0947b103 183 fh2PtTrueDeltaMRel[i] = 0;
f4b02da3 184 fhnMassResponse[i] = 0;
185 fh2PtTrueSubFacV1[i] = 0;
186 fh2PtRawSubFacV1[i] = 0;
187 fh2PtCorrSubFacV1[i] = 0;
188 fh2NConstSubFacV1[i] = 0;
189 fh2PtTrueSubFacV2[i] = 0;
190 fh2PtRawSubFacV2[i] = 0;
191 fh2PtCorrSubFacV2[i] = 0;
192 fh2NConstSubFacV2[i] = 0;
193 }
194
195 SetMakeGeneralHistograms(kTRUE);
196 DefineOutput(2, TTree::Class());
197}
198
199//________________________________________________________________________
200AliAnalysisTaskJetShapeDeriv::~AliAnalysisTaskJetShapeDeriv()
201{
202 // Destructor.
203}
204
205//________________________________________________________________________
206void AliAnalysisTaskJetShapeDeriv::UserCreateOutputObjects()
207{
208 // Create user output.
209
210 AliAnalysisTaskEmcalJet::UserCreateOutputObjects();
211
212 Bool_t oldStatus = TH1::AddDirectoryStatus();
213 TH1::AddDirectory(kFALSE);
214
215 const Int_t nBinsPt = 200;
216 const Double_t minPt = -50.;
217 const Double_t maxPt = 150.;
218
219 const Int_t nBinsM = 150;
220 const Double_t minM = -50.;
221 const Double_t maxM = 100.;
222
8d95aff1 223 const Int_t nBinsMT = 50;
224 const Double_t minMT = 0.;
225 const Double_t maxMT = 50.;
226
f4b02da3 227 const Int_t nBinsV1 = 60;
228 const Double_t minV1 = -60.;
229 const Double_t maxV1 = 0.;
230
231 const Int_t nBinsV2 = 60;
232 const Double_t minV2 = -30.;
233 const Double_t maxV2 = 0.;
234
235 //Binning for THnSparse
8d95aff1 236 const Int_t nBinsSparse0 = 5;
237 const Int_t nBins0[nBinsSparse0] = {nBinsM,nBinsM,nBinsPt,nBinsPt,nBinsMT};
238 const Double_t xmin0[nBinsSparse0] = { minM, minM, minPt, minPt, minMT};
239 const Double_t xmax0[nBinsSparse0] = { maxM, maxM, maxPt, maxPt, maxMT};
f4b02da3 240
241 TString histName = "";
242 TString histTitle = "";
243
244 for (Int_t i = 0; i < fNcentBins; i++) {
245 histName = Form("fh2MSubMatch_%d",i);
246 histTitle = Form("fh2MSubMatch_%d;#it{M}_{sub};match",i);
247 fh2MSubMatch[i] = new TH2F(histName.Data(),histTitle.Data(),nBinsM,minM,maxM,2,-0.5,1.5);
248 fOutput->Add(fh2MSubMatch[i]);
249
250 histName = Form("fh2MSubPtRawAll_%d",i);
251 histTitle = Form("fh2MSubPtRawAll_%d;#it{M}_{sub};#it{p}_{T}",i);
252 fh2MSubPtRawAll[i] = new TH2F(histName.Data(),histTitle.Data(),nBinsM,minM,maxM,nBinsPt,minPt,maxPt);
253 fOutput->Add(fh2MSubPtRawAll[i]);
254
255 histName = Form("fh2MSubPtRawMatch_%d",i);
256 histTitle = Form("fh2MSubPtRawMatch_%d;#it{M}_{sub};#it{p}_{T}",i);
257 fh2MSubPtRawMatch[i] = new TH2F(histName.Data(),histTitle.Data(),nBinsM,minM,maxM,nBinsPt,minPt,maxPt);
258 fOutput->Add(fh2MSubPtRawMatch[i]);
259
260 histName = Form("fh2MSubPtTrue_%d",i);
261 histTitle = Form("fh2MSubPtTrue_%d;#it{M}_{sub};#it{p}_{T}",i);
262 fh2MSubPtTrue[i] = new TH2F(histName.Data(),histTitle.Data(),nBinsM,minM,maxM,nBinsPt,minPt,maxPt);
263 fOutput->Add(fh2MSubPtTrue[i]);
264
265 histName = Form("fh2MTruePtTrue_%d",i);
266 histTitle = Form("fh2MTruePtTrue_%d;#it{M}_{sub};#it{p}_{T}",i);
267 fh2MTruePtTrue[i] = new TH2F(histName.Data(),histTitle.Data(),nBinsM,minM,maxM,nBinsPt,minPt,maxPt);
268 fOutput->Add(fh2MTruePtTrue[i]);
269
270 histName = Form("fh2PtTrueDeltaM_%d",i);
271 histTitle = Form("fh2PtTrueDeltaM_%d;#it{p}_{T,true};#it{M}_{sub}-#it{M}_{true}",i);
272 fh2PtTrueDeltaM[i] = new TH2F(histName.Data(),histTitle.Data(),nBinsPt,minPt,maxPt,100,-50.,50.);
273 fOutput->Add(fh2PtTrueDeltaM[i]);
274
0947b103 275 histName = Form("fh2PtTrueDeltaMRel_%d",i);
276 histTitle = Form("fh2PtTrueDeltaMRel_%d;#it{p}_{T,true};(#it{M}_{sub}-#it{M}_{true})/#it{M}_{true}",i);
277 fh2PtTrueDeltaMRel[i] = new TH2F(histName.Data(),histTitle.Data(),nBinsPt,minPt,maxPt,200,-1.,1.);
278 fOutput->Add(fh2PtTrueDeltaMRel[i]);
279
f4b02da3 280 histName = Form("fhnMassResponse_%d",i);
8d95aff1 281 histTitle = Form("fhnMassResponse_%d;#it{M}_{sub};#it{M}_{true};#it{p}_{T,sub};#it{p}_{T,true};#it{M}_{sub}^{tagged}",i);
f4b02da3 282 fhnMassResponse[i] = new THnSparseF(histName.Data(),histTitle.Data(),nBinsSparse0,nBins0,xmin0,xmax0);
283 fOutput->Add(fhnMassResponse[i]);
284
285 //derivative histograms
286 histName = Form("fh2PtTrueSubFacV1_%d",i);
287 histTitle = Form("fh2PtTrueSubFacV1_%d;#it{p}_{T,true};-(#rho+#rho_{m})V_{1}",i);
288 fh2PtTrueSubFacV1[i] = new TH2F(histName.Data(),histTitle.Data(),nBinsPt,minPt,maxPt,nBinsV1,minV1,maxV1);
289
290 histName = Form("fh2PtRawSubFacV1_%d",i);
291 histTitle = Form("fh2PtRawSubFacV1_%d;#it{p}_{T,raw};-(#rho+#rho_{m})V_{1}",i);
292 fh2PtRawSubFacV1[i] = new TH2F(histName.Data(),histTitle.Data(),nBinsPt,minPt,maxPt,nBinsV1,minV1,maxV1);
293
294 histName = Form("fh2PtCorrSubFacV1_%d",i);
295 histTitle = Form("fh2PtCorrSubFacV1_%d;#it{p}_{T,corr};-(#rho+#rho_{m})V_{1}",i);
296 fh2PtCorrSubFacV1[i] = new TH2F(histName.Data(),histTitle.Data(),nBinsPt,minPt,maxPt,nBinsV1,minV1,maxV1);
297
298 histName = Form("fh2NConstSubFacV1_%d",i);
299 histTitle = Form("fh2NConstSubFacV1_%d;#it{N}_{const};-(#rho+#rho_{m})V_{1}",i);
300 fh2NConstSubFacV1[i] = new TH2F(histName.Data(),histTitle.Data(),nBinsPt,minPt,maxPt,100,0.,200.);
301
302 histName = Form("fh2PtTrueSubFacV2_%d",i);
303 histTitle = Form("fh2PtTrueSubFacV2_%d;#it{p}_{T,true};0.5(#rho+#rho_{m})^{2}V_{2}",i);
304 fh2PtTrueSubFacV2[i] = new TH2F(histName.Data(),histTitle.Data(),nBinsPt,minPt,maxPt,nBinsV2,minV2,maxV2);
305
306 histName = Form("fh2PtRawSubFacV2_%d",i);
307 histTitle = Form("fh2PtRawSubFacV2_%d;#it{p}_{T,raw};0.5(#rho+#rho_{m})^{2}V_{2}",i);
308 fh2PtRawSubFacV2[i] = new TH2F(histName.Data(),histTitle.Data(),nBinsPt,minPt,maxPt,nBinsV2,minV2,maxV2);
309
310 histName = Form("fh2PtCorrSubFacV2_%d",i);
311 histTitle = Form("fh2PtCorrSubFacV2_%d;#it{p}_{T,corr};0.5(#rho+#rho_{m})^{2}V_{2}",i);
312 fh2PtCorrSubFacV2[i] = new TH2F(histName.Data(),histTitle.Data(),nBinsPt,minPt,maxPt,nBinsV2,minV2,maxV2);
313
314 histName = Form("fh2NConstSubFacV2_%d",i);
315 histTitle = Form("fh2NConstSubFacV2_%d;#it{N}_{const};0.5(#rho+#rho_{m})^{2}V_{2}",i);
316 fh2NConstSubFacV2[i] = new TH2F(histName.Data(),histTitle.Data(),nBinsPt,minPt,maxPt,100,0.,200.);
317
318 }
319
320 // =========== Switch on Sumw2 for all histos ===========
321 for (Int_t i=0; i<fOutput->GetEntries(); ++i) {
322 TH1 *h1 = dynamic_cast<TH1*>(fOutput->At(i));
323 if (h1){
324 h1->Sumw2();
325 continue;
326 }
327 THnSparse *hn = dynamic_cast<THnSparse*>(fOutput->At(i));
328 if(hn)hn->Sumw2();
329 }
330
331 TH1::AddDirectory(oldStatus);
332
333 // Create a tree.
334 if(fCreateTree) {
335 fTreeJetBkg = new TTree("fTreeJetBkg", "fTreeJetBkg");
336 fTreeJetBkg->Branch("fJet1Vec","TLorentzVector",&fJet1Vec);
337 fTreeJetBkg->Branch("fJet2Vec","TLorentzVector",&fJet2Vec);
338 fTreeJetBkg->Branch("fArea",&fArea,"fArea/F");
339 fTreeJetBkg->Branch("fAreaPhi",&fAreaPhi,"fAreaPhi/F");
340 fTreeJetBkg->Branch("fAreaEta",&fAreaEta,"fAreaEta/F");
341 fTreeJetBkg->Branch("fRho",&fRho,"fRho/F");
342 fTreeJetBkg->Branch("fRhoM",&fRhoM,"fRhoM/F");
343 fTreeJetBkg->Branch("fNConst",&fNConst,"fNConst/I");
344 fTreeJetBkg->Branch("fM1st",&fM1st,"fM1st/F");
345 fTreeJetBkg->Branch("fM2nd",&fM2nd,"fM2nd/F");
346 fTreeJetBkg->Branch("fDeriv1st",&fDeriv1st,"fDeriv1st/F");
347 fTreeJetBkg->Branch("fDeriv2nd",&fDeriv2nd,"fDeriv2nd/F");
348 fTreeJetBkg->Branch("fMatch",&fMatch,"fMatch/I");
349 }
350
351 PostData(1, fOutput); // Post data for ALL output slots > 0 here.
352 if(fCreateTree) PostData(2, fTreeJetBkg);
353}
354
355//________________________________________________________________________
356Bool_t AliAnalysisTaskJetShapeDeriv::Run()
357{
358 // Run analysis code here, if needed. It will be executed before FillHistograms().
359
360 return kTRUE;
361}
362
363//________________________________________________________________________
364Bool_t AliAnalysisTaskJetShapeDeriv::FillHistograms()
365{
366 // Fill histograms.
367
8d95aff1 368 AliEmcalJet* jet1 = NULL; //AA jet
369 AliEmcalJet *jet2 = NULL; //Embedded Pythia jet
370 AliEmcalJet *jet1T = NULL; //tagged AA jet
371 // AliEmcalJet *jet2T = NULL; //tagged Pythia jet
f4b02da3 372 AliJetContainer *jetCont = GetJetContainer(fContainerBase);
373 fRho = (Float_t)jetCont->GetRhoVal();
374 fRhoM = (Float_t)jetCont->GetRhoMassVal();
375 if(jetCont) {
376 jetCont->ResetCurrentID();
377 while((jet1 = jetCont->GetNextAcceptJet())) {
8d95aff1 378 jet2 = NULL;
379 jet1T = NULL;
380 // jet2T = NULL;
381 if(jet1->GetTagStatus()>0) jet1T = jet1->GetTaggedJet();
f4b02da3 382 //Fill histograms for all AA jets
383 fh2MSubPtRawAll[fCentBin]->Fill(jet1->GetSecondOrderSubtracted(),jet1->Pt()-jetCont->GetRhoVal()*jet1->Area());
384 fh2PtRawSubFacV1[fCentBin]->Fill(jet1->Pt(),-1.*(fRho+fRhoM)*jet1->GetFirstDerivative());
385 fh2PtCorrSubFacV1[fCentBin]->Fill(jet1->Pt()-fRho*jet1->Area(),-1.*(fRho+fRhoM)*jet1->GetFirstDerivative());
386 fh2NConstSubFacV1[fCentBin]->Fill(jet1->GetNumberOfTracks(),-1.*(fRho+fRhoM)*jet1->GetFirstDerivative());
387 fh2PtRawSubFacV2[fCentBin]->Fill(jet1->Pt(),0.5*(fRho+fRhoM)*(fRho+fRhoM)*jet1->GetSecondDerivative());
388 fh2PtCorrSubFacV2[fCentBin]->Fill(jet1->Pt()-fRho*jet1->Area(),0.5*(fRho+fRhoM)*(fRho+fRhoM)*jet1->GetSecondDerivative());
389 fh2NConstSubFacV2[fCentBin]->Fill(jet1->GetNumberOfTracks(),0.5*(fRho+fRhoM)*(fRho+fRhoM)*jet1->GetSecondDerivative());
390
391 Double_t fraction = 0.;
392 fMatch = 0;
393 fJet2Vec->SetPtEtaPhiM(0.,0.,0.,0.);
394 if(fSingleTrackEmb) {
395 AliVParticle *vp = GetEmbeddedConstituent(jet1);
396 if(vp) {
397 fJet2Vec->SetPxPyPzE(vp->Px(),vp->Py(),vp->Pz(),vp->E());
398 fMatch = 1;
399 }
400 } else {
401 jet2 = jet1->ClosestJet();
402 fraction = jetCont->GetFractionSharedPt(jet1);
403 fMatch = 1;
404 if(fMinFractionShared>0.) {
405 if(fraction>fMinFractionShared) {
406 fJet2Vec->SetPxPyPzE(jet2->Px(),jet2->Py(),jet2->Pz(),jet2->E());
407 fMatch = 1;
408 } else
409 fMatch = 0;
410 }
411 }
412
8d95aff1 413 // if(fMatch==1 && jet2->GetTagStatus()>0) jet2T = jet2->GetTaggedJet();
414
f4b02da3 415 //Fill histograms for matched jets
416 fh2MSubMatch[fCentBin]->Fill(jet1->GetSecondOrderSubtracted(),fMatch);
417 if(fMatch==1) {
418 fh2MSubPtRawMatch[fCentBin]->Fill(jet1->GetSecondOrderSubtracted(),jet1->Pt()-jetCont->GetRhoVal()*jet1->Area());
419 if(jet2) {
420 fh2MSubPtTrue[fCentBin]->Fill(jet1->GetSecondOrderSubtracted(),jet2->Pt());
421 fh2MTruePtTrue[fCentBin]->Fill(jet2->M(),jet2->Pt());
422 fh2PtTrueDeltaM[fCentBin]->Fill(jet2->Pt(),jet1->GetSecondOrderSubtracted()-jet2->M());
0947b103 423 if(jet2->M()>0.) fh2PtTrueDeltaMRel[fCentBin]->Fill(jet2->Pt(),(jet1->GetSecondOrderSubtracted()-jet2->M())/jet2->M());
8d95aff1 424 Double_t mJet1Tagged = -1.;
425 if(jet1T) mJet1Tagged = jet1T->M();
426 Double_t var[5] = {jet1->GetSecondOrderSubtracted(),jet2->M(),jet1->Pt()-jetCont->GetRhoVal()*jet1->Area(),jet2->Pt(),mJet1Tagged};
f4b02da3 427 fhnMassResponse[fCentBin]->Fill(var);
428 }
429 }
430
431 // if(jet2) Printf("unsubtracted: %f pt: %f true: %f pt: %f",jet1->M(),jet1->Pt()-jetCont->GetRhoVal()*jet1->Area(),jet2->M(),jet2->Pt());
432 // Printf("1st derivative: %f",jet1->GetFirstDerivative());
433 // Printf("2nd derivative: %f",jet1->GetSecondDerivative());
434 // Printf("1st order subtracted: %f",jet1->GetFirstOrderSubtracted());
435 // Printf("2nd order subtracted: %f",jet1->GetSecondOrderSubtracted());
436
437 if(fCreateTree) {
438 fJet1Vec->SetPxPyPzE(jet1->Px(),jet1->Py(),jet1->Pz(),jet1->E());
439 fArea = (Float_t)jet1->Area();
440 fAreaPhi = (Float_t)jet1->AreaPhi();
441 fAreaEta = (Float_t)jet1->AreaEta();
442 fNConst = (Int_t)jet1->GetNumberOfTracks();
443 fM1st = (Float_t)jet1->GetFirstOrderSubtracted();
444 fM2nd = (Float_t)jet1->GetSecondOrderSubtracted();
445 fDeriv1st = (Float_t)jet1->GetFirstDerivative();
446 fDeriv2nd = (Float_t)jet1->GetSecondDerivative();
447 fTreeJetBkg->Fill();
448 }
449 }
450 }
451 return kTRUE;
452}
453
454//________________________________________________________________________
455AliVParticle* AliAnalysisTaskJetShapeDeriv::GetEmbeddedConstituent(AliEmcalJet *jet) {
456
457 AliJetContainer *jetCont = GetJetContainer(fContainerBase);
458 AliVParticle *vp = 0x0;
459 AliVParticle *vpe = 0x0; //embedded particle
460 Int_t nc = 0;
461 for(Int_t i=0; i<jet->GetNumberOfTracks(); i++) {
462 vp = static_cast<AliVParticle*>(jet->TrackAt(i, jetCont->GetParticleContainer()->GetArray())); //check if fTracks is the correct track branch
463 if (vp->TestBits(TObject::kBitMask) != (Int_t)(TObject::kBitMask) ) continue;
464 if(!vpe) vpe = vp;
465 else if(vp->Pt()>vpe->Pt()) vpe = vp;
466 nc++;
467 }
468 AliDebug(11,Form("Found %d embedded particles",nc));
469 return vpe;
470}
471
472
473//________________________________________________________________________
474Bool_t AliAnalysisTaskJetShapeDeriv::RetrieveEventObjects() {
475 //
476 // retrieve event objects
477 //
478
479 if (!AliAnalysisTaskEmcalJet::RetrieveEventObjects())
480 return kFALSE;
481
482 AliJetContainer *jetCont = GetJetContainer(fContainerBase);
483 jetCont->LoadRhoMass(InputEvent());
484
485 return kTRUE;
486}
487
488//_______________________________________________________________________
489void AliAnalysisTaskJetShapeDeriv::Terminate(Option_t *)
490{
491 // Called once at the end of the analysis.
492}
493