How to generate particles in random positions - Optical photons physics list

Dear Geant4 users,

I am trying to simulate optical photons and I would like that each particle generated was in a random position inside my box volume. I am currently using G4GeneralParticleSource instead of G4ParticleGun and I have the following configuration in my macro:

/control/verbose 0
/tracking/verbose 0
/event/verbose 0
/gps/verbose 2
/gps/particle opticalphoton
/gps/ene/type Mono
/gps/ene/mono 2.883 eV
/gps/polarization 1 1 1
/gps/pos/type Plane
/gps/pos/shape Square
/gps/pos/centre 0 0 0 cm
/gps/pos/confine DetectorVolume
/gps/pos/halfx 5 cm
/gps/pos/halfy 5 cm
/gps/ang/type iso
/gps/number 99

So what I get is that in every single run all of my photons are generated in a random position together, but I need that all of them were generated in separately positions. Does anyone know how to do that?

Also, I was wondering about the physics list for this case, I am only interested in physical process for optical photons, so using G4OpticalPhysics is enough, or is it necessary to activate something else? Because I was currently using like this:

auto *factory = new G4PhysListFactory();
auto *physicsList = factory->GetReferencePhysList("FTFP_BERT");
physicsList->RegisterPhysics( new G4OpticalPhysics() );
runManager->SetUserInitialization( physicsList );

Any help would be very appreciated! Thanks in advance!
Best regards.

GPS works by generating a ranom vertex, then it puts /gps/number particles at that vertex in a single event. Modify your macro to have

/gps/number 1
/run/beamOn 99

or however many different photons you want to simulate. You’ll get one photon per event (which makes understanding the connection between the thrown particle and the hits much simpler), and each event will have a different, random position.

1 Like

Thanks for your answer!
Although when I do as you suggested I only get a single photon per time. I guess for some reason it is not working my /run/beamOn 99. I have checked and I have both
#include "G4ParticleGun.hh" #include "G4GeneralParticleSource.hh"
in my primaryGeneratorAction.hh. Any idea what it could be?

Yes, you should get one primary photon per event, but if you’re collecting statistics, you should get one output file per run, not per event. Check how you’ve set up your RunAction and EventAction. Are you accumulating data per run, or just for the first event?

I am just a beginner in Geant4, so I am trying to follow the basic example B4a for RunAction and EventAction, as follows:

MyRunAction::MyRunAction()
: G4UserRunAction()
{
// set printing event number per each event
G4RunManager::GetRunManager()->SetPrintProgress(1);
// Create analysis manager
auto analysisManager = G4AnalysisManager::Instance();
G4cout << "Using " << analysisManager->GetType() << G4endl;
// Create directories
analysisManager->SetHistoDirectoryName(“histograms”);
analysisManager->SetNtupleDirectoryName(“ntuple”);
analysisManager->SetVerboseLevel(1);
analysisManager->SetFirstHistoId(1);
analysisManager->SetNtupleMerging(true);
// Creating histograms
analysisManager->CreateH1(“Eabs”,“Edep in absorber”, 100, 0., 800MeV);
// Creating ntuple
analysisManager->CreateNtuple(“B4”, “Edep and TrackL”);
analysisManager->CreateNtupleDColumn(“Eabs”);
analysisManager->FinishNtuple();
}
MyRunAction::~MyRunAction()
{
delete G4AnalysisManager::Instance();
}
void MyRunAction::BeginOfRunAction(const G4Run
aRun)
{
//inform the runManager to save random number seed
G4RunManager::GetRunManager()->SetRandomNumberStore(true);
// Get analysis manager
auto analysisManager = G4AnalysisManager::Instance();
// Open an output file
G4String fileName = “B4”;
analysisManager->OpenFile(fileName);
analysisManager->SetVerboseLevel(1);
// Indicates ID of first histogram created – default = 0
analysisManager->SetFirstHistoId(1); // starts couting on 1
// Creates the histogram. H1: 1 dimension
analysisManager->CreateH1(“Eabs”, “Edep in absorber”, 10, 0., 800MeV);
}
void MyRunAction::EndOfRunAction(const G4Run
aRun)
{
// save histograms & ntuple
auto analysisManager = G4AnalysisManager::Instance();
analysisManager->Write();
analysisManager->CloseFile();
delete G4AnalysisManager::Instance();
}

and

MyEventAction::MyEventAction()
: G4UserEventAction(),
fEnergyAbs(0.),
{}
MyEventAction::~MyEventAction()
{}
void MyEventAction::BeginOfEventAction(const G4Event* aEvent)
{
// initialisation per event
fEnergyAbs = 0.;
}
void MyEventAction::EndOfEventAction(const G4Event* event)
{
// Accumulate statistics
// get analysis manager
auto analysisManager = G4AnalysisManager::Instance();
// fill histograms
analysisManager->FillH1(0, fEnergyAbs);
// fill ntuple
//analysisManager->FillNtupleDColumn(0, fEnergyAbs);
//analysisManager->AddNtupleRow();
// Print per event (modulo n)
auto eventID = event->GetEventID();
auto printModulo = G4RunManager::GetRunManager()->GetPrintProgress();
if ( ( printModulo > 0 ) && ( eventID % printModulo == 0 ) ) {
G4cout << "—> End of event: " << eventID << G4endl;
G4cout
<< " Absorber: total energy: " << std::setw(7)
<< G4BestUnit(fEnergyAbs,“Energy”)
<< " total track length: " << std::setw(7)
<< G4BestUnit(fTrackLAbs,“Length”)
<< G4endl
<< " Gap: total energy: " << std::setw(7)
<< G4BestUnit(fEnergyGap,“Energy”)
<< " total track length: " << std::setw(7)
<< G4BestUnit(fTrackLGap,“Length”)
<< G4endl;
}
}

It seems like it’s generating data for each run, but I’m not sure if it’s accumulating data in the output .root file, see:

Your code is doing the right thing with the ROOT file, but I don’t see where you actually collect any useful information to put into the histogram. Look at your MyEventAction class: in the BeginOfEventAction() function, you initialize the data member fEnergyAbs = 0.; in the EndOfEventAction() function, you use the data member fEnergyAbs to add to (Fill) a histogram. But where do you do the work to put values into fEnergyAbs?

If you open the ROOT file, I predict that the histogram in that file will have 100 entries, all in the zero bin.

Indeed, that’s exactly what is happening, I’m not putting values in fEnergyAbs so i’m actually getting an empty file with 0 entries in my ROOT file. Thank you!