hej
is hosted by
Hepforge
,
IPPP Durham
HEJ
2.1.4
High energy resummation for hadron colliders
Loading...
Searching...
No Matches
include
HEJ
Event.hh
Go to the documentation of this file.
1
8
#pragma once
9
10
#include <array>
11
#include <cstddef>
12
#include <iosfwd>
13
#include <iterator>
14
#include <unordered_map>
15
#include <utility>
16
#include <vector>
17
18
#include "boost/iterator/filter_iterator.hpp"
19
20
#include "fastjet/ClusterSequence.hh"
21
#include "fastjet/PseudoJet.hh"
22
23
#include "
HEJ/Constants.hh
"
24
#include "
HEJ/Parameters.hh
"
25
#include "
HEJ/Particle.hh
"
26
#include "
HEJ/event_types.hh
"
27
28
namespace
LHEF
{
29
class
HEPEUP;
30
class
HEPRUP;
31
}
32
33
namespace
fastjet
{
34
class
JetDefinition;
35
}
36
37
namespace
HEJ
{
38
struct
RNG;
39
struct
UnclusteredEvent;
40
47
class
Event
{
48
public
:
49
class
EventData
;
50
52
using
ConstPartonIterator
= boost::filter_iterator<
53
bool (*)(
Particle
const
&),
54
std::vector<Particle>::const_iterator
55
>;
57
using
ConstReversePartonIterator
= std::reverse_iterator<
58
ConstPartonIterator
>;
60
Event
() =
delete
;
64
[[deprecated(
"UnclusteredEvent got superseded by EventData"
)]]
65
Event
(
66
UnclusteredEvent
const
& ev,
67
fastjet::JetDefinition
const
&
jet_def
,
double
min_jet_pt
68
);
69
72
74
std::array<Particle, 2>
const
&
incoming
()
const
{
75
return
incoming_;
76
}
78
std::vector<Particle>
const
&
outgoing
()
const
{
79
return
outgoing_;
80
}
82
ConstPartonIterator
begin_partons
()
const
;
84
ConstPartonIterator
cbegin_partons
()
const
;
85
87
ConstPartonIterator
end_partons
()
const
;
89
ConstPartonIterator
cend_partons
()
const
;
90
92
ConstReversePartonIterator
rbegin_partons
()
const
;
94
ConstReversePartonIterator
crbegin_partons
()
const
;
96
ConstReversePartonIterator
rend_partons
()
const
;
98
ConstReversePartonIterator
crend_partons
()
const
;
99
101
105
std::unordered_map<std::size_t, std::vector<Particle>>
const
&
decays
()
106
const
{
107
return
decays_;
108
}
110
std::vector<fastjet::PseudoJet>
const
&
jets
()
const
{
111
return
jets_;
112
}
114
117
119
Parameters<EventParameters>
const
&
parameters
()
const
{
120
return
parameters_;
121
}
123
Parameters<EventParameters>
&
parameters
(){
124
return
parameters_;
125
}
126
128
EventParameters
const
&
central
()
const
{
129
return
parameters_.central;
130
}
132
EventParameters
&
central
(){
133
return
parameters_.central;
134
}
135
137
std::vector<EventParameters>
const
&
variations
()
const
{
138
return
parameters_.variations;
139
}
141
std::vector<EventParameters> &
variations
(){
142
return
parameters_.variations;
143
}
144
146
149
EventParameters
const
&
variations
(std::size_t i)
const
{
150
return
parameters_.variations.at(i);
151
}
153
156
EventParameters
&
variations
(std::size_t i){
157
return
parameters_.variations.at(i);
158
}
160
162
169
std::vector<int>
particle_jet_indices
(
170
std::vector<fastjet::PseudoJet>
const
&
jets
171
)
const
{
172
return
cs_.particle_jet_indices(
jets
);
173
}
175
std::vector<int>
particle_jet_indices
()
const
{
176
return
particle_jet_indices
(
jets
());
177
}
178
180
fastjet::JetDefinition
const
&
jet_def
()
const
{
181
return
cs_.jet_def();
182
}
183
185
double
min_jet_pt
()
const
{
186
return
min_jet_pt_;
187
}
188
190
event_type::EventType
type
()
const
{
191
return
type_;
192
}
193
195
203
bool
generate_colours
(
RNG
&
/*ran*/
);
204
206
216
bool
is_leading_colour
()
const
;
217
235
bool
valid_hej_state
(
236
double
soft_pt_regulator =
DEFAULT_SOFT_PT_REGULATOR
,
237
double
min_pt = 0.
238
)
const
;
239
240
private
:
243
248
Event
(
249
std::array<Particle, 2> &&
incoming
,
250
std::vector<Particle> &&
outgoing
,
251
std::unordered_map<std::size_t, std::vector<Particle>> &&
decays
,
252
Parameters<EventParameters>
&&
parameters
,
253
fastjet::JetDefinition
const
&
jet_def
,
254
double
min_jet_pt
255
);
256
258
using
PartonIterator = boost::filter_iterator<
259
bool (*)(
Particle
const
&),
260
std::vector<Particle>::iterator
261
>;
263
using
ReversePartonIterator = std::reverse_iterator<PartonIterator>;
264
266
PartonIterator
begin_partons
();
268
PartonIterator
end_partons
();
269
271
ReversePartonIterator
rbegin_partons
();
273
ReversePartonIterator
rend_partons
();
274
275
std::array<Particle, 2> incoming_;
276
std::vector<Particle> outgoing_;
277
std::unordered_map<std::size_t, std::vector<Particle>> decays_;
278
std::vector<fastjet::PseudoJet> jets_;
279
Parameters<EventParameters>
parameters_;
280
fastjet::ClusterSequence cs_;
281
double
min_jet_pt_;
282
event_type::EventType
type_;
283
};
// end class Event
284
286
class
Event::EventData
{
287
public
:
289
EventData
() =
default
;
291
EventData
(LHEF::HEPEUP
const
& hepeup);
293
EventData
(
294
std::array<Particle, 2>
incoming
,
295
std::vector<Particle>
outgoing
,
296
std::unordered_map<std::size_t, std::vector<Particle>>
decays
,
297
Parameters<EventParameters>
parameters
298
):
299
incoming
(std::move(
incoming
)),
outgoing
(std::move(
outgoing
)),
300
decays
(std::move(
decays
)),
parameters
(std::move(
parameters
))
301
{}
302
304
315
Event
cluster
(
316
fastjet::JetDefinition
const
&
jet_def
,
double
min_jet_pt
);
317
319
Event
operator()
(
320
fastjet::JetDefinition
const
&
jet_def
,
double
const
min_jet_pt
){
321
return
cluster
(
jet_def
,
min_jet_pt
);
322
}
323
325
void
sort
();
326
328
333
void
reconstruct_intermediate
();
334
336
std::array<Particle, 2>
incoming
;
338
std::vector<Particle>
outgoing
;
340
std::unordered_map<std::size_t, std::vector<Particle>>
decays
;
342
Parameters<EventParameters>
parameters
;
343
};
// end class EventData
344
346
std::ostream&
operator<<
(std::ostream & os,
Event
const
& ev);
347
349
double
shat
(
Event
const
& ev);
350
352
LHEF::HEPEUP
to_HEPEUP
(
Event
const
& event, LHEF::HEPRUP *
/*heprup*/
);
353
354
// put deprecated warning at the end, so don't get the warning inside Event.hh,
355
// additionally doxygen can not identify [[deprecated]] correctly
356
struct
[[deprecated(
"UnclusteredEvent will be replaced by EventData"
)]]
357
UnclusteredEvent
;
361
struct
UnclusteredEvent
{
363
UnclusteredEvent
() =
default
;
365
UnclusteredEvent
(LHEF::HEPEUP
const
& hepeup);
366
367
std::array<Particle, 2>
incoming
;
368
std::vector<Particle>
outgoing
;
370
std::unordered_map<std::size_t, std::vector<Particle>>
decays
;
372
EventParameters
central
;
373
std::vector<EventParameters>
variations
;
374
};
375
376
}
// namespace HEJ
Constants.hh
Header file defining all global constants used for HEJ.
Parameters.hh
Containers for Parameter variations, e.g. different Weights.
Particle.hh
Contains the particle struct.
HEJ::Event::EventData
Class to store general Event setup, i.e. Phase space and weights.
Definition:
Event.hh:286
HEJ::Event::EventData::parameters
Parameters< EventParameters > parameters
Parameters, e.g. scale or inital weight.
Definition:
Event.hh:342
HEJ::Event::EventData::operator()
Event operator()(fastjet::JetDefinition const &jet_def, double const min_jet_pt)
Alias for cluster()
Definition:
Event.hh:319
HEJ::Event::EventData::EventData
EventData()=default
Default Constructor.
HEJ::Event::EventData::EventData
EventData(std::array< Particle, 2 > incoming, std::vector< Particle > outgoing, std::unordered_map< std::size_t, std::vector< Particle > > decays, Parameters< EventParameters > parameters)
Constructor with all values given.
Definition:
Event.hh:293
HEJ::Event::EventData::outgoing
std::vector< Particle > outgoing
Outcoing particles.
Definition:
Event.hh:338
HEJ::Event::EventData::decays
std::unordered_map< std::size_t, std::vector< Particle > > decays
Particle decays in the format {outgoing index, decay products}.
Definition:
Event.hh:340
HEJ::Event::EventData::cluster
Event cluster(fastjet::JetDefinition const &jet_def, double min_jet_pt)
Generate an Event from the stored EventData.
HEJ::Event::EventData::incoming
std::array< Particle, 2 > incoming
Incoming particles.
Definition:
Event.hh:336
HEJ::Event::EventData::sort
void sort()
Sort particles in rapidity.
HEJ::Event::EventData::reconstruct_intermediate
void reconstruct_intermediate()
Reconstruct intermediate particles from final-state leptons.
HEJ::Event::EventData::EventData
EventData(LHEF::HEPEUP const &hepeup)
Constructor from LesHouches event information.
HEJ::Event
An event with clustered jets.
Definition:
Event.hh:47
HEJ::Event::particle_jet_indices
std::vector< int > particle_jet_indices(std::vector< fastjet::PseudoJet > const &jets) const
Indices of the jets the outgoing partons belong to.
Definition:
Event.hh:169
HEJ::Event::end_partons
ConstPartonIterator end_partons() const
Iterator to the end of the outgoing partons.
HEJ::Event::cbegin_partons
ConstPartonIterator cbegin_partons() const
Iterator to the first outgoing parton.
HEJ::Event::central
EventParameters const & central() const
Central parameter choice (const version)
Definition:
Event.hh:128
HEJ::Event::variations
std::vector< EventParameters > & variations()
Parameter (scale) variations.
Definition:
Event.hh:141
HEJ::Event::incoming
std::array< Particle, 2 > const & incoming() const
Incoming particles.
Definition:
Event.hh:74
HEJ::Event::begin_partons
ConstPartonIterator begin_partons() const
Iterator to the first outgoing parton.
HEJ::Event::ConstPartonIterator
boost::filter_iterator< bool(*)(Particle const &), std::vector< Particle >::const_iterator > ConstPartonIterator
Iterator over partons.
Definition:
Event.hh:55
HEJ::Event::variations
EventParameters const & variations(std::size_t i) const
Parameter (scale) variation (const version)
Definition:
Event.hh:149
HEJ::Event::variations
std::vector< EventParameters > const & variations() const
Parameter (scale) variations (const version)
Definition:
Event.hh:137
HEJ::Event::is_leading_colour
bool is_leading_colour() const
Check that current colours are leading in the high energy limit.
HEJ::Event::decays
std::unordered_map< std::size_t, std::vector< Particle > > const & decays() const
Particle decays.
Definition:
Event.hh:105
HEJ::Event::parameters
Parameters< EventParameters > const & parameters() const
All chosen parameter, i.e. scale choices (const version)
Definition:
Event.hh:119
HEJ::Event::Event
Event(UnclusteredEvent const &ev, fastjet::JetDefinition const &jet_def, double min_jet_pt)
HEJ::Event::min_jet_pt
double min_jet_pt() const
Minimum jet transverse momentum.
Definition:
Event.hh:185
HEJ::Event::cend_partons
ConstPartonIterator cend_partons() const
Iterator to the end of the outgoing partons.
HEJ::Event::valid_hej_state
bool valid_hej_state(double soft_pt_regulator=DEFAULT_SOFT_PT_REGULATOR, double min_pt=0.) const
Check if given event could have been produced by HEJ.
HEJ::Event::crend_partons
ConstReversePartonIterator crend_partons() const
Reverse Iterator to the first outgoing parton.
HEJ::Event::rend_partons
ConstReversePartonIterator rend_partons() const
Reverse Iterator to the first outgoing parton.
HEJ::Event::central
EventParameters & central()
Central parameter choice.
Definition:
Event.hh:132
HEJ::Event::rbegin_partons
ConstReversePartonIterator rbegin_partons() const
Reverse Iterator to the first outgoing parton.
HEJ::Event::jet_def
fastjet::JetDefinition const & jet_def() const
Jet definition used for clustering.
Definition:
Event.hh:180
HEJ::Event::crbegin_partons
ConstReversePartonIterator crbegin_partons() const
Reverse Iterator to the first outgoing parton.
HEJ::Event::type
event_type::EventType type() const
Event type.
Definition:
Event.hh:190
HEJ::Event::jets
std::vector< fastjet::PseudoJet > const & jets() const
The jets formed by the outgoing partons, sorted in rapidity.
Definition:
Event.hh:110
HEJ::Event::ConstReversePartonIterator
std::reverse_iterator< ConstPartonIterator > ConstReversePartonIterator
Reverse Iterator over partons.
Definition:
Event.hh:58
HEJ::Event::generate_colours
bool generate_colours(RNG &)
Give colours to each particle.
HEJ::Event::parameters
Parameters< EventParameters > & parameters()
All chosen parameter, i.e. scale choices.
Definition:
Event.hh:123
HEJ::Event::Event
Event()=delete
No default Constructor.
HEJ::Event::particle_jet_indices
std::vector< int > particle_jet_indices() const
particle_jet_indices() of the Event jets()
Definition:
Event.hh:175
HEJ::Event::variations
EventParameters & variations(std::size_t i)
Parameter (scale) variation.
Definition:
Event.hh:156
HEJ::Event::outgoing
std::vector< Particle > const & outgoing() const
Outgoing particles.
Definition:
Event.hh:78
event_types.hh
Define different types of events.
HEJ::event_type::EventType
EventType
Possible event types.
Definition:
event_types.hh:19
HEJ
Main HEJ 2 Namespace.
Definition:
mainpage.dox:1
HEJ::to_HEPEUP
LHEF::HEPEUP to_HEPEUP(Event const &event, LHEF::HEPRUP *)
Convert an event to a LHEF::HEPEUP.
HEJ::shat
double shat(Event const &ev)
Square of the partonic centre-of-mass energy .
HEJ::operator<<
std::ostream & operator<<(std::ostream &os, const CrossSectionAccumulator &xs)
Print CrossSectionAccumulator to stream.
HEJ::DEFAULT_SOFT_PT_REGULATOR
constexpr double DEFAULT_SOFT_PT_REGULATOR
default value for the maximal pt fraction of soft radiation in any tagging jets This cut is needed to...
Definition:
Constants.hh:29
LHEF
Definition:
Analysis.hh:14
fastjet
Definition:
Event.hh:33
HEJ::EventParameters
Event parameters.
Definition:
Parameters.hh:107
HEJ::Parameters
Collection of parameters, e.g. Weights, assigned to a single event.
Definition:
Parameters.hh:26
HEJ::Particle
Class representing a particle.
Definition:
Particle.hh:24
HEJ::RNG
Interface for random number generator.
Definition:
RNG.hh:19
HEJ::UnclusteredEvent
Definition:
Event.hh:361
HEJ::UnclusteredEvent::outgoing
std::vector< Particle > outgoing
Definition:
Event.hh:368
HEJ::UnclusteredEvent::decays
std::unordered_map< std::size_t, std::vector< Particle > > decays
Particle decays in the format {outgoing index, decay products}.
Definition:
Event.hh:370
HEJ::UnclusteredEvent::UnclusteredEvent
UnclusteredEvent(LHEF::HEPEUP const &hepeup)
Constructor from LesHouches event information.
HEJ::UnclusteredEvent::incoming
std::array< Particle, 2 > incoming
Definition:
Event.hh:367
HEJ::UnclusteredEvent::central
EventParameters central
Central parameter (e.g. scale) choice.
Definition:
Event.hh:372
HEJ::UnclusteredEvent::UnclusteredEvent
UnclusteredEvent()=default
Default Constructor.
HEJ::UnclusteredEvent::variations
std::vector< EventParameters > variations
Definition:
Event.hh:373
Generated by
1.9.5