Abstract class. Provides a data caching mechanism for for concrete implementations of the EventRecordVisitorI interface, generating kinematics and wishing to cache maximum differential xsecs. More...

#include <KineGeneratorWithCache.h>

Inheritance diagram for genie::KineGeneratorWithCache:

[legend]

Collaboration diagram for genie::KineGeneratorWithCache:

[legend]

Protected Member Functions
	KineGeneratorWithCache ()

	KineGeneratorWithCache (string name)

	KineGeneratorWithCache (string name, string config)

	~KineGeneratorWithCache ()

virtual double	ComputeMaxXSec (const Interaction *in) const =0

virtual double	ComputeMaxXSec (const Interaction *in, const int nkey) const

virtual double	MaxXSec (GHepRecord *evrec, const int nkey=0) const

virtual double	FindMaxXSec (const Interaction *in, const int nkey=0) const

virtual void	CacheMaxXSec (const Interaction *in, double xsec, const int nkey=0) const

virtual double	Energy (const Interaction *in) const

virtual CacheBranchFx *	AccessCacheBranch (const Interaction *in, const int nkey=0) const

virtual void	AssertXSecLimits (const Interaction *in, double xsec, double xsec_max) const

Protected Member Functions inherited from genie::EventRecordVisitorI
	EventRecordVisitorI ()

	EventRecordVisitorI (string name)

	EventRecordVisitorI (string name, string config)

Protected Member Functions inherited from genie::Algorithm
	Algorithm ()

	Algorithm (string name)

	Algorithm (string name, string config)

void	Initialize (void)

void	DeleteConfig (void)

void	DeleteSubstructure (void)

Registry *	ExtractLocalConfig (const Registry &in) const

Registry *	ExtractLowerConfig (const Registry &in, const string &alg_key) const
	Split an incoming configuration Registry into a block valid for the sub-algo identified by alg_key. More...

template<class T >
bool	GetParam (const RgKey &name, T &p, bool is_top_call=true) const

template<class T >
bool	GetParamDef (const RgKey &name, T &p, const T &def) const

template<class T >
int	GetParamVect (const std::string &comm_name, std::vector< T > &v, bool is_top_call=true) const
	Handle to load vectors of parameters. More...

int	GetParamVectKeys (const std::string &comm_name, std::vector< RgKey > &k, bool is_top_call=true) const

template<class T >
int	GetParamMat (const std::string &comm_name, TMatrixT< T > &mat, bool is_top_call=true) const
	Handle to load matrix of parameters. More...

template<class T >
int	GetParamMatSym (const std::string &comm_name, TMatrixTSym< T > &mat, bool is_top_call=true) const

int	GetParamMatKeys (const std::string &comm_name, std::vector< RgKey > &k, bool is_top_call=true) const

int	AddTopRegistry (Registry *rp, bool owns=true)
	add registry with top priority, also update ownership More...

int	AddLowRegistry (Registry *rp, bool owns=true)
	add registry with lowest priority, also update ownership More...

int	MergeTopRegistry (const Registry &r)

int	AddTopRegisties (const vector< Registry * > &rs, bool owns=false)
	Add registries with top priority, also udated Ownerships. More...

Protected Attributes
const XSecAlgorithmI *	fXSecModel

double	fSafetyFactor
	ComputeMaxXSec -> ComputeMaxXSec * fSafetyFactor. More...

std::vector< double >	vSafetyFactors
	MaxXSec -> MaxXSec * fSafetyFactors[nkey]. More...

int	fNumOfSafetyFactors
	Number of given safety factors. More...

std::vector< string >	vInterpolatorTypes
	Type of interpolator for each key in a branch. More...

int	fNumOfInterpolatorTypes
	Number of given interpolators types. More...

double	fMaxXSecDiffTolerance
	max{100(xsec-maxxsec)/.5(xsec+maxxsec)} if xsec>maxxsec More...

double	fEMin
	min E for which maxxsec is cached - forcing explicit calc. More...

bool	fGenerateUniformly
	uniform over allowed phase space + event weight? More...

Protected Attributes inherited from genie::Algorithm
bool	fAllowReconfig

bool	fOwnsSubstruc
	true if it owns its substructure (sub-algs,...) More...

AlgId	fID
	algorithm name and configuration set More...

vector< Registry * >	fConfVect

vector< bool >	fOwnerships
	ownership for every registry in fConfVect More...

AlgStatus_t	fStatus
	algorithm execution status More...

AlgMap *	fOwnedSubAlgMp
	local pool for owned sub-algs (taken out of the factory pool) More...

