GAMOS an easy and flexible framework for GEANT4 simulations Pedro Arce Pedro Rato Mendes

1. GAMOS an easy and flexible framework for GEANT4 simulations Pedro Arce Pedro Rato Mendes Juan Ignacio Lagares CIEMAT, Madrid 2008 NSS-MIC-RTSD workshop Dresde (D), 19-25 October 2008

2. Outline • Debugging and optimisation • Extracting information • Time studies • Variance reduction • Setting cuts automatically • Applications • PET • Radiotherapy • Accelerator • Dose in phantom • Documentation • Summary • Introduction • GEANT4-based frameworks • GAMOS objectives • GAMOS functionality • Geometry • Movements • Generator • Physics • User actions • Sensitive detector and hits • Histograms • Parameter management • Scoring

3. Introduction

4. GEANT4-based frameworks • GEANT4 is a very powerful and flexible toolkit • Well-proven physics for medical applications • Powerful geometry and visualisation tools • Modern language and object-oriented technologies • Easy to extract information and modify running conditions • But: it is not easy to use • Everything requires writing C++, and a good knowledge of GEANT4 details • Several tools try to make GEANT4 easy to use in an specific domain: • Providing an scripting language so that no C++ is needed • Trying to cover with it all the requirements of a user in the field

5. GEANT4-based frameworks • But often these tools become a quite rigid framework • Users can only do what the authors thought they would want to do Most of the GEANT4 users nowadays are researchers • They always want to have a deep understanding of what happens in the simulation • Length travelled by electrons produced by compton interactions in a crystal... • Angle of bremsstrahlung photons when original electron > 1 MeV... • Time spent in each volume… • … • They want to try new things, that the framework authors have not imagined before • Non-standard geometries… • Uncommon source position/angle/energy/time distributions… • Use a different sensitive detector with some peculiarities… • Particular trigger or DAQ conditions… • …. • A framework should be easy to use and flexible

6. GAMOS objectives The first objective of GAMOS is to allow the user to: • Simulate fully a medical physics project with a minimal knowledge of GEANT4 and no need of C++ (only text scripts) The second objective is to allow the user to: • Easilyadd new functionality and combine it with the existing functionality in GAMOS • And also reuse any GEANT4 C++ code (from the user, from GEANT4 examples, …) • Transform dinamically C++ code into user commands + detailed documentation + step-by-step tutorials + templates on how to do your own plug-in  GAMOS is based on the plug-in technology

7. GAMOS plug-in´s • Most users will be able to do their simulation with simple user commands • If they want to have something not foreseen by the framework  use plug-in’s • GAMOS has no predefined components • At run-time user selects which components to load through user commands • They are loaded through the plug-in mechanism • User has full freedom in choosing components • User can define a component not foreseen by GAMOS • Take a C++ class and use it through an user command • Mix it with any other component For example: instead of GAMOS generator use the one from G4 example underground_physics DEFINE_GAMOS_GENERATOR(DMXPrimaryGeneratorAction); and then you can select it in your script /gamos/generator DMXPrimaryGeneratorAction For the plug-in's implementation in GAMOS it has been chosen the CERN library: SEAL

8. GAMOS functionality

9. Geometry Three different ways to define: C++ code: • The usual GEANT4 way • Add one line to transform your class in a plug-in DEFINE_GAMOS_GEOMETRY (MyGeometry); so that you can select it in your input macro /gamos/geometry MyGeometry Define it in ASCII (text) files • The easiest way to define a geometry • Based on simple tags • Same order of parameters as corresponding GEANT4 classes Using one of the GAMOS examples for simple cases • Simple PET can be defined through an 8-parameters file (n_crystals, crystal_x/y/z, radius, …) • Simple phantom with a few user commands

