PFTauSecondaryVertexProducer.cc

Go to the documentation of this file.

/* class PFTauSecondaryVertexProducer
 * EDProducer of the 
 * authors: Ian M. Nugent
 * This work is based on the impact parameter work by Rosamaria Venditti and reconstructing the 3 prong taus.
 * The idea of the fully reconstructing the tau using a kinematic fit comes from
 * Lars Perchalla and Philip Sauerland Theses under Achim Stahl supervision. This
 * work was continued by Ian M. Nugent and Vladimir Cherepanov.
 * Thanks goes to Christian Veelken and Evan Klose Friis for their help and suggestions.
 */

#include "DataFormats/TauReco/interface/PFTauTagInfo.h"
#include "DataFormats/VertexReco/interface/Vertex.h"
#include "DataFormats/VertexReco/interface/VertexFwd.h"

#include "FWCore/Framework/interface/EventSetup.h"
#include "FWCore/Framework/interface/ESHandle.h"
#include "FWCore/Framework/interface/Event.h"
#include "FWCore/Framework/interface/MakerMacros.h"
#include "FWCore/Framework/interface/Frameworkfwd.h"
#include "FWCore/Framework/interface/global/EDProducer.h"
#include "FWCore/ParameterSet/interface/ParameterSet.h"
#include "FWCore/MessageLogger/interface/MessageLogger.h"
#include "FWCore/Utilities/interface/Exception.h"

#include <FWCore/ParameterSet/interface/ConfigurationDescriptions.h>
#include <FWCore/ParameterSet/interface/ParameterSetDescription.h>

#include "TrackingTools/TransientTrack/interface/TransientTrackBuilder.h"
#include "TrackingTools/Records/interface/TransientTrackRecord.h"
#include "RecoVertex/VertexPrimitives/interface/TransientVertex.h"
#include "RecoVertex/KalmanVertexFit/interface/KalmanVertexFitter.h"

#include "DataFormats/TauReco/interface/PFTau.h"
#include "DataFormats/TauReco/interface/PFTauFwd.h"
#include "DataFormats/VertexReco/interface/Vertex.h"
#include "DataFormats/VertexReco/interface/VertexFwd.h"
#include "DataFormats/TrackReco/interface/Track.h"
#include "DataFormats/TrackReco/interface/TrackFwd.h"
#include "DataFormats/GsfTrackReco/interface/GsfTrack.h"
#include "DataFormats/GsfTrackReco/interface/GsfTrackFwd.h"
#include "DataFormats/Common/interface/RefToBase.h"
#include "DataFormats/PatCandidates/interface/PackedCandidate.h"

#include "DataFormats/Common/interface/Association.h"
#include "DataFormats/Common/interface/AssociationVector.h"
#include "DataFormats/Common/interface/RefProd.h"

#include <memory>

using namespace reco;
using namespace edm;
using namespace std;

class PFTauSecondaryVertexProducer : public edm::global::EDProducer<> {
public:
  enum Alg { useInputPV = 0, usePVwithMaxSumPt, useTauPV };

  explicit PFTauSecondaryVertexProducer(const edm::ParameterSet& iConfig);
  ~PFTauSecondaryVertexProducer() override;
  void produce(edm::StreamID, edm::Event&, const edm::EventSetup&) const override;
  static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);

private:
  const edm::InputTag PFTauTag_;
  const edm::EDGetTokenT<std::vector<reco::PFTau>> PFTauToken_;
};

PFTauSecondaryVertexProducer::PFTauSecondaryVertexProducer(const edm::ParameterSet& iConfig)
    : PFTauTag_(iConfig.getParameter<edm::InputTag>("PFTauTag")),
      PFTauToken_(consumes<std::vector<reco::PFTau>>(iConfig.getParameter<edm::InputTag>("PFTauTag"))) {
  produces<edm::AssociationVector<PFTauRefProd, std::vector<std::vector<reco::VertexRef>>>>();
  produces<VertexCollection>("PFTauSecondaryVertices");
}

