MCTruth/test/testReader.cc

0001 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0002 #include "SimMuon/MCTruth/test/testReader.h"
0003
0004 testReader::testReader(const edm::ParameterSet& parset)
0005     : tracksTag_(parset.getParameter<edm::InputTag>("tracksTag")),
0006       tpTag_(parset.getParameter<edm::InputTag>("tpTag")),
0007       assoMapsTag_(parset.getParameter<edm::InputTag>("assoMapsTag")),
0008       tracksToken_(consumes<edm::View<reco::Track>>(tracksTag_)),
0009       tpToken_(consumes<TrackingParticleCollection>(tpTag_)),
0010       recoToSimToken_(consumes<reco::RecoToSimCollection>(assoMapsTag_)),
0011       simToRecoToken_(consumes<reco::SimToRecoCollection>(assoMapsTag_)) {}
0012
0013 void testReader::analyze(const edm::Event& event, const edm::EventSetup& setup) {
0014   LogTrace("testReader") << "testReader::analyze : getting reco::Track collection, " << tracksTag_;
0015   edm::View<reco::Track> trackCollection;
0016   const edm::Handle<edm::View<reco::Track>>& trackCollectionH = event.getHandle(tracksToken_);
0017   if (trackCollectionH.isValid()) {
0018     trackCollection = *(trackCollectionH.product());
0019     LogTrace("testReader") << "... size = " << trackCollection.size();
0020   } else
0021     LogTrace("testReader") << "... NOT FOUND.";
0022
0023   LogTrace("testReader") << "testReader::analyze : getting TrackingParticle collection, " << tpTag_;
0024   TrackingParticleCollection tPC;
0025   const edm::Handle<TrackingParticleCollection>& TPCollectionH = event.getHandle(tpToken_);
0026   if (TPCollectionH.isValid()) {
0027     tPC = *(TPCollectionH.product());
0028     LogTrace("testReader") << "... size = " << tPC.size();
0029   } else
0030     LogTrace("testReader") << "... NOT FOUND.";
0031
0032   LogTrace("testReader") << "testReader::analyze : getting  RecoToSimCollection - " << assoMapsTag_;
0033   reco::RecoToSimCollection recSimColl;
0034   const edm::Handle<reco::RecoToSimCollection>& recSimH = event.getHandle(recoToSimToken_);
0035   if (recSimH.isValid()) {
0036     recSimColl = *(recSimH.product());
0037     LogTrace("testReader") << "... size = " << recSimColl.size();
0038   } else {
0039     LogTrace("testReader") << "... NOT FOUND.";
0040   }
0041
0042   LogTrace("testReader") << "testReader::analyze : getting  SimToRecoCollection - " << assoMapsTag_;
0043   reco::SimToRecoCollection simRecColl;
0044   const edm::Handle<reco::SimToRecoCollection>& simRecH = event.getHandle(simToRecoToken_);
0045   if (simRecH.isValid()) {
0046     simRecColl = *(simRecH.product());
0047     LogTrace("testReader") << "... size = " << simRecColl.size();
0048   } else {
0049     LogTrace("testReader") << "... NOT FOUND.";
0050   }
0051
0052   edm::LogVerbatim("testReader") << "\n === Event ID = " << event.id() << " ===";
0053
0054   // RECOTOSIM
0055   edm::LogVerbatim("testReader") << "\n                      ****************** Reco To Sim "
0056                                     "****************** ";
0057   if (recSimH.isValid()) {
0058     edm::LogVerbatim("testReader") << "\n There are " << trackCollection.size() << " reco::Tracks "
0059                                    << "(" << recSimColl.size() << " matched) \n";
0060
0061     for (edm::View<reco::Track>::size_type i = 0; i < trackCollection.size(); ++i) {
0062       edm::RefToBase<reco::Track> track(trackCollectionH, i);
0063
0064       if (recSimColl.find(track) != recSimColl.end()) {
0065         std::vector<std::pair<TrackingParticleRef, double>> recSimAsso = recSimColl[track];
0066
0067         for (std::vector<std::pair<TrackingParticleRef, double>>::const_iterator IT = recSimAsso.begin();
0068              IT != recSimAsso.end();
0069              ++IT) {
0070           TrackingParticleRef trpart = IT->first;
0071           double purity = IT->second;
0072           edm::LogVerbatim("testReader") << "reco::Track #" << int(i) << " with pt = " << track->pt()
0073                                          << " associated to TrackingParticle #" << trpart.key()
0074                                          << " (pt = " << trpart->pt() << ") with Quality = " << purity;
0075         }
0076       } else {
0077         edm::LogVerbatim("testReader") << "reco::Track #" << int(i) << " with pt = " << track->pt()
0078                                        << " NOT associated to any TrackingParticle"
0079                                        << "\n";
0080       }
0081     }
0082
0083   } else
0084     edm::LogVerbatim("testReader") << "\n RtS map not found in the Event.";
0085
0086   // SIMTORECO
0087   edm::LogVerbatim("testReader") << "\n                      ****************** Sim To Reco "
0088                                     "****************** ";
0089   if (simRecH.isValid()) {
0090     edm::LogVerbatim("testReader") << "\n There are " << tPC.size() << " TrackingParticles "
0091                                    << "(" << simRecColl.size() << " matched) \n";
0092     bool any_trackingParticle_matched = false;
0093
0094     for (TrackingParticleCollection::size_type i = 0; i < tPC.size(); i++) {
0095       TrackingParticleRef trpart(TPCollectionH, i);
0096
0097       std::vector<std::pair<edm::RefToBase<reco::Track>, double>> simRecAsso;
0098
0099       if (simRecColl.find(trpart) != simRecColl.end()) {
0100         simRecAsso = (std::vector<std::pair<edm::RefToBase<reco::Track>, double>>)simRecColl[trpart];
0101
0102         for (std::vector<std::pair<edm::RefToBase<reco::Track>, double>>::const_iterator IT = simRecAsso.begin();
0103              IT != simRecAsso.end();
0104              ++IT) {
0105           edm::RefToBase<reco::Track> track = IT->first;
0106           double quality = IT->second;
0107           any_trackingParticle_matched = true;
0108
0109           // find the purity from RecoToSim association (set purity = -1 for
0110           // unmatched recoToSim)
0111           double purity = -1.;
0112           if (recSimColl.find(track) != recSimColl.end()) {
0113             std::vector<std::pair<TrackingParticleRef, double>> recSimAsso = recSimColl[track];
0114             for (std::vector<std::pair<TrackingParticleRef, double>>::const_iterator ITS = recSimAsso.begin();
0115                  ITS != recSimAsso.end();
0116                  ++ITS) {
0117               TrackingParticleRef tp = ITS->first;
0118               if (tp == trpart)
0119                 purity = ITS->second;
0120             }
0121           }
0122
0123           edm::LogVerbatim("testReader") << "TrackingParticle #" << int(i) << " with pt = " << trpart->pt()
0124                                          << " associated to reco::Track #" << track.key() << " (pt = " << track->pt()
0125                                          << ") with Quality = " << quality << " and Purity = " << purity;
0126         }
0127       }
0128
0129       //    else
0130       //      {
0131       //    LogTrace("testReader") << "TrackingParticle #" << int(i)<< "
0132       // with pt = " << trpart->pt()
0133       //                   << " NOT associated to any reco::Track" ;
0134       //      }
0135     }
0136
0137     if (!any_trackingParticle_matched) {
0138       edm::LogVerbatim("testReader") << "NO TrackingParticle associated to ANY input reco::Track !"
0139                                      << "\n";
0140     }
0141
0142   } else
0143     edm::LogVerbatim("testReader") << "\n StR map not found in the Event.";
0144 }
0145 #include "FWCore/Framework/interface/MakerMacros.h"
0146 DEFINE_FWK_MODULE(testReader);