#include <G4FastPathHadronicCrossSection.hh>

Public Member Functions
	fastPathEntry (const G4ParticleDefinition par, const G4Material mat, G4double min_cutoff)

	~fastPathEntry ()

G4double	GetCrossSection (G4double ene) const

void	Initialize (G4CrossSectionDataStore *)

Public Attributes
const G4ParticleDefinition *const	particle

const G4Material *const	material

const G4double	min_cutoff

XSParam *	physicsVector

Detailed Description

Definition at line 54 of file G4FastPathHadronicCrossSection.hh.

Constructor & Destructor Documentation

◆ fastPathEntry()

fastPathEntry::fastPathEntry	(	const G4ParticleDefinition *	par,
		const G4Material *	mat,
		G4double	min_cutoff
	)

Definition at line 65 of file G4FastPathHadronicCrossSection.cc.

                                                                                                  :
        particle(part),material(mat),min_cutoff(min),physicsVector(nullptr)
    {
    DBG("Initializing a fastPathEntry");
#ifdef FPDEBUG
    count = 0;
    slowpath_sum=0.;
    max_delta=0.;
    min_delta=0.;
    sum_delta=0.;
    sum_delta_square=0.;
#endif
}

◆ ~fastPathEntry()

fastPathEntry::~fastPathEntry ( )

Definition at line 79 of file G4FastPathHadronicCrossSection.cc.

{
    DBG("Deleting fastPathEntry");
    DBG("Dumping status for: "<<(particle?particle->GetParticleName():"PART_NONE")<<" "\
          <<(material?material->GetName():"MAT_NONE")<<" min_cutoff:"<<min_cutoff<<" "\
          <<" count:"<<count<<" slowpath_sum:"<<slowpath_sum<<" max_delta:"<<max_delta\
          <<" min_delta"<<min_delta<<" sum_delta"<<sum_delta<<" sum_delta_squared:"<<sum_delta_square);
    delete physicsVector;
}

Member Function Documentation

◆ GetCrossSection()

G4double G4FastPathHadronicCrossSection::fastPathEntry::GetCrossSection ( G4double ene ) const

inline

Definition at line 58 of file G4FastPathHadronicCrossSection.hh.

58{ return physicsVector->Value(ene); }

G4PhysicsVector::Value

G4double Value(G4double theEnergy, std::size_t &lastidx) const

Definition: G4PhysicsVector.cc:519

◆ Initialize()

void fastPathEntry::Initialize ( G4CrossSectionDataStore * xsds )

Definition at line 95 of file G4FastPathHadronicCrossSection.cc.

{
    //Check this method is called when G4CrossSectionDataStore is in the correct state:
    //  FastPath is enabled and we are indeed initializing
    assert( xsds->GetFastPathControlFlags().useFastPathIfAvailable &&
            xsds->GetFastPathControlFlags().initializationPhase );
    using std::log10;
    std::vector<Point_t> data_in;
    const fastPathParameters& params = xsds->GetFastPathParameters();
    G4double xs;
 
      //G4double max_query = params.queryMax;
      //G4int count = sampleCount;
      //G4double tol = dpTol;
 
    //Shift so max and min are >= 1.
    //Don't forget to shift back before computing XS
    G4double min = params.sampleMin;
    G4double max = params.sampleMax;
    G4double shift = 0.0;
    if(min < 1.0){
        shift = 1.0 - min;
      }
    min += shift;
    max += shift;
 
    G4double log_max = std::log10(params.sampleMax);
    G4double log_min = std::log10(params.sampleMin);
    G4double log_step = (log_max-log_min)/(1.0*params.sampleCount);
 
    G4double max_xs = 0.0;
 
    //Utility particle to calculate XS, with 0 kin energy by default
    static const G4ThreeVector constDirection(0.,0.,1.);
    G4DynamicParticle* probingParticle = new G4DynamicParticle( particle , constDirection , 0 );
 
    //add the cutoff energy
    probingParticle->SetKineticEnergy(min_cutoff);
    //Sample cross-section
    xs = xsds->GetCrossSection(probingParticle,material);
    data_in.push_back({min_cutoff,xs});
 
    G4double currEnergy = 0.0;
    //log results
    auto exp10 = [](G4double x){ return std::exp( M_LN10*x); };
    for(G4double log_currEnergy = log_min; log_currEnergy < log_max; log_currEnergy += log_step){
        currEnergy = exp10(log_currEnergy) - shift;
        if (currEnergy <  min_cutoff) continue;
        probingParticle->SetKineticEnergy(currEnergy);
        xs=xsds->GetCrossSection(probingParticle,material);
        //G4cout << "PRUTH: energy value " << currEnergy << ", XS value " << xs << G4endl;
        if (xs > max_xs) max_xs = xs;
        data_in.push_back({currEnergy,xs});
    } // --- end of loop i
    probingParticle->SetKineticEnergy(max-shift);
    xs = xsds->GetCrossSection(probingParticle,material);
    data_in.push_back({max-shift,xs});
 
    G4double tol = max_xs * 0.01;
    std::vector<Point_t> decimated_data;
    simplify_function(tol,  data_in,  decimated_data);
    std::vector<Point_t> debiased_data;
    RemoveBias( data_in,  decimated_data,  debiased_data);
    if ( physicsVector != nullptr ) delete physicsVector;
    physicsVector = new XSParam(decimated_data.size());
    G4int physicsVectorIndex = 0;
    for(size_t i = 0; i < decimated_data.size(); i++){
        physicsVector->PutValue(physicsVectorIndex++, decimated_data[i].e, decimated_data[i].xs);
    }
    //xsds->DumpFastPath(particle,material,G4cout);
}