StripDesign.icc

Go to the documentation of this file.

/*
  Copyright (C) 2002-2023 CERN for the benefit of the ATLAS collaboration
*/
#ifndef MUONREADOUTGEOMETRYR4_STRIPDESIGN_ICC
#define MUONREADOUTGEOMETRYR4_STRIPDESIGN_ICC
 
#define SET_DIRECTION(METHOD, DIRVEC, FROM, TO)                             \
     inline const Amg::Vector2D& StripDesign::METHOD() const {              \
        if (!DIRVEC) {                                                      \
             DIRVEC.set(std::make_unique<Amg::Vector2D>((TO-FROM).unit())); \
        }                                                                   \
        return (*DIRVEC);                                                   \
     }
 
namespace MuonGMR4 {
    using CheckVector2D = StripDesign::CheckVector2D;
    inline const Amg::Vector2D& StripDesign::cornerBotLeft() const { return m_bottomLeft; }
    inline const Amg::Vector2D& StripDesign::cornerBotRight() const { return m_bottomRight; }
    inline const Amg::Vector2D& StripDesign::cornerTopLeft() const { return m_topLeft; }
    inline const Amg::Vector2D& StripDesign::cornerTopRight() const { return m_topRight; }
    inline const Amg::Vector2D& StripDesign::stripNormal() const { return m_stripNormal; }
    inline const Amg::Vector2D& StripDesign::stripDir() const { return m_stripDir; }
    inline const Amg::Vector2D& StripDesign::firstStripPos() const { return m_firstStripPos; }
 
/// Additional Functions for Diamonds
    inline double StripDesign::yCutout() const { return m_yCutout; }
/// ==============================================================================================
    SET_DIRECTION(edgeDirTop, m_dirTopEdge, m_topLeft, m_topRight)
    SET_DIRECTION(edgeDirBottom, m_dirBotEdge, m_bottomLeft, m_bottomRight)
    SET_DIRECTION(edgeDirLeft, m_dirLeftEdge, m_topLeft, m_bottomLeft)
    SET_DIRECTION(edgeDirRight, m_dirRightEdge, m_topRight, m_bottomRight)
 
    inline double StripDesign::lenTopEdge() const {
        return isFlipped() ? 2.*(longHalfHeight()  + m_cutLongEdge):  m_lenSlopEdge;
    }
    inline double StripDesign::lenLeftEdge() const {return isFlipped() ? m_lenSlopEdge : 2.* shortHalfHeight();}
    inline double StripDesign::lenBottomEdge() const { return isFlipped() ? 2.* shortHalfHeight() : m_lenSlopEdge; }
    inline double StripDesign::lenRightEdge() const{ return isFlipped() ? m_lenSlopEdge  : 2.* longHalfHeight(); }
 
    inline double StripDesign::stripPitch() const { return m_stripPitch; }
    inline double StripDesign::stripWidth() const { return m_stripWidth; }
    inline double StripDesign::halfWidth() const { return m_halfX; }
    inline double StripDesign::shortHalfHeight() const { return m_shortHalfY; }
    inline double StripDesign::longHalfHeight() const { return m_longHalfY; }
    inline int StripDesign::numStrips() const { return m_numStrips; }
    inline bool StripDesign::hasStereoAngle() const { return m_hasStereo; }
    inline bool StripDesign::isFlipped() const { return m_isFlipped; }
    inline double StripDesign::stereoAngle() const {return m_stereoAngle; }
    inline int StripDesign::firstStripNumber() const { return m_channelShift; }
    inline Amg::Vector2D StripDesign::stripPosition(int stripCh) const {
        return m_firstStripPos + (1.*stripCh) * stripPitch() *  Amg::Vector2D::UnitX();
    }
    
