INukeOsetTable Class Reference

Table-based implementation of Oset model. More...

#include <INukeOsetTable.h>

Inheritance diagram for INukeOsetTable:

Collaboration diagram for INukeOsetTable:

Classes
struct	PointHandler
	handle table's index and weights for given density and energy More...

Public Member Functions
	INukeOsetTable (const char *filename)
	constructor
void	setupOset (const double &density, const double &pionTk, const int &pionPDG, const double &protonFraction)
	use to set up Oset class (assign pion Tk, nuclear density etc)
Public Member Functions inherited from INukeOset
	INukeOset ()
	contructor
double	getTotalCrossSection () const
	return total = (qel+cex+abs) cross section
double	getCexCrossSection () const
	return cex cross section
double	getAbsorptionCrossSection () const
	return absorption cross section
double	getCexFraction () const
	return fraction of cex events
double	getAbsorptionFraction () const
	return fraction of absorption events

Private Member Functions
double	interpolate (const std::vector< double > &data) const
	interpolate cross section (method fixed for Oset tables)
int	processLine (const std::string &line)
	process single line from table file, push values to proper vector (method fixed for Oset tables)
int	checkIntegrity (const double &densityValue, const double &energyValue)
	check if data in file is consistent (method fixed for Oset tables)
void	badFile (const char *file, const int &errorCode, const int &line=0) const
	stop program and through an error if input file is corrupted (method fixed for Oset tables)
void	setCrossSections ()
	calculalte cross sections for each channel

Private Attributes
std::vector< double >	fQelCrossSectionTable [fNChannels]
	quasi-elastic piN cross section
std::vector< double >	fCexCrossSectionTable [fNChannels]
	charge-exchange piN cross section
std::vector< double >	fAbsorptionCrossSectionTable
	pi absorption cross section
unsigned int	fNDensityBins
	number of denisty bins
unsigned int	fNEnergyBins
	number of energy bins
double	fDensityBinWidth
	density step (must be fixed)
double	fEnergyBinWidth
	energy step (must be fixed)
PointHandler	fDensityHandler
	nuclear density handler
PointHandler	fEnergyHandler
	pion kinetic energy handler

Additional Inherited Members
Protected Member Functions inherited from INukeOset
void	setCrossSections (const int &pionPDG, const double &protonFraction)
	calculate avg cross sections according to proton / neutron fraction
Protected Attributes inherited from INukeOset
double	fNuclearDensity
	nuclear density in fm-3
double	fPionKineticEnergy
	pion kinetic energy in MeV
double	fTotalCrossSection
	el+cex+abs cross section (averaged over proton / neutron fraction)
double	fCexCrossSection
	cex cross section (averaged over proton / neutron fraction)
double	fAbsorptionCrossSection
	absorption cross section (averaged over proton / neutron fraction)
double	fQelCrossSections [fNChannels]
	total qel (el+cex) cross section for each channel
double	fCexCrossSections [fNChannels]
	cex cross section for each channel
Static Protected Attributes inherited from INukeOset
static const unsigned int	fNChannels = 3
	number of possible channels: pi+n, pi+p, pi0

Detailed Description

Table-based implementation of Oset model.

Author

Tomasz Golan

Date

2015

Warning

Applicable for pion with Tk < 350 MeV

Remarks

Tables taken from NuWro's implementation of: E. Oset et al., Nucl. Phys. A484 (1988) 557-592

This implementation is kept from historical reasons. Due to different approach to cascade in GENIE and NuWro there are some normalization issues. Default Oset model can be found in INukeOsetFormula

Definition at line 25 of file INukeOsetTable.h.

Constructor & Destructor Documentation

INukeOsetTable()

INukeOsetTable::INukeOsetTable

(

const char *

filename

)

constructor

load tables from file and check its integrity

