Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-02-14 13:33:27

0001 #ifndef PhysicsTools_Utilities_Integral_h
0002 #define PhysicsTools_Utilities_Integral_h
0003 #include "PhysicsTools/Utilities/interface/Primitive.h"
0004 #include "PhysicsTools/Utilities/interface/NumericalIntegration.h"
0005 
0006 namespace funct {
0007 
0008   struct no_var;
0009 
0010   template <typename F, typename X = no_var>
0011   struct IntegralStruct {
0012     IntegralStruct(const F& f) : p(primitive<X>(f)) {}
0013     double operator()(double min, double max) const {
0014       X::set(min);
0015       double pMin = p();
0016       X::set(max);
0017       double pMax = p();
0018       return pMax - pMin;
0019     }
0020 
0021   private:
0022     typename Primitive<F, X>::type p;
0023   };
0024 
0025   template <typename F>
0026   struct IntegralStruct<F> {
0027     IntegralStruct(const F& f) : p(primitive(f)) {}
0028     double operator()(double min, double max) const { return p(max) - p(min); }
0029 
0030   private:
0031     typename Primitive<F>::type p;
0032   };
0033 
0034   template <typename Integrator, typename F, typename X = no_var>
0035   struct NumericalIntegral {
0036     NumericalIntegral(const F& f, const Integrator& integrator) : f_(f), integrator_(integrator) {}
0037     inline double operator()(double min, double max) const { return integrator_(f_, min, max); }
0038 
0039   private:
0040     struct function {
0041       function(const F& f) : f_(f) {}
0042       inline double operator()(double x) const {
0043         X::set(x);
0044         return f_();
0045       }
0046 
0047     private:
0048       F f_;
0049     };
0050     function f_;
0051     Integrator integrator_;
0052   };
0053 
0054   template <typename Integrator, typename F>
0055   struct NumericalIntegral<Integrator, F, no_var> {
0056     NumericalIntegral(const F& f, const Integrator& integrator) : f_(f), integrator_(integrator) {}
0057     double operator()(double min, double max) const { return integrator_(f_, min, max); }
0058 
0059   private:
0060     F f_;
0061     Integrator integrator_;
0062   };
0063 
0064   template <typename F, typename X = no_var>
0065   struct Integral {
0066     typedef IntegralStruct<F, X> type;
0067   };
0068 
0069   template <typename X, typename F>
0070   typename Integral<F, X>::type integral(const F& f) {
0071     return typename Integral<F, X>::type(f);
0072   }
0073 
0074   template <typename X, typename F, typename Integrator>
0075   typename Integral<F, X>::type integral(const F& f, const Integrator& integrator) {
0076     return typename Integral<F, X>::type(f, integrator);
0077   }
0078 
0079   template <typename F, typename Integrator>
0080   typename Integral<F>::type integral_f(const F& f, const Integrator& integrator) {
0081     return typename Integral<F>::type(f, integrator);
0082   }
0083 
0084   template <typename X, typename F>
0085   double integral(const F& f, double min, double max) {
0086     return integral<X>(f)(min, max);
0087   }
0088 
0089   template <typename X, typename F, typename Integrator>
0090   double integral(const F& f, double min, double max, const Integrator& integrator) {
0091     return integral<X>(f, integrator)(min, max);
0092   }
0093 
0094   template <typename F>
0095   typename Integral<F>::type integral_f(const F& f) {
0096     return typename Integral<F>::type(f);
0097   }
0098 
0099   template <typename F>
0100   double integral_f(const F& f, double min, double max) {
0101     return integral_f(f)(min, max);
0102   }
0103 
0104   template <typename F, typename Integrator>
0105   double integral_f(const F& f, double min, double max, const Integrator& integrator) {
0106     return integral_f(f, integrator)(min, max);
0107   }
0108 
0109   template <typename F, typename MIN, typename MAX, typename Integrator = no_var, typename X = no_var>
0110   struct DefIntegral {
0111     DefIntegral(const F& f, const MIN& min, const MAX& max, const Integrator& integrator)
0112         : f_(f), min_(min), max_(max), integrator_(integrator) {}
0113     double operator()() const { return integral<X>(f_, min_(), max_(), integrator_); }
0114 
0115   private:
0116     F f_;
0117     MIN min_;
0118     MAX max_;
0119     Integrator integrator_;
0120   };
0121 
0122   template <typename F, typename MIN, typename MAX, typename X>
0123   struct DefIntegral<F, MIN, MAX, no_var, X> {
0124     DefIntegral(const F& f, const MIN& min, const MAX& max) : f_(f), min_(min), max_(max) {}
0125     double operator()(double x) const { return integral<X>(f_, min_(x), max_(x)); }
0126 
0127   private:
0128     F f_;
0129     MIN min_;
0130     MAX max_;
0131   };
0132 
0133   template <typename F, typename MIN, typename MAX, typename Integrator>
0134   struct DefIntegral<F, MIN, MAX, Integrator, no_var> {
0135     DefIntegral(const F& f, const MIN& min, const MAX& max, const Integrator& integrator)
0136         : f_(f), min_(min), max_(max), integrator_(integrator) {}
0137     double operator()(double x) const { return integral_f(f_, min_(x), max_(x), integrator_); }
0138 
0139   private:
0140     F f_;
0141     MIN min_;
0142     MAX max_;
0143     Integrator integrator_;
0144   };
0145 
0146   template <typename F, typename MIN, typename MAX>
0147   struct DefIntegral<F, MIN, MAX, no_var, no_var> {
0148     DefIntegral(const F& f, const MIN& min, const MAX& max) : f_(f), min_(min), max_(max) {}
0149     double operator()(double x) const { return integral_f(f_, min_(x), max_(x)); }
0150 
0151   private:
0152     F f_;
0153     MIN min_;
0154     MAX max_;
0155   };
0156 
0157 }  // namespace funct
0158 
0159 #define NUMERICAL_INTEGRAL(X, F, INTEGRATOR)            \
0160   namespace funct {                                     \
0161     template <typename X>                               \
0162     struct Integral<F, X> {                             \
0163       typedef NumericalIntegral<INTEGRATOR, F, X> type; \
0164     };                                                  \
0165   }                                                     \
0166   struct __useless_ignoreme
0167 
0168 #define NUMERICAL_FUNCT_INTEGRAL(F, INTEGRATOR)      \
0169   namespace funct {                                  \
0170     template <>                                      \
0171     struct Integral<F, no_var> {                     \
0172       typedef NumericalIntegral<INTEGRATOR, F> type; \
0173     };                                               \
0174   }                                                  \
0175   struct __useless_ignoreme
0176 
0177 #endif