trf++/include/trfbase
diff -N TrfVector_LA.h
--- TrfVector_LA.h 7 Jul 2010 16:10:38 -0000 1.1
+++ /dev/null 1 Jan 1970 00:00:00 -0000
@@ -1,429 +0,0 @@
-// TrfVector_LA.h
-
-#ifndef TrfVector_LA_H
-#define TrfVector_LA_H
-
-// Version of TrfVector using FNAL ZOOM linearAlgebra package.
-
-// Contains headers for underlying TrfVector, TrfMatrix and TrfSMatrix
-// classes and corresponding typedefs.
-
-#include <cassert>
-#include <cmath>
-#include <vector>
-#include "LinearAlgebra/ColumnVector.h"
-#include "LinearAlgebra/Matrix.h"
-#include "trfutil/TRFMath.h"
-
-//**********************************************************************
-
-class TrfVector : public ColumnVector {
-
-public: // methods
-
- // Constructor.
- explicit TrfVector(int len) : ColumnVector(len) { }
-
- // Constructor from subclass.
- TrfVector(const ColumnVector& vec)
- : ColumnVector(vec) { }
-
- // Constructor from subclass.
- TrfVector(const MatrixD& mtx) : ColumnVector(mtx.rows()) {
- assert( mtx.columns() == 1 );
- MatrixD::operator=(mtx);
- }
-
- // Return the length.
- int length() const { return rows(); }
-
- // Return the length.
- int nrow() const { return rows(); }
-
- // Fill the vector.
- void fill(double val) {
- for ( int i=0; i<rows(); ++i ) operator()(i) = val;
- }
-
- // Return the minimum value in the vector.
- double min() const {
- double val = operator()(0);
- for ( int i=1; i<rows(); ++i ) {
- double newval = operator()(i);
- if ( newval < val ) val = newval;
- }
- return val;
- }
-
- // Return the maximum value in the vector.
- double max() const {
- double val = operator()(0);
- for ( int i=1; i<rows(); ++i ) {
- double newval = operator()(i);
- if ( newval > val ) val = newval;
- }
- return val;
- }
-
- // Return the minimum absolute value in the vector.
- double amin() const {
- double val = std::abs(operator()(0));
- for ( int i=1; i<rows(); ++i ) {
- double newval = std::abs(operator()(i));
- if ( newval < val ) val = newval;
- }
- return val;
- }
-
- // Return the minimum absolute value in the vector.
- double amax() const {
- double val = std::abs(operator()(0));
- for ( int i=1; i<rows(); ++i ) {
- double newval = std::abs(operator()(i));
- if ( newval > val ) val = newval;
- }
- return val;
- }
-
- // Addition.
- TrfVector& operator+=(const TrfVector& rhs) {
- TrfVector& lhs = *this;
- assert( lhs.nrow() == rhs.nrow() );
- for ( int i=0; i<lhs.nrow(); ++i ) {
- lhs(i) += rhs(i);
- }
- return lhs;
- }
-
- // Subtraction.
- TrfVector& operator-=(const TrfVector& rhs) {
- TrfVector& lhs = *this;
- assert( lhs.nrow() == rhs.nrow() );
- for ( int i=0; i<lhs.nrow(); ++i ) {
- lhs(i) -= rhs(i);
- }
- return lhs;
- }
-
-};
-
-// Addition.
-inline
-TrfVector operator+(const TrfVector& rhs1, const TrfVector& rhs2) {
- TrfVector lhs = rhs1;
- lhs += rhs2;
- return lhs;
-}
-
-// Subtraction.
-inline
-TrfVector operator-(const TrfVector& rhs1, const TrfVector& rhs2) {
- TrfVector lhs = rhs1;
- lhs -= rhs2;
- return lhs;
-}
-
-// Equality.
-inline
-bool operator==(const TrfVector lhs, const TrfVector& rhs) {
- if ( lhs.nrow() != rhs.nrow() ) {
- return false;
- }
- for ( int i=0; i<lhs.nrow(); ++i ) {
- if ( lhs(i) != rhs(i) ) {
- return false;
- }
- }
- return true;
-}
-
-// Near equality.
-inline
-bool is_equal(const TrfVector lhs, const TrfVector& rhs) {
- if ( lhs.nrow() != rhs.nrow() ) {
- return false;
- }
- for ( int i=0; i<lhs.nrow(); ++i ) {
- if ( ! is_equal(lhs(i), rhs(i)) ) {
- return false;
- }
- }
- return true;
-}
-
-//**********************************************************************
-
-class TrfMatrix : public MatrixD {
-
-public: // methods
-
- // constructor
- explicit TrfMatrix(int ncol, int nrow)
- : MatrixD(ncol,nrow) { }
-
- // Constructor from subclass.
- TrfMatrix(const MatrixD& mtx)
- : MatrixD(mtx) { }
-
- // Return the length.
- int nrow() const { return rows(); }
- int ncol() const { return columns(); }
-
- // Fill the matrix.
- void fill(double val) {
- for ( int i=0; i<rows(); ++i ) {
- for ( int j=0; j<columns(); ++j ) {
- operator()(i,j) = val;
- }
- }
- }
-
- // Return the minimum value in the matrix.
- double min() const {
- double val = operator()(0,0);
- for ( int i=0; i<rows(); ++i ) {
- for ( int j=0; j<columns(); ++j ) {
- double newval = operator()(i,j);
- if ( newval < val ) val = newval;
- }
- }
- return val;
- }
-
- // Return the maximum value in the matrix.
- double max() const {
- double val = operator()(0,0);
- for ( int i=0; i<rows(); ++i ) {
- for ( int j=0; j<columns(); ++j ) {
- double newval = operator()(i,j);
- if ( newval > val ) val = newval;
- }
- }
- return val;
- }
-
- // Return the minimum absolute value in the matrix.
- double amin() const {
- double val = std::abs(operator()(0,0));
- for ( int i=0; i<rows(); ++i ) {
- for ( int j=0; j<columns(); ++j ) {
- double newval = std::abs(operator()(i,j));
- if ( newval < val ) val = newval;
- }
- }
- return val;
- }
-
- // Return the maximum absolute value in the matrix.
- double amax() const {
- double val = std::abs(operator()(0,0));
- for ( int i=0; i<rows(); ++i ) {
- for ( int j=0; j<columns(); ++j ) {
- double newval = std::abs(operator()(i,j));
- if ( newval > val ) val = newval;
- }
- }
- return val;
- }
-
- // Transpose.
- TrfMatrix transpose() const {
- return MatrixD::transpose();
- }
-
-};
-
-//**********************************************************************
-
-class TrfSMatrix : public MatrixD {
-
-public: // methods
-
- // constructor
- explicit TrfSMatrix(int dim) : MatrixD(dim,dim) {
- declareSymmetric();
- }
-
- // Constructor from subclass.
- TrfSMatrix(const MatrixD& mtx)
- : MatrixD(mtx) { }
-
- // Return the length.
- int nrow() const { return rows(); }
- int ncol() const { return columns(); }
-
- // Fill the matrix.
- void fill(double val) {
- for ( int i=0; i<rows(); ++i ) {
- for ( int j=0; j<=i; ++j ) {
- operator()(i,j) = val;
- }
- }
- }
-
- // Return the minimum value in the matrix.
- double min() const {
- double val = operator()(0,0);
- for ( int i=1; i<rows(); ++i ) {
- for ( int j=0; j<=i; ++j ) {
- double newval = operator()(i,j);
- if ( newval < val ) val = newval;
- }
- }
- return val;
- }
-
- // Return the maximum value in the matrix.
- double max() const {
- double val = operator()(0,0);
- for ( int i=1; i<rows(); ++i ) {
- for ( int j=0; j<i; ++j ) {
- double newval = operator()(i,j);
- if ( newval > val ) val = newval;
- }
- }
- return val;
- }
-
- // Return the minimum absolute value in the matrix.
- double amin() const {
- double val = std::abs(operator()(0,0));
- for ( int i=1; i<rows(); ++i ) {
- for ( int j=0; j<=i; ++j ) {
- double newval = std::abs(operator()(i,j));
- if ( newval > val ) val = newval;
- }
- }
- return val;
- }
-
- // Return the maximum absolute value in the matrix.
- double amax() const {
- double val = std::abs(operator()(0,0));
- for ( int i=1; i<rows(); ++i ) {
- for ( int j=0; j<=i; ++j ) {
- double newval = std::abs(operator()(i,j));
- if ( newval < val ) val = newval;
- }
- }
- return val;
- }
-
- // Xform S --> M * S * MT
- TrfSMatrix xform(const TrfMatrix& mtx) const {
- return xsxt(mtx);
- }
-
- // Xform S --> MT * S * M
- TrfSMatrix xformt(const TrfMatrix& mtx) const {
- return xtsx(mtx);
- }
-
-};
-
-// Normalize a matrix.
-// Transformation which makes all diagonal elements unity.
-inline
-void normalize(TrfSMatrix& sma) {
- std::vector<double> vec(sma.nrow());
- for ( int irow=0; irow<sma.nrow(); irow++ ) {
- vec[irow] = 1.0/std::sqrt( sma(irow,irow) );
- double rowfac = sqrt( vec[irow] );
- for ( int icol=0; icol<=irow; icol++ ) {
- sma(irow,icol) *= vec[irow]*vec[icol];
- }
- }
-}
-
-//**********************************************************************
-
-// Chi-square difference.
-inline
-double chisq_diff(const TrfVector& vecdiff, const TrfSMatrix& inverr) {
- MatrixD mtx = inverr.xtsx(vecdiff);
- assert( mtx.rows() == 1 );
- assert( mtx.columns() == 1 );
- return mtx(0,0);
-}
-
-// Matrix Equality.
-inline
-bool operator==(const TrfMatrix& lhs, const TrfMatrix& rhs) {
- if ( lhs.nrow() != rhs.nrow() ) {
- return false;
- }
- if ( lhs.ncol() != rhs.ncol() ) {
- return false;
- }
- for ( int i=0; i<lhs.nrow(); ++i ) {
- for ( int j=0; j<=lhs.ncol(); ++j ) {
- if ( lhs(i,j) != rhs(i,j) ) {
- return false;
- }
- }
- }
- return true;
-}
-
-// Symmetric matrix Equality.
-inline
-bool operator==(const TrfSMatrix& lhs, const TrfSMatrix& rhs) {
- if ( lhs.nrow() != rhs.nrow() ) {
- return false;
- }
- for ( int i=0; i<lhs.nrow(); ++i ) {
- for ( int j=0; j<=i; ++j ) {
- if ( lhs(i,j) != rhs(i,j) ) {
- return false;
- }
- }
- }
- return true;
-}
-
-// Matrix near equality.
-inline
-bool is_equal(const TrfMatrix& lhs, const TrfMatrix& rhs) {
- if ( lhs.nrow() != rhs.nrow() ) {
- return false;
- }
- if ( lhs.ncol() != rhs.ncol() ) {
- return false;
- }
- for ( int i=0; i<lhs.nrow(); ++i ) {
- for ( int j=0; j<=lhs.ncol(); ++j ) {
- if ( ! is_equal(lhs(i,j), rhs(i,j)) ) {
- return false;
- }
- }
- }
- return true;
-}
-
-// Near equality.
-inline
-bool is_equal(const TrfSMatrix& lhs, const TrfSMatrix& rhs) {
- if ( lhs.nrow() != rhs.nrow() ) {
- return false;
- }
- for ( int i=0; i<lhs.nrow(); ++i ) {
- for ( int j=0; j<=i; ++j ) {
- if ( ! is_equal(lhs(i,j), rhs(i,j)) ) {
- return false;
- }
- }
- }
- return true;
-}
-
-// Invert.
-// Return 0 for success.
-inline
-int invert(TrfSMatrix& mtx) {
- mtx = mtx.inverse();
- return 0;
-}
-
-//**********************************************************************
-
-#endif