include/LHAPDF/Reweighting.h

0001 // -*- C++ -*-
0002 //
0003 // This file is part of LHAPDF
0004 // Copyright (C) 2012-2022 The LHAPDF collaboration (see AUTHORS for details)
0005 //
0006 #pragma once
0007 #ifndef LHAPDF_Reweighting_H
0008 #define LHAPDF_Reweighting_H
0009
0010 #include "LHAPDF/PDF.h"
0011 #include "LHAPDF/PDFSet.h"
0012
0013 namespace LHAPDF {
0014
0015
0016   /// @defgroup reweight PDF reweighting
0017   /// @{
0018
0019   namespace {
0020     inline bool _checkAlphasQ2(double Q2, const PDF& pdfa, const PDF& pdfb, double aschk) {
0021       if (aschk < 0) return true;
0022       const double as_a = pdfa.alphasQ2(Q2);
0023       const double as_b = pdfb.alphasQ2(Q2);
0024       if (2 * std::abs(as_a - as_b) / (std::abs(as_a) + std::abs(as_b)) < aschk) return true;
0025       std::cerr << "WARNING: alpha_s(Q2) mismatch in PDF reweighting "
0026                 << "at Q2 = " << Q2 << " GeV2:\n  "
0027                 << as_a << " for " << pdfa.set().name() << "/" << pdfa.memberID() << " vs. "
0028                 << as_b << " for " << pdfb.set().name() << "/" << pdfb.memberID()
0029                 << std::endl;
0030       return false;
0031     }
0032   }
0033
0034
0035   /// @defgroup reweight_single Single-beam reweighting
0036   /// @{
0037
0038   /// Get the PDF reweighting factor for a single beam with id,x,Q parameters, from basepdf to newpdf
0039   ///
0040   /// @note For NLO calculations, in general different PDF values enter for each counterterm: be careful.
0041   inline double weightxQ2(int id, double x, double Q2, const PDF& basepdf, const PDF& newpdf, double aschk=5e-2) {
0042     if (aschk >= 0) _checkAlphasQ2(Q2, basepdf, newpdf, aschk);
0043     const double xf_base = basepdf.xfxQ2(id, x, Q2);
0044     const double xf_new = newpdf.xfxQ2(id, x, Q2);
0045     return xf_new / xf_base;
0046   }
0047
0048   /// Get the PDF reweighting factor for a single beam with id,x,Q parameters, from basepdf to newpdf
0049   ///
0050   /// @note For NLO calculations, in general different PDF values enter for each counterterm: be careful.
0051   template <typename PDFPTR>
0052   inline double weightxQ2(int id, double x, double Q2, const PDFPTR basepdf, const PDFPTR newpdf, double aschk=5e-2) {
0053     return weightxQ2(id, x, Q2, *basepdf, *newpdf, aschk);
0054   }
0055
0056   /// Get the PDF reweighting factor for a single beam with id,x,Q parameters, from basepdf to newpdf
0057   ///
0058   /// @note For NLO calculations, in general different PDF values enter for each counterterm: be careful.
0059   inline double weightxQ(int id, double x, double Q, const PDF& basepdf, const PDF& newpdf, double aschk=5e-2) {
0060     return weightxQ2(id, x, sqr(Q), basepdf, newpdf, aschk);
0061   }
0062
0063   /// Get the PDF reweighting factor for a single beam with id,x,Q parameters, from basepdf to newpdf
0064   ///
0065   /// @note For NLO calculations, in general different PDF values enter for each counterterm: be careful.
0066   template <typename PDFPTR>
0067   inline double weightxQ(int id, double x, double Q, const PDFPTR basepdf, const PDFPTR newpdf, double aschk=5e-2) {
0068     return weightxQ(id, x, Q, *basepdf, *newpdf, aschk);
0069   }
0070
0071   /// @}
0072
0073
0074   /// @defgroup reweight_double Two-beam reweighting
0075   /// @{
0076
0077   /// Get the PDF reweighting factor for two beams, one with id1,x1 and the other with id2,x2, from basepdf to newpdf
0078   ///
0079   /// @note For NLO calculations, in general different PDF values enter for each counterterm: be careful.
0080   inline double weightxxQ2(int id1, int id2, double x1, double x2, double Q2, const PDF& basepdf, const PDF& newpdf, double aschk=5e-2) {
0081     if (aschk >= 0) _checkAlphasQ2(Q2, basepdf, newpdf, aschk);
0082     const double w1 = weightxQ2(id1, x1, Q2, basepdf, newpdf, -1);
0083     const double w2 = weightxQ2(id2, x2, Q2, basepdf, newpdf, -1);
0084     return w1 * w2;
0085   }
0086
0087   /// Get the PDF reweighting factor for two beams, one with id1,x1 and the other with id2,x2, from basepdf to newpdf
0088   ///
0089   /// @note For NLO calculations, in general different PDF values enter for each counterterm: be careful.
0090   template <typename PDFPTR>
0091   inline double weightxxQ2(int id1, int id2, double x1, double x2, double Q2, const PDFPTR basepdf, const PDFPTR newpdf, double aschk=5e-2) {
0092     return weightxxQ2(id1, id2, x1, x2, Q2, *basepdf, *newpdf, aschk);
0093   }
0094
0095   /// Get the PDF reweighting factor for two beams, one with id1,x1 and the other with id2,x2, from basepdf to newpdf
0096   ///
0097   /// @note For NLO calculations, in general different PDF values enter for each counterterm: be careful.
0098   inline double weightxxQ(int id1, int id2, double x1, double x2, double Q, const PDF& basepdf, const PDF& newpdf, double aschk=5e-2) {
0099     return weightxxQ2(id1, id2, x1, x2, sqr(Q), basepdf, newpdf, aschk);
0100   }
0101
0102   /// Get the PDF reweighting factor for two beams, one with id1,x1 and the other with id2,x2, from basepdf to newpdf
0103   ///
0104   /// @note For NLO calculations, in general different PDF values enter for each counterterm: be careful.
0105   template <typename PDFPTR>
0106   inline double weightxxQ(int id1, int id2, double x1, double x2, double Q, const PDFPTR basepdf, const PDFPTR newpdf, double aschk=5e-2) {
0107     return weightxxQ(id1, id2, x1, x2, Q, *basepdf, *newpdf, aschk);
0108   }
0109
0110   /// @}
0111
0112   /// @}
0113
0114 }
0115 #endif