72 std::ostringstream ost;
89 if(!std::getenv(
"G4NEUTRONHPDATA"))
90 throw G4HadronicException(__FILE__, __LINE__,
"Please setenv G4NEUTRONHPDATA to point to the neutron cross-section files where Inelastic/Gammas data is found.");
91 G4String tBase = std::getenv(
"G4NEUTRONHPDATA");
98 if( std::getenv(
"G4ParticleHPDebug") )
G4cout <<
" G4ParticleHPInelasticCompFS::Init FILE " << filename <<
G4endl;
110 if(std::getenv(
"G4ParticleHPDebug_NamesLogging"))
G4cout <<
"Skipped = "<< filename <<
" "<<A<<
" "<<Z<<
G4endl;
131 G4int infoType, dataType, dummy;
134 while (theData >> infoType)
138 theData >> sfType >> dummy;
140 if(sfType>=600||(sfType<100&&sfType>=50)) it = sfType%50;
148 theData >> dqi >> ilr;
174 else if(dataType==12)
179 else if(dataType==13)
184 else if(dataType==14)
188 else if(dataType==15)
194 throw G4HadronicException(__FILE__, __LINE__,
"Data-type unknown to G4ParticleHPInelasticCompFS");
209 if(i!=0) running[i]=running[i-1];
221 for(i0=0; i0<50; i0++)
225 if(random < running[i0]/sum)
break;
284 boosted.
Lorentz(incidReactionProduct, theTarget);
318 (targetMass - residualMass);
320 if ( availableEnergy < 0 )
325 G4int nothingWasKnownOnHadron = 0;
337 (aHadron.
GetMass()+residualMass));
346 if ( p2 > 0.0 ) p = std::sqrt( p );
369 if (
QI[it] < 0 || 849 <
QI[it] ) dqi =
QI[
it];
374 G4int icounter_max=1024;
377 if ( icounter > icounter_max ) {
378 G4cout <<
"Loop-counter exceeded the threshold value at " << __LINE__ <<
"th line of " << __FILE__ <<
"." <<
G4endl;
382 }
while(eSecN>MaxEne);
385 eGamm = eKinetic-eSecN;
391 iLevel+=
G4int(random);
398 while (eKinetic-eExcitation < 0 && iLevel>0)
407 if ( dqi < 0 || 849 < dqi ) useQI =
true;
413 eExcitation =
QI[0] -
QI[
it];
414 if(eExcitation < 20*
CLHEP::keV) eExcitation = 0;
441 if (iLevel == -1) iLevel = 0;
449 if (!find) iLevel = imaxEx;
452 if(std::getenv(
"G4ParticleHPDebug") && eKinetic-eExcitation < 0)
454 throw G4HadronicException(__FILE__, __LINE__,
"SEVERE: InelasticCompFS: Consistency of data not good enough, please file report");
456 if(eKinetic-eExcitation < 0) eExcitation = 0;
475 theRestEnergy->
SetMomentum(eGamm*std::sin(std::acos(costh))*std::cos(phi),
476 eGamm*std::sin(std::acos(costh))*std::sin(phi),
479 thePhotons->push_back(theRestEnergy);
490 if ( theParticles != NULL ) {
495 for (
G4int j = 0 ; j != (
G4int)theParticles->size() ; j++ ) {
496 if ( theParticles->at(j)->GetDefinition()->GetBaryonNumber() >
maxA ) {
497 maxA = theParticles->at(j)->GetDefinition()->GetBaryonNumber();
500 sumA += theParticles->at(j)->GetDefinition()->GetBaryonNumber();
501 sumZ +=
G4int( theParticles->at(j)->GetDefinition()->GetPDGCharge() +
eps );
505 if ( dA < 0 || dZ < 0 ) {
506 G4int newA = theParticles->at(jAtMaxA)->GetDefinition()->GetBaryonNumber() + dA ;
507 G4int newZ =
G4int( theParticles->at(jAtMaxA)->GetDefinition()->GetPDGCharge() +
eps ) + dZ;
509 theParticles->at( jAtMaxA )->SetDefinition( pd );
518 nothingWasKnownOnHadron = 1;
532 boosted_tmp.
Lorentz(incidReactionProduct, theTarget);
537 if(thePhotons!=0 && thePhotons->size()!=0)
538 { aBaseEnergy-=thePhotons->operator[](0)->GetTotalEnergy(); }
544 G4bool foundMatchingLevel =
false;
562 if ( thePhotons != NULL ) thePhotons->push_back(theNext->operator[](0));
563 aBaseEnergy = testEnergy-theNext->operator[](0)->GetTotalEnergy();
565 foundMatchingLevel =
true;
574 if(!foundMatchingLevel)
578 if ( thePhotons != NULL ) thePhotons->push_back(theNext->operator[](0));
579 aBaseEnergy = aBaseEnergy-theNext->operator[](0)->GetTotalEnergy();
588 for(i0=0; i0<thePhotons->size(); i0++)
591 thePhotons->operator[](i0)->Lorentz(*(thePhotons->operator[](i0)), -1.*theTarget);
595 if(nothingWasKnownOnHadron)
604 if(thePhotons==0 && eExcitation>0){
657 G4int nSecondaries = 2;
658 G4bool needsSeparateRecoil =
false;
659 G4int totalBaryonNumber = 0;
660 G4int totalCharge = 0;
662 if(theParticles != 0)
664 nSecondaries = theParticles->size();
667 for(ii0=0; ii0<theParticles->size(); ii0++)
669 aDef = theParticles->operator[](ii0)->GetDefinition();
672 totalMomentum += theParticles->operator[](ii0)->GetMomentum();
676 needsSeparateRecoil =
true;
686 if(thePhotons!=0) { nPhotons = thePhotons->size(); }
687 nSecondaries += nPhotons;
691 if( theParticles==0 )
701 aHadron.
Lorentz(aHadron, theTarget);
704 ->GetIon(static_cast<G4int>(residualZ), static_cast<G4int>(residualA), 0));
712 theResidual.
Lorentz(theResidual, -1.*theTarget);
716 for(i=0; i<nPhotons; i++)
718 totalPhotonMomentum += thePhotons->operator[](i)->GetMomentum();
731 for(i0=0; i0<theParticles->size(); i0++)
734 theSec->
SetDefinition(theParticles->operator[](i0)->GetDefinition());
735 theSec->
SetMomentum(theParticles->operator[](i0)->GetMomentum());
740 delete theParticles->operator[](i0);
743 if(needsSeparateRecoil && residualZ!=0)
747 ->GetIon(static_cast<G4int>(residualZ), static_cast<G4int>(residualA), 0));
749 resiualKineticEnergy += totalMomentum*totalMomentum;
750 resiualKineticEnergy = std::sqrt(resiualKineticEnergy) - theResidual.
GetMass();
773 for(i=0; i<nPhotons; i++)
778 theSec->
SetDefinition( thePhotons->operator[](i)->GetDefinition() );
780 theSec->
SetMomentum(thePhotons->operator[](i)->GetMomentum());
786 delete thePhotons->operator[](i);
828 theCMSproj.
Lorentz(*proj,theCMS);
829 theCMStarg.
Lorentz(*targ,theCMS);
834 G4double fmomsquared=1./4./totE/totE*(totE*totE-(prodmass-resmass)*(prodmass-resmass))*(totE*totE-(prodmass+resmass)*(prodmass+resmass));
837 fmom=std::sqrt(fmomsquared);
845 product->
SetMomentum(fmom*sinth*std::cos(phi),fmom*sinth*std::sin(phi),fmom*cosTh);
848 product->
Lorentz(*product,-1.*theCMS);
890 for (
G4int j=0; j<4; j++ )
892 theProds[j].
Lorentz(theProds[j], -1.*theTarget);
930 for (
G4int j=0; j<2; j++ )
933 theProds[j].
Lorentz(theProds[j], -1.*theTarget);
944 G4Exception(
"G4ParticleHPInelasticCompFS::CompositeApply()",
"G4ParticleInelasticCompFS.cc",
FatalException,
"Alpha production with LR!=0.");