Energy deposited by individual particles in the sensitive volume

[gokulgklh]

Hi all,
I was trying to simulate the energy deposited by the radiations from a radioactive source in a water column. My radioactive source is 252Cf, which gives out neutrons as result of spontaneous fission. The source(yellow colour cube) is enclosed in a stainless steel casing and is placed inside a water body ( blue colour cube).

There are different particles coming out of the source including gamma and neutron. I am using GPS to create the radioactive source 252Cf and I have successfully calculated the total energy deposited in the water column by all the radiations combined. But I want to calculate the energy deposited by each of the primary particles separately(neutrons and gammas separately). Also I want identify the different particles in the water column. I am using a modified version of example B4a and I have made the following changes in my stepping action.

void SteppingAction::UserSteppingAction(const G4Step* step)
{

  auto volume = step->GetPreStepPoint()->GetTouchableHandle()->GetVolume();


  G4double neutronEdep=0;
  G4double gammaEdep=0;
  auto edep = step->GetTotalEnergyDeposit();
  if(step->GetTrack()->GetDefinition()==G4Neutron::Neutron()){
      neutronEdep=step->GetTotalEnergyDeposit();
    }
    else if (step->GetTrack()->GetDefinition()==G4Gamma::Gamma()){
      gammaEdep=step->GetTotalEnergyDeposit();
    }
  G4Track *track = (G4Track*)(step->GetTrack());
  // step length
  G4double stepLength = 0.;

  if ( step->GetTrack()->GetDefinition()->GetPDGCharge() != 0. ) {
    stepLength = step->GetStepLength();
  } 

  if ( volume == fDetConstruction->GetAbsorberPV() ) {
    fEventAction->AddAbs(edep,stepLength);
  }
  if ( volume == fDetConstruction->GetGapPV() ) {
    G4int pdgcharge=step->GetTrack()->GetDefinition()->GetPDGCharge();
    fEventAction->AddGap(edep,neutronEdep,gammaEdep,stepLength,pdgcharge);
  }
}

//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......

}

Here gap is my sensitive volume. The function AddGap is defined in the EventAction.hh as follows.

class EventAction : public G4UserEventAction
{
  public:
    EventAction() = default;
    ~EventAction() override = default;

    void  BeginOfEventAction(const G4Event* event) override;
    void    EndOfEventAction(const G4Event* event) override;

    void AddAbs(G4double de, G4double dl);
    void AddGap(G4double de, G4double den, G4double deg, G4double dl, G4int i);


  private:
    G4double  fEnergyAbs = 0.;
    G4double  fEnergyGap = 0.;
    G4double  fEnergyGapGamma=0.;
    G4double  fEnergyGapNeutron=0.;
    G4double  fTrackLAbs = 0.;
    G4double  fTrackLGap = 0.;
    G4int fpdg=0;
};

// inline functions

inline void EventAction::AddAbs(G4double de, G4double dl) {
  fEnergyAbs += de;
  fTrackLAbs += dl;
}



inline void EventAction::AddGap(G4double de, G4double den, G4double deg, G4double dl, G4int i) {
  fEnergyGap += de;
  fEnergyGapNeutron+=den;
  fEnergyGapGamma+=deg;
  fTrackLGap += dl;
  fpdg=i;
}

}

//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......

#endif

Also my EventAction.cc looks like this:

void EventAction::BeginOfEventAction(const G4Event* /*event*/)
{
  // initialisation per event
  fEnergyAbs = 0.;
  fEnergyGap = 0.;
  fTrackLAbs = 0.;
  fTrackLGap = 0.;
  fEnergyGapGamma=0;
  fEnergyGapNeutron=0;
  fpdg=0;
}

//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......

void EventAction::EndOfEventAction(const G4Event* event)
{
  // Accumulate statistics
  //

  // get analysis manager
  auto analysisManager = G4AnalysisManager::Instance();

  // fill histograms
  analysisManager->FillH1(0, fEnergyAbs);
  analysisManager->FillH1(1, fEnergyGap);
  analysisManager->FillH1(2, fTrackLAbs);
  analysisManager->FillH1(3, fTrackLGap);
  analysisManager->FillH1(4, fEnergyGapNeutron);
  analysisManager->FillH1(5, fEnergyGapGamma); 

  // fill ntuple
  analysisManager->FillNtupleDColumn(0, fEnergyAbs);
  analysisManager->FillNtupleDColumn(1, fEnergyGap);
  analysisManager->FillNtupleDColumn(2, fTrackLAbs);
  analysisManager->FillNtupleDColumn(3, fTrackLGap);
  analysisManager->FillNtupleDColumn(4, fEnergyGapNeutron);
  analysisManager->FillNtupleDColumn(5, fEnergyGapGamma);
  analysisManager->FillNtupleDColumn(6, fpdg);


  analysisManager->AddNtupleRow();

I have saved these ntuples into a root file named B4.root in my runaction.cc, but I think I am not getting the correct values of energy deposited by the neutrons and gammas. The results I got for the total energy deposited, energy deposited by neutrons, energy deposited by gammas and pdg code are as follows

Total

Neutrons

Gammas

PDG code

I am attaching my relevant codes with this post.
RunAction.cc (5.8 KB)
SteppingAction.cc (3.6 KB)
EventAction.cc (4.1 KB)
EventAction.hh (3.4 KB)

Any help with this issue will be greatly appreciated.
Thanks a lot in advance

_Geant4 Version:_11.1.2
_Operating System:_Ubuntu 22.04

---