Save data of geant in tree

Hello,
I need to save the beam information in the tree

In the analysis messenger I fill a vector with beam information

void AnalysisManager::FillhPrimaries(G4ThreeVector fPos,G4double fkinE,G4ThreeVector fmom)
{
  BeamX.push_back(fPos.x());
  BeamY.push_back(fPos.y());
  BeamZ.push_back(fPos.z());
  BeamKinE.push_back(fkinE);
  BeamTheta.push_back(fmom.theta()/deg);  
}

moreover, to clean the vector at each primary particle, I’ve

void AnalysisManager::ClearPrimaries(){
  BeamX.clear();
  BeamY.clear();
  BeamZ.clear();
  BeamKinE.clear();
  BeamTheta.clear();
}

the ClearPrimaries function is called in the event action, because I want that for each event, it saves the information in the tree

void EventAction::BeginOfEventAction(const G4Event* evt)
{
 ....
  analysismanager->ClearPrimaries();

but when running the simulation, the number of entries isn’t correct. For example, I simulated 1000 primary particles, but I don’t get entries in the TTree

if I call the Clear Primaries in

void EventAction::EndOfEventAction(const G4Event* evt)
{

I get 57 entries (instead of 1000)
Then I don’t know how to save these information in the TTree

AnalysisManager.cc (10.5 KB)
EventAction.cc (3.7 KB)

PS. FillPrimaryParticles function is called in the primarygeneratoraction

PrimaryGeneratorAction.cc (1.9 KB)

I also tried to call the clearprimaries in the primarygeneratoraction, but I get 57 entries instead of 1000. Moreover, it would be better to not clean the variables in primarygeneratoraction (if possible) because it is a shared files with other experiments.

That sounds to me like you’re running your job with many worker threads. Each thread does it’s own collection. Do you have G4AnalysisManager set up to merge the tuples across threads at the end of the run? Otherwise, you’ll only get one thread’s worth recorded.

Hi @mkelsey thank you.
If I correctly understood what you mean, I added in AnalysisManager.cc

	AnalysisManager::AnalysisManager():
  fileName("Default.root")
{
	 auto analysismanager = G4AnalysisManager::Instance();
	// Default settings
  #ifdef G4MULTITHREADED
  analysismanager->SetNtupleMerging(true);
  #endif

to merge the Ntuple and I call the cleanPrimaries in

void EventAction::EndOfEventAction(const G4Event* evt)
{

...

  analysismanager->ClearPrimaries();

but I still get only 57 entries

AnalysisManager.cc (10.8 KB)
EventAction.cc (3.7 KB)

That is what I meant. It looks like AnalysisManager is your own class wrapped around G4AnalysisManager? This code is very confusing. You set up the G4AnalysisManager, but then you don’t ever use it for anything. You open your own ROOT file, you create your own TTrees and ROOT histograms. If you don’t make use of G4AnalysisManager to manage your data, then it can’t do anything about managing it for you.

None of this seems to be properly thread safe or protected by mutexes.

Do you have a separate AnalysisManager instance on every worker thread? If so, then each of those instances is opening the same ROOT file, and each one is therefore stomping on all of the other worker threads.

If you are comfortable and familiar with writing thread-safe code, global singletons and mutexes, then you should bring your code up to spec. This forum is not the place for you to learn how to do that.

If you are not comfortable with that (and it’s not trivial), then you should consider using AnalysisManager and its functionality directly, instead of using raw ROOT yourself.

Hi @mkelsey

How can I know it?
I’m not the author of the simulation. I’m using it and I made some changes to adapt it to my need. I.e. I added the storage of secondary particles in the Tree (because the author just provided the storage in histograms.
I created the tree for secondary particles and it works…but I’m having problems with primary particles.

Do you mean saving data in histograms?
As you can see in the code, it fills the histograms

void AnalysisManager::FillhPrimaries(G4ThreeVector fPos,G4double fkinE,G4ThreeVector fmom)
{
  hPrimaryXvsY->Fill(fPos.x(),fPos.y());
  hPrimaryZ->Fill(fPos.z());
  hPrimarykinE->Fill(fkinE);
  hPrimaryDiv->Fill(fmom.theta()/deg);
}

but they are not so comfortable as the Tree for the off line analysis…

Hello, I solved myself