AliPhysics  3aa38c9 (3aa38c9)
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 <array>
28 #include <iostream>
29 #include <string>
30 #include <set>
31 #include <sstream>
32 #include <vector>
33 
34 #include <fastjet/ClusterSequence.hh>
35 #include <fastjet/contrib/Nsubjettiness.hh>
36 #include <fastjet/contrib/SoftDrop.hh>
37 
38 #include <THistManager.h>
39 #include <TLinearBinning.h>
40 #include <TLorentzVector.h>
41 #include <TMath.h>
42 #include <TObjString.h>
43 #include <TString.h>
44 #include <TVector3.h>
45 
46 #include "AliAODEvent.h"
47 #include "AliAODInputHandler.h"
48 #include "AliAnalysisManager.h"
49 #include "AliAnalysisDataSlot.h"
50 #include "AliAnalysisDataContainer.h"
52 #include "AliCDBEntry.h"
53 #include "AliCDBManager.h"
54 #include "AliClusterContainer.h"
55 #include "AliJetContainer.h"
58 #include "AliEmcalJet.h"
59 #include "AliEmcalList.h"
61 #include "AliLog.h"
62 #include "AliParticleContainer.h"
63 #include "AliRhoParameter.h"
64 #include "AliTrackContainer.h"
65 #include "AliTriggerCluster.h"
66 #include "AliTriggerConfiguration.h"
67 #include "AliVCluster.h"
68 #include "AliVParticle.h"
69 
70 #ifdef EXPERIMENTAL_JETCONSTITUENTS
73 #endif
74 
75 
79 
80 namespace EmcalTriggerJets {
81 
84  fJetSubstructureTree(nullptr),
85  fQAHistos(nullptr),
86  fLumiMonitor(nullptr),
87  fSDZCut(0.1),
88  fSDBetaCut(0),
89  fReclusterizer(kCAAlgo),
90  fTriggerSelectionBits(AliVEvent::kAny),
91  fTriggerSelectionString(""),
92  fNameTriggerDecisionContainer("EmcalTriggerDecision"),
93  fUseTriggerSelectionForData(false),
94  fUseDownscaleWeight(false),
95  fUseChargedConstituents(true),
96  fUseNeutralConstituents(true),
97  fFillPart(true),
98  fFillAcceptance(true),
99  fFillRho(true),
100  fFillMass(true),
101  fFillSoftDrop(true),
102  fFillNSub(true),
103  fFillStructGlob(true)
104 {
105  memset(fJetTreeData, 0, sizeof(Double_t) * kTNVar);
106 }
107 
109  AliAnalysisTaskEmcalJet(name, kTRUE),
113  fSDZCut(0.1),
114  fSDBetaCut(0),
116  fTriggerSelectionBits(AliVEvent::kAny),
118  fNameTriggerDecisionContainer("EmcalTriggerDecision"),
120  fUseDownscaleWeight(false),
123  fFillPart(true),
124  fFillAcceptance(true),
125  fFillRho(true),
126  fFillMass(true),
127  fFillSoftDrop(true),
128  fFillNSub(true),
129  fFillStructGlob(true)
130 {
131  memset(fJetTreeData, 0, sizeof(Double_t) * kTNVar);
132  DefineOutput(2, TTree::Class());
133 }
134 
136 
137 }
138 
141 
142  // Make QA for constituent clusters
143  TLinearBinning jetptbinning(9, 20, 200),
144  clusterenergybinning(200, 0., 200),
145  timebinning(1000, -500., 500.),
146  m02binning(100, 0., 1.),
147  ncellbinning(101, -0.5, 100.5);
148  fQAHistos = new THistManager("QAhistos");
149  fQAHistos->CreateTH1("hEventCounter", "Event counter", 1, 0.5, 1.5);
150  fQAHistos->CreateTH2("hClusterConstE", "EMCAL cluster energy vs jet pt; p_{t, jet} (GeV/c); E_{cl} (GeV)", jetptbinning, clusterenergybinning);
151  fQAHistos->CreateTH2("hClusterConstTime", "EMCAL cluster time vs. jet pt; p_{t, jet} (GeV/c); t_{cl} (ns)", jetptbinning, timebinning);
152  fQAHistos->CreateTH2("hClusterConstM02", "EMCAL cluster M02 vs. jet pt; p{t, jet} (GeV/c); M02", jetptbinning, m02binning);
153  fQAHistos->CreateTH2("hClusterConstNcell", "EMCAL cluster ncell vs. jet pt; p{t, jet} (GeV/c); Number of cells", jetptbinning, ncellbinning);
154 
155  // Test of constituent QA
156 #ifdef EXPERIMENTAL_JETCONSTITUENTS
157  fQAHistos->CreateTH2("hChargedConstituentPt", "charged constituent pt vs jet pt (via constituent map); p_{t,jet} (GeV/c); p_{t,ch} (GeV/c)", jetptbinning, clusterenergybinning);
158  fQAHistos->CreateTH2("hChargedIndexPt", "charged constituent pt vs jet pt (via index map); p_{t, jet} (GeV/c); p_{t, ch} (GeV/c)", jetptbinning, clusterenergybinning);
159 
160  fQAHistos->CreateTH2("hClusterConstituentEDefault", "cluster constituent default energy vs. jet pt (va constituent map); p_{t, jet} (GeV/c); E_{cl} (GeV)", jetptbinning, clusterenergybinning);
161  fQAHistos->CreateTH2("hClusterConstituentENLC", "cluster constituent non-linearity-corrected energy vs. jet pt (va constituent map); p_{t, jet} (GeV/c); E_{cl} (GeV)", jetptbinning, clusterenergybinning);
162  fQAHistos->CreateTH2("hClusterConstituentEHC", "cluster constituent hadronic-corrected energy vs. jet pt (va constituent map); p_{t, jet} (GeV/c); E_{cl} (GeV)", jetptbinning, clusterenergybinning);
163  fQAHistos->CreateTH2("hClusterIndexENLC", "cluster constituent non-linearity-corrected energy vs. jet pt (via index map); p_{t, jet} (GeV/c); E_{cl} (GeV)", jetptbinning, clusterenergybinning);
164  fQAHistos->CreateTH2("hClusterIndexEHC", "cluster constituent hadronic-corrected energy vs. jet pt (via index map); p_{t, jet} (GeV/c); E_{cl} (GeV)", jetptbinning, clusterenergybinning);
165 #endif
166  for(auto h : *(fQAHistos->GetListOfHistograms())) fOutput->Add(h);
167 
168  OpenFile(2);
169  std::string treename = this->GetOutputSlot(2)->GetContainer()->GetName();
170  fJetSubstructureTree = new TTree(treename.data(), "Tree with jet substructure information");
171  std::vector<std::string> varnames = {
172  "Radius", "EventWeight", "PtJetRec", "PtJetSim", "EJetRec",
173  "EJetSim", "EtaRec", "EtaSim", "PhiRec", "PhiSim",
174  "RhoPtRec", "RhoPtSim", "RhoMassRec", "RhoMassSim", "AreaRec",
175  "AreaSim", "NEFRec", "NEFSim", "MassRec", "MassSim",
176  "ZgMeasured", "ZgTrue", "RgMeasured", "RgTrue", "MgMeasured",
177  "MgTrue", "PtgMeasured", "PtgTrue", "MugMeasured", "MugTrue",
178  "OneSubjettinessMeasured", "OneSubjettinessTrue", "TwoSubjettinessMeasured", "TwoSubjettinessTrue", "AngularityMeasured",
179  "AngularityTrue", "PtDMeasured", "PtDTrue", "NCharged", "NNeutral",
180  "NConstTrue", "NDroppedMeasured", "NDroppedTrue" };
181 
182  for(int ib = 0; ib < kTNVar; ib++){
183  LinkOutputBranch(varnames[ib], fJetTreeData + ib);
184  }
185  PostData(1, fOutput);
186  PostData(2, fJetSubstructureTree);
187 }
188 
189 void AliAnalysisTaskEmcalJetSubstructureTree::LinkOutputBranch(const std::string &branchname, Double_t *datalocation) {
190  // Check whether branch is rejected
191  if(!fFillPart && IsPartBranch(branchname)) return;
192  if(!fFillAcceptance && IsAcceptanceBranch(branchname)) return;
193  if(!fFillRho && IsRhoBranch(branchname)) return;
194  if(!fFillMass && IsMassBranch(branchname)) return;
195  if(!fFillSoftDrop && IsSoftdropBranch(branchname)) return;
196  if(!fFillNSub && IsNSubjettinessBranch(branchname)) return;
197  if(!fFillStructGlob && IsStructbranch(branchname)) return;
198 
199  std::cout << "Adding branch " << branchname << std::endl;
200  fJetSubstructureTree->Branch(branchname.data(), datalocation, Form("%s/D", branchname.data()));
201 }
202 
206  }
207 }
208 
212  AliParticleContainer *particles = GetParticleContainer("mcparticles");
213 
214  AliJetContainer *mcjets = GetJetContainer("mcjets");
215  AliJetContainer *datajets = GetJetContainer("datajets");
216 
217  FillLuminosity(); // Makes only sense in data
218 
219  // for(auto e : *(fInputEvent->GetList())) std::cout << e->GetName() << std::endl;
220 
221  std::stringstream rhoTagData, rhoTagMC;
222  if(datajets) rhoTagData << "R" << std::setw(2) << std::setfill('0') << static_cast<Int_t>(datajets->GetJetRadius() * 10.);
223  if(mcjets) rhoTagMC << "R" << std::setw(2) << std::setfill('0') << static_cast<Int_t>(mcjets->GetJetRadius() * 10.);
224 
225  std::string rhoSparseData = "RhoSparse_Full_" + rhoTagData.str(), rhoSparseMC = "RhoSparse_Full_" + rhoTagMC.str(),
226  rhoMassData = "RhoMassSparse_Full_" + rhoTagData.str(), rhoMassMC = "RhoMassSparse_Full_" + rhoTagMC.str();
227  AliRhoParameter *rhoPtRec = GetRhoFromEvent(rhoSparseData.data()),
228  *rhoMassRec = GetRhoFromEvent(rhoMassData.data()),
229  *rhoPtSim = GetRhoFromEvent(rhoSparseMC.data()),
230  *rhoMassSim = GetRhoFromEvent(rhoMassMC.data());
231  AliDebugStream(2) << "Found rho parameter for reconstructed pt: " << (rhoPtRec ? "yes" : "no") << ", value: " << (rhoPtRec ? rhoPtRec->GetVal() : 0.) << std::endl;
232  AliDebugStream(2) << "Found rho parameter for sim pt: " << (rhoPtSim ? "yes" : "no") << ", value: " << (rhoPtSim ? rhoPtSim->GetVal() : 0.) << std::endl;
233  AliDebugStream(2) << "Found rho parameter for reconstructed Mass: " << (rhoMassRec ? "yes" : "no") << ", value: " << (rhoMassRec ? rhoMassRec->GetVal() : 0.) << std::endl;
234  AliDebugStream(2) << "Found rho parameter for sim Mass: " << (rhoMassSim ? "yes" : "no") << ", value: " << (rhoMassSim ? rhoMassSim->GetVal() : 0.) << std::endl;
235 
236  AliDebugStream(1) << "Inspecting jet radius " << (datajets ? datajets->GetJetRadius() : mcjets->GetJetRadius()) << std::endl;
237 
238  // Run trigger selection (not on pure MCgen train)
239  if(datajets){
240  if(!(fInputHandler->IsEventSelected() & fTriggerSelectionBits)) return false;
241  if(!mcjets){
242  // Pure data - do EMCAL trigger selection from selection string
243  if(fTriggerSelectionString.Length()) {
244  if(!fInputEvent->GetFiredTriggerClasses().Contains(fTriggerSelectionString)) return false;
246  auto trgselresult = static_cast<PWG::EMCAL::AliEmcalTriggerDecisionContainer *>(fInputEvent->FindListObject(fNameTriggerDecisionContainer));
247  AliDebugStream(1) << "Found trigger decision object: " << (trgselresult ? "yes" : "no") << std::endl;
248  if(!trgselresult){
249  AliErrorStream() << "Trigger decision container with name " << fNameTriggerDecisionContainer << " not found in event - not possible to select EMCAL triggers" << std::endl;
250  return false;
251  }
252  if(!trgselresult->IsEventSelected(fTriggerSelectionString)) return false;
253  }
254  }
255  } else {
257  // Simulation - do EMCAL trigger selection from trigger selection object
258  auto mctrigger = static_cast<PWG::EMCAL::AliEmcalTriggerDecisionContainer *>(fInputEvent->FindListObject(fNameTriggerDecisionContainer));
259  AliDebugStream(1) << "Found trigger decision object: " << (mctrigger ? "yes" : "no") << std::endl;
260  if(!mctrigger){
261  AliErrorStream() << "Trigger decision container with name " << fNameTriggerDecisionContainer << " not found in event - not possible to select EMCAL triggers" << std::endl;
262  return false;
263  }
264  if(!mctrigger->IsEventSelected(fTriggerSelectionString)) return false;
265  }
266  }
267  }
268 
269  double weight = 1.;
271  AliDebugStream(2) << "Trigger selection string: " << fTriggerSelectionString << std::endl;
272  TString selectionString = (fTriggerSelectionBits & AliVEvent::kINT7) ? "INT7" : fTriggerSelectionString;
273  auto triggerstring = MatchTrigger(selectionString.Data());
274  AliDebugStream(2) << "Getting downscale correction factor for trigger string " << triggerstring << std::endl;
276  }
277  AliDebugStream(1) << "Using downscale weight " << weight << std::endl;
278 
279 
280  // Count events (for spectrum analysis)
281  fQAHistos->FillTH1("hEventCounter", 1);
282 
283  Double_t rhoparameters[4]; memset(rhoparameters, 0, sizeof(Double_t) * 4);
284  if(rhoPtRec) rhoparameters[0] = rhoPtRec->GetVal();
285  if(rhoPtSim) rhoparameters[1] = rhoPtSim->GetVal();
286  if(rhoMassRec) rhoparameters[2] = rhoMassRec->GetVal();
287  if(rhoMassSim) rhoparameters[3] = rhoMassSim->GetVal();
288 
289  AliSoftdropDefinition softdropSettings;
290  softdropSettings.fBeta = fSDBetaCut;
291  softdropSettings.fZ = fSDZCut;
292  switch(fReclusterizer) {
293  case kCAAlgo: softdropSettings.fRecluserAlgo = fastjet::cambridge_aachen_algorithm; break;
294  case kKTAlgo: softdropSettings.fRecluserAlgo = fastjet::kt_algorithm; break;
295  case kAKTAlgo: softdropSettings.fRecluserAlgo = fastjet::antikt_algorithm; break;
296  };
297 
298  AliNSubjettinessDefinition nsubjettinessSettings;
299  nsubjettinessSettings.fBeta = 1.;
300  nsubjettinessSettings.fRadius = 0.4;
301 
302  if(datajets) {
303  AliDebugStream(1) << "In data jets branch: found " << datajets->GetNJets() << " jets, " << datajets->GetNAcceptedJets() << " were accepted\n";
304  AliDebugStream(1) << "Having MC information: " << (mcjets ? TString::Format("yes, with %d jets", mcjets->GetNJets()) : "no") << std::endl;
305  if(mcjets) {
306  AliDebugStream(1) << "In MC jets branch: found " << mcjets->GetNJets() << " jets, " << mcjets->GetNAcceptedJets() << " were accepted\n";
307  }
308  for(auto jet : datajets->accepted()) {
309  double pt = jet->Pt(), pz = jet->Pz(), E = jet->E(), M = TMath::Sqrt(E*E - pt*pt - pz*pz);
310  AliDebugStream(2) << "Next jet: pt:" << jet->Pt() << ", E: " << E << ", pz: " << pz << ", M(self): " << M << "M(fj)" << jet->M() << std::endl;
311  AliEmcalJet *associatedJet = jet->ClosestJet();
312 
313  if(mcjets) {
314  if(!associatedJet) {
315  AliDebugStream(2) << "Not found associated jet" << std::endl;
316  continue;
317  }
318  if(!(SelectJet(*jet, tracks) && SelectJet(*associatedJet, particles))) continue;
319  try {
320  DoConstituentQA(jet, tracks, clusters);
321  AliJetSubstructureData structureData = MakeJetSubstructure(*jet, datajets->GetJetRadius() * 2., tracks, clusters, {softdropSettings, nsubjettinessSettings}),
322  structureMC = fFillPart ? MakeJetSubstructure(*associatedJet, mcjets->GetJetRadius() * 2, particles, nullptr, {softdropSettings, nsubjettinessSettings}) : AliJetSubstructureData();
323  Double_t angularity[2] = {fFillStructGlob ? MakeAngularity(*jet, tracks, clusters) : 0., (fFillStructGlob && fFillPart) ? MakeAngularity(*associatedJet, particles, nullptr) : 0.},
324  ptd[2] = {fFillStructGlob ? MakePtD(*jet, tracks, clusters) : 0., (fFillStructGlob && fFillPart) ? MakePtD(*associatedJet, particles, nullptr) : 0};
325  FillTree(datajets->GetJetRadius(), weight, jet, associatedJet, &(structureData.fSoftDrop), &(structureMC.fSoftDrop), &(structureData.fNsubjettiness), &(structureMC.fNsubjettiness), angularity, ptd, rhoparameters);
326  } catch(ReclusterizerException &e) {
327  AliErrorStream() << "Error in reclusterization - skipping jet" << std::endl;
328  } catch(SubstructureException &e) {
329  AliErrorStream() << "Error in substructure observable - skipping jet" << std::endl;
330  }
331  } else {
332  if(!SelectJet(*jet, tracks)) continue;
333  try {
334  DoConstituentQA(jet, tracks, clusters);
335  AliJetSubstructureData structure = MakeJetSubstructure(*jet, 0.4, tracks, clusters, {softdropSettings, nsubjettinessSettings});
336  Double_t angularity[2] = {fFillStructGlob ? MakeAngularity(*jet, tracks, clusters): 0., 0.},
337  ptd[2] = {fFillStructGlob ? MakePtD(*jet, tracks, clusters) : 0., 0.};
338  FillTree(datajets->GetJetRadius(), weight, jet, nullptr, &(structure.fSoftDrop), nullptr, &(structure.fNsubjettiness), nullptr, angularity, ptd, rhoparameters);
339  } catch(ReclusterizerException &e) {
340  AliErrorStream() << "Error in reclusterization - skipping jet" << std::endl;
341  } catch(SubstructureException &e) {
342  AliErrorStream() << "Error in substructure observable - skipping jet" << std::endl;
343  }
344  }
345  }
346  } else {
347  if(mcjets) {
348  // for MCgen train
349  AliDebugStream(1) << "In MC pure jet branch: found " << mcjets->GetNJets() << " jets, " << mcjets->GetNAcceptedJets() << " were accepted\n";
350  for(auto j : mcjets->accepted()){
351  AliEmcalJet *mcjet = static_cast<AliEmcalJet *>(j);
352  try {
353  AliJetSubstructureData structure = MakeJetSubstructure(*mcjet, mcjets->GetJetRadius() * 2., particles, nullptr,{softdropSettings, nsubjettinessSettings});
354  Double_t angularity[2] = {0., MakeAngularity(*mcjet, particles, nullptr)},
355  ptd[2] = {0., MakePtD(*mcjet, particles, nullptr)};
356  FillTree(mcjets->GetJetRadius(), weight, nullptr, mcjet, nullptr, &(structure.fSoftDrop), nullptr, &(structure.fNsubjettiness), angularity, ptd, rhoparameters);
357  } catch (ReclusterizerException &e) {
358  AliErrorStream() << "Error in reclusterization - skipping jet" << std::endl;
359  } catch (SubstructureException &e) {
360  AliErrorStream() << "Error in substructure observable - skipping jet" << std::endl;
361  }
362  }
363  }
364  }
365 
366  return true;
367 }
368 
370  AliCDBManager * cdb = AliCDBManager::Instance();
371  if(!fMCEvent && cdb){
372  // Get List of trigger clusters
373  AliCDBEntry *en = cdb->Get("GRP/CTP/Config");
374  AliTriggerConfiguration *trg = static_cast<AliTriggerConfiguration *>(en->GetObject());
375  std::vector<std::string> clusternames;
376  for(auto c : trg->GetClusters()) {
377  AliTriggerCluster *clust = static_cast<AliTriggerCluster *>(c);
378  std::string clustname = clust->GetName();
379  auto iscent = clustname.find("CENT") != std::string::npos, iscalo = clustname.find("CALO") != std::string::npos;
380  if(!(iscalo || iscent)) continue;
381  clusternames.emplace_back(clustname);
382  }
383 
384  // Set the x-axis of the luminosity monitor histogram
385  fLumiMonitor = new TH1F("hLumiMonitor", "Luminosity monitor", clusternames.size(), 0, clusternames.size());
386  int currentbin(1);
387  for(auto c : clusternames) {
388  fLumiMonitor->GetXaxis()->SetBinLabel(currentbin++, c.data());
389  }
390  fOutput->Add(fLumiMonitor);
391  }
392 }
393 
395  const AliEmcalJet *datajet, const AliEmcalJet *mcjet,
396  AliSoftDropParameters *dataSoftdrop, AliSoftDropParameters *mcSoftdrop,
397  AliNSubjettinessParameters *dataSubjettiness, AliNSubjettinessParameters *mcSubjettiness,
398  Double_t *angularity, Double_t *ptd, Double_t *rhoparameters){
399 
400  memset(fJetTreeData,0, sizeof(Double_t) * kTNVar);
401  fJetTreeData[kTRadius] = r;
402  fJetTreeData[kTWeight] = weight;
403  if(fFillRho){
404  fJetTreeData[kTRhoPtRec] = rhoparameters[0];
405  if(fFillMass) fJetTreeData[kTRhoMassRec] = rhoparameters[2];
406  if(fFillPart){
407  fJetTreeData[kTRhoPtSim] = rhoparameters[1];
408  if(fFillMass) fJetTreeData[kTRhoMassSim] = rhoparameters[3];
409  }
410  }
411  if(datajet) {
412  fJetTreeData[kTPtJetRec] = TMath::Abs(datajet->Pt());
414  fJetTreeData[kTNNeutral] = fUseNeutralConstituents ? datajet->GetNumberOfClusters() : 0.;
415  fJetTreeData[kTAreaRec] = datajet->Area();
416  fJetTreeData[kTNEFRec] = datajet->NEF();
417  if(fFillMass) fJetTreeData[kTMassRec] = datajet->M();
418  fJetTreeData[kTEJetRec] = datajet->E();
419  if(fFillAcceptance) {
420  fJetTreeData[kTEtaRec] = datajet->Eta();
421  fJetTreeData[kTPhiRec] = datajet->Phi();
422  }
423  }
424 
425  if(fFillPart && mcjet){
426  fJetTreeData[kTPtJetSim] = TMath::Abs(mcjet->Pt());
428  fJetTreeData[kTAreaSim] = mcjet->Area();
429  fJetTreeData[kTNEFSim] = mcjet->NEF();
430  if(fFillMass) fJetTreeData[kTMassSim] = mcjet->M();
431  fJetTreeData[kTEJetSim] = mcjet->E();
432  if(fFillAcceptance){
433  fJetTreeData[kTEtaSim] = mcjet->Eta();
434  fJetTreeData[kTPhiSim] = mcjet->Phi();
435  }
436  }
437 
438  if(fFillSoftDrop){
439  if(dataSoftdrop) {
440  fJetTreeData[kTZgMeasured] = dataSoftdrop->fZg;
441  fJetTreeData[kTRgMeasured] = dataSoftdrop->fRg;
442  fJetTreeData[kTMgMeasured] = dataSoftdrop->fMg;
443  fJetTreeData[kTPtgMeasured] = dataSoftdrop->fPtg;
444  fJetTreeData[kTMugMeasured] = dataSoftdrop->fMug;
445  fJetTreeData[kTNDroppedMeasured] = dataSoftdrop->fNDropped;
446  }
447 
448  if(fFillPart && mcSoftdrop) {
449  fJetTreeData[kTZgTrue] = mcSoftdrop->fZg;
450  fJetTreeData[kTRgTrue] = mcSoftdrop->fRg;
451  fJetTreeData[kTMgTrue] = mcSoftdrop->fMg;
452  fJetTreeData[kTPtgTrue] = mcSoftdrop->fPtg;
453  fJetTreeData[kTMugTrue] = mcSoftdrop->fMug;
454  fJetTreeData[kTNDroppedTrue] = mcSoftdrop->fNDropped;
455  }
456  }
457 
458  if(fFillNSub){
459  if(dataSubjettiness) {
462  }
463 
464  if(fFillPart && mcSubjettiness) {
467  }
468  }
469 
470  if(fFillStructGlob){
471  fJetTreeData[kTAngularityMeasured] = angularity[0];
472  fJetTreeData[kTPtDMeasured] = ptd[0];
473 
474  if(fFillPart){
475  fJetTreeData[kTAngularityTrue] = angularity[1];
476  fJetTreeData[kTPtDTrue] = ptd[1];
477  }
478  }
479 
480  fJetSubstructureTree->Fill();
481 }
482 
486  if(fInputEvent->GetFiredTriggerClasses().Contains("INT7")) {
487  for(auto trigger : DecodeTriggerString(fInputEvent->GetFiredTriggerClasses().Data())){
488  auto int7trigger = trigger.IsTriggerClass("INT7");
489  auto bunchcrossing = trigger.fBunchCrossing == "B";
490  auto nopf = trigger.fPastFutureProtection == "NOPF";
491  AliDebugStream(4) << "Full name: " << trigger.ExpandClassName() << ", INT7 trigger: " << (int7trigger ? "Yes" : "No") << ", bunch crossing: " << (bunchcrossing ? "Yes" : "No") << ", no past-future protection: " << (nopf ? "Yes" : "No") << ", Cluster: " << trigger.fTriggerCluster << std::endl;
492  if(int7trigger && bunchcrossing && nopf) {
493  double downscale = downscalefactors->GetDownscaleFactorForTriggerClass(trigger.ExpandClassName());
494  AliDebugStream(5) << "Using downscale " << downscale << std::endl;
495  fLumiMonitor->Fill(trigger.fTriggerCluster.data(), 1./downscale);
496  }
497  }
498  }
499  }
500 }
501 
503  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
504  std::vector<fastjet::PseudoJet> constituents;
505  bool isMC = dynamic_cast<const AliMCParticleContainer *>(tracks);
506  AliDebugStream(2) << "Make new jet substrucutre for " << (isMC ? "MC" : "data") << " jet: Number of tracks " << jet.GetNumberOfTracks() << ", clusters " << jet.GetNumberOfClusters() << std::endl;
507  if(tracks && (fUseChargedConstituents || isMC)){ // Neutral particles part of particle container in case of MC
508  for(int itrk = 0; itrk < jet.GetNumberOfTracks(); itrk++){
509  auto track = jet.TrackAt(itrk, tracks->GetArray());
510  if(!track->Charge() && !fUseNeutralConstituents) continue; // Reject neutral constituents in case of using only charged consituents
511  if(track->Charge() && !fUseChargedConstituents) continue; // Reject charged constituents in case of using only neutral consituents
512  fastjet::PseudoJet constituentTrack(track->Px(), track->Py(), track->Pz(), track->E());
513  constituentTrack.set_user_index(jet.TrackAt(itrk));
514  constituents.push_back(constituentTrack);
515  }
516  }
517 
518  if(clusters && fUseNeutralConstituents){
519  for(int icl = 0; icl < jet.GetNumberOfClusters(); icl++) {
520  auto cluster = jet.ClusterAt(icl, clusters->GetArray());
521  TLorentzVector clustervec;
522  cluster->GetMomentum(clustervec, fVertex, (AliVCluster::VCluUserDefEnergy_t)clusters->GetDefaultClusterEnergy());
523  fastjet::PseudoJet constituentCluster(clustervec.Px(), clustervec.Py(), clustervec.Pz(), cluster->GetHadCorrEnergy());
524  constituentCluster.set_user_index(jet.ClusterAt(icl) + kClusterOffset);
525  constituents.push_back(constituentCluster);
526  }
527  }
528 
529  AliDebugStream(3) << "Found " << constituents.size() << " constituents for jet with pt=" << jet.Pt() << " GeV/c" << std::endl;
530  if(!constituents.size()){
531  AliErrorStream() << "Jet has 0 constituents." << std::endl;
532  throw ReclusterizerException();
533  }
534  // Redo jet finding on constituents with a
535  fastjet::JetDefinition jetdef(fastjet::antikt_algorithm, jetradius*2, static_cast<fastjet::RecombinationScheme>(0), fastjet::BestFJ30 );
536  std::vector<fastjet::PseudoJet> outputjets;
537  try {
538  fastjet::ClusterSequence jetfinder(constituents, jetdef);
539  outputjets = jetfinder.inclusive_jets(0);
541  return result;
542  } catch (fastjet::Error &e) {
543  AliErrorStream() << " FJ Exception caught: " << e.message() << std::endl;
544  throw ReclusterizerException();
545  } catch (SoftDropException &e) {
546  AliErrorStream() << "Softdrop exception caught: " << e.what() << std::endl;
547  throw ReclusterizerException();
548  }
549 }
550 
552  fastjet::contrib::SoftDrop softdropAlgorithm(cutparameters.fBeta, cutparameters.fZ);
553  softdropAlgorithm.set_verbose_structure(kTRUE);
554  std::unique_ptr<fastjet::contrib::Recluster> reclusterizer(new fastjet::contrib::Recluster(cutparameters.fRecluserAlgo, 1, true));
555  softdropAlgorithm.set_reclustering(kTRUE, reclusterizer.get());
556  AliDebugStream(4) << "Jet has " << jet.constituents().size() << " constituents" << std::endl;
557  auto groomed = softdropAlgorithm(jet);
558  try {
559  auto softdropstruct = groomed.structure_of<fastjet::contrib::SoftDrop>();
560 
561  AliSoftDropParameters result({softdropstruct.symmetry(),
562  groomed.m(),
563  softdropstruct.delta_R(),
564  groomed.perp(),
565  softdropstruct.mu(),
566  softdropstruct.dropped_count()});
567  return result;
568  } catch(std::bad_cast &e) {
569  throw SoftDropException();
570  }
571 }
572 
575  fastjet::contrib::Nsubjettiness (1,fastjet::contrib::KT_Axes(),fastjet::contrib::NormalizedCutoffMeasure(cut.fBeta, cut.fRadius, 1e100)).result(jet),
576  fastjet::contrib::Nsubjettiness (2,fastjet::contrib::KT_Axes(),fastjet::contrib::NormalizedCutoffMeasure(cut.fBeta, cut.fRadius, 1e100)).result(jet)
577  });
578  return result;
579 }
580 
582  if(!(jet.GetNumberOfTracks() || jet.GetNumberOfClusters()))
583  throw SubstructureException();
584  TVector3 jetvec(jet.Px(), jet.Py(), jet.Pz());
585  Double_t den(0.), num(0.);
586  bool isMC = dynamic_cast<const AliMCParticleContainer *>(tracks);
587  if(tracks && (fUseChargedConstituents || isMC)){
588  for(UInt_t itrk = 0; itrk < jet.GetNumberOfTracks(); itrk++) {
589  auto track = jet.TrackAt(itrk, tracks->GetArray());
590  if(!track){
591  AliErrorStream() << "Associated constituent particle / track not found\n";
592  continue;
593  }
594  if(!track->Charge() && !fUseNeutralConstituents) continue; // Reject neutral constituents in case of using only charged consituents
595  if(track->Charge() && !fUseChargedConstituents) continue; // Reject charged constituents in case of using only neutral consituents
596  TVector3 trackvec(track->Px(), track->Py(), track->Pz());
597 
598  num += track->Pt() * trackvec.DrEtaPhi(jetvec);
599  den += +track->Pt();
600  }
601  }
602  if(clusters && fUseNeutralConstituents) {
603  for(UInt_t icl = 0; icl < jet.GetNumberOfClusters(); icl++){
604  auto clust = jet.ClusterAt(icl, clusters->GetArray());
605  if(!clust) {
606  AliErrorStream() << "Associated constituent cluster not found\n";
607  continue;
608  }
609  TLorentzVector clusterp;
610  clust->GetMomentum(clusterp, fVertex, (AliVCluster::VCluUserDefEnergy_t)clusters->GetDefaultClusterEnergy());
611 
612  num += clusterp.Pt() * clusterp.Vect().DrEtaPhi(jetvec);
613  den += clusterp.Pt();
614  }
615  }
616  return num/den;
617 }
618 
620  if (!(jet.GetNumberOfTracks() || jet.GetNumberOfClusters()))
621  throw SubstructureException();
622  Double_t den(0.), num(0.);
623  bool isMC = dynamic_cast<const AliMCParticleContainer *>(particles);
624  if(particles && (fUseChargedConstituents || isMC)){
625  for(UInt_t itrk = 0; itrk < jet.GetNumberOfTracks(); itrk++) {
626  auto trk = jet.TrackAt(itrk, particles->GetArray());
627  if(!trk){
628  AliErrorStream() << "Associated constituent particle / track not found\n";
629  continue;
630  }
631  if(!trk->Charge() && !fUseNeutralConstituents) continue; // Reject neutral constituents in case of using only charged consituents
632  if(trk->Charge() && !fUseChargedConstituents) continue; // Reject charged constituents in case of using only neutral consituents
633  num += trk->Pt() * trk->Pt();
634  den += trk->Pt();
635  }
636  }
637  if(clusters && fUseNeutralConstituents){
638  for(UInt_t icl = 0; icl < jet.GetNumberOfClusters(); icl++){
639  auto clust = jet.ClusterAt(icl, clusters->GetArray());
640  if(!clust) {
641  AliErrorStream() << "Associated constituent cluster not found\n";
642  continue;
643  }
644  TLorentzVector clusterp;
645  clust->GetMomentum(clusterp, fVertex, (AliVCluster::VCluUserDefEnergy_t)clusters->GetDefaultClusterEnergy());
646  num += clusterp.Pt() * clusterp.Pt();
647  den += clusterp.Pt();
648  }
649  }
650  return TMath::Sqrt(num)/den;
651 }
652 
654  for(int icl = 0; icl < jet->GetNumberOfClusters(); icl++){
655  auto clust = jet->ClusterAt(icl, clusters->GetArray());
656  AliDebugStream(3) << "cluster time " << clust->GetTOF() << std::endl;
657  fQAHistos->FillTH2("hClusterConstE", jet->Pt(),clust->GetUserDefEnergy(clusters->GetDefaultClusterEnergy()));
658  fQAHistos->FillTH2("hClusterConstTime", jet->Pt(), clust->GetTOF());
659  fQAHistos->FillTH2("hClusterConstM02", jet->Pt(), clust->GetM02());
660  fQAHistos->FillTH2("hClusterConstNcell", jet->Pt(), clust->GetNCells());
661 
662 #ifdef EXPERIMENTAL_JETCONSTITUENTS
663  fQAHistos->FillTH2("hClusterIndexENLC", jet->Pt(), clust->GetNonLinCorrEnergy());
664  fQAHistos->FillTH2("hClusterIndexEHC", jet->Pt(), clust->GetHadCorrEnergy());
665 #endif
666  }
667 
668 #ifdef EXPERIMENTAL_JETCONSTITUENTS
669  // Loop over charged particles - fill test histogram
670  for(int itrk = 0; itrk < jet->GetNumberOfTracks(); itrk++){
671  auto part = jet->TrackAt(itrk, cont->GetArray());
672  fQAHistos->FillTH2("hChargedIndexPt", jet->Pt(), part->Pt());
673  }
674 
675  // Look over charged constituents
676  AliDebugStream(2) << "Jet: Number of particle constituents: " << jet->GetParticleConstituents().GetEntriesFast() << std::endl;
677  for(auto pconst : jet->GetParticleConstituents()) {
678  auto part = static_cast<PWG::JETFW::AliEmcalParticleJetConstituent *>(pconst);
679  AliDebugStream(3) << "Found particle constituent with pt " << part->Pt() << ", from VParticle " << part->GetParticle()->Pt() << std::endl;
680  fQAHistos->FillTH2("hChargedConstituentPt", jet->Pt(), part->Pt());
681  }
682 
683  // Loop over neutral constituents
684  AliDebugStream(2) << "Jet: Number of cluster constituents: " << jet->GetClusterConstituents().GetEntriesFast() << std::endl;
685  for(auto cconst : jet->GetClusterConstituents()){
686  auto clust = static_cast<PWG::JETFW::AliEmcalClusterJetConstituent *>(cconst);
687  AliDebugStream(3) << "Found cluster constituent with energy " << clust->E() << " using energy definition " << static_cast<int>(clust->GetDefaultEnergyType()) << std::endl;
688  fQAHistos->FillTH2("hClusterConstituentEDefault", jet->Pt(), clust->E());
689  fQAHistos->FillTH2("hClusterConstituentENLC", jet->Pt(), clust->GetCluster()->GetNonLinCorrEnergy());
690  fQAHistos->FillTH2("hClusterConstituentEHC", jet->Pt(), clust->GetCluster()->GetHadCorrEnergy());
691  }
692 #endif
693 }
694 
696  int ncharged = 0, nneutral = jet.GetNumberOfClusters();
697  if(particles) {
698  for(decltype(jet.GetNumberOfTracks()) ipart = 0; ipart < jet.GetNumberOfTracks(); ipart++){
699  auto part = jet.TrackAt(ipart, particles->GetArray());
700  if(!part) continue;
701  if(part->Charge()) ncharged++;
702  else nneutral++;
703  }
704  }
705  // check if the jet has at least one consituent for jet substructure
706  int nallowed = 0;
707  nallowed += fUseChargedConstituents ? ncharged : 0;
708  nallowed += fUseNeutralConstituents ? nneutral : 0;
709  return nallowed > 0;
710 }
711 
712 std::vector<Triggerinfo> AliAnalysisTaskEmcalJetSubstructureTree::DecodeTriggerString(const std::string &triggerstring) const {
713  std::vector<Triggerinfo> result;
714  std::stringstream triggerparser(triggerstring);
715  std::string currenttrigger;
716  while(std::getline(triggerparser, currenttrigger, ' ')){
717  if(!currenttrigger.length()) continue;
718  std::vector<std::string> tokens;
719  std::stringstream triggerdecoder(currenttrigger);
720  std::string token;
721  while(std::getline(triggerdecoder, token, '-')) tokens.emplace_back(token);
722  result.emplace_back(Triggerinfo({tokens[0], tokens[1], tokens[2], tokens[3]}));
723  }
724  return result;
725 }
726 
727 std::string AliAnalysisTaskEmcalJetSubstructureTree::MatchTrigger(const std::string &triggertoken) const {
728  std::vector<std::string> tokens;
729  std::string result;
730  std::stringstream decoder(fInputEvent->GetFiredTriggerClasses().Data());
731  while(std::getline(decoder, result, ' ')) tokens.emplace_back(result);
732  result.clear();
733  for(auto t : tokens) {
734  if(t.find(triggertoken) != std::string::npos) {
735  // take first occurrence - downscale factor should normally be the same
736  result = t;
737  break;
738  }
739  }
740  return result;
741 }
742 
743 bool AliAnalysisTaskEmcalJetSubstructureTree::IsSelectEmcalTriggers(const std::string &triggerstring) const {
744  const std::array<std::string, 8> kEMCALTriggers = {
745  "EJ1", "EJ2", "DJ1", "DJ2", "EG1", "EG2", "DG1", "DG2"
746  };
747  bool isEMCAL = false;
748  for(auto emcaltrg : kEMCALTriggers) {
749  if(triggerstring.find(emcaltrg) != std::string::npos) {
750  isEMCAL = true;
751  break;
752  }
753  }
754  return isEMCAL;
755 }
756 
757 bool AliAnalysisTaskEmcalJetSubstructureTree::IsPartBranch(const std::string &branchname) const{
758  return (branchname.find("Sim") != std::string::npos) || (branchname.find("True") != std::string::npos);
759 }
760 
761 bool AliAnalysisTaskEmcalJetSubstructureTree::IsAcceptanceBranch(const std::string &branchname) const {
762  return (branchname.find("Eta") != std::string::npos) || (branchname.find("Phi") != std::string::npos);
763 }
764 
765 bool AliAnalysisTaskEmcalJetSubstructureTree::IsRhoBranch(const std::string &branchname) const{
766  return (branchname.find("Rho") != std::string::npos);
767 }
768 
769 bool AliAnalysisTaskEmcalJetSubstructureTree::IsMassBranch(const std::string &branchname) const{
770  return (branchname.find("Mass") != std::string::npos); // also disable rho mass branch
771 }
772 
773 bool AliAnalysisTaskEmcalJetSubstructureTree::IsSoftdropBranch(const std::string &branchname) const{
774  return (branchname.find("gMeasured") != std::string::npos) || (branchname.find("gTrue") != std::string::npos) || (branchname.find("NDropped") != std::string::npos);
775 }
776 
777 bool AliAnalysisTaskEmcalJetSubstructureTree::IsNSubjettinessBranch(const std::string &branchname) const{
778  return (branchname.find("Subjettiness") != std::string::npos);
779 }
780 
781 bool AliAnalysisTaskEmcalJetSubstructureTree::IsStructbranch(const std::string &branchname) const{
782  return (branchname.find("Angularity") != std::string::npos) || (branchname.find("PtD") != std::string::npos);
783 }
784 
786  AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
787 
788  Bool_t isAOD(kFALSE);
789  AliInputEventHandler *inputhandler = static_cast<AliInputEventHandler *>(mgr->GetInputEventHandler());
790  if(inputhandler) {
791  if(inputhandler->IsA() == AliAODInputHandler::Class()){
792  std::cout << "Analysing AOD events\n";
793  isAOD = kTRUE;
794  } else {
795  std::cout << "Analysing ESD events\n";
796  }
797  }
798 
799  std::stringstream taskname;
800  taskname << "JetSubstructureTreemaker_R" << std::setw(2) << std::setfill('0') << int(jetradius*10) << trigger;
802  mgr->AddTask(treemaker);
803  treemaker->SetMakeGeneralHistograms(kTRUE);
804 
805  // Adding containers
806  if(isMC) {
807  AliParticleContainer *particles = treemaker->AddMCParticleContainer("mcparticles");
808  particles->SetMinPt(0.);
809 
810  AliJetContainer *mcjets = treemaker->AddJetContainer(
811  jettype,
813  recombinationScheme,
814  jetradius,
816  particles, nullptr);
817  mcjets->SetName("mcjets");
818  mcjets->SetJetPtCut(20.);
819  }
820 
821  if(isData) {
822  AliTrackContainer *tracks(nullptr);
823  if((jettype == AliJetContainer::kChargedJet) || (jettype == AliJetContainer::kFullJet)){
825  std::cout << "Track container name: " << tracks->GetName() << std::endl;
826  tracks->SetMinPt(0.15);
827  }
828  AliClusterContainer *clusters(nullptr);
830  std::cout << "Using full or neutral jets ..." << std::endl;
832  std::cout << "Cluster container name: " << clusters->GetName() << std::endl;
833  clusters->SetClusHadCorrEnergyCut(0.3); // 300 MeV E-cut
834  clusters->SetDefaultClusterEnergy(AliVCluster::kHadCorr);
835  } else {
836  std::cout << "Using charged jets ... " << std::endl;
837  }
838 
839  AliJetContainer *datajets = treemaker->AddJetContainer(
840  jettype,
842  recombinationScheme,
843  jetradius,
845  tracks, clusters);
846  datajets->SetName("datajets");
847  datajets->SetJetPtCut(20.);
848 
849  treemaker->SetUseAliAnaUtils(true, true);
850  treemaker->SetVzRange(-10., 10);
851 
852  // configure trigger selection
853  std::string triggerstring(trigger);
854  if(triggerstring.find("INT7") != std::string::npos) {
855  treemaker->SetTriggerBits(AliVEvent::kINT7);
856  } else if(triggerstring.find("EJ1") != std::string::npos) {
857  treemaker->SetTriggerBits(AliVEvent::kEMCEJE);
858  treemaker->SetTriggerString("EJ1");
859  } else if(triggerstring.find("EJ2") != std::string::npos) {
860  treemaker->SetTriggerBits(AliVEvent::kEMCEJE);
861  treemaker->SetTriggerString("EJ2");
862  }
863  }
864 
865  std::string jettypestring;
866  switch(jettype) {
867  case AliJetContainer::kFullJet: jettypestring = "FullJets"; break;
868  case AliJetContainer::kChargedJet: jettypestring = "ChargedJets"; break;
869  case AliJetContainer::kNeutralJet: jettypestring = "NeutralJets"; break;
870  default: jettypestring = "Undef";
871  };
872 
873  // Connecting containers
874  std::stringstream outputfile, histname, treename;
875  outputfile << mgr->GetCommonFileName() << ":JetSubstructure_" << jettypestring << "_R" << std::setw(2) << std::setfill('0') << int(jetradius * 10.) << "_" << trigger;
876  histname << "JetSubstructureHistos_" << jettypestring << "_R" << std::setw(2) << std::setfill('0') << int(jetradius * 10.) << "_" << trigger;
877  treename << "JetSubstructureTree_" << jettypestring << "_R" << std::setw(2) << std::setfill('0') << int(jetradius * 10.) << "_" << trigger;
878  mgr->ConnectInput(treemaker, 0, mgr->GetCommonInputContainer());
879  mgr->ConnectOutput(treemaker, 1, mgr->CreateContainer(histname.str().data(), AliEmcalList::Class(), AliAnalysisManager::kOutputContainer, outputfile.str().data()));
880  mgr->ConnectOutput(treemaker, 2, mgr->CreateContainer(treename.str().data(), TTree::Class(), AliAnalysisManager::kOutputContainer, mgr->GetCommonFileName()));
881 
882  return treemaker;
883 }
884 
885 std::string Triggerinfo::ExpandClassName() const {
886  std::string result = fTriggerClass + "-" + fBunchCrossing + "-" + fPastFutureProtection + "-" + fTriggerCluster;
887  return result;
888 }
889 
890 bool Triggerinfo::IsTriggerClass(const std::string &triggerclass) const {
891  return fTriggerClass.substr(1) == triggerclass; // remove C from trigger class part
892 }
893 
894 } /* namespace EmcalTriggerJets */
Double_t Area() const
Definition: AliEmcalJet.h:130
double Double_t
Definition: External.C:58
static AliAnalysisTaskEmcalJetSubstructureTree * AddEmcalJetSubstructureTreeMaker(Bool_t isMC, Bool_t isData, Double_t jetradius, AliJetContainer::EJetType_t jettype, AliJetContainer::ERecoScheme_t recombinationScheme, const char *name)
Class creating a linear binning, used in the histogram manager.
Double_t MakePtD(const AliEmcalJet &jet, const AliParticleContainer *const particles, const AliClusterContainer *const clusters) const
AliEmcalJet * ClosestJet() const
Definition: AliEmcalJet.h:327
AliJetContainer * GetJetContainer(Int_t i=0) const
virtual double E() const
Access to constituent energy.
Double_t Eta() const
Definition: AliEmcalJet.h:121
static AliEmcalDownscaleFactorsOCDB * Instance()
Double_t Py() const
Definition: AliEmcalJet.h:107
Double_t Phi() const
Definition: AliEmcalJet.h:117
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
void FillTH2(const char *hname, double x, double y, double weight=1., Option_t *opt="")
Fill a 2D histogram within the container.
Int_t ClusterAt(Int_t idx) const
Definition: AliEmcalJet.h:137
std::vector< Triggerinfo > DecodeTriggerString(const std::string &triggerstring) const
TCanvas * c
Definition: TestFitELoss.C:172
AliJetContainer * AddJetContainer(const char *n, TString defaultCutType, Float_t jetRadius=0.4)
const std::vector< PWG::JETFW::AliEmcalParticleJetConstituent > & GetParticleConstituents() const
Get container with particle (track / MC particle) constituents.
Definition: AliEmcalJet.h:184
Structure for results from the soft drop algorithm.
virtual void UserExecOnce()
Task initializations handled in user tasks.
Double_t E() const
Definition: AliEmcalJet.h:119
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.
Double_t GetDownscaleFactorForTriggerClass(const TString &trigger) const
void SetVzRange(Double_t min, Double_t max)
AliClusterContainer * AddClusterContainer(const char *n)
Create new cluster container and attach it to the task.
TString part
use mixed event to constrain combinatorial background
Definition: InvMassFit.C:52
fastjet::JetAlgorithm fRecluserAlgo
Reclusterization algorithm.
UShort_t GetNumberOfConstituents() const
Definition: AliEmcalJet.h:140
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 Px() const
Definition: AliEmcalJet.h:106
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
AliRhoParameter * GetRhoFromEvent(const char *name)
TH2 * CreateTH2(const char *name, const char *title, int nbinsx, double xmin, double xmax, int nbinsy, double ymin, double ymax, Option_t *opt="")
Create a new TH2 within the container.
int Int_t
Definition: External.C:63
UShort_t GetNumberOfClusters() const
Definition: AliEmcalJet.h:138
void SetJetPtCut(Float_t cut)
TPC acceptance.
Definition: AliEmcalJet.h:67
unsigned int UInt_t
Definition: External.C:33
THashList * GetListOfHistograms() const
Get the list of histograms.
Definition: THistManager.h:671
Int_t GetNJets() const
AliNSubjettinessParameters MakeNsubjettinessParameters(const fastjet::PseudoJet &jet, const AliNSubjettinessDefinition &cut) const
TH1 * CreateTH1(const char *name, const char *title, int nbins, double xmin, double xmax, Option_t *opt="")
Create a new TH1 within the container.
Definition for the algorithm obtaining the softdrop parameters.
AliClusterContainer * GetClusterContainer(Int_t i=0) const
Get cluster container attached to this task.
virtual void RunChanged(Int_t newrun)
Process tasks relevant when a file with a different run number is processed.
AliMCParticleContainer * AddMCParticleContainer(const char *n)
Create new container for MC particles and attach it to the task.
void FillTH1(const char *hname, double x, double weight=1., Option_t *opt="")
Fill a 1D histogram within the container.
Double_t MakeAngularity(const AliEmcalJet &jet, const AliParticleContainer *tracks, const AliClusterContainer *clusters) const
Implementation of a jet constituent for constituent clusters.
bool IsTriggerClass(const std::string &triggerclass) const
Double_t Pt() const
Definition: AliEmcalJet.h:109
Bool_t isMC
void DoConstituentQA(const AliEmcalJet *jet, const AliParticleContainer *tracks, const AliClusterContainer *clusters)
Float_t GetJetRadius() const
const std::vector< PWG::JETFW::AliEmcalClusterJetConstituent > & GetClusterConstituents() const
Get container with cluster constituents.
Definition: AliEmcalJet.h:190
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.
Handler for downscale factors for various triggers obtained from the OCDB.
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:51
Container class for histograms.
Definition: THistManager.h:99
Double_t Pz() const
Definition: AliEmcalJet.h:108
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
void LinkOutputBranch(const std::string &branchname, Double_t *datalocation)
Double_t NEF() const
Definition: AliEmcalJet.h:148
void SetDefaultClusterEnergy(Int_t d)
bool SelectJet(const AliEmcalJet &jet, const AliParticleContainer *particles) const
Container structure for EMCAL clusters.
Container for MC-true particles within the EMCAL framework.
Double_t M() const
Definition: AliEmcalJet.h:120
EMCal fiducial acceptance (each eta, phi edge narrowed by jet R)
Definition: AliEmcalJet.h:70
Container for jet within the EMCAL jet framework.
Bool_t fUseTriggerSelectionForData
Use trigger selection on data (require trigger patch in addition to trigger selection string) ...
TList * OpenFile(const char *fname)
Definition: DrawAnaELoss.C:65
void SetClusHadCorrEnergyCut(Double_t cut)
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.