PandoraAnalysis  01.02.01
PfoAnalysis.h
1 
9 #ifndef PFO_ANALYSIS_H
10 #define PFO_ANALYSIS_H 1
11 
12 #include "EVENT/LCStrVec.h"
13 
14 #include "marlin/Processor.h"
15 
16 #include "CalibrationHelper.h"
17 
18 #include <set>
19 #include <vector>
20 
21 class TFile;
22 class TH1F;
23 class TTree;
24 
25 //------------------------------------------------------------------------------------------------------------------------------------------
26 
30 class PfoAnalysis : public marlin::Processor
31 {
32 public:
36  PfoAnalysis();
37 
41  virtual Processor *newProcessor();
42 
46  virtual void init();
47 
53  virtual void processRunHeader(lcio::LCRunHeader *pLCRunHeader);
54 
60  virtual void processEvent(EVENT::LCEvent *pLCEvent);
61 
67  virtual void check(EVENT::LCEvent *pLCEvent);
68 
72  virtual void end();
73 
74 private:
75  typedef std::vector<const EVENT::ReconstructedParticle *> ParticleVector;
76  typedef std::vector<const EVENT::MCParticle*> MCParticleVector;
77  typedef std::set<const EVENT::MCParticle*> MCParticleList;
78  typedef std::vector<std::string> StringVector;
79 
83  void Clear();
84 
90  void ExtractCollections(EVENT::LCEvent *pLCEvent);
91 
98  void ApplyPfoSelectionRules(const EVENT::MCParticle *pMCParticle, MCParticleList &mcPfoCandidates) const;
99 
105  void MakeQuarkVariables(EVENT::LCEvent *pLCEvent);
106 
110  void PerformPfoAnalysis();
111 
118  static bool SortPfoTargetsByEnergy(const EVENT::MCParticle *const pLhs, const EVENT::MCParticle *const pRhs);
119 
120  int m_nRun;
121  int m_nEvt;
122 
123  int m_nRunSum;
124  int m_nEvtSum;
125 
126  std::string m_inputPfoCollection;
127  std::string m_mcParticleCollection;
128 
129  int m_printing;
130  std::string m_rootFile;
131 
132  int m_lookForQuarksWithMotherZ;
133 
134  float m_mcPfoSelectionRadius;
135  float m_mcPfoSelectionMomentum;
136  float m_mcPfoSelectionLowEnergyNPCutOff;
137 
138  ParticleVector m_pfoVector;
139  MCParticleVector m_pfoTargetVector;
140 
141  int m_nPfosTotal;
142  int m_nPfosNeutralHadrons;
143  int m_nPfosPhotons;
144  int m_nPfosTracks;
145  float m_pfoEnergyTotal;
146  float m_pfoEnergyNeutralHadrons;
147  float m_pfoEnergyPhotons;
148 
149  float m_pfoEnergyTracks;
150  float m_pfoECalToEmEnergy;
151  float m_pfoECalToHadEnergy;
152  float m_pfoHCalToEmEnergy;
153  float m_pfoHCalToHadEnergy;
154  float m_pfoMuonToEnergy;
155  float m_pfoOtherEnergy;
156 
157  float m_pfoMassTotal;
158 
159  typedef std::vector<float> FloatVector;
160  FloatVector m_pfoEnergies;
161  FloatVector m_pfoPx;
162  FloatVector m_pfoPy;
163  FloatVector m_pfoPz;
164  FloatVector m_pfoCosTheta;
165 
166  FloatVector m_pfoTargetEnergies;
167  FloatVector m_pfoTargetPx;
168  FloatVector m_pfoTargetPy;
169  FloatVector m_pfoTargetPz;
170  FloatVector m_pfoTargetCosTheta;
171 
172  typedef std::vector<int> IntVector;
173  IntVector m_pfoPdgCodes;
174  IntVector m_pfoTargetPdgCodes;
175 
176  int m_nPfoTargetsTotal;
177  int m_nPfoTargetsNeutralHadrons;
178  int m_nPfoTargetsPhotons;
179  int m_nPfoTargetsTracks;
180 
181  float m_pfoTargetsEnergyTotal;
182  float m_pfoTargetsEnergyNeutralHadrons;
183  float m_pfoTargetsEnergyPhotons;
184  float m_pfoTargetsEnergyTracks;
185 
186  float m_mcEnergyENu;
187  float m_mcEnergyFwd;
188  float m_eQQ;
189  float m_eQ1;
190  float m_eQ2;
191  float m_costQQ;
192  float m_costQ1;
193  float m_costQ2;
194  float m_mQQ;
195  float m_thrust;
196  int m_qPdg;
197 
198  TFile *m_pTFile;
199  TTree *m_pTTree;
200  TH1F *m_hPfoEnergySum;
201  TH1F *m_hPfoEnergySumL7A;
202 
203  int m_collectCalibrationDetails;
204  pandora_analysis::CalibrationHelper *m_pCalibrationHelper;
205  pandora_analysis::CalibrationHelper::Settings m_calibrationHelperSettings;
206 };
207 
208 //------------------------------------------------------------------------------------------------------------------------------------------
209 
210 inline marlin::Processor *PfoAnalysis::newProcessor()
211 {
212  return new PfoAnalysis;
213 }
214 
215 #endif // #ifndef PFO_ANALYSIS_H
virtual void end()
End, called at shutdown.
Definition: PfoAnalysis.cc:320
PfoAnalysis()
Default constructor.
Definition: PfoAnalysis.cc:29
CalibrationHelper class.
Definition: CalibrationHelper.h:29
virtual void processEvent(EVENT::LCEvent *pLCEvent)
Process event, main entry point.
Definition: PfoAnalysis.cc:292
virtual void processRunHeader(lcio::LCRunHeader *pLCRunHeader)
Process run header.
Definition: PfoAnalysis.cc:283
Settings class.
Definition: CalibrationHelper.h:35
virtual void check(EVENT::LCEvent *pLCEvent)
Checks for event.
Definition: PfoAnalysis.cc:314
virtual void init()
Initialize, called at startup.
Definition: PfoAnalysis.cc:204
virtual Processor * newProcessor()
Create new processor.
Definition: PfoAnalysis.h:210
PfoAnalysis class.
Definition: PfoAnalysis.h:30