Measuring mass attenuation coefficient of CsI

kevin.islami · June 3, 2025, 9:22am

Please fill out the following information to help in answering your question, and also see tips for posting code snippets. If you don’t provide this information it will take more time to help with your problem!

_Geant4 Version:_11.3.1
_Operating System:_Windows 10, WSL 2, Ubuntu 24.04
_Compiler/Version:_13.3
_CMake Version:_3.28.3

Hi,

I am currently trying to recreate the mass attenuation coefficient vs energy plot of CsI that is found in the NIST database.

For this I shoot a photon of different energies (1 to 511 keV in steps of 10) towards a CsI crystal for 1000 times per energy and compute the mass attenuation coefficient out of the energy the photon loses in the process with the Beer-Lamberdt law.

However, the values I obtain in the simulation in the range of ~ 400 - 511 keV are less apart than below 400 keV where they start to differ immensly the closer I get to 1 keV. In general all simulated values tend to be lower than the ones given in the NIST database.

Do you have an idea why this might occur? Maybe it is because of the choice of physics lists? (I am using the G4EmLivermorePhysics list).

Cheers,
Kevin

jrellin · June 3, 2025, 3:30pm

TestEM13 I think does what you want although more generally to any particle type. It also lets you change EM lists with macro commands to do quick comparisons.

atolosad · June 4, 2025, 4:17pm

Hi @kevin.islami,

Welcome to the Geant4 forum!

To get the mass attenuation coefficient (MAC) more directly and accurately, I recommend using Geant4’s EM utilities instead. Here’s an example code snippet that builds a custom CsI material and computes the MAC using G4EmCalculator:

auto* nist = G4NistManager::Instance();

// Get Cs and I elements from NIST
G4Element* elCs = nist->FindOrBuildElement("Cs");
G4Element* elI  = nist->FindOrBuildElement("I");

// Define CsI material with density (approx. 4.51 g/cm³)
G4double density = 4.51 * g/cm3;
G4Material* csI = new G4Material("Custom_CsI", density, 2);
csI->AddElement(elCs, 1); // 1 atom of Cs
csI->AddElement(elI, 1);  // 1 atom of I

G4double energy = 100 * keV;

G4EmCalculator emCalc;
G4double lambda = emCalc.ComputeGammaAttenuationLength(energy, csI);
G4double mu_over_rho = 1.0 / (lambda * csI->GetDensity());

You are correct that the MAC values can vary somewhat depending on the physics list used, because:

The Livermore physics list relies on updated evaluated data (EPICS2014),
The Penelope physics list uses model-specific parameterizations,
The Standard physics list uses analytical models rather than tabulated data.

Additionally:

Differences with the NIST XCOM database may arise due to varying data sources or smoothing algorithms
At energies below ~10 keV, atomic binding and subshell effects become important, which can cause further differences depending on model

If you can provide a comparison plot of MAC vs energy showing NIST data alongside Geant4 results from different physics lists (Standard, Livermore, Penelope), it will be very helpful for investigating the sources of discrepancy.

Thanks for your time and interest!

Best regards,
Alvaro

[1] https://github.com/Geant4/geant4/blob/20a218bbe1d8b4a78f819e396f8be615e6a38358/source/processes/electromagnetic/utils/src/G4EmCalculator.cc#L685

kevin.islami · June 5, 2025, 8:49am

Hi Alvaro,

Thank you very much for the quick and detailed response!

I tried implementing your suggestion and calculating the MAC in MyEventAction, where I fetched the initial energy of every gamma of every run, however, the MAC I get for e.g 511 keV is of the order e-21 which (I think) seems off.

Is this maybe not the right environment to compute this in?

Do you maybe also know how the peaks in the NIST data arise and whether they are simulatable in Geant4 since they seem to shift the “ideal” logarithmic curve which the Beer-Lambert law follows sightly upwards?

(Apologies if some of these questions seem trivial since it’s part of my first ever project).

Best regards,
Kevin

Plot is taken from: NIST: X-Ray Mass Attenuation Coefficients - Cesium Iodide

atolosad · June 5, 2025, 3:47pm

Dear @kevin.islami

Thank you for your quick response.

Could you please try to modify the Geant4 application distributed as an example here [1] to extract the mass attenuation coefficients?

To define a new material, you could for example add the material definition in the DetectorConstructor [2], for example after the line 165 you could add these piece of code:

    auto* nist = G4NistManager::Instance();
    G4Element* elCs = nist->FindOrBuildElement("Cs");
    G4Element* elI  = nist->FindOrBuildElement("I");

    // Define CsI material with density (approx. 4.51 g/cm³)
    density = 4.51 * g/cm3;
    G4Material* csI = new G4Material("Custom_CsI", density, 2);
    csI->AddElement(elCs, 1); // 1 atom of Cs
    csI->AddElement(elI, 1);  // 1 atom of I

After compiling, you can execute it passing as argument the following macro file:

# Macro file for "TestEm0.cc"
#
/control/verbose 0
/run/verbose 1
#
#/testem/phys/addPhysics local
/testem/phys/addPhysics emstandard_opt4
#
/process/eLoss/verbose 0
/process/em/printParameters
#
/run/initialize
#
/testem/det/setMat Custom_CsI
/gun/particle gamma
/gun/energy 100 keV
/run/beamOn
/gun/energy 200 keV
/run/beamOn
/gun/energy 300 keV
/run/beamOn
/gun/energy 500 keV
/run/beamOn
/gun/energy 1000 keV
/run/beamOn

These are the results for few energies, values from NIST and Geant4 electromagnetic physics standard 4 (EM4) agree within 0.5%

	CsI mu / (cm2/g)
Energy / MeV	NIST	EM4
0.1	2.035	2.03644
0.2	0.3805	0.382286
0.3	0.1818	0.181533
0.5	0.09809	0.0977995
1	0.05848	0.0583394

Other physics list may give slightly different numbers, but not far away from these. If you evolve the TestEM0 to print these coefficients for a range of energies and different physics list, it would be interesting to see how big the differences are

Regarding your method for extracting the attenuation coefficient by using Beer-Lamberdt law, maybe it is worth to try out with more photons?

Best,
Alvaro

[1] geant4/examples/extended/electromagnetic/TestEm0 at 20a218bbe1d8b4a78f819e396f8be615e6a38358 · Geant4/geant4 · GitHub
[2] geant4/examples/extended/electromagnetic/TestEm0/src/DetectorConstruction.cc at 20a218bbe1d8b4a78f819e396f8be615e6a38358 · Geant4/geant4 · GitHub

kevin.islami · June 11, 2025, 11:05am

Dear @atolosad,

Apologies for the late response. I modified the Geant4 application and extracted the MAC for the EM4, Livermore and Penelope lists.

In my compuation EM4 and Livermore gave the same results. Penelope had closer ones to NIST in the area of around 1 keV.

Here is a Plot I made of the three lists and the NIST data combined:

In case you want to have a look at the data yourself, here are the .txt files of the energies I used and the values I received:

AttenuationTestEM4.txt (901 Bytes)
AttenuationTestLiv.txt (838 Bytes)
AttenuationTestPen.txt (840 Bytes)
NISTDataMAC.txt (4.0 MB)

As for the Beer-Lamerdt law, I tried to use 1000 photons but somehow the result still seemed to be the same. I assume that using the built in Geant4 calculators is more accurate, however, it still gives a value that is too low. Maybe there has to be some pyhsics to be added to get the right MAC values?

Best,
Kevin