G4NeutronHPVector Class Reference

#include <G4NeutronHPVector.hh>

Public Member Functions
	G4NeutronHPVector ()
	G4NeutronHPVector (G4int n)
	~G4NeutronHPVector ()
G4NeutronHPVector &	operator= (const G4NeutronHPVector &right)
void	SetVerbose (G4int ff)
void	Times (G4double factor)
void	SetPoint (G4int i, const G4NeutronHPDataPoint &it)
void	SetData (G4int i, G4double x, G4double y)
void	SetX (G4int i, G4double e)
void	SetEnergy (G4int i, G4double e)
void	SetY (G4int i, G4double x)
void	SetXsec (G4int i, G4double x)
G4double	GetEnergy (G4int i) const
G4double	GetXsec (G4int i)
G4double	GetX (G4int i) const
const G4NeutronHPDataPoint &	GetPoint (G4int i) const
void	Hash ()
void	ReHash ()
G4double	GetXsec (G4double e)
G4double	GetXsec (G4double e, G4int min)
G4double	GetY (G4double x)
G4int	GetVectorLength () const
G4double	GetY (G4int i)
G4double	GetY (G4int i) const
void	Dump ()
void	InitInterpolation (std::istream &aDataFile)
void	Init (std::istream &aDataFile, G4int total, G4double ux=1., G4double uy=1.)
void	Init (std::istream &aDataFile, G4double ux=1., G4double uy=1.)
void	ThinOut (G4double precision)
void	SetLabel (G4double aLabel)
G4double	GetLabel ()
void	CleanUp ()
void	Merge (G4NeutronHPVector *active, G4NeutronHPVector *passive)
void	Merge (G4InterpolationScheme aScheme, G4double aValue, G4NeutronHPVector *active, G4NeutronHPVector *passive)
G4double	SampleLin ()
G4double	Sample ()
G4double *	Debug ()
void	IntegrateAndNormalise ()
void	Integrate ()
G4double	GetIntegral ()
void	SetInterpolationManager (const G4InterpolationManager &aManager)
const G4InterpolationManager &	GetInterpolationManager () const
void	SetInterpolationManager (G4InterpolationManager &aMan)
void	SetScheme (G4int aPoint, const G4InterpolationScheme &aScheme)
G4InterpolationScheme	GetScheme (G4int anIndex)
G4double	GetMeanX ()
std::vector< G4double >	GetBlocked ()
std::vector< G4double >	GetBuffered ()
G4double	Get15percentBorder ()
G4double	Get50percentBorder ()

Friends
G4NeutronHPVector &	operator+ (G4NeutronHPVector &left, G4NeutronHPVector &right)

Detailed Description

Definition at line 50 of file G4NeutronHPVector.hh.

Constructor & Destructor Documentation

G4NeutronHPVector::G4NeutronHPVector

(

)

Definition at line 80 of file G4NeutronHPVector.cc.

References DBL_MAX.

  {
    theData = new G4NeutronHPDataPoint[20]; 
    nPoints=20;
    nEntries=0;
    Verbose=0;
    theIntegral=0;
    totalIntegral=-1;
    isFreed = 0;
    maxValue = -DBL_MAX;
    the15percentBorderCash = -DBL_MAX;
    the50percentBorderCash = -DBL_MAX;
    label = -DBL_MAX;

  }

G4NeutronHPVector::G4NeutronHPVector

(

G4int

n

)

Definition at line 96 of file G4NeutronHPVector.cc.

References DBL_MAX, and G4INCL::Math::max().

  {
    nPoints=std::max(n, 20);
    theData = new G4NeutronHPDataPoint[nPoints]; 
    nEntries=0;
    Verbose=0;
    theIntegral=0;
    totalIntegral=-1;
    isFreed = 0;
    maxValue = -DBL_MAX;
    the15percentBorderCash = -DBL_MAX;
    the50percentBorderCash = -DBL_MAX;
  }

G4NeutronHPVector::~G4NeutronHPVector

(

)

Definition at line 110 of file G4NeutronHPVector.cc.

  {
//    if(Verbose==1)G4cout <<"G4NeutronHPVector::~G4NeutronHPVector"<<G4endl;
      delete [] theData;
//    if(Verbose==1)G4cout <<"Vector: delete theData"<<G4endl;
      delete [] theIntegral;
//    if(Verbose==1)G4cout <<"Vector: delete theIntegral"<<G4endl;
    isFreed = 1;
  }

Member Function Documentation

void G4NeutronHPVector::CleanUp ( )

inline

Definition at line 255 of file G4NeutronHPVector.hh.

References G4InterpolationManager::CleanUp(), G4NeutronHPHash::Clear(), and DBL_MAX.

Referenced by Merge().

  {
    nEntries=0;   
    theManager.CleanUp();
    maxValue = -DBL_MAX;
    theHash.Clear();
//080811 TK DB 
    delete[] theIntegral;
    theIntegral = NULL;
  }

G4double* G4NeutronHPVector::Debug ( )

inline

Definition at line 362 of file G4NeutronHPVector.hh.

  {
    return theIntegral;
  }

void G4NeutronHPVector::Dump

(

)

Definition at line 194 of file G4NeutronHPVector.cc.

References G4cout, G4endl, G4NeutronHPDataPoint::GetX(), and G4NeutronHPDataPoint::GetY().

  {
    G4cout << nEntries<<G4endl;
    for(G4int i=0; i<nEntries; i++)
    {
      G4cout << theData[i].GetX()<<" ";
      G4cout << theData[i].GetY()<<" ";
//      if (i!=1&&i==5*(i/5)) G4cout << G4endl;
      G4cout << G4endl;
    }
    G4cout << G4endl;
  }

G4double G4NeutronHPVector::Get15percentBorder

(

)

Definition at line 441 of file G4NeutronHPVector.cc.

