Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2024-08-30 02:10:35

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