Back to home page

Project CMSSW displayed by LXR
 
 

    

File indexing completed on 2025-05-08 02:18:41 

0001 /////////////////////////////// BToTrkLLBuilder ///////////////////////////////
0002 /// original authors: G Karathanasis (CERN),  G Melachroinos (NKUA)
0003 // takes selected track collection and a dilepton collection and produces B
0004 // moth// - ers using a three-track vertex
0005 
0006 #include <algorithm>
0007 #include <limits>
0008 #include <map>
0009 #include <memory>
0010 #include <string>
0011 #include <vector>
0012 
0013 #include "CommonTools/Statistics/interface/ChiSquaredProbability.h"
0014 #include "CommonTools/Utils/interface/StringCutObjectSelector.h"
0015 #include "DataFormats/BeamSpot/interface/BeamSpot.h"
0016 #include "DataFormats/Math/interface/deltaR.h"
0017 #include "DataFormats/PatCandidates/interface/CompositeCandidate.h"
0018 #include "DataFormats/PatCandidates/interface/PackedCandidate.h"
0019 #include "FWCore/Framework/interface/Event.h"
0020 #include "FWCore/Framework/interface/global/EDProducer.h"
0021 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0022 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0023 #include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
0024 #include "FWCore/Utilities/interface/InputTag.h"
0025 #include "KinVtxFitter.h"
0026 #include "MagneticField/Engine/interface/MagneticField.h"
0027 #include "MagneticField/Records/interface/IdealMagneticFieldRecord.h"
0028 #include "TrackingTools/Records/interface/TransientTrackRecord.h"
0029 #include "TrackingTools/TransientTrack/interface/TransientTrack.h"
0030 #include "TrackingTools/TransientTrack/interface/TransientTrackBuilder.h"
0031 #include "helper.h"
0032 
0033 class BToTrkLLBuilder : public edm::global::EDProducer<> {
0034   // perhaps we need better structure here (begin run etc)
0035 public:
0036   typedef std::vector<reco::TransientTrack> TransientTrackCollection;
0037 
0038   explicit BToTrkLLBuilder(const edm::ParameterSet &cfg)
0039       : bFieldToken_{esConsumes<MagneticField, IdealMagneticFieldRecord>()},
0040         pre_vtx_selection_{cfg.getParameter<std::string>("preVtxSelection")},
0041         post_vtx_selection_{cfg.getParameter<std::string>("postVtxSelection")},
0042         dileptons_{consumes<pat::CompositeCandidateCollection>(cfg.getParameter<edm::InputTag>("dileptons"))},
0043         leptons_ttracks_{consumes<TransientTrackCollection>(cfg.getParameter<edm::InputTag>("leptonTransientTracks"))},
0044         kaons_{consumes<pat::CompositeCandidateCollection>(cfg.getParameter<edm::InputTag>("kaons"))},
0045         kaons_ttracks_{consumes<TransientTrackCollection>(cfg.getParameter<edm::InputTag>("kaonsTransientTracks"))},
0046         track_mass_{cfg.getParameter<double>("trackMass")},
0047         pu_tracks_(consumes<pat::CompositeCandidateCollection>(cfg.getParameter<edm::InputTag>("PUtracks"))),
0048         beamspot_{consumes<reco::BeamSpot>(cfg.getParameter<edm::InputTag>("beamSpot"))},
0049         dilepton_constraint_{cfg.getParameter<bool>("dileptonMassContraint")} {
0050     produces<pat::CompositeCandidateCollection>();
0051   }
0052 
0053   ~BToTrkLLBuilder() override {}
0054 
0055   void produce(edm::StreamID, edm::Event &, const edm::EventSetup &) const override;
0056 
0057 private:
0058   const edm::ESGetToken<MagneticField, IdealMagneticFieldRecord> bFieldToken_;
0059 
0060   // selections
0061   const StringCutObjectSelector<pat::CompositeCandidate> pre_vtx_selection_;
0062   const StringCutObjectSelector<pat::CompositeCandidate> post_vtx_selection_;
0063 
0064   // inputs
0065   const edm::EDGetTokenT<pat::CompositeCandidateCollection> dileptons_;
0066   const edm::EDGetTokenT<TransientTrackCollection> leptons_ttracks_;
0067   const edm::EDGetTokenT<pat::CompositeCandidateCollection> kaons_;
0068   const edm::EDGetTokenT<TransientTrackCollection> kaons_ttracks_;
0069   const double track_mass_;
0070   const edm::EDGetTokenT<pat::CompositeCandidateCollection> pu_tracks_;
0071   const edm::EDGetTokenT<reco::BeamSpot> beamspot_;
0072   const bool dilepton_constraint_;
0073 };
0074 
0075 void BToTrkLLBuilder::produce(edm::StreamID, edm::Event &evt, edm::EventSetup const &iSetup) const {
0076   // input
0077   edm::Handle<pat::CompositeCandidateCollection> dileptons;
0078   evt.getByToken(dileptons_, dileptons);
0079   edm::Handle<TransientTrackCollection> leptons_ttracks;
0080   evt.getByToken(leptons_ttracks_, leptons_ttracks);
0081 
0082   edm::Handle<pat::CompositeCandidateCollection> kaons;
0083   evt.getByToken(kaons_, kaons);
0084   edm::Handle<TransientTrackCollection> kaons_ttracks;
0085   evt.getByToken(kaons_ttracks_, kaons_ttracks);
0086 
0087   edm::Handle<pat::CompositeCandidateCollection> pu_tracks;
0088   evt.getByToken(pu_tracks_, pu_tracks);
0089 
0090   edm::Handle<reco::BeamSpot> beamspot;
0091   evt.getByToken(beamspot_, beamspot);
0092 
0093   edm::ESHandle<MagneticField> fieldHandle;
0094   const auto &bField = iSetup.getData(bFieldToken_);
0095   AnalyticalImpactPointExtrapolator extrapolator(&bField);
0096 
0097   // output
0098   std::unique_ptr<pat::CompositeCandidateCollection> ret_val(new pat::CompositeCandidateCollection());
0099 
0100   for (size_t k_idx = 0; k_idx < kaons->size(); ++k_idx) {
0101     edm::Ptr<pat::CompositeCandidate> k_ptr(kaons, k_idx);
0102 
0103     math::PtEtaPhiMLorentzVector k_p4(k_ptr->pt(), k_ptr->eta(), k_ptr->phi(), bph::K_MASS);
0104 
0105     for (size_t ll_idx = 0; ll_idx < dileptons->size(); ++ll_idx) {
0106       edm::Ptr<pat::CompositeCandidate> ll_prt(dileptons, ll_idx);
0107       edm::Ptr<reco::Candidate> l1_ptr = ll_prt->userCand("l1");
0108       edm::Ptr<reco::Candidate> l2_ptr = ll_prt->userCand("l2");
0109       int l1_idx = ll_prt->userInt("l1_idx");
0110       int l2_idx = ll_prt->userInt("l2_idx");
0111 
0112       pat::CompositeCandidate cand;
0113       cand.setP4(ll_prt->p4() + k_p4);
0114       cand.setCharge(ll_prt->charge() + k_ptr->charge());
0115       // Use UserCands as they should not use memory but keep the Ptr itself
0116       // Put the lepton passing the corresponding selection
0117       cand.addUserCand("l1", l1_ptr);
0118       cand.addUserCand("l2", l2_ptr);
0119       cand.addUserCand("trk", k_ptr);
0120       cand.addUserCand("dilepton", ll_prt);
0121 
0122       cand.addUserInt("l1_idx", l1_idx);
0123       cand.addUserInt("l2_idx", l2_idx);
0124       cand.addUserInt("trk_idx", k_idx);
0125 
0126       auto dr_info = bph::min_max_dr({l1_ptr, l2_ptr, k_ptr});
0127       cand.addUserFloat("min_dr", dr_info.first);
0128       cand.addUserFloat("max_dr", dr_info.second);
0129 
0130       if (!pre_vtx_selection_(cand))
0131         continue;
0132 
0133       KinVtxFitter fitter({leptons_ttracks->at(l1_idx), leptons_ttracks->at(l2_idx), kaons_ttracks->at(k_idx)},
0134                           {l1_ptr->mass(), l2_ptr->mass(), bph::K_MASS},
0135                           {bph::LEP_SIGMA, bph::LEP_SIGMA, bph::K_SIGMA}  // some small sigma for the lepton mass
0136       );
0137 
0138       if (!fitter.success())
0139         continue;  // hardcoded, but do we need otherwise?
0140       cand.setVertex(reco::Candidate::Point(fitter.fitted_vtx().x(), fitter.fitted_vtx().y(), fitter.fitted_vtx().z()));
0141       cand.addUserInt("sv_OK", fitter.success());
0142       cand.addUserFloat("sv_chi2", fitter.chi2());
0143       cand.addUserFloat("sv_ndof", fitter.dof());  // float??
0144       cand.addUserFloat("sv_prob", fitter.prob());
0145       cand.addUserFloat("fitted_mll", (fitter.daughter_p4(0) + fitter.daughter_p4(1)).mass());
0146       auto fit_p4 = fitter.fitted_p4();
0147       cand.addUserFloat("fitted_pt", fit_p4.pt());
0148       cand.addUserFloat("fitted_eta", fit_p4.eta());
0149       cand.addUserFloat("fitted_phi", fit_p4.phi());
0150       cand.addUserFloat("fitted_mass", fitter.fitted_candidate().mass());
0151       cand.addUserFloat("fitted_massErr", sqrt(fitter.fitted_candidate().kinematicParametersError().matrix()(6, 6)));
0152 
0153       cand.addUserFloat("cos_theta_2D", bph::cos_theta_2D(fitter, *beamspot, cand.p4()));
0154       cand.addUserFloat("fitted_cos_theta_2D", bph::cos_theta_2D(fitter, *beamspot, fit_p4));
0155 
0156       auto lxy = bph::l_xy(fitter, *beamspot);
0157       cand.addUserFloat("l_xy", lxy.value());
0158       cand.addUserFloat("l_xy_unc", lxy.error());
0159       // track impact parameter from SV
0160       TrajectoryStateOnSurface tsos =
0161           extrapolator.extrapolate(kaons_ttracks->at(k_idx).impactPointState(), fitter.fitted_vtx());
0162       std::pair<bool, Measurement1D> cur2DIP =
0163           bph::signedTransverseImpactParameter(tsos, fitter.fitted_refvtx(), *beamspot);
0164       cand.addUserFloat("k_svip2d", cur2DIP.second.value());
0165       cand.addUserFloat("k_svip2d_err", cur2DIP.second.error());
0166 
0167       if (!post_vtx_selection_(cand))
0168         continue;
0169 
0170       cand.addUserFloat("vtx_x", cand.vx());
0171       cand.addUserFloat("vtx_y", cand.vy());
0172       cand.addUserFloat("vtx_z", cand.vz());
0173 
0174       const auto &covMatrix = fitter.fitted_vtx_uncertainty();
0175       cand.addUserFloat("vtx_cxx", covMatrix.cxx());
0176       cand.addUserFloat("vtx_cyy", covMatrix.cyy());
0177       cand.addUserFloat("vtx_czz", covMatrix.czz());
0178       cand.addUserFloat("vtx_cyx", covMatrix.cyx());
0179       cand.addUserFloat("vtx_czx", covMatrix.czx());
0180       cand.addUserFloat("vtx_czy", covMatrix.czy());
0181 
0182       // refitted daughters (leptons/tracks)
0183       std::vector<std::string> dnames{"l1", "l2", "trk"};
0184 
0185       for (size_t idaughter = 0; idaughter < dnames.size(); idaughter++) {
0186         cand.addUserFloat("fitted_" + dnames[idaughter] + "_pt", fitter.daughter_p4(idaughter).pt());
0187 
0188         cand.addUserFloat("fitted_" + dnames[idaughter] + "_eta", fitter.daughter_p4(idaughter).eta());
0189 
0190         cand.addUserFloat("fitted_" + dnames[idaughter] + "_phi", fitter.daughter_p4(idaughter).phi());
0191       }
0192 
0193       // compute isolation
0194       std::vector<float> isos = bph::TrackerIsolation(pu_tracks, cand, dnames);
0195       for (size_t idaughter = 0; idaughter < dnames.size(); idaughter++) {
0196         cand.addUserFloat(dnames[idaughter] + "_iso04", isos[idaughter]);
0197       }
0198 
0199       float constraint_sv_prob = -9;
0200       float constraint_pt = -9;
0201       float constraint_eta = -9;
0202       float constraint_phi = -9;
0203       float constraint_mass = -9;
0204       float constraint_massErr = -9;
0205       float constraint_mll = -9;
0206 
0207       const double dilepton_mass = ll_prt->userFloat("fitted_mass");
0208       const double jpsi_bin[2] = {2.8, 3.35};
0209       const double psi2s_bin[2] = {3.45, 3.85};
0210 
0211       if (dilepton_constraint_ && ((dilepton_mass > jpsi_bin[0] && dilepton_mass < jpsi_bin[1]) ||
0212                                    (dilepton_mass > psi2s_bin[0] && dilepton_mass < psi2s_bin[1]))) {
0213         ParticleMass JPsi_mass = 3.0969;   // Jpsi mass 3.096900±0.000006
0214         ParticleMass Psi2S_mass = 3.6861;  // Psi2S mass 3.6861093±0.0000034
0215         ParticleMass mass_constraint = (dilepton_mass < jpsi_bin[1]) ? JPsi_mass : Psi2S_mass;
0216 
0217         // Mass constraint is applied to the first two particles in the
0218         // "particles" vector Make sure that the first two particles are the
0219         // ones you want to constrain
0220         KinVtxFitter constraint_fitter(
0221             {leptons_ttracks->at(l1_idx), leptons_ttracks->at(l2_idx), kaons_ttracks->at(k_idx)},
0222             {l1_ptr->mass(), l2_ptr->mass(), bph::K_MASS},
0223             {bph::LEP_SIGMA, bph::LEP_SIGMA, bph::K_SIGMA},
0224             mass_constraint);
0225         if (constraint_fitter.success()) {
0226           auto constraint_p4 = constraint_fitter.fitted_p4();
0227           constraint_sv_prob = constraint_fitter.prob();
0228           constraint_pt = constraint_p4.pt();
0229           constraint_eta = constraint_p4.eta();
0230           constraint_phi = constraint_p4.phi();
0231           constraint_mass = constraint_fitter.fitted_candidate().mass();
0232           constraint_massErr = sqrt(constraint_fitter.fitted_candidate().kinematicParametersError().matrix()(6, 6));
0233           constraint_mll = (constraint_fitter.daughter_p4(0) + constraint_fitter.daughter_p4(1)).mass();
0234         }
0235       }
0236       cand.addUserFloat("constraint_sv_prob", constraint_sv_prob);
0237       cand.addUserFloat("constraint_pt", constraint_pt);
0238       cand.addUserFloat("constraint_eta", constraint_eta);
0239       cand.addUserFloat("constraint_phi", constraint_phi);
0240       cand.addUserFloat("constraint_mass", constraint_mass);
0241       cand.addUserFloat("constraint_massErr", constraint_massErr);
0242       cand.addUserFloat("constraint_mll", constraint_mll);
0243 
0244       ret_val->push_back(cand);
0245     }  // for(size_t ll_idx = 0; ll_idx < dileptons->size(); ++ll_idx) {
0246   }  // for(size_t k_idx = 0; k_idx < kaons->size(); ++k_idx)
0247 
0248   evt.put(std::move(ret_val));
0249 }
0250 
0251 #include "FWCore/Framework/interface/MakerMacros.h"
0252 DEFINE_FWK_MODULE(BToTrkLLBuilder);
[ Source navigation ] [ Diff markup ] [ Identifier search ] [ general search ]

This page was automatically generated by the 2.2.1 LXR engine.

The LXR team
Valid CSS 2.1! Valid HTML 4.01!