Gamma and neutron simulation

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.01
_Operating System:_ubuntu 20.04.6
_CMake Version:_3.16.3
iam a new learner of geant4, i want to simulate a gamma and neutron through a shielding material from a cesium-137 source which can be detected by a NaI detector and later analyse the result using root. how can i record the deposited energy on the detector? what is the best way to analyse the efficiency of my shielding? please i need some advice on how to do this. below are some of the codes i developed in trying to develop the Geometries, including, detectorconstruction, main/sim

#include "MyDetectorconstruction.hh"



G4VPhysicalVolume *MyDetectorConstruction::Construct()
	G4NistManager *nist = G4NistManager::Instance();
	G4Material *worldMat = nist->FindOrBuildMaterial("G4_AIR");
	G4Box *solidWorld = new G4Box("solidWorld", 1.0*m, 1.0*m, 1.5*m);
	G4LogicalVolume *logicWorld = new G4LogicalVolume(solidWorld, worldMat, "logicWorld");
	G4VPhysicalVolume *physWorld = new G4PVPlacement(0, G4ThreeVector(0., 0., 0.), logicWorld, "physWorld", 0, false, 0, true);
	std::vector<G4double> energy = {1.239841939*eV/0.9, 1.239841939*eV/0.2};
	//G4double startAngle = 0.0 * degree;
	//G4double spanningAngle = 360.0 * degree;
	// Define the shielding material

	G4Material *Bi2O3 = new G4Material("Bi2O3", 8.9*g/cm3, 2);
	Bi2O3->AddElement(nist->FindOrBuildElement("Bi"), 2);
	Bi2O3->AddElement(nist->FindOrBuildElement("O"), 3);
	G4Material *Polystyrene = new G4Material("Polystyrene", 1.05*g/cm3, 2);
	Polystyrene->AddElement(nist->FindOrBuildElement("C"), 8);
	Polystyrene->AddElement(nist->FindOrBuildElement("H"), 8);
	// Define the polystyreneBi2O3 material
	G4Material* PolystyreneBi2O3 = new G4Material("PolystyreneBi2O3", 1.4425 * g/cm3, 2);
	PolystyreneBi2O3->AddMaterial(Polystyrene, 95*perCent); // 95% polystyrene
	PolystyreneBi2O3->AddMaterial(Bi2O3, 5*perCent); // 5% bismuth oxide

		// Defining the shielding geometry
	G4double shielding_innerRadius = 0.0 * cm;
	G4double shielding_outerRadius = 40.0 * cm;
	G4double shielding_halfLength = 5.0 * cm;
	G4double startAngle = 0.0 * degree;
	G4double spanningAngle = 360.0 * degree;
	G4VSolid* solidshielding = new G4Tubs("shielding", shielding_innerRadius, shielding_outerRadius, shielding_halfLength, startAngle, spanningAngle);
	G4LogicalVolume *logicshielding = new G4LogicalVolume(solidshielding, PolystyreneBi2O3, "logicshielding");
	// Set the visualization attributes for the shielding
	G4VisAttributes* shieldingVisAtt = new G4VisAttributes(G4Colour(1.0, 1.0, 1.0)); // Yellow colour
	G4VPhysicalVolume *physShielding = new G4PVPlacement(0, G4ThreeVector(0., 0., 0.5*m), logicshielding, "physShielding", logicWorld, false, 0, true);
	//for detector
	G4Material* NaI = new G4Material("NaI", 3.67 * g/cm3, 2);
	NaI->AddElement(nist->FindOrBuildElement("Na"), 1);
	NaI->AddElement(nist->FindOrBuildElement("I"), 1);
	// Define the NaI detector
	G4double naI_innerRadius = 0.0 * cm;
	G4double naI_outerRadius = 40.0 * cm;
	G4double naI_halfLength = 10.0 * cm;
	G4VSolid* solidNaI = new G4Tubs("NaIDetector", naI_innerRadius, naI_outerRadius, naI_halfLength, startAngle, spanningAngle);
	logicNaI = new G4LogicalVolume(solidNaI, NaI, "logicNaI");
	// Set the visualization attributes for the NaI Detector
	G4VisAttributes* NaIVisAtt = new G4VisAttributes(G4Colour(1.0, 1.0, 0.0)); // Yellow colour
	G4VPhysicalVolume *physNaI = new G4PVPlacement(0, G4ThreeVector(0., 0., 1.2*m), logicNaI, "physNaI", logicWorld, false, 0, true);
	// Define the lead material
	G4Material* leadMaterial = new G4Material("Lead", 82, 207.2 * g/mole, 11.35 * g/cm3);
	// Define the lead collimator
	G4double collimator_innerRadius = 0.0 * cm;
	G4double collimator_outerRadius = 50.0 * cm;
	G4double collimator_halfLength = 15.0 * cm;
	//G4double startAngle = 0.0 * //degree; since it is define above no need to define again
	//G4double spanningAngle = 360.0 * degree; //degree; since it is define above no need to define again
	G4VSolid* leadCollimatorSolid = new G4Tubs("LeadCollimator", collimator_innerRadius, collimator_outerRadius, collimator_halfLength, startAngle, spanningAngle);
	G4LogicalVolume* leadCollimatorLogical = new G4LogicalVolume(leadCollimatorSolid, leadMaterial, "leadCollimatorLogical");
	// Set the visualization attributes for the lead collimator
	G4VisAttributes* leadCollimatorVisAtt = new G4VisAttributes(G4Colour(0.0, 1.0, 0.0)); // Green colour
	G4VPhysicalVolume* leadCollimatorPhysical = new G4PVPlacement(0, G4ThreeVector(0., 0., 1.2*m), leadCollimatorLogical, "leadCollimatorPhysical", logicWorld, false, 0, true);
	// Define the Cesium-137 material
	G4Material* CesiumMaterial = new G4Material("Cesium", 55, 132.9054519 * g/mole, 1.85 * g/cm3);
	// Define the particle source geometry
	G4double source_innerRadius = 0.0 * cm;
	G4double source_outerRadius = 5.0 * cm;
	G4double source_halfLength = 1.0 * cm;
	//G4double startAngle = 0.0 * degree;
	//G4double spanningAngle = 360.0 * degree;
	G4VSolid* solidCesiumsource = new G4Tubs("CesiumSource", source_innerRadius, source_outerRadius, source_halfLength, startAngle, spanningAngle);
	G4LogicalVolume *logicCesiumsource = new G4LogicalVolume(solidCesiumsource, CesiumMaterial, "logicCesiumsource");
	// Set the visualization attributes for the NaI Detector
	G4VisAttributes* CesiumVisAtt = new G4VisAttributes(G4Colour(1.0, 0.0, 0.0)); // Red colour
	G4VPhysicalVolume *physCesiumSource = new G4PVPlacement(0, G4ThreeVector(0., 0., 0), logicCesiumsource, "physCesiumSource", logicWorld, false, 0, true);
	return physWorld;

