24                                                        const std::string& 
name,
 
   36    ATH_CHECK( m_backExtrapolator.retrieve() );
 
   38    ATH_CHECK( m_alphaBetaEstimate.retrieve() );
 
   41    ATH_CHECK( m_ptFromAlphaBeta.retrieve() );
 
   44    ATH_CHECK( m_nswStationFitter.retrieve(DisableTool{m_nswStationFitter.empty()}) );
 
   46    return StatusCode::SUCCESS;
 
   53   m_use_mcLUT = use_mcLUT;
 
   57      if ( ptEndcapLUTSvc.retrieve().isFailure() ) {
 
   59         return StatusCode::FAILURE;
 
   62      m_alphaBetaEstimate->setMCFlag(m_use_mcLUT, &*ptEndcapLUTSvc);
 
   64      m_ptFromAlphaBeta->setMCFlag(m_use_mcLUT, &*ptEndcapLUTSvc);
 
   67      if ( ptEndcapLUTSvc.retrieve().isFailure() ) {
 
   69         return StatusCode::FAILURE;
 
   72      m_alphaBetaEstimate->setMCFlag(m_use_mcLUT, &*ptEndcapLUTSvc);
 
   74      m_ptFromAlphaBeta->setMCFlag(m_use_mcLUT, &*ptEndcapLUTSvc);
 
   77   ATH_MSG_DEBUG( 
"Completed tp set " << (m_use_mcLUT?
"MC":
"not MC") << 
" flag" );
 
   79   return StatusCode::SUCCESS;
 
   87                                                               std::vector<TrigL2MuonSA::TrackPattern>& v_trackPatterns)
 const 
   99         itTrack.isRpcFailure = !rpcFitResult.
isSuccess;
 
  105     return StatusCode::SUCCESS;
 
  114                                     std::vector<TrigL2MuonSA::TrackPattern>& v_trackPatterns,
 
  122       itTrack.phiMSDir = tgcFitResult.
phiDir;
 
  129       itTrack.isTgcFailure = 
true;
 
  134     if(!m_nswStationFitter.empty())
 
  135       ATH_CHECK( m_nswStationFitter->superPointFitter(p_roids, itTrack, stgcHits, mmHits) );
 
  138   return StatusCode::SUCCESS;
 
  147                                                               std::vector<TrigL2MuonSA::TrackPattern>& v_trackPatterns,
 
  155       itTrack.phiMSDir = tgcFitResult.
phiDir;
 
  158       itTrack.isTgcFailure = 
true;
 
  161     ATH_CHECK( superPointFitter(itTrack, muonRoad) );
 
  163     makeReferenceLine(itTrack, muonRoad);
 
  164     ATH_CHECK( m_alphaBetaEstimate->setAlphaBeta(p_roids, tgcFitResult, itTrack, muonRoad) );
 
  166     if ( itTrack.etaBin < -1 ) {
 
  168       if(itTrack.etaBin <= -1) itTrack.etaBin = 0;
 
  171     ATH_CHECK( m_ptFromAlphaBeta->setPt(itTrack,tgcFitResult) );
 
  173     double exInnerA = fromAlphaPtToInn(tgcFitResult,itTrack);
 
  174     double bw = muonRoad.
bw[3][0];
 
  175     double aw = muonRoad.
aw[3][0];
 
  176     if(std::abs(exInnerA) > 
ZERO_LIMIT) updateInnSP(itTrack, exInnerA, aw,bw);
 
  178     if(!m_nswStationFitter.empty())
 
  179       ATH_CHECK( m_nswStationFitter->superPointFitter(p_roids, itTrack, stgcHits, mmHits) );
 
  183   return StatusCode::SUCCESS;
 
  192     const unsigned int MAX_STATION = 10; 
 
  193     const float SIGMA              = 0.0080;
 
  194     const float DRIFTSPACE_LIMIT   = 16.;
 
  195     const int   MIN_MDT_FOR_FIT    = 3;
 
  201         if (mdtSegment.size()==0) 
continue;
 
  204         float sigma {0.}, phim {0.}, Xor {0.}, Yor {0.}, Ymid {0.}; 
 
  210             if (itMdtHit.isOutlier) 
continue;
 
  213                 Ymid = itMdtHit.cYmid;
 
  220             phim  = itMdtHit.cPhip;
 
  221             sigma = (std::abs(itMdtHit.DriftSigma) > 
ZERO_LIMIT)? itMdtHit.DriftSigma: SIGMA;
 
  224                 std::abs(itMdtHit.DriftSpace) < DRIFTSPACE_LIMIT &&
 
  229                 pbFitResult.
RILIN[
count] = itMdtHit.DriftSpace;
 
  240             for(
int i=0;
i<pbFitResult.
NPOI;
i++) {
 
  242                         << 
i << 
"/" << pbFitResult.
XILIN[
i] << 
"/" << pbFitResult.
YILIN[
i]
 
  243                         << 
"/" << pbFitResult.
RILIN[
i] << 
"/" << pbFitResult.
WILIN[
i]);
 
  249         if (
count >= MIN_MDT_FOR_FIT) {
 
  253             superPoint.Npoint = pbFitResult.
NPOI;
 
  260                     superPoint.R = (Ymid - Yor - pbFitResult.
BLIN)/pbFitResult.
ALIN  + Xor;
 
  261                     superPoint.Alin = 1./pbFitResult.
ALIN;
 
  262                     superPoint.Blin = -pbFitResult.
BLIN/pbFitResult.
ALIN;
 
  267                 superPoint.Z      = pbFitResult.
ALIN*(Ymid - Xor) + pbFitResult.
BLIN + Yor;
 
  268                 superPoint.Alin   = pbFitResult.
ALIN;
 
  269                 superPoint.Blin   = pbFitResult.
BLIN;
 
  272             superPoint.Phim   = phim;
 
  273             superPoint.Xor    = Xor;
 
  274             superPoint.Yor    = Yor;
 
  275             superPoint.Chi2   = pbFitResult.
CHI2;
 
  276             superPoint.PChi2  = pbFitResult.
PCHI2;
 
  277             for(
int i=0;
i<pbFitResult.
NPOI;
i++) superPoint.Residual[
i] =  pbFitResult.
RESI[
i];
 
  280         ATH_MSG_DEBUG(
"... Superpoint chamber/s_address/count/R/Z/Alin/Blin/Phim/Xor/Yor/Chi2/PChi2=" 
  282             << superPoint.R << 
"/" << superPoint.Z << 
"/" << superPoint.Alin << 
"/" 
  283             << superPoint.Blin << 
"/" << superPoint.Phim << 
"/" << superPoint.Xor << 
"/" 
  284             << superPoint.Yor << 
"/" << superPoint.Chi2 << 
"/" << superPoint.PChi2);
 
  287     return StatusCode::SUCCESS;
 
  295     const unsigned int MAX_STATION = 10;
 
  296     const float SIGMA              = 0.0080;
 
  297     const float DRIFTSPACE_LIMIT   = 16.;
 
  298     const int   MIN_MDT_FOR_FIT    = 3;
 
  306         if (mdtSegment.size()==0) 
continue;
 
  313             float sigma {0.}, phim {0.}, Xor {0.}, Yor {0.}, Ymid {0.}; 
 
  318                 if (itMdtHit.isOutlier) 
continue;
 
  322                     Ymid = itMdtHit.cYmid;
 
  329                 phim  = itMdtHit.cPhip;
 
  330                 sigma = (std::abs(itMdtHit.DriftSigma) > 
ZERO_LIMIT)? itMdtHit.DriftSigma: SIGMA;
 
  333                         std::abs(itMdtHit.DriftSpace) < DRIFTSPACE_LIMIT &&
 
  336                     pbFitResult.XILIN[
count] = itMdtHit.R - Xor;
 
  337                     pbFitResult.YILIN[
count] = itMdtHit.Z - Yor;
 
  338                     pbFitResult.RILIN[
count] = itMdtHit.DriftSpace;
 
  340                     pbFitResult.RESI[
count] = 0.;
 
  342                     pbFitResult.NPOI = 
count;
 
  348             ATH_MSG_DEBUG(
"... MDT hit used in fit #=" << pbFitResult.NPOI);
 
  349             for(
int i=0;
i<pbFitResult.NPOI;
i++) {
 
  351                     << 
i << 
"/" << pbFitResult.XILIN[
i] << 
"/" << pbFitResult.YILIN[
i]
 
  352                     << 
"/" << pbFitResult.RILIN[
i] << 
"/" << pbFitResult.WILIN[
i]);
 
  354             if (
count >= MIN_MDT_FOR_FIT) {
 
  355                 Evlfit( pbFitResult);
 
  357                 float bc = (Ymid - Xor);
 
  358                 float X = (pbFitResult.ALIN*bc)+pbFitResult.BLIN ;
 
  360                 superPoint.Npoint = pbFitResult.NPOI;
 
  362                     if (std::abs(pbFitResult.ALIN) > 
ZERO_LIMIT) {
 
  364                             superPoint.R = (Ymid-Yor)/pbFitResult.ALIN - pbFitResult.BLIN/pbFitResult.ALIN + Xor;
 
  365                             superPoint.Alin = 1./pbFitResult.ALIN;
 
  366                             superPoint.Blin = -pbFitResult.BLIN/pbFitResult.ALIN;
 
  368                             superPoint.R      = bc + Xor;
 
  369                             superPoint.Z      = 
X + Yor;
 
  370                             superPoint.Alin   = pbFitResult.ALIN;
 
  371                             superPoint.Blin   = pbFitResult.BLIN;
 
  375                 superPoint.Phim   = phim;
 
  376                 superPoint.Xor    = Xor;
 
  377                 superPoint.Yor    = Yor;
 
  378                 superPoint.Chi2   = pbFitResult.CHI2;
 
  379                 superPoint.PChi2  = pbFitResult.PCHI2;
 
  380                 for(
int i=0;
i<pbFitResult.NPOI;
i++) superPoint.Residual[
i] =  pbFitResult.RESI[
i];
 
  383             ATH_MSG_DEBUG(
"...Special Superpoint chamber/s_address/count/R/Z/Alin/Blin/Phim/Xor/Yor/Chi2/PChi2=" 
  385                     << superPoint.R << 
"/" << superPoint.Z << 
"/" << superPoint.Alin << 
"/" 
  386                     << superPoint.Blin << 
"/" << superPoint.Phim << 
"/" << superPoint.Xor << 
"/" 
  387                     << superPoint.Yor << 
"/" << superPoint.Chi2 << 
"/" << superPoint.PChi2);
 
  394         unsigned int  sumN = 0;
 
  396         if(
chamber==3) {nrWidth = m_rwidth_Endcapinn_first;}
 
  397         if(
chamber==4) {nrWidth = m_rwidth_Endcapmid_first;}
 
  398         if(
chamber==5) {nrWidth = m_rwidth_Endcapout_first;}
 
  401             if (std::abs(itMdtHit.DriftSpace) < m_mdt_driftspace_downlimit ||
 
  402                 std::abs(itMdtHit.DriftSpace) > m_mdt_driftspace_uplimit){
 
  403                 itMdtHit.isOutlier = 2;
 
  407             if(itMdtHit.isOutlier > 1)
continue;
 
  408             double Z = itMdtHit.Z;
 
  409             double R = itMdtHit.R;
 
  410             double nbw = aw*Z + bw;
 
  411             if (R>(nbw-nrWidth) && R<(nbw+nrWidth)){
 
  412                 itMdtHit.isOutlier = 0;
 
  416                 itMdtHit.isOutlier = 2;
 
  420         if(sumN==0) 
continue;
 
  422         stationSPFit(mdtSegment, superPoint,pbFitResult, trackPattern.
s_address,
chamber,aw);
 
  426     return StatusCode::SUCCESS;
 
  436     unsigned int i_layer_max = 0;
 
  438     float Xor, Yor, 
sigma, phim=0;
 
  440     const unsigned int MAX_STATION = 8;
 
  441     const float SIGMA              = 0.0080;
 
  442     const float DRIFTSPACE_LIMIT   = 16.;
 
  443     const int   MIN_MDT_FOR_FIT    = 3;
 
  445     const unsigned int MAX_LAYER = 12;
 
  448     std:: vector<unsigned int> MdtLayerHits_index[MAX_STATION][MAX_LAYER];
 
  450     for (
unsigned int i_st=0; i_st<MAX_STATION; i_st++) {
 
  451         for (
unsigned int i_ly=0; i_ly<MAX_LAYER; i_ly++) {
 
  453             MdtLayerHits_index[i_st][i_ly].clear();
 
  457     float Chbest = 1.e25;
 
  462     unsigned int sumN[8];
 
  465     for (
unsigned int i_st=0; i_st<8; i_st++) {
 
  473     for (
unsigned int i_hit=0; i_hit<mdtSegment.size(); i_hit++){
 
  475         if (mdtSegment.at(i_hit).isOutlier>1) 
continue;
 
  477         double Z = mdtSegment.at(i_hit).Z;
 
  478         double R = mdtSegment.at(i_hit).R;
 
  480         sumZ[i_station] = sumZ[i_station] + Z;
 
  481         sumR[i_station] = sumR[i_station] + R;
 
  484         if (sumN[i_station]!=0) {
 
  485             avZ[i_station] = sumZ[i_station]/sumN[i_station];
 
  486             avR[i_station] = sumR[i_station]/sumN[i_station];
 
  490      if (sumN[i_station]==0) 
return;
 
  492     for (
unsigned int i_hit=0; i_hit<mdtSegment.size(); i_hit++) {
 
  494         if (mdtSegment.at(i_hit).isOutlier>1) 
continue;
 
  496         double Z = mdtSegment.at(i_hit).Z;
 
  497         double R = mdtSegment.at(i_hit).R;
 
  499         if (i_station==3) nsWidth = m_rwidth_Endcapinn_second;
 
  500         if (i_station==4) nsWidth = m_rwidth_Endcapmid_second;
 
  501         if (i_station==5) nsWidth = m_rwidth_Endcapout_second;
 
  503         double nbw = aw*Z+(avR[i_station]-aw*avZ[i_station]);
 
  505         if ( R>(nbw-nsWidth) && R<(nbw+nsWidth) ) {
 
  506             mdtSegment.at(i_hit).isOutlier = 0;
 
  508             mdtSegment.at(i_hit).isOutlier = 2;
 
  513     for (
unsigned int i_hit=0; i_hit<mdtSegment.size(); i_hit++) {
 
  515         unsigned int i_layer =mdtSegment.at(i_hit).Layer;
 
  516         if (mdtSegment.at(i_hit).isOutlier>1) 
continue;
 
  517         if (i_layer > i_layer_max) i_layer_max = i_layer;
 
  520         MdtLayerHits_index[i_station][i_layer].push_back(i_hit);
 
  523     std::vector<unsigned int> Ly_1st{};
 
  524     std::vector<float> Line_A{};
 
  525     std::vector<float> Line_B{};
 
  526     std::vector<float> Line_Chi2{};
 
  527     std::vector<unsigned int> Line_count{};
 
  528     std::vector<float>  Line_Xor{};
 
  529     std::vector<float>Line_Yor{};
 
  530     std::vector<float>  Line_sum_Z{};
 
  531     std::vector<float>  Line_phim{};
 
  532     std::vector<float>  Maxlayers_A{};
 
  533     std::vector<float>  Maxlayers_B{};
 
  534     std::vector<float>  Maxlayers_Chi2{};
 
  535     std::vector<float>  Maxlayers_RESI{};
 
  536     float Maxlayers_Phim  = 0;
 
  537     float Maxlayers_R     = 0;
 
  538     float Maxlayers_Z     = 0;
 
  539     float Maxlayers_Xor   = 0;
 
  540     float Maxlayers_Yor   = 0;
 
  541     float Maxlayers_PChi2 = 0;
 
  544     for (
unsigned int i_layer=0; i_layer<=i_layer_max; i_layer++) {
 
  548         Ly_1st.push_back(i_layer);
 
  551     const int real_layer= Ly_1st.size();
 
  552     std::vector<std::vector<unsigned int> > Ly_flg{};
 
  554     for (
int pr=real_layer; 
pr>=3; 
pr--) {
 
  567         findLayerCombination(Ly_1st, real_layer, 
pr,Ly_flg, total_cp);
 
  569         for (
unsigned int i=0;
i<Ly_flg.size(); 
i++) {
 
  571             std::vector<std::vector<int> >tID{};
 
  572             std::vector<std::vector<int> >tIndex{};
 
  577             for (
int ti=0; ti<8; ti++) {
 
  578                 for (
int tj=0; tj<2; tj++) {
 
  580                     tubeindex[ti][tj] = 0;
 
  584             for (
unsigned int j=0; j<Ly_flg[
i].size(); j++) {
 
  586                 int i_layer = Ly_flg[
i][j];
 
  587                 std::vector<int> tid{};
 
  588                 std::vector<int> tindex{}; 
 
  592                 float tube_1st = 999999.;
 
  593                 float tube_2nd = 999999.;
 
  595                 int layer_2nd = 9999;
 
  597                 for (
unsigned int i_hit=0; i_hit< 
MdtLayerHits[i_station][i_layer]; i_hit++) {
 
  599                     unsigned int i_index = MdtLayerHits_index[i_station][i_layer].at(i_hit);
 
  601                     if (mdtSegment.at(i_index).isOutlier>1) 
continue;
 
  603                     float nbw3 = (mdtSegment.at(i_index).Z)*(aw) + (avR[i_station]-(aw)*avZ[i_station]) ;
 
  604                     float dis_tube = std::abs(std::abs(nbw3-mdtSegment.at(i_index).R)- mdtSegment.at(i_index).DriftSpace);
 
  606                     if (dis_tube<tube_1st) {
 
  608                         layer_2nd = layer_1st;
 
  611                     } 
else if (dis_tube<tube_2nd) {
 
  617                 if ( layer_1st != 9999 ) {
 
  618                     mdtSegment.at(layer_1st).isOutlier = 0;
 
  620                     tindex.push_back(layer_1st);
 
  623                 if ( layer_2nd != 9999 ) {
 
  624                     mdtSegment.at(layer_2nd).isOutlier = 1;
 
  626                     tindex.push_back(layer_2nd);
 
  630                 tIndex.push_back(tindex);
 
  633             for (
unsigned int ti=0; ti<tID.size();ti++) {
 
  634                 for (
unsigned int tj=0; tj<tID[ti].size();tj++) {
 
  635                 tube_ID[ti][tj]   = tID[ti][tj];
 
  636                 tubeindex[ti][tj] = tIndex[ti][tj];
 
  640             std::vector<int> isg;
 
  641             std::vector<int> hitarray;
 
  643             int ntry = (
int)floor(
pow(2.,
pr))-1;
 
  645             for (
int ntryi=0; ntryi<=ntry; ntryi++) {
 
  651                 for (
int ntryj=1; ntryj<=
pr; ntryj++) {
 
  653                     int Isg = (ntryi&(
int)
pow(2.,ntryj-1))? 1 : 0;
 
  655                     if (tube_ID[ntryj-1][Isg] != 0)  yhit = 1;
 
  656                     sumid = sumid * yhit;
 
  660                     for (
unsigned int tt=0;
tt<isg.size(); 
tt++) {
 
  661                         int tindex = tubeindex[
tt][isg[
tt]];
 
  662                         hitarray.push_back(tindex);
 
  670                 float sum_Z_used = 0.;
 
  671                 float sum_R_used = 0.;
 
  673                 if (hitarray.size()==0) 
continue;
 
  675                 for (itMdtHit=mdtSegment.begin(); itMdtHit!=mdtSegment.end(); ++itMdtHit) { 
 
  679                     if(mdtSegment.at(hit_index).isOutlier>1) 
continue;
 
  685                     for (
unsigned int j=0; j<hitarray.size(); j++) {
 
  687                         if (hitarray[j]==hit_index) {
 
  702                     phim  = itMdtHit->cPhip;
 
  703                     sigma = (std::abs(itMdtHit->DriftSigma) > 
ZERO_LIMIT)? itMdtHit->DriftSigma: SIGMA;
 
  706                         std::abs(itMdtHit->DriftSpace) < DRIFTSPACE_LIMIT &&
 
  707                         std::abs(itMdtHit->DriftTime) > 
ZERO_LIMIT ) {
 
  711                         pbFitResult.
RILIN[
count] = itMdtHit->DriftSpace;
 
  718                         sum_Z_used = sum_Z_used + itMdtHit->Z;
 
  719                         sum_R_used = sum_R_used + itMdtHit->R;
 
  726                 for (
int i=0;
i<pbFitResult.
NPOI;
i++) {
 
  728                     << 
i << 
"/" << pbFitResult.
XILIN[
i] << 
"/" << pbFitResult.
YILIN[
i]
 
  729                     << 
"/" << pbFitResult.
RILIN[
i] << 
"/" << pbFitResult.
WILIN[
i]);
 
  732                 if (
count >= MIN_MDT_FOR_FIT) {
 
  738                     for (
int cand=0; cand<6; cand++) {
 
  741                         Line_A.push_back(1/pbFitResult.
SlopeCand[cand]);
 
  743                         Line_Chi2.push_back(pbFitResult.
Chi2Cand[cand]);
 
  744                         Line_count.push_back(pbFitResult.
NPOI);
 
  745                         Line_Xor.push_back(Xor);
 
  746                         Line_Yor .push_back(Yor);
 
  747                         Line_sum_Z.push_back(sum_Z_used/
count);
 
  748                         Line_phim.push_back(phim);
 
  755         if (Line_Chi2.size()==0) 
continue;
 
  757         std::multimap<float, int>chi_map{};
 
  758         std::vector<float> t_A;
 
  759         std::vector<float> t_B;
 
  760         std::vector<float> t_Chi2;
 
  761         std::vector<float> t_count;
 
  762         std::vector<float> t_Xor;
 
  763         std::vector<float> t_Yor;
 
  764         std::vector<float> t_sum_Z;
 
  765         std::vector<float> t_phim;
 
  776         for (
unsigned int ir=0; 
ir<Line_Chi2.size(); 
ir++) chi_map.insert(std::make_pair(Line_Chi2.at(
ir), 
ir));
 
  779             t_A.push_back(Line_A.at(jt->second));
 
  780             t_B.push_back(Line_B.at(jt->second));
 
  781             t_Chi2.push_back(Line_Chi2.at(jt->second));
 
  782             t_count.push_back(Line_count.at(jt->second));
 
  783             t_Xor.push_back(Line_Xor.at(jt->second));
 
  784             t_Yor.push_back(Line_Yor.at(jt->second));
 
  785             t_sum_Z.push_back(Line_sum_Z.at(jt->second));
 
  786             t_phim.push_back(Line_phim.at(jt->second));
 
  788                 Maxlayers_A.push_back(Line_A.at(jt->second));
 
  789                 Maxlayers_B.push_back(Line_B.at(jt->second));
 
  790                 Maxlayers_Chi2.push_back(Line_Chi2.at(jt->second));
 
  794         superPoint.
Npoint = t_count[0];
 
  795         if(i_station==4 && 
pr==real_layer){
 
  796             Maxlayers_Z     = t_sum_Z[0];
 
  797             Maxlayers_R     = t_A[0]*t_sum_Z[0]+t_B[0];
 
  798             Maxlayers_Phim  = t_phim[0];
 
  799             Maxlayers_Xor   = t_Xor[0];
 
  800             Maxlayers_Yor   = t_Yor[0];
 
  801             Maxlayers_PChi2 = pbFitResult.
PCHI2;
 
  802             Maxlayers_N     = t_count[0];
 
  805         if (s_address == -1) { 
 
  808                 superPoint.
Z = t_sum_Z[0];
 
  809                 superPoint.
R = t_A[0]*t_sum_Z[0]+t_B[0];
 
  810                 superPoint.
Alin =t_A[0];
 
  811                 superPoint.
Blin =t_B[0];
 
  814             superPoint.
Phim   = t_phim[0];
 
  815             superPoint.
Xor    = t_Xor[0];
 
  816             superPoint.
Yor    = t_Yor[0];
 
  817             superPoint.
Chi2   = t_Chi2[0];
 
  822             for (
int cand=0; cand<6; cand++) {
 
  826                 superPoint.
Chi2Cand[cand]      = t_Chi2[cand];
 
  834             if ((i_station == 3 || i_station == 5) && 
pr==4 && Chbest > m_endcapmid_mdt_chi2_limit) {
 
  841                 for (
int cand=0; cand<6; cand++) {
 
  850         if (Chbest<m_endcapmid_mdt_chi2_limit){ 
 
  852             ATH_MSG_DEBUG(
"... Superpoint chamber/s_address/count/R/Z/Alin/Blin/Phim/Xor/Yor/Chi2/PChi2=" 
  853                     << i_station << 
"/" << s_address << 
"/" << 
count << 
"/" 
  854                     << superPoint.
R << 
"/" << superPoint.
Z << 
"/" << superPoint.
Alin << 
"/" 
  855                     << superPoint.
Blin << 
"/" << superPoint.
Phim << 
"/" << superPoint.
Xor << 
"/" 
  856                     << superPoint.
Yor << 
"/" << superPoint.
Chi2 << 
"/" << superPoint.
PChi2);
 
  860             if(i_station==4 && Maxlayers_A.size()>0){
 
  861                 superPoint.
Npoint = Maxlayers_N;
 
  862                 superPoint.
Z = Maxlayers_Z;
 
  863                 superPoint.
R = Maxlayers_R;
 
  864                 superPoint.
Alin =Maxlayers_A[0];
 
  865                 superPoint.
Blin =Maxlayers_B[0];
 
  866                 superPoint.
Phim   = Maxlayers_Phim;
 
  867                 superPoint.
Xor    = Maxlayers_Xor;
 
  868                 superPoint.
Yor    = Maxlayers_Yor;
 
  869                 superPoint.
Chi2   = Maxlayers_Chi2[0];
 
  870                 superPoint.
PChi2  = Maxlayers_PChi2;
 
  871                 for (
int cand=0; cand<6; cand++) {
 
  872                     if (std::abs(Maxlayers_A[cand]) > 
ZERO_LIMIT ) {
 
  873                         superPoint.
SlopeCand[cand]     = Maxlayers_A[cand];
 
  875                         superPoint.
Chi2Cand[cand]      = Maxlayers_Chi2[cand];
 
  890     const unsigned int MAX_STATION = 8;
 
  896     float Chi2_cand[8][6];
 
  900     const int middle = 4;
 
  903     for (
unsigned int i_station=4; i_station<MAX_STATION; i_station++) {
 
  910         spChi2[i_station]=0.;
 
  912         for (
unsigned int ci=0; ci<
NCAND; ci++) {
 
  913             A_cand[i_station][ci]    =0.;
 
  914             B_cand[i_station][ci]    =0.;
 
  915             Chi2_cand[i_station][ci] =0.;
 
  918         if (i_station<4 || i_station>5) 
continue; 
 
  921         aw[i_station] = muonRoad.
aw[i_station][0];
 
  922         spZ[i_station]  = superPoint.Z;
 
  923         spR[i_station]  = superPoint.R;
 
  924         spA[i_station]  = superPoint.Alin;
 
  925         spB[i_station]  = superPoint.Blin;
 
  927         for (
unsigned int cand=0; cand<
NCAND; cand++) {
 
  928             A_cand[i_station][cand]    = superPoint.SlopeCand[cand];
 
  929             B_cand[i_station][cand]    = superPoint.InterceptCand[cand];
 
  930             Chi2_cand[i_station][cand] = superPoint.Chi2Cand[cand];
 
  934     float test_diff  = 1.e25;
 
  935     float best_diff  = 1.e25;
 
  937     float match_midA = 0.;
 
  938     float match_outA = 0.;
 
  957         for (
int m=0; 
m<6; 
m++) {
 
  959             test_diff = std::abs(A_cand[middle][
m] - aw[middle]);
 
  961             if (test_diff<best_diff) {
 
  962                 best_diff      = test_diff;
 
  963                 rmatched       = A_cand[middle][
m];
 
  964                 spB[middle]    = B_cand[middle][
m];
 
  965                 spChi2[middle] = Chi2_cand[middle][
m];
 
  966                 spR[middle]    = rmatched * spZ[middle] + spB[middle];
 
  976         if(std::abs(spZ[outer]-spZ[middle]) > 
ZERO_LIMIT) sp_line = (spR[outer]-spR[middle])/(spZ[outer]-spZ[middle]);
 
  978         for (
int t=0; 
t<2; ++
t) {
 
  982                 for (
int i=0; 
i<6; ++
i) {
 
  983                     if (
t==0) test_diff = std::abs(A_cand[middle][
i] - sp_line);
 
  984                     else if (
t==1) test_diff = std::abs(A_cand[outer][
i] - sp_line);
 
  985                     if (test_diff<best_diff) {
 
  986                         best_diff = test_diff;
 
  988                             match_midA      = A_cand[middle][
i];
 
  989                             spB[middle]     = B_cand[middle][
i];
 
  990                             spChi2[middle]  = Chi2_cand[middle][
i];
 
  991                             spR[middle]     = match_midA * spZ[middle] + spB[middle];
 
  993                             match_outA     = A_cand[outer][
i];
 
  994                             spB[outer]     = B_cand[outer][
i];
 
  995                             spChi2[outer]  = Chi2_cand[outer][
i];
 
  996                             spR[outer]     = match_outA * spZ[outer] + spB[outer];
 
 1006             spA[middle] = match_midA;
 
 1007         } 
else if (std::abs(rmatched) > 
ZERO_LIMIT) {
 
 1008             spA[middle] = rmatched;
 
 1011         if (std::abs(match_outA) > 
ZERO_LIMIT) spA[outer] = match_outA;
 
 1014     for (
unsigned int i_station=4; i_station<MAX_STATION; i_station++) {
 
 1015         if (i_station<4 || i_station>5) 
continue; 
 
 1018             superPoint.
Alin =spA[i_station];
 
 1019             superPoint.Blin =spB[i_station];
 
 1020             superPoint.Chi2 =spChi2[i_station];
 
 1031     float MiddleSlope     = 0;
 
 1032     float OuterSlope      = 0;
 
 1034     for (
int i_station=4; i_station<6; i_station++) {
 
 1043         if ( superPoint.Npoint > 2 && superPoint.R > 0.) {
 
 1044             if ( i_station==4 ) { 
 
 1045                 MiddleSlope = superPoint.
Alin;
 
 1047             if ( i_station==5 ) {
 
 1048                 OuterSlope  = superPoint.Alin;
 
 1053     double mdtpT    = std::abs(tgcFitResult.
tgcPT);
 
 1058     } 
else if (std::abs(tgcFitResult.
tgcPT)>=8.0 && std::abs(MiddleSlope) > 
ZERO_LIMIT) {
 
 1062     mdtpT = (std::abs(tgcFitResult.
tgcPT)>1
e-5)? mdtpT*(tgcFitResult.
tgcPT/std::abs(tgcFitResult.
tgcPT)) : 0;
 
 1067     double extrInnerEta = 0;
 
 1072     muonSA->
setPt(mdtpT);
 
 1080     StatusCode sc = m_backExtrapolator->give_eta_phi_at_vertex(muonSA, 
eta,sigma_eta,
phi,sigma_phi,0.);
 
 1082     if (
sc.isSuccess() ){
 
 1085         extrInnerEta = etaMiddle;
 
 1092         naw = 
std::tan(
theta)*(std::abs(extrInnerEta)/extrInnerEta);
 
 1105     double nrWidth = m_rwidth_Endcapinn_first;
 
 1106     unsigned int  sumN[8];
 
 1108     for (
unsigned int i_st=0; i_st<8;i_st++) {
 
 1113     const int i_station = 3;
 
 1114     int chamberID = i_station;
 
 1119     if (mdtSegment.size()==0) 
return;
 
 1122         if (std::abs(itMdtHit.DriftSpace) < m_mdt_driftspace_downlimit ||
 
 1123             std::abs(itMdtHit.DriftSpace) > m_mdt_driftspace_uplimit){
 
 1124             itMdtHit.isOutlier = 2;
 
 1128         if (itMdtHit.isOutlier > 1) 
continue;
 
 1130         double Z = itMdtHit.Z;
 
 1131         double R = itMdtHit.R;
 
 1132         double nbw = tgc_aw*Z + bw;
 
 1134         if (R>(nbw-nrWidth) && R<(nbw+nrWidth)) {
 
 1135             itMdtHit.isOutlier = 0;
 
 1138             itMdtHit.isOutlier = 2;
 
 1143     if (sumN[i_station]==0) 
return;
 
 1145     stationSPFit(mdtSegment, superPoint,pbFitResult, trackPattern.
s_address,i_station, aw);
 
 1149     for (
int cand=0; cand<
NCAND; cand++) {
 
 1150         float ds=std::abs(superPoint.SlopeCand[cand]-aw);
 
 1153             superPoint.Alin = superPoint.SlopeCand[cand];
 
 1154             superPoint.Blin = superPoint.InterceptCand[cand];
 
 1155             superPoint.Chi2 = superPoint.Chi2Cand[cand];
 
 1166                                                              std::vector<std::vector<unsigned int> > &
c,
 
 1169   std::vector<unsigned int> 
b(
r,0);
 
 1171   findSubLayerCombination(
a,
n,
r,
b,0,
r,
c,nr);
 
 1182                                                                 std::vector<unsigned int> &
b,
 
 1185                                                                 std::vector<std::vector<unsigned int> > &
c,
 
 1191     std::vector<unsigned int> 
t;
 
 1196       for (
int j = 0;j<
r; j++) 
t.push_back(
b[j]);
 
 1203       findSubLayerCombination(
a,
n,
r,
b,
i+1,
num-1,
c,nr);
 
 1229     float Xnor = 1. / std::sqrt(1. + pbFitResult.
ALIN * pbFitResult.
ALIN);
 
 1231     for(
int j=0;j<pbFitResult.
NPOI;j++) {
 
 1233         float distj = (pbFitResult.
ALIN * pbFitResult.
XILIN[j] + pbFitResult.
BLIN - pbFitResult.
YILIN[j]) * Xnor;
 
 1234         float rlin = (distj>=0.) ? pbFitResult.
RILIN[j] : -pbFitResult.
RILIN[j];
 
 1235         pbFitResult.
RESI[j] = distj - rlin;
 
 1269                                                  const std::array<float, NMEAMX>& XI,
 
 1270                                                  const std::array<float, NMEAMX>& YI,
 
 1271                                                  const std::array<float, NMEAMX>& RI,
 
 1272                                                  const std::array<float, NMEAMX>& WI,
 
 1273                                                  float& 
A,
float& B, 
float& Chi2,
float& Pchi2)
 const 
 1275     if (Nmeas<=2||Nmeas>=
NMEAMX+1) 
return;
 
 1279     Xline(XI,YI,WI,Nmeas,
A,B,SAA,SBB,SAB,Chi2);
 
 1282     const float WIlim = 0.1 * (*std::max_element(WI.begin(), WI.begin() + Nmeas));
 
 1283     const int Ngood = std::count_if(WI.begin(), WI.begin() + Nmeas, [WIlim](
const float& 
w) {
 
 1287     std::vector<int> bestIdx{};
 
 1288     for (
int j=0;j<Nmeas;j++) {    
 
 1289         if (WI[j]>=WIlim || Ngood<=3) {    
 
 1290             if (bestIdx.size() < 4){
 
 1291                 bestIdx.push_back(j);
 
 1294                 bestIdx.at(bestIdx.size()-2) = bestIdx.back();
 
 1299     const int Nbest = bestIdx.size();   
 
 1301     const int Ntry = (1 << Nbest) - 1;    
 
 1302     std::array<float,NMEAMX> RRi{};
 
 1303     float Abest{0.}, Bbest {0.}, CHbest{1.e25};
 
 1304     for (
int j=0;j<=Ntry;j++) {   
 
 1305         for (
int k=0;
k<Nbest;
k++) {   
 
 1306             const int Isig = (j & (1<<
k)) ? 1 : 0;   
 
 1307             RRi[bestIdx[
k]] = (Isig==1) ? -RI[bestIdx[
k]] : RI[bestIdx[
k]];
 
 1312         Circfit(Nmeas,XI,YI,RRi,WI,Aj,Bj,chi2j,&bestIdx);
 
 1314         if (chi2j>=0.0&&chi2j<=CHbest) {
 
 1325     Circfit(Nmeas,XI,YI,RI,WI,
A,B,Chi2);
 
 1328         Pchi2 = TMath::Prob(Chi2, Nmeas - 2);
 
 1345                                                  const std::array<float, NMEAMX>& XI,
 
 1346                                                  const std::array<float, NMEAMX>& YI,
 
 1347                                                  const std::array<float, NMEAMX>& RI,
 
 1348                                                  const std::array<float, NMEAMX>& WI,
 
 1349                                                  float& 
A,
float& B,
float& Chi2, std::vector<int>* idx_vec)
 const 
 1351     const bool use_all {!idx_vec};     
 
 1352     std::vector<int> temp{};
 
 1355         std::iota(temp.begin(), temp.end(), 0);   
 
 1359     std::array<float,NMEAMX> XX{},YY{};
 
 1361     float Test, Toll {.1};
 
 1366         float Xnor  = 1. / std::sqrt(1. + 
A * 
A);
 
 1371             for(
const int& j : *idx_vec) {
 
 1372                 const int Epsi {(
A * XI[j] + B - YI[j]>=0.) ? 1 : -1};
 
 1373                 XX[j] = XI[j] - Epsi * Xnor * std::abs(RI[j]) * 
A;
 
 1374                 YY[j] = YI[j] + Epsi * Xnor * std::abs(RI[j]);
 
 1378             for(
const int& j : *idx_vec){
 
 1379                 XX[j] = XI[j] - Xnor * RI[j] * 
A;
 
 1380                 YY[j] = YI[j] + Xnor * RI[j];
 
 1384         Xline(XX,YY,WI,Nmeas,
A,B,SAA,SBB,SAB,Chi2, idx_vec);
 
 1386         Test = ((Aold-
A)*(Aold-
A))/ SAA + ((Bold-B)*(Bold-B))/ SBB;
 
 1388     } 
while(
Test>=Toll&&Niter<=20);
 
 1401                                                const std::array<float, NMEAMX>& 
Y,
 
 1402                                                const std::array<float, NMEAMX>& 
W,
 
 1404                                                float& 
A,
float& B,
float& SAA,
float& SBB,
float& SAB,
float& Square, std::vector<int>* idx_vec)
 const 
 1406     std::vector<int> temp{};
 
 1409         std::iota(temp.begin(), temp.end(), 0);
 
 1413     float S1{0.},SX{0.},SY{0.},SXX{0.},SXY{0.},SYY{0.};
 
 1415     for(
const int& j : *idx_vec) {
 
 1417         SX  = SX  + 
W[j] * 
X[j];
 
 1418         SY  = SY  + 
W[j] * 
Y[j];
 
 1419         SXX = SXX + 
W[j] * 
X[j] * 
X[j];
 
 1420         SXY = SXY + 
W[j] * 
X[j] * 
Y[j];
 
 1421         SYY = SYY + 
W[j] * 
Y[j] * 
Y[j];
 
 1424     float Deter  = 
S1 * SXX - SX * SX;
 
 1427         A      = (
S1 * SXY - SX * SY)  / Deter;
 
 1428         B      = (SY * SXX - SX * SXY) / Deter;    
 
 1434         A= (
S1 * SXY - SX * SY > 0.) ? 9.e+5 : -9.e+5;
 
 1435         B      = SY/
S1 - SX/
S1 * 
A;
 
 1441     for(
const int& j : *idx_vec) {
 
 1442         float DY =(
Y[j] - 
A * 
X[j] - B)/std::sqrt(1 + 
A * 
A);
 
 1443         Square = Square + 
W[j] * DY * DY;
 
 1451                                                   const std::array<float, NMEAMX>& XI,
 
 1452                                                   const std::array<float, NMEAMX>& YI,
 
 1453                                                   const std::array<float, NMEAMX>& RI,
 
 1454                                                   const std::array<float, NMEAMX>& WI,
 
 1455                                                   float& 
A,
float& B,
float& Chi2,
float& Pchi2,
 
 1456                                                   std::array<float, NCAND>& SlopeCand, 
 
 1457                                                   std::array<float, NCAND>& InterceptCand,
 
 1458                                                   std::array<float, NCAND>& Chi2Cand)
 const  
 1460     if (Nmeas<=2||Nmeas>=
NMEAMX+1) 
return;
 
 1463     float A0,B0, SAA,SBB,SAB,Square;  
 
 1464     Xline(XI,YI,WI,Nmeas,A0,B0,SAA,SBB,SAB,Square);
 
 1467     const float WIlim = 0.1 * (*std::max_element(WI.begin(), WI.begin() + Nmeas));
 
 1468     const int Ngood = std::count_if(WI.begin(), WI.begin() + Nmeas, [WIlim](
const float& 
w) {
 
 1472     std::vector<int> bestIdx{};
 
 1473     for (
int j=0;j<Nmeas;j++) {    
 
 1474         if (WI[j]>=WIlim || Ngood<=3) {    
 
 1475             if (bestIdx.size() < 4){
 
 1476                 bestIdx.push_back(j);
 
 1479                 bestIdx.at(bestIdx.size()-2) = bestIdx.back();
 
 1484     const int Nbest = bestIdx.size();   
 
 1486     std::vector<float> st_chi2{};
 
 1487     std::vector<float> st_A{};
 
 1488     std::vector<float> st_B{};
 
 1492         InterceptCand[
i] = 0.;
 
 1496     const int Ntry = (1 << Nbest) - 1;    
 
 1497     std::array<float,NMEAMX> RRi{}; 
 
 1498     float CHbest{1.e25}, Abest{0.}, Bbest{0.};
 
 1500     for (
int j=0;j<=Ntry;j++) {   
 
 1501         for (
int k=0;
k<Nbest;
k++) {   
 
 1502             int Isig = (j & (1<<
k)) ? 1 : 0;   
 
 1503             RRi[bestIdx[
k]] = (Isig==1) ? -RI[bestIdx[
k]] : RI[bestIdx[
k]];
 
 1508         Circfit(Nmeas,XI,YI,RRi,WI,Aj,Bj,Chi2,&bestIdx);
 
 1509         Circfit(Nmeas,XI,YI,RI,WI,Aj,Bj,Chi2);
 
 1510         st_A.push_back(Aj); st_B.push_back(Bj); st_chi2.push_back(Chi2);
 
 1512         if (Chi2>=0.0&&Chi2<=CHbest) {
 
 1519     std::multimap<float, int>chi_map {};
 
 1520     std::vector<float> t_A {};
 
 1521     std::vector<float> t_B {};
 
 1522     std::vector<float> t_chi2 {};
 
 1524     for (
size_t ir=0; 
ir<st_chi2.size(); 
ir++) chi_map.insert(std::make_pair(st_chi2.at(
ir), 
ir));
 
 1527         t_A.push_back(st_A.at(jt->second));
 
 1528         t_B.push_back(st_B.at(jt->second));
 
 1529         t_chi2.push_back(st_chi2.at(jt->second));
 
 1532     for (
int nv=0; nv<6; nv++) {
 
 1533         SlopeCand[nv]     = t_A[nv];
 
 1534         InterceptCand[nv] = t_B[nv];
 
 1535         Chi2Cand[nv]      = t_chi2[nv];
 
 1541     Circfit(Nmeas,XI,YI,RI,WI,
A,B,Chi2);
 
 1544         Pchi2 = TMath::Prob(Chi2, Nmeas - 2);