diff --git a/Common/MathUtils/include/MathUtils/CartesianGPU.h b/Common/MathUtils/include/MathUtils/CartesianGPU.h
index 99ae8bb4b85aa..905f2f660f9c0 100644
--- a/Common/MathUtils/include/MathUtils/CartesianGPU.h
+++ b/Common/MathUtils/include/MathUtils/CartesianGPU.h
@@ -25,7 +25,7 @@ namespace detail
 template <typename T, int I>
 struct GPUPoint2D {
   GPUdDefault() GPUPoint2D() = default;
-  GPUd() GPUPoint2D(T a, T b) : xx(a), yy(b) {}
+  GPUhd() GPUPoint2D(T a, T b) : xx(a), yy(b) {}
   GPUhd() float X() const { return xx; }
   GPUhd() float Y() const { return yy; }
   GPUd() float R() const { return o2::gpu::CAMath::Sqrt(xx * xx + yy * yy); }
@@ -38,7 +38,7 @@ struct GPUPoint2D {
 template <typename T, int I>
 struct GPUPoint3D : public GPUPoint2D<T, I> {
   GPUdDefault() GPUPoint3D() = default;
-  GPUd() GPUPoint3D(T a, T b, T c) : GPUPoint2D<T, I>(a, b), zz(c) {}
+  GPUhd() GPUPoint3D(T a, T b, T c) : GPUPoint2D<T, I>(a, b), zz(c) {}
   GPUhd() float Z() const { return zz; }
   GPUd() float R() const { return o2::gpu::CAMath::Sqrt(GPUPoint2D<T, I>::xx * GPUPoint2D<T, I>::xx + GPUPoint2D<T, I>::yy * GPUPoint2D<T, I>::yy + zz * zz); }
   GPUd() void SetZ(float v) { zz = v; }
diff --git a/DataFormats/Reconstruction/include/ReconstructionDataFormats/Vertex.h b/DataFormats/Reconstruction/include/ReconstructionDataFormats/Vertex.h
index a1663a256612c..da68048420f4d 100644
--- a/DataFormats/Reconstruction/include/ReconstructionDataFormats/Vertex.h
+++ b/DataFormats/Reconstruction/include/ReconstructionDataFormats/Vertex.h
@@ -42,7 +42,7 @@ class VertexBase
   static constexpr int kNCov = 6;
   GPUhdDefault() VertexBase() = default;
   GPUhdDefault() ~VertexBase() = default;
-  GPUd() VertexBase(const math_utils::Point3D<float>& pos, const gpu::gpustd::array<float, kNCov>& cov) : mPos(pos), mCov(cov)
+  GPUhd() VertexBase(const math_utils::Point3D<float>& pos, const gpu::gpustd::array<float, kNCov>& cov) : mPos(pos), mCov(cov)
   {
   }
 
@@ -121,7 +121,7 @@ class Vertex : public VertexBase
 
   GPUhdDefault() Vertex() = default;
   GPUhdDefault() ~Vertex() = default;
-  GPUd() Vertex(const math_utils::Point3D<float>& pos, const gpu::gpustd::array<float, kNCov>& cov, ushort nCont, float chi2)
+  GPUhd() Vertex(const math_utils::Point3D<float>& pos, const gpu::gpustd::array<float, kNCov>& cov, ushort nCont, float chi2)
     : VertexBase(pos, cov), mChi2(chi2), mNContributors(nCont)
   {
   }
diff --git a/Detectors/ITSMFT/ITS/macros/test/CMakeLists.txt b/Detectors/ITSMFT/ITS/macros/test/CMakeLists.txt
index befc776289e25..ee9479edb17ec 100644
--- a/Detectors/ITSMFT/ITS/macros/test/CMakeLists.txt
+++ b/Detectors/ITSMFT/ITS/macros/test/CMakeLists.txt
@@ -66,10 +66,6 @@ o2_add_test_root_macro(CheckTracksCA.C
                                              O2::DataFormatsITSMFT
                        LABELS its)
 
-o2_add_test_root_macro(RunGPUTracking.C
-                       PUBLIC_LINK_LIBRARIES O2::GPUTracking
-                       LABELS its)
-
 o2_add_test_root_macro(DisplayTrack.C
                        PUBLIC_LINK_LIBRARIES O2::ITSBase
                                              O2::DataFormatsITSMFT
@@ -95,3 +91,6 @@ o2_add_test_root_macro(ITSMisaligner.C
                        PUBLIC_LINK_LIBRARIES O2::CCDB
                                              O2::ITSReconstruction
                        LABELS its)
+
+o2_add_test_root_macro(CompareArtefacts.C
+                       LABELS its)
diff --git a/Detectors/ITSMFT/ITS/macros/test/CompareArtefacts.C b/Detectors/ITSMFT/ITS/macros/test/CompareArtefacts.C
new file mode 100644
index 0000000000000..26f0b718774c9
--- /dev/null
+++ b/Detectors/ITSMFT/ITS/macros/test/CompareArtefacts.C
@@ -0,0 +1,164 @@
+#if !defined(__CLING__) || defined(__ROOTCLING__)
+#include "TFile.h"
+#include "TTree.h"
+#include "TCanvas.h"
+#include "TH1F.h"
+#include "TLegend.h"
+#include "TMath.h"
+#include "TString.h"
+#include "TStyle.h"
+
+#include <string>
+#include <vector>
+#endif
+
+using std::string;
+using std::vector;
+
+void CompareArtefacts(const string cpu_file = "artefacts_tf.root", const string gpu_file = "artefacts_tf_gpu.root")
+{
+  gStyle->SetOptStat(0);
+
+  auto f_cpu = TFile::Open(cpu_file.data(), "r");
+  auto f_gpu = TFile::Open(gpu_file.data(), "r");
+
+  auto tree_cpu_tracklets = (TTree*)f_cpu->Get("tracklets");
+  auto tree_gpu_tracklets = (TTree*)f_gpu->Get("tracklets");
+
+  auto hist_cpu_tracklets0_phi = new TH1F("hist_cpu_tracklets0_phi", "hist_cpu_tracklets0_phi", 100, -TMath::Pi() - 1, TMath::Pi() + 1);
+  auto hist_gpu_tracklets0_phi = new TH1F("hist_gpu_tracklets0_phi", "hist_gpu_tracklets0_phi", 100, -TMath::Pi() - 1, TMath::Pi() + 1);
+  auto hist_cpu_tracklets0_tanL = new TH1F("hist_cpu_tracklets0_tanL", "hist_cpu_tracklets0_tanL", 100, -80, 80);
+  auto hist_gpu_tracklets0_tanL = new TH1F("hist_gpu_tracklets0_tanL", "hist_gpu_tracklets0_tanL", 100, -80, 80);
+  auto hist_cpu_tracklets0_firstClusterIndex = new TH1F("hist_cpu_tracklets0_firstClusterIndex", "hist_cpu_tracklets0_firstClusterIndex", 200, 0, 600);
+  auto hist_gpu_tracklets0_firstClusterIndex = new TH1F("hist_gpu_tracklets0_firstClusterIndex", "hist_gpu_tracklets0_firstClusterIndex", 200, 0, 600);
+
+  auto hist_cpu_tracklets1_phi = new TH1F("hist_cpu_tracklets1_phi", "hist_cpu_tracklets1_phi", 100, -TMath::Pi() - 1, TMath::Pi() + 1);
+  auto hist_gpu_tracklets1_phi = new TH1F("hist_gpu_tracklets1_phi", "hist_gpu_tracklets1_phi", 100, -TMath::Pi() - 1, TMath::Pi() + 1);
+  auto hist_cpu_tracklets1_tanL = new TH1F("hist_cpu_tracklets1_tanL", "hist_cpu_tracklets1_tanL", 100, -80, 80);
+  auto hist_gpu_tracklets1_tanL = new TH1F("hist_gpu_tracklets1_tanL", "hist_gpu_tracklets1_tanL", 100, -80, 80);
+  auto hist_cpu_tracklets1_firstClusterIndex = new TH1F("hist_cpu_tracklets1_firstClusterIndex", "hist_cpu_tracklets1_firstClusterIndex", 200, 0, 600);
+  auto hist_gpu_tracklets1_firstClusterIndex = new TH1F("hist_gpu_tracklets1_firstClusterIndex", "hist_gpu_tracklets1_firstClusterIndex", 200, 0, 600);
+
+  hist_cpu_tracklets0_phi->SetLineColor(kRed);
+  hist_gpu_tracklets0_phi->SetLineColor(kBlue);
+  hist_cpu_tracklets0_tanL->SetLineColor(kRed);
+  hist_gpu_tracklets0_tanL->SetLineColor(kBlue);
+  hist_cpu_tracklets0_firstClusterIndex->SetLineColor(kRed);
+  hist_gpu_tracklets0_firstClusterIndex->SetLineColor(kBlue);
+
+  hist_cpu_tracklets1_phi->SetLineColor(kRed);
+  hist_gpu_tracklets1_phi->SetLineColor(kBlue);
+  hist_cpu_tracklets1_tanL->SetLineColor(kRed);
+  hist_gpu_tracklets1_tanL->SetLineColor(kBlue);
+  hist_cpu_tracklets1_firstClusterIndex->SetLineColor(kRed);
+  hist_gpu_tracklets1_firstClusterIndex->SetLineColor(kBlue);
+
+  auto c1 = new TCanvas("c1", "c1", 800, 800);
+  c1->cd();
+
+  tree_cpu_tracklets->Draw("Tracklets0.phi >> hist_cpu_tracklets0_phi");
+  tree_gpu_tracklets->Draw("Tracklets0.phi >> hist_gpu_tracklets0_phi");
+
+  hist_gpu_tracklets0_phi->Draw();
+  hist_cpu_tracklets0_phi->Draw("same");
+
+  auto legend = new TLegend(0.4, 0.4, 0.2, 0.2);
+  legend->SetHeader("Tracklets0 #varphi", "C");
+  legend->AddEntry(hist_cpu_tracklets0_phi, Form("CPU: %1.f", hist_cpu_tracklets0_phi->GetEntries()), "l");
+  legend->AddEntry(hist_gpu_tracklets0_phi, Form("GPU: %1.f", hist_gpu_tracklets0_phi->GetEntries()), "l");
+
+  legend->Draw();
+
+  // // Tan(L)
+  auto c2 = new TCanvas("c2", "c2", 800, 800);
+  c2->cd();
+
+  tree_cpu_tracklets->Draw("Tracklets0.tanLambda >> hist_cpu_tracklets0_tanL");
+  tree_gpu_tracklets->Draw("Tracklets0.tanLambda >> hist_gpu_tracklets0_tanL");
+
+  hist_gpu_tracklets0_tanL->Draw();
+  hist_cpu_tracklets0_tanL->Draw("same");
+
+  auto legend2 = new TLegend(0.4, 0.4, 0.2, 0.2);
+  legend2->SetHeader("Tracklets0 tan(#lambda)", "C");
+  legend2->AddEntry(hist_cpu_tracklets0_tanL, Form("CPU: %1.f", hist_cpu_tracklets0_tanL->GetEntries()), "l");
+  legend2->AddEntry(hist_gpu_tracklets0_tanL, Form("GPU: %1.f", hist_gpu_tracklets0_tanL->GetEntries()), "l");
+
+  legend2->Draw();
+
+  // first cluster index
+  auto c3 = new TCanvas("c3", "c3", 800, 800);
+  c3->cd();
+
+  tree_cpu_tracklets->Draw("Tracklets0.firstClusterIndex >> hist_cpu_tracklets0_firstClusterIndex");
+  tree_gpu_tracklets->Draw("Tracklets0.firstClusterIndex >> hist_gpu_tracklets0_firstClusterIndex");
+  auto hist_cpu_tracklets0_firstClusterIndex_clone = (TH1F*)hist_cpu_tracklets0_firstClusterIndex->Clone("hist_cpu_tracklets0_firstClusterIndex_clone");
+  hist_cpu_tracklets0_firstClusterIndex_clone->Add(hist_gpu_tracklets0_firstClusterIndex, -1);
+  hist_cpu_tracklets0_firstClusterIndex_clone->SetLineColor(kBlack);
+
+  hist_gpu_tracklets0_firstClusterIndex->Draw();
+  hist_cpu_tracklets0_firstClusterIndex->Draw("same");
+  hist_cpu_tracklets0_firstClusterIndex_clone->Draw("same");
+
+  auto legend3 = new TLegend(0.4, 0.4, 0.2, 0.2);
+  legend3->SetHeader("Tracklets0 first cluster index", "C");
+  legend3->AddEntry(hist_cpu_tracklets0_firstClusterIndex, Form("CPU: %1.f", hist_cpu_tracklets0_firstClusterIndex->GetEntries()), "l");
+  legend3->AddEntry(hist_gpu_tracklets0_firstClusterIndex, Form("GPU: %1.f", hist_gpu_tracklets0_firstClusterIndex->GetEntries()), "l");
+
+  legend3->Draw();
+
+  // // Tracklets1
+  auto c4 = new TCanvas("c4", "c4", 800, 800);
+  c4->cd();
+
+  tree_cpu_tracklets->Draw("Tracklets1.phi >> hist_cpu_tracklets1_phi");
+  tree_gpu_tracklets->Draw("Tracklets1.phi >> hist_gpu_tracklets1_phi");
+
+  hist_cpu_tracklets1_phi->Draw();
+  hist_gpu_tracklets1_phi->Draw("same");
+
+  auto legend4 = new TLegend(0.4, 0.4, 0.2, 0.2);
+  legend4->SetHeader("Tracklets1 #varphi", "C");
+  legend4->AddEntry(hist_cpu_tracklets1_phi, Form("CPU: %1.f", hist_cpu_tracklets1_phi->GetEntries()), "l");
+  legend4->AddEntry(hist_gpu_tracklets1_phi, Form("GPU: %1.f", hist_gpu_tracklets1_phi->GetEntries()), "l");
+
+  legend4->Draw();
+
+  // // Tan(L)
+  auto c5 = new TCanvas("c5", "c5", 800, 800);
+  c5->cd();
+
+  tree_cpu_tracklets->Draw("Tracklets1.tanLambda >> hist_cpu_tracklets1_tanL");
+  tree_gpu_tracklets->Draw("Tracklets1.tanLambda >> hist_gpu_tracklets1_tanL");
+
+  hist_cpu_tracklets1_tanL->Draw();
+  hist_gpu_tracklets1_tanL->Draw("same");
+
+  auto legend5 = new TLegend(0.4, 0.4, 0.2, 0.2);
+  legend5->SetHeader("Tracklets1 tan(#lambda)", "C");
+  legend5->AddEntry(hist_cpu_tracklets1_tanL, Form("CPU: %1.f", hist_cpu_tracklets1_tanL->GetEntries()), "l");
+  legend5->AddEntry(hist_gpu_tracklets1_tanL, Form("GPU: %1.f", hist_gpu_tracklets1_tanL->GetEntries()), "l");
+
+  legend5->Draw();
+
+  // first cluster index
+  auto c6 = new TCanvas("c6", "c6", 800, 800);
+  c6->cd();
+
+  tree_cpu_tracklets->Draw("Tracklets1.firstClusterIndex >> hist_cpu_tracklets1_firstClusterIndex");
+  tree_gpu_tracklets->Draw("Tracklets1.firstClusterIndex >> hist_gpu_tracklets1_firstClusterIndex");
+  auto hist_cpu_tracklets1_firstClusterIndex_clone = (TH1F*)hist_cpu_tracklets1_firstClusterIndex->Clone("hist_cpu_tracklets1_firstClusterIndex_clone");
+  hist_cpu_tracklets1_firstClusterIndex_clone->Add(hist_gpu_tracklets1_firstClusterIndex, -1);
+  hist_cpu_tracklets1_firstClusterIndex_clone->SetLineColor(kBlack);
+
+  hist_cpu_tracklets1_firstClusterIndex->Draw();
+  hist_gpu_tracklets1_firstClusterIndex->Draw("same");
+  hist_cpu_tracklets1_firstClusterIndex_clone->Draw("same");
+  auto legend6 = new TLegend(0.4, 0.4, 0.2, 0.2);
+
+  legend6->SetHeader("Tracklets1 first cluster index", "C");
+  legend6->AddEntry(hist_cpu_tracklets1_firstClusterIndex, Form("CPU: %1.f", hist_cpu_tracklets1_firstClusterIndex->GetEntries()), "l");
+  legend6->AddEntry(hist_gpu_tracklets1_firstClusterIndex, Form("GPU: %1.f", hist_gpu_tracklets1_firstClusterIndex->GetEntries()), "l");
+
+  legend6->Draw();
+}
\ No newline at end of file
diff --git a/Detectors/ITSMFT/ITS/macros/test/RunGPUTracking.C b/Detectors/ITSMFT/ITS/macros/test/RunGPUTracking.C
deleted file mode 100644
index 6eab73649a68a..0000000000000
--- a/Detectors/ITSMFT/ITS/macros/test/RunGPUTracking.C
+++ /dev/null
@@ -1,186 +0,0 @@
-#if !defined(__CLING__) || defined(__ROOTCLING__)
-#include "ITStracking/Vertexer.h"
-#include "GPUO2Interface.h"
-#include "GPUReconstruction.h"
-#include "GPUChainITS.h"
-
-#include <TChain.h>
-#include <TFile.h>
-#include <TTree.h>
-#include <TGeoGlobalMagField.h>
-#include <FairEventHeader.h>
-#include <FairGeoParSet.h>
-#include <FairLogger.h>
-#include "DetectorsCommonDataFormats/DetectorNameConf.h"
-#include "SimulationDataFormat/MCEventHeader.h"
-
-#include "DataFormatsITSMFT/TopologyDictionary.h"
-#include "DataFormatsITSMFT/CompCluster.h"
-#include "DetectorsCommonDataFormats/DetID.h"
-#include "DataFormatsParameters/GRPObject.h"
-#include "DetectorsBase/GeometryManager.h"
-#include "Field/MagneticField.h"
-#include "ITSBase/GeometryTGeo.h"
-#include "ITStracking/ROframe.h"
-#include "ITStracking/IOUtils.h"
-#include "ITStracking/Vertexer.h"
-#include "ITStracking/Tracker.h"
-#include "ITStrackingGPU/TimeFrameGPU.h"
-#include "MathUtils/Utils.h"
-#include "DetectorsBase/Propagator.h"
-#include "CCDB/BasicCCDBManager.h"
-#include "CCDB/CCDBTimeStampUtils.h"
-
-#include "SimulationDataFormat/MCCompLabel.h"
-#include "SimulationDataFormat/MCTruthContainer.h"
-#include "ITStracking/Configuration.h"
-
-using MCLabCont = o2::dataformats::MCTruthContainer<o2::MCCompLabel>;
-using Vertex = o2::dataformats::Vertex<o2::dataformats::TimeStamp<int>>;
-
-void RunGPUTracking(bool useLUT = true,
-                    std::string path = "./",
-                    std::string outputfile = "o2trac_its.root",
-                    std::string inputClustersITS = "o2clus_its.root",
-                    std::string matLUTFile = "matbud.root",
-                    std::string inputGRP = "o2sim_grp.root",
-                    long timestamp = 0)
-{
-  o2::its::ROframe event(0, 7);
-
-  if (path.back() != '/') {
-    path += '/';
-  }
-
-  //-------- init geometry and field --------//
-  const auto grp = o2::parameters::GRPObject::loadFrom(path + inputGRP);
-  if (!grp) {
-    LOG(fatal) << "Cannot run w/o GRP object";
-  }
-
-  o2::base::GeometryManager::loadGeometry(path);
-  o2::base::Propagator::initFieldFromGRP(grp);
-  auto field = static_cast<o2::field::MagneticField*>(TGeoGlobalMagField::Instance()->GetField());
-  if (!field) {
-    LOG(fatal) << "Failed to load ma";
-  }
-  double origD[3] = {0., 0., 0.};
-  // tracker.setBz(field->getBz(origD));
-
-  //-------- init lookuptable --------//
-  if (useLUT) {
-    auto* lut = o2::base::MatLayerCylSet::loadFromFile(matLUTFile);
-    o2::base::Propagator::Instance()->setMatLUT(lut);
-  } else {
-    // tracker.setCorrType(o2::base::PropagatorImpl<float>::MatCorrType::USEMatCorrTGeo);
-  }
-
-  // if (tracker.isMatLUT()) {
-  //   LOG(info) << "Loaded material LUT from " << matLUTFile;
-  // } else {
-  //   LOG(info) << "Material LUT " << matLUTFile << " file is absent, only TGeo can be used";
-  // }
-
-  auto gman = o2::its::GeometryTGeo::Instance();
-  gman->fillMatrixCache(o2::math_utils::bit2Mask(o2::math_utils::TransformType::T2L, o2::math_utils::TransformType::T2GRot,
-                                                 o2::math_utils::TransformType::L2G)); // request cached transforms
-
-  auto& mgr = o2::ccdb::BasicCCDBManager::instance();
-  mgr.setURL("http://alice-ccdb.cern.ch");
-  mgr.setTimestamp(timestamp ? timestamp : o2::ccdb::getCurrentTimestamp());
-  const o2::itsmft::TopologyDictionary* dict = mgr.get<o2::itsmft::TopologyDictionary>("ITS/Calib/ClusterDictionary");
-
-  //>>>---------- attach input data --------------->>>
-  TChain itsClusters("o2sim");
-  itsClusters.AddFile((path + inputClustersITS).data());
-
-  if (!itsClusters.GetBranch("ITSClusterComp")) {
-    LOG(fatal) << "Did not find ITS clusters branch ITSClusterComp in the input tree";
-  }
-  std::vector<o2::itsmft::CompClusterExt>* cclusters = nullptr;
-  itsClusters.SetBranchAddress("ITSClusterComp", &cclusters);
-
-  if (!itsClusters.GetBranch("ITSClusterPatt")) {
-    LOG(fatal) << "Did not find ITS cluster patterns branch ITSClusterPatt in the input tree";
-  }
-  std::vector<unsigned char>* patterns = nullptr;
-  itsClusters.SetBranchAddress("ITSClusterPatt", &patterns);
-
-  MCLabCont* labels = nullptr;
-  if (!itsClusters.GetBranch("ITSClusterMCTruth")) {
-    LOG(warning) << "Did not find ITS clusters branch ITSClusterMCTruth in the input tree";
-  } else {
-    itsClusters.SetBranchAddress("ITSClusterMCTruth", &labels);
-  }
-
-  if (!itsClusters.GetBranch("ITSClustersROF")) {
-    LOG(fatal) << "Did not find ITS clusters branch ITSClustersROF in the input tree";
-  }
-
-  std::vector<o2::itsmft::MC2ROFRecord>* mc2rofs = nullptr;
-  if (!itsClusters.GetBranch("ITSClustersMC2ROF")) {
-    LOG(fatal) << "Did not find ITS clusters branch ITSClustersROF in the input tree";
-  }
-  itsClusters.SetBranchAddress("ITSClustersMC2ROF", &mc2rofs);
-
-  std::vector<o2::itsmft::ROFRecord>* rofs = nullptr;
-  itsClusters.SetBranchAddress("ITSClustersROF", &rofs);
-  itsClusters.GetEntry(0);
-
-  //-------------------------------------------------
-  std::unique_ptr<o2::gpu::GPUReconstruction> recCUDA(o2::gpu::GPUReconstruction::CreateInstance(o2::gpu::GPUDataTypes::DeviceType::CUDA, true));
-  auto* chainITSCUDA = recCUDA->AddChain<o2::gpu::GPUChainITS>();
-  std::unique_ptr<o2::its::Vertexer> vertexerCUDA = std::make_unique<o2::its::Vertexer>(chainITSCUDA->GetITSVertexerTraits());
-  std::unique_ptr<o2::its::Tracker> trackerCUDA = std::make_unique<o2::its::Tracker>(chainITSCUDA->GetITSTrackerTraits());
-
-  std::unique_ptr<o2::gpu::GPUReconstruction> recHIP(o2::gpu::GPUReconstruction::CreateInstance(o2::gpu::GPUDataTypes::DeviceType::HIP, true));
-  auto* chainITSHIP = recHIP->AddChain<o2::gpu::GPUChainITS>();
-  std::unique_ptr<o2::its::Vertexer> vertexerHIP = std::make_unique<o2::its::Vertexer>(chainITSHIP->GetITSVertexerTraits());
-  std::unique_ptr<o2::its::Tracker> trackerHIP = std::make_unique<o2::its::Tracker>(chainITSHIP->GetITSTrackerTraits());
-
-  o2::its::VertexingParameters parameters;
-  parameters.phiCut = 0.005f;
-  parameters.tanLambdaCut = 0.002f;
-
-  vertexerCUDA->setParameters(parameters);
-  vertexerHIP->setParameters(parameters);
-
-  gsl::span<const unsigned char> patt(patterns->data(), patterns->size());
-  auto pattIt = patt.begin();
-  auto clSpan = gsl::span(cclusters->data(), cclusters->size());
-
-  for (auto& rof : *rofs) {
-
-    auto it = pattIt;
-    o2::its::ioutils::loadROFrameData(rof, event, clSpan, pattIt, dict, labels);
-
-    //   // CUDA
-    vertexerCUDA->initialiseVertexer(&event);
-    vertexerCUDA->findTracklets();
-    // vertexerCUDA.filterMCTracklets(); // to use MC check
-    vertexerCUDA->validateTracklets();
-    vertexerCUDA->findVertices();
-    std::vector<Vertex> vertITSCU = vertexerCUDA->exportVertices();
-    if (!vertITSCU.empty()) {
-      std::cout << " - Reconstructed vertex: x = " << vertITSCU[0].getX() << " y = " << vertITSCU[0].getY() << " x = " << vertITSCU[0].getZ() << std::endl;
-      event.addPrimaryVertex(vertITSCU[0].getX(), vertITSCU[0].getY(), vertITSCU[0].getZ());
-    } else {
-      std::cout << " - Vertex not reconstructed" << std::endl;
-    }
-
-    //   // HIP
-    vertexerHIP->initialiseVertexer(&event);
-    vertexerHIP->findTracklets();
-    //   // vertexerHIP.filterMCTracklets(); // to use MC check
-    vertexerHIP->validateTracklets();
-    vertexerHIP->findVertices();
-    std::vector<Vertex> vertITSHIP = vertexerHIP->exportVertices();
-    if (!vertITSHIP.empty()) {
-      std::cout << " - Reconstructed vertex: x = " << vertITSHIP[0].getX() << " y = " << vertITSHIP[0].getY() << " x = " << vertITSHIP[0].getZ() << std::endl;
-      event.addPrimaryVertex(vertITSHIP[0].getX(), vertITSHIP[0].getY(), vertITSHIP[0].getZ());
-    } else {
-      std::cout << " - Vertex not reconstructed" << std::endl;
-    }
-  }
-}
-#endif
diff --git a/Detectors/ITSMFT/ITS/tracking/GPU/ITStrackingGPU/DeviceStoreGPU.h b/Detectors/ITSMFT/ITS/tracking/GPU/ITStrackingGPU/DeviceStoreGPU.h
deleted file mode 100644
index 107d94e8e3e34..0000000000000
--- a/Detectors/ITSMFT/ITS/tracking/GPU/ITStrackingGPU/DeviceStoreGPU.h
+++ /dev/null
@@ -1,132 +0,0 @@
-// 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 DeviceStoreNV.h
-/// \brief
-///
-
-#ifndef ITSTRACKINGGPU_DEVICESTOREGPU_H_
-#define ITSTRACKINGGPU_DEVICESTOREGPU_H_
-
-#ifndef GPUCA_GPUCODE_GENRTC
-#include <cub/cub.cuh>
-#include <cstdint>
-#endif
-
-#include "ITStracking/Cell.h"
-#include "ITStracking/Configuration.h"
-#include "ITStracking/Cluster.h"
-#include "ITStracking/Constants.h"
-#include "ITStracking/Tracklet.h"
-
-#include "Array.h"
-#include "UniquePointer.h"
-#include "Vector.h"
-
-namespace o2
-{
-namespace its
-{
-
-namespace gpu
-{
-
-class DeviceStoreNV final
-{
- public:
-  DeviceStoreNV();
-
-  UniquePointer<DeviceStoreNV> initialise(const float3&,
-                                          const std::array<std::vector<Cluster>, constants::its2::LayersNumber>&,
-                                          const std::array<std::vector<Tracklet>, constants::its2::TrackletsPerRoad>&,
-                                          const std::array<std::vector<Cell>, constants::its2::CellsPerRoad>&,
-                                          const std::array<std::vector<int>, constants::its2::CellsPerRoad - 1>&,
-                                          const std::array<float, constants::its2::LayersNumber>&,
-                                          const std::array<float, constants::its2::LayersNumber>&);
-  GPUd() const float3& getPrimaryVertex() { return *mPrimaryVertex; };
-  GPUhd() Array<Vector<Cluster>, constants::its2::LayersNumber>& getClusters();
-  GPUd() Array<Array<int, constants::its2::ZBins * constants::its2::PhiBins + 1>,
-               constants::its2::TrackletsPerRoad>& getIndexTables() { return mIndexTables; };
-  GPUhd() Array<Vector<Tracklet>, constants::its2::TrackletsPerRoad>& getTracklets();
-  GPUhd() Array<Vector<int>, constants::its2::CellsPerRoad>& getTrackletsLookupTable();
-  GPUhd() Array<Vector<int>, constants::its2::CellsPerRoad>& getTrackletsPerClusterTable();
-  GPUhd() Array<Vector<Cell>, constants::its2::CellsPerRoad>& getCells();
-  GPUhd() Array<Vector<int>, constants::its2::CellsPerRoad - 1>& getCellsLookupTable();
-  GPUhd() Array<Vector<int>, constants::its2::CellsPerRoad - 1>& getCellsPerTrackletTable();
-  Array<Vector<int>, constants::its2::CellsPerRoad>& getTempTableArray();
-
-  GPUhd() float getRmin(int layer);
-  GPUhd() float getRmax(int layer);
-
- private:
-  UniquePointer<float3> mPrimaryVertex;
-  Array<Vector<Cluster>, constants::its2::LayersNumber> mClusters;
-  Array<float, constants::its2::LayersNumber> mRmin;
-  Array<float, constants::its2::LayersNumber> mRmax;
-  Array<Array<int, constants::its2::ZBins * constants::its2::PhiBins + 1>, constants::its2::TrackletsPerRoad>
-    mIndexTables;
-  Array<Vector<Tracklet>, constants::its2::TrackletsPerRoad> mTracklets;
-  Array<Vector<int>, constants::its2::CellsPerRoad> mTrackletsLookupTable;
-  Array<Vector<int>, constants::its2::CellsPerRoad> mTrackletsPerClusterTable;
-  Array<Vector<Cell>, constants::its2::CellsPerRoad> mCells;
-  Array<Vector<int>, constants::its2::CellsPerRoad - 1> mCellsLookupTable;
-  Array<Vector<int>, constants::its2::CellsPerRoad - 1> mCellsPerTrackletTable;
-};
-
-GPUhd() Array<Vector<Cluster>, constants::its2::LayersNumber>& DeviceStoreNV::getClusters()
-{
-  return mClusters;
-}
-
-GPUd() Array<Vector<Tracklet>, constants::its2::TrackletsPerRoad>& DeviceStoreNV::getTracklets()
-{
-  return mTracklets;
-}
-
-GPUd() Array<Vector<int>, constants::its2::CellsPerRoad>& DeviceStoreNV::getTrackletsLookupTable()
-{
-  return mTrackletsLookupTable;
-}
-
-GPUd() Array<Vector<int>, constants::its2::CellsPerRoad>& DeviceStoreNV::getTrackletsPerClusterTable()
-{
-  return mTrackletsPerClusterTable;
-}
-
-GPUhd() Array<Vector<Cell>, constants::its2::CellsPerRoad>& DeviceStoreNV::getCells()
-{
-  return mCells;
-}
-
-GPUhd() Array<Vector<int>, constants::its2::CellsPerRoad - 1>& DeviceStoreNV::getCellsLookupTable()
-{
-  return mCellsLookupTable;
-}
-
-GPUhd() Array<Vector<int>, constants::its2::CellsPerRoad - 1>& DeviceStoreNV::getCellsPerTrackletTable()
-{
-  return mCellsPerTrackletTable;
-}
-
-GPUhd() float DeviceStoreNV::getRmin(int layer)
-{
-  return mRmin[layer];
-}
-
-GPUhd() float DeviceStoreNV::getRmax(int layer)
-{
-  return mRmax[layer];
-}
-
-} // namespace gpu
-} // namespace its
-} // namespace o2
-#endif
diff --git a/Detectors/ITSMFT/ITS/tracking/GPU/ITStrackingGPU/TimeFrameGPU.h b/Detectors/ITSMFT/ITS/tracking/GPU/ITStrackingGPU/TimeFrameGPU.h
index 5abf07d4cbd9d..60ed5cab65ef7 100644
--- a/Detectors/ITSMFT/ITS/tracking/GPU/ITStrackingGPU/TimeFrameGPU.h
+++ b/Detectors/ITSMFT/ITS/tracking/GPU/ITStrackingGPU/TimeFrameGPU.h
@@ -21,12 +21,14 @@
 #include "ITStracking/Configuration.h"
 
 #include "ITStrackingGPU/ClusterLinesGPU.h"
+#include "ITStrackingGPU/Stream.h"
 
 #include "Array.h"
 #include "Vector.h"
 
 #include "GPUCommonDef.h"
 #include "GPUCommonMath.h"
+#include "GPUCommonLogger.h"
 
 namespace o2
 {
@@ -38,6 +40,50 @@ using namespace constants::its2;
 class TimeFrameGPUConfig;
 namespace gpu
 {
+
+template <int NLayers>
+struct StaticTrackingParameters {
+  StaticTrackingParameters<NLayers>& operator=(const StaticTrackingParameters<NLayers>& t) = default;
+  void set(const TrackingParameters& pars);
+
+  /// General parameters
+  int ClusterSharing = 0;
+  int MinTrackLength = NLayers;
+  float NSigmaCut = 5;
+  float PVres = 1.e-2f;
+  int DeltaROF = 0;
+  int ZBins{256};
+  int PhiBins{128};
+
+  /// Cell finding cuts
+  float CellDeltaTanLambdaSigma = 0.007f;
+};
+
+template <int NLayers>
+void StaticTrackingParameters<NLayers>::set(const TrackingParameters& pars)
+{
+  ClusterSharing = pars.ClusterSharing;
+  MinTrackLength = pars.MinTrackLength;
+  NSigmaCut = pars.NSigmaCut;
+  PVres = pars.PVres;
+  DeltaROF = pars.DeltaROF;
+  ZBins = pars.ZBins;
+  PhiBins = pars.PhiBins;
+  CellDeltaTanLambdaSigma = pars.CellDeltaTanLambdaSigma;
+}
+
+template <class T>
+GPUhd() T* getPtrFromRuler(int index, T* src, const int* ruler, const int stride = 1)
+{
+  return src + ruler[index] * stride;
+}
+
+template <class T>
+GPUhd() const T* getPtrFromRuler(int index, const T* src, const int* ruler, const int stride = 1)
+{
+  return src + ruler[index] * stride;
+}
+
 template <int NLayers = 7>
 class TimeFrameGPU : public TimeFrame
 {
@@ -45,33 +91,49 @@ class TimeFrameGPU : public TimeFrame
  public:
   TimeFrameGPU();
   ~TimeFrameGPU();
-  void loadToDevice(const int maxLayers);
 
   void checkBufferSizes();
   void initialise(const int iteration,
                   const TrackingParameters& trkParam,
                   const int maxLayers);
+  template <unsigned char isTracker = false>
+  void initialiseDevice(const TrackingParameters&);
   /// Getters
   float getDeviceMemory();
   Cluster* getDeviceClustersOnLayer(const int rofId, const int layerId) const;
+  unsigned char* getDeviceUsedClustersOnLayer(const int rofId, const int layerId);
+  int* getDeviceROframesClustersOnLayer(const int layerId) const { return mROframesClustersD[layerId].get(); }
   int getNClustersLayer(const int rofId, const int layerId) const;
   TimeFrameGPUConfig& getConfig() { return mConfig; }
+  gpu::Stream& getStream(const int iLayer) { return mStreamArray[iLayer]; }
+  std::vector<int>& getTrackletSizeHost() { return mTrackletSizeHost; }
+  std::vector<int>& getCellSizeHost() { return mCellSizeHost; }
 
   // Vertexer only
   int* getDeviceNTrackletsCluster(int rofId, int combId);
-  int* getDeviceIndexTableL0(const int rofId) { return mIndexTablesLayer0D.get() + rofId * (ZBins * PhiBins + 1); }
-  int* getDeviceIndexTableL2(const int rofId) { return mIndexTablesLayer2D.get() + rofId * (ZBins * PhiBins + 1); }
+  int* getDeviceIndexTables(const int layerId) { return mIndexTablesD[layerId].get(); }
+  int* getDeviceIndexTableAtRof(const int layerId, const int rofId) { return mIndexTablesD[layerId].get() + rofId * (ZBins * PhiBins + 1); }
   unsigned char* getDeviceUsedTracklets(const int rofId);
   Line* getDeviceLines(const int rofId);
+  Tracklet* getDeviceTrackletsVertexerOnly(const int rofId, const int layerId); // this method uses the cluster table for layer 1 for any layer. It is used for the vertexer only.
   Tracklet* getDeviceTracklets(const int rofId, const int layerId);
+  Tracklet* getDeviceTrackletsAll(const int layerId);
+  Cell* getDeviceCells(const int layerId);
+  int* getDeviceTrackletsLookupTable(const int rofId, const int layerId);
+  int* getDeviceCellsLookupTable(const int layerId);
   int* getDeviceNFoundLines(const int rofId);
   int* getDeviceExclusiveNFoundLines(const int rofId);
   int* getDeviceCUBBuffer(const size_t rofId);
+  int* getDeviceNFoundTracklets() const { return mDeviceFoundTracklets; };
+  int* getDeviceNFoundCells() const { return mDeviceFoundCells; };
   float* getDeviceXYCentroids(const int rofId);
   float* getDeviceZCentroids(const int rofId);
   int* getDeviceXHistograms(const int rofId);
   int* getDeviceYHistograms(const int rofId);
   int* getDeviceZHistograms(const int rofId);
+  gpu::StaticTrackingParameters<NLayers>* getDeviceTrackingParameters() const { return mDeviceTrackingParams; }
+  IndexTableUtils* getDeviceIndexTableUtils() const { return mDeviceIndexTableUtils; }
+
 #ifdef __HIPCC__
   hipcub::KeyValuePair<int, int>* getTmpVertexPositionBins(const int rofId);
 #else
@@ -82,15 +144,25 @@ class TimeFrameGPU : public TimeFrame
 
  private:
   TimeFrameGPUConfig mConfig;
-
+  std::array<gpu::Stream, NLayers + 1> mStreamArray;
+  std::vector<int> mTrackletSizeHost;
+  std::vector<int> mCellSizeHost;
   // Per-layer information, do not expand at runtime
   std::array<Vector<Cluster>, NLayers> mClustersD;
+  std::array<Vector<unsigned char>, NLayers> mUsedClustersD;
   std::array<Vector<TrackingFrameInfo>, NLayers> mTrackingFrameInfoD;
   std::array<Vector<int>, NLayers - 1> mIndexTablesD;
   std::array<Vector<int>, NLayers> mClusterExternalIndicesD;
-  std::array<Vector<int>, NLayers> mROframesClustersD;
   std::array<Vector<Tracklet>, NLayers - 1> mTrackletsD;
-  int* mCUBTmpBuffers; // don't know whether will be used by the tracker
+  std::array<Vector<int>, NLayers - 1> mTrackletsLookupTablesD;
+  std::array<Vector<Cell>, NLayers - 2> mCellsD;
+  std::array<Vector<int>, NLayers - 2> mCellsLookupTablesD;
+  std::array<Vector<int>, NLayers> mROframesClustersD; // layers x roframes
+  int* mCUBTmpBuffers;
+  int* mDeviceFoundTracklets;
+  int* mDeviceFoundCells;
+  gpu::StaticTrackingParameters<NLayers>* mDeviceTrackingParams;
+  IndexTableUtils* mDeviceIndexTableUtils;
 
   // Vertexer only
   Vector<Line> mLines;
@@ -115,139 +187,120 @@ class TimeFrameGPU : public TimeFrame
 template <int NLayers>
 inline Cluster* TimeFrameGPU<NLayers>::getDeviceClustersOnLayer(const int rofId, const int layerId) const
 {
-  if (rofId < 0 || rofId >= mNrof) {
-    LOG(error) << "Invalid rofId: " << rofId << "/" << mNrof << ", returning nullptr";
-    return nullptr;
-  }
-  return mClustersD[layerId].get() + mROframesClusters[layerId][rofId];
+  return getPtrFromRuler<Cluster>(rofId, mClustersD[layerId].get(), mROframesClusters[layerId].data());
+}
+
+template <int NLayers>
+inline unsigned char* TimeFrameGPU<NLayers>::getDeviceUsedClustersOnLayer(const int rofId, const int layerId)
+{
+  return getPtrFromRuler<unsigned char>(rofId, mUsedClustersD[layerId].get(), mROframesClusters[layerId].data());
 }
 
 template <int NLayers>
 inline int TimeFrameGPU<NLayers>::getNClustersLayer(const int rofId, const int layerId) const
 {
-  if (rofId < 0 || rofId >= mNrof) {
-    LOG(error) << "Invalid rofId: " << rofId << "/" << mNrof << ", returning 0 as value";
-    return 0;
-  }
   return static_cast<int>(mROframesClusters[layerId][rofId + 1] - mROframesClusters[layerId][rofId]);
 }
 
 template <int NLayers>
 inline int* TimeFrameGPU<NLayers>::getDeviceNTrackletsCluster(int rofId, int combId)
 {
-  if (rofId < 0 || rofId >= mNrof) {
-    LOG(error) << "Invalid rofId: " << rofId << "/" << mNrof << ", returning nullptr";
-    return nullptr;
-  }
-  return mNTrackletsPerClusterD[combId].get() + mROframesClusters[1][rofId];
+  return getPtrFromRuler<int>(rofId, mNTrackletsPerClusterD[combId].get(), mROframesClusters[1].data());
 }
 
 template <int NLayers>
 inline unsigned char* TimeFrameGPU<NLayers>::getDeviceUsedTracklets(const int rofId)
 {
-  if (rofId < 0 || rofId >= mNrof) {
-    LOG(error) << "Invalid rofId: " << rofId << "/" << mNrof << ", returning nullptr";
-    return nullptr;
-  }
-  return mUsedTracklets.get() + mROframesClusters[1][rofId] * mConfig.maxTrackletsPerCluster;
+  return getPtrFromRuler<unsigned char>(rofId, mUsedTracklets.get(), mROframesClusters[1].data(), mConfig.maxTrackletsPerCluster);
 }
 
 template <int NLayers>
 inline Line* TimeFrameGPU<NLayers>::getDeviceLines(const int rofId)
 {
-  if (rofId < 0 || rofId >= mNrof) {
-    LOG(error) << "Invalid rofId: " << rofId << "/" << mNrof << ", returning nullptr";
-    return nullptr;
-  }
-  return mLines.get() + mROframesClusters[1][rofId];
+  return getPtrFromRuler(rofId, mLines.get(), mROframesClusters[1].data());
+}
+
+template <int NLayers>
+inline Tracklet* TimeFrameGPU<NLayers>::getDeviceTrackletsVertexerOnly(const int rofId, const int layerId)
+{
+  return getPtrFromRuler(rofId, mTrackletsD[layerId].get(), mROframesClusters[1].data(), mConfig.maxTrackletsPerCluster);
 }
 
 template <int NLayers>
 inline Tracklet* TimeFrameGPU<NLayers>::getDeviceTracklets(const int rofId, const int layerId)
 {
-  if (rofId < 0 || rofId >= mNrof) {
-    LOG(error) << "Invalid rofId: " << rofId << "/" << mNrof << ", returning nullptr";
-    return nullptr;
-  }
-  return mTrackletsD[layerId].get() + mROframesClusters[1][rofId] * mConfig.maxTrackletsPerCluster;
+  return getPtrFromRuler(rofId, mTrackletsD[layerId].get(), mROframesClusters[layerId].data(), mConfig.maxTrackletsPerCluster);
+}
+
+template <int NLayers>
+inline Tracklet* TimeFrameGPU<NLayers>::getDeviceTrackletsAll(const int layerId)
+{
+  return mTrackletsD[layerId].get();
+}
+
+template <int NLayers>
+inline Cell* TimeFrameGPU<NLayers>::getDeviceCells(const int layerId)
+{
+  return mCellsD[layerId].get();
+}
+
+template <int NLayers>
+inline int* TimeFrameGPU<NLayers>::getDeviceTrackletsLookupTable(const int rofId, const int layerId)
+{
+  return getPtrFromRuler(rofId, mTrackletsLookupTablesD[layerId].get(), mROframesClusters[layerId].data());
+}
+
+template <int NLayers>
+inline int* TimeFrameGPU<NLayers>::getDeviceCellsLookupTable(const int layerId)
+{
+  return mCellsLookupTablesD[layerId].get();
 }
 
 template <int NLayers>
 inline int* TimeFrameGPU<NLayers>::getDeviceNFoundLines(const int rofId)
 {
-  if (rofId < 0 || rofId >= mNrof) {
-    LOG(error) << "Invalid rofId: " << rofId << "/" << mNrof << ", returning nullptr";
-    return nullptr;
-  }
-  return mNFoundLines.get() + mROframesClusters[1][rofId];
+  return getPtrFromRuler<int>(rofId, mNFoundLines.get(), mROframesClusters[1].data());
 }
 
 template <int NLayers>
 inline int* TimeFrameGPU<NLayers>::getDeviceExclusiveNFoundLines(const int rofId)
 {
-  if (rofId < 0 || rofId >= mNrof) {
-    LOG(error) << "Invalid rofId: " << rofId << "/" << mNrof << ", returning nullptr";
-    return nullptr;
-  }
-  return mNExclusiveFoundLines.get() + mROframesClusters[1][rofId];
+  return getPtrFromRuler<int>(rofId, mNExclusiveFoundLines.get(), mROframesClusters[1].data());
 }
 
 template <int NLayers>
 inline int* TimeFrameGPU<NLayers>::getDeviceCUBBuffer(const size_t rofId)
 {
-  if (rofId >= mNrof) {
-    LOG(error) << "Invalid rofId: " << rofId << "/" << mNrof << ", returning nullptr";
-  }
   return reinterpret_cast<int*>(reinterpret_cast<char*>(mCUBTmpBuffers) + (static_cast<size_t>(rofId * mConfig.tmpCUBBufferSize) & 0xFFFFFFFFFFFFF000));
 }
 
 template <int NLayers>
 inline float* TimeFrameGPU<NLayers>::getDeviceXYCentroids(const int rofId)
 {
-  if (rofId < 0 || rofId >= mNrof) {
-    LOG(error) << "Invalid rofId: " << rofId << "/" << mNrof << ", returning nullptr";
-    return nullptr;
-  }
   return mXYCentroids.get() + 2 * rofId * mConfig.maxCentroidsXYCapacity;
 }
 
 template <int NLayers>
 inline float* TimeFrameGPU<NLayers>::getDeviceZCentroids(const int rofId)
 {
-  if (rofId < 0 || rofId >= mNrof) {
-    LOG(error) << "Invalid rofId: " << rofId << "/" << mNrof << ", returning nullptr";
-    return nullptr;
-  }
   return mZCentroids.get() + rofId * mConfig.maxLinesCapacity;
 }
 
 template <int NLayers>
 inline int* TimeFrameGPU<NLayers>::getDeviceXHistograms(const int rofId)
 {
-  if (rofId < 0 || rofId >= mNrof) {
-    LOG(error) << "Invalid rofId: " << rofId << "/" << mNrof << ", returning nullptr";
-    return nullptr;
-  }
   return mXYZHistograms[0].get() + rofId * mConfig.histConf.nBinsXYZ[0];
 }
 
 template <int NLayers>
 inline int* TimeFrameGPU<NLayers>::getDeviceYHistograms(const int rofId)
 {
-  if (rofId < 0 || rofId >= mNrof) {
-    LOG(error) << "Invalid rofId: " << rofId << "/" << mNrof << ", returning nullptr";
-    return nullptr;
-  }
   return mXYZHistograms[1].get() + rofId * mConfig.histConf.nBinsXYZ[1];
 }
 
 template <int NLayers>
 inline int* TimeFrameGPU<NLayers>::getDeviceZHistograms(const int rofId)
 {
-  if (rofId < 0 || rofId >= mNrof) {
-    LOG(error) << "Invalid rofId: " << rofId << "/" << mNrof << ", returning nullptr";
-    return nullptr;
-  }
   return mXYZHistograms[2].get() + rofId * mConfig.histConf.nBinsXYZ[2];
 }
 
@@ -258,30 +311,18 @@ inline hipcub::KeyValuePair<int, int>* TimeFrameGPU<NLayers>::getTmpVertexPositi
 inline cub::KeyValuePair<int, int>* TimeFrameGPU<NLayers>::getTmpVertexPositionBins(const int rofId)
 #endif
 {
-  if (rofId < 0 || rofId >= mNrof) {
-    LOG(error) << "Invalid rofId: " << rofId << "/" << mNrof << ", returning nullptr";
-    return nullptr;
-  }
   return mTmpVertexPositionBins.get() + 3 * rofId;
 }
 
 template <int NLayers>
 inline float* TimeFrameGPU<NLayers>::getDeviceBeamPosition(const int rofId)
 {
-  if (rofId < 0 || rofId >= mNrof) {
-    LOG(error) << "Invalid rofId: " << rofId << "/" << mNrof << ", returning nullptr";
-    return nullptr;
-  }
   return mBeamPosition.get() + 2 * rofId;
 }
 
 template <int NLayers>
 inline Vertex* TimeFrameGPU<NLayers>::getDeviceVertices(const int rofId)
 {
-  if (rofId < 0 || rofId >= mNrof) {
-    LOG(error) << "Invalid rofId: " << rofId << "/" << mNrof << ", returning nullptr";
-    return nullptr;
-  }
   return mGPUVertices.get() + rofId * mConfig.maxVerticesCapacity;
 }
 
diff --git a/Detectors/ITSMFT/ITS/tracking/GPU/ITStrackingGPU/TrackerTraitsGPU.h b/Detectors/ITSMFT/ITS/tracking/GPU/ITStrackingGPU/TrackerTraitsGPU.h
index 01c5675115621..1f0e198836188 100644
--- a/Detectors/ITSMFT/ITS/tracking/GPU/ITStrackingGPU/TrackerTraitsGPU.h
+++ b/Detectors/ITSMFT/ITS/tracking/GPU/ITStrackingGPU/TrackerTraitsGPU.h
@@ -22,11 +22,6 @@ namespace o2
 {
 namespace its
 {
-namespace gpu
-{
-template <int NLayers>
-struct StaticTrackingParameters;
-}
 
 template <int NLayers = 7>
 class TrackerTraitsGPU : public TrackerTraits
@@ -36,13 +31,32 @@ class TrackerTraitsGPU : public TrackerTraits
   ~TrackerTraitsGPU() override = default;
 
   // void computeLayerCells() final;
+  void adoptTimeFrame(TimeFrame* tf) override;
+  void initialiseTimeFrame(const int iteration) override;
   void computeLayerTracklets(const int iteration) final;
+  void computeLayerCells(const int iteration) override;
+  void setBz(float) override;
+  void findCellsNeighbours(const int iteration) override;
+  void findRoads(const int iteration) override;
+  void findTracks() override;
+  void extendTracks(const int iteration) override;
   // void refitTracks(const std::vector<std::vector<TrackingFrameInfo>>& tf, std::vector<TrackITSExt>& tracks) override;
 
+  // TimeFrameGPU information forwarding
+  int getTFNumberOfClusters() const override;
+  int getTFNumberOfTracklets() const override;
+  int getTFNumberOfCells() const override;
+
  private:
-  gpu::TimeFrameGPU<NLayers> mTimeFrameGPU;
+  gpu::TimeFrameGPU<7>* mTimeFrameGPU;
   gpu::StaticTrackingParameters<NLayers>* mStaticTrkPars;
 };
+
+template <int NLayers>
+inline void TrackerTraitsGPU<NLayers>::adoptTimeFrame(TimeFrame* tf)
+{
+  mTimeFrameGPU = static_cast<gpu::TimeFrameGPU<NLayers>*>(tf);
+}
 } // namespace its
 } // namespace o2
 
diff --git a/Detectors/ITSMFT/ITS/tracking/GPU/ITStrackingGPU/Utils.h b/Detectors/ITSMFT/ITS/tracking/GPU/ITStrackingGPU/Utils.h
index 71ef7dad82691..acd28de061690 100644
--- a/Detectors/ITSMFT/ITS/tracking/GPU/ITStrackingGPU/Utils.h
+++ b/Detectors/ITSMFT/ITS/tracking/GPU/ITStrackingGPU/Utils.h
@@ -25,10 +25,8 @@ namespace its
 {
 namespace gpu
 {
-
 namespace utils
 {
-
 namespace host
 {
 
diff --git a/Detectors/ITSMFT/ITS/tracking/GPU/ITStrackingGPU/Vector.h b/Detectors/ITSMFT/ITS/tracking/GPU/ITStrackingGPU/Vector.h
index 64b7c2dababc2..38694156e0299 100644
--- a/Detectors/ITSMFT/ITS/tracking/GPU/ITStrackingGPU/Vector.h
+++ b/Detectors/ITSMFT/ITS/tracking/GPU/ITStrackingGPU/Vector.h
@@ -55,6 +55,8 @@ class Vector final
   void resize(const size_t);
   void reset(const size_t, const size_t = 0);
   void reset(const T* const, const size_t, const size_t = 0);
+
+  void resetWithInt(const size_t, const int value = 0);
   void copyIntoSizedVector(std::vector<T>&);
 
   GPUhd() T* get() const;
@@ -74,8 +76,8 @@ class Vector final
  private:
   GPUhd() Vector(const Vector&, const bool);
 
-  T* mArrayPointer = nullptr;
-  size_t* mDeviceSize = nullptr;
+  T* mArrayPtr = nullptr;
+  size_t* mDeviceSizePtr = nullptr;
   size_t mCapacity;
   bool mIsWeak;
 };
@@ -98,17 +100,17 @@ Vector<T>::Vector(const T* const source, const size_t size, const size_t initial
   if (size > 0) {
     try {
 
-      utils::host::gpuMalloc(reinterpret_cast<void**>(&mArrayPointer), size * sizeof(T));
-      utils::host::gpuMalloc(reinterpret_cast<void**>(&mDeviceSize), sizeof(size_t));
+      utils::host::gpuMalloc(reinterpret_cast<void**>(&mArrayPtr), size * sizeof(T));
+      utils::host::gpuMalloc(reinterpret_cast<void**>(&mDeviceSizePtr), sizeof(size_t));
 
       if (source != nullptr) {
 
-        utils::host::gpuMemcpyHostToDevice(mArrayPointer, source, size * sizeof(T));
-        utils::host::gpuMemcpyHostToDevice(mDeviceSize, &size, sizeof(size_t));
+        utils::host::gpuMemcpyHostToDevice(mArrayPtr, source, size * sizeof(T));
+        utils::host::gpuMemcpyHostToDevice(mDeviceSizePtr, &size, sizeof(size_t));
 
       } else {
 
-        utils::host::gpuMemcpyHostToDevice(mDeviceSize, &initialSize, sizeof(size_t));
+        utils::host::gpuMemcpyHostToDevice(mDeviceSizePtr, &initialSize, sizeof(size_t));
       }
 
     } catch (...) {
@@ -122,8 +124,8 @@ Vector<T>::Vector(const T* const source, const size_t size, const size_t initial
 
 template <typename T>
 GPUhd() Vector<T>::Vector(const Vector& other, const bool isWeak)
-  : mArrayPointer{other.mArrayPointer},
-    mDeviceSize{other.mDeviceSize},
+  : mArrayPtr{other.mArrayPtr},
+    mDeviceSizePtr{other.mDeviceSizePtr},
     mCapacity{other.mCapacity},
     mIsWeak{isWeak}
 {
@@ -146,13 +148,13 @@ GPUhd() Vector<T>::~Vector()
 
 template <typename T>
 GPUhd() Vector<T>::Vector(Vector<T>&& other)
-  : mArrayPointer{other.mArrayPointer},
-    mDeviceSize{other.mDeviceSize},
+  : mArrayPtr{other.mArrayPtr},
+    mDeviceSizePtr{other.mDeviceSizePtr},
     mCapacity{other.mCapacity},
     mIsWeak{other.mIsWeak}
 {
-  other.mArrayPointer = nullptr;
-  other.mDeviceSize = nullptr;
+  other.mArrayPtr = nullptr;
+  other.mDeviceSizePtr = nullptr;
 }
 
 template <typename T>
@@ -160,13 +162,13 @@ Vector<T>& Vector<T>::operator=(Vector<T>&& other)
 {
   destroy();
 
-  mArrayPointer = other.mArrayPointer;
-  mDeviceSize = other.mDeviceSize;
+  mArrayPtr = other.mArrayPtr;
+  mDeviceSizePtr = other.mDeviceSizePtr;
   mCapacity = other.mCapacity;
   mIsWeak = other.mIsWeak;
 
-  other.mArrayPointer = nullptr;
-  other.mDeviceSize = nullptr;
+  other.mArrayPtr = nullptr;
+  other.mDeviceSizePtr = nullptr;
 
   return *this;
 }
@@ -175,7 +177,7 @@ template <typename T>
 size_t Vector<T>::getSizeFromDevice() const
 {
   size_t size;
-  utils::host::gpuMemcpyDeviceToHost(&size, mDeviceSize, sizeof(size_t));
+  utils::host::gpuMemcpyDeviceToHost(&size, mDeviceSizePtr, sizeof(size_t));
 
   return size;
 }
@@ -183,7 +185,7 @@ size_t Vector<T>::getSizeFromDevice() const
 template <typename T>
 void Vector<T>::resize(const size_t size)
 {
-  utils::host::gpuMemcpyHostToDevice(mDeviceSize, &size, sizeof(size_t));
+  utils::host::gpuMemcpyHostToDevice(mDeviceSizePtr, &size, sizeof(size_t));
 }
 
 template <typename T>
@@ -196,49 +198,63 @@ template <typename T>
 void Vector<T>::reset(const T* const source, const size_t size, const size_t initialSize)
 {
   if (size > mCapacity) {
-    if (mArrayPointer != nullptr) {
-      utils::host::gpuFree(mArrayPointer);
+    if (mArrayPtr != nullptr) {
+      utils::host::gpuFree(mArrayPtr);
     }
-    utils::host::gpuMalloc(reinterpret_cast<void**>(&mArrayPointer), size * sizeof(T));
+    utils::host::gpuMalloc(reinterpret_cast<void**>(&mArrayPtr), size * sizeof(T));
     mCapacity = size;
   }
+  if (mDeviceSizePtr == nullptr) {
+    utils::host::gpuMalloc(reinterpret_cast<void**>(&mDeviceSizePtr), sizeof(size_t));
+  }
 
   if (source != nullptr) {
-    utils::host::gpuMemcpyHostToDevice(mArrayPointer, source, size * sizeof(T));
-    utils::host::gpuMemcpyHostToDevice(mDeviceSize, &size, sizeof(size_t));
-
+    utils::host::gpuMemcpyHostToDevice(mArrayPtr, source, size * sizeof(T));
+    utils::host::gpuMemcpyHostToDevice(mDeviceSizePtr, &size, sizeof(size_t));
   } else {
-    if (mDeviceSize == nullptr) {
-      utils::host::gpuMalloc(reinterpret_cast<void**>(&mDeviceSize), sizeof(size_t));
+    utils::host::gpuMemcpyHostToDevice(mDeviceSizePtr, &initialSize, sizeof(size_t));
+  }
+}
+
+template <typename T>
+void Vector<T>::resetWithInt(const size_t size, const int value)
+{
+  if (size > mCapacity) {
+    if (mArrayPtr != nullptr) {
+      utils::host::gpuFree(mArrayPtr);
     }
-    utils::host::gpuMemcpyHostToDevice(mDeviceSize, &initialSize, sizeof(size_t));
+    utils::host::gpuMalloc(reinterpret_cast<void**>(&mArrayPtr), size * sizeof(int));
+    mCapacity = size;
   }
+  if (mDeviceSizePtr == nullptr) {
+    utils::host::gpuMalloc(reinterpret_cast<void**>(&mDeviceSizePtr), sizeof(int));
+  }
+
+  utils::host::gpuMemset(mArrayPtr, value, size * sizeof(int));
+  utils::host::gpuMemcpyHostToDevice(mDeviceSizePtr, &size, sizeof(int));
 }
 
 template <typename T>
 void Vector<T>::copyIntoSizedVector(std::vector<T>& destinationVector)
 {
-  utils::host::gpuMemcpyDeviceToHost(destinationVector.data(), mArrayPointer, destinationVector.size() * sizeof(T));
+  utils::host::gpuMemcpyDeviceToHost(destinationVector.data(), mArrayPtr, destinationVector.size() * sizeof(T));
 }
 
 template <typename T>
 inline void Vector<T>::destroy()
 {
-  if (mArrayPointer != nullptr) {
-
-    utils::host::gpuFree(mArrayPointer);
+  if (mArrayPtr != nullptr) {
+    utils::host::gpuFree(mArrayPtr);
   }
-
-  if (mDeviceSize != nullptr) {
-
-    utils::host::gpuFree(mDeviceSize);
+  if (mDeviceSizePtr != nullptr) {
+    utils::host::gpuFree(mDeviceSizePtr);
   }
 }
 
 template <typename T>
 GPUhd() T* Vector<T>::get() const
 {
-  return mArrayPointer;
+  return mArrayPtr;
 }
 
 template <typename T>
@@ -256,14 +272,14 @@ GPUhd() Vector<T> Vector<T>::getWeakCopy() const
 template <typename T>
 GPUd() T& Vector<T>::operator[](const size_t index) const
 {
-  return mArrayPointer[index];
+  return mArrayPtr[index];
 }
 
 template <typename T>
 GPUh() T Vector<T>::getElementFromDevice(const size_t index) const
 {
   T element;
-  utils::host::gpuMemcpyDeviceToHost(&element, mArrayPointer + index, sizeof(T));
+  utils::host::gpuMemcpyDeviceToHost(&element, mArrayPtr + index, sizeof(T));
 
   return element;
 }
@@ -271,21 +287,21 @@ GPUh() T Vector<T>::getElementFromDevice(const size_t index) const
 template <typename T>
 GPUd() size_t Vector<T>::size() const
 {
-  return *mDeviceSize;
+  return *mDeviceSizePtr;
 }
 
 template <typename T>
 template <typename... Args>
 GPUd() void Vector<T>::emplace(const size_t index, Args&&... arguments)
 {
-  new (mArrayPointer + index) T(std::forward<Args>(arguments)...);
+  new (mArrayPtr + index) T(std::forward<Args>(arguments)...);
 }
 
 template <typename T>
 GPUhd() void Vector<T>::dump()
 {
-  printf("mArrayPointer = %p\nmDeviceSize   = %p\nmCapacity     = %d\nmIsWeak       = %s\n",
-         mArrayPointer, mDeviceSize, mCapacity, mIsWeak ? "true" : "false");
+  printf("mArrayPtr = %p\nmDeviceSize   = %p\nmCapacity     = %d\nmIsWeak       = %s\n",
+         mArrayPtr, mDeviceSizePtr, mCapacity, mIsWeak ? "true" : "false");
 }
 } // namespace gpu
 } // namespace its
diff --git a/Detectors/ITSMFT/ITS/tracking/GPU/ITStrackingGPU/VertexerTraitsGPU.h b/Detectors/ITSMFT/ITS/tracking/GPU/ITStrackingGPU/VertexerTraitsGPU.h
index a08a2bad61a41..c9d923191d31f 100644
--- a/Detectors/ITSMFT/ITS/tracking/GPU/ITStrackingGPU/VertexerTraitsGPU.h
+++ b/Detectors/ITSMFT/ITS/tracking/GPU/ITStrackingGPU/VertexerTraitsGPU.h
@@ -13,6 +13,7 @@
 /// \brief
 /// \author matteo.concas@cern.ch
 
+// #define VTX_DEBUG
 #ifndef ITSTRACKINGGPU_VERTEXERTRAITSGPU_H_
 #define ITSTRACKINGGPU_VERTEXERTRAITSGPU_H_
 
@@ -45,25 +46,16 @@ class VertexerTraitsGPU : public VertexerTraits
   void computeTracklets() override;
   void computeTrackletMatching() override;
   void computeVertices() override;
+  void updateVertexingParameters(const VertexingParameters& vrtPar) override;
   // void computeMCFiltering() override;
 
-  // GPU-specific getters
-  GPUd() static const int2 getBinsPhiRectWindow(const Cluster&, float maxdeltaphi);
+  void computeVerticesHist();
 
  protected:
   IndexTableUtils* mDeviceIndexTableUtils;
   gpu::TimeFrameGPU<7>* mTimeFrameGPU;
 };
 
-inline GPUd() const int2 VertexerTraitsGPU::getBinsPhiRectWindow(const Cluster& currentCluster, float phiCut)
-{
-  // This function returns the lowest PhiBin and the number of phi bins to be spanned, In the form int2{phiBinLow, PhiBinSpan}
-  const int phiBinMin{constants::its2::getPhiBinIndex(
-    math_utils::getNormalizedPhi(currentCluster.phi - phiCut))};
-  const int phiBinSpan{static_cast<int>(MATH_CEIL(phiCut * InversePhiBinSize))};
-  return int2{phiBinMin, phiBinSpan};
-}
-
 inline void VertexerTraitsGPU::adoptTimeFrame(TimeFrame* tf) { mTimeFrameGPU = static_cast<gpu::TimeFrameGPU<7>*>(tf); }
 
 } // namespace its
diff --git a/Detectors/ITSMFT/ITS/tracking/GPU/cuda/DeviceStoreGPU.cu b/Detectors/ITSMFT/ITS/tracking/GPU/cuda/DeviceStoreGPU.cu
deleted file mode 100644
index 3a3fd910226ea..0000000000000
--- a/Detectors/ITSMFT/ITS/tracking/GPU/cuda/DeviceStoreGPU.cu
+++ /dev/null
@@ -1,184 +0,0 @@
-// 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 DeviceStoreNV.cxx
-/// \brief
-///
-
-#include "ITStrackingGPU/DeviceStoreGPU.h"
-#include "ITStrackingGPU/Stream.h"
-#include <sstream>
-
-namespace
-{
-
-using namespace o2::its;
-
-__device__ void fillIndexTables(o2::its::gpu::DeviceStoreNV& primaryVertexContext, const int layerIndex)
-{
-
-  const int currentClusterIndex{static_cast<int>(blockDim.x * blockIdx.x + threadIdx.x)};
-  const int nextLayerClustersNum{static_cast<int>(primaryVertexContext.getClusters()[layerIndex + 1].size())};
-
-  if (currentClusterIndex < nextLayerClustersNum) {
-
-    const int currentBinIndex{
-      primaryVertexContext.getClusters()[layerIndex + 1][currentClusterIndex].indexTableBinIndex};
-    int previousBinIndex;
-
-    if (currentClusterIndex == 0) {
-
-      primaryVertexContext.getIndexTables()[layerIndex][0] = 0;
-      previousBinIndex = 0;
-
-    } else {
-
-      previousBinIndex = primaryVertexContext.getClusters()[layerIndex + 1][currentClusterIndex - 1].indexTableBinIndex;
-    }
-
-    if (currentBinIndex > previousBinIndex) {
-
-      for (int iBin{previousBinIndex + 1}; iBin <= currentBinIndex; ++iBin) {
-
-        primaryVertexContext.getIndexTables()[layerIndex][iBin] = currentClusterIndex;
-      }
-
-      previousBinIndex = currentBinIndex;
-    }
-
-    if (currentClusterIndex == nextLayerClustersNum - 1) {
-
-      for (int iBin{currentBinIndex + 1}; iBin <= o2::its::constants::its2::ZBins * o2::its::constants::its2::PhiBins;
-           iBin++) {
-
-        primaryVertexContext.getIndexTables()[layerIndex][iBin] = nextLayerClustersNum;
-      }
-    }
-  }
-}
-
-__device__ void fillTrackletsPerClusterTables(o2::its::gpu::DeviceStoreNV& primaryVertexContext, const int layerIndex)
-{
-  const int currentClusterIndex{static_cast<int>(blockDim.x * blockIdx.x + threadIdx.x)};
-  const int clustersSize{static_cast<int>(primaryVertexContext.getClusters()[layerIndex + 1].size())};
-
-  if (currentClusterIndex < clustersSize) {
-
-    primaryVertexContext.getTrackletsPerClusterTable()[layerIndex][currentClusterIndex] = 0;
-  }
-}
-
-__device__ void fillCellsPerClusterTables(o2::its::gpu::DeviceStoreNV& primaryVertexContext, const int layerIndex)
-{
-  const int totalThreadNum{static_cast<int>(primaryVertexContext.getClusters()[layerIndex + 1].size())};
-  const int trackletsSize{static_cast<int>(primaryVertexContext.getTracklets()[layerIndex + 1].capacity())};
-  const int trackletsPerThread{1 + (trackletsSize - 1) / totalThreadNum};
-  const int firstTrackletIndex{static_cast<int>(blockDim.x * blockIdx.x + threadIdx.x) * trackletsPerThread};
-
-  if (firstTrackletIndex < trackletsSize) {
-
-    const int trackletsToSet{min(trackletsSize, firstTrackletIndex + trackletsPerThread) - firstTrackletIndex};
-    memset(&primaryVertexContext.getCellsPerTrackletTable()[layerIndex][firstTrackletIndex], 0,
-           trackletsToSet * sizeof(int));
-  }
-}
-
-__global__ void fillDeviceStructures(o2::its::gpu::DeviceStoreNV& primaryVertexContext, const int layerIndex)
-{
-  fillIndexTables(primaryVertexContext, layerIndex);
-
-  if (layerIndex < o2::its::constants::its2::CellsPerRoad) {
-
-    fillTrackletsPerClusterTables(primaryVertexContext, layerIndex);
-  }
-
-  if (layerIndex < o2::its::constants::its2::CellsPerRoad - 1) {
-
-    fillCellsPerClusterTables(primaryVertexContext, layerIndex);
-  }
-}
-} // namespace
-
-namespace o2
-{
-namespace its
-{
-namespace gpu
-{
-
-DeviceStoreNV::DeviceStoreNV() = default;
-
-UniquePointer<DeviceStoreNV> DeviceStoreNV::initialise(const float3& primaryVertex,
-                                                       const std::array<std::vector<Cluster>, constants::its2::LayersNumber>& clusters,
-                                                       const std::array<std::vector<Tracklet>, constants::its2::TrackletsPerRoad>& tracklets,
-                                                       const std::array<std::vector<Cell>, constants::its2::CellsPerRoad>& cells,
-                                                       const std::array<std::vector<int>, constants::its2::CellsPerRoad - 1>& cellsLookupTable,
-                                                       const std::array<float, constants::its2::LayersNumber>& rmin,
-                                                       const std::array<float, constants::its2::LayersNumber>& rmax)
-{
-  mPrimaryVertex = UniquePointer<float3>{primaryVertex};
-
-  for (int iLayer{0}; iLayer < constants::its2::LayersNumber; ++iLayer) {
-    this->mRmin[iLayer] = rmin[iLayer];
-    this->mRmax[iLayer] = rmax[iLayer];
-
-    this->mClusters[iLayer] =
-      Vector<Cluster>{&clusters[iLayer][0], static_cast<int>(clusters[iLayer].size())};
-
-    if (iLayer < constants::its2::TrackletsPerRoad) {
-      this->mTracklets[iLayer].reset(tracklets[iLayer].capacity());
-    }
-
-    if (iLayer < constants::its2::CellsPerRoad) {
-
-      this->mTrackletsLookupTable[iLayer].reset(static_cast<int>(clusters[iLayer + 1].size()));
-      this->mTrackletsPerClusterTable[iLayer].reset(static_cast<int>(clusters[iLayer + 1].size()));
-      this->mCells[iLayer].reset(static_cast<int>(cells[iLayer].capacity()));
-    }
-
-    if (iLayer < constants::its2::CellsPerRoad - 1) {
-
-      this->mCellsLookupTable[iLayer].reset(static_cast<int>(cellsLookupTable[iLayer].size()));
-      this->mCellsPerTrackletTable[iLayer].reset(static_cast<int>(cellsLookupTable[iLayer].size()));
-    }
-  }
-
-  UniquePointer<DeviceStoreNV> gpuContextDevicePointer{*this};
-
-  std::array<Stream, constants::its2::LayersNumber> streamArray;
-
-  for (int iLayer{0}; iLayer < constants::its2::TrackletsPerRoad; ++iLayer) {
-
-    const int nextLayerClustersNum = static_cast<int>(clusters[iLayer + 1].size());
-
-    dim3 threadsPerBlock{utils::host::getBlockSize(nextLayerClustersNum)};
-    dim3 blocksGrid{utils::host::getBlocksGrid(threadsPerBlock, nextLayerClustersNum)};
-
-    fillDeviceStructures<<<blocksGrid, threadsPerBlock, 0, streamArray[iLayer].get()>>>(*gpuContextDevicePointer, iLayer);
-
-    cudaError_t error = cudaGetLastError();
-
-    if (error != cudaSuccess) {
-
-      std::ostringstream errorString{};
-      errorString << __FILE__ << ":" << __LINE__ << " CUDA API returned error [" << cudaGetErrorString(error)
-                  << "] (code " << error << ")" << std::endl;
-
-      throw std::runtime_error{errorString.str()};
-    }
-  }
-
-  return gpuContextDevicePointer;
-}
-
-} // namespace gpu
-} // namespace its
-} // namespace o2
diff --git a/Detectors/ITSMFT/ITS/tracking/GPU/cuda/Stream.cu b/Detectors/ITSMFT/ITS/tracking/GPU/cuda/Stream.cu
index ae2a7b2a0263a..f5a2706fe5b8e 100644
--- a/Detectors/ITSMFT/ITS/tracking/GPU/cuda/Stream.cu
+++ b/Detectors/ITSMFT/ITS/tracking/GPU/cuda/Stream.cu
@@ -21,7 +21,7 @@ namespace gpu
 
 Stream::Stream()
 {
-  discardResult(cudaStreamCreateWithFlags(&mStream, cudaStreamNonBlocking));
+  discardResult(cudaStreamCreateWithFlags(&mStream, cudaStreamDefault));
 }
 // usles
 Stream::~Stream()
diff --git a/Detectors/ITSMFT/ITS/tracking/GPU/cuda/TimeFrameGPU.cu b/Detectors/ITSMFT/ITS/tracking/GPU/cuda/TimeFrameGPU.cu
index 74668d006c4fe..bf680252011fa 100644
--- a/Detectors/ITSMFT/ITS/tracking/GPU/cuda/TimeFrameGPU.cu
+++ b/Detectors/ITSMFT/ITS/tracking/GPU/cuda/TimeFrameGPU.cu
@@ -13,6 +13,8 @@
 #include <fmt/format.h>
 #include <sstream>
 
+#include <thrust/fill.h>
+
 #include "ITStracking/Constants.h"
 
 #include "ITStrackingGPU/Utils.h"
@@ -25,7 +27,7 @@ namespace its
 using constants::MB;
 namespace gpu
 {
-
+using utils::host::checkGPUError;
 GPUh() void gpuThrowOnError()
 {
   cudaError_t error = cudaGetLastError();
@@ -40,48 +42,22 @@ GPUh() void gpuThrowOnError()
 template <int NLayers>
 TimeFrameGPU<NLayers>::TimeFrameGPU()
 {
-  getDeviceMemory(); // We don't check if we can store the data in the GPU for the moment, only log it.
-
-  for (int iLayer{0}; iLayer < NLayers; ++iLayer) { // Tracker and vertexer
-    mClustersD[iLayer] = Vector<Cluster>{mConfig.clustersPerLayerCapacity, mConfig.clustersPerLayerCapacity};
-    mTrackingFrameInfoD[iLayer] = Vector<TrackingFrameInfo>{mConfig.clustersPerLayerCapacity, mConfig.clustersPerLayerCapacity};
-    mClusterExternalIndicesD[iLayer] = Vector<int>{mConfig.clustersPerLayerCapacity, mConfig.clustersPerLayerCapacity};
-    mROframesClustersD[iLayer] = Vector<int>{mConfig.clustersPerROfCapacity, mConfig.clustersPerROfCapacity};
-    if (iLayer < NLayers - 1) {
-      mTrackletsD[iLayer] = Vector<Tracklet>{mConfig.trackletsCapacity,
-                                             mConfig.trackletsCapacity};
-    }
-  }
-
-  for (auto iComb{0}; iComb < 2; ++iComb) { // Vertexer only
-    mNTrackletsPerClusterD[iComb] = Vector<int>{mConfig.clustersPerLayerCapacity, mConfig.clustersPerLayerCapacity};
-  }
-  mIndexTablesLayer0D = Vector<int>{mConfig.nMaxROFs * (ZBins * PhiBins + 1), mConfig.nMaxROFs * (ZBins * PhiBins + 1)};
-  mIndexTablesLayer2D = Vector<int>{mConfig.nMaxROFs * (ZBins * PhiBins + 1), mConfig.nMaxROFs * (ZBins * PhiBins + 1)};
-  mLines = Vector<Line>{mConfig.trackletsCapacity, mConfig.trackletsCapacity};
-  mNFoundLines = Vector<int>{mConfig.clustersPerLayerCapacity, mConfig.clustersPerLayerCapacity};
-  mNExclusiveFoundLines = Vector<int>{mConfig.clustersPerLayerCapacity, mConfig.clustersPerLayerCapacity};
-  mUsedTracklets = Vector<unsigned char>{mConfig.trackletsCapacity, mConfig.trackletsCapacity};
-  discardResult(cudaMalloc(&mCUBTmpBuffers, mConfig.nMaxROFs * mConfig.tmpCUBBufferSize));
-  mXYCentroids = Vector<float>{2 * mConfig.nMaxROFs * mConfig.maxCentroidsXYCapacity, 2 * mConfig.nMaxROFs * mConfig.maxCentroidsXYCapacity};
-  mZCentroids = Vector<float>{mConfig.nMaxROFs * mConfig.maxLinesCapacity, mConfig.nMaxROFs * mConfig.maxLinesCapacity};
-  for (size_t i{0}; i < 3; ++i) {
-    mXYZHistograms[i] = Vector<int>{mConfig.nMaxROFs * mConfig.histConf.nBinsXYZ[i], mConfig.nMaxROFs * mConfig.histConf.nBinsXYZ[i]};
-  }
-  mTmpVertexPositionBins = Vector<cub::KeyValuePair<int, int>>{3 * mConfig.nMaxROFs, 3 * mConfig.nMaxROFs};
-  mBeamPosition = Vector<float>{2 * mConfig.nMaxROFs, 2 * mConfig.nMaxROFs};
-  mGPUVertices = Vector<Vertex>{mConfig.nMaxROFs * mConfig.maxVerticesCapacity, mConfig.nMaxROFs * mConfig.maxVerticesCapacity};
+  mIsGPU = true;
+  // getDeviceMemory(); To be updated
 }
 
 template <int NLayers>
 float TimeFrameGPU<NLayers>::getDeviceMemory()
 {
+  // We don't check if we can store the data in the GPU for the moment, only log it.
   float totalMemory{0};
   totalMemory += NLayers * mConfig.clustersPerLayerCapacity * sizeof(Cluster);
+  totalMemory += NLayers * mConfig.clustersPerLayerCapacity * sizeof(unsigned char);
   totalMemory += NLayers * mConfig.clustersPerLayerCapacity * sizeof(TrackingFrameInfo);
   totalMemory += NLayers * mConfig.clustersPerLayerCapacity * sizeof(int);
   totalMemory += NLayers * mConfig.clustersPerROfCapacity * sizeof(int);
   totalMemory += (NLayers - 1) * mConfig.trackletsCapacity * sizeof(Tracklet);
+  totalMemory += (NLayers - 1) * mConfig.nMaxROFs * (256 * 128 + 1) * sizeof(int);
   totalMemory += 2 * mConfig.clustersPerLayerCapacity * sizeof(int);
   totalMemory += 2 * mConfig.nMaxROFs * (ZBins * PhiBins + 1) * sizeof(int);
   totalMemory += mConfig.trackletsCapacity * sizeof(Line);
@@ -98,47 +74,98 @@ float TimeFrameGPU<NLayers>::getDeviceMemory()
   totalMemory += 2 * mConfig.nMaxROFs * sizeof(float);
   totalMemory += mConfig.nMaxROFs * mConfig.maxVerticesCapacity * sizeof(Vertex);
 
-  LOGP(debug, "Total requested memory for GPU: {:.2f} MB", totalMemory / MB);
-  LOGP(debug, "\t- Clusters: {:.2f} MB", NLayers * mConfig.clustersPerLayerCapacity * sizeof(Cluster) / MB);
-  LOGP(debug, "\t- Tracking frame info: {:.2f} MB", NLayers * mConfig.clustersPerLayerCapacity * sizeof(TrackingFrameInfo) / MB);
-  LOGP(debug, "\t- Cluster external indices: {:.2f} MB", NLayers * mConfig.clustersPerLayerCapacity * sizeof(int) / MB);
-  LOGP(debug, "\t- Clusters per ROf: {:.2f} MB", NLayers * mConfig.clustersPerROfCapacity * sizeof(int) / MB);
-  LOGP(debug, "\t- Tracklets: {:.2f} MB", (NLayers - 1) * mConfig.trackletsCapacity * sizeof(Tracklet) / MB);
-  LOGP(debug, "\t- N tracklets per cluster: {:.2f} MB", 2 * mConfig.clustersPerLayerCapacity * sizeof(int) / MB);
-  LOGP(debug, "\t- Index tables: {:.2f} MB", 2 * mConfig.nMaxROFs * (ZBins * PhiBins + 1) * sizeof(int) / MB);
-  LOGP(debug, "\t- Lines: {:.2f} MB", mConfig.trackletsCapacity * sizeof(Line) / MB);
-  LOGP(debug, "\t- N found lines: {:.2f} MB", mConfig.clustersPerLayerCapacity * sizeof(int) / MB);
-  LOGP(debug, "\t- N exclusive-scan found lines: {:.2f} MB", mConfig.clustersPerLayerCapacity * sizeof(int) / MB);
-  LOGP(debug, "\t- Used tracklets: {:.2f} MB", mConfig.trackletsCapacity * sizeof(unsigned char) / MB);
-  LOGP(debug, "\t- CUB tmp buffers: {:.2f} MB", mConfig.nMaxROFs * mConfig.tmpCUBBufferSize / MB);
-  LOGP(debug, "\t- XY centroids: {:.2f} MB", 2 * mConfig.nMaxROFs * mConfig.maxCentroidsXYCapacity * sizeof(float) / MB);
-  LOGP(debug, "\t- Z centroids: {:.2f} MB", mConfig.nMaxROFs * mConfig.maxLinesCapacity * sizeof(float) / MB);
-  LOGP(debug, "\t- XY histograms: {:.2f} MB", 2 * mConfig.nMaxROFs * mConfig.histConf.nBinsXYZ[0] * sizeof(int) / MB);
-  LOGP(debug, "\t- Z histograms: {:.2f} MB", mConfig.nMaxROFs * mConfig.histConf.nBinsXYZ[2] * sizeof(int) / MB);
-  LOGP(debug, "\t- TMP Vertex position bins: {:.2f} MB", 3 * mConfig.nMaxROFs * sizeof(cub::KeyValuePair<int, int>) / MB);
-  LOGP(debug, "\t- Beam positions: {:.2f} MB", 2 * mConfig.nMaxROFs * sizeof(float) / MB);
-  LOGP(debug, "\t- Vertices: {:.2f} MB", mConfig.nMaxROFs * mConfig.maxVerticesCapacity * sizeof(Vertex) / MB);
+  LOG(info) << fmt::format("Total requested memory for GPU: {:.2f} MB", totalMemory / MB);
+  LOG(info) << fmt::format("\t- Clusters: {:.2f} MB", NLayers * mConfig.clustersPerLayerCapacity * sizeof(Cluster) / MB);
+  LOG(info) << fmt::format("\t- Used clusters: {:.2f} MB", NLayers * mConfig.clustersPerLayerCapacity * sizeof(unsigned char) / MB);
+  LOG(info) << fmt::format("\t- Tracking frame info: {:.2f} MB", NLayers * mConfig.clustersPerLayerCapacity * sizeof(TrackingFrameInfo) / MB);
+  LOG(info) << fmt::format("\t- Cluster external indices: {:.2f} MB", NLayers * mConfig.clustersPerLayerCapacity * sizeof(int) / MB);
+  LOG(info) << fmt::format("\t- Clusters per ROf: {:.2f} MB", NLayers * mConfig.clustersPerROfCapacity * sizeof(int) / MB);
+  LOG(info) << fmt::format("\t- Tracklets: {:.2f} MB", (NLayers - 1) * mConfig.trackletsCapacity * sizeof(Tracklet) / MB);
+  LOG(info) << fmt::format("\t- Tracklet index tables: {:.2f} MB", (NLayers - 1) * mConfig.nMaxROFs * (256 * 128 + 1) * sizeof(int) / MB);
+  LOG(info) << fmt::format("\t- N tracklets per cluster: {:.2f} MB", 2 * mConfig.clustersPerLayerCapacity * sizeof(int) / MB);
+  LOG(info) << fmt::format("\t- Index tables: {:.2f} MB", 2 * mConfig.nMaxROFs * (ZBins * PhiBins + 1) * sizeof(int) / MB);
+  LOG(info) << fmt::format("\t- Lines: {:.2f} MB", mConfig.trackletsCapacity * sizeof(Line) / MB);
+  LOG(info) << fmt::format("\t- N found lines: {:.2f} MB", mConfig.clustersPerLayerCapacity * sizeof(int) / MB);
+  LOG(info) << fmt::format("\t- N exclusive-scan found lines: {:.2f} MB", mConfig.clustersPerLayerCapacity * sizeof(int) / MB);
+  LOG(info) << fmt::format("\t- Used tracklets: {:.2f} MB", mConfig.trackletsCapacity * sizeof(unsigned char) / MB);
+  LOG(info) << fmt::format("\t- CUB tmp buffers: {:.2f} MB", mConfig.nMaxROFs * mConfig.tmpCUBBufferSize / MB);
+  LOG(info) << fmt::format("\t- XY centroids: {:.2f} MB", 2 * mConfig.nMaxROFs * mConfig.maxCentroidsXYCapacity * sizeof(float) / MB);
+  LOG(info) << fmt::format("\t- Z centroids: {:.2f} MB", mConfig.nMaxROFs * mConfig.maxLinesCapacity * sizeof(float) / MB);
+  LOG(info) << fmt::format("\t- XY histograms: {:.2f} MB", 2 * mConfig.nMaxROFs * mConfig.histConf.nBinsXYZ[0] * sizeof(int) / MB);
+  LOG(info) << fmt::format("\t- Z histograms: {:.2f} MB", mConfig.nMaxROFs * mConfig.histConf.nBinsXYZ[2] * sizeof(int) / MB);
+  LOG(info) << fmt::format("\t- TMP Vertex position bins: {:.2f} MB", 3 * mConfig.nMaxROFs * sizeof(cub::KeyValuePair<int, int>) / MB);
+  LOG(info) << fmt::format("\t- Beam positions: {:.2f} MB", 2 * mConfig.nMaxROFs * sizeof(float) / MB);
+  LOG(info) << fmt::format("\t- Vertices: {:.2f} MB", mConfig.nMaxROFs * mConfig.maxVerticesCapacity * sizeof(Vertex) / MB);
 
   return totalMemory;
 }
 
 template <int NLayers>
-void TimeFrameGPU<NLayers>::loadToDevice(const int maxLayers)
+template <unsigned char isTracker>
+void TimeFrameGPU<NLayers>::initialiseDevice(const TrackingParameters& trkParam)
 {
-  for (int iLayer{0}; iLayer < maxLayers; ++iLayer) {
+  mTrackletSizeHost.resize(NLayers - 1, 0);
+  mCellSizeHost.resize(NLayers - 2, 0);
+  for (int iLayer{0}; iLayer < NLayers - 1; ++iLayer) { // Tracker and vertexer
+    mTrackletsD[iLayer] = Vector<Tracklet>{mConfig.trackletsCapacity, mConfig.trackletsCapacity};
+    auto thrustTrackletsBegin = thrust::device_ptr<Tracklet>(mTrackletsD[iLayer].get());
+    auto thrustTrackletsEnd = thrustTrackletsBegin + mConfig.trackletsCapacity;
+    thrust::fill(thrustTrackletsBegin, thrustTrackletsEnd, Tracklet{});
+    mTrackletsLookupTablesD[iLayer].resetWithInt(mClusters[iLayer].size());
+    if (iLayer < NLayers - 2) {
+      mCellsD[iLayer] = Vector<Cell>{mConfig.validatedTrackletsCapacity, mConfig.validatedTrackletsCapacity};
+      mCellsLookupTablesD[iLayer] = Vector<int>{mConfig.cellsLUTsize, mConfig.cellsLUTsize};
+      mCellsLookupTablesD[iLayer].resetWithInt(mConfig.cellsLUTsize);
+    }
+  }
+
+  for (auto iComb{0}; iComb < 2; ++iComb) { // Vertexer only
+    mNTrackletsPerClusterD[iComb] = Vector<int>{mConfig.clustersPerLayerCapacity, mConfig.clustersPerLayerCapacity};
+  }
+  mLines = Vector<Line>{mConfig.trackletsCapacity, mConfig.trackletsCapacity};
+  mNFoundLines = Vector<int>{mConfig.clustersPerLayerCapacity, mConfig.clustersPerLayerCapacity};
+  mNFoundLines.resetWithInt(mConfig.clustersPerLayerCapacity);
+  mNExclusiveFoundLines = Vector<int>{mConfig.clustersPerLayerCapacity, mConfig.clustersPerLayerCapacity};
+  mNExclusiveFoundLines.resetWithInt(mConfig.clustersPerLayerCapacity);
+  mUsedTracklets = Vector<unsigned char>{mConfig.trackletsCapacity, mConfig.trackletsCapacity};
+  discardResult(cudaMalloc(&mCUBTmpBuffers, mConfig.nMaxROFs * mConfig.tmpCUBBufferSize));
+  discardResult(cudaMalloc(&mDeviceFoundTracklets, (NLayers - 1) * sizeof(int)));
+  discardResult(cudaMemset(mDeviceFoundTracklets, 0, (NLayers - 1) * sizeof(int)));
+  discardResult(cudaMalloc(&mDeviceFoundCells, (NLayers - 2) * sizeof(int)));
+  discardResult(cudaMemset(mDeviceFoundCells, 0, (NLayers - 2) * sizeof(int)));
+  mXYCentroids = Vector<float>{2 * mConfig.nMaxROFs * mConfig.maxCentroidsXYCapacity, 2 * mConfig.nMaxROFs * mConfig.maxCentroidsXYCapacity};
+  mZCentroids = Vector<float>{mConfig.nMaxROFs * mConfig.maxLinesCapacity, mConfig.nMaxROFs * mConfig.maxLinesCapacity};
+  for (size_t i{0}; i < 3; ++i) {
+    mXYZHistograms[i] = Vector<int>{mConfig.nMaxROFs * mConfig.histConf.nBinsXYZ[i], mConfig.nMaxROFs * mConfig.histConf.nBinsXYZ[i]};
+  }
+  mTmpVertexPositionBins = Vector<cub::KeyValuePair<int, int>>{3 * mConfig.nMaxROFs, 3 * mConfig.nMaxROFs};
+  mBeamPosition = Vector<float>{2 * mConfig.nMaxROFs, 2 * mConfig.nMaxROFs};
+  mGPUVertices = Vector<Vertex>{mConfig.nMaxROFs * mConfig.maxVerticesCapacity, mConfig.nMaxROFs * mConfig.maxVerticesCapacity};
+  //////////////////////////////////////////////////////////////////////////////
+  constexpr int layers = isTracker ? NLayers : 3;
+  for (int iLayer{0}; iLayer < layers; ++iLayer) {
     mClustersD[iLayer].reset(mClusters[iLayer].data(), static_cast<int>(mClusters[iLayer].size()));
-    mROframesClustersD[iLayer].reset(mROframesClusters[iLayer].data(), static_cast<int>(mROframesClusters[iLayer].size()));
   }
-  if (maxLayers == NLayers) {
+  if constexpr (isTracker) {
+    StaticTrackingParameters<NLayers> pars;
+    pars.set(trkParam);
+    checkGPUError(cudaMalloc(reinterpret_cast<void**>(&mDeviceTrackingParams), sizeof(gpu::StaticTrackingParameters<NLayers>)), __FILE__, __LINE__);
+    checkGPUError(cudaMalloc(reinterpret_cast<void**>(&mDeviceIndexTableUtils), sizeof(IndexTableUtils)), __FILE__, __LINE__);
+    checkGPUError(cudaMemcpy(mDeviceTrackingParams, &pars, sizeof(gpu::StaticTrackingParameters<NLayers>), cudaMemcpyHostToDevice), __FILE__, __LINE__);
+    checkGPUError(cudaMemcpy(mDeviceIndexTableUtils, &mIndexTableUtils, sizeof(IndexTableUtils), cudaMemcpyHostToDevice), __FILE__, __LINE__);
     // Tracker-only: we don't need to copy data in vertexer
-    for (int iLayer{0}; iLayer < maxLayers; ++iLayer) {
+    for (int iLayer{0}; iLayer < NLayers; ++iLayer) {
+      mUsedClustersD[iLayer].reset(mUsedClusters[iLayer].data(), static_cast<int>(mUsedClusters[iLayer].size()));
       mTrackingFrameInfoD[iLayer].reset(mTrackingFrameInfo[iLayer].data(), static_cast<int>(mTrackingFrameInfo[iLayer].size()));
       mClusterExternalIndicesD[iLayer].reset(mClusterExternalIndices[iLayer].data(), static_cast<int>(mClusterExternalIndices[iLayer].size()));
+      mROframesClustersD[iLayer].reset(mROframesClusters[iLayer].data(), static_cast<int>(mROframesClusters[iLayer].size()));
+      mIndexTablesD[iLayer].reset(mIndexTables[iLayer].data(), static_cast<int>(mIndexTables[iLayer].size()));
     }
   } else {
-    mIndexTablesLayer0D.reset(getIndexTableWhole(0).data(), static_cast<int>(getIndexTableWhole(0).size()));
-    mIndexTablesLayer2D.reset(getIndexTableWhole(2).data(), static_cast<int>(getIndexTableWhole(2).size()));
+    mIndexTablesD[0].reset(getIndexTableWhole(0).data(), static_cast<int>(getIndexTableWhole(0).size()));
+    mIndexTablesD[2].reset(getIndexTableWhole(2).data(), static_cast<int>(getIndexTableWhole(2).size()));
   }
+
   gpuThrowOnError();
 }
 
@@ -149,11 +176,22 @@ void TimeFrameGPU<NLayers>::initialise(const int iteration,
 {
   o2::its::TimeFrame::initialise(iteration, trkParam, maxLayers);
   checkBufferSizes();
-  loadToDevice(maxLayers);
+  if (maxLayers < NLayers) {
+    initialiseDevice<false>(trkParam); // vertexer
+  } else {
+    initialiseDevice<true>(trkParam); // tracker
+  }
 }
 
 template <int NLayers>
-TimeFrameGPU<NLayers>::~TimeFrameGPU() = default;
+TimeFrameGPU<NLayers>::~TimeFrameGPU()
+{
+  discardResult(cudaFree(mCUBTmpBuffers));
+  discardResult(cudaFree(mDeviceFoundTracklets));
+  discardResult(cudaFree(mDeviceTrackingParams));
+  discardResult(cudaFree(mDeviceIndexTableUtils));
+  discardResult(cudaFree(mDeviceFoundCells));
+}
 
 template <int NLayers>
 void TimeFrameGPU<NLayers>::checkBufferSizes()
diff --git a/Detectors/ITSMFT/ITS/tracking/GPU/cuda/TrackerTraitsGPU.cu b/Detectors/ITSMFT/ITS/tracking/GPU/cuda/TrackerTraitsGPU.cu
index a04e51dfcd69f..8132e96399271 100644
--- a/Detectors/ITSMFT/ITS/tracking/GPU/cuda/TrackerTraitsGPU.cu
+++ b/Detectors/ITSMFT/ITS/tracking/GPU/cuda/TrackerTraitsGPU.cu
@@ -13,23 +13,24 @@
 #include <array>
 #include <sstream>
 #include <iostream>
+#include <unistd.h>
 
-// #ifndef GPUCA_GPUCODE_GENRTC
-// #include <cooperative_groups.h>
-// #include "cub/cub.cuh"
-// #endif
+#include <thrust/execution_policy.h>
+#include <thrust/device_ptr.h>
+#include <thrust/sort.h>
+#include <thrust/reduce.h>
+#include <thrust/functional.h>
+#include <thrust/unique.h>
 
 #include "ITStracking/Constants.h"
 #include "ITStracking/Configuration.h"
 #include "ITStracking/IndexTableUtils.h"
 #include "ITStracking/MathUtils.h"
-#include "ITStracking/TimeFrame.h"
 
-#include "ITStrackingGPU/Context.h"
-#include "ITStrackingGPU/Stream.h"
-#include "ITStrackingGPU/Vector.h"
 #include "ITStrackingGPU/TrackerTraitsGPU.h"
 
+#include "GPUCommonLogger.h"
+#include "GPUCommonAlgorithmThrust.h"
 namespace o2
 {
 namespace its
@@ -37,7 +38,11 @@ namespace its
 using gpu::utils::host::checkGPUError;
 using namespace constants::its2;
 
+namespace gpu
+{
+
 GPUd() const int4 getBinsRect(const Cluster& currentCluster, const int layerIndex,
+                              const o2::its::IndexTableUtils& utils,
                               const float z1, const float z2, float maxdeltaz, float maxdeltaphi)
 {
   const float zRangeMin = o2::gpu::GPUCommonMath::Min(z1, z2) - maxdeltaz;
@@ -52,483 +57,408 @@ GPUd() const int4 getBinsRect(const Cluster& currentCluster, const int layerInde
   }
 
   return int4{o2::gpu::GPUCommonMath::Max(0, getZBinIndex(layerIndex + 1, zRangeMin)),
-              getPhiBinIndex(phiRangeMin),
+              utils.getPhiBinIndex(math_utils::getNormalizedPhi(phiRangeMin)),
               o2::gpu::GPUCommonMath::Min(ZBins - 1, getZBinIndex(layerIndex + 1, zRangeMax)),
-              getPhiBinIndex(phiRangeMax)};
+              utils.getPhiBinIndex(math_utils::getNormalizedPhi(phiRangeMax))};
 }
 
-// template <int NLayers>
-// void TrackerTraitsGPU<NLayers>::loadToDevice()
-// {
-//   mTimeFrameGPU.loadToDevice();
-// }
-
-namespace gpu
+GPUhd() float Sq(float q)
 {
+  return q * q;
+}
 
-template <int NLayers>
-struct StaticTrackingParameters {
-  // StaticTrackingParameters<NLayers>& operator=(const StaticTrackingParameters<NLayers>& t);
-  // int CellMinimumLevel();
-  /// General parameters
-  int ClusterSharing = 0;
-  int MinTrackLength = NLayers;
-  /// Trackleting cuts
-  float TrackletMaxDeltaPhi = 0.3f;
-  float TrackletMaxDeltaZ[NLayers - 1] = {0.1f, 0.1f, 0.3f, 0.3f, 0.3f, 0.3f};
-  /// Cell finding cuts
-  // float CellMaxDeltaTanLambda = 0.025f;
-  // float CellMaxDCA[NLayers - 2] = {0.05f, 0.04f, 0.05f, 0.2f, 0.4f};
-  // float CellMaxDeltaPhi = 0.14f;
-  // float CellMaxDeltaZ[NLayers - 2] = {0.2f, 0.4f, 0.5f, 0.6f, 3.0f};
-  // /// Neighbour finding cuts
-  // float NeighbourMaxDeltaCurvature[NLayers - 3] = {0.008f, 0.0025f, 0.003f, 0.0035f};
-  // float NeighbourMaxDeltaN[NLayers - 3] = {0.002f, 0.0090f, 0.002f, 0.005f};
+template <typename T>
+struct trackletSortEmptyFunctor : public thrust::binary_function<T, T, bool> {
+  GPUhd() bool operator()(const T& lhs, const T& rhs) const
+  {
+    return lhs.firstClusterIndex > rhs.firstClusterIndex;
+  }
 };
 
-template struct gpu::StaticTrackingParameters<7>;
-__constant__ StaticTrackingParameters<7> kTrkPar;
-
-// GPUd() void computeLayerTracklets(DeviceStoreNV& devStore, const int layerIndex,
-//                                   Vector<Tracklet>& trackletsVector)
-// {
-//   const int currentClusterIndex = static_cast<int>(blockDim.x * blockIdx.x + threadIdx.x);
-//   int clusterTrackletsNum = 0;
-
-//   if (currentClusterIndex < devStore.getClusters()[layerIndex].size()) {
-
-//     Vector<Cluster> nextLayerClusters{devStore.getClusters()[layerIndex + 1].getWeakCopy()};
-//     const Cluster currentCluster{devStore.getClusters()[layerIndex][currentClusterIndex]};
-
-//     /*if (mUsedClustersTable[currentCluster.clusterId] != constants::its::UnusedIndex) {
-
-//      continue;
-//      }*/
-
-//     const float tanLambda{(currentCluster.zCoordinate - devStore.getPrimaryVertex().z) / currentCluster.radius};
-//     const float zAtRmin{tanLambda * (devStore.getRmin(layerIndex + 1) - currentCluster.radius) + currentCluster.zCoordinate};
-//     const float zAtRmax{tanLambda * (devStore.getRmax(layerIndex + 1) - currentCluster.radius) + currentCluster.zCoordinate};
-
-//     const int4 selectedBinsRect{getBinsRect(currentCluster, layerIndex, zAtRmin, zAtRmax,
-//                                             kTrkPar.TrackletMaxDeltaZ[layerIndex], kTrkPar.TrackletMaxDeltaPhi)};
-
-//     if (selectedBinsRect.x != 0 || selectedBinsRect.y != 0 || selectedBinsRect.z != 0 || selectedBinsRect.w != 0) {
-
-//       const int nextLayerClustersNum{static_cast<int>(nextLayerClusters.size())};
-//       int phiBinsNum{selectedBinsRect.w - selectedBinsRect.y + 1};
-
-//       if (phiBinsNum < 0) {
-
-//         phiBinsNum += constants::its2::PhiBins;
-//       }
-
-//       for (int iPhiBin{selectedBinsRect.y}, iPhiCount{0}; iPhiCount < phiBinsNum;
-//            iPhiBin = ++iPhiBin == constants::its2::PhiBins ? 0 : iPhiBin, iPhiCount++) {
-
-//         const int firstBinIndex{constants::its2::getBinIndex(selectedBinsRect.x, iPhiBin)};
-//         const int firstRowClusterIndex = devStore.getIndexTables()[layerIndex][firstBinIndex];
-//         const int maxRowClusterIndex = devStore.getIndexTables()[layerIndex][{firstBinIndex + selectedBinsRect.z - selectedBinsRect.x + 1}];
-
-//         for (int iNextLayerCluster{firstRowClusterIndex};
-//              iNextLayerCluster <= maxRowClusterIndex && iNextLayerCluster < nextLayerClustersNum; ++iNextLayerCluster) {
-
-//           const Cluster& nextCluster{nextLayerClusters[iNextLayerCluster]};
-
-//           const float deltaZ{o2::gpu::GPUCommonMath::Abs(
-//             tanLambda * (nextCluster.radius - currentCluster.radius) + currentCluster.zCoordinate - nextCluster.zCoordinate)};
-//           const float deltaPhi{o2::gpu::GPUCommonMath::Abs(currentCluster.phi - nextCluster.phi)};
-
-//           if (deltaZ < kTrkPar.TrackletMaxDeltaZ[layerIndex] && (deltaPhi < kTrkPar.TrackletMaxDeltaPhi || o2::gpu::GPUCommonMath::Abs(deltaPhi - constants::math::TwoPi) < kTrkPar.TrackletMaxDeltaPhi)) {
-
-//             cooperative_groups::coalesced_group threadGroup = cooperative_groups::coalesced_threads();
-//             int currentIndex{};
-
-//             if (threadGroup.thread_rank() == 0) {
-
-//               currentIndex = trackletsVector.extend(threadGroup.size());
-//             }
-
-//             currentIndex = threadGroup.shfl(currentIndex, 0) + threadGroup.thread_rank();
-
-//             trackletsVector.emplace(currentIndex, currentClusterIndex, iNextLayerCluster, currentCluster, nextCluster);
-//             ++clusterTrackletsNum;
-//           }
-//         }
-//       }
-
-//       if (layerIndex > 0) {
-
-//         devStore.getTrackletsPerClusterTable()[layerIndex - 1][currentClusterIndex] = clusterTrackletsNum;
-//       }
-//     }
-//   }
-// }
-
-// GPUd() void computeLayerCells(DeviceStoreNV& devStore, const int layerIndex,
-//                               Vector<Cell>& cellsVector)
-// {
-//   const int currentTrackletIndex = static_cast<int>(blockDim.x * blockIdx.x + threadIdx.x);
-//   const float3& primaryVertex = devStore.getPrimaryVertex();
-//   int trackletCellsNum = 0;
-
-//   if (currentTrackletIndex < devStore.getTracklets()[layerIndex].size()) {
-
-//     const Tracklet& currentTracklet{devStore.getTracklets()[layerIndex][currentTrackletIndex]};
-//     const int nextLayerClusterIndex{currentTracklet.secondClusterIndex};
-//     const int nextLayerFirstTrackletIndex{
-//       devStore.getTrackletsLookupTable()[layerIndex][nextLayerClusterIndex]};
-//     const int nextLayerTrackletsNum{static_cast<int>(devStore.getTracklets()[layerIndex + 1].size())};
-
-//     if (devStore.getTracklets()[layerIndex + 1][nextLayerFirstTrackletIndex].firstClusterIndex == nextLayerClusterIndex) {
-
-//       const Cluster& firstCellCluster{
-//         devStore.getClusters()[layerIndex][currentTracklet.firstClusterIndex]};
-//       const Cluster& secondCellCluster{
-//         devStore.getClusters()[layerIndex + 1][currentTracklet.secondClusterIndex]};
-//       const float firstCellClusterQuadraticRCoordinate{firstCellCluster.radius * firstCellCluster.radius};
-//       const float secondCellClusterQuadraticRCoordinate{secondCellCluster.radius * secondCellCluster.radius};
-//       const float3 firstDeltaVector{secondCellCluster.xCoordinate - firstCellCluster.xCoordinate,
-//                                     secondCellCluster.yCoordinate - firstCellCluster.yCoordinate, secondCellClusterQuadraticRCoordinate - firstCellClusterQuadraticRCoordinate};
-
-//       for (int iNextLayerTracklet{nextLayerFirstTrackletIndex};
-//            iNextLayerTracklet < nextLayerTrackletsNum && devStore.getTracklets()[layerIndex + 1][iNextLayerTracklet].firstClusterIndex == nextLayerClusterIndex; ++iNextLayerTracklet) {
-
-//         const Tracklet& nextTracklet{devStore.getTracklets()[layerIndex + 1][iNextLayerTracklet]};
-//         const float deltaTanLambda{o2::gpu::GPUCommonMath::Abs(currentTracklet.tanLambda - nextTracklet.tanLambda)};
-//         const float deltaPhi{o2::gpu::GPUCommonMath::Abs(currentTracklet.phi - nextTracklet.phi)};
-
-//         if (deltaTanLambda < kTrkPar.CellMaxDeltaTanLambda && (deltaPhi < kTrkPar.CellMaxDeltaPhi || o2::gpu::GPUCommonMath::Abs(deltaPhi - constants::math::TwoPi) < kTrkPar.CellMaxDeltaPhi)) {
-
-//           const float averageTanLambda{0.5f * (currentTracklet.tanLambda + nextTracklet.tanLambda)};
-//           const float directionZIntersection{-averageTanLambda * firstCellCluster.radius + firstCellCluster.zCoordinate};
-//           const float deltaZ{o2::gpu::GPUCommonMath::Abs(directionZIntersection - primaryVertex.z)};
-
-//           if (deltaZ < kTrkPar.CellMaxDeltaZ[layerIndex]) {
-
-//             const Cluster& thirdCellCluster{
-//               devStore.getClusters()[layerIndex + 2][nextTracklet.secondClusterIndex]};
-
-//             const float thirdCellClusterQuadraticRCoordinate{thirdCellCluster.radius * thirdCellCluster.radius};
-
-//             const float3 secondDeltaVector{thirdCellCluster.xCoordinate - firstCellCluster.xCoordinate,
-//                                            thirdCellCluster.yCoordinate - firstCellCluster.yCoordinate, thirdCellClusterQuadraticRCoordinate - firstCellClusterQuadraticRCoordinate};
-
-//             float3 cellPlaneNormalVector{math_utils::crossProduct(firstDeltaVector, secondDeltaVector)};
-
-//             const float vectorNorm{o2::gpu::GPUCommonMath::Sqrt(
-//               cellPlaneNormalVector.x * cellPlaneNormalVector.x + cellPlaneNormalVector.y * cellPlaneNormalVector.y + cellPlaneNormalVector.z * cellPlaneNormalVector.z)};
-
-//             if (!(vectorNorm < constants::math::FloatMinThreshold || o2::gpu::GPUCommonMath::Abs(cellPlaneNormalVector.z) < constants::math::FloatMinThreshold)) {
-
-//               const float inverseVectorNorm{1.0f / vectorNorm};
-//               const float3 normalizedPlaneVector{cellPlaneNormalVector.x * inverseVectorNorm, cellPlaneNormalVector.y * inverseVectorNorm, cellPlaneNormalVector.z * inverseVectorNorm};
-//               const float planeDistance{-normalizedPlaneVector.x * (secondCellCluster.xCoordinate - primaryVertex.x) - (normalizedPlaneVector.y * secondCellCluster.yCoordinate - primaryVertex.y) - normalizedPlaneVector.z * secondCellClusterQuadraticRCoordinate};
-//               const float normalizedPlaneVectorQuadraticZCoordinate{normalizedPlaneVector.z * normalizedPlaneVector.z};
-//               const float cellTrajectoryRadius{o2::gpu::GPUCommonMath::Sqrt(
-//                 (1.0f - normalizedPlaneVectorQuadraticZCoordinate - 4.0f * planeDistance * normalizedPlaneVector.z) / (4.0f * normalizedPlaneVectorQuadraticZCoordinate))};
-//               const float2 circleCenter{-0.5f * normalizedPlaneVector.x / normalizedPlaneVector.z, -0.5f * normalizedPlaneVector.y / normalizedPlaneVector.z};
-//               const float distanceOfClosestApproach{o2::gpu::GPUCommonMath::Abs(
-//                 cellTrajectoryRadius - o2::gpu::GPUCommonMath::Sqrt(circleCenter.x * circleCenter.x + circleCenter.y * circleCenter.y))};
-
-//               if (distanceOfClosestApproach <= kTrkPar.CellMaxDCA[layerIndex]) {
-
-//                 cooperative_groups::coalesced_group threadGroup = cooperative_groups::coalesced_threads();
-//                 int currentIndex{};
-
-//                 if (threadGroup.thread_rank() == 0) {
-
-//                   currentIndex = cellsVector.extend(threadGroup.size());
-//                 }
-
-//                 currentIndex = threadGroup.shfl(currentIndex, 0) + threadGroup.thread_rank();
-
-//                 cellsVector.emplace(currentIndex, currentTracklet.firstClusterIndex,
-//                                     nextTracklet.firstClusterIndex, nextTracklet.secondClusterIndex, currentTrackletIndex,
-//                                     iNextLayerTracklet, averageTanLambda);
-//                 ++trackletCellsNum;
-//               }
-//             }
-//           }
-//         }
-//       }
-
-//       if (layerIndex > 0) {
-
-//         devStore.getCellsPerTrackletTable()[layerIndex - 1][currentTrackletIndex] = trackletCellsNum;
-//       }
-//     }
-//   }
-// }
-
-// GPUg() void layerTrackletsKernel(DeviceStoreNV& devStore, const int layerIndex,
-//                                  Vector<Tracklet> trackletsVector)
-// {
-//   computeLayerTracklets(devStore, layerIndex, trackletsVector);
-// }
-
-// GPUg() void sortTrackletsKernel(DeviceStoreNV& devStore, const int layerIndex,
-//                                 Vector<Tracklet> tempTrackletArray)
-// {
-//   const int currentTrackletIndex{static_cast<int>(blockDim.x * blockIdx.x + threadIdx.x)};
-
-//   if (currentTrackletIndex < tempTrackletArray.size()) {
-
-//     const int firstClusterIndex = tempTrackletArray[currentTrackletIndex].firstClusterIndex;
-//     const int offset = atomicAdd(&devStore.getTrackletsPerClusterTable()[layerIndex - 1][firstClusterIndex],
-//                                  -1) -
-//                        1;
-//     const int startIndex = devStore.getTrackletsLookupTable()[layerIndex - 1][firstClusterIndex];
-
-//     memcpy(&devStore.getTracklets()[layerIndex][startIndex + offset],
-//            &tempTrackletArray[currentTrackletIndex], sizeof(Tracklet));
-//   }
-// }
-
-// GPUg() void layerCellsKernel(DeviceStoreNV& devStore, const int layerIndex,
-//                              Vector<Cell> cellsVector)
-// {
-//   computeLayerCells(devStore, layerIndex, cellsVector);
-// }
-
-// GPUg() void sortCellsKernel(DeviceStoreNV& devStore, const int layerIndex,
-//                             Vector<Cell> tempCellsArray)
-// {
-//   const int currentCellIndex = static_cast<int>(blockDim.x * blockIdx.x + threadIdx.x);
-
-//   if (currentCellIndex < tempCellsArray.size()) {
+template <typename T>
+struct trackletSortIndexFunctor : public thrust::binary_function<T, T, bool> {
+  GPUhd() bool operator()(const T& lhs, const T& rhs) const
+  {
+    return lhs.firstClusterIndex < rhs.firstClusterIndex || (lhs.firstClusterIndex == rhs.firstClusterIndex && lhs.secondClusterIndex < rhs.secondClusterIndex);
+  }
+};
 
-//     const int firstTrackletIndex = tempCellsArray[currentCellIndex].getFirstTrackletIndex();
-//     const int offset = atomicAdd(&devStore.getCellsPerTrackletTable()[layerIndex - 1][firstTrackletIndex],
-//                                  -1) -
-//                        1;
-//     const int startIndex = devStore.getCellsLookupTable()[layerIndex - 1][firstTrackletIndex];
+// Compute the tracklets for a given layer
+template <int NLayers = 7>
+GPUg() void computeLayerTrackletsKernel(
+  const int rof0,
+  const int maxRofs,
+  const int layerIndex,
+  const Cluster* clustersCurrentLayer,        // input data rof0
+  const Cluster* clustersNextLayer,           // input data rof0-delta <rof0< rof0+delta (up to 3 rofs)
+  const int* indexTable,                      // input data rof0-delta <rof0< rof0+delta (up to 3 rofs)
+  const int* roFrameClusters,                 // input data O(1)
+  const int* roFrameClustersNext,             // input data O(1)
+  const unsigned char* usedClustersLayer,     // input data rof0
+  const unsigned char* usedClustersNextLayer, // input data rof1
+  const Vertex* vertices,                     // input data
+  int* trackletsLookUpTable,                  // output data
+  Tracklet* tracklets,                        // output data
+  const int nVertices,
+  const int currentLayerClustersSize,
+  const float phiCut,
+  const float minR,
+  const float maxR,
+  const float meanDeltaR,
+  const float positionResolution,
+  const float mSAngle,
+  const StaticTrackingParameters<NLayers>* trkPars,
+  const IndexTableUtils* utils,
+  const unsigned int maxTrackletsPerCluster = 50)
+{
+  for (int currentClusterIndex = blockIdx.x * blockDim.x + threadIdx.x; currentClusterIndex < currentLayerClustersSize; currentClusterIndex += blockDim.x * gridDim.x) {
+    unsigned int storedTracklets{0};
+    const Cluster& currentCluster{clustersCurrentLayer[currentClusterIndex]};
+    const int currentSortedIndex{roFrameClusters[rof0] + currentClusterIndex};
+    if (usedClustersLayer[currentSortedIndex]) {
+      continue;
+    }
+    int minRof = (rof0 >= trkPars->DeltaROF) ? rof0 - trkPars->DeltaROF : 0;
+    int maxRof = (rof0 == maxRofs - trkPars->DeltaROF) ? rof0 : rof0 + trkPars->DeltaROF;
+    const float inverseR0{1.f / currentCluster.radius};
+    for (int iPrimaryVertex{0}; iPrimaryVertex < nVertices; iPrimaryVertex++) {
+      const auto& primaryVertex{vertices[iPrimaryVertex]};
+      if (primaryVertex.getX() == 0.f && primaryVertex.getY() == 0.f && primaryVertex.getZ() == 0.f) {
+        continue;
+      }
+      const float resolution{o2::gpu::GPUCommonMath::Sqrt(Sq(trkPars->PVres) / primaryVertex.getNContributors() + Sq(positionResolution))};
+      const float tanLambda{(currentCluster.zCoordinate - primaryVertex.getZ()) * inverseR0};
+      const float zAtRmin{tanLambda * (minR - currentCluster.radius) + currentCluster.zCoordinate};
+      const float zAtRmax{tanLambda * (maxR - currentCluster.radius) + currentCluster.zCoordinate};
+      const float sqInverseDeltaZ0{1.f / (Sq(currentCluster.zCoordinate - primaryVertex.getZ()) + 2.e-8f)}; /// protecting from overflows adding the detector resolution
+      const float sigmaZ{std::sqrt(Sq(resolution) * Sq(tanLambda) * ((Sq(inverseR0) + sqInverseDeltaZ0) * Sq(meanDeltaR) + 1.f) + Sq(meanDeltaR * mSAngle))};
+
+      const int4 selectedBinsRect{getBinsRect(currentCluster, layerIndex, *utils, zAtRmin, zAtRmax, sigmaZ * trkPars->NSigmaCut, phiCut)};
+      if (selectedBinsRect.x == 0 && selectedBinsRect.y == 0 && selectedBinsRect.z == 0 && selectedBinsRect.w == 0) {
+        continue;
+      }
+      int phiBinsNum{selectedBinsRect.w - selectedBinsRect.y + 1};
+      if (phiBinsNum < 0) {
+        phiBinsNum += trkPars->PhiBins;
+      }
+      constexpr int tableSize{256 * 128 + 1}; // hardcoded for the time being
+
+      for (int rof1{minRof}; rof1 <= maxRof; ++rof1) {
+        if (!(roFrameClustersNext[rof1 + 1] - roFrameClustersNext[rof1])) { // number of clusters on next layer > 0
+          continue;
+        }
+        for (int iPhiCount{0}; iPhiCount < phiBinsNum; iPhiCount++) {
+          int iPhiBin = (selectedBinsRect.y + iPhiCount) % trkPars->PhiBins;
+          const int firstBinIndex{utils->getBinIndex(selectedBinsRect.x, iPhiBin)};
+          const int maxBinIndex{firstBinIndex + selectedBinsRect.z - selectedBinsRect.x + 1};
+          const int firstRowClusterIndex = indexTable[rof1 * tableSize + firstBinIndex];
+          const int maxRowClusterIndex = indexTable[rof1 * tableSize + maxBinIndex];
+          for (int iNextCluster{firstRowClusterIndex}; iNextCluster < maxRowClusterIndex; ++iNextCluster) {
+            if (iNextCluster >= (roFrameClustersNext[rof1 + 1] - roFrameClustersNext[rof1])) {
+              break;
+            }
+            const Cluster& nextCluster{getPtrFromRuler<Cluster>(rof1, clustersNextLayer, roFrameClustersNext)[iNextCluster]};
+            if (usedClustersNextLayer[nextCluster.clusterId]) {
+              continue;
+            }
+            const float deltaPhi{o2::gpu::GPUCommonMath::Abs(currentCluster.phi - nextCluster.phi)};
+            const float deltaZ{o2::gpu::GPUCommonMath::Abs(tanLambda * (nextCluster.radius - currentCluster.radius) + currentCluster.zCoordinate - nextCluster.zCoordinate)};
+
+            if (deltaZ / sigmaZ < trkPars->NSigmaCut && (deltaPhi < phiCut || o2::gpu::GPUCommonMath::Abs(deltaPhi - constants::math::TwoPi) < phiCut)) {
+              trackletsLookUpTable[currentSortedIndex]++; // Race-condition safe
+              const float phi{o2::gpu::GPUCommonMath::ATan2(currentCluster.yCoordinate - nextCluster.yCoordinate, currentCluster.xCoordinate - nextCluster.xCoordinate)};
+              const float tanL{(currentCluster.zCoordinate - nextCluster.zCoordinate) / (currentCluster.radius - nextCluster.radius)};
+              const size_t stride{currentClusterIndex * maxTrackletsPerCluster};
+              new (tracklets + stride + storedTracklets) Tracklet{currentSortedIndex, roFrameClustersNext[rof1] + iNextCluster, tanL, phi, rof0, rof1};
+              ++storedTracklets;
+            }
+          }
+        }
+      }
+    }
+    if (storedTracklets > maxTrackletsPerCluster) {
+      printf("its-gpu-tracklet finder: found more tracklets per clusters (%d) than maximum set (%d), check the configuration!\n", maxTrackletsPerCluster, storedTracklets);
+    }
+  }
+}
 
-//     memcpy(&devStore.getCells()[layerIndex][startIndex + offset], &tempCellsArray[currentCellIndex],
-//            sizeof(Cell));
-//   }
-// }
+// Decrease LUT entries corresponding to duplicated tracklets. NB: duplicate tracklets are removed separately (see const Tracklets*).
+GPUg() void removeDuplicateTrackletsEntriesLUTKernel(
+  int* trackletsLookUpTable,
+  const Tracklet* tracklets,
+  const int* nTracklets,
+  const int layerIndex)
+{
+  int id0{-1}, id1{-1};
+  for (int iTracklet{0}; iTracklet < nTracklets[layerIndex]; ++iTracklet) {
+    auto& trk = tracklets[iTracklet];
+    if (trk.firstClusterIndex == id0 && trk.secondClusterIndex == id1) {
+      trackletsLookUpTable[id0]--;
+    } else {
+      id0 = trk.firstClusterIndex;
+      id1 = trk.secondClusterIndex;
+    }
+  }
+}
 
+// Compute cells kernel
+template <bool initRun, int NLayers = 7>
+GPUg() void computeLayerCellsKernel(
+  const Tracklet* trackletsCurrentLayer,
+  const Tracklet* trackletsNextLayer,
+  const int* trackletsCurrentLayerLUT,
+  const int nTrackletsCurrent,
+  Cell* cells,
+  int* cellsLUT,
+  const StaticTrackingParameters<NLayers>* trkPars)
+{
+  for (int iCurrentTrackletIndex = blockIdx.x * blockDim.x + threadIdx.x; iCurrentTrackletIndex < nTrackletsCurrent; iCurrentTrackletIndex += blockDim.x * gridDim.x) {
+    const Tracklet& currentTracklet = trackletsCurrentLayer[iCurrentTrackletIndex];
+    const int nextLayerClusterIndex{currentTracklet.secondClusterIndex};
+    const int nextLayerFirstTrackletIndex{trackletsCurrentLayerLUT[nextLayerClusterIndex]};
+    const int nextLayerLastTrackletIndex{trackletsCurrentLayerLUT[nextLayerClusterIndex + 1]};
+    if (nextLayerFirstTrackletIndex == nextLayerLastTrackletIndex) {
+      continue;
+    }
+    int foundCells{0};
+    for (int iNextTrackletIndex{nextLayerFirstTrackletIndex}; iNextTrackletIndex < nextLayerLastTrackletIndex; ++iNextTrackletIndex) {
+      if (trackletsNextLayer[iNextTrackletIndex].firstClusterIndex != nextLayerClusterIndex) {
+        break;
+      }
+      const Tracklet& nextTracklet = trackletsNextLayer[iNextTrackletIndex];
+      const float deltaTanLambda{o2::gpu::GPUCommonMath::Abs(currentTracklet.tanLambda - nextTracklet.tanLambda)};
+      const float tanLambda{(currentTracklet.tanLambda + nextTracklet.tanLambda) * 0.5f};
+
+      if (deltaTanLambda / trkPars->CellDeltaTanLambdaSigma < trkPars->NSigmaCut) {
+        if constexpr (!initRun) {
+          new (cells + cellsLUT[iCurrentTrackletIndex] + foundCells) Cell{currentTracklet.firstClusterIndex, nextTracklet.firstClusterIndex,
+                                                                          nextTracklet.secondClusterIndex,
+                                                                          iCurrentTrackletIndex,
+                                                                          iNextTrackletIndex,
+                                                                          tanLambda};
+        }
+        ++foundCells;
+      }
+    }
+    if constexpr (initRun) {
+      // Fill cell Lookup table
+      cellsLUT[iCurrentTrackletIndex] = foundCells;
+    }
+  }
+}
 } // namespace gpu
 
-// void TrackeTraitsGPU::adoptTimeFrame(TimeFrame* tf)
-// {
-//   mTimeFrameGPU = tf;
-// }
-// TrackerTraits* createTrackerTraitsGPU()
-// {
-//   return new TrackerTraitsGPU;
-// }
+template <int NLayers>
+void TrackerTraitsGPU<NLayers>::initialiseTimeFrame(const int iteration)
+{
+  mTimeFrameGPU->initialise(iteration, mTrkParams[iteration], NLayers);
+  setIsGPU(true);
+}
 
 template <int NLayers>
 void TrackerTraitsGPU<NLayers>::computeLayerTracklets(const int iteration)
 {
-  //   PrimaryVertexContextNV* primaryVertexContext = static_cast<PrimaryVertexContextNV*>(nullptr); //TODO: FIX THIS with Time Frames
-
-  checkGPUError(cudaMemcpyToSymbol(gpu::kTrkPar, &mTrkParams, sizeof(gpu::StaticTrackingParameters<7>)), __FILE__, __LINE__);
-  //   std::array<size_t, constants::its2::CellsPerRoad> tempSize;
-  //   std::array<int, constants::its2::CellsPerRoad> trackletsNum;
-  //   std::array<gpu::Stream, constants::its2::TrackletsPerRoad> streamArray;
-
-  //   for (int iLayer{0}; iLayer < constants::its2::CellsPerRoad; ++iLayer) {
-
-  //     tempSize[iLayer] = 0;
-  //     primaryVertexContext->getTempTrackletArray()[iLayer].reset(
-  //       static_cast<int>(primaryVertexContext->getDeviceTracklets()[iLayer + 1].capacity()));
-
-  //     cub::DeviceScan::ExclusiveSum(static_cast<void*>(NULL), tempSize[iLayer],
-  //                                   primaryVertexContext->getDeviceTrackletsPerClustersTable()[iLayer].get(),
-  //                                   primaryVertexContext->getDeviceTrackletsLookupTable()[iLayer].get(),
-  //                                   primaryVertexContext->getClusters()[iLayer + 1].size());
-
-  //     primaryVertexContext->getTempTableArray()[iLayer].reset(static_cast<int>(tempSize[iLayer]));
-  //   }
-
-  //   cudaDeviceSynchronize();
-
-  //   for (int iLayer{0}; iLayer < constants::its2::TrackletsPerRoad; ++iLayer) {
-
-  //     const gpu::DeviceProperties& deviceProperties = gpu::Context::getInstance().getDeviceProperties();
-  //     const int clustersNum{static_cast<int>(primaryVertexContext->getClusters()[iLayer].size())};
-  //     dim3 threadsPerBlock{gpu::utils::host::getBlockSize(clustersNum, 1, 192)};
-  //     dim3 blocksGrid{gpu::utils::host::getBlocksGrid(threadsPerBlock, clustersNum)};
-
-  //     if (iLayer == 0) {
-
-  //       gpu::layerTrackletsKernel<<<blocksGrid, threadsPerBlock, 0, streamArray[iLayer].get()>>>(primaryVertexContext->getDeviceContext(),
-  //                                                                                                iLayer, primaryVertexContext->getDeviceTracklets()[iLayer].getWeakCopy());
-
-  //     } else {
-
-  //       gpu::layerTrackletsKernel<<<blocksGrid, threadsPerBlock, 0, streamArray[iLayer].get()>>>(primaryVertexContext->getDeviceContext(),
-  //                                                                                                iLayer, primaryVertexContext->getTempTrackletArray()[iLayer - 1].getWeakCopy());
-  //     }
-
-  //     cudaError_t error = cudaGetLastError();
-
-  //     if (error != cudaSuccess) {
-
-  //       std::ostringstream errorString{};
-  //       errorString << "CUDA API returned error [" << cudaGetErrorString(error) << "] (code " << error << ")"
-  //                   << std::endl;
-
-  //       throw std::runtime_error{errorString.str()};
-  //     }
-  //   }
-
-  //   cudaDeviceSynchronize();
-
-  //   for (int iLayer{0}; iLayer < constants::its2::CellsPerRoad; ++iLayer) {
-
-  //     trackletsNum[iLayer] = primaryVertexContext->getTempTrackletArray()[iLayer].getSizeFromDevice();
-  //     if (trackletsNum[iLayer] == 0) {
-  //       continue;
-  //     }
-  //     primaryVertexContext->getDeviceTracklets()[iLayer + 1].resize(trackletsNum[iLayer]);
-
-  //     cub::DeviceScan::ExclusiveSum(static_cast<void*>(primaryVertexContext->getTempTableArray()[iLayer].get()), tempSize[iLayer],
-  //                                   primaryVertexContext->getDeviceTrackletsPerClustersTable()[iLayer].get(),
-  //                                   primaryVertexContext->getDeviceTrackletsLookupTable()[iLayer].get(),
-  //                                   primaryVertexContext->getClusters()[iLayer + 1].size(), streamArray[iLayer + 1].get());
-
-  //     dim3 threadsPerBlock{gpu::utils::host::getBlockSize(trackletsNum[iLayer])};
-  //     dim3 blocksGrid{gpu::utils::host::getBlocksGrid(threadsPerBlock, trackletsNum[iLayer])};
-
-  //     gpu::sortTrackletsKernel<<<blocksGrid, threadsPerBlock, 0, streamArray[iLayer + 1].get()>>>(primaryVertexContext->getDeviceContext(),
-  //                                                                                                 iLayer + 1, primaryVertexContext->getTempTrackletArray()[iLayer].getWeakCopy());
-
-  //     cudaError_t error = cudaGetLastError();
-
-  //     if (error != cudaSuccess) {
+  const Vertex diamondVert({mTrkParams[iteration].Diamond[0], mTrkParams[iteration].Diamond[1], mTrkParams[iteration].Diamond[2]}, {25.e-6f, 0.f, 0.f, 25.e-6f, 0.f, 36.f}, 1, 1.f);
+  gsl::span<const Vertex> diamondSpan(&diamondVert, 1);
+  size_t bufferSize = mTimeFrameGPU->getConfig().tmpCUBBufferSize;
+
+  for (int rof0{0}; rof0 < mTimeFrameGPU->getNrof(); ++rof0) {
+    gsl::span<const Vertex> primaryVertices = mTrkParams[iteration].UseDiamond ? diamondSpan : mTimeFrameGPU->getPrimaryVertices(rof0); // replace with GPU one
+    std::vector<Vertex> paddedVertices;
+    for (int iVertex{0}; iVertex < mTimeFrameGPU->getConfig().maxVerticesCapacity; ++iVertex) {
+      if (iVertex < primaryVertices.size()) {
+        paddedVertices.emplace_back(primaryVertices[iVertex]);
+      } else {
+        paddedVertices.emplace_back(Vertex());
+      }
+    }
+    checkGPUError(cudaMemcpy(mTimeFrameGPU->getDeviceVertices(rof0), paddedVertices.data(), mTimeFrameGPU->getConfig().maxVerticesCapacity * sizeof(Vertex), cudaMemcpyHostToDevice), __FILE__, __LINE__);
+    for (int iLayer{0}; iLayer < NLayers - 1; ++iLayer) {
+      const dim3 threadsPerBlock{gpu::utils::host::getBlockSize(mTimeFrameGPU->getNClustersLayer(rof0, iLayer))};
+      const dim3 blocksGrid{gpu::utils::host::getBlocksGrid(threadsPerBlock, mTimeFrameGPU->getNClustersLayer(rof0, iLayer))};
+      const float meanDeltaR{mTrkParams[iteration].LayerRadii[iLayer + 1] - mTrkParams[iteration].LayerRadii[iLayer]};
+
+      if (!mTimeFrameGPU->getClustersOnLayer(rof0, iLayer).size()) {
+        LOGP(info, "Skipping ROF0: {}, no clusters found on layer {}", rof0, iLayer);
+        continue;
+      }
+
+      gpu::computeLayerTrackletsKernel<<<blocksGrid, threadsPerBlock, 0, mTimeFrameGPU->getStream(iLayer).get()>>>(
+        rof0,
+        mTimeFrameGPU->getNrof(),
+        iLayer,
+        mTimeFrameGPU->getDeviceClustersOnLayer(rof0, iLayer),       // :checked:
+        mTimeFrameGPU->getDeviceClustersOnLayer(0, iLayer + 1),      // :checked:
+        mTimeFrameGPU->getDeviceIndexTables(iLayer + 1),             // :checked:
+        mTimeFrameGPU->getDeviceROframesClustersOnLayer(iLayer),     // :checked:
+        mTimeFrameGPU->getDeviceROframesClustersOnLayer(iLayer + 1), // :checked:
+        mTimeFrameGPU->getDeviceUsedClustersOnLayer(0, iLayer),      // :checked:
+        mTimeFrameGPU->getDeviceUsedClustersOnLayer(0, iLayer + 1),  // :checked:
+        mTimeFrameGPU->getDeviceVertices(rof0),
+        mTimeFrameGPU->getDeviceTrackletsLookupTable(0, iLayer),
+        mTimeFrameGPU->getDeviceTracklets(rof0, iLayer),
+        mTimeFrameGPU->getConfig().maxVerticesCapacity,
+        mTimeFrameGPU->getNClustersLayer(rof0, iLayer),
+        mTimeFrameGPU->getPhiCut(iLayer),
+        mTimeFrameGPU->getMinR(iLayer + 1),
+        mTimeFrameGPU->getMaxR(iLayer + 1),
+        meanDeltaR,
+        mTimeFrameGPU->getPositionResolution(iLayer),
+        mTimeFrameGPU->getMSangle(iLayer),
+        mTimeFrameGPU->getDeviceTrackingParameters(),
+        mTimeFrameGPU->getDeviceIndexTableUtils(),
+        mTimeFrameGPU->getConfig().maxTrackletsPerCluster);
+    }
+  }
 
-  //       std::ostringstream errorString{};
-  //       errorString << "CUDA API returned error [" << cudaGetErrorString(error) << "] (code " << error << ")"
-  //                   << std::endl;
+  for (int iLayer{0}; iLayer < NLayers - 1; ++iLayer) {
+    // Sort tracklets to put empty ones on the right side of the array.
+    auto thrustTrackletsBegin = thrust::device_ptr<o2::its::Tracklet>(mTimeFrameGPU->getDeviceTrackletsAll(iLayer));
+    auto thrustTrackletsEnd = thrust::device_ptr<o2::its::Tracklet>(mTimeFrameGPU->getDeviceTrackletsAll(iLayer) + mTimeFrameGPU->getConfig().trackletsCapacity);
+    thrust::sort(thrustTrackletsBegin, thrustTrackletsEnd, gpu::trackletSortEmptyFunctor<o2::its::Tracklet>());
+
+    // Get number of found tracklets
+    // With thrust:
+    // auto thrustTrackletsLUTbegin = thrust::device_ptr<int>(mTimeFrameGPU->getDeviceTrackletsLookupTable(0, iLayer));
+    // auto thrustTrackletsLUTend = thrust::device_ptr<int>(mTimeFrameGPU->getDeviceTrackletsLookupTable(0, iLayer) + mTimeFrameGPU->mClusters[iLayer].size());
+    // mTimeFrameGPU->getTrackletSizeHost()[iLayer] = thrust::reduce(thrustTrackletsLUTbegin, thrustTrackletsLUTend, 0);
+
+    discardResult(cub::DeviceReduce::Sum(reinterpret_cast<void*>(mTimeFrameGPU->getDeviceCUBBuffer(iLayer)), // d_temp_storage
+                                         bufferSize,                                                         // temp_storage_bytes
+                                         mTimeFrameGPU->getDeviceTrackletsLookupTable(0, iLayer),            // d_in
+                                         mTimeFrameGPU->getDeviceNFoundTracklets() + iLayer,                 // d_out
+                                         mTimeFrameGPU->mClusters[iLayer].size(),                            // num_items
+                                         mTimeFrameGPU->getStream(iLayer).get()));
+  }
+  discardResult(cudaDeviceSynchronize());
+  checkGPUError(cudaMemcpy(mTimeFrameGPU->getTrackletSizeHost().data(), mTimeFrameGPU->getDeviceNFoundTracklets(), (NLayers - 1) * sizeof(int), cudaMemcpyDeviceToHost), __FILE__, __LINE__);
+  for (int iLayer{0}; iLayer < NLayers - 1; ++iLayer) {
+
+    // Sort tracklets according to cluster ids
+    auto thrustTrackletsBegin = thrust::device_ptr<o2::its::Tracklet>(mTimeFrameGPU->getDeviceTrackletsAll(iLayer));
+    auto thrustTrackletsEnd = thrust::device_ptr<o2::its::Tracklet>(mTimeFrameGPU->getDeviceTrackletsAll(iLayer) + mTimeFrameGPU->getTrackletSizeHost()[iLayer]);
+    thrust::sort(thrustTrackletsBegin, thrustTrackletsEnd, gpu::trackletSortIndexFunctor<o2::its::Tracklet>());
+  }
+  discardResult(cudaDeviceSynchronize());
+  for (int iLayer{0}; iLayer < NLayers - 1; ++iLayer) {
+    // Remove duplicate entries in LUTs, done by single thread so far
+    gpu::removeDuplicateTrackletsEntriesLUTKernel<<<1, 1, 0, mTimeFrameGPU->getStream(iLayer).get()>>>(
+      mTimeFrameGPU->getDeviceTrackletsLookupTable(0, iLayer),
+      mTimeFrameGPU->getDeviceTrackletsAll(iLayer),
+      mTimeFrameGPU->getDeviceNFoundTracklets(),
+      iLayer);
+  }
+  // Remove actual tracklet duplicates
+  for (int iLayer{0}; iLayer < NLayers - 1; ++iLayer) {
+    auto begin = thrust::device_ptr<o2::its::Tracklet>(mTimeFrameGPU->getDeviceTrackletsAll(iLayer));
+    auto end = thrust::device_ptr<o2::its::Tracklet>(mTimeFrameGPU->getDeviceTrackletsAll(iLayer) + mTimeFrameGPU->getTrackletSizeHost()[iLayer]);
 
-  //       throw std::runtime_error{errorString.str()};
-  //     }
-  //   }
+    auto new_end = thrust::unique(begin, end);
+    mTimeFrameGPU->getTrackletSizeHost()[iLayer] = new_end - begin;
+  }
+  discardResult(cudaDeviceSynchronize());
+  for (int iLayer{0}; iLayer < NLayers - 1; ++iLayer) {
+    // Compute LUT
+    discardResult(cub::DeviceScan::ExclusiveSum(reinterpret_cast<void*>(mTimeFrameGPU->getDeviceCUBBuffer(iLayer)), // d_temp_storage
+                                                bufferSize,                                                         // temp_storage_bytes
+                                                mTimeFrameGPU->getDeviceTrackletsLookupTable(0, iLayer),            // d_in
+                                                mTimeFrameGPU->getDeviceTrackletsLookupTable(0, iLayer),            // d_out
+                                                mTimeFrameGPU->mClusters[iLayer].size(),                            // num_items
+                                                mTimeFrameGPU->getStream(iLayer).get()));
+  }
+  discardResult(cudaDeviceSynchronize());
+
+  // Create tracklets labels, at the moment on the host
+  if (mTimeFrameGPU->hasMCinformation()) {
+    for (int iLayer{0}; iLayer < mTrkParams[iteration].TrackletsPerRoad(); ++iLayer) {
+      std::vector<o2::its::Tracklet> tracklets(mTimeFrameGPU->getTrackletSizeHost()[iLayer]);
+      checkGPUError(cudaMemcpy(tracklets.data(), mTimeFrameGPU->getDeviceTrackletsAll(iLayer), mTimeFrameGPU->getTrackletSizeHost()[iLayer] * sizeof(o2::its::Tracklet), cudaMemcpyDeviceToHost), __FILE__, __LINE__);
+      for (auto& trk : tracklets) {
+        MCCompLabel label;
+        int currentId{mTimeFrameGPU->mClusters[iLayer][trk.firstClusterIndex].clusterId};
+        int nextId{mTimeFrameGPU->mClusters[iLayer + 1][trk.secondClusterIndex].clusterId};
+        for (auto& lab1 : mTimeFrameGPU->getClusterLabels(iLayer, currentId)) {
+          for (auto& lab2 : mTimeFrameGPU->getClusterLabels(iLayer + 1, nextId)) {
+            if (lab1 == lab2 && lab1.isValid()) {
+              label = lab1;
+              break;
+            }
+          }
+          if (label.isValid()) {
+            break;
+          }
+        }
+        mTimeFrameGPU->getTrackletsLabel(iLayer).emplace_back(label);
+      }
+    }
+  }
 }
 
-// void TrackerTraitsGPU::computeLayerCells()
-// {
-
-//   PrimaryVertexContextNV* primaryVertexContext = static_cast<PrimaryVertexContextNV*>(nullptr); //TODO: FIX THIS with Time Frames
-//   std::array<size_t, constants::its2::CellsPerRoad - 1> tempSize;
-//   std::array<int, constants::its2::CellsPerRoad - 1> trackletsNum;
-//   std::array<int, constants::its2::CellsPerRoad - 1> cellsNum;
-//   std::array<gpu::Stream, constants::its2::CellsPerRoad> streamArray;
-
-//   for (int iLayer{0}; iLayer < constants::its2::CellsPerRoad - 1; ++iLayer) {
-
-//     tempSize[iLayer] = 0;
-//     trackletsNum[iLayer] = primaryVertexContext->getDeviceTracklets()[iLayer + 1].getSizeFromDevice();
-//     primaryVertexContext->getTempCellArray()[iLayer].reset(
-//       static_cast<int>(primaryVertexContext->getDeviceCells()[iLayer + 1].capacity()));
-//     if (trackletsNum[iLayer] == 0) {
-//       continue;
-//     }
-//     cub::DeviceScan::ExclusiveSum(static_cast<void*>(NULL), tempSize[iLayer],
-//                                   primaryVertexContext->getDeviceCellsPerTrackletTable()[iLayer].get(),
-//                                   primaryVertexContext->getDeviceCellsLookupTable()[iLayer].get(), trackletsNum[iLayer]);
-
-//     primaryVertexContext->getTempTableArray()[iLayer].reset(static_cast<int>(tempSize[iLayer]));
-//   }
-
-//   cudaDeviceSynchronize();
-
-//   for (int iLayer{0}; iLayer < constants::its2::CellsPerRoad; ++iLayer) {
-//     const gpu::DeviceProperties& deviceProperties = gpu::Context::getInstance().getDeviceProperties();
-//     const int trackletsSize = primaryVertexContext->getDeviceTracklets()[iLayer].getSizeFromDevice();
-//     if (trackletsSize == 0) {
-//       continue;
-//     }
-//     dim3 threadsPerBlock{gpu::utils::host::getBlockSize(trackletsSize)};
-//     dim3 blocksGrid{gpu::utils::host::getBlocksGrid(threadsPerBlock, trackletsSize)};
-
-//     if (iLayer == 0) {
-
-//       gpu::layerCellsKernel<<<blocksGrid, threadsPerBlock, 0, streamArray[iLayer].get()>>>(primaryVertexContext->getDeviceContext(),
-//                                                                                            iLayer, primaryVertexContext->getDeviceCells()[iLayer].getWeakCopy());
-
-//     } else {
-
-//       gpu::layerCellsKernel<<<blocksGrid, threadsPerBlock, 0, streamArray[iLayer].get()>>>(primaryVertexContext->getDeviceContext(),
-//                                                                                            iLayer, primaryVertexContext->getTempCellArray()[iLayer - 1].getWeakCopy());
-//     }
-
-//     cudaError_t error = cudaGetLastError();
-
-//     if (error != cudaSuccess) {
-
-//       std::ostringstream errorString{};
-//       errorString << "CUDA API returned error [" << cudaGetErrorString(error) << "] (code " << error << ")"
-//                   << std::endl;
-
-//       throw std::runtime_error{errorString.str()};
-//     }
-//   }
-
-//   cudaDeviceSynchronize();
-
-//   for (int iLayer{0}; iLayer < constants::its2::CellsPerRoad - 1; ++iLayer) {
-//     cellsNum[iLayer] = primaryVertexContext->getTempCellArray()[iLayer].getSizeFromDevice();
-//     if (cellsNum[iLayer] == 0) {
-//       continue;
-//     }
-//     primaryVertexContext->getDeviceCells()[iLayer + 1].resize(cellsNum[iLayer]);
-
-//     cub::DeviceScan::ExclusiveSum(static_cast<void*>(primaryVertexContext->getTempTableArray()[iLayer].get()), tempSize[iLayer],
-//                                   primaryVertexContext->getDeviceCellsPerTrackletTable()[iLayer].get(),
-//                                   primaryVertexContext->getDeviceCellsLookupTable()[iLayer].get(), trackletsNum[iLayer],
-//                                   streamArray[iLayer + 1].get());
-
-//     dim3 threadsPerBlock{gpu::utils::host::getBlockSize(trackletsNum[iLayer])};
-//     dim3 blocksGrid{gpu::utils::host::getBlocksGrid(threadsPerBlock, trackletsNum[iLayer])};
-
-//     gpu::sortCellsKernel<<<blocksGrid, threadsPerBlock, 0, streamArray[iLayer + 1].get()>>>(primaryVertexContext->getDeviceContext(),
-//                                                                                             iLayer + 1, primaryVertexContext->getTempCellArray()[iLayer].getWeakCopy());
-
-//     cudaError_t error = cudaGetLastError();
-
-//     if (error != cudaSuccess) {
-
-//       std::ostringstream errorString{};
-//       errorString << "CUDA API returned error [" << cudaGetErrorString(error) << "] (code " << error << ")"
-//                   << std::endl;
-
-//       throw std::runtime_error{errorString.str()};
-//     }
-//   }
-
-//   cudaDeviceSynchronize();
-
-//   for (int iLayer{0}; iLayer < constants::its2::CellsPerRoad; ++iLayer) {
-
-//     int cellsSize = 0;
-//     if (iLayer == 0) {
-
-//       cellsSize = primaryVertexContext->getDeviceCells()[iLayer].getSizeFromDevice();
-//       if (cellsSize == 0) {
-//         continue;
-//       }
-//     } else {
-
-//       cellsSize = cellsNum[iLayer - 1];
-//       if (cellsSize == 0) {
-//         continue;
-//       }
-//       primaryVertexContext->getDeviceCellsLookupTable()[iLayer - 1].copyIntoVector(
-//         primaryVertexContext->getCellsLookupTable()[iLayer - 1], trackletsNum[iLayer - 1]);
-//     }
+template <int NLayers>
+void TrackerTraitsGPU<NLayers>::computeLayerCells(const int iteration)
+{
+  size_t bufferSize = mTimeFrameGPU->getConfig().tmpCUBBufferSize;
+  for (int iLayer{0}; iLayer < NLayers - 2; ++iLayer) {
+    if (!mTimeFrameGPU->getTrackletSizeHost()[iLayer + 1] ||
+        !mTimeFrameGPU->getTrackletSizeHost()[iLayer]) {
+      continue;
+    }
+    const dim3 threadsPerBlock{gpu::utils::host::getBlockSize(mTimeFrameGPU->getTrackletSizeHost()[iLayer])};
+    const dim3 blocksGrid{gpu::utils::host::getBlocksGrid(threadsPerBlock, mTimeFrameGPU->getTrackletSizeHost()[iLayer])};
+    gpu::computeLayerCellsKernel<true><<<blocksGrid, threadsPerBlock, 0, mTimeFrameGPU->getStream(iLayer).get()>>>(
+      mTimeFrameGPU->getDeviceTrackletsAll(iLayer),
+      mTimeFrameGPU->getDeviceTrackletsAll(iLayer + 1),
+      mTimeFrameGPU->getDeviceTrackletsLookupTable(0, iLayer + 1),
+      mTimeFrameGPU->getTrackletSizeHost()[iLayer],
+      nullptr,
+      mTimeFrameGPU->getDeviceCellsLookupTable(iLayer),
+      mTimeFrameGPU->getDeviceTrackingParameters());
+
+    // Compute number of found Cells
+    discardResult(cub::DeviceReduce::Sum(reinterpret_cast<void*>(mTimeFrameGPU->getDeviceCUBBuffer(iLayer)), // d_temp_storage
+                                         bufferSize,                                                         // temp_storage_bytes
+                                         mTimeFrameGPU->getDeviceCellsLookupTable(iLayer),                   // d_in
+                                         mTimeFrameGPU->getDeviceNFoundCells() + iLayer,                     // d_out
+                                         mTimeFrameGPU->getTrackletSizeHost()[iLayer],                       // num_items
+                                         mTimeFrameGPU->getStream(iLayer).get()));
+
+    // Compute LUT
+    discardResult(cub::DeviceScan::ExclusiveSum(reinterpret_cast<void*>(mTimeFrameGPU->getDeviceCUBBuffer(iLayer)), // d_temp_storage
+                                                bufferSize,                                                         // temp_storage_bytes
+                                                mTimeFrameGPU->getDeviceCellsLookupTable(iLayer),                   // d_in
+                                                mTimeFrameGPU->getDeviceCellsLookupTable(iLayer),                   // d_out
+                                                mTimeFrameGPU->getTrackletSizeHost()[iLayer],                       // num_items
+                                                mTimeFrameGPU->getStream(iLayer).get()));
+    discardResult(cudaStreamSynchronize(mTimeFrameGPU->getStream(iLayer).get()));
+
+    gpu::computeLayerCellsKernel<false><<<blocksGrid, threadsPerBlock, 0, mTimeFrameGPU->getStream(iLayer).get()>>>(
+      mTimeFrameGPU->getDeviceTrackletsAll(iLayer),
+      mTimeFrameGPU->getDeviceTrackletsAll(iLayer + 1),
+      mTimeFrameGPU->getDeviceTrackletsLookupTable(0, iLayer + 1),
+      mTimeFrameGPU->getTrackletSizeHost()[iLayer],
+      mTimeFrameGPU->getDeviceCells(iLayer),
+      mTimeFrameGPU->getDeviceCellsLookupTable(iLayer),
+      mTimeFrameGPU->getDeviceTrackingParameters());
+  }
 
-//     primaryVertexContext->getDeviceCells()[iLayer].copyIntoVector(primaryVertexContext->getCells()[iLayer], cellsSize);
-//   }
-// }
+  checkGPUError(cudaMemcpy(mTimeFrameGPU->getCellSizeHost().data(), mTimeFrameGPU->getDeviceNFoundCells(), (NLayers - 2) * sizeof(int), cudaMemcpyDeviceToHost), __FILE__, __LINE__);
+  /// Create cells labels
+  if (mTimeFrameGPU->hasMCinformation()) {
+    for (int iLayer{0}; iLayer < NLayers - 2; ++iLayer) {
+      std::vector<o2::its::Cell> cells(mTimeFrameGPU->getCellSizeHost()[iLayer]);
+      checkGPUError(cudaMemcpy(cells.data(), mTimeFrameGPU->getDeviceCells(iLayer), mTimeFrameGPU->getCellSizeHost()[iLayer] * sizeof(o2::its::Cell), cudaMemcpyDeviceToHost), __FILE__, __LINE__);
+      for (auto& cell : cells) {
+        MCCompLabel currentLab{mTimeFrameGPU->getTrackletsLabel(iLayer)[cell.getFirstTrackletIndex()]};
+        MCCompLabel nextLab{mTimeFrameGPU->getTrackletsLabel(iLayer + 1)[cell.getSecondTrackletIndex()]};
+        mTimeFrameGPU->getCellsLabel(iLayer).emplace_back(currentLab == nextLab ? currentLab : MCCompLabel());
+      }
+    }
+  }
+}
 
 // void TrackerTraitsGPU::refitTracks(const std::vector<std::vector<TrackingFrameInfo>>& tf, std::vector<TrackITSExt>& tracks)
 // {
 //   PrimaryVertexContextNV* pvctx = static_cast<PrimaryVertexContextNV*>(nullptr); //TODO: FIX THIS with Time Frames
-
 //   std::array<const Cell*, 5> cells;
 //   for (int iLayer = 0; iLayer < 5; iLayer++) {
 //     cells[iLayer] = pvctx->getDeviceCells()[iLayer].get();
@@ -540,6 +470,44 @@ void TrackerTraitsGPU<NLayers>::computeLayerTracklets(const int iteration)
 //   //TODO: restore this
 //   // mChainRunITSTrackFit(*mChain, mPrimaryVertexContext->getRoads(), clusters, cells, tf, tracks);
 // }
+
+template <int NLayers>
+void TrackerTraitsGPU<NLayers>::findCellsNeighbours(const int iteration){};
+
+template <int NLayers>
+void TrackerTraitsGPU<NLayers>::findRoads(const int iteration){};
+
+template <int NLayers>
+void TrackerTraitsGPU<NLayers>::findTracks(){};
+
+template <int NLayers>
+void TrackerTraitsGPU<NLayers>::extendTracks(const int iteration){};
+
+template <int NLayers>
+void TrackerTraitsGPU<NLayers>::setBz(float bz)
+{
+  mBz = bz;
+  mTimeFrameGPU->setBz(bz);
+}
+
+template <int NLayers>
+int TrackerTraitsGPU<NLayers>::getTFNumberOfClusters() const
+{
+  return mTimeFrameGPU->getNumberOfClusters();
+}
+
+template <int NLayers>
+int TrackerTraitsGPU<NLayers>::getTFNumberOfTracklets() const
+{
+  return mTimeFrameGPU->getNumberOfTracklets();
+}
+
+template <int NLayers>
+int TrackerTraitsGPU<NLayers>::getTFNumberOfCells() const
+{
+  return mTimeFrameGPU->getNumberOfCells();
+}
+
 template class TrackerTraitsGPU<7>;
 } // namespace its
 } // namespace o2
diff --git a/Detectors/ITSMFT/ITS/tracking/GPU/cuda/Utils.cu b/Detectors/ITSMFT/ITS/tracking/GPU/cuda/Utils.cu
index de858d876f71c..e6b07965c3172 100644
--- a/Detectors/ITSMFT/ITS/tracking/GPU/cuda/Utils.cu
+++ b/Detectors/ITSMFT/ITS/tracking/GPU/cuda/Utils.cu
@@ -60,7 +60,8 @@ void utils::host::checkGPUError(const cudaError_t error, const char* file, const
 {
   if (error != cudaSuccess) {
     std::ostringstream errorString{};
-    errorString << file << ":" << line << GPU_ARCH << " API returned error [" << cudaGetErrorString(error) << "] (code "
+    errorString << file << ":" << line << std::endl
+                << GPU_ARCH << " API returned error [" << cudaGetErrorString(error) << "] (code "
                 << error << ")" << std::endl;
     throw std::runtime_error{errorString.str()};
   }
diff --git a/Detectors/ITSMFT/ITS/tracking/GPU/cuda/VertexerTraitsGPU.cu b/Detectors/ITSMFT/ITS/tracking/GPU/cuda/VertexerTraitsGPU.cu
index 68cac94936658..92d49c53e2ea7 100644
--- a/Detectors/ITSMFT/ITS/tracking/GPU/cuda/VertexerTraitsGPU.cu
+++ b/Detectors/ITSMFT/ITS/tracking/GPU/cuda/VertexerTraitsGPU.cu
@@ -132,6 +132,24 @@ GPUg() void printBufferOnThreadF(const float* v, size_t size, const unsigned int
   }
 }
 
+GPUg() void resetTrackletsKernel(Tracklet* tracklets, const int nTracklets)
+{
+  for (int iCurrentLayerClusterIndex = blockIdx.x * blockDim.x + threadIdx.x; iCurrentLayerClusterIndex < nTracklets; iCurrentLayerClusterIndex += blockDim.x * gridDim.x) {
+    new (tracklets + iCurrentLayerClusterIndex) Tracklet{};
+  }
+}
+
+GPUg() void dumpFoundTrackletsKernel(const Tracklet* tracklets, const int* nTracklet, const size_t nClustersMiddleLayer, const int maxTrackletsPerCluster)
+{
+  for (int iCurrentLayerClusterIndex = blockIdx.x * blockDim.x + threadIdx.x; iCurrentLayerClusterIndex < nClustersMiddleLayer; iCurrentLayerClusterIndex += blockDim.x * gridDim.x) {
+    const int stride{iCurrentLayerClusterIndex * maxTrackletsPerCluster};
+    for (int iTracklet{0}; iTracklet < nTracklet[iCurrentLayerClusterIndex]; ++iTracklet) {
+      auto& t = tracklets[stride + iTracklet];
+      t.dump();
+    }
+  }
+}
+
 GPUg() void dumpMaximaKernel(const cub::KeyValuePair<int, int>* tmpVertexBins, const int threadId)
 {
   if (blockIdx.x * blockDim.x + threadIdx.x == threadId) {
@@ -153,7 +171,8 @@ GPUg() void trackleterKernel(
   Tracklet* Tracklets,
   int* foundTracklets,
   const IndexTableUtils* utils,
-  const size_t maxTrackletsPerCluster = 10)
+  const int rofId,
+  const size_t maxTrackletsPerCluster = 1e2)
 {
   const int phiBins{utils->getNphiBins()};
   const int zBins{utils->getNzBins()};
@@ -180,9 +199,9 @@ GPUg() void trackleterKernel(
             if (o2::gpu::GPUCommonMath::Abs(currentCluster.phi - nextCluster.phi) < phiCut) {
               if (storedTracklets < maxTrackletsPerCluster) {
                 if constexpr (Mode == TrackletMode::Layer0Layer1) {
-                  new (Tracklets + stride + storedTracklets) Tracklet{iNextLayerClusterIndex, iCurrentLayerClusterIndex, nextCluster, currentCluster};
+                  new (Tracklets + stride + storedTracklets) Tracklet{iNextLayerClusterIndex, iCurrentLayerClusterIndex, nextCluster, currentCluster, rofId, rofId};
                 } else {
-                  new (Tracklets + stride + storedTracklets) Tracklet{iCurrentLayerClusterIndex, iNextLayerClusterIndex, currentCluster, nextCluster};
+                  new (Tracklets + stride + storedTracklets) Tracklet{iCurrentLayerClusterIndex, iNextLayerClusterIndex, currentCluster, nextCluster, rofId, rofId};
                 }
                 ++storedTracklets;
               }
@@ -191,6 +210,9 @@ GPUg() void trackleterKernel(
         }
       }
       foundTracklets[iCurrentLayerClusterIndex] = storedTracklets;
+      if (storedTracklets >= maxTrackletsPerCluster) {
+        printf("gpu tracklet finder: some lines will be left behind for cluster %d. valid: %u max: %zu\n", iCurrentLayerClusterIndex, storedTracklets, maxTrackletsPerCluster);
+      }
     }
   }
 }
@@ -208,7 +230,7 @@ GPUg() void trackletSelectionKernel(
   Line* lines,
   int* nFoundLines,
   int* nExclusiveFoundLines,
-  const int maxTrackletsPerCluster = 10,
+  const int maxTrackletsPerCluster = 1e2,
   const float tanLambdaCut = 0.025f,
   const float phiCut = 0.002f)
 {
@@ -237,16 +259,17 @@ GPUg() void trackletSelectionKernel(
   }
 }
 
-GPUg() void computeCentroidsKernel(Line* lines,
-                                   int* nFoundLines,
-                                   int* nExclusiveFoundLines,
-                                   const size_t nClustersMiddleLayer,
-                                   float* centroids,
-                                   const float lowHistX,
-                                   const float highHistX,
-                                   const float lowHistY,
-                                   const float highHistY,
-                                   const float pairCut)
+GPUg() void computeCentroidsKernel(
+  Line* lines,
+  int* nFoundLines,
+  int* nExclusiveFoundLines,
+  const size_t nClustersMiddleLayer,
+  float* centroids,
+  const float lowHistX,
+  const float highHistX,
+  const float lowHistY,
+  const float highHistY,
+  const float pairCut)
 {
   const int nLines = nExclusiveFoundLines[nClustersMiddleLayer - 1] + nFoundLines[nClustersMiddleLayer - 1];
   const int maxIterations{nLines * (nLines - 1) / 2};
@@ -281,23 +304,24 @@ GPUg() void computeCentroidsKernel(Line* lines,
   }
 }
 
-GPUg() void computeZCentroidsKernel(const int nLines,
-                                    const cub::KeyValuePair<int, int>* tmpVtX,
-                                    float* beamPosition,
-                                    Line* lines,
-                                    float* centroids,
-                                    const int* histX, // X
-                                    const float lowHistX,
-                                    const float binSizeHistX,
-                                    const int nBinsHistX,
-                                    const int* histY, // Y
-                                    const float lowHistY,
-                                    const float binSizeHistY,
-                                    const int nBinsHistY,
-                                    const float lowHistZ, // Z
-                                    const float pairCut,
-                                    const int binOpeningX,
-                                    const int binOpeningY)
+GPUg() void computeZCentroidsKernel(
+  const int nLines,
+  const cub::KeyValuePair<int, int>* tmpVtX,
+  float* beamPosition,
+  Line* lines,
+  float* centroids,
+  const int* histX, // X
+  const float lowHistX,
+  const float binSizeHistX,
+  const int nBinsHistX,
+  const int* histY, // Y
+  const float lowHistY,
+  const float binSizeHistY,
+  const int nBinsHistY,
+  const float lowHistZ, // Z
+  const float pairCut,
+  const int binOpeningX,
+  const int binOpeningY)
 {
   for (size_t currentThreadIndex = blockIdx.x * blockDim.x + threadIdx.x; currentThreadIndex < nLines; currentThreadIndex += blockDim.x * gridDim.x) {
     if (tmpVtX[0].value || tmpVtX[1].value) {
@@ -334,16 +358,17 @@ GPUg() void computeZCentroidsKernel(const int nLines,
   }
 }
 
-GPUg() void computeVertexKernel(cub::KeyValuePair<int, int>* tmpVertexBins,
-                                int* histZ, // Z
-                                const float lowHistZ,
-                                const float binSizeHistZ,
-                                const int nBinsHistZ,
-                                Vertex* vertices,
-                                float* beamPosition,
-                                const int vertIndex,
-                                const int minContributors,
-                                const int binOpeningZ)
+GPUg() void computeVertexKernel(
+  cub::KeyValuePair<int, int>* tmpVertexBins,
+  int* histZ, // Z
+  const float lowHistZ,
+  const float binSizeHistZ,
+  const int nBinsHistZ,
+  Vertex* vertices,
+  float* beamPosition,
+  const int vertIndex,
+  const int minContributors,
+  const int binOpeningZ)
 {
   for (size_t currentThreadIndex = blockIdx.x * blockDim.x + threadIdx.x; currentThreadIndex < binOpeningZ; currentThreadIndex += blockDim.x * gridDim.x) {
     if (currentThreadIndex == 0) {
@@ -374,6 +399,15 @@ GPUg() void computeVertexKernel(cub::KeyValuePair<int, int>* tmpVertexBins,
 }
 } // namespace gpu
 
+void VertexerTraitsGPU::updateVertexingParameters(const VertexingParameters& vrtPar)
+{
+  mVrtParams = vrtPar;
+  mIndexTableUtils.setTrackingParameters(vrtPar);
+  mVrtParams.phiSpan = static_cast<int>(std::ceil(mIndexTableUtils.getNphiBins() * mVrtParams.phiCut /
+                                                  constants::math::TwoPi));
+  mVrtParams.zSpan = static_cast<int>(std::ceil(mVrtParams.zCut * mIndexTableUtils.getInverseZCoordinate(0)));
+}
+
 void VertexerTraitsGPU::computeTracklets()
 {
   if (!mTimeFrameGPU->getClusters().size()) {
@@ -382,39 +416,49 @@ void VertexerTraitsGPU::computeTracklets()
   for (int rofId{0}; rofId < mTimeFrameGPU->getNrof(); ++rofId) {
     const dim3 threadsPerBlock{gpu::utils::host::getBlockSize(mTimeFrameGPU->getNClustersLayer(rofId, 1))};
     const dim3 blocksGrid{gpu::utils::host::getBlocksGrid(threadsPerBlock, mTimeFrameGPU->getNClustersLayer(rofId, 1))};
-    gpu::trackleterKernel<TrackletMode::Layer0Layer1><<<blocksGrid, threadsPerBlock>>>(
+    gpu::trackleterKernel<TrackletMode::Layer0Layer1><<<blocksGrid, threadsPerBlock, 0, mTimeFrameGPU->getStream(0).get()>>>(
       mTimeFrameGPU->getDeviceClustersOnLayer(rofId, 0),
       mTimeFrameGPU->getDeviceClustersOnLayer(rofId, 1),
       mTimeFrameGPU->getNClustersLayer(rofId, 0),
       mTimeFrameGPU->getNClustersLayer(rofId, 1),
-      mTimeFrameGPU->getDeviceIndexTableL0(rofId),
+      mTimeFrameGPU->getDeviceIndexTableAtRof(0, rofId),
       mVrtParams.phiCut,
-      mTimeFrameGPU->getDeviceTracklets(rofId, 0),
+      mTimeFrameGPU->getDeviceTrackletsVertexerOnly(rofId, 0),
       mTimeFrameGPU->getDeviceNTrackletsCluster(rofId, 0),
       mDeviceIndexTableUtils,
-      mTimeFrameGPU->getConfig().maxTrackletsPerCluster);
-
-    gpu::trackleterKernel<TrackletMode::Layer1Layer2><<<blocksGrid, threadsPerBlock>>>(
+      rofId,
+      mVrtParams.maxTrackletsPerCluster);
+    gpu::trackleterKernel<TrackletMode::Layer1Layer2><<<blocksGrid, threadsPerBlock, 0, mTimeFrameGPU->getStream(1).get()>>>(
       mTimeFrameGPU->getDeviceClustersOnLayer(rofId, 2),
       mTimeFrameGPU->getDeviceClustersOnLayer(rofId, 1),
       mTimeFrameGPU->getNClustersLayer(rofId, 2),
       mTimeFrameGPU->getNClustersLayer(rofId, 1),
-      mTimeFrameGPU->getDeviceIndexTableL2(rofId),
+      mTimeFrameGPU->getDeviceIndexTableAtRof(2, rofId),
       mVrtParams.phiCut,
-      mTimeFrameGPU->getDeviceTracklets(rofId, 1),
+      mTimeFrameGPU->getDeviceTrackletsVertexerOnly(rofId, 1),
       mTimeFrameGPU->getDeviceNTrackletsCluster(rofId, 1),
       mDeviceIndexTableUtils,
-      mTimeFrameGPU->getConfig().maxTrackletsPerCluster);
+      rofId,
+      mVrtParams.maxTrackletsPerCluster);
   }
 #ifdef VTX_DEBUG
+  discardResult(cudaDeviceSynchronize());
   std::ofstream out01("NTC01.txt"), out12("NTC12.txt");
   std::vector<std::vector<int>> NtrackletsClusters01(mTimeFrameGPU->getNrof());
   std::vector<std::vector<int>> NtrackletsClusters12(mTimeFrameGPU->getNrof());
   for (int iRof{0}; iRof < mTimeFrameGPU->getNrof(); ++iRof) {
-    NtrackletsClusters01[iRof].resize(mTimeFrameGPU->getClustersOnLayer(iRof, 1).size());
-    NtrackletsClusters12[iRof].resize(mTimeFrameGPU->getClustersOnLayer(iRof, 1).size());
-    checkGPUError(cudaMemcpy(NtrackletsClusters01[iRof].data(), mTimeFrameGPU->getDeviceNTrackletsCluster(iRof, 0), sizeof(int) * mTimeFrameGPU->getClustersOnLayer(iRof, 1).size(), cudaMemcpyDeviceToHost), __FILE__, __LINE__);
-    checkGPUError(cudaMemcpy(NtrackletsClusters12[iRof].data(), mTimeFrameGPU->getDeviceNTrackletsCluster(iRof, 1), sizeof(int) * mTimeFrameGPU->getClustersOnLayer(iRof, 1).size(), cudaMemcpyDeviceToHost), __FILE__, __LINE__);
+    NtrackletsClusters01[iRof].resize(mTimeFrameGPU->getNClustersLayer(iRof, 1));
+    NtrackletsClusters12[iRof].resize(mTimeFrameGPU->getNClustersLayer(iRof, 1));
+    checkGPUError(cudaMemcpy(NtrackletsClusters01[iRof].data(),
+                             mTimeFrameGPU->getDeviceNTrackletsCluster(iRof, 0),
+                             sizeof(int) * mTimeFrameGPU->getNClustersLayer(iRof, 1),
+                             cudaMemcpyDeviceToHost),
+                  __FILE__, __LINE__);
+    checkGPUError(cudaMemcpy(NtrackletsClusters12[iRof].data(),
+                             mTimeFrameGPU->getDeviceNTrackletsCluster(iRof, 1),
+                             sizeof(int) * mTimeFrameGPU->getNClustersLayer(iRof, 1),
+                             cudaMemcpyDeviceToHost),
+                  __FILE__, __LINE__);
 
     std::copy(NtrackletsClusters01[iRof].begin(), NtrackletsClusters01[iRof].end(), std::ostream_iterator<double>(out01, "\t"));
     std::copy(NtrackletsClusters12[iRof].begin(), NtrackletsClusters12[iRof].end(), std::ostream_iterator<double>(out12, "\t"));
@@ -426,37 +470,33 @@ void VertexerTraitsGPU::computeTracklets()
   // Dump lines on root file
   TFile* trackletFile = TFile::Open("artefacts_tf_gpu.root", "recreate");
   TTree* tr_tre = new TTree("tracklets", "tf");
-  std::vector<o2::its::Tracklet> tracklets_vec01(0);
-  std::vector<o2::its::Tracklet> tracklets_vec12(0);
-  std::vector<o2::its::Tracklet> tracklets_clean01(0);
-  std::vector<o2::its::Tracklet> tracklets_clean12(0);
-  tr_tre->Branch("Tracklets0", &tracklets_clean01);
-  tr_tre->Branch("Tracklets1", &tracklets_clean12);
+  std::vector<o2::its::Tracklet> tracklets01;
+  std::vector<o2::its::Tracklet> tracklets12;
+  tr_tre->Branch("Tracklets0", &tracklets01);
+  tr_tre->Branch("Tracklets1", &tracklets12);
   for (int iRof{0}; iRof < mTimeFrameGPU->getNrof(); ++iRof) {
-    tracklets_clean01.clear();
-    tracklets_clean12.clear();
-    tracklets_vec01.resize(mTimeFrameGPU->getClustersOnLayer(iRof, 1).size() * mTimeFrameGPU->getConfig().maxTrackletsPerCluster); // Nclusters * TrackletsPerCluster
-    tracklets_vec12.resize(mTimeFrameGPU->getClustersOnLayer(iRof, 1).size() * mTimeFrameGPU->getConfig().maxTrackletsPerCluster);
-    checkGPUError(cudaMemcpy(tracklets_vec01.data(),
-                             mTimeFrameGPU->getDeviceTracklets(iRof, 0),
-                             mTimeFrameGPU->getClustersOnLayer(iRof, 1).size() * sizeof(Tracklet) * mTimeFrameGPU->getConfig().maxTrackletsPerCluster,
+    tracklets01.clear();
+    tracklets12.clear();
+    std::vector<o2::its::Tracklet> rawTracklets0(mTimeFrameGPU->getNClustersLayer(iRof, 1) * mVrtParams.maxTrackletsPerCluster);
+    std::vector<o2::its::Tracklet> rawTracklets1(mTimeFrameGPU->getNClustersLayer(iRof, 1) * mVrtParams.maxTrackletsPerCluster);
+    checkGPUError(cudaMemcpy(rawTracklets0.data(),
+                             mTimeFrameGPU->getDeviceTrackletsVertexerOnly(iRof, 0),
+                             mTimeFrameGPU->getNClustersLayer(iRof, 1) * mVrtParams.maxTrackletsPerCluster * sizeof(o2::its::Tracklet),
                              cudaMemcpyDeviceToHost),
                   __FILE__, __LINE__);
-    checkGPUError(cudaMemcpy(tracklets_vec12.data(),
-                             mTimeFrameGPU->getDeviceTracklets(iRof, 1),
-                             mTimeFrameGPU->getClustersOnLayer(iRof, 1).size() * sizeof(Tracklet) * mTimeFrameGPU->getConfig().maxTrackletsPerCluster,
+    checkGPUError(cudaMemcpy(rawTracklets1.data(),
+                             mTimeFrameGPU->getDeviceTrackletsVertexerOnly(iRof, 1),
+                             mTimeFrameGPU->getNClustersLayer(iRof, 1) * mVrtParams.maxTrackletsPerCluster * sizeof(o2::its::Tracklet),
                              cudaMemcpyDeviceToHost),
                   __FILE__, __LINE__);
-    for (auto iCluster{0}; iCluster < NtrackletsClusters01[iRof].size(); ++iCluster) {
-      auto nTracklets{NtrackletsClusters01[iRof][iCluster]};
-      for (auto iTracklet{0}; iTracklet < nTracklets; ++iTracklet) {
-        tracklets_clean01.push_back(tracklets_vec01[iCluster * mTimeFrameGPU->getConfig().maxTrackletsPerCluster + iTracklet]);
+    for (int iCluster{0}; iCluster < mTimeFrameGPU->getNClustersLayer(iRof, 1); ++iCluster) {
+      for (int iTracklet{0}; iTracklet < NtrackletsClusters01[iRof][iCluster]; ++iTracklet) {
+        tracklets01.push_back(rawTracklets0[iCluster * mVrtParams.maxTrackletsPerCluster + iTracklet]);
       }
     }
-    for (auto iCluster{0}; iCluster < NtrackletsClusters12[iRof].size(); ++iCluster) {
-      auto nTracklets{NtrackletsClusters12[iRof][iCluster]};
-      for (auto iTracklet{0}; iTracklet < nTracklets; ++iTracklet) {
-        tracklets_clean12.push_back(tracklets_vec12[iCluster * mTimeFrameGPU->getConfig().maxTrackletsPerCluster + iTracklet]);
+    for (int iCluster{0}; iCluster < mTimeFrameGPU->getNClustersLayer(iRof, 1); ++iCluster) {
+      for (int iTracklet{0}; iTracklet < NtrackletsClusters12[iRof][iCluster]; ++iTracklet) {
+        tracklets12.push_back(rawTracklets1[iCluster * mVrtParams.maxTrackletsPerCluster + iTracklet]);
       }
     }
     tr_tre->Fill();
@@ -479,21 +519,21 @@ void VertexerTraitsGPU::computeTrackletMatching()
     size_t bufferSize = mTimeFrameGPU->getConfig().tmpCUBBufferSize;
 
     // Reset used tracklets
-    checkGPUError(cudaMemset(mTimeFrameGPU->getDeviceUsedTracklets(rofId), false, sizeof(unsigned char) * mTimeFrameGPU->getConfig().maxTrackletsPerCluster * mTimeFrameGPU->getNClustersLayer(rofId, 1)), __FILE__, __LINE__);
+    checkGPUError(cudaMemset(mTimeFrameGPU->getDeviceUsedTracklets(rofId), false, sizeof(unsigned char) * mVrtParams.maxTrackletsPerCluster * mTimeFrameGPU->getNClustersLayer(rofId, 1)), __FILE__, __LINE__);
 
     gpu::trackletSelectionKernel<true><<<blocksGrid, threadsPerBlock>>>(
       mTimeFrameGPU->getDeviceClustersOnLayer(rofId, 0),
       mTimeFrameGPU->getDeviceClustersOnLayer(rofId, 1),
       mTimeFrameGPU->getNClustersLayer(rofId, 1),
-      mTimeFrameGPU->getDeviceTracklets(rofId, 0),
-      mTimeFrameGPU->getDeviceTracklets(rofId, 1),
+      mTimeFrameGPU->getDeviceTrackletsVertexerOnly(rofId, 0),
+      mTimeFrameGPU->getDeviceTrackletsVertexerOnly(rofId, 1),
       mTimeFrameGPU->getDeviceNTrackletsCluster(rofId, 0),
       mTimeFrameGPU->getDeviceNTrackletsCluster(rofId, 1),
       mTimeFrameGPU->getDeviceUsedTracklets(rofId),
       mTimeFrameGPU->getDeviceLines(rofId),
       mTimeFrameGPU->getDeviceNFoundLines(rofId),
       mTimeFrameGPU->getDeviceExclusiveNFoundLines(rofId),
-      mTimeFrameGPU->getConfig().maxTrackletsPerCluster,
+      mVrtParams.maxTrackletsPerCluster,
       mVrtParams.tanLambdaCut,
       mVrtParams.phiCut);
 
@@ -503,26 +543,36 @@ void VertexerTraitsGPU::computeTrackletMatching()
                                                 mTimeFrameGPU->getDeviceExclusiveNFoundLines(rofId),
                                                 mTimeFrameGPU->getNClustersLayer(rofId, 1)));
     // Reset used tracklets
-    checkGPUError(cudaMemset(mTimeFrameGPU->getDeviceUsedTracklets(rofId), false, sizeof(unsigned char) * mTimeFrameGPU->getConfig().maxTrackletsPerCluster * mTimeFrameGPU->getNClustersLayer(rofId, 1)), __FILE__, __LINE__);
+    checkGPUError(cudaMemset(mTimeFrameGPU->getDeviceUsedTracklets(rofId), false, sizeof(unsigned char) * mVrtParams.maxTrackletsPerCluster * mTimeFrameGPU->getNClustersLayer(rofId, 1)), __FILE__, __LINE__);
 
     gpu::trackletSelectionKernel<false><<<blocksGrid, threadsPerBlock>>>(
       mTimeFrameGPU->getDeviceClustersOnLayer(rofId, 0),
       mTimeFrameGPU->getDeviceClustersOnLayer(rofId, 1),
       mTimeFrameGPU->getNClustersLayer(rofId, 1),
-      mTimeFrameGPU->getDeviceTracklets(rofId, 0),
-      mTimeFrameGPU->getDeviceTracklets(rofId, 1),
+      mTimeFrameGPU->getDeviceTrackletsVertexerOnly(rofId, 0),
+      mTimeFrameGPU->getDeviceTrackletsVertexerOnly(rofId, 1),
       mTimeFrameGPU->getDeviceNTrackletsCluster(rofId, 0),
       mTimeFrameGPU->getDeviceNTrackletsCluster(rofId, 1),
       mTimeFrameGPU->getDeviceUsedTracklets(rofId),
       mTimeFrameGPU->getDeviceLines(rofId),
       mTimeFrameGPU->getDeviceNFoundLines(rofId),
       mTimeFrameGPU->getDeviceExclusiveNFoundLines(rofId),
-      mTimeFrameGPU->getConfig().maxTrackletsPerCluster,
+      mVrtParams.maxTrackletsPerCluster,
       mVrtParams.tanLambdaCut,
       mVrtParams.phiCut);
     gpuThrowOnError();
+
+    // Code to be removed in case of GPU vertex computing
+    int excLas, nLas, nLines;
+    checkGPUError(cudaMemcpy(&excLas, mTimeFrameGPU->getDeviceExclusiveNFoundLines(rofId) + mTimeFrameGPU->getNClustersLayer(rofId, 1) - 1, sizeof(int), cudaMemcpyDeviceToHost), __FILE__, __LINE__);
+    checkGPUError(cudaMemcpy(&nLas, mTimeFrameGPU->getDeviceNFoundLines(rofId) + mTimeFrameGPU->getNClustersLayer(rofId, 1) - 1, sizeof(int), cudaMemcpyDeviceToHost), __FILE__, __LINE__);
+    nLines = excLas + nLas;
+    mTimeFrameGPU->getLines(rofId).resize(nLines);
+    checkGPUError(cudaMemcpy(mTimeFrameGPU->getLines(rofId).data(), mTimeFrameGPU->getDeviceLines(rofId), sizeof(Line) * nLines, cudaMemcpyDeviceToHost), __FILE__, __LINE__);
   }
+
 #ifdef VTX_DEBUG
+  discardResult(cudaDeviceSynchronize());
   std::vector<std::vector<int>> NFoundLines(mTimeFrameGPU->getNrof()), ExcNFoundLines(mTimeFrameGPU->getNrof());
   std::ofstream nlines_out("N_lines_gpu.txt");
   for (size_t rofId{0}; rofId < mTimeFrameGPU->getNrof(); ++rofId) {
@@ -557,6 +607,125 @@ void VertexerTraitsGPU::computeTrackletMatching()
 }
 
 void VertexerTraitsGPU::computeVertices()
+{
+#ifdef VTX_DEBUG
+  std::vector<std::vector<ClusterLines>> dbg_clusLines(mTimeFrameGPU->getNrof());
+#endif
+  std::vector<int> noClustersVec(mTimeFrameGPU->getNrof(), 0);
+  for (int rofId{0}; rofId < mTimeFrameGPU->getNrof(); ++rofId) {
+    const int numTracklets{static_cast<int>(mTimeFrameGPU->getLines(rofId).size())};
+    std::vector<bool> usedTracklets(numTracklets, false);
+    for (int tracklet1{0}; tracklet1 < numTracklets; ++tracklet1) {
+      if (usedTracklets[tracklet1]) {
+        continue;
+      }
+      for (int tracklet2{tracklet1 + 1}; tracklet2 < numTracklets; ++tracklet2) {
+        if (usedTracklets[tracklet2]) {
+          continue;
+        }
+        if (Line::getDCA(mTimeFrameGPU->getLines(rofId)[tracklet1], mTimeFrameGPU->getLines(rofId)[tracklet2]) < mVrtParams.pairCut) {
+          mTimeFrameGPU->getTrackletClusters(rofId).emplace_back(tracklet1, mTimeFrameGPU->getLines(rofId)[tracklet1], tracklet2, mTimeFrameGPU->getLines(rofId)[tracklet2]);
+          std::array<float, 3> tmpVertex{mTimeFrameGPU->getTrackletClusters(rofId).back().getVertex()};
+          if (tmpVertex[0] * tmpVertex[0] + tmpVertex[1] * tmpVertex[1] > 4.f) {
+            mTimeFrameGPU->getTrackletClusters(rofId).pop_back();
+            break;
+          }
+          usedTracklets[tracklet1] = true;
+          usedTracklets[tracklet2] = true;
+          for (int tracklet3{0}; tracklet3 < numTracklets; ++tracklet3) {
+            if (usedTracklets[tracklet3]) {
+              continue;
+            }
+            if (Line::getDistanceFromPoint(mTimeFrameGPU->getLines(rofId)[tracklet3], tmpVertex) < mVrtParams.pairCut) {
+              mTimeFrameGPU->getTrackletClusters(rofId).back().add(tracklet3, mTimeFrameGPU->getLines(rofId)[tracklet3]);
+              usedTracklets[tracklet3] = true;
+              tmpVertex = mTimeFrameGPU->getTrackletClusters(rofId).back().getVertex();
+            }
+          }
+          break;
+        }
+      }
+    }
+    std::sort(mTimeFrameGPU->getTrackletClusters(rofId).begin(), mTimeFrameGPU->getTrackletClusters(rofId).end(),
+              [](ClusterLines& cluster1, ClusterLines& cluster2) { return cluster1.getSize() > cluster2.getSize(); });
+    noClustersVec[rofId] = static_cast<int>(mTimeFrameGPU->getTrackletClusters(rofId).size());
+    for (int iCluster1{0}; iCluster1 < noClustersVec[rofId]; ++iCluster1) {
+      std::array<float, 3> vertex1{mTimeFrameGPU->getTrackletClusters(rofId)[iCluster1].getVertex()};
+      std::array<float, 3> vertex2{};
+      for (int iCluster2{iCluster1 + 1}; iCluster2 < noClustersVec[rofId]; ++iCluster2) {
+        vertex2 = mTimeFrameGPU->getTrackletClusters(rofId)[iCluster2].getVertex();
+        if (std::abs(vertex1[2] - vertex2[2]) < mVrtParams.clusterCut) {
+          float distance{(vertex1[0] - vertex2[0]) * (vertex1[0] - vertex2[0]) +
+                         (vertex1[1] - vertex2[1]) * (vertex1[1] - vertex2[1]) +
+                         (vertex1[2] - vertex2[2]) * (vertex1[2] - vertex2[2])};
+          if (distance < mVrtParams.pairCut * mVrtParams.pairCut) {
+            for (auto label : mTimeFrameGPU->getTrackletClusters(rofId)[iCluster2].getLabels()) {
+              mTimeFrameGPU->getTrackletClusters(rofId)[iCluster1].add(label, mTimeFrameGPU->getLines(rofId)[label]);
+              vertex1 = mTimeFrameGPU->getTrackletClusters(rofId)[iCluster1].getVertex();
+            }
+          }
+          mTimeFrameGPU->getTrackletClusters(rofId).erase(mTimeFrameGPU->getTrackletClusters(rofId).begin() + iCluster2);
+          --iCluster2;
+          --noClustersVec[rofId];
+        }
+      }
+    }
+  }
+
+  for (int rofId{0}; rofId < mTimeFrameGPU->getNrof(); ++rofId) {
+#ifdef VTX_DEBUG
+    for (auto& cl : mTimeFrameGPU->getTrackletClusters(rofId)) {
+      dbg_clusLines[rofId].push_back(cl);
+    }
+#endif
+    for (int iCluster{0}; iCluster < noClustersVec[rofId]; ++iCluster) {
+      if (mTimeFrameGPU->getTrackletClusters(rofId)[iCluster].getSize() < mVrtParams.clusterContributorsCut && noClustersVec[rofId] > 1) {
+        mTimeFrameGPU->getTrackletClusters(rofId).erase(mTimeFrameGPU->getTrackletClusters(rofId).begin() + iCluster);
+        noClustersVec[rofId]--;
+        continue;
+      }
+      if (mTimeFrameGPU->getTrackletClusters(rofId)[iCluster].getVertex()[0] * mTimeFrameGPU->getTrackletClusters(rofId)[iCluster].getVertex()[0] +
+            mTimeFrameGPU->getTrackletClusters(rofId)[iCluster].getVertex()[1] * mTimeFrameGPU->getTrackletClusters(rofId)[iCluster].getVertex()[1] <
+          1.98 * 1.98) {
+        mVertices.emplace_back(mTimeFrameGPU->getTrackletClusters(rofId)[iCluster].getVertex()[0],
+                               mTimeFrameGPU->getTrackletClusters(rofId)[iCluster].getVertex()[1],
+                               mTimeFrameGPU->getTrackletClusters(rofId)[iCluster].getVertex()[2],
+                               mTimeFrameGPU->getTrackletClusters(rofId)[iCluster].getRMS2(),         // Symm matrix. Diagonal: RMS2 components,
+                                                                                                      // off-diagonal: square mean of projections on planes.
+                               mTimeFrameGPU->getTrackletClusters(rofId)[iCluster].getSize(),         // Contributors
+                               mTimeFrameGPU->getTrackletClusters(rofId)[iCluster].getAvgDistance2(), // In place of chi2
+                               rofId);
+      }
+    }
+
+    mTimeFrameGPU->addPrimaryVertices(mVertices);
+    mVertices.clear();
+  }
+#ifdef VTX_DEBUG
+  TFile* dbg_file = TFile::Open("artefacts_tf.root", "update");
+  TTree* ln_clus_lines_tree = new TTree("clusterlines", "tf");
+  std::vector<o2::its::ClusterLines> cl_lines_vec_pre(0);
+  std::vector<o2::its::ClusterLines> cl_lines_vec_post(0);
+  ln_clus_lines_tree->Branch("cllines_pre", &cl_lines_vec_pre);
+  ln_clus_lines_tree->Branch("cllines_post", &cl_lines_vec_post);
+  for (auto rofId{0}; rofId < mTimeFrameGPU->getNrof(); ++rofId) {
+    cl_lines_vec_pre.clear();
+    cl_lines_vec_post.clear();
+    for (auto& clln : mTimeFrameGPU->getTrackletClusters(rofId)) {
+      cl_lines_vec_post.push_back(clln);
+    }
+    for (auto& cl : dbg_clusLines[rofId]) {
+      cl_lines_vec_pre.push_back(cl);
+    }
+    ln_clus_lines_tree->Fill();
+  }
+  dbg_file->cd();
+  ln_clus_lines_tree->Write();
+  dbg_file->Close();
+#endif
+}
+
+void VertexerTraitsGPU::computeVerticesHist()
 {
   if (!mTimeFrameGPU->getClusters().size()) {
     return;
diff --git a/Detectors/ITSMFT/ITS/tracking/GPU/hip/CMakeLists.txt b/Detectors/ITSMFT/ITS/tracking/GPU/hip/CMakeLists.txt
index 3c1c83f71a9b2..7f9643b89ccb4 100644
--- a/Detectors/ITSMFT/ITS/tracking/GPU/hip/CMakeLists.txt
+++ b/Detectors/ITSMFT/ITS/tracking/GPU/hip/CMakeLists.txt
@@ -31,7 +31,7 @@ if(HIP_ENABLED)
   set(CMAKE_CXX_EXTENSIONS OFF)
   set(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} ${O2_HIP_CMAKE_CXX_FLAGS} -fgpu-rdc")
 
-  message(STATUS "Building ITS HIP vertexer")
+  message(STATUS "Building ITS HIP tracker")
   o2_add_library(ITStrackingHIP
                  SOURCES ClusterLinesGPU.hip.cxx
                          Context.hip.cxx
diff --git a/Detectors/ITSMFT/ITS/tracking/include/ITStracking/Cell.h b/Detectors/ITSMFT/ITS/tracking/include/ITStracking/Cell.h
index bcb449e261ae0..adf5bf095143c 100644
--- a/Detectors/ITSMFT/ITS/tracking/include/ITStracking/Cell.h
+++ b/Detectors/ITSMFT/ITS/tracking/include/ITStracking/Cell.h
@@ -31,16 +31,17 @@ namespace its
 class Cell final
 {
  public:
+  GPUhd() Cell();
   GPUd() Cell(const int, const int, const int, const int, const int, const float);
 
-  GPUhdni() int getFirstClusterIndex() const { return mFirstClusterIndex; };
-  GPUhdni() int getSecondClusterIndex() const { return mSecondClusterIndex; };
-  GPUhdni() int getThirdClusterIndex() const { return mThirdClusterIndex; };
-  GPUhdni() int getFirstTrackletIndex() const { return mFirstTrackletIndex; };
-  GPUhdni() int getSecondTrackletIndex() const { return mSecondTrackletIndex; };
-  GPUhdni() float getTanLambda() const { return mTanLambda; };
-  int getLevel() const { return mLevel; };
-  void setLevel(const int level) { mLevel = level; };
+  GPUhd() int getFirstClusterIndex() const { return mFirstClusterIndex; };
+  GPUhd() int getSecondClusterIndex() const { return mSecondClusterIndex; };
+  GPUhd() int getThirdClusterIndex() const { return mThirdClusterIndex; };
+  GPUhd() int getFirstTrackletIndex() const { return mFirstTrackletIndex; };
+  GPUhd() int getSecondTrackletIndex() const { return mSecondTrackletIndex; };
+  GPUhd() float getTanLambda() const { return mTanLambda; };
+  GPUhd() int getLevel() const { return mLevel; };
+  GPUhd() void setLevel(const int level) { mLevel = level; };
 
  private:
   const int mFirstClusterIndex;
@@ -52,13 +53,25 @@ class Cell final
   int mLevel;
 };
 
+GPUhdi() Cell::Cell()
+  : mFirstClusterIndex{0},
+    mSecondClusterIndex{0},
+    mThirdClusterIndex{0},
+    mFirstTrackletIndex{0},
+    mSecondTrackletIndex{0},
+    mTanLambda{0},
+    mLevel{0}
+{
+  // Nothing to do
+}
+
 GPUdi() Cell::Cell(const int firstClusterIndex, const int secondClusterIndex, const int thirdClusterIndex,
                    const int firstTrackletIndex, const int secondTrackletIndex, const float tanL)
   : mFirstClusterIndex{firstClusterIndex},
     mSecondClusterIndex{secondClusterIndex},
     mThirdClusterIndex{thirdClusterIndex},
-    mFirstTrackletIndex(firstTrackletIndex),
-    mSecondTrackletIndex(secondTrackletIndex),
+    mFirstTrackletIndex{firstTrackletIndex},
+    mSecondTrackletIndex{secondTrackletIndex},
     mTanLambda{tanL},
     mLevel{1}
 {
diff --git a/Detectors/ITSMFT/ITS/tracking/include/ITStracking/Configuration.h b/Detectors/ITSMFT/ITS/tracking/include/ITStracking/Configuration.h
index 96cdd862ec31f..681647a35bd7c 100644
--- a/Detectors/ITSMFT/ITS/tracking/include/ITStracking/Configuration.h
+++ b/Detectors/ITSMFT/ITS/tracking/include/ITStracking/Configuration.h
@@ -147,13 +147,16 @@ struct TimeFrameGPUConfig {
                      size_t cluLayCap,
                      size_t cluROfCap,
                      size_t maxTrkCap,
-                     size_t maxVertCap);
+                     size_t maxVertCap,
+                     size_t maxROFs);
 
   size_t tmpCUBBufferSize = 1e5; // In average in pp events there are required 4096 bytes
-  size_t maxTrackletsPerCluster = 50;
-  size_t clustersPerLayerCapacity = 5e5;
+  size_t maxTrackletsPerCluster = 1e2;
+  size_t clustersPerLayerCapacity = 2.5e5;
   size_t clustersPerROfCapacity = 1e4;
   size_t trackletsCapacity = maxTrackletsPerCluster * clustersPerLayerCapacity;
+  size_t validatedTrackletsCapacity = 1e5;
+  size_t cellsLUTsize = validatedTrackletsCapacity;
   size_t maxLinesCapacity = 1e2;
   size_t maxCentroidsXYCapacity = std::ceil(maxLinesCapacity * (maxLinesCapacity - 1) / (float)2);
   size_t maxVerticesCapacity = 10;
@@ -167,12 +170,14 @@ inline TimeFrameGPUConfig::TimeFrameGPUConfig(size_t cubBufferSize,
                                               size_t cluLayCap,
                                               size_t cluROfCap,
                                               size_t maxTrkCap,
-                                              size_t maxVertCap) : tmpCUBBufferSize{cubBufferSize},
-                                                                   maxTrackletsPerCluster{maxTrkClu},
-                                                                   clustersPerLayerCapacity{cluLayCap},
-                                                                   clustersPerROfCapacity{cluROfCap},
-                                                                   maxLinesCapacity{maxTrkCap},
-                                                                   maxVerticesCapacity{maxVertCap}
+                                              size_t maxVertCap,
+                                              size_t maxROFs) : tmpCUBBufferSize{cubBufferSize},
+                                                                maxTrackletsPerCluster{maxTrkClu},
+                                                                clustersPerLayerCapacity{cluLayCap},
+                                                                clustersPerROfCapacity{cluROfCap},
+                                                                maxLinesCapacity{maxTrkCap},
+                                                                maxVerticesCapacity{maxVertCap},
+                                                                nMaxROFs{maxROFs}
 {
   maxCentroidsXYCapacity = std::ceil(maxLinesCapacity * (maxLinesCapacity - 1) / 2);
   trackletsCapacity = maxTrackletsPerCluster * clustersPerLayerCapacity;
diff --git a/Detectors/ITSMFT/ITS/tracking/include/ITStracking/Definitions.h b/Detectors/ITSMFT/ITS/tracking/include/ITStracking/Definitions.h
index a6506f03b82bb..de67d50344c62 100644
--- a/Detectors/ITSMFT/ITS/tracking/include/ITStracking/Definitions.h
+++ b/Detectors/ITSMFT/ITS/tracking/include/ITStracking/Definitions.h
@@ -16,6 +16,7 @@
 #define TRACKINGITS_DEFINITIONS_H_
 
 // #define CA_DEBUG
+// #define VTX_DEBUG
 
 template <typename T>
 void discardResult(const T&)
diff --git a/Detectors/ITSMFT/ITS/tracking/include/ITStracking/IndexTableUtils.h b/Detectors/ITSMFT/ITS/tracking/include/ITStracking/IndexTableUtils.h
index e91deb6d72127..66de759f943e3 100644
--- a/Detectors/ITSMFT/ITS/tracking/include/ITStracking/IndexTableUtils.h
+++ b/Detectors/ITSMFT/ITS/tracking/include/ITStracking/IndexTableUtils.h
@@ -16,18 +16,11 @@
 #ifndef TRACKINGITSU_INCLUDE_INDEXTABLEUTILS_H_
 #define TRACKINGITSU_INCLUDE_INDEXTABLEUTILS_H_
 
-#ifndef GPUCA_GPUCODE_DEVICE
-#include <array>
-#include <utility>
-#include <vector>
-#endif
-
 #include "ITStracking/Constants.h"
 #include "ITStracking/Configuration.h"
 #include "ITStracking/Definitions.h"
 #include "GPUCommonMath.h"
 #include "GPUCommonDef.h"
-#include "GPUCommonLogger.h"
 
 namespace o2
 {
@@ -52,8 +45,8 @@ class IndexTableUtils
   int mNzBins = 0;
   int mNphiBins = 0;
   float mInversePhiBinSize = 0.f;
-  float mLayerZ[7];
-  float mInverseZBinSize[7];
+  float mLayerZ[7] = {0.f, 0.f, 0.f, 0.f, 0.f, 0.f, 0.f};
+  float mInverseZBinSize[7] = {0.f, 0.f, 0.f, 0.f, 0.f, 0.f, 0.f};
 };
 
 template <class T>
diff --git a/Detectors/ITSMFT/ITS/tracking/include/ITStracking/TimeFrame.h b/Detectors/ITSMFT/ITS/tracking/include/ITStracking/TimeFrame.h
index 729af6d8b5924..c24580d2a681a 100644
--- a/Detectors/ITSMFT/ITS/tracking/include/ITStracking/TimeFrame.h
+++ b/Detectors/ITSMFT/ITS/tracking/include/ITStracking/TimeFrame.h
@@ -58,6 +58,17 @@ namespace its
 {
 using Vertex = o2::dataformats::Vertex<o2::dataformats::TimeStamp<int>>;
 
+struct lightVertex {
+  lightVertex(float x, float y, float z, std::array<float, 6> rms2, int cont, float avgdis2, int stamp);
+  float mX;
+  float mY;
+  float mZ;
+  std::array<float, 6> mRMS2;
+  float mAvgDistance2;
+  int mContributors;
+  int mTimeStamp;
+};
+
 class TimeFrame
 {
  public:
@@ -70,6 +81,7 @@ class TimeFrame
   int getPrimaryVerticesNum(int rofID = -1) const;
   void addPrimaryVertices(const std::vector<Vertex>& vertices);
   void addPrimaryVertices(const gsl::span<const Vertex>& vertices);
+  void addPrimaryVertices(const std::vector<lightVertex>&);
   void removePrimaryVerticesInROf(const int rofId);
   int loadROFrameData(const o2::itsmft::ROFRecord& rof, gsl::span<const itsmft::Cluster> clusters,
                       const dataformats::MCTruthContainer<MCCompLabel>* mcLabels = nullptr);
@@ -84,6 +96,7 @@ class TimeFrame
   bool empty() const;
 
   int getSortedIndex(int rof, int layer, int i) const;
+  int getSortedStartIndex(const int, const int) const;
   int getNrof() const;
 
   void resetBeamXY(const float x, const float y, const float w = 0);
@@ -172,11 +185,14 @@ class TimeFrame
 
   IndexTableUtils mIndexTableUtils;
 
+  bool mIsGPU = false;
   std::vector<std::vector<Cluster>> mClusters;
   std::vector<std::vector<TrackingFrameInfo>> mTrackingFrameInfo;
   std::vector<std::vector<int>> mClusterExternalIndices;
   std::vector<std::vector<int>> mROframesClusters;
   std::vector<std::vector<int>> mIndexTables;
+  std::vector<std::vector<int>> mTrackletsLookupTable;
+  std::vector<std::vector<unsigned char>> mUsedClusters;
   int mNrof = 0;
 
  private:
@@ -199,7 +215,6 @@ class TimeFrame
   std::vector<Vertex> mPrimaryVertices;
   std::vector<std::array<float, 2>> mPValphaX; /// PV x and alpha for track propagation
   std::vector<std::vector<Cluster>> mUnsortedClusters;
-  std::vector<std::vector<bool>> mUsedClusters;
   const dataformats::MCTruthContainer<MCCompLabel>* mClusterLabels = nullptr;
   std::vector<std::vector<MCCompLabel>> mTrackletLabels;
   std::vector<std::vector<MCCompLabel>> mCellLabels;
@@ -212,7 +227,6 @@ class TimeFrame
   std::vector<int> mBogusClusters; /// keep track of clusters with wild coordinates
 
   std::vector<std::vector<Tracklet>> mTracklets;
-  std::vector<std::vector<int>> mTrackletsLookupTable;
 
   std::vector<std::pair<unsigned long long, bool>> mRoadLabels;
   int mCutClusterMult;
@@ -259,6 +273,8 @@ inline bool TimeFrame::empty() const { return getTotalClusters() == 0; }
 
 inline int TimeFrame::getSortedIndex(int rof, int layer, int index) const { return mROframesClusters[layer][rof] + index; }
 
+inline int TimeFrame::getSortedStartIndex(const int rof, const int layer) const { return mROframesClusters[layer][rof]; }
+
 inline int TimeFrame::getNrof() const { return mNrof; }
 
 inline void TimeFrame::resetBeamXY(const float x, const float y, const float w)
diff --git a/Detectors/ITSMFT/ITS/tracking/include/ITStracking/Tracker.h b/Detectors/ITSMFT/ITS/tracking/include/ITStracking/Tracker.h
index 067de88839170..7f08904c2d231 100644
--- a/Detectors/ITSMFT/ITS/tracking/include/ITStracking/Tracker.h
+++ b/Detectors/ITSMFT/ITS/tracking/include/ITStracking/Tracker.h
@@ -105,11 +105,6 @@ inline void Tracker::setParameters(const std::vector<TrackingParameters>& trkPar
   mTrkParams = trkPars;
 }
 
-inline void Tracker::initialiseTimeFrame(int& iteration)
-{
-  mTimeFrame->initialise(iteration, mTrkParams[iteration]);
-}
-
 template <typename... T>
 float Tracker::evaluateTask(void (Tracker::*task)(T...), const char* taskName, std::function<void(std::string s)> logger,
                             T&&... args)
diff --git a/Detectors/ITSMFT/ITS/tracking/include/ITStracking/TrackerTraits.h b/Detectors/ITSMFT/ITS/tracking/include/ITStracking/TrackerTraits.h
index 9605a4d193de2..fb0dfa198bbd9 100644
--- a/Detectors/ITSMFT/ITS/tracking/include/ITStracking/TrackerTraits.h
+++ b/Detectors/ITSMFT/ITS/tracking/include/ITStracking/TrackerTraits.h
@@ -52,6 +52,8 @@ class TrackerTraits
 {
  public:
   virtual ~TrackerTraits() = default;
+  virtual void adoptTimeFrame(TimeFrame* tf);
+  virtual void initialiseTimeFrame(const int iteration);
   virtual void computeLayerTracklets(const int iteration);
   virtual void computeLayerCells(const int iteration);
   virtual void findCellsNeighbours(const int iteration);
@@ -60,12 +62,13 @@ class TrackerTraits
   virtual void extendTracks(const int iteration);
   virtual void findShortPrimaries();
   virtual void refitTracks(const int iteration, const std::vector<std::vector<TrackingFrameInfo>>&, std::vector<TrackITSExt>&);
+  virtual void setBz(float bz);
   virtual bool trackFollowing(TrackITSExt* track, int rof, bool outward, const int iteration);
 
   void UpdateTrackingParameters(const std::vector<TrackingParameters>& trkPars);
   TimeFrame* getTimeFrame() { return mTimeFrame; }
-  void adoptTimeFrame(TimeFrame* tf) { mTimeFrame = tf; }
-  void setBz(float bz);
+
+  void setIsGPU(const unsigned char isgpu) { mIsGPU = isgpu; };
   float getBz() const;
   void setCorrType(const o2::base::PropagatorImpl<float>::MatCorrType type) { mCorrType = type; }
   bool isMatLUT() const;
@@ -85,6 +88,13 @@ class TrackerTraits
   void setNThreads(int n);
   int getNThreads() const { return mNThreads; }
 
+  // TimeFrame information forwarding
+  virtual int getTFNumberOfClusters() const;
+  virtual int getTFNumberOfTracklets() const;
+  virtual int getTFNumberOfCells() const;
+
+  float mBz = 5.f;
+
  private:
   void traverseCellsTree(const int, const int);
   track::TrackParCov buildTrackSeed(const Cluster& cluster1, const Cluster& cluster2, const Cluster& cluster3, const TrackingFrameInfo& tf3, float resolution);
@@ -93,12 +103,13 @@ class TrackerTraits
   int mNThreads = 1;
   bool mApplySmoothing = false;
   o2::base::PropagatorImpl<float>::MatCorrType mCorrType = o2::base::PropagatorImpl<float>::MatCorrType::USEMatCorrNONE;
-  float mBz = 5.f;
+
+  // virtual bool checkTFMemory(const int iteration);
 
  protected:
   TimeFrame* mTimeFrame;
   std::vector<TrackingParameters> mTrkParams;
-
+  bool mIsGPU = false;
   o2::gpu::GPUChainITS* mChain = nullptr;
   FuncRunITSTrackFit_t mChainRunITSTrackFit;
 };
@@ -123,6 +134,12 @@ inline const int4 TrackerTraits::getBinsRect(const Cluster& currentCluster, int
   return getBinsRect(layerIndex, currentCluster.phi, maxdeltaphi, z1, z2, maxdeltaz);
 }
 
+inline void TrackerTraits::initialiseTimeFrame(const int iteration)
+{
+  mTimeFrame->initialise(iteration, mTrkParams[iteration], 7);
+  setIsGPU(false);
+}
+
 inline const int4 TrackerTraits::getBinsRect(const int layerIndex, float phi, float maxdeltaphi,
                                              float z1, float z2, float maxdeltaz)
 {
diff --git a/Detectors/ITSMFT/ITS/tracking/include/ITStracking/TrackingConfigParam.h b/Detectors/ITSMFT/ITS/tracking/include/ITStracking/TrackingConfigParam.h
index 910439094e8ab..8d6ce5d55887c 100644
--- a/Detectors/ITSMFT/ITS/tracking/include/ITStracking/TrackingConfigParam.h
+++ b/Detectors/ITSMFT/ITS/tracking/include/ITStracking/TrackingConfigParam.h
@@ -30,7 +30,7 @@ struct VertexerParamConfig : public o2::conf::ConfigurableParamHelper<VertexerPa
   float histPairCut = 0.04f;
   float tanLambdaCut = 0.002f; // tanLambda = deltaZ/deltaR
   int clusterContributorsCut = 16;
-  int maxTrackletsPerCluster = 2e3;
+  int maxTrackletsPerCluster = 1e2;
   int phiSpan = -1;
   int zSpan = -1;
 
diff --git a/Detectors/ITSMFT/ITS/tracking/include/ITStracking/Tracklet.h b/Detectors/ITSMFT/ITS/tracking/include/ITStracking/Tracklet.h
index 877e348d80d81..3508c541caa1c 100644
--- a/Detectors/ITSMFT/ITS/tracking/include/ITStracking/Tracklet.h
+++ b/Detectors/ITSMFT/ITS/tracking/include/ITStracking/Tracklet.h
@@ -27,17 +27,18 @@ namespace its
 {
 
 struct Tracklet final {
-  Tracklet();
-  GPUdi() Tracklet(const int, const int, const Cluster&, const Cluster&, int rof0, int rof1);
-  GPUdi() Tracklet(const int, const int, float tanL, float phi, int rof0, int rof1);
-  bool operator==(const Tracklet&) const;
-  bool operator!=(const Tracklet&) const;
-  GPUhd() unsigned char isEmpty() const
+  GPUhdi() Tracklet();
+  GPUhdi() Tracklet(const int, const int, const Cluster&, const Cluster&, int rof0, int rof1);
+  GPUhdi() Tracklet(const int, const int, float tanL, float phi, int rof0, int rof1);
+  GPUhdi() bool operator==(const Tracklet&) const;
+  GPUhdi() bool operator!=(const Tracklet&) const;
+  GPUhdi() unsigned char isEmpty() const
   {
-    return !firstClusterIndex && !secondClusterIndex && !tanLambda && !phi;
+    return firstClusterIndex < 0 || secondClusterIndex < 0;
   }
-  void dump();
-  unsigned char operator<(const Tracklet&) const;
+  GPUhdi() void dump();
+  GPUhdi() void dump() const;
+  GPUhdi() unsigned char operator<(const Tracklet&) const;
 
   int firstClusterIndex;
   int secondClusterIndex;
@@ -46,13 +47,12 @@ struct Tracklet final {
   unsigned short rof[2];
 };
 
-inline Tracklet::Tracklet() : firstClusterIndex{0}, secondClusterIndex{0}, tanLambda{0.0f}, phi{0.0f}
+GPUhdi() Tracklet::Tracklet() : firstClusterIndex{-1}, secondClusterIndex{-1}, tanLambda{0.0f}, phi{0.0f}
 {
-  // Nothing to do
 }
 
-GPUdi() Tracklet::Tracklet(const int firstClusterOrderingIndex, const int secondClusterOrderingIndex,
-                           const Cluster& firstCluster, const Cluster& secondCluster, int rof0 = -1, int rof1 = -1)
+GPUhdi() Tracklet::Tracklet(const int firstClusterOrderingIndex, const int secondClusterOrderingIndex,
+                            const Cluster& firstCluster, const Cluster& secondCluster, int rof0 = -1, int rof1 = -1)
   : firstClusterIndex{firstClusterOrderingIndex},
     secondClusterIndex{secondClusterOrderingIndex},
     tanLambda{(firstCluster.zCoordinate - secondCluster.zCoordinate) /
@@ -64,7 +64,7 @@ GPUdi() Tracklet::Tracklet(const int firstClusterOrderingIndex, const int second
   // Nothing to do
 }
 
-GPUdi() Tracklet::Tracklet(const int idx0, const int idx1, float tanL, float phi, int rof0, int rof1)
+GPUhdi() Tracklet::Tracklet(const int idx0, const int idx1, float tanL, float phi, int rof0, int rof1)
   : firstClusterIndex{idx0},
     secondClusterIndex{idx1},
     tanLambda{tanL},
@@ -74,15 +74,17 @@ GPUdi() Tracklet::Tracklet(const int idx0, const int idx1, float tanL, float phi
   // Nothing to do
 }
 
-inline bool Tracklet::operator==(const Tracklet& rhs) const
+GPUhdi() bool Tracklet::operator==(const Tracklet& rhs) const
 {
   return this->firstClusterIndex == rhs.firstClusterIndex &&
          this->secondClusterIndex == rhs.secondClusterIndex &&
          this->tanLambda == rhs.tanLambda &&
-         this->phi == rhs.phi;
+         this->phi == rhs.phi &&
+         this->rof[0] == rhs.rof[0] &&
+         this->rof[1] == rhs.rof[1];
 }
 
-inline bool Tracklet::operator!=(const Tracklet& rhs) const
+GPUhdi() bool Tracklet::operator!=(const Tracklet& rhs) const
 {
   return this->firstClusterIndex != rhs.firstClusterIndex ||
          this->secondClusterIndex != rhs.secondClusterIndex ||
@@ -90,25 +92,25 @@ inline bool Tracklet::operator!=(const Tracklet& rhs) const
          this->phi != rhs.phi;
 }
 
-inline unsigned char Tracklet::operator<(const Tracklet& t) const
+GPUhdi() unsigned char Tracklet::operator<(const Tracklet& t) const
 {
-  if (isEmpty() && t.isEmpty()) {
+  if (isEmpty()) {
     return false;
-  } else {
-    if (isEmpty()) {
-      return false;
-    }
   }
   return true;
 }
 
-inline void Tracklet::dump()
+GPUhdi() void Tracklet::dump()
+{
+  printf("fClIdx: %d sClIdx: %d  rof1: %hu rof2: %hu\n", firstClusterIndex, secondClusterIndex, rof[0], rof[1]);
+}
+
+GPUhdi() void Tracklet::dump() const
 {
-  std::cout << "firstClusterIndex: " << firstClusterIndex << std::endl;
-  std::cout << "secondClusterIndex: " << secondClusterIndex << std::endl;
-  std::cout << "tanLambda: " << tanLambda << std::endl;
-  std::cout << "phi: " << phi << std::endl;
+  printf("fClIdx: %d sClIdx: %d  rof1: %hu rof2: %hu\n", firstClusterIndex, secondClusterIndex, rof[0], rof[1]);
 }
+// tanL: %f phi: %f
+// tanLambda, phi,
 
 } // namespace its
 } // namespace o2
diff --git a/Detectors/ITSMFT/ITS/tracking/include/ITStracking/Vertexer.h b/Detectors/ITSMFT/ITS/tracking/include/ITStracking/Vertexer.h
index a4573c2f5d649..73e19a9dd55b5 100644
--- a/Detectors/ITSMFT/ITS/tracking/include/ITStracking/Vertexer.h
+++ b/Detectors/ITSMFT/ITS/tracking/include/ITStracking/Vertexer.h
@@ -62,8 +62,7 @@ class Vertexer
   std::vector<Vertex> exportVertices();
   VertexerTraits* getTraits() const { return mTraits; };
 
-  float clustersToVertices(
-    const bool useMc = false, std::function<void(std::string s)> = [](std::string s) { std::cout << s << std::endl; });
+  float clustersToVertices(std::function<void(std::string s)> = [](std::string s) { std::cout << s << std::endl; });
   void filterMCTracklets();
   void validateTracklets();
 
@@ -82,15 +81,7 @@ class Vertexer
   // Utils
   void dumpTraits();
   template <typename... T>
-  float evaluateTask(void (Vertexer::*)(T...), bool, const char*, std::function<void(std::string s)> logger, T&&... args);
-
-  // debug
-  void setDebugCombinatorics();
-  void setDebugTrackletSelection();
-  void setDebugLines();
-  void setDebugSummaryLines();
-  void setDebugCentroidsHistograms();
-  // \debug
+  float evaluateTask(void (Vertexer::*)(T...), const char*, std::function<void(std::string s)> logger, T&&... args);
 
  private:
   std::uint32_t mROframe = 0;
@@ -116,11 +107,6 @@ void Vertexer::findTracklets(T&&... args)
   mTraits->computeTracklets(std::forward<T>(args)...);
 }
 
-// inline void Vertexer::findTrivialMCTracklets()
-// {
-//   mTraits->computeTrackletsPureMontecarlo();
-// }
-
 inline VertexingParameters Vertexer::getVertParameters() const
 {
   return mTraits->getVertexingParameters();
@@ -152,7 +138,7 @@ inline std::vector<Vertex> Vertexer::exportVertices()
 }
 
 template <typename... T>
-float Vertexer::evaluateTask(void (Vertexer::*task)(T...), bool verbose, const char* taskName, std::function<void(std::string s)> logger,
+float Vertexer::evaluateTask(void (Vertexer::*task)(T...), const char* taskName, std::function<void(std::string s)> logger,
                              T&&... args)
 {
   float diff{0.f};
@@ -171,9 +157,7 @@ float Vertexer::evaluateTask(void (Vertexer::*task)(T...), bool verbose, const c
     } else {
       sstream << std::setw(2) << " - " << taskName << " completed in: " << diff << " ms";
     }
-    if (verbose) {
-      logger(sstream.str());
-    }
+    logger(sstream.str());
   } else {
     (this->*task)(std::forward<T>(args)...);
   }
@@ -181,31 +165,6 @@ float Vertexer::evaluateTask(void (Vertexer::*task)(T...), bool verbose, const c
   return diff;
 }
 
-inline void Vertexer::setDebugCombinatorics()
-{
-  mTraits->setDebugFlag(VertexerDebug::CombinatoricsTreeAll);
-}
-
-inline void Vertexer::setDebugTrackletSelection()
-{
-  mTraits->setDebugFlag(VertexerDebug::TrackletTreeAll);
-}
-
-inline void Vertexer::setDebugLines()
-{
-  mTraits->setDebugFlag(VertexerDebug::LineTreeAll);
-}
-
-inline void Vertexer::setDebugSummaryLines()
-{
-  mTraits->setDebugFlag(VertexerDebug::LineSummaryAll);
-}
-
-inline void Vertexer::setDebugCentroidsHistograms()
-{
-  mTraits->setDebugFlag(VertexerDebug::HistCentroids);
-}
-
 } // namespace its
 } // namespace o2
 #endif
diff --git a/Detectors/ITSMFT/ITS/tracking/include/ITStracking/VertexerTraits.h b/Detectors/ITSMFT/ITS/tracking/include/ITStracking/VertexerTraits.h
index 56a56115e78be..b7af712e3d7f4 100644
--- a/Detectors/ITSMFT/ITS/tracking/include/ITStracking/VertexerTraits.h
+++ b/Detectors/ITSMFT/ITS/tracking/include/ITStracking/VertexerTraits.h
@@ -38,44 +38,13 @@ class MCCompLabel;
 namespace its
 {
 class ROframe;
-
 using constants::its::LayersNumberVertexer;
 
-struct lightVertex {
-  lightVertex(float x, float y, float z, std::array<float, 6> rms2, int cont, float avgdis2, int stamp);
-  float mX;
-  float mY;
-  float mZ;
-  std::array<float, 6> mRMS2;
-  float mAvgDistance2;
-  int mContributors;
-  int mTimeStamp;
-};
-
-struct ClusterMCLabelInfo {
-  int TrackId;
-  int MotherId;
-  int EventId;
-  float Pt;
-};
-
-enum class VertexerDebug : unsigned int {
-  TrackletTreeAll = 0x1 << 1,
-  LineTreeAll = 0x1 << 2,
-  CombinatoricsTreeAll = 0x1 << 3,
-  LineSummaryAll = 0x1 << 4,
-  HistCentroids = 0x1 << 5
-};
-
 enum class TrackletMode {
   Layer0Layer1 = 0,
   Layer1Layer2 = 2
 };
 
-inline lightVertex::lightVertex(float x, float y, float z, std::array<float, 6> rms2, int cont, float avgdis2, int stamp) : mX{x}, mY{y}, mZ{z}, mRMS2{rms2}, mAvgDistance2{avgdis2}, mContributors{cont}, mTimeStamp{stamp}
-{
-}
-
 class VertexerTraits
 {
  public:
@@ -96,18 +65,11 @@ class VertexerTraits
   virtual void initialise(const TrackingParameters& trackingParams);
   virtual void computeTracklets();
   virtual void computeTrackletMatching();
-  // virtual void computeMCFiltering();
-  // virtual void filterTrackletsWithMC(std::vector<Tracklet>&,
-  //                                    std::vector<Tracklet>&,
-  //                                    std::vector<int>&,
-  //                                    std::vector<int>&,
-  //                                    const int);
-
-  // virtual void computeTrackletsPureMontecarlo();
   virtual void computeVertices();
+  virtual void adoptTimeFrame(TimeFrame* tf);
+  virtual void updateVertexingParameters(const VertexingParameters& vrtPar);
   // virtual void computeHistVertices();
 
-  void updateVertexingParameters(const VertexingParameters& vrtPar);
   VertexingParameters getVertexingParameters() const { return mVrtParams; }
   static const std::vector<std::pair<int, int>> selectClusters(const int* indexTable,
                                                                const std::array<int, 4>& selectedBinsRect,
@@ -115,18 +77,13 @@ class VertexerTraits
   std::vector<lightVertex> getVertices() const { return mVertices; }
 
   // utils
-  virtual void adoptTimeFrame(TimeFrame* tf);
+
   void setIsGPU(const unsigned char isgpu) { mIsGPU = isgpu; };
   unsigned char getIsGPU() const { return mIsGPU; };
   void dumpVertexerTraits();
 
-  void setDebugFlag(VertexerDebug flag, const unsigned char on);
-  unsigned char isDebugFlag(const VertexerDebug& flags) const;
-  unsigned int getDebugFlags() const { return static_cast<unsigned int>(mDBGFlags); }
-
  protected:
   unsigned char mIsGPU;
-  unsigned int mDBGFlags = 0;
 
   VertexingParameters mVrtParams;
   IndexTableUtils mIndexTableUtils;
@@ -145,15 +102,6 @@ inline void VertexerTraits::initialise(const TrackingParameters& trackingParams)
   setIsGPU(false);
 }
 
-inline void VertexerTraits::updateVertexingParameters(const VertexingParameters& vrtPar)
-{
-  mVrtParams = vrtPar;
-  mIndexTableUtils.setTrackingParameters(vrtPar);
-  mVrtParams.phiSpan = static_cast<int>(std::ceil(mIndexTableUtils.getNphiBins() * mVrtParams.phiCut /
-                                                  constants::math::TwoPi));
-  mVrtParams.zSpan = static_cast<int>(std::ceil(mVrtParams.zCut * mIndexTableUtils.getInverseZCoordinate(0)));
-}
-
 GPUhdi() const int2 VertexerTraits::getPhiBins(float phi, float dPhi)
 {
   return VertexerTraits::getPhiBins(phi, dPhi, mIndexTableUtils);
@@ -191,21 +139,6 @@ GPUhdi() const int4 VertexerTraits::getBinsRect(const Cluster& currentCluster, c
   return VertexerTraits::getBinsRect(currentCluster, layerIndex, directionZIntersection, maxdeltaz, maxdeltaphi, mIndexTableUtils);
 }
 
-// debug
-inline void VertexerTraits::setDebugFlag(VertexerDebug flag, const unsigned char on = true)
-{
-  if (on) {
-    mDBGFlags |= static_cast<unsigned int>(flag);
-  } else {
-    mDBGFlags &= ~static_cast<unsigned int>(flag);
-  }
-}
-
-inline unsigned char VertexerTraits::isDebugFlag(const VertexerDebug& flags) const
-{
-  return mDBGFlags & static_cast<unsigned int>(flags);
-}
-
 inline void VertexerTraits::adoptTimeFrame(TimeFrame* tf) { mTimeFrame = tf; }
 
 } // namespace its
diff --git a/Detectors/ITSMFT/ITS/tracking/src/TimeFrame.cxx b/Detectors/ITSMFT/ITS/tracking/src/TimeFrame.cxx
index b61f73afa062a..69fe1aca5ce5a 100644
--- a/Detectors/ITSMFT/ITS/tracking/src/TimeFrame.cxx
+++ b/Detectors/ITSMFT/ITS/tracking/src/TimeFrame.cxx
@@ -51,6 +51,10 @@ namespace o2
 namespace its
 {
 
+lightVertex::lightVertex(float x, float y, float z, std::array<float, 6> rms2, int cont, float avgdis2, int stamp) : mX{x}, mY{y}, mZ{z}, mRMS2{rms2}, mAvgDistance2{avgdis2}, mContributors{cont}, mTimeStamp{stamp}
+{
+}
+
 constexpr float DefClusErrorRow = o2::itsmft::SegmentationAlpide::PitchRow * 0.5;
 constexpr float DefClusErrorCol = o2::itsmft::SegmentationAlpide::PitchCol * 0.5;
 constexpr float DefClusError2Row = DefClusErrorRow * DefClusErrorRow;
@@ -93,6 +97,16 @@ void TimeFrame::addPrimaryVertices(const gsl::span<const Vertex>& vertices)
   mROframesPV.push_back(mPrimaryVertices.size());
 }
 
+void TimeFrame::addPrimaryVertices(const std::vector<lightVertex>& lVertices)
+{
+  std::vector<Vertex> vertices;
+  for (auto& vertex : lVertices) {
+    vertices.emplace_back(o2::math_utils::Point3D<float>(vertex.mX, vertex.mY, vertex.mZ), vertex.mRMS2, vertex.mContributors, vertex.mAvgDistance2);
+    vertices.back().setTimeStamp(vertex.mTimeStamp);
+  }
+  addPrimaryVertices(vertices);
+}
+
 int TimeFrame::loadROFrameData(const o2::itsmft::ROFRecord& rof, gsl::span<const itsmft::Cluster> clusters,
                                const dataformats::MCTruthContainer<MCCompLabel>* mcLabels)
 {
@@ -216,7 +230,6 @@ void TimeFrame::initialise(const int iteration, const TrackingParameters& trkPar
     mTracklets.resize(std::min(trkParam.TrackletsPerRoad(), maxLayers - 1));
     mTrackletLabels.resize(trkParam.TrackletsPerRoad());
     mTrackletsLookupTable.resize(trkParam.CellsPerRoad());
-    mIndexTables.clear();
     mIndexTableUtils.setTrackingParameters(trkParam);
     mPositionResolution.resize(trkParam.NLayers);
     mBogusClusters.resize(trkParam.NLayers, 0);
diff --git a/Detectors/ITSMFT/ITS/tracking/src/Tracker.cxx b/Detectors/ITSMFT/ITS/tracking/src/Tracker.cxx
index b603469d8cb2a..e21fd29b778c5 100644
--- a/Detectors/ITSMFT/ITS/tracking/src/Tracker.cxx
+++ b/Detectors/ITSMFT/ITS/tracking/src/Tracker.cxx
@@ -52,24 +52,24 @@ void Tracker::clustersToTracks(std::function<void(std::string s)> logger, std::f
   for (int iteration = 0; iteration < (int)mTrkParams.size(); ++iteration) {
     total += evaluateTask(&Tracker::initialiseTimeFrame, "Timeframe initialisation", logger, iteration);
     total += evaluateTask(&Tracker::computeTracklets, "Tracklet finding", logger, iteration);
-    logger(fmt::format("\t- Number of tracklets: {}", mTimeFrame->getNumberOfTracklets()));
+    logger(fmt::format("\t- Number of tracklets: {}", mTraits->getTFNumberOfTracklets()));
     if (!mTimeFrame->checkMemory(mTrkParams[iteration].MaxMemory)) {
       error("Too much memory used during trackleting, check the detector status and/or the selections.");
       break;
     }
-    float trackletsPerCluster = mTimeFrame->getNumberOfClusters() > 0 ? float(mTimeFrame->getNumberOfTracklets()) / mTimeFrame->getNumberOfClusters() : 0.f;
+    float trackletsPerCluster = mTraits->getTFNumberOfClusters() > 0 ? float(mTraits->getTFNumberOfTracklets()) / mTraits->getTFNumberOfClusters() : 0.f;
     if (trackletsPerCluster > mTrkParams[iteration].TrackletsPerClusterLimit) {
       error(fmt::format("Too many tracklets per cluster ({}), check the detector status and/or the selections.", trackletsPerCluster));
       break;
     }
 
     total += evaluateTask(&Tracker::computeCells, "Cell finding", logger, iteration);
-    logger(fmt::format("\t- Number of Cells: {}", mTimeFrame->getNumberOfCells()));
+    logger(fmt::format("\t- Number of Cells: {}", mTraits->getTFNumberOfCells()));
     if (!mTimeFrame->checkMemory(mTrkParams[iteration].MaxMemory)) {
       error("Too much memory used during cell finding, check the detector status and/or the selections.");
       break;
     }
-    float cellsPerCluster = mTimeFrame->getNumberOfClusters() > 0 ? float(mTimeFrame->getNumberOfCells()) / mTimeFrame->getNumberOfClusters() : 0.f;
+    float cellsPerCluster = mTraits->getTFNumberOfClusters() > 0 ? float(mTraits->getTFNumberOfCells()) / mTraits->getTFNumberOfClusters() : 0.f;
     if (cellsPerCluster > mTrkParams[iteration].CellsPerClusterLimit) {
       error(fmt::format("Too many cells per cluster ({}), check the detector status and/or the selections.", cellsPerCluster));
       break;
@@ -99,6 +99,11 @@ void Tracker::clustersToTracks(std::function<void(std::string s)> logger, std::f
   mNumberOfRuns++;
 }
 
+void Tracker::initialiseTimeFrame(int& iteration)
+{
+  mTraits->initialiseTimeFrame(iteration);
+}
+
 void Tracker::computeTracklets(int& iteration)
 {
   mTraits->computeLayerTracklets(iteration);
diff --git a/Detectors/ITSMFT/ITS/tracking/src/TrackerTraits.cxx b/Detectors/ITSMFT/ITS/tracking/src/TrackerTraits.cxx
index 8b89f37e28896..c805c9e4ba11f 100644
--- a/Detectors/ITSMFT/ITS/tracking/src/TrackerTraits.cxx
+++ b/Detectors/ITSMFT/ITS/tracking/src/TrackerTraits.cxx
@@ -130,11 +130,12 @@ void TrackerTraits::computeLayerTracklets(const int iteration)
               const int maxRowClusterIndex = tf->getIndexTable(rof1, iLayer + 1)[maxBinIndex];
 
               for (int iNextCluster{firstRowClusterIndex}; iNextCluster < maxRowClusterIndex; ++iNextCluster) {
+
                 if (iNextCluster >= (int)layer1.size()) {
                   break;
                 }
-                const Cluster& nextCluster{layer1[iNextCluster]};
 
+                const Cluster& nextCluster{layer1[iNextCluster]};
                 if (tf->isClusterUsed(iLayer + 1, nextCluster.clusterId)) {
                   continue;
                 }
@@ -183,7 +184,6 @@ void TrackerTraits::computeLayerTracklets(const int iteration)
       }
     }
   }
-  /// Cold code, fixups
 
   for (int iLayer{0}; iLayer < mTrkParams[iteration].CellsPerRoad(); ++iLayer) {
     /// Sort tracklets
@@ -267,9 +267,10 @@ void TrackerTraits::computeLayerCells(const int iteration)
       continue;
     }
 
+#ifdef OPTIMISATION_OUTPUT
     float resolution{std::sqrt(Sq(mTrkParams[iteration].LayerMisalignment[iLayer]) + Sq(mTrkParams[iteration].LayerMisalignment[iLayer + 1]) + Sq(mTrkParams[iteration].LayerMisalignment[iLayer + 2])) / mTrkParams[iteration].LayerResolution[iLayer]};
     resolution = resolution > 1.e-12 ? resolution : 1.f;
-
+#endif
     const int currentLayerTrackletsNum{static_cast<int>(tf->getTracklets()[iLayer].size())};
 
     for (int iTracklet{0}; iTracklet < currentLayerTrackletsNum; ++iTracklet) {
@@ -996,5 +997,30 @@ void TrackerTraits::setNThreads(int n)
   mNThreads = 1;
 #endif
 }
+int TrackerTraits::getTFNumberOfClusters() const
+{
+  return mTimeFrame->getNumberOfClusters();
+}
+
+int TrackerTraits::getTFNumberOfTracklets() const
+{
+  return mTimeFrame->getNumberOfTracklets();
+}
+
+int TrackerTraits::getTFNumberOfCells() const
+{
+  return mTimeFrame->getNumberOfCells();
+}
+
+void TrackerTraits::adoptTimeFrame(TimeFrame* tf)
+{
+  mTimeFrame = tf;
+}
+
+// bool TrackerTraits::checkTFMemory(const int iteration)
+// {
+//   return mTimeFrame->checkMemory(mTrkParams[iteration].MaxMemory);
+// }
+
 } // namespace its
 } // namespace o2
diff --git a/Detectors/ITSMFT/ITS/tracking/src/Vertexer.cxx b/Detectors/ITSMFT/ITS/tracking/src/Vertexer.cxx
index c530f676a0dd8..f00a5b78bbd11 100644
--- a/Detectors/ITSMFT/ITS/tracking/src/Vertexer.cxx
+++ b/Detectors/ITSMFT/ITS/tracking/src/Vertexer.cxx
@@ -36,14 +36,14 @@ Vertexer::Vertexer(VertexerTraits* traits)
   mTraits = traits;
 }
 
-float Vertexer::clustersToVertices(const bool useMc, std::function<void(std::string s)> logger)
+float Vertexer::clustersToVertices(std::function<void(std::string s)> logger)
 {
   float total{0.f};
   TrackingParameters trkPars;
-  total += evaluateTask(&Vertexer::initialiseVertexer, false, "Vertexer initialisation", logger, trkPars);
-  total += evaluateTask(&Vertexer::findTracklets, false, "Tracklet finding", logger);
-  total += evaluateTask(&Vertexer::validateTracklets, false, "Adjacent tracklets validation", logger);
-  total += evaluateTask(&Vertexer::findVertices, false, "Vertex finding", logger);
+  total += evaluateTask(&Vertexer::initialiseVertexer, "Vertexer initialisation", logger, trkPars);
+  total += evaluateTask(&Vertexer::findTracklets, "Tracklet finding", logger);
+  total += evaluateTask(&Vertexer::validateTracklets, "Adjacent tracklets validation", logger);
+  total += evaluateTask(&Vertexer::findVertices, "Vertex finding", logger);
 
   return total;
 }
diff --git a/Detectors/ITSMFT/ITS/tracking/src/VertexerTraits.cxx b/Detectors/ITSMFT/ITS/tracking/src/VertexerTraits.cxx
index 967daef6f5e63..859414545a15a 100644
--- a/Detectors/ITSMFT/ITS/tracking/src/VertexerTraits.cxx
+++ b/Detectors/ITSMFT/ITS/tracking/src/VertexerTraits.cxx
@@ -24,8 +24,6 @@
 #include <ostream>
 #endif
 
-// #define VTX_DEBUG
-
 namespace o2
 {
 namespace its
@@ -133,6 +131,15 @@ const std::vector<std::pair<int, int>> VertexerTraits::selectClusters(const int*
   return filteredBins;
 }
 
+void VertexerTraits::updateVertexingParameters(const VertexingParameters& vrtPar)
+{
+  mVrtParams = vrtPar;
+  mIndexTableUtils.setTrackingParameters(vrtPar);
+  mVrtParams.phiSpan = static_cast<int>(std::ceil(mIndexTableUtils.getNphiBins() * mVrtParams.phiCut /
+                                                  constants::math::TwoPi));
+  mVrtParams.zSpan = static_cast<int>(std::ceil(mVrtParams.zCut * mIndexTableUtils.getInverseZCoordinate(0)));
+}
+
 void VertexerTraits::computeTracklets()
 {
   for (int rofId{0}; rofId < mTimeFrame->getNrof(); ++rofId) {
diff --git a/Detectors/ITSMFT/ITS/workflow/src/CookedTrackerSpec.cxx b/Detectors/ITSMFT/ITS/workflow/src/CookedTrackerSpec.cxx
index 1195f02afffc6..b221f01fcde1d 100644
--- a/Detectors/ITSMFT/ITS/workflow/src/CookedTrackerSpec.cxx
+++ b/Detectors/ITSMFT/ITS/workflow/src/CookedTrackerSpec.cxx
@@ -139,7 +139,7 @@ void CookedTrackerDPL::run(ProcessingContext& pc)
   mTimeFrame.setMultiplicityCutMask(processingMask);
   float vertexerElapsedTime;
   if (mRunVertexer) {
-    vertexerElapsedTime = mVertexerPtr->clustersToVertices(false, [&](std::string s) { LOG(info) << s; });
+    vertexerElapsedTime = mVertexerPtr->clustersToVertices([&](std::string s) { LOG(info) << s; });
   }
   LOG(info) << fmt::format(" - Vertex seeding total elapsed time: {} ms in {} ROFs", vertexerElapsedTime, rofspan.size());
   for (size_t iRof{0}; iRof < rofspan.size(); ++iRof) {
diff --git a/Detectors/ITSMFT/ITS/workflow/src/TrackerSpec.cxx b/Detectors/ITSMFT/ITS/workflow/src/TrackerSpec.cxx
index 7e7fb817142a2..1f6e23cacd5dd 100644
--- a/Detectors/ITSMFT/ITS/workflow/src/TrackerSpec.cxx
+++ b/Detectors/ITSMFT/ITS/workflow/src/TrackerSpec.cxx
@@ -209,11 +209,10 @@ void TrackerDPL::run(ProcessingContext& pc)
   std::vector<bool> processingMask;
   int cutVertexMult{0}, cutRandomMult = int(rofs.size()) - multEst.selectROFs(rofs, compClusters, physTriggers, processingMask);
   timeFrame->setMultiplicityCutMask(processingMask);
-
   float vertexerElapsedTime{0.f};
   if (mRunVertexer) {
     // Run seeding vertexer
-    vertexerElapsedTime = mVertexer->clustersToVertices(false, logger);
+    vertexerElapsedTime = mVertexer->clustersToVertices(logger);
   }
   const auto& multEstConf = FastMultEstConfig::Instance(); // parameters for mult estimation and cuts
   for (auto iRof{0}; iRof < rofspan.size(); ++iRof) {
@@ -258,13 +257,13 @@ void TrackerDPL::run(ProcessingContext& pc)
   } else {
 
     timeFrame->setMultiplicityCutMask(processingMask);
+    // Run CA tracker
     mTracker->clustersToTracks(logger, errorLogger);
     if (timeFrame->hasBogusClusters()) {
       LOG(warning) << fmt::format(" - The processed timeframe had {} clusters with wild z coordinates, check the dictionaries", timeFrame->hasBogusClusters());
     }
 
     for (unsigned int iROF{0}; iROF < rofs.size(); ++iROF) {
-
       auto& rof{rofs[iROF]};
       tracks = timeFrame->getTracks(iROF);
       trackLabels = timeFrame->getTracksLabel(iROF);
diff --git a/GPU/Common/GPUCommonMath.h b/GPU/Common/GPUCommonMath.h
index 18c9d695ce242..995d4eb503e8b 100644
--- a/GPU/Common/GPUCommonMath.h
+++ b/GPU/Common/GPUCommonMath.h
@@ -60,7 +60,7 @@ class GPUCommonMath
   GPUd() static float ASin(float x);
   GPUd() static float ACos(float x);
   GPUd() static float ATan(float x);
-  GPUd() static float ATan2(float y, float x);
+  GPUhd() static float ATan2(float y, float x);
   GPUd() static float Sin(float x);
   GPUd() static float Cos(float x);
   GPUhdni() static void SinCos(float x, float& s, float& c);
@@ -243,7 +243,7 @@ GPUdi() bool GPUCommonMath::Finite(float x) { return CHOICE(std::isfinite(x), tr
 
 GPUdi() float GPUCommonMath::ATan(float x) { return CHOICE(atanf(x), atanf(x), atan(x)); }
 
-GPUdi() float GPUCommonMath::ATan2(float y, float x) { return CHOICE(atan2f(y, x), atan2f(y, x), atan2(y, x)); }
+GPUhdi() float GPUCommonMath::ATan2(float y, float x) { return CHOICE(atan2f(y, x), atan2f(y, x), atan2(y, x)); }
 
 GPUdi() float GPUCommonMath::Sin(float x) { return CHOICE(sinf(x), sinf(x), sin(x)); }
 
diff --git a/macro/run_trac_ca_its.C b/macro/run_trac_ca_its.C
index c135f81899732..2bac28c38fbc1 100644
--- a/macro/run_trac_ca_its.C
+++ b/macro/run_trac_ca_its.C
@@ -222,7 +222,7 @@ void run_trac_ca_its(bool cosmics = false,
 
   int rofId{0};
   vertexer.adoptTimeFrame(tf);
-  vertexer.clustersToVertices(false);
+  vertexer.clustersToVertices();
 
   tf.printVertices();
 
diff --git a/macro/run_trac_its.C b/macro/run_trac_its.C
index c2ce2cd438bc5..67bc075c5754d 100644
--- a/macro/run_trac_its.C
+++ b/macro/run_trac_its.C
@@ -169,7 +169,7 @@ void run_trac_its(std::string path = "./", std::string outputfile = "o2trac_its.
     tf.loadROFrameData(rofspan, clSpan, pattIt_vertexer, dict, labels);
     tf.setMultiplicityCutMask(processingMask);
     vertexer.adoptTimeFrame(tf);
-    vertexer.clustersToVertices(mcTruth);
+    vertexer.clustersToVertices();
     int iRof = 0;
     for (auto& rof : *rofs) {
       auto it = pattIt;