diff --git a/CMakeLists.txt b/CMakeLists.txt
index b37a529..31886c8 100644
--- a/CMakeLists.txt
+++ b/CMakeLists.txt
@@ -35,6 +35,17 @@ if(WITH_ATLTileCalTB_PulseOutput)
   add_compile_definitions(ATLTileCalTB_PulseOutput)
 endif()
 
+#----------------------------------------------------------------------------
+#Option to enable leakage particle spectrum analysis
+#
+option(WITH_LEAKAGEANALYSIS "enable SpectrumAnalysis on leakage" OFF)
+if(WITH_LEAKAGEANALYSIS)
+    add_compile_definitions(ATLTileCalTB_LEAKANALYSIS)
+    message(STATUS "WITH_LEAKAGEANALYSIS = ON : building with ATLTileCalTB_LEAKANALYSIS compiler definition.")
+else()
+    message(STATUS "WITH_LEAKAGEANALYSIS = OFF : building without ATLTileCalTB_LEAKANALYSIS compiler definition.")
+endif()
+
 #----------------------------------------------------------------------------
 # Option to disable noise
 #
diff --git a/README.md b/README.md
index fe32ce1..62a4477 100644
--- a/README.md
+++ b/README.md
@@ -287,7 +287,8 @@ Custom options:
 -  `WITH_ATLTileCalTB_NoNoise`: if set to `ON`, the simulation will not put electronic noise on the
    signal (per cell) and disable the 2 sigma noise cut. Only relevant for noise calibration.
 -  `WITH_GEANT4_UIVIS`: if set to `ON` (default), build with UI and visualization drivers.
--   `G4_USE_FLUKA`: if set to `ON` build against the Fluka.Cern interface (default `OFF`).
+-  `G4_USE_FLUKA`: if set to `ON` build against the Fluka.Cern interface (default `OFF`).
+-  `WITH_LEAKAGEANALYSIS`: if set to `ON` build with leakage spectrum analyzer (default `OFF`).
 
 Relevant built-in options:
 -  `CMAKE_BUILD_TYPE`: set to `Debug` for debugging and to `Release` for production (faster).
diff --git a/include/SpectrumAnalyzer.hh b/include/SpectrumAnalyzer.hh
new file mode 100644
index 0000000..4bef476
--- /dev/null
+++ b/include/SpectrumAnalyzer.hh
@@ -0,0 +1,90 @@
+//**************************************************
+// \file SpectrumAnalyzer.hh
+// \brief: Declaration of SpectrumAnalyzer class
+// \author: Lorenzo Pezzotti (CERN EP-SFT-sim)
+//          @lopezzot
+// \start date: 28 August 2023
+//**************************************************
+
+// A portable Geant4-based particle spectrum analyzer
+// to be used within a Geant4 simulation without affecting it.
+// Instead of coding it in the simulation, create a singleton
+// and manage its usage with (#ifdef) compiler definition.
+
+#ifndef SpectrumAnalyzer_h
+#  define SpectrumAnalyzer_h
+
+// Includers from Geant4
+//
+#  include "G4Step.hh"
+#  include "G4ThreadLocalSingleton.hh"
+
+// Includers from C++
+//
+#  include <functional>
+
+class SpectrumAnalyzer
+{
+    friend class G4ThreadLocalSingleton<SpectrumAnalyzer>;
+
+  public:
+    // Return pointer to class instance
+    static SpectrumAnalyzer* GetInstance()
+    {
+      static G4ThreadLocalSingleton<SpectrumAnalyzer> instance{};
+      return instance.Instance();
+    }
+
+    // Methods
+    //
+    // Run-wise methods
+    void CreateNtupleAndScorer(const G4String scName = "te");
+    inline void ClearNtupleID() { ntupleID = 99; }
+    // Event-wise methods
+    inline void ClearEventFields()
+    {
+      neutronScore = 0.;
+      protonScore = 0., pionScore = 0., gammaScore = 0., electronScore = 0., othersScore = 0.;
+    }
+    void FillEventFields() const;
+    // Step-wise methods
+    void Analyze(const G4Step* step);
+
+  private:
+    // Members
+    //
+    // Run-wise members
+    G4int ntupleID;
+    std::function<G4double(const G4Step* step)> scorer;
+    G4String scorerName{};
+    // Event-wise members
+    G4double neutronScore;
+    G4double protonScore;
+    G4double pionScore;
+    G4double gammaScore;
+    G4double electronScore;
+    G4double othersScore;
+
+    // Scoring quantities
+    inline static G4double GetMomentum(const G4Step* step)
+    {
+      return step->GetTrack()->GetMomentum().mag();
+    };
+    inline static G4double GetKE(const G4Step* step)
+    {
+      return step->GetTrack()->GetKineticEnergy();
+    };
+    inline static G4double GetTE(const G4Step* step) { return step->GetTrack()->GetTotalEnergy(); };
+
+  private:
+    // Private constructor
+    SpectrumAnalyzer() = default;
+
+  public:
+    SpectrumAnalyzer(SpectrumAnalyzer const&) = delete;
+    void operator=(SpectrumAnalyzer const&) = delete;
+};
+
+#endif  // SpectrumAnalyzer_h
+
+//**************************************************
diff --git a/src/ATLTileCalTBEventAction.cc b/src/ATLTileCalTBEventAction.cc
index f0cad78..727cbdc 100644
--- a/src/ATLTileCalTBEventAction.cc
+++ b/src/ATLTileCalTBEventAction.cc
@@ -13,6 +13,9 @@
 #include "ATLTileCalTBGeometry.hh"
 #include "ATLTileCalTBConstants.hh"
 #include "ATLTileCalTBPrimaryGenAction.hh"
