Alpha particles producing unexpected isotopes in 13C

I am simulating Am-241 alpha particles incident on C-13. I expect to see O-16 produced via the (alpha,n) reaction, but I am also seeing unexpected production of O-17 at the same order of magnitude as O-16. I get the same result using QGSP_BIC_HP and QGSP_BERT_HP physics lists.

Can anyone provide some insight into how Geant4 models alpha particle interactions? If this is purely a cross section issue, is it possible to view the tables that are being used? I wonder if it could be related to how the (alpha,n) reaction is modeled, and somehow an inelastic reaction occurs without the neutron being emitted.

I do not know exactly how nuclear reactions are modeled in Geant4. I wish only illustrate what can be learned from G4 public examples.

1- Energy spectrum of alphas from Am-241 decay chain is ~4 MeV → 9 MeV
See attached macro for example rdecay01.
rdecay01.mac.txt (213 Bytes)

2- (alpha, C-13) reaction. For alpha below or equal to 4 MeV, the only opened channel is O-17 production. For alpha = 9 MeV the main channel is O-16 production. The transition is in 5 MeV → 5.5 MeV. See macro for example Hadr03.
hadr03.mac.txt (305 Bytes)

3- in a realistic simulation, eg. with electromagnetic processes registered, alphas loose their energy by ionisation before to, eventually, interact.
Here a macro for example Hadr06, and the printout which confirms your observation.
hadr06.mac.txt (388 Bytes)
hadr06.out.txt (841 Bytes)

It would be good to print or plot the kinetic energy of alpha just before they interact. In SteppingAction, do something like :

if (process->GetProcessName() == “alphaInelastic”) {
G4double Ekin1 = step->GetPreStepPoint()->GetKineticEnergy();
G4double dEstep = step->GetTotalEnergyDeposit();
G4double Ekin2 = Ekin1 - dEstep;
print or plot Ekin2

Thanks for the reply, @maire. I have done a little more digging into this, as I recently ran into another similar issue.

When I run Am-241 alphas incident on C-13 in the Hadr03 example, the reaction numbers match up with what I expect, with the O-17 production being greatly reduced relative to O-16 production.

Also, in another simulation I am working on, I have thermal neutrons incident on BF3 gas. In this simulation, the number of neutron captures in the gas makes sense, but the branching ratio for secondary products and their energies (i.e. the gamma-producing state vs ground state) is wrong. Again, when I test neutrons on B-10 in Hadr03, the branching ratios are correct.

So my question has evolved a bit. Most of my simulations are built from the exampleB1 as a starting framework. It seems there is some physics process that is invoked for Hadr03 which leads to the correct secondary products in each case, which is absent from B1. I have tried using the same PhysicsList setup from Hadr03 in my simulation, but that removes electromagnetic processes (which I want to keep to get energy deposition by the secondaries), and seems to add a lot of extraneous detail to my output files, making them quite a bit larger.

Does anyone know which physics list is the key element that lets Hadr03 get secondary production right where B1 fails?

the answer is the third part of my previous mail.
Hadr03 computes hadronic cross sections and final states at nominal incident energy.
In realistic simulation, hadronic interactions do not occur at that energy.