ATLAS Offline Software
Public Member Functions | Protected Member Functions | Private Types | Private Member Functions | Private Attributes | List of all members
TrigL2MuonSA::NswStationFitter Class Reference

#include <NswStationFitter.h>

Inheritance diagram for TrigL2MuonSA::NswStationFitter:
Collaboration diagram for TrigL2MuonSA::NswStationFitter:

Public Member Functions

 NswStationFitter (const std::string &type, const std::string &name, const IInterface *parent)
 
StatusCode superPointFitter (const TrigRoiDescriptor *p_roids, TrigL2MuonSA::TrackPattern &trackPattern, TrigL2MuonSA::StgcHits &stgcHits, TrigL2MuonSA::MmHits &mmHits) const
 
StatusCode selectStgcHits (const TrigRoiDescriptor *p_roids, TrigL2MuonSA::StgcHits &stgcHits) const
 
StatusCode selectMmHits (const TrigRoiDescriptor *p_roids, TrigL2MuonSA::MmHits &mmHits) const
 
StatusCode calcWeightedSumHit (TrigL2MuonSA::TrackPattern &trackPattern) const
 
StatusCode findStgcHitsInSegment (TrigL2MuonSA::StgcHits &stgcHits) const
 
void findSetOfStgcHitIds (TrigL2MuonSA::StgcHits &stgcHits, std::array< std::vector< int >, 8 > hitIdByLayer, std::vector< std::array< int, 8 >> &hitIdsCandidate) const
 
StatusCode findMmHitsInSegment (TrigL2MuonSA::MmHits &mmHits) const
 
void findSetOfMmHitIds (TrigL2MuonSA::MmHits &mmHits, std::array< std::vector< int >, 8 > hitIdByLayer, std::vector< std::array< int, 8 >> &hitIdsCandidate) const
 
StatusCode MakeSegment (TrigL2MuonSA::TrackPattern &trackPattern, TrigL2MuonSA::StgcHits &stgcHits) const
 
StatusCode MakeSegment (TrigL2MuonSA::TrackPattern &trackPattern, TrigL2MuonSA::MmHits &mmHits) const
 
StatusCode calcMergedHit (TrigL2MuonSA::TrackPattern &trackPattern) const
 
void LinearFit (std::vector< double > &x, std::vector< double > &y, double *slope, double *intercept, double *mse) const
 
void LinearFitWeight (std::vector< double > &x, std::vector< double > &y, std::vector< bool > &isStgc, double *slope, double *intercept, double *mse, double eta) const
 
void getNswResolution (double *stgcDeltaR, double *mmDeltaR, unsigned int size) const
 
ServiceHandle< StoreGateSvc > & evtStore ()
 The standard StoreGateSvc (event store) Returns (kind of) a pointer to the StoreGateSvc. More...
 
const ServiceHandle< StoreGateSvc > & evtStore () const
 The standard StoreGateSvc (event store) Returns (kind of) a pointer to the StoreGateSvc. More...
 
const ServiceHandle< StoreGateSvc > & detStore () const
 The standard StoreGateSvc/DetectorStore Returns (kind of) a pointer to the StoreGateSvc. More...
 
virtual StatusCode sysInitialize () override
 Perform system initialization for an algorithm. More...
 
virtual StatusCode sysStart () override
 Handle START transition. More...
 
virtual std::vector< Gaudi::DataHandle * > inputHandles () const override
 Return this algorithm's input handles. More...
 
virtual std::vector< Gaudi::DataHandle * > outputHandles () const override
 Return this algorithm's output handles. More...
 
Gaudi::Details::PropertyBase & declareProperty (Gaudi::Property< T > &t)
 
Gaudi::Details::PropertyBase * declareProperty (const std::string &name, SG::VarHandleKey &hndl, const std::string &doc, const SG::VarHandleKeyType &)
 Declare a new Gaudi property. More...
 
Gaudi::Details::PropertyBase * declareProperty (const std::string &name, SG::VarHandleBase &hndl, const std::string &doc, const SG::VarHandleType &)
 Declare a new Gaudi property. More...
 
Gaudi::Details::PropertyBase * declareProperty (const std::string &name, SG::VarHandleKeyArray &hndArr, const std::string &doc, const SG::VarHandleKeyArrayType &)
 
Gaudi::Details::PropertyBase * declareProperty (const std::string &name, T &property, const std::string &doc, const SG::NotHandleType &)
 Declare a new Gaudi property. More...
 
Gaudi::Details::PropertyBase * declareProperty (const std::string &name, T &property, const std::string &doc="none")
 Declare a new Gaudi property. More...
 
void updateVHKA (Gaudi::Details::PropertyBase &)
 
MsgStream & msg () const
 
MsgStream & msg (const MSG::Level lvl) const
 
bool msgLvl (const MSG::Level lvl) const
 

Protected Member Functions

void renounceArray (SG::VarHandleKeyArray &handlesArray)
 remove all handles from I/O resolution More...
 
std::enable_if_t< std::is_void_v< std::result_of_t< decltype(&T::renounce)(T)> > &&!std::is_base_of_v< SG::VarHandleKeyArray, T > &&std::is_base_of_v< Gaudi::DataHandle, T >, void > renounce (T &h)
 
void extraDeps_update_handler (Gaudi::Details::PropertyBase &ExtraDeps)
 Add StoreName to extra input/output deps as needed. More...
 

Private Types

typedef ServiceHandle< StoreGateSvc > StoreGateSvc_t
 

Private Member Functions

Gaudi::Details::PropertyBase & declareGaudiProperty (Gaudi::Property< T > &hndl, const SG::VarHandleKeyType &)
 specialization for handling Gaudi::Property<SG::VarHandleKey> More...
 
Gaudi::Details::PropertyBase & declareGaudiProperty (Gaudi::Property< T > &hndl, const SG::VarHandleKeyArrayType &)
 specialization for handling Gaudi::Property<SG::VarHandleKeyArray> More...
 
Gaudi::Details::PropertyBase & declareGaudiProperty (Gaudi::Property< T > &hndl, const SG::VarHandleType &)
 specialization for handling Gaudi::Property<SG::VarHandleBase> More...
 
Gaudi::Details::PropertyBase & declareGaudiProperty (Gaudi::Property< T > &t, const SG::NotHandleType &)
 specialization for handling everything that's not a Gaudi::Property<SG::VarHandleKey> or a <SG::VarHandleKeyArray> More...
 

Private Attributes

StoreGateSvc_t m_evtStore
 Pointer to StoreGate (event store by default) More...
 
StoreGateSvc_t m_detStore
 Pointer to StoreGate (detector store by default) More...
 
std::vector< SG::VarHandleKeyArray * > m_vhka
 
bool m_varHandleArraysDeclared
 

Detailed Description

Definition at line 20 of file NswStationFitter.h.

Member Typedef Documentation

◆ StoreGateSvc_t

typedef ServiceHandle<StoreGateSvc> AthCommonDataStore< AthCommonMsg< AlgTool > >::StoreGateSvc_t
privateinherited

Definition at line 388 of file AthCommonDataStore.h.

Constructor & Destructor Documentation

◆ NswStationFitter()

TrigL2MuonSA::NswStationFitter::NswStationFitter ( const std::string &  type,
const std::string &  name,
const IInterface *  parent 
)

Definition at line 26 of file NswStationFitter.cxx.

28  :
29  AthAlgTool(type,name,parent)
30 {
31 }

Member Function Documentation

◆ calcMergedHit()

StatusCode TrigL2MuonSA::NswStationFitter::calcMergedHit ( TrigL2MuonSA::TrackPattern &  trackPattern) const

Definition at line 924 of file NswStationFitter.cxx.

