ProHEPaM/G4ava
3
0
Fork 0
Code Issues Pull requests Projects Releases Packages Wiki Activity Actions

Compare commits

base: ProHEPaM:f5e1a55e29c3bc619dafb62fd359f9a686ff3226
Branches Tags
ProHEPaM:main
..
compare: ProHEPaM:01fa1fbc555fe7f4cbbd7226f9aeafad2b23e428
Branches Tags
ProHEPaM:main

No commits in common. "f5e1a55e29c3bc619dafb62fd359f9a686ff3226" and "01fa1fbc555fe7f4cbbd7226f9aeafad2b23e428" have entirely different histories.
f5e1a55e29 ... 01fa1fbc55

20 changed files with 259 additions and 1244 deletions
Download patch file Download diff file Expand all files Collapse all files

1
.gitignore vendored
View file

@ -1 +0,0 @@
build

20
CMakeLists.txt
View file

@ -2,30 +2,18 @@ cmake_minimum_required(VERSION 3.2 FATAL_ERROR)
project(Tutorial) project(Tutorial)
# Nur einmal find_package, mit allen nötigen Features find_package(Geant4 REQUIRED ui_all vis_all)
find_package(Geant4 REQUIRED ui_all vis_all gdml)
# Setze Include-Verzeichnisse
include(${Geant4_USE_FILE}) include(${Geant4_USE_FILE})
include_directories(${PROJECT_SOURCE_DIR}/include) include_directories(${PROJECT_SOURCE_DIR}/include)
# Sammle alle .cc Dateien
file(GLOB sources ${PROJECT_SOURCE_DIR}/src/*.cc) file(GLOB sources ${PROJECT_SOURCE_DIR}/src/*.cc)
# Entferne unerwünschte Dateien
list(REMOVE_ITEM sources
${PROJECT_SOURCE_DIR}/src/ADetectorConstructionv2.cc
${PROJECT_SOURCE_DIR}/src/ASensitiveDetectorv2.cc
${PROJECT_SOURCE_DIR}/src/ARunMessenger.cc
${PROJECT_SOURCE_DIR}/src/ARunActionv2.cc
)
# Makros kopieren
file(GLOB MACRO_FILES "macros/*.mac") file(GLOB MACRO_FILES "macros/*.mac")
file(COPY ${MACRO_FILES} DESTINATION ${CMAKE_BINARY_DIR}/) file(COPY ${MACRO_FILES} DESTINATION ${CMAKE_BINARY_DIR}/)
# Executable erstellen
add_executable(sim sim.cc ${sources}) add_executable(sim sim.cc ${sources})
target_link_libraries(sim ${Geant4_LIBRARIES})
# Geant4 Libraries verlinken add_custom_target(Tutorial DEPENDS sim)
target_link_libraries(sim PUBLIC ${Geant4_LIBRARIES})

80
gdml/bgo.gdml
View file

@ -1,80 +0,0 @@
<?xml version="1.0" encoding="UTF-8"?>
<gdml xmlns:xsi="http://www.w3.org/2001/XMLSchema-instance"
xsi:noNamespaceSchemaLocation="http://service-spi.web.cern.ch/service-spi/app/releases/GDML/schema/gdml.xsd">
<define>
<position name="origin" unit="cm" x="0" y="0" z="0" />
</define>
<materials>
<element Z="1." formula="H" name="H">
<atom type="A" value="1.0079"/>
</element>
<element Z="6." formula="C" name="C">
<atom type="A" value="12.011"/>
</element>
<element Z="7." formula="N" name="N">
<atom type="A" value="14.007"/>
</element>
<element Z="8." formula="O" name="O">
<atom type="A" value="15.999"/>
</element>
<element Z="32." formula="Ge" name="Ge">
<atom type="A" value="72.63"/>
</element>
<element Z="83." formula="Bi" name="Bi">
<atom type="A" value="208.98"/>
</element>
<material Z="14." formula="Si" name="Silicon">
<D value="2.3296" unit="g/cm3"/>
<atom type="A" value="28.0855"/>
</material>
<material name="HiVacuum" formula="Vacuum">
<D value="1.29e-16" unit="g/cm3"/>
<fraction n="0.7" ref="N"/>
<fraction n="0.3" ref="O"/>
</material>
<material name="BGO" formula="BGO">
<D value="7.13" unit="g/cm3"/>
<composite n="12" ref="O"/>
<composite n="3" ref="Ge"/>
<composite n="4" ref="Bi"/>
</material>
<material name="BC430" formula="BC430">
<D value="1.032" unit="g/cm3"/>
<composite n="523" ref="H"/>
<composite n="472" ref="C"/>
</material>
</materials>
<solids>
<sphere name="world" rmin="0" rmax="25" deltaphi="360" deltatheta="180" aunit="deg" lunit="cm"/>
<tube name="bgo_solid" rmin="0" rmax="8" deltaphi="2*pi" z="10" aunit="radian" lunit="cm"/>
</solids>
<structure>
<volume name="BGO_vol">
<materialref ref="BGO"/>
<solidref ref="bgo_solid"/>
<auxiliary auxtype="sensi" auxvalue="BGO" auxunit="name"/>
</volume>
<volume name="World">
<materialref ref="HiVacuum"/>
<solidref ref="world"/>
<physvol>
<volumeref ref="BGO_vol"/>
<positionref ref="origin"/>
</physvol>
</volume>
</structure>
<setup name="Default" version="1.0">
<world ref="World"/>
</setup>
</gdml>

94
gdml/bgo_zylinder_2cm.gdml
View file

@ -1,94 +0,0 @@
<?xml version="1.0" encoding="UTF-8"?>
<gdml xmlns:xsi="http://www.w3.org/2001/XMLSchema-instance"
xsi:noNamespaceSchemaLocation="http://service-spi.web.cern.ch/service-spi/app/releases/GDML/schema/gdml.xsd">
<define>
<position name="origin" unit="cm" x="0" y="0" z="0"/>
</define>
<materials>
<element Z="1" name="H">
<atom value="1.0079"/>
</element>
<element Z="6" name="C">
<atom value="12.011"/>
</element>
<element Z="7" name="N">
<atom value="14.007"/>
</element>
<element Z="8" name="O">
<atom value="15.999"/>
</element>
<element Z="32" name="Ge">
<atom value="72.63"/>
</element>
<element Z="83" name="Bi">
<atom value="208.98"/>
</element>
<material name="HiVacuum">
<D value="1.29e-16" unit="g/cm3"/>
<fraction n="0.7" ref="N"/>
<fraction n="0.3" ref="O"/>
</material>
<material name="BGO">
<D value="7.13" unit="g/cm3"/>
<composite n="12" ref="O"/>
<composite n="3" ref="Ge"/>
<composite n="4" ref="Bi"/>
</material>
</materials>
<solids>
<sphere name="world" rmin="0" rmax="25"
deltaphi="360"
deltatheta="180"
aunit="deg"
lunit="cm"/>
<tube name="bgo_cylinder"
rmin="0"
rmax="2.5"
z="2"
deltaphi="360"
aunit="deg"
lunit="cm"/>
</solids>
<structure>
<volume name="BGO_vol">
<materialref ref="BGO"/>
<solidref ref="bgo_cylinder"/>
<auxiliary auxtype="sensi" auxvalue="BGO"/>
</volume>
<volume name="World">
<materialref ref="HiVacuum"/>
<solidref ref="world"/>
<physvol>
<volumeref ref="BGO_vol"/>
<positionref ref="origin"/>
</physvol>
</volume>
</structure>
<setup name="Default" version="1.0">
<world ref="World"/>
</setup>
</gdml>

115
gdml/shower_setup1.gdml
View file

@ -1,115 +0,0 @@
<?xml version="1.0" encoding="UTF-8"?>
<gdml xmlns:xsi="http://www.w3.org/2001/XMLSchema-instance"
xsi:noNamespaceSchemaLocation="http://service-spi.web.cern.ch/service-spi/app/releases/GDML/schema/gdml.xsd">
<define>
<position name="origin" unit="cm" x="0" y="0" z="0"/>
</define>
<materials>
<!-- Elemente -->
<element Z="1" formula="H" name="H"><atom type="A" value="1.0079"/></element>
<element Z="6" formula="C" name="C"><atom type="A" value="12.011"/></element>
<element Z="8" formula="O" name="O"><atom type="A" value="15.999"/></element>
<element Z="14" formula="Si" name="Si"><atom type="A" value="28.0855"/></element>
<element Z="32" formula="Ge" name="Ge"><atom type="A" value="72.63"/></element>
<element Z="83" formula="Bi" name="Bi"><atom type="A" value="208.98"/></element>
<!-- Materialien -->
<material name="HiVacuum" formula="Vacuum">
<D value="1.29e-16" unit="g/cm3"/>
<fraction n="0.7" ref="H"/>
<fraction n="0.3" ref="O"/>
</material>
<material name="BGO" formula="BGO">
<D value="7.13" unit="g/cm3"/>
<composite n="12" ref="O"/>
<composite n="3" ref="Ge"/>
<composite n="4" ref="Bi"/>
</material>
<material name="Silicon">
<D value="2.3296" unit="g/cm3"/>
<fraction n="1" ref="Si"/>
</material>
</materials>
<solids>
<box name="World" x="50" y="50" z="100" lunit="cm"/>
<box name="SSD_solid" x="5" y="5" z="0.025" lunit="cm"/>
<box name="BGO_solid" x="5" y="5" z="2.0" lunit="cm"/>
</solids>
<structure>
<!-- Layer Volumes mit Transparenz -->
<volume name="SSD1">
<materialref ref="Silicon"/>
<solidref ref="SSD_solid"/>
<auxiliary auxtype="sensi" auxvalue="shower_setup"/>
<auxiliary auxtype="VisAttributes" auxvalue="orange 0.6"/>
</volume>
<volume name="BGO1">
<materialref ref="BGO"/>
<solidref ref="BGO_solid"/>
<auxiliary auxtype="sensi" auxvalue="shower_setup"/>
<auxiliary auxtype="VisAttributes" auxvalue="green 0.6"/>
</volume>
<volume name="SSD2">
<materialref ref="Silicon"/>
<solidref ref="SSD_solid"/>
<auxiliary auxtype="sensi" auxvalue="shower_setup"/>
<auxiliary auxtype="VisAttributes" auxvalue="orange 0.6"/>
</volume>
<volume name="BGO2">
<materialref ref="BGO"/>
<solidref ref="BGO_solid"/>
<auxiliary auxtype="sensi" auxvalue="shower_setup"/>
<auxiliary auxtype="VisAttributes" auxvalue="green 0.6"/>
</volume>
<volume name="SSD3">
<materialref ref="Silicon"/>
<solidref ref="SSD_solid"/>
<auxiliary auxtype="sensi" auxvalue="shower_setup"/>
<auxiliary auxtype="VisAttributes" auxvalue="orange 0.6"/>
</volume>
<volume name="BGO3">
<materialref ref="BGO"/>
<solidref ref="BGO_solid"/>
<auxiliary auxtype="sensi" auxvalue="shower_setup"/>
<auxiliary auxtype="VisAttributes" auxvalue="green 0.6"/>
</volume>
<volume name="SSD4">
<materialref ref="Silicon"/>
<solidref ref="SSD_solid"/>
<auxiliary auxtype="sensi" auxvalue="shower_setup"/>
<auxiliary auxtype="VisAttributes" auxvalue="orange 0.6"/>
</volume>
<!-- World Volume mit allen Layern physisch -->
<volume name="World_vol">
<materialref ref="HiVacuum"/>
<solidref ref="World"/>
<physvol><volumeref ref="SSD1"/><position unit="cm" x="0" y="0" z="0.025"/></physvol>
<physvol><volumeref ref="BGO1"/><position unit="cm" x="0" y="0" z="1.55"/></physvol>
<physvol><volumeref ref="SSD2"/><position unit="cm" x="0" y="0" z="3.075"/></physvol>
<physvol><volumeref ref="BGO2"/><position unit="cm" x="0" y="0" z="4.6"/></physvol>
<physvol><volumeref ref="SSD3"/><position unit="cm" x="0" y="0" z="6.125"/></physvol>
<physvol><volumeref ref="BGO3"/><position unit="cm" x="0" y="0" z="7.65"/></physvol>
<physvol><volumeref ref="SSD4"/><position unit="cm" x="0" y="0" z="9.175"/></physvol>
</volume>
</structure>
<setup name="Default" version="1.0">
<world ref="World_vol"/>
</setup>
</gdml>

24
include/AActionInitialization.hh
View file

@ -1,24 +1,18 @@
#ifndef AACTIONINITIALIZATION_H #ifndef AACTIONINITIALIZATION_HH
#define AACTIONINITIALIZATION_H #define AACTIONINITIALIZATION_HH
#pragma once
#include "G4VUserActionInitialization.hh" #include "G4VUserActionInitialization.hh"
#include "APrimaryGenerator.hh"
class ARunAction; #include "ARunAction.hh"
class APrimaryGenerator;
class ASensitiveDetector;
class AActionInitialization : public G4VUserActionInitialization class AActionInitialization : public G4VUserActionInitialization
{ {
public: public:
AActionInitialization(ASensitiveDetector* sd); AActionInitialization();
virtual ~AActionInitialization(); ~AActionInitialization();
virtual void Build() const override;
virtual void BuildForMaster() const override;
private: virtual void BuildForMaster() const;
ASensitiveDetector* fSD; // Pointer auf den Sensitive Detector virtual void Build() const;
}; };
#endif #endif

48
include/ADetectorConstruction.hh
View file

@ -2,34 +2,48 @@
#define ADETECTORCONSTRUCTION_HH #define ADETECTORCONSTRUCTION_HH
#include "G4VUserDetectorConstruction.hh" #include "G4VUserDetectorConstruction.hh"
#include "G4GDMLParser.hh" #include "G4UImessenger.hh"
#include "G4String.hh" #include "G4Cache.hh"
#include "G4Box.hh"
#include "G4LogicalVolume.hh"
#include "G4VPhysicalVolume.hh"
#include "G4PVPlacement.hh"
#include "G4Material.hh"
#include "G4NistManager.hh"
#include "G4SystemOfUnits.hh"
#include "G4UnitsTable.hh"
#include "G4VisAttributes.hh"
#include "G4Color.hh"
#include "G4SDManager.hh"
#include "ASensitiveDetector.hh"
#include "G4UImessenger.hh"
#include "G4Cache.hh"
#include "G4UImessenger.hh"
#include "G4Cache.hh"
#include <vector> #include <vector>
class G4LogicalVolume; //#include "G4GDMLParser.hh"
class G4VPhysicalVolume;
class ASensitiveDetector;
class ADetectorConstruction : public G4VUserDetectorConstruction
class ADetectorConstruction : public G4VUserDetectorConstruction, public G4UImessenger
{ {
public: public:
ADetectorConstruction(); ADetectorConstruction();
virtual ~ADetectorConstruction(); virtual ~ADetectorConstruction();
virtual G4VPhysicalVolume* Construct(); // GDML laden virtual G4VPhysicalVolume *Construct();
virtual void ConstructSDandField(); // SD setzen
void ConstructSDandField();
void SetGDMLFile(const G4String& file); // GDML-Datei setzen
ASensitiveDetector* GetSensitiveDetector() const { return fSD; }
private: private:
G4GDMLParser fParser; // GDML Parser std::vector<G4LogicalVolume*> detectorVolumes;
G4String fGDMLFile; // GDML Datei
ASensitiveDetector* fSD; // Sensitive Detector
std::vector<G4LogicalVolume*> detectorVolumes; // Alle Volumes mit SD
}; };
#endif #endif

34
include/AMessenger.hh
View file

@ -1,34 +0,0 @@
#pragma once
#include "G4UImessenger.hh"
#include "G4UIcmdWithAString.hh"
#include "G4UIcmdWithADoubleAndUnit.hh"
class ASensitiveDetector;
class APrimaryGenerator;
class ADetectorConstruction;
class AMessenger : public G4UImessenger
{
public:
// Konstruktor: SD, optional Primärgenerator, optional DetectorConstruction
AMessenger(ASensitiveDetector* sd,
APrimaryGenerator* generator = nullptr,
ADetectorConstruction* detector = nullptr);
~AMessenger();
// UI-Kommandos verarbeiten
void SetNewValue(G4UIcommand* cmd, G4String value) override;
private:
ASensitiveDetector* fSD; // Referenz zum SensitiveDetector
APrimaryGenerator* fGenerator; // Referenz zum Primärgenerator
ADetectorConstruction* fDetector; // Referenz zur DetectorConstruction für GDML
// UI-Kommandos
G4UIcmdWithAString* fOutputFileCmd;
G4UIcmdWithAString* fParticleTypeCmd;
G4UIcmdWithADoubleAndUnit* fEnergyCmd;
G4UIcmdWithAString* fOutputColumnsCmd;
G4UIcmdWithAString* fGDMLFileCmd; // Neu: GDML-Datei setzen
};

7
include/APrimaryGenerator.hh
View file

@ -13,15 +13,10 @@ public:
APrimaryGenerator(); APrimaryGenerator();
~APrimaryGenerator(); ~APrimaryGenerator();
virtual void GeneratePrimaries(G4Event * anEvent) override; virtual void GeneratePrimaries(G4Event *);
void SetParticleEnergy(G4double energy) { fParticleEnergy = energy; }
void SetParticleName(const G4String& name) { fParticleName = name; }
private: private:
G4ParticleGun *fParticleGun; G4ParticleGun *fParticleGun;
G4double fParticleEnergy;
G4String fParticleName;
}; };
#endif #endif

143
include/ASensitiveDetector.hh
View file

@ -1,149 +1,28 @@
// #pragma once
// #include "G4VSensitiveDetector.hh"
// #include "G4THitsCollection.hh"
// #include "G4Step.hh"
// #include <fstream>
// #include <vector>
// #include <string>
// class ASensitiveDetector : public G4VSensitiveDetector
// {
// public:
// ASensitiveDetector(const G4String& name);
// ~ASensitiveDetector() override;
// void Initialize(G4HCofThisEvent*) override;
// G4bool ProcessHits(G4Step* step, G4TouchableHistory*) override;
// void EndOfEvent(G4HCofThisEvent*) override;
// void SetOutputFilename(const G4String& filename) { fOutputFilename = filename; }
// void SetOutputColumns(const std::vector<G4String>& cols) { fOutputColumns = cols; }
// private:
// std::ofstream fOutputFile;
// G4String fOutputFilename;
// std::vector<G4String> fOutputColumns;
// };
#ifndef ASENSITIVEDETECTOR_HH #ifndef ASENSITIVEDETECTOR_HH
#define ASENSITIVEDETECTOR_HH #define ASENSITIVEDETECTOR_HH
#include "G4VSensitiveDetector.hh" #include "G4VSensitiveDetector.hh"
#include "G4HCofThisEvent.hh" #include "G4HCofThisEvent.hh"
#include "G4TouchableHistory.hh" #include "G4TouchableHistory.hh"
#include "G4Step.hh" #include "G4Step.hh"
#include "G4Track.hh"
#include "G4SystemOfUnits.hh" #include "G4SystemOfUnits.hh"
#include "G4UnitsTable.hh" #include "G4UnitsTable.hh"
#include <unordered_map>
#include <fstream>
#include <functional>
#include <map>
#include <vector>
class AMessenger;
//
// Makro für die OutputMap – gleiche Logik wie beim Kollegen
//
#define OUTMAP [](const HitInfo& hit, std::ostream& out)
//#define OUTMAP [&, safeString, safeDouble, safeVec3](const HitInfo& hit, std::ostream& out)
class ASensitiveDetector : public G4VSensitiveDetector class ASensitiveDetector : public G4VSensitiveDetector
{ {
public: public:
ASensitiveDetector(G4String name); ASensitiveDetector(G4String);
~ASensitiveDetector(); ~ASensitiveDetector();
virtual void Initialize(G4HCofThisEvent* hce) override;
virtual G4bool ProcessHits(G4Step* step, G4TouchableHistory* touch) override;
virtual void EndOfEvent(G4HCofThisEvent* hce) override;
void SetOutputFilename(const G4String& f) { fOutputFilename = f; }
//void SetOutputColumns(const std::vector<G4String>& cols) { fOutputColumns = cols; }
void SetOutputColumns(const std::vector<G4String>& cols);
void SetMessenger(AMessenger* m) { fMessenger = m; }
private: private:
AMessenger* fMessenger = nullptr; G4double fTotalEnergyDeposited;
G4String fOutputFilename;
G4int fCurrentEventID = -1; virtual void Initialize(G4HCofThisEvent *) override; //hit collection
virtual void EndOfEvent(G4HCofThisEvent *) override;
public:
// virtual G4bool ProcessHits(G4Step *, G4TouchableHistory *);
// HitInfo enthält ALLE üblichen Parameter
//
struct HitInfo
{
// --- Basis-identifikation ---
G4int eventID; // Nummer des Events
G4int trackID; // Track ID
G4int parentID; // Parent Track ID
// --- Teilchen ---
G4String particleName; // Name des Teilchens
G4int pdg; // PDG-Code
G4double particleEnergy;// kinetische Energie am Vertex
// --- Primärinformation ---
G4String primaryName; // Name des Primärteilchens
G4double primaryEnergy; // Primärenergie
// --- Energien ---
G4double edep; // Deposited Energy
G4double ekin; // aktuelle kinetische Energie
// --- Positionen ---
G4ThreeVector prePos; // Position vor dem Step
G4ThreeVector postPos; // Position nach dem Step
// --- Richtung ---
G4ThreeVector preDir; // Richtung vor dem Step
G4ThreeVector postDir; // Richtung nach dem Step
// --- Zeiten ---
G4double globalTime; // absolute Zeit
G4double localTime; // lokale Zeit
// --- Step-Information ---
G4double stepLength; // Step-Länge
G4String processName; // Prozess, der diesen Step erzeugt hat
// --- Detektor ---
G4String detectorName; // voller Name
G4String cleanName; // aufgeräumter Name ohne _phys
HitInfo() :
eventID(0), trackID(0), parentID(0),
pdg(0),
particleEnergy(0.), primaryEnergy(0.),
edep(0.), ekin(0.),
globalTime(0.), localTime(0.),
stepLength(0.)
{}
};
private:
//
// OutputMap: string → Lambda(hit, out)
//
using OutputLambda = std::function<void(const HitInfo&, std::ostream&)>;
std::map<std::string, OutputLambda> outputMap;
//
// interne Daten
//
std::vector<G4String> fOutputColumns;
std::unordered_map<G4int, HitInfo> fTrackHitMap;
//std::vector<StepHit> fStepHits;
//
// interne Hilfsfunktionen
//
void InitializeOutputMap();
}; };
#endif #endif

17
macros/defrun.mac
View file

@ -1,17 +0,0 @@
/run/numberOfThreads 2
# --- Output Datei setzen ---
/my/outputFilename outfiles/test1000_100MeV_e
# --- Particle Type setzen ---
/my/particleType e-
# --- Particle Energy setzen ---
/my/particleEnergy 100 MeV
# --- Spalten setzen (ohne Kommata) ---
/my/outputColumns event track parent part pdg ekin edep x y z dir step proc det
/run/initialize
/run/beamOn 10

2
macros/run.mac
View file

@ -2,5 +2,5 @@
/run/initialize /run/initialize
/run/beamOn 1 /run/beamOn 1000

20
macros/test_gdml.mac
View file

@ -1,20 +0,0 @@
# Anzahl Threads
/run/numberOfThreads 1
# --- Output-Datei ---
/my/outputFilename outfiles/gdml_test_zylinder.hits
# --- Partikeltyp ---
/my/particleType e-
# --- Partikelenergie ---
/my/particleEnergy 100 MeV
# --- Output-Spalten ---
/my/outputColumns event track parent part pdg edep ekin proc det
# --- Initialisierung Geant4 ---
/run/initialize
# --- Simulation starten (z.B. 50.000 Hits) ---
/run/beamOn 1000

2
macros/vis.mac
View file

@ -1,6 +1,6 @@
/run/initialize /run/initialize
/vis/open OGLI 600x600-0+0 # OGLI = OpenGL Immediate /vis/open OGL 600x600-0+x600
/vis/viewer/set/viewpointVector 1 1 1 /vis/viewer/set/viewpointVector 1 1 1
/vis/viewer/set/autoRefresh true /vis/viewer/set/autoRefresh true

156
sim.cc
View file

@ -1,4 +1,7 @@
//Including header files
#include <iostream> #include <iostream>
// Including Geant4 stuff, G4MT is for multithreading
#include "G4RunManager.hh" #include "G4RunManager.hh"
#include "G4MTRunManager.hh" #include "G4MTRunManager.hh"
#include "G4UImanager.hh" #include "G4UImanager.hh"
@ -12,121 +15,52 @@
int main(int argc, char** argv) int main(int argc, char** argv)
{ {
// --- UI: Interaktiv, falls keine Argumente --- //G4UIExecutive *ui = new G4UIExecutive(argc, argv);
G4UIExecutive* ui = (argc == 1) ? new G4UIExecutive(argc, argv) : nullptr; G4UIExecutive *ui = 0;
if (argc == 1)
#ifdef G4MULTITHREADED {
G4MTRunManager* runManager = new G4MTRunManager; ui = new G4UIExecutive(argc, argv);
#else
G4RunManager* runManager = new G4RunManager;
#endif
// --- Physics ---
runManager->SetUserInitialization(new APhysicsList());
// --- Detector Construction ---
ADetectorConstruction* detector = new ADetectorConstruction();
runManager->SetUserInitialization(detector);
// GDML-Datei setzen (1. Argument)
if(argc > 1) {
detector->SetGDMLFile(argv[1]);
} }
#ifdef G4MULTITHREADED
G4MTRunManager *runManager = new G4MTRunManager;
#else
G4RunManager *runManager = new G4RunManager;
#endif
// Physics list
runManager->SetUserInitialization(static_cast<G4VUserPhysicsList*>(new APhysicsList()));
// Detector Construction
runManager->SetUserInitialization(static_cast<G4VUserDetectorConstruction*>(new ADetectorConstruction()));
// Action initialization
runManager->SetUserInitialization(static_cast<G4VUserActionInitialization*>(new AActionInitialization()));
// --- User Actions ---
runManager->SetUserInitialization(new AActionInitialization(detector->GetSensitiveDetector())); G4VisManager *visManager = new G4VisExecutive();
// --- Visualization ---
G4VisManager* visManager = new G4VisExecutive();
visManager->Initialize(); visManager->Initialize();
// --- UI Manager --- G4UImanager *UImanager = G4UImanager::GetUIpointer();
G4UImanager* UImanager = G4UImanager::GetUIpointer();
if (ui)
if(ui) { {
// Interaktiv G4cout << "ui exists" << G4endl;
UImanager->ApplyCommand("/control/execute ../macros/vis_shower.mac"); UImanager->ApplyCommand("/control/execute vis.mac");
ui->SessionStart(); ui->SessionStart();
delete ui;
}
else if(argc > 2) {
// Batch mit Macro
G4String macroFile = argv[2];
UImanager->ApplyCommand("/control/execute " + macroFile);
}
else {
// Batch ohne Macro
UImanager->ApplyCommand("/run/initialize");
} }
else
// --- Clean up --- {
G4cout << "ui does not exists" << G4endl;
G4String command = "/control/execute ";
G4String fileName = argv[1];
UImanager->ApplyCommand(command + fileName);
}
//clean up
delete runManager; delete runManager;
delete visManager; delete visManager;
if (ui) delete ui;
return 0; return 0;
}
}
// #include "G4RunManager.hh"
// #include "G4MTRunManager.hh"
// #include "G4UImanager.hh"
// #include "G4VisManager.hh"
// #include "G4VisExecutive.hh"
// #include "G4UIExecutive.hh"
// #include "APhysicsList.hh"
// #include "ADetectorConstruction.hh"
// #include "AActionInitialization.hh"
// int main(int argc, char** argv)
// {
// // UI Executive für interaktive Sitzung
// G4UIExecutive* ui = nullptr;
// if(argc == 1) {
// ui = new G4UIExecutive(argc, argv);
// }
// // Run Manager (MT oder Single Thread)
// #ifdef G4MULTITHREADED
// G4MTRunManager* runManager = new G4MTRunManager;
// #else
// G4RunManager* runManager = new G4RunManager;
// #endif
// // Physics
// runManager->SetUserInitialization(new APhysicsList());
// // Detector
// G4String gdmlFile = (argc > 1) ? argv[1] : "";
// ADetectorConstruction* detector = new ADetectorConstruction(gdmlFile);
// runManager->SetUserInitialization(detector);
// // Actions
// runManager->SetUserInitialization(new AActionInitialization(detector->GetSensitiveDetector()));
// // Initialize RunManager
// runManager->Initialize();
// // Visualization
// G4VisManager* visManager = new G4VisExecutive();
// visManager->Initialize();
// // UI Manager
// G4UImanager* UImanager = G4UImanager::GetUIpointer();
// if(ui) {
// // Interaktive Sitzung: vis.mac automatisch ausführen
// UImanager->ApplyCommand("/control/execute ../macros/vis.mac");
// ui->SessionStart();
// } else if(argc > 2) {
// // Batch: Macro ausführen
// G4String macroFile = argv[2];
// UImanager->ApplyCommand("/control/execute " + macroFile);
// }
// // Clean up
// delete runManager;
// delete visManager;
// if(ui) delete ui;
// return 0;
// }

46
src/AActionInitialization.cc
View file

@ -1,28 +1,40 @@
#include "AActionInitialization.hh" #include "AActionInitialization.hh"
#include "AMessenger.hh" //#include "PrimaryGeneratorAction.hh"
#include "APrimaryGenerator.hh"
#include "ARunAction.hh"
AActionInitialization::AActionInitialization(ASensitiveDetector* sd) AActionInitialization::AActionInitialization()
: G4VUserActionInitialization(), fSD(sd)
{ {
} }
AActionInitialization::~AActionInitialization() {} AActionInitialization::~AActionInitialization()
void AActionInitialization::Build() const
{ {
// Primärgenerator
APrimaryGenerator* primaryGen = new APrimaryGenerator();
SetUserAction(primaryGen);
AMessenger* messenger = new AMessenger(fSD, primaryGen);
// RunAction
SetUserAction(new ARunAction());
} }
void AActionInitialization::BuildForMaster() const void AActionInitialization::BuildForMaster() const
{ {
SetUserAction(new ARunAction()); // This is needed to avoid problems due to multithreading when writing the data.
ARunAction *runAction = new ARunAction();
SetUserAction(runAction);
} }
void AActionInitialization::Build() const
{
// PrimaryGenerator and SteppingAction do not have to be set in BuildForMaster because they only
// regard individual threads.
APrimaryGenerator *generator = new APrimaryGenerator();
SetUserAction(generator);
ARunAction *runAction = new ARunAction();
SetUserAction(runAction);
// SteppingAction *steppingAction = new SteppingAction();
// SetUserAction(steppingAction);
}
// void UserActionInitialization::Build() const {
// SetUserAction(new PrimaryGeneratorAction);
// }
//
// void UserActionInitialization::BuildForMaster() const {
// }

165
src/ADetectorConstruction.cc
View file

@ -1,101 +1,114 @@
#include "ADetectorConstruction.hh" #include "ADetectorConstruction.hh"
#include "ASensitiveDetector.hh"
#include "AMessenger.hh"
#include "G4LogicalVolumeStore.hh"
#include "G4SDManager.hh"
#include "G4GDMLAuxStructType.hh"
#include "G4SystemOfUnits.hh"
#include "G4Tubs.hh"
#include "G4Box.hh"
#include "G4VisAttributes.hh"
#include "G4Color.hh"
#include "G4PVPlacement.hh"
#include "G4NistManager.hh"
#include "G4ios.hh"
ADetectorConstruction::ADetectorConstruction() ADetectorConstruction::ADetectorConstruction()
: fSD(nullptr)
{ {
// Ein einziges SensitiveDetector erzeugen
fSD = new ASensitiveDetector("shower_setup");
// Messenger für den SD erzeugen
AMessenger* messenger = new AMessenger(fSD);
fSD->SetMessenger(messenger);
} }
ADetectorConstruction::~ADetectorConstruction() ADetectorConstruction::~ADetectorConstruction()
{ {
delete fSD; // SD + Messenger aufräumen
} }
void ADetectorConstruction::SetGDMLFile(const G4String& file) G4VPhysicalVolume *ADetectorConstruction::Construct()
{ {
fGDMLFile = file; G4bool checkOverlaps = true;
}
G4VPhysicalVolume* ADetectorConstruction::Construct() // --- Materials ---
{ G4NistManager *nist = G4NistManager::Instance();
if(fGDMLFile.empty()) //G4Material *worldMat = nist->FindOrBuildMaterial("G4_AIR");
G4Material *worldMat = nist->FindOrBuildMaterial("G4_Galactic");
G4Material *ssdMat = nist->FindOrBuildMaterial("G4_Si");
G4Material *bgoMat = nist->FindOrBuildMaterial("G4_BGO");
// --- World volume ---
G4double worldSize = 1.0 * m;
G4Box *solidWorld = new G4Box("solidWorld", 0.5*worldSize, 0.5*worldSize, 0.5*worldSize);
G4LogicalVolume *logicWorld = new G4LogicalVolume(solidWorld, worldMat, "logicalWorld");
G4VPhysicalVolume *physWorld =
new G4PVPlacement(0, G4ThreeVector(), logicWorld,
"physWorld", 0, false, 0, checkOverlaps);
// --- Layer sizes ---
G4double sizeXY = 5.0 * cm;
G4double thSSD = 0.5 * mm;
G4double thBGO = 20.0 * mm;
G4double gap = 5.0 * mm;
// --- Startposition ---
G4double totalThickness = 4*thSSD + 3*thBGO + 6*gap;
G4double zPos = 0.0;
// --- Visualization ---
G4VisAttributes *ssdVis = new G4VisAttributes(G4Color(0.0, 0.0, 1.0, 0.6));
ssdVis->SetForceSolid(true);
G4VisAttributes *bgoVis = new G4VisAttributes(G4Color(0.0, 1.0, 0.0, 0.6));
bgoVis->SetForceSolid(true);
// --- Build SSD – BGO – SSD – BGO – SSD – BGO – SSD ---
for(int i = 0; i < 7; i++)
{ {
G4Exception("ADetectorConstruction::Construct()", G4bool isSSD = (i % 2 == 0);
"NoGDMLFile", G4double thickness = isSSD ? thSSD : thBGO;
FatalException, G4Material *mat = isSSD ? ssdMat : bgoMat;
"No GDML file specified.");
}
G4cout << "Loading GDML geometry: " << fGDMLFile << G4endl; // G4String baseName = isSSD ? "SSD_" : "BGO_";
// baseName += std::to_string(i);
int ssdIndex = 1;
int bgoIndex = 1;
G4String baseName;
// GDML-Datei einlesen if(isSSD) {
fParser.Read(fGDMLFile, false); baseName = "SSD" + std::to_string(ssdIndex);
ssdIndex++;
// Weltvolumen aus GDML holen
G4VPhysicalVolume* world = fParser.GetWorldVolume();
if(!world)
{
G4Exception("ADetectorConstruction::Construct()",
"InvalidGDML",
FatalException,
"GDML world volume not found.");
}
G4cout << "GDML geometry successfully loaded." << G4endl;
// --- LogicalVolumes iterieren und SD zuweisen ---
detectorVolumes.clear();
auto lvStore = G4LogicalVolumeStore::GetInstance();
for(auto lv : *lvStore)
{
// AuxList vom GDML Parser holen
G4GDMLAuxListType auxList = fParser.GetVolumeAuxiliaryInformation(lv);
for(const auto& aux : auxList)
{
if(aux.type == "sensi")
{
lv->SetSensitiveDetector(fSD);
detectorVolumes.push_back(lv); // merken für SD-Registrierung
}
} }
else {
baseName = "BGO" + std::to_string(bgoIndex);
bgoIndex++;
}
G4Box *solid = new G4Box(baseName,
0.5*sizeXY,
0.5*sizeXY,
0.5*thickness);
G4LogicalVolume *lv =
new G4LogicalVolume(solid, mat, baseName + "_logic");
new G4PVPlacement(0,
G4ThreeVector(0. ,0. , zPos + 0.5*thickness),
lv,
baseName + "_phys",
logicWorld,
false,
i,
checkOverlaps);
// Farbe
if(isSSD)
lv->SetVisAttributes(ssdVis);
else
lv->SetVisAttributes(bgoVis);
// Für SD merken
detectorVolumes.push_back(lv);
// Update z
zPos += thickness + gap;
} }
G4cout << "Assigned SensitiveDetector to " << detectorVolumes.size() << " volumes." << G4endl; return physWorld;
return world;
} }
void ADetectorConstruction::ConstructSDandField() void ADetectorConstruction::ConstructSDandField()
{ {
// SD bei Geant4 registrieren ASensitiveDetector *sensDet = new ASensitiveDetector("MySensitiveDetector");
G4SDManager::GetSDMpointer()->AddNewDetector(fSD); G4SDManager::GetSDMpointer()->AddNewDetector(sensDet);
// Alle Layer bekommen denselben SD (bereits in Construct gesetzt) // alle SSDs und BGOs sensitive machen
for(auto lv : detectorVolumes) for(size_t i = 0; i < detectorVolumes.size(); i++)
{ {
lv->SetSensitiveDetector(fSD); detectorVolumes[i]->SetSensitiveDetector(sensDet);
} }
}
G4cout << "SensitiveDetector registered for " << detectorVolumes.size() << " volumes." << G4endl;
}

85
src/AMessenger.cc
View file

@ -1,85 +0,0 @@
#include "AMessenger.hh"
#include "ASensitiveDetector.hh"
#include "APrimaryGenerator.hh"
#include "ADetectorConstruction.hh"
#include "G4SystemOfUnits.hh"
#include "G4UIcmdWithAString.hh"
#include "G4UIcmdWithADoubleAndUnit.hh"
#include "G4UImanager.hh"
#include <sstream>
AMessenger::AMessenger(ASensitiveDetector* sd, APrimaryGenerator* generator, ADetectorConstruction* detector)
: fSD(sd), fGenerator(generator), fDetector(detector)
{
// Kommando: Output-Datei für SD
fOutputFileCmd = new G4UIcmdWithAString("/my/outputFilename", this);
fOutputFileCmd->SetGuidance("Set output filename for SD text output.");
fOutputFileCmd->SetParameterName("filename", false);
// Kommando: Partikeltyp
fParticleTypeCmd = new G4UIcmdWithAString("/my/particleType", this);
fParticleTypeCmd->SetGuidance("Set particle type (e-, e+, gamma, proton, neutron).");
fParticleTypeCmd->SetParameterName("particle", false);
fParticleTypeCmd->AvailableForStates(G4State_PreInit, G4State_Idle);
// Kommando: Partikelenergie
fEnergyCmd = new G4UIcmdWithADoubleAndUnit("/my/particleEnergy", this);
fEnergyCmd->SetGuidance("Set particle energy.");
fEnergyCmd->SetParameterName("energy", false);
fEnergyCmd->SetUnitCategory("Energy");
fEnergyCmd->AvailableForStates(G4State_PreInit, G4State_Idle);
// Kommando: Output-Spalten
fOutputColumnsCmd = new G4UIcmdWithAString("/my/outputColumns", this);
fOutputColumnsCmd->SetGuidance("Set output columns and order (space separated). E.g.: event track parent part edep");
fOutputColumnsCmd->SetParameterName("columns", false);
fOutputColumnsCmd->AvailableForStates(G4State_PreInit, G4State_Idle);
// Kommando: GDML-Datei
fGDMLFileCmd = new G4UIcmdWithAString("/my/gdmlFile", this);
fGDMLFileCmd->SetGuidance("Set GDML input file");
fGDMLFileCmd->SetParameterName("filename", false);
fGDMLFileCmd->AvailableForStates(G4State_PreInit);
}
AMessenger::~AMessenger()
{
delete fOutputFileCmd;
delete fParticleTypeCmd;
delete fEnergyCmd;
delete fOutputColumnsCmd;
delete fGDMLFileCmd;
}
void AMessenger::SetNewValue(G4UIcommand* cmd, G4String value)
{
if (cmd == fOutputFileCmd && fSD)
{
fSD->SetOutputFilename(value);
}
else if (cmd == fEnergyCmd && fGenerator)
{
fGenerator->SetParticleEnergy(fEnergyCmd->GetNewDoubleValue(value));
}
else if (cmd == fParticleTypeCmd && fGenerator)
{
fGenerator->SetParticleName(value);
}
else if (cmd == fOutputColumnsCmd && fSD)
{
std::vector<G4String> cols;
std::istringstream iss(value);
G4String token;
while (iss >> token)
{
cols.push_back(token);
}
fSD->SetOutputColumns(cols);
}
else if (cmd == fGDMLFileCmd && fDetector)
{
fDetector->SetGDMLFile(value);
}
}

38
src/APrimaryGenerator.cc
View file

@ -3,16 +3,29 @@
APrimaryGenerator::APrimaryGenerator() APrimaryGenerator::APrimaryGenerator()
{ {
fParticleGun = new G4ParticleGun(1); // 1 particles per event fParticleGun = new G4ParticleGun(1); // 1 particles per event
fParticleEnergy = 100.*MeV; // Default-Energie
fParticleName = "e-"; //default
G4ThreeVector pos(0.,0.,-1.*cm); //Particle position
G4ThreeVector mom(0.,0.,1.); G4double x = 0. * m;
G4ParticleTable* particleTable = G4ParticleTable::GetParticleTable(); G4double y = 0. * m;
G4ParticleDefinition* particle = particleTable->FindParticle("e-"); G4double z = -1. * cm;
G4ThreeVector pos(x, y, z);
//Particle direction
G4double px = 0.;
G4double py = 0.;
G4double pz = 1.;
G4ThreeVector mom(px,py,pz);
//Particle type
G4ParticleTable *particleTable = G4ParticleTable::GetParticleTable();
// possibilities: "e+","e-","gamma", ,"proton", "neutron"
G4ParticleDefinition *particle = particleTable->FindParticle("e-");
fParticleGun->SetParticlePosition(pos); fParticleGun->SetParticlePosition(pos);
fParticleGun->SetParticleMomentumDirection(mom); fParticleGun->SetParticleMomentumDirection(mom);
fParticleGun->SetParticleEnergy(100 * MeV); //or GeV
fParticleGun->SetParticleDefinition(particle); fParticleGun->SetParticleDefinition(particle);
} }
@ -23,15 +36,6 @@ APrimaryGenerator::~APrimaryGenerator()
void APrimaryGenerator::GeneratePrimaries(G4Event *anEvent) void APrimaryGenerator::GeneratePrimaries(G4Event *anEvent)
{ {
G4ParticleTable* particleTable = G4ParticleTable::GetParticleTable(); //Create vertex
G4ParticleDefinition* particle = particleTable->FindParticle(fParticleName);
if (!particle)
{
G4cerr << "Unknown particle: " << fParticleName << ". Using e- as default." << G4endl;
particle = particleTable->FindParticle("e-");
}
fParticleGun->SetParticleDefinition(particle);
fParticleGun->SetParticleEnergy(fParticleEnergy); // dynamisch vom Macro
fParticleGun->GeneratePrimaryVertex(anEvent); fParticleGun->GeneratePrimaryVertex(anEvent);
} }

400
src/ASensitiveDetector.cc
View file

@ -1,404 +1,32 @@
#include "ASensitiveDetector.hh" #include "ASensitiveDetector.hh"
#include <cmath>
// #include "G4Step.hh"
// #include "G4SystemOfUnits.hh"
// #include <iomanip>
// ASensitiveDetector::ASensitiveDetector(const G4String& name) ASensitiveDetector::ASensitiveDetector(G4String name) : G4VSensitiveDetector(name)
// : G4VSensitiveDetector(name)
// {
// }
// ASensitiveDetector::~ASensitiveDetector()
// {
// if(fOutputFile.is_open())
// fOutputFile.close();
// }
// void ASensitiveDetector::Initialize(G4HCofThisEvent*)
// {
// if(!fOutputFile.is_open() && !fOutputFilename.empty())
// fOutputFile.open(fOutputFilename);
// }
// G4bool ASensitiveDetector::ProcessHits(G4Step* step, G4TouchableHistory*)
// {
// // Beispiel: speichere Position und Energieverlust
// auto pos = step->GetPreStepPoint()->GetPosition();
// auto edep = step->GetTotalEnergyDeposit();
// if(fOutputFile.is_open())
// fOutputFile << std::fixed << std::setprecision(3)
// << pos.z()/cm << " " << pos.perp()/cm << " " << edep/MeV << G4endl;
// return true;
// }
// void ASensitiveDetector::EndOfEvent(G4HCofThisEvent*)
// {
// fOutputFile.flush();
// }
#include "ASensitiveDetector.hh"
#include "AMessenger.hh"
#include "G4Event.hh"
#include "G4RunManager.hh"
#include "G4Threading.hh"
#include <filesystem>
#include <iostream>
ASensitiveDetector::ASensitiveDetector(G4String name)
: G4VSensitiveDetector(name), fOutputFilename("outfiles/shower_setup")
{ {
std::filesystem::create_directories("outfiles"); fTotalEnergyDeposited = 0.;
InitializeOutputMap();
if(fOutputColumns.empty()) {
fOutputColumns = {"event", "track", "parent", "part", "edep", "det"};
}
} }
ASensitiveDetector::~ASensitiveDetector() {} ASensitiveDetector::~ASensitiveDetector()
void ASensitiveDetector::Initialize(G4HCofThisEvent*)
{ {
fTrackHitMap.clear();
//fStepHits.clear();
} }
// void ASensitiveDetector::InitializeOutputMap() void ASensitiveDetector::Initialize(G4HCofThisEvent *)
// {
// // Lambda für sicheren Zugriff auf hit-Felder
// auto safeString = [](const G4String& s) -> G4String {
// return s.empty() ? "undefined" : s;
// };
// auto safeDouble = [](G4double value, G4double scale=1.0) -> G4double {
// return value; // optional: prüfen auf NaN oder negative Werte
// };
// auto safeVec3 = [](const G4ThreeVector& v, G4double scale=1.0) -> G4ThreeVector {
// return v / scale;
// };
// outputMap = {
// {"event", OUTMAP { out << hit.eventID; }},
// {"track", OUTMAP { out << hit.trackID; }},
// {"parent", OUTMAP { out << hit.parentID; }},
// {"part", OUTMAP { out << safeString(hit.particleName); }},
// {"pdg", OUTMAP { out << hit.pdg; }},
// {"ekin", OUTMAP { out << safeDouble(hit.ekin/MeV); }},
// {"edep", OUTMAP { out << safeDouble(hit.edep/MeV); }},
// {"x", OUTMAP { out << safeDouble(hit.postPos.x()/mm); }},
// {"y", OUTMAP { out << safeDouble(hit.postPos.y()/mm); }},
// {"z", OUTMAP { out << safeDouble(hit.postPos.z()/mm); }},
// {"pos", OUTMAP {
// auto v = safeVec3(hit.postPos, mm);
// out << v.x() << " " << v.y() << " " << v.z();
// }},
// {"r", OUTMAP {
// G4double r = std::sqrt(hit.postPos.x()*hit.postPos.x() + hit.postPos.y()*hit.postPos.y());
// out << r;
// }},
// {"dx", OUTMAP { out << hit.postDir.x(); }},
// {"dy", OUTMAP { out << hit.postDir.y(); }},
// {"dz", OUTMAP { out << hit.postDir.z(); }},
// {"dir", OUTMAP { out << hit.postDir.x() << " " << hit.postDir.y() << " " << hit.postDir.z(); }},
// {"global_t", OUTMAP { out << hit.globalTime; }},
// {"local_t", OUTMAP { out << hit.localTime; }},
// {"step", OUTMAP { out << hit.stepLength / mm; }},
// {"proc", OUTMAP { out << safeString(hit.processName); }},
// {"det", OUTMAP { out << safeString(hit.cleanName); }}
// };
// }
void ASensitiveDetector::InitializeOutputMap()
{ {
outputMap = { fTotalEnergyDeposited = 0.;
{"event", OUTMAP { out<<hit.eventID; }},
{"track", OUTMAP { out<<hit.trackID; }},
{"parent",OUTMAP { out<<hit.parentID; }},
{"part", OUTMAP { out<<hit.particleName; }},
{"pdg", OUTMAP { out<<hit.pdg; }},
{"primary", OUTMAP { out<<hit.primaryName; }},
{"primaryE",OUTMAP { out<<hit.primaryEnergy; }},
{"edep", OUTMAP { out<<hit.edep; }},
{"ekin", OUTMAP { out<<hit.ekin; }},
{"x", OUTMAP { out<<hit.postPos.x(); }},
{"y", OUTMAP { out<<hit.postPos.y(); }},
{"z", OUTMAP { out<<hit.postPos.z(); }},
{"r", OUTMAP {
out<<std::sqrt(hit.postPos.x()*hit.postPos.x()
+ hit.postPos.y()*hit.postPos.y());
}},
{"dx", OUTMAP { out<<hit.postDir.x(); }},
{"dy", OUTMAP { out<<hit.postDir.y(); }},
{"dz", OUTMAP { out<<hit.postDir.z(); }},
{"global_t", OUTMAP { out<<hit.globalTime; }},
{"local_t", OUTMAP { out<<hit.localTime; }},
{"step", OUTMAP { out<<hit.stepLength; }},
{"proc", OUTMAP { out<<hit.processName; }},
{"det", OUTMAP { out<<hit.cleanName; }}
};
} }
void ASensitiveDetector::SetOutputColumns( void ASensitiveDetector::EndOfEvent(G4HCofThisEvent *)
const std::vector<G4String>& cols)
{ {
for(const auto& c : cols) G4cout << "Deposited energy: " << fTotalEnergyDeposited << G4endl;
{
if(!outputMap.count(c))
{
G4Exception("InvalidColumn",
"BadColumn",
FatalException,
("Unknown column: " + c).c_str());
}
}
fOutputColumns = cols;
} }
// G4bool ASensitiveDetector::ProcessHits(G4Step* step, G4TouchableHistory*) G4bool ASensitiveDetector::ProcessHits(G4Step *aStep, G4TouchableHistory *touchHist)
// {
// if(!step) return false;
// G4double edep = step->GetTotalEnergyDeposit();
// if(edep <= 0) return false;
// auto track = step->GetTrack();
// auto post = step->GetPostStepPoint();
// if(!track || !post) return false;
// auto pos = post->GetPosition();
// StepHit hit;
// hit.eventID = G4RunManager::GetRunManager()
// ->GetCurrentEvent()->GetEventID();
// hit.trackID = track->GetTrackID();
// hit.parentID = track->GetParentID();
// hit.particleName = track->GetDefinition()->GetParticleName();
// hit.edep = edep / MeV;
// hit.z = pos.z() / cm;
// hit.r = std::sqrt(pos.x()*pos.x() + pos.y()*pos.y()) / cm;
// fStepHits.push_back(hit);
// return true;
// }
G4bool ASensitiveDetector::ProcessHits(G4Step* step, G4TouchableHistory*)
{ {
if(!step) return false; G4double fEnergyDeposited = aStep->GetTotalEnergyDeposit();
G4Track* track = step->GetTrack(); if(fEnergyDeposited > 0)
if(!track) return false;
const G4Event* event =
G4RunManager::GetRunManager()->GetCurrentEvent();
if(!event) return false;
G4int trackID = track->GetTrackID();
G4int eventID = event->GetEventID();
// Referenz auf Track-Container
HitInfo& hit = fTrackHitMap[trackID];
// ===== Basis =====
hit.eventID = eventID;
hit.trackID = trackID;
hit.parentID = track->GetParentID();
if(track->GetDefinition())
{ {
hit.particleName = fTotalEnergyDeposited += fEnergyDeposited;
track->GetDefinition()->GetParticleName();
hit.pdg =
track->GetDefinition()->GetPDGEncoding();
} }
hit.ekin = track->GetKineticEnergy()/MeV;
// ===== Primärinfo =====
if(auto vertex = event->GetPrimaryVertex())
{
if(auto prim = vertex->GetPrimary())
{
if(auto code = prim->GetG4code())
{
hit.primaryName =
code->GetParticleName();
hit.primaryEnergy =
prim->GetKineticEnergy()/MeV;
}
}
}
// ===== Energie aufsummieren =====
hit.edep += step->GetTotalEnergyDeposit()/MeV;
// ===== Step Info =====
auto pre = step->GetPreStepPoint();
auto post = step->GetPostStepPoint();
if(pre && post)
{
hit.prePos = pre->GetPosition()/mm;
hit.postPos = post->GetPosition()/mm;
hit.preDir = pre->GetMomentumDirection();
hit.postDir = post->GetMomentumDirection();
hit.globalTime = post->GetGlobalTime()/ns;
hit.localTime = post->GetLocalTime()/ns;
hit.stepLength += step->GetStepLength()/mm;
hit.processName =
post->GetProcessDefinedStep()
? post->GetProcessDefinedStep()->GetProcessName()
: "none";
if(pre->GetTouchable())
{
auto vol = pre->GetTouchable()->GetVolume();
if(vol)
{
hit.detectorName = vol->GetName();
auto pos = hit.detectorName.find("_phys");
hit.cleanName =
(pos == std::string::npos)
? hit.detectorName
: hit.detectorName.substr(0,pos);
}
}
}
return true; return true;
} }
// void ASensitiveDetector::EndOfEvent(G4HCofThisEvent*)
// {
// // Thread-ID für Datei
// G4int tid = G4Threading::G4GetThreadId();
// // Thread-local Datei
// thread_local std::ofstream out;
// // Datei öffnen, wenn noch nicht offen
// if(!out.is_open()) {
// std::filesystem::create_directories("outfiles");
// std::string threadFile = fOutputFilename + "_" + std::to_string(tid) + ".hits";
// out.open(threadFile, std::ios::app);
// if(!out.is_open()) {
// G4cerr << "Cannot open output file: " << threadFile << G4endl;
// return;
// }
// // Header einmal pro Datei
// for(size_t i=0; i<fOutputColumns.size(); ++i) {
// out << fOutputColumns[i];
// if(i+1 < fOutputColumns.size()) out << "\t";
// }
// out << "\n";
// }
// // Hits schreiben
// // for(auto& [id, hit] : fTrackHitMap) {
// // for(size_t i=0; i<fOutputColumns.size(); ++i) {
// // const auto& col = fOutputColumns[i];
// // if(outputMap.count(col))
// // outputMap[col](hit, out);
// // else
// // out << "NA";
// // if(i+1 < fOutputColumns.size()) out << "\t";
// // }
// // out << "\n";
// // }
// // // Map zurücksetzen
// // fTrackHitMap.clear();
// for(const auto& hit : fStepHits) {
// out << hit.eventID << "\t"
// << hit.trackID << "\t"
// << hit.edep << "\t"
// << hit.z << "\t"
// << hit.r << "\n";
// }
// fStepHits.clear();
// }
void ASensitiveDetector::EndOfEvent(G4HCofThisEvent*)
{
G4int tid = G4Threading::G4GetThreadId();
thread_local std::ofstream out;
if(!out.is_open())
{
std::filesystem::create_directories("outfiles");
std::string filename =
fOutputFilename + "_" +
std::to_string(tid) + ".hits";
//out.open(filename, std::ios::app); //anhängen
out.open(filename, std::ios::out); //datei überschreieben
if(!out.is_open())
{
G4cerr << "Cannot open file "
<< filename << G4endl;
return;
}
// Header schreiben
for(size_t i=0;i<fOutputColumns.size();++i)
{
out << fOutputColumns[i];
if(i+1<fOutputColumns.size())
out << "\t";
}
out << "\n";
}
// ===== Daten schreiben =====
for(const auto& [id, hit] : fTrackHitMap)
{
for(size_t i=0;i<fOutputColumns.size();++i)
{
const auto& col = fOutputColumns[i];
if(outputMap.count(col))
outputMap[col](hit,out);
else
out << "NA";
if(i+1<fOutputColumns.size())
out << "\t";
}
out << "\n";
}
fTrackHitMap.clear();
}