-
Notifications
You must be signed in to change notification settings - Fork 511
ITS: new CPU + GPU seeding vertexer #15733
New issue
Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.
By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.
Already on GitHub? Sign in to your account
base: dev
Are you sure you want to change the base?
Changes from all commits
File filter
Filter by extension
Conversations
Jump to
Diff view
Diff view
There are no files selected for viewing
| Original file line number | Diff line number | Diff line change |
|---|---|---|
| @@ -0,0 +1,194 @@ | ||
| // Copyright 2019-2020 CERN and copyright holders of ALICE O2. | ||
| // See https://alice-o2.web.cern.ch/copyright for details of the copyright holders. | ||
| // All rights not expressly granted are reserved. | ||
| // | ||
| // This software is distributed under the terms of the GNU General Public | ||
| // License v3 (GPL Version 3), copied verbatim in the file "COPYING". | ||
| // | ||
| // 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 ClusterLinesGPU.h | ||
| /// \brief device-side line + N-line vertex fit for the GPU seeding vertexer. | ||
|
|
||
| #ifndef O2_ITS_CLUSTERLINES_GPU_H | ||
| #define O2_ITS_CLUSTERLINES_GPU_H | ||
|
|
||
| #include "DataFormatsITS/TimeEstBC.h" | ||
| #include "GPUCommonDef.h" | ||
| #include "GPUCommonMath.h" | ||
| #include "ITStracking/LineProjection.h" | ||
|
|
||
| namespace o2::its::gpu | ||
| { | ||
|
|
||
| using LineTime = o2::its::LineTime; | ||
| using LineWindow = o2::its::LineWindow; | ||
|
|
||
| struct LineProjSoA { | ||
| float* z{nullptr}; // projected z at the beamline; sort key and binary-search key, kept dense | ||
| LineTime* t{nullptr}; // time centre + half-width | ||
| int* idx{nullptr}; // sorted slot -> original line index | ||
| int* rof{nullptr}; // ROF of the line | ||
| }; | ||
|
|
||
| struct VertexCand { | ||
| float x, y, z; | ||
| float rms2[6]; | ||
| float avgDist2; | ||
| int nGood; | ||
| float seed[3]; | ||
| o2::its::TimeEstBC time; | ||
| int size; | ||
| uint8_t ok; // 1 if the candidate passed the fit cuts | ||
| uint8_t keep; // 1 if it survived duplicate suppression (subset of ok) | ||
| uint8_t fine; | ||
| }; | ||
|
|
||
| // Device-side line: origin point + unit direction, with a time stamp | ||
| struct GPULine { | ||
| GPUhdDefault() GPULine() = default; | ||
|
|
||
| GPUhdi() GPULine(const float origin[3], const float direction[3], const o2::its::TimeEstBC& t) : mTime(t) | ||
| { | ||
| const float norm = o2::gpu::GPUCommonMath::Sqrt(direction[0] * direction[0] + | ||
| direction[1] * direction[1] + | ||
| direction[2] * direction[2]); | ||
| const float inv = norm > 0.f ? 1.f / norm : 0.f; | ||
| for (int i = 0; i < 3; ++i) { | ||
| originPoint[i] = origin[i]; | ||
| cosinesDirector[i] = direction[i] * inv; | ||
| } | ||
| } | ||
|
|
||
| // Squared distance from a point to the (infinite) line: |delta - (delta.u) u|^2 | ||
| GPUhdi() static float getDistance2FromPoint(const GPULine& line, const float point[3]) | ||
| { | ||
| float delta[3]; | ||
| float proj = 0.f; | ||
| for (int i = 0; i < 3; ++i) { | ||
| delta[i] = point[i] - line.originPoint[i]; | ||
| proj += delta[i] * line.cosinesDirector[i]; | ||
| } | ||
| float d2 = 0.f; | ||
| for (int i = 0; i < 3; ++i) { | ||
| const float residual = delta[i] - proj * line.cosinesDirector[i]; | ||
| d2 += residual * residual; | ||
| } | ||
| return d2; | ||
| } | ||
|
|
||
| GPUhdi() static void getDCAComponents(const GPULine& line, const float point[3], float out[6]) | ||
| { | ||
| float delta[3]; | ||
| float proj = 0.f; | ||
| for (int i = 0; i < 3; ++i) { | ||
| delta[i] = line.originPoint[i] - point[i]; | ||
| proj += delta[i] * line.cosinesDirector[i]; | ||
| } | ||
| float r[3]; | ||
| for (int i = 0; i < 3; ++i) { | ||
| r[i] = delta[i] - proj * line.cosinesDirector[i]; | ||
| } | ||
| out[0] = r[0]; // (0,0) XX | ||
| out[1] = o2::gpu::GPUCommonMath::Hypot(r[0], r[1]); // (0,1) XY | ||
| out[2] = r[1]; // (1,1) YY | ||
| out[3] = o2::gpu::GPUCommonMath::Hypot(r[0], r[2]); // (0,2) XZ | ||
| out[4] = o2::gpu::GPUCommonMath::Hypot(r[1], r[2]); // (1,2) YZ | ||
| out[5] = r[2]; // (2,2) ZZ | ||
| } | ||
|
|
||
| float originPoint[3] = {0.f, 0.f, 0.f}; | ||
| float cosinesDirector[3] = {0.f, 0.f, 0.f}; | ||
| o2::its::TimeEstBC mTime; | ||
| }; | ||
|
|
||
| class GPUClusterLinesFit | ||
| { | ||
| public: | ||
| GPUhdDefault() GPUClusterLinesFit() = default; | ||
|
|
||
| // Add one line's contribution: A_ij += (delta_ij*|d|^2 - d_i*d_j)/|d|^2, | ||
| // b_i += (d_i*(d.o) - |d|^2*o_i)/|d|^2. For a unit director |d|^2 == 1. | ||
| GPUhdi() void add(const GPULine& line) | ||
| { | ||
| const double d0 = line.cosinesDirector[0], d1 = line.cosinesDirector[1], d2 = line.cosinesDirector[2]; | ||
| const double o0 = line.originPoint[0], o1 = line.originPoint[1], o2 = line.originPoint[2]; | ||
| const double det = d0 * d0 + d1 * d1 + d2 * d2; // == 1 for a normalised director | ||
| if (det <= 0.) { | ||
| return; | ||
| } | ||
| if (mNContributors <= 0) { | ||
| mTime = line.mTime; | ||
| } else { | ||
| mTime += line.mTime; | ||
| } | ||
| mA[0] += (det - d0 * d0) / det; | ||
| mA[1] += (-d0 * d1) / det; | ||
| mA[2] += (-d0 * d2) / det; | ||
| mA[3] += (det - d1 * d1) / det; | ||
| mA[4] += (-d1 * d2) / det; | ||
| mA[5] += (det - d2 * d2) / det; | ||
| const double dDotO = d0 * o0 + d1 * o1 + d2 * o2; | ||
| mB[0] += (d0 * dDotO - det * o0) / det; | ||
| mB[1] += (d1 * dDotO - det * o1) / det; | ||
| mB[2] += (d2 * dDotO - det * o2) / det; | ||
| ++mNContributors; | ||
| } | ||
|
|
||
| // Solve the symmetric system and write the vertex (= -A^-1 B) | ||
| GPUhdi() bool solve(float vertex[3]) const | ||
| { | ||
| const double a = mA[0], b = mA[1], c = mA[2], d = mA[3], e = mA[4], f = mA[5]; | ||
| const double c00 = d * f - e * e; | ||
| const double c01 = c * e - b * f; | ||
| const double c02 = b * e - c * d; | ||
| const double c11 = a * f - c * c; | ||
| const double c12 = b * c - a * e; | ||
| const double c22 = a * d - b * b; | ||
| const double det = a * c00 + b * c01 + c * c02; | ||
| if (o2::gpu::GPUCommonMath::Abs(det) < 1.e-12) { | ||
| return false; | ||
| } | ||
| const double invDet = 1. / det; | ||
| const double x0 = (c00 * mB[0] + c01 * mB[1] + c02 * mB[2]) * invDet; | ||
| const double x1 = (c01 * mB[0] + c11 * mB[1] + c12 * mB[2]) * invDet; | ||
| const double x2 = (c02 * mB[0] + c12 * mB[1] + c22 * mB[2]) * invDet; | ||
| vertex[0] = static_cast<float>(-x0); | ||
| vertex[1] = static_cast<float>(-x1); | ||
| vertex[2] = static_cast<float>(-x2); | ||
| return true; | ||
| } | ||
|
|
||
| GPUhdi() void addResidual(const GPULine& line, const float vertex[3]) | ||
| { | ||
| float dca[6]; | ||
| GPULine::getDCAComponents(line, vertex, dca); | ||
| const float d2 = GPULine::getDistance2FromPoint(line, vertex); | ||
| ++mResidualCount; | ||
| const float inv = 1.f / static_cast<float>(mResidualCount); | ||
| for (int i = 0; i < 6; ++i) { | ||
| mRMS2[i] += (dca[i] - mRMS2[i]) * inv; | ||
| } | ||
| mAvgDistance2 += (d2 - mAvgDistance2) * inv; | ||
| } | ||
|
|
||
| GPUhdi() int getNContributors() const { return mNContributors; } | ||
| GPUhdi() const float* getRMS2() const { return mRMS2; } // Packed symmetric covariance in {XX, XY, YY, XZ, YZ, ZZ} order | ||
| GPUhdi() float getAvgDistance2() const { return mAvgDistance2; } | ||
| GPUhdi() const o2::its::TimeEstBC& getTimeStamp() const { return mTime; } | ||
|
|
||
| private: | ||
| double mA[6] = {0., 0., 0., 0., 0., 0.}; | ||
| double mB[3] = {0., 0., 0.}; | ||
| int mNContributors = 0; | ||
| float mRMS2[6] = {0.f, 0.f, 0.f, 0.f, 0.f, 0.f}; | ||
| float mAvgDistance2 = 0.f; | ||
| int mResidualCount = 0; | ||
| o2::its::TimeEstBC mTime; | ||
| }; | ||
|
|
||
| } // namespace o2::its::gpu | ||
|
|
||
| #endif /* O2_ITS_CLUSTERLINES_GPU_H */ |
| Original file line number | Diff line number | Diff line change |
|---|---|---|
|
|
@@ -21,6 +21,7 @@ | |
| #include "ITStracking/Configuration.h" | ||
| #include "ITStracking/TrackExtensionHypothesis.h" | ||
| #include "ITStrackingGPU/Utils.h" | ||
| #include "ITStrackingGPU/ClusterLinesGPU.h" | ||
|
|
||
| namespace o2::its::gpu | ||
| { | ||
|
|
@@ -54,10 +55,14 @@ class TimeFrameGPU : public TimeFrame<NLayers> | |
| void createTrackingFrameInfoDeviceArray(const int = NLayers); | ||
| void loadUnsortedClustersDevice(const int); | ||
| void createUnsortedClustersDeviceArray(const int = NLayers); | ||
| void loadClustersDevice(const int); | ||
| void createClustersDeviceArray(const int = NLayers); | ||
| void loadClustersIndexTables(const int); | ||
| void createClustersIndexTablesArray(const int = NLayers); | ||
| void createClustersDevice(const int); | ||
| void createClustersIndexTables(const int); | ||
| void createClusterRadiiDevice(); | ||
| void uploadClusterRadii(); | ||
| void sortClustersDevice(const int layer, const TrackingParameters& trkParam); | ||
| void createUsedClustersDevice(const int); | ||
| void createUsedClustersDeviceArray(const int = NLayers); | ||
| void loadUsedClustersDevice(); | ||
|
|
@@ -87,6 +92,35 @@ class TimeFrameGPU : public TimeFrame<NLayers> | |
| void createTrackExtensionScratchDevice(const int nThreads, const int maxHypotheses); | ||
| void downloadTrackITSExtDevice(); | ||
|
|
||
| // Seeding-vertexer | ||
| void createClusterOwnersDeviceArray(); | ||
| void createClusterOwnersDevice(); | ||
| void resetClusterOwnersDevice(); | ||
| void createClusterSortScratchDevice(const int layer); | ||
|
|
||
| protected: | ||
| void prepareClusters(const TrackingParameters& trkParam, const int maxLayers) override | ||
| { | ||
| if (maxLayers < NLayers) { // only if former seeding vertexer is run | ||
| TimeFrame<NLayers>::prepareClusters(trkParam, maxLayers); | ||
| } | ||
| } | ||
| void allocateClusterSortStorage(const TrackingParameters& trkParam, const int maxLayers) override | ||
| { | ||
| if (maxLayers < NLayers) { // only if former seeding vertexer is run | ||
| TimeFrame<NLayers>::allocateClusterSortStorage(trkParam, maxLayers); | ||
| } | ||
| } | ||
|
Comment on lines
+102
to
+113
Collaborator
There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. why defined in the header? |
||
|
|
||
| public: | ||
| void createLinesDevice(const int nCells); | ||
| void createDiamondDevice(const Vertex& diamond); | ||
| unsigned int downloadLinesDevice(); | ||
| unsigned int getNLines(); | ||
| const auto& getHostLines() const { return mLinesHost; } | ||
| const auto& getHostLineRof() const { return mLineRofHost; } | ||
| const auto& getHostLineClusters() const { return mLineClustersHost; } | ||
|
|
||
| /// synchronization | ||
| auto& getStream(const size_t stream) { return mGpuStreams[stream]; } | ||
| auto& getStreams() { return mGpuStreams; } | ||
|
|
@@ -111,6 +145,19 @@ class TimeFrameGPU : public TimeFrame<NLayers> | |
| auto& getTrackITSExt() { return mTrackITSExt; } | ||
| auto& getTrackIndices() { return mTrackIndices; } | ||
| Vertex* getDeviceVertices() { return mPrimaryVerticesDevice; } | ||
| int* getDeviceROFramesClusters(const int layer) { return mROFramesClustersDevice[layer]; } | ||
| int* getDeviceClusterSortKeys(const int layer) { return mClusterSortKeysDevice[layer]; } | ||
| int* getDeviceClusterSortPerm(const int layer) { return mClusterSortPermDevice[layer]; } | ||
| Cluster* getDeviceUnsortedClusters(const int layer) { return mUnsortedClustersDevice[layer]; } | ||
| Cluster* getDeviceClusters(const int layer) { return mClustersDevice[layer]; } | ||
| int* getDeviceClustersIndexTable(const int layer) { return mClustersIndexTablesDevice[layer]; } | ||
| const float* getDeviceMinRs() const { return mClusterMinRDevice; } | ||
| const float* getDeviceMaxRs() const { return mClusterMaxRDevice; } | ||
| int* getDeviceROFramesPV() { return mROFramesPVDevice; } | ||
| unsigned char* getDeviceUsedClusters(const int); | ||
| const o2::base::Propagator* getChainPropagator(); | ||
| bool arePersistentTablesLoaded() { return mPersistentTablesLoaded; } | ||
| void setPersistentTablesLoaded(bool setValue) { mPersistentTablesLoaded = setValue; } | ||
|
|
||
| // Hybrid | ||
| TrackITSExt* getDeviceTrackITSExt() { return mTrackITSExtDevice; } | ||
|
|
@@ -119,6 +166,45 @@ class TimeFrameGPU : public TimeFrame<NLayers> | |
| TrackExtensionHypothesis<NLayers>* getDeviceNextTrackExtensionHypotheses() { return mNextTrackExtensionHypothesesDevice; } | ||
| int* getDeviceNeighboursLUT(const int layer) { return mNeighboursLUTDevice[layer]; } | ||
| CellNeighbour** getDeviceArrayNeighbours() { return mNeighboursDeviceArray; } | ||
| unsigned long long** getDeviceArrayClusterOwners() { return mClusterOwnersDeviceArray; } | ||
| GPULine* getDeviceLines() { return mLinesDevice; } | ||
| int* getDeviceLineSlots() { return mLineSlotsDevice; } | ||
| int* getDeviceLineRof() { return mLineRofDevice; } | ||
| int* getDeviceLineClusters() { return mLineClustersDevice; } | ||
| float* getDeviceLineChi2() { return mLineChi2Device; } | ||
| float* getDeviceLinePt() { return mLinePtDevice; } | ||
| float* getDeviceLineZs() { return mLineZsDevice; } | ||
| gpu::LineTime* getDeviceLineTimes() { return mLineTimesDevice; } | ||
| int* getDeviceLineSortedIdx() { return mLinesSortedIdx; } | ||
| LineProjSoA getLineProjSoA() { return {mLineZsDevice, mLineTimesDevice, mLinesSortedIdx, mLineRofDevice}; } | ||
| LineProjSoA getLineProjSortedSoA() { return {mLineZsSortedDevice, mLineTimesSortedDevice, mLinesSortedIdx, mLineRofSortedDevice}; } | ||
| int* getDeviceRofLineOffsets() { return mRofLineOffsetsDevice; } | ||
| int* getDeviceLineDensity() { return mLineDensityDevice; } | ||
| gpu::LineWindow* getDeviceLineWin() { return mLineWinDevice; } | ||
| uint8_t* getDeviceLineIsPeak() { return mLineIsPeakDevice; } | ||
| int* getDeviceLineDensityFine() { return mLineDensityFineDevice; } | ||
| gpu::LineWindow* getDeviceLineWinFine() { return mLineWinFineDevice; } | ||
| uint8_t* getDeviceLineIsPeakFine() { return mLineIsPeakFineDevice; } | ||
| int* getDevicePeakScan() { return mPeakScanDevice; } | ||
| int* getDevicePeakLineIdx() { return mPeakLineIdxDevice; } | ||
| int* getDevicePeakOffsets() { return mPeakOffsetsDevice; } | ||
| const int* getDeviceNPeaks() { return mPeakOffsetsDevice + this->getNrof(1); } | ||
| VertexCand* getDeviceVertexCands() { return mVertexCandsDevice; } | ||
| int* getDeviceMemberOffsets() { return mMemberOffsetsDevice; } | ||
| int* getDeviceMemberLines() { return mMemberLinesDevice; } | ||
| int downloadVertexCandsDevice(); | ||
| int getNMembers() const { return mNMembers; } | ||
| void downloadMemberOffsetsDevice(); // (MC only) | ||
| void createMemberLinesMCDevice(const int nMembers); // (MC only) | ||
| void downloadMemberLinesDevice(); // (MC only) | ||
| const auto& getHostVertexCands() const { return mVertexCandsHost; } | ||
| const auto& getHostPeakOffsets() const { return mPeakOffsetsHost; } | ||
| const auto& getHostMemberOffsets() const { return mMemberOffsetsHost; } | ||
| const auto& getHostMemberLines() const { return mMemberLinesHost; } | ||
| std::vector<o2::MCCompLabel>& getLineLabelFlat() { return mLineLabelFlatHost; } | ||
| const std::vector<o2::MCCompLabel>& getLineLabelFlat() const { return mLineLabelFlatHost; } | ||
| Vertex* getDeviceDiamond() { return mDiamondDevice; } | ||
| std::array<CellNeighbour*, MaxCells>& getDeviceNeighboursAll() { return mNeighboursDevice; } | ||
| CellNeighbour* getDeviceNeighbours(const int layer) { return mNeighboursDevice[layer]; } | ||
| const TrackingFrameInfo** getDeviceArrayTrackingFrameInfo() const { return mTrackingFrameInfoDeviceArray; } | ||
| const Cluster** getDeviceArrayClusters() const { return mClustersDeviceArray; } | ||
|
|
@@ -215,6 +301,11 @@ class TimeFrameGPU : public TimeFrame<NLayers> | |
| const int** mClustersIndexTablesDeviceArray{nullptr}; | ||
| uint8_t** mUsedClustersDeviceArray{nullptr}; | ||
| const int** mROFramesClustersDeviceArray{nullptr}; | ||
| int* mROFramesPVDevice; | ||
| std::array<int*, NLayers> mClusterSortKeysDevice{}; | ||
| std::array<int*, NLayers> mClusterSortPermDevice{}; | ||
| float* mClusterMinRDevice{nullptr}; | ||
| float* mClusterMaxRDevice{nullptr}; | ||
| std::array<Tracklet*, MaxLinks> mTrackletsDevice{}; | ||
| std::array<int*, MaxLinks> mTrackletsLUTDevice{}; | ||
| std::array<int*, MaxCells> mCellsLUTDevice{}; | ||
|
|
@@ -239,6 +330,47 @@ class TimeFrameGPU : public TimeFrame<NLayers> | |
| CellNeighbour** mNeighboursDeviceArray{nullptr}; | ||
| std::array<TrackingFrameInfo*, NLayers> mTrackingFrameInfoDevice{}; | ||
| const TrackingFrameInfo** mTrackingFrameInfoDeviceArray{nullptr}; | ||
| std::array<unsigned long long*, 3> mClusterOwnersDevice{}; | ||
| unsigned long long** mClusterOwnersDeviceArray{nullptr}; | ||
| int* mLineSlotsDevice{nullptr}; | ||
| GPULine* mLinesDevice{nullptr}; | ||
| int* mLineRofDevice{nullptr}; | ||
| int* mLineClustersDevice{nullptr}; | ||
| float* mLineChi2Device{nullptr}; | ||
| float* mLinePtDevice{nullptr}; | ||
| float* mLineZsDevice{nullptr}; | ||
| gpu::LineTime* mLineTimesDevice{nullptr}; | ||
| float* mLineZsSortedDevice{nullptr}; | ||
| gpu::LineTime* mLineTimesSortedDevice{nullptr}; | ||
| int* mLinesSortedIdx{nullptr}; | ||
| int* mLineRofSortedDevice{nullptr}; // per (sorted) line's ROF | ||
| int* mRofLineOffsetsDevice{nullptr}; // CSR offsets into the (rof,z)-sorted lines, size nRofs+1 | ||
| int* mLineDensityDevice{nullptr}; // per (sorted) line: count of time-compatible neighbours in its z-window | ||
| gpu::LineWindow* mLineWinDevice{nullptr}; // per (sorted) line: [lo,hi) bounds of its z-window (sorted coords) | ||
| uint8_t* mLineIsPeakDevice{nullptr}; // per (sorted) line: 1 if it is a local density peak (vertex candidate) | ||
| int* mLineDensityFineDevice{nullptr}; | ||
| gpu::LineWindow* mLineWinFineDevice{nullptr}; | ||
| uint8_t* mLineIsPeakFineDevice{nullptr}; | ||
| int* mPeakScanDevice{nullptr}; // per (sorted) line: number of peaks strictly before it | ||
| int* mPeakLineIdxDevice{nullptr}; // per peak slot: the sorted line index it came from | ||
| int* mPeakOffsetsDevice{nullptr}; // CSR offsets into the compacted peaks | ||
| VertexCand* mVertexCandsDevice{nullptr}; | ||
| int* mMemberOffsetsDevice{nullptr}; | ||
| int* mMemberLinesDevice{nullptr}; | ||
| int mNLinesCapacity{0}; // = nCells the line buffers were sized for | ||
| std::vector<GPULine> mLinesHost; | ||
| std::vector<int> mLineRofHost; | ||
| std::vector<int> mLineClustersHost; | ||
| std::vector<VertexCand> mVertexCandsHost; | ||
| std::vector<int> mPeakOffsetsHost; | ||
| std::vector<int> mMemberOffsetsHost; | ||
| std::vector<int> mMemberLinesHost; | ||
| std::vector<o2::MCCompLabel> mLineLabelFlatHost; | ||
| int mNMembers{0}; | ||
| Vertex* mDiamondDevice{nullptr}; | ||
| bool mPersistentTablesLoaded{false}; | ||
|
Collaborator
There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. Do we actually need this flag, can we not just add to the configuration the equivalent step? |
||
| std::bitset<NLayers> mUnsortedClustersUploaded{}; | ||
| std::bitset<NLayers> mTrackingFrameInfoUploaded{}; | ||
|
Comment on lines
+372
to
+373
Collaborator
There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. These are not used anywhere? |
||
|
|
||
| // State | ||
| Streams mGpuStreams; | ||
|
|
||
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Do we need to duplicate the already existing ClusterLines math, can we not have only one of these classes?