File indexing completed on 2023-03-17 11:17:10
0001 #ifndef RecoBTag_SecondaryVertex_CombinedSVSoftLeptonComputer_h
0002 #define RecoBTag_SecondaryVertex_CombinedSVSoftLeptonComputer_h
0003
0004 #include "DataFormats/BTauReco/interface/CandSoftLeptonTagInfo.h"
0005 #include "DataFormats/BTauReco/interface/TaggingVariable.h"
0006 #include "DataFormats/BTauReco/interface/TemplatedSoftLeptonTagInfo.h"
0007 #include "DataFormats/BTauReco/interface/VertexTypes.h"
0008 #include "DataFormats/TrackReco/interface/TrackFwd.h"
0009
0010 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0011
0012 #include "RecoBTag/SecondaryVertex/interface/CombinedSVComputer.h"
0013
0014 class CombinedSVSoftLeptonComputer : public CombinedSVComputer {
0015 public:
0016 explicit CombinedSVSoftLeptonComputer(const edm::ParameterSet ¶ms);
0017 ~CombinedSVSoftLeptonComputer() override = default;
0018 double flipSoftLeptonValue(double value) const;
0019
0020 template <class IPTI, class SVTI>
0021 reco::TaggingVariableList operator()(const IPTI &ipInfo,
0022 const SVTI &svInfo,
0023 const reco::CandSoftLeptonTagInfo &muonInfo,
0024 const reco::CandSoftLeptonTagInfo &elecInfo) const;
0025
0026 private:
0027 bool SoftLeptonFlip;
0028 };
0029
0030 inline double CombinedSVSoftLeptonComputer::flipSoftLeptonValue(double value) const {
0031 return SoftLeptonFlip ? -value : value;
0032 }
0033
0034 template <class IPTI, class SVTI>
0035 reco::TaggingVariableList CombinedSVSoftLeptonComputer::operator()(const IPTI &ipInfo,
0036 const SVTI &svInfo,
0037 const reco::CandSoftLeptonTagInfo &muonInfo,
0038 const reco::CandSoftLeptonTagInfo &elecInfo) const {
0039 using namespace reco;
0040
0041
0042 TaggingVariableList vars = CombinedSVComputer::operator()(ipInfo, svInfo);
0043
0044
0045 unsigned int vtxType =
0046 (vars.checkTag(reco::btau::vertexCategory) ? (unsigned int)(vars.get(reco::btau::vertexCategory)) : 99);
0047 if (vtxType == 99)
0048 return vars;
0049
0050
0051 int leptonCategory = 0;
0052
0053 for (unsigned int i = 0; i < muonInfo.leptons();
0054 ++i)
0055 {
0056 leptonCategory = 1;
0057 const SoftLeptonProperties &propertiesMuon = muonInfo.properties(i);
0058 vars.insert(btau::leptonPtRel, propertiesMuon.ptRel, true);
0059 vars.insert(btau::leptonSip3d, flipSoftLeptonValue(propertiesMuon.sip3d), true);
0060 vars.insert(btau::leptonDeltaR, propertiesMuon.deltaR, true);
0061 vars.insert(btau::leptonRatioRel, propertiesMuon.ratioRel, true);
0062 vars.insert(btau::leptonEtaRel, propertiesMuon.etaRel, true);
0063 vars.insert(btau::leptonRatio, propertiesMuon.ratio, true);
0064 }
0065
0066 if (leptonCategory != 1)
0067 {
0068 for (unsigned int i = 0; i < elecInfo.leptons();
0069 ++i)
0070 {
0071 leptonCategory = 2;
0072 const SoftLeptonProperties &propertiesElec = elecInfo.properties(i);
0073 vars.insert(btau::leptonPtRel, propertiesElec.ptRel, true);
0074 vars.insert(btau::leptonSip3d, flipSoftLeptonValue(propertiesElec.sip3d), true);
0075 vars.insert(btau::leptonDeltaR, propertiesElec.deltaR, true);
0076 vars.insert(btau::leptonRatioRel, propertiesElec.ratioRel, true);
0077 vars.insert(btau::leptonEtaRel, propertiesElec.etaRel, true);
0078 vars.insert(btau::leptonRatio, propertiesElec.ratio, true);
0079 }
0080 }
0081
0082
0083 int vertexLepCat = 2;
0084
0085 if (leptonCategory == 0)
0086 {
0087 if (vtxType == (unsigned int)(btag::Vertices::RecoVertex))
0088 vertexLepCat = 0;
0089 else if (vtxType == (unsigned int)(btag::Vertices::PseudoVertex))
0090 vertexLepCat = 1;
0091 else
0092 vertexLepCat = 2;
0093 } else if (leptonCategory == 1)
0094 {
0095 if (vtxType == (unsigned int)(btag::Vertices::RecoVertex))
0096 vertexLepCat = 3;
0097 else if (vtxType == (unsigned int)(btag::Vertices::PseudoVertex))
0098 vertexLepCat = 4;
0099 else
0100 vertexLepCat = 5;
0101 } else if (leptonCategory == 2)
0102 {
0103 if (vtxType == (unsigned int)(btag::Vertices::RecoVertex))
0104 vertexLepCat = 6;
0105 else if (vtxType == (unsigned int)(btag::Vertices::PseudoVertex))
0106 vertexLepCat = 7;
0107 else
0108 vertexLepCat = 8;
0109 }
0110 vars.insert(btau::vertexLeptonCategory, vertexLepCat, true);
0111
0112 vars.finalize();
0113 return vars;
0114 }
0115
0116 #endif