VecGeom/volumes/Plane.h

0001 /// \file Plane.h
0002 /// \author mgheata
0003
0004 #ifndef VECGEOM_VOLUMES_PLANE_H_
0005 #define VECGEOM_VOLUMES_PLANE_H_
0006
0007 #include "VecGeom/base/Global.h"
0008
0009 #include "VecGeom/base/AlignedBase.h"
0010 #include "VecGeom/volumes/kernel/GenericKernels.h"
0011
0012 namespace vecgeom {
0013
0014 VECGEOM_DEVICE_FORWARD_DECLARE(class Plane;);
0015 VECGEOM_DEVICE_DECLARE_CONV(class, Plane);
0016
0017 inline namespace VECGEOM_IMPL_NAMESPACE {
0018
0019 /// A plane defined by the distance to origin and a normal vector. The normal
0020 /// direction points to the outside halfspace.
0021
0022 class Plane : public AlignedBase {
0023
0024 private:
0025   Vector3D<Precision> fNormal; ///< Normalized normal of the plane.
0026   Precision fDistance;         ///< Distance from plane to origin (0, 0, 0).
0027
0028 public:
0029   VECCORE_ATT_HOST_DEVICE
0030   Plane() : fNormal(), fDistance(0.) {}
0031
0032   VECCORE_ATT_HOST_DEVICE
0033   ~Plane() {}
0034
0035   VECCORE_ATT_HOST_DEVICE
0036   VECGEOM_FORCE_INLINE
0037   Vector3D<Precision> const &GetNormal() const { return fNormal; }
0038
0039   VECCORE_ATT_HOST_DEVICE
0040   VECGEOM_FORCE_INLINE
0041   Precision GetDistance() const { return fDistance; }
0042
0043   VECCORE_ATT_HOST_DEVICE
0044   void Set(Vector3D<Precision> const &normal, Vector3D<Precision> const &origin)
0045   {
0046     Vector3D<Precision> fixedNormal(normal);
0047     fixedNormal.FixZeroes();
0048     Precision inverseLength = 1. / fixedNormal.Mag();
0049     fNormal                 = inverseLength * fixedNormal;
0050     fDistance               = inverseLength * -fixedNormal.Dot(origin);
0051   }
0052
0053   VECCORE_ATT_HOST_DEVICE
0054   void Set(Vector3D<Precision> const &normal, Precision distance)
0055   {
0056     fNormal   = normal;
0057     fDistance = distance;
0058   }
0059
0060   template <typename Real_v>
0061   VECGEOM_FORCE_INLINE
0062   VECCORE_ATT_HOST_DEVICE
0063   Real_v DistPlane(Vector3D<Real_v> const &point) const;
0064
0065   template <typename Real_v, typename Bool_v>
0066   VECGEOM_FORCE_INLINE
0067   VECCORE_ATT_HOST_DEVICE
0068   void Contains(Vector3D<Real_v> const &point, Bool_v &inside) const;
0069
0070   template <typename Real_v, typename Inside_v>
0071   VECGEOM_FORCE_INLINE
0072   VECCORE_ATT_HOST_DEVICE
0073   void Inside(Vector3D<Real_v> const &point, Inside_v &inside) const;
0074
0075   template <typename Real_v>
0076   VECGEOM_FORCE_INLINE
0077   VECCORE_ATT_HOST_DEVICE
0078   void DistanceToIn(Vector3D<Real_v> const &point, Vector3D<Real_v> const &direction, Real_v &distance) const;
0079
0080   template <typename Real_v>
0081   VECGEOM_FORCE_INLINE
0082   VECCORE_ATT_HOST_DEVICE
0083   void DistanceToOut(Vector3D<Real_v> const &point, Vector3D<Real_v> const &direction, Real_v &distance) const;
0084
0085   template <typename Real_v>
0086   VECGEOM_FORCE_INLINE
0087   VECCORE_ATT_HOST_DEVICE
0088   void SafetyToIn(Vector3D<Real_v> const &point, Real_v &distance) const;
0089
0090   template <typename Real_v>
0091   VECGEOM_FORCE_INLINE
0092   VECCORE_ATT_HOST_DEVICE
0093   void SafetyToOut(Vector3D<Real_v> const &point, Real_v &distance) const;
0094
0095 }; // class Plane
0096
0097 std::ostream &operator<<(std::ostream &os, Plane const &plane);
0098
0099 template <typename Real_v>
0100 VECGEOM_FORCE_INLINE
0101 VECCORE_ATT_HOST_DEVICE
0102 Real_v Plane::DistPlane(Vector3D<Real_v> const &point) const
0103 {
0104   // Returns distance from point to plane. This is positive if the point is on
0105   // the outside halfspace, negative otherwise.
0106   return (point.Dot(fNormal) + fDistance);
0107 }
0108
0109 template <typename Real_v, typename Bool_v>
0110 VECGEOM_FORCE_INLINE
0111 VECCORE_ATT_HOST_DEVICE
0112 void Plane::Contains(Vector3D<Real_v> const &point, Bool_v &inside) const
0113 {
0114   // Returns true if the point is in the halfspace behind the normal.
0115   inside = DistPlane(point) < Real_v(0.);
0116 }
0117
0118 template <typename Real_v, typename Inside_v>
0119 VECGEOM_FORCE_INLINE
0120 VECCORE_ATT_HOST_DEVICE
0121 void Plane::Inside(Vector3D<Real_v> const &point, Inside_v &inside) const
0122 {
0123   Real_v dplane = DistPlane(point);
0124   inside        = vecCore::Blend(dplane < Real_v(0.0), Inside_v(EInside::kInside), Inside_v(EInside::kOutside));
0125   vecCore::MaskedAssign(inside, vecCore::math::Abs(dplane) < Real_v(kTolerance), Inside_v(EInside::kSurface));
0126 }
0127
0128 template <typename Real_v>
0129 VECGEOM_FORCE_INLINE
0130 VECCORE_ATT_HOST_DEVICE
0131 void Plane::DistanceToIn(Vector3D<Real_v> const &point, Vector3D<Real_v> const &direction, Real_v &distance) const
0132 {
0133   // The function returns a negative distance for points already inside or
0134   // direction going outwards (along the normal)
0135   using Bool_v = vecCore::Mask_v<Real_v>;
0136   distance     = -InfinityLength<Real_v>();
0137
0138   Real_v ndd   = NonZero(direction.Dot(fNormal));
0139   Real_v saf   = DistPlane(point);
0140   Bool_v valid = ndd < Real_v(0.) && saf > Real_v(-kTolerance);
0141
0142   if (vecCore::EarlyReturnAllowed()) {
0143     if (vecCore::MaskEmpty(valid)) return;
0144   }
0145   // If competing with other planes, the maximum distance is winning
0146   vecCore__MaskedAssignFunc(distance, valid, -saf / ndd);
0147 }
0148
0149 template <typename Real_v>
0150 VECGEOM_FORCE_INLINE
0151 VECCORE_ATT_HOST_DEVICE
0152 void Plane::DistanceToOut(Vector3D<Real_v> const &point, Vector3D<Real_v> const &direction, Real_v &distance) const
0153 {
0154   // The function returns infinity if the plane is not hit from inside, negative
0155   // if the point is outside
0156   using Bool_v = vecCore::Mask_v<Real_v>;
0157   distance     = InfinityLength<Real_v>();
0158
0159   Real_v ndd   = NonZero(direction.Dot(fNormal));
0160   Real_v saf   = DistPlane(point);
0161   Bool_v valid = ndd > Real_v(0.) && saf < Real_v(kTolerance);
0162
0163   vecCore__MaskedAssignFunc(distance, saf > Real_v(kTolerance), -InfinityLength<Real_v>());
0164
0165   if (vecCore::EarlyReturnAllowed()) {
0166     if (vecCore::MaskEmpty(valid)) return;
0167   }
0168
0169   // If competing with other planes, the minimum distance is winning
0170   vecCore__MaskedAssignFunc(distance, valid, -saf / ndd);
0171 }
0172
0173 template <typename Real_v>
0174 VECGEOM_FORCE_INLINE
0175 VECCORE_ATT_HOST_DEVICE
0176 void Plane::SafetyToIn(Vector3D<Real_v> const &point, Real_v &distance) const
0177 {
0178   // The safety contains the sign, i.e. if point is inside it is negative.
0179   // If competing with other planes, the maximum distance is winning
0180   distance = DistPlane(point);
0181 }
0182
0183 template <typename Real_v>
0184 VECGEOM_FORCE_INLINE
0185 VECCORE_ATT_HOST_DEVICE
0186 void Plane::SafetyToOut(Vector3D<Real_v> const &point, Real_v &distance) const
0187 {
0188   // The safety contains the sign, i.e. if point is outside it is negative.
0189   // If competing with other planes, the minimum distance is winning
0190   distance = -DistPlane(point);
0191 }
0192
0193 } // namespace VECGEOM_IMPL_NAMESPACE
0194
0195 } // namespace vecgeom
0196
0197 #endif // VECGEOM_VOLUMES_PLANE_H_