References DBL_MAX, GetVectorLength(), G4NeutronHPDataPoint::GetX(), IntegrateAndNormalise(), and G4INCL::Math::min().

  {    
    if(the15percentBorderCash>-DBL_MAX/2.) return the15percentBorderCash;
    G4double result;
    if(GetVectorLength()==1)
    {
      result = theData[0].GetX();
      the15percentBorderCash = result;
    }
    else
    {
      if(theIntegral==0) { IntegrateAndNormalise(); }
      G4int i;
      result = theData[GetVectorLength()-1].GetX();
      for(i=0;i<GetVectorLength();i++)
      {
    if(theIntegral[i]/theIntegral[GetVectorLength()-1]>0.15) 
    {
      result = theData[std::min(i+1, GetVectorLength()-1)].GetX();
          the15percentBorderCash = result;
      break;
    }
      }
      the15percentBorderCash = result;
    }
    return result;
  }

G4double G4NeutronHPVector::Get50percentBorder

(

)

Definition at line 469 of file G4NeutronHPVector.cc.

References DBL_MAX, GetVectorLength(), G4NeutronHPDataPoint::GetX(), IntegrateAndNormalise(), G4NeutronHPInterpolator::Lin(), and test::x.

  {    
    if(the50percentBorderCash>-DBL_MAX/2.) return the50percentBorderCash;
    G4double result;
    if(GetVectorLength()==1)
    {
      result = theData[0].GetX();
      the50percentBorderCash = result;
    }
    else
    {
      if(theIntegral==0) { IntegrateAndNormalise(); }
      G4int i;
      G4double x = 0.5;
      result = theData[GetVectorLength()-1].GetX();
      for(i=0;i<GetVectorLength();i++)
      {
    if(theIntegral[i]/theIntegral[GetVectorLength()-1]>x) 
    {
      G4int it;
      it = i;
      if(it == GetVectorLength()-1)
      {
        result = theData[GetVectorLength()-1].GetX();
      }
      else
      {
        G4double x1, x2, y1, y2;
        x1 = theIntegral[i-1]/theIntegral[GetVectorLength()-1];
        x2 = theIntegral[i]/theIntegral[GetVectorLength()-1];
        y1 = theData[i-1].GetX();
        y2 = theData[i].GetX();
        result = theLin.Lin(x, x1, x2, y1, y2);
      }
          the50percentBorderCash = result;
      break;
    }
      }
      the50percentBorderCash = result;
    }
    return result;
  }

std::vector<G4double> G4NeutronHPVector::GetBlocked ( )

inline

Definition at line 532 of file G4NeutronHPVector.hh.

{return theBlocked;}

std::vector<G4double> G4NeutronHPVector::GetBuffered ( )

inline

Definition at line 533 of file G4NeutronHPVector.hh.

{return theBuffered;}

G4double G4NeutronHPVector::GetEnergy ( G4int  i ) const

inline

Definition at line 119 of file G4NeutronHPVector.hh.

References G4NeutronHPDataPoint::GetX().

Referenced by G4NeutronHPElementData::Harmonise(), G4NeutronHPChannel::Harmonise(), G4NeutronHPArbitaryTab::Init(), and Merge().

{ return theData[i].GetX(); }

G4double G4NeutronHPVector::GetIntegral ( )

inline

Definition at line 471 of file G4NeutronHPVector.hh.

References Integrate().

Referenced by G4NeutronHPLabAngularEnergy::Sample().

  { 
    if(totalIntegral<-0.5) Integrate();
    return totalIntegral; 
  }

const G4InterpolationManager& G4NeutronHPVector::GetInterpolationManager ( ) const

inline

Definition at line 482 of file G4NeutronHPVector.hh.

Referenced by G4NeutronHPLabAngularEnergy::Sample().

  {
    return theManager;
  }

G4double G4NeutronHPVector::GetLabel ( )

inline

Definition at line 250 of file G4NeutronHPVector.hh.

Referenced by Merge(), G4NeutronHPLabAngularEnergy::Sample(), and G4NeutronHPArbitaryTab::Sample().

  {
    return label;
  }

G4double G4NeutronHPVector::GetMeanX ( )

inline

Definition at line 502 of file G4NeutronHPVector.hh.

References G4NeutronHPInterpolator::GetBinIntegral(), G4InterpolationManager::GetScheme(), G4NeutronHPInterpolator::GetWeightedBinIntegral(), G4NeutronHPDataPoint::GetX(), and G4NeutronHPDataPoint::GetY().

Referenced by G4NeutronHPLabAngularEnergy::Sample().

  {
    G4double result;
    G4double running = 0;
    G4double weighted = 0;
    for(G4int i=1; i<nEntries; i++)
    {
      running += theInt.GetBinIntegral(theManager.GetScheme(i-1),
                           theData[i-1].GetX(), theData[i].GetX(),
                           theData[i-1].GetY(), theData[i].GetY());
      weighted += theInt.GetWeightedBinIntegral(theManager.GetScheme(i-1),
                           theData[i-1].GetX(), theData[i].GetX(),
                           theData[i-1].GetY(), theData[i].GetY());
    }  
    result = weighted / running;  
    return result;
  }

const G4NeutronHPDataPoint& G4NeutronHPVector::GetPoint ( G4int  i ) const

inline

Definition at line 127 of file G4NeutronHPVector.hh.

Referenced by G4NeutronHPIsoData::FillChannelData(), and operator=().

{ return theData[i]; }

G4InterpolationScheme G4NeutronHPVector::GetScheme ( G4int  anIndex )

inline

Definition at line 497 of file G4NeutronHPVector.hh.

References G4InterpolationManager::GetScheme().

Referenced by Merge().

  {
    return theManager.GetScheme(anIndex);
  }

G4int G4NeutronHPVector::GetVectorLength ( ) const

inline

Definition at line 193 of file G4NeutronHPVector.hh.

