ATLAS Offline Software
Loading...
Searching...
No Matches
Trk::MatrixTool Class Reference

#include <MatrixTool.h>

Inheritance diagram for Trk::MatrixTool:
Collaboration diagram for Trk::MatrixTool:

Public Types

enum  SolveOption {
  NONE = 0 , SOLVE = 1 , SOLVE_FAST = 2 , DIRECT_SOLVE = 3 ,
  DIRECT_SOLVE_FAST = 4 , DIRECT_SOLVE_CLUSTER = 5 , SOLVE_ROOT = 6 , SOLVE_CLHEP = 7
}

Public Member Functions

 MatrixTool (const std::string &type, const std::string &name, const IInterface *parent)
 Constructor.
virtual ~MatrixTool ()
 Virtual destructor.
StatusCode initialize ()
 initialize
StatusCode finalize ()
 initialize
StatusCode allocateMatrix (int nDoF=0)
 allocates memory for big matrix and big vector
void prepareBinaryFiles (int solveOption)
 reads/writes matrix entries from/to binary files as necessary
void addFirstDerivatives (AlVec *vector)
 adds first derivative to vector
void addFirstDerivatives (std::list< int, double > &derivatives)
 adds first derivative to vector for only some entries
void addFirstDerivative (int irow, double firstderiv)
void addSecondDerivatives (AlSymMatBase *matrix)
 adds second derivatives to matrix
void addSecondDerivatives (std::list< std::pair< int, int >, double > &derivatives)
 adds first derivative to vector for only some entries
void addSecondDerivative (int irow, int icol, double secondderiv)
bool accumulateFromFiles ()
 accumulates derivates from files.
bool accumulateFromBinaries ()
 accumulates derivates from binary files
int solve ()
 solves for alignment parameters
void storeInTFile (const TString &filename)
 Store Files in a tfile.
bool accumulateFromTFiles ()
 Store Files in a tfile.
void printModuleSolution (std::ostream &os, const AlignModule *module, const CLHEP::HepSymMatrix *cov) const
 namespace { class RestoreIOSFlags { public: RestoreIOSFlags (std::ostream &os) : m_os(&os), m_flags(m_os->flags()), m_precision(m_os->precision()) {} ~RestoreIOSFlags() { m_os->flags(m_flags); m_os->precision(m_precision); } private: std::ostream *m_os; std::ios_base::fmtflags m_flags; std::streamsize m_precision; }; }
void printGlobalSolution (std::ostream &os, const CLHEP::HepSymMatrix *cov)
void printGlobalSolution (std::ostream &os, const TMatrixDSym *cov)
ServiceHandle< StoreGateSvc > & evtStore ()
 The standard StoreGateSvc (event store) Returns (kind of) a pointer to the StoreGateSvc.
const ServiceHandle< StoreGateSvc > & detStore () const
 The standard StoreGateSvc/DetectorStore Returns (kind of) a pointer to the StoreGateSvc.
virtual StatusCode sysInitialize () override
 Perform system initialization for an algorithm.
virtual StatusCode sysStart () override
 Handle START transition.
virtual std::vector< Gaudi::DataHandle * > inputHandles () const override
 Return this algorithm's input handles.
virtual std::vector< Gaudi::DataHandle * > outputHandles () const override
 Return this algorithm's output handles.
Gaudi::Details::PropertyBase & declareProperty (Gaudi::Property< T, V, H > &t)
void updateVHKA (Gaudi::Details::PropertyBase &)
MsgStream & msg () const
bool msgLvl (const MSG::Level lvl) const
void addModule (int alignModuleIndex, int nAlignParam)
int entryNumber (int alignModuleIndex)
void setNHits (int n)
 set module identifier
int nHits () const
void setNTracks (int n)
 set number of tracks
int nTracks () const
void setNMeasurements (int n)
 set number of measurements
int nMeasurements () const
virtual void setLogStream (std::ostream *os)
 sets the output stream for the logfile

Static Public Member Functions

static const InterfaceID & interfaceID ()
 Retrieve interface ID.

Protected Member Functions

void renounceArray (SG::VarHandleKeyArray &handlesArray)
 remove all handles from I/O resolution
std::enable_if_t< std::is_void_v< std::result_of_t< decltype(&T::renounce)(T)> > &&!std::is_base_of_v< SG::VarHandleKeyArray, T > &&std::is_base_of_v< Gaudi::DataHandle, T >, void > renounce (T &h)
void extraDeps_update_handler (Gaudi::Details::PropertyBase &ExtraDeps)
 Add StoreName to extra input/output deps as needed.

Protected Attributes

std::ostream * m_logStream
 logfile output stream
int m_nHits
int m_nTracks
int m_nMeasurements

Private Types

typedef ServiceHandle< StoreGateSvcStoreGateSvc_t

Private Member Functions

int solveROOT ()
int solveCLHEP ()
int solveLapack ()
int solveSparseEigen ()
int solveLocal ()
StatusCode spuriousRemoval ()
void postSolvingLapack (AlVec *dChi2, AlSymMat *d2Chi2, AlVec &w, AlMat &z, int size)
void writeHitmap ()
void readHitmaps ()
Gaudi::Details::PropertyBase & declareGaudiProperty (Gaudi::Property< T, V, H > &hndl, const SG::VarHandleKeyType &)
 specialization for handling Gaudi::Property<SG::VarHandleKey>

Static Private Member Functions

static int fillVecMods ()

Private Attributes

PublicToolHandle< IAlignModuleToolm_alignModuleTool
AlSymMatBasem_bigmatrix = nullptr
 matrix to contain second derivative terms to be used for alignment
AlVecm_bigvector = nullptr
 vector to contain first derivative terms to be used for alignment
Gaudi::Property< bool > m_useSparse {this, "UseSparse", false}
 flag to use AlSpaMat for the big matrix (default is AlSymMat)
Gaudi::Property< bool > m_diagonalize
Gaudi::Property< double > m_eigenvaluethreshold {this, "EigenvalueThreshold", 0., "cut on the minimum eigenvalue"}
Gaudi::Property< int > m_solveOption {this, "SolveOption", NONE, "solving option"}
Gaudi::Property< int > m_modcut
Gaudi::Property< int > m_minNumHits
Gaudi::Property< int > m_minNumTrks
Gaudi::Property< float > m_pullcut
Gaudi::Property< float > m_eigenvalueStep
Gaudi::Property< float > m_Align_db_step
Gaudi::Property< bool > m_calDet
Gaudi::Property< bool > m_wSqMatrix
Gaudi::Property< bool > m_writeMat
Gaudi::Property< bool > m_writeMatTxt
Gaudi::Property< bool > m_writeEigenMat
Gaudi::Property< bool > m_writeEigenMatTxt
Gaudi::Property< bool > m_writeModuleNames
Gaudi::Property< bool > m_writeHitmap
Gaudi::Property< bool > m_writeHitmapTxt
Gaudi::Property< bool > m_readHitmaps
Gaudi::Property< bool > m_writeTFile
Gaudi::Property< bool > m_readTFiles
Gaudi::Property< bool > m_runLocal
double m_scale = -1.
 scale for big matrix and vector normalization
Gaudi::Property< bool > m_scaleMatrix
Gaudi::Property< double > m_softEigenmodeCut
Gaudi::Property< double > m_removeSpurious
Gaudi::Property< double > m_calculateFullCovariance
Gaudi::Property< std::string > m_pathbin
Gaudi::Property< std::string > m_pathtxt
Gaudi::Property< std::string > m_prefixName
Gaudi::Property< std::string > m_tfileName {this, "TFileName", "AlignmentTFile.root", "prefix string to filenames"}
Gaudi::Property< std::string > m_scalaMatName {this, "ScalapackMatrixName", "eigenvectors.bin", "Scalapack matrix name"}
Gaudi::Property< std::string > m_scalaVecName {this, "ScalapackVectorName", "eigenvalues.bin", "Scalapack vector name"}
Gaudi::Property< std::vector< std::string > > m_inputMatrixFiles
Gaudi::Property< std::vector< std::string > > m_inputVectorFiles
Gaudi::Property< std::vector< std::string > > m_inputHitmapFiles
Gaudi::Property< std::vector< std::string > > m_inputTFiles
std::vector< int > m_activeIndices {}
 vector of indices which pass the min-hits cut
int m_aNDoF = 0
 number of active DoF (size of m_activeIndices)
Gaudi::Property< int > m_maxReadErrors
Gaudi::Property< bool > m_AlignIBLbutNotPixel {this, "AlignIBLbutNotPixel", false}
Gaudi::Property< bool > m_AlignPixelbutNotIBL {this, "AlignPixelbutNotIBL", false}
Gaudi::Property< bool > m_DeactivateSCT_ECA_LastDisk {this, "DeactivateSCT_ECA_LastDisk", false}
Gaudi::Property< bool > m_Remove_Pixel_Tx {this, "Remove_Pixel_Tx", false}
Gaudi::Property< bool > m_Remove_Pixel_Ty {this, "Remove_Pixel_Ty", false}
Gaudi::Property< bool > m_Remove_Pixel_Tz {this, "Remove_Pixel_Tz", false}
Gaudi::Property< bool > m_Remove_Pixel_Rx {this, "Remove_Pixel_Rx", false}
Gaudi::Property< bool > m_Remove_Pixel_Ry {this, "Remove_Pixel_Ry", false}
Gaudi::Property< bool > m_Remove_Pixel_Rz {this, "Remove_Pixel_Rz", false}
Gaudi::Property< bool > m_Remove_IBL_Tx {this, "Remove_IBL_Tx", false}
Gaudi::Property< bool > m_Remove_IBL_Ty {this, "Remove_IBL_Ty", false}
Gaudi::Property< bool > m_Remove_IBL_Tz {this, "Remove_IBL_Tz", false}
Gaudi::Property< bool > m_Remove_IBL_Rx {this, "Remove_IBL_Rx", false}
Gaudi::Property< bool > m_Remove_IBL_Ry {this, "Remove_IBL_Ry", false}
Gaudi::Property< bool > m_Remove_IBL_Rz {this, "Remove_IBL_Rz", false}
StoreGateSvc_t m_evtStore
 Pointer to StoreGate (event store by default)
StoreGateSvc_t m_detStore
 Pointer to StoreGate (detector store by default)
std::vector< SG::VarHandleKeyArray * > m_vhka
bool m_varHandleArraysDeclared
std::map< int, int > m_alignModuleMap
int m_nentries

Detailed Description

Definition at line 55 of file MatrixTool.h.

Member Typedef Documentation

◆ StoreGateSvc_t

typedef ServiceHandle<StoreGateSvc> AthCommonDataStore< AthCommonMsg< AlgTool > >::StoreGateSvc_t
privateinherited

Definition at line 388 of file AthCommonDataStore.h.

Member Enumeration Documentation

◆ SolveOption

Enumerator
NONE 

not solve in any case (to be used when ipc)

SOLVE 

solving after data accumulation (LAPACK)

SOLVE_FAST 

Fast (Eigen method) solving after data accumulation.

DIRECT_SOLVE 

direct solving (LAPACK), already available matrix & vector

DIRECT_SOLVE_FAST 

direct Fast (Eigen method) solving, already available matrix & vector

DIRECT_SOLVE_CLUSTER 

computation of alignment parameters from SCALAPAK already solved matrix

SOLVE_ROOT 

computation using ROOT

SOLVE_CLHEP 

computation using CLHEP

Definition at line 58 of file MatrixTool.h.

58 {
59 NONE = 0,
60 SOLVE = 1,
61 SOLVE_FAST = 2,
62 DIRECT_SOLVE = 3,
65 SOLVE_ROOT = 6,
66 SOLVE_CLHEP = 7
67 };
@ SOLVE_FAST
Fast (Eigen method) solving after data accumulation.
Definition MatrixTool.h:61
@ SOLVE_ROOT
computation using ROOT
Definition MatrixTool.h:65
@ SOLVE
solving after data accumulation (LAPACK)
Definition MatrixTool.h:60
@ DIRECT_SOLVE_FAST
direct Fast (Eigen method) solving, already available matrix & vector
Definition MatrixTool.h:63
@ DIRECT_SOLVE
direct solving (LAPACK), already available matrix & vector
Definition MatrixTool.h:62
@ DIRECT_SOLVE_CLUSTER
computation of alignment parameters from SCALAPAK already solved matrix
Definition MatrixTool.h:64
@ NONE
not solve in any case (to be used when ipc)
Definition MatrixTool.h:59
@ SOLVE_CLHEP
computation using CLHEP
Definition MatrixTool.h:66

Constructor & Destructor Documentation

◆ MatrixTool()

MatrixTool::MatrixTool ( const std::string & type,
const std::string & name,
const IInterface * parent )

Constructor.

Definition at line 58 of file MatrixTool.cxx.

60 : IMatrixTool()
61 , AthAlgTool(type,name,parent)
62 {
63 declareInterface<IMatrixTool>(this);
64 }
AthAlgTool()
Default constructor:
IMatrixTool()
constructor

◆ ~MatrixTool()

MatrixTool::~MatrixTool ( )
virtual

Virtual destructor.

Definition at line 67 of file MatrixTool.cxx.

68 {
69 delete m_bigmatrix;
70 delete m_bigvector;
71 }
AlVec * m_bigvector
vector to contain first derivative terms to be used for alignment
Definition MatrixTool.h:149
AlSymMatBase * m_bigmatrix
matrix to contain second derivative terms to be used for alignment
Definition MatrixTool.h:146

Member Function Documentation

◆ accumulateFromBinaries()

bool MatrixTool::accumulateFromBinaries ( )

accumulates derivates from binary files

Definition at line 635 of file MatrixTool.cxx.

636 {
637
638
639 DataVector<AlignPar>* alignParList = m_alignModuleTool->alignParList1D();
640 int nDoF=alignParList->size();
641
642 std::map<int,unsigned long long> modIndexMap;
643 float dummyVersion(0.);
644 double totalscale=0.;
645 for (int ivec=0;ivec<(int)m_inputVectorFiles.size();ivec++) {
646
647 ATH_MSG_DEBUG("Reading vector "<<ivec<<" from file "<<m_inputVectorFiles[ivec]);
648
649 AlVec newVector(nDoF);
650 std::map<int,unsigned long long> newModIndexMap;
651 newVector.SetPathBin(m_pathbin.value()+m_prefixName.value());
652 newVector.SetPathTxt(m_pathtxt.value()+m_prefixName.value());
653 double scale=0;
654 StatusCode sc = newVector.ReadPartial(m_inputVectorFiles[ivec],scale,newModIndexMap,dummyVersion);
655 totalscale += scale;
656 if (sc==StatusCode::FAILURE) {
657 msg(MSG::FATAL)<<"Problem reading vector from "<<m_inputVectorFiles[ivec]<<endmsg;
658 return false;
659 }
660 if (newVector.size()!=m_bigvector->size()) {
661 msg(MSG::FATAL) <<"vector wrong size! newVector size "<<newVector.size()
662 <<", bigvector size "<<m_bigvector->size()<<endmsg;
663 return false;
664 }
665
666 // check modIndexMaps to make sure they are the same
667 if (ivec==0)
668 modIndexMap = newModIndexMap;
669 else if (modIndexMap!=newModIndexMap) {
670 msg(MSG::FATAL)<<"module index maps don't agree!"<<endmsg;
671 return false;
672 }
673 if (ivec>0)
674 *m_bigvector += newVector;
675 else
676 *m_bigvector = newVector;
677 }
678
679 m_scale = totalscale;
680
681 AlSymMat * symBigMatrix=dynamic_cast<AlSymMat*>(m_bigmatrix);
682 AlSpaMat * spaBigMatrix=dynamic_cast<AlSpaMat*>(m_bigmatrix);
683
684
685 for (int imat=0;imat<(int)m_inputMatrixFiles.size();imat++) {
686 ATH_MSG_DEBUG("Reading matrix "<<imat<<" from file "<<m_inputMatrixFiles[imat]);
687
688 // create new matrix to read data from current file
689 int nDoF=modIndexMap.size();
690 bool triang;
692 if (symBigMatrix) {
693 AlSymMat newMatrix(nDoF);
694 sc = newMatrix.Read(m_inputMatrixFiles[imat],nDoF,triang,dummyVersion);
695 if (sc==StatusCode::SUCCESS)
696 *symBigMatrix += newMatrix;
697 }
698 else {
699 if (!spaBigMatrix) {
700 throw std::logic_error("Unhandled matrix type");
701 }
702
703 AlSpaMat newMatrix(nDoF);
704 sc = newMatrix.Read(m_inputMatrixFiles[imat],nDoF,triang,dummyVersion);
705
706 if (sc==StatusCode::SUCCESS) {
707 if (imat>0)
708 *spaBigMatrix += newMatrix;
709 else
710 *spaBigMatrix = newMatrix;
711 }
712 }
713
714 if (sc==StatusCode::FAILURE) {
715 msg(MSG::FATAL)<<"problem reading matrix from "<<m_inputMatrixFiles[imat]<<endmsg;
716 return false;
717 }
718
719 if (!m_useSparse && triang==m_wSqMatrix) {
720 ATH_MSG_WARNING("matrix not expected format! Changing m_wSqMatrix to "<<!triang);
721 m_wSqMatrix=!triang;
722 }
723
724 }
725
726 // accumulate hitmap from hitmap files
727 if(m_readHitmaps)
728 readHitmaps();
729
730 return true;
731 }
#define endmsg
#define ATH_MSG_WARNING(x)
#define ATH_MSG_DEBUG(x)
static Double_t sc
MsgStream & msg() const
size_type size() const noexcept
Returns the number of elements in the collection.
Gaudi::Property< std::string > m_prefixName
Definition MatrixTool.h:222
Gaudi::Property< bool > m_readHitmaps
Definition MatrixTool.h:193
Gaudi::Property< bool > m_wSqMatrix
Definition MatrixTool.h:176
double m_scale
scale for big matrix and vector normalization
Definition MatrixTool.h:204
Gaudi::Property< std::vector< std::string > > m_inputVectorFiles
Definition MatrixTool.h:237
Gaudi::Property< std::vector< std::string > > m_inputMatrixFiles
Definition MatrixTool.h:234
PublicToolHandle< IAlignModuleTool > m_alignModuleTool
Definition MatrixTool.h:142
Gaudi::Property< std::string > m_pathbin
Definition MatrixTool.h:218
Gaudi::Property< bool > m_useSparse
flag to use AlSpaMat for the big matrix (default is AlSymMat)
Definition MatrixTool.h:152
Gaudi::Property< std::string > m_pathtxt
Definition MatrixTool.h:220
vec_fb< typename boost::int_t< sizeof(T) *8 >::exact, N > ivec
Definition vec_fb.h:53
::StatusCode StatusCode
StatusCode definition for legacy code.

◆ accumulateFromFiles()

bool MatrixTool::accumulateFromFiles ( )
virtual

