From 08baf58c4e5c6e2aeccead9e1752fe5a4c64c9db Mon Sep 17 00:00:00 2001 From: David Rohr Date: Wed, 9 Dec 2020 20:24:13 +0100 Subject: [PATCH 1/9] GPU: Fix macro name --- .../include/ReconstructionDataFormats/DCA.h | 2 +- .../include/ReconstructionDataFormats/PrimaryVertex.h | 4 ++-- .../include/ReconstructionDataFormats/Vertex.h | 9 ++++++--- DataFormats/Reconstruction/src/DCA.cxx | 4 ++-- DataFormats/Reconstruction/src/PrimaryVertex.cxx | 2 +- DataFormats/Reconstruction/src/Vertex.cxx | 2 +- .../common/include/CommonDataFormat/InteractionRecord.h | 6 +++--- DataFormats/common/src/InteractionRecord.cxx | 2 +- GPU/Common/GPUCommonDef.h | 2 +- 9 files changed, 18 insertions(+), 15 deletions(-) diff --git a/DataFormats/Reconstruction/include/ReconstructionDataFormats/DCA.h b/DataFormats/Reconstruction/include/ReconstructionDataFormats/DCA.h index 2b827d5bb740d..ff6c357d4ca9a 100644 --- a/DataFormats/Reconstruction/include/ReconstructionDataFormats/DCA.h +++ b/DataFormats/Reconstruction/include/ReconstructionDataFormats/DCA.h @@ -16,7 +16,7 @@ #ifndef __OPENCL__ #include #endif -#ifndef ALIGPU_GPUCODE +#ifndef GPUCA_ALIGPUCODE #include #endif diff --git a/DataFormats/Reconstruction/include/ReconstructionDataFormats/PrimaryVertex.h b/DataFormats/Reconstruction/include/ReconstructionDataFormats/PrimaryVertex.h index 8d319333a24da..872e7517e8ff6 100644 --- a/DataFormats/Reconstruction/include/ReconstructionDataFormats/PrimaryVertex.h +++ b/DataFormats/Reconstruction/include/ReconstructionDataFormats/PrimaryVertex.h @@ -36,7 +36,7 @@ class PrimaryVertex : public Vertex> void setIR(const InteractionRecord& ir) { mIRMin = mIRMax = ir; } bool hasUniqueIR() const { return !mIRMin.isDummy() && (mIRMin == mIRMax); } -#ifndef ALIGPU_GPUCODE +#ifndef GPUCA_ALIGPUCODE void print() const; std::string asString() const; #endif @@ -48,7 +48,7 @@ class PrimaryVertex : public Vertex> ClassDefNV(PrimaryVertex, 1); }; -#ifndef ALIGPU_GPUCODE +#ifndef GPUCA_ALIGPUCODE std::ostream& operator<<(std::ostream& os, const o2::dataformats::PrimaryVertex& v); #endif diff --git a/DataFormats/Reconstruction/include/ReconstructionDataFormats/Vertex.h b/DataFormats/Reconstruction/include/ReconstructionDataFormats/Vertex.h index 7abcae8423447..f4d798f863868 100644 --- a/DataFormats/Reconstruction/include/ReconstructionDataFormats/Vertex.h +++ b/DataFormats/Reconstruction/include/ReconstructionDataFormats/Vertex.h @@ -11,12 +11,15 @@ #ifndef ALICEO2_VERTEX_H #define ALICEO2_VERTEX_H +#include "GPUCommonDef.h" +#include "GPUCommonMath.h" #include + #include "CommonDataFormat/TimeStamp.h" #ifndef __OPENCL__ #include #endif -#ifndef ALIGPU_GPUCODE +#ifndef GPUCA_ALIGPUCODE #include #endif @@ -42,7 +45,7 @@ class VertexBase { } -#ifndef ALIGPU_GPUCODE +#ifndef GPUCA_ALIGPUCODE void print() const; std::string asString() const; #endif @@ -144,7 +147,7 @@ class Vertex : public VertexBase ClassDefNV(Vertex, 3); }; -#ifndef ALIGPU_GPUCODE +#ifndef GPUCA_ALIGPUCODE std::ostream& operator<<(std::ostream& os, const o2::dataformats::VertexBase& v); #endif diff --git a/DataFormats/Reconstruction/src/DCA.cxx b/DataFormats/Reconstruction/src/DCA.cxx index 925948b2bc127..d840cefc0ccd9 100644 --- a/DataFormats/Reconstruction/src/DCA.cxx +++ b/DataFormats/Reconstruction/src/DCA.cxx @@ -19,7 +19,7 @@ namespace o2 namespace dataformats { -#ifndef ALIGPU_GPUCODE +#ifndef GPUCA_ALIGPUCODE std::ostream& operator<<(std::ostream& os, const o2::dataformats::DCA& d) { // stream itself @@ -30,7 +30,7 @@ std::ostream& operator<<(std::ostream& os, const o2::dataformats::DCA& d) void DCA::print() const { -#ifndef ALIGPU_GPUCODE +#ifndef GPUCA_ALIGPUCODE std::cout << *this << '\n'; #endif } diff --git a/DataFormats/Reconstruction/src/PrimaryVertex.cxx b/DataFormats/Reconstruction/src/PrimaryVertex.cxx index f94606c2e203a..fb4f7d3b652e5 100644 --- a/DataFormats/Reconstruction/src/PrimaryVertex.cxx +++ b/DataFormats/Reconstruction/src/PrimaryVertex.cxx @@ -18,7 +18,7 @@ namespace o2 namespace dataformats { -#ifndef ALIGPU_GPUCODE +#ifndef GPUCA_ALIGPUCODE std::string PrimaryVertex::asString() const { diff --git a/DataFormats/Reconstruction/src/Vertex.cxx b/DataFormats/Reconstruction/src/Vertex.cxx index b6165566e8aa1..f4a98e4223992 100644 --- a/DataFormats/Reconstruction/src/Vertex.cxx +++ b/DataFormats/Reconstruction/src/Vertex.cxx @@ -17,7 +17,7 @@ namespace o2 namespace dataformats { -#ifndef ALIGPU_GPUCODE +#ifndef GPUCA_ALIGPUCODE std::string VertexBase::asString() const { diff --git a/DataFormats/common/include/CommonDataFormat/InteractionRecord.h b/DataFormats/common/include/CommonDataFormat/InteractionRecord.h index 841197c2af4ee..7ed2b72b0ccb2 100644 --- a/DataFormats/common/include/CommonDataFormat/InteractionRecord.h +++ b/DataFormats/common/include/CommonDataFormat/InteractionRecord.h @@ -14,7 +14,7 @@ #define ALICEO2_INTERACTIONRECORD_H #include "GPUCommonRtypes.h" -#ifndef ALIGPU_GPUCODE +#ifndef GPUCA_ALIGPUCODE #include #include #endif @@ -246,7 +246,7 @@ struct InteractionRecord { return InteractionRecord(l % o2::constants::lhc::LHCMaxBunches, l / o2::constants::lhc::LHCMaxBunches); } -#ifndef ALIGPU_GPUCODE +#ifndef GPUCA_ALIGPUCODE void print() const; std::string asString() const; friend std::ostream& operator<<(std::ostream& stream, InteractionRecord const& ir); @@ -324,7 +324,7 @@ struct InteractionTimeRecord : public InteractionRecord { return !((*this) > other); } -#ifndef ALIGPU_GPUCODE +#ifndef GPUCA_ALIGPUCODE void print() const; std::string asString() const; friend std::ostream& operator<<(std::ostream& stream, InteractionTimeRecord const& ir); diff --git a/DataFormats/common/src/InteractionRecord.cxx b/DataFormats/common/src/InteractionRecord.cxx index 4ceb8f957a45c..59962961548b1 100644 --- a/DataFormats/common/src/InteractionRecord.cxx +++ b/DataFormats/common/src/InteractionRecord.cxx @@ -14,7 +14,7 @@ namespace o2 { -#ifndef ALIGPU_GPUCODE +#ifndef GPUCA_ALIGPUCODE std::string InteractionRecord::asString() const { diff --git a/GPU/Common/GPUCommonDef.h b/GPU/Common/GPUCommonDef.h index 1eabbfc6c4b58..3adf1a81c21a8 100644 --- a/GPU/Common/GPUCommonDef.h +++ b/GPU/Common/GPUCommonDef.h @@ -74,7 +74,7 @@ #endif //Set AliRoot / O2 namespace -#if defined(GPUCA_STANDALONE) || (defined(GPUCA_O2_LIB) && !defined(GPUCA_O2_INTERFACE)) || defined(GPUCA_ALIROOT_LIB) || defined(GPUCA_GPULIBRARY) +#if defined(GPUCA_STANDALONE) || (defined(GPUCA_O2_LIB) && !defined(GPUCA_O2_INTERFACE)) || defined(GPUCA_ALIROOT_LIB) || defined(GPUCA_GPULIBRARY) || defined (GPUCA_GPUCODE) #define GPUCA_ALIGPUCODE #endif #ifdef GPUCA_ALIROOT_LIB From ef43df8dfd314cd236bb448ab7c717da922cca99 Mon Sep 17 00:00:00 2001 From: David Rohr Date: Wed, 9 Dec 2020 20:31:26 +0100 Subject: [PATCH 2/9] GPU: Was missing the GPU annotations for the GPU versions of CartesianND functions... --- .../include/MathUtils/CartesianGPU.h | 26 ++++++++++--------- 1 file changed, 14 insertions(+), 12 deletions(-) diff --git a/Common/MathUtils/include/MathUtils/CartesianGPU.h b/Common/MathUtils/include/MathUtils/CartesianGPU.h index 93cc00633dd2d..4ec08e51bc6c2 100644 --- a/Common/MathUtils/include/MathUtils/CartesianGPU.h +++ b/Common/MathUtils/include/MathUtils/CartesianGPU.h @@ -14,6 +14,8 @@ #ifndef ALICEO2_CARTESIANGPU_H #define ALICEO2_CARTESIANGPU_H +#include "GPUCommonDef.h" + namespace o2::math_utils { @@ -21,25 +23,25 @@ namespace detail { template struct GPUPoint2D { - GPUPoint2D() = default; - GPUPoint2D(T a, T b) : xx(a), yy(b) {} + GPUdDefault() GPUPoint2D() = default; + GPUd() GPUPoint2D(T a, T b) : xx(a), yy(b) {} + GPUd() float X() const { return xx; } + GPUd() float Y() const { return yy; } + GPUd() float R() const { return o2::gpu::CAMath::Sqrt(xx * xx + yy * yy); } + GPUd() void SetX(float v) { xx = v; } + GPUd() void SetY(float v) { yy = v; } T xx; T yy; - float X() const { return xx; } - float Y() const { return yy; } - float R() const { return o2::gpu::CAMath::Sqrt(xx * xx + yy * yy); } - void SetX(float v) { xx = v; } - void SetY(float v) { yy = v; } }; template struct GPUPoint3D : public GPUPoint2D { - GPUPoint3D() = default; - GPUPoint3D(T a, T b, T c) : GPUPoint2D(a, b), zz(c) {} + GPUdDefault() GPUPoint3D() = default; + GPUd() GPUPoint3D(T a, T b, T c) : GPUPoint2D(a, b), zz(c) {} + GPUd() float Z() const { return zz; } + GPUd() float R() const { return o2::gpu::CAMath::Sqrt(GPUPoint2D::xx * GPUPoint2D::xx + GPUPoint2D::yy * GPUPoint2D::yy + zz * zz); } + GPUd() void SetZ(float v) { zz = v; } T zz; - float Z() const { return zz; } - float R() const { return o2::gpu::CAMath::Sqrt(GPUPoint2D::xx * GPUPoint2D::xx + GPUPoint2D::yy * GPUPoint2D::yy + zz * zz); } - void SetZ(float v) { zz = v; } }; } // namespace detail From a75e6fb5550bda54e2db15befb95127c2f7d72bd Mon Sep 17 00:00:00 2001 From: David Rohr Date: Wed, 9 Dec 2020 20:31:58 +0100 Subject: [PATCH 3/9] Make rotateZ / sincos / sector2angle compatible to GPU --- Common/MathUtils/include/MathUtils/Utils.h | 6 +- .../include/MathUtils/detail/CircleXY.h | 9 +- .../include/MathUtils/detail/IntervalXY.h | 109 ++++++++++-------- .../include/MathUtils/detail/trigonometric.h | 46 ++++---- 4 files changed, 91 insertions(+), 79 deletions(-) diff --git a/Common/MathUtils/include/MathUtils/Utils.h b/Common/MathUtils/include/MathUtils/Utils.h index 4b656580d5a5a..0f7b2e16185c8 100644 --- a/Common/MathUtils/include/MathUtils/Utils.h +++ b/Common/MathUtils/include/MathUtils/Utils.h @@ -114,17 +114,17 @@ GPUdi() void sincosd(double ang, double& s, double& c) detail::sincos(ang, s, c); } -#ifndef GPUCA_GPUCODE_DEVICE -inline void rotateZ(float xL, float yL, float& xG, float& yG, float snAlp, float csAlp) +GPUdi() void rotateZ(float xL, float yL, float& xG, float& yG, float snAlp, float csAlp) { return detail::rotateZ(xL, yL, xG, yG, snAlp, csAlp); } -inline void rotateZd(double xL, double yL, double& xG, double& yG, double snAlp, double csAlp) +GPUdi() void rotateZd(double xL, double yL, double& xG, double& yG, double snAlp, double csAlp) { return detail::rotateZ(xL, yL, xG, yG, snAlp, csAlp); } +#ifndef GPUCA_GPUCODE_DEVICE inline void rotateZInv(float xG, float yG, float& xL, float& yL, float snAlp, float csAlp) { detail::rotateZInv(xG, yG, xL, yL, snAlp, csAlp); diff --git a/Common/MathUtils/include/MathUtils/detail/CircleXY.h b/Common/MathUtils/include/MathUtils/detail/CircleXY.h index f33f167df68e5..cf740b455441c 100644 --- a/Common/MathUtils/include/MathUtils/detail/CircleXY.h +++ b/Common/MathUtils/include/MathUtils/detail/CircleXY.h @@ -16,6 +16,7 @@ #ifndef MATHUTILS_INCLUDE_MATHUTILS_DETAIL_CIRCLEXY_H_ #define MATHUTILS_INCLUDE_MATHUTILS_DETAIL_CIRCLEXY_H_ +#include "GPUCommonDef.h" #include "GPUCommonRtypes.h" namespace o2 @@ -32,18 +33,18 @@ struct CircleXY { T rC; // circle radius T xC; // x-center T yC; // y-center - CircleXY(T r = T(), T x = T(), T y = T()); - T getCenterD2() const; + GPUd() CircleXY(T r = T(), T x = T(), T y = T()); + GPUd() T getCenterD2() const; ClassDefNV(CircleXY, 2); }; template -CircleXY::CircleXY(T r, T x, T y) : rC(std::move(r)), xC(std::move(x)), yC(std::move(y)) +GPUdi() CircleXY::CircleXY(T r, T x, T y) : rC(r), xC(x), yC(y) { } template -inline T CircleXY::getCenterD2() const +GPUdi() T CircleXY::getCenterD2() const { return xC * xC + yC * yC; } diff --git a/Common/MathUtils/include/MathUtils/detail/IntervalXY.h b/Common/MathUtils/include/MathUtils/detail/IntervalXY.h index 67d8719a9d3aa..2d77c57f2bb21 100644 --- a/Common/MathUtils/include/MathUtils/detail/IntervalXY.h +++ b/Common/MathUtils/include/MathUtils/detail/IntervalXY.h @@ -16,9 +16,12 @@ #ifndef MATHUTILS_INCLUDE_MATHUTILS_DETAIL_INTERVALXY_H_ #define MATHUTILS_INCLUDE_MATHUTILS_DETAIL_INTERVALXY_H_ +#include "GPUCommonDef.h" #include "GPUCommonRtypes.h" +#ifndef GPUCA_GPUCODE_DEVICE #include #include +#endif #include "MathUtils/detail/CircleXY.h" @@ -36,30 +39,32 @@ class IntervalXY using value_t = T; ///< 2D interval in lab frame defined by its one edge and signed projection lengths on X,Y axes - IntervalXY(T x = T(), T y = T(), T dx = T(), T dy = T()); - T getX0() const; - T getY0() const; - T& getX0(); - T& getY0(); - T getX1() const; - T getY1() const; - T getDX() const; - T getDY() const; - T& getDX(); - T& getDY(); - - void setX0(T x0); - void setY0(T y0); - void setX1(T x1); - void setY1(T y1); - void setDX(T dx); - void setDY(T dy); - void setEdges(T x0, T y0, T x1, T y1); - - void eval(T t, T& x, T& y) const; + GPUd() IntervalXY(T x = T(), T y = T(), T dx = T(), T dy = T()); + GPUd() T getX0() const; + GPUd() T getY0() const; + GPUd() T& getX0(); + GPUd() T& getY0(); + GPUd() T getX1() const; + GPUd() T getY1() const; + GPUd() T getDX() const; + GPUd() T getDY() const; + GPUd() T& getDX(); + GPUd() T& getDY(); + + GPUd() void setX0(T x0); + GPUd() void setY0(T y0); + GPUd() void setX1(T x1); + GPUd() void setY1(T y1); + GPUd() void setDX(T dx); + GPUd() void setDY(T dy); + GPUd() void setEdges(T x0, T y0, T x1, T y1); + + GPUd() void eval(T t, T& x, T& y) const; +#ifndef GPUCA_GPUCODE_DEVICE std::tuple eval(T t) const; +#endif - void getLineCoefs(T& a, T& b, T& c) const; + GPUd() void getLineCoefs(T& a, T& b, T& c) const; /** check if XY interval is seen by the circle. * The tolerance parameter eps is interpreted as a fraction of the interval @@ -69,9 +74,9 @@ class IntervalXY * y = yc + dy*t * with 0& circle, T eps) const; + GPUd() bool seenByCircle(const CircleXY& circle, T eps) const; - bool circleCrossParam(const CircleXY& circle, T& t) const; + GPUd() bool circleCrossParam(const CircleXY& circle, T& t) const; /** * check if XY interval is seen by the line defined by other interval @@ -82,12 +87,12 @@ class IntervalXY * y = yc + dy*t * with 0& other, T eps) const; + GPUd() bool seenByLine(const IntervalXY& other, T eps) const; /** * get crossing parameter of 2 intervals */ - bool lineCrossParam(const IntervalXY& other, T& t) const; + GPUd() bool lineCrossParam(const IntervalXY& other, T& t) const; private: T mX, mY; ///< one of edges @@ -97,100 +102,100 @@ class IntervalXY }; template -IntervalXY::IntervalXY(T x, T y, T dx, T dy) : mX(std::move(x)), mY(std::move(y)), mDx(std::move(dx)), mDY(std::move(dy)) +GPUdi() IntervalXY::IntervalXY(T x, T y, T dx, T dy) : mX(x), mY(y), mDx(dx), mDY(dy) { } template -inline T IntervalXY::getX0() const +GPUdi() T IntervalXY::getX0() const { return mX; } template -inline T IntervalXY::getY0() const +GPUdi() T IntervalXY::getY0() const { return mY; } template -inline T& IntervalXY::getX0() +GPUdi() T& IntervalXY::getX0() { return mX; } template -inline T& IntervalXY::getY0() +GPUdi() T& IntervalXY::getY0() { return mY; } template -inline T IntervalXY::getX1() const +GPUdi() T IntervalXY::getX1() const { return mX + mDx; } template -inline T IntervalXY::getY1() const +GPUdi() T IntervalXY::getY1() const { return mY + mDY; } template -inline T IntervalXY::getDX() const +GPUdi() T IntervalXY::getDX() const { return mDx; } template -inline T IntervalXY::getDY() const +GPUdi() T IntervalXY::getDY() const { return mDY; } template -inline T& IntervalXY::getDX() +GPUdi() T& IntervalXY::getDX() { return mDx; } template -inline T& IntervalXY::getDY() +GPUdi() T& IntervalXY::getDY() { return mDY; } template -inline void IntervalXY::setX0(T x0) +GPUdi() void IntervalXY::setX0(T x0) { mX = x0; } template -inline void IntervalXY::setY0(T y0) +GPUdi() void IntervalXY::setY0(T y0) { mY = y0; } template -inline void IntervalXY::setX1(T x1) +GPUdi() void IntervalXY::setX1(T x1) { mDx = x1 - mX; } template -inline void IntervalXY::setY1(T y1) +GPUdi() void IntervalXY::setY1(T y1) { mDY = y1 - mY; } template -inline void IntervalXY::setDX(T dx) +GPUdi() void IntervalXY::setDX(T dx) { mDx = dx; } template -inline void IntervalXY::setDY(T dy) +GPUdi() void IntervalXY::setDY(T dy) { mDY = dy; } template -inline void IntervalXY::setEdges(T x0, T y0, T x1, T y1) +GPUdi() void IntervalXY::setEdges(T x0, T y0, T x1, T y1) { mX = x0; mY = y0; @@ -198,20 +203,22 @@ inline void IntervalXY::setEdges(T x0, T y0, T x1, T y1) mDY = y1 - y0; } +#ifndef GPUCA_GPUCODE_DEVICE template -inline std::tuple IntervalXY::eval(T t) const +GPUdi() std::tuple IntervalXY::eval(T t) const { return {mX + t * mDx, mY + t * mDY}; } +#endif template -inline void IntervalXY::eval(T t, T& x, T& y) const +GPUdi() void IntervalXY::eval(T t, T& x, T& y) const { std::tie(x, y) = eval(t); } template -void IntervalXY::getLineCoefs(T& a, T& b, T& c) const +GPUdi() void IntervalXY::getLineCoefs(T& a, T& b, T& c) const { // convert to line parameters in canonical form: a*x+b*y+c = 0 c = mX * mDY - mY * mDx; @@ -230,7 +237,7 @@ void IntervalXY::getLineCoefs(T& a, T& b, T& c) const } template -bool IntervalXY::seenByCircle(const CircleXY& circle, T eps) const +GPUdi() bool IntervalXY::seenByCircle(const CircleXY& circle, T eps) const { T dx0 = mX - circle.xC; T dy0 = mY - circle.yC; @@ -251,7 +258,7 @@ bool IntervalXY::seenByCircle(const CircleXY& circle, T eps) const return (d02 - rC2) * (d12 - rC2) < 0; } template -bool IntervalXY::circleCrossParam(const CircleXY& circle, T& t) const +GPUdi() bool IntervalXY::circleCrossParam(const CircleXY& circle, T& t) const { const T dx = mX - circle.xC; const T dy = mY - circle.yC; @@ -271,7 +278,7 @@ bool IntervalXY::circleCrossParam(const CircleXY& circle, T& t) const } template -bool IntervalXY::seenByLine(const IntervalXY& other, T eps) const +GPUdi() bool IntervalXY::seenByLine(const IntervalXY& other, T eps) const { T a, b, c; // find equation of the line a*x+b*y+c = 0 other.getLineCoefs(a, b, c); @@ -282,7 +289,7 @@ bool IntervalXY::seenByLine(const IntervalXY& other, T eps) const } template -bool IntervalXY::lineCrossParam(const IntervalXY& other, T& t) const +GPUdi() bool IntervalXY::lineCrossParam(const IntervalXY& other, T& t) const { // tolerance constexpr float eps = 1.e-9f; diff --git a/Common/MathUtils/include/MathUtils/detail/trigonometric.h b/Common/MathUtils/include/MathUtils/detail/trigonometric.h index aa8a0dc891bb7..36c8eee7e90e4 100644 --- a/Common/MathUtils/include/MathUtils/detail/trigonometric.h +++ b/Common/MathUtils/include/MathUtils/detail/trigonometric.h @@ -16,10 +16,10 @@ #define MATHUTILS_INCLUDE_MATHUTILS_DETAIL_TRIGONOMETRIC_H_ #ifndef GPUCA_GPUCODE_DEVICE -#include #include #include #endif +#include "GPUCommonArray.h" #include "GPUCommonDef.h" #include "GPUCommonMath.h" #include "CommonConstants/MathConstants.h" @@ -32,7 +32,7 @@ namespace detail { template -GPUdi() T to02Pi(T phi) +GPUhdi() T to02Pi(T phi) { T res = phi; // ensure angle in [0:2pi] for the input in [-pi:pi] or [0:pi] @@ -43,7 +43,7 @@ GPUdi() T to02Pi(T phi) } template -GPUdi() void bringTo02Pi(T& phi) +GPUhdi() void bringTo02Pi(T& phi) { phi = to02Pi(phi); } @@ -69,7 +69,7 @@ inline void bringTo02PiGen(T& phi) } template -inline T toPMPi(T phi) +GPUhdi() T toPMPi(T phi) { T res = phi; // ensure angle in [-pi:pi] for the input in [-pi:pi] or [0:pi] @@ -82,7 +82,7 @@ inline T toPMPi(T phi) } template -inline void bringToPMPi(T& phi) +GPUhdi() void bringToPMPi(T& phi) { phi = toPMPi(phi); } @@ -108,13 +108,13 @@ inline void bringToPMPiGen(T& phi) } template -GPUdi() void sincos(T ang, GPUgeneric() T& s, GPUgeneric() T& c) +GPUhdi() void sincos(T ang, GPUgeneric() T& s, GPUgeneric() T& c) { return o2::gpu::GPUCommonMath::SinCos(ang, s, c); } template <> -GPUdi() void sincos(double ang, GPUgeneric() double& s, GPUgeneric() double& c) +GPUhdi() void sincos(double ang, GPUgeneric() double& s, GPUgeneric() double& c) { return o2::gpu::GPUCommonMath::SinCosd(ang, s, c); } @@ -138,39 +138,43 @@ inline std::tuple rotateZ(T xL, T yL, T snAlp, T csAlp) } template -inline std::tuple rotateZInv(T xG, T yG, T snAlp, T csAlp) +GPUhdi() std::tuple rotateZInv(T xG, T yG, T snAlp, T csAlp) { // inverse 2D rotation of the point by angle alpha (global to local) return rotateZ(xG, yG, -snAlp, csAlp); } -template -inline void rotateZ(T xL, T yL, T& xG, T& yG, T snAlp, T csAlp) -{ - std::tie(xG, yG) = rotateZ(xL, yL, snAlp, csAlp); -} +#endif template -inline void rotateZ(std::array& xy, T alpha) +GPUhdi() void rotateZ(gpu::gpustd::array& xy, T alpha) { // transforms vector in tracking frame alpha to global frame - auto [sin, cos] = sincos(alpha); + T sin, cos; + sincos(alpha, sin, cos); const T x = xy[0]; xy[0] = x * cos - xy[1] * sin; xy[1] = x * sin + xy[1] * cos; } template -inline void rotateZInv(T xG, T yG, T& xL, T& yL, T snAlp, T csAlp) +GPUhdi() void rotateZ(T xL, T yL, T& xG, T& yG, T snAlp, T csAlp) +{ + xG = xL * csAlp - yL * snAlp; + yG = xL * snAlp + yL * csAlp; + ; + // std::tie(xG, yG) = rotateZ(xL, yL, snAlp, csAlp); // must not use std:: - versions on GPU +} + +template +GPUhdi() void rotateZInv(T xG, T yG, T& xL, T& yL, T snAlp, T csAlp) { // inverse 2D rotation of the point by angle alpha (global to local) rotateZ(xG, yG, xL, yL, -snAlp, csAlp); } -#endif - template -inline int angle2Sector(T phi) +GPUhdi() int angle2Sector(T phi) { // convert angle to sector ID, phi can be either in 0:2pi or -pi:pi convention int sect = phi * o2::constants::math::Rad2Deg / o2::constants::math::SectorSpanDeg; @@ -181,7 +185,7 @@ inline int angle2Sector(T phi) } template -inline T sector2Angle(int sect) +GPUhdi() T sector2Angle(int sect) { // convert sector to its angle center, in -pi:pi convention T ang = o2::constants::math::SectorSpanRad * (0.5f + sect); @@ -190,7 +194,7 @@ inline T sector2Angle(int sect) } template -inline T angle2Alpha(T phi) +GPUhdi() T angle2Alpha(T phi) { // convert angle to its sector alpha return sector2Angle(angle2Sector(phi)); From 469a8510cbf271f983d07672e61d5aefe33e0e97 Mon Sep 17 00:00:00 2001 From: David Rohr Date: Wed, 9 Dec 2020 20:33:44 +0100 Subject: [PATCH 4/9] GPU: Add a dummy stream logger implementation for GPU in order not to fail on LOG(...) << ... --- GPU/Common/GPUCommonLogger.h | 14 +++++++++++++- 1 file changed, 13 insertions(+), 1 deletion(-) diff --git a/GPU/Common/GPUCommonLogger.h b/GPU/Common/GPUCommonLogger.h index 6c0c851a1c82c..103d0638f0c53 100644 --- a/GPU/Common/GPUCommonLogger.h +++ b/GPU/Common/GPUCommonLogger.h @@ -20,7 +20,19 @@ #define LOGP(...) #elif defined(GPUCA_GPUCODE_DEVICE) -#define LOG(...) static_assert("LOG(...) << ... unsupported in GPU code"); +#include "GPUCommonDef.h" +namespace o2::gpu::detail +{ +struct DummyLogger { + template + GPUd() DummyLogger& operator<<(Args... args) + { + return *this; + } +}; +} // namespace o2::gpu::detail +#define LOG(...) o2::gpu::detail::DummyLogger() +//#define LOG(...) static_assert(false, "LOG(...) << ... unsupported in GPU code"); #define LOGF(type, string, ...) \ { \ printf(string "\n", ##__VA_ARGS__); \ From 0098fd80f8c8b23859c0a27344d4c613e8ac1689 Mon Sep 17 00:00:00 2001 From: David Rohr Date: Wed, 9 Dec 2020 20:34:16 +0100 Subject: [PATCH 5/9] GPU: Add gpu::gpustd::array to forward std::array on host / mimic as array on GPU --- .../TrackParametrization.h | 14 +++---- .../TrackParametrizationWithError.h | 8 ++-- .../ReconstructionDataFormats/TrackUtils.h | 8 ++-- .../src/TrackParametrization.cxx | 18 ++++----- .../src/TrackParametrizationWithError.cxx | 20 +++++----- GPU/Common/CMakeLists.txt | 1 + GPU/Common/GPUCommonArray.h | 37 +++++++++++++++++++ 7 files changed, 71 insertions(+), 35 deletions(-) create mode 100644 GPU/Common/GPUCommonArray.h diff --git a/DataFormats/Reconstruction/include/ReconstructionDataFormats/TrackParametrization.h b/DataFormats/Reconstruction/include/ReconstructionDataFormats/TrackParametrization.h index d9ffc0be528ce..a1184146c0430 100644 --- a/DataFormats/Reconstruction/include/ReconstructionDataFormats/TrackParametrization.h +++ b/DataFormats/Reconstruction/include/ReconstructionDataFormats/TrackParametrization.h @@ -28,10 +28,11 @@ #include "GPUCommonDef.h" #include "GPUCommonRtypes.h" #include "GPUCommonMath.h" +#include "GPUCommonArray.h" +#include "GPUROOTCartesianFwd.h" #ifndef __OPENCL__ #include -#include #include #include #include @@ -50,9 +51,6 @@ #include "ReconstructionDataFormats/TrackUtils.h" -//Forward declarations, since we cannot include the headers if we eventually want to use track.h on GPU -#include "GPUROOTCartesianFwd.h" - namespace o2 { template @@ -122,9 +120,9 @@ class TrackParametrization public: using value_t = value_T; - using dim2_t = std::array; - using dim3_t = std::array; - using params_t = std::array; + using dim2_t = gpu::gpustd::array; + using dim3_t = gpu::gpustd::array; + using params_t = gpu::gpustd::array; static_assert(std::is_floating_point_v); @@ -188,7 +186,7 @@ class TrackParametrization math_utils::Point3D getXYZGlo() const; void getXYZGlo(dim3_t& xyz) const; bool getPxPyPzGlo(dim3_t& pxyz) const; - bool getPosDirGlo(std::array& posdirp) const; + bool getPosDirGlo(gpu::gpustd::array& posdirp) const; // methods for track params estimate at other point bool getYZAt(value_t xk, value_t b, value_t& y, value_t& z) const; diff --git a/DataFormats/Reconstruction/include/ReconstructionDataFormats/TrackParametrizationWithError.h b/DataFormats/Reconstruction/include/ReconstructionDataFormats/TrackParametrizationWithError.h index 9910623f25f24..01af5b2cf0cf9 100644 --- a/DataFormats/Reconstruction/include/ReconstructionDataFormats/TrackParametrizationWithError.h +++ b/DataFormats/Reconstruction/include/ReconstructionDataFormats/TrackParametrizationWithError.h @@ -37,12 +37,12 @@ class TrackParametrizationWithError : public TrackParametrization static_assert(std::is_floating_point_v); public: - using covMat_t = std::array; + using covMat_t = gpu::gpustd::array; TrackParametrizationWithError(); - TrackParametrizationWithError(value_t x, value_t alpha, const params_t& par, const std::array& cov, int charge = 1); + TrackParametrizationWithError(value_t x, value_t alpha, const params_t& par, const covMat_t& cov, int charge = 1); TrackParametrizationWithError(const dim3_t& xyz, const dim3_t& pxpypz, - const std::array& cv, int sign, bool sectorAlpha = true); + const gpu::gpustd::array& cv, int sign, bool sectorAlpha = true); TrackParametrizationWithError(const TrackParametrizationWithError& src) = default; TrackParametrizationWithError(TrackParametrizationWithError&& src) = default; @@ -70,7 +70,7 @@ class TrackParametrizationWithError : public TrackParametrization value_t getCovarElem(int i, int j) const; value_t getDiagError2(int i) const; - bool getCovXYZPxPyPzGlo(std::array& c) const; + bool getCovXYZPxPyPzGlo(gpu::gpustd::array& c) const; void print() const; #ifndef GPUCA_ALIGPUCODE diff --git a/DataFormats/Reconstruction/include/ReconstructionDataFormats/TrackUtils.h b/DataFormats/Reconstruction/include/ReconstructionDataFormats/TrackUtils.h index 1ff1f87489aac..737687b2257bb 100644 --- a/DataFormats/Reconstruction/include/ReconstructionDataFormats/TrackUtils.h +++ b/DataFormats/Reconstruction/include/ReconstructionDataFormats/TrackUtils.h @@ -17,9 +17,9 @@ #define INCLUDE_RECONSTRUCTIONDATAFORMATS_TRACKUTILS_H_ #include "GPUCommonRtypes.h" +#include "GPUCommonArray.h" #ifndef __OPENCL__ -#include #include #endif @@ -33,13 +33,13 @@ namespace track // helper function template value_T BetheBlochSolid(value_T bg, value_T rho = 2.33, value_T kp1 = 0.20, value_T kp2 = 3.00, value_T meanI = 173e-9, - value_T meanZA = 0.49848); + value_T meanZA = 0.49848); template -void g3helx3(value_T qfield, value_T step, std::array& vect); +void g3helx3(value_T qfield, value_T step, gpu::gpustd::array& vect); //____________________________________________________ template -void g3helx3(value_T qfield, value_T step, std::array& vect) +void g3helx3(value_T qfield, value_T step, gpu::gpustd::array& vect) { /****************************************************************** * * diff --git a/DataFormats/Reconstruction/src/TrackParametrization.cxx b/DataFormats/Reconstruction/src/TrackParametrization.cxx index 17fd10838210f..c3a663cc0e84d 100644 --- a/DataFormats/Reconstruction/src/TrackParametrization.cxx +++ b/DataFormats/Reconstruction/src/TrackParametrization.cxx @@ -130,7 +130,7 @@ bool TrackParametrization::getPxPyPzGlo(dim3_t& pxyz) const //____________________________________________________ template -bool TrackParametrization::getPosDirGlo(std::array& posdirp) const +bool TrackParametrization::getPosDirGlo(gpu::gpustd::array& posdirp) const { // fill vector with lab x,y,z,px/p,py/p,pz/p,p,sinAlpha,cosAlpha value_t ptI = std::fabs(getQ2Pt()); @@ -229,7 +229,7 @@ bool TrackParametrization::propagateParamTo(value_t xk, const dim3_t& b step *= std::sqrt(1.f + getTgl() * getTgl()); // // get the track x,y,z,px/p,py/p,pz/p,p,sinAlpha,cosAlpha in the Global System - std::array vecLab{0.f}; + gpu::gpustd::array vecLab{0.f}; if (!getPosDirGlo(vecLab)) { return false; } @@ -248,13 +248,13 @@ bool TrackParametrization::propagateParamTo(value_t xk, const dim3_t& b costet = b[2] / bb; sintet = bt / bb; } - std::array vect{costet * cosphi * vecLab[0] + costet * sinphi * vecLab[1] - sintet * vecLab[2], - -sinphi * vecLab[0] + cosphi * vecLab[1], - sintet * cosphi * vecLab[0] + sintet * sinphi * vecLab[1] + costet * vecLab[2], - costet * cosphi * vecLab[3] + costet * sinphi * vecLab[4] - sintet * vecLab[5], - -sinphi * vecLab[3] + cosphi * vecLab[4], - sintet * cosphi * vecLab[3] + sintet * sinphi * vecLab[4] + costet * vecLab[5], - vecLab[6]}; + gpu::gpustd::array vect{costet * cosphi * vecLab[0] + costet * sinphi * vecLab[1] - sintet * vecLab[2], + -sinphi * vecLab[0] + cosphi * vecLab[1], + sintet * cosphi * vecLab[0] + sintet * sinphi * vecLab[1] + costet * vecLab[2], + costet * cosphi * vecLab[3] + costet * sinphi * vecLab[4] - sintet * vecLab[5], + -sinphi * vecLab[3] + cosphi * vecLab[4], + sintet * cosphi * vecLab[3] + sintet * sinphi * vecLab[4] + costet * vecLab[5], + vecLab[6]}; // Do the helix step value_t q = getCharge(); diff --git a/DataFormats/Reconstruction/src/TrackParametrizationWithError.cxx b/DataFormats/Reconstruction/src/TrackParametrizationWithError.cxx index f02b78ff2b05a..a52b773b476d4 100644 --- a/DataFormats/Reconstruction/src/TrackParametrizationWithError.cxx +++ b/DataFormats/Reconstruction/src/TrackParametrizationWithError.cxx @@ -262,7 +262,7 @@ bool TrackParametrizationWithError::propagateToDCA(const o2::dataformat //______________________________________________________________ template TrackParametrizationWithError::TrackParametrizationWithError(const dim3_t& xyz, const dim3_t& pxpypz, - const std::array& cv, int charge, bool sectorAlpha) + const gpu::gpustd::array& cv, int charge, bool sectorAlpha) { // construct track param and covariance from kinematics and lab errors @@ -454,7 +454,7 @@ bool TrackParametrizationWithError::propagateTo(value_t xk, const dim3_ step *= std::sqrt(1.f + this->getTgl() * this->getTgl()); // // get the track x,y,z,px/p,py/p,pz/p,p,sinAlpha,cosAlpha in the Global System - std::array vecLab{0.f}; + gpu::gpustd::array vecLab{0.f}; if (!this->getPosDirGlo(vecLab)) { return false; } @@ -521,13 +521,13 @@ bool TrackParametrizationWithError::propagateTo(value_t xk, const dim3_ costet = b[2] / bb; sintet = bt / bb; } - std::array vect{costet * cosphi * vecLab[0] + costet * sinphi * vecLab[1] - sintet * vecLab[2], - -sinphi * vecLab[0] + cosphi * vecLab[1], - sintet * cosphi * vecLab[0] + sintet * sinphi * vecLab[1] + costet * vecLab[2], - costet * cosphi * vecLab[3] + costet * sinphi * vecLab[4] - sintet * vecLab[5], - -sinphi * vecLab[3] + cosphi * vecLab[4], - sintet * cosphi * vecLab[3] + sintet * sinphi * vecLab[4] + costet * vecLab[5], - vecLab[6]}; + gpu::gpustd::array vect{costet * cosphi * vecLab[0] + costet * sinphi * vecLab[1] - sintet * vecLab[2], + -sinphi * vecLab[0] + cosphi * vecLab[1], + sintet * cosphi * vecLab[0] + sintet * sinphi * vecLab[1] + costet * vecLab[2], + costet * cosphi * vecLab[3] + costet * sinphi * vecLab[4] - sintet * vecLab[5], + -sinphi * vecLab[3] + cosphi * vecLab[4], + sintet * cosphi * vecLab[3] + sintet * sinphi * vecLab[4] + costet * vecLab[5], + vecLab[6]}; // Do the helix step value_t sgn = this->getSign(); @@ -1004,7 +1004,7 @@ bool TrackParametrizationWithError::correctForMaterial(value_t x2x0, va //______________________________________________________________ template -bool TrackParametrizationWithError::getCovXYZPxPyPzGlo(std::array& cv) const +bool TrackParametrizationWithError::getCovXYZPxPyPzGlo(gpu::gpustd::array& cv) const { //--------------------------------------------------------------------- // This function returns the global covariance matrix of the track params diff --git a/GPU/Common/CMakeLists.txt b/GPU/Common/CMakeLists.txt index f3a8502796a19..57f0f083e9f8f 100644 --- a/GPU/Common/CMakeLists.txt +++ b/GPU/Common/CMakeLists.txt @@ -18,6 +18,7 @@ set(HDRS_INSTALL GPUCommonLogger.h GPUCommonMath.h GPUCommonRtypes.h + GPUCommonArray.h GPUCommonTransform3D.h GPUDef.h GPUDefConstantsAndSettings.h diff --git a/GPU/Common/GPUCommonArray.h b/GPU/Common/GPUCommonArray.h new file mode 100644 index 0000000000000..c70443bfc13a9 --- /dev/null +++ b/GPU/Common/GPUCommonArray.h @@ -0,0 +1,37 @@ +// Copyright CERN and copyright holders of ALICE O2. This software is +// distributed under the terms of the GNU General Public License v3 (GPL +// Version 3), copied verbatim in the file "COPYING". +// +// See http://alice-o2.web.cern.ch/license for full licensing information. +// +// In applying this license CERN does not waive the privileges and immunities +// granted to it by virtue of its status as an Intergovernmental Organization +// or submit itself to any jurisdiction. + +/// \file GPUCommonArray.h +/// \author David Rohr + +#ifndef GPUCOMMONFAIRARRAY_H +#define GPUCOMMONFAIRARRAY_H + +#ifndef GPUCA_GPUCODE_DEVICE +#include +#endif + +#include "GPUCommonDef.h" +namespace o2::gpu::gpustd +{ +#ifdef GPUCA_GPUCODE_DEVICE +template +struct array { + GPUd() T& operator[](size_t i) { return m_internal_V__[i]; } + GPUd() const T& operator[](size_t i) const { return m_internal_V__[i]; } + T m_internal_V__[N]; +}; +#else +template +using array = std::array; +#endif +} // namespace o2::gpu::gpustd + +#endif From 7de784435cb9b71a2b923a9522c7c7da6d6a5136 Mon Sep 17 00:00:00 2001 From: David Rohr Date: Wed, 9 Dec 2020 21:22:22 +0100 Subject: [PATCH 6/9] Hide all stuff that uses ROOT from GPU compilation --- .../Reconstruction/src/TrackParametrizationWithError.cxx | 6 +++++- 1 file changed, 5 insertions(+), 1 deletion(-) diff --git a/DataFormats/Reconstruction/src/TrackParametrizationWithError.cxx b/DataFormats/Reconstruction/src/TrackParametrizationWithError.cxx index a52b773b476d4..23e26b1298122 100644 --- a/DataFormats/Reconstruction/src/TrackParametrizationWithError.cxx +++ b/DataFormats/Reconstruction/src/TrackParametrizationWithError.cxx @@ -27,10 +27,10 @@ #include "ReconstructionDataFormats/Vertex.h" #include "ReconstructionDataFormats/DCA.h" #include -#include "Math/SMatrix.h" #ifndef GPUCA_GPUCODE_DEVICE #include +#include "Math/SMatrix.h" #endif #ifndef GPUCA_ALIGPUCODE @@ -686,6 +686,8 @@ typename TrackParametrizationWithError::value_t TrackParametrizationWit return (d * (szz * d - sdz * z) + z * (sdd * z - d * sdz)) / det; } +#ifndef GPUCA_GPUCODE_DEVICE // Disable function relying on ROOT SMatrix on GPU + //______________________________________________ template void TrackParametrizationWithError::buildCombinedCovMatrix(const TrackParametrizationWithError& rhs, MatrixDSym5& cov) const @@ -831,6 +833,8 @@ bool TrackParametrizationWithError::update(const TrackParametrizationWi return update(rhs, covI); } +#endif + //______________________________________________ template bool TrackParametrizationWithError::update(const dim2_t& p, const dim3_t& cov) From a71928f127be1e239eb5b89461d0ff2c7f268c2b Mon Sep 17 00:00:00 2001 From: David Rohr Date: Wed, 9 Dec 2020 21:25:39 +0100 Subject: [PATCH 7/9] Move PID constants from static class member to a namespace, since static constexpr is not supported by GPU --- .../include/ReconstructionDataFormats/PID.h | 105 ++++++++++-------- 1 file changed, 58 insertions(+), 47 deletions(-) diff --git a/DataFormats/Reconstruction/include/ReconstructionDataFormats/PID.h b/DataFormats/Reconstruction/include/ReconstructionDataFormats/PID.h index 25f010b20d2e4..ace4dcec33077 100644 --- a/DataFormats/Reconstruction/include/ReconstructionDataFormats/PID.h +++ b/DataFormats/Reconstruction/include/ReconstructionDataFormats/PID.h @@ -24,11 +24,55 @@ namespace track { namespace o2cp = o2::constants::physics; +namespace pid_constants // GPUs currently cannot have static constexpr array members +{ +typedef uint8_t ID; +static constexpr ID NIDsTot = 14; +GPUconstexpr() const char* sNames[NIDsTot + 1] = ///< defined particle names + {"Electron", "Muon", "Pion", "Kaon", "Proton", "Deuteron", "Triton", "He3", "Alpha", + "Pion0", "Photon", "K0", "Lambda", "HyperTriton", + nullptr}; + +GPUconstexpr() const float sMasses[NIDsTot] = ///< defined particle masses + {o2cp::MassElectron, o2cp::MassMuon, o2cp::MassPionCharged, o2cp::MassKaonCharged, + o2cp::MassProton, o2cp::MassDeuteron, o2cp::MassTriton, o2cp::MassHelium3, + o2cp::MassAlpha, o2cp::MassPionNeutral, o2cp::MassPhoton, + o2cp::MassKaonNeutral, o2cp::MassLambda, o2cp::MassHyperTriton}; + +GPUconstexpr() const float sMasses2[NIDsTot] = ///< defined particle masses^2 + {o2cp::MassElectron * o2cp::MassElectron, + o2cp::MassMuon* o2cp::MassMuon, + o2cp::MassPionCharged* o2cp::MassPionCharged, + o2cp::MassKaonCharged* o2cp::MassKaonCharged, + o2cp::MassProton* o2cp::MassProton, + o2cp::MassDeuteron* o2cp::MassDeuteron, + o2cp::MassTriton* o2cp::MassTriton, + o2cp::MassHelium3* o2cp::MassHelium3, + o2cp::MassAlpha* o2cp::MassAlpha, + o2cp::MassPionNeutral* o2cp::MassPionNeutral, + o2cp::MassPhoton* o2cp::MassPhoton, + o2cp::MassKaonNeutral* o2cp::MassKaonNeutral, + o2cp::MassLambda* o2cp::MassLambda, + o2cp::MassHyperTriton* o2cp::MassHyperTriton}; + +GPUconstexpr() const float sMasses2Z[NIDsTot] = ///< defined particle masses / Z + {o2cp::MassElectron, o2cp::MassMuon, + o2cp::MassPionCharged, o2cp::MassKaonCharged, + o2cp::MassProton, o2cp::MassDeuteron, + o2cp::MassTriton, o2cp::MassHelium3 / 2., + o2cp::MassAlpha / 2., + 0, 0, 0, 0, o2cp::MassHyperTriton}; + +GPUconstexpr() const int sCharges[NIDsTot] = ///< defined particle charges + {1, 1, 1, 1, 1, 1, 1, 2, 2, + 0, 0, 0, 0, 1}; +} // namespace pid_constants + class PID { public: // particle identifiers, continuos starting from 0 - typedef uint8_t ID; + typedef pid_constants::ID ID; static constexpr ID Electron = 0; static constexpr ID Muon = 1; @@ -52,7 +96,8 @@ class PID static constexpr ID HyperTriton = 13; static constexpr ID FirstExt = PI0; static constexpr ID LastExt = HyperTriton; - static constexpr ID NIDsTot = LastExt + 1; ///< total number of defined IDs + static constexpr ID NIDsTot = pid_constants::NIDsTot; ///< total number of defined IDs + static_assert(NIDsTot == LastExt + 1, "Incorrect NIDsTot, please update!"); PID() = default; PID(ID id) : mID(id) {} @@ -66,13 +111,18 @@ class PID float getMass() const { return getMass(mID); } float getMass2Z() const { return getMass2Z(mID); } int getCharge() const { return getCharge(mID); } - const char* getName() const { return getName(mID); } - static constexpr const char* getName(ID id) { return sNames[id]; } - static constexpr float getMass(ID id) { return sMasses[id]; } - static constexpr float getMass2(ID id) { return sMasses2[id]; } - static constexpr float getMass2Z(ID id) { return sMasses2Z[id]; } - static constexpr int getCharge(ID id) { return sCharges[id]; } + static float getMass(ID id) { return pid_constants::sMasses[id]; } + static float getMass2(ID id) { return pid_constants::sMasses2[id]; } + static float getMass2Z(ID id) { return pid_constants::sMasses2Z[id]; } + static int getCharge(ID id) { return pid_constants::sCharges[id]; } +#ifndef GPUCA_GPUCODE_DEVICE + const char* getName() const + { + return getName(mID); + } + static const char* getName(ID id) { return pid_constants::sNames[id]; } +#endif private: ID mID = Pion; @@ -88,45 +138,6 @@ class PID return id > LastExt ? id : sameStr(name, sNames[id]) ? id : nameToID(name, id + 1); } - static constexpr const char* sNames[NIDsTot + 1] = ///< defined particle names - {"Electron", "Muon", "Pion", "Kaon", "Proton", "Deuteron", "Triton", "He3", "Alpha", - "Pion0", "Photon", "K0", "Lambda", "HyperTriton", - nullptr}; - - static constexpr const float sMasses[NIDsTot] = ///< defined particle masses - {o2cp::MassElectron, o2cp::MassMuon, o2cp::MassPionCharged, o2cp::MassKaonCharged, - o2cp::MassProton, o2cp::MassDeuteron, o2cp::MassTriton, o2cp::MassHelium3, - o2cp::MassAlpha, o2cp::MassPionNeutral, o2cp::MassPhoton, - o2cp::MassKaonNeutral, o2cp::MassLambda, o2cp::MassHyperTriton}; - - static constexpr const float sMasses2[NIDsTot] = ///< defined particle masses^2 - {o2cp::MassElectron * o2cp::MassElectron, - o2cp::MassMuon* o2cp::MassMuon, - o2cp::MassPionCharged* o2cp::MassPionCharged, - o2cp::MassKaonCharged* o2cp::MassKaonCharged, - o2cp::MassProton* o2cp::MassProton, - o2cp::MassDeuteron* o2cp::MassDeuteron, - o2cp::MassTriton* o2cp::MassTriton, - o2cp::MassHelium3* o2cp::MassHelium3, - o2cp::MassAlpha* o2cp::MassAlpha, - o2cp::MassPionNeutral* o2cp::MassPionNeutral, - o2cp::MassPhoton* o2cp::MassPhoton, - o2cp::MassKaonNeutral* o2cp::MassKaonNeutral, - o2cp::MassLambda* o2cp::MassLambda, - o2cp::MassHyperTriton* o2cp::MassHyperTriton}; - - static constexpr const float sMasses2Z[NIDsTot] = ///< defined particle masses / Z - {o2cp::MassElectron, o2cp::MassMuon, - o2cp::MassPionCharged, o2cp::MassKaonCharged, - o2cp::MassProton, o2cp::MassDeuteron, - o2cp::MassTriton, o2cp::MassHelium3 / 2., - o2cp::MassAlpha / 2., - 0, 0, 0, 0, o2cp::MassHyperTriton}; - - static constexpr const int sCharges[NIDsTot] = ///< defined particle charges - {1, 1, 1, 1, 1, 1, 1, 2, 2, - 0, 0, 0, 0, 1}; - ClassDefNV(PID, 2); }; } // namespace track From 4e37ae8c6faf40ebfcbc2a7e848a716736e601ce Mon Sep 17 00:00:00 2001 From: David Rohr Date: Wed, 9 Dec 2020 22:56:48 +0100 Subject: [PATCH 8/9] GPU: Make TrackPar / TrackParCov compile for GPU --- Common/MathUtils/include/MathUtils/Utils.h | 3 +- .../include/MathUtils/detail/IntervalXY.h | 17 +- .../include/MathUtils/detail/trigonometric.h | 13 +- .../include/ReconstructionDataFormats/DCA.h | 31 +- .../include/ReconstructionDataFormats/PID.h | 41 +-- .../TrackParametrization.h | 298 +++++++++--------- .../TrackParametrizationWithError.h | 154 ++++----- .../ReconstructionDataFormats/TrackUtils.h | 27 +- .../ReconstructionDataFormats/Vertex.h | 106 +++---- DataFormats/Reconstruction/src/DCA.cxx | 4 +- DataFormats/Reconstruction/src/PID.cxx | 5 - .../src/TrackParametrization.cxx | 192 +++++------ .../src/TrackParametrizationWithError.cxx | 206 ++++++------ DataFormats/Reconstruction/src/Vertex.cxx | 2 +- Detectors/Base/src/Propagator.cxx | 2 +- GPU/Common/GPUCommonLogger.h | 16 +- .../Base/GPUReconstructionIncludesDevice.h | 4 + GPU/GPUTracking/Base/cuda/CMakeLists.txt | 1 + GPU/GPUTracking/Base/hip/CMakeLists.txt | 1 + GPU/GPUTracking/Base/opencl2/CMakeLists.txt | 1 + GPU/GPUTracking/Standalone/CMakeLists.txt | 3 + 21 files changed, 580 insertions(+), 547 deletions(-) diff --git a/Common/MathUtils/include/MathUtils/Utils.h b/Common/MathUtils/include/MathUtils/Utils.h index 0f7b2e16185c8..e621226a6614f 100644 --- a/Common/MathUtils/include/MathUtils/Utils.h +++ b/Common/MathUtils/include/MathUtils/Utils.h @@ -108,11 +108,12 @@ GPUdi() void sincos(float ang, float& s, float& c) { detail::sincos(ang, s, c); } - +#ifndef __OPENCL__ GPUdi() void sincosd(double ang, double& s, double& c) { detail::sincos(ang, s, c); } +#endif GPUdi() void rotateZ(float xL, float yL, float& xG, float& yG, float snAlp, float csAlp) { diff --git a/Common/MathUtils/include/MathUtils/detail/IntervalXY.h b/Common/MathUtils/include/MathUtils/detail/IntervalXY.h index 2d77c57f2bb21..2a65634ff0778 100644 --- a/Common/MathUtils/include/MathUtils/detail/IntervalXY.h +++ b/Common/MathUtils/include/MathUtils/detail/IntervalXY.h @@ -18,6 +18,7 @@ #include "GPUCommonDef.h" #include "GPUCommonRtypes.h" +#include "GPUCommonMath.h" #ifndef GPUCA_GPUCODE_DEVICE #include #include @@ -205,7 +206,7 @@ GPUdi() void IntervalXY::setEdges(T x0, T y0, T x1, T y1) #ifndef GPUCA_GPUCODE_DEVICE template -GPUdi() std::tuple IntervalXY::eval(T t) const +std::tuple IntervalXY::eval(T t) const { return {mX + t * mDx, mY + t * mDY}; } @@ -214,7 +215,8 @@ GPUdi() std::tuple IntervalXY::eval(T t) const template GPUdi() void IntervalXY::eval(T t, T& x, T& y) const { - std::tie(x, y) = eval(t); + x = mX + t * mDx; + y = mY + t * mDY; } template @@ -269,11 +271,11 @@ GPUdi() bool IntervalXY::circleCrossParam(const CircleXY& circle, T& t) co if (det < 0.f) { return false; } - det = std::sqrt(det); + det = gpu::CAMath::Sqrt(det); const T t0 = -b - det; const T t1 = -b + det; // select the one closer to [0:1] interval - t = (std::fabs(t0 - 0.5f) < std::fabs(t1 - 0.5f)) ? t0 : t1; + t = (gpu::CAMath::Abs(t0 - 0.5f) < gpu::CAMath::Abs(t1 - 0.5f)) ? t0 : t1; return true; } @@ -282,8 +284,9 @@ GPUdi() bool IntervalXY::seenByLine(const IntervalXY& other, T eps) const { T a, b, c; // find equation of the line a*x+b*y+c = 0 other.getLineCoefs(a, b, c); - const auto [x0, y0] = eval(-eps); - const auto [x1, y1] = eval(1.f + eps); + T x0, y0, x1, y1; + eval(-eps, x0, y0); + eval(1.f + eps, x0, y0); return (a * x0 + b * y0 + c) * (a * x1 + b * y1 + c) < 0; } @@ -296,7 +299,7 @@ GPUdi() bool IntervalXY::lineCrossParam(const IntervalXY& other, T& t) con const T dx = other.getX0() - mX; const T dy = other.getY0() - mY; const T det = -mDx * other.getY0() + mDY * other.getX0(); - if (std::fabs(det) < eps) { + if (gpu::CAMath::Abs(det) < eps) { return false; // parallel } t = (-dx * other.getY0() + dy * other.getX0()) / det; diff --git a/Common/MathUtils/include/MathUtils/detail/trigonometric.h b/Common/MathUtils/include/MathUtils/detail/trigonometric.h index 36c8eee7e90e4..af6eb69274783 100644 --- a/Common/MathUtils/include/MathUtils/detail/trigonometric.h +++ b/Common/MathUtils/include/MathUtils/detail/trigonometric.h @@ -107,17 +107,24 @@ inline void bringToPMPiGen(T& phi) phi = toPMPiGen(phi); } +#ifdef __OPENCL__ // TODO: get rid of that stupid workaround for OpenCL template address spaces template -GPUhdi() void sincos(T ang, GPUgeneric() T& s, GPUgeneric() T& c) +GPUhdi() void sincos(T ang, float& s, float& c) +{ + return o2::gpu::GPUCommonMath::SinCos(ang, s, c); +} +#else +template +GPUhdi() void sincos(T ang, T& s, T& c) { return o2::gpu::GPUCommonMath::SinCos(ang, s, c); } - template <> -GPUhdi() void sincos(double ang, GPUgeneric() double& s, GPUgeneric() double& c) +GPUhdi() void sincos(double ang, double& s, double& c) { return o2::gpu::GPUCommonMath::SinCosd(ang, s, c); } +#endif #ifndef GPUCA_GPUCODE_DEVICE diff --git a/DataFormats/Reconstruction/include/ReconstructionDataFormats/DCA.h b/DataFormats/Reconstruction/include/ReconstructionDataFormats/DCA.h index ff6c357d4ca9a..ff4167b59fdbb 100644 --- a/DataFormats/Reconstruction/include/ReconstructionDataFormats/DCA.h +++ b/DataFormats/Reconstruction/include/ReconstructionDataFormats/DCA.h @@ -11,12 +11,11 @@ #ifndef ALICEO2_DCA_H #define ALICEO2_DCA_H +#include "GPUCommonDef.h" #include "GPUCommonRtypes.h" +#include "GPUCommonArray.h" -#ifndef __OPENCL__ -#include -#endif -#ifndef GPUCA_ALIGPUCODE +#ifndef GPUCA_GPUCODE_DEVICE #include #endif @@ -32,14 +31,14 @@ class DCA { public: - DCA() = default; + GPUdDefault() DCA() = default; - DCA(float y, float z, float syy = 0.f, float syz = 0.f, float szz = 0.f) + GPUd() DCA(float y, float z, float syy = 0.f, float syz = 0.f, float szz = 0.f) { set(y, z, syy, syz, szz); } - void set(float y, float z, float syy, float syz, float szz) + GPUd() void set(float y, float z, float syy, float syz, float szz) { mY = y; mZ = z; @@ -48,30 +47,32 @@ class DCA mCov[2] = szz; } - void set(float y, float z) + GPUd() void set(float y, float z) { mY = y; mZ = z; } - auto getY() const { return mY; } - auto getZ() const { return mZ; } - auto getSigmaY2() const { return mCov[0]; } - auto getSigmaYZ() const { return mCov[1]; } - auto getSigmaZ2() const { return mCov[2]; } - const auto& getCovariance() const { return mCov; } + GPUd() auto getY() const { return mY; } + GPUd() auto getZ() const { return mZ; } + GPUd() auto getSigmaY2() const { return mCov[0]; } + GPUd() auto getSigmaYZ() const { return mCov[1]; } + GPUd() auto getSigmaZ2() const { return mCov[2]; } + GPUd() const auto& getCovariance() const { return mCov; } void print() const; private: float mY = 0.f; float mZ = 0.f; - std::array mCov; ///< s2y, syz, s2z + gpu::gpustd::array mCov; ///< s2y, syz, s2z ClassDefNV(DCA, 1); }; +#ifndef GPUCA_GPUCODE_DEVICE std::ostream& operator<<(std::ostream& os, const DCA& d); +#endif } // namespace dataformats } // namespace o2 diff --git a/DataFormats/Reconstruction/include/ReconstructionDataFormats/PID.h b/DataFormats/Reconstruction/include/ReconstructionDataFormats/PID.h index ace4dcec33077..da25dc910e778 100644 --- a/DataFormats/Reconstruction/include/ReconstructionDataFormats/PID.h +++ b/DataFormats/Reconstruction/include/ReconstructionDataFormats/PID.h @@ -15,6 +15,7 @@ #ifndef ALICEO2_track_PID_H_ #define ALICEO2_track_PID_H_ +#include "GPUCommonDef.h" #include "GPUCommonRtypes.h" #include "CommonConstants/PhysicsConstants.h" @@ -99,44 +100,46 @@ class PID static constexpr ID NIDsTot = pid_constants::NIDsTot; ///< total number of defined IDs static_assert(NIDsTot == LastExt + 1, "Incorrect NIDsTot, please update!"); - PID() = default; - PID(ID id) : mID(id) {} - PID(const char* name); - PID(const PID& src) = default; - PID& operator=(const PID& src) = default; + GPUdDefault() PID() = default; + GPUd() PID(ID id) : mID(id) {} + GPUd() PID(const char* name); + GPUdDefault() PID(const PID& src) = default; + GPUdDefault() PID& operator=(const PID& src) = default; - ID getID() const { return mID; } - operator ID() const { return getID(); } + GPUd() ID getID() const { return mID; } + GPUd() operator ID() const { return getID(); } - float getMass() const { return getMass(mID); } - float getMass2Z() const { return getMass2Z(mID); } - int getCharge() const { return getCharge(mID); } + GPUd() float getMass() const { return getMass(mID); } + GPUd() float getMass2Z() const { return getMass2Z(mID); } + GPUd() int getCharge() const { return getCharge(mID); } - static float getMass(ID id) { return pid_constants::sMasses[id]; } - static float getMass2(ID id) { return pid_constants::sMasses2[id]; } - static float getMass2Z(ID id) { return pid_constants::sMasses2Z[id]; } - static int getCharge(ID id) { return pid_constants::sCharges[id]; } + GPUd() static float getMass(ID id) { return pid_constants::sMasses[id]; } + GPUd() static float getMass2(ID id) { return pid_constants::sMasses2[id]; } + GPUd() static float getMass2Z(ID id) { return pid_constants::sMasses2Z[id]; } + GPUd() static int getCharge(ID id) { return pid_constants::sCharges[id]; } #ifndef GPUCA_GPUCODE_DEVICE - const char* getName() const + GPUd() const char* getName() const { return getName(mID); } - static const char* getName(ID id) { return pid_constants::sNames[id]; } + GPUd() static const char* getName(ID id) { return pid_constants::sNames[id]; } #endif private: ID mID = Pion; // are 2 strings equal ? (trick from Giulio) - inline static constexpr bool sameStr(char const* x, char const* y) + GPUdi() static constexpr bool sameStr(char const* x, char const* y) { return !*x && !*y ? true : /* default */ (*x == *y && sameStr(x + 1, y + 1)); } - inline static constexpr ID nameToID(char const* name, ID id) +#ifndef GPUCA_GPUCODE_DEVICE + GPUdi() static constexpr ID nameToID(char const* name, ID id) { - return id > LastExt ? id : sameStr(name, sNames[id]) ? id : nameToID(name, id + 1); + return id > LastExt ? id : sameStr(name, pid_constants::sNames[id]) ? id : nameToID(name, id + 1); } +#endif ClassDefNV(PID, 2); }; diff --git a/DataFormats/Reconstruction/include/ReconstructionDataFormats/TrackParametrization.h b/DataFormats/Reconstruction/include/ReconstructionDataFormats/TrackParametrization.h index a1184146c0430..77bd1c372ed80 100644 --- a/DataFormats/Reconstruction/include/ReconstructionDataFormats/TrackParametrization.h +++ b/DataFormats/Reconstruction/include/ReconstructionDataFormats/TrackParametrization.h @@ -103,14 +103,14 @@ constexpr float kCY2max = 100 * 100, // SigmaY<=100cm kCalcdEdxAuto = -999.f; // value indicating request for dedx calculation // access to covariance matrix by row and column -constexpr int CovarMap[kNParams][kNParams] = {{0, 1, 3, 6, 10}, - {1, 2, 4, 7, 11}, - {3, 4, 5, 8, 12}, - {6, 7, 8, 9, 13}, - {10, 11, 12, 13, 14}}; +GPUconstexpr() int CovarMap[kNParams][kNParams] = {{0, 1, 3, 6, 10}, + {1, 2, 4, 7, 11}, + {3, 4, 5, 8, 12}, + {6, 7, 8, 9, 13}, + {10, 11, 12, 13, 14}}; // access to covariance matrix diagonal elements -constexpr int DiagMap[kNParams] = {0, 2, 5, 9, 14}; +GPUconstexpr() int DiagMap[kNParams] = {0, 2, 5, 9, 14}; constexpr float HugeF = o2::constants::math::VeryBig; @@ -124,100 +124,102 @@ class TrackParametrization using dim3_t = gpu::gpustd::array; using params_t = gpu::gpustd::array; +#ifndef GPUCA_GPUCODE_DEVICE static_assert(std::is_floating_point_v); +#endif - TrackParametrization() = default; - TrackParametrization(value_t x, value_t alpha, const params_t& par, int charge = 1); - TrackParametrization(const dim3_t& xyz, const dim3_t& pxpypz, int charge, bool sectorAlpha = true); - TrackParametrization(const TrackParametrization&) = default; - TrackParametrization(TrackParametrization&&) = default; - TrackParametrization& operator=(const TrackParametrization& src) = default; - TrackParametrization& operator=(TrackParametrization&& src) = default; - ~TrackParametrization() = default; - - const value_t* getParams() const; - value_t getParam(int i) const; - value_t getX() const; - value_t getAlpha() const; - value_t getY() const; - value_t getZ() const; - value_t getSnp() const; - value_t getTgl() const; - value_t getQ2Pt() const; - value_t getCharge2Pt() const; - int getAbsCharge() const; - PID getPID() const; - void setPID(const PID pid); + GPUdDefault() TrackParametrization() = default; + GPUd() TrackParametrization(value_t x, value_t alpha, const params_t& par, int charge = 1); + GPUd() TrackParametrization(const dim3_t& xyz, const dim3_t& pxpypz, int charge, bool sectorAlpha = true); + GPUdDefault() TrackParametrization(const TrackParametrization&) = default; + GPUdDefault() TrackParametrization(TrackParametrization&&) = default; + GPUdDefault() TrackParametrization& operator=(const TrackParametrization& src) = default; + GPUdDefault() TrackParametrization& operator=(TrackParametrization&& src) = default; + GPUdDefault() ~TrackParametrization() = default; + + GPUd() const value_t* getParams() const; + GPUd() value_t getParam(int i) const; + GPUd() value_t getX() const; + GPUd() value_t getAlpha() const; + GPUd() value_t getY() const; + GPUd() value_t getZ() const; + GPUd() value_t getSnp() const; + GPUd() value_t getTgl() const; + GPUd() value_t getQ2Pt() const; + GPUd() value_t getCharge2Pt() const; + GPUd() int getAbsCharge() const; + GPUd() PID getPID() const; + GPUd() void setPID(const PID pid); /// calculate cos^2 and cos of track direction in rphi-tracking - value_t getCsp2() const; - value_t getCsp() const; - - void setX(value_t v); - void setParam(value_t v, int i); - void setAlpha(value_t v); - void setY(value_t v); - void setZ(value_t v); - void setSnp(value_t v); - void setTgl(value_t v); - void setQ2Pt(value_t v); - void setAbsCharge(int q); + GPUd() value_t getCsp2() const; + GPUd() value_t getCsp() const; + + GPUd() void setX(value_t v); + GPUd() void setParam(value_t v, int i); + GPUd() void setAlpha(value_t v); + GPUd() void setY(value_t v); + GPUd() void setZ(value_t v); + GPUd() void setSnp(value_t v); + GPUd() void setTgl(value_t v); + GPUd() void setQ2Pt(value_t v); + GPUd() void setAbsCharge(int q); // derived getters - bool getXatLabR(value_t r, value_t& x, value_t bz, DirType dir = DirAuto) const; - void getCircleParamsLoc(value_t bz, o2::math_utils::CircleXY& circle) const; - void getCircleParams(value_t bz, o2::math_utils::CircleXY& circle, value_t& sna, value_t& csa) const; - void getLineParams(o2::math_utils::IntervalXY& line, value_t& sna, value_t& csa) const; - value_t getCurvature(value_t b) const; - int getCharge() const; - int getSign() const; - value_t getPhi() const; - value_t getPhiPos() const; - - value_t getPtInv() const; - value_t getP2Inv() const; - value_t getP2() const; - value_t getPInv() const; - value_t getP() const; - value_t getPt() const; - - value_t getTheta() const; - value_t getEta() const; - math_utils::Point3D getXYZGlo() const; - void getXYZGlo(dim3_t& xyz) const; - bool getPxPyPzGlo(dim3_t& pxyz) const; - bool getPosDirGlo(gpu::gpustd::array& posdirp) const; + GPUd() bool getXatLabR(value_t r, value_t& x, value_t bz, DirType dir = DirAuto) const; + GPUd() void getCircleParamsLoc(value_t bz, o2::math_utils::CircleXY& circle) const; + GPUd() void getCircleParams(value_t bz, o2::math_utils::CircleXY& circle, value_t& sna, value_t& csa) const; + GPUd() void getLineParams(o2::math_utils::IntervalXY& line, value_t& sna, value_t& csa) const; + GPUd() value_t getCurvature(value_t b) const; + GPUd() int getCharge() const; + GPUd() int getSign() const; + GPUd() value_t getPhi() const; + GPUd() value_t getPhiPos() const; + + GPUd() value_t getPtInv() const; + GPUd() value_t getP2Inv() const; + GPUd() value_t getP2() const; + GPUd() value_t getPInv() const; + GPUd() value_t getP() const; + GPUd() value_t getPt() const; + + GPUd() value_t getTheta() const; + GPUd() value_t getEta() const; + GPUd() math_utils::Point3D getXYZGlo() const; + GPUd() void getXYZGlo(dim3_t& xyz) const; + GPUd() bool getPxPyPzGlo(dim3_t& pxyz) const; + GPUd() bool getPosDirGlo(gpu::gpustd::array& posdirp) const; // methods for track params estimate at other point - bool getYZAt(value_t xk, value_t b, value_t& y, value_t& z) const; - value_t getZAt(value_t xk, value_t b) const; - value_t getYAt(value_t xk, value_t b) const; - math_utils::Point3D getXYZGloAt(value_t xk, value_t b, bool& ok) const; + GPUd() bool getYZAt(value_t xk, value_t b, value_t& y, value_t& z) const; + GPUd() value_t getZAt(value_t xk, value_t b) const; + GPUd() value_t getYAt(value_t xk, value_t b) const; + GPUd() math_utils::Point3D getXYZGloAt(value_t xk, value_t b, bool& ok) const; // parameters manipulation - bool correctForELoss(value_t xrho, value_t mass, bool anglecorr = false, value_t dedx = kCalcdEdxAuto); - bool rotateParam(value_t alpha); - bool propagateParamTo(value_t xk, value_t b); - bool propagateParamTo(value_t xk, const dim3_t& b); + GPUd() bool correctForELoss(value_t xrho, value_t mass, bool anglecorr = false, value_t dedx = kCalcdEdxAuto); + GPUd() bool rotateParam(value_t alpha); + GPUd() bool propagateParamTo(value_t xk, value_t b); + GPUd() bool propagateParamTo(value_t xk, const dim3_t& b); - bool propagateParamToDCA(const math_utils::Point3D& vtx, value_t b, dim2_t* dca = nullptr, value_t maxD = 999.f); + GPUd() bool propagateParamToDCA(const math_utils::Point3D& vtx, value_t b, dim2_t* dca = nullptr, value_t maxD = 999.f); - void invertParam(); + GPUd() void invertParam(); - bool isValid() const; - void invalidate(); + GPUd() bool isValid() const; + GPUd() void invalidate(); - uint16_t getUserField() const; - void setUserField(uint16_t v); + GPUd() uint16_t getUserField() const; + GPUd() void setUserField(uint16_t v); - void printParam() const; + GPUd() void printParam() const; #ifndef GPUCA_ALIGPUCODE std::string asString() const; #endif protected: - void updateParam(value_t delta, int i); - void updateParams(const value_t delta[kNParams]); + GPUd() void updateParam(value_t delta, int i); + GPUd() void updateParams(const value_t delta[kNParams]); private: // @@ -234,100 +236,102 @@ class TrackParametrization //____________________________________________________________ template -inline TrackParametrization::TrackParametrization(value_t x, value_t alpha, const params_t& par, int charge) - : mX{x}, mAlpha{alpha}, mAbsCharge{char(std::abs(charge))} +GPUdi() TrackParametrization::TrackParametrization(value_t x, value_t alpha, const params_t& par, int charge) + : mX{x}, mAlpha{alpha}, mAbsCharge{char(gpu::CAMath::Abs(charge))} { // explicit constructor - std::copy(par.begin(), par.end(), mP); + for (int i = 0; i < kNParams; i++) { + mP[i] = par[i]; + } } //____________________________________________________________ template -inline const typename TrackParametrization::value_t* TrackParametrization::getParams() const +GPUdi() const typename TrackParametrization::value_t* TrackParametrization::getParams() const { return mP; } //____________________________________________________________ template -inline typename TrackParametrization::value_t TrackParametrization::getParam(int i) const +GPUdi() typename TrackParametrization::value_t TrackParametrization::getParam(int i) const { return mP[i]; } //____________________________________________________________ template -inline typename TrackParametrization::value_t TrackParametrization::getX() const +GPUdi() typename TrackParametrization::value_t TrackParametrization::getX() const { return mX; } //____________________________________________________________ template -inline typename TrackParametrization::value_t TrackParametrization::getAlpha() const +GPUdi() typename TrackParametrization::value_t TrackParametrization::getAlpha() const { return mAlpha; } //____________________________________________________________ template -inline typename TrackParametrization::value_t TrackParametrization::getY() const +GPUdi() typename TrackParametrization::value_t TrackParametrization::getY() const { return mP[kY]; } //____________________________________________________________ template -inline typename TrackParametrization::value_t TrackParametrization::getZ() const +GPUdi() typename TrackParametrization::value_t TrackParametrization::getZ() const { return mP[kZ]; } //____________________________________________________________ template -inline typename TrackParametrization::value_t TrackParametrization::getSnp() const +GPUdi() typename TrackParametrization::value_t TrackParametrization::getSnp() const { return mP[kSnp]; } //____________________________________________________________ template -inline typename TrackParametrization::value_t TrackParametrization::getTgl() const +GPUdi() typename TrackParametrization::value_t TrackParametrization::getTgl() const { return mP[kTgl]; } //____________________________________________________________ template -inline typename TrackParametrization::value_t TrackParametrization::getQ2Pt() const +GPUdi() typename TrackParametrization::value_t TrackParametrization::getQ2Pt() const { return mP[kQ2Pt]; } //____________________________________________________________ template -inline typename TrackParametrization::value_t TrackParametrization::getCharge2Pt() const +GPUdi() typename TrackParametrization::value_t TrackParametrization::getCharge2Pt() const { return mAbsCharge ? mP[kQ2Pt] : 0.f; } //____________________________________________________________ template -inline int TrackParametrization::getAbsCharge() const +GPUdi() int TrackParametrization::getAbsCharge() const { return mAbsCharge; } //____________________________________________________________ template -inline PID TrackParametrization::getPID() const +GPUdi() PID TrackParametrization::getPID() const { return mPID; } //____________________________________________________________ template -inline void TrackParametrization::setPID(const PID pid) +GPUdi() void TrackParametrization::setPID(const PID pid) { mPID = pid; setAbsCharge(pid.getCharge()); @@ -335,7 +339,7 @@ inline void TrackParametrization::setPID(const PID pid) //____________________________________________________________ template -inline typename TrackParametrization::value_t TrackParametrization::getCsp2() const +GPUdi() typename TrackParametrization::value_t TrackParametrization::getCsp2() const { const value_t csp2 = (1.f - mP[kSnp]) * (1.f + mP[kSnp]); return csp2 > o2::constants::math::Almost0 ? csp2 : o2::constants::math::Almost0; @@ -343,88 +347,88 @@ inline typename TrackParametrization::value_t TrackParametrization -inline typename TrackParametrization::value_t TrackParametrization::getCsp() const +GPUdi() typename TrackParametrization::value_t TrackParametrization::getCsp() const { - return std::sqrt(getCsp2()); + return gpu::CAMath::Sqrt(getCsp2()); } //____________________________________________________________ template -inline void TrackParametrization::setX(value_t v) +GPUdi() void TrackParametrization::setX(value_t v) { mX = v; } //____________________________________________________________ template -inline void TrackParametrization::setParam(value_t v, int i) +GPUdi() void TrackParametrization::setParam(value_t v, int i) { mP[i] = v; } //____________________________________________________________ template -inline void TrackParametrization::setAlpha(value_t v) +GPUdi() void TrackParametrization::setAlpha(value_t v) { mAlpha = v; } //____________________________________________________________ template -inline void TrackParametrization::setY(value_t v) +GPUdi() void TrackParametrization::setY(value_t v) { mP[kY] = v; } //____________________________________________________________ template -inline void TrackParametrization::setZ(value_t v) +GPUdi() void TrackParametrization::setZ(value_t v) { mP[kZ] = v; } //____________________________________________________________ template -inline void TrackParametrization::setSnp(value_t v) +GPUdi() void TrackParametrization::setSnp(value_t v) { mP[kSnp] = v; } //____________________________________________________________ template -inline void TrackParametrization::setTgl(value_t v) +GPUdi() void TrackParametrization::setTgl(value_t v) { mP[kTgl] = v; } //____________________________________________________________ template -inline void TrackParametrization::setQ2Pt(value_t v) +GPUdi() void TrackParametrization::setQ2Pt(value_t v) { mP[kQ2Pt] = v; } //____________________________________________________________ template -inline void TrackParametrization::setAbsCharge(int q) +GPUdi() void TrackParametrization::setAbsCharge(int q) { - mAbsCharge = std::abs(q); + mAbsCharge = gpu::CAMath::Abs(q); } //_______________________________________________________ template -inline void TrackParametrization::getCircleParamsLoc(value_t bz, o2::math_utils::CircleXY& c) const +GPUdi() void TrackParametrization::getCircleParamsLoc(value_t bz, o2::math_utils::CircleXY& c) const { // get circle params in track local frame, for straight line just set to local coordinates c.rC = getCurvature(bz); // treat as straight track if sagitta between the vertex and middle of TPC is below 0.01 cm constexpr value_t MinSagitta = 0.01f, TPCMidR = 160.f, MinCurv = 8 * MinSagitta / (TPCMidR * TPCMidR); - if (std::abs(c.rC) > MinCurv) { + if (gpu::CAMath::Abs(c.rC) > MinCurv) { c.rC = 1.f / getCurvature(bz); - value_t sn = getSnp(), cs = std::sqrt((1.f - sn) * (1.f + sn)); + value_t sn = getSnp(), cs = gpu::CAMath::Sqrt((1.f - sn) * (1.f + sn)); c.xC = getX() - sn * c.rC; // center in tracking c.yC = getY() + cs * c.rC; // frame. Note: r is signed!!! - c.rC = std::abs(c.rC); + c.rC = gpu::CAMath::Abs(c.rC); } else { c.rC = 0.f; // signal straight line c.xC = getX(); @@ -434,7 +438,7 @@ inline void TrackParametrization::getCircleParamsLoc(value_t bz, o2::ma //_______________________________________________________ template -inline void TrackParametrization::getCircleParams(value_t bz, o2::math_utils::CircleXY& c, value_t& sna, value_t& csa) const +GPUdi() void TrackParametrization::getCircleParams(value_t bz, o2::math_utils::CircleXY& c, value_t& sna, value_t& csa) const { // get circle params in loc and lab frame, for straight line just set to global coordinates getCircleParamsLoc(bz, c); @@ -444,69 +448,69 @@ inline void TrackParametrization::getCircleParams(value_t bz, o2::math_ //_______________________________________________________ template -inline void TrackParametrization::getLineParams(o2::math_utils::IntervalXY& ln, value_t& sna, value_t& csa) const +GPUdi() void TrackParametrization::getLineParams(o2::math_utils::IntervalXY& ln, value_t& sna, value_t& csa) const { // get line parameterization as { x = x0 + xSlp*t, y = y0 + ySlp*t } o2::math_utils::detail::sincos(getAlpha(), sna, csa); o2::math_utils::detail::rotateZ(getX(), getY(), ln.getX0(), ln.getY0(), sna, csa); // reference point in global frame - value_t snp = getSnp(), csp = std::sqrt((1.f - snp) * (1.f + snp)); + value_t snp = getSnp(), csp = gpu::CAMath::Sqrt((1.f - snp) * (1.f + snp)); ln.setDX(csp * csa - snp * sna); ln.setDY(snp * csa + csp * sna); } //____________________________________________________________ template -inline typename TrackParametrization::value_t TrackParametrization::getCurvature(value_t b) const +GPUdi() typename TrackParametrization::value_t TrackParametrization::getCurvature(value_t b) const { return mAbsCharge ? mP[kQ2Pt] * b * o2::constants::math::B2C : 0.; } //____________________________________________________________ template -inline int TrackParametrization::getCharge() const +GPUdi() int TrackParametrization::getCharge() const { return getSign() > 0 ? mAbsCharge : -mAbsCharge; } //____________________________________________________________ template -inline int TrackParametrization::getSign() const +GPUdi() int TrackParametrization::getSign() const { return mAbsCharge ? (mP[kQ2Pt] > 0.f ? 1 : -1) : 0; } //_______________________________________________________ template -inline typename TrackParametrization::value_t TrackParametrization::getPhi() const +GPUdi() typename TrackParametrization::value_t TrackParametrization::getPhi() const { // track pt direction phi (in 0:2pi range) - value_t phi = std::asin(getSnp()) + getAlpha(); + value_t phi = gpu::CAMath::ASin(getSnp()) + getAlpha(); math_utils::detail::bringTo02Pi(phi); return phi; } //_______________________________________________________ template -inline typename TrackParametrization::value_t TrackParametrization::getPhiPos() const +GPUdi() typename TrackParametrization::value_t TrackParametrization::getPhiPos() const { // angle of track position (in -pi:pi range) - value_t phi = std::atan2(getY(), getX()) + getAlpha(); + value_t phi = gpu::CAMath::ATan2(getY(), getX()) + getAlpha(); math_utils::detail::bringTo02Pi(phi); return phi; } //____________________________________________________________ template -inline typename TrackParametrization::value_t TrackParametrization::getPtInv() const +GPUdi() typename TrackParametrization::value_t TrackParametrization::getPtInv() const { // return the inverted track pT - const value_t ptInv = std::fabs(mP[kQ2Pt]); + const value_t ptInv = gpu::CAMath::Abs(mP[kQ2Pt]); return (mAbsCharge > 1) ? ptInv / mAbsCharge : ptInv; } //____________________________________________________________ template -inline typename TrackParametrization::value_t TrackParametrization::getP2Inv() const +GPUdi() typename TrackParametrization::value_t TrackParametrization::getP2Inv() const { // return the inverted track momentum^2 const value_t p2 = mP[kQ2Pt] * mP[kQ2Pt] / (1.f + getTgl() * getTgl()); @@ -515,7 +519,7 @@ inline typename TrackParametrization::value_t TrackParametrization -inline typename TrackParametrization::value_t TrackParametrization::getP2() const +GPUdi() typename TrackParametrization::value_t TrackParametrization::getP2() const { // return the track momentum^2 const value_t p2inv = getP2Inv(); @@ -524,16 +528,16 @@ inline typename TrackParametrization::value_t TrackParametrization -inline typename TrackParametrization::value_t TrackParametrization::getPInv() const +GPUdi() typename TrackParametrization::value_t TrackParametrization::getPInv() const { // return the inverted track momentum^2 - const value_t pInv = std::fabs(mP[kQ2Pt]) / std::sqrt(1.f + getTgl() * getTgl()); + const value_t pInv = gpu::CAMath::Abs(mP[kQ2Pt]) / gpu::CAMath::Sqrt(1.f + getTgl() * getTgl()); return (mAbsCharge > 1) ? pInv / mAbsCharge : pInv; } //____________________________________________________________ template -inline typename TrackParametrization::value_t TrackParametrization::getP() const +GPUdi() typename TrackParametrization::value_t TrackParametrization::getP() const { // return the track momentum const value_t pInv = getPInv(); @@ -542,10 +546,10 @@ inline typename TrackParametrization::value_t TrackParametrization -inline typename TrackParametrization::value_t TrackParametrization::getPt() const +GPUdi() typename TrackParametrization::value_t TrackParametrization::getPt() const { // return the track transverse momentum - value_t ptI = std::fabs(mP[kQ2Pt]); + value_t ptI = gpu::CAMath::Abs(mP[kQ2Pt]); if (mAbsCharge > 1) { ptI /= mAbsCharge; } @@ -554,34 +558,34 @@ inline typename TrackParametrization::value_t TrackParametrization -inline typename TrackParametrization::value_t TrackParametrization::getTheta() const +GPUdi() typename TrackParametrization::value_t TrackParametrization::getTheta() const { - return constants::math::PIHalf - std::atan(mP[3]); + return constants::math::PIHalf - gpu::CAMath::ATan(mP[3]); } //____________________________________________________________ template -inline typename TrackParametrization::value_t TrackParametrization::getEta() const +GPUdi() typename TrackParametrization::value_t TrackParametrization::getEta() const { - return -std::log(std::tan(0.5f * getTheta())); + return -gpu::CAMath::Log(gpu::CAMath::Tan(0.5f * getTheta())); } //_______________________________________________________ template -inline math_utils::Point3D::value_t> TrackParametrization::getXYZGlo() const +GPUdi() math_utils::Point3D::value_t> TrackParametrization::getXYZGlo() const { #ifndef GPUCA_ALIGPUCODE return math_utils::Rotation2D(getAlpha())(math_utils::Point3D(getX(), getY(), getZ())); #else // mockup on GPU without ROOT float sina, cosa; - o2::gpu::CAMath::SinCos(getAlpha(), sina, cosa); + gpu::CAMath::SinCos(getAlpha(), sina, cosa); return math_utils::Point3D(cosa * getX() + sina * getY(), cosa * getY() - sina * getX(), getZ()); #endif } //_______________________________________________________ template -inline void TrackParametrization::getXYZGlo(dim3_t& xyz) const +GPUdi() void TrackParametrization::getXYZGlo(dim3_t& xyz) const { // track coordinates in lab frame xyz[0] = getX(); @@ -592,7 +596,7 @@ inline void TrackParametrization::getXYZGlo(dim3_t& xyz) const //_______________________________________________________ template -inline math_utils::Point3D::value_t> TrackParametrization::getXYZGloAt(value_t xk, value_t b, bool& ok) const +GPUdi() math_utils::Point3D::value_t> TrackParametrization::getXYZGloAt(value_t xk, value_t b, bool& ok) const { //---------------------------------------------------------------- // estimate global X,Y,Z in global frame at given X @@ -604,7 +608,7 @@ inline math_utils::Point3D::value_t> Trac return math_utils::Rotation2D(getAlpha())(math_utils::Point3D(xk, y, z)); #else // mockup on GPU without ROOT float sina, cosa; - o2::gpu::CAMath::SinCos(getAlpha(), sina, cosa); + gpu::CAMath::SinCos(getAlpha(), sina, cosa); return math_utils::Point3D(cosa * xk + sina * y, cosa * y - sina * xk, z); #endif } else { @@ -614,40 +618,40 @@ inline math_utils::Point3D::value_t> Trac //____________________________________________________________ template -inline bool TrackParametrization::isValid() const +GPUdi() bool TrackParametrization::isValid() const { return mX != InvalidX; } //____________________________________________________________ template -inline void TrackParametrization::invalidate() +GPUdi() void TrackParametrization::invalidate() { mX = InvalidX; } template -inline uint16_t TrackParametrization::getUserField() const +GPUdi() uint16_t TrackParametrization::getUserField() const { return mUserField; } template -inline void TrackParametrization::setUserField(uint16_t v) +GPUdi() void TrackParametrization::setUserField(uint16_t v) { mUserField = v; } //____________________________________________________________ template -inline void TrackParametrization::updateParam(value_t delta, int i) +GPUdi() void TrackParametrization::updateParam(value_t delta, int i) { mP[i] += delta; } //____________________________________________________________ template -inline void TrackParametrization::updateParams(const value_t delta[kNParams]) +GPUdi() void TrackParametrization::updateParams(const value_t delta[kNParams]) { for (int i = kNParams; i--;) { mP[i] += delta[i]; diff --git a/DataFormats/Reconstruction/include/ReconstructionDataFormats/TrackParametrizationWithError.h b/DataFormats/Reconstruction/include/ReconstructionDataFormats/TrackParametrizationWithError.h index 01af5b2cf0cf9..d6c78eae4f50b 100644 --- a/DataFormats/Reconstruction/include/ReconstructionDataFormats/TrackParametrizationWithError.h +++ b/DataFormats/Reconstruction/include/ReconstructionDataFormats/TrackParametrizationWithError.h @@ -34,82 +34,82 @@ class TrackParametrizationWithError : public TrackParametrization using typename TrackParametrization::dim2_t; using typename TrackParametrization::params_t; +#ifndef GPUCA_GPUCODE_DEVICE static_assert(std::is_floating_point_v); +#endif public: using covMat_t = gpu::gpustd::array; - TrackParametrizationWithError(); - TrackParametrizationWithError(value_t x, value_t alpha, const params_t& par, const covMat_t& cov, int charge = 1); - TrackParametrizationWithError(const dim3_t& xyz, const dim3_t& pxpypz, + GPUd() TrackParametrizationWithError(); + GPUd() TrackParametrizationWithError(value_t x, value_t alpha, const params_t& par, const covMat_t& cov, int charge = 1); + GPUd() TrackParametrizationWithError(const dim3_t& xyz, const dim3_t& pxpypz, const gpu::gpustd::array& cv, int sign, bool sectorAlpha = true); - TrackParametrizationWithError(const TrackParametrizationWithError& src) = default; - TrackParametrizationWithError(TrackParametrizationWithError&& src) = default; - TrackParametrizationWithError& operator=(const TrackParametrizationWithError& src) = default; - TrackParametrizationWithError& operator=(TrackParametrizationWithError&& src) = default; - ~TrackParametrizationWithError() = default; + GPUdDefault() TrackParametrizationWithError(const TrackParametrizationWithError& src) = default; + GPUdDefault() TrackParametrizationWithError(TrackParametrizationWithError&& src) = default; + GPUdDefault() TrackParametrizationWithError& operator=(const TrackParametrizationWithError& src) = default; + GPUdDefault() TrackParametrizationWithError& operator=(TrackParametrizationWithError&& src) = default; + GPUdDefault() ~TrackParametrizationWithError() = default; using TrackParametrization::TrackParametrization; - const value_t* getCov() const; - value_t getSigmaY2() const; - value_t getSigmaZY() const; - value_t getSigmaZ2() const; - value_t getSigmaSnpY() const; - value_t getSigmaSnpZ() const; - value_t getSigmaSnp2() const; - value_t getSigmaTglY() const; - value_t getSigmaTglZ() const; - value_t getSigmaTglSnp() const; - value_t getSigmaTgl2() const; - value_t getSigma1PtY() const; - value_t getSigma1PtZ() const; - value_t getSigma1PtSnp() const; - value_t getSigma1PtTgl() const; - value_t getSigma1Pt2() const; - value_t getCovarElem(int i, int j) const; - value_t getDiagError2(int i) const; - - bool getCovXYZPxPyPzGlo(gpu::gpustd::array& c) const; - - void print() const; + GPUd() const value_t* getCov() const; + GPUd() value_t getSigmaY2() const; + GPUd() value_t getSigmaZY() const; + GPUd() value_t getSigmaZ2() const; + GPUd() value_t getSigmaSnpY() const; + GPUd() value_t getSigmaSnpZ() const; + GPUd() value_t getSigmaSnp2() const; + GPUd() value_t getSigmaTglY() const; + GPUd() value_t getSigmaTglZ() const; + GPUd() value_t getSigmaTglSnp() const; + GPUd() value_t getSigmaTgl2() const; + GPUd() value_t getSigma1PtY() const; + GPUd() value_t getSigma1PtZ() const; + GPUd() value_t getSigma1PtSnp() const; + GPUd() value_t getSigma1PtTgl() const; + GPUd() value_t getSigma1Pt2() const; + GPUd() value_t getCovarElem(int i, int j) const; + GPUd() value_t getDiagError2(int i) const; + + GPUd() bool getCovXYZPxPyPzGlo(gpu::gpustd::array& c) const; + + GPUd() void print() const; #ifndef GPUCA_ALIGPUCODE std::string asString() const; #endif // parameters + covmat manipulation - bool rotate(value_t alpha); - bool propagateTo(value_t xk, value_t b); - bool propagateTo(value_t xk, const dim3_t& b); - bool propagateToDCA(const o2::dataformats::VertexBase& vtx, value_t b, o2::dataformats::DCA* dca = nullptr, value_t maxD = 999.f); - void invert(); + GPUd() bool rotate(value_t alpha); + GPUd() bool propagateTo(value_t xk, value_t b); + GPUd() bool propagateTo(value_t xk, const dim3_t& b); + GPUd() bool propagateToDCA(const o2::dataformats::VertexBase& vtx, value_t b, o2::dataformats::DCA* dca = nullptr, value_t maxD = 999.f); + GPUd() void invert(); - value_t getPredictedChi2(const dim2_t& p, const dim3_t& cov) const; + GPUd() value_t getPredictedChi2(const dim2_t& p, const dim3_t& cov) const; template - value_t getPredictedChi2(const BaseCluster& p) const; - - value_t getPredictedChi2(const TrackParametrizationWithError& rhs) const; + GPUd() value_t getPredictedChi2(const BaseCluster& p) const; void buildCombinedCovMatrix(const TrackParametrizationWithError& rhs, MatrixDSym5& cov) const; value_t getPredictedChi2(const TrackParametrizationWithError& rhs, MatrixDSym5& covToSet) const; + value_t getPredictedChi2(const TrackParametrizationWithError& rhs) const; bool update(const TrackParametrizationWithError& rhs, const MatrixDSym5& covInv); + bool update(const TrackParametrizationWithError& rhs); - bool update(const dim2_t& p, const dim3_t& cov); + GPUd() bool update(const dim2_t& p, const dim3_t& cov); template - bool update(const BaseCluster& p); - - bool update(const TrackParametrizationWithError& rhs); + GPUd() bool update(const BaseCluster& p); - bool correctForMaterial(value_t x2x0, value_t xrho, value_t mass, bool anglecorr = false, value_t dedx = kCalcdEdxAuto); + GPUd() bool correctForMaterial(value_t x2x0, value_t xrho, value_t mass, bool anglecorr = false, value_t dedx = kCalcdEdxAuto); - void resetCovariance(value_t s2 = 0); - void checkCovariance(); - void setCov(value_t v, int i); + GPUd() void resetCovariance(value_t s2 = 0); + GPUd() void checkCovariance(); + GPUd() void setCov(value_t v, int i); - void updateCov(const value_t delta[kCovMatSize]); - void updateCov(value_t delta, int i); + GPUd() void updateCov(const value_t delta[kCovMatSize]); + GPUd() void updateCov(value_t delta, int i); protected: value_t mC[kCovMatSize] = {0.f}; // 15 covariance matrix elements @@ -119,142 +119,144 @@ class TrackParametrizationWithError : public TrackParametrization //__________________________________________________________________________ template -inline TrackParametrizationWithError::TrackParametrizationWithError() : TrackParametrization{} +GPUdi() TrackParametrizationWithError::TrackParametrizationWithError() : TrackParametrization{} { } //__________________________________________________________________________ template -inline TrackParametrizationWithError::TrackParametrizationWithError(value_t x, value_t alpha, const params_t& par, - const std::array& cov, int charge) +GPUdi() TrackParametrizationWithError::TrackParametrizationWithError(value_t x, value_t alpha, const params_t& par, + const covMat_t& cov, int charge) : TrackParametrization{x, alpha, par, charge} { // explicit constructor - std::copy(cov.begin(), cov.end(), mC); + for (int i = 0; i < kCovMatSize; i++) { + mC[i] = cov[i]; + } } //__________________________________________________________________________ template -inline const typename TrackParametrizationWithError::value_t* TrackParametrizationWithError::getCov() const +GPUdi() const typename TrackParametrizationWithError::value_t* TrackParametrizationWithError::getCov() const { return mC; } //__________________________________________________________________________ template -inline typename TrackParametrizationWithError::value_t TrackParametrizationWithError::getSigmaY2() const +GPUdi() typename TrackParametrizationWithError::value_t TrackParametrizationWithError::getSigmaY2() const { return mC[kSigY2]; } //__________________________________________________________________________ template -inline typename TrackParametrizationWithError::value_t TrackParametrizationWithError::getSigmaZY() const +GPUdi() typename TrackParametrizationWithError::value_t TrackParametrizationWithError::getSigmaZY() const { return mC[kSigZY]; } //__________________________________________________________________________ template -inline typename TrackParametrizationWithError::value_t TrackParametrizationWithError::getSigmaZ2() const +GPUdi() typename TrackParametrizationWithError::value_t TrackParametrizationWithError::getSigmaZ2() const { return mC[kSigZ2]; } //__________________________________________________________________________ template -inline typename TrackParametrizationWithError::value_t TrackParametrizationWithError::getSigmaSnpY() const +GPUdi() typename TrackParametrizationWithError::value_t TrackParametrizationWithError::getSigmaSnpY() const { return mC[kSigSnpY]; } //__________________________________________________________________________ template -inline typename TrackParametrizationWithError::value_t TrackParametrizationWithError::getSigmaSnpZ() const +GPUdi() typename TrackParametrizationWithError::value_t TrackParametrizationWithError::getSigmaSnpZ() const { return mC[kSigSnpZ]; } //__________________________________________________________________________ template -inline typename TrackParametrizationWithError::value_t TrackParametrizationWithError::getSigmaSnp2() const +GPUdi() typename TrackParametrizationWithError::value_t TrackParametrizationWithError::getSigmaSnp2() const { return mC[kSigSnp2]; } //__________________________________________________________________________ template -inline typename TrackParametrizationWithError::value_t TrackParametrizationWithError::getSigmaTglY() const +GPUdi() typename TrackParametrizationWithError::value_t TrackParametrizationWithError::getSigmaTglY() const { return mC[kSigTglY]; } //__________________________________________________________________________ template -inline typename TrackParametrizationWithError::value_t TrackParametrizationWithError::getSigmaTglZ() const +GPUdi() typename TrackParametrizationWithError::value_t TrackParametrizationWithError::getSigmaTglZ() const { return mC[kSigTglZ]; } //__________________________________________________________________________ template -inline typename TrackParametrizationWithError::value_t TrackParametrizationWithError::getSigmaTglSnp() const +GPUdi() typename TrackParametrizationWithError::value_t TrackParametrizationWithError::getSigmaTglSnp() const { return mC[kSigTglSnp]; } //__________________________________________________________________________ template -inline typename TrackParametrizationWithError::value_t TrackParametrizationWithError::getSigmaTgl2() const +GPUdi() typename TrackParametrizationWithError::value_t TrackParametrizationWithError::getSigmaTgl2() const { return mC[kSigTgl2]; } //__________________________________________________________________________ template -inline typename TrackParametrizationWithError::value_t TrackParametrizationWithError::getSigma1PtY() const +GPUdi() typename TrackParametrizationWithError::value_t TrackParametrizationWithError::getSigma1PtY() const { return mC[kSigQ2PtY]; } //__________________________________________________________________________ template -inline typename TrackParametrizationWithError::value_t TrackParametrizationWithError::getSigma1PtZ() const +GPUdi() typename TrackParametrizationWithError::value_t TrackParametrizationWithError::getSigma1PtZ() const { return mC[kSigQ2PtZ]; } //__________________________________________________________________________ template -inline typename TrackParametrizationWithError::value_t TrackParametrizationWithError::getSigma1PtSnp() const +GPUdi() typename TrackParametrizationWithError::value_t TrackParametrizationWithError::getSigma1PtSnp() const { return mC[kSigQ2PtSnp]; } //__________________________________________________________________________ template -inline typename TrackParametrizationWithError::value_t TrackParametrizationWithError::getSigma1PtTgl() const +GPUdi() typename TrackParametrizationWithError::value_t TrackParametrizationWithError::getSigma1PtTgl() const { return mC[kSigQ2PtTgl]; } //__________________________________________________________________________ template -inline typename TrackParametrizationWithError::value_t TrackParametrizationWithError::getSigma1Pt2() const +GPUdi() typename TrackParametrizationWithError::value_t TrackParametrizationWithError::getSigma1Pt2() const { return mC[kSigQ2Pt2]; } //__________________________________________________________________________ template -inline typename TrackParametrizationWithError::value_t TrackParametrizationWithError::getCovarElem(int i, int j) const +GPUdi() typename TrackParametrizationWithError::value_t TrackParametrizationWithError::getCovarElem(int i, int j) const { return mC[CovarMap[i][j]]; } //__________________________________________________________________________ template -inline typename TrackParametrizationWithError::value_t TrackParametrizationWithError::getDiagError2(int i) const +GPUdi() typename TrackParametrizationWithError::value_t TrackParametrizationWithError::getDiagError2(int i) const { return mC[DiagMap[i]]; } @@ -262,7 +264,7 @@ inline typename TrackParametrizationWithError::value_t TrackParametriza //__________________________________________________________________________ template template -typename TrackParametrizationWithError::value_t TrackParametrizationWithError::getPredictedChi2(const BaseCluster& p) const +GPUdi() typename TrackParametrizationWithError::value_t TrackParametrizationWithError::getPredictedChi2(const BaseCluster& p) const { const dim2_t pyz = {p.getY(), p.getZ()}; const dim3_t cov = {p.getSigmaY2(), p.getSigmaYZ(), p.getSigmaZ2()}; @@ -272,7 +274,7 @@ typename TrackParametrizationWithError::value_t TrackParametrizationWit //__________________________________________________________________________ template template -bool TrackParametrizationWithError::update(const BaseCluster& p) +GPUdi() bool TrackParametrizationWithError::update(const BaseCluster& p) { const dim2_t pyz = {p.getY(), p.getZ()}; const dim3_t cov = {p.getSigmaY2(), p.getSigmaYZ(), p.getSigmaZ2()}; @@ -281,21 +283,21 @@ bool TrackParametrizationWithError::update(const BaseCluster& p) //__________________________________________________________________________ template -inline void TrackParametrizationWithError::setCov(value_t v, int i) +GPUdi() void TrackParametrizationWithError::setCov(value_t v, int i) { mC[i] = v; } //__________________________________________________________________________ template -inline void TrackParametrizationWithError::updateCov(value_t delta, int i) +GPUdi() void TrackParametrizationWithError::updateCov(value_t delta, int i) { mC[i] += delta; } //__________________________________________________________________________ template -inline void TrackParametrizationWithError::updateCov(const value_t delta[kCovMatSize]) +GPUdi() void TrackParametrizationWithError::updateCov(const value_t delta[kCovMatSize]) { for (int i = kCovMatSize; i--;) { mC[i] += delta[i]; diff --git a/DataFormats/Reconstruction/include/ReconstructionDataFormats/TrackUtils.h b/DataFormats/Reconstruction/include/ReconstructionDataFormats/TrackUtils.h index 737687b2257bb..bfd65642e1b79 100644 --- a/DataFormats/Reconstruction/include/ReconstructionDataFormats/TrackUtils.h +++ b/DataFormats/Reconstruction/include/ReconstructionDataFormats/TrackUtils.h @@ -32,14 +32,14 @@ namespace track { // helper function template -value_T BetheBlochSolid(value_T bg, value_T rho = 2.33, value_T kp1 = 0.20, value_T kp2 = 3.00, value_T meanI = 173e-9, +GPUd() value_T BetheBlochSolid(value_T bg, value_T rho = 2.33, value_T kp1 = 0.20, value_T kp2 = 3.00, value_T meanI = 173e-9, value_T meanZA = 0.49848); template -void g3helx3(value_T qfield, value_T step, gpu::gpustd::array& vect); +GPUd() void g3helx3(value_T qfield, value_T step, gpu::gpustd::array& vect); //____________________________________________________ template -void g3helx3(value_T qfield, value_T step, gpu::gpustd::array& vect) +GPUd() void g3helx3(value_T qfield, value_T step, gpu::gpustd::array& vect) { /****************************************************************** * * @@ -58,8 +58,9 @@ void g3helx3(value_T qfield, value_T step, gpu::gpustd::array& vect) * vect[7](cm,GeV/c) - input/output x, y, z, px/p, py/p ,pz/p, p * * * ******************************************************************/ - +#ifndef GPUCA_GPUCODE_DEVICE static_assert(std::is_floating_point_v); +#endif const int ix = 0, iy = 1, iz = 2, ipx = 3, ipy = 4, ipz = 5, ipp = 6; constexpr value_T kOvSqSix = 0.408248f; // std::sqrt(1./6.); @@ -70,11 +71,11 @@ void g3helx3(value_T qfield, value_T step, gpu::gpustd::array& vect) value_T tet = rho * step; value_T tsint, sintt, sint, cos1t; - if (std::fabs(tet) > 0.03f) { - sint = std::sin(tet); + if (gpu::CAMath::Abs(tet) > 0.03f) { + sint = gpu::CAMath::Sin(tet); sintt = sint / tet; tsint = (tet - sint) / tet; - value_T t = std::sin(0.5f * tet); + value_T t = gpu::CAMath::Sin(0.5f * tet); cos1t = 2 * t * t / tet; } else { tsint = tet * tet / 6.f; @@ -99,8 +100,8 @@ void g3helx3(value_T qfield, value_T step, gpu::gpustd::array& vect) //____________________________________________________ template -value_T BetheBlochSolid(value_T bg, value_T rho, value_T kp1, value_T kp2, value_T meanI, - value_T meanZA) +GPUd() value_T BetheBlochSolid(value_T bg, value_T rho, value_T kp1, value_T kp2, value_T meanI, + value_T meanZA) { // // This is the parameterization of the Bethe-Bloch formula inspired by Geant. @@ -115,7 +116,9 @@ value_T BetheBlochSolid(value_T bg, value_T rho, value_T kp1, value_T kp2, value // The default values for the kp* parameters are for silicon. // The returned value is in [GeV/(g/cm^2)]. // +#ifndef GPUCA_GPUCODE_DEVICE static_assert(std::is_floating_point_v); +#endif constexpr value_T mK = 0.307075e-3f; // [GeV*cm^2/g] constexpr value_T me = 0.511e-3f; // [GeV/c^2] @@ -126,15 +129,15 @@ value_T BetheBlochSolid(value_T bg, value_T rho, value_T kp1, value_T kp2, value //*** Density effect value_T d2 = 0.; - const value_T x = std::log(bg); - const value_T lhwI = std::log(28.816f * 1e-9f * std::sqrt(rho * meanZA) / meanI); + const value_T x = gpu::CAMath::Log(bg); + const value_T lhwI = gpu::CAMath::Log(28.816f * 1e-9f * gpu::CAMath::Sqrt(rho * meanZA) / meanI); if (x > kp2) { d2 = lhwI + x - 0.5f; } else if (x > kp1) { double r = (kp2 - x) / (kp2 - kp1); d2 = lhwI + x - 0.5f + (0.5f - lhwI - kp1) * r * r * r; } - return mK * meanZA * (1 + bg2) / bg2 * (0.5f * std::log(2 * me * bg2 * maxT / (meanI * meanI)) - bg2 / (1 + bg2) - d2); + return mK * meanZA * (1 + bg2) / bg2 * (0.5f * gpu::CAMath::Log(2 * me * bg2 * maxT / (meanI * meanI)) - bg2 / (1 + bg2) - d2); } } // namespace track diff --git a/DataFormats/Reconstruction/include/ReconstructionDataFormats/Vertex.h b/DataFormats/Reconstruction/include/ReconstructionDataFormats/Vertex.h index f4d798f863868..eb8bf67d5d474 100644 --- a/DataFormats/Reconstruction/include/ReconstructionDataFormats/Vertex.h +++ b/DataFormats/Reconstruction/include/ReconstructionDataFormats/Vertex.h @@ -13,13 +13,11 @@ #include "GPUCommonDef.h" #include "GPUCommonMath.h" +#include "GPUCommonArray.h" #include #include "CommonDataFormat/TimeStamp.h" -#ifndef __OPENCL__ -#include -#endif -#ifndef GPUCA_ALIGPUCODE +#ifndef GPUCA_GPUCODE_DEVICE #include #endif @@ -39,51 +37,51 @@ class VertexBase kCovYZ, kCovZZ }; static constexpr int kNCov = 6; - VertexBase() = default; - ~VertexBase() = default; - VertexBase(const math_utils::Point3D& pos, const std::array& cov) : mPos(pos), mCov(cov) + GPUdDefault() VertexBase() = default; + GPUdDefault() ~VertexBase() = default; + GPUd() VertexBase(const math_utils::Point3D& pos, const gpu::gpustd::array& cov) : mPos(pos), mCov(cov) { } -#ifndef GPUCA_ALIGPUCODE +#ifndef GPUCA_GPUCODE_DEVICE void print() const; std::string asString() const; #endif // getting the cartesian coordinates and errors - float getX() const { return mPos.X(); } - float getY() const { return mPos.Y(); } - float getZ() const { return mPos.Z(); } - float getSigmaX2() const { return mCov[kCovXX]; } - float getSigmaY2() const { return mCov[kCovYY]; } - float getSigmaZ2() const { return mCov[kCovZZ]; } - float getSigmaXY() const { return mCov[kCovXY]; } - float getSigmaXZ() const { return mCov[kCovXZ]; } - float getSigmaYZ() const { return mCov[kCovYZ]; } - const std::array& getCov() const { return mCov; } - - math_utils::Point3D getXYZ() const { return mPos; } - math_utils::Point3D& getXYZ() { return mPos; } - - void setX(float x) { mPos.SetX(x); } - void setY(float y) { mPos.SetY(y); } - void setZ(float z) { mPos.SetZ(z); } - - void setXYZ(float x, float y, float z) + GPUd() float getX() const { return mPos.X(); } + GPUd() float getY() const { return mPos.Y(); } + GPUd() float getZ() const { return mPos.Z(); } + GPUd() float getSigmaX2() const { return mCov[kCovXX]; } + GPUd() float getSigmaY2() const { return mCov[kCovYY]; } + GPUd() float getSigmaZ2() const { return mCov[kCovZZ]; } + GPUd() float getSigmaXY() const { return mCov[kCovXY]; } + GPUd() float getSigmaXZ() const { return mCov[kCovXZ]; } + GPUd() float getSigmaYZ() const { return mCov[kCovYZ]; } + GPUd() const gpu::gpustd::array& getCov() const { return mCov; } + + GPUd() math_utils::Point3D getXYZ() const { return mPos; } + GPUd() math_utils::Point3D& getXYZ() { return mPos; } + + GPUd() void setX(float x) { mPos.SetX(x); } + GPUd() void setY(float y) { mPos.SetY(y); } + GPUd() void setZ(float z) { mPos.SetZ(z); } + + GPUd() void setXYZ(float x, float y, float z) { setX(x); setY(y); setZ(z); } - void setPos(const math_utils::Point3D& p) { mPos = p; } - - void setSigmaX2(float v) { mCov[kCovXX] = v; } - void setSigmaY2(float v) { mCov[kCovYY] = v; } - void setSigmaZ2(float v) { mCov[kCovZZ] = v; } - void setSigmaXY(float v) { mCov[kCovXY] = v; } - void setSigmaXZ(float v) { mCov[kCovXZ] = v; } - void setSigmaYZ(float v) { mCov[kCovYZ] = v; } - void setCov(float sxx, float sxy, float syy, float sxz, float syz, float szz) + GPUd() void setPos(const math_utils::Point3D& p) { mPos = p; } + + GPUd() void setSigmaX2(float v) { mCov[kCovXX] = v; } + GPUd() void setSigmaY2(float v) { mCov[kCovYY] = v; } + GPUd() void setSigmaZ2(float v) { mCov[kCovZZ] = v; } + GPUd() void setSigmaXY(float v) { mCov[kCovXY] = v; } + GPUd() void setSigmaXZ(float v) { mCov[kCovXZ] = v; } + GPUd() void setSigmaYZ(float v) { mCov[kCovYZ] = v; } + GPUd() void setCov(float sxx, float sxy, float syy, float sxz, float syz, float szz) { setSigmaX2(sxx); setSigmaY2(syy); @@ -92,11 +90,11 @@ class VertexBase setSigmaXZ(sxz); setSigmaYZ(syz); } - void setCov(const std::array& cov) { mCov = cov; } + GPUd() void setCov(const gpu::gpustd::array& cov) { mCov = cov; } protected: math_utils::Point3D mPos{0., 0., 0.}; ///< cartesian position - std::array mCov{}; ///< errors, see CovElems enum + gpu::gpustd::array mCov{}; ///< errors, see CovElems enum ClassDefNV(VertexBase, 1); }; @@ -115,28 +113,28 @@ class Vertex : public VertexBase FlagsMask = 0xffff }; - Vertex() = default; - ~Vertex() = default; - Vertex(const math_utils::Point3D& pos, const std::array& cov, ushort nCont, float chi2) + GPUdDefault() Vertex() = default; + GPUdDefault() ~Vertex() = default; + GPUd() Vertex(const math_utils::Point3D& pos, const gpu::gpustd::array& cov, ushort nCont, float chi2) : VertexBase(pos, cov), mNContributors(nCont), mChi2(chi2) { } - ushort getNContributors() const { return mNContributors; } - void setNContributors(ushort v) { mNContributors = v; } - void addContributor() { mNContributors++; } + GPUd() ushort getNContributors() const { return mNContributors; } + GPUd() void setNContributors(ushort v) { mNContributors = v; } + GPUd() void addContributor() { mNContributors++; } - ushort getFlags() const { return mBits; } - bool isFlagSet(uint f) const { return mBits & (FlagsMask & f); } - void setFlags(ushort f) { mBits |= FlagsMask & f; } - void resetFrags(ushort f = FlagsMask) { mBits &= ~(FlagsMask & f); } + GPUd() ushort getFlags() const { return mBits; } + GPUd() bool isFlagSet(uint f) const { return mBits & (FlagsMask & f); } + GPUd() void setFlags(ushort f) { mBits |= FlagsMask & f; } + GPUd() void resetFrags(ushort f = FlagsMask) { mBits &= ~(FlagsMask & f); } - void setChi2(float v) { mChi2 = v; } - float getChi2() const { return mChi2; } + GPUd() void setChi2(float v) { mChi2 = v; } + GPUd() float getChi2() const { return mChi2; } - const Stamp& getTimeStamp() const { return mTimeStamp; } - Stamp& getTimeStamp() { return mTimeStamp; } - void setTimeStamp(const Stamp& v) { mTimeStamp = v; } + GPUd() const Stamp& getTimeStamp() const { return mTimeStamp; } + GPUd() Stamp& getTimeStamp() { return mTimeStamp; } + GPUd() void setTimeStamp(const Stamp& v) { mTimeStamp = v; } protected: float mChi2 = 0; ///< chi2 or quality of tracks to vertex attachment @@ -147,7 +145,7 @@ class Vertex : public VertexBase ClassDefNV(Vertex, 3); }; -#ifndef GPUCA_ALIGPUCODE +#ifndef GPUCA_GPUCODE_DEVICE std::ostream& operator<<(std::ostream& os, const o2::dataformats::VertexBase& v); #endif diff --git a/DataFormats/Reconstruction/src/DCA.cxx b/DataFormats/Reconstruction/src/DCA.cxx index d840cefc0ccd9..34d6f714d0feb 100644 --- a/DataFormats/Reconstruction/src/DCA.cxx +++ b/DataFormats/Reconstruction/src/DCA.cxx @@ -19,7 +19,7 @@ namespace o2 namespace dataformats { -#ifndef GPUCA_ALIGPUCODE +#ifndef GPUCA_GPUCODE_DEVICE std::ostream& operator<<(std::ostream& os, const o2::dataformats::DCA& d) { // stream itself @@ -30,7 +30,7 @@ std::ostream& operator<<(std::ostream& os, const o2::dataformats::DCA& d) void DCA::print() const { -#ifndef GPUCA_ALIGPUCODE +#ifndef GPUCA_GPUCODE_DEVICE std::cout << *this << '\n'; #endif } diff --git a/DataFormats/Reconstruction/src/PID.cxx b/DataFormats/Reconstruction/src/PID.cxx index 6b3986de9c7ab..ffe6dab6c778d 100644 --- a/DataFormats/Reconstruction/src/PID.cxx +++ b/DataFormats/Reconstruction/src/PID.cxx @@ -18,11 +18,6 @@ using namespace o2::track; -constexpr const char* PID::sNames[NIDsTot + 1]; -constexpr const float PID::sMasses[NIDsTot]; -constexpr const float PID::sMasses2Z[NIDsTot]; -constexpr const int PID::sCharges[NIDsTot]; - //_______________________________ PID::PID(const char* name) : mID(nameToID(name, First)) { diff --git a/DataFormats/Reconstruction/src/TrackParametrization.cxx b/DataFormats/Reconstruction/src/TrackParametrization.cxx index c3a663cc0e84d..d4a9c8ac2da2e 100644 --- a/DataFormats/Reconstruction/src/TrackParametrization.cxx +++ b/DataFormats/Reconstruction/src/TrackParametrization.cxx @@ -26,8 +26,8 @@ #include "ReconstructionDataFormats/TrackParametrization.h" #include "ReconstructionDataFormats/Vertex.h" #include "ReconstructionDataFormats/DCA.h" +#include #include -#include "Math/SMatrix.h" #ifndef GPUCA_GPUCODE_DEVICE #include @@ -37,14 +37,12 @@ #include #endif -namespace o2 -{ -namespace track -{ +using namespace o2::gpu; +using namespace o2::track; //______________________________________________________________ template -TrackParametrization::TrackParametrization(const dim3_t& xyz, const dim3_t& pxpypz, int charge, bool sectorAlpha) +GPUd() TrackParametrization::TrackParametrization(const dim3_t& xyz, const dim3_t& pxpypz, int charge, bool sectorAlpha) : mX{0.f}, mAlpha{0.f}, mP{0.f} { // construct track param from kinematics @@ -59,9 +57,9 @@ TrackParametrization::TrackParametrization(const dim3_t& xyz, const dim value_t radPos2 = xyz[0] * xyz[0] + xyz[1] * xyz[1]; value_t alp = 0; if (sectorAlpha || radPos2 < 1) { - alp = std::atan2(pxpypz[1], pxpypz[0]); + alp = gpu::CAMath::ATan2(pxpypz[1], pxpypz[0]); } else { - alp = std::atan2(xyz[1], xyz[0]); + alp = gpu::CAMath::ATan2(xyz[1], xyz[0]); } if (sectorAlpha) { alp = math_utils::detail::angle2Alpha(alp); @@ -70,14 +68,14 @@ TrackParametrization::TrackParametrization(const dim3_t& xyz, const dim value_t sn, cs; math_utils::detail::sincos(alp, sn, cs); // protection: avoid alpha being too close to 0 or +-pi/2 - if (std::fabs(sn) < 2 * kSafe) { + if (gpu::CAMath::Abs(sn) < 2 * kSafe) { if (alp > 0) { alp += alp < constants::math::PIHalf ? 2 * kSafe : -2 * kSafe; } else { alp += alp > -constants::math::PIHalf ? -2 * kSafe : 2 * kSafe; } math_utils::detail::sincos(alp, sn, cs); - } else if (std::fabs(cs) < 2 * kSafe) { + } else if (gpu::CAMath::Abs(cs) < 2 * kSafe) { if (alp > 0) { alp += alp > constants::math::PIHalf ? 2 * kSafe : -2 * kSafe; } else { @@ -100,12 +98,12 @@ TrackParametrization::TrackParametrization(const dim3_t& xyz, const dim mP[kZ] = ver[2]; mP[kSnp] = mom[1] * ptI; mP[kTgl] = mom[2] * ptI; - mAbsCharge = std::abs(charge); + mAbsCharge = gpu::CAMath::Abs(charge); mP[kQ2Pt] = charge ? ptI * charge : ptI; // - if (std::fabs(1 - getSnp()) < kSafe) { + if (gpu::CAMath::Abs(1 - getSnp()) < kSafe) { mP[kSnp] = 1.f - kSafe; // Protection - } else if (std::fabs(-1 - getSnp()) < kSafe) { + } else if (gpu::CAMath::Abs(-1 - getSnp()) < kSafe) { mP[kSnp] = -1.f + kSafe; // Protection } // @@ -113,14 +111,14 @@ TrackParametrization::TrackParametrization(const dim3_t& xyz, const dim //_______________________________________________________ template -bool TrackParametrization::getPxPyPzGlo(dim3_t& pxyz) const +GPUd() bool TrackParametrization::getPxPyPzGlo(dim3_t& pxyz) const { // track momentum - if (std::fabs(getQ2Pt()) < constants::math::Almost0 || std::fabs(getSnp()) > constants::math::Almost1) { + if (gpu::CAMath::Abs(getQ2Pt()) < constants::math::Almost0 || gpu::CAMath::Abs(getSnp()) > constants::math::Almost1) { return false; } value_t cs, sn, pt = getPt(); - value_t r = std::sqrt((1.f - getSnp()) * (1.f + getSnp())); + value_t r = gpu::CAMath::Sqrt((1.f - getSnp()) * (1.f + getSnp())); math_utils::detail::sincos(getAlpha(), sn, cs); pxyz[0] = pt * (r * cs - getSnp() * sn); pxyz[1] = pt * (getSnp() * cs + r * sn); @@ -130,17 +128,17 @@ bool TrackParametrization::getPxPyPzGlo(dim3_t& pxyz) const //____________________________________________________ template -bool TrackParametrization::getPosDirGlo(gpu::gpustd::array& posdirp) const +GPUd() bool TrackParametrization::getPosDirGlo(gpu::gpustd::array& posdirp) const { // fill vector with lab x,y,z,px/p,py/p,pz/p,p,sinAlpha,cosAlpha - value_t ptI = std::fabs(getQ2Pt()); + value_t ptI = gpu::CAMath::Abs(getQ2Pt()); value_t snp = getSnp(); - if (ptI < constants::math::Almost0 || std::fabs(snp) > constants::math::Almost1) { + if (ptI < constants::math::Almost0 || gpu::CAMath::Abs(snp) > constants::math::Almost1) { return false; } value_t &sn = posdirp[7], &cs = posdirp[8]; - value_t csp = std::sqrt((1.f - snp) * (1.f + snp)); - value_t cstht = std::sqrt(1.f + getTgl() * getTgl()); + value_t csp = gpu::CAMath::Sqrt((1.f - snp) * (1.f + snp)); + value_t cstht = gpu::CAMath::Sqrt(1.f + getTgl() * getTgl()); value_t csthti = 1.f / cstht; math_utils::detail::sincos(getAlpha(), sn, cs); posdirp[0] = getX() * cs - getY() * sn; @@ -155,10 +153,10 @@ bool TrackParametrization::getPosDirGlo(gpu::gpustd::array& //______________________________________________________________ template -bool TrackParametrization::rotateParam(value_t alpha) +GPUd() bool TrackParametrization::rotateParam(value_t alpha) { // rotate to alpha frame - if (std::fabs(getSnp()) > constants::math::Almost1) { + if (gpu::CAMath::Abs(getSnp()) > constants::math::Almost1) { LOGP(WARNING, "Precondition is not satisfied: |sin(phi)|>1 ! {:f}", getSnp()); return false; } @@ -167,7 +165,7 @@ bool TrackParametrization::rotateParam(value_t alpha) // value_t ca = 0, sa = 0; math_utils::detail::sincos(alpha - getAlpha(), sa, ca); - value_t snp = getSnp(), csp = std::sqrt((1.f - snp) * (1.f + snp)); // Improve precision + value_t snp = getSnp(), csp = gpu::CAMath::Sqrt((1.f - snp) * (1.f + snp)); // Improve precision // RS: check if rotation does no invalidate track model (cos(local_phi)>=0, i.e. particle // direction in local frame is along the X axis if ((csp * ca + snp * sa) < 0) { @@ -176,7 +174,7 @@ bool TrackParametrization::rotateParam(value_t alpha) } // value_t tmp = snp * ca - csp * sa; - if (std::fabs(tmp) > constants::math::Almost1) { + if (gpu::CAMath::Abs(tmp) > constants::math::Almost1) { LOGP(WARNING, "Rotation failed: new snp {:.2f}", tmp); return false; } @@ -190,7 +188,7 @@ bool TrackParametrization::rotateParam(value_t alpha) //____________________________________________________________ template -bool TrackParametrization::propagateParamTo(value_t xk, const dim3_t& b) +GPUd() bool TrackParametrization::propagateParamTo(value_t xk, const dim3_t& b) { //---------------------------------------------------------------- // Extrapolate this track params (w/o cov matrix) to the plane X=xk in the field b[]. @@ -200,11 +198,11 @@ bool TrackParametrization::propagateParamTo(value_t xk, const dim3_t& b //---------------------------------------------------------------- value_t dx = xk - getX(); - if (std::fabs(dx) < constants::math::Almost0) { + if (gpu::CAMath::Abs(dx) < constants::math::Almost0) { return true; } // Do not propagate tracks outside the ALICE detector - if (std::fabs(dx) > 1e5 || std::fabs(getY()) > 1e5 || std::fabs(getZ()) > 1e5) { + if (gpu::CAMath::Abs(dx) > 1e5 || gpu::CAMath::Abs(getY()) > 1e5 || gpu::CAMath::Abs(getZ()) > 1e5) { LOGP(WARNING, "Anomalous track, target X:{:f}", xk); // print(); return false; @@ -212,21 +210,21 @@ bool TrackParametrization::propagateParamTo(value_t xk, const dim3_t& b value_t crv = getCurvature(b[2]); value_t x2r = crv * dx; value_t f1 = getSnp(), f2 = f1 + x2r; - if (std::fabs(f1) > constants::math::Almost1 || std::fabs(f2) > constants::math::Almost1) { + if (gpu::CAMath::Abs(f1) > constants::math::Almost1 || gpu::CAMath::Abs(f2) > constants::math::Almost1) { return false; } - value_t r1 = std::sqrt((1.f - f1) * (1.f + f1)); - if (std::fabs(r1) < constants::math::Almost0) { + value_t r1 = gpu::CAMath::Sqrt((1.f - f1) * (1.f + f1)); + if (gpu::CAMath::Abs(r1) < constants::math::Almost0) { return false; } - value_t r2 = std::sqrt((1.f - f2) * (1.f + f2)); - if (std::fabs(r2) < constants::math::Almost0) { + value_t r2 = gpu::CAMath::Sqrt((1.f - f2) * (1.f + f2)); + if (gpu::CAMath::Abs(r2) < constants::math::Almost0) { return false; } value_t dy2dx = (f1 + f2) / (r1 + r2); - value_t step = (std::fabs(x2r) < 0.05f) ? dx * std::fabs(r2 + f2 * dy2dx) // chord - : 2.f * asinf(0.5f * dx * std::sqrt(1.f + dy2dx * dy2dx) * crv) / crv; // arc - step *= std::sqrt(1.f + getTgl() * getTgl()); + value_t step = (gpu::CAMath::Abs(x2r) < 0.05f) ? dx * gpu::CAMath::Abs(r2 + f2 * dy2dx) // chord + : 2.f * CAMath::ASin(0.5f * dx * gpu::CAMath::Sqrt(1.f + dy2dx * dy2dx) * crv) / crv; // arc + step *= gpu::CAMath::Sqrt(1.f + getTgl() * getTgl()); // // get the track x,y,z,px/p,py/p,pz/p,p,sinAlpha,cosAlpha in the Global System gpu::gpustd::array vecLab{0.f}; @@ -236,13 +234,13 @@ bool TrackParametrization::propagateParamTo(value_t xk, const dim3_t& b // rotate to the system where Bx=By=0. value_t bxy2 = b[0] * b[0] + b[1] * b[1]; - value_t bt = std::sqrt(bxy2); + value_t bt = gpu::CAMath::Sqrt(bxy2); value_t cosphi = 1.f, sinphi = 0.f; if (bt > constants::math::Almost0) { cosphi = b[0] / bt; sinphi = b[1] / bt; } - value_t bb = std::sqrt(bxy2 + b[2] * b[2]); + value_t bb = gpu::CAMath::Sqrt(bxy2 + b[2] * b[2]); value_t costet = 1.f, sintet = 0.f; if (bb > constants::math::Almost0) { costet = b[2] / bb; @@ -280,8 +278,8 @@ bool TrackParametrization::propagateParamTo(value_t xk, const dim3_t& b // Do the final correcting step to the target plane (linear approximation) value_t x = vecLab[0], y = vecLab[1], z = vecLab[2]; - if (std::fabs(dx) > constants::math::Almost0) { - if (std::fabs(vecLab[3]) < constants::math::Almost0) { + if (gpu::CAMath::Abs(dx) > constants::math::Almost0) { + if (gpu::CAMath::Abs(vecLab[3]) < constants::math::Almost0) { return false; } dx = xk - vecLab[0]; @@ -291,7 +289,7 @@ bool TrackParametrization::propagateParamTo(value_t xk, const dim3_t& b } // Calculate the track parameters - t = 1.f / std::sqrt(vecLab[3] * vecLab[3] + vecLab[4] * vecLab[4]); + t = 1.f / gpu::CAMath::Sqrt(vecLab[3] * vecLab[3] + vecLab[4] * vecLab[4]); mX = x; mP[kY] = y; mP[kZ] = z; @@ -304,7 +302,7 @@ bool TrackParametrization::propagateParamTo(value_t xk, const dim3_t& b //____________________________________________________________ template -bool TrackParametrization::propagateParamTo(value_t xk, value_t b) +GPUd() bool TrackParametrization::propagateParamTo(value_t xk, value_t b) { //---------------------------------------------------------------- // propagate this track to the plane X=xk (cm) in the field "b" (kG) @@ -312,28 +310,28 @@ bool TrackParametrization::propagateParamTo(value_t xk, value_t b) // distances only ( constants::math::Almost1) || (std::fabs(f2) > constants::math::Almost1)) { + if ((gpu::CAMath::Abs(f1) > constants::math::Almost1) || (gpu::CAMath::Abs(f2) > constants::math::Almost1)) { return false; } - value_t r1 = std::sqrt((1.f - f1) * (1.f + f1)); - if (std::fabs(r1) < constants::math::Almost0) { + value_t r1 = gpu::CAMath::Sqrt((1.f - f1) * (1.f + f1)); + if (gpu::CAMath::Abs(r1) < constants::math::Almost0) { return false; } - value_t r2 = std::sqrt((1.f - f2) * (1.f + f2)); - if (std::fabs(r2) < constants::math::Almost0) { + value_t r2 = gpu::CAMath::Sqrt((1.f - f2) * (1.f + f2)); + if (gpu::CAMath::Abs(r2) < constants::math::Almost0) { return false; } mX = xk; double dy2dx = (f1 + f2) / (r1 + r2); mP[kY] += dx * dy2dx; mP[kSnp] += x2r; - if (std::fabs(x2r) < 0.05f) { + if (gpu::CAMath::Abs(x2r) < 0.05f) { mP[kZ] += dx * (r2 + f2 * dy2dx) * getTgl(); } else { // for small dx/R the linear apporximation of the arc by the segment is OK, @@ -344,7 +342,7 @@ bool TrackParametrization::propagateParamTo(value_t xk, value_t b) // double rot = 2*TMath::ASin(0.5*chord*crv); // angular difference seen from the circle center // track1 += rot/crv*track3; // - value_t rot = asinf(r1 * f2 - r2 * f1); // more economic version from Yura. + value_t rot = CAMath::ASin(r1 * f2 - r2 * f1); // more economic version from Yura. if (f1 * f1 + f2 * f2 > 1.f && f1 * f2 < 0.f) { // special cases of large rotations or large abs angles if (f2 > 0.f) { rot = constants::math::PI - rot; // @@ -359,32 +357,32 @@ bool TrackParametrization::propagateParamTo(value_t xk, value_t b) //_______________________________________________________________________ template -bool TrackParametrization::propagateParamToDCA(const math_utils::Point3D& vtx, value_t b, dim2_t* dca, value_t maxD) +GPUd() bool TrackParametrization::propagateParamToDCA(const math_utils::Point3D& vtx, value_t b, dim2_t* dca, value_t maxD) { // propagate track to DCA to the vertex value_t sn, cs, alp = getAlpha(); math_utils::detail::sincos(alp, sn, cs); - value_t x = getX(), y = getY(), snp = getSnp(), csp = std::sqrt((1.f - snp) * (1.f + snp)); + value_t x = getX(), y = getY(), snp = getSnp(), csp = gpu::CAMath::Sqrt((1.f - snp) * (1.f + snp)); value_t xv = vtx.X() * cs + vtx.Y() * sn, yv = -vtx.X() * sn + vtx.Y() * cs, zv = vtx.Z(); x -= xv; y -= yv; //Estimate the impact parameter neglecting the track curvature - value_t d = std::abs(x * snp - y * csp); + value_t d = gpu::CAMath::Abs(x * snp - y * csp); if (d > maxD) { return false; } value_t crv = getCurvature(b); value_t tgfv = -(crv * x - snp) / (crv * y + csp); - sn = tgfv / std::sqrt(1.f + tgfv * tgfv); - cs = std::sqrt((1.f - sn) * (1.f + sn)); - cs = (std::abs(tgfv) > constants::math::Almost0) ? sn / tgfv : constants::math::Almost1; + sn = tgfv / gpu::CAMath::Sqrt(1.f + tgfv * tgfv); + cs = gpu::CAMath::Sqrt((1.f - sn) * (1.f + sn)); + cs = (gpu::CAMath::Abs(tgfv) > constants::math::Almost0) ? sn / tgfv : constants::math::Almost1; x = xv * cs + yv * sn; yv = -xv * sn + yv * cs; xv = x; auto tmpT(*this); // operate on the copy to recover after the failure - alp += std::asin(sn); + alp += gpu::CAMath::ASin(sn); if (!tmpT.rotateParam(alp) || !tmpT.propagateParamTo(xv, b)) { LOG(WARNING) << "failed to propagate to alpha=" << alp << " X=" << xv << " for vertex " << vtx.X() << ' ' << vtx.Y() << ' ' << vtx.Z() << " | Track is: "; @@ -401,33 +399,33 @@ bool TrackParametrization::propagateParamToDCA(const math_utils::Point3 //____________________________________________________________ template -bool TrackParametrization::getYZAt(value_t xk, value_t b, value_t& y, value_t& z) const +GPUd() bool TrackParametrization::getYZAt(value_t xk, value_t b, value_t& y, value_t& z) const { //---------------------------------------------------------------- // estimate Y,Z in tracking frame at given X //---------------------------------------------------------------- value_t dx = xk - getX(); - if (std::fabs(dx) < constants::math::Almost0) { + if (gpu::CAMath::Abs(dx) < constants::math::Almost0) { return true; } value_t crv = getCurvature(b); value_t x2r = crv * dx; value_t f1 = getSnp(), f2 = f1 + x2r; - if ((std::fabs(f1) > constants::math::Almost1) || (std::fabs(f2) > constants::math::Almost1)) { + if ((gpu::CAMath::Abs(f1) > constants::math::Almost1) || (gpu::CAMath::Abs(f2) > constants::math::Almost1)) { return false; } - value_t r1 = std::sqrt((1.f - f1) * (1.f + f1)); - if (std::fabs(r1) < constants::math::Almost0) { + value_t r1 = gpu::CAMath::Sqrt((1.f - f1) * (1.f + f1)); + if (gpu::CAMath::Abs(r1) < constants::math::Almost0) { return false; } - value_t r2 = std::sqrt((1.f - f2) * (1.f + f2)); - if (std::fabs(r2) < constants::math::Almost0) { + value_t r2 = gpu::CAMath::Sqrt((1.f - f2) * (1.f + f2)); + if (gpu::CAMath::Abs(r2) < constants::math::Almost0) { return false; } double dy2dx = (f1 + f2) / (r1 + r2); y = mP[kY] + dx * dy2dx; z = mP[kZ]; - if (std::fabs(x2r) < 0.05f) { + if (gpu::CAMath::Abs(x2r) < 0.05f) { z += dx * (r2 + f2 * dy2dx) * getTgl(); } else { // for small dx/R the linear apporximation of the arc by the segment is OK, @@ -438,7 +436,7 @@ bool TrackParametrization::getYZAt(value_t xk, value_t b, value_t& y, v // double rot = 2*TMath::ASin(0.5*chord*crv); // angular difference seen from the circle center // track1 += rot/crv*track3; // - value_t rot = asinf(r1 * f2 - r2 * f1); // more economic version from Yura. + value_t rot = CAMath::ASin(r1 * f2 - r2 * f1); // more economic version from Yura. if (f1 * f1 + f2 * f2 > 1.f && f1 * f2 < 0.f) { // special cases of large rotations or large abs angles if (f2 > 0.f) { rot = constants::math::PI - rot; // @@ -453,7 +451,7 @@ bool TrackParametrization::getYZAt(value_t xk, value_t b, value_t& y, v //______________________________________________________________ template -void TrackParametrization::invertParam() +GPUd() void TrackParametrization::invertParam() { // Transform this track to the local coord. system rotated by 180 deg. mX = -mX; @@ -468,7 +466,7 @@ void TrackParametrization::invertParam() //______________________________________________________________ template -typename TrackParametrization::value_t TrackParametrization::getZAt(value_t xk, value_t b) const +GPUd() typename TrackParametrization::value_t TrackParametrization::getZAt(value_t xk, value_t b) const { ///< this method is just an alias for obtaining Z @ X in the tree->Draw() value_t y, z; @@ -477,7 +475,7 @@ typename TrackParametrization::value_t TrackParametrization::g //______________________________________________________________ template -typename TrackParametrization::value_t TrackParametrization::getYAt(value_t xk, value_t b) const +GPUd() typename TrackParametrization::value_t TrackParametrization::getYAt(value_t xk, value_t b) const { ///< this method is just an alias for obtaining Z @ X in the tree->Draw() value_t y, z; @@ -497,7 +495,7 @@ std::string TrackParametrization::asString() const //______________________________________________________________ template -void TrackParametrization::printParam() const +GPUd() void TrackParametrization::printParam() const { // print parameters #ifndef GPUCA_ALIGPUCODE @@ -507,7 +505,7 @@ void TrackParametrization::printParam() const //______________________________________________________________ template -bool TrackParametrization::getXatLabR(value_t r, value_t& x, value_t bz, track::DirType dir) const +GPUd() bool TrackParametrization::getXatLabR(value_t r, value_t& x, value_t bz, track::DirType dir) const { // Get local X of the track position estimated at the radius lab radius r. // The track curvature is accounted exactly @@ -521,16 +519,16 @@ bool TrackParametrization::getXatLabR(value_t r, value_t& x, value_t bz const value_t kEps = 1.e-6; // auto crv = getCurvature(bz); - if (std::fabs(crv) > constants::math::Almost0) { // helix + if (gpu::CAMath::Abs(crv) > constants::math::Almost0) { // helix // get center of the track circle math_utils::CircleXY circle; getCircleParamsLoc(bz, circle); - value_t r0 = std::sqrt(circle.getCenterD2()); + value_t r0 = gpu::CAMath::Sqrt(circle.getCenterD2()); if (r0 <= constants::math::Almost0) { return false; // the track is concentric to circle } value_t tR2r0 = 1.f, g = 0.f, tmp = 0.f; - if (std::fabs(circle.rC - r0) > kEps) { + if (gpu::CAMath::Abs(circle.rC - r0) > kEps) { tR2r0 = circle.rC / r0; g = 0.5f * (r * r / (r0 * circle.rC) - tR2r0 - 1.f / tR2r0); tmp = 1.f + g * tR2r0; @@ -543,7 +541,7 @@ bool TrackParametrization::getXatLabR(value_t r, value_t& x, value_t bz if (det < 0.f) { return false; // does not reach raduis r } - det = std::sqrt(det); + det = gpu::CAMath::Sqrt(det); // // the intersection happens in 2 points: {circle.xC+tR*C,circle.yC+tR*S} // with C=f*c0+-|s0|*det and S=f*s0-+c0 sign(s0)*det @@ -551,8 +549,8 @@ bool TrackParametrization::getXatLabR(value_t r, value_t& x, value_t bz // x = circle.xC * tmp; value_t y = circle.yC * tmp; - if (std::fabs(circle.yC) > constants::math::Almost0) { // when circle.yC==0 the x,y is unique - value_t dfx = tR2r0 * std::fabs(circle.yC) * det; + if (gpu::CAMath::Abs(circle.yC) > constants::math::Almost0) { // when circle.yC==0 the x,y is unique + value_t dfx = tR2r0 * gpu::CAMath::Abs(circle.yC) * det; value_t dfy = tR2r0 * circle.xC * (circle.yC > 0.f ? det : -det); if (dir == DirAuto) { // chose the one which corresponds to smallest step value_t delta = (x - mX) * dfx - (y - fy) * dfy; // the choice of + in C will lead to smaller step if delta<0 @@ -563,7 +561,7 @@ bool TrackParametrization::getXatLabR(value_t r, value_t& x, value_t bz if (dfeps < -kEps) { return true; } - if (std::fabs(dfeps) < kEps && std::fabs(mX * mX + fy * fy - r * r) < kEps) { // are we already in right r? + if (gpu::CAMath::Abs(dfeps) < kEps && gpu::CAMath::Abs(mX * mX + fy * fy - r * r) < kEps) { // are we already in right r? return mX; } x += dfx + dfx; @@ -580,7 +578,7 @@ bool TrackParametrization::getXatLabR(value_t r, value_t& x, value_t bz if (dfeps < -kEps) { return true; } - if (std::fabs(dfeps) < kEps && std::fabs(mX * mX + fy * fy - r * r) < kEps) { // are we already in right r? + if (gpu::CAMath::Abs(dfeps) < kEps && gpu::CAMath::Abs(mX * mX + fy * fy - r * r) < kEps) { // are we already in right r? return mX; } x -= dfx + dfx; @@ -598,8 +596,8 @@ bool TrackParametrization::getXatLabR(value_t r, value_t& x, value_t bz return false; } } - } else { // this is a straight track - if (std::fabs(sn) >= constants::math::Almost1) { // || to Y axis + } else { // this is a straight track + if (gpu::CAMath::Abs(sn) >= constants::math::Almost1) { // || to Y axis value_t det = (r - mX) * (r + mX); if (det < 0.f) { return false; // does not reach raduis r @@ -608,7 +606,7 @@ bool TrackParametrization::getXatLabR(value_t r, value_t& x, value_t bz if (dir == DirAuto) { return true; } - det = std::sqrt(det); + det = gpu::CAMath::Sqrt(det); if (dir == DirOutward) { // along the track direction if (sn > 0.f) { if (fy > det) { @@ -628,12 +626,12 @@ bool TrackParametrization::getXatLabR(value_t r, value_t& x, value_t bz return false; // track is against Y axis } } - } else if (std::fabs(sn) <= constants::math::Almost0) { // || to X axis + } else if (gpu::CAMath::Abs(sn) <= constants::math::Almost0) { // || to X axis value_t det = (r - fy) * (r + fy); if (det < 0.f) { return false; // does not reach raduis r } - det = std::sqrt(det); + det = gpu::CAMath::Sqrt(det); if (dir == DirAuto) { x = mX > 0.f ? det : -det; // choose the solution requiring the smalest step return true; @@ -651,13 +649,13 @@ bool TrackParametrization::getXatLabR(value_t r, value_t& x, value_t bz } } } else { // general case of straight line - value_t cs = std::sqrt((1.f - sn) * (1.f + sn)); + value_t cs = gpu::CAMath::Sqrt((1.f - sn) * (1.f + sn)); value_t xsyc = mX * sn - fy * cs; value_t det = (r - xsyc) * (r + xsyc); if (det < 0.f) { return false; // does not reach raduis r } - det = std::sqrt(det); + det = gpu::CAMath::Sqrt(det); value_t xcys = mX * cs + fy * sn; value_t t = -xcys; if (dir == DirAuto) { @@ -684,7 +682,7 @@ bool TrackParametrization::getXatLabR(value_t r, value_t& x, value_t bz //______________________________________________ template -bool TrackParametrization::correctForELoss(value_t xrho, value_t mass, bool anglecorr, value_t dedx) +GPUd() bool TrackParametrization::correctForELoss(value_t xrho, value_t mass, bool anglecorr, value_t dedx) { //------------------------------------------------------------------ // This function corrects the track parameters for the energy loss in crossed material. @@ -702,7 +700,7 @@ bool TrackParametrization::correctForELoss(value_t xrho, value_t mass, if (anglecorr) { value_t csp2 = (1.f - getSnp()) * (1.f + getSnp()); // cos(phi)^2 value_t cst2I = (1.f + getTgl() * getTgl()); // 1/cos(lambda)^2 - value_t angle = std::sqrt(cst2I / (csp2)); + value_t angle = gpu::CAMath::Sqrt(cst2I / (csp2)); xrho *= angle; } @@ -717,15 +715,15 @@ bool TrackParametrization::correctForELoss(value_t xrho, value_t mass, // Calculating the energy loss corrections************************ if ((xrho != 0.f) && (beta2 < 1.f)) { if (dedx < kCalcdEdxAuto + constants::math::Almost1) { // request to calculate dedx on the fly - dedx = BetheBlochSolid(p / std::fabs(mass)); + dedx = BetheBlochSolid(p / gpu::CAMath::Abs(mass)); if (mAbsCharge != 1) { dedx *= mAbsCharge * mAbsCharge; } } value_t dE = dedx * xrho; - value_t e = std::sqrt(e2); - if (std::fabs(dE) > kMaxELossFrac * e) { + value_t e = gpu::CAMath::Sqrt(e2); + if (gpu::CAMath::Abs(dE) > kMaxELossFrac * e) { return false; // 30% energy loss is too much! } value_t eupd = e + dE; @@ -733,14 +731,16 @@ bool TrackParametrization::correctForELoss(value_t xrho, value_t mass, if (pupd2 < kMinP * kMinP) { return false; } - setQ2Pt(getQ2Pt() * p / std::sqrt(pupd2)); + setQ2Pt(getQ2Pt() * p / gpu::CAMath::Sqrt(pupd2)); } return true; } +namespace o2::track +{ template class TrackParametrization; +#ifndef GPUCA_GPUCODE_DEVICE template class TrackParametrization; - -} // namespace track -} // namespace o2 +#endif +} // namespace o2::track diff --git a/DataFormats/Reconstruction/src/TrackParametrizationWithError.cxx b/DataFormats/Reconstruction/src/TrackParametrizationWithError.cxx index 23e26b1298122..fd0e20c0f9f01 100644 --- a/DataFormats/Reconstruction/src/TrackParametrizationWithError.cxx +++ b/DataFormats/Reconstruction/src/TrackParametrizationWithError.cxx @@ -37,14 +37,12 @@ #include #endif -namespace o2 -{ -namespace track -{ +using namespace o2::track; +using namespace o2::gpu; //______________________________________________________________ template -void TrackParametrizationWithError::invert() +GPUd() void TrackParametrizationWithError::invert() { // Transform this track to the local coord. system rotated by 180 deg. this->invertParam(); @@ -59,27 +57,27 @@ void TrackParametrizationWithError::invert() //______________________________________________________________ template -bool TrackParametrizationWithError::propagateTo(value_t xk, value_t b) +GPUd() bool TrackParametrizationWithError::propagateTo(value_t xk, value_t b) { //---------------------------------------------------------------- // propagate this track to the plane X=xk (cm) in the field "b" (kG) //---------------------------------------------------------------- value_t dx = xk - this->getX(); - if (std::fabs(dx) < constants::math::Almost0) { + if (gpu::CAMath::Abs(dx) < constants::math::Almost0) { return true; } value_t crv = this->getCurvature(b); value_t x2r = crv * dx; value_t f1 = this->getSnp(), f2 = f1 + x2r; - if ((std::fabs(f1) > constants::math::Almost1) || (std::fabs(f2) > constants::math::Almost1)) { + if ((gpu::CAMath::Abs(f1) > constants::math::Almost1) || (gpu::CAMath::Abs(f2) > constants::math::Almost1)) { return false; } - value_t r1 = std::sqrt((1.f - f1) * (1.f + f1)); - if (std::fabs(r1) < constants::math::Almost0) { + value_t r1 = gpu::CAMath::Sqrt((1.f - f1) * (1.f + f1)); + if (gpu::CAMath::Abs(r1) < constants::math::Almost0) { return false; } - value_t r2 = std::sqrt((1.f - f2) * (1.f + f2)); - if (std::fabs(r2) < constants::math::Almost0) { + value_t r2 = gpu::CAMath::Sqrt((1.f - f2) * (1.f + f2)); + if (gpu::CAMath::Abs(r2) < constants::math::Almost0) { return false; } this->setX(xk); @@ -87,7 +85,7 @@ bool TrackParametrizationWithError::propagateTo(value_t xk, value_t b) value_t dP[kNParams] = {0.f}; dP[kY] = dx * dy2dx; dP[kSnp] = x2r; - if (std::fabs(x2r) < 0.05f) { + if (gpu::CAMath::Abs(x2r) < 0.05f) { dP[kZ] = dx * (r2 + f2 * dy2dx) * this->getTgl(); } else { // for small dx/R the linear apporximation of the arc by the segment is OK, @@ -98,7 +96,7 @@ bool TrackParametrizationWithError::propagateTo(value_t xk, value_t b) // double rot = 2*TMath::ASin(0.5*chord*crv); // angular difference seen from the circle center // mP1 += rot/crv*mP3; // - value_t rot = std::asin(r1 * f2 - r2 * f1); // more economic version from Yura. + value_t rot = gpu::CAMath::ASin(r1 * f2 - r2 * f1); // more economic version from Yura. if (f1 * f1 + f2 * f2 > 1.f && f1 * f2 < 0.f) { // special cases of large rotations or large abs angles if (f2 > 0.f) { rot = constants::math::PI - rot; // @@ -164,10 +162,10 @@ bool TrackParametrizationWithError::propagateTo(value_t xk, value_t b) //______________________________________________________________ template -bool TrackParametrizationWithError::rotate(value_t alpha) +GPUd() bool TrackParametrizationWithError::rotate(value_t alpha) { // rotate to alpha frame - if (std::fabs(this->getSnp()) > constants::math::Almost1) { + if (gpu::CAMath::Abs(this->getSnp()) > constants::math::Almost1) { LOGP(WARNING, "Precondition is not satisfied: |sin(phi)|>1 ! {:f}", this->getSnp()); return false; } @@ -176,7 +174,7 @@ bool TrackParametrizationWithError::rotate(value_t alpha) // value_t ca = 0, sa = 0; math_utils::detail::sincos(alpha - this->getAlpha(), sa, ca); - value_t snp = this->getSnp(), csp = std::sqrt((1.f - snp) * (1.f + snp)); // Improve precision + value_t snp = this->getSnp(), csp = gpu::CAMath::Sqrt((1.f - snp) * (1.f + snp)); // Improve precision // RS: check if rotation does no invalidate track model (cos(local_phi)>=0, i.e. particle // direction in local frame is along the X axis if ((csp * ca + snp * sa) < 0) { @@ -186,7 +184,7 @@ bool TrackParametrizationWithError::rotate(value_t alpha) // value_t updSnp = snp * ca - csp * sa; - if (std::fabs(updSnp) > constants::math::Almost1) { + if (gpu::CAMath::Abs(updSnp) > constants::math::Almost1) { LOGP(WARNING, "Rotation failed: new snp {:.2f}", updSnp); return false; } @@ -196,7 +194,7 @@ bool TrackParametrizationWithError::rotate(value_t alpha) this->setY(-xold * sa + yold * ca); this->setSnp(updSnp); - if (std::fabs(csp) < constants::math::Almost0) { + if (gpu::CAMath::Abs(csp) < constants::math::Almost0) { LOGP(WARNING, "Too small cosine value {:f}", csp); csp = constants::math::Almost0; } @@ -219,32 +217,32 @@ bool TrackParametrizationWithError::rotate(value_t alpha) //_______________________________________________________________________ template -bool TrackParametrizationWithError::propagateToDCA(const o2::dataformats::VertexBase& vtx, value_t b, o2::dataformats::DCA* dca, value_t maxD) +GPUd() bool TrackParametrizationWithError::propagateToDCA(const o2::dataformats::VertexBase& vtx, value_t b, o2::dataformats::DCA* dca, value_t maxD) { // propagate track to DCA to the vertex value_t sn, cs, alp = this->getAlpha(); o2::math_utils::detail::sincos(alp, sn, cs); - value_t x = this->getX(), y = this->getY(), snp = this->getSnp(), csp = std::sqrt((1.f - snp) * (1.f + snp)); + value_t x = this->getX(), y = this->getY(), snp = this->getSnp(), csp = gpu::CAMath::Sqrt((1.f - snp) * (1.f + snp)); value_t xv = vtx.getX() * cs + vtx.getY() * sn, yv = -vtx.getX() * sn + vtx.getY() * cs, zv = vtx.getZ(); x -= xv; y -= yv; //Estimate the impact parameter neglecting the track curvature - value_t d = std::abs(x * snp - y * csp); + value_t d = gpu::CAMath::Abs(x * snp - y * csp); if (d > maxD) { return false; } value_t crv = this->getCurvature(b); value_t tgfv = -(crv * x - snp) / (crv * y + csp); - sn = tgfv / std::sqrt(1.f + tgfv * tgfv); - cs = std::sqrt((1.f - sn) * (1.f + sn)); - cs = (std::abs(tgfv) > constants::math::Almost0) ? sn / tgfv : constants::math::Almost1; + sn = tgfv / gpu::CAMath::Sqrt(1.f + tgfv * tgfv); + cs = gpu::CAMath::Sqrt((1.f - sn) * (1.f + sn)); + cs = (gpu::CAMath::Abs(tgfv) > constants::math::Almost0) ? sn / tgfv : constants::math::Almost1; x = xv * cs + yv * sn; yv = -xv * sn + yv * cs; xv = x; auto tmpT(*this); // operate on the copy to recover after the failure - alp += std::asin(sn); + alp += gpu::CAMath::ASin(sn); if (!tmpT.rotate(alp) || !tmpT.propagateTo(xv, b)) { LOG(WARNING) << "failed to propagate to alpha=" << alp << " X=" << xv << vtx << " | Track is: "; tmpT.print(); @@ -261,7 +259,7 @@ bool TrackParametrizationWithError::propagateToDCA(const o2::dataformat //______________________________________________________________ template -TrackParametrizationWithError::TrackParametrizationWithError(const dim3_t& xyz, const dim3_t& pxpypz, +GPUd() TrackParametrizationWithError::TrackParametrizationWithError(const dim3_t& xyz, const dim3_t& pxpypz, const gpu::gpustd::array& cv, int charge, bool sectorAlpha) { // construct track param and covariance from kinematics and lab errors @@ -276,9 +274,9 @@ TrackParametrizationWithError::TrackParametrizationWithError(const dim3 value_t radPos2 = xyz[0] * xyz[0] + xyz[1] * xyz[1]; value_t alp = 0; if (sectorAlpha || radPos2 < 1) { - alp = std::atan2(pxpypz[1], pxpypz[0]); + alp = gpu::CAMath::ATan2(pxpypz[1], pxpypz[0]); } else { - alp = std::atan2(xyz[1], xyz[0]); + alp = gpu::CAMath::ATan2(xyz[1], xyz[0]); } if (sectorAlpha) { alp = math_utils::detail::angle2Alpha(alp); @@ -287,14 +285,14 @@ TrackParametrizationWithError::TrackParametrizationWithError(const dim3 value_t sn, cs; math_utils::detail::sincos(alp, sn, cs); // protection: avoid alpha being too close to 0 or +-pi/2 - if (std::fabs(sn) < 2.f * kSafe) { + if (gpu::CAMath::Abs(sn) < 2.f * kSafe) { if (alp > 0) { alp += alp < constants::math::PIHalf ? 2.f * kSafe : -2.f * kSafe; } else { alp += alp > -constants::math::PIHalf ? -2.f * kSafe : 2.f * kSafe; } math_utils::detail::sincos(alp, sn, cs); - } else if (std::fabs(cs) < 2.f * kSafe) { + } else if (gpu::CAMath::Abs(cs) < 2.f * kSafe) { if (alp > 0) { alp += alp > constants::math::PIHalf ? 2.f * kSafe : -2.f * kSafe; } else { @@ -310,7 +308,7 @@ TrackParametrizationWithError::TrackParametrizationWithError(const dim3 math_utils::detail::rotateZ(ver, -alp); math_utils::detail::rotateZ(mom, -alp); // - value_t pt = std::sqrt(mom[0] * mom[0] + mom[1] * mom[1]); + value_t pt = gpu::CAMath::Sqrt(mom[0] * mom[0] + mom[1] * mom[1]); value_t ptI = 1.f / pt; this->setX(ver[0]); this->setAlpha(alp); @@ -318,25 +316,25 @@ TrackParametrizationWithError::TrackParametrizationWithError(const dim3 this->setZ(ver[2]); this->setSnp(mom[1] * ptI); // cos(phi) this->setTgl(mom[2] * ptI); // tg(lambda) - this->setAbsCharge(std::abs(charge)); + this->setAbsCharge(gpu::CAMath::Abs(charge)); this->setQ2Pt(charge ? ptI * charge : ptI); // - if (std::fabs(1.f - this->getSnp()) < kSafe) { + if (gpu::CAMath::Abs(1.f - this->getSnp()) < kSafe) { this->setSnp(1.f - kSafe); // Protection - } else if (std::fabs(-1.f - this->getSnp()) < kSafe) { + } else if (gpu::CAMath::Abs(-1.f - this->getSnp()) < kSafe) { this->setSnp(-1.f + kSafe); // Protection } // // Covariance matrix (formulas to be simplified) value_t r = mom[0] * ptI; // cos(phi) - value_t cv34 = std::sqrt(cv[3] * cv[3] + cv[4] * cv[4]); + value_t cv34 = gpu::CAMath::Sqrt(cv[3] * cv[3] + cv[4] * cv[4]); // int special = 0; value_t sgcheck = r * sn + this->getSnp() * cs; - if (std::fabs(sgcheck) > 1 - kSafe) { // special case: lab phi is +-pi/2 + if (gpu::CAMath::Abs(sgcheck) > 1 - kSafe) { // special case: lab phi is +-pi/2 special = 1; sgcheck = sgcheck < 0 ? -1.f : 1.f; - } else if (std::fabs(sgcheck) < kSafe) { + } else if (gpu::CAMath::Abs(sgcheck) < kSafe) { sgcheck = cs < 0 ? -1.0f : 1.0f; special = 2; // special case: lab phi is 0 } @@ -350,29 +348,29 @@ TrackParametrizationWithError::TrackParametrizationWithError(const dim3 if (special == 1) { mC[kSigSnpY] = cv[6] * ptI; mC[kSigSnpZ] = -sgcheck * cv[8] * r * ptI; - mC[kSigSnp2] = std::fabs(cv[9] * r * r * ptI2); + mC[kSigSnp2] = gpu::CAMath::Abs(cv[9] * r * r * ptI2); mC[kSigTglY] = (cv[10] * this->getTgl() - sgcheck * cv[15]) * ptI / r; mC[kSigTglZ] = (cv[17] - sgcheck * cv[12] * this->getTgl()) * ptI; mC[kSigTglSnp] = (-sgcheck * cv[18] + cv[13] * this->getTgl()) * r * ptI2; - mC[kSigTgl2] = std::fabs(cv[20] - 2 * sgcheck * cv[19] * mC[4] + cv[14] * tgl2) * ptI2; + mC[kSigTgl2] = gpu::CAMath::Abs(cv[20] - 2 * sgcheck * cv[19] * mC[4] + cv[14] * tgl2) * ptI2; mC[kSigQ2PtY] = cv[10] * ptI2 / r * charge; mC[kSigQ2PtZ] = -sgcheck * cv[12] * ptI2 * charge; mC[kSigQ2PtSnp] = cv[13] * r * ptI * ptI2 * charge; mC[kSigQ2PtTgl] = (-sgcheck * cv[19] + cv[14] * this->getTgl()) * r * ptI2 * ptI; - mC[kSigQ2Pt2] = std::fabs(cv[14] * ptI2 * ptI2); + mC[kSigQ2Pt2] = gpu::CAMath::Abs(cv[14] * ptI2 * ptI2); } else if (special == 2) { mC[kSigSnpY] = -cv[10] * ptI * cs / sn; mC[kSigSnpZ] = cv[12] * cs * ptI; - mC[kSigSnp2] = std::fabs(cv[14] * cs * cs * ptI2); + mC[kSigSnp2] = gpu::CAMath::Abs(cv[14] * cs * cs * ptI2); mC[kSigTglY] = (sgcheck * cv[6] * this->getTgl() - cv[15]) * ptI / sn; mC[kSigTglZ] = (cv[17] - sgcheck * cv[8] * this->getTgl()) * ptI; mC[kSigTglSnp] = (cv[19] - sgcheck * cv[13] * this->getTgl()) * cs * ptI2; - mC[kSigTgl2] = std::fabs(cv[20] - 2 * sgcheck * cv[18] * this->getTgl() + cv[9] * tgl2) * ptI2; + mC[kSigTgl2] = gpu::CAMath::Abs(cv[20] - 2 * sgcheck * cv[18] * this->getTgl() + cv[9] * tgl2) * ptI2; mC[kSigQ2PtY] = sgcheck * cv[6] * ptI2 / sn * charge; mC[kSigQ2PtZ] = -sgcheck * cv[8] * ptI2 * charge; mC[kSigQ2PtSnp] = -sgcheck * cv[13] * cs * ptI * ptI2 * charge; mC[kSigQ2PtTgl] = (-sgcheck * cv[18] + cv[9] * this->getTgl()) * ptI2 * ptI * charge; - mC[kSigQ2Pt2] = std::fabs(cv[9] * ptI2 * ptI2); + mC[kSigQ2Pt2] = gpu::CAMath::Abs(cv[9] * ptI2 * ptI2); } else { double m00 = -sn; // m10=cs; double m23 = -pt * (sn + this->getSnp() * cs / r), m43 = -pt * pt * (r * cs - this->getSnp() * sn); @@ -396,8 +394,8 @@ TrackParametrizationWithError::TrackParametrizationWithError(const dim3 mC[kSigSnpZ] = (cv[12] * m43 - cv[8] * m44) / (m24 * m43 - m23 * m44); mC[kSigQ2PtZ] = (cv[8] - mC[kSigSnpZ] * m23) / m43; mC[kSigTglZ] = cv[17] / m35 - mC[kSigQ2PtZ] * m45 / m35; - mC[kSigSnp2] = std::fabs((a4 * a3 - a6 * a1) / (a5 * a3 - a6 * a2)); - mC[kSigQ2Pt2] = std::fabs((a1 - a2 * mC[kSigSnp2]) / a3); + mC[kSigSnp2] = gpu::CAMath::Abs((a4 * a3 - a6 * a1) / (a5 * a3 - a6 * a2)); + mC[kSigQ2Pt2] = gpu::CAMath::Abs((a1 - a2 * mC[kSigSnp2]) / a3); mC[kSigQ2PtSnp] = (cv[9] - mC[kSigSnp2] * m23 * m23 - mC[kSigQ2Pt2] * m43 * m43) / m23 / m43; double b1 = cv[18] - mC[kSigQ2PtSnp] * m23 * m45 - mC[kSigQ2Pt2] * m43 * m45; double b2 = m23 * m35; @@ -407,14 +405,14 @@ TrackParametrizationWithError::TrackParametrizationWithError(const dim3 double b6 = m44 * m35; mC[kSigTglSnp] = (b4 - b6 * b1 / b3) / (b5 - b6 * b2 / b3); mC[kSigQ2PtTgl] = b1 / b3 - b2 * mC[kSigTglSnp] / b3; - mC[kSigTgl2] = std::fabs((cv[20] - mC[kSigQ2Pt2] * (m45 * m45) - mC[kSigQ2PtTgl] * 2.f * m35 * m45) / (m35 * m35)); + mC[kSigTgl2] = gpu::CAMath::Abs((cv[20] - mC[kSigQ2Pt2] * (m45 * m45) - mC[kSigQ2PtTgl] * 2.f * m35 * m45) / (m35 * m35)); } checkCovariance(); } //____________________________________________________________ template -bool TrackParametrizationWithError::propagateTo(value_t xk, const dim3_t& b) +GPUd() bool TrackParametrizationWithError::propagateTo(value_t xk, const dim3_t& b) { //---------------------------------------------------------------- // Extrapolate this track to the plane X=xk in the field b[]. @@ -424,34 +422,34 @@ bool TrackParametrizationWithError::propagateTo(value_t xk, const dim3_ //---------------------------------------------------------------- value_t dx = xk - this->getX(); - if (std::fabs(dx) < constants::math::Almost0) { + if (gpu::CAMath::Abs(dx) < constants::math::Almost0) { return true; } // Do not propagate tracks outside the ALICE detector - if (std::fabs(dx) > 1e5 || std::fabs(this->getY()) > 1e5 || std::fabs(this->getZ()) > 1e5) { + if (gpu::CAMath::Abs(dx) > 1e5 || gpu::CAMath::Abs(this->getY()) > 1e5 || gpu::CAMath::Abs(this->getZ()) > 1e5) { LOGP(WARNING, "Anomalous track, target X:{:f}", xk); // print(); return false; } - value_t crv = (std::fabs(b[2]) < constants::math::Almost0) ? 0.f : this->getCurvature(b[2]); + value_t crv = (gpu::CAMath::Abs(b[2]) < constants::math::Almost0) ? 0.f : this->getCurvature(b[2]); value_t x2r = crv * dx; value_t f1 = this->getSnp(), f2 = f1 + x2r; - if ((std::fabs(f1) > constants::math::Almost1) || (std::fabs(f2) > constants::math::Almost1)) { + if ((gpu::CAMath::Abs(f1) > constants::math::Almost1) || (gpu::CAMath::Abs(f2) > constants::math::Almost1)) { return false; } - value_t r1 = std::sqrt((1.f - f1) * (1.f + f1)); - if (std::fabs(r1) < constants::math::Almost0) { + value_t r1 = gpu::CAMath::Sqrt((1.f - f1) * (1.f + f1)); + if (gpu::CAMath::Abs(r1) < constants::math::Almost0) { return false; } - value_t r2 = std::sqrt((1.f - f2) * (1.f + f2)); - if (std::fabs(r2) < constants::math::Almost0) { + value_t r2 = gpu::CAMath::Sqrt((1.f - f2) * (1.f + f2)); + if (gpu::CAMath::Abs(r2) < constants::math::Almost0) { return false; } value_t dy2dx = (f1 + f2) / (r1 + r2); - value_t step = (std::fabs(x2r) < 0.05f) ? dx * std::fabs(r2 + f2 * dy2dx) // chord - : 2.f * std::asin(0.5f * dx * std::sqrt(1.f + dy2dx * dy2dx) * crv) / crv; // arc - step *= std::sqrt(1.f + this->getTgl() * this->getTgl()); + value_t step = (gpu::CAMath::Abs(x2r) < 0.05f) ? dx * gpu::CAMath::Abs(r2 + f2 * dy2dx) // chord + : 2.f * gpu::CAMath::ASin(0.5f * dx * gpu::CAMath::Sqrt(1.f + dy2dx * dy2dx) * crv) / crv; // arc + step *= gpu::CAMath::Sqrt(1.f + this->getTgl() * this->getTgl()); // // get the track x,y,z,px/p,py/p,pz/p,p,sinAlpha,cosAlpha in the Global System gpu::gpustd::array vecLab{0.f}; @@ -509,13 +507,13 @@ bool TrackParametrizationWithError::propagateTo(value_t xk, const dim3_ // Rotate to the system where Bx=By=0. value_t bxy2 = b[0] * b[0] + b[1] * b[1]; - value_t bt = std::sqrt(bxy2); + value_t bt = gpu::CAMath::Sqrt(bxy2); value_t cosphi = 1.f, sinphi = 0.f; if (bt > constants::math::Almost0) { cosphi = b[0] / bt; sinphi = b[1] / bt; } - value_t bb = std::sqrt(bxy2 + b[2] * b[2]); + value_t bb = gpu::CAMath::Sqrt(bxy2 + b[2] * b[2]); value_t costet = 1., sintet = 0.; if (bb > constants::math::Almost0) { costet = b[2] / bb; @@ -553,8 +551,8 @@ bool TrackParametrizationWithError::propagateTo(value_t xk, const dim3_ // Do the final correcting step to the target plane (linear approximation) value_t x = vecLab[0], y = vecLab[1], z = vecLab[2]; - if (std::fabs(dx) > constants::math::Almost0) { - if (std::fabs(vecLab[3]) < constants::math::Almost0) { + if (gpu::CAMath::Abs(dx) > constants::math::Almost0) { + if (gpu::CAMath::Abs(vecLab[3]) < constants::math::Almost0) { return false; } dx = xk - vecLab[0]; @@ -564,7 +562,7 @@ bool TrackParametrizationWithError::propagateTo(value_t xk, const dim3_ } // Calculate the track parameters - t = 1.f / std::sqrt(vecLab[3] * vecLab[3] + vecLab[4] * vecLab[4]); + t = 1.f / gpu::CAMath::Sqrt(vecLab[3] * vecLab[3] + vecLab[4] * vecLab[4]); this->setX(x); this->setY(y); this->setZ(z); @@ -577,51 +575,51 @@ bool TrackParametrizationWithError::propagateTo(value_t xk, const dim3_ //______________________________________________ template -void TrackParametrizationWithError::checkCovariance() +GPUd() void TrackParametrizationWithError::checkCovariance() { // This function forces the diagonal elements of the covariance matrix to be positive. // In case the diagonal element is bigger than the maximal allowed value, it is set to // the limit and the off-diagonal elements that correspond to it are set to zero. - mC[kSigY2] = std::fabs(mC[kSigY2]); + mC[kSigY2] = gpu::CAMath::Abs(mC[kSigY2]); if (mC[kSigY2] > kCY2max) { - value_t scl = std::sqrt(kCY2max / mC[kSigY2]); + value_t scl = gpu::CAMath::Sqrt(kCY2max / mC[kSigY2]); mC[kSigY2] = kCY2max; mC[kSigZY] *= scl; mC[kSigSnpY] *= scl; mC[kSigTglY] *= scl; mC[kSigQ2PtY] *= scl; } - mC[kSigZ2] = std::fabs(mC[kSigZ2]); + mC[kSigZ2] = gpu::CAMath::Abs(mC[kSigZ2]); if (mC[kSigZ2] > kCZ2max) { - value_t scl = std::sqrt(kCZ2max / mC[kSigZ2]); + value_t scl = gpu::CAMath::Sqrt(kCZ2max / mC[kSigZ2]); mC[kSigZ2] = kCZ2max; mC[kSigZY] *= scl; mC[kSigSnpZ] *= scl; mC[kSigTglZ] *= scl; mC[kSigQ2PtZ] *= scl; } - mC[kSigSnp2] = std::fabs(mC[kSigSnp2]); + mC[kSigSnp2] = gpu::CAMath::Abs(mC[kSigSnp2]); if (mC[kSigSnp2] > kCSnp2max) { - value_t scl = std::sqrt(kCSnp2max / mC[kSigSnp2]); + value_t scl = gpu::CAMath::Sqrt(kCSnp2max / mC[kSigSnp2]); mC[kSigSnp2] = kCSnp2max; mC[kSigSnpY] *= scl; mC[kSigSnpZ] *= scl; mC[kSigTglSnp] *= scl; mC[kSigQ2PtSnp] *= scl; } - mC[kSigTgl2] = std::fabs(mC[kSigTgl2]); + mC[kSigTgl2] = gpu::CAMath::Abs(mC[kSigTgl2]); if (mC[kSigTgl2] > kCTgl2max) { - value_t scl = std::sqrt(kCTgl2max / mC[kSigTgl2]); + value_t scl = gpu::CAMath::Sqrt(kCTgl2max / mC[kSigTgl2]); mC[kSigTgl2] = kCTgl2max; mC[kSigTglY] *= scl; mC[kSigTglZ] *= scl; mC[kSigTglSnp] *= scl; mC[kSigQ2PtTgl] *= scl; } - mC[kSigQ2Pt2] = std::fabs(mC[kSigQ2Pt2]); + mC[kSigQ2Pt2] = gpu::CAMath::Abs(mC[kSigQ2Pt2]); if (mC[kSigQ2Pt2] > kC1Pt2max) { - value_t scl = std::sqrt(kC1Pt2max / mC[kSigQ2Pt2]); + value_t scl = gpu::CAMath::Sqrt(kC1Pt2max / mC[kSigQ2Pt2]); mC[kSigQ2Pt2] = kC1Pt2max; mC[kSigQ2PtY] *= scl; mC[kSigQ2PtZ] *= scl; @@ -632,7 +630,7 @@ void TrackParametrizationWithError::checkCovariance() //______________________________________________ template -void TrackParametrizationWithError::resetCovariance(value_t s2) +GPUd() void TrackParametrizationWithError::resetCovariance(value_t s2) { // Reset the covarince matrix to "something big" double d0(kCY2max), d1(kCZ2max), d2(kCSnp2max), d3(kCTgl2max), d4(kC1Pt2max); @@ -658,7 +656,9 @@ void TrackParametrizationWithError::resetCovariance(value_t s2) d4 = kC1Pt2max; } } - memset(mC, 0, kCovMatSize * sizeof(value_t)); + for (int i = 0; i < kCovMatSize; i++) { + mC[i] = 0; + } mC[kSigY2] = d0; mC[kSigZ2] = d1; mC[kSigSnp2] = d2; @@ -668,7 +668,7 @@ void TrackParametrizationWithError::resetCovariance(value_t s2) //______________________________________________ template -typename TrackParametrizationWithError::value_t TrackParametrizationWithError::getPredictedChi2(const dim2_t& p, const dim3_t& cov) const +GPUd() typename TrackParametrizationWithError::value_t TrackParametrizationWithError::getPredictedChi2(const dim2_t& p, const dim3_t& cov) const { // Estimate the chi2 of the space point "p" with the cov. matrix "cov" auto sdd = static_cast(getSigmaY2()) + static_cast(cov[0]); @@ -676,7 +676,7 @@ typename TrackParametrizationWithError::value_t TrackParametrizationWit auto szz = static_cast(getSigmaZ2()) + static_cast(cov[2]); auto det = sdd * szz - sdz * sdz; - if (std::fabs(det) < constants::math::Almost0) { + if (gpu::CAMath::Abs(det) < constants::math::Almost0) { return constants::math::VeryBig; } @@ -686,7 +686,7 @@ typename TrackParametrizationWithError::value_t TrackParametrizationWit return (d * (szz * d - sdz * z) + z * (sdd * z - d * sdz)) / det; } -#ifndef GPUCA_GPUCODE_DEVICE // Disable function relying on ROOT SMatrix on GPU +#ifndef GPUCA_GPUCODE // Disable function relying on ROOT SMatrix on GPU //______________________________________________ template @@ -725,11 +725,11 @@ typename TrackParametrizationWithError::value_t TrackParametrizationWit // get chi2 wrt other track, which must be defined at the same parameters X,alpha // Supplied non-initialized covToSet matrix is filled by inverse combined matrix for further use - if (std::abs(this->getAlpha() - rhs.getAlpha()) > FLT_EPSILON) { + if (gpu::CAMath::Abs(this->getAlpha() - rhs.getAlpha()) > FLT_EPSILON) { LOG(ERROR) << "The reference Alpha of the tracks differ: " << this->getAlpha() << " : " << rhs.getAlpha(); return 2.f * HugeF; } - if (std::abs(this->getX() - rhs.getX()) > FLT_EPSILON) { + if (gpu::CAMath::Abs(this->getX() - rhs.getX()) > FLT_EPSILON) { LOG(ERROR) << "The reference X of the tracks differ: " << this->getX() << " : " << rhs.getX(); return 2.f * HugeF; } @@ -758,11 +758,11 @@ bool TrackParametrizationWithError::update(const TrackParametrizationWi // update track with other track, the inverted combined cov matrix should be supplied // consider skipping this check, since it is usually already done upstream - if (std::abs(this->getAlpha() - rhs.getAlpha()) > FLT_EPSILON) { + if (gpu::CAMath::Abs(this->getAlpha() - rhs.getAlpha()) > FLT_EPSILON) { LOG(ERROR) << "The reference Alpha of the tracks differ: " << this->getAlpha() << " : " << rhs.getAlpha(); return false; } - if (std::abs(this->getX() - rhs.getX()) > FLT_EPSILON) { + if (gpu::CAMath::Abs(this->getX() - rhs.getX()) > FLT_EPSILON) { LOG(ERROR) << "The reference X of the tracks differ: " << this->getX() << " : " << rhs.getX(); return false; } @@ -837,7 +837,7 @@ bool TrackParametrizationWithError::update(const TrackParametrizationWi //______________________________________________ template -bool TrackParametrizationWithError::update(const dim2_t& p, const dim3_t& cov) +GPUd() bool TrackParametrizationWithError::update(const dim2_t& p, const dim3_t& cov) { // Update the track parameters with the space point "p" having // the covariance matrix "cov" @@ -853,7 +853,7 @@ bool TrackParametrizationWithError::update(const dim2_t& p, const dim3_ double r11 = static_cast(cov[2]) + static_cast(cm11); double det = r00 * r11 - r01 * r01; - if (std::fabs(det) < constants::math::Almost0) { + if (gpu::CAMath::Abs(det) < constants::math::Almost0) { return false; } double detI = 1. / det; @@ -870,7 +870,7 @@ bool TrackParametrizationWithError::update(const dim2_t& p, const dim3_ value_t dy = p[kY] - this->getY(), dz = p[kZ] - this->getZ(); value_t dsnp = k20 * dy + k21 * dz; - if (std::fabs(this->getSnp() + dsnp) > constants::math::Almost1) { + if (gpu::CAMath::Abs(this->getSnp() + dsnp) > constants::math::Almost1) { return false; } @@ -908,7 +908,7 @@ bool TrackParametrizationWithError::update(const dim2_t& p, const dim3_ //______________________________________________ template -bool TrackParametrizationWithError::correctForMaterial(value_t x2x0, value_t xrho, value_t mass, bool anglecorr, value_t dedx) +GPUd() bool TrackParametrizationWithError::correctForMaterial(value_t x2x0, value_t xrho, value_t mass, bool anglecorr, value_t dedx) { //------------------------------------------------------------------ // This function corrects the track parameters for the crossed material. @@ -933,7 +933,7 @@ bool TrackParametrizationWithError::correctForMaterial(value_t x2x0, va value_t cst2I = (1.f + this->getTgl() * this->getTgl()); // 1/cos(lambda)^2 // Apply angle correction, if requested if (anglecorr) { - value_t angle = std::sqrt(cst2I / csp2); + value_t angle = gpu::CAMath::Sqrt(cst2I / csp2); x2x0 *= angle; xrho *= angle; } @@ -946,7 +946,7 @@ bool TrackParametrizationWithError::correctForMaterial(value_t x2x0, va // Calculating the multiple scattering corrections****************** value_t cC22(0.f), cC33(0.f), cC43(0.f), cC44(0.f); if (x2x0 != 0.f) { - value_t theta2 = kMSConst2 / (beta2 * p2) * std::fabs(x2x0); + value_t theta2 = kMSConst2 / (beta2 * p2) * gpu::CAMath::Abs(x2x0); if (this->getAbsCharge() != 1) { theta2 *= this->getAbsCharge() * this->getAbsCharge(); } @@ -970,15 +970,15 @@ bool TrackParametrizationWithError::correctForMaterial(value_t x2x0, va value_t cP4 = 1.f; if ((xrho != 0.f) && (beta2 < 1.f)) { if (dedx < kCalcdEdxAuto + constants::math::Almost1) { // request to calculate dedx on the fly - dedx = BetheBlochSolid(p / std::fabs(mass)); + dedx = BetheBlochSolid(p / gpu::CAMath::Abs(mass)); if (this->getAbsCharge() != 1) { dedx *= this->getAbsCharge() * this->getAbsCharge(); } } value_t dE = dedx * xrho; - value_t e = std::sqrt(e2); - if (std::fabs(dE) > kMaxELossFrac * e) { + value_t e = gpu::CAMath::Sqrt(e2); + if (gpu::CAMath::Abs(dE) > kMaxELossFrac * e) { return false; // 30% energy loss is too much! } value_t eupd = e + dE; @@ -986,11 +986,11 @@ bool TrackParametrizationWithError::correctForMaterial(value_t x2x0, va if (pupd2 < kMinP * kMinP) { return false; } - cP4 = p / std::sqrt(pupd2); + cP4 = p / gpu::CAMath::Sqrt(pupd2); // // Approximate energy loss fluctuation (M.Ivanov) constexpr value_t knst = 0.07f; // To be tuned. - value_t sigmadE = knst * std::sqrt(std::fabs(dE)) * e / p2 * this->getCharge2Pt(); + value_t sigmadE = knst * gpu::CAMath::Sqrt(gpu::CAMath::Abs(dE)) * e / p2 * this->getCharge2Pt(); cC44 += sigmadE * sigmadE; } @@ -1008,7 +1008,7 @@ bool TrackParametrizationWithError::correctForMaterial(value_t x2x0, va //______________________________________________________________ template -bool TrackParametrizationWithError::getCovXYZPxPyPzGlo(gpu::gpustd::array& cv) const +GPUd() bool TrackParametrizationWithError::getCovXYZPxPyPzGlo(gpu::gpustd::array& cv) const { //--------------------------------------------------------------------- // This function returns the global covariance matrix of the track params @@ -1022,7 +1022,7 @@ bool TrackParametrizationWithError::getCovXYZPxPyPzGlo(gpu::gpustd::arr // // Results for (nearly) straight tracks are meaningless ! //--------------------------------------------------------------------- - if (std::abs(this->getQ2Pt()) <= constants::math::Almost0 || std::abs(this->getSnp()) > constants::math::Almost1) { + if (gpu::CAMath::Abs(this->getQ2Pt()) <= constants::math::Almost0 || gpu::CAMath::Abs(this->getSnp()) > constants::math::Almost1) { for (int i = 0; i < 21; i++) { cv[i] = 0.; } @@ -1032,7 +1032,7 @@ bool TrackParametrizationWithError::getCovXYZPxPyPzGlo(gpu::gpustd::arr auto pt = this->getPt(); value_t sn, cs; o2::math_utils::detail::sincos(this->getAlpha(), sn, cs); - auto r = std::sqrt((1. - this->getSnp()) * (1. + this->getSnp())); + auto r = gpu::CAMath::Sqrt((1. - this->getSnp()) * (1. + this->getSnp())); auto m00 = -sn, m10 = cs; auto m23 = -pt * (sn + this->getSnp() * cs / r), m43 = -pt * pt * (r * cs - this->getSnp() * sn); auto m24 = pt * (cs - this->getSnp() * sn / r), m44 = -pt * pt * (r * sn + this->getSnp() * cs); @@ -1089,7 +1089,7 @@ std::string TrackParametrizationWithError::asString() const //______________________________________________________________ template -void TrackParametrizationWithError::print() const +GPUd() void TrackParametrizationWithError::print() const { // print parameters #ifndef GPUCA_ALIGPUCODE @@ -1097,8 +1097,10 @@ void TrackParametrizationWithError::print() const #endif } +namespace o2::track +{ template class TrackParametrizationWithError; +#ifndef GPUCA_GPUCODE_DEVICE template class TrackParametrizationWithError; - -} // namespace track -} // namespace o2 +#endif +} // namespace o2::track diff --git a/DataFormats/Reconstruction/src/Vertex.cxx b/DataFormats/Reconstruction/src/Vertex.cxx index f4a98e4223992..4c85cbdae24ca 100644 --- a/DataFormats/Reconstruction/src/Vertex.cxx +++ b/DataFormats/Reconstruction/src/Vertex.cxx @@ -17,7 +17,7 @@ namespace o2 namespace dataformats { -#ifndef GPUCA_ALIGPUCODE +#ifndef GPUCA_GPUCODE_DEVICE std::string VertexBase::asString() const { diff --git a/Detectors/Base/src/Propagator.cxx b/Detectors/Base/src/Propagator.cxx index c9e83f16436b4..e778898792b7a 100644 --- a/Detectors/Base/src/Propagator.cxx +++ b/Detectors/Base/src/Propagator.cxx @@ -9,7 +9,7 @@ // or submit itself to any jurisdiction. #include "DetectorsBase/Propagator.h" -#include +#include "GPUCommonLogger.h" #include "Field/MagFieldFast.h" #include "MathUtils/Utils.h" #include "ReconstructionDataFormats/Vertex.h" diff --git a/GPU/Common/GPUCommonLogger.h b/GPU/Common/GPUCommonLogger.h index 103d0638f0c53..fc5fa6aac505d 100644 --- a/GPU/Common/GPUCommonLogger.h +++ b/GPU/Common/GPUCommonLogger.h @@ -14,13 +14,9 @@ #ifndef GPUCOMMONFAIRLOGGER_H #define GPUCOMMONFAIRLOGGER_H -#if defined(__OPENCL__) -#define LOG(...) -#define LOGF(...) -#define LOGP(...) - -#elif defined(GPUCA_GPUCODE_DEVICE) #include "GPUCommonDef.h" + +#if defined(GPUCA_GPUCODE_DEVICE) namespace o2::gpu::detail { struct DummyLogger { @@ -31,6 +27,14 @@ struct DummyLogger { } }; } // namespace o2::gpu::detail +#endif + +#if defined(__OPENCL__) +#define LOG(...) o2::gpu::detail::DummyLogger() +#define LOGF(...) +#define LOGP(...) + +#elif defined(GPUCA_GPUCODE_DEVICE) #define LOG(...) o2::gpu::detail::DummyLogger() //#define LOG(...) static_assert(false, "LOG(...) << ... unsupported in GPU code"); #define LOGF(type, string, ...) \ diff --git a/GPU/GPUTracking/Base/GPUReconstructionIncludesDevice.h b/GPU/GPUTracking/Base/GPUReconstructionIncludesDevice.h index 2a38c3f5a2ef9..ff17a30c19bd6 100644 --- a/GPU/GPUTracking/Base/GPUReconstructionIncludesDevice.h +++ b/GPU/GPUTracking/Base/GPUReconstructionIncludesDevice.h @@ -58,6 +58,10 @@ using namespace GPUCA_NAMESPACE::gpu; #include "MatLayerCyl.cxx" #include "Ray.cxx" +// O2 track model +#include "TrackParametrization.cxx" +#include "TrackParametrizationWithError.cxx" + // Files for GPU dEdx #include "GPUdEdx.cxx" diff --git a/GPU/GPUTracking/Base/cuda/CMakeLists.txt b/GPU/GPUTracking/Base/cuda/CMakeLists.txt index 8ec7ef7ce6b87..36b776e059cd5 100644 --- a/GPU/GPUTracking/Base/cuda/CMakeLists.txt +++ b/GPU/GPUTracking/Base/cuda/CMakeLists.txt @@ -115,6 +115,7 @@ if(ALIGPU_BUILD_TYPE STREQUAL "O2") PRIVATE_INCLUDE_DIRECTORIES ${CMAKE_SOURCE_DIR}/Detectors/Base/src ${CMAKE_SOURCE_DIR}/Detectors/TRD/base/src + ${CMAKE_SOURCE_DIR}/DataFormats/Reconstruction/src PUBLIC_LINK_LIBRARIES O2::GPUTracking O2::ITStrackingCUDA TARGETVARNAME targetName) diff --git a/GPU/GPUTracking/Base/hip/CMakeLists.txt b/GPU/GPUTracking/Base/hip/CMakeLists.txt index 62a171ef7cf28..775b58d6d9198 100644 --- a/GPU/GPUTracking/Base/hip/CMakeLists.txt +++ b/GPU/GPUTracking/Base/hip/CMakeLists.txt @@ -31,6 +31,7 @@ if(ALIGPU_BUILD_TYPE STREQUAL "O2") PUBLIC_LINK_LIBRARIES O2::GPUTracking O2::ITStrackingHIP hip::host hip::device hip::hipcub roc::rocthrust PUBLIC_INCLUDE_DIRECTORIES ${CMAKE_SOURCE_DIR}/Detectors/TRD/base/src ${CMAKE_SOURCE_DIR}/Detectors/Base/src + ${CMAKE_SOURCE_DIR}/DataFormats/Reconstruction/src TARGETVARNAME targetName) target_compile_definitions( diff --git a/GPU/GPUTracking/Base/opencl2/CMakeLists.txt b/GPU/GPUTracking/Base/opencl2/CMakeLists.txt index da8f89a2113e2..7e9593ff45cbe 100644 --- a/GPU/GPUTracking/Base/opencl2/CMakeLists.txt +++ b/GPU/GPUTracking/Base/opencl2/CMakeLists.txt @@ -30,6 +30,7 @@ set(OCL_DEFINECL "-D$,$-I>" -I${CMAKE_SOURCE_DIR}/Detectors/TRD/base/src -I${CMAKE_SOURCE_DIR}/Detectors/Base/src + -I${CMAKE_SOURCE_DIR}/DataFormats/Reconstruction/src -I${CMAKE_SOURCE_DIR}/Detectors/ITSMFT/ITS/tracking/cuda/include -DGPUCA_GPULIBRARY=OCL2 -D__OPENCLCPP__ diff --git a/GPU/GPUTracking/Standalone/CMakeLists.txt b/GPU/GPUTracking/Standalone/CMakeLists.txt index 04361760182b6..5a3ad4dd18f9e 100644 --- a/GPU/GPUTracking/Standalone/CMakeLists.txt +++ b/GPU/GPUTracking/Standalone/CMakeLists.txt @@ -147,6 +147,7 @@ include_directories(${CMAKE_SOURCE_DIR}/TPCClusterFinder ${CMAKE_SOURCE_DIR}/../../../DataFormats/Headers/include ${CMAKE_SOURCE_DIR}/../../../DataFormats/MemoryResources/include ${CMAKE_SOURCE_DIR}/../../../DataFormats/Reconstruction/include + ${CMAKE_SOURCE_DIR}/../../../DataFormats/Reconstruction/src ${CMAKE_SOURCE_DIR}/../../../DataFormats/simulation/include ${CMAKE_SOURCE_DIR}/../../../Detectors/Base/include ${CMAKE_SOURCE_DIR}/../../../Detectors/Base/src @@ -183,6 +184,8 @@ if(CONFIG_O2_EXTENSIONS) ../../../DataFormats/Reconstruction/src/TrackParametrizationWithError.cxx ../../../DataFormats/Reconstruction/src/Vertex.cxx ../../../DataFormats/Reconstruction/src/TrackLTIntegral.cxx + ../../../DataFormats/Reconstruction/src/TrackParametrization.cxx + ../../../DataFormats/Reconstruction/src/TrackParametrizationWithError.cxx ../../../Detectors/TRD/base/src/GeometryBase.cxx ../../../Detectors/Base/src/MatLayerCylSet.cxx ../../../Detectors/Base/src/MatLayerCyl.cxx From ba9650255e76f6725437727cb000f004f9536692 Mon Sep 17 00:00:00 2001 From: David Rohr Date: Fri, 11 Dec 2020 21:42:30 +0100 Subject: [PATCH 9/9] GPU: Workaround for reoccuring HIP problem parsing host-code during device compilation --- GPU/Common/GPUCommonLogger.h | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/GPU/Common/GPUCommonLogger.h b/GPU/Common/GPUCommonLogger.h index fc5fa6aac505d..5cabc470ad865 100644 --- a/GPU/Common/GPUCommonLogger.h +++ b/GPU/Common/GPUCommonLogger.h @@ -16,12 +16,12 @@ #include "GPUCommonDef.h" -#if defined(GPUCA_GPUCODE_DEVICE) +#if defined(GPUCA_GPUCODE_DEVICE) || defined(__HIPCC__) namespace o2::gpu::detail { struct DummyLogger { template - GPUd() DummyLogger& operator<<(Args... args) + GPUhd() DummyLogger& operator<<(Args... args) { return *this; } @@ -34,7 +34,7 @@ struct DummyLogger { #define LOGF(...) #define LOGP(...) -#elif defined(GPUCA_GPUCODE_DEVICE) +#elif defined(GPUCA_GPUCODE_DEVICE) || defined(__HIPCC__) #define LOG(...) o2::gpu::detail::DummyLogger() //#define LOG(...) static_assert(false, "LOG(...) << ... unsupported in GPU code"); #define LOGF(type, string, ...) \