Skip to content

Bugfix for HA2018 and HA2025 Pion Absorption #503

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Open
mking99456 wants to merge 2 commits into
base: master
Choose a base branch
from
Open

Bugfix for HA2018 and HA2025 Pion Absorption #503

Show file tree
Hide file tree
Changes from 1 commit
Commits
Show all changes
2 commits
Select commit Hold shift + click to select a range
2982370
Bugfix: For HA Pion Absorption, adding in the pion mass to the final …
mking99456 Jun 17, 2026
7816948
Updating bugfix code, fixing indentation and utilizing IsPion() function
mking99456 Jun 17, 2026
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
17 changes: 11 additions & 6 deletions src/Physics/HadronTransport/HAIntranuke2018.cxx
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -2,7 +2,7 @@
/*
Copyright (c) 2003-2025, The GENIE Collaboration
For the full text of the license visit http://copyright.genie-mc.org


Author: Steve Dytman <dytman+@pitt.edu>, Pittsburgh Univ.
Aaron Meyer <asm58@pitt.edu>, Pittsburgh Univ.
Expand Down Expand Up @@ -1241,10 +1241,13 @@ void HAIntranuke2018::Inelastic(
// simple for now, each (of 5) in hadron cluster has 1/5 of mom and KE

double probM = pLib->Find(pdgc) ->Mass();
probM -= .025; // BE correction
double BE_correction_per_nucleon = .005; //GeV
probM -= BE_correction_per_nucleon*5; // BE correction
TVector3 pP3 = p->P4()->Vect() * (1./5.);
double probKE = p->P4()->E() -probM;
double clusKE = probKE * (1./5.);
double probEAvail;
if ( pdgc==kPdgPiP || pdgc==kPdgPi0 || pdgc==kPdgPiM) probEAvail = p->P4()->E() + BE_correction_per_nucleon*5; //pion probe -- includes mass, BE correction for 4 nucleons

@nusense nusense Jun 17, 2026

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

You might consider the use of genie::pdg::IsPion(pgdc) from Framework/ParticleData/PDGUtils.h here and elsewhere.

@mking99456 mking99456 Jun 17, 2026

Copy link
Copy Markdown
Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Good call -- I'll update that

else probEAvail = p->P4()->E() -probM; //proton or neutron probe -- does not include mass
double clusKE = probEAvail * (1./5.);
TLorentzVector clusP4(pP3,clusKE); //no mass
LOG("HAIntranuke2018",pINFO) << "probM = " << probM << " ;clusKE= " << clusKE;
TLorentzVector X4(*p->X4());
Expand Down Expand Up @@ -1391,8 +1394,10 @@ void HAIntranuke2018::Inelastic(
double probM = pLib->Find(pdgc) ->Mass();
double probBE = (np+nn)*.005; // BE correction
TVector3 pP3 = p->P4()->Vect();
double probKE = p->P4()->E() - (probM - probBE);
double clusKE = probKE; // + np*0.9383 + nn*.9396;
double probEAvail;
if ( pdgc==kPdgPiP || pdgc==kPdgPi0 || pdgc==kPdgPiM) probEAvail = p->P4()->E() + probBE; //pion probe -- includes mass, BE correction for 4 nucleons
else probEAvail = p->P4()->E() - (probM - probBE); //proton or neutron probe -- does not include mass
double clusKE = probEAvail; // + np*0.9383 + nn*.9396;
TLorentzVector clusP4(pP3,clusKE); //no mass is correct
LOG("HAIntranuke2018",pINFO) << "probM = " << probM << " ;clusKE= " << clusKE;
TLorentzVector X4(*p->X4());
Expand Down
21 changes: 13 additions & 8 deletions src/Physics/HadronTransport/HAIntranuke2025.cxx
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -1258,11 +1258,14 @@ void HAIntranuke2025::Inelastic(
// simple for now, each (of 5) in hadron cluster has 1/5 of mom and KE

double probM = pLib->Find(pdgc) ->Mass();
probM -= .025; // BE correction
TVector3 pP3 = p->P4()->Vect() * (1./5.);
double probKE = p->P4()->E() -probM;
double clusKE = probKE * (1./5.);
TLorentzVector clusP4(pP3,clusKE); //no mass
double BE_correction_per_nucleon = .005; //GeV
probM -= BE_correction_per_nucleon*5; // BE correction
TVector3 pP3 = p->P4()->Vect() * (1./5.);

@nusense nusense Jun 17, 2026

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

The indent for this line seems borked. Potentially due to use of tabs perhaps.

double probEAvail;
if ( pdgc==kPdgPiP || pdgc==kPdgPi0 || pdgc==kPdgPiM) probEAvail = p->P4()->E() + BE_correction_per_nucleon*5; //pion probe -- includes mass, BE correction for 4 nucleons
else probEAvail = p->P4()->E() -probM; //proton or neutron probe -- does not include mass
double clusKE = probEAvail * (1./5.);
TLorentzVector clusP4(pP3,clusKE); //no mass
LOG("HAIntranuke2025",pINFO) << "probM = " << probM << " ;clusKE= " << clusKE;
TLorentzVector X4(*p->X4());
GHepStatus_t ist = kIStNucleonClusterTarget;
Expand Down Expand Up @@ -1408,9 +1411,11 @@ void HAIntranuke2025::Inelastic(
double probM = pLib->Find(pdgc) ->Mass();
double probBE = (np+nn)*.005; // BE correction
TVector3 pP3 = p->P4()->Vect();
double probKE = p->P4()->E() - (probM - probBE);
double clusKE = probKE; // + np*0.9383 + nn*.9396;
TLorentzVector clusP4(pP3,clusKE); //no mass is correct
double probEAvail;
if ( pdgc==kPdgPiP || pdgc==kPdgPi0 || pdgc==kPdgPiM) probEAvail = p->P4()->E() + probBE; //pion probe -- includes mass, BE correction for 4 nucleons
else probEAvail = p->P4()->E() - (probM - probBE); //proton or neutron probe -- does not include mass
double clusKE = probEAvail; // + np*0.9383 + nn*.9396;
TLorentzVector clusP4(pP3,clusKE); //no mass is correct

@nusense nusense Jun 17, 2026

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Another weird indent. Also you might want to re-word "no mass is correct" -- I assume it means "using no mass here is correct"

@mking99456 mking99456 Jun 17, 2026

Copy link
Copy Markdown
Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I'll go in and fix the indents. I can update the comment, but I will note that for my bugfix, I did not change the line with the comment as written (I think I only changed the indentation on that line by accident).

LOG("HAIntranuke2025",pINFO) << "probM = " << probM << " ;clusKE= " << clusKE;
TLorentzVector X4(*p->X4());
GHepStatus_t ist = kIStNucleonClusterTarget;
Expand Down