|
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
|
// TODO: change file name to SiPixelCompareVerticesSoA.cc when CUDA code is removed
// -*- C++ -*-
// Package: SiPixelCompareVertices
// Class: SiPixelCompareVertices
//
/**\class SiPixelCompareVertices SiPixelCompareVertices.cc
*/
//
// Author: Suvankar Roy Chowdhury
//
#include "DQMServices/Core/interface/DQMEDAnalyzer.h"
#include "DQMServices/Core/interface/DQMStore.h"
#include "DQMServices/Core/interface/MonitorElement.h"
#include "DataFormats/BeamSpot/interface/BeamSpot.h"
#include "DataFormats/Common/interface/Handle.h"
#include "DataFormats/VertexSoA/interface/ZVertexHost.h"
#include "FWCore/Framework/interface/Event.h"
#include "FWCore/Framework/interface/Frameworkfwd.h"
#include "FWCore/Framework/interface/MakerMacros.h"
#include "FWCore/MessageLogger/interface/MessageLogger.h"
#include "FWCore/ParameterSet/interface/ParameterSet.h"
#include "FWCore/Utilities/interface/InputTag.h"
// TODO: change class name to SiPixelCompareVerticesSoA when CUDA code is removed
class SiPixelCompareVertices : public DQMEDAnalyzer {
public:
using IndToEdm = std::vector<uint16_t>;
explicit SiPixelCompareVertices(const edm::ParameterSet&);
~SiPixelCompareVertices() override = default;
void bookHistograms(DQMStore::IBooker& ibooker, edm::Run const& iRun, edm::EventSetup const& iSetup) override;
void analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) override;
// analyzeSeparate is templated to accept distinct types of SoAs
// The default use case is to use vertices from Alpaka reconstructed on CPU and GPU;
template <typename U, typename V>
void analyzeSeparate(U tokenRef, V tokenTar, const edm::Event& iEvent);
static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
private:
// these two are both on Host but originally they have been produced on Host or on Device
const edm::EDGetTokenT<ZVertexHost> tokenSoAVertexReferenceSoA_;
const edm::EDGetTokenT<ZVertexHost> tokenSoAVertexTargetSoA_;
const edm::EDGetTokenT<reco::BeamSpot> tokenBeamSpot_;
const std::string topFolderName_;
const float dzCut_;
MonitorElement* hnVertex_;
MonitorElement* hx_;
MonitorElement* hy_;
MonitorElement* hz_;
MonitorElement* hchi2_;
MonitorElement* hchi2oNdof_;
MonitorElement* hptv2_;
MonitorElement* hntrks_;
MonitorElement* hxdiff_;
MonitorElement* hydiff_;
MonitorElement* hzdiff_;
};
//
// constructors
//
SiPixelCompareVertices::SiPixelCompareVertices(const edm::ParameterSet& iConfig)
: tokenSoAVertexReferenceSoA_(
consumes<ZVertexHost>(iConfig.getParameter<edm::InputTag>("pixelVertexReferenceSoA"))),
tokenSoAVertexTargetSoA_(consumes<ZVertexHost>(iConfig.getParameter<edm::InputTag>("pixelVertexTargetSoA"))),
tokenBeamSpot_(consumes<reco::BeamSpot>(iConfig.getParameter<edm::InputTag>("beamSpotSrc"))),
topFolderName_(iConfig.getParameter<std::string>("topFolderName")),
dzCut_(iConfig.getParameter<double>("dzCut")) {}
template <typename U, typename V>
void SiPixelCompareVertices::analyzeSeparate(U tokenRef, V tokenTar, const edm::Event& iEvent) {
const auto& vsoaHandleRef = iEvent.getHandle(tokenRef);
const auto& vsoaHandleTar = iEvent.getHandle(tokenTar);
if (not vsoaHandleRef or not vsoaHandleTar) {
edm::LogWarning out("SiPixelCompareVertices");
if (not vsoaHandleRef) {
out << "reference vertices not found; ";
}
if (not vsoaHandleTar) {
out << "target vertices not found; ";
}
out << "the comparison will not run.";
return;
}
auto const& vsoaRef = *vsoaHandleRef;
int nVerticesRef = vsoaRef.view().nvFinal();
auto const& vsoaTar = *vsoaHandleTar;
int nVerticesTar = vsoaTar.view().nvFinal();
auto bsHandle = iEvent.getHandle(tokenBeamSpot_);
float x0 = 0., y0 = 0., z0 = 0., dxdz = 0., dydz = 0.;
if (not bsHandle.isValid()) {
edm::LogWarning("SiPixelCompareVertices") << "No beamspot found, returning vertexes with (0,0,Z).";
} else {
const reco::BeamSpot& bs = *bsHandle;
x0 = bs.x0();
y0 = bs.y0();
z0 = bs.z0();
dxdz = bs.dxdz();
dydz = bs.dydz();
}
for (int ivc = 0; ivc < nVerticesRef; ivc++) {
auto sic = vsoaRef.view()[ivc].sortInd();
auto zc = vsoaRef.view()[sic].zv();
auto xc = x0 + dxdz * zc;
auto yc = y0 + dydz * zc;
zc += z0;
auto ndofRef = vsoaRef.template view<reco::ZVertexTracksSoA>()[sic].ndof();
auto chi2Ref = vsoaRef.view()[sic].chi2();
const int32_t notFound = -1;
int32_t closestVtxidx = notFound;
float mindz = dzCut_;
for (int ivg = 0; ivg < nVerticesTar; ivg++) {
auto sig = vsoaTar.view()[ivg].sortInd();
auto zgc = vsoaTar.view()[sig].zv() + z0;
auto zDist = std::abs(zc - zgc);
//insert some matching condition
if (zDist > dzCut_)
continue;
if (mindz > zDist) {
mindz = zDist;
closestVtxidx = sig;
}
}
if (closestVtxidx == notFound)
continue;
auto zg = vsoaTar.view()[closestVtxidx].zv();
auto xg = x0 + dxdz * zg;
auto yg = y0 + dydz * zg;
zg += z0;
auto ndofTar = vsoaTar.template view<reco::ZVertexTracksSoA>()[closestVtxidx].ndof();
auto chi2Tar = vsoaTar.view()[closestVtxidx].chi2();
hx_->Fill(xc - x0, xg - x0);
hy_->Fill(yc - y0, yg - y0);
hz_->Fill(zc, zg);
hxdiff_->Fill(xc - xg);
hydiff_->Fill(yc - yg);
hzdiff_->Fill(zc - zg);
hchi2_->Fill(chi2Ref, chi2Tar);
hchi2oNdof_->Fill(chi2Ref / ndofRef, chi2Tar / ndofTar);
hptv2_->Fill(vsoaRef.view()[sic].ptv2(), vsoaTar.view()[closestVtxidx].ptv2());
hntrks_->Fill(ndofRef + 1, ndofTar + 1);
}
hnVertex_->Fill(nVerticesRef, nVerticesTar);
}
//
// -- Analyze
//
void SiPixelCompareVertices::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) {
// The default use case is to use vertices from Alpaka reconstructed on CPU and GPU;
// The function is left templated if any other cases need to be added
analyzeSeparate(tokenSoAVertexReferenceSoA_, tokenSoAVertexTargetSoA_, iEvent);
}
//
// -- Book Histograms
//
void SiPixelCompareVertices::bookHistograms(DQMStore::IBooker& ibooker,
edm::Run const& iRun,
edm::EventSetup const& iSetup) {
ibooker.cd();
ibooker.setCurrentFolder(topFolderName_);
// 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
// these should be moved to a less resource consuming format
hnVertex_ = ibooker.book2I("nVertex", "# of Vertices;Reference;Target", 101, -0.5, 100.5, 101, -0.5, 100.5);
hx_ = ibooker.book2I("vx", "Vertez x - Beamspot x;Reference;Target", 50, -0.1, 0.1, 50, -0.1, 0.1);
hy_ = ibooker.book2I("vy", "Vertez y - Beamspot y;Reference;Target", 50, -0.1, 0.1, 50, -0.1, 0.1);
hz_ = ibooker.book2I("vz", "Vertez z;Reference;Target", 30, -30., 30., 30, -30., 30.);
hchi2_ = ibooker.book2I("chi2", "Vertex chi-squared;Reference;Target", 40, 0., 20., 40, 0., 20.);
hchi2oNdof_ = ibooker.book2I("chi2oNdof", "Vertex chi-squared/Ndof;Reference;Target", 40, 0., 20., 40, 0., 20.);
hptv2_ = ibooker.book2I("ptsq", "Vertex #sum (p_{T})^{2};Reference;Target", 200, 0., 200., 200, 0., 200.);
hntrks_ = ibooker.book2I("ntrk", "#tracks associated;Reference;Target", 100, -0.5, 99.5, 100, -0.5, 99.5);
hntrks_ = ibooker.book2I("ntrk", "#tracks associated;Reference;Target", 100, -0.5, 99.5, 100, -0.5, 99.5);
hxdiff_ = ibooker.book1D("vxdiff", ";Vertex x difference (Reference - Target);#entries", 100, -0.001, 0.001);
hydiff_ = ibooker.book1D("vydiff", ";Vertex y difference (Reference - Target);#entries", 100, -0.001, 0.001);
hzdiff_ = ibooker.book1D("vzdiff", ";Vertex z difference (Reference - Target);#entries", 100, -2.5, 2.5);
}
void SiPixelCompareVertices::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
// monitorpixelVertexSoA
edm::ParameterSetDescription desc;
desc.add<edm::InputTag>("pixelVertexReferenceSoA", edm::InputTag("pixelVerticesAlpakaSerial"));
desc.add<edm::InputTag>("pixelVertexTargetSoA", edm::InputTag("pixelVerticesAlpaka"));
desc.add<edm::InputTag>("beamSpotSrc", edm::InputTag("offlineBeamSpot"));
desc.add<std::string>("topFolderName", "SiPixelHeterogeneous/PixelVertexCompareSoADeviceVSHost");
desc.add<double>("dzCut", 1.);
descriptions.addWithDefaultLabel(desc);
}
// TODO: change module name to SiPixelCompareVerticesSoA when CUDA code is removed
DEFINE_FWK_MODULE(SiPixelCompareVertices);
|