LHAPDF is hosted by Hepforge, IPPP Durham
LHAPDF  6.1.6
PDF.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_PDF_H
8 #define LHAPDF_PDF_H
9 
10 #include "LHAPDF/PDFInfo.h"
11 #include "LHAPDF/PDFIndex.h"
12 #include "LHAPDF/Factories.h"
13 #include "LHAPDF/AlphaS.h"
14 #include "LHAPDF/Utils.h"
15 #include "LHAPDF/Paths.h"
16 #include "LHAPDF/Exceptions.h"
17 #include "LHAPDF/Version.h"
18 #include "LHAPDF/Config.h"
19 
20 namespace LHAPDF {
21 
22 
26  class PDF {
27  protected: //< These constructors should only be called by subclasses
28 
30  typedef unique_ptr<AlphaS> AlphaSPtr;
31 
33  PDF() : _forcePos(0) { }
34 
35 
36  public:
37 
39  virtual ~PDF() { }
40 
42 
43 
44  protected:
45 
46 
48 
49 
50  void _loadInfo(const std::string& mempath) {
51  if (mempath.empty())
52  throw UserError("Tried to initialize a PDF with a null data file path... oops");
53  _mempath = mempath;
54  _info = PDFInfo(mempath);
55  //_info = PDFInfo(_setname(), memberID());
57  if (_info.has_key("MinLHAPDFVersion")) {
58  if (_info.get_entry_as<int>("MinLHAPDFVersion") > LHAPDF_VERSION_CODE) {
59  throw VersionError("Current LHAPDF version " + to_str(LHAPDF_VERSION_CODE)
60  + " less than required " + _info.get_entry("MinLHAPDFVersion"));
61  }
62  }
64  const int v = verbosity();
65  if (v > 0) {
66  std::cout << "LHAPDF " << version() << " loading " << mempath << std::endl;
67  print(std::cout, v);
68  }
70  if (_info.get_entry_as<int>("DataVersion", -1) <= 0) {
71  std::cerr << "WARNING: this PDF is preliminary, unvalidated, and not for production use!" << std::endl;
72  }
73  }
74 
75  void _loadInfo(const std::string& setname, int member) {
76  const string searchpath = findpdfmempath(setname, member);
77  if (searchpath.empty())
78  throw UserError("Can't find a valid PDF " + setname + "/" + to_str(member));
79  _loadInfo(searchpath);
80  }
81 
82  void _loadInfo(int lhaid) {
83  const pair<string,int> setname_memid = lookupPDF(lhaid);
84  if (setname_memid.second == -1)
85  throw IndexError("Can't find a PDF with LHAPDF ID = " + to_str(lhaid));
86  _loadInfo(setname_memid.first, setname_memid.second);
87  }
88 
90 
91 
92  public:
93 
95 
96 
106  double xfxQ2(int id, double x, double q2) const {
107  // Physical x range check
108  if (!inPhysicalRangeX(x)) {
109  throw RangeError("Unphysical x given: " + to_str(x));
110  }
111  // Physical Q2 range check
112  if (!inPhysicalRangeQ2(q2)) {
113  throw RangeError("Unphysical Q2 given: " + to_str(q2));
114  }
115  // Treat PID = 0 as always equivalent to a gluon: query as PID = 21
116  const int id2 = (id != 0) ? id : 21; //< @note Treat 0 as an alias for 21
117  // Undefined PIDs
118  if (!hasFlavor(id2)) return 0.0;
119  // Call the delegated method in the concrete PDF object to calculate the in-range value
120  double xfx = _xfxQ2(id2, x, q2);
121  // Apply positivity forcing at the enabled level
122  switch (forcePositive()) {
123  case 0: break;
124  case 1: if (xfx < 0) xfx = 0; break;
125  case 2: if (xfx < 1e-10) xfx = 1e-10; break;
126  default: throw LogicError("ForcePositive value not in expected range!");
127  }
128  // Return
129  return xfx;
130  }
131 
132 
143  double xfxQ(int id, double x, double q) const {
144  return xfxQ2(id, x, q*q);
145  }
146 
147 
156  void xfxQ2(double x, double q2, std::map<int, double>& rtn) const {
157  rtn.clear();
158  for (int id : flavors()) rtn[id] = xfxQ2(id, x, q2);
159  }
160 
161 
170  void xfxQ(double x, double q, std::map<int, double>& rtn) const {
171  xfxQ2(x, q*q, rtn);
172  }
173 
174 
187  void xfxQ2(double x, double q2, std::vector<double>& rtn) const {
188  rtn.clear();
189  rtn.resize(13);
190  for (int i = 0; i < 13; ++i) {
191  const int id = i-6; // PID = 0 is automatically treated as PID = 21
192  rtn[i] = xfxQ2(id, x, q2);
193  }
194  }
195 
196 
209  void xfxQ(double x, double q, std::vector<double>& rtn) const {
210  xfxQ2(x, q*q, rtn);
211  }
212 
213 
222  std::map<int, double> xfxQ2(double x, double q2) const {
223  std::map<int, double> rtn;
224  xfxQ2(x, q2, rtn);
225  return rtn;
226  }
227 
240  std::map<int, double> xfxQ(double x, double q) const {
241  return xfxQ2(x, q*q);
242  }
243 
244 
245  protected:
246 
259  virtual double _xfxQ2(int id, double x, double q2) const = 0;
260 
262 
263 
264  public:
265 
267 
268 
270  virtual double xMin() {
271  if (info().has_key("XMin"))
272  return info().get_entry_as<double>("XMin");
273  return numeric_limits<double>::epsilon();
274  }
275 
277  virtual double xMax() {
278  if (info().has_key("XMax"))
279  return info().get_entry_as<double>("XMax");
280  return 1.0;
281  }
282 
285  virtual double qMin() {
286  return info().get_entry_as<double>("QMin", 0);
287  }
288 
291  virtual double qMax() {
292  return info().get_entry_as<double>("QMax", numeric_limits<double>::max());
293  }
294 
296  virtual double q2Min() {
297  return sqr(this->qMin());
298  }
299 
301  virtual double q2Max() {
302  // Explicitly re-access this from the info, to avoid an overflow from squaring double_max
303  return (info().has_key("QMax")) ? sqr(info().get_entry_as<double>("QMax")) : numeric_limits<double>::max();
304  }
305 
311  int forcePositive() const {
312  if (_forcePos < 0) //< Caching
313  _forcePos = info().get_entry_as<unsigned int>("ForcePositive", 0);
314  return _forcePos;
315  }
316 
321  bool inPhysicalRangeX(double x) const {
322  return x >= 0.0 && x <= 1.0;
323  }
324 
328  bool inPhysicalRangeQ2(double q2) const {
329  return q2 >= 0.0;
330  }
331 
335  bool inPhysicalRangeQ(double q) const {
336  return inPhysicalRangeQ2(q*q);
337  }
338 
340  bool inPhysicalRangeXQ2(double x, double q2) const {
341  return inPhysicalRangeX(x) && inPhysicalRangeQ2(q2);
342  }
343 
345  bool inPhysicalRangeXQ(double x, double q) const {
346  return inPhysicalRangeX(x) && inPhysicalRangeQ(q);
347  }
348 
356  virtual bool inRangeQ(double q) const {
357  return inRangeQ2(q*q);
358  }
359 
366  virtual bool inRangeQ2(double q2) const = 0;
367 
374  virtual bool inRangeX(double x) const = 0;
375 
377  virtual bool inRangeXQ(double x, double q) const {
378  return inRangeX(x) && inRangeQ(q);
379  }
380 
382  bool inRangeXQ2(double x, double q2) const {
383  return inRangeX(x) && inRangeQ2(q2);
384  }
385 
387 
388 
390 
391 
393  PDFInfo& info() { return _info; }
394 
396  const PDFInfo& info() const { return _info; }
397 
401  PDFSet& set() const {
402  return getPDFSet(_setname());
403  }
404 
406 
407 
409 
410 
414  int memberID() const {
415  const string memname = file_stem(_mempath);
416  assert(memname.length() > 5); // There must be more to the filename stem than just the _nnnn suffix
417  const int memid = lexical_cast<int>(memname.substr(memname.length()-4)); //< Last 4 chars should be the member number
418  return memid;
419  }
420 
424  int lhapdfID() const;
425 
427  std::string description() const {
428  return info().get_entry("PdfDesc", "");
429  }
430 
432  int dataversion() const {
433  return info().get_entry_as<int>("DataVersion", -1);
434  }
435 
437  std::string type() const {
438  return to_lower(info().get_entry("PdfType"));
439  }
440 
442 
443 
445  void print(std::ostream& os=std::cout, int verbosity=1) const;
446 
447 
449 
450 
459  virtual const std::vector<int>& flavors() const {
460  if (_flavors.empty()) {
461  _flavors = info().get_entry_as< vector<int> >("Flavors");
462  sort(_flavors.begin(), _flavors.end());
463  }
464  return _flavors;
465  }
466 
468  bool hasFlavor(int id) const {
469  const int id2 = (id != 0) ? id : 21; //< @note Treat 0 as an alias for 21
470  const vector<int>& ids = flavors();
471  return find(ids.begin(), ids.end(), id2) != ids.end();
472  }
473 
479  int orderQCD() const {
480  return info().get_entry_as<int>("OrderQCD");
481  }
483  int qcdOrder() const { return orderQCD(); }
484 
489  double quarkMass(int id) const;
490 
494  double quarkThreshold(int id) const;
495 
497 
498 
500 
501 
507  void setAlphaS(AlphaS* alphas) {
508  _alphas.reset(alphas);
509  }
510 
512  bool hasAlphaS() const {
513  return _alphas.get() != 0;
514  }
515 
518  return *_alphas;
519  }
520 
522  const AlphaS& alphaS() const {
523  return *_alphas;
524  }
525 
530  double alphasQ(double q) const {
531  return alphasQ2(q*q);
532  }
533 
538  double alphasQ2(double q2) const {
539  if (!hasAlphaS())
540  _alphas.reset( mkAlphaS(info()) );
541  return _alphas->alphasQ2(q2);
542  }
543 
545 
546 
547  protected:
548 
550  std::string _setname() const {
551  return basename(dirname(_mempath));
552  }
553 
555  std::string _mempath;
556 
559 
561  mutable vector<int> _flavors;
562 
564  mutable AlphaSPtr _alphas;
565 
571  mutable int _forcePos;
572 
573  };
574 
575 
576 }
577 #endif
std::string version()
Get the LHAPDF library version code (as a string)
Definition: Version.h:33
int memberID() const
PDF member local ID number.
Definition: PDF.h:414
Error for places where it should not have been possible to get to!
Definition: Exceptions.h:46
PDF is the general interface for access to parton density information.
Definition: PDF.h:26
virtual double xMin()
Minimum valid x value for this PDF.
Definition: PDF.h:270
int forcePositive() const
Check whether PDF is set to only return positive (definite) values or not.
Definition: PDF.h:311
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:143
std::string to_str(const T &val)
Make a string representation of val.
Definition: Utils.h:61
int lhapdfID() const
PDF member global LHAPDF ID number.
virtual double qMin()
Definition: PDF.h:285
double quarkMass(int id) const
Get a quark mass in GeV by PDG code (|PID| = 1-6 only)
bool inPhysicalRangeXQ2(double x, double q2) const
Check whether the given (x,Q2) is physically valid.
Definition: PDF.h:340
Error to be thrown when out of the valid range of a PDF.
Definition: Exceptions.h:38
virtual ~PDF()
Virtual destructor, to allow unfettered inheritance.
Definition: PDF.h:39
virtual double _xfxQ2(int id, double x, double q2) const =0
Calculate the PDF xf(x) value at (x,q2) for the given PID.
double xfxQ2(int id, double x, double q2) const
Get the PDF xf(x) value at (x,q2) for the given PID.
Definition: PDF.h:106
unique_ptr< AlphaS > AlphaSPtr
Internal convenience typedef for the AlphaS object handle.
Definition: PDF.h:30
virtual bool inRangeQ2(double q2) const =0
Grid range check for Q2.
Error to be raised when a newer LHAPDF version is needed.
Definition: Exceptions.h:102
std::string _setname() const
Get the set name from the member data file path (for internal use only)
Definition: PDF.h:550
virtual double q2Max()
Maximum valid Q2 value for this PDF (in GeV2).
Definition: PDF.h:301
int dataversion() const
Version of this PDF&#39;s data file.
Definition: PDF.h:432
vector< int > _flavors
Locally cached list of supported PIDs.
Definition: PDF.h:561
bool inRangeXQ2(double x, double q2) const
Combined range check for x and Q2.
Definition: PDF.h:382
double alphasQ(double q) const
Value of alpha_s(Q2) used by this PDF.
Definition: PDF.h:530
virtual double xMax()
Maximum valid x value for this PDF.
Definition: PDF.h:277
bool hasAlphaS() const
Check if an AlphaS calculator is set.
Definition: PDF.h:512
bool inPhysicalRangeQ2(double q2) const
Check whether the given Q2 is physically valid.
Definition: PDF.h:328
std::string description() const
Description of this PDF member.
Definition: PDF.h:427
PDFSet & getPDFSet(const std::string &setname)
T lexical_cast(const U &in)
Convert between any types via stringstream.
Definition: Utils.h:47
std::string basename(const std::string &p)
Get the basename (i.e. terminal file name) from a path p.
Definition: Utils.h:184
virtual bool inRangeXQ(double x, double q) const
Combined range check for x and Q.
Definition: PDF.h:377
virtual bool inRangeX(double x) const =0
Grid range check for x.
bool hasFlavor(int id) const
Checks whether id is a valid parton for this PDF.
Definition: PDF.h:468
AlphaS & alphaS()
Retrieve the AlphaS object for this PDF.
Definition: PDF.h:517
virtual const std::vector< int > & flavors() const
List of flavours defined by this PDF set.
Definition: PDF.h:459
PDFInfo & info()
Get the info class that actually stores and handles the metadata.
Definition: PDF.h:393
void xfxQ2(double x, double q2, std::map< int, double > &rtn) const
Get the PDF xf(x) value at (x,q2) for all supported PIDs.
Definition: PDF.h:156
virtual bool inRangeQ(double q) const
Grid range check for Q.
Definition: PDF.h:356
const PDFInfo & info() const
Get the info class that actually stores and handles the metadata (const version)
Definition: PDF.h:396
N sqr(const N &x)
Convenience function for squaring (of any type)
Definition: Utils.h:217
AlphaS * mkAlphaS(const Info &info)
Make an AlphaS object from an Info object.
T get_entry_as(const std::string &key) const
Definition: Info.h:128
void xfxQ2(double x, double q2, std::vector< double > &rtn) const
Get the PDF xf(x) value at (x,q2) for "standard" PIDs.
Definition: PDF.h:187
std::string to_lower(const std::string &s)
Convert a string to lower-case (not in-place)
Definition: Utils.h:138
bool inPhysicalRangeQ(double q) const
Check whether the given Q is physically valid.
Definition: PDF.h:335
std::string dirname(const std::string &p)
Get the dirname (i.e. path to the penultimate directory) from a path p.
Definition: Utils.h:190
double alphasQ2(double q2) const
Value of alpha_s(Q2) used by this PDF.
Definition: PDF.h:538
Class for PDF set metadata and manipulation.
Definition: PDFSet.h:39
PDF()
Force initialization of the only non-class member.
Definition: PDF.h:33
int qcdOrder() const
Definition: PDF.h:483
int verbosity()
Definition: Config.h:65
Metadata class for PDF members.
Definition: PDFInfo.h:18
PDFInfo _info
Metadata container.
Definition: PDF.h:558
std::string type() const
Get the type of PDF member that this object represents (central, error)
Definition: PDF.h:437
Namespace for all LHAPDF functions and classes.
Definition: AlphaS.h:14
Error to be raised when a LHAPDF ID indexing fails.
Definition: Exceptions.h:86
int orderQCD() const
Order of QCD at which this PDF has been constructed.
Definition: PDF.h:479
Calculator interface for computing alpha_s(Q2) in various ways.
Definition: AlphaS.h:23
std::string file_stem(const std::string &f)
Get the stem (i.e. part without a file extension) from a filename f.
Definition: Utils.h:196
bool has_key(const std::map< K, T > &container, const K &key)
Does the map<K,T> container have a key K key?
Definition: Utils.h:260
void print(std::ostream &os=std::cout, int verbosity=1) const
Summary printout.
pair< std::string, int > lookupPDF(int lhaid)
Definition: PDFIndex.h:51
int _forcePos
Cached flag for whether to return only positive (or postive definite) PDF values. ...
Definition: PDF.h:571
void setAlphaS(AlphaS *alphas)
Set the AlphaS calculator by pointer.
Definition: PDF.h:507
void xfxQ(double x, double q, std::vector< double > &rtn) const
Get the PDF xf(x) value at (x,q) for "standard" PIDs.
Definition: PDF.h:209
std::string _mempath
Member data file path.
Definition: PDF.h:555
bool inPhysicalRangeX(double x) const
Check whether the given x is physically valid.
Definition: PDF.h:321
std::map< int, double > xfxQ(double x, double q) const
Get the PDF xf(x) value at (x,q) for all supported PIDs.
Definition: PDF.h:240
std::map< int, double > xfxQ2(double x, double q2) const
Get the PDF xf(x) value at (x,q2) for all supported PIDs.
Definition: PDF.h:222
double quarkThreshold(int id) const
Get a flavor scale threshold in GeV by PDG code (|PID| = 1-6 only) Convenience interface to the Mass*...
bool inPhysicalRangeXQ(double x, double q) const
Check whether the given (x,Q) is physically valid.
Definition: PDF.h:345
const AlphaS & alphaS() const
Retrieve the AlphaS object for this PDF (const)
Definition: PDF.h:522
virtual double q2Min()
Minimum valid Q2 value for this PDF (in GeV2).
Definition: PDF.h:296
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 double qMax()
Maximum valid Q value for this PDF (in GeV).
Definition: PDF.h:291
Problem exists between keyboard and chair.
Definition: Exceptions.h:110
void _loadInfo(const std::string &mempath)
Definition: PDF.h:50
const std::string & get_entry(const std::string &key) const
Retrieve a metadata string by key name.
void xfxQ(double x, double q, std::map< int, double > &rtn) const
Get the PDF xf(x) value at (x,q) for all supported PIDs.
Definition: PDF.h:170
AlphaSPtr _alphas
Optionally loaded AlphaS object.
Definition: PDF.h:564