accumulates derivates from files.

Flag decides if it is binary or TFiles

Implements Trk::IMatrixTool.

Definition at line 622 of file MatrixTool.cxx.

623 {
624
625 if(m_readTFiles){
626 ATH_MSG_INFO("Info to obtained from from TFiles");
627 return accumulateFromTFiles();
628 }else{
629 ATH_MSG_INFO("Info to obtained from from Binary files");
630 return accumulateFromBinaries();
631 }
632 }
#define ATH_MSG_INFO(x)
bool accumulateFromTFiles()
Store Files in a tfile.
bool accumulateFromBinaries()
accumulates derivates from binary files
Gaudi::Property< bool > m_readTFiles
Definition MatrixTool.h:198

◆ accumulateFromTFiles()

bool MatrixTool::accumulateFromTFiles ( )

Store Files in a tfile.

Definition at line 829 of file MatrixTool.cxx.

830 {
831 DataVector<AlignPar>* alignParList = m_alignModuleTool->alignParList1D();
832 int nDoF=alignParList->size();
833 ATH_MSG_DEBUG("OPENING TFILES");
834
835 std::map<int,unsigned long long> modIndexMap;
836 std::map<int,unsigned long long> DoFMap;
837 double totalscale=0.;
838
839 AlSymMat * symBigMatrix=dynamic_cast<AlSymMat*>(m_bigmatrix);
840 AlSpaMat * spaBigMatrix=dynamic_cast<AlSpaMat*>(m_bigmatrix);
841 //TMatrixDSparse *accumMatrix(0);
842 AlSpaMat *accumMatrix = nullptr;
843
844 const AlignModuleList * moduleList = m_alignModuleTool->alignModules1D();
845 int nModules = moduleList->size();
846
847 TVectorD TotalHits(nModules);
848 TVectorD TotalTracks(nModules);
849
850 int numberOfReadErrors = 0;
851
852 struct rusage myusage{};
853 int itworked = getrusage(RUSAGE_SELF,&myusage);
854 if(itworked)//note: rusage returns zero if it succeeds!
855 ATH_MSG_DEBUG("ItWorked");
856
857 long intialMemUse = myusage.ru_maxrss;
858
859 for (int ifile = 0; ifile < (int)m_inputTFiles.size(); ifile++) {
860 if (numberOfReadErrors > m_maxReadErrors){
861 msg(MSG::FATAL) << " number of errors when reading the TFiles already exceed " << m_maxReadErrors << endmsg;
862 return false;
863 }
864
865 ATH_MSG_DEBUG("Reading File number " << ifile << ", " << m_inputTFiles[ifile]);
866
867 itworked = getrusage(RUSAGE_SELF,&myusage);
868 ATH_MSG_DEBUG("Memory usage [MB], total " << myusage.ru_maxrss/1024 << ", increase " << (myusage.ru_maxrss-intialMemUse)/1024);
869
870 TFile* myFile = TFile::Open(m_inputTFiles[ifile].c_str());
871
872 if ( myFile->IsZombie() || !(myFile->IsOpen()) ) {
873 ++numberOfReadErrors;
874 ATH_MSG_ERROR( " Problem reading TFile " << m_inputTFiles[ifile] );
875 continue;
876 }
877
878 std::map<int,unsigned long long> newModIndexMap;
879
880 TVectorD* myModuleIDs;
881 myModuleIDs = (TVectorD*)myFile->Get("ModuleID");
882 if( !myModuleIDs ){
883 ++numberOfReadErrors;
884 ATH_MSG_ERROR("Modules ID not read!!!");
885 continue;
886 }
887
888 for (int i(0); i<myModuleIDs->GetNrows(); ++i){
889 //Coverting back from a double to a unvi signed long long
890 double source = (*myModuleIDs)(i);
892 memcpy(&target, &source, sizeof(target));
893 newModIndexMap[i]=target;
894 //std::cout << i<< " " <<target <<std::endl;
895 }
896
897 delete myModuleIDs;
898
899 std::map<int,unsigned long long> newDoFMap;
900
901 TVectorD* myDoFs;
902 myDoFs = (TVectorD*)myFile->Get("dof");
903 if( !myDoFs ){
904 ++numberOfReadErrors;
905 ATH_MSG_ERROR("DoFs not read!!!");
906 continue;
907 }
908
909 for (int i(0); i<myDoFs->GetNrows(); ++i){
910 //Coverting back from a double to a unsigned long long
911 double source = (*myDoFs)(i);
913 memcpy(&target, &source, sizeof(target));
914 newDoFMap[i]=target;
915 }
916 delete myDoFs;
917
918
919 TVectorD* Scale;
920 Scale = (TVectorD*)myFile->Get("Scale");
921 if( !Scale ){
922 ++numberOfReadErrors;
923 ATH_MSG_ERROR("Scale not read!!!");
924 continue;
925 }
926
927 double scale=(*Scale)(0);
928 totalscale += scale;
929 delete Scale;
930
931
932 ATH_MSG_DEBUG("Reading Vector");
933 TVectorD* vector = (TVectorD*)myFile->Get("Vector");
934 if( !vector ){
935 ++numberOfReadErrors;
936 ATH_MSG_ERROR("Vector not read!!!");
937 continue;
938 }
939
940 AlVec* newVector = new AlVec(nDoF);
941 newVector->SetPathBin(m_pathbin.value()+m_prefixName.value());
942 newVector->SetPathTxt(m_pathtxt.value()+m_prefixName.value());
943
944 if (newVector->size() != m_bigvector->size() ) {
945 msg(MSG::FATAL) << "vector wrong size! newVector size " << newVector->size()
946 << ", bigvector size " << m_bigvector->size()<<endmsg;
947 delete newVector;
948 delete vector;
949 return false;
950 }
951
952 if (m_bigvector->size() != vector->GetNrows() ) {
953 msg(MSG::FATAL) << "File vector wrong size! File Vector size " << vector->GetNrows()
954 << ", bigvector size " << m_bigvector->size()<<endmsg;
955 delete newVector;
956 delete vector;
957 return false;
958 }
959
960
961 for (int i=0;i<nDoF;i++) {
962 (*newVector)[i] = (*vector)(i);
963 }
964 delete vector;
965
966 // check modIndexMaps to make sure they are the same
967 if (ifile == 0){
968 DoFMap = newDoFMap;
969 } else if (DoFMap!=newDoFMap) {
970 msg(MSG::FATAL) << "module dofs don't agree!" << endmsg;
971 return false;
972 }
973
974 if (ifile == 0){
975 modIndexMap = newModIndexMap;
976 } else if (modIndexMap!=newModIndexMap) {
977 msg(MSG::FATAL) << "module index maps don't agree!" << endmsg;
978 return false;
979 }
980
981 if (ifile>0){
982 *m_bigvector += *newVector;
983 delete newVector;
984 } else {
985 delete m_bigvector;
986 m_bigvector = newVector;
987 }
988
989
990 ATH_MSG_DEBUG("Reading matrix ");
991 TMatrixDSparse* matrix = (TMatrixDSparse*)myFile->Get("Matrix");
992
993 if( !matrix ){
994 ++numberOfReadErrors;
995 ATH_MSG_ERROR("Matrix not read!!!");
996 continue;
997 }
998
999
1000 if (ifile == 0 ){
1001
1002 accumMatrix = new AlSpaMat(nDoF);
1003 ATH_MSG_DEBUG("Matrix size b4 "<< accumMatrix->ptrMap()->size() );
1004
1005 //This method is ok for large matrix files... really only access the non zero elements
1006 for (int ii=0;ii<nDoF;ii++) {
1007 const TMatrixTSparseRow_const<double> myRow = (*matrix)[ii];
1008 int i = myRow.GetRowIndex();
1009 for (int jj=0;jj<myRow.GetNindex();jj++) {
1010 int j = (myRow.GetColPtr())[jj];
1011 const double myElement= (myRow.GetDataPtr())[jj];
1012 if (i<j){
1013 ATH_MSG_DEBUG("i < j " );
1014 j = i;
1015 i = (myRow.GetColPtr())[jj];
1016 }
1017 (*accumMatrix)[i][j] = myElement;
1018 }
1019 }
1020 ATH_MSG_DEBUG("Matrix size AF "<< accumMatrix->ptrMap()->size() );
1021
1022 } else if ( accumMatrix) {
1023 ATH_MSG_DEBUG("Matrix size b4 "<< accumMatrix->ptrMap()->size() );
1024
1025 for (int ii=0;ii<nDoF;ii++) {
1026 const TMatrixTSparseRow_const<double> myRow = (*matrix)[ii];
1027 int i = myRow.GetRowIndex();
1028 for (int jj=0;jj<myRow.GetNindex();jj++) {
1029 int j = (myRow.GetColPtr())[jj];
1030 const double myElement= (myRow.GetDataPtr())[jj];
1031 if (i<j){
1032 ATH_MSG_DEBUG("i < j " );
1033 j = i;
1034 i = (myRow.GetColPtr())[jj];
1035 }
1036 (*accumMatrix)[i][j] += myElement;
1037 }
1038 }
1039 ATH_MSG_DEBUG("Matrix size AF "<< accumMatrix->ptrMap()->size() );
1040
1041 } else {
1042 delete matrix;
1043 ++numberOfReadErrors;
1044 ATH_MSG_ERROR("Matrix allocation error!!!");
1045 continue;
1046 }
1047
1048 delete matrix;
1049
1050 TVectorD* hits;
1051 TVectorD* tracks;
1052
1053 ATH_MSG_DEBUG("Reading hitmap ");
1054 hits = (TVectorD*)myFile->Get("Hits");
1055 if( !hits ){
1056 ++numberOfReadErrors;
1057 ATH_MSG_ERROR("Hitmap 1 not read!!!");
1058 continue;
1059 }
1060
1061 tracks = (TVectorD*)myFile->Get("Tracks");
1062 if( !tracks ){
1063 delete hits;
1064 ++numberOfReadErrors;
1065 ATH_MSG_ERROR("Hitmap 2 not read!!!");
1066 continue;
1067 }
1068
1069 if(hits->GetNrows() != TotalHits.GetNrows() ){
1070 delete hits;
1071 delete tracks;
1072 ++numberOfReadErrors;
1073 ATH_MSG_ERROR("Hitmap size incorrect!!!");
1074 continue;
1075 }
1076
1077 TotalHits += (*hits);
1078 TotalTracks += (*tracks);
1079
1080 delete hits;
1081 delete tracks;
1082
1083 myFile->Close("R");
1084 delete myFile;
1085
1086 itworked = getrusage(RUSAGE_SELF,&myusage);
1087 ATH_MSG_DEBUG("Memory usage [MB], total " << myusage.ru_maxrss/1024 << ", increase " << (myusage.ru_maxrss-intialMemUse)/1024);
1088
1089 }
1090
1091
1092
1093 // create new matrix to read data from current file
1094 if(accumMatrix){
1095 if (symBigMatrix) {
1096 AlSymMat newMatrix(nDoF);
1097 //This method is ok for small matrix files
1098 for (int i=0;i<nDoF;i++) {
1099 for (int j=0;j<=i;j++) {
1100 newMatrix[i][j] = (*accumMatrix)[i][j];
1101 }
1102 }
1103
1104 *symBigMatrix += newMatrix;
1105 delete accumMatrix;
1106 } else if (spaBigMatrix) {
1107 ATH_MSG_DEBUG( "should reassign matrix "<< spaBigMatrix->ptrMap()->size() );
1108 *spaBigMatrix += *accumMatrix;
1109 ATH_MSG_DEBUG( "?????? "<< spaBigMatrix->ptrMap()->size() );
1110 delete accumMatrix;
1111 }
1112 }
1113
1114 ATH_MSG_DEBUG( "?????? "<< m_bigmatrix->ptrMap()->size() );
1115
1116 AlignModuleList::const_iterator imod = moduleList->begin();
1117 AlignModuleList::const_iterator imod_end = moduleList->end();
1118 int index = 0;
1119 int totalhits = 0;
1120 for(; imod != imod_end; ++imod, ++index ) {
1121 AlignModule * module = *imod;
1122 module->setNHits((int)TotalHits(index));
1123 module->setNTracks((int)TotalTracks(index));
1124 totalhits += (int)TotalHits(index);
1125 }
1126
1127
1128 m_nHits = totalhits;
1129 m_nTracks = 0;
1130 m_nMeasurements = 0;
1131 m_scale = totalscale;
1132
1133 return true;
1134 }
class TMatrixTSparse< double > TMatrixDSparse
Definition AlSpaMat.h:14
#define ATH_MSG_ERROR(x)
Gaudi::Property< int > m_maxReadErrors
Definition MatrixTool.h:251
Gaudi::Property< std::vector< std::string > > m_inputTFiles
Definition MatrixTool.h:245
void Scale(TH1 *h, double d=1)
str index
Definition DeMoScan.py:362
std::vector< AlignModule * > AlignModuleList

◆ addFirstDerivative()

void MatrixTool::addFirstDerivative ( int irow,
double firstderiv )
virtual

Implements Trk::IMatrixTool.

Definition at line 1481 of file MatrixTool.cxx.

1482 {
1483 (*m_bigvector)[irow] += firstderiv;
1484 }

◆ addFirstDerivatives() [1/2]

void MatrixTool::addFirstDerivatives ( AlVec * vector)
virtual

adds first derivative to vector

Implements Trk::IMatrixTool.

Definition at line 1461 of file MatrixTool.cxx.

1462 {
1463 }

◆ addFirstDerivatives() [2/2]

void MatrixTool::addFirstDerivatives ( std::list< int, double > & derivatives)
virtual

adds first derivative to vector for only some entries

Implements Trk::IMatrixTool.

Definition at line 1471 of file MatrixTool.cxx.

1472 {
1473 }

◆ addModule()

void IMatrixTool::addModule ( int alignModuleIndex,
int nAlignParam )
inlineinherited

Definition at line 116 of file IMatrixTool.h.

116{ m_alignModuleMap[alignModuleIndex]=m_nentries; m_nentries += nAlignParam; }
std::map< int, int > m_alignModuleMap

◆ addSecondDerivative()

void MatrixTool::addSecondDerivative ( int irow,
int icol,
double secondderiv )
virtual

Implements Trk::IMatrixTool.

Definition at line 1487 of file MatrixTool.cxx.

1488 {
1489 (*m_bigmatrix)[irow][icol] += secondderiv;
1490 }

◆ addSecondDerivatives() [1/2]

void MatrixTool::addSecondDerivatives ( AlSymMatBase * matrix)
virtual

adds second derivatives to matrix

Implements Trk::IMatrixTool.

Definition at line 1466 of file MatrixTool.cxx.

1467 {
1468 }

◆ addSecondDerivatives() [2/2]

void MatrixTool::addSecondDerivatives ( std::list< std::pair< int, int >, double > & derivatives)
virtual

adds first derivative to vector for only some entries

Implements Trk::IMatrixTool.

Definition at line 1476 of file MatrixTool.cxx.

1477 {
1478 }

◆ allocateMatrix()

StatusCode MatrixTool::allocateMatrix ( int nDoF = 0)
virtual

allocates memory for big matrix and big vector

Implements Trk::IMatrixTool.

Definition at line 104 of file MatrixTool.cxx.

105 {
106 ATH_MSG_INFO("allocating matrix and vector with nDoF = "<<nDoF);
107
108 if (nullptr!=m_bigmatrix || nullptr!=m_bigvector)
109 ATH_MSG_ERROR("big matrix already allocated!");
110
111 // Decide upon the big matrix representation:
112 if( m_useSparse )
113 m_bigmatrix = new AlSpaMat(nDoF);
114 else
115 m_bigmatrix = new AlSymMat(nDoF);
116
117 m_bigvector = new AlVec(nDoF);
118
119 ATH_MSG_INFO(" After Matrix and Vector allocation");
120
121 // set paths for matrix and vector output
122 m_bigmatrix->SetPathBin(m_pathbin.value()+m_prefixName);
123 m_bigmatrix->SetPathTxt(m_pathtxt.value()+m_prefixName);
124 m_bigvector->SetPathBin(m_pathbin.value()+m_prefixName);
125 m_bigvector->SetPathTxt(m_pathtxt.value()+m_prefixName);
126
127 ATH_MSG_INFO("set path to "<<m_pathbin.value()+m_prefixName.value());
128 return StatusCode::SUCCESS;
129 }

◆ declareGaudiProperty()

Gaudi::Details::PropertyBase & AthCommonDataStore< AthCommonMsg< AlgTool > >::declareGaudiProperty ( Gaudi::Property< T, V, H > & hndl,
const SG::VarHandleKeyType &  )
inlineprivateinherited

specialization for handling Gaudi::Property<SG::VarHandleKey>

Definition at line 156 of file AthCommonDataStore.h.

158 {
160 hndl.value(),
161 hndl.documentation());
162
163 }
Gaudi::Details::PropertyBase & declareProperty(Gaudi::Property< T, V, H > &t)

◆ declareProperty()

Gaudi::Details::PropertyBase & AthCommonDataStore< AthCommonMsg< AlgTool > >::declareProperty ( Gaudi::Property< T, V, H > & t)
inlineinherited

Definition at line 145 of file AthCommonDataStore.h.

145 {
146 typedef typename SG::HandleClassifier<T>::type htype;
148 }
Gaudi::Details::PropertyBase & declareGaudiProperty(Gaudi::Property< T, V, H > &hndl, const SG::VarHandleKeyType &)
specialization for handling Gaudi::Property<SG::VarHandleKey>

◆ detStore()

const ServiceHandle< StoreGateSvc > & AthCommonDataStore< AthCommonMsg< AlgTool > >::detStore ( ) const
inlineinherited

The standard StoreGateSvc/DetectorStore Returns (kind of) a pointer to the StoreGateSvc.

Definition at line 95 of file AthCommonDataStore.h.

◆ entryNumber()

int IMatrixTool::entryNumber ( int alignModuleIndex)
inlineinherited

Definition at line 118 of file IMatrixTool.h.

118{ if ( m_alignModuleMap.find(alignModuleIndex) == m_alignModuleMap.end()) return -1; else return m_alignModuleMap[alignModuleIndex]; }

◆ evtStore()

ServiceHandle< StoreGateSvc > & AthCommonDataStore< AthCommonMsg< AlgTool > >::evtStore ( )
inlineinherited

The standard StoreGateSvc (event store) Returns (kind of) a pointer to the StoreGateSvc.

Definition at line 85 of file AthCommonDataStore.h.

◆ extraDeps_update_handler()

void AthCommonDataStore< AthCommonMsg< AlgTool > >::extraDeps_update_handler ( Gaudi::Details::PropertyBase & ExtraDeps)
protectedinherited

Add StoreName to extra input/output deps as needed.

use the logic of the VarHandleKey to parse the DataObjID keys supplied via the ExtraInputs and ExtraOuputs Properties to add the StoreName if it's not explicitly given

