74 {
75 bool doMdtGasGainCorrectionForMc = false;
76
77
78 SG::ReadHandle<xAOD::EventInfo> eventInfo(
m_eventInfo);
79
80
82 doMdtGasGainCorrectionForMc = true;
83 }
84
85
86
88 if (!states) {
89 ATH_MSG_INFO(
"Cannot get track states on surface for TrackParticle");
90 return -9999.;
91 }
92
95
96 int nhitsadc = 0;
97 double absR = 0;
98 double datfit = 0;
99 double mcfit = 0;
101 double testEta = 0;
102 double maxhit = -999.;
103 float PhiFit = 0;
104
105 double meandADC = 0.;
106 std::vector<double> dADCvec;
107 dADCvec.clear();
108 double meanMDTdADC = 0.;
109
111 double track_phi =
track.perigeeParameters()->parameters()[
Trk::phi];
112
113 for (int nhits = 0; tsit != tsit_end; ++tsit, ++nhits) {
114
116
117 const Trk::MeasurementBase* measurement = (*tsit)->measurementOnTrack();
118 if (!measurement) { continue; }
121 continue;
122 }
123 if (!id.is_valid()) { continue; }
124
126 const Muon::MdtDriftCircleOnTrack* mdtcirc = dynamic_cast<const Muon::MdtDriftCircleOnTrack*>(measurement);
127 if (!mdtcirc) {
128 ATH_MSG_WARNING(
"cannot cast Trk::MeasurementBase to Muon::MdtDriftCircleOnTrack");
129 continue;
130 }
131 const Muon::MdtPrepData* rawdata = mdtcirc->
prepRawData();
132 if (rawdata) {
133
134
136
137 bool isInBme = fabs(testEta) >= 0.644 && fabs(testEta) <= 0.772 && track_phi >= -1.72 && track_phi <= -1.42;
138
139 if ((!isInBme && (absR <= 0.4 || absR >= 14.)) || (isInBme && (absR <= 0.4 || absR >= 6.5))) continue;
140 nhitsadc++;
141
142 if (fabs(testEta) < 1.) {
143 datfit = 105.088 + 27.0638 *
pow(absR, 1) - 4.72089 *
pow(absR, 2) + 0.110274 *
pow(absR, 3) +
144 0.041508 *
pow(absR, 4) - 0.00403678 *
pow(absR, 5) + 0.000111072 *
pow(absR, 6);
145 mcfit = 73.8974 + 12.0642 *
pow(absR, 1) + 0.975372 *
pow(absR, 2) - 0.922337 *
pow(absR, 3) +
146 0.140759 *
pow(absR, 4) - 0.00881345 *
pow(absR, 5) + 0.000202078 *
pow(absR, 6);
147 } else {
148 datfit = 106.329 + 26.5296 *
pow(absR, 1) - 4.07423 *
pow(absR, 2) - 0.0594686 *
pow(absR, 3) +
149 0.0608916 *
pow(absR, 4) - 0.00506114 *
pow(absR, 5) + 0.000131493 *
pow(absR, 6);
150 mcfit = 73.3596 + 12.9939 *
pow(absR, 1) + 0.494472 *
pow(absR, 2) - 0.812082 *
pow(absR, 3) +
151 0.128345 *
pow(absR, 4) - 0.00814268 *
pow(absR, 5) + 0.000188203 *
pow(absR, 6);
152 }
153
154 if (doMdtGasGainCorrectionForMc) {
157 dADCvec.push_back(correction * rawdata->
adc() - datfit);
158
160 } else {
161 meandADC += (rawdata->
adc() - datfit);
162 dADCvec.push_back(rawdata->
adc() - datfit);
163
164 maxhit = ((rawdata->
adc() - datfit) >= maxhit) ? (rawdata->
adc() - datfit) : maxhit;
165 }
166 }
167 }
168 }
169
170 if (nhitsadc == 0) { meandADC = -9999; }
171
172 if (nhitsadc == 1) { meandADC = meandADC /
double(nhitsadc); }
173 if (nhitsadc >= 2) { meandADC =
double(meandADC - maxhit) /
double(nhitsadc - 1); }
174
175 if (doMdtGasGainCorrectionForMc) {
176 if (track_phi > -3.2 && track_phi <= -2.87) PhiFit = -13.5471 *
pow(track_phi, 2) - 39.0001 * track_phi + 4.23613;
177 if (track_phi > -2.87 && track_phi <= -2.13) PhiFit = 53.209 *
pow(track_phi, 2) + 272.502 * track_phi + 342.867;
178 if (track_phi > -2.13 && track_phi <= -1.37) PhiFit = 72.4707 *
pow(track_phi, 2) + 254.911 * track_phi + 216.811;
179 if (track_phi > -1.37 && track_phi <= -0.57) PhiFit = 77.218 *
pow(track_phi, 2) + 149.676 * track_phi + 65.371;
180 if (track_phi > -0.57 && track_phi < 0.21) PhiFit = 86.977 *
pow(track_phi, 2) + 29.0558 * track_phi - 6.09313;
181 if (track_phi >= 0.21 && track_phi <= 1.00) PhiFit = 79.203 *
pow(track_phi, 2) - 92.9512 * track_phi + 21.6361;
182 if (track_phi > 1.00 && track_phi <= 1.79) PhiFit = 85.5711 *
pow(track_phi, 2) - 239.068 * track_phi + 161.918;
183 if (track_phi > 1.79 && track_phi <= 2.60) PhiFit = 82.8996 *
pow(track_phi, 2) - 362.665 * track_phi + 391.419;
184 if (track_phi > 2.60 && track_phi < 3.20) PhiFit = 73.8744 *
pow(track_phi, 2) - 443.274 * track_phi + 656.926;
185
186 meandADC += PhiFit - 5.01 + 9.46;
187 }
188
189 dADCvec.clear();
190
191 meanMDTdADC = meandADC;
192
193 return meanMDTdADC;
194 }
#define ATH_MSG_WARNING(x)
constexpr int pow(int base, int exp) noexcept
DataModel_detail::const_iterator< DataVector > const_iterator
double driftRadius() const
Returns the value of the drift radius.
virtual const MdtPrepData * prepRawData() const override final
Returns the PrepRawData used to create this corrected measurement.
int adc() const
Returns the ADC (typically range is 0 to 250)
@ Measurement
This is a measurement, and will at least contain a Trk::MeasurementBase.
@ Outlier
This TSoS contains an outlier, that is, it contains a MeasurementBase/RIO_OnTrack which was not used ...
@ IS_SIMULATION
true: simulation, false: data
DataVector< const Trk::TrackStateOnSurface > TrackStates
asinh(x)
helper methods ---------------------------------------------------------—