AliPhysics  a3be53f (a3be53f)
AliAnalysisTaskEmcalSoftDropResponse.cxx
Go to the documentation of this file.
1 /************************************************************************************
2  * Copyright (C) 2019, 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 <vector>
29 
30 #include <TArrayD.h>
31 #include <TBinning.h>
32 #include <TCustomBinning.h>
33 #include <TH1.h>
34 #include <TH2.h>
35 #include <TLorentzVector.h>
36 #include <TRandom.h>
37 
38 #include <fastjet/ClusterSequence.hh>
39 #include <fastjet/PseudoJet.hh>
40 #include <fastjet/contrib/SoftDrop.hh>
41 
42 #include <RooUnfoldResponse.h>
43 
44 #include <AliAODEvent.h>
45 #include <AliAODInputHandler.h>
46 #include <AliAnalysisManager.h>
48 #include <AliClusterContainer.h>
49 #include <AliEmcalJet.h>
51 #include <AliJetContainer.h>
52 #include <AliMCParticleContainer.h>
53 #include <AliTrackContainer.h>
54 #include <AliVCluster.h>
55 #include <AliVParticle.h>
56 #include <AliVTrack.h>
57 
59 
60 using namespace PWGJE::EMCALJetTasks;
61 
64  fBinningMode(kSDModeINT7),
65  fFractionResponseClosure(0.5),
66  fZcut(0.1),
67  fBeta(0.),
68  fReclusterizer(kCAAlgo),
69  fUseChargedConstituents(true),
70  fUseNeutralConstituents(true),
71  fSampleSplitter(nullptr),
72  fPartLevelPtBinning(nullptr),
73  fDetLevelPtBinning(nullptr),
74  fZgResponse(nullptr),
75  fZgResponseClosure(nullptr),
76  fZgPartLevel(nullptr),
77  fZgDetLevel(nullptr),
78  fZgPartLevelTruncated(nullptr),
79  fZgPartLevelClosureNoResp(nullptr),
80  fZgDetLevelClosureNoResp(nullptr),
81  fZgPartLevelClosureResp(nullptr),
82  fZgDetLevelClosureResp(nullptr)
83 {
84 
85 }
86 
87 AliAnalysisTaskEmcalSoftDropResponse::AliAnalysisTaskEmcalSoftDropResponse(const char *name):
88  AliAnalysisTaskEmcalJet(name, kTRUE),
89  fBinningMode(kSDModeINT7),
90  fFractionResponseClosure(0.5),
91  fZcut(0.1),
92  fBeta(0.),
93  fReclusterizer(kCAAlgo),
94  fUseChargedConstituents(true),
95  fUseNeutralConstituents(true),
96  fSampleSplitter(nullptr),
97  fPartLevelPtBinning(nullptr),
98  fDetLevelPtBinning(nullptr),
99  fZgResponse(nullptr),
100  fZgResponseClosure(nullptr),
101  fZgPartLevel(nullptr),
102  fZgDetLevel(nullptr),
103  fZgPartLevelTruncated(nullptr),
104  fZgPartLevelClosureNoResp(nullptr),
105  fZgDetLevelClosureNoResp(nullptr),
106  fZgPartLevelClosureResp(nullptr),
107  fZgDetLevelClosureResp(nullptr)
108 {
110 }
111 
115  if(fSampleSplitter) delete fSampleSplitter;
116 }
117 
120 
121  fSampleSplitter = new TRandom;
122 
125  auto zgbinning = GetZgBinning();
126  TArrayD binEdgesZg, binEdgesPtPart, binEdgesPtDet;
127  zgbinning->CreateBinEdges(binEdgesZg);
128  fPartLevelPtBinning->CreateBinEdges(binEdgesPtPart);
129  fDetLevelPtBinning->CreateBinEdges(binEdgesPtDet);
130 
131  fZgDetLevel = new TH2D("hZgDetLevel", "Zg response at detector level", binEdgesZg.GetSize() - 1, binEdgesZg.GetArray(), binEdgesPtDet.GetSize() - 1, binEdgesPtDet.GetArray());
132  fZgPartLevel = new TH2D("hZgPartLevel", "Zg response at particle level", binEdgesZg.GetSize() - 1, binEdgesZg.GetArray(), binEdgesPtDet.GetSize() - 1, binEdgesPtDet.GetArray());
133  fZgPartLevelTruncated = new TH2D("hZgPartLevelTruncated", "Zg response at particle level (truncated)", binEdgesZg.GetSize() - 1, binEdgesZg.GetArray(), binEdgesPtDet.GetSize() - 1, binEdgesPtDet.GetArray());
134 
135  // For closure test
136  fZgPartLevelClosureNoResp = new TH2D("hZgPartLevelClosureNoResp", "Zg response at particle level (closure test, jets not used for the response matrix)", binEdgesZg.GetSize() - 1, binEdgesZg.GetArray(), binEdgesPtDet.GetSize() - 1, binEdgesPtDet.GetArray());
137  fZgDetLevelClosureNoResp = new TH2D("hZgDetLevelClosureNoResp", "Zg response at detector level (closure test, jets not used for the response matrix)", binEdgesZg.GetSize() - 1, binEdgesZg.GetArray(), binEdgesPtDet.GetSize() - 1, binEdgesPtDet.GetArray());
138  fZgPartLevelClosureResp = new TH2D("hZgPartLevelClosureResp", "Zg response at particle level (closure test, jets used for the response matrix)", binEdgesZg.GetSize() - 1, binEdgesZg.GetArray(), binEdgesPtDet.GetSize() - 1, binEdgesPtDet.GetArray());
139  fZgDetLevelClosureResp = new TH2D("hZgDetLevelClosureResp", "Zg response at detector level (closure test, jets used for the response matrix)", binEdgesZg.GetSize() - 1, binEdgesZg.GetArray(), binEdgesPtDet.GetSize() - 1, binEdgesPtDet.GetArray());
140 
141  fZgResponse = new RooUnfoldResponse("hZgResponse", "z_{g} response matrix");
142  fZgResponse->Setup(fZgDetLevel, fZgPartLevel);
143  fZgResponseClosure = new RooUnfoldResponse("hZgResponseClosure", "z_{g} response matrix for the closure test");
144  fZgResponseClosure->Setup(fZgDetLevel, fZgPartLevel);
145 
146  fOutput->Add(fZgResponse);
148  fOutput->Add(fZgDetLevel);
149  fOutput->Add(fZgPartLevel);
150  fOutput->Add(fZgPartLevelTruncated);
151  fOutput->Add(fZgPartLevelClosureNoResp);
152  fOutput->Add(fZgDetLevelClosureNoResp);
153  fOutput->Add(fZgPartLevelClosureResp);
154  fOutput->Add(fZgDetLevelClosureResp);
155 
156  PostData(1, fOutput);
157 }
158 
160  AliJetContainer *partLevelJets = this->GetJetContainer("partLevel"),
161  *detLevelJets = GetJetContainer("detLevel");
164  AliParticleContainer *particles = GetParticleContainer("mcparticles");
165  if(!(partLevelJets || detLevelJets)) {
166  AliErrorStream() << "Either of the jet containers not found" << std::endl;
167  return kFALSE;
168  }
169 
170  // get truncations at detector level
171  auto ptmindet = fZgDetLevel->GetYaxis()->GetBinLowEdge(1),
172  ptmaxdet = fZgDetLevel->GetYaxis()->GetBinUpEdge(fZgDetLevel->GetYaxis()->GetNbins());
173 
174  for(auto detjet : detLevelJets->accepted()){
175  auto partjet = detjet->ClosestJet();
176  if(!partjet) continue;
177 
178  // sample splitting (for closure test)
179  bool closureUseResponse = (fSampleSplitter->Uniform() < fFractionResponseClosure);
180 
181  // Get the softdrop response
182  std::vector<double> softdropDet, softdropPart;
183  try {
184  softdropDet = MakeSoftdrop(*detjet, detLevelJets->GetJetRadius(), tracks, clusters);
185  softdropPart = MakeSoftdrop(*partjet, partLevelJets->GetJetRadius(), particles, nullptr);
186  fZgPartLevel->Fill(softdropPart[0], partjet->Pt());
187  if(detjet->Pt() >= ptmindet && detjet->Pt() <= ptmaxdet) {
188  fZgPartLevelTruncated->Fill(softdropPart[0], partjet->Pt());
189  fZgDetLevel->Fill(softdropDet[0], detjet->Pt());
190  fZgResponse->Fill(softdropDet[0], detjet->Pt(), softdropPart[0], partjet->Pt());
191  if(closureUseResponse) {
192  fZgResponseClosure->Fill(softdropDet[0], detjet->Pt(), softdropPart[0], partjet->Pt());
193  fZgDetLevelClosureResp->Fill(softdropDet[0], detjet->Pt());
194  fZgPartLevelClosureResp->Fill(softdropPart[0], partjet->Pt());
195  } else {
196  fZgDetLevelClosureNoResp->Fill(softdropDet[0], detjet->Pt());
197  fZgPartLevelClosureNoResp->Fill(softdropPart[0], partjet->Pt());
198  }
199  }
200  } catch(...){
201  AliErrorStream() << "Error in softdrop evaluation - jet will be ignored" << std::endl;
202  continue;
203  }
204 
205  }
206  return kTRUE;
207 }
208 
209 std::vector<double> AliAnalysisTaskEmcalSoftDropResponse::MakeSoftdrop(const AliEmcalJet &jet, double jetradius, const AliParticleContainer *tracks, const AliClusterContainer *clusters) const {
210  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
211  std::vector<fastjet::PseudoJet> constituents;
212  bool isMC = dynamic_cast<const AliMCParticleContainer *>(tracks);
213  AliDebugStream(2) << "Make new jet substrucutre for " << (isMC ? "MC" : "data") << " jet: Number of tracks " << jet.GetNumberOfTracks() << ", clusters " << jet.GetNumberOfClusters() << std::endl;
214  if(tracks && (fUseChargedConstituents || isMC)){ // Neutral particles part of particle container in case of MC
215  AliDebugStream(1) << "Jet substructure: Using charged constituents" << std::endl;
216  for(int itrk = 0; itrk < jet.GetNumberOfTracks(); itrk++){
217  auto track = jet.TrackAt(itrk, tracks->GetArray());
218  if(!track->Charge() && !fUseNeutralConstituents) continue; // Reject neutral constituents in case of using only charged consituents
219  if(track->Charge() && !fUseChargedConstituents) continue; // Reject charged constituents in case of using only neutral consituents
220  fastjet::PseudoJet constituentTrack(track->Px(), track->Py(), track->Pz(), track->E());
221  constituentTrack.set_user_index(jet.TrackAt(itrk));
222  constituents.push_back(constituentTrack);
223  }
224  }
225 
226  if(clusters && fUseNeutralConstituents){
227  AliDebugStream(1) << "Jet substructure: Using neutral constituents" << std::endl;
228  for(int icl = 0; icl < jet.GetNumberOfClusters(); icl++) {
229  auto cluster = jet.ClusterAt(icl, clusters->GetArray());
230  TLorentzVector clustervec;
231  cluster->GetMomentum(clustervec, fVertex, (AliVCluster::VCluUserDefEnergy_t)clusters->GetDefaultClusterEnergy());
232  fastjet::PseudoJet constituentCluster(clustervec.Px(), clustervec.Py(), clustervec.Pz(), cluster->GetHadCorrEnergy());
233  constituentCluster.set_user_index(jet.ClusterAt(icl) + kClusterOffset);
234  constituents.push_back(constituentCluster);
235  }
236  }
237 
238  AliDebugStream(3) << "Found " << constituents.size() << " constituents for jet with pt=" << jet.Pt() << " GeV/c" << std::endl;
239  if(!constituents.size()){
240  AliErrorStream() << "Jet has 0 constituents." << std::endl;
241  throw 1;
242  }
243  // Redo jet finding on constituents with a
244  fastjet::JetDefinition jetdef(fastjet::antikt_algorithm, jetradius*2, static_cast<fastjet::RecombinationScheme>(0), fastjet::BestFJ30 );
245  fastjet::ClusterSequence jetfinder(constituents, jetdef);
246  std::vector<fastjet::PseudoJet> outputjets = jetfinder.inclusive_jets(0);
247  auto sdjet = outputjets[0];
248  fastjet::contrib::SoftDrop softdropAlgorithm(fBeta, fZcut);
249  softdropAlgorithm.set_verbose_structure(kTRUE);
250  fastjet::JetAlgorithm reclusterizingAlgorithm;
251  switch(fReclusterizer) {
252  case kCAAlgo: reclusterizingAlgorithm = fastjet::cambridge_aachen_algorithm; break;
253  case kKTAlgo: reclusterizingAlgorithm = fastjet::kt_algorithm; break;
254  case kAKTAlgo: reclusterizingAlgorithm = fastjet::antikt_algorithm; break;
255  };
256  std::unique_ptr<fastjet::contrib::Recluster> reclusterizer(new fastjet::contrib::Recluster(reclusterizingAlgorithm, 1, true));
257  softdropAlgorithm.set_reclustering(kTRUE, reclusterizer.get());
258  AliDebugStream(4) << "Jet has " << sdjet.constituents().size() << " constituents" << std::endl;
259  auto groomed = softdropAlgorithm(sdjet);
260  auto softdropstruct = groomed.structure_of<fastjet::contrib::SoftDrop>();
261 
262  std::vector<double> result ={softdropstruct.symmetry(),
263  groomed.m(),
264  softdropstruct.delta_R(),
265  groomed.perp(),
266  softdropstruct.mu(),
267  static_cast<double>(softdropstruct.dropped_count())};
268  return result;
269 }
270 
272  auto binning = new TCustomBinning;
273  binning->SetMinimum(0);
274  switch(fBinningMode) {
275  case kSDModeINT7: {
276  binning->AddStep(20., 20.);
277  binning->AddStep(40., 10.);
278  binning->AddStep(80., 20.);
279  binning->AddStep(120., 40.);
280  binning->AddStep(240., 120.);
281  break;
282  }
283  case kSDModeEJ1: {
284  binning->AddStep(80., 80.);
285  binning->AddStep(140., 10.);
286  binning->AddStep(200., 20.);
287  binning->AddStep(240., 40.);
288  binning->AddStep(400., 160.);
289  break;
290  }
291  case kSDModeEJ2: {
292  binning->AddStep(70., 70.);
293  binning->AddStep(100., 10.);
294  binning->AddStep(140., 20.);
295  binning->AddStep(400., 260.);
296  break;
297  }
298  };
299  return binning;
300 }
301 
303  auto binning = new TCustomBinning;
304  switch(fBinningMode) {
305  case kSDModeINT7: {
306  binning->SetMinimum(20);
307  binning->AddStep(40., 5.);
308  binning->AddStep(60., 10.);
309  binning->AddStep(80., 20.);
310  binning->AddStep(120., 40.);
311  break;
312  }
313  case kSDModeEJ1: {
314  binning->SetMinimum(80.);
315  binning->AddStep(120., 5.);
316  binning->AddStep(160., 10.);
317  binning->AddStep(200., 20.);
318  break;
319  }
320  case kSDModeEJ2: {
321  binning->SetMinimum(70.);
322  binning->AddStep(100., 5.);
323  binning->AddStep(120., 10.);
324  binning->AddStep(140., 20.);
325  break;
326  }
327  };
328  return binning;
329 }
330 
332  auto binning = new TCustomBinning;
333  binning->SetMinimum(0.);
334  binning->AddStep(0.1, 0.1);
335  binning->AddStep(0.5, 0.05);
336  return binning;
337 }
338 
340  AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
341 
342  Bool_t isAOD(kFALSE);
343  AliInputEventHandler *inputhandler = static_cast<AliInputEventHandler *>(mgr->GetInputEventHandler());
344  if(inputhandler) {
345  if(inputhandler->IsA() == AliAODInputHandler::Class()){
346  std::cout << "Analysing AOD events\n";
347  isAOD = kTRUE;
348  } else {
349  std::cout << "Analysing ESD events\n";
350  }
351  }
352 
353  std::stringstream taskname;
354  taskname << "SoftdropResponsemaker_R" << std::setw(2) << std::setfill('0') << int(jetradius*10) << trigger;
355  AliAnalysisTaskEmcalSoftDropResponse *responsemaker = new AliAnalysisTaskEmcalSoftDropResponse(taskname.str().data());
356  responsemaker->SelectCollisionCandidates(AliVEvent::kINT7);
357  mgr->AddTask(responsemaker);
358 
359  AliParticleContainer *particles = responsemaker->AddMCParticleContainer("mcparticles");
360  particles->SetMinPt(0.);
361 
362  AliJetContainer *mcjets = responsemaker->AddJetContainer(
363  jettype,
365  recombinationScheme,
366  jetradius,
368  particles, nullptr);
369  mcjets->SetName("partLevel");
370  mcjets->SetJetPtCut(0.);
371  mcjets->SetMaxTrackPt(1000.);
372 
373  AliTrackContainer *tracks(nullptr);
374  if((jettype == AliJetContainer::kChargedJet) || (jettype == AliJetContainer::kFullJet)){
376  std::cout << "Track container name: " << tracks->GetName() << std::endl;
377  tracks->SetMinPt(0.15);
378  }
379  AliClusterContainer *clusters(nullptr);
380  if((jettype == AliJetContainer::kFullJet) || (jettype == AliJetContainer::kNeutralJet)){
381  std::cout << "Using full or neutral jets ..." << std::endl;
383  std::cout << "Cluster container name: " << clusters->GetName() << std::endl;
384  clusters->SetClusHadCorrEnergyCut(0.3); // 300 MeV E-cut
385  clusters->SetDefaultClusterEnergy(AliVCluster::kHadCorr);
386  } else {
387  std::cout << "Using charged jets ... " << std::endl;
388  }
389 
390  AliJetContainer *datajets = responsemaker->AddJetContainer(
391  jettype,
393  recombinationScheme,
394  jetradius,
396  tracks, clusters);
397  datajets->SetName("detLevel");
398  datajets->SetJetPtCut(0.);
399  datajets->SetMaxTrackPt(1000.);
400 
401  std::string jettypestring;
402  switch(jettype) {
403  case AliJetContainer::kFullJet: jettypestring = "FullJets"; break;
404  case AliJetContainer::kChargedJet: jettypestring = "ChargedJets"; break;
405  case AliJetContainer::kNeutralJet: jettypestring = "NeutralJets"; break;
406  default: jettypestring = "Undef";
407  };
408 
409  EBinningMode_t binmode(kSDModeINT7);
410  std::string triggerstring(trigger);
411  if(triggerstring == "EJ1") binmode = kSDModeEJ1;
412  else if(triggerstring == "EJ2") binmode = kSDModeEJ2;
413  responsemaker->SetBinningMode(binmode);
414 
415  // Connecting containers
416  std::stringstream outputfile, histname;
417  outputfile << mgr->GetCommonFileName() << ":SoftDropResponse_" << jettypestring << "_R" << std::setw(2) << std::setfill('0') << int(jetradius * 10.) << "_" << trigger;
418  histname << "SoftDropResponseHistos_" << jettypestring << "_R" << std::setw(2) << std::setfill('0') << int(jetradius * 10.) << "_" << trigger;
419  mgr->ConnectInput(responsemaker, 0, mgr->GetCommonInputContainer());
420  mgr->ConnectOutput(responsemaker, 1, mgr->CreateContainer(histname.str().data(), AliEmcalList::Class(), AliAnalysisManager::kOutputContainer, outputfile.str().data()));
421 
422  return responsemaker;
423 }
double Double_t
Definition: External.C:58
TH2 * fZgPartLevelClosureNoResp
! Zg vs. pt at particle level for closure test (jets not used for response)
AliJetContainer * GetJetContainer(Int_t i=0) const
TH2 * fZgPartLevelClosureResp
! Zg vs. pt at particle level for closure test (jets used for response)
Container with name, TClonesArray and cuts for particles.
Int_t ClusterAt(Int_t idx) const
Definition: AliEmcalJet.h:137
TH2 * fZgDetLevelClosureNoResp
! Zg vs. pt at detector level for closure test (jets not used for response)
AliJetContainer * AddJetContainer(const char *n, TString defaultCutType, Float_t jetRadius=0.4)
virtual void CreateBinEdges(TArrayD &binedges) const =0
TH2 * fZgPartLevelTruncated
! Zg vs. pt at particle level after truncation at
TH2 * fZgDetLevelClosureResp
! Zg vs. pt at detector level for closure test (jets used for response)
Interface for binnings used by the histogram handler.
Definition: TBinning.h:23
AliClusterContainer * AddClusterContainer(const char *n)
Create new cluster container and attach it to the task.
Container for particles within the EMCAL framework.
Int_t GetDefaultClusterEnergy() const
Int_t TrackAt(Int_t idx) const
Definition: AliEmcalJet.h:160
UShort_t GetNumberOfTracks() const
Definition: AliEmcalJet.h:139
Double_t fFractionResponseClosure
Fraction of events used for the response matrix in the closure test.
AliParticleContainer * GetParticleContainer(Int_t i=0) const
Get particle container attached to this task.
TPC fiducial acceptance (each eta edge narrowed by jet R)
Definition: AliEmcalJet.h:68
UShort_t GetNumberOfClusters() const
Definition: AliEmcalJet.h:138
void SetJetPtCut(Float_t cut)
TPC acceptance.
Definition: AliEmcalJet.h:67
Definition: External.C:228
AliClusterContainer * GetClusterContainer(Int_t i=0) const
Get cluster container attached to this task.
Helper class creating user defined custom binning.
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.
RooUnfoldResponse * fZgResponseClosure
! RooUnfold response for the closure test
AliMCParticleContainer * AddMCParticleContainer(const char *n)
Create new container for MC particles and attach it to the task.
Double_t Pt() const
Definition: AliEmcalJet.h:109
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.
Represent a jet reconstructed using the EMCal jet framework.
Definition: AliEmcalJet.h:51
void UserCreateOutputObjects()
Main initialization function on the worker.
bool Bool_t
Definition: External.C:53
static AliAnalysisTaskEmcalSoftDropResponse * AddTaskEmcalSoftDropResponse(Double_t jetradius, AliJetContainer::EJetType_t jettype, AliJetContainer::ERecoScheme_t recombinationScheme, const char *trigger)
void SetDefaultClusterEnergy(Int_t d)
void SetMaxTrackPt(Float_t b)
Container structure for EMCAL clusters.
Container for MC-true particles within the EMCAL framework.
EMCal fiducial acceptance (each eta, phi edge narrowed by jet R)
Definition: AliEmcalJet.h:70
Container for jet within the EMCAL jet framework.
std::vector< double > MakeSoftdrop(const AliEmcalJet &jet, double jetradius, const AliParticleContainer *tracks, const AliClusterContainer *clusters) const
void SetClusHadCorrEnergyCut(Double_t cut)
void SetMinimum(Double_t min)
static TString TrackContainerNameFactory(Bool_t isAOD)
Get name of the default track container.
static TString ClusterContainerNameFactory(Bool_t isAOD)
Get name of the default cluster container.