    inline CheckVector2D StripDesign::leftInterSect(int stripCh, bool uncapped) const {
        /// Nominal strip position position
        return leftInterSect(stripPosition(stripCh), uncapped);
    }
    inline CheckVector2D StripDesign::leftInterSect(const Amg::Vector2D& stripPos, bool uncapped) const {
        /// Flipped diamonds may not exceed the halfLength along the x-coordinate
        if (yCutout() && isFlipped() && std::abs(stripPos.x()) > longHalfHeight()) {
            return std::nullopt;
        }        
        /// We expect lambda to be positive
        const std::optional<double> lambda = Amg::intersect<2>(stripPos, m_stripDir, cornerTopLeft(), edgeDirTop());
        if (!lambda) {
            ATH_MSG_WARNING("The strip "<<Amg::toString(stripPos, 2)
                          <<" does not intersect the top edge "<<Amg::toString(cornerTopLeft(), 2)
                          <<" + lambda "<<Amg::toString(edgeDirTop(), 2));
            return std::nullopt;
        }
        /// If the channel is a stereo channel && lamda is either smaller 0 or longer
        /// then the bottom edge, then it's a routed strip
        if (!uncapped && ((*lambda) < 0  || (*lambda) > lenTopEdge())) {
            const bool pickLeft =  (*lambda) < 0;
            const Amg::Vector2D& edgePos{pickLeft ?  cornerTopLeft() : cornerTopRight()};
            const Amg::Vector2D& edgeDir{pickLeft ?  edgeDirLeft() : edgeDirRight()};
            const double edgeLen = pickLeft ? lenLeftEdge(): lenRightEdge();
            const std::optional<double> top_line = Amg::intersect<2>(stripPos, m_stripDir, edgePos, edgeDir);
            if (!top_line || (*top_line) <0. || (*top_line) > edgeLen) {
                ATH_MSG_DEBUG(__func__<<"() The strip "<<Amg::toString(stripPos, 2)
                                <<" does not intersect with the edge "<<Amg::toString(edgePos));
                return std::nullopt;
            }
            return std::make_optional<Amg::Vector2D>(edgePos + (*top_line) * edgeDir);
       }
       Amg::Vector2D isect = cornerTopLeft() + (*lambda) * edgeDirTop();
        /// In case of diamonds, ensure that the trapezoid is cut on the upper half
       if (!uncapped && yCutout() && !isFlipped()) {
            isect.y() = std::min(isect.y(), longHalfHeight());
       }
       return std::make_optional<Amg::Vector2D>(std::move(isect));
    }
    inline CheckVector2D StripDesign::rightInterSect(const Amg::Vector2D& stripPos, bool uncapped) const {
        std::optional<double> lambda = Amg::intersect<2>(stripPos, m_stripDir, cornerBotLeft(), edgeDirBottom());
        if (!lambda) {
            ATH_MSG_WARNING("The strip "<<Amg::toString(stripPos,2)
                          <<" does not intersect the left edge "<<Amg::toString(cornerBotLeft(), 2)
                          <<" + lambda "<<Amg::toString(edgeDirBottom(), 2));
            return std::nullopt;
        }
        /// If the channel is a stereo channel && lamda is either smaller 0 or longer
        /// then the bottom edge, then it's a routed strip
        if (!uncapped  && ( (*lambda) < 0 || (*lambda) > lenBottomEdge())) {
            const bool selLeftEdge = (*lambda) < 0; 
            const Amg::Vector2D& edgePos{selLeftEdge ? cornerBotLeft() : cornerBotRight()};
            const Amg::Vector2D& edgeDir{selLeftEdge ? edgeDirLeft() : edgeDirRight()};
            const double edgeLen{selLeftEdge ? lenLeftEdge() : lenRightEdge()};
            std::optional<double> bottom_line = Amg::intersect<2>(stripPos, m_stripDir, edgePos, edgeDir);
            if (!bottom_line || (*bottom_line) > 0 || (*bottom_line) < -edgeLen) {
                ATH_MSG_DEBUG("The strip "<<Amg::toString(stripPos, 2)<<" does not intersect with the left/right edges ");
                return std::nullopt;
            }
            return std::make_optional<Amg::Vector2D>(edgePos + (*bottom_line) * edgeDir);
        }
        Amg::Vector2D isect = cornerBotLeft() + (*lambda) * edgeDirBottom();
        /// In case of diamonds, ensure that the trapezoid is cut on the upper half
        if (!uncapped && yCutout() && !isFlipped()) {
            isect.y() = std::max(isect.y(), -longHalfHeight());
        }
        return std::make_optional<Amg::Vector2D>(std::move(isect));
    }
    //============================================================================
    inline double StripDesign::stripLength(int stripNumber) const {
        const int stripCh = (stripNumber - firstStripNumber());
        if (stripCh < 0 ||  stripCh > numStrips()) {
            ATH_MSG_WARNING("center() -- Invalid strip number given "<<stripNumber<<" allowed range ["
                         <<firstStripNumber()<<";"<<numStrips()<<"]");
            return 0.;
        }
        const Amg::Vector2D strip{stripPosition(stripCh)};
        CheckVector2D bottomIsect{}, topIsect{};
        if (!(bottomIsect = rightInterSect(strip)) || !(topIsect = leftInterSect(strip))) {
            ATH_MSG_WARNING("The strip "<<Amg::toString(strip)<<" is outside of the trapezoid");
            return 0;
        }
        const Amg::Vector2D delta = (*bottomIsect) - (*topIsect);
        double stripLen = std::hypot(delta.x(), delta.y());
        ATH_MSG_VERBOSE("Strip "<<stripNumber<<" right edge: "<<Amg::toString((*bottomIsect), 2)
                     <<" left edge: "<<Amg::toString( (*topIsect),2)<<" length: "<<stripLen);
        return stripLen;
    }
    inline CheckVector2D StripDesign::rightInterSect(int stripCh, bool uncapped) const {
        return rightInterSect(stripPosition(stripCh), uncapped);
    }
    
