File indexing completed on 2024-09-13 22:52:25
0001
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011 #include "FWCore/Framework/interface/Frameworkfwd.h"
0012 #include "FWCore/Framework/interface/Event.h"
0013 #include "FWCore/Framework/interface/ESHandle.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 #include "DataFormats/Common/interface/Handle.h"
0020
0021 #include "DQMServices/Core/interface/MonitorElement.h"
0022 #include "DQMServices/Core/interface/DQMEDAnalyzer.h"
0023 #include "DQMServices/Core/interface/DQMStore.h"
0024 #include "DataFormats/VertexSoA/interface/ZVertexHost.h"
0025 #include "DataFormats/BeamSpot/interface/BeamSpot.h"
0026
0027 class SiPixelMonitorVertexSoAAlpaka : public DQMEDAnalyzer {
0028 public:
0029 using IndToEdm = std::vector<uint16_t>;
0030 explicit SiPixelMonitorVertexSoAAlpaka(const edm::ParameterSet&);
0031 ~SiPixelMonitorVertexSoAAlpaka() 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 const edm::EDGetTokenT<ZVertexHost> tokenSoAVertex_;
0038 const edm::EDGetTokenT<reco::BeamSpot> tokenBeamSpot_;
0039 std::string topFolderName_;
0040 MonitorElement* hnVertex;
0041 MonitorElement* hx;
0042 MonitorElement* hy;
0043 MonitorElement* hz;
0044 MonitorElement* hchi2;
0045 MonitorElement* hchi2oNdof;
0046 MonitorElement* hptv2;
0047 MonitorElement* hntrks;
0048 };
0049
0050
0051
0052
0053
0054 SiPixelMonitorVertexSoAAlpaka::SiPixelMonitorVertexSoAAlpaka(const edm::ParameterSet& iConfig)
0055 : tokenSoAVertex_(consumes<ZVertexHost>(iConfig.getParameter<edm::InputTag>("pixelVertexSrc"))),
0056 tokenBeamSpot_(consumes<reco::BeamSpot>(iConfig.getParameter<edm::InputTag>("beamSpotSrc"))),
0057 topFolderName_(iConfig.getParameter<std::string>("topFolderName")) {}
0058
0059
0060
0061
0062 void SiPixelMonitorVertexSoAAlpaka::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) {
0063 const auto& vsoaHandle = iEvent.getHandle(tokenSoAVertex_);
0064 if (!vsoaHandle.isValid()) {
0065 edm::LogWarning("SiPixelMonitorVertexSoAAlpaka") << "No Vertex SoA found \n returning!" << std::endl;
0066 return;
0067 }
0068
0069 auto const& vsoa = *vsoaHandle;
0070 auto vtx_view = vsoa.view<reco::ZVertexSoA>();
0071 auto trk_view = vsoa.view<reco::ZVertexTracksSoA>();
0072 int nVertices = vtx_view.nvFinal();
0073 auto bsHandle = iEvent.getHandle(tokenBeamSpot_);
0074 float x0 = 0., y0 = 0., z0 = 0., dxdz = 0., dydz = 0.;
0075 if (!bsHandle.isValid()) {
0076 edm::LogWarning("SiPixelMonitorVertexSoAAlpaka") << "No beamspot found. returning vertexes with (0,0,Z) ";
0077 } else {
0078 const reco::BeamSpot& bs = *bsHandle;
0079 x0 = bs.x0();
0080 y0 = bs.y0();
0081 z0 = bs.z0();
0082 dxdz = bs.dxdz();
0083 dydz = bs.dydz();
0084 }
0085
0086 for (int iv = 0; iv < nVertices; iv++) {
0087 auto si = vtx_view[iv].sortInd();
0088 auto z = vtx_view[si].zv();
0089 auto x = x0 + dxdz * z;
0090 auto y = y0 + dydz * z;
0091
0092 z += z0;
0093 hx->Fill(x);
0094 hy->Fill(y);
0095 hz->Fill(z);
0096 auto ndof = trk_view[si].ndof();
0097 hchi2->Fill(vtx_view[si].chi2());
0098 hchi2oNdof->Fill(vtx_view[si].chi2() / ndof);
0099 hptv2->Fill(vtx_view[si].ptv2());
0100 hntrks->Fill(ndof + 1);
0101 }
0102 hnVertex->Fill(nVertices);
0103 }
0104
0105
0106
0107
0108 void SiPixelMonitorVertexSoAAlpaka::bookHistograms(DQMStore::IBooker& ibooker,
0109 edm::Run const& iRun,
0110 edm::EventSetup const& iSetup) {
0111
0112 ibooker.cd();
0113 ibooker.setCurrentFolder(topFolderName_);
0114 hnVertex = ibooker.book1D("nVertex", ";# of Vertices;#entries", 101, -0.5, 100.5);
0115 hx = ibooker.book1D("vx", ";Vertex x;#entries", 10, -5., 5.);
0116 hy = ibooker.book1D("vy", ";Vertex y;#entries", 10, -5., 5.);
0117 hz = ibooker.book1D("vz", ";Vertex z;#entries", 30, -30., 30);
0118 hchi2 = ibooker.book1D("chi2", ";Vertex chi-squared;#entries", 40, 0., 20.);
0119 hchi2oNdof = ibooker.book1D("chi2oNdof", ";Vertex chi-squared/Ndof;#entries", 40, 0., 20.);
0120 hptv2 = ibooker.book1D("ptsq", ";Vertex #sum (p_{T})^{2};#entries", 200, 0., 200.);
0121 hntrks = ibooker.book1D("ntrk", ";#tracks associated;#entries", 100, -0.5, 99.5);
0122 }
0123
0124 void SiPixelMonitorVertexSoAAlpaka::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0125
0126 edm::ParameterSetDescription desc;
0127 desc.add<edm::InputTag>("pixelVertexSrc", edm::InputTag("pixelVerticesAlpaka"));
0128 desc.add<edm::InputTag>("beamSpotSrc", edm::InputTag("offlineBeamSpot"));
0129 desc.add<std::string>("topFolderName", "SiPixelHeterogeneous/PixelVertexAlpaka");
0130 descriptions.addWithDefaultLabel(desc);
0131 }
0132
0133 DEFINE_FWK_MODULE(SiPixelMonitorVertexSoAAlpaka);