HLTDisplacedtktkVtxProducer.cc

Go to the documentation of this file.

#include <iostream>

#include "FWCore/Framework/interface/ESHandle.h"
#include "FWCore/MessageLogger/interface/MessageLogger.h"
#include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"

#include "DataFormats/Candidate/interface/CandidateFwd.h"
#include "DataFormats/RecoCandidate/interface/RecoCandidate.h"
#include "DataFormats/RecoCandidate/interface/RecoChargedCandidate.h"
#include "DataFormats/RecoCandidate/interface/RecoChargedCandidateFwd.h"
#include "DataFormats/Candidate/interface/Candidate.h"
#include "DataFormats/HLTReco/interface/TriggerFilterObjectWithRefs.h"
#include "DataFormats/HLTReco/interface/TriggerRefsCollections.h"

#include "DataFormats/TrackReco/interface/TrackFwd.h"
#include "DataFormats/TrackReco/interface/Track.h"
#include "TrackingTools/TransientTrack/interface/TransientTrackBuilder.h"
#include "TrackingTools/Records/interface/TransientTrackRecord.h"

#include "RecoVertex/KalmanVertexFit/interface/KalmanVertexFitter.h"
#include "RecoVertex/VertexPrimitives/interface/TransientVertex.h"
#include "DataFormats/VertexReco/interface/Vertex.h"
#include "DataFormats/VertexReco/interface/VertexFwd.h"
#include "DataFormats/Common/interface/RefToBase.h"

#include "HLTDisplacedtktkVtxProducer.h"

using namespace edm;
using namespace reco;
using namespace std;
using namespace trigger;
//
// constructors and destructor
//
HLTDisplacedtktkVtxProducer::HLTDisplacedtktkVtxProducer(const edm::ParameterSet& iConfig)
    : srcTag_(iConfig.getParameter<edm::InputTag>("Src")),
      srcToken_(consumes<reco::RecoChargedCandidateCollection>(srcTag_)),
      previousCandTag_(iConfig.getParameter<edm::InputTag>("PreviousCandTag")),
      previousCandToken_(consumes<trigger::TriggerFilterObjectWithRefs>(previousCandTag_)),
      maxEta_(iConfig.getParameter<double>("MaxEta")),
      minPt_(iConfig.getParameter<double>("MinPt")),
      minPtPair_(iConfig.getParameter<double>("MinPtPair")),
      minInvMass_(iConfig.getParameter<double>("MinInvMass")),
      maxInvMass_(iConfig.getParameter<double>("MaxInvMass")),
      massParticle1_(iConfig.getParameter<double>("massParticle1")),
      massParticle2_(iConfig.getParameter<double>("massParticle2")),
      chargeOpt_(iConfig.getParameter<int>("ChargeOpt")),
      triggerTypeDaughters_(iConfig.getParameter<int>("triggerTypeDaughters"))

{
  produces<VertexCollection>();
}

HLTDisplacedtktkVtxProducer::~HLTDisplacedtktkVtxProducer() = default;

void HLTDisplacedtktkVtxProducer::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
  edm::ParameterSetDescription desc;
  desc.add<edm::InputTag>("Src", edm::InputTag("hltL3MuonCandidates"));
  desc.add<edm::InputTag>("PreviousCandTag", edm::InputTag(""));
  desc.add<double>("MaxEta", 2.5);
  desc.add<double>("MinPt", 0.0);
  desc.add<double>("MinPtPair", 0.0);
  desc.add<double>("MinInvMass", 1.0);
  desc.add<double>("MaxInvMass", 20.0);
  desc.add<double>("massParticle1", 0.1396);
  desc.add<double>("massParticle2", 0.4937);
  desc.add<int>("ChargeOpt", -1);
  desc.add<int>("triggerTypeDaughters", 0);

  descriptions.add("hltDisplacedtktkVtxProducer", desc);
}