    inline CheckVector2D StripDesign::stripCenter(int stripCh) const {
       CheckVector2D botIsect{}, topIsect{};
       if (!(botIsect = leftInterSect(stripCh)) || !(topIsect = rightInterSect(stripCh))){
           return std::nullopt;
       }
       return std::make_optional<Amg::Vector2D>(0.5 *(std::move(*botIsect) + std::move(*topIsect)));
    }
 
    inline CheckVector2D StripDesign::center(int stripNumber) const {
        const int stripCh = (stripNumber - firstStripNumber());
        if (stripCh < 0 ||  stripCh > numStrips()) {
            ATH_MSG_WARNING("center() -- Invalid strip number given "<<stripNumber<<" allowed range ["
                         <<firstStripNumber()<<";"<<numStrips()<<"]");
            return std::nullopt;
        }
        CheckVector2D cen{stripCenter(stripCh)};
        if (!cen) return std::nullopt;
        return std::make_optional<Amg::Vector2D>(m_etaToStereo * std::move(*cen));
    }    
    inline CheckVector2D StripDesign::leftEdge(int stripNumber) const {
       const int stripCh = (stripNumber - firstStripNumber());
       if (stripCh < 0 ||  stripCh > numStrips()) {
            ATH_MSG_WARNING("leftEdge() -- Invalid strip number given "<<stripNumber<<" allowed range ["
                          <<firstStripNumber()<<";"<<numStrips()<<"]");
            return std::nullopt;
       }
       CheckVector2D iSect{leftInterSect(stripCh)};
       if (!iSect) return std::nullopt;
       return std::make_optional<Amg::Vector2D>(m_etaToStereo * std::move(*iSect));
    }
    inline CheckVector2D StripDesign::rightEdge(int stripNumber) const {
        const int stripCh = (stripNumber - firstStripNumber());
        if (stripCh < 0 ||  stripCh > numStrips()) {
            ATH_MSG_WARNING("rightEdge() -- Invalid strip number given "<<stripNumber<<" allowed range ["
                          <<firstStripNumber()<<";"<<numStrips()<<"]");
            return std::nullopt;
        }
        CheckVector2D iSect{rightInterSect(stripCh)};
        if (!iSect) return std::nullopt;
        return std::make_optional<Amg::Vector2D>(m_etaToStereo * std::move(*iSect));
    }
    inline double StripDesign::distanceToStrip(const Amg::Vector2D& pos, int stripNumber) const {
        CheckVector2D stripCent = center(stripNumber);
        if (!stripCent) return std::numeric_limits<double>::max();
        return (pos - (*stripCent)).x();
    }
    inline int StripDesign::stripNumber(const Amg::Vector2D& extPos) const {
        const Amg::Vector2D pos = m_stereoToEta * extPos;
        std::optional<double> stripDist = Amg::intersect<2>(m_firstStripPos, stripDir(),
                                                            pos, -Amg::Vector2D::UnitX());
        if (!stripDist) {
            return -1;
        }
        if ((*stripDist) < -0.5*stripPitch()) {
            ATH_MSG_VERBOSE("Object "<<Amg::toString(pos)<<" is way before the first strip "<<Amg::toString(m_firstStripPos)
                        <<" distance: "<<(*stripDist)<<", half-pitch: "<<0.5*stripPitch()<<".");
            return -1;
        }
        const double stripPitches = std::abs((*stripDist) / stripPitch());
        if (stripPitches > (1.*numStrips() + 0.5)) {
            ATH_MSG_VERBOSE("Object "<<Amg::toString(pos)<<" is not covered by any strip."
                          <<" Virtual strip pitches: "<<stripPitches<<" max strips: "<<numStrips());
             return -1;
        }
        const int chNum = std::min(static_cast<int>(std::round(stripPitches)), numStrips());       
        if (!insideBoundaries(pos)){
            ATH_MSG_VERBOSE("Object "<<Amg::toString(pos)<<" is outside the active area of the trapezoid. "
                          <<"Distance to strip: "<<stripNormal().dot(pos - (*stripCenter(chNum)))
                          <<","<<stripDir().dot(pos - (*stripCenter(chNum)))
                          <<" stripLength: "<<0.5*stripLength(chNum + firstStripNumber()));
            return -1;
        }
        return chNum + firstStripNumber();
    }
 
    inline bool StripDesign::insideTrapezoid(const Amg::Vector2D& extPos) const {
        return insideBoundaries(m_stereoToEta * extPos);
    }
    inline bool StripDesign::insideBoundaries(const Amg::Vector2D& pos) const {
        /// One needs to go backward to reach the bottom edge
        CheckVector2D botIsect = rightInterSect(pos);
        if (!botIsect || stripDir().dot(pos - (*botIsect)) < 0) return false;
        CheckVector2D topIsect = leftInterSect(pos);
        return topIsect && stripDir().dot(pos - (*topIsect)) < 0;       
    }
}
#undef SET_DIRECTION
#endif