925 {
926  TrigL2MuonSA::StgcHits stgcHits = trackPattern.stgcSegment;
927  TrigL2MuonSA::MmHits mmHits = trackPattern.mmSegment;
928 
929  double side_mm = 0;
930  std::vector<double> r, z;
931  std::vector<bool> isStgc;
932  for(unsigned int iHit = 0; iHit < stgcHits.size(); ++iHit) {
933  if (stgcHits.at(iHit).channelType == 1) {
934  r.push_back(stgcHits.at(iHit).r);
935  z.push_back(stgcHits.at(iHit).z);
936  isStgc.push_back(true);
937  }
938  }
939  for(unsigned int iHit = 0; iHit < mmHits.size(); ++iHit) {
940  if (mmHits.at(iHit).layerNumber < 2 || mmHits.at(iHit).layerNumber > 5) {
941  r.push_back(mmHits.at(iHit).r);
942  z.push_back(mmHits.at(iHit).z);
943  isStgc.push_back(false);
944  side_mm = std::abs(mmHits.at(iHit).z)/mmHits.at(iHit).z;
945  }
946  }
947  double slopefit=0., interceptfit=99999., mse=-1.;
948  LinearFit(z,r,&slopefit,&interceptfit,&mse);
949  ATH_MSG_DEBUG("@@Merge@@ stgc_mmX_fit slope= " << slopefit);
950  ATH_MSG_DEBUG("@@Merge@@ stgc_mmX_fit intercept= " << interceptfit);
951  ATH_MSG_DEBUG("@@Merge@@ stgc_mmX_fit mse= " << mse);
952 
953  std::vector<double> phiLocal;
954  double localPhiCenter = 3.*M_PI;
955  for (unsigned int iHit = 0; iHit < stgcHits.size(); ++iHit) {
956  if (stgcHits.at(iHit).channelType != 2) {continue;}
957  if (stgcHits.at(iHit).stationPhi<=5) {
958  localPhiCenter = 0.25 * M_PI * ((double)stgcHits.at(iHit).stationPhi-1.);
959  } else {
960  localPhiCenter = 0.25 * M_PI * ((double)stgcHits.at(iHit).stationPhi-9.);
961  }
962 
963  if (stgcHits.at(iHit).stationName == 57){
964  localPhiCenter += M_PI/8.; // small sTGC sectors
965  if (stgcHits.at(iHit).stationPhi == 5) localPhiCenter -= 2. * M_PI;
966  }
967 
968  double phiProj = stgcHits.at(iHit).phi - localPhiCenter;
969  if (phiProj > M_PI) phiProj -= 2.0*M_PI;
970  if (phiProj < -1.*M_PI) phiProj += 2.0*M_PI;
971  double rInterpolate = slopefit * stgcHits.at(iHit).z + interceptfit;
972  double rProj = stgcHits.at(iHit).r;
973  phiLocal.push_back( std::atan(rProj/rInterpolate*std::tan(phiProj)) );
974  ATH_MSG_DEBUG("@@Merge@@ philocalwire " << stgcHits.at(iHit).stationPhi << " " << stgcHits.at(iHit).stationName << " "
975  << localPhiCenter << " " << stgcHits.at(iHit).phi << " " << stgcHits.at(iHit).r << " " << stgcHits.at(iHit).z);
976  ATH_MSG_DEBUG("@@Merge@@ philocalwire " << rProj << " " << rInterpolate << " " << std::tan(phiProj) );
977  ATH_MSG_DEBUG("@@Merge@@ philocalwire " << std::atan(rProj/rInterpolate*std::tan(phiProj)) );
978  }
979  double tanTiltAngleU = 0,
980  tanTiltAngleV = 0;
981  double cosTiltAngleU = 0,
982  cosTiltAngleV = 0;
983  double sinTiltAngleU = 0,
984  sinTiltAngleV = 0;
985  if(side_mm > ZERO_LIMIT){
986  tanTiltAngleU = TanM1p5,
987  tanTiltAngleV = TanP1p5;
988  cosTiltAngleU = CosM1p5,
989  cosTiltAngleV = CosP1p5;
990  sinTiltAngleU = SinM1p5,
991  sinTiltAngleV = SinP1p5;
992  } else if(side_mm < -1.*ZERO_LIMIT){
993  tanTiltAngleU = TanP1p5,
994  tanTiltAngleV = TanM1p5;
995  cosTiltAngleU = CosP1p5,
996  cosTiltAngleV = CosM1p5;
997  sinTiltAngleU = SinP1p5,
998  sinTiltAngleV = SinM1p5;
999  } else {
1000  ATH_MSG_DEBUG("@@Merge@@ no U, V layer hits -> not consider tilt of U/V layers");
1001  }
1002  for (unsigned int iHit = 0; iHit < mmHits.size(); ++iHit) {
1003  if (localPhiCenter > 2.*M_PI) {
1004  if (mmHits.at(iHit).stationPhi<=5) {
1005  localPhiCenter = 0.25 * M_PI * ((double)mmHits.at(iHit).stationPhi-1.);
1006  } else {
1007  localPhiCenter = 0.25 * M_PI * ((double)mmHits.at(iHit).stationPhi-9.);
1008  }
1009  if (mmHits.at(iHit).stationName == 55){
1010  localPhiCenter += M_PI/8.; // small MM sectors
1011  if (mmHits.at(iHit).stationPhi == 5) localPhiCenter -= 2 * M_PI;
1012  }
1013  }
1014  if (mmHits.at(iHit).layerNumber >1 && mmHits.at(iHit).layerNumber < 6){
1015  double rInterpolate = slopefit * mmHits.at(iHit).z + interceptfit;
1016  double rProj = mmHits.at(iHit).r;
1017  if(std::abs(side_mm) < ZERO_LIMIT) {
1018  phiLocal.push_back(0);
1019  ATH_MSG_DEBUG("@@Merge@@ philocalmm 0");
1020  }
1021  else if ((mmHits.at(iHit).layerNumber)%2 == 0) { // layer U
1022  phiLocal.push_back( std::atan((rProj-rInterpolate)/tanTiltAngleU/rInterpolate));
1023  ATH_MSG_DEBUG("@@Merge@@ philocalmmU " << std::atan((rProj-rInterpolate)/tanTiltAngleU/rInterpolate));
1024  } else { // layer V
1025  phiLocal.push_back( std::atan((rProj-rInterpolate)/tanTiltAngleV/rInterpolate));
1026  ATH_MSG_DEBUG("@@Merge@@ philocalmmV " << std::atan((rProj-rInterpolate)/tanTiltAngleV/rInterpolate));
1027  }
1028  }
1029  }
1030  double sumPhiLocal = 0;
1031  for (unsigned int iHit = 0; iHit < phiLocal.size(); ++iHit ) {
1032  sumPhiLocal += phiLocal.at(iHit);
1033  }
1034  double phiLocalAvg = sumPhiLocal/phiLocal.size();
1035  ATH_MSG_DEBUG("@@Merge@@ philocalAvg " << phiLocalAvg);
1036 
1037  r.clear();
1038  z.clear();
1039  isStgc.clear();
1040  std::vector<double> r_stgc, z_stgc, r_mm, z_mm;
1041  std::vector<bool> isStgc_stgc, isStgc_mm;
1042  double side_stgc = 0;
1043  for(unsigned int iHit = 0; iHit < stgcHits.size(); ++iHit) {
1044  if (stgcHits.at(iHit).stationPhi<=5) localPhiCenter = 0.25 * M_PI * ((double)stgcHits.at(iHit).stationPhi-1.);
1045  if (stgcHits.at(iHit).stationPhi> 5) localPhiCenter = 0.25 * M_PI * ((double)stgcHits.at(iHit).stationPhi-9.);
1046  if (stgcHits.at(iHit).channelType == 1) {
1047  r_stgc.push_back(stgcHits.at(iHit).r/std::cos(phiLocalAvg));
1048  z_stgc.push_back(stgcHits.at(iHit).z);
1049  isStgc_stgc.push_back(true);
1050  side_stgc = std::abs(stgcHits.at(iHit).z)/stgcHits.at(iHit).z;
1051 
1052  ATH_MSG_DEBUG("@@Merge@@ stgc strip_r " << phiLocalAvg << " " << stgcHits.at(iHit).z << " " << stgcHits.at(iHit).r/std::cos(phiLocalAvg));
1053  }
1054 
1055  }
1056  double slopefit_stgc=0., interceptfit_stgc=99999., mse_stgc=1.e20;
1057  if(r_stgc.size() == 0) {
1058  ATH_MSG_DEBUG("No STGC hit to calculate superoint");
1059  } else {
1060  LinearFit(z_stgc,r_stgc,&slopefit_stgc,&interceptfit_stgc,&mse_stgc);
1061  }
1062 
1063  for(unsigned int iHit = 0; iHit < mmHits.size(); ++iHit) {
1064  if (mmHits.at(iHit).layerNumber < 2 || mmHits.at(iHit).layerNumber > 5) {
1065  r_mm.push_back(mmHits.at(iHit).r/std::cos(phiLocalAvg));
1066  z_mm.push_back(mmHits.at(iHit).z);
1067  isStgc_mm.push_back(false);
1068  } else {
1069  z_mm.push_back(mmHits.at(iHit).z);
1070  isStgc_mm.push_back(false);
1071 
1072  double rProj = mmHits.at(iHit).r;
1073  if(std::abs(side_mm) < ZERO_LIMIT) { // no layer U/V
1074  r_mm.push_back(rProj);
1075  }
1076  else if ((mmHits.at(iHit).layerNumber)%2 == 0) { // layer U
1077  double rPrime = (rProj * cosTiltAngleU)/(cos(phiLocalAvg)*cosTiltAngleU + sin(phiLocalAvg)*sinTiltAngleU);
1078  r_mm.push_back(rPrime);
1079  } else { //layer V
1080  double rPrime = (rProj * cosTiltAngleV)/(cos(phiLocalAvg)*cosTiltAngleV + sin(phiLocalAvg)*sinTiltAngleV);
1081  r_mm.push_back(rPrime);
1082  }
1083  }
1084  }
1085  double slopefit_mm=0., interceptfit_mm=99999., mse_mm=1.e20;
1086  if(r_mm.size() == 0) {
1087  ATH_MSG_WARNING("No MM hit to calculate superoint");
1088  } else {
1089  LinearFit(z_mm,r_mm,&slopefit_mm,&interceptfit_mm,&mse_mm);
1090  }
1091 
1092  unsigned int fmerge = 0;
1093  slopefit=0., interceptfit=99999., mse=1.e20;
1094  double side = 0;
1095  double StgcSegZ = 7526.329;
1096  double StgcSegR = 0;
1097  double MmSegZ = 7526.329;
1098  double MmSegR = 0;
1099  if (mse_stgc < 1.e7 && mse_mm < 1.e7) {
1100  r = r_stgc;
1101  copy(r_mm.begin(), r_mm.end(), back_inserter(r));
1102  z = z_stgc;
1103  copy(z_mm.begin(), z_mm.end(), back_inserter(z));
1104  isStgc = isStgc_stgc;
1105  copy(isStgc_mm.begin(), isStgc_mm.end(), back_inserter(isStgc));
1106 
1107  if(side_stgc < -1.*ZERO_LIMIT){
1108  StgcSegZ = -7526.329;
1109  }
1110  StgcSegR = slopefit_stgc * StgcSegZ + interceptfit_stgc;
1111  double StgcSegOriginTheta = std::atan(StgcSegR / StgcSegZ);
1112  double StgcSegEta = side_stgc * (- std::log(std::abs(std::tan(StgcSegOriginTheta / 2))));
1113  if(side_mm < -1.*ZERO_LIMIT){
1114  MmSegZ = -7526.329;
1115  }
1116  MmSegR = slopefit_mm * MmSegZ + interceptfit_mm;
1117  double MmSegOriginTheta = std::atan(MmSegR / MmSegZ);
1118  double MmSegEta = side_stgc * (- std::log(std::abs(std::tan(MmSegOriginTheta / 2))));
1119 
1120  double SegEtaAve = 0;
1121  if(std::abs(side_stgc) > ZERO_LIMIT || std::abs(side_mm) > ZERO_LIMIT){
1122  if(side_stgc*side_mm > 0){
1123  side = side_stgc;
1124  SegEtaAve = (StgcSegEta + MmSegEta)/2;
1125  } else if(std::abs(side_stgc) < ZERO_LIMIT) {
1126  side = side_mm;
1127  SegEtaAve = MmSegEta;
1128  } else if(std::abs(side_mm) < ZERO_LIMIT) {
1129  side = side_stgc;
1130  SegEtaAve = StgcSegEta;
1131  }
1132  }
1133 
1134  LinearFitWeight(z,r,isStgc,&slopefit,&interceptfit,&mse,SegEtaAve);
1135  fmerge = 1;
1136  }
1137  if (mse > 1.e7) {
1138  if (mse_stgc < mse_mm) {
1139  slopefit = slopefit_stgc;
1140  interceptfit = interceptfit_stgc;
1141  mse = mse_stgc;
1142  z = z_stgc;
1143  side = side_stgc;
1144  fmerge = 2;
1145  } else {
1146  slopefit = slopefit_mm;
1147  interceptfit = interceptfit_mm;
1148  mse = mse_mm;
1149  z = z_mm;
1150  side = side_mm;
1151  fmerge = 3;
1152  }
1153  }
1154  if (mse > 1.e19) {
1155  ATH_MSG_WARNING("No sTGC and MM hit to calculate superoint");
1156  return StatusCode::SUCCESS;
1157  }
1158 
1159  // store superpoint info in TrackData
1160  xAOD::L2MuonParameters::Chamber inner = xAOD::L2MuonParameters::Chamber::EndcapInner;
1161  TrigL2MuonSA::SuperPoint* superPoint = &(trackPattern.superPoints[inner]);
1162  double NSWCenterZ = 7526.329;
1163  if(side < -1.*ZERO_LIMIT){
1164  NSWCenterZ = -7526.329;
1165  }
1166  superPoint->R = slopefit * NSWCenterZ + interceptfit;
1167  superPoint->Phim = phiLocalAvg+localPhiCenter;
1168  superPoint->Z = NSWCenterZ;
1169  superPoint->Npoint = z.size();
1170 
1171  if (NSWCenterZ != 0) superPoint->Alin = slopefit;
1172  superPoint->Blin = interceptfit;
1173 
1174  ATH_MSG_DEBUG("Nsw Super Point r/phi/z/slope = "<<superPoint->R<<"/"<<superPoint->Phim<<"/"<<superPoint->Z<<"/"<<superPoint->Alin);
1175 
1176  ATH_MSG_DEBUG("@@Merge@@ Nsw Super Point r/phi/z/slope = "<<superPoint->R<<"/"<<superPoint->Phim<<"/"<<superPoint->Z<<"/"<<superPoint->Alin);
1177  ATH_MSG_DEBUG("@@Merge@@ fit slope= " << slopefit << " " << slopefit_stgc << " " << slopefit_mm);
1178  ATH_MSG_DEBUG("@@Merge@@ fit intercept= " << interceptfit << " " << interceptfit_stgc << " " << interceptfit_mm);
1179  ATH_MSG_DEBUG("@@Merge@@ fit mse= " << mse << " " << mse_stgc << " " << mse_mm);
1180  ATH_MSG_DEBUG("@@Merge@@ fit tech= " << fmerge);
1181 
1182 
1183  return StatusCode::SUCCESS;
1184 
1185 }

◆ calcWeightedSumHit()

StatusCode TrigL2MuonSA::NswStationFitter::calcWeightedSumHit ( TrigL2MuonSA::TrackPattern &  trackPattern) const

Definition at line 184 of file NswStationFitter.cxx.

185 {
186 
187  TrigL2MuonSA::StgcHits stgcHits = trackPattern.stgcSegment;
188  TrigL2MuonSA::MmHits mmHits = trackPattern.mmSegment;
189 
190  if(stgcHits.size()==0 && mmHits.size()==0)
191  return StatusCode::SUCCESS;
192 
193  //superpoint information
194  float rWeightedCenter=0.,zWeightedCenter=0.,phiWeightedCenter=0.;
195  int totNumRWeightedCenter=0.,totNumZWeightedCenter=0.,totNumPhiWeightedCenter=0.;
196  float localPhiCenter=0., localPhi;
197 
198  // loop over sTGC digits.
199  unsigned int iHit;
200  if (stgcHits.size()>0) {
201  for (iHit = 0; iHit < stgcHits.size(); iHit++){
202 
203  // Get the digit point.
204  TrigL2MuonSA::StgcHitData& hit = stgcHits[iHit];
205 
206  if (iHit==0 && hit.stationPhi<=5) localPhiCenter = 0.25 * M_PI * ((float)hit.stationPhi-1.);
207  if (iHit==0 && hit.stationPhi> 5) localPhiCenter = 0.25 * M_PI * ((float)hit.stationPhi-9.);
208  if (hit.stationName == 57){
209  localPhiCenter += M_PI/8.; // small sTGC sectors
210  if (hit.stationPhi == 5) localPhiCenter -= 2 * M_PI;
211  }
212 
213  localPhi = hit.phi - localPhiCenter;
214  if (localPhi > M_PI) localPhi -= 2.0*M_PI;
215  if (localPhi < -1.*M_PI) localPhi += 2.0*M_PI;
216  // calculate Weighted Center
217  if (hit.channelType == 1) { // strip
218  rWeightedCenter += hit.r;
219  totNumRWeightedCenter += 1;
220  }
221  if (hit.channelType == 2) { // wire
222  phiWeightedCenter += localPhi;
223  totNumPhiWeightedCenter += 1;
224  }
225  zWeightedCenter += hit.z;
226  totNumZWeightedCenter += 1;
227  }
228  }
229 
230  // loop over MM digits.
231  if (mmHits.size()>0) {
232  for (iHit = 0; iHit < mmHits.size(); iHit++){
233 
234  // Get the digit point.
235  TrigL2MuonSA::MmHitData& hit = mmHits[iHit];
236 
237  if (hit.stationName == 55){
238  localPhiCenter += M_PI/8.; // small MM sectors
239  if (hit.stationPhi == 5) localPhiCenter -= 2 * M_PI;
240  }
241 
242  localPhi = hit.phi - localPhiCenter;
243  if (localPhi > M_PI) localPhi -= 2.0*M_PI;
244  if (localPhi < -1.*M_PI) localPhi += 2.0*M_PI;
245  }
246  }
247 
248  // calculate Weighted Center
249  if (totNumRWeightedCenter!=0) rWeightedCenter /= totNumRWeightedCenter;
250  if (totNumPhiWeightedCenter!=0) phiWeightedCenter /= totNumPhiWeightedCenter;
251  if (totNumZWeightedCenter!=0) zWeightedCenter /= totNumZWeightedCenter;
252 
253  // store superpoint info in TrackData
254  xAOD::L2MuonParameters::Chamber inner = xAOD::L2MuonParameters::Chamber::EndcapInner;
255  TrigL2MuonSA::SuperPoint* superPoint = &(trackPattern.superPoints[inner]);
256  superPoint->R = rWeightedCenter;
257  superPoint->Phim = phiWeightedCenter+localPhiCenter;
258  superPoint->Z = zWeightedCenter;
259  superPoint->Npoint = totNumZWeightedCenter;
260  if (zWeightedCenter!=0) superPoint->Alin = rWeightedCenter/zWeightedCenter;
261  superPoint->Blin = 0.;
262 
263  ATH_MSG_DEBUG("Nsw Super Point r/phi/z/slope = "<<superPoint->R<<"/"<<superPoint->Phim<<"/"<<superPoint->Z<<"/"<<superPoint->Alin);
264 
265  return StatusCode::SUCCESS;
266 
267 }

◆ declareGaudiProperty() [1/4]

Gaudi::Details::PropertyBase& AthCommonDataStore< AthCommonMsg< AlgTool > >::declareGaudiProperty ( Gaudi::Property< T > &  hndl,
const SG::VarHandleKeyArrayType &   
)
inlineprivateinherited

specialization for handling Gaudi::Property<SG::VarHandleKeyArray>

Definition at line 170 of file AthCommonDataStore.h.

172  {
173  return *AthCommonDataStore<PBASE>::declareProperty(hndl.name(),
174  hndl.value(),
175  hndl.documentation());
176 
177  }

◆ declareGaudiProperty() [2/4]

Gaudi::Details::PropertyBase& AthCommonDataStore< AthCommonMsg< AlgTool > >::declareGaudiProperty ( Gaudi::Property< T > &  hndl,
const SG::VarHandleKeyType &   
)
inlineprivateinherited

specialization for handling Gaudi::Property<SG::VarHandleKey>

Definition at line 156 of file AthCommonDataStore.h.

158  {
159  return *AthCommonDataStore<PBASE>::declareProperty(hndl.name(),
160  hndl.value(),
161  hndl.documentation());
162 
163  }

◆ declareGaudiProperty() [3/4]

Gaudi::Details::PropertyBase& AthCommonDataStore< AthCommonMsg< AlgTool > >::declareGaudiProperty ( Gaudi::Property< T > &  hndl,
const SG::VarHandleType &   
)
inlineprivateinherited

specialization for handling Gaudi::Property<SG::VarHandleBase>

Definition at line 184 of file AthCommonDataStore.h.

186  {
187  return *AthCommonDataStore<PBASE>::declareProperty(hndl.name(),
188  hndl.value(),
189  hndl.documentation());
190  }

◆ declareGaudiProperty() [4/4]

Gaudi::Details::PropertyBase& AthCommonDataStore< AthCommonMsg< AlgTool > >::declareGaudiProperty ( Gaudi::Property< T > &  t,
const SG::NotHandleType &   
)
inlineprivateinherited

specialization for handling everything that's not a Gaudi::Property<SG::VarHandleKey> or a <SG::VarHandleKeyArray>

Definition at line 199 of file AthCommonDataStore.h.

200  {
201  return PBASE::declareProperty(t);
202  }

◆ declareProperty() [1/6]