Referenced by G4NeutronHPData::DoPhysicsVector(), G4NeutronHPIsoData::FillChannelData(), Get15percentBorder(), Get50percentBorder(), G4NeutronHPPhotonDist::GetPhotons(), GetX(), GetY(), G4NeutronHPElementData::Harmonise(), G4NeutronHPChannel::Harmonise(), G4NeutronHPArbitaryTab::Init(), Integrate(), IntegrateAndNormalise(), Merge(), operator+(), G4NeutronHPContAngularPar::Sample(), G4NeutronHPLabAngularEnergy::Sample(), G4NeutronHPDiscreteTwoBody::Sample(), Sample(), SampleLin(), and ThinOut().

{return nEntries;}

G4double G4NeutronHPVector::GetX ( G4int  i ) const

inline

Definition at line 121 of file G4NeutronHPVector.hh.

References GetVectorLength(), and G4NeutronHPDataPoint::GetX().

Referenced by G4NeutronHPData::DoPhysicsVector(), G4NeutronHPPhotonDist::GetPhotons(), GetXsec(), Hash(), Integrate(), operator+(), G4NeutronHPContAngularPar::Sample(), G4NeutronHPLabAngularEnergy::Sample(), G4NeutronHPLegendreStore::Sample(), G4NeutronHPDiscreteTwoBody::Sample(), and Sample().

  { 
    if (i<0) i=0;
    if(i>=GetVectorLength()) i=GetVectorLength()-1;
    return theData[i].GetX();
  }

G4double G4NeutronHPVector::GetXsec ( G4int  i )

inline

Definition at line 120 of file G4NeutronHPVector.hh.

References G4NeutronHPDataPoint::GetY().

Referenced by G4NeutronHPPhotonDist::GetPhotons(), G4NeutronHPIsoData::GetXsec(), G4NeutronHPChannel::GetXsec(), GetY(), G4NeutronHPElementData::Harmonise(), G4NeutronHPChannel::Harmonise(), and Merge().

{ return theData[i].GetY(); }

G4double G4NeutronHPVector::GetXsec

(

G4double

e

)

Definition at line 145 of file G4NeutronHPVector.cc.

References G4NeutronHPHash::GetMinIndex(), G4InterpolationManager::GetScheme(), G4NeutronHPDataPoint::GetX(), GetX(), G4NeutronHPDataPoint::GetY(), Hash(), G4NeutronHPInterpolator::Interpolate(), G4INCL::Math::min(), and G4NeutronHPHash::Prepared().

  {
    if(nEntries == 0) return 0;
    if(!theHash.Prepared()) Hash();
    G4int min = theHash.GetMinIndex(e);
    G4int i;
    for(i=min ; i<nEntries; i++)
    {
      //if(theData[i].GetX()>e) break;
      if(theData[i].GetX() >= e) break;
    }
    G4int low = i-1;
    G4int high = i;
    if(i==0)
    {
      low = 0;
      high = 1;
    }
    else if(i==nEntries)
    {
      low = nEntries-2;
      high = nEntries-1;
    }
    G4double y;
    if(e<theData[nEntries-1].GetX()) 
    {
      // Protect against doubled-up x values
      //if( (theData[high].GetX()-theData[low].GetX())/theData[high].GetX() < 0.000001)
      if ( theData[high].GetX() !=0 
       //080808 TKDB
       //&&( theData[high].GetX()-theData[low].GetX())/theData[high].GetX() < 0.000001)
       &&( std::abs( (theData[high].GetX()-theData[low].GetX())/theData[high].GetX() ) < 0.000001 ) )
      {
        y = theData[low].GetY();
      }
      else
      {
        y = theInt.Interpolate(theManager.GetScheme(high), e, 
                               theData[low].GetX(), theData[high].GetX(),
                   theData[low].GetY(), theData[high].GetY());
      }
    }
    else
    {
      y=theData[nEntries-1].GetY();
    }
    return y;
  }

G4double G4NeutronHPVector::GetXsec ( G4double  e, G4int  min  )

inline

Definition at line 151 of file G4NeutronHPVector.hh.

References G4InterpolationManager::GetScheme(), G4NeutronHPDataPoint::GetX(), GetX(), G4NeutronHPDataPoint::GetY(), and G4NeutronHPInterpolator::Interpolate().

  {
    G4int i;
    for(i=min ; i<nEntries; i++)
    {
      if(theData[i].GetX()>e) break;
    }
    G4int low = i-1;
    G4int high = i;
    if(i==0)
    {
      low = 0;
      high = 1;
    }
    else if(i==nEntries)
    {
      low = nEntries-2;
      high = nEntries-1;
    }
    G4double y;
    if(e<theData[nEntries-1].GetX()) 
    {
      // Protect against doubled-up x values
      if( (theData[high].GetX()-theData[low].GetX())/theData[high].GetX() < 0.000001)
      {
        y = theData[low].GetY();
      }
      else
      {
        y = theInt.Interpolate(theManager.GetScheme(high), e, 
                               theData[low].GetX(), theData[high].GetX(),
                   theData[low].GetY(), theData[high].GetY());
      }
    }
    else
    {
      y=theData[nEntries-1].GetY();
    }
    return y;
  }

G4double G4NeutronHPVector::GetY ( G4double  x )

inline

Definition at line 192 of file G4NeutronHPVector.hh.

References GetXsec().

