cepgen
is hosted by
Hepforge
,
IPPP Durham
CepGen
1.2.5
Central exclusive processes event generator
Loading...
Searching...
No Matches
KTIntegratedFlux.cpp
Go to the documentation of this file.
1
/*
2
* CepGen: a central exclusive processes event generator
3
* Copyright (C) 2018-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/CollinearFluxes/CollinearFlux.h
"
22
#include "
CepGen/Core/Exception.h
"
23
#include "
CepGen/Integration/AnalyticIntegrator.h
"
24
#include "
CepGen/KTFluxes/KTFlux.h
"
25
#include "
CepGen/Modules/AnalyticIntegratorFactory.h
"
26
#include "
CepGen/Modules/PartonFluxFactory.h
"
27
#include "
CepGen/Physics/PDG.h
"
28
#include "
CepGen/Utils/FunctionsWrappers.h
"
29
#include "
CepGen/Utils/Limits.h
"
30
31
namespace
cepgen
{
32
class
KTIntegratedFlux
:
public
CollinearFlux
{
33
public
:
34
explicit
KTIntegratedFlux
(
const
ParametersList
& params)
35
:
CollinearFlux
(params),
36
integr_(AnalyticIntegratorFactory::get().build(
steer
<
ParametersList
>(
"integrator"
))),
37
flux_(KTFluxFactory::get().build(
steer
<
ParametersList
>(
"ktFlux"
))),
38
kt2_range_(
steer
<
Limits
>(
"kt2range"
)),
39
func_q2_([this](double kt2, void* pars) {
40
const
auto
& args = *
static_cast<
std::pair<double, double>*
>
(pars);
41
return
flux_->fluxQ2(args.first, kt2, args.second);
42
}),
43
func_mx2_([
this
](
double
kt2,
void
* pars) {
44
const
auto
& args = *
static_cast<
std::pair<double, double>*
>
(pars);
45
return
flux_->fluxMX2(args.first, kt2, args.second);
46
}) {
47
if
(!flux_->ktFactorised())
48
throw
CG_FATAL
(
"GammaIntegrated"
) <<
"Input flux has to be unintegrated."
;
49
// initialise the functions to integrate
50
CG_INFO
(
"KTIntegratedFlux"
) <<
"kt flux-integrated collinear flux evaluator initialised.\n\t"
51
<<
"Analytical integrator: "
<< integr_->name() <<
"\n\t"
52
<<
"Q^2 integration range: "
<< kt2_range_ <<
" GeV^2\n\t"
53
<<
"Unintegrated flux: "
<< flux_->name() <<
"."
;
54
}
55
56
bool
fragmenting
() const override final {
return
flux_->fragmenting(); }
57
pdgid_t
partonPdgId
() const override final {
return
flux_->partonPdgId(); }
58
double
mass2
() const override final {
return
flux_->mass2(); }
59
60
static
ParametersDescription
description
() {
61
auto
desc =
CollinearFlux::description
();
62
desc.setDescription(
"kt-integr. coll.flux"
);
63
desc.add<
ParametersDescription
>(
"integrator"
, AnalyticIntegratorFactory::get().describeParameters(
"gsl"
))
64
.setDescription(
"Steering parameters for the analytical integrator"
);
65
desc.add<
ParametersDescription
>(
"ktFlux"
,
PartonFluxFactory::get
().
describeParameters
(
"BudnevElastic"
))
66
.setDescription(
"Type of unintegrated kT-dependent parton flux"
);
67
desc.add<
Limits
>(
"kt2range"
, {0., 1.e4})
68
.setDescription(
"kinematic range for the parton transverse virtuality, in GeV^2"
);
69
return
desc;
70
}
71
72
double
fluxQ2
(
double
x,
double
q2)
const override
{
73
if
(!
x_range_
.
contains
(x,
true
))
74
return
0.;
75
return
integr_->integrate(func_q2_, std::make_pair(x, q2), kt2_range_) / x;
76
}
77
78
double
fluxMX2
(
double
x,
double
mx2)
const override
{
79
if
(!
x_range_
.
contains
(x,
true
))
80
return
0.;
81
return
integr_->integrate(func_mx2_, std::make_pair(x, mx2), kt2_range_) / x;
82
}
83
84
private
:
85
const
std::unique_ptr<AnalyticIntegrator> integr_;
86
const
std::unique_ptr<KTFlux> flux_;
87
const
Limits
kt2_range_;
88
const
utils::Function1D
func_q2_, func_mx2_;
89
};
90
}
// namespace cepgen
91
REGISTER_COLLINEAR_FLUX
(
"KTIntegrated"
, KTIntegratedFlux);
AnalyticIntegratorFactory.h
AnalyticIntegrator.h
CollinearFlux.h
Exception.h
CG_FATAL
#define CG_FATAL(mod)
Definition
Exception.h:61
FunctionsWrappers.h
KTFlux.h
Limits.h
CG_INFO
#define CG_INFO(mod)
Definition
Message.h:216
PDG.h
PartonFluxFactory.h
REGISTER_COLLINEAR_FLUX
#define REGISTER_COLLINEAR_FLUX(name, obj)
Add a generic collinear parton flux evaluator builder definition.
Definition
PartonFluxFactory.h:25
cepgen::CollinearFlux
Definition
CollinearFlux.h:25
cepgen::CollinearFlux::description
static ParametersDescription description()
Definition
CollinearFlux.cpp:25
cepgen::KTIntegratedFlux
Definition
KTIntegratedFlux.cpp:32
cepgen::KTIntegratedFlux::fragmenting
bool fragmenting() const override final
Is initiator particle fragmenting after parton emission?
Definition
KTIntegratedFlux.cpp:56
cepgen::KTIntegratedFlux::fluxQ2
double fluxQ2(double x, double q2) const override
Compute the collinear flux for this x value and virtuality.
Definition
KTIntegratedFlux.cpp:72
cepgen::KTIntegratedFlux::partonPdgId
pdgid_t partonPdgId() const override final
Parton PDG identifier.
Definition
KTIntegratedFlux.cpp:57
cepgen::KTIntegratedFlux::KTIntegratedFlux
KTIntegratedFlux(const ParametersList ¶ms)
Definition
KTIntegratedFlux.cpp:34
cepgen::KTIntegratedFlux::fluxMX2
double fluxMX2(double x, double mx2) const override
Compute the collinear flux for this x value and remnant mass.
Definition
KTIntegratedFlux.cpp:78
cepgen::KTIntegratedFlux::description
static ParametersDescription description()
Definition
KTIntegratedFlux.cpp:60
cepgen::KTIntegratedFlux::mass2
double mass2() const override final
Initiator particle squared mass (in )
Definition
KTIntegratedFlux.cpp:58
cepgen::Limits
Validity interval for a variable.
Definition
Limits.h:28
cepgen::Limits::contains
bool contains(double val, bool exclude_boundaries=false) const
Check if value is inside limits' boundaries.
Definition
Limits.cpp:77
cepgen::ParametersDescription
A description object for parameters collection.
Definition
ParametersDescription.h:26
cepgen::ParametersList
Definition
ParametersList.h:52
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::Function1D
Wrapper to a 1-dimensional function with optional parameters.
Definition
FunctionsWrappers.h:29
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::PartonFluxFactory::get
static PartonFluxFactory & get()
Definition
PartonFluxFactory.h:53
cepgen::PartonFluxFactory::describeParameters
ParametersDescription describeParameters(const std::string &name, const ParametersList ¶ms=ParametersList()) const
Definition
PartonFluxFactory.cpp:27
CepGen
CollinearFluxes
KTIntegratedFlux.cpp
Generated on Mon Jul 29 2024 for CepGen by
1.9.7