Hi everyone,
I'm currently working on a Geant4 simulation where I have a low-energy visible light beam source that immediately gets absorbed by any material it encounters, and it seems to only move in a vacuum. When I increase the beam energy to around 1.9 MeV, it can pass through materials without any issues; however, with a lower energy (such as 1.99975 eV), the beam doesn't penetrate any materials at all.
I've posted the relevant code below, where I set up a simple environment with a beam and some scoring, and I suspect the issue might be related to the physics list or Geant4’s capabilities to run simulations on low energy beams. I have enabled optical properties for Air and Water, and I am fairly sure there is an issue with the physics of low energy beams.
Does anyone have insights on why this is happening at low energy levels and how I might adjust the setup to allow the beam to interact more realistically with materials at lower energies? I have tried implementing multiple different materials and all of them have interacted the same way. Below I have posted the code and an image of what is happening to the beam. Please ignore the artifacts on the axis. 🙂
`
d:Ge/World/HLX = 20 cm
d:Ge/World/HLY = 20 cm
d:Ge/World/HLZ = 20 cm
s:Ge/World/Material = "Air"
i:Ts/NumberOfThreads = 4
##############################
Graphics
##############################
s:Gr/view/Parent = "World"
s:Gr/view/Type = "OpenGl"
d:Gr/view/Theta = 45 deg
d:Gr/view/Phi = 45 deg
u:Gr/view/Zoom = 2
b:Gr/view/IncludeAxes = "True"
d:Gr/view/AxesSize = 3.2 cm
b:Gr/view/Active = "True" # defaults to "True"
s:Gr/view/CheckForOverlaps = "False"
b:Ge/QuitIfOverlapDetected = "False"
#Color Values
iv:Gr/Color/White1 = 4 255 255 255 25
iv:Gr/Color/DarkBlue = 4 0 0 255 80
b:Ts/QuitIfManyHistoriesSeemAnomalous = "False"
##############################
To End
##############################
b:Ts/PauseBeforeQuit = "True"
##############################
Materials
##############################
b:Ma/Air/EnableOpticalProperties = "True"
dv:Ma/Air/RefractiveIndex/Energies = 7 0.5 1.0 1.5 2.0 2.5 3.0 3.5 eV
uv:Ma/Air/RefractiveIndex/Values = 7 1.00 1.00 1.00 1.00 1.00 1.00 1.00
b:Ma/G4_WATER/EnableOpticalProperties = "True"
##############################
Physics List
##############################
s😛h/ListName = "Optical"
s😛h/Optical/Type = "Geant4_Modular"
sv😛h/Optical/Modules = 2 "g4em-standard_opt3" "g4optical"
sv😛h/Default/Modules = 1 "g4em-standard_opt0"
###################################
Sequence
###################################
i:Ts/SequenceVerbosity = 1
b:Ts/ShowCPUTime = "True"
i:Ts/ShowHistoryCountAtInterval = 10
##############################
Beam Source
##############################
s:Ge/BeamPosition/Type = "TsBox"
s:Ge/BeamPosition/Parent = "World"
s:Ge/BeamPosition/Material = "Vacuum"
d:Ge/BeamPosition/TransX = 0.0 cm
d:Ge/BeamPosition/TransY = 9.0 cm
d:Ge/BeamPosition/TransZ = 0.0 cm
d:Ge/BeamPosition/RotX = -90 deg
d:Ge/BeamPosition/RotY = 0.0 deg
d:Ge/BeamPosition/RotZ = 0.0 deg
d:Ge/BeamPosition/HLX = 1. cm # Half Length
d:Ge/BeamPosition/HLY = 1. cm
d:Ge/BeamPosition/HLZ = 1. cm
b:Ge/BeamPosition/Invisible = "False"
s:So/redLight/Type = "Beam" # Beam, Isotropic, Emittance or PhaseSpace
s:So/redLight/Component = "BeamPosition"
#s:So/redLight/Parent = "World"
s:So/redLight/BeamParticle = "22"
d:So/redLight/BeamEnergy = 1.99975 eV
#u:So/redLight/BeamEnergySpread = 0.757504 # dont want beam energy spread
s:So/redLight/BeamPositionDistribution = "Flat" # None, Flat or Gaussian
s:So/redLight/BeamPositionCutoffShape = "Rectangle" # Rectangle or Ellipse (if Flat or Gaussian)
d:So/redLight/BeamPositionCutoffX = 1. cm # X extent of position (if Flat or Gaussian)
d:So/redLight/BeamPositionCutoffY = 1. cm # Y extent of position (if Flat or Gaussian)
d:So/redLight/BeamPositionSpreadX = 0.65 cm # distribution (if Gaussian) # no smoothing need exact
d:So/redLight/BeamPositionSpreadY = 0.65 cm # distribution (if Gaussian)
s:So/redLight/BeamAngularDistribution = "Flat" # None, Flat or Gaussian
d:So/redLight/BeamAngularCutoffX = 22.57 deg # X cutoff of angular distrib (if Flat or Gaussian)
d:So/redLight/BeamAngularCutoffZ = 22.57 deg
d:So/redLight/BeamAngularCutoffY = 22.57 deg # Y cutoff of angular distrib (if Flat or Gaussian)
#maybe need beam polarization
d:So/redLight/BeamAngularSpreadX = 0.0032 rad # X angular distribution (if Gaussian)
d:So/redLight/BeamAngularSpreadY = 0.0032 rad # Y angular distribution (if Gaussian)
i:So/redLight/NumberOfHistoriesInRun = 1000
s:Ge/Volume/Type = "TsBox"
s:Ge/Volume/Parent = "World"
s:Ge/Volume/Material = "G4_WATER"
d:Ge/Volume/TransX = 0. mm
d:Ge/Volume/TransY = 0. mm
d:Ge/Volume/TransZ = 0. mm
d:Ge/Volume/RotX = 0. deg
d:Ge/Volume/RotY = 0. deg
d:Ge/Volume/RotZ = 0. deg
d:Ge/Volume/HLX = 5. cm
d:Ge/Volume/HLY = 5. cm
d:Ge/Volume/HLZ = 5. cm
s:Ge/Volume/DrawingStyle = "Wireframe"
s:Ge/Volume/Color = "DarkBlue"
##############################
Scoring
##############################
s:Sc/Scorer1/Quantity = "OpticalPhotonCount"
s:Sc/Scorer1/Component = "Volume"
s:Sc/Scorer1/OutputFile = "Volume_OPC"
s:Sc/Scorer1/OutputType = "ASCII"
b:Sc/Scorer1/OutputtoConsole = "True"
s:Sc/Scorer1/IfOutputFileAlreadyExists = "Overwrite"
#i:Sc/Scorer1/BounceLimit = 1000000
s:Sc/Scorer2/Quantity = "EnergyDeposit"
s:Sc/Scorer2/Component = "Volume"
i:Ma/Verbosity = 1
b:Sc/Scorer2/OutputtoConsole = "True"
s:Sc/Scorer2/Material = "G4_WATER"
i:Sc/Scorer2/EBins = 20 # defaults to 1, that is, un-binned
d:Sc/Scorer2/EBinMin = 0.0 eV # defaults to zero
d:Sc/Scorer2/EBinMax = 5.0 eV # must be specified if EBins is greater than 1
s:Sc/Scorer2/EBinEnergy = "PreStep" # "IncidentTrack", "PreStep" or "DepositedInStep"
s:Sc/Scorer2/OutputFile = "DoseToWater_ED"
s:Sc/Scorer2/OutputType = "CSV"
s:Sc/Scorer2/IfOutputFileAlreadyExists = "Overwrite"
sv:Sc/MyScorer/Report = 3 "Sum" "Mean" "Count_In_Bin"# One or more of Sum, Mean, Histories, Count_In_Bin, Second_Moment, Variance, Standard_Deviation, Min, Max
s:Sc/Scorer3/Quantity = "DoseToMedium"
s:Sc/Scorer3/Component = "Volume"
b:Sc/Scorer3/OutputtoConsole = "True"
s:Sc/Scorer3/IfOutputFileAlreadyExists = "Overwrite"
s:Sc/Scorer4/Quantity = "PhaseSpace"
s:Sc/Scorer4/Surface = "Volume/AnySurface"
s:Sc/Scorer4/OutputType = "ASCII"
s:Sc/Scorer4/OutputFile = "DoseToVolume"
s:Sc/Scorer4/IfOutputFileAlreadyExists = "Overwrite"
b:Sc/Scorer4/OutputtoConsole = "True"
b:Sc/Scorer4/IncludeTimeOfFlight = "True"
b:Sc/Scorer4/UsePDGEncoding = "True"
`
