#include "GaussianEliminationAlg.h"

Public Member Functions
	GaussianEliminationAlg (float, float)

double	GetDistance (float) const

const std::vector< float > &	SolveEquations (const std::vector< float > &meanVector, const std::vector< float > &sigmaVector, const std::vector< float > &heightVector)

void	FillAugmentedMatrix (const std::vector< float > &meanVector, const std::vector< float > &sigmaVector, const std::vector< float > &heightVector)

void	GaussianElimination ()

const std::vector< float > &	GetSolutions ()

void	Print ()

Private Member Functions
void	FillDistanceLookupTable ()

Private Attributes
float	fDistanceStepSize

float	fDistanceMax

std::vector< double >	fDistanceLookupTable

std::vector< std::vector< double > >	fMatrix

std::vector< float >	fSolutions

Detailed Description

Definition at line 20 of file GaussianEliminationAlg.h.

Constructor & Destructor Documentation

util::GaussianEliminationAlg::GaussianEliminationAlg	(	float	step,
		float	max
	)

Title: GaussianEliminationAlg Class Author: Wes Ketchum (wketc.nosp@m.hum@.nosp@m.lanl..nosp@m.gov)

Description: Class that solves system of linear equations via Gaussian Elimination. Intended for use with RFFHitFitter

Definition at line 17 of file GaussianEliminationAlg.cxx.

References fDistanceMax, fDistanceStepSize, FillDistanceLookupTable(), and max.

 {
   fDistanceStepSize = step;
   fDistanceMax = max;
 
   if(fDistanceStepSize<0 || fDistanceMax<0)
     throw std::runtime_error("Error in GaussianEliminationAlg: Cannot construct with negative step or max.");
 
   FillDistanceLookupTable();
 }

Member Function Documentation

void util::GaussianEliminationAlg::FillAugmentedMatrix	(	const std::vector< float > &	meanVector,
		const std::vector< float > &	sigmaVector,
		const std::vector< float > &	heightVector
	)

Definition at line 72 of file GaussianEliminationAlg.cxx.

References fMatrix, and GetDistance().

Referenced by SolveEquations().

 {
 
   fMatrix.resize(meanVector.size());
   for(size_t i=0; i<meanVector.size(); i++){
     
     fMatrix[i].resize(meanVector.size()+1);
     for(size_t j=0; j<meanVector.size(); j++){
       if(sigmaVector[j] < std::numeric_limits<float>::epsilon()){
         if(i==j)
           fMatrix[i][j] = 1.0;
         else
           fMatrix[i][j] = 0.0;
       }
       else
         fMatrix[i][j] = GetDistance( (meanVector[i]-meanVector[j])/sigmaVector[j] );
     }
     fMatrix[i][meanVector.size()] = heightVector[i];
     
   }
 
 }

void util::GaussianEliminationAlg::FillDistanceLookupTable ( )

private

Definition at line 28 of file GaussianEliminationAlg.cxx.

References fDistanceLookupTable, fDistanceMax, and fDistanceStepSize.

Referenced by GaussianEliminationAlg().

 {
 
   fDistanceLookupTable.clear();
   fDistanceLookupTable.reserve(std::ceil(fDistanceMax/fDistanceStepSize)+2);
 
   //let's be extra safe, and make sure we have more sampling points than we will need
   double x_val=0.0;
   while(x_val <= fDistanceMax+fDistanceStepSize){
     fDistanceLookupTable.push_back( std::exp(x_val*x_val*0.5*-1) );
     x_val += fDistanceStepSize;
   }
   //do one more to be sure to push beyond...
   fDistanceLookupTable.push_back( std::exp(x_val*x_val*0.5*-1) ); 
 }

void util::GaussianEliminationAlg::GaussianElimination ( )

Definition at line 97 of file GaussianEliminationAlg.cxx.

References fMatrix, and fSolutions.

