AliPhysics  5364b50 (5364b50)
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  fGlobalTreeParams(nullptr),
86  fSoftDropMeasured(nullptr),
87  fSoftDropTrue(nullptr),
88  fNSubMeasured(nullptr),
89  fNSubTrue(nullptr),
90  fKineRec(nullptr),
91  fKineSim(nullptr),
92  fJetStructureMeasured(nullptr),
93  fJetStructureTrue(nullptr),
94  fQAHistos(nullptr),
95  fLumiMonitor(nullptr),
96  fSDZCut(0.1),
97  fSDBetaCut(0),
98  fReclusterizer(kCAAlgo),
99  fTriggerSelectionBits(AliVEvent::kAny),
100  fTriggerSelectionString(""),
101  fNameTriggerDecisionContainer("EmcalTriggerDecision"),
102  fUseTriggerSelectionForData(false),
103  fUseDownscaleWeight(false),
104  fUseChargedConstituents(true),
105  fUseNeutralConstituents(true),
106  fFillPart(true),
107  fFillRho(true),
108  fFillSoftDrop(true),
109  fFillNSub(true),
110  fFillStructGlob(true)
111 {
112 }
113 
115  AliAnalysisTaskEmcalJet(name, kTRUE),
122  fKineRec(nullptr),
123  fKineSim(nullptr),
128  fSDZCut(0.1),
129  fSDBetaCut(0),
131  fTriggerSelectionBits(AliVEvent::kAny),
133  fNameTriggerDecisionContainer("EmcalTriggerDecision"),
135  fUseDownscaleWeight(false),
138  fFillPart(true),
139  fFillRho(true),
140  fFillSoftDrop(true),
141  fFillNSub(true),
142  fFillStructGlob(true)
143 {
144  DefineOutput(2, TTree::Class());
145 }
146 
150  if(fSoftDropTrue) delete fSoftDropTrue;
151  if(fNSubMeasured) delete fNSubMeasured;
152  if(fNSubTrue) delete fNSubTrue;
153  if(fKineRec) delete fKineRec;
154  if(fKineSim) delete fKineSim;
157 }
158 
161 
162  // Make QA for constituent clusters
163  TLinearBinning jetptbinning(9, 20, 200),
164  clusterenergybinning(200, 0., 200),
165  timebinning(1000, -500., 500.),
166  m02binning(100, 0., 1.),
167  ncellbinning(101, -0.5, 100.5);
168  fQAHistos = new THistManager("QAhistos");
169  fQAHistos->CreateTH1("hEventCounter", "Event counter", 1, 0.5, 1.5);
170  fQAHistos->CreateTH2("hClusterConstE", "EMCAL cluster energy vs jet pt; p_{t, jet} (GeV/c); E_{cl} (GeV)", jetptbinning, clusterenergybinning);
171  fQAHistos->CreateTH2("hClusterConstTime", "EMCAL cluster time vs. jet pt; p_{t, jet} (GeV/c); t_{cl} (ns)", jetptbinning, timebinning);
172  fQAHistos->CreateTH2("hClusterConstM02", "EMCAL cluster M02 vs. jet pt; p{t, jet} (GeV/c); M02", jetptbinning, m02binning);
173  fQAHistos->CreateTH2("hClusterConstNcell", "EMCAL cluster ncell vs. jet pt; p{t, jet} (GeV/c); Number of cells", jetptbinning, ncellbinning);
174 
175  // Test of constituent QA
176 #ifdef EXPERIMENTAL_JETCONSTITUENTS
177  fQAHistos->CreateTH2("hChargedConstituentPt", "charged constituent pt vs jet pt (via constituent map); p_{t,jet} (GeV/c); p_{t,ch} (GeV/c)", jetptbinning, clusterenergybinning);
178  fQAHistos->CreateTH2("hChargedIndexPt", "charged constituent pt vs jet pt (via index map); p_{t, jet} (GeV/c); p_{t, ch} (GeV/c)", jetptbinning, clusterenergybinning);
179 
180  fQAHistos->CreateTH2("hClusterConstituentEDefault", "cluster constituent default energy vs. jet pt (va constituent map); p_{t, jet} (GeV/c); E_{cl} (GeV)", jetptbinning, clusterenergybinning);
181  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);
182  fQAHistos->CreateTH2("hClusterConstituentEHC", "cluster constituent hadronic-corrected energy vs. jet pt (va constituent map); p_{t, jet} (GeV/c); E_{cl} (GeV)", jetptbinning, clusterenergybinning);
183  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);
184  fQAHistos->CreateTH2("hClusterIndexEHC", "cluster constituent hadronic-corrected energy vs. jet pt (via index map); p_{t, jet} (GeV/c); E_{cl} (GeV)", jetptbinning, clusterenergybinning);
185 #endif
186  for(auto h : *(fQAHistos->GetListOfHistograms())) fOutput->Add(h);
187 
188  OpenFile(2);
189  std::string treename = this->GetOutputSlot(2)->GetContainer()->GetName();
190  fJetSubstructureTree = new TTree(treename.data(), "Tree with jet substructure information");
191 
193  fGlobalTreeParams->LinkJetTreeBranches(fJetSubstructureTree, fFillRho);
196  if(fFillPart) {
199  }
200  if(fFillSoftDrop) {
203  if(fFillPart) {
206  }
207  }
208  if(fFillNSub) {
211  if(fFillPart) {
214  }
215  }
216 
217  if(fFillStructGlob){
220  if(fFillPart){
223  }
224  }
225 
226  PostData(1, fOutput);
227  PostData(2, fJetSubstructureTree);
228 }
229 
233  }
234 }
235 
239  AliParticleContainer *particles = GetParticleContainer("mcparticles");
240 
241  AliJetContainer *mcjets = GetJetContainer("mcjets");
242  AliJetContainer *datajets = GetJetContainer("datajets");
243 
244  FillLuminosity(); // Makes only sense in data
245 
246  // for(auto e : *(fInputEvent->GetList())) std::cout << e->GetName() << std::endl;
247 
248  std::stringstream rhoTagData, rhoTagMC;
249  if(datajets) rhoTagData << "R" << std::setw(2) << std::setfill('0') << static_cast<Int_t>(datajets->GetJetRadius() * 10.);
250  if(mcjets) rhoTagMC << "R" << std::setw(2) << std::setfill('0') << static_cast<Int_t>(mcjets->GetJetRadius() * 10.);
251 
252  if(fFillRho){
253  std::string rhoSparseData = "RhoSparse_Full_" + rhoTagData.str(), rhoSparseMC = "RhoSparse_Full_" + rhoTagMC.str(),
254  rhoMassData = "RhoMassSparse_Full_" + rhoTagData.str(), rhoMassMC = "RhoMassSparse_Full_" + rhoTagMC.str();
255  AliRhoParameter *rhoPtRec = GetRhoFromEvent(rhoSparseData.data()),
256  *rhoMassRec = GetRhoFromEvent(rhoMassData.data()),
257  *rhoPtSim = GetRhoFromEvent(rhoSparseMC.data()),
258  *rhoMassSim = GetRhoFromEvent(rhoMassMC.data());
259  AliDebugStream(2) << "Found rho parameter for reconstructed pt: " << (rhoPtRec ? "yes" : "no") << ", value: " << (rhoPtRec ? rhoPtRec->GetVal() : 0.) << std::endl;
260  AliDebugStream(2) << "Found rho parameter for sim pt: " << (rhoPtSim ? "yes" : "no") << ", value: " << (rhoPtSim ? rhoPtSim->GetVal() : 0.) << std::endl;
261  AliDebugStream(2) << "Found rho parameter for reconstructed Mass: " << (rhoMassRec ? "yes" : "no") << ", value: " << (rhoMassRec ? rhoMassRec->GetVal() : 0.) << std::endl;
262  AliDebugStream(2) << "Found rho parameter for sim Mass: " << (rhoMassSim ? "yes" : "no") << ", value: " << (rhoMassSim ? rhoMassSim->GetVal() : 0.) << std::endl;
263  Double_t rhopars[4] = {
264  rhoPtRec ? rhoPtRec->GetVal() : 0.,
265  rhoPtSim ? rhoPtSim->GetVal() : 0.,
266  rhoMassRec ? rhoMassRec->GetVal() : 0.,
267  rhoMassSim ? rhoMassSim->GetVal() : 0.
268  };
269  memcpy(this->fGlobalTreeParams->fRhoParamters, rhopars, sizeof(Double_t) * 4);
270  }
271 
272  AliDebugStream(1) << "Inspecting jet radius " << (datajets ? datajets->GetJetRadius() : mcjets->GetJetRadius()) << std::endl;
273  this->fGlobalTreeParams->fJetRadius = (datajets ? datajets->GetJetRadius() : mcjets->GetJetRadius());
274 
275  // Run trigger selection (not on pure MCgen train)
276  if(datajets){
277  if(!(fInputHandler->IsEventSelected() & fTriggerSelectionBits)) return false;
278  if(!mcjets){
279  // Pure data - do EMCAL trigger selection from selection string
280  if(fTriggerSelectionString.Length()) {
281  if(!fInputEvent->GetFiredTriggerClasses().Contains(fTriggerSelectionString)) return false;
283  auto trgselresult = static_cast<PWG::EMCAL::AliEmcalTriggerDecisionContainer *>(fInputEvent->FindListObject(fNameTriggerDecisionContainer));
284  AliDebugStream(1) << "Found trigger decision object: " << (trgselresult ? "yes" : "no") << std::endl;
285  if(!trgselresult){
286  AliErrorStream() << "Trigger decision container with name " << fNameTriggerDecisionContainer << " not found in event - not possible to select EMCAL triggers" << std::endl;
287  return false;
288  }
289  if(!trgselresult->IsEventSelected(fTriggerSelectionString)) return false;
290  }
291  }
292  } else {
294  // Simulation - do EMCAL trigger selection from trigger selection object
295  auto mctrigger = static_cast<PWG::EMCAL::AliEmcalTriggerDecisionContainer *>(fInputEvent->FindListObject(fNameTriggerDecisionContainer));
296  AliDebugStream(1) << "Found trigger decision object: " << (mctrigger ? "yes" : "no") << std::endl;
297  if(!mctrigger){
298  AliErrorStream() << "Trigger decision container with name " << fNameTriggerDecisionContainer << " not found in event - not possible to select EMCAL triggers" << std::endl;
299  return false;
300  }
301  if(!mctrigger->IsEventSelected(fTriggerSelectionString)) return false;
302  }
303  }
304  }
305 
306  double weight = 1.;
308  AliDebugStream(2) << "Trigger selection string: " << fTriggerSelectionString << std::endl;
309  TString selectionString = (fTriggerSelectionBits & AliVEvent::kINT7) ? "INT7" : fTriggerSelectionString;
310  auto triggerstring = MatchTrigger(selectionString.Data());
311  AliDebugStream(2) << "Getting downscale correction factor for trigger string " << triggerstring << std::endl;
313  }
314  AliDebugStream(1) << "Using downscale weight " << weight << std::endl;
315  this->fGlobalTreeParams->fEventWeight = weight;
316 
317 
318  // Count events (for spectrum analysis)
319  fQAHistos->FillTH1("hEventCounter", 1);
320 
321  AliSoftdropDefinition softdropSettings;
322  softdropSettings.fBeta = fSDBetaCut;
323  softdropSettings.fZ = fSDZCut;
324  switch(fReclusterizer) {
325  case kCAAlgo: softdropSettings.fRecluserAlgo = fastjet::cambridge_aachen_algorithm; break;
326  case kKTAlgo: softdropSettings.fRecluserAlgo = fastjet::kt_algorithm; break;
327  case kAKTAlgo: softdropSettings.fRecluserAlgo = fastjet::antikt_algorithm; break;
328  };
329 
330  AliNSubjettinessDefinition nsubjettinessSettings;
331  nsubjettinessSettings.fBeta = 1.;
332  nsubjettinessSettings.fRadius = 0.4;
333 
334  if(datajets) {
335  AliDebugStream(1) << "In data jets branch: found " << datajets->GetNJets() << " jets, " << datajets->GetNAcceptedJets() << " were accepted\n";
336  AliDebugStream(1) << "Having MC information: " << (mcjets ? TString::Format("yes, with %d jets", mcjets->GetNJets()) : "no") << std::endl;
337  if(mcjets) {
338  AliDebugStream(1) << "In MC jets branch: found " << mcjets->GetNJets() << " jets, " << mcjets->GetNAcceptedJets() << " were accepted\n";
339  }
340  for(auto jet : datajets->accepted()) {
341  double pt = jet->Pt(), pz = jet->Pz(), E = jet->E(), M = TMath::Sqrt(E*E - pt*pt - pz*pz);
342  AliDebugStream(2) << "Next jet: pt:" << jet->Pt() << ", E: " << E << ", pz: " << pz << ", M(self): " << M << "M(fj)" << jet->M() << std::endl;
343  AliEmcalJet *associatedJet = jet->ClosestJet();
344 
345  if(mcjets) {
346  if(!associatedJet) {
347  AliDebugStream(2) << "Not found associated jet" << std::endl;
348  continue;
349  }
350  if(!(SelectJet(*jet, tracks) && SelectJet(*associatedJet, particles))) continue;
351  try {
352  DoConstituentQA(jet, tracks, clusters);
353  AliJetSubstructureData structureData = MakeJetSubstructure(*jet, datajets->GetJetRadius() * 2., tracks, clusters, {softdropSettings, nsubjettinessSettings}),
354  structureMC = fFillPart ? MakeJetSubstructure(*associatedJet, mcjets->GetJetRadius() * 2, particles, nullptr, {softdropSettings, nsubjettinessSettings}) : AliJetSubstructureData();
356  if(fKineSim) *fKineSim = MakeJetKineParameters(*associatedJet);
357  if(fSoftDropMeasured) *fSoftDropMeasured = structureData.fSoftDrop;
358  if(fSoftDropTrue) *fSoftDropTrue = structureMC.fSoftDrop;
359  if(fNSubMeasured) *fNSubMeasured = structureData.fNsubjettiness;
360  if(fNSubTrue) *fNSubTrue = structureMC.fNsubjettiness;
361  if(fJetStructureMeasured) *fJetStructureMeasured = {MakeAngularity(*jet, tracks, clusters), MakePtD(*jet, tracks, clusters)};
362  if(fJetStructureTrue) *fJetStructureTrue = {MakeAngularity(*associatedJet, particles, nullptr), MakePtD(*associatedJet, particles, nullptr)};
363  fJetSubstructureTree->Fill();
364  } catch(ReclusterizerException &e) {
365  AliErrorStream() << "Error in reclusterization - skipping jet" << std::endl;
366  } catch(SubstructureException &e) {
367  AliErrorStream() << "Error in substructure observable - skipping jet" << std::endl;
368  }
369  } else {
370  if(!SelectJet(*jet, tracks)) continue;
371  try {
372  DoConstituentQA(jet, tracks, clusters);
373  AliJetSubstructureData structure = MakeJetSubstructure(*jet, 0.4, tracks, clusters, {softdropSettings, nsubjettinessSettings});
376  if(fNSubMeasured) *fNSubMeasured = structure.fNsubjettiness;
377  if(fJetStructureMeasured) *fJetStructureMeasured = {MakeAngularity(*jet, tracks, clusters), MakePtD(*jet, tracks, clusters)};
378  fJetSubstructureTree->Fill();
379  } catch(ReclusterizerException &e) {
380  AliErrorStream() << "Error in reclusterization - skipping jet" << std::endl;
381  } catch(SubstructureException &e) {
382  AliErrorStream() << "Error in substructure observable - skipping jet" << std::endl;
383  }
384  }
385  }
386  } else {
387  if(mcjets) {
388  // for MCgen train
389  AliDebugStream(1) << "In MC pure jet branch: found " << mcjets->GetNJets() << " jets, " << mcjets->GetNAcceptedJets() << " were accepted\n";
390  for(auto j : mcjets->accepted()){
391  AliEmcalJet *mcjet = static_cast<AliEmcalJet *>(j);
392  try {
393  AliJetSubstructureData structure = MakeJetSubstructure(*mcjet, mcjets->GetJetRadius() * 2., particles, nullptr,{softdropSettings, nsubjettinessSettings});
394  if(this->fKineSim) *fKineSim = MakeJetKineParameters(*mcjet);
395  if(fSoftDropTrue) *fSoftDropTrue = structure.fSoftDrop;
396  if(fNSubTrue) *fNSubTrue = structure.fNsubjettiness;
397  if(fJetStructureTrue) *fJetStructureTrue = {MakeAngularity(*mcjet, particles, nullptr), MakePtD(*mcjet, particles, nullptr)};
398  fJetSubstructureTree->Fill();
399  } catch (ReclusterizerException &e) {
400  AliErrorStream() << "Error in reclusterization - skipping jet" << std::endl;
401  } catch (SubstructureException &e) {
402  AliErrorStream() << "Error in substructure observable - skipping jet" << std::endl;
403  }
404  }
405  }
406  }
407 
408  return true;
409 }
410 
412  AliCDBManager * cdb = AliCDBManager::Instance();
413  if(!fMCEvent && cdb){
414  // Get List of trigger clusters
415  AliCDBEntry *en = cdb->Get("GRP/CTP/Config");
416  AliTriggerConfiguration *trg = static_cast<AliTriggerConfiguration *>(en->GetObject());
417  std::vector<std::string> clusternames;
418  for(auto c : trg->GetClusters()) {
419  AliTriggerCluster *clust = static_cast<AliTriggerCluster *>(c);
420  std::string clustname = clust->GetName();
421  auto iscent = clustname.find("CENT") != std::string::npos, iscalo = clustname.find("CALO") != std::string::npos;
422  if(!(iscalo || iscent)) continue;
423  clusternames.emplace_back(clustname);
424  }
425 
426  // Set the x-axis of the luminosity monitor histogram
427  fLumiMonitor = new TH1F("hLumiMonitor", "Luminosity monitor", clusternames.size(), 0, clusternames.size());
428  int currentbin(1);
429  for(auto c : clusternames) {
430  fLumiMonitor->GetXaxis()->SetBinLabel(currentbin++, c.data());
431  }
432  fOutput->Add(fLumiMonitor);
433  }
434 }
435 
438  auto downscalefactors = PWG::EMCAL::AliEmcalDownscaleFactorsOCDB::Instance();
439  if(fInputEvent->GetFiredTriggerClasses().Contains("INT7")) {
440  for(auto trigger : DecodeTriggerString(fInputEvent->GetFiredTriggerClasses().Data())){
441  auto int7trigger = trigger.IsTriggerClass("INT7");
442  auto bunchcrossing = trigger.fBunchCrossing == "B";
443  auto nopf = trigger.fPastFutureProtection == "NOPF";
444  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;
445  if(int7trigger && bunchcrossing && nopf) {
446  double downscale = downscalefactors->GetDownscaleFactorForTriggerClass(trigger.ExpandClassName());
447  AliDebugStream(5) << "Using downscale " << downscale << std::endl;
448  fLumiMonitor->Fill(trigger.fTriggerCluster.data(), 1./downscale);
449  }
450  }
451  }
452  }
453 }
454 
456  AliJetKineParameters result;
457  result.fPt = TMath::Abs(jet.Pt());
458  result.fE = jet.E();
459  result.fEta = jet.Eta();
460  result.fPhi = jet.Phi();
461  result.fArea = jet.Area();
462  result.fMass = jet.M();
463  result.fNEF = jet.NEF();
464  result.fNCharged = jet.GetNumberOfTracks();
465  result.fNNeutral = jet.GetNumberOfClusters();
466  return result;
467 }
468 
470  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
471  std::vector<fastjet::PseudoJet> constituents;
472  bool isMC = dynamic_cast<const AliMCParticleContainer *>(tracks);
473  AliDebugStream(2) << "Make new jet substrucutre for " << (isMC ? "MC" : "data") << " jet: Number of tracks " << jet.GetNumberOfTracks() << ", clusters " << jet.GetNumberOfClusters() << std::endl;
474  if(tracks && (fUseChargedConstituents || isMC)){ // Neutral particles part of particle container in case of MC
475  for(int itrk = 0; itrk < jet.GetNumberOfTracks(); itrk++){
476  auto track = jet.TrackAt(itrk, tracks->GetArray());
477  if(!track->Charge() && !fUseNeutralConstituents) continue; // Reject neutral constituents in case of using only charged consituents
478  if(track->Charge() && !fUseChargedConstituents) continue; // Reject charged constituents in case of using only neutral consituents
479  fastjet::PseudoJet constituentTrack(track->Px(), track->Py(), track->Pz(), track->E());
480  constituentTrack.set_user_index(jet.TrackAt(itrk));
481  constituents.push_back(constituentTrack);
482  }
483  }
484 
485  if(clusters && fUseNeutralConstituents){
486  for(int icl = 0; icl < jet.GetNumberOfClusters(); icl++) {
487  auto cluster = jet.ClusterAt(icl, clusters->GetArray());
488  TLorentzVector clustervec;
489  cluster->GetMomentum(clustervec, fVertex, (AliVCluster::VCluUserDefEnergy_t)clusters->GetDefaultClusterEnergy());
490  fastjet::PseudoJet constituentCluster(clustervec.Px(), clustervec.Py(), clustervec.Pz(), cluster->GetHadCorrEnergy());
491  constituentCluster.set_user_index(jet.ClusterAt(icl) + kClusterOffset);
492  constituents.push_back(constituentCluster);
493  }
494  }
495 
496  AliDebugStream(3) << "Found " << constituents.size() << " constituents for jet with pt=" << jet.Pt() << " GeV/c" << std::endl;
497  if(!constituents.size()){
498  AliErrorStream() << "Jet has 0 constituents." << std::endl;
499  throw ReclusterizerException();
500  }
501  // Redo jet finding on constituents with a
502  fastjet::JetDefinition jetdef(fastjet::antikt_algorithm, jetradius*2, static_cast<fastjet::RecombinationScheme>(0), fastjet::BestFJ30 );
503  std::vector<fastjet::PseudoJet> outputjets;
504  try {
505  fastjet::ClusterSequence jetfinder(constituents, jetdef);
506  outputjets = jetfinder.inclusive_jets(0);
508  return result;
509  } catch (fastjet::Error &e) {
510  AliErrorStream() << " FJ Exception caught: " << e.message() << std::endl;
511  throw ReclusterizerException();
512  } catch (SoftDropException &e) {
513  AliErrorStream() << "Softdrop exception caught: " << e.what() << std::endl;
514  throw ReclusterizerException();
515  }
516 }
517 
519  fastjet::contrib::SoftDrop softdropAlgorithm(cutparameters.fBeta, cutparameters.fZ);
520  softdropAlgorithm.set_verbose_structure(kTRUE);
521  std::unique_ptr<fastjet::contrib::Recluster> reclusterizer(new fastjet::contrib::Recluster(cutparameters.fRecluserAlgo, 1, true));
522  softdropAlgorithm.set_reclustering(kTRUE, reclusterizer.get());
523  AliDebugStream(4) << "Jet has " << jet.constituents().size() << " constituents" << std::endl;
524  auto groomed = softdropAlgorithm(jet);
525  try {
526  auto softdropstruct = groomed.structure_of<fastjet::contrib::SoftDrop>();
527 
528  AliSoftDropParameters result({softdropstruct.symmetry(),
529  groomed.m(),
530  softdropstruct.delta_R(),
531  groomed.perp(),
532  softdropstruct.delta_R(),
533  softdropstruct.mu(),
534  softdropstruct.dropped_count()});
535  return result;
536  } catch(std::bad_cast &e) {
537  throw SoftDropException();
538  }
539 }
540 
543  fastjet::contrib::Nsubjettiness (1,fastjet::contrib::KT_Axes(),fastjet::contrib::NormalizedCutoffMeasure(cut.fBeta, cut.fRadius, 1e100)).result(jet),
544  fastjet::contrib::Nsubjettiness (2,fastjet::contrib::KT_Axes(),fastjet::contrib::NormalizedCutoffMeasure(cut.fBeta, cut.fRadius, 1e100)).result(jet)
545  });
546  return result;
547 }
548 
550  if(!(jet.GetNumberOfTracks() || jet.GetNumberOfClusters()))
551  throw SubstructureException();
552  TVector3 jetvec(jet.Px(), jet.Py(), jet.Pz());
553  Double_t den(0.), num(0.);
554  bool isMC = dynamic_cast<const AliMCParticleContainer *>(tracks);
555  if(tracks && (fUseChargedConstituents || isMC)){
556  for(UInt_t itrk = 0; itrk < jet.GetNumberOfTracks(); itrk++) {
557  auto track = jet.TrackAt(itrk, tracks->GetArray());
558  if(!track){
559  AliErrorStream() << "Associated constituent particle / track not found\n";
560  continue;
561  }
562  if(!track->Charge() && !fUseNeutralConstituents) continue; // Reject neutral constituents in case of using only charged consituents
563  if(track->Charge() && !fUseChargedConstituents) continue; // Reject charged constituents in case of using only neutral consituents
564  TVector3 trackvec(track->Px(), track->Py(), track->Pz());
565 
566  num += track->Pt() * trackvec.DrEtaPhi(jetvec);
567  den += +track->Pt();
568  }
569  }
570  if(clusters && fUseNeutralConstituents) {
571  for(UInt_t icl = 0; icl < jet.GetNumberOfClusters(); icl++){
572  auto clust = jet.ClusterAt(icl, clusters->GetArray());
573  if(!clust) {
574  AliErrorStream() << "Associated constituent cluster not found\n";
575  continue;
576  }
577  TLorentzVector clusterp;
578  clust->GetMomentum(clusterp, fVertex, (AliVCluster::VCluUserDefEnergy_t)clusters->GetDefaultClusterEnergy());
579 
580  num += clusterp.Pt() * clusterp.Vect().DrEtaPhi(jetvec);
581  den += clusterp.Pt();
582  }
583  }
584  return num/den;
585 }
586 
588  if (!(jet.GetNumberOfTracks() || jet.GetNumberOfClusters()))
589  throw SubstructureException();
590  Double_t den(0.), num(0.);
591  bool isMC = dynamic_cast<const AliMCParticleContainer *>(particles);
592  if(particles && (fUseChargedConstituents || isMC)){
593  for(UInt_t itrk = 0; itrk < jet.GetNumberOfTracks(); itrk++) {
594  auto trk = jet.TrackAt(itrk, particles->GetArray());
595  if(!trk){
596  AliErrorStream() << "Associated constituent particle / track not found\n";
597  continue;
598  }
599  if(!trk->Charge() && !fUseNeutralConstituents) continue; // Reject neutral constituents in case of using only charged consituents
600  if(trk->Charge() && !fUseChargedConstituents) continue; // Reject charged constituents in case of using only neutral consituents
601  num += trk->Pt() * trk->Pt();
602  den += trk->Pt();
603  }
604  }
605  if(clusters && fUseNeutralConstituents){
606  for(UInt_t icl = 0; icl < jet.GetNumberOfClusters(); icl++){
607  auto clust = jet.ClusterAt(icl, clusters->GetArray());
608  if(!clust) {
609  AliErrorStream() << "Associated constituent cluster not found\n";
610  continue;
611  }
612  TLorentzVector clusterp;
613  clust->GetMomentum(clusterp, fVertex, (AliVCluster::VCluUserDefEnergy_t)clusters->GetDefaultClusterEnergy());
614  num += clusterp.Pt() * clusterp.Pt();
615  den += clusterp.Pt();
616  }
617  }
618  return TMath::Sqrt(num)/den;
619 }
620 
622  for(int icl = 0; icl < jet->GetNumberOfClusters(); icl++){
623  auto clust = jet->ClusterAt(icl, clusters->GetArray());
624  AliDebugStream(3) << "cluster time " << clust->GetTOF() << std::endl;
625  fQAHistos->FillTH2("hClusterConstE", jet->Pt(),clust->GetUserDefEnergy(clusters->GetDefaultClusterEnergy()));
626  fQAHistos->FillTH2("hClusterConstTime", jet->Pt(), clust->GetTOF());
627  fQAHistos->FillTH2("hClusterConstM02", jet->Pt(), clust->GetM02());
628  fQAHistos->FillTH2("hClusterConstNcell", jet->Pt(), clust->GetNCells());
629 
630 #ifdef EXPERIMENTAL_JETCONSTITUENTS
631  fQAHistos->FillTH2("hClusterIndexENLC", jet->Pt(), clust->GetNonLinCorrEnergy());
632  fQAHistos->FillTH2("hClusterIndexEHC", jet->Pt(), clust->GetHadCorrEnergy());
633 #endif
634  }
635 
636 #ifdef EXPERIMENTAL_JETCONSTITUENTS
637  // Loop over charged particles - fill test histogram
638  for(int itrk = 0; itrk < jet->GetNumberOfTracks(); itrk++){
639  auto part = jet->TrackAt(itrk, cont->GetArray());
640  fQAHistos->FillTH2("hChargedIndexPt", jet->Pt(), part->Pt());
641  }
642 
643  // Look over charged constituents
644  AliDebugStream(2) << "Jet: Number of particle constituents: " << jet->GetParticleConstituents().GetEntriesFast() << std::endl;
645  for(auto pconst : jet->GetParticleConstituents()) {
646  auto part = static_cast<PWG::JETFW::AliEmcalParticleJetConstituent *>(pconst);
647  AliDebugStream(3) << "Found particle constituent with pt " << part->Pt() << ", from VParticle " << part->GetParticle()->Pt() << std::endl;
648  fQAHistos->FillTH2("hChargedConstituentPt", jet->Pt(), part->Pt());
649  }
650 
651  // Loop over neutral constituents
652  AliDebugStream(2) << "Jet: Number of cluster constituents: " << jet->GetClusterConstituents().GetEntriesFast() << std::endl;
653  for(auto cconst : jet->GetClusterConstituents()){
654  auto clust = static_cast<PWG::JETFW::AliEmcalClusterJetConstituent *>(cconst);
655  AliDebugStream(3) << "Found cluster constituent with energy " << clust->E() << " using energy definition " << static_cast<int>(clust->GetDefaultEnergyType()) << std::endl;
656  fQAHistos->FillTH2("hClusterConstituentEDefault", jet->Pt(), clust->E());
657  fQAHistos->FillTH2("hClusterConstituentENLC", jet->Pt(), clust->GetCluster()->GetNonLinCorrEnergy());
658  fQAHistos->FillTH2("hClusterConstituentEHC", jet->Pt(), clust->GetCluster()->GetHadCorrEnergy());
659  }
660 #endif
661 }
662 
664  int ncharged = 0, nneutral = jet.GetNumberOfClusters();
665  if(particles) {
666  for(decltype(jet.GetNumberOfTracks()) ipart = 0; ipart < jet.GetNumberOfTracks(); ipart++){
667  auto part = jet.TrackAt(ipart, particles->GetArray());
668  if(!part) continue;
669  if(part->Charge()) ncharged++;
670  else nneutral++;
671  }
672  }
673  // check if the jet has at least one consituent for jet substructure
674  int nallowed = 0;
675  nallowed += fUseChargedConstituents ? ncharged : 0;
676  nallowed += fUseNeutralConstituents ? nneutral : 0;
677  return nallowed > 0;
678 }
679 
680 std::vector<Triggerinfo> AliAnalysisTaskEmcalJetSubstructureTree::DecodeTriggerString(const std::string &triggerstring) const {
681  std::vector<Triggerinfo> result;
682  std::stringstream triggerparser(triggerstring);
683  std::string currenttrigger;
684  while(std::getline(triggerparser, currenttrigger, ' ')){
685  if(!currenttrigger.length()) continue;
686  std::vector<std::string> tokens;
687  std::stringstream triggerdecoder(currenttrigger);
688  std::string token;
689  while(std::getline(triggerdecoder, token, '-')) tokens.emplace_back(token);
690  result.emplace_back(Triggerinfo({tokens[0], tokens[1], tokens[2], tokens[3]}));
691  }
692  return result;
693 }
694 
695 std::string AliAnalysisTaskEmcalJetSubstructureTree::MatchTrigger(const std::string &triggertoken) const {
696  std::vector<std::string> tokens;
697  std::string result;
698  std::stringstream decoder(fInputEvent->GetFiredTriggerClasses().Data());
699  while(std::getline(decoder, result, ' ')) tokens.emplace_back(result);
700  result.clear();
701  for(auto t : tokens) {
702  if(t.find(triggertoken) != std::string::npos) {
703  // take first occurrence - downscale factor should normally be the same
704  result = t;
705  break;
706  }
707  }
708  return result;
709 }
710 
711 bool AliAnalysisTaskEmcalJetSubstructureTree::IsSelectEmcalTriggers(const std::string &triggerstring) const {
712  const std::array<std::string, 8> kEMCALTriggers = {
713  "EJ1", "EJ2", "DJ1", "DJ2", "EG1", "EG2", "DG1", "DG2"
714  };
715  bool isEMCAL = false;
716  for(auto emcaltrg : kEMCALTriggers) {
717  if(triggerstring.find(emcaltrg) != std::string::npos) {
718  isEMCAL = true;
719  break;
720  }
721  }
722  return isEMCAL;
723 }
724 
726  AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
727 
728  Bool_t isAOD(kFALSE);
729  AliInputEventHandler *inputhandler = static_cast<AliInputEventHandler *>(mgr->GetInputEventHandler());
730  if(inputhandler) {
731  if(inputhandler->IsA() == AliAODInputHandler::Class()){
732  std::cout << "Analysing AOD events\n";
733  isAOD = kTRUE;
734  } else {
735  std::cout << "Analysing ESD events\n";
736  }
737  }
738 
739  std::stringstream taskname;
740  taskname << "JetSubstructureTreemaker_R" << std::setw(2) << std::setfill('0') << int(jetradius*10) << trigger;
742  mgr->AddTask(treemaker);
743  treemaker->SetMakeGeneralHistograms(kTRUE);
744 
745  // Adding containers
746  if(isMC) {
747  AliParticleContainer *particles = treemaker->AddMCParticleContainer("mcparticles");
748  particles->SetMinPt(0.);
749 
750  AliJetContainer *mcjets = treemaker->AddJetContainer(
751  jettype,
753  recombinationScheme,
754  jetradius,
756  particles, nullptr);
757  mcjets->SetName("mcjets");
758  mcjets->SetJetPtCut(20.);
759  }
760 
761  if(isData) {
762  AliTrackContainer *tracks(nullptr);
763  if((jettype == AliJetContainer::kChargedJet) || (jettype == AliJetContainer::kFullJet)){
765  std::cout << "Track container name: " << tracks->GetName() << std::endl;
766  tracks->SetMinPt(0.15);
767  }
768  AliClusterContainer *clusters(nullptr);
770  std::cout << "Using full or neutral jets ..." << std::endl;
772  std::cout << "Cluster container name: " << clusters->GetName() << std::endl;
773  clusters->SetClusHadCorrEnergyCut(0.3); // 300 MeV E-cut
774  clusters->SetDefaultClusterEnergy(AliVCluster::kHadCorr);
775  } else {
776  std::cout << "Using charged jets ... " << std::endl;
777  }
778 
779  AliJetContainer *datajets = treemaker->AddJetContainer(
780  jettype,
782  recombinationScheme,
783  jetradius,
785  tracks, clusters);
786  datajets->SetName("datajets");
787  datajets->SetJetPtCut(20.);
788 
789  treemaker->SetUseAliAnaUtils(true, true);
790  treemaker->SetVzRange(-10., 10);
791 
792  // configure trigger selection
793  std::string triggerstring(trigger);
794  if(triggerstring.find("INT7") != std::string::npos) {
795  treemaker->SetTriggerBits(AliVEvent::kINT7);
796  } else if(triggerstring.find("EJ1") != std::string::npos) {
797  treemaker->SetTriggerBits(AliVEvent::kEMCEJE);
798  treemaker->SetTriggerString("EJ1");
799  } else if(triggerstring.find("EJ2") != std::string::npos) {
800  treemaker->SetTriggerBits(AliVEvent::kEMCEJE);
801  treemaker->SetTriggerString("EJ2");
802  }
803  }
804 
805  std::string jettypestring;
806  switch(jettype) {
807  case AliJetContainer::kFullJet: jettypestring = "FullJets"; break;
808  case AliJetContainer::kChargedJet: jettypestring = "ChargedJets"; break;
809  case AliJetContainer::kNeutralJet: jettypestring = "NeutralJets"; break;
810  default: jettypestring = "Undef";
811  };
812 
813  // Connecting containers
814  std::stringstream outputfile, histname, treename;
815  outputfile << mgr->GetCommonFileName() << ":JetSubstructure_" << jettypestring << "_R" << std::setw(2) << std::setfill('0') << int(jetradius * 10.) << "_" << trigger;
816  histname << "JetSubstructureHistos_" << jettypestring << "_R" << std::setw(2) << std::setfill('0') << int(jetradius * 10.) << "_" << trigger;
817  treename << "JetSubstructureTree_" << jettypestring << "_R" << std::setw(2) << std::setfill('0') << int(jetradius * 10.) << "_" << trigger;
818  mgr->ConnectInput(treemaker, 0, mgr->GetCommonInputContainer());
819  mgr->ConnectOutput(treemaker, 1, mgr->CreateContainer(histname.str().data(), AliEmcalList::Class(), AliAnalysisManager::kOutputContainer, outputfile.str().data()));
820  mgr->ConnectOutput(treemaker, 2, mgr->CreateContainer(treename.str().data(), TTree::Class(), AliAnalysisManager::kOutputContainer, mgr->GetCommonFileName()));
821 
822  return treemaker;
823 }
824 
825 std::string Triggerinfo::ExpandClassName() const {
826  std::string result = fTriggerClass + "-" + fBunchCrossing + "-" + fPastFutureProtection + "-" + fTriggerCluster;
827  return result;
828 }
829 
830 bool Triggerinfo::IsTriggerClass(const std::string &triggerclass) const {
831  return fTriggerClass.substr(1) == triggerclass; // remove C from trigger class part
832 }
833 
834 void AliSoftDropParameters::LinkJetTreeBranches(TTree *jettree, const char *tag) {
835  LinkBranch(jettree, &fZg, Form("Zg%s", tag), "D");
836  LinkBranch(jettree, &fRg, Form("Rg%s", tag), "D");
837  LinkBranch(jettree, &fMg, Form("Mg%s", tag), "D");
838  LinkBranch(jettree, &fPtg, Form("Ptg%s", tag), "D");
839  LinkBranch(jettree, &fMug, Form("Mug%s", tag), "D");
840  LinkBranch(jettree, &fDeltaR, Form("DeltaRg%s", tag), "D");
841  LinkBranch(jettree, &fNDropped, Form("NDropped%s", tag), "I");
842 };
843 
845  LinkBranch(jettree, &fOneSubjettiness, Form("OneSubjettiness%s", tag), "D");
846  LinkBranch(jettree, &fTwoSubjettiness, Form("TwoSubjettiness%s", tag), "D");
847 }
848 
850  LinkBranch(jettree, &fAngularity, Form("Angularity%s", tag), "D");
851  LinkBranch(jettree, &fPtD, Form("PtD%s", tag), "D");
852 }
853 
854 void AliJetKineParameters::LinkJetTreeBranches(TTree *jettree, const char *tag){
855  LinkBranch(jettree, &fPt, Form("PtJet%s", tag), "D");
856  LinkBranch(jettree, &fE, Form("EJet%s", tag), "D");
857  LinkBranch(jettree, &fEta, Form("Eta%s", tag), "D");
858  LinkBranch(jettree, &fPhi, Form("Phi%s", tag), "D");
859  LinkBranch(jettree, &fArea, Form("Area%s", tag), "D");
860  LinkBranch(jettree, &fMass, Form("Mass%s", tag), "D");
861  LinkBranch(jettree, &fNEF, Form("NEFt%s", tag), "D");
862  LinkBranch(jettree, &fPt, Form("PtJet%s", tag), "D");
863  LinkBranch(jettree, &fNCharged, Form("NCharged%s", tag), "I");
864  LinkBranch(jettree, &fNNeutral, Form("NNeutral%s", tag), "I");
865 }
866 
868  LinkBranch(jettree, &fJetRadius, "Radius", "D");
869  LinkBranch(jettree, &fEventWeight, "EventWeight", "D");
870  if(fillRho) {
871  std::string varnames[] = {"RhoPtRec", "RhoPtSim", "RhoMassRec", "RhoMassSim"};
872  for(int i = 0; i < 4; i++){
873  LinkBranch(jettree, fRhoParamters + i, varnames[i].data(), "D");
874  }
875  }
876 }
877 
878 void LinkBranch(TTree *jettree, void *data, const char *branchname, const char *type) {
879  jettree->Branch(branchname, data, Form("%s/%s", branchname, type));
880 }
881 } /* 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
AliNSubjettinessParameters * fNSubTrue
! Data field for true n-subjettiness parameters in jet tree
virtual double E() const
Access to constituent energy.
Double_t Eta() const
Definition: AliEmcalJet.h:121
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
static AliEmcalDownscaleFactorsOCDB * Instance()
AliSoftDropParameters * fSoftDropMeasured
! Data field for measured soft drop parameters in jet tree
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.
Double_t GetDownscaleFactorForTriggerClass(const TString &trigger) const
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.
void SetVzRange(Double_t min, Double_t max)
AliNSubjettinessParameters * fNSubMeasured
! Data field for measured n-subjettiness parameters in jet tree
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.
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
AliJetKineParameters MakeJetKineParameters(const AliEmcalJet &jet) const
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
void LinkBranch(TTree *jettree, void *data, const char *branchname, const char *type)
Helper function linking struct members to branches in the jet substructure tree.
AliJetStructureParameters * fJetStructureTrue
! True jet substructure paramteres
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
AliJetStructureParameters * fJetStructureMeasured
! Measured jet substructure parameters
Double_t Pt() const
Definition: AliEmcalJet.h:109
Bool_t isMC
AliJetKineParameters * fKineRec
! Detector level jet kinematics
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
AliJetKineParameters * fKineSim
! Particle level jet kinematics
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
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
AliJetTreeGlobalParameters * fGlobalTreeParams
! Global jet tree parameters (same for all jets in event)
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:148
void SetDefaultClusterEnergy(Int_t d)
bool SelectJet(const AliEmcalJet &jet, const AliParticleContainer *particles) const
AliSoftDropParameters * fSoftDropTrue
! Data field for true soft drop parameters in jet tree
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.