File indexing completed on 2023-03-17 11:15:53
0001 #include "PhysicsTools/IsolationUtils/interface/IntegrandThetaFunction.h"
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011
0012
0013
0014
0015
0016
0017
0018
0019
0020
0021
0022
0023 #include <iostream>
0024 #include <iomanip>
0025 #include <vector>
0026
0027
0028 #include <TMath.h>
0029
0030
0031 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0032
0033 #include "PhysicsTools/IsolationUtils/interface/IntegralOverPhiFunction.h"
0034 #include "DataFormats/Math/interface/normalizedPhi.h"
0035
0036
0037
0038
0039
0040 IntegrandThetaFunction::IntegrandThetaFunction()
0041 : ROOT::Math::ParamFunction<ROOT::Math::IParametricGradFunctionOneDim>(3) {
0042 theta0_ = 0.;
0043 phi0_ = 0.;
0044 alpha_ = 0.;
0045
0046 fPhi_ = new IntegralOverPhiFunction();
0047 }
0048
0049 IntegrandThetaFunction::IntegrandThetaFunction(const IntegrandThetaFunction& bluePrint)
0050 : ROOT::Math::ParamFunction<ROOT::Math::IParametricGradFunctionOneDim>(bluePrint) {
0051 theta0_ = bluePrint.theta0_;
0052 phi0_ = bluePrint.phi0_;
0053 alpha_ = bluePrint.alpha_;
0054
0055 fPhi_ = new IntegralOverPhiFunction(*bluePrint.fPhi_);
0056 }
0057
0058 IntegrandThetaFunction::~IntegrandThetaFunction() { delete fPhi_; }
0059
0060
0061
0062
0063
0064 IntegrandThetaFunction& IntegrandThetaFunction::operator=(const IntegrandThetaFunction& bluePrint) {
0065 theta0_ = bluePrint.theta0_;
0066 phi0_ = bluePrint.phi0_;
0067 alpha_ = bluePrint.alpha_;
0068
0069 (*fPhi_) = (*bluePrint.fPhi_);
0070
0071 return (*this);
0072 }
0073
0074
0075
0076
0077
0078 void IntegrandThetaFunction::SetParameterTheta0(double theta0) { theta0_ = theta0; }
0079
0080 void IntegrandThetaFunction::SetParameterPhi0(double phi0) {
0081 phi0_ = normalizedPhi(phi0);
0082 }
0083
0084 void IntegrandThetaFunction::SetParameterAlpha(double alpha) { alpha_ = alpha; }
0085
0086 void IntegrandThetaFunction::SetParameters(const double* param) {
0087 theta0_ = param[0];
0088 phi0_ = param[1];
0089 alpha_ = param[2];
0090 }
0091
0092 double IntegrandThetaFunction::DoEvalPar(double x, const double* param) const
0093 {
0094 theta0_ = param[0];
0095 phi0_ = param[1];
0096 alpha_ = param[2];
0097
0098 return DoEval(x);
0099 }
0100
0101 double IntegrandThetaFunction::DoEval(double x) const {
0102
0103
0104 const double epsilon = 1.e-3;
0105 if (x < epsilon || x > (TMath::Pi() - epsilon))
0106 return 0.;
0107
0108
0109
0110
0111 double sinTheta = TMath::Sin(x);
0112 double cscTheta = 1. / sinTheta;
0113
0114 double detJacobi =
0115 -cscTheta;
0116
0117
0118
0119
0120 fPhi_->SetParameterTheta0(theta0_);
0121 fPhi_->SetParameterPhi0(phi0_);
0122 fPhi_->SetParameterAlpha(alpha_);
0123
0124 double integralOverPhi = (*fPhi_)(x);
0125
0126 if (debugLevel_ > 0) {
0127 edm::LogVerbatim("") << "integralOverPhi = " << integralOverPhi << std::endl
0128 << " theta0 = " << theta0_ << std::endl
0129 << " phi0 = " << phi0_ << std::endl
0130 << " alpha = " << alpha_ << std::endl
0131 << " theta = " << x << std::endl
0132 << std::endl;
0133 }
0134
0135
0136
0137
0138
0139
0140
0141
0142
0143
0144
0145 return TMath::Abs(detJacobi) * integralOverPhi;
0146 }
0147
0148 double IntegrandThetaFunction::DoDerivative(double x) const {
0149
0150
0151 edm::LogWarning("") << "Function not implemented yet !" << std::endl;
0152
0153 return 0.;
0154 }
0155
0156 double IntegrandThetaFunction::DoParameterDerivative(double, const double*, unsigned int) const {
0157
0158
0159 edm::LogWarning("") << "Function not implemented yet !" << std::endl;
0160
0161 return 0.;
0162 }
0163
0164 void IntegrandThetaFunction::DoParameterGradient(double x, double* paramGradient) const {
0165
0166
0167 edm::LogWarning("") << "Function not implemented yet !" << std::endl;
0168 }