SiPixelPhase1HitsV/src/SiPixelPhase1HitsV.cc

0001 // -*- C++ -*-
0002 //
0003 // Package:     SiPixelPhase1HitsV
0004 // Class:       SiPixelPhase1HitsV
0005 //
0006
0007 // Original Author: Marcel Schneider
0008 // Additional Authors: Alexander Morton - modifying code for validation use
0009
0010 #include "FWCore/Framework/interface/MakerMacros.h"
0011 #include "FWCore/Utilities/interface/Exception.h"
0012 #include "Validation/SiPixelPhase1HitsV/interface/SiPixelPhase1HitsV.h"
0013
0014 #include "DataFormats/GeometryVector/interface/LocalPoint.h"
0015 #include "FWCore/Framework/interface/ESHandle.h"
0016
0017 #include "Geometry/CommonDetUnit/interface/TrackingGeometry.h"
0018
0019 #include "SimDataFormats/TrackingHit/interface/PSimHit.h"
0020 #include "SimDataFormats/TrackingHit/interface/PSimHitContainer.h"
0021
0022 // class TrackAssociatorByHits;
0023
0024 SiPixelPhase1HitsV::SiPixelPhase1HitsV(const edm::ParameterSet &iConfig)
0025     : SiPixelPhase1Base(iConfig),
0026       pixelBarrelLowToken_(consumes<edm::PSimHitContainer>(iConfig.getParameter<edm::InputTag>("pixBarrelLowSrc"))),
0027       pixelBarrelHighToken_(consumes<edm::PSimHitContainer>(iConfig.getParameter<edm::InputTag>("pixBarrelHighSrc"))),
0028       pixelForwardLowToken_(consumes<edm::PSimHitContainer>(iConfig.getParameter<edm::InputTag>("pixForwardLowSrc"))),
0029       pixelForwardHighToken_(consumes<edm::PSimHitContainer>(iConfig.getParameter<edm::InputTag>("pixForwardHighSrc"))),
0030
0031       tracksToken_(consumes<edm::View<reco::Track>>(iConfig.getParameter<edm::InputTag>("tracksTag"))),
0032       tpToken_(consumes<TrackingParticleCollection>(iConfig.getParameter<edm::InputTag>("tpTag"))),
0033       trackAssociatorByHitsToken_(consumes<reco::TrackToTrackingParticleAssociator>(
0034           iConfig.getParameter<edm::InputTag>("trackAssociatorByHitsTag"))),
0035
0036       trackerGeomToken_(esConsumes<TrackerGeometry, TrackerDigiGeometryRecord>()) {}
0037
0038 void SiPixelPhase1HitsV::analyze(const edm::Event &iEvent, const edm::EventSetup &iSetup) {
0039   edm::Handle<edm::PSimHitContainer> barrelLowInput;
0040   iEvent.getByToken(pixelBarrelLowToken_, barrelLowInput);
0041   if (!barrelLowInput.isValid())
0042     return;
0043
0044   edm::Handle<edm::PSimHitContainer> barrelHighInput;
0045   iEvent.getByToken(pixelBarrelHighToken_, barrelHighInput);
0046   if (!barrelHighInput.isValid())
0047     return;
0048
0049   edm::Handle<edm::PSimHitContainer> forwardLowInput;
0050   iEvent.getByToken(pixelForwardLowToken_, forwardLowInput);
0051   if (!forwardLowInput.isValid())
0052     return;
0053
0054   edm::Handle<edm::PSimHitContainer> forwardHighInput;
0055   iEvent.getByToken(pixelForwardHighToken_, forwardHighInput);
0056   if (!forwardHighInput.isValid())
0057     return;
0058
0059   edm::PSimHitContainer::const_iterator it;
0060
0061   // Get geometry information
0062
0063   edm::ESHandle<TrackerGeometry> tracker = iSetup.getHandle(trackerGeomToken_);
0064
0065   // get low barrel info
0066   for (it = barrelLowInput->begin(); it != barrelLowInput->end(); ++it) {
0067     auto id = DetId(it->detUnitId());
0068     const GeomDetUnit *det = (const GeomDetUnit *)tracker->idToDetUnit(id);
0069     GlobalPoint gpos = det->toGlobal(it->localPosition());
0070
0071     float tof = it->timeOfFlight();
0072     float globalR = gpos.mag();
0073
0074     float energyLoss = it->energyLoss();
0075
0076     float entryExitX = (it->entryPoint().x() - it->exitPoint().x());
0077     float entryExitY = (it->entryPoint().y() - it->exitPoint().y());
0078     float entryExitZ = std::abs(it->entryPoint().z() - it->exitPoint().z());
0079
0080     float localX = it->localPosition().x();
0081     float localY = it->localPosition().y();
0082     float localZ = it->localPosition().z();
0083     float localPhi = it->localPosition().phi();
0084     float localEta = it->localPosition().eta();
0085
0086     histo[TOF_R].fill(globalR, tof, id, &iEvent);
0087     histo[ELOSS].fill(energyLoss, id, &iEvent);
0088     histo[ENTRY_EXIT_X].fill(entryExitX, id, &iEvent);
0089     histo[ENTRY_EXIT_Y].fill(entryExitY, id, &iEvent);
0090     histo[ENTRY_EXIT_Z].fill(entryExitZ, id, &iEvent);
0091     histo[LOCAL_X].fill(localX, id, &iEvent);
0092     histo[LOCAL_Y].fill(localY, id, &iEvent);
0093     histo[LOCAL_Z].fill(localZ, id, &iEvent);
0094     histo[LOCAL_PHI].fill(localPhi, id, &iEvent);
0095     histo[LOCAL_ETA].fill(localEta, id, &iEvent);
0096   }
0097   // get high barrel info
0098   for (it = barrelHighInput->begin(); it != barrelHighInput->end(); ++it) {
0099     auto id = DetId(it->detUnitId());
0100     const GeomDetUnit *det = (const GeomDetUnit *)tracker->idToDetUnit(id);
0101     GlobalPoint gpos = det->toGlobal(it->localPosition());
0102
0103     float tof = it->timeOfFlight();
0104     float globalR = gpos.mag();
0105
0106     float energyLoss = it->energyLoss();
0107
0108     float entryExitX = (it->entryPoint().x() - it->exitPoint().x());
0109     float entryExitY = (it->entryPoint().y() - it->exitPoint().y());
0110     float entryExitZ = std::abs(it->entryPoint().z() - it->exitPoint().z());
0111
0112     float localX = it->localPosition().x();
0113     float localY = it->localPosition().y();
0114     float localZ = it->localPosition().z();
0115     float localPhi = it->localPosition().phi();
0116     float localEta = it->localPosition().eta();
0117
0118     histo[TOF_R].fill(globalR, tof, id, &iEvent);
0119     histo[ELOSS].fill(energyLoss, id, &iEvent);
0120     histo[ENTRY_EXIT_X].fill(entryExitX, id, &iEvent);
0121     histo[ENTRY_EXIT_Y].fill(entryExitY, id, &iEvent);
0122     histo[ENTRY_EXIT_Z].fill(entryExitZ, id, &iEvent);
0123     histo[LOCAL_X].fill(localX, id, &iEvent);
0124     histo[LOCAL_Y].fill(localY, id, &iEvent);
0125     histo[LOCAL_Z].fill(localZ, id, &iEvent);
0126     histo[LOCAL_PHI].fill(localPhi, id, &iEvent);
0127     histo[LOCAL_ETA].fill(localEta, id, &iEvent);
0128   }
0129
0130   // get low forward info
0131   for (it = forwardLowInput->begin(); it != forwardLowInput->end(); ++it) {
0132     auto id = DetId(it->detUnitId());
0133     const GeomDetUnit *det = (const GeomDetUnit *)tracker->idToDetUnit(id);
0134     GlobalPoint gpos = det->toGlobal(it->localPosition());
0135
0136     float tof = it->timeOfFlight();
0137     float globalR = gpos.mag();
0138
0139     float energyLoss = it->energyLoss();
0140
0141     float entryExitX = (it->entryPoint().x() - it->exitPoint().x());
0142     float entryExitY = (it->entryPoint().y() - it->exitPoint().y());
0143     float entryExitZ = std::abs(it->entryPoint().z() - it->exitPoint().z());
0144
0145     float localX = it->localPosition().x();
0146     float localY = it->localPosition().y();
0147     float localZ = it->localPosition().z();
0148     float localPhi = it->localPosition().phi();
0149     float localEta = it->localPosition().eta();
0150
0151     histo[TOF_R].fill(globalR, tof, id, &iEvent);
0152     histo[ELOSS].fill(energyLoss, id, &iEvent);
0153     histo[ENTRY_EXIT_X].fill(entryExitX, id, &iEvent);
0154     histo[ENTRY_EXIT_Y].fill(entryExitY, id, &iEvent);
0155     histo[ENTRY_EXIT_Z].fill(entryExitZ, id, &iEvent);
0156     histo[LOCAL_X].fill(localX, id, &iEvent);
0157     histo[LOCAL_Y].fill(localY, id, &iEvent);
0158     histo[LOCAL_Z].fill(localZ, id, &iEvent);
0159     histo[LOCAL_PHI].fill(localPhi, id, &iEvent);
0160     histo[LOCAL_ETA].fill(localEta, id, &iEvent);
0161   }
0162
0163   // get high forward info
0164   for (it = forwardHighInput->begin(); it != forwardHighInput->end(); ++it) {
0165     auto id = DetId(it->detUnitId());
0166     const GeomDetUnit *det = (const GeomDetUnit *)tracker->idToDetUnit(id);
0167     GlobalPoint gpos = det->toGlobal(it->localPosition());
0168
0169     float tof = it->timeOfFlight();
0170     float globalR = gpos.mag();
0171
0172     float energyLoss = it->energyLoss();
0173
0174     float entryExitX = (it->entryPoint().x() - it->exitPoint().x());
0175     float entryExitY = (it->entryPoint().y() - it->exitPoint().y());
0176     float entryExitZ = std::abs(it->entryPoint().z() - it->exitPoint().z());
0177
0178     float localX = it->localPosition().x();
0179     float localY = it->localPosition().y();
0180     float localZ = it->localPosition().z();
0181     float localPhi = it->localPosition().phi();
0182     float localEta = it->localPosition().eta();
0183
0184     histo[TOF_R].fill(globalR, tof, id, &iEvent);
0185     histo[ELOSS].fill(energyLoss, id, &iEvent);
0186     histo[ENTRY_EXIT_X].fill(entryExitX, id, &iEvent);
0187     histo[ENTRY_EXIT_Y].fill(entryExitY, id, &iEvent);
0188     histo[ENTRY_EXIT_Z].fill(entryExitZ, id, &iEvent);
0189     histo[LOCAL_X].fill(localX, id, &iEvent);
0190     histo[LOCAL_Y].fill(localY, id, &iEvent);
0191     histo[LOCAL_Z].fill(localZ, id, &iEvent);
0192     histo[LOCAL_PHI].fill(localPhi, id, &iEvent);
0193     histo[LOCAL_ETA].fill(localEta, id, &iEvent);
0194   }
0195
0196   // Sim Hit efficiency info
0197   edm::Handle<edm::View<reco::Track>> trackCollectionH;
0198   iEvent.getByToken(tracksToken_, trackCollectionH);
0199   const edm::View<reco::Track> &tC = *(trackCollectionH.product());
0200
0201   edm::Handle<TrackingParticleCollection> TPCollectionH;
0202   iEvent.getByToken(tpToken_, TPCollectionH);
0203
0204   edm::Handle<reco::TrackToTrackingParticleAssociator> theHitsAssociator;
0205   iEvent.getByToken(trackAssociatorByHitsToken_, theHitsAssociator);
0206   if (!theHitsAssociator.isValid()) {
0207     throw cms::Exception("NO VALID HIT ASSOCIATOR");
0208   }
0209   reco::TrackToTrackingParticleAssociator const *associatorByHits = theHitsAssociator.product();
0210
0211   if (TPCollectionH.isValid() && trackCollectionH.isValid()) {
0212     reco::RecoToSimCollection const &p = associatorByHits->associateRecoToSim(trackCollectionH, TPCollectionH);
0213
0214     for (edm::View<reco::Track>::size_type i = 0; i < tC.size(); ++i) {
0215       edm::RefToBase<reco::Track> track(trackCollectionH, i);
0216       //      const reco::Track& t = *track;
0217       auto id = DetId(track->innerDetId());  // histo manager requires a det ID,
0218                                              // use innermost ID for ease
0219
0220       auto iter = p.find(track);
0221       histo[EFFICIENCY_TRACK].fill(iter != p.end() ? 1 : 0, id, &iEvent);
0222     }
0223   }
0224 }
0225
0226 DEFINE_FWK_MODULE(SiPixelPhase1HitsV);