Referenced by SolveEquations().

 {
 
   fSolutions.resize(fMatrix.size(),0.0);
 
   for(size_t i=0; i<fMatrix.size(); i++){
 
     for(size_t j=i+1; j<(fMatrix[i].size()-1); j++){
       float scale_value = fMatrix[j][i] / fMatrix[i][i];
       for(size_t k=i; k<fMatrix[i].size(); k++)
         fMatrix[j][k] -= fMatrix[i][k]*scale_value;      
     }//end column loop
   
   }//end row loop
 
   for(int i=fMatrix.size()-1; i>=0; i--){
     fSolutions[i] = fMatrix[i].back();
 
     for(size_t j=i+1; j<fMatrix.size(); j++)
       fSolutions[i] -= fMatrix[i][j]*fSolutions[j];
 
     fSolutions[i] /= fMatrix[i][i];
   }
 
 }

double util::GaussianEliminationAlg::GetDistance ( float d ) const

Definition at line 44 of file GaussianEliminationAlg.cxx.

References fDistanceLookupTable, fDistanceMax, and fDistanceStepSize.

Referenced by FillAugmentedMatrix().

 {
   double d_abs = std::abs(d);
   if(d_abs > fDistanceMax)
     return 0.0;
 
   size_t low_bin = std::floor(d_abs/fDistanceStepSize);
   return fDistanceLookupTable[low_bin] -
     (d_abs/fDistanceStepSize-(double)low_bin)*(fDistanceLookupTable[low_bin]-fDistanceLookupTable[low_bin+1]);
   
 }

const std::vector<float>& util::GaussianEliminationAlg::GetSolutions ( )

inline

Definition at line 34 of file GaussianEliminationAlg.h.

References fSolutions, and Print().

34 { return fSolutions; }

util::GaussianEliminationAlg::fSolutions

std::vector< float > fSolutions

Definition: GaussianEliminationAlg.h:46

void util::GaussianEliminationAlg::Print ( )

Definition at line 123 of file GaussianEliminationAlg.cxx.

References fDistanceLookupTable, fDistanceMax, fDistanceStepSize, fMatrix, and fSolutions.

Referenced by GetSolutions().

 {
   std::cout << "GaussianEliminationAlg." << std::endl;
 
   std::cout << "\tLookup table (step=" << fDistanceStepSize << ", max=" << fDistanceMax << ")" << std::endl;
   for(size_t i=0; i<fDistanceLookupTable.size(); i++)
     std::cout << "\t\tGaussian(" << fDistanceStepSize*i << ") = " << fDistanceLookupTable[i] << std::endl;
 
 
   std::cout << "\tAugmented matrix " << std::endl;
   for(size_t i=0; i<fMatrix.size(); i++){
     std::cout << "\t\t | ";
     for(size_t j=0; j<fMatrix[i].size()-1; j++)
       std::cout << fMatrix[i][j] << " ";
     std::cout << " | " << fMatrix[i][fMatrix[i].size()-1] << " |" << std::endl;
   }
 
   std::cout << "\tSolutions" << std::endl;
   for(size_t i=0; i<fSolutions.size(); i++)
     std::cout << "\t\t" << i << " " << fSolutions[i] << std::endl;
 }

const std::vector< float > & util::GaussianEliminationAlg::SolveEquations	(	const std::vector< float > &	meanVector,
		const std::vector< float > &	sigmaVector,
		const std::vector< float > &	heightVector
	)

Definition at line 56 of file GaussianEliminationAlg.cxx.

References FillAugmentedMatrix(), fSolutions, and GaussianElimination().

Referenced by hit::RFFHitFitter::CalculateAmplitudes().

 {
 
   if(meanVector.size()!=sigmaVector.size() ||
      meanVector.size()!=heightVector.size())
     throw std::runtime_error("Error in GaussianEliminationAlg: Vector sizes not equal.");
 
   FillAugmentedMatrix(meanVector,sigmaVector,heightVector);
   GaussianElimination();
 
   return fSolutions;
 
 }