AliPhysics  fe039ad (fe039ad)
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
AliAnalysisTaskEmcalJetSubstructureTree.h
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 #ifndef ALIANALYSISTASKEMCALJETSUBSTRUCTURETREE_H
28 #define ALIANALYSISTASKEMCALJETSUBSTRUCTURETREE_H
29 
31 #include <exception>
32 #include <TString.h>
33 #include <fastjet/PseudoJet.hh>
34 #include <fastjet/JetDefinition.hh>
35 
36 class THistManager;
37 class TTree;
39 class AliEmcalJet;
41 class AliTrackContainer;
42 
43 namespace EmcalTriggerJets {
44 
53 };
54 
63 };
64 
77 };
78 
87  fastjet::JetAlgorithm fRecluserAlgo;
88 };
89 
93 };
94 
98 };
99 
107 public:
108  class ReclusterizerException : public std::exception {
109  public:
110  ReclusterizerException() : std::exception() {}
111  virtual ~ReclusterizerException() throw() {}
112 
113  virtual const char *what() const throw() { return "Error in reclusterizing in fastjet"; }
114  };
115  class SubstructureException : public std::exception {
116  public:
117  SubstructureException() : std::exception() {}
118  virtual ~SubstructureException() throw() {}
119 
120  virtual const char *what() const throw() { return "Error in builing substructure observable"; }
121  };
122  class SoftDropException : public std::exception {
123  public:
124  SoftDropException() : std::exception() {}
125  virtual ~SoftDropException() throw() {}
126 
127  virtual const char *what() const throw() { return "No associated softdrop structure found for jet - softdrop algorithm failing"; }
128  };
130  kCAAlgo = 0,
131  kKTAlgo = 1,
133  };
135  kTRadius = 0,
136  kTWeight = 1,
141  kTEtaRec = 6,
142  kTEtaSim = 7,
143  kTPhiRec = 8,
144  kTPhiSim = 9,
149  kTAreaRec = 14,
150  kTAreaSim = 15,
151  kTNEFRec = 16,
152  kTNEFSim = 17,
153  kTMassRec = 18,
154  kTMassSim = 19,
156  kTZgTrue = 21,
158  kTRgTrue = 23,
160  kTMgTrue = 25,
162  kTPtgTrue = 27,
164  kTMugTrue = 29,
172  kTPtDTrue = 37,
178  kTNVar = 43
179  };
180 
182  AliAnalysisTaskEmcalJetSubstructureTree(const char *name);
184 
185  void SetTriggerBits(UInt_t triggersel) { fTriggerSelectionBits = triggersel; }
186  void SetTriggerString(TString triggerstring) { fTriggerSelectionString = triggerstring; }
187  void SetUseDownscaleWeight(Bool_t usedownscale) { fUseDownscaleWeight = usedownscale; }
188 
189  void SetSoftdropDefiniion(Double_t zcut, Double_t betacut, Reclusterizer_t reclusterizer) {
190  fSDZCut = zcut;
191  fSDBetaCut = betacut;
192  fReclusterizer = reclusterizer;
193  }
194 
195  void SetFillPartLevelBranches(Bool_t doFill) { fFillPart = doFill; }
196  void SetFillAcceptance(Bool_t doFill) { fFillAcceptance = doFill; }
197  void SetFillRhoBranches(Bool_t doFill) { fFillRho = doFill; }
198  void SetFillMassBranches(Bool_t doFill) { fFillMass = doFill; }
199  void SetFillSoftdropBranches(Bool_t doFill) { fFillSoftDrop = doFill; }
200  void SetFillNSubjettinessBranches(Bool_t doFill) { fFillNSub = doFill; }
202 
204 
205 protected:
206  virtual void UserCreateOutputObjects();
207  virtual bool Run();
208  virtual void RunChanged(Int_t newrun);
209 
210  AliJetSubstructureData MakeJetSubstructure(const AliEmcalJet &jet, double jetradius, const AliParticleContainer *tracks, const AliClusterContainer *clusters, const AliJetSubstructureSettings &settings) const;
211 
212  AliSoftDropParameters MakeSoftDropParameters(const fastjet::PseudoJet &jet, const AliSoftdropDefinition &cut) const;
213 
214  AliNSubjettinessParameters MakeNsubjettinessParameters(const fastjet::PseudoJet &jet, const AliNSubjettinessDefinition &cut) const;
215 
216  Double_t MakeAngularity(const AliEmcalJet &jet, const AliParticleContainer *tracks, const AliClusterContainer *clusters) const;
217 
218  Double_t MakePtD(const AliEmcalJet &jet, const AliParticleContainer *const particles, const AliClusterContainer *const clusters) const;
219 
220  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);
221 
222  void DoConstituentQA(const AliEmcalJet *jet, const AliParticleContainer *tracks, const AliClusterContainer *clusters);
223 
224  void LinkOutputBranch(const TString &branchname, Double_t *datalocation);
225 
226  bool IsPartBranch(const TString &branchname) const;
227  bool IsAcceptanceBranch(const TString &branchname) const;
228  bool IsRhoBranch(const TString &branchname) const;
229  bool IsMassBranch(const TString &branchname) const;
230  bool IsSoftdropBranch(const TString &branchname) const;
231  bool IsNSubjettinessBranch(const TString &branchname) const;
232  bool IsStructbranch(const TString &branchname) const;
233 
234 private:
238 
242 
246 
247  // Fill levels for tree (save disk space when information is not needed)
255 
259 };
260 
261 } /* namespace EmcalTriggerJets */
262 
263 #endif /* ALIANALYSISTASKEMCALJETSUBSTRUCTURETREE_H */
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)
Double_t MakePtD(const AliEmcalJet &jet, const AliParticleContainer *const particles, const AliClusterContainer *const clusters) const
void SetSoftdropDefiniion(Double_t zcut, Double_t betacut, Reclusterizer_t reclusterizer)
Container with name, TClonesArray and cuts for particles.
AliJetSubstructureData MakeJetSubstructure(const AliEmcalJet &jet, double jetradius, const AliParticleContainer *tracks, const AliClusterContainer *clusters, const AliJetSubstructureSettings &settings) const
Structure for results from the soft drop algorithm.
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.
fastjet::JetAlgorithm fRecluserAlgo
Reclusterization algorithm.
Container for particles within the EMCAL framework.
int Int_t
Definition: External.C:63
unsigned int UInt_t
Definition: External.C:33
AliNSubjettinessParameters MakeNsubjettinessParameters(const fastjet::PseudoJet &jet, const AliNSubjettinessDefinition &cut) const
Definition for the algorithm obtaining the softdrop parameters.
virtual void RunChanged(Int_t newrun)
Process tasks relevant when a file with a different run number is processed.
Double_t MakeAngularity(const AliEmcalJet &jet, const AliParticleContainer *tracks, const AliClusterContainer *clusters) const
Bool_t isMC
void DoConstituentQA(const AliEmcalJet *jet, const AliParticleContainer *tracks, const AliClusterContainer *clusters)
void LinkOutputBranch(const TString &branchname, Double_t *datalocation)
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
AliSoftDropParameters MakeSoftDropParameters(const fastjet::PseudoJet &jet, const AliSoftdropDefinition &cut) const
bool Bool_t
Definition: External.C:53
Container structure for EMCAL clusters.