lhapdf is hosted by Hepforge, IPPP Durham
LHAPDF 6.5.5
Loading...
Searching...
No Matches
Interpolator.h
1// -*- C++ -*-
2//
3// This file is part of LHAPDF
4// Copyright (C) 2012-2024 The LHAPDF collaboration (see AUTHORS for details)
5//
6#pragma once
7#ifndef LHAPDF_Interpolator_H
8#define LHAPDF_Interpolator_H
9
10#include "LHAPDF/Utils.h"
11#include "LHAPDF/KnotArray.h"
12
13namespace LHAPDF {
14
15
16 // Forward declaration
17 class GridPDF;
18
19
20 /// The general interface for interpolating between grid points
22 public:
23
24 /// Destructor to allow inheritance
25 virtual ~Interpolator() { }
26
27
28 /// @name Binding to a PDF object
29 ///@{
30
31 /// Bind to a GridPDF
32 void bind(const GridPDF* pdf) { _pdf = pdf; }
33
34 /// Unbind from GridPDF
35 void unbind() { _pdf = 0; }
36
37 /// Identify whether this Interpolator has an associated PDF
38 bool hasPDF() { return _pdf != 0; }
39
40 /// Get the associated GridPDF
41 const GridPDF& pdf() const { return *_pdf; }
42
43 ///@}
44
45
46 /// @name Interpolation methods
47 ///@{
48
49 /// Interpolate a single-point in (x,Q)
50 double interpolateXQ(int id, double x, double q) const {
51 return interpolateXQ2(id, x, q*q);
52 }
53
54 /// Interpolate a single-point in (x,Q2)
55 double interpolateXQ2(int id, double x, double q2) const;
56
57 void interpolateXQ2(double x, double q2, std::vector<double>& ret) const;
58
59 /// @todo Make an all-PID version of interpolateQ and Q2?
60
61 ///@}
62
63 /// The name of this type of interpolator
64 ///
65 /// @todo Would name() or scheme() be a better name? "Type" maybe confuses with the language type-system
66 std::string type() const {
67 return _type;
68 }
69
70 /// Set the interpolation type
71 void setType(std::string t) {
72 _type = t;
73 }
74
75
76 protected:
77
78 /// @brief Interpolate a single-point in (x,Q2), given x/Q2 values and subgrid indices.
79 ///
80 /// The key function to be overridden in derived classes: the subgrid and
81 /// x/Q2 index lookup (and their caching) are done centrally in the
82 /// Interpolator base class so do not need to be re-implemented in each
83 /// flavour of interpolator.
84 virtual double _interpolateXQ2(const KnotArray& grid, double x, size_t ix, double q2, size_t iq2, int id) const = 0;
85
86 virtual void _interpolateXQ2(const KnotArray& grid, double x, size_t ix, double q2, size_t iq2, std::vector<double>& ret) const = 0;
87
88
89 private:
90 const GridPDF* _pdf;
91 std::string _type;
92 };
93
94
95}
96
97#endif
A PDF defined via an interpolation grid.
Definition GridPDF.h:19
The general interface for interpolating between grid points.
Definition Interpolator.h:21
double interpolateXQ2(int id, double x, double q2) const
Interpolate a single-point in (x,Q2)
virtual double _interpolateXQ2(const KnotArray &grid, double x, size_t ix, double q2, size_t iq2, int id) const =0
Interpolate a single-point in (x,Q2), given x/Q2 values and subgrid indices.
void setType(std::string t)
Set the interpolation type.
Definition Interpolator.h:71
std::string type() const
Definition Interpolator.h:66
void unbind()
Unbind from GridPDF.
Definition Interpolator.h:35
bool hasPDF()
Identify whether this Interpolator has an associated PDF.
Definition Interpolator.h:38
void bind(const GridPDF *pdf)
Bind to a GridPDF.
Definition Interpolator.h:32
virtual ~Interpolator()
Destructor to allow inheritance.
Definition Interpolator.h:25
double interpolateXQ(int id, double x, double q) const
Interpolate a single-point in (x,Q)
Definition Interpolator.h:50
const GridPDF & pdf() const
Get the associated GridPDF.
Definition Interpolator.h:41
Namespace for all LHAPDF functions and classes.
Definition AlphaS.h:14