Gaudi::Details::PropertyBase* AthCommonDataStore< AthCommonMsg< AlgTool > >::declareProperty ( const std::string &  name,
SG::VarHandleBase &  hndl,
const std::string &  doc,
const SG::VarHandleType &   
)
inlineinherited

Declare a new Gaudi property.

Parameters
nameName of the property.
hndlObject holding the property value.
docDocumentation string for the property.

This is the version for types that derive from SG::VarHandleBase. The property value object is put on the input and output lists as appropriate; then we forward to the base class.

Definition at line 245 of file AthCommonDataStore.h.

249  {
250  this->declare(hndl.vhKey());
251  hndl.vhKey().setOwner(this);
252 
253  return PBASE::declareProperty(name,hndl,doc);
254  }

◆ declareProperty() [2/6]

Gaudi::Details::PropertyBase* AthCommonDataStore< AthCommonMsg< AlgTool > >::declareProperty ( const std::string &  name,
SG::VarHandleKey &  hndl,
const std::string &  doc,
const SG::VarHandleKeyType &   
)
inlineinherited

Declare a new Gaudi property.

Parameters
nameName of the property.
hndlObject holding the property value.
docDocumentation string for the property.

This is the version for types that derive from SG::VarHandleKey. The property value object is put on the input and output lists as appropriate; then we forward to the base class.

Definition at line 221 of file AthCommonDataStore.h.

225  {
226  this->declare(hndl);
227  hndl.setOwner(this);
228 
229  return PBASE::declareProperty(name,hndl,doc);
230  }

◆ declareProperty() [3/6]

Gaudi::Details::PropertyBase* AthCommonDataStore< AthCommonMsg< AlgTool > >::declareProperty ( const std::string &  name,
SG::VarHandleKeyArray &  hndArr,
const std::string &  doc,
const SG::VarHandleKeyArrayType &   
)
inlineinherited

Definition at line 259 of file AthCommonDataStore.h.

263  {
264 
265  // std::ostringstream ost;
266  // ost << Algorithm::name() << " VHKA declareProp: " << name
267  // << " size: " << hndArr.keys().size()
268  // << " mode: " << hndArr.mode()
269  // << " vhka size: " << m_vhka.size()
270  // << "\n";
271  // debug() << ost.str() << endmsg;
272 
273  hndArr.setOwner(this);
274  m_vhka.push_back(&hndArr);
275 
276  Gaudi::Details::PropertyBase* p = PBASE::declareProperty(name, hndArr, doc);
277  if (p != 0) {
278  p->declareUpdateHandler(&AthCommonDataStore<PBASE>::updateVHKA, this);
279  } else {
280  ATH_MSG_ERROR("unable to call declareProperty on VarHandleKeyArray "
281  << name);
282  }
283 
284  return p;
285 
286  }

◆ declareProperty() [4/6]

Gaudi::Details::PropertyBase* AthCommonDataStore< AthCommonMsg< AlgTool > >::declareProperty ( const std::string &  name,
T &  property,
const std::string &  doc,
const SG::NotHandleType &   
)
inlineinherited

Declare a new Gaudi property.

Parameters
nameName of the property.
propertyObject holding the property value.
docDocumentation string for the property.

This is the generic version, for types that do not derive from SG::VarHandleKey. It just forwards to the base class version of declareProperty.

Definition at line 333 of file AthCommonDataStore.h.

337  {
338  return PBASE::declareProperty(name, property, doc);
339  }

◆ declareProperty() [5/6]

Gaudi::Details::PropertyBase* AthCommonDataStore< AthCommonMsg< AlgTool > >::declareProperty ( const std::string &  name,
T &  property,
const std::string &  doc = "none" 
)
inlineinherited

Declare a new Gaudi property.

Parameters
nameName of the property.
propertyObject holding the property value.
docDocumentation string for the property.

This dispatches to either the generic declareProperty or the one for VarHandle/Key/KeyArray.

Definition at line 352 of file AthCommonDataStore.h.

355  {
356  typedef typename SG::HandleClassifier<T>::type htype;
357  return declareProperty (name, property, doc, htype());
358  }

◆ declareProperty() [6/6]

Gaudi::Details::PropertyBase& AthCommonDataStore< AthCommonMsg< AlgTool > >::declareProperty ( Gaudi::Property< T > &  t)
inlineinherited

Definition at line 145 of file AthCommonDataStore.h.

145  {
146  typedef typename SG::HandleClassifier<T>::type htype;
147  return AthCommonDataStore<PBASE>::declareGaudiProperty(t, htype());
148  }

◆ detStore()

const ServiceHandle<StoreGateSvc>& AthCommonDataStore< AthCommonMsg< AlgTool > >::detStore ( ) const
inlineinherited

The standard StoreGateSvc/DetectorStore Returns (kind of) a pointer to the StoreGateSvc.

Definition at line 95 of file AthCommonDataStore.h.

95 { return m_detStore; }

◆ evtStore() [1/2]

ServiceHandle<StoreGateSvc>& AthCommonDataStore< AthCommonMsg< AlgTool > >::evtStore ( )
inlineinherited

The standard StoreGateSvc (event store) Returns (kind of) a pointer to the StoreGateSvc.

Definition at line 85 of file AthCommonDataStore.h.

85 { return m_evtStore; }

◆ evtStore() [2/2]

const ServiceHandle<StoreGateSvc>& AthCommonDataStore< AthCommonMsg< AlgTool > >::evtStore ( ) const
inlineinherited

The standard StoreGateSvc (event store) Returns (kind of) a pointer to the StoreGateSvc.

Definition at line 90 of file AthCommonDataStore.h.

90 { return m_evtStore; }

◆ extraDeps_update_handler()

void AthCommonDataStore< AthCommonMsg< AlgTool > >::extraDeps_update_handler ( Gaudi::Details::PropertyBase &  ExtraDeps)
protectedinherited

Add StoreName to extra input/output deps as needed.

use the logic of the VarHandleKey to parse the DataObjID keys supplied via the ExtraInputs and ExtraOuputs Properties to add the StoreName if it's not explicitly given

◆ findMmHitsInSegment()

StatusCode TrigL2MuonSA::NswStationFitter::findMmHitsInSegment ( TrigL2MuonSA::MmHits &  mmHits) const

Definition at line 1187 of file NswStationFitter.cxx.

1188 {
1189  if(mmHits.size() == 0) return StatusCode::SUCCESS;
1190  int hitsInRoad = 0;
1191  for(unsigned int iHit = 0; iHit < mmHits.size(); iHit++){
1192  if(mmHits.at(iHit).isOutlier == 0){
1193  hitsInRoad++;
1194  mmHits.at(iHit).isOutlier = 1;
1195  }
1196  }
1197  if(hitsInRoad == 0) return StatusCode::SUCCESS;
1198 
1199  if(hitsInRoad < 6) {
1200  ATH_MSG_DEBUG("Number of MM hits is too small, at least 6 hits required : "<<hitsInRoad<<" hits");
1201  return StatusCode::SUCCESS;
1202  } else if(hitsInRoad > 100) {
1203  ATH_MSG_WARNING("Number of MM hits is too large, at most (2^16 - 1) hits allowed : "<<hitsInRoad<<" hits");
1204  return StatusCode::SUCCESS;
1205  }
1206 
1207  std::array< std::vector<int>, 8 > hitIdByLayer;
1208  for(unsigned int iHit = 0; iHit < mmHits.size(); ++iHit){
1209  if(mmHits.at(iHit).isOutlier != 1) continue;
1210  int layerNumber = mmHits.at(iHit).layerNumber;
1211  if (layerNumber > 7) {
1212  ATH_MSG_WARNING("MM hit layer number > 7");
1213  continue;
1214  }
1215  hitIdByLayer[layerNumber].push_back(iHit);
1216  }
1217  ATH_MSG_DEBUG("@@MM@@ Nhits " << hitIdByLayer[0].size()
1218  << " " << hitIdByLayer[1].size()
1219  << " " << hitIdByLayer[2].size()
1220  << " " << hitIdByLayer[3].size()
1221  << " " << hitIdByLayer[4].size()
1222  << " " << hitIdByLayer[5].size()
1223  << " " << hitIdByLayer[6].size()
1224  << " " << hitIdByLayer[7].size());
1225 
1226  std::vector< std::array<int, 8> > mmHitIds;
1227  findSetOfMmHitIds(mmHits, hitIdByLayer, mmHitIds);
1228  for (unsigned int iHit = 0; iHit < mmHitIds.size(); ++iHit) {
1229  std::array<int, 8> hitIds = mmHitIds.at(iHit);
1230  for (unsigned int iLayer = 0; iLayer < 8; ++iLayer) {
1231  if (hitIds[iLayer] != -1) {
1232  mmHits.at(hitIds[iLayer]).isOutlier = 0;
1233  }
1234  }
1235  }
1236  return StatusCode::SUCCESS;
1237 }

◆ findSetOfMmHitIds()

void TrigL2MuonSA::NswStationFitter::findSetOfMmHitIds ( TrigL2MuonSA::MmHits &  mmHits,
std::array< std::vector< int >, 8 >  hitIdByLayer,
std::vector< std::array< int, 8 >> &  hitIdsCandidate 
) const

Definition at line 1239 of file NswStationFitter.cxx.