Additional Inherited Members
Public Member Functions inherited from genie::EventRecordVisitorI
virtual	~EventRecordVisitorI ()

virtual void	ProcessEventRecord (GHepRecord *event_rec) const =0

Public Member Functions inherited from genie::Algorithm
virtual	~Algorithm ()

virtual void	Configure (const Registry &config)

virtual void	Configure (string config)

virtual void	FindConfig (void)

virtual const Registry &	GetConfig (void) const

Registry *	GetOwnedConfig (void)

virtual const AlgId &	Id (void) const
	Get algorithm ID. More...

virtual AlgStatus_t	GetStatus (void) const
	Get algorithm status. More...

virtual bool	AllowReconfig (void) const

virtual AlgCmp_t	Compare (const Algorithm *alg) const
	Compare with input algorithm. More...

virtual void	SetId (const AlgId &id)
	Set algorithm ID. More...

virtual void	SetId (string name, string config)

const Algorithm *	SubAlg (const RgKey &registry_key) const

void	AdoptConfig (void)

void	AdoptSubstructure (void)

virtual void	Print (ostream &stream) const
	Print algorithm info. More...

Static Public Member Functions inherited from genie::Algorithm
static string	BuildParamVectKey (const std::string &comm_name, unsigned int i)

static string	BuildParamVectSizeKey (const std::string &comm_name)

static string	BuildParamMatKey (const std::string &comm_name, unsigned int i, unsigned int j)

static string	BuildParamMatRowSizeKey (const std::string &comm_name)

static string	BuildParamMatColSizeKey (const std::string &comm_name)

Detailed Description

Abstract class. Provides a data caching mechanism for for concrete implementations of the EventRecordVisitorI interface, generating kinematics and wishing to cache maximum differential xsecs.

      This class provides some common implementation for handling
      (retrieving, creating, searching, adding to) the cache.
      The various super-classes should implement the ComputeMaxXSec(...)
      method for computing the maximum xsec in case it has not already
      being pushed into the cache at a previous iteration. \n

      Update May 15, 2022 IK: 
      It makes possible to cache several values having different keys. 
      The example of using this opportunity see in 
      the class QELEventGeneratorSM.

Author: Costas Andreopoulos <c.andreopoulos cern.ch> University of Liverpool
Igor Kakorin kakor.nosp@m.in@j.nosp@m.inr.r.nosp@m.u Joint Institute for Nuclear Research

Created:: December 15, 2004

License:: Copyright (c) 2003-2024, The GENIE Collaboration For the full text of the license visit http://copyright.genie-mc.org

Definition at line 49 of file KineGeneratorWithCache.h.

Constructor & Destructor Documentation

KineGeneratorWithCache::KineGeneratorWithCache ( )

protected

Definition at line 34 of file KineGeneratorWithCache.cxx.

                                                : 
 EventRecordVisitorI(), fSafetyFactor(1.), fNumOfSafetyFactors(-1), fNumOfInterpolatorTypes(-1)
 {
 
 }

KineGeneratorWithCache::KineGeneratorWithCache ( string name )

protected

Definition at line 40 of file KineGeneratorWithCache.cxx.

                                                           : 
 EventRecordVisitorI(name), fSafetyFactor(1.), fNumOfSafetyFactors(-1), fNumOfInterpolatorTypes(-1)
 {
 
 }

KineGeneratorWithCache::KineGeneratorWithCache	(	string	name,
		string	config
	)

protected

Definition at line 46 of file KineGeneratorWithCache.cxx.

                                                                          : 
 EventRecordVisitorI(name, config), fSafetyFactor(1.), fNumOfSafetyFactors(-1), fNumOfInterpolatorTypes(-1)
 {
 
 }

KineGeneratorWithCache::~KineGeneratorWithCache ( )

protected

Definition at line 52 of file KineGeneratorWithCache.cxx.

 {
 
 }

Member Function Documentation

CacheBranchFx * KineGeneratorWithCache::AccessCacheBranch	(	const Interaction *	in,
		const int	nkey = `0`
	)		const

protectedvirtual

Definition at line 224 of file KineGeneratorWithCache.cxx.

References genie::Cache::AddCacheBranch(), genie::Interaction::AsString(), genie::Cache::CacheBranchKey(), genie::Cache::FindCacheBranch(), genie::Algorithm::Id(), genie::Cache::Instance(), genie::AlgId::Key(), LOG, and pINFO.

