Skip to content

Commit 2f5a51a

Browse files
perematoperemato
committed
Added JuliaAction example
1 parent 125aaa9 commit 2f5a51a

File tree

6 files changed

+308
-3
lines changed

6 files changed

+308
-3
lines changed

‎docs/Project.toml

+1
Original file line numberDiff line numberDiff line change
@@ -5,6 +5,7 @@ DisplayAs = "0b91fe84-8a4c-11e9-3e1d-67c38462b6d6"
55
Documenter = "e30172f5-a6a5-5a46-863b-614d45cd2de4"
66
FHist = "68837c9b-b678-4cd5-9925-8a54edc8f695"
77
Geant4 = "559df036-b7a0-42fd-85df-7d5dd9d70f44"
8+
Geant4_jll = "872b6946-528a-5ac7-9145-d37eec569368"
89
GeometryBasics = "5c1252a2-5f33-56bf-86c9-59e7332b4326"
910
IGLWrap_jll = "283677c1-8365-580c-84e5-ef4b5d190868"
1011
IJulia = "7073ff75-c697-5162-941a-fcdaad2a7d2a"

‎docs/make.jl

+3-3
Original file line numberDiff line numberDiff line change
@@ -38,10 +38,10 @@ end
3838

3939
basic_mds = process_literate("B1", "B2a", "B2aVis", "B3a")
4040
extend_mds = process_literate("GPS", "RE03", "TestEm3", "Solids")
41-
advanced_mds = process_literate("TPCSim", "HBC30", "WaterPhantom", "UserLib")
41+
advanced_mds = process_literate("TPCSim", "HBC30", "WaterPhantom", "UserLib", "JuliaAction")
4242
extra_mds = create_extras("B2aDetector.jl", "B2aVisSettings.jl", "B3Detector.jl", "GPSDetector.jl",
43-
"RE03Detector.jl", "TestEm3Detector.jl", "HBC30Detector.jl", "UserLibrary.cpp")
44-
43+
"RE03Detector.jl", "TestEm3Detector.jl", "HBC30Detector.jl", "UserLibrary.cpp",
44+
"MyCode.jl", "G4example.cpp")
4545
examples_mds = []
4646
append!(examples_mds, basic_mds, extend_mds, advanced_mds)
4747

‎docs/src/examples/G4example.cpp

