1012 {
1014
1015
1016
1018 std::vector<Tracklet> myTracks = tracks;
1019 tracks.clear();
1020 for (const Tracklet &track : myTracks) {
1021 Identifier id1 =
track.getML1seg().mdtHitsOnTrack().at(0)->identify();
1022 Identifier id2 =
track.getML2seg().mdtHitsOnTrack().at(0)->identify();
1023 int nLayerML1 =
m_idHelperSvc->mdtIdHelper().tubeLayerMax(id1);
1024 int nLayerML2 =
m_idHelperSvc->mdtIdHelper().tubeLayerMax(id2);
1025 double ratio = (
double)(
track.mdtHitsOnTrack().size()) / (nLayerML1 + nLayerML2);
1026 if (ratio > 0.75) tracks.push_back(track);
1027 }
1028 }
1029
1030 std::vector<Tracklet> UniqueTracks;
1031 std::vector<unsigned int> AmbigTrks;
1032 for (unsigned int tk1 = 0; tk1 < tracks.size(); ++tk1) {
1033 int nShared = 0;
1034
1035 bool isResolved = false;
1036 for (unsigned int AmbigTrksIdx : AmbigTrks) {
1037 if (tk1 == AmbigTrksIdx) {
1038 isResolved = true;
1039 break;
1040 }
1041 }
1042 if (isResolved) continue;
1043 std::vector<Tracklet> AmbigTracks;
1044 AmbigTracks.push_back(tracks.at(tk1));
1045
1046 double Trk1ML1R = tracks.at(tk1).getML1seg().globalPosition().perp();
1047 double Trk1ML1Z = tracks.at(tk1).getML1seg().globalPosition().z();
1048 double Trk1ML2R = tracks.at(tk1).getML2seg().globalPosition().perp();
1049 double Trk1ML2Z = tracks.at(tk1).getML2seg().globalPosition().z();
1050
1051 Identifier tk1ID = tracks.at(tk1).muonIdentifier();
1052 bool tk1_isBarrel =
m_idHelperSvc->mdtIdHelper().isBarrel(tk1ID);
1053 bool tk1_isEndcap =
m_idHelperSvc->mdtIdHelper().isEndcap(tk1ID);
1054
1055
1056 for (unsigned int tk2 = (tk1 + 1); tk2 < tracks.size(); ++tk2) {
1057 if (tracks.at(tk1).mdtChamber() == tracks.at(tk2).mdtChamber() && tracks.at(tk1).mdtChPhi() == tracks.at(tk2).mdtChPhi() &&
1058 (tracks.at(tk1).mdtChEta()) * (tracks.at(tk2).mdtChEta()) > 0) {
1059
1060 for (unsigned int AmbigTrksIdx : AmbigTrks) {
1061 if (tk2 == AmbigTrksIdx) {
1062 isResolved = true;
1063 break;
1064 }
1065 }
1066 if (isResolved) continue;
1067
1068 double Trk2ML1R = tracks.at(tk2).getML1seg().globalPosition().perp();
1069 double Trk2ML1Z = tracks.at(tk2).getML1seg().globalPosition().z();
1070 double Trk2ML2R = tracks.at(tk2).getML2seg().globalPosition().perp();
1071 double Trk2ML2Z = tracks.at(tk2).getML2seg().globalPosition().z();
1072
1073
1074 double DistML1(1000), DistML2(1000);
1075 if (tk1_isBarrel) {
1076 DistML1 = std::abs(Trk1ML1Z - Trk2ML1Z);
1077 DistML2 = std::abs(Trk1ML2Z - Trk2ML2Z);
1078 } else if (tk1_isEndcap) {
1079 DistML1 = std::abs(Trk1ML1R - Trk2ML1R);
1080 DistML2 = std::abs(Trk1ML2R - Trk2ML2R);
1081 }
1082 if (DistML1 < 40 || DistML2 < 40) {
1083
1084 std::vector<const Muon::MdtPrepData*> mdts1 = tracks.at(tk1).mdtHitsOnTrack();
1085 std::vector<const Muon::MdtPrepData*> mdts2 = tracks.at(tk2).mdtHitsOnTrack();
1086 nShared = 0;
1087 for (const Muon::MdtPrepData *mdt1 : mdts1) {
1088 for (const Muon::MdtPrepData *mdt2 : mdts2) {
1089 if (mdt1->identify() == mdt2->identify()) {
1090 ++nShared;
1091 break;
1092 }
1093 }
1094 }
1095
1096 if (nShared <= 1) continue;
1097
1098 AmbigTracks.push_back(tracks.at(tk2));
1099 AmbigTrks.push_back(tk2);
1100 }
1101 }
1102 }
1103
1104 if (AmbigTracks.size() == 1) {
1105 UniqueTracks.push_back(tracks.at(tk1));
1106 continue;
1107 }
1108
1109
1110 if (tk1_isBarrel) {
1111 bool hasMomentum = tracks.at(tk1).charge() != 0;
1112 double aveX(0), aveY(0), aveZ(0), aveAlpha(0);
1113 double aveP(0), nAmbigP(0), TrkCharge(tracks.at(tk1).charge());
1114 bool allSameSign(true);
1115
1116 for (const Tracklet &AmbigTrack : AmbigTracks) {
1117 if (!hasMomentum) {
1118 aveX += AmbigTrack.globalPosition().x();
1119 aveY += AmbigTrack.globalPosition().y();
1120 aveZ += AmbigTrack.globalPosition().z();
1121 aveAlpha += AmbigTrack.getML1seg().alpha();
1122 } else {
1123
1124 if (std::abs(AmbigTrack.charge() - TrkCharge) > 0.1) allSameSign = false;
1125
1126 aveP += AmbigTrack.momentum().mag();
1127 ++nAmbigP;
1128 aveAlpha += AmbigTrack.alpha();
1129 aveX += AmbigTrack.globalPosition().x();
1130 aveY += AmbigTrack.globalPosition().y();
1131 aveZ += AmbigTrack.globalPosition().z();
1132 }
1133 }
1134 if (!hasMomentum) {
1135 aveX = aveX / (
double)AmbigTracks.size();
1136 aveY = aveY / (
double)AmbigTracks.size();
1137 aveZ = aveZ / (
double)AmbigTracks.size();
1139 aveAlpha = aveAlpha / (
double)AmbigTracks.size();
1140 double alphaErr = tracks.at(tk1).getML1seg().alphaError();
1141 double rErr = tracks.at(tk1).getML1seg().rError();
1142 double zErr = tracks.at(tk1).getML1seg().zError();
1143
1144 TrackletSegment aveSegML1(
m_idHelperSvc.get(), tracks.at(tk1).getML1seg().mdtHitsOnTrack(), gpos, aveAlpha, alphaErr, rErr, zErr, 0);
1145 double pT =
m_maxpTot * std::sin(aveSegML1.alpha());
1146 double pz =
m_maxpTot * std::cos(aveSegML1.alpha());
1148 pT * std::sin(aveSegML1.globalPosition().phi()),
1149 pz);
1151 matrix.setIdentity();
1152 matrix(0, 0) = std::pow(tracks.at(tk1).getML1seg().rError(),2);
1153 matrix(1, 1) = std::pow(tracks.at(tk1).getML1seg().zError(),2);
1154 matrix(2, 2) = std::pow(0.0000001,2);
1155 matrix(3, 3) = std::pow(tracks.at(tk1).getML1seg().alphaError(),2);
1157 Tracklet aveTrack(aveSegML1, momentum, matrix, 0);
1158 UniqueTracks.push_back(aveTrack);
1159 } else
if (allSameSign) {
1160 aveP = aveP / nAmbigP;
1161 double pT = aveP * std::sin(tracks.at(tk1).getML1seg().alpha());
1162 double pz = aveP * std::cos(tracks.at(tk1).getML1seg().alpha());
1164 pT * std::sin(tracks.at(tk1).globalPosition().phi()),
1165 pz);
1166 Tracklet MyTrack = tracks.at(tk1);
1168 MyTrack.
charge(tracks.at(tk1).charge());
1169 UniqueTracks.push_back(MyTrack);
1170 } else {
1171 aveX = aveX / (
double)AmbigTracks.size();
1172 aveY = aveY / (
double)AmbigTracks.size();
1173 aveZ = aveZ / (
double)AmbigTracks.size();
1175 aveAlpha = aveAlpha / (
double)AmbigTracks.size();
1176 double alphaErr = tracks.at(tk1).getML1seg().alphaError();
1177 double rErr = tracks.at(tk1).getML1seg().rError();
1178 double zErr = tracks.at(tk1).getML1seg().zError();
1179
1180 TrackletSegment aveSegML1(
m_idHelperSvc.get(), tracks.at(tk1).getML1seg().mdtHitsOnTrack(), gpos, aveAlpha, alphaErr, rErr, zErr, 0);
1181 double pT =
m_maxpTot * std::sin(aveSegML1.alpha());
1182 double pz =
m_maxpTot * std::cos(aveSegML1.alpha());
1184 pT * std::sin(aveSegML1.globalPosition().phi()),
1185 pz);
1187 matrix.setIdentity();
1188 matrix(0, 0) = std::pow(tracks.at(tk1).getML1seg().rError(),2);
1189 matrix(1, 1) = std::pow(tracks.at(tk1).getML1seg().zError(),2);
1190 matrix(2, 2) = std::pow(0.0000001,2);
1191 matrix(3, 3) = std::pow(tracks.at(tk1).getML1seg().alphaError(),2);
1193 Tracklet aveTrack(aveSegML1, momentum, matrix, 0);
1194 UniqueTracks.push_back(aveTrack);
1195 }
1196 }
1197
1198
1199 else
if (tk1_isEndcap) {
1200 std::vector<const Muon::MdtPrepData*> AllMdts;
1201 for (Tracklet const &AmbigTrack : AmbigTracks) {
1202 std::vector<const Muon::MdtPrepData*> mdts = AmbigTrack.mdtHitsOnTrack();
1203 std::vector<const Muon::MdtPrepData*> tmpAllMdt = AllMdts;
1204 for (const Muon::MdtPrepData *mdt : mdts) {
1205 bool isNewHit = true;
1206 for (const Muon::MdtPrepData *tmpmdt : tmpAllMdt) {
1207 if (mdt->identify() == tmpmdt->identify()) {
1208 isNewHit = false;
1209 break;
1210 }
1211 }
1212 if (isNewHit) AllMdts.push_back(mdt);
1213 }
1214 }
1215
1217 if (!MyECsegs.empty()) {
1218 TrackletSegment ECseg = MyECsegs.at(0);
1224 pz);
1226 matrix.setIdentity();
1227 matrix(0, 0) = std::pow(ECseg.rError(),2);
1228 matrix(1, 1) = std::pow(ECseg.zError(),2);
1229 matrix(2, 2) = std::pow(0.0000001,2);
1230 matrix(3, 3) = std::pow(ECseg.alphaError(),2);
1232 Tracklet MyCombTrack(MyECsegs.at(0), ECseg, momentum, matrix, 0);
1233 UniqueTracks.push_back(MyCombTrack);
1234 } else
1235 UniqueTracks.push_back(tracks.at(tk1));
1236 }
1237
1238 }
1239
1240 return UniqueTracks;
1241 }
void momentum(const Amg::Vector3D &p)
void charge(double charge)