Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2022-09-16 00:05:54

0001 
0002 #include "RecoTauTag/HLTProducers/interface/L2TauPixelTrackMatch.h"
0003 #include "FWCore/Utilities/interface/Exception.h"
0004 #include "DataFormats/Common/interface/Handle.h"
0005 #include "DataFormats/Common/interface/RefToBase.h"
0006 #include "DataFormats/JetReco/interface/CaloJet.h"
0007 #include "DataFormats/JetReco/interface/CaloJetCollection.h"
0008 #include "DataFormats/HLTReco/interface/TriggerTypeDefs.h"
0009 #include "DataFormats/Math/interface/deltaPhi.h"
0010 
0011 // all these debug printouts will need to be removed at some point
0012 //#define DBG_PRINT(arg) (arg)
0013 #define DBG_PRINT(arg)
0014 
0015 L2TauPixelTrackMatch::L2TauPixelTrackMatch(const edm::ParameterSet& conf) {
0016   m_jetSrc = consumes<trigger::TriggerFilterObjectWithRefs>(conf.getParameter<edm::InputTag>("JetSrc"));
0017   m_jetMinPt = conf.getParameter<double>("JetMinPt");
0018   m_jetMaxEta = conf.getParameter<double>("JetMaxEta");
0019   //m_jetMinN        = conf.getParameter<int>("JetMinN");
0020   m_trackSrc = consumes<reco::TrackCollection>(conf.getParameter<edm::InputTag>("TrackSrc"));
0021   m_trackMinPt = conf.getParameter<double>("TrackMinPt");
0022   m_deltaR = conf.getParameter<double>("deltaR");
0023   m_beamSpotTag = consumes<reco::BeamSpot>(conf.getParameter<edm::InputTag>("BeamSpotSrc"));
0024 
0025   produces<reco::CaloJetCollection>();
0026 }
0027 
0028 L2TauPixelTrackMatch::~L2TauPixelTrackMatch() {}
0029 
0030 void L2TauPixelTrackMatch::produce(edm::StreamID, edm::Event& ev, const edm::EventSetup& es) const {
0031   using namespace std;
0032   using namespace reco;
0033 
0034   // *** Pick up beam spot ***
0035 
0036   // use beam spot for vertex x,y
0037   edm::Handle<BeamSpot> bsHandle;
0038   ev.getByToken(m_beamSpotTag, bsHandle);
0039   const reco::BeamSpot& bs = *bsHandle;
0040   math::XYZPoint beam_spot(bs.x0(), bs.y0(), bs.z0());
0041   DBG_PRINT(cout << endl << "beamspot " << beam_spot << endl);
0042 
0043   // *** Pick up pixel tracks ***
0044 
0045   edm::Handle<TrackCollection> tracksHandle;
0046   ev.getByToken(m_trackSrc, tracksHandle);
0047 
0048   // *** Pick up L2 tau jets that were previously selected by some other filter ***
0049 
0050   // first, get L2 object refs by the label
0051   edm::Handle<trigger::TriggerFilterObjectWithRefs> jetsHandle;
0052   ev.getByToken(m_jetSrc, jetsHandle);
0053 
0054   // now we can get pre-selected L2 tau jets
0055   std::vector<CaloJetRef> tau_jets;
0056   jetsHandle->getObjects(trigger::TriggerTau, tau_jets);
0057   const size_t n_jets = tau_jets.size();
0058 
0059   // *** Selects interesting tracks ***
0060 
0061   vector<TinyTrack> good_tracks;
0062   for (TrackCollection::const_iterator itrk = tracksHandle->begin(); itrk != tracksHandle->end(); ++itrk) {
0063     if (itrk->pt() < m_trackMinPt)
0064       continue;
0065     if (std::abs(itrk->eta()) > m_jetMaxEta + m_deltaR)
0066       continue;
0067 
0068     TinyTrack trk;
0069     trk.pt = itrk->pt();
0070     trk.phi = itrk->phi();
0071     trk.eta = itrk->eta();
0072     double dz = itrk->dz(beam_spot);
0073     trk.vtx = math::XYZPoint(bs.x(dz), bs.y(dz), dz);
0074     good_tracks.push_back(trk);
0075   }
0076   DBG_PRINT(cout << "got " << good_tracks.size() << " good tracks;   " << n_jets << " tau jets" << endl);
0077 
0078   // *** Match tau jets to intertesting tracks  and assign them new vertices ***
0079 
0080   // the new product
0081   std::unique_ptr<CaloJetCollection> new_tau_jets(new CaloJetCollection);
0082   int n_uniq = 0;
0083   if (!good_tracks.empty())
0084     for (size_t i = 0; i < n_jets; ++i) {
0085       reco::CaloJetRef jet = tau_jets[i];
0086       if (jet->pt() < m_jetMinPt || std::abs(jet->eta()) > m_jetMaxEta)
0087         continue;
0088 
0089       DBG_PRINT(cout << i << " :" << endl);
0090 
0091       size_t n0 = new_tau_jets->size();
0092 
0093       for (vector<TinyTrack>::const_iterator itrk = good_tracks.begin(); itrk != good_tracks.end(); ++itrk) {
0094         DBG_PRINT(cout << "  trk pt,eta,phi,z: " << itrk->pt << " " << itrk->eta << " " << itrk->phi << " "
0095                        << itrk->vtx.z() << " \t\t ");
0096 
0097         math::XYZTLorentzVector new_jet_dir = Jet::physicsP4(itrk->vtx, *jet, itrk->vtx);
0098         float dphi = reco::deltaPhi(new_jet_dir.phi(), itrk->phi);
0099         float deta = new_jet_dir.eta() - itrk->eta;
0100 
0101         DBG_PRINT(cout << " jet pt,deta,dphi,dr: " << jet->pt() << " " << deta << " " << dphi << " "
0102                        << sqrt(dphi * dphi + deta * deta) << endl);
0103 
0104         if (dphi * dphi + deta * deta > m_deltaR * m_deltaR)
0105           continue;
0106 
0107         DBG_PRINT(cout << "  jet-trk match!" << endl);
0108 
0109         // create a jet copy and assign a new vertex to it
0110         CaloJet new_jet = *jet;
0111         new_jet.setVertex(itrk->vtx);
0112 
0113         new_tau_jets->push_back(new_jet);
0114       }
0115       DBG_PRINT(cout << "  nmatchedjets " << new_tau_jets->size() - n0 << endl);
0116       if (new_tau_jets->size() - n0 > 0)
0117         n_uniq++;
0118 
0119       ///if (jet_with_vertices.size()) new_tau_jets->push_back(jet_with_vertices);
0120     }
0121   DBG_PRINT(cout << "n_uniq_matched_jets " << n_uniq << endl << "storing njets " << new_tau_jets->size() << endl);
0122 
0123   // store the result
0124   ev.put(std::move(new_tau_jets));
0125 }