TRD calibration subdirectory added (#5894)

* TRD calibration subdirectory added

- template for vDrift calibration
- after global tracking angular differences between TRD-only tracks and tracklets are collected
- small correction to always update to the correct pad plane during track following, since the pad row size depends also on the stack, not only on the layer

* Remove copy-paste artifact (wrong comment)

DataFormats/Detectors/TRD/include/DataFormatsTRD/Constants.h

@@ -56,6 +56,7 @@ constexpr float GRANULARITYTRKLSLOPE = 1.f / PADGRANULARITYTRKLSLOPE; // granula
 
 // OS: Should this not be flexible for example in case of Kr calib?
 constexpr int TIMEBINS = 30; // the number of time bins
+constexpr int NBINSANGLEDIFF = 26; // the number of bins for the track angle used for the vDrift and ExB calibration based on the tracking (last bin is for under-/overflow entries)
 
 // Trigger parameters
 constexpr double READOUT_TIME = 3000;                  // the time the readout takes, as 30 TB = 3 micro-s.

Detectors/TRD/CMakeLists.txt

@@ -9,6 +9,7 @@
 # submit itself to any jurisdiction.
 
 add_subdirectory(base)
+add_subdirectory(calibration)
 add_subdirectory(simulation)
 add_subdirectory(reconstruction)
 add_subdirectory(macros)

Detectors/TRD/calibration/CMakeLists.txt

@@ -0,0 +1,17 @@
+#Copyright CERN and copyright holders of ALICE O2.This software is distributed
+#under the terms of the GNU General Public License v3(GPL Version 3), copied
+#verbatim in the file "COPYING".
+#
+#See http: //alice-o2.web.cern.ch/license for full licensing information.
+#
+#In applying this license CERN does not waive the privileges and immunities
+#granted to it by virtue of its status as an Intergovernmental Organization or
+#submit itself to any jurisdiction.
+
+o2_add_library(TRDCalibration
+               SOURCES src/CalibVDrift.cxx
+               PUBLIC_LINK_LIBRARIES O2::TRDBase
+                                     O2::DataFormatsTRD)
+
+ o2_target_root_dictionary(TRDCalibration
+                           HEADERS include/TRDCalibration/CalibVDrift.h)

Detectors/TRD/calibration/include/TRDCalibration/CalibVDrift.h

@@ -0,0 +1,68 @@
+// Copyright CERN and copyright holders of ALICE O2. This software is
+// distributed under the terms of the GNU General Public License v3 (GPL
+// Version 3), copied verbatim in the file "COPYING".
+//
+// See http://alice-o2.web.cern.ch/license for full licensing information.
+//
+// In applying this license CERN does not waive the privileges and immunities
+// granted to it by virtue of its status as an Intergovernmental Organization
+// or submit itself to any jurisdiction.
+
+#ifndef ALICEO2_TRD_CALIBVDRIFT_H_
+#define ALICEO2_TRD_CALIBVDRIFT_H_
+
+/// \file   CalibVDrift.h
+/// \author Ole Schmidt, ole.schmidt@cern.ch
+
+#include "DataFormatsTRD/Constants.h"
+#include <array>
+
+namespace o2
+{
+namespace trd
+{
+
+/// \brief VDrift calibration class
+///
+/// This class is used to determine chamber-wise vDrift values
+///
+/// origin: TRD
+/// \author Ole Schmidt, ole.schmidt@cern.ch
+
+class CalibVDrift
+{
+ public:
+  /// default constructor
+  CalibVDrift() = default;
+
+  /// default destructor
+  ~CalibVDrift() = default;
+
+  /// set input angular difference sums
+  void setAngleDiffSums(float* input)
+  {
+    for (int i = 0; i < constants::MAXCHAMBER * constants::NBINSANGLEDIFF; ++i) {
+      mAngleDiffSums[i] = input[i];
+    }
+  }
+
+  /// set input angular difference bin counters
+  void setAngleDiffCounters(short* input)
+  {
+    for (int i = 0; i < constants::MAXCHAMBER * constants::NBINSANGLEDIFF; ++i) {
+      mAngleDiffCounters[i] = input[i];
+    }
+  }
+
+  /// main processing function
+  void process();
+
+ private:
+  std::array<float, constants::MAXCHAMBER * constants::NBINSANGLEDIFF> mAngleDiffSums{};     ///< input TRD track to tracklet angular difference sums per bin
+  std::array<short, constants::MAXCHAMBER * constants::NBINSANGLEDIFF> mAngleDiffCounters{}; ///< input bin counters
+};
+
+} // namespace trd
+
+} // namespace o2
+#endif

Detectors/TRD/calibration/src/CalibVDrift.cxx

@@ -0,0 +1,50 @@
+// Copyright CERN and copyright holders of ALICE O2. This software is
+// distributed under the terms of the GNU General Public License v3 (GPL
+// Version 3), copied verbatim in the file "COPYING".
+//
+// See http://alice-o2.web.cern.ch/license for full licensing information.
+//
+// In applying this license CERN does not waive the privileges and immunities
+// granted to it by virtue of its status as an Intergovernmental Organization
+// or submit itself to any jurisdiction.
+
+/// \file   CalibVDrift.cxx
+/// \author Ole Schmidt, ole.schmidt@cern.ch
+
+#include "TFile.h"
+#include "TH2F.h"
+
+#include <fairlogger/Logger.h>
+
+#include "TRDCalibration/CalibVDrift.h"
+
+using namespace o2::trd;
+using namespace o2::trd::constants;
+
+void CalibVDrift::process()
+{
+  LOG(info) << "Started processing for vDrift calibration";
+
+  for (int iDet = 0; iDet < constants::MAXCHAMBER; ++iDet) {
+    for (int iBin = 0; iBin < constants::NBINSANGLEDIFF; ++iBin) { // note: iBin = constants::NBINSANGLEDIFF - 1 is under-/overflow bin
+      short nEntries = mAngleDiffCounters[iDet * constants::NBINSANGLEDIFF + iBin];
+      float angleDiffSum = mAngleDiffSums[iDet * constants::NBINSANGLEDIFF + iBin];
+      if (nEntries > 0) {
+        LOGF(INFO, "Found %i entrie(s) in chamber %i, bin %i. Average angular deviation: %f", nEntries, iDet, iBin, angleDiffSum / nEntries);
+      }
+    }
+  }
+
+  /*
+  // as an example I loop over the input, create a histogram and write it to a file
+  auto fOut = TFile::Open("trdcalibdummy.root", "recreate");
+  auto hXY = std::make_unique<TH2F>("histDummy", "foo", 100, -60, 60, 100, 250, 400); // xy distribution of TRD space points
+  for (int i = 0; i < mAngulerDeviationProf.size(); ++i) {
+    hXY->Fill(mAngulerDeviationProf[i].mX[0], mAngulerDeviationProf[i].mR);
+  }
+  fOut->cd();
+  hXY->Write();
+  hXY.reset(); // delete the histogram before closing the output file
+  fOut->Close();
+  */
+}

Detectors/TRD/calibration/src/TRDCalibrationLinkDef.h

@@ -0,0 +1,19 @@
+// Copyright CERN and copyright holders of ALICE O2. This software is
+// distributed under the terms of the GNU General Public License v3 (GPL
+// Version 3), copied verbatim in the file "COPYING".
+//
+// See http://alice-o2.web.cern.ch/license for full licensing information.
+//
+// In applying this license CERN does not waive the privileges and immunities
+// granted to it by virtue of its status as an Intergovernmental Organization
+// or submit itself to any jurisdiction.
+
+#ifdef __CLING__
+
+#pragma link off all globals;
+#pragma link off all classes;
+#pragma link off all functions;
+
+#pragma link C++ class o2::trd::CalibVDrift + ;
+
+#endif

Detectors/TRD/workflow/CMakeLists.txt

@@ -25,8 +25,8 @@ o2_add_library(TRDWorkflow
                        src/TRDTrackWriterSpec.cxx
                        src/TRDTrackingWorkflow.cxx
                        src/EntropyDecoderSpec.cxx
-                       src/EntropyEncoderSpec.cxx                       
-                       PUBLIC_LINK_LIBRARIES O2::Framework O2::DPLUtils O2::Steer O2::Algorithm O2::DataFormatsTRD O2::TRDSimulation O2::TRDReconstruction O2::DetectorsBase O2::SimulationDataFormat O2::TRDBase O2::GPUTracking O2::GlobalTrackingWorkflow)
+                       src/EntropyEncoderSpec.cxx
+                       PUBLIC_LINK_LIBRARIES O2::Framework O2::DPLUtils O2::Steer O2::Algorithm O2::DataFormatsTRD O2::TRDSimulation O2::TRDReconstruction O2::DetectorsBase O2::SimulationDataFormat O2::TRDBase O2::TRDCalibration O2::GPUTracking O2::GlobalTrackingWorkflow)
 
                    #o2_target_root_dictionary(TRDWorkflow
                    # HEADERS include/TRDWorkflow/TRDTrapSimulatorSpec.h)