Referenced by CacheMaxXSec(), and FindMaxXSec().

 {
 // Returns the cache branch for this algorithm and this interaction. If no
 // branch is found then one is created.
 
   Cache * cache = Cache::Instance();
 
   // build the cache branch key as: namespace::algorithm/config/interaction/nkey
   string algkey = this->Id().Key();
   string intkey = interaction->AsString();
   string key    = cache->CacheBranchKey(algkey, intkey, std::to_string(nkey));
 
   CacheBranchFx * cache_branch =
               dynamic_cast<CacheBranchFx *> (cache->FindCacheBranch(key));
   if(!cache_branch) {
     //-- create the cache branch at the first pass
     LOG("Kinematics", pINFO) << "No cache branch found";
     LOG("Kinematics", pINFO) << "Creating cache branch - key = " << key;
 
     cache_branch = new CacheBranchFx("Max over phase space");
     cache->AddCacheBranch(key, cache_branch);
   }
   assert(cache_branch);
 
   return cache_branch;
 }

void KineGeneratorWithCache::AssertXSecLimits	(	const Interaction *	in,
		double	xsec,
		double	xsec_max
	)		const

protectedvirtual

Definition at line 252 of file KineGeneratorWithCache.cxx.

References fMaxXSecDiffTolerance, LOG, pERROR, pFATAL, and pWARN.

 {
   // check the computed cross section for the current kinematics against the
   // maximum cross section used in the rejection MC method for the current
   // interaction at the current energy.
   if(xsec>xsec_max) {
     double f = 200*(xsec-xsec_max)/(xsec_max+xsec);
     if(f>fMaxXSecDiffTolerance) {
        LOG("Kinematics", pFATAL)
           << "xsec: (curr) = " << xsec
                  << " > (max) = " << xsec_max << "\n for " << *interaction;
        LOG("Kinematics", pFATAL)
           << "*** Exceeding estimated maximum differential cross section";
        std::terminate();
     } else {
        LOG("Kinematics", pWARN)
           << "xsec: (curr) = " << xsec
                  << " > (max) = " << xsec_max << "\n for " << *interaction;
        LOG("Kinematics", pWARN)
             << "*** The fractional deviation of " << f << " % was allowed";
     }
   }
 
   // this should never happen - print an error mesg just in case...
   if(xsec<0) {
     LOG("Kinematics", pERROR)
      << "Negative cross section for current kinematics!! \n" << *interaction;
   }
 }

void KineGeneratorWithCache::CacheMaxXSec	(	const Interaction *	in,
		double	xsec,
		const int	nkey = `0`
	)		const

protectedvirtual

Definition at line 190 of file KineGeneratorWithCache.cxx.

References AccessCacheBranch(), genie::CacheBranchFx::AddValues(), genie::CacheBranchFx::CreateSpline(), Energy(), fEMin, fNumOfInterpolatorTypes, LOG, genie::CacheBranchFx::Map(), pINFO, genie::CacheBranchFx::Spl(), vInterpolatorTypes, and genie::Spline::XMax().

Referenced by MaxXSec().

 {
   LOG("Kinematics", pINFO)
                        << "Adding the computed max value to cache";
   CacheBranchFx * cb = this->AccessCacheBranch(interaction, nkey);
 
   double E = this->Energy(interaction);
   if (E<fEMin) return;
   if(max_xsec>0) cb->AddValues(E,max_xsec);
 
   if(! cb->Spl() ) {
     if( cb->Map().size() > 40 ) 
       cb->CreateSpline(nkey<=fNumOfInterpolatorTypes-1?vInterpolatorTypes[nkey]:"");
   }
 
   if( cb->Spl() ) {
      if( E < cb->Spl()->XMin() || E > cb->Spl()->XMax() ) {
         cb->CreateSpline(nkey<=fNumOfInterpolatorTypes-1?vInterpolatorTypes[nkey]:"");
      }
   }
 }

virtual double genie::KineGeneratorWithCache::ComputeMaxXSec ( const Interaction * in ) const

protectedpure virtual

double KineGeneratorWithCache::ComputeMaxXSec	(	const Interaction *	in,
		const int	nkey
	)		const

protectedvirtual

Reimplemented in genie::QELEventGeneratorSM.

Definition at line 283 of file KineGeneratorWithCache.cxx.

References ComputeMaxXSec().

 {
   if (nkey == 0)
   {
      return this->ComputeMaxXSec(in);
   }
   else
   {
     return -1;
   }
 }

double KineGeneratorWithCache::Energy ( const Interaction * in ) const

protectedvirtual

