AliPhysics  e6c8d43 (e6c8d43)
AliAnalysisTaskEmcalJetShapesMC.cxx
Go to the documentation of this file.
1 //
2 // Jet QG tagging analysis task.
3 //
4 // Author: D. Caffarri, L. Cunqueiro
5 
6 #include <TClonesArray.h>
7 #include <TH1F.h>
8 #include <TH2F.h>
9 #include <TH3F.h>
10 #include <THnSparse.h>
11 #include <TTree.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 <AliAnalysisDataSlot.h>
20 #include <AliAnalysisDataContainer.h>
21 #include "TMatrixD.h"
22 #include "TMatrixDSym.h"
23 #include "TMatrixDSymEigen.h"
24 #include "TVector3.h"
25 #include "TVector2.h"
26 #include "AliVCluster.h"
27 #include "AliVTrack.h"
28 #include "AliEmcalJet.h"
29 #include "AliRhoParameter.h"
30 #include "AliLog.h"
31 #include "AliEmcalParticle.h"
32 #include "AliMCEvent.h"
33 #include "AliGenPythiaEventHeader.h"
34 #include "AliAODMCHeader.h"
35 #include "AliMCEvent.h"
36 #include "AliAnalysisManager.h"
37 #include "AliJetContainer.h"
38 #include "AliParticleContainer.h"
39 #include "AliEmcalPythiaInfo.h"
40 #include "TRandom3.h"
41 #include "TF1.h"
42 #include "AliEmcalJetFinder.h"
43 #include "AliAODEvent.h"
45 
46 using std::cout;
47 using std::endl;
48 
50 
51 //________________________________________________________________________
54  fContainer(0),
55  fMinFractionShared(0),
56  fJetShapeType(kGenShapes),
57  fJetShapeSub(kNoSub),
58  fJetSelection(kInclusive),
59  fPtThreshold(-9999.),
60  fRMatching(0.2),
61  fJetRadius(0.4),
62  fSubjetRadius(0.2),
63  fSelectedShapes(0),
64  fSwitchKtNSub(0),
65  fSwitchMinNSub(0),
66  fSwitchAktNSub(0),
67  fSwitchSDKtNSub(0),
68  fSwitchSDMinNSub(0),
69  fAdditionalTracks(0),
70  fHardCutoff(0),
71  fminpTTrig(20.),
72  fmaxpTTrig(50.),
73  fangWindowRecoil(0.6),
74  fSemigoodCorrect(0),
75  fHolePos(0),
76  fHoleWidth(0),
77  fRandom(0),
78  fqhat(1),
79  fxlength(2),
80  fCentSelectOn(kTRUE),
81  fCentMin(0),
82  fCentMax(10),
83  fOneConstSelectOn(kFALSE),
84  fDerivSubtrOrder(0),
85  fPhiJetCorr6(0x0),
86  fPhiJetCorr7(0x0),
87  fEtaJetCorr6(0x0),
88  fEtaJetCorr7(0x0),
89  fPtJetCorr(0x0),
90  fPtJet(0x0),
91  fhpTjetpT(0x0),
92  fhPt(0x0),
93  fhPhi(0x0),
94  fHLundIterative(0x0),
95  fHLundIterativeInject(0x0),
96  fNbOfConstvspT(0x0),
97  fTreeObservableTagging(0x0),
98  fTf1Omega(0x0),
99  fTf1Kt(0x0),
100  fScaleELoss(kFALSE),
101  xfraction(1),
102  fAddMedScat(kFALSE),
103  fAddMedScatPtFrac(1),
104  fAddMedScatN(100)
105 
106 
107 {
108  for(Int_t i=0;i<11;i++){
109  fShapesVar[i]=0;}
110  SetMakeGeneralHistograms(kTRUE);
111  DefineOutput(1, TList::Class());
112  DefineOutput(2, TTree::Class());
113 }
114 
115 //________________________________________________________________________
117  AliAnalysisTaskEmcalJet(name, kTRUE),
118  fContainer(0),
119  fMinFractionShared(0),
120  fJetShapeType(kGenShapes),
121  fJetShapeSub(kNoSub),
122  fJetSelection(kInclusive),
123  fPtThreshold(-9999.),
124  fRMatching(0.2),
125  fSelectedShapes(0),
126  fSwitchKtNSub(0),
127  fSwitchMinNSub(0),
128  fSwitchAktNSub(0),
129  fSwitchSDKtNSub(0),
130  fSwitchSDMinNSub(0),
131  fAdditionalTracks(0),
132  fHardCutoff(0),
133  fminpTTrig(20.),
134  fmaxpTTrig(50.),
135  fangWindowRecoil(0.6),
136  fSemigoodCorrect(0),
137  fHolePos(0),
138  fHoleWidth(0),
139  fRandom(0),
140  fqhat(1),
141  fxlength(2),
142  fCentSelectOn(kTRUE),
143  fCentMin(0),
144  fCentMax(10),
145  fOneConstSelectOn(kFALSE),
146  fDerivSubtrOrder(0),
147  fPhiJetCorr6(0x0),
148  fPhiJetCorr7(0x0),
149  fEtaJetCorr6(0x0),
150  fEtaJetCorr7(0x0),
151  fPtJetCorr(0x0),
152  fPtJet(0x0),
153  fhpTjetpT(0x0),
154  fhPt(0x0),
155  fhPhi(0x0),
156  fHLundIterative(0x0),
157  fHLundIterativeInject(0x0),
158  fNbOfConstvspT(0x0),
159  fTreeObservableTagging(0x0),
160  fTf1Omega(0x0),
161  fTf1Kt(0x0),
162  fScaleELoss(kFALSE),
163  xfraction(1),
164  fAddMedScat(kFALSE),
165  fAddMedScatPtFrac(1),
166  fAddMedScatN(100)
167 {
168  // Standard constructor.
169 
170 
171  for(Int_t i=0;i<11;i++){
172  fShapesVar[i]=0;}
173 
175 
176  DefineOutput(1, TList::Class());
177  DefineOutput(2, TTree::Class());
178 
179 
180 }
181 
182 //________________________________________________________________________
184 {
186  delete fTreeObservableTagging;
188  }
189 
190  if(fRandom) delete fRandom;
191  if(fTf1Omega) delete fTf1Omega;
192  if(fTf1Kt) delete fTf1Kt;
193 
194 }
195 
196 //________________________________________________________________________
198 {
199  // Create user output.
201 
202  Bool_t oldStatus = TH1::AddDirectoryStatus();
203  TH1::AddDirectory(oldStatus);
204 
205  //fTreeObservableTagging = new TTree("fTreeJetShape", "fTreeJetShape");
206 
207  //TH1::AddDirectory(oldStatus);
208 
209  const char* nameoutput = GetOutputSlot(2)->GetContainer()->GetName();
210  fTreeObservableTagging = new TTree(nameoutput, nameoutput);
211 
212  const Int_t nVar = 11;
213 
214  TString *fShapesVarNames = new TString [nVar];
215 
216  fShapesVarNames[0] = "partonCode";
217  fShapesVarNames[1] = "ptJet";
218  fShapesVarNames[2] = "ptDJet";
219  fShapesVarNames[3] = "mJet";
220  fShapesVarNames[4] = "nbOfConst";
221  fShapesVarNames[5] = "angularity";
222  fShapesVarNames[6] = "nitersd";
223  fShapesVarNames[7] = "niterall";
224  fShapesVarNames[8] = "weightPythia";
225  //fShapesVarNames[6] = "Nsubjet1kt";
226  // fShapesVarNames[7] = "Nsubjet2kt";
227  // fShapesVarNames[8] = "Nsubjet1Min";
228  // fShapesVarNames[9] = "Nsubjet2Min";
229  // fShapesVarNames[10] = "DeltaRkt";
230  // fShapesVarNames[11] = "DeltaRMin";
231  fShapesVarNames[9] = "SDSymm";
232  fShapesVarNames[10] = "scaledptJet";
233  // fShapesVarNames[13] = "SDDeltaR";
234  // fShapesVarNames[14] = "SDGroomedFrac";
235  // fShapesVarNames[15] = "SDGroomedN";
236  // fShapesVarNames[16] = "SDMass";
237  // fShapesVarNames[17] = "SDSymmkt";
238  // fShapesVarNames[18] = "SDDeltaRkt";
239  // fShapesVarNames[19] = "SDGroomedFrackt";
240  // fShapesVarNames[20] = "SDGroomedNkt";
241  // fShapesVarNames[21] = "SDMasskt";
242  // fShapesVarNames[22] = "SDSymmAkt";
243  // fShapesVarNames[23] = "SDDeltaRAkt";
244  // fShapesVarNames[24] = "SDGroomedFracAkt";
245  // fShapesVarNames[25] = "SDGroomedNAkt";
246  // fShapesVarNames[26] = "SDMassAkt";
247  // fShapesVarNames[27] = "SDSymmktForm";
248  // fShapesVarNames[28] = "SDDeltaRktForm";
249  // fShapesVarNames[29] = "SDGroomedFracktForm";
250  // fShapesVarNames[30] = "SDGroomedNktForm";
251  // fShapesVarNames[31] = "SDMassktForm";
252  // fShapesVarNames[32] = "SDSymmDemo";
253  // fShapesVarNames[33] = "SDDeltaRDemo";
254  // fShapesVarNames[34] = "SDGroomedFracDemo";
255  // fShapesVarNames[35] = "SDGroomedNDemo";
256  // fShapesVarNames[36] = "SDMassDemo";
257  // fShapesVarNames[42] = "SDSymmForm";
258  // fShapesVarNames[43] = "SDDeltaRForm";
259  // fShapesVarNames[44] = "SDGroomedFracForm";
260  // fShapesVarNames[45] = "SDGroomedNForm";
261  // fShapesVarNames[46] = "SDMassForm";
262  // fShapesVarNames[47] = "weightPythia";
263  // fShapesVarNames[42] = "SDSymmNoCut";
264  // fShapesVarNames[43] = "SDDeltaRNoCut";
265  // fShapesVarNames[44] = "SDGroomedFracNoCut";
266  // fShapesVarNames[45] = "SDGroomedNNoCut";
267  // fShapesVarNames[46] = "SDMassNoCut";
268 
269 
270  //fShapesVarNames[7] = "lesub";
271  //fShapesVarNames[8] = "CoreFraction";
272  //fShapesVarNames[9] = "Nsubjet1";
273  //fShapesVarNames[10] = "Nsubjet2";
274  //fShapesVarNames[11] = "DeltaR";
275  //fShapesVarNames[12] = "OpenAngle";
276  //fShapesVarNames[13] = "weightPythia";
277 
278  //fShapesVarNames[14] = "NT70";
279  //fShapesVarNames[15] = "nConstNT70";
280  //fShapesVarNames[16] = "NT80";
281  //fShapesVarNames[17] = "nConstNT80";
282  //fShapesVarNames[18] = "NT90";
283  //fShapesVarNames[19] = "nConstNT90";
284  //fShapesVarNames[20] = "NT95";
285  //fShapesVarNames[21] = "nConstNT95";
286 
287  //fShapesVarNames[22] = "SubjetFraction";
288 
289 
290  for(Int_t ivar=0; ivar < nVar; ivar++){
291  cout<<"looping over variables"<<endl;
292  fTreeObservableTagging->Branch(fShapesVarNames[ivar].Data(), &fShapesVar[ivar], Form("%s/F", fShapesVarNames[ivar].Data()));}
293 
294 
295 
296  fPhiJetCorr6= new TH2F("fPhiJetCorr6", "fPhiJetCorr6", 50, 0, 2*TMath::Pi(), 50, 0, 2*TMath::Pi());
297  fOutput->Add(fPhiJetCorr6);
298  fEtaJetCorr6= new TH2F("fEtaJetCorr6", "fEtaJetCorr6", 50, -1.5, 1.5, 50, -1.5, 1.5);
299  fOutput->Add(fEtaJetCorr6);
300 
301  fPhiJetCorr7= new TH2F("fPhiJetCorr7", "fPhiJetCorr7", 50, 0, 2*TMath::Pi(), 50, 0, 2*TMath::Pi());
302  fOutput->Add(fPhiJetCorr7);
303  fEtaJetCorr7= new TH2F("fEtaJetCorr7", "fEtaJetCorr7", 50, -1.5, 1.5, 50, -1.5, 1.5);
304  fOutput->Add(fEtaJetCorr7);
305 
306  fPtJetCorr= new TH2F("fPtJetCorr", "fPtJetCorr", 100, 0, 200, 100, 0, 200);
307  fOutput->Add(fPtJetCorr);
308  fPtJet= new TH1F("fPtJet", "fPtJet", 100, 0, 200);
309  fOutput->Add(fPtJet);
310 
311  fhpTjetpT= new TH2F("fhpTjetpT", "fhpTjetpT", 200, 0, 200, 200, 0, 200);
312  fOutput->Add(fhpTjetpT);
313  fhPt= new TH1F("fhPt", "fhPt", 200, 0, 200);
314  fOutput->Add(fhPt);
315  fhPhi= new TH1F("fhPhi", "fhPhi", 100, -TMath::Pi(), TMath::Pi());
316  fOutput->Add(fhPhi);
317 
318 
319 
320  //log(1/theta),log(z*theta),jetpT,algo//
321  const Int_t dimSpec = 6;
322  const Int_t nBinsSpec[6] = {50,50,10,3,22,10};
323  const Double_t lowBinSpec[6] = {0.0,-10, 0,0,0,0};
324  const Double_t hiBinSpec[6] = {5.0, 0,200,3,22,10};
325  fHLundIterative = new THnSparseF("fHLundIterative",
326  "LundIterativePlot [log(1/theta),log(z*theta),pTjet,algo,partonFlavor,depth]",
327  dimSpec,nBinsSpec,lowBinSpec,hiBinSpec);
328  fOutput->Add(fHLundIterative);
329 
330  if(fAdditionalTracks>0){
331  //log(1/theta),log(z*theta),jetpT of added tracks
332  const Int_t dimSpecb = 4;
333  const Int_t nBinsSpecb[4] = {50,50,20,20};
334  const Double_t lowBinSpecb[4] = {0.0,-10, 0,0};
335  const Double_t hiBinSpecb[4] = {5.0, 0,200,200};
336  fHLundIterativeInject = new THnSparseF("fHLundIterativeInject",
337  "LundIterativePlotInject [log(1/theta),log(z*theta),pTjet,algo]",
338  dimSpecb,nBinsSpecb,lowBinSpecb,hiBinSpecb);
340 
341 
342  fNbOfConstvspT=new TH2F("fNbOfConstvspT", "fNbOfConstvspT", 100, 0, 100, 200, 0, 200);
343  fOutput->Add(fNbOfConstvspT);
344 
345  //fOutput->Add(fTreeObservableTagging);
346 
347  fRandom = new TRandom3(0);
348  PostData(1, fOutput); // Post data for ALL output slots > 0 here
349  PostData(2, fTreeObservableTagging);
350 
351  delete [] fShapesVarNames;
352 
353  }
354 
355 //________________________________________________________________________
357 {
358  // Run analysis code here, if needed. It will be executed before FillHistograms().
359  // if (gRandom) delete gRandom;
360  // gRandom = new TRandom3(0);
361  return kTRUE;
362 }
363 
364 //________________________________________________________________________
366 {
367  // Fill histograms.
368  //cout<<"IntoFillHistograms"<<endl;
369  AliEmcalJet* jet1 = NULL;
370  AliJetContainer *jetCont = GetJetContainer(0);
371 
372  Float_t kWeight=1;
373  if (fCentSelectOn)
374  if ((fCent>fCentMax) || (fCent<fCentMin)) return 0;
375 
376  AliAODTrack *triggerHadron = 0x0;
377 
378  if (fJetSelection == kRecoil) {
379  //Printf("Recoil jets!!!, fminpTTrig = %f, fmaxpTTrig = %f", fminpTTrig, fmaxpTTrig);
380  Int_t triggerHadronLabel = SelectTrigger(fminpTTrig, fmaxpTTrig);
381 
382 
383  if (triggerHadronLabel==-99999) {
384  //Printf ("Trigger Hadron not found, return");
385  return 0;}
386 
387 
389  TClonesArray *trackArrayAn = partContAn->GetArray();
390  triggerHadron = static_cast<AliAODTrack*>(trackArrayAn->At(triggerHadronLabel));
391 
392  if (!triggerHadron) {
393  //Printf("No Trigger hadron with the found label!!");
394  return 0;
395  }
396 
397  if(fSemigoodCorrect){
398  Double_t disthole=RelativePhi(triggerHadron->Phi(),fHolePos);
399  if(TMath::Abs(disthole)+fHoleWidth>TMath::Pi()-fangWindowRecoil){
400  return 0;}
401  }
402 
403  fhPt->Fill(triggerHadron->Pt());
404 
405  }
406 
407  if(jetCont) {
408  jetCont->ResetCurrentID();
409  while((jet1 = jetCont->GetNextAcceptJet())) {
410  //Printf("jet1=%p", jet1);
411  if (!jet1) continue;
412 
413  fPtJet->Fill(jet1->Pt());
414 
415 
416 
417 
419  Double_t disthole=RelativePhi(jet1->Phi(),fHolePos);
420  if(TMath::Abs(disthole)<fHoleWidth){
421  continue;
422  }
423  }
424 
425  Float_t dphiRecoil = 0.;
426  if (fJetSelection == kRecoil){
427  dphiRecoil = RelativePhi(triggerHadron->Phi(), jet1->Phi());
428  if (TMath::Abs(dphiRecoil) < (TMath::Pi() - fangWindowRecoil)) {
429  // Printf("Recoil jets back to back not found! continuing");
430  continue;
431  }
432 
433  fhpTjetpT->Fill(triggerHadron->Pt(), jet1->Pt());
434  //Printf(" ************ FILLING HISTOS****** shapeSub = %d, triggerHadron = %f, jet1 = %f", fJetShapeSub, triggerHadron->Pt(), jet1->Pt());
435  fhPhi->Fill(RelativePhi(triggerHadron->Phi(), jet1->Phi()));
436 
437  }
438 
439 
440  fShapesVar[0] = 0.;
441 
442  if (fJetShapeType == kGenShapes){
443  const AliEmcalPythiaInfo *partonsInfo = 0x0;
444  partonsInfo = GetPythiaInfo();
445  //Printf("partonsInfo=%p", partonsInfo);
446  Double_t jp1=RelativePhi(jet1->Phi(),partonsInfo->GetPartonPhi6());
447  Double_t detap1=(jet1->Eta())-(partonsInfo->GetPartonEta6());
448  kWeight=partonsInfo->GetPythiaEventWeight();
449  //Printf("kWeight=%f", kWeight);
450  fShapesVar[8] = kWeight;
451 
452  Float_t dRp1 = TMath::Sqrt(jp1 * jp1 + detap1 * detap1);
453  fEtaJetCorr6->Fill(jet1->Eta(), partonsInfo->GetPartonEta6());
454  fPhiJetCorr6->Fill(jet1->Phi(), partonsInfo->GetPartonPhi6());
455  if(dRp1 < fRMatching) {
456  fShapesVar[0] = partonsInfo->GetPartonFlag6();
457  fPtJetCorr ->Fill(partonsInfo->GetPartonPt6(), jet1->Pt());
458  }
459  else {
460  jp1=RelativePhi(jet1->Phi(),partonsInfo->GetPartonPhi7());
461  detap1=(jet1->Eta())-(partonsInfo->GetPartonEta7());
462  dRp1 = TMath::Sqrt(jp1 * jp1 + detap1 * detap1);
463  fEtaJetCorr7->Fill(jet1->Eta(), partonsInfo->GetPartonEta7());
464  fPhiJetCorr7->Fill(jet1->Phi(), partonsInfo->GetPartonPhi7());
465  if(dRp1 < fRMatching) {
466  fShapesVar[0] = partonsInfo->GetPartonFlag7();
467  fPtJetCorr->Fill(partonsInfo->GetPartonPt7(), jet1->Pt());
468  }
469  else fShapesVar[0]=0;
470  }
471  }
472 
473  Double_t ptSubtracted = 0;
474  ptSubtracted= jet1->Pt();
475  //Printf("ptSubtracted=%f", ptSubtracted);
476 
477 
478  if (ptSubtracted < fPtThreshold) continue;
479 
480  if ((fCentSelectOn == kFALSE) && (jet1->GetNumberOfTracks() < 2)) continue;
481 
482  //AliEmcalJetFinder *Reclusterer1; //Object containg Subjets from Subtracted Hybrid Jets
483  //Reclusterer1 = Recluster(jet1, 0, fJetRadius, fSubjetRadius, 1, 0, "SubJetFinder_1");
484 
485 
486 
487 
488 
489 
490 
491  fShapesVar[1] = ptSubtracted;
492  fShapesVar[2] = GetJetpTD(jet1,0);
493  fShapesVar[3] = GetJetMass(jet1,0);
494  fShapesVar[4] = 1.*GetJetNumberOfConstituents(jet1,0);
495  fShapesVar[5] = GetJetAngularity(jet1,0);
496 
497 
498  RecursiveParents(jet1,jetCont,0,fShapesVar[0]);
499  RecursiveParents(jet1,jetCont,1,fShapesVar[0]);
500  RecursiveParents(jet1,jetCont,2,fShapesVar[0]);
501 
502 
503  fTreeObservableTagging->Fill();
504 
505  }
506 
507  }
508 
509  return kTRUE;
510 }
511 
512 //________________________________________________________________________
514  //calc subtracted jet mass
515  if((fJetShapeSub==kDerivSub)&&(jetContNb==0))
517  else return jet->GetShapeProperties()->GetSecondOrderSubtracted();
518  else
519  return jet->M();
520 }
521 
522 //________________________________________________________________________
524 
525  AliJetContainer *jetCont = GetJetContainer(jetContNb);
526  if (!jet->GetNumberOfTracks())
527  return 0;
528  Double_t den=0.;
529  Double_t num = 0.;
530  AliVParticle *vp1 = 0x0;
531  for(UInt_t i = 0; i < jet->GetNumberOfTracks(); i++) {
532  vp1 = static_cast<AliVParticle*>(jet->TrackAt(i, jetCont->GetParticleContainer()->GetArray()));
533 
534  if (!vp1){
535  Printf("AliVParticle associated to constituent not found");
536  continue;
537  }
538 
539  Double_t dphi = RelativePhi(vp1->Phi(),jet->Phi());
540  Double_t dr2 = (vp1->Eta()-jet->Eta())*(vp1->Eta()-jet->Eta()) + dphi*dphi;
541  Double_t dr = TMath::Sqrt(dr2);
542  num=num+vp1->Pt()*dr;
543  den=den+vp1->Pt();
544  }
545  return num/den;
546 }
547 
548 //________________________________________________________________________
550 
551  if((fJetShapeSub==kDerivSub) && (jetContNb==0))
554  else
555  return Angularity(jet, jetContNb);
556 
557 }
558 
559 
560 //________________________________________________________________________
562 
563  AliJetContainer *jetCont = GetJetContainer(jetContNb);
564  if (!jet->GetNumberOfTracks())
565  return 0;
566  Double_t den=0.;
567  Double_t num = 0.;
568  AliVParticle *vp1 = 0x0;
569  for(UInt_t i = 0; i < jet->GetNumberOfTracks(); i++) {
570  vp1 = static_cast<AliVParticle*>(jet->TrackAt(i, jetCont->GetParticleContainer()->GetArray()));
571 
572  if (!vp1){
573  Printf("AliVParticle associated to constituent not found");
574  continue;
575  }
576 
577  num=num+vp1->Pt()*vp1->Pt();
578  den=den+vp1->Pt();
579  }
580  return TMath::Sqrt(num)/den;
581 }
582 
583 //________________________________________________________________________
585  //calc subtracted jet mass
586  if((fJetShapeSub==kDerivSub)&&(jetContNb==0))
588  else return jet->GetShapeProperties()->GetSecondOrderSubtractedpTD();
589  else
590  return PTD(jet, jetContNb);
591 
592 }
593 
594 //_____________________________________________________________________________
596 
597  AliJetContainer *jetCont = GetJetContainer(jetContNb);
598  if (!jet->GetNumberOfTracks())
599  return 0;
600  Double_t mxx = 0.;
601  Double_t myy = 0.;
602  Double_t mxy = 0.;
603  int nc = 0;
604  Double_t sump2 = 0.;
605  Double_t pxjet=jet->Px();
606  Double_t pyjet=jet->Py();
607  Double_t pzjet=jet->Pz();
608 
609 
610  //2 general normalized vectors perpendicular to the jet
611  TVector3 ppJ1(pxjet, pyjet, pzjet);
612  TVector3 ppJ3(- pxjet* pzjet, - pyjet * pzjet, pxjet * pxjet + pyjet * pyjet);
613  ppJ3.SetMag(1.);
614  TVector3 ppJ2(-pyjet, pxjet, 0);
615  ppJ2.SetMag(1.);
616  AliVParticle *vp1 = 0x0;
617  for(UInt_t i = 0; i < jet->GetNumberOfTracks(); i++) {
618  vp1 = static_cast<AliVParticle*>(jet->TrackAt(i, jetCont->GetParticleContainer()->GetArray()));
619 
620  if (!vp1){
621  Printf("AliVParticle associated to constituent not found");
622  continue;
623  }
624 
625  TVector3 pp(vp1->Px(), vp1->Py(), vp1->Pz());
626 
627  //local frame
628  TVector3 pLong = pp.Dot(ppJ1) / ppJ1.Mag2() * ppJ1;
629  TVector3 pPerp = pp - pLong;
630  //projection onto the two perpendicular vectors defined above
631 
632  Float_t ppjX = pPerp.Dot(ppJ2);
633  Float_t ppjY = pPerp.Dot(ppJ3);
634  Float_t ppjT = TMath::Sqrt(ppjX * ppjX + ppjY * ppjY);
635  if(ppjT<=0) return 0;
636 
637  mxx += (ppjX * ppjX / ppjT);
638  myy += (ppjY * ppjY / ppjT);
639  mxy += (ppjX * ppjY / ppjT);
640  nc++;
641  sump2 += ppjT;}
642 
643  if(nc<2) return 0;
644  if(sump2==0) return 0;
645  // Sphericity Matrix
646  Double_t ele[4] = {mxx / sump2, mxy / sump2, mxy / sump2, myy / sump2};
647  TMatrixDSym m0(2,ele);
648 
649  // Find eigenvectors
650  TMatrixDSymEigen m(m0);
651  TVectorD eval(2);
652  TMatrixD evecm = m.GetEigenVectors();
653  eval = m.GetEigenValues();
654  // Largest eigenvector
655  int jev = 0;
656  // cout<<eval[0]<<" "<<eval[1]<<endl;
657  if (eval[0] < eval[1]) jev = 1;
658  TVectorD evec0(2);
659  // Principle axis
660  evec0 = TMatrixDColumn(evecm, jev);
661  Double_t compx=evec0[0];
662  Double_t compy=evec0[1];
663  TVector2 evec(compx, compy);
664  Double_t circ=0;
665  if(jev==1) circ=2*eval[0];
666  if(jev==0) circ=2*eval[1];
667 
668  return circ;
669 
670 
671 
672 }
673 
674 
675 
676 
677 //________________________________________________________________________
679  //calc subtracted jet mass
680 
681  if((fJetShapeSub==kDerivSub)&&(jetContNb==0))
684  else
685  return Circularity(jet, jetContNb);
686 
687 }
688 
689 //________________________________________________________________________
691 
692  AliJetContainer *jetCont = GetJetContainer(jetContNb);
693  if (!jet->GetNumberOfTracks())
694  return 0;
695  Double_t den=0.;
696  Double_t num = 0.;
697  AliVParticle *vp1 = 0x0;
698  AliVParticle *vp2 = 0x0;
699  std::vector<int> ordindex;
700  ordindex=jet->GetPtSortedTrackConstituentIndexes(jetCont->GetParticleContainer()->GetArray());
701  //Printf("Nbof const = %d", jet->GetNumberOfTracks());
702  //Printf("ordindex[0] = %d, ordindex[1] = %d", ordindex[0], ordindex[1]);
703 
704  if(ordindex.size()<2) return -1;
705 
706  vp1 = static_cast<AliVParticle*>(jet->TrackAt(ordindex[0], jetCont->GetParticleContainer()->GetArray()));
707  if (!vp1){
708  Printf("AliVParticle associated to Leading constituent not found");
709  return -1;
710  }
711 
712  vp2 = static_cast<AliVParticle*>(jet->TrackAt(ordindex[1], jetCont->GetParticleContainer()->GetArray()));
713  if (!vp2){
714  Printf("AliVParticle associated to Subleading constituent not found");
715  return -1;
716  }
717 
718 
719  num=vp1->Pt();
720  den=vp2->Pt();
721  //Printf("vp1->Pt() =%f, vp2->Pt() =%f", vp1->Pt(), vp2->Pt());
722 
723 return num-den;
724 }
725 
726 //________________________________________________________________________
728  //calc subtracted jet mass
729 
730  if((fJetShapeSub==kDerivSub)&&(jetContNb==0))
733  else
734  return LeSub(jet, jetContNb);
735 
736 }
737 
738 //________________________________________________________________________
740  //calc subtracted jet mass
741 
742  if((fJetShapeSub==kDerivSub)&&(jetContNb==0))
745  else
746  return jet->GetNumberOfTracks();
747 
748  }
749 
750 
751 //________________________________________________________________________
752 AliEmcalJetFinder *AliAnalysisTaskEmcalJetShapesMC::Recluster(AliEmcalJet *Jet, Int_t JetContNb, Double_t JetRadius, Double_t SubJetRadius, Double_t SubJetMinPt, Int_t Algorithm, const char* Name){
753 
754  AliJetContainer *JetCont = GetJetContainer(JetContNb);
755  AliEmcalJetFinder *Reclusterer = new AliEmcalJetFinder(Name); //JetFinder Object for reclustered jets
756  Reclusterer->SetRadius(SubJetRadius);
757  Reclusterer->SetJetMinPt(SubJetMinPt);
758  Reclusterer->SetJetAlgorithm(Algorithm); //0 for anti-kt 1 for kt
759  Reclusterer->SetJetMaxEta(0.9-JetRadius);
760  Reclusterer->SetRecombSheme(0);
761 
762 
763  //Double_t dVtx[3]={vert->GetX(),vert->GetY(),vert->GetZ()};
764  Double_t dVtx[3]={0.,0.,0.};
765  if(Reclusterer->AliEmcalJetFinder::Filter(Jet, JetCont, dVtx)){;} //reclustering jet1 using the jetfinderobject Reclusterer
766  return Reclusterer;
767 }
768 
769 
770 
771 
772 //----------------------------------------------------------------------
774  AliJetContainer *JetCont = GetJetContainer(JetContNb);
775 
776 
777 
778  Double_t SubJetiness_Numerator = 0;
779  Double_t SubJetiness_Denominator = 0;
780  Double_t Index=-2;
781  if (Reclusterer->GetNumberOfJets() < 1) return -2;
782  Index=SubJetOrdering(Jet,Reclusterer,1,0,kTRUE);
783  if(Index==-999) return -2;
784  SubJetiness_Numerator=(Reclusterer->GetJet(Index)->Pt());
785  SubJetiness_Denominator=Jet->Pt();
786  return SubJetiness_Numerator/SubJetiness_Denominator;
787 
788 
789 }
790 //__________________________________________________________________________________
792 
793  AliJetContainer *jetCont = GetJetContainer(jetContNb);
794  if (!jet->GetNumberOfTracks())
795  return 0;
796  Double_t den=0.;
797  Double_t num = 0.;
798  AliVParticle *vp1 = 0x0;
799  for(UInt_t i = 0; i < jet->GetNumberOfTracks(); i++) {
800  vp1 = static_cast<AliVParticle*>(jet->TrackAt(i, jetCont->GetParticleContainer()->GetArray()));
801 
802  if (!vp1){
803  Printf("AliVParticle associated to constituent not found");
804  continue;
805  }
806 
807  Double_t dphi = RelativePhi(vp1->Phi(),jet->Phi());
808  Double_t dr2 = (vp1->Eta()-jet->Eta())*(vp1->Eta()-jet->Eta()) + dphi*dphi;
809  Double_t dr = TMath::Sqrt(dr2);
810  if(dr<=fSubjetRadius) num=num+vp1->Pt();
811 
812  }
813  return num/jet->Pt();
814 }
815 
816 
817 
818 
819 //________________________________________________________________________
821  //calc subtracted jet mass
822 
823  if((fJetShapeSub==kDerivSub)&&(jetContNb==0))
826  else
827  return CoreFrac(jet, jetContNb);
828 
829 }
830 
831 
832 
833 
834 //----------------------------------------------------------------------
836  AliEmcalJet *SubJet=NULL;
837  Double_t SortingVariable;
838  Int_t ArraySize =N+1;
839  TArrayD *JetSorter = new TArrayD(ArraySize);
840  TArrayD *JetIndexSorter = new TArrayD(ArraySize);
841  for (Int_t i=0; i<ArraySize; i++){
842  JetSorter->SetAt(0,i);
843  }
844  for (Int_t i=0; i<ArraySize; i++){
845  JetIndexSorter->SetAt(0,i);
846  }
847  if(Reclusterer->GetNumberOfJets()<N) return -999;
848  for (Int_t i=0; i<Reclusterer->GetNumberOfJets(); i++){
849  SubJet=Reclusterer->GetJet(i);
850  if (Type==0) SortingVariable=SubJet->Pt();
851  else if (Type==1) SortingVariable=SubJet->E();
852  else if (Type==2) SortingVariable=SubJet->M();
853  for (Int_t j=0; j<N; j++){
854  if (SortingVariable>JetSorter->GetAt(j)){
855  for (Int_t k=N-1; k>=j; k--){
856  JetSorter->SetAt(JetSorter->GetAt(k),k+1);
857  JetIndexSorter->SetAt(JetIndexSorter->GetAt(k),k+1);
858  }
859  JetSorter->SetAt(SortingVariable,j);
860  JetIndexSorter->SetAt(i,j);
861  break;
862  }
863  }
864  }
865  if (!Index) return JetSorter->GetAt(N-1);
866  else return JetIndexSorter->GetAt(N-1);
867 }
868 
869 
870 
871 //returns -1 if the Nth hardest jet is requested where N>number of available jets
872 //type: 0=Pt 1=E 2=M
873 //Index TRUE=returns index FALSE=returns value of quantatiy in question
874 
875 
876 
877 
878 
879 //______________________________________________________________________________
881 
882  AliJetContainer *jetCont = GetJetContainer(jetContNb);
883  if (!jet->GetNumberOfTracks())
884  return 0;
885  Double_t mxx = 0.;
886  Double_t myy = 0.;
887  Double_t mxy = 0.;
888  int nc = 0;
889  Double_t sump2 = 0.;
890 
891  AliVParticle *vp1 = 0x0;
892  for(UInt_t i = 0; i < jet->GetNumberOfTracks(); i++) {
893  vp1 = static_cast<AliVParticle*>(jet->TrackAt(i, jetCont->GetParticleContainer()->GetArray()));
894 
895  if (!vp1){
896  Printf("AliVParticle associated to constituent not found");
897  continue;
898  }
899 
900  Double_t ppt=vp1->Pt();
901  Double_t dphi = RelativePhi(vp1->Phi(),jet->Phi());
902 
903  Double_t deta = vp1->Eta()-jet->Eta();
904  mxx += ppt*ppt*deta*deta;
905  myy += ppt*ppt*dphi*dphi;
906  mxy -= ppt*ppt*deta*TMath::Abs(dphi);
907  nc++;
908  sump2 += ppt*ppt;
909 
910  }
911  if(nc<2) return 0;
912  if(sump2==0) return 0;
913  // Sphericity Matrix
914  Double_t ele[4] = {mxx , mxy , mxy , myy };
915  TMatrixDSym m0(2,ele);
916 
917  // Find eigenvectors
918  TMatrixDSymEigen m(m0);
919  TVectorD eval(2);
920  TMatrixD evecm = m.GetEigenVectors();
921  eval = m.GetEigenValues();
922  // Largest eigenvector
923  int jev = 0;
924  // cout<<eval[0]<<" "<<eval[1]<<endl;
925  if (eval[0] < eval[1]) jev = 1;
926  TVectorD evec0(2);
927  // Principle axis
928  evec0 = TMatrixDColumn(evecm, jev);
929  Double_t compx=evec0[0];
930  Double_t compy=evec0[1];
931  TVector2 evec(compx, compy);
932  Double_t sig=0;
933  if(jev==1) sig=TMath::Sqrt(TMath::Abs(eval[0])/sump2);
934  if(jev==0) sig=TMath::Sqrt(TMath::Abs(eval[1])/sump2);
935 
936  return sig;
937 
938 }
939 
940 //________________________________________________________________________
942  //calc subtracted jet mass
943 
944  if((fJetShapeSub==kDerivSub)&&(jetContNb==0))
947  else
948  return Sigma2(jet, jetContNb);
949 
950 }
951 
952 
953 //_________________________________________________________________________
955 
956  AliJetContainer *jetCont = GetJetContainer(jetContNb);
957  if (!jet->GetNumberOfTracks())
958  return;
959 
960  Double_t ptJet = jet->Pt();
961 
962  AliVParticle *vp1 = 0x0;
963  std::vector<int> ordindex;
964  ordindex=jet->GetPtSortedTrackConstituentIndexes(jetCont->GetParticleContainer()->GetArray());
965  //Printf("Nbof const = %d", jet->GetNumberOfTracks());
966  //Printf("ordindex[0] = %d, ordindex[1] = %d", ordindex[0], ordindex[1]);
967  //if(ordindex.size()<2) return -1;
968 
969  for(Int_t iconst =0; iconst<jet->GetNumberOfTracks(); iconst++){
970 
971  vp1 = static_cast<AliVParticle*>(jet->TrackAt(ordindex[iconst], jetCont->GetParticleContainer()->GetArray()));
972  if (!vp1){
973  Printf("AliVParticle associated to Leading constituent not found");
974  return;
975  }
976 
977  if (nTFractions[0] <= 0.7*ptJet){
978  nTFractions[0] += vp1->Pt();
979  nTFractions[1] +=1;
980  }
981 
982  if (nTFractions[2] <= 0.8*ptJet){
983  nTFractions[2] += vp1->Pt();
984  nTFractions[3] +=1;
985  }
986 
987  if (nTFractions[4] <= 0.9*ptJet){
988  nTFractions[4] += vp1->Pt();
989  nTFractions[5] +=1;
990  }
991 
992  if (nTFractions[6] <= 0.95*ptJet){
993  nTFractions[6] += vp1->Pt();
994  nTFractions[7] +=1;
995  }
996  }
997 }
998 //_________________________________________________________________________________________________
1000 
1001  //WARNING!!! Only works for parent jets that are clustered with Anti-Kt! To change go to AliEmcalJetFinder.cxx and look at the Nsubjettiness() function
1002 
1003  //Algorithm==0 -> kt_axes;
1004  // Algorithm==1 -> ca_axes;
1005  //Algorithm==2 -> antikt_0p2_axes;
1006  //Algorithm==3 -> wta_kt_axes;
1007  //Algorithm==4 -> wta_ca_axes;
1008  //Algorithm==5 -> onepass_kt_axes;
1009  //Algorithm==6 -> onepass_ca_axes;
1010  //Algorithm==7 -> onepass_antikt_0p2_axes;
1011  //Algorithm==8 -> onepass_wta_kt_axes;
1012  //Algorithm==9 -> onepass_wta_ca_axes;
1013  //Algorithm==10 -> min_axes;
1014 
1015 
1016  //Option==0 returns Nsubjettiness Value
1017  //Option==1 && N==2 returns opening angle between two subjet axes(Delta R?)
1018  //Option==2 && N==2 returns Delta R
1019 
1020  if (Jet->GetNumberOfTracks()>=N){
1021  AliJetContainer *JetCont = GetJetContainer(JetContNb);
1022  AliEmcalJetFinder *JetFinder=new AliEmcalJetFinder("Nsubjettiness");
1023  JetFinder->SetJetMaxEta(0.9-fJetRadius);
1024  JetFinder->SetRadius(fJetRadius);
1025  JetFinder->SetJetAlgorithm(0); //0 for anti-kt 1 for kt //this is for the JET!!!!!!!!!! Not the SubJets
1026  JetFinder->SetRecombSheme(0);
1027  JetFinder->SetJetMinPt(Jet->Pt());
1028  //Double_t dVtx[3]={vert->GetX(),vert->GetY(),vert->GetZ()};
1029  Double_t dVtx[3]={1,1,1};
1030  //Printf("JetFinder->Nsubjettiness =%f", JetFinder->Nsubjettiness(Jet,JetCont,dVtx,N,Algorithm,fSubjetRadius,Beta,Option));
1031  return JetFinder->Nsubjettiness(Jet,JetCont,dVtx,N,Algorithm,0.2,Beta,Option,0,Beta_SD,ZCut,SoftDropOn);
1032 
1033  }
1034  else return -2;
1035 }
1036 
1037 
1038 
1039 //________________________________________________________________________
1041 
1043  TClonesArray *tracksArray = partCont->GetArray();
1044 
1045  if(!partCont || !tracksArray) return -99999;
1046  AliAODTrack *track = 0x0;
1047  AliEmcalParticle *emcPart = 0x0;
1048 
1049 
1050  TList *trackList = new TList();
1051  Int_t triggers[100];
1052  for (Int_t iTrigger=0; iTrigger<100; iTrigger++) triggers[iTrigger] = 0;
1053  Int_t iTT = 0;
1054 
1055  for(Int_t iTrack=0; iTrack <= tracksArray->GetEntriesFast(); iTrack++){
1056 
1057 
1058  if (fJetShapeSub == kConstSub){
1059  emcPart = static_cast<AliEmcalParticle*>(tracksArray->At(iTrack));
1060  if (!emcPart) continue;
1061  if(TMath::Abs(emcPart->Eta())>0.9) continue;
1062  if (emcPart->Pt()<0.15) continue;
1063 
1064  if ((emcPart->Pt() >= minpT) && (emcPart->Pt()< maxpT)) {
1065  trackList->Add(emcPart);
1066  triggers[iTT] = iTrack;
1067  iTT++;
1068  }
1069  }
1070  else{
1071  track = static_cast<AliAODTrack*>(tracksArray->At(iTrack));
1072  if (!track) continue;
1073  if(TMath::Abs(track->Eta())>0.9) continue;
1074  if (track->Pt()<0.15) continue;
1075  if (!(track->TestFilterBit(768))) continue;
1076 
1077  if ((track->Pt() >= minpT) && (track->Pt()< maxpT)) {
1078  trackList->Add(track);
1079  triggers[iTT] = iTrack;
1080  iTT++;
1081 
1082  }
1083  }
1084  }
1085 
1086  if (iTT == 0) return -99999;
1087  Int_t nbRn = 0, index = 0 ;
1088  TRandom3* random = new TRandom3(0);
1089  nbRn = random->Integer(iTT);
1090 
1091  index = triggers[nbRn];
1092  //Printf("iTT Total= %d, nbRn = %d, Index = %d",iTT, nbRn, index );
1093  return index;
1094 
1095 }
1096 
1097 //__________________________________________________________________________________
1099 
1100  if (vphi < -1*TMath::Pi()) vphi += (2*TMath::Pi());
1101  else if (vphi > TMath::Pi()) vphi -= (2*TMath::Pi());
1102  if (mphi < -1*TMath::Pi()) mphi += (2*TMath::Pi());
1103  else if (mphi > TMath::Pi()) mphi -= (2*TMath::Pi());
1104  double dphi = mphi-vphi;
1105  if (dphi < -1*TMath::Pi()) dphi += (2*TMath::Pi());
1106  else if (dphi > TMath::Pi()) dphi -= (2*TMath::Pi());
1107  return dphi;//dphi in [-Pi, Pi]
1108 }
1109 
1110 
1111 //________________________________________________________________________
1113  //
1114  // retrieve event objects
1115  //
1117  return kFALSE;
1118 
1119  return kTRUE;
1120 }
1121 
1122 //_______________________________________________________________________
1124 {
1125  // Called once at the end of the analysis.
1126 
1127  AliInfo("Terminate");
1128  AliAnalysisTaskSE::Terminate();
1129 
1130  fOutput = dynamic_cast<AliEmcalList*> (GetOutputData(1));
1131  if (!fOutput) {
1132  AliError("fOutput not available");
1133  return;
1134  }
1135 
1136  fTreeObservableTagging = dynamic_cast<TTree*>(GetOutputData(2));
1137  if (!fTreeObservableTagging){
1138  Printf("ERROR: fTreeObservableTagging not available");
1139  return;
1140  }
1141 
1142 }
1143 
1144 //_________________________________________________________________________
1146 
1147  std::vector<fastjet::PseudoJet> fInputVectors;
1148  fInputVectors.clear();
1149  Double_t JetInvMass=0, PseudJetInvMass=0, TrackMom = 0, TrackEnergy = 0;
1150  fastjet::PseudoJet PseudoTracks;
1151  fastjet::PseudoJet MyJet;
1152  fastjet::PseudoJet PseudoTracksCMS;
1153  AliParticleContainer *fTrackCont = fJetCont->GetParticleContainer();
1154  //cout<<"CALL TO SOFTDROP"<<endl;
1155 
1156  Double_t zeta=0;
1157  Double_t angle=0;
1158  if (fTrackCont) for (Int_t i=0; i<fJet->GetNumberOfTracks(); i++) {
1159  AliVParticle *fTrk = fJet->TrackAt(i, fTrackCont->GetArray());
1160  if (!fTrk) continue;
1161  JetInvMass += fTrk->M();
1162 
1163  PseudoTracks.reset(fTrk->Px(), fTrk->Py(), fTrk->Pz(),fTrk->E());
1164  TrackMom += TMath::Sqrt(TMath::Power(fTrk->Px(),2)+TMath::Power(fTrk->Py(),2)+TMath::Power(fTrk->Pz(),2));
1165  TrackEnergy += fTrk->E();
1166  PseudoTracks.set_user_index(fJet->TrackAt(i)+100);
1167  PseudJetInvMass += PseudoTracks.m();
1168  fInputVectors.push_back(PseudoTracks);
1169 
1170  }
1171 
1172  //fRandom->SetSeed(0);
1173  //here add N tracks with random phi and eta and theta according to bdmps distrib.
1174 
1175  MyJet.reset(fJet->Px(),fJet->Py(),fJet->Pz(),fJet->E());
1176  Double_t omegac=0.5*fqhat*fxlength*fxlength/0.2;
1177  Double_t thetac=TMath::Sqrt(12*0.2/(fqhat*TMath::Power(fxlength,3)));
1178  Double_t xQs=TMath::Sqrt(fqhat*fxlength);
1179  //cout<<"medium parameters "<<omegac<<" "<<thetac<<" "<<xQs<<endl;
1180 
1181  for(Int_t i=0;i<fAdditionalTracks;i++){
1182 
1183  Double_t ppx,ppy,ppz,kTscale,lim2o,lim1o;
1184  Double_t lim2=xQs;
1185  Double_t lim1=10000;
1186 
1187  //generation of kT according to 1/kT^4, with minimum QS=2 GeV and maximum ~sqrt(ptjet*T)
1188  fTf1Kt= new TF1("fTf1Kt","1/(x*x*x*x)",lim2,lim1);
1189  kTscale=fTf1Kt->GetRandom();
1190  //generation within the jet cone
1191 
1192  //generation of w according to 1/w, with minimum wc
1193  //omega needs to be larger than kT so to have well defined angles
1194  lim2o=kTscale;
1195  lim1o=kTscale/TMath::Sin(0.1);
1196  fTf1Omega= new TF1("fTf1Omega","1/x",lim2o,lim1o);
1197  Double_t omega=fTf1Omega->GetRandom();
1198 
1199  Double_t sinpptheta=kTscale/omega;
1200  Double_t pptheta=TMath::ASin(sinpptheta);
1201  //cout<<"angle_omega_kt"<<pptheta<<" "<<omega<<" "<<kTscale<<endl;
1202  if(pptheta>fJetRadius) continue;
1203 
1204  PseudoTracksCMS.reset(kTscale/TMath::Sqrt(2),kTscale/TMath::Sqrt(2),omega*TMath::Cos(pptheta),omega);
1205  //boost the particle in the rest frame of the jet to the lab frame
1206  fastjet::PseudoJet PseudoTracksLab=PseudoTracksCMS.boost(MyJet);
1207  PseudoTracksLab.set_user_index(i+fJet->GetNumberOfTracks()+100);
1208  fInputVectors.push_back(PseudoTracksLab);
1209  //in the frame of the jet
1210  zeta=omega/fJet->E();
1211  angle=pptheta;
1212 
1213 
1214  }
1215 
1216 
1217 
1218 
1219 
1220 
1221 
1222 
1223 
1224  fastjet::JetDefinition fJetDef(fastjet::antikt_algorithm, fJetRadius*2, static_cast<fastjet::RecombinationScheme>(0), fastjet::BestFJ30 );
1225  fastjet::contrib::Recluster *recluster;
1226  try {
1227  fastjet::ClusterSequence fClustSeqSA(fInputVectors, fJetDef);
1228  std::vector<fastjet::PseudoJet> fOutputJets;
1229  fOutputJets.clear();
1230  fOutputJets=fClustSeqSA.inclusive_jets(0);
1231 
1232  //cout<<fOutputJets[0].perp()<<" "<<fJet->Pt()<<endl;
1233 
1234  fastjet::contrib::SoftDrop softdrop(beta, zcut);
1235 
1236  softdrop.set_verbose_structure(kTRUE);
1237  fastjet::JetAlgorithm jetalgo(fastjet::cambridge_algorithm);
1238  if(ReclusterAlgo==2) jetalgo=fastjet::antikt_algorithm;
1239  if(ReclusterAlgo==1) jetalgo=fastjet::kt_algorithm;
1240  if(ReclusterAlgo==0) jetalgo=fastjet::cambridge_algorithm;
1241 
1242  recluster = new fastjet::contrib::Recluster(jetalgo,1,true);
1243  softdrop.set_reclustering(true,recluster);
1244  fastjet::PseudoJet finaljet = softdrop(fOutputJets[0]);
1245  Int_t NDroppedTracks = fJet->GetNumberOfTracks()-finaljet.constituents().size();
1246 
1247  Double_t SymParam, Mu, DeltaR, GroomedPt,GroomedMass;
1248  Int_t NGroomedBranches;
1249  SymParam=(finaljet.structure_of<fastjet::contrib::SoftDrop>().symmetry());
1250  Mu=(finaljet.structure_of<fastjet::contrib::SoftDrop>().mu());
1251  DeltaR=(finaljet.structure_of<fastjet::contrib::SoftDrop>().delta_R());
1252  NGroomedBranches=finaljet.structure_of<fastjet::contrib::SoftDrop>().dropped_count();
1253  GroomedPt=finaljet.perp();
1254  GroomedMass=finaljet.m();
1255  // if(beta==0){
1256  // if(ReclusterAlgo==0){
1257  // fShapesVar[12]=SymParam;
1258  // fShapesVar[13]=DeltaR;
1259  // fShapesVar[14]=zeta;
1260  // fShapesVar[15]=angle;
1261  // fShapesVar[16]=GroomedMass;}
1262  // if(ReclusterAlgo==1){
1263  // fShapesVar[17]=SymParam;
1264  // fShapesVar[18]=DeltaR;
1265  // fShapesVar[19]=zeta;
1266  // fShapesVar[20]=angle;
1267  // fShapesVar[21]=GroomedMass; }
1268 
1269  // if(ReclusterAlgo==2){
1270  // fShapesVar[22]=SymParam;
1271  // fShapesVar[23]=DeltaR;
1272  // fShapesVar[24]=zeta;
1273  // fShapesVar[25]=angle;
1274  // fShapesVar[26]=GroomedMass;
1275  // }}
1276  // if(beta==1){
1277  // fShapesVar[27]=SymParam;
1278  // fShapesVar[28]=DeltaR;
1279  // fShapesVar[29]=zeta;
1280  // fShapesVar[30]=angle;
1281  // fShapesVar[31]=GroomedMass;
1282  // }
1283  // //this one kills soft and large angle radiation
1284  // if((beta==1.5) && (zcut==0.5)){
1285  // fShapesVar[32]=SymParam;
1286  // fShapesVar[33]=DeltaR;
1287  // fShapesVar[34]=zeta;
1288  // fShapesVar[35]=angle;
1289  // fShapesVar[36]=GroomedMass; }
1290  // //this option favour democratic branches at large kt
1291  // if((beta==-1) && (zcut==0.005)){
1292  // fShapesVar[37]=SymParam;
1293  // fShapesVar[38]=DeltaR;
1294  // fShapesVar[39]=zeta;
1295  // fShapesVar[40]=angle;
1296  // fShapesVar[41]=GroomedMass; }
1297 
1298  // if((beta==-2) && (zcut==0.005)){
1299  // fShapesVar[42]=SymParam;
1300  // fShapesVar[43]=DeltaR;
1301  // fShapesVar[44]=zeta;
1302  // fShapesVar[45]=angle;
1303  // fShapesVar[46]=GroomedMass; }
1304 
1305 
1306 
1307 
1308 
1309 
1310 
1311 
1312 
1313  } catch (fastjet::Error) {
1314  AliError(" [w] FJ Exception caught.");
1315  //return -1;
1316  }
1317 
1318 
1319 
1320 
1321  if(recluster) { delete recluster; recluster = NULL; }
1322 
1323  if(fTf1Kt){ delete fTf1Kt;}
1324  if(fTf1Omega){ delete fTf1Omega;}
1325 
1326 
1327 
1328  return;
1329 
1330 
1331 }
1332 
1333 
1334 //_________________________________________________________________________
1336 
1337  std::vector<fastjet::PseudoJet> fInputVectors;
1338  fInputVectors.clear();
1339  fastjet::PseudoJet PseudoTracks;
1340  fastjet::PseudoJet PseudoTracksLab;
1341  double xflagalgo=0;
1342  double lnpt_relinject=0;
1343  double yinject=0;
1344  int xflagAdded=0;
1345  double zinject,angleinject,pptheta,sinpptheta,omega,omega2,angle2;
1346  AliParticleContainer *fTrackCont = fJetCont->GetParticleContainer();
1347  Float_t pTscale=0., phiscale=0., thetascale=0., pXscale=0., pYscale=0., pZscale=0., pscale=0.;
1348 
1349  if (fTrackCont) for (Int_t i=0; i<fJet->GetNumberOfTracks(); i++) {
1350  AliVParticle *fTrk = fJet->TrackAt(i, fTrackCont->GetArray());
1351  if (!fTrk) continue;
1352  if (fScaleELoss){
1353  pTscale = xfraction*sqrt(pow(fTrk->Px(),2)+pow(fTrk->Py(),2));
1354  phiscale = fTrk->Phi();
1355  thetascale = 2.*TMath::ATan(TMath::Exp(-1.*(fTrk->Eta())));
1356  pXscale = pTscale * TMath::Cos(phiscale);
1357  pYscale = pTscale * TMath::Sin(phiscale);
1358  pZscale = pTscale/TMath::Tan(thetascale);
1359  pscale = TMath::Sqrt(pTscale*pTscale+pZscale*pZscale);
1360  PseudoTracks.reset(pXscale, pYscale, pZscale, pscale);
1361  }
1362  else PseudoTracks.reset(fTrk->Px(), fTrk->Py(), fTrk->Pz(),fTrk->E());
1363  PseudoTracks.set_user_index(fJet->TrackAt(i)+100);
1364  fInputVectors.push_back(PseudoTracks);
1365 
1366  }
1367  if(fAddMedScat){
1368  for(int i = 0; i < fAddMedScatN; i++){
1369  TRandom3 rand1(0),rand2(0),rand3(0); //set range +- jet R
1370  Double_t randN1 = 0.4*0.4*rand1.Rndm();
1371  Double_t randN2 = 2*TMath::Pi()*rand2.Rndm();
1372  Double_t phi_rand = (fJet->Phi())+TMath::Sqrt(randN1)*TMath::Sin(randN2);
1373  Double_t eta_rand = (fJet->Eta())+TMath::Sqrt(randN1)*TMath::Cos(randN2);
1374  Double_t fAddMedScatPt = (fAddMedScatPtFrac*fJet->Pt())/fAddMedScatN;
1375  PseudoTracks.reset(fAddMedScatPt*TMath::Cos(phi_rand),fAddMedScatPt*TMath::Sin(phi_rand),fAddMedScatPt/TMath::Tan(eta_rand),fAddMedScatPt);
1376  PseudoTracks.set_user_index(i+fJet->GetNumberOfTracks()+100);
1377  fInputVectors.push_back(PseudoTracks);
1378  }
1379  }
1380 
1381 
1382 
1383  //add tracks to the jet prior to the reclusterer in case of iterative mapping of splittings
1384 
1385  Double_t omegac=0.5*fqhat*fxlength*fxlength/0.2;
1386  Double_t thetac=TMath::Sqrt(12*0.2/(fqhat*TMath::Power(fxlength,3)));
1387  Double_t xQs=TMath::Sqrt(fqhat*fxlength);
1388 
1389 
1390  for(Int_t i=0;i<fAdditionalTracks;i++){
1391 
1392  Double_t ppx,ppy,ppz,kTscale,lim2o,lim1o;
1393  Double_t lim2=xQs;
1394  Double_t lim1=10000;
1395 
1396  //generation of kT according to 1/kT^4, with minimum QS=2 GeV and maximum ~sqrt(ptjet*T)
1397  fTf1Kt= new TF1("fTf1Kt","1/(x*x*x*x)",lim2,lim1);
1398  kTscale=fTf1Kt->GetRandom();
1399  //generation within the jet cone
1400 
1401  //generation of w according to 1/w, with minimum wc
1402  //omega needs to be larger than kT so to have well defined angles
1403  lim2o=kTscale;
1404  lim1o=kTscale/TMath::Sin(0.1);
1405  fTf1Omega= new TF1("fTf1Omega","1/x",lim2o,lim1o);
1406  omega=fTf1Omega->GetRandom();
1407  sinpptheta=kTscale/omega;
1408  pptheta=TMath::ASin(sinpptheta);
1409  if(pptheta>fJetRadius) continue;
1410 
1411  //Lorentz vector in the frame where the jet moves along z axis
1412  TLorentzVector pTrackCMS(kTscale/TMath::Sqrt(2),kTscale/TMath::Sqrt(2),omega*TMath::Cos(pptheta),omega);
1413  TVector3 MyJet(fJet->Px(),fJet->Py(),fJet->Pz());
1414  TVector3 direction = MyJet.Unit();
1415  //rotate the track to the jet frame
1416  pTrackCMS.RotateUz(direction);
1417 
1418  //add the rotated track to the jet
1419  PseudoTracksLab.reset(pTrackCMS.Px(),pTrackCMS.Py(),pTrackCMS.Pz(),pTrackCMS.E());
1420 
1421  PseudoTracksLab.set_user_index(i+fJet->GetNumberOfTracks()+100);
1422 
1423  omega2=PseudoTracksLab.perp();
1424  angle2=pTrackCMS.Angle(MyJet);
1425 
1426  fInputVectors.push_back(PseudoTracksLab);
1427  xflagAdded=1;
1428  }
1429 
1430 
1431 
1432 
1433  fastjet::JetAlgorithm jetalgo(fastjet::antikt_algorithm);
1434  if(ReclusterAlgo==0){ xflagalgo=0.5;
1435  jetalgo=fastjet::kt_algorithm ;}
1436 
1437  if(ReclusterAlgo==1){ xflagalgo=1.5;
1438  jetalgo=fastjet::cambridge_algorithm;}
1439  if(ReclusterAlgo==2){ xflagalgo=2.5;
1440  jetalgo=fastjet::antikt_algorithm;}
1441 
1442  fastjet::JetDefinition fJetDef(jetalgo, 1., static_cast<fastjet::RecombinationScheme>(0), fastjet::BestFJ30 );
1443 
1444  try {
1445  fastjet::ClusterSequence fClustSeqSA(fInputVectors, fJetDef);
1446  std::vector<fastjet::PseudoJet> fOutputJets;
1447  fOutputJets.clear();
1448  fOutputJets=fClustSeqSA.inclusive_jets(0);
1449 
1450  fastjet::PseudoJet jj;
1451  fastjet::PseudoJet j1;
1452  fastjet::PseudoJet j2;
1453  jj=fOutputJets[0];
1454  double ndepth=0;
1455  double nsd=0;
1456  double nall=0;
1457  fShapesVar[10]=fOutputJets[0].perp();
1458  while(jj.has_parents(j1,j2)){
1459 
1460  if(j1.perp() < j2.perp()) swap(j1,j2);
1461  double delta_R=j1.delta_R(j2);
1462  double z=j2.perp()/(j1.perp()+j2.perp());
1463  double y =log(1.0/delta_R);
1464  double lnpt_rel=log(z*delta_R);
1465  nall=nall+1;
1466  if(z>fHardCutoff){
1467  ndepth=ndepth+1;
1468  nsd=nsd+1;
1469  if(nsd==1 && ReclusterAlgo==1) fShapesVar[9]=z;
1470  Double_t LundEntries[6] = {y,lnpt_rel,fOutputJets[0].perp(),xflagalgo,partonFlavor,ndepth};
1471  fHLundIterative->Fill(LundEntries);}
1472  jj=j1;}
1473 
1474  if(ReclusterAlgo==1){
1475  fShapesVar[6]=nsd;
1476  fShapesVar[7]=nall;}
1477 
1478  if(fAdditionalTracks>0 && xflagAdded>0){
1479  zinject=omega2/fOutputJets[0].perp();
1480  angleinject=angle2;
1481 
1482  yinject =log(1.0/angleinject);
1483  lnpt_relinject=log(zinject*angleinject);
1484  Double_t LundEntriesInject[4] = {yinject,lnpt_relinject,fOutputJets[0].perp(),fJet->Pt()};
1485  fHLundIterativeInject->Fill(LundEntriesInject);}
1486 
1487 
1488  } catch (fastjet::Error) {
1489  AliError(" [w] FJ Exception caught.");
1490  //return -1;
1491  }
1492 
1493 
1494 
1495  if(fTf1Kt){ delete fTf1Kt;}
1496  if(fTf1Omega){ delete fTf1Omega;}
1497 
1498 
1499  return;
1500 
1501 
1502 }
1503 
1504 
1506  const Double_t jetradius,
1507  const Double_t subjetradius,
1508  const char *ntracksPartLevel,
1509  const char *type,
1510  const char *CentEst,
1511  Int_t pSel,
1512  TString trigClass,
1513  TString kEmcalTriggers,
1514  TString tag,
1515  const char *rhoName,
1519  Float_t minpTHTrigger, Float_t maxpTHTrigger,
1521 
1522 
1523 
1524  AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
1525  if (!mgr)
1526  {
1527  ::Error("AliAnalysisTaskEmcalJetShapesMC","No analysis manager found.");
1528  return 0;
1529  }
1530  Bool_t ismc=kFALSE;
1531  ismc = (mgr->GetMCtruthEventHandler())?kTRUE:kFALSE;
1532 
1533  // Check the analysis type using the event handlers connected to the analysis manager.
1534  //==============================================================================
1535  if (!mgr->GetInputEventHandler())
1536  {
1537  ::Error("AddTaskJetShapesMC", "This task requires an input event handler");
1538  return NULL;
1539  }
1540 
1541  TString wagonName1 = Form("JetShapesMC_%s_Histos%s%s",njetsBase,trigClass.Data(),tag.Data());
1542  TString wagonName2 = Form("JetShapesMC_%s_Tree%s%s",njetsBase,trigClass.Data(),tag.Data());
1543 
1544  //Configure jet tagger task
1546 
1547  //task->SetNCentBins(4);
1548  task->SetJetShapeType(jetShapeType);
1549  task->SetJetShapeSub(jetShapeSub);
1550  task->SetJetSelection(jetSelection);
1551  task->SetDerivativeSubtractionOrder(derivSubtrOrder);
1552  task->SetJetRadius(jetradius);
1553  task->SetSubjetRadius(subjetradius);
1554 
1555  if (jetSelection == AliAnalysisTaskEmcalJetShapesMC::kRecoil) task->SetPtTriggerSelections(minpTHTrigger, maxpTHTrigger);
1556 
1557  TString thename(njetsBase);
1558 
1559  AliParticleContainer *trackContPartLevel = task->AddMCParticleContainer(ntracksPartLevel);
1560 
1561  AliJetContainer *jetContBase=0x0;
1562  TString strType(type);
1563 
1564  jetContBase = task->AddJetContainer(njetsBase,strType,jetradius);
1565  if(jetContBase) {
1566  jetContBase->SetRhoName(rhoName);
1567  jetContBase->ConnectParticleContainer(trackContPartLevel);
1568  //jetContBase->ConnectClusterContainer(clusterCont);
1569  jetContBase->SetPercAreaCut(0.6);
1570  }
1571 
1572  task->SetCaloTriggerPatchInfoName(kEmcalTriggers.Data());
1573  task->SetCentralityEstimator(CentEst);
1574  task->SelectCollisionCandidates(pSel);
1575  task->SetUseAliAnaUtils(kFALSE);
1576 
1577  mgr->AddTask(task);
1578 
1579  //Connnect input
1580  mgr->ConnectInput (task, 0, mgr->GetCommonInputContainer() );
1581 
1582  //Connect output
1583  TString contName1(wagonName1);
1584  TString contName2(wagonName2);
1585 
1586  if (jetShapeType == AliAnalysisTaskEmcalJetShapesMC::kGenShapes) {
1587  contName1 += "_GenShapes";
1588  contName2 += "_GenShapes";
1589  }
1590 
1591  switch (jetShapeSub) {
1592 
1594  contName1 += "_NoSub";
1595  contName2 += "_NoSub";
1596  break;
1597 
1599  contName1 += "_ConstSub";
1600  contName2 += "_ConstSub";
1601  break;
1602 
1604  contName1 += "_DerivSub";
1605  contName2 += "_DerivSub";
1606  break;
1607 
1608 
1609  }
1610 
1611  switch (jetSelection) {
1613  contName1 += "_Incl";
1614  contName2 += "_Incl";
1615  break;
1616 
1617 
1619  TString recoilTriggerString = Form("_Recoil_%.0f_%0.f", minpTHTrigger, maxpTHTrigger);
1620  contName1 += recoilTriggerString;
1621  contName2 += recoilTriggerString;
1622 
1623  break;
1624 
1625  }
1626 
1627 
1628  TString outputfile = Form("%s",AliAnalysisManager::GetCommonFileName());
1629 
1630 
1631  AliAnalysisDataContainer *coutput1 = mgr->CreateContainer(contName1.Data(),TList::Class(),AliAnalysisManager::kOutputContainer,outputfile.Data());
1632  AliAnalysisDataContainer *coutput2 = mgr->CreateContainer(contName2.Data(),TTree::Class(),AliAnalysisManager::kOutputContainer,outputfile.Data());
1633 
1634  mgr->ConnectOutput(task,1,coutput1);
1635  mgr->ConnectOutput(task,2,coutput2);
1636 
1637  return task;
1638 
1639 }
1640 
1641 
1642 
void SetRadius(Double_t val)
Float_t GetSigma2(AliEmcalJet *jet, Int_t jetContNb=0)
Double_t GetFirstOrderSubtractedAngularity() const
Float_t GetJetMass(AliEmcalJet *jet, Int_t jetContNb=0)
double Double_t
Definition: External.C:58
void SetJetMinPt(Double_t val)
Double_t GetSecondOrderSubtractedSigma2() const
Definition: External.C:236
Float_t GetPartonEta6() const
AliJetContainer * GetJetContainer(Int_t i=0) const
Float_t GetPartonEta7() const
Double_t GetSecondOrderSubtractedConstituent() const
Double_t Eta() const
Definition: AliEmcalJet.h:121
Double_t FjNSubJettiness(AliEmcalJet *Jet, Int_t JetContNb, Int_t N, Int_t Algorithm, Double_t Beta, Int_t Option, Double_t Beta_SD=0.0, Double_t ZCut=0.1, Int_t SoftDropOn=0)
Double_t Py() const
Definition: AliEmcalJet.h:107
Double_t Phi() const
Definition: AliEmcalJet.h:117
Float_t GetPythiaEventWeight() const
void RecursiveParents(AliEmcalJet *fJet, AliJetContainer *fJetCont, Int_t ReclusterAlgo, Float_t PartonFlavor)
Float_t GetPartonPhi7() const
void SetPercAreaCut(Float_t p)
Float_t CoreFrac(AliEmcalJet *jet, Int_t jetContNb=0)
Float_t GetJetCircularity(AliEmcalJet *jet, Int_t jetContNb=0)
Double_t E() const
Definition: AliEmcalJet.h:119
Float_t GetJetpTD(AliEmcalJet *jet, Int_t jetContNb=0)
Int_t GetPartonFlag7() const
Container for particles within the EMCAL framework.
Float_t GetPartonPhi6() const
Int_t TrackAt(Int_t idx) const
Definition: AliEmcalJet.h:160
UShort_t GetNumberOfTracks() const
Definition: AliEmcalJet.h:139
Double_t Px() const
Definition: AliEmcalJet.h:106
Float_t GetJetCoreFrac(AliEmcalJet *jet, Int_t jetContNb=0)
Double_t Nsubjettiness(AliEmcalJet *pJet, AliJetContainer *pContJets, Double_t dVtx[3], Int_t N, Int_t Algorithm, Double_t Radius, Double_t Beta, Int_t Option=0, Int_t Measure=0, Double_t Beta_SD=0, Double_t ZCut=0.1, Int_t SoftDropOn=0)
AliParticleContainer * GetParticleContainer(Int_t i=0) const
Get particle container attached to this task.
Int_t GetPartonFlag6() const
void SetRhoName(const char *n)
AliParticleContainer * GetParticleContainer() const
Double_t GetFirstOrderSubtractedConstituent() const
void SetRecombSheme(Int_t val)
static AliAnalysisTaskEmcalJetShapesMC * AddTaskJetShapesMC(const char *njetsBase, const Double_t jetradius, const Double_t subjetradius, const char *ntracksPartLevel, const char *type, const char *CentEst, Int_t pSel, TString trigClass="", TString kEmcalTriggers="", TString tag="", const char *rhoName="", AliAnalysisTaskEmcalJetShapesMC::JetShapeType jetShapeType=AliAnalysisTaskEmcalJetShapesMC::kGenShapes, AliAnalysisTaskEmcalJetShapesMC::JetShapeSub jetShapeSub=AliAnalysisTaskEmcalJetShapesMC::kNoSub, AliAnalysisTaskEmcalJetShapesMC::JetSelectionType jetSelection=AliAnalysisTaskEmcalJetShapesMC::kInclusive, Float_t minpTHTrigger=0., Float_t maxpTHTrigger=0., AliAnalysisTaskEmcalJetShapesMC::DerivSubtrOrder derivSubtrOrder=AliAnalysisTaskEmcalJetShapesMC::kSecondOrder)
void SetJetAlgorithm(Int_t val)
int Int_t
Definition: External.C:63
Double_t Eta() const
Double_t Pt() const
unsigned int UInt_t
Definition: External.C:33
Float_t GetJetNumberOfConstituents(AliEmcalJet *jet, Int_t jetContNb=0)
float Float_t
Definition: External.C:68
Float_t fqhat
Random number generator.
Task to store and correlate the MC shapes.
TTree * fTreeObservableTagging
! Tree with tagging variables subtracted MC or true MC or raw
Double_t RelativePhi(Double_t mphi, Double_t vphi)
Bool_t Run()
Run function. This is the core function of the analysis and contains the user code. Therefore users have to implement this function.
Double_t GetSecondOrderSubtractedAngularity() const
Double_t GetSubjetFraction(AliEmcalJet *Jet, Int_t JetContNb, Double_t JetRadius, AliEmcalJetFinder *Reclusterer)
Definition: Option.C:68
AliEmcalJetFinder * Recluster(AliEmcalJet *Jet, Int_t JetContNb, Double_t JetRadius, Double_t SubJetRadius, Double_t SubJetMinPt, Int_t Algorithm, const char *Name)
Double_t fCent
!event centrality
Float_t PTD(AliEmcalJet *jet, Int_t jetContNb=0)
void SetJetMaxEta(Double_t val)
AliEmcalJet * GetNextAcceptJet()
TF1 * fTf1Kt
to generate omega according to BDMPS tail
Bool_t Data(TH1F *h, Double_t *rangefit, Bool_t writefit, Double_t &sgn, Double_t &errsgn, Double_t &bkg, Double_t &errbkg, Double_t &sgnf, Double_t &errsgnf, Double_t &sigmafit, Int_t &status)
Float_t GetJetAngularity(AliEmcalJet *jet, Int_t jetContNb=0)
Int_t SelectTrigger(Float_t minpT, Float_t maxpT)
AliEmcalJet * GetJet(Int_t index)
void ConnectParticleContainer(AliParticleContainer *c)
Double_t Pt() const
Definition: AliEmcalJet.h:109
Enhanced TList-derived class that implements correct merging for pt_hard binned production.
Definition: AliEmcalList.h:25
Float_t Angularity(AliEmcalJet *jet, Int_t jetContNb=0)
Bool_t FillHistograms()
Function filling histograms.
Double_t GetFirstOrderSubtractedCircularity() const
Bool_t RetrieveEventObjects()
Retrieve common objects from event.
AliEmcalList * fOutput
!output list
Double_t SubJetOrdering(AliEmcalJet *Jet, AliEmcalJetFinder *Reclusterer, Int_t N, Int_t Type, Bool_t Index)
void NTValues(AliEmcalJet *jet, Int_t jetContNb, Float_t *nTFractions)
const Int_t nVar
Double_t GetSecondOrderSubtractedCircularity() const
Double_t DeltaR(const AliVParticle *part1, const AliVParticle *part2)
Float_t GetJetLeSub(AliEmcalJet *jet, Int_t jetContNb=0)
Store some informaion about a Pythia eventThis class is used to store some information about a Pythia...
void SetMakeGeneralHistograms(Bool_t g)
Base task in the EMCAL jet framework.
Represent a jet reconstructed using the EMCal jet framework.
Definition: AliEmcalJet.h:51
const AliEmcalPythiaInfo * GetPythiaInfo() const
Double_t Pz() const
Definition: AliEmcalJet.h:108
Float_t Circularity(AliEmcalJet *jet, Int_t jetContNb=0)
void SoftDrop(AliEmcalJet *fJet, AliJetContainer *fJetCont, Double_t zcut, Double_t beta, Int_t ReclusterAlgo)
Declaration of class AliEmcalPythiaInfo.
const char Option_t
Definition: External.C:48
void UserCreateOutputObjects()
Main initialization function on the worker.
Double_t GetFirstOrderSubtractedSigma2() const
Float_t Sigma2(AliEmcalJet *jet, Int_t jetContNb=0)
bool Bool_t
Definition: External.C:53
AliEmcalJetShapeProperties * GetShapeProperties() const
Definition: AliEmcalJet.h:361
Float_t GetPartonPt7() const
Float_t LeSub(AliEmcalJet *jet, Int_t jetContNb=0)
void swap(PWGJE::EMCALJetTasks::AliAnalysisTaskEmcalJetHPerformance &first, PWGJE::EMCALJetTasks::AliAnalysisTaskEmcalJetHPerformance &second)
Double_t M() const
Definition: AliEmcalJet.h:120
Container for jet within the EMCAL jet framework.
std::vector< int > GetPtSortedTrackConstituentIndexes(TClonesArray *tracks) const
Float_t GetPartonPt6() const