Root/src/RootPatternRecognitionPerformanceWriter.cpp

0001 // This file is part of the ACTS project.
0002 //
0003 // Copyright (C) 2016 CERN for the benefit of the ACTS project
0004 //
0005 // This Source Code Form is subject to the terms of the Mozilla Public
0006 // License, v. 2.0. If a copy of the MPL was not distributed with this
0007 // file, You can obtain one at https://mozilla.org/MPL/2.0/.
0008
0009 #include "ActsExamples/Io/Root/RootPatternRecognitionPerformanceWriter.hpp"
0010
0011 #include "Acts/Utilities/VectorHelpers.hpp"
0012 #include "ActsPlugins/Root/HistogramConverter.hpp"
0013
0014 #include <format>
0015 #include <stdexcept>
0016
0017 #include <TEfficiency.h>
0018 #include <TFile.h>
0019 #include <TH1.h>
0020 #include <TH2.h>
0021 #include <TProfile.h>
0022 #include <TTree.h>
0023 #include <TVectorFfwd.h>
0024 #include <TVectorT.h>
0025
0026 using ActsPlugins::toRoot;
0027
0028 namespace ActsExamples {
0029
0030 namespace {
0031
0032 void writeTrackSummaryPlots(const TrackSummaryPlotTool& tool) {
0033   for (const auto& [name, prof] : tool.profiles()) {
0034     toRoot(prof)->Write();
0035   }
0036 }
0037
0038 }  // namespace
0039
0040 RootPatternRecognitionPerformanceWriter::
0041     RootPatternRecognitionPerformanceWriter(
0042         RootPatternRecognitionPerformanceWriter::Config cfg,
0043         Acts::Logging::Level lvl)
0044     : WriterT(cfg.inputTracks, "RootPatternRecognitionPerformanceWriter", lvl),
0045       m_cfg(std::move(cfg)),
0046       m_collector(
0047           PatternRecognitionPerformanceCollector::Config{
0048               m_cfg.label, m_cfg.effPlotToolConfig, m_cfg.fakePlotToolConfig,
0049               m_cfg.duplicationPlotToolConfig, m_cfg.trackSummaryPlotToolConfig,
0050               m_cfg.trackQualityPlotToolConfig,
0051               m_cfg.subDetectorTrackSummaryVolumes},
0052           logger().clone()) {
0053   // tracks collection name is already checked by base ctor
0054   if (m_cfg.inputParticles.empty()) {
0055     throw std::invalid_argument("Missing particles input collection");
0056   }
0057   if (m_cfg.inputTrackParticleMatching.empty()) {
0058     throw std::invalid_argument("Missing input track particles matching");
0059   }
0060   if (m_cfg.inputParticleTrackMatching.empty()) {
0061     throw std::invalid_argument("Missing input particle track matching");
0062   }
0063   if (m_cfg.inputParticleMeasurementsMap.empty()) {
0064     throw std::invalid_argument("Missing input measurement particles map");
0065   }
0066   if (m_cfg.filePath.empty()) {
0067     throw std::invalid_argument("Missing output filename");
0068   }
0069
0070   m_inputParticles.initialize(m_cfg.inputParticles);
0071   m_inputTrackParticleMatching.initialize(m_cfg.inputTrackParticleMatching);
0072   m_inputParticleTrackMatching.initialize(m_cfg.inputParticleTrackMatching);
0073   m_inputParticleMeasurementsMap.initialize(m_cfg.inputParticleMeasurementsMap);
0074
0075   // the output file can not be given externally since TFile accesses to the
0076   // same file from multiple threads are unsafe.
0077   // must always be opened internally
0078   m_outputFile = TFile::Open(m_cfg.filePath.c_str(), m_cfg.fileMode.c_str());
0079   if (m_outputFile == nullptr) {
0080     throw std::invalid_argument("Could not open '" + m_cfg.filePath + "'");
0081   }
0082
0083   if (m_cfg.writeMatchingDetails) {
0084     m_matchingTree = new TTree("matchingdetails", "matchingdetails");
0085
0086     m_matchingTree->Branch("event_nr", &m_treeEventNr);
0087     m_matchingTree->Branch("particle_id_vertex_primary",
0088                            &m_treeParticleVertexPrimary);
0089     m_matchingTree->Branch("particle_id_vertex_secondary",
0090                            &m_treeParticleVertexSecondary);
0091     m_matchingTree->Branch("particle_id_particle", &m_treeParticleParticle);
0092     m_matchingTree->Branch("particle_id_generation", &m_treeParticleGeneration);
0093     m_matchingTree->Branch("particle_id_sub_particle",
0094                            &m_treeParticleSubParticle);
0095     m_matchingTree->Branch("matched", &m_treeIsMatched);
0096   }
0097 }
0098
0099 RootPatternRecognitionPerformanceWriter::
0100     ~RootPatternRecognitionPerformanceWriter() {
0101   if (m_outputFile != nullptr) {
0102     m_outputFile->Close();
0103   }
0104 }
0105
0106 ProcessCode RootPatternRecognitionPerformanceWriter::finalize() {
0107   m_collector.logSummary();
0108
0109   auto writeFloat = [&](float f, const char* name) {
0110     TVectorF v(1);
0111     v[0] = f;
0112     m_outputFile->WriteObject(&v, name);
0113   };
0114
0115   std::string labelPlural = std::format("{}s", m_cfg.label);
0116
0117   if (m_outputFile != nullptr) {
0118     m_outputFile->cd();
0119
0120     // Write efficiency histograms
0121     for (const auto& [name, eff] : m_collector.effPlotTool().efficiencies1D()) {
0122       toRoot(eff)->Write();
0123     }
0124     for (const auto& [name, eff] : m_collector.effPlotTool().efficiencies2D()) {
0125       toRoot(eff)->Write();
0126     }
0127     for (const auto& eff :
0128          m_collector.effPlotTool().trackEffVsEtaInPtRanges()) {
0129       toRoot(eff)->Write();
0130     }
0131     for (const auto& eff :
0132          m_collector.effPlotTool().trackEffVsPtInAbsEtaRanges()) {
0133       toRoot(eff)->Write();
0134     }
0135
0136     // Write fake ratio histograms
0137     for (const auto& [name, hist] : m_collector.fakePlotTool().histograms()) {
0138       toRoot(hist)->Write();
0139     }
0140     for (const auto& [name, eff] : m_collector.fakePlotTool().efficiencies()) {
0141       toRoot(eff)->Write();
0142     }
0143
0144     // Write duplication ratio histograms
0145     for (const auto& [name, prof] :
0146          m_collector.duplicationPlotTool().profiles()) {
0147       toRoot(prof)->Write();
0148     }
0149     for (const auto& [name, eff] :
0150          m_collector.duplicationPlotTool().efficiencies()) {
0151       toRoot(eff)->Write();
0152     }
0153
0154     // Write track summary histograms
0155     writeTrackSummaryPlots(m_collector.trackSummaryPlotTool());
0156     for (const auto& [key, tool] : m_collector.subDetectorSummaryTools()) {
0157       writeTrackSummaryPlots(tool);
0158     }
0159
0160     // Write track quality histograms
0161     for (const auto& [name, prof] :
0162          m_collector.trackQualityPlotTool().profiles()) {
0163       toRoot(prof)->Write();
0164     }
0165
0166     // Write summary scalars derived from the collector's accumulated counts.
0167     const auto s = m_collector.stats();
0168     float eff_val = static_cast<float>(s.nTotalMatchedTracks) / s.nTotalTracks;
0169     float fakeRatio_val =
0170         static_cast<float>(s.nTotalFakeTracks) / s.nTotalTracks;
0171     float duplicationRatio_val =
0172         static_cast<float>(s.nTotalDuplicateTracks) / s.nTotalTracks;
0173     float eff_particle =
0174         static_cast<float>(s.nTotalMatchedParticles) / s.nTotalParticles;
0175     float fakeRatio_particle =
0176         static_cast<float>(s.nTotalFakeParticles) / s.nTotalParticles;
0177     float duplicationRatio_particle =
0178         static_cast<float>(s.nTotalDuplicateParticles) / s.nTotalParticles;
0179
0180     writeFloat(eff_val, std::format("eff_{}", labelPlural).c_str());
0181     writeFloat(fakeRatio_val, std::format("fakeratio_{}", labelPlural).c_str());
0182     writeFloat(duplicationRatio_val,
0183                std::format("duplicateratio_{}", labelPlural).c_str());
0184     writeFloat(eff_particle, "eff_particles");
0185     writeFloat(fakeRatio_particle, "fakeratio_particles");
0186     writeFloat(duplicationRatio_particle, "duplicateratio_particles");
0187
0188     if (m_matchingTree != nullptr) {
0189       m_matchingTree->Write();
0190     }
0191
0192     ACTS_INFO("Wrote performance plots to '" << m_outputFile->GetPath() << "'");
0193
0194     m_outputFile->Close();
0195     m_outputFile = nullptr;
0196   }
0197   return ProcessCode::SUCCESS;
0198 }
0199
0200 ProcessCode RootPatternRecognitionPerformanceWriter::writeT(
0201     const AlgorithmContext& ctx, const ConstTrackContainer& tracks) {
0202   // Read truth input collections
0203   const auto& particles = m_inputParticles(ctx);
0204   const auto& trackParticleMatching = m_inputTrackParticleMatching(ctx);
0205   const auto& particleTrackMatching = m_inputParticleTrackMatching(ctx);
0206   const auto& particleMeasurementsMap = m_inputParticleMeasurementsMap(ctx);
0207
0208   // Exclusive access to the histograms while filling
0209   std::lock_guard<std::mutex> lock(m_writeMutex);
0210
0211   m_collector.fill(ctx.geoContext, tracks, particles, trackParticleMatching,
0212                    particleTrackMatching, particleMeasurementsMap);
0213
0214   // Write additional matching details to TTree
0215   if (m_cfg.writeMatchingDetails && m_matchingTree != nullptr) {
0216     for (const auto& particle : particles) {
0217       auto particleId = particle.particleId();
0218
0219       m_treeEventNr = ctx.eventNumber;
0220       m_treeParticleVertexPrimary = particleId.vertexPrimary();
0221       m_treeParticleVertexSecondary = particleId.vertexSecondary();
0222       m_treeParticleParticle = particleId.particle();
0223       m_treeParticleGeneration = particleId.generation();
0224       m_treeParticleSubParticle = particleId.subParticle();
0225
0226       m_treeIsMatched = false;
0227       if (auto imatched = particleTrackMatching.find(particleId);
0228           imatched != particleTrackMatching.end()) {
0229         m_treeIsMatched = imatched->second.track.has_value();
0230       }
0231
0232       m_matchingTree->Fill();
0233     }
0234   }
0235
0236   return ProcessCode::SUCCESS;
0237 }
0238
0239 }  // namespace ActsExamples