SiPixelHeterogeneous/plugins/SiPixelCompareVertexSoAAlpaka.cc

0001 // -*- C++ -*-
0002 // Package:    SiPixelCompareVertexSoAAlpaka
0003 // Class:      SiPixelCompareVertexSoAAlpaka
0004 //
0005 /**\class SiPixelCompareVertexSoAAlpaka SiPixelCompareVertexSoAAlpaka.cc
0006 */
0007 //
0008 // Author: Suvankar Roy Chowdhury
0009 //
0010 #include "FWCore/Framework/interface/Frameworkfwd.h"
0011 #include "FWCore/Framework/interface/Event.h"
0012 #include "FWCore/Framework/interface/MakerMacros.h"
0013 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0014 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0015 #include "FWCore/Utilities/interface/InputTag.h"
0016 #include "DataFormats/Common/interface/Handle.h"
0017 // DQM Histograming
0018 #include "DQMServices/Core/interface/MonitorElement.h"
0019 #include "DQMServices/Core/interface/DQMEDAnalyzer.h"
0020 #include "DQMServices/Core/interface/DQMStore.h"
0021 #include "DataFormats/VertexSoA/interface/ZVertexHost.h"
0022 #include "DataFormats/BeamSpot/interface/BeamSpot.h"
0023
0024 class SiPixelCompareVertexSoAAlpaka : public DQMEDAnalyzer {
0025 public:
0026   using IndToEdm = std::vector<uint16_t>;
0027   explicit SiPixelCompareVertexSoAAlpaka(const edm::ParameterSet&);
0028   ~SiPixelCompareVertexSoAAlpaka() override = default;
0029   void bookHistograms(DQMStore::IBooker& ibooker, edm::Run const& iRun, edm::EventSetup const& iSetup) override;
0030   void analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) override;
0031   static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0032
0033 private:
0034   const edm::EDGetTokenT<ZVertexHost> tokenSoAVertexHost_;
0035   const edm::EDGetTokenT<ZVertexHost> tokenSoAVertexDevice_;
0036   const edm::EDGetTokenT<reco::BeamSpot> tokenBeamSpot_;
0037   const std::string topFolderName_;
0038   const float dzCut_;
0039   MonitorElement* hnVertex_;
0040   MonitorElement* hx_;
0041   MonitorElement* hy_;
0042   MonitorElement* hz_;
0043   MonitorElement* hchi2_;
0044   MonitorElement* hchi2oNdof_;
0045   MonitorElement* hptv2_;
0046   MonitorElement* hntrks_;
0047   MonitorElement* hxdiff_;
0048   MonitorElement* hydiff_;
0049   MonitorElement* hzdiff_;
0050 };
0051
0052 //
0053 // constructors
0054 //
0055
0056 // Note tokenSoAVertexDevice_ contains data copied from device to host, hence is a HostCollection
0057 SiPixelCompareVertexSoAAlpaka::SiPixelCompareVertexSoAAlpaka(const edm::ParameterSet& iConfig)
0058     : tokenSoAVertexHost_(consumes<ZVertexHost>(iConfig.getParameter<edm::InputTag>("pixelVertexSrcHost"))),
0059       tokenSoAVertexDevice_(consumes<ZVertexHost>(iConfig.getParameter<edm::InputTag>("pixelVertexSrcDevice"))),
0060       tokenBeamSpot_(consumes<reco::BeamSpot>(iConfig.getParameter<edm::InputTag>("beamSpotSrc"))),
0061       topFolderName_(iConfig.getParameter<std::string>("topFolderName")),
0062       dzCut_(iConfig.getParameter<double>("dzCut")) {}
0063
0064 //
0065 // -- Analyze
0066 //
0067 void SiPixelCompareVertexSoAAlpaka::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) {
0068   const auto& vsoaHandleHost = iEvent.getHandle(tokenSoAVertexHost_);
0069   const auto& vsoaHandleDevice = iEvent.getHandle(tokenSoAVertexDevice_);
0070   if (not vsoaHandleHost or not vsoaHandleDevice) {
0071     edm::LogWarning out("SiPixelCompareVertexSoAAlpaka");
0072     if (not vsoaHandleHost) {
0073       out << "reference (cpu) tracks not found; ";
0074     }
0075     if (not vsoaHandleDevice) {
0076       out << "target (gpu) tracks not found; ";
0077     }
0078     out << "the comparison will not run.";
0079     return;
0080   }
0081
0082   auto const& vsoaHost = *vsoaHandleHost;
0083   int nVerticesHost = vsoaHost.view().nvFinal();
0084   auto const& vsoaDevice = *vsoaHandleDevice;
0085   int nVerticesDevice = vsoaDevice.view().nvFinal();
0086
0087   auto bsHandle = iEvent.getHandle(tokenBeamSpot_);
0088   float x0 = 0., y0 = 0., z0 = 0., dxdz = 0., dydz = 0.;
0089   if (!bsHandle.isValid()) {
0090     edm::LogWarning("SiPixelCompareVertexSoAAlpaka") << "No beamspot found. returning vertexes with (0,0,Z) ";
0091   } else {
0092     const reco::BeamSpot& bs = *bsHandle;
0093     x0 = bs.x0();
0094     y0 = bs.y0();
0095     z0 = bs.z0();
0096     dxdz = bs.dxdz();
0097     dydz = bs.dydz();
0098   }
0099
0100   for (int ivc = 0; ivc < nVerticesHost; ivc++) {
0101     auto sic = vsoaHost.view()[ivc].sortInd();
0102     auto zc = vsoaHost.view()[sic].zv();
0103     auto xc = x0 + dxdz * zc;
0104     auto yc = y0 + dydz * zc;
0105     zc += z0;
0106
0107     auto ndofHost = vsoaHost.view()[sic].ndof();
0108     auto chi2Host = vsoaHost.view()[sic].chi2();
0109
0110     const int32_t notFound = -1;
0111     int32_t closestVtxidx = notFound;
0112     float mindz = dzCut_;
0113
0114     for (int ivg = 0; ivg < nVerticesDevice; ivg++) {
0115       auto sig = vsoaDevice.view()[ivg].sortInd();
0116       auto zgc = vsoaDevice.view()[sig].zv() + z0;
0117       auto zDist = std::abs(zc - zgc);
0118       //insert some matching condition
0119       if (zDist > dzCut_)
0120         continue;
0121       if (mindz > zDist) {
0122         mindz = zDist;
0123         closestVtxidx = sig;
0124       }
0125     }
0126     if (closestVtxidx == notFound)
0127       continue;
0128
0129     auto zg = vsoaDevice.view()[closestVtxidx].zv();
0130     auto xg = x0 + dxdz * zg;
0131     auto yg = y0 + dydz * zg;
0132     zg += z0;
0133     auto ndofDevice = vsoaDevice.view()[closestVtxidx].ndof();
0134     auto chi2Device = vsoaDevice.view()[closestVtxidx].chi2();
0135
0136     hx_->Fill(xc - x0, xg - x0);
0137     hy_->Fill(yc - y0, yg - y0);
0138     hz_->Fill(zc, zg);
0139     hxdiff_->Fill(xc - xg);
0140     hydiff_->Fill(yc - yg);
0141     hzdiff_->Fill(zc - zg);
0142     hchi2_->Fill(chi2Host, chi2Device);
0143     hchi2oNdof_->Fill(chi2Host / ndofHost, chi2Device / ndofDevice);
0144     hptv2_->Fill(vsoaHost.view()[sic].ptv2(), vsoaDevice.view()[closestVtxidx].ptv2());
0145     hntrks_->Fill(ndofHost + 1, ndofDevice + 1);
0146   }
0147   hnVertex_->Fill(nVerticesHost, nVerticesDevice);
0148 }
0149
0150 //
0151 // -- Book Histograms
0152 //
0153 void SiPixelCompareVertexSoAAlpaka::bookHistograms(DQMStore::IBooker& ibooker,
0154                                                    edm::Run const& iRun,
0155                                                    edm::EventSetup const& iSetup) {
0156   ibooker.cd();
0157   ibooker.setCurrentFolder(topFolderName_);
0158
0159   // FIXME: all the 2D correlation plots are quite heavy in terms of memory consumption, so a as soon as DQM supports either TH2I or THnSparse
0160   // these should be moved to a less resource consuming format
0161   hnVertex_ = ibooker.book2I("nVertex", "# of Vertices;Host;Device", 101, -0.5, 100.5, 101, -0.5, 100.5);
0162   hx_ = ibooker.book2I("vx", "Vertez x - Beamspot x;Host;Device", 50, -0.1, 0.1, 50, -0.1, 0.1);
0163   hy_ = ibooker.book2I("vy", "Vertez y - Beamspot y;Host;Device", 50, -0.1, 0.1, 50, -0.1, 0.1);
0164   hz_ = ibooker.book2I("vz", "Vertez z;Host;Device", 30, -30., 30., 30, -30., 30.);
0165   hchi2_ = ibooker.book2I("chi2", "Vertex chi-squared;Host;Device", 40, 0., 20., 40, 0., 20.);
0166   hchi2oNdof_ = ibooker.book2I("chi2oNdof", "Vertex chi-squared/Ndof;Host;Device", 40, 0., 20., 40, 0., 20.);
0167   hptv2_ = ibooker.book2I("ptsq", "Vertex #sum (p_{T})^{2};Host;Device", 200, 0., 200., 200, 0., 200.);
0168   hntrks_ = ibooker.book2I("ntrk", "#tracks associated;Host;Device", 100, -0.5, 99.5, 100, -0.5, 99.5);
0169   hntrks_ = ibooker.book2I("ntrk", "#tracks associated;Host;Device", 100, -0.5, 99.5, 100, -0.5, 99.5);
0170   hxdiff_ = ibooker.book1D("vxdiff", ";Vertex x difference (Host - Device);#entries", 100, -0.001, 0.001);
0171   hydiff_ = ibooker.book1D("vydiff", ";Vertex y difference (Host - Device);#entries", 100, -0.001, 0.001);
0172   hzdiff_ = ibooker.book1D("vzdiff", ";Vertex z difference (Host - Device);#entries", 100, -2.5, 2.5);
0173 }
0174
0175 void SiPixelCompareVertexSoAAlpaka::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0176   // monitorpixelVertexSoA
0177   edm::ParameterSetDescription desc;
0178   desc.add<edm::InputTag>("pixelVertexSrcHost", edm::InputTag("pixelVerticesAlpakaSerial"));
0179   desc.add<edm::InputTag>("pixelVertexSrcDevice", edm::InputTag("pixelVerticesAlpaka"));
0180   desc.add<edm::InputTag>("beamSpotSrc", edm::InputTag("offlineBeamSpot"));
0181   desc.add<std::string>("topFolderName", "SiPixelHeterogeneous/PixelVertexCompareSoADeviceVSHost");
0182   desc.add<double>("dzCut", 1.);
0183   descriptions.addWithDefaultLabel(desc);
0184 }
0185
0186 DEFINE_FWK_MODULE(SiPixelCompareVertexSoAAlpaka);