ITS: speedup vertexer

Signed-off-by: Felix Schlepper <felix.schlepper@cern.ch>

Author: f3sch

Detectors/ITSMFT/ITS/tracking/include/ITStracking/ClusterLines.h

@@ -31,12 +31,15 @@ struct Line final {
 
   Line() = default;
   Line(const Tracklet&, const Cluster*, const Cluster*);
+  bool operator==(const Line&) const = default;
 
+  static float getDistance2FromPoint(const Line& line, const std::array<float, 3>& point);
   static float getDistanceFromPoint(const Line& line, const std::array<float, 3>& point);
   static SMatrix3f getDCAComponents(const Line& line, const std::array<float, 3>& point);
+  static float getDCA2(const Line&, const Line&, const float precision = constants::Tolerance);
   static float getDCA(const Line&, const Line&, const float precision = constants::Tolerance);
   bool isEmpty() const noexcept;
-  bool operator==(const Line&) const = default;
+  void print() const;
 
   SVector3f originPoint;
   SVector3f cosinesDirector;

Detectors/ITSMFT/ITS/tracking/include/ITStracking/TimeFrame.h

@@ -108,6 +108,8 @@ struct TimeFrame {
 
   float getBeamX() const { return mBeamPos[0]; }
   float getBeamY() const { return mBeamPos[1]; }
+  std::array<float, 2>& getBeamXY() { return mBeamPos; }
+
   auto& getMinRs() { return mMinR; }
   auto& getMaxRs() { return mMaxR; }
   float getMinR(int layer) const { return mMinR[layer]; }
@@ -230,7 +232,6 @@ struct TimeFrame {
   gsl::span<int> getExclusiveNTrackletsCluster(int rofId, int combId);
   uint32_t getTotalTrackletsTF(const int iLayer) { return mTotalTracklets[iLayer]; }
   int getTotalClustersPerROFrange(int rofMin, int range, int layerId) const;
-  std::array<float, 2>& getBeamXY() { return mBeamPos; }
   // \Vertexer
 
   int hasBogusClusters() const { return std::accumulate(mBogusClusters.begin(), mBogusClusters.end(), 0); }

Detectors/ITSMFT/ITS/tracking/src/ClusterLines.cxx

@@ -10,6 +10,7 @@
 // or submit itself to any jurisdiction.
 
 #include <cmath>
+#include "Framework/Logger.h"
 #include "ITStracking/ClusterLines.h"
 
 namespace o2::its
@@ -27,16 +28,21 @@ Line::Line(const Tracklet& tracklet, const Cluster* innerClusters, const Cluster
   cosinesDirector /= std::sqrt(ROOT::Math::Dot(cosinesDirector, cosinesDirector));
 }
 
-float Line::getDistanceFromPoint(const Line& line, const std::array<float, 3>& point)
+float Line::getDistance2FromPoint(const Line& line, const std::array<float, 3>& point)
 {
   const SVector3f p(point.data(), 3);
   const SVector3f delta = p - line.originPoint;
   const float proj = ROOT::Math::Dot(delta, line.cosinesDirector);
   const SVector3f residual = delta - proj * line.cosinesDirector;
-  return std::sqrt(ROOT::Math::Dot(residual, residual));
+  return ROOT::Math::Dot(residual, residual);
 }
 
-float Line::getDCA(const Line& firstLine, const Line& secondLine, const float precision)
+float Line::getDistanceFromPoint(const Line& line, const std::array<float, 3>& point)
+{
+  return std::sqrt(getDistance2FromPoint(line, point));
+}
+
+float Line::getDCA2(const Line& firstLine, const Line& secondLine, const float precision)
 {
   const SVector3f n = ROOT::Math::Cross(firstLine.cosinesDirector, secondLine.cosinesDirector);
   const float norm2 = ROOT::Math::Dot(n, n);
@@ -46,11 +52,17 @@ float Line::getDCA(const Line& firstLine, const Line& secondLine, const float pr
     const SVector3f d = secondLine.originPoint - firstLine.originPoint;
     const float proj = ROOT::Math::Dot(d, firstLine.cosinesDirector);
     const SVector3f residual = d - proj * firstLine.cosinesDirector;
-    return std::sqrt(ROOT::Math::Dot(residual, residual));
+    return ROOT::Math::Dot(residual, residual);
   }
 
   const SVector3f delta = secondLine.originPoint - firstLine.originPoint;
-  return std::abs(ROOT::Math::Dot(delta, n)) / std::sqrt(norm2);
+  const float numerator = ROOT::Math::Dot(delta, n);
+  return (numerator * numerator) / norm2;
+}
+
+float Line::getDCA(const Line& firstLine, const Line& secondLine, const float precision)
+{
+  return std::sqrt(getDCA2(firstLine, secondLine, precision));
 }
 
 Line::SMatrix3f Line::getDCAComponents(const Line& line, const std::array<float, 3>& point)
@@ -77,6 +89,14 @@ bool Line::isEmpty() const noexcept
          ROOT::Math::Dot(cosinesDirector, cosinesDirector) == 0.f;
 }
 
+void Line::print() const
+{
+  LOGP(info, "\tLine: originPoint = ({}, {}, {}), cosinesDirector = ({}, {}, {}) ts={}+-{}",
+       originPoint(0), originPoint(1), originPoint(2),
+       cosinesDirector(0), cosinesDirector(1), cosinesDirector(2),
+       mTime.getTimeStamp(), mTime.getTimeStampError());
+}
+
 // Accumulate the weighted normal equation contributions (A matrix and B vector)
 // from a single line into the running sums. The covariance is assumed to be
 // diagonal and uniform ({1,1,1}) so the weights simplify accordingly.
@@ -124,8 +144,8 @@ ClusterLines::ClusterLines(const int firstLabel, const Line& firstLine, const in
   mRMS2 += (tmpRMS2 - mRMS2) * (1.f / static_cast<float>(getSize()));
 
   // AvgDistance2
-  mAvgDistance2 = Line::getDistanceFromPoint(firstLine, mVertex) * Line::getDistanceFromPoint(firstLine, mVertex);
-  mAvgDistance2 += (Line::getDistanceFromPoint(secondLine, mVertex) * Line::getDistanceFromPoint(secondLine, mVertex) - mAvgDistance2) / (float)getSize();
+  mAvgDistance2 = Line::getDistance2FromPoint(firstLine, mVertex);
+  mAvgDistance2 += (Line::getDistance2FromPoint(secondLine, mVertex) - mAvgDistance2) / (float)getSize();
 }
 
 void ClusterLines::add(const int lineLabel, const Line& line)
@@ -135,7 +155,7 @@ void ClusterLines::add(const int lineLabel, const Line& line)
 
   accumulate(line);
   computeClusterCentroid();
-  mAvgDistance2 += (Line::getDistanceFromPoint(line, mVertex) * Line::getDistanceFromPoint(line, mVertex) - mAvgDistance2) / (float)getSize();
+  mAvgDistance2 += (Line::getDistance2FromPoint(line, mVertex) - mAvgDistance2) / (float)getSize();
 }
 
 void ClusterLines::computeClusterCentroid()

Detectors/ITSMFT/ITS/tracking/src/TrackerTraits.cxx

@@ -981,7 +981,6 @@ void TrackerTraits<NLayers>::setNThreads(int n, std::shared_ptr<tbb::task_arena>
     LOGP(info, "Setting tracker with {} threads.", n);
   } else {
     mTaskArena = arena;
-    LOGP(info, "Attaching tracker to calling thread's arena");
   }
 #endif
 }