TrTrack.cpp

Go to the documentation of this file.

#include "fmt/format.h"
#include <bitset>

#include "Containers/TrTrack.h"

namespace NAIA {
const std::string TrTrackBase::BranchName = "TrTrackBaseData";
const std::string TrTrackPlus::BranchName = "TrTrackPlusData";

using namespace TrTrack;

float ApplyGeomCorrection(const TrackFitVariable<float> &Rigidity, Fit fit, Span span, TrGeometry trGeom) {
  if (!ContainsKeys(Rigidity, fit, span)) {
    return 0;
  }

  static constexpr float corr = 0.0993e-3;

  switch (trGeom) {
  case TrGeometry::Default:
    if (span == Span::InnerOnly || span == Span::InnerL1)
      break;
    else
      return Rigidity.at(fit).at(span);
  case TrGeometry::InnerOnly:
  case TrGeometry::InnerL1:
    break;
  default:
    return Rigidity.at(fit).at(span);
  }

  return 1.0 / (1.0 / Rigidity.at(fit).at(span) + corr);
}

float TrTrackBaseData::GetGeomCorrectedRigidityCorr(Fit fit, Span span, TrGeometry trGeom) const {
  return ApplyGeomCorrection(RigidityCorr, fit, span, trGeom);
}

float TrTrackPlusData::GetGeomCorrectedRigidity(Fit fit, Span span, TrGeometry trGeom) const {
  return ApplyGeomCorrection(Rigidity, fit, span, trGeom);
}

float TrTrackPlusData::GetGeomCorrectedRigidityTOI(Fit fit, Span span, TrGeometry trGeom) const {
  return ApplyGeomCorrection(RigidityTOI, fit, span, trGeom);
}

bool TrTrackBaseData::FitIDExists(Fit fit, Span span) const {
  return (RigidityCorr.find(fit) != end(RigidityCorr) && RigidityCorr.at(fit).find(span) != end(RigidityCorr.at(fit)));
}

Span TrTrackBaseData::GetBestSpan(Fit fit) const {
  Span bestspan = FitIDExists(fit, Span::FullSpan)
                      ? Span::FullSpan
                      : (FitIDExists(fit, Span::InnerL1) > 0
                             ? Span::InnerL1
                             : (FitIDExists(fit, Span::InnerL9) > 0 ? Span::InnerL9 : Span::InnerOnly));
  return bestspan;
}

void TrTrackBaseData::Clear() {
  RigidityCorr.clear();
  TrChiSq.clear();

  Charge.clear();
  InnerCharge.clear();
  LayerChargeXY.clear();
  LayerChargeStatus.clear();

  _beta = 1.0f;
  _fit_ID.clear();
}

void TrTrackPlusData::Clear() {
  TrTrackHitPos.clear();
  TrackFeetDistance.clear();

  Rigidity.clear();
  RigidityTOI.clear();
  InvRigErr.clear();
  PartialRigidity.clear();
  PartialInvRigErr.clear();
  PartialTrChiSq.clear();

  TrTrackResidual.clear();
  TrTrackFitPos.clear();

  UnbiasedCharge.clear();
  ChargeRMS.clear();
  InnerChargeRMS.clear();

  LayerCharge.clear();

  LayerEdep.clear();

  NClusters.clear();
  ClustersEdep.clear();
  MaxClusterEdep.clear();
  MaxClusterDistance.clear();
  ClusterSignalRatio.clear();

  _fit_ID.clear();
}

void TrTrackBaseData::Dump() const {
  fmt::print("{:=^120} \n", " TrTrackBaseData ");
  if (!RigidityCorr.empty()) {
    for (const auto &fitIter : RigidityCorr) {
      fmt::print("Fit {}\n", TrTrack::fitNames[fitIter.first]);
      fmt::print("  {:<30}   {:>15} {:>15} {:>15}\n", "Span", "RigidityCorr", "TrChiSqX", "TrChiSqY");
      for (auto spanIter : fitIter.second) {
        fmt::print("  - {:<28} = {:> 12.4f} GV {:> 15.4f} {:> 15.4f}\n", TrTrack::spanNames[spanIter.first],
                   spanIter.second, TrChiSq.at(fitIter.first).at(spanIter.first).at(Side::X),
                   TrChiSq.at(fitIter.first).at(spanIter.first).at(Side::Y));
      }
    }
  }

  if (!Charge.empty()) {
    fmt::print("{:-^120}\n", "");
    fmt::print("{:<32}   {:>15} {:>15} {:>15} {:>32}\n", "Charge type", "Standard", "Hu Liu", "Yi Jia", "status");

    fmt::print("  - {:<28} = ", "Global charge");
    for (auto chType : TrTrack::recoTypes) {
      fmt::print("{:> 15.4f} ", ContainsKeys(Charge, chType) ? Charge.at(chType) : 0);
    }
    fmt::print("\n");

    fmt::print("  - {:<28} = ", "Inner tracker charge");
    for (auto chType : TrTrack::recoTypes) {
      fmt::print("{:> 15.4f} ", ContainsKeys(InnerCharge, chType) ? InnerCharge.at(chType) : 0);
    }
    fmt::print("\n");

    for (size_t ijl = 0; ijl < 9; ijl++) {
      fmt::print("  - {:<28} = ", fmt::format("Layer {}", ijl));
      for (auto chType : TrTrack::recoTypes) {
        if (ContainsKeys(LayerChargeXY, ijl, chType))
          fmt::print("{:> 15.4f} ", LayerChargeXY.at(ijl).at(chType));
        else
          fmt::print("{:>15} ", "");
      }
      std::bitset<32> dummy(ContainsKeys(LayerChargeStatus, ijl) ? LayerChargeStatus.at(ijl) : 0xFFFF);
      fmt::print("{:>32} ", dummy.to_string());
      fmt::print("\n");
    }
  }

  fmt::print("{:-^120}\n", "");
}

void TrTrackPlusData::Dump() const {
  fmt::print("{:=^120} \n", " TrTrackPlusData ");
  if (!Rigidity.empty()) {
    for (const auto &fitIter : Rigidity) {
      fmt::print("Fit {}\n", TrTrack::fitNames[fitIter.first]);
      if (RigidityTOI.find(fitIter.first) == end(RigidityTOI)) {
        fmt::print("  {:<30}   {:>15} {:>15} {:>15}\n", "Span", "Rigidity", "", "InvRigErr");
        for (auto spanIter : fitIter.second) {
          fmt::print("  - {:<28} = {:> 12.4f} GV {:>15} {:>15}\n", TrTrack::spanNames[spanIter.first], spanIter.second,
                     "", InvRigErr.at(fitIter.first).at(spanIter.first));
        }
      } else {
        fmt::print("  {:<30}   {:>15} {:>15} {:>15}\n", "Span", "Rigidity", "RigidityTOI", "InvRigErr");
        for (auto spanIter : fitIter.second) {
          fmt::print("  - {:<28} = {:> 12.4f} GV {:> 12.4f} GV {:>15}\n", TrTrack::spanNames[spanIter.first],
                     spanIter.second, RigidityTOI.at(fitIter.first).at(spanIter.first),
                     InvRigErr.at(fitIter.first).at(spanIter.first));
        }
      }
    }
  }

  fmt::print("{:-^120}\n", "");
  fmt::print("{:<32}   {:>15} {:>15} {:>15} {:>15}\n", "\"Partial\" variables", "Rigidity", "InvRigErr", "TrChiSqX",
             "TrChiSqY");
  for (unsigned int ijl = 0; ijl < 9; ijl++) {
    if (!ContainsKeys(PartialRigidity, ijl))
      continue;

    for (auto &fitIter : PartialRigidity.at(ijl)) {
      auto fit = fitIter.first;

      fmt::print("  - {:<28} = {:>15} {:>15} {:>15} {:>15}\n",
                 fmt::format("Layer {} - Fit {}", ijl + 1, TrTrack::fitNames[fit]), PartialRigidity.at(ijl).at(fit),
                 PartialInvRigErr.at(ijl).at(fit), PartialTrChiSq.at(ijl).at(fit).at(TrTrack::Side::X),
                 PartialTrChiSq.at(ijl).at(fit).at(TrTrack::Side::Y));
    }
  }

  if (!ChargeRMS.empty()) {
    fmt::print("{:-^120}\n", "");
    fmt::print("{:<32}   {:>15} {:>15} {:>15}\n", "Charge type", "Standard", "Hu Liu", "Yi Jia");

    fmt::print("  - {:<28} = ", "Global charge RMS");
    for (auto chType : TrTrack::recoTypes) {
      fmt::print("{:> 15.4f} ", ContainsKeys(ChargeRMS, chType) ? ChargeRMS.at(chType) : 0);
    }
    fmt::print("\n");

    fmt::print("  - {:<28} = ", "Inner tracker charge RMS");
    for (auto chType : TrTrack::recoTypes) {
      fmt::print("{:> 15.4f} ", ContainsKeys(InnerChargeRMS, chType) ? InnerChargeRMS.at(chType) : 0);
    }
    fmt::print("\n");

    for (size_t ijl = 0; ijl < 9; ijl++) {
      fmt::print("  - {:<28} = ", fmt::format("Layer {} X", ijl));
      for (auto chType : TrTrack::recoTypes) {
        if (ContainsKeys(LayerCharge, ijl, chType, TrTrack::Side::X))
          fmt::print("{:> 15.4f} ", LayerCharge.at(ijl).at(chType).at(TrTrack::Side::X));
        else
          fmt::print("{:>15} ", "");
      }
      fmt::print("\n");
      fmt::print("  - {:<28} = ", fmt::format("Layer {} Y", ijl));
      for (auto chType : TrTrack::recoTypes) {
        if (ContainsKeys(LayerCharge, ijl, chType, TrTrack::Side::Y))
          fmt::print("{:> 15.4f} ", LayerCharge.at(ijl).at(chType).at(TrTrack::Side::Y));
        else
          fmt::print("{:>15} ", "");
      }
      fmt::print("\n");
    }
  }
  fmt::print("{:-^120}\n", "");
}
} // namespace NAIA