Referenced by G4NeutronHPData::DoPhysicsVector(), G4NeutronHPNeutronYield::GetDelayed(), G4NeutronHPEvapSpectrum::GetFractionalProbability(), G4NeutronHPFissionSpectrum::GetFractionalProbability(), G4NeutronHPSimpleEvapSpectrum::GetFractionalProbability(), G4NeutronHPWattSpectrum::GetFractionalProbability(), G4NeutronHPMadlandNixSpectrum::GetFractionalProbability(), G4NeutronHPArbitaryTab::GetFractionalProbability(), G4NeutronHPNeutronYield::GetMean(), G4NeutronHPProduct::GetMeanYield(), G4NeutronHPPhotonDist::GetPhotons(), G4NeutronHPNeutronYield::GetPrompt(), G4NeutronHPFissionBaseFS::GetXsec(), G4NeutronHPInelasticBaseFS::GetXsec(), Hash(), operator+(), G4NeutronHPContAngularPar::Sample(), G4NeutronHPFissionSpectrum::Sample(), G4NeutronHPEvapSpectrum::Sample(), G4NeutronHPSimpleEvapSpectrum::Sample(), G4NeutronHPWattSpectrum::Sample(), G4NeutronHPLabAngularEnergy::Sample(), G4NeutronHPMadlandNixSpectrum::Sample(), G4NeutronHPLegendreStore::Sample(), G4NeutronHPDiscreteTwoBody::Sample(), G4NeutronHPProduct::Sample(), Sample(), ThinOut(), and Times().

{return GetXsec(x);}

G4double G4NeutronHPVector::GetY ( G4int  i )

inline

Definition at line 195 of file G4NeutronHPVector.hh.

References GetVectorLength(), and G4NeutronHPDataPoint::GetY().

  { 
    if (i<0) i=0;
    if(i>=GetVectorLength()) i=GetVectorLength()-1;
    return theData[i].GetY(); 
  }

G4double G4NeutronHPVector::GetY ( G4int  i ) const

inline

Definition at line 202 of file G4NeutronHPVector.hh.

References GetVectorLength(), and G4NeutronHPDataPoint::GetY().

  {
    if (i<0) i=0;
    if(i>=GetVectorLength()) i=GetVectorLength()-1;
    return theData[i].GetY(); 
  }

void G4NeutronHPVector::Hash ( )

inline

Definition at line 129 of file G4NeutronHPVector.hh.

References GetX(), GetY(), G4NeutronHPHash::SetData(), and test::x.

Referenced by GetXsec(), and ReHash().

  {
    G4int i;
    G4double x, y;
    for(i=0 ; i<nEntries; i++)
    {
      if(0 == (i+1)%10)
      {
        x = GetX(i);
    y = GetY(i);
    theHash.SetData(i, x, y);
      }
    }
  }

void G4NeutronHPVector::Init ( std::istream &  aDataFile, G4int  total, G4double  ux = 1., G4double  uy = 1.  )

inline

Definition at line 215 of file G4NeutronHPVector.hh.

References SetData(), G4NeutronHPHash::SetData(), G4INCL::CrossSections::total(), and test::x.

Referenced by G4NeutronHPEvapSpectrum::Init(), G4NeutronHPFissionSpectrum::Init(), G4NeutronHPSimpleEvapSpectrum::Init(), G4NeutronHPWattSpectrum::Init(), G4NeutronHPFissionBaseFS::Init(), G4NeutronHPArbitaryTab::Init(), G4NeutronHPMadlandNixSpectrum::Init(), G4NeutronHPPhotonXSection::Init(), G4NeutronHPInelasticBaseFS::Init(), G4NeutronHPIsoData::Init(), G4NeutronHPProduct::Init(), G4NeutronHPLabAngularEnergy::Init(), G4NeutronHPInelasticCompFS::Init(), Init(), G4NeutronHPNeutronYield::InitDelayed(), G4NeutronHPPhotonDist::InitEnergies(), G4NeutronHPNeutronYield::InitMean(), G4NeutronHPPhotonDist::InitMean(), G4NeutronHPPhotonDist::InitPartials(), and G4NeutronHPNeutronYield::InitPrompt().

  {
    G4double x,y;
    for (G4int i=0;i<total;i++)
    {
      aDataFile >> x >> y;
      x*=ux;
      y*=uy;
      SetData(i,x,y);
      if(0 == nEntries%10)
      {
        theHash.SetData(nEntries-1, x, y);
      }
    }
  }

void G4NeutronHPVector::Init ( std::istream &  aDataFile, G4double  ux = 1., G4double  uy = 1.  )

inline

Definition at line 231 of file G4NeutronHPVector.hh.

References G4InterpolationManager::Init(), Init(), and G4INCL::CrossSections::total().

  {
    G4int total;
    aDataFile >> total;
    if(theData!=0) delete [] theData;
    theData = new G4NeutronHPDataPoint[total]; 
    nPoints=total;
    nEntries=0;    
    theManager.Init(aDataFile);
    Init(aDataFile, total, ux, uy);
  }

void G4NeutronHPVector::InitInterpolation ( std::istream &  aDataFile )

inline

Definition at line 210 of file G4NeutronHPVector.hh.

References G4InterpolationManager::Init().

  {
    theManager.Init(aDataFile);
  }

void G4NeutronHPVector::Integrate ( )

inline

Definition at line 417 of file G4NeutronHPVector.hh.

References test::a, test::b, CHISTO, CLINLIN, CLINLOG, CLOGLIN, CLOGLOG, G4InterpolationManager::GetScheme(), GetVectorLength(), G4NeutronHPDataPoint::GetX(), GetX(), G4NeutronHPDataPoint::GetY(), HISTO, LINLIN, LINLOG, LOGLIN, LOGLOG, UHISTO, ULINLIN, ULINLOG, ULOGLIN, and ULOGLOG.

