Back to home page

EIC code displayed by LXR

 
 

    


Warning, file /include/VecGeom/volumes/kernel/OrbImplementation.h was not indexed or was modified since last indexation (in which case cross-reference links may be missing, inaccurate or erroneous).

0001 // This file is part of VecGeom and is distributed under the
0002 // conditions in the file LICENSE.txt in the top directory.
0003 // For the full list of authors see CONTRIBUTORS.txt and `git log`.
0004 
0005 /// \brief This file implements the algorithms for Orb
0006 /// \file volumes/kernel/orbImplementation.h
0007 /// \author Raman Sehgal
0008 
0009 /// History notes:
0010 /// 2014 - 2015: original development (abstracted kernels); Raman Sehgal
0011 /// July 2016: revision + moving to new backend structure (Raman Sehgal)
0012 
0013 #ifndef VECGEOM_VOLUMES_KERNEL_ORBIMPLEMENTATION_H_
0014 #define VECGEOM_VOLUMES_KERNEL_ORBIMPLEMENTATION_H_
0015 
0016 #include "VecGeom/base/Vector3D.h"
0017 #include "VecGeom/volumes/OrbStruct.h"
0018 #include "VecGeom/volumes/kernel/GenericKernels.h"
0019 #include <VecCore/VecCore>
0020 
0021 #include <cstdio>
0022 
0023 namespace vecgeom {
0024 
0025 VECGEOM_DEVICE_FORWARD_DECLARE(struct OrbImplementation;);
0026 VECGEOM_DEVICE_DECLARE_CONV(struct, OrbImplementation);
0027 
0028 inline namespace VECGEOM_IMPL_NAMESPACE {
0029 
0030 class PlacedOrb;
0031 template <typename T>
0032 struct OrbStruct;
0033 class UnplacedOrb;
0034 
0035 struct OrbImplementation {
0036 
0037   using PlacedShape_t    = PlacedOrb;
0038   using UnplacedStruct_t = OrbStruct<Precision>;
0039   using UnplacedVolume_t = UnplacedOrb;
0040 
0041   template <typename Real_v, typename Bool_v>
0042   VECGEOM_FORCE_INLINE VECCORE_ATT_HOST_DEVICE static void Contains(UnplacedStruct_t const &orb,
0043                                                                     Vector3D<Real_v> const &point, Bool_v &inside)
0044   {
0045     Bool_v unused(false), outside(false);
0046     GenericKernelForContainsAndInside<Real_v, Bool_v, false>(orb, point, unused, outside);
0047     inside = !outside;
0048   }
0049 
0050   // BIG QUESTION: DO WE WANT TO GIVE ALL 3 TEMPLATE PARAMETERS
0051   // -- OR -- DO WE WANT TO DEDUCE Bool_v, Index_t from Real_v???
0052   template <typename Real_v, typename Inside_t>
0053   VECGEOM_FORCE_INLINE VECCORE_ATT_HOST_DEVICE static void Inside(UnplacedStruct_t const &orb,
0054                                                                   Vector3D<Real_v> const &point, Inside_t &inside)
0055   {
0056 
0057     using Bool_v       = vecCore::Mask_v<Real_v>;
0058     using InsideBool_v = vecCore::Mask_v<Inside_t>;
0059     Bool_v completelyinside, completelyoutside;
0060     GenericKernelForContainsAndInside<Real_v, Bool_v, true>(orb, point, completelyinside, completelyoutside);
0061     inside = EInside::kSurface;
0062     vecCore::MaskedAssign(inside, (InsideBool_v)completelyoutside, Inside_t(EInside::kOutside));
0063     vecCore::MaskedAssign(inside, (InsideBool_v)completelyinside, Inside_t(EInside::kInside));
0064   }
0065 
0066   template <typename Real_v, typename Bool_v, bool ForInside>
0067   VECGEOM_FORCE_INLINE VECCORE_ATT_HOST_DEVICE static void GenericKernelForContainsAndInside(
0068       UnplacedStruct_t const &orb, Vector3D<Real_v> const &localPoint, Bool_v &completelyinside,
0069       Bool_v &completelyoutside)
0070   {
0071     Precision fR = orb.fR;
0072     Real_v rad2  = localPoint.Mag2();
0073     Real_v tolR  = fR - Real_v(kTolerance);
0074     if (ForInside) completelyinside = (rad2 <= tolR * tolR);
0075     tolR              = fR + Real_v(kTolerance);
0076     completelyoutside = (rad2 >= tolR * tolR);
0077     return;
0078   }
0079 
0080   template <typename Real_v>
0081   VECGEOM_FORCE_INLINE VECCORE_ATT_HOST_DEVICE static void DistanceToIn(UnplacedStruct_t const &orb,
0082                                                                         Vector3D<Real_v> const &point,
0083                                                                         Vector3D<Real_v> const &direction,
0084                                                                         Real_v const & /*stepMax*/, Real_v &distance)
0085   {
0086     using Bool_v         = vecCore::Mask_v<Real_v>;
0087     distance             = kInfLength;
0088     Real_v rad           = point.Mag();
0089     Bool_v isPointInside = (rad < Real_v(orb.fR - kTolerance));
0090     vecCore__MaskedAssignFunc(distance, isPointInside, Real_v(-1.));
0091     Bool_v done = isPointInside;
0092     if (vecCore::MaskFull(done)) return;
0093 
0094     Real_v pDotV3D          = point.Dot(direction);
0095     Bool_v isPointOnSurface = (rad >= Real_v(orb.fR - kTolerance)) && (rad <= Real_v(orb.fR + kTolerance));
0096     Bool_v cond             = (isPointOnSurface && (pDotV3D < Real_v(0.)));
0097     vecCore__MaskedAssignFunc(distance, !done && cond, Real_v(0.));
0098     done |= cond;
0099     if (vecCore::MaskFull(done)) return;
0100     Real_v dist(kInfLength);
0101     vecCore::MaskedAssign(
0102         distance, !done && DetectIntersectionAndCalculateDistance<Real_v, true>(orb, point, direction, dist), dist);
0103   }
0104 
0105   template <typename Real_v>
0106   VECGEOM_FORCE_INLINE VECCORE_ATT_HOST_DEVICE static void DistanceToOut(UnplacedStruct_t const &orb,
0107                                                                          Vector3D<Real_v> const &point,
0108                                                                          Vector3D<Real_v> const &direction,
0109                                                                          Real_v const & /* stepMax */, Real_v &distance)
0110   {
0111     using Bool_v = vecCore::Mask_v<Real_v>;
0112 
0113     distance = kInfLength;
0114 
0115     Real_v rad            = point.Mag();
0116     Bool_v isPointOutside = (rad > Real_v(orb.fR + kTolerance));
0117     vecCore__MaskedAssignFunc(distance, isPointOutside, Real_v(-1.));
0118     Bool_v done = isPointOutside;
0119     if (vecCore::MaskFull(done)) return;
0120 
0121     Real_v pDotV3D          = point.Dot(direction);
0122     Bool_v isPointOnSurface = (rad >= Real_v(orb.fR - kTolerance)) && (rad <= Real_v(orb.fR + kTolerance));
0123     Bool_v cond             = (isPointOnSurface && (pDotV3D > Real_v(0.)));
0124     vecCore__MaskedAssignFunc(distance, !done && cond, Real_v(0.));
0125     done |= cond;
0126     if (vecCore::MaskFull(done)) return;
0127     Real_v dist(kInfLength);
0128     vecCore::MaskedAssign(
0129         distance, !done && DetectIntersectionAndCalculateDistance<Real_v, false>(orb, point, direction, dist), dist);
0130 
0131     return;
0132   }
0133 
0134   template <typename Real_v>
0135   VECGEOM_FORCE_INLINE VECCORE_ATT_HOST_DEVICE static void SafetyToIn(UnplacedStruct_t const &orb,
0136                                                                       Vector3D<Real_v> const &point, Real_v &safety)
0137   {
0138     using Bool_v         = vecCore::Mask_v<Real_v>;
0139     Real_v rad           = point.Mag();
0140     safety               = rad - Real_v(orb.fR);
0141     Bool_v isPointInside = (rad < Real_v(orb.fR - kTolerance));
0142     vecCore__MaskedAssignFunc(safety, isPointInside, Real_v(-1.));
0143     if (vecCore::MaskFull(isPointInside)) return;
0144 
0145     Bool_v isPointOnSurface = (rad > Real_v(orb.fR - kTolerance)) && (rad < Real_v(orb.fR + kTolerance));
0146     vecCore__MaskedAssignFunc(safety, isPointOnSurface, Real_v(0.));
0147   }
0148 
0149   template <typename Real_v>
0150   VECGEOM_FORCE_INLINE VECCORE_ATT_HOST_DEVICE static void SafetyToOut(UnplacedStruct_t const &orb,
0151                                                                        Vector3D<Real_v> const &point, Real_v &safety)
0152   {
0153     using Bool_v = vecCore::Mask_v<Real_v>;
0154 
0155     Real_v rad = point.Mag();
0156     safety     = Real_v(orb.fR) - rad;
0157 
0158     Bool_v isPointOutside = (rad > Real_v(orb.fR + kTolerance));
0159     vecCore__MaskedAssignFunc(safety, isPointOutside, Real_v(-1.));
0160     if (vecCore::MaskFull(isPointOutside)) return;
0161 
0162     Bool_v isPointOnSurface = (rad > Real_v(orb.fR - kTolerance)) && (rad < Real_v(orb.fR + kTolerance));
0163     vecCore__MaskedAssignFunc(safety, isPointOnSurface, Real_v(0.));
0164   }
0165 
0166   template <typename Real_v, bool ForDistanceToIn>
0167   VECGEOM_FORCE_INLINE VECCORE_ATT_HOST_DEVICE static typename vecCore::Mask_v<Real_v>
0168   DetectIntersectionAndCalculateDistance(UnplacedStruct_t const &orb, Vector3D<Real_v> const &point,
0169                                          Vector3D<Real_v> const &direction, Real_v &distance)
0170   {
0171 
0172     using Bool_v   = vecCore::Mask_v<Real_v>;
0173     Real_v rad2    = point.Mag2();
0174     Real_v pDotV3D = point.Dot(direction);
0175     Precision fR   = orb.fR;
0176     Real_v c       = rad2 - fR * fR;
0177     Real_v d2      = (pDotV3D * pDotV3D - c);
0178 
0179     if (ForDistanceToIn) {
0180       Bool_v cond = ((d2 >= Real_v(0.)) && (pDotV3D <= Real_v(0.)));
0181       vecCore__MaskedAssignFunc(distance, cond, (-pDotV3D - Sqrt(vecCore::math::Abs(d2))));
0182       return cond;
0183     } else {
0184       vecCore__MaskedAssignFunc(distance, (d2 >= Real_v(0.)), (-pDotV3D + Sqrt(vecCore::math::Abs(d2))));
0185       return (d2 >= Real_v(0.));
0186     }
0187   }
0188 
0189   template <typename Real_v>
0190   VECGEOM_FORCE_INLINE VECCORE_ATT_HOST_DEVICE static Vector3D<Real_v> NormalKernel(
0191       UnplacedStruct_t const &orb, Vector3D<Real_v> const &point, typename vecCore::Mask_v<Real_v> &valid)
0192   {
0193     Real_v rad2             = point.Mag2();
0194     Real_v invRadius        = Real_v(1.) / Sqrt(rad2);
0195     Vector3D<Real_v> normal = point * invRadius;
0196 
0197     Real_v tolRMaxO = orb.fR + kTolerance;
0198     Real_v tolRMaxI = orb.fR - kTolerance;
0199 
0200     // Check radial surface
0201     valid = ((rad2 <= tolRMaxO * tolRMaxO) && (rad2 >= tolRMaxI * tolRMaxI)); // means we are on surface
0202     return normal;
0203   }
0204 };
0205 } // namespace VECGEOM_IMPL_NAMESPACE
0206 } // namespace vecgeom
0207 
0208 #endif // VECGEOM_VOLUMES_KERNEL_orbIMPLEMENTATION_H_