From 01be99b17d0114525c7fad5cf64bdcc0602d4f9d Mon Sep 17 00:00:00 2001 From: "se.gorbunov" <se.gorbunov@gsi.de> Date: Mon, 20 Nov 2023 23:40:46 +0000 Subject: [PATCH] update of Kf track fit; Bba alignment --- algo/ca/core/tracking/CaTrackFit.cxx | 4 +- macro/mcbm/mcbm_reco_event.C | 1 + macro/run/run_reco.C | 1 + reco/KF/CbmKfFitTracksTask.cxx | 165 ++++++++++++ reco/KF/CbmKfFitTracksTask.h | 62 +++++ ...KFTrackFitter.cxx => CbmKfTrackFitter.cxx} | 239 +++++++++--------- ...{CbmKFTrackFitter.h => CbmKfTrackFitter.h} | 14 +- reco/KF/KF.cmake | 3 +- reco/KF/KFLinkDef.h | 3 +- reco/L1/CbmL1Util.h | 52 ++++ reco/L1/qa/CbmCaTrackFitQa.cxx | 20 +- reco/alignment/CbmBbaAlignmentTask.cxx | 28 +- reco/alignment/CbmBbaAlignmentTask.h | 10 +- 13 files changed, 438 insertions(+), 164 deletions(-) create mode 100644 reco/KF/CbmKfFitTracksTask.cxx create mode 100644 reco/KF/CbmKfFitTracksTask.h rename reco/KF/{CbmKFTrackFitter.cxx => CbmKfTrackFitter.cxx} (72%) rename reco/KF/{CbmKFTrackFitter.h => CbmKfTrackFitter.h} (94%) diff --git a/algo/ca/core/tracking/CaTrackFit.cxx b/algo/ca/core/tracking/CaTrackFit.cxx index 508b7e0fff..20a00699cd 100644 --- a/algo/ca/core/tracking/CaTrackFit.cxx +++ b/algo/ca/core/tracking/CaTrackFit.cxx @@ -194,14 +194,14 @@ namespace cbm::algo::ca mx.SetSinPhi(fvec::Zero()); mx.SetU(mxy.X()); mx.SetDu2(mxy.Dx2()); - mx.SetNdf(fvec::One()); + mx.SetNdf(mxy.NdfX()); ca::MeasurementU<fvec> mu; mu.SetCosPhi(-mxy.Dxy() / mxy.Dx2()); mu.SetSinPhi(fvec::One()); mu.SetU(mu.CosPhi() * mxy.X() + mxy.Y()); mu.SetDu2(mxy.Dy2() - mxy.Dxy() * mxy.Dxy() / mxy.Dx2()); - mu.SetNdf(fvec::One()); + mu.SetNdf(mxy.NdfY()); Filter1d(mx); Filter1d(mu); diff --git a/macro/mcbm/mcbm_reco_event.C b/macro/mcbm/mcbm_reco_event.C index bbb8a21af1..e2000fcfff 100644 --- a/macro/mcbm/mcbm_reco_event.C +++ b/macro/mcbm/mcbm_reco_event.C @@ -395,6 +395,7 @@ void mcbm_reco_event(Int_t nEvents = 10, TString dataset = "data/test", FairTask* globalFindTracks = new CbmL1GlobalFindTracksEvents(globalTrackFinder); run->AddTask(globalFindTracks); + //run->AddTask(new CbmKfFitTracksTask(CbmKfFitTracksTask::FitMode::kMcbm)); // ----- Parameter database -------------------------------------------- std::cout << std::endl << std::endl; diff --git a/macro/run/run_reco.C b/macro/run/run_reco.C index 0f3fcecacb..9caf00b532 100644 --- a/macro/run/run_reco.C +++ b/macro/run/run_reco.C @@ -430,6 +430,7 @@ void run_reco(TString input = "", Int_t nTimeSlices = -1, Int_t firstTimeSlice = run->AddTask(stsFindTracks); std::cout << "-I- " << myName << ": Added task " << stsFindTracks->GetName() << std::endl; } + //run->AddTask(new CbmKfFitTracksTask(CbmKfFitTracksTask::FitMode::kSts)); } // ------------------------------------------------------------------------ diff --git a/reco/KF/CbmKfFitTracksTask.cxx b/reco/KF/CbmKfFitTracksTask.cxx new file mode 100644 index 0000000000..fc66ab460e --- /dev/null +++ b/reco/KF/CbmKfFitTracksTask.cxx @@ -0,0 +1,165 @@ +/* Copyright (C) 2023-2023 GSI Helmholtzzentrum fuer Schwerionenforschung, Darmstadt + SPDX-License-Identifier: GPL-3.0-only + Authors: S.Gorbunov[committer] */ + +/// @file CbmKfFitTracksTask.cxx +/// @author Sergey Gorbunov +/// @date 15.11.2023 +/// @brief Task class for refitting global or sts tracks + + +#include "CbmKfFitTracksTask.h" + +#include "CbmCaTimeSliceReader.h" +#include "CbmGlobalTrack.h" +#include "CbmL1.h" +#include "CbmL1Util.h" +#include "CbmMuchTrack.h" +#include "CbmMuchTrackingInterface.h" +#include "CbmMvdHit.h" +#include "CbmMvdTrackingInterface.h" +#include "CbmStsHit.h" +#include "CbmStsSetup.h" +#include "CbmStsTrack.h" +#include "CbmStsTrackingInterface.h" +#include "CbmTofTrack.h" +#include "CbmTofTrackingInterface.h" +#include "CbmTrdTrack.h" +#include "CbmTrdTrackingInterface.h" + +#include "FairRootManager.h" + +#include "TClonesArray.h" + +#include <iostream> + +#include <cmath> + + +ClassImp(CbmKfFitTracksTask); + +namespace +{ + using namespace cbm::algo; +} + + +CbmKfFitTracksTask::CbmKfFitTracksTask(FitMode mode, Int_t iVerbose) + : FairTask("CbmKfFitTracksTask", iVerbose) + , fFitMode(mode) +{ +} + +CbmKfFitTracksTask::~CbmKfFitTracksTask() {} + + +InitStatus CbmKfFitTracksTask::Init() +{ + + fFitter.Init(); + + //Get ROOT Manager + FairRootManager* ioman = FairRootManager::Instance(); + + if (!ioman) { + LOG(error) << "CbmKfFitTracksTask::Init :: RootManager not instantiated!"; + return kERROR; + } + + // Get global tracks + + fGlobalTracks = dynamic_cast<TClonesArray*>(ioman->GetObject("GlobalTrack")); + + // Get detector tracks + fStsTracks = dynamic_cast<TClonesArray*>(ioman->GetObject("StsTrack")); + fMuchTracks = dynamic_cast<TClonesArray*>(ioman->GetObject("MuchTrack")); + fTrdTracks = dynamic_cast<TClonesArray*>(ioman->GetObject("TrdTrack")); + fTofTracks = dynamic_cast<TClonesArray*>(ioman->GetObject("TofTrack")); + + if (FitMode::kSts != fFitMode && !fGlobalTracks) { + LOG(error) << "CbmKfFitTracksTask::Init: Global track array not found!"; + return kERROR; + } + + if (FitMode::kSts == fFitMode && !fStsTracks) { + LOG(error) << "CbmKfFitTracksTask::Init: Sts track array not found!"; + return kERROR; + } + + return kSUCCESS; +} + +void CbmKfFitTracksTask::Exec(Option_t* /*opt*/) +{ + + LOG(info) << "CbmKfFitTracksTask: exec event N " << fNeventsProcessed++; + + // select tracks for alignment and store them + if (FitMode::kSts == fFitMode && fStsTracks) { + + for (int iTr = 0; iTr < fStsTracks->GetEntriesFast(); iTr++) { + CbmStsTrack* stsTrack = dynamic_cast<CbmStsTrack*>(fStsTracks->At(iTr)); + if (!stsTrack) { + LOG(fatal) << "CbmKfFitTracksTask: null pointer to the sts track!"; + return; + } + CbmKfTrackFitter::Track t; + if (!fFitter.CreateMvdStsTrack(t, *stsTrack)) { + LOG(fatal) << "CbmKfFitTracksTask: can not create the sts track for the fit! "; + return; + } + fFitter.FitTrack(t); + { + const auto& parV = t.fNodes[t.fFirstHitNode].fTrack; + cbm::algo::ca::TrackParamD parD; + parD.Set(parV, 0); + FairTrackParam trackFirst = cbm::L1Util::ConvertTrackParam(parD); + stsTrack->SetParamFirst(&trackFirst); + } + { + const auto& parV = t.fNodes[t.fLastHitNode].fTrack; + cbm::algo::ca::TrackParamD parD; + parD.Set(parV, 0); + FairTrackParam trackLast = cbm::L1Util::ConvertTrackParam(parD); + stsTrack->SetParamLast(&trackLast); + } + } + } + + if (FitMode::kMcbm == fFitMode && fGlobalTracks) { + + for (int iTr = 0; iTr < fGlobalTracks->GetEntriesFast(); iTr++) { + CbmGlobalTrack* globalTrack = dynamic_cast<CbmGlobalTrack*>(fGlobalTracks->At(iTr)); + if (!globalTrack) { + LOG(fatal) << "CbmKfFitTracksTask: null pointer to the global track!"; + return; + } + CbmKfTrackFitter::Track t; + if (!fFitter.CreateGlobalTrack(t, *globalTrack)) { + LOG(fatal) << "CbmKfFitTracksTask: can not create the global track for the fit! "; + return; + } + t.fMsQp0 = 1. / 0.1; + t.fIsMsQp0Set = true; + t.fNodes[t.fFirstHitNode].fTrack.Qp() = 0.; + fFitter.FitTrack(t); + { + const auto& parV = t.fNodes[t.fFirstHitNode].fTrack; + cbm::algo::ca::TrackParamD parD; + parD.Set(parV, 0); + FairTrackParam trackFirst = cbm::L1Util::ConvertTrackParam(parD); + globalTrack->SetParamFirst(&trackFirst); + } + { + const auto& parV = t.fNodes[t.fLastHitNode].fTrack; + cbm::algo::ca::TrackParamD parD; + parD.Set(parV, 0); + FairTrackParam trackLast = cbm::L1Util::ConvertTrackParam(parD); + globalTrack->SetParamLast(&trackLast); + } + } + } +} + + +void CbmKfFitTracksTask::Finish() {} diff --git a/reco/KF/CbmKfFitTracksTask.h b/reco/KF/CbmKfFitTracksTask.h new file mode 100644 index 0000000000..69de5445fb --- /dev/null +++ b/reco/KF/CbmKfFitTracksTask.h @@ -0,0 +1,62 @@ +/* Copyright (C) 2023-2023 GSI Helmholtzzentrum fuer Schwerionenforschung, Darmstadt + SPDX-License-Identifier: GPL-3.0-only + Authors: S.Gorbunov[committer] */ + +/// @file CbmKfFitTracksTask.h +/// @author Sergey Gorbunov +/// @date 15.11.2023 +/// @brief Task class for refitting global or sts tracks + +#pragma once + +#include "CbmKfTrackFitter.h" + +#include "FairTask.h" + +class TClonesArray; + +/// @class CbmKfFitTracksTask +/// @brief Task class for refitting global or sts tracks +/// +class CbmKfFitTracksTask : public FairTask { +public: + enum FitMode + { + kSts, + kMcbm, + kGlobal + }; + + // Constructors/Destructors --------- + CbmKfFitTracksTask(FitMode mode = FitMode::kSts, Int_t iVerbose = 0); + + const CbmKfFitTracksTask& operator=(const CbmKfFitTracksTask&) = delete; + CbmKfFitTracksTask(const CbmKfFitTracksTask&) = delete; + + virtual ~CbmKfFitTracksTask(); + + InitStatus Init() override; + void Exec(Option_t* opt) override; + void Finish() override; + + void SetFitGlobalTracks() { fFitMode = FitMode::kGlobal; } + void SetFitStsTracks() { fFitMode = FitMode::kSts; } + void SetFitMcbmTracks() { fFitMode = FitMode::kMcbm; } + +private: + FitMode fFitMode {FitMode::kGlobal}; ///> fit mode + + /// input data arrays + + TClonesArray* fGlobalTracks {nullptr}; ///< global tracks + TClonesArray* fStsTracks {nullptr}; ///< sts tracks + TClonesArray* fMuchTracks {nullptr}; ///< much tracks + TClonesArray* fTrdTracks {nullptr}; ///< trd tracks + TClonesArray* fTofTracks {nullptr}; ///< tof tracks + + CbmKfTrackFitter fFitter; ///< track fitter + + Int_t fNeventsProcessed {0}; ///< number of processed events + + ClassDefOverride(CbmKfFitTracksTask, 0); +}; diff --git a/reco/KF/CbmKFTrackFitter.cxx b/reco/KF/CbmKfTrackFitter.cxx similarity index 72% rename from reco/KF/CbmKFTrackFitter.cxx rename to reco/KF/CbmKfTrackFitter.cxx index 281f08cd3c..37700a9965 100644 --- a/reco/KF/CbmKFTrackFitter.cxx +++ b/reco/KF/CbmKfTrackFitter.cxx @@ -2,7 +2,7 @@ SPDX-License-Identifier: GPL-3.0-only Authors: Sergey Gorbunov [committer] */ -#include "CbmKFTrackFitter.h" +#include "CbmKfTrackFitter.h" #include "CbmGlobalTrack.h" #include "CbmL1.h" @@ -47,7 +47,7 @@ namespace using namespace cbm::algo; } -void CbmKFTrackFitter::Track::MakeConsistent() +void CbmKfTrackFitter::Track::MakeConsistent() { // sort the nodes in z std::sort(fNodes.begin(), fNodes.end(), [](const FitNode& a, const FitNode& b) { return a.fZ < b.fZ; }); @@ -64,21 +64,21 @@ void CbmKFTrackFitter::Track::MakeConsistent() } -CbmKFTrackFitter::CbmKFTrackFitter() {} +CbmKfTrackFitter::CbmKfTrackFitter() {} -CbmKFTrackFitter::~CbmKFTrackFitter() {} +CbmKfTrackFitter::~CbmKfTrackFitter() {} -void CbmKFTrackFitter::Init() +void CbmKfTrackFitter::Init() { if (fIsInitialized) return; if (!CbmL1::Instance() || !CbmL1::Instance()->fpAlgo) { - LOG(fatal) << "CbmKFTrackFitter: no CbmL1 task initialized "; + LOG(fatal) << "CbmKfTrackFitter: no CbmL1 task initialized "; } FairRootManager* ioman = FairRootManager::Instance(); - if (!ioman) { LOG(fatal) << "CbmKFTrackFitter: no FairRootManager"; } + if (!ioman) { LOG(fatal) << "CbmKfTrackFitter: no FairRootManager"; } // Get hits @@ -100,28 +100,28 @@ void CbmKFTrackFitter::Init() fIsInitialized = true; } -void CbmKFTrackFitter::SetParticleHypothesis(int pdg) +void CbmKfTrackFitter::SetParticleHypothesis(int pdg) { TParticlePDG* particlePDG = TDatabasePDG::Instance()->GetParticle(pdg); if (!particlePDG) { - LOG(fatal) << "CbmKFTrackFitter: particle PDG " << pdg << " is not in the data base, please set the mass manually"; + LOG(fatal) << "CbmKfTrackFitter: particle PDG " << pdg << " is not in the data base, please set the mass manually"; return; } fMass = particlePDG->Mass(); fIsElectron = (abs(pdg) == 11); } -void CbmKFTrackFitter::SetMassHypothesis(double mass) +void CbmKfTrackFitter::SetMassHypothesis(double mass) { assert(mass >= 0.); fMass = mass; } -void CbmKFTrackFitter::SetElectronFlag(bool isElectron) { fIsElectron = isElectron; } +void CbmKfTrackFitter::SetElectronFlag(bool isElectron) { fIsElectron = isElectron; } -bool CbmKFTrackFitter::CreateGlobalTrack(CbmKFTrackFitter::Track& kfTrack, const CbmGlobalTrack& globalTrack) +bool CbmKfTrackFitter::CreateGlobalTrack(CbmKfTrackFitter::Track& kfTrack, const CbmGlobalTrack& globalTrack) { Init(); if (!fIsInitialized) return false; @@ -141,16 +141,16 @@ bool CbmKFTrackFitter::CreateGlobalTrack(CbmKFTrackFitter::Track& kfTrack, const int stsTrackIndex = globalTrack.GetStsTrackIndex(); if (!fInputStsTracks) { - LOG(error) << "CbmKFTrackFitter: Sts track array not found!"; + LOG(error) << "CbmKfTrackFitter: Sts track array not found!"; return false; } if (stsTrackIndex >= fInputStsTracks->GetEntriesFast()) { - LOG(error) << "CbmKFTrackFitter: Sts track index " << stsTrackIndex << " is out of range!"; + LOG(error) << "CbmKfTrackFitter: Sts track index " << stsTrackIndex << " is out of range!"; return false; } auto* stsTrack = dynamic_cast<const CbmStsTrack*>(fInputStsTracks->At(stsTrackIndex)); if (!stsTrack) { - LOG(error) << "CbmKFTrackFitter: Sts track is null!"; + LOG(error) << "CbmKfTrackFitter: Sts track is null!"; return false; } @@ -159,14 +159,14 @@ bool CbmKFTrackFitter::CreateGlobalTrack(CbmKFTrackFitter::Track& kfTrack, const int nMvdHits = stsTrack->GetNofMvdHits(); if (nMvdHits > 0) { if (!fInputMvdHits) { - LOG(error) << "CbmKFTrackFitter: Mvd hit array not found!"; + LOG(error) << "CbmKfTrackFitter: Mvd hit array not found!"; return false; } mvdHits.reserve(nMvdHits); for (int ih = 0; ih < nMvdHits; ih++) { int hitIndex = stsTrack->GetMvdHitIndex(ih); if (hitIndex >= fInputMvdHits->GetEntriesFast()) { - LOG(error) << "CbmKFTrackFitter: Mvd hit index " << hitIndex << " is out of range!"; + LOG(error) << "CbmKfTrackFitter: Mvd hit index " << hitIndex << " is out of range!"; return false; } mvdHits.push_back(*dynamic_cast<const CbmMvdHit*>(fInputMvdHits->At(hitIndex))); @@ -178,14 +178,14 @@ bool CbmKFTrackFitter::CreateGlobalTrack(CbmKFTrackFitter::Track& kfTrack, const int nStsHits = stsTrack->GetNofStsHits(); if (nStsHits > 0) { if (!fInputStsHits) { - LOG(error) << "CbmKFTrackFitter: Sts hit array not found!"; + LOG(error) << "CbmKfTrackFitter: Sts hit array not found!"; return false; } stsHits.reserve(nStsHits); for (int ih = 0; ih < nStsHits; ih++) { int hitIndex = stsTrack->GetStsHitIndex(ih); if (hitIndex >= fInputStsHits->GetEntriesFast()) { - LOG(error) << "CbmKFTrackFitter: Sts hit index " << hitIndex << " is out of range!"; + LOG(error) << "CbmKfTrackFitter: Sts hit index " << hitIndex << " is out of range!"; return false; } stsHits.push_back(*dynamic_cast<const CbmStsHit*>(fInputStsHits->At(hitIndex))); @@ -199,29 +199,29 @@ bool CbmKFTrackFitter::CreateGlobalTrack(CbmKFTrackFitter::Track& kfTrack, const if (globalTrack.GetMuchTrackIndex() >= 0) { int muchTrackIndex = globalTrack.GetMuchTrackIndex(); if (!fInputMuchTracks) { - LOG(error) << "CbmKFTrackFitter: Much track array not found!"; + LOG(error) << "CbmKfTrackFitter: Much track array not found!"; return false; } if (muchTrackIndex >= fInputMuchTracks->GetEntriesFast()) { - LOG(error) << "CbmKFTrackFitter: Much track index " << muchTrackIndex << " is out of range!"; + LOG(error) << "CbmKfTrackFitter: Much track index " << muchTrackIndex << " is out of range!"; return false; } auto* track = dynamic_cast<const CbmMuchTrack*>(fInputMuchTracks->At(muchTrackIndex)); if (!track) { - LOG(error) << "CbmKFTrackFitter: Much track is null!"; + LOG(error) << "CbmKfTrackFitter: Much track is null!"; return false; } int nHits = track->GetNofHits(); if (nHits > 0) { if (!fInputMuchHits) { - LOG(error) << "CbmKFTrackFitter: Much hit array not found!"; + LOG(error) << "CbmKfTrackFitter: Much hit array not found!"; return false; } muchHits.reserve(nHits); for (int ih = 0; ih < nHits; ih++) { int hitIndex = track->GetHitIndex(ih); if (hitIndex >= fInputMuchHits->GetEntriesFast()) { - LOG(error) << "CbmKFTrackFitter: Much hit index " << hitIndex << " is out of range!"; + LOG(error) << "CbmKfTrackFitter: Much hit index " << hitIndex << " is out of range!"; return false; } muchHits.push_back(*dynamic_cast<const CbmMuchPixelHit*>(fInputMuchHits->At(hitIndex))); @@ -234,29 +234,29 @@ bool CbmKFTrackFitter::CreateGlobalTrack(CbmKFTrackFitter::Track& kfTrack, const if (globalTrack.GetTrdTrackIndex() >= 0) { int trdTrackIndex = globalTrack.GetTrdTrackIndex(); if (!fInputTrdTracks) { - LOG(error) << "CbmKFTrackFitter: Trd track array not found!"; + LOG(error) << "CbmKfTrackFitter: Trd track array not found!"; return false; } if (trdTrackIndex >= fInputTrdTracks->GetEntriesFast()) { - LOG(error) << "CbmKFTrackFitter: Trd track index " << trdTrackIndex << " is out of range!"; + LOG(error) << "CbmKfTrackFitter: Trd track index " << trdTrackIndex << " is out of range!"; return false; } auto* track = dynamic_cast<const CbmTrdTrack*>(fInputTrdTracks->At(trdTrackIndex)); if (!track) { - LOG(error) << "CbmKFTrackFitter: Trd track is null!"; + LOG(error) << "CbmKfTrackFitter: Trd track is null!"; return false; } int nHits = track->GetNofHits(); if (nHits > 0) { if (!fInputTrdHits) { - LOG(error) << "CbmKFTrackFitter: Trd hit array not found!"; + LOG(error) << "CbmKfTrackFitter: Trd hit array not found!"; return false; } trdHits.reserve(nHits); for (int ih = 0; ih < nHits; ih++) { int hitIndex = track->GetHitIndex(ih); if (hitIndex >= fInputTrdHits->GetEntriesFast()) { - LOG(error) << "CbmKFTrackFitter: Trd hit index " << hitIndex << " is out of range!"; + LOG(error) << "CbmKfTrackFitter: Trd hit index " << hitIndex << " is out of range!"; return false; } trdHits.push_back(*dynamic_cast<const CbmTrdHit*>(fInputTrdHits->At(hitIndex))); @@ -270,30 +270,30 @@ bool CbmKFTrackFitter::CreateGlobalTrack(CbmKFTrackFitter::Track& kfTrack, const if (globalTrack.GetTofTrackIndex() >= 0) { int tofTrackIndex = globalTrack.GetTofTrackIndex(); if (!fInputTofTracks) { - LOG(error) << "CbmKFTrackFitter: Trd track array not found!"; + LOG(error) << "CbmKfTrackFitter: Trd track array not found!"; return false; } if (tofTrackIndex >= fInputTofTracks->GetEntriesFast()) { - LOG(error) << "CbmKFTrackFitter: Trd track index " << tofTrackIndex << " is out of range!"; + LOG(error) << "CbmKfTrackFitter: Trd track index " << tofTrackIndex << " is out of range!"; return false; } auto* track = dynamic_cast<const CbmTofTrack*>(fInputTofTracks->At(tofTrackIndex)); if (!track) { - LOG(error) << "CbmKFTrackFitter: Tof track is null!"; + LOG(error) << "CbmKfTrackFitter: Tof track is null!"; return false; } int nHits = track->GetNofHits(); if (nHits > 0) { if (!fInputTofHits) { - LOG(error) << "CbmKFTrackFitter: Tof hit array not found!"; + LOG(error) << "CbmKfTrackFitter: Tof hit array not found!"; return false; } tofHits.reserve(nHits); for (int ih = 0; ih < nHits; ih++) { int hitIndex = track->GetHitIndex(ih); if (hitIndex >= fInputTofHits->GetEntriesFast()) { - LOG(error) << "CbmKFTrackFitter: Tof hit index " << hitIndex << " is out of range!"; + LOG(error) << "CbmKfTrackFitter: Tof hit index " << hitIndex << " is out of range!"; return false; } tofHits.push_back(*dynamic_cast<const CbmTofHit*>(fInputTofHits->At(hitIndex))); @@ -305,7 +305,7 @@ bool CbmKFTrackFitter::CreateGlobalTrack(CbmKFTrackFitter::Track& kfTrack, const } -bool CbmKFTrackFitter::CreateMvdStsTrack(CbmKFTrackFitter::Track& kfTrack, const CbmStsTrack& stsTrack) +bool CbmKfTrackFitter::CreateMvdStsTrack(CbmKfTrackFitter::Track& kfTrack, const CbmStsTrack& stsTrack) { Init(); if (!fIsInitialized) return false; @@ -323,7 +323,7 @@ bool CbmKFTrackFitter::CreateMvdStsTrack(CbmKFTrackFitter::Track& kfTrack, const int nMvdHits = stsTrack.GetNofMvdHits(); if (nMvdHits > 0) { if (!fInputMvdHits) { - LOG(error) << "CbmKFTrackFitter: Mvd hit array not found!"; + LOG(error) << "CbmKfTrackFitter: Mvd hit array not found!"; return false; } mvdHits.reserve(nMvdHits); @@ -338,7 +338,7 @@ bool CbmKFTrackFitter::CreateMvdStsTrack(CbmKFTrackFitter::Track& kfTrack, const int nStsHits = stsTrack.GetNofStsHits(); if (nStsHits > 0) { if (!fInputStsHits) { - LOG(error) << "CbmKFTrackFitter: Sts hit array not found!"; + LOG(error) << "CbmKfTrackFitter: Sts hit array not found!"; return false; } stsHits.reserve(nStsHits); @@ -352,7 +352,7 @@ bool CbmKFTrackFitter::CreateMvdStsTrack(CbmKFTrackFitter::Track& kfTrack, const } -bool CbmKFTrackFitter::CreateTrack(CbmKFTrackFitter::Track& kfTrack, const FairTrackParam& trackFirst, +bool CbmKfTrackFitter::CreateTrack(CbmKfTrackFitter::Track& kfTrack, const FairTrackParam& trackFirst, const std::vector<CbmMvdHit>& mvdHits, const std::vector<CbmStsHit>& stsHits, const std::vector<CbmMuchPixelHit>& muchHits, const std::vector<CbmTrdHit>& trdHits, const std::vector<CbmTofHit>& tofHits @@ -363,72 +363,36 @@ bool CbmKFTrackFitter::CreateTrack(CbmKFTrackFitter::Track& kfTrack, const FairT Init(); if (!fIsInitialized) return false; - std::vector<const CbmPixelHit*> hits; - - std::vector<int> hitStations; - const ca::Parameters& caPar = CbmL1::Instance()->fpAlgo->GetParameters(); - for (auto& h : mvdHits) { - hits.push_back(dynamic_cast<const CbmPixelHit*>(&h)); - int stIdx = CbmMvdTrackingInterface::Instance()->GetTrackingStationIndex(&h); - hitStations.push_back(caPar.GetStationIndexActive(stIdx, ca::EDetectorID::kMvd)); - } - - for (auto& h : stsHits) { - hits.push_back(dynamic_cast<const CbmPixelHit*>(&h)); - int stIdx = CbmStsTrackingInterface::Instance()->GetTrackingStationIndex(&h); - hitStations.push_back(caPar.GetStationIndexActive(stIdx, ca::EDetectorID::kSts)); - } - - for (auto& h : muchHits) { - hits.push_back(dynamic_cast<const CbmPixelHit*>(&h)); - int stIdx = CbmMuchTrackingInterface::Instance()->GetTrackingStationIndex(&h); - hitStations.push_back(caPar.GetStationIndexActive(stIdx, ca::EDetectorID::kMuch)); - } - - for (auto& h : trdHits) { - hits.push_back(dynamic_cast<const CbmPixelHit*>(&h)); - int stIdx = CbmTrdTrackingInterface::Instance()->GetTrackingStationIndex(&h); - hitStations.push_back(caPar.GetStationIndexActive(stIdx, ca::EDetectorID::kTrd)); - } - - for (auto& h : tofHits) { - hits.push_back(dynamic_cast<const CbmPixelHit*>(&h)); - int stIdx = CbmTofTrackingInterface::Instance()->GetTrackingStationIndex(&h); - hitStations.push_back(caPar.GetStationIndexActive(stIdx, ca::EDetectorID::kTof)); - } - - CbmKFTrackFitter::Track t; - int nStations = caPar.GetNstationsActive(); - t.fNodes.resize(nStations); - for (int i = 0; i < nStations; i++) { - t.fNodes[i].fMaterialLayer = i; - t.fNodes[i].fZ = caPar.GetStation(i).GetZScal(); - t.fNodes[i].fRadThick = 0.; - t.fNodes[i].fIsRadThickSet = false; - t.fNodes[i].fIsFitted = false; - } - - t.fFirstHitNode = nStations - 1; - t.fLastHitNode = 0; + CbmKfTrackFitter::Track t; - for (unsigned int i = 0; i < hits.size(); i++) { - - assert(hits[i]); - const CbmPixelHit& h = *hits[i]; - - int ista = hitStations[i]; + { + t.fNodes.resize(nStations); + for (int i = 0; i < nStations; i++) { + CbmKfTrackFitter::FitNode& n = t.fNodes[i]; + n = {}; + n.fMaterialLayer = i; + n.fZ = caPar.GetStation(i).GetZScal(); + n.fRadThick = 0.; + n.fIsRadThickSet = false; + n.fIsFitted = false; + n.fIsXySet = false; + n.fIsTimeSet = false; + } + t.fFirstHitNode = nStations - 1; + t.fLastHitNode = 0; + } - if (ista < 0) continue; + // lambda to set the node from the pixel hit + auto setNode = [&](const CbmPixelHit& h, int stIdx, ca::EDetectorID detId) { + int ista = caPar.GetStationIndexActive(stIdx, detId); + if (ista < 0) return; assert(ista < nStations); - - CbmKFTrackFitter::FitNode& n = t.fNodes[ista]; - - n.fZ = h.GetZ(); - + CbmKfTrackFitter::FitNode& n = t.fNodes[ista]; + n.fZ = h.GetZ(); n.fMxy.SetX(h.GetX()); n.fMxy.SetY(h.GetY()); n.fMxy.SetDx2(h.GetDx() * h.GetDx()); @@ -436,17 +400,39 @@ bool CbmKFTrackFitter::CreateTrack(CbmKFTrackFitter::Track& kfTrack, const FairT n.fMxy.SetDxy(h.GetDxy()); n.fMxy.SetNdfX(1); n.fMxy.SetNdfY(1); - + n.fIsXySet = true; n.fMt.SetT(h.GetTime()); n.fMt.SetDt2(h.GetTimeError() * h.GetTimeError()); n.fMt.SetNdfT(1); + n.fIsTimeSet = (detId != ca::EDetectorID::kMvd); + n.fRadThick = 0.; + n.fIsRadThickSet = false; + n.fSystemId = cbm::L1Util::ConvertDetSystemId(detId); + n.fAddress = h.GetAddress(); + }; - n.fRadThick = 0.; + for (auto& h : mvdHits) { + setNode(h, CbmMvdTrackingInterface::Instance()->GetTrackingStationIndex(&h), ca::EDetectorID::kMvd); + } + + for (auto& h : stsHits) { + setNode(h, CbmStsTrackingInterface::Instance()->GetTrackingStationIndex(&h), ca::EDetectorID::kSts); + } + + for (auto& h : muchHits) { + setNode(h, CbmMuchTrackingInterface::Instance()->GetTrackingStationIndex(&h), ca::EDetectorID::kMuch); + } + + for (auto& h : trdHits) { + setNode(h, CbmTrdTrackingInterface::Instance()->GetTrackingStationIndex(&h), ca::EDetectorID::kTrd); + } - t.fFirstHitNode = std::min(t.fFirstHitNode, ista); - t.fLastHitNode = std::max(t.fLastHitNode, ista); + for (auto& h : tofHits) { + setNode(h, CbmTofTrackingInterface::Instance()->GetTrackingStationIndex(&h), ca::EDetectorID::kTof); } + t.MakeConsistent(); + ca::TrackParamD tmp = cbm::L1Util::ConvertTrackParam(trackFirst); t.fNodes[t.fFirstHitNode].fTrack.Set(tmp); @@ -457,7 +443,7 @@ bool CbmKFTrackFitter::CreateTrack(CbmKFTrackFitter::Track& kfTrack, const FairT } -void CbmKFTrackFitter::FilterFirstMeasurement(const FitNode& n) +void CbmKfTrackFitter::FilterFirstMeasurement(const FitNode& n) { // a special routine to filter the first measurement. // the measurement errors are simply copied to the track covariance matrix @@ -473,7 +459,7 @@ void CbmKFTrackFitter::FilterFirstMeasurement(const FitNode& n) tr.SetChiSq(0.); tr.SetChiSqTime(0.); tr.SetNdf(-5 + mxy.NdfX() + mxy.NdfY()); - if (mt.NdfT()[0] > 0) { + if (n.fIsTimeSet) { tr.SetTime(mt.T()); tr.SetC55(mt.Dt2()); tr.SetNdfTime(-2 + 1); @@ -485,30 +471,27 @@ void CbmKFTrackFitter::FilterFirstMeasurement(const FitNode& n) } -void CbmKFTrackFitter::AddMaterialEffects(const CbmKFTrackFitter::Track& t, CbmKFTrackFitter::FitNode& n, +void CbmKfTrackFitter::AddMaterialEffects(const CbmKfTrackFitter::Track& t, CbmKfTrackFitter::FitNode& n, bool upstreamDirection) { // add material effects if (n.fMaterialLayer < 0) { return; } + auto& tr = n.fIsFitted ? n.fTrack : fFit.Tr(); + // calculate the radiation thickness from the current track - if (!n.fIsRadThickSet && !n.fIsFitted) { - n.fRadThick = CbmL1::Instance()->fpAlgo->GetParameters().GetMaterialThicknessScal( - n.fMaterialLayer, fFit.Tr().GetX()[0], fFit.Tr().GetY()[0]); + if (!n.fIsRadThickSet) { + n.fRadThick = + CbmL1::Instance()->fpAlgo->GetParameters().GetMaterialThicknessScal(n.fMaterialLayer, tr.GetX()[0], tr.GetY()[0]); } - fvec msQp0 = t.fMsQp0; - if (!t.fIsMsQp0Set) { - if (n.fIsFitted) { msQp0 = n.fTrack.GetQp(); } - else { - msQp0 = fFit.Tr().GetQp(); - } - } + fvec msQp0 = t.fIsMsQp0Set ? t.fMsQp0 : tr.GetQp(); + fFit.MultipleScattering(n.fRadThick, msQp0); fFit.EnergyLossCorrection(n.fRadThick, upstreamDirection ? fvec::One() : fvec::Zero()); } -void CbmKFTrackFitter::FitTrack(CbmKFTrackFitter::Track& t) +void CbmKfTrackFitter::FitTrack(CbmKfTrackFitter::Track& t) { // fit the track @@ -519,6 +502,7 @@ void CbmKFTrackFitter::FitTrack(CbmKFTrackFitter::Track& t) fFit.SetMask(fmask::One()); fFit.SetParticleMass(fMass); + t.fIsMsQp0Set = false; ca::FieldRegion field _fvecalignment; field.SetUseOriginalField(); @@ -539,15 +523,18 @@ void CbmKFTrackFitter::FitTrack(CbmKFTrackFitter::Track& t) fFit.Extrapolate(n.fZ, field); if (n.fIsFitted) { fFit.SetQp0(n.fTrack.GetQp()); } AddMaterialEffects(t, n, false); + + // store track fitted at the previous node and extrapolated to this node n.fTrack = fFit.Tr(); n.fIsFitted = false; - fFit.FilterXY(n.fMxy); - fFit.FilterTime(n.fMt); - if (iNode == t.fLastHitNode) { n.fTrack = fFit.Tr(); } - if (iNode >= t.fLastHitNode) { n.fIsFitted = true; } + if (n.fIsXySet) { fFit.FilterXY(n.fMxy); } + if (n.fIsTimeSet) { fFit.FilterTime(n.fMt); } + if (iNode >= t.fLastHitNode) { // store track fitted at the last node + n.fTrack = fFit.Tr(); + n.fIsFitted = true; + } } - // fit upstream { FitNode& n = t.fNodes[t.fLastHitNode]; @@ -559,8 +546,9 @@ void CbmKFTrackFitter::FitTrack(CbmKFTrackFitter::Track& t) for (int iNode = t.fLastHitNode - 1; iNode >= 0; iNode--) { FitNode& n = t.fNodes[iNode]; fFit.Extrapolate(n.fZ, field); - fFit.FilterXY(n.fMxy); - fFit.FilterTime(n.fMt); + + if (n.fIsXySet) { fFit.FilterXY(n.fMxy); } + if (n.fIsTimeSet) { fFit.FilterTime(n.fMt); } // combine partially fitted downstream and upstream tracks if (iNode > t.fFirstHitNode) { Smooth(n.fTrack, fFit.Tr()); } @@ -570,7 +558,6 @@ void CbmKFTrackFitter::FitTrack(CbmKFTrackFitter::Track& t) n.fIsFitted = true; fFit.SetQp0(n.fTrack.GetQp()); AddMaterialEffects(t, n, true); - if (iNode == t.fFirstHitNode) { n.fTrack = fFit.Tr(); } } // distribute the final chi2, ndf to all nodes @@ -585,7 +572,7 @@ void CbmKFTrackFitter::FitTrack(CbmKFTrackFitter::Track& t) } -void CbmKFTrackFitter::Smooth(ca::TrackParamV& t1, const ca::TrackParamV& t2) +void CbmKfTrackFitter::Smooth(ca::TrackParamV& t1, const ca::TrackParamV& t2) { // combine two tracks std::tie(t1.X(), t1.Y(), t1.C00(), t1.C10(), t1.C11()) = @@ -619,7 +606,7 @@ void CbmKFTrackFitter::Smooth(ca::TrackParamV& t1, const ca::TrackParamV& t2) t1.C64() = 0.; } -std::tuple<ca::fvec, ca::fvec> CbmKFTrackFitter::Smooth1D(ca::fvec x1, ca::fvec Cxx1, ca::fvec x2, ca::fvec Cxx2) +std::tuple<ca::fvec, ca::fvec> CbmKfTrackFitter::Smooth1D(ca::fvec x1, ca::fvec Cxx1, ca::fvec x2, ca::fvec Cxx2) { // combine two 1D values ca::fvec x = (x1 * Cxx1 + x2 * Cxx2) / (Cxx1 + Cxx2); @@ -628,7 +615,7 @@ std::tuple<ca::fvec, ca::fvec> CbmKFTrackFitter::Smooth1D(ca::fvec x1, ca::fvec } std::tuple<ca::fvec, ca::fvec, ca::fvec, ca::fvec, ca::fvec> -CbmKFTrackFitter::Smooth2D(ca::fvec x1, ca::fvec y1, ca::fvec Cxx1, ca::fvec /*Cxy1*/, ca::fvec Cyy1, ca::fvec x2, +CbmKfTrackFitter::Smooth2D(ca::fvec x1, ca::fvec y1, ca::fvec Cxx1, ca::fvec /*Cxy1*/, ca::fvec Cyy1, ca::fvec x2, ca::fvec y2, ca::fvec Cxx2, ca::fvec /*Cxy2*/, ca::fvec Cyy2) { // combine two 2D values diff --git a/reco/KF/CbmKFTrackFitter.h b/reco/KF/CbmKfTrackFitter.h similarity index 94% rename from reco/KF/CbmKFTrackFitter.h rename to reco/KF/CbmKfTrackFitter.h index 25d88910bb..660c8a4989 100644 --- a/reco/KF/CbmKFTrackFitter.h +++ b/reco/KF/CbmKfTrackFitter.h @@ -40,7 +40,7 @@ namespace /// A fitter for the Cbm tracks /// -class CbmKFTrackFitter { +class CbmKfTrackFitter { public: /// A node on the trajectory where the track parameters are: /// a) measured and / or @@ -67,11 +67,11 @@ public: ca::MeasurementTime<ca::fvec> fMt {}; ///< time measurement at fZ ECbmModuleId fSystemId {ECbmModuleId::kNotExist}; ///< detector system ID of the hit - int fDetectorId {-1}; ///< detector ID of the hit - int fHitId {-1}; ///< hit ID + int fAddress {-1}; ///< detector ID of the hit /// == Flags etc - + bool fIsXySet {false}; ///< true if the XY measurement is set + bool fIsTimeSet {false}; ///< true if the time measurement is set bool fIsFitted {false}; ///< true if the node is fitted, false if the fit failed bool fIsRadThickSet {false}; ///< true if the radiation thickness is set int fReference1 {-1}; ///< some reference can be set by the user @@ -94,8 +94,8 @@ public: void MakeConsistent(); // make the structure fields consistent }; - CbmKFTrackFitter(); - ~CbmKFTrackFitter(); + CbmKfTrackFitter(); + ~CbmKfTrackFitter(); /// initialize the fitter. It must be called in the Init() of the user task. /// when called later, it sees track branches always empty for whatever reason @@ -118,7 +118,7 @@ public: bool CreateGlobalTrack(Track& kfTrack, const CbmGlobalTrack& globalTrack); /// fit the track - void FitTrack(CbmKFTrackFitter::Track& t); + void FitTrack(CbmKfTrackFitter::Track& t); /// fit sts tracks // void FitStsTracks(vector<CbmStsTrack>& Tracks, const vector<int>& pidHypo); diff --git a/reco/KF/KF.cmake b/reco/KF/KF.cmake index af1bd8f834..5b1da2183f 100644 --- a/reco/KF/KF.cmake +++ b/reco/KF/KF.cmake @@ -21,7 +21,8 @@ set(SRCS CbmKFTrackInterface.cxx CbmKFUMeasurement.cxx CbmKFVertexInterface.cxx - CbmKFTrackFitter.cxx + CbmKfTrackFitter.cxx + CbmKfFitTracksTask.cxx #Interface/CbmEcalTrackExtrapolationKF.cxx Interface/CbmKFStsHit.cxx diff --git a/reco/KF/KFLinkDef.h b/reco/KF/KFLinkDef.h index d34ed89840..d2fe9ead38 100644 --- a/reco/KF/KFLinkDef.h +++ b/reco/KF/KFLinkDef.h @@ -40,7 +40,8 @@ #pragma link C++ class CbmL1TofMerger + ; #pragma link C++ class CbmL1TrdTracklet + ; #pragma link C++ class CbmL1TrdTracklet4 + ; -#pragma link C++ class CbmKFTrackFitter + ; +#pragma link C++ class CbmKfTrackFitter + ; +#pragma link C++ class CbmKfFitTracksTask + ; //#pragma link C++ class CbmL1TrdTrackFinderSts+; //#pragma link C++ class CbmL1CATrdTrackFinderSA+; diff --git a/reco/L1/CbmL1Util.h b/reco/L1/CbmL1Util.h index 1e6f88c61b..6d62b6e9ac 100644 --- a/reco/L1/CbmL1Util.h +++ b/reco/L1/CbmL1Util.h @@ -5,6 +5,14 @@ #ifndef CbmL1Util_H #define CbmL1Util_H 1 +#include "CbmDefs.h" +#include "CbmL1DetectorID.h" +#include "CbmMuchTrackingInterface.h" +#include "CbmMvdTrackingInterface.h" +#include "CbmStsTrackingInterface.h" +#include "CbmTofTrackingInterface.h" +#include "CbmTrdTrackingInterface.h" + #include "Rtypes.h" #include "CaTrackParam.h" @@ -23,6 +31,50 @@ namespace cbm::L1Util /// copy Ca track param to fair track param FairTrackParam ConvertTrackParam(const cbm::algo::ca::TrackParamD& t); + /// convert Ca detector Id to Cbm detector Id + inline ECbmModuleId ConvertDetSystemId(const cbm::algo::ca::EDetectorID caDetId) + { + switch (caDetId) { + case cbm::algo::ca::EDetectorID::kMvd: return ECbmModuleId::kMvd; + case cbm::algo::ca::EDetectorID::kSts: return ECbmModuleId::kSts; + case cbm::algo::ca::EDetectorID::kMuch: return ECbmModuleId::kMuch; + case cbm::algo::ca::EDetectorID::kTrd: return ECbmModuleId::kTrd; + case cbm::algo::ca::EDetectorID::kTof: return ECbmModuleId::kTof; + default: return ECbmModuleId::kNotExist; + } + } + + /// convert Cbm detector Id to Ca detector Id + inline cbm::algo::ca::EDetectorID ConvertDetSystemId(const ECbmModuleId cbmDetId) + { + switch (cbmDetId) { + case ECbmModuleId::kMvd: return cbm::algo::ca::EDetectorID::kMvd; + case ECbmModuleId::kSts: return cbm::algo::ca::EDetectorID::kSts; + case ECbmModuleId::kMuch: return cbm::algo::ca::EDetectorID::kMuch; + case ECbmModuleId::kTrd: return cbm::algo::ca::EDetectorID::kTrd; + case ECbmModuleId::kTof: return cbm::algo::ca::EDetectorID::kTof; + default: return cbm::algo::ca::EDetectorID::kEND; + } + } + + inline const CbmTrackingDetectorInterfaceBase* GetTrackingInterface(const cbm::algo::ca::EDetectorID caDetId) + { + switch (caDetId) { + case cbm::algo::ca::EDetectorID::kMvd: return CbmMvdTrackingInterface::Instance(); + case cbm::algo::ca::EDetectorID::kSts: return CbmStsTrackingInterface::Instance(); + case cbm::algo::ca::EDetectorID::kMuch: return CbmMuchTrackingInterface::Instance(); + case cbm::algo::ca::EDetectorID::kTrd: return CbmTrdTrackingInterface::Instance(); + case cbm::algo::ca::EDetectorID::kTof: return CbmTofTrackingInterface::Instance(); + default: return nullptr; + } + } + + inline const CbmTrackingDetectorInterfaceBase* GetTrackingInterface(const ECbmModuleId cbmDetId) + { + cbm::algo::ca::EDetectorID caDetId = ConvertDetSystemId(cbmDetId); + return GetTrackingInterface(caDetId); + } + } // namespace cbm::L1Util #endif diff --git a/reco/L1/qa/CbmCaTrackFitQa.cxx b/reco/L1/qa/CbmCaTrackFitQa.cxx index 715a95fe81..3359040e84 100644 --- a/reco/L1/qa/CbmCaTrackFitQa.cxx +++ b/reco/L1/qa/CbmCaTrackFitQa.cxx @@ -210,17 +210,17 @@ void TrackFitQa::SetDefaultProperties() // ** Pulls distribution parameters ** fvPBins[ETrackParType::kX] = 200; ///< Number of bins, pull of x - fvPLo[ETrackParType::kX] = -4.; ///< Lower boundary, pull of x [cm] - fvPUp[ETrackParType::kX] = +4.; ///< Upper boundary, pull of x [cm] + fvPLo[ETrackParType::kX] = -10.; ///< Lower boundary, pull of x [cm] + fvPUp[ETrackParType::kX] = +10.; ///< Upper boundary, pull of x [cm] fvPBins[ETrackParType::kY] = 200; ///< Number of bins, pull of y - fvPLo[ETrackParType::kY] = -4.; ///< Lower boundary, pull of y [cm] - fvPUp[ETrackParType::kY] = +4.; ///< Upper boundary, pull of y [cm] + fvPLo[ETrackParType::kY] = -10.; ///< Lower boundary, pull of y [cm] + fvPUp[ETrackParType::kY] = +10.; ///< Upper boundary, pull of y [cm] fvPBins[ETrackParType::kTX] = 200; ///< Number of bins, pull of slope along x-axis - fvPLo[ETrackParType::kTX] = -4.; ///< Lower boundary, pull of slope along x-axis - fvPUp[ETrackParType::kTX] = +4.; ///< Upper boundary, pull of slope along x-axis + fvPLo[ETrackParType::kTX] = -10.; ///< Lower boundary, pull of slope along x-axis + fvPUp[ETrackParType::kTX] = +10.; ///< Upper boundary, pull of slope along x-axis fvPBins[ETrackParType::kTY] = 200; ///< Number of bins, pull of slope along y-axis - fvPLo[ETrackParType::kTY] = -4.; ///< Lower boundary, pull of slope along y-axis - fvPUp[ETrackParType::kTY] = +4.; ///< Upper boundary, pull of slope along y-axis + fvPLo[ETrackParType::kTY] = -10.; ///< Lower boundary, pull of slope along y-axis + fvPUp[ETrackParType::kTY] = +10.; ///< Upper boundary, pull of slope along y-axis fvPBins[ETrackParType::kQP] = 200; ///< Number of bins, pull of q/p fvPLo[ETrackParType::kQP] = -10.; ///< Lower boundary, pull of q/p [ec/GeV] fvPUp[ETrackParType::kQP] = +10.; ///< Upper boundary, pull of q/p [ec/GeV] @@ -228,8 +228,8 @@ void TrackFitQa::SetDefaultProperties() fvPLo[ETrackParType::kTIME] = -10.; ///< Lower boundary, pull of time [ns] fvPUp[ETrackParType::kTIME] = +10.; ///< Upper boundary, pull of time [ns] fvPBins[ETrackParType::kVI] = 200; ///< Number of bins, pull of inverse speed - fvPLo[ETrackParType::kVI] = -2.; ///< Lower boundary, pull of inverse speed [1/c] - fvPUp[ETrackParType::kVI] = +2.; ///< Upper boundary, pull of inverse speed [1/c] + fvPLo[ETrackParType::kVI] = -10.; ///< Lower boundary, pull of inverse speed [1/c] + fvPUp[ETrackParType::kVI] = +10.; ///< Upper boundary, pull of inverse speed [1/c] } // --------------------------------------------------------------------------------------------------------------------- diff --git a/reco/alignment/CbmBbaAlignmentTask.cxx b/reco/alignment/CbmBbaAlignmentTask.cxx index aa1197cc50..c7e143158b 100644 --- a/reco/alignment/CbmBbaAlignmentTask.cxx +++ b/reco/alignment/CbmBbaAlignmentTask.cxx @@ -202,6 +202,10 @@ void CbmBbaAlignmentTask::Exec(Option_t* /*opt*/) break; } t.MakeConsistent(); + t.fUnalignedTrack.fMsQp0 = 1. / 0.1; + t.fUnalignedTrack.fIsMsQp0Set = true; + t.fUnalignedTrack.fNodes[t.fUnalignedTrack.fFirstHitNode].fTrack.Qp() = 0.; + fFitter.FitTrack(t.fUnalignedTrack); t.fAlignedTrack = t.fUnalignedTrack; //if (t.fNstsHits < 1) continue; @@ -222,8 +226,8 @@ void CbmBbaAlignmentTask::ApplyAlignment(const std::vector<double>& par) for (TrackContainer& t : fTracks) { for (unsigned int in = 0; in < t.fUnalignedTrack.fNodes.size(); in++) { - CbmKFTrackFitter::FitNode& node = t.fUnalignedTrack.fNodes[in]; - CbmKFTrackFitter::FitNode& nodeAligned = t.fAlignedTrack.fNodes[in]; + CbmKfTrackFitter::FitNode& node = t.fUnalignedTrack.fNodes[in]; + CbmKfTrackFitter::FitNode& nodeAligned = t.fAlignedTrack.fNodes[in]; int iAlignmentModule = node.fReference1; if (iAlignmentModule < 0) continue; assert(nodeAligned.fReference1 == iAlignmentModule); @@ -309,7 +313,7 @@ void CbmBbaAlignmentTask::Finish() p.SetBoundaryMax(2.); p.SetStepMin(1.e-6); p.SetStepMax(0.1); - p.SetStep(10.e-4); + p.SetStep(100.e-4); } par[0].SetActive(0); // fix the first station @@ -328,7 +332,11 @@ void CbmBbaAlignmentTask::Finish() par[3 * i + 2].SetActive(0); } - par[3 * 1 + 0].SetActive(1); + for (int i = 1; i < nStations - 1; i++) { + par[3 * i + 0].SetActive(1); + par[3 * i + 1].SetActive(1); + par[3 * i + 2].SetActive(1); + } gRandom->SetSeed(1); @@ -336,14 +344,9 @@ void CbmBbaAlignmentTask::Finish() bba::Parameter& px = par[3 * is + 0]; bba::Parameter& py = par[3 * is + 1]; bba::Parameter& pz = par[3 * is + 2]; - //py.SetActive(1); - //pz.SetStepMin(1.e-4); - // +- 0.5 cm misalignment - //if (px.IsActive()) { px.SetValue(.1); } //gRandom->Uniform(2 * d) - d); } - //if (py.IsActive()) { py.SetValue(.1); } //gRandom->Uniform(2. * d) - d); } - //if (pz.IsActive()) { pz.SetValue(.1); } //gRandom->Uniform(2. * d) - d); } - double d = 0.5; + // +- 0.5 cm misalignment + double d = 10.; if (px.IsActive()) { px.SetValue(gRandom->Uniform(2. * d) - d); } if (py.IsActive()) { py.SetValue(gRandom->Uniform(2. * d) - d); } if (pz.IsActive()) { pz.SetValue(gRandom->Uniform(2. * d) - d); } @@ -432,8 +435,9 @@ void CbmBbaAlignmentTask::WriteHistosCurFile(TObject* obj) sub->cd(); TList* listSub = (static_cast<TDirectory*>(obj))->GetList(); TIter it(listSub); - while (TObject* obj1 = it()) + while (TObject* obj1 = it()) { WriteHistosCurFile(obj1); + } cur->cd(); gFile = currentFile; gDirectory = cur; diff --git a/reco/alignment/CbmBbaAlignmentTask.h b/reco/alignment/CbmBbaAlignmentTask.h index f2e91a5ddd..dd03c17f61 100644 --- a/reco/alignment/CbmBbaAlignmentTask.h +++ b/reco/alignment/CbmBbaAlignmentTask.h @@ -13,7 +13,7 @@ #define CbmBbaAlignmentTask_H -#include "CbmKFTrackFitter.h" +#include "CbmKfTrackFitter.h" #include "FairTask.h" @@ -61,11 +61,11 @@ public: }; /// information about a track - /// aligned and unaligned hits are stored in the corresponding CbmKFTrackFitter::Track objects + /// aligned and unaligned hits are stored in the corresponding CbmKfTrackFitter::Track objects struct TrackContainer { - CbmKFTrackFitter::Track fUnalignedTrack; // track before alignment - CbmKFTrackFitter::Track fAlignedTrack; // track after alignment + CbmKfTrackFitter::Track fUnalignedTrack; // track before alignment + CbmKfTrackFitter::Track fAlignedTrack; // track after alignment int fNmvdHits {0}; // number of MVD hits int fNstsHits {0}; // number of STS hits @@ -98,7 +98,7 @@ private: TClonesArray* fInputGlobalTrackMatches {nullptr}; // Mc info for debugging TClonesArray* fInputStsTrackMatches {nullptr}; // Mc info for debugging - CbmKFTrackFitter fFitter; + CbmKfTrackFitter fFitter; // collection of selected tracks and hits std::vector<TrackContainer> fTracks; -- GitLab