Referenced by GetIntegral().

  {
    G4int i;
    if(nEntries == 1)
    {
      totalIntegral = 0;
      return;
    }
    G4double sum = 0;
    for(i=1;i<GetVectorLength();i++)
    {
      if(std::abs((theData[i].GetX()-theData[i-1].GetX())/theData[i].GetX())>0.0000001)
      {
        G4double x1 = theData[i-1].GetX();
        G4double x2 = theData[i].GetX();
        G4double y1 = theData[i-1].GetY();
        G4double y2 = theData[i].GetY();
        G4InterpolationScheme aScheme = theManager.GetScheme(i);
        if(aScheme==LINLIN||aScheme==CLINLIN||aScheme==ULINLIN)
        {
          sum+= 0.5*(y2+y1)*(x2-x1);
        }
        else if(aScheme==LINLOG||aScheme==CLINLOG||aScheme==ULINLOG)
        {
          G4double a = y1;
          G4double b = (y2-y1)/(std::log(x2)-std::log(x1));
          sum+= (a-b)*(x2-x1) + b*(x2*std::log(x2)-x1*std::log(x1));
        }
        else if(aScheme==LOGLIN||aScheme==CLOGLIN||aScheme==ULOGLIN)
        {
          G4double a = std::log(y1);
          G4double b = (std::log(y2)-std::log(y1))/(x2-x1);
          sum += (std::exp(a)/b)*(std::exp(b*x2)-std::exp(b*x1));
        }
        else if(aScheme==HISTO||aScheme==CHISTO||aScheme==UHISTO)
        {
          sum+= y1*(x2-x1);
        }
        else if(aScheme==LOGLOG||aScheme==CLOGLOG||aScheme==ULOGLOG)
        {
          G4double a = std::log(y1);
          G4double b = (std::log(y2)-std::log(y1))/(std::log(x2)-std::log(x1));
          sum += (std::exp(a)/(b+1))*(std::pow(x2,b+1)-std::pow(x1,b+1));
        }
        else
        {
          throw G4HadronicException(__FILE__, __LINE__, "Unknown interpolation scheme in G4NeutronHPVector::Integrate");
        }
          
      }
    }
    totalIntegral = sum;
  }

void G4NeutronHPVector::IntegrateAndNormalise ( )

inline

Definition at line 367 of file G4NeutronHPVector.hh.

References G4NeutronHPInterpolator::GetBinIntegral(), G4InterpolationManager::GetScheme(), GetVectorLength(), G4NeutronHPDataPoint::GetX(), G4NeutronHPDataPoint::GetY(), and G4INCL::CrossSections::total().

Referenced by Get15percentBorder(), Get50percentBorder(), Sample(), and SampleLin().

  {
    G4int i;
    if(theIntegral!=0) return;
    theIntegral = new G4double[nEntries];
    if(nEntries == 1)
    {
      theIntegral[0] = 1;
      return;
    }
    theIntegral[0] = 0;
    G4double sum = 0;
    G4double x1 = 0;
    G4double x0 = 0;
    for(i=1;i<GetVectorLength();i++)
    {
      x1 = theData[i].GetX();
      x0 = theData[i-1].GetX();
      if (std::abs(x1-x0) > std::abs(x1*0.0000001) )
      {
    //********************************************************************
    //EMendoza -> the interpolation scheme is not always lin-lin
    /*
        sum+= 0.5*(theData[i].GetY()+theData[i-1].GetY())*
                  (x1-x0);
    */
    //********************************************************************
        G4InterpolationScheme aScheme = theManager.GetScheme(i);
        G4double y0 = theData[i-1].GetY();
        G4double y1 = theData[i].GetY();
    G4double integ=theInt.GetBinIntegral(aScheme,x0,x1,y0,y1);
#if defined WIN32-VC
    if(!_finite(integ)){integ=0;}
#elif defined __IBMCPP__
    if(isinf(integ)||isnan(integ)){integ=0;}
#else
    if(std::isinf(integ)||std::isnan(integ)){integ=0;}
#endif
    sum+=integ;
    //********************************************************************
      }
      theIntegral[i] = sum;
    }
    G4double total = theIntegral[GetVectorLength()-1];
    for(i=1;i<GetVectorLength();i++)
    {
      theIntegral[i]/=total;
    }
  }

void G4NeutronHPVector::Merge ( G4NeutronHPVector *  active, G4NeutronHPVector *  passive  )

inline

Definition at line 267 of file G4NeutronHPVector.hh.

References test::a, active, G4InterpolationManager::AppendScheme(), CleanUp(), GetEnergy(), GetScheme(), GetVectorLength(), GetXsec(), n, and SetData().

Referenced by G4NeutronHPLabAngularEnergy::Sample(), and G4NeutronHPDiscreteTwoBody::Sample().

  {
    CleanUp();
    G4int s_tmp = 0, n=0, m_tmp=0;
    G4NeutronHPVector * tmp;
    G4int a = s_tmp, p = n, t;
    while (a<active->GetVectorLength()&&p<passive->GetVectorLength())
    {
      if(active->GetEnergy(a) <= passive->GetEnergy(p))
      {
        G4double xa = active->GetEnergy(a);
        G4double yy = active->GetXsec(a);
        SetData(m_tmp, xa, yy);
        theManager.AppendScheme(m_tmp, active->GetScheme(a));
        m_tmp++;
        a++;
        G4double xp = passive->GetEnergy(p);

//080409 TKDB 
        //if( std::abs(std::abs(xp-xa)/xa)<0.001 ) p++;
        if ( !( xa == 0 ) && std::abs(std::abs(xp-xa)/xa)<0.001 ) p++;
      } else {
        tmp = active; 
        t=a;
        active = passive; 
        a=p;
        passive = tmp; 
        p=t;
      }
    }
    while (a!=active->GetVectorLength())
    {
      SetData(m_tmp, active->GetEnergy(a), active->GetXsec(a));
      theManager.AppendScheme(m_tmp++, active->GetScheme(a));
      a++;
    }
    while (p!=passive->GetVectorLength())
    {
      if(std::abs(GetEnergy(m_tmp-1)-passive->GetEnergy(p))/passive->GetEnergy(p)>0.001)
      //if(std::abs(GetEnergy(m)-passive->GetEnergy(p))/passive->GetEnergy(p)>0.001)
      {
        SetData(m_tmp, passive->GetEnergy(p), passive->GetXsec(p));
        theManager.AppendScheme(m_tmp++, active->GetScheme(p));
      }
      p++;
    }
  }    