void MyDetectorConstruction::ConstructSDandField()
	MySensitiveDetector *sensDet = new MySensitiveDetector("SensitiveDetector");

#include "G4VUserDetectorConstruction.hh"
#include "G4VPhysicalVolume.hh"
#include "G4LogicalVolume.hh"
#include "G4Box.hh"
#include "G4PVPlacement.hh"
#include "G4NistManager.hh"
#include "G4SystemOfUnits.hh"
#include "G4Tubs.hh"
#include "G4VSolid.hh"
#include "G4Colour.hh"
#include "G4RotationMatrix.hh"
#include "G4ThreeVector.hh"
#include "G4VisAttributes.hh"

#include "MySensitiveDetector.hh"

class MyDetectorConstruction : public G4VUserDetectorConstruction
	virtual G4VPhysicalVolume *Construct();
	G4LogicalVolume *logicNaI;
	virtual void ConstructSDandField();




#include <iostream>

#include "G4RunManager.hh"
#include "G4UImanager.hh"
#include "G4VisManager.hh"
#include "G4VisExecutive.hh"
#include "G4UIExecutive.hh"

#include "MyDetectorconstruction.hh"
#include "physics.hh"
#include "ActionInitialization.hh"

int main(int argc, char** argv)
	G4RunManager *runManager = new G4RunManager();
	runManager->SetUserInitialization(new MyDetectorConstruction());
	runManager->SetUserInitialization(new MyPhysicsList());
	runManager->SetUserInitialization(new MyActionInitialization());

	G4UIExecutive *ui = new G4UIExecutive(argc, argv);
	G4VisManager *visManager = new G4VisExecutive();
	G4UImanager *UImanager = G4UImanager::GetUIpointer();
	UImanager->ApplyCommand("/vis/open OGL ");
	UImanager->ApplyCommand("/vis/viewer/set/viewpointVector 1 1 1");
	UImanager->ApplyCommand("/vis/viewer/set/autoRefresh true");
	UImanager->ApplyCommand("/vis/scene/add/trajectories smooth");
	UImanager->ApplyCommand("/vis/scene/endOfEventAction accumulate");


	return 0;

