CbmLitTrackFinderNN.cxx

Go to the documentation of this file.

/* Copyright (C) 2008-2017 GSI/JINR-LIT, Darmstadt/Dubna
   SPDX-License-Identifier: GPL-3.0-only
   Authors: Andrey Lebedev [committer] */
 
#include "finder/CbmLitTrackFinderNN.h"
 
#include "data/CbmLitHit.h"
#include "data/CbmLitPixelHit.h"
#include "data/CbmLitTrack.h"
#include "data/CbmLitTrackParam.h"
#include "interface/CbmLitTrackPropagator.h"
#include "interface/CbmLitTrackSelection.h"
#include "interface/CbmLitTrackUpdate.h"
#include "utils/CbmLitMath.h"
#include "utils/CbmLitMemoryManagment.h"
 
#include <algorithm>
#include <cmath>
#include <iostream>
#include <limits>
#include <map>
#include <set>
 
using std::make_pair;
using std::map;
using std::set;
using std::sqrt;
 
CbmLitTrackFinderNN::CbmLitTrackFinderNN()
  : fTracks()
  , fFilter()
  , fSeedSelection()
  , fFinalSelection()
  , fPropagator()
  , fNofStations(0)
  , fNofIterations(0)
  , fIteration(0)
  , fHitData()
  , fUsedHitsSet()
  , fUsedSeedsSet()
  , fMaxNofMissingHits()
  , fPDG()
  , fChiSqStripHitCut()
  , fChiSqPixelHitCut()
  , fSigmaCoef()
{
}
 
CbmLitTrackFinderNN::~CbmLitTrackFinderNN() {}
 
LitStatus CbmLitTrackFinderNN::DoFind(HitPtrVector& hits, TrackPtrVector& trackSeeds, TrackPtrVector& tracks)
{
  fTracks.clear();
  fUsedSeedsSet.clear();
  fUsedHitsSet.clear();
  fHitData.SetNofStations(fNofStations);
 
  for (fIteration = 0; fIteration < fNofIterations; fIteration++) {
    // std::cout << "CbmLitTrackFinderNN::DoFind: iteration=" << fIteration << std::endl;
    ArrangeHits(hits.begin(), hits.end());
    // std::cout << fHitData.ToString();
 
    InitTrackSeeds(trackSeeds.begin(), trackSeeds.end());
    FollowTracks(fTracks.begin(), fTracks.end());
 
    fFinalSelection->DoSelect(fTracks.begin(), fTracks.end());
    RemoveHits(fTracks.begin(), fTracks.end());
    CopyToOutput(fTracks.begin(), fTracks.end(), tracks);
 
    for_each(fTracks.begin(), fTracks.end(), DeleteObject());
    fTracks.clear();
    fHitData.Clear();
  }
 
  //std::cout << "tracks.size()=" << tracks.size() << std::endl;
  //for (Int_t i = 0; i < tracks.size(); i++) {
  //   std::cout << i << " " << tracks[i]->ToString();
  //}
 
  static Int_t eventNo = 0;
  std::cout << "CbmLitTrackFinderNN::DoFind: " << eventNo++ << " events processed" << std::endl;
  return kLITSUCCESS;
}
 
void CbmLitTrackFinderNN::ArrangeHits(HitPtrIterator itBegin, HitPtrIterator itEnd)
{
  for (HitPtrIterator it = itBegin; it != itEnd; it++) {
    CbmLitHit* hit = *it;
    if (fUsedHitsSet.find(hit->GetRefId()) != fUsedHitsSet.end()) {
      continue;
    }  // FIXME Check if it works correctly for MUCH-TRD
    fHitData.AddHit(hit);
  }
  fHitData.Arrange();
}
 
void CbmLitTrackFinderNN::InitTrackSeeds(TrackPtrIterator itBegin, TrackPtrIterator itEnd)
{
  for (TrackPtrIterator it = itBegin; it != itEnd; it++) {
    (*it)->SetQuality(kLITGOOD);
  }
 
  fSeedSelection->DoSelect(itBegin, itEnd);
 
  for (TrackPtrIterator it = itBegin; it != itEnd; it++) {
    CbmLitTrack* track = *it;
    if (track->GetQuality() == kLITBAD) {
      continue;
    }
    if (fUsedSeedsSet.find(track->GetPreviousTrackId()) != fUsedSeedsSet.end()) {
      continue;
    }
    track->SetPDG(fPDG[fIteration]);
    track->SetChi2(0.);
    fTracks.push_back(new CbmLitTrack(*track));
  }
}
 