1242 {
1243 
1244  double NSWCenterZ = 7526.329;
1245  int side = 0;
1246  for (unsigned int iLayer = 0; iLayer < 8; ++iLayer) {
1247  if ( hitIdByLayer[iLayer].size() > 0) {
1248  side = std::abs(mmHits.at(hitIdByLayer[iLayer].at(0)).z)/mmHits.at(hitIdByLayer[iLayer].at(0)).z;
1249  break;
1250  }
1251  }
1252  NSWCenterZ = NSWCenterZ * side;
1253 
1254  std::array<std::vector<unsigned long int>,4> hitIdsInTwo;
1255  std::array<std::vector<double>,4> slopeInTwo;
1256  std::array<std::vector<double>,4> interceptInTwo;
1257  // Loop over pairs of the i-th and the (i+6)-th layers
1258  for(unsigned int iPair = 0; iPair < 4; ++iPair){
1259 
1260  unsigned int nHitsInInner = hitIdByLayer[iPair].size();
1261  unsigned int nHitsInOuter;
1262  if (iPair < 2) { // paris of X layers
1263  nHitsInOuter = hitIdByLayer[iPair+6].size();
1264  } else { // pairs of U or V layers
1265  nHitsInOuter = hitIdByLayer[iPair+2].size();
1266  }
1267 
1268  if ( nHitsInInner > 0xffff-1 || nHitsInOuter > 0xffff-1) {
1269  ATH_MSG_WARNING("Number of Mm hits in layers exceeds the limit of (2^16 - 1) : Number of Mm hits in "<<iPair<<"th layer = "<< nHitsInInner
1270  <<", Number of Mm hits in "<<iPair+4<<"th layer = "<<nHitsInOuter);
1271  ATH_MSG_WARNING("Number of Mm hits is limitted to (2^16 - 1) and hits with id more than (2^16 -1) will be trancated.");
1272  if (nHitsInInner > 0xffff-1) {nHitsInInner = 0xffff-1;}
1273  if (nHitsInOuter > 0xffff-1) {nHitsInOuter = 0xffff-1;}
1274  }
1275 
1276  bool foundCounterparts[0xffff] = {};
1277  // Loop over hits in the i-th layer
1278  for(unsigned int iHit = 0; iHit < nHitsInInner; ++iHit){
1279 
1280  bool foundCounterpart = 0;
1281 
1282  double z[2] = {};
1283  double r[2] = {};
1284 
1285  int iHitId = hitIdByLayer[iPair].at(iHit);
1286  r[0] = mmHits.at(iHitId).r;
1287  z[0] = mmHits.at(iHitId).z;
1288 
1289  // Loop over hits in the (i+6)-th layer
1290  for(unsigned int jHit = 0; jHit < nHitsInOuter; ++jHit){
1291 
1292  int jHitId;
1293  if (iPair < 2) {
1294  jHitId = hitIdByLayer[iPair+6].at(jHit);
1295  } else {
1296  jHitId = hitIdByLayer[iPair+2].at(jHit);
1297  }
1298  r[1] = mmHits.at(jHitId).r;
1299  z[1] = mmHits.at(jHitId).z;
1300 
1301  double slope = (r[1] - r[0]) / (z[1] - z[0]);
1302  double intercept = slope*(0. - z[0]) + r[0];
1303  // select pairs whose slops in limited regions
1304  if(std::abs(slope) < 0.1 || std::abs(slope) > 0.7 || std::abs(intercept) > 500.) continue;
1305 
1306  int encodedIds = (iHitId<<16) + jHitId;
1307  hitIdsInTwo[iPair].push_back(encodedIds);
1308  slopeInTwo[iPair].push_back(slope);
1309  interceptInTwo[iPair].push_back(intercept);
1310 
1311  foundCounterpart = 1;
1312  foundCounterparts[jHit] = 1;
1313  }//end of jHit in the (i+6)-th layer
1314  if(!foundCounterpart){ // in case of no counterpart in the (i+4)-th layer
1315  int encodedIds = (iHitId<<16) + 0xffff; // fill all bits with 1 for hit id for the layer with no hit
1316  hitIdsInTwo[iPair].push_back(encodedIds);
1317  slopeInTwo[iPair].push_back(r[0]/z[0]);
1318  interceptInTwo[iPair].push_back(0.);
1319  }
1320  }//end of iHit in the i-th layer
1321  // Loop over hits in the (i+4)-th layer
1322  for(unsigned int jHit = 0; jHit < nHitsInOuter; ++jHit){
1323  if (!foundCounterparts[jHit]) {
1324  int jHitId;
1325  if (iPair < 2) {
1326  jHitId = hitIdByLayer[iPair+6].at(jHit);
1327  } else {
1328  jHitId = hitIdByLayer[iPair+2].at(jHit);
1329  }
1330  int encodedIds = (0xFFFFu<<16) + jHitId; // fill all bits with 1 for hit id for the layer with no hit
1331  hitIdsInTwo[iPair].push_back(encodedIds);
1332  slopeInTwo[iPair].push_back(mmHits.at(jHitId).r/mmHits.at(jHitId).z);
1333  interceptInTwo[iPair].push_back(0.);
1334  }
1335  }
1336  }//end of pair loop
1337  ATH_MSG_DEBUG("@@MM@@ Npairs " << hitIdsInTwo[0].size() << " " << hitIdsInTwo[1].size() << " " << hitIdsInTwo[2].size() << " " << hitIdsInTwo[3].size());
1338  for (unsigned int iLayer = 0; iLayer < 4; ++ iLayer) {
1339  for (unsigned int iPair = 0; iPair < slopeInTwo[iLayer].size(); ++iPair) {
1340  ATH_MSG_DEBUG("@@MM@@ pair fit slope= " << slopeInTwo[iLayer].at(iPair) << " intercept= " << interceptInTwo[iLayer].at(iPair));
1341  }
1342  }
1343 
1344  std::vector<std::array<int, 4>> hitIdsInFourX;
1345  std::vector<double> slopeInFourX;
1346  std::vector<double> interceptInFourX;
1347  std::vector<double> mseInFourX;
1348 
1349  unsigned int nPairsInInnerX = hitIdsInTwo[0].size();
1350  unsigned int nPairsInOuterX = hitIdsInTwo[1].size();
1351 
1352  for(unsigned int iPairX = 0; iPairX < nPairsInInnerX; ++iPairX){
1353 
1354  double slope[2];
1355  double intercept[2];
1356  double spR[2];
1357 
1358  slope[0] = slopeInTwo[0].at(iPairX);
1359  intercept[0] = interceptInTwo[0].at(iPairX);
1360  spR[0] = slope[0] * NSWCenterZ + intercept[0];
1361  for(unsigned int jPairX = 0; jPairX < nPairsInOuterX; ++jPairX){
1362  int ihitIds = hitIdsInTwo[0].at(iPairX);
1363  int jhitIds = hitIdsInTwo[1].at(jPairX);
1364  if ( ((ihitIds>>16 & 0xffff) == 0xffff || (ihitIds & 0xffff) == 0xffff) &&
1365  ((jhitIds>>16 & 0xffff) == 0xffff || (jhitIds & 0xffff) == 0xffff )) continue; // require at least 3 hits in 4 layers
1366 
1367  slope[1] = slopeInTwo[1].at(jPairX);
1368  intercept[1] = interceptInTwo[1].at(jPairX);
1369  spR[1] = slope[1] * NSWCenterZ + intercept[1];
1370 
1371  if(std::abs(spR[1] - spR[0]) > 50. ||
1372  std::abs((intercept[1] + intercept[0]) / 2) > 200.) continue;
1373 
1374  std::array<int, 4> setOfHitIds{};
1375  setOfHitIds[0] = (ihitIds>>16 & 0xffff);
1376  setOfHitIds[1] = (ihitIds & 0xffff);
1377  setOfHitIds[2] = (jhitIds>>16 & 0xffff);
1378  setOfHitIds[3] = (jhitIds & 0xffff);
1379  std::vector<double> r;
1380  std::vector<double> z;
1381  for(unsigned int iLayer = 0; iLayer < 4; ++iLayer){
1382  if(setOfHitIds[iLayer] == 0xffff) {
1383  continue;
1384  }
1385  double rhit = mmHits.at(setOfHitIds[iLayer]).r;
1386  double zhit = mmHits.at(setOfHitIds[iLayer]).z;
1387  r.push_back(rhit);
1388  z.push_back(zhit);
1389  }
1390  double slopefit=0., interceptfit=99999., mse=-1.;
1391  LinearFit(z,r,&slopefit, &interceptfit, &mse);
1392 
1393  hitIdsInFourX.push_back(setOfHitIds);
1394  slopeInFourX.push_back(slopefit);
1395  interceptInFourX.push_back(interceptfit);
1396  mseInFourX.push_back(mse);
1397 
1398  }// end of iPair of the i-th and (i+4)-th layers, i=0,2
1399  }// end of jPair of the j-th and (j+4)-th layers, j=1,3
1400  ATH_MSG_DEBUG("@@MM@@ X Nquads " << hitIdsInFourX.size());
1401  for (unsigned int iQuad = 0; iQuad < slopeInFourX.size(); ++iQuad) {
1402  ATH_MSG_DEBUG("@@MM@@ X quad fit slope= " << slopeInFourX.at(iQuad) << " intercept= " << interceptInFourX.at(iQuad) << " mse= " << mseInFourX.at(iQuad));
1403  }
1404 
1405  if(!hitIdsInFourX.size()){
1406  ATH_MSG_WARNING("No candidate segment found in MM X layers");
1407  return;
1408  }
1409 
1410  double tanTiltAngleU = 0,
1411  tanTiltAngleV = 0;
1412  if(side > ZERO_LIMIT){
1413  tanTiltAngleU = tan( 1.5/360.*2.*M_PI),
1414  tanTiltAngleV = tan(-1.5/360.*2.*M_PI);
1415  } else if(side < -1.*ZERO_LIMIT){
1416  tanTiltAngleU = tan(-1.5/360.*2.*M_PI),
1417  tanTiltAngleV = tan(1.5/360.*2.*M_PI);
1418  }
1419 
1420  std::vector< std::array<int, 8> > hitIdsInEight;
1421  std::vector<double> mseInEight;
1422 
1423  for(unsigned int iQuadX = 0; iQuadX < hitIdsInFourX.size(); ++iQuadX){
1424  if(mseInFourX.at(iQuadX) > 10) continue;
1425 
1426  double slopeX = slopeInFourX.at(iQuadX);
1427  double interceptX = interceptInFourX.at(iQuadX);
1428  std::array<int,4> hitIdsX{};
1429  hitIdsX = hitIdsInFourX.at(iQuadX);
1430 
1431  for (unsigned int iPairU = 0; iPairU < hitIdsInTwo[2].size(); ++iPairU) {
1432 
1433  int hitIdsU[2];
1434  hitIdsU[0] = hitIdsInTwo[2].at(iPairU)>>16 & 0xffff;
1435  hitIdsU[1] = hitIdsInTwo[2].at(iPairU) & 0xffff;
1436  double phiLocalU[2]={-99999,-99999};
1437  for(unsigned int iLayer = 0; iLayer < 2; ++iLayer) {
1438  if (hitIdsU[iLayer] == 0xffff) continue;
1439  if (hitIdsU[iLayer] < 0) {
1440  ATH_MSG_DEBUG("@@MM@@ hitIdsU[iLayer] iLayer= " << iLayer << " hitIdsU[iLayer]= " << hitIdsU[iLayer]);
1441  }
1442  double rInterpolate = slopeX * mmHits.at(hitIdsU[iLayer]).z + interceptX;
1443  double rProj = mmHits.at(hitIdsU[iLayer]).r;
1444  if(std::abs(tanTiltAngleU) < ZERO_LIMIT)
1445  phiLocalU[iLayer] = 0;
1446  else
1447  phiLocalU[iLayer] = std::atan((rProj-rInterpolate)/tanTiltAngleU/rInterpolate);
1448  }
1449 
1450  for(unsigned int iPairV = 0; iPairV < hitIdsInTwo[3].size(); ++iPairV) {
1451 
1452  int hitIdsV[2];
1453  hitIdsV[0] = hitIdsInTwo[3].at(iPairV)>>16 & 0xffff;
1454  hitIdsV[1] = hitIdsInTwo[3].at(iPairV) & 0xffff;
1455 
1456  if( (hitIdsU[0] == 0xffff || hitIdsU[1] == 0xffff) &&
1457  (hitIdsV[0] == 0xffff || hitIdsV[1] == 0xffff) ) continue; // require at least 3 UV layers having hits
1458 
1459  double phiLocalV[2]={-99999,-99999};
1460  for(unsigned int iLayer = 0; iLayer < 2; ++iLayer) {
1461  if (hitIdsV[iLayer] == 0xffff) continue;
1462  if (hitIdsV[iLayer] < 0) {
1463  ATH_MSG_DEBUG("@@MM@@ hitIdsV[iLayer] iLayer= " << iLayer << " hitIdsV[iLayer]= " << hitIdsV[iLayer]);
1464  }
1465  double rInterpolate = slopeX * mmHits.at(hitIdsV[iLayer]).z + interceptX;
1466  double rProj = mmHits.at(hitIdsV[iLayer]).r;
1467  if(std::abs(tanTiltAngleV) < ZERO_LIMIT)
1468  phiLocalV[iLayer] = 0;
1469  else
1470  phiLocalV[iLayer] = std::atan((rProj-rInterpolate)/tanTiltAngleV/rInterpolate);
1471  }
1472 
1473  if ( std::abs(phiLocalU[0]-phiLocalV[0]) > 0.05 &&
1474  std::abs(phiLocalU[1]-phiLocalV[1]) > 0.05) continue;
1475 
1476  // average of phis in 4 UV layers
1477  double phiLocalUV = 0;
1478  int nPhi = 0;
1479  for(unsigned int iLayer = 0; iLayer < 2; ++iLayer) {
1480  if(phiLocalU[iLayer] > -99999.) {
1481  phiLocalUV += phiLocalU[iLayer];
1482  ++nPhi;
1483  }
1484  if(phiLocalV[iLayer] > -99999.) {
1485  phiLocalUV += phiLocalV[iLayer];
1486  ++nPhi;
1487  }
1488  }
1489  phiLocalUV /= nPhi;
1490 
1491  std::array<int, 8> setOfHitIds = {-1,-1,-1,-1,-1,-1,-1,-1};
1492  std::vector<double> r, z;
1493  for (unsigned int iLayer = 0; iLayer < 4; ++iLayer) {
1494  if ( hitIdsX[iLayer] != 0xffff) {
1495  if (hitIdsX[iLayer] < 0) {
1496  ATH_MSG_DEBUG("@@MM@@ hitIdsX[iLayer] iLayer= " << iLayer << " hitIdsX[iLayer]= " << hitIdsX[iLayer]);
1497  }
1498  z.push_back(mmHits.at(hitIdsX[iLayer]).z);
1499  r.push_back(mmHits.at(hitIdsX[iLayer]).r / std::cos(phiLocalUV));
1500  setOfHitIds[iLayer] = hitIdsX[iLayer];
1501  }
1502  }
1503  for (unsigned int iLayer = 0; iLayer < 2; ++iLayer) {
1504  if ( hitIdsU[iLayer] != 0xffff) {
1505  if (hitIdsU[iLayer] < 0) {
1506  ATH_MSG_DEBUG("@@MM@@ 2 hitIdsU[iLayer] iLayer= " << iLayer << " hitIdsU[iLayer]= " << hitIdsU[iLayer]);
1507  }
1508  z.push_back(mmHits.at(hitIdsU[iLayer]).z);
1509  r.push_back(mmHits.at(hitIdsU[iLayer]).r*(std::cos(phiLocalUV) + 1/std::cos(phiLocalUV))/2.);
1510  setOfHitIds[iLayer+4] = hitIdsU[iLayer];
1511  }
1512  if ( hitIdsV[iLayer] != 0xffff) {
1513  if (hitIdsV[iLayer] < 0) {
1514  ATH_MSG_DEBUG("@@MM@@ 2 hitIdsV[iLayer] iLayer= " << iLayer << " hitIdsV[iLayer]= " << hitIdsV[iLayer]);
1515  }
1516  z.push_back(mmHits.at(hitIdsV[iLayer]).z);
1517  r.push_back(mmHits.at(hitIdsV[iLayer]).r*(std::cos(phiLocalUV) + 1/std::cos(phiLocalUV))/2.);
1518  setOfHitIds[iLayer+6] = hitIdsV[iLayer];
1519  }
1520  }
1521  double slopefit=0., interceptfit=99999., mse=-1.;
1522  LinearFit(z,r,&slopefit,&interceptfit,&mse);
1523 
1524  hitIdsInEight.push_back(setOfHitIds);
1525  mseInEight.push_back(mse);
1526  } // end of Pair loop of V layers
1527  } // end of Pair loop of U Layers
1528  }// end of Quad loop of X layers
1529  ATH_MSG_DEBUG("@@MM@@ Noctets " << hitIdsInEight.size());
1530 
1531  std::vector<int> nOctetSegments;
1532  std::vector<int> patternStationName;
1533  for (unsigned int iOctet = 0; iOctet < hitIdsInEight.size(); ++iOctet) {
1534  bool isFirstHit = true;
1535  int hitStationName = 0;
1536  int nOctetSegment = 0;
1537  ATH_MSG_DEBUG("@@MM@@ octet fit mse " << mseInEight.at(iOctet));
1538  std::array<int, 8> tmpOctet = hitIdsInEight.at(iOctet);
1539  for (unsigned int iLayer = 0; iLayer < 8; ++iLayer) {
1540  if (tmpOctet[iLayer] != -1) {
1541 
1542  if(isFirstHit){
1543  hitStationName = mmHits.at(tmpOctet[iLayer]).stationName;
1544  isFirstHit = false;
1545 
1546  ATH_MSG_DEBUG("@@MM@@ octet pos r= " << mmHits.at(tmpOctet[iLayer]).r << " phi= " << mmHits.at(tmpOctet[iLayer]).phi << " z= " << mmHits.at(tmpOctet[iLayer]).z);
1547  nOctetSegment++;
1548  }
1549  else if(mmHits.at(tmpOctet[iLayer]).stationName == hitStationName){
1550  ATH_MSG_DEBUG("@@MM@@ octet pos r= " << mmHits.at(tmpOctet[iLayer]).r << " phi= " << mmHits.at(tmpOctet[iLayer]).phi << " z= " << mmHits.at(tmpOctet[iLayer]).z);
1551  nOctetSegment++;
1552  }
1553 
1554  }
1555  }
1556  nOctetSegments.push_back(nOctetSegment);
1557  patternStationName.push_back(hitStationName);
1558  }
1559 
1560  double mseminL = 100000.;
1561  double mseminS = 100000.;
1562  std::vector<int> octetIds(2,-1);
1563  for(unsigned int iOctet = 0; iOctet < hitIdsInEight.size(); ++iOctet){
1564  if(patternStationName.at(iOctet) == 56){
1565  if( mseInEight.at(iOctet) < mseminL) {
1566  mseminL = mseInEight.at(iOctet);
1567  }
1568  }
1569  else if(patternStationName.at(iOctet) == 55){
1570  if( mseInEight.at(iOctet) < mseminS) {
1571  mseminS = mseInEight.at(iOctet);
1572  }
1573  }
1574  }// end of Octet loop
1575 
1576  for(unsigned int iOctet = 0; iOctet < hitIdsInEight.size(); ++iOctet){
1577  if(patternStationName.at(iOctet) == 56){
1578  if( mseInEight.at(iOctet) != mseminL) {
1579  continue;
1580  }
1581  }
1582  else if(patternStationName.at(iOctet) == 55){
1583  if( mseInEight.at(iOctet) != mseminS) {
1584  continue;
1585  }
1586  }
1587  octetIds.push_back(iOctet);
1588  }
1589 
1590  for(unsigned int ids = 0; ids < octetIds.size(); ids++){
1591  if (octetIds.at(ids) != -1) {
1592  hitIdsCandidate.push_back(hitIdsInEight.at(octetIds.at(ids)));
1593  }
1594  }
1595 }

◆ findSetOfStgcHitIds()

void TrigL2MuonSA::NswStationFitter::findSetOfStgcHitIds ( TrigL2MuonSA::StgcHits &  stgcHits,
std::array< std::vector< int >, 8 >  hitIdByLayer,
std::vector< std::array< int, 8 >> &  hitIdsCandidate 
) const

Definition at line 361 of file NswStationFitter.cxx.

