Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2025-07-01 00:07:41

0001 // -*- C++ -*-
0002 //
0003 // Package:    PhysicsTools/NanoAOD
0004 // Class:      PVertexBPHTable
0005 //
0006 /*
0007  Simple primary vertex table  
0008  Description: [one line class summary]
0009 
0010  Implementation:
0011      [Notes on implementation]
0012 */
0013 //
0014 // Original Author:  gkaratha
0015 //         Created:  Mon, 28 Aug 2024 09:26:39 GMT
0016 //
0017 //
0018 
0019 // system include files
0020 #include <memory>
0021 
0022 // user include files
0023 #include "FWCore/Framework/interface/Frameworkfwd.h"
0024 #include "FWCore/Framework/interface/stream/EDProducer.h"
0025 
0026 #include "FWCore/Framework/interface/Event.h"
0027 #include "FWCore/Framework/interface/MakerMacros.h"
0028 
0029 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0030 #include "FWCore/Utilities/interface/StreamID.h"
0031 
0032 #include "DataFormats/VertexReco/interface/Vertex.h"
0033 #include "DataFormats/Candidate/interface/VertexCompositePtrCandidate.h"
0034 
0035 #include "CommonTools/Utils/interface/StringCutObjectSelector.h"
0036 
0037 #include "DataFormats/NanoAOD/interface/FlatTable.h"
0038 #include "RecoVertex/VertexTools/interface/VertexDistance3D.h"
0039 #include "RecoVertex/VertexTools/interface/VertexDistanceXY.h"
0040 #include "RecoVertex/VertexPrimitives/interface/ConvertToFromReco.h"
0041 #include "RecoVertex/VertexPrimitives/interface/VertexState.h"
0042 #include "DataFormats/Common/interface/ValueMap.h"
0043 #include "DataFormats/PatCandidates/interface/CompositeCandidate.h"
0044 
0045 //
0046 // class declaration
0047 //
0048 
0049 class PVertexBPHTable : public edm::stream::EDProducer<> {
0050 public:
0051   explicit PVertexBPHTable(const edm::ParameterSet&);
0052   ~PVertexBPHTable() override;
0053 
0054   static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0055 
0056 private:
0057   void beginStream(edm::StreamID) override;
0058   void produce(edm::Event&, const edm::EventSetup&) override;
0059   void endStream() override;
0060 
0061   // ----------member data ---------------------------
0062 
0063   const edm::EDGetTokenT<std::vector<reco::Vertex>> pvs_;
0064   const edm::EDGetTokenT<edm::ValueMap<float>> pvsScore_;
0065   const edm::EDGetTokenT<pat::CompositeCandidateCollection> dileptonToken_;
0066   const double maxDzDilep_;
0067   const StringCutObjectSelector<reco::Vertex> goodPvCut_;
0068   const std::string pvName_;
0069 };
0070 
0071 //
0072 // constructors and destructor
0073 //
0074 PVertexBPHTable::PVertexBPHTable(const edm::ParameterSet& params)
0075     : pvs_(consumes<std::vector<reco::Vertex>>(params.getParameter<edm::InputTag>("pvSrc"))),
0076       pvsScore_(consumes<edm::ValueMap<float>>(params.getParameter<edm::InputTag>("pvSrc"))),
0077       dileptonToken_(consumes<pat::CompositeCandidateCollection>(params.getParameter<edm::InputTag>("dileptons"))),
0078       maxDzDilep_(params.getParameter<double>("maxDzDilep")),
0079       goodPvCut_(params.getParameter<std::string>("goodPvCut"), true),
0080       pvName_(params.getParameter<std::string>("pvName"))
0081 
0082 {
0083   produces<nanoaod::FlatTable>("pv");
0084   produces<edm::PtrVector<reco::Candidate>>();
0085 }
0086 
0087 PVertexBPHTable::~PVertexBPHTable() {
0088   // do anything here that needs to be done at destruction time
0089   // (e.g. close files, deallocate resources etc.)
0090 }
0091 
0092 //
0093 // member functions
0094 //
0095 
0096 // ------------ method called to produce the data  ------------
0097 
0098 void PVertexBPHTable::produce(edm::Event& iEvent, const edm::EventSetup& iSetup) {
0099   using namespace edm;
0100   const auto& pvsScoreProd = iEvent.get(pvsScore_);
0101   auto pvsCol = iEvent.getHandle(pvs_);
0102 
0103   edm::Handle<pat::CompositeCandidateCollection> dileptons;
0104   iEvent.getByToken(dileptonToken_, dileptons);
0105 
0106   auto selCandPv = std::make_unique<PtrVector<reco::Candidate>>();
0107   std::vector<float> pvscore, chi2, covXX, covYY, covZZ, covXY, covXZ, covYZ, vx, vy, vz, pt, eta, phi, mass, ndof;
0108   std::vector<int> charge, ntracks;
0109 
0110   size_t i = 0;
0111   for (const auto& pv : *pvsCol) {
0112     if (!goodPvCut_(pv) and i != 0) {
0113       i++;
0114       continue;
0115     }
0116     bool within_dz = false;
0117     for (const pat::CompositeCandidate& dilep : *dileptons) {
0118       if (fabs(pv.z() - dilep.vz()) < maxDzDilep_ && maxDzDilep_ > 0)
0119         within_dz = true;
0120     }
0121     if (!within_dz)
0122       continue;
0123 
0124     // int sum_charge = 0;
0125     pvscore.push_back(pvsScoreProd.get(pvsCol.id(), i));
0126     ntracks.push_back(pv.tracksSize());
0127     chi2.push_back(pv.chi2());
0128     covXX.push_back(pv.covariance(0, 0));
0129     covYY.push_back(pv.covariance(1, 1));
0130     covZZ.push_back(pv.covariance(2, 2));
0131     covXY.push_back(pv.covariance(0, 1));
0132     covXZ.push_back(pv.covariance(0, 2));
0133     covYZ.push_back(pv.covariance(1, 2));
0134     vx.push_back(pv.x());
0135     vy.push_back(pv.y());
0136     vz.push_back(pv.z());
0137     pt.push_back(pv.p4().pt());
0138     eta.push_back(pv.p4().eta());
0139     phi.push_back(pv.p4().phi());
0140     mass.push_back(pv.p4().M());
0141     ndof.push_back(pv.ndof());
0142     i++;
0143   }
0144   auto table = std::make_unique<nanoaod::FlatTable>(pvscore.size(), pvName_, false, false);
0145   table->addColumn<float>("score", pvscore, "", 12);
0146   table->addColumn<float>("vx", vx, "", 16);
0147   table->addColumn<float>("vy", vy, "", 16);
0148   table->addColumn<float>("vz", vz, "", 16);
0149   table->addColumn<float>("pt", pt, "", 12);
0150   table->addColumn<float>("eta", eta, "", 12);
0151   table->addColumn<float>("phi", phi, "", 12);
0152   table->addColumn<float>("mass", mass, "", 12);
0153   table->addColumn<float>("chi2", chi2, "", 12);
0154   table->addColumn<float>("ndof", ndof, "", 12);
0155   table->addColumn<float>("covXX", covXX, "", 12);
0156   table->addColumn<float>("covYY", covYY, "", 12);
0157   table->addColumn<float>("covZZ", covZZ, "", 12);
0158   table->addColumn<float>("covXY", covXY, "", 12);
0159   table->addColumn<float>("covXZ", covXZ, "", 12);
0160   table->addColumn<float>("covYZ", covYZ, "", 12);
0161   table->addColumn<uint8_t>("ntracks", ntracks, "");
0162 
0163   iEvent.put(std::move(table), "pv");
0164 }
0165 
0166 // ------------ method called once each stream before processing any runs, lumis or events  ------------
0167 void PVertexBPHTable::beginStream(edm::StreamID) {}
0168 
0169 // ------------ method called once each stream after processing all runs, lumis and events  ------------
0170 void PVertexBPHTable::endStream() {}
0171 
0172 // ------------ method fills 'descriptions' with the allowed parameters for the module  ------------
0173 void PVertexBPHTable::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0174   edm::ParameterSetDescription desc;
0175 
0176   desc.add<edm::InputTag>("pvSrc")->setComment(
0177       "std::vector<reco::Vertex> and ValueMap<float> primary vertex input collections");
0178 
0179   desc.add<edm::InputTag>("dileptons")->setComment("input dilepton collection");
0180 
0181   desc.add<double>("maxDzDilep")->setComment("maxDz cut wrt dilepton vertices to keep only useful PVs"),
0182 
0183       desc.add<std::string>("goodPvCut")->setComment("selection on the primary vertex");
0184 
0185   desc.add<std::string>("pvName")->setComment("name of the flat table ouput");
0186 
0187   descriptions.addWithDefaultLabel(desc);
0188 }
0189 
0190 //define this as a plug-in
0191 DEFINE_FWK_MODULE(PVertexBPHTable);