#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, and FillDistanceLookupTable().

 {
   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 71 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 29 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 94 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 45 of file GaussianEliminationAlg.cxx.

References util::abs(), 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:45

void util::GaussianEliminationAlg::Print ( )

Definition at line 119 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;
 }