635 {
636
637
638
639 MagField::AtlasFieldCache fieldCache;
640
644 return nullptr;
645 }
646
647 fieldCondObj->getInitializedCache (fieldCache);
648
649
650
651 SG::ReadHandle<InDet::PixelClusterContainer> pixcontainer(
m_pixcontainerkey, ctx);
652
653 if (!pixcontainer.isValid()) {
655 return nullptr;
656 }
657
659
660 SG::ReadHandle<InDet::SCT_ClusterContainer> sctcontainer(
m_sctcontainerkey, ctx);
661
662 if (!sctcontainer.isValid()) {
664 return nullptr;
665 }
666
668
669
670
672
673 std::vector<const Trk::PrepRawData*> assignedHits(nModules, nullptr);
674
675 std::vector<int> moduleStatus(nModules, 0);
676
677 std::vector<Identifier> seedIdents;
678 std::vector<const Trk::PrepRawData*> seedHits;
679
680 for(
const auto&
sp : seed) {
681 const Trk::PrepRawData* prd =
sp->clusterList().first;
683 seedHits.push_back(prd);
684 prd =
sp->clusterList().second;
685 if(prd == nullptr) continue;
687 seedHits.push_back(prd);
688 }
689
690
691
692 unsigned int seedSize = seedIdents.size();
693 int nUnassigned = seedSize;
694
695 int startModuleIdx = -1;
696
697 for(
int moduleIdx = 0;moduleIdx<
nModules;moduleIdx++) {
698
699 Identifier
ident = road.at(moduleIdx)->identify();
700
701 for(unsigned int clIdx=0;clIdx<seedSize;clIdx++) {
702
703 if(seedIdents[clIdx] != ident) continue;
704
705 assignedHits[moduleIdx] = seedHits[clIdx];
706 moduleStatus[moduleIdx] = 1;
707
708 startModuleIdx = moduleIdx;
709 --nUnassigned;
710 break;
711 }
712
713 if(nUnassigned == 0) break;
714 }
715
716 if(nUnassigned > 0) return nullptr;
717
718 std::unique_ptr<TrigFTF_ExtendedTrackState> initialState =
fitTheSeed(seed, fieldCache);
719
720 if(initialState == nullptr) return nullptr;
721
722 initialState->correctAngles();
723
724 TrigFTF_ExtendedTrackState& theState = *initialState;
725
726
727
728 std::vector<int> layerSequence[2];
729 std::unordered_map<int, std::vector<int> > layerMap[2];
730
733
734 for(int passIdx=0;passIdx<2;passIdx++) {
735
736 int start = passIdx==0 ? startModuleIdx + 1 : startModuleIdx;
738 int step = passIdx==0 ? 1 : -1;
739
740 for(
int moduleIdx = start;moduleIdx!=
end;moduleIdx+=
step) {
741 const InDetDD::SiDetectorElement* de = road.at(moduleIdx);
743 int l = de->
isPixel() ? vPixelL.at(h) : vStripL.at(h);
746 if(h % 2) {
748 }
749 }
750 auto it = layerMap[passIdx].find(l);
751 if(it != layerMap[passIdx].
end()) (*it).second.push_back(moduleIdx);
752 else {
753 layerSequence[passIdx].push_back(l);
754 std::vector<int>
v = {moduleIdx};
755 layerMap[passIdx].insert(std::make_pair(l,v));
756 }
757 }
758 }
759
760
761
762 for(int passIdx=0;passIdx<2;passIdx++) {
763
764 for(const auto& lkey : layerSequence[passIdx]) {
765
767 return nullptr;
768 }
769
770 const auto& lp = (*layerMap[passIdx].find(lkey));
771
772
773
774 std::vector<std::tuple<double, const Trk::PrepRawData*, int> > hitLinks;
775
776 for(const auto& moduleIdx : lp.second) {
777
778 if(moduleIdx == startModuleIdx) {
780 }
781
782 if(moduleStatus[moduleIdx] < 0) continue;
783
784 if (assignedHits[moduleIdx] != nullptr) {
785 hitLinks.emplace_back(std::make_tuple(-1.0,assignedHits[moduleIdx],moduleIdx));
786 continue;
787 }
788
789
790
791 const InDetDD::SiDetectorElement* de = road.at(moduleIdx);
792
793
794
796
797 bool noHits = false;
798
800 noHits = (p_pixcontainer == nullptr) || ((*p_pixcontainer).indexFindPtr(moduleHash) == nullptr);
801 } else {
802 noHits = (p_sctcontainer == nullptr) || ((*p_sctcontainer).indexFindPtr(moduleHash) == nullptr);
803 }
804
805 if (noHits) {
806 moduleStatus[moduleIdx] = -3;
807 continue;
808 }
809
810
811
812 const Trk::PlaneSurface* plane =
static_cast<const Trk::PlaneSurface*
>(&de->
surface());
813
814 double trackParams[2];
815
817 moduleStatus[moduleIdx] = -2;
818 continue;
819 }
820
821 const double bound_tol = 0.2;
822
824
826 moduleStatus[moduleIdx] = -2;
827 continue;
828 }
829
830
831
832 moduleStatus[moduleIdx] = -4;
833
835 findNearestHit(moduleIdx, (*p_pixcontainer).indexFindPtr(moduleHash), trackParams, hitLinks);
836 }
837 else {
839 }
840 }
841
842
843
844 if(hitLinks.empty()) {
846 continue;
847 }
848
849
850
851 std::sort(hitLinks.begin(), hitLinks.end());
852
853
854
855 for(const auto& hl : hitLinks) {
856
857 int moduleIdx =
get<2>(hl);
858
859 const Trk::PrepRawData* pPRD =
get<1>(hl);
860
861 const InDetDD::SiDetectorElement* de = road.at(moduleIdx);
862
863 const Trk::PlaneSurface* plane =
static_cast<const Trk::PlaneSurface*
>(&de->
surface());
864
865
866
868
869 if(rkCode!=0) {
870 moduleStatus[moduleIdx] = -2;
871 continue;
872 }
873
874 const Trk::PrepRawData* acceptedHit = nullptr;
875
877 const InDet::PixelCluster* pPixelHit = dynamic_cast<const InDet::PixelCluster*>(pPRD);
879 }
880 else {
881 const InDet::SCT_Cluster* pStripHit = dynamic_cast<const InDet::SCT_Cluster*>(pPRD);
883 }
884
885 if(acceptedHit == nullptr) {
886
887 }
888 else {
889 assignedHits[moduleIdx] = acceptedHit;
890 moduleStatus[moduleIdx] = de->
isPixel() ? 2 : 3;
891 }
892 }
893 }
894 }
895
896
897
900 double chi2tot = theState.
m_chi2;
901
903 return nullptr;
904 }
905
907 return nullptr;
908 }
909
911
912 if(rkCode !=0) {
913 return nullptr;
914 }
915
916 int ndoftot = theState.
m_ndof;
921
922 double z0 = theState.
m_Xk[1];
923 double d0 = theState.
m_Xk[0];
924
925 bool bad_cov = false;
926
928
929 for(
int i=0;
i<5;
i++) {
930 for(
int j=0;j<=
i;j++) {
931 double c = theState.
m_Gk[
i][j];
932 if (i == j && c < 0) {
933 bad_cov = true;
934 ATH_MSG_DEBUG(
"REGTEST: cov(" << i <<
"," << i <<
") =" << c <<
" < 0, reject track");
935 break;
936 }
937 cov.fillSymmetric(i,j, c);
938 }
939 }
940
941 if((ndoftot<0) || (fabs(pt)<100.0) || (std::isnan(pt)) || bad_cov) {
942 ATH_MSG_DEBUG(
"Track following failed - possibly floating point problem");
943 return nullptr;
944 }
945
946 Trk::PerigeeSurface perigeeSurface;
947
949
950 std::bitset<Trk::TrackStateOnSurface::NumberOfTrackStateOnSurfaceTypes> perType;
952
953 auto pParVec = std::make_unique<Trk::TrackStates>();
954
955 pParVec->reserve(50);
956 pParVec->push_back(new Trk::TrackStateOnSurface(nullptr, std::move(pPP), nullptr, perType));
957
958 std::bitset<Trk::TrackStateOnSurface::NumberOfTrackStateOnSurfaceTypes> rioType(0);
961
962 for(
const auto& ha : theState.
m_track) {
963
964 const Trk::PrepRawData* pPRD =
ha.m_pPRD;
965
966 if(pPRD == nullptr) continue;
967
968
969
971
972 const Trk::PlaneSurface* pPS =
dynamic_cast<const Trk::PlaneSurface*
>(&pPRD->
detectorElement()->surface());
973
974 if(pPS==nullptr) continue;
975
976 const InDet::SiCluster* pCL = dynamic_cast<const InDet::SiCluster*>(pPRD);
977
978 Trk::LocalParameters locPos = Trk::LocalParameters(pCL->
localPosition());
979
981
983
984 std::unique_ptr<Trk::MeasurementBase> pRIO{};
985
986 if(ndof == 2) {
987 const InDet::PixelCluster* pPixel = static_cast<const InDet::PixelCluster*>(pCL);
988 if(pPixel) {
989 pRIO = std::make_unique<InDet::PixelClusterOnTrack>(pPixel, std::move(locPos),
Amg::MatrixX(cov),
991 }
992 }
993 else {
994 const InDet::SCT_Cluster* pStrip = static_cast<const InDet::SCT_Cluster*>(pCL);
995 if(pStrip) {
996 pRIO = std::make_unique<InDet::SCT_ClusterOnTrack>(pStrip, std::move(locPos),
Amg::MatrixX(cov),
998 }
999 }
1000
1002
1004 for(
int i=0;
i<5;
i++) {
1005 for(
int j=0;j<=
i;j++) {
1006 pM.fillSymmetric(i,j,
ha.m_Ck[idx++]);
1007 }
1008 }
1009
1011
1012 Trk::FitQualityOnSurface FQ = Trk::FitQualityOnSurface(
ha.m_dchi2, ndof);
1013
1014 Trk::TrackStateOnSurface* pTSS = new Trk::TrackStateOnSurface(FQ, std::move(pRIO), std::move(pTP), nullptr, rioType);
1015
1016 pParVec->push_back(pTSS);
1017 }
1018
1019 auto pFQ = std::make_unique<Trk::FitQuality>(chi2tot, ndoftot);
1020
1021 Trk::TrackInfo
info{};
1022
1025
1026 Trk::Track* foundTrack = new Trk::Track(info, std::move(pParVec), std::move(pFQ));
1027
1028 return foundTrack;
1029}
Scalar theta() const
theta method
#define AmgSymMatrix(dim)
virtual DetectorShape shape() const
Shape of element.
virtual IdentifierHash identifyHash() const override final
identifier hash (inline)
SiIntersect inDetector(const Amg::Vector2D &localPosition, double phiTol, double etaTol) const
Test that it is in the active region.
Trk::Surface & surface()
Element Surface.
const Amg::Vector3D & globalPosition() const
return global position reference
bool gangedPixel() const
return the flag of this cluster containing a gangedPixel
Trk::PrepRawDataContainer< SCT_ClusterCollection > SCT_ClusterContainer
Trk::PrepRawDataContainer< PixelClusterCollection > PixelClusterContainer
virtual Surface::ChargedTrackParametersUniquePtr createUniqueTrackParameters(double l1, double l2, double phi, double theta, double qop, std::optional< AmgSymMatrix(5)> cov=std::nullopt) const override final
Use the Surface as a ParametersBase constructor, from local parameters - charged.
virtual Surface::ChargedTrackParametersUniquePtr createUniqueTrackParameters(double l1, double l2, double phi, double theta, double qop, std::optional< AmgSymMatrix(5)> cov=std::nullopt) const override final
Use the Surface as a ParametersBase constructor, from local parameters - charged.
const Amg::MatrixX & localCovariance() const
return const ref to the error matrix
@ Measurement
This is a measurement, and will at least contain a Trk::MeasurementBase.
@ Perigee
This represents a perigee, and so will contain a Perigee object only.
@ Scatterer
This represents a scattering point on the track, and so will contain TrackParameters and MaterialEffe...
virtual IdentifierHash identifyHash() const =0
Identifier hash.
virtual Identifier identify() const =0
Identifier.
std::vector< std::string > intersection(std::vector< std::string > &v1, std::vector< std::string > &v2)
T * get(TKey *tobj)
get a TObject* from a TKey* (why can't a TObject be a TKey?)
Eigen::Matrix< double, Eigen::Dynamic, Eigen::Dynamic > MatrixX
Dynamic Matrix - dynamic allocation.
Eigen::Matrix< double, 2, 1 > Vector2D
l
Printing final latex table to .tex output file.
void sort(typename DataModel_detail::iterator< DVL > beg, typename DataModel_detail::iterator< DVL > end)
Specialization of sort for DataVector/List.
std::list< TrigFTF_HitAssignment > m_track