+#ifdef ATLTileCalTB_LEAKANALYSIS
+#include "SpectrumAnalyzer.hh"
+#endif
 
 //Includers from Geant4
 //
@@ -67,7 +70,10 @@ void ATLTileCalTBEventAction::BeginOfEventAction([[maybe_unused]] const G4Event*
     pulse_event_path = std::filesystem::path("ATLTileCalTBpulse_Run" + std::to_string(runNumber) + "/Ev" + std::to_string(eventNumber));
     std::filesystem::create_directory(pulse_event_path);
     #endif
-
+    
+    #ifdef ATLTileCalTB_LEAKANALYSIS
+    SpectrumAnalyzer::GetInstance()->ClearEventFields();
+    #endif
 }
 
 //GetHitsCollection method()
@@ -200,7 +206,10 @@ void ATLTileCalTBEventAction::EndOfEventAction( const G4Event* event ) {
     analysisManager->FillNtupleFColumn(7, fPrimaryGenAction->GetParticlenGun()->GetParticleEnergy());
 
     analysisManager->AddNtupleRow();
-
+    
+    #ifdef ATLTileCalTB_LEAKANALYSIS
+    SpectrumAnalyzer::GetInstance()->FillEventFields();
+    #endif
 } 
 
 //**************************************************
diff --git a/src/ATLTileCalTBRunAction.cc b/src/ATLTileCalTBRunAction.cc
index 177ddcc..061dc69 100644
--- a/src/ATLTileCalTBRunAction.cc
+++ b/src/ATLTileCalTBRunAction.cc
@@ -11,6 +11,9 @@
 //
 #include "ATLTileCalTBRunAction.hh"
 #include "ATLTileCalTBEventAction.hh"
+#ifdef ATLTileCalTB_LEAKANALYSIS
+#include "SpectrumAnalyzer.hh"
+#endif
 
 //Includers from Geant4
 //
