/* Copyright (C) 2022 GSI Helmholtzzentrum fuer Schwerionenforschung, Darmstadt SPDX-License-Identifier: GPL-3.0-only Authors: Sergey Gorbunov, Sergei Zharko [committer] */ /// \file L1InitManager.cxx /// \brief Input parameters management class for L1Algo /// \since 19.01.2022 #include "L1InitManager.h" #include <boost/archive/binary_iarchive.hpp> #include <boost/archive/binary_oarchive.hpp> #include <boost/archive/text_iarchive.hpp> #include <algorithm> #include <fstream> #include <sstream> #include "L1Algo.h" #include "L1Assert.h" // ---------------------------------------------------------------------------------------------------------------------- // void L1InitManager::AddStation(const L1BaseStationInfo& inStation) { // Check if other fields were defined already // Active detector IDs L1MASSERT(0, fInitController.GetFlag(EInitKey::kActiveDetectorIDs), "Attempt to add a station info before the active detectors set had been initialized"); // Number of stations check L1MASSERT(0, fInitController.GetFlag(EInitKey::kStationsNumberCrosscheck), "Attempt to add a station info before the numbers of stations for each detector had been initialized"); // Field function L1MASSERT(0, fInitController.GetFlag(EInitKey::kFieldFunction), "Attempt to add a station info before the magnetic field function had been initialized"); // Check activeness of this station type bool isStationActive = inStation.GetTrackingStatus() && fActiveDetectorIDs.find(inStation.GetDetectorID()) != fActiveDetectorIDs.end(); if (isStationActive) { // initialize magnetic field slice L1BaseStationInfo inStationCopy = L1BaseStationInfo(inStation); // make a copy of station so it can be initialized inStationCopy.SetFieldFunction(fFieldFunction); // check, if material map is set if (!inStationCopy.GetInitController().GetFlag(L1BaseStationInfo::EInitKey::kThicknessMap)) { LOG(fatal) << "Station material map was not set for detectorID = " << static_cast<int>(inStationCopy.GetDetectorID()) << ", stationID = " << inStationCopy.GetStationID(); } // Check station init LOG(debug) << "L1InitManager::AddStation:(original) L1BaseStationInfo " << inStation.GetInitController().ToString(); LOG(debug) << "L1InitManager::AddStation:(copy) L1BaseStationInfo " << inStation.GetInitController().ToString(); if (!inStationCopy.GetInitController().IsFinalized()) { LOG(fatal) << "Attempt to add an incompletely initialized station info object (detectorID = " << static_cast<int>(inStationCopy.GetDetectorID()) << ", stationID = " << inStationCopy.GetStationID() << ")"; } // insert the station in a set auto insertionResult = fStationsInfo.insert(std::move(inStationCopy)); if (!insertionResult.second) { LOG(fatal) << "Attempt to add a duplicated station info object (detectorID = " << static_cast<int>(inStationCopy.GetDetectorID()) << ", stationID = " << inStationCopy.GetStationID() << ")"; } int index = fStationsInfo.size() - 1 + fParameters.fNstationsGeometryTotal - fParameters.fNstationsActiveTotal; fParameters.fActiveStationGlobalIDs[index] = fStationsInfo.size() - 1; } else { int index = fStationsInfo.size() + fParameters.fNstationsGeometryTotal - fParameters.fNstationsActiveTotal; fParameters.fActiveStationGlobalIDs[index] = -1; fParameters.fNstationsActive[static_cast<L1DetectorID_t>(inStation.GetDetectorID())]--; fParameters.fNstationsActiveTotal--; } LOG(debug) << "L1InitManager: adding a station with stationID = " << inStation.GetStationID() << " and detectorID = " << static_cast<int>(inStation.GetDetectorID()) << ". Is active: " << isStationActive; } // ---------------------------------------------------------------------------------------------------------------------- // void L1InitManager::CheckInit() { this->CheckCAIterationsInit(); this->CheckStationsInfoInit(); } // ---------------------------------------------------------------------------------------------------------------------- // void L1InitManager::ClearSetupInfo() { // Clear stations set and a thickness map fStationsInfo.clear(); std::fill(fParameters.fThickMap.begin(), fParameters.fThickMap.end(), L1Material()); fInitController.SetFlag(EInitKey::kStationsInfo, false); // Set number of stations do default values std::fill(fParameters.fNstationsGeometry.begin(), fParameters.fNstationsGeometry.end(), 0); std::fill(fParameters.fNstationsActive.begin(), fParameters.fNstationsActive.end(), 0); std::fill(fParameters.fActiveStationGlobalIDs.begin(), fParameters.fActiveStationGlobalIDs.end(), 0); fParameters.fNstationsGeometryTotal = -1; fParameters.fNstationsActiveTotal = -1; fInitController.SetFlag(EInitKey::kStationsNumberCrosscheck, false); // Clear active detectors fActiveDetectorIDs.clear(); fInitController.SetFlag(EInitKey::kActiveDetectorIDs, false); // Clear field info fParameters.fVertexFieldRegion = L1FieldRegion(); fParameters.fVertexFieldValue = L1FieldValue(); fInitController.SetFlag(EInitKey::kPrimaryVertexField, false); // Clear target position std::fill(fParameters.fTargetPos.begin(), fParameters.fTargetPos.end(), L1Utils::kNaN); fTargetZ = 0.; fInitController.SetFlag(EInitKey::kTargetPos, false); // Clear field function fFieldFunction = L1FieldFunction_t([](const double(&)[3], double(&)[3]) {}); fInitController.SetFlag(EInitKey::kFieldFunction, false); // Clear other flags fParameters.fTrackingLevel = 0; fParameters.fGhostSuppression = 0; fParameters.fMomentumCutOff = 0; fInitController.SetFlag(EInitKey::kTrackingLevel, false); fInitController.SetFlag(EInitKey::kGhostSuppression, false); fInitController.SetFlag(EInitKey::kMomentumCutOff, false); fParameters.fDevIsIgnoreHitSearchAreas = false; fParameters.fDevIsUseOfOriginalField = false; fParameters.fDevIsMatchDoubletsViaMc = false; fParameters.fDevIsMatchTripletsViaMc = false; fParameters.fDevIsExtendTracksViaMc = false; fParameters.fDevIsSuppressOverlapHitsViaMc = false; } // ---------------------------------------------------------------------------------------------------------------------- // void L1InitManager::ClearCAIterations() { fParameters.fCAIterations.clear(); fCAIterationsNumberCrosscheck = -1; fInitController.SetFlag(EInitKey::kCAIterations, false); fInitController.SetFlag(EInitKey::kCAIterationsNumberCrosscheck, false); } // ---------------------------------------------------------------------------------------------------------------------- // NOTE: this function should be called once in the SendParameters bool L1InitManager::FormParametersContainer() { // Read configuration file // NOTE: We consider parameters from the configuration file as ones with a higher priority, so all the defined // variables there would be rewritten by the configuration if (fConfigInputName != "") { try { fConfigRW.ReadYaml(fConfigInputName); LOG(info) << "L1InitManager: parameters configuration file \"" << fConfigInputName << "\" was loaded"; } catch (const std::runtime_error& err) { LOG(error) << "L1InitManager: parameters configuration file \"" << fConfigInputName << "\" was defined, " << "but the loading failed. Reason: " << err.what(); return false; } } if (!fParameters.fDevIsParSearchWUsed) { fInitController.SetFlag(EInitKey::kSearchWindows, true); } // Check initialization this->CheckInit(); if (!fInitController.IsFinalized()) { LOG(error) << "L1InitManager: Attempt to form parameters container before all necessary fields were initialized" << fInitController.ToString(); return false; } // Form array of stations auto destinationArrayIterator = fParameters.fStations.begin(); for (const auto& item : fStationsInfo) { *destinationArrayIterator = item.GetL1Station(); ++destinationArrayIterator; } // Form array of material map auto thickMapIt = fParameters.fThickMap.begin(); for (auto it = fStationsInfo.begin(); it != fStationsInfo.end(); ++it) { auto node = fStationsInfo.extract(it); *thickMapIt = std::move(node.value().TakeMaterialMap()); fStationsInfo.insert(std::move(node)); ++thickMapIt; } // Check the consistency of the parameters object. If object inconsistent, it throws std::logic_error try { fParameters.CheckConsistency(); } catch (const std::logic_error& err) { LOG(error) << "L1InitManager: parameters container consistency check failed. Reason: " << err.what(); return false; } return true; } // ---------------------------------------------------------------------------------------------------------------------- // void L1InitManager::InitTargetField(double zStep) { if (fInitController.GetFlag(EInitKey::kPrimaryVertexField)) { LOG(warn) << "L1InitManager::InitTargetField: attempt to reinitialize the field value and field region " << "near target. Ignore"; return; } // Check for field function L1MASSERT(0, fInitController.GetFlag(EInitKey::kFieldFunction), "Attempt to initialize the field value and field region near target before initializing field function"); // Check for target defined L1MASSERT( 0, fInitController.GetFlag(EInitKey::kTargetPos), "Attempt to initialize the field value and field region near target before the target position initialization"); constexpr int nDimensions {3}; constexpr int nPointsNodal {3}; std::array<double, nPointsNodal> inputNodalZ {fTargetZ, fTargetZ + zStep, fTargetZ + 2. * zStep}; std::array<L1FieldValue, nPointsNodal> B {}; std::array<fvec, nPointsNodal> z {}; // loop over nodal points for (int idx = 0; idx < nPointsNodal; ++idx) { double point[nDimensions] {0., 0., inputNodalZ[idx]}; double field[nDimensions] {}; fFieldFunction(point, field); z[idx] = inputNodalZ[idx]; B[idx].x = field[0]; B[idx].y = field[1]; B[idx].z = field[2]; } // loop over nodal points: end fParameters.fVertexFieldRegion.Set(B[0], z[0], B[1], z[1], B[2], z[2]); fParameters.fVertexFieldValue = B[0]; fInitController.SetFlag(EInitKey::kPrimaryVertexField); } // ---------------------------------------------------------------------------------------------------------------------- // void L1InitManager::PushBackCAIteration(const L1CAIteration& iteration) { // TODO: probably some checks must be inserted here (S.Zharko) bool control = fInitController.GetFlag(EInitKey::kCAIterationsNumberCrosscheck); L1MASSERT(0, control, //fInitController.GetFlag(EInitKey::kCAIterationsNumberCrosscheck), "Attempt to push back a CA track finder iteration before the number of iterations was defined"); fParameters.fCAIterations.push_back(iteration); } // --------------------------------------------------------------------------------------------------------------------- // void L1InitManager::ReadParametersObject(const std::string& fileName) { // Open input binary file std::ifstream ifs(fileName); if (!ifs) { LOG(fatal) << "L1InitManager: parameters data file \"" << fileName << "\" was not found"; } // Get L1InputData object try { boost::archive::text_iarchive ia(ifs); ia >> fParameters; } catch (const std::exception&) { LOG(fatal) << "L1InitManager: parameters file \"" << fileName << "\" has incorrect data format or was corrupted"; } } // --------------------------------------------------------------------------------------------------------------------- // void L1InitManager::ReadSearchWindows(const std::string& fileName) { // Open input binary file std::ifstream ifs(fileName); if (!ifs) { LOG(fatal) << "L1InitManager: search window file \"" << fileName << "\" was not found"; } try { boost::archive::text_iarchive ia(ifs); int nPars = -1; int nWindows = -1; ia >> nPars; assert(nPars == 1); // Currently only the constant windows are available ia >> nWindows; std::stringstream errMsg; for (int iW = 0; iW < nWindows; ++iW) { L1SearchWindow swBuffer; ia >> swBuffer; int iStationID = swBuffer.GetStationID(); int iTrackGrID = swBuffer.GetTrackGroupID(); if (iStationID < 0 || iStationID > L1Constants::size::kMaxNstations) { errMsg << "\t- wrong station id for entry " << iW << ": " << iStationID << " (should be between 0 and " << L1Constants::size::kMaxNstations << ")\n"; } if (iTrackGrID < 0 || iTrackGrID > L1Constants::size::kMaxNtrackGroups) { errMsg << "\t- wrong track group id for entry " << iW << ": " << iTrackGrID << " (should be between 0 and " << L1Constants::size::kMaxNtrackGroups << ")\n"; } fParameters.fSearchWindows[iTrackGrID * L1Constants::size::kMaxNstations + iStationID] = swBuffer; } if (errMsg.str().size()) { LOG(fatal) << "L1InitManager: some errors occurred while reading search windows: " << errMsg.str(); } } catch (const std::exception&) { LOG(fatal) << "L1InitManager: search windows file \"" << fileName << "\" has incorrect data format or was corrupted"; } fInitController.SetFlag(EInitKey::kSearchWindows, true); } // --------------------------------------------------------------------------------------------------------------------- // bool L1InitManager::SendParameters(L1Algo* pAlgo) { assert(pAlgo); pAlgo->ReceiveParameters(std::move(fParameters)); return true; } // ---------------------------------------------------------------------------------------------------------------------- // void L1InitManager::SetActiveDetectorIDs(const L1DetectorIDSet_t& detectorIDs) { // TODO: To think about redefinition possibilities: should it be allowed or not? (S.Zh.) fActiveDetectorIDs = detectorIDs; fInitController.SetFlag(EInitKey::kActiveDetectorIDs); } // ---------------------------------------------------------------------------------------------------------------------- // void L1InitManager::SetCAIterationsNumberCrosscheck(int nIterations) { fCAIterationsNumberCrosscheck = nIterations; L1Vector<L1CAIteration>& iterationsContainer = fParameters.fCAIterations; // NOTE: should be called to prevent multiple copies of objects between the memory reallocations iterationsContainer.reserve(nIterations); fInitController.SetFlag(EInitKey::kCAIterationsNumberCrosscheck); } // ---------------------------------------------------------------------------------------------------------------------- // void L1InitManager::SetFieldFunction(const L1FieldFunction_t& fieldFunction) { if (!fInitController.GetFlag(EInitKey::kFieldFunction)) { fFieldFunction = fieldFunction; fInitController.SetFlag(EInitKey::kFieldFunction); } else { LOG(warn) << "L1InitManager::SetFieldFunction: attempt to reinitialize the field function. Ignored"; } } // ---------------------------------------------------------------------------------------------------------------------- // void L1InitManager::SetGhostSuppression(int ghostSuppression) { if (fInitController.GetFlag(EInitKey::kGhostSuppression)) { LOG(warn) << "L1InitManager::SetGhostSuppression: attempt of reinitializating the ghost suppresion flag. Ignore"; return; } fParameters.fGhostSuppression = ghostSuppression; fInitController.SetFlag(EInitKey::kGhostSuppression); } // ---------------------------------------------------------------------------------------------------------------------- // void L1InitManager::SetMomentumCutOff(float momentumCutOff) { if (fInitController.GetFlag(EInitKey::kMomentumCutOff)) { LOG(warn) << "L1InitManager::SetMomentumCutOff: attempt of reinitializating the momentum cutoff value. Ignore"; return; } fParameters.fMomentumCutOff = momentumCutOff; fInitController.SetFlag(EInitKey::kMomentumCutOff); } // ---------------------------------------------------------------------------------------------------------------------- // void L1InitManager::SetNstations(L1DetectorID detectorID, int nStations) { L1MASSERT(0, fInitController.GetFlag(EInitKey::kActiveDetectorIDs), "Attempt to set crosscheck number of stations before the active detetors set had been initialized"); // NOTE: We add and check only those detectors which will be active (?) // For INACTIVE detectors the initialization code for it inside CbmL1/BmnL1 can (and must) be still in, // but it will be ignored inside L1InitManager. if (fActiveDetectorIDs.find(detectorID) != fActiveDetectorIDs.end()) { if (nStations) { fParameters.fNstationsGeometry[static_cast<L1DetectorID_t>(detectorID)] = nStations; fParameters.fNstationsActive[static_cast<L1DetectorID_t>(detectorID)] = nStations; } else { // TODO: Probably it is better to replace fatal with warn and remove the detectorID from active detectors (S.Zharko) LOG(fatal) << "L1InitManager::SetNstations: attempt to initialize zero stations for active detector: " << static_cast<L1DetectorID_t>(detectorID); } } // Check if all the station numbers for active detectors are initialized now: LOG(debug) << "L1InitManager::SetNstations called for detectorID = " << static_cast<int>(detectorID); if (!fInitController.GetFlag(EInitKey::kStationsNumberCrosscheck)) { bool ifInitialized = true; for (auto item : fActiveDetectorIDs) { if (fParameters.fNstationsGeometry[static_cast<L1DetectorID_t>(item)] == 0) { ifInitialized = false; break; } } fInitController.SetFlag(EInitKey::kStationsNumberCrosscheck, ifInitialized); } if (fInitController.GetFlag(EInitKey::kStationsNumberCrosscheck)) { fParameters.fNstationsGeometryTotal = std::accumulate(fParameters.fNstationsGeometry.begin(), fParameters.fNstationsGeometry.end(), 0); fParameters.fNstationsActiveTotal = fParameters.fNstationsGeometryTotal; } } // ---------------------------------------------------------------------------------------------------------------------- // void L1InitManager::SetTargetPosition(double x, double y, double z) { if (fInitController.GetFlag(EInitKey::kTargetPos)) { LOG(warn) << "L1InitManager::SetTargetPosition: attempt to reinitialize the target position. Ignore"; return; } /// Fill fvec target fields fParameters.fTargetPos[0] = x; fParameters.fTargetPos[1] = y; fParameters.fTargetPos[2] = z; /// Set additional field for z component in double precision fTargetZ = z; fInitController.SetFlag(EInitKey::kTargetPos); } // ---------------------------------------------------------------------------------------------------------------------- // void L1InitManager::SetTrackingLevel(int trackingLevel) { if (fInitController.GetFlag(EInitKey::kTrackingLevel)) { LOG(warn) << "L1InitManager::SetTrackingLevel: attempt of reinitialization the tracking level. Ignore"; return; } fParameters.fTrackingLevel = trackingLevel; fInitController.SetFlag(EInitKey::kTrackingLevel); } // ---------------------------------------------------------------------------------------------------------------------- // void L1InitManager::WriteParametersObject(const std::string& fileName) const { // Open output binary file std::ofstream ofs(fileName, std::ios::binary); if (!ofs) { LOG(error) << "L1InitManager: failed opening file \"" << fileName << " for writing parameters object\""; return; } // Serialize L1Parameters object and write boost::archive::binary_oarchive oa(ofs); oa << fParameters; } // // INIT CHECKERS // // ---------------------------------------------------------------------------------------------------------------------- // void L1InitManager::CheckCAIterationsInit() { // // 1) Check number of iterations // bool ifInitPassed = true; if (!fInitController.GetFlag(EInitKey::kCAIterations)) { int nIterationsActual = fParameters.fCAIterations.size(); int nIterationsExpected = fCAIterationsNumberCrosscheck; if (nIterationsActual != nIterationsExpected) { LOG(warn) << "L1InitManager::CheckCAIterations: incorrect number of iterations registered: " << nIterationsActual << " of " << nIterationsExpected << " expected"; ifInitPassed = false; } } fInitController.SetFlag(EInitKey::kCAIterations, ifInitPassed); } // ---------------------------------------------------------------------------------------------------------------------- // TODO: REWRITE! and add const qualifier (S.Zharko) void L1InitManager::CheckStationsInfoInit() { bool ifInitPassed = true; if (!fInitController.GetFlag(EInitKey::kStationsInfo)) { // // 1) Check numbers of stations passed // // loop over active detectors for (auto itemDetector : fActiveDetectorIDs) { auto selectDetector = [&itemDetector](const L1BaseStationInfo& station) { return station.GetDetectorID() == itemDetector; }; int nStationsExpected = GetNstationsActive(itemDetector); int nStations = std::count_if(fStationsInfo.begin(), fStationsInfo.end(), selectDetector); if (nStations != nStationsExpected) { LOG(error) << "L1InitManager::CheckStationsInfoInit: Incorrect number of L1BaseStationInfo objects passed" << " to the L1Manager for L1DetectorID = " << static_cast<int>(itemDetector) << ": " << nStations << " of " << nStationsExpected << " expected"; ifInitPassed = false; } } // loop over active detectors: end L1MASSERT(0, ifInitPassed, "Station info initialization failed"); // // 2) Check for maximum allowed number of stations // int nStationsTotal = fParameters.fNstationsGeometryTotal; if (nStationsTotal > L1Constants::size::kMaxNstations) { LOG(fatal) << "Actual total number of registered stations in geometry (" << nStationsTotal << ") is larger then possible (" << L1Constants::size::kMaxNstations << "). Please, select another set of active tracking detectors or recompile the code with enlarged" << " L1Constants::size::kMaxNstations value"; } } fInitController.SetFlag(EInitKey::kStationsInfo, ifInitPassed); }