diff --git a/src/Physics/HadronTransport/HAIntranuke2018.cxx b/src/Physics/HadronTransport/HAIntranuke2018.cxx index 682470dde..47dcc7865 100644 --- a/src/Physics/HadronTransport/HAIntranuke2018.cxx +++ b/src/Physics/HadronTransport/HAIntranuke2018.cxx @@ -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 , Pittsburgh Univ. Aaron Meyer , Pittsburgh Univ. @@ -1241,11 +1241,14 @@ 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.); - TLorentzVector clusP4(pP3,clusKE); //no mass + double probEAvail; + if ( genie::pdg::IsPion(pgdc) ) 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); //using no mass here is correct LOG("HAIntranuke2018",pINFO) << "probM = " << probM << " ;clusKE= " << clusKE; TLorentzVector X4(*p->X4()); GHepStatus_t ist = kIStNucleonClusterTarget; @@ -1391,9 +1394,11 @@ 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; - TLorentzVector clusP4(pP3,clusKE); //no mass is correct + double probEAvail; + if ( genie::pdg::IsPion(pgdc) ) 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); //using no mass here is correct LOG("HAIntranuke2018",pINFO) << "probM = " << probM << " ;clusKE= " << clusKE; TLorentzVector X4(*p->X4()); GHepStatus_t ist = kIStNucleonClusterTarget; diff --git a/src/Physics/HadronTransport/HAIntranuke2025.cxx b/src/Physics/HadronTransport/HAIntranuke2025.cxx index cb5027f69..5885a98d4 100644 --- a/src/Physics/HadronTransport/HAIntranuke2025.cxx +++ b/src/Physics/HadronTransport/HAIntranuke2025.cxx @@ -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.); + double probEAvail; + if ( genie::pdg::IsPion(pgdc) ) 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); //using no mass here is correct LOG("HAIntranuke2025",pINFO) << "probM = " << probM << " ;clusKE= " << clusKE; TLorentzVector X4(*p->X4()); GHepStatus_t ist = kIStNucleonClusterTarget; @@ -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 ( genie::pdg::IsPion(pgdc) ) 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); //using no mass here is correct LOG("HAIntranuke2025",pINFO) << "probM = " << probM << " ;clusKE= " << clusKE; TLorentzVector X4(*p->X4()); GHepStatus_t ist = kIStNucleonClusterTarget;