◆ fillVecMods()

int MatrixTool::fillVecMods ( )
staticprivate

Definition at line 1631 of file MatrixTool.cxx.

1632 {
1633 return 0;
1634 }

◆ finalize()

StatusCode MatrixTool::finalize ( )

initialize

Definition at line 95 of file MatrixTool.cxx.

96 {
97 ATH_MSG_DEBUG("finalize() of MatrixTool");
98
99 return StatusCode::SUCCESS;
100 }

◆ initialize()

StatusCode MatrixTool::initialize ( )

initialize

Definition at line 74 of file MatrixTool.cxx.

75 {
76 ATH_MSG_DEBUG("initialize() of MatrixTool");
77
78 // get AlignModuleTool
79 if (m_alignModuleTool.retrieve().isSuccess())
80 ATH_MSG_INFO("Retrieved " << m_alignModuleTool);
81 else{
82 msg(MSG::FATAL) << "Could not get " << m_alignModuleTool << endmsg;
83 return StatusCode::FAILURE;
84 }
85
86 ATH_MSG_INFO("Retrieving data from the following files: ");
87 for (auto & inputVectorFile : m_inputVectorFiles) {
88 ATH_MSG_INFO(m_pathbin+inputVectorFile);
89 }
90
91 return StatusCode::SUCCESS;
92 }

◆ inputHandles()

virtual std::vector< Gaudi::DataHandle * > AthCommonDataStore< AthCommonMsg< AlgTool > >::inputHandles ( ) const
overridevirtualinherited

Return this algorithm's input handles.

We override this to include handle instances from key arrays if they have not yet been declared. See comments on updateVHKA.

◆ interfaceID()

const InterfaceID & IMatrixTool::interfaceID ( )
inlinestaticinherited

Retrieve interface ID.

Definition at line 110 of file IMatrixTool.h.

110 {
112 }
static const InterfaceID IID_TRKALIGNINTERFACES_IMatrixTool("IMatrixTool", 1, 0)

◆ msg()

MsgStream & AthCommonMsg< AlgTool >::msg ( ) const
inlineinherited

Definition at line 24 of file AthCommonMsg.h.

24 {
25 return this->msgStream();
26 }

◆ msgLvl()

bool AthCommonMsg< AlgTool >::msgLvl ( const MSG::Level lvl) const
inlineinherited

Definition at line 30 of file AthCommonMsg.h.

30 {
31 return this->msgLevel(lvl);
32 }

◆ nHits()

int Trk::IMatrixTool::nHits ( ) const
inlineinherited

Definition at line 84 of file IMatrixTool.h.

84{ return m_nHits; }

◆ nMeasurements()

int Trk::IMatrixTool::nMeasurements ( ) const
inlineinherited

Definition at line 92 of file IMatrixTool.h.

92{ return m_nMeasurements; }

◆ nTracks()

int Trk::IMatrixTool::nTracks ( ) const
inlineinherited

Definition at line 88 of file IMatrixTool.h.

88{ return m_nTracks; }

◆ outputHandles()

virtual std::vector< Gaudi::DataHandle * > AthCommonDataStore< AthCommonMsg< AlgTool > >::outputHandles ( ) const
overridevirtualinherited

Return this algorithm's output handles.

We override this to include handle instances from key arrays if they have not yet been declared. See comments on updateVHKA.

◆ postSolvingLapack()

void MatrixTool::postSolvingLapack ( AlVec * dChi2,
AlSymMat * d2Chi2,
AlVec & w,
AlMat & z,
int size )
private

Definition at line 1830 of file MatrixTool.cxx.

1831 {
1832 ATH_MSG_DEBUG("in postSolvinglapack()");
1833
1834 if( z.ncol() != size) {
1835 msg(MSG::ERROR)<<"Eigenvector matrix has incorrect size : "<<z.ncol()<<" != "<<size<<endmsg;
1836 return;
1837 }
1838
1839 if( (int)m_activeIndices.size() != size) {
1840 msg(MSG::ERROR)<<"Number of active parameters is incorrect : "<<m_activeIndices.size()<<" != "<<size<<endmsg;
1841 return;
1842 }
1843
1844 // Compute bigvector in diagonal basis (Vb = Ut * bigvector)
1845 AlVec D(size);
1846 D = z*(*dChi2);
1847
1848 if (m_writeEigenMat) {
1849
1850 ATH_MSG_INFO("writing the eigenvectors in a matrix: "<< z.nrow() << "x" << z.ncol());
1851
1852 // Set Path for the z matrix (eigenvector matrix)
1853 z.SetPathBin(m_pathbin.value()+m_prefixName.value());
1854 z.SetPathTxt(m_pathtxt.value()+m_prefixName.value());
1855
1856 ATH_MSG_INFO("writing the eigenvector matrix: "<< m_scalaMatName);
1857 ATH_MSG_DEBUG("matrix will be in: "<< m_pathbin.value()+m_prefixName.value()+m_scalaMatName.value());
1858
1859 StatusCode sc = z.Write(m_scalaMatName, true); // write the eigenvector matrix
1860
1861 if (sc!=StatusCode::SUCCESS)
1862 msg(MSG::ERROR)<<"Problem writing eigenvector matrix"<<endmsg;
1863
1864 // Set Path for the w matrix (eigenvalues matrix - diagonal bigmatrix)
1865 w.SetPathBin(m_pathbin.value()+m_prefixName.value());
1866 w.SetPathTxt(m_pathtxt.value()+m_prefixName.value());
1867
1868 ATH_MSG_INFO("writing the eigenvectors in a vector: "<< w.size());
1869 ATH_MSG_INFO("writing the eigenvalues vector (diagonal bigmatrix): "<< m_scalaVecName);
1870 ATH_MSG_DEBUG("vector will be in: "<< m_pathbin.value()+m_prefixName.value()+m_scalaVecName.value());
1871
1872 sc = w.WriteEigenvalueVec(m_scalaVecName, true); // write the eigenvalues vecor
1873
1874 if (sc!=StatusCode::SUCCESS)
1875 msg(MSG::ERROR)<<"Problem writing eigenvector matrix"<<endmsg;
1876
1877 if (m_writeEigenMatTxt) {
1878 sc = z.Write("eigenvectors.txt", false);
1879 if (sc!=StatusCode::SUCCESS)
1880 msg(MSG::ERROR)<<"Problem writing eigenvector matrix to text file"<<endmsg;
1881 sc = w.WriteEigenvalueVec("eigenvalues.txt", false);
1882 if (sc!=StatusCode::SUCCESS)
1883 msg(MSG::ERROR)<<"Problem writing eigenvalue vector to text file"<<endmsg;
1884 }
1885
1886 }
1887
1888 // Set eigenvalue thresholds
1889 const double eigenvalue_threshold = 1e-19;
1890
1891 // weak mode removal
1892 if (m_modcut == -1) {
1893
1894 ATH_MSG_INFO(" Starting the automatic Weak Mode Removal method");
1895
1896 // create a pull vector for the alignment corrections in diagonal basis (db_pulls)
1897 //int nDoF=m_alignModuleTool->nAlignParameters();
1898 AlVec* Align_db = new AlVec(size);
1899 AlVec* Align_error_db = new AlVec(size);
1900 AlVec* AlignPull = new AlVec(size);
1901 ATH_MSG_DEBUG("AlignPull vector size is: "<< (*AlignPull).size());
1902
1903 m_modcut = 0;
1904 bool wm_stop = false;
1905
1906 // -----------------------------------------------------------------------
1907 // First pass: removing weak modes because of large pull
1908 // compute alignment pulls for corrections in diagonal basis (db)
1909 for(int i=0; i<size; i++) {
1910
1911 (*Align_db)[i] = (-D[i]/w[i]);
1912 (*Align_error_db)[i] = sqrt(1.0/w[i]/m_scale);
1913
1914 if (w[i]<eigenvalue_threshold) {
1915 ATH_MSG_INFO(" + EigenMode " << i
1916 << " removed as eigenvalue lower than the threshold " << eigenvalue_threshold
1917 << ": " << w[i]);
1918 (*AlignPull)[i] = 0.0;
1919 m_modcut++;
1920 }
1921 else
1922 (*AlignPull)[i] = (*Align_db)[i] / (*Align_error_db)[i];
1923
1924 ATH_MSG_DEBUG(i << ". AlignPar: " << (*Align_db)[i] << " +- " << (*Align_error_db)[i]
1925 << " (pull: " << (*AlignPull)[i] << ") ; w[i]: " << w[i]);
1926 }
1927 ATH_MSG_INFO(" +++ Weak Mode removal ++ stop after mode "<< m_modcut << " (First pass)");
1928 // -----------------------------------------------------------------------
1929
1930 // -----------------------------------------------------------------------
1931 // Second pass
1932 // if the error is greater than the correction -> cut this mode
1933 for(int i=m_modcut; (i<size && !wm_stop); i++) {
1934
1935 // if the error is greater than the correction -> cut this mode
1936 if (fabs((*AlignPull)[i])<m_pullcut) {
1937 ATH_MSG_INFO(" + EigenMode " << i
1938 << " removed as pull is lower than " << m_pullcut << ": "
1939 << (*AlignPull)[i]);
1940 m_modcut++;
1941 }
1942 else
1943 wm_stop = true;
1944
1945 }
1946 ATH_MSG_INFO(" +++ Weak Mode removal ++ stop after mode "<< m_modcut << " (Second pass)");
1947 // -----------------------------------------------------------------------
1948
1949 wm_stop = false;
1950
1951 // ----------------------------------------------------------------------
1952 // Third pass
1953 // Check if the next eigenvalues is far away. If it is two orders of
1954 // magnitude bigger remove also this mode and allow the search
1955 for(int i=m_modcut; (i<size && !wm_stop); i++) {
1956
1957 // if the next eigenvalues is far away -> cut this mode
1958 if (m_eigenvalueStep*w[i]<w[i+1]) {
1959 ATH_MSG_INFO(" + EigenMode " << i
1960 << " removed as diff between eigenvalues, " << w[i] << " and " << w[i+1]
1961 << ", is greater than " << m_eigenvalueStep);
1962 m_modcut++;
1963 }
1964 else
1965 wm_stop = true;
1966
1967 }
1968 ATH_MSG_INFO(" +++ Weak Mode removal ++ stop after mode "<< m_modcut << " (Third pass)");
1969 // -----------------------------------------------------------------------
1970
1971 wm_stop = false;
1972
1973 // -----------------------------------------------------------------------
1974 // Fourth pass
1975 // Check if the next eigenvalues is far away. If it is two orders of
1976 // magnitude bigger remove also this mode and allow the search
1977 for(int i=m_modcut; (i<size && !wm_stop); i++) {
1978
1979 // if the next eigenvalues is far away -> cut this mode
1980 if ( fabs((*Align_db)[i]) > m_Align_db_step*fabs((*Align_db)[i+1]) ) {
1981 ATH_MSG_INFO(" + EigenMode " << i
1982 << " removed as diff between corrections, " << w[i] << " and " << w[i+1]
1983 << ", is greater than "
1984 << m_Align_db_step);
1985 m_modcut++;
1986 }
1987 else
1988 wm_stop = true;
1989
1990 }
1991 ATH_MSG_INFO(" +++ Weak Mode removal ++ stop after mode "<< m_modcut << " (Fourth pass)");
1992 // -----------------------------------------------------------------------------------------------------
1993
1994 // Free memory and clear the pointer to
1995 // prevent using invalid memory reference
1996 delete Align_db;
1997 delete Align_error_db;
1998 delete AlignPull;
1999 Align_db = nullptr;
2000 Align_error_db = nullptr;
2001 AlignPull = nullptr;
2002
2003 } // end of if(m_modcut == -1)
2004
2005 // Save some stuff to debug purposes
2006 /*
2007 if (m_storeDia) {
2008 std::string path = m_pathtxt+m_prefixName+"align_dia";
2009 std::fstream orthogon(path.c_str(), std::ios::out);
2010 orthogon.setf(std::ios::fixed);
2011 orthogon.setf(std::ios::showpoint);
2012 orthogon.precision(6);
2013
2014 orthogon << std::setw(10)
2015 << "--------------------------------------------------------------------------------"
2016 << std::endl;
2017 orthogon << std::setw(10) << " ModeCut = " << m_modcut << std::endl;
2018 orthogon << std::setw(10)
2019 << "--------------------------------------------------------------------------------"
2020 << std::endl;
2021
2022 orthogon << std::setw(10) << "mode"
2023 << std::setw(20) << "eigenmode"
2024 << std::setw(20) << "eigenmode error"
2025 << std::setw(25) << "eigenvalue"
2026 << std::endl;
2027
2028 for( int m=0; m<size; m++) {
2029
2030 // mode
2031 orthogon << std::setw(10) << m;
2032
2033 // eigenmode (db)
2034 if( w[m]>1.0e-15) orthogon << std::setw(20) << -D[m]/w[m];
2035 else orthogon << std::setw(20) << 0.0;
2036
2037 // error eigenmode (error_db)
2038 if( w[m]>1.0e-15) orthogon << std::setw(20) << sqrt(1.0/w[m]/m_scale);
2039 else orthogon << std::setw(20) << 0.0;
2040
2041 // eigenvalues
2042 orthogon << std::setw(25) << w[m] << std::endl;
2043 }
2044 orthogon.close();
2045 } // end store align_dia.txt
2046 */
2047
2048 AlVec delta(size);
2049 AlVec deltafull(size);
2050 AlVec errSq(size);
2051
2052 // full covariance matrix
2053 CLHEP::HepSymMatrix * cov = nullptr;
2055 // Warning ! The matrix can be huge!
2056 // This can lead to memory problems
2057 cov = new CLHEP::HepSymMatrix(size,0);
2058
2059 if(m_logStream)
2060 *m_logStream<<"/------ The Eigenvalue Spectrum -------"<<std::endl;
2061
2062 for (int i=0;i<size;i++) {
2063 AlVec thisdelta(size);
2064 for(int j=0;j<size;j++)
2065 thisdelta[j] = z[i][j] * (-D[i]/w[i]);
2066 deltafull += thisdelta;
2067
2068 ATH_MSG_DEBUG("eigenvalue "<<w[i]);
2069 if( i<m_modcut ) {
2070 ATH_MSG_INFO("skipping eigenvalue "<<w[i]<<" , modcut is "<<m_modcut);
2071 if(m_logStream)
2072 *m_logStream<<"| skipping eigenvalue "<<w[i]<<std::endl;
2073 }
2074 else if( w[i] < m_eigenvaluethreshold ) {
2075 ATH_MSG_INFO("skipping eigenvalue "<<w[i]<<" , cut is "<<m_eigenvaluethreshold);
2076 if(m_logStream)
2077 *m_logStream<<"| skipping eigenvalue "<<w[i]<<std::endl;
2078 }
2079 else {
2080 if(m_logStream)
2081 *m_logStream<<"| "<<w[i]<<std::endl;
2082
2083 delta += thisdelta;
2084 for(int j=0;j<size;j++) {
2085 errSq[j] += z[i][j] * z[i][j] / w[i];
2087 for(int k=0;k<=j;k++)
2088 (*cov)[j][k] += z[i][j] * z[i][k] / w[i];
2089 }
2090 }
2091 }
2092 }
2093
2094 if(m_logStream)
2095 *m_logStream<<"\\----- End of Eigenvalue Spectrum -----"<<std::endl;
2096
2097 ATH_MSG_DEBUG("Alignment constants:");
2098
2099 DataVector<AlignPar> * alignParList = m_alignModuleTool->alignParList1D();
2100
2101 // compute alignment corrections (translations in mm and rotations in rad) and their variances
2102 for(int i=0; i<size; i++) {
2103
2104 double param = delta[i];
2105 double err = sqrt(2.*std::fabs(errSq[i]));
2106
2107 int idof = m_activeIndices[i];
2108 AlignPar * alignPar=(*alignParList)[idof];
2109
2110 // undo the sigma scaling
2111 double sigma = alignPar->sigma();
2112
2113 param *= sigma;
2114 err *= sigma;
2115
2116 // undo normalization scaling of error
2117 if(m_scaleMatrix && m_scale>0.)
2118 err /= sqrt(m_scale);
2119
2120 ATH_MSG_DEBUG(i <<" : "<< param << " +/- "<< err);
2121 ATH_MSG_DEBUG("cov("<<i<<")="<<errSq[i]<<", sigma: "<<sigma);
2122 ATH_MSG_DEBUG("init par: "<<alignPar->initPar());
2123 alignPar->setPar(param, err);
2124 ATH_MSG_DEBUG("set param to "<<param<<" for alignPar "<<alignPar);
2125 ATH_MSG_DEBUG(*(*alignParList)[idof]);
2126 }
2127
2128 if(m_logStream) {
2130
2131 // norm of first derivative
2132 double norm1st = dChi2->norm();
2133 if(m_scaleMatrix && m_scale>0.) // undo normalization scaling
2134 norm1st *= m_scale;
2135 *m_logStream<<"norm of first derivative : "<<norm1st<<std::endl;
2136
2137 if(d2Chi2) {
2138 // distance to solution
2139 double dist = ( (*d2Chi2) * deltafull + (*dChi2) ).norm();
2140 if(m_scaleMatrix && m_scale>0.) // undo normalization scaling
2141 dist *= m_scale;
2142 *m_logStream<<"distance to solution : "<<dist<<std::endl;
2143
2144 // calculate chi2 of the alignment change
2145 double chi2 = delta * (*d2Chi2) * delta * .5;
2146 if(m_scaleMatrix && m_scale>0.) // undo normalization scaling
2147 chi2 *= m_scale;
2148 *m_logStream<<"delta(chi2) of the alignment change : "<<chi2<<" / "<<size<<std::endl;
2149 }
2150 }
2151
2152 delete cov;
2153 }
std::ostream * m_logStream
logfile output stream
Definition IMatrixTool.h:98
Gaudi::Property< bool > m_writeEigenMatTxt
Definition MatrixTool.h:184
std::vector< int > m_activeIndices
vector of indices which pass the min-hits cut
Definition MatrixTool.h:248
Gaudi::Property< double > m_calculateFullCovariance
Definition MatrixTool.h:215
Gaudi::Property< std::string > m_scalaVecName
Definition MatrixTool.h:231
Gaudi::Property< float > m_eigenvalueStep
Definition MatrixTool.h:169
Gaudi::Property< float > m_Align_db_step
Definition MatrixTool.h:171
Gaudi::Property< bool > m_scaleMatrix
Definition MatrixTool.h:205
Gaudi::Property< bool > m_writeEigenMat
Definition MatrixTool.h:182
Gaudi::Property< int > m_modcut
Definition MatrixTool.h:161
Gaudi::Property< double > m_eigenvaluethreshold
Definition MatrixTool.h:157
Gaudi::Property< std::string > m_scalaMatName
Definition MatrixTool.h:229
Gaudi::Property< float > m_pullcut
Definition MatrixTool.h:167
void printGlobalSolution(std::ostream &os, const CLHEP::HepSymMatrix *cov)
double chi2(TH1 *h0, TH1 *h1)
@ z
global position (cartesian)
Definition ParamDefs.h:57

