From dd9dc87657596a29ae1c6c45b1be7636709a053e Mon Sep 17 00:00:00 2001 From: Michal Bluj Date: Thu, 7 Dec 2023 13:47:06 +0100 Subject: [PATCH] Produce lepton time-life info table without using intermediate lepton collection --- DataFormats/NanoAOD/BuildFile.xml | 1 + .../NanoAOD/interface/LeptonTimeLifeInfo.h | 73 ++++ DataFormats/NanoAOD/src/LeptonTimeLifeInfo.cc | 15 + DataFormats/NanoAOD/src/classes.h | 1 + DataFormats/NanoAOD/src/classes_def.xml | 1 + PhysicsTools/NanoAOD/plugins/BuildFile.xml | 1 + .../LeptonTimeLifeInfoTableProducer.cc | 340 ++++++++++++++++++ PhysicsTools/NanoAOD/python/electrons_cff.py | 17 +- .../python/leptonTimeLifeInfo_common_cff.py | 56 +-- PhysicsTools/NanoAOD/python/muons_cff.py | 17 +- PhysicsTools/NanoAOD/python/taus_cff.py | 16 +- PhysicsTools/NanoAOD/python/vertices_cff.py | 14 +- .../plugins/PATLeptonTimeLifeInfoUpdater.cc | 331 ----------------- 13 files changed, 487 insertions(+), 396 deletions(-) create mode 100644 DataFormats/NanoAOD/interface/LeptonTimeLifeInfo.h create mode 100644 DataFormats/NanoAOD/src/LeptonTimeLifeInfo.cc create mode 100644 PhysicsTools/NanoAOD/plugins/LeptonTimeLifeInfoTableProducer.cc delete mode 100644 PhysicsTools/PatAlgos/plugins/PATLeptonTimeLifeInfoUpdater.cc diff --git a/DataFormats/NanoAOD/BuildFile.xml b/DataFormats/NanoAOD/BuildFile.xml index 048178ece8470..59a9c3e543991 100644 --- a/DataFormats/NanoAOD/BuildFile.xml +++ b/DataFormats/NanoAOD/BuildFile.xml @@ -1,4 +1,5 @@ + diff --git a/DataFormats/NanoAOD/interface/LeptonTimeLifeInfo.h b/DataFormats/NanoAOD/interface/LeptonTimeLifeInfo.h new file mode 100644 index 0000000000000..5e1a8abdb0149 --- /dev/null +++ b/DataFormats/NanoAOD/interface/LeptonTimeLifeInfo.h @@ -0,0 +1,73 @@ +#ifndef DataFormats_NanoAOD_LeptonTimeLifeInfo_h +#define DataFormats_NanoAOD_LeptonTimeLifeInfo_h + +/** + \class LeptonTimeLifeInfo + \brief Structure to hold lepton life-time information + + \author Michal Bluj, NCBJ, Warsaw +*/ + +#include "DataFormats/GeometryCommonDetAlgo/interface/GlobalError.h" +#include "DataFormats/GeometryCommonDetAlgo/interface/Measurement1D.h" +#include "DataFormats/GeometryVector/interface/GlobalPoint.h" +#include "DataFormats/GeometryVector/interface/GlobalVector.h" +#include "DataFormats/TrackReco/interface/Track.h" +#include "DataFormats/VertexReco/interface/Vertex.h" + +class LeptonTimeLifeInfo { +public: + LeptonTimeLifeInfo(); + ~LeptonTimeLifeInfo() {} + + // Secondary vertex + void setSV(reco::Vertex sv) { + sv_ = sv; + hasSV_ = true; + } + const reco::Vertex& sv() const { return sv_; } + bool hasSV() const { return hasSV_; } + void setFlightVector(GlobalVector flight_vec, const GlobalError& flight_cov) { + flight_vec_ = flight_vec; + flight_cov_ = flight_cov; + } + // Flight-path + const GlobalVector& flightVector() const { return flight_vec_; } + const GlobalError& flightCovariance() const { return flight_cov_; } + void setFlightLength(Measurement1D flightLength) { flightLength_ = flightLength; } + const Measurement1D& flightLength() const { return flightLength_; } + // Point of closest approach + void setPCA(GlobalPoint pca, const GlobalError& pca_cov) { + pca_ = pca; + pca_cov_ = pca_cov_; + } + const GlobalPoint& pca() const { return pca_; } + const GlobalError& pcaCovariance() const { return pca_cov_; } + // Impact parameter + void setIP(GlobalVector ip_vec, const GlobalError& ip_cov) { + ip_vec_ = ip_vec; + ip_cov_ = ip_cov; + } + const GlobalVector& ipVector() const { return ip_vec_; } + const GlobalError& ipCovariance() const { return ip_cov_; } + void setIPLength(Measurement1D ipLength) { ipLength_ = ipLength; } + const Measurement1D& ipLength() const { return ipLength_; } + // Track + void setTrack(const reco::Track* track) { track_ = track; } + const reco::Track* track() const { return track_; } + bool hasTrack() const { return track_ != nullptr; } + void setBFiled_z(float bField_z) { bField_z_ = bField_z; } + float bField_z() const { return bField_z_; } + +private: + bool hasSV_; + reco::Vertex sv_; + GlobalVector flight_vec_, ip_vec_; + GlobalPoint pca_; + GlobalError flight_cov_, pca_cov_, ip_cov_; + Measurement1D flightLength_, ipLength_; + const reco::Track* track_; + float bField_z_; +}; + +#endif diff --git a/DataFormats/NanoAOD/src/LeptonTimeLifeInfo.cc b/DataFormats/NanoAOD/src/LeptonTimeLifeInfo.cc new file mode 100644 index 0000000000000..b3c17fff8912a --- /dev/null +++ b/DataFormats/NanoAOD/src/LeptonTimeLifeInfo.cc @@ -0,0 +1,15 @@ +#include "DataFormats/NanoAOD/interface/LeptonTimeLifeInfo.h" + +LeptonTimeLifeInfo::LeptonTimeLifeInfo() + : hasSV_(false), + sv_(reco::Vertex()), + flight_vec_(GlobalVector()), + ip_vec_(GlobalVector()), + pca_(GlobalPoint()), + flight_cov_(GlobalError()), + pca_cov_(GlobalError()), + ip_cov_(GlobalError()), + flightLength_(Measurement1D()), + ipLength_(Measurement1D()), + track_(nullptr), + bField_z_(0.){}; diff --git a/DataFormats/NanoAOD/src/classes.h b/DataFormats/NanoAOD/src/classes.h index acdf6f161f7cb..905d954a31cf3 100644 --- a/DataFormats/NanoAOD/src/classes.h +++ b/DataFormats/NanoAOD/src/classes.h @@ -3,4 +3,5 @@ #include "DataFormats/NanoAOD/interface/FlatTable.h" #include "DataFormats/NanoAOD/interface/MergeableCounterTable.h" #include "DataFormats/NanoAOD/interface/UniqueString.h" +#include "DataFormats/NanoAOD/interface/LeptonTimeLifeInfo.h" #include "DataFormats/Common/interface/Wrapper.h" diff --git a/DataFormats/NanoAOD/src/classes_def.xml b/DataFormats/NanoAOD/src/classes_def.xml index 582d353813d4d..f71aa91d9d513 100644 --- a/DataFormats/NanoAOD/src/classes_def.xml +++ b/DataFormats/NanoAOD/src/classes_def.xml @@ -42,5 +42,6 @@ + diff --git a/PhysicsTools/NanoAOD/plugins/BuildFile.xml b/PhysicsTools/NanoAOD/plugins/BuildFile.xml index 2d06b9ec6bfd6..037f2983eb364 100644 --- a/PhysicsTools/NanoAOD/plugins/BuildFile.xml +++ b/PhysicsTools/NanoAOD/plugins/BuildFile.xml @@ -25,6 +25,7 @@ + diff --git a/PhysicsTools/NanoAOD/plugins/LeptonTimeLifeInfoTableProducer.cc b/PhysicsTools/NanoAOD/plugins/LeptonTimeLifeInfoTableProducer.cc new file mode 100644 index 0000000000000..568c969cf2d0b --- /dev/null +++ b/PhysicsTools/NanoAOD/plugins/LeptonTimeLifeInfoTableProducer.cc @@ -0,0 +1,340 @@ +/** + \class LeptonTimeLifeInfoTableProducer + \brief Produces FlatTable with lepton life-time information + + \author Michal Bluj, NCBJ, Warsaw +*/ + +#include "FWCore/Framework/interface/stream/EDProducer.h" +#include "FWCore/Framework/interface/Event.h" +#include "FWCore/Framework/interface/EventSetup.h" +#include "FWCore/Utilities/interface/InputTag.h" +#include "FWCore/ParameterSet/interface/ParameterSet.h" +#include "FWCore/ParameterSet/interface/ParameterSetDescription.h" +#include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h" + +#include "DataFormats/NanoAOD/interface/FlatTable.h" +#include "DataFormats/PatCandidates/interface/Electron.h" +#include "DataFormats/PatCandidates/interface/Muon.h" +#include "DataFormats/PatCandidates/interface/Tau.h" +#include "DataFormats/PatCandidates/interface/PackedCandidate.h" +#include "DataFormats/VertexReco/interface/Vertex.h" +#include "DataFormats/TrackReco/interface/Track.h" + +#include "MagneticField/Engine/interface/MagneticField.h" +#include "RecoVertex/VertexPrimitives/interface/TransientVertex.h" +#include "RecoVertex/VertexTools/interface/VertexDistance3D.h" +#include "RecoVertex/VertexPrimitives/interface/ConvertToFromReco.h" +#include "RecoVertex/KalmanVertexFit/interface/KalmanVertexFitter.h" +#include "TrackingTools/TransientTrack/interface/TransientTrackBuilder.h" +#include "TrackingTools/Records/interface/TransientTrackRecord.h" +#include "TrackingTools/GeomPropagators/interface/AnalyticalTrajectoryExtrapolatorToLine.h" +#include "TrackingTools/GeomPropagators/interface/AnalyticalImpactPointExtrapolator.h" +#include "CommonTools/UtilAlgos/interface/StringCutObjectSelector.h" +#include "Utilities/General/interface/ClassName.h" +#include "PhysicsTools/NanoAOD/interface/SimpleFlatTableProducer.h" +#include "DataFormats/NanoAOD/interface/LeptonTimeLifeInfo.h" + +#include + +namespace { + typedef FuncVariable, int32_t> IntLeptonTimeLifeInfoVar; + typedef FuncVariable, uint32_t> UIntLeptonTimeLifeInfoVar; + typedef FuncVariable, float> FloatLeptonTimeLifeInfoVar; + typedef FuncVariable, double> DoubleLeptonTimeLifeInfoVar; + typedef FuncVariable, int8_t> Int8LeptonTimeLifeInfoVar; + typedef FuncVariable, uint8_t> UInt8LeptonTimeLifeInfoVar; + typedef FuncVariable, int16_t> Int16LeptonTimeLifeInfoVar; + typedef FuncVariable, uint16_t> + UInt16LeptonTimeLifeInfoVar; + typedef FuncVariable, bool> BoolLeptonTimeLifeInfoVar; +} // namespace + +template +class LeptonTimeLifeInfoTableProducer : public edm::stream::EDProducer<> { +public: + explicit LeptonTimeLifeInfoTableProducer(const edm::ParameterSet&); + ~LeptonTimeLifeInfoTableProducer() override{}; + + static void fillDescriptions(edm::ConfigurationDescriptions& descriptions); + void produce(edm::Event&, const edm::EventSetup&) override; + +private: + //--- private utility methods + const reco::Track* getTrack(const T&); + void produceAndFillIPInfo(const T&, const TransientTrackBuilder&, const reco::Vertex&, LeptonTimeLifeInfo&); + void produceAndFillSVInfo(const T&, const TransientTrackBuilder&, const reco::Vertex&, LeptonTimeLifeInfo&); + static bool fitVertex(const std::vector& transTrk, TransientVertex& transVtx) { + if (transTrk.size() < 2) + return false; + KalmanVertexFitter kvf(true); + transVtx = kvf.vertex(transTrk); + return transVtx.hasRefittedTracks() && transVtx.refittedTracks().size() == transTrk.size(); + } + + //--- configuration parameters + edm::EDGetTokenT> leptonsToken_; + edm::EDGetTokenT pvToken_; + edm::ESGetToken transTrackBuilderToken_; + std::string name_, doc_; + bool extension_; + const StringCutObjectSelector selector_; + int pvChoice_; + std::vector>> vars_; + + enum PVChoice { useFront = 0, useClosestInDz }; +}; + +template +LeptonTimeLifeInfoTableProducer::LeptonTimeLifeInfoTableProducer(const edm::ParameterSet& cfg) + : leptonsToken_(consumes>(cfg.getParameter("src"))), + pvToken_(consumes(cfg.getParameter("pvSource"))), + transTrackBuilderToken_(esConsumes(edm::ESInputTag("", "TransientTrackBuilder"))), + name_(cfg.getParameter("name")), + doc_(cfg.getParameter("doc")), + extension_(cfg.getParameter("extension")), + selector_(cfg.getParameter("selection")), + pvChoice_(cfg.getParameter("pvChoice")) { + //Time-life info variables + edm::ParameterSet const& varsPSet = cfg.getParameter("variables"); + for (const std::string& vname : varsPSet.getParameterNamesForType()) { + const auto& varPSet = varsPSet.getParameter(vname); + const std::string& type = varPSet.getParameter("type"); + if (type == "int") + vars_.push_back(std::make_unique(vname, varPSet)); + else if (type == "uint") + vars_.push_back(std::make_unique(vname, varPSet)); + else if (type == "float") + vars_.push_back(std::make_unique(vname, varPSet)); + else if (type == "double") + vars_.push_back(std::make_unique(vname, varPSet)); + else if (type == "int8") + vars_.push_back(std::make_unique(vname, varPSet)); + else if (type == "uint8") + vars_.push_back(std::make_unique(vname, varPSet)); + else if (type == "int16") + vars_.push_back(std::make_unique(vname, varPSet)); + else if (type == "uint16") + vars_.push_back(std::make_unique(vname, varPSet)); + else if (type == "bool") + vars_.push_back(std::make_unique(vname, varPSet)); + else + throw cms::Exception("Configuration", "unsupported type " + type + " for variable " + vname); + } + produces(); +} + +template +void LeptonTimeLifeInfoTableProducer::produce(edm::Event& evt, const edm::EventSetup& es) { + // Get leptons + edm::Handle> leptons; + evt.getByToken(leptonsToken_, leptons); + + // Get the vertices + edm::Handle vertices; + evt.getByToken(pvToken_, vertices); + + // Get transient track builder + const TransientTrackBuilder& transTrackBuilder = es.getData(transTrackBuilderToken_); + + std::vector infos; + infos.reserve(leptons->size()); + + for (const auto& lepton : *leptons) { + LeptonTimeLifeInfo* info = new LeptonTimeLifeInfo(); + + // Do nothing for lepton not passing selection + if (!selector_(lepton)) { + infos.push_back(info); + continue; + } + size_t pv_idx = 0; + if (pvChoice_ == useClosestInDz && getTrack(lepton) != nullptr) { + float dz_min = 999; + size_t vtx_idx = 0; + for (const auto& vtx : *vertices) { + float dz_tmp = std::abs(getTrack(lepton)->dz(vtx.position())); + if (dz_tmp < dz_min) { + dz_min = dz_tmp; + pv_idx = vtx_idx; + } + vtx_idx++; + } + } + const reco::Vertex& pv = !vertices->empty() ? (*vertices)[pv_idx] : reco::Vertex(); + + // Obtain IP vector and set related info into lepton + produceAndFillIPInfo(lepton, transTrackBuilder, pv, *info); + + // Fit SV and set related info for taus or do nothing for other lepton types + produceAndFillSVInfo(lepton, transTrackBuilder, pv, *info); + infos.push_back(info); + } // end of lepton loop + + // Define and fill table + auto timelifeTable = std::make_unique(leptons->size(), name_, false, extension_); + timelifeTable->setDoc(doc_); + for (const auto& var : vars_) + var->fill(infos, *timelifeTable); + for (auto info : infos) + delete info; + + // Store table in the event + evt.put(std::move(timelifeTable)); +} + +template <> +const reco::Track* LeptonTimeLifeInfoTableProducer::getTrack(const pat::Electron& electron) { + return electron.gsfTrack().isNonnull() ? electron.gsfTrack().get() : nullptr; +} + +template <> +const reco::Track* LeptonTimeLifeInfoTableProducer::getTrack(const pat::Muon& muon) { + return muon.innerTrack().isNonnull() ? muon.innerTrack().get() : nullptr; +} + +template <> +const reco::Track* LeptonTimeLifeInfoTableProducer::getTrack(const pat::Tau& tau) { + const reco::Track* track = nullptr; + if (tau.leadChargedHadrCand().isNonnull()) + track = tau.leadChargedHadrCand()->bestTrack(); + return track; +} + +template +void LeptonTimeLifeInfoTableProducer::produceAndFillIPInfo(const T& lepton, + const TransientTrackBuilder& transTrackBuilder, + const reco::Vertex& pv, + LeptonTimeLifeInfo& info) { + const reco::Track* track = getTrack(lepton); + if (track != nullptr) { + info.setTrack(track); + info.setBFiled_z(transTrackBuilder.field()->inInverseGeV(GlobalPoint(track->vx(), track->vy(), track->vz())).z()); + + // Extrapolate track to the point closest to PV + reco::TransientTrack transTrack = transTrackBuilder.build(track); + AnalyticalImpactPointExtrapolator extrapolator(transTrack.field()); + TrajectoryStateOnSurface closestState = + extrapolator.extrapolate(transTrack.impactPointState(), RecoVertex::convertPos(pv.position())); + GlobalPoint pca = closestState.globalPosition(); + GlobalError pca_cov = closestState.cartesianError().position(); + GlobalVector ip_vec = GlobalVector(pca.x() - pv.x(), pca.y() - pv.y(), pca.z() - pv.z()); + GlobalError ip_cov = pca_cov + GlobalError(pv.covariance()); + VertexDistance3D pca_dist; + Measurement1D ip_mes = pca_dist.distance(pv, VertexState(pca, pca_cov)); + if (ip_vec.dot(GlobalVector(lepton.px(), lepton.py(), lepton.pz())) < 0) + ip_mes = Measurement1D(-1. * ip_mes.value(), ip_mes.error()); + + // Store PCA info + info.setPCA(pca, pca_cov); + info.setIP(ip_vec, ip_cov); + info.setIPLength(ip_mes); + } +} + +template +void LeptonTimeLifeInfoTableProducer::produceAndFillSVInfo(const T& lepton, + const TransientTrackBuilder& transTrackBuilder, + const reco::Vertex& pv, + LeptonTimeLifeInfo& info) {} + +template <> +void LeptonTimeLifeInfoTableProducer::produceAndFillSVInfo(const pat::Tau& tau, + const TransientTrackBuilder& transTrackBuilder, + const reco::Vertex& pv, + LeptonTimeLifeInfo& info) { + // Fit SV with tracks of charged tau decay products + int fitOK = 0; + if (tau.signalChargedHadrCands().size() + tau.signalLostTracks().size() > 1) { + // Get tracks from tau signal charged candidates + std::vector transTrks; + TransientVertex transVtx; + for (const auto& cand : tau.signalChargedHadrCands()) { + if (cand.isNull()) + continue; + const reco::Track* track = cand->bestTrack(); + if (track != nullptr) + transTrks.push_back(transTrackBuilder.build(track)); + } + for (const auto& cand : tau.signalLostTracks()) { + if (cand.isNull()) + continue; + const reco::Track* track = cand->bestTrack(); + if (track != nullptr) + transTrks.push_back(transTrackBuilder.build(track)); + } + // Fit SV with KalmanVertexFitter + fitOK = fitVertex(transTrks, transVtx) ? 1 : -1; + if (fitOK > 0) { + reco::Vertex sv = transVtx; + // Get flight-length + // Full PV->SV flight vector with its covariance + GlobalVector flight_vec = GlobalVector(sv.x() - pv.x(), sv.y() - pv.y(), sv.z() - pv.z()); + GlobalError flight_cov = transVtx.positionError() + GlobalError(pv.covariance()); + //MB: can be taken from tau itself (but with different fit of PV and SV) as follows: + //tau.flightLength().mag2()); + //tau.flightLengthSig(); + VertexDistance3D sv_dist; + Measurement1D flightLength_mes = sv_dist.signedDistance(pv, sv, GlobalVector(tau.px(), tau.py(), tau.pz())); + + // Store SV info + info.setSV(sv); + info.setFlightVector(flight_vec, flight_cov); + info.setFlightLength(flightLength_mes); + } + } +} + +template +void LeptonTimeLifeInfoTableProducer::fillDescriptions(edm::ConfigurationDescriptions& descriptions) { + // pat{Electron,Muon,Tau}TimeLifeInfoUpdater + edm::ParameterSetDescription desc; + + std::string lepCollName; + if (typeid(T) == typeid(pat::Electron)) + lepCollName = "slimmedElectrons"; + else if (typeid(T) == typeid(pat::Muon)) + lepCollName = "slimmedMuons"; + else if (typeid(T) == typeid(pat::Tau)) + lepCollName = "slimmedTaus"; + desc.add("src", edm::InputTag(lepCollName)); + desc.add("pvSource", edm::InputTag("offlineSlimmedPrimaryVertices")); + desc.add("selection", "")->setComment("Selection required to produce and store time-life information"); + desc.add("name")->setComment("Name of the branch in the time-life info table for " + + ClassName::name()); + desc.add("doc", "")->setComment("Few words of self documentation"); + desc.add("extension", true)->setComment("Whether or not to extend an existing same table (default: true)"); + desc.add("pvChoice", useFront) + ->setComment( + "Define PV to compute IP: 0: first PV, 1: PV with the smallest dz of the tau leading track (default: " + + std::to_string(useFront) + ")"); + + // variables + edm::ParameterSetDescription variable; + variable.add("expr")->setComment("a function to define the content of the branch in the flat table"); + variable.add("doc")->setComment("few words description of the branch content"); + variable.ifValue( + edm::ParameterDescription( + "type", "int", true, edm::Comment("the c++ type of the branch in the flat table")), + edm::allowedValues("int", "uint", "float", "double", "int8", "uint8", "int16", "uint16", "bool")); + variable.addOptionalNode( + edm::ParameterDescription( + "precision", true, edm::Comment("the precision with which to store the value in the flat table")) xor + edm::ParameterDescription( + "precision", true, edm::Comment("the precision with which to store the value in the flat table")), + false); + edm::ParameterSetDescription variables; + variables.setComment("a parameters set to define additional variables describing main vertex"); + variables.addNode(edm::ParameterWildcard("*", edm::RequireZeroOrMore, true, variable)); + desc.add("variables", variables)->setComment("Time-life info variables"); + + descriptions.addWithDefaultLabel(desc); +} + +#include "FWCore/Framework/interface/MakerMacros.h" +typedef LeptonTimeLifeInfoTableProducer ElectronTimeLifeInfoTableProducer; +DEFINE_FWK_MODULE(ElectronTimeLifeInfoTableProducer); +typedef LeptonTimeLifeInfoTableProducer MuonTimeLifeInfoTableProducer; +DEFINE_FWK_MODULE(MuonTimeLifeInfoTableProducer); +typedef LeptonTimeLifeInfoTableProducer TauTimeLifeInfoTableProducer; +DEFINE_FWK_MODULE(TauTimeLifeInfoTableProducer); diff --git a/PhysicsTools/NanoAOD/python/electrons_cff.py b/PhysicsTools/NanoAOD/python/electrons_cff.py index 021092bbb7062..4297a67da4f15 100644 --- a/PhysicsTools/NanoAOD/python/electrons_cff.py +++ b/PhysicsTools/NanoAOD/python/electrons_cff.py @@ -394,20 +394,17 @@ def _get_bitmapVIDForEle_docstring(modules,WorkingPoints): ########## electronTable extension defn ########## # Produce and add time-life info (e.g. for electrons from tau decays) -from PhysicsTools.PatAlgos.patElectronTimeLifeInfoUpdater_cfi import patElectronTimeLifeInfoUpdater +from PhysicsTools.NanoAOD.electronTimeLifeInfoTableProducer_cfi import electronTimeLifeInfoTableProducer import PhysicsTools.NanoAOD.leptonTimeLifeInfo_common_cff as tli from TrackingTools.TransientTrack.TransientTrackBuilder_cfi import * -electronsWithTimeLifeInfo = patElectronTimeLifeInfoUpdater.clone( - src = electronTable.src, - pvSource = tli.prod_common.pvSource, - pvChoice = tli.prod_common.pvChoice -) - -electronTimeLifeInfoTable = simpleCandidateFlatTableProducer.clone( - src = "electronsWithTimeLifeInfo", +electronTimeLifeInfoTable = electronTimeLifeInfoTableProducer.clone( name = electronTable.name, + src = electronTable.src, + selection = 'pt > 20', doc = cms.string("Additional time-life info for non-prompt electrons"), extension = True, + pvSource = tli.prod_common.pvSource, + pvChoice = tli.prod_common.pvChoice, variables = cms.PSet( tli.ipVars, tli.trackVars @@ -465,7 +462,7 @@ def _get_bitmapVIDForEle_docstring(modules,WorkingPoints): electronTask = cms.Task(bitmapVIDForEle,bitmapVIDForEleFall17V2,bitmapVIDForEleHEEP,isoForEle,isoForEleFall17V2,ptRatioRelForEle,seedGainEle,calibratedPatElectronsNano,slimmedElectronsWithUserData,finalElectrons) electronTablesTask = cms.Task(electronMVATTH, electronTable) electronMCTask = cms.Task(tautaggerForMatching, matchingElecPhoton, electronsMCMatchForTable, electronsMCMatchForTableAlt, electronMCTable) -electronTimeLifeInfoTask = cms.Task(electronsWithTimeLifeInfo,electronTimeLifeInfoTable) +electronTimeLifeInfoTask = cms.Task(electronTimeLifeInfoTable) electronTablesTask.add(electronTimeLifeInfoTask) _electronTask_Run2 = electronTask.copy() diff --git a/PhysicsTools/NanoAOD/python/leptonTimeLifeInfo_common_cff.py b/PhysicsTools/NanoAOD/python/leptonTimeLifeInfo_common_cff.py index 35c85bcf43482..dd3696bb17b3d 100644 --- a/PhysicsTools/NanoAOD/python/leptonTimeLifeInfo_common_cff.py +++ b/PhysicsTools/NanoAOD/python/leptonTimeLifeInfo_common_cff.py @@ -1,6 +1,6 @@ # -# Common definition of time-life variables for pat-leptons updated -# with pat{Electron,Muon,Tau}TimeLifeInfoUpdater +# Common definition of time-life variables for pat-leptons produced +# with {Electron,Muon,Tau}TimeLifeInfoTableProducer # import FWCore.ParameterSet.Config as cms from PhysicsTools.NanoAOD.common_cff import * @@ -15,47 +15,47 @@ # impact parameter ipVars = cms.PSet( - ipLength = Var("?hasUserFloat('IP')?userFloat('IP'):0", float, doc="lenght of impact parameter (3d)", precision=16), - ipLengthSig = Var("?hasUserFloat('IP_sig')?userFloat('IP_sig'):0", float, doc="Significance of impact parameter", precision=16), - IPx = Var("?hasUserFloat('IP_x')?userFloat('IP_x'):0", float, doc="x coordinate of impact parameter vector", precision=16), - IPy = Var("?hasUserFloat('IP_x')?userFloat('IP_x'):0", float, doc="x coordinate of impact parameter vector", precision=16), - IPz = Var("?hasUserFloat('IP_x')?userFloat('IP_x'):0", float, doc="x coordinate of impact parameter vector", precision=16) + ipLength = Var("ipLength().value()", float, doc="lenght of impact parameter (3d)", precision=10), + ipLengthSig = Var("ipLength().significance()", float, doc="significance of impact parameter", precision=10), + IPx = Var("ipVector().x()", float, doc="x coordinate of impact parameter vector", precision=10), + IPy = Var("ipVector().y()", float, doc="y coordinate of impact parameter vector", precision=10), + IPz = Var("ipVector().z()", float, doc="z coordinate of impact parameter vector", precision=10) ) # track parameters and covariance at ref. point trackVars = cms.PSet( - track_qoverp = Var("?hasUserFloat('track_qoverp')?userFloat('track_qoverp'):0", float, doc="track q/p", precision=16), - track_lambda = Var("?hasUserFloat('track_lambda')?userFloat('track_lambda'):0", float, doc="track lambda", precision=16), - track_phi = Var("?hasUserFloat('track_phi')?userFloat('track_phi'):0", float, doc="track phi", precision=16), - track_dxy = Var("?hasUserFloat('track_dxy')?userFloat('track_dxy'):0", float, doc="track dxy", precision=16), - track_dsz = Var("?hasUserFloat('track_dsz')?userFloat('track_dsz'):0", float, doc="track dsz", precision=16), - bField_z = Var("?hasUserFloat('bField_z')?userFloat('bField_z'):0", float, doc="z coordinate of magnetic field at track ref. point", precision=16), + track_qoverp = Var("?hasTrack()?track().parameter(0):0", float, doc="track q/p", precision=10), + track_lambda = Var("?hasTrack()?track().parameter(1):0", float, doc="track lambda", precision=10), + track_phi = Var("?hasTrack()?track().parameter(2):0", float, doc="track phi", precision=10), + #track_deltaPhi = Var("?hasTrack()?deltaPhi(track().parameter(2), phi):0", float, doc="track phi minus lepton phi", precision=10), + track_dxy = Var("?hasTrack()?track().parameter(3):0", float, doc="track dxy", precision=10), + track_dsz = Var("?hasTrack()?track().parameter(4):0", float, doc="track dsz", precision=10), + bField_z = Var("?hasTrack()?bField_z:0", float, doc="z coordinate of magnetic field at track ref. point", precision=10), ) # track covariance elements (adding to trackVars) for i in range(0,5): for j in range(i,5): jistr = str(j)+str(i) - setattr(trackVars, 'track_cov'+jistr,Var("?hasUserFloat('trackCov_"+jistr+"')?userFloat('trackCov_"+jistr+"'):0", float, doc="track covariance element ("+str(j)+","+str(i)+")", precision=16)) + setattr(trackVars, 'track_cov'+jistr, Var("?hasTrack()?track().covariance("+str(j)+","+str(i)+"):0", float, doc="track covariance element ("+str(j)+","+str(i)+")", precision=10)) # secondary vertex svVars = cms.PSet( # SV - SVstatus = Var("?hasUserInt('refitSV_status')?userInt('refitSV_status'):-1", int, doc="Status of SV refit: 0: not fitted or not present, 1: OK, -1: failed"), - SVx = Var("?hasUserFloat('refitSV_x')?userFloat('refitSV_x'):0", float, doc="x coordinate of SV", precision=16), - SVy = Var("?hasUserFloat('refitSV_y')?userFloat('refitSV_y'):0", float, doc="y coordinate of SV", precision=16), - SVz = Var("?hasUserFloat('refitSV_z')?userFloat('refitSV_z'):0", float, doc="z coordinate of SV", precision=16), - SVchi2 = Var("?hasUserFloat('refitSV_chi2')?userFloat('refitSV_chi2'):0", float, doc="chi2 of SV fit", precision=10), - SVndof = Var("?hasUserFloat('refitSV_ndof')?userFloat('refitSV_ndof'):0", float, doc="ndof of SV fit", precision=10), - SVcxx = Var("?hasUserFloat('refitSV_cxx')?userFloat('refitSV_cxx'):0", float, doc="Covariance of SV (0,0)", precision=16), - SVcyx = Var("?hasUserFloat('refitSV_cyx')?userFloat('refitSV_cyx'):0", float, doc="Covariance of SV (1,0)", precision=16), - SVczx = Var("?hasUserFloat('refitSV_czx')?userFloat('refitSV_czx'):0", float, doc="Covariance of SV (2,0)", precision=16), - SVcyy = Var("?hasUserFloat('refitSV_cyy')?userFloat('refitSV_cyy'):0", float, doc="Covariance of SV (1,1)", precision=16), - SVczy = Var("?hasUserFloat('refitSV_czy')?userFloat('refitSV_czy'):0", float, doc="Covariance of SV (2,1)", precision=16), - SVczz = Var("?hasUserFloat('refitSV_czz')?userFloat('refitSV_czz'):0", float, doc="Covariance of SV (2,2)", precision=16), + hasRefitSV = Var("hasSV()", bool, doc="has SV refit using miniAOD quantities"), + refitSVx = Var("?hasSV()?sv().x():0", float, doc="x coordinate of SV", precision=10), + refitSVy = Var("?hasSV()?sv().y():0", float, doc="y coordinate of SV", precision=10), + refitSVz = Var("?hasSV()?sv().z():0", float, doc="z coordinate of SV", precision=10), + refitSVchi2 = Var("?hasSV()?sv().chi2():0", float, doc="chi2 of SV fit", precision=8), + refitSVndof = Var("?hasSV()?sv().ndof():0", float, doc="ndof of SV fit", precision=8), # flight-length - flightLength = Var("?hasUserFloat('flightLength')?userFloat('flightLength'):0", float, doc="flight-length,i.e. the PV to SV distance", precision=16), - flightLengthSig = Var("?hasUserFloat('flightLength_sig')?userFloat('flightLength_sig'):0", float, doc="Significance of flight-length", precision=16) + #refitFlightLength = Var("?hasSV()?flightLength().value():0", float, doc="flight-length,i.e. the PV to SV distance", precision=10), + #refitFlightLengthSig = Var("?hasSV()?flightLength().significance():0", float, doc="Significance of flight-length", precision=10) ) +# secondary vertex covariance elements (adding to svVars) +for i in range(0,3): + for j in range(i,3): + jistr = str(j)+str(i) + setattr(svVars, 'refitSVcov'+jistr, Var("?hasSV()?sv().covariance("+str(j)+","+str(i)+"):0", float, doc="Covariance of SV ("+str(j)+","+str(i)+")", precision=10)) # Module to refit PV with beam-spot constraint that is not present in Run-2 samples refittedPV = miniAODRefitVertexProducer.clone( diff --git a/PhysicsTools/NanoAOD/python/muons_cff.py b/PhysicsTools/NanoAOD/python/muons_cff.py index 70084a93e0017..8a066f1a93cb5 100644 --- a/PhysicsTools/NanoAOD/python/muons_cff.py +++ b/PhysicsTools/NanoAOD/python/muons_cff.py @@ -5,7 +5,7 @@ from PhysicsTools.NanoAOD.simpleCandidateFlatTableProducer_cfi import simpleCandidateFlatTableProducer import PhysicsTools.PatAlgos.producersLayer1.muonProducer_cfi -from PhysicsTools.PatAlgos.patMuonTimeLifeInfoUpdater_cfi import patMuonTimeLifeInfoUpdater +from PhysicsTools.NanoAOD.muonTimeLifeInfoTableProducer_cfi import muonTimeLifeInfoTableProducer import PhysicsTools.NanoAOD.leptonTimeLifeInfo_common_cff as tli # this below is used only in some eras @@ -189,17 +189,14 @@ ) # Produce and add time-life info (e.g. for muons from tau decays) -muonsWithTimeLifeInfo = patMuonTimeLifeInfoUpdater.clone( - src = muonTable.src, - pvSource = tli.prod_common.pvSource, - pvChoice = tli.prod_common.pvChoice -) - -muonTimeLifeInfoTable = simpleCandidateFlatTableProducer.clone( - src = "muonsWithTimeLifeInfo", +muonTimeLifeInfoTable = muonTimeLifeInfoTableProducer.clone( name = muonTable.name, + src = muonTable.src, + selection = 'pt > 15', doc = cms.string("Additional time-life info for non-prompt muons"), extension = True, + pvSource = tli.prod_common.pvSource, + pvChoice = tli.prod_common.pvChoice, variables = cms.PSet( tli.ipVars, tli.trackVars @@ -244,7 +241,7 @@ muonTask = cms.Task(slimmedMuonsUpdated,isoForMu,ptRatioRelForMu,slimmedMuonsWithUserData,finalMuons,finalLooseMuons ) muonMCTask = cms.Task(muonsMCMatchForTable,muonMCTable) muonTablesTask = cms.Task(muonMVATTH,muonMVALowPt,muonBSConstrain,muonTable,muonMVAID) -muonTimeLifeInfoTask = cms.Task(muonsWithTimeLifeInfo,muonTimeLifeInfoTable) +muonTimeLifeInfoTask = cms.Task(muonTimeLifeInfoTable) # Refit PV with beam-spot constraint that is not present in Run-2 samples _muonTimeLifeInfoTaskRun2 = muonTimeLifeInfoTask.copy() _muonTimeLifeInfoTaskRun2.add(tli.refittedPV) diff --git a/PhysicsTools/NanoAOD/python/taus_cff.py b/PhysicsTools/NanoAOD/python/taus_cff.py index c2a7da94a5daa..33621940ffe0e 100644 --- a/PhysicsTools/NanoAOD/python/taus_cff.py +++ b/PhysicsTools/NanoAOD/python/taus_cff.py @@ -7,7 +7,7 @@ from PhysicsTools.JetMCAlgos.TauGenJetsDecayModeSelectorAllHadrons_cfi import tauGenJetsSelectorAllHadrons from PhysicsTools.PatAlgos.patTauSignalCandidatesProducer_cfi import patTauSignalCandidatesProducer -from PhysicsTools.PatAlgos.patTauTimeLifeInfoUpdater_cfi import patTauTimeLifeInfoUpdater +from PhysicsTools.NanoAOD.tauTimeLifeInfoTableProducer_cfi import tauTimeLifeInfoTableProducer import PhysicsTools.NanoAOD.leptonTimeLifeInfo_common_cff as tli ##################### Updated tau collection with MVA-based tau-Ids rerun ####### @@ -176,17 +176,13 @@ def _tauIdWPMask(pattern, choices, doc="", from_raw=False, wp_thrs=None): # Produce and add time-life info from TrackingTools.TransientTrack.TransientTrackBuilder_cfi import * -tausWithTimeLifeInfo = patTauTimeLifeInfoUpdater.clone( - src = tauTable.src, - pvSource = tli.prod_common.pvSource, - pvChoice = tli.prod_common.pvChoice -) - -tauTimeLifeInfoTable = simpleCandidateFlatTableProducer.clone( - src = "tausWithTimeLifeInfo", +tauTimeLifeInfoTable = tauTimeLifeInfoTableProducer.clone( name = tauTable.name, + src = tauTable.src, doc = cms.string("Additional tau time-life info"), extension = True, + pvSource = tli.prod_common.pvSource, + pvChoice = tli.prod_common.pvChoice, variables = cms.PSet( tli.svVars, tli.ipVars, @@ -263,7 +259,7 @@ def _tauIdWPMask(pattern, choices, doc="", from_raw=False, wp_thrs=None): tauTablesTask = cms.Task(tauTable) tauSignalCandsTask = cms.Task(tauSignalCands,tauSignalCandsTable) tauTablesTask.add(tauSignalCandsTask) -tauTimeLifeInfoTask = cms.Task(tausWithTimeLifeInfo,tauTimeLifeInfoTable) +tauTimeLifeInfoTask = cms.Task(tauTimeLifeInfoTable) # Refit PV with beam-spot constraint that is not present in Run-2 samples _tauTimeLifeInfoTaskRun2 = tauTimeLifeInfoTask.copy() _tauTimeLifeInfoTaskRun2.add(tli.refittedPV) diff --git a/PhysicsTools/NanoAOD/python/vertices_cff.py b/PhysicsTools/NanoAOD/python/vertices_cff.py index b0e5965ed915f..54569f753a32d 100644 --- a/PhysicsTools/NanoAOD/python/vertices_cff.py +++ b/PhysicsTools/NanoAOD/python/vertices_cff.py @@ -39,15 +39,15 @@ svCandidateTable.variables.phi.precision=12 covarianceVars = cms.PSet( - cxx = Var("covariance(0,0)", float, doc="vertex covariance (0,0)", precision = 16), - cyx = Var("covariance(1,0)", float, doc="vertex covariance (1,0)", precision = 16), - czx = Var("covariance(2,0)", float, doc="vertex covariance (2,0)", precision = 16), - cyy = Var("covariance(1,1)", float, doc="vertex covariance (1,1)", precision = 16), - czy = Var("covariance(2,1)", float, doc="vertex covariance (2,1)", precision = 16), - czz = Var("covariance(2,2)", float, doc="vertex covariance (2,2)", precision = 16) + cxx = Var("covariance(0,0)", float, doc="vertex covariance (0,0)", precision = 10), + cyx = Var("covariance(1,0)", float, doc="vertex covariance (1,0)", precision = 10), + czx = Var("covariance(2,0)", float, doc="vertex covariance (2,0)", precision = 10), + cyy = Var("covariance(1,1)", float, doc="vertex covariance (1,1)", precision = 10), + czy = Var("covariance(2,1)", float, doc="vertex covariance (2,1)", precision = 10), + czz = Var("covariance(2,2)", float, doc="vertex covariance (2,2)", precision = 10) ) -vertexTable.optionalPvVariables = covarianceVars +#vertexTable.optionalPvVariables = covarianceVars verticesWithBS = prod_common.pvSource pvbsTable = simpleVertexFlatTableProducer.clone( diff --git a/PhysicsTools/PatAlgos/plugins/PATLeptonTimeLifeInfoUpdater.cc b/PhysicsTools/PatAlgos/plugins/PATLeptonTimeLifeInfoUpdater.cc deleted file mode 100644 index 1bb901c2cdb36..0000000000000 --- a/PhysicsTools/PatAlgos/plugins/PATLeptonTimeLifeInfoUpdater.cc +++ /dev/null @@ -1,331 +0,0 @@ -/** - \class PATLeptonTimeLifeInfoUpdater - \brief Produces updated collection of pat::Electron's or pat::Muon's - - The PATLeptonInfoUpdater produces analysis-level pat::Electron's or - pat::Muon's with updated time-life information. - - \author Michal Bluj, NCBJ, Warsaw -*/ - -#include "FWCore/Framework/interface/stream/EDProducer.h" -#include "FWCore/Framework/interface/Event.h" -#include "FWCore/Framework/interface/EventSetup.h" -#include "FWCore/Utilities/interface/InputTag.h" -#include "FWCore/ParameterSet/interface/ParameterSet.h" -#include "FWCore/ParameterSet/interface/ParameterSetDescription.h" -#include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h" - -#include "DataFormats/PatCandidates/interface/Electron.h" -#include "DataFormats/PatCandidates/interface/Muon.h" -#include "DataFormats/PatCandidates/interface/Tau.h" -#include "DataFormats/PatCandidates/interface/PackedCandidate.h" -#include "DataFormats/VertexReco/interface/Vertex.h" -#include "DataFormats/TrackReco/interface/Track.h" - -#include "MagneticField/Engine/interface/MagneticField.h" -#include "RecoVertex/VertexPrimitives/interface/TransientVertex.h" -#include "RecoVertex/VertexTools/interface/VertexDistance3D.h" -#include "RecoVertex/VertexPrimitives/interface/ConvertToFromReco.h" -#include "RecoVertex/KalmanVertexFit/interface/KalmanVertexFitter.h" -#include "TrackingTools/TransientTrack/interface/TransientTrackBuilder.h" -#include "TrackingTools/Records/interface/TransientTrackRecord.h" -#include "TrackingTools/GeomPropagators/interface/AnalyticalTrajectoryExtrapolatorToLine.h" -#include "TrackingTools/GeomPropagators/interface/AnalyticalImpactPointExtrapolator.h" - -#include - -template -class PATLeptonTimeLifeInfoUpdater : public edm::stream::EDProducer<> { -public: - explicit PATLeptonTimeLifeInfoUpdater(const edm::ParameterSet&); - ~PATLeptonTimeLifeInfoUpdater() override{}; - - static void fillDescriptions(edm::ConfigurationDescriptions& descriptions); - void produce(edm::Event&, const edm::EventSetup&) override; - -private: - //--- private utility methods - const reco::Track* getTrack(const T&); - void produceAndFillIPInfo(T&, const TransientTrackBuilder&, const reco::Vertex&); - void produceAndFillSVInfo(T&, const TransientTrackBuilder&, const reco::Vertex&); - - //--- configuration parameters - edm::EDGetTokenT> leptonsToken_; - edm::EDGetTokenT pvToken_; - edm::ESGetToken transTrackBuilderToken_; - int pvChoice_; - - enum PVChoice { useFront = 0, useClosestInDz }; -}; - -template -PATLeptonTimeLifeInfoUpdater::PATLeptonTimeLifeInfoUpdater(const edm::ParameterSet& cfg) - : leptonsToken_(consumes>(cfg.getParameter("src"))), - pvToken_(consumes(cfg.getParameter("pvSource"))), - transTrackBuilderToken_(esConsumes(edm::ESInputTag("", "TransientTrackBuilder"))), - pvChoice_(cfg.getParameter("pvChoice")) { - produces>(); -} - -template -void PATLeptonTimeLifeInfoUpdater::produce(edm::Event& evt, const edm::EventSetup& es) { - // Get leptons - edm::Handle> inputLeptons; - evt.getByToken(leptonsToken_, inputLeptons); - - auto outputLeptons = std::make_unique>(); - outputLeptons->reserve(inputLeptons->size()); - - // Get the vertices - edm::Handle vertices; - evt.getByToken(pvToken_, vertices); - - // Get transient track builder - const TransientTrackBuilder& transTrackBuilder = es.getData(transTrackBuilderToken_); - - for (const auto& inputLepton : *inputLeptons) { - T outputLepton(inputLepton); - size_t pv_idx = 0; - if (pvChoice_ == useClosestInDz && getTrack(inputLepton) != nullptr) { - float dz_min = 999; - size_t vtx_idx = 0; - for (const auto& vtx : *vertices) { - float dz_tmp = std::abs(getTrack(inputLepton)->dz(vtx.position())); - if (dz_tmp < dz_min) { - dz_min = dz_tmp; - pv_idx = vtx_idx; - } - vtx_idx++; - } - } - const reco::Vertex& pv = !vertices->empty() ? (*vertices)[pv_idx] : reco::Vertex(); - // Obtain IP vector and set related info into lepton - produceAndFillIPInfo(outputLepton, transTrackBuilder, pv); - - // Fit SV and set related info in taus or do nothing for other lepton types - produceAndFillSVInfo(outputLepton, transTrackBuilder, pv); - - // Set updated lepton to output collection - outputLeptons->push_back(outputLepton); - } // end of lepton loop - - evt.put(std::move(outputLeptons)); -} - -template <> -const reco::Track* PATLeptonTimeLifeInfoUpdater::getTrack(const pat::Electron& electron) { - return electron.gsfTrack().isNonnull() ? electron.gsfTrack().get() : nullptr; -} - -template <> -const reco::Track* PATLeptonTimeLifeInfoUpdater::getTrack(const pat::Muon& muon) { - return muon.innerTrack().isNonnull() ? muon.innerTrack().get() : nullptr; -} - -template <> -const reco::Track* PATLeptonTimeLifeInfoUpdater::getTrack(const pat::Tau& tau) { - const reco::Track* track = nullptr; - if (tau.leadChargedHadrCand().isNonnull()) - track = tau.leadChargedHadrCand()->bestTrack(); - return track; -} - -template -void PATLeptonTimeLifeInfoUpdater::produceAndFillIPInfo(T& lepton, - const TransientTrackBuilder& transTrackBuilder, - const reco::Vertex& pv) { - GlobalPoint pca; - GlobalError pca_cov; - GlobalVector ip_vec; - GlobalError ip_cov; - float ipLength = -1; - float ipLength_sig = -1; - float bField_z = 0; - const reco::Track* track = getTrack(lepton); - if (track != nullptr) { - bField_z = transTrackBuilder.field()->inInverseGeV(GlobalPoint(track->vx(), track->vy(), track->vz())).z(); - // Extrapolate track to the point closest to PV - reco::TransientTrack transTrack = transTrackBuilder.build(track); - AnalyticalImpactPointExtrapolator extrapolator(transTrack.field()); - TrajectoryStateOnSurface closestState = - extrapolator.extrapolate(transTrack.impactPointState(), RecoVertex::convertPos(pv.position())); - pca = closestState.globalPosition(); - pca_cov = closestState.cartesianError().position(); - ip_vec = GlobalVector(pca.x() - pv.x(), pca.y() - pv.y(), pca.z() - pv.z()); - ip_cov = pca_cov + GlobalError(pv.covariance()); - VertexDistance3D pca_dist; - Measurement1D ip_mes = pca_dist.distance(pv, VertexState(pca, pca_cov)); - if (ip_vec.dot(GlobalVector(lepton.px(), lepton.py(), lepton.pz())) < 0) - ip_mes = Measurement1D(-1. * ip_mes.value(), ip_mes.error()); - - ipLength = ip_mes.value(); - ipLength_sig = ip_mes.significance(); - } - // Store PCA info - lepton.addUserFloat("refitPCA_x", pca.x()); - lepton.addUserFloat("refitPCA_y", pca.y()); - lepton.addUserFloat("refitPCA_z", pca.z()); - lepton.addUserFloat("refitPCACov_xx", pca_cov.cxx()); - lepton.addUserFloat("refitPCACov_yx", pca_cov.cyx()); - lepton.addUserFloat("refitPCACov_zx", pca_cov.czx()); - lepton.addUserFloat("refitPCACov_yy", pca_cov.cyy()); - lepton.addUserFloat("refitPCACov_zy", pca_cov.czy()); - lepton.addUserFloat("refitPCACov_zz", pca_cov.czz()); - // Store IP info - lepton.addUserFloat("IP", ipLength); - lepton.addUserFloat("IP_sig", ipLength_sig); - lepton.addUserFloat("IP_x", ip_vec.x()); - lepton.addUserFloat("IP_y", ip_vec.y()); - lepton.addUserFloat("IP_z", ip_vec.z()); - lepton.addUserFloat("IPCov_xx", ip_cov.cxx()); - lepton.addUserFloat("IPCov_yx", ip_cov.cyx()); - lepton.addUserFloat("IPCov_zx", ip_cov.czx()); - lepton.addUserFloat("IPCov_yy", ip_cov.cyy()); - lepton.addUserFloat("IPCov_zy", ip_cov.czy()); - lepton.addUserFloat("IPCov_zz", ip_cov.czz()); - // Store leading track parameters and covariance at reference point - if (track != nullptr) { - lepton.addUserFloat("track_qoverp", track->parameter(reco::TrackBase::i_qoverp)); - lepton.addUserFloat("track_lambda", track->parameter(reco::TrackBase::i_lambda)); - lepton.addUserFloat("track_phi", track->parameter(reco::TrackBase::i_phi)); - lepton.addUserFloat("track_dxy", track->parameter(reco::TrackBase::i_dxy)); - lepton.addUserFloat("track_dsz", track->parameter(reco::TrackBase::i_dsz)); - for (size_t i = 0; i < reco::TrackBase::dimension; ++i) - for (size_t j = i; j < reco::TrackBase::dimension; ++j) - lepton.addUserFloat("trackCov_" + std::to_string(j + i), track->covariance(j, i)); - } else { - lepton.addUserFloat("track_qoverp", 0); - lepton.addUserFloat("track_lambda", 0); - lepton.addUserFloat("track_phi", 0); - lepton.addUserFloat("track_dxy", 0); - lepton.addUserFloat("track_dsz", 0); - for (size_t i = 0; i < reco::TrackBase::dimension; ++i) - for (size_t j = i; j < reco::TrackBase::dimension; ++j) - lepton.addUserFloat("trackCov_" + std::to_string(j + i), 0); - } - // ... and z coordinate of magnetic field at track reference point - lepton.addUserFloat("bField_z", bField_z); -} - -template -void PATLeptonTimeLifeInfoUpdater::produceAndFillSVInfo(T& lepton, - const TransientTrackBuilder& transTrackBuilder, - const reco::Vertex& pv) {} - -template <> -void PATLeptonTimeLifeInfoUpdater::produceAndFillSVInfo(pat::Tau& tau, - const TransientTrackBuilder& transTrackBuilder, - const reco::Vertex& pv) { - //MB: Set initially SV to the tau vertex (tau.vertex()) or the PV choosen? - reco::Vertex sv; - GlobalVector flight_vec; - GlobalError flight_cov; - float flightLength = -1; - float flightLength_sig = -1; - - // Fit SV with tracks of charged tau decay products - int fitOK = 0; - if (tau.signalChargedHadrCands().size() + tau.signalLostTracks().size() > 1) { - // Get tracks form tau signal charged candidates - std::vector transTrks; - TransientVertex transVtx; - for (const auto& cand : tau.signalChargedHadrCands()) { - if (cand.isNull()) - continue; - const reco::Track* track = cand->bestTrack(); - if (track != nullptr) - transTrks.push_back(transTrackBuilder.build(track)); - } - for (const auto& cand : tau.signalLostTracks()) { - if (cand.isNull()) - continue; - const reco::Track* track = cand->bestTrack(); - if (track != nullptr) - transTrks.push_back(transTrackBuilder.build(track)); - } - // Fit SV with KalmanVertexFitter - fitOK = 1; - KalmanVertexFitter kvf(true); - if (transTrks.size() > 1) { - try { - transVtx = kvf.vertex(transTrks); - } catch (...) { - fitOK = -1; - } - } else { - fitOK = -1; - } - if (!transVtx.hasRefittedTracks() || transVtx.refittedTracks().size() != transTrks.size()) - fitOK = -1; - if (fitOK > 0) { - sv = transVtx; - // Get flight-length - // Full PV->SV flight vector with its covariance - flight_vec = GlobalVector(sv.x() - pv.x(), sv.y() - pv.y(), sv.z() - pv.z()); - flight_cov = transVtx.positionError() + GlobalError(pv.covariance()); - //MB: can be taken from tau itself (but with different fit of PV and SV) as follows: - //flightLength = sqrt(tau.flightLength().mag2()); - //flightLength_sig = tau.flightLengthSig(); - VertexDistance3D sv_dist; - Measurement1D flightLength_mes = sv_dist.signedDistance(pv, sv, GlobalVector(tau.px(), tau.py(), tau.pz())); - flightLength = flightLength_mes.value(); - flightLength_sig = flightLength_mes.significance(); - } - } - // Embed SV info in tau - // Store SV info - tau.addUserInt("refitSV_status", fitOK); - tau.addUserFloat("refitSV_x", sv.x()); - tau.addUserFloat("refitSV_y", sv.y()); - tau.addUserFloat("refitSV_z", sv.z()); - tau.addUserFloat("refitSV_chi2", sv.chi2()); - tau.addUserFloat("refitSV_ndof", sv.ndof()); - tau.addUserFloat("refitSVCov_xx", sv.covariance(0, 0)); - tau.addUserFloat("refitSVCov_yx", sv.covariance(1, 0)); - tau.addUserFloat("refitSVCov_zx", sv.covariance(2, 0)); - tau.addUserFloat("refitSVCov_yy", sv.covariance(1, 1)); - tau.addUserFloat("refitSVCov_zy", sv.covariance(2, 1)); - tau.addUserFloat("refitSVCov_zz", sv.covariance(2, 2)); - // Store flight-length info - tau.addUserFloat("flightLength", flightLength); - tau.addUserFloat("flightLength_sig", flightLength_sig); - tau.addUserFloat("flight_x", flight_vec.x()); - tau.addUserFloat("flight_y", flight_vec.y()); - tau.addUserFloat("flight_z", flight_vec.z()); - tau.addUserFloat("flightCov_xx", flight_cov.cxx()); - tau.addUserFloat("flightCov_yx", flight_cov.cyx()); - tau.addUserFloat("flightCov_zx", flight_cov.czx()); - tau.addUserFloat("flightCov_yy", flight_cov.cyy()); - tau.addUserFloat("flightCov_zy", flight_cov.czy()); - tau.addUserFloat("flightCov_zz", flight_cov.czz()); -} - -template -void PATLeptonTimeLifeInfoUpdater::fillDescriptions(edm::ConfigurationDescriptions& descriptions) { - // pat{Electron,Muon,Tau}TimeLifeInfoUpdater - edm::ParameterSetDescription desc; - - std::string lepCollName; - if (typeid(T) == typeid(pat::Electron)) - lepCollName = "slimmedElectrons"; - else if (typeid(T) == typeid(pat::Muon)) - lepCollName = "slimmedMuons"; - else if (typeid(T) == typeid(pat::Tau)) - lepCollName = "slimmedTaus"; - desc.add("src", edm::InputTag(lepCollName)); - desc.add("pvSource", edm::InputTag("offlineSlimmedPrimaryVertices")); - desc.add("pvChoice", useFront) - ->setComment( - "Define PV to compute IP: 0: first PV, 1: PV with the smallest dz of the tau leading track (default: 0)"); - - descriptions.addWithDefaultLabel(desc); -} - -#include "FWCore/Framework/interface/MakerMacros.h" -typedef PATLeptonTimeLifeInfoUpdater PATElectronTimeLifeInfoUpdater; -DEFINE_FWK_MODULE(PATElectronTimeLifeInfoUpdater); -typedef PATLeptonTimeLifeInfoUpdater PATMuonTimeLifeInfoUpdater; -DEFINE_FWK_MODULE(PATMuonTimeLifeInfoUpdater); -typedef PATLeptonTimeLifeInfoUpdater PATTauTimeLifeInfoUpdater; -DEFINE_FWK_MODULE(PATTauTimeLifeInfoUpdater);