103{
104
105
106 v_trackPatterns.clear();
107 TrigL2MuonSA::TrackPattern trackPattern;
109
110
111 trackPattern.
s_address = (muonRoad.
isEndcap)? -1: muonRoad.Special + 2*muonRoad.LargeSmall;
112
113 const unsigned int MAX_STATION = 11;
114 const unsigned int MAX_LAYER = 12;
115
116 TrigL2MuonSA::MdtLayerHits v_mdtLayerHits[MAX_STATION][MAX_LAYER];
117
118 unsigned int i_layer_max = 0;
119 unsigned int chamber_max = 0;
120
122
123 for(unsigned int i_hit=0; i_hit<mdtHits.size(); i_hit++) {
124
125 unsigned int chamber = mdtHits[i_hit].Chamber;
126 if( chamber > chamber_max ) chamber_max =
chamber;
128
130 std::string
name =
m_idHelperSvc->mdtIdHelper().stationNameString(mdtHits[i_hit].name);
131 int chamber_this = 99;
132 int sector_this = 99;
133 bool isEndcap;
136 mdtHits[i_hit].isOutlier = 3;
137 continue;
138 }
139
142 double Z = mdtHits[i_hit].Z;
143 double R = mdtHits[i_hit].R;
145
147 mdtHits[i_hit].isOutlier = 0;
148 unsigned int i_layer = mdtHits[i_hit].Layer;
150
152 << chamber << "/" << Z << "/" << R << "/" << aw << "/" << bw << "/" << residual << "/" << rWidth);
153
154 if( std::abs(residual) > rWidth ) {
155 mdtHits[i_hit].isOutlier = 2;
156 continue;
157 }
158
160 if( mdtHits[i_hit].TubeLayer < 1 || i_layer > (MAX_LAYER-1) ) {
162 return StatusCode::FAILURE;
163 }
164 if( i_layer > i_layer_max ) i_layer_max = i_layer;
165
172 }
173
174 const double DeltaMin = 0.025;
175
178
180 while(1) {
181 if (chamber==9) break;
182 if (chamber==10) break;
183 unsigned int layer = 999999;
184 double DistMax = 0.;
185 double Residual = 0.;
186 unsigned int i_hit_max = 999999;
187 double ResMed = (v_mdtLayerHits[
chamber][0].
ntot!=0) ? v_mdtLayerHits[chamber][0].ResSum/v_mdtLayerHits[chamber][0].ntot : 0.;
188
189 for(unsigned int i_layer=0; i_layer<=i_layer_max; i_layer++) {
190 for(
unsigned int idigi=0; idigi<v_mdtLayerHits[
chamber][i_layer].
ndigi_all;idigi++) {
191
192 unsigned int i_hit = v_mdtLayerHits[
chamber][i_layer].
indexes[idigi];
193 if(mdtHits[i_hit].isOutlier > 0) continue;
194 double DistMed = std::abs(mdtHits[i_hit].Residual - ResMed);
195 if(DistMed>=DistMax) {
196 DistMax = DistMed;
197 Residual = mdtHits[i_hit].Residual;
199 i_hit_max = i_hit;
200 }
201 }
202 }
203 ATH_MSG_DEBUG(
"ResMed=" << ResMed <<
": DistMax/layer/i_hit_max/ntot="
204 << DistMax << "/" << layer << "/" << i_hit_max << "/" << v_mdtLayerHits[chamber][0].ntot);
205
206
207 if(layer == 999999) break;
208 if(v_mdtLayerHits[chamber][layer].ndigi==1) break;
209
210 double Mednew = (v_mdtLayerHits[
chamber][0].
ResSum - Residual)/(v_mdtLayerHits[chamber][0].ntot - 1);
211 double Delta = 2.*std::abs((ResMed - Mednew)/(ResMed + Mednew));
212 ATH_MSG_DEBUG(
"Mednew/Delta/DeltaMin=" << Mednew <<
"/" << Delta <<
"/" << DeltaMin);
213 if(Delta<=DeltaMin) break;
214
215
219 mdtHits[i_hit_max].isOutlier = 2;
220 }
221
222 ATH_MSG_DEBUG(
"choosing one at each layer, and record it in segment...");
224 mdtSegment.clear();
225
227 for(unsigned int i_layer=0; i_layer<=i_layer_max; i_layer++) {
228
229 ATH_MSG_DEBUG(
"i_layer=" << i_layer <<
": ndigi=" << v_mdtLayerHits[chamber][i_layer].ndigi);
230
231
232 while( v_mdtLayerHits[chamber][i_layer].ndigi>=2 ) {
233 double ResMax = 0.;
234 unsigned int i_hit_max = 999999;
235 for(
unsigned int idigi=0;idigi<v_mdtLayerHits[
chamber][i_layer].
ndigi_all;idigi++) {
236 unsigned int i_hit = v_mdtLayerHits[
chamber][i_layer].
indexes[idigi];
237 if( mdtHits[i_hit].isOutlier > 0 ) continue;
238 if( std::abs(mdtHits[i_hit].Residual) >= ResMax ) {
239 ResMax = std::abs(mdtHits[i_hit].Residual);
240 i_hit_max = i_hit;
241 }
242 }
243 ATH_MSG_DEBUG(
"ResMax=" << ResMax <<
": i_hit_max=" << i_hit_max);
244 if( i_hit_max == 999999 ) break;
248 mdtHits[i_hit_max].isOutlier = 1;
249 }
250
251
252 for(
unsigned int i_digi=0; i_digi< v_mdtLayerHits[
chamber][i_layer].
ndigi_all; i_digi++) {
253 unsigned int i_hit = v_mdtLayerHits[
chamber][i_layer].
indexes[i_digi];
254 if (i_layer==0)
phi0 = mdtHits[i_hit].cPhi0;
256 if (mdtHits[i_hit].isOutlier > 1) continue;
257 mdtSegment.push_back(mdtHits[i_hit]);
258 }
259 }
260 ATH_MSG_DEBUG(
"nr of hits in segment=" << mdtSegment.size());
261 }
262
263 v_trackPatterns.push_back(trackPattern);
264
265 return StatusCode::SUCCESS;
266}
#define ATH_MSG_WARNING(x)
static void find_station_sector(const std::string &name, int phi, bool &endcap, int &chamber, int §or)
double calc_residual(double aw, double bw, double x, double y) const
void doMdtCalibration(const EventContext &ctx, TrigL2MuonSA::MdtHitData &mdtHit, double track_phi, double phi0, bool isEndcap) const
double aw[N_STATION][N_SECTOR]
double bw[N_STATION][N_SECTOR]
double rWidth[N_STATION][N_LAYER]
TrigL2MuonSA::MdtHits mdtSegments[s_NCHAMBER]
constexpr uint8_t stationPhi
station Phi 1 to 8
std::vector< MdtHitData > MdtHits
@ MaxChamber
Number of measurement point definitions.
std::vector< unsigned int > indexes