1 #ifndef PhysicsTools_Heppy_MuScleFitCorrector_h 2 #define PhysicsTools_Heppy_MuScleFitCorrector_h 16 #include "TLorentzVector.h" 55 double pt = lorentzVector.Pt();
65 return (sigmaPtData*sigmaPtData - sigmaPtMC*sigmaPtMC);
73 double pt = lorentzVector.Pt();
74 double eta = lorentzVector.Eta();
76 if (squaredDiff < 0)
return pt;
78 if (
fileName_.Contains(
"smearPrompt")) Cfact = 0.85;
79 else if (
fileName_.Contains(
"smearReReco")) Cfact = 0.75;
80 else if (
fileName_.Contains(
"2011_MC")) Cfact = 0.8;
82 std::cout<<
"Are you sure you want to smear data??"<<std::endl;
85 double sPar = Cfact*
sqrt(squaredDiff);
86 double curv = ((double)chg/pt);
87 double normSmearFact =
gRandom_->Gaus(0,sPar);
88 if (fake) normSmearFact = sPar;
89 double smearedCurv = curv + fabs(curv)*normSmearFact;
91 return 1./((double)chg*smearedCurv);
98 double eta = lorentzVector.Eta();
99 double phi = lorentzVector.Phi();
100 double m = lorentzVector.M();
102 lorentzVector.SetPtEtaPhiM(corrPt,eta,phi,m);
108 double eta = lorentzVector.Eta();
109 double phi = lorentzVector.Phi();
110 double m = lorentzVector.M();
111 double smearedPt =
getSmearedPt(lorentzVector, chg, fake);
112 lorentzVector.SetPtEtaPhiM(smearedPt,eta,phi,m);
143 int resolDataParVecSize = resolDataParVec_.size();
144 int resolMCParVecSize = resolMCParVec_.size();
145 int scaleParVecSize = scaleParVec_.size();
148 if (resolDataParVecSize!=0 && resolMCParVecSize!=0)
useResol_ =
true;
150 if( resolParNum != resolDataParVecSize &&
useResol_) {
151 std::cout <<
"Error: inconsistent number of parameters: resolutionFunction has "<<resolParNum<<
" parameters but "<<resolDataParVecSize<<
" have been read from file" << std::endl;
155 if( resolParNum != resolMCParVecSize &&
useResol_) {
156 std::cout <<
"Error: inconsistent number of parameters: resolutionFunction has "<<resolParNum<<
" parameters but "<<resolMCParVecSize<<
" have been read from file" << std::endl;
160 if( scaleParNum != scaleParVecSize ) {
161 std::cout <<
"Error: inconsistent number of parameters: scaleFunction has "<<scaleParNum<<
" parameters but "<<scaleParVecSize<<
" have been read from file" << std::endl;
166 std::vector<double>::const_iterator scaleParIt = scaleParVec_.begin();
167 std::vector<double>::const_iterator resolDataParIt = resolDataParVec_.begin();
168 std::vector<double>::const_iterator resolMCParIt = resolMCParVec_.begin();
170 scaleParArray_ =
new double[scaleParNum];
171 resolDataParArray_ =
new double[resolParNum];
172 resolMCParArray_ =
new double[resolParNum];
174 for (
int iPar=0; iPar<scaleParNum; ++iPar) {
181 for (
int iPar=0; iPar<resolParNum; ++iPar) {
187 for (
int iPar=0; iPar<resolParNum; ++iPar) {
201 resolDataParArray_ = 0;
202 resolMCParArray_ = 0;
205 std::ifstream parametersFile(fileName.Data());
207 if( !parametersFile.is_open() ) {
208 std::cout <<
"Error: file " << fileName <<
" not found. Aborting." << std::endl;
216 while (parametersFile) {
218 getline( parametersFile, line );
219 size_t lineInt = line.find(
"value");
220 size_t iterationSubStr = line.find(iteration);
222 if( iterationSubStr != std::string::npos ) {
224 std::stringstream sLine(line);
226 int scaleFunctionNum = 0;
227 int resolutionFunctionNum = 0;
231 while( sLine >> num ) {
233 if( wordCounter == 8 ) {
234 std::stringstream
in(num);
235 in >> resolutionFunctionNum;
238 if( wordCounter == 9 ) {
239 std::stringstream
in(num);
240 in >> scaleFunctionNum;
244 scaleFunctionId_ = scaleFunctionNum;
246 resolutionFunctionId_ = resolutionFunctionNum;
258 if ( (lineInt != std::string::npos) ) {
259 size_t subStr1 = line.find(
"value");
260 std::stringstream paramStr;
262 paramStr << line.substr(subStr1+5);
265 if (nReadPar<nScalePar) {
266 scaleParVec_.push_back(param);
268 else if (nReadPar < (nScalePar+nResolPar) ) {
269 resolDataParVec_.push_back(param);
271 else if (nReadPar < (nScalePar+2*nResolPar) ) {
272 resolMCParVec_.push_back(param);
void readParameters(const TString &fileName)
scaleFunctionBase< double * > * scaleFunctionService(const int identifier)
Service to build the scale functor corresponding to the passed identifier.
MuScleFitCorrector(const TString &identifier)
std::vector< double > resolDataParVec_
virtual double scale(const double &pt, const double &eta, const double &phi, const int chg, const T &parScale) const =0
double getSmearedPt(const TLorentzVector &lorentzVector, const int &chg, bool fake)
reco::Particle::LorentzVector lorentzVector
std::vector< double > scaleParVec_
virtual double sigmaPt(const double &pt, const double &eta, const T &parval)=0
std::vector< double > resolMCParVec_
int resolutionFunctionId_
resolutionFunctionBase< double * > * getResolFunction()
scaleFunctionBase< double * > * scaleFunction_
resolutionFunctionBase< double * > * resolutionFunction_
virtual int parNum() const
std::vector< double > getResolMCParVec()
TAKEN FROM http://cmssw.cvs.cern.ch/cgi-bin/cmssw.cgi/CMSSW/ElectroWeakAnalysis/Utilities/src/PdfWeig...
double * resolMCParArray_
void applyPtSmearing(TLorentzVector &lorentzVector, const int &chg, bool fake=false)
double * resolDataParArray_
double getSigmaPtDiffSquared(const double &pt, const double &eta)
resolutionFunctionBase< double * > * resolutionFunctionService(const int identifier)
Service to build the resolution functor corresponding to the passed identifier.
virtual int parNum() const
double getCorrectPt(const TLorentzVector &lorentzVector, const int &chg)
void applyPtCorrection(TLorentzVector &lorentzVector, const int &chg)
scaleFunctionBase< double * > * getScaleFunction()