+170
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,170 @@
1+
#include "G4VUserDetectorConstruction.hh"
2+
#include "G4VUserPrimaryGeneratorAction.hh"
3+
#include "G4UserEventAction.hh"
4+
#include "G4UserRunAction.hh"
5+
#include "G4ParticleGun.hh"
6+
#include "G4RunManager.hh"
7+
#include "G4UImanager.hh"
8+
#include "G4RunManagerFactory.hh"
9+
#include "QBBC.hh"
10+
#include "G4NistManager.hh"
11+
#include "G4Box.hh"
12+
#include "G4LogicalVolume.hh"
13+
#include "G4PVPlacement.hh"
14+
//#include "G4SystemOfUnits.hh"
15+
#include "globals.hh"
16+
#include "G4VUserActionInitialization.hh"
17+
#include <julia.h>
18+
#include <iostream>
19+
20+
//---Detector construction class-------------------------------------------------------------------
21+
class DetectorConstruction : public G4VUserDetectorConstruction
22+
{
23+
public:
24+
DetectorConstruction() = default;
25+
~DetectorConstruction() override = default;
26+
G4VPhysicalVolume* Construct() override {
27+
auto nist = G4NistManager::Instance();
28+
auto world_mat = nist->FindOrBuildMaterial("G4_AIR");
29+
auto core_mat = nist->FindOrBuildMaterial("G4_WATER");
30+
auto world_size = 1.0*CLHEP::m;
31+
auto solidWorld = new G4Box("World", world_size, world_size, world_size);
32+
auto logicWorld = new G4LogicalVolume(solidWorld, world_mat, "World");
33+
auto physWorld = new G4PVPlacement(0, G4ThreeVector(), logicWorld, "World", 0, false, 0);
34+
auto core_size = 0.5*CLHEP::m;
35+
auto solidCore = new G4Box("Core", core_size, core_size, core_size);
36+
auto logicCore = new G4LogicalVolume(solidCore, core_mat, "Core");
37+
new G4PVPlacement(0, G4ThreeVector(), logicCore, "Core", logicWorld, false, 0);
38+
return physWorld;
39+
}
40+
};
41+
42+
//---Primary generator action class----------------------------------------------------------------
43+
class PrimaryGeneratorAction : public G4VUserPrimaryGeneratorAction
44+
{
45+
public:
46+
PrimaryGeneratorAction() { fParticleGun = new G4ParticleGun(); }
47+
~PrimaryGeneratorAction() { delete fParticleGun; }
48+
void GeneratePrimaries(G4Event* anEvent) override {
49+
fPrimaryParticleName = fParticleGun->GetParticleDefinition()->GetParticleName();
50+
fParticleGun->SetParticleMomentumDirection(G4ThreeVector(0.,0.,1.));
51+
fParticleGun->SetParticlePosition(G4ThreeVector(0.,0.,-1.*CLHEP::m));
52+
fParticleGun->GeneratePrimaryVertex(anEvent);
53+
}
54+
G4ParticleGun* GetParticleGun() {return fParticleGun;};
55+
const G4String& GetParticleName() { return fPrimaryParticleName;}
56+
private:
57+
G4String fPrimaryParticleName;
58+
G4ParticleGun* fParticleGun;
59+
};
60+
61+
typedef void (*stepaction_f)(const G4Step*);
62+
class SteppingAction : public G4UserSteppingAction
63+
{
64+
public:
65+
SteppingAction() {
66+
action_jl = (stepaction_f)(jl_unbox_voidpointer(jl_eval_string("@cfunction(stepping_action, Nothing, (CxxPtr{G4Step},))")));
67+
std::cout << "=====> " << action_jl << std::endl;
68+
}
69+
~SteppingAction() {}
70+
void UserSteppingAction(const G4Step* step) override { action_jl(step); }
71+
private:
72+
stepaction_f action_jl;
73+
};
74+
75+
typedef void (*eventaction_f)(const G4Event*);
76+
class EventAction : public G4UserEventAction
77+
{
78+
public:
79+
EventAction() {
80+
beginevent_jl = (eventaction_f)(jl_unbox_voidpointer(jl_eval_string("@cfunction(begin_of_event_action, Nothing, (CxxPtr{G4Event},))")));
81+
endevent_jl = (eventaction_f)(jl_unbox_voidpointer(jl_eval_string("@cfunction(end_of_event_action, Nothing, (CxxPtr{G4Event},))")));
82+
}
83+
~EventAction() override = default;
84+
85+
void BeginOfEventAction(const G4Event* event) override { beginevent_jl(event); }
86+
void EndOfEventAction(const G4Event* event) override { endevent_jl(event); }
87+
private:
88+
eventaction_f beginevent_jl;
89+
eventaction_f endevent_jl;
90+
};
91+
92+
typedef void (*runaction_f)(const G4Run*);
93+
class RunAction : public G4UserRunAction
94+
{
95+
public:
96+
RunAction() {
97+
beginrun_jl = (runaction_f)(jl_unbox_voidpointer(jl_eval_string("@cfunction(begin_of_run_action, Nothing, (CxxPtr{G4Run},))")));
98+
endrun_jl = (runaction_f)(jl_unbox_voidpointer(jl_eval_string("@cfunction(end_of_run_action, Nothing, (CxxPtr{G4Run},))")));
99+
}
100+
~RunAction() override = default;
101+
102+
void BeginOfRunAction(const G4Run* run) override { beginrun_jl(run); }
103+
void EndOfRunAction(const G4Run* run) override { endrun_jl(run); }
104+
105+
private:
106+
runaction_f beginrun_jl;
107+
runaction_f endrun_jl;
108+
};
109+
110+
//---Action initialization class-------------------------------------------------------------------
111+
class ActionInitialization : public G4VUserActionInitialization
112+
{
113+
public:
114+
ActionInitialization() = default;
115+
~ActionInitialization() override = default;
116+
void BuildForMaster() const override {}
117+
void Build() const override {
118+
SetUserAction(new PrimaryGeneratorAction);
119+
SetUserAction(new RunAction);
120+
SetUserAction(new EventAction);
121+
SetUserAction(new SteppingAction);
122+
}
123+
};
124+
125+
JULIA_DEFINE_FAST_TLS // only define this once, in an executable (not in a shared library) if you want fast code.
126+
127+
//----Main program---------------------------------------------------------------------------------
128+
int main(int, char**)
129+
{
130+
//--- Required to setup the Julia context
131+
jl_init();
132+
/* run Julia commands */
133+
jl_eval_string("include(\"MyCode.jl\")");
134+
if (jl_exception_occurred()) {
135+
std::cout << "=====> " << jl_typeof_str(jl_exception_occurred()) << std::endl;
136+
}
137+
138+
//---Construct the default run manager
139+
auto runManager = G4RunManagerFactory::CreateRunManager(G4RunManagerType::Serial);
140+
141+
//---Set mandatory initialization classes
142+
// Detector construction
143+
runManager->SetUserInitialization(new DetectorConstruction());
144+
145+
// Physics list
146+
runManager->SetUserInitialization(new QBBC(0));
147+
148+
// User action initialization
149+
runManager->SetUserInitialization(new ActionInitialization());
150+
151+
// Initialize G4 kernel
152+
runManager->Initialize();
153+
154+
// Get the pointer to the User Interface manager
155+
auto UImanager = G4UImanager::GetUIpointer();
156+
UImanager->ApplyCommand("/control/verbose 1");
157+
UImanager->ApplyCommand("/run/verbose 1");
158+
//UImanager->ApplyCommand("/event/verbose 0");
159+
//UImanager->ApplyCommand("/tracking/verbose 1");
160+
UImanager->ApplyCommand("/gun/particle e+");
161+
UImanager->ApplyCommand("/gun/energy 100 MeV");
162+
UImanager->ApplyCommand("/run/beamOn 100000");
163+
164+
// Job termination
165+
delete runManager;
166+
167+
// strongly recommended: notify Julia that the program is about to terminate.
168+
// this allows Julia time to cleanup pending write requests and run all finalizers
169+
jl_atexit_hook(0);
170+
}