Definition at line 10 of file INukeOsetTable.cxx.

                                                      : fNDensityBins (0), fNEnergyBins (0), 
                                                          fDensityBinWidth (0.0), fEnergyBinWidth (0.0)
{
  // open file with Oset table
  std::ifstream fileWithTables (filename);
  // check if file is open
  if (not fileWithTables.is_open()) badFile (filename, 1);
  // temporary string for getline
  std::string tempLine;
  // line counter (used in the case of error)
  int lineCounter = 0;
 
  // skip comments lines at the beginning of the file
  while (fileWithTables.get() == '#' and ++lineCounter)
    getline (fileWithTables, tempLine);
  // process other lines
  while (getline (fileWithTables, tempLine) and ++lineCounter)
    if (int error = processLine (tempLine))   // get exit code
      badFile (filename, error, lineCounter); // stop if exit code != 0
 
  // reset counter (in the case other file will be loaded)
  checkIntegrity (-1.0, -1.0);
 
  fDensityHandler.setHandler (fDensityBinWidth, fNDensityBins);
  fEnergyHandler.setHandler  (fEnergyBinWidth,  fNEnergyBins);
}

References badFile(), checkIntegrity(), fDensityBinWidth, fDensityHandler, fEnergyBinWidth, fEnergyHandler, fNDensityBins, fNEnergyBins, and processLine().

Member Function Documentation

badFile()

void INukeOsetTable::badFile ( const char * file, const int & errorCode, const int & line = 0 ) const

private

stop program and through an error if input file is corrupted (method fixed for Oset tables)

possible errors:

• the file with tables could not be loaded
• density bin width is not constant
• there is different number of energy point per density value
• energy bin width is not constant

Definition at line 141 of file INukeOsetTable.cxx.

{
  std::cerr << "\nERROR: " << file << " is corrupted (";
  
  switch (errorCode)
  {
    case 1: std::cerr << "could not be opened).\n"; exit(errorCode);
    case 2: std::cerr << "wrong density"; break;
    case 3: std::cerr << "different number of energy bins per density"; break;
    case 4: std::cerr << "wrong energy"; break;
    default: std::cerr << "undefined error";
  }
 
  std::cerr << " in line " << line << ").\n";
  exit(errorCode);
}

Referenced by INukeOsetTable().

checkIntegrity()

int INukeOsetTable::checkIntegrity ( const double & densityValue, const double & energyValue )

private

check if data in file is consistent (method fixed for Oset tables)

return 0 if all bins widths are constistent or return error code

Definition at line 80 of file INukeOsetTable.cxx.

{
  static unsigned int energyBinCounter = 0; // #energy bins for current density
  
  if (densityValue < 0) // checkIntegrity(-1.0, -1.0) is called to reset counter
    energyBinCounter = -1;
  else if (not isEqual (fNuclearDensity, densityValue)) // on new density value
  {
    fNDensityBins++; // update nof density bins
 
    if (fNuclearDensity >= 0) // skip first entry (density = -1.0)
    {
      // calculate current bin width
      double currentDensityBinWidth = densityValue - fNuclearDensity;
 
      // assign first density bin width
      // or check if next bins are the same as first
      // stop program with corresponding error message (exit code 2) if not ==
      if (fDensityBinWidth == 0.0) fDensityBinWidth = currentDensityBinWidth;
      else if (not isEqual (fDensityBinWidth, currentDensityBinWidth)) return 2;
      
      // assign energy bins for first density value
      // or check if next bins are the same as first
      if (fNEnergyBins == 0) fNEnergyBins = energyBinCounter;
      else if (fNEnergyBins != energyBinCounter) return 3; // exit code 3
 
      // reset counter and last energy
      energyBinCounter = 0;
      fPionKineticEnergy  = -1;
    }
  }
  else if (fPionKineticEnergy >= 0) // check energy width consistency (skip first)
  {
    // calculate current energy bin width
    double currentEnergyBinWidth = energyValue - fPionKineticEnergy;
 
    // assign first energy bin width
    // or check if next bins are the same as first
    // stop program with corresponding error message (exit code 4) if not ==
    if (fEnergyBinWidth == 0.0) fEnergyBinWidth = currentEnergyBinWidth;
    else if (not isEqual(fEnergyBinWidth, currentEnergyBinWidth)) return 4;
  }
 
  // increase energy bin counter (for current density)
  energyBinCounter++;
  // update last values
  fPionKineticEnergy  = energyValue;
  fNuclearDensity = densityValue;
 
  return 0; // no errors
}

References fDensityBinWidth, fEnergyBinWidth, fNDensityBins, fNEnergyBins, INukeOset::fNuclearDensity, INukeOset::fPionKineticEnergy, and osetUtils::isEqual().

Referenced by INukeOsetTable(), and processLine().

interpolate()

double INukeOsetTable::interpolate ( const std::vector< double > & data ) const

private