Detectors/TRD/workflow/include/TRDWorkflow/TRDGlobalTrackingSpec.h

@@ -17,6 +17,7 @@
 #include "Framework/Task.h"
 #include "TStopwatch.h"
 #include "TRDBase/GeometryFlat.h"
+#include "TRDCalibration/CalibVDrift.h"
 #include "GPUO2Interface.h"
 #include "GPUTRDTracker.h"
 
@@ -41,6 +42,7 @@ class TRDGlobalTracking : public o2::framework::Task
   std::unique_ptr<GeometryFlat> mFlatGeo{nullptr};    ///< flat TRD geometry
   bool mUseMC{false};                                 ///< MC flag
   bool mUseTrackletTransform{false};                  ///< if true, output from TrackletTransformer is used instead of uncalibrated Tracklet64 directly
+  CalibVDrift mCalibVDrift{};                         ///< steers the vDrift calibration
   TStopwatch mTimer;
 };
 

Detectors/TRD/workflow/src/TRDGlobalTrackingSpec.cxx

@@ -22,6 +22,7 @@
 #include "DataFormatsTRD/Tracklet64.h"
 #include "DataFormatsTRD/CalibratedTracklet.h"
 #include "DataFormatsTRD/TriggerRecord.h"
+#include "DataFormatsTRD/Constants.h"
 
 // GPU header
 #include "GPUReconstruction.h"
