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
93 int nhitsadc = 0;
94 double absR = 0;
95 double datfit = 0;
96 double mcfit = 0;
98 double testEta = 0;
99 double maxhit = -999.;
100 float PhiFit = 0;
101
102 double meandADC = 0.;
103 std::vector<double> dADCvec;
104 dADCvec.clear();
105 double meanMDTdADC = 0.;
106
108 double track_phi =
track.perigeeParameters()->parameters()[
Trk::phi];
109
110 for (
const Trk::TrackStateOnSurface*
ts : *states) {
111
113
114 const Trk::MeasurementBase* measurement =
ts->measurementOnTrack();
115 if (!measurement) { continue; }
118 continue;
119 }
120 if (!id.is_valid()) { continue; }
121
123 const Muon::MdtDriftCircleOnTrack* mdtcirc = dynamic_cast<const Muon::MdtDriftCircleOnTrack*>(measurement);
124 if (!mdtcirc) {
125 ATH_MSG_WARNING(
"cannot cast Trk::MeasurementBase to Muon::MdtDriftCircleOnTrack");
126 continue;
127 }
128 const Muon::MdtPrepData* rawdata = mdtcirc->
prepRawData();
129 if (rawdata) {
130
131
133
134 bool isInBme = fabs(testEta) >= 0.644 && fabs(testEta) <= 0.772 && track_phi >= -1.72 && track_phi <= -1.42;
135
136 if ((!isInBme && (absR <= 0.4 || absR >= 14.)) || (isInBme && (absR <= 0.4 || absR >= 6.5))) continue;
137 nhitsadc++;
138
139 if (fabs(testEta) < 1.) {
140 datfit = 105.088 + 27.0638 *
pow(absR, 1) - 4.72089 *
pow(absR, 2) + 0.110274 *
pow(absR, 3) +
141 0.041508 *
pow(absR, 4) - 0.00403678 *
pow(absR, 5) + 0.000111072 *
pow(absR, 6);
142 mcfit = 73.8974 + 12.0642 *
pow(absR, 1) + 0.975372 *
pow(absR, 2) - 0.922337 *
pow(absR, 3) +
143 0.140759 *
pow(absR, 4) - 0.00881345 *
pow(absR, 5) + 0.000202078 *
pow(absR, 6);
144 } else {
145 datfit = 106.329 + 26.5296 *
pow(absR, 1) - 4.07423 *
pow(absR, 2) - 0.0594686 *
pow(absR, 3) +
146 0.0608916 *
pow(absR, 4) - 0.00506114 *
pow(absR, 5) + 0.000131493 *
pow(absR, 6);
147 mcfit = 73.3596 + 12.9939 *
pow(absR, 1) + 0.494472 *
pow(absR, 2) - 0.812082 *
pow(absR, 3) +
148 0.128345 *
pow(absR, 4) - 0.00814268 *
pow(absR, 5) + 0.000188203 *
pow(absR, 6);
149 }
150
151 if (doMdtGasGainCorrectionForMc) {
154 dADCvec.push_back(correction * rawdata->
adc() - datfit);
155
157 } else {
158 meandADC += (rawdata->
adc() - datfit);
159 dADCvec.push_back(rawdata->
adc() - datfit);
160
161 maxhit = ((rawdata->
adc() - datfit) >= maxhit) ? (rawdata->
adc() - datfit) : maxhit;
162 }
163 }
164 }
165 }
166
167 if (nhitsadc == 0) { meandADC = -9999; }
168
169 if (nhitsadc == 1) { meandADC = meandADC /
double(nhitsadc); }
170 if (nhitsadc >= 2) { meandADC =
double(meandADC - maxhit) /
double(nhitsadc - 1); }
171
172 if (doMdtGasGainCorrectionForMc) {
173 if (track_phi > -3.2 && track_phi <= -2.87) PhiFit = -13.5471 *
pow(track_phi, 2) - 39.0001 * track_phi + 4.23613;
174 if (track_phi > -2.87 && track_phi <= -2.13) PhiFit = 53.209 *
pow(track_phi, 2) + 272.502 * track_phi + 342.867;
175 if (track_phi > -2.13 && track_phi <= -1.37) PhiFit = 72.4707 *
pow(track_phi, 2) + 254.911 * track_phi + 216.811;
176 if (track_phi > -1.37 && track_phi <= -0.57) PhiFit = 77.218 *
pow(track_phi, 2) + 149.676 * track_phi + 65.371;
177 if (track_phi > -0.57 && track_phi < 0.21) PhiFit = 86.977 *
pow(track_phi, 2) + 29.0558 * track_phi - 6.09313;
178 if (track_phi >= 0.21 && track_phi <= 1.00) PhiFit = 79.203 *
pow(track_phi, 2) - 92.9512 * track_phi + 21.6361;
179 if (track_phi > 1.00 && track_phi <= 1.79) PhiFit = 85.5711 *
pow(track_phi, 2) - 239.068 * track_phi + 161.918;
180 if (track_phi > 1.79 && track_phi <= 2.60) PhiFit = 82.8996 *
pow(track_phi, 2) - 362.665 * track_phi + 391.419;
181 if (track_phi > 2.60 && track_phi < 3.20) PhiFit = 73.8744 *
pow(track_phi, 2) - 443.274 * track_phi + 656.926;
182
183 meandADC += PhiFit - 5.01 + 9.46;
184 }
185
186 dADCvec.clear();
187
188 meanMDTdADC = meandADC;
189
190 return meanMDTdADC;
191 }
#define ATH_MSG_WARNING(x)
constexpr int pow(int base, int exp) noexcept
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 ---------------------------------------------------------—