Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-03-26 00:00:19

0001 // -*- C++ -*-
0002 //
0003 // Package:    V0Validator
0004 // Class:      V0Validator
0005 //
0006 /**\class V0Validator V0Validator.cc Validation/RecoVertex/src/V0Validator.cc
0007 
0008  Description: Creates validation histograms for RecoVertex/V0Producer
0009 
0010  Implementation:
0011      <Notes on implementation>
0012 */
0013 //
0014 // Original Author:  Brian Drell
0015 //         Created:  Wed Feb 18 17:21:04 MST 2009
0016 //
0017 //
0018 
0019 #include "Validation/RecoVertex/interface/V0Validator.h"
0020 #include "DataFormats/VertexReco/interface/VertexFwd.h"
0021 
0022 typedef std::vector<TrackingVertex> TrackingVertexCollection;
0023 typedef edm::Ref<TrackingVertexCollection> TrackingVertexRef;
0024 typedef edm::RefVector<edm::HepMCProduct, HepMC::GenVertex> GenVertexRefVector;
0025 typedef edm::RefVector<edm::HepMCProduct, HepMC::GenParticle> GenParticleRefVector;
0026 
0027 V0Validator::V0Validator(const edm::ParameterSet& iConfig)
0028     : theDQMRootFileName(iConfig.getUntrackedParameter<std::string>("DQMRootFileName")),
0029       dirName(iConfig.getUntrackedParameter<std::string>("dirName")),
0030       recoRecoToSimCollectionToken_(
0031           consumes<reco::RecoToSimCollection>(iConfig.getUntrackedParameter<edm::InputTag>("trackAssociatorMap"))),
0032       recoSimToRecoCollectionToken_(
0033           consumes<reco::SimToRecoCollection>(iConfig.getUntrackedParameter<edm::InputTag>("trackAssociatorMap"))),
0034       trackingVertexCollection_Token_(
0035           consumes<TrackingVertexCollection>(iConfig.getUntrackedParameter<edm::InputTag>("trackingVertexCollection"))),
0036       vec_recoVertex_Token_(
0037           consumes<std::vector<reco::Vertex> >(iConfig.getUntrackedParameter<edm::InputTag>("vertexCollection"))),
0038       recoVertexCompositeCandidateCollection_k0s_Token_(consumes<reco::VertexCompositeCandidateCollection>(
0039           iConfig.getUntrackedParameter<edm::InputTag>("kShortCollection"))),
0040       recoVertexCompositeCandidateCollection_lambda_Token_(consumes<reco::VertexCompositeCandidateCollection>(
0041           iConfig.getUntrackedParameter<edm::InputTag>("lambdaCollection"))) {}
0042 
0043 V0Validator::~V0Validator() {}
0044 
0045 void V0Validator::bookHistograms(DQMStore::IBooker& ibooker, edm::Run const&, edm::EventSetup const&) {
0046   double minKsMass = 0.49767 - 0.07;
0047   double maxKsMass = 0.49767 + 0.07;
0048   double minLamMass = 1.1156 - 0.05;
0049   double maxLamMass = 1.1156 + 0.05;
0050   int ksMassNbins = 100;
0051   double ksMassXmin = minKsMass;
0052   double ksMassXmax = maxKsMass;
0053   int lamMassNbins = 100;
0054   double lamMassXmin = minLamMass;
0055   double lamMassXmax = maxLamMass;
0056 
0057   ibooker.cd();
0058   std::string subDirName = V0Validator::dirName + "/K0";
0059   ibooker.setCurrentFolder(subDirName);
0060 
0061   candidateEffVsR_num_[V0Validator::KSHORT] =
0062       ibooker.book1D("K0sEffVsR_num", "K^{0}_{S} Efficiency vs #rho", 80, 0., 40.);
0063   candidateEffVsEta_num_[V0Validator::KSHORT] =
0064       ibooker.book1D("K0sEffVsEta_num", "K^{0}_{S} Efficiency vs #eta", 40, -2.5, 2.5);
0065   candidateEffVsPt_num_[V0Validator::KSHORT] =
0066       ibooker.book1D("K0sEffVsPt_num", "K^{0}_{S} Efficiency vs p_{T}", 70, 0., 20.);
0067 
0068   candidateTkEffVsR_num_[V0Validator::KSHORT] =
0069       ibooker.book1D("K0sTkEffVsR_num", "K^{0}_{S} Tracking Efficiency vs #rho", 80, 0., 40.);
0070   candidateTkEffVsEta_num_[V0Validator::KSHORT] =
0071       ibooker.book1D("K0sTkEffVsEta_num", "K^{0}_{S} Tracking Efficiency vs #eta", 40, -2.5, 2.5);
0072   candidateTkEffVsPt_num_[V0Validator::KSHORT] =
0073       ibooker.book1D("K0sTkEffVsPt_num", "K^{0}_{S} Tracking Efficiency vs p_{T}", 70, 0., 20.);
0074 
0075   candidateEffVsR_denom_[V0Validator::KSHORT] =
0076       ibooker.book1D("K0sEffVsR_denom", "K^{0}_{S} Efficiency vs #rho", 80, 0., 40.);
0077   candidateEffVsEta_denom_[V0Validator::KSHORT] =
0078       ibooker.book1D("K0sEffVsEta_denom", "K^{0}_{S} Efficiency vs #eta", 40, -2.5, 2.5);
0079   candidateEffVsPt_denom_[V0Validator::KSHORT] =
0080       ibooker.book1D("K0sEffVsPt_denom", "K^{0}_{S} Efficiency vs p_{T}", 70, 0., 20.);
0081 
0082   candidateFakeVsR_num_[V0Validator::KSHORT] =
0083       ibooker.book1D("K0sFakeVsR_num", "K^{0}_{S} Fake Rate vs #rho", 80, 0., 40.);
0084   candidateFakeVsEta_num_[V0Validator::KSHORT] =
0085       ibooker.book1D("K0sFakeVsEta_num", "K^{0}_{S} Fake Rate vs #eta", 40, -2.5, 2.5);
0086   candidateFakeVsPt_num_[V0Validator::KSHORT] =
0087       ibooker.book1D("K0sFakeVsPt_num", "K^{0}_{S} Fake Rate vs p_{T}", 70, 0., 20.);
0088   candidateTkFakeVsR_num_[V0Validator::KSHORT] =
0089       ibooker.book1D("K0sTkFakeVsR_num", "K^{0}_{S} Tracking Fake Rate vs #rho", 80, 0., 80.);
0090   candidateTkFakeVsEta_num_[V0Validator::KSHORT] =
0091       ibooker.book1D("K0sTkFakeVsEta_num", "K^{0}_{S} Tracking Fake Rate vs #eta", 40, -2.5, 2.5);
0092   candidateTkFakeVsPt_num_[V0Validator::KSHORT] =
0093       ibooker.book1D("K0sTkFakeVsPt_num", "K^{0}_{S} Tracking Fake Rate vs p_{T}", 70, 0., 20.);
0094 
0095   candidateFakeVsR_denom_[V0Validator::KSHORT] =
0096       ibooker.book1D("K0sFakeVsR_denom", "K^{0}_{S} Fake Rate vs #rho", 80, 0., 40.);
0097   candidateFakeVsEta_denom_[V0Validator::KSHORT] =
0098       ibooker.book1D("K0sFakeVsEta_denom", "K^{0}_{S} Fake Rate vs #eta", 40, -2.5, 2.5);
0099   candidateFakeVsPt_denom_[V0Validator::KSHORT] =
0100       ibooker.book1D("K0sFakeVsPt_denom", "K^{0}_{S} Fake Rate vs p_{T}", 70, 0., 20.);
0101   nCandidates_[V0Validator::KSHORT] = ibooker.book1D("nK0s", "Number of K^{0}_{S} found per event", 60, 0., 60.);
0102   fakeCandidateMass_[V0Validator::KSHORT] =
0103       ibooker.book1D("ksMassFake", "Mass of fake K0S", ksMassNbins, minKsMass, maxKsMass);
0104   goodCandidateMass[V0Validator::KSHORT] =
0105       ibooker.book1D("ksMassGood", "Mass of good reco K0S", ksMassNbins, minKsMass, maxKsMass);
0106   candidateMassAll[V0Validator::KSHORT] =
0107       ibooker.book1D("ksMassAll", "Invariant mass of all K0S", ksMassNbins, ksMassXmin, ksMassXmax);
0108   candidateFakeDauRadDist_[V0Validator::KSHORT] =
0109       ibooker.book1D("radDistFakeKs", "Production radius of daughter particle of Ks fake", 100, 0., 15.);
0110   candidateStatus_[V0Validator::KSHORT] = ibooker.book1D("ksCandStatus", "Fake type by cand status", 10, 0., 10.);
0111 
0112   // Lambda Plots follow
0113 
0114   subDirName = V0Validator::dirName + "/Lambda";
0115   ibooker.setCurrentFolder(subDirName);
0116 
0117   candidateEffVsR_num_[V0Validator::LAMBDA] =
0118       ibooker.book1D("LamEffVsR_num", "#Lambda^{0} Efficiency vs #rho", 80, 0., 40.);
0119   candidateEffVsEta_num_[V0Validator::LAMBDA] =
0120       ibooker.book1D("LamEffVsEta_num", "#Lambda^{0} Efficiency vs #eta", 40, -2.5, 2.5);
0121   candidateEffVsPt_num_[V0Validator::LAMBDA] =
0122       ibooker.book1D("LamEffVsPt_num", "#Lambda^{0} Efficiency vs p_{T}", 70, 0., 20.);
0123 
0124   candidateTkEffVsR_num_[V0Validator::LAMBDA] =
0125       ibooker.book1D("LamTkEffVsR_num", "#Lambda^{0} TrackingEfficiency vs #rho", 80, 0., 40.);
0126   candidateTkEffVsEta_num_[V0Validator::LAMBDA] =
0127       ibooker.book1D("LamTkEffVsEta_num", "#Lambda^{0} Tracking Efficiency vs #eta", 40, -2.5, 2.5);
0128   candidateTkEffVsPt_num_[V0Validator::LAMBDA] =
0129       ibooker.book1D("LamTkEffVsPt_num", "#Lambda^{0} Tracking Efficiency vs p_{T}", 70, 0., 20.);
0130 
0131   candidateEffVsR_denom_[V0Validator::LAMBDA] =
0132       ibooker.book1D("LamEffVsR_denom", "#Lambda^{0} Efficiency vs #rho", 80, 0., 40.);
0133   candidateEffVsEta_denom_[V0Validator::LAMBDA] =
0134       ibooker.book1D("LamEffVsEta_denom", "#Lambda^{0} Efficiency vs #eta", 40, -2.5, 2.5);
0135   candidateEffVsPt_denom_[V0Validator::LAMBDA] =
0136       ibooker.book1D("LamEffVsPt_denom", "#Lambda^{0} Efficiency vs p_{T}", 70, 0., 20.);
0137 
0138   candidateFakeVsR_num_[V0Validator::LAMBDA] =
0139       ibooker.book1D("LamFakeVsR_num", "#Lambda^{0} Fake Rate vs #rho", 80, 0., 40.);
0140   candidateFakeVsEta_num_[V0Validator::LAMBDA] =
0141       ibooker.book1D("LamFakeVsEta_num", "#Lambda^{0} Fake Rate vs #eta", 40, -2.5, 2.5);
0142   candidateFakeVsPt_num_[V0Validator::LAMBDA] =
0143       ibooker.book1D("LamFakeVsPt_num", "#Lambda^{0} Fake Rate vs p_{T}", 70, 0., 20.);
0144   candidateTkFakeVsR_num_[V0Validator::LAMBDA] =
0145       ibooker.book1D("LamTkFakeVsR_num", "#Lambda^{0} Tracking Fake Rate vs #rho", 80, 0., 40.);
0146   candidateTkFakeVsEta_num_[V0Validator::LAMBDA] =
0147       ibooker.book1D("LamTkFakeVsEta_num", "#Lambda^{0} Tracking Fake Rate vs #eta", 40, -2.5, 2.5);
0148   candidateTkFakeVsPt_num_[V0Validator::LAMBDA] =
0149       ibooker.book1D("LamTkFakeVsPt_num", "#Lambda^{0} Tracking Fake Rate vs p_{T}", 70, 0., 20.);
0150 
0151   candidateFakeVsR_denom_[V0Validator::LAMBDA] =
0152       ibooker.book1D("LamFakeVsR_denom", "#Lambda^{0} Fake Rate vs #rho", 80, 0., 40.);
0153   candidateFakeVsEta_denom_[V0Validator::LAMBDA] =
0154       ibooker.book1D("LamFakeVsEta_denom", "#Lambda^{0} Fake Rate vs #eta", 40, -2.5, 2.5);
0155   candidateFakeVsPt_denom_[V0Validator::LAMBDA] =
0156       ibooker.book1D("LamFakeVsPt_denom", "#Lambda^{0} Fake Rate vs p_{T}", 70, 0., 20.);
0157 
0158   nCandidates_[V0Validator::LAMBDA] = ibooker.book1D("nLam", "Number of #Lambda^{0} found per event", 60, 0., 60.);
0159   fakeCandidateMass_[V0Validator::LAMBDA] =
0160       ibooker.book1D("lamMassFake", "Mass of fake Lambda", lamMassNbins, minLamMass, maxLamMass);
0161   goodCandidateMass[V0Validator::LAMBDA] =
0162       ibooker.book1D("lamMassGood", "Mass of good Lambda", lamMassNbins, minLamMass, maxLamMass);
0163 
0164   candidateMassAll[V0Validator::LAMBDA] =
0165       ibooker.book1D("lamMassAll", "Invariant mass of all #Lambda^{0}", lamMassNbins, lamMassXmin, lamMassXmax);
0166   candidateFakeDauRadDist_[V0Validator::LAMBDA] =
0167       ibooker.book1D("radDistFakeLam", "Production radius of daughter particle of Lam fake", 100, 0., 15.);
0168 
0169   candidateStatus_[V0Validator::LAMBDA] = ibooker.book1D("ksCandStatus", "Fake type by cand status", 10, 0., 10.);
0170 }
0171 
0172 void V0Validator::doFakeRates(const reco::VertexCompositeCandidateCollection& collection,
0173                               const reco::RecoToSimCollection& recotosimCollection,
0174                               V0Type v0_type,
0175                               int particle_pdgid,
0176                               int misreconstructed_particle_pdgid) {
0177   using namespace edm;
0178 
0179   int numCandidateFound = 0;
0180   int realCandidateFound = 0;
0181   double mass = 0.;
0182   float CandidatepT = 0.;
0183   float CandidateEta = 0.;
0184   float CandidateR = 0.;
0185   int CandidateStatus = 0;
0186   const unsigned int NUM_DAUGHTERS = 2;
0187   if (!collection.empty()) {
0188     for (reco::VertexCompositeCandidateCollection::const_iterator iCandidate = collection.begin();
0189          iCandidate != collection.end();
0190          iCandidate++) {
0191       // Fill values to be histogrammed
0192       mass = iCandidate->mass();
0193       CandidatepT = (sqrt(iCandidate->momentum().perp2()));
0194       CandidateEta = iCandidate->momentum().eta();
0195       CandidateR = (sqrt(iCandidate->vertex().perp2()));
0196       candidateMassAll[v0_type]->Fill(mass);
0197       CandidateStatus = 0;
0198 
0199       std::array<reco::TrackRef, NUM_DAUGHTERS> theDaughterTracks = {
0200           {(*(dynamic_cast<const reco::RecoChargedCandidate*>(iCandidate->daughter(0)))).track(),
0201            (*(dynamic_cast<const reco::RecoChargedCandidate*>(iCandidate->daughter(1)))).track()}};
0202 
0203       TrackingParticleRef tpref;
0204       TrackingParticleRef firstDauTP;
0205       TrackingVertexRef candidateVtx;
0206 
0207       std::array<double, NUM_DAUGHTERS> radDist;
0208       // Loop through candidate's daugher tracks
0209       for (View<reco::Track>::size_type i = 0; i < theDaughterTracks.size(); ++i) {
0210         radDist = {{-1., -1.}};
0211         // Found track from theDaughterTracks
0212         RefToBase<reco::Track> track(theDaughterTracks.at(i));
0213 
0214         if (recotosimCollection.find(track) != recotosimCollection.end()) {
0215           const std::vector<std::pair<TrackingParticleRef, double> >& tp = recotosimCollection[track];
0216           if (!tp.empty()) {
0217             tpref = tp.begin()->first;
0218 
0219             TrackingVertexRef parentVertex = tpref->parentVertex();
0220             if (parentVertex.isNonnull()) {
0221               radDist[i] = parentVertex->position().R();
0222               if (candidateVtx.isNonnull()) {
0223                 if (candidateVtx->position() == parentVertex->position()) {
0224                   if (parentVertex->nDaughterTracks() == 2) {
0225                     if (parentVertex->nSourceTracks() == 0) {
0226                       // No source tracks found for candidate's
0227                       // vertex: it shouldn't happen, but does for
0228                       // evtGen events
0229                       CandidateStatus = 6;
0230                     }
0231 
0232                     for (TrackingVertex::tp_iterator iTP = parentVertex->sourceTracks_begin();
0233                          iTP != parentVertex->sourceTracks_end();
0234                          iTP++) {
0235                       if (abs((*iTP)->pdgId()) == particle_pdgid) {
0236                         CandidateStatus = 1;
0237                         realCandidateFound++;
0238                         numCandidateFound += 1.;
0239                         goodCandidateMass[v0_type]->Fill(mass);
0240                       } else {
0241                         CandidateStatus = 2;
0242                         if (abs((*iTP)->pdgId()) == misreconstructed_particle_pdgid) {
0243                           CandidateStatus = 7;
0244                         }
0245                       }
0246                     }
0247                   } else {
0248                     // Found a bad match because the mother has too
0249                     // many daughters
0250                     CandidateStatus = 3;
0251                   }
0252                 } else {
0253                   // Found a bad match because the parent vertices
0254                   // from the two tracks are different
0255                   CandidateStatus = 4;
0256                 }
0257               } else {
0258                 // if candidateVtx is null, fill it with parentVertex
0259                 // to compare to the parentVertex from the second
0260                 // track
0261                 candidateVtx = parentVertex;
0262                 firstDauTP = tpref;
0263               }
0264             }  // parent vertex is null
0265           }    // check on associated tp size zero
0266         } else {
0267           CandidateStatus = 5;
0268         }
0269       }  // Loop on candidate's daughter tracks
0270 
0271       // fill the fake rate histograms
0272       if (CandidateStatus > 1) {
0273         candidateFakeVsR_num_[v0_type]->Fill(CandidateR);
0274         candidateFakeVsEta_num_[v0_type]->Fill(CandidateEta);
0275         candidateFakeVsPt_num_[v0_type]->Fill(CandidatepT);
0276         candidateStatus_[v0_type]->Fill((float)CandidateStatus);
0277         fakeCandidateMass_[v0_type]->Fill(mass);
0278         for (auto distance : radDist) {
0279           if (distance > 0)
0280             candidateFakeDauRadDist_[v0_type]->Fill(distance);
0281         }
0282       }
0283       if (CandidateStatus == 5) {
0284         candidateTkFakeVsR_num_[v0_type]->Fill(CandidateR);
0285         candidateTkFakeVsEta_num_[v0_type]->Fill(CandidateEta);
0286         candidateTkFakeVsPt_num_[v0_type]->Fill(CandidatepT);
0287       }
0288       candidateFakeVsR_denom_[v0_type]->Fill(CandidateR);
0289       candidateFakeVsEta_denom_[v0_type]->Fill(CandidateEta);
0290       candidateFakeVsPt_denom_[v0_type]->Fill(CandidatepT);
0291     }  // Loop on candidates
0292   }    // check on presence of candidate's collection in the event
0293   nCandidates_[v0_type]->Fill((float)numCandidateFound);
0294 }
0295 
0296 void V0Validator::doEfficiencies(const TrackingVertexCollection& gen_vertices,
0297                                  V0Type v0_type,
0298                                  int parent_particle_id,
0299                                  int first_daughter_id,  /* give only positive charge */
0300                                  int second_daughter_id, /* give only positive charge */
0301                                  const reco::VertexCompositeCandidateCollection& collection,
0302                                  const reco::SimToRecoCollection& simtorecoCollection) {
0303   /* We store the TrackRef of the tracks that have been used to
0304    * produce the V0 under consideration here. This is used later to
0305    * check if a specific V0 has been really reconstructed or not. The
0306    * ordering is based on the key_index of the reference, since it
0307    * indeed does not matter that much. */
0308 
0309   std::set<V0Couple> reconstructed_V0_couples;
0310   if (!collection.empty()) {
0311     for (reco::VertexCompositeCandidateCollection::const_iterator iCandidate = collection.begin();
0312          iCandidate != collection.end();
0313          iCandidate++) {
0314       reconstructed_V0_couples.insert(
0315           V0Couple((dynamic_cast<const reco::RecoChargedCandidate*>(iCandidate->daughter(0)))->track(),
0316                    (dynamic_cast<const reco::RecoChargedCandidate*>(iCandidate->daughter(1)))->track()));
0317     }
0318   }
0319 
0320   /* PSEUDO CODE
0321      for v in gen_vertices
0322        if v.eventId().BX() !=0 continue
0323        if v.nDaughterTracks != 2 continue
0324        for source in v.sourceTracks_begin
0325          if source is parent_particle_id
0326           for daughter in v.daughterTracks_begin
0327            if daughter in region_and_kine_cuts
0328              decay_found
0329    */
0330   unsigned int candidateEff[2] = {0, 0};
0331   for (auto const& gen_vertex : gen_vertices) {
0332     if (gen_vertex.eventId().bunchCrossing() != 0)
0333       continue;  // Consider only in-time events
0334     if (gen_vertex.nDaughterTracks() != 2)
0335       continue;  // Keep only V0 vertices
0336     for (TrackingVertex::tp_iterator source = gen_vertex.sourceTracks_begin(); source != gen_vertex.sourceTracks_end();
0337          ++source) {
0338       if (std::abs((*source)->pdgId()) == parent_particle_id) {
0339         if ((std::abs((gen_vertex.daughterTracks().at(0))->pdgId()) == first_daughter_id &&
0340              std::abs((gen_vertex.daughterTracks().at(1))->pdgId()) == second_daughter_id) ||
0341             (std::abs((gen_vertex.daughterTracks().at(0))->pdgId()) == second_daughter_id &&
0342              std::abs((gen_vertex.daughterTracks().at(1))->pdgId()) == first_daughter_id)) {
0343           if ((std::abs((gen_vertex.daughterTracks().at(0))->momentum().eta()) < 2.4 &&
0344                gen_vertex.daughterTracks().at(0)->pt() > 0.9) &&
0345               (std::abs((gen_vertex.daughterTracks().at(1))->momentum().eta()) < 2.4 &&
0346                gen_vertex.daughterTracks().at(1)->pt() > 0.9)) {
0347             // found desired generated Candidate
0348             float candidateGenpT = sqrt((*source)->momentum().perp2());
0349             float candidateGenEta = (*source)->momentum().eta();
0350             float candidateGenR = sqrt((*source)->vertex().perp2());
0351             candidateEffVsPt_denom_[v0_type]->Fill(candidateGenpT);
0352             candidateEffVsEta_denom_[v0_type]->Fill(candidateGenEta);
0353             candidateEffVsR_denom_[v0_type]->Fill(candidateGenR);
0354 
0355             std::array<reco::TrackRef, 2> reco_daughter;
0356 
0357             for (unsigned int daughter = 0; daughter < 2; ++daughter) {
0358               if (simtorecoCollection.find(gen_vertex.daughterTracks()[daughter]) != simtorecoCollection.end()) {
0359                 if (!simtorecoCollection[gen_vertex.daughterTracks()[daughter]].empty()) {
0360                   candidateEff[daughter] = 1;  // Found a daughter track
0361                   reco_daughter[daughter] =
0362                       simtorecoCollection[gen_vertex.daughterTracks()[daughter]].begin()->first.castTo<reco::TrackRef>();
0363                 }
0364               } else {
0365                 candidateEff[daughter] = 2;  // First daughter not found
0366               }
0367             }
0368             if ((candidateEff[0] == 1 && candidateEff[1] == 1) && (reco_daughter[0].key() != reco_daughter[1].key()) &&
0369                 (reconstructed_V0_couples.find(V0Couple(reco_daughter[0], reco_daughter[1])) !=
0370                  reconstructed_V0_couples.end())) {
0371               candidateEffVsPt_num_[v0_type]->Fill(candidateGenpT);
0372               candidateEffVsEta_num_[v0_type]->Fill(candidateGenEta);
0373               candidateEffVsR_num_[v0_type]->Fill(candidateGenR);
0374             }
0375           }  // Check that daughters are inside the desired kinematic region
0376         }    // Check decay products of the current generatex vertex
0377       }      // Check pdgId of the source of the current generated vertex
0378     }        // Loop over all sources of the current generated vertex
0379   }          // Loop over all generated vertices
0380 }
0381 
0382 void V0Validator::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) {
0383   using std::cout;
0384   using std::endl;
0385   using namespace edm;
0386   using namespace std;
0387 
0388   // Make matching collections
0389   Handle<reco::RecoToSimCollection> recotosimCollectionH;
0390   iEvent.getByToken(recoRecoToSimCollectionToken_, recotosimCollectionH);
0391 
0392   Handle<reco::SimToRecoCollection> simtorecoCollectionH;
0393   iEvent.getByToken(recoSimToRecoCollectionToken_, simtorecoCollectionH);
0394 
0395   // Get Monte Carlo information
0396   edm::Handle<TrackingVertexCollection> TVCollectionH;
0397   iEvent.getByToken(trackingVertexCollection_Token_, TVCollectionH);
0398 
0399   // Select the primary vertex, create a new reco::Vertex to hold it
0400   edm::Handle<std::vector<reco::Vertex> > primaryVtxCollectionH;
0401   iEvent.getByToken(vec_recoVertex_Token_, primaryVtxCollectionH);
0402 
0403   std::vector<reco::Vertex>::const_iterator iVtxPH = primaryVtxCollectionH->begin();
0404   for (std::vector<reco::Vertex>::const_iterator iVtx = primaryVtxCollectionH->begin();
0405        iVtx < primaryVtxCollectionH->end();
0406        iVtx++) {
0407     if (primaryVtxCollectionH->size() > 1) {
0408       if (iVtx->tracksSize() > iVtxPH->tracksSize()) {
0409         iVtxPH = iVtx;
0410       }
0411     } else
0412       iVtxPH = iVtx;
0413   }
0414 
0415   // get the V0s;
0416   edm::Handle<reco::VertexCompositeCandidateCollection> k0sCollection;
0417   edm::Handle<reco::VertexCompositeCandidateCollection> lambdaCollection;
0418   iEvent.getByToken(recoVertexCompositeCandidateCollection_k0s_Token_, k0sCollection);
0419   iEvent.getByToken(recoVertexCompositeCandidateCollection_lambda_Token_, lambdaCollection);
0420 
0421   // Do fake rate and efficiency calculation
0422 
0423   // Get gen vertex collection out of the event, as done in the Vertex
0424   // validation package!!!
0425   if (k0sCollection.isValid()) {
0426     doFakeRates(*k0sCollection.product(), *recotosimCollectionH.product(), V0Type::KSHORT, 310, 3122);
0427     doEfficiencies(*TVCollectionH.product(),
0428                    V0Type::KSHORT,
0429                    310,
0430                    211,
0431                    211,
0432                    *k0sCollection.product(),
0433                    *simtorecoCollectionH.product());
0434   }
0435   if (lambdaCollection.isValid()) {
0436     doFakeRates(*lambdaCollection.product(), *recotosimCollectionH.product(), V0Type::LAMBDA, 3122, 310);
0437     doEfficiencies(*TVCollectionH.product(),
0438                    V0Type::LAMBDA,
0439                    3122,
0440                    211,
0441                    2212,
0442                    *lambdaCollection.product(),
0443                    *simtorecoCollectionH.product());
0444   }
0445 }
0446 
0447 // define this as a plug-in
0448 // DEFINE_FWK_MODULE(V0Validator);