File indexing completed on 2024-04-06 12:24:18
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 }
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