LHAPDF
is hosted by
Hepforge
,
IPPP Durham
LHAPDF
6.2.1
include
LHAPDF
GridPDF.h
1
// -*- C++ -*-
2
//
3
// This file is part of LHAPDF
4
// Copyright (C) 2012-2016 The LHAPDF collaboration (see AUTHORS for details)
5
//
6
#pragma once
7
#ifndef LHAPDF_GridPDF_H
8
#define LHAPDF_GridPDF_H
9
10
#include "LHAPDF/PDF.h"
11
#include "LHAPDF/Interpolator.h"
12
#include "LHAPDF/Extrapolator.h"
13
#include "LHAPDF/KnotArray.h"
14
15
namespace
LHAPDF
{
16
17
19
class
GridPDF
:
public
PDF
{
20
public
:
21
23
24
26
GridPDF
() {
27
_mempath
=
""
;
28
_info
=
PDFInfo
();
29
_forcePos
= -1;
30
}
31
39
GridPDF
(
const
std::string& path) {
40
_loadInfo(path);
// Sets _mempath
41
_loadPlugins
();
42
_loadData
(
_mempath
);
43
_forcePos
= -1;
44
}
45
47
GridPDF
(
const
std::string& setname,
int
member) {
48
_loadInfo(setname, member);
// Sets _mempath
49
_loadPlugins
();
50
_loadData
(
_mempath
);
51
_forcePos
= -1;
52
}
53
55
GridPDF
(
int
lhaid) {
56
_loadInfo(lhaid);
// Sets _mempath
57
_loadPlugins
();
58
_loadData
(
_mempath
);
59
_forcePos
= -1;
60
}
61
63
virtual
~GridPDF
() { }
64
66
67
68
protected
:
69
71
void
_loadInterpolator
();
72
74
void
_loadExtrapolator
();
75
77
void
_loadPlugins
() {
78
_loadAlphaS();
79
_loadInterpolator
();
80
_loadExtrapolator
();
81
}
82
84
void
_loadData
(
const
std::string& mempath);
85
86
87
public
:
88
90
91
97
void
setInterpolator
(
Interpolator
* ipol);
98
105
template
<
typename
INTERPOLATOR>
106
void
setInterpolator
(INTERPOLATOR ipol) {
107
setInterpolator
(
new
INTERPOLATOR(ipol));
108
}
109
114
void
setInterpolator
(
const
std::string& ipolname);
115
117
bool
hasInterpolator
()
const
{
return
bool(
_interpolator
); }
118
120
const
Interpolator
&
interpolator
()
const
;
121
122
128
void
setExtrapolator
(
Extrapolator
* xpol);
129
136
template
<
typename
EXTRAPOLATOR>
137
void
setExtrapolator
(EXTRAPOLATOR xpol) {
138
setExtrapolator
(
new
EXTRAPOLATOR(xpol));
139
}
140
145
void
setExtrapolator
(
const
std::string& xpolname);
146
148
bool
hasExtrapolator
()
const
{
return
bool(
_extrapolator
); }
149
151
const
Extrapolator
&
extrapolator
()
const
;
152
154
155
156
protected
:
157
159
double
_xfxQ2
(
int
id
,
double
x,
double
q2)
const
;
160
161
162
public
:
163
165
166
168
std::map<double, KnotArrayNF>&
knotarrays
() {
169
return
_knotarrays
;
170
}
171
173
const
KnotArrayNF
&
subgrid
(
double
q2)
const
;
174
176
const
KnotArray1F
&
subgrid
(
int
id
,
double
q2)
const
{
177
return
subgrid
(q2).
get_pid
(
id
);
178
}
179
183
const
vector<double>&
xKnots
()
const
{
184
const
KnotArrayNF
& subgrid1 =
_knotarrays
.begin()->second;
185
const
KnotArray1F
& grid1 = subgrid1.
get_first
();
186
return
grid1.
xs
();
187
}
188
192
const
vector<double>&
q2Knots
()
const
;
193
194
public
:
195
197
bool
inRangeX
(
double
x)
const
{
198
assert(!
xKnots
().empty());
199
if
(x <
xKnots
().front())
return
false
;
200
if
(x >
xKnots
().back())
return
false
;
201
return
true
;
202
}
203
205
bool
inRangeQ2
(
double
q2)
const
{
206
assert(!
q2Knots
().empty());
207
if
(q2 <
q2Knots
().front())
return
false
;
208
if
(q2 >
q2Knots
().back())
return
false
;
209
return
true
;
210
}
211
213
214
215
private
:
216
218
std::map<double, KnotArrayNF>
_knotarrays
;
219
220
// /// Caching vector of x knot values
221
// mutable std::vector<double> _xknots;
222
224
mutable
std::vector<double>
_q2knots
;
225
227
typedef
unique_ptr<Interpolator>
InterpolatorPtr
;
228
230
typedef
unique_ptr<Extrapolator>
ExtrapolatorPtr
;
231
233
mutable
InterpolatorPtr
_interpolator
;
234
236
mutable
ExtrapolatorPtr
_extrapolator
;
237
238
};
239
240
241
}
242
#endif
LHAPDF::PDF
PDF is the general interface for access to parton density information.
Definition:
PDF.h:26
LHAPDF::KnotArrayNF::get_pid
const KnotArray1F & get_pid(int id) const
Get the KnotArray1F for PID code id.
Definition:
KnotArray.h:183
LHAPDF::GridPDF::_loadExtrapolator
void _loadExtrapolator()
Load the PDF grid data block, based on current metadata.
LHAPDF::GridPDF::inRangeX
bool inRangeX(double x) const
Check if x is in the grid range.
Definition:
GridPDF.h:197
LHAPDF::GridPDF::ExtrapolatorPtr
unique_ptr< Extrapolator > ExtrapolatorPtr
Typedef of smart pointer for xpol memory handling.
Definition:
GridPDF.h:230
LHAPDF::GridPDF::_interpolator
InterpolatorPtr _interpolator
Associated interpolator (mutable to allow laziness)
Definition:
GridPDF.h:233
LHAPDF::KnotArrayNF
A collection of {KnotArray1F}s accessed by PID code.
Definition:
KnotArray.h:168
LHAPDF::GridPDF::_loadData
void _loadData(const std::string &mempath)
Load the PDF grid data block (not the metadata) from the given PDF member file.
LHAPDF::GridPDF::interpolator
const Interpolator & interpolator() const
Get the current interpolator.
LHAPDF::GridPDF::_loadPlugins
void _loadPlugins()
Load the alphaS, interpolator, and extrapolator based on current metadata.
Definition:
GridPDF.h:77
LHAPDF::GridPDF
A PDF defined via an interpolation grid.
Definition:
GridPDF.h:19
LHAPDF::GridPDF::setInterpolator
void setInterpolator(Interpolator *ipol)
Set the interpolator by pointer.
LHAPDF::GridPDF::GridPDF
GridPDF(int lhaid)
Constructor from an LHAPDF ID.
Definition:
GridPDF.h:55
LHAPDF::GridPDF::_q2knots
std::vector< double > _q2knots
Caching vector of Q2 knot values.
Definition:
GridPDF.h:224
LHAPDF::GridPDF::subgrid
const KnotArrayNF & subgrid(double q2) const
Get the N-flavour subgrid containing Q2 = q2.
LHAPDF::GridPDF::GridPDF
GridPDF(const std::string &setname, int member)
Constructor from a set name and member ID.
Definition:
GridPDF.h:47
LHAPDF::GridPDF::_loadInterpolator
void _loadInterpolator()
Load the interpolator, based on current metadata.
LHAPDF::GridPDF::setExtrapolator
void setExtrapolator(Extrapolator *xpol)
Set the extrapolator by pointer.
LHAPDF::GridPDF::inRangeQ2
bool inRangeQ2(double q2) const
Check if q2 is in the grid range.
Definition:
GridPDF.h:205
LHAPDF::GridPDF::setExtrapolator
void setExtrapolator(EXTRAPOLATOR xpol)
Set the extrapolator by value.
Definition:
GridPDF.h:137
LHAPDF::GridPDF::_xfxQ2
double _xfxQ2(int id, double x, double q2) const
Get PDF xf(x,Q2) value (via grid inter/extrapolators)
LHAPDF::GridPDF::~GridPDF
virtual ~GridPDF()
Virtual destructor to allow inheritance.
Definition:
GridPDF.h:63
LHAPDF::GridPDF::_knotarrays
std::map< double, KnotArrayNF > _knotarrays
Map of multi-flavour KnotArrays "binned" for lookup by low edge in Q2.
Definition:
GridPDF.h:218
LHAPDF::GridPDF::extrapolator
const Extrapolator & extrapolator() const
Get the current extrapolator.
LHAPDF::KnotArray1F::xs
const std::vector< double > & xs() const
x knot accessor
Definition:
KnotArray.h:63
LHAPDF::PDFInfo
Metadata class for PDF members.
Definition:
PDFInfo.h:18
LHAPDF::PDF::_info
PDFInfo _info
Metadata container.
Definition:
PDF.h:503
LHAPDF
Namespace for all LHAPDF functions and classes.
Definition:
AlphaS.h:14
LHAPDF::GridPDF::hasInterpolator
bool hasInterpolator() const
Find whether an extrapolator has been set on this PDF.
Definition:
GridPDF.h:117
LHAPDF::GridPDF::GridPDF
GridPDF(const std::string &path)
Constructor from a file path.
Definition:
GridPDF.h:39
LHAPDF::Interpolator
The general interface for interpolating between grid points.
Definition:
Interpolator.h:21
LHAPDF::GridPDF::knotarrays
std::map< double, KnotArrayNF > & knotarrays()
Directly access the knot arrays in non-const mode, for programmatic filling.
Definition:
GridPDF.h:168
LHAPDF::GridPDF::setInterpolator
void setInterpolator(INTERPOLATOR ipol)
Set the interpolator by value.
Definition:
GridPDF.h:106
LHAPDF::PDF::_forcePos
int _forcePos
Cached flag for whether to return only positive (or postive definite) PDF values. ...
Definition:
PDF.h:516
LHAPDF::KnotArrayNF::get_first
const KnotArray1F & get_first() const
Convenience accessor for any valid subgrid, to get access to the x/Q2/etc. arrays.
Definition:
KnotArray.h:189
LHAPDF::PDF::_mempath
std::string _mempath
Member data file path.
Definition:
PDF.h:500
LHAPDF::GridPDF::subgrid
const KnotArray1F & subgrid(int id, double q2) const
Get the 1-flavour subgrid for PID=id containing Q2 = q2.
Definition:
GridPDF.h:176
LHAPDF::GridPDF::xKnots
const vector< double > & xKnots() const
Return a representative list of interpolation knots in x.
Definition:
GridPDF.h:183
LHAPDF::GridPDF::_extrapolator
ExtrapolatorPtr _extrapolator
Associated extrapolator (mutable to allow laziness)
Definition:
GridPDF.h:236
LHAPDF::KnotArray1F
Internal storage class for PDF data point grids.
Definition:
KnotArray.h:19
LHAPDF::GridPDF::InterpolatorPtr
unique_ptr< Interpolator > InterpolatorPtr
Typedef of smart pointer for ipol memory handling.
Definition:
GridPDF.h:227
LHAPDF::Extrapolator
The general interface for extrapolating beyond grid boundaries.
Definition:
Extrapolator.h:20
LHAPDF::GridPDF::q2Knots
const vector< double > & q2Knots() const
Return a representative list of interpolation knots in Q2.
LHAPDF::GridPDF::hasExtrapolator
bool hasExtrapolator() const
Find whether an extrapolator has been set on this PDF.
Definition:
GridPDF.h:148
LHAPDF::GridPDF::GridPDF
GridPDF()
Default constructor, making an empty PDF to be populated by hand.
Definition:
GridPDF.h:26
Generated on Wed Apr 18 2018 09:35:57 for LHAPDF by
1.8.13