Back to home page

Project CMSSW displayed by LXR

 
 

    


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

0001 #include "FWCore/Framework/interface/global/EDProducer.h"
0002 #include "FWCore/Framework/interface/Event.h"
0003 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0004 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0005 #include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
0006 #include "FWCore/Utilities/interface/InputTag.h"
0007 #include "DataFormats/BeamSpot/interface/BeamSpot.h"
0008 #include "TrackingTools/TransientTrack/interface/TransientTrack.h"
0009 #include "TrackingTools/TransientTrack/interface/TransientTrackBuilder.h"
0010 #include "TrackingTools/Records/interface/TransientTrackRecord.h"
0011 
0012 #include "MagneticField/Engine/interface/MagneticField.h"
0013 #include "MagneticField/Records/interface/IdealMagneticFieldRecord.h"
0014 
0015 #include <vector>
0016 #include <memory>
0017 #include <map>
0018 #include <string>
0019 #include "DataFormats/PatCandidates/interface/PackedCandidate.h"
0020 #include "CommonTools/Utils/interface/StringCutObjectSelector.h"
0021 #include "DataFormats/PatCandidates/interface/CompositeCandidate.h"
0022 #include "DataFormats/Candidate/interface/Candidate.h"
0023 #include "DataFormats/Math/interface/deltaR.h"
0024 #include "CommonTools/Statistics/interface/ChiSquaredProbability.h"
0025 #include "helper.h"
0026 #include <limits>
0027 #include <algorithm>
0028 #include "KinVtxFitter.h"
0029 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0030 
0031 using namespace std;
0032 
0033 class BToV0TrkLLBuilder : public edm::global::EDProducer<> {
0034   // perhaps we need better structure here (begin run etc)
0035 public:
0036   typedef std::vector<reco::TransientTrack> TransientTrackCollection;
0037   explicit BToV0TrkLLBuilder(const edm::ParameterSet &cfg)
0038       :  // selections
0039         pre_vtx_selection_{cfg.getParameter<std::string>("preVtxSelection")},
0040         post_vtx_selection_{cfg.getParameter<std::string>("postVtxSelection")},
0041         //inputs
0042         dileptons_{consumes<pat::CompositeCandidateCollection>(cfg.getParameter<edm::InputTag>("dileptons"))},
0043         leptons_ttracks_{consumes<TransientTrackCollection>(cfg.getParameter<edm::InputTag>("leptonTransientTracks"))},
0044         V0s_ttracks_{consumes<TransientTrackCollection>(cfg.getParameter<edm::InputTag>("V0s_ttracks"))},
0045         V0s_{consumes<pat::CompositeCandidateCollection>(cfg.getParameter<edm::InputTag>("V0s"))},
0046         pions_{consumes<pat::CompositeCandidateCollection>(cfg.getParameter<edm::InputTag>("pions"))},
0047         pions_ttracks_{consumes<TransientTrackCollection>(cfg.getParameter<edm::InputTag>("pionsTransientTracks"))},
0048         beamspot_{consumes<reco::BeamSpot>(cfg.getParameter<edm::InputTag>("beamSpot"))},
0049         vertex_src_{consumes<reco::VertexCollection>(cfg.getParameter<edm::InputTag>("offlinePrimaryVertexSrc"))} {
0050     //output
0051     produces<pat::CompositeCandidateCollection>();
0052   }
0053 
0054   ~BToV0TrkLLBuilder() override {}
0055 
0056   void produce(edm::StreamID, edm::Event &, const edm::EventSetup &) const override;
0057 
0058   static void fillDescriptions(edm::ConfigurationDescriptions &descriptions) {}
0059 
0060 private:
0061   // selections
0062   const StringCutObjectSelector<pat::CompositeCandidate> pre_vtx_selection_;   // cut on the di-lepton before the SV fit
0063   const StringCutObjectSelector<pat::CompositeCandidate> post_vtx_selection_;  // cut on the di-lepton after the SV fit
0064 
0065   const edm::EDGetTokenT<pat::CompositeCandidateCollection> dileptons_;
0066   const edm::EDGetTokenT<TransientTrackCollection> leptons_ttracks_;
0067   const edm::EDGetTokenT<TransientTrackCollection> V0s_ttracks_;
0068   const edm::EDGetTokenT<pat::CompositeCandidateCollection> V0s_;
0069   const edm::EDGetTokenT<pat::CompositeCandidateCollection> pions_;
0070   const edm::EDGetTokenT<TransientTrackCollection> pions_ttracks_;
0071 
0072   const edm::EDGetTokenT<reco::BeamSpot> beamspot_;
0073   const edm::EDGetTokenT<reco::VertexCollection> vertex_src_;
0074 };
0075 
0076 void BToV0TrkLLBuilder::produce(edm::StreamID, edm::Event &evt, edm::EventSetup const &iSetup) const {
0077   //input
0078   edm::Handle<pat::CompositeCandidateCollection> dileptons;
0079   evt.getByToken(dileptons_, dileptons);
0080 
0081   edm::Handle<TransientTrackCollection> leptons_ttracks;
0082   evt.getByToken(leptons_ttracks_, leptons_ttracks);
0083   edm::Handle<pat::CompositeCandidateCollection> pions;
0084   evt.getByToken(pions_, pions);
0085 
0086   edm::Handle<TransientTrackCollection> pions_ttracks;
0087   evt.getByToken(pions_ttracks_, pions_ttracks);
0088 
0089   edm::Handle<pat::CompositeCandidateCollection> V0s;
0090   evt.getByToken(V0s_, V0s);
0091 
0092   edm::Handle<TransientTrackCollection> V0s_ttracks;
0093   evt.getByToken(V0s_ttracks_, V0s_ttracks);
0094 
0095   edm::Handle<reco::BeamSpot> beamspot;
0096   evt.getByToken(beamspot_, beamspot);
0097 
0098   edm::Handle<reco::VertexCollection> pvtxs;
0099   evt.getByToken(vertex_src_, pvtxs);
0100 
0101   edm::ESHandle<MagneticField> fieldHandle;
0102   const MagneticField *fMagneticField = fieldHandle.product();
0103   AnalyticalImpactPointExtrapolator extrapolator(fMagneticField);
0104 
0105   std::vector<int> used_lep1_id, used_lep2_id, used_pi_id, used_V0_id;
0106 
0107   // output
0108   std::unique_ptr<pat::CompositeCandidateCollection> ret_val(new pat::CompositeCandidateCollection());
0109   for (size_t V0_idx = 0; V0_idx < V0s->size(); ++V0_idx) {
0110     edm::Ptr<pat::CompositeCandidate> V0_ptr(V0s, V0_idx);
0111     math::PtEtaPhiMLorentzVector V0_p4(V0_ptr->userFloat("fitted_pt"),
0112                                        V0_ptr->userFloat("fitted_eta"),
0113                                        V0_ptr->userFloat("fitted_phi"),
0114                                        V0_ptr->userFloat("fitted_mass"));
0115     edm::Ptr<reco::Candidate> pi1_ptr = V0_ptr->userCand("trk1");
0116     edm::Ptr<reco::Candidate> pi2_ptr = V0_ptr->userCand("trk2");
0117     unsigned int pi1_idx = V0_ptr->userInt("trk1_idx");
0118     unsigned int pi2_idx = V0_ptr->userInt("trk2_idx");
0119     //    float pi1_dr = V0_ptr->userFloat("trk1_dr");
0120     //    float pi2_dr = V0_ptr->userFloat("trk2_dr");
0121     for (size_t pi_idx = 0; pi_idx < pions->size(); ++pi_idx) {
0122       edm::Ptr<pat::CompositeCandidate> pi_ptr(pions, pi_idx);
0123       if (pi1_idx == pi_idx || pi2_idx == pi_idx)
0124         continue;
0125       edm::Ptr<reco::Candidate> pi1_ptr(pions, pi1_idx);
0126       edm::Ptr<reco::Candidate> pi2_ptr(pions, pi2_idx);
0127       math::PtEtaPhiMLorentzVector pi_p4(pi_ptr->pt(), pi_ptr->eta(), pi_ptr->phi(), bph::PI_MASS);
0128       pat::CompositeCandidate cand;
0129       cand.setP4(pi_ptr->p4() + V0_p4);
0130       //      cand.setCharge(V0_ptr->userInt("fit_trk1_charge") + V0_ptr->userInt("fit_trk2_charge") + pi_ptr->charge());
0131       //      cand.addUserInt("fitted_charge",V0_ptr->userInt("fit_trk1_charge") + V0_ptr->userInt("fit_trk2_charge") + pi_ptr->charge());
0132       cand.addUserInt("pi_idx", pi_idx);
0133       cand.addUserInt("V0_idx", V0_idx);
0134       cand.addUserCand("pi", pi_ptr);
0135       cand.addUserCand("V0", V0_ptr);
0136       float dr = deltaR(pi_ptr->eta(), pi_ptr->phi(), V0_ptr->userFloat("fitted_eta"), V0_ptr->userFloat("fitted_phi"));
0137       cand.addUserFloat("V0pi_dr", dr);
0138       for (size_t ll_idx = 0; ll_idx < dileptons->size(); ++ll_idx) {
0139         edm::Ptr<pat::CompositeCandidate> ll_ptr(dileptons, ll_idx);
0140         edm::Ptr<reco::Candidate> l1_ptr = ll_ptr->userCand("l1");
0141         edm::Ptr<reco::Candidate> l2_ptr = ll_ptr->userCand("l2");
0142         int l1_idx = ll_ptr->userInt("l1_idx");
0143         int l2_idx = ll_ptr->userInt("l2_idx");
0144         cand.addUserCand("l1", l1_ptr);
0145         cand.addUserCand("l2", l2_ptr);
0146         cand.addUserCand("dilepton", ll_ptr);
0147         cand.addUserInt("l1_idx", l1_idx);
0148         cand.addUserInt("l2_idx", l2_idx);
0149         cand.addUserInt("ll_idx", ll_idx);
0150         cand.addUserCand("l1_ptr", l1_ptr);
0151         cand.addUserCand("l2_ptr", l2_ptr);
0152         auto lep1_p4 = l1_ptr->polarP4();
0153         auto lep2_p4 = l2_ptr->polarP4();
0154         lep1_p4.SetM(l1_ptr->mass());
0155         lep2_p4.SetM(l2_ptr->mass());
0156         auto pi_p4 = pi_ptr->polarP4();
0157         auto V0_p4 = V0_ptr->polarP4();
0158         pi_p4.SetM(bph::PI_MASS);
0159         V0_p4.SetM(V0_ptr->mass());
0160         cand.setP4(ll_ptr->p4() + pi_p4 + V0_p4);
0161         cand.setCharge(ll_ptr->charge() + V0_ptr->userInt("fit_trk1_charge") + V0_ptr->userInt("fit_trk2_charge") +
0162                        pi_ptr->charge());
0163         cand.addUserFloat("ll_V0_deltaR", reco::deltaR(*ll_ptr, *V0_ptr));
0164         cand.addUserFloat("ll_pi_deltaR", reco::deltaR(*ll_ptr, *pi_ptr));
0165         auto V0_dr_info = bph::min_max_dr({l1_ptr, l2_ptr, V0_ptr});
0166         cand.addUserFloat("V0_min_dr", V0_dr_info.first);
0167         cand.addUserFloat("V0_max_dr", V0_dr_info.second);
0168         auto pi_dr_info = bph::min_max_dr({l1_ptr, l2_ptr, pi_ptr});
0169         cand.addUserFloat("pi_min_dr", pi_dr_info.first);
0170         cand.addUserFloat("pi_max_dr", pi_dr_info.second);
0171         cand.addUserFloat("mIntermediate_unfitted", (pi_p4 + V0_p4).M());
0172         if (!pre_vtx_selection_(cand))
0173           continue;
0174 
0175         KinVtxFitter fitter;
0176         try {
0177           fitter = KinVtxFitter({leptons_ttracks->at(l1_idx),
0178                                  leptons_ttracks->at(l2_idx),
0179                                  pions_ttracks->at(pi_idx),
0180                                  V0s_ttracks->at(V0_idx)},
0181                                 {l1_ptr->mass(), l2_ptr->mass(), bph::PI_MASS, V0_ptr->mass()},
0182                                 {bph::LEP_SIGMA, bph::LEP_SIGMA, bph::PI_SIGMA, V0_ptr->userFloat("massErr")});
0183 
0184         } catch (const VertexException &e) {
0185           edm::LogWarning("KinematicFit") << "V0TrkLL Builder: Skipping candidate due to fit failure: " << e.what();
0186           continue;
0187         }
0188         if (!fitter.success())
0189           continue;
0190 
0191         cand.setVertex(
0192             reco::Candidate::Point(fitter.fitted_vtx().x(), fitter.fitted_vtx().y(), fitter.fitted_vtx().z()));
0193 
0194         TrajectoryStateOnSurface V0tsos =
0195             extrapolator.extrapolate(V0s_ttracks->at(V0_idx).impactPointState(), fitter.fitted_vtx());
0196         cand.addUserFloat("V0_dz", V0tsos.globalPosition().z() - fitter.fitted_vtx().z());
0197         cand.addUserFloat("V0_x", V0tsos.globalPosition().x());
0198         cand.addUserFloat("V0_y", V0tsos.globalPosition().y());
0199         cand.addUserFloat("V0_z", V0tsos.globalPosition().z());
0200         TrajectoryStateOnSurface pitsos =
0201             extrapolator.extrapolate(pions_ttracks->at(pi_idx).impactPointState(), fitter.fitted_vtx());
0202         cand.addUserFloat("pi_dz", pions_ttracks->at(pi_idx).track().vz() - fitter.fitted_vtx().z());  //
0203         cand.addUserFloat("pi_x", pions_ttracks->at(pi_idx).track().vx());  //pitsos.globalPosition().x());
0204         cand.addUserFloat("pi_y", pions_ttracks->at(pi_idx).track().vy());  //pitsos.globalPosition().y());
0205         cand.addUserFloat("pi_z", pions_ttracks->at(pi_idx).track().vz());  //pitsos.globalPosition().z());
0206         cand.addUserFloat(
0207             "V0trk_dz",
0208             V0tsos.globalPosition().z() - pions_ttracks->at(pi_idx).track().vz());  //pitsos.globalPosition().z());
0209 
0210         // vertex vars
0211         cand.addUserFloat("sv_chi2", fitter.chi2());
0212         cand.addUserFloat("sv_ndof", fitter.dof());
0213         cand.addUserFloat("sv_prob", fitter.prob());
0214         // refitted kinematic vars
0215         cand.addUserFloat("fitted_KstarPlus_mass", (fitter.daughter_p4(2) + fitter.daughter_p4(3)).mass());
0216         cand.addUserFloat("fitted_KstarPlus_pt", (fitter.daughter_p4(2) + fitter.daughter_p4(3)).pt());
0217         cand.addUserFloat("fitted_KstarPlus_eta", (fitter.daughter_p4(2) + fitter.daughter_p4(3)).eta());
0218         cand.addUserFloat("fitted_KstarPlus_phi", (fitter.daughter_p4(2) + fitter.daughter_p4(3)).phi());
0219         auto fit_p4 = fitter.fitted_p4();
0220         cand.addUserFloat("fitted_mass", fit_p4.mass());
0221         cand.addUserFloat("fitted_massErr", sqrt(fitter.fitted_candidate().kinematicParametersError().matrix()(6, 6)));
0222         cand.addUserFloat("fitted_mll_mass", (fitter.daughter_p4(0) + fitter.daughter_p4(1)).mass());
0223         cand.addUserFloat("fitted_mll_pt", (fitter.daughter_p4(0) + fitter.daughter_p4(1)).pt());
0224         cand.addUserFloat("fitted_mll_eta", (fitter.daughter_p4(0) + fitter.daughter_p4(1)).eta());
0225         cand.addUserFloat("fitted_mll_phi", (fitter.daughter_p4(0) + fitter.daughter_p4(1)).phi());
0226         cand.addUserFloat("fitted_pt", fit_p4.pt());
0227         cand.addUserFloat("fitted_eta", fit_p4.eta());
0228         cand.addUserFloat("fitted_phi", fit_p4.phi());
0229 
0230         const reco::BeamSpot &beamSpot = *beamspot;
0231         TrajectoryStateClosestToPoint theDCAXBS = fitter.fitted_candidate_ttrk().trajectoryStateClosestToPoint(
0232             GlobalPoint(beamSpot.position().x(), beamSpot.position().y(), beamSpot.position().z()));
0233         double DCAB0BS = -99.;
0234         double DCAB0BSErr = -99.;
0235 
0236         if (theDCAXBS.isValid() == true) {
0237           DCAB0BS = theDCAXBS.perigeeParameters().transverseImpactParameter();
0238           DCAB0BSErr = theDCAXBS.perigeeError().transverseImpactParameterError();
0239         }
0240         cand.addUserFloat("dca", DCAB0BS);
0241         cand.addUserFloat("dcaErr", DCAB0BSErr);
0242 
0243         cand.addUserFloat("vtx_x", cand.vx());
0244         cand.addUserFloat("vtx_y", cand.vy());
0245         cand.addUserFloat("vtx_z", cand.vz());
0246         cand.addUserFloat("vtx_ex", sqrt(fitter.fitted_vtx_uncertainty().cxx()));
0247         cand.addUserFloat("vtx_ey", sqrt(fitter.fitted_vtx_uncertainty().cyy()));
0248         cand.addUserFloat("vtx_ez", sqrt(fitter.fitted_vtx_uncertainty().czz()));
0249         // refitted daughters (leptons/tracks)
0250         std::vector<std::string> dnames{"l1", "l2", "pi", "V0"};
0251         for (size_t idaughter = 0; idaughter < dnames.size(); idaughter++) {
0252           cand.addUserFloat("fitted_" + dnames[idaughter] + "_pt", fitter.daughter_p4(idaughter).pt());
0253           cand.addUserFloat("fitted_" + dnames[idaughter] + "_eta", fitter.daughter_p4(idaughter).eta());
0254           cand.addUserFloat("fitted_" + dnames[idaughter] + "_phi", fitter.daughter_p4(idaughter).phi());
0255         }
0256         // other vars
0257         cand.addUserFloat("cos_theta_2D", bph::cos_theta_2D(fitter, *beamspot, cand.p4()));
0258         cand.addUserFloat("fitted_cos_theta_2D", bph::cos_theta_2D(fitter, *beamspot, fit_p4));
0259         auto lxy = bph::l_xy(fitter, *beamspot);
0260         cand.addUserFloat("l_xy", lxy.value());
0261         cand.addUserFloat("l_xy_unc", lxy.error());
0262 
0263         std::pair<bool, Measurement1D> cur2DIP_V0 =
0264             bph::signedTransverseImpactParameter(V0tsos, fitter.fitted_refvtx(), *beamspot);
0265         std::pair<bool, Measurement1D> cur3DIP_V0 =
0266             bph::signedImpactParameter3D(V0tsos, fitter.fitted_refvtx(), *beamspot, (*pvtxs)[0].position().z());
0267         std::pair<bool, Measurement1D> cur2DIP_pi =
0268             bph::signedTransverseImpactParameter(pitsos, fitter.fitted_refvtx(), *beamspot);
0269         std::pair<bool, Measurement1D> cur3DIP_pi =
0270             bph::signedImpactParameter3D(pitsos, fitter.fitted_refvtx(), *beamspot, (*pvtxs)[0].position().z());
0271         cand.addUserFloat("pi_svip2d", cur2DIP_pi.second.value());
0272         cand.addUserFloat("pi_svip2d_err", cur2DIP_pi.second.error());
0273         cand.addUserFloat("pi_svip3d", cur3DIP_pi.second.value());
0274         cand.addUserFloat("pi_svip3d_err", cur3DIP_pi.second.error());
0275         cand.addUserFloat("V0_svip2d", cur2DIP_V0.second.value());
0276         cand.addUserFloat("V0_svip2d_err", cur2DIP_V0.second.error());
0277         cand.addUserFloat("V0_svip3d", cur3DIP_V0.second.value());
0278         cand.addUserFloat("V0_svip3d_err", cur3DIP_V0.second.error());
0279         if (!post_vtx_selection_(cand))
0280           continue;
0281 
0282         // /////////////////////////////////////////////////////////////////
0283         // ///     Mass constrained fit START                            ///
0284         // /////////////////////////////////////////////////////////////////
0285 
0286         // Define variables
0287         bool sv_OK_withMC = false;
0288         float sv_chi2_withMC = cand.userFloat("sv_chi2");
0289         float sv_ndof_withMC = cand.userFloat("sv_ndof");
0290         float sv_prob_withMC = cand.userFloat("sv_prob");
0291         float fitted_mll_withMC = cand.userFloat("fitted_mll_mass");
0292         float fitted_pt_withMC = cand.userFloat("fitted_pt");
0293         float fitted_eta_withMC = cand.userFloat("fitted_eta");
0294         float fitted_phi_withMC = cand.userFloat("fitted_phi");
0295         float fitted_mass_withMC = cand.userFloat("fitted_mass");
0296         float fitted_massErr_withMC = cand.userFloat("fitted_massErr");
0297         float fitted_cos_theta_2D_withMC = cand.userFloat("fitted_cos_theta_2D");
0298         float l_xy_withMC = cand.userFloat("l_xy");
0299         float l_xy_unc_withMC = cand.userFloat("l_xy_unc");
0300         float vtx_x_withMC = cand.userFloat("vtx_x");
0301         float vtx_y_withMC = cand.userFloat("vtx_y");
0302         float vtx_z_withMC = cand.userFloat("vtx_z");
0303         float vtx_ex_withMC = cand.userFloat("vtx_ex");
0304         float vtx_ey_withMC = cand.userFloat("vtx_ey");
0305         float vtx_ez_withMC = cand.userFloat("vtx_ez");
0306         float fitted_l1_pt_withMC = cand.userFloat("fitted_l1_pt");
0307         float fitted_l1_eta_withMC = cand.userFloat("fitted_l1_eta");
0308         float fitted_l1_phi_withMC = cand.userFloat("fitted_l1_phi");
0309         float fitted_l2_pt_withMC = cand.userFloat("fitted_l2_pt");
0310         float fitted_l2_eta_withMC = cand.userFloat("fitted_l2_eta");
0311         float fitted_l2_phi_withMC = cand.userFloat("fitted_l2_phi");
0312         float fitted_V0_pt_withMC = cand.userFloat("fitted_V0_pt");
0313         float fitted_V0_eta_withMC = cand.userFloat("fitted_V0_eta");
0314         float fitted_V0_phi_withMC = cand.userFloat("fitted_V0_phi");
0315         float fitted_pi_pt_withMC = cand.userFloat("fitted_pi_pt");
0316         float fitted_pi_eta_withMC = cand.userFloat("fitted_pi_eta");
0317         float fitted_pi_phi_withMC = cand.userFloat("fitted_pi_phi");
0318 
0319         // Check dilepton mass from Bparticles to be in the jpsi bin
0320         const double dilepton_mass = ll_ptr->userFloat("fitted_mass");
0321         // const double dilepton_mass = (fitter.daughter_p4(0) + fitter.daughter_p4(1)).mass();
0322         const double jpsi_bin[2] = {2.8,
0323                                     3.2};  // {2.9, 3.2}; Start bin from 2.8 to be able to measure systematics later
0324         const double psi2s_bin[2] = {3.55, 3.8};
0325         if ((dilepton_mass > jpsi_bin[0] && dilepton_mass < jpsi_bin[1]) ||
0326             (dilepton_mass > psi2s_bin[0] && dilepton_mass < psi2s_bin[1])) {
0327           // JPsi  mass constrait
0328           // do mass constrained vertex fit
0329 
0330           ParticleMass JPsi_mass = 3.0969;   // Jpsi mass 3.096900±0.000006
0331           ParticleMass Psi2S_mass = 3.6861;  // Psi2S mass 3.6861093±0.0000034
0332           ParticleMass mass_constraint = (dilepton_mass < jpsi_bin[1]) ? JPsi_mass : Psi2S_mass;
0333           // Mass constraint is applied to the first two particles in the "particles" vector
0334           // Make sure that the first two particles are the ones you want to constrain
0335           KinVtxFitter constrained_fitter;
0336           try {
0337             constrained_fitter =
0338                 KinVtxFitter({leptons_ttracks->at(l1_idx),
0339                               leptons_ttracks->at(l2_idx),
0340                               pions_ttracks->at(pi_idx),
0341                               V0s_ttracks->at(V0_idx)},
0342                              {l1_ptr->mass(), l2_ptr->mass(), bph::PI_MASS, V0_ptr->mass()},
0343                              {bph::LEP_SIGMA, bph::LEP_SIGMA, bph::PI_SIGMA, V0_ptr->userFloat("massErr")},
0344                              mass_constraint);
0345           } catch (const VertexException &e) {
0346             edm::LogWarning("KinematicFit")
0347                 << "V0TrkLL Builder constraint: Skipping candidate due to fit failure: " << e.what();
0348             continue;
0349           }
0350           if (!constrained_fitter.success()) {
0351             // Save default values and continue
0352             cand.addUserInt("sv_OK_withMC", sv_OK_withMC);
0353             cand.addUserFloat("sv_chi2_withMC", sv_chi2_withMC);
0354             cand.addUserFloat("sv_ndof_withMC", sv_ndof_withMC);
0355             cand.addUserFloat("sv_prob_withMC", sv_prob_withMC);
0356             cand.addUserFloat("fitted_mll_withMC", fitted_mll_withMC);
0357             cand.addUserFloat("fitted_pt_withMC", fitted_pt_withMC);
0358             cand.addUserFloat("fitted_eta_withMC", fitted_eta_withMC);
0359             cand.addUserFloat("fitted_phi_withMC", fitted_phi_withMC);
0360             cand.addUserFloat("fitted_mass_withMC", fitted_mass_withMC);
0361             cand.addUserFloat("fitted_massErr_withMC", fitted_massErr_withMC);
0362             cand.addUserFloat("fitted_cos_theta_2D_withMC", fitted_cos_theta_2D_withMC);
0363             cand.addUserFloat("l_xy_withMC", l_xy_withMC);
0364             cand.addUserFloat("l_xy_unc_withMC", l_xy_unc_withMC);
0365             cand.addUserFloat("vtx_x_withMC", vtx_x_withMC);
0366             cand.addUserFloat("vtx_y_withMC", vtx_y_withMC);
0367             cand.addUserFloat("vtx_z_withMC", vtx_z_withMC);
0368             cand.addUserFloat("vtx_ex_withMC", vtx_ex_withMC);
0369             cand.addUserFloat("vtx_ey_withMC", vtx_ey_withMC);
0370             cand.addUserFloat("vtx_ez_withMC", vtx_ez_withMC);
0371             cand.addUserFloat("fitted_l1_pt_withMC", fitted_l1_pt_withMC);
0372             cand.addUserFloat("fitted_l1_eta_withMC", fitted_l1_eta_withMC);
0373             cand.addUserFloat("fitted_l1_phi_withMC", fitted_l1_phi_withMC);
0374             cand.addUserFloat("fitted_l2_pt_withMC", fitted_l2_pt_withMC);
0375             cand.addUserFloat("fitted_l2_eta_withMC", fitted_l2_eta_withMC);
0376             cand.addUserFloat("fitted_l2_phi_withMC", fitted_l2_phi_withMC);
0377             cand.addUserFloat("fitted_V0_pt_withMC", fitted_V0_pt_withMC);
0378             cand.addUserFloat("fitted_V0_eta_withMC", fitted_V0_eta_withMC);
0379             cand.addUserFloat("fitted_V0_phi_withMC", fitted_V0_phi_withMC);
0380             cand.addUserFloat("fitted_pi_pt_withMC", fitted_pi_pt_withMC);
0381             cand.addUserFloat("fitted_pi_eta_withMC", fitted_pi_eta_withMC);
0382             cand.addUserFloat("fitted_pi_phi_withMC", fitted_pi_phi_withMC);
0383             ret_val->push_back(cand);
0384             continue;
0385           }
0386           auto fit_p4_withMC = constrained_fitter.fitted_p4();
0387           sv_OK_withMC = constrained_fitter.success();
0388           sv_chi2_withMC = constrained_fitter.chi2();
0389           sv_ndof_withMC = constrained_fitter.dof();
0390           sv_prob_withMC = constrained_fitter.prob();
0391           fitted_mll_withMC = (constrained_fitter.daughter_p4(0) + constrained_fitter.daughter_p4(1)).mass();
0392           fitted_pt_withMC = fit_p4_withMC.pt();
0393           fitted_eta_withMC = fit_p4_withMC.eta();
0394           fitted_phi_withMC = fit_p4_withMC.phi();
0395           fitted_mass_withMC = constrained_fitter.fitted_candidate().mass();
0396           fitted_massErr_withMC = sqrt(constrained_fitter.fitted_candidate().kinematicParametersError().matrix()(6, 6));
0397           fitted_cos_theta_2D_withMC = bph::cos_theta_2D(constrained_fitter, *beamspot, fit_p4_withMC);
0398           auto lxy_withMC = bph::l_xy(constrained_fitter, *beamspot);
0399           l_xy_withMC = lxy_withMC.value();
0400           l_xy_unc_withMC = lxy_withMC.error();
0401           vtx_x_withMC = cand.vx();
0402           vtx_y_withMC = cand.vy();
0403           vtx_z_withMC = cand.vz();
0404           vtx_ex_withMC = sqrt(constrained_fitter.fitted_vtx_uncertainty().cxx());
0405           vtx_ey_withMC = sqrt(constrained_fitter.fitted_vtx_uncertainty().cyy());
0406           vtx_ez_withMC = sqrt(constrained_fitter.fitted_vtx_uncertainty().czz());
0407           fitted_l1_pt_withMC = constrained_fitter.daughter_p4(0).pt();
0408           fitted_l1_eta_withMC = constrained_fitter.daughter_p4(0).eta();
0409           fitted_l1_phi_withMC = constrained_fitter.daughter_p4(0).phi();
0410           fitted_l2_pt_withMC = constrained_fitter.daughter_p4(1).pt();
0411           fitted_l2_eta_withMC = constrained_fitter.daughter_p4(1).eta();
0412           fitted_l2_phi_withMC = constrained_fitter.daughter_p4(1).phi();
0413           fitted_V0_pt_withMC = constrained_fitter.daughter_p4(3).pt();
0414           fitted_V0_eta_withMC = constrained_fitter.daughter_p4(3).eta();
0415           fitted_V0_phi_withMC = constrained_fitter.daughter_p4(3).phi();
0416           fitted_pi_pt_withMC = constrained_fitter.daughter_p4(2).pt();
0417           fitted_pi_eta_withMC = constrained_fitter.daughter_p4(2).eta();
0418           fitted_pi_phi_withMC = constrained_fitter.daughter_p4(2).phi();
0419         }
0420         cand.addUserInt("sv_OK_withMC", sv_OK_withMC);
0421         cand.addUserFloat("sv_chi2_withMC", sv_chi2_withMC);
0422         cand.addUserFloat("sv_ndof_withMC", sv_ndof_withMC);
0423         cand.addUserFloat("sv_prob_withMC", sv_prob_withMC);
0424         cand.addUserFloat("fitted_mll_withMC", fitted_mll_withMC);
0425         cand.addUserFloat("fitted_pt_withMC", fitted_pt_withMC);
0426         cand.addUserFloat("fitted_eta_withMC", fitted_eta_withMC);
0427         cand.addUserFloat("fitted_phi_withMC", fitted_phi_withMC);
0428         cand.addUserFloat("fitted_mass_withMC", fitted_mass_withMC);
0429         cand.addUserFloat("fitted_massErr_withMC", fitted_massErr_withMC);
0430         cand.addUserFloat("fitted_cos_theta_2D_withMC", fitted_cos_theta_2D_withMC);
0431         cand.addUserFloat("l_xy_withMC", l_xy_withMC);
0432         cand.addUserFloat("l_xy_unc_withMC", l_xy_unc_withMC);
0433         cand.addUserFloat("vtx_x_withMC", vtx_x_withMC);
0434         cand.addUserFloat("vtx_y_withMC", vtx_y_withMC);
0435         cand.addUserFloat("vtx_z_withMC", vtx_z_withMC);
0436         cand.addUserFloat("vtx_ex_withMC", vtx_ex_withMC);
0437         cand.addUserFloat("vtx_ey_withMC", vtx_ey_withMC);
0438         cand.addUserFloat("vtx_ez_withMC", vtx_ez_withMC);
0439         cand.addUserFloat("fitted_l1_pt_withMC", fitted_l1_pt_withMC);
0440         cand.addUserFloat("fitted_l1_eta_withMC", fitted_l1_eta_withMC);
0441         cand.addUserFloat("fitted_l1_phi_withMC", fitted_l1_phi_withMC);
0442         cand.addUserFloat("fitted_l2_pt_withMC", fitted_l2_pt_withMC);
0443         cand.addUserFloat("fitted_l2_eta_withMC", fitted_l2_eta_withMC);
0444         cand.addUserFloat("fitted_l2_phi_withMC", fitted_l2_phi_withMC);
0445         cand.addUserFloat("fitted_V0_pt_withMC", fitted_V0_pt_withMC);
0446         cand.addUserFloat("fitted_V0_eta_withMC", fitted_V0_eta_withMC);
0447         cand.addUserFloat("fitted_V0_phi_withMC", fitted_V0_phi_withMC);
0448         cand.addUserFloat("fitted_pi_pt_withMC", fitted_pi_pt_withMC);
0449         cand.addUserFloat("fitted_pi_eta_withMC", fitted_pi_eta_withMC);
0450         cand.addUserFloat("fitted_pi_phi_withMC", fitted_pi_phi_withMC);
0451 
0452         // /////////////////////////////////////////////////////////////////
0453         // ///     Mass constrained fit END                              ///
0454         // /////////////////////////////////////////////////////////////////
0455         ret_val->push_back(cand);
0456         //  std::cout << "BGIKE" << endl;
0457       }  //     for(size_t ll_idx = 0; ll_idx < dileptons->size(); ++ll_idx)
0458     }  //   for(size_t pi_idx = 0; pi_idx < pions->size(); ++V0_idx)
0459   }  // for(size_t V0_idx = 0; V0_idx < V0s->size(); ++V0_idx)
0460   for (auto &cand : *ret_val) {
0461     cand.addUserInt("n_pi_used", std::count(used_pi_id.begin(), used_pi_id.end(), cand.userInt("pi_idx")));
0462     cand.addUserInt("n_V0_used", std::count(used_V0_id.begin(), used_V0_id.end(), cand.userInt("V0_idx")));
0463     cand.addUserInt("n_l1_used",
0464                     std::count(used_lep1_id.begin(), used_lep1_id.end(), cand.userInt("l1_idx")) +
0465                         std::count(used_lep2_id.begin(), used_lep2_id.end(), cand.userInt("l1_idx")));
0466     cand.addUserInt("n_l2_used",
0467                     std::count(used_lep1_id.begin(), used_lep1_id.end(), cand.userInt("l2_idx")) +
0468                         std::count(used_lep2_id.begin(), used_lep2_id.end(), cand.userInt("l2_idx")));
0469   }
0470   evt.put(std::move(ret_val));
0471 }
0472 #include "FWCore/Framework/interface/MakerMacros.h"
0473 DEFINE_FWK_MODULE(BToV0TrkLLBuilder);