Reimplemented in genie::COHKinematicsGenerator, genie::SKKinematicsGenerator, genie::DMEKinematicsGenerator, genie::NuEKinematicsGenerator, and genie::HELeptonKinematicsGenerator.

Definition at line 213 of file KineGeneratorWithCache.cxx.

References genie::Interaction::InitState(), genie::kRfHitNucRest, and genie::InitialState::ProbeE().

Referenced by CacheMaxXSec(), and FindMaxXSec().

 {
 // Returns the neutrino energy at the struck nucleon rest frame. Kinematic
 // generators should override this method if they need to cache the max-xsec
 // values for another energy value (eg kinematic generators for IMD or COH)
 
   const InitialState & init_state = interaction->InitState();
   double E = init_state.ProbeE(kRfHitNucRest);
   return E;
 }

double KineGeneratorWithCache::FindMaxXSec	(	const Interaction *	in,
		const int	nkey = `0`
	)		const

protectedvirtual

Definition at line 105 of file KineGeneratorWithCache.cxx.

References AccessCacheBranch(), Energy(), genie::Spline::Evaluate(), fEMin, LOG, genie::CacheBranchFx::Map(), pINFO, genie::CacheBranchFx::Spl(), genie::Spline::XMax(), and genie::Spline::XMin().

Referenced by MaxXSec().

 {
 // Find a cached max xsec for the specified xsec algorithm & interaction and
 // close to the specified energy
 
   // get neutrino energy
   double E = this->Energy(interaction);
   LOG("Kinematics", pINFO) << "E = " << E;
 
   if(E < fEMin) {
      LOG("Kinematics", pINFO)
          << "Below minimum energy - Forcing explicit calculation";
      return -1.;
   }
 
   // access the the cache branch
   CacheBranchFx * cb = this->AccessCacheBranch(interaction, nkey);
 
   // if there are enough points stored in the cache buffer to build a
   // spline, then intepolate
   if( cb->Spl() ) {
      if( E >= cb->Spl()->XMin() && E <= cb->Spl()->XMax()) {
        double spl_max_xsec = cb->Spl()->Evaluate(E);
        LOG("Kinematics", pINFO)
           << "\nInterpolated: max (E=" << E << ") = " << spl_max_xsec;
        return spl_max_xsec;
      }
      LOG("Kinematics", pINFO)
           << "Outside spline boundaries - Forcing explicit calculation";
      return -1.;
   }
 
   // if there are not enough points at the cache buffer to have a spline,
   // look whether there is another point that is sufficiently close
   double dE = TMath::Min(0.25, 0.05*E);
   const map<double,double> & fmap = cb->Map();
   map<double,double>::const_iterator iter = fmap.lower_bound(E);
   if(iter != fmap.end()) {
      if(TMath::Abs(E - iter->first) < dE) return iter->second;
   }
 
   return -1;
 
 /*
   // build the search rule
   double dE = TMath::Min(0.25, 0.05*E);
   ostringstream search;
   search << "(x-" << E << " < " << dE << ") && (x>=" << E << ")";
 
   // query for all the entries at a window around the current energy
   TSQLResult * result = cb->Ntuple()->Query("x:y", search.str().c_str());
   int nrows = result->GetRowCount();
   LOG("Kinematics", pDEBUG)
             << "Found " << nrows << " rows with " << search.str();
   if(nrows <= 0) {
      delete result;
      return -1;
   }
 
   // and now select the entry with the closest energy
   double max_xsec = -1.0;
   double Ep       = 0;
   double dEmin    = 999;
   TSQLRow * row   = 0;
   while( (row = result->Next()) ) {
      double cE    = atof( row->GetField(0) );
      double cxsec = atof( row->GetField(1) );
      double dE    = TMath::Abs(E-cE);
      if(dE < dEmin) {
         max_xsec = cxsec;
         Ep       = cE;
         dEmin    = TMath::Min(dE,dEmin);
      }
      delete row;
   }
   delete result;
 
   LOG("Kinematics", pINFO)
      << "\nRetrieved: max xsec = " << max_xsec << " cached at E = " << Ep;
 
   return max_xsec;
 */
 }

double KineGeneratorWithCache::MaxXSec	(	GHepRecord *	evrec,
		const int	nkey = `0`
	)		const

protectedvirtual

Definition at line 57 of file KineGeneratorWithCache.cxx.