PFTauSecondaryVertexProducer::~PFTauSecondaryVertexProducer() {}

namespace {
  const reco::Track* getTrack(const reco::Candidate& cand) {
    const reco::PFCandidate* pfCand = dynamic_cast<const reco::PFCandidate*>(&cand);
    if (pfCand != nullptr) {
      if (pfCand->trackRef().isNonnull())
        return &*pfCand->trackRef();
      else if (pfCand->gsfTrackRef().isNonnull())
        return &*pfCand->gsfTrackRef();
    }
    const pat::PackedCandidate* pCand = dynamic_cast<const pat::PackedCandidate*>(&cand);
    if (pCand != nullptr && pCand->hasTrackDetails())
      return &pCand->pseudoTrack();
    return nullptr;
  }
}  // namespace
void PFTauSecondaryVertexProducer::produce(edm::StreamID, edm::Event& iEvent, const edm::EventSetup& iSetup) const {
  // Obtain
  edm::ESHandle<TransientTrackBuilder> transTrackBuilder;
  iSetup.get<TransientTrackRecord>().get("TransientTrackBuilder", transTrackBuilder);

  edm::Handle<std::vector<reco::PFTau>> Tau;
  iEvent.getByToken(PFTauToken_, Tau);

  // Set Association Map
  auto AVPFTauSV = std::make_unique<edm::AssociationVector<PFTauRefProd, std::vector<std::vector<reco::VertexRef>>>>(
      PFTauRefProd(Tau));
  auto VertexCollection_out = std::make_unique<VertexCollection>();
  reco::VertexRefProd VertexRefProd_out = iEvent.getRefBeforePut<reco::VertexCollection>("PFTauSecondaryVertices");

  // For each Tau Run Algorithim
  if (Tau.isValid()) {
    for (reco::PFTauCollection::size_type iPFTau = 0; iPFTau < Tau->size(); iPFTau++) {
      reco::PFTauRef RefPFTau(Tau, iPFTau);
      std::vector<reco::VertexRef> SV;
      if (RefPFTau->decayMode() >= 5) {
        // Get tracks form PFTau daugthers
        std::vector<reco::TransientTrack> transTrk;
        TransientVertex transVtx;
        const std::vector<edm::Ptr<reco::Candidate>> cands = RefPFTau->signalChargedHadrCands();
        for (const auto& cand : cands) {
          if (cand.isNull())
            continue;
          const reco::Track* track = getTrack(*cand);
          if (track != nullptr)
            transTrk.push_back(transTrackBuilder->build(*track));
        }
        // Fit the secondary vertex
        bool FitOk(true);
        KalmanVertexFitter kvf(true);
        if (transTrk.size() > 1) {
          transVtx = kvf.vertex(transTrk);  //KalmanVertexFitter
        } else {
          FitOk = false;
        }
        if (!transVtx.hasRefittedTracks())
          FitOk = false;
        if (transVtx.refittedTracks().size() != transTrk.size())
          FitOk = false;
        if (FitOk) {
          SV.push_back(reco::VertexRef(VertexRefProd_out, VertexCollection_out->size()));
          VertexCollection_out->push_back(transVtx);
        }
      }
      AVPFTauSV->setValue(iPFTau, SV);
    }
  }
  iEvent.put(std::move(VertexCollection_out), "PFTauSecondaryVertices");
  iEvent.put(std::move(AVPFTauSV));
}

void PFTauSecondaryVertexProducer::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
  // PFTauSecondaryVertexProducer
  edm::ParameterSetDescription desc;
  desc.add<edm::InputTag>("PFTauTag", edm::InputTag("hpsPFTauProducer"));
  descriptions.add("PFTauSecondaryVertexProducer", desc);
}

DEFINE_FWK_MODULE(PFTauSecondaryVertexProducer);