CMS 3D CMS Logo

All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Modules Pages
SiPixelCompareRecHits.cc
Go to the documentation of this file.
1 // TODO: change file name to SiPixelCompareRecHitsSoA.cc when CUDA code is removed
2 
21 
22 // TODO: change class name to SiPixelCompareRecHitsSoA when CUDA code is removed
23 template <typename T>
25 public:
27 
29  ~SiPixelCompareRecHits() override = default;
30  void dqmBeginRun(const edm::Run&, const edm::EventSetup&) override;
31  void bookHistograms(DQMStore::IBooker& ibooker, edm::Run const& iRun, edm::EventSetup const& iSetup) override;
32  void analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) override;
33  // analyzeSeparate is templated to accept distinct types of SoAs
34  // The default use case is to use rechits from Alpaka reconstructed on CPU and GPU;
35  template <typename U, typename V>
36  void analyzeSeparate(U tokenRef, V tokenTar, const edm::Event& iEvent);
37  static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
38 
39 private:
42  // these two are both on Host but originally they have been produced on Host or on Device
46  const float mind2cut_;
48  static constexpr float micron_ = 10000.;
49  const TrackerGeometry* tkGeom_ = nullptr;
50  const TrackerTopology* tTopo_ = nullptr;
52  MonitorElement* hBchargeL_[4]; // max 4 barrel hits
57  MonitorElement* hFchargeD_[2][12]; // max 12 endcap disks
62  //differences
73 };
74 
75 //
76 // constructors
77 //
78 template <typename T>
81  topoToken_(esConsumes<TrackerTopology, TrackerTopologyRcd, edm::Transition::BeginRun>()),
82  tokenSoAHitsReference_(consumes(iConfig.getParameter<edm::InputTag>("pixelHitsReferenceSoA"))),
83  tokenSoAHitsTarget_(consumes(iConfig.getParameter<edm::InputTag>("pixelHitsTargetSoA"))),
84  topFolderName_(iConfig.getParameter<std::string>("topFolderName")),
85  mind2cut_(iConfig.getParameter<double>("minD2cut")) {}
86 
87 //
88 // Begin Run
89 //
90 template <typename T>
92  tkGeom_ = &iSetup.getData(geomToken_);
93  tTopo_ = &iSetup.getData(topoToken_);
94 }
95 
96 template <typename T>
97 template <typename U, typename V>
99  const auto& rhsoaHandleRef = iEvent.getHandle(tokenRef);
100  const auto& rhsoaHandleTar = iEvent.getHandle(tokenTar);
101 
102  // Exit early if any handle is invalid
103  if (!rhsoaHandleRef || !rhsoaHandleTar) {
104  edm::LogWarning out("SiPixelCompareRecHits");
105  if (!rhsoaHandleRef)
106  out << "reference rechits not found; ";
107  if (!rhsoaHandleTar)
108  out << "target rechits not found; ";
109  out << "the comparison will not run.";
110  return;
111  }
112 
113  const auto& rhsoaRef = *rhsoaHandleRef;
114  const auto& rhsoaTar = *rhsoaHandleTar;
115 
116  auto const& soa2dRef = rhsoaRef.const_view();
117  auto const& soa2dTar = rhsoaTar.const_view();
118 
119  uint32_t nHitsRef = soa2dRef.metadata().size();
120  uint32_t nHitsTar = soa2dTar.metadata().size();
121 
122  hnHits_->Fill(nHitsRef, nHitsTar);
123 
124  // Map detector indices to target hits for quick access
125  std::unordered_map<uint16_t, std::vector<size_t>> detectorIndexMap;
126  detectorIndexMap.reserve(nHitsTar);
127  for (size_t j = 0; j < nHitsTar; ++j) {
128  detectorIndexMap[soa2dTar[j].detectorIndex()].push_back(j);
129  }
130 
131  auto detIds = tkGeom_->detUnitIds();
132 
133  // Loop through reference hits
134  for (uint32_t i = 0; i < nHitsRef; i++) {
135  float minD = mind2cut_;
136  uint32_t matchedHit = invalidHit_;
137  uint16_t indRef = soa2dRef[i].detectorIndex();
138  float xLocalRef = soa2dRef[i].xLocal();
139  float yLocalRef = soa2dRef[i].yLocal();
140 
141  // Look up hits in target with matching detector index
142  auto it = detectorIndexMap.find(indRef);
143  if (it != detectorIndexMap.end()) {
144  for (auto j : it->second) {
145  float dx = xLocalRef - soa2dTar[j].xLocal();
146  float dy = yLocalRef - soa2dTar[j].yLocal();
147  float distance = dx * dx + dy * dy;
148  if (distance < minD) {
149  minD = distance;
150  matchedHit = j;
151  }
152  }
153  }
154 
155  // Gather reference hit properties
156  DetId id = detIds[indRef];
157  uint32_t chargeRef = soa2dRef[i].chargeAndStatus().charge;
158  int16_t sizeXRef = (soa2dRef[i].clusterSizeX() + 7) / 8;
159  int16_t sizeYRef = (soa2dRef[i].clusterSizeY() + 7) / 8;
160 
161  // Initialize target hit properties
162  uint32_t chargeTar = 0;
163  int16_t sizeXTar = -99;
164  int16_t sizeYTar = -99;
165  float xLocalTar = -999.;
166  float yLocalTar = -999.;
167 
168  if (matchedHit != invalidHit_) {
169  chargeTar = soa2dTar[matchedHit].chargeAndStatus().charge;
170  sizeXTar = (soa2dTar[matchedHit].clusterSizeX() + 7) / 8;
171  sizeYTar = (soa2dTar[matchedHit].clusterSizeY() + 7) / 8;
172  xLocalTar = soa2dTar[matchedHit].xLocal();
173  yLocalTar = soa2dTar[matchedHit].yLocal();
174  }
175 
176  // Populate histograms based on subdetector type
177  switch (id.subdetId()) {
179  hBchargeL_[tTopo_->pxbLayer(id) - 1]->Fill(chargeRef, chargeTar);
180  hBsizexL_[tTopo_->pxbLayer(id) - 1]->Fill(sizeXRef, sizeXTar);
181  hBsizeyL_[tTopo_->pxbLayer(id) - 1]->Fill(sizeYRef, sizeYTar);
182  hBposxL_[tTopo_->pxbLayer(id) - 1]->Fill(xLocalRef, xLocalTar);
183  hBposyL_[tTopo_->pxbLayer(id) - 1]->Fill(yLocalRef, yLocalTar);
184  hBchargeDiff_->Fill(chargeRef - chargeTar);
185  hBsizeXDiff_->Fill(sizeXRef - sizeXTar);
186  hBsizeYDiff_->Fill(sizeYRef - sizeYTar);
187  hBposXDiff_->Fill(micron_ * (xLocalRef - xLocalTar));
188  hBposYDiff_->Fill(micron_ * (yLocalRef - yLocalTar));
189  break;
191  hFchargeD_[tTopo_->pxfSide(id) - 1][tTopo_->pxfDisk(id) - 1]->Fill(chargeRef, chargeTar);
192  hFsizexD_[tTopo_->pxfSide(id) - 1][tTopo_->pxfDisk(id) - 1]->Fill(sizeXRef, sizeXTar);
193  hFsizeyD_[tTopo_->pxfSide(id) - 1][tTopo_->pxfDisk(id) - 1]->Fill(sizeYRef, sizeYTar);
194  hFposxD_[tTopo_->pxfSide(id) - 1][tTopo_->pxfDisk(id) - 1]->Fill(xLocalRef, xLocalTar);
195  hFposyD_[tTopo_->pxfSide(id) - 1][tTopo_->pxfDisk(id) - 1]->Fill(yLocalRef, yLocalTar);
196  hFchargeDiff_->Fill(chargeRef - chargeTar);
197  hFsizeXDiff_->Fill(sizeXRef - sizeXTar);
198  hFsizeYDiff_->Fill(sizeYRef - sizeYTar);
199  hFposXDiff_->Fill(micron_ * (xLocalRef - xLocalTar));
200  hFposYDiff_->Fill(micron_ * (yLocalRef - yLocalTar));
201  break;
202  }
203  }
204 }
205 
206 //
207 // -- Analyze
208 //
209 template <typename T>
211  // The default use case is to use vertices from Alpaka reconstructed on CPU and GPU;
212  // The function is left templated if any other cases need to be added
213  analyzeSeparate(tokenSoAHitsReference_, tokenSoAHitsTarget_, iEvent);
214 }
215 
216 //
217 // -- Book Histograms
218 //
219 template <typename T>
221  edm::Run const& iRun,
222  edm::EventSetup const& iSetup) {
223  iBook.cd();
224  iBook.setCurrentFolder(topFolderName_);
225 
226  // clang-format off
227  //Global
228  hnHits_ = iBook.book2I("nHits", "ReferencevsTarget RecHits per event;#Reference RecHits;#Target RecHits", 200, 0, 5000,200, 0, 5000);
229  //Barrel Layer
230  for(unsigned int il=0;il<tkGeom_->numberOfLayers(PixelSubdetector::PixelBarrel);il++){
231  hBchargeL_[il] = iBook.book2I(Form("recHitsBLay%dCharge",il+1), Form("ReferencevsTarget RecHits Charge Barrel Layer%d;Reference Charge;Target Charge",il+1), 250, 0, 100000, 250, 0, 100000);
232  hBsizexL_[il] = iBook.book2I(Form("recHitsBLay%dSizex",il+1), Form("ReferencevsTarget RecHits SizeX Barrel Layer%d;Reference SizeX;Target SizeX",il+1), 30, 0, 30, 30, 0, 30);
233  hBsizeyL_[il] = iBook.book2I(Form("recHitsBLay%dSizey",il+1), Form("ReferencevsTarget RecHits SizeY Barrel Layer%d;Reference SizeY;Target SizeY",il+1), 30, 0, 30, 30, 0, 30);
234  hBposxL_[il] = iBook.book2D(Form("recHitsBLay%dPosx",il+1), Form("ReferencevsTarget RecHits x-pos in Barrel Layer%d;Reference pos x;Target pos x",il+1), 200, -5, 5, 200,-5,5);
235  hBposyL_[il] = iBook.book2D(Form("recHitsBLay%dPosy",il+1), Form("ReferencevsTarget RecHits y-pos in Barrel Layer%d;Reference pos y;Target pos y",il+1), 200, -5, 5, 200,-5,5);
236  }
237  //Endcaps
238  //Endcaps Disk
239  for(int is=0;is<2;is++){
240  int sign=is==0? -1:1;
241  for(unsigned int id=0;id<tkGeom_->numberOfLayers(PixelSubdetector::PixelEndcap);id++){
242  hFchargeD_[is][id] = iBook.book2I(Form("recHitsFDisk%+dCharge",id*sign+sign), Form("ReferencevsTarget RecHits Charge Endcaps Disk%+d;Reference Charge;Target Charge",id*sign+sign), 250, 0, 100000, 250, 0, 100000);
243  hFsizexD_[is][id] = iBook.book2I(Form("recHitsFDisk%+dSizex",id*sign+sign), Form("ReferencevsTarget RecHits SizeX Endcaps Disk%+d;Reference SizeX;Target SizeX",id*sign+sign), 30, 0, 30, 30, 0, 30);
244  hFsizeyD_[is][id] = iBook.book2I(Form("recHitsFDisk%+dSizey",id*sign+sign), Form("ReferencevsTarget RecHits SizeY Endcaps Disk%+d;Reference SizeY;Target SizeY",id*sign+sign), 30, 0, 30, 30, 0, 30);
245  hFposxD_[is][id] = iBook.book2D(Form("recHitsFDisk%+dPosx",id*sign+sign), Form("ReferencevsTarget RecHits x-pos Endcaps Disk%+d;Reference pos x;Target pos x",id*sign+sign), 200, -5, 5, 200, -5, 5);
246  hFposyD_[is][id] = iBook.book2D(Form("recHitsFDisk%+dPosy",id*sign+sign), Form("ReferencevsTarget RecHits y-pos Endcaps Disk%+d;Reference pos y;Target pos y",id*sign+sign), 200, -5, 5, 200, -5, 5);
247  }
248  }
249  //1D differences
250  hBchargeDiff_ = iBook.book1D("rechitChargeDiffBpix","Charge differnce of rechits in BPix; rechit charge difference (Reference - Target)", 101, -50.5, 50.5);
251  hFchargeDiff_ = iBook.book1D("rechitChargeDiffFpix","Charge differnce of rechits in FPix; rechit charge difference (Reference - Target)", 101, -50.5, 50.5);
252  hBsizeXDiff_ = iBook.book1D("rechitsizeXDiffBpix","SizeX difference of rechits in BPix; rechit sizex difference (Reference - Target)", 21, -10.5, 10.5);
253  hFsizeXDiff_ = iBook.book1D("rechitsizeXDiffFpix","SizeX difference of rechits in FPix; rechit sizex difference (Reference - Target)", 21, -10.5, 10.5);
254  hBsizeYDiff_ = iBook.book1D("rechitsizeYDiffBpix","SizeY difference of rechits in BPix; rechit sizey difference (Reference - Target)", 21, -10.5, 10.5);
255  hFsizeYDiff_ = iBook.book1D("rechitsizeYDiffFpix","SizeY difference of rechits in FPix; rechit sizey difference (Reference - Target)", 21, -10.5, 10.5);
256  hBposXDiff_ = iBook.book1D("rechitsposXDiffBpix","x-position difference of rechits in BPix; rechit x-pos difference (Reference - Target)", 1000, -10, 10);
257  hFposXDiff_ = iBook.book1D("rechitsposXDiffFpix","x-position difference of rechits in FPix; rechit x-pos difference (Reference - Target)", 1000, -10, 10);
258  hBposYDiff_ = iBook.book1D("rechitsposYDiffBpix","y-position difference of rechits in BPix; rechit y-pos difference (Reference - Target)", 1000, -10, 10);
259  hFposYDiff_ = iBook.book1D("rechitsposYDiffFpix","y-position difference of rechits in FPix; rechit y-pos difference (Reference - Target)", 1000, -10, 10);
260 }
261 
262 template<typename T>
264  // monitorpixelRecHitsSoAAlpaka
266  desc.add<edm::InputTag>("pixelHitsReferenceSoA", edm::InputTag("siPixelRecHitsPreSplittingAlpakaSerial"));
267  desc.add<edm::InputTag>("pixelHitsTargetSoA", edm::InputTag("siPixelRecHitsPreSplittingAlpaka"));
268  desc.add<std::string>("topFolderName", "SiPixelHeterogeneous/PixelRecHitsCompareDeviceVSHost");
269  desc.add<double>("minD2cut", 0.0001);
270  descriptions.addWithDefaultLabel(desc);
271 }
272 
276 
278 // TODO: change module names to SiPixel*CompareRecHitsSoA when CUDA code is removed
282 
const edm::ESGetToken< TrackerTopology, TrackerTopologyRcd > topoToken_
void addWithDefaultLabel(ParameterSetDescription const &psetDescription)
ESGetTokenH3DDVariant esConsumes(std::string const &Record, edm::ConsumesCollector &)
Definition: DeDxTools.cc:283
MonitorElement * hFsizexD_[2][12]
T const & getData(const ESGetToken< T, R > &iToken) const noexcept(false)
Definition: EventSetup.h:119
virtual void setCurrentFolder(std::string const &fullpath)
Definition: DQMStore.cc:36
static void fillDescriptions(edm::ConfigurationDescriptions &descriptions)
MonitorElement * hBchargeDiff_
MonitorElement * hBchargeL_[4]
MonitorElement * hBposxL_[4]
const std::string topFolderName_
const edm::EDGetTokenT< HitsSoA > tokenSoAHitsReference_
MonitorElement * hBsizeyL_[4]
const edm::ESGetToken< TrackerGeometry, TrackerDigiGeometryRecord > geomToken_
void analyzeSeparate(U tokenRef, V tokenTar, const edm::Event &iEvent)
void bookHistograms(DQMStore::IBooker &ibooker, edm::Run const &iRun, edm::EventSetup const &iSetup) override
const edm::EDGetTokenT< HitsSoA > tokenSoAHitsTarget_
MonitorElement * hFchargeD_[2][12]
~SiPixelCompareRecHits() override=default
int iEvent
Definition: GenABIO.cc:224
uint32_t T const *__restrict__ uint32_t const *__restrict__ int32_t int Histo::index_type cudaStream_t V
MonitorElement * hFchargeDiff_
MonitorElement * hFsizeyD_[2][12]
void analyze(const edm::Event &iEvent, const edm::EventSetup &iSetup) override
const TrackerTopology * tTopo_
Transition
Definition: Transition.h:12
#define DEFINE_FWK_MODULE(type)
Definition: MakerMacros.h:16
MonitorElement * hBposyL_[4]
MonitorElement * hFposxD_[2][12]
Definition: DetId.h:17
MonitorElement * hBsizexL_[4]
SiPixelCompareRecHits(const edm::ParameterSet &)
static constexpr uint32_t invalidHit_
MonitorElement * book2D(TString const &name, TString const &title, int nchX, double lowX, double highX, int nchY, double lowY, double highY, FUNC onbooking=NOOP())
Definition: DQMStore.h:221
HLT enums.
const TrackerGeometry * tkGeom_
MonitorElement * book1D(TString const &name, TString const &title, int const nchX, double const lowX, double const highX, FUNC onbooking=NOOP())
Definition: DQMStore.h:98
MonitorElement * book2I(TString const &name, TString const &title, int nchX, double lowX, double highX, int nchY, double lowY, double highY, FUNC onbooking=NOOP())
Definition: DQMStore.h:305
static constexpr float micron_
void dqmBeginRun(const edm::Run &, const edm::EventSetup &) override
Definition: Run.h:45
MonitorElement * hFposyD_[2][12]