References CacheMaxXSec(), ComputeMaxXSec(), genie::GHepRecord::EventFlags(), FindMaxXSec(), fNumOfSafetyFactors, genie::kISkipKinematicChk, genie::kISkipProcessChk, genie::kKineGenErr, genie::kPSNull, LOG, pINFO, pNOTICE, genie::GHepRecord::SetDiffXSec(), genie::exceptions::EVGThreadException::SetReason(), genie::GHepRecord::Summary(), genie::exceptions::EVGThreadException::SwitchOnFastForward(), and vSafetyFactors.

 {
   LOG("Kinematics", pINFO)
                 << "Getting max. for the rejection method";
 
   double xsec_max = -1;
   Interaction * interaction = event_rec->Summary();
 
   LOG("Kinematics", pINFO)
                   << "Attempting to find a cached max value";
   xsec_max = this->FindMaxXSec(interaction, nkey);
   if(xsec_max>0) return nkey<=fNumOfSafetyFactors-1?vSafetyFactors[nkey]*xsec_max:xsec_max;
 
   LOG("Kinematics", pINFO)
                   << "Attempting to compute the max value";
   if (nkey == 0)
   {
     xsec_max = this->ComputeMaxXSec(interaction);
   }
   else
   {
     xsec_max = this->ComputeMaxXSec(interaction, nkey);
   }
   
   if(xsec_max>0) {
      LOG("Kinematics", pINFO) << "max = " << xsec_max;
      this->CacheMaxXSec(interaction, xsec_max, nkey);
      return nkey<=fNumOfSafetyFactors-1?vSafetyFactors[nkey]*xsec_max:xsec_max;
   }
 
   LOG("Kinematics", pNOTICE)
             << "Can not generate event kinematics max_xsec<=0)";
   // xsec for selected kinematics = 0
   event_rec->SetDiffXSec(0,kPSNull);
   // switch on error flag
   event_rec->EventFlags()->SetBitNumber(kKineGenErr, true);
   // reset 'trust' bits
   interaction->ResetBit(kISkipProcessChk);
   interaction->ResetBit(kISkipKinematicChk);
   // throw exception
   genie::exceptions::EVGThreadException exception;
   exception.SetReason("kinematics generation: max_xsec<=0");
   exception.SwitchOnFastForward();
   throw exception;
 
   return 0;
 }

Member Data Documentation

double genie::KineGeneratorWithCache::fEMin

protected

min E for which maxxsec is cached - forcing explicit calc.

Definition at line 76 of file KineGeneratorWithCache.h.

bool genie::KineGeneratorWithCache::fGenerateUniformly

protected

uniform over allowed phase space + event weight?

Definition at line 77 of file KineGeneratorWithCache.h.

double genie::KineGeneratorWithCache::fMaxXSecDiffTolerance

protected

max{100*(xsec-maxxsec)/.5*(xsec+maxxsec)} if xsec>maxxsec

Definition at line 75 of file KineGeneratorWithCache.h.

int genie::KineGeneratorWithCache::fNumOfInterpolatorTypes

protected

Number of given interpolators types.

Definition at line 74 of file KineGeneratorWithCache.h.

Referenced by CacheMaxXSec(), and genie::QELEventGeneratorSM::LoadConfig().

int genie::KineGeneratorWithCache::fNumOfSafetyFactors

protected

Number of given safety factors.

Definition at line 72 of file KineGeneratorWithCache.h.

Referenced by genie::QELEventGeneratorSM::LoadConfig(), and MaxXSec().

double genie::KineGeneratorWithCache::fSafetyFactor

protected

ComputeMaxXSec -> ComputeMaxXSec * fSafetyFactor.

Definition at line 70 of file KineGeneratorWithCache.h.

const XSecAlgorithmI* genie::KineGeneratorWithCache::fXSecModel

mutableprotected

Definition at line 68 of file KineGeneratorWithCache.h.

std::vector<string> genie::KineGeneratorWithCache::vInterpolatorTypes

protected

Type of interpolator for each key in a branch.

Definition at line 73 of file KineGeneratorWithCache.h.

Referenced by CacheMaxXSec(), and genie::QELEventGeneratorSM::LoadConfig().

std::vector<double> genie::KineGeneratorWithCache::vSafetyFactors

protected

MaxXSec -> MaxXSec * fSafetyFactors[nkey].

Definition at line 71 of file KineGeneratorWithCache.h.

Referenced by genie::QELEventGeneratorSM::LoadConfig(), and MaxXSec().

The documentation for this class was generated from the following files:

Protected Member Functions

Protected Attributes

Additional Inherited Members

Detailed Description

Constructor & Destructor Documentation

Member Function Documentation

Member Data Documentation