Fork me on GitHub

source: git/classes/DelphesClasses.h @ 0e2f49b

ImprovedOutputFileTimingdual_readoutllp
Last change on this file since 0e2f49b was 0e2f49b, checked in by Michele Selvaggi <michele.selvaggi@…>, 4 years ago

added classes for Vertex Finding

  • Property mode set to 100644
File size: 17.3 KB
Line 
1/*
2 *  Delphes: a framework for fast simulation of a generic collider experiment
3 *  Copyright (C) 2012-2014  Universite catholique de Louvain (UCL), Belgium
4 *
5 *  This program is free software: you can redistribute it and/or modify
6 *  it under the terms of the GNU General Public License as published by
7 *  the Free Software Foundation, either version 3 of the License, or
8 *  (at your option) any later version.
9 *
10 *  This program is distributed in the hope that it will be useful,
11 *  but WITHOUT ANY WARRANTY; without even the implied warranty of
12 *  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
13 *  GNU General Public License for more details.
14 *
15 *  You should have received a copy of the GNU General Public License
16 *  along with this program.  If not, see <http://www.gnu.org/licenses/>.
17 */
18
19#ifndef DelphesClasses_h
20#define DelphesClasses_h
21
22/**
23 *
24 *  Definition of classes to be stored in the root tree.
25 *  Function CompareXYZ sorts objects by the variable XYZ that MUST be
26 *  present in the data members of the root tree class of the branch.
27 *
28 *  \author P. Demin - UCL, Louvain-la-Neuve
29 *
30 */
31
32// Dependencies (#includes)
33
34#include "TRef.h"
35#include "TObject.h"
36#include "TRefArray.h"
37#include "TLorentzVector.h"
38
39#include "classes/SortableObject.h"
40
41class DelphesFactory;
42
43//---------------------------------------------------------------------------
44
45class Event: public TObject
46{
47public:
48
49  Long64_t Number; // event number
50
51  Float_t ReadTime;
52  Float_t ProcTime;
53
54  ClassDef(Event, 1)
55};
56
57//---------------------------------------------------------------------------
58
59class LHCOEvent: public Event
60{
61public:
62
63  Int_t Trigger; // trigger word
64
65  ClassDef(LHCOEvent, 1)
66};
67
68//---------------------------------------------------------------------------
69
70class LHEFEvent: public Event
71{
72public:
73
74  Int_t ProcessID; // subprocess code for the event | hepup.IDPRUP
75
76  Float_t Weight; // weight for the event | hepup.XWGTUP
77  Float_t ScalePDF; // scale in GeV used in the calculation of the PDFs in the event | hepup.SCALUP
78  Float_t AlphaQED; // value of the QED coupling used in the event | hepup.AQEDUP
79  Float_t AlphaQCD; // value of the QCD coupling used in the event | hepup.AQCDUP
80
81  ClassDef(LHEFEvent, 2)
82};
83
84//---------------------------------------------------------------------------
85
86class LHEFWeight: public TObject
87{
88public:
89  Int_t ID; // weight ID
90  Float_t Weight; // weight value
91
92  ClassDef(LHEFWeight, 1)
93};
94
95//---------------------------------------------------------------------------
96
97class HepMCEvent: public Event
98{
99public:
100
101  Int_t ProcessID; // unique signal process id | signal_process_id()
102  Int_t MPI; // number of multi parton interactions | mpi ()
103
104  Float_t Weight; // weight for the event
105
106  Float_t Scale; // energy scale, see hep-ph/0109068 | event_scale()
107  Float_t AlphaQED; // QED coupling, see hep-ph/0109068 | alphaQED()
108  Float_t AlphaQCD; // QCD coupling, see hep-ph/0109068 | alphaQCD()
109
110  Int_t ID1; // flavour code of first parton | pdf_info()->id1()
111  Int_t ID2; // flavour code of second parton | pdf_info()->id2()
112
113  Float_t X1; // fraction of beam momentum carried by first parton ("beam side") | pdf_info()->x1()
114  Float_t X2; // fraction of beam momentum carried by second parton ("target side") | pdf_info()->x2()
115
116  Float_t ScalePDF; // Q-scale used in evaluation of PDF's (in GeV) | pdf_info()->scalePDF()
117
118  Float_t PDF1; // PDF (id1, x1, Q) | pdf_info()->pdf1()
119  Float_t PDF2; // PDF (id2, x2, Q) | pdf_info()->pdf2()
120
121  ClassDef(HepMCEvent, 2)
122};
123
124//---------------------------------------------------------------------------
125
126class GenParticle: public SortableObject
127{
128public:
129  Int_t PID; // particle HEP ID number | hepevt.idhep[number]
130
131  Int_t Status; // particle status | hepevt.isthep[number]
132  Int_t IsPU; // 0 or 1 for particles from pile-up interactions
133
134  Int_t M1; // particle 1st mother | hepevt.jmohep[number][0] - 1
135  Int_t M2; // particle 2nd mother | hepevt.jmohep[number][1] - 1
136
137  Int_t D1; // particle 1st daughter | hepevt.jdahep[number][0] - 1
138  Int_t D2; // particle last daughter | hepevt.jdahep[number][1] - 1
139
140  Int_t Charge; // particle charge
141
142  Float_t Mass; // particle mass
143
144  Float_t E; // particle energy | hepevt.phep[number][3]
145  Float_t Px; // particle momentum vector (x component) | hepevt.phep[number][0]
146  Float_t Py; // particle momentum vector (y component) | hepevt.phep[number][1]
147  Float_t Pz; // particle momentum vector (z component) | hepevt.phep[number][2]
148
149  Float_t D0;
150  Float_t DZ;
151  Float_t P;
152  Float_t PT;
153  Float_t CtgTheta;
154  Float_t Phi;
155  Float_t Eta; // particle pseudorapidity
156  Float_t Rapidity; // particle rapidity
157
158  Float_t T; // particle vertex position (t component) | hepevt.vhep[number][3]
159  Float_t X; // particle vertex position (x component) | hepevt.vhep[number][0]
160  Float_t Y; // particle vertex position (y component) | hepevt.vhep[number][1]
161  Float_t Z; // particle vertex position (z component) | hepevt.vhep[number][2]
162
163  static CompBase *fgCompare; //!
164  const CompBase *GetCompare() const { return fgCompare; }
165
166  TLorentzVector P4() const;
167
168  ClassDef(GenParticle, 1)
169};
170
171//---------------------------------------------------------------------------
172
173class Vertex: public TObject
174{
175public:
176  Float_t T; // vertex position (t component)
177  Float_t X; // vertex position (x component)
178  Float_t Y; // vertex position (y component)
179  Float_t Z; // vertex position (z component)
180 
181  Double_t ErrorX;
182  Double_t ErrorY;
183  Double_t ErrorZ;
184
185  Int_t Index;
186  Int_t NDF;
187  Double_t Sigma;
188  Double_t SumPT2;
189  Double_t BTVSumPT2;
190  Double_t GenDeltaZ;
191  Double_t GenSumPT2;
192
193  ClassDef(Vertex, 2)
194};
195
196//---------------------------------------------------------------------------
197
198class MissingET: public TObject
199{
200public:
201  Float_t MET; // mising transverse energy
202  Float_t Eta; // mising energy pseudorapidity
203  Float_t Phi; // mising energy azimuthal angle
204
205  TLorentzVector P4() const;
206
207  ClassDef(MissingET, 1)
208};
209
210//---------------------------------------------------------------------------
211
212class ScalarHT: public TObject
213{
214public:
215  Float_t HT; // scalar sum of transverse momenta
216
217  ClassDef(ScalarHT, 1)
218};
219
220//---------------------------------------------------------------------------
221
222class Rho: public TObject
223{
224public:
225  Float_t Rho; // rho energy density
226  Float_t Edges[2]; // pseudorapidity range edges
227
228  ClassDef(Rho, 1)
229};
230
231//---------------------------------------------------------------------------
232
233class Weight: public TObject
234{
235public:
236  Float_t Weight; // weight for the event
237
238  ClassDef(Weight, 1)
239};
240
241//---------------------------------------------------------------------------
242
243class Photon: public SortableObject
244{
245public:
246
247  Float_t PT; // photon transverse momentum
248  Float_t Eta; // photon pseudorapidity
249  Float_t Phi; // photon azimuthal angle
250
251  Float_t E; // photon energy
252
253  Float_t T; //particle arrival time of flight
254
255  Float_t EhadOverEem; // ratio of the hadronic versus electromagnetic energy deposited in the calorimeter
256
257  TRefArray Particles; // references to generated particles
258
259  // Isolation variables
260
261  Float_t IsolationVar;
262  Float_t IsolationVarRhoCorr;
263  Float_t SumPtCharged;
264  Float_t SumPtNeutral;
265  Float_t SumPtChargedPU;
266  Float_t SumPt;
267
268  static CompBase *fgCompare; //!
269  const CompBase *GetCompare() const { return fgCompare; }
270
271  TLorentzVector P4() const;
272
273  ClassDef(Photon, 3)
274};
275
276//---------------------------------------------------------------------------
277
278class Electron: public SortableObject
279{
280public:
281
282  Float_t PT; // electron transverse momentum
283  Float_t Eta; // electron pseudorapidity
284  Float_t Phi; // electron azimuthal angle
285
286  Float_t T; //particle arrival time of flight
287
288  Int_t Charge; // electron charge
289
290  Float_t EhadOverEem; // ratio of the hadronic versus electromagnetic energy deposited in the calorimeter
291
292  TRef Particle; // reference to generated particle
293
294  // Isolation variables
295
296  Float_t IsolationVar;
297  Float_t IsolationVarRhoCorr;
298  Float_t SumPtCharged;
299  Float_t SumPtNeutral;
300  Float_t SumPtChargedPU;
301  Float_t SumPt;
302
303  static CompBase *fgCompare; //!
304  const CompBase *GetCompare() const { return fgCompare; }
305
306  TLorentzVector P4() const;
307
308  ClassDef(Electron, 3)
309};
310
311//---------------------------------------------------------------------------
312
313class Muon: public SortableObject
314{
315public:
316
317  Float_t PT; // muon transverse momentum
318  Float_t Eta; // muon pseudorapidity
319  Float_t Phi; // muon azimuthal angle
320
321  Float_t T; //particle arrival time of flight
322
323  Int_t Charge; // muon charge
324
325  TRef Particle; // reference to generated particle
326
327   // Isolation variables
328
329  Float_t IsolationVar;
330  Float_t IsolationVarRhoCorr;
331  Float_t SumPtCharged;
332  Float_t SumPtNeutral;
333  Float_t SumPtChargedPU;
334  Float_t SumPt;
335
336  static CompBase *fgCompare; //!
337  const CompBase *GetCompare() const { return fgCompare; }
338
339  TLorentzVector P4() const;
340
341  ClassDef(Muon, 3)
342};
343
344//---------------------------------------------------------------------------
345
346class Jet: public SortableObject
347{
348public:
349
350  Float_t PT; // jet transverse momentum
351  Float_t Eta; // jet pseudorapidity
352  Float_t Phi; // jet azimuthal angle
353
354  Float_t T; //particle arrival time of flight
355
356  Float_t Mass; // jet invariant mass
357
358  Float_t DeltaEta;  // jet radius in pseudorapidity
359  Float_t DeltaPhi;  // jet radius in azimuthal angle
360
361  UInt_t Flavor;
362  UInt_t FlavorAlgo;
363  UInt_t FlavorPhys;
364
365  UInt_t BTag; // 0 or 1 for a jet that has been tagged as containing a heavy quark
366  UInt_t BTagAlgo;
367  UInt_t BTagPhys;
368
369  UInt_t TauTag; // 0 or 1 for a jet that has been tagged as a tau
370
371  Int_t Charge; // tau charge
372
373  Float_t EhadOverEem; // ratio of the hadronic versus electromagnetic energy deposited in the calorimeter
374
375  Int_t NCharged; // number of charged constituents
376  Int_t NNeutrals; // number of neutral constituents
377  Float_t Beta; // (sum pt of charged pile-up constituents)/(sum pt of charged constituents)
378  Float_t BetaStar; // (sum pt of charged constituents coming from hard interaction)/(sum pt of charged constituents)
379  Float_t MeanSqDeltaR; // average distance (squared) between constituent and jet weighted by pt (squared) of constituent
380  Float_t PTD; // average pt between constituent and jet weighted by pt of constituent
381  Float_t FracPt[5]; // (sum pt of constituents within a ring 0.1*i < DeltaR < 0.1*(i+1))/(sum pt of constituents)
382
383  Float_t Tau[5]; // N-subjettiness
384
385  TLorentzVector TrimmedP4[5]; // first entry (i = 0) is the total Trimmed Jet 4-momenta and from i = 1 to 4 are the trimmed subjets 4-momenta
386  TLorentzVector PrunedP4[5]; // first entry (i = 0) is the total Pruned Jet 4-momenta and from i = 1 to 4 are the pruned subjets 4-momenta
387  TLorentzVector SoftDroppedP4[5]; // first entry (i = 0) is the total SoftDropped Jet 4-momenta and from i = 1 to 4 are the pruned subjets 4-momenta
388
389  Int_t NSubJetsTrimmed; // number of subjets trimmed
390  Int_t NSubJetsPruned; // number of subjets pruned
391  Int_t NSubJetsSoftDropped; // number of subjets soft-dropped
392
393  TRefArray Constituents; // references to constituents
394  TRefArray Particles; // references to generated particles
395
396  static CompBase *fgCompare; //!
397  const CompBase *GetCompare() const { return fgCompare; }
398
399  TLorentzVector P4() const;
400  TLorentzVector Area;
401
402  ClassDef(Jet, 3)
403};
404
405//---------------------------------------------------------------------------
406
407class Track: public SortableObject
408{
409public:
410  Int_t PID; // HEP ID number
411
412  Int_t Charge; // track charge
413
414  Float_t Eta; // track pseudorapidity
415 
416  Float_t EtaOuter; // track pseudorapidity at the tracker edge
417  Float_t PhiOuter; // track azimuthal angle at the tracker edge
418
419  Float_t X; // track vertex position (x component)
420  Float_t Y; // track vertex position (y component)
421  Float_t Z; // track vertex position (z component)
422  Float_t T; // track vertex position (z component)
423
424  Float_t XOuter; // track position (x component) at the tracker edge
425  Float_t YOuter; // track position (y component) at the tracker edge
426  Float_t ZOuter; // track position (z component) at the tracker edge
427  Float_t TOuter; // track position (z component) at the tracker edge
428
429  Float_t L; // track path length
430   
431  Float_t D0;     // track signed transverse impact parameter
432  Float_t ErrorD0;    // signed error on the track signed transverse impact parameter
433 
434  Float_t DZ; // track transverse momentum
435  Float_t ErrorDZ; // track transverse momentum error
436 
437  Float_t P; // track transverse momentum
438  Float_t ErrorP; // track transverse momentum error
439 
440  Float_t PT; // track transverse momentum
441  Float_t ErrorPT; // track transverse momentum error
442 
443  Float_t CtgTheta; // track transverse momentum
444  Float_t ErrorCtgTheta; // track transverse momentum error
445 
446  Float_t Phi; // track azimuthal angle
447  Float_t ErrorPhi; // track azimuthal angle
448 
449  Float_t Xd;      // X coordinate of point of closest approach to vertex
450  Float_t Yd;      // Y coordinate of point of closest approach to vertex
451  Float_t Zd;      // Z coordinate of point of closest approach to vertex
452
453  TRef Particle; // reference to generated particle
454
455  static CompBase *fgCompare; //!
456  const CompBase *GetCompare() const { return fgCompare; }
457
458  TLorentzVector P4() const;
459
460  ClassDef(Track, 2)
461};
462
463//---------------------------------------------------------------------------
464
465class Tower: public SortableObject
466{
467public:
468  Float_t ET; // calorimeter tower transverse energy
469  Float_t Eta; // calorimeter tower pseudorapidity
470  Float_t Phi; // calorimeter tower azimuthal angle
471
472  Float_t E; // calorimeter tower energy
473
474  Float_t T; // ecal deposit time, averaged by sqrt(EM energy) over all particles, not smeared
475  Int_t NTimeHits; // number of hits contributing to time measurement
476
477  Float_t Eem; // calorimeter tower electromagnetic energy
478  Float_t Ehad; // calorimeter tower hadronic energy
479
480  Float_t Edges[4]; // calorimeter tower edges
481
482  TRefArray Particles; // references to generated particles
483
484  static CompBase *fgCompare; //!
485  const CompBase *GetCompare() const { return fgCompare; }
486
487  TLorentzVector P4() const;
488
489  ClassDef(Tower, 2)
490};
491
492//---------------------------------------------------------------------------
493
494class HectorHit: public SortableObject
495{
496public:
497  Float_t E; // reconstructed energy [GeV]
498
499  Float_t Tx; // angle of the momentum in the horizontal (x,z) plane [urad]
500  Float_t Ty; // angle of the momentum in the verical (y,z) plane [urad]
501
502  Float_t T; // time of flight to the detector [s]
503
504  Float_t X; // horizontal distance to the beam [um]
505  Float_t Y; // vertical distance to the beam [um]
506  Float_t S; // distance to the interaction point [m]
507
508  TRef Particle; // reference to generated particle
509
510  static CompBase *fgCompare; //!
511  const CompBase *GetCompare() const { return fgCompare; }
512
513  ClassDef(HectorHit, 1)
514};
515
516//---------------------------------------------------------------------------
517
518class Candidate: public SortableObject
519{
520  friend class DelphesFactory;
521
522public:
523  Candidate();
524
525  Int_t PID;
526
527  Int_t Status;
528  Int_t M1, M2, D1, D2;
529
530  Int_t Charge;
531
532  Float_t Mass;
533
534  Int_t IsPU;
535  Int_t IsRecoPU;
536
537  Int_t IsConstituent;
538
539  Int_t IsFromConversion;
540
541  UInt_t Flavor;
542  UInt_t FlavorAlgo;
543  UInt_t FlavorPhys;
544
545  UInt_t BTag;
546  UInt_t BTagAlgo;
547  UInt_t BTagPhys;
548
549  UInt_t TauTag;
550
551  Float_t Eem;
552  Float_t Ehad;
553
554  Float_t Edges[4];
555  Float_t DeltaEta;
556  Float_t DeltaPhi;
557
558  TLorentzVector Momentum, Position, InitialPosition, PositionError, Area;
559
560  Float_t L; // path length
561  Float_t D0;
562  Float_t ErrorD0;
563  Float_t DZ;
564  Float_t ErrorDZ;
565  Float_t P;
566  Float_t ErrorP;
567  Float_t PT;
568  Float_t ErrorPT;
569  Float_t CtgTheta;
570  Float_t ErrorCtgTheta;
571  Float_t Phi;
572  Float_t ErrorPhi;
573
574  Float_t Xd;
575  Float_t Yd;
576  Float_t Zd;
577
578  // tracking resolution
579 
580  Float_t TrackResolution;
581
582  // PileUpJetID variables
583
584  Int_t NCharged;
585  Int_t NNeutrals;
586  Float_t Beta;
587  Float_t BetaStar;
588  Float_t MeanSqDeltaR;
589  Float_t PTD;
590  Float_t FracPt[5];
591
592  // Timing information
593
594  Int_t NTimeHits;
595  std::vector< std::pair< Float_t, Float_t > > ECalEnergyTimePairs;
596
597  // Isolation variables
598
599  Float_t IsolationVar;
600  Float_t IsolationVarRhoCorr;
601  Float_t SumPtCharged;
602  Float_t SumPtNeutral;
603  Float_t SumPtChargedPU;
604  Float_t SumPt;
605
606  // vertex variables
607 
608  Int_t ClusterIndex;
609  Int_t ClusterNDF;
610  Double_t ClusterSigma;
611  Double_t SumPT2;
612  Double_t BTVSumPT2;
613  Double_t GenDeltaZ;
614  Double_t GenSumPT2;
615
616  // N-subjettiness variables
617
618  Float_t Tau[5];
619
620  // Other Substructure variables
621
622  TLorentzVector TrimmedP4[5]; // first entry (i = 0) is the total Trimmed Jet 4-momenta and from i = 1 to 4 are the trimmed subjets 4-momenta
623  TLorentzVector PrunedP4[5]; // first entry (i = 0) is the total Pruned Jet 4-momenta and from i = 1 to 4 are the pruned subjets 4-momenta
624  TLorentzVector SoftDroppedP4[5]; // first entry (i = 0) is the total SoftDropped Jet 4-momenta and from i = 1 to 4 are the pruned subjets 4-momenta
625
626  Int_t NSubJetsTrimmed; // number of subjets trimmed
627  Int_t NSubJetsPruned; // number of subjets pruned
628  Int_t NSubJetsSoftDropped; // number of subjets soft-dropped
629
630
631  static CompBase *fgCompare; //!
632  const CompBase *GetCompare() const { return fgCompare; }
633
634  void AddCandidate(Candidate *object);
635  TObjArray *GetCandidates();
636
637  Bool_t Overlaps(const Candidate *object) const;
638
639  virtual void Copy(TObject &object) const;
640  virtual TObject *Clone(const char *newname = "") const;
641  virtual void Clear(Option_t* option = "");
642
643private:
644  DelphesFactory *fFactory; //!
645  TObjArray *fArray; //!
646
647  void SetFactory(DelphesFactory *factory) { fFactory = factory; }
648
649  ClassDef(Candidate, 4)
650};
651
652#endif // DelphesClasses_h
653
654
Note: See TracBrowser for help on using the repository browser.