interpolate cross section (method fixed for Oset tables)

make bilinear interpolation between four points around (density, energy)

Definition at line 159 of file INukeOsetTable.cxx.

{
  // take four points adjacent to (density, energy) = (d,E):
  // (d0, E0), (d1, E0), (d0, E1), (d1, E1)
  // where d0 < d < d1, E0 < E < E1
  // each point goes in with weight = proper distance
 
  // total bin width for normalization
  static const double totalBinWidth = fDensityBinWidth * fEnergyBinWidth;
  
  // extract low boundary value from table
  double lowValue = data[fEnergyHandler.index +
                         fDensityHandler.index * fNEnergyBins]; // (d0, E0)
  // high boundaries = low boundary if on edge
  double highDensityValue = lowValue; // (d1, E0)
  double highEnergyValue  = lowValue; // (d0, E1)
  double highValue        = lowValue; // (d1, E1)
 
  // extract high boundaries values from table (if not on edge)
  if (not fDensityHandler.isEdge)
  {
    highDensityValue = data[fEnergyHandler.index +
                            (fDensityHandler.index + 1) * fNEnergyBins];
 
    if (not fEnergyHandler.isEdge)
    {
      highEnergyValue = data[fEnergyHandler.index + 1 +
                             fDensityHandler.index * fNEnergyBins];
      highValue       = data[fEnergyHandler.index + 1 +
                             (fDensityHandler.index + 1) * fNEnergyBins];
    }
  }
  else if (not fEnergyHandler.isEdge)
    highEnergyValue = data[fEnergyHandler.index + 1 +
                           fDensityHandler.index * fNEnergyBins];
 
  lowValue         *= fDensityHandler.lowWeight  * fEnergyHandler.lowWeight;
  highDensityValue *= fDensityHandler.highWeight * fEnergyHandler.lowWeight;
  highEnergyValue  *= fDensityHandler.lowWeight  * fEnergyHandler.highWeight;
  highValue        *= fDensityHandler.highWeight * fEnergyHandler.highWeight;
  
  return (lowValue + highDensityValue + highEnergyValue + highValue) /
          totalBinWidth;
}

References fDensityBinWidth, fDensityHandler, fEnergyBinWidth, fEnergyHandler, and fNEnergyBins.

Referenced by setCrossSections().

processLine()

int INukeOsetTable::processLine ( const std::string & line )

private

process single line from table file, push values to proper vector (method fixed for Oset tables)

• split line read from the file; the following structure is used:
  
  d; Tk; pi+n: tot fabs fcex; pi+p: tot fabs fcex; pi0: tot fabs fcex
  
  
• check its integrity
• save values into proper variables

Definition at line 47 of file INukeOsetTable.cxx.

{
  std::istringstream splitLine (line); // use istringstream to split string
  // line = d; KE; pi+n: tot fabs fcex; pi+p: tot fabs fcex; pi0: tot fabs fcex
 
  double currentDensityValue, currentEnergyValue; // will be extracted from line
 
  // get current density and energy value
  splitLine >> currentDensityValue >> currentEnergyValue;
 
  if (int exitCode = checkIntegrity (currentDensityValue, currentEnergyValue))
    return exitCode; // stop in the case of error (exit code != 0)
 
  for (unsigned int i = 0; i < fNChannels; i++) // channel loop
  {
    // get qel and cex cross sections for i-th channel
    double xsecQel, xsecCex;
    splitLine >> xsecQel >> xsecCex;
    // save them in proper vector
    fQelCrossSectionTable[i].push_back (xsecQel);
    fCexCrossSectionTable[i].push_back (xsecCex);
  } // channel loop
 
  // get absorption cross section
  double absorption;
  splitLine >> absorption;
  // save it in proper vector
  fAbsorptionCrossSectionTable.push_back (absorption);
 
  return 0; // no errors
}

References checkIntegrity(), fAbsorptionCrossSectionTable, fCexCrossSectionTable, INukeOset::fNChannels, and fQelCrossSectionTable.

Referenced by INukeOsetTable().

setCrossSections()

void INukeOsetTable::setCrossSections ( )

privatevirtual

calculalte cross sections for each channel

assign cross sections values to proper variables using bilinear interpolation

Implements INukeOset.

Definition at line 243 of file INukeOsetTable.cxx.

