cepgen is hosted by Hepforge, IPPP Durham
CepGen N/A
Central exclusive processes event generator
CepGenEvent.h
Go to the documentation of this file.
1/*
2 * CepGen: a central exclusive processes event generator
3 * Copyright (C) 2016-2025 Laurent Forthomme
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 * 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 CepGenPythia8_CepGenEvent_h
20#define CepGenPythia8_CepGenEvent_h
21
22#include <Pythia8/Pythia.h>
23
24#include <unordered_map>
25
26namespace cepgen {
27 class RunParameters;
28 class Event;
29 class Particle;
30} // namespace cepgen
31
32namespace Pythia8 {
34 class CepGenEvent final : public LHAup {
35 public:
36 explicit CepGenEvent();
38
47 void feedEvent(const cepgen::Event& ev, const Type& type);
52 void setCrossSection(int id, double cross_section, double cross_section_err);
59 void setProcess(int id, double cross_section, double q2_scale, double alpha_qed, double alpha_qcd);
60
61 void addComments(const std::string& comments);
62
66 unsigned short cepgenId(unsigned short pythia_id) const;
70 unsigned short pythiaId(unsigned short cepgen_id) const;
73 int status = INVALID_ID,
74 const std::pair<int, int>& mothers = {0, 0},
75 const std::pair<int, int>& colours = {0, 0});
79 void addCorresp(unsigned short pythia_id, unsigned short cepgen_id);
80 void dumpCorresp() const;
81
82 static constexpr unsigned short INVALID_ID = 999;
83 static constexpr unsigned short MIN_COLOUR_INDEX = 501;
84
85 inline bool setInit() override { return true; }
86#if defined(PYTHIA_VERSION_INTEGER) && PYTHIA_VERSION_INTEGER >= 8200
87 bool setEvent(int) override { return true; }
88#else
89 bool setEvent(int, double) override { return true; }
90#endif
91
92 private:
93 std::pair<int, int> findMothers(const cepgen::Event& cepgen_event, const cepgen::Particle& cepgen_particle) const;
94 const double mp_, mp2_;
95 bool inel1_{false}, inel2_{false};
96 std::unordered_map<unsigned short, unsigned short> py_cg_corresp_;
97 const cepgen::RunParameters* params_{nullptr}; // borrowed
98 };
99} // namespace Pythia8
100#endif
Interfacing between CepGen and Pythia8 event definitions.
Definition CepGenEvent.h:34
bool setInit() override
Definition CepGenEvent.h:85
bool setEvent(int, double) override
Definition CepGenEvent.h:89
void feedEvent(const cepgen::Event &ev, const Type &type)
@ centralAndPartons
only include initiators and central system
@ centralAndFullBeamRemnants
include dissociated beam remnants and central system
@ centralAndBeamRemnants
include non-dissociated beam remnants and central system
void setCrossSection(int id, double cross_section, double cross_section_err)
Set the cross-section for a given process.
static constexpr unsigned short INVALID_ID
Invalid id association.
Definition CepGenEvent.h:82
void addCorresp(unsigned short pythia_id, unsigned short cepgen_id)
Register a new Pythia8 / CepGen particle mapping.
unsigned short pythiaId(unsigned short cepgen_id) const
Retrieve the Pythia8 particle index given its CepGen event id.
void initialise(const cepgen::RunParameters &)
Initialise this conversion object with CepGen parameters.
unsigned short cepgenId(unsigned short pythia_id) const
Retrieve the CepGen particle index given its Pythia8 event id.
void addComments(const std::string &comments)
Feed comments to the LHEF block.
void addCepGenParticle(const cepgen::Particle &part, int status=INVALID_ID, const std::pair< int, int > &mothers={0, 0}, const std::pair< int, int > &colours={0, 0})
Add a CepGen particle to the event content.
static constexpr unsigned short MIN_COLOUR_INDEX
Minimal colour indexing number.
Definition CepGenEvent.h:83
void dumpCorresp() const
Print all Pythia8 / CepGen Particles correspondences.
void setProcess(int id, double cross_section, double q2_scale, double alpha_qed, double alpha_qcd)
Specify new process attributes.
Container for the information on the in- and outgoing particles' kinematics.
Definition Event.h:26
Kinematic information for one particle.
Definition Particle.h:32
List of parameters used to start and run the simulation job.
Common namespace for this Monte Carlo generator.
Definition Handler.h:26