cepgen is hosted by Hepforge, IPPP Durham
CepGen 1.2.5
Central exclusive processes event generator
Loading...
Searching...
No Matches
EPACollinearFlux.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 * 2009-2012 Nicolas Schul, Jerome de Favereau de Jeneret
5 *
6 * This program is free software: you can redistribute it and/or modify
7 * it under the terms of the GNU General Public License as published by
8 * the Free Software Foundation, either version 3 of the License, or
9 * any later version.
10 *
11 * This program is distributed in the hope that it will be useful,
12 * but WITHOUT ANY WARRANTY; without even the implied warranty of
13 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14 * GNU General Public License for more details.
15 *
16 * You should have received a copy of the GNU General Public License
17 * along with this program. If not, see <http://www.gnu.org/licenses/>.
18 */
19
26#include "CepGen/Physics/PDG.h"
28#include "CepGen/Utils/Limits.h"
29
30namespace cepgen {
31 class EPACollinearFlux final : public CollinearFlux {
32 public:
33 explicit EPACollinearFlux(const ParametersList& params)
34 : CollinearFlux(params), ff_(FormFactorsFactory::get().build(steer<ParametersList>("formFactors"))) {}
35
37 auto desc = CollinearFlux::description();
38 desc.setDescription("EPA FF-dependent flux");
39 desc.add<ParametersDescription>("formFactors", FormFactorsFactory::get().describeParameters("StandardDipole"));
40 return desc;
41 }
42
43 bool fragmenting() const override { return ff_->fragmenting(); }
44 pdgid_t partonPdgId() const override { return PDG::photon; }
45
46 double fluxQ2(double x, double q2) const override {
47 if (!x_range_.contains(x, true))
48 return 0.;
49 const auto q2min = utils::kt::q2(x, 0., mass2());
50 if (q2min == 0. || q2 < q2min)
51 return 0.;
52 const auto form_factors = (*ff_)(q2);
53 return prefactor_ * ((1. - x) * (1. - q2min / q2) * form_factors.FE + 0.5 * x * x * form_factors.FM) / x;
54 }
55
56 private:
57 double mass2() const override { return mp2_; }
58 const std::unique_ptr<formfac::Parameterisation> ff_;
59 };
60} // namespace cepgen
61REGISTER_COLLINEAR_FLUX("EPAFlux", EPACollinearFlux);
#define REGISTER_COLLINEAR_FLUX(name, obj)
Add a generic collinear parton flux evaluator builder definition.
static ParametersDescription description()
bool fragmenting() const override
Is initiator particle fragmenting after parton emission?
double fluxQ2(double x, double q2) const override
Compute the collinear flux for this x value and virtuality.
static ParametersDescription description()
EPACollinearFlux(const ParametersList &params)
pdgid_t partonPdgId() const override
Parton PDG identifier.
bool contains(double val, bool exclude_boundaries=false) const
Check if value is inside limits' boundaries.
Definition Limits.cpp:77
@ photon
Definition PDG.h:41
A description object for parameters collection.
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.