File indexing completed on 2023-10-25 10:00:57
0001 #include "RecoMET/METAlgorithms/interface/significanceAlgo.h"
0002 #include <iostream>
0003 #include <string>
0004 #include "TMath.h"
0005 #include <cmath>
0006 #include "TROOT.h"
0007 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0008
0009
0010
0011
0012
0013
0014
0015
0016
0017
0018
0019
0020
0021
0022
0023
0024
0025
0026
0027 metsig::significanceAlgo::significanceAlgo()
0028 :
0029 set_worker_(0),
0030 xmet_(0),
0031 ymet_(0)
0032
0033 {
0034
0035 signifmatrix_(0, 0) = signifmatrix_(1, 0) = signifmatrix_(0, 1) = signifmatrix_(1, 1) = 0;
0036 }
0037
0038
0039 const void metsig::significanceAlgo::addSignifMatrix(const reco::METCovMatrix &input) {
0040 signifmatrix_ += input;
0041 return;
0042 }
0043
0044
0045
0046 const void metsig::significanceAlgo::setSignifMatrix(const reco::METCovMatrix &input,
0047 const double &met_r,
0048 const double &met_phi,
0049 const double &met_set) {
0050 signifmatrix_ = input;
0051 set_worker_ = met_set;
0052 xmet_ = met_r * cos(met_phi);
0053 ymet_ = met_r * sin(met_phi);
0054
0055 return;
0056 }
0057
0058
0059
0060 metsig::significanceAlgo::~significanceAlgo() {}
0061
0062
0063
0064
0065
0066
0067
0068
0069
0070
0071
0072
0073
0074
0075
0076
0077
0078
0079 const void metsig::significanceAlgo::subtractObjects(const std::vector<metsig::SigInputObj> &eventVec) {
0080 reco::METCovMatrix v_tot = signifmatrix_;
0081
0082
0083 for (std::vector<SigInputObj>::const_iterator obj = eventVec.begin(); obj != eventVec.end(); ++obj) {
0084 double et_tmp = obj->get_energy();
0085 double phi_tmp = obj->get_phi();
0086 double sigma_et = obj->get_sigma_e();
0087 double sigma_tan = obj->get_sigma_tan();
0088
0089 double cosphi = cos(phi_tmp);
0090 double sinphi = sin(phi_tmp);
0091 double xval = et_tmp * cosphi;
0092 double yval = et_tmp * sinphi;
0093 xmet_ += xval;
0094 ymet_ += yval;
0095 set_worker_ -= et_tmp;
0096
0097 double sigma0_2 = sigma_et * sigma_et;
0098 double sigma1_2 = sigma_tan * sigma_tan;
0099
0100 v_tot(0, 0) -= sigma0_2 * cosphi * cosphi + sigma1_2 * sinphi * sinphi;
0101 v_tot(0, 1) -= cosphi * sinphi * (sigma0_2 - sigma1_2);
0102 v_tot(1, 0) -= cosphi * sinphi * (sigma0_2 - sigma1_2);
0103 v_tot(1, 1) -= sigma1_2 * cosphi * cosphi + sigma0_2 * sinphi * sinphi;
0104 }
0105 signifmatrix_ = v_tot;
0106 }
0107
0108
0109 const void metsig::significanceAlgo::addObjects(const std::vector<metsig::SigInputObj> &eventVec) {
0110 reco::METCovMatrix v_tot = signifmatrix_;
0111
0112
0113 for (std::vector<SigInputObj>::const_iterator obj = eventVec.begin(); obj != eventVec.end(); ++obj) {
0114 double et_tmp = obj->get_energy();
0115 double phi_tmp = obj->get_phi();
0116 double sigma_et = obj->get_sigma_e();
0117 double sigma_tan = obj->get_sigma_tan();
0118
0119 double cosphi = cos(phi_tmp);
0120 double sinphi = sin(phi_tmp);
0121 double xval = et_tmp * cosphi;
0122 double yval = et_tmp * sinphi;
0123 xmet_ -= xval;
0124 ymet_ -= yval;
0125 set_worker_ += et_tmp;
0126
0127 double sigma0_2 = sigma_et * sigma_et;
0128 double sigma1_2 = sigma_tan * sigma_tan;
0129
0130 v_tot(0, 0) += sigma0_2 * cosphi * cosphi + sigma1_2 * sinphi * sinphi;
0131 v_tot(0, 1) += cosphi * sinphi * (sigma0_2 - sigma1_2);
0132 v_tot(1, 0) += cosphi * sinphi * (sigma0_2 - sigma1_2);
0133 v_tot(1, 1) += sigma1_2 * cosphi * cosphi + sigma0_2 * sinphi * sinphi;
0134 }
0135 signifmatrix_ = v_tot;
0136 }
0137
0138
0139 const double metsig::significanceAlgo::significance(double &met_r, double &met_phi, double &met_set) {
0140 if (signifmatrix_(0, 0) == 0 && signifmatrix_(1, 1) == 0 && signifmatrix_(1, 0) == 0 && signifmatrix_(0, 1) == 0) {
0141
0142 return (-1);
0143 }
0144
0145
0146
0147 reco::METCovMatrix v_tot;
0148 ROOT::Math::SVector<double, 2> metvec;
0149
0150
0151 v_tot = signifmatrix_;
0152
0153
0154
0155 met_r = sqrt(xmet_ * xmet_ + ymet_ * ymet_);
0156 met_set = set_worker_;
0157 met_phi = TMath::ATan2(ymet_, xmet_);
0158
0159
0160
0161 double det = 0;
0162 v_tot.Det(det);
0163 if (fabs(det) < 0.000001)
0164 return -1;
0165
0166
0167
0168 v_tot.Invert();
0169
0170
0171 metvec(0) = xmet_;
0172 metvec(1) = ymet_;
0173 double lnSignificance = ROOT::Math::Dot(metvec, (v_tot * metvec));
0174
0175
0176
0177 return lnSignificance;
0178 }
0179