◆ prepareBinaryFiles()

void MatrixTool::prepareBinaryFiles ( int solveOption)
virtual

reads/writes matrix entries from/to binary files as necessary

Implements Trk::IMatrixTool.

Definition at line 132 of file MatrixTool.cxx.

133 {
134 }

◆ printGlobalSolution() [1/2]

void MatrixTool::printGlobalSolution ( std::ostream & os,
const CLHEP::HepSymMatrix * cov )

Definition at line 1493 of file MatrixTool.cxx.

1494 {
1495 const AlignModuleList * alignModules = m_alignModuleTool->alignModules1D();
1496
1497 AlignModuleList::const_iterator imod = alignModules->begin();
1498 AlignModuleList::const_iterator imod_end = alignModules->end();
1499 for( ; imod!=imod_end; ++imod) {
1500 AlignModule * module = *imod;
1501
1502 DataVector<AlignPar> * alignPars = m_alignModuleTool->getAlignPars(module);
1503 int thisNDoF = alignPars->size();
1504
1505 // fill local covariance matrix
1506 CLHEP::HepSymMatrix * covsub = nullptr;
1507 if(cov && module->nHits() >= m_minNumHits && module->nTracks() >= m_minNumTrks) {
1508 covsub = new CLHEP::HepSymMatrix(thisNDoF,0);
1509 for (int i=0;i<thisNDoF;++i) {
1510 int ipar = alignPars->at(i)->index();
1511 double sigma_i = alignPars->at(i)->sigma();
1512
1513 std::vector<int>::iterator itActive = std::find(m_activeIndices.begin(),m_activeIndices.end(),ipar);
1514 if( itActive == m_activeIndices.end() )
1515 continue;
1516 int iActive = std::distance(m_activeIndices.begin(),itActive);
1517
1518 for (int j=0;j<=i;++j) {
1519 int jpar = alignPars->at(j)->index();
1520 double sigma_j = alignPars->at(j)->sigma();
1521
1522 std::vector<int>::iterator jtActive = std::find(m_activeIndices.begin(),m_activeIndices.end(),jpar);
1523 if( jtActive == m_activeIndices.end() )
1524 continue;
1525 int jActive = std::distance(m_activeIndices.begin(),jtActive);
1526
1527 (*covsub)[i][j] = (*cov)[iActive][jActive] * sigma_i * sigma_j;
1528 }
1529 }
1530 }
1531
1532 printModuleSolution(os,module,covsub);
1533
1534 delete covsub;
1535 }
1536 os << "--------------------------------------------------------------------------------" << std::endl;
1537 }
const T * at(size_type n) const
Access an element, as an rvalue.
Gaudi::Property< int > m_minNumHits
Definition MatrixTool.h:163
void printModuleSolution(std::ostream &os, const AlignModule *module, const CLHEP::HepSymMatrix *cov) const
namespace { class RestoreIOSFlags { public: RestoreIOSFlags (std::ostream &os) : m_os(&os),...
Gaudi::Property< int > m_minNumTrks
Definition MatrixTool.h:165

◆ printGlobalSolution() [2/2]

void MatrixTool::printGlobalSolution ( std::ostream & os,
const TMatrixDSym * cov )

Definition at line 1540 of file MatrixTool.cxx.

1541 {
1542 CLHEP::HepSymMatrix * cov = nullptr;
1543 if(cov0) {
1544 int nsize = cov0->GetNrows();
1545 cov = new CLHEP::HepSymMatrix(nsize,0);
1546
1547 for(int i=0; i<nsize; i++)
1548 for(int j=0; j<=i; j++)
1549 (*cov)[i][j] = (*cov0)[i][j];
1550 }
1551
1552 printGlobalSolution(os,cov);
1553
1554 delete cov;
1555 }

◆ printModuleSolution()

void MatrixTool::printModuleSolution ( std::ostream & os,
const AlignModule * module,
const CLHEP::HepSymMatrix * cov ) const

namespace { class RestoreIOSFlags { public: RestoreIOSFlags (std::ostream &os) : m_os(&os), m_flags(m_os->flags()), m_precision(m_os->precision()) {} ~RestoreIOSFlags() { m_os->flags(m_flags); m_os->precision(m_precision); } private: std::ostream *m_os; std::ios_base::fmtflags m_flags; std::streamsize m_precision; }; }

Definition at line 1578 of file MatrixTool.cxx.

1579 {
1580 boost::io::ios_all_saver ias( os ); //save the stream state
1581
1582 os << "--------------------------------------------------------------------------------" << std::endl;
1583 os << "Alignment parameters for module: " << module->name() << std::endl;
1584 os << "Number of tracks passing: " << module->nTracks() << std::endl;
1585 if(m_minNumHits>0 && module->nHits()<m_minNumHits) {
1586 os << "Number of hits too small: "<<module->nHits()<<" < "<<m_minNumHits<<" Skipping the module"<<std::endl;
1587 return;
1588 }
1589 if(m_minNumTrks>0 && module->nTracks()<m_minNumTrks) {
1590 os << "Number of tracks too small: "<<module->nTracks()<<" < "<<m_minNumTrks<<" Skipping the module"<<std::endl;
1591 return;
1592 }
1593 os << "Number of hits seen: " << module->nHits() << std::endl;
1594 os << "Number of tracks seen: " << module->nTracks() << std::endl;
1595
1596 DataVector<AlignPar> * alignPars = m_alignModuleTool->getAlignPars(module);
1597 int thisNDoF = alignPars->size();
1598
1599 if(alignPars->empty())
1600 os << "No active parameters" << std::endl;
1601 else
1602 {
1603 //RestoreIOSFlags restore_flags(os);
1604
1605 os.unsetf(std::ios_base::floatfield);
1606 os << std::setiosflags(std::ios_base::left) << std::setprecision(5);
1607
1608 // output alignment parameters and errors
1609 DataVector<AlignPar>::const_iterator ipar = alignPars->begin();
1610 DataVector<AlignPar>::const_iterator ipar_end = alignPars->end();
1611 for ( ; ipar != ipar_end; ++ipar) {
1612 const AlignPar * par = *ipar;
1613 os << std::setw(10) << par->dumpType()
1614 << std::setw(12) << par->par() << " +/- " << std::setw(12) << par->err()
1615 << std::endl;
1616 }
1617
1618 if(cov) {
1619 // calculate local correlation matrix
1620 CLHEP::HepSymMatrix corrsub(thisNDoF,0);
1621 for(int irow=0; irow<thisNDoF; ++irow)
1622 for(int icol=0; icol<=irow; ++icol)
1623 corrsub[irow][icol] = (*cov)[irow][icol] / sqrt((*cov)[irow][irow] * (*cov)[icol][icol]);
1624 os << "Local correlation matrix: " << corrsub << std::flush;
1625 }
1626 }
1627 ias.restore(); //restore the stream state
1628 }
DataModel_detail::const_iterator< DataVector > const_iterator
Standard const_iterator.
Definition DataVector.h:838
const_iterator end() const noexcept
Return a const_iterator pointing past the end of the collection.
const_iterator begin() const noexcept
Return a const_iterator pointing at the beginning of the collection.
bool empty() const noexcept
Returns true if the collection is empty.

◆ readHitmaps()

void MatrixTool::readHitmaps ( )
private

Definition at line 2327 of file MatrixTool.cxx.

2328 {
2329 ATH_MSG_INFO("read hitmaps from files");
2330
2331 const AlignModuleList * moduleList = m_alignModuleTool->alignModules1D();
2332 int nModules = moduleList->size();
2333
2334 AlMat hitmap(nModules,2);
2335 int nFiles = (int)m_inputHitmapFiles.size();
2336 for(int imap=0;imap<nFiles; imap++) {
2337 AlMat nextmap(nModules,2);
2338
2339 ATH_MSG_INFO("Reading hitmap "<<imap<<" from file "<<m_inputHitmapFiles[imap]);
2340
2341 if(nextmap.ReadScalaPack(m_inputHitmapFiles[imap]).isFailure()) {
2342 ATH_MSG_WARNING("Problem reading hitmap from \'"<<m_inputHitmapFiles[imap]<<"\'. Skipping.");
2343 continue;
2344 }
2345
2346 if(nextmap.nrow()!=nModules || nextmap.ncol()!=2) {
2347 ATH_MSG_WARNING("Matrix in file \'"<<m_inputHitmapFiles[imap]<<"\' has wrong size ("
2348 <<nextmap.nrow()<<" x "<<nextmap.ncol()<<"), should be ("<<nModules<<" x 2). Skipping.");
2349 continue;
2350 }
2351
2352 hitmap += nextmap;
2353 }
2354
2355 AlignModuleList::const_iterator imod = moduleList->begin();
2356 AlignModuleList::const_iterator imod_end = moduleList->end();
2357 int index = 0;
2358 int totalhits = 0;
2359 for(; imod != imod_end; ++imod) {
2360 AlignModule * module = *imod;
2361 //if ((int)hitmap[index][0]!=(int)module->identify().get_identifier32().get_compact()) ATH_MSG_ERROR("bad module identifier");
2362 //module->setIdentifier((Identifier)hitmap[index][0]);
2363 module->setNHits((int)hitmap[index][0]);
2364 module->setNTracks((int)hitmap[index][1]);
2365 totalhits += (int)hitmap[index][0];
2366 index++;
2367 }
2368
2369 m_nHits = totalhits;
2370 m_nTracks = 0;
2371 m_nMeasurements = 0;
2372
2373 ATH_MSG_INFO("Hitmap accumulated from "<<nFiles<<" files with total of "<<totalhits<<" hits.");
2374 }
Gaudi::Property< std::vector< std::string > > m_inputHitmapFiles
Definition MatrixTool.h:241

◆ renounce()

std::enable_if_t< std::is_void_v< std::result_of_t< decltype(&T::renounce)(T)> > &&!std::is_base_of_v< SG::VarHandleKeyArray, T > &&std::is_base_of_v< Gaudi::DataHandle, T >, void > AthCommonDataStore< AthCommonMsg< AlgTool > >::renounce ( T & h)
inlineprotectedinherited

Definition at line 380 of file AthCommonDataStore.h.

381 {
382 h.renounce();
384 }
std::enable_if_t< std::is_void_v< std::result_of_t< decltype(&T::renounce)(T)> > &&!std::is_base_of_v< SG::VarHandleKeyArray, T > &&std::is_base_of_v< Gaudi::DataHandle, T >, void > renounce(T &h)

◆ renounceArray()

void AthCommonDataStore< AthCommonMsg< AlgTool > >::renounceArray ( SG::VarHandleKeyArray & handlesArray)
inlineprotectedinherited

remove all handles from I/O resolution

Definition at line 364 of file AthCommonDataStore.h.

364 {
366 }

◆ setLogStream()

virtual void Trk::IMatrixTool::setLogStream ( std::ostream * os)
inlinevirtualinherited

sets the output stream for the logfile

Definition at line 95 of file IMatrixTool.h.

95{ m_logStream = os; }

◆ setNHits()

void Trk::IMatrixTool::setNHits ( int n)
inlineinherited

set module identifier

set number of hits

Definition at line 83 of file IMatrixTool.h.

◆ setNMeasurements()

void Trk::IMatrixTool::setNMeasurements ( int n)
inlineinherited

set number of measurements

Definition at line 91 of file IMatrixTool.h.

91{ m_nMeasurements = n; }

◆ setNTracks()

void Trk::IMatrixTool::setNTracks ( int n)
inlineinherited

set number of tracks

Definition at line 87 of file IMatrixTool.h.

87{ m_nTracks = n; }

◆ solve()

int MatrixTool::solve ( )
virtual

solves for alignment parameters

Implements Trk::IMatrixTool.

Definition at line 1140 of file MatrixTool.cxx.

1141 {
1142 // ============
1143 // solve
1144 // ============
1145 ATH_MSG_DEBUG("in MatrixTool::solve()");
1146
1147 // set normalization scale to number of hits for now
1148 if(m_scale<0)
1149 m_scale = m_nHits;
1150
1151 //-------------------------------------------------------
1152 // write matrix and vector to file
1153 if (m_writeMat) {
1154 // version has to be 2 to for reading matrices and vectors back in to work properly
1155 double dummyVersion(2.);
1156
1157 // make map of matrix entry to module index (set by geometry manager tool)
1158 std::map<int,unsigned long long> modIndexMap;
1159 std::map<int,std::string> modNameMap;
1160 DataVector<AlignPar>* alignPars = m_alignModuleTool->alignParList1D();
1161 for (int i=0;i<(int)alignPars->size();i++) {
1162 modIndexMap[i]=(*alignPars)[i]->alignModule()->identify().get_compact();
1163 modNameMap [i]=(*alignPars)[i]->alignModule()->name();
1164 }
1165
1166 // binary files
1167 ATH_MSG_DEBUG("writing binary files");
1168 StatusCode sc1 = m_bigmatrix->Write("matrix.bin",true,m_wSqMatrix,m_scale,dummyVersion);
1169 StatusCode sc2 = m_bigvector->WritePartial("vector.bin",true,m_scale,modIndexMap,dummyVersion);
1170 if (!sc1.isSuccess() || !sc2.isSuccess()) {
1171 msg(MSG::ERROR)<<"problem writing matrix or vector"<<endmsg;
1172 return -1;
1173 }
1174
1175 if (m_writeMatTxt) {
1176
1177 // text files
1178 ATH_MSG_DEBUG("writing text files");
1179 sc1 = m_bigmatrix->Write("matrix.txt",false,m_wSqMatrix,m_scale,dummyVersion);
1180 sc2 = m_writeModuleNames ?
1181 m_bigvector->WritePartial("vector.txt",false,m_scale,modNameMap,dummyVersion) :
1182 m_bigvector->WritePartial("vector.txt",false,m_scale,modIndexMap,dummyVersion);
1183
1184 if (!sc1.isSuccess() || !sc2.isSuccess()) {
1185 msg(MSG::ERROR)<<"problem writing matrix or vector"<<endmsg;
1186 return -1;
1187 }
1188 }
1189
1190 ATH_MSG_DEBUG("matrix and vector written to: "<<m_pathbin.value()+m_prefixName.value()<<"matrix.bin (.txt) and "<<m_pathbin.value()+m_prefixName.value()<<"vector.bin (.txt)");
1191 }
1192
1193 //-------------------------------------------------------
1194 // write hitmap to file
1195 if (m_writeHitmap)
1196 writeHitmap();
1197
1198 if(m_writeTFile)
1199 storeInTFile(m_pathbin.value()+m_prefixName.value()+m_tfileName.value());
1200
1201
1202 if(!m_runLocal && m_solveOption==0) {
1203 ATH_MSG_DEBUG("No solving requested.");
1204 return 1;
1205 }
1206
1207 //-------------------------------------------------------
1208 // rescale the vector and the matrix according to sigmas
1209 // and apply soft mode cut
1210
1211 ATH_MSG_DEBUG("rescaling the matrix/vector and applying the soft-mode-cut");
1212
1213 DataVector<AlignPar>* alignParList = m_alignModuleTool->alignParList1D();
1214 int nDoF = alignParList->size();
1215
1216
1217 const AlSymMat * chkMatrix = dynamic_cast<const AlSymMat*>(m_bigmatrix);
1218 if(chkMatrix){
1219 // Method when using the dense matrix
1220 for (int i=0;i<nDoF;i++) {
1221 // scale the vector
1222 double sigma_i = (*alignParList)[i]->sigma();
1223 double softCut = 2 * pow( (*alignParList)[i]->softCut() , -2 );
1224 (*m_bigvector)[i] *= sigma_i;
1225
1226 for (int j=0;j<=i;j++) {
1227 // scale the matrix
1228 if ((*chkMatrix)[i][j] != 0.) {
1229 double sigma_j = (*alignParList)[j]->sigma();
1230 (*m_bigmatrix)[i][j] *= sigma_i * sigma_j;
1231 }
1232 // apply soft-mode-cut
1233 if (i==j && m_softEigenmodeCut>0.){
1234 (*m_bigmatrix)[i][j] += m_softEigenmodeCut * softCut;
1235
1236 }
1237
1238 // set first and second derivatives on AlignPar
1239 if (i==j) {
1240 (*alignParList)[i]->setFirstDeriv((*m_bigvector)[i]/sigma_i);
1241 (*alignParList)[i]->setSecndDeriv((*m_bigmatrix)[i][i]/sigma_i/sigma_i);
1242 }
1243 }
1244 }
1245 } else {
1246 // Method when using the sparse matrix
1247 for (const datamap::value_type& p : *m_bigmatrix->ptrMap()) {
1248 int i = p.first.first;
1249 int j = p.first.second;
1250
1251 // Scale matrix
1252 double sigma_i = (*alignParList)[i]->sigma();
1253 double sigma_j = (*alignParList)[j]->sigma();
1254
1255 (*m_bigmatrix)[i][j] *= sigma_i * sigma_j;
1256
1257 }
1258
1259
1260 for (int i=0;i<nDoF;i++) {
1261 // scale the vector
1262 double sigma_i = (*alignParList)[i]->sigma();
1263 (*m_bigvector)[i] *= sigma_i;
1264 if (m_softEigenmodeCut >0. ){
1265 double softCut = 2 * pow( (*alignParList)[i]->softCut() , -2 );
1266 (*m_bigmatrix)[i][i] += m_softEigenmodeCut * softCut;
1267 ATH_MSG_DEBUG( "DOF "<< i <<" Nhits "<< (*alignParList)[i]->alignModule()->nHits() << " soft-mode-cut "<< (*alignParList)[i]->softCut() <<" -> " << m_softEigenmodeCut * softCut << " sigma_i "<< sigma_i << " Matrix: " << (*m_bigmatrix)[i][i] << " Vector: " << (*m_bigvector)[i]);
1268 }
1269 (*alignParList)[i]->setFirstDeriv((*m_bigvector)[i]/sigma_i);
1270 (*alignParList)[i]->setSecndDeriv((*m_bigmatrix)[i][i]/sigma_i/sigma_i);
1271 }
1272 }
1273
1274 unsigned long long OldPixelIdentifier = 37769216; //Identifier for the Pixel Detector
1275 unsigned long long IBLIdentifier = 33574912; //Identifier for the Pixel Detector
1276
1277 unsigned long long SCT_ECA_8_Identifier = 218116096; //Identifier for the SCT ECA Last Disk
1278 std::string SCT_ECA_8_Name = "SCT/EndcapA/Disk_8";
1279
1280
1281
1282 ATH_MSG_INFO("rescaling done");
1283 ATH_MSG_INFO("Javi: Printing (*alignParList)[i]->alignModule()->identify32()");
1284
1285 // select modules with non-zero tracks
1286 for(int i=0;i<nDoF;i++)
1287 {
1288 ATH_MSG_DEBUG(i);
1289 ATH_MSG_DEBUG((*alignParList)[i]->alignModule()->identify32());
1290 ATH_MSG_DEBUG((*alignParList)[i]->alignModule());
1291 ATH_MSG_DEBUG((*alignParList)[i]->alignModule()->name());
1292 ATH_MSG_DEBUG((*alignParList)[i]->paramType());
1293
1294 //Skip solving for Pixel or IBL:
1295 const auto & theParameterList = *alignParList;
1296 const auto & thisIdentifier = theParameterList[i]->alignModule()->identify32();
1297 const auto & thisName = theParameterList[i]->alignModule()->name();
1298 const auto & thisParameterType = theParameterList[i]->paramType();
1299 const bool oldPixel = (thisIdentifier == OldPixelIdentifier);
1300 const bool ibl = (thisIdentifier == IBLIdentifier);
1301 const bool SCTECA8 = (thisIdentifier == SCT_ECA_8_Identifier);
1302 const bool SCTECA8_n = (thisName.find(SCT_ECA_8_Name)!= std::string::npos);
1303
1305 if (SCTECA8)
1306 {ATH_MSG_INFO( "SCT ECA Last Disk DoF have been skipped in the solving because DeactivateSCT_ECA_LastDisk is set to True");
1307 continue;}
1308 if (SCTECA8_n)
1309 {ATH_MSG_INFO( "SCT ECA Last Disk DoF have been skipped in the solving because DeactivateSCT_ECA_LastDisk is set to True");
1310 continue;}
1311 }
1312
1313 if (m_AlignIBLbutNotPixel) //If m_AlignIBLbutNotPixel is set to True, Pixel will be skipped in the solving.
1314 if (oldPixel)
1315 {ATH_MSG_INFO( "Pixel DoF have been skipped in the solving because AlignIBLbutNotPixel is set to True");
1316 continue;}
1317
1318 if (m_AlignPixelbutNotIBL) //If m_AlignPixelbutNotIBL is set to True, IBL will be skipped in the solving.
1319 if (ibl)
1320 {ATH_MSG_INFO( "IBL DoF have been skipped in the solving because AlignPixelbutNotIBL is set to True");
1321 continue;}
1322
1323 //For specific DoF: (*alignParList)[i]->paramType() = 0,1,2,3,4,5,6 for Tx,Ty,Tz,Rx,Ry,Rz,Bx
1324 //Pixel Dofs:
1325 if (m_Remove_Pixel_Tx) //If m_Remove_Pixel_Tx is set to True, Pixel Tx will be skipped in the solving.
1326 if ( oldPixel and (thisParameterType == 0))
1327 {ATH_MSG_INFO( "Pixel Tx DoF has been skipped in the solving because Remove_Pixel_Tx is set to True");
1328 continue;}
1329
1330 if (m_Remove_Pixel_Ty) //If m_Remove_Pixel_Ty is set to True, Pixel Ty will be skipped in the solving.
1331 if ( oldPixel and (thisParameterType == 1))
1332 {ATH_MSG_INFO( "Pixel Ty DoF has been skipped in the solving because Remove_Pixel_Ty is set to True");
1333 continue;}
1334
1335 if (m_Remove_Pixel_Tz) //If m_Remove_Pixel_Tz is set to True, Pixel Tz will be skipped in the solving.
1336 if (oldPixel and (thisParameterType == 2))
1337 {ATH_MSG_INFO( "Pixel Tz DoF has been skipped in the solving because Remove_Pixel_Tz is set to True");
1338 continue;}
1339
1340 if (m_Remove_Pixel_Rx) //If m_Remove_Pixel_Rx is set to True, Pixel Rx will be skipped in the solving.
1341 if (oldPixel and (thisParameterType == 3))
1342 {ATH_MSG_INFO( "Pixel Rx DoF has been skipped in the solving because Remove_Pixel_Rx is set to True");
1343 continue;}
1344
1345 if (m_Remove_Pixel_Ry) //If m_Remove_Pixel_Ry is set to True, Pixel Ry will be skipped in the solving.
1346 if (oldPixel and (thisParameterType == 4))
1347 {ATH_MSG_INFO( "Pixel Ry DoF has been skipped in the solving because Remove_Pixel_Ry is set to True");
1348 continue;}
1349
1350 if (m_Remove_Pixel_Rz) //If m_Remove_Pixel_Rz is set to True, Pixel Rz will be skipped in the solving.
1351 if (oldPixel and (thisParameterType == 5))
1352 {ATH_MSG_INFO( "Pixel Rz DoF has been skipped in the solving because Remove_Pixel_Rz is set to True");
1353 continue;}
1354
1355 //IBL Dofs:
1356 if (m_Remove_IBL_Tx) //If m_Remove_IBL_Tx is set to True, IBL Tx will be skipped in the solving.
1357 if (ibl and (thisParameterType == 0))
1358 {ATH_MSG_INFO( "IBL Tx DoF has been skipped in the solving because Remove_IBL_Tx is set to True");
1359 continue;}
1360
1361 if (m_Remove_IBL_Ty) //If m_Remove_IBL_Ty is set to True, IBL Ty will be skipped in the solving.
1362 if (ibl and (thisParameterType == 1))
1363 {ATH_MSG_INFO( "IBL Ty DoF has been skipped in the solving because Remove_IBL_Ty is set to True");
1364 continue;}
1365
1366 if (m_Remove_IBL_Tz) //If m_Remove_IBL_Tz is set to True, IBL Tz will be skipped in the solving.
1367 if (ibl and (thisParameterType == 2))
1368 {ATH_MSG_INFO( "IBL Tz DoF has been skipped in the solving because Remove_IBL_Tz is set to True");
1369 continue;}
1370
1371 if (m_Remove_IBL_Rx) //If m_Remove_IBL_Rx is set to True, IBL Rx will be skipped in the solving.
1372 if (ibl and (thisParameterType == 3))
1373 {ATH_MSG_INFO( "IBL Rx DoF has been skipped in the solving because Remove_IBL_Rx is set to True");
1374 continue;}
1375
1376 if (m_Remove_IBL_Ry) //If m_Remove_IBL_Ry is set to True, IBL Ry will be skipped in the solving.
1377 if (ibl and (thisParameterType == 4))
1378 {ATH_MSG_INFO( "IBL Ry DoF has been skipped in the solving because Remove_IBL_Ry is set to True");
1379 continue;}
1380
1381 if (m_Remove_IBL_Rz) //If m_Remove_IBL_Rz is set to True, IBL Rz will be skipped in the solving.
1382 if (ibl and (thisParameterType == 5))
1383 {ATH_MSG_INFO( "IBL Rz DoF has been skipped in the solving because Remove_IBL_Rz is set to True");
1384 continue;}
1385
1386 if(theParameterList[i]->alignModule()->nHits() >= m_minNumHits && theParameterList[i]->alignModule()->nTracks() >= m_minNumTrks)
1387 m_activeIndices.push_back(i);
1388 }
1389 m_aNDoF = m_activeIndices.size();
1390 ATH_MSG_DEBUG("aNDoF/nDoF: "<<m_aNDoF<<"/"<<nDoF);
1391
1392 // --------------------
1393 // now do the SOLVING
1394 // --------------------
1395
1396 int info = 0;
1397
1398 // first Local solving
1399 if (m_runLocal)
1400 info = solveLocal();
1401
1402 // remove spurious modules and resize
1403 if (m_removeSpurious) {
1404
1405 ATH_MSG_INFO("Spurious removal not implemented at the moment.");
1406/* if (StatusCode::SUCCESS != spuriousRemoval()) {
1407 ATH_MSG_ERROR("Problem while trying to remove spurious. Stopping solving");
1408 return -1;
1409 }
1410
1411 // if nDoF=0, bad job...
1412 int NDoF = m_alignModuleTool->nAlignParameters();
1413 if(NDoF==0) {
1414 ATH_MSG_WARNING("Removal removed everything: NDoF=" << NDoF << " !!!");
1415 return 1;
1416 }
1417 ATH_MSG_DEBUG("NDoF: " << NDoF);
1418*/
1419 }
1420
1421 // --------------------
1422 // now Global solving
1423 switch(m_solveOption) {
1424
1425 case NONE:
1426 ATH_MSG_DEBUG("No global solving requested.");
1427 break;
1428
1429 case SOLVE_ROOT:
1430 info = solveROOT();
1431 break;
1432
1433 case SOLVE_CLHEP:
1434 info = solveCLHEP();
1435 break;
1436
1437 case SOLVE:
1438 case DIRECT_SOLVE:
1440 info = solveLapack();
1441 break;
1442
1443 case SOLVE_FAST:
1444 case DIRECT_SOLVE_FAST:
1446 break;
1447
1448 default:
1449 ATH_MSG_INFO("Unknown solving option.");
1450 info = 0;
1451 break;
1452 }
1453
1454 ATH_MSG_INFO("Return value from solving: "<<info);
1455
1456 return info;
1457 }
constexpr int pow(int base, int exp) noexcept
int nTracks() const
Definition IMatrixTool.h:88
int nHits() const
Definition IMatrixTool.h:84
Gaudi::Property< std::string > m_tfileName
Definition MatrixTool.h:226
Gaudi::Property< bool > m_AlignIBLbutNotPixel
Definition MatrixTool.h:256
Gaudi::Property< bool > m_Remove_Pixel_Ry
Definition MatrixTool.h:268
Gaudi::Property< bool > m_writeHitmap
Definition MatrixTool.h:189
Gaudi::Property< bool > m_Remove_Pixel_Rz
Definition MatrixTool.h:269
Gaudi::Property< bool > m_Remove_IBL_Rx
Definition MatrixTool.h:275
Gaudi::Property< bool > m_runLocal
Definition MatrixTool.h:201
Gaudi::Property< double > m_removeSpurious
Definition MatrixTool.h:211
Gaudi::Property< bool > m_Remove_Pixel_Tz
Definition MatrixTool.h:266
Gaudi::Property< bool > m_Remove_IBL_Rz
Definition MatrixTool.h:277
Gaudi::Property< bool > m_Remove_Pixel_Ty
Definition MatrixTool.h:265
Gaudi::Property< bool > m_Remove_IBL_Tx
Definition MatrixTool.h:272
int m_aNDoF
number of active DoF (size of m_activeIndices)
Definition MatrixTool.h:249
Gaudi::Property< bool > m_Remove_Pixel_Tx
Definition MatrixTool.h:264
Gaudi::Property< bool > m_Remove_IBL_Tz
Definition MatrixTool.h:274
Gaudi::Property< bool > m_writeMatTxt
Definition MatrixTool.h:180
Gaudi::Property< int > m_solveOption
Definition MatrixTool.h:160
Gaudi::Property< bool > m_DeactivateSCT_ECA_LastDisk
Definition MatrixTool.h:261
Gaudi::Property< bool > m_writeMat
Definition MatrixTool.h:178
void storeInTFile(const TString &filename)
Store Files in a tfile.
Gaudi::Property< double > m_softEigenmodeCut
Definition MatrixTool.h:208
Gaudi::Property< bool > m_writeTFile
Definition MatrixTool.h:196
Gaudi::Property< bool > m_Remove_Pixel_Rx
Definition MatrixTool.h:267
Gaudi::Property< bool > m_AlignPixelbutNotIBL
Definition MatrixTool.h:258
Gaudi::Property< bool > m_Remove_IBL_Ry
Definition MatrixTool.h:276
Gaudi::Property< bool > m_Remove_IBL_Ty
Definition MatrixTool.h:273
Gaudi::Property< bool > m_writeModuleNames
Definition MatrixTool.h:186

◆ solveCLHEP()

int MatrixTool::solveCLHEP ( )
private

Definition at line 263 of file MatrixTool.cxx.

264 {
265 ATH_MSG_INFO("solving Global using CLHEP");
266 if(m_logStream) {
267 *m_logStream<<"*************************************************************"<<std::endl;
268 *m_logStream<<"************** solving using Global method ****************"<<std::endl;
269 *m_logStream<<"************** using CLHEP ****************"<<std::endl;
270 *m_logStream<<"*************************************************************"<<std::endl;
271 }
272
273 // start measuring time
274 clock_t starttime = clock();
275
276 DataVector<AlignPar> * alignParList = m_alignModuleTool->alignParList1D();
277 //const AlignModuleList* alignModules = m_alignModuleTool->alignModules1D();
278 for (int i=0;i<(int)alignParList->size();i++)
279 ATH_MSG_DEBUG("ap["<<i<<"]="<<(*alignParList)[i]);
280
281 int nDoF = m_alignModuleTool->nAlignParameters();
282
283 // some debugging output
284 if (msgLvl(MSG::DEBUG)) {
285 msg(MSG::DEBUG)<<"dumping matrix and vector to screen"<<endmsg;
286 for (int i=0;i<nDoF;i++)
287 for (int j=0;j<nDoF;j++)
288 //if (std::fabs((*m_bigmatrix)[i][j])>.0001)
289 msg(MSG::DEBUG)<<i<<", "<<j<<" : "<<(*m_bigmatrix)[i][j] <<endmsg;
290
291 for (int i=0;i<nDoF;i++)
292 msg(MSG::DEBUG)<<i <<" : "<<(*m_bigvector)[i]<<endmsg;
293 }
294
295 // get rescaled first and second derivatives
296 CLHEP::HepSymMatrix * d2Chi2 = new CLHEP::HepSymMatrix(m_aNDoF,0);
297 CLHEP::HepVector * dChi2 = new CLHEP::HepVector(m_aNDoF,0);
298 for (int iActive=0;iActive<m_aNDoF;iActive++) {
299 int i = m_activeIndices[iActive];
300 (*dChi2)[iActive] = (*m_bigvector)[i];
301 for (int jActive=0;jActive<m_aNDoF;jActive++) {
302 int j = m_activeIndices[jActive];
303 (*d2Chi2)[iActive][jActive] = (*m_bigmatrix)[i][j];
304 }
305 }
306
307 ATH_MSG_DEBUG("First derivatives:" << (*dChi2));
308 ATH_MSG_DEBUG("Second derivatives:" << (*d2Chi2));
309
310 CLHEP::HepSymMatrix cov(m_aNDoF,0);
311 CLHEP::HepVector delta(m_aNDoF,0);
312 CLHEP::HepVector deltafull(m_aNDoF,0);
313
314 bool status=true;
315 int ierr(0);
316 if(!m_diagonalize) {
317 // ==========================================================
318 // Run Matrix Inversion
319 ATH_MSG_INFO("Running matrix inversion");
320 if(m_logStream)
321 *m_logStream<<"Running matrix inversion"<<std::endl;
322
323 cov = *d2Chi2;
324 cov.invert(ierr);
325 if(ierr>0)
326 msg(MSG::ERROR)<<"CLHEP inversion status flag = "<<ierr<<endmsg;
327 else
328 ATH_MSG_INFO("CLHEP inversion OK");
329 if(m_logStream)
330 *m_logStream<<"CLHEP inversion status flag = "<<ierr<<std::endl;
331
332 // calculate corrections
333 delta = cov * (*dChi2);
334
335 // the covariance matrix is actually defined as 2 * d2Chi2^-1
336 ATH_MSG_DEBUG("Result: "<<delta);
337 ATH_MSG_DEBUG("cov: "<<cov*2);
338
339 // -----------------------
340 // calculate also for matrix and vector multiplied by factor 0.5
341 // this should make no difference if everything is correct but
342 // it can go wrong if insensitive DoF is included
343 CLHEP::HepSymMatrix cov2 = *d2Chi2 * .5;
344 // invert the matrix
345 int ierr2 = 0;
346 cov2.invert(ierr2);
347 if(ierr2>0)
348 msg(MSG::WARNING)<<"Second CLHEP inversion status flag = "<<ierr2<<endmsg;
349
350 CLHEP::HepVector delta2 = cov2 * (*dChi2) * .5;
351 for (int i=0;i<delta.num_row(); ++i)
352 if ( fabs((delta[i] - delta2[i])/delta[i]) > 1e-5 ) {
353 msg(MSG::WARNING)<<"Something's wrong with the matrix inversion: delta["<<i<<"] = "<<delta[i]<<" delta2["<<i<<"] = "<<delta2[i]<<endmsg;
354 status=false;
355 break;
356 }
357
358 if(m_logStream && (ierr2>0 || !status)) {
359 *m_logStream<<"CLHEP inversion status flag for halfed matrix = "<<ierr2<<std::endl;
360 *m_logStream<<"Matrix inversion check failed"<<std::endl;
361 *m_logStream<<std::endl;
362 }
363 // -- end of check of matrix inversion
364 }
365 else {
366 // ==========================================================
367 // Run Diagonalization
368 ATH_MSG_INFO("Running diagonalization");
369 if(m_logStream)
370 *m_logStream<<"Running diagonalization"<<std::endl;
371
372 CLHEP::HepSymMatrix D = *d2Chi2;
373 CLHEP::HepMatrix U = CLHEP::diagonalize( &D );
374
375 ATH_MSG_INFO("Diagonalization done");
376 //sold = U*sdiag*U.T.
377
378 // reorder eigenvalues ascending
379 // eigenvectors need to be reordered consistently
380 ATH_MSG_DEBUG(" Reordering eigenvalues ascending ");
381 for (int i=0; i<m_aNDoF-1; i++)
382 for (int j=i+1; j<m_aNDoF; j++)
383 if(D[j][j] < D[i][i]) {
384 // swap eigenvalues
385 double ei = D[i][i];
386 D[i][i] = D[j][j];
387 D[j][j] = ei;
388 // swap eigenvectors
389 for(int k=0;k<m_aNDoF; k++) {
390 double ev = U[k][i];
391 U[k][i] = U[k][j];
392 U[k][j] = ev;
393 }
394 }
395
396 // how do I now get the eigenvalues? this cannot be the most
397 // efficient way ... CLHEP::HepSymMatrix D = d2Chi2->similarityT( U );
398 CLHEP::HepVector eigenvector(m_aNDoF);
399
400 if(m_logStream)
401 *m_logStream<<"/------ The Eigenvalue Spectrum -------"<<std::endl;
402
403 ATH_MSG_DEBUG("Calculating eigenvalues");
404 for(int imode=0; imode<m_aNDoF; ++imode) {
405
406 // get the relevant eigenvector
407 for(int irow=0; irow<m_aNDoF; ++irow)
408 eigenvector[irow] = U[irow][imode];
409
410 // calculate the eigenvalue
411 //double eigenvalue = d2Chi2->similarity( eigenvector );
412 double eigenvalue = D[imode][imode];
413 ATH_MSG_DEBUG("eigenvalue "<<eigenvalue);
414
415 double evdotb = dot(*dChi2,eigenvector);
416 CLHEP::HepVector thisdelta = evdotb/eigenvalue * eigenvector;
417 deltafull += thisdelta;
418
419 if(imode<m_modcut) {
420 ATH_MSG_INFO("skipping eigenvalue "<<imode<<" : "<<eigenvalue<<" , modcut is "<<m_modcut);
421 if(m_logStream)
422 *m_logStream<<"| skipping eigenvalue "<<eigenvalue<<std::endl;
423 }
424 else if( eigenvalue < m_eigenvaluethreshold ) {
425 ATH_MSG_INFO("skipping eigenvalue "<<eigenvalue<<" , cut is "<<m_eigenvaluethreshold);
426 if(m_logStream)
427 *m_logStream<<"| skipping eigenvalue "<<eigenvalue<<std::endl;
428 }
429 else {
430 if(m_logStream)
431 *m_logStream<<"| "<<eigenvalue<<std::endl;
432
433 delta += thisdelta;
434
435 // this is the time consuming part
436 for(int irow=0; irow<m_aNDoF; ++irow)
437 for(int icol=0; icol<=irow; ++icol)
438 cov[irow][icol] += eigenvector[irow] * eigenvector[icol] / eigenvalue;
439 }
440 }
441
442 // the covariance matrix is actually defined as 2 * d2Chi2^-1
443
444 ATH_MSG_DEBUG("Result: "<<delta);
445 ATH_MSG_DEBUG("cov: "<<cov);
446
447 if(m_logStream)
448 *m_logStream<<"\\----- End of Eigenvalue Spectrum -----"<<std::endl;
449
450 // end of diagonalization
451 // ==========================================================
452 }
453
454 // stop measuring time
455 clock_t stoptime = clock();
456 double totaltime = (stoptime-starttime)/double(CLOCKS_PER_SEC);
457 ATH_MSG_INFO("Time spent in solveCLHEP: "<<totaltime<<" s");
458
459 if(ierr==0 && status)
460 {
461 ATH_MSG_DEBUG("Alignment constants:");
462 for (int iAdof=0;iAdof<m_aNDoF;iAdof++) {
463
464 int idof = m_activeIndices[iAdof];
465 AlignPar * alignPar=(*alignParList)[idof];
466
467 double sigma = alignPar->sigma();
468 double param = -delta[iAdof] * sigma;
469 double err = std::sqrt(2.*std::fabs(cov[iAdof][iAdof])) * sigma;
470
471 ATH_MSG_DEBUG(iAdof <<" : "<< param << " +/- "<< err);
472 ATH_MSG_DEBUG("cov("<<iAdof<<")="<<cov[iAdof][iAdof]<<", sigma: "<<sigma);
473 ATH_MSG_DEBUG("init par: "<<alignPar->initPar());
474 alignPar->setPar(param,err);
475 ATH_MSG_DEBUG("set param to "<<param<<" for alignPar "<<alignPar);
476 ATH_MSG_DEBUG(*(*alignParList)[idof]);
477 }
478
479 if(m_logStream) {
481
482 // norm of first derivative
483 *m_logStream<<"norm of first derivative : "<<dChi2->norm()<<std::endl;
484
485 // distance to solution
486 double dist = ( - (*d2Chi2) * deltafull + (*dChi2) ).norm();
487 *m_logStream<<"distance to solution : "<<dist<<std::endl;
488
489 // calculate chi2 of the alignment change
490 double chi2 = d2Chi2->similarity(delta) * .5;
491 *m_logStream<<"delta(chi2) of the alignment change : "<<chi2<<" / "<<m_aNDoF<<std::endl;
492
493 // time spent here
494 *m_logStream<<"time spent in solve : "<<totaltime<<" s"<<std::endl;
495 }
496 }
497
498 delete d2Chi2;
499 delete dChi2;
500
501 return 1;
502 }
bool msgLvl(const MSG::Level lvl) const
Gaudi::Property< bool > m_diagonalize
Definition MatrixTool.h:154
int ev
Definition globals.cxx:25
dot(G, fn, nodesToHighlight=[])
Definition dot.py:5
status
Definition merge.py:16

◆ solveLapack()

int MatrixTool::solveLapack ( )
private

Definition at line 1637 of file MatrixTool.cxx.

1638 {
1639 ATH_MSG_INFO("solving Global using Lapack");
1640 if(m_logStream) {
1641 *m_logStream<<"*************************************************************"<<std::endl;
1642 *m_logStream<<"************** solving using Global method ****************"<<std::endl;
1643 *m_logStream<<"************** using LAPACK ****************"<<std::endl;
1644 *m_logStream<<"*************************************************************"<<std::endl;
1645 }
1646
1647 // get rescaled first and second derivatives
1648 AlSymMat* aBetterMat = new AlSymMat(m_aNDoF);
1649 AlVec* aBetterVec = new AlVec(m_aNDoF);
1650 for (int iActive=0;iActive<m_aNDoF;iActive++) {
1651 int i = m_activeIndices[iActive];
1652 (*aBetterVec)[iActive] = (*m_bigvector)[i];
1653 for (int jActive=0;jActive<m_aNDoF;jActive++) {
1654 int j = m_activeIndices[jActive];
1655 (*aBetterMat)[iActive][jActive] = (*m_bigmatrix)[i][j];
1656 }
1657 }
1658
1659 // normalize bigmatrix and bigvector
1660 if(m_scaleMatrix) {
1661 if(m_scale<=0.)
1662 ATH_MSG_WARNING("Scaling requested but scale not set. Not scaling matrix and vector.");
1663 else {
1664 (*aBetterVec) *= 1./m_scale;
1665 (*aBetterMat) *= 1./m_scale;
1666 }
1667 }
1668
1669 ATH_MSG_DEBUG("Now Solving alignment using lapack diagonalization routine dspev...");
1670
1671 if (m_calDet) {
1672 const double tol = 1.e-20;
1673 // compute final determinant
1674 double determ = (*aBetterMat).determinant();
1675 ATH_MSG_INFO("Determinant: " << determ);
1676 if (fabs(determ) < tol)
1677 ATH_MSG_WARNING("Matrix is singular!");
1678 }
1679
1680 // store the original matrix for checks
1681 AlSymMat * d2Chi2 = nullptr;
1683 d2Chi2 = new AlSymMat(*aBetterMat);
1684
1685 clock_t starttime = clock();
1686
1687 // declare transition matrix + vector to store eigenvalues
1688 AlMat z(m_aNDoF,m_aNDoF);
1689 AlVec w(m_aNDoF); // vector to store the eigenvalues
1690 ATH_MSG_DEBUG("MatrixTool::after z/w allocation");
1691
1692 char jobz = 'V';
1693 int info = (*aBetterMat).diagonalize(jobz,w,z);
1694 ATH_MSG_DEBUG(" info: " << info);
1695 ATH_MSG_INFO("MatrixTool::after diagonalization");
1696
1697 // stop time calculation
1698 clock_t stoptime = clock();
1699 double time_diag = (stoptime-starttime)/double(CLOCKS_PER_SEC);
1700 ATH_MSG_INFO(" - time spent diagonalizing the matrix: "<<time_diag<<" s");
1701
1702 double time_solve = 0.;
1703 if (info==0) {
1704 starttime = clock();
1705 postSolvingLapack(aBetterVec,d2Chi2,w,z,m_aNDoF);
1706 stoptime = clock();
1707 time_solve = (stoptime-starttime)/double(CLOCKS_PER_SEC);
1708 ATH_MSG_INFO(" - time spent solving the system: "<<time_solve<<" s");
1709 if(m_logStream) {
1710 *m_logStream<<"time spent for diagonalization: "<<time_diag<<" s"<<std::endl;
1711 *m_logStream<<"time spent for post-solving: "<<time_solve<<" s"<<std::endl;
1712 }
1713 }
1714 else {
1715 ATH_MSG_ERROR("Problem in diagonalization. Solving skipped.");
1716 if(m_logStream)
1717 *m_logStream<<"time spent for diagonalization: "<<time_diag<<" s"<<std::endl;
1718 }
1719
1720 if(m_logStream) {
1721 *m_logStream<<"total time spent in solve: "<<time_diag+time_solve<<" s"<<std::endl;
1722 *m_logStream<<std::endl;
1723 }
1724
1725 delete d2Chi2;
1726 delete aBetterMat;
1727 delete aBetterVec;
1728
1729 // need to do this since success return value from Lapack is 0
1730 // and from solveLapack() it is 1
1731 if (info==0)
1732 info = 1;
1733
1734 return info;
1735 }
void postSolvingLapack(AlVec *dChi2, AlSymMat *d2Chi2, AlVec &w, AlMat &z, int size)
Gaudi::Property< bool > m_calDet
Definition MatrixTool.h:174

◆ solveLocal()

int MatrixTool::solveLocal ( )
private

Definition at line 505 of file MatrixTool.cxx.

506 {
507 ATH_MSG_INFO("solving using Local method");
508 if(m_logStream) {
509 *m_logStream<<"*************************************************************"<<std::endl;
510 *m_logStream<<"************** solving using Local method *****************"<<std::endl;
511 *m_logStream<<"*************************************************************"<<std::endl;
512 }
513
514 int totalNDoF(0);
515 double totalChi2(0.);
516
517 const AlignModuleList * alignModules = m_alignModuleTool->alignModules1D();
518
519 AlignModuleList::const_iterator imod = alignModules->begin();
520 AlignModuleList::const_iterator imod_end = alignModules->end();
521 for( ; imod!=imod_end; ++imod) {
522
523 AlignModule * module = *imod;
524
525 ATH_MSG_INFO("Solving for module: "<<module->name());
526
527 DataVector<AlignPar> * alignPars = m_alignModuleTool->getAlignPars(module);
528
529 int thisNDoF = alignPars->size();
530
531 CLHEP::HepSymMatrix d2Chi2(thisNDoF,0);
532 CLHEP::HepVector dChi2(thisNDoF,0);
533 for (int i=0;i<thisNDoF;++i) {
534 int ipar = alignPars->at(i)->index();
535 dChi2[i] = (*m_bigvector)[ipar];
536 for (int j=0;j<thisNDoF;++j) {
537 int jpar = alignPars->at(j)->index();
538 d2Chi2[i][j] = (*m_bigmatrix)[ipar][jpar];
539 }
540 }
541
542 ATH_MSG_DEBUG("First derivatives:" << dChi2);
543 ATH_MSG_DEBUG("Second derivatives:" << d2Chi2);
544
545 if(module->nHits() < m_minNumHits || module->nTracks() < m_minNumTrks) {
546 ATH_MSG_INFO("Not enough hits in module \'"<<module->name()<<"\': "
547 <<module->nHits()<<" < "<<m_minNumHits<<" or "
548 <<module->nTracks()<<" < "<<m_minNumTrks
549 <<". Skipping...");
550
551 // print module summary even when solving is not done
552 if(m_logStream)
553 printModuleSolution(*m_logStream,module,nullptr);
554
555 continue;
556 }
557
558 ATH_MSG_DEBUG("First derivatives:" << dChi2);
559 ATH_MSG_DEBUG("Second derivatives:" << d2Chi2);
560
561
562 totalNDoF += thisNDoF;
563
564 CLHEP::HepSymMatrix cov(d2Chi2);
565
566 // invert the matrix
567 int ierr = 0;
568 cov.invert(ierr);
569 if(ierr>0)
570 ATH_MSG_WARNING("CLHEP inversion status flag = "<<ierr);
571 else
572 ATH_MSG_DEBUG("CLHEP inversion status flag = "<<ierr);
573
574 // calculate corrections
575 CLHEP::HepVector delta = cov * dChi2;
576 //ATH_MSG_DEBUG("d2Chi2: "<<d2Chi2);
577 //ATH_MSG_DEBUG("cov: "<<cov);
578 //ATH_MSG_DEBUG("d2Chi2*cov: "<<d2Chi2*cov);
579 //ATH_MSG_DEBUG("dChi2: "<<dChi2);
580 ATH_MSG_DEBUG("Result: "<<delta);
581
582 ATH_MSG_DEBUG("Alignment constants:");
583 for (int idof=0;idof<thisNDoF;++idof) {
584 AlignPar * alignPar = alignPars->at(idof);
585
586 double sigma = alignPar->sigma();
587 double param = -delta[idof] * sigma;
588 double err = std::sqrt(2.*std::fabs(cov[idof][idof])) * sigma;
589
590 ATH_MSG_DEBUG(idof <<" : "<< param << " +/- "<< err);
591 ATH_MSG_DEBUG("cov("<<idof<<")="<<cov[idof][idof]<<", sigma: "<<sigma);
592 ATH_MSG_DEBUG("init par: "<<alignPar->initPar());
593
594 ATH_MSG_DEBUG("Filling constants obtained using Local method");
595 alignPar->setPar(param,err);
596 ATH_MSG_DEBUG("set param to "<<param<<" for alignPar "<<alignPar);
597 ATH_MSG_DEBUG(*alignPar);
598 }
599
600 if(m_logStream) {
601 printModuleSolution(*m_logStream,module,&cov);
602
603 *m_logStream<<"CLHEP inversion status flag = "<<ierr<<std::endl;
604
605 // calculate chi2 of the alignment change
606 double chi2 = d2Chi2.similarity(delta) * .5;
607 totalChi2 += chi2;
608 *m_logStream<<"delta(chi2) of the alignment change : "<<chi2<<" / "<<thisNDoF<<std::endl;
609 }
610 }
611
612 if(m_logStream) {
613 *m_logStream<<"--------------------------------------------------------------------------------"<<std::endl;
614 *m_logStream<<"Total delta(chi2) of the alignment change from the local method : "<<totalChi2<<" / "<<totalNDoF<<std::endl;
615 *m_logStream<<std::endl;
616 }
617
618 return 1;
619 }

◆ solveROOT()

int MatrixTool::solveROOT ( )
private

Definition at line 137 of file MatrixTool.cxx.

138 {
139 ATH_MSG_INFO("solving Global using ROOT");
140 if(m_logStream) {
141 *m_logStream<<"*************************************************************"<<std::endl;
142 *m_logStream<<"************** solving using Global method ****************"<<std::endl;
143 *m_logStream<<"************** using ROOT ****************"<<std::endl;
144 *m_logStream<<"*************************************************************"<<std::endl;
145 }
146
147 // start measuring time
148 clock_t starttime = clock();
149
150 DataVector<AlignPar>* alignParList = m_alignModuleTool->alignParList1D();
151 //const AlignModuleList* alignModules = m_alignModuleTool->alignModules1D();
152
153 int nDoF=m_alignModuleTool->nAlignParameters();
154
155 // some debugging output
156 if (msgLvl(MSG::VERBOSE)) {
157 msg(MSG::VERBOSE)<<"dumping matrix and vector to screen"<<endmsg;
158 for (int i=0;i<nDoF;i++)
159 for (int j=0;j<nDoF;j++)
160 //if (std::fabs((*m_bigmatrix)[i][j])>.0001)
161 msg(MSG::VERBOSE)<<i<<", "<<j<<" : "<<(*m_bigmatrix)[i][j] <<endmsg;
162
163 for (int i=0;i<nDoF;i++)
164 msg(MSG::VERBOSE)<<i <<" : "<<(*m_bigvector)[i]<<endmsg;
165 }
166
167 // get rescaled first and second derivatives
168 double * secderiv = new double[m_aNDoF*m_aNDoF];
169 double * firstderiv = new double[m_aNDoF];
170 for (int iActive=0;iActive<m_aNDoF;iActive++) {
171 int i = m_activeIndices[iActive];
172 firstderiv[iActive] = (*m_bigvector)[i];
173 for (int jActive=0;jActive<m_aNDoF;jActive++) {
174 int j = m_activeIndices[jActive];
175 secderiv[iActive*m_aNDoF+jActive] = (*m_bigmatrix)[i][j];
176 }
177 }
178
179 // attention, the dimension of matrix a and b is m_aNDoF not nDoF,
180 // this means some alignment parameters have not been calculated
181 // if the corresponding modules did not satify the select cut
182
183 TMatrixDSym a(m_aNDoF,secderiv);
184 TVectorD b(m_aNDoF,firstderiv);
185
186 if(msgLvl(MSG::DEBUG)) {
187 msg(MSG::DEBUG)<<"First derivatives:"<<endmsg;
188 b.Print();
189 msg(MSG::DEBUG)<<"Second derivatives:"<<endmsg;
190 a.Print();
191 }
192
193 TDecompBK c(a);
194 Bool_t status;
195 TMatrixDSym ainv(c.Invert(status));
196
197 TVectorD r(b.GetNrows());
198 if(status)
199 r = c.Solve(b,status);
200
201 // stop measuring time
202 clock_t stoptime = clock();
203 double totaltime = (stoptime-starttime)/double(CLOCKS_PER_SEC);
204 ATH_MSG_INFO("Time spent in solveROOT: "<<totaltime<<" s");
205
206 if(!status) {
207 msg(MSG::ERROR)<<"ROOT inversion failed"<<endmsg;
208 if(m_logStream) {
209 *m_logStream<<"ROOT inversion failed"<<std::endl;
210 *m_logStream<<std::endl;
211 }
212 }
213 else {
214 ATH_MSG_INFO("ROOT inversion ok");
215
216 ATH_MSG_DEBUG("Alignment constants:");
217 for (int iAdof=0;iAdof<m_aNDoF;iAdof++) {
218
219 int idof = m_activeIndices[iAdof];
220 AlignPar * alignPar=(*alignParList)[idof];
221
222 double sigma = alignPar->sigma();
223 double param = -r[iAdof] * sigma;
224 double err = std::sqrt(2.*std::fabs(ainv(iAdof,iAdof))) * sigma;
225
226 ATH_MSG_DEBUG(iAdof <<" : "<< param << " +/- "<< err);
227 ATH_MSG_DEBUG("ainv("<<iAdof<<")="<<ainv(iAdof,iAdof)<<", sigma: "<<sigma);
228 ATH_MSG_DEBUG("init par: "<<alignPar->initPar());
229 alignPar->setPar(param,err);
230 ATH_MSG_DEBUG("set param to "<<param<<" for alignPar "<<alignPar);
231 ATH_MSG_DEBUG(*(*alignParList)[idof]);
232 }
233
234 if(m_logStream)
235 {
236 *m_logStream<<"ROOT inversion ok"<<std::endl;
237
239
240 // norm of first derivative
241 *m_logStream<<"norm of first derivative : "<<sqrt(b.Norm2Sqr())<<std::endl;
242
243 // distance to solution
244 double dist = sqrt( ( b - (a * r) ).Norm2Sqr() );
245 *m_logStream<<"distance to solution : "<<dist<<std::endl;
246
247 // calculate chi2 of the alignment change
248 double chi2 = a.Similarity(r) * .5;
249 *m_logStream<<"delta(chi2) of the alignment change : "<<chi2<<" / "<<m_aNDoF<<std::endl;
250
251 // time spent here
252 *m_logStream<<"time spent in solve : "<<totaltime<<" s"<<std::endl;
253 }
254 }
255
256 delete [] secderiv;
257 delete [] firstderiv;
258
259 return 1;
260 }
static Double_t a
int r
Definition globals.cxx:22

◆ solveSparseEigen()

int MatrixTool::solveSparseEigen ( )
private

Definition at line 2156 of file MatrixTool.cxx.

2157 {
2158 ATH_MSG_INFO("solving Global using SparseEigen");
2159 if(m_logStream) {
2160 *m_logStream<<"*************************************************************"<<std::endl;
2161 *m_logStream<<"************** solving using Global method ****************"<<std::endl;
2162 *m_logStream<<"************** using SparseEigen ****************"<<std::endl;
2163 *m_logStream<<"*************************************************************"<<std::endl;
2164 }
2165
2166 // start measuring time
2167 clock_t starttime = clock();
2168
2169 DataVector<AlignPar> * alignParList = m_alignModuleTool->alignParList1D();
2170
2171 AlSpaMat * ABetterMat = nullptr;
2172 bool isCopy = false;
2173 if ( dynamic_cast<AlSymMat*>(m_bigmatrix) ) {
2174 ATH_MSG_INFO("Converting Matrix Format for fast solving");
2175 ABetterMat = new AlSpaMat(*(dynamic_cast<AlSymMat*>(m_bigmatrix)));
2176 isCopy = true;
2177 }
2178 else if ( dynamic_cast<AlSpaMat*>(m_bigmatrix) ) {
2179 ATH_MSG_INFO("Matrix format native to the fast solving");
2180 ABetterMat = (dynamic_cast<AlSpaMat*>(m_bigmatrix));
2181 }
2182 else {
2183 ATH_MSG_ERROR("Cannot cast to neither AlSymMat nor AlSpaMat");
2184 return 0;
2185 }
2186
2187 ATH_MSG_DEBUG("checking active indices");
2188
2189 // use const matrix when checking for non-zero elements to avoid
2190 // filling of the whole matrix
2191 const AlSpaMat * chkMatrix = ABetterMat;
2192
2193 AlSpaMat * aBetterMat = new AlSpaMat(m_aNDoF);
2194 AlVec * aBetterVec = new AlVec(m_aNDoF);
2195
2196
2197 for (int iActive=0;iActive<m_aNDoF;iActive++) {
2198 int i = m_activeIndices[iActive];
2199 (*aBetterVec)[iActive] = (*m_bigvector)[i];
2200 for (int jActive=0;jActive<m_aNDoF;jActive++) {
2201 int j = m_activeIndices[jActive];
2202 // only fill if non-zero !!!
2203 if ( (*chkMatrix)[iActive][jActive] != 0. )
2204 (*aBetterMat)[iActive][jActive]=(*ABetterMat)[i][j];
2205 }
2206 }
2207
2208 // store original vector for cross-checks
2209 AlVec origVec(*aBetterVec);
2210
2211 ATH_MSG_DEBUG("running the solving");
2212
2213 // solve
2214 int info = (*aBetterMat).SolveWithEigen(*aBetterVec);
2215
2216 if(info == 0) {
2217 ATH_MSG_INFO("SolveWithEigen solving OK");
2218 if(m_logStream)
2219 *m_logStream<<"SolveWithEigen solving OK."<<std::endl;
2220 }
2221 else {
2222 ATH_MSG_ERROR( "SolveWithEigen error code (0 if OK) = "<<info );
2223 if(m_logStream)
2224 *m_logStream<<"SolveWithEigen error code (0 if OK) = "<<info<<std::endl;
2225 }
2226
2227 if( isCopy )
2228 delete ABetterMat;
2229 ABetterMat = nullptr;
2230
2231 // stop measuring time
2232 clock_t stoptime = clock();
2233 double totaltime = (stoptime-starttime)/double(CLOCKS_PER_SEC);
2234 ATH_MSG_INFO("Time spent in SolveWithEigen: "<<totaltime<<" s");
2235
2236 ATH_MSG_DEBUG("Alignment constants:");
2237 // compute alignment corrections (translations in mm and rotations in rad)
2238 // for solveSparseEigen variances are not calculated
2239 for(int i=0; i<m_aNDoF; i++) {
2240
2241 double param = -(*aBetterVec)[i];
2242 double err = 0.;
2243
2244 int idof = m_activeIndices[i];
2245 AlignPar * alignPar=(*alignParList)[idof];
2246
2247 // undo the sigma scaling
2248 double sigma = alignPar->sigma();
2249 param *= sigma;
2250
2251 ATH_MSG_DEBUG(i <<" : "<< param << " +/- "<< err);
2252 ATH_MSG_DEBUG("sigma: "<<sigma);
2253 ATH_MSG_DEBUG("init par: "<<alignPar->initPar());
2254 alignPar->setPar(param, err);
2255 ATH_MSG_DEBUG("set param to "<<param<<" for alignPar "<<alignPar);
2256 ATH_MSG_DEBUG(*(*alignParList)[idof]);
2257 }
2258
2259 if(m_logStream) {
2260 CLHEP::HepSymMatrix * cov = nullptr;
2262
2263 // norm of first derivative
2264 *m_logStream<<"norm of first derivative : "<<origVec.norm()<<std::endl;
2265
2266 // distance to solution
2267 double dist = ( (*aBetterMat) * (*aBetterVec) - origVec ).norm();
2268 *m_logStream<<"distance to solution : "<<dist<<std::endl;
2269
2270 // calculate chi2 of the alignment change
2271 double chi2 = (*aBetterVec) * (*aBetterMat) * (*aBetterVec) * .5;
2272 *m_logStream<<"delta(chi2) of the alignment change : "<<chi2<<" / "<<m_aNDoF<<std::endl;
2273
2274 // time spent here
2275 *m_logStream<<"time spent in solve : "<<totaltime<<" s"<<std::endl;
2276 }
2277
2278 delete aBetterMat;
2279 delete aBetterVec;
2280
2281 // need to do this since success return value from Lapack is 0
2282 // and from SolveWithEigen() it is 1
2283 if (info==0)
2284 info = 1;
2285
2286 return info;
2287 }

◆ spuriousRemoval()

StatusCode MatrixTool::spuriousRemoval ( )
private

Definition at line 1738 of file MatrixTool.cxx.

1739 {
1740
1741 // copied from SiGlobalChi2Algs
1742 ATH_MSG_DEBUG("in spuriousRemoval");
1743
1744 // compute determinant before resizing
1745 if (m_calDet) {
1746 const double tol = 1.e-20;
1747 double determ = m_bigmatrix->determinant();
1748 ATH_MSG_INFO("Determinant: " << determ);
1749 if (std::fabs(determ) < tol)
1750 ATH_MSG_WARNING("Matrix is singular!");
1751 }
1752
1753 // fill vector with modules that need to be removed
1754 int fillvecmods=fillVecMods();
1755 if (fillvecmods==0) {
1756
1757 //ATH_MSG_INFO(" No resize needed (NhitsCut = "
1758 // << m_hitscut << ")");
1759
1760 if (msgLvl(MSG::DEBUG)) {
1761 //m_bigmatrix->Write("bigmatrix.txt", false, m_wSqMatrix, m_scale,
1762 // MatVersion);
1763 //m_bigvector->Write("bigvector.txt", false, m_scale, m_modcodemap,
1764 // VecVersion, m_alignProcessLevel, m_fitParam);
1765 }
1766
1767 return StatusCode::SUCCESS;
1768 }
1769 else if (fillvecmods==2)
1770 return StatusCode::FAILURE;
1771
1772 /* this is a bit difficult to implement for now....
1773
1774 // remove matrix/vector elements
1775 int cont=0;
1776 ModuleIndexMap::iterator itcode;
1777 ATH_MSG_INFO("Eliminating module...");
1778
1779
1780
1781 for (std::vector<int>::const_iterator it=m_dropmods.begin();
1782 it!= m_dropmods.end(); ++it) {
1783
1784 itcode = m_modcodemap.find((*it));
1785
1786 if (itcode == m_modcodemap.end()) {
1787 ATH_MSG_WARNING("Could not find module " << *it << " in map.");
1788 return StatusCode::FAILURE;
1789 }
1790
1791 ATH_MSG_INFO(" - Removing mcode: " << (itcode->second)
1792 << " (old index: " << (*it) << " -> new index: " << (*it)-cont << ")");
1793
1794 m_bigmatrix->RemoveModule((*it)-cont);
1795 m_bigvector->RemoveModule((*it)-cont);
1796 cont++;
1797 }
1798 ATH_MSG_INFO("Modules removed from the Matrix and from the Vector Successfully!");
1799 ATH_MSG_DEBUG("(DoF: " << nDoF << ")");
1800
1801 // resizing...
1802 ATH_MSG_INFO("Resizing the bigvector in memory...");
1803 m_bigvector->reSize(nDoF-6*m_dropmods.size());
1804 ATH_MSG_INFO("Resizing the bigmatrix in memory...");
1805 m_bigmatrix->reSize(nDoF-6*m_dropmods.size());
1806
1807 ATH_MSG_INFO(m_dropmods.size() << " modules eliminated from the matrix, i.e, "
1808 << 6*m_dropmods.size() << " DoFs");
1809 ATH_MSG_INFO(nDoF/6 - m_dropmods.size() << " modules to align (" << nDoF-6*m_dropmods.size() << " DoFs)");
1810 ATH_MSG_INFO("New bigmatrix size is: " << m_bigmatrix->size());
1811 ATH_MSG_INFO("New bigvector size is: " << (*m_bigvector).size());
1812
1813 NDoF = nDoF-6*(m_dropmods.size());
1814
1815 // resize (update) nDoF variable
1816 nDoF = NDoF;
1817
1818 // Resizing vectors to store results
1819 m_alignPar->reSize(nDoF);
1820 m_alignSqErr->reSize(nDoF);
1821
1822 // Fill the mapping module map and updating m_modcodemap
1823 UpdateModMap();
1824 */
1825
1826 return StatusCode::SUCCESS;
1827 }
static int fillVecMods()

◆ storeInTFile()

void MatrixTool::storeInTFile ( const TString & filename)

Store Files in a tfile.

Definition at line 734 of file MatrixTool.cxx.

735 {
736 //Store reults in a single TFile....
737 //Including Matrix Vector Hitmap.. Soluton EVs etc.
738
739 ATH_MSG_DEBUG("Writing Results to a TFile");
740
741 DataVector<AlignPar> * alignParList = m_alignModuleTool->alignParList1D();
742 int nDoF = alignParList->size();
743
744 TMatrixDSparse* myTMatrix = m_bigmatrix->makeTMatrix();
745
746 ATH_MSG_DEBUG( "Created TMatrixDSparse" );
747
748
749 double *val = new double[nDoF];
750 for (int i=0;i<nDoF;i++) {
751 val[i] = (*m_bigvector)[i];
752 }
753
754 TVectorD myTVector(nDoF, val);
755 delete [] val;
756
757 ATH_MSG_DEBUG( "Created TVectorD" );
758
759
760 const AlignModuleList * moduleList = m_alignModuleTool->alignModules1D();
761 int nModules = moduleList->size();
762
763 double *hitmapA = new double[nModules];
764 double *hitmapB = new double[nModules];
765 AlignModuleList::const_iterator imod = moduleList->begin();
766 AlignModuleList::const_iterator imod_end = moduleList->end();
767 int index(0);
768 for(; imod != imod_end; ++imod) {
769 AlignModule * module = *imod;
770 hitmapA[index] = (double)module->nHits();
771 hitmapB[index] = (double)module->nTracks();
772 index++;
773 }
774
775 TVectorD hitmapHits(nModules, hitmapA);
776 TVectorD hitmapTracks(nModules, hitmapB);
777
778 delete [] hitmapA;
779 delete [] hitmapB;
780
781
782
783 TFile myFile(filename,"recreate");
784 hitmapHits.Write("Hits");
785 hitmapTracks.Write("Tracks");
786 myTMatrix->Write("Matrix");
787 myTVector.Write("Vector");
788
789 TVectorD scale(1, &m_scale) ;
790 scale.Write("Scale");
791
792
793 double *moduleInfoA = new double[nDoF];//unsigned long long
794 double *dofInfoA = new double[nDoF];//int
795
796 if (sizeof(unsigned long long) != sizeof(double))
797 ATH_MSG_ERROR("Module Identifiers will not be saved. sizeof(double)!=sizeof(ulonglong)");
798 else{
799
800 DataVector<AlignPar>* alignPars = m_alignModuleTool->alignParList1D();
801 for (int i=0;i<(int)alignPars->size();i++) {
802 //Do a direct memory copy to store unsigned long long in the momory of a double
803 double target;
804 uint64_t id = (*alignPars)[i]->alignModule()->identify().get_compact();
805 memcpy(&target, &id, sizeof(target));
806 moduleInfoA[i]=target;
807 //moduleInfoB[i]=(*alignPars)[i]->alignModule()->name();
808 uint64_t dof = (*alignPars)[i]->paramType();
809 memcpy(&target, &dof, sizeof(target));
810 dofInfoA[i]=target;
811 }
812
813 TVectorD moduleIDs(nDoF, moduleInfoA) ;
814 TVectorD moduleDoFs(nDoF,dofInfoA);
815 delete [] moduleInfoA;
816 moduleIDs.Write("ModuleID");
817 moduleDoFs.Write("dof");
818 }
819
820 myFile.Write();
821 myFile.Close();
822
823 delete myTMatrix;
824 ATH_MSG_DEBUG("Finshed writing TFILE");
825
826 }

◆ sysInitialize()

virtual StatusCode AthCommonDataStore< AthCommonMsg< AlgTool > >::sysInitialize ( )
overridevirtualinherited

Perform system initialization for an algorithm.

We override this to declare all the elements of handle key arrays at the end of initialization. See comments on updateVHKA.

Reimplemented in asg::AsgMetadataTool, AthCheckedComponent< AthAlgTool >, AthCheckedComponent<::AthAlgTool >, and DerivationFramework::CfAthAlgTool.

◆ sysStart()

virtual StatusCode AthCommonDataStore< AthCommonMsg< AlgTool > >::sysStart ( )
overridevirtualinherited

Handle START transition.

We override this in order to make sure that conditions handle keys can cache a pointer to the conditions container.

◆ updateVHKA()

void AthCommonDataStore< AthCommonMsg< AlgTool > >::updateVHKA ( Gaudi::Details::PropertyBase & )
inlineinherited

Definition at line 308 of file AthCommonDataStore.h.

308 {
309 // debug() << "updateVHKA for property " << p.name() << " " << p.toString()
310 // << " size: " << m_vhka.size() << endmsg;
311 for (auto &a : m_vhka) {
313 for (auto k : keys) {
314 k->setOwner(this);
315 }
316 }
317 }
std::vector< SG::VarHandleKeyArray * > m_vhka

◆ writeHitmap()

void MatrixTool::writeHitmap ( )
private

Definition at line 2290 of file MatrixTool.cxx.

2291 {
2292 ATH_MSG_INFO("writing the hitmap to file");
2293
2294 const AlignModuleList * moduleList = m_alignModuleTool->alignModules1D();
2295 int nModules = moduleList->size();
2296
2297 AlMat hitmap(nModules,2);
2298 AlignModuleList::const_iterator imod = moduleList->begin();
2299 AlignModuleList::const_iterator imod_end = moduleList->end();
2300 int index(0);
2301 for(; imod != imod_end; ++imod) {
2302 AlignModule * module = *imod;
2303 hitmap[index][0] = module->nHits();
2304 hitmap[index][1] = module->nTracks();
2305 index++;
2306 }
2307
2308 // Set Path for the hitmap matrix
2309 hitmap.SetPathBin(m_pathbin.value()+m_prefixName.value());
2310 hitmap.SetPathTxt(m_pathtxt.value()+m_prefixName.value());
2311
2312 StatusCode sc = hitmap.Write("hitmap.bin",true); // write the hitmap matrix
2313
2314 if (sc!=StatusCode::SUCCESS)
2315 ATH_MSG_ERROR("Problem writing hitmap matrix");
2316
2317 if (m_writeHitmapTxt) {
2318 sc = hitmap.Write("hitmap.txt", false, 0);
2319 if (sc!=StatusCode::SUCCESS)
2320 ATH_MSG_ERROR("Problem writing hitmap matrix to text file");
2321 }
2322
2323 ATH_MSG_DEBUG("hitmap written to: "<< m_pathbin.value()+m_prefixName.value() <<"hitmap.bin (.txt)");
2324 }
Gaudi::Property< bool > m_writeHitmapTxt
Definition MatrixTool.h:191

Member Data Documentation

◆ m_activeIndices

std::vector<int> Trk::MatrixTool::m_activeIndices {}
private

vector of indices which pass the min-hits cut

Definition at line 248 of file MatrixTool.h.

248{};

◆ m_Align_db_step

Gaudi::Property<float> Trk::MatrixTool::m_Align_db_step
private
Initial value:
{this, "AlignCorrDBStep", 10.,
"corr in the diagonal basis step for the third pass in the auto weak mode removal method"}

Definition at line 171 of file MatrixTool.h.

171 {this, "AlignCorrDBStep", 10.,
172 "corr in the diagonal basis step for the third pass in the auto weak mode removal method"};

◆ m_AlignIBLbutNotPixel

Gaudi::Property<bool> Trk::MatrixTool::m_AlignIBLbutNotPixel {this, "AlignIBLbutNotPixel", false}
private

Definition at line 255 of file MatrixTool.h.

256{this, "AlignIBLbutNotPixel", false};

◆ m_alignModuleMap

std::map<int,int> Trk::IMatrixTool::m_alignModuleMap
privateinherited

Definition at line 105 of file IMatrixTool.h.

◆ m_alignModuleTool

PublicToolHandle<IAlignModuleTool> Trk::MatrixTool::m_alignModuleTool
private
Initial value:
{
this, "AlignModuleTool", "Trk::AlignModuleTool/AlignModuleTool"}

Definition at line 142 of file MatrixTool.h.

142 {
143 this, "AlignModuleTool", "Trk::AlignModuleTool/AlignModuleTool"};

◆ m_AlignPixelbutNotIBL

Gaudi::Property<bool> Trk::MatrixTool::m_AlignPixelbutNotIBL {this, "AlignPixelbutNotIBL", false}
private

Definition at line 257 of file MatrixTool.h.

258{this, "AlignPixelbutNotIBL", false};

◆ m_aNDoF

int Trk::MatrixTool::m_aNDoF = 0
private

number of active DoF (size of m_activeIndices)

Definition at line 249 of file MatrixTool.h.

◆ m_bigmatrix

AlSymMatBase* Trk::MatrixTool::m_bigmatrix = nullptr
private

matrix to contain second derivative terms to be used for alignment

Definition at line 146 of file MatrixTool.h.

◆ m_bigvector

AlVec* Trk::MatrixTool::m_bigvector = nullptr
private

vector to contain first derivative terms to be used for alignment

Definition at line 149 of file MatrixTool.h.

◆ m_calculateFullCovariance

Gaudi::Property<double> Trk::MatrixTool::m_calculateFullCovariance
private
Initial value:
{this, "CalculateFullCovariance", true,
"calculate full covariance matrix for Lapack"}

Definition at line 214 of file MatrixTool.h.

215 {this, "CalculateFullCovariance", true,
216 "calculate full covariance matrix for Lapack"};

◆ m_calDet

Gaudi::Property<bool> Trk::MatrixTool::m_calDet
private
Initial value:
{this, "MatrixDet", false,
"compute bigmatrix's determinant ?"}

Definition at line 174 of file MatrixTool.h.

174 {this, "MatrixDet", false,
175 "compute bigmatrix's determinant ?"};

◆ m_DeactivateSCT_ECA_LastDisk

Gaudi::Property<bool> Trk::MatrixTool::m_DeactivateSCT_ECA_LastDisk {this, "DeactivateSCT_ECA_LastDisk", false}
private

Definition at line 260 of file MatrixTool.h.

261{this, "DeactivateSCT_ECA_LastDisk", false};

◆ m_detStore

StoreGateSvc_t AthCommonDataStore< AthCommonMsg< AlgTool > >::m_detStore
privateinherited

Pointer to StoreGate (detector store by default)

Definition at line 393 of file AthCommonDataStore.h.

◆ m_diagonalize

Gaudi::Property<bool> Trk::MatrixTool::m_diagonalize
private
Initial value:
{this, "Diagonalize", true,
"run diagonalization instead of inversion"}

Definition at line 154 of file MatrixTool.h.

154 {this, "Diagonalize", true,
155 "run diagonalization instead of inversion"};

◆ m_eigenvalueStep

Gaudi::Property<float> Trk::MatrixTool::m_eigenvalueStep
private
Initial value:
{this, "EigenvalueStep", 1e3,
"eigenvalue step for the second pass in the automatic weak mode removal method"}

Definition at line 169 of file MatrixTool.h.

169 {this, "EigenvalueStep", 1e3,
170 "eigenvalue step for the second pass in the automatic weak mode removal method"};

◆ m_eigenvaluethreshold

Gaudi::Property<double> Trk::MatrixTool::m_eigenvaluethreshold {this, "EigenvalueThreshold", 0., "cut on the minimum eigenvalue"}
private

Definition at line 156 of file MatrixTool.h.

157{this, "EigenvalueThreshold", 0., "cut on the minimum eigenvalue"};

◆ m_evtStore

StoreGateSvc_t AthCommonDataStore< AthCommonMsg< AlgTool > >::m_evtStore
privateinherited

Pointer to StoreGate (event store by default)

Definition at line 390 of file AthCommonDataStore.h.

◆ m_inputHitmapFiles

Gaudi::Property<std::vector<std::string> > Trk::MatrixTool::m_inputHitmapFiles
private
Initial value:
{this, "InputHitmapFiles", {"hitmap.bin"},
"input binary files containing the hitmaps"}

Definition at line 240 of file MatrixTool.h.

241 {this, "InputHitmapFiles", {"hitmap.bin"},
242 "input binary files containing the hitmaps"};

◆ m_inputMatrixFiles

Gaudi::Property<std::vector<std::string> > Trk::MatrixTool::m_inputMatrixFiles
private
Initial value:
{this, "InputMatrixFiles", {"matrix.bin"},
"input binary files containing matrix terms"}

Definition at line 233 of file MatrixTool.h.

234 {this, "InputMatrixFiles", {"matrix.bin"},
235 "input binary files containing matrix terms"};

◆ m_inputTFiles

Gaudi::Property<std::vector<std::string> > Trk::MatrixTool::m_inputTFiles
private
Initial value:
{this, "InputTFiles", {"AlignmentTFile.root"},
"input binary files containing matrix terms"}

Definition at line 244 of file MatrixTool.h.

245 {this, "InputTFiles", {"AlignmentTFile.root"},
246 "input binary files containing matrix terms"};

◆ m_inputVectorFiles

Gaudi::Property<std::vector<std::string> > Trk::MatrixTool::m_inputVectorFiles
private
Initial value:
{this, "InputVectorFiles", {"vector.bin"},
"input binary files containing vector terms"}

Definition at line 236 of file MatrixTool.h.

237 {this, "InputVectorFiles", {"vector.bin"},
238 "input binary files containing vector terms"};

◆ m_logStream

std::ostream* Trk::IMatrixTool::m_logStream
protectedinherited

logfile output stream

Definition at line 98 of file IMatrixTool.h.

◆ m_maxReadErrors

Gaudi::Property<int> Trk::MatrixTool::m_maxReadErrors
private
Initial value:
{this, "MaxReadErrors", 10,
"maximum number of reading TFile errors"}

Definition at line 251 of file MatrixTool.h.

251 {this, "MaxReadErrors", 10,
252 "maximum number of reading TFile errors"};

◆ m_minNumHits

Gaudi::Property<int> Trk::MatrixTool::m_minNumHits
private
Initial value:
{this, "MinNumHitsPerModule", 0,
"cut on the minimum number of hits per module"}

Definition at line 163 of file MatrixTool.h.

163 {this, "MinNumHitsPerModule", 0,
164 "cut on the minimum number of hits per module"};

◆ m_minNumTrks

Gaudi::Property<int> Trk::MatrixTool::m_minNumTrks
private
Initial value:
{this, "MinNumTrksPerModule", 0,
"cut on the minimum number of tracks per module"}

Definition at line 165 of file MatrixTool.h.

165 {this, "MinNumTrksPerModule", 0,
166 "cut on the minimum number of tracks per module"};

◆ m_modcut

Gaudi::Property<int> Trk::MatrixTool::m_modcut
private
Initial value:
{this, "ModCut", 0,
"cut on the weak modes which number is <par_modcut"}

Definition at line 161 of file MatrixTool.h.

161 {this, "ModCut", 0,
162 "cut on the weak modes which number is <par_modcut"};

◆ m_nentries

int Trk::IMatrixTool::m_nentries
privateinherited

Definition at line 106 of file IMatrixTool.h.

◆ m_nHits

int Trk::IMatrixTool::m_nHits
protectedinherited

Definition at line 100 of file IMatrixTool.h.

◆ m_nMeasurements

int Trk::IMatrixTool::m_nMeasurements
protectedinherited

Definition at line 102 of file IMatrixTool.h.

◆ m_nTracks

int Trk::IMatrixTool::m_nTracks
protectedinherited

Definition at line 101 of file IMatrixTool.h.

◆ m_pathbin

Gaudi::Property<std::string> Trk::MatrixTool::m_pathbin
private
Initial value:
{this, "PathBinName", "./",
"path binary files (in/out)"}

Definition at line 218 of file MatrixTool.h.

218 {this, "PathBinName", "./",
219 "path binary files (in/out)"};

◆ m_pathtxt

Gaudi::Property<std::string> Trk::MatrixTool::m_pathtxt
private
Initial value:
{this, "PathTxtName", "./",
"path ascii files (in/out)"}

Definition at line 220 of file MatrixTool.h.

220 {this, "PathTxtName", "./",
221 "path ascii files (in/out)"};

◆ m_prefixName

Gaudi::Property<std::string> Trk::MatrixTool::m_prefixName
private
Initial value:
{this, "PrefixName", "",
"prefix string to filenames"}

Definition at line 222 of file MatrixTool.h.

222 {this, "PrefixName", "",
223 "prefix string to filenames"};

◆ m_pullcut

Gaudi::Property<float> Trk::MatrixTool::m_pullcut
private
Initial value:
{this, "PullCut", 1.0,
"pull cut for the automatic weak mode removal method"}

Definition at line 167 of file MatrixTool.h.

167 {this, "PullCut", 1.0,
168 "pull cut for the automatic weak mode removal method"};

◆ m_readHitmaps

Gaudi::Property<bool> Trk::MatrixTool::m_readHitmaps
private
Initial value:
{this, "ReadHitmaps", false,
"accumulate hitymap from files"}

Definition at line 193 of file MatrixTool.h.

193 {this, "ReadHitmaps", false,
194 "accumulate hitymap from files"};

◆ m_readTFiles

Gaudi::Property<bool> Trk::MatrixTool::m_readTFiles
private
Initial value:
{this, "ReadTFile", false,
"if True then files will be read from TFiles instead of Binary files"}

Definition at line 198 of file MatrixTool.h.

198 {this, "ReadTFile", false,
199 "if True then files will be read from TFiles instead of Binary files"};

◆ m_Remove_IBL_Rx

Gaudi::Property<bool> Trk::MatrixTool::m_Remove_IBL_Rx {this, "Remove_IBL_Rx", false}
private

Definition at line 275 of file MatrixTool.h.

275{this, "Remove_IBL_Rx", false};

◆ m_Remove_IBL_Ry

Gaudi::Property<bool> Trk::MatrixTool::m_Remove_IBL_Ry {this, "Remove_IBL_Ry", false}
private

Definition at line 276 of file MatrixTool.h.

276{this, "Remove_IBL_Ry", false};

◆ m_Remove_IBL_Rz

Gaudi::Property<bool> Trk::MatrixTool::m_Remove_IBL_Rz {this, "Remove_IBL_Rz", false}
private

Definition at line 277 of file MatrixTool.h.

277{this, "Remove_IBL_Rz", false};

◆ m_Remove_IBL_Tx

Gaudi::Property<bool> Trk::MatrixTool::m_Remove_IBL_Tx {this, "Remove_IBL_Tx", false}
private

Definition at line 272 of file MatrixTool.h.

272{this, "Remove_IBL_Tx", false};

◆ m_Remove_IBL_Ty

Gaudi::Property<bool> Trk::MatrixTool::m_Remove_IBL_Ty {this, "Remove_IBL_Ty", false}
private

Definition at line 273 of file MatrixTool.h.

273{this, "Remove_IBL_Ty", false};

◆ m_Remove_IBL_Tz

Gaudi::Property<bool> Trk::MatrixTool::m_Remove_IBL_Tz {this, "Remove_IBL_Tz", false}
private

Definition at line 274 of file MatrixTool.h.

274{this, "Remove_IBL_Tz", false};

◆ m_Remove_Pixel_Rx

Gaudi::Property<bool> Trk::MatrixTool::m_Remove_Pixel_Rx {this, "Remove_Pixel_Rx", false}
private

Definition at line 267 of file MatrixTool.h.

267{this, "Remove_Pixel_Rx", false};

◆ m_Remove_Pixel_Ry

Gaudi::Property<bool> Trk::MatrixTool::m_Remove_Pixel_Ry {this, "Remove_Pixel_Ry", false}
private

Definition at line 268 of file MatrixTool.h.

268{this, "Remove_Pixel_Ry", false};

◆ m_Remove_Pixel_Rz

Gaudi::Property<bool> Trk::MatrixTool::m_Remove_Pixel_Rz {this, "Remove_Pixel_Rz", false}
private

Definition at line 269 of file MatrixTool.h.

269{this, "Remove_Pixel_Rz", false};

◆ m_Remove_Pixel_Tx

Gaudi::Property<bool> Trk::MatrixTool::m_Remove_Pixel_Tx {this, "Remove_Pixel_Tx", false}
private

Definition at line 264 of file MatrixTool.h.

264{this, "Remove_Pixel_Tx", false};

◆ m_Remove_Pixel_Ty

Gaudi::Property<bool> Trk::MatrixTool::m_Remove_Pixel_Ty {this, "Remove_Pixel_Ty", false}
private

Definition at line 265 of file MatrixTool.h.

265{this, "Remove_Pixel_Ty", false};

◆ m_Remove_Pixel_Tz

Gaudi::Property<bool> Trk::MatrixTool::m_Remove_Pixel_Tz {this, "Remove_Pixel_Tz", false}
private

Definition at line 266 of file MatrixTool.h.

266{this, "Remove_Pixel_Tz", false};

◆ m_removeSpurious

Gaudi::Property<double> Trk::MatrixTool::m_removeSpurious
private
Initial value:
{this, "RemoveSpurious", false,
"run spurious removal"}

Definition at line 211 of file MatrixTool.h.

211 {this, "RemoveSpurious", false,
212 "run spurious removal"};

◆ m_runLocal

Gaudi::Property<bool> Trk::MatrixTool::m_runLocal
private
Initial value:
{this, "RunLocalMethod", true,
"Run solving using Local method"}

Definition at line 201 of file MatrixTool.h.

201 {this, "RunLocalMethod", true,
202 "Run solving using Local method"};

◆ m_scalaMatName

Gaudi::Property<std::string> Trk::MatrixTool::m_scalaMatName {this, "ScalapackMatrixName", "eigenvectors.bin", "Scalapack matrix name"}
private

Definition at line 228 of file MatrixTool.h.

229{this, "ScalapackMatrixName", "eigenvectors.bin", "Scalapack matrix name"};

◆ m_scalaVecName

Gaudi::Property<std::string> Trk::MatrixTool::m_scalaVecName {this, "ScalapackVectorName", "eigenvalues.bin", "Scalapack vector name"}
private

Definition at line 230 of file MatrixTool.h.

231{this, "ScalapackVectorName", "eigenvalues.bin", "Scalapack vector name"};

◆ m_scale

double Trk::MatrixTool::m_scale = -1.
private

scale for big matrix and vector normalization

Definition at line 204 of file MatrixTool.h.

◆ m_scaleMatrix

Gaudi::Property<bool> Trk::MatrixTool::m_scaleMatrix
private
Initial value:
{this, "ScaleMatrix", false,
"scale matrix by number of hits before solving"}

Definition at line 205 of file MatrixTool.h.

205 {this, "ScaleMatrix", false,
206 "scale matrix by number of hits before solving"};

◆ m_softEigenmodeCut

Gaudi::Property<double> Trk::MatrixTool::m_softEigenmodeCut
private
Initial value:
{this, "SoftEigenmodeCut", 0.,
"add constant to diagonal to effectively cut on weak eigenmodes"}

Definition at line 208 of file MatrixTool.h.

208 {this, "SoftEigenmodeCut", 0.,
209 "add constant to diagonal to effectively cut on weak eigenmodes"};

◆ m_solveOption

Gaudi::Property<int> Trk::MatrixTool::m_solveOption {this, "SolveOption", NONE, "solving option"}
private

Definition at line 159 of file MatrixTool.h.

160{this, "SolveOption", NONE, "solving option"};

◆ m_tfileName

Gaudi::Property<std::string> Trk::MatrixTool::m_tfileName {this, "TFileName", "AlignmentTFile.root", "prefix string to filenames"}
private

Definition at line 225 of file MatrixTool.h.

226{this, "TFileName", "AlignmentTFile.root", "prefix string to filenames"};

◆ m_useSparse

Gaudi::Property<bool> Trk::MatrixTool::m_useSparse {this, "UseSparse", false}
private

flag to use AlSpaMat for the big matrix (default is AlSymMat)

Definition at line 152 of file MatrixTool.h.

152{this, "UseSparse", false};

◆ m_varHandleArraysDeclared

bool AthCommonDataStore< AthCommonMsg< AlgTool > >::m_varHandleArraysDeclared
privateinherited

Definition at line 399 of file AthCommonDataStore.h.

◆ m_vhka

std::vector<SG::VarHandleKeyArray*> AthCommonDataStore< AthCommonMsg< AlgTool > >::m_vhka
privateinherited

Definition at line 398 of file AthCommonDataStore.h.

◆ m_writeEigenMat

Gaudi::Property<bool> Trk::MatrixTool::m_writeEigenMat
private
Initial value:
{this, "WriteEigenMat", true,
"write eigenvalues and eigenvectors into files ?"}

Definition at line 182 of file MatrixTool.h.

182 {this, "WriteEigenMat", true,
183 "write eigenvalues and eigenvectors into files ?"};

◆ m_writeEigenMatTxt

Gaudi::Property<bool> Trk::MatrixTool::m_writeEigenMatTxt
private
Initial value:
{this, "WriteEigenMatTxt", true,
"also write eigenvalues and eigenvectors into txt files ?"}

Definition at line 184 of file MatrixTool.h.

184 {this, "WriteEigenMatTxt", true,
185 "also write eigenvalues and eigenvectors into txt files ?"};

◆ m_writeHitmap

Gaudi::Property<bool> Trk::MatrixTool::m_writeHitmap
private
Initial value:
{this, "WriteHitmap", false,
"write hitmap into file"}

Definition at line 189 of file MatrixTool.h.

189 {this, "WriteHitmap", false,
190 "write hitmap into file"};

◆ m_writeHitmapTxt

Gaudi::Property<bool> Trk::MatrixTool::m_writeHitmapTxt
private
Initial value:
{this, "WriteHitmapTxt", false,
"write hitmap into text file"}

Definition at line 191 of file MatrixTool.h.

191 {this, "WriteHitmapTxt", false,
192 "write hitmap into text file"};

◆ m_writeMat

Gaudi::Property<bool> Trk::MatrixTool::m_writeMat
private
Initial value:
{this, "WriteMat", true,
"write big matrix and vector into files ?"}

Definition at line 178 of file MatrixTool.h.

178 {this, "WriteMat", true,
179 "write big matrix and vector into files ?"};

◆ m_writeMatTxt

Gaudi::Property<bool> Trk::MatrixTool::m_writeMatTxt
private
Initial value:
{this, "WriteMatTxt", true,
"also write big matrix and vector into txt files ?"}

Definition at line 180 of file MatrixTool.h.

180 {this, "WriteMatTxt", true,
181 "also write big matrix and vector into txt files ?"};

◆ m_writeModuleNames

Gaudi::Property<bool> Trk::MatrixTool::m_writeModuleNames
private
Initial value:
{this, "WriteModuleNames", false,
"write module name instead of Identifier to vector file"}

Definition at line 186 of file MatrixTool.h.

186 {this, "WriteModuleNames", false,
187 "write module name instead of Identifier to vector file"};

◆ m_writeTFile

Gaudi::Property<bool> Trk::MatrixTool::m_writeTFile
private
Initial value:
{this, "WriteTFile", false,
"write out files to a root file"}

Definition at line 196 of file MatrixTool.h.

196 {this, "WriteTFile", false,
197 "write out files to a root file"};

◆ m_wSqMatrix

Gaudi::Property<bool> Trk::MatrixTool::m_wSqMatrix
private
Initial value:
{this, "WriteSquareMatrix", false,
"write a triangular matrix by default (true: square format) ?"}

Definition at line 176 of file MatrixTool.h.

176 {this, "WriteSquareMatrix", false,
177 "write a triangular matrix by default (true: square format) ?"};

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