Skip to content

Commit 257f2ee

Browse files
shahor02shahor02
committed
Revert "MatLUT: faster phi-sector determination; reduction of divisions (#12030)"
This reverts commit 030e396.
1 parent 8188378 commit 257f2ee

2 files changed

Lines changed: 7 additions & 58 deletions

File tree

‎Detectors/Base/include/DetectorsBase/Ray.h‎

Lines changed: 1 addition & 55 deletions
Original file line numberDiff line numberDiff line change
@@ -46,7 +46,7 @@ class Ray
4646
static constexpr float InvalidT = -1e9;
4747
static constexpr float Tiny = 1e-9;
4848

49-
GPUd() Ray() : mP{0.f}, mD{0.f}, mDistXY2(0.f), mDistXY2i(0.f), mDistXYZ(0.f), mXDxPlusYDy(0.f), mXDxPlusYDyRed(0.f), mXDxPlusYDy2(0.f), mR02(0.f), mR12(0.f), mInvDz(Tiny)
49+
GPUd() Ray() : mP{0.f}, mD{0.f}, mDistXY2(0.f), mDistXY2i(0.f), mDistXYZ(0.f), mXDxPlusYDy(0.f), mXDxPlusYDyRed(0.f), mXDxPlusYDy2(0.f), mR02(0.f), mR12(0.f)
5050
{
5151
}
5252
GPUdDefault() ~Ray() CON_DEFAULT;
@@ -61,7 +61,6 @@ class Ray
6161
GPUd() float crossRadial(const MatLayerCyl& lr, int sliceID) const;
6262
GPUd() float crossRadial(float cs, float sn) const;
6363
GPUd() float crossZ(float z) const;
64-
GPUd() float crossZ_fast(float z) const; // without check
6564

6665
GPUd() void getCrossParams(int i, float& par1, float& par2) const
6766
{
@@ -85,49 +84,6 @@ class Ray
8584
return p;
8685
}
8786

88-
// fast, approximate arctan calculation - sufficient for phi sector determination in material budget layers
89-
// source of inspirations:
90-
// - https://mazzo.li/posts/vectorized-atan2.html
91-
// - https://www.coranac.com/documents/arctangent/
92-
GPUd() float atan_approx_5terms(float x) const
93-
{
94-
// fast atan implementation
95-
// function has a max error of ~1.16825e-05f
96-
float a1 = 0.9998660f;
97-
float a3 = -0.3302995f;
98-
float a5 = 0.1801410f;
99-
float a7 = -0.0851330f;
100-
float a9 = 0.0208351f;
101-
102-
float x_sq = x * x;
103-
// todo: force fuse-multiply add
104-
return x * (a1 + x_sq * (a3 + x_sq * (a5 + x_sq * (a7 + x_sq * a9))));
105-
}
106-
107-
GPUd() float atan2_approx(float y, float x) const
108-
{
109-
const bool swap = fabs(x) < fabs(y);
110-
const float atan_input = (swap ? x : y) / (swap ? y : x);
111-
const float pi2 = M_PI_2;
112-
const float pi = M_PI;
113-
float res = atan_approx_5terms(atan_input);
114-
115-
// if swapped, adjust atan output
116-
res = swap ? (atan_input >= 0.0f ? pi2 : -pi2) - res : res;
117-
// adjust quadrants
118-
if (x < 0.0f) {
119-
res = ((y < 0.0f) ? -pi : pi) + res;
120-
}
121-
return res;
122-
}
123-
124-
GPUd() float getPhi_approx(float t) const
125-
{
126-
float p = atan2_approx(mP[1] + t * mD[1], mP[0] + t * mD[0]);
127-
o2::math_utils::bringTo02Pi(p);
128-
return p;
129-
}
130-
13187
GPUd() float getZ(float t) const { return mP[2] + t * mD[2]; }
13288