@@ -73,7 +76,10 @@ ATLTileCalTBRunAction::ATLTileCalTBRunAction( ATLTileCalTBEventAction* eventActi
     analysisManager->CreateNtupleIColumn("PDGID");
     analysisManager->CreateNtupleFColumn("EBeam");
     analysisManager->FinishNtuple();
-
+    
+    #ifdef ATLTileCalTB_LEAKANALYSIS
+    SpectrumAnalyzer::GetInstance()->CreateNtupleAndScorer("ke");
+    #endif
 }
 
 ATLTileCalTBRunAction::~ATLTileCalTBRunAction() {
@@ -109,7 +115,10 @@ void ATLTileCalTBRunAction::EndOfRunAction(const G4Run* /*run*/) {
     auto analysisManager = G4AnalysisManager::Instance();
     analysisManager->Write();
     analysisManager->CloseFile();
-
+    
+    #ifdef ATLTileCalTB_LEAKANALYSIS
+    SpectrumAnalyzer::GetInstance()->ClearNtupleID();
+    #endif
 }
 
 //**************************************************
diff --git a/src/ATLTileCalTBStepAction.cc b/src/ATLTileCalTBStepAction.cc
index 501e443..0c44a42 100644
--- a/src/ATLTileCalTBStepAction.cc
+++ b/src/ATLTileCalTBStepAction.cc
@@ -10,6 +10,9 @@
 //Includers from project files
 //
 #include "ATLTileCalTBStepAction.hh"
+#ifdef ATLTileCalTB_LEAKANALYSIS
+#include "SpectrumAnalyzer.hh"
+#endif
 
 //Constructor and de-constructor
 //
@@ -27,6 +30,9 @@ void ATLTileCalTBStepAction::UserSteppingAction( const G4Step* aStep ) {
     //
     if ( !aStep->GetTrack()->GetNextVolume() ){
         fEventAction->Add( 0, aStep->GetTrack()->GetKineticEnergy() ); 
+        #ifdef ATLTileCalTB_LEAKANALYSIS
+        SpectrumAnalyzer::GetInstance()->Analyze(aStep);
+        #endif
     }
 
     if ( aStep->GetTrack()->GetTouchableHandle()->GetVolume()->GetName() != "CALO::CALO" ||
diff --git a/src/SpectrumAnalyzer.cc b/src/SpectrumAnalyzer.cc
new file mode 100644
index 0000000..546d6ca
--- /dev/null
+++ b/src/SpectrumAnalyzer.cc
@@ -0,0 +1,97 @@
+//**************************************************
+// \file SpectrumAnalyzer.cc
+// \brief: Definition of SpectrumAnalyzer class
+// \author: Lorenzo Pezzotti (CERN EP-SFT-sim)
+//          @lopezzot
+// \start date: 28 August 2023
+//**************************************************
+
+// Includers from project files
+//
+#include "SpectrumAnalyzer.hh"
+
+// Includers from Geant4
+//
+#include "G4Version.hh"
+#if G4VERSION_NUMBER < 1100
+#  include "g4root.hh"
+#else
+#  include "G4AnalysisManager.hh"
+#endif
+#include "G4Step.hh"
+
+// #define DEBUG
+
+void SpectrumAnalyzer::CreateNtupleAndScorer(const G4String scName)
+{
+  auto AM = G4AnalysisManager::Instance();
+
+  ntupleID = AM->CreateNtuple("Spectrum", "Spectrum");
+  AM->CreateNtupleDColumn("neutronScore");
+  AM->CreateNtupleDColumn("protonScore");
+  AM->CreateNtupleDColumn("pionScore");
+  AM->CreateNtupleDColumn("gammaScore");
+  AM->CreateNtupleDColumn("electronScore");
+  AM->CreateNtupleDColumn("othersScore");
+  AM->FinishNtuple();
+
+  // Define scorer type
+  scorerName = scName;
+  if (scorerName == "te") {
+    scorer = GetTE;
+  }
+  if (scorerName == "momentum") {
+    scorer = GetMomentum;
+  }
+  else if (scorerName == "ke") {
+    scorer = GetKE;
+  }
+  else {
+    scorer = GetTE;
+  }  // default case
+}
+
+void SpectrumAnalyzer::FillEventFields() const
+{
+  auto AM = G4AnalysisManager::Instance();
+  AM->FillNtupleDColumn(ntupleID, 0, neutronScore);
+  AM->FillNtupleDColumn(ntupleID, 1, protonScore);
+  AM->FillNtupleDColumn(ntupleID, 2, pionScore);
+  AM->FillNtupleDColumn(ntupleID, 3, gammaScore);
+  AM->FillNtupleDColumn(ntupleID, 4, electronScore);
+  AM->FillNtupleDColumn(ntupleID, 5, othersScore);
+  AM->AddNtupleRow(ntupleID);
+}
+
+void SpectrumAnalyzer::Analyze(const G4Step* step)
+{
+  auto PDGID = step->GetTrack()->GetParticleDefinition()->GetPDGEncoding();
+  auto val = scorer(step);
+  if (PDGID == 2112 || PDGID == -2112) {
+    neutronScore += val;
+  }
+  else if (PDGID == 2212 || PDGID == 2212) {
+    protonScore += val;
+  }
+  else if (PDGID == 211 || PDGID == 211) {
+    pionScore += val;
+  }
+  else if (PDGID == 22) {
+    gammaScore += val;
+  }
+  else if (PDGID == -11 || PDGID == 11) {
+    electronScore += val;
+  }
+  else {
+    othersScore += val;
+  }
+
+#ifdef DEBUG
+  G4cout << "-->SpectrumAnalyzer::Analyze, scorer name " << scorerName << " " << PDGID << " "
+         << step->GetTrack()->GetParticleDefinition()->GetParticleName() << " Total Energy "
+         << GetTE(step) << " Momentum " << GetMomentum(step) << " Kinetic Energy " << GetKE(step)
+         << G4endl;
+#endif
+}
+
+//**************************************************