17 declareInterface<IMDT_DigitizationTool>(
this);
24 return StatusCode::SUCCESS;
29 CLHEP::HepRandomEngine *rndmEngine)
const {
32 const Identifier DigitId =
input.getHitID();
33 const double maxTubeRadius{detMgr->getMdtReadoutElement(DigitId)->innerTubeRadius()};
48 double innerTubeRadius,
49 const Identifier& DigitId,
50 CLHEP::HepRandomEngine *rndmEngine)
const {
54 if (!calibConstants.isValid()) {
56 throw std::runtime_error(
"No Mdt calibration constants");
59 const RtRelationPtr&
data{calibConstants->getCalibData(DigitId, msgStream())->rtRelation};
63 bool outOfBound =
false;
72 time = trRelation->
tFromR(std::abs(measRadius), outOfBound);
86 ATH_MSG_WARNING(
"Drift velocity <=0 ! Time will not be smeared with resolution but will take the default measRadius-t value");
89 double timeWidth = radiusWidth / velocity;
92 double tUp = rtRelation->
tUpper();
95 double tmin =
time - 3.4 * timeWidth;
96 double tmax =
time + 3.4 * timeWidth;
98 bool outOfBound2 =
false;
99 if (tmin < 0.0) tmin = 0.0;
101 tmax = trRelation->
tFromR(innerTubeRadius, outOfBound2);
104 constexpr
double sqrt_one_over_two_pi = 0.39894228;
105 double p1r = 0.8480 *
std::exp(-0.5879 * measRadius);
111 t = CLHEP::RandFlat::shoot(rndmEngine, tmin, tmax);
113 gaussian = (1 - p1r) * sqrt_one_over_two_pi *
exp(-(
t -
time) * (
t -
time) / (2.0 * timeWidth * timeWidth));
114 if (gaussian >= CLHEP::RandFlat::shoot(rndmEngine, 0.0, 1.0) || cutoff > 200) {
flag = 1; }
118 ATH_MSG_DEBUG(
"t from r = " <<
time <<
" t resolution = " << timeWidth <<
"\nr resolution = " << radiusWidth
119 <<
" driftvelocity = " << velocity <<
"\nsmeared t = " <<
t <<
" cutoff = " << cutoff);
121 ATH_MSG_ERROR(
"Null pointer returned from CalibDBSvc. Returning 0.0 as drift time");
130 constexpr
double p0 = 57.38141;
131 constexpr
double p1 = 8.616943;
132 constexpr
double p2 = 2.497827;
133 constexpr
double p3 = -1.625900;
134 constexpr
double p4 = 0.3125281;
135 constexpr
double p5 = -0.02929554;
136 constexpr
double p6 = 0.001367115;
137 constexpr
double p7 = -0.00002541936;
143 constexpr
double g0 = 10.27808;
144 constexpr
double g1 = -0.3774593;
145 constexpr
double g2 = 0.02751001;
146 constexpr
double g3 = -0.0005994742;
151 double adc = CLHEP::RandGaussZiggurat::shoot(rndmEngine, adcfunc, adcWidth);