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
21
#include "
CepGen/Core/Exception.h
"
22
#include "
CepGen/FormFactors/Parameterisation.h
"
23
#include "
CepGen/KTFluxes/KTFlux.h
"
24
#include "
CepGen/Modules/FormFactorsFactory.h
"
25
#include "
CepGen/Modules/PartonFluxFactory.h
"
26
#include "
CepGen/Physics/PDG.h
"
27
#include "
CepGen/Physics/Utils.h
"
28
29
namespace
cepgen
{
30
class
ElasticNucleonKTFlux
:
public
KTFlux
{
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
39
static
ParametersDescription
description
() {
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
63
struct
BudnevElasticNucleonKTFlux
:
public
ElasticNucleonKTFlux
{
64
using
ElasticNucleonKTFlux::ElasticNucleonKTFlux
;
65
static
ParametersDescription
description
() {
66
auto
desc =
ElasticNucleonKTFlux::description
();
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
82
class
BudnevElasticLeptonKTFlux
final :
public
BudnevElasticNucleonKTFlux
{
83
public
:
84
explicit
BudnevElasticLeptonKTFlux
(
const
ParametersList
& params)
85
:
BudnevElasticNucleonKTFlux
(params), ml2_(std::pow(
PDG
::get().mass(
steer
<
pdgid_t
>(
"pdgId"
)), 2)) {}
86
static
ParametersDescription
description
() {
87
auto
desc =
BudnevElasticNucleonKTFlux::description
();
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
99
class
ElasticHeavyIonKTFlux
final :
public
ElasticNucleonKTFlux
{
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
110
static
ParametersDescription
description
() {
111
auto
desc =
ElasticNucleonKTFlux::description
();
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
131
REGISTER_KT_FLUX
(
"Elastic"
, 0, ElasticNucleonKTFlux);
132
REGISTER_KT_FLUX
(
"BudnevElastic"
, 10, BudnevElasticNucleonKTFlux);
133
REGISTER_KT_FLUX
(
"BudnevElasticLepton"
, 12, BudnevElasticLeptonKTFlux);
134
REGISTER_KT_FLUX
(
"ElasticHeavyIon"
, 100, ElasticHeavyIonKTFlux);
Exception.h
CG_FATAL
#define CG_FATAL(mod)
Definition
Exception.h:61
FormFactorsFactory.h
Parameterisation.h
KTFlux.h
CG_DEBUG
#define CG_DEBUG(mod)
Definition
Message.h:220
PDG.h
PartonFluxFactory.h
REGISTER_KT_FLUX
#define REGISTER_KT_FLUX(name, id, obj)
Add a generic KT-factorised flux evaluator builder definition.
Definition
PartonFluxFactory.h:34
Utils.h
cepgen::BudnevElasticLeptonKTFlux
Definition
ElasticKTFluxes.cpp:82
cepgen::BudnevElasticLeptonKTFlux::BudnevElasticLeptonKTFlux
BudnevElasticLeptonKTFlux(const ParametersList ¶ms)
Definition
ElasticKTFluxes.cpp:84
cepgen::BudnevElasticLeptonKTFlux::description
static ParametersDescription description()
Definition
ElasticKTFluxes.cpp:86
cepgen::BudnevElasticLeptonKTFlux::mass2
double mass2() const override
Initiator particle squared mass (in )
Definition
ElasticKTFluxes.cpp:93
cepgen::ElasticHeavyIonKTFlux
Definition
ElasticKTFluxes.cpp:99
cepgen::ElasticHeavyIonKTFlux::ElasticHeavyIonKTFlux
ElasticHeavyIonKTFlux(const ParametersList ¶ms)
Definition
ElasticKTFluxes.cpp:101
cepgen::ElasticHeavyIonKTFlux::description
static ParametersDescription description()
Definition
ElasticKTFluxes.cpp:110
cepgen::ElasticHeavyIonKTFlux::fluxMX2
double fluxMX2(double x, double kt2, double mx2) const override
Compute the kt-dependent flux for this x value and remnant mass.
Definition
ElasticKTFluxes.cpp:120
cepgen::ElasticHeavyIonKTFlux::mass2
double mass2() const override
Initiator particle squared mass (in )
Definition
ElasticKTFluxes.cpp:118
cepgen::ElasticNucleonKTFlux
Definition
ElasticKTFluxes.cpp:30
cepgen::ElasticNucleonKTFlux::fragmenting
bool fragmenting() const override final
Is initiator particle fragmenting after parton emission?
Definition
ElasticKTFluxes.cpp:45
cepgen::ElasticNucleonKTFlux::ElasticNucleonKTFlux
ElasticNucleonKTFlux(const ParametersList ¶ms)
Definition
ElasticKTFluxes.cpp:32
cepgen::ElasticNucleonKTFlux::fluxMX2
double fluxMX2(double x, double kt2, double) const override
Compute the kt-dependent flux for this x value and remnant mass.
Definition
ElasticKTFluxes.cpp:49
cepgen::ElasticNucleonKTFlux::description
static ParametersDescription description()
Definition
ElasticKTFluxes.cpp:39
cepgen::ElasticNucleonKTFlux::ff_
const std::unique_ptr< formfac::Parameterisation > ff_
Elastic form factors computation.
Definition
ElasticKTFluxes.cpp:60
cepgen::ElasticNucleonKTFlux::partonPdgId
pdgid_t partonPdgId() const override
Parton PDG identifier.
Definition
ElasticKTFluxes.cpp:47
cepgen::ElasticNucleonKTFlux::mass2
double mass2() const override
Initiator particle squared mass (in )
Definition
ElasticKTFluxes.cpp:46
cepgen::KTFlux
Definition
KTFlux.h:25
cepgen::KTFlux::description
static ParametersDescription description()
Definition
KTFlux.cpp:26
cepgen::Limits::contains
bool contains(double val, bool exclude_boundaries=false) const
Check if value is inside limits' boundaries.
Definition
Limits.cpp:77
cepgen::PDG
A singleton holding all physics constants associated to particles.
Definition
PDG.h:28
cepgen::PDG::electron
@ electron
Definition
PDG.h:37
cepgen::PDG::photon
@ photon
Definition
PDG.h:41
cepgen::ParametersDescription
A description object for parameters collection.
Definition
ParametersDescription.h:26
cepgen::ParametersDescription::add
ParametersDescription & add(const std::string &name, const T &def)
Add the description to a new parameter.
Definition
ParametersDescription.h:59
cepgen::ParametersList
Definition
ParametersList.h:52
cepgen::PartonFlux::prefactor_
const double prefactor_
Definition
PartonFlux.h:38
cepgen::PartonFlux::mp2_
const double mp2_
Definition
PartonFlux.h:39
cepgen::PartonFlux::x_range_
const Limits x_range_
Definition
PartonFlux.h:40
cepgen::Steerable::steer
T steer(const std::string &key) const
Retrieve a parameters as previously steered.
Definition
Steerable.h:39
cepgen::utils::kt::q2
double q2(double x, double kt2, double mi2, double mx2)
Compute the virtuality from longitudinal loss/transverse virtuality/diffractive mass.
Definition
Utils.cpp:60
cepgen
Common namespace for this Monte Carlo generator.
Definition
CommandLineHandler.cpp:36
cepgen::pdgid_t
unsigned long long pdgid_t
Alias for the integer-like particle PDG id.
Definition
ParticleProperties.h:26
cepgen::BudnevElasticNucleonKTFlux
Definition
ElasticKTFluxes.cpp:63
cepgen::BudnevElasticNucleonKTFlux::fluxMX2
double fluxMX2(double x, double kt2, double) const override final
Compute the kt-dependent flux for this x value and remnant mass.
Definition
ElasticKTFluxes.cpp:70
cepgen::BudnevElasticNucleonKTFlux::description
static ParametersDescription description()
Definition
ElasticKTFluxes.cpp:65
cepgen::HeavyIon
Heavy ion container (Z+A)
Definition
HeavyIon.h:44
cepgen::HeavyIon::Pb
static HeavyIon Pb()
Standard lead.
Definition
HeavyIon.h:77
cepgen::HeavyIon::Z
Element Z
Atomic number.
Definition
HeavyIon.h:85
cepgen::HeavyIon::mass
double mass() const
Total heavy ion mass, in GeV/c2.
Definition
HeavyIon.h:56
CepGen
KTFluxes
ElasticKTFluxes.cpp
Generated on Mon Jul 29 2024 for CepGen by
1.9.7