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