364 {
365  double NSWCenterZ = 7526.329;
366  int side = 0;
367 
368  bool isStrip = 0;
369  for (unsigned int iLayer = 0; iLayer < 8; ++iLayer) {
370  if ( hitIdByLayer[iLayer].size() > 0) {
371  side = std::abs(stgcHits.at(hitIdByLayer[iLayer].at(0)).z)/stgcHits.at(hitIdByLayer[iLayer].at(0)).z;
372  if ( stgcHits.at(hitIdByLayer[iLayer].at(0)).channelType == 1 ){
373  isStrip = 1;
374  break;
375  }
376  }
377  }
378  NSWCenterZ = NSWCenterZ * side;
379 
380  std::array<std::vector<unsigned long int>,4> hitIdsInTwo;
381  std::array<std::vector<double>,4> slopeInTwo;
382  std::array<std::vector<double>,4> interceptInTwo;
383 
384  // Loop over pairs of the i-th and the (i+4)-th layers
385  for(unsigned int iPair = 0; iPair < 4; ++iPair){
386  unsigned int nHitsInInner = hitIdByLayer[iPair].size();
387  unsigned int nHitsInOuter = hitIdByLayer[iPair+4].size();
388  if ( nHitsInInner > 0xffff-1 || nHitsInOuter > 0xffff-1) {
389  ATH_MSG_WARNING("Number of Stgc hits in layers exceeds the limit of (2^16 - 1) : Number of Stgc hits in "<<iPair<<"th layer = "<< nHitsInInner
390  <<", Number of Stgc hits in "<<iPair+4<<"th layer = "<<nHitsInOuter);
391  ATH_MSG_WARNING("Number of Stgc hits is limitted to (2^16 - 1) and hits with id more than (2^16 -1) will be trancated.");
392  if (nHitsInInner > 0xffff-1) {nHitsInInner = 0xffff-1;}
393  if (nHitsInOuter > 0xffff-1) {nHitsInOuter = 0xffff-1;}
394  }
395 
396  bool foundCounterparts[256] = {};
397  // Loop over hits in the i-th layer
398  for(unsigned int iHit = 0; iHit < nHitsInInner; ++iHit){
399  bool foundCounterpart = 0;
400 
401  double z[2] = {};
402  double r[2] = {};
403 
404  int iHitId = hitIdByLayer[iPair].at(iHit);
405  if(isStrip){
406  r[0] = stgcHits.at(iHitId).r;
407  z[0] = stgcHits.at(iHitId).z;
408  } else {
409  double localPhiCenter;
410  if (stgcHits.at(iHitId).stationPhi<=5) {
411  localPhiCenter = 0.25 * M_PI * ((double)stgcHits.at(iHitId).stationPhi-1.);
412  } else {
413  localPhiCenter = 0.25 * M_PI * ((double)stgcHits.at(iHitId).stationPhi-9.);
414  }
415 
416  if (stgcHits.at(iHitId).stationName == 57){
417  localPhiCenter += M_PI/8.; // small sTGC sectors
418  if (stgcHits.at(iHitId).stationPhi == 5) localPhiCenter -= 2. * M_PI;
419  }
420 
421  double phiProj = stgcHits.at(iHitId).phi - localPhiCenter;
422  if (phiProj > M_PI) phiProj -= 2.0*M_PI;
423  if (phiProj < -1.*M_PI) phiProj += 2.0*M_PI;
424  r[0] = stgcHits.at(iHitId).r;
425  z[0] = phiProj;
426  }
427 
428  // Loop over hits in the (i+4)-th layer
429  for(unsigned int jHit = 0; jHit < nHitsInOuter; ++jHit){
430  int jHitId = hitIdByLayer[iPair+4].at(jHit);
431  double slope, intercept;
432  if(isStrip) {
433  r[1] = stgcHits.at(jHitId).r;
434  z[1] = stgcHits.at(jHitId).z;
435  slope = (r[1] - r[0]) / (z[1] - z[0]);
436  intercept = slope*(0. - z[0]) + r[0];
437  // select pairs whose slops in limited regions
438  if(std::abs(slope) < 0.14 || std::abs(slope) > 0.6 || std::abs(intercept) > 300.) continue;
439  } else {
440  double localPhiCenter;
441  if (stgcHits.at(jHitId).stationPhi<=5) {
442  localPhiCenter = 0.25 * M_PI * ((double)stgcHits.at(jHitId).stationPhi-1.);
443  } else {
444  localPhiCenter = 0.25 * M_PI * ((double)stgcHits.at(jHitId).stationPhi-9.);
445  }
446 
447  if (stgcHits.at(jHitId).stationName == 57){
448  localPhiCenter += M_PI/8.; // small sTGC sectors
449  if (stgcHits.at(jHitId).stationPhi == 5) localPhiCenter -= 2 * M_PI;
450  }
451 
452  double phiProj = stgcHits.at(jHitId).phi - localPhiCenter;
453  if (phiProj > M_PI) phiProj -= 2.0*M_PI;
454  if (phiProj < -1.*M_PI) phiProj += 2.0*M_PI;
455  r[1] = stgcHits.at(jHitId).r;
456  z[1] = phiProj;
457  slope = (z[0]+z[1])/2.;
458  intercept = (r[0]+r[1])/2.;
459  if(std::abs(r[0]*std::sin(z[0]) - r[1]*std::sin(z[1])) > 300.) continue;
460  }
461 
462  unsigned int encodedIds = (iHitId<<16) + jHitId;
463  hitIdsInTwo[iPair].push_back(encodedIds);
464  slopeInTwo[iPair].push_back(slope);
465  interceptInTwo[iPair].push_back(intercept);
466 
467  foundCounterpart = 1;
468  foundCounterparts[jHit] = 1;
469  }//end of jHit in the (i+4)-th layer
470  if(!foundCounterpart){ // in case of no counterpart in the (i+4)-th layer
471  unsigned int encodedIds = (iHitId<<16) + 0xffff; // fill all bits with 1 for hit id for the layer with no hit
472  hitIdsInTwo[iPair].push_back(encodedIds);
473  if(isStrip) {
474  slopeInTwo[iPair].push_back(r[0]/z[0]);
475  interceptInTwo[iPair].push_back(0.);
476  } else {
477  slopeInTwo[iPair].push_back(z[0]);
478  interceptInTwo[iPair].push_back(r[0]);
479  }
480  }
481  }//end of iHit in the i-th layer
482  // Loop over hits in the (i+4)-th layer
483  for(unsigned int jHit = 0; jHit < nHitsInOuter; ++jHit){
484  if (!foundCounterparts[jHit]) {
485  int jHitId = hitIdByLayer[iPair+4].at(jHit);
486  unsigned int encodedIds = 0xffff0000 + jHitId; // fill all bits with 1 for hit id for the layer with no hit
487  hitIdsInTwo[iPair].push_back(encodedIds);
488  if(isStrip) {
489  slopeInTwo[iPair].push_back(stgcHits.at(jHitId).r/stgcHits.at(jHitId).z);
490  interceptInTwo[iPair].push_back(0.);
491  } else {
492  double localPhiCenter;
493  if (stgcHits.at(jHitId).stationPhi<=5) {
494  localPhiCenter = 0.25 * M_PI * ((double)stgcHits.at(jHitId).stationPhi-1.);
495  } else {
496  localPhiCenter = 0.25 * M_PI * ((double)stgcHits.at(jHitId).stationPhi-9.);
497  }
498  if (stgcHits.at(jHitId).stationName == 57){
499  localPhiCenter += M_PI/8.; // small sTGC sectors
500  if (stgcHits.at(jHitId).stationPhi == 5) localPhiCenter -= 2 * M_PI;
501  }
502 
503  double phiProj = stgcHits.at(jHitId).phi - localPhiCenter;
504  if (phiProj > M_PI) phiProj -= 2.0*M_PI;
505  if (phiProj < -1.*M_PI) phiProj += 2.0*M_PI;
506  slopeInTwo[iPair].push_back(phiProj);
507  interceptInTwo[iPair].push_back(stgcHits.at(jHitId).r);
508  }
509  }
510  }
511  }//end of pair loop
512 
513  ATH_MSG_DEBUG("@@STGC@@ isStrip= " << isStrip << " Npairs " << hitIdsInTwo[0].size() << " " << hitIdsInTwo[1].size() << " " << hitIdsInTwo[2].size() << " " << hitIdsInTwo[3].size());
514  for (unsigned int iLayer = 0; iLayer < 4; ++ iLayer) {
515  for (unsigned int iPair = 0; iPair < slopeInTwo[iLayer].size(); ++iPair) {
516  ATH_MSG_DEBUG("@@STGC@@ pair fit isStrip= " << isStrip << " slope= " << slopeInTwo[iLayer].at(iPair) << " intercept= " << interceptInTwo[iLayer].at(iPair));
517  }
518  }
519 
520  std::array<std::vector<unsigned long int>,2> hitIdsInFour;
521  std::array<std::vector<double>,2> slopeInFour;
522  std::array<std::vector<double>,2> interceptInFour;
523  for(unsigned int iQuad = 0; iQuad < 2; ++iQuad){
524  unsigned int nPairsInInner = hitIdsInTwo[iQuad].size();
525  unsigned int nPairsInOuter = hitIdsInTwo[iQuad+2].size();
526 
527  bool foundCounterparts[0xffff] = {};
528  for(unsigned int iPair = 0; iPair < nPairsInInner; ++iPair){
529  bool foundCounterpart = 0;
530  double slope[2];
531  double intercept[2];
532  slope[0] = slopeInTwo[iQuad].at(iPair);
533  intercept[0] = interceptInTwo[iQuad].at(iPair);
534 
535  for(unsigned int jPair = 0; jPair < nPairsInOuter; ++jPair){
536  unsigned int ihitIds = hitIdsInTwo[iQuad].at(iPair);
537  unsigned int jhitIds = hitIdsInTwo[iQuad+2].at(jPair);
538  if ( !(((ihitIds>>16 & 0xffff) != 0xffff || (ihitIds & 0xffff) != 0xffff) &&
539  ((jhitIds>>16 & 0xffff) != 0xffff || (jhitIds & 0xffff) != 0xffff )) ) continue; // require at least 2 hits in 4 layers
540 
541  slope[1] = slopeInTwo[iQuad+2].at(jPair);
542  intercept[1] = interceptInTwo[iQuad+2].at(jPair);
543 
544  if (isStrip) {
545  double spR0 = slope[0] * NSWCenterZ + intercept[0];
546  double spR1 = slope[1] * NSWCenterZ + intercept[1];
547  if(std::abs(spR1 - spR0) > 50.) continue; //OPTIMIZE ME!!!
548  } else {
549  if(std::abs(intercept[0]*std::sin(slope[0]) - intercept[1]*std::sin(slope[1])) > 100.) continue;
550  }
551 
552  foundCounterpart = 1;
553  foundCounterparts[jPair] = 1;
554 
555  unsigned long int encodedIds = (hitIdsInTwo[iQuad].at(iPair) << 32 ) + hitIdsInTwo[iQuad+2].at(jPair);
556  hitIdsInFour[iQuad].push_back(encodedIds);
557  slopeInFour[iQuad].push_back((slope[1] + slope[0])/2.);
558  interceptInFour[iQuad].push_back((intercept[1] + intercept[0])/2.);
559  }// end of iPair of the i-th and (i+4)-th layers, i=0,2
560 
561  if(foundCounterpart) continue; // in case of no counterpart of pairs in the inner layers
562  if((hitIdsInTwo[iQuad].at(iPair)>>16 & 0xffff) == 0xffff || (hitIdsInTwo[iQuad].at(iPair) & 0xffff) == 0xffff) continue;
563 
564  unsigned long int encodedIds = (hitIdsInTwo[iQuad].at(iPair) << 32 ) + 0xffffffff;
565  hitIdsInFour[iQuad].push_back(encodedIds);
566  slopeInFour[iQuad].push_back(slope[0]);
567  interceptInFour[iQuad].push_back(intercept[0]);
568  }// end of jPair of the j-th and (j+4)-th layers, j=1,3
569  for (unsigned int jPair = 0; jPair < nPairsInOuter; ++jPair) {
570  if(foundCounterparts[jPair]) continue; // in case of no counterpart of pairs in the outner layers
571  if((hitIdsInTwo[iQuad+2].at(jPair)>>16 & 0xffff) == 0xffff || (hitIdsInTwo[iQuad+2].at(jPair) & 0xffff) == 0xffff) continue;
572 // unsigned long int encodedIds = (0xffffffff << 32) + hitIdsInTwo[iQuad+2].at(jPair);
573  unsigned long int encodedIds = (0xffffffff00000000 ) + hitIdsInTwo[iQuad+2].at(jPair);
574  hitIdsInFour[iQuad].push_back(encodedIds);
575  slopeInFour[iQuad].push_back(slopeInTwo[iQuad+2].at(jPair));
576  interceptInFour[iQuad].push_back(interceptInTwo[iQuad+2].at(jPair));
577  }
578  }// end of quad loop
579 
580  ATH_MSG_DEBUG("@@STGC@@ isStrip= " << isStrip << " Nquads " << hitIdsInFour[0].size() << " " << hitIdsInFour[1].size());
581  for (unsigned int iLayer = 0; iLayer < 2; ++ iLayer) {
582  for (unsigned int iQuad = 0; iQuad < slopeInFour[iLayer].size(); ++iQuad) {
583  ATH_MSG_DEBUG("@@STGC@@ quad fit isStrip= " << isStrip << " slope= " << slopeInFour[iLayer].at(iQuad) << " intercept= " << interceptInFour[iLayer].at(iQuad));
584  }
585  }
586 
587  std::vector< std::array<int, 8> > hitIdsInEight;
588  std::vector<double> mseInEight;
589 
590  unsigned int nQuadInInner = hitIdsInFour[0].size();
591  unsigned int nQuadInOuter = hitIdsInFour[1].size();
592 
593  for(unsigned int iQuad = 0; iQuad < nQuadInInner; ++iQuad){
594  double slope[2];
595  double intercept[2];
596  slope[0] = slopeInFour[0].at(iQuad);
597  intercept[0] = interceptInFour[0].at(iQuad);
598 
599  for(unsigned int jQuad = 0; jQuad < nQuadInOuter; ++jQuad){
600  unsigned long int ihitIds = hitIdsInFour[0].at(iQuad);
601  unsigned long int jhitIds = hitIdsInFour[1].at(jQuad);
602  int nOfLayersWithNoHit = 0;
603  for (unsigned int iLayer = 0; iLayer < 4; ++iLayer) {
604  if ( (ihitIds>>(3-iLayer)*16 & 0xffff) == 0xffff ) {++nOfLayersWithNoHit;}
605  if ( (jhitIds>>(3-iLayer)*16 & 0xffff) == 0xffff ) {++nOfLayersWithNoHit;}
606  }
607  if (nOfLayersWithNoHit > 4) continue; // require at least 4 hits in 8 layers
608 
609  slope[1] = slopeInFour[1].at(jQuad);
610  intercept[1] = interceptInFour[1].at(jQuad);
611 
612  if(isStrip) {
613  double spR0 = slope[0] * NSWCenterZ + intercept[0];
614  double spR1 = slope[1] * NSWCenterZ + intercept[1];
615  if(std::abs(spR1 - spR0) > 10. ||
616  std::abs(intercept[1] + intercept[0]) / 2 > 100.) continue; // OPTIMIZE ME!!!
617  } else {
618  if(std::abs(intercept[0]*std::sin(slope[0]) - intercept[1]*std::sin(slope[1])) > 100. ) continue;
619  }
620 
621  std::array<int,8> setOfHitIds = {-1,-1,-1,-1,-1,-1,-1,-1};
622  std::vector<double> r, z;
623  for(unsigned int i = 0; i < 8; ++i) {
624  unsigned int iHitId, iLayer = 0;
625  if (i <= 3) {
626  iHitId = (unsigned int) ((ihitIds>>(3-i)*16) & 0xffff);
627  iLayer = i+3*(i%2);
628  } else {
629  iHitId = (unsigned int) ((jhitIds>>(3-i%4)*16) & 0xffff);
630  iLayer = (i-4)+3*(i%2)+1;
631  }
632  if ( iHitId != 0xffff ) {
633  if (isStrip) {
634  r.push_back(stgcHits.at(iHitId).r);
635  z.push_back(stgcHits.at(iHitId).z);
636  }
637  else {
638  double localPhiCenter;
639  if (stgcHits.at(iHitId).stationPhi<=5) {
640  localPhiCenter = 0.25 * M_PI * ((double)stgcHits.at(iHitId).stationPhi-1.);
641  } else {
642  localPhiCenter = 0.25 * M_PI * ((double)stgcHits.at(iHitId).stationPhi-9.);
643  }
644  if (stgcHits.at(iHitId).stationName == 57){
645  localPhiCenter += M_PI/8.; // small sTGC sectors
646  if (stgcHits.at(iHitId).stationPhi == 5) localPhiCenter -= 2 * M_PI;
647  }
648 
649  double phiProj = stgcHits.at(iHitId).phi - localPhiCenter;
650  if (phiProj > M_PI) phiProj -= 2.0*M_PI;
651  if (phiProj < -1.*M_PI) phiProj += 2.0*M_PI;
652  r.push_back(phiProj);
653  z.push_back(stgcHits.at(iHitId).r);
654  }
655  setOfHitIds[iLayer] = iHitId;
656  ATH_MSG_DEBUG("@@STGC@@ strip_pos iHitId " << iLayer << " " << iHitId);
657  }
658  }
659  double slopefit=0., interceptfit=99999., mse =-1.;
660  if (isStrip) {
661  LinearFit(z,r,&slopefit,&interceptfit,&mse);
662  } else {
663  double phiavg = 0.;
664  for (unsigned int iHit = 0; iHit < r.size(); ++iHit){
665  phiavg += r.at(iHit);
666  }
667  phiavg /= r.size();
668  mse = 0.;
669  for (unsigned int iHit = 0; iHit < r.size(); ++iHit){
670  mse += std::pow(r.at(iHit) - phiavg,2);
671  }
672  }
673  hitIdsInEight.push_back(setOfHitIds);
674  mseInEight.push_back(mse);
675  } // end of quad loop in the outer layers
676  } // end of quad loop in the inner layers
677  if(!hitIdsInEight.size()){
678  ATH_MSG_DEBUG("No candidate segment found in STGC");
679  return;
680  }
681 
682  ATH_MSG_DEBUG("@@STGC@@ isStrip= " << isStrip << " Noctets " << hitIdsInEight.size());
683  std::vector<int> nOctetSegments;
684  std::vector<int> patternStationName;
685 
686  for (unsigned int iOctet = 0; iOctet < hitIdsInEight.size(); ++iOctet) {
687 
688  bool isFirstHit = true;
689  int hitStationName = 0;
690 
691  int nOctetSegment = 0;
692  ATH_MSG_DEBUG("@@STGC@@ octet fit isStrip= " << isStrip << " mse " << mseInEight.at(iOctet));
693  std::array<int, 8> tmpOctet = hitIdsInEight.at(iOctet);
694  for (unsigned int iLayer = 0; iLayer < 8; ++iLayer) {
695  if (tmpOctet[iLayer] != -1) {
696 
697  if(isFirstHit){
698  hitStationName = stgcHits.at(tmpOctet[iLayer]).stationName;
699  isFirstHit = false;
700 
701  ATH_MSG_DEBUG("@@STGC@@ octet pos isStrip= " << isStrip << " r= " << stgcHits.at(tmpOctet[iLayer]).r << " phi= " << stgcHits.at(tmpOctet[iLayer]).phi << " z= " << stgcHits.at(tmpOctet[iLayer]).z);
702  nOctetSegment++;
703  }
704  else if(stgcHits.at(tmpOctet[iLayer]).stationName == hitStationName){
705  ATH_MSG_DEBUG("@@STGC@@ octet pos isStrip= " << isStrip << " r= " << stgcHits.at(tmpOctet[iLayer]).r << " phi= " << stgcHits.at(tmpOctet[iLayer]).phi << " z= " << stgcHits.at(tmpOctet[iLayer]).z);
706  nOctetSegment++;
707  }
708  }
709  }
710  nOctetSegments.push_back(nOctetSegment);
711  patternStationName.push_back(hitStationName);
712  }
713  double nOcSegMax = 0;
714  for(unsigned int iOctet = 0; iOctet < hitIdsInEight.size(); ++iOctet){
715  if(nOctetSegments.at(iOctet) > nOcSegMax){
716  nOcSegMax = nOctetSegments.at(iOctet);
717  }
718  }
719 
720  double msemin = 1000000.;
721  double mseminWireL = 1000000.; // Large sector
722  double mseminWireS = 1000000.; // Small sector
723 
724  std::vector<int> octetIds(2,-1);
725 
726  if(isStrip){
727  for(unsigned int iOctet = 0; iOctet < hitIdsInEight.size(); ++iOctet){
728  if(nOctetSegments.at(iOctet) != nOcSegMax){
729  continue;
730  }
731  if( mseInEight.at(iOctet) < msemin) {
732  msemin = mseInEight.at(iOctet);
733  }
734  }// end of Octet loop
735  for(unsigned int iOctet = 0; iOctet < hitIdsInEight.size(); ++iOctet){
736  if(nOctetSegments.at(iOctet) != nOcSegMax) continue;
737  if(mseInEight.at(iOctet) != msemin){
738  continue;
739  }
740  octetIds.push_back(iOctet);
741  }
742  } else {
743  for(unsigned int iOctet = 0; iOctet < hitIdsInEight.size(); ++iOctet){
744  if(patternStationName.at(iOctet) == 58){
745  if( mseInEight.at(iOctet) < mseminWireL) {
746  mseminWireL = mseInEight.at(iOctet);
747  }
748  }
749  else if(patternStationName.at(iOctet) == 57){
750  if( mseInEight.at(iOctet) < mseminWireS) {
751  mseminWireS = mseInEight.at(iOctet);
752  }
753  }
754  }// end of Octet loop
755  for(unsigned int iOctet = 0; iOctet < hitIdsInEight.size(); ++iOctet){
756  if(patternStationName.at(iOctet) == 58){
757  if(mseInEight.at(iOctet) != mseminWireL){
758  continue;
759  }
760  }
761  else if(patternStationName.at(iOctet) == 57){
762  if(mseInEight.at(iOctet) != mseminWireS){
763  continue;
764  }
765  }
766  octetIds.push_back(iOctet);
767  }
768  }
769  for(unsigned int ids = 0; ids < octetIds.size(); ids++){
770  if (octetIds.at(ids) != -1) {
771  hitIdsCandidate.push_back(hitIdsInEight.at(octetIds.at(ids)));
772  }
773  }
774 }

