cepgen is hosted by Hepforge, IPPP Durham
CepGen 1.2.5
Central exclusive processes event generator
Loading...
Searching...
No Matches
ElasticKTFluxes.cpp
Go to the documentation of this file.
1/*
2 * CepGen: a central exclusive processes event generator
3 * Copyright (C) 2023-2024 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#include <cmath>
20
26#include "CepGen/Physics/PDG.h"
28
29namespace cepgen {
31 public:
32 explicit ElasticNucleonKTFlux(const ParametersList& params)
33 : KTFlux(params), ff_(FormFactorsFactory::get().build(steer<ParametersList>("formFactors"))) {
34 if (!ff_)
35 throw CG_FATAL("ElasticNucleonKTFlux")
36 << "Elastic kT flux requires a modelling of electromagnetic form factors!";
37 }
38
40 auto desc = KTFlux::description();
41 desc.setDescription("Nucl. el. photon emission");
42 desc.add<ParametersDescription>("formFactors", FormFactorsFactory::get().describeParameters("StandardDipole"));
43 return desc;
44 }
45 bool fragmenting() const override final { return false; }
46 double mass2() const override { return mp2_; }
47 pdgid_t partonPdgId() const override { return PDG::photon; }
48
49 double fluxMX2(double x, double kt2, double) const override {
50 if (!x_range_.contains(x))
51 return 0.;
52 const auto q2 = utils::kt::q2(x, kt2, mass2()), q2min = q2 - kt2 / (1. - x);
53 const double qnorm = 1. - q2min / q2;
54 const auto& formfac = (*ff_)(q2);
55 return prefactor_ * formfac.FE * qnorm * qnorm / q2;
56 }
57
58 protected:
60 const std::unique_ptr<formfac::Parameterisation> ff_;
61 };
62
67 desc.setDescription("Nucl. el. photon emission (Budnev flux)");
68 return desc;
69 }
70 double fluxMX2(double x, double kt2, double) const override final {
71 if (!x_range_.contains(x))
72 return 0.;
73 const auto q2 = utils::kt::q2(x, kt2, mass2()), q2min = q2 - kt2 / (1. - x);
74 const double qnorm = 1. - q2min / q2;
75 const auto& formfac = (*ff_)(q2);
76 const double f_D = formfac.FE * (1. - x) * qnorm;
77 const double f_C = formfac.FM;
78 return prefactor_ * (f_D + 0.5 * x * x * f_C) * (1. - x) / q2;
79 }
80 };
81
83 public:
85 : BudnevElasticNucleonKTFlux(params), ml2_(std::pow(PDG::get().mass(steer<pdgid_t>("pdgId")), 2)) {}
88 desc.setDescription("Lepton el. photon emission (Budnev flux)");
89 desc.add<ParametersDescription>("formFactors", FormFactorsFactory::get().describeParameters("PointLikeFermion"));
90 desc.add<pdgid_t>("pdgId", PDG::electron).setDescription("lepton flavour");
91 return desc;
92 }
93 double mass2() const override { return ml2_; }
94
95 private:
96 const double ml2_;
97 };
98
100 public:
101 explicit ElasticHeavyIonKTFlux(const ParametersList& params)
102 : ElasticNucleonKTFlux(params),
103 hi_(HeavyIon::fromPdgId(steer<pdgid_t>("heavyIon"))),
104 mass2_(hi_.mass() * hi_.mass()) {
105 CG_DEBUG("ElasticHeavyIonKTFlux") << "KT-factorised elastic photon-from-HI flux evaluator built for HI=" << hi_
106 << ", (mass=" << hi_.mass()
107 << "), electromagnetic form factors: " << ff_->parameters() << ".";
108 }
109
112 desc.setDescription("HI el. photon emission");
113 desc.addAs<pdgid_t, HeavyIon>("heavyIon", HeavyIon::Pb());
114 desc.add<ParametersDescription>("formFactors", FormFactorsFactory::get().describeParameters("HeavyIonDipole"));
115 return desc;
116 }
117
118 double mass2() const override { return mass2_; }
119
120 double fluxMX2(double x, double kt2, double mx2) const override {
121 const auto z = (unsigned short)hi_.Z;
122 return z * z * ElasticNucleonKTFlux::fluxMX2(x, kt2, mx2);
123 }
124
125 private:
126 const HeavyIon hi_;
127 const double mass2_;
128 };
129} // namespace cepgen
130
131REGISTER_KT_FLUX("Elastic", 0, ElasticNucleonKTFlux);
132REGISTER_KT_FLUX("BudnevElastic", 10, BudnevElasticNucleonKTFlux);
133REGISTER_KT_FLUX("BudnevElasticLepton", 12, BudnevElasticLeptonKTFlux);
134REGISTER_KT_FLUX("ElasticHeavyIon", 100, ElasticHeavyIonKTFlux);
#define CG_FATAL(mod)
Definition Exception.h:61
#define CG_DEBUG(mod)
Definition Message.h:220
#define REGISTER_KT_FLUX(name, id, obj)
Add a generic KT-factorised flux evaluator builder definition.
BudnevElasticLeptonKTFlux(const ParametersList &params)
static ParametersDescription description()
double mass2() const override
Initiator particle squared mass (in )
ElasticHeavyIonKTFlux(const ParametersList &params)
static ParametersDescription description()
double fluxMX2(double x, double kt2, double mx2) const override
Compute the kt-dependent flux for this x value and remnant mass.
double mass2() const override
Initiator particle squared mass (in )
bool fragmenting() const override final
Is initiator particle fragmenting after parton emission?
ElasticNucleonKTFlux(const ParametersList &params)
double fluxMX2(double x, double kt2, double) const override
Compute the kt-dependent flux for this x value and remnant mass.
static ParametersDescription description()
const std::unique_ptr< formfac::Parameterisation > ff_
Elastic form factors computation.
pdgid_t partonPdgId() const override
Parton PDG identifier.
double mass2() const override
Initiator particle squared mass (in )
static ParametersDescription description()
Definition KTFlux.cpp:26
bool contains(double val, bool exclude_boundaries=false) const
Check if value is inside limits' boundaries.
Definition Limits.cpp:77
A singleton holding all physics constants associated to particles.
Definition PDG.h:28
@ electron
Definition PDG.h:37
@ photon
Definition PDG.h:41
A description object for parameters collection.
ParametersDescription & add(const std::string &name, const T &def)
Add the description to a new parameter.
const double prefactor_
Definition PartonFlux.h:38
const double mp2_
Definition PartonFlux.h:39
const Limits x_range_
Definition PartonFlux.h:40
T steer(const std::string &key) const
Retrieve a parameters as previously steered.
Definition Steerable.h:39
double q2(double x, double kt2, double mi2, double mx2)
Compute the virtuality from longitudinal loss/transverse virtuality/diffractive mass.
Definition Utils.cpp:60
Common namespace for this Monte Carlo generator.
unsigned long long pdgid_t
Alias for the integer-like particle PDG id.
double fluxMX2(double x, double kt2, double) const override final
Compute the kt-dependent flux for this x value and remnant mass.
static ParametersDescription description()
Heavy ion container (Z+A)
Definition HeavyIon.h:44
static HeavyIon Pb()
Standard lead.
Definition HeavyIon.h:77
Element Z
Atomic number.
Definition HeavyIon.h:85
double mass() const
Total heavy ion mass, in GeV/c2.
Definition HeavyIon.h:56