cepgen is hosted by Hepforge, IPPP Durham
CepGen N/A
Central exclusive processes event generator
GridHandler.h
Go to the documentation of this file.
1/*
2 * CepGen: a central exclusive processes event generator
3 * Copyright (C) 2018-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 CepGen_Utils_GridHandler_h
20#define CepGen_Utils_GridHandler_h
21
22#include <gsl/gsl_interp.h>
23#include <gsl/gsl_spline.h>
24#include <gsl/gsl_version.h>
25#if defined(GSL_MAJOR_VERSION) && (GSL_MAJOR_VERSION > 2 || (GSL_MAJOR_VERSION == 2 && GSL_MINOR_VERSION >= 1))
26#include <gsl/gsl_spline2d.h>
27#define GSL_VERSION_ABOVE_2_1 1
28#endif
29
30#include <array>
31#include <map>
32#include <memory>
33
34#include "CepGen/Utils/Limits.h"
35
36namespace cepgen {
38 enum struct GridType { linear, logarithmic, square };
42 template <size_t D, size_t N = 1>
44 public:
45 explicit GridHandler(const GridType& grid_type);
46 virtual ~GridHandler() = default;
47
48 using coord_t = std::vector<double>;
49 using values_t = std::array<double, N>;
50
51 values_t eval(const coord_t& in_coords) const;
52
53 void insert(const coord_t& coord, values_t value);
54 inline std::map<coord_t, values_t> values() const { return values_raw_; }
55
56 void initialise();
57 std::array<Limits, D> boundaries() const;
58 std::array<double, D> min() const;
59 std::array<double, D> max() const;
60
61 protected:
63 std::map<coord_t, values_t> values_raw_;
65 std::vector<std::unique_ptr<gsl_interp_accel, void (*)(gsl_interp_accel*)> > accelerators_;
66 std::vector<std::unique_ptr<gsl_spline, void (*)(gsl_spline*)> > splines_1d_;
67#ifdef GSL_VERSION_ABOVE_2_1
69 std::vector<std::unique_ptr<gsl_spline2d, void (*)(gsl_spline2d*)> > splines_2d_;
70#endif
71 std::array<coord_t, D> coordinates_;
72 std::array<std::unique_ptr<double[]>, N> values_;
73
74 private:
75 void findIndices(const coord_t& coord, coord_t& min, coord_t& max) const;
77 struct grid_point_t : values_t {
78 grid_point_t(const values_t& arr) : values_t(arr) {}
79 grid_point_t operator*(double c) const;
80 grid_point_t operator+(const grid_point_t& rhs) const;
81 };
82 bool initialised_{false};
83 };
84} // namespace cepgen
85
86#endif
A generic class for grid interpolation.
Definition GridHandler.h:43
std::vector< std::unique_ptr< gsl_spline, void(*)(gsl_spline *)> > splines_1d_
Splines for linear interpolations.
Definition GridHandler.h:66
values_t eval(const coord_t &in_coords) const
Interpolate a point to a given coordinate.
void initialise()
Initialise the grid and all useful interpolators/accelerators.
std::map< coord_t, values_t > values_raw_
List of coordinates and associated value(s) in the grid Grid interpolation accelerator.
Definition GridHandler.h:63
std::array< std::unique_ptr< double[]>, N > values_
Values for all points in the grid.
Definition GridHandler.h:72
std::array< double, N > values_t
Value(s) at a given coordinate.
Definition GridHandler.h:49
void insert(const coord_t &coord, values_t value)
Insert a new value in the grid.
std::vector< std::unique_ptr< gsl_interp_accel, void(*)(gsl_interp_accel *)> > accelerators_
Definition GridHandler.h:65
std::array< coord_t, D > coordinates_
Coordinates building up the grid.
Definition GridHandler.h:71
std::array< Limits, D > boundaries() const
Grid boundaries (collection of (min,max))
virtual ~GridHandler()=default
std::vector< double > coord_t
Coordinates container.
Definition GridHandler.h:48
std::array< double, D > max() const
Highest bound of the grid coordinates.
const GridType grid_type_
Type of interpolation for the grid members.
Definition GridHandler.h:62
std::array< double, D > min() const
Lowest bound of the grid coordinates.
GridHandler(const GridType &grid_type)
Build a grid interpolator from a grid type.
std::map< coord_t, values_t > values() const
List of values in the grid.
Definition GridHandler.h:54
Common namespace for this Monte Carlo generator.
Definition Handler.h:26
GridType
Interpolation type for the grid coordinates.
Definition GridHandler.h:38