◆ findStgcHitsInSegment()

StatusCode TrigL2MuonSA::NswStationFitter::findStgcHitsInSegment ( TrigL2MuonSA::StgcHits &  stgcHits) const

Definition at line 269 of file NswStationFitter.cxx.

270 {
271  if(stgcHits.size() == 0) return StatusCode::SUCCESS;
272  int hitsInRoad = 0;
273  for(unsigned int iHit = 0; iHit < stgcHits.size(); iHit++){
274  if(stgcHits.at(iHit).isOutlier == 0){
275  hitsInRoad++;
276  stgcHits.at(iHit).isOutlier = 1;
277  }
278  }
279  if(hitsInRoad == 0) return StatusCode::SUCCESS;
280 
281  if(hitsInRoad < 9) {
282  ATH_MSG_DEBUG("Number of STGC hits is too small, at least 9 hits required : "<<hitsInRoad<<" hits");
283  return StatusCode::SUCCESS;
284  } else if(hitsInRoad > 100) {
285  ATH_MSG_WARNING("Number of STGC hits is too large, at most 100 hits allowed : "<<hitsInRoad<<" hits");
286  return StatusCode::SUCCESS;
287  }
288 
289  // cassifyDataEachLayer
290  std::array<std::vector<int>, 8> strHitIdByLayer;
291  std::array<std::vector<int>, 8> wireHitIdByLayer;
292  for(unsigned int iHit = 0; iHit < stgcHits.size(); ++iHit){
293  if(stgcHits.at(iHit).isOutlier != 1) continue;
294  int layerNumber = stgcHits.at(iHit).layerNumber;
295  if (layerNumber > 7) {
296  ATH_MSG_WARNING("STGC hit layer number > 7");
297  continue;
298  }
299  if(stgcHits.at(iHit).channelType == 1){
300  strHitIdByLayer[layerNumber].push_back(iHit);
301  }else if(stgcHits.at(iHit).channelType == 2){
302  wireHitIdByLayer[layerNumber].push_back(iHit);
303  }
304  }
305  ATH_MSG_DEBUG("@@STGC@@ strip Nhits " << strHitIdByLayer[0].size()
306  << " " << strHitIdByLayer[1].size()
307  << " " << strHitIdByLayer[2].size()
308  << " " << strHitIdByLayer[3].size()
309  << " " << strHitIdByLayer[4].size()
310  << " " << strHitIdByLayer[5].size()
311  << " " << strHitIdByLayer[6].size()
312  << " " << strHitIdByLayer[7].size());
313  ATH_MSG_DEBUG("@@STGC@@ wire Nhits " << wireHitIdByLayer[0].size()
314  << " " << wireHitIdByLayer[1].size()
315  << " " << wireHitIdByLayer[2].size()
316  << " " << wireHitIdByLayer[3].size()
317  << " " << wireHitIdByLayer[4].size()
318  << " " << wireHitIdByLayer[5].size()
319  << " " << wireHitIdByLayer[6].size()
320  << " " << wireHitIdByLayer[7].size());
321 
322  std::vector< std::array<int, 8> > strHitIds; // Candidates' sets of strip hit ids in 8 layers
323  findSetOfStgcHitIds(stgcHits, strHitIdByLayer, strHitIds);
324  std::vector< std::array<int, 8> > wireHitIds; // Candidates' sets of wire hit ids in 8 layers
325  findSetOfStgcHitIds(stgcHits, wireHitIdByLayer, wireHitIds);
326  ATH_MSG_DEBUG("@@STGC@@ strip wire " << strHitIds.size() << " " << wireHitIds.size());
327 
328  bool isLargeStrip = false;
329  bool isSmallStrip = false;
330  for (unsigned int iHit = 0; iHit < strHitIds.size(); ++iHit) {
331  std::array<int, 8> hitIds = strHitIds.at(iHit);
332  for (unsigned int iLayer = 0; iLayer < 8; ++iLayer) {
333  if (hitIds[iLayer] != -1) {
334  if(stgcHits.at(hitIds[iLayer]).stationName == 58) isLargeStrip = true;
335  else if(stgcHits.at(hitIds[iLayer]).stationName == 57) isSmallStrip = true;
336  stgcHits.at(hitIds[iLayer]).isOutlier = 0;
337  }
338  }
339  }
340  for (unsigned int iHit = 0; iHit < wireHitIds.size(); ++iHit) {
341  std::array<int, 8> hitIds = wireHitIds.at(iHit);
342  for (unsigned int iLayer = 0; iLayer < 8; ++iLayer) {
343  if (hitIds[iLayer] != -1) {
344  if(isLargeStrip){
345  if(stgcHits.at(hitIds[iLayer]).stationName == 58){
346  stgcHits.at(hitIds[iLayer]).isOutlier = 0;
347  }
348  }
349  else if(isSmallStrip){
350  if(stgcHits.at(hitIds[iLayer]).stationName == 57){
351  stgcHits.at(hitIds[iLayer]).isOutlier = 0;
352  }
353  }
354  }
355  }
356  }
357  return StatusCode::SUCCESS;
358 
359 }

◆ getNswResolution()

void TrigL2MuonSA::NswStationFitter::getNswResolution ( double *  stgcDeltaR,
double *  mmDeltaR,
unsigned int  size 
) const

Definition at line 914 of file NswStationFitter.cxx.

915 {
916  double RmsDeltarEtaStgc[12] = {1.43,1.53,1.53,1.56,1.59,1.54,1.70,1.69,1.76,1.81,1.83,1.84};
917  double RmsDeltarEtaMm[12] = {0.49,0.46,0.48,0.40,0.39,0.39,0.38,0.35,0.36,0.33,0.33,0.40};
918  for(unsigned int bin=0; bin < size; bin++){
919  stgcDeltaR[bin] = RmsDeltarEtaStgc[bin];
920  mmDeltaR[bin] = RmsDeltarEtaMm[bin];
921  }
922 }

◆ inputHandles()

virtual std::vector<Gaudi::DataHandle*> AthCommonDataStore< AthCommonMsg< AlgTool > >::inputHandles ( ) const
overridevirtualinherited

Return this algorithm's input handles.

We override this to include handle instances from key arrays if they have not yet been declared. See comments on updateVHKA.

◆ LinearFit()

void TrigL2MuonSA::NswStationFitter::LinearFit ( std::vector< double > &  x,
std::vector< double > &  y,
double *  slope,
double *  intercept,
double *  mse 
) const

Definition at line 806 of file NswStationFitter.cxx.

