diff --git a/InnerDetector/InDetRecAlgs/SiSpacePointFormation/SiSpacePointFormation/SiTrackerSpacePointFinder.h b/InnerDetector/InDetRecAlgs/SiSpacePointFormation/SiSpacePointFormation/SiTrackerSpacePointFinder.h index 1b3800fdbe47f1a7f89963634a649650164b59a1..88a38ac85da89705881cbb02d2f95d644d3deaa8 100755 --- a/InnerDetector/InDetRecAlgs/SiSpacePointFormation/SiSpacePointFormation/SiTrackerSpacePointFinder.h +++ b/InnerDetector/InDetRecAlgs/SiSpacePointFormation/SiSpacePointFormation/SiTrackerSpacePointFinder.h @@ -130,21 +130,7 @@ namespace InDet { (const SCT_ClusterCollection* next, const SiElementPropertiesTable* properties, const InDetDD::SiDetectorElementCollection* elements, - SpacePointCollection* spacepointCollection, SpacePointOverlapCollection* spacepointOverlapCollection, SPFCache&) const; - - void checkForSCT_Points - (const SCT_ClusterCollection* clusters1, - const IdentifierHash id2, - const InDetDD::SiDetectorElementCollection* elements, - double minDiff, double maxDiff, - SpacePointCollection* spacepointCollection, bool overlapColl, SpacePointOverlapCollection* spacepointOverlapCollection, SPFCache&) const; - - void checkForSCT_Points - (const SCT_ClusterCollection* clusters1, - const IdentifierHash id2, - const InDetDD::SiDetectorElementCollection* elements, - double min1, double max1, - double min2, double max2, SpacePointOverlapCollection* spacepointOverlapCollection, SPFCache&) const; + SpacePointCollection* spacepointCollection, SpacePointOverlapCollection* spacepointOverlapCollection, SPFCache&) const; //@} // data members @@ -253,7 +239,7 @@ namespace InDet { mutable std::atomic<int> m_sctCacheHits{0}; mutable std::atomic<int> m_pixCacheHits{0}; //@} - + }; } diff --git a/InnerDetector/InDetRecAlgs/SiSpacePointFormation/src/SiTrackerSpacePointFinder.cxx b/InnerDetector/InDetRecAlgs/SiSpacePointFormation/src/SiTrackerSpacePointFinder.cxx index 3f0f2723034e3d4ef7862fdf9cde2a027ebe3848..6692a055a91300ffa946b9f037dc432f119fb5fc 100755 --- a/InnerDetector/InDetRecAlgs/SiSpacePointFormation/src/SiTrackerSpacePointFinder.cxx +++ b/InnerDetector/InDetRecAlgs/SiSpacePointFormation/src/SiTrackerSpacePointFinder.cxx @@ -30,6 +30,7 @@ ATLAS Collaboration #include "StoreGate/ReadCondHandle.h" #include "AthenaMonitoringKernel/Monitored.h" + namespace InDet { //------------------------------------------------------------------------ @@ -52,7 +53,6 @@ namespace InDet { declareProperty("SpacePointCacheSCT", m_SpacePointCache_SCTKey=""); declareProperty("SpacePointCachePix", m_SpacePointCache_PixKey=""); - } //----------------------------------------------------------------------- @@ -115,6 +115,7 @@ StatusCode SiTrackerSpacePointFinder::initialize() if (!m_monTool.empty()) CHECK(m_monTool.retrieve()); ATH_MSG_INFO( "SiTrackerSpacePointFinder::initialized for package version " << PACKAGE_VERSION ); + return StatusCode::SUCCESS; } @@ -122,8 +123,6 @@ StatusCode SiTrackerSpacePointFinder::initialize() StatusCode SiTrackerSpacePointFinder::execute (const EventContext& ctx) const { - - ++m_numberOfEvents; const InDetDD::SiDetectorElementCollection* elements = nullptr; const SiElementPropertiesTable* properties = nullptr; @@ -232,12 +231,13 @@ StatusCode SiTrackerSpacePointFinder::execute (const EventContext& ctx) const auto spacepointCollection = std::make_unique<SpacePointCollection>(idHash); spacepointCollection->setIdentifier(elementID); - if ( colNext->size() != 0){ + if(!colNext->empty() && !elements->getDetectorElement(colNext->identifyHash())->isStereo()) { addSCT_SpacePoints(colNext, properties, elements, spacepointCollection.get(), spacepointoverlapCollection.ptr(), r_cache); } else { ATH_MSG_DEBUG( "Empty SCT cluster collection" ); } + size_t size = spacepointCollection->size(); nSCTspacePoints = size; if (size == 0){ @@ -338,9 +338,9 @@ StatusCode SiTrackerSpacePointFinder::execute (const EventContext& ctx) const m_sctCacheHits += sctCacheCount; m_pixCacheHits += pixCacheCount; } - - + return StatusCode::SUCCESS; + } //--------------------------------------------------------------------------- @@ -361,149 +361,147 @@ StatusCode SiTrackerSpacePointFinder::finalize() //-------------------------------------------------------------------------- -void SiTrackerSpacePointFinder:: -addSCT_SpacePoints(const SCT_ClusterCollection* next, - const SiElementPropertiesTable* properties, - const InDetDD::SiDetectorElementCollection* elements, - SpacePointCollection* spacepointCollection, SpacePointOverlapCollection* spacepointOverlapCollection, SPFCache &r_cache) const +void SiTrackerSpacePointFinder::addSCT_SpacePoints(const SCT_ClusterCollection* sctClusters, + const SiElementPropertiesTable* properties, + const InDetDD::SiDetectorElementCollection* elements, + SpacePointCollection* spacepointCollection, + SpacePointOverlapCollection* spacepointOverlapCollection, + SPFCache &r_cache) const { + // For each trigger element, first evaluate and collect the quantities you need to build the space points. + // The detector element array has capacity 6: + // [0] is the trigger element + // [1] is the opposite element + // [2]-[3] are the elements tested for eta overlaps + // [4]-[5] are the elements tested for phi overlaps + // For each element you save the correspoding cluster collections and the + // space point compatibility range as described below. + // + // For the opposite element and the ones tested for eta overlaps, you have to check + // if clusters are compatible with the local position of the trigger cluster + // requiring that the distance between the two clusters in phi is withing a specified range. + // - For the clusters on the opposite element: [-m_overlapLimitOpposite, m_overlapLimitOpposite] + // + // - For the eta overlapping clusters : you use m_overlapLimitEtaMin and m_overlapLimitEtaMax + // in different combination depending if you are checking a stereo module or not + // + // For each element, the extremes of these ranges are saved in the overlapExtents array + // which is later on used in SiSpacePointMakerTool::fillSCT_SpacePointCollection + // - overlapExtents[0], overlapExtents[1] are filled for the opposite element + // - overlapExtents[2], overlapExtents[3], overlapExtents[4], overlapExtents[5] are filled for the eta overlapping elements + // + // For the elements tested for phi overlaps, you have to check + // if clusters are compatible with the local position of the trigger cluster. + // This needs that the trigger cluster is at the edge of the trigger module: + // e.g. -hwidth < locX_trigger < -hwidth+m_overlapLimitPhi (or hwidth-m_overlapLimitPhi < locX_trigger < hwidth) + // and that the other cluster is on the compatible edge of its module: + // e.g. hwidth-m_overlapLimitPhi < locX_other < hwidth (or -hwidth < locX_other < -hwidth+m_overlapLimitPhi) + // + // For each element, the extremes of these ranges are saved in the overlapExtents array + // which is later on used in SiSpacePointMakerTool::fillSCT_SpacePointCollection + // - overlapExtents[6], overlapExtents[7], overlapExtents[10], overlapExtents[11] + // overlapExtents[8], overlapExtents[9], overlapExtents[12], overlapExtents[13] are filled for the phi overlapping elements + + enum NeighbourIndices{ThisOne, Opposite, PhiPlus, PhiMinus, EtaPlus, EtaMinus, nNeighbours}; + + std::array<const SCT_ClusterCollection*, nNeighbours> neighbourClusters{}; + std::array<const InDetDD::SiDetectorElement*, nNeighbours> neighbourElements{}; + std::array<double, 14> overlapExtents{}; + + // Get the detector element (and its Identifier) correspoding to the cluster collection + IdentifierHash triggerIdHash(sctClusters->identifyHash()); + const InDetDD::SiDetectorElement *triggerElement = elements->getDetectorElement(triggerIdHash); + Identifier thisID = triggerElement->identify(); + + // Retrieve the neighbours of the detector element + const std::vector<IdentifierHash>* others(properties->neighbours(triggerIdHash)); + if (others==0 || others->empty() ) return; + + // Save the current detector element and clusters + neighbourElements[0] = triggerElement; + neighbourClusters[0] = sctClusters; + int n = 0 ; + + // Default case: you test the opposite element and the overlapping in phi (total 3 elements) + int Nmax = 4 ; + + // In the barrel, test the eta overlaps as well (total 5 elements) + if (m_idHelper->is_barrel(thisID)) Nmax = 6; + + // You can remove all the overlaps if requrested. Here you test only the opposite element + if(!m_overlap) Nmax = 2; - // Do nothing unless this is a side 1 detector (strips of const phi). - IdentifierHash thisHash(next->identifyHash()); - - // if it is not the stereo side - const InDetDD::SiDetectorElement *element = elements->getDetectorElement(thisHash); - if (element && !(element->isStereo())){ - //if (m_idHelper->side(thisID)==1){ - // Space points are created from clusters for all possibly - // overlapping pairs in the given range if m_overlap=true. - // Otherwise just the opposite pair are processed. - // Pick up the five neighbours of detector, and check ranges - // for which overlap is possible. - // "check1" is used for opposite and eta overlaps. - // check2 for phi overlaps - - const std::vector<IdentifierHash>* - others(properties->neighbours(thisHash)); - if (others==0 || others->empty() ) return; - std::vector<IdentifierHash>::const_iterator otherHash = others->begin(); - - bool overlapColl = false; - // check opposite wafer - checkForSCT_Points(next, *otherHash, - elements, - -m_overlapLimitOpposite, +m_overlapLimitOpposite, - spacepointCollection,overlapColl,spacepointOverlapCollection, r_cache); - - if (! m_overlap) return; - - // check phi overlaps - overlapColl = true; - ++otherHash; - if (otherHash == others->end() ) return; - float hwidth(properties->halfWidth(thisHash)); - // half-width of wafer - - checkForSCT_Points(next, *otherHash, - elements, - -hwidth, -hwidth+m_overlapLimitPhi, - +hwidth-m_overlapLimitPhi, +hwidth,spacepointOverlapCollection, r_cache); - ++otherHash; - if (otherHash == others->end() ) return; - checkForSCT_Points(next, *otherHash, - elements, - +hwidth-m_overlapLimitPhi, +hwidth, - -hwidth, -hwidth+m_overlapLimitPhi,spacepointOverlapCollection, r_cache); - - // if barrel, check the eta overlaps too. - // In this case, action depends on whether layer is even or odd - // Also, whether we look at "lower in eta", which comes first, - // or "higher". - ++otherHash; - if (otherHash == others->end() ) return; - // get Identifier without hashing or so - Identifier thisID = element->identify(); - - //if (m_idHelper->barrel_ec(thisID)==0) - if (m_idHelper->is_barrel(thisID)) - { - //if (m_idHelper->layer_disk(thisID)==0 || m_idHelper->layer_disk(thisID)==2) - if (m_idHelper->layer_disk(thisID)%2 == 0) - { - checkForSCT_Points(next, *otherHash, - elements, - +m_overlapLimitEtaMin, - +m_overlapLimitEtaMax, - spacepointCollection,overlapColl,spacepointOverlapCollection, r_cache); - ++otherHash; - if (otherHash == others->end() ) return; - - checkForSCT_Points(next, *otherHash, - elements, - -m_overlapLimitEtaMax, - -m_overlapLimitEtaMin, - spacepointCollection,overlapColl,spacepointOverlapCollection, r_cache); - }else{ - checkForSCT_Points(next, *otherHash, - elements, - -m_overlapLimitEtaMax, - -m_overlapLimitEtaMin, - spacepointCollection,overlapColl,spacepointOverlapCollection, r_cache); - ++otherHash; - if (otherHash == others->end() ) return; - - - checkForSCT_Points(next, *otherHash, - elements, - +m_overlapLimitEtaMin, - +m_overlapLimitEtaMax, - spacepointCollection,overlapColl,spacepointOverlapCollection, r_cache); + float hwidth(properties->halfWidth(triggerIdHash)); + + // This flag is use to trigger if the search should be performed. + // In case there are no other detector elements and/or cluster collection + // this flag stays false. + bool search = false; + + // The order of the elements in others is such that you first get the opposite element, + // the overlapping in phi and then the overlapping in eta + // For this reason you need to re-order the indices, since the SiSpacePointMakerTool will process + // first the eta overlaps and then the phi ones + const std::array<size_t, nNeighbours> neigbourIndices{ThisOne, Opposite, EtaPlus, EtaMinus, PhiPlus, PhiMinus}; + + for (const auto& otherHash : *others) { + + if(++n==Nmax) break; + const SCT_ClusterCollection* otherClusters = r_cache.SCTCContainer->indexFindPtr (otherHash); + const InDetDD::SiDetectorElement* otherElement = elements ->getDetectorElement(otherHash); + + if(!otherElement || !otherClusters) continue; + + neighbourElements[neigbourIndices[n]] = otherElement; + neighbourClusters[neigbourIndices[n]] = otherClusters; + search = true ; + + switch (n) { + case 1: { + overlapExtents[ 0] = -m_overlapLimitOpposite; + overlapExtents[ 1] = m_overlapLimitOpposite; + break; + } + case 2: { + overlapExtents[ 6] =-hwidth; + overlapExtents[ 7] =-hwidth+m_overlapLimitPhi; + overlapExtents[ 8] = hwidth-m_overlapLimitPhi; + overlapExtents[ 9] = hwidth; + break; + } + case 3: { + overlapExtents[10] = hwidth-m_overlapLimitPhi; + overlapExtents[11] = hwidth; + overlapExtents[12] =-hwidth; + overlapExtents[13] =-hwidth+m_overlapLimitPhi; + break; + } + case 4: { + if (m_idHelper->layer_disk(thisID)%2 == 0) { + overlapExtents[ 2] = m_overlapLimitEtaMin; + overlapExtents[ 3] = m_overlapLimitEtaMax; + } else { + overlapExtents[ 2] =-m_overlapLimitEtaMax; + overlapExtents[ 3] =-m_overlapLimitEtaMin; + } + break; + } + default: { + if (m_idHelper->layer_disk(thisID)%2 == 0) { + overlapExtents[ 4] = -m_overlapLimitEtaMax; + overlapExtents[ 5] = -m_overlapLimitEtaMin; + } else { + overlapExtents[ 4] = m_overlapLimitEtaMin; + overlapExtents[ 5] = m_overlapLimitEtaMax; + } + break; } } } -} - -//-------------------------------------------------------------------------- - -void SiTrackerSpacePointFinder:: -checkForSCT_Points(const SCT_ClusterCollection* clusters1, - const IdentifierHash id2, - const InDetDD::SiDetectorElementCollection* elements, - double min, double max, - SpacePointCollection* spacepointCollection, bool overlapColl, SpacePointOverlapCollection* spacepointOverlapCollection, SPFCache &r_cache) const -{ - - - // Get the cluster collections for these two detectors. - // Require that (xPhi2 - xPhi1) lie in the range specified. - // Used for opposite and eta modules - - //indexFindPtr is faster in the MT implementation - const SCT_ClusterCollection * clusters2 = r_cache.SCTCContainer->indexFindPtr(id2); - if (clusters2==nullptr) return; - - if (!overlapColl) { - m_SiSpacePointMakerTool->fillSCT_SpacePointCollection(clusters1, clusters2, min, max, m_allClusters, r_cache.vertex, elements, spacepointCollection); - } - else { - m_SiSpacePointMakerTool->fillSCT_SpacePointEtaOverlapCollection(clusters1, clusters2, min, max, m_allClusters, r_cache.vertex, elements, spacepointOverlapCollection); + + if(search) { + m_SiSpacePointMakerTool->fillSCT_SpacePointCollection(neighbourElements,neighbourClusters,overlapExtents,m_allClusters,r_cache.vertex,spacepointCollection,spacepointOverlapCollection); } } -//-------------------------------------------------------------------------- -void SiTrackerSpacePointFinder:: -checkForSCT_Points(const SCT_ClusterCollection* clusters1, - const IdentifierHash id2, - const InDetDD::SiDetectorElementCollection* elements, - double min1, double max1, double min2, double max2, SpacePointOverlapCollection* spacepointOverlapCollection, SPFCache &r_cache) const -{ - - // get the cluster collections for these two detectors. Clus1 must lie - // within min1 and max1 and clus between min2 and max2. Used for phi - // overlaps. - const SCT_ClusterCollection * clusters2 (r_cache.SCTCContainer->indexFindPtr(id2)); - if (clusters2==nullptr) return; - - m_SiSpacePointMakerTool->fillSCT_SpacePointPhiOverlapCollection(clusters1, clusters2, min1, max1, min2, max2, m_allClusters, r_cache.vertex, elements, spacepointOverlapCollection); -} } //namespace diff --git a/InnerDetector/InDetRecTools/SiSpacePointTool/SiSpacePointTool/SCTinformation.h b/InnerDetector/InDetRecTools/SiSpacePointTool/SiSpacePointTool/SCTinformation.h index c77e4ae7a76c7faf550644bfb5755986512d443f..f3e00981d1430831ffcbfb973d3eb2b62e495e10 100644 --- a/InnerDetector/InDetRecTools/SiSpacePointTool/SiSpacePointTool/SCTinformation.h +++ b/InnerDetector/InDetRecTools/SiSpacePointTool/SiSpacePointTool/SCTinformation.h @@ -1,7 +1,7 @@ // -*- C++ -*- /* - Copyright (C) 2002-2019 CERN for the benefit of the ATLAS collaboration + Copyright (C) 2002-2020 CERN for the benefit of the ATLAS collaboration */ @@ -27,85 +27,63 @@ namespace InDet { ///////////////////////////////////////////////////////////////////////////////// public: - - SCTinformation(const InDet::SiCluster* CL, const Amg::Vector3D& a, const Amg::Vector3D& b, const Amg::Vector3D& vec); + SCTinformation(); + SCTinformation(const InDet::SiCluster* cluster, + const Amg::Vector3D& strip_start, + const Amg::Vector3D& strip_end, + const Amg::Vector3D& vec, + const double& locx); SCTinformation(const SCTinformation&) = default; virtual ~SCTinformation() = default; SCTinformation& operator = (const SCTinformation&) = default; - const InDet::SiCluster* cluster() const; - const Amg::Vector3D& a() const; - const Amg::Vector3D& q() const; - const Amg::Vector3D& s() const; - const Amg::Vector3D& qs() const; - double qm() const; + const InDet::SiCluster* cluster() const {return m_cluster;} + const Amg::Vector3D& strip_center () const {return m_center ;} + const Amg::Vector3D& strip_direction () const {return m_sdir ;} + const Amg::Vector3D& traj_direction () const {return m_tdir ;} + const Amg::Vector3D& normal() const {return m_normal;} + const double& oneOverStrip() const {return m_oneOverStrip;} + const double& locX() const {return m_locX;} double step(SCTinformation&) const; - Amg::Vector3D position(const double s) const; + Amg::Vector3D position(const double& s) const; + void set(const InDet::SiCluster* CL, + const Amg::Vector3D& strip_start, + const Amg::Vector3D& strip_end, + const Amg::Vector3D& vec, + const double& locx); + private: - const InDet::SiCluster* m_CL; - Amg::Vector3D m_A; - Amg::Vector3D m_Q; - Amg::Vector3D m_S; - Amg::Vector3D m_QS; - double m_QM; + const InDet::SiCluster* m_cluster; // SCT cluster + Amg::Vector3D m_center; // Center of strip + Amg::Vector3D m_sdir; // Direction of strip (strip_start-strip_end) + Amg::Vector3D m_tdir; // Direction of trajectory (strip_start+strip_end-2*vertexVec) + Amg::Vector3D m_normal; // Normal to strip diretion and trjectory direction plane + double m_oneOverStrip; // Invert length of the strip + double m_locX; // cluster local X }; ///////////////////////////////////////////////////////////////////////////////// // Inline methods ///////////////////////////////////////////////////////////////////////////////// - + inline SCTinformation::SCTinformation() : m_cluster(nullptr) {} + ///////////////////////////////////////////////////////////////////////////////// // Constructor with initialisation ///////////////////////////////////////////////////////////////////////////////// - inline SCTinformation::SCTinformation - (const InDet::SiCluster* CL, const Amg::Vector3D& a, const Amg::Vector3D& b, const Amg::Vector3D& vec) - { - m_CL = CL ; - m_A = 0.5*(a+b) ; // Center of strip - m_Q = a-b ; // Direction of strip a-b - m_S = 2.*(m_A-vec) ; // Direction of trajectory a+b-2*vertexVec - m_QS = m_Q.cross(m_S); // Normal to strip diretion and trjectory direction plane - m_QM = 1./m_Q.mag() ; // Invert length of the strip + inline SCTinformation::SCTinformation (const InDet::SiCluster* cluster, + const Amg::Vector3D& strip_start, + const Amg::Vector3D& strip_end, + const Amg::Vector3D& vec, + const double& locx) { + this->set(cluster, strip_start, strip_end, vec, locx); } - ///////////////////////////////////////////////////////////////////////////////// - // Getter methods - ///////////////////////////////////////////////////////////////////////////////// - inline const InDet::SiCluster* SCTinformation::cluster() const - { - return m_CL; - } - - inline const Amg::Vector3D& SCTinformation::a() const - { - return m_A; - } - - inline const Amg::Vector3D& SCTinformation::q() const - { - return m_Q; - } - - inline const Amg::Vector3D& SCTinformation::s() const - { - return m_S; - } - - inline const Amg::Vector3D& SCTinformation::qs() const - { - return m_QS; - } - - inline double SCTinformation::qm() const - { - return m_QM; - } ///////////////////////////////////////////////////////////////////////////////// // Step along strip @@ -113,16 +91,38 @@ namespace InDet { inline double SCTinformation::step(SCTinformation& In) const { - double a = m_Q.dot(In.m_QS); - if (a!=0.) return (-m_S.dot(In.m_QS)/a); + // this function returns the ratio between the components of the + // In normal vector over the strip direction and trajectory direction. + + // first evaluate the denominator: + // a is component of the In normal vector on the strip direction. + double a = m_sdir.dot(In.m_normal); + + // a = 0 if the In.m_normal and m_sdir are perpendicular. + // If so, return a very big number. + // Otherwise, evaluate the numerator and return + if (a!=0.) return (-m_tdir.dot(In.m_normal)/a); else return 10000.; } - inline Amg::Vector3D SCTinformation::position(const double s) const + inline Amg::Vector3D SCTinformation::position(const double& s) const { - return (m_A+(0.5*s)*m_Q); + return (m_center+(0.5*s)*m_sdir); } - + inline void SCTinformation::set (const InDet::SiCluster* cluster, + const Amg::Vector3D& strip_start, + const Amg::Vector3D& strip_end, + const Amg::Vector3D& vec, + const double& locx) { + m_cluster = cluster; + m_center = 0.5*(strip_start+strip_end); + m_sdir = strip_start-strip_end; + m_tdir = 2.*(m_center-vec); + m_normal = m_sdir.cross(m_tdir); + m_oneOverStrip = 1./m_sdir.mag(); + m_locX = locx; + } + } // end of name space #endif // SCTinformation_h diff --git a/InnerDetector/InDetRecTools/SiSpacePointTool/SiSpacePointTool/SiSpacePointMakerTool.h b/InnerDetector/InDetRecTools/SiSpacePointTool/SiSpacePointTool/SiSpacePointMakerTool.h index 76b5ef5de79d615e947b66c800b64ef97ea1632b..7dc072a7da2374c69ebe18778f65dfa8ccbec320 100644 --- a/InnerDetector/InDetRecTools/SiSpacePointTool/SiSpacePointTool/SiSpacePointMakerTool.h +++ b/InnerDetector/InDetRecTools/SiSpacePointTool/SiSpacePointTool/SiSpacePointMakerTool.h @@ -29,6 +29,9 @@ namespace InDet { } namespace InDet { + + constexpr size_t nNeighbours{6}; + /** * @class SiSpacePointMakerTool * Used by SiTrackerSpacePointFinder. @@ -62,6 +65,8 @@ namespace InDet { Trk::SpacePoint* makeSCT_SpacePoint(const InDet::SiCluster& cluster1, const InDet::SiCluster& cluster2, const Amg::Vector3D& vertexVec, const InDetDD::SiDetectorElement* element1, const InDetDD::SiDetectorElement* element2, double stripLengthGapTolerance) const; + + Trk::SpacePoint* makeSCT_SpacePoint(InDet::SCTinformation&,InDet::SCTinformation&,IdentifierHash,IdentifierHash,double,double) const; /// Convert clusters to space points: SCT_Clusters -> SCT_SpacePoints void fillSCT_SpacePointCollection(const InDet::SCT_ClusterCollection* clusters1, @@ -69,6 +74,11 @@ namespace InDet { const Amg::Vector3D& vertexVec, const InDetDD::SiDetectorElementCollection* elements, SpacePointCollection* spacepointCollection) const; + void fillSCT_SpacePointCollection(std::array<const InDetDD::SiDetectorElement*, nNeighbours>&, + std::array<const SCT_ClusterCollection*, nNeighbours>&, + std::array<double, 14>&,bool,const Amg::Vector3D&, + SpacePointCollection*,SpacePointOverlapCollection*) const; + /// Convert clusters to space points: PixelClusters -> PixelSpacePoints void fillPixelSpacePointCollection(const InDet::PixelClusterCollection* clusters, SpacePointCollection* spacepointCollection) const; @@ -86,6 +96,8 @@ namespace InDet { const InDetDD::SiDetectorElementCollection* elements, SpacePointOverlapCollection* spacepointOverlapCollection) const; + + private: /// @name Cut parameters //@{ @@ -138,19 +150,22 @@ namespace InDet { /// Guarded by m_mutex in const methods. mutable SG::SlotSpecificObj<CacheEntry> m_cache ATLAS_THREAD_SAFE; + /// update range accordingly to the gap between the stereo modules + void updateRange(const InDetDD::SiDetectorElement* element1, + const InDetDD::SiDetectorElement* element2, + double& stripLengthGapTolerance, double& min, double& max) const; + /// Get stripLengthGapTolerance and return offset value for two SiDetectorElement's - double offset(const InDetDD::SiDetectorElement* element1, const InDetDD::SiDetectorElement* element2, double& stripLengthGapTolerance) const; - - /// Get stripLengthGapTolerance for two SiDetectorElement's - void offset(double& stripLengthGapTolerance, const InDetDD::SiDetectorElement* element1, const InDetDD::SiDetectorElement* element2) const; + double offset(const InDetDD::SiDetectorElement* element1, + const InDetDD::SiDetectorElement* element2, + double& stripLengthGapTolerance) const; /// Not implemented yet - bool fillSCT_Information(const InDet::SCT_ClusterCollection* clusters1, const InDet::SCT_ClusterCollection* clusters2, + bool fillSCT_Information(const InDet::SCT_ClusterCollection* clusters1, + const InDet::SCT_ClusterCollection* clusters2, const Amg::Vector3D& vertexVec, const InDetDD::SiDetectorElementCollection* elements) const; - - /// Convert clusters to space points using CacheEntry - void makeSCT_SpacePoints(const double stripLengthGapTolerance) const; + }; } diff --git a/InnerDetector/InDetRecTools/SiSpacePointTool/src/SiSpacePointMakerTool.cxx b/InnerDetector/InDetRecTools/SiSpacePointTool/src/SiSpacePointMakerTool.cxx index 7f0680ba638fdb62cc4111578c04d76ef7e0579f..1c685db9bb13f54467eb86ae63db0b4dbed8c366 100644 --- a/InnerDetector/InDetRecTools/SiSpacePointTool/src/SiSpacePointMakerTool.cxx +++ b/InnerDetector/InDetRecTools/SiSpacePointTool/src/SiSpacePointMakerTool.cxx @@ -226,9 +226,7 @@ namespace InDet { } if (m_SCTgapParameter != 0.) { - double dm = offset(element1, element2, stripLengthGapTolerance); - min -= dm; - max += dm; + updateRange(element1, element2, stripLengthGapTolerance, min, max); } for (; clusters2Next!=clusters2Finish; ++clusters2Next){ @@ -257,7 +255,7 @@ namespace InDet { InDet::PixelClusterCollection::const_iterator clusStart = clusters->begin(); InDet::PixelClusterCollection::const_iterator clusFinish = clusters->end(); - + if ((*clusStart)->detectorElement()) { IdentifierHash idHash = clusters->identifyHash(); const Amg::Transform3D& T = (*clusStart)->detectorElement()->surface().transform(); @@ -337,9 +335,7 @@ namespace InDet { break; } if (m_SCTgapParameter != 0.) { - double dm = offset(element1, element2, stripLengthGapTolerance); - min -= dm; - max += dm; + updateRange(element1, element2, stripLengthGapTolerance, min, max); } for (; clusters2Next!=clusters2Finish; ++clusters2Next){ @@ -409,9 +405,7 @@ namespace InDet { } if (m_SCTgapParameter != 0.) { - double dm = offset(element1, element2, stripLengthGapTolerance); - min2 -= dm; - max2 += dm; + updateRange(element1, element2, stripLengthGapTolerance, min2, max2); } for (; clusters2Next!=clusters2Finish; ++clusters2Next) { @@ -450,120 +444,301 @@ namespace InDet { if (fabs(T1(2,2)) > 0.7) d*=(r/fabs(T1(2,3))); // endcap d = d*R/Z stripLengthGapTolerance = d; + return dm; } + + /////////////////////////////////////////////////////////////////// + // Updating the range used to search for overlaps + /////////////////////////////////////////////////////////////////// + + void SiSpacePointMakerTool::updateRange(const InDetDD::SiDetectorElement* element1, + const InDetDD::SiDetectorElement* element2, + double& stripLengthGapTolerance, + double& min, double& max) const { + double dm = offset(element1, element2, stripLengthGapTolerance); + min -= dm; + max += dm; + } + - void SiSpacePointMakerTool::offset(double& stripLengthGapTolerance, - const InDetDD::SiDetectorElement* element1, - const InDetDD::SiDetectorElement* element2) const { - const Amg::Transform3D& T1 = element1->transform(); - const Amg::Transform3D& T2 = element2->transform(); - Amg::Vector3D C = element1->center() ; + /////////////////////////////////////////////////////////////////// + // New methods for sct space points production + /////////////////////////////////////////////////////////////////// - double r = sqrt(C[0]*C[0]+C[1]*C[1]) ; - double x12 = T1(0,0)*T2(0,0) + T1(1,0)*T2(1,0) + T1(2,0)*T2(2,0) ; - double s = (T1(0,3)-T2(0,3))*T1(0,2) + (T1(1,3)-T2(1,3))*T1(1,2) + (T1(2,3)-T2(2,3))*T1(2,2); + void SiSpacePointMakerTool::fillSCT_SpacePointCollection (std::array<const InDetDD::SiDetectorElement*, nNeighbours>& elements, + std::array<const SCT_ClusterCollection*, nNeighbours>& clusters, + std::array<double, 14>& overlapExtents, + bool allClusters, + const Amg::Vector3D& vertexVec, + SpacePointCollection* spacepointCollection, + SpacePointOverlapCollection* spacepointoverlapCollection) const + { + + // This function is called once all the needed quantities are collected. + // It is used to build space points checking the compatibility of clusters on pairs of detector elements. + // Detector elements and cluster collections are elements and clusters, respectively. + // [0] is the trigger element + // [1] is the opposite element + // [2]-[3] are the elements tested for eta overlaps + // [4]-[5] are the elements tested for phi overlaps + // + // To build space points: + // - For the opposite element and the ones tested for eta overlaps, you have to check + // if clusters are compatible with the local position of the trigger cluster + // requiring that the distance between the two clusters in phi is withing a specified range. + // - overlapExtents[0], overlapExtents[1] are filled for the opposite element + // - overlapExtents[2], overlapExtents[3], overlapExtents[4], overlapExtents[5] are filled for the eta overlapping elements + // - For the elements tested for phi overlaps, you have to check + // if clusters are compatible with the local position of the trigger cluster. + // This needs that the trigger cluster is at the edge of the trigger module + // and that the other cluster is on the compatible edge of its module + // - overlapExtents[6], overlapExtents[7], overlapExtents[10], overlapExtents[11] + // overlapExtents[8], overlapExtents[9], overlapExtents[12], overlapExtents[13] are filled for the phi overlapping elements + + constexpr int otherSideIndex{1}; + constexpr int maxEtaIndex{3}; + std::array<int, nNeighbours-1> elementIndex{0}; + int nElements = 0; + + // For the nNeighbours sides, fill elementIndex with the indices of the existing elements. + // Same the number of elements in nElements to loop on the later on + for(int n=1; n!=nNeighbours; ++n) { + if(elements[n]) { + elementIndex[nElements++] = n; + } + } + // return if all detector elements are nullptr + if(!nElements) return; - double dm = (m_SCTgapParameter*r)*fabs(s*x12); + // trigger element and clusters + const InDetDD::SiDetectorElement* element = elements[0]; + IdentifierHash Id = clusters[0]->identifyHash(); - double d = 0.; - if (element1->design().shape() == InDetDD::Trapezoid || element1->design().shape() == InDetDD::Annulus) { - d = dm*(1./0.04); - } else { - d = dm/sqrt((1.-x12)*(1.+x12)); + std::vector<SCTinformation> sctInfos; + sctInfos.reserve(clusters[0]->size()); + + // loop on all clusters on the trigger detector element and save the related information + for (const auto& cluster : *clusters[0]) { + Amg::Vector2D locpos = cluster->localPosition(); + std::pair<Amg::Vector3D, Amg::Vector3D > ends(element->endsOfStrip(InDetDD::SiLocalPosition(locpos.y(),locpos.x(),0.))); + InDet::SCTinformation sct(cluster,ends.first, ends.second, vertexVec,locpos.x()); + sctInfos.push_back(sct); } - if (fabs(T1(2,2)) > 0.7) d *= (r/fabs(T1(2,3))); // endcap d = d*R/Z - stripLengthGapTolerance = d; - } + double limit = 1. + m_stripLengthTolerance; + double slimit = 0. ; + std::vector<Trk::SpacePoint*> tmpSpacePoints; + tmpSpacePoints.reserve(sctInfos.size()); - /////////////////////////////////////////////////////////////////// - // Compare SCT strips and space points production - /////////////////////////////////////////////////////////////////// + if(!allClusters) { - void SiSpacePointMakerTool::makeSCT_SpacePoints(const double stripLengthGapTolerance) const { - // Find intersection of a line through a cluster on one sct detector and - // a line through a cluster on its stereo pair. Return zero if lines - // don't intersect. + // Start processing the opposite side and the eta overlapping elements + int n = 0; + for(; n < nElements; ++n) { + + int currentIndex = elementIndex[n]; + + if(currentIndex > maxEtaIndex) break; - // A general point on the line joining point a to point b is - // x, where 2*x=(1+m)*a + (1-m)*b. Similarly for 2*y=(1+n)*c + (1-n)*d. - // Suppose that v is the vertex, which we take as (0,0,0); this could - // an input parameter later, if required. Requiring that the two 'general - // points' lie on a straight through v means that the vector x-v is a - // multiple of y-v. This condition fixes the parameters m and n. - // We then return the 'space-point' x, supposed to be the phi-layer. - // We require that -1<m<1, otherwise x lies - // outside the segment a to b; and similarly for n. + // get the detector element and the IdentifierHash + const InDetDD::SiDetectorElement* currentElement = elements[currentIndex]; + IdentifierHash currentId = clusters[currentIndex]->identifyHash(); - const EventContext& ctx{Gaudi::Hive::currentContext()}; - std::lock_guard<std::mutex> lock{m_mutex}; - CacheEntry* ent{m_cache.get(ctx)}; - if (ent->m_evt!=ctx.evt()) { // New event in this slot - ent->clear(); - ent->m_evt = ctx.evt(); - } else { - if (ent->m_tmpSpacePoints.size()) { - for (Trk::SpacePoint* sp : ent->m_tmpSpacePoints) { - delete sp; + // retrieve the range + double min = overlapExtents[currentIndex*2-2]; + double max = overlapExtents[currentIndex*2-1]; + + if (m_SCTgapParameter != 0.) { + updateRange(element, currentElement, slimit, min, max); + } + + InDet::SCTinformation sctInfo; + + for (const auto& cluster : *clusters[currentIndex]) { + + bool processed = false; + const Amg::Vector2D& locpos = cluster->localPosition(); + double lx1 = locpos.x(); + + for(auto& sct : sctInfos) { + + double diff = lx1-sct.locX(); + + if(diff < min || diff > max) continue; + + if (not processed) { + processed = true; + std::pair<Amg::Vector3D, Amg::Vector3D > ends(currentElement->endsOfStrip(InDetDD::SiLocalPosition(locpos.y(),locpos.x(),0.))); + sctInfo.set(cluster,ends.first, ends.second, vertexVec,lx1); + } + + Trk::SpacePoint* sp = makeSCT_SpacePoint(sct,sctInfo,Id,currentId,limit,slimit); + if(sp) tmpSpacePoints.push_back(sp); + } + } + // If you are processing the opposite element, save the space points into + // the spacepointCollection and clear the temporary collection + if( currentIndex==otherSideIndex && !tmpSpacePoints.empty() ) { + spacepointCollection->reserve(tmpSpacePoints.size()+spacepointCollection->size()); + for (Trk::SpacePoint* sp: tmpSpacePoints) { + spacepointCollection->push_back(sp); + } + tmpSpacePoints.clear(); + } + } + + // process the phi overlapping elements + // if possible n starts from 4 + for(; n < nElements; ++n) { + + int currentIndex = elementIndex[n]; + const InDetDD::SiDetectorElement* currentElement = elements[currentIndex]; + + double min = overlapExtents[4*currentIndex-10]; + double max = overlapExtents[4*currentIndex- 9]; + + if (m_SCTgapParameter != 0.) { + updateRange(element, currentElement, slimit, min, max); + } + + std::vector<SCTinformation*> sctPhiInfos; + sctPhiInfos.reserve(sctInfos.size()); + + for(auto& sct : sctInfos) { + double lx0 = sct.locX(); + if (min <= lx0 && lx0 <= max) { + sctPhiInfos.push_back(&sct); + } + } + // continue if you have no cluster from the phi overlapping region of the trigger element + if(sctPhiInfos.empty()) continue; + + IdentifierHash currentId = clusters[currentIndex]->identifyHash(); + + min = overlapExtents[4*currentIndex-8]; + max = overlapExtents[4*currentIndex-7]; + + if (m_SCTgapParameter != 0.) { + updateRange(element, currentElement, slimit, min, max); + } + + for (const auto& cluster : *clusters[currentIndex]) { + + const Amg::Vector2D& locpos = cluster->localPosition(); + double lx1 = locpos.x(); + if(lx1 < min || lx1 > max ) continue; + + std::pair<Amg::Vector3D, Amg::Vector3D > ends(currentElement->endsOfStrip(InDetDD::SiLocalPosition(locpos.y(),locpos.x(),0.))); + InDet::SCTinformation sctInfo(cluster,ends.first, ends.second, vertexVec,lx1); + + for(auto& sct : sctPhiInfos) { + Trk::SpacePoint* sp = makeSCT_SpacePoint(*sct,sctInfo,Id,currentId,limit,slimit); + if(sp) tmpSpacePoints.push_back(sp); + } + } + } + + // fill the space point collection for eta/phi overlapping clusters + if(!tmpSpacePoints.empty()) { + spacepointoverlapCollection->reserve(tmpSpacePoints.size()+spacepointoverlapCollection->size()); + for (Trk::SpacePoint* sp: tmpSpacePoints) { + spacepointoverlapCollection->push_back(sp); } } - ent->m_tmpSpacePoints.clear(); + return; } - std::vector<SCTinformation>::iterator I = ent->m_SCT0.begin(), IE = ent->m_SCT0.end(); - std::vector<SCTinformation>::iterator J,JB = ent->m_SCT1.begin(), JE = ent->m_SCT1.end(); - - double limit = 1. + m_stripLengthTolerance; + // the following code is used to create spacepoints processing all clusters without limits + + for(int n=0; n!=nElements; ++n) { + + int currentIndex = elementIndex[n]; + const InDetDD::SiDetectorElement* currentElement = elements[currentIndex]; + IdentifierHash currentId = clusters[currentIndex]->identifyHash(); - for (; I!=IE; ++I) { - double qm = (*I).qm(); - for (J=JB; J!=JE; ++J) { - double limitm = limit+(stripLengthGapTolerance*qm); + if (m_SCTgapParameter != 0.) { + offset(element, currentElement, slimit); + } + + for (const auto& cluster : *clusters[currentIndex]) { - double a =-(*I).s().dot((*J).qs()); - double b = (*I).q().dot((*J).qs()); - if (fabs(a) > fabs(b)*limitm) continue; + const Amg::Vector2D& locpos = cluster->localPosition(); + std::pair<Amg::Vector3D, Amg::Vector3D > ends(currentElement->endsOfStrip(InDetDD::SiLocalPosition(locpos.y(),locpos.x(),0.))); + InDet::SCTinformation sctInfo(cluster,ends.first, ends.second,vertexVec,locpos.x()); + + for(auto& sct : sctInfos) { + Trk::SpacePoint* sp = makeSCT_SpacePoint(sct,sctInfo,Id,currentId,limit,slimit); + if(sp) tmpSpacePoints.push_back(sp); + } + } + // If you are processing the opposite element, save the space points into + // the spacepointCollection and clear the temporary collection + if( currentIndex==otherSideIndex && !tmpSpacePoints.empty() ) { + spacepointCollection->reserve(tmpSpacePoints.size()+spacepointCollection->size()); + for (Trk::SpacePoint* sp: tmpSpacePoints) { + spacepointCollection->push_back(sp); + } + tmpSpacePoints.clear(); + } + } + // fill the space point collection for eta/phi overlapping clusters + if(!tmpSpacePoints.empty()) { + spacepointoverlapCollection->reserve(tmpSpacePoints.size()+spacepointoverlapCollection->size()); + for (Trk::SpacePoint* sp: tmpSpacePoints) { + spacepointoverlapCollection->push_back(sp); + } + } + } - double qn = (*J).qm(); - double limitn = limit+(stripLengthGapTolerance*qn); + //-------------------------------------------------------------------------- + // + Trk::SpacePoint* SiSpacePointMakerTool::makeSCT_SpacePoint + (InDet::SCTinformation& In0,InDet::SCTinformation& In1,IdentifierHash ID0,IdentifierHash ID1,double limit,double slimit) const { - double c =-(*J).s().dot((*I).qs()); - double d = (*J).q().dot((*I).qs()); - if (fabs(c) > fabs(d)*limitn) continue; - double m = a/b; - double n = c/d; + double a =-In0.traj_direction().dot(In1.normal()); + double b = In0.strip_direction().dot(In1.normal()); + double l0 = In0.oneOverStrip()*slimit+limit ; - if (m > limit or n > limit) { + if(fabs(a) > (fabs(b)*l0)) return nullptr; - double cs = (*I).q().dot((*J).q())*(qm*qm); - double dm = (m-1.); - double dmn = (n-1.)*cs; - if (dmn > dm) dm = dmn; - m-=dm; n-=(dm/cs); + double c =-In1.traj_direction().dot(In0.normal()); + double d = In1.strip_direction().dot(In0.normal()); + double l1 = In1.oneOverStrip()*slimit+limit ; - } else if (m < -limit or n < -limit) { + if(fabs(c) > (fabs(d)*l1)) return nullptr; - double cs = (*I).q().dot((*J).q())*(qm*qm); - double dm = -(1.+m); - double dmn = -(1.+n)*cs; - if (dmn > dm) dm = dmn; - m+=dm; n+=(dm/cs); + double m = a/b; - } + if(slimit!=0.) { - if (fabs(m) > limit or fabs(n) > limit) continue; + double n = c/d; - Amg::Vector3D point((*I).position(m)); + if (m > limit || n > limit) { - ATH_MSG_VERBOSE("SpacePoint generated at: ( " << point.x() << " , " << point.y() << " , " << point.z() << " ) "); - const std::pair<IdentifierHash, IdentifierHash> elementIdList(ent->m_element0->identifyHash(), ent->m_element1->identifyHash()); - const std::pair<const Trk::PrepRawData*, const Trk::PrepRawData*>* - clusList = new std::pair<const Trk::PrepRawData*, const Trk::PrepRawData*>((*I).cluster(), (*J).cluster()); - ent->m_tmpSpacePoints.push_back(new InDet::SCT_SpacePoint(elementIdList, new Amg::Vector3D(point), clusList)); + double cs = In0.strip_direction().dot(In1.strip_direction())*(In0.oneOverStrip()*In0.oneOverStrip()); + double dm = (m-1); + double dmn = (n-1.)*cs; + if(dmn > dm) dm = dmn; + m-=dm; n-=(dm/cs); + if(fabs(m) > limit || fabs(n) > limit) return nullptr; + } else if(m < -limit || n < -limit) { + + double cs = In0.strip_direction().dot(In1.strip_direction())*(In0.oneOverStrip()*In0.oneOverStrip()); + double dm = -(1.+m); + double dmn = -(1.+n)*cs; + if(dmn > dm) dm = dmn; + m+=dm; n+=(dm/cs); + if(fabs(m) > limit || fabs(n) > limit) return nullptr; } } + Amg::Vector3D point(In0.position(m)); + + const std::pair<IdentifierHash,IdentifierHash> elementIdList(ID0,ID1); + const std::pair<const Trk::PrepRawData*,const Trk::PrepRawData*>* + clusList = new std::pair<const Trk::PrepRawData*,const Trk::PrepRawData*>(In0.cluster(),In1.cluster()); + return new InDet::SCT_SpacePoint(elementIdList,new Amg::Vector3D(point),clusList); } - + }