void G4NeutronHPVector::Merge	(	G4InterpolationScheme	aScheme,
		G4double	aValue,
		G4NeutronHPVector *	active,
		G4NeutronHPVector *	passive
	)

Definition at line 222 of file G4NeutronHPVector.cc.

References test::a, active, G4InterpolationManager::AppendScheme(), CleanUp(), GetEnergy(), GetLabel(), GetScheme(), GetVectorLength(), GetXsec(), G4NeutronHPInterpolator::Interpolate(), n, G4NeutronHPHash::Prepared(), ReHash(), and SetData().

  { 
    // interpolate between labels according to aScheme, cut at aValue, 
    // continue in unknown areas by substraction of the last difference.
    
    CleanUp();
    G4int s_tmp = 0, n=0, m_tmp=0;
    G4NeutronHPVector * tmp;
    G4int a = s_tmp, p = n, t;
    while ( a<active->GetVectorLength() )
    {
      if(active->GetEnergy(a) <= passive->GetEnergy(p))
      {
        G4double xa  = active->GetEnergy(a);
        G4double yy = theInt.Interpolate(aScheme, aValue, active->GetLabel(), passive->GetLabel(),
                                                          active->GetXsec(a), passive->GetXsec(xa));
        SetData(m_tmp, xa, yy);
        theManager.AppendScheme(m_tmp, active->GetScheme(a));
        m_tmp++;
        a++;
        G4double xp = passive->GetEnergy(p);
        //if( std::abs(std::abs(xp-xa)/xa)<0.0000001&&a<active->GetVectorLength() ) 
        if ( xa != 0 
          && std::abs(std::abs(xp-xa)/xa) < 0.0000001 
          && a < active->GetVectorLength() )
        {
          p++;
          tmp = active; t=a;
          active = passive; a=p;
          passive = tmp; p=t;
        }
      } else {
        tmp = active; t=a;
        active = passive; a=p;
        passive = tmp; p=t;
      }
    }
    
    G4double deltaX = passive->GetXsec(GetEnergy(m_tmp-1)) - GetXsec(m_tmp-1);
    while (p!=passive->GetVectorLength()&&passive->GetEnergy(p)<=aValue)
    {
      G4double anX;
      anX = passive->GetXsec(p)-deltaX;
      if(anX>0)
      {
        //if(std::abs(GetEnergy(m-1)-passive->GetEnergy(p))/passive->GetEnergy(p)>0.0000001)
        if ( passive->GetEnergy(p) == 0 
          || std::abs(GetEnergy(m_tmp-1)-passive->GetEnergy(p))/passive->GetEnergy(p) > 0.0000001 )
        {
          SetData(m_tmp, passive->GetEnergy(p), anX);
          theManager.AppendScheme(m_tmp++, passive->GetScheme(p));
        }
      }
      p++;
    }
    // Rebuild the Hash;
    if(theHash.Prepared()) 
    {
      ReHash();
    }
  }

G4NeutronHPVector & G4NeutronHPVector::operator=

(

const G4NeutronHPVector &

right

)

Definition at line 121 of file G4NeutronHPVector.cc.

References GetPoint(), and SetPoint().

  {
    if(&right == this) return *this;
    
    G4int i;
    
    totalIntegral = right.totalIntegral;
    if(right.theIntegral!=0) theIntegral = new G4double[right.nEntries];
    for(i=0; i<right.nEntries; i++)
    {
      SetPoint(i, right.GetPoint(i)); // copy theData
      if(right.theIntegral!=0) theIntegral[i] = right.theIntegral[i];
    }
    theManager = right.theManager; 
    label = right.label;
  
    Verbose = right.Verbose;
    the15percentBorderCash = right.the15percentBorderCash;
    the50percentBorderCash = right.the50percentBorderCash;
    theHash = right.theHash;
   return *this;
  }

void G4NeutronHPVector::ReHash ( )

inline

Definition at line 144 of file G4NeutronHPVector.hh.

References G4NeutronHPHash::Clear(), and Hash().

Referenced by Merge(), and ThinOut().

  {
    theHash.Clear();
    Hash();
  }

G4double G4NeutronHPVector::Sample

(

)

Definition at line 348 of file G4NeutronHPVector.cc.

References G4cout, G4endl, G4UniformRand, GetVectorLength(), G4NeutronHPDataPoint::GetX(), GetX(), G4NeutronHPDataPoint::GetY(), GetY(), IntegrateAndNormalise(), G4INCL::Math::max(), SetY(), and mcscore::test().

