Belle II Software development
PXDDigitizerModule.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/pxdSimulation/PXDDigitizerModule.h>
10#include <vxd/geometry/GeoCache.h>
11#include <vxd/geometry/GeoTools.h>
12#include <pxd/reconstruction/PXDGainCalibrator.h>
13#include <pxd/reconstruction/PXDPixelMasker.h>
14
15#include <framework/logging/Logger.h>
16#include <framework/gearbox/Unit.h>
17#include <framework/gearbox/Const.h>
18#include <framework/datastore/DataStore.h>
19#include <framework/datastore/StoreArray.h>
20#include <framework/datastore/RelationArray.h>
21#include <framework/datastore/RelationIndex.h>
22
23#include <mdst/dataobjects/MCParticle.h>
24#include <pxd/dataobjects/PXDSimHit.h>
25#include <pxd/dataobjects/PXDTrueHit.h>
26#include <pxd/dataobjects/PXDDigit.h>
27#include <pxd/dataobjects/PXDInjectionBGTiming.h>
28#include <cmath>
29
30#include <TRandom.h>
31
32using namespace std;
33using namespace Belle2;
34using namespace Belle2::PXD;
35
36//-----------------------------------------------------------------
37// Register the Module
38//-----------------------------------------------------------------
39REG_MODULE(PXDDigitizer);
40
41//-----------------------------------------------------------------
42// Implementation
43//-----------------------------------------------------------------
44
49 , m_triggerGate(0), m_gated(false)
50{
51 //Set module properties
52 setDescription("Digitize PXDSimHits");
54 addParam("ElectronicEffects", m_applyNoise, "Apply electronic effects?",
55 true);
56 addParam("ElectronicNoise", m_elNoise,
57 "Noise added by the electronics, set in ENC", 150.0);
58
59 addParam("MCParticles", m_storeMCParticlesName,
60 "MCParticle collection name", string(""));
61 addParam("Digits", m_storeDigitsName, "Digits collection name", string(""));
62 addParam("SimHits", m_storeSimHitsName, "SimHit collection name",
63 string(""));
64 addParam("TrueHits", m_storeTrueHitsName, "TrueHit collection name",
65 string(""));
66
67 addParam("PoissonSmearing", m_applyPoisson,
68 "Apply Poisson smearing of electrons collected on pixels?", true);
69 addParam("IntegrationWindow", m_applyWindow, "Use integration window?",
70 true);
71
72 addParam("SegmentLength", m_segmentLength, "Maximum segment length (in mm)",
73 0.005);
74 addParam("ElectronGroupSize", m_elGroupSize,
75 "Split Signalpoints in smaller groups of N electrons (in e)", 100);
76 addParam("ElectronStepTime", m_elStepTime,
77 "Time step for tracking electron groups in readout plane (in ns)",
78 1.0);
79 addParam("ElectronMaxSteps", m_elMaxSteps,
80 "Maximum number of steps when propagating electrons", 200);
81 addParam("Gq", m_gq, "Gq of a pixel in nA/electron", 0.6);
82 addParam("ADCUnit", m_ADCUnit, "Slope of the linear ADC transfer curve in nA/ADU", 130.0);
83 addParam("PedestalMean", m_pedestalMean, "Mean of pedestals in ADU", 100.0);
84 addParam("PedestalRMS", m_pedestalRMS, "RMS of pedestals in ADU", 30.0);
85 addParam("GatingTime", m_gatingTime,
86 "Time window during which the PXD is not collecting charge in nano seconds", 400.0);
87 addParam("GatingWithoutReadout", m_gatingWithoutReadout,
88 "Digits from gated rows not sent to DHH ", true);
89 addParam("GatingWithoutReadoutTime", m_gatingWithoutReadoutTime,
90 "Time window during which digits from gated rows are not sent to DHH in nano seconds", 1400.0);
91 addParam("HardwareDelay", m_hwdelay,
92 "Constant time delay between bunch crossing and switching on triggergate in nano seconds", 0.0);
93
94}
95
97{
98 //Register output collections
100 storeDigits.registerInDataStore();
102 storeDigits.registerRelationTo(storeParticles);
104 storeDigits.registerRelationTo(storeTrueHits);
105
106 m_storePXDTiming.isOptional();
107
108 //Set default names for the relations
121
122 //Convert parameters to correct units
127
128 // Get number of PXD gates
129 auto gTools = VXD::GeoCache::getInstance().getGeoTools();
130 m_nGates = gTools->getNumberOfPXDReadoutGates();
131
132 // Active integration time (per gate) in ns
133 m_pxdIntegrationTime = 20000.0;
134 // Readout time per gate (sample - clear - cycle of rolling shutter) in ns
136
137 B2DEBUG(20, "PXDDigitizer Parameters (in system units, *=calculated +=set in xml):");
138 B2DEBUG(20, " --> ElectronicEffects: " << (m_applyNoise ? "true" : "false"));
139 B2DEBUG(20, " --> ElectronicNoise: " << m_elNoise << " e-");
140 B2DEBUG(20, " --> +ChargeThreshold: " << "set in xml by sensor, nominal 4 ADU");
141 B2DEBUG(20, " --> *NoiseFraction: " << "set in xml by sensor, nominal 1.0e-5");
142 B2DEBUG(20,
144 B2DEBUG(20,
146 B2DEBUG(20,
148 B2DEBUG(20,
150 B2DEBUG(20, " --> MCSimHitRel: " << m_relMCParticleSimHitName);
151 B2DEBUG(20, " --> DigitMCRel: " << m_relDigitMCParticleName);
152 B2DEBUG(20, " --> TrueSimRel: " << m_relTrueHitSimHitName);
153 B2DEBUG(20, " --> DigitTrueRel: " << m_relDigitTrueHitName);
154 B2DEBUG(20, " --> PoissonSmearing: " << (m_applyPoisson ? "true" : "false"));
155 B2DEBUG(20, " --> +IntegrationWindow: " << (m_applyWindow ? "true" : "false") << ", size defined in xml");
156 B2DEBUG(20, " --> SegmentLength: " << m_segmentLength << " cm");
157 B2DEBUG(20, " --> ElectronGroupSize: " << m_elGroupSize << " e-");
158 B2DEBUG(20, " --> ElectronStepTime: " << m_elStepTime << " ns");
159 B2DEBUG(20, " --> ElectronMaxSteps: " << m_elMaxSteps);
160 B2DEBUG(20, " --> ADU unit: " << m_eToADU << " e-/ADU");
161}
162
164{
165 //Fill map with all possible sensors This is too slow to be done every event so
166 //we fill it once and only clear the content of the sensors per event, not
167 //the whole map.
168 // NOTE: Some VXDIDs will be added to the map on the way if we have multiple
169 // frames, but they will stay and will be cleared appropriately, so this is not
170 // a problem, and after a few events the performance will stabilize.
171 m_sensors.clear();
173 for (VxdID layer : geo.getLayers(SensorInfo::PXD)) {
174 for (VxdID ladder : geo.getLadders(layer)) {
175 for (VxdID sensor : geo.getSensors(ladder)) {
176 m_sensors[sensor] = Sensor();
177 }
178 }
179 }
180}
181
183{
184
185 //Check for injection vetos
186 if (m_storePXDTiming.isValid()) {
187 m_gated = (*m_storePXDTiming).isGated();
188 m_triggerGate = (*m_storePXDTiming).getTriggerGate();
189 m_gatingStartTimes = (*m_storePXDTiming).getGatingStartTimes();
191
192 B2DEBUG(20, "Reading from PXDInjectionBGTiming: Event has TriggerGate " << m_triggerGate << " and is gated " << m_gated);
193
195 for (auto gatingStartTime : m_gatingStartTimes) {
196 B2DEBUG(20, "Gating start time " << gatingStartTime << "ns ");
197
198 // Gated readout can only come from gating in the time period from t=0us to t=20us.
199 if (gatingStartTime >= 0) {
200 // Compute the gate where gating started. Note the wrap around of the rolling shutter
201 int firstGated = (m_triggerGate + int(gatingStartTime / m_timePerGate)) % m_nGates ;
202 // Compute the gate where gating stopped. Note the wrap around of the rolling shutter
203 int lastGated = (m_triggerGate + int((gatingStartTime + m_gatingWithoutReadoutTime) / m_timePerGate)) % m_nGates ;
204
205 if (lastGated >= firstGated) {
206 // Gated channels in the same frame
207 B2DEBUG(20, "Gated channel interval from " << firstGated << " to " << lastGated);
208 m_gatedChannelIntervals.push_back(pair<int, int>(firstGated, lastGated));
209 } else {
210 // Gated channels in two consecutive frames
211 B2DEBUG(20, "Gated channel interval from " << firstGated << " to " << m_nGates);
212 B2DEBUG(20, "Gated channel interval from " << 0 << " to " << lastGated);
213 m_gatedChannelIntervals.push_back(pair<int, int>(firstGated, m_nGates));
214 m_gatedChannelIntervals.push_back(pair<int, int>(0, lastGated));
215 }
216 }
217 }
218 }
219 } else {
220 m_gated = false;
221 m_triggerGate = gRandom->Integer(m_nGates);
222 m_gatingStartTimes.clear();
224 }
225
226 //Clear sensors and process SimHits
227 for (Sensors::value_type& sensor : m_sensors) {
228 sensor.second.clear();
229 }
230 m_currentSensor = 0;
232 //Clear return value
234
238
239 RelationArray mcParticlesToSimHits(storeMCParticles, storeSimHits,
241 RelationArray mcParticlesToTrueHits(storeMCParticles, storeTrueHits); // not used here
242 RelationArray trueHitsToSimHits(storeTrueHits, storeSimHits,
244
246 storeDigits.clear();
247
248 RelationArray relDigitMCParticle(storeDigits, storeMCParticles,
250 relDigitMCParticle.clear();
251
252 RelationArray relDigitTrueHit(storeDigits, storeTrueHits,
254 relDigitTrueHit.clear();
255
256 unsigned int nSimHits = storeSimHits.getEntries();
257 if (nSimHits == 0)
258 return;
259
260 RelationIndex<MCParticle, PXDSimHit> relMCParticleSimHit(storeMCParticles,
261 storeSimHits, m_relMCParticleSimHitName);
262 RelationIndex<PXDTrueHit, PXDSimHit> relTrueHitSimHit(storeTrueHits,
263 storeSimHits, m_relTrueHitSimHitName);
264
265 //Check sensor info and set pointers to current sensor
266 for (unsigned int i = 0; i < nSimHits; ++i) {
267 m_currentHit = storeSimHits[i];
268 if (!m_currentHit->getElectrons()) continue;
270 relMCParticleSimHit.getFirstElementTo(m_currentHit);
271 if (mcRel) {
273 if (mcRel->weight < 0) {
274 //This simhit is from a particle which was not saved by the simulation
275 //so we do not take it into account for relations. Otherwise we might
276 //end up adding positive and negative weights together
278 }
279 } else {
280 // Don't warn if this is a background SimHit
281 if (m_currentHit->getBackgroundTag() == BackgroundMetaData::bg_none)
282 B2WARNING(
283 "Could not find MCParticle which produced PXDSimhit " << i);
285 }
287 relTrueHitSimHit.getFirstElementTo(m_currentHit);
288 //We only care about true hits from particles which have not been ignored
289 if (trueRel && trueRel->weight > 0) {
290 m_currentTrueHit = trueRel->indexFrom;
291 } else {
292 m_currentTrueHit = -1;
293 }
294
295 VxdID sensorID = m_currentHit->getSensorID();
296 if (!m_currentSensorInfo || sensorID != m_currentSensorInfo->getID()) {
297 m_currentSensorInfo = dynamic_cast<const SensorInfo*>(&VXD::GeoCache::getInstance().getSensorInfo(sensorID));
299 B2FATAL(
300 "SensorInformation for Sensor " << sensorID << " not found, make sure that the geometry is set up correctly");
301 m_currentSensor = &m_sensors[sensorID];
302 B2DEBUG(20, "Sensor Parameters for Sensor " << sensorID << ": " << endl
303 << " --> Width: " << m_currentSensorInfo->getWidth() << endl
304 << " --> Length: " << m_currentSensorInfo->getLength() << endl
305 << " --> uPitch: " << m_currentSensorInfo->getUPitch() << endl
306 << " --> vPitch " << m_currentSensorInfo->getVPitch(-m_currentSensorInfo->getLength() / 2.0)
307 << ", " << m_currentSensorInfo->getVPitch(m_currentSensorInfo->getLength() / 2.0) << endl
308 << " --> Thickness: " << m_currentSensorInfo->getThickness() << endl
309 << " --> BulkDoping: " << m_currentSensorInfo->getBulkDoping() << endl
310 << " --> BackVoltage: " << m_currentSensorInfo->getBackVoltage() << endl
311 << " --> TopVoltage: " << m_currentSensorInfo->getTopVoltage() << endl
312 << " --> SourceBorder: " << m_currentSensorInfo->getSourceBorder(-0.4 * m_currentSensorInfo->getLength())
313 << ", " << m_currentSensorInfo->getSourceBorder(+0.4 * m_currentSensorInfo->getLength()) << endl
314 << " --> ClearBorder: " << m_currentSensorInfo->getClearBorder(-0.4 * m_currentSensorInfo->getLength())
315 << ", " << m_currentSensorInfo->getClearBorder(+0.4 * m_currentSensorInfo->getLength()) << endl
316 << " --> DrainBorder: " << m_currentSensorInfo->getDrainBorder(-0.4 * m_currentSensorInfo->getLength())
317 << ", " << m_currentSensorInfo->getDrainBorder(+0.4 * m_currentSensorInfo->getLength()) << endl
318 << " --> GateDepth: " << m_currentSensorInfo->getGateDepth() << endl
319 << " --> DoublePixel: " << m_currentSensorInfo->getDoublePixel() << endl);
320
321 }
322 B2DEBUG(10,
323 "Processing hit " << i << " in Sensor " << sensorID << ", related to MCParticle " << m_currentParticle);
324 processHit();
325 }
326
328 saveDigits();
329 //Return number of created digits
330 setReturnValue(storeDigits.getEntries());
331}
332
334{
335 if (m_applyWindow) {
336 //Ignore a hit which is outside of the active frame time of the hit gate
337 float currentHitTime = m_currentHit->getGlobalTime();
338
339 // First we have to now the current hit gate
340 const ROOT::Math::XYZVector startPoint = m_currentHit->getPosIn();
341 const ROOT::Math::XYZVector stopPoint = m_currentHit->getPosOut();
342 auto row = m_currentSensorInfo->getVCellID(0.5 * (stopPoint.Y() + startPoint.Y()));
343
344 if (m_currentSensorInfo->getID().getLayerNumber() == 1)
345 row = 767 - row;
346 int gate = row / 4;
347
348 // Compute the distance in unit of gates to the trigger gate. Do it in
349 // direction of PXD rolling shutter.
350 int distanceToTriggerGate = 0;
351 if (gate >= m_triggerGate) distanceToTriggerGate = gate - m_triggerGate;
352 else distanceToTriggerGate = 192 - m_triggerGate + gate;
353
354 // The triggergate is switched ON at t0=0ns. Compute the integration time window for the current hit gate.
355 // Time intervals of subsequent readout gates are shifted by timePerGate.
356 float gateMinTime = -m_pxdIntegrationTime + m_timePerGate + m_hwdelay + distanceToTriggerGate * m_timePerGate ;
357 float gateMaxTime = gateMinTime + m_pxdIntegrationTime;
358
359 B2DEBUG(30,
360 "Checking if hit is in timeframe " << gateMinTime << " <= " << currentHitTime <<
361 " <= " << gateMaxTime);
362
363 if (currentHitTime
364 < gateMinTime
365 || currentHitTime
366 > gateMaxTime)
367 return;
368
369 if (m_gated) {
370 //Ignore simHits when the PXD is gated, i.e. PXD does not collect any charge.
371 for (auto gatingStartTime : m_gatingStartTimes) {
372 B2DEBUG(30,
373 "Checking if hit is in gated timeframe " << gatingStartTime << " <= " << currentHitTime <<
374 " <= " << gatingStartTime + m_gatingTime);
375
376 if (currentHitTime
377 > gatingStartTime
378 && currentHitTime
379 < gatingStartTime + m_gatingTime)
380 return;
381 }
382 }
383 }
384
385 //Get step length and direction
386 const ROOT::Math::XYZVector startPoint = m_currentHit->getPosIn();
387 ROOT::Math::XYZVector stopPoint = m_currentHit->getPosOut();
388 double dx = stopPoint.X() - startPoint.X();
389 double dy = stopPoint.Y() - startPoint.Y();
390 double dz = stopPoint.Z() - startPoint.Z();
391 double trackLength2 = dx * dx + dy * dy + dz * dz;
392
393 // Set magnetic field to save calls to getBField()
394 m_currentBField = m_currentSensorInfo->getBField(0.5 * (startPoint + stopPoint));
395
396 if (m_currentHit->getPDGcode() == Const::photon.getPDGCode() || trackLength2 <= 0.01 * Unit::um * Unit::um) {
397 //Photons deposit the energy at the end of their step
398 driftCharge(stopPoint, m_currentHit->getElectrons());
399 } else {
400 //Otherwise, split into segments of (default) max. 5µm and
401 //drift the charges from the center of each segment
402 auto segments = m_currentHit->getElectronsConstantDistance(m_segmentLength);
403 double lastFraction {0};
404 double lastElectrons {0};
405
406 for (auto& segment : segments) {
407 //Simhit returns step fraction and cumulative electrons. We want the
408 //center of these steps and electrons in this step
409 const double f = (segment.first + lastFraction) * 0.5;
410 const double e = segment.second - lastElectrons;
411 //Update last values
412 std::tie(lastFraction, lastElectrons) = segment;
413
414 //And drift charge from that position
415 stopPoint.SetXYZ(startPoint.X() + f * dx, startPoint.Y() + f * dy, startPoint.Z() + f * dz);
416 driftCharge(stopPoint, e);
417 }
418 }
419}
420
421inline double pol3(double x, const double* c)
422{
423 return c[0] + x * (c[1] + x * (c[2] + x * c[3]));
424};
425
426void PXDDigitizerModule::driftCharge(const ROOT::Math::XYZVector& r, double electrons)
427{
428 //Get references to current sensor/info for ease of use
429 const SensorInfo& info = *m_currentSensorInfo;
430 Sensor& sensor = *m_currentSensor;
431
432 static const double cx[] = {0, 2.611995e-01, -1.948316e+01, 9.167064e+02};
433 static const double cu[] = {4.223817e-04, -1.041434e-03, 5.139596e-02, -1.180229e+00};
434 static const double cv[] = {4.085681e-04, -8.615459e-04, 5.706439e-02, -1.774319e+00};
435 static const double ct[] = {2.823743e-04, -1.138627e-06, 4.310888e-09, -1.559542e-11}; // temp range [250, 350] degrees with m_elStepTime = 1 ns
436
437 double dz = 0.5 * info.getThickness() - info.getGateDepth() - r.Z(), adz = fabs(dz);
438
439 double dx = fabs(m_currentBField.Y()) > Unit::T ? copysign(pol3(adz, cx),
440 dz) : 0; // two configurations: with default B field (B~1.5T) and B=0T
441 double sigmaDrift_u = pol3(adz, cu);
442 double sigmaDrift_v = pol3(adz, cv);
443 double sigmaDiffus = pol3(info.getTemperature() - 300.0,
444 ct); // sqrt(2 * Const::uTherm * info.getElectronMobility(0) * m_elStepTime);
445
446 int nGroups = (int)(electrons / m_elGroupSize) + 1;
447 double groupCharge = electrons / nGroups;
448 while (nGroups--) {
449 double du = gRandom->Gaus(0.0, sigmaDrift_u), dv = gRandom->Gaus(0.0, sigmaDrift_v);
450 double uPos = r.X() + du + dx;
451 double vPos = r.Y() + dv;
452 int step = 0;
453 for (; step < m_elMaxSteps; ++step) {
454 int id = info.getTrappedID(uPos, vPos);
455 //Check if cloud inside of IG region
456 if (id >= 0) {
457 // Trapped in the internal gate region
458 sensor[Digit(id % 250, id / 250)].add(groupCharge, m_currentParticle, m_currentTrueHit);
459 break;
460 }
461 //Random walk with drift
462 dv = gRandom->Gaus(0.0, sigmaDiffus), du = gRandom->Gaus(0.0, sigmaDiffus);
463 uPos += du;
464 vPos += dv;
465 }
466 if (step == m_elMaxSteps) {
467 // cloud is not trapped but save it anyway into the closest pixel
468 int iu = info.getUCellID(uPos, vPos, 1), iv = info.getVCellID(vPos, 1);
469 sensor[Digit(iu, iv)].add(groupCharge, m_currentParticle, m_currentTrueHit);
470 }
471 }
472}
473
474
475double PXDDigitizerModule::addNoise(double charge)
476{
477 if (charge <= 0) {
478 //Noise Pixel, add chargeThreshold electrons;
479 charge = 1.1 * m_chargeThresholdElectrons;
480 } else {
481 if (m_applyPoisson) {
482 // For big charge assume Gaussian distr.
483 if (charge > (1000. * Unit::e))
484 charge = gRandom->Gaus(charge, sqrt(charge));
485 else
486 // Otherwise Poisson distr.
487 charge = gRandom->PoissonD(charge);
488 }
489 if (m_applyNoise) {
490 charge += gRandom->Gaus(0., m_elNoise);
491 }
492 }
493 return charge;
494}
495
497{
498 if (!m_applyNoise)
499 return;
500
501 for (Sensors::value_type& sensor : m_sensors) {
502 Sensor& s = sensor.second;
503
504 //Calculate the number of pixels on an empty sensor which will exceed the noise cut
505 const SensorInfo& info =
506 dynamic_cast<const SensorInfo&>(VXD::GeoCache::getInstance().getSensorInfo(sensor.first));
507 int nU = info.getUCells();
508 int nV = info.getVCells();
509 int nPixels = gRandom->Poisson(info.getNoiseFraction() * nU * nV);
510
511 // Existing digits will get noise in PXDDIgitizer::AddNoise.
512 // Here we add zero charge to nPixels digits.
513 // In part, this will create some new zero-charge digits.
514 // If we happen to select an existing (fired) digit, it remains unchanged.
515 for (int i = 0; i < nPixels; ++i) {
516 Digit d(gRandom->Integer(nU), gRandom->Integer(nV));
517 s[d].add(0.0);
518 }
519 }
520}
521
523{
524 return std::any_of(m_gatedChannelIntervals.begin(), m_gatedChannelIntervals.end(),
525 [gate](auto interval) {return gate >= interval.first && gate < interval.second;});
526}
527
529{
533 RelationArray relDigitMCParticle(storeDigits, storeMCParticles,
535 RelationArray relDigitTrueHit(storeDigits, storeTrueHits,
537
538
539 for (Sensors::value_type& sensor : m_sensors) {
540 auto sensorID = sensor.first;
541 const SensorInfo& info =
542 dynamic_cast<const SensorInfo&>(VXD::GeoCache::getInstance().getSensorInfo(sensorID));
543 m_chargeThreshold = info.getChargeThreshold();
545 for (Sensor::value_type& digitAndValue : sensor.second) {
546 const Digit& d = digitAndValue.first;
547 const DigitValue& v = digitAndValue.second;
548
549 //Add Noise where applicable
550 double charge = addNoise(v.charge());
551
552 // Draw a pedestal value
553 double pedestal = std::max(gRandom->Gaus(m_pedestalMean, m_pedestalRMS), 0.0);
554
555 // Get gain correction
556 double gain = PXDGainCalibrator::getInstance().getGainCorrection(sensorID, d.u(), d.v());
557
558 // Add pedestal to charge
559 charge = round(charge * gain / m_eToADU + pedestal);
560
561 // Clipping of ADC codes at 255
562 charge = std::min(charge, 255.0);
563
564 // Subtraction of integer pedestal in DHP
565 charge = charge - round(pedestal);
566
567 // Zero Suppression in DHP
568 if (charge <= m_chargeThreshold)
569 continue;
570
571 // Check if the readout digits is coming from a gated row
573 int row = d.v();
574 if (sensorID.getLayerNumber() == 1)
575 row = 767 - row;
576 int gate = row / 4;
577 if (checkIfGated(gate)) continue;
578 }
579
580 // Check if the readout digit is coming from a masked or dead area
581 if (PXD::PXDPixelMasker::getInstance().pixelDead(sensorID, d.u(), d.v())
582 || !PXD::PXDPixelMasker::getInstance().pixelOK(sensorID, d.u(), d.v())) {
583 continue;
584 }
585
586 //Add the digit to datastore
587 int digIndex = storeDigits.getEntries();
588 storeDigits.appendNew(
589 PXDDigit(sensorID, d.u(), d.v(), charge));
590
591 //If the digit has any relations to MCParticles, add the Relation
592 if (v.particles().size() > 0) {
593 relDigitMCParticle.add(digIndex, v.particles().begin(),
594 v.particles().end());
595 }
596 //If the digit has any truehits to TrueHit, add the Relation
597 if (v.truehits().size() > 0) {
598 relDigitTrueHit.add(digIndex, v.truehits().begin(),
599 v.truehits().end());
600 }
601 }
602 }
603}
604
605
static const ParticleType photon
photon particle
Definition Const.h:673
static std::string relationName(const std::string &fromName, const std::string &toName, std::string const &namedRelation="")
Return storage name for a relation between two arrays of the given names.
Definition DataStore.h:180
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
void setReturnValue(int value)
Sets the return value for this module as integer.
Definition Module.cc:220
@ 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 PXD digit class.
Definition PXDDigit.h:27
Class representing the charge and particle contributions for one pixel.
const relations_map & truehits() const
Return the map containing all truehit contributions to the pixel charge.
const relations_map & particles() const
Return the map containing all particle contributions to the pixel charge.
double charge() const
Return the charge collected in the pixel.
Class to represent the coordinates of one pixel.
double m_pedestalRMS
RMS pedestal in ADU.
void processHit()
Process one PXDSimHit by dividing the step in smaller steps and drifting the charge.
double m_pxdIntegrationTime
Integration time for each gate of the PXD in ns.
double m_pedestalMean
Mean pedestal in ADU.
bool m_gatingWithoutReadout
Digits from gated rows not sent to DHH.
void initialize() override final
Initialize the module and check the parameters.
std::string m_relTrueHitSimHitName
Name of the relation between PXDTrueHits and PXDSimHits.
void saveDigits()
Save all digits to the datastore.
int m_currentParticle
Index of the particle which caused the current hit.
double m_ADCUnit
Slope of the linear ADC transfer curve in nA/ADU.
Sensors m_sensors
Structure containing all existing sensors.
double m_elStepTime
Timeframe for one random walk step.
std::string m_relDigitMCParticleName
Name of the relation between PXDDigits and MCParticles.
int m_elMaxSteps
Maximum number of random walks before abort.
bool m_applyWindow
Whether or not to apply a time window cut.
std::string m_storeTrueHitsName
Name of the collection for the PXDTrueHits.
double m_timePerGate
Time needed to sample and clear a readout gate.
double m_gatingWithoutReadoutTime
Time window during which digits from gated rows are not sent to DHH.
void driftCharge(const ROOT::Math::XYZVector &position, double electrons)
Drift the charge inside the silicon.
const PXDSimHit * m_currentHit
Pointer to the PXDSimhit currently digitized.
std::vector< float > m_gatingStartTimes
Vector of start times for gating.
std::string m_storeMCParticlesName
Name of the collection for the MCParticles.
double m_noiseFraction
Fraction of noisy pixels per sensor.
bool m_applyNoise
Whether or not to apply noise.
void addNoiseDigits()
Add pure noise digits to the Sensors.
void event() override final
Digitize one event.
int m_nGates
Number of readout gates (or total number of Switcher channels)
double m_chargeThreshold
Zero-suppression threshold in ADU.
bool m_applyPoisson
Whether or not to apply poission fluctuation of charge.
Sensor * m_currentSensor
Pointer to the sensor in which the current hit occurred.
double m_chargeThresholdElectrons
... and its equivalent in electrons
std::vector< std::pair< int, int > > m_gatedChannelIntervals
Vector of gated readout channels.
ROOT::Math::XYZVector m_currentBField
Current magnetic field.
double m_elNoise
Amount of noise to apply.
std::string m_storeDigitsName
Name of the collection for the PXDDigits.
const SensorInfo * m_currentSensorInfo
Pointer to the SensorInfo of the current sensor.
std::string m_relDigitTrueHitName
Name of the relation between PXDDigits and PXDTrueHits.
double m_gatingTime
Time window during which the PXD is not collecting charge.
void beginRun() override final
Initialize the list of existing PXD Sensors.
double m_eToADU
ENC equivalent of 1 ADU.
bool checkIfGated(int gate)
Check if gate was read while in gated mode.
double m_gq
g_q of a pixel in nA/electrons.
StoreObjPtr< PXDInjectionBGTiming > m_storePXDTiming
Input array for timings.
double m_hwdelay
Hardware delay between time of bunch crossing and switching on triggergate in ns.
int m_elGroupSize
Max number of electrons per random walk.
int m_currentTrueHit
Index of the TrueHit the current hit belongs to.
std::string m_relMCParticleSimHitName
Name of the relation between MCParticles and PXDSimHits.
double addNoise(double charge)
Calculate the noise contribution to one pixel with given charge.
std::string m_storeSimHitsName
Name of the collection for the PXDSimhits.
float getGainCorrection(VxdID id, unsigned int uid, unsigned int vid) const
Get gain correction.
static PXDGainCalibrator & getInstance()
Main (and only) way to access the PXDGainCalibrator.
static PXDPixelMasker & getInstance()
Main (and only) way to access the PXDPixelMasker.
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.
Provides access to fast ( O(log n) ) bi-directional lookups on a specified relation.
const Element * getFirstElementTo(const TO &to) const
Return a pointer to the first relation Element of the given object.
RelationIndexContainer< FROM, TO >::Element Element
Struct representing a single element in the index.
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
bool registerRelationTo(const StoreArray< TO > &toArray, DataStore::EDurability durability=DataStore::c_Event, DataStore::EStoreFlags storeFlags=DataStore::c_WriteOut, const std::string &namedRelation="") const
Register a relation to the given StoreArray.
Definition StoreArray.h:140
static const double mm
[millimeters]
Definition Unit.h:70
static const double e
Standard of [electric charge].
Definition Unit.h:53
static const double um
[micrometers]
Definition Unit.h:71
static const double ns
Standard of [time].
Definition Unit.h:48
static const double T
[tesla]
Definition Unit.h:120
Class to facilitate easy access to sensor information of the VXD like coordinate transformations or p...
Definition GeoCache.h:38
const std::set< Belle2::VxdID > getLayers(SensorInfoBase::SensorType sensortype=SensorInfoBase::VXD)
Return a set of all known Layers.
Definition GeoCache.cc:176
const SensorInfoBase & getSensorInfo(Belle2::VxdID id) const
Return a reference to the SensorInfo of a given SensorID.
Definition GeoCache.cc:67
const std::set< Belle2::VxdID > & getSensors(Belle2::VxdID ladder) const
Return a set of all sensor IDs belonging to a given ladder.
Definition GeoCache.cc:204
static GeoCache & getInstance()
Return a reference to the singleton instance.
Definition GeoCache.cc:214
const GeoTools * getGeoTools()
Return a raw pointer to a GeoTools object.
Definition GeoCache.h:141
const std::set< Belle2::VxdID > & getLadders(Belle2::VxdID layer) const
Return a set of all ladder IDs belonging to a given layer.
Definition GeoCache.cc:193
int getUCells() const
Return number of pixel/strips in u direction.
Class to uniquely identify a any structure of the PXD and SVD.
Definition VxdID.h:32
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.
std::map< Digit, DigitValue > Sensor
Map of all hits in one Sensor.
Abstract base class for different kinds of events.
STL namespace.
RelationElement::index_type indexFrom
index of the element from which the relation points.
RelationElement::weight_type weight
weight of the relation.