‎docs/src/examples/JuliaAction_lit.jl

+94
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,94 @@
1+
# # Calling G4 actions in Julia
2+
#
3+
# This is a very simple example of calling user actions in Julia from a C++
4+
# Geant4 application.
5+
# We define the user actions in Julia language in the file [`MyCode.jl`](@ref)
6+
# and call them from the C++ application. The name and signatures of the functions
7+
# are important since the C++ will associate them in the corresponding inherited
8+
# classes.
9+
#
10+
# The C++ code is a single file [`G4example.cpp`](@ref) that defines the Geant4 the minimal
11+
# set of classes to run a simulation.
12+
# - The main program is responsible of initializing Julia by calling `julia_init` and
13+
# loading the Julia code executing.
14+
# ```cpp
15+
# jl_init()
16+
# jl_eval_string("include(\"MyCode.jl\")");
17+
# ```
18+
# - Each constructor of a user action class needs to initialize a C++ function pointer to the
19+
# corresponding Julia function. This is done in the constructor to avoid any dynamic dispatch
20+
# at runtime. For example, for the `EventAction` class:
21+
# ```cpp
22+
# typedef void (*eventaction_f)(const G4Event*);
23+
# class EventAction : public G4UserEventAction {
24+
# public:
25+
# EventAction() {
26+
# beginevent_jl = (eventaction_f)(jl_unbox_voidpointer(jl_eval_string("@cfunction(begin_of_event_action, Nothing, (CxxPtr{G4Event},))")));
27+
# endevent_jl = (eventaction_f)(jl_unbox_voidpointer(jl_eval_string("@cfunction(end_of_event_action, Nothing, (CxxPtr{G4Event},))")));
28+
# }
29+
# ...
30+
# private:
31+
# eventaction_f beginevent_jl;
32+
# eventaction_f endevent_jl;
33+
# };
34+
# ```
35+
# - Finally the actions are called in the corresponding Geant4 classes. For example in the `EventAction` class:
36+
# ```cpp
37+
# void EventAction::BeginOfEventAction(const G4Event* event) {
38+
# beginevent_jl(event);
39+
# }
40+
# ...
41+
# ```
42+
43+
#md # !!! note "Note that"
44+
#md # You can also download this example as a
45+
#md # [Jupyter notebook](JuliaAction.ipynb) and a plain
46+
#md # [Julia source file](JuliaAction.jl).
47+
#md #
48+
#md # The C++ code is available as a [source file](G4example.cpp) and the
49+
#md # Julia code is available as a [source file](MyCode.jl).
50+
#md #
51+
#md # #### Table of contents
52+
#md # ```@contents
53+
#md # Pages = ["JuliaAction.md"]
54+
#md # Depth = 2:3
55+
#md # ```
56+
57+
# ## Loading the necessary Julia modules
58+
using Geant4_jll # Needed to locate the Geant4 installation directory
59+
#md import DisplayAs: PNG #hide
60+
61+
# ## Building G4Example Application
62+
# The custom library is defined in the C++ file `G4example.cpp`. It is a single file to
63+
# facilitate the building of the executable.
64+
#
65+
# The attribute `Geant4_jll.artifact_dir` provides the path to the Geant4 installation directory.
66+
# Sources are in the same location as this script.
67+
cd(@__DIR__)
68+
g4prefix = Geant4_jll.artifact_dir
69+
jlprefix = dirname(Sys.BINDIR);
70+
71+
# We use the executables `geant4-config` and `julia-config.jl` to get the needed
72+
# libraries and compiler/linker flags.
73+
74+
g4libs = read(`$g4prefix/bin/geant4-config --libs`, String) |> split
75+
filter!(x -> x != "-lG4gdml", g4libs)
76+
jllibs = read(`$jlprefix/share/julia/julia-config.jl --ldlibs`, String) |> split
77+
append!(jllibs, ["-L$jlprefix/lib"])
78+
cflags = read(`$g4prefix/bin/geant4-config --cflags`, String) |> split
79+
ldflags = ["-Wl,-rpath,$g4prefix/lib", "-Wl,-rpath,$jlprefix/lib"];
80+
Sys.KERNEL == :Linux && append!(ldflags, ["-Wl,--no-as-needed"]);
81+
82+
# Run the compilation and link command
83+
Base.run(`c++ -O2 -fPIC $cflags -I$jlprefix/include/julia $ldflags $g4libs $jllibs
84+
-o G4example $(@__DIR__)/G4example.cpp`).exitcode == 0 || error("Compilation failed");
85+
86+
# ## Run the application
87+
# We need to set the variable `JULIA_PROJECT` pointing to correctly setup Julia environment.
88+
withenv("JULIA_PROJECT" => "@.") do
89+
Base.run(`./G4example`).exitcode == 0 || error("Execution failed");
90+
end
91+
92+
# ## Display the results
93+
println("=====> The file edepHist.png should have been saved")
94+
#md # ![](edepHist.png)

