Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2022-06-30 23:16:54

0001 //
0002 //
0003 
0004 #include "DataFormats/PatCandidates/interface/MET.h"
0005 
0006 using namespace pat;
0007 
0008 /// default constructor
0009 MET::MET() { initCorMap(); }
0010 
0011 /// constructor from reco::MET
0012 MET::MET(const reco::MET &aMET) : PATObject<reco::MET>(aMET) {
0013   const reco::CaloMET *calo = dynamic_cast<const reco::CaloMET *>(&aMET);
0014   if (calo != nullptr)
0015     caloMET_.push_back(calo->getSpecific());
0016   const reco::PFMET *pf = dynamic_cast<const reco::PFMET *>(&aMET);
0017   if (pf != nullptr)
0018     pfMET_.push_back(pf->getSpecific());
0019   const pat::MET *pm = dynamic_cast<const pat::MET *>(&aMET);
0020   if (pm != nullptr)
0021     this->operator=(*pm);
0022 
0023   metSig_ = 0.;
0024   sumPtUnclustered_ = 0.;
0025   initCorMap();
0026 }
0027 
0028 /// constructor from ref to reco::MET
0029 MET::MET(const edm::RefToBase<reco::MET> &aMETRef) : PATObject<reco::MET>(aMETRef) {
0030   const reco::CaloMET *calo = dynamic_cast<const reco::CaloMET *>(aMETRef.get());
0031   if (calo != nullptr)
0032     caloMET_.push_back(calo->getSpecific());
0033   const reco::PFMET *pf = dynamic_cast<const reco::PFMET *>(aMETRef.get());
0034   if (pf != nullptr)
0035     pfMET_.push_back(pf->getSpecific());
0036   const pat::MET *pm = dynamic_cast<const pat::MET *>(aMETRef.get());
0037   if (pm != nullptr)
0038     this->operator=(*pm);
0039 
0040   metSig_ = 0.;
0041   sumPtUnclustered_ = 0.;
0042   initCorMap();
0043 }
0044 
0045 /// constructor from ref to reco::MET
0046 MET::MET(const edm::Ptr<reco::MET> &aMETRef) : PATObject<reco::MET>(aMETRef) {
0047   const reco::CaloMET *calo = dynamic_cast<const reco::CaloMET *>(aMETRef.get());
0048   if (calo != nullptr)
0049     caloMET_.push_back(calo->getSpecific());
0050   const reco::PFMET *pf = dynamic_cast<const reco::PFMET *>(aMETRef.get());
0051   if (pf != nullptr)
0052     pfMET_.push_back(pf->getSpecific());
0053   const pat::MET *pm = dynamic_cast<const pat::MET *>(aMETRef.get());
0054   if (pm != nullptr)
0055     this->operator=(*pm);
0056 
0057   metSig_ = 0.;
0058   sumPtUnclustered_ = 0.;
0059   initCorMap();
0060 }
0061 
0062 /// copy constructor
0063 MET::MET(MET const &iOther)
0064     : PATObject<reco::MET>(iOther),
0065       genMET_(iOther.genMET_),
0066       caloMET_(iOther.caloMET_),
0067       pfMET_(iOther.pfMET_),
0068       metSig_(iOther.metSig_),
0069       sumPtUnclustered_(iOther.sumPtUnclustered_),
0070       uncertaintiesRaw_(iOther.uncertaintiesRaw_),          //74X reading compatibility
0071       uncertaintiesType1_(iOther.uncertaintiesType1_),      //74X compatibility
0072       uncertaintiesType1p2_(iOther.uncertaintiesType1p2_),  //74X compatibility
0073       uncertainties_(iOther.uncertainties_),
0074       corrections_(iOther.corrections_),
0075       caloPackedMet_(iOther.caloPackedMet_) {
0076   initCorMap();
0077 }
0078 
0079 /// constructor for corrected mets, keeping track of srcMET informations,
0080 // old uncertainties discarded on purpose to avoid confusion
0081 MET::MET(const reco::MET &corMET, const MET &srcMET)
0082     : PATObject<reco::MET>(corMET),
0083       genMET_(srcMET.genMET_),
0084       caloMET_(srcMET.caloMET_),
0085       pfMET_(srcMET.pfMET_),
0086       metSig_(srcMET.metSig_),
0087       sumPtUnclustered_(srcMET.sumPtUnclustered_),
0088       caloPackedMet_(srcMET.caloPackedMet_) {
0089   setSignificanceMatrix(srcMET.getSignificanceMatrix());
0090 
0091   initCorMap();
0092 }
0093 
0094 /// destructor
0095 MET::~MET() {}
0096 
0097 MET &MET::operator=(MET const &iOther) {
0098   PATObject<reco::MET>::operator=(iOther);
0099   genMET_ = iOther.genMET_;
0100   caloMET_ = iOther.caloMET_;
0101   pfMET_ = iOther.pfMET_;
0102   uncertaintiesRaw_ = iOther.uncertaintiesRaw_;  //74X compatibility
0103   uncertaintiesType1_ = iOther.uncertaintiesType1_;
0104   uncertaintiesType1p2_ = iOther.uncertaintiesType1p2_;
0105   uncertainties_ = iOther.uncertainties_;
0106   corrections_ = iOther.corrections_;
0107   metSig_ = iOther.metSig_;
0108   sumPtUnclustered_ = iOther.sumPtUnclustered_;
0109   caloPackedMet_ = iOther.caloPackedMet_;
0110 
0111   return *this;
0112 }
0113 
0114 /// return the generated MET from neutrinos
0115 const reco::GenMET *MET::genMET() const { return (!genMET_.empty() ? &genMET_.front() : nullptr); }
0116 
0117 /// method to set the generated MET
0118 void MET::setGenMET(const reco::GenMET &gm) {
0119   genMET_.clear();
0120   genMET_.push_back(gm);
0121 }
0122 
0123 //Method to set the MET significance
0124 void MET::setMETSignificance(const double &metSig) { metSig_ = metSig; }
0125 
0126 double MET::metSignificance() const { return metSig_; }
0127 
0128 void MET::setMETSumPtUnclustered(const double &sumPtUnclustered) { sumPtUnclustered_ = sumPtUnclustered; }
0129 
0130 double MET::metSumPtUnclustered() const { return sumPtUnclustered_; }
0131 
0132 void MET::initCorMap() {
0133   std::vector<MET::METCorrectionType> tmpRaw;
0134   std::vector<MET::METCorrectionType> tmpType1;
0135   std::vector<MET::METCorrectionType> tmpType01;
0136   std::vector<MET::METCorrectionType> tmpTypeXY;
0137   std::vector<MET::METCorrectionType> tmpType1XY;
0138   std::vector<MET::METCorrectionType> tmpType01XY;
0139   std::vector<MET::METCorrectionType> tmpType1Smear;
0140   std::vector<MET::METCorrectionType> tmpType01Smear;
0141   std::vector<MET::METCorrectionType> tmpType1SmearXY;
0142   std::vector<MET::METCorrectionType> tmpType01SmearXY;
0143 
0144   tmpRaw.push_back(MET::None);
0145 
0146   tmpType1.push_back(MET::T1);
0147   tmpType01.push_back(MET::T1);
0148   tmpType1XY.push_back(MET::T1);
0149   tmpType01XY.push_back(MET::T1);
0150   tmpType1Smear.push_back(MET::T1);
0151   tmpType01Smear.push_back(MET::T1);
0152   tmpType1SmearXY.push_back(MET::T1);
0153   tmpType01SmearXY.push_back(MET::T1);
0154 
0155   tmpType01.push_back(MET::T0);
0156   tmpType01XY.push_back(MET::T0);
0157   tmpType01Smear.push_back(MET::T0);
0158   tmpType01SmearXY.push_back(MET::T0);
0159 
0160   tmpType1Smear.push_back(MET::Smear);
0161   tmpType01Smear.push_back(MET::Smear);
0162   tmpType1SmearXY.push_back(MET::Smear);
0163   tmpType01SmearXY.push_back(MET::Smear);
0164 
0165   tmpTypeXY.push_back(MET::TXYForRaw);
0166   tmpType1XY.push_back(MET::TXY);
0167   tmpType01XY.push_back(MET::TXYForT01);
0168   tmpType1SmearXY.push_back(MET::TXYForT1Smear);
0169   tmpType01SmearXY.push_back(MET::TXYForT01Smear);
0170 
0171   corMap_[MET::Raw] = tmpRaw;
0172   corMap_[MET::Type1] = tmpType1;
0173   corMap_[MET::Type01] = tmpType01;
0174   corMap_[MET::TypeXY] = tmpTypeXY;
0175   corMap_[MET::Type1XY] = tmpType1XY;
0176   corMap_[MET::Type01XY] = tmpType01XY;
0177   corMap_[MET::Type1Smear] = tmpType1Smear;
0178   corMap_[MET::Type01Smear] = tmpType01Smear;
0179   corMap_[MET::Type1SmearXY] = tmpType1SmearXY;
0180   corMap_[MET::Type01SmearXY] = tmpType01SmearXY;
0181 
0182   //specific calo case
0183   std::vector<MET::METCorrectionType> tmpRawCalo;
0184   tmpRawCalo.push_back(MET::Calo);
0185   corMap_[MET::RawCalo] = tmpRawCalo;
0186 
0187   //specific chs case
0188   std::vector<MET::METCorrectionType> tmpRawChs;
0189   tmpRawChs.push_back(MET::Chs);
0190   corMap_[MET::RawChs] = tmpRawChs;
0191 
0192   //specific trk case
0193   std::vector<MET::METCorrectionType> tmpRawTrk;
0194   tmpRawTrk.push_back(MET::Trk);
0195   corMap_[MET::RawTrk] = tmpRawTrk;
0196 
0197   //specific deep response tune case
0198   std::vector<MET::METCorrectionType> tmpDeepResponse;
0199   tmpDeepResponse.push_back(MET::DeepResponseTune);
0200   corMap_[MET::RawDeepResponseTune] = tmpDeepResponse;
0201 
0202   //specific deep resolution tune case
0203   std::vector<MET::METCorrectionType> tmpDeepResolution;
0204   tmpDeepResolution.push_back(MET::DeepResolutionTune);
0205   corMap_[MET::RawDeepResolutionTune] = tmpDeepResolution;
0206 }
0207 
0208 MET::UnpackedMETUncertainty MET::findMETTotalShift(MET::METCorrectionLevel cor, MET::METUncertainty shift) const {
0209   //find corrections shifts =============================
0210   std::map<MET::METCorrectionLevel, std::vector<MET::METCorrectionType> >::const_iterator itCor_ = corMap_.find(cor);
0211   if (itCor_ == corMap_.end())
0212     throw cms::Exception("Unsupported", "Specified MET correction scheme does not exist");
0213 
0214   bool isSmeared = false;
0215   MET::UnpackedMETUncertainty totShift;
0216   unsigned int scor = itCor_->second.size();
0217   for (unsigned int i = 0; i < scor; i++) {
0218     auto up = corrections_[itCor_->second[i]].unpack();
0219     totShift.add(up.dpx(), up.dpy(), up.dsumEt());
0220 
0221     if (itCor_->first >= MET::Type1Smear)
0222       isSmeared = true;
0223   }
0224 
0225   //find uncertainty shift =============================
0226 
0227   if (uncertainties_.empty())
0228     return totShift;
0229 
0230   if (shift >= MET::METUncertaintySize)
0231     throw cms::Exception("Unsupported", "MET uncertainty does not exist");
0232   if (isSmeared && shift <= MET::JetResDown)
0233     shift = (MET::METUncertainty)(MET::METUncertaintySize + shift + 1);
0234 
0235   auto up = uncertainties_[shift].unpack();
0236   totShift.add(up.dpx(), up.dpy(), up.dsumEt());
0237 
0238   return totShift;
0239 }
0240 
0241 MET::Vector2 MET::shiftedP2(MET::METUncertainty shift, MET::METCorrectionLevel cor) const {
0242   Vector2 vo;
0243 
0244   //backward compatibility with 74X samples -> the only one
0245   // with uncertaintiesType1_/uncertaintiesRaw_ not empty
0246   //will be removed once 74X is not used anymore
0247   if (!uncertaintiesType1_.empty() || !uncertaintiesRaw_.empty()) {
0248     if (cor != MET::METCorrectionLevel::RawCalo) {
0249       vo = shiftedP2_74x(shift, cor);
0250     } else {
0251       Vector2 ret{caloPackedMet_.unpackDpx(), caloPackedMet_.unpackDpy()};
0252       vo = ret;
0253     }
0254   } else {
0255     auto v = findMETTotalShift(cor, shift);
0256     Vector2 ret{(px() + v.dpx()), (py() + v.dpy())};
0257     //return ret;
0258     vo = ret;
0259   }
0260   return vo;
0261 }
0262 MET::Vector MET::shiftedP3(MET::METUncertainty shift, MET::METCorrectionLevel cor) const {
0263   Vector vo;
0264 
0265   //backward compatibility with 74X samples -> the only one
0266   // with uncertaintiesType1_/uncertaintiesRaw_ not empty
0267   //will be removed once 74X is not used anymore
0268   if (!uncertaintiesType1_.empty() || !uncertaintiesRaw_.empty()) {
0269     if (cor != MET::METCorrectionLevel::RawCalo) {
0270       vo = shiftedP3_74x(shift, cor);
0271     } else {
0272       Vector tmp(caloPackedMet_.unpackDpx(), caloPackedMet_.unpackDpy(), 0);
0273       vo = tmp;
0274     }
0275   } else {
0276     const MET::UnpackedMETUncertainty &v = findMETTotalShift(cor, shift);
0277     //return Vector(px() + v.dpx(), py() + v.dpy(), 0);
0278     Vector tmp(px() + v.dpx(), py() + v.dpy(), 0);
0279     vo = tmp;
0280   }
0281   return vo;
0282 }
0283 MET::LorentzVector MET::shiftedP4(METUncertainty shift, MET::METCorrectionLevel cor) const {
0284   LorentzVector vo;
0285 
0286   //backward compatibility with 74X samples -> the only one
0287   // with uncertaintiesType1_/uncertaintiesRaw_ not empty
0288   //will be removed once 74X is not used anymore
0289   if (!uncertaintiesType1_.empty() || !uncertaintiesRaw_.empty()) {
0290     if (cor != MET::METCorrectionLevel::RawCalo) {
0291       vo = shiftedP4_74x(shift, cor);
0292     } else {
0293       double x = caloPackedMet_.unpackDpx(), y = caloPackedMet_.unpackDpy();
0294       LorentzVector tmp(x, y, 0, std::hypot(x, y));
0295       vo = tmp;
0296     }
0297   } else {
0298     const auto v = findMETTotalShift(cor, shift);
0299     double x = px() + v.dpx(), y = py() + v.dpy();
0300     //return LorentzVector(x, y, 0, std::hypot(x,y));
0301     LorentzVector tmp(x, y, 0, std::hypot(x, y));
0302     vo = tmp;
0303   }
0304   return vo;
0305 }
0306 double MET::shiftedSumEt(MET::METUncertainty shift, MET::METCorrectionLevel cor) const {
0307   double sumEto;
0308 
0309   //backward compatibility with 74X samples -> the only one
0310   // with uncertaintiesType1_/uncertaintiesRaw_ not empty
0311   //will be removed once 74X is not used anymore
0312   if (!uncertaintiesType1_.empty() || !uncertaintiesRaw_.empty()) {
0313     if (cor != MET::METCorrectionLevel::RawCalo) {
0314       sumEto = shiftedSumEt_74x(shift, cor);
0315     } else {
0316       sumEto = caloPackedMet_.unpackDSumEt();
0317     }
0318   } else {
0319     const auto v = findMETTotalShift(cor, shift);
0320     //return sumEt() + v.dsumEt();
0321     sumEto = sumEt() + v.dsumEt();
0322   }
0323   return sumEto;
0324 }
0325 
0326 MET::Vector2 MET::corP2(MET::METCorrectionLevel cor) const { return shiftedP2(MET::NoShift, cor); }
0327 MET::Vector MET::corP3(MET::METCorrectionLevel cor) const { return shiftedP3(MET::NoShift, cor); }
0328 MET::LorentzVector MET::corP4(MET::METCorrectionLevel cor) const { return shiftedP4(MET::NoShift, cor); }
0329 double MET::corSumEt(MET::METCorrectionLevel cor) const { return shiftedSumEt(MET::NoShift, cor); }
0330 
0331 MET::Vector2 MET::uncorP2() const { return shiftedP2(MET::NoShift, MET::Raw); }
0332 MET::Vector MET::uncorP3() const { return shiftedP3(MET::NoShift, MET::Raw); }
0333 MET::LorentzVector MET::uncorP4() const { return shiftedP4(MET::NoShift, MET::Raw); }
0334 double MET::uncorSumEt() const { return shiftedSumEt(MET::NoShift, MET::Raw); }
0335 
0336 void MET::setUncShift(double px, double py, double sumEt, METUncertainty shift, bool isSmeared) {
0337   if (uncertainties_.empty()) {
0338     uncertainties_.resize(METUncertainty::METFullUncertaintySize);
0339   }
0340 
0341   if (isSmeared && shift <= MET::JetResDown) {
0342     //changing reference to only get the uncertainty shift and not the smeared one
0343     // which is performed independently
0344     shift = (MET::METUncertainty)(METUncertainty::METUncertaintySize + shift + 1);
0345     const PackedMETUncertainty &ref = corrections_[METCorrectionType::Smear];
0346     uncertainties_[shift].set(px - ref.unpackDpx() - this->px(),
0347                               py - ref.unpackDpy() - this->py(),
0348                               sumEt - ref.unpackDSumEt() - this->sumEt());
0349   } else
0350     uncertainties_[shift].set(px - this->px(), py - this->py(), sumEt - this->sumEt());
0351 }
0352 
0353 void MET::setCorShift(double px, double py, double sumEt, MET::METCorrectionType level) {
0354   if (corrections_.empty()) {
0355     corrections_.resize(MET::METCorrectionType::METCorrectionTypeSize);
0356   }
0357 
0358   corrections_[level].set(px - this->px(), py - this->py(), sumEt - this->sumEt());
0359 }
0360 
0361 MET::Vector2 MET::caloMETP2() const {
0362   return shiftedP2(MET::METUncertainty::NoShift, MET::METCorrectionLevel::RawCalo);
0363 }
0364 
0365 double MET::caloMETPt() const { return caloMETP2().pt(); }
0366 
0367 double MET::caloMETPhi() const { return caloMETP2().phi(); }
0368 
0369 double MET::caloMETSumEt() const { return shiftedSumEt(MET::NoShift, MET::RawCalo); }
0370 
0371 // functions to access to 74X samples ========================================================
0372 MET::Vector2 MET::shiftedP2_74x(MET::METUncertainty shift, MET::METCorrectionLevel level) const {
0373   if (level != Type1 && level != Raw)
0374     throw cms::Exception("Unsupported", "MET uncertainties only supported for Raw and Type1 in 74X samples \n");
0375   const std::vector<PackedMETUncertainty> &v = (level == Type1 ? uncertaintiesType1_ : uncertaintiesRaw_);
0376   if (v.empty())
0377     throw cms::Exception("Unsupported", "MET uncertainties not available for the specified correction type\n");
0378   if (v.size() == 1) {
0379     if (shift != MET::METUncertainty::NoShift)
0380       throw cms::Exception(
0381           "Unsupported",
0382           "MET uncertainties not available for the specified correction type (only central value available)\n");
0383     auto const &p = v.front();
0384     return Vector2{(px() + p.unpackDpx()), (py() + p.unpackDpy())};
0385   }
0386   auto const &p = v[shift];
0387   Vector2 ret{(px() + p.unpackDpx()), (py() + p.unpackDpy())};
0388   return ret;
0389 }
0390 
0391 MET::Vector MET::shiftedP3_74x(MET::METUncertainty shift, MET::METCorrectionLevel level) const {
0392   if (level != Type1 && level != Raw)
0393     throw cms::Exception("Unsupported", "MET uncertainties only supported for Raw and Type1 in 74X samples \n");
0394   const std::vector<PackedMETUncertainty> &v = (level == Type1 ? uncertaintiesType1_ : uncertaintiesRaw_);
0395   if (v.empty())
0396     throw cms::Exception("Unsupported", "MET uncertainties not available for the specified correction type\n");
0397   if (v.size() == 1) {
0398     if (shift != MET::METUncertainty::NoShift)
0399       throw cms::Exception(
0400           "Unsupported",
0401           "MET uncertainties not available for the specified correction type (only central value available)\n");
0402     auto const &p = v.front();
0403     return Vector(px() + p.unpackDpx(), py() + p.unpackDpy(), 0);
0404   }
0405   auto const &p = v[shift];
0406   return Vector(px() + p.unpackDpx(), py() + p.unpackDpy(), 0);
0407 }
0408 
0409 MET::LorentzVector MET::shiftedP4_74x(METUncertainty shift, MET::METCorrectionLevel level) const {
0410   if (level != Type1 && level != Raw)
0411     throw cms::Exception("Unsupported", "MET uncertainties only supported for Raw and Type1 in 74X samples\n");
0412   const std::vector<PackedMETUncertainty> &v = (level == Type1 ? uncertaintiesType1_ : uncertaintiesRaw_);
0413   if (v.empty())
0414     throw cms::Exception("Unsupported", "MET uncertainties not available for the specified correction type\n");
0415   if (v.size() == 1) {
0416     if (shift != MET::METUncertainty::NoShift)
0417       throw cms::Exception(
0418           "Unsupported",
0419           "MET uncertainties not available for the specified correction type (only central value available)\n");
0420     auto const &p = v.front();
0421     double x = px() + p.unpackDpx(), y = py() + p.unpackDpy();
0422     return LorentzVector(x, y, 0, std::hypot(x, y));
0423   }
0424   auto const &p = v[shift];
0425   double x = px() + p.unpackDpx(), y = py() + p.unpackDpy();
0426   return LorentzVector(x, y, 0, std::hypot(x, y));
0427 }
0428 
0429 double MET::shiftedSumEt_74x(MET::METUncertainty shift, MET::METCorrectionLevel level) const {
0430   if (level != Type1 && level != Raw)
0431     throw cms::Exception("Unsupported", "MET uncertainties only supported for Raw and Type1 in 74X samples\n");
0432   const std::vector<PackedMETUncertainty> &v = (level == Type1 ? uncertaintiesType1_ : uncertaintiesRaw_);
0433   if (v.empty())
0434     throw cms::Exception("Unsupported", "MET uncertainties not available for the specified correction type\n");
0435   if (v.size() == 1) {
0436     if (shift != MET::METUncertainty::NoShift)
0437       throw cms::Exception(
0438           "Unsupported",
0439           "MET uncertainties not available for the specified correction type (only central value available)\n");
0440     return sumEt() + v.front().unpackDSumEt();
0441   }
0442   return sumEt() + v[shift].unpackDSumEt();
0443 }
0444 
0445 #include "DataFormats/Math/interface/libminifloat.h"
0446 
0447 MET::UnpackedMETUncertainty MET::PackedMETUncertainty::unpack() const {
0448   auto dpx = MiniFloatConverter::float16to32(packedDpx_);
0449   auto dpy = MiniFloatConverter::float16to32(packedDpy_);
0450   auto dsumEt = MiniFloatConverter::float16to32(packedDSumEt_);
0451   return UnpackedMETUncertainty(dpx, dpy, dsumEt);
0452 }
0453 
0454 float MET::PackedMETUncertainty::unpackDpx() const { return MiniFloatConverter::float16to32(packedDpx_); }
0455 
0456 float MET::PackedMETUncertainty::unpackDpy() const { return MiniFloatConverter::float16to32(packedDpy_); }
0457 
0458 float MET::PackedMETUncertainty::unpackDSumEt() const { return MiniFloatConverter::float16to32(packedDSumEt_); }
0459 
0460 void MET::PackedMETUncertainty::pack(float dpx, float dpy, float dsumEt) {
0461   packedDpx_ = MiniFloatConverter::float32to16(dpx);
0462   packedDpy_ = MiniFloatConverter::float32to16(dpy);
0463   packedDSumEt_ = MiniFloatConverter::float32to16(dsumEt);
0464 }