Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2022-06-04 01:25:53

0001 // -*- C++ -*-
0002 ///bookLayer
0003 // Package:    SiPixelPhase1MonitorTrackSoA
0004 // Class:      SiPixelPhase1MonitorTrackSoA
0005 //
0006 /**\class SiPixelPhase1MonitorTrackSoA SiPixelPhase1MonitorTrackSoA.cc 
0007 */
0008 //
0009 // Author: Suvankar Roy Chowdhury
0010 //
0011 #include "DataFormats/Common/interface/Handle.h"
0012 #include "FWCore/Framework/interface/ESHandle.h"
0013 #include "FWCore/Framework/interface/Event.h"
0014 #include "FWCore/Framework/interface/Frameworkfwd.h"
0015 #include "FWCore/Framework/interface/MakerMacros.h"
0016 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0017 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0018 #include "FWCore/ServiceRegistry/interface/Service.h"
0019 #include "FWCore/Utilities/interface/InputTag.h"
0020 // DQM Histograming
0021 #include "DQMServices/Core/interface/MonitorElement.h"
0022 #include "DQMServices/Core/interface/DQMEDAnalyzer.h"
0023 #include "DQMServices/Core/interface/DQMStore.h"
0024 #include "CUDADataFormats/Track/interface/PixelTrackHeterogeneous.h"
0025 // for string manipulations
0026 #include <fmt/printf.h>
0027 
0028 class SiPixelPhase1MonitorTrackSoA : public DQMEDAnalyzer {
0029 public:
0030   explicit SiPixelPhase1MonitorTrackSoA(const edm::ParameterSet&);
0031   ~SiPixelPhase1MonitorTrackSoA() override = default;
0032   void bookHistograms(DQMStore::IBooker& ibooker, edm::Run const& iRun, edm::EventSetup const& iSetup) override;
0033   void analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) override;
0034   static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0035 
0036 private:
0037   edm::EDGetTokenT<PixelTrackHeterogeneous> tokenSoATrack_;
0038   std::string topFolderName_;
0039   bool useQualityCut_;
0040   pixelTrack::Quality minQuality_;
0041   MonitorElement* hnTracks;
0042   MonitorElement* hnLooseAndAboveTracks;
0043   MonitorElement* hnHits;
0044   MonitorElement* hnHitsVsPhi;
0045   MonitorElement* hnHitsVsEta;
0046   MonitorElement* hnLayers;
0047   MonitorElement* hnLayersVsPhi;
0048   MonitorElement* hnLayersVsEta;
0049   MonitorElement* hchi2;
0050   MonitorElement* hChi2VsPhi;
0051   MonitorElement* hChi2VsEta;
0052   MonitorElement* hpt;
0053   MonitorElement* heta;
0054   MonitorElement* hphi;
0055   MonitorElement* hz;
0056   MonitorElement* htip;
0057   MonitorElement* hquality;
0058 };
0059 
0060 //
0061 // constructors
0062 //
0063 
0064 SiPixelPhase1MonitorTrackSoA::SiPixelPhase1MonitorTrackSoA(const edm::ParameterSet& iConfig) {
0065   tokenSoATrack_ = consumes<PixelTrackHeterogeneous>(iConfig.getParameter<edm::InputTag>("pixelTrackSrc"));
0066   topFolderName_ = iConfig.getParameter<std::string>("topFolderName");  //"SiPixelHeterogeneous/PixelTrackSoA";
0067   useQualityCut_ = iConfig.getParameter<bool>("useQualityCut");
0068   minQuality_ = pixelTrack::qualityByName(iConfig.getParameter<std::string>("minQuality"));
0069 }
0070 
0071 //
0072 // -- Analyze
0073 //
0074 void SiPixelPhase1MonitorTrackSoA::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) {
0075   const auto& tsoaHandle = iEvent.getHandle(tokenSoATrack_);
0076   if (!tsoaHandle.isValid()) {
0077     edm::LogWarning("SiPixelPhase1MonitorTrackSoA") << "No Track SoA found \n returning!" << std::endl;
0078     return;
0079   }
0080 
0081   auto const& tsoa = *((tsoaHandle.product())->get());
0082   auto maxTracks = tsoa.stride();
0083   auto const* quality = tsoa.qualityData();
0084   int32_t nTracks = 0;
0085   int32_t nLooseAndAboveTracks = 0;
0086 
0087   for (int32_t it = 0; it < maxTracks; ++it) {
0088     auto nHits = tsoa.nHits(it);
0089     auto nLayers = tsoa.nLayers(it);
0090     if (nHits == 0)
0091       break;  // this is a guard
0092     float pt = tsoa.pt(it);
0093     if (!(pt > 0.))
0094       continue;
0095 
0096     // fill the quality for all tracks
0097     pixelTrack::Quality qual = tsoa.quality(it);
0098     hquality->Fill(int(qual));
0099     nTracks++;
0100 
0101     if (useQualityCut_ && quality[it] < minQuality_)
0102       continue;
0103 
0104     // fill parameters only for quality >= loose
0105     float chi2 = tsoa.chi2(it);
0106     float phi = tsoa.phi(it);
0107     float zip = tsoa.zip(it);
0108     float eta = tsoa.eta(it);
0109     float tip = tsoa.tip(it);
0110 
0111     hchi2->Fill(chi2);
0112     hChi2VsPhi->Fill(phi, chi2);
0113     hChi2VsEta->Fill(eta, chi2);
0114     hnHits->Fill(nHits);
0115     hnLayers->Fill(nLayers);
0116     hnHitsVsPhi->Fill(phi, nHits);
0117     hnHitsVsEta->Fill(eta, nHits);
0118     hnLayersVsPhi->Fill(phi, nLayers);
0119     hnLayersVsEta->Fill(eta, nLayers);
0120     hpt->Fill(pt);
0121     heta->Fill(eta);
0122     hphi->Fill(phi);
0123     hz->Fill(zip);
0124     htip->Fill(tip);
0125     nLooseAndAboveTracks++;
0126   }
0127   hnTracks->Fill(nTracks);
0128   hnLooseAndAboveTracks->Fill(nLooseAndAboveTracks);
0129 }
0130 
0131 //
0132 // -- Book Histograms
0133 //
0134 void SiPixelPhase1MonitorTrackSoA::bookHistograms(DQMStore::IBooker& iBook,
0135                                                   edm::Run const& iRun,
0136                                                   edm::EventSetup const& iSetup) {
0137   iBook.cd();
0138   iBook.setCurrentFolder(topFolderName_);
0139 
0140   // clang-format off
0141   std::string toRep = "Number of tracks";
0142   hnTracks = iBook.book1D("nTracks", fmt::sprintf(";%s per event;#events",toRep), 1001, -0.5, 1000.5);
0143   hnLooseAndAboveTracks = iBook.book1D("nLooseAndAboveTracks", fmt::sprintf(";%s (quality #geq loose) per event;#events",toRep), 1001, -0.5, 1000.5);
0144 
0145   toRep = "Number of all RecHits per track (quality #geq loose)";
0146   hnHits = iBook.book1D("nRecHits", fmt::sprintf(";%s;#tracks",toRep), 15, -0.5, 14.5);
0147   hnHitsVsPhi = iBook.bookProfile("nHitsPerTrackVsPhi", fmt::sprintf("%s vs track #phi;Track #phi;%s",toRep,toRep), 30, -M_PI, M_PI,0., 15.);
0148   hnHitsVsEta = iBook.bookProfile("nHitsPerTrackVsEta", fmt::sprintf("%s vs track #eta;Track #eta;%s",toRep,toRep), 30, -3., 3., 0., 15.);
0149 
0150   toRep = "Number of all layers per track (quality #geq loose)";
0151   hnLayers = iBook.book1D("nLayers", fmt::sprintf(";%s;#tracks",toRep), 15, -0.5, 14.5);
0152   hnLayersVsPhi = iBook.bookProfile("nLayersPerTrackVsPhi", fmt::sprintf("%s vs track #phi;Track #phi;%s",toRep,toRep), 30, -M_PI, M_PI,0., 15.);
0153   hnLayersVsEta = iBook.bookProfile("nLayersPerTrackVsEta", fmt::sprintf("%s vs track #eta;Track #eta;%s",toRep,toRep), 30, -3., 3., 0., 15.);
0154 
0155   toRep = "Track (quality #geq loose) #chi^{2}/ndof";
0156   hchi2 = iBook.book1D("nChi2ndof", fmt::sprintf(";%s;#tracks",toRep), 40, 0., 20.);
0157   hChi2VsPhi = iBook.bookProfile("nChi2ndofVsPhi", fmt::sprintf("%s vs track #phi;Track #phi;%s",toRep,toRep), 30, -M_PI, M_PI, 0., 20.);
0158   hChi2VsEta = iBook.bookProfile("nChi2ndofVsEta", fmt::sprintf("%s vs track #eta;Track #eta;%s",toRep,toRep), 30, -3., 3., 0., 20.);
0159   // clang-format on
0160 
0161   hpt = iBook.book1D("pt", ";Track (quality #geq loose) p_{T} [GeV];#tracks", 200, 0., 200.);
0162   heta = iBook.book1D("eta", ";Track (quality #geq loose) #eta;#tracks", 30, -3., 3.);
0163   hphi = iBook.book1D("phi", ";Track (quality #geq loose) #phi;#tracks", 30, -M_PI, M_PI);
0164   hz = iBook.book1D("z", ";Track (quality #geq loose) z [cm];#tracks", 30, -30., 30.);
0165   htip = iBook.book1D("tip", ";Track (quality #geq loose) TIP [cm];#tracks", 100, -0.5, 0.5);
0166   hquality = iBook.book1D("quality", ";Track Quality;#tracks", 7, -0.5, 6.5);
0167   uint i = 1;
0168   for (const auto& q : pixelTrack::qualityName) {
0169     hquality->setBinLabel(i, q, 1);
0170     i++;
0171   }
0172 }
0173 
0174 void SiPixelPhase1MonitorTrackSoA::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0175   // monitorpixelTrackSoA
0176   edm::ParameterSetDescription desc;
0177   desc.add<edm::InputTag>("pixelTrackSrc", edm::InputTag("pixelTracksSoA"));
0178   desc.add<std::string>("topFolderName", "SiPixelHeterogeneous/PixelTrackSoA");
0179   desc.add<bool>("useQualityCut", true);
0180   desc.add<std::string>("minQuality", "loose");
0181   descriptions.addWithDefaultLabel(desc);
0182 }
0183 DEFINE_FWK_MODULE(SiPixelPhase1MonitorTrackSoA);