From 90c146b3b5ffbbcd19e63622fb569201506300bd Mon Sep 17 00:00:00 2001 From: Victor Schwan Date: Mon, 23 Jun 2025 17:35:48 +0200 Subject: [PATCH 01/38] init --- Tracking/components/TrackD0Printer.cpp | 69 ++++++++++++++++++++++++++ Tracking/options/printD0.py | 23 +++++++++ 2 files changed, 92 insertions(+) create mode 100644 Tracking/components/TrackD0Printer.cpp create mode 100644 Tracking/options/printD0.py diff --git a/Tracking/components/TrackD0Printer.cpp b/Tracking/components/TrackD0Printer.cpp new file mode 100644 index 00000000..ef208297 --- /dev/null +++ b/Tracking/components/TrackD0Printer.cpp @@ -0,0 +1,69 @@ +#include "Gaudi/Property.h" +#include "edm4hep/Track.h" +#include "edm4hep/TrackCollection.h" +#include "edm4hep/TrackState.h" +#include "k4FWCore/Consumer.h" +#include "podio/RelationRange.h" + +#include +#include +#include +#include + +// Type alias for TrackCollection to improve readability +using TrackColl = edm4hep::TrackCollection; + +// Consumer that processes track collections and prints D0 values +struct TrackD0Printer final : k4FWCore::Consumer { + // Constructor: define the input collections (ClupatraTracks and SiTracks) + TrackD0Printer(const std::string& name, ISvcLocator* svcLoc) + : Consumer(name, svcLoc, + { + KeyValues("InputSiTracks", {"SiTracks"}), + KeyValues("InputCluTracks", {"ClupatraTracks"}), + }) {} + + // This function will be called to process the data + void operator()(const TrackColl& inSiTracks, const TrackColl& inCluTracks) const override { + + debug() << "Received SiTracks collection with " << inSiTracks.size() << " tracks" << endmsg; + debug() << "Received ClupatraTracks collection with " << inCluTracks.size() << " tracks" << endmsg; + + // Check that both collections are of the same size + if (inSiTracks.size() != inCluTracks.size()) { + fatal() << "Track collections have different sizes: SiTracks (" << inSiTracks.size() << ") and ClupatraTracks (" + << inCluTracks.size() << "). Exiting!" << endmsg; + return; + } + + // Process both collections (SiTracks and ClupatraTracks) at the same time + for (size_t i = 0; i < inSiTracks.size(); ++i) { + + // Process SiTrack + processTrack(inSiTracks[i], "SiTrack"); + + // Process CluTrack + processTrack(inCluTracks[i], "CluTrack"); + } + } + +private: + // Optional: You can add a property to filter or adjust behavior if needed (e.g., track state index). + Gaudi::Property m_trackStateIndex{this, "TrackStateIndex", 2, "Index of track state to print (default 2)"}; + + // Helper function to process each track (either SiTrack or CluTrack) + void processTrack(const edm4hep::Track& track, const std::string& trackType) const { + auto trackStates = track.getTrackStates(); // RelationRange + + // Check if there are enough track states (e.g., third track state) + if (trackStates.size() > 2) { + const edm4hep::TrackState& state = trackStates[2]; // Get the third track state + std::cout << trackType << " D0: " << state.D0 << std::endl; // Print D0 value + } else { + warning() << trackType << " has less than 3 track states, skipping D0 print" << endmsg; + } + } +}; + +// Declare the consumer component +DECLARE_COMPONENT(TrackD0Printer) \ No newline at end of file diff --git a/Tracking/options/printD0.py b/Tracking/options/printD0.py new file mode 100644 index 00000000..aff33370 --- /dev/null +++ b/Tracking/options/printD0.py @@ -0,0 +1,23 @@ +from pathlib import Path + +from Configurables import TrackD0Printer +from Gaudi.Configuration import INFO +from k4FWCore import ApplicationMgr, IOSvc + +iosvc = IOSvc() +iosvc.Input = str( + Path.home() + / "promotion/code/ILDConfig/StandardConfig/production/data/test_tracking_3_detmods_V02_REC.edm4hep.root" +) + +iosvc.CollectionNames = ["SiTracks", "ClupatraTracks"] + +printer = TrackD0Printer("TrackD0Printer") + +ApplicationMgr( + TopAlg=[printer], + EvtSel="NONE", + EvtMax=10, + ExtSvc=[iosvc], + OutputLevel=INFO, +) From a1fa01b68ad0f8755e52a53cd039c973467375ff Mon Sep 17 00:00:00 2001 From: Victor Schwan Date: Tue, 24 Jun 2025 14:28:28 +0200 Subject: [PATCH 02/38] switch from print D0 to phi --- Tracking/components/TrackD0Printer.cpp | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/Tracking/components/TrackD0Printer.cpp b/Tracking/components/TrackD0Printer.cpp index ef208297..4bfa03cb 100644 --- a/Tracking/components/TrackD0Printer.cpp +++ b/Tracking/components/TrackD0Printer.cpp @@ -57,10 +57,10 @@ struct TrackD0Printer final : k4FWCore::Consumer 2) { - const edm4hep::TrackState& state = trackStates[2]; // Get the third track state - std::cout << trackType << " D0: " << state.D0 << std::endl; // Print D0 value + const edm4hep::TrackState& state = trackStates[2]; // Get the third track state + std::cout << trackType << " phi: " << state.phi << std::endl; // Print phi value } else { - warning() << trackType << " has less than 3 track states, skipping D0 print" << endmsg; + warning() << trackType << " has less than 3 track states, skipping phi print" << endmsg; } } }; From 35fdd04b43a739cec3044c2c0fb04767524a4506 Mon Sep 17 00:00:00 2001 From: Victor Schwan Date: Tue, 24 Jun 2025 16:13:28 +0200 Subject: [PATCH 03/38] print error on phi, info() instead of std::cout --- Tracking/components/TrackD0Printer.cpp | 16 +++++++++++++--- 1 file changed, 13 insertions(+), 3 deletions(-) diff --git a/Tracking/components/TrackD0Printer.cpp b/Tracking/components/TrackD0Printer.cpp index 4bfa03cb..f3c90870 100644 --- a/Tracking/components/TrackD0Printer.cpp +++ b/Tracking/components/TrackD0Printer.cpp @@ -1,10 +1,12 @@ #include "Gaudi/Property.h" +#include "GaudiKernel/MsgStream.h" #include "edm4hep/Track.h" #include "edm4hep/TrackCollection.h" #include "edm4hep/TrackState.h" #include "k4FWCore/Consumer.h" #include "podio/RelationRange.h" +#include #include #include #include @@ -12,6 +14,7 @@ // Type alias for TrackCollection to improve readability using TrackColl = edm4hep::TrackCollection; +using TP = edm4hep::TrackParams; // Consumer that processes track collections and prints D0 values struct TrackD0Printer final : k4FWCore::Consumer { @@ -25,6 +28,9 @@ struct TrackD0Printer final : k4FWCore::Consumer m_trackStateIndex{this, "TrackStateIndex", 2, "Index of track state to print (default 2)"}; + float getSigmaPhi(const edm4hep::TrackState& ts) const { return std::sqrt(ts.getCovMatrix(TP::phi, TP::phi)); } + // Helper function to process each track (either SiTrack or CluTrack) void processTrack(const edm4hep::Track& track, const std::string& trackType) const { auto trackStates = track.getTrackStates(); // RelationRange // Check if there are enough track states (e.g., third track state) + std::string varName = "phi"; if (trackStates.size() > 2) { - const edm4hep::TrackState& state = trackStates[2]; // Get the third track state - std::cout << trackType << " phi: " << state.phi << std::endl; // Print phi value + const edm4hep::TrackState& state = trackStates[2]; // Get the third track state + info() << trackType << std::string(7, ' ') + varName + ": " << state.phi << endmsg; // Print phi value + info() << trackType << " sigma " + varName + ": " << getSigmaPhi(state) << endmsg; } else { - warning() << trackType << " has less than 3 track states, skipping phi print" << endmsg; + warning() << trackType << " has less than 3 track states, skipping " + varName + " print" << endmsg; } } }; From 49a1c8fb6a4f96e4ea67a35e2bf73b45a232d5bf Mon Sep 17 00:00:00 2001 From: Victor Schwan Date: Tue, 24 Jun 2025 18:59:09 +0200 Subject: [PATCH 04/38] start switch to fmt::format and outsource value + uncertainty printing into func --- Tracking/components/TrackD0Printer.cpp | 17 ++++++++++++++--- 1 file changed, 14 insertions(+), 3 deletions(-) diff --git a/Tracking/components/TrackD0Printer.cpp b/Tracking/components/TrackD0Printer.cpp index f3c90870..54d5d429 100644 --- a/Tracking/components/TrackD0Printer.cpp +++ b/Tracking/components/TrackD0Printer.cpp @@ -5,7 +5,9 @@ #include "edm4hep/TrackState.h" #include "k4FWCore/Consumer.h" #include "podio/RelationRange.h" +#include +#include #include #include #include @@ -59,16 +61,25 @@ struct TrackD0Printer final : k4FWCore::Consumer{}}{:>10.5f}", strTrType, vName, maxLabelWidth, varValue); + } + // Helper function to process each track (either SiTrack or CluTrack) void processTrack(const edm4hep::Track& track, const std::string& trackType) const { auto trackStates = track.getTrackStates(); // RelationRange // Check if there are enough track states (e.g., third track state) std::string varName = "phi"; + std::string sigmaVarName = "sigma " + varName; + int maxVarWidth = std::max(varName.size(), sigmaVarName.size()) + 2; + + // Check if there are enough track states (e.g., third track state) if (trackStates.size() > 2) { - const edm4hep::TrackState& state = trackStates[2]; // Get the third track state - info() << trackType << std::string(7, ' ') + varName + ": " << state.phi << endmsg; // Print phi value - info() << trackType << " sigma " + varName + ": " << getSigmaPhi(state) << endmsg; + const edm4hep::TrackState& state = trackStates[2]; // Get the third track state + info() << printValueUnc(trackType, varName, maxVarWidth, state.phi) << endmsg; + info() << printValueUnc(trackType, sigmaVarName, maxVarWidth, getSigmaPhi(state)) << endmsg; } else { warning() << trackType << " has less than 3 track states, skipping " + varName + " print" << endmsg; } From 097f4a4ea4167c4dd70213ebe6f0491269fe9d9b Mon Sep 17 00:00:00 2001 From: Victor Schwan Date: Tue, 24 Jun 2025 18:59:37 +0200 Subject: [PATCH 05/38] avoid auto type --- Tracking/components/TrackD0Printer.cpp | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/Tracking/components/TrackD0Printer.cpp b/Tracking/components/TrackD0Printer.cpp index 54d5d429..9105b309 100644 --- a/Tracking/components/TrackD0Printer.cpp +++ b/Tracking/components/TrackD0Printer.cpp @@ -68,9 +68,9 @@ struct TrackD0Printer final : k4FWCore::Consumer - // Check if there are enough track states (e.g., third track state) + podio::RelationRange trackStates = track.getTrackStates(); + std::string varName = "phi"; std::string sigmaVarName = "sigma " + varName; int maxVarWidth = std::max(varName.size(), sigmaVarName.size()) + 2; From 965a104ca4db10ce3c3535d42b543a10cdc94817 Mon Sep 17 00:00:00 2001 From: Victor Schwan Date: Wed, 25 Jun 2025 11:00:38 +0200 Subject: [PATCH 06/38] new func printInStars && move to fmt::format --- Tracking/components/TrackD0Printer.cpp | 11 ++++++++--- 1 file changed, 8 insertions(+), 3 deletions(-) diff --git a/Tracking/components/TrackD0Printer.cpp b/Tracking/components/TrackD0Printer.cpp index 9105b309..0c898f6f 100644 --- a/Tracking/components/TrackD0Printer.cpp +++ b/Tracking/components/TrackD0Printer.cpp @@ -30,9 +30,8 @@ struct TrackD0Printer final : k4FWCore::Consumer m_trackStateIndex{this, "TrackStateIndex", 2, "Index of track state to print (default 2)"}; + void printInStars(const std::string& msg, const int lineWidth) const { + info() << fmt::format("{:*^{}}", "", lineWidth) << endmsg; + info() << fmt::format("{:*^{}}", msg, lineWidth) << endmsg; + info() << fmt::format("{:*^{}}", "", lineWidth) << endmsg; + } + float getSigmaPhi(const edm4hep::TrackState& ts) const { return std::sqrt(ts.getCovMatrix(TP::phi, TP::phi)); } std::string printValueUnc(const std::string& strTrType, const std::string& vName, const int maxLabelWidth, From c7fd24a00fa8f34212c9979702db4944f66c3e69 Mon Sep 17 00:00:00 2001 From: Victor Schwan Date: Wed, 25 Jun 2025 11:06:09 +0200 Subject: [PATCH 07/38] fix comments --- Tracking/components/TrackD0Printer.cpp | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/Tracking/components/TrackD0Printer.cpp b/Tracking/components/TrackD0Printer.cpp index 0c898f6f..4f80e7d0 100644 --- a/Tracking/components/TrackD0Printer.cpp +++ b/Tracking/components/TrackD0Printer.cpp @@ -14,11 +14,11 @@ #include #include -// Type alias for TrackCollection to improve readability +// Type aliases for improved readability using TrackColl = edm4hep::TrackCollection; using TP = edm4hep::TrackParams; -// Consumer that processes track collections and prints D0 values +// Consumer that processes track collections and prints phi values struct TrackD0Printer final : k4FWCore::Consumer { // Constructor: define the input collections (ClupatraTracks and SiTracks) TrackD0Printer(const std::string& name, ISvcLocator* svcLoc) From 528336bacdd3184947749c954c0e90e03bf46a5d Mon Sep 17 00:00:00 2001 From: Victor Schwan Date: Wed, 25 Jun 2025 14:31:08 +0200 Subject: [PATCH 08/38] python param for nStars in print stars in box --- Tracking/components/TrackD0Printer.cpp | 10 +++++----- Tracking/options/printD0.py | 2 +- 2 files changed, 6 insertions(+), 6 deletions(-) diff --git a/Tracking/components/TrackD0Printer.cpp b/Tracking/components/TrackD0Printer.cpp index 4f80e7d0..578b889c 100644 --- a/Tracking/components/TrackD0Printer.cpp +++ b/Tracking/components/TrackD0Printer.cpp @@ -31,7 +31,7 @@ struct TrackD0Printer final : k4FWCore::Consumer m_trackStateIndex{this, "TrackStateIndex", 2, "Index of track state to print (default 2)"}; + Gaudi::Property n_stars{this, "nStars", 20, "line-width of message in star box"}; void printInStars(const std::string& msg, const int lineWidth) const { info() << fmt::format("{:*^{}}", "", lineWidth) << endmsg; @@ -81,8 +80,9 @@ struct TrackD0Printer final : k4FWCore::Consumer 2) { - const edm4hep::TrackState& state = trackStates[2]; // Get the third track state + const size_t trackStateIndex = 2; + if (trackStates.size() > trackStateIndex) { + const edm4hep::TrackState& state = trackStates[trackStateIndex]; // Get the third track state info() << printValueUnc(trackType, varName, maxVarWidth, state.phi) << endmsg; info() << printValueUnc(trackType, sigmaVarName, maxVarWidth, getSigmaPhi(state)) << endmsg; } else { diff --git a/Tracking/options/printD0.py b/Tracking/options/printD0.py index aff33370..da8d4e63 100644 --- a/Tracking/options/printD0.py +++ b/Tracking/options/printD0.py @@ -12,7 +12,7 @@ iosvc.CollectionNames = ["SiTracks", "ClupatraTracks"] -printer = TrackD0Printer("TrackD0Printer") +printer = TrackD0Printer("TrackD0Printer", nStars=40) ApplicationMgr( TopAlg=[printer], From 9968891f761265114de6a0202a9825bad5c0141f Mon Sep 17 00:00:00 2001 From: Victor Schwan Date: Wed, 25 Jun 2025 16:18:10 +0200 Subject: [PATCH 09/38] Switch to track state AtIP instead of 2 --- Tracking/components/TrackD0Printer.cpp | 16 +++++++--------- 1 file changed, 7 insertions(+), 9 deletions(-) diff --git a/Tracking/components/TrackD0Printer.cpp b/Tracking/components/TrackD0Printer.cpp index 578b889c..abc87dbc 100644 --- a/Tracking/components/TrackD0Printer.cpp +++ b/Tracking/components/TrackD0Printer.cpp @@ -17,6 +17,7 @@ // Type aliases for improved readability using TrackColl = edm4hep::TrackCollection; using TP = edm4hep::TrackParams; +using TS = edm4hep::TrackState; // Consumer that processes track collections and prints phi values struct TrackD0Printer final : k4FWCore::Consumer { @@ -73,20 +74,17 @@ struct TrackD0Printer final : k4FWCore::Consumer trackStates = track.getTrackStates(); - std::string varName = "phi"; std::string sigmaVarName = "sigma " + varName; int maxVarWidth = std::max(varName.size(), sigmaVarName.size()) + 2; - // Check if there are enough track states (e.g., third track state) - const size_t trackStateIndex = 2; - if (trackStates.size() > trackStateIndex) { - const edm4hep::TrackState& state = trackStates[trackStateIndex]; // Get the third track state - info() << printValueUnc(trackType, varName, maxVarWidth, state.phi) << endmsg; - info() << printValueUnc(trackType, sigmaVarName, maxVarWidth, getSigmaPhi(state)) << endmsg; + // assuming there is only one Track State at IP + if (auto trackAtIP = std::ranges::find(track.getTrackStates(), TS::AtIP, &TS::location); + trackAtIP != track.getTrackStates().end()) { + info() << printValueUnc(trackType, varName, maxVarWidth, trackAtIP->phi) << endmsg; + info() << printValueUnc(trackType, sigmaVarName, maxVarWidth, getSigmaPhi(*trackAtIP)) << endmsg; } else { - warning() << trackType << " has less than 3 track states, skipping " + varName + " print" << endmsg; + fatal() << fmt::format("No track at IP found for {}!", trackType) << endmsg; } } }; From 86f5b565920757022c4ea5865cd24f7f8bd89e6a Mon Sep 17 00:00:00 2001 From: Victor Schwan Date: Wed, 25 Jun 2025 17:58:04 +0200 Subject: [PATCH 10/38] Switch to DEBUG output lvl for TrackD0Printer --- Tracking/options/printD0.py | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/Tracking/options/printD0.py b/Tracking/options/printD0.py index da8d4e63..4033aa55 100644 --- a/Tracking/options/printD0.py +++ b/Tracking/options/printD0.py @@ -1,7 +1,7 @@ from pathlib import Path from Configurables import TrackD0Printer -from Gaudi.Configuration import INFO +from Gaudi.Configuration import INFO, DEBUG from k4FWCore import ApplicationMgr, IOSvc iosvc = IOSvc() @@ -13,6 +13,7 @@ iosvc.CollectionNames = ["SiTracks", "ClupatraTracks"] printer = TrackD0Printer("TrackD0Printer", nStars=40) +printer.OutputLevel = DEBUG ApplicationMgr( TopAlg=[printer], From 296feecb3f009d5d86fe4d44b205ad7aee224887 Mon Sep 17 00:00:00 2001 From: Victor Schwan Date: Wed, 25 Jun 2025 18:47:44 +0200 Subject: [PATCH 11/38] Outsource printStars --- Tracking/components/TrackD0Printer.cpp | 9 ++------- Tracking/components/printStars.cpp | 10 ++++++++++ Tracking/components/printStars.h | 4 ++++ 3 files changed, 16 insertions(+), 7 deletions(-) create mode 100644 Tracking/components/printStars.cpp create mode 100644 Tracking/components/printStars.h diff --git a/Tracking/components/TrackD0Printer.cpp b/Tracking/components/TrackD0Printer.cpp index abc87dbc..257afd03 100644 --- a/Tracking/components/TrackD0Printer.cpp +++ b/Tracking/components/TrackD0Printer.cpp @@ -5,6 +5,7 @@ #include "edm4hep/TrackState.h" #include "k4FWCore/Consumer.h" #include "podio/RelationRange.h" +#include "printStars.h" #include #include @@ -32,7 +33,7 @@ struct TrackD0Printer final : k4FWCore::Consumer n_stars{this, "nStars", 20, "line-width of message in star box"}; - void printInStars(const std::string& msg, const int lineWidth) const { - info() << fmt::format("{:*^{}}", "", lineWidth) << endmsg; - info() << fmt::format("{:*^{}}", msg, lineWidth) << endmsg; - info() << fmt::format("{:*^{}}", "", lineWidth) << endmsg; - } - float getSigmaPhi(const edm4hep::TrackState& ts) const { return std::sqrt(ts.getCovMatrix(TP::phi, TP::phi)); } std::string printValueUnc(const std::string& strTrType, const std::string& vName, const int maxLabelWidth, diff --git a/Tracking/components/printStars.cpp b/Tracking/components/printStars.cpp new file mode 100644 index 00000000..4aa37376 --- /dev/null +++ b/Tracking/components/printStars.cpp @@ -0,0 +1,10 @@ +#include "printStars.h" +#include +#include +#include + +void printInStars(const Gaudi::Algorithm* thisAlg, const std::string& msg, const int lineWidth) { + thisAlg->debug() << fmt::format("{:*^{}}", "", lineWidth) << endmsg; + thisAlg->debug() << fmt::format("{:*^{}}", msg, lineWidth) << endmsg; + thisAlg->debug() << fmt::format("{:*^{}}", "", lineWidth) << endmsg; +} diff --git a/Tracking/components/printStars.h b/Tracking/components/printStars.h new file mode 100644 index 00000000..fc621318 --- /dev/null +++ b/Tracking/components/printStars.h @@ -0,0 +1,4 @@ +#include +#include + +void printInStars(const Gaudi::Algorithm* thisAlg, const std::string& msg, const int lineWidth); From 6e3e502572c0c581a4d9a78a3e9dfed28beec256 Mon Sep 17 00:00:00 2001 From: Victor Schwan Date: Fri, 27 Jun 2025 10:30:51 +0200 Subject: [PATCH 12/38] experimental --- Tracking/components/ExtractTrackParams.cpp | 99 ++++++++++++++++++++++ Tracking/options/ExtractTrackParams.py | 33 ++++++++ 2 files changed, 132 insertions(+) create mode 100644 Tracking/components/ExtractTrackParams.cpp create mode 100644 Tracking/options/ExtractTrackParams.py diff --git a/Tracking/components/ExtractTrackParams.cpp b/Tracking/components/ExtractTrackParams.cpp new file mode 100644 index 00000000..7ffbfd19 --- /dev/null +++ b/Tracking/components/ExtractTrackParams.cpp @@ -0,0 +1,99 @@ +#include "Gaudi/Property.h" +#include "edm4hep/TrackCollection.h" +#include "k4FWCore/Transformer.h" +#include "podio/UserDataCollection.h" +#include "printStars.h" + +#include +#include + +// Which type of collection we are reading +using FloatColl = podio::UserDataCollection; +using TrackColl = edm4hep::TrackCollection; +using TS = edm4hep::TrackState; + +struct TrackParamExtractor final + : k4FWCore::MultiTransformer(const TrackColl&, const TrackColl&)> { + TrackParamExtractor(const std::string& name, ISvcLocator* svcLoc) + : MultiTransformer(name, svcLoc, + { + KeyValues("InputSiTracks", {"SiTracks"}), + KeyValues("InputCluTracks", {"ClupatraTracks"}), + }, + { + KeyValues("OutCollSiPhi", {"SiTrackPhi"}), + KeyValues("OutCollSiOmega", {"SiTrackOmega"}), + KeyValues("OutCollSiD0", {"SiTrackD0"}), + KeyValues("OutCollSiTanL", {"SiTrackTanL"}), + KeyValues("OutCollSiZ0", {"SiTrackZ0"}), + KeyValues("OutCollCluPhi", {"CluTrackPhi"}), + KeyValues("OutCollCluOmega", {"CluTrackOmega"}), + KeyValues("OutCollCluD0", {"CluTrackD0"}), + KeyValues("OutCollCluTanL", {"CluTrackTanL"}), + KeyValues("OutCollCluZ0", {"CluTrackZ0"}), + + }) {} + + // This is the function that will be called to transform the data + // Note that the function has to be const, as well as the collections + // we get from the input + std::tuple + operator()(const TrackColl& inSiTracks, const TrackColl& inCluTracks) const override { + + printInStars(this, "New Event", n_stars); + + debug() << "Received SiTracks collection with " << inSiTracks.size() << " tracks" << endmsg; + debug() << "Received ClupatraTracks collection with " << inCluTracks.size() << " tracks" << endmsg; + + // Check that both collections are of the same size + if (inSiTracks.size() != inCluTracks.size()) { + fatal() << "Track collections have different sizes: SiTracks (" << inSiTracks.size() << ") and ClupatraTracks (" + << inCluTracks.size() << ")." << endmsg; + } + + // Process both collections (SiTracks and ClupatraTracks) at the same time + std::tuple siColls; + std::tuple cluColls; + for (size_t i = 0; i < inSiTracks.size(); ++i) { + + // Process SiTrack + siColls = processTrack(inSiTracks[i], "SiTrack"); + + // Process CluTrack + cluColls = processTrack(inCluTracks[i], "CluTrack"); + } + + return std::tuple_cat(std::move(siColls), std::move(cluColls)); + }; + +private: + Gaudi::Property n_stars{this, "nStars", 20, "line-width of message in star box"}; + + std::tuple processTrack(const edm4hep::Track& track, + const std::string& trackType) const { + + // assuming there is only one Track State at IP + auto trackAtIP = std::ranges::find(track.getTrackStates(), TS::AtIP, &TS::location); + if (trackAtIP != track.getTrackStates().end()) { + verbose() << fmt::format("Track at IP found for {}.", trackType) << endmsg; + } else { + fatal() << fmt::format("No track at IP found for {}!", trackType) << endmsg; + } + FloatColl trackPhi; + FloatColl trackOmega; + FloatColl trackD0; + FloatColl trackTanL; + FloatColl trackZ0; + trackPhi.push_back(trackAtIP->phi); + trackOmega.push_back(trackAtIP->omega); + trackD0.push_back(trackAtIP->D0); + trackTanL.push_back(trackAtIP->tanLambda); + trackZ0.push_back(trackAtIP->Z0); + + return std::make_tuple(std::move(trackPhi), std::move(trackOmega), std::move(trackD0), std::move(trackTanL), + std::move(trackZ0)); + } +}; +DECLARE_COMPONENT(TrackParamExtractor) \ No newline at end of file diff --git a/Tracking/options/ExtractTrackParams.py b/Tracking/options/ExtractTrackParams.py new file mode 100644 index 00000000..1fafba8a --- /dev/null +++ b/Tracking/options/ExtractTrackParams.py @@ -0,0 +1,33 @@ +from os import getenv +from pathlib import Path + +from Configurables import TrackParamExtractor +from Gaudi.Configuration import INFO, VERBOSE +from k4FWCore import ApplicationMgr, IOSvc + +fileSuffix = ".edm4hep.root" +versionName = "v0" +procName = "TrackParamExtractor" +basePath = Path(getenv("prmDir", Path.home() / "promotion")) + +iosvc = IOSvc() +iosvc.Input = str( + ( + basePath + / "code/ILDConfig/StandardConfig/production/data/test_tracking_3_detmods_V02_REC" + ).with_suffix(fileSuffix) +) +iosvc.Output = str((basePath / "data" / procName / versionName).with_suffix(fileSuffix)) +iosvc.CollectionNames = ["SiTracks", "ClupatraTracks"] +# iosvc.outputCommands = ["drop *", "keep SiTrackPhi"] + +printer = TrackParamExtractor(procName, nStars=40) +printer.OutputLevel = VERBOSE + +ApplicationMgr( + TopAlg=[printer], + EvtSel="NONE", + EvtMax=10, + ExtSvc=[iosvc], + OutputLevel=INFO, +) From 83bc354a79b430b3bb0628e78e848d528ba39fea Mon Sep 17 00:00:00 2001 From: Victor Schwan Date: Fri, 27 Jun 2025 21:41:16 +0200 Subject: [PATCH 13/38] use argparsing --- Tracking/options/ExtractTrackParams.py | 31 ++++++++++++++++++++++---- 1 file changed, 27 insertions(+), 4 deletions(-) diff --git a/Tracking/options/ExtractTrackParams.py b/Tracking/options/ExtractTrackParams.py index 1fafba8a..ca637868 100644 --- a/Tracking/options/ExtractTrackParams.py +++ b/Tracking/options/ExtractTrackParams.py @@ -4,20 +4,43 @@ from Configurables import TrackParamExtractor from Gaudi.Configuration import INFO, VERBOSE from k4FWCore import ApplicationMgr, IOSvc +from k4FWCore.parseArgs import parser + +detModOpts = ["v02", "if1", "if2"] +class CaseInsensitiveDict(dict): + def __getitem__(self, key): + return super().__getitem__(key.lower()) + + def __setitem__(self, key, value): + super().__setitem__(key.lower(), value) +detModNames = CaseInsensitiveDict({el: el for el in detModOpts}) +parser.add_argument( + "--detectorModel", + "-m", + help="Which detector model to run reconstruction for", + choices=detModOpts + [el.upper() for el in detModOpts], + type=str, + default="V02", +) +parser.add_argument( + "--version", type=str, help="str to identify a run through the pipeline" +) +args = parser.parse_known_args()[0] + + fileSuffix = ".edm4hep.root" -versionName = "v0" procName = "TrackParamExtractor" basePath = Path(getenv("prmDir", Path.home() / "promotion")) +corePath = f"{args.version}_{detModNames[args.detectorModel]}" iosvc = IOSvc() iosvc.Input = str( ( - basePath - / "code/ILDConfig/StandardConfig/production/data/test_tracking_3_detmods_V02_REC" + basePath / "code/ILDConfig/StandardConfig/production/data" / f"{corePath}_REC" ).with_suffix(fileSuffix) ) -iosvc.Output = str((basePath / "data" / procName / versionName).with_suffix(fileSuffix)) +iosvc.Output = str((basePath / "data" / procName / corePath).with_suffix(fileSuffix)) iosvc.CollectionNames = ["SiTracks", "ClupatraTracks"] # iosvc.outputCommands = ["drop *", "keep SiTrackPhi"] From 5c61ba3057106e695586b99027af8119ff965f4f Mon Sep 17 00:00:00 2001 From: Victor Schwan Date: Fri, 27 Jun 2025 21:43:58 +0200 Subject: [PATCH 14/38] add script to hist diff in track params --- Tracking/options/histTrackParams.py | 84 +++++++++++++++++++++++++++++ 1 file changed, 84 insertions(+) create mode 100644 Tracking/options/histTrackParams.py diff --git a/Tracking/options/histTrackParams.py b/Tracking/options/histTrackParams.py new file mode 100644 index 00000000..36d5a318 --- /dev/null +++ b/Tracking/options/histTrackParams.py @@ -0,0 +1,84 @@ +from itertools import product +from os import getenv +from pathlib import Path +from argparse import ArgumentParser + +import matplotlib as mpl +import matplotlib.pyplot as plt +import uproot + +# general plotting options +labelsize = 20 +linewidth = 1.5 +majorTickSize = 10 +plot_grid_alpha = .7 +params = { + "xtick.direction": "in", + "ytick.direction": "in", + "xtick.top": True, + "ytick.right": True, + "xtick.major.size": majorTickSize, # Tick length + "ytick.major.size": majorTickSize, + "xtick.major.width": linewidth, # Tick line width + "ytick.major.width": linewidth, + "axes.linewidth": linewidth, + "legend.fontsize": labelsize, # "x-large", + "axes.labelsize": labelsize, # "x-large", + "axes.titlesize": labelsize, # "x-large", + "xtick.labelsize": labelsize, # "x-large", + "ytick.labelsize": labelsize, # "x-large", + "figure.autolayout": True, +} #'figure.figsize': (15, 5), +mpl.rcParams.update(params) + +# argparse +detModOpts = ["v02", "if1", "if2"] +class CaseInsensitiveDict(dict): + def __getitem__(self, key): + return super().__getitem__(key.lower()) + + def __setitem__(self, key, value): + super().__setitem__(key.lower(), value) +detModNames = CaseInsensitiveDict({el: el for el in detModOpts}) +parser = ArgumentParser() +parser.add_argument( + "--detectorModel", + "-m", + help="Which detector model to run reconstruction for", + choices=detModOpts + [el.upper() for el in detModOpts], + type=str, + default="V02", +) +parser.add_argument( + "--version", type=str, help="str to identify a run through the pipeline" +) +args = parser.parse_known_args()[0] +corePath = Path(f"{args.version}_{detModNames[args.detectorModel]}") + +# strings to build path +processor = "TrackParamExtractor" +basePath = Path(getenv("dtDir", str(Path.home() / "promotion" / "data"))) + +# Lists to build branch names to be analyzed +trackNames = ["SiTrack", "CluTrack"] +varNames = ["D0", "Omega", "Phi", "TanL", "Z0"] + +# build vars based on above vars +keys = [f"{trackName}{varName}" for trackName, varName in product(trackNames, varNames)] +in_file = basePath / processor / corePath.with_suffix(".edm4hep.root") + +print(in_file) +with uproot.open(str(in_file) + ":events") as events: + regex = f"/^({'|'.join(trackNames)})({'|'.join(varNames)})$/" + data = events.arrays(filter_name=regex, library="pd") + for var in varNames: + data[f"d_{var}"] = data[f"{trackNames[0]}{var}"] - data[f"{trackNames[1]}{var}"] + +for varName in varNames: + plt.figure() + plt.grid(True, which="both", linestyle="--", linewidth=linewidth, alpha=plot_grid_alpha) + plt.hist(data[f"d_{varName}"], bins=30) + plt.xlabel(rf"$\Delta$ {varName}") + plt.ylabel("Frequency") + plt.title(args.detectorModel) + plt.show() From 6d8441131cda303f204e9180ea9a05c6a09ba1c4 Mon Sep 17 00:00:00 2001 From: Victor Schwan Date: Fri, 27 Jun 2025 22:29:45 +0200 Subject: [PATCH 15/38] plot multiple vars in shared figure --- Tracking/options/histTrackParams.py | 38 ++++++++++++++++++++++------- 1 file changed, 29 insertions(+), 9 deletions(-) diff --git a/Tracking/options/histTrackParams.py b/Tracking/options/histTrackParams.py index 36d5a318..5811466a 100644 --- a/Tracking/options/histTrackParams.py +++ b/Tracking/options/histTrackParams.py @@ -1,17 +1,18 @@ +from argparse import ArgumentParser from itertools import product from os import getenv from pathlib import Path -from argparse import ArgumentParser +import awkward as ak import matplotlib as mpl import matplotlib.pyplot as plt import uproot # general plotting options -labelsize = 20 +labelsize = 22 linewidth = 1.5 majorTickSize = 10 -plot_grid_alpha = .7 +plot_grid_alpha = 0.7 params = { "xtick.direction": "in", "ytick.direction": "in", @@ -61,7 +62,9 @@ def __setitem__(self, key, value): # Lists to build branch names to be analyzed trackNames = ["SiTrack", "CluTrack"] -varNames = ["D0", "Omega", "Phi", "TanL", "Z0"] +varSimilar = ["Phi", "Omega", "TanL"] +varSpread = ["D0", "Z0"] +varNames = varSimilar + varSpread # build vars based on above vars keys = [f"{trackName}{varName}" for trackName, varName in product(trackNames, varNames)] @@ -74,11 +77,28 @@ def __setitem__(self, key, value): for var in varNames: data[f"d_{var}"] = data[f"{trackNames[0]}{var}"] - data[f"{trackNames[1]}{var}"] -for varName in varNames: +# for varName in varNames: +# plt.figure() +# plt.grid(True, which="both", linestyle="--", linewidth=linewidth, alpha=plot_grid_alpha) +# plt.hist(data[f"d_{varName}"], bins=30) +# plt.xlabel(rf"$\Delta$ {varName}") +# plt.ylabel("Frequency") +# plt.title(args.detectorModel) +# plt.show() + +for group, xlim in zip([varSpread, varSimilar], [1.5, 0.0015]): plt.figure() - plt.grid(True, which="both", linestyle="--", linewidth=linewidth, alpha=plot_grid_alpha) - plt.hist(data[f"d_{varName}"], bins=30) - plt.xlabel(rf"$\Delta$ {varName}") + plt.grid( + True, which="both", linestyle="--", linewidth=linewidth, alpha=plot_grid_alpha + ) + plt.hist( + x=[ak.to_numpy(ak.flatten(data[f"d_{varName}"])) for varName in group], + bins=30, + label=group, + range=(-xlim, xlim) + ) + plt.xlabel(rf"$\Delta$") plt.ylabel("Frequency") - plt.title(args.detectorModel) + plt.title(f"{args.detectorModel}: {','.join(group)}") + plt.legend() plt.show() From c8a05f4f53f8c6e8c66875ad96830e8a70cc5c21 Mon Sep 17 00:00:00 2001 From: Victor Schwan Date: Tue, 1 Jul 2025 15:32:50 +0200 Subject: [PATCH 16/38] Make TrackState nullable w/ std::optional, use TrackParam enum to assign params to tuple indices --- Tracking/components/ExtractTrackParams.cpp | 48 ++++++++++++---------- 1 file changed, 27 insertions(+), 21 deletions(-) diff --git a/Tracking/components/ExtractTrackParams.cpp b/Tracking/components/ExtractTrackParams.cpp index 7ffbfd19..4fbb827a 100644 --- a/Tracking/components/ExtractTrackParams.cpp +++ b/Tracking/components/ExtractTrackParams.cpp @@ -3,7 +3,9 @@ #include "k4FWCore/Transformer.h" #include "podio/UserDataCollection.h" #include "printStars.h" +#include +#include #include #include @@ -22,16 +24,16 @@ struct TrackParamExtractor final KeyValues("InputCluTracks", {"ClupatraTracks"}), }, { + KeyValues("OutCollSiD0", {"SiTrackD0"}), KeyValues("OutCollSiPhi", {"SiTrackPhi"}), KeyValues("OutCollSiOmega", {"SiTrackOmega"}), - KeyValues("OutCollSiD0", {"SiTrackD0"}), - KeyValues("OutCollSiTanL", {"SiTrackTanL"}), KeyValues("OutCollSiZ0", {"SiTrackZ0"}), + KeyValues("OutCollSiTanL", {"SiTrackTanL"}), + KeyValues("OutCollCluD0", {"CluTrackD0"}), KeyValues("OutCollCluPhi", {"CluTrackPhi"}), KeyValues("OutCollCluOmega", {"CluTrackOmega"}), - KeyValues("OutCollCluD0", {"CluTrackD0"}), - KeyValues("OutCollCluTanL", {"CluTrackTanL"}), KeyValues("OutCollCluZ0", {"CluTrackZ0"}), + KeyValues("OutCollCluTanL", {"CluTrackTanL"}), }) {} @@ -59,10 +61,26 @@ struct TrackParamExtractor final for (size_t i = 0; i < inSiTracks.size(); ++i) { // Process SiTrack - siColls = processTrack(inSiTracks[i], "SiTrack"); + const auto oSiTrackStateIP = getOTrackAtIP(inSiTracks[i], "SiTrack"); + if (oSiTrackStateIP.has_value()) { + std::get(edm4hep::TrackParams::d0)>(siColls).push_back(oSiTrackStateIP->D0); + std::get(edm4hep::TrackParams::phi)>(siColls).push_back(oSiTrackStateIP->phi); + std::get(edm4hep::TrackParams::omega)>(siColls).push_back(oSiTrackStateIP->omega); + std::get(edm4hep::TrackParams::z0)>(siColls).push_back(oSiTrackStateIP->Z0); + std::get(edm4hep::TrackParams::tanLambda)>(siColls).push_back( + oSiTrackStateIP->tanLambda); + } // Process CluTrack - cluColls = processTrack(inCluTracks[i], "CluTrack"); + const auto oCluTrackStateIP = getOTrackAtIP(inCluTracks[i], "CluTrack"); + if (oCluTrackStateIP.has_value()) { + std::get(edm4hep::TrackParams::d0)>(cluColls).push_back(oCluTrackStateIP->D0); + std::get(edm4hep::TrackParams::phi)>(cluColls).push_back(oCluTrackStateIP->phi); + std::get(edm4hep::TrackParams::omega)>(cluColls).push_back(oCluTrackStateIP->omega); + std::get(edm4hep::TrackParams::z0)>(cluColls).push_back(oCluTrackStateIP->Z0); + std::get(edm4hep::TrackParams::tanLambda)>(cluColls).push_back( + oCluTrackStateIP->tanLambda); + } } return std::tuple_cat(std::move(siColls), std::move(cluColls)); @@ -71,29 +89,17 @@ struct TrackParamExtractor final private: Gaudi::Property n_stars{this, "nStars", 20, "line-width of message in star box"}; - std::tuple processTrack(const edm4hep::Track& track, - const std::string& trackType) const { + std::optional getOTrackAtIP(const edm4hep::Track& track, const std::string& trackType) const { // assuming there is only one Track State at IP auto trackAtIP = std::ranges::find(track.getTrackStates(), TS::AtIP, &TS::location); if (trackAtIP != track.getTrackStates().end()) { verbose() << fmt::format("Track at IP found for {}.", trackType) << endmsg; + return *trackAtIP; } else { fatal() << fmt::format("No track at IP found for {}!", trackType) << endmsg; + return std::nullopt; } - FloatColl trackPhi; - FloatColl trackOmega; - FloatColl trackD0; - FloatColl trackTanL; - FloatColl trackZ0; - trackPhi.push_back(trackAtIP->phi); - trackOmega.push_back(trackAtIP->omega); - trackD0.push_back(trackAtIP->D0); - trackTanL.push_back(trackAtIP->tanLambda); - trackZ0.push_back(trackAtIP->Z0); - - return std::make_tuple(std::move(trackPhi), std::move(trackOmega), std::move(trackD0), std::move(trackTanL), - std::move(trackZ0)); } }; DECLARE_COMPONENT(TrackParamExtractor) \ No newline at end of file From d933db607b7d31bc3ac6a22a7612c5e4344cd13c Mon Sep 17 00:00:00 2001 From: Victor Schwan Date: Tue, 1 Jul 2025 15:34:31 +0200 Subject: [PATCH 17/38] make ExtractTrackParams compatible with FCC models --- Tracking/options/ExtractTrackParams.py | 16 ++++++++++++++-- 1 file changed, 14 insertions(+), 2 deletions(-) diff --git a/Tracking/options/ExtractTrackParams.py b/Tracking/options/ExtractTrackParams.py index ca637868..10ad2071 100644 --- a/Tracking/options/ExtractTrackParams.py +++ b/Tracking/options/ExtractTrackParams.py @@ -7,13 +7,18 @@ from k4FWCore.parseArgs import parser detModOpts = ["v02", "if1", "if2"] + + class CaseInsensitiveDict(dict): def __getitem__(self, key): return super().__getitem__(key.lower()) def __setitem__(self, key, value): super().__setitem__(key.lower(), value) + + detModNames = CaseInsensitiveDict({el: el for el in detModOpts}) + parser.add_argument( "--detectorModel", "-m", @@ -28,7 +33,6 @@ def __setitem__(self, key, value): args = parser.parse_known_args()[0] - fileSuffix = ".edm4hep.root" procName = "TrackParamExtractor" basePath = Path(getenv("prmDir", Path.home() / "promotion")) @@ -41,12 +45,20 @@ def __setitem__(self, key, value): ).with_suffix(fileSuffix) ) iosvc.Output = str((basePath / "data" / procName / corePath).with_suffix(fileSuffix)) -iosvc.CollectionNames = ["SiTracks", "ClupatraTracks"] # iosvc.outputCommands = ["drop *", "keep SiTrackPhi"] printer = TrackParamExtractor(procName, nStars=40) printer.OutputLevel = VERBOSE +# the collection name with the SiTracks differs between ILC and FCC models +if "IF" in args.detectorModel: + printer.InputSiTracks = ["SiTracksCT"] + siTrackCollName = "SiTracksCT" +else: + siTrackCollName = "SiTracks" +iosvc.CollectionNames = ["ClupatraTracks", siTrackCollName] + + ApplicationMgr( TopAlg=[printer], EvtSel="NONE", From 6fa7d895ae7dd82278ab272b113e5238cb4078f2 Mon Sep 17 00:00:00 2001 From: Victor Schwan Date: Wed, 2 Jul 2025 20:15:23 +0200 Subject: [PATCH 18/38] outsource commonArgParsing, hists for several detMods --- Tracking/options/ExtractTrackParams.py | 51 ++++------- Tracking/options/histTrackParams.py | 115 ++++++++++++++----------- 2 files changed, 82 insertions(+), 84 deletions(-) diff --git a/Tracking/options/ExtractTrackParams.py b/Tracking/options/ExtractTrackParams.py index 10ad2071..0fd88ae6 100644 --- a/Tracking/options/ExtractTrackParams.py +++ b/Tracking/options/ExtractTrackParams.py @@ -1,4 +1,5 @@ -from os import getenv +import os +import sys from pathlib import Path from Configurables import TrackParamExtractor @@ -6,44 +7,28 @@ from k4FWCore import ApplicationMgr, IOSvc from k4FWCore.parseArgs import parser -detModOpts = ["v02", "if1", "if2"] +sys.path.append(os.getenv("trckOptDir")) +from commonArgParsing import addCommonArgs, detModNames - -class CaseInsensitiveDict(dict): - def __getitem__(self, key): - return super().__getitem__(key.lower()) - - def __setitem__(self, key, value): - super().__setitem__(key.lower(), value) - - -detModNames = CaseInsensitiveDict({el: el for el in detModOpts}) - -parser.add_argument( - "--detectorModel", - "-m", - help="Which detector model to run reconstruction for", - choices=detModOpts + [el.upper() for el in detModOpts], - type=str, - default="V02", +args = addCommonArgs(parser).parse_known_args()[0] +assert len(args.detectorModels) == 1, ( + f"Only provide one detector model! You provided {args.detectorModels}" ) -parser.add_argument( - "--version", type=str, help="str to identify a run through the pipeline" -) -args = parser.parse_known_args()[0] - +args.detectorModels = args.detectorModels[0] fileSuffix = ".edm4hep.root" procName = "TrackParamExtractor" -basePath = Path(getenv("prmDir", Path.home() / "promotion")) -corePath = f"{args.version}_{detModNames[args.detectorModel]}" +basePath = Path(os.getenv("prmDir", Path.home() / "promotion")) +corePath = f"{args.version}_{detModNames[args.detectorModels]}" + +# assert that the input path exists +inputPath = ( + basePath / "code/ILDConfig/StandardConfig/production/data" / f"{corePath}_REC" +).with_suffix(fileSuffix) +assert inputPath.exists(), f"ERROR: The input path ({inputPath}) does not exist!" iosvc = IOSvc() -iosvc.Input = str( - ( - basePath / "code/ILDConfig/StandardConfig/production/data" / f"{corePath}_REC" - ).with_suffix(fileSuffix) -) +iosvc.Input = str(inputPath) iosvc.Output = str((basePath / "data" / procName / corePath).with_suffix(fileSuffix)) # iosvc.outputCommands = ["drop *", "keep SiTrackPhi"] @@ -51,7 +36,7 @@ def __setitem__(self, key, value): printer.OutputLevel = VERBOSE # the collection name with the SiTracks differs between ILC and FCC models -if "IF" in args.detectorModel: +if "IF" in args.detectorModels: printer.InputSiTracks = ["SiTracksCT"] siTrackCollName = "SiTracksCT" else: diff --git a/Tracking/options/histTrackParams.py b/Tracking/options/histTrackParams.py index 5811466a..521d2802 100644 --- a/Tracking/options/histTrackParams.py +++ b/Tracking/options/histTrackParams.py @@ -8,11 +8,13 @@ import matplotlib.pyplot as plt import uproot +from commonArgParsing import addCommonArgs, detModNames + # general plotting options labelsize = 22 linewidth = 1.5 majorTickSize = 10 -plot_grid_alpha = 0.7 +plotGridAlpha = 0.7 params = { "xtick.direction": "in", "ytick.direction": "in", @@ -32,33 +34,15 @@ } #'figure.figsize': (15, 5), mpl.rcParams.update(params) -# argparse -detModOpts = ["v02", "if1", "if2"] -class CaseInsensitiveDict(dict): - def __getitem__(self, key): - return super().__getitem__(key.lower()) - - def __setitem__(self, key, value): - super().__setitem__(key.lower(), value) -detModNames = CaseInsensitiveDict({el: el for el in detModOpts}) -parser = ArgumentParser() -parser.add_argument( - "--detectorModel", - "-m", - help="Which detector model to run reconstruction for", - choices=detModOpts + [el.upper() for el in detModOpts], - type=str, - default="V02", -) -parser.add_argument( - "--version", type=str, help="str to identify a run through the pipeline" -) -args = parser.parse_known_args()[0] -corePath = Path(f"{args.version}_{detModNames[args.detectorModel]}") +args = addCommonArgs(ArgumentParser()).parse_known_args()[0] -# strings to build path -processor = "TrackParamExtractor" -basePath = Path(getenv("dtDir", str(Path.home() / "promotion" / "data"))) +################################# +# commands to access cov Matrix +################################# +# import ROOT +# ROOT.gInterpreter.LoadFile("edm4hep/utils/cov_matrix_utils.h") +# # ... +# edm4hep.utils.get_cov_value(cov_m, edm4hep.TrackParams.d0, edm4hep.TrackParams.d0) # Lists to build branch names to be analyzed trackNames = ["SiTrack", "CluTrack"] @@ -66,39 +50,68 @@ def __setitem__(self, key, value): varSpread = ["D0", "Z0"] varNames = varSimilar + varSpread -# build vars based on above vars -keys = [f"{trackName}{varName}" for trackName, varName in product(trackNames, varNames)] -in_file = basePath / processor / corePath.with_suffix(".edm4hep.root") +data = {} + +# extract data +for detMod in args.detectorModels: + corePath = Path(f"{args.version}_{detModNames[detMod]}") + + # strings to build path + processor = "TrackParamExtractor" + basePath = Path(getenv("dtDir", str(Path.home() / "promotion" / "data"))) + + # build vars based on above vars + keys = [ + f"{trackName}{varName}" for trackName, varName in product(trackNames, varNames) + ] + in_file = basePath / processor / corePath.with_suffix(".edm4hep.root") -print(in_file) -with uproot.open(str(in_file) + ":events") as events: - regex = f"/^({'|'.join(trackNames)})({'|'.join(varNames)})$/" - data = events.arrays(filter_name=regex, library="pd") - for var in varNames: - data[f"d_{var}"] = data[f"{trackNames[0]}{var}"] - data[f"{trackNames[1]}{var}"] + with uproot.open(str(in_file) + ":events") as events: + regex = f"/^({'|'.join(trackNames)})({'|'.join(varNames)})$/" + data[detMod] = events.arrays(filter_name=regex, library="pd") + for var in varNames: + data[detMod][f"d_{var}"] = ( + data[detMod][f"{trackNames[0]}{var}"] + - data[detMod][f"{trackNames[1]}{var}"] + ) -# for varName in varNames: -# plt.figure() -# plt.grid(True, which="both", linestyle="--", linewidth=linewidth, alpha=plot_grid_alpha) -# plt.hist(data[f"d_{varName}"], bins=30) -# plt.xlabel(rf"$\Delta$ {varName}") -# plt.ylabel("Frequency") -# plt.title(args.detectorModel) -# plt.show() +# # plot data +# for detMod in args.detectorModels: +# for group, xlim in zip([varSpread, varSimilar], [1.5, 0.0015]): +# plt.figure() +# plt.grid( +# True, which="both", linestyle="--", linewidth=linewidth, alpha=plotGridAlpha +# ) +# plt.hist( +# x=[ +# ak.to_numpy(ak.flatten(data[detMod][f"d_{varName}"])) +# for varName in group +# ], +# bins=30, +# label=group, +# range=(-xlim, xlim), +# ) +# plt.xlabel(rf"$\Delta$ Si-Clu") +# plt.ylabel("Frequency") +# plt.title(f"{detMod}: {','.join(group)}") +# plt.legend() +# plt.show() -for group, xlim in zip([varSpread, varSimilar], [1.5, 0.0015]): +for var in ["D0"]: plt.figure() plt.grid( - True, which="both", linestyle="--", linewidth=linewidth, alpha=plot_grid_alpha + True, which="both", linestyle="--", linewidth=linewidth, alpha=plotGridAlpha ) plt.hist( - x=[ak.to_numpy(ak.flatten(data[f"d_{varName}"])) for varName in group], + x=[ + ak.to_numpy(ak.flatten(data[detMod][f"{trackNames[0]}{var}"])) + for detMod in args.detectorModels + ], bins=30, - label=group, - range=(-xlim, xlim) + label=args.detectorModels, + range=(-0.03, 0.03), ) - plt.xlabel(rf"$\Delta$") plt.ylabel("Frequency") - plt.title(f"{args.detectorModel}: {','.join(group)}") + plt.title(f"{var}") plt.legend() plt.show() From 8d15e1bfdc122f6166b348edf736d7e2886d9fb2 Mon Sep 17 00:00:00 2001 From: Victor Schwan Date: Wed, 2 Jul 2025 20:16:41 +0200 Subject: [PATCH 19/38] printD0: reduce output level to INFO --- Tracking/options/printD0.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/Tracking/options/printD0.py b/Tracking/options/printD0.py index 4033aa55..e6c5dd3c 100644 --- a/Tracking/options/printD0.py +++ b/Tracking/options/printD0.py @@ -1,7 +1,7 @@ from pathlib import Path from Configurables import TrackD0Printer -from Gaudi.Configuration import INFO, DEBUG +from Gaudi.Configuration import INFO from k4FWCore import ApplicationMgr, IOSvc iosvc = IOSvc() @@ -13,7 +13,7 @@ iosvc.CollectionNames = ["SiTracks", "ClupatraTracks"] printer = TrackD0Printer("TrackD0Printer", nStars=40) -printer.OutputLevel = DEBUG +printer.OutputLevel = INFO ApplicationMgr( TopAlg=[printer], From ef8b0d596e74d2f560286764f121d437d76687b3 Mon Sep 17 00:00:00 2001 From: Victor Schwan Date: Wed, 2 Jul 2025 22:53:19 +0200 Subject: [PATCH 20/38] alignment of det mod handling with BS repo --- Tracking/options/ExtractTrackParams.py | 4 +- Tracking/options/histTrackParams.py | 51 +++++++++++++------------- 2 files changed, 28 insertions(+), 27 deletions(-) diff --git a/Tracking/options/ExtractTrackParams.py b/Tracking/options/ExtractTrackParams.py index 0fd88ae6..593e4b1b 100644 --- a/Tracking/options/ExtractTrackParams.py +++ b/Tracking/options/ExtractTrackParams.py @@ -8,9 +8,9 @@ from k4FWCore.parseArgs import parser sys.path.append(os.getenv("trckOptDir")) -from commonArgParsing import addCommonArgs, detModNames +from commonArgParsing import add_common_args, detModNames -args = addCommonArgs(parser).parse_known_args()[0] +args = add_common_args(parser).parse_known_args()[0] assert len(args.detectorModels) == 1, ( f"Only provide one detector model! You provided {args.detectorModels}" ) diff --git a/Tracking/options/histTrackParams.py b/Tracking/options/histTrackParams.py index 521d2802..7af41c8d 100644 --- a/Tracking/options/histTrackParams.py +++ b/Tracking/options/histTrackParams.py @@ -8,7 +8,7 @@ import matplotlib.pyplot as plt import uproot -from commonArgParsing import addCommonArgs, detModNames +from commonArgParsing import add_common_args, detModNames, registry # general plotting options labelsize = 22 @@ -34,7 +34,7 @@ } #'figure.figsize': (15, 5), mpl.rcParams.update(params) -args = addCommonArgs(ArgumentParser()).parse_known_args()[0] +args = add_common_args(ArgumentParser()).parse_known_args()[0] ################################# # commands to access cov Matrix @@ -45,7 +45,7 @@ # edm4hep.utils.get_cov_value(cov_m, edm4hep.TrackParams.d0, edm4hep.TrackParams.d0) # Lists to build branch names to be analyzed -trackNames = ["SiTrack", "CluTrack"] +trackType = ["SiTrack", "CluTrack"] varSimilar = ["Phi", "Omega", "TanL"] varSpread = ["D0", "Z0"] varNames = varSimilar + varSpread @@ -62,17 +62,17 @@ # build vars based on above vars keys = [ - f"{trackName}{varName}" for trackName, varName in product(trackNames, varNames) + f"{trackName}{varName}" for trackName, varName in product(trackType, varNames) ] in_file = basePath / processor / corePath.with_suffix(".edm4hep.root") with uproot.open(str(in_file) + ":events") as events: - regex = f"/^({'|'.join(trackNames)})({'|'.join(varNames)})$/" + regex = f"/^({'|'.join(trackType)})({'|'.join(varNames)})$/" data[detMod] = events.arrays(filter_name=regex, library="pd") for var in varNames: data[detMod][f"d_{var}"] = ( - data[detMod][f"{trackNames[0]}{var}"] - - data[detMod][f"{trackNames[1]}{var}"] + data[detMod][f"{trackType[0]}{var}"] + - data[detMod][f"{trackType[1]}{var}"] ) # # plot data @@ -97,21 +97,22 @@ # plt.legend() # plt.show() -for var in ["D0"]: - plt.figure() - plt.grid( - True, which="both", linestyle="--", linewidth=linewidth, alpha=plotGridAlpha - ) - plt.hist( - x=[ - ak.to_numpy(ak.flatten(data[detMod][f"{trackNames[0]}{var}"])) - for detMod in args.detectorModels - ], - bins=30, - label=args.detectorModels, - range=(-0.03, 0.03), - ) - plt.ylabel("Frequency") - plt.title(f"{var}") - plt.legend() - plt.show() +for type, xlim in zip(trackType, [.03,1]): + for var in ["D0"]: + plt.figure() + plt.grid( + True, which="both", linestyle="--", linewidth=linewidth, alpha=plotGridAlpha + ) + plt.hist( + x=[ + ak.to_numpy(ak.flatten(data[detMod][f"{type}{var}"])) + for detMod in args.detectorModels + ], + bins=30, + label=[registry.get(detMod).get_name(args) for detMod in args.detectorModels], + range=(-xlim, xlim), + ) + plt.ylabel("Frequency") + plt.title(f"{type}: {var}") + plt.legend() + plt.show() From 2524f735ebdc7490b6c798a2f8c0bbaf20447c42 Mon Sep 17 00:00:00 2001 From: Victor Schwan Date: Thu, 3 Jul 2025 14:46:29 +0200 Subject: [PATCH 21/38] Introduce DetMod property at_fcc/ilc, use enum to avoid to hard code accelerator names --- Tracking/options/ExtractTrackParams.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/Tracking/options/ExtractTrackParams.py b/Tracking/options/ExtractTrackParams.py index 593e4b1b..a1262ca1 100644 --- a/Tracking/options/ExtractTrackParams.py +++ b/Tracking/options/ExtractTrackParams.py @@ -8,7 +8,7 @@ from k4FWCore.parseArgs import parser sys.path.append(os.getenv("trckOptDir")) -from commonArgParsing import add_common_args, detModNames +from commonArgParsing import add_common_args, detModNames, registry args = add_common_args(parser).parse_known_args()[0] assert len(args.detectorModels) == 1, ( @@ -36,7 +36,7 @@ printer.OutputLevel = VERBOSE # the collection name with the SiTracks differs between ILC and FCC models -if "IF" in args.detectorModels: +if registry.get(args.detectorModels).at_fcc: printer.InputSiTracks = ["SiTracksCT"] siTrackCollName = "SiTracksCT" else: From 6aee4d61186ce7a6578a51814e927221674bf01e Mon Sep 17 00:00:00 2001 From: Victor Schwan Date: Thu, 3 Jul 2025 15:08:09 +0200 Subject: [PATCH 22/38] Outsource hist xlims in proper dict --- Tracking/options/histTrackParams.py | 13 +++++++++---- 1 file changed, 9 insertions(+), 4 deletions(-) diff --git a/Tracking/options/histTrackParams.py b/Tracking/options/histTrackParams.py index 7af41c8d..42a2f221 100644 --- a/Tracking/options/histTrackParams.py +++ b/Tracking/options/histTrackParams.py @@ -97,8 +97,11 @@ # plt.legend() # plt.show() -for type, xlim in zip(trackType, [.03,1]): - for var in ["D0"]: +xlims = {type: None for type in trackType} +xlims["SiTrack"] = {"D0": 0.03, "Omega": 0.0002} +xlims["CluTrack"] = {"D0": .8, "Omega": 0.00025} +for type in trackType: + for var in ["D0", "Omega"]: plt.figure() plt.grid( True, which="both", linestyle="--", linewidth=linewidth, alpha=plotGridAlpha @@ -109,8 +112,10 @@ for detMod in args.detectorModels ], bins=30, - label=[registry.get(detMod).get_name(args) for detMod in args.detectorModels], - range=(-xlim, xlim), + label=[ + registry.get(detMod).get_name(args) for detMod in args.detectorModels + ], + range=(-xlims[type][var], xlims[type][var]) if var in xlims[type] else None, ) plt.ylabel("Frequency") plt.title(f"{type}: {var}") From b76c92a9a3b01a09d522ae357c62d1731461c759 Mon Sep 17 00:00:00 2001 From: Victor Schwan Date: Thu, 3 Jul 2025 16:13:07 +0200 Subject: [PATCH 23/38] HistTrackParams: use typewriter font in the title --- Tracking/options/histTrackParams.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/Tracking/options/histTrackParams.py b/Tracking/options/histTrackParams.py index 42a2f221..884f863a 100644 --- a/Tracking/options/histTrackParams.py +++ b/Tracking/options/histTrackParams.py @@ -93,7 +93,7 @@ # ) # plt.xlabel(rf"$\Delta$ Si-Clu") # plt.ylabel("Frequency") -# plt.title(f"{detMod}: {','.join(group)}") +# plt.title(rf"Diff Si-Clu in $\mathtt{{{registry.get(detMod).get_name(args)}}}$: {','.join(group)}") # plt.legend() # plt.show() @@ -118,6 +118,6 @@ range=(-xlims[type][var], xlims[type][var]) if var in xlims[type] else None, ) plt.ylabel("Frequency") - plt.title(f"{type}: {var}") + plt.title(f"Diff DetMods $\mathtt{{{type}}}$: {var}") plt.legend() plt.show() From 13f3f73bbd04bb203ad743c33f08f44f3d3918ec Mon Sep 17 00:00:00 2001 From: Victor Schwan Date: Fri, 4 Jul 2025 13:12:14 +0200 Subject: [PATCH 24/38] histTrackParams: add argparse opts to choose which hists to show --- Tracking/options/histTrackParams.py | 98 +++++++++++++++-------------- 1 file changed, 52 insertions(+), 46 deletions(-) diff --git a/Tracking/options/histTrackParams.py b/Tracking/options/histTrackParams.py index 884f863a..ec6b2d14 100644 --- a/Tracking/options/histTrackParams.py +++ b/Tracking/options/histTrackParams.py @@ -34,7 +34,11 @@ } #'figure.figsize': (15, 5), mpl.rcParams.update(params) -args = add_common_args(ArgumentParser()).parse_known_args()[0] +parser = add_common_args(ArgumentParser()) +plot_opts = parser.add_argument_group("Plotting opts", "which plots should be shown") +plot_opts.add_argument("--track", action="store_true", help="Show difference between Silicon and Clupatra tracks") +plot_opts.add_argument("--detmods", action="store_true", help="Show difference between detector models") +args = parser.parse_known_args()[0] ################################# # commands to access cov Matrix @@ -75,49 +79,51 @@ - data[detMod][f"{trackType[1]}{var}"] ) -# # plot data -# for detMod in args.detectorModels: -# for group, xlim in zip([varSpread, varSimilar], [1.5, 0.0015]): -# plt.figure() -# plt.grid( -# True, which="both", linestyle="--", linewidth=linewidth, alpha=plotGridAlpha -# ) -# plt.hist( -# x=[ -# ak.to_numpy(ak.flatten(data[detMod][f"d_{varName}"])) -# for varName in group -# ], -# bins=30, -# label=group, -# range=(-xlim, xlim), -# ) -# plt.xlabel(rf"$\Delta$ Si-Clu") -# plt.ylabel("Frequency") -# plt.title(rf"Diff Si-Clu in $\mathtt{{{registry.get(detMod).get_name(args)}}}$: {','.join(group)}") -# plt.legend() -# plt.show() +if args.track: + # plot data + for detMod in args.detectorModels: + for group, xlim in zip([varSpread, varSimilar], [1.5, 0.0015]): + plt.figure() + plt.grid( + True, which="both", linestyle="--", linewidth=linewidth, alpha=plotGridAlpha + ) + plt.hist( + x=[ + ak.to_numpy(ak.flatten(data[detMod][f"d_{varName}"])) + for varName in group + ], + bins=30, + label=group, + range=(-xlim, xlim), + ) + plt.xlabel(rf"$\Delta$ Si-Clu") + plt.ylabel("Frequency") + plt.title(rf"Diff Si-Clu in $\mathtt{{{registry.get(detMod).get_name(args)}}}$: {','.join(group)}") + plt.legend() + plt.show() -xlims = {type: None for type in trackType} -xlims["SiTrack"] = {"D0": 0.03, "Omega": 0.0002} -xlims["CluTrack"] = {"D0": .8, "Omega": 0.00025} -for type in trackType: - for var in ["D0", "Omega"]: - plt.figure() - plt.grid( - True, which="both", linestyle="--", linewidth=linewidth, alpha=plotGridAlpha - ) - plt.hist( - x=[ - ak.to_numpy(ak.flatten(data[detMod][f"{type}{var}"])) - for detMod in args.detectorModels - ], - bins=30, - label=[ - registry.get(detMod).get_name(args) for detMod in args.detectorModels - ], - range=(-xlims[type][var], xlims[type][var]) if var in xlims[type] else None, - ) - plt.ylabel("Frequency") - plt.title(f"Diff DetMods $\mathtt{{{type}}}$: {var}") - plt.legend() - plt.show() +if args.detmods: + xlims = {type: None for type in trackType} + xlims["SiTrack"] = {"D0": 0.03, "Omega": 0.0002} + xlims["CluTrack"] = {"D0": .8, "Omega": 0.00025} + for type in trackType: + for var in ["D0", "Omega"]: + plt.figure() + plt.grid( + True, which="both", linestyle="--", linewidth=linewidth, alpha=plotGridAlpha + ) + plt.hist( + x=[ + ak.to_numpy(ak.flatten(data[detMod][f"{type}{var}"])) + for detMod in args.detectorModels + ], + bins=30, + label=[ + rf"$\mathtt{{{registry.get(detMod).get_name(args)}}}$" for detMod in args.detectorModels + ], + range=(-xlims[type][var], xlims[type][var]) if var in xlims[type] else None, + ) + plt.ylabel("Frequency") + plt.title(f"Diff DetMods $\mathtt{{{type}}}$: {var}") + plt.legend() + plt.show() From 0353f97d7bb4447967da238ead1f9071d6b3eae2 Mon Sep 17 00:00:00 2001 From: Victor Schwan Date: Fri, 4 Jul 2025 13:36:19 +0200 Subject: [PATCH 25/38] Update plot params: e.g. colorblind colors --- Tracking/options/histTrackParams.py | 5 +++-- 1 file changed, 3 insertions(+), 2 deletions(-) diff --git a/Tracking/options/histTrackParams.py b/Tracking/options/histTrackParams.py index ec6b2d14..dfeee417 100644 --- a/Tracking/options/histTrackParams.py +++ b/Tracking/options/histTrackParams.py @@ -9,12 +9,13 @@ import uproot from commonArgParsing import add_common_args, detModNames, registry +plt.style.use("seaborn-v0_8-colorblind") # general plotting options -labelsize = 22 +labelsize = 24 linewidth = 1.5 majorTickSize = 10 -plotGridAlpha = 0.7 +plotGridAlpha = .7 params = { "xtick.direction": "in", "ytick.direction": "in", From 0ea8168f724d84123667d7cc15345ed703d49743 Mon Sep 17 00:00:00 2001 From: Victor Schwan Date: Fri, 4 Jul 2025 13:37:03 +0200 Subject: [PATCH 26/38] uniform format for floats --- Tracking/options/histTrackParams.py | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/Tracking/options/histTrackParams.py b/Tracking/options/histTrackParams.py index dfeee417..d77bc461 100644 --- a/Tracking/options/histTrackParams.py +++ b/Tracking/options/histTrackParams.py @@ -83,7 +83,7 @@ if args.track: # plot data for detMod in args.detectorModels: - for group, xlim in zip([varSpread, varSimilar], [1.5, 0.0015]): + for group, xlim in zip([varSpread, varSimilar], [1.5, .0015]): plt.figure() plt.grid( True, which="both", linestyle="--", linewidth=linewidth, alpha=plotGridAlpha @@ -105,8 +105,8 @@ if args.detmods: xlims = {type: None for type in trackType} - xlims["SiTrack"] = {"D0": 0.03, "Omega": 0.0002} - xlims["CluTrack"] = {"D0": .8, "Omega": 0.00025} + xlims["SiTrack"] = {"D0": .03, "Omega": .0002} + xlims["CluTrack"] = {"D0": .8, "Omega": .00025} for type in trackType: for var in ["D0", "Omega"]: plt.figure() From 898a8e9f3ca8c231ff1ba09b4a3542f3a545f1f7 Mon Sep 17 00:00:00 2001 From: Victor Schwan Date: Fri, 4 Jul 2025 13:43:47 +0200 Subject: [PATCH 27/38] Outsource plotting params dict --- Tracking/options/histTrackParams.py | 24 ++---------------------- Tracking/options/plotting.py | 22 ++++++++++++++++++++++ 2 files changed, 24 insertions(+), 22 deletions(-) create mode 100644 Tracking/options/plotting.py diff --git a/Tracking/options/histTrackParams.py b/Tracking/options/histTrackParams.py index d77bc461..5884ee8c 100644 --- a/Tracking/options/histTrackParams.py +++ b/Tracking/options/histTrackParams.py @@ -9,30 +9,10 @@ import uproot from commonArgParsing import add_common_args, detModNames, registry +from plotting import linewidth, params, plotGridAlpha + plt.style.use("seaborn-v0_8-colorblind") -# general plotting options -labelsize = 24 -linewidth = 1.5 -majorTickSize = 10 -plotGridAlpha = .7 -params = { - "xtick.direction": "in", - "ytick.direction": "in", - "xtick.top": True, - "ytick.right": True, - "xtick.major.size": majorTickSize, # Tick length - "ytick.major.size": majorTickSize, - "xtick.major.width": linewidth, # Tick line width - "ytick.major.width": linewidth, - "axes.linewidth": linewidth, - "legend.fontsize": labelsize, # "x-large", - "axes.labelsize": labelsize, # "x-large", - "axes.titlesize": labelsize, # "x-large", - "xtick.labelsize": labelsize, # "x-large", - "ytick.labelsize": labelsize, # "x-large", - "figure.autolayout": True, -} #'figure.figsize': (15, 5), mpl.rcParams.update(params) parser = add_common_args(ArgumentParser()) diff --git a/Tracking/options/plotting.py b/Tracking/options/plotting.py new file mode 100644 index 00000000..ecaf7cf2 --- /dev/null +++ b/Tracking/options/plotting.py @@ -0,0 +1,22 @@ +# general plotting options +labelsize = 24 +linewidth = 1.5 +majorTickSize = 10 +plotGridAlpha = .7 +params = { + "xtick.direction": "in", + "ytick.direction": "in", + "xtick.top": True, + "ytick.right": True, + "xtick.major.size": majorTickSize, # Tick length + "ytick.major.size": majorTickSize, + "xtick.major.width": linewidth, # Tick line width + "ytick.major.width": linewidth, + "axes.linewidth": linewidth, + "legend.fontsize": labelsize, # "x-large", + "axes.labelsize": labelsize, # "x-large", + "axes.titlesize": labelsize, # "x-large", + "xtick.labelsize": labelsize, # "x-large", + "ytick.labelsize": labelsize, # "x-large", + "figure.autolayout": True, +} #'figure.figsize': (15, 5), From 373558d25b4c53e6bd3d4e535d151adb66a0e6b5 Mon Sep 17 00:00:00 2001 From: Victor Schwan Date: Fri, 4 Jul 2025 13:45:33 +0200 Subject: [PATCH 28/38] make ruff happy --- Tracking/options/histTrackParams.py | 39 +++++++++++++++++++++-------- Tracking/options/plotting.py | 2 +- 2 files changed, 30 insertions(+), 11 deletions(-) diff --git a/Tracking/options/histTrackParams.py b/Tracking/options/histTrackParams.py index 5884ee8c..b292a087 100644 --- a/Tracking/options/histTrackParams.py +++ b/Tracking/options/histTrackParams.py @@ -17,8 +17,14 @@ parser = add_common_args(ArgumentParser()) plot_opts = parser.add_argument_group("Plotting opts", "which plots should be shown") -plot_opts.add_argument("--track", action="store_true", help="Show difference between Silicon and Clupatra tracks") -plot_opts.add_argument("--detmods", action="store_true", help="Show difference between detector models") +plot_opts.add_argument( + "--track", + action="store_true", + help="Show difference between Silicon and Clupatra tracks", +) +plot_opts.add_argument( + "--detmods", action="store_true", help="Show difference between detector models" +) args = parser.parse_known_args()[0] ################################# @@ -63,10 +69,14 @@ if args.track: # plot data for detMod in args.detectorModels: - for group, xlim in zip([varSpread, varSimilar], [1.5, .0015]): + for group, xlim in zip([varSpread, varSimilar], [1.5, 0.0015]): plt.figure() plt.grid( - True, which="both", linestyle="--", linewidth=linewidth, alpha=plotGridAlpha + True, + which="both", + linestyle="--", + linewidth=linewidth, + alpha=plotGridAlpha, ) plt.hist( x=[ @@ -79,19 +89,25 @@ ) plt.xlabel(rf"$\Delta$ Si-Clu") plt.ylabel("Frequency") - plt.title(rf"Diff Si-Clu in $\mathtt{{{registry.get(detMod).get_name(args)}}}$: {','.join(group)}") + plt.title( + rf"Diff Si-Clu in $\mathtt{{{registry.get(detMod).get_name(args)}}}$: {','.join(group)}" + ) plt.legend() plt.show() if args.detmods: xlims = {type: None for type in trackType} - xlims["SiTrack"] = {"D0": .03, "Omega": .0002} - xlims["CluTrack"] = {"D0": .8, "Omega": .00025} + xlims["SiTrack"] = {"D0": 0.03, "Omega": 0.0002} + xlims["CluTrack"] = {"D0": 0.8, "Omega": 0.00025} for type in trackType: for var in ["D0", "Omega"]: plt.figure() plt.grid( - True, which="both", linestyle="--", linewidth=linewidth, alpha=plotGridAlpha + True, + which="both", + linestyle="--", + linewidth=linewidth, + alpha=plotGridAlpha, ) plt.hist( x=[ @@ -100,9 +116,12 @@ ], bins=30, label=[ - rf"$\mathtt{{{registry.get(detMod).get_name(args)}}}$" for detMod in args.detectorModels + rf"$\mathtt{{{registry.get(detMod).get_name(args)}}}$" + for detMod in args.detectorModels ], - range=(-xlims[type][var], xlims[type][var]) if var in xlims[type] else None, + range=(-xlims[type][var], xlims[type][var]) + if var in xlims[type] + else None, ) plt.ylabel("Frequency") plt.title(f"Diff DetMods $\mathtt{{{type}}}$: {var}") diff --git a/Tracking/options/plotting.py b/Tracking/options/plotting.py index ecaf7cf2..44514509 100644 --- a/Tracking/options/plotting.py +++ b/Tracking/options/plotting.py @@ -2,7 +2,7 @@ labelsize = 24 linewidth = 1.5 majorTickSize = 10 -plotGridAlpha = .7 +plotGridAlpha = 0.7 params = { "xtick.direction": "in", "ytick.direction": "in", From 55d0cfceb11d6a921629f1ce9147f7c31bc38ae2 Mon Sep 17 00:00:00 2001 From: Victor Schwan Date: Fri, 4 Jul 2025 13:49:48 +0200 Subject: [PATCH 29/38] minor: add comments --- Tracking/options/histTrackParams.py | 22 ++++++++++++++++++++-- 1 file changed, 20 insertions(+), 2 deletions(-) diff --git a/Tracking/options/histTrackParams.py b/Tracking/options/histTrackParams.py index b292a087..d72cf2a9 100644 --- a/Tracking/options/histTrackParams.py +++ b/Tracking/options/histTrackParams.py @@ -15,7 +15,15 @@ mpl.rcParams.update(params) + +############################################# +# arg parsing +############################################# + +# import common args parser = add_common_args(ArgumentParser()) + +# add plotting options plot_opts = parser.add_argument_group("Plotting opts", "which plots should be shown") plot_opts.add_argument( "--track", @@ -25,17 +33,23 @@ plot_opts.add_argument( "--detmods", action="store_true", help="Show difference between detector models" ) + +# parse args args = parser.parse_known_args()[0] -################################# +############################################# # commands to access cov Matrix -################################# +############################################# # import ROOT # ROOT.gInterpreter.LoadFile("edm4hep/utils/cov_matrix_utils.h") # # ... # edm4hep.utils.get_cov_value(cov_m, edm4hep.TrackParams.d0, edm4hep.TrackParams.d0) + +############################################# # Lists to build branch names to be analyzed +############################################# + trackType = ["SiTrack", "CluTrack"] varSimilar = ["Phi", "Omega", "TanL"] varSpread = ["D0", "Z0"] @@ -66,6 +80,10 @@ - data[detMod][f"{trackType[1]}{var}"] ) +############################################# +# plotting +############################################# + if args.track: # plot data for detMod in args.detectorModels: From ed2efce3c6a2be2f206dcb30f0c4fa734ef46522 Mon Sep 17 00:00:00 2001 From: Victor Schwan Date: Fri, 4 Jul 2025 14:00:15 +0200 Subject: [PATCH 30/38] minor: update env var name --- Tracking/options/ExtractTrackParams.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/Tracking/options/ExtractTrackParams.py b/Tracking/options/ExtractTrackParams.py index a1262ca1..308628ae 100644 --- a/Tracking/options/ExtractTrackParams.py +++ b/Tracking/options/ExtractTrackParams.py @@ -7,7 +7,7 @@ from k4FWCore import ApplicationMgr, IOSvc from k4FWCore.parseArgs import parser -sys.path.append(os.getenv("trckOptDir")) +sys.path.append(os.getenv("pytrkDir")) from commonArgParsing import add_common_args, detModNames, registry args = add_common_args(parser).parse_known_args()[0] From 712f1f51d804badd36e4fbf29e79c0aa3dffe6f8 Mon Sep 17 00:00:00 2001 From: Victor Schwan Date: Fri, 4 Jul 2025 14:10:39 +0200 Subject: [PATCH 31/38] use type alias for edm4hep::TrackParams --- Tracking/components/ExtractTrackParams.cpp | 23 +++++++++++----------- 1 file changed, 11 insertions(+), 12 deletions(-) diff --git a/Tracking/components/ExtractTrackParams.cpp b/Tracking/components/ExtractTrackParams.cpp index 4fbb827a..57491ca5 100644 --- a/Tracking/components/ExtractTrackParams.cpp +++ b/Tracking/components/ExtractTrackParams.cpp @@ -13,6 +13,7 @@ using FloatColl = podio::UserDataCollection; using TrackColl = edm4hep::TrackCollection; using TS = edm4hep::TrackState; +using TP = edm4hep::TrackParams; struct TrackParamExtractor final : k4FWCore::MultiTransformer(edm4hep::TrackParams::d0)>(siColls).push_back(oSiTrackStateIP->D0); - std::get(edm4hep::TrackParams::phi)>(siColls).push_back(oSiTrackStateIP->phi); - std::get(edm4hep::TrackParams::omega)>(siColls).push_back(oSiTrackStateIP->omega); - std::get(edm4hep::TrackParams::z0)>(siColls).push_back(oSiTrackStateIP->Z0); - std::get(edm4hep::TrackParams::tanLambda)>(siColls).push_back( - oSiTrackStateIP->tanLambda); + std::get(TP::d0)>(siColls).push_back(oSiTrackStateIP->D0); + std::get(TP::phi)>(siColls).push_back(oSiTrackStateIP->phi); + std::get(TP::omega)>(siColls).push_back(oSiTrackStateIP->omega); + std::get(TP::z0)>(siColls).push_back(oSiTrackStateIP->Z0); + std::get(TP::tanLambda)>(siColls).push_back(oSiTrackStateIP->tanLambda); } // Process CluTrack const auto oCluTrackStateIP = getOTrackAtIP(inCluTracks[i], "CluTrack"); if (oCluTrackStateIP.has_value()) { - std::get(edm4hep::TrackParams::d0)>(cluColls).push_back(oCluTrackStateIP->D0); - std::get(edm4hep::TrackParams::phi)>(cluColls).push_back(oCluTrackStateIP->phi); - std::get(edm4hep::TrackParams::omega)>(cluColls).push_back(oCluTrackStateIP->omega); - std::get(edm4hep::TrackParams::z0)>(cluColls).push_back(oCluTrackStateIP->Z0); - std::get(edm4hep::TrackParams::tanLambda)>(cluColls).push_back( - oCluTrackStateIP->tanLambda); + std::get(TP::d0)>(cluColls).push_back(oCluTrackStateIP->D0); + std::get(TP::phi)>(cluColls).push_back(oCluTrackStateIP->phi); + std::get(TP::omega)>(cluColls).push_back(oCluTrackStateIP->omega); + std::get(TP::z0)>(cluColls).push_back(oCluTrackStateIP->Z0); + std::get(TP::tanLambda)>(cluColls).push_back(oCluTrackStateIP->tanLambda); } } From c2cbe4bef84590691e0478ca244ab2719bd91d78 Mon Sep 17 00:00:00 2001 From: Victor Schwan Date: Fri, 15 Aug 2025 16:17:25 +0200 Subject: [PATCH 32/38] Minor: avoid None return type --- Tracking/options/ExtractTrackParams.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/Tracking/options/ExtractTrackParams.py b/Tracking/options/ExtractTrackParams.py index 308628ae..8af45bd9 100644 --- a/Tracking/options/ExtractTrackParams.py +++ b/Tracking/options/ExtractTrackParams.py @@ -7,7 +7,7 @@ from k4FWCore import ApplicationMgr, IOSvc from k4FWCore.parseArgs import parser -sys.path.append(os.getenv("pytrkDir")) +sys.path.append(os.environ["pytrkDir"]) from commonArgParsing import add_common_args, detModNames, registry args = add_common_args(parser).parse_known_args()[0] From 94d0b2ec2fc86991c16f2b26fa017237e53062f4 Mon Sep 17 00:00:00 2001 From: Victor Schwan Date: Fri, 15 Aug 2025 16:18:08 +0200 Subject: [PATCH 33/38] move commonArgParsing out of this repo to shared --- Tracking/options/ExtractTrackParams.py | 5 +---- 1 file changed, 1 insertion(+), 4 deletions(-) diff --git a/Tracking/options/ExtractTrackParams.py b/Tracking/options/ExtractTrackParams.py index 8af45bd9..97f5a65b 100644 --- a/Tracking/options/ExtractTrackParams.py +++ b/Tracking/options/ExtractTrackParams.py @@ -1,15 +1,12 @@ import os -import sys from pathlib import Path +from commonArgParsing import add_common_args, detModNames, registry from Configurables import TrackParamExtractor from Gaudi.Configuration import INFO, VERBOSE from k4FWCore import ApplicationMgr, IOSvc from k4FWCore.parseArgs import parser -sys.path.append(os.environ["pytrkDir"]) -from commonArgParsing import add_common_args, detModNames, registry - args = add_common_args(parser).parse_known_args()[0] assert len(args.detectorModels) == 1, ( f"Only provide one detector model! You provided {args.detectorModels}" From 21f0abc0ea6f569716fb2973b8a29a16f04160dd Mon Sep 17 00:00:00 2001 From: Victor Schwan Date: Fri, 15 Aug 2025 17:46:02 +0200 Subject: [PATCH 34/38] new structure of data dir with input /output subdirs --- Tracking/options/ExtractTrackParams.py | 13 +++++++++---- 1 file changed, 9 insertions(+), 4 deletions(-) diff --git a/Tracking/options/ExtractTrackParams.py b/Tracking/options/ExtractTrackParams.py index 97f5a65b..4fa6409d 100644 --- a/Tracking/options/ExtractTrackParams.py +++ b/Tracking/options/ExtractTrackParams.py @@ -16,17 +16,22 @@ fileSuffix = ".edm4hep.root" procName = "TrackParamExtractor" basePath = Path(os.getenv("prmDir", Path.home() / "promotion")) +in_out_base_path = basePath / "data" / procName corePath = f"{args.version}_{detModNames[args.detectorModels]}" # assert that the input path exists -inputPath = ( - basePath / "code/ILDConfig/StandardConfig/production/data" / f"{corePath}_REC" -).with_suffix(fileSuffix) +inputPath = (in_out_base_path / "input_data" / f"{corePath}_REC").with_suffix( + fileSuffix +) assert inputPath.exists(), f"ERROR: The input path ({inputPath}) does not exist!" iosvc = IOSvc() iosvc.Input = str(inputPath) -iosvc.Output = str((basePath / "data" / procName / corePath).with_suffix(fileSuffix)) +iosvc.Output = str( + (in_out_base_path / "out_track_params" / f"{corePath}_track_params").with_suffix( + fileSuffix + ) +) # iosvc.outputCommands = ["drop *", "keep SiTrackPhi"] printer = TrackParamExtractor(procName, nStars=40) From 145372e47acc94cd72f53f77130ddaafcb4f8496 Mon Sep 17 00:00:00 2001 From: Victor Schwan Date: Fri, 15 Aug 2025 17:58:28 +0200 Subject: [PATCH 35/38] FORMAT/REFAC: make pylint happy and obey naming convention --- Tracking/options/ExtractTrackParams.py | 40 +++++++++++++------------- 1 file changed, 20 insertions(+), 20 deletions(-) diff --git a/Tracking/options/ExtractTrackParams.py b/Tracking/options/ExtractTrackParams.py index 4fa6409d..c6163a45 100644 --- a/Tracking/options/ExtractTrackParams.py +++ b/Tracking/options/ExtractTrackParams.py @@ -7,43 +7,43 @@ from k4FWCore import ApplicationMgr, IOSvc from k4FWCore.parseArgs import parser -args = add_common_args(parser).parse_known_args()[0] -assert len(args.detectorModels) == 1, ( - f"Only provide one detector model! You provided {args.detectorModels}" +ARGS = add_common_args(parser).parse_known_args()[0] +assert len(ARGS.detectorModels) == 1, ( + f"Only provide one detector model! You provided {ARGS.detectorModels}" ) -args.detectorModels = args.detectorModels[0] +ARGS.detectorModels = ARGS.detectorModels[0] -fileSuffix = ".edm4hep.root" -procName = "TrackParamExtractor" -basePath = Path(os.getenv("prmDir", Path.home() / "promotion")) -in_out_base_path = basePath / "data" / procName -corePath = f"{args.version}_{detModNames[args.detectorModels]}" +FILE_SUFFIX = ".edm4hep.root" +PROCESSOR_NAME = "TrackParamExtractor" +BASE_PATH = Path(os.getenv("prmDir", Path.home() / "promotion")) +IN_OUT_BASE_PATH = BASE_PATH / "data" / PROCESSOR_NAME +CORE_PATH = f"{ARGS.version}_{detModNames[ARGS.detectorModels]}" # assert that the input path exists -inputPath = (in_out_base_path / "input_data" / f"{corePath}_REC").with_suffix( - fileSuffix +INPUT_PATH = (IN_OUT_BASE_PATH / "input_data" / f"{CORE_PATH}_REC").with_suffix( + FILE_SUFFIX ) -assert inputPath.exists(), f"ERROR: The input path ({inputPath}) does not exist!" +assert INPUT_PATH.exists(), f"ERROR: The input path ({INPUT_PATH}) does not exist!" iosvc = IOSvc() -iosvc.Input = str(inputPath) +iosvc.Input = str(INPUT_PATH) iosvc.Output = str( - (in_out_base_path / "out_track_params" / f"{corePath}_track_params").with_suffix( - fileSuffix + (IN_OUT_BASE_PATH / "out_track_params" / f"{CORE_PATH}_track_params").with_suffix( + FILE_SUFFIX ) ) # iosvc.outputCommands = ["drop *", "keep SiTrackPhi"] -printer = TrackParamExtractor(procName, nStars=40) +printer = TrackParamExtractor(PROCESSOR_NAME, nStars=40) printer.OutputLevel = VERBOSE # the collection name with the SiTracks differs between ILC and FCC models -if registry.get(args.detectorModels).at_fcc: +if registry.get(ARGS.detectorModels).at_fcc: printer.InputSiTracks = ["SiTracksCT"] - siTrackCollName = "SiTracksCT" + SI_TRACK_COLL_NAME = "SiTracksCT" else: - siTrackCollName = "SiTracks" -iosvc.CollectionNames = ["ClupatraTracks", siTrackCollName] + SI_TRACK_COLL_NAME = "SiTracks" +iosvc.CollectionNames = ["ClupatraTracks", SI_TRACK_COLL_NAME] ApplicationMgr( From b91a64ea240937c3a310dadc4ca8ab4f7e792bb2 Mon Sep 17 00:00:00 2001 From: Victor Schwan Date: Fri, 4 Jul 2025 17:45:35 +0200 Subject: [PATCH 36/38] Extract Uncertainties of Track Params as well --- Tracking/components/ExtractTrackParams.cpp | 60 +++- Tracking/options/histTrackParams.py | 313 +++++++++++++++------ Tracking/options/plotting.py | 28 +- Tracking/options/utils.py | 35 +++ 4 files changed, 321 insertions(+), 115 deletions(-) create mode 100644 Tracking/options/utils.py diff --git a/Tracking/components/ExtractTrackParams.cpp b/Tracking/components/ExtractTrackParams.cpp index 57491ca5..03652ee9 100644 --- a/Tracking/components/ExtractTrackParams.cpp +++ b/Tracking/components/ExtractTrackParams.cpp @@ -14,10 +14,13 @@ using FloatColl = podio::UserDataCollection; using TrackColl = edm4hep::TrackCollection; using TS = edm4hep::TrackState; using TP = edm4hep::TrackParams; +using edm4hep::utils::detail::to_index; struct TrackParamExtractor final : k4FWCore::MultiTransformer(const TrackColl&, const TrackColl&)> { + FloatColl, FloatColl, FloatColl, FloatColl, FloatColl, FloatColl, FloatColl, + FloatColl, FloatColl, FloatColl, FloatColl, FloatColl, FloatColl>( + const TrackColl&, const TrackColl&)> { TrackParamExtractor(const std::string& name, ISvcLocator* svcLoc) : MultiTransformer(name, svcLoc, { @@ -35,14 +38,25 @@ struct TrackParamExtractor final KeyValues("OutCollCluOmega", {"CluTrackOmega"}), KeyValues("OutCollCluZ0", {"CluTrackZ0"}), KeyValues("OutCollCluTanL", {"CluTrackTanL"}), - + // uncertainties + KeyValues("OutCollSiUncD0", {"SiTrackUncD0"}), + KeyValues("OutCollSiUncPhi", {"SiTrackUncPhi"}), + KeyValues("OutCollSiUncOmega", {"SiTrackUncOmega"}), + KeyValues("OutCollSiUncZ0", {"SiTrackUncZ0"}), + KeyValues("OutCollSiUncTanL", {"SiTrackUncTanL"}), + KeyValues("OutCollCluUncD0", {"CluTrackUncD0"}), + KeyValues("OutCollCluUncPhi", {"CluTrackUncPhi"}), + KeyValues("OutCollCluUncOmega", {"CluTrackUncOmega"}), + KeyValues("OutCollCluUncZ0", {"CluTrackUncZ0"}), + KeyValues("OutCollCluUncTanL", {"CluTrackUncTanL"}), }) {} // This is the function that will be called to transform the data // Note that the function has to be const, as well as the collections // we get from the input std::tuple + FloatColl, FloatColl, FloatColl, FloatColl, FloatColl, FloatColl, FloatColl, FloatColl, FloatColl, + FloatColl, FloatColl> operator()(const TrackColl& inSiTracks, const TrackColl& inCluTracks) const override { printInStars(this, "New Event", n_stars); @@ -59,30 +73,46 @@ struct TrackParamExtractor final // Process both collections (SiTracks and ClupatraTracks) at the same time std::tuple siColls; std::tuple cluColls; + std::tuple siUncColls; + std::tuple cluUncColls; + for (size_t i = 0; i < inSiTracks.size(); ++i) { // Process SiTrack const auto oSiTrackStateIP = getOTrackAtIP(inSiTracks[i], "SiTrack"); if (oSiTrackStateIP.has_value()) { - std::get(TP::d0)>(siColls).push_back(oSiTrackStateIP->D0); - std::get(TP::phi)>(siColls).push_back(oSiTrackStateIP->phi); - std::get(TP::omega)>(siColls).push_back(oSiTrackStateIP->omega); - std::get(TP::z0)>(siColls).push_back(oSiTrackStateIP->Z0); - std::get(TP::tanLambda)>(siColls).push_back(oSiTrackStateIP->tanLambda); + // values + std::get(siColls).push_back(oSiTrackStateIP->D0); + std::get(siColls).push_back(oSiTrackStateIP->phi); + std::get(siColls).push_back(oSiTrackStateIP->omega); + std::get(siColls).push_back(oSiTrackStateIP->Z0); + std::get(siColls).push_back(oSiTrackStateIP->tanLambda); + // uncertainties + std::get(siUncColls).push_back(getSigmaVar(*oSiTrackStateIP, TP::d0)); + std::get(siUncColls).push_back(getSigmaVar(*oSiTrackStateIP, TP::phi)); + std::get(siUncColls).push_back(getSigmaVar(*oSiTrackStateIP, TP::omega)); + std::get(siUncColls).push_back(getSigmaVar(*oSiTrackStateIP, TP::z0)); + std::get(siUncColls).push_back(getSigmaVar(*oSiTrackStateIP, TP::tanLambda)); } // Process CluTrack const auto oCluTrackStateIP = getOTrackAtIP(inCluTracks[i], "CluTrack"); if (oCluTrackStateIP.has_value()) { - std::get(TP::d0)>(cluColls).push_back(oCluTrackStateIP->D0); - std::get(TP::phi)>(cluColls).push_back(oCluTrackStateIP->phi); - std::get(TP::omega)>(cluColls).push_back(oCluTrackStateIP->omega); - std::get(TP::z0)>(cluColls).push_back(oCluTrackStateIP->Z0); - std::get(TP::tanLambda)>(cluColls).push_back(oCluTrackStateIP->tanLambda); + std::get(cluColls).push_back(oCluTrackStateIP->D0); + std::get(cluColls).push_back(oCluTrackStateIP->phi); + std::get(cluColls).push_back(oCluTrackStateIP->omega); + std::get(cluColls).push_back(oCluTrackStateIP->Z0); + std::get(cluColls).push_back(oCluTrackStateIP->tanLambda); + // uncertainties + std::get(cluUncColls).push_back(getSigmaVar(*oCluTrackStateIP, TP::d0)); + std::get(cluUncColls).push_back(getSigmaVar(*oCluTrackStateIP, TP::phi)); + std::get(cluUncColls).push_back(getSigmaVar(*oCluTrackStateIP, TP::omega)); + std::get(cluUncColls).push_back(getSigmaVar(*oCluTrackStateIP, TP::z0)); + std::get(cluUncColls).push_back(getSigmaVar(*oCluTrackStateIP, TP::tanLambda)); } } - return std::tuple_cat(std::move(siColls), std::move(cluColls)); + return std::tuple_cat(std::move(siColls), std::move(cluColls), std::move(siUncColls), std::move(cluUncColls)); }; private: @@ -100,5 +130,7 @@ struct TrackParamExtractor final return std::nullopt; } } + + float getSigmaVar(const edm4hep::TrackState& ts, const TP var) const { return std::sqrt(ts.getCovMatrix(var, var)); } }; DECLARE_COMPONENT(TrackParamExtractor) \ No newline at end of file diff --git a/Tracking/options/histTrackParams.py b/Tracking/options/histTrackParams.py index d72cf2a9..5158c8bd 100644 --- a/Tracking/options/histTrackParams.py +++ b/Tracking/options/histTrackParams.py @@ -1,3 +1,7 @@ +############################################# +# call with `python3` NOT `k4run` +############################################# + from argparse import ArgumentParser from itertools import product from os import getenv @@ -6,15 +10,20 @@ import awkward as ak import matplotlib as mpl import matplotlib.pyplot as plt +import numpy as np import uproot from commonArgParsing import add_common_args, detModNames, registry -from plotting import linewidth, params, plotGridAlpha +from plotting import my_line_styles +from utils import is_outlier -plt.style.use("seaborn-v0_8-colorblind") +plt.style.use(["seaborn-v0_8-colorblind", "vics_basic"]) -mpl.rcParams.update(params) +threshold_outlier_detection = 4 +my_hist_type = "step" +my_line_width = 2.5 +my_n_bins = 30 ############################################# # arg parsing @@ -22,7 +31,14 @@ # import common args parser = add_common_args(ArgumentParser()) - +parser.add_argument( + "--mode", + choices=["nominal", "uncertainty", "both"], + default="nominal", + help="Choose what to process: 'nominal' (default), 'uncertainty', or 'both'.", +) +parser.add_argument("--debug", action="store_true") +parser.add_argument("--rm-outliers", action="store_true") # add plotting options plot_opts = parser.add_argument_group("Plotting opts", "which plots should be shown") plot_opts.add_argument( @@ -36,28 +52,24 @@ # parse args args = parser.parse_known_args()[0] - -############################################# -# commands to access cov Matrix -############################################# -# import ROOT -# ROOT.gInterpreter.LoadFile("edm4hep/utils/cov_matrix_utils.h") -# # ... -# edm4hep.utils.get_cov_value(cov_m, edm4hep.TrackParams.d0, edm4hep.TrackParams.d0) +assert len(parser.parse_known_args()[1]) == 0, ( + f"Unknown args provided: {parser.parse_known_args()[1]}" +) ############################################# # Lists to build branch names to be analyzed ############################################# -trackType = ["SiTrack", "CluTrack"] -varSimilar = ["Phi", "Omega", "TanL"] -varSpread = ["D0", "Z0"] -varNames = varSimilar + varSpread +track_types = ["SiTrack", "CluTrack"] +var_similar = ["Phi", "Omega", "TanL"] +var_spread = ["D0", "Z0"] +var_names = var_similar + var_spread +############################################# +# extract data from root file +############################################# data = {} - -# extract data for detMod in args.detectorModels: corePath = Path(f"{args.version}_{detModNames[detMod]}") @@ -67,81 +79,222 @@ # build vars based on above vars keys = [ - f"{trackName}{varName}" for trackName, varName in product(trackType, varNames) + f"{trackName}{varName}" + for trackName, varName in product(track_types, var_names) ] in_file = basePath / processor / corePath.with_suffix(".edm4hep.root") with uproot.open(str(in_file) + ":events") as events: - regex = f"/^({'|'.join(trackType)})({'|'.join(varNames)})$/" + # regex to match desired branch names + regex = ( + f"/^({'|'.join(track_types)}){'(Unc)?' if args.mode != 'nominal' else ''}" + + f"({'|'.join(var_names)})$/" + ) + if args.debug: + print(f"Regex to match branches we are interested in: {regex}") data[detMod] = events.arrays(filter_name=regex, library="pd") - for var in varNames: + for var in var_names: data[detMod][f"d_{var}"] = ( - data[detMod][f"{trackType[0]}{var}"] - - data[detMod][f"{trackType[1]}{var}"] + data[detMod][f"{track_types[0]}{var}"] + - data[detMod][f"{track_types[1]}{var}"] ) + if args.debug: + print("Matched branches are:") + print(data[detMod].columns) + ############################################# -# plotting +# plotting funcs ############################################# +def process_data_for_hist( + data, det_mod, var_column_name, rm_outliers, thresh_outlier_detection +): + # TODO: remove comment? + """ + Process the data for plotting, handling outliers if required. + + Args: + data: The full data dictionary. + det_mod: The detector model name. + var_column_name: The variable name. + rm_outliers: Boolean flag to remove outliers. + thresh_outlier_detection: Threshold for outlier detection. + + Returns: + numpy array: Processed data ready for plotting. + """ + data_array = ak.to_numpy(ak.flatten(data[det_mod][f"{var_column_name}"])) + + if rm_outliers: + # Remove outliers + return data_array[~is_outlier(data_array, thresh=thresh_outlier_detection)] + + return data_array + -if args.track: - # plot data - for detMod in args.detectorModels: - for group, xlim in zip([varSpread, varSimilar], [1.5, 0.0015]): - plt.figure() - plt.grid( - True, - which="both", - linestyle="--", - linewidth=linewidth, - alpha=plotGridAlpha, +def debug_info_outlier_removal( + data, thresh, args, var_list, fixed_data_key, det_mod_variable: bool +): + if args.rm_outliers and args.debug: + star_string = "*" * 50 + print(f"\n\n{star_string}\nOutlier threshold value is: {thresh}\n{star_string}") + for current_var in var_list: + data_to_clean = ( + data[current_var][fixed_data_key] + if det_mod_variable + else data[fixed_data_key][current_var] ) - plt.hist( - x=[ - ak.to_numpy(ak.flatten(data[detMod][f"d_{varName}"])) - for varName in group - ], - bins=30, - label=group, - range=(-xlim, xlim), + pre_processed_data = ak.to_numpy(ak.flatten(data_to_clean)) + print( + f"{np.sum(is_outlier(pre_processed_data))} outliers will be removed out of {len(pre_processed_data)} values for {current_var}" ) - plt.xlabel(rf"$\Delta$ Si-Clu") - plt.ylabel("Frequency") - plt.title( - rf"Diff Si-Clu in $\mathtt{{{registry.get(detMod).get_name(args)}}}$: {','.join(group)}" + + +# general plotting func +def plot_track_param_hist(data, thresh_outlier_detection, args, var, labels, hist_args): + if var: # TODO: support for var_group case + debug_info_outlier_removal( + data, + thresh_outlier_detection, + args, + args.detectorModels, + var, + det_mod_variable=True, + ) + with mpl.rc_context({"axes.titlesize": 18}): + plt.figure() + plt.hist( + **hist_args, + bins=my_n_bins, + histtype=my_hist_type, + linewidth=my_line_width, + linestyle=my_line_styles, + ) + plt.xlabel(labels["xlabel"] if "xlabel" in labels else None) + plt.ylabel("Frequency") + plt.suptitle(labels["suptitle"]) + plt.title(labels["title"]) + plt.legend() + plt.show() + + +# plotting func for collective plot of group of vars +def plot_track_param_hist_var_groups( + data, thresh_outlier_detection, args, det_mod, labels, group +): + hist_args_diff_detmods = { + "x": [ + process_data_for_hist( + data, + det_mod, + f"d_{var_name}", + args.rm_outliers, + thresh_outlier_detection, ) - plt.legend() - plt.show() - -if args.detmods: - xlims = {type: None for type in trackType} - xlims["SiTrack"] = {"D0": 0.03, "Omega": 0.0002} - xlims["CluTrack"] = {"D0": 0.8, "Omega": 0.00025} - for type in trackType: - for var in ["D0", "Omega"]: - plt.figure() - plt.grid( - True, - which="both", - linestyle="--", - linewidth=linewidth, - alpha=plotGridAlpha, + for var_name in group + ], + "label": group, + } + plot_track_param_hist( + data, + thresh_outlier_detection, + args, + None, + labels, + hist_args_diff_detmods, + ) + + +# plotting func for collective plot of all det mods +def plot_track_param_hist_diff_detmods( + data, thresh_outlier_detection, args, var, labels +): + hist_args_diff_detmods = { + "x": [ + process_data_for_hist( + data, det_mod, var, args.rm_outliers, thresh_outlier_detection ) - plt.hist( - x=[ - ak.to_numpy(ak.flatten(data[detMod][f"{type}{var}"])) - for detMod in args.detectorModels - ], - bins=30, - label=[ - rf"$\mathtt{{{registry.get(detMod).get_name(args)}}}$" - for detMod in args.detectorModels - ], - range=(-xlims[type][var], xlims[type][var]) - if var in xlims[type] - else None, + for det_mod in args.detectorModels + ], + "label": [ + registry.get(det_mod).get_name(args.detname) + for det_mod in args.detectorModels + ], + } + plot_track_param_hist( + data, + thresh_outlier_detection, + args, + var, + labels, + hist_args_diff_detmods, + ) + + +############################################# +# actual plotting (calling the funcs) +############################################# + +if args.mode != "nominal": # uncertainty or both + for track_type in [track_types[0]]: + for var in var_names: + # define labels + uncertainty_labels = { + "xlabel": r"$\sigma$", + "suptitle": rf"$\sigma$({var})", + "title": rf"single $\mu$, {track_type} {' (no outliers)' if args.rm_outliers else ''}", + } + # plotting ;) + plot_track_param_hist_diff_detmods( + data, + threshold_outlier_detection, + args, + f"{track_type}Unc{var}", + uncertainty_labels, ) - plt.ylabel("Frequency") - plt.title(f"Diff DetMods $\mathtt{{{type}}}$: {var}") - plt.legend() - plt.show() + + +if args.mode != "uncertainty": # nominal or both + if args.track: + for det_mod in args.detectorModels: + for group in [var_spread, var_similar]: + # define labels + nominal_var_groups_labels = { + "xlabel": r"$\Delta$ Si-Clu", + "suptitle": r"$\Delta$ Si-Clu:", + "title": rf"single $\mu$, {registry.get(det_mod).get_name(args.detname)}" + f" {' (no outliers)' if args.rm_outliers else ''}", + } + plot_track_param_hist_var_groups( + data, + threshold_outlier_detection, + args, + det_mod, + nominal_var_groups_labels, + group, + ) + + if args.detmods: + for track_type in track_types: + for var in ["D0", "Omega"]: + nominal_diff_det_mods_labels = { + "suptitle": f"Diff DetMods: {var}", + "title": rf"single $\mu$, {track_type}" + f" {' (no outliers)' if args.rm_outliers else ''}", + } + plot_track_param_hist_diff_detmods( + data, + threshold_outlier_detection, + args, + f"{track_type}{var}", + nominal_diff_det_mods_labels, + ) + + +############################################# +# commands to access cov Matrix +############################################# +# import ROOT +# ROOT.gInterpreter.LoadFile("edm4hep/utils/cov_matrix_utils.h") +# # ... +# edm4hep.utils.get_cov_value(cov_m, edm4hep.TrackParams.d0, edm4hep.TrackParams.d0) diff --git a/Tracking/options/plotting.py b/Tracking/options/plotting.py index 44514509..39ddb2bc 100644 --- a/Tracking/options/plotting.py +++ b/Tracking/options/plotting.py @@ -1,22 +1,8 @@ # general plotting options -labelsize = 24 -linewidth = 1.5 -majorTickSize = 10 -plotGridAlpha = 0.7 -params = { - "xtick.direction": "in", - "ytick.direction": "in", - "xtick.top": True, - "ytick.right": True, - "xtick.major.size": majorTickSize, # Tick length - "ytick.major.size": majorTickSize, - "xtick.major.width": linewidth, # Tick line width - "ytick.major.width": linewidth, - "axes.linewidth": linewidth, - "legend.fontsize": labelsize, # "x-large", - "axes.labelsize": labelsize, # "x-large", - "axes.titlesize": labelsize, # "x-large", - "xtick.labelsize": labelsize, # "x-large", - "ytick.labelsize": labelsize, # "x-large", - "figure.autolayout": True, -} #'figure.figsize': (15, 5), +my_line_styles = [ + "solid", + "dashed", + (0, (3, 1, 1, 1)), # "densely dashdotted" + (0, (3, 1, 1, 1, 1, 1)), # "densely dashdotdotted", + "dashdot", +] diff --git a/Tracking/options/utils.py b/Tracking/options/utils.py new file mode 100644 index 00000000..8bfb7f9a --- /dev/null +++ b/Tracking/options/utils.py @@ -0,0 +1,35 @@ +import numpy as np + + +def is_outlier(points, thresh=3.5): + """ + Returns a boolean array with True if points are outliers and False + otherwise. + + Parameters: + ----------- + points : An numobservations by numdimensions array of observations + thresh : The modified z-score to use as a threshold. Observations with + a modified z-score (based on the median absolute deviation) greater + than this value will be classified as outliers. + + Returns: + -------- + mask : A numobservations-length boolean array. + + References: + ---------- + Boris Iglewicz and David Hoaglin (1993), "Volume 16: How to Detect and + Handle Outliers", The ASQC Basic References in Quality Control: + Statistical Techniques, Edward F. Mykytka, Ph.D., Editor. + """ + if len(points.shape) == 1: + points = points[:, None] + median = np.median(points, axis=0) + diff = np.sum((points - median) ** 2, axis=-1) + diff = np.sqrt(diff) + med_abs_deviation = np.median(diff) + + modified_z_score = 0.6745 * diff / med_abs_deviation + + return modified_z_score > thresh From 1d43c6fb1b9c70027e4737be517d3abd59751ee6 Mon Sep 17 00:00:00 2001 From: Thomas Madlener Date: Wed, 5 Nov 2025 16:41:38 +0100 Subject: [PATCH 37/38] Make python files adhere to ruff formatting rules --- Tracking/options/ExtractTrackParams.py | 8 ++----- Tracking/options/histTrackParams.py | 29 +++++++------------------- 2 files changed, 9 insertions(+), 28 deletions(-) diff --git a/Tracking/options/ExtractTrackParams.py b/Tracking/options/ExtractTrackParams.py index c6163a45..6196beae 100644 --- a/Tracking/options/ExtractTrackParams.py +++ b/Tracking/options/ExtractTrackParams.py @@ -20,17 +20,13 @@ CORE_PATH = f"{ARGS.version}_{detModNames[ARGS.detectorModels]}" # assert that the input path exists -INPUT_PATH = (IN_OUT_BASE_PATH / "input_data" / f"{CORE_PATH}_REC").with_suffix( - FILE_SUFFIX -) +INPUT_PATH = (IN_OUT_BASE_PATH / "input_data" / f"{CORE_PATH}_REC").with_suffix(FILE_SUFFIX) assert INPUT_PATH.exists(), f"ERROR: The input path ({INPUT_PATH}) does not exist!" iosvc = IOSvc() iosvc.Input = str(INPUT_PATH) iosvc.Output = str( - (IN_OUT_BASE_PATH / "out_track_params" / f"{CORE_PATH}_track_params").with_suffix( - FILE_SUFFIX - ) + (IN_OUT_BASE_PATH / "out_track_params" / f"{CORE_PATH}_track_params").with_suffix(FILE_SUFFIX) ) # iosvc.outputCommands = ["drop *", "keep SiTrackPhi"] diff --git a/Tracking/options/histTrackParams.py b/Tracking/options/histTrackParams.py index 5158c8bd..47f59aea 100644 --- a/Tracking/options/histTrackParams.py +++ b/Tracking/options/histTrackParams.py @@ -78,10 +78,7 @@ basePath = Path(getenv("dtDir", str(Path.home() / "promotion" / "data"))) # build vars based on above vars - keys = [ - f"{trackName}{varName}" - for trackName, varName in product(track_types, var_names) - ] + keys = [f"{trackName}{varName}" for trackName, varName in product(track_types, var_names)] in_file = basePath / processor / corePath.with_suffix(".edm4hep.root") with uproot.open(str(in_file) + ":events") as events: @@ -95,8 +92,7 @@ data[detMod] = events.arrays(filter_name=regex, library="pd") for var in var_names: data[detMod][f"d_{var}"] = ( - data[detMod][f"{track_types[0]}{var}"] - - data[detMod][f"{track_types[1]}{var}"] + data[detMod][f"{track_types[0]}{var}"] - data[detMod][f"{track_types[1]}{var}"] ) if args.debug: print("Matched branches are:") @@ -106,9 +102,7 @@ ############################################# # plotting funcs ############################################# -def process_data_for_hist( - data, det_mod, var_column_name, rm_outliers, thresh_outlier_detection -): +def process_data_for_hist(data, det_mod, var_column_name, rm_outliers, thresh_outlier_detection): # TODO: remove comment? """ Process the data for plotting, handling outliers if required. @@ -179,9 +173,7 @@ def plot_track_param_hist(data, thresh_outlier_detection, args, var, labels, his # plotting func for collective plot of group of vars -def plot_track_param_hist_var_groups( - data, thresh_outlier_detection, args, det_mod, labels, group -): +def plot_track_param_hist_var_groups(data, thresh_outlier_detection, args, det_mod, labels, group): hist_args_diff_detmods = { "x": [ process_data_for_hist( @@ -206,20 +198,13 @@ def plot_track_param_hist_var_groups( # plotting func for collective plot of all det mods -def plot_track_param_hist_diff_detmods( - data, thresh_outlier_detection, args, var, labels -): +def plot_track_param_hist_diff_detmods(data, thresh_outlier_detection, args, var, labels): hist_args_diff_detmods = { "x": [ - process_data_for_hist( - data, det_mod, var, args.rm_outliers, thresh_outlier_detection - ) - for det_mod in args.detectorModels - ], - "label": [ - registry.get(det_mod).get_name(args.detname) + process_data_for_hist(data, det_mod, var, args.rm_outliers, thresh_outlier_detection) for det_mod in args.detectorModels ], + "label": [registry.get(det_mod).get_name(args.detname) for det_mod in args.detectorModels], } plot_track_param_hist( data, From e0a07ebdb7f077e54bd1d1853fe930f6054b3381 Mon Sep 17 00:00:00 2001 From: Thomas Madlener Date: Wed, 5 Nov 2025 21:06:29 +0100 Subject: [PATCH 38/38] Merge utility functionality into existing structure --- Tracking/components/ExtractTrackParams.cpp | 16 ++++++++++------ Tracking/components/TrackD0Printer.cpp | 14 +++++++++----- Tracking/components/printStars.cpp | 10 ---------- Tracking/components/printStars.h | 4 ---- Tracking/include/utils.hpp | 9 ++++++++- Tracking/src/utils.cpp | 8 ++++++++ 6 files changed, 35 insertions(+), 26 deletions(-) delete mode 100644 Tracking/components/printStars.cpp delete mode 100644 Tracking/components/printStars.h diff --git a/Tracking/components/ExtractTrackParams.cpp b/Tracking/components/ExtractTrackParams.cpp index 03652ee9..5d41e9c9 100644 --- a/Tracking/components/ExtractTrackParams.cpp +++ b/Tracking/components/ExtractTrackParams.cpp @@ -1,9 +1,13 @@ -#include "Gaudi/Property.h" -#include "edm4hep/TrackCollection.h" -#include "k4FWCore/Transformer.h" -#include "podio/UserDataCollection.h" -#include "printStars.h" +#include "utils.hpp" + #include +#include + +#include + +#include "k4FWCore/Transformer.h" + +#include "Gaudi/Property.h" #include #include @@ -133,4 +137,4 @@ struct TrackParamExtractor final float getSigmaVar(const edm4hep::TrackState& ts, const TP var) const { return std::sqrt(ts.getCovMatrix(var, var)); } }; -DECLARE_COMPONENT(TrackParamExtractor) \ No newline at end of file +DECLARE_COMPONENT(TrackParamExtractor) diff --git a/Tracking/components/TrackD0Printer.cpp b/Tracking/components/TrackD0Printer.cpp index 257afd03..d4a9104a 100644 --- a/Tracking/components/TrackD0Printer.cpp +++ b/Tracking/components/TrackD0Printer.cpp @@ -1,11 +1,15 @@ -#include "Gaudi/Property.h" -#include "GaudiKernel/MsgStream.h" +#include "utils.hpp" + #include "edm4hep/Track.h" #include "edm4hep/TrackCollection.h" #include "edm4hep/TrackState.h" -#include "k4FWCore/Consumer.h" #include "podio/RelationRange.h" -#include "printStars.h" + +#include "Gaudi/Property.h" +#include "GaudiKernel/MsgStream.h" + +#include "k4FWCore/Consumer.h" + #include #include @@ -85,4 +89,4 @@ struct TrackD0Printer final : k4FWCore::Consumer -#include -#include - -void printInStars(const Gaudi::Algorithm* thisAlg, const std::string& msg, const int lineWidth) { - thisAlg->debug() << fmt::format("{:*^{}}", "", lineWidth) << endmsg; - thisAlg->debug() << fmt::format("{:*^{}}", msg, lineWidth) << endmsg; - thisAlg->debug() << fmt::format("{:*^{}}", "", lineWidth) << endmsg; -} diff --git a/Tracking/components/printStars.h b/Tracking/components/printStars.h deleted file mode 100644 index fc621318..00000000 --- a/Tracking/components/printStars.h +++ /dev/null @@ -1,4 +0,0 @@ -#include -#include - -void printInStars(const Gaudi::Algorithm* thisAlg, const std::string& msg, const int lineWidth); diff --git a/Tracking/include/utils.hpp b/Tracking/include/utils.hpp index 9a731685..d6c3743b 100644 --- a/Tracking/include/utils.hpp +++ b/Tracking/include/utils.hpp @@ -20,6 +20,8 @@ #include "extension/MutableTrack.h" #include "extension/TrackCollection.h" +#include + //=== Others === #include #include @@ -189,4 +191,9 @@ int getHypotesisCharge(int pdg); TMatrixDSym computeTrackStateCovMatrix(TVectorD stateTrack, TVectorD params, TVector3 referencePoint, double timeError, TMatrixDSym statecovMatrix); -#endif // UTILS_HPP \ No newline at end of file +/** + * @brief Print the message in a block of stars ('*') at level DEBUG + */ +void printInStars(const Gaudi::Algorithm* thisAlg, const std::string& msg, const int lineWidth); + +#endif // UTILS_HPP diff --git a/Tracking/src/utils.cpp b/Tracking/src/utils.cpp index fe17829e..9644173d 100644 --- a/Tracking/src/utils.cpp +++ b/Tracking/src/utils.cpp @@ -1,5 +1,7 @@ #include "utils.hpp" +#include + dd4hep::rec::LayeredCalorimeterData* getExtension(unsigned int includeFlag, unsigned int excludeFlag) { dd4hep::rec::LayeredCalorimeterData* theExtension = 0; @@ -353,3 +355,9 @@ TMatrixDSym computeTrackStateCovMatrix(TVectorD stateTrack, TVectorD params, TVe return covarianceTrackState; } + +void printInStars(const Gaudi::Algorithm* thisAlg, const std::string& msg, const int lineWidth) { + thisAlg->debug() << fmt::format("{:*^{}}", "", lineWidth) << endmsg; + thisAlg->debug() << fmt::format("{:*^{}}", msg, lineWidth) << endmsg; + thisAlg->debug() << fmt::format("{:*^{}}", "", lineWidth) << endmsg; +}