ImpactParameter/interface/TemplatedTrackCountingComputer.h

0001 #ifndef ImpactParameter_TemplatedTrackCountingComputer_h
0002 #define ImpactParameter_TemplatedTrackCountingComputer_h
0003
0004 #include "DataFormats/TrackReco/interface/Track.h"
0005 #include "DataFormats/VertexReco/interface/Vertex.h"
0006 #include "DataFormats/BTauReco/interface/TrackCountingTagInfo.h"
0007 #include "DataFormats/BTauReco/interface/IPTagInfo.h"
0008 #include "Math/GenVector/VectorUtil.h"
0009 #include "RecoBTau/JetTagComputer/interface/JetTagComputer.h"
0010
0011 template <class Container, class Base>
0012 class TemplatedTrackCountingComputer : public JetTagComputer {
0013 public:
0014   using Tokens = void;
0015
0016   typedef reco::IPTagInfo<Container, Base> TagInfo;
0017
0018   TemplatedTrackCountingComputer(const edm::ParameterSet& parameters) {
0019     m_minIP = parameters.existsAs<double>("minimumImpactParameter")
0020                   ? parameters.getParameter<double>("minimumImpactParameter")
0021                   : -1;
0022     m_useSignedIPSig = parameters.existsAs<bool>("useSignedImpactParameterSig")
0023                            ? parameters.getParameter<bool>("useSignedImpactParameterSig")
0024                            : true;
0025     m_nthTrack = parameters.getParameter<int>("nthTrack");
0026     m_ipType = parameters.getParameter<int>("impactParameterType");
0027     m_deltaR = parameters.getParameter<double>("deltaR");
0028     m_cutMaxDecayLen = parameters.getParameter<double>("maximumDecayLength");          //used
0029     m_cutMaxDistToAxis = parameters.getParameter<double>("maximumDistanceToJetAxis");  //used
0030     //
0031     // access track quality class; "any" takes everything
0032     //
0033     std::string trackQualityType = parameters.getParameter<std::string>("trackQualityClass");  //used
0034     m_trackQuality = reco::TrackBase::qualityByName(trackQualityType);
0035     m_useAllQualities = false;
0036     if (trackQualityType == "any" || trackQualityType == "Any" || trackQualityType == "ANY")
0037       m_useAllQualities = true;
0038
0039     uses("ipTagInfos");
0040
0041     useVariableJTA_ =
0042         parameters.existsAs<bool>("useVariableJTA") ? parameters.getParameter<bool>("useVariableJTA") : false;
0043     if (useVariableJTA_) {
0044       varJTApars = {parameters.getParameter<double>("a_dR"),
0045                     parameters.getParameter<double>("b_dR"),
0046                     parameters.getParameter<double>("a_pT"),
0047                     parameters.getParameter<double>("b_pT"),
0048                     parameters.getParameter<double>("min_pT"),
0049                     parameters.getParameter<double>("max_pT"),
0050                     parameters.getParameter<double>("min_pT_dRcut"),
0051                     parameters.getParameter<double>("max_pT_dRcut"),
0052                     parameters.getParameter<double>("max_pT_trackPTcut")};
0053     }
0054   }
0055
0056   float discriminator(const TagInfoHelper& ti) const override {
0057     const TagInfo& tkip = ti.get<TagInfo>();
0058     std::multiset<float> significances = orderedSignificances(tkip);
0059     std::multiset<float>::reverse_iterator nth = significances.rbegin();
0060     for (int i = 0; i < m_nthTrack - 1 && nth != significances.rend(); i++)
0061       nth++;
0062     if (nth != significances.rend())
0063       return *nth;
0064     else
0065       return -100.;
0066   }
0067
0068 protected:
0069   std::multiset<float> orderedSignificances(const TagInfo& tkip) const {
0070     const std::vector<reco::btag::TrackIPData>& impactParameters((tkip.impactParameterData()));
0071     const Container& tracks(tkip.selectedTracks());
0072     std::multiset<float> significances;
0073     int i = 0;
0074     if (tkip.primaryVertex().isNull()) {
0075       return std::multiset<float>();
0076     }
0077
0078     GlobalPoint pv(tkip.primaryVertex()->position().x(),
0079                    tkip.primaryVertex()->position().y(),
0080                    tkip.primaryVertex()->position().z());
0081
0082     for (std::vector<reco::btag::TrackIPData>::const_iterator it = impactParameters.begin();
0083          it != impactParameters.end();
0084          ++it, i++) {
0085       if (fabs(impactParameters[i].distanceToJetAxis.value()) < m_cutMaxDistToAxis &&  // distance to JetAxis
0086           (impactParameters[i].closestToJetAxis - pv).mag() < m_cutMaxDecayLen &&      // max decay len
0087           (m_useAllQualities == true ||
0088            reco::btag::toTrack(tracks[i])->quality(m_trackQuality)) &&       // use selected track qualities
0089           (fabs(((m_ipType == 0) ? it->ip3d : it->ip2d).value()) > m_minIP)  // minimum impact parameter
0090       ) {
0091         //calculate the signed or un-signed significance
0092         float signed_sig = ((m_ipType == 0) ? it->ip3d : it->ip2d).significance();
0093         float unsigned_sig = fabs(signed_sig);
0094         float significance = (m_useSignedIPSig) ? signed_sig : unsigned_sig;
0095
0096         if (useVariableJTA_) {
0097           if (tkip.variableJTA(varJTApars)[i])
0098             significances.insert(significance);
0099         } else  // no using variable JTA, use the default method
0100           if (m_deltaR <= 0 ||
0101               ROOT::Math::VectorUtil::DeltaR((*tkip.jet()).p4().Vect(), (*tracks[i]).momentum()) < m_deltaR)
0102             significances.insert(significance);
0103       }
0104     }
0105
0106     return significances;
0107   }
0108
0109   bool useVariableJTA_;
0110   reco::btag::variableJTAParameters varJTApars;
0111
0112   double m_minIP;
0113   bool m_useSignedIPSig;
0114
0115   int m_nthTrack;
0116   int m_ipType;
0117   double m_deltaR;
0118   double m_cutMaxDecayLen;
0119   double m_cutMaxDistToAxis;
0120   reco::TrackBase::TrackQuality m_trackQuality;
0121   bool m_useAllQualities;
0122 };
0123
0124 #endif  // ImpactParameter_TemplatedTrackCountingComputer_h