808 {
809  double sumX=0, sumY=0, sumXY=0, sumX2=0;
810  int nHits = x.size();
811  *mse = 0.;
812  for (unsigned int iHit = 0; iHit < x.size(); ++iHit){
813  sumX += x.at(iHit);
814  sumY += y.at(iHit);
815  sumXY += x.at(iHit)*y.at(iHit);
816  sumX2 += x.at(iHit)*x.at(iHit);
817  }
818 
819  if(nHits > 1) {
820  if(nHits * sumX2 - (sumX * sumX) > ZERO_LIMIT) {
821  *slope = (nHits * sumXY - sumX * sumY) / (nHits * sumX2 - (sumX * sumX));
822  *intercept = (sumX2 * sumY - sumXY * sumX) / (nHits * sumX2 - sumX * sumX);
823  } else {
824  *slope = 0.;
825  *intercept = sumY/nHits;
826  }
827  }
828  else if(nHits == 1) {
829  *slope = sumY/sumX;
830  *intercept = 0.;
831  }
832 
833  if(nHits > 2) {
834  for(unsigned int iHit = 0; iHit< x.size(); ++iHit){
835  *mse += std::pow(y.at(iHit) - (*slope * x.at(iHit) + *intercept), 2.0);
836  }
837  *mse = *mse / (nHits - 2);
838  }
839  else {
840  *mse = 1000.;
841  }
842 }

◆ LinearFitWeight()

void TrigL2MuonSA::NswStationFitter::LinearFitWeight ( std::vector< double > &  x,
std::vector< double > &  y,
std::vector< bool > &  isStgc,
double *  slope,
double *  intercept,
double *  mse,
double  eta 
) const

Definition at line 844 of file NswStationFitter.cxx.

846 {
847  double RmsDeltarEtaStgc[12] = {};
848  double RmsDeltarEtaMm[12]= {};
849  getNswResolution(RmsDeltarEtaStgc, RmsDeltarEtaMm, 12);
850 
851  double weightStgc[12] = {};
852  double weightMm[12] = {};
853  for(int i_weight=0; i_weight<12; i_weight++){
854  weightStgc[i_weight] = 1/std::pow(RmsDeltarEtaStgc[i_weight],2);
855  weightMm[i_weight] = 1/std::pow(RmsDeltarEtaMm[i_weight],2);
856  }
857  int weightBin = 0;
858  double minEta = 1.3;
859  double maxEta = 1.4;
860  for(int iBin=0; iBin<12; iBin++){
861  if(std::abs(eta) >= minEta && std::abs(eta) < maxEta){
862  weightBin = iBin;
863  break;
864  } else {
865  minEta += 0.1;
866  maxEta += 0.1;
867  }
868  }
869 
870  double sumX=0, sumY=0, sumXY=0, sumX2=0, sumW=0;
871  int nHits = x.size();
872  *mse = 0.;
873  for (unsigned int iHit = 0; iHit < x.size(); ++iHit){
874  if(isStgc.at(iHit)){
875  sumX += weightStgc[weightBin] * x.at(iHit);
876  sumY += weightStgc[weightBin] * y.at(iHit);
877  sumXY += weightStgc[weightBin] * x.at(iHit) * y.at(iHit);
878  sumX2 += weightStgc[weightBin] * x.at(iHit) * x.at(iHit);
879  sumW += weightStgc[weightBin];
880  } else {
881  sumX += weightMm[weightBin] * x.at(iHit);
882  sumY += weightMm[weightBin] * y.at(iHit);
883  sumXY += weightMm[weightBin] * x.at(iHit) * y.at(iHit);
884  sumX2 += weightMm[weightBin] * x.at(iHit) * x.at(iHit);
885  sumW += weightMm[weightBin];
886  }
887  }
888 
889  if(nHits > 1) {
890  if(nHits * sumX2 - (sumX * sumX) > ZERO_LIMIT) {
891  *slope = (sumW * sumXY - sumX * sumY) / (sumW * sumX2 - (sumX * sumX));
892  *intercept = (sumX2 * sumY - sumXY * sumX) / (sumW * sumX2 - sumX * sumX);
893  } else {
894  *slope = 0.;
895  *intercept = sumY/nHits;
896  }
897  }
898  else if(nHits == 1) {
899  *slope = sumY/sumX;
900  *intercept = 0.;
901  }
902 
903  if(nHits > 2) {
904  for(unsigned int iHit = 0; iHit< x.size(); ++iHit){
905  *mse += std::pow((y.at(iHit) - (*slope * x.at(iHit) + *intercept)), 2.0);
906  }
907  *mse = *mse / (nHits - 2);
908  }
909  else {
910  *mse = 1000.;
911  }
912 }

◆ MakeSegment() [1/2]

StatusCode TrigL2MuonSA::NswStationFitter::MakeSegment ( TrigL2MuonSA::TrackPattern &  trackPattern,
TrigL2MuonSA::MmHits &  mmHits 
) const

Definition at line 790 of file NswStationFitter.cxx.

792 {
793  TrigL2MuonSA::MmHits selectedMmHits;
794  selectedMmHits.clear();
795  if(mmHits.size() == 0) return StatusCode::SUCCESS;
796  for(unsigned int iHit = 0; iHit < mmHits.size(); iHit++){
797  if(mmHits.at(iHit).isOutlier != 0) continue;
798  selectedMmHits.push_back(mmHits.at(iHit));
799  }
800  trackPattern.mmSegment.clear();
801  trackPattern.mmSegment = selectedMmHits;
802  return StatusCode::SUCCESS;
803 }

◆ MakeSegment() [2/2]

StatusCode TrigL2MuonSA::NswStationFitter::MakeSegment ( TrigL2MuonSA::TrackPattern &  trackPattern,
TrigL2MuonSA::StgcHits &  stgcHits 
) const

Definition at line 776 of file NswStationFitter.cxx.

778 {
779  TrigL2MuonSA::StgcHits selectedStgcHits;
780  selectedStgcHits.clear();
781  if(stgcHits.size() == 0) return StatusCode::SUCCESS;
782  for(unsigned int iHit = 0; iHit < stgcHits.size(); iHit++){
783  if(stgcHits.at(iHit).isOutlier != 0) continue;
784  selectedStgcHits.push_back(stgcHits.at(iHit));
785  }
786  trackPattern.stgcSegment.clear();
787  trackPattern.stgcSegment = selectedStgcHits;
788  return StatusCode::SUCCESS;
789 }

◆ msg() [1/2]

MsgStream& AthCommonMsg< AlgTool >::msg ( ) const
inlineinherited

Definition at line 24 of file AthCommonMsg.h.

24  {
25  return this->msgStream();
26  }

◆ msg() [2/2]

MsgStream& AthCommonMsg< AlgTool >::msg ( const MSG::Level  lvl) const
inlineinherited

Definition at line 27 of file AthCommonMsg.h.

27  {
28  return this->msgStream(lvl);
29  }

◆ msgLvl()

bool AthCommonMsg< AlgTool >::msgLvl ( const MSG::Level  lvl) const
inlineinherited

Definition at line 30 of file AthCommonMsg.h.

30  {
31  return this->msgLevel(lvl);
32  }

◆ outputHandles()

virtual std::vector<Gaudi::DataHandle*> AthCommonDataStore< AthCommonMsg< AlgTool > >::outputHandles ( ) const
overridevirtualinherited

Return this algorithm's output handles.

We override this to include handle instances from key arrays if they have not yet been declared. See comments on updateVHKA.

◆ renounce()

std::enable_if_t<std::is_void_v<std::result_of_t<decltype(&T::renounce)(T)> > && !std::is_base_of_v<SG::VarHandleKeyArray, T> && std::is_base_of_v<Gaudi::DataHandle, T>, void> AthCommonDataStore< AthCommonMsg< AlgTool > >::renounce ( T &  h)
inlineprotectedinherited

Definition at line 380 of file AthCommonDataStore.h.

381  {
382  h.renounce();
383  PBASE::renounce (h);
384  }

◆ renounceArray()

void AthCommonDataStore< AthCommonMsg< AlgTool > >::renounceArray ( SG::VarHandleKeyArray &  handlesArray)
inlineprotectedinherited

remove all handles from I/O resolution

Definition at line 364 of file AthCommonDataStore.h.

364  {
365  handlesArray.renounce();
366  }

◆ selectMmHits()

StatusCode TrigL2MuonSA::NswStationFitter::selectMmHits ( const TrigRoiDescriptor *  p_roids,
TrigL2MuonSA::MmHits &  mmHits 
) const

Definition at line 134 of file NswStationFitter.cxx.

136 {
137 
138  TrigL2MuonSA::MmHits selectedMmHits;
139  selectedMmHits.clear();
140 
141  // define region where RoI is near
142  double etaMin = p_roids->eta() - 1.;
143  double etaMax = p_roids->eta() + 1.;
144  double phiMin = p_roids->phi() - 1.;
145  double phiMax = p_roids->phi() + 1.;
146  if( phiMin > M_PI ) phiMin -= 2*M_PI;
147  if( phiMax > M_PI ) phiMax -= 2*M_PI;
148  if( phiMin < -1.*M_PI ) phiMin += 2*M_PI;
149  if( phiMax < -1.*M_PI ) phiMax += 2*M_PI;
150  // boolian to check if sTGC hits are near RoI
151  bool inRoiEta, inRoiPhi, inRoi;
152 
153  // loop over MM digits.
154  unsigned int iHit;
155  if (mmHits.size()>0) {
156  for (iHit = 0; iHit < mmHits.size(); iHit++){
157 
158  // Get the digit point.
159  TrigL2MuonSA::MmHitData& hit = mmHits[iHit];
160 
161  // check if MM hits are near RoI
162  inRoiEta = false; inRoiPhi = false; inRoi = false;
163  if ( etaMin<=hit.eta && hit.eta<=etaMax ) inRoiEta = true;
164  if ( phiMin<=phiMax && phiMin<=hit.phi && hit.phi<=phiMax ) inRoiPhi= true;
165  if ( phiMin>phiMax && (phiMin<=hit.phi || hit.phi<=phiMax)) inRoiPhi= true;
166  if ( inRoiEta && inRoiPhi ) inRoi = true;
167  ATH_MSG_DEBUG("MM hits eta = "<<hit.eta<<", phi = "<<hit.phi<<", r = "<<hit.r<<", z = "<<hit.z<<", stationEta = "<<hit.stationEta<<", stationPhi = "<<hit.stationPhi<<", stationName = "<<hit.stationName<<", matched RoI? = "<<inRoi);
168 
169  // pushback if the hit is near RoI
170  if (!inRoi){
171  continue;
172  }
173  selectedMmHits.push_back(hit);
174  }
175  }
176 
177  mmHits.clear();
178  mmHits = selectedMmHits;
179 
180  return StatusCode::SUCCESS;
181 
182 }

◆ selectStgcHits()

StatusCode TrigL2MuonSA::NswStationFitter::selectStgcHits ( const TrigRoiDescriptor *  p_roids,
TrigL2MuonSA::StgcHits &  stgcHits 
) const

Definition at line 86 of file NswStationFitter.cxx.

88 {
89 
90  TrigL2MuonSA::StgcHits selectedStgcHits;
91  selectedStgcHits.clear();
92 
93  // define region where RoI is near
94  double etaMin = p_roids->eta() - 1.;
95  double etaMax = p_roids->eta() + 1.;
96  double phiMin = p_roids->phi() - 1.;
97  double phiMax = p_roids->phi() + 1.;
98  if( phiMin > M_PI ) phiMin -= 2*M_PI;
99  if( phiMax > M_PI ) phiMax -= 2*M_PI;
100  if( phiMin < -1.*M_PI ) phiMin += 2*M_PI;
101  if( phiMax < -1.*M_PI ) phiMax += 2*M_PI;
102  // boolian to check if sTGC hits are near RoI
103  bool inRoiEta, inRoiPhi, inRoi;
104 
105  // loop over sTGC digits.
106  unsigned int iHit;
107  if (stgcHits.size()>0) {
108  for (iHit = 0; iHit < stgcHits.size(); iHit++){
109  // Get the digit point.
110  TrigL2MuonSA::StgcHitData& hit = stgcHits[iHit];
111 
112  // check if sTGC hits are near RoI
113  inRoiEta = false; inRoiPhi = false; inRoi = false;
114  if ( etaMin<=hit.eta && hit.eta<=etaMax ) inRoiEta = true;
115  if ( phiMin<=phiMax && phiMin<=hit.phi && hit.phi<=phiMax ) inRoiPhi= true;
116  if ( phiMin>phiMax && (phiMin<=hit.phi || hit.phi<=phiMax)) inRoiPhi= true;
117  if ( inRoiEta && inRoiPhi ) inRoi = true;
118  ATH_MSG_DEBUG("sTGC hits eta = "<<hit.eta<<", phi = "<<hit.phi<<", r = "<<hit.r<<", z = "<<hit.z<<", stationEta = "<<hit.stationEta<<", stationPhi = "<<hit.stationPhi<<", channelType = "<<hit.channelType<<", stationName = "<<hit.stationName<<",layerNumber ="<<hit.layerNumber<<", matched RoI? = "<<inRoi);
119  // pushback if the hit is near RoI
120  if (!inRoi){
121  continue;
122  }
123  selectedStgcHits.push_back(hit);
124  }
125  }
126 
127  stgcHits.clear();
128  stgcHits = selectedStgcHits;
129 
130  return StatusCode::SUCCESS;
131 
132 }

◆ superPointFitter()

StatusCode TrigL2MuonSA::NswStationFitter::superPointFitter ( const TrigRoiDescriptor *  p_roids,
TrigL2MuonSA::TrackPattern &  trackPattern,
TrigL2MuonSA::StgcHits &  stgcHits,
TrigL2MuonSA::MmHits &  mmHits 
) const

Definition at line 33 of file NswStationFitter.cxx.

37 {
38 
39  ATH_MSG_DEBUG("NswStationFitter::findSuperPoints() was called.");
40 
41  // selection for sTGC hits, RoI matching based
42  ATH_CHECK( selectStgcHits(p_roids,trackPattern.stgcSegment) );
43  ATH_CHECK( selectMmHits(p_roids,trackPattern.mmSegment) );
44 
45  ATH_CHECK( findStgcHitsInSegment(stgcHits) );
46  ATH_CHECK( findMmHitsInSegment(mmHits) );
47 
48  bool isLargeStgc = false;
49  bool isSmallStgc = false;
50  if(stgcHits.size() != 0){
51  for(unsigned int iHit = 0; iHit < stgcHits.size(); iHit++){
52  if(stgcHits.at(iHit).isOutlier == 0){
53  if(stgcHits.at(iHit).stationName == 58) isLargeStgc = true;
54  else if(stgcHits.at(iHit).stationName == 57) isSmallStgc = true;
55  if(isLargeStgc && isSmallStgc) continue;
56  }
57  }
58  }
59  if(mmHits.size() != 0){
60  for(unsigned int iHit = 0; iHit < mmHits.size(); iHit++){
61  if(mmHits.at(iHit).isOutlier == 0){
62  if(isLargeStgc){
63  if(mmHits.at(iHit).stationName == 55) mmHits.at(iHit).isOutlier = 1;
64  } else if(isSmallStgc) {
65  if(mmHits.at(iHit).stationName == 56) mmHits.at(iHit).isOutlier = 1;
66  }
67  }
68  }
69  }
70  ATH_CHECK( MakeSegment(trackPattern,stgcHits) );
71  ATH_CHECK( MakeSegment(trackPattern,mmHits) );
72 
73  ATH_MSG_DEBUG("Number of sTGC and MM hits for SPs " << trackPattern.stgcSegment.size() << " " << trackPattern.mmSegment.size());
74  if (trackPattern.stgcSegment.size() < 9 && trackPattern.mmSegment.size() < 6) {
75  ATH_MSG_DEBUG("Number of sTGC and MM hits for SPs is small " << trackPattern.stgcSegment.size() + trackPattern.mmSegment.size());
76  }
77  else {
78  ATH_CHECK( calcMergedHit(trackPattern) );
79  }
80 
81 
82  return StatusCode::SUCCESS;
83 
84 }