this is what i developed, but iam not sure if the detector is correctly defined to detect the deposited energy so also i dont know if the source is correctly initiated, even the shielding materials too.

You can use your sensitive detector function and an analysis manager histogram to track energy deposited in each material you are interested in. You need to generate your particles (gammas, neutrons) using a particle gun to define the location, momenum, and energy of your source (you can randomize the energy over a range as well in your primary generator action class).

For efficiency, you just need to look at number of absorptions/number produced particles. You will be analyzing your energy deposited spectrum and number of interactions as a function of particle and shielding geometer, I assume

1 Like

Today I will try to give an attempt with files you shared.


To start with,
Some of my suggestions are as follows;

In DetectorConstruction, replace source material with G4_Galactic,

In PrimaryGeneratorAction replace particleGun with gps (as it is much versatile)
e.g. as

#### This will be Cs137 source.
/gps/particle ion
/gps/ion 55 137

Use Shielding Physics [#include “Shielding.hh”] as it is based on best-guess selection of electromagnetic and hadronic physics processes.

Record events/hits in the end of RunAction. Therefore, create a function that use concrete class G4Run.
Look e.g. rdecay01/


Thank you Sir, i will surely try it this way.

Also @BashirAminu ,

when you use gps
from the macro you confine it with the source.
/gps/pos/type Volume
/gps/pos/confine CesiumSource #This is logical name.
/gps/pos/shape Sphere
/gps/pos/centre 0 0 0 mm
/gps/pos/radius 5 cm

Thank you Sir, besides can you share a link to access this rdecay01 example? Also the main target of my program is to determine a new shielding material(using my target material) that is more reliable in shielding a gamma and a Neutron. see the the visualized picture attached. thanks again for your prompt reply.

iam also very new to Geant4. if there is any simple method to simulate what i explain above please guide me thanks. from the above setup, i shoot a gamma from a point source passing through my shielding including the material and i want to detect it on the NaI detector which is enclosed/collimated by lead. i got all these setup from the above codes i shared.

Hi @BashirAminu

The ground principle is attenuation, DE, CounterParticles.

The concept to keep source materiel “G4_Galactic” is to fill it with Cs137 source from the GPS.
For that reason I confined with source logical.

Geometry has no problem. It is what you want.

Are you able to confine gun or not ? If not, this is how you define gps in PGA. (3.5 KB)
PrimaryGeneratorAction.hh (500 Bytes)

For the rdecay01,
check examples/extended/radioactivedecay

Please go through some basic example too, examples/basic/B4a. It will also help you to design your model.


thank you so much for your kind guidance, it is not yet done i will go through it, and i will give you the feedback.

Hi @BashirAminu

For sure.
Make efforts and ask key point questions. You will develop the model soon.
Practice it.

This topic was automatically closed 30 days after the last reply. New replies are no longer allowed.