Referenced by G4NeutronHPPhotonDist::GetPhotons(), G4NeutronHPContAngularPar::Sample(), G4NeutronHPEvapSpectrum::Sample(), G4NeutronHPLabAngularEnergy::Sample(), G4NeutronHPDiscreteTwoBody::Sample(), G4NeutronHPArbitaryTab::Sample(), and G4NeutronHPPartial::Sample().

  {
    G4double result;
    G4int j;
    for(j=0; j<GetVectorLength(); j++)
    {
      if(GetY(j)<0) SetY(j, 0);
    }
    
    if(theBuffered.size() !=0 && G4UniformRand()<0.5) 
    {
      result = theBuffered[0];
      theBuffered.erase(theBuffered.begin());
      if(result < GetX(GetVectorLength()-1) ) return result;
    }
    if(GetVectorLength()==1)
    {
      result = theData[0].GetX();
    }
    else
    {
      if(theIntegral==0) { IntegrateAndNormalise(); }
      do
      {
//080808
/*
        G4double rand;
        G4double value, test, baseline;
        baseline = theData[GetVectorLength()-1].GetX()-theData[0].GetX();
        do
        {
          value = baseline*G4UniformRand();
          value += theData[0].GetX();
          test = GetY(value)/maxValue;
          rand = G4UniformRand();
        }
        //while(test<rand);
        while( test < rand && test > 0 );
        result = value;
*/
        G4double rand;
        G4double value, test;
        do 
        {
           rand = G4UniformRand();
           G4int ibin = -1;
           for ( G4int i = 0 ; i < GetVectorLength() ; i++ )
           {
              if ( rand < theIntegral[i] ) 
              {
                 ibin = i; 
                 break;
              }
           }
           if ( ibin < 0 ) G4cout << "TKDB 080807 " << rand << G4endl; 
           // result 
           rand = G4UniformRand();
           G4double x1, x2; 
           if ( ibin == 0 ) 
           {
              x1 = theData[ ibin ].GetX(); 
              value = x1; 
              break;
           }
           else 
           {
              x1 = theData[ ibin-1 ].GetX();
           }
           
           x2 = theData[ ibin ].GetX();
           value = rand * ( x2 - x1 ) + x1;
       //***********************************************************************
       /*
           test = GetY ( value ) / std::max ( GetY( ibin-1 ) , GetY ( ibin ) ); 
       */
       //***********************************************************************
       //EMendoza - Always linear interpolation:
       G4double y1=theData[ ibin-1 ].GetY();
           G4double y2=theData[ ibin ].GetY();
       G4double mval=(y2-y1)/(x2-x1);
       G4double bval=y1-mval*x1;
       test =(mval*value+bval)/std::max ( GetY( ibin-1 ) , GetY ( ibin ) ); 
       //***********************************************************************
        }
        while ( G4UniformRand() > test );
        result = value;
//080807
      }
      while(IsBlocked(result));
    }
    return result;
  }

G4double G4NeutronHPVector::SampleLin ( )

inline

Definition at line 318 of file G4NeutronHPVector.hh.

References G4UniformRand, GetVectorLength(), G4NeutronHPDataPoint::GetX(), IntegrateAndNormalise(), and G4NeutronHPInterpolator::Lin().

  {
    G4double result;
    if(theIntegral==0) IntegrateAndNormalise();
    if(GetVectorLength()==1)
    {
      result = theData[0].GetX();
    }
    else
    {
      G4int i;
      G4double rand = G4UniformRand();
      
      // this was replaced 
//      for(i=1;i<GetVectorLength();i++)
//      {
//  if(rand<theIntegral[i]/theIntegral[GetVectorLength()-1]) break;
//      }

// by this (begin)
      for(i=GetVectorLength()-1; i>=0 ;i--)
      {
    if(rand>theIntegral[i]/theIntegral[GetVectorLength()-1]) break;
      }
      if(i!=GetVectorLength()-1) i++;
// until this (end)
      
      G4double x1, x2, y1, y2;
      y1 = theData[i-1].GetX();
      x1 = theIntegral[i-1];
      y2 = theData[i].GetX();
      x2 = theIntegral[i];
      if(std::abs((y2-y1)/y2)<0.0000001) // not really necessary, since the case is excluded by construction
      {
    y1 = theData[i-2].GetX();
    x1 = theIntegral[i-2];
      }
      result = theLin.Lin(rand, x1, x2, y1, y2);
    }
    return result;
  }

void G4NeutronHPVector::SetData ( G4int  i, G4double  x, G4double  y  )

inline

Definition at line 90 of file G4NeutronHPVector.hh.

References G4NeutronHPDataPoint::SetData().

Referenced by G4NeutronHPPartial::GetY(), G4NeutronHPElementData::Harmonise(), G4NeutronHPChannel::Harmonise(), Init(), Merge(), operator+(), G4NeutronHPLabAngularEnergy::Sample(), G4NeutronHPLegendreStore::Sample(), G4NeutronHPDiscreteTwoBody::Sample(), and SetPoint().

  { 
//    G4cout <<"G4NeutronHPVector::SetData called"<<nPoints<<" "<<nEntries<<G4endl;
    Check(i);
    if(y>maxValue) maxValue=y;
    theData[i].SetData(x, y);
  }

void G4NeutronHPVector::SetEnergy ( G4int  i, G4double  e  )

inline

Definition at line 102 of file G4NeutronHPVector.hh.

References G4NeutronHPDataPoint::SetX().

  {
    Check(i);
    theData[i].SetX(e);
  }

void G4NeutronHPVector::SetInterpolationManager ( const G4InterpolationManager &  aManager )

inline

Definition at line 477 of file G4NeutronHPVector.hh.

Referenced by G4NeutronHPPartial::GetY(), G4NeutronHPContAngularPar::Sample(), G4NeutronHPLabAngularEnergy::Sample(), G4NeutronHPDiscreteTwoBody::Sample(), and G4NeutronHPPartial::Sample().

  {
    theManager = aManager;
  }

void G4NeutronHPVector::SetInterpolationManager ( G4InterpolationManager &  aMan )

inline

Definition at line 487 of file G4NeutronHPVector.hh.

  {
    theManager = aMan;
  }

void G4NeutronHPVector::SetLabel ( G4double  aLabel )

inline

Definition at line 245 of file G4NeutronHPVector.hh.

Referenced by Lesson1Wx.MyText::ChangeColor(), Lesson2Wx.MyText::ChangeColor(), G4NeutronHPArbitaryTab::Init(), and G4NeutronHPLabAngularEnergy::Init().

  {
    label = aLabel;
  }

void G4NeutronHPVector::SetPoint ( G4int  i, const G4NeutronHPDataPoint &  it  )

inline

Definition at line 83 of file G4NeutronHPVector.hh.

References G4NeutronHPDataPoint::GetX(), G4NeutronHPDataPoint::GetY(), SetData(), and test::x.

Referenced by G4NeutronHPIsoData::FillChannelData(), and operator=().

  {
    G4double x = it.GetX();
    G4double y = it.GetY();
    SetData(i, x, y);
  }