13389
GPUd() bool validateZRange(float& cpar1, float& cpar2, const MatLayerCyl& lr) const;
@@ -145,7 +101,6 @@ class Ray
145101
float mR12; ///< radius^2 of mP1
146102
float mCrossParams1[2]; ///< parameters of crossing the layer (first parameter)
147103
float mCrossParams2[2]; ///< parameters of crossing the layer (second parameter)
148-
float mInvDz;
149104

150105
ClassDefNV(Ray, 1);
151106
};
@@ -164,7 +119,6 @@ inline Ray::Ray(const math_utils::Point3D<float> point0, const math_utils::Point
164119
mXDxPlusYDy2 = mXDxPlusYDy * mXDxPlusYDy;
165120
mR02 = point0.Perp2();
166121
mR12 = point1.Perp2();
167-
mInvDz = 1.f / mD[2];
168122
}
169123
#endif // !GPUCA_ALIGPUCODE
170124

@@ -180,7 +134,6 @@ GPUdi() Ray::Ray(float x0, float y0, float z0, float x1, float y1, float z1)
180134
mXDxPlusYDy2 = mXDxPlusYDy * mXDxPlusYDy;
181135
mR02 = x0 * x0 + y0 * y0;
182136
mR12 = x1 * x1 + y1 * y1;
183-
mInvDz = 1.f / mD[2];
184137
}
185138

186139
//______________________________________________________
@@ -222,13 +175,6 @@ GPUdi() float Ray::crossZ(float z) const
222175
return mD[2] != 0. ? (z - mP[2]) / mD[2] : InvalidT;
223176
}
224177

225-
//______________________________________________________
226-
GPUdi() float Ray::crossZ_fast(float z) const
227-
{
228-
// calculate t of crossing XY plane at Z
229-
return (z - mP[2]) * mInvDz;
230-
}
231-
232178
//______________________________________________________
233179
GPUdi() bool Ray::validateZRange(float& cpar1, float& cpar2, const MatLayerCyl& lr) const
234180
{

‎Detectors/Base/src/MatLayerCylSet.cxx‎

Lines changed: 6 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -294,8 +294,8 @@ GPUd() MatBudget MatLayerCylSet::getMatBudget(float x0, float y0, float z0, floa
294294
for (int ic = nc; ic--;) {
295295
float cross1, cross2;
296296
ray.getCrossParams(ic, cross1, cross2); // tmax,tmin of crossing the layer
297-
auto phi0 = ray.getPhi_approx(cross1), phi1 = ray.getPhi_approx(cross2), dPhi = phi0 - phi1;
298-
// auto phi0 = ray.getPhi(cross1), phi1 = ray.getPhi(cross2), dPhi = phi0 - phi1;
297+
298+
auto phi0 = ray.getPhi(cross1), phi1 = ray.getPhi(cross2), dPhi = phi0 - phi1;
299299
auto phiID = lr.getPhiSliceID(phi0), phiIDLast = lr.getPhiSliceID(phi1);
300300
// account for eventual wrapping around 0
301301
if (dPhi > 0.f) {
@@ -339,7 +339,10 @@ GPUd() MatBudget MatLayerCylSet::getMatBudget(float x0, float y0, float z0, floa
339339
tEndZ = tEndPhi;
340340
checkMoreZ = false;
341341
} else {
342-
tEndZ = ray.crossZ_fast(lr.getZBinMin(stepZID > 0 ? zID + 1 : zID));
342+
tEndZ = ray.crossZ(lr.getZBinMin(stepZID > 0 ? zID + 1 : zID));
343+
if (tEndZ == Ray::InvalidT) { // track normal to Z axis, abandon Zbin change test
344+
break;
345+
}
343346
}
344347
// account materials of this step
345348
float step = tEndZ > tStartZ ? tEndZ - tStartZ : tStartZ - tEndZ; // the real step is ray.getDist(tEnd-tStart), will rescale all later

0 commit comments

Comments
 (0)