10. Geometry from text files • Based on simple tags, with the same order of parameters as corresponding GEANT4 classes :MATE Cu 29 63.54 8.9333*g/cm3 :VOLU yoke TUBE  62.*cm 820. 1.27*m Cu :PLACE  yoke  1   expHall  R0   0.0   0.0   370.*cm MATERIALS: • Isotopes • Elements • Simple materials • Material mixtures by weight, volume or number of atoms • GEANT4 intrinsic materials SOLIDS: • All GEANT4 “CSG” and “specific” solids • Twisted solids • Tesellated solids • Boolean solids

11. Geometry from text files ROTATIONMATRICES: • 3 rotation angles around X,Y,Z • 6 theta and phi angles of X,Y,Z axis • 9 matrix values PLACEMENTS: • Simple placements • Divisions • Replicas • Parameterisations • Linear, circular or square • For complicated parameterisations example of how to mix the C++ parameterisation with the ASCII geometry file COLOUR VISUALISATION ON/OFF

12. Geometry from text files PARAMETERS: • Can be defined to use them later :P InnerR 12. :VOLU yoke :TUBS  Iron   3  $InnerR 820. 1270. ARITHMETIC EXPRESSIONS: • Elaborated expressions can be used :SOLID yoke TUBE sin($ANGX)*2+4*exp(1.5) 820. 1270. UNITS: • Default units for each parameter • Each vaule can be overridden by user INCLUDE OTHER FILES (hierarchical approach): #include mygeom2.txt.

