Belle II Software development
PXDClusterizerModule.cc
1/**************************************************************************
2 * basf2 (Belle II Analysis Software Framework) *
3 * Author: The Belle II Collaboration *
4 * *
5 * See git log for contributors and copyright holders. *
6 * This file is licensed under LGPL-3.0, see LICENSE.md. *
7 **************************************************************************/
8
9#include <pxd/modules/pxdReconstruction/PXDClusterizerModule.h>
10#include <pxd/dbobjects/PXDClusterPositionErrorPar.h>
11#include <pxd/reconstruction/ClusterCache.h>
12#include <pxd/reconstruction/ClusterProjection.h>
13
14#include <framework/datastore/StoreArray.h>
15#include <framework/datastore/RelationArray.h>
16#include <framework/logging/Logger.h>
17
18#include <vxd/dataobjects/VxdID.h>
19#include <vxd/geometry/GeoCache.h>
20
21#include <mdst/dataobjects/MCParticle.h>
22#include <mdst/dataobjects/EventLevelTrackingInfo.h>
23#include <pxd/dataobjects/PXDDigit.h>
24#include <pxd/dataobjects/PXDCluster.h>
25#include <pxd/dataobjects/PXDTrueHit.h>
26#include <pxd/geometry/SensorInfo.h>
27
28#include <pxd/reconstruction/PXDClusterPositionEstimator.h>
29
30#include <pxd/utilities/PXDUtilities.h>
31
32using namespace std;
33using namespace Belle2;
34using namespace Belle2::PXD;
35
36//-----------------------------------------------------------------
37// Register the Module
38//-----------------------------------------------------------------
39REG_MODULE(PXDClusterizer);
40
41//-----------------------------------------------------------------
42// Implementation
43//-----------------------------------------------------------------
44
46 , m_elNoise(0.7), m_cutSeed(5.0), m_cutAdjacent(3.0), m_cutCluster(8.0)
48{
49 //Set module properties
50 setDescription("Cluster PXDHits");
52
53 addParam("ElectronicNoise", m_elNoise,
54 "Noise added by the electronics, set in ADU", m_elNoise);
55 addParam("NoiseSN", m_cutAdjacent,
56 "SN for digits to be considered for clustering", m_cutAdjacent);
57 addParam("SeedSN", m_cutSeed, "SN for digits to be considered as seed",
58 m_cutSeed);
59 addParam("ClusterSN", m_cutCluster, "Minimum SN for clusters", m_cutCluster);
60 addParam("ClusterCacheSize", m_clusterCacheSize,
61 "Maximum desired number of sensor rows", 0);
62 addParam("HeadTailSize", m_sizeHeadTail,
63 "Minimum cluster size to switch to Analog head tail algorithm for cluster center",
65 addParam("Digits", m_storeDigitsName, "Digits collection name", string(""));
66 addParam("Clusters", m_storeClustersName, "Cluster collection name",
67 string(""));
68 addParam("TrueHits", m_storeTrueHitsName, "TrueHit collection name",
69 string(""));
70 addParam("MCParticles", m_storeMCParticlesName, "MCParticles collection name",
71 string(""));
72 addParam("ErrorFromDB", m_errorFromDB, "Assign cluster position error from DB", true);
73 addParam("PositionErrorUPayloadName", m_positionErrorUName, "Payload name for cluster position error in U",
74 string("PXDClusterPositionErrorUPar"));
75 addParam("PositionErrorVPayloadName", m_positionErrorVName, "Payload name for cluster position error in V",
76 string("PXDClusterPositionErrorVPar"));
77 addParam("createPXDClustersForAbortedTrackingEvents", m_createPXDClustersForAbortedTrackingEvents,
78 "Create PXDClusters for events where either the SVDSpacePointCreator abort flag or the VXDTF2 and SVDCKF abort flags are set.",
80
81}
82
84{
85 //Register collections
90
92 storeDigits.isRequired();
93 storeTrueHits.isOptional();
94 storeMCParticles.isOptional();
95
96 RelationArray relClusterDigits(storeClusters, storeDigits);
97 RelationArray relClusterTrueHits(storeClusters, storeTrueHits);
98 RelationArray relClusterMCParticles(storeClusters, storeMCParticles);
99 RelationArray relDigitMCParticles(storeDigits, storeMCParticles);
100 RelationArray relDigitTrueHits(storeDigits, storeTrueHits);
101
102 relClusterDigits.registerInDataStore();
103 //Relations to Simulation objects are only needed if the parent one already
104 //exists
105 if (relDigitTrueHits.isOptional()) {
106 relClusterTrueHits.registerInDataStore();
107 }
108 if (relDigitMCParticles.isOptional()) {
109 relClusterMCParticles.registerInDataStore();
110 }
111
112 //Now save all names to speed creation later, mainly if they had default names
113 m_storeClustersName = storeClusters.getName();
114 m_storeDigitsName = storeDigits.getName();
115 m_storeTrueHitsName = storeTrueHits.getName();
116 m_storeMCParticlesName = storeMCParticles.getName();
117
118 m_relClusterDigitName = relClusterDigits.getName();
119 m_relClusterTrueHitName = relClusterTrueHits.getName();
120 m_relClusterMCParticleName = relClusterMCParticles.getName();
121 m_relDigitTrueHitName = relDigitTrueHits.getName();
122 m_relDigitMCParticleName = relDigitMCParticles.getName();
123
124 B2DEBUG(20,
125 "PXDClusterizer Parameters (in default system units, *=cannot be set directly):");
126 B2DEBUG(20, " --> ElectronicNoise: " << m_elNoise);
127 B2DEBUG(20, " --> NoiseSN: " << m_cutAdjacent);
128 B2DEBUG(20, " --> SeedSN: " << m_cutSeed);
129 B2DEBUG(20, " --> ClusterSN: " << m_cutCluster);
130 B2DEBUG(20,
132 B2DEBUG(20,
134 B2DEBUG(20,
136 B2DEBUG(20,
138 B2DEBUG(20, " --> DigitMCRel: " << m_relDigitMCParticleName);
139 B2DEBUG(20, " --> ClusterMCRel: " << m_relClusterMCParticleName);
140 B2DEBUG(20, " --> ClusterDigitRel: " << m_relClusterDigitName);
141 B2DEBUG(20, " --> DigitTrueRel: " << m_relDigitTrueHitName);
142 B2DEBUG(20, " --> ClusterTrueRel: " << m_relClusterTrueHitName);
143
144
145 m_noiseMap.setNoiseLevel(m_elNoise);
147 if (m_clusterCacheSize > 0)
148 m_cache = std::unique_ptr<ClusterCache>(new ClusterCache(m_clusterCacheSize));
149 else
150 m_cache = std::unique_ptr<ClusterCache>(new ClusterCache());
151
152 // Cluster position error from DB
153 if (m_errorFromDB) {
154 if (!m_positionErrorUName.size() || !m_positionErrorVName.size()) {
155 B2WARNING("You chose to use cluster position errors from DB but did not provide PositionErrorUPayloadName ("
156 << m_positionErrorUName << ") and/or PositionErrorVPayloadName (" << m_positionErrorVName
157 << "). Disabling DB option.");
158 m_errorFromDB = false;
159 }
160 m_clusterPositionErrorUPar = std::make_unique<DBObjPtr<PXDClusterPositionErrorPar>>(m_positionErrorUName);
161 m_clusterPositionErrorVPar = std::make_unique<DBObjPtr<PXDClusterPositionErrorPar>>(m_positionErrorVName);
162 if (m_clusterPositionErrorUPar == nullptr || m_clusterPositionErrorVPar == nullptr) {
163 B2FATAL("DB objects for ClusterPositionError not valid");
164 }
165 }
166
167 m_eventLevelTrackingInfo.isOptional();
168
169}
170
172{
173 // Need to check if payload should be re-loaded at run change (currently just one revision for the entire period)
174 if (m_errorFromDB && (!m_clusterPositionErrorUPar->isValid() || !m_clusterPositionErrorVPar.get()->isValid())) {
175 B2FATAL("DB objects for ClusterPositionError not valid for this run");
176 }
177}
178
180{
181 // Abort in case SVDSpacePointCreator was aborted (high occupancy events) or if both VXDTF2 and SVDCKF were aborted
182 // as we do not add PXD hits to CDC standalone tracks, so we don't need to run the PXDClusterizer.
183 // This veto can be overwritten by setting m_createPXDClustersForAbortedTrackingEvents to true to create them regardless.
184 if (m_eventLevelTrackingInfo.isValid()) {
186 (m_eventLevelTrackingInfo->hasSVDSpacePointCreatorAbortionFlag() or
187 (m_eventLevelTrackingInfo->hasSVDCKFAbortionFlag() and m_eventLevelTrackingInfo->hasVXDTF2AbortionFlag()))) {
188 return;
189 }
190 }
191
192 const StoreArray<MCParticle> storeMCParticles(m_storeMCParticlesName);
193 const StoreArray<PXDTrueHit> storeTrueHits(m_storeTrueHitsName);
194 const StoreArray<PXDDigit> storeDigits(m_storeDigitsName);
196
197 storeClusters.clear();
198
199 RelationArray relClusterMCParticle(storeClusters, storeMCParticles,
201 if (relClusterMCParticle) relClusterMCParticle.clear();
202
203 RelationArray relClusterDigit(storeClusters, storeDigits,
205 if (relClusterDigit) relClusterDigit.clear();
206
207 RelationArray relClusterTrueHit(storeClusters, storeTrueHits,
209 if (relClusterTrueHit) relClusterTrueHit.clear();
210
211 int nPixels = storeDigits.getEntries();
212 if (nPixels == 0)
213 return;
214
215 //Build lookup tables for relations
216 RelationArray relDigitMCParticle(storeDigits, storeMCParticles, m_relDigitMCParticleName);
217 RelationArray relDigitTrueHit(storeDigits, storeTrueHits, m_relDigitTrueHitName);
218 createRelationLookup(relDigitMCParticle, m_mcRelation, storeDigits.getEntries());
219 createRelationLookup(relDigitTrueHit, m_trueRelation, storeDigits.getEntries());
220
221 m_cache->clear();
222
223 //We require all pixels are already sorted and directly cluster them. Once
224 //the sensorID changes, we write out all existing clusters and continue.
225 VxdID sensorID(0);
226 //To check sorting
227 Pixel lastPixel;
228 for (int i = 0; i < nPixels; i++) {
229 const PXDDigit* const storeDigit = storeDigits[i];
230 Pixel px(storeDigit, i);
231 //New sensor, write clusters
232 if (sensorID != storeDigit->getSensorID()) {
233 writeClusters(sensorID);
234 sensorID = storeDigit->getSensorID();
235 //Load the correct noise map for the new sensor
236 m_noiseMap.setSensorID(sensorID);
237 } else if (px <= lastPixel) {
238 //Check for sorting as precaution
239 B2FATAL("Pixels are not sorted correctly, please include the "
240 "PXDDigitSorter module before running the Clusterizer or fix "
241 "the input to be ordered by v,u in ascending order");
242 }
243 //Remember last pixel to check sorting
244 lastPixel = px;
245
246 //Ignore digits with not enough signal. Has to be done after the check of
247 //the SensorID to make sure we compare to the correct noise level
248 if (!m_noiseMap(px, m_cutAdjacent)) continue;
249
250 // TODO: If we would like to cluster across dead channels we would need a
251 // sorted list of dead pixels. Assuming we have such a list m_deadChannels
252 // containing Pixel instances with all dead channels of this sensor and
253 // m_currentDeadChannel is an iterator to that list (initialized to
254 // begin(m_deadChannels) when the sensorID changes), we would do
255 //
256 // while(m_currentDeadChannel != end(m_deadChannels) && *m_currentDeadChannel < px){
257 // m_cache->findCluster(m_currentDeadChannel->getU(), m_currentDeadChannel->getV());
258 // m_currentDeadChannel++;
259 // }
260 //
261 // This would have the effect of marking the pixel address as belonging
262 // to a cluster but would not modify the clusters itself (except for
263 // possible merging of clusters) so we would end up with clusters that
264 // contain holes or are disconnected.
265
266 // Find correct cluster and add pixel to cluster
267 try {
268 m_cache->findCluster(px.getU(), px.getV()).add(px);
269 } catch (std::out_of_range& e) {
270 B2WARNING("PXD clustering: Ignoring pixel " << px.getU() << "," << px.getV() << ": " << e.what());
271 }
272 }
273 writeClusters(sensorID);
274}
275
277{
278 lookup.clear();
279 //If we don't have a relation we don't build a lookuptable
280 if (!relation) return;
281 //Resize to number of digits and set all values
282 lookup.resize(digits);
283 for (const RelationElement& element : relation) {
284 lookup[element.getFromIndex()] = &element;
285 }
286}
287
288void PXDClusterizerModule::fillRelationMap(const RelationLookup& lookup, std::map<unsigned int, float>& relation,
289 unsigned int index)
290{
291 //If the lookup table is not empty and the element is set
292 if (!lookup.empty() && lookup[index]) {
293 const RelationElement& element = *lookup[index];
294 const unsigned int size = element.getSize();
295 //Add all Relations to the map
296 for (unsigned int i = 0; i < size; ++i) {
297 //negative weights are from ignored particles, we don't like them and
298 //thus ignore them :D
299 if (element.getWeight(i) < 0) continue;
300 relation[element.getToIndex(i)] += element.getWeight(i);
301 }
302 }
303}
304
306{
307 if (m_cache->empty())
308 return;
309
310 //Get all datastore elements
311 const StoreArray<MCParticle> storeMCParticles(m_storeMCParticlesName);
312 const StoreArray<PXDDigit> storeDigits(m_storeDigitsName);
313 const StoreArray<PXDTrueHit> storeTrueHits(m_storeTrueHitsName);
315 RelationArray relClusterMCParticle(storeClusters, storeMCParticles,
317 RelationArray relClusterDigit(storeClusters, storeDigits,
319 RelationArray relClusterTrueHit(storeClusters, storeTrueHits,
321
322 //Get Geometry information
323 const SensorInfo& info = dynamic_cast<const SensorInfo&>(VXD::GeoCache::getInstance().getSensorInfo(
324 sensorID));
325
326 map<unsigned int, float> mc_relations;
327 map<unsigned int, float> truehit_relations;
328 vector<pair<unsigned int, float> > digit_weights;
329
330
331 for (ClusterCandidate& cls : *m_cache) {
332 //Check for noise cuts
333 if (!(cls.size() > 0 && m_noiseMap(cls.getCharge(), m_cutCluster) && m_noiseMap(cls.getSeed(), m_cutSeed))) continue;
334
335 double rho(0);
336 ClusterProjection projU, projV;
337 mc_relations.clear();
338 truehit_relations.clear();
339 digit_weights.clear();
340 digit_weights.reserve(cls.size());
341
342 const Pixel& seed = cls.getSeed();
343
344 for (const PXD::Pixel& px : cls.pixels()) {
345 //Add the Pixel information to the two projections
346 projU.add(px.getU(), info.getUCellPosition(px.getU()), px.getCharge());
347 projV.add(px.getV(), info.getVCellPosition(px.getV()), px.getCharge());
348
349 //Obtain relations from MCParticles
350 fillRelationMap(m_mcRelation, mc_relations, px.getIndex());
351 //Obtain relations from PXDTrueHits
352 fillRelationMap(m_trueRelation, truehit_relations, px.getIndex());
353 //Save the weight of the digits for the Cluster->Digit relation
354 digit_weights.emplace_back(px.getIndex(), px.getCharge());
355 }
356 projU.finalize();
357 projV.finalize();
358
359 const double pitchU = info.getUPitch();
360 const double pitchV = info.getVPitch(projV.getPos());
361 // Calculate shape correlation coefficient: only for non-trivial shapes
362 if (projU.getSize() > 1 && projV.getSize() > 1) {
363 // Add in-pixel position noise to smear the correlation
364 double posUU = cls.getCharge() * pitchU * pitchU / 12.0;
365 double posVV = cls.getCharge() * pitchV * pitchV / 12.0;
366 double posUV(0);
367 for (const Pixel& px : cls.pixels()) {
368 const double du = info.getUCellPosition(px.getU()) - projU.getPos();
369 const double dv = info.getVCellPosition(px.getV()) - projV.getPos();
370 posUU += px.getCharge() * du * du;
371 posVV += px.getCharge() * dv * dv;
372 posUV += px.getCharge() * du * dv;
373 }
374 rho = posUV / sqrt(posUU * posVV);
375 }
376
377 //Calculate position and error with u as primary axis, fixed pitch size
378 calculatePositionError(cls, projU, projV, pitchU, pitchU, pitchU);
379 //Calculate position and error with v as primary axis, possibly different pitch sizes
380 calculatePositionError(cls, projV, projU, info.getVPitch(projV.getMinPos()), pitchV, info.getVPitch(projV.getMaxPos()));
381
382 if (m_errorFromDB) { // Overwrite cluster position error with value from DB (keep the above calculation untouched for now)
383 unsigned int uID = info.getUCellID(projU.getPos());
384 unsigned int vID = info.getVCellID(projV.getPos());
385 bool isUedge = PXD::isClusterAtUEdge(sensorID, projU.getMinCell(), projU.getMaxCell());
386 bool isVedge = (PXD::isClusterAtVEdge(sensorID, projV.getMinCell(), projV.getMaxCell())
387 || PXD::isClusterAtLadderJoint(sensorID, projV.getMinCell(), projV.getMaxCell()));
388 assignPositionErrorFromDB(projU, **m_clusterPositionErrorUPar, sensorID, uID, vID, pitchU, isUedge, isVedge);
389 assignPositionErrorFromDB(projV, **m_clusterPositionErrorVPar, sensorID, uID, vID, pitchV, isUedge, isVedge);
390 }
391
392 ROOT::Math::XYZVector lorentzShift = info.getLorentzShift(projU.getPos(), projV.getPos());
393 projU.setPos(projU.getPos() - lorentzShift.X());
394 projV.setPos(projV.getPos() - lorentzShift.Y());
395 B2DEBUG(20, "Lorentz shift: " << lorentzShift.X() << " " << lorentzShift.Y());
396
397 // Pre classification of cluster looking at pitch type of pixels and if they touch sensor edges
398 int clusterkind = PXDClusterPositionEstimator::getInstance().getClusterkind(cls.pixels(), sensorID);
399
400 // Compute sorted set of pixel
401 // FIXME: I am not 100% sure if cls.pixels() are sorted
402 set<Pixel> pixelSet(cls.pixels().begin(), cls.pixels().end());
403
404 // Compute classifier variables needed for later retrieval of position correction in PXD CKF
405 vector<float> sectorEtaValues = {0, 0, 0, 0};
406 sectorEtaValues[0] = PXDClusterPositionEstimator::getInstance().computeEta(pixelSet, projV.getMinCell(), projV.getSize(), +1.0,
407 +1.0);
408 sectorEtaValues[1] = PXDClusterPositionEstimator::getInstance().computeEta(pixelSet, projV.getMinCell(), projV.getSize(), -1.0,
409 +1.0);
410 sectorEtaValues[2] = PXDClusterPositionEstimator::getInstance().computeEta(pixelSet, projV.getMinCell(), projV.getSize(), -1.0,
411 -1.0);
412 sectorEtaValues[3] = PXDClusterPositionEstimator::getInstance().computeEta(pixelSet, projV.getMinCell(), projV.getSize(), +1.0,
413 -1.0);
414
415 vector<int> sectorShapeIndices = { -1, -1, -1, -1};
416 sectorShapeIndices[0] = PXDClusterPositionEstimator::getInstance().computeShapeIndex(pixelSet, projU.getMinCell(),
417 projV.getMinCell(), projV.getSize(), +1.0, +1.0);
418 sectorShapeIndices[1] = PXDClusterPositionEstimator::getInstance().computeShapeIndex(pixelSet, projU.getMinCell(),
419 projV.getMinCell(), projV.getSize(), -1.0, +1.0);
420 sectorShapeIndices[2] = PXDClusterPositionEstimator::getInstance().computeShapeIndex(pixelSet, projU.getMinCell(),
421 projV.getMinCell(), projV.getSize(), -1.0, -1.0);
422 sectorShapeIndices[3] = PXDClusterPositionEstimator::getInstance().computeShapeIndex(pixelSet, projU.getMinCell(),
423 projV.getMinCell(), projV.getSize(), +1.0, -1.0);
424
425 //Store Cluster into Datastore ...
426 int clsIndex = storeClusters.getEntries();
427 storeClusters.appendNew(sensorID, projU.getPos(), projV.getPos(), projU.getError(), projV.getError(),
428 rho, cls.getCharge(), seed.getCharge(),
429 cls.size(), projU.getSize(), projV.getSize(), projU.getMinCell(), projV.getMinCell(), clusterkind,
430 sectorEtaValues, sectorShapeIndices
431 );
432
433 //Create Relations to this Digit
434 if (!mc_relations.empty()) relClusterMCParticle.add(clsIndex, mc_relations.begin(), mc_relations.end());
435 if (!truehit_relations.empty()) relClusterTrueHit.add(clsIndex, truehit_relations.begin(), truehit_relations.end());
436 relClusterDigit.add(clsIndex, digit_weights.begin(), digit_weights.end());
437 }
438
439 m_cache->clear();
440}
441
443 const ClusterProjection& secondary, double minPitch, double centerPitch, double maxPitch)
444{
445 if (primary.getSize() >= static_cast<unsigned int>(m_sizeHeadTail)) { //Analog head tail
446 //Average charge in the central area
447 const double centerCharge = primary.getCenterCharge() / (primary.getSize() - 2);
448 const double minCharge = (primary.getMinCharge() < centerCharge) ? primary.getMinCharge() : centerCharge;
449 const double maxCharge = (primary.getMaxCharge() < centerCharge) ? primary.getMaxCharge() : centerCharge;
450 primary.setPos(0.5 * (primary.getMinPos() + primary.getMaxPos()
451 + (maxCharge * maxPitch - minCharge * minPitch) / centerCharge));
452 const double snHead = centerCharge / m_cutAdjacentSignal / minPitch;
453 const double snTail = centerCharge / m_cutAdjacentSignal / maxPitch;
454 const double landauHead = minCharge / centerCharge * minPitch;
455 const double landauTail = maxCharge / centerCharge * maxPitch;
456 primary.setError(0.5 * sqrt(1.0 / snHead / snHead + 1.0 / snTail / snTail
457 + 0.5 * landauHead * landauHead + 0.5 * landauTail * landauTail));
458 } else if (primary.getSize() <= 2) { // Add a phantom charge to second strip
459 primary.setError(centerPitch * (secondary.getSize() + 2) * m_cutAdjacentSignal / (primary.getCharge() +
460 (secondary.getSize() + 3) * m_cutAdjacentSignal));
461 } else {
462 const double sn = cls.getSeedCharge() / m_elNoise;
463 primary.setError(2.0 * centerPitch / sn);
464 }
465}
466
468 VxdID sensorID, unsigned int uCell, unsigned int vCell, double centerPitch,
469 bool isAtUEdge, bool isAtVEdge, bool isAdjacentDead)
470{
471 // Get bins from cell ID
472 unsigned short uBin = PXD::getBinU(sensorID, uCell, vCell, errorPar.getBinsU());
473 unsigned short vBin = PXD::getBinV(sensorID, vCell, errorPar.getBinsV());
474 // Get error from DB [in units of pix]
475 double error = errorPar.getContent(sensorID.getID(), uBin, vBin, primary.getSize());
476 double sf = 1.;
477 // Apply additional factor if at sensor edges or adjacent to dead rows/columns
478 if (isAtUEdge) sf *= errorPar.getSensorUEdgeFactor(sensorID.getID(), uBin, vBin, primary.getSize());
479 if (isAtVEdge) sf *= errorPar.getSensorVEdgeFactor(sensorID.getID(), uBin, vBin, primary.getSize());
480 if (isAdjacentDead) sf *= errorPar.getDeadNeighbourFactor(sensorID.getID(), uBin, vBin, primary.getSize());
481 // Set error (convert to [um])
482 if (error) primary.setError(sf * error * centerPitch); // zero means default values to use analytic error
483
484}
@ c_ErrorIfAlreadyRegistered
If the object/array was already registered, produce an error (aborting initialisation).
Definition DataStore.h:72
static std::string arrayName(const TClass *t, const std::string &name)
Return the storage name for an object of the given TClass and name.
Definition DataStore.cc:163
void setDescription(const std::string &description)
Sets the description of the module.
Definition Module.cc:214
void setPropertyFlags(unsigned int propertyFlags)
Sets the flags for the module properties.
Definition Module.cc:208
Module()
Constructor.
Definition Module.cc:30
@ c_ParallelProcessingCertified
This module can be run in parallel processing mode safely (All I/O must be done through the data stor...
Definition Module.h:80
The payload class for PXD cluster position error.
float getContent(unsigned short sensorID, unsigned short globalID) const
Get content.
unsigned short getBinsV() const
Get number of bins along sensor v side.
float getDeadNeighbourFactor(unsigned short sensorID, unsigned short globalID) const
Get scaling factor when neighbouring dead rows/column.
unsigned short getBinsU() const
Get number of bins along sensor u side.
float getSensorVEdgeFactor(unsigned short sensorID, unsigned short globalID) const
Get scaling factor at sensor edge in V.
float getSensorUEdgeFactor(unsigned short sensorID, unsigned short globalID) const
Get scaling factor at sensor edge in U.
The PXD digit class.
Definition PXDDigit.h:27
VxdID getSensorID() const
Get the sensor ID.
Definition PXDDigit.h:78
Class to remember recently assigned clusters This class will remember the current and the last pixel ...
Class representing a possible cluster during clustering of the PXD It supports merging of different c...
float getSeedCharge() const
get the seed charge of the cluster
Helper struct to collect information about the 1D projection of a Pixel cluster.
double getMinPos() const
Return the position of the minimum cell of the cluster.
void add(unsigned int cell, float position, float charge)
Add Pixel information to the projection.
double getMaxCharge() const
Return the charge in the maximum cell of the cluster.
void finalize()
Finish calculation of center of gravity and set correct cluster size.
unsigned int getMinCell() const
Return the minimum cell part of the cluster.
double getCharge() const
Return the total charge of the cluster.
unsigned int getMaxCell() const
Return the maximum cell part of the cluster.
void setPos(double pos)
Set the position of the cluster.
double getMinCharge() const
Return the charge in the minimum cell of the cluster.
void setError(double error)
Set the error of the cluster.
double getError() const
Return the error of the cluster.
double getMaxPos() const
Return the position of the maximum cell of the cluster.
double getPos() const
Return the projected position of the cluster.
double getCenterCharge() const
Return the center charge of the cluster, that is total charge minus minimum and maximum cell charge.
unsigned int getSize() const
Return the projected size of the cluster.
static PXDClusterPositionEstimator & getInstance()
Main (and only) way to access the PXDClusterPositionEstimator.
float computeEta(const std::set< Pixel > &pixels, int vStart, int vSize, double thetaU, double thetaV) const
Return the normed charge ratio between head and tail pixels (size>=2) or the charge of the seed (size...
int computeShapeIndex(const std::set< Pixel > &pixels, int uStart, int vStart, int vSize, double thetaU, double thetaV) const
Return the shape index of the pixels.
int getClusterkind(const PXDCluster &cluster) const
Return kind of cluster needed to find cluster position correction.
NoiseMap m_noiseMap
Noise map for the currently active sensor.
int m_clusterCacheSize
Size of cluster Cache (0 = default)
std::unique_ptr< DBObjPtr< PXDClusterPositionErrorPar > > m_clusterPositionErrorVPar
DB object for cluster position errors in V.
double m_cutAdjacentSignal
Signal in ADU for Adjacent cut, basically m_elNoise*m_cutAdjacent.
virtual void initialize() override
Initialize the module.
std::string m_positionErrorVName
Name of the DB payload containing cluster position errors in V.
StoreObjPtr< EventLevelTrackingInfo > m_eventLevelTrackingInfo
StoreObject to access the event level tracking information.
std::string m_relDigitMCParticleName
Name of the relation between PXDDigits and MCParticles.
virtual void event() override
do the clustering
void calculatePositionError(const ClusterCandidate &cls, ClusterProjection &primary, const ClusterProjection &secondary, double minPitch, double centerPitch, double maxPitch)
Calculate position and error for a given cluster.
std::vector< const RelationElement * > RelationLookup
Container for a RelationArray Lookup table.
std::string m_relClusterDigitName
Name of the relation between PXDClusters and PXDDigits.
bool m_createPXDClustersForAbortedTrackingEvents
bool to override the EventLevelTrackingInfo abort flag decision
double m_cutCluster
Cluster cut in sigma.
std::string m_storeTrueHitsName
Name of the collection to use for the PXDTrueHits.
PXDClusterizerModule()
Constructor defining the parameters.
void fillRelationMap(const RelationLookup &lookup, std::map< unsigned int, float > &relation, unsigned int index)
Add the relation from a given PXDDigit index to a map.
std::string m_storeMCParticlesName
Name of the collection to use for the MCParticles.
RelationLookup m_trueRelation
Lookup table for PXDDigit->PXDTrueHit relation.
std::string m_positionErrorUName
Name of the DB payload containing cluster position errors in U.
virtual void beginRun() override
do at every run change
std::unique_ptr< DBObjPtr< PXDClusterPositionErrorPar > > m_clusterPositionErrorUPar
DB object for cluster position errors in U.
void assignPositionErrorFromDB(ClusterProjection &primary, PXDClusterPositionErrorPar errorPar, VxdID sensorID, unsigned int uCell, unsigned int vCell, double centerPitch, bool isAtUEdge=false, bool isAtVEdge=false, bool isAdjacentDead=false)
Assign position error for a given cluster from DB.
int m_sizeHeadTail
Size of the cluster at which we switch from Center of Gravity to Analog Head Tail.
std::string m_storeDigitsName
Name of the collection to use for the PXDDigits.
void createRelationLookup(const RelationArray &relation, RelationLookup &lookup, size_t digits)
Create lookup maps for Relations We do not use the RelationIndex as we know much more about the relat...
std::string m_storeClustersName
Name of the collection to use for the PXDClusters.
void writeClusters(VxdID sensorID)
Write clusters to collection.
std::string m_relDigitTrueHitName
Name of the relation between PXDDigits and PXDTrueHits.
std::string m_relClusterMCParticleName
Name of the relation between PXDClusters and MCParticles.
RelationLookup m_mcRelation
Lookup table for PXDDigit->MCParticle relation.
std::unique_ptr< ClusterCache > m_cache
cache of the last seen clusters to speed up clustering
bool m_errorFromDB
Flag to set cluster position error from DB (default = true)
double m_cutAdjacent
Noise cut in sigma.
std::string m_relClusterTrueHitName
Name of the relation between PXDClusters and PXDTrueHits.
Class to represent one pixel, used in clustering for fast access.
Definition Pixel.h:36
unsigned short getU() const
Return the CellID in u.
Definition Pixel.h:68
unsigned short getV() const
Return the CellID in v.
Definition Pixel.h:70
Specific implementation of SensorInfo for PXD Sensors which provides additional pixel specific inform...
Definition SensorInfo.h:23
Low-level class to create/modify relations between StoreArrays.
void add(index_type from, index_type to, weight_type weight=1.0)
Add a new element to the relation.
void clear() override
Clear all elements from the relation.
Class to store a single element of a relation.
const std::string & getName() const
Return name under which the object is saved in the DataStore.
bool isRequired(const std::string &name="")
Ensure this array/object has been registered previously.
bool isOptional(const std::string &name="")
Tell the DataStore about an optional input.
bool registerInDataStore(DataStore::EStoreFlags storeFlags=DataStore::c_WriteOut)
Register the object/array in the DataStore.
Accessor to arrays stored in the data store.
Definition StoreArray.h:113
T * appendNew()
Construct a new T object at the end of the array.
Definition StoreArray.h:246
int getEntries() const
Get the number of objects in the array.
Definition StoreArray.h:216
void clear() override
Delete all entries in this array.
Definition StoreArray.h:207
const SensorInfoBase & getSensorInfo(Belle2::VxdID id) const
Return a reference to the SensorInfo of a given SensorID.
Definition GeoCache.cc:67
static GeoCache & getInstance()
Return a reference to the singleton instance.
Definition GeoCache.cc:214
Class to uniquely identify a any structure of the PXD and SVD.
Definition VxdID.h:32
baseType getID() const
Get the unique id.
Definition VxdID.h:93
void addParam(const std::string &name, T &paramVariable, const std::string &description, const T &defaultValue)
Adds a new parameter to the module.
Definition Module.h:559
#define REG_MODULE(moduleName)
Register the given module (without 'Module' suffix) with the framework.
Definition Module.h:649
double sqrt(double a)
sqrt for double
Definition beamHelpers.h:28
Namespace to encapsulate code needed for simulation and reconstrucion of the PXD.
unsigned short getBinU(VxdID id, unsigned int uid, unsigned int vid, unsigned short nBinsU)
Function to return a bin number for equal sized binning in U.
bool isClusterAtUEdge(VxdID id, unsigned int umin, unsigned int umax)
Helper function to check if one of the end pixels are at the edge of the sensor.
bool isClusterAtVEdge(VxdID id, unsigned int vmin, unsigned int vmax)
Helper function to check if one of the end pixels are at the edge of the sensor.
bool isClusterAtLadderJoint(VxdID id, unsigned int vmin, unsigned int vmax)
Helper function to check if one of the end pixels are at the ladder joint.
unsigned short getBinV(VxdID id, unsigned int vid, unsigned short nBinsV)
Function to return a bin number for equal sized binning in V.
Abstract base class for different kinds of events.
STL namespace.