void G4NeutronHPVector::SetScheme ( G4int  aPoint, const G4InterpolationScheme &  aScheme  )

inline

Definition at line 492 of file G4NeutronHPVector.hh.

References G4InterpolationManager::AppendScheme().

Referenced by G4NeutronHPPartial::GetY(), and G4NeutronHPPartial::Sample().

  {
    theManager.AppendScheme(aPoint, aScheme);
  }

void G4NeutronHPVector::SetVerbose ( G4int  ff )

inline

Definition at line 65 of file G4NeutronHPVector.hh.

  {
    Verbose = ff;
  }

void G4NeutronHPVector::SetX ( G4int  i, G4double  e  )

inline

Definition at line 97 of file G4NeutronHPVector.hh.

References G4NeutronHPDataPoint::SetX().

Referenced by G4NeutronHPContAngularPar::Sample(), G4NeutronHPLabAngularEnergy::Sample(), G4NeutronHPDiscreteTwoBody::Sample(), and G4NeutronHPPartial::Sample().

  {
    Check(i);
    theData[i].SetX(e);
  }

void G4NeutronHPVector::SetXsec ( G4int  i, G4double  x  )

inline

Definition at line 113 of file G4NeutronHPVector.hh.

References G4NeutronHPDataPoint::SetY(), and test::x.

  {
    Check(i);
    if(x>maxValue) maxValue=x;
    theData[i].SetY(x);
  }

void G4NeutronHPVector::SetY ( G4int  i, G4double  x  )

inline

Definition at line 107 of file G4NeutronHPVector.hh.

References G4NeutronHPDataPoint::SetY(), and test::x.

Referenced by G4NeutronHPContAngularPar::Sample(), G4NeutronHPLabAngularEnergy::Sample(), G4NeutronHPDiscreteTwoBody::Sample(), G4NeutronHPPartial::Sample(), and Sample().

  {
    Check(i);
    if(x>maxValue) maxValue=x;
    theData[i].SetY(x);
  }

void G4NeutronHPVector::ThinOut

(

G4double

precision

)

Definition at line 285 of file G4NeutronHPVector.cc.

References GetVectorLength(), G4NeutronHPDataPoint::GetX(), G4NeutronHPDataPoint::GetY(), GetY(), G4NeutronHPInterpolator::Lin(), G4NeutronHPHash::Prepared(), ReHash(), and test::x.

Referenced by G4NeutronHPElementData::Init(), operator+(), and G4NeutronHPIsoData::ThinOut().

  {
    // anything in there?
    if(GetVectorLength()==0) return;
    // make the new vector
    G4NeutronHPDataPoint * aBuff = new G4NeutronHPDataPoint[nPoints];
    G4double x, x1, x2, y, y1, y2;
    G4int count = 0, current = 2, start = 1;
    
    // First element always goes and is never tested.
    aBuff[0] = theData[0];
    
    // Find the rest
    while(current < GetVectorLength())
    {
      x1=aBuff[count].GetX();
      y1=aBuff[count].GetY();
      x2=theData[current].GetX();
      y2=theData[current].GetY();
      for(G4int j=start; j<current; j++)
      {
    x = theData[j].GetX();
    if(x1-x2 == 0) y = (y2+y1)/2.;
    else y = theInt.Lin(x, x1, x2, y1, y2);
    if (std::abs(y-theData[j].GetY())>precision*y)
    {
      aBuff[++count] = theData[current-1]; // for this one, everything was fine
          start = current; // the next candidate
      break;
    }
      }
      current++ ;
    }
    // The last one also always goes, and is never tested.
    aBuff[++count] = theData[GetVectorLength()-1];
    delete [] theData;
    theData = aBuff;
    nEntries = count+1;
    
    // Rebuild the Hash;
    if(theHash.Prepared()) 
    {
      ReHash();
    }
  }

void G4NeutronHPVector::Times ( G4double  factor )

inline

Definition at line 70 of file G4NeutronHPVector.hh.

References GetY(), and G4NeutronHPDataPoint::SetY().

Referenced by G4NeutronHPChannel::UpdateData().

  {
    G4int i;
    for(i=0; i<nEntries; i++)
    {
      theData[i].SetY(theData[i].GetY()*factor);
    }
    if(theIntegral!=0)
    {
      theIntegral[i] *= factor;
    }
  }

Friends And Related Function Documentation

G4NeutronHPVector& operator+ ( G4NeutronHPVector &  left, G4NeutronHPVector &  right  )

friend

Definition at line 37 of file G4NeutronHPVector.cc.

  {
    G4NeutronHPVector * result = new G4NeutronHPVector;
    G4int j=0;
    G4double x;
    G4double y;
    G4int running = 0;
    for(G4int i=0; i<left.GetVectorLength(); i++)
    {
      while(j<right.GetVectorLength())
      {
        if(right.GetX(j)<left.GetX(i)*1.001)
        {
          x = right.GetX(j);
          y = right.GetY(j)+left.GetY(x);
          result->SetData(running++, x, y);
          j++;
        }
        //else if(std::abs((right.GetX(j)-left.GetX(i))/(left.GetX(i)+right.GetX(j)))>0.001)
        else if( left.GetX(i)+right.GetX(j) == 0 
              || std::abs((right.GetX(j)-left.GetX(i))/(left.GetX(i)+right.GetX(j))) > 0.001 )
        {
          x = left.GetX(i);
          y = left.GetY(i)+right.GetY(x);
          result->SetData(running++, x, y);
          break;
        }
        else
        {
          break;
        }
      }
      if(j==right.GetVectorLength())
      {
        x = left.GetX(i);
        y = left.GetY(i)+right.GetY(x);
        result->SetData(running++, x, y);     
      }
    }
    result->ThinOut(0.02);
    return *result;
  }

---

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

• G4NeutronHPVector.hh
• G4NeutronHPVector.cc