void CbmLitTrackFinderNN::FollowTracks(TrackPtrIterator itBegin, TrackPtrIterator itEnd)
{
  for (TrackPtrIterator itTrack = itBegin; itTrack != itEnd; itTrack++) {
    CbmLitTrack* track = *itTrack;
 
    for (Int_t iStation = 0; iStation < fNofStations; iStation++) {
      CbmLitTrackParam par(*track->GetParamLast());
 
      litfloat zMin = fHitData.GetMinZPos(iStation);
      if (fPropagator->Propagate(&par, zMin, fPDG[fIteration]) == kLITERROR) {
        break;
      }
 
      const vector<Int_t>& bins = fHitData.GetZPosBins(iStation);
      // map<bin index, pair<track parameter for the bin, true if track was propagated correctly >>
      map<Int_t, pair<CbmLitTrackParam, Bool_t>> binParamMap;
      vector<Int_t>::const_iterator itBins;
      for (itBins = bins.begin(); itBins != bins.end(); itBins++) {
        binParamMap[*itBins] = make_pair<CbmLitTrackParam, Bool_t>(CbmLitTrackParam(), true);
      }
 
      // Extrapolate track parameters to each Z position in the map.
      // This is done to improve calculation speed.
      // In case of planar station only 1 track extrapolation is required,
      // since all hits located at the same Z.
      map<Int_t, pair<CbmLitTrackParam, Bool_t>>::iterator itMap;
      for (itMap = binParamMap.begin(); itMap != binParamMap.end(); itMap++) {
        (*itMap).second.first = par;
        litfloat z            = fHitData.GetZPosByBin(iStation, (*itMap).first);
        if (fPropagator->Propagate(&(*itMap).second.first, z, fPDG[fIteration]) == kLITERROR) {
          (*itMap).second.second = false;
        }
      }
 
      // Loop over hits
      litfloat minChiSq       = std::numeric_limits<litfloat>::max();  // minimum chi-square of hit
      const CbmLitHit* minHit = NULL;                                  // Pointer to hit with minimum chi-square
      CbmLitTrackParam minPar;                                         // Track parameters for closest hit
      const HitPtrVector& hits = fHitData.GetHits(iStation);
      for (HitPtrConstIterator itHit = hits.begin(); itHit != hits.end(); itHit++) {
        const CbmLitHit* hit = *itHit;
        Int_t bin            = fHitData.GetBinByZPos(iStation, hit->GetZ());
        assert(binParamMap.find(bin) != binParamMap.end());
        if (!binParamMap[bin].second) continue;  // Track parameters are wrong for this propagation
        CbmLitTrackParam tpar(binParamMap[bin].first);
 
        // Check preliminary if hit is in the validation gate.
        // This is done in order to speed up the algorithm.
        // Based on the predicted track position (w/o KF update step)
        // and maximum hit position error.
        if (hit->GetType() == kLITPIXELHIT) {
          const CbmLitPixelHit* pixelHit = static_cast<const CbmLitPixelHit*>(hit);
          litfloat maxErrX               = fHitData.GetMaxErrX(iStation);
          litfloat devX                  = fSigmaCoef[fIteration] * (sqrt(tpar.GetCovariance(0) + maxErrX * maxErrX));
          litfloat maxErrY               = fHitData.GetMaxErrY(iStation);
          litfloat devY                  = fSigmaCoef[fIteration] * (sqrt(tpar.GetCovariance(6) + maxErrY * maxErrY));
          litfloat maxErrT               = fHitData.GetMaxErrT(iStation);
          litfloat devT                  = fSigmaCoef[fIteration] * (sqrt(tpar.GetCovariance(20) + maxErrT * maxErrT));
          bool hitInside = (pixelHit->GetX() < (tpar.GetX() + devX)) && (pixelHit->GetX() > (tpar.GetX() - devX))
                           && (pixelHit->GetY() < (tpar.GetY() + devY)) && (pixelHit->GetY() > (tpar.GetY() - devY))
                           && (pixelHit->GetT() < (tpar.GetTime() + devT))
                           && (pixelHit->GetT() > (tpar.GetTime() - devT));
          if (!hitInside) continue;
        }
 
        litfloat chi = std::numeric_limits<litfloat>::max();
        fFilter->Update(&tpar, hit, chi);
        bool hitInValidationGate = (hit->GetType() == kLITPIXELHIT && chi < fChiSqPixelHitCut[fIteration])
                                   || (hit->GetType() == kLITSTRIPHIT && chi < fChiSqStripHitCut[fIteration]);
        if (hitInValidationGate && chi < minChiSq) {  // Check if hit is inside validation gate and closer to the track.
          minChiSq = chi;
          minHit   = hit;
          minPar   = tpar;
        }
      }
 
      if (minHit != NULL) {  // Check if hit was added
        track->AddHit(minHit);
        track->SetLastStationId(iStation);
        track->SetParamLast(&minPar);
        track->SetChi2(track->GetChi2() + minChiSq);
        track->SetNDF(lit::NDF(track));
      }
      else {  // Missing hit
        track->SetNofMissingHits(track->GetNofMissingHits() + 1);
        if (track->GetNofMissingHits() > fMaxNofMissingHits[fIteration]) {
          break;
        }
      }
    }
  }
}
 
void CbmLitTrackFinderNN::RemoveHits(TrackPtrIterator itBegin, TrackPtrIterator itEnd)
{
  for (TrackPtrIterator it = itBegin; it != itEnd; it++) {
    CbmLitTrack* track = *it;
    if (track->GetQuality() == kLITBAD) {
      continue;
    }
    for (Int_t hit = 0; hit < track->GetNofHits(); hit++) {
      fUsedHitsSet.insert(track->GetHit(hit)->GetRefId());
    }
  }
}
 
void CbmLitTrackFinderNN::CopyToOutput(TrackPtrIterator itBegin, TrackPtrIterator itEnd, TrackPtrVector& tracks)
{
  for (TrackPtrIterator it = itBegin; it != itEnd; it++) {
    CbmLitTrack* track = *it;
    if (track->GetQuality() == kLITBAD) {
      continue;
    }
    if (!track->CheckParams()) {
      continue;
    }
    fUsedSeedsSet.insert(track->GetPreviousTrackId());
    tracks.push_back(new CbmLitTrack(*track));
  }
}