1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
|
/// \file Plane.h
/// \author mgheata
#ifndef VECGEOM_VOLUMES_PLANE_H_
#define VECGEOM_VOLUMES_PLANE_H_
#include "VecGeom/base/Global.h"
#include "VecGeom/base/AlignedBase.h"
#include "VecGeom/volumes/kernel/GenericKernels.h"
namespace vecgeom {
VECGEOM_DEVICE_FORWARD_DECLARE(class Plane;);
VECGEOM_DEVICE_DECLARE_CONV(class, Plane);
inline namespace VECGEOM_IMPL_NAMESPACE {
/// A plane defined by the distance to origin and a normal vector. The normal
/// direction points to the outside halfspace.
class Plane : public AlignedBase {
private:
Vector3D<Precision> fNormal; ///< Normalized normal of the plane.
Precision fDistance; ///< Distance from plane to origin (0, 0, 0).
public:
VECCORE_ATT_HOST_DEVICE
Plane() : fNormal(), fDistance(0.) {}
VECCORE_ATT_HOST_DEVICE
~Plane() {}
VECCORE_ATT_HOST_DEVICE
VECGEOM_FORCE_INLINE
Vector3D<Precision> const &GetNormal() const { return fNormal; }
VECCORE_ATT_HOST_DEVICE
VECGEOM_FORCE_INLINE
Precision GetDistance() const { return fDistance; }
VECCORE_ATT_HOST_DEVICE
void Set(Vector3D<Precision> const &normal, Vector3D<Precision> const &origin)
{
Vector3D<Precision> fixedNormal(normal);
fixedNormal.FixZeroes();
Precision inverseLength = 1. / fixedNormal.Mag();
fNormal = inverseLength * fixedNormal;
fDistance = inverseLength * -fixedNormal.Dot(origin);
}
VECCORE_ATT_HOST_DEVICE
void Set(Vector3D<Precision> const &normal, Precision distance)
{
fNormal = normal;
fDistance = distance;
}
template <typename Real_v>
VECGEOM_FORCE_INLINE
VECCORE_ATT_HOST_DEVICE
Real_v DistPlane(Vector3D<Real_v> const &point) const;
template <typename Real_v, typename Bool_v>
VECGEOM_FORCE_INLINE
VECCORE_ATT_HOST_DEVICE
void Contains(Vector3D<Real_v> const &point, Bool_v &inside) const;
template <typename Real_v, typename Inside_v>
VECGEOM_FORCE_INLINE
VECCORE_ATT_HOST_DEVICE
void Inside(Vector3D<Real_v> const &point, Inside_v &inside) const;
template <typename Real_v>
VECGEOM_FORCE_INLINE
VECCORE_ATT_HOST_DEVICE
void DistanceToIn(Vector3D<Real_v> const &point, Vector3D<Real_v> const &direction, Real_v &distance) const;
template <typename Real_v>
VECGEOM_FORCE_INLINE
VECCORE_ATT_HOST_DEVICE
void DistanceToOut(Vector3D<Real_v> const &point, Vector3D<Real_v> const &direction, Real_v &distance) const;
template <typename Real_v>
VECGEOM_FORCE_INLINE
VECCORE_ATT_HOST_DEVICE
void SafetyToIn(Vector3D<Real_v> const &point, Real_v &distance) const;
template <typename Real_v>
VECGEOM_FORCE_INLINE
VECCORE_ATT_HOST_DEVICE
void SafetyToOut(Vector3D<Real_v> const &point, Real_v &distance) const;
}; // class Plane
std::ostream &operator<<(std::ostream &os, Plane const &plane);
template <typename Real_v>
VECGEOM_FORCE_INLINE
VECCORE_ATT_HOST_DEVICE
Real_v Plane::DistPlane(Vector3D<Real_v> const &point) const
{
// Returns distance from point to plane. This is positive if the point is on
// the outside halfspace, negative otherwise.
return (point.Dot(fNormal) + fDistance);
}
template <typename Real_v, typename Bool_v>
VECGEOM_FORCE_INLINE
VECCORE_ATT_HOST_DEVICE
void Plane::Contains(Vector3D<Real_v> const &point, Bool_v &inside) const
{
// Returns true if the point is in the halfspace behind the normal.
inside = DistPlane(point) < Real_v(0.);
}
template <typename Real_v, typename Inside_v>
VECGEOM_FORCE_INLINE
VECCORE_ATT_HOST_DEVICE
void Plane::Inside(Vector3D<Real_v> const &point, Inside_v &inside) const
{
Real_v dplane = DistPlane(point);
inside = vecCore::Blend(dplane < Real_v(0.0), Inside_v(EInside::kInside), Inside_v(EInside::kOutside));
vecCore::MaskedAssign(inside, vecCore::math::Abs(dplane) < Real_v(kTolerance), Inside_v(EInside::kSurface));
}
template <typename Real_v>
VECGEOM_FORCE_INLINE
VECCORE_ATT_HOST_DEVICE
void Plane::DistanceToIn(Vector3D<Real_v> const &point, Vector3D<Real_v> const &direction, Real_v &distance) const
{
// The function returns a negative distance for points already inside or
// direction going outwards (along the normal)
using Bool_v = vecCore::Mask_v<Real_v>;
distance = -InfinityLength<Real_v>();
Real_v ndd = NonZero(direction.Dot(fNormal));
Real_v saf = DistPlane(point);
Bool_v valid = ndd < Real_v(0.) && saf > Real_v(-kTolerance);
if (vecCore::EarlyReturnAllowed()) {
if (vecCore::MaskEmpty(valid)) return;
}
// If competing with other planes, the maximum distance is winning
vecCore__MaskedAssignFunc(distance, valid, -saf / ndd);
}
template <typename Real_v>
VECGEOM_FORCE_INLINE
VECCORE_ATT_HOST_DEVICE
void Plane::DistanceToOut(Vector3D<Real_v> const &point, Vector3D<Real_v> const &direction, Real_v &distance) const
{
// The function returns infinity if the plane is not hit from inside, negative
// if the point is outside
using Bool_v = vecCore::Mask_v<Real_v>;
distance = InfinityLength<Real_v>();
Real_v ndd = NonZero(direction.Dot(fNormal));
Real_v saf = DistPlane(point);
Bool_v valid = ndd > Real_v(0.) && saf < Real_v(kTolerance);
vecCore__MaskedAssignFunc(distance, saf > Real_v(kTolerance), -InfinityLength<Real_v>());
if (vecCore::EarlyReturnAllowed()) {
if (vecCore::MaskEmpty(valid)) return;
}
// If competing with other planes, the minimum distance is winning
vecCore__MaskedAssignFunc(distance, valid, -saf / ndd);
}
template <typename Real_v>
VECGEOM_FORCE_INLINE
VECCORE_ATT_HOST_DEVICE
void Plane::SafetyToIn(Vector3D<Real_v> const &point, Real_v &distance) const
{
// The safety contains the sign, i.e. if point is inside it is negative.
// If competing with other planes, the maximum distance is winning
distance = DistPlane(point);
}
template <typename Real_v>
VECGEOM_FORCE_INLINE
VECCORE_ATT_HOST_DEVICE
void Plane::SafetyToOut(Vector3D<Real_v> const &point, Real_v &distance) const
{
// The safety contains the sign, i.e. if point is outside it is negative.
// If competing with other planes, the minimum distance is winning
distance = -DistPlane(point);
}
} // namespace VECGEOM_IMPL_NAMESPACE
} // namespace vecgeom
#endif // VECGEOM_VOLUMES_PLANE_H_
|
