AliPhysics  3337bb0 (3337bb0)
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
AliAnalysisTaskEmcalJetSubstructureTree.cxx
Go to the documentation of this file.
1 /************************************************************************************
2  * Copyright (C) 2017, Copyright Holders of the ALICE Collaboration *
3  * All rights reserved. *
4  * *
5  * Redistribution and use in source and binary forms, with or without *
6  * modification, are permitted provided that the following conditions are met: *
7  * * Redistributions of source code must retain the above copyright *
8  * notice, this list of conditions and the following disclaimer. *
9  * * Redistributions in binary form must reproduce the above copyright *
10  * notice, this list of conditions and the following disclaimer in the *
11  * documentation and/or other materials provided with the distribution. *
12  * * Neither the name of the <organization> nor the *
13  * names of its contributors may be used to endorse or promote products *
14  * derived from this software without specific prior written permission. *
15  * *
16  * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND *
17  * ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED *
18  * WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE *
19  * DISCLAIMED. IN NO EVENT SHALL ALICE COLLABORATION BE LIABLE FOR ANY *
20  * DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES *
21  * (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; *
22  * LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND *
23  * ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT *
24  * (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS *
25  * SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. *
26  ************************************************************************************/
27 #include <iostream>
28 #include <string>
29 #include <sstream>
30 #include <vector>
31 
32 #include <fastjet/ClusterSequence.hh>
33 #include <fastjet/contrib/Nsubjettiness.hh>
34 #include <fastjet/contrib/SoftDrop.hh>
35 
36 #include <TLorentzVector.h>
37 #include <TMath.h>
38 #include <TString.h>
39 #include <TVector3.h>
40 
41 #include "AliAODInputHandler.h"
42 #include "AliAnalysisManager.h"
44 #include "AliClusterContainer.h"
45 #include "AliJetContainer.h"
47 #include "AliEmcalJet.h"
48 #include "AliEmcalList.h"
49 #include "AliLog.h"
50 #include "AliParticleContainer.h"
51 #include "AliTrackContainer.h"
52 #include "AliRhoParameter.h"
53 #include "AliVCluster.h"
54 #include "AliVParticle.h"
55 
59 
60 namespace EmcalTriggerJets {
61 
64  fJetSubstructureTree(nullptr),
65  fSDZCut(0.1),
66  fSDBetaCut(0),
67  fReclusterizer(kCAAlgo),
68  fTriggerSelectionBits(AliVEvent::kAny),
69  fTriggerSelectionString("")
70 {
71  memset(fJetTreeData, 0, sizeof(Double_t) * kTNVar);
72 }
73 
75  AliAnalysisTaskEmcalJet(name, kTRUE),
76  fJetSubstructureTree(nullptr),
77  fSDZCut(0.1),
78  fSDBetaCut(0),
79  fReclusterizer(kCAAlgo),
80  fTriggerSelectionBits(AliVEvent::kAny),
81  fTriggerSelectionString("")
82 {
83  memset(fJetTreeData, 0, sizeof(Double_t) * kTNVar);
84  DefineOutput(2, TTree::Class());
85 }
86 
88 
89 }
90 
93 
94  OpenFile(2);
95  fJetSubstructureTree = new TTree("jetSubstructure", "Tree with jet substructure information");
96  TString varnames[kTNVar];
97  varnames[0] = "Radius";
98  varnames[1] = "EventWeight";
99  varnames[2] = "PtJetRec";
100  varnames[3] = "PtJetSim";
101  varnames[4] = "EJetRec";
102  varnames[5] = "EJetSim";
103  varnames[6] = "RhoPtRec";
104  varnames[7] = "RhoPtSim";
105  varnames[8] = "RhoMassRec";
106  varnames[9] = "RhoMassSim";
107  varnames[10] = "AreaRec";
108  varnames[11] = "AreaSim";
109  varnames[12] = "NEFRec";
110  varnames[13] = "NEFSim";
111  varnames[14] = "MassRec";
112  varnames[15] = "MassSim";
113  varnames[16] = "ZgMeasured";
114  varnames[17] = "ZgTrue";
115  varnames[18] = "RgMeasured";
116  varnames[19] = "RgTrue";
117  varnames[20] = "MgMeasured";
118  varnames[21] = "MgTrue";
119  varnames[22] = "PtgMeasured";
120  varnames[23] = "PtgTrue";
121  varnames[24] = "MugMeasured";
122  varnames[25] = "MugTrue";
123  varnames[26] = "OneSubjettinessMeasured";
124  varnames[27] = "OneSubjettinessTrue";
125  varnames[28] = "TwoSubjettinessMeasured";
126  varnames[29] = "TwoSubjettinessTrue";
127  varnames[30] = "AngularityMeasured";
128  varnames[31] = "AngularityTrue";
129  varnames[32] = "PtDMeasured";
130  varnames[33] = "PtDTrue";
131  varnames[34] = "NCharged";
132  varnames[35] = "NNeutral";
133  varnames[36] = "NConstTrue";
134  varnames[37] = "NDroppedMeasured";
135  varnames[38] = "NDroppedTrue";
136 
137  for(int ib = 0; ib < kTNVar; ib++){
138  fJetSubstructureTree->Branch(varnames[ib], fJetTreeData + ib, Form("%s/D", varnames[ib].Data()));
139  }
140  PostData(1, fOutput);
141  PostData(2, fJetSubstructureTree);
142 }
143 
145  AliClusterContainer *clusters = GetClusterContainer("caloClusters");
146  AliTrackContainer *tracks = GetTrackContainer("tracks");
147  AliParticleContainer *particles = GetParticleContainer("mcparticles");
148 
149  AliJetContainer *mcjets = GetJetContainer("mcjets");
150  AliJetContainer *datajets = GetJetContainer("datajets");
151 
152  TString rhoTagData = datajets ? TString::Format("R%02d", static_cast<Int_t>(datajets->GetJetRadius() * 10.)) : "",
153  rhoTagMC = mcjets ? TString::Format("R%02d", static_cast<Int_t>(mcjets->GetJetRadius() * 10.)) : "";
154 
155  AliRhoParameter *rhoPtRec = GetRhoFromEvent("RhoSparse_Full_" + rhoTagData),
156  *rhoMassRec = GetRhoFromEvent("RhoMassSparse_Full_" + rhoTagData),
157  *rhoPtSim = GetRhoFromEvent("RhoSparse_Full_" + rhoTagMC),
158  *rhoMassSim = GetRhoFromEvent("RhoMassSparse_Full_" + rhoTagMC);
159  AliDebugStream(2) << "Found rho parameter for reconstructed pt: " << (rhoPtRec ? "yes" : "no") << ", value: " << (rhoPtRec ? rhoPtRec->GetVal() : 0.) << std::endl;
160  AliDebugStream(2) << "Found rho parameter for sim pt: " << (rhoPtSim ? "yes" : "no") << ", value: " << (rhoPtSim ? rhoPtSim->GetVal() : 0.) << std::endl;
161  AliDebugStream(2) << "Found rho parameter for reconstructed Mass: " << (rhoMassRec ? "yes" : "no") << ", value: " << (rhoMassRec ? rhoMassRec->GetVal() : 0.) << std::endl;
162  AliDebugStream(2) << "Found rho parameter for sim Mass: " << (rhoMassSim ? "yes" : "no") << ", value: " << (rhoMassSim ? rhoMassSim->GetVal() : 0.) << std::endl;
163 
164  // Run trigger selection (not on pure MCgen train)
165  if(datajets){
166  if(!(fInputHandler->IsEventSelected() & fTriggerSelectionBits)) return false;
167  if(fTriggerSelectionString.Length()) {
168  if(!fInputEvent->GetFiredTriggerClasses().Contains(fTriggerSelectionString)) return false;
169  }
170  }
171 
172  Double_t weight = 1;
173  Double_t rhoparameters[4]; memset(rhoparameters, 0, sizeof(Double_t) * 4);
174  if(rhoPtRec) rhoparameters[0] = rhoPtRec->GetVal();
175  if(rhoPtSim) rhoparameters[1] = rhoPtSim->GetVal();
176  if(rhoMassRec) rhoparameters[2] = rhoMassRec->GetVal();
177  if(rhoMassSim) rhoparameters[3] = rhoMassSim->GetVal();
178 
179  AliSoftdropDefinition softdropSettings;
180  softdropSettings.fBeta = fSDBetaCut;
181  softdropSettings.fZ = fSDZCut;
182  switch(fReclusterizer) {
183  case kCAAlgo: softdropSettings.fRecluserAlgo = fastjet::cambridge_aachen_algorithm; break;
184  case kKTAlgo: softdropSettings.fRecluserAlgo = fastjet::kt_algorithm; break;
185  case kAKTAlgo: softdropSettings.fRecluserAlgo = fastjet::antikt_algorithm; break;
186  };
187 
188  AliNSubjettinessDefinition nsubjettinessSettings;
189  nsubjettinessSettings.fBeta = 1.;
190  nsubjettinessSettings.fRadius = 0.4;
191 
192  if(mcjets && !datajets) {
193  // pure MC (gen) train - run over MC jets
194  AliDebugStream(1) << "In MC pure jet branch: found " << mcjets->GetNJets() << " jets, " << mcjets->GetNAcceptedJets() << " were accepted\n";
195  for(auto jet : mcjets->accepted()) {
196  try {
197  AliJetSubstructureData structure = MakeJetSubstructure(*jet, mcjets->GetJetRadius() * 2., particles, nullptr,{softdropSettings, nsubjettinessSettings});
198  Double_t angularity[2] = {0., MakeAngularity(*jet, particles, nullptr)},
199  ptd[2] = {0., MakePtD(*jet, particles, nullptr)};
200  FillTree(mcjets->GetJetRadius(), weight, nullptr, jet, nullptr, &(structure.fSoftDrop), nullptr, &(structure.fNsubjettiness), angularity, ptd, rhoparameters);
201  } catch (ReclusterizerException &e) {
202  AliErrorStream() << "Error in reclusterization - skipping jet" << std::endl;
203  }
204  }
205  }
206 
207 
208  if(datajets) {
209  AliDebugStream(1) << "In data jets branch: found" << datajets->GetNJets() << " jets, " << datajets->GetNAcceptedJets() << " were accepted\n";
210  AliDebugStream(1) << "Having MC information: " << (mcjets ? TString::Format("yes, with %d jets", mcjets->GetNJets()) : "no") << std::endl;
211  for(auto jet : datajets->accepted()) {
212  double pt = jet->Pt(), pz = jet->Pz(), E = jet->E(), M = TMath::Sqrt(E*E - pt*pt - pz*pz);
213  AliDebugStream(2) << "Next jet: pt:" << jet->Pt() << ", E: " << E << ", pz: " << pz << ", M(self): " << M << "M(fj)" << jet->M() << std::endl;
214  AliEmcalJet *associatedJet = jet->ClosestJet();
215  if(mcjets) {
216  if(!associatedJet) continue;
217  try {
218  AliJetSubstructureData structureData = MakeJetSubstructure(*jet, datajets->GetJetRadius() * 2., tracks, clusters, {softdropSettings, nsubjettinessSettings}),
219  structureMC = MakeJetSubstructure(*associatedJet, mcjets->GetJetRadius() * 2, particles, nullptr, {softdropSettings, nsubjettinessSettings});
220  Double_t angularity[2] = {MakeAngularity(*jet, tracks, clusters), MakeAngularity(*associatedJet, particles, nullptr)},
221  ptd[2] = {MakePtD(*jet, tracks, clusters), MakePtD(*associatedJet, particles, nullptr)};
222  FillTree(datajets->GetJetRadius(), weight, jet, associatedJet, &(structureData.fSoftDrop), &(structureMC.fSoftDrop), &(structureData.fNsubjettiness), &(structureMC.fNsubjettiness), angularity, ptd, rhoparameters);
223  } catch(ReclusterizerException &e) {
224  AliErrorStream() << "Error in reclusterization - skipping jet" << std::endl;
225  }
226  } else {
227  try {
228  AliJetSubstructureData structure = MakeJetSubstructure(*jet, 0.4, tracks, clusters, {softdropSettings, nsubjettinessSettings});
229  Double_t angularity[2] = {MakeAngularity(*jet, tracks, clusters), 0.},
230  ptd[2] = {MakePtD(*jet, tracks, clusters), 0.};
231  FillTree(datajets->GetJetRadius(), weight, jet, nullptr, &(structure.fSoftDrop), nullptr, &(structure.fNsubjettiness), nullptr, angularity, ptd, rhoparameters);
232  } catch(ReclusterizerException &e) {
233  AliErrorStream() << "Error in reclusterization - skipping jet" << std::endl;
234  }
235  }
236  }
237  }
238 
239  return true;
240 }
241 
243  const AliEmcalJet *datajet, const AliEmcalJet *mcjet,
244  AliSoftDropParameters *dataSoftdrop, AliSoftDropParameters *mcSoftdrop,
245  AliNSubjettinessParameters *dataSubjettiness, AliNSubjettinessParameters *mcSubjettiness,
246  Double_t *angularity, Double_t *ptd, Double_t *rhoparameters){
247  fJetTreeData[kTRadius] = r;
248  fJetTreeData[kTWeight] = weight;
249  fJetTreeData[kTRhoPtRec] = rhoparameters[0];
250  fJetTreeData[kTRhoPtSim] = rhoparameters[1];
251  fJetTreeData[kTRhoMassRec] = rhoparameters[2];
252  fJetTreeData[kTRhoMassSim] = rhoparameters[3];
253  if(datajet) {
254  fJetTreeData[kTPtJetRec] = TMath::Abs(datajet->Pt());
256  fJetTreeData[kTNNeutral] = datajet->GetNumberOfClusters();
257  fJetTreeData[kTAreaRec] = datajet->Area();
258  fJetTreeData[kTNEFRec] = datajet->NEF();
259  fJetTreeData[kTMassRec] = datajet->M();
260  fJetTreeData[kTEJetRec] = datajet->E();
261  } else {
262  fJetTreeData[kTPtJetRec] = 0.;
265  fJetTreeData[kTAreaRec] = 0.;
266  fJetTreeData[kTNEFRec] = 0.;
267  fJetTreeData[kTMassRec] = 0.;
268  fJetTreeData[kTEJetRec] = 0.;
269  }
270 
271  if(mcjet) {
272  fJetTreeData[kTPtJetSim] = TMath::Abs(mcjet->Pt());
274  fJetTreeData[kTAreaSim] = mcjet->Area();
275  fJetTreeData[kTNEFSim] = mcjet->NEF();
276  fJetTreeData[kTMassSim] = mcjet->M();
277  fJetTreeData[kTEJetSim] = mcjet->E();
278  } else {
279  fJetTreeData[kTPtJetSim] = 0.;
281  fJetTreeData[kTAreaSim] = 0.;
282  fJetTreeData[kTNEFSim] = 0.;
283  fJetTreeData[kTMassSim] = 0.;
284  fJetTreeData[kTEJetSim] = 0.;
285  }
286 
287  if(dataSoftdrop) {
288  fJetTreeData[kTZgMeasured] = dataSoftdrop->fZg;
289  fJetTreeData[kTRgMeasured] = dataSoftdrop->fRg;
290  fJetTreeData[kTMgMeasured] = dataSoftdrop->fMg;
291  fJetTreeData[kTPtgMeasured] = dataSoftdrop->fPtg;
292  fJetTreeData[kTMugMeasured] = dataSoftdrop->fMug;
293  fJetTreeData[kTNDroppedMeasured] = dataSoftdrop->fNDropped;
294  } else {
301  }
302 
303  if(mcSoftdrop) {
304  fJetTreeData[kTZgTrue] = mcSoftdrop->fZg;
305  fJetTreeData[kTRgTrue] = mcSoftdrop->fRg;
306  fJetTreeData[kTMgTrue] = mcSoftdrop->fMg;
307  fJetTreeData[kTPtgTrue] = mcSoftdrop->fPtg;
308  fJetTreeData[kTMugTrue] = mcSoftdrop->fMug;
309  fJetTreeData[kTNDroppedTrue] = mcSoftdrop->fNDropped;
310  } else {
311  fJetTreeData[kTZgTrue] = 0.;
312  fJetTreeData[kTRgTrue] = 0.;
313  fJetTreeData[kTMgTrue] = 0.;
314  fJetTreeData[kTPtgTrue] = 0.;
315  fJetTreeData[kTMugTrue] = 0.;
317  }
318 
319  if(dataSubjettiness) {
322  } else {
325  }
326 
327  if(mcSubjettiness) {
330  } else {
333  }
334 
335  fJetTreeData[kTAngularityMeasured] = angularity[0];
336  fJetTreeData[kTAngularityTrue] = angularity[1];
337  fJetTreeData[kTPtDMeasured] = ptd[0];
338  fJetTreeData[kTPtDTrue] = ptd[1];
339 
340  fJetSubstructureTree->Fill();
341 }
342 
343 
345  const int kClusterOffset = 30000; // In order to handle tracks and clusters in the same index space the cluster index needs and offset, large enough so that there is no overlap with track indices
346  std::vector<fastjet::PseudoJet> constituents;
347  for(int itrk = 0; itrk < jet.GetNumberOfTracks(); itrk++){
348  AliVTrack *track = static_cast<AliVTrack *>(jet.TrackAt(itrk, tracks->GetArray()));
349  fastjet::PseudoJet constituentTrack(track->Px(), track->Py(), track->Pz(), track->E());
350  constituentTrack.set_user_index(jet.TrackAt(itrk));
351  constituents.push_back(constituentTrack);
352  }
353 
354  if(clusters){
355  for(int icl = 0; icl < jet.GetNumberOfClusters(); icl++) {
356  AliVCluster *cluster = jet.ClusterAt(icl, clusters->GetArray());
357  TLorentzVector clustervec;
358  cluster->GetMomentum(clustervec, fVertex, (AliVCluster::VCluUserDefEnergy_t)clusters->GetDefaultClusterEnergy());
359  fastjet::PseudoJet constituentCluster(clustervec.Px(), clustervec.Py(), clustervec.Pz(), cluster->GetHadCorrEnergy());
360  constituentCluster.set_user_index(jet.ClusterAt(icl) + kClusterOffset);
361  constituents.push_back(constituentCluster);
362  }
363  }
364 
365  // Redo jet finding on constituents with a
366  fastjet::JetDefinition jetdef(fastjet::antikt_algorithm, jetradius*2, static_cast<fastjet::RecombinationScheme>(0), fastjet::BestFJ30 );
367  std::vector<fastjet::PseudoJet> outputjets;
368  try {
369  fastjet::ClusterSequence jetfinder(constituents, jetdef);
370  outputjets = jetfinder.inclusive_jets(0);
371  AliJetSubstructureData result({MakeSoftDropParameters(outputjets[0], settings.fSoftdropSettings), MakeNsubjettinessParameters(outputjets[0], settings.fSubjettinessSettings)});
372  return result;
373  } catch (fastjet::Error &e) {
374  AliErrorStream() << " FJ Exception caught: " << e.message() << std::endl;
375  throw ReclusterizerException();
376  }
377 }
378 
380  fastjet::contrib::SoftDrop softdropAlgorithm(cutparameters.fBeta, cutparameters.fZ);
381  softdropAlgorithm.set_verbose_structure(kTRUE);
382  std::unique_ptr<fastjet::contrib::Recluster> reclusterizer(new fastjet::contrib::Recluster(cutparameters.fRecluserAlgo, 1, true));
383  softdropAlgorithm.set_reclustering(kTRUE, reclusterizer.get());
384  fastjet::PseudoJet groomed = softdropAlgorithm(jet);
385 
386  AliSoftDropParameters result({groomed.structure_of<fastjet::contrib::SoftDrop>().symmetry(),
387  groomed.m(),
388  groomed.structure_of<fastjet::contrib::SoftDrop>().delta_R(),
389  groomed.perp(),
390  groomed.structure_of<fastjet::contrib::SoftDrop>().mu(),
391  groomed.structure_of<fastjet::contrib::SoftDrop>().dropped_count()});
392  return result;
393 }
394 
397  fastjet::contrib::Nsubjettiness (1,fastjet::contrib::KT_Axes(),fastjet::contrib::NormalizedCutoffMeasure(cut.fBeta, cut.fRadius, 1e100)).result(jet),
398  fastjet::contrib::Nsubjettiness (2,fastjet::contrib::KT_Axes(),fastjet::contrib::NormalizedCutoffMeasure(cut.fBeta, cut.fRadius, 1e100)).result(jet)
399  });
400  return result;
401 }
402 
404  if(!jet.GetNumberOfTracks()) return 0;
405  TVector3 jetvec(jet.Px(), jet.Py(), jet.Pz());
406  Double_t den(0.), num(0.);
407  if(tracks){
408  for(UInt_t itrk = 0; itrk < jet.GetNumberOfTracks(); itrk++) {
409  AliVParticle *track = jet.TrackAt(itrk, tracks->GetArray());
410  if(!track){
411  AliErrorStream() << "Associated constituent particle / track not found\n";
412  continue;
413  }
414  TVector3 trackvec(track->Px(), track->Py(), track->Pz());
415 
416  num += track->Pt() * trackvec.DrEtaPhi(jetvec);
417  den += +track->Pt();
418  }
419  }
420  if(clusters) {
421  for(UInt_t icl = 0; icl < jet.GetNumberOfClusters(); icl++){
422  AliVCluster *clust = jet.ClusterAt(icl, clusters->GetArray());
423  if(!clust) {
424  AliErrorStream() << "Associated constituent cluster not found\n";
425  continue;
426  }
427  TLorentzVector clusterp;
428  clust->GetMomentum(clusterp, fVertex, (AliVCluster::VCluUserDefEnergy_t)clusters->GetDefaultClusterEnergy());
429 
430  num += clusterp.Pt() * clusterp.Vect().DrEtaPhi(jetvec);
431  den += clusterp.Pt();
432  }
433  }
434  return num/den;
435 }
436 
438  if (!jet.GetNumberOfTracks()) return 0;
439  Double_t den(0.), num(0.);
440  if(particles){
441  for(UInt_t itrk = 0; itrk < jet.GetNumberOfTracks(); itrk++) {
442  AliVParticle *trk = jet.TrackAt(itrk, particles->GetArray());
443  if(!trk){
444  AliErrorStream() << "Associated constituent particle / track not found\n";
445  continue;
446  }
447  num += trk->Pt() * trk->Pt();
448  den += trk->Pt();
449  }
450  }
451  if(clusters){
452  for(UInt_t icl = 0; icl < jet.GetNumberOfClusters(); icl++){
453  AliVCluster *clust = jet.ClusterAt(icl, clusters->GetArray());
454  if(!clust) {
455  AliErrorStream() << "Associated constituent cluster not found\n";
456  continue;
457  }
458  TLorentzVector clusterp;
459  clust->GetMomentum(clusterp, fVertex, (AliVCluster::VCluUserDefEnergy_t)clusters->GetDefaultClusterEnergy());
460  num += clusterp.Pt() * clusterp.Pt();
461  den += clusterp.Pt();
462  }
463  }
464  return TMath::Sqrt(num)/den;
465 }
466 
468  AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
469 
470  Bool_t isAOD(kFALSE);
471  AliInputEventHandler *inputhandler = static_cast<AliInputEventHandler *>(mgr->GetInputEventHandler());
472  if(inputhandler) {
473  if(inputhandler->IsA() == AliAODInputHandler::Class()) isAOD = kTRUE;
474  }
475 
476  AliAnalysisTaskEmcalJetSubstructureTree *treemaker = new AliAnalysisTaskEmcalJetSubstructureTree("JetSubstructureTreemaker_" + TString::Format("R%02d_", int(jetradius * 10.)) + trigger);
477  mgr->AddTask(treemaker);
478  treemaker->SetMakeGeneralHistograms(kTRUE);
479 
480  // Adding containers
481  if(isMC) {
482  AliParticleContainer *particles = treemaker->AddMCParticleContainer("mcparticles");
483  particles->SetMinPt(0.);
484 
485  AliJetContainer *mcjets = treemaker->AddJetContainer(
488  recombinationScheme,
489  jetradius,
491  particles, nullptr);
492  mcjets->SetName("mcjets");
493  mcjets->SetJetPtCut(20.);
494  }
495 
496  if(isData) {
498  tracks->SetMinPt(0.15);
499  AliClusterContainer *clusters(nullptr);
500  if(jettype == kFull){
501  std::cout << "Using full jets ..." << std::endl;
503  clusters->SetClusHadCorrEnergyCut(0.3); // 300 MeV E-cut
504  clusters->SetDefaultClusterEnergy(AliVCluster::kHadCorr);
505  } else {
506  std::cout << "Using charged jets ... " << std::endl;
507  }
508 
509  AliJetContainer *datajets = treemaker->AddJetContainer(
512  recombinationScheme,
513  jetradius,
515  tracks, clusters);
516  datajets->SetName("datajets");
517  datajets->SetJetPtCut(20.);
518 
519  treemaker->SetUseAliAnaUtils(true, true);
520  treemaker->SetVzRange(-10., 10);
521 
522  // configure trigger selection
523  TString triggerstring(trigger);
524  if(triggerstring.Contains("INT7")) {
525  treemaker->SetTriggerBits(AliVEvent::kINT7);
526  } else if(triggerstring.Contains("EJ1")) {
527  treemaker->SetTriggerBits(AliVEvent::kEMCEJE);
528  treemaker->SetTriggerString("EJ1");
529  } else if(triggerstring.Contains("EJ2")) {
530  treemaker->SetTriggerBits(AliVEvent::kEMCEJE);
531  treemaker->SetTriggerString("EJ2");
532  }
533  }
534 
535  // Connecting containers
536  TString outputfile = mgr->GetCommonFileName();
537  outputfile += TString::Format(":JetSubstructure_R%02d_%s", int(jetradius * 10.), trigger);
538  mgr->ConnectInput(treemaker, 0, mgr->GetCommonInputContainer());
539  mgr->ConnectOutput(treemaker, 1, mgr->CreateContainer("JetSubstructureHistos_" + TString::Format("R%0d_", int(jetradius * 10.)) + trigger, AliEmcalList::Class(), AliAnalysisManager::kOutputContainer, outputfile));
540  mgr->ConnectOutput(treemaker, 2, mgr->CreateContainer("JetSubstuctureTree_" + TString::Format("R%0d_", int(jetradius * 10.)) + trigger, TTree::Class(), AliAnalysisManager::kOutputContainer, Form("JetSubstructureTree_R%02d_%s.root", static_cast<int>(jetradius*10.), trigger)));
541 
542  return treemaker;
543 }
544 
545 } /* namespace EmcalTriggerJets */
Double_t Area() const
Definition: AliEmcalJet.h:123
double Double_t
Definition: External.C:58
Double_t MakePtD(const AliEmcalJet &jet, const AliParticleContainer *const particles, const AliClusterContainer *const clusters) const
AliEmcalJet * ClosestJet() const
Definition: AliEmcalJet.h:222
AliJetContainer * GetJetContainer(Int_t i=0) const
Double_t Py() const
Definition: AliEmcalJet.h:100
Container with name, TClonesArray and cuts for particles.
void SetUseAliAnaUtils(Bool_t b, Bool_t bRejPilup=kTRUE)
AliJetSubstructureData MakeJetSubstructure(const AliEmcalJet &jet, double jetradius, const AliParticleContainer *tracks, const AliClusterContainer *clusters, const AliJetSubstructureSettings &settings) const
Int_t ClusterAt(Int_t idx) const
Definition: AliEmcalJet.h:130
AliJetContainer * AddJetContainer(const char *n, TString defaultCutType, Float_t jetRadius=0.4)
Structure for results from the soft drop algorithm.
Double_t E() const
Definition: AliEmcalJet.h:112
virtual bool Run()
Run function. This is the core function of the analysis and contains the user code. Therefore users have to implement this function.
void SetVzRange(Double_t min, Double_t max)
AliClusterContainer * AddClusterContainer(const char *n)
Create new cluster container and attach it to the task.
fastjet::JetAlgorithm fRecluserAlgo
Reclusterization algorithm.
UShort_t GetNumberOfConstituents() const
Definition: AliEmcalJet.h:133
Container for particles within the EMCAL framework.
Int_t GetDefaultClusterEnergy() const
Int_t TrackAt(Int_t idx) const
Definition: AliEmcalJet.h:153
UShort_t GetNumberOfTracks() const
Definition: AliEmcalJet.h:132
Double_t Px() const
Definition: AliEmcalJet.h:99
AliParticleContainer * GetParticleContainer(Int_t i=0) const
Get particle container attached to this task.
static AliAnalysisTaskEmcalJetSubstructureTree * AddEmcalJetSubstructureTreeMaker(Bool_t isMC, Bool_t isData, Double_t jetradius, JetType_t jettype, AliJetContainer::ERecoScheme_t recombinationScheme, const char *name)
AliRhoParameter * GetRhoFromEvent(const char *name)
UShort_t GetNumberOfClusters() const
Definition: AliEmcalJet.h:131
void SetJetPtCut(Float_t cut)
TPC acceptance.
Definition: AliEmcalJet.h:60
unsigned int UInt_t
Definition: External.C:33
Int_t GetNJets() const
AliNSubjettinessParameters MakeNsubjettinessParameters(const fastjet::PseudoJet &jet, const AliNSubjettinessDefinition &cut) const
Definition for the algorithm obtaining the softdrop parameters.
AliClusterContainer * GetClusterContainer(Int_t i=0) const
Get cluster container attached to this task.
AliMCParticleContainer * AddMCParticleContainer(const char *n)
Create new container for MC particles and attach it to the task.
Double_t MakeAngularity(const AliEmcalJet &jet, const AliParticleContainer *tracks, const AliClusterContainer *clusters) const
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)
Double_t Pt() const
Definition: AliEmcalJet.h:102
Bool_t isMC
Float_t GetJetRadius() const
AliEmcalList * fOutput
!output list
Double_t fVertex[3]
!event vertex
AliTrackContainer * AddTrackContainer(const char *n)
Create new track container and attach it to the task.
AliTrackContainer * GetTrackContainer(Int_t i=0) const
void SetMakeGeneralHistograms(Bool_t g)
Base task in the EMCAL jet framework.
void FillTree(double r, double weight, const AliEmcalJet *datajet, const AliEmcalJet *mcjet, AliSoftDropParameters *dataSoftdrop, AliSoftDropParameters *mcsoftdrop, AliNSubjettinessParameters *dataSubjettiness, AliNSubjettinessParameters *mcSubjettiness, Double_t *angularity, Double_t *ptd, Double_t *rhoparameters)
Represent a jet reconstructed using the EMCal jet framework.
Definition: AliEmcalJet.h:44
Double_t Pz() const
Definition: AliEmcalJet.h:101
AliSoftDropParameters MakeSoftDropParameters(const fastjet::PseudoJet &jet, const AliSoftdropDefinition &cut) const
void UserCreateOutputObjects()
Main initialization function on the worker.
const AliJetIterableContainer accepted() const
bool Bool_t
Definition: External.C:53
Double_t NEF() const
Definition: AliEmcalJet.h:141
void SetDefaultClusterEnergy(Int_t d)
Container structure for EMCAL clusters.
Double_t M() const
Definition: AliEmcalJet.h:113
EMCal fiducial acceptance (each eta, phi edge narrowed by jet R)
Definition: AliEmcalJet.h:63
Container for jet within the EMCAL jet framework.
TList * OpenFile(const char *fname)
Definition: DrawAnaELoss.C:65
void SetClusHadCorrEnergyCut(Double_t cut)
static TString ClusterContainerNameFactory(Bool_t isAOD)
Get name of the default cluster container.
static TString TrackContainerNameFactory(Bool_t isAOD)
Get name of the default track container.