‎docs/src/examples/MyCode.jl

+39
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,39 @@
1+
#---Load the needed Julia modules------------------------------------------------------------------
2+
3+
using Geant4
4+
using GeometryBasics
5+
using FHist
6+
using Plots
7+
8+
println("=====> Loading MyCode.jl")
9+
10+
edepHist = H1D("Event total Edep distribution", 100, 0., 110.)
11+
edep = 0.0
12+
13+
function end_of_event_action(event)
14+
push!(edepHist, edep)
15+
return # This is mandatory to force to return nothing
16+
end
17+
18+
function begin_of_event_action(event)
19+
global edep = 0.0
20+
return # This is mandatory to force to return nothing
21+
end
22+
23+
function stepping_action(step)
24+
global edep += GetTotalEnergyDeposit(step)
25+
return # This is mandatory to force to return nothing
26+
end
27+
28+
function begin_of_run_action(run)
29+
return # This is mandatory to force to return nothing
30+
end
31+
32+
function end_of_run_action(run)
33+
println("=====> End of run")
34+
h = edepHist
35+
img = plot(h.hist, title=h.title)
36+
savefig(img, "edepHist.png")
37+
println("=====> edepHist.png saved")
38+
return # This is mandatory to force to return nothing
39+
end

‎docs/src/releases.md

+1
Original file line numberDiff line numberDiff line change
@@ -2,6 +2,7 @@
22
### 0.2.x
33
- New Examples:
44
- UserLib: how to build and call a user custom library providing additional Geant4 functionally that is not provided by the set of wrapped classes
5+
- JuliaAction: emmbeding Julia in a C++ application. In this example we call user actions implemented in Julia
56

67
### 0.2.0
78
- New Features

0 commit comments

Comments
 (0)