63 fCache.Get()->currentMeanEnergy = -2;
64 fCache.Get()->fresh =
true;
66 if ( getenv(
"G4PHP_DO_NOT_ADJUST_FINAL_STATE" ) )
adjustResult =
false;
76 if ( getenv(
"G4PHP_DO_NOT_ADJUST_FINAL_STATE" ) )
adjustResult =
false;
100 if( getenv(
"G4PHPTEST") )
G4cout <<
" G4ParticleHPContAngularPar::Sample " << anEnergy <<
" " << massCode <<
" " << angularRep <<
G4endl;
131 if(Z!=2)
throw G4HadronicException(__FILE__, __LINE__,
"G4ParticleHPContAngularPar: Unknown ion case 1");
142 if( angularRep == 1 )
153 if (
fCache.Get()->fresh == true )
158 fCache.Get()->fresh =
false;
175 if (
theAngular[j].GetValue(0) != 0.0 )
break;
189 running[j+1] = running[j] + delta;
207 if ( j == nDiscreteEnergies ) {
218 if ( j != nEnergies-1 ) {
223 throw G4HadronicException(__FILE__, __LINE__,
"G4ParticleHPContAngularPar: Unexpected non zero value of theAngular[nEnergies-1].GetValue(0)");
228 running[j+1] = running[j] + ( ( e_high - e_low ) * delta );
244 random *= (tot_prob_DIS + tot_prob_CON);
246 if ( random <= ( tot_prob_DIS / ( tot_prob_DIS + tot_prob_CON ) ) || nDiscreteEnergies ==
nEnergies )
252 if ( random < running[ j+1 ] )
276 if ( random < running[ j ] )
287 if ( it != nDiscreteEnergies )
301 if ( it == nDiscreteEnergies ) itt = it+1;
330 if (
fCache.Get()->fresh == true )
333 fCache.Get()->fresh =
false;
355 running[i]=running[i-1];
368 if ( nEnergies == 1 || running[nEnergies-1] == 0 )
369 fCache.Get()->currentMeanEnergy = 0.0;
372 fCache.Get()->currentMeanEnergy = weighted/running[nEnergies-1];
376 if ( nEnergies == 1 ) it = 0;
379 if ( running[nEnergies-1] != 0 )
384 if ( random < running [ i ] / running [ nEnergies-1 ] )
break;
389 if ( running [ nEnergies-1 ] == 0 ) it = 0;
413 running[it-1]/running[nEnergies-1],
414 running[it]/running[nEnergies-1],
432 G4double x1 = running[it-1]/running[nEnergies-1];
433 G4double x2 = running[it]/running[nEnergies-1];
457 else if(angularRep==2)
464 if( getenv(
"G4PHPTEST") )
G4cout <<
" G4ParticleHPContAngularPar::Sample nEnergies " <<
nEnergies <<
G4endl;
467 if(j!=0) running[j]=running[j-1];
474 if( getenv(
"G4PHPTEST") )
G4cout <<
" G4ParticleHPContAngularPar::Sample " << j <<
" running " << running[j]
480 if ( nEnergies == 1 )
481 fCache.Get()->currentMeanEnergy = 0.0;
483 fCache.Get()->currentMeanEnergy = weighted/running[nEnergies-1];
492 if(randkal<running[j]/running[nEnergies-1])
break;
498 x = randkal*running[nEnergies-1];
507 if( getenv(
"G4PHPTEST") )
G4cout << itt <<
" G4particleHPContAngularPar fsEnergy " << fsEnergy <<
" " <<
theManager.
GetInverseScheme(itt-1) <<
" x " << x <<
" " << x1 <<
" " << x2 <<
" y " << y1 <<
" " << y2 <<
G4endl;
512 fsEnergy, y1, y2, cLow,cHigh);
513 if( getenv(
"G4PHPTEST") )
G4cout << itt <<
" G4particleHPContAngularPar compoundFraction " << compoundFraction <<
" E " << fsEnergy <<
" " <<
theManager.
GetScheme(itt) <<
" ener " << fsEnergy <<
" y " << y1 <<
" " << y2 <<
" cLH " << cLow <<
" " << cHigh <<
G4endl;
525 targetA =
G4int (
fCache.Get()->theTarget->GetMass()/amu_c2 + 0.5 );
527 G4int residualA = targetA+1-
A;
528 G4int residualZ = targetZ-Z;
533 incidentEnergy, incidentMass,
534 productEnergy, productMass,
535 residualMass, residualA, residualZ,
536 targetMass, targetA, targetZ);
537 cosTh = theKallbach.
Sample(anEnergy);
538 if( getenv(
"G4PHPTEST") )
G4cout <<
" G4ParticleHPKallbachMannSyst::Sample resulttest " << cosTh <<
G4endl;
540 else if(angularRep>10&&angularRep<16)
548 if(i!=0) running[i]=running[i-1];
558 if ( nEnergies == 1 )
559 fCache.Get()->currentMeanEnergy = 0.0;
561 fCache.Get()->currentMeanEnergy = weighted/running[nEnergies-1];
564 if ( nEnergies == 1 ) it = 0;
569 if(random<running[i]/running[nEnergies-1])
break;
585 aMan.
Init(angularRep-10, nAngularParameters-1);
587 cosTh = theStore.
Sample();
603 running[it-1]/running[nEnergies-1],
604 running[it]/running[nEnergies-1],
609 cosTh = theStore.
Sample();
614 G4double x1 = running[it-1]/running[nEnergies-1];
615 G4double x2 = running[it]/running[nEnergies-1];
653 x = theBuff1.
GetX(i);
654 y1 = theBuff1.
GetY(i);
655 y2 = theBuff2.
GetY(i);
662 cosTh = theStore.
Sample();
668 throw G4HadronicException(__FILE__, __LINE__,
"G4ParticleHPContAngularPar::Sample: Unknown angular representation");
675 G4ThreeVector tempVector(mtot*sinth*std::cos(phi), mtot*sinth*std::sin(phi), mtot*std::cos(theta) );
692 if( !angParPrev )
return;
699 for(ie=0; ie<nDiscreteEnergiesPrev; ie++) {
704 for(ie=nDiscreteEnergies; ie<
nEnergies; ie++) {
713 for(ie=nDiscreteEnergiesPrev; ie<nEnergiesPrev; ie++) {
715 G4double enerT = (ener-minEnerPrev)/(maxEnerPrev-minEnerPrev);
728 G4int ie,ie1,ie2, ie1Prev, ie2Prev;
734 std::set<G4double>::const_iterator itede;
737 std::map<G4double,G4int>::const_iterator itedeo;
741 itedeo = discEnerOwn1.find(discEner);
742 if( itedeo == discEnerOwn1.end() ) {
747 itedeo = discEnerOwn2.find(discEner);
748 if( itedeo == discEnerOwn2.end() ) {
772 if( getenv(
"G4PHPTEST2") )
G4cout << ie <<
" " << ip <<
" G4ParticleHPContAngularPar::Merge DiscreteEnergies val1 " << val1 <<
" val2 " << val2 <<
" value " << value <<
G4endl;
790 std::set<G4double>::const_iterator iteet = energiesTransformed.begin();
803 G4double e1 = (maxEner1-minEner1) * eT + minEner1;
805 for(ie1=nDiscreteEnergies1; ie1<nEnergies1; ie1++) {
809 if( ie1 == 0 ) ie1Prev++;
810 if( ie1 == nEnergies1 ) {
815 G4double e2 = (maxEner2-minEner2) * eT + minEner2;
817 for(ie2=nDiscreteEnergies2; ie2<nEnergies2; ie2++) {
822 if( ie2 == 0 ) ie2Prev++;
823 if( ie2 == nEnergies2 ) {
830 if( getenv(
"G4PHPTEST2") )
G4cout << ie <<
" " << ie1 <<
" " << ie2 <<
" G4ParticleHPContAngularPar::loop eT " << eT <<
" -> eN " << eN <<
" e1 " << e1 <<
" e2 " << e2 <<
G4endl;
853 if( getenv(
"G4PHPTEST2") )
G4cout << ie <<
" " << ip <<
" G4ParticleHPContAngularPar::Merge val1 " << val1 <<
" val2 " << val2 <<
" value " << value <<
G4endl;
862 if( getenv(
"G4PHPTEST2") ) {
863 G4cout <<
" G4ParticleHPContAngularPar::Merge ANGPAR1 " <<
G4endl;
865 G4cout <<
" G4ParticleHPContAngularPar::Merge ANGPAR2 " <<
G4endl;
867 G4cout <<
" G4ParticleHPContAngularPar::Merge ANGPARNEW " <<
G4endl;
std::map< G4double, G4int > GetDiscreteEnergiesOwn() const
void Init(G4int i, G4double e, G4int n)
G4int GetVectorLength() const
G4double GetTotalMomentum() const
CLHEP::Hep3Vector G4ThreeVector
void Init(std::istream &aDataFile, G4int nPar, G4double unit=1.)
G4ReactionProduct * Sample(G4double anEnergy, G4double massCode, G4double mass, G4int angularRep, G4int interpol)
void BuildByInterpolation(G4double anEnergy, G4InterpolationScheme aScheme, G4ParticleHPContAngularPar &store1, G4ParticleHPContAngularPar &store2)
G4double Interpolate2(G4InterpolationScheme aScheme, G4double x, G4double x1, G4double x2, G4double y1, G4double y2) const
void Init(G4int aScheme, G4int aRange)
void SetKineticEnergy(const G4double en)
G4double GetMinEner() const
void SetMomentum(const G4double x, const G4double y, const G4double z)
G4double SampleMax(G4double energy)
G4Cache< toBeCached * > fCache
G4ParticleDefinition * theProjectile
G4InterpolationManager theManager
void SetInterpolationManager(const G4InterpolationManager &aManager)
void SetDefinition(const G4ParticleDefinition *aParticleDefinition)
G4ParticleHPContAngularPar()
std::map< G4double, G4int > theDiscreteEnergiesOwn
std::set< G4double > GetEnergiesTransformed() const
void SetY(G4int i, G4double x)
G4GLOB_DLL std::ostream G4cout
double A(double temperature)
void Init(std::istream &aDataFile, G4ParticleDefinition *projectile)
G4int GetNEnergiesTransformed() const
void SetValue(G4int i, G4double y)
G4InterpolationScheme GetScheme(G4int index) const
static const double twopi
static G4Triton * Triton()
static G4Proton * Proton()
G4double GetX(G4int i) const
static G4Neutron * Neutron()
G4double GetY(G4double x)
G4double GetMaxEner() const
G4double GetWeightedBinIntegral(const G4InterpolationScheme &aScheme, const G4double x1, const G4double x2, const G4double y1, const G4double y2)
static G4Deuteron * Deuteron()
G4double Interpolate(G4InterpolationScheme aScheme, G4double x, G4double x1, G4double x2, G4double y1, G4double y2) const
void SetLabel(G4double aLabel)
static G4IonTable * GetIonTable()
void SetManager(G4InterpolationManager &aManager)
static G4Positron * Positron()
G4int GetNDiscreteEnergies() const
G4ParticleHPInterpolator theInt
G4double GetPDGMass() const
static G4double GetBindingEnergy(const G4int A, const G4int Z)
T max(const T t1, const T t2)
brief Return the largest of the two arguments
G4ParticleHPList * theAngular
void SetCoeff(G4int i, G4int l, G4double coeff)
std::set< G4double > theDiscreteEnergies
const G4double x[NPOINTSGL]
std::set< G4double > theEnergiesTransformed
T min(const T t1, const T t2)
brief Return the smallest of the two arguments
G4InterpolationScheme GetInverseScheme(G4int index)
static G4Electron * Electron()
G4double GetValue(G4int i)
G4double GetBinIntegral(const G4InterpolationScheme &aScheme, const G4double x1, const G4double x2, const G4double y1, const G4double y2)
void PrepareTableInterpolation(const G4ParticleHPContAngularPar *angularPrev)
void SetX(G4int i, G4double e)
G4double Sample(G4double anEnergy)
G4int GetNEnergies() const