{
    for (unsigned int i = 0; i < fNChannels; i++) // channel loop
    {
      fQelCrossSections[i] = interpolate (fQelCrossSectionTable[i]);
      fCexCrossSections[i] = interpolate (fCexCrossSectionTable[i]);
    }
 
    fAbsorptionCrossSection = interpolate (fAbsorptionCrossSectionTable); 
}

References INukeOset::fAbsorptionCrossSection, fAbsorptionCrossSectionTable, INukeOset::fCexCrossSections, fCexCrossSectionTable, INukeOset::fNChannels, INukeOset::fQelCrossSections, fQelCrossSectionTable, and interpolate().

Referenced by setupOset().

setupOset()

void INukeOsetTable::setupOset ( const double & density, const double & pionTk, const int & pionPDG, const double & protonFraction )

virtual

use to set up Oset class (assign pion Tk, nuclear density etc)

• set up density, pion Tk and their handlers (for interpolation)
• get interpolated cross sections
• set up proper cross section variables

Implements INukeOset.

Definition at line 229 of file INukeOsetTable.cxx.

{
  fNuclearDensity    = density;
  fPionKineticEnergy = pionTk;
  fDensityHandler.update (fNuclearDensity);    // update density index / weights
  fEnergyHandler.update  (fPionKineticEnergy); // update energy index / weights
  setCrossSections();
  INukeOset::setCrossSections (pionPDG, protonFraction);  
}

References fDensityHandler, fEnergyHandler, INukeOset::fNuclearDensity, INukeOset::fPionKineticEnergy, INukeOset::setCrossSections(), and setCrossSections().

Member Data Documentation

fAbsorptionCrossSectionTable

std::vector<double> INukeOsetTable::fAbsorptionCrossSectionTable

private

pi absorption cross section

vector contains values in the following order: d0 e0, d0 e1, ... , d0 en, d1 e0 ...

Definition at line 55 of file INukeOsetTable.h.

Referenced by processLine(), and setCrossSections().

fCexCrossSectionTable

std::vector<double> INukeOsetTable::fCexCrossSectionTable[fNChannels]

private

charge-exchange piN cross section

vector contains values in the following order: d0 e0, d0 e1, ... , d0 en, d1 e0 ...
channel = 0 -> pi+n or pi-p, 1 -> pi+p or pi-n, 2 -> pi0

Definition at line 49 of file INukeOsetTable.h.

Referenced by processLine(), and setCrossSections().

fDensityBinWidth

double INukeOsetTable::fDensityBinWidth

private

density step (must be fixed)

Definition at line 59 of file INukeOsetTable.h.

Referenced by checkIntegrity(), interpolate(), and INukeOsetTable().

fDensityHandler

PointHandler INukeOsetTable::fDensityHandler

private

nuclear density handler

Definition at line 100 of file INukeOsetTable.h.

Referenced by interpolate(), INukeOsetTable(), and setupOset().

fEnergyBinWidth

double INukeOsetTable::fEnergyBinWidth

private

energy step (must be fixed)

Definition at line 60 of file INukeOsetTable.h.

Referenced by checkIntegrity(), interpolate(), and INukeOsetTable().

fEnergyHandler

PointHandler INukeOsetTable::fEnergyHandler

private

pion kinetic energy handler

Definition at line 101 of file INukeOsetTable.h.

Referenced by interpolate(), INukeOsetTable(), and setupOset().

fNDensityBins

unsigned int INukeOsetTable::fNDensityBins

private

number of denisty bins

Definition at line 57 of file INukeOsetTable.h.

Referenced by checkIntegrity(), and INukeOsetTable().

fNEnergyBins

unsigned int INukeOsetTable::fNEnergyBins

private

number of energy bins

Definition at line 58 of file INukeOsetTable.h.

Referenced by checkIntegrity(), interpolate(), and INukeOsetTable().

fQelCrossSectionTable

std::vector<double> INukeOsetTable::fQelCrossSectionTable[fNChannels]

private

quasi-elastic piN cross section

vector contains values in the following order: d0 e0, d0 e1, ... , d0 en, d1 e0 ...
channel = 0 -> pi+n or pi-p, 1 -> pi+p or pi-n, 2 -> pi0

Definition at line 42 of file INukeOsetTable.h.

Referenced by processLine(), and setCrossSections().

---

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

• INukeOsetTable.h
• INukeOsetTable.cxx