diff --git a/InnerDetector/InDetExample/InDetAlignExample/CalibrationLoopScripts/Solve_trf.py b/InnerDetector/InDetExample/InDetAlignExample/CalibrationLoopScripts/Solve_trf.py
index bc2a3f9447092c14d063b9c4e83b7100c0c8ffb4..970f470275224d8b74cbe2c8ff10e224387cc28c 100644
--- a/InnerDetector/InDetExample/InDetAlignExample/CalibrationLoopScripts/Solve_trf.py
+++ b/InnerDetector/InDetExample/InDetAlignExample/CalibrationLoopScripts/Solve_trf.py
@@ -194,7 +194,7 @@ solvingOption = 1
## solving options
## 0 - No global solving
## 1 - Lapack
-## 2 - ***REMOVED***
+## 2 - Eigen
## 6 - ROOT
## 7 - CLHEP
@@ -342,7 +342,7 @@ solvingOption = 1
## solving options
## 0 - No global solving
## 1 - Lapack
-## 2 - ***REMOVED***
+## 2 - Eigen
## 6 - ROOT
## 7 - CLHEP
diff --git a/InnerDetector/InDetExample/InDetAlignExample/python/IteratorGridClasses.py b/InnerDetector/InDetExample/InDetAlignExample/python/IteratorGridClasses.py
index c7e37a97b9884d9d32b00f97b536286ad7379a00..c4a8c776d30628a80426c388bdb08f52103ddada 100644
--- a/InnerDetector/InDetExample/InDetAlignExample/python/IteratorGridClasses.py
+++ b/InnerDetector/InDetExample/InDetAlignExample/python/IteratorGridClasses.py
@@ -418,7 +418,7 @@ class writeScriptGridForTFile :
("pixelAlignmentLevel" in self.AlignmentOptions and self.AlignmentOptions["pixelAlignmentLevel"] == 3) or
("pixelAlignmentLevelBarrel" in self.AlignmentOptions and self.AlignmentOptions["pixelAlignmentLevelBarrel"] == 3) or
("pixelAlignmentLevelEndcaps" in self.AlignmentOptions and self.AlignmentOptions["pixelAlignmentLevelEndcaps"] == 3) ) :
- print "hmn, you are going to run L3 alignment, ***REMOVED*** is using, so no eigen informations!!! "
+ print "hmn, you are going to run L3 alignment, Eigen is being used, so no eigen value information!!! "
else :
@@ -1103,7 +1103,7 @@ class writeScriptGrid :
("pixelAlignmentLevel" in self.AlignmentOptions and self.AlignmentOptions["pixelAlignmentLevel"] == 3) or
("pixelAlignmentLevelBarrel" in self.AlignmentOptions and self.AlignmentOptions["pixelAlignmentLevelBarrel"] == 3) or
("pixelAlignmentLevelEndcaps" in self.AlignmentOptions and self.AlignmentOptions["pixelAlignmentLevelEndcaps"] == 3) ) :
- print "hmn, you are going to run L3 alignment, ***REMOVED*** is using, so no eigen informations!!! "
+ print "hmn, you are going to run L3 alignment, Eigen is going to be used, so no eigen value information!!! "
else :
if ( "writeEigenMat" in self.GridOptions and self.GridOptions["writeEigenMat"] == True ) and ( self.AlignmentOptions["runLocal"] == False ):
diff --git a/InnerDetector/InDetExample/InDetAlignExample/share/AlignmentATNSimple.py b/InnerDetector/InDetExample/InDetAlignExample/share/AlignmentATNSimple.py
index c3bf2fb0c60f8ce004e3a7cd80b7cab018bc031e..936b0d37d2089aaaceb87dd7cb9b7c73c9cd05f0 100644
--- a/InnerDetector/InDetExample/InDetAlignExample/share/AlignmentATNSimple.py
+++ b/InnerDetector/InDetExample/InDetAlignExample/share/AlignmentATNSimple.py
@@ -256,7 +256,7 @@ runLocal = False
## solving option
## 0 - No global solving
## 1 - Lapack
-## 2 - ***REMOVED***
+## 2 - Eigen
## 6 - ROOT
## 7 - CLHEP
solvingOption = 1
diff --git a/InnerDetector/InDetExample/InDetAlignExample/share/InDetAlignAlgSetup_DBM.py b/InnerDetector/InDetExample/InDetAlignExample/share/InDetAlignAlgSetup_DBM.py
index 7a9ba8861b2ebddc3eed186744314ca9627a9e78..74cd5bc82ebcc07fcb3275031d6810b104a387df 100644
--- a/InnerDetector/InDetExample/InDetAlignExample/share/InDetAlignAlgSetup_DBM.py
+++ b/InnerDetector/InDetExample/InDetAlignExample/share/InDetAlignAlgSetup_DBM.py
@@ -35,7 +35,7 @@ newInDetAlignAlg_Options = {
,"inputMatrixFiles" : [ "matrix.bin" ] # list of matrix files when solving only
,"inputVectorFiles" : [ "vector.bin" ] # list of vector files when solving only
,"inputTFiles" : ["AlignmentTFile.root"] # list of the Alignment TFiles, used only if WriteTFile is True
- ,"solvingOption" : 3 # which global solver to run 0-none, 1-Lapack, 2-***REMOVED***, 6-ROOT, 7-CLHEP
+ ,"solvingOption" : 3 # which global solver to run 0-none, 1-Lapack, 2-Eigen, 6-ROOT, 7-CLHEP
,"solveLocal" : True # whether to run locaal solving
,"writeMatrixFile" : True # whether to write matrix to file before solving
,"writeMatrixFileTxt" : True # whether to write matrix to text file
@@ -264,7 +264,7 @@ if not newInDetAlignAlg_Options["runLocal"]:
):
matrixTool.WriteEigenMat = False
matrixTool.UseSparse = True
- matrixTool.SolveOption = 2 # run ***REMOVED*** for L3 by default
+ matrixTool.SolveOption = 2 # run Eigen for L3 by default
matrixTool.WriteMatTxt = False
matrixTool.WriteHitmapTxt = False
matrixTool.CalculateFullCovariance = False
diff --git a/InnerDetector/InDetExample/InDetAlignExample/share/NewInDetAlignAlgSetup.py b/InnerDetector/InDetExample/InDetAlignExample/share/NewInDetAlignAlgSetup.py
index f021c34f418dae94f74e9c9e039d484ddc279bd7..64ef5afcffc4529e74f2d8f9ac100aadbc688575 100644
--- a/InnerDetector/InDetExample/InDetAlignExample/share/NewInDetAlignAlgSetup.py
+++ b/InnerDetector/InDetExample/InDetAlignExample/share/NewInDetAlignAlgSetup.py
@@ -41,7 +41,7 @@ newInDetAlignAlg_Options = {
,"inputMatrixFiles" : [ "matrix.bin" ] # list of matrix files when solving only
,"inputVectorFiles" : [ "vector.bin" ] # list of vector files when solving only
,"inputTFiles" : ["AlignmentTFile.root"] # list of the Alignment TFiles, used only if WriteTFile is True
- ,"solvingOption" : 3 # which global solver to run 0-none, 1-Lapack, 2-***REMOVED***, 6-ROOT, 7-CLHEP
+ ,"solvingOption" : 3 # which global solver to run 0-none, 1-Lapack, 2-Eigen, 6-ROOT, 7-CLHEP
,"solveLocal" : True # whether to run locaal solving
,"writeMatrixFile" : True # whether to write matrix to file before solving
,"writeMatrixFileTxt" : True # whether to write matrix to text file
@@ -307,7 +307,7 @@ if not newInDetAlignAlg_Options["runLocal"]:
):
matrixTool.WriteEigenMat = False
matrixTool.UseSparse = True
- matrixTool.SolveOption = 2 # run ***REMOVED*** for L3 by default
+ matrixTool.SolveOption = 2 # run Eigen for L3 by default
matrixTool.WriteMatTxt = False
matrixTool.WriteHitmapTxt = False
matrixTool.CalculateFullCovariance = False
diff --git a/InnerDetector/InDetExample/InDetAlignExample/share/NewInDetAlignment_DBMRel17.py b/InnerDetector/InDetExample/InDetAlignExample/share/NewInDetAlignment_DBMRel17.py
index ec431f65f8f8dd7b6a61b9ed1dc65670231de2ff..c5e33ca483e6aefcc121a4b0c597c2498ea57351 100644
--- a/InnerDetector/InDetExample/InDetAlignExample/share/NewInDetAlignment_DBMRel17.py
+++ b/InnerDetector/InDetExample/InDetAlignExample/share/NewInDetAlignment_DBMRel17.py
@@ -156,7 +156,7 @@ runLocal = True
## solving option
## 0 - No global solving
## 1 - Lapack
-## 2 - ***REMOVED***
+## 2 - Eigen
## 6 - ROOT
## 7 - CLHEP
solvingOption = 1
diff --git a/InnerDetector/InDetExample/InDetAlignExample/share/NewInDetAlignment_Run2Rel17.py b/InnerDetector/InDetExample/InDetAlignExample/share/NewInDetAlignment_Run2Rel17.py
index 99c1b14041ede80ce48bcc7c477bb6df7cbd1a9b..c9574e5b9dd6c8a7568d8819791aab2e5fb7ea6b 100644
--- a/InnerDetector/InDetExample/InDetAlignExample/share/NewInDetAlignment_Run2Rel17.py
+++ b/InnerDetector/InDetExample/InDetAlignExample/share/NewInDetAlignment_Run2Rel17.py
@@ -326,7 +326,7 @@ runLocal = False
## solving option
## 0 - No global solving
## 1 - Lapack
-## 2 - ***REMOVED***
+## 2 - Eigen
## 6 - ROOT
## 7 - CLHEP
solvingOption = 1
diff --git a/InnerDetector/InDetExample/InDetAlignExample/share/NewInDetAlignment_Run2Rel19.py b/InnerDetector/InDetExample/InDetAlignExample/share/NewInDetAlignment_Run2Rel19.py
index 1d6d95783e0871c6a97679e2e584d00d52d47f3a..b1b7c21021336b32251290dfcf8de0a78be23bac 100644
--- a/InnerDetector/InDetExample/InDetAlignExample/share/NewInDetAlignment_Run2Rel19.py
+++ b/InnerDetector/InDetExample/InDetAlignExample/share/NewInDetAlignment_Run2Rel19.py
@@ -328,7 +328,7 @@ runLocal = True
## solving option
## 0 - No global solving
## 1 - Lapack
-## 2 - ***REMOVED***
+## 2 - Eigen
## 6 - ROOT
## 7 - CLHEP
solvingOption = 1
diff --git a/Tracking/TrkAlignment/TrkAlgebraUtils/CMakeLists.txt b/Tracking/TrkAlignment/TrkAlgebraUtils/CMakeLists.txt
index e280a0e1f104ccffd355c8afe016c0a0fe04e618..cb2f982db3fbebf1532a67f81c001c2fbd65de02 100644
--- a/Tracking/TrkAlignment/TrkAlgebraUtils/CMakeLists.txt
+++ b/Tracking/TrkAlignment/TrkAlgebraUtils/CMakeLists.txt
@@ -11,14 +11,16 @@ atlas_depends_on_subdirs( PRIVATE
TestPolicy )
# External dependencies:
+find_package( BLAS )
find_package( LAPACK )
+find_package( Eigen )
find_package( ROOT COMPONENTS Core Matrix Tree MathCore Hist RIO pthread )
# Component(s) in the package:
atlas_add_library( TrkAlgebraUtils
src/*.cxx
- ***REMOVED***/*.f
PUBLIC_HEADERS TrkAlgebraUtils
- PRIVATE_INCLUDE_DIRS ${ROOT_INCLUDE_DIRS} ${LAPACK_INCLUDE_DIRS}
- PRIVATE_LINK_LIBRARIES ${ROOT_LIBRARIES} ${LAPACK_LIBRARIES} GaudiKernel )
+ PRIVATE_INCLUDE_DIRS ${ROOT_INCLUDE_DIRS} ${LAPACK_INCLUDE_DIRS} ${BLAS_INCLUDE_DIRS} ${EIGEN_INCLUDE_DIRS}
+ PRIVATE_LINK_LIBRARIES ${ROOT_LIBRARIES} ${LAPACK_LIBRARIES} ${BLAS_LIBRARIES} ${EIGEN_LIBRARIES} GaudiKernel ${CMAKE_Fortran_IMPLICIT_LINK_LIBRARIES})
+
diff --git a/Tracking/TrkAlignment/TrkAlgebraUtils/TrkAlgebraUtils/AlMat.h b/Tracking/TrkAlignment/TrkAlgebraUtils/TrkAlgebraUtils/AlMat.h
index 1aa61aa11399d99dc22fffe9b3296993b87f77bc..3c495e3e75f87411789a710879dc9d4a0310022d 100644
--- a/Tracking/TrkAlignment/TrkAlgebraUtils/TrkAlgebraUtils/AlMat.h
+++ b/Tracking/TrkAlignment/TrkAlgebraUtils/TrkAlgebraUtils/AlMat.h
@@ -7,6 +7,7 @@
#include <iostream>
#include <string>
+#include <exception>
class StatusCode;
@@ -37,8 +38,8 @@ class AlMat {
inline AlMat_row(AlMat&,int);
double & operator[](int);
private:
- AlMat& _a;
- int _r;
+ AlMat& m_a;
+ int m_r;
};
class AlMat_row_const {
@@ -46,8 +47,8 @@ class AlMat {
inline AlMat_row_const(const AlMat&,int);
const double & operator[](int) const;
private:
- const AlMat& _a;
- int _r;
+ const AlMat& m_a;
+ int m_r;
};
// helper class to implement m[i][j]
@@ -104,12 +105,12 @@ class AlMat {
void copy(const AlSymMatBase& m);
protected:
- int _ncol, _nrow;
- int _nele;
+ int m_ncol, m_nrow;
+ int m_nele;
std::string m_pathbin;
std::string m_pathtxt;
- double* ptr_data;
- bool transpose;
+ double* m_ptr_data;
+ bool m_transpose;
virtual long int elem(long int, long int) const;
@@ -128,41 +129,41 @@ inline AlMat::AlMat_row_const AlMat::operator[] (int r) const{
}
inline double &AlMat::AlMat_row::operator[](int c){
- if(_r<0||_r>=_a.nrow() || c<0||c>=_a.ncol())
- std::cerr << "Range error in AlMat::operator[][]" << std::endl;
+ if(m_r<0||m_r>=m_a.nrow() || c<0||c>=m_a.ncol())
+ throw std::out_of_range( "Range error in AlMat::operator[][]" );
- return *(_a.ptr_data+_r*_a.ncol()+c);
+ return *(m_a.m_ptr_data+m_r*m_a.ncol()+c);
}
inline const double & AlMat::AlMat_row_const::operator[](int c) const {
- if(_r<0||_r>=_a.nrow() || c<0||c>=_a.ncol())
- std::cerr << "Range error in AlMat::operator[][]" << std::endl;
+ if(m_r<0||m_r>=m_a.nrow() || c<0||c>=m_a.ncol())
+ throw std::out_of_range( "Range error in AlMat::operator[][]" );
- return *(_a.ptr_data+_r*_a.ncol()+c);
+ return *(m_a.m_ptr_data+m_r*m_a.ncol()+c);
}
inline AlMat::AlMat_row::AlMat_row(AlMat &a,int r)
- : _a(a), _r(r)
+ : m_a(a), m_r(r)
{}
inline AlMat::AlMat_row_const::AlMat_row_const
(const AlMat&a,int r)
- : _a(a), _r(r)
+ : m_a(a), m_r(r)
{}
inline long int AlMat::nrow() const {
- if( transpose ) return _ncol;
- return _nrow;
+ if( m_transpose ) return m_ncol;
+ return m_nrow;
}
inline long int AlMat::ncol() const {
- if( transpose ) return _nrow;
- return _ncol;
+ if( m_transpose ) return m_nrow;
+ return m_ncol;
}
inline double* AlMat::ptrData() const {
- return ptr_data;
+ return m_ptr_data;
}
inline std::string AlMat::pathBin() const {
diff --git a/Tracking/TrkAlignment/TrkAlgebraUtils/TrkAlgebraUtils/AlSpaMat.h b/Tracking/TrkAlignment/TrkAlgebraUtils/TrkAlgebraUtils/AlSpaMat.h
index 4a771249ca13eb24c13765eef926ba0fd028f7c7..6ae73cdfd775aa61b811b06bf7748ce00078e523 100644
--- a/Tracking/TrkAlignment/TrkAlgebraUtils/TrkAlgebraUtils/AlSpaMat.h
+++ b/Tracking/TrkAlignment/TrkAlgebraUtils/TrkAlgebraUtils/AlSpaMat.h
@@ -7,7 +7,7 @@
#include "TrkAlgebraUtils/AlSymMatBase.h"
-#include <iostream>
+#include <exception>
class StatusCode;
template<class Element> class TMatrixTSparse;
@@ -50,7 +50,7 @@ class AlSpaMat : public AlSymMatBase {
virtual AlSpaMat& operator*=(const double&);
// ADVANCED:
- int ***REMOVED***Solve(AlVec& RHS);
+ int SolveWithEigen(AlVec& RHS);
virtual void RemoveModule( int);
virtual void RemoveAlignPar(int, int);
virtual void RemoveDoF( int, int nelem=1);
@@ -85,8 +85,8 @@ class AlSpaMat : public AlSymMatBase {
std::string m_pathbin;
std::string m_pathtxt;
- int* ptr_row;
- int* ptr_col;
+ int* m_ptr_row;
+ int* m_ptr_col;
};
@@ -95,15 +95,15 @@ class AlSpaMat : public AlSymMatBase {
inline void AlSpaMat::elem(const indices& key, long int& i, long int& j) const {
i = key.first;
j = key.second;
- if(j>i) std::cerr << "AlSpaMat::elem: i<j ! " << std::endl;
+ if(j>i) throw std::out_of_range( "AlSpaMat::elem: i<j ! " );
return;
}
inline long int AlSpaMat::nele() {
// first check for intrinsic consistency:
- if( int(ptr_map.size()) != _nele ) std::cerr << "AlSpaMat::_nele has been corrupted!" << std::endl;
- _nele = ptr_map.size();
- return _nele;
+ if( int(m_ptr_map.size()) != m_nele ) throw std::range_error( "AlSpaMat::m_nele has been corrupted!" );
+ m_nele = m_ptr_map.size();
+ return m_nele;
}
inline std::string AlSpaMat::pathBin() const {
diff --git a/Tracking/TrkAlignment/TrkAlgebraUtils/TrkAlgebraUtils/AlSymMat.h b/Tracking/TrkAlignment/TrkAlgebraUtils/TrkAlgebraUtils/AlSymMat.h
index f70a5c9d929055e761960f1caf636ed7fbfc87f2..7c3e717c85ecc7b011b36d55f1e2ec9c698b8a09 100644
--- a/Tracking/TrkAlignment/TrkAlgebraUtils/TrkAlgebraUtils/AlSymMat.h
+++ b/Tracking/TrkAlignment/TrkAlgebraUtils/TrkAlgebraUtils/AlSymMat.h
@@ -83,7 +83,7 @@ class AlSymMat : public AlSymMatBase {
virtual void copy(const AlSpaMat& m);
virtual void copy(const AlMat& m);
- double* ptr_data;
+ double* m_ptr_data;
std::string m_pathbin;
std::string m_pathtxt;
@@ -94,10 +94,10 @@ class AlSymMat : public AlSymMatBase {
inline long int AlSymMat::elem(long int i,long int j) const {
#ifdef _DEBUG
- if( i<0 ) { std::cerr << "AlSymMat::elem: Index 1 < zero! " << i << std::endl; return *(ptr_data); };
- if( i>=size() ) { std::cerr << "AlSymMat::elem: Index 1 too large! " << i << std::endl; return *(ptr_data); };
- if( j<0 ) { std::cerr << "AlSymMat::elem: Index 2 < zero! " << j << std::endl; return *(ptr_data); };
- if( j>=size() ) { std::cerr << "AlSymMat::elem: Index 2 too large! " << j << std::endl; return *(ptr_data); };
+ if( i<0 ) { throw std::out_of_range( "AlSymMat::elem: Index 1 < zero! " ); };
+ if( i>=size() ) { throw std::out_of_range( "AlSymMat::elem: Index 1 too large! " ); };
+ if( j<0 ) { throw std::out_of_range( "AlSymMat::elem: Index 2 < zero! " ); };
+ if( j>=size() ) { throw std::out_of_range( "AlSymMat::elem: Index 2 too large! " ); };
#endif
if( j<=i ) {
return ((i+1)*i/2+j);
@@ -108,7 +108,7 @@ inline long int AlSymMat::elem(long int i,long int j) const {
}
inline double* AlSymMat::ptrData() const {
- return ptr_data;
+ return m_ptr_data;
}
inline std::string AlSymMat::pathBin() const {
diff --git a/Tracking/TrkAlignment/TrkAlgebraUtils/TrkAlgebraUtils/AlSymMatBase.h b/Tracking/TrkAlignment/TrkAlgebraUtils/TrkAlgebraUtils/AlSymMatBase.h
index 93ee825e8c9e980b62b155209001842809184c08..47353fe1a6554bd92b84b491ade316cc375589cc 100644
--- a/Tracking/TrkAlignment/TrkAlgebraUtils/TrkAlgebraUtils/AlSymMatBase.h
+++ b/Tracking/TrkAlignment/TrkAlgebraUtils/TrkAlgebraUtils/AlSymMatBase.h
@@ -9,7 +9,7 @@
// modified from SiGlobalChi2AlgebrUtils::AlSymMatBase, modified to make
// pure virtual
-#include <iostream>
+#include <exception>
#include <map>
#include <vector>
@@ -46,8 +46,8 @@ class AlSymMatBase {
inline AlSymMatBase_row(AlSymMatBase&,long int);
double & operator[](long int);
private:
- AlSymMatBase& _a;
- long int _r;
+ AlSymMatBase& m_a;
+ long int m_r;
};
class AlSymMatBase_row_const {
@@ -55,8 +55,8 @@ class AlSymMatBase {
inline AlSymMatBase_row_const(const AlSymMatBase&,long int);
double operator[](long int) const;
private:
- const AlSymMatBase& _a;
- long int _r;
+ const AlSymMatBase& m_a;
+ long int m_r;
};
// helper class to implement m[i][j]
@@ -97,11 +97,11 @@ class AlSymMatBase {
AlSymMatBase(const AlSymMatBase&);
AlSymMatBase& operator=(const AlSymMatBase&);
- int _matrix_type;
- datamap ptr_map;
+ int m_matrix_type;
+ datamap m_ptr_map;
- long int _size;
- long int _nele;
+ long int m_size;
+ long int m_nele;
};
@@ -118,38 +118,36 @@ inline AlSymMatBase::AlSymMatBase_row_const AlSymMatBase::operator[] (long int r
}
inline double & AlSymMatBase::AlSymMatBase_row::operator[](long int c) {
- if(_r<0||_r>=_a.nrow() || c<0||c>=_a.ncol()) {
- std::cerr << "Range error in AlSymMatBase::operator[][]" << std::endl;
- return _a.elemr(0,0);
+ if(m_r<0||m_r>=m_a.nrow() || c<0||c>=m_a.ncol()) {
+ throw std::out_of_range( "Range error in AlSymMatBase::operator[][]" );
} else {
- return _a.elemr(_r,c);
+ return m_a.elemr(m_r,c);
}
}
inline double AlSymMatBase::AlSymMatBase_row_const::operator[](long int c) const {
- if(_r<0||_r>=_a.nrow() || c<0||c>=_a.ncol()) {
- std::cerr << "Range error in AlSymMatBase::operator[][]" << std::endl;
- return _a.elemc(0,0);
+ if(m_r<0||m_r>=m_a.nrow() || c<0||c>=m_a.ncol()) {
+ throw std::out_of_range( "Range error in AlSymMatBase::operator[][]" );
} else {
- return _a.elemc(_r,c);
+ return m_a.elemc(m_r,c);
}
}
inline AlSymMatBase::AlSymMatBase_row::AlSymMatBase_row(AlSymMatBase &a,long int r)
- : _a(a), _r(r)
+ : m_a(a), m_r(r)
{}
inline AlSymMatBase::AlSymMatBase_row_const::AlSymMatBase_row_const
(const AlSymMatBase&a,long int r)
- : _a(a), _r(r)
+ : m_a(a), m_r(r)
{}
-inline long int AlSymMatBase::size() const { return _size; }
-inline long int AlSymMatBase::nrow() const { return _size; }
-inline long int AlSymMatBase::ncol() const { return _size; }
-inline int AlSymMatBase::matrix_type() const { return _matrix_type; }
-inline const datamap * AlSymMatBase::ptrMap() const { return &ptr_map; }
+inline long int AlSymMatBase::size() const { return m_size; }
+inline long int AlSymMatBase::nrow() const { return m_size; }
+inline long int AlSymMatBase::ncol() const { return m_size; }
+inline int AlSymMatBase::matrix_type() const { return m_matrix_type; }
+inline const datamap * AlSymMatBase::ptrMap() const { return &m_ptr_map; }
} // end namespace Trk
diff --git a/Tracking/TrkAlignment/TrkAlgebraUtils/TrkAlgebraUtils/AlVec.h b/Tracking/TrkAlignment/TrkAlgebraUtils/TrkAlgebraUtils/AlVec.h
index 18ba7c3f7999d6b0809f322aa4cfe6d5deeca9b1..e32c8ca5c359460968523476f2cce2a04ee70d99 100644
--- a/Tracking/TrkAlignment/TrkAlgebraUtils/TrkAlgebraUtils/AlVec.h
+++ b/Tracking/TrkAlignment/TrkAlgebraUtils/TrkAlgebraUtils/AlVec.h
@@ -7,7 +7,7 @@
// PBdR (17Apr2007)
-#include <iostream>
+#include <exception>
#include <map>
#include <vector>
@@ -66,8 +66,8 @@ class AlVec {
inline double* ptrData() const;
protected:
- int _size;
- double* ptr_data;
+ int m_size;
+ double* m_ptr_data;
std::string m_pathbin;
std::string m_pathtxt;
@@ -82,39 +82,35 @@ class AlVec {
// inline methods:
inline int AlVec::size() const {
- return _size;
+ return m_size;
}
inline double* AlVec::ptrData() const {
- return ptr_data;
+ return m_ptr_data;
}
inline double& AlVec::operator[](int i) {
if( i < 0 ) {
- std::cerr << "AlVec: Index < zero! " << std::endl;
- return ptr_data[0];
+ throw std::out_of_range( "AlVec: Index < zero! " );
}
- if( i >= _size ) {
- std::cerr << "AlVec: Index too large! " << std::endl;
- return ptr_data[0];
+ if( i >= m_size ) {
+ throw std::out_of_range( "AlVec: Index too large! ");
}
- return *(ptr_data+i);
+ return *(m_ptr_data+i);
}
inline const double& AlVec::operator[](int i) const {
if( i < 0 ) {
- std::cerr << "AlVec: Index < zero! " << std::endl;
- return ptr_data[0];
+ throw std::out_of_range( "AlVec: Index < zero! " );
}
- if( i >= _size ) {
- std::cerr << "AlVec: Index too large! " << std::endl;
- return ptr_data[0];
+ if( i >= m_size ) {
+ throw std::out_of_range( "AlVec: Index too large! " );
}
- return *(ptr_data+i);
+ return *(m_ptr_data+i);
}
} // end namespace Trk
diff --git a/Tracking/TrkAlignment/TrkAlgebraUtils/TrkAlgebraUtils/IntVec.h b/Tracking/TrkAlignment/TrkAlgebraUtils/TrkAlgebraUtils/IntVec.h
index 3d41b0cba467e79706a06df463db46c546449106..ef2ebec8caed37c04bedf5f559adb1aa5f35ace5 100644
--- a/Tracking/TrkAlignment/TrkAlgebraUtils/TrkAlgebraUtils/IntVec.h
+++ b/Tracking/TrkAlignment/TrkAlgebraUtils/TrkAlgebraUtils/IntVec.h
@@ -24,11 +24,11 @@ class IntVec {
IntVec& operator-=(const IntVec&);
void reSize(int);
- int n_elem() const { return Nele; }
+ int n_elem() const { return m_Nele; }
private:
- int Nele;
- int* ptr_data;
+ int m_Nele;
+ int* m_ptr_data;
};
} // end namespace Trk
diff --git a/Tracking/TrkAlignment/TrkAlgebraUtils/cmt/requirements b/Tracking/TrkAlignment/TrkAlgebraUtils/cmt/requirements
index c9abc1d7b0d019f6251d0186cd179f86763a49eb..37ec1df5871343209b18ecb64beaeb017f13167c 100644
--- a/Tracking/TrkAlignment/TrkAlgebraUtils/cmt/requirements
+++ b/Tracking/TrkAlignment/TrkAlgebraUtils/cmt/requirements
@@ -19,5 +19,5 @@ public
macro_append ROOT_linkopts " -lMatrix"
-library TrkAlgebraUtils *.cxx -s=../***REMOVED*** *.f
+library TrkAlgebraUtils *.cxx
apply_pattern installed_library
diff --git a/Tracking/TrkAlignment/TrkAlgebraUtils/src/AlMat.cxx b/Tracking/TrkAlignment/TrkAlgebraUtils/src/AlMat.cxx
index 3f0e87d02b046b9494c79a146c264594016efb24..27a4f310eee9dd1170de72aec27ecf1f002ab528 100644
--- a/Tracking/TrkAlignment/TrkAlgebraUtils/src/AlMat.cxx
+++ b/Tracking/TrkAlignment/TrkAlgebraUtils/src/AlMat.cxx
@@ -14,49 +14,50 @@
#include <fstream>
#include <sstream>
#include <stdint.h>
+#include <exception>
namespace Trk {
AlMat::AlMat() {
- _ncol = 0; _nrow = 0; _nele = 0;
- ptr_data = 0; // set pointer to null
- transpose = false;
+ m_ncol = 0; m_nrow = 0; m_nele = 0;
+ m_ptr_data = 0; // set pointer to null
+ m_transpose = false;
m_pathbin="./";
m_pathtxt="./";
}
AlMat::AlMat(int N, int M) {
- _nrow = N;
- _ncol = M;
- _nele = N*M;
- ptr_data = new double[_nele];
- transpose = false;
+ m_nrow = N;
+ m_ncol = M;
+ m_nele = N*M;
+ m_ptr_data = new double[m_nele];
+ m_transpose = false;
m_pathbin="./";
m_pathtxt="./";
- double* p = ptr_data + _nele;
- while (p > ptr_data) *(--p) = 0.;
+ double* p = m_ptr_data + m_nele;
+ while (p > m_ptr_data) *(--p) = 0.;
}
AlMat::AlMat(const AlMat& m) {
- _nrow = m.nrow(); _ncol = m.ncol(); _nele = _nrow*_ncol;
- ptr_data = new double[_nele];
- transpose = false;
+ m_nrow = m.nrow(); m_ncol = m.ncol(); m_nele = m_nrow*m_ncol;
+ m_ptr_data = new double[m_nele];
+ m_transpose = false;
m_pathbin = m.m_pathbin;
m_pathtxt = m.m_pathtxt;
copy(m);
}
AlMat::~AlMat()
-{if( ptr_data != NULL) delete [] ptr_data;}
+{if( m_ptr_data != nullptr ) delete [] m_ptr_data;}
void AlMat::copy(const AlMat& m) {
int nr = nrow();
int nc = ncol();
if( nr != m.nrow() || nc != m.ncol() ) {
- std::cerr << "AlMat::copy: sizes do not match!" << std::endl;
- return; }
+ throw std::range_error( "AlMat::copy: sizes do not match!" );
+ }
for( int i=0; i<nr; i++ ) {
for( int j=0; j<nc; j++ ) {
@@ -68,57 +69,57 @@ void AlMat::copy(const AlMat& m) {
void AlMat::copy(const AlSymMatBase& m) {
long int n = m.size();
if( nrow() != n || ncol() != n ) {
- std::cerr << "AlMat::copy: sizes do not match!" << std::endl;
- return; }
+ throw std::range_error( "AlMat::copy: sizes do not match!" );
+ }
for( int i=0; i<n; i++ ) {
for( int j=0; j<n; j++ ) {
- *(ptr_data+i*_ncol+j) = m.elemc(i, j);
+ *(m_ptr_data+i*m_ncol+j) = m.elemc(i, j);
}
}
}
double& AlMat::elemr(long int i,long int j) {
#ifdef _DEBUG
- if( i<0 ) { std::cerr << "AlMat::elemr: Index 1 < zero! " << i << std::endl; return *(ptr_data); };
- if( i>=size ) { std::cerr << "AlMat::elemr: Index 1 too large! " << i << std::endl; return *(ptr_data); };
- if( j<0 ) { std::cerr << "AlMat::elemr: Index 2 < zero! " << j << std::endl; return *(ptr_data); };
- if( j>=size ) { std::cerr << "AlMat::elemr: Index 2 too large! " << j << std::endl; return *(ptr_data); };
+ if( i<0 ) { throw std::out_of_range( "AlMat::elemr: Index 1 < zero! " ); };
+ if( i>=size ) { throw std::out_of_range( "AlMat::elemr: Index 1 too large!" ); };
+ if( j<0 ) { throw std::out_of_range( "AlMat::elemr: Index 2 < zero! " ); };
+ if( j>=size ) { throw std::out_of_range( "AlMat::elemr: Index 2 too large!" ); };
#endif
- if( transpose ) return *(ptr_data+j*_ncol+i);
- return *(ptr_data+i*_ncol+j);
+ if( m_transpose ) return *(m_ptr_data+j*m_ncol+i);
+ return *(m_ptr_data+i*m_ncol+j);
}
double AlMat::elemc(long int i,long int j) const {
#ifdef _DEBUG
- if( i<0 ) { std::cerr << "AlMat::elemc: Index 1 < zero! " << i << std::endl; return *(ptr_data); };
- if( i>=size ) { std::cerr << "AlMat::elemc: Index 1 too large! " << i << std::endl; return *(ptr_data); };
- if( j<0 ) { std::cerr << "AlMat::elemc: Index 2 < zero! " << j << std::endl; return *(ptr_data); };
- if( j>=size ) { std::cerr << "AlMat::elemc: Index 2 too large! " << j << std::endl; return *(ptr_data); };
+ if( i<0 ) { throw std::out_of_range( "AlMat::elemr: Index 1 < zero! " ); };
+ if( i>=size ) { throw std::out_of_range( "AlMat::elemr: Index 1 too large!" ); };
+ if( j<0 ) { throw std::out_of_range( "AlMat::elemr: Index 2 < zero! " ); };
+ if( j>=size ) { throw std::out_of_range( "AlMat::elemr: Index 2 too large!" ); };
#endif
- if( transpose ) return *(ptr_data+j*_ncol+i);
- return *(ptr_data+i*_ncol+j);
+ if( m_transpose ) return *(m_ptr_data+j*m_ncol+i);
+ return *(m_ptr_data+i*m_ncol+j);
}
long int AlMat::elem(long int i,long int j) const {
#ifdef _DEBUG
- if( i<0 ) { std::cerr << "AlMat::elem: Index 1 < zero! " << i << std::endl; return *(ptr_data); };
- if( i>=size ) { std::cerr << "AlMat::elem: Index 1 too large! " << i << std::endl; return *(ptr_data); };
- if( j<0 ) { std::cerr << "AlMat::elem: Index 2 < zero! " << j << std::endl; return *(ptr_data); };
- if( j>=size ) { std::cerr << "AlMat::elem: Index 2 too large! " << j << std::endl; return *(ptr_data); };
+ if( i<0 ) { throw std::out_of_range( "AlMat::elemr: Index 1 < zero! " ); };
+ if( i>=size ) { throw std::out_of_range( "AlMat::elemr: Index 1 too large!" ); };
+ if( j<0 ) { throw std::out_of_range( "AlMat::elemr: Index 2 < zero! " ); };
+ if( j>=size ) { throw std::out_of_range( "AlMat::elemr: Index 2 too large!" ); };
#endif
- if( transpose ) return (j*_ncol+i);
- return (i*_ncol+j);
+ if( m_transpose ) return (j*m_ncol+i);
+ return (i*m_ncol+j);
}
AlMat& AlMat::operator=(const double& d) {
- double* p = ptr_data + _nele;
- while (p > ptr_data) *(--p) = d;
+ double* p = m_ptr_data + m_nele;
+ while (p > m_ptr_data) *(--p) = d;
return *this;
}
@@ -126,10 +127,12 @@ AlMat& AlMat::operator=(const double& d) {
AlMat& AlMat::operator=(const AlMat& m) {
if (this==&m) return *this;
- if(( _nrow!=0 || _ncol!=0) && (_nrow != m.nrow() || _ncol != m.ncol() ))
- { std::cerr << "AlMat=AlMat Assignment: size do not match!" << std::endl; return *this; }
+ if(( m_nrow!=0 || m_ncol!=0) && (m_nrow != m.nrow() || m_ncol != m.ncol() ))
+ {
+ throw std::range_error( "AlMat=AlMat Assignment: size do not match!" );
+ }
- if ( ptr_data != m.ptr_data ) {
+ if ( m_ptr_data != m.m_ptr_data ) {
reSize(m.nrow(), m.ncol());
copy(m);
};
@@ -137,10 +140,12 @@ AlMat& AlMat::operator=(const AlMat& m) {
}
AlMat& AlMat::operator=(const AlSymMat& m) {
- if( ( _nrow!=0 || _ncol!=0) && (_nrow != m.size() || _ncol != m.size() ))
- { std::cerr << "AlMat=AlSymMatBase Assignment: size do not match!" << std::endl; return *this; }
+ if( ( m_nrow!=0 || m_ncol!=0) && (m_nrow != m.size() || m_ncol != m.size() ))
+ {
+ throw std::range_error( "AlMat=AlSymMatBase Assignment: size do not match!" );
+ }
- if ( ptr_data != m.ptrData() ) {
+ if ( m_ptr_data != m.ptrData() ) {
reSize(m.size(), m.size());
copy(m);
}
@@ -148,59 +153,69 @@ AlMat& AlMat::operator=(const AlSymMat& m) {
}
AlMat AlMat::operator+(const AlMat& m) const {
- if( _nrow != m._nrow || _ncol != m._ncol )
- { std::cerr << "AlMat: operator+: size do not match!" << std::endl; return *this; }
+ if( m_nrow != m.m_nrow || m_ncol != m.m_ncol )
+ {
+ throw std::range_error( "AlMat: operator+: size do not match!" );
+ }
- AlMat b( _nrow, _ncol );
+ AlMat b( m_nrow, m_ncol );
- double* p = ptr_data + _nele;
- double* q = m.ptr_data + _nele;
- double* r = b.ptr_data + _nele;
- while (p > ptr_data) *(--r) = (*(--p))+(*(--q));
+ double* p = m_ptr_data + m_nele;
+ double* q = m.m_ptr_data + m_nele;
+ double* r = b.m_ptr_data + m_nele;
+ while (p > m_ptr_data) *(--r) = (*(--p))+(*(--q));
return b;
}
AlMat& AlMat::operator+=(const AlMat& m) {
- if( _nrow != m._nrow || _ncol != m._ncol )
- { std::cerr << "AlMat: operator+=: size do not match!" << std::endl; return *this; }
+ if( m_nrow != m.m_nrow || m_ncol != m.m_ncol )
+ {
+ throw std::range_error( "AlMat: operator+=: size do not match!" );
+ }
- double* p = ptr_data + _nele;
- double* q = m.ptr_data + _nele;
- while (p > ptr_data) *(--p) += *(--q);
+ double* p = m_ptr_data + m_nele;
+ double* q = m.m_ptr_data + m_nele;
+ while (p > m_ptr_data) *(--p) += *(--q);
return *this;
}
AlMat AlMat::operator-(const AlMat& m) const {
- if( _nrow != m._nrow || _ncol != m._ncol )
- { std::cerr << "AlMat: operator-: size do not match!" << std::endl; return *this; }
+ if( m_nrow != m.m_nrow || m_ncol != m.m_ncol )
+ {
+ throw std::range_error( "AlMat: operator-: size do not match!" );
+ }
- AlMat b( _nrow, _ncol );
+ AlMat b( m_nrow, m_ncol );
- double* p = ptr_data + _nele;
- double* q = m.ptr_data + _nele;
- double* r = b.ptr_data + _nele;
- while (p > ptr_data) *(--r) = (*(--p))-(*(--q));
+ double* p = m_ptr_data + m_nele;
+ double* q = m.m_ptr_data + m_nele;
+ double* r = b.m_ptr_data + m_nele;
+ while (p > m_ptr_data) *(--r) = (*(--p))-(*(--q));
return b;
}
AlMat& AlMat::operator-=(const AlMat& m) {
- if( _nrow != m._nrow || _ncol != m._ncol )
- { std::cerr << "AlMat: operator-=: size do not match!" << std::endl; return *this; }
+ if( m_nrow != m.m_nrow || m_ncol != m.m_ncol )
+ {
+ throw std::range_error( "AlMat: operator-=: size do not match!");
+ }
- double* p = ptr_data + _nele;
- double* q = m.ptr_data + _nele;
- while (p > ptr_data) *(--p) -= *(--q);
+ double* p = m_ptr_data + m_nele;
+ double* q = m.m_ptr_data + m_nele;
+ while (p > m_ptr_data) *(--p) -= *(--q);
return *this;
}
AlMat AlMat::operator*(const AlMat& m) const {
- if( ncol() != m.nrow() ) { std::cerr << "AlMat: operator*: size do not match!" << std::endl;
- return m; }
+ if( ncol() != m.nrow() )
+ {
+ throw std::range_error( "AlMat: operator*: size do not match!" );
+ }
int k(nrow());
int l(m.ncol());
@@ -216,8 +231,10 @@ AlMat AlMat::operator*(const AlMat& m) const {
}
AlMat AlMat::operator*(const AlSymMatBase& m) const {
- if( ncol() != m.size()) { std::cerr << "AlMat: operator*: size do not match!" << std::endl;
- return *this; }
+ if( ncol() != m.size())
+ {
+ throw std::range_error( "AlMat: operator*: size do not match!" );
+ }
int k(nrow());
int l(m.size());
@@ -234,7 +251,9 @@ AlMat AlMat::operator*(const AlSymMatBase& m) const {
AlVec AlMat::operator*(const AlVec& v) const {
if( ncol() != v.size() )
- { std::cerr << "AlMat: operator*: size do not match! " << std::endl; return v; }
+ {
+ throw std::range_error( "AlMat: operator*: size do not match! " );
+ }
int k(nrow());
int l(ncol());
@@ -244,7 +263,7 @@ AlVec AlMat::operator*(const AlVec& v) const {
AlVec b(k);
for( int i=0; i<k; i++ ) {
p = b.ptrData()+i;
- q = ptr_data+i*l;
+ q = m_ptr_data+i*l;
for( int j=0; j<l; j++ ) *p += (*(q+j))*(*(v.ptrData()+j));
}
@@ -252,8 +271,8 @@ AlVec AlMat::operator*(const AlVec& v) const {
}
AlMat& AlMat::operator*=(const double& d) {
- double* p = ptr_data+_nele;
- while (p > ptr_data)
+ double* p = m_ptr_data+m_nele;
+ while (p > m_ptr_data)
*(--p) *= d;
return *this;
@@ -261,10 +280,10 @@ AlMat& AlMat::operator*=(const double& d) {
// transposition
AlMat AlMat::T() const {
- AlMat b(_ncol,_nrow);
- for( int i=0; i<b._nrow; i++ ) {
- for( int j=0; j<b._ncol; j++ ) {
- b[i][j]= *(ptr_data+j*_ncol+i);
+ AlMat b(m_ncol,m_nrow);
+ for( int i=0; i<b.m_nrow; i++ ) {
+ for( int j=0; j<b.m_ncol; j++ ) {
+ b[i][j]= *(m_ptr_data+j*m_ncol+i);
}
}
return b;
@@ -273,7 +292,7 @@ AlMat AlMat::T() const {
// transposition
AlMat& AlMat::Transpose() {
- transpose = true;
+ m_transpose = true;
return *this;
}
@@ -281,7 +300,7 @@ AlMat& AlMat::Transpose() {
// normal position
AlMat& AlMat::Normal() {
- transpose = false;
+ m_transpose = false;
return *this;
}
@@ -289,12 +308,11 @@ AlMat& AlMat::Normal() {
//invert sym matrix declared as non-symetric for convenience
void AlMat::invertS(int& ierr, double Norm = 1.) {
- if(_nrow!=_ncol) {
- std::cerr << "AlMat invertS: non-square matrix!" << std::endl;
- return;
+ if(m_nrow!=m_ncol) {
+ throw std::range_error( "AlMat invertS: non-square matrix!" );
}
- AlSymMat b(_nrow);
+ AlSymMat b(m_nrow);
for( int i=0; i<b.size(); i++ ) {
for( int j=0; j<=i; j++ ) {
b.elemr(i, j) = elemc(i, j);
@@ -311,23 +329,23 @@ void AlMat::invertS(int& ierr, double Norm = 1.) {
// reSize
void AlMat::reSize(int Nnew, int Mnew) {
- if ( Nnew != _nrow || Mnew != _ncol ) {
+ if ( Nnew != m_nrow || Mnew != m_ncol ) {
- double* p = ptr_data;
- _nele = Nnew*Mnew;
- ptr_data = new double[_nele];
- int _nrow_old = _nrow;
- int _ncol_old = _ncol;
- transpose = false;
+ double* p = m_ptr_data;
+ m_nele = Nnew*Mnew;
+ m_ptr_data = new double[m_nele];
+ int m_nrow_old = m_nrow;
+ int m_ncol_old = m_ncol;
+ m_transpose = false;
- _nrow = Nnew;
- _ncol = Mnew;
- int k = _nrow <= _nrow_old ? _nrow : _nrow_old;
- int l = _ncol <= _ncol_old ? _ncol : _ncol_old;
+ m_nrow = Nnew;
+ m_ncol = Mnew;
+ int k = m_nrow <= m_nrow_old ? m_nrow : m_nrow_old;
+ int l = m_ncol <= m_ncol_old ? m_ncol : m_ncol_old;
for( int i=0; i<k; i++ ) {
for( int j=0; j<l; j++ ) {
- *(ptr_data+i*_ncol+j) = *(p+i*_ncol_old+j);
+ *(m_ptr_data+i*m_ncol+j) = *(p+i*m_ncol_old+j);
}
}
@@ -341,25 +359,25 @@ StatusCode AlMat::ReadScalaPack(const std::string &filename){
if(inmat.fail())
return StatusCode::FAILURE;
- int32_t mNrow = _nrow;
+ int32_t mNrow = m_nrow;
inmat.read((char*)&mNrow, sizeof (mNrow));
- int32_t mNcol = _ncol;
+ int32_t mNcol = m_ncol;
inmat.read((char*)&mNcol, sizeof (mNcol));
- _nrow=abs(mNrow);
- _ncol=abs(mNcol);
- _nele=_ncol*_nrow;
- transpose = false;
+ m_nrow=abs(mNrow);
+ m_ncol=abs(mNcol);
+ m_nele=m_ncol*m_nrow;
+ m_transpose = false;
- // printf("ALMat::nrow: %d \n",_nrow);
- // printf("ALMat::ncol: %d \n",_ncol);
- // printf("ALMat::nele: %d \n",_nele);
+ // printf("ALMat::nrow: %d \n",m_nrow);
+ // printf("ALMat::ncol: %d \n",m_ncol);
+ // printf("ALMat::nele: %d \n",m_nele);
double melem=0;
- for(int i=0; i<_nrow; i++) {
- for(int j=0; j<_ncol; j++) {
+ for(int i=0; i<m_nrow; i++) {
+ for(int j=0; j<m_ncol; j++) {
inmat.read((char*)&melem, sizeof (melem));
- *(ptr_data+i*_ncol+j) = melem;
+ *(m_ptr_data+i*m_ncol+j) = melem;
// printf("(%d,%d) = %.16lf \n",i,j,melem);
}
}
@@ -390,9 +408,9 @@ StatusCode AlMat::Write(const std::string &filename, bool binary, unsigned int p
if(outmat.fail())
return StatusCode::FAILURE;
- int32_t mNrow = _nrow;
+ int32_t mNrow = m_nrow;
outmat.write((char*)&mNrow, sizeof (mNrow));
- int32_t mNcol = _ncol;
+ int32_t mNcol = m_ncol;
outmat.write((char*)&mNcol, sizeof (mNcol));
}
else {
@@ -407,9 +425,9 @@ StatusCode AlMat::Write(const std::string &filename, bool binary, unsigned int p
double melem=0;
- for( int i=0; i<_nrow; i++) {
- for( int j=0; j<_ncol; j++) {
- melem = *(ptr_data+i*_ncol+j);
+ for( int i=0; i<m_nrow; i++) {
+ for( int j=0; j<m_ncol; j++) {
+ melem = *(m_ptr_data+i*m_ncol+j);
if(binary)
outmat.write((char*)&(melem), sizeof (melem));
else
@@ -455,7 +473,7 @@ const std::string AlMat::Print(const int NColsPerSet)
for (int row=0; row<nrow(); row++) {
textmatrix << " |" << std::setw(6) << row << " | ";
for (int col=FirstCol; col<LastCol; col++) {
- double melem = *(ptr_data+row*ncol()+col);
+ double melem = *(m_ptr_data+row*ncol()+col);
textmatrix << std::setprecision(5) << std:: setw(6) << melem << " | ";
}
textmatrix << std::endl;
diff --git a/Tracking/TrkAlignment/TrkAlgebraUtils/src/AlSpaMat.cxx b/Tracking/TrkAlignment/TrkAlgebraUtils/src/AlSpaMat.cxx
index 84b7bde47b1bf796bdc28daddeff7ce742afd8b7..528a06926685e515c74805a424e967a162f25f11 100644
--- a/Tracking/TrkAlignment/TrkAlgebraUtils/src/AlSpaMat.cxx
+++ b/Tracking/TrkAlignment/TrkAlgebraUtils/src/AlSpaMat.cxx
@@ -19,88 +19,76 @@
#include <math.h>
#include <float.h>
#include <stdint.h>
+#include <stdexcept>
#include <TMatrixDSparse.h>
-extern"C" {
- void ma27id_(int ICNTL[],double CNTL[]);
- void ma27ad_(int *N,int *NZ,int IRN[],int ICN[],
- int IW[],int *LIW,int IKEEP[],int IW1p[],
- int *NSTEPS,int *IFLAG,int ICNTL[],double CNTL[],
- int INFO[],double *OPS);
+#include <Eigen/Core>
+#include <Eigen/IterativeLinearSolvers>
+#include <Eigen/SparseCholesky>
- void ma27bd_(int *N,int *NZ,int IRN[],int ICN[],
- double A[],int *LA,int IW[],int *LIW,int IKEEP[],
- int *NSTEPS,int *MAXFRT, int IW1[],int ICNTL[],
- double CNTL[], int INFO[]);
-
- void ma27cd_(int *N, double A[],int *LA,int IW[],int *LIW,
- double W[],int *MAXFRT,double RHS[],int IW1[],int *NSTEPS,
- int ICNTL[], int INFO[]);
-}
namespace Trk {
//______________________________________________________________________________
AlSpaMat::AlSpaMat()
{
- _matrix_type = 2;
- _size = 0;
- _nele = 0;
- ptr_row = NULL;
- ptr_col = NULL; // set pointer to null
+ m_matrix_type = 2;
+ m_size = 0;
+ m_nele = 0;
+ m_ptr_row = nullptr;
+ m_ptr_col = nullptr; // set pointer to null
}
//______________________________________________________________________________
AlSpaMat::AlSpaMat(long int N)
{
- _matrix_type = 2;
- _size = N;
- _nele = 0;
- ptr_row = NULL;
- ptr_col = NULL; // set pointer to null
+ m_matrix_type = 2;
+ m_size = N;
+ m_nele = 0;
+ m_ptr_row = nullptr;
+ m_ptr_col = nullptr; // set pointer to null
}
//______________________________________________________________________________
AlSpaMat::AlSpaMat(const AlSpaMat& m)
: AlSymMatBase(m)
{
- _matrix_type = 2;
- _size = m.size();
- ptr_row = NULL;
- ptr_col = NULL; // set pointer to null
+ m_matrix_type = 2;
+ m_size = m.size();
+ m_ptr_row = nullptr;
+ m_ptr_col = nullptr; // set pointer to null
copy(m);
}
//______________________________________________________________________________
AlSpaMat::AlSpaMat(const AlSymMat& m)
{
- _matrix_type = 2;
- _size = m.size();
- ptr_row = NULL;
- ptr_col = NULL; // set pointer to null
+ m_matrix_type = 2;
+ m_size = m.size();
+ m_ptr_row = nullptr;
+ m_ptr_col = nullptr; // set pointer to null
copy(m);
}
//______________________________________________________________________________
AlSpaMat::~AlSpaMat()
{
- ptr_map.clear();
- if( ptr_row != NULL ) delete [] ptr_row;
- if( ptr_col != NULL ) delete [] ptr_col;
+ m_ptr_map.clear();
+ if( m_ptr_row != nullptr ) delete [] m_ptr_row;
+ if( m_ptr_col != nullptr ) delete [] m_ptr_col;
}
//______________________________________________________________________________
void AlSpaMat::copy(const AlSpaMat& m)
{
if( size() != m.size()) {
- std::cerr << "AlSpaMat::copy: size do not match!" << std::endl;
- return;
+ throw std::range_error( "AlSpaMat::copy: size do not match!" );
}
- ptr_map.clear();
- _nele=m._nele;
- ptr_map = m.ptr_map;
+ m_ptr_map.clear();
+ m_nele=m.m_nele;
+ m_ptr_map = m.m_ptr_map;
return;
}
@@ -108,18 +96,17 @@ void AlSpaMat::copy(const AlSpaMat& m)
void AlSpaMat::copy(const AlSymMat& m)
{
if( size() != m.size()) {
- std::cerr << "AlSpaMat::copy: size do not match!" << std::endl;
- return;
+ throw std::range_error( "AlSpaMat::copy: size do not match!" );
}
- ptr_map.clear();
- _nele=0;
+ m_ptr_map.clear();
+ m_nele=0;
// Convert Storage System
- for (int i = 0; i<_size; i++)
+ for (int i = 0; i<m_size; i++)
for (int j = 0; j<=i; j++)
if (m.elemc(i,j) != 0.) {
- ptr_map.insert(std::make_pair(std::make_pair(i,j), m.elemc(i,j)));
- _nele++;
+ m_ptr_map.insert(std::make_pair(std::make_pair(i,j), m.elemc(i,j)));
+ m_nele++;
}
return;
@@ -129,20 +116,19 @@ void AlSpaMat::copy(const AlSymMat& m)
void AlSpaMat::copy(const AlMat& m)
{
if( size() != m.nrow() || size() != m.ncol() ) {
- std::cerr << "AlSpaMat::copy: size do not match!" << std::endl;
- return;
+ throw std::range_error( "AlSpaMat::copy: size do not match!" );
}
// copy just the lower triangle:
- ptr_map.clear();
- _nele=0;
+ m_ptr_map.clear();
+ m_nele=0;
//Convert Storage System
- for (int i = 0; i<_size; i++)
+ for (int i = 0; i<m_size; i++)
for (int j = 0; j<=i; j++)
if (m.elemc(i,j) != 0.) {
- ptr_map.insert(std::make_pair(std::make_pair(i,j), m.elemc(i,j)));
- _nele++;
+ m_ptr_map.insert(std::make_pair(std::make_pair(i,j), m.elemc(i,j)));
+ m_nele++;
}
return;
@@ -153,38 +139,34 @@ double& AlSpaMat::elemr(long int i,long int j)
{
#ifdef _DEBUG
if( i<0 ) {
- std::cerr << "AlSpaMat::elemr: Index 1 < zero! " << i << std::endl;
- return ptr_map.begin()->second;
+ throw std::out_of_range( "AlSpaMat::elemr: Index 1 < zero! " );
}
if( i>=size() ) {
- std::cerr << "AlSpaMat::elemr: Index 1 too large! " << i << std::endl;
- return ptr_map.begin()->second;
+ throw std::out_of_range( "AlSpaMat::elemr: Index 1 too large! " );
}
if( j<0 ) {
- std::cerr << "AlSpaMat::elemr: Index 2 < zero! " << j << std::endl;
- return ptr_map.begin()->second;
+ throw std::out_of_range( "AlSpaMat::elemr: Index 2 < zero! " );
}
if( j>=size() ) {
- std::cerr << "AlSpaMat::elemr: Index 2 too large! " << j << std::endl;
- return ptr_map.begin()->second;
+ throw std::out_of_range("AlSpaMat::elemr: Index 2 too large! " );;
}
#endif
// try fast referencing:
mapiterator pos;
indices key = j < i ? std::make_pair(i,j) : std::make_pair(j,i);
- pos=ptr_map.find(key);
- if( pos!=ptr_map.end() ) // already exists
+ pos=m_ptr_map.find(key);
+ if( pos!=m_ptr_map.end() ) // already exists
return pos->second;
else { // does not yet exist
// create it with the null value:
- _nele++;
- return (ptr_map.insert(std::make_pair(key, 0.))).first->second;
+ m_nele++;
+ return (m_ptr_map.insert(std::make_pair(key, 0.))).first->second;
}
// simpler implementation of the above:
/*
- ptr_map[key]=0;
- _nele = ptr_map.size();
- return ptr_map[key];
+ m_ptr_map[key]=0;
+ m_nele = m_ptr_map.size();
+ return m_ptr_map[key];
*/
}
@@ -193,28 +175,24 @@ double AlSpaMat::elemc(long int i,long int j) const
{
#ifdef _DEBUG
if( i<0 ) {
- std::cerr << "AlSpaMat::elemc: Index 1 < zero! " << i << std::endl;
- return 0.0;
+ throw std::out_of_range( "AlSpaMat::elemc: Index 1 < zero! " );
}
if( i>=size() ) {
- std::cerr << "AlSpaMat::elemc: Index 1 too large! " << i << std::endl;
- return 0.0;
+ throw std::out_of_range( "AlSpaMat::elemc: Index 1 too large! " );
}
if( j<0 ) {
- std::cerr << "AlSpaMat::elemc: Index 2 < zero! " << j << std::endl;
- return 0.0;
+ throw std::out_of_range( "AlSpaMat::elemc: Index 2 < zero! " );
}
if( j>=size() ) {
- std::cerr << "AlSpaMat::elemc: Index 2 too large! " << j << std::endl;
- return 0.0;
+ throw std::out_of_range( "AlSpaMat::elemc: Index 2 too large! " );
}
#endif
// try fast referencing:
const_mapiterator pos;
indices key = j < i ? std::make_pair(i,j) : std::make_pair(j,i);
- pos=ptr_map.find(key);
+ pos=m_ptr_map.find(key);
- if( pos!=ptr_map.end() ) // already exists
+ if( pos!=m_ptr_map.end() ) // already exists
return pos->second;
// does not yet exist
@@ -228,20 +206,16 @@ indices AlSpaMat::elem(long int i,long int j) const
// ATTENTION! the key value is returned:
#ifdef _DEBUG
if( i<0 ) {
- std::cerr << "AlSpaMat::elem: Index 1 < zero! " << i << std::endl;
- return std::make_pair(0,0);
+ throw std::out_of_range( "AlSpaMat::elem: Index 1 < zero! " );
}
if( i>=size() ) {
- std::cerr << "AlSpaMat::elem: Index 1 too large! " << i << std::endl;
- return std::make_pair(0,0);
+ throw std::out_of_range( "AlSpaMat::elem: Index 1 too large! " );
}
if( j<0 ) {
- std::cerr << "AlSpaMat::elem: Index 2 < zero! " << j << std::endl;
- return std::make_pair(0,0);
+ throw std::out_of_range( "AlSpaMat::elem: Index 2 < zero! " );
}
if( j>=size() ) {
- std::cerr << "AlSpaMat::elem: Index 2 too large! " << j << std::endl;
- return std::make_pair(0,0);
+ throw std::out_of_range( "AlSpaMat::elem: Index 2 too large! " );
}
#endif
@@ -254,7 +228,7 @@ AlSpaMat& AlSpaMat::operator=(const AlSpaMat& m)
{
if ( this != &m ) {
- _size=m.size();
+ m_size=m.size();
copy(m);
}
return *this;
@@ -263,7 +237,7 @@ AlSpaMat& AlSpaMat::operator=(const AlSpaMat& m)
//______________________________________________________________________________
AlSpaMat& AlSpaMat::operator=(const AlSymMat& m)
{
- _size=m.size();
+ m_size=m.size();
copy(m);
return *this;
}
@@ -272,11 +246,10 @@ AlSpaMat& AlSpaMat::operator=(const AlSymMat& m)
AlSpaMat& AlSpaMat::operator=(const AlMat& m)
{
if( m.nrow() != m.ncol() ) {
- std::cerr << "AlSpaMat::=operator: allowed for square matrices only!" << std::endl;
- return *this;
+ throw std::range_error( "AlSpaMat::=operator: allowed for square matrices only!" );
}
- _size=m.nrow();
+ m_size=m.nrow();
copy(m);
return *this;
}
@@ -285,7 +258,7 @@ AlSpaMat& AlSpaMat::operator=(const AlMat& m)
AlSpaMat& AlSpaMat::operator=(const double& d)
{
mapiterator pos;
- for (pos = ptr_map.begin(); pos!=ptr_map.end(); pos++)
+ for (pos = m_ptr_map.begin(); pos!=m_ptr_map.end(); pos++)
pos->second = d;
return *this;
@@ -295,15 +268,14 @@ AlSpaMat& AlSpaMat::operator=(const double& d)
AlSpaMat AlSpaMat::operator+(const AlSpaMat& m) const
{
if( size() != m.size()) {
- std::cerr << "AlSpaMat::operator+: size do not match!" << std::endl;
- return *this;
+ throw std::range_error( "AlSpaMat::operator+: size do not match!" );
}
AlSpaMat b(m);
const_mapiterator pos;
- for (pos = ptr_map.begin(); pos!=ptr_map.end(); pos++)
- b.ptr_map[pos->first] += pos->second;
- b._nele = b.ptr_map.size();
+ for (pos = m_ptr_map.begin(); pos!=m_ptr_map.end(); pos++)
+ b.m_ptr_map[pos->first] += pos->second;
+ b.m_nele = b.m_ptr_map.size();
return b;
}
@@ -312,14 +284,13 @@ AlSpaMat AlSpaMat::operator+(const AlSpaMat& m) const
AlSpaMat& AlSpaMat::operator+=(const AlSpaMat& m)
{
if( size() != m.size()) {
- std::cerr << "AlSpaMat::operator+=: size do not match!" << std::endl;
- return *this;
+ throw std::range_error( "AlSpaMat::operator+=: size do not match!" );
}
const_mapiterator pos;
- for (pos = m.ptr_map.begin(); pos!=m.ptr_map.end(); pos++)
- (*this).ptr_map[pos->first] += pos->second;
- _nele = ptr_map.size();
+ for (pos = m.m_ptr_map.begin(); pos!=m.m_ptr_map.end(); pos++)
+ (*this).m_ptr_map[pos->first] += pos->second;
+ m_nele = m_ptr_map.size();
return *this;
}
@@ -328,15 +299,14 @@ AlSpaMat& AlSpaMat::operator+=(const AlSpaMat& m)
AlSpaMat AlSpaMat::operator-(const AlSpaMat& m) const
{
if( size() != m.size()) {
- std::cerr << "AlSpaMat::operator-: size do not match!" << std::endl;
- return *this;
+ throw std::range_error( "AlSpaMat::operator-: size do not match!" );
}
AlSpaMat b(m);
const_mapiterator pos;
- for (pos = ptr_map.begin(); pos!=ptr_map.end(); pos++)
- b.ptr_map[pos->first] -= pos->second;
- b._nele = b.ptr_map.size();
+ for (pos = m_ptr_map.begin(); pos!=m_ptr_map.end(); pos++)
+ b.m_ptr_map[pos->first] -= pos->second;
+ b.m_nele = b.m_ptr_map.size();
return b;
}
@@ -345,14 +315,13 @@ AlSpaMat AlSpaMat::operator-(const AlSpaMat& m) const
AlSpaMat& AlSpaMat::operator-=(const AlSpaMat& m)
{
if( size() != m.size()) {
- std::cerr << "AlSpaMat::operator-=: size do not match!" << std::endl;
- return *this;
+ throw std::range_error( "AlSpaMat::operator-=: size do not match!" );
}
const_mapiterator pos;
- for (pos = m.ptr_map.begin(); pos!=m.ptr_map.end(); pos++)
- (*this).ptr_map[pos->first] -= pos->second;
- _nele = ptr_map.size();
+ for (pos = m.m_ptr_map.begin(); pos!=m.m_ptr_map.end(); pos++)
+ (*this).m_ptr_map[pos->first] -= pos->second;
+ m_nele = m_ptr_map.size();
return *this;
}
@@ -361,9 +330,7 @@ AlSpaMat& AlSpaMat::operator-=(const AlSpaMat& m)
AlMat AlSpaMat::operator*(const AlSymMatBase& m) const
{
if( size() != m.size() ) {
- std::cerr << "AlSpaMat::operator*: size do not match!" << std::endl;
- AlMat b( size(), m.size());
- return b;
+ throw std::range_error( "AlSpaMat::operator*: size do not match!" );
}
long int isiz(size());
@@ -380,8 +347,7 @@ AlMat AlSpaMat::operator*(const AlSymMatBase& m) const
//______________________________________________________________________________
AlMat AlSpaMat::operator*(const AlMat& m) const {
if( size() != m.nrow() ) {
- std::cerr << "AlSpaMat::operator*: size do not match!" << std::endl;
- return m;
+ throw std::range_error( "AlSpaMat::operator*: size do not match!" );
}
long int isiz(size());
@@ -400,8 +366,7 @@ AlMat AlSpaMat::operator*(const AlMat& m) const {
AlVec AlSpaMat::operator*(const AlVec& v) const
{
if( size() != v.size() ) {
- std::cerr << "AlSpaMat::operator*: size do not match! " << std::endl;
- return v;
+ throw std::range_error( "AlSpaMat::operator*: size do not match! " );
}
long int isiz(size());
@@ -418,7 +383,7 @@ AlVec AlSpaMat::operator*(const AlVec& v) const
AlSpaMat& AlSpaMat::operator*=(const double& d)
{
mapiterator pos;
- for (pos = ptr_map.begin(); pos!=ptr_map.end(); pos++)
+ for (pos = m_ptr_map.begin(); pos!=m_ptr_map.end(); pos++)
pos->second *= d;
return *this;
@@ -429,8 +394,8 @@ AlSpaMat AlSpaMat::operator*(const double& d) const
{
AlSpaMat a(size());
const_mapiterator pos;
- for (pos = ptr_map.begin(); pos!=ptr_map.end(); pos++)
- a.ptr_map.insert(std::make_pair(pos->first, (pos->second)*d));
+ for (pos = m_ptr_map.begin(); pos!=m_ptr_map.end(); pos++)
+ a.m_ptr_map.insert(std::make_pair(pos->first, (pos->second)*d));
return a;
}
@@ -438,118 +403,94 @@ AlSpaMat AlSpaMat::operator*(const double& d) const
//______________________________________________________________________________
double AlSpaMat::determinant()
{
- double deter = 1.;
-
- std::cerr << "AlSpaMat::determinant: not implemented!" << std::endl;
-
- return deter;
+ throw std::invalid_argument( "AlSpaMat::determinant: not implemented!" );
}
-//______________________________________________________________________________
-int AlSpaMat::***REMOVED***Solve(AlVec& RHS)
-{
- int ICNTL[30];
- double CNTL[5];
- ma27id_(ICNTL,CNTL);
-
- double* ptr_data = new double[_nele];
- ptr_row = new int[_nele];
- ptr_col = new int[_nele];
-
- int Size(_size);
- int Nele(_nele);
-
- //Convert Storage System
+int AlSpaMat::SolveWithEigen(AlVec& RHS){
+
+ if(RHS.size() != size() ){
+ throw std::range_error( "AlSpaMat::SolveWithEigen vector size is incorrect" );
+ }
+
+ Eigen::VectorXd eigenBigVector( RHS.size() );
+ for(int i(0); i< RHS.size(); ++i ){
+ eigenBigVector[i] = RHS[i];
+ }
+ Eigen::SparseMatrix<double> eigenBigMatrix( m_size, m_size );
+
+ typedef Eigen::Triplet<double> Triplet;
+ std::vector<Triplet> tripletList;
+ tripletList.reserve(m_nele);
long int i, j;
long int counter(0);
mapiterator pos;
- for (pos = ptr_map.begin(); pos!=ptr_map.end(); pos++){
+ for (pos = m_ptr_map.begin(); pos!=m_ptr_map.end(); pos++){
elem(pos->first, i, j);
- *(ptr_data+counter)=pos->second;
- *(ptr_row+counter)= i+1;
- *(ptr_col+counter)= j+1;
+ tripletList.push_back(Triplet(i,j,pos->second));
+ if(i!=j) tripletList.push_back(Triplet(j,i,pos->second));
counter++;
}
+ eigenBigMatrix.setFromTriplets(tripletList.begin(), tripletList.end());
- int LIW =3/2*(2*Nele+3*Size+1);
- int* IW = new int[LIW];
- int* IKEEP = new int[3*Size];
- int* IW1 = new int[2*Size];
- int NSTEPS;
- int IFLAG = 0;
- int INFO[20];
- double OPS;
-
- ma27ad_(&Size, &Nele,ptr_row,ptr_col,
- IW,&LIW,IKEEP,IW1,
- &NSTEPS,&IFLAG,ICNTL,CNTL,
- INFO,&OPS);
-
- int MAXFRT;
- int LA =2*INFO[4];
+ // Eigen::CholmodSupernodalLLT is much much quicker (x50) so it would be great to move to that in the future
+ // requires an external package SuiteSparse
+ // BiCGSTAB was the fastest iterative solver in Eigen from a quick test
+ // SimplicialLDLT was the fastest direct solver in Eigen (x2 slower than BiCGSTAB )
- double* TempA = new double[Nele];
- for (int i=0; i<Nele ;i++)
- *(TempA+i) = *(ptr_data+i);
+ // Eigen::BiCGSTAB<Eigen::SparseMatrix<double> > solver;
+ Eigen::SimplicialLDLT<Eigen::SparseMatrix<double> > solver;
- delete [] ptr_data;
- ptr_data = new double[LA];
-
- for (int i=0; i<LA ;i++)
- *(ptr_data+i)=0;
-
- for (int i=0; i<Nele ;i++)
- *(ptr_data+i)= *(TempA+i);
-
- delete [] TempA;
-
- ma27bd_(&Size,&Nele,ptr_row,ptr_col,
- ptr_data,&LA,IW,&LIW,IKEEP,
- &NSTEPS,&MAXFRT,IW1,ICNTL,
- CNTL,INFO);
+ solver.compute(eigenBigMatrix);
+ if(solver.info()!=Eigen::Success) {
+ // decomposition failed
+ throw std::domain_error("AlSpaMat::SolveWithEigen: failed to compute: -- is the input matrix singular?" );
+ return 1;
+ }
- double* W =new double[MAXFRT];
+ Eigen::VectorXd x = solver.solve( eigenBigVector );
+ if(solver.info()!=Eigen::Success) {
+ // solving failed
+ throw std::domain_error("AlSpaMat::SolveWithEigen: Failed to solve: -- is your matrix singular? ");
+ return 2;
+ }
+
+ //Copy results into vector
+ for(int i(0); i< RHS.size(); ++i ){
+ RHS[i] = x[i];
+ }
- double* RHStemp =new double[Size];
- for (int i=0; i<Size; i++)
- RHStemp[i] = RHS[i];
+ //Check accuracy
+ Eigen::VectorXd residual = eigenBigMatrix * x - eigenBigVector;
+ double sumresidual = 0;
+ for( int i=0; i<residual.size(); ++i){
+ sumresidual += fabs(residual[i]);
+ }
+ sumresidual /= (double) residual.size();
+
+ if( sumresidual > 1e-3 ){
+ throw std::overflow_error( "AlSpaMat::SolveWithEigen: your solution is no good! ");
+ return 3;
+ }
- ma27cd_(&Size, ptr_data,&LA,IW,&LIW,
- W,&MAXFRT,RHStemp,IW1,&NSTEPS,
- ICNTL,INFO);
+ return 0;
- for (int i=0; i<Size; i++)
- RHS[i] = RHStemp[i];
- delete [] RHStemp;
- delete [] IW;
- delete [] IKEEP;
- delete [] IW1;
- delete [] W;
+}
- delete [] ptr_data; ptr_data=NULL;
- delete [] ptr_row; ptr_row=NULL;
- delete [] ptr_col; ptr_col=NULL;
- return INFO[0];
-}
//______________________________________________________________________________
//jlove int AlSpaMat::diagonalize(char jobz, AlVec& w, AlMat& z) {
int AlSpaMat::diagonalize(char, AlVec&, AlMat&)
{
- std::cerr << "AlSpaMat::diagonalize: not implemented!" << std::endl;
- int ierr = -1;
- return ierr;
+ throw std::invalid_argument( "AlSpaMat::diagonalize: not implemented!" );
}
//______________________________________________________________________________
int AlSpaMat::invert()
{
- std::cerr << "AlSpaMat::invert: not implemented!" << std::endl;
- int ierr = -1;
- return ierr;
+ throw std::invalid_argument( "AlSpaMat::invert: not implemented!" );
}
//______________________________________________________________________________
@@ -557,16 +498,16 @@ void AlSpaMat::RemoveDoF(int index, int nelem)
{
// index here is the first alignment parameter to remove, nelem number of consecutive ones
- _size-=nelem;
+ m_size-=nelem;
long int n=index+nelem-1;
mapiterator pos;
- mapiterator pos_obs=ptr_map.end();
- for (pos = ptr_map.begin(); pos!=ptr_map.end(); pos++) {
+ mapiterator pos_obs=m_ptr_map.end();
+ for (pos = m_ptr_map.begin(); pos!=m_ptr_map.end(); pos++) {
- if( pos_obs!=ptr_map.end() )
- ptr_map.erase(pos_obs);
+ if( pos_obs!=m_ptr_map.end() )
+ m_ptr_map.erase(pos_obs);
- pos_obs = ptr_map.end();
+ pos_obs = m_ptr_map.end();
long int k, l;
elem(pos->first, k, l);
@@ -585,9 +526,9 @@ void AlSpaMat::RemoveDoF(int index, int nelem)
}
}
- if( pos_obs!=ptr_map.end() )
- ptr_map.erase(pos_obs);
- _nele = ptr_map.size();
+ if( pos_obs!=m_ptr_map.end() )
+ m_ptr_map.erase(pos_obs);
+ m_nele = m_ptr_map.size();
}
//______________________________________________________________________________
@@ -595,12 +536,11 @@ int AlSpaMat::RemoveCollsRows(std::vector<int> indices)
{
int n = indices.size();
if (n==0) {
- std::cerr<<"Vector of indices to remove is empty."<<std::endl;
- return _size;
+ return m_size;
}
- if (n>_size) {
- std::cerr<<"Vector of indices larger than matrix size."<<std::endl;
- return _size;
+ if (n>m_size) {
+ throw std::invalid_argument( "AlSpaMat::RemoveCollsRows: Vector of indices larger than matrix size." );
+ return m_size;
}
// first sort the list of indices descending
@@ -615,14 +555,13 @@ int AlSpaMat::RemoveCollsRows(std::vector<int> indices)
// remove rows and columns starting from largest indices
for (int i=0;i<n;i++) {
- if (indices[i] > _size-1) {
- std::cerr<<"Index "<<indices[i]<<" goes beyond matrix (size "<<_size<<")."<<std::endl;
- continue;
+ if (indices[i] > m_size-1) {
+ throw std::out_of_range( "AlSpaMat::RemoveCollsRows: Index goes beyond matrix." );
}
RemoveDoF(indices[i]);
}
- return _size;
+ return m_size;
}
//______________________________________________________________________________
@@ -652,15 +591,15 @@ void AlSpaMat::RemoveAlignPar(int index, int control)
counterCol=0;
}
- pos = ptr_map.find(std::make_pair(row+shiftRow,col+shiftCol));
- if( pos!=ptr_map.end() ) { // exists
+ pos = m_ptr_map.find(std::make_pair(row+shiftRow,col+shiftCol));
+ if( pos!=m_ptr_map.end() ) { // exists
elemr(row,col) = pos->second;
- ptr_map.erase(pos);
+ m_ptr_map.erase(pos);
}
else {
- pos = ptr_map.find(std::make_pair(row,col));
- if( pos!=ptr_map.end() )
- ptr_map.erase(pos);
+ pos = m_ptr_map.find(std::make_pair(row,col));
+ if( pos!=m_ptr_map.end() )
+ m_ptr_map.erase(pos);
}
counterCol++;
if (counterCol==5-control) {
@@ -676,7 +615,7 @@ void AlSpaMat::RemoveAlignPar(int index, int control)
}
}
- _nele = ptr_map.size();
+ m_nele = m_ptr_map.size();
}
//______________________________________________________________________________
@@ -694,21 +633,21 @@ void AlSpaMat::RemoveModule(int index)
void AlSpaMat::reSize(long int n)
{
// balanced binary tree!
- if( n!=_size ) {
+ if( n!=m_size ) {
// make a temporary copy:
AlSpaMat m(*this);
- ptr_map.clear();
- _size = n;
+ m_ptr_map.clear();
+ m_size = n;
long int i, j;
mapiterator pos;
- for (pos = m.ptr_map.begin(); pos!=m.ptr_map.end(); pos++) {
+ for (pos = m.m_ptr_map.begin(); pos!=m.m_ptr_map.end(); pos++) {
m.elem(pos->first, i, j);
if( i<n && j<n )
- ptr_map.insert(*pos);
+ m_ptr_map.insert(*pos);
}
}
- _nele = ptr_map.size();
+ m_nele = m_ptr_map.size();
return;
}
@@ -728,12 +667,12 @@ void AlSpaMat::SetPathTxt(const std::string &path)
//jlove StatusCode AlSpaMat::Write(const std::string &filename, bool binary,
// bool square, double scale, float version){
StatusCode AlSpaMat::Write(const std::string &filename, bool binary,
- bool square, double, float version)
+ bool /*square*/, double, float version)
{
std::ofstream outmat;
int32_t msizz = 1000000+size();
- int32_t nelem = ptr_map.size();
+ int32_t nelem = m_ptr_map.size();
if(binary) {
outmat.open((m_pathbin+filename).c_str(), std::ios::binary);
@@ -744,9 +683,6 @@ StatusCode AlSpaMat::Write(const std::string &filename, bool binary,
outmat.write((char*)&msizz, sizeof (msizz));
outmat.write((char*)&version, sizeof (version));
outmat.write((char*)&nelem, sizeof (nelem));
- // std::cout << "AlSpaMat::Write: msizz = " << msizz << std::endl;
- // std::cout << "AlSpaMat::Write: version = " << version << std::endl;
- // std::cout << "AlSpaMat::Write: nelem = " << nelem << std::endl;
}
else {
outmat.open((m_pathtxt+filename).c_str());
@@ -765,7 +701,7 @@ StatusCode AlSpaMat::Write(const std::string &filename, bool binary,
int32_t i, j;
mapiterator pos;
- for (pos = ptr_map.begin(); pos!=ptr_map.end(); pos++) {
+ for (pos = m_ptr_map.begin(); pos!=m_ptr_map.end(); pos++) {
melem = pos->second;
elem(pos->first, ii, jj); i=ii; j=jj; // just a type conversion
if(binary) {
@@ -777,8 +713,6 @@ StatusCode AlSpaMat::Write(const std::string &filename, bool binary,
outmat << std::setw(6) << i << std::setw(6) << j << std::setw(18) << melem << std::endl;
}
outmat.close();
- if(square)
- std::cout << "AlSpaMat::Write: square flag has been ignored!" << std::endl;
return StatusCode::SUCCESS;
}
@@ -790,7 +724,7 @@ StatusCode AlSpaMat::CheckMatVersion(const std::string filename, bool &StdUnits)
if(inmat.fail())
return StatusCode::FAILURE;
- ptr_map.clear();
+ m_ptr_map.clear();
int32_t msiz=0;
inmat.read((char*)&msiz, sizeof (msiz));
@@ -805,8 +739,6 @@ StatusCode AlSpaMat::CheckMatVersion(const std::string filename, bool &StdUnits)
inmat.close();
- // std::cout << "AlSpaMat::StdUnits: " << StdUnits << std::endl;
-
return StatusCode::SUCCESS;
}
@@ -817,13 +749,11 @@ StatusCode AlSpaMat::Read(const std::string &filename, int &dofs, bool &triang,
if (StatusCode::SUCCESS != CheckMatVersion(m_pathbin+filename, m_StdUnits))
return StatusCode::FAILURE;
- // std::cout << "AlSpaMat::StdUnits: " << m_StdUnits << std::endl;
-
std::ifstream inmat((m_pathbin+filename).c_str(), std::ios::binary);
if(inmat.fail())
return StatusCode::FAILURE;
- ptr_map.clear();
+ m_ptr_map.clear();
int32_t msiz=0;
int32_t nelem;
@@ -832,11 +762,10 @@ StatusCode AlSpaMat::Read(const std::string &filename, int &dofs, bool &triang,
dofs = msiz-1000000;
else
dofs = abs(msiz);
- _size = dofs;
+ m_size = dofs;
if (m_StdUnits)
inmat.read((char*)&version, sizeof (version));
- // std::cout << "AlSpaMat::Write: version = " << version << std::endl;
double melem=0;
int32_t i, j;
@@ -844,13 +773,12 @@ StatusCode AlSpaMat::Read(const std::string &filename, int &dofs, bool &triang,
if(msiz>999999) { // sparse format
triang=false;
inmat.read((char*)&nelem, sizeof (nelem));
- _nele=nelem;
- // std::cout << "AlSpaMat::Write: nelem = " << nelem << std::endl;
+ m_nele=nelem;
for(int k=0; k<nelem; k++) {
inmat.read((char*)&i, sizeof (i));
inmat.read((char*)&j, sizeof (j));
inmat.read((char*)&melem, sizeof (melem));
- ptr_map.insert(std::make_pair(std::make_pair(i,j), melem));
+ m_ptr_map.insert(std::make_pair(std::make_pair(i,j), melem));
}
}
else if(msiz>0) { // square format
@@ -858,26 +786,24 @@ StatusCode AlSpaMat::Read(const std::string &filename, int &dofs, bool &triang,
for(int32_t i=0; i<msiz; i++) {
for(int32_t j=0; j<msiz; j++) {
inmat.read((char*)&melem, sizeof (melem));
- // std::cout << "AlSpaMat::Write: nelem = " << nelem << std::endl;
if( i>=j && melem!=0. )
- ptr_map.insert(std::make_pair(std::make_pair(i,j), melem));
+ m_ptr_map.insert(std::make_pair(std::make_pair(i,j), melem));
}
}
}
else { // triangular format
triang=true;
msiz = (-1)*msiz;
- // std::cout << "msiz=" << msiz << std::endl;
for( int32_t i=0; i<msiz; i++) {
for( int32_t j=0; j<=i; j++) {
inmat.read((char*)&melem, sizeof (melem));
if( melem!=0. )
- ptr_map.insert(std::make_pair(std::make_pair(i,j), melem));
+ m_ptr_map.insert(std::make_pair(std::make_pair(i,j), melem));
}
}
}
- _nele = ptr_map.size();
+ m_nele = m_ptr_map.size();
inmat.close();
return StatusCode::SUCCESS;
}
@@ -890,7 +816,7 @@ StatusCode AlSpaMat::ReadProjected(const std::string &filename, int &dofs,
if(inmat.fail())
return StatusCode::FAILURE;
- ptr_map.clear();
+ m_ptr_map.clear();
int32_t msiz=0;
int32_t nelem;
@@ -899,7 +825,7 @@ StatusCode AlSpaMat::ReadProjected(const std::string &filename, int &dofs,
dofs = msiz-1000000;
else
dofs = abs(msiz);
- _size = dofs;
+ m_size = dofs;
inmat.read((char*)&version, sizeof (version));
@@ -909,13 +835,12 @@ StatusCode AlSpaMat::ReadProjected(const std::string &filename, int &dofs,
if(msiz>999999) { // sparse format
triang=false;
inmat.read((char*)&nelem, sizeof (nelem));
- _nele=nelem;
-// std::cout << "AlSpaMat::Write: nelem = " << nelem << std::endl;
+ m_nele=nelem;
for(int k=0; k<nelem; k++) {
inmat.read((char*)&i, sizeof (i));
inmat.read((char*)&j, sizeof (j));
inmat.read((char*)&melem, sizeof (melem));
- ptr_map.insert(std::make_pair(std::make_pair(i,j), melem));
+ m_ptr_map.insert(std::make_pair(std::make_pair(i,j), melem));
}
}
else if(msiz>0) { // square format
@@ -924,24 +849,23 @@ StatusCode AlSpaMat::ReadProjected(const std::string &filename, int &dofs,
for(int32_t j=0; j<msiz; j++) {
inmat.read((char*)&melem, sizeof (melem));
if( i>=j && melem!=0. )
- ptr_map.insert(std::make_pair(std::make_pair(i,j), melem));
+ m_ptr_map.insert(std::make_pair(std::make_pair(i,j), melem));
}
}
}
else { // triangular format
triang=true;
msiz = (-1)*msiz;
-// std::cout << "msiz=" << msiz << std::endl;
for( int32_t i=0; i<msiz; i++) {
for( int32_t j=0; j<=i; j++) {
inmat.read((char*)&melem, sizeof (melem));
if( melem!=0. )
- ptr_map.insert(std::make_pair(std::make_pair(i,j), melem));
+ m_ptr_map.insert(std::make_pair(std::make_pair(i,j), melem));
}
}
}
- _nele = ptr_map.size();
+ m_nele = m_ptr_map.size();
inmat.close();
return StatusCode::SUCCESS;
}
@@ -949,15 +873,15 @@ StatusCode AlSpaMat::ReadProjected(const std::string &filename, int &dofs,
TMatrixDSparse* AlSpaMat::makeTMatrix()
{
- long int nonZeroElements = ptr_map.size();
+ long int nonZeroElements = m_ptr_map.size();
int *irow = new int[nonZeroElements];
int *icol = new int[nonZeroElements];
double *val = new double[nonZeroElements];
long int i, j;
long int counter(0);
- const_mapiterator pos = ptr_map.begin();
- for(pos=ptr_map.begin(); pos!=ptr_map.end(); pos++){
+ const_mapiterator pos = m_ptr_map.begin();
+ for(pos=m_ptr_map.begin(); pos!=m_ptr_map.end(); pos++){
i = pos->first.first;
j = pos->first.second;
*(val+counter)=pos->second;
@@ -966,9 +890,7 @@ TMatrixDSparse* AlSpaMat::makeTMatrix()
counter++;
}
- std::cout << counter << " " << nonZeroElements << " size " << _size <<std::endl;
-
- TMatrixDSparse* myTMatrix = new TMatrixDSparse(0,_size-1,0,_size-1);
+ TMatrixDSparse* myTMatrix = new TMatrixDSparse(0,m_size-1,0,m_size-1);
myTMatrix->SetMatrixArray(nonZeroElements,irow,icol,val);
delete [] irow;
diff --git a/Tracking/TrkAlignment/TrkAlgebraUtils/src/AlSymMat.cxx b/Tracking/TrkAlignment/TrkAlgebraUtils/src/AlSymMat.cxx
index c31c41e2ea5372400ad5c25c40e18f131f46bdfa..b63210ce40c59c064edbfb918a6327d407cc922f 100644
--- a/Tracking/TrkAlignment/TrkAlgebraUtils/src/AlSymMat.cxx
+++ b/Tracking/TrkAlignment/TrkAlgebraUtils/src/AlSymMat.cxx
@@ -16,6 +16,7 @@
#include <math.h>
#include <float.h> //for DBL_EPSILON
#include <stdint.h>
+#include <stdexcept>
#include <TMatrixDSparse.h>
@@ -38,10 +39,10 @@ namespace Trk {
//______________________________________________________________________________
AlSymMat::AlSymMat()
{
- _matrix_type = 1;
- _size = 0;
- _nele = 0;
- ptr_data = NULL; // set pointer to null
+ m_matrix_type = 1;
+ m_size = 0;
+ m_nele = 0;
+ m_ptr_data = nullptr; // set pointer to null
m_pathbin="./";
m_pathtxt="./";
}
@@ -50,15 +51,15 @@ AlSymMat::AlSymMat()
//______________________________________________________________________________
AlSymMat::AlSymMat(long int N)
{
- _matrix_type = 1;
- _size = N;
- _nele = N*(N+1)/2;
- ptr_data = new double[_nele];
+ m_matrix_type = 1;
+ m_size = N;
+ m_nele = N*(N+1)/2;
+ m_ptr_data = new double[m_nele];
m_pathbin="./";
m_pathtxt="./";
- double* p = ptr_data + _nele;
- while (p > ptr_data) *(--p) = 0.;
+ double* p = m_ptr_data + m_nele;
+ while (p > m_ptr_data) *(--p) = 0.;
}
@@ -66,12 +67,12 @@ AlSymMat::AlSymMat(long int N)
AlSymMat::AlSymMat(const AlSymMat& m)
: AlSymMatBase(m)
{
- _matrix_type = 1;
- _size = m.size();
- _nele = m._nele;
+ m_matrix_type = 1;
+ m_size = m.size();
+ m_nele = m.m_nele;
m_pathbin = m.m_pathbin;
m_pathtxt = m.m_pathtxt;
- ptr_data = new double[_nele];
+ m_ptr_data = new double[m_nele];
copy(m);
}
@@ -79,25 +80,25 @@ AlSymMat::AlSymMat(const AlSymMat& m)
//______________________________________________________________________________
AlSymMat::AlSymMat(const AlSpaMat& m)
{
- _matrix_type = 1;
- _size = m.size();
- _nele = _size*(_size+1)/2;
+ m_matrix_type = 1;
+ m_size = m.size();
+ m_nele = m_size*(m_size+1)/2;
m_pathbin = m.pathBin();
m_pathtxt = m.pathTxt();
- ptr_data = new double[_nele];
+ m_ptr_data = new double[m_nele];
copy(m);
}
//______________________________________________________________________________
AlSymMat& AlSymMat::operator=(const AlSpaMat& m)
{
- _matrix_type = 1;
- _size = m.size();
- _nele = _size*(_size+1)/2;
+ m_matrix_type = 1;
+ m_size = m.size();
+ m_nele = m_size*(m_size+1)/2;
m_pathbin = m.pathBin();
m_pathtxt = m.pathTxt();
- ptr_data = new double[_nele];
+ m_ptr_data = new double[m_nele];
copy(m);
return *this;
}
@@ -105,7 +106,7 @@ AlSymMat& AlSymMat::operator=(const AlSpaMat& m)
//______________________________________________________________________________
AlSymMat::~AlSymMat()
{
- if( ptr_data != NULL ) delete [] ptr_data;
+ if( m_ptr_data != nullptr ) delete [] m_ptr_data;
//ptr_map.clear();
}
@@ -114,13 +115,12 @@ void AlSymMat::copy(const AlSymMat& m)
{
// this one implements the fast copy alg.
if( size() != m.size()) {
- std::cerr << "AlSymMat::copy: size do not match!" << std::endl;
- return;
+ throw std::range_error( "AlSymMat::copy: size do not match!" );
}
- double * p = ptr_data + _nele;
- double* r = m.ptr_data + _nele;
- while (p > ptr_data) *(--p) = *(--r);
+ double * p = m_ptr_data + m_nele;
+ double* r = m.m_ptr_data + m_nele;
+ while (p > m_ptr_data) *(--p) = *(--r);
return;
}
@@ -130,8 +130,7 @@ void AlSymMat::copy(const AlSpaMat& m)
{
//
if( size() != m.size()) {
- std::cerr << "AlSymMat::copy: size do not match!" << std::endl;
- return;
+ throw std::range_error( "AlSymMat::copy: size do not match!" );
}
(*this) = 0.;
@@ -152,8 +151,7 @@ void AlSymMat::copy(const AlMat& m)
// copy the lower triangle only!
int si = size();
if( si != m.nrow() || si != m.ncol() ) {
- std::cerr << "AlSymMat::copy: sizes do not match!" << std::endl;
- return;
+ throw std::range_error("AlSymMat::copy: sizes do not match!" );
}
for( int i=0; i<si; i++ )
@@ -167,7 +165,7 @@ AlSymMat& AlSymMat::operator=(const AlSymMat& m)
if (this==&m)
return *this;
- if ( ptr_data != m.ptr_data ) {
+ if ( m_ptr_data != m.m_ptr_data ) {
reSize(m.size());
copy(m);
}
@@ -178,8 +176,7 @@ AlSymMat& AlSymMat::operator=(const AlSymMat& m)
AlSymMat& AlSymMat::operator=(const AlMat& m)
{
if( m.nrow() != m.ncol() ) {
- std::cerr << "AlSymMat::operator=: a squared matrix is required!" << std::endl;
- return *this;
+ throw std::range_error( "AlSymMat::operator=: a squared matrix is required!" );
}
reSize(m.nrow());
@@ -191,8 +188,8 @@ AlSymMat& AlSymMat::operator=(const AlMat& m)
//______________________________________________________________________________
AlSymMat& AlSymMat::operator=(const double& d)
{
- double * p = ptr_data + _nele;
- while (p > ptr_data) *(--p) = d;
+ double * p = m_ptr_data + m_nele;
+ while (p > m_ptr_data) *(--p) = d;
return *this;
}
@@ -201,15 +198,14 @@ AlSymMat& AlSymMat::operator=(const double& d)
AlSymMat AlSymMat::operator+(const AlSymMat& m) const
{
if( size() != m.size()) {
- std::cerr << "AlSymMat::operator+: size do not match!" << std::endl;
- return *this;
+ throw std::range_error( "AlSymMat::operator+: size do not match!" );
}
AlSymMat b(size());
- double * p = ptr_data + _nele;
- double * q = m.ptr_data + _nele;
- double * r = b.ptr_data + _nele;
- while (p > ptr_data) *(--r) = (*(--p))+(*(--q));
+ double * p = m_ptr_data + m_nele;
+ double * q = m.m_ptr_data + m_nele;
+ double * r = b.m_ptr_data + m_nele;
+ while (p > m_ptr_data) *(--r) = (*(--p))+(*(--q));
return b;
}
@@ -218,13 +214,12 @@ AlSymMat AlSymMat::operator+(const AlSymMat& m) const
AlSymMat& AlSymMat::operator+=(const AlSymMat& m)
{
if( size() != m.size()){
- std::cerr << "AlSymMat::operator+=: size do not match!" << std::endl;
- return *this;
+ throw std::range_error( "AlSymMat::operator+=: size do not match!" );
}
- double * p = ptr_data + _nele;
- double * q = m.ptr_data + _nele;
- while (p > ptr_data) *(--p) += *(--q);
+ double * p = m_ptr_data + m_nele;
+ double * q = m.m_ptr_data + m_nele;
+ while (p > m_ptr_data) *(--p) += *(--q);
return *this;
}
@@ -233,15 +228,14 @@ AlSymMat& AlSymMat::operator+=(const AlSymMat& m)
AlSymMat AlSymMat::operator-(const AlSymMat& m) const
{
if( size() != m.size()) {
- std::cerr << "AlSymMat::operator-: size do not match!" << std::endl;
- return *this;
+ throw std::range_error( "AlSymMat::operator-: size do not match!" );
}
AlSymMat b(size());
- double * p = ptr_data + _nele;
- double * q = m.ptr_data + _nele;
- double * r = b.ptr_data + _nele;
- while (p > ptr_data) *(--r) = (*(--p))-(*(--q));
+ double * p = m_ptr_data + m_nele;
+ double * q = m.m_ptr_data + m_nele;
+ double * r = b.m_ptr_data + m_nele;
+ while (p > m_ptr_data) *(--r) = (*(--p))-(*(--q));
return b;
}
@@ -250,13 +244,12 @@ AlSymMat AlSymMat::operator-(const AlSymMat& m) const
AlSymMat& AlSymMat::operator-=(const AlSymMat& m)
{
if( size() != m.size()) {
- std::cerr << "AlSymMat::operator-=: size do not match!" << std::endl;
- return *this;
+ throw std::range_error( "AlSymMat::operator-=: size do not match!" );
}
- double * p = ptr_data + _nele;
- double * q = m.ptr_data + _nele;
- while (p > ptr_data) *(--p) -= *(--q);
+ double * p = m_ptr_data + m_nele;
+ double * q = m.m_ptr_data + m_nele;
+ while (p > m_ptr_data) *(--p) -= *(--q);
return *this;
}
@@ -265,9 +258,7 @@ AlSymMat& AlSymMat::operator-=(const AlSymMat& m)
AlMat AlSymMat::operator*(const AlSymMat& m) const
{
if( size() != m.size() ) {
- std::cerr << "AlSymMat::operator*: size do not match!" << std::endl;
- AlMat b( size(), m.size());
- return b;
+ throw std::range_error( "AlSymMat::operator*: size do not match!" );
}
AlMat b( size(), size());
@@ -284,17 +275,17 @@ AlMat AlSymMat::operator*(const AlSymMat& m) const
ii = (i+1)*i/2;
jj = (j+1)*j/2;
for( k=0; k<l; k++)
- x += (*(ptr_data+ii+k))*(*(m.ptr_data+jj+k));
+ x += (*(m_ptr_data+ii+k))*(*(m.m_ptr_data+jj+k));
if( i<j ) {
for( k=l; k<n; k++)
- x += (*(ptr_data+(k+1)*k/2+i))*(*(m.ptr_data+jj+k));
+ x += (*(m_ptr_data+(k+1)*k/2+i))*(*(m.m_ptr_data+jj+k));
}
else {
for( k=l; k<n; k++)
- x += (*(ptr_data+ii+k))*(*(m.ptr_data+(k+1)*k/2+j));
+ x += (*(m_ptr_data+ii+k))*(*(m.m_ptr_data+(k+1)*k/2+j));
}
for( k=n; k<siz; k++)
- x += (*(ptr_data+(k+1)*k/2+i))*(*(m.ptr_data+(k+1)*k/2+j));
+ x += (*(m_ptr_data+(k+1)*k/2+i))*(*(m.m_ptr_data+(k+1)*k/2+j));
b.elemr(i,j) = x;
}
@@ -305,8 +296,7 @@ AlMat AlSymMat::operator*(const AlSymMat& m) const
AlMat AlSymMat::operator*(const AlMat& m) const
{
if( size() != m.nrow() ) {
- std::cerr << "AlSymMat::operator*: size do not match!" << std::endl;
- return m;
+ throw std::range_error("AlSymMat::operator*: size do not match!" );
}
AlMat b( size(), m.ncol());
@@ -321,9 +311,9 @@ AlMat AlSymMat::operator*(const AlMat& m) const
x = 0.;
ii = (i+1)*i/2;
for( k=0; k<i; k++)
- x += (*(ptr_data+ii+k))*m.elemc(k,j);
+ x += (*(m_ptr_data+ii+k))*m.elemc(k,j);
for( k=i; k<siz; k++)
- x += (*(ptr_data+(k+1)*k/2+i))*m.elemc(k,j);
+ x += (*(m_ptr_data+(k+1)*k/2+i))*m.elemc(k,j);
b.elemr(i,j) = x;
}
@@ -334,8 +324,7 @@ AlMat AlSymMat::operator*(const AlMat& m) const
AlVec AlSymMat::operator*(const AlVec& v) const
{
if( size() != v.size() ) {
- std::cerr << "AlSymMat::operator*: size do not match! " << std::endl;
- return v;
+ throw std::range_error( "AlSymMat::operator*: size do not match! " );
}
AlVec b(size());
@@ -345,9 +334,9 @@ AlVec AlSymMat::operator*(const AlVec& v) const
for( long int i=0; i<size(); i++ ) {
ii = (i+1)*i/2;
for( long int j=0; j<i; j++ )
- b[i] += (*(ptr_data+ii+j))*v[j];
+ b[i] += (*(m_ptr_data+ii+j))*v[j];
for( long int j=i; j<size(); j++ )
- b[i] += (*(ptr_data+(j+1)*j/2+i))*v[j];
+ b[i] += (*(m_ptr_data+(j+1)*j/2+i))*v[j];
}
return b;
@@ -356,8 +345,8 @@ AlVec AlSymMat::operator*(const AlVec& v) const
//______________________________________________________________________________
AlSymMat& AlSymMat::operator*=(const double& d)
{
- double * p = ptr_data + _nele;
- while (p > ptr_data) *(--p) *= d;
+ double * p = m_ptr_data + m_nele;
+ while (p > m_ptr_data) *(--p) *= d;
return *this;
}
@@ -366,9 +355,9 @@ AlSymMat& AlSymMat::operator*=(const double& d)
AlSymMat AlSymMat::operator*(const double& d) const
{
AlSymMat a(size());
- double * p = ptr_data + _nele;
- double * q = a.ptr_data + _nele;
- while (p > ptr_data) *(--q) = (*(--p))*d;
+ double * p = m_ptr_data + m_nele;
+ double * q = a.m_ptr_data + m_nele;
+ while (p > m_ptr_data) *(--q) = (*(--p))*d;
return a;
}
@@ -387,7 +376,7 @@ int AlSymMat::invert()
int N = size();
int * ipiv = new int[N];
double * work = new double[N];
- double * mx = ptr_data;
+ double * mx = m_ptr_data;
dsptrf_(&uplo, &N, mx, ipiv, &ierr);
@@ -475,7 +464,7 @@ int AlSymMat::diagonalize(char jobz, AlVec &w, AlMat &z)
int N = size();
double * work = new double[3*N];
- double * ap = ptr_data;
+ double * ap = m_ptr_data;
dspev_(&jobz, &uplo, &N, ap, w.ptrData(), z.ptrData(), &N, work, &info);
@@ -494,12 +483,12 @@ void AlSymMat::RemoveModule(int index)
for(int row=(shift*index); row<(size()-shift); row+=shift) {
for(int k=0; k<shift; k++) {
for(int col=0; col<row+k+1; col++) {
- // cout << "(" << row+k << "," << col << ") -> " << *(ptr_data+elem(row+k,col)) << endl;
+ // cout << "(" << row+k << "," << col << ") -> " << *(m_ptr_data+elem(row+k,col)) << endl;
if(col<shift*index)
- *(ptr_data+elem(row+k,col))=*(ptr_data+elem(row+k+shift,col));
+ *(m_ptr_data+elem(row+k,col))=*(m_ptr_data+elem(row+k+shift,col));
else
- *(ptr_data+elem(row+k,col))=*(ptr_data+elem(row+k+shift,col+shift));
+ *(m_ptr_data+elem(row+k,col))=*(m_ptr_data+elem(row+k+shift,col+shift));
}
}
}
@@ -512,12 +501,10 @@ int AlSymMat::RemoveCollsRows(std::vector<int> indices)
{
int n = indices.size();
if (n==0) {
- std::cerr<<"Vector of indices to remove is empty."<<std::endl;
- return _size;
+ return m_size;
}
- if (n>_size) {
- std::cerr<<"Vector of indices larger than matrix size."<<std::endl;
- return _size;
+ if (n>m_size) {
+ throw std::range_error("Vector of indices larger than matrix size.");
}
// first sort the list of indices descending
@@ -533,22 +520,21 @@ int AlSymMat::RemoveCollsRows(std::vector<int> indices)
// remove rows and columns starting from largest indices
for (int i=0;i<n;i++) {
int index = indices[i];
- if (index > _size-1) {
- std::cerr<<"Index "<<index<<" goes beyond matrix (size "<<_size<<")."<<std::endl;
- continue;
+ if (index > m_size-1) {
+ throw std::out_of_range("AlSymMat::RemoveCollsRows: Index goes beyond matrix.");
}
- for (int j=index; j<_size-1; j++)
+ for (int j=index; j<m_size-1; j++)
for (int col=0; col<=j; col++)
if(col<index)
- *(ptr_data+elem(j,col)) = *(ptr_data+elem(j+1,col));
+ *(m_ptr_data+elem(j,col)) = *(m_ptr_data+elem(j+1,col));
else
- *(ptr_data+elem(j,col)) = *(ptr_data+elem(j+1,col+1));
+ *(m_ptr_data+elem(j,col)) = *(m_ptr_data+elem(j+1,col+1));
- _size--;
+ m_size--;
}
- return _size;
+ return m_size;
}
//______________________________________________________________________________
@@ -577,7 +563,7 @@ void AlSymMat::RemoveAlignPar(int index, int control)
counterCol=0;
}
- *(ptr_data+elem(row,col))=*(ptr_data+elem(row+shiftRow,col+shiftCol));
+ *(m_ptr_data+elem(row,col))=*(m_ptr_data+elem(row+shiftRow,col+shiftCol));
counterCol++;
if (counterCol==5-control) {
@@ -639,7 +625,7 @@ StatusCode AlSymMat::Write(const std::string &filename, bool binary,
if (!square) {
for( int i=0; i<size(); i++) {
for( int j=0; j<=i; j++) {
- melem = *(ptr_data+(i+1)*i/2+j);
+ melem = *(m_ptr_data+(i+1)*i/2+j);
if(binary)
outmat.write((char*)&(melem), sizeof (melem));
else
@@ -652,8 +638,8 @@ StatusCode AlSymMat::Write(const std::string &filename, bool binary,
else {
for( int i=0; i<size(); i++) {
for( int j=0; j<size(); j++) {
- if(i>=j) melem = *(ptr_data+(i+1)*i/2+j);
- else melem = *(ptr_data+(j+1)*j/2+i);
+ if(i>=j) melem = *(m_ptr_data+(i+1)*i/2+j);
+ else melem = *(m_ptr_data+(j+1)*j/2+i);
if(binary)
outmat.write((char*)&(melem), sizeof (melem));
else
@@ -709,7 +695,7 @@ StatusCode AlSymMat::Read(const std::string &filename, int &dofs,
int32_t msiz=0;
inmat.read((char*)&msiz, sizeof (msiz));
dofs = abs(msiz);
- _size = abs(msiz);
+ m_size = abs(msiz);
if (m_StdUnits)
inmat.read((char*)&version, sizeof (version));
@@ -722,7 +708,7 @@ StatusCode AlSymMat::Read(const std::string &filename, int &dofs,
for(int j=0; j<msiz; j++) {
inmat.read((char*)&melem, sizeof (melem));
if( i>=j )
- *(ptr_data+(i+1)*i/2+j) = melem;
+ *(m_ptr_data+(i+1)*i/2+j) = melem;
}
}
}
@@ -733,7 +719,7 @@ StatusCode AlSymMat::Read(const std::string &filename, int &dofs,
for( int i=0; i<msiz; i++) {
for( int j=0; j<=i; j++) {
inmat.read((char*)&melem, sizeof (melem));
- *(ptr_data+(i+1)*i/2+j) = melem;
+ *(m_ptr_data+(i+1)*i/2+j) = melem;
}
}
}
@@ -753,7 +739,7 @@ StatusCode AlSymMat::ReadProjected(const std::string &filename, int &dofs,
int32_t msiz=0;
inmat.read((char*)&msiz, sizeof (msiz));
dofs = abs(msiz);
- _size = abs(msiz);
+ m_size = abs(msiz);
inmat.read((char*)&version, sizeof (version));
@@ -765,7 +751,7 @@ StatusCode AlSymMat::ReadProjected(const std::string &filename, int &dofs,
for(int j=0; j<msiz; j++) {
inmat.read((char*)&melem, sizeof (melem));
if( i>=j )
- *(ptr_data+(i+1)*i/2+j) = melem;
+ *(m_ptr_data+(i+1)*i/2+j) = melem;
}
}
}
@@ -776,7 +762,7 @@ StatusCode AlSymMat::ReadProjected(const std::string &filename, int &dofs,
for( int i=0; i<msiz; i++) {
for( int j=0; j<=i; j++) {
inmat.read((char*)&melem, sizeof (melem));
- *(ptr_data+(i+1)*i/2+j) = melem;
+ *(m_ptr_data+(i+1)*i/2+j) = melem;
}
}
}
@@ -789,18 +775,18 @@ StatusCode AlSymMat::ReadProjected(const std::string &filename, int &dofs,
void AlSymMat::reSize(long int Nnew)
{
if ( Nnew != size() ) {
- double * p = ptr_data;
+ double * p = m_ptr_data;
long int ii;
- _nele = Nnew*(Nnew+1)/2;
- ptr_data = new double[_nele];
+ m_nele = Nnew*(Nnew+1)/2;
+ m_ptr_data = new double[m_nele];
long int Size_old = size();
- _size = Nnew;
+ m_size = Nnew;
long int k = size() <= Size_old ? size() : Size_old;
for( long int i=0; i<k; i++ ) {
ii = (i+1)*i/2;
for( long int j=0; j<=i; j++ )
- *(ptr_data+ii+j) = *(p+ii+j);
+ *(m_ptr_data+ii+j) = *(p+ii+j);
}
delete [] p;
}
@@ -811,26 +797,22 @@ double& AlSymMat::elemr(long int i,long int j)
{
#ifdef _DEBUG
if( i<0 ) {
- std::cerr << "AlSymMat::elemr: Index 1 < zero! " << i << std::endl;
- return *(ptr_data);
+ throw std::underflow_error( "AlSymMat::elemr: Index 1 < zero!" );
}
if( i>=size() ) {
- std::cerr << "AlSymMat::elemr: Index 1 too large! " << i << std::endl;
- return *(ptr_data);
+ throw std::overflow_error( "AlSymMat::elemr: Index 1 too large!" );
}
if( j<0 ) {
- std::cerr << "AlSymMat::elemr: Index 2 < zero! " << j << std::endl;
- return *(ptr_data);
+ throw std::underflow_error( "AlSymMat::elemr: Index 2 < zero!" );
}
if( j>=size() ) {
- std::cerr << "AlSymMat::elemr: Index 2 too large! " << j << std::endl;
- return *(ptr_data);
+ throw std::overflow_error( "AlSymMat::elemr: Index 2 too large!" );
}
#endif
if( j<=i )
- return *(ptr_data+(i+1)*i/2+j);
+ return *(m_ptr_data+(i+1)*i/2+j);
- return *(ptr_data+(j+1)*j/2+i);
+ return *(m_ptr_data+(j+1)*j/2+i);
}
//______________________________________________________________________________
@@ -838,26 +820,22 @@ double AlSymMat::elemc(long int i,long int j) const
{
#ifdef _DEBUG
if( i<0 ) {
- std::cerr << "AlSymMat::elemc: Index 1 < zero! " << i << std::endl;
- return *(ptr_data);
+ throw std::underflow_error( "AlSymMat::elemc: Index 1 < zero!" );
}
if( i>=size() ) {
- std::cerr << "AlSymMat::elemc: Index 1 too large! " << i << std::endl;
- return *(ptr_data);
+ throw std::overflow_error( "AlSymMat::elemc: Index 1 too large!" );
}
if( j<0 ) {
- std::cerr << "AlSymMat::elemc: Index 2 < zero! " << j << std::endl;
- return *(ptr_data);
+ throw std::underflow_error( "AlSymMat::elemc: Index 2 < zero!" );
}
if( j>=size() ) {
- std::cerr << "AlSymMat::elemc: Index 2 too large! " << j << std::endl;
- return *(ptr_data);
+ throw std::overflow_error( "AlSymMat::elemc: Index 2 too large!" );
}
#endif
if( j<=i )
- return *(ptr_data+(i+1)*i/2+j);
+ return *(m_ptr_data+(i+1)*i/2+j);
- return *(ptr_data+(j+1)*j/2+i);
+ return *(m_ptr_data+(j+1)*j/2+i);
}
@@ -869,7 +847,7 @@ TMatrixDSparse* AlSymMat::makeTMatrix(){
for( int i=0; i<si; i++ )
for( int j=0; j<=i; j++ )
- if( *(ptr_data+(i+1)*i/2+j) != 0 )
+ if( *(m_ptr_data+(i+1)*i/2+j) != 0 )
++nnz;
@@ -880,10 +858,10 @@ TMatrixDSparse* AlSymMat::makeTMatrix(){
for (int i=0;i<si;i++) {
for (int j=0;j<=i;j++) {
- if ( *(ptr_data+(i+1)*i/2+j) != 0.) {
+ if ( *(m_ptr_data+(i+1)*i/2+j) != 0.) {
irow[n] = i;
icol[n] = j;
- val[n] = *(ptr_data+(i+1)*i/2+j);
+ val[n] = *(m_ptr_data+(i+1)*i/2+j);
++n;
}
}
diff --git a/Tracking/TrkAlignment/TrkAlgebraUtils/src/AlSymMatBase.cxx b/Tracking/TrkAlignment/TrkAlgebraUtils/src/AlSymMatBase.cxx
index 2da1a90f68af883201cdc8cd59f6a50fc27f5d59..f5c61b7a228f743ebcffbfe7b0e67476896af7da 100644
--- a/Tracking/TrkAlignment/TrkAlgebraUtils/src/AlSymMatBase.cxx
+++ b/Tracking/TrkAlignment/TrkAlgebraUtils/src/AlSymMatBase.cxx
@@ -19,26 +19,26 @@ namespace Trk {
//______________________________________________________________________________
AlSymMatBase::AlSymMatBase()
{
- _matrix_type = 0;
- _size = 0;
- _nele = 0;
+ m_matrix_type = 0;
+ m_size = 0;
+ m_nele = 0;
}
//______________________________________________________________________________
AlSymMatBase::AlSymMatBase(long int N)
{
- _matrix_type = 0;
- _size = N;
- _nele = N*(N+1)/2;
+ m_matrix_type = 0;
+ m_size = N;
+ m_nele = N*(N+1)/2;
}
//______________________________________________________________________________
AlSymMatBase::AlSymMatBase(const AlSymMatBase & m)
{
- _matrix_type = 0;
- _size = m.size();
- _nele = m._nele;
+ m_matrix_type = 0;
+ m_size = m.size();
+ m_nele = m.m_nele;
}
//______________________________________________________________________________
@@ -46,8 +46,8 @@ AlSymMatBase & AlSymMatBase::operator=(const AlSymMatBase & m)
{
if (this==&m)
return *this;
- _matrix_type=m._matrix_type;
- _nele=m._nele;
+ m_matrix_type=m.m_matrix_type;
+ m_nele=m.m_nele;
return *this;
}
diff --git a/Tracking/TrkAlignment/TrkAlgebraUtils/src/AlVec.cxx b/Tracking/TrkAlignment/TrkAlgebraUtils/src/AlVec.cxx
index 35b3893858a21e1719309d26efb7aa5d74167514..ccdafd5be070a1508e4b60be39622e63f0544520 100644
--- a/Tracking/TrkAlignment/TrkAlgebraUtils/src/AlVec.cxx
+++ b/Tracking/TrkAlignment/TrkAlgebraUtils/src/AlVec.cxx
@@ -12,34 +12,35 @@
#include <fstream>
#include <stdint.h>
#include <math.h>
+#include <exception>
namespace Trk {
//______________________________________________________________________________
AlVec::AlVec() {
- _size = 0;
- ptr_data = 0; // set pointer to null
+ m_size = 0;
+ m_ptr_data = 0; // set pointer to null
m_pathbin="./";
m_pathtxt="./";
}
//______________________________________________________________________________
AlVec::AlVec(int N) {
- _size = N;
- ptr_data = new double[_size];
+ m_size = N;
+ m_ptr_data = new double[m_size];
m_pathbin="./";
m_pathtxt="./";
- double* p = ptr_data + _size;
- while (p > ptr_data) *(--p) = 0.;
+ double* p = m_ptr_data + m_size;
+ while (p > m_ptr_data) *(--p) = 0.;
}
//______________________________________________________________________________
AlVec::AlVec(const AlVec& v) {
- _size = v._size;
+ m_size = v.m_size;
m_pathbin = v.m_pathbin;
m_pathtxt = v.m_pathtxt;
- ptr_data = new double[_size];
+ m_ptr_data = new double[m_size];
copy(v);
}
@@ -47,26 +48,25 @@ AlVec::AlVec(const AlVec& v) {
//______________________________________________________________________________
AlVec::~AlVec()
{
- delete [] ptr_data;
+ delete [] m_ptr_data;
}
//______________________________________________________________________________
void AlVec::copy(const AlVec& v) {
- if(_size!=v._size) {
- std::cerr << "AlVec Assignment: size does not match!" << std::endl;
- return;
+ if(m_size!=v.m_size) {
+ throw std::range_error( "AlVec Assignment: size does not match!" );
}
- double* p = ptr_data + _size;
- double* q = v.ptr_data + _size;
- while (p > ptr_data) *(--p) = (*(--q));
+ double* p = m_ptr_data + m_size;
+ double* q = v.m_ptr_data + m_size;
+ while (p > m_ptr_data) *(--p) = (*(--q));
}
//______________________________________________________________________________
AlVec& AlVec::operator=(const double& d) {
- double* p = ptr_data + _size;
- while (p > ptr_data) *(--p) = d;
+ double* p = m_ptr_data + m_size;
+ while (p > m_ptr_data) *(--p) = d;
return *this;
}
@@ -75,74 +75,69 @@ AlVec& AlVec::operator=(const double& d) {
AlVec& AlVec::operator=(const AlVec& v) {
if (this==&v) return *this;
- if(_size!=0 && _size!=v._size) {
- std::cerr << "AlVec Assignment: size does not match!" << std::endl;
- return *this;
+ if(m_size!=0 && m_size!=v.m_size) {
+ throw std::range_error( "AlVec Assignment: size does not match!" );
}
- if ( ptr_data != v.ptr_data ) copy(v);
+ if ( m_ptr_data != v.m_ptr_data ) copy(v);
return *this;
}
//______________________________________________________________________________
AlVec AlVec::operator+(const AlVec& v) const {
- if( _size != v._size ) {
- std::cerr << "operator+: vectors size does not match!" << std::endl;
- return *this;
+ if( m_size != v.m_size ) {
+ throw std::range_error( "operator+: vectors size does not match!" );
}
- AlVec b(_size);
+ AlVec b(m_size);
- double* p = ptr_data + _size;
- double* q = v.ptr_data + _size;
- double* r = b.ptr_data + _size;
- while (p > ptr_data) *(--r) = (*(--p))+(*(--q));
+ double* p = m_ptr_data + m_size;
+ double* q = v.m_ptr_data + m_size;
+ double* r = b.m_ptr_data + m_size;
+ while (p > m_ptr_data) *(--r) = (*(--p))+(*(--q));
return b;
}
//______________________________________________________________________________
AlVec& AlVec::operator+=(const AlVec& v) {
- if( _size != v._size ) {
- std::cerr << "operator+=: vectors size does not match!" << std::endl;
- return *this;
+ if( m_size != v.m_size ) {
+ throw std::range_error( "operator+=: vectors size does not match!" );
}
- double* p = ptr_data + _size;
- double* q = v.ptr_data + _size;
- while (p > ptr_data) *(--p) += (*(--q));
+ double* p = m_ptr_data + m_size;
+ double* q = v.m_ptr_data + m_size;
+ while (p > m_ptr_data) *(--p) += (*(--q));
return *this;
}
//______________________________________________________________________________
AlVec AlVec::operator-(const AlVec& v) const {
- if( _size != v._size ) {
- std::cerr << "operator-: vectors size does not match!" << std::endl;
- return *this;
+ if( m_size != v.m_size ) {
+ throw std::range_error( "operator-: vectors size does not match!" );
}
- AlVec b(_size);
+ AlVec b(m_size);
- double* p = ptr_data + _size;
- double* q = v.ptr_data + _size;
- double* r = b.ptr_data + _size;
- while (p > ptr_data) *(--r) = (*(--p))-(*(--q));
+ double* p = m_ptr_data + m_size;
+ double* q = v.m_ptr_data + m_size;
+ double* r = b.m_ptr_data + m_size;
+ while (p > m_ptr_data) *(--r) = (*(--p))-(*(--q));
return b;
}
//______________________________________________________________________________
AlVec& AlVec::operator-=(const AlVec& v) {
- if( _size != v._size ) {
- std::cerr << "operator+=: vectors size does not match!" << std::endl;
- return *this;
+ if( m_size != v.m_size ) {
+ throw std::range_error( "operator+=: vectors size does not match!" );
}
- double* p = ptr_data + _size;
- double* q = v.ptr_data + _size;
- while (p > ptr_data) *(--p) -= (*(--q));
+ double* p = m_ptr_data + m_size;
+ double* q = v.m_ptr_data + m_size;
+ while (p > m_ptr_data) *(--p) -= (*(--q));
return *this;
}
@@ -150,44 +145,41 @@ AlVec& AlVec::operator-=(const AlVec& v) {
//______________________________________________________________________________
double AlVec::operator*(const AlVec& v) const {
double b=0.;
- if( _size != v._size ) {
- std::cerr << "scalar product: vectors size does not match!" << std::endl;
- return b;
+ if( m_size != v.m_size ) {
+ throw std::range_error( "scalar product: vectors size does not match!" );
}
- double* p = ptr_data + _size;
- double* q = v.ptr_data + _size;
- while (p > ptr_data) b += *(--p)*(*(--q));
+ double* p = m_ptr_data + m_size;
+ double* q = v.m_ptr_data + m_size;
+ while (p > m_ptr_data) b += *(--p)*(*(--q));
return b;
}
//______________________________________________________________________________
AlVec AlVec::operator*(const AlMat& m) const {
- if (_size != m.nrow()) {
- std::cerr << "Left hand vector-matrix multiplication: size does not match!" << std::endl;
- return *this;
+ if (m_size != m.nrow()) {
+ throw std::range_error( "Left hand vector-matrix multiplication: size does not match!" );
}
AlVec b(m.ncol());
for (int i=0;i<m.ncol();i++) {
- for (int j=0;j<_size;j++) *(b.ptr_data+i) += *(ptr_data+j)*m.elemc(j, i);
+ for (int j=0;j<m_size;j++) *(b.m_ptr_data+i) += *(m_ptr_data+j)*m.elemc(j, i);
}
return b;
}
//______________________________________________________________________________
AlVec AlVec::operator*(const AlSymMatBase& m) const {
- if (_size != m.size()) {
- std::cerr << "Left hand vector-matrix multiplication: size does not match!" << std::endl;
- return *this;
+ if (m_size != m.size()) {
+ throw std::range_error( "Left hand vector-matrix multiplication: size does not match!" );
}
- AlVec b(_size);
+ AlVec b(m_size);
- for (int i=0;i<_size;i++) {
- for (int j=0;j<_size;j++) *(b.ptr_data+i) += *(ptr_data+j)*m.elemc(j, i);
+ for (int i=0;i<m_size;i++) {
+ for (int j=0;j<m_size;j++) *(b.m_ptr_data+i) += *(m_ptr_data+j)*m.elemc(j, i);
}
return b;
}
@@ -195,8 +187,8 @@ AlVec AlVec::operator*(const AlSymMatBase& m) const {
//______________________________________________________________________________
AlVec& AlVec::operator*=(const double& d) {
- double* p = ptr_data + _size;
- while (p > ptr_data) *(--p) *= d;
+ double* p = m_ptr_data + m_size;
+ while (p > m_ptr_data) *(--p) *= d;
return *this;
}
@@ -207,8 +199,8 @@ AlVec& AlVec::operator*=(const double& d) {
double AlVec::norm() const {
double result(0.);
- double * p = ptr_data + _size;
- while (p > ptr_data) {
+ double * p = m_ptr_data + m_size;
+ while (p > m_ptr_data) {
--p;
result += *p * *p ;
}
@@ -218,16 +210,16 @@ double AlVec::norm() const {
//______________________________________________________________________________
void AlVec::reSize(int Nnew) {
- if ( Nnew>=0 && Nnew != _size ) {
- double* p = ptr_data;
- int _size_old = _size;
- ptr_data = new double[Nnew];
- _size = Nnew;
- int k = _size <= _size_old ? _size : _size_old;
+ if ( Nnew>=0 && Nnew != m_size ) {
+ double* p = m_ptr_data;
+ int m_size_old = m_size;
+ m_ptr_data = new double[Nnew];
+ m_size = Nnew;
+ int k = m_size <= m_size_old ? m_size : m_size_old;
p += k;
- double* q = ptr_data + k;
- while (q > ptr_data) *(--q) = *(--p);
+ double* q = m_ptr_data + k;
+ while (q > m_ptr_data) *(--q) = *(--p);
delete [] p;
}
@@ -239,9 +231,9 @@ void AlVec::RemoveModule(int index) {
const int shift=6;
const int start = shift*index;
- if(start<_size && index>=0){
- for(int row=start; row<_size; row++) {
- *(ptr_data+row)=*(ptr_data+row+shift);
+ if(start<m_size && index>=0){
+ for(int row=start; row<m_size; row++) {
+ *(m_ptr_data+row)=*(m_ptr_data+row+shift);
}
}
}
@@ -256,8 +248,8 @@ void AlVec::RemoveAlignPar(int index, int control)
index = index-control;
// std::cout << "index: " << index << " - control: " << control << std::endl;
- for(int row=index; row<_size; row++) {
- *(ptr_data+row) = *(ptr_data+row+shift);
+ for(int row=index; row<m_size; row++) {
+ *(m_ptr_data+row) = *(m_ptr_data+row+shift);
counter++;
if (counter==5-control) { counter=0; shift++; }
}
@@ -268,12 +260,10 @@ int AlVec::RemoveElements(std::vector<int> indices)
{
int n = indices.size();
if (n==0) {
- std::cerr<<"Vector of indices to remove is empty."<<std::endl;
- return _size;
+ return m_size;
}
- if (n>_size) {
- std::cerr<<"Vector of indices larger than matrix size."<<std::endl;
- return _size;
+ if (n>m_size) {
+ throw std::range_error( "Vector of indices larger than matrix size." );
}
// first sort the list of indices descending
@@ -289,18 +279,18 @@ int AlVec::RemoveElements(std::vector<int> indices)
// remove elements starting from largest indices
for (int i=0;i<n;i++) {
int index = indices[i];
- if (index > _size-1) {
- std::cerr<<"Index "<<index<<" goes beyond matrix (size "<<_size<<")."<<std::endl;
+ if (index > m_size-1) {
+ throw std::out_of_range( "AlVec::RemoveElements: Index goes beyond matrix " );
continue;
}
- for (int j=index; j<_size-1; j++)
- *(ptr_data+j) = *(ptr_data+j+1);
+ for (int j=index; j<m_size-1; j++)
+ *(m_ptr_data+j) = *(m_ptr_data+j+1);
- _size--;
+ m_size--;
}
- return _size;
+ return m_size;
}
//______________________________________________________________________________
@@ -320,7 +310,7 @@ StatusCode AlVec::Write(const std::string &filename, bool binary, double scale,
std::map<int,unsigned long long> moduleIndexMap, float version)
{
std::ofstream outvec;
- int32_t io_size=_size;
+ int32_t io_size=m_size;
// int32_t io_scale=scale;
if(binary) {
@@ -338,7 +328,7 @@ StatusCode AlVec::Write(const std::string &filename, bool binary, double scale,
outvec.setf(std::ios::fixed);
outvec.setf(std::ios::showpoint);
outvec.precision(6);
- outvec << "DoF: " << std::setw(6) << _size << std::endl;
+ outvec << "DoF: " << std::setw(6) << m_size << std::endl;
outvec << "scale: " << std::setw(18) << scale << std::endl;
outvec << "AlVec version: " << std::setw(6) << version << std::endl;
}
@@ -347,13 +337,13 @@ StatusCode AlVec::Write(const std::string &filename, bool binary, double scale,
double velem=0;
std::map<int,unsigned long long>::iterator itcod;
- for( int i=0; i<_size; i++) {
+ for( int i=0; i<m_size; i++) {
itcod=moduleIndexMap.find(i/6);
if (itcod != moduleIndexMap.end()) {
ielem = (itcod->second);
}
- velem = *(ptr_data+i);
+ velem = *(m_ptr_data+i);
if(binary){
outvec.write((char*)&(ielem), sizeof (ielem));
@@ -373,7 +363,7 @@ StatusCode AlVec::WritePartial(const std::string &filename, bool binary, double
std::ofstream outvec;
/*
- int32_t io_size=_size;
+ int32_t io_size=m_size;
// int32_t io_scale=scale;
if(binary) {
@@ -392,7 +382,7 @@ StatusCode AlVec::WritePartial(const std::string &filename, bool binary, double
outvec.setf(std::ios::fixed);
outvec.setf(std::ios::showpoint);
outvec.precision(6);
- outvec << "DoF: " << std::setw(6) << _size << std::endl;
+ outvec << "DoF: " << std::setw(6) << m_size << std::endl;
outvec << "scale: " << std::setw(18) << scale << std::endl;
outvec << "AlVec version: " << std::setw(6) << version << std::endl;
}
@@ -405,12 +395,12 @@ StatusCode AlVec::WritePartial(const std::string &filename, bool binary, double
double velem=0;
std::map<int,unsigned long long>::iterator itcod;
- for( int i=0; i<_size; i++) {
+ for( int i=0; i<m_size; i++) {
itcod=moduleIndexMap.find(i);
if (itcod != moduleIndexMap.end())
ielem = (itcod->second);
- velem = *(ptr_data+i);
+ velem = *(m_ptr_data+i);
if(binary){
outvec.write((char*)&(ielem), sizeof (ielem));
@@ -436,13 +426,13 @@ StatusCode AlVec::WritePartial(const std::string &filename, bool binary, double
double velem=0;
std::map<int,std::string>::iterator itcod;
- for( int i=0; i<_size; i++) {
+ for( int i=0; i<m_size; i++) {
itcod=moduleNameMap.find(i);
if (itcod != moduleNameMap.end()) {
elem = (itcod->second);
}
- velem = *(ptr_data+i);
+ velem = *(m_ptr_data+i);
if(binary){
std::cout<<"can't write module name in binary output!"<<std::endl;
@@ -459,7 +449,7 @@ StatusCode AlVec::WritePartial(const std::string &filename, bool binary, double
StatusCode AlVec::InitializeOutputVector(const std::string& filename, bool binary, double scale,
float version, std::ofstream& outvec)
{
- int32_t io_size=_size;
+ int32_t io_size=m_size;
if(binary) {
outvec.open((m_pathbin+filename).c_str(), std::ios::binary);
@@ -477,7 +467,7 @@ StatusCode AlVec::InitializeOutputVector(const std::string& filename, bool binar
outvec.setf(std::ios::fixed);
outvec.setf(std::ios::showpoint);
outvec.precision(6);
- outvec << "DoF: " << std::setw(6) << _size << std::endl;
+ outvec << "DoF: " << std::setw(6) << m_size << std::endl;
outvec << "scale: " << std::setw(18) << scale << std::endl;
outvec << "AlVec version: " << std::setw(6) << version << std::endl;
}
@@ -490,19 +480,19 @@ StatusCode AlVec::ReadPartial(const std::string &filename, double &scale,
{
bool m_StdUnits = true;
if (StatusCode::SUCCESS != CheckVecVersion(m_pathbin+filename, m_StdUnits)) {
- std::cout<<"CheckVecVersion failed"<<std::endl;
+ //std::cout<<"CheckVecVersion failed"<<std::endl;
return StatusCode::FAILURE;
}
std::ifstream invec((m_pathbin+filename).c_str(), std::ios::binary);
if(invec.fail()) {
- std::cout<<"ifstream failed"<<std::endl;
+ //std::cout<<"ifstream failed"<<std::endl;
return StatusCode::FAILURE;
}
int32_t vsiz=0;
invec.read((char*)&vsiz, sizeof (vsiz));
- _size = vsiz;
+ m_size = vsiz;
// int32_t io_scale;
invec.read((char*)&scale, sizeof(scale));
@@ -517,12 +507,12 @@ StatusCode AlVec::ReadPartial(const std::string &filename, double &scale,
int64_t ielem=0;
double velem=0.0;
- for( int i=0; i<_size; i++) {
+ for( int i=0; i<m_size; i++) {
invec.read((char*)&ielem, sizeof (ielem));
modmap[i] = ielem;
invec.read((char*)&velem, sizeof (velem));
- *(ptr_data+i) = velem;
+ *(m_ptr_data+i) = velem;
// std::cout << "AlVec (" << i << "):: " << ielem << ", " << velem << std::endl;
@@ -568,7 +558,7 @@ StatusCode AlVec::Read(const std::string &filename, double &scale,
bool m_StdUnits = true;
if (StatusCode::SUCCESS != CheckVecVersion(m_pathbin+filename, m_StdUnits)) {
- std::cout<<"CheckVecVersion failed"<<std::endl;
+ //std::cout<<"CheckVecVersion failed"<<std::endl;
return StatusCode::FAILURE;
}
@@ -578,8 +568,8 @@ StatusCode AlVec::Read(const std::string &filename, double &scale,
int32_t vsiz=0;
invec.read((char*)&vsiz, sizeof (vsiz));
- _size = vsiz;
-// std::cout<<"size="<<_size<<std::endl;
+ m_size = vsiz;
+// std::cout<<"size="<<m_size<<std::endl;
// int32_t io_scale;
invec.read((char*)&scale, sizeof (scale));
@@ -595,12 +585,12 @@ StatusCode AlVec::Read(const std::string &filename, double &scale,
int64_t ielem=0;
double velem=0.0;
- for( int i=0; i<_size; i++) {
+ for( int i=0; i<m_size; i++) {
invec.read((char*)&ielem, sizeof (ielem));
modmap[i/6] = ielem;
invec.read((char*)&velem, sizeof (velem));
- *(ptr_data+i) = velem;
+ *(m_ptr_data+i) = velem;
// std::cout << "AlVec (" << i << "):: " << ielem << ", " << velem << std::endl;
@@ -620,7 +610,7 @@ StatusCode AlVec::ReadProjected(const std::string &filename, double &scale,
int32_t vsiz=0;
invec.read((char*)&vsiz, sizeof (vsiz));
- _size = vsiz;
+ m_size = vsiz;
// int32_t io_scale;
invec.read((char*)&scale, sizeof (scale));
@@ -633,12 +623,12 @@ StatusCode AlVec::ReadProjected(const std::string &filename, double &scale,
int64_t ielem=0;
double velem=0.0;
- for( int i=0; i<_size; i++) {
+ for( int i=0; i<m_size; i++) {
invec.read((char*)&ielem, sizeof (ielem));
modmap[i/6] = ielem;
invec.read((char*)&velem, sizeof (velem));
- *(ptr_data+i) = velem;
+ *(m_ptr_data+i) = velem;
// std::cout << "AlVec (" << i << "):: " << ielem << ", " << velem << std::endl;
@@ -656,13 +646,13 @@ StatusCode AlVec::ReadScalaPack(const std::string &filename){
int32_t vsiz=0;
eigenvec.read((char*)&vsiz, sizeof (vsiz));
- _size=vsiz;
+ m_size=vsiz;
double velem=0;
- for( int i=0; i<_size; i++) {
+ for( int i=0; i<m_size; i++) {
eigenvec.read((char*)&velem, sizeof (velem));
// printf("v[%d] = %.16lf \n",i,velem);
- *(ptr_data+i) = velem;
+ *(m_ptr_data+i) = velem;
}
eigenvec.close();
@@ -673,7 +663,7 @@ StatusCode AlVec::ReadScalaPack(const std::string &filename){
StatusCode AlVec::WriteEigenvalueVec(const std::string &filename, bool binary){
std::ofstream outvec;
- int32_t io_size=_size;
+ int32_t io_size=m_size;
if(binary) {
outvec.open((m_pathbin+filename).c_str(), std::ios::binary);
@@ -692,16 +682,16 @@ StatusCode AlVec::WriteEigenvalueVec(const std::string &filename, bool binary){
outvec.setf(std::ios::fixed);
outvec.setf(std::ios::showpoint);
outvec.precision(6);
- outvec << "AlVec::DoF: " << std::setw(6) << _size << std::endl;
+ outvec << "AlVec::DoF: " << std::setw(6) << m_size << std::endl;
}
//jlove int32_t ielem=0;
double velem=0;
// std::map<int,int>::iterator itcod;
- for( int i=0; i<_size; i++) {
+ for( int i=0; i<m_size; i++) {
- velem = *(ptr_data+i);
+ velem = *(m_ptr_data+i);
if(binary)
outvec.write((char*)&(velem), sizeof (velem));
diff --git a/Tracking/TrkAlignment/TrkAlgebraUtils/src/IntVec.cxx b/Tracking/TrkAlignment/TrkAlgebraUtils/src/IntVec.cxx
index d52767e7b634075cb02dccb4af6c4029dcaaba60..b9dbfde12b1b14e0993986a83f0721366c82460a 100644
--- a/Tracking/TrkAlignment/TrkAlgebraUtils/src/IntVec.cxx
+++ b/Tracking/TrkAlignment/TrkAlgebraUtils/src/IntVec.cxx
@@ -6,48 +6,48 @@
#include "TrkAlgebraUtils/IntVec.h"
#include <iostream>
#include <stdint.h>
+#include <exception>
namespace Trk {
IntVec::IntVec()
- : Nele(0), ptr_data(0)
+ : m_Nele(0), m_ptr_data(0)
{}
IntVec::IntVec(int N) {
- Nele = N;
- ptr_data = new int[Nele];
- for(int i=0;i<Nele;i++)
- *(ptr_data+i)=0;
+ m_Nele = N;
+ m_ptr_data = new int[m_Nele];
+ for(int i=0;i<m_Nele;i++)
+ *(m_ptr_data+i)=0;
}
IntVec::IntVec(int N, int init) {
- Nele = N;
- ptr_data = new int[Nele];
- for(int i=0;i<Nele;i++)
- *(ptr_data+i)=init;
+ m_Nele = N;
+ m_ptr_data = new int[m_Nele];
+ for(int i=0;i<m_Nele;i++)
+ *(m_ptr_data+i)=init;
}
IntVec::IntVec(const IntVec& v) {
- Nele = v.Nele;
- ptr_data = new int[Nele];
- for(int i=0;i<Nele;i++)
- *(ptr_data+i)=v[i];
+ m_Nele = v.m_Nele;
+ m_ptr_data = new int[m_Nele];
+ for(int i=0;i<m_Nele;i++)
+ *(m_ptr_data+i)=v[i];
}
IntVec::~IntVec(){
- delete [] ptr_data;
+ delete [] m_ptr_data;
}
IntVec& IntVec::operator=(const IntVec& v) {
- if(Nele!=0 && Nele!=v.Nele) {
- std::cerr << "IntVec Assignment: size does not match!" << std::endl;
- return *this;
+ if(m_Nele!=0 && m_Nele!=v.m_Nele) {
+ throw std::range_error( "IntVec Assignment: size does not match!" );
}
- if ( ptr_data != v.ptr_data ) {
- reSize(v.Nele);
- for(int i=0;i<Nele;i++)
- *(ptr_data+i)=v[i];
+ if ( m_ptr_data != v.m_ptr_data ) {
+ reSize(v.m_Nele);
+ for(int i=0;i<m_Nele;i++)
+ *(m_ptr_data+i)=v[i];
}
return *this;
@@ -55,92 +55,84 @@ IntVec& IntVec::operator=(const IntVec& v) {
int& IntVec::operator[](int i) {
if(i<0) {
- std::cerr << "IntVec: Index < zero! " << std::endl;
- return ptr_data[0];
+ throw std::out_of_range( "IntVec: Index < zero! " );
}
- else if(i>=Nele) {
- std::cerr << "IntVec: Index too large! " << std::endl;
- return ptr_data[0];
+ else if(i>=m_Nele) {
+ throw std::out_of_range( "IntVec: Index too large! " );
}
- return *(ptr_data+i);
+ return *(m_ptr_data+i);
}
const int& IntVec::operator[](int i) const {
if(i<0) {
- std::cerr << "IntVec: Index < zero! " << std::endl;
- return ptr_data[0];
+ throw std::out_of_range( "IntVec: Index < zero! " );
}
- else if(i>=Nele) {
- std::cerr << "IntVec: Index too large! " << std::endl;
- return ptr_data[0];
+ else if(i>=m_Nele) {
+ throw std::out_of_range( "IntVec: Index too large! " );
}
- return *(ptr_data+i);
+ return *(m_ptr_data+i);
}
IntVec IntVec::operator+(const IntVec& v) {
- if( Nele != v.Nele ) {
- std::cerr << "operator+: vectors size does not match!" << std::endl;
- return *this;
+ if( m_Nele != v.m_Nele ) {
+ throw std::range_error( "operator+: vectors size does not match!" );
}
- IntVec b(Nele);
- for (int i=0;i<Nele;i++)
- b[i] = *(ptr_data+i) + v[i];
+ IntVec b(m_Nele);
+ for (int i=0;i<m_Nele;i++)
+ b[i] = *(m_ptr_data+i) + v[i];
return b;
}
IntVec& IntVec::operator+=(const IntVec& v) {
- if( Nele != v.Nele ) {
- std::cerr << "operator+=: vectors size does not match!" << std::endl;
- return *this;
+ if( m_Nele != v.m_Nele ) {
+ throw std::range_error( "operator+=: vectors size does not match!" );
}
- for (int i=0;i<Nele;i++)
- *(ptr_data+i)+=v[i];
+ for (int i=0;i<m_Nele;i++)
+ *(m_ptr_data+i)+=v[i];
return *this;
}
IntVec IntVec::operator-(const IntVec& v) {
- if( Nele != v.Nele ) {
- std::cerr << "operator+: vectors size does not match!" << std::endl;
- return *this;
+ if( m_Nele != v.m_Nele ) {
+ throw std::range_error( "operator+: vectors size does not match!" );
}
- IntVec b(Nele);
- for (int i=0;i<Nele;i++)
- b[i] = *(ptr_data+i) - v[i];
+ IntVec b(m_Nele);
+ for (int i=0;i<m_Nele;i++)
+ b[i] = *(m_ptr_data+i) - v[i];
return b;
}
IntVec& IntVec::operator-=(const IntVec& v) {
- if( Nele != v.Nele ) {
- std::cerr << "operator+=: vectors size does not match!" << std::endl;
- return *this;
+ if( m_Nele != v.m_Nele ) {
+ throw std::range_error( "operator+=: vectors size does not match!" );
}
- for (int i=0;i<Nele;i++)
- *(ptr_data+i)-=v[i];
+ for (int i=0;i<m_Nele;i++)
+ *(m_ptr_data+i)-=v[i];
return *this;
}
void IntVec::reSize(int Nnew) {
- if ( Nnew>=0 && Nnew != Nele ) {
- int* p = ptr_data;
- int Nele_old = Nele;
- ptr_data = new int[Nnew];
- Nele = Nnew;
- int k = Nele <= Nele_old ? Nele : Nele_old;
+ if ( Nnew>=0 && Nnew != m_Nele ) {
+ int* p = m_ptr_data;
+ int Nele_old = m_Nele;
+ m_ptr_data = new int[Nnew];
+ m_Nele = Nnew;
+ int k = m_Nele <= Nele_old ? m_Nele : Nele_old;
p += k;
- int* q = ptr_data + k;
- while (q > ptr_data)
+ int* q = m_ptr_data + k;
+ while (q > m_ptr_data)
*(--q) = *(--p);
delete [] p;
diff --git a/Tracking/TrkAlignment/TrkAlignGenTools/TrkAlignGenTools/MatrixTool.h b/Tracking/TrkAlignment/TrkAlignGenTools/TrkAlignGenTools/MatrixTool.h
index b3cf31e85c98192c568ea10e1ba3dd183a411cfb..26a6bdb38b7322ec26b043c38ebfe56805f8f7b2 100644
--- a/Tracking/TrkAlignment/TrkAlignGenTools/TrkAlignGenTools/MatrixTool.h
+++ b/Tracking/TrkAlignment/TrkAlignGenTools/TrkAlignGenTools/MatrixTool.h
@@ -58,9 +58,9 @@ namespace Trk {
enum SolveOption {
NONE = 0, //!< not solve in any case (to be used when ipc)
SOLVE = 1, //!< solving after data accumulation (LAPACK)
- SOLVE_FAST = 2, //!< Fast (***REMOVED*** method) solving after data accumulation
+ SOLVE_FAST = 2, //!< Fast (Eigen method) solving after data accumulation
DIRECT_SOLVE = 3, //!< direct solving (LAPACK), already available matrix & vector
- DIRECT_SOLVE_FAST = 4, //!< direct Fast (***REMOVED*** method) solving, already available matrix & vector
+ DIRECT_SOLVE_FAST = 4, //!< direct Fast (Eigen method) solving, already available matrix & vector
DIRECT_SOLVE_CLUSTER = 5, //!< computation of alignment parameters from SCALAPAK already solved matrix
SOLVE_ROOT = 6, //!< computation using ROOT
SOLVE_CLHEP = 7 //!< computation using CLHEP
@@ -128,7 +128,7 @@ namespace Trk {
int solveROOT();
int solveCLHEP();
int solveLapack();
- int solve***REMOVED***();
+ int solveSparseEigen();
int solveLocal();
StatusCode spuriousRemoval();
diff --git a/Tracking/TrkAlignment/TrkAlignGenTools/src/MatrixTool.cxx b/Tracking/TrkAlignment/TrkAlignGenTools/src/MatrixTool.cxx
index a477f5f25402349c89a067cb8cadf1867f9ceec1..b7c9846b05cf126d70934b04570800c968897e85 100644
--- a/Tracking/TrkAlignment/TrkAlignGenTools/src/MatrixTool.cxx
+++ b/Tracking/TrkAlignment/TrkAlignGenTools/src/MatrixTool.cxx
@@ -1426,7 +1426,7 @@ namespace Trk {
case SOLVE_FAST:
case DIRECT_SOLVE_FAST:
- info = solve***REMOVED***();
+ info = solveSparseEigen();
break;
default:
@@ -2147,13 +2147,13 @@ namespace Trk {
}
//________________________________________________________________________
- int MatrixTool::solve***REMOVED***()
+ int MatrixTool::solveSparseEigen()
{
- ATH_MSG_INFO("solving Global using ***REMOVED***");
+ ATH_MSG_INFO("solving Global using SparseEigen");
if(m_logStream) {
*m_logStream<<"*************************************************************"<<std::endl;
*m_logStream<<"************** solving using Global method ****************"<<std::endl;
- *m_logStream<<"************** using ***REMOVED*** ****************"<<std::endl;
+ *m_logStream<<"************** using SparseEigen ****************"<<std::endl;
*m_logStream<<"*************************************************************"<<std::endl;
}
@@ -2205,17 +2205,17 @@ namespace Trk {
ATH_MSG_DEBUG("running the solving");
// solve
- int info = (*aBetterMat).***REMOVED***Solve(*aBetterVec);
+ int info = (*aBetterMat).SolveWithEigen(*aBetterVec);
if(info == 0) {
- ATH_MSG_INFO("***REMOVED*** solving OK");
+ ATH_MSG_INFO("SolveWithEigen solving OK");
if(m_logStream)
- *m_logStream<<"***REMOVED*** solving OK."<<std::endl;
+ *m_logStream<<"SolveWithEigen solving OK."<<std::endl;
}
else {
- msg(MSG::ERROR)<<"***REMOVED*** error code (0 if OK) = "<<info<<endreq;
+ msg(MSG::ERROR)<<"SolveWithEigen error code (0 if OK) = "<<info<<endreq;
if(m_logStream)
- *m_logStream<<"***REMOVED*** error code (0 if OK) = "<<info<<std::endl;
+ *m_logStream<<"SolveWithEigen error code (0 if OK) = "<<info<<std::endl;
}
if( isCopy )
@@ -2225,11 +2225,11 @@ namespace Trk {
// stop measuring time
clock_t stoptime = clock();
double totaltime = (stoptime-starttime)/double(CLOCKS_PER_SEC);
- ATH_MSG_INFO("Time spent in solve***REMOVED***: "<<totaltime<<" s");
+ ATH_MSG_INFO("Time spent in SolveWithEigen: "<<totaltime<<" s");
ATH_MSG_DEBUG("Alignment constants:");
// compute alignment corrections (translations in mm and rotations in rad)
- // for ***REMOVED*** variances are not calculated
+ // for solveSparseEigen variances are not calculated
for(int i=0; i<m_aNDoF; i++) {
double param = -(*aBetterVec)[i];
@@ -2273,7 +2273,7 @@ namespace Trk {
delete aBetterVec;
// need to do this since success return value from Lapack is 0
- // and from solve***REMOVED***() it is 1
+ // and from SolveWithEigen() it is 1
if (info==0)
info = 1;
diff --git a/Tracking/TrkAlignment/TrkAlignInterfaces/TrkAlignInterfaces/IDerivCalcTool.h b/Tracking/TrkAlignment/TrkAlignInterfaces/TrkAlignInterfaces/IDerivCalcTool.h
index d5bb42b75811f6f0c532f4f39220c763c9bc9713..5d80d02c73f0d4a26aec5c574e0f781093c0232d 100644
--- a/Tracking/TrkAlignment/TrkAlignInterfaces/TrkAlignInterfaces/IDerivCalcTool.h
+++ b/Tracking/TrkAlignment/TrkAlignInterfaces/TrkAlignInterfaces/IDerivCalcTool.h
@@ -41,9 +41,9 @@ namespace Trk {
enum SolveOption {
NONE = 0, //!< not solve in any case (to be used when ipc)
SOLVE = 1, //!< solving after data accumulation (LAPACK)
- SOLVE_FAST = 2, //!< Fast (***REMOVED*** method) solving after data accumulation
+ SOLVE_FAST = 2, //!< Fast (Eigen method) solving after data accumulation
DIRECT_SOLVE = 3, //!< direct solving (LAPACK), already available matrix & vector
- DIRECT_SOLVE_FAST = 4, //!< direct Fast (***REMOVED*** method) solving, already available matrix & vector
+ DIRECT_SOLVE_FAST = 4, //!< direct Fast (Eigen method) solving, already available matrix & vector
DIRECT_SOLVE_CLUSTER = 5 //!< computation of alignment parameters from SCALAPAK already solved matrix
}; // this is also defined in TrkGlobAlign class