@@ -68,6 +69,7 @@ void TRDGlobalTracking::init(InitContext& ic)
   mTracker->SetProcessPerTimeFrame();
   mTracker->SetNMaxCollisions(mRec->GetProcessingSettings().trdNMaxCollisions);
   mTracker->SetTrkltTransformNeeded(!mUseTrackletTransform);
+  //mTracker->SetDoImpactAngleHistograms(true);
 
   mRec->RegisterGPUProcessor(mTracker, false);
   mChainTracking->SetTRDGeometry(std::move(mFlatGeo));
@@ -99,7 +101,7 @@ void TRDGlobalTracking::run(ProcessingContext& pc)
   int nTrackletsCal = 0;
 
   if (mUseTrackletTransform) {
-    cTrklts.emplace(pc.inputs().get<gsl::span<CalibratedTracklet>>("trdctracklets")); // MC labels associated to the input digits
+    cTrklts.emplace(pc.inputs().get<gsl::span<CalibratedTracklet>>("trdctracklets"));
     cTrkltsPtr = &cTrklts.value();
     nTrackletsCal = cTrkltsPtr->size();
     LOGF(INFO, "Got %i calibrated tracklets as input", nTrackletsCal);
@@ -171,11 +173,18 @@ void TRDGlobalTracking::run(ProcessingContext& pc)
   mTracker->SetTriggerRecordIndices(&(trdTriggerIndices[0]));
   mTracker->SetNCollisions(nCollisions);
   //mTracker->DumpTracks();
+  mTracker->ResetImpactAngleHistograms();
   mTracker->DoTracking(mChainTracking);
   //mTracker->DumpTracks();
 
   std::vector<GPUTRDTrack> tracksOut(mTracker->NTracks());
   std::copy(mTracker->Tracks(), mTracker->Tracks() + mTracker->NTracks(), tracksOut.begin());
+
+  // Temporary until it is transferred to its own DPL device for calibrations
+  mCalibVDrift.setAngleDiffSums(mTracker->AngleDiffSums());
+  mCalibVDrift.setAngleDiffCounters(mTracker->AngleDiffCounters());
+  mCalibVDrift.process();
+
   pc.outputs().snapshot(Output{o2::header::gDataOriginTRD, "MATCHTRD", 0, Lifetime::Timeframe}, tracksOut);
 
   mTimer.Stop();

GPU/GPUTracking/TRDTracking/GPUTRDInterfaces.h

@@ -301,7 +301,7 @@ class trackInterface<GPUTPCGMTrackParam> : public GPUTPCGMTrackParam
   GPUd() const float* getPar() const { return GetPar(); }
   GPUd() const float* getCov() const { return GetCov(); }
   GPUd() float getTime() const { return -1.f; }
-
+  GPUd() void resetCovariance(float s) { ResetCovariance(); }
   GPUd() void setAlpha(float alpha) { mAlpha = alpha; }
   GPUd() void set(float x, float alpha, const float param[5], const float cov[15])
   {