13. Geometry from text files CMS detector • User can extend it: add new tags and process it without touching base code • Can mix a C++ geometry with a text geometry • GEANT4 in memory geometry  text files • Install and use it as another GEANT4 library G4VPhysicalVolume* MyDetectorConstruction::Construct(){ G4tgbVolumeMgr* volmgr = G4tgbVolumeMgr::GetInstance(); volmgr->AddTextFile(filename); // SEVERAL FILES CAN BE ADDED return = volmgr->ReadAndConstructDetector(); HISTORY: • In use to build GEANT4 geometries since 10 years ago • An evolving code… • It will be part of the new GEANT4 release in December ‘08

14. Movements User can move a volume by using commands • Displacements or rotations • Every N events or every interval of time • N times or forever • Offset can be defined • Several movements of the same volume or different ones can be set

15. Generator • C++ code • The usual GEANT4 way • Add one line to transform your class in a plug-in • DEFINE_GAMOS_GENERATOR(MyGenerator); • so that you can select it in your input macro • /gamos/generator MyGenerator • GAMOS generator • Combine any number of single particles or isotopes decaying to e+, e-, g • For each particle or isotope user may select by user commands any combination of time, energy, position and direction distributions • Or create its own and select it by a user command (transforming it into a plug-in)

16. Generator distributions • POSITION: • Point Square Rectangle • Disc DiscGaussian LineSteps • InG4Volumes (sphere/box/tube_section/ellipsoid) • InG4VolumeSurfaces (sphere/box/tube_section/ellipsoid) • InG4VolumesGeneral (any solid) VoxelPhantomMaterials • DIRECTION: • Random Const Cone • ENERGY: • Constant Gaussian RandomFlat • BetaDecay ConstantIsotopeDecay • TIME: • Constant Decay • POSITION & DIRECTION: • InVolumeSurfaceTowardsCentre TowardsBox • A template provided for each type: user only has to write a few C++ lines to create her/his own distribution F18 decay

17. Physics • C++ code • The usual GEANT4 way • Add one line to transform your class in a plug-in • DEFINE_GAMOS_PHYSICS_LIST (ExN02PhysicsList); • so that you can select it in your input macro • /gamos/physicsList ExN02PhysicsList • Any GEANT4 physics list can be easily selected through a user command • GAMOS physics list • Based on hadrotherapy advanced example • User can combine different physics lists for photons, electrons, positrons, muons, protons and ions • A dummy one for visualisation

18. User actions • They are the main way for a user to extract information of what is happening and modify the running conditions: flexible framework is needed • User can have as many user actions of any type as she/he wants • Only 1 of each type in GEANT4 • One user action can be of several types (run/event/stacking/tracking/stepping) • Only 1 class - 1 type in GEANT4 • User actions can use filters and indexers(= classifiers) /gamos/userAction GmTrackHistosUA GmGammaFilter • User can activate any user action by a user command • User can create new user actions • Just adding a line to transform it into a plug-in DEFINE_GAMOS_USER_ACTION(MyUserAction); that can be selected with a user command /gamos/userAction MyUserAction

19. Sensitive Detectors • To produce hits in GEANT4 a user has to (with C++): • Define a class inheriting from G4VSensitiveDetector • Associate it to a G4LogicalVolume • Create hits in the ProcessHits method • Clean the hits at EndOfEvent • In GAMOS you can do all this with a user command • /gamos/assocSD2LogVol SD_CLASS SD_TYPE LOGVOL_NAME • SD_CLASS: the G4VSensitiveDetector class • SD_TYPE: an identifier string, so that different SD/hits can have different treatment • User can create his/her own SD class and make it a plugin • DEFINE_GAMOS_SENSDET(MySD);

20. Hits • A GAMOS hit has the following information • G4int theDetUnitID; ID of the sensitive volume copy • G4int theEventID; • G4double theEnergy; • G4double theTimeMin; time of the first E deposit • G4double theTimeMax; time of the last E deposit • G4ThreeVector thePosition; • std::set<G4int> theTrackIDs; list of all tracks that contributed • std::set<G4int> thePrimaryTrackIDs; list of all ‘primary´tracks that contributed • std::vector<GamosEDepo*> theEDepos; list of all deposited energies • G4String theSDType; • User can create his/her own hit class • Hits can be stored in a job and read in another job • Format compatible with STIR (PET reconstruction)

21. Hits Bi fluorescence peak 511 keV photoelectric (PE) peak backscattering peak escape peak E > 511 keV photoelectric peaks /gamos/userAction GmHitsHistosUA Hits energy spectra in BGO crystals

22. Hits reconstruction • Framework to transform simulated hits digital signals reconstructed hits (energy/time/position) • Digitization and hits reconstruction are very detector specific : it is not possible to provide a general solution • GAMOS just provide a few simple digitizers and hits reconstructor s • 1 hit  1 digit • Merge hits close enough • Same set of sensitive volumes • Closer than a given distance • … and a basic structure • Hits compatible in time • spanning various events • Trigger • Pulse simulation • Sampling • Noise -- hits -- clusters of hits

23. Detector effects • Measuring time • - A detector is not able to separate signals from different events if they come close in time • Dead time • When a detector is triggered, this detector (or even the whole group it belongs to) is not able to take data during some time • Paralizable or non-paralizable • Both can be set by the user in the input macro • A different time for each SD_TYPE • /gamos/setParam SD:MeasuringTime:MySD_1 1000.*ns • /gamos/setParam SD:MeasuringTime:MySD_2 250.*ns

24. Histograms • Own format: can write histograms in text files (CSV) without any external dependency • Read them in Excel, Matlab, Origin, … • Also ROOT histograms supported • AIDA will be supported soon /gamos/userAction GmGenerHistosUA /gamos/analysis/fileFormat CSV MS Excel graphics of energy of e+ from F18 decay

25. Histograms • Same code to create and fill histograms independent of the format • GAMOS takes care of writing the file in the chosen format at the end of job • GmAnalysisMgrkeeps a list of histograms so that they can be accessed from any part of the code, by number or name • GmHitsEventMgr::GetInstance(“pet”)->GetHisto1(1234)->Fill(ener); • GmHitsEventMgr::GetInstance(“pet”)->GetHisto1(“CalorSD: hits energy”)->Fill(ener); • There can be several files, each one with its own histograms • When creating an histogram, user chooses file name

26. Parameter management • Many classes change their behaviour depending on the value of some parameters that the user can set • GAMOS utilities use parameters extensively to maximise the flexibility • Always a default value in case the user does not care • Parameters can be numbers, strings, list of numbers or list of strings • GAMOS provides a unique command to manage all these parameters • A parameter is defined in the input macro • /gamos/setParam SD:EnergyResol 0.1 • Any class can use this value in any part of the code • float enerResol = GmParameterMgr::GetInstance() • ->GetNumericValue(“SD:EnergyResol”,0.);

27. Scoring • Scoring is an important part of your simulation  powerful and flexible framework developed: • Many possible quantities can be scored in one or several volumes • For each scored quantity one of several filters can be used • only electrons, only particles in a given volume, … • Several ways to classify the different scores • One different score for each different volume copy, or volume name, or energy bin, … • Results can be printed in one or several formats for each scored quantity • All scored quantities can be calculated with/without errors • All scored quantities can be calculated per event or per job • Taking into account correlations from particles from same event

28. Scoring Associate one/several scorer(s) to any logical volume(s): CellCharge CellFlux DoseDeposit EnergyDeposit FlatSurfaceCurrent FlatSurfaceFlux MinKinEAtGeneration NofCollision NofSecondary NofStep PassageCellCurrent PassageCellFlux PassageTrackLength Population SphereSurfaceCurrent SphereSurfaceFlux TrackCounter TrackLength Associate one/several filter(s) to each scorer Gamma Electron Positron ElectronOrPositron EMParticle Particle Charged Neutral Primary Secondary KineticEnergy Region LogicalVolume PhysicalVolume Associate one/several printer(s) to each scorer PrinterDefault Printer3ddose PrinterSqdose RTDoseHistos Associate one/several indexer(s) to each scorer By1Ancestor ByAncestors ByKineticEnergy ByLogicalVolume ByPhysicalVolume ByRegion

29. Debugging and Optimisation

30. Debugging and optimisation • We have made a substantial effort to provide tools to help users in • Knowing all the details of what is happening in the simulation • Optimising GEANT4 for her/his application • Flexible user actions • Flexible and powerful scoring • Many control histograms • Detailed verbosity control • Detailed CPU time study • Profiling your simulation • Setting cuts by regions through user commands • Automatic optimisation of cuts • Variance reduction techniques

31. Extracting information • Flexible user actions • Flexible and powerful scoring • Many control histograms available through user commands • Source particles, hits, reconstructed hits, track information, PET classification, phase space, dose, … • Controlled by filters and indexers • Verbosity controlled through user commands • Different verbosities for different simulation parts (geometry, physics, generator, …) • 6 levels of verbosity (silent, error, warning, info, debug, test) • Verbosity managers are plug-ins: user can easily create its own one • GEANT4 tracking verbosity can be controlled by event or by track /gamos/setParam GmTrackingVerboseUA:EventMin 14632 /gamos/setparam GmTrackingVerboseUA:TrackMin 10 /gamos/setparam GmTrackingVerboseUA:TrackMax 20

32. Optimisation: time studies • User commands to get CPU time study • By particle, energy bins, volume, region (or combination of them) Just add a user command: /gamos/userAction GmTimeStudyUA GmIndexerByParticle GmIndexerByEnergy • gamma/0.001-0.01: User=0.01 Real=0 Sys=0 • gamma/0.01-0.1: User=2.01 Real=2.45 Sys=0.27 • gamma/0.1-1: User=19.12 Real=22.05 Sys=1.51 • gamma/1-10: User=4.25 Real=5.4 Sys=0.3 • e-/0.0001-0.001: User=0.07 Real=0.1 Sys=0 • e-/0.001-0.01: User=0.54 Real=0.69 Sys=0.06 • e-/0.01-0.1: User=4.71 Real=5.41 Sys=0.38 • e-/0.1-1: User=15.59 Real=18.19 Sys=1.79 • e-/1-10: User=82.83 Real=98.62 Sys=7.45 • Example to get detailed gprof profiling about where (in which methods) the time is spent • Time spent in a method and integrated time in all the methods called by it

33. Optimisation: setting cuts • All GEANT4 commands are available, for example to control cuts and electromagnetic physics parameters • Cuts by region can be defined through user commands or in the TEXT geometry file • Automatic conversion energy cut  range cut Just add two user commands: /gamos/userAction GmCutsEnergy2RangeUA /gamos/ECuts2RangeCuts * 10.*keV gamma • MATERIAL: G4_AIR PART: gamma ENERGY CUT: 0.01 (MeV) = RANGE CUT: 17407.7 (mm) • MATERIAL: G4_Cu PART: gamma ENERGY CUT: 0.01 (MeV) = RANGE CUT: 0.185575 (mm) • MATERIAL: G4_W PART: gamma ENERGY CUT: 0.01 (MeV) = RANGE CUT: 0.0149411 (mm) • MATERIAL: G4_KAPTON PART: gamma ENERGY CUT: 0.01 (MeV) = RANGE CUT: 21.5443 (mm) • User limits can be attached to one/several logical volume(s) through user commands • MaxStep, MaxTrkLength, MaxTOF, MinKinE, MinRange

34. How to get best production cuts/user limits? • Normally you are worried that cuts do not change your results more than few % (or less) • need to run very big statistics • Cuts can be very different for different particle/regions, and a cut in one particle/region may affect another • need to run many sets of values RESULT: usually people think it is too complicated • optimise cuts quickly in an approximate way • or take somebody’s cuts, not really optimised for her/his setup In GAMOS you can get the optimal production cuts and user limits in a single job with limited statistics • We use an inverse reasoning: • Information of all tracks/steps that contribute to your results is stored • At the end of job results are analysed to get automatically best cut values • It is a indeed a complicated technique: details in GAMOS user’s guide

35. RESULTS: Optimising cuts 1 single job to optimise production cuts for a radiotherapy accelerator simulation: Number of particles that would not reach the patient for each cut group 1 single job to optimise minimum range cuts for a dose in phantom simulation: Dose in patient lost for each cut group (Lost dose distribution histograms also provided)

36. Optimisation: variance reduction • Two bremsstrahlung splitting techniques are available through user commands • Uniform bremsstrahlung splitting • Z-plane direction bremsstrahlung splitting • A third one (inspired on EGSnrc DBS) under development • Other techniques are under study

37. Applications

38. Applications: PET • Any PET detector can be simulated with simple text format • Several sensitive detectors types selectable by user commands • User can create it own one and make it a plug-in • Automatic creation of hits with detailed information • Writing hits into text or binary file • Retrieving them in next job • Format compatible with STIR software • Digitizer and reconstructed hits framework and examples • Histos of hits or reconstructed hits with a user command

39. Applications: PET • Detector effects • Energy resolution • Position resolution • Time resolution • Dead time • Measuring time • CTI ECAT Exact (922) & GE Discovery ST (PET/CT) detectors have been simulated • Results agree with data (it is GEANT4…) • ClearPET is being simulated with GAMOS • BrainPET is being simulated with GAMOS

40. Applications: radiotherapy Accelerator simulation • Any accelerator (with any kind of MLCs) can be simulated with simple text format • Write phase space files in IAEA format • Several Z planes in the same job • Stop after last Z plane or not • Save header after N events (not to lose everything if job is aborted) • Save info of regions particle traversed • Kill particles at big XY • Automatic of user-defined limits • Optional histograms by particle type at each Z plane

41. Applications: radiotherapy Dose in phantom simulation: • Reading IAEA phase space files • Displace or rotate phase space particles • Reuse phase space particles • Optional automatic calculation of reuse number • Optional mirroring in X & Y • Recycle phase space files • Optional skipping of first N events • Optional histograms of particles read

42. Applications: radiotherapy Dose in phantom simulation: • Building phantom geometry • Use DICOM files • Read GEANT4 format (example advanced/medical/DICOM: DICOM  ASCII) • Read EGS format • Split one material in several if voxel densities are different • User defines an interval X and a material is created for the densities in each X interval • Displace or rotate read-in phantom geometry • PTV/CTV/GTV management • Simple phantoms can be created with a few commands • Define voxel limits • Define number of voxels • Define material for each Z plane

43. Applications: radiotherapy Dose in phantom simulation: • Dose in phantom • Fast navigation (see POSTER N02-134) • Different printing formats in each run • Text in standard output • Dose histograms (X, Y, Z, XY, XZ, YZ, dose, dose-volume) • User can define several XY, XZ, YZ planes (= simple dose visualisation) • File with dose and error at each voxel • File with dose and dose2 at each voxel (to properly take into account correlations) • Dose in total or by event • Proper management of original number of events when reading phase space files

44. Applications: radiotherapy gMocren visualisation of DICOM geometry and tracks gMocren view of GEANT4 geometry and tracks transported with the new algorithm gMocren view of GEANT4 geometry and tracks transported with the new algorithm gMocren view of GEANT4 geometry and tracks transported with the new algorithm gMocren view of GEANT4 geometry and tracks transported with the new algorithm gMocren view of GEANT4 geometry and tracks transported with the new algorithm

45. Applications: radiotherapy • Analysis of results: A set of easy-to-use executables and ROOT scripts Phase space files • Sum phase space files • Analyse phase space files • Statistics • Histograms by particle type at each Z plane (= when writing PS) • Compare two phase space files Dose files • Sum dose files • Analyse dose files • Statistics • Dose histograms (X, Y, Z, XY, XZ, YZ, dose, dose-volume) (= when writing dose) • Compare two dose files • GUI under development (MIRAS project)

46. Applications: radiotherapy optimisation • GEANT4 electromagnetic parameters are not optimised for radiotherapy • They cover a wide range of energy and applications  they have to be conservative • We have started to optimise GEANT4 electromagnetic parameters for radiotherapy: http://fismed.ciemat.es/GAMOS/RToptim VARIAN 2100 gamma accelerator: 106 events on Pentium Dual-Core 3 GHz * Same parameters as for 6MeV Dose in 104 5x5x3 mm water phantom: 106 events on Pentium Dual-Core 3 GHz Dose in 4.5 106 patient, 23 materials: 106 events on Pentium Dual-Core 3 GHz

47. Documentation

48. Documentation • User’s Guide: • Installation • Documentation of all available functionality • Documentation on how to provide new functionality by creating a plug-in • Examples: • A simple one and a few more complicated ones • /gamos/setParam GmGeometryFromText:FileName mygeom.txt • /gamos/geometry GmGeometryFromText • /gamos/physics GmEMPhysics • /gamos/generator GmGenerator • /run/initialize • /gamos/generator/addIsotopeSource myF18 F18 1.E6*becquerel • /run/beamOn 10

49. Documentation • Tutorials: • Three tutorials • PET tutorial • Radiotherapy tutorial • plug-in tutorial • Propose about 10 exercises each • Increasing in difficulty • Reference output provided • Solutions provided • 4 GAMOS tutorial courses have been given • About 70 attendees

50. Summary • GAMOS is GEANT4-based framework • The computations are made by GEANT4, it only provides a simple user interface • GAMOS user-friendly and flexible • many utilities that allow to do full GEANT4 simulation through user commands • plug-in’s allow to extend functionality by converting C++ classes into user commands by adding one line • flexible and powerful framework to extract detailed information • several tools to optimise CPU performance • GAMOS core is application independent • Several full medical applications are being built on top of GAMOS core http://fismed.ciemat.es/GAMOS