// ------------ method called on each new Event  ------------
void HLTDisplacedtktkVtxProducer::produce(edm::Event& iEvent, const edm::EventSetup& iSetup) {
  double const firstTrackMass = massParticle1_;
  double const firstTrackMass2 = firstTrackMass * firstTrackMass;
  double const secondTrackMass = massParticle2_;
  double const secondTrackMass2 = secondTrackMass * secondTrackMass;

  // get hold of track trks
  Handle<RecoChargedCandidateCollection> trackcands;
  iEvent.getByToken(srcToken_, trackcands);

  //get the transient track builder:
  edm::ESHandle<TransientTrackBuilder> theB;
  iSetup.get<TransientTrackRecord>().get("TransientTrackBuilder", theB);

  std::unique_ptr<VertexCollection> vertexCollection(new VertexCollection());

  // look at all trackcands,  check cuts and make vertices
  double e1, e2;
  Particle::LorentzVector p, p1, p2;

  RecoChargedCandidateCollection::const_iterator cand1;
  RecoChargedCandidateCollection::const_iterator cand2;

  // get the objects passing the previous filter
  Handle<TriggerFilterObjectWithRefs> previousCands;
  iEvent.getByToken(previousCandToken_, previousCands);

  vector<RecoChargedCandidateRef> vPrevCands;
  previousCands->getObjects(triggerTypeDaughters_, vPrevCands);

  std::vector<bool> candComp;
  for (cand1 = trackcands->begin(); cand1 != trackcands->end(); cand1++)
    candComp.push_back(checkPreviousCand(cand1->get<TrackRef>(), vPrevCands));

  for (cand1 = trackcands->begin(); cand1 != trackcands->end(); cand1++) {
    TrackRef tk1 = cand1->get<TrackRef>();
    LogDebug("HLTDisplacedtktkVtxProducer") << " 1st track in loop: q*pt= " << cand1->charge() * cand1->pt()
                                            << ", eta= " << cand1->eta() << ", hits= " << tk1->numberOfValidHits();

    //first check if this track passed the previous filter
    if (!candComp[cand1 - trackcands->begin()])
      continue;
    // if( ! checkPreviousCand( tk1, vPrevCands) ) continue;

    // cuts
    if (abs(cand1->eta()) > maxEta_)
      continue;
    if (cand1->pt() < minPt_)
      continue;

    cand2 = trackcands->begin();
    if (massParticle1_ == massParticle2_) {
      cand2 = cand1 + 1;
    }

    for (; cand2 != trackcands->end(); cand2++) {
      TrackRef tk2 = cand2->get<TrackRef>();
      if (tk1 == tk2)
        continue;

      // eta cut
      LogDebug("HLTDisplacedtktkVtxProducer")
          << " 2nd track in loop: q*pt= " << cand2->charge() * cand2->pt() << ", eta= " << cand2->eta()
          << ", hits= " << tk2->numberOfValidHits() << ", d0= " << tk2->d0();
      //first check if this track passed the previous filter
      if (!candComp[cand2 - trackcands->begin()])
        continue;
      // if( ! checkPreviousCand( tk2, vPrevCands) ) continue;

      // cuts
      if (abs(cand2->eta()) > maxEta_)
        continue;
      if (cand2->pt() < minPt_)
        continue;

      // opposite sign or same sign
      if (chargeOpt_ < 0) {
        if (cand1->charge() * cand2->charge() > 0)
          continue;
      } else if (chargeOpt_ > 0) {
        if (cand1->charge() * cand2->charge() < 0)
          continue;
      }

      // Combined ditrack system
      e1 = sqrt(cand1->momentum().Mag2() + firstTrackMass2);
      e2 = sqrt(cand2->momentum().Mag2() + secondTrackMass2);
      p1 = Particle::LorentzVector(cand1->px(), cand1->py(), cand1->pz(), e1);
      p2 = Particle::LorentzVector(cand2->px(), cand2->py(), cand2->pz(), e2);
      p = p1 + p2;

      if (p.pt() < minPtPair_)
        continue;

      double invmass = abs(p.mass());
      LogDebug("HLTDisplacedtktkVtxProducer") << " ... 1-2 invmass= " << invmass;

      if (invmass < minInvMass_)
        continue;
      if (invmass > maxInvMass_)
        continue;

      // do the vertex fit
      vector<TransientTrack> t_tks;
      TransientTrack ttkp1 = (*theB).build(&tk1);
      TransientTrack ttkp2 = (*theB).build(&tk2);
      t_tks.push_back(ttkp1);
      t_tks.push_back(ttkp2);

      if (t_tks.size() != 2)
        continue;

      KalmanVertexFitter kvf;
      TransientVertex tv = kvf.vertex(t_tks);

      if (!tv.isValid())
        continue;

      Vertex vertex = tv;

      // put vertex in the event
      vertexCollection->push_back(vertex);
    }
  }
  iEvent.put(std::move(vertexCollection));
}

bool HLTDisplacedtktkVtxProducer::checkPreviousCand(const TrackRef& trackref,
                                                    const vector<RecoChargedCandidateRef>& refVect) const {
  bool ok = false;
  for (auto& i : refVect) {
    if (i->get<TrackRef>() == trackref) {
      ok = true;
      break;
    }
  }
  return ok;
}