File indexing completed on 2023-03-17 11:20:55
0001 #include "RecoMuon/TrackerSeedGenerator/interface/SeedMvaEstimator.h"
0002
0003 #include "DataFormats/Math/interface/deltaPhi.h"
0004 #include "DataFormats/Math/interface/deltaR.h"
0005
0006 #include "FWCore/ParameterSet/interface/FileInPath.h"
0007 #include "CommonTools/MVAUtils/interface/GBRForestTools.h"
0008
0009 SeedMvaEstimator::SeedMvaEstimator(const edm::FileInPath& weightsfile,
0010 const std::vector<double>& scale_mean,
0011 const std::vector<double>& scale_std,
0012 const bool isFromL1,
0013 const int minL1Qual)
0014 : scale_mean_(scale_mean), scale_std_(scale_std), isFromL1_(isFromL1), minL1Qual_(minL1Qual) {
0015 gbrForest_ = createGBRForest(weightsfile);
0016 }
0017
0018 SeedMvaEstimator::~SeedMvaEstimator() {}
0019
0020 namespace {
0021 enum inputIndexes {
0022 kTsosErr0,
0023 kTsosErr2,
0024 kTsosErr5,
0025 kTsosDxdz,
0026 kTsosDydz,
0027 kTsosQbp,
0028 kDRdRL1SeedP,
0029 kDPhidRL1SeedP,
0030 kLastL1,
0031
0032 kDRdRL2SeedP = 8,
0033 kDPhidRL2SeedP,
0034 kLastL2,
0035 };
0036 }
0037
0038 void SeedMvaEstimator::getL1MuonVariables(const GlobalVector& global_p,
0039 const l1t::MuonBxCollection& l1Muons,
0040 float& dR2dRL1SeedP,
0041 float& dPhidRL1SeedP) const {
0042 for (int ibx = l1Muons.getFirstBX(); ibx <= l1Muons.getLastBX(); ++ibx) {
0043 if (ibx != 0)
0044 continue;
0045
0046 for (auto it = l1Muons.begin(ibx); it != l1Muons.end(ibx); it++) {
0047 if (it->hwQual() < minL1Qual_)
0048 continue;
0049
0050 float dR2tmp = reco::deltaR2(it->etaAtVtx(), it->phiAtVtx(), global_p.eta(), global_p.phi());
0051 if (dR2tmp < dR2dRL1SeedP) {
0052 dR2dRL1SeedP = dR2tmp;
0053 dPhidRL1SeedP = reco::deltaPhi(it->phiAtVtx(), global_p.phi());
0054 }
0055 }
0056 }
0057 }
0058
0059 void SeedMvaEstimator::getL2MuonVariables(const GlobalVector& global_p,
0060 const reco::RecoChargedCandidateCollection& l2Muons,
0061 float& dR2dRL2SeedP,
0062 float& dPhidRL2SeedP) const {
0063 for (auto it = l2Muons.begin(); it != l2Muons.end(); it++) {
0064 float dR2tmp = reco::deltaR2(*it, global_p);
0065 if (dR2tmp < dR2dRL2SeedP) {
0066 dR2dRL2SeedP = dR2tmp;
0067 dPhidRL2SeedP = reco::deltaPhi(it->phi(), global_p.phi());
0068 }
0069 }
0070 }
0071
0072 double SeedMvaEstimator::computeMva(const TrajectorySeed& seed,
0073 const GlobalVector& global_p,
0074 const l1t::MuonBxCollection& l1Muons,
0075 const reco::RecoChargedCandidateCollection& l2Muons) const {
0076 static constexpr float initDRdPhi(99999.);
0077 auto kLast = isFromL1_ ? kLastL1 : kLastL2;
0078 float var[kLast];
0079
0080 var[kTsosErr0] = seed.startingState().error(0);
0081 var[kTsosErr2] = seed.startingState().error(2);
0082 var[kTsosErr5] = seed.startingState().error(5);
0083 var[kTsosDxdz] = seed.startingState().parameters().dxdz();
0084 var[kTsosDydz] = seed.startingState().parameters().dydz();
0085 var[kTsosQbp] = seed.startingState().parameters().qbp();
0086
0087 float dR2dRL1SeedP = initDRdPhi;
0088 float dPhidRL1SeedP = initDRdPhi;
0089 getL1MuonVariables(global_p, l1Muons, dR2dRL1SeedP, dPhidRL1SeedP);
0090
0091 var[kDRdRL1SeedP] = std::sqrt(dR2dRL1SeedP);
0092 var[kDPhidRL1SeedP] = dPhidRL1SeedP;
0093
0094 if (!isFromL1_) {
0095 float dR2dRL2SeedP = initDRdPhi;
0096 float dPhidRL2SeedP = initDRdPhi;
0097 getL2MuonVariables(global_p, l2Muons, dR2dRL2SeedP, dPhidRL2SeedP);
0098
0099 var[kDRdRL2SeedP] = std::sqrt(dR2dRL2SeedP);
0100 var[kDPhidRL2SeedP] = dPhidRL2SeedP;
0101 }
0102
0103 for (int iv = 0; iv < kLast; ++iv) {
0104 var[iv] = (var[iv] - scale_mean_.at(iv)) / scale_std_.at(iv);
0105 }
0106
0107 return gbrForest_->GetResponse(var);
0108 }