The following examples should help you to get to grips with using LHAPDF from C++:
Using and testing PDFs in C++
#include <iostream>
#include <fstream>
using namespace std;
int main(int argc, char* argv[]) {
if (argc < 3) {
cerr << "You must specify a PDF set and member number" << endl;
return 1;
}
const string setname = argv[1];
const string smem = argv[2];
const int imem = lexical_cast<int>(smem);
const PDF* pdf =
mkPDF(setname, imem);
vector<int> pids = pdf->flavors();
const double MINLOGX = -10;
const double MAXLOGX = 0;
const double DX = 0.01;
const int NX = (int) floor((MAXLOGX - MINLOGX)/DX) + 1;
const double MINLOGQ2 = 1;
const double MAXLOGQ2 = 8;
const double DQ2 = 0.01;
const int NQ2 = (int) floor((MAXLOGQ2 - MINLOGQ2)/DQ2) + 1;
for (int pid : pids) {
const string spid = lexical_cast<string>(pid);
const string filename = setname + "_" + smem + "_" + spid + ".dat";
ofstream f(filename.c_str());
for (int ix = 0; ix < NX; ++ix) {
const double log10x = (MINLOGX + ix*DX < -1e-3) ? MINLOGX + ix*DX : 0;
const double x = pow(10, log10x);
for (int iq2 = 0; iq2 < NQ2; ++iq2) {
const double log10q2 = MINLOGQ2 + iq2*DQ2;
const double q2 = pow(10, log10q2);
const double xf = pdf->xfxQ2(pid, x, q2);
f << x << " " << q2 << " " << xf << endl;
}
}
f.close();
}
for (double log10q2 = MINLOGQ2; log10q2 <= MAXLOGQ2; log10q2 += 0.2) {
const double q2 = pow(10, log10q2);
cout << "alpha_s(" << setprecision(1) << fixed << sqrt(q2) << " GeV) = "
<< setprecision(5) << pdf->alphasQ2(q2) << endl;
}
delete pdf;
return 0;
}
PDF * mkPDF(const std::string &setname, size_t member)
Namespace for all LHAPDF functions and classes.
Definition: AlphaS.h:14
#include <iostream>
#include <cmath>
#include <ctime>
#ifdef HAVE_MPI
#include <mpi.h>
#endif
using namespace std;
int main(int argc, char* argv[]) {
#ifdef HAVE_MPI
MPI_Init(&argc, &argv);
#endif
const double MINLOGX = -7.5;
const double MINLOGQ = 1;
const double MAXLOGQ = 3;
const double dx = 0.005;
const double dq = 0.005;
const clock_t start = clock();
#if LHAPDF_MAJOR_VERSION > 5
#define XFS(X, Q) pdf->xfxQ(X, Q, xfs)
#else
LHAPDF::initPDFSetByName("CT10nlo.LHgrid");
LHAPDF::initPDF(0);
#define XFS(X, Q) LHAPDF::xfx(X, Q, &xfs[0])
#endif
const clock_t init = clock();
vector<double> xfs; xfs.resize(13);
for (double log10x = MINLOGX; log10x <= 0.0; log10x += dx) {
for (double log10q = MINLOGQ; log10q <= MAXLOGQ; log10q += dq) {
XFS(pow(10, log10x), pow(10, log10q));
}
}
const clock_t end = clock();
std::cout << "Init = " << (init - start) << std::endl;
std::cout << "Query = " << (end - init) << std::endl;
std::cout << "Total = " << (end - start) << std::endl;
#if LHAPDF_MAJOR_VERSION > 5
delete pdf;
#endif
#ifdef HAVE_MPI
MPI_Finalize();
#endif
return 0;
}
PDF is the general interface for access to parton density information.
Definition: PDF.h:40
Using PDFs in Fortran
program example1
implicit double precision (a-h,o-z)
character name*64
double precision f(-6:6)
character*20 lparm
logical has_photon
dimension z(10), xx(10)
Data (z(i), i=1,10) /.05, .1, .2, .3, .4, .5, .6, .7, .8, .9/
Do i = 1, 10
xx(i) = z(i) **3
EndDo
name='CT10nlo.LHgrid'
call initpdfsetbyname(name)
qmz=91.18d0
write(*,*)
call numberpdf(n)
print *,'There are ',n,' PDF sets'
do i=0,n
write(*,*) '---------------------------------------------'
call initpdf(i)
write(*,*) 'PDF set ',i
call getxmin(i,xmin)
call getxmax(i,xmax)
call getq2min(i,q2min)
call getq2max(i,q2max)
print *,'xmin=',xmin,' xmax=',xmax,' Q2min=',q2min,' Q2max=',q2max
call getminmax(i,xmin,xmax,q2min,q2max)
print *,'xmin=',xmin,' xmax=',xmax,' Q2min=',q2min,' Q2max=',q2max
call setlhaparm('EXTRAPOLATE')
call getlhaparm(18,lparm)
print *,'lhaparm(18)=',lparm
write(*,*)
a=alphaspdf(qmz)
write(*,*) 'alpha_S(M_Z) = ',a
call getlam4m(1,i,xlam4)
call getlam5m(1,i,xlam5)
print *,' lambda5: ',xlam5, ' lambda4: ',xlam4
write(*,*)
write(*,*) 'x*up'
write(*,*) ' x Q=10 GeV Q=100 GeV Q=1000 GeV'
q = 50.0d0
print *,q
do ix=1,10
x = xx(ix)
if(has_photon()) then
print *,"This set has a photon"
call evolvepdfphoton(x,q,f,photon)
else
call evolvepdf(x,q,f)
endif
g = f(0)
u = f(2)
d = f(1)
s = f(3)
c = f(4)
b = f(5)
ubar = f(-2)
dbar = f(-1)
sbar = f(-3)
cbar = f(-4)
bbar = f(-5)
write(*,'(F7.4,13(1pE10.3))') x,u,d,ubar,dbar,s,sbar,c,cbar,b,bbar,g,photon
enddo
enddo
end program example1
Using and testing PDFs in Python
import lhapdf
p = lhapdf.mkPDF("CT10nlo", 0)
p = lhapdf.mkPDF("CT10nlo/0")
print(p.xfxQ2(21, 1e-3, 1e4))
for pid in p.flavors():
print(p.xfxQ(pid, 0.01, 91.2))
pset = lhapdf.getPDFSet("CT10nlo")
print(pset.description)
pcentral = pset.mkPDF(0)
pdfs1 = pset.mkPDFs()
pdfs2 = lhapdf.mkPDFs("CT10nlo")
import numpy as np
xs = [x for x in np.logspace(-7, 0, 5)]
qs = [q for q in np.logspace(1, 4, 4)]
gluon_xfs = np.empty([len(xs), len(qs)])
for ix, x in enumerate(xs):
for iq, q in enumerate(qs):
gluon_xfs[ix,iq] = p.xfxQ(21, x, q)
print(gluon_xfs)
print(lhapdf.version())
print(lhapdf.__version__)
lhapdf.pathsPrepend("/path/to/extra/pdfsets")
print(lhapdf.paths())
Handling LHAPDF v5/v6 compatibility in C++
#include <iostream>
int main() {
const double x = 1e-3, Q = 200;
#if LHAPDF_MAJOR_VERSION == 6
std::cout <<
"xf_g = " << pdf->
xfxQ(21, x, Q) << std::endl;
delete pdf;
#else
LHAPDF::initPDFSet("CT10nlo", LHAPDF::LHGRID, 0);
std::cout << "xf_g = " << LHAPDF::xfx(x, Q, 0) << std::endl;
#endif
return 0;
}
double xfxQ(int id, double x, double q) const
Get the PDF xf(x) value at (x,q) for the given PID.
Definition: PDF.h:108
Using and testing PDF grids in C++
#include "LHAPDF/GridPDF.h"
#include <iostream>
#include <fstream>
#ifdef HAVE_MPI
#include <mpi.h>
#endif
using namespace std;
void safeprint(
const LHAPDF::PDF& pdf,
const string& key) {
}
int main(int argc, char* argv[]) {
#ifdef HAVE_MPI
MPI_Init(&argc, &argv);
#endif
if (argc < 2) {
cout << "Usage: testgrid <PDFNAME=CT10nlo>" << endl;
}
const string setname = (argc < 2) ? "CT10nlo" : argv[1];
cout << endl;
safeprint(pdf, "Verbosity");
safeprint(pdf, "PdfDesc");
safeprint(pdf, "SetDesc");
cout <<
"Flavors (str) = " << pdf.
info().
get_entry(
"Flavors") << endl;
cout << "Flavors (ints) = ";
for (int f : pids) cout << f << " ";
cout << endl;
cout <<
"x0, Q0 = " << pdf.knotarray().
xf(21, 0, 0) << endl;
cout <<
"x1, Q0 = " << pdf.knotarray().
xf(21, 1, 0) << endl;
cout <<
"x0, Q1 = " << pdf.knotarray().
xf(21, 0, 1) << endl;
cout <<
"x1, Q1 = " << pdf.knotarray().
xf(21, 1, 1) << endl;
cout << pdf.
xfxQ(21, 0.7, 10.0) << endl;
cout << pdf.
xfxQ(21, 0.2, 126) << endl;
cout << pid <<
" = " << pdf.
xfxQ(pid, 0.2, 124) << endl;
}
ofstream f("pdf.dat");
for (double x = 0; x <= 1; x += 0.02) {
for (double log10q2 = 1; log10q2 < 5; log10q2 += 0.05) {
f << x <<
" " << log10q2 <<
" " << pdf.
xfxQ2(21, x, pow(10, log10q2)) << endl;
}
}
f.close();
cout << endl;
#ifdef HAVE_MPI
MPI_Finalize();
#endif
return 0;
}
A PDF defined via an interpolation grid.
Definition: GridPDF.h:19
T get_entry_as(const std::string &key) const
Definition: Info.h:139
double xf(int ix, int iq2, int ipid) const
convenient accessor to the grid values
Definition: KnotArray.h:68
const std::string & get_entry(const std::string &key) const
Retrieve a metadata string by key name.
bool has_key(const std::string &key) const
Can this Info object return a value for the given key? (it may be defined non-locally)
virtual const std::vector< int > & flavors() const
List of flavours defined by this PDF set.
Definition: PDF.h:418
double xfxQ2(int id, double x, double q2) const
Get the PDF xf(x) value at (x,q2) for the given PID.
PDFInfo & info()
Get the info class that actually stores and handles the metadata.
Definition: PDF.h:352
std::vector< std::string > paths()
Get the ordered list of search paths, from $LHAPDF_DATA_PATH and the install location.
std::string to_str(const T &val)
Make a string representation of val.
Definition: Utils.h:61
Testing the Info system in C++
#include "LHAPDF/Info.h"
#include "LHAPDF/Config.h"
#include "LHAPDF/PDFInfo.h"
#include "LHAPDF/PDFSet.h"
#include "LHAPDF/Factories.h"
#include <iostream>
#ifdef HAVE_MPI
#include <mpi.h>
#endif
using namespace std;
int main(int argc, char* argv[]) {
#ifdef HAVE_MPI
MPI_Init(&argc, &argv);
#endif
cout <<
"Verbosity: " << cfg.
get_entry(
"Verbosity") << endl;
cout <<
"New Verbosity from second Config: " << cfg2.
get_entry(
"Verbosity") << endl;
cout << "SetDesc: " << set.get_entry("SetDesc") << endl;
cout << "Verbosity from set: " << set.get_entry("Verbosity") << endl;
if (info.has_key("PdfDesc")) cout << "PdfDesc: " << info.get_entry("PdfDesc") << endl;
cout << "PdfType: " << info.get_entry("PdfType") << endl;
cout << "Verbosity from PDF: " << info.get_entry("Verbosity") << endl;
vector<int> pids = info.get_entry_as< vector<int> >("Flavors");
cout << "PIDs (1): "; for (int f : pids) { cout << f << " "; } cout << endl;
cout << "Testing PDFInfo for " << name << endl;
delete i;
}
#ifdef HAVE_MPI
MPI_Finalize();
#endif
return 0;
}
Metadata base class for PDFs, PDF sets, or global configuration.
Definition: Info.h:29
void set_entry(const std::string &key, const T &val)
Set a keyed value entry.
Definition: Info.h:158
virtual const std::string & get_entry(const std::string &key) const
Definition: Info.h:119
Metadata class for PDF members.
Definition: PDFInfo.h:18
Class for PDF-set metadata and manipulation.
Definition: PDFSet.h:105
PDFInfo * mkPDFInfo(const std::string &setname, size_t member)
const std::vector< std::string > & availablePDFSets()
Get the names of all available PDF sets in the search path.
Testing the LHAPDF ID indexing system in C++
#include "LHAPDF/PDFIndex.h"
#include <iostream>
#ifdef HAVE_MPI
#include <mpi.h>
#endif
using namespace std;
void lookup(int id) {
cout << "ID=" << id << " -> set=" << set_id.first << ", mem=" << set_id.second << endl;
}
int main(int argc, char* argv[]) {
#ifdef HAVE_MPI
MPI_Init(&argc, &argv);
#endif
lookup(10800);
lookup(10801);
lookup(10042);
lookup(10041);
lookup(10799);
lookup(12346);
#ifdef HAVE_MPI
MPI_Finalize();
#endif
return 0;
}
std::pair< std::string, int > lookupPDF(int lhaid)
Testing the path searching system in C++
#include "LHAPDF/Paths.h"
#include <iostream>
using namespace std;
#ifdef HAVE_MPI
#include <mpi.h>
#endif
int main(int argc, char* argv[]) {
#ifdef HAVE_MPI
MPI_Init(&argc, &argv);
#endif
cout << p << endl;
cout << "List of available PDFs:" << endl;
cout << " " << s << endl;
#ifdef HAVE_MPI
MPI_Finalize();
#endif
return 0;
}
std::string findFile(const std::string &target)
Testing the AlphaS classes in C++
#include <iostream>
#include <fstream>
#include <limits>
#ifdef HAVE_MPI
#include <mpi.h>
#endif
using namespace std;
int main(int argc, char* argv[]) {
#ifdef HAVE_MPI
MPI_Init(&argc, &argv);
#endif
AlphaS_Analytic as_ana;
as_ana.setOrderQCD(4);
as_ana.setQuarkMass(1, 0.0017);
as_ana.setQuarkMass(2, 0.0041);
as_ana.setQuarkMass(3, 0.1);
as_ana.setQuarkMass(4, 1.29);
as_ana.setQuarkMass(5, 4.1);
as_ana.setQuarkMass(6, 172.5);
as_ana.setLambda(3, 0.339);
as_ana.setLambda(4, 0.296);
as_ana.setLambda(5, 0.213);
AlphaS_ODE as_ode;
as_ode.setOrderQCD(5);
as_ode.setMZ(91);
as_ode.setAlphaSMZ(0.118);
as_ode.setQuarkMass(1, 0.0017);
as_ode.setQuarkMass(2, 0.0041);
as_ode.setQuarkMass(3, 0.1);
as_ode.setQuarkMass(4, 1.29);
as_ode.setQuarkMass(5, 4.1);
as_ode.setQuarkMass(6, 172.5);
AlphaS_Ipol as_ipol;
vector<double> qs = { 1.300000e+00, 1.300000e+00, 1.560453e+00, 1.873087e+00, 2.248357e+00, 2.698811e+00, 3.239513e+00, 3.888544e+00, 4.667607e+00, 5.602754e+00, 6.725257e+00, 8.072650e+00, 9.689992e+00, 1.163137e+01, 1.396169e+01, 1.675889e+01, 2.011651e+01, 2.414681e+01, 2.898459e+01, 3.479160e+01, 4.176203e+01, 5.012899e+01, 6.017224e+01, 7.222765e+01, 8.669834e+01, 1.040682e+02, 1.249181e+02, 1.499452e+02, 1.799865e+02, 2.160465e+02, 2.593310e+02, 3.112875e+02, 3.736534e+02, 4.485143e+02, 5.383733e+02, 6.462355e+02, 7.757077e+02, 9.311194e+02, 1.117668e+03, 1.341590e+03, 1.610376e+03, 1.933012e+03, 2.320287e+03, 2.785153e+03, 3.343154e+03, 4.012949e+03, 4.816936e+03, 4.816936e+03 };
vector<double> alphas = { 4.189466e-01, 4.189466e-01, 3.803532e-01, 3.482705e-01, 3.211791e-01, 2.979983e-01, 2.779383e-01, 2.604087e-01, 2.451922e-01, 2.324903e-01, 2.210397e-01, 2.106640e-01, 2.012187e-01, 1.925841e-01, 1.846600e-01, 1.773623e-01, 1.706194e-01, 1.643704e-01, 1.585630e-01, 1.531520e-01, 1.480981e-01, 1.433671e-01, 1.389290e-01, 1.347574e-01, 1.308290e-01, 1.271232e-01, 1.236215e-01, 1.203076e-01, 1.171667e-01, 1.141857e-01, 1.113525e-01, 1.086566e-01, 1.060881e-01, 1.036382e-01, 1.012990e-01, 9.906295e-02, 9.692353e-02, 9.487457e-02, 9.291044e-02, 9.102599e-02, 8.921646e-02, 8.747747e-02, 8.580498e-02, 8.419524e-02, 8.264479e-02, 8.115040e-02, 7.970910e-02, 7.970910e-02 };
as_ipol.setQValues(qs);
as_ipol.setAlphaSValues(alphas);
PDF* pdf =
mkPDF(
"CT10", 0);
const double inf = numeric_limits<double>::infinity();
ofstream fa("alphas_ana.dat"), fo("alphas_ode.dat"), fi("alphas_ipol.dat"), fc("alphas_ct10nlo.dat");
cout << endl;
for (double log10q = -0.5; log10q < 3; log10q += 0.05) {
const double q = pow(10, log10q);
const double as_ana_q = as_ana.alphasQ(q);
cout << fixed;
cout << "Q = " << setprecision(3) << q << " GeV" << endl;
cout << "Analytical solution: " << setprecision(3) << setw(6) << ( (as_ana_q > 2) ? inf : as_ana_q )
<< " num flavs = " << as_ana.numFlavorsQ(q) << endl;
fa << q << " " << as_ana_q << endl;
const double as_ode_q = as_ode.alphasQ(q);
cout << "ODE solution: " << setprecision(3) << setw(6) << ( (as_ode_q > 2) ? inf : as_ode_q )
<< " num flavs = " << as_ode.numFlavorsQ(q) << endl;
fo << q << " " << as_ode_q << endl;
}
fa.close(); fo.close(); fi.close(); fc.close();
delete pdf;
#ifdef HAVE_MPI
MPI_Finalize();
#endif
return 0;
}
Making a new PDF subclass in C++
#include "LHAPDF/PDF.h"
#include <iostream>
using namespace std;
AnalyticPDF() {
info().set_entry("Flavors", "-5,-4,-3,-2,-1,21,1,2,3,4,5");
}
double _xfxQ2(int id, double x, double q2) const {
if (abs(id) > 5 && id != 21) return 0;
return 0.15 * sin(20.0*x) * sin(20.0*q2);
}
void _xfxQ2(double x, double q2, std::vector<double>& ret) const {
for (int id(-5); id<5; ++id)
_xfxQ2(id,x,q2);
}
bool inRangeX(double x) const { return true; }
bool inRangeQ2(double q2) const { return true; }
};
int main(int argc, const char* argv[]) {
AnalyticPDF apdf;
for (double x = 0; x < 1.0; x += 0.1) {
for (double logq2 = 1; logq2 < 6; logq2 += 0.5) {
const double q2 = pow(10, logq2);
cout << x <<
" " << q2 <<
" " << pdf.
xfxQ2(21, x, q2) << endl;
}
}
return 0;
}
Calculating PDF uncertainties in C++
#include <random>
using namespace std;
void printUncs(
const LHAPDF::PDFSet& set,
const vector<double>& vals,
double cl,
const string& varname,
bool alt=
false) {
cout <<
"PDF uncertainties on " << varname <<
" computed with " << set.
name() <<
" at CL=" << cl <<
"%" << endl;
if (cl >= 0) cout << "Scaled PDF uncertainties using scale = " << err.scale << endl;
cout.precision(5);
cout << varname <<
" = " << err.
central <<
" +" << err.errplus <<
" -" << err.errminus <<
" (+-" << err.errsymm <<
")" << endl;
for (
size_t i = 0; i < errinfo.
qparts.size(); ++i) {
cout <<
" " << setw(12) << err.
errparts[i].first << setw(12) << err.
errparts[i].second <<
" " << errinfo.
qpartName(i) << endl;
}
}
int main(int argc, char* argv[]) {
if (argc < 2) {
cerr << "You must specify a PDF set: ./testpdfunc setname" << endl;
return 1;
}
const string setname = argv[1];
const size_t nmem = set.
size()-1;
double x = 0.1;
double Q = 100.0;
const vector<LHAPDF::PDF*> pdfs = set.
mkPDFs();
vector<double> xgAll, xuAll;
vector<string> pdftypes;
for (size_t imem = 0; imem <= nmem; imem++) {
xgAll.push_back(pdfs[imem]->xfxQ(21,x,Q));
xuAll.push_back(pdfs[imem]->xfxQ(2,x,Q));
pdftypes.push_back(pdfs[imem]->type());
}
cout << endl;
cout <<
"ErrorType: " << errinfo.
errtype << endl;
cout <<
"ErrorConfLevel: " << errinfo.
conflevel << endl;
const size_t npar = errinfo.
nmemPar();
if (npar > 0) cout << "Last " << npar << " members are parameter variations" << endl;
cout << endl;
cout << "Gluon distribution at Q = " << Q << " GeV (normal uncertainties)" << endl;
printUncs(set, xgAll, -1, "xg");
cout << endl;
cout << "Up-quark distribution at Q = " << Q << " GeV (normal uncertainties)" << endl;
printUncs(set, xuAll, -1, "xu");
cout << endl;
cout << "Self-correlation of xg = " << autocorr << endl;
cout << endl;
cout << "Correlation between xg and xu = " << corr << endl;
cout << endl;
cout << "Gluon distribution at Q = " << Q << " GeV (scaled uncertainties)" << endl;
printUncs(set, xgAll, 90, "xg");
cout << endl;
cout << "Up-quark distribution at Q = " << Q << " GeV (scaled uncertainties)" << endl;
printUncs(set, xuAll, 0, "xu");
cout << endl;
cout << "Gluon distribution at Q = " << Q << " GeV (median and 90% C.L.)" << endl;
printUncs(set, xgAll, 90, "xg", true);
cout << endl;
cout << "Up-quark distribution at Q = " << Q << " GeV (median and 68% C.L.)" << endl;
printUncs(set, xuAll, 68, "xu", true);
cout << endl;
const unsigned seed = 1234;
default_random_engine generator(seed);
normal_distribution<double> distribution;
const int nrand = 5;
for (int irand = 1; irand <= nrand; irand++) {
vector<double> randoms;
for (int ieigen=1; ieigen <= neigen; ieigen++) {
double r = distribution(generator);
randoms.push_back(r);
}
const bool symmetrise = true;
cout << "Random " << irand << ": xg = " << xgrand << ", xu = " << xurand << endl;
}
cout << endl;
}
return 0;
}
size_t size() const
Number of members in this set.
Definition: PDFSet.h:160
PDFErrInfo errorInfo() const
Get the structured decomposition of the error-type string.
double randomValueFromHessian(const std::vector< double > &values, const std::vector< double > &randoms, bool symmetrise=true) const
Generate a random value from Hessian values and Gaussian random numbers.
void _checkPdfType(const std::vector< string > &pdftypes) const
double correlation(const std::vector< double > &valuesA, const std::vector< double > &valuesB) const
Calculate the PDF correlation between valuesA and valuesB using appropriate formulae for this set.
PDFUncertainty uncertainty(const std::vector< double > &values, double cl=CL1SIGMA, bool alternative=false) const
Calculate the central value and PDF uncertainty on an observable.
Definition: PDFSet.h:330
std::string name() const
PDF set name.
Definition: PDFSet.h:123
std::string errorType() const
Get the type of PDF errors in this set (replicas, symmhessian, hessian, custom, etc....
Definition: PDFSet.h:143
void mkPDFs(std::vector< PTR > &pdfs) const
Definition: PDFSet.h:220
bool startswith(const std::string &s, const std::string &sub)
Does a string s start with the sub substring?
Definition: Utils.h:115
const double CL1SIGMA
CL percentage for a Gaussian 1-sigma.
Definition: PDFSet.h:20
Structure encoding the structure of the PDF error-set.
Definition: PDFSet.h:63
QuadParts qparts
Error-set quadrature parts.
Definition: PDFSet.h:77
size_t nmemCore() const
Number of core-set members.
double conflevel
Default confidence-level.
Definition: PDFSet.h:80
size_t nmemPar() const
Number of par-set members.
std::string errtype
Error-type annotation.
Definition: PDFSet.h:83
std::string qpartName(size_t iq) const
Calculated name of a quadrature part.
Structure for storage of uncertainty info calculated over a PDF error set.
Definition: PDFSet.h:36
double central
Variables for the central value, +ve, -ve & symmetrised errors, and a CL scalefactor.
Definition: PDFSet.h:49
ErrPairs errparts
Full error-breakdown of all quadrature uncertainty components, as (+,-) pairs.
Definition: PDFSet.h:57