◆ sysInitialize()

virtual StatusCode AthCommonDataStore< AthCommonMsg< AlgTool > >::sysInitialize ( )
overridevirtualinherited

Perform system initialization for an algorithm.

We override this to declare all the elements of handle key arrays at the end of initialization. See comments on updateVHKA.

Reimplemented in DerivationFramework::CfAthAlgTool, AthCheckedComponent< AthAlgTool >, AthCheckedComponent<::AthAlgTool >, and asg::AsgMetadataTool.

◆ sysStart()

virtual StatusCode AthCommonDataStore< AthCommonMsg< AlgTool > >::sysStart ( )
overridevirtualinherited

Handle START transition.

We override this in order to make sure that conditions handle keys can cache a pointer to the conditions container.

◆ updateVHKA()

void AthCommonDataStore< AthCommonMsg< AlgTool > >::updateVHKA ( Gaudi::Details::PropertyBase &  )
inlineinherited

Definition at line 308 of file AthCommonDataStore.h.

308  {
309  // debug() << "updateVHKA for property " << p.name() << " " << p.toString()
310  // << " size: " << m_vhka.size() << endmsg;
311  for (auto &a : m_vhka) {
312  std::vector<SG::VarHandleKey*> keys = a->keys();
313  for (auto k : keys) {
314  k->setOwner(this);
315  }
316  }
317  }

Member Data Documentation

◆ m_detStore

StoreGateSvc_t AthCommonDataStore< AthCommonMsg< AlgTool > >::m_detStore
privateinherited

Pointer to StoreGate (detector store by default)

Definition at line 393 of file AthCommonDataStore.h.

◆ m_evtStore

StoreGateSvc_t AthCommonDataStore< AthCommonMsg< AlgTool > >::m_evtStore
privateinherited

Pointer to StoreGate (event store by default)

Definition at line 390 of file AthCommonDataStore.h.

◆ m_varHandleArraysDeclared

bool AthCommonDataStore< AthCommonMsg< AlgTool > >::m_varHandleArraysDeclared
privateinherited

Definition at line 399 of file AthCommonDataStore.h.

◆ m_vhka

std::vector<SG::VarHandleKeyArray*> AthCommonDataStore< AthCommonMsg< AlgTool > >::m_vhka
privateinherited

Definition at line 398 of file AthCommonDataStore.h.

The documentation for this class was generated from the following files:
beamspotman.r
def r
Definition: beamspotman.py:676
TrigL2MuonSA::TrackPattern::stgcSegment
TrigL2MuonSA::StgcHits stgcSegment
Definition: TrackData.h:59
TrigL2MuonSA::TrackPattern::superPoints
TrigL2MuonSA::SuperPoint superPoints[s_NCHAMBER]
Definition: TrackData.h:60
xAOD::L2MuonParameters::Chamber
Chamber
Define chamber types and locations.
Definition: TrigMuonDefs.h:15
CaloCellPos2Ntuple.int
int
Definition: CaloCellPos2Ntuple.py:24
TrigL2MuonSA::StgcHitData::stationEta
int stationEta
Definition: StgcData.h:36
AthCommonDataStore::declareProperty
Gaudi::Details::PropertyBase & declareProperty(Gaudi::Property< T > &t)
Definition: AthCommonDataStore.h:145
TrigL2MuonSA::TrackPattern::mmSegment
TrigL2MuonSA::MmHits mmSegment
Definition: TrackData.h:58
TrigL2MuonSA::MmHitData::eta
double eta
Definition: MmData.h:28
TrigL2MuonSA::NswStationFitter::MakeSegment
StatusCode MakeSegment(TrigL2MuonSA::TrackPattern &trackPattern, TrigL2MuonSA::StgcHits &stgcHits) const
Definition: NswStationFitter.cxx:776
TrigL2MuonSA::MmHitData::phi
double phi
Definition: MmData.h:32
conifer::pow
constexpr int pow(int x)
Definition: conifer.h:20
TrigL2MuonSA::MmHitData::stationName
int stationName
Definition: MmData.h:37
M_PI
#define M_PI
Definition: ActiveFraction.h:11
bin
Definition: BinsDiffFromStripMedian.h:43
AthCommonDataStore< AthCommonMsg< AlgTool > >::m_evtStore
StoreGateSvc_t m_evtStore
Pointer to StoreGate (event store by default)
Definition: AthCommonDataStore.h:390
AthCommonDataStore< AthCommonMsg< AlgTool > >::m_vhka
std::vector< SG::VarHandleKeyArray * > m_vhka
Definition: AthCommonDataStore.h:398
xAOD::etaMax
etaMax
Definition: HIEventShape_v2.cxx:46
TrigL2MuonSA::NswStationFitter::findStgcHitsInSegment
StatusCode findStgcHitsInSegment(TrigL2MuonSA::StgcHits &stgcHits) const
Definition: NswStationFitter.cxx:269
read_hist_ntuple.t
t
Definition: read_hist_ntuple.py:5
drawFromPickle.cos
cos
Definition: drawFromPickle.py:36
TrigL2MuonSA::StgcHitData::phi
double phi
Definition: StgcData.h:33
x
#define x
TrigL2MuonSA::NswStationFitter::LinearFitWeight
void LinearFitWeight(std::vector< double > &x, std::vector< double > &y, std::vector< bool > &isStgc, double *slope, double *intercept, double *mse, double eta) const
Definition: NswStationFitter.cxx:844
TrigL2MuonSA::SuperPoint::Blin
float Blin
Definition: SuperPointData.h:106
TrigVSI::AlgConsts::nPhi
constexpr int nPhi
Default bin number of phi for vertex map.
Definition: Trigger/TrigTools/TrigVrtSecInclusive/TrigVrtSecInclusive/Constants.h:27
Trk::u
@ u
Enums for curvilinear frames.
Definition: ParamDefs.h:77
TrigL2MuonSA::SuperPoint::R
float R
Definition: SuperPointData.h:102
drawFromPickle.atan
atan
Definition: drawFromPickle.py:36
PUfitVar::maxEta
constexpr float maxEta
Definition: GepMETPufitAlg.cxx:13
SG::VarHandleKeyArray::setOwner
virtual void setOwner(IDataHandleHolder *o)=0
IDTPMcnv.htype
htype
Definition: IDTPMcnv.py:27
TRT::Hit::side
@ side
Definition: HitInfo.h:83
python.setupRTTAlg.size
int size
Definition: setupRTTAlg.py:39
TrigL2MuonSA::SuperPoint::Phim
float Phim
Definition: SuperPointData.h:104
TrigL2MuonSA::MmHitData::r
double r
Definition: MmData.h:33
python.utils.AtlRunQueryDQUtils.p
p
Definition: AtlRunQueryDQUtils.py:210
AthCommonDataStore
Definition: AthCommonDataStore.h:52
ATH_MSG_ERROR
#define ATH_MSG_ERROR(x)
Definition: AthMsgStreamMacros.h:33
TrigL2MuonSA::MmHitData::stationPhi
int stationPhi
Definition: MmData.h:36
lumiFormat.i
int i
Definition: lumiFormat.py:85
TrigL2MuonSA::StgcHitData::r
double r
Definition: StgcData.h:34
z
#define z
TrigL2MuonSA::StgcHitData
Definition: StgcData.h:14
TrigL2MuonSA::NswStationFitter::findSetOfStgcHitIds
void findSetOfStgcHitIds(TrigL2MuonSA::StgcHits &stgcHits, std::array< std::vector< int >, 8 > hitIdByLayer, std::vector< std::array< int, 8 >> &hitIdsCandidate) const
Definition: NswStationFitter.cxx:361
TrigL2MuonSA::NswStationFitter::findSetOfMmHitIds
void findSetOfMmHitIds(TrigL2MuonSA::MmHits &mmHits, std::array< std::vector< int >, 8 > hitIdByLayer, std::vector< std::array< int, 8 >> &hitIdsCandidate) const
Definition: NswStationFitter.cxx:1239
ATH_MSG_DEBUG
#define ATH_MSG_DEBUG(x)
Definition: AthMsgStreamMacros.h:29
test_pyathena.parent
parent
Definition: test_pyathena.py:15
ATH_CHECK
#define ATH_CHECK
Definition: AthCheckMacros.h:40
TrigL2MuonSA::NswStationFitter::selectMmHits
StatusCode selectMmHits(const TrigRoiDescriptor *p_roids, TrigL2MuonSA::MmHits &mmHits) const
Definition: NswStationFitter.cxx:134
drawFromPickle.tan
tan
Definition: drawFromPickle.py:36
AthCommonDataStore< AthCommonMsg< AlgTool > >::m_detStore
StoreGateSvc_t m_detStore
Pointer to StoreGate (detector store by default)
Definition: AthCommonDataStore.h:393
xAOD::double
double
Definition: CompositeParticle_v1.cxx:159
TrigL2MuonSA::StgcHitData::channelType
int channelType
Definition: StgcData.h:43
AthAlgTool::AthAlgTool
AthAlgTool()
Default constructor:
ZERO_LIMIT
const float ZERO_LIMIT
Definition: VP1TriggerHandleL2.cxx:37
SG::VarHandleKeyArray::renounce
virtual void renounce()=0
SG::HandleClassifier::type
std::conditional< std::is_base_of< SG::VarHandleKeyArray, T >::value, VarHandleKeyArrayType, type2 >::type type
Definition: HandleClassifier.h:54
dumpTgcDigiThreshold.isStrip
list isStrip
Definition: dumpTgcDigiThreshold.py:33
TrigL2MuonSA::StgcHitData::stationPhi
int stationPhi
Definition: StgcData.h:37
merge_scale_histograms.doc
string doc
Definition: merge_scale_histograms.py:9
name
std::string name
Definition: Control/AthContainers/Root/debug.cxx:221
python.subdetectors.mmg.ids
ids
Definition: mmg.py:8
TrigL2MuonSA::MmHitData
Definition: MmData.h:14
plotBeamSpotVxVal.bin
int bin
Definition: plotBeamSpotVxVal.py:83
TrigL2MuonSA::StgcHitData::stationName
int stationName
Definition: StgcData.h:38
TrigL2MuonSA::NswStationFitter::getNswResolution
void getNswResolution(double *stgcDeltaR, double *mmDeltaR, unsigned int size) const
Definition: NswStationFitter.cxx:914
TrigL2MuonSA::StgcHitData::layerNumber
unsigned int layerNumber
Definition: StgcData.h:42
TrigL2MuonSA::NswStationFitter::calcMergedHit
StatusCode calcMergedHit(TrigL2MuonSA::TrackPattern &trackPattern) const
Definition: NswStationFitter.cxx:924
a
TList * a
Definition: liststreamerinfos.cxx:10
TrigL2MuonSA::MmHits
std::vector< MmHitData > MmHits
Definition: MmData.h:47
LArCellBinning.etaMin
etaMin
Definition: LArCellBinning.py:84
y
#define y
h
ATH_MSG_WARNING
#define ATH_MSG_WARNING(x)
Definition: AthMsgStreamMacros.h:32
RoiDescriptor::phi
virtual double phi() const override final
Methods to retrieve data members.
Definition: RoiDescriptor.h:100
xAOD::L2MuonParameters::EndcapInner
@ EndcapInner
Inner station in the endcap spectrometer.
Definition: TrigMuonDefs.h:19
TrigL2MuonSA::SuperPoint::Npoint
int Npoint
Definition: SuperPointData.h:97
python.CaloScaleNoiseConfig.type
type
Definition: CaloScaleNoiseConfig.py:78
TrigL2MuonSA::StgcHits
std::vector< StgcHitData > StgcHits
Definition: StgcData.h:49
python.CaloCondTools.log
log
Definition: CaloCondTools.py:20
SG::VarHandleBase::vhKey
SG::VarHandleKey & vhKey()
Return a non-const reference to the HandleKey.
Definition: StoreGate/src/VarHandleBase.cxx:623
RoiDescriptor::eta
virtual double eta() const override final
Definition: RoiDescriptor.h:101
TrigL2MuonSA::NswStationFitter::findMmHitsInSegment
StatusCode findMmHitsInSegment(TrigL2MuonSA::MmHits &mmHits) const
Definition: NswStationFitter.cxx:1187
TrigL2MuonSA::SuperPoint::Z
float Z
Definition: SuperPointData.h:103
TrigL2MuonSA::StgcHitData::z
double z
Definition: StgcData.h:35
TrigL2MuonSA::SuperPoint
Definition: SuperPointData.h:74
python.Bindings.keys
keys
Definition: Control/AthenaPython/python/Bindings.py:798
calibdata.copy
bool copy
Definition: calibdata.py:27
TrigL2MuonSA::StgcHitData::eta
double eta
Definition: StgcData.h:29
drawFromPickle.sin
sin
Definition: drawFromPickle.py:36
TrigL2MuonSA::MmHitData::z
double z
Definition: MmData.h:34
TrigL2MuonSA::NswStationFitter::selectStgcHits
StatusCode selectStgcHits(const TrigRoiDescriptor *p_roids, TrigL2MuonSA::StgcHits &stgcHits) const
Definition: NswStationFitter.cxx:86
AthCommonDataStore::declareGaudiProperty
Gaudi::Details::PropertyBase & declareGaudiProperty(Gaudi::Property< T > &hndl, const SG::VarHandleKeyType &)
specialization for handling Gaudi::Property<SG::VarHandleKey>
Definition: AthCommonDataStore.h:156
TrigL2MuonSA::SuperPoint::Alin
float Alin
Definition: SuperPointData.h:105
readCCLHist.float
float
Definition: readCCLHist.py:83
fitman.k
k
Definition: fitman.py:528
TrigL2MuonSA::NswStationFitter::LinearFit
void LinearFit(std::vector< double > &x, std::vector< double > &y, double *slope, double *intercept, double *mse) const
Definition: NswStationFitter.cxx:806
TrigL2MuonSA::MmHitData::stationEta
int stationEta
Definition: MmData.h:35
Generated on Thu Nov 7 2024 22:27:47 for ATLAS Offline Software by   doxygen 1.8.18