I am trying to build a model of a plastic scintillator with scintillation by particle type (after poring over the application developers book), but I’m running into an issue where I have defined the ELECTRONSCINTILLATIONYIELD in my material properties table, yet somehow Geant4 is kicking an error that this property is not defined.
This is the error:
-------- EEEE ------- G4Exception-START -------- EEEE -------
*** G4Exception : Scint01
issued by : G4Scintillation::PostStepDoIt
G4Scintillation::PostStepDoIt(): Request for scintillation yield for energy deposit and particle
type without correct entry in MaterialPropertiesTable.
ScintillationByParticleType requires at minimum that
ELECTRONSCINTILLATIONYIELD is set by the user
Missing MaterialPropertiesTable entry - No correct entry in MaterialPropertiesTable
*** Fatal Exception *** core dump ***
G4WT0 > G4Track (0x164c592d0) - track ID = 4, parent ID = 1
G4WT0 > Particle type : C12 - creator process : conv, creator model : model_AugerElectron
G4WT0 > Kinetic energy : 0 eV - Momentum direction : (-0.114823,0.779454,-0.615847)
G4WT0 > Step length : 10.1354 nm - total energy deposit : 81.697 eV
G4WT0 > Pre-step point : (9.64211,-1.06817,5.49988) - Physical volume : CherenkovDet (BC422)
G4WT0 > - defined by : not available
G4WT0 > Post-step point : (9.64211,-1.06817,5.49987) - Physical volume : CherenkovDet (BC422)
G4WT0 > - defined by : ionIoni - step status : 3
G4WT0 > *** Note: Step information might not be properly updated.
G4WT0 >
-------- EEEE -------- G4Exception-END --------- EEEE -------
Here is my material definition:
//BC-422 Plastic Scintillator---------------------------------------------------------------------
matBC422 = new G4Material("BC422", density=1.032*g/cm3, ncomponents=2);
matBC422->AddElement(C, natoms=10);
matBC422->AddElement(H, natoms=11);
G4MaterialPropertiesTable* BC422_MPT = new G4MaterialPropertiesTable();
G4double BC422_rise_time = 0.35*ns; //Source: St. Gobain data sheet
G4double BC422_decay_time = 1.6*ns;
vector<G4double> BC422_rindex = {1.573,1.575,1.58,1.595,1.605,1.63,1.644}; //Source: Nakamura 2015
vector<G4double> BC422_energy_r = {WavelengthToEnergy(700)*eV,WavelengthToEnergy(656)*eV,WavelengthToEnergy(590)*eV,WavelengthToEnergy(487)*eV,
WavelengthToEnergy(436)*eV,WavelengthToEnergy(365)*eV,WavelengthToEnergy(340)*eV};
vector<G4double> BC422_abs_length = {8.0*cm,8.0*cm,8.0*cm,8.0*cm,8.0*cm,8.0*cm,8.0*cm}; //Source: St. Gobain
vector<G4double> BC422_emission = {2.37,7.19,11.09,15.17,19.81,25.28,30.95,37.63,43.38,50.44,57.30,64.92,71.22,79.02,86.64,96.10,100.,98.89,87.0,66.59,46.73,27.61,9.42,2.92};
vector<G4double> BC422_energy_em = {WavelengthToEnergy(458)*eV,WavelengthToEnergy(450)*eV,WavelengthToEnergy(445)*eV,WavelengthToEnergy(440)*eV,
WavelengthToEnergy(435)*eV,WavelengthToEnergy(430)*eV,WavelengthToEnergy(425)*eV,WavelengthToEnergy(420)*eV,WavelengthToEnergy(415)*eV,
WavelengthToEnergy(410)*eV,WavelengthToEnergy(405)*eV,WavelengthToEnergy(400)*eV,WavelengthToEnergy(395)*eV,WavelengthToEnergy(390)*eV,
WavelengthToEnergy(385)*eV,WavelengthToEnergy(380)*eV,WavelengthToEnergy(376.5)*eV,WavelengthToEnergy(375)*eV,WavelengthToEnergy(370)*eV,
WavelengthToEnergy(365)*eV,WavelengthToEnergy(360)*eV,WavelengthToEnergy(355)*eV,WavelengthToEnergy(350)*eV,WavelengthToEnergy(347)*eV};
//energy in MeV
vector<G4double> BC422_LY_energy = {0.01*MeV,1.*MeV,5.*MeV,10.*MeV,20.*MeV,30.*MeV,40.*MeV,50.*MeV,60.*MeV,70.*MeV,80.*MeV,90.*MeV,100.*MeV,110.*MeV,120.*MeV,130.*MeV};
//Light yield in number of photons, Source: St. Gobain BC422 relative light yield plot
vector<G4double> BC422_electron_LY = {110.,11000.,55000.,110000.,220000.,330000.,440000.,550000.,660000.,770000.,880000.,990000.,1100000.,1210000.,1320000.,1430000.};
vector<G4double> BC422_proton_LY = {90.,9041.,28630.,63288.,141644.,226027.,317945.,408356.,504795.,607260.,703699.,800137.,902603.,1005068.,1109041.,1211507.};
vector<G4double> BC422_alpha_LY = {10.,1506.,10548.,24110.,51233.,81370.,120548.,159726.,204932.,250137.,301370.,355616.,414384.,474658.,537945.,602740.};
//apply properties to material table
BC422_MPT->AddProperty("RINDEX",BC422_energy_r,BC422_rindex,false,true);
BC422_MPT->AddProperty("ABSLENGTH",BC422_energy_r,BC422_abs_length,false,true);
BC422_MPT->AddProperty("SCINTILLATIONCOMPONENT1",BC422_energy_em,BC422_emission,false,true);
BC422_MPT->AddProperty("ELECTRONSCINTILLATIONYIELD",BC422_LY_energy,BC422_electron_LY,false,true);
BC422_MPT->AddProperty("PROTONSCINTILLATIONYIELD",BC422_LY_energy,BC422_proton_LY,false,true);
BC422_MPT->AddProperty("ALPHASCINTILLATIONYIELD",BC422_LY_energy,BC422_alpha_LY,false,true);
BC422_MPT->AddConstProperty("SCINTILLATIONTIMECONSTANT1",BC422_decay_time,false);
BC422_MPT->AddConstProperty("SCINTILLATIONRISETIME1",BC422_rise_time,false);
//BC422_MPT->AddConstProperty("ELECTRONSCINTILLATIONYIELD1",1.0,false);
BC422_MPT->AddConstProperty("RESOLUTIONSCALE",1.0,false);
matBC422->SetMaterialPropertiesTable(BC422_MPT);
BC422_MPT->DumpTable();
Interestingly, the DumpTable function seems to show that ELECTRONSCINTILLATIONYIELD is set, yet somehow it cannot be accessed when the particle is being tracked? DumpTable output below:
0: RINDEX
1.7712e-06 1.573
1.89e-06 1.575
2.10143e-06 1.58
2.54588e-06 1.595
2.84367e-06 1.605
3.39683e-06 1.63
3.64659e-06 1.644
9: GROUPVEL
1.7712e-06 186.925
1.8306e-06 186.809
1.99571e-06 184.527
2.32365e-06 179.982
2.69478e-06 177.35
3.12025e-06 170.516
3.64659e-06 162.814
16: ABSLENGTH
1.7712e-06 80
1.89e-06 80
2.10143e-06 80
2.54588e-06 80
2.84367e-06 80
3.39683e-06 80
3.64659e-06 80
17: PROTONSCINTILLATIONYIELD
0.01 90
1 9041
5 28630
10 63288
20 141644
30 226027
40 317945
50 408356
60 504795
70 607260
80 703699
90 800137
100 902603
110 1.00507e+06
120 1.10904e+06
130 1.21151e+06
20: ALPHASCINTILLATIONYIELD
0.01 10
1 1506
5 10548
10 24110
20 51233
30 81370
40 120548
50 159726
60 204932
70 250137
80 301370
90 355616
100 414384
110 474658
120 537945
130 602740
22: ELECTRONSCINTILLATIONYIELD
0.01 110
1 11000
5 55000
10 110000
20 220000
30 330000
40 440000
50 550000
60 660000
70 770000
80 880000
90 990000
100 1.1e+06
110 1.21e+06
120 1.32e+06
130 1.43e+06
23: SCINTILLATIONCOMPONENT1
2.70708e-06 2.37
2.7552e-06 7.19
2.78616e-06 11.09
2.81782e-06 15.17
2.85021e-06 19.81
2.88335e-06 25.28
2.91727e-06 30.95
2.952e-06 37.63
2.98757e-06 43.38
3.024e-06 50.44
3.06134e-06 57.3
3.0996e-06 64.92
3.13884e-06 71.22
3.17908e-06 79.02
3.22037e-06 86.64
3.26274e-06 96.1
3.29307e-06 100
3.30624e-06 98.89
3.35092e-06 87
3.39683e-06 66.59
3.444e-06 46.73
3.49251e-06 27.61
3.54241e-06 9.42
3.57303e-06 2.92
11: RESOLUTIONSCALE 1
30: SCINTILLATIONTIMECONSTANT1 1.6
33: SCINTILLATIONRISETIME1 0.35
I am able to run the extended optical examples LXe and OpNovice2 without problems, so I’m sure I’ve done something wrong, but I just can’t spot it. Also, if I turn off scintillation by particle and just specify a constant property for the scintillation yield, it works fine. Problem is I don’t have a Birks’ constant for the material, so scintillation by particle is the way to get the response to neutrons, etc.
As an aside, it was not clear to me from the application book that RESOLUTIONSCALE is a required parameter and does not just default to 1. I was previously getting errors of:
-------- EEEE ------- G4Exception-START -------- EEEE -------
*** G4Exception : mat202
issued by : G4MaterialPropertiesTable::GetConstProperty()
Constant Material Property Index 11 not found.
*** Fatal Exception *** core dump ***
and it turned out that Index 11 is RESOLUTIONSCALE. Just a side note, but maybe that information will help someone.
Geant4 Version: 11.0.2
Operating System: MacOS 12.6.9
Compiler/Version: Apple clang 14.0.0
CMake Version: 3.24.1