From b36198abdc05fb7768a530507e6f2ce5684ae0b8 Mon Sep 17 00:00:00 2001 From: Marc Granado Gonzalez Date: Wed, 17 Dec 2025 15:37:34 +0100 Subject: [PATCH 01/12] First changes with Copilot for spatial resconstruction uncertainty --- source/digits_hits/include/GateDigi.hh | 12 +++++ source/digits_hits/src/GateDigi.cc | 3 ++ .../digits_hits/src/GateSpatialResolution.cc | 7 ++- source/general/include/GateRootDefs.hh | 9 ++++ source/general/src/GateRootDefs.cc | 49 +++++++++++++++++-- 5 files changed, 76 insertions(+), 4 deletions(-) diff --git a/source/digits_hits/include/GateDigi.hh b/source/digits_hits/include/GateDigi.hh index cf513cb44..95908f6df 100644 --- a/source/digits_hits/include/GateDigi.hh +++ b/source/digits_hits/include/GateDigi.hh @@ -173,6 +173,13 @@ public: inline void SetParentID(G4int parentID) { m_parentID = parentID; } inline G4int GetParentID() const { return m_parentID; } + inline void SetSpatialRes2DStdDevX(G4double v) { m_spatialRes2DStdDevX = v; } + inline G4double GetSpatialRes2DStdDevX() const { return m_spatialRes2DStdDevX; } + inline void SetSpatialRes2DStdDevY(G4double v) { m_spatialRes2DStdDevY = v; } + inline G4double GetSpatialRes2DStdDevY() const { return m_spatialRes2DStdDevY; } + inline void SetSpatialRes2DStdDevZ(G4double v) { m_spatialRes2DStdDevZ = v; } + inline G4double GetSpatialRes2DStdDevZ() const { return m_spatialRes2DStdDevZ; } + //AE inline G4double GetEnergyError() const { return m_energyError; } inline void SetEnergyError(G4double value) { m_energyError = value; } @@ -252,6 +259,11 @@ public: G4int m_nCrystalConv; //-------------------- + // Spatial resolution info set by spatial resolution digitizer when using 2D distributions + G4double m_spatialRes2DStdDevX; + G4double m_spatialRes2DStdDevY; + G4double m_spatialRes2DStdDevZ; + //! Pointer to the original crystal hit if known const void* m_mother; diff --git a/source/digits_hits/src/GateDigi.cc b/source/digits_hits/src/GateDigi.cc index d6bea9e2f..78a09ad5d 100644 --- a/source/digits_hits/src/GateDigi.cc +++ b/source/digits_hits/src/GateDigi.cc @@ -42,6 +42,9 @@ GateDigi::GateDigi(const void* itsMother): m_energyError(0.0), m_globalPosError(0.0), m_localPosError(0.0), + m_spatialRes2DStdDevX(0.0), + m_spatialRes2DStdDevY(0.0), + m_spatialRes2DStdDevZ(0.0), m_mother(itsMother) { } diff --git a/source/digits_hits/src/GateSpatialResolution.cc b/source/digits_hits/src/GateSpatialResolution.cc index f5e795ff0..a63ee9b1b 100644 --- a/source/digits_hits/src/GateSpatialResolution.cc +++ b/source/digits_hits/src/GateSpatialResolution.cc @@ -196,7 +196,7 @@ void GateSpatialResolution::Digitize(){ G4double Px = P.x(); G4double Py = P.y(); G4double Pz = P.z(); - G4double stddevX, stddevY, stddevZ; + G4double stddevX = 0., stddevY = 0., stddevZ = 0.; if (m_fwhmDistrib2D) { if (m_nameAxis.size() != 2) { @@ -240,6 +240,11 @@ void GateSpatialResolution::Digitize(){ G4double PyNew ; G4double PzNew ; +// store the computed stddevs into the digi for later ROOT output + m_outputDigi->SetSpatialRes2DStdDevX(stddevX); + m_outputDigi->SetSpatialRes2DStdDevY(stddevY); + m_outputDigi->SetSpatialRes2DStdDevZ(stddevZ); + if (m_IsConfined) diff --git a/source/general/include/GateRootDefs.hh b/source/general/include/GateRootDefs.hh index 747095c08..14727a544 100644 --- a/source/general/include/GateRootDefs.hh +++ b/source/general/include/GateRootDefs.hh @@ -338,6 +338,10 @@ class GateRootSingleBuffer Float_t energyIni; Int_t septalNb; //!< HDS : septal penetration + // Spatial resolution standard deviations (from 2D FWHM distribution) in mm + Float_t spatialRes2DStdDevX; + Float_t spatialRes2DStdDevY; + Float_t spatialRes2DStdDevZ; //@} }; @@ -434,6 +438,11 @@ class GateRootCoincBuffer Char_t comptonVolumeName2[40]; Char_t RayleighVolumeName2[40]; + // Spatial resolution standard deviations (averaged from both digis) in mm + Float_t spatialRes2DStdDevX; + Float_t spatialRes2DStdDevY; + Float_t spatialRes2DStdDevZ; + Float_t sinogramTheta; Float_t sinogramS; //@} diff --git a/source/general/src/GateRootDefs.cc b/source/general/src/GateRootDefs.cc index 140a8cfd3..6935e5e17 100644 --- a/source/general/src/GateRootDefs.cc +++ b/source/general/src/GateRootDefs.cc @@ -462,6 +462,11 @@ void GateRootSingleBuffer::Clear() // HDS : septal septalNb = 0; + + // initialize spatial resolution stddev fields + spatialRes2DStdDevX = 0.; + spatialRes2DStdDevY = 0.; + spatialRes2DStdDevZ = 0.; for ( d = 0 ; d < ROOT_VOLUMEIDSIZE ; ++d ) volumeID[d] = -1; @@ -490,6 +495,10 @@ void GateRootSingleBuffer::Clear() nCrystalCompt=0; nCrystalRayl=0; + spatialRes2DStdDevX = 0.; + spatialRes2DStdDevY = 0.; + spatialRes2DStdDevZ = 0.; + } @@ -551,6 +560,11 @@ void GateRootSingleBuffer::Fill(GateDigi* aDigi) } + // spatial resolution stddevs stored by digitizer (in internal units), convert to mm for output + spatialRes2DStdDevX = static_cast(aDigi->GetSpatialRes2DStdDevX()/mm); + spatialRes2DStdDevY = static_cast(aDigi->GetSpatialRes2DStdDevY()/mm); + spatialRes2DStdDevZ = static_cast(aDigi->GetSpatialRes2DStdDevZ()/mm); + aDigi->GetVolumeID().StoreDaughterIDs(volumeID,ROOT_VOLUMEIDSIZE); } @@ -630,8 +644,13 @@ void GateSingleTree::Init(GateRootSingleBuffer& buffer) } - //Initialized by default.TO DO: Mask option should be included or a flag - Branch("volumeID", (void *)buffer.volumeID,"volumeID[10]/I"); + // spatial resolution stddevs (2D distribution) - in mm + Branch("spatialRes2DStdDevX", &buffer.spatialRes2DStdDevX, "spatialRes2DStdDevX/F"); + Branch("spatialRes2DStdDevY", &buffer.spatialRes2DStdDevY, "spatialRes2DStdDevY/F"); + Branch("spatialRes2DStdDevZ", &buffer.spatialRes2DStdDevZ, "spatialRes2DStdDevZ/F"); + + //Initialized by default.TO DO: Mask option should be included or a flag + Branch("volumeID", (void *)buffer.volumeID,"volumeID[10]/I"); } @@ -682,7 +701,11 @@ void GateRootCoincBuffer::Clear() RayleighCrystal2 = -1; strcpy (comptonVolumeName2," "); strcpy (RayleighVolumeName2," "); -} + + // initialize spatial resolution stddev fields for coincidences + spatialRes2DStdDevX = 0.; + spatialRes2DStdDevY = 0.; + spatialRes2DStdDevZ = 0.; @@ -736,6 +759,21 @@ void GateRootCoincBuffer::Fill(GateCoincidenceDigi* aDigi) strcpy (comptonVolumeName2,((aDigi->GetDigi(1))->GetComptonVolumeName()).c_str()); strcpy (RayleighVolumeName2,((aDigi->GetDigi(1))->GetRayleighVolumeName()).c_str()); + // spatial resolution stddevs: average the two constituent digis and convert to mm for output + { + G4double sX1 = (aDigi->GetDigi(0))->GetSpatialRes2DStdDevX(); + G4double sX2 = (aDigi->GetDigi(1))->GetSpatialRes2DStdDevX(); + spatialRes2DStdDevX = static_cast(((sX1 + sX2) / 2.0) / mm); + + G4double sY1 = (aDigi->GetDigi(0))->GetSpatialRes2DStdDevY(); + G4double sY2 = (aDigi->GetDigi(1))->GetSpatialRes2DStdDevY(); + spatialRes2DStdDevY = static_cast(((sY1 + sY2) / 2.0) / mm); + + G4double sZ1 = (aDigi->GetDigi(0))->GetSpatialRes2DStdDevZ(); + G4double sZ2 = (aDigi->GetDigi(1))->GetSpatialRes2DStdDevZ(); + spatialRes2DStdDevZ = static_cast(((sZ1 + sZ2) / 2.0) / mm); + } + sinogramTheta = ComputeSinogramTheta(); sinogramS = ComputeSinogramS(); } @@ -864,6 +902,11 @@ void GateCoincTree::Init(GateRootCoincBuffer& buffer) if ( GateCoincidenceDigi::GetCoincidenceASCIIMask(19) ) Branch("sinogramS", &buffer.sinogramS,"sinogramS/F"); + // spatial resolution stddevs (averaged per coincidence) - in mm + Branch("spatialRes2DStdDevX", &buffer.spatialRes2DStdDevX, "spatialRes2DStdDevX/F"); + Branch("spatialRes2DStdDevY", &buffer.spatialRes2DStdDevY, "spatialRes2DStdDevY/F"); + Branch("spatialRes2DStdDevZ", &buffer.spatialRes2DStdDevZ, "spatialRes2DStdDevZ/F"); + if ( GateCoincidenceDigi::GetCoincidenceASCIIMask(20) ) Branch("comptVolName1", (void *)buffer.comptonVolumeName1,"comptVolName1/C"); if ( GateCoincidenceDigi::GetCoincidenceASCIIMask(20) ) From 1b361d3f2446861b57e5694dfb4e1e8b60b6c632 Mon Sep 17 00:00:00 2001 From: granadogmarc Date: Thu, 18 Dec 2025 14:53:26 +0100 Subject: [PATCH 02/12] New generalisation of the 2D distribution allowing the 3 axis to be changed with the same or different distributions chosing what are the variables to use to determine the position within the distribution --- CMakeLists.txt | 11 +++ .../include/GateSpatialResolution.hh | 17 +++- .../include/GateSpatialResolutionMessenger.hh | 3 + .../digits_hits/src/GateSpatialResolution.cc | 91 +++++++++---------- .../src/GateSpatialResolutionMessenger.cc | 59 +++++++++--- .../src/GateVirtualSegmentationSD.cc | 2 +- .../digits_hits/test/TestSpatialResStdDev.cc | 38 ++++++++ source/general/src/GateRootDefs.cc | 1 + 8 files changed, 160 insertions(+), 62 deletions(-) create mode 100644 source/digits_hits/test/TestSpatialResStdDev.cc diff --git a/CMakeLists.txt b/CMakeLists.txt index e009bc4b5..7f41e332f 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -300,6 +300,17 @@ IF(GATE_COMPILE_GATEDIGIT) INSTALL(TARGETS Convert_CCMod2PETCoinc DESTINATION bin) ENDIF(GATE_COMPILE_GATEDIGIT) +# Option to build small unit tests for Gate (requires Geant4/ROOT as configured above) +OPTION(GATE_BUILD_UNIT_TESTS "Build unit tests for Gate" OFF) +IF(GATE_BUILD_UNIT_TESTS) + ADD_EXECUTABLE(TestSpatialResStdDev ${PROJECT_SOURCE_DIR}/source/digits_hits/test/TestSpatialResStdDev.cc $) + TARGET_LINK_LIBRARIES(TestSpatialResStdDev GateLib) + INSTALL(TARGETS TestSpatialResStdDev DESTINATION bin) + IF(BUILD_TESTING) + add_test(NAME TestSpatialResStdDev COMMAND TestSpatialResStdDev) + ENDIF(BUILD_TESTING) +ENDIF(GATE_BUILD_UNIT_TESTS) + #========================================================= # We remove the warning option "shadow", because there are tons of # such warning related to clhep/g4 system of units. diff --git a/source/digits_hits/include/GateSpatialResolution.hh b/source/digits_hits/include/GateSpatialResolution.hh index 8fd30748a..7296ec953 100644 --- a/source/digits_hits/include/GateSpatialResolution.hh +++ b/source/digits_hits/include/GateSpatialResolution.hh @@ -51,7 +51,11 @@ public: GateVDistribution* GetFWHMydistrib() { return m_fwhmYdistrib; } GateVDistribution* GetFWHMzdistrib() { return m_fwhmZdistrib; } G4String GetNameAxis() { return m_nameAxis; } - GateVDistribution* GetFWHMDistrib2D() { return m_fwhmDistrib2D; } + GateVDistribution* GetFWHMDistrib2D() { return m_fwhmXDistrib2D ? m_fwhmXDistrib2D : (m_fwhmYDistrib2D? m_fwhmYDistrib2D : m_fwhmZDistrib2D); } + + GateVDistribution* GetFWHMXDistrib2D() { return m_fwhmXDistrib2D; } + GateVDistribution* GetFWHMYDistrib2D() { return m_fwhmYDistrib2D; } + GateVDistribution* GetFWHMZDistrib2D() { return m_fwhmZDistrib2D; } G4double GetFWHMx() { return m_fwhmX; } G4double GetFWHMy() { return m_fwhmY; } @@ -68,7 +72,12 @@ public: void SetNameAxis(const G4String& name) {m_nameAxis=name;} - void SetFWHMDistrib2D(GateVDistribution* dist) { m_fwhmDistrib2D= dist;} + // Backwards-compatible: sets all 3 axis 2D distributions to the same distribution + void SetFWHMDistrib2D(GateVDistribution* dist) { m_fwhmXDistrib2D = m_fwhmYDistrib2D = m_fwhmZDistrib2D = dist; } + + void SetFWHMXDistrib2D(GateVDistribution* dist) { m_fwhmXDistrib2D = dist; } + void SetFWHMYDistrib2D(GateVDistribution* dist) { m_fwhmYDistrib2D = dist; } + void SetFWHMZDistrib2D(GateVDistribution* dist) { m_fwhmZDistrib2D = dist; } void SetFWHMx(G4double val) { m_fwhmX = val; } void SetFWHMy(G4double val) { m_fwhmY = val; } @@ -101,7 +110,9 @@ protected: GateVDistribution* m_fwhmYdistrib; GateVDistribution* m_fwhmZdistrib; - GateVDistribution* m_fwhmDistrib2D; + GateVDistribution* m_fwhmXDistrib2D; + GateVDistribution* m_fwhmYDistrib2D; + GateVDistribution* m_fwhmZDistrib2D; G4String m_nameAxis; diff --git a/source/digits_hits/include/GateSpatialResolutionMessenger.hh b/source/digits_hits/include/GateSpatialResolutionMessenger.hh index 33d33213a..8672820ee 100644 --- a/source/digits_hits/include/GateSpatialResolutionMessenger.hh +++ b/source/digits_hits/include/GateSpatialResolutionMessenger.hh @@ -49,6 +49,9 @@ private: G4UIcmdWithAString *spresolutionZdistribCmd;// Command declaration for 1D Y-resolution distribution G4UIcmdWithAString *nameAxisCmd; G4UIcmdWithAString *spresolutionDistrib2DCmd; // Command declaration for 2D resolution distribution will be used with nameAxisCmd + G4UIcmdWithAString *spresolutionDistrib2DXCmd; // per-axis 2D distribution commands + G4UIcmdWithAString *spresolutionDistrib2DYCmd; + G4UIcmdWithAString *spresolutionDistrib2DZCmd; G4UIcmdWithABool* confineCmd; G4UIcmdWithABool* useTruncatedGaussianCmd; diff --git a/source/digits_hits/src/GateSpatialResolution.cc b/source/digits_hits/src/GateSpatialResolution.cc index a63ee9b1b..57f170f3f 100644 --- a/source/digits_hits/src/GateSpatialResolution.cc +++ b/source/digits_hits/src/GateSpatialResolution.cc @@ -55,11 +55,13 @@ GateSpatialResolution::GateSpatialResolution(GateSinglesDigitizer *digitizer, G4 m_fwhmX(0), m_fwhmY(0), m_fwhmZ(0), - m_fwhmXdistrib(0), - m_fwhmYdistrib(0), - m_fwhmZdistrib(0), - m_nameAxis("XY"), - m_fwhmDistrib2D(0), + m_fwhmXdistrib(0), + m_fwhmYdistrib(0), + m_fwhmZdistrib(0), + m_nameAxis("YZ"), + m_fwhmXDistrib2D(0), + m_fwhmYDistrib2D(0), + m_fwhmZDistrib2D(0), m_IsConfined(true), m_UseTruncatedGaussian(true), m_Navigator(0), @@ -83,28 +85,29 @@ GateSpatialResolution::~GateSpatialResolution() } void GateSpatialResolution::SetSpatialResolutionParameters() { // Check FWHM parameters - if (m_fwhm != 0 && (m_fwhmX != 0 || m_fwhmY != 0 || m_fwhmZ != 0 || m_fwhmDistrib2D != 0 || m_fwhmXdistrib !=0 || m_fwhmYdistrib !=0 || m_fwhmZdistrib !=0 )) { + if (m_fwhm != 0 && (m_fwhmX != 0 || m_fwhmY != 0 || m_fwhmZ != 0 || m_fwhmXDistrib2D != 0 || m_fwhmYDistrib2D != 0 || m_fwhmZDistrib2D != 0 || m_fwhmXdistrib !=0 || m_fwhmYdistrib !=0 || m_fwhmZdistrib !=0 )) { G4cout << "***ERROR*** Spatial Resolution is ambiguous: you can set a unique FWHM for all 3 axes OR set FWHM for X, Y, Z individually." << G4endl; abort(); } - if (m_fwhmDistrib2D) - { - if ((m_fwhmX != 0 || m_fwhmXdistrib !=0) && (m_nameAxis.find('X') != std::string::npos)) { - G4cout << "***ERROR*** Spatial Resolution is ambiguous: you can set FWHM for X OR set FWHM for X distribution." << G4endl; - abort(); - } - - if ((m_fwhmY != 0|| m_fwhmYdistrib !=0) && (m_nameAxis.find('Y') != std::string::npos)) { - G4cout << "***ERROR*** Spatial Resolution is ambiguous: you can set FWHM for Y OR set FWHM for Y distribution." << G4endl; - abort(); - } - - if ((m_fwhmZ != 0 || m_fwhmZdistrib !=0) && (m_nameAxis.find('Z') != std::string::npos)) { - G4cout << "***ERROR*** Spatial Resolution is ambiguous: you can set FWHM for Z OR set FWHM for Z distribution." << G4endl; - abort(); - } - } + if (m_fwhmXDistrib2D || m_fwhmYDistrib2D || m_fwhmZDistrib2D) + { + // If a per-axis 2D distribution is provided for an axis, it is + // ambiguous to also provide a scalar FWHM or a 1D distribution for + // the same axis. Check per-axis instead of using nameAxis combinatorics. + if (m_fwhmXDistrib2D && (m_fwhmX != 0 || m_fwhmXdistrib != 0)) { + G4cout << "***ERROR*** Spatial Resolution is ambiguous: you can set FWHM for X OR set FWHM distribution for X." << G4endl; + abort(); + } + if (m_fwhmYDistrib2D && (m_fwhmY != 0 || m_fwhmYdistrib != 0)) { + G4cout << "***ERROR*** Spatial Resolution is ambiguous: you can set FWHM for Y OR set FWHM distribution for Y." << G4endl; + abort(); + } + if (m_fwhmZDistrib2D && (m_fwhmZ != 0 || m_fwhmZdistrib != 0)) { + G4cout << "***ERROR*** Spatial Resolution is ambiguous: you can set FWHM for Z OR set FWHM distribution for Z." << G4endl; + abort(); + } + } if (m_fwhmX != 0 && m_fwhmXdistrib !=0){ G4cout << "***ERROR*** Spatial Resolution is ambiguous: you can set FWHM for X OR set FWHM for Z distribution." << G4endl; @@ -198,28 +201,18 @@ void GateSpatialResolution::Digitize(){ G4double Pz = P.z(); G4double stddevX = 0., stddevY = 0., stddevZ = 0.; - if (m_fwhmDistrib2D) { - if (m_nameAxis.size() != 2) { - GateError(" *** ERROR*** GateSpatialResolution::Digitize. " - "Attempt to use fwhmDistrib2D but the length of the named axis is not 2!\n"); - } else { - if (m_nameAxis == "XY") { - stddevX = m_fwhmDistrib2D->Value2D(P.x() * mm, P.y() * mm); - stddevY = m_fwhmDistrib2D->Value2D(P.x() * mm, P.y() * mm); - if (fwhmZ) stddevZ = fwhmZ / GateConstants::fwhm_to_sigma; - } else if (m_nameAxis == "XZ") { - stddevX = m_fwhmDistrib2D->Value2D(P.x() * mm, P.z() * mm); - stddevZ = m_fwhmDistrib2D->Value2D(P.x() * mm, P.z() * mm); - if (fwhmY) stddevY = fwhmY / GateConstants::fwhm_to_sigma; - } else if (m_nameAxis == "YZ") { - stddevY = m_fwhmDistrib2D->Value2D(P.y() * mm, P.z() * mm); - stddevZ = m_fwhmDistrib2D->Value2D(P.y() * mm, P.z() * mm); - if (fwhmX) stddevX = fwhmX / GateConstants::fwhm_to_sigma; - } else { - GateError(" *** ERROR*** GateSpatialResolution::Digitize. " - "Unrecognized axis configuration: " + m_nameAxis + "\n"); - } - } + // Use the configured axis pair (m_nameAxis) to evaluate Value2D. + // Allowed pairs for PET: "XZ" or "YZ" (default "YZ"). + if (m_fwhmXDistrib2D || m_fwhmYDistrib2D || m_fwhmZDistrib2D) { + if (m_nameAxis == "XZ") { + if (m_fwhmXDistrib2D) stddevX = m_fwhmXDistrib2D->Value2D(P.x() * mm, P.z() * mm); + if (m_fwhmYDistrib2D) stddevY = m_fwhmYDistrib2D->Value2D(P.x() * mm, P.z() * mm); + if (m_fwhmZDistrib2D) stddevZ = m_fwhmZDistrib2D->Value2D(P.x() * mm, P.z() * mm); + } else { // YZ + if (m_fwhmXDistrib2D) stddevX = m_fwhmXDistrib2D->Value2D(P.y() * mm, P.z() * mm); + if (m_fwhmYDistrib2D) stddevY = m_fwhmYDistrib2D->Value2D(P.y() * mm, P.z() * mm); + if (m_fwhmZDistrib2D) stddevZ = m_fwhmZDistrib2D->Value2D(P.y() * mm, P.z() * mm); + } } else { @@ -345,7 +338,13 @@ void GateSpatialResolution::Digitize(){ if (nVerboseLevel>1) G4cout << "[GateSpatialResolution::Digitize]: input digi collection is null -> nothing to do\n\n"; return; - } + // Ensure the chosen axis configuration is allowed for PET scanners + if (!(m_nameAxis == "XZ" || m_nameAxis == "YZ")) { + G4cout << "***ERROR*** GateSpatialResolution::SetSpatialResolutionParameters: " + "Only 'XZ' and 'YZ' are allowed as nameAxis values for 2D spatial resolution distributions.\n"; + abort(); + } + } StoreDigiCollection(m_OutputDigiCollection); } diff --git a/source/digits_hits/src/GateSpatialResolutionMessenger.cc b/source/digits_hits/src/GateSpatialResolutionMessenger.cc index 38e2711d7..f3595e92d 100644 --- a/source/digits_hits/src/GateSpatialResolutionMessenger.cc +++ b/source/digits_hits/src/GateSpatialResolutionMessenger.cc @@ -59,13 +59,23 @@ GateSpatialResolutionMessenger::GateSpatialResolutionMessenger (GateSpatialResol cmdName = GetDirectoryName() + "fwhmDistrib2D"; spresolutionDistrib2DCmd = new G4UIcmdWithAString(cmdName,this); - spresolutionDistrib2DCmd->SetGuidance("Set the distribution 2D of spatial resolution in position for gaussian spblurring the name of the axis must be defined. default=XY"); + spresolutionDistrib2DCmd->SetGuidance("Set the distribution 2D of spatial resolution in position for gaussian blurring. The name of the axis must be defined. default=XY. This legacy command sets the same distribution for X, Y and Z 2D distributors."); + cmdName = GetDirectoryName() + "fwhmDistrib2DX"; + spresolutionDistrib2DXCmd = new G4UIcmdWithAString(cmdName,this); + spresolutionDistrib2DXCmd->SetGuidance("Set the 2D distribution for X axis (expects a 2D distribution object)"); + cmdName = GetDirectoryName() + "fwhmDistrib2DY"; + spresolutionDistrib2DYCmd = new G4UIcmdWithAString(cmdName,this); + spresolutionDistrib2DYCmd->SetGuidance("Set the 2D distribution for Y axis (expects a 2D distribution object)"); + cmdName = GetDirectoryName() + "fwhmDistrib2DZ"; + spresolutionDistrib2DZCmd = new G4UIcmdWithAString(cmdName,this); + spresolutionDistrib2DZCmd->SetGuidance("Set the 2D distribution for Z axis (expects a 2D distribution object)"); - cmdName = GetDirectoryName()+"nameAxis"; - nameAxisCmd = new G4UIcmdWithAString(cmdName,this); - nameAxisCmd ->SetGuidance("Provide the name of the axis that will follow the spatial resolution distribution, X, Y, Z, XY, XZ, or YZ"); - nameAxisCmd ->SetCandidates("X Y Z XY XZ YZ"); + + cmdName = GetDirectoryName()+"nameAxis"; + nameAxisCmd = new G4UIcmdWithAString(cmdName,this); + nameAxisCmd ->SetGuidance("Provide the coordinate pair used by 2D distributions: 'XZ' or 'YZ'. Default is 'YZ'."); + nameAxisCmd ->SetCandidates("XZ YZ"); cmdName = GetDirectoryName() + "confineInsideOfSmallestElement"; @@ -93,6 +103,10 @@ GateSpatialResolutionMessenger::~GateSpatialResolutionMessenger() delete spresolutionYdistribCmd; delete spresolutionZdistribCmd; + delete spresolutionDistrib2DXCmd; + delete spresolutionDistrib2DYCmd; + delete spresolutionDistrib2DZCmd; + delete nameAxisCmd; delete spresolutionDistrib2DCmd; @@ -123,13 +137,34 @@ void GateSpatialResolutionMessenger::SetNewValue(G4UIcommand * aCommand,G4String // Handle command for 2D-distribution resolution if (aCommand == nameAxisCmd) - { - m_SpatialResolution->SetNameAxis(newValue); - } - else if (aCommand == spresolutionDistrib2DCmd) - {GateVDistribution* distrib = (GateVDistribution*)GateDistributionListManager::GetInstance()->FindElementByBaseName(newValue); - if (distrib) m_SpatialResolution->SetFWHMDistrib2D(distrib); - } + { + // Only accept the PET-relevant options + if (newValue == "XZ" || newValue == "YZ") { + m_SpatialResolution->SetNameAxis(newValue); + } else { + G4cout << "***ERROR*** GateSpatialResolution::SetNewValue: nameAxis must be 'XZ' or 'YZ'\n"<< G4endl; + } + } + else if (aCommand == spresolutionDistrib2DCmd) + { + GateVDistribution* distrib = (GateVDistribution*)GateDistributionListManager::GetInstance()->FindElementByBaseName(newValue); + if (distrib) m_SpatialResolution->SetFWHMDistrib2D(distrib); + } + else if (aCommand == spresolutionDistrib2DXCmd) + { + GateVDistribution* distrib = (GateVDistribution*)GateDistributionListManager::GetInstance()->FindElementByBaseName(newValue); + if (distrib) m_SpatialResolution->SetFWHMXDistrib2D(distrib); + } + else if (aCommand == spresolutionDistrib2DYCmd) + { + GateVDistribution* distrib = (GateVDistribution*)GateDistributionListManager::GetInstance()->FindElementByBaseName(newValue); + if (distrib) m_SpatialResolution->SetFWHMYDistrib2D(distrib); + } + else if (aCommand == spresolutionDistrib2DZCmd) + { + GateVDistribution* distrib = (GateVDistribution*)GateDistributionListManager::GetInstance()->FindElementByBaseName(newValue); + if (distrib) m_SpatialResolution->SetFWHMZDistrib2D(distrib); + } else if ( aCommand==spresolutionXCmd ) { m_SpatialResolution->SetFWHMx(spresolutionXCmd->GetNewDoubleValue(newValue)); } diff --git a/source/digits_hits/src/GateVirtualSegmentationSD.cc b/source/digits_hits/src/GateVirtualSegmentationSD.cc index 73107e648..125be63c5 100644 --- a/source/digits_hits/src/GateVirtualSegmentationSD.cc +++ b/source/digits_hits/src/GateVirtualSegmentationSD.cc @@ -426,7 +426,7 @@ void GateVirtualSegmentationSD::SetParameters() - if(digi_SpatialResolution->GetFWHMxdistrib()||digi_SpatialResolution->GetFWHMydistrib()||digi_SpatialResolution->GetFWHMDistrib2D()) + if(digi_SpatialResolution->GetFWHMxdistrib()||digi_SpatialResolution->GetFWHMydistrib()||digi_SpatialResolution->GetFWHMXDistrib2D()||digi_SpatialResolution->GetFWHMYDistrib2D()||digi_SpatialResolution->GetFWHMZDistrib2D()) { GateError("***ERROR*** No value of the target pitch has been provided and no value can be obtained from the spatial resolution distribution. /n Please provide a value for the pitch that is at least half of the minimum value of the distribution. "); } diff --git a/source/digits_hits/test/TestSpatialResStdDev.cc b/source/digits_hits/test/TestSpatialResStdDev.cc new file mode 100644 index 000000000..220cd958b --- /dev/null +++ b/source/digits_hits/test/TestSpatialResStdDev.cc @@ -0,0 +1,38 @@ +/* + Simple unit test to verify GateDigi stores spatial-resolution stddevs + This executable is optional and requires a Geant4-enabled build. +*/ +#include "GateDigi.hh" +#include +#include + +int main() +{ + GateDigi digi; + + const double x = 1.2345; + const double y = 2.3456; + const double z = 3.4567; + + digi.SetSpatialRes2DStdDevX(x); + digi.SetSpatialRes2DStdDevY(y); + digi.SetSpatialRes2DStdDevZ(z); + + const double eps = 1e-12; + + if (std::fabs(digi.GetSpatialRes2DStdDevX() - x) > eps) { + std::cerr << "Get/Set mismatch for StdDevX: got " << digi.GetSpatialRes2DStdDevX() << " expected " << x << "\n"; + return 1; + } + if (std::fabs(digi.GetSpatialRes2DStdDevY() - y) > eps) { + std::cerr << "Get/Set mismatch for StdDevY: got " << digi.GetSpatialRes2DStdDevY() << " expected " << y << "\n"; + return 2; + } + if (std::fabs(digi.GetSpatialRes2DStdDevZ() - z) > eps) { + std::cerr << "Get/Set mismatch for StdDevZ: got " << digi.GetSpatialRes2DStdDevZ() << " expected " << z << "\n"; + return 3; + } + + std::cout << "TestSpatialResStdDev: OK\n"; + return 0; +} diff --git a/source/general/src/GateRootDefs.cc b/source/general/src/GateRootDefs.cc index 6935e5e17..bfc90339d 100644 --- a/source/general/src/GateRootDefs.cc +++ b/source/general/src/GateRootDefs.cc @@ -706,6 +706,7 @@ void GateRootCoincBuffer::Clear() spatialRes2DStdDevX = 0.; spatialRes2DStdDevY = 0.; spatialRes2DStdDevZ = 0.; +} From 61a76fb2ffd13fc613e7851d3a55fb5d014d19fa Mon Sep 17 00:00:00 2001 From: granadogmarc Date: Tue, 27 Jan 2026 16:20:58 +0100 Subject: [PATCH 03/12] Adding authorship for Spatial Resolution Changes --- AUTHORS | 1 + source/digits_hits/src/GateSpatialResolution.cc | 11 +++++++++-- 2 files changed, 10 insertions(+), 2 deletions(-) diff --git a/AUTHORS b/AUTHORS index 102b936ab..87495aa9a 100644 --- a/AUTHORS +++ b/AUTHORS @@ -25,6 +25,7 @@ Garcia Marie-Paule Garpebring Anders Germano Russo Gil Alex Vergara +Granado-Gonzalez Marc Grevillot Loic Groiselle Corinne Gueth Pierre diff --git a/source/digits_hits/src/GateSpatialResolution.cc b/source/digits_hits/src/GateSpatialResolution.cc index 57f170f3f..4b3d2e544 100644 --- a/source/digits_hits/src/GateSpatialResolution.cc +++ b/source/digits_hits/src/GateSpatialResolution.cc @@ -23,9 +23,16 @@ This blurring has been validated up to a given FWHM of 10mm. At higher FWHM, the number of "relocated" digis is no longer negligible. The blurring effect is then so compensated that resolution will improve compared to lower values of FWHM. -modified by Radia Oudihat 06/2024 - Added support for 1D and 2D FWHM distributions for X and Y, and applied Gaussian blurring. + Added support for 1D FWHM distribution for X and Y, and applied Gaussian blurring. Implemented logic to determine standard deviations (stddevX, stddevY) based on defined 1D and 2D FWHM distributions for the X and Y axes. -*/ +-modified by Marc Granado-Gonzalez 2025 + - Added 2D FWHM distribution for X, Y and Z, and applied Gaussian blurring. + Implemented logic to determine standard deviations (stddevX, stddevY, stddevZ) based on defined 2D FWHM distributions for the X, Y and Z axes. + - Added option to choose the axis pair for 2D distributions (nameAxis): "XZ" or "YZ" (default "YZ"). + - Added Truncated Gaussian option for confined and non-confined cases. + */ + + #include "GateSpatialResolution.hh" #include "GateSpatialResolutionMessenger.hh" #include "GateDigi.hh" From 323f19fd9ffdd223a2a06b95ab2900f498543467 Mon Sep 17 00:00:00 2001 From: granadogmarc Date: Thu, 29 Jan 2026 17:06:14 +0100 Subject: [PATCH 04/12] Changes for consistency in spatialResolutionDistribution2D --- .../include/GateSpatialResolution.hh | 18 ++++---- .../include/GateSpatialResolutionMessenger.hh | 7 ++-- .../digits_hits/src/GateSpatialResolution.cc | 26 ++++++------ .../src/GateSpatialResolutionMessenger.cc | 41 ++++++++----------- 4 files changed, 41 insertions(+), 51 deletions(-) diff --git a/source/digits_hits/include/GateSpatialResolution.hh b/source/digits_hits/include/GateSpatialResolution.hh index 7296ec953..549c19f44 100644 --- a/source/digits_hits/include/GateSpatialResolution.hh +++ b/source/digits_hits/include/GateSpatialResolution.hh @@ -47,9 +47,9 @@ public: //! These functions return the resolution in use. G4double GetFWHM() { return m_fwhm; } - GateVDistribution* GetFWHMxdistrib() { return m_fwhmXdistrib; } - GateVDistribution* GetFWHMydistrib() { return m_fwhmYdistrib; } - GateVDistribution* GetFWHMzdistrib() { return m_fwhmZdistrib; } + GateVDistribution* GetFWHMxdistrib() { return m_fwhmXDistrib; } + GateVDistribution* GetFWHMydistrib() { return m_fwhmYDistrib; } + GateVDistribution* GetFWHMzdistrib() { return m_fwhmZDistrib; } G4String GetNameAxis() { return m_nameAxis; } GateVDistribution* GetFWHMDistrib2D() { return m_fwhmXDistrib2D ? m_fwhmXDistrib2D : (m_fwhmYDistrib2D? m_fwhmYDistrib2D : m_fwhmZDistrib2D); } @@ -66,9 +66,9 @@ public: If you want a resolution of 10%, SetSpresolution(0.1) */ void SetFWHM(G4double val) { m_fwhm = val; } - void SetFWHMxdistrib(GateVDistribution* dist) { m_fwhmXdistrib= dist; } - void SetFWHMydistrib(GateVDistribution* dist) { m_fwhmYdistrib = dist; } - void SetFWHMzdistrib(GateVDistribution* dist) { m_fwhmZdistrib = dist; } + void SetFWHMxdistrib(GateVDistribution* dist) { m_fwhmXDistrib= dist; } + void SetFWHMydistrib(GateVDistribution* dist) { m_fwhmYDistrib = dist; } + void SetFWHMzdistrib(GateVDistribution* dist) { m_fwhmZDistrib = dist; } void SetNameAxis(const G4String& name) {m_nameAxis=name;} @@ -106,9 +106,9 @@ protected: G4double m_fwhmY; G4double m_fwhmZ; - GateVDistribution* m_fwhmXdistrib; - GateVDistribution* m_fwhmYdistrib; - GateVDistribution* m_fwhmZdistrib; + GateVDistribution* m_fwhmXDistrib; + GateVDistribution* m_fwhmYDistrib; + GateVDistribution* m_fwhmZDistrib; GateVDistribution* m_fwhmXDistrib2D; GateVDistribution* m_fwhmYDistrib2D; diff --git a/source/digits_hits/include/GateSpatialResolutionMessenger.hh b/source/digits_hits/include/GateSpatialResolutionMessenger.hh index 8672820ee..bfe534c3a 100644 --- a/source/digits_hits/include/GateSpatialResolutionMessenger.hh +++ b/source/digits_hits/include/GateSpatialResolutionMessenger.hh @@ -48,10 +48,9 @@ private: G4UIcmdWithAString *spresolutionYdistribCmd;// Command declaration for 1D Y-resolution distribution G4UIcmdWithAString *spresolutionZdistribCmd;// Command declaration for 1D Y-resolution distribution G4UIcmdWithAString *nameAxisCmd; - G4UIcmdWithAString *spresolutionDistrib2DCmd; // Command declaration for 2D resolution distribution will be used with nameAxisCmd - G4UIcmdWithAString *spresolutionDistrib2DXCmd; // per-axis 2D distribution commands - G4UIcmdWithAString *spresolutionDistrib2DYCmd; - G4UIcmdWithAString *spresolutionDistrib2DZCmd; + G4UIcmdWithAString *spresolutionXDistrib2DCmd; // per-axis 2D distribution commands will be used nameAxis to determine which axes to apply the distribution to + G4UIcmdWithAString *spresolutionYDistrib2DCmd; + G4UIcmdWithAString *spresolutionZDistrib2DCmd; G4UIcmdWithABool* confineCmd; G4UIcmdWithABool* useTruncatedGaussianCmd; diff --git a/source/digits_hits/src/GateSpatialResolution.cc b/source/digits_hits/src/GateSpatialResolution.cc index 4b3d2e544..6c0b3573e 100644 --- a/source/digits_hits/src/GateSpatialResolution.cc +++ b/source/digits_hits/src/GateSpatialResolution.cc @@ -62,9 +62,9 @@ GateSpatialResolution::GateSpatialResolution(GateSinglesDigitizer *digitizer, G4 m_fwhmX(0), m_fwhmY(0), m_fwhmZ(0), - m_fwhmXdistrib(0), - m_fwhmYdistrib(0), - m_fwhmZdistrib(0), + m_fwhmXDistrib(0), + m_fwhmYDistrib(0), + m_fwhmZDistrib(0), m_nameAxis("YZ"), m_fwhmXDistrib2D(0), m_fwhmYDistrib2D(0), @@ -92,7 +92,7 @@ GateSpatialResolution::~GateSpatialResolution() } void GateSpatialResolution::SetSpatialResolutionParameters() { // Check FWHM parameters - if (m_fwhm != 0 && (m_fwhmX != 0 || m_fwhmY != 0 || m_fwhmZ != 0 || m_fwhmXDistrib2D != 0 || m_fwhmYDistrib2D != 0 || m_fwhmZDistrib2D != 0 || m_fwhmXdistrib !=0 || m_fwhmYdistrib !=0 || m_fwhmZdistrib !=0 )) { + if (m_fwhm != 0 && (m_fwhmX != 0 || m_fwhmY != 0 || m_fwhmZ != 0 || m_fwhmXDistrib2D != 0 || m_fwhmYDistrib2D != 0 || m_fwhmZDistrib2D != 0 || m_fwhmXDistrib !=0 || m_fwhmYDistrib !=0 || m_fwhmZDistrib !=0 )) { G4cout << "***ERROR*** Spatial Resolution is ambiguous: you can set a unique FWHM for all 3 axes OR set FWHM for X, Y, Z individually." << G4endl; abort(); } @@ -102,29 +102,29 @@ void GateSpatialResolution::SetSpatialResolutionParameters() { // If a per-axis 2D distribution is provided for an axis, it is // ambiguous to also provide a scalar FWHM or a 1D distribution for // the same axis. Check per-axis instead of using nameAxis combinatorics. - if (m_fwhmXDistrib2D && (m_fwhmX != 0 || m_fwhmXdistrib != 0)) { + if (m_fwhmXDistrib2D && (m_fwhmX != 0 || m_fwhmXDistrib != 0)) { G4cout << "***ERROR*** Spatial Resolution is ambiguous: you can set FWHM for X OR set FWHM distribution for X." << G4endl; abort(); } - if (m_fwhmYDistrib2D && (m_fwhmY != 0 || m_fwhmYdistrib != 0)) { + if (m_fwhmYDistrib2D && (m_fwhmY != 0 || m_fwhmYDistrib != 0)) { G4cout << "***ERROR*** Spatial Resolution is ambiguous: you can set FWHM for Y OR set FWHM distribution for Y." << G4endl; abort(); } - if (m_fwhmZDistrib2D && (m_fwhmZ != 0 || m_fwhmZdistrib != 0)) { + if (m_fwhmZDistrib2D && (m_fwhmZ != 0 || m_fwhmZDistrib != 0)) { G4cout << "***ERROR*** Spatial Resolution is ambiguous: you can set FWHM for Z OR set FWHM distribution for Z." << G4endl; abort(); } } - if (m_fwhmX != 0 && m_fwhmXdistrib !=0){ + if (m_fwhmX != 0 && m_fwhmXDistrib !=0){ G4cout << "***ERROR*** Spatial Resolution is ambiguous: you can set FWHM for X OR set FWHM for Z distribution." << G4endl; abort(); } - if (m_fwhmY != 0 && m_fwhmYdistrib !=0){ + if (m_fwhmY != 0 && m_fwhmYDistrib !=0){ G4cout << "***ERROR*** Spatial Resolution is ambiguous: you can set FWHM for Y OR set FWHM for Z distribution." << G4endl; abort(); } - if (m_fwhmZ != 0 && m_fwhmZdistrib !=0){ + if (m_fwhmZ != 0 && m_fwhmZDistrib !=0){ G4cout << "***ERROR*** Spatial Resolution is ambiguous: you can set FWHM for Z OR set FWHM for Z distribution." << G4endl; abort(); } @@ -223,13 +223,13 @@ void GateSpatialResolution::Digitize(){ } else { - if (m_fwhmXdistrib) stddevX = m_fwhmXdistrib->Value(P.x() * mm); + if (m_fwhmXDistrib) stddevX = m_fwhmXDistrib->Value(P.x() * mm); else if (fwhmX) stddevX = fwhmX / GateConstants::fwhm_to_sigma; - if (m_fwhmYdistrib) stddevY = m_fwhmYdistrib->Value(P.y() * mm); + if (m_fwhmYDistrib) stddevY = m_fwhmYDistrib->Value(P.y() * mm); else if (fwhmY) stddevY = fwhmY / GateConstants::fwhm_to_sigma; - if (m_fwhmZdistrib) stddevZ = m_fwhmZdistrib->Value(P.z() * mm); + if (m_fwhmZDistrib) stddevZ = m_fwhmZDistrib->Value(P.z() * mm); else if (fwhmZ) stddevZ = fwhmZ / GateConstants::fwhm_to_sigma; } diff --git a/source/digits_hits/src/GateSpatialResolutionMessenger.cc b/source/digits_hits/src/GateSpatialResolutionMessenger.cc index f3595e92d..413b271e2 100644 --- a/source/digits_hits/src/GateSpatialResolutionMessenger.cc +++ b/source/digits_hits/src/GateSpatialResolutionMessenger.cc @@ -57,19 +57,15 @@ GateSpatialResolutionMessenger::GateSpatialResolutionMessenger (GateSpatialResol spresolutionZdistribCmd = new G4UIcmdWithAString(cmdName,this); spresolutionZdistribCmd->SetGuidance("Set the distribution resolution in position for gaussian spblurring"); - cmdName = GetDirectoryName() + "fwhmDistrib2D"; - spresolutionDistrib2DCmd = new G4UIcmdWithAString(cmdName,this); - spresolutionDistrib2DCmd->SetGuidance("Set the distribution 2D of spatial resolution in position for gaussian blurring. The name of the axis must be defined. default=XY. This legacy command sets the same distribution for X, Y and Z 2D distributors."); - - cmdName = GetDirectoryName() + "fwhmDistrib2DX"; - spresolutionDistrib2DXCmd = new G4UIcmdWithAString(cmdName,this); - spresolutionDistrib2DXCmd->SetGuidance("Set the 2D distribution for X axis (expects a 2D distribution object)"); - cmdName = GetDirectoryName() + "fwhmDistrib2DY"; - spresolutionDistrib2DYCmd = new G4UIcmdWithAString(cmdName,this); - spresolutionDistrib2DYCmd->SetGuidance("Set the 2D distribution for Y axis (expects a 2D distribution object)"); - cmdName = GetDirectoryName() + "fwhmDistrib2DZ"; - spresolutionDistrib2DZCmd = new G4UIcmdWithAString(cmdName,this); - spresolutionDistrib2DZCmd->SetGuidance("Set the 2D distribution for Z axis (expects a 2D distribution object)"); + cmdName = GetDirectoryName() + "fwhmXDistrib2D"; + spresolutionXDistrib2DCmd = new G4UIcmdWithAString(cmdName,this); + spresolutionXDistrib2DCmd->SetGuidance("Set the 2D distribution for X axis (expects a 2D distribution object)"); + cmdName = GetDirectoryName() + "fwhmYDistrib2D"; + spresolutionYDistrib2DCmd = new G4UIcmdWithAString(cmdName,this); + spresolutionYDistrib2DCmd->SetGuidance("Set the 2D distribution for Y axis (expects a 2D distribution object)"); + cmdName = GetDirectoryName() + "fwhmZDistrib2D"; + spresolutionZDistrib2DCmd = new G4UIcmdWithAString(cmdName,this); + spresolutionZDistrib2DCmd->SetGuidance("Set the 2D distribution for Z axis (expects a 2D distribution object)"); cmdName = GetDirectoryName()+"nameAxis"; @@ -103,12 +99,12 @@ GateSpatialResolutionMessenger::~GateSpatialResolutionMessenger() delete spresolutionYdistribCmd; delete spresolutionZdistribCmd; - delete spresolutionDistrib2DXCmd; - delete spresolutionDistrib2DYCmd; - delete spresolutionDistrib2DZCmd; + delete spresolutionXDistrib2DCmd; + delete spresolutionYDistrib2DCmd; + delete spresolutionZDistrib2DCmd; delete nameAxisCmd; - delete spresolutionDistrib2DCmd; + delete confineCmd; delete useTruncatedGaussianCmd; @@ -145,22 +141,17 @@ void GateSpatialResolutionMessenger::SetNewValue(G4UIcommand * aCommand,G4String G4cout << "***ERROR*** GateSpatialResolution::SetNewValue: nameAxis must be 'XZ' or 'YZ'\n"<< G4endl; } } - else if (aCommand == spresolutionDistrib2DCmd) - { - GateVDistribution* distrib = (GateVDistribution*)GateDistributionListManager::GetInstance()->FindElementByBaseName(newValue); - if (distrib) m_SpatialResolution->SetFWHMDistrib2D(distrib); - } - else if (aCommand == spresolutionDistrib2DXCmd) + else if (aCommand == spresolutionXDistrib2DCmd) { GateVDistribution* distrib = (GateVDistribution*)GateDistributionListManager::GetInstance()->FindElementByBaseName(newValue); if (distrib) m_SpatialResolution->SetFWHMXDistrib2D(distrib); } - else if (aCommand == spresolutionDistrib2DYCmd) + else if (aCommand == spresolutionYDistrib2DCmd) { GateVDistribution* distrib = (GateVDistribution*)GateDistributionListManager::GetInstance()->FindElementByBaseName(newValue); if (distrib) m_SpatialResolution->SetFWHMYDistrib2D(distrib); } - else if (aCommand == spresolutionDistrib2DZCmd) + else if (aCommand == spresolutionZDistrib2DCmd) { GateVDistribution* distrib = (GateVDistribution*)GateDistributionListManager::GetInstance()->FindElementByBaseName(newValue); if (distrib) m_SpatialResolution->SetFWHMZDistrib2D(distrib); From 8f4b7159677e89a441d08abe0ec3a4fcf508098d Mon Sep 17 00:00:00 2001 From: granadogmarc Date: Tue, 3 Feb 2026 17:22:43 +0100 Subject: [PATCH 05/12] Update digitizer_and_detector_modeling.rst Changes in documentation for spatial resolution 2D distribution --- docs/digitizer_and_detector_modeling.rst | 15 +++++++++++---- 1 file changed, 11 insertions(+), 4 deletions(-) diff --git a/docs/digitizer_and_detector_modeling.rst b/docs/digitizer_and_detector_modeling.rst index 2bcae1a4f..052ed47e1 100644 --- a/docs/digitizer_and_detector_modeling.rst +++ b/docs/digitizer_and_detector_modeling.rst @@ -526,7 +526,10 @@ BEWARE: This relocation procedure is validated only for the first group level of **Configuring Spatial Resolution with 1D and 2D Distributions**:: -This approach is particularly essential for monolithic crystal detectors, where factors like edge effects and interaction positions significantly may influence spatial resolution. +This approach is particularly essential for monolithic crystal detectors, where factors such as edge effects and interaction positions may significantly influence spatial resolution. + +In order to map the 2D distribution with the crystal, the command "nameAxis" is used. The value of the axis provided will be the ones attributed to the 2 dimensions (columns and rows) of the distribution file. Assuming the crystal will be placed along the Z direction the options are "XZ" for a ring starting on top, and "YZ" for a ring starting on the sides. The default value for "nameAxis" is "YZ". + Here is an example of how to configure this in a macro file: **Example for 2D distribution**:: @@ -534,15 +537,18 @@ Here is an example of how to configure this in a macro file: /gate/distributions/name my_distrib2D /gate/distributions/insert File - /gate/distributions/my_distrib2D/setFileName Lut(X,Y).txt + /gate/distributions/my_distrib2D/setFileName Lut_XY.txt /gate/distributions/my_distrib2D/readMatrix2d /gate/digitizerMgr/crystalUnit/SinglesDigitizer/Singles/insert spatialResolution - /gate/digitizerMgr/crystalUnit/SinglesDigitizer/Singles/spatialResolution/fwhmXYdistrib2D my_distrib2D + /gate/digitizerMgr/crystalUnit/SinglesDigitizer/Singles/spatialResolution/nameAxis YZ + /gate/digitizerMgr/crystalUnit/SinglesDigitizer/Singles/spatialResolution/fwhmYdistrib2D my_distrib2D + /gate/digitizerMgr/crystalUnit/SinglesDigitizer/Singles/spatialResolution/fwhmZdistrib2D my_distrib2D + **Example for 1D distribution**:: /gate/distributions/name my_distrib1D /gate/distributions/insert File - /gate/distributions/my_distrib1D/setFileName macros/LutY.txt + /gate/distributions/my_distrib1D/setFileName macros/Lut_Y.txt /gate/distributions/my_distrib1D/read /gate/digitizerMgr/crystalUnit/SinglesDigitizer/Singles/insert spatialResolution /gate/digitizerMgr/crystalUnit/SinglesDigitizer/Singles/spatialResolution/fwhmYdistrib my_distrib1D @@ -566,6 +572,7 @@ BEWARE : The file for 2D Distribution should be structured such that: -28.50 11.38 11.18 10.23 -27.50 12.82 10.43 9.70 + Energy Framing ^^^^^^^^^^^^^^ *Previously Thresholder and Upholder* From 8d89e81b75de61f014bda0b4147fc3d2ce83c19e Mon Sep 17 00:00:00 2001 From: granadogmarc Date: Tue, 3 Feb 2026 17:25:46 +0100 Subject: [PATCH 06/12] Update digitizer_and_detector_modeling.rst --- docs/digitizer_and_detector_modeling.rst | 9 +++++---- 1 file changed, 5 insertions(+), 4 deletions(-) diff --git a/docs/digitizer_and_detector_modeling.rst b/docs/digitizer_and_detector_modeling.rst index 052ed47e1..35c856554 100644 --- a/docs/digitizer_and_detector_modeling.rst +++ b/docs/digitizer_and_detector_modeling.rst @@ -567,10 +567,11 @@ BEWARE : The file for 2D Distribution should be structured such that: **Example**:: --29.50 -28.50 -27.50 --29.50 9.62 13.66 10.22 --28.50 11.38 11.18 10.23 --27.50 12.82 10.43 9.70 +-30 -15 0 15 30 +-15 9.31 7.25 6.22 7.31 9.72 + 0 9.42 6.25 3.25 6.22 9.74 + 15 9.42 6.53 3.15 6.32 9.73 + 30 9.42 7.45 6.25 7.32 9.74 Energy Framing From e05587930f20d8d7b6bf50b14138d7b75a61d16a Mon Sep 17 00:00:00 2001 From: granadogmarc Date: Tue, 3 Feb 2026 17:27:23 +0100 Subject: [PATCH 07/12] Update digitizer_and_detector_modeling.rst --- docs/digitizer_and_detector_modeling.rst | 10 +++++----- 1 file changed, 5 insertions(+), 5 deletions(-) diff --git a/docs/digitizer_and_detector_modeling.rst b/docs/digitizer_and_detector_modeling.rst index 35c856554..d6ef8a8de 100644 --- a/docs/digitizer_and_detector_modeling.rst +++ b/docs/digitizer_and_detector_modeling.rst @@ -567,11 +567,11 @@ BEWARE : The file for 2D Distribution should be structured such that: **Example**:: --30 -15 0 15 30 --15 9.31 7.25 6.22 7.31 9.72 - 0 9.42 6.25 3.25 6.22 9.74 - 15 9.42 6.53 3.15 6.32 9.73 - 30 9.42 7.45 6.25 7.32 9.74 +-30 -15 0 15 30 +-15 9.31 7.25 6.22 7.31 9.73 +0 9.42 6.25 3.25 6.22 9.72 +15 9.42 6.53 3.15 6.32 9.71 +30 9.42 7.45 6.25 7.32 9.74 Energy Framing From 18f026bc1646f72c78ad3b3f902cbeb9b7281be9 Mon Sep 17 00:00:00 2001 From: granadogmarc Date: Tue, 3 Feb 2026 17:30:30 +0100 Subject: [PATCH 08/12] Update digitizer_and_detector_modeling.rst Correction in the distribution example --- docs/digitizer_and_detector_modeling.rst | 12 +++++++----- 1 file changed, 7 insertions(+), 5 deletions(-) diff --git a/docs/digitizer_and_detector_modeling.rst b/docs/digitizer_and_detector_modeling.rst index d6ef8a8de..227120224 100644 --- a/docs/digitizer_and_detector_modeling.rst +++ b/docs/digitizer_and_detector_modeling.rst @@ -567,11 +567,13 @@ BEWARE : The file for 2D Distribution should be structured such that: **Example**:: --30 -15 0 15 30 --15 9.31 7.25 6.22 7.31 9.73 -0 9.42 6.25 3.25 6.22 9.72 -15 9.42 6.53 3.15 6.32 9.71 -30 9.42 7.45 6.25 7.32 9.74 + + -30 -15 0 15 30 + -15 9.31 7.25 6.22 7.31 9.73 + 0 9.42 6.25 3.25 6.22 9.72 + 15 9.42 6.53 3.15 6.32 9.71 + 30 9.42 7.45 6.25 7.32 9.74 + Energy Framing From d0941a5b9173fa8da9d9863c2017a5b2154b6b90 Mon Sep 17 00:00:00 2001 From: granadogmarc Date: Tue, 3 Feb 2026 17:36:09 +0100 Subject: [PATCH 09/12] Update digitizer_and_detector_modeling.rst Changes in documentation for spatial resolution 2D distribution --- docs/digitizer_and_detector_modeling.rst | 46 +++++++++++------------- 1 file changed, 20 insertions(+), 26 deletions(-) diff --git a/docs/digitizer_and_detector_modeling.rst b/docs/digitizer_and_detector_modeling.rst index 9b6f03e8a..d4d254131 100644 --- a/docs/digitizer_and_detector_modeling.rst +++ b/docs/digitizer_and_detector_modeling.rst @@ -533,32 +533,31 @@ BEWARE: The confined and the use of the truncated Gaussian are default options. /gate/digitizerMgr/crystal/SinglesDigitizer/Singles/spatialResolution/confineInsideOfSmallestElement true /gate/digitizerMgr/pseudoCrystal/SinglesDigitizer/Singles/spatialResolution/useTruncatedGaussian true - **Configuring Spatial Resolution with 1D and 2D Distributions**:: -This approach is particularly essential for monolithic crystal detectors, where factors like edge effects and interaction positions significantly may influence spatial resolution. The spatial distribution resolution allows to select what axis will be using such distribution. Currently only one distribution is allowed for all axis. +This approach is particularly essential for monolithic crystal detectors, where factors such as edge effects and interaction positions may significantly influence spatial resolution. + +In order to map the 2D distribution with the crystal, the command "nameAxis" is used. The value of the axis provided will be the ones attributed to the 2 dimensions (columns and rows) of the distribution file. Assuming the crystal will be placed along the Z direction the options are "XZ" for a ring starting on top, and "YZ" for a ring starting on the sides. The default value for "nameAxis" is "YZ". Here is an example of how to configure this in a macro file: **Example for 2D distribution**:: -The axis are selected with a "nameAxis". IMPORTANT NOTE: Only the "XY" axis are available: - -/gate/distributions/name my_distrib2D -/gate/distributions/insert File -/gate/distributions/my_distrib2D/setFileName data/my_stddev_distribuion_file.txt -/gate/distributions/my_distrib2D/readMatrix2d -/gate/digitizerMgr/pseudoCrystal/SinglesDigitizer/Singles/insert spatialResolution -/gate/digitizerMgr/pseudoCrystal/SinglesDigitizer/Singles/spatialResolution/nameAxis XY -/gate/digitizerMgr/pseudoCrystal/SinglesDigitizer/Singles/spatialResolution/fwhmDistrib2D my_distrib2D - + /gate/distributions/name my_distrib2D + /gate/distributions/insert File + /gate/distributions/my_distrib2D/setFileName Lut_XY.txt + /gate/distributions/my_distrib2D/readMatrix2d + /gate/digitizerMgr/crystalUnit/SinglesDigitizer/Singles/insert spatialResolution + /gate/digitizerMgr/crystalUnit/SinglesDigitizer/Singles/spatialResolution/nameAxis YZ + /gate/digitizerMgr/crystalUnit/SinglesDigitizer/Singles/spatialResolution/fwhmYdistrib2D my_distrib2D + /gate/digitizerMgr/crystalUnit/SinglesDigitizer/Singles/spatialResolution/fwhmZdistrib2D my_distrib2D **Example for 1D distribution**:: /gate/distributions/name my_distrib1D /gate/distributions/insert File - /gate/distributions/my_distrib1D/setFileName macros/LutY.txt + /gate/distributions/my_distrib1D/setFileName macros/Lut_Y.txt /gate/distributions/my_distrib1D/read /gate/digitizerMgr/crystalUnit/SinglesDigitizer/Singles/insert spatialResolution /gate/digitizerMgr/crystalUnit/SinglesDigitizer/Singles/spatialResolution/fwhmYdistrib my_distrib1D @@ -569,23 +568,18 @@ These commands allow for more precise control over the spatial resolution by usi BEWARE : The file for 2D Distribution should be structured such that: --The first line contains the values of the x position. See the example below where the first line has values from -29.5 to 29.5 in this example for a 59 mm crystal. +-The first line contains the x values. --Each subsequent line begins with the value of the y position followed by the standard deviation (stddev) values corresponding to each x value and y value pair. +-Each subsequent line begins with a y value followed by the standard deviation (stddev) values corresponding to each x value and y value pair. **Example**:: --29.50 -22.94 -16.39 -9.83 -3.28 3.28 9.83 16.39 22.94 29.50 --29.50 8.05 5.1 4.24 4.23 4.23 4.19 4.56 4.6 5.02 7.68 --22.94 4.76 2.39 2.31 2.35 2.4 2.35 2.36 2.32 2.44 4.52 --16.39 4.53 2.45 2.28 2.39 2.35 2.36 2.46 2.46 2.55 4.52 --9.83 4.38 2.35 2.13 2.09 2.13 2.07 2.14 2.11 2.33 4.44 --3.28 4.18 2.18 1.97 1.96 1.97 1.97 2.03 2.02 2.19 4.2 -3.28 4.2 2.27 2.04 2.01 2.03 2.06 1.98 1.97 2.24 4.13 -9.83 4.24 2.4 2.23 2.24 2.21 2.28 2.22 2.21 2.35 4.36 -16.39 4.32 2.59 2.47 2.55 2.56 2.52 2.45 2.42 2.55 4.51 -22.94 4.07 2.32 2.33 2.34 2.34 2.27 2.36 2.29 2.4 4.3 -29.50 7.2 4.58 4.07 3.9 3.89 3.83 4.09 3.99 4.47 7.08 + + -30 -15 0 15 30 + -15 9.31 7.25 6.22 7.31 9.73 + 0 9.42 6.25 3.25 6.22 9.72 + 15 9.42 6.53 3.15 6.32 9.71 + 30 9.42 7.45 6.25 7.32 9.74 Energy Framing ^^^^^^^^^^^^^^ From 170cde907ad6c8034ec6bc6e2d7a3ae8dc40529d Mon Sep 17 00:00:00 2001 From: granadogmarc Date: Tue, 10 Feb 2026 12:50:32 +0100 Subject: [PATCH 10/12] addition of sigma correction as done in GATE10 for TruncatedGaussian --- .../src/GateDistributionTruncatedGaussian.cc | 39 ++++++++++++++++++- 1 file changed, 37 insertions(+), 2 deletions(-) diff --git a/source/digits_hits/src/GateDistributionTruncatedGaussian.cc b/source/digits_hits/src/GateDistributionTruncatedGaussian.cc index 03d22c17e..fae837e99 100644 --- a/source/digits_hits/src/GateDistributionTruncatedGaussian.cc +++ b/source/digits_hits/src/GateDistributionTruncatedGaussian.cc @@ -21,6 +21,7 @@ See LICENSE.md for further details //#include "GateDistributionTruncatedGaussianMessenger.hh" #include +#include //#include #include "Randomize.hh" #include "GateTools.hh" @@ -144,7 +145,24 @@ G4double GateDistributionTruncatedGaussian::computeTruncatedSigmaStatic(G4double double mean_shift = (phi_lowLim - phi_hiLim) / Z; double variance_correction = 1 - (lowLim_std * phi_lowLim - hiLim_std * phi_hiLim) / Z - mean_shift * mean_shift; - return sigma * sqrt(variance_correction); + double truncatedSigma = sigma * sqrt(variance_correction); + + // Edge correction: empirical fit to compensate for edge effects + constexpr double A = 3.61070188; + constexpr double B = -0.86538264; + + double delta = std::min(mu - lowLimit, highLimit - mu); + delta = std::max(0., delta); + + double delta_norm = delta / sigma; + double correction = 1; + + // Cuts on delta_norm > 0.95 and delta_norm < 4 comply with the safe limits of the fit + delta_norm = std::max(.95, delta_norm); + if (delta_norm < 4) + correction = A * std::exp(B * delta_norm); + + return truncatedSigma * correction; } //___________________________________________________________________ @@ -162,5 +180,22 @@ G4double GateDistributionTruncatedGaussian::computeTruncatedSigma() const{ double mean_shift = (phi_lowLim - phi_hiLim) / Z; double variance_correction = 1 - (lowLim_std * phi_lowLim - hiLim_std * phi_hiLim) / Z - mean_shift * mean_shift; - return m_Sigma * sqrt(variance_correction); + double truncatedSigma = m_Sigma * sqrt(variance_correction); + + // Edge correction: empirical fit to compensate for edge effects + constexpr double A = 3.61070188; + constexpr double B = -0.86538264; + + double delta = std::min(m_Mu - m_lowLimit, m_highLimit - m_Mu); + delta = std::max(0., delta); + + double delta_norm = delta / m_Sigma; + double correction = 1; + + // Cuts on delta_norm > 0.95 and delta_norm < 4 comply with the safe limits of the fit + delta_norm = std::max(.95, delta_norm); + if (delta_norm < 4) + correction = A * std::exp(B * delta_norm); + + return truncatedSigma * correction; } \ No newline at end of file From 85fb9e3d6311816ceb4c7d8e6c3392956b3e7562 Mon Sep 17 00:00:00 2001 From: kochebina Date: Wed, 11 Feb 2026 10:13:26 +0100 Subject: [PATCH 11/12] main.yml --- .github/workflows/main.yml | 57 ++++++++++++++++++++++---------------- 1 file changed, 33 insertions(+), 24 deletions(-) diff --git a/.github/workflows/main.yml b/.github/workflows/main.yml index be995ba51..975acea62 100644 --- a/.github/workflows/main.yml +++ b/.github/workflows/main.yml @@ -23,9 +23,9 @@ jobs: options: torch env: - ROOT_VERSION: 'v6-32-02' + ROOT_VERSION: 'v6-38-00' ITK_VERSION: 'v5.4.2' - GEANT4_VERSION: 'v11.3.0' + GEANT4_VERSION: 'v11.4.0' ROOT_DIR: $(HOME)/software/root GEANT4_DIR: $(HOME)/software/geant4 @@ -41,8 +41,8 @@ jobs: uses: actions/cache@v3 with: path: ~/software - key: ${{ matrix.os }}-geant4-${{ env.GEANT4_VERSION }}-root-${{ env.ROOT_VERSION }}-build1 - restore-keys: ${{ matrix.os }}-geant4-${{ env.GEANT4_VERSION }}-root-${{ env.ROOT_VERSION }}-build1 + key: ${{ matrix.os }}-geant4-${{ env.GEANT4_VERSION }}-root-${{ env.ROOT_VERSION }}-qt6-root3800 + restore-keys: ${{ matrix.os }}-geant4-${{ env.GEANT4_VERSION }}-root-${{ env.ROOT_VERSION }}-qt6-root3800 - name: Install dependencies run: | if [ "${{ matrix.os }}" == "ubuntu-latest" ]; then @@ -56,8 +56,10 @@ jobs: libxpm-dev \ libxft-dev \ libxext-dev \ - qtbase5-dev \ - qt5-qmake \ + qt6-base-dev \ + qt6-base-dev-tools \ + qt6-tools-dev \ + qt6-tools-dev-tools \ git \ cmake \ libssl-dev \ @@ -65,24 +67,28 @@ jobs: ccache \ libfftw3-dev gcc -v - elif [ "${{ matrix.os }}" == "macos-latest" ]; then + elif [ "${{ matrix.os }}" == "macos-latest" ]; then + brew update || true brew install python@3.12 || true - brew link --overwrite python@3.12 - #brew update - brew cleanup - brew config - #rm -rf /usr/local/bin/python3.11-config /usr/local/bin/2to3-3.11 /usr/local/bin/idle3.11 /usr/local/bin/pydoc3.11 /usr/local/bin/python3.11 - #rm -rf /usr/local/bin/python3-config /usr/local/bin/2to3 /usr/local/bin/idle3 /usr/local/bin/pydoc3 /usr/local/bin/python3 - brew install --force --verbose --overwrite --debug qt5 \ + + brew link --overwrite --force python@3.12 + export Python3_EXECUTABLE="$(brew --prefix python@3.12)/bin/python3" + brew install --force --verbose --overwrite --debug cmake \ + qt \ ccache \ tbb \ xrootd \ fftw - brew link qt5 --force && sudo ln -s /usr/local/opt/qt/mkspecs /usr/local/mkspecs && sudo ln -s /usr/local/opt/qt/plugins /usr/local/plugins - export PATH=/usr/local/opt/qt/bin:$PATH - export LDFLAGS="-L/usr/local/opt/qt/lib -L/usr/local/opt/llvm/lib" - export CPPFLAGS="-I/usr/local/opt/qt/include -I/usr/local/opt/llvm/include -fopenmp" - fi + + brew cleanup + #LLVM_PREFIX="$(brew --prefix llvm || true)" + #export SDKROOT="$(xcrun --sdk macosx --show-sdk-path)" + export CC=/usr/bin/clang + export CXX=/usr/bin/clang++ + + cmake --version + which cmake + fi cd $HOME/ mkdir -p software option_dependencies cmake --version @@ -95,13 +101,14 @@ jobs: mkdir src bin install git clone --branch $ROOT_VERSION https://github.com/root-project/root.git --depth 1 src cd bin + export CC=/usr/bin/clang + export CXX=/usr/bin/clang++ cmake ../src -DCMAKE_CXX_STANDARD=17 \ - -Dpython=OFF \ - -Dpyroot=OFF \ - -Dclad=OFF \ - -Dxrootd=OFF \ + -DCMAKE_C_COMPILER=$CC \ + -DCMAKE_CXX_COMPILER=$CXX \ + -DCMAKE_OSX_ARCHITECTURES=arm64 \ -DCMAKE_INSTALL_PREFIX=$HOME/software/root/install - make -j4 install + make install -j 4 cd .. rm -rf bin src fi @@ -116,8 +123,10 @@ jobs: cd bin cmake -DGEANT4_INSTALL_DATA=ON \ -DGEANT4_BUILD_MULTITHREADED=OFF \ + -DGEANT4_USE_QT=ON \ -DGEANT4_INSTALL_DATADIR=$HOME/software/geant4/data \ -DCMAKE_INSTALL_PREFIX=$HOME/software/geant4/install \ + -DQt6_DIR="$Qt6_DIR" \ ../src make -j4 install cd .. From 7c3715a1861cb778512d0e9ffe1f657dfff21f8f Mon Sep 17 00:00:00 2001 From: kochebina <42149373+kochebina@users.noreply.github.com> Date: Wed, 11 Feb 2026 10:59:04 +0100 Subject: [PATCH 12/12] Remove unit test build option for Gate Removed the option to build unit tests for Gate, including associated executable and installation commands. --- CMakeLists.txt | 10 ---------- 1 file changed, 10 deletions(-) diff --git a/CMakeLists.txt b/CMakeLists.txt index 7f41e332f..395a2798c 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -300,16 +300,6 @@ IF(GATE_COMPILE_GATEDIGIT) INSTALL(TARGETS Convert_CCMod2PETCoinc DESTINATION bin) ENDIF(GATE_COMPILE_GATEDIGIT) -# Option to build small unit tests for Gate (requires Geant4/ROOT as configured above) -OPTION(GATE_BUILD_UNIT_TESTS "Build unit tests for Gate" OFF) -IF(GATE_BUILD_UNIT_TESTS) - ADD_EXECUTABLE(TestSpatialResStdDev ${PROJECT_SOURCE_DIR}/source/digits_hits/test/TestSpatialResStdDev.cc $) - TARGET_LINK_LIBRARIES(TestSpatialResStdDev GateLib) - INSTALL(TARGETS TestSpatialResStdDev DESTINATION bin) - IF(BUILD_TESTING) - add_test(NAME TestSpatialResStdDev COMMAND TestSpatialResStdDev) - ENDIF(BUILD_TESTING) -ENDIF(GATE_BUILD_UNIT_TESTS) #========================================================= # We remove the warning option "shadow", because there are tons of