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< StoreGateSvc >	StoreGateSvc_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< IAlignModuleTool >	m_alignModuleTool
AlSymMatBase *	m_bigmatrix = nullptr
	matrix to contain second derivative terms to be used for alignment
AlVec *	m_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

enum Trk::MatrixTool::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.

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

Constructor & Destructor Documentation

MatrixTool()

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

Constructor.

Definition at line 54 of file MatrixTool.cxx.

    : IMatrixTool()
    , AthAlgTool(type,name,parent)
  {
    declareInterface<IMatrixTool>(this);
  }

~MatrixTool()

MatrixTool::~MatrixTool ( )

virtual

Virtual destructor.

Definition at line 63 of file MatrixTool.cxx.

  {
    delete m_bigmatrix;
    delete m_bigvector;
  }

Member Function Documentation

accumulateFromBinaries()

bool MatrixTool::accumulateFromBinaries

(

)

accumulates derivates from binary files

Definition at line 631 of file MatrixTool.cxx.

  {
    
    
    DataVector<AlignPar>* alignParList = m_alignModuleTool->alignParList1D();
    int nDoF=alignParList->size();
 
    std::map<int,unsigned long long> modIndexMap;
    float dummyVersion(0.);
    double totalscale=0.;
    for (int ivec=0;ivec<(int)m_inputVectorFiles.size();ivec++) {
 
      ATH_MSG_DEBUG("Reading vector "<<ivec<<" from file "<<m_inputVectorFiles[ivec]);
 
      AlVec newVector(nDoF);
      std::map<int,unsigned long long> newModIndexMap;
      newVector.SetPathBin(m_pathbin.value()+m_prefixName.value());
      newVector.SetPathTxt(m_pathtxt.value()+m_prefixName.value());
      double scale=0;
      StatusCode sc = newVector.ReadPartial(m_inputVectorFiles[ivec],scale,newModIndexMap,dummyVersion);
      totalscale += scale;
      if (sc==StatusCode::FAILURE) {
        msg(MSG::FATAL)<<"Problem reading vector from "<<m_inputVectorFiles[ivec]<<endmsg;
        return false;
      }
      if (newVector.size()!=m_bigvector->size()) {
        msg(MSG::FATAL) <<"vector wrong size!  newVector size "<<newVector.size()
                        <<", bigvector size "<<m_bigvector->size()<<endmsg;
        return false;
      }
 
      // check modIndexMaps to make sure they are the same
      if (ivec==0)
        modIndexMap = newModIndexMap;
      else if (modIndexMap!=newModIndexMap) {
        msg(MSG::FATAL)<<"module index maps don't agree!"<<endmsg;
        return false;
      }
      if (ivec>0)
        *m_bigvector += newVector;
      else
        *m_bigvector  = newVector;
    }
 
    m_scale = totalscale;
 
    AlSymMat * symBigMatrix=dynamic_cast<AlSymMat*>(m_bigmatrix);
    AlSpaMat * spaBigMatrix=dynamic_cast<AlSpaMat*>(m_bigmatrix);
 
 
    for (int imat=0;imat<(int)m_inputMatrixFiles.size();imat++) {
      ATH_MSG_DEBUG("Reading matrix "<<imat<<" from file "<<m_inputMatrixFiles[imat]);
 
      // create new matrix to read data from current file
      int nDoF=modIndexMap.size();
      bool triang;
      StatusCode sc;
      if (symBigMatrix) {
        AlSymMat newMatrix(nDoF);
        sc = newMatrix.Read(m_inputMatrixFiles[imat],nDoF,triang,dummyVersion);
        if (sc==StatusCode::SUCCESS)
          *symBigMatrix += newMatrix;
      }
      else {
        if (!spaBigMatrix) {
          throw std::logic_error("Unhandled matrix type");
        }
 
        AlSpaMat newMatrix(nDoF);
        sc = newMatrix.Read(m_inputMatrixFiles[imat],nDoF,triang,dummyVersion);
 
        if (sc==StatusCode::SUCCESS) {
          if (imat>0)
            *spaBigMatrix += newMatrix;
          else
            *spaBigMatrix  = newMatrix;
        }
      }
 
      if (sc==StatusCode::FAILURE) {
        msg(MSG::FATAL)<<"problem reading matrix from "<<m_inputMatrixFiles[imat]<<endmsg;
        return false;
      }
 
      if (!m_useSparse && triang==m_wSqMatrix) {
        ATH_MSG_WARNING("matrix not expected format!  Changing m_wSqMatrix to "<<!triang);
        m_wSqMatrix=!triang;
      }
 
    }
 
    // accumulate hitmap from hitmap files
    if(m_readHitmaps)
      readHitmaps();
 
    return true;
  }

accumulateFromFiles()

bool MatrixTool::accumulateFromFiles ( )

virtual

accumulates derivates from files.

Flag decides if it is binary or TFiles

Implements Trk::IMatrixTool.

Definition at line 618 of file MatrixTool.cxx.

  {
    
    if(m_readTFiles){ 
      ATH_MSG_INFO("Info to obtained from from TFiles");
      return accumulateFromTFiles();
    }else{ 
      ATH_MSG_INFO("Info to obtained from from Binary files");
      return accumulateFromBinaries();
    }
  }

accumulateFromTFiles()

bool MatrixTool::accumulateFromTFiles

(

)

Store Files in a tfile.

Definition at line 825 of file MatrixTool.cxx.

  {
    DataVector<AlignPar>* alignParList = m_alignModuleTool->alignParList1D();
    int nDoF=alignParList->size();
    ATH_MSG_DEBUG("OPENING TFILES");
 
    std::map<int,unsigned long long> modIndexMap;
    std::map<int,unsigned long long> DoFMap;
    double totalscale=0.;
 
    AlSymMat * symBigMatrix=dynamic_cast<AlSymMat*>(m_bigmatrix);
    AlSpaMat * spaBigMatrix=dynamic_cast<AlSpaMat*>(m_bigmatrix);
    //TMatrixDSparse *accumMatrix(0);
    AlSpaMat *accumMatrix = nullptr;
    
    const AlignModuleList * moduleList = m_alignModuleTool->alignModules1D();
    int nModules = moduleList->size();
    
    TVectorD TotalHits(nModules);
    TVectorD TotalTracks(nModules);
   
    int numberOfReadErrors = 0; 
   
    struct rusage myusage{};
    int itworked =  getrusage(RUSAGE_SELF,&myusage);
    if(itworked)//note: rusage returns zero if it succeeds!
          ATH_MSG_DEBUG("ItWorked");
 
    long intialMemUse = myusage.ru_maxrss;
 
    for (int ifile = 0; ifile < (int)m_inputTFiles.size(); ifile++) {
      if (numberOfReadErrors > m_maxReadErrors){
        msg(MSG::FATAL) << " number of errors when reading the TFiles already exceed " << m_maxReadErrors << endmsg;
        return false;
      }
   
      ATH_MSG_DEBUG("Reading File number " << ifile << ",  " << m_inputTFiles[ifile]);
      
      itworked =  getrusage(RUSAGE_SELF,&myusage);
      ATH_MSG_DEBUG("Memory usage [MB], total " << myusage.ru_maxrss/1024 << ", increase  " << (myusage.ru_maxrss-intialMemUse)/1024);
       
      TFile* myFile = TFile::Open(m_inputTFiles[ifile].c_str()); 
      
      if ( myFile->IsZombie() || !(myFile->IsOpen()) ) {
        ++numberOfReadErrors;
        ATH_MSG_ERROR( " Problem reading TFile " << m_inputTFiles[ifile] );
        continue;
      }
      
      std::map<int,unsigned long long> newModIndexMap;
      
      TVectorD* myModuleIDs;
      myModuleIDs = (TVectorD*)myFile->Get("ModuleID");
      if( !myModuleIDs ){
        ++numberOfReadErrors;
        ATH_MSG_ERROR("Modules ID not read!!!");
        continue;
      }
      
      for (int i(0); i<myModuleIDs->GetNrows(); ++i){
        //Coverting back from a double to a unvi signed long long
        double source = (*myModuleIDs)(i);
        uint64_t target;
        memcpy(&target, &source, sizeof(target));
        newModIndexMap[i]=target;
        //std::cout << i<< " " <<target <<"\n"; 
      }
      
      delete myModuleIDs;
      
      std::map<int,unsigned long long> newDoFMap;
      
      TVectorD* myDoFs;
      myDoFs = (TVectorD*)myFile->Get("dof");
      if( !myDoFs ){
        ++numberOfReadErrors;
        ATH_MSG_ERROR("DoFs not read!!!");
        continue;
      }
      
      for (int i(0); i<myDoFs->GetNrows(); ++i){
        //Coverting back from a double to a unsigned long long
        double source = (*myDoFs)(i);
        uint64_t target;
        memcpy(&target, &source, sizeof(target));
        newDoFMap[i]=target;
      }
      delete myDoFs;
      
      
      TVectorD* Scale;
      Scale = (TVectorD*)myFile->Get("Scale");
      if( !Scale ){
        ++numberOfReadErrors;
        ATH_MSG_ERROR("Scale not read!!!");
        continue;
      }
      
      double scale=(*Scale)(0);
      totalscale += scale;
      delete Scale;
      
            
      ATH_MSG_DEBUG("Reading Vector");
      TVectorD* vector = (TVectorD*)myFile->Get("Vector");
      if( !vector ){
        ++numberOfReadErrors;
        ATH_MSG_ERROR("Vector not read!!!");
        continue;
      }
      
      AlVec* newVector = new AlVec(nDoF);
      newVector->SetPathBin(m_pathbin.value()+m_prefixName.value());
      newVector->SetPathTxt(m_pathtxt.value()+m_prefixName.value());
      
      if (newVector->size()  != m_bigvector->size() ) {
         msg(MSG::FATAL) << "vector wrong size!  newVector size " << newVector->size()
                         << ", bigvector size " << m_bigvector->size()<<endmsg;
         delete newVector;
         delete vector;
         return false;
      }
      
      if (m_bigvector->size() != vector->GetNrows() ) {
         msg(MSG::FATAL) << "File vector wrong size!  File Vector size " << vector->GetNrows()
                         << ", bigvector size " << m_bigvector->size()<<endmsg;
         delete newVector;
         delete vector;
         return false;
      }
      
        
      for (int i=0;i<nDoF;i++) {
        (*newVector)[i] = (*vector)(i);
      }
      delete vector;
      
      // check modIndexMaps to make sure they are the same
      if (ifile == 0){
        DoFMap = newDoFMap;
      } else if (DoFMap!=newDoFMap) {
        msg(MSG::FATAL) << "module dofs don't agree!" << endmsg;
        return false;
      }
      
      if (ifile == 0){
         modIndexMap = newModIndexMap;
      } else if (modIndexMap!=newModIndexMap) {
         msg(MSG::FATAL) << "module index maps don't agree!" << endmsg;
         return false;
      }
      
      if (ifile>0){
        *m_bigvector += *newVector;
        delete newVector;  
      } else {
        delete m_bigvector;
        m_bigvector  = newVector;
      }
      
      
      ATH_MSG_DEBUG("Reading matrix ");
      TMatrixDSparse* matrix = (TMatrixDSparse*)myFile->Get("Matrix");
      
      if( !matrix ){
        ++numberOfReadErrors;
        ATH_MSG_ERROR("Matrix not read!!!");
        continue;
      }
      
      
      if (ifile == 0 ){ 
      
        accumMatrix = new AlSpaMat(nDoF);
        ATH_MSG_DEBUG("Matrix size b4 "<< accumMatrix->ptrMap()->size()  );
        
        //This method is ok for large matrix files... really only access the non zero elements
        for (int ii=0;ii<nDoF;ii++) {
          const TMatrixTSparseRow_const<double> myRow = (*matrix)[ii];
          int i =  myRow.GetRowIndex();
          for (int jj=0;jj<myRow.GetNindex();jj++) {
            int j = (myRow.GetColPtr())[jj];
            const double  myElement= (myRow.GetDataPtr())[jj];
            if (i<j){
              ATH_MSG_DEBUG("i < j " );
              j = i;
              i = (myRow.GetColPtr())[jj];
            }
            (*accumMatrix)[i][j] = myElement;
          }
        }         
        ATH_MSG_DEBUG("Matrix size AF "<< accumMatrix->ptrMap()->size()  );
  
      } else if ( accumMatrix) {        
        ATH_MSG_DEBUG("Matrix size b4 "<< accumMatrix->ptrMap()->size()  );
  
        for (int ii=0;ii<nDoF;ii++) {
          const TMatrixTSparseRow_const<double> myRow = (*matrix)[ii];
          int i =  myRow.GetRowIndex();
          for (int jj=0;jj<myRow.GetNindex();jj++) {
            int j = (myRow.GetColPtr())[jj];
            const double  myElement= (myRow.GetDataPtr())[jj];
            if (i<j){
              ATH_MSG_DEBUG("i < j " );
              j = i;
              i = (myRow.GetColPtr())[jj];
            }
            (*accumMatrix)[i][j] += myElement;
          }
        }
        ATH_MSG_DEBUG("Matrix size AF "<< accumMatrix->ptrMap()->size()  );
 
      } else {
        delete matrix;
        ++numberOfReadErrors;
        ATH_MSG_ERROR("Matrix allocation error!!!");
        continue;
      }
            
      delete matrix;
      
      TVectorD* hits;
      TVectorD* tracks;
      
      ATH_MSG_DEBUG("Reading hitmap ");
      hits = (TVectorD*)myFile->Get("Hits");
      if( !hits ){
        ++numberOfReadErrors;
        ATH_MSG_ERROR("Hitmap 1 not read!!!");
        continue;
      }
      
      tracks = (TVectorD*)myFile->Get("Tracks");
      if( !tracks ){
        delete hits;
        ++numberOfReadErrors;
        ATH_MSG_ERROR("Hitmap 2 not read!!!");
        continue;
      }
      
      if(hits->GetNrows() != TotalHits.GetNrows() ){
        delete hits;
        delete tracks; 
        ++numberOfReadErrors;
        ATH_MSG_ERROR("Hitmap size incorrect!!!");
        continue;
      }
      
      TotalHits   += (*hits);    
      TotalTracks += (*tracks);
      
      delete hits;
      delete tracks;
      
      myFile->Close("R");
      delete myFile;
 
      itworked =  getrusage(RUSAGE_SELF,&myusage);
      ATH_MSG_DEBUG("Memory usage [MB], total " << myusage.ru_maxrss/1024 << ", increase  " << (myusage.ru_maxrss-intialMemUse)/1024);
 
    }
 
 
 
    // create new matrix to read data from current file
    if(accumMatrix){
      if (symBigMatrix) {
        AlSymMat newMatrix(nDoF);
        //This method is ok for small matrix files
        for (int i=0;i<nDoF;i++) {
          for (int j=0;j<=i;j++) {
            newMatrix[i][j] = (*accumMatrix)[i][j];
          }
        }       
 
        *symBigMatrix += newMatrix;
        delete accumMatrix;
      } else if (spaBigMatrix) { 
        ATH_MSG_DEBUG( "should reassign matrix "<< spaBigMatrix->ptrMap()->size() );        
        *spaBigMatrix  += *accumMatrix;     
        ATH_MSG_DEBUG( "??????  "<< spaBigMatrix->ptrMap()->size() );
        delete accumMatrix;
      }
    }
 
    ATH_MSG_DEBUG( "??????  "<< m_bigmatrix->ptrMap()->size() );
 
    AlignModuleList::const_iterator imod     = moduleList->begin();
    AlignModuleList::const_iterator imod_end = moduleList->end();
    int index = 0;
    int totalhits = 0;
    for(; imod != imod_end; ++imod, ++index ) {
          AlignModule * module = *imod;
          module->setNHits((int)TotalHits(index));
          module->setNTracks((int)TotalTracks(index));
          totalhits += (int)TotalHits(index);
    }
 
 
    m_nHits = totalhits;
    m_nTracks = 0;
    m_nMeasurements = 0;
    m_scale = totalscale;
 
    return true;
  }

addFirstDerivative()

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

virtual

Implements Trk::IMatrixTool.

Definition at line 1477 of file MatrixTool.cxx.

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

addFirstDerivatives() [1/2]

void MatrixTool::addFirstDerivatives ( AlVec * vector )

virtual

adds first derivative to vector

Implements Trk::IMatrixTool.

Definition at line 1457 of file MatrixTool.cxx.

  {
     }

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 1467 of file MatrixTool.cxx.

  {
     }

addModule()

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

inlineinherited

Definition at line 116 of file IMatrixTool.h.

{ m_alignModuleMap[alignModuleIndex]=m_nentries; m_nentries += nAlignParam; }

addSecondDerivative()

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

virtual

Implements Trk::IMatrixTool.

Definition at line 1483 of file MatrixTool.cxx.

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

addSecondDerivatives() [1/2]

void MatrixTool::addSecondDerivatives ( AlSymMatBase * matrix )

virtual

adds second derivatives to matrix

Implements Trk::IMatrixTool.

Definition at line 1462 of file MatrixTool.cxx.

  {
     }

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 1472 of file MatrixTool.cxx.

  {
     }

allocateMatrix()

StatusCode MatrixTool::allocateMatrix ( int nDoF = 0 )

virtual

allocates memory for big matrix and big vector

Implements Trk::IMatrixTool.

Definition at line 100 of file MatrixTool.cxx.

  {
    ATH_MSG_INFO("allocating matrix and vector with nDoF = "<<nDoF);
 
    if (nullptr!=m_bigmatrix || nullptr!=m_bigvector)
      ATH_MSG_ERROR("big matrix already allocated!");
 
    // Decide upon the big matrix representation:
    if( m_useSparse )
      m_bigmatrix = new AlSpaMat(nDoF);
    else
      m_bigmatrix = new AlSymMat(nDoF);
 
    m_bigvector       = new AlVec(nDoF);
 
    ATH_MSG_INFO(" After Matrix and Vector allocation");
 
    // set paths for matrix and vector output
    m_bigmatrix->SetPathBin(m_pathbin.value()+m_prefixName);
    m_bigmatrix->SetPathTxt(m_pathtxt.value()+m_prefixName);
    m_bigvector->SetPathBin(m_pathbin.value()+m_prefixName);
    m_bigvector->SetPathTxt(m_pathtxt.value()+m_prefixName);
 
    ATH_MSG_INFO("set path to "<<m_pathbin.value()+m_prefixName.value());
    return StatusCode::SUCCESS;
  }

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.

  {
    return *AthCommonDataStore<PBASE>::declareProperty(hndl.name(),
                                                       hndl.value(), 
                                                       hndl.documentation());
 
  }

declareProperty()

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

inlineinherited

Definition at line 145 of file AthCommonDataStore.h.

                                                                     {
    typedef typename SG::HandleClassifier<T>::type htype;
    return AthCommonDataStore<PBASE>::declareGaudiProperty(t, htype());
  }

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.

{ return m_detStore; }

entryNumber()

int IMatrixTool::entryNumber ( int alignModuleIndex )

inlineinherited

Definition at line 118 of file IMatrixTool.h.

{ 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.

{ return m_evtStore; }

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 1620 of file MatrixTool.cxx.

  {
    return 0;
  }

finalize()

StatusCode MatrixTool::finalize

(

)

initialize

Definition at line 91 of file MatrixTool.cxx.

  {
    ATH_MSG_DEBUG("finalize() of MatrixTool");
 
    return StatusCode::SUCCESS;
  }

initialize()

StatusCode MatrixTool::initialize

(

)

initialize

Definition at line 70 of file MatrixTool.cxx.

  {
    ATH_MSG_DEBUG("initialize() of MatrixTool");
 
    // get AlignModuleTool
    if (m_alignModuleTool.retrieve().isSuccess())
      ATH_MSG_INFO("Retrieved " << m_alignModuleTool);
    else{
      msg(MSG::FATAL) << "Could not get " << m_alignModuleTool << endmsg;
      return StatusCode::FAILURE;
    }
 
    ATH_MSG_INFO("Retrieving data from the following files: ");
    for (auto & inputVectorFile : m_inputVectorFiles) {
      ATH_MSG_INFO(m_pathbin+inputVectorFile);
    }    
 
    return StatusCode::SUCCESS;
  }

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.

                                                     {
    return IID_TRKALIGNINTERFACES_IMatrixTool;
  }

msg()

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

inlineinherited

Definition at line 24 of file AthCommonMsg.h.

                                {
    return this->msgStream();
  }

msgLvl()

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

inlineinherited

Definition at line 30 of file AthCommonMsg.h.

                                               {
    return this->msgLevel(lvl);
  }

nHits()

int Trk::IMatrixTool::nHits ( ) const

inlineinherited

Definition at line 84 of file IMatrixTool.h.

{ return m_nHits; }

nMeasurements()

int Trk::IMatrixTool::nMeasurements ( ) const

inlineinherited

Definition at line 92 of file IMatrixTool.h.

{ return m_nMeasurements; }

nTracks()

int Trk::IMatrixTool::nTracks ( ) const

inlineinherited

Definition at line 88 of file IMatrixTool.h.

{ 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 1819 of file MatrixTool.cxx.

  {
    ATH_MSG_DEBUG("in postSolvinglapack()");
 
    if( z.ncol() != size) {
      msg(MSG::ERROR)<<"Eigenvector matrix has incorrect size : "<<z.ncol()<<" != "<<size<<endmsg;
      return;
    }
 
    if( (int)m_activeIndices.size() != size) {
      msg(MSG::ERROR)<<"Number of active parameters is incorrect : "<<m_activeIndices.size()<<" != "<<size<<endmsg;
      return;
    }
 
    // Compute bigvector in diagonal basis (Vb = Ut * bigvector)
    AlVec D(size);
    D = z*(*dChi2);
 
    if (m_writeEigenMat) {
 
      ATH_MSG_INFO("writing the eigenvectors in a matrix: "<< z.nrow() << "x" << z.ncol());
 
      // Set Path for the z matrix (eigenvector matrix)
      z.SetPathBin(m_pathbin.value()+m_prefixName.value());
      z.SetPathTxt(m_pathtxt.value()+m_prefixName.value());
 
      ATH_MSG_INFO("writing the eigenvector matrix: "<< m_scalaMatName);
      ATH_MSG_DEBUG("matrix will be in: "<< m_pathbin.value()+m_prefixName.value()+m_scalaMatName.value());
 
      StatusCode sc = z.Write(m_scalaMatName, true); // write the eigenvector matrix
 
      if (sc!=StatusCode::SUCCESS)
        msg(MSG::ERROR)<<"Problem writing eigenvector matrix"<<endmsg;
 
      // Set Path for the w matrix (eigenvalues matrix - diagonal bigmatrix)
      w.SetPathBin(m_pathbin.value()+m_prefixName.value());
      w.SetPathTxt(m_pathtxt.value()+m_prefixName.value());
 
      ATH_MSG_INFO("writing the eigenvectors in a vector: "<< w.size());
      ATH_MSG_INFO("writing the eigenvalues vector (diagonal bigmatrix): "<< m_scalaVecName);
      ATH_MSG_DEBUG("vector will be in: "<< m_pathbin.value()+m_prefixName.value()+m_scalaVecName.value());
 
      sc = w.WriteEigenvalueVec(m_scalaVecName, true); // write the eigenvalues vecor
 
      if (sc!=StatusCode::SUCCESS)
        msg(MSG::ERROR)<<"Problem writing eigenvector matrix"<<endmsg;
 
      if (m_writeEigenMatTxt) {
        sc = z.Write("eigenvectors.txt", false);
        if (sc!=StatusCode::SUCCESS)
          msg(MSG::ERROR)<<"Problem writing eigenvector matrix to text file"<<endmsg;
        sc = w.WriteEigenvalueVec("eigenvalues.txt", false);
        if (sc!=StatusCode::SUCCESS)
          msg(MSG::ERROR)<<"Problem writing eigenvalue vector to text file"<<endmsg;
      }
 
    }
 
    // Set eigenvalue thresholds
    const double eigenvalue_threshold = 1e-19;
 
    // weak mode removal
    if (m_modcut == -1) {
 
      ATH_MSG_INFO(" Starting the automatic Weak Mode Removal method");
 
      // create a pull vector for the alignment corrections in diagonal basis (db_pulls)
      //int nDoF=m_alignModuleTool->nAlignParameters();
      AlVec* Align_db = new AlVec(size);
      AlVec* Align_error_db = new AlVec(size);
      AlVec* AlignPull = new AlVec(size);
      ATH_MSG_DEBUG("AlignPull vector size is: "<< (*AlignPull).size());
 
      m_modcut = 0;
      bool wm_stop = false;
 
      // -----------------------------------------------------------------------
      // First pass: removing weak modes because of large pull
      // compute alignment pulls for corrections in diagonal basis (db)
      for(int i=0; i<size; i++) {
 
        (*Align_db)[i] = (-D[i]/w[i]);
        if(m_scale<=0.)[[unlikely]]{
          ATH_MSG_WARNING("postSolvingLapack: Scaling requested but scale not set. Not scaling matrix and vector.");
        } else {
          (*Align_error_db)[i] = sqrt(1.0/w[i]/m_scale);
        }
        
        if (w[i]<eigenvalue_threshold) {
          ATH_MSG_INFO("  + EigenMode " << i
                         << " removed as eigenvalue lower than the threshold " << eigenvalue_threshold
                         << ": " << w[i]);
          (*AlignPull)[i] = 0.0;
          ++m_modcut;
        }
        else
          (*AlignPull)[i] = (*Align_db)[i] / (*Align_error_db)[i];
 
        ATH_MSG_DEBUG(i << ". AlignPar: " << (*Align_db)[i] << " +- " << (*Align_error_db)[i]
                        << " (pull: " << (*AlignPull)[i]  << ") ; w[i]: " << w[i]);
      }
      ATH_MSG_INFO(" +++ Weak Mode removal ++ stop after mode "<< m_modcut << " (First pass)");
      // -----------------------------------------------------------------------
 
      // -----------------------------------------------------------------------
      // Second pass
      // if the error is greater than the correction -> cut this mode
      for(int i=m_modcut; (i<size && !wm_stop); i++) {
 
        // if the error is greater than the correction -> cut this mode
        if (fabs((*AlignPull)[i])<m_pullcut) {
          ATH_MSG_INFO("  + EigenMode " << i
                         << " removed as pull is lower than " << m_pullcut << ": "
                         << (*AlignPull)[i]);
          ++m_modcut;
        }
        else
          wm_stop = true;
 
      }
      ATH_MSG_INFO(" +++ Weak Mode removal ++ stop after mode "<< m_modcut << " (Second pass)");
      // -----------------------------------------------------------------------
 
      wm_stop = false;
 
      // ----------------------------------------------------------------------
      // Third pass
      // Check if the next eigenvalues is far away. If it is two orders of
      // magnitude bigger remove also this mode and allow the search
      for(int i=m_modcut; (i<size && !wm_stop); i++) {
 
        // if the next eigenvalues is far away -> cut this mode
        if (m_eigenvalueStep*w[i]<w[i+1]) {
          ATH_MSG_INFO("  + EigenMode " << i
                         << " removed as diff between eigenvalues, " << w[i] << " and " << w[i+1]
                         << ", is greater than " << m_eigenvalueStep);
          ++m_modcut;
        }
        else
          wm_stop = true;
 
      }
      ATH_MSG_INFO(" +++ Weak Mode removal ++ stop after mode "<< m_modcut << " (Third pass)");
      // -----------------------------------------------------------------------
 
      wm_stop = false;
 
      // -----------------------------------------------------------------------
      // Fourth pass
      // Check if the next eigenvalues is far away. If it is two orders of
      // magnitude bigger remove also this mode and allow the search
      for(int i=m_modcut; (i<size && !wm_stop); i++) {
 
        // if the next eigenvalues is far away -> cut this mode
        if ( fabs((*Align_db)[i]) > m_Align_db_step*fabs((*Align_db)[i+1]) ) {
          ATH_MSG_INFO("  + EigenMode " << i
                         << " removed as diff between corrections, " << w[i] << " and " << w[i+1]
                         << ", is greater than "
                         << m_Align_db_step);
          ++m_modcut;
        }
        else
          wm_stop = true;
 
      }
      ATH_MSG_INFO(" +++ Weak Mode removal ++ stop after mode "<< m_modcut << " (Fourth pass)");
      // -----------------------------------------------------------------------------------------------------
 
      // Free memory and clear the pointer to
      // prevent using invalid memory reference
      delete Align_db;
      delete Align_error_db;
      delete AlignPull;
      Align_db = nullptr;
      Align_error_db = nullptr;
      AlignPull = nullptr;
 
    } // end of if(m_modcut == -1)
 
    // Save some stuff to debug purposes
    /*
      if (m_storeDia) {
      std::string path = m_pathtxt+m_prefixName+"align_dia";
      std::fstream orthogon(path.c_str(), std::ios::out);
      orthogon.setf(std::ios::fixed);
      orthogon.setf(std::ios::showpoint);
      orthogon.precision(6);
 
      orthogon << std::setw(10)
      << "--------------------------------------------------------------------------------"
      << std::endl;
      orthogon << std::setw(10) << " ModeCut = " << m_modcut << std::endl;
      orthogon << std::setw(10)
      << "--------------------------------------------------------------------------------"
      << std::endl;
 
      orthogon << std::setw(10) << "mode"
      << std::setw(20) << "eigenmode"
      << std::setw(20) << "eigenmode error"
      << std::setw(25) << "eigenvalue"
      << std::endl;
 
      for( int m=0; m<size; m++) {
 
      // mode
      orthogon << std::setw(10) << m;
 
      // eigenmode (db)
      if( w[m]>1.0e-15) orthogon << std::setw(20) << -D[m]/w[m];
      else orthogon << std::setw(20) << 0.0;
 
      // error eigenmode (error_db)
      if( w[m]>1.0e-15) orthogon << std::setw(20) << sqrt(1.0/w[m]/m_scale);
      else orthogon << std::setw(20) << 0.0;
 
      // eigenvalues
      orthogon << std::setw(25) << w[m] << std::endl;
      }
      orthogon.close();
      } // end store align_dia.txt
    */
 
    AlVec delta(size);
    AlVec deltafull(size);
    AlVec errSq(size);
 
    // full covariance matrix
    CLHEP::HepSymMatrix * cov = nullptr;
    if(m_calculateFullCovariance)
      // Warning ! The matrix can be huge!
      // This can lead to memory problems
      cov = new CLHEP::HepSymMatrix(size,0);
 
    if(m_logStream)
      *m_logStream<<"/------ The Eigenvalue Spectrum -------\n";
 
    for (int i=0;i<size;i++) {
      AlVec thisdelta(size);
      for(int j=0;j<size;j++)
        thisdelta[j] = z[i][j] * (-D[i]/w[i]);
      deltafull += thisdelta;
 
      ATH_MSG_DEBUG("eigenvalue "<<w[i]);
      if( i<m_modcut ) {
        ATH_MSG_INFO("skipping eigenvalue "<<w[i]<<" , modcut is "<<m_modcut);
        if(m_logStream)
           *m_logStream<<"| skipping eigenvalue "<<w[i]<<"\n";
      }
      else if( w[i] < m_eigenvaluethreshold ) {
        ATH_MSG_INFO("skipping eigenvalue "<<w[i]<<" , cut is "<<m_eigenvaluethreshold);
        if(m_logStream)
           *m_logStream<<"| skipping eigenvalue "<<w[i]<<"\n";
      }
      else {
        if(m_logStream)
          *m_logStream<<"| "<<w[i]<<"\n";
 
        delta += thisdelta;
        for(int j=0;j<size;j++) {
          errSq[j] += z[i][j] * z[i][j] / w[i];
          if(m_calculateFullCovariance) {
            for(int k=0;k<=j;k++)
              (*cov)[j][k] += z[i][j] * z[i][k] / w[i];
          }
        }
      }
    }
 
    if(m_logStream)
      *m_logStream<<"\\----- End of Eigenvalue Spectrum -----\n";
 
    ATH_MSG_DEBUG("Alignment constants:");
 
    DataVector<AlignPar> * alignParList = m_alignModuleTool->alignParList1D();
 
    // compute alignment corrections (translations in mm and rotations in rad) and their variances
    for(int i=0; i<size; i++) {
 
      double param = delta[i];
      double err = sqrt(2.*std::fabs(errSq[i]));
 
      int idof = m_activeIndices[i];
      AlignPar * alignPar=(*alignParList)[idof];
 
      // undo the sigma scaling
      double sigma = alignPar->sigma();
 
      param *= sigma;
      err *= sigma;
 
      // undo normalization scaling of error
      if(m_scaleMatrix && m_scale>0.)
        err /= sqrt(m_scale);
 
      ATH_MSG_DEBUG(i <<" : "<< param << " +/- "<< err);
      ATH_MSG_DEBUG("cov("<<i<<")="<<errSq[i]<<", sigma: "<<sigma);
      ATH_MSG_DEBUG("init par: "<<alignPar->initPar());
      alignPar->setPar(param, err);
      ATH_MSG_DEBUG("set param to "<<param<<" for alignPar "<<alignPar);
      ATH_MSG_DEBUG(*(*alignParList)[idof]);
    }
 
    if(m_logStream) {
      printGlobalSolution(*m_logStream, cov);
 
      // norm of first derivative
      double norm1st = dChi2->norm();
      if(m_scaleMatrix && m_scale>0.) // undo normalization scaling
        norm1st *= m_scale;
      *m_logStream<<"norm of first derivative :            "<<norm1st<<"\n";
 
      if(d2Chi2) {
        // distance to solution
        double dist = ( (*d2Chi2) * deltafull + (*dChi2) ).norm();
        if(m_scaleMatrix && m_scale>0.) // undo normalization scaling
          dist *= m_scale;
        *m_logStream<<"distance to solution :                "<<dist<<"\n";
 
        // calculate chi2 of the alignment change
        double chi2 = delta * (*d2Chi2) * delta * .5;
        if(m_scaleMatrix && m_scale>0.) // undo normalization scaling
          chi2 *= m_scale;
        *m_logStream<<"delta(chi2) of the alignment change : "<<chi2<<" / "<<size<<"\n";
      }
    }
 
    delete cov;
  }

prepareBinaryFiles()

void MatrixTool::prepareBinaryFiles ( int solveOption )

virtual

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

Implements Trk::IMatrixTool.

Definition at line 128 of file MatrixTool.cxx.

  {
     }

printGlobalSolution() [1/2]

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

Definition at line 1489 of file MatrixTool.cxx.

  {
    const AlignModuleList * alignModules = m_alignModuleTool->alignModules1D();
 
    AlignModuleList::const_iterator imod     = alignModules->begin();
    AlignModuleList::const_iterator imod_end = alignModules->end();
    for( ; imod!=imod_end; ++imod) {
      AlignModule * module = *imod;
 
      DataVector<AlignPar> * alignPars = m_alignModuleTool->getAlignPars(module);
      int thisNDoF = alignPars->size();
 
      // fill local covariance matrix
      CLHEP::HepSymMatrix * covsub = nullptr;
      if(cov && module->nHits() >= m_minNumHits && module->nTracks() >= m_minNumTrks) {
        covsub = new CLHEP::HepSymMatrix(thisNDoF,0);
        for (int i=0;i<thisNDoF;++i) {
          int ipar = alignPars->at(i)->index();
          double sigma_i = alignPars->at(i)->sigma();
 
          std::vector<int>::iterator itActive = std::find(m_activeIndices.begin(),m_activeIndices.end(),ipar);
          if( itActive == m_activeIndices.end() )
            continue;
          int iActive = std::distance(m_activeIndices.begin(),itActive);
 
          for (int j=0;j<=i;++j) {
            int jpar = alignPars->at(j)->index();
            double sigma_j = alignPars->at(j)->sigma();
 
            std::vector<int>::iterator jtActive = std::find(m_activeIndices.begin(),m_activeIndices.end(),jpar);
            if( jtActive == m_activeIndices.end() )
              continue;
            int jActive = std::distance(m_activeIndices.begin(),jtActive);
 
            (*covsub)[i][j] = (*cov)[iActive][jActive] * sigma_i * sigma_j;
          }
        }
      }
 
      printModuleSolution(os,module,covsub);
 
      delete covsub;
    }
    os << "--------------------------------------------------------------------------------" << std::endl;
  }

printGlobalSolution() [2/2]

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

Definition at line 1536 of file MatrixTool.cxx.

  {
    CLHEP::HepSymMatrix * cov = nullptr;
    if(cov0) {
      int nsize = cov0->GetNrows();
      cov = new CLHEP::HepSymMatrix(nsize,0);
 
      for(int i=0; i<nsize; i++)
        for(int j=0; j<=i; j++)
          (*cov)[i][j] = (*cov0)[i][j];
    }
 
    printGlobalSolution(os,cov);
 
    delete cov;
  }

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 1574 of file MatrixTool.cxx.

  {
    os << "--------------------------------------------------------------------------------" << std::endl;
    os << "Alignment parameters for module: " << module->name() << std::endl;
    os << "Number of tracks passing: " << module->nTracks() << std::endl;
    if(m_minNumHits>0 && module->nHits()<m_minNumHits) {
      os << "Number of hits too small: "<<module->nHits()<<" < "<<m_minNumHits<<"  Skipping the module\n";
      return;
    }
    if(m_minNumTrks>0 && module->nTracks()<m_minNumTrks) {
      os << "Number of tracks too small: "<<module->nTracks()<<" < "<<m_minNumTrks<<"  Skipping the module\n";
      return;
    }
    os << "Number of hits seen:      " << module->nHits() << std::endl;
    os << "Number of tracks seen:    " << module->nTracks() << std::endl;
 
    DataVector<AlignPar> * alignPars = m_alignModuleTool->getAlignPars(module);
    int thisNDoF = alignPars->size();
 
    if(alignPars->empty())
      os << "No active parameters" << std::endl;
    else
    {
      // output alignment parameters and errors
      DataVector<AlignPar>::const_iterator ipar     = alignPars->begin();
      DataVector<AlignPar>::const_iterator ipar_end = alignPars->end();
      for ( ; ipar != ipar_end; ++ipar) {
        const AlignPar * par = *ipar;
        os << std::format("{:<10}{:<12.5g} +/- {:<12.5g}\n",
                          par->dumpType(),
                          par->par(),
                          par->err());
      }
 
      if(cov) {
        // calculate local correlation matrix
        CLHEP::HepSymMatrix corrsub(thisNDoF,0);
        for(int irow=0; irow<thisNDoF; ++irow)
          for(int icol=0; icol<=irow; ++icol)
            corrsub[irow][icol] = (*cov)[irow][icol] / sqrt((*cov)[irow][irow] * (*cov)[icol][icol]);
        os << "Local correlation matrix: " << corrsub << std::flush;
      }
    }
  }

readHitmaps()

void MatrixTool::readHitmaps ( )

private

Definition at line 2320 of file MatrixTool.cxx.

  {
    ATH_MSG_INFO("read hitmaps from files");
 
    const AlignModuleList * moduleList = m_alignModuleTool->alignModules1D();
    int nModules = moduleList->size();
 
    AlMat hitmap(nModules,2);
    int nFiles = (int)m_inputHitmapFiles.size();
    for(int imap=0;imap<nFiles; imap++) {
      AlMat nextmap(nModules,2);
 
      ATH_MSG_INFO("Reading hitmap "<<imap<<" from file "<<m_inputHitmapFiles[imap]);
 
      if(nextmap.ReadScalaPack(m_inputHitmapFiles[imap]).isFailure()) {
        ATH_MSG_WARNING("Problem reading hitmap from \'"<<m_inputHitmapFiles[imap]<<"\'. Skipping.");
        continue;
      }
 
      if(nextmap.nrow()!=nModules || nextmap.ncol()!=2) {
        ATH_MSG_WARNING("Matrix in file \'"<<m_inputHitmapFiles[imap]<<"\' has wrong size ("
            <<nextmap.nrow()<<" x "<<nextmap.ncol()<<"), should be ("<<nModules<<" x 2). Skipping.");
        continue;
      }
 
      hitmap += nextmap;
    }
 
    AlignModuleList::const_iterator imod     = moduleList->begin();
    AlignModuleList::const_iterator imod_end = moduleList->end();
    int index = 0;
    int totalhits = 0;
    for(; imod != imod_end; ++imod) {
      AlignModule * module = *imod;
      //if ((int)hitmap[index][0]!=(int)module->identify().get_identifier32().get_compact()) ATH_MSG_ERROR("bad module identifier");
      //module->setIdentifier((Identifier)hitmap[index][0]);
      module->setNHits((int)hitmap[index][0]);
      module->setNTracks((int)hitmap[index][1]);
      totalhits += (int)hitmap[index][0];
      index++;
    }
 
    m_nHits = totalhits;
    m_nTracks = 0;
    m_nMeasurements = 0;
 
    ATH_MSG_INFO("Hitmap accumulated from "<<nFiles<<" files with total of "<<totalhits<<" hits.");
  }

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.

  {
    h.renounce();
    PBASE::renounce (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.

                                                          {
    handlesArray.renounce();
  }

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.

{ 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.

{ m_nHits = n; }

setNMeasurements()

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

inlineinherited

set number of measurements

Definition at line 91 of file IMatrixTool.h.

{ m_nMeasurements = n; }

setNTracks()

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

inlineinherited

set number of tracks

Definition at line 87 of file IMatrixTool.h.

{ m_nTracks = n; }

solve()

int MatrixTool::solve ( )

virtual

solves for alignment parameters

Implements Trk::IMatrixTool.

Definition at line 1136 of file MatrixTool.cxx.

  {
    // ============
    // solve
    // ============
    ATH_MSG_DEBUG("in MatrixTool::solve()");    
 
    // set normalization scale to number of hits for now
    if(m_scale<0)
      m_scale = m_nHits;
 
    //-------------------------------------------------------
    // write matrix and vector to file
    if (m_writeMat) {
      // version has to be 2 to for reading matrices and vectors back in to work properly
      double dummyVersion(2.);
 
      // make map of matrix entry to module index (set by geometry manager tool)
      std::map<int,unsigned long long> modIndexMap;
      std::map<int,std::string>        modNameMap;
      DataVector<AlignPar>* alignPars = m_alignModuleTool->alignParList1D();
      for (int i=0;i<(int)alignPars->size();i++) {
        modIndexMap[i]=(*alignPars)[i]->alignModule()->identify().get_compact();
        modNameMap [i]=(*alignPars)[i]->alignModule()->name();
      }
 
      // binary files
      ATH_MSG_DEBUG("writing binary files");
      StatusCode sc1 = m_bigmatrix->Write("matrix.bin",true,m_wSqMatrix,m_scale,dummyVersion);
      StatusCode sc2 = m_bigvector->WritePartial("vector.bin",true,m_scale,modIndexMap,dummyVersion);
      if (!sc1.isSuccess() || !sc2.isSuccess()) {
        msg(MSG::ERROR)<<"problem writing matrix or vector"<<endmsg;
        return -1;
      }
 
      if (m_writeMatTxt) {
 
        // text files
        ATH_MSG_DEBUG("writing text files");
        sc1 = m_bigmatrix->Write("matrix.txt",false,m_wSqMatrix,m_scale,dummyVersion);
        sc2 = m_writeModuleNames ? 
          m_bigvector->WritePartial("vector.txt",false,m_scale,modNameMap,dummyVersion) :
          m_bigvector->WritePartial("vector.txt",false,m_scale,modIndexMap,dummyVersion);
  
        if (!sc1.isSuccess() || !sc2.isSuccess()) {
          msg(MSG::ERROR)<<"problem writing matrix or vector"<<endmsg;
          return -1;
        }
      }
 
      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)");
    }
 
    //-------------------------------------------------------
    // write hitmap to file
    if (m_writeHitmap)
      writeHitmap();
 
    if(m_writeTFile)
      storeInTFile(m_pathbin.value()+m_prefixName.value()+m_tfileName.value());
 
 
    if(!m_runLocal && m_solveOption==0) {
      ATH_MSG_DEBUG("No solving requested.");
      return 1;
    }
 
    //-------------------------------------------------------
    // rescale the vector and the matrix according to sigmas
    // and apply soft mode cut
 
    ATH_MSG_DEBUG("rescaling the matrix/vector and applying the soft-mode-cut");
 
    DataVector<AlignPar>* alignParList = m_alignModuleTool->alignParList1D();
    int nDoF = alignParList->size();
 
 
    const AlSymMat * chkMatrix = dynamic_cast<const AlSymMat*>(m_bigmatrix);
    if(chkMatrix){
      // Method when using the dense matrix
      for (int i=0;i<nDoF;i++) {
        // scale the vector
        double sigma_i = (*alignParList)[i]->sigma();
        double softCut = 2 * pow( (*alignParList)[i]->softCut() , -2 );
        (*m_bigvector)[i] *= sigma_i;
  
        for (int j=0;j<=i;j++) {
          // scale the matrix
          if ((*chkMatrix)[i][j] != 0.) {
            double sigma_j = (*alignParList)[j]->sigma();
            (*m_bigmatrix)[i][j] *= sigma_i * sigma_j;
          }
          // apply soft-mode-cut
          if (i==j && m_softEigenmodeCut>0.){
            (*m_bigmatrix)[i][j] += m_softEigenmodeCut * softCut;
         
          }
  
          // set first and second derivatives on AlignPar
          if (i==j) {
            (*alignParList)[i]->setFirstDeriv((*m_bigvector)[i]/sigma_i);
            (*alignParList)[i]->setSecndDeriv((*m_bigmatrix)[i][i]/sigma_i/sigma_i);
          }
        }
      }    
    } else {
      // Method when using the sparse matrix
      for (const datamap::value_type& p : *m_bigmatrix->ptrMap()) {
        int i = p.first.first;
        int j = p.first.second;
  
        // Scale matrix
        double sigma_i = (*alignParList)[i]->sigma();
        double sigma_j = (*alignParList)[j]->sigma();
        
        (*m_bigmatrix)[i][j] *= sigma_i * sigma_j;
  
      }
  
  
      for (int i=0;i<nDoF;i++) {
        // scale the vector
        double sigma_i = (*alignParList)[i]->sigma();
        (*m_bigvector)[i] *= sigma_i;
        if (m_softEigenmodeCut >0. ){
          double softCut = 2 * pow( (*alignParList)[i]->softCut() , -2 );
          (*m_bigmatrix)[i][i] += m_softEigenmodeCut * softCut;
          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]);
        }
        (*alignParList)[i]->setFirstDeriv((*m_bigvector)[i]/sigma_i);
        (*alignParList)[i]->setSecndDeriv((*m_bigmatrix)[i][i]/sigma_i/sigma_i);
      }
    }
    
    unsigned long long OldPixelIdentifier = 37769216; //Identifier for the Pixel Detector
    unsigned long long IBLIdentifier      = 33574912; //Identifier for the Pixel Detector
    
    unsigned long long SCT_ECA_8_Identifier = 218116096; //Identifier for the SCT ECA Last Disk
    std::string SCT_ECA_8_Name = "SCT/EndcapA/Disk_8";
 
    
 
    ATH_MSG_INFO("rescaling done");
    ATH_MSG_INFO("Javi: Printing (*alignParList)[i]->alignModule()->identify32()");
 
    // select modules with non-zero tracks
    for(int i=0;i<nDoF;i++)
      {
        ATH_MSG_DEBUG(i);
        ATH_MSG_DEBUG((*alignParList)[i]->alignModule()->identify32());
        ATH_MSG_DEBUG((*alignParList)[i]->alignModule());
        ATH_MSG_DEBUG((*alignParList)[i]->alignModule()->name());
        ATH_MSG_DEBUG((*alignParList)[i]->paramType());
 
        //Skip solving for Pixel or IBL:
    const auto & theParameterList = *alignParList;      
    const auto & thisIdentifier = theParameterList[i]->alignModule()->identify32();
    const auto & thisName       = theParameterList[i]->alignModule()->name();
    const auto & thisParameterType = theParameterList[i]->paramType();
    const bool oldPixel = (thisIdentifier == OldPixelIdentifier);
    const bool ibl = (thisIdentifier == IBLIdentifier);
    const bool SCTECA8   = (thisIdentifier == SCT_ECA_8_Identifier);
    const bool SCTECA8_n = (thisName.find(SCT_ECA_8_Name)!= std::string::npos);
 
    if (m_DeactivateSCT_ECA_LastDisk){ 
      if (SCTECA8)
        {ATH_MSG_INFO( "SCT ECA Last Disk DoF have been skipped in the solving because DeactivateSCT_ECA_LastDisk is set to True");
          continue;}
      if (SCTECA8_n)
        {ATH_MSG_INFO( "SCT ECA Last Disk DoF have been skipped in the solving because DeactivateSCT_ECA_LastDisk is set to True");
          continue;}
    }
 
    if (m_AlignIBLbutNotPixel) //If m_AlignIBLbutNotPixel is set to True, Pixel will be skipped in the solving.
      if  (oldPixel)  
        {ATH_MSG_INFO( "Pixel DoF have been skipped in the solving because AlignIBLbutNotPixel is set to True");
          continue;} 
    
    if (m_AlignPixelbutNotIBL) //If m_AlignPixelbutNotIBL is set to True, IBL will be skipped in the solving.
      if  (ibl)  
        {ATH_MSG_INFO( "IBL DoF have been skipped in the solving because AlignPixelbutNotIBL is set to True");
          continue;}
 
        //For specific DoF: (*alignParList)[i]->paramType() = 0,1,2,3,4,5,6 for Tx,Ty,Tz,Rx,Ry,Rz,Bx
    //Pixel Dofs:
    if (m_Remove_Pixel_Tx) //If m_Remove_Pixel_Tx is set to True, Pixel Tx will be skipped in the solving.
      if  ( oldPixel and (thisParameterType == 0))  
        {ATH_MSG_INFO( "Pixel Tx DoF has been skipped in the solving because Remove_Pixel_Tx is set to True");
          continue;} 
        
    if (m_Remove_Pixel_Ty) //If m_Remove_Pixel_Ty is set to True, Pixel Ty will be skipped in the solving.
      if  ( oldPixel and (thisParameterType == 1))  
        {ATH_MSG_INFO( "Pixel Ty DoF has been skipped in the solving because Remove_Pixel_Ty is set to True");
          continue;} 
        
    if (m_Remove_Pixel_Tz) //If m_Remove_Pixel_Tz is set to True, Pixel Tz will be skipped in the solving.
      if  (oldPixel and (thisParameterType == 2))  
        {ATH_MSG_INFO( "Pixel Tz DoF has been skipped in the solving because Remove_Pixel_Tz is set to True");
          continue;} 
        
    if (m_Remove_Pixel_Rx) //If m_Remove_Pixel_Rx is set to True, Pixel Rx will be skipped in the solving.
      if  (oldPixel and (thisParameterType == 3))
        {ATH_MSG_INFO( "Pixel Rx DoF has been skipped in the solving because Remove_Pixel_Rx is set to True");
          continue;} 
        
    if (m_Remove_Pixel_Ry) //If m_Remove_Pixel_Ry is set to True, Pixel Ry will be skipped in the solving.
      if  (oldPixel and (thisParameterType == 4)) 
        {ATH_MSG_INFO( "Pixel Ry DoF has been skipped in the solving because Remove_Pixel_Ry is set to True");
          continue;} 
        
    if (m_Remove_Pixel_Rz) //If m_Remove_Pixel_Rz is set to True, Pixel Rz will be skipped in the solving.
      if  (oldPixel and (thisParameterType == 5)) 
        {ATH_MSG_INFO( "Pixel Rz DoF has been skipped in the solving because Remove_Pixel_Rz is set to True");
          continue;} 
        
        //IBL Dofs:
    if (m_Remove_IBL_Tx) //If m_Remove_IBL_Tx is set to True, IBL Tx will be skipped in the solving.
      if  (ibl and (thisParameterType == 0))  
        {ATH_MSG_INFO( "IBL Tx DoF has been skipped in the solving because Remove_IBL_Tx is set to True");
          continue;} 
        
    if (m_Remove_IBL_Ty) //If m_Remove_IBL_Ty is set to True, IBL Ty will be skipped in the solving.
      if  (ibl and (thisParameterType == 1)) 
        {ATH_MSG_INFO( "IBL Ty DoF has been skipped in the solving because Remove_IBL_Ty is set to True");
          continue;} 
        
    if (m_Remove_IBL_Tz) //If m_Remove_IBL_Tz is set to True, IBL Tz will be skipped in the solving.
      if  (ibl and (thisParameterType == 2))  
        {ATH_MSG_INFO( "IBL Tz DoF has been skipped in the solving because Remove_IBL_Tz is set to True");
          continue;} 
        
    if (m_Remove_IBL_Rx) //If m_Remove_IBL_Rx is set to True, IBL Rx will be skipped in the solving.
      if  (ibl and (thisParameterType == 3))
        {ATH_MSG_INFO( "IBL Rx DoF has been skipped in the solving because Remove_IBL_Rx is set to True");
          continue;} 
        
    if (m_Remove_IBL_Ry) //If m_Remove_IBL_Ry is set to True, IBL Ry will be skipped in the solving.
      if  (ibl and (thisParameterType == 4))  
        {ATH_MSG_INFO( "IBL Ry DoF has been skipped in the solving because Remove_IBL_Ry is set to True");
          continue;} 
 
    if (m_Remove_IBL_Rz) //If m_Remove_IBL_Rz is set to True, IBL Rz will be skipped in the solving.
      if  (ibl and (thisParameterType == 5))  
        {ATH_MSG_INFO( "IBL Rz DoF has been skipped in the solving because Remove_IBL_Rz is set to True");
          continue;}
 
    if(theParameterList[i]->alignModule()->nHits() >= m_minNumHits && theParameterList[i]->alignModule()->nTracks() >= m_minNumTrks)
      m_activeIndices.push_back(i);
      }
    m_aNDoF = m_activeIndices.size();
    ATH_MSG_DEBUG("aNDoF/nDoF: "<<m_aNDoF<<"/"<<nDoF);
 
    // --------------------
    // now do the SOLVING
    // --------------------
 
    int info = 0;
 
    // first Local solving
    if (m_runLocal)
      info = solveLocal();
 
    // remove spurious modules and resize
    if (m_removeSpurious) {
 
      ATH_MSG_INFO("Spurious removal not implemented at the moment.");
/*      if (StatusCode::SUCCESS != spuriousRemoval()) {
        ATH_MSG_ERROR("Problem while trying to remove spurious. Stopping solving");
        return -1;
      }
 
      // if nDoF=0, bad job...
      int NDoF = m_alignModuleTool->nAlignParameters();
      if(NDoF==0) {
        ATH_MSG_WARNING("Removal removed everything: NDoF=" << NDoF << " !!!");
        return 1;
      }
      ATH_MSG_DEBUG("NDoF: " << NDoF);
*/
    }
 
    // --------------------
    // now Global solving
    switch(m_solveOption) {
 
      case NONE:
        ATH_MSG_DEBUG("No global solving requested.");
        break;
 
      case SOLVE_ROOT:
        info = solveROOT();
        break;
 
      case SOLVE_CLHEP:
        info = solveCLHEP();
        break;
 
      case SOLVE:
      case DIRECT_SOLVE:
      case DIRECT_SOLVE_CLUSTER:
        info = solveLapack();
        break;
 
      case SOLVE_FAST:
      case DIRECT_SOLVE_FAST:
        info = solveSparseEigen();
        break;
 
      default:
        ATH_MSG_INFO("Unknown solving option.");
        info = 0;
        break;
    }
 
    ATH_MSG_INFO("Return value from solving: "<<info);
 
    return info;
  }

solveCLHEP()

int MatrixTool::solveCLHEP ( )

private

Definition at line 259 of file MatrixTool.cxx.

  {
    ATH_MSG_INFO("solving Global using CLHEP");
    if(m_logStream) {
      *m_logStream<<"*************************************************************\n";
      *m_logStream<<"**************  solving using Global method  ****************\n";
      *m_logStream<<"**************          using CLHEP          ****************\n";
      *m_logStream<<"*************************************************************\n";
    }
 
    // start measuring time
    clock_t starttime = clock();
 
    DataVector<AlignPar> * alignParList = m_alignModuleTool->alignParList1D();
    //const AlignModuleList* alignModules = m_alignModuleTool->alignModules1D();
    for (int i=0;i<(int)alignParList->size();i++)
      ATH_MSG_DEBUG("ap["<<i<<"]="<<(*alignParList)[i]);
 
    int nDoF = m_alignModuleTool->nAlignParameters();
 
    // some debugging output
    if (msgLvl(MSG::DEBUG)) {
      msg(MSG::DEBUG)<<"dumping matrix and vector to screen"<<endmsg;
      for (int i=0;i<nDoF;i++)
        for (int j=0;j<nDoF;j++)
        //if (std::fabs((*m_bigmatrix)[i][j])>.0001)
          msg(MSG::DEBUG)<<i<<", "<<j<<" : "<<(*m_bigmatrix)[i][j] <<endmsg;
 
      for (int i=0;i<nDoF;i++)
        msg(MSG::DEBUG)<<i <<" : "<<(*m_bigvector)[i]<<endmsg;
    }
 
    // get rescaled first and second derivatives
    CLHEP::HepSymMatrix * d2Chi2 = new CLHEP::HepSymMatrix(m_aNDoF,0);
    CLHEP::HepVector    * dChi2  = new CLHEP::HepVector(m_aNDoF,0);
    for (int iActive=0;iActive<m_aNDoF;iActive++) {
      int i = m_activeIndices[iActive];
      (*dChi2)[iActive] = (*m_bigvector)[i];
      for (int jActive=0;jActive<m_aNDoF;jActive++) {
        int j = m_activeIndices[jActive];
        (*d2Chi2)[iActive][jActive] = (*m_bigmatrix)[i][j];
      }
    }
 
    ATH_MSG_DEBUG("First derivatives:" << (*dChi2));
    ATH_MSG_DEBUG("Second derivatives:" << (*d2Chi2));
 
    CLHEP::HepSymMatrix cov(m_aNDoF,0);
    CLHEP::HepVector delta(m_aNDoF,0);
    CLHEP::HepVector deltafull(m_aNDoF,0);
 
    bool status=true;
    int ierr(0);
    if(!m_diagonalize) {
      // ==========================================================
      //  Run Matrix Inversion
      ATH_MSG_INFO("Running matrix inversion");
      if(m_logStream)
        *m_logStream<<"Running matrix inversion\n";
 
      cov = *d2Chi2;
      cov.invert(ierr);
      if(ierr>0)
        msg(MSG::ERROR)<<"CLHEP inversion status flag = "<<ierr<<endmsg;
      else
        ATH_MSG_INFO("CLHEP inversion OK");
      if(m_logStream)
        *m_logStream<<"CLHEP inversion status flag = "<<ierr<<"\n";
 
      // calculate corrections
      delta = cov * (*dChi2);
 
      // the covariance matrix is actually defined as 2 * d2Chi2^-1
      ATH_MSG_DEBUG("Result: "<<delta);
      ATH_MSG_DEBUG("cov: "<<cov*2);
 
      // -----------------------
      // calculate also for matrix and vector multiplied by factor 0.5
      // this should make no difference if everything is correct but
      // it can go wrong if insensitive DoF is included
      CLHEP::HepSymMatrix cov2 = *d2Chi2 * .5;
      // invert the matrix
      int ierr2 = 0;
      cov2.invert(ierr2);
      if(ierr2>0)
        msg(MSG::WARNING)<<"Second CLHEP inversion status flag = "<<ierr2<<endmsg;
 
      CLHEP::HepVector delta2 = cov2 * (*dChi2) * .5;
      for (int i=0;i<delta.num_row(); ++i)
        if ( fabs((delta[i] - delta2[i])/delta[i]) > 1e-5 ) {
          msg(MSG::WARNING)<<"Something's wrong with the matrix inversion: delta["<<i<<"] = "<<delta[i]<<"    delta2["<<i<<"] = "<<delta2[i]<<endmsg;
          status=false;
          break;
        }
 
      if(m_logStream && (ierr2>0 || !status)) {
        *m_logStream<<"CLHEP inversion status flag for halfed matrix = "<<ierr2<<"\n";
        *m_logStream<<"Matrix inversion check failed\n";
        *m_logStream<<"\n";
      }
      // -- end of check of matrix inversion
    }
    else {
      // ==========================================================
      // Run Diagonalization
      ATH_MSG_INFO("Running diagonalization");
      if(m_logStream)
        *m_logStream<<"Running diagonalization\n";
 
      CLHEP::HepSymMatrix D = *d2Chi2;
      CLHEP::HepMatrix U = CLHEP::diagonalize( &D );
 
      ATH_MSG_INFO("Diagonalization done");
      //sold = U*sdiag*U.T.
 
      // reorder eigenvalues ascending
      // eigenvectors need to be reordered consistently
      ATH_MSG_DEBUG(" Reordering eigenvalues ascending ");
      for (int i=0; i<m_aNDoF-1; i++)
        for (int j=i+1; j<m_aNDoF; j++)
          if(D[j][j] < D[i][i]) {
            // swap eigenvalues
            double ei = D[i][i];
            D[i][i] = D[j][j];
            D[j][j] = ei;
            // swap eigenvectors
            for(int k=0;k<m_aNDoF; k++) {
              double ev = U[k][i];
              U[k][i] = U[k][j];
              U[k][j] = ev;
            }
          }
 
      // how do I now get the eigenvalues? this cannot be the most
      // efficient way ...  CLHEP::HepSymMatrix D = d2Chi2->similarityT( U );
      CLHEP::HepVector eigenvector(m_aNDoF);
 
      if(m_logStream)
        *m_logStream<<"/------ The Eigenvalue Spectrum -------\n";
 
      ATH_MSG_DEBUG("Calculating eigenvalues");
      for(int imode=0; imode<m_aNDoF; ++imode) {
 
        // get the relevant eigenvector
        for(int irow=0; irow<m_aNDoF; ++irow)
          eigenvector[irow] = U[irow][imode];
 
        // calculate the eigenvalue
        //double eigenvalue = d2Chi2->similarity( eigenvector );
        double eigenvalue = D[imode][imode];
        ATH_MSG_DEBUG("eigenvalue "<<eigenvalue);
 
        double evdotb = dot(*dChi2,eigenvector);
        CLHEP::HepVector thisdelta = evdotb/eigenvalue * eigenvector;
        deltafull += thisdelta;
 
        if(imode<m_modcut) {
          ATH_MSG_INFO("skipping eigenvalue "<<imode<<" : "<<eigenvalue<<" , modcut is "<<m_modcut);
          if(m_logStream)
             *m_logStream<<"| skipping eigenvalue "<<eigenvalue<<"\n";
        }
        else if( eigenvalue < m_eigenvaluethreshold ) {
          ATH_MSG_INFO("skipping eigenvalue "<<eigenvalue<<" , cut is "<<m_eigenvaluethreshold);
          if(m_logStream)
             *m_logStream<<"| skipping eigenvalue "<<eigenvalue<<"\n";
        }
        else {
          if(m_logStream)
            *m_logStream<<"| "<<eigenvalue<<"\n";
 
          delta += thisdelta;
 
          // this is the time consuming part
          for(int irow=0; irow<m_aNDoF; ++irow)
            for(int icol=0; icol<=irow; ++icol)
              cov[irow][icol] += eigenvector[irow] * eigenvector[icol] / eigenvalue;
        }
      }
 
      // the covariance matrix is actually defined as 2 * d2Chi2^-1
 
      ATH_MSG_DEBUG("Result: "<<delta);
      ATH_MSG_DEBUG("cov: "<<cov);
 
      if(m_logStream)
        *m_logStream<<"\\----- End of Eigenvalue Spectrum -----\n";
 
      // end of diagonalization
      // ==========================================================
    }
 
    // stop measuring time
    clock_t stoptime = clock();
    double totaltime = (stoptime-starttime)/double(CLOCKS_PER_SEC);
    ATH_MSG_INFO("Time spent in solveCLHEP: "<<totaltime<<" s");
 
    if(ierr==0 && status)
    {
      ATH_MSG_DEBUG("Alignment constants:");
      for (int iAdof=0;iAdof<m_aNDoF;iAdof++) {
 
        int idof = m_activeIndices[iAdof];
        AlignPar * alignPar=(*alignParList)[idof];
 
        double sigma = alignPar->sigma();
        double param = -delta[iAdof] * sigma;
        double err = std::sqrt(2.*std::fabs(cov[iAdof][iAdof])) * sigma;
 
        ATH_MSG_DEBUG(iAdof <<" : "<< param << " +/- "<< err);
        ATH_MSG_DEBUG("cov("<<iAdof<<")="<<cov[iAdof][iAdof]<<", sigma: "<<sigma);
        ATH_MSG_DEBUG("init par: "<<alignPar->initPar());
        alignPar->setPar(param,err);
        ATH_MSG_DEBUG("set param to "<<param<<" for alignPar "<<alignPar);
        ATH_MSG_DEBUG(*(*alignParList)[idof]);
      }
 
      if(m_logStream) {
        printGlobalSolution(*m_logStream,&cov);
 
        // norm of first derivative
        *m_logStream<<"norm of first derivative :            "<<dChi2->norm()<<"\n";
 
        // distance to solution
        double dist = ( - (*d2Chi2) * deltafull + (*dChi2) ).norm();
        *m_logStream<<"distance to solution :                "<<dist<<"\n";
 
        // calculate chi2 of the alignment change
        double chi2 = d2Chi2->similarity(delta) * .5;
        *m_logStream<<"delta(chi2) of the alignment change : "<<chi2<<" / "<<m_aNDoF<<"\n";
 
        // time spent here
        *m_logStream<<"time spent in solve :                 "<<totaltime<<" s\n";
      }
    }
 
    delete d2Chi2;
    delete dChi2;
 
    return 1;
  }

solveLapack()

int MatrixTool::solveLapack ( )

private

Definition at line 1626 of file MatrixTool.cxx.

  {
    ATH_MSG_INFO("solving Global using Lapack");
    if(m_logStream) {
      *m_logStream<<"*************************************************************\n";
      *m_logStream<<"**************  solving using Global method  ****************\n";
      *m_logStream<<"**************         using LAPACK          ****************\n";
      *m_logStream<<"*************************************************************\n";
    }
 
    // get rescaled first and second derivatives
    AlSymMat* aBetterMat = new AlSymMat(m_aNDoF);
    AlVec*    aBetterVec = new AlVec(m_aNDoF);
    for (int iActive=0;iActive<m_aNDoF;iActive++) {
      int i = m_activeIndices[iActive];
      (*aBetterVec)[iActive] = (*m_bigvector)[i];
      for (int jActive=0;jActive<m_aNDoF;jActive++) {
        int j = m_activeIndices[jActive];
        (*aBetterMat)[iActive][jActive] = (*m_bigmatrix)[i][j];
      }
    }
 
    // normalize bigmatrix and bigvector
    if(m_scaleMatrix) {
      if(m_scale<=0.)
        ATH_MSG_WARNING("Scaling requested but scale not set. Not scaling matrix and vector.");
      else {
        (*aBetterVec) *= 1./m_scale;
        (*aBetterMat) *= 1./m_scale;
      }
    }
 
    ATH_MSG_DEBUG("Now Solving alignment using lapack diagonalization routine dspev...");
 
    if (m_calDet) {
      const double tol = 1.e-20;
      // compute final determinant
      double determ = (*aBetterMat).determinant();
      ATH_MSG_INFO("Determinant: " << determ);
      if (fabs(determ) < tol)
        ATH_MSG_WARNING("Matrix is singular!");
    }
 
    // store the original matrix for checks
    AlSymMat * d2Chi2 = nullptr;
    if (m_calculateFullCovariance)
      d2Chi2 = new AlSymMat(*aBetterMat);
 
    clock_t starttime = clock();
 
    // declare transition matrix + vector to store eigenvalues
    AlMat z(m_aNDoF,m_aNDoF);
    AlVec w(m_aNDoF); // vector to store the eigenvalues
    ATH_MSG_DEBUG("MatrixTool::after z/w allocation");
 
    char jobz = 'V';
    int info = (*aBetterMat).diagonalize(jobz,w,z);
    ATH_MSG_DEBUG(" info: " << info);
    ATH_MSG_INFO("MatrixTool::after diagonalization");
 
    // stop time calculation
    clock_t stoptime = clock();
    double time_diag = (stoptime-starttime)/double(CLOCKS_PER_SEC);
    ATH_MSG_INFO(" - time spent diagonalizing the matrix: "<<time_diag<<" s");
 
    double time_solve = 0.;
    if (info==0) {
      starttime = clock();
      postSolvingLapack(aBetterVec,d2Chi2,w,z,m_aNDoF);
      stoptime = clock();
      time_solve = (stoptime-starttime)/double(CLOCKS_PER_SEC);
      ATH_MSG_INFO(" - time spent solving the system: "<<time_solve<<" s");
      if(m_logStream) {
        *m_logStream<<"time spent for diagonalization: "<<time_diag<<" s\n";
        *m_logStream<<"time spent for post-solving: "<<time_solve<<" s\n";
      }
    }
    else {
      ATH_MSG_ERROR("Problem in diagonalization. Solving skipped.");
      if(m_logStream)
        *m_logStream<<"time spent for diagonalization: "<<time_diag<<" s\n";
    }
 
    if(m_logStream) {
      *m_logStream<<"total time spent in solve: "<<time_diag+time_solve<<" s\n";
      *m_logStream<<"\n";
    }
 
    delete d2Chi2;
    delete aBetterMat;
    delete aBetterVec;
 
    // need to do this since success return value from Lapack is 0
    // and from solveLapack() it is 1
    if (info==0)
      info = 1;
 
    return info;
  }

solveLocal()

int MatrixTool::solveLocal ( )

private

Definition at line 501 of file MatrixTool.cxx.

  {
    ATH_MSG_INFO("solving using Local method");
    if(m_logStream) {
      *m_logStream<<"*************************************************************\n";
      *m_logStream<<"**************  solving using Local method  *****************\n";
      *m_logStream<<"*************************************************************\n";
    }
 
    int    totalNDoF(0);
    double totalChi2(0.);
 
    const AlignModuleList * alignModules = m_alignModuleTool->alignModules1D();
 
    AlignModuleList::const_iterator imod     = alignModules->begin();
    AlignModuleList::const_iterator imod_end = alignModules->end();
    for( ; imod!=imod_end; ++imod) {
 
      AlignModule * module = *imod;
 
      ATH_MSG_INFO("Solving for module: "<<module->name());
 
      DataVector<AlignPar> * alignPars = m_alignModuleTool->getAlignPars(module);
 
      int thisNDoF = alignPars->size();
 
      CLHEP::HepSymMatrix d2Chi2(thisNDoF,0);
      CLHEP::HepVector dChi2(thisNDoF,0);
      for (int i=0;i<thisNDoF;++i) {
        int ipar = alignPars->at(i)->index();
        dChi2[i] = (*m_bigvector)[ipar];
        for (int j=0;j<thisNDoF;++j) {
          int jpar = alignPars->at(j)->index();
          d2Chi2[i][j] = (*m_bigmatrix)[ipar][jpar];
        }
      }
 
      ATH_MSG_DEBUG("First derivatives:" << dChi2);
      ATH_MSG_DEBUG("Second derivatives:" << d2Chi2);
 
      if(module->nHits() < m_minNumHits || module->nTracks() < m_minNumTrks) {
        ATH_MSG_INFO("Not enough hits in module \'"<<module->name()<<"\':  "
         <<module->nHits()<<" < "<<m_minNumHits<<" or "
         <<module->nTracks()<<" <  "<<m_minNumTrks
         <<". Skipping...");
 
        // print module summary even when solving is not done
        if(m_logStream)
          printModuleSolution(*m_logStream,module,nullptr);
 
        continue;
      }
 
      ATH_MSG_DEBUG("First derivatives:" << dChi2);
      ATH_MSG_DEBUG("Second derivatives:" << d2Chi2);
 
 
      totalNDoF += thisNDoF;
 
      CLHEP::HepSymMatrix cov(d2Chi2);
 
      // invert the matrix
      int ierr = 0;
      cov.invert(ierr);
      if(ierr>0)
        ATH_MSG_WARNING("CLHEP inversion status flag = "<<ierr);
      else
        ATH_MSG_DEBUG("CLHEP inversion status flag = "<<ierr);
 
      // calculate corrections
      CLHEP::HepVector delta = cov * dChi2;
      //ATH_MSG_DEBUG("d2Chi2: "<<d2Chi2);
      //ATH_MSG_DEBUG("cov: "<<cov);
      //ATH_MSG_DEBUG("d2Chi2*cov: "<<d2Chi2*cov);
      //ATH_MSG_DEBUG("dChi2: "<<dChi2);
      ATH_MSG_DEBUG("Result: "<<delta);
 
      ATH_MSG_DEBUG("Alignment constants:");
      for (int idof=0;idof<thisNDoF;++idof) {
        AlignPar * alignPar = alignPars->at(idof);
 
        double sigma = alignPar->sigma();
        double param = -delta[idof] * sigma;
        double err = std::sqrt(2.*std::fabs(cov[idof][idof])) * sigma;
 
        ATH_MSG_DEBUG(idof <<" : "<< param << " +/- "<< err);
        ATH_MSG_DEBUG("cov("<<idof<<")="<<cov[idof][idof]<<", sigma: "<<sigma);
        ATH_MSG_DEBUG("init par: "<<alignPar->initPar());
 
        ATH_MSG_DEBUG("Filling constants obtained using Local method");
        alignPar->setPar(param,err);
        ATH_MSG_DEBUG("set param to "<<param<<" for alignPar "<<alignPar);
        ATH_MSG_DEBUG(*alignPar);
      }
 
      if(m_logStream) {
        printModuleSolution(*m_logStream,module,&cov);
 
        *m_logStream<<"CLHEP inversion status flag = "<<ierr<<"\n";
 
        // calculate chi2 of the alignment change
        double chi2 = d2Chi2.similarity(delta) * .5;
        totalChi2 += chi2;
        *m_logStream<<"delta(chi2) of the alignment change : "<<chi2<<" / "<<thisNDoF<<"\n";
      }
    }
 
    if(m_logStream) {
      *m_logStream<<"--------------------------------------------------------------------------------\n";
      *m_logStream<<"Total delta(chi2) of the alignment change from the local method : "<<totalChi2<<" / "<<totalNDoF<<"\n";
      *m_logStream<<"\n";
    }
 
    return 1;
  }

solveROOT()

int MatrixTool::solveROOT ( )

private

Definition at line 133 of file MatrixTool.cxx.

  {
    ATH_MSG_INFO("solving Global using ROOT");
    if(m_logStream) {
      *m_logStream<<"*************************************************************\n";
      *m_logStream<<"**************  solving using Global method  ****************\n";
      *m_logStream<<"**************          using ROOT           ****************\n";
      *m_logStream<<"*************************************************************\n";
    }
 
    // start measuring time
    clock_t starttime = clock();
 
    DataVector<AlignPar>*  alignParList = m_alignModuleTool->alignParList1D();
    //const AlignModuleList* alignModules = m_alignModuleTool->alignModules1D();
 
    int nDoF=m_alignModuleTool->nAlignParameters();
 
    // some debugging output
    if (msgLvl(MSG::VERBOSE)) {
      msg(MSG::VERBOSE)<<"dumping matrix and vector to screen"<<endmsg;
      for (int i=0;i<nDoF;i++)
        for (int j=0;j<nDoF;j++)
        //if (std::fabs((*m_bigmatrix)[i][j])>.0001)
          msg(MSG::VERBOSE)<<i<<", "<<j<<" : "<<(*m_bigmatrix)[i][j] <<endmsg;
 
      for (int i=0;i<nDoF;i++)
        msg(MSG::VERBOSE)<<i <<" : "<<(*m_bigvector)[i]<<endmsg;
    }
 
    // get rescaled first and second derivatives
    double * secderiv   = new double[m_aNDoF*m_aNDoF];
    double * firstderiv = new double[m_aNDoF];
    for (int iActive=0;iActive<m_aNDoF;iActive++) {
      int i = m_activeIndices[iActive];
      firstderiv[iActive] = (*m_bigvector)[i];
      for (int jActive=0;jActive<m_aNDoF;jActive++) {
        int j = m_activeIndices[jActive];
        secderiv[iActive*m_aNDoF+jActive] = (*m_bigmatrix)[i][j];
      }
    }
 
    // attention, the dimension of matrix a and b is m_aNDoF not nDoF,
    // this means some alignment parameters have not been calculated
    // if the corresponding modules did not satify the select cut
 
    TMatrixDSym a(m_aNDoF,secderiv);
    TVectorD    b(m_aNDoF,firstderiv);
 
    if(msgLvl(MSG::DEBUG)) {
      msg(MSG::DEBUG)<<"First derivatives:"<<endmsg;
      b.Print();
      msg(MSG::DEBUG)<<"Second derivatives:"<<endmsg;
      a.Print();
    }
 
    TDecompBK c(a);
    Bool_t status;
    TMatrixDSym ainv(c.Invert(status));
 
    TVectorD r(b.GetNrows());
    if(status)
      r = c.Solve(b,status);
 
    // stop measuring time
    clock_t stoptime = clock();
    double totaltime = (stoptime-starttime)/double(CLOCKS_PER_SEC);
    ATH_MSG_INFO("Time spent in solveROOT: "<<totaltime<<" s");
 
    if(!status) {
      msg(MSG::ERROR)<<"ROOT inversion failed"<<endmsg;
      if(m_logStream) {
        *m_logStream<<"ROOT inversion failed\n";
        *m_logStream<<"\n";
      }
    }
    else {
      ATH_MSG_INFO("ROOT inversion ok");
 
      ATH_MSG_DEBUG("Alignment constants:");
      for (int iAdof=0;iAdof<m_aNDoF;iAdof++) {
 
        int idof = m_activeIndices[iAdof];
        AlignPar * alignPar=(*alignParList)[idof];
 
        double sigma = alignPar->sigma();
        double param = -r[iAdof] * sigma;
        double err = std::sqrt(2.*std::fabs(ainv(iAdof,iAdof))) * sigma;
 
        ATH_MSG_DEBUG(iAdof <<" : "<< param << " +/- "<< err);
        ATH_MSG_DEBUG("ainv("<<iAdof<<")="<<ainv(iAdof,iAdof)<<", sigma: "<<sigma);
        ATH_MSG_DEBUG("init par: "<<alignPar->initPar());
        alignPar->setPar(param,err);
        ATH_MSG_DEBUG("set param to "<<param<<" for alignPar "<<alignPar);
        ATH_MSG_DEBUG(*(*alignParList)[idof]);
      }
 
      if(m_logStream)
      {
        *m_logStream<<"ROOT inversion ok\n";
 
        printGlobalSolution(*m_logStream,&ainv);
 
        // norm of first derivative
        *m_logStream<<"norm of first derivative :            "<<sqrt(b.Norm2Sqr())<<"\n";
 
        // distance to solution
        double dist = sqrt( ( b - (a * r) ).Norm2Sqr() );
        *m_logStream<<"distance to solution :                "<<dist<<"\n";
 
        // calculate chi2 of the alignment change
        double chi2 = a.Similarity(r) * .5;
        *m_logStream<<"delta(chi2) of the alignment change : "<<chi2<<" / "<<m_aNDoF<<"\n";
 
        // time spent here
        *m_logStream<<"time spent in solve :                 "<<totaltime<<" s\n";
      }
    }
 
    delete [] secderiv;
    delete [] firstderiv;
 
    return 1;
  }

solveSparseEigen()

int MatrixTool::solveSparseEigen ( )

private

Definition at line 2149 of file MatrixTool.cxx.

  {
    ATH_MSG_INFO("solving Global using SparseEigen");
    if(m_logStream) {
      *m_logStream<<"*************************************************************\n";
      *m_logStream<<"**************  solving using Global method  ****************\n";
      *m_logStream<<"**************      using SparseEigen        ****************\n";
      *m_logStream<<"*************************************************************\n";
    }
 
    // start measuring time
    clock_t starttime = clock();
 
    DataVector<AlignPar> * alignParList = m_alignModuleTool->alignParList1D();
 
    AlSpaMat * ABetterMat = nullptr;
    bool isCopy = false;
    if ( dynamic_cast<AlSymMat*>(m_bigmatrix) ) {
      ATH_MSG_INFO("Converting Matrix Format for fast solving");
      ABetterMat = new AlSpaMat(*(dynamic_cast<AlSymMat*>(m_bigmatrix)));
      isCopy = true;
    }
    else if ( dynamic_cast<AlSpaMat*>(m_bigmatrix) ) {
      ATH_MSG_INFO("Matrix format native to the fast solving");
      ABetterMat = (dynamic_cast<AlSpaMat*>(m_bigmatrix));
    }
    else {
      ATH_MSG_ERROR("Cannot cast to neither AlSymMat nor AlSpaMat");
      return 0;
    }
 
    ATH_MSG_DEBUG("checking active indices");
 
    // use const matrix when checking for non-zero elements to avoid
    // filling of the whole matrix
    const AlSpaMat * chkMatrix = ABetterMat;
 
    AlSpaMat * aBetterMat = new AlSpaMat(m_aNDoF);
    AlVec    * aBetterVec = new AlVec(m_aNDoF);
 
      
    for (int iActive=0;iActive<m_aNDoF;iActive++) {
      int i = m_activeIndices[iActive];
      (*aBetterVec)[iActive] = (*m_bigvector)[i];
      for (int jActive=0;jActive<m_aNDoF;jActive++) {
        int j = m_activeIndices[jActive];
        // only fill if non-zero !!!
        if ( (*chkMatrix)[iActive][jActive] != 0. )
          (*aBetterMat)[iActive][jActive]=(*ABetterMat)[i][j];
      }
    }
 
    // store original vector for cross-checks
    AlVec origVec(*aBetterVec);
 
    ATH_MSG_DEBUG("running the solving");
 
    // solve
    int info = (*aBetterMat).SolveWithEigen(*aBetterVec);
 
    if(info == 0) {
      ATH_MSG_INFO("SolveWithEigen solving OK");
      if(m_logStream)
        *m_logStream<<"SolveWithEigen solving OK.\n";
    }
    else {
      ATH_MSG_ERROR( "SolveWithEigen error code (0 if OK) = "<<info );
      if(m_logStream)
        *m_logStream<<"SolveWithEigen error code (0 if OK) = "<<info<<"\n";
    }
 
    if( isCopy )
      delete ABetterMat;
    ABetterMat = nullptr;
 
    // stop measuring time
    clock_t stoptime = clock();
    double totaltime = (stoptime-starttime)/double(CLOCKS_PER_SEC);
    ATH_MSG_INFO("Time spent in SolveWithEigen: "<<totaltime<<" s");
 
    ATH_MSG_DEBUG("Alignment constants:");
    // compute alignment corrections (translations in mm and rotations in rad)
    // for solveSparseEigen variances are not calculated
    for(int i=0; i<m_aNDoF; i++) {
 
      double param = -(*aBetterVec)[i];
      double err = 0.;
 
      int idof = m_activeIndices[i];
      AlignPar * alignPar=(*alignParList)[idof];
 
      // undo the sigma scaling
      double sigma = alignPar->sigma();
      param *= sigma;
 
      ATH_MSG_DEBUG(i <<" : "<< param << " +/- "<< err);
      ATH_MSG_DEBUG("sigma: "<<sigma);
      ATH_MSG_DEBUG("init par: "<<alignPar->initPar());
      alignPar->setPar(param, err);
      ATH_MSG_DEBUG("set param to "<<param<<" for alignPar "<<alignPar);
      ATH_MSG_DEBUG(*(*alignParList)[idof]);
    }
 
    if(m_logStream) {
      CLHEP::HepSymMatrix * cov = nullptr;
      printGlobalSolution(*m_logStream, cov);
 
      // norm of first derivative
      *m_logStream<<"norm of first derivative :            "<<origVec.norm()<<"\n";
 
      // distance to solution
      double dist = ( (*aBetterMat) * (*aBetterVec) - origVec ).norm();
      *m_logStream<<"distance to solution :                "<<dist<<"\n";
 
      // calculate chi2 of the alignment change
      double chi2 = (*aBetterVec) * (*aBetterMat) * (*aBetterVec) * .5;
      *m_logStream<<"delta(chi2) of the alignment change : "<<chi2<<" / "<<m_aNDoF<<"\n";
 
      // time spent here
      *m_logStream<<"time spent in solve :                 "<<totaltime<<" s\n";
    }
 
    delete aBetterMat;
    delete aBetterVec;
 
    // need to do this since success return value from Lapack is 0
    // and from SolveWithEigen() it is 1
    if (info==0)
      info = 1;
 
    return info;
  }

spuriousRemoval()

StatusCode MatrixTool::spuriousRemoval ( )

private

Definition at line 1727 of file MatrixTool.cxx.

  {
 
    // copied from SiGlobalChi2Algs
    ATH_MSG_DEBUG("in spuriousRemoval");
 
    // compute determinant before resizing
    if (m_calDet) {
      const double tol = 1.e-20;
      double determ = m_bigmatrix->determinant();
      ATH_MSG_INFO("Determinant: " << determ);
      if (std::fabs(determ) < tol)
        ATH_MSG_WARNING("Matrix is singular!");
    }
 
    // fill vector with modules that need to be removed
    int fillvecmods=fillVecMods();
    if (fillvecmods==0) {
 
      //ATH_MSG_INFO(" No resize needed (NhitsCut = "
      //     << m_hitscut << ")");
 
      if (msgLvl(MSG::DEBUG)) {
        //m_bigmatrix->Write("bigmatrix.txt", false, m_wSqMatrix, m_scale,
        //           MatVersion);
        //m_bigvector->Write("bigvector.txt", false, m_scale, m_modcodemap,
        //           VecVersion, m_alignProcessLevel, m_fitParam);
      }
 
      return StatusCode::SUCCESS;
    }
    else if (fillvecmods==2)
      return StatusCode::FAILURE;
 
    /* this is a bit difficult to implement for now....
 
    // remove matrix/vector elements
    int cont=0;
    ModuleIndexMap::iterator itcode;
    ATH_MSG_INFO("Eliminating module...");
 
 
 
    for (std::vector<int>::const_iterator it=m_dropmods.begin();
    it!= m_dropmods.end(); ++it) {
 
    itcode = m_modcodemap.find((*it));
 
    if (itcode == m_modcodemap.end()) {
    ATH_MSG_WARNING("Could not find module " << *it << " in map.");
    return StatusCode::FAILURE;
    }
 
    ATH_MSG_INFO(" - Removing mcode: " << (itcode->second)
    << " (old index: " << (*it) << " -> new index: " << (*it)-cont << ")");
 
    m_bigmatrix->RemoveModule((*it)-cont);
    m_bigvector->RemoveModule((*it)-cont);
    cont++;
    }
    ATH_MSG_INFO("Modules removed from the Matrix and from the Vector Successfully!");
    ATH_MSG_DEBUG("(DoF: " << nDoF << ")");
 
    // resizing...
    ATH_MSG_INFO("Resizing the bigvector in memory...");
    m_bigvector->reSize(nDoF-6*m_dropmods.size());
    ATH_MSG_INFO("Resizing the bigmatrix in memory...");
    m_bigmatrix->reSize(nDoF-6*m_dropmods.size());
 
    ATH_MSG_INFO(m_dropmods.size() << " modules eliminated from the matrix, i.e, "
    << 6*m_dropmods.size() << " DoFs");
    ATH_MSG_INFO(nDoF/6 - m_dropmods.size() << " modules to align (" << nDoF-6*m_dropmods.size() << " DoFs)");
    ATH_MSG_INFO("New bigmatrix size is: " << m_bigmatrix->size());
    ATH_MSG_INFO("New bigvector size is: " << (*m_bigvector).size());
 
    NDoF = nDoF-6*(m_dropmods.size());
 
    // resize (update) nDoF variable
    nDoF = NDoF;
 
    // Resizing vectors to store results
    m_alignPar->reSize(nDoF);
    m_alignSqErr->reSize(nDoF);
 
    // Fill the mapping module map and updating m_modcodemap
    UpdateModMap();
    */
 
    return StatusCode::SUCCESS;
  }

storeInTFile()

void MatrixTool::storeInTFile

(

const TString &

filename

)

Store Files in a tfile.

Definition at line 730 of file MatrixTool.cxx.

  {
    //Store reults in a single TFile....
    //Including Matrix Vector Hitmap.. Soluton EVs etc.
    
    ATH_MSG_DEBUG("Writing Results to a TFile");
    
    DataVector<AlignPar> * alignParList = m_alignModuleTool->alignParList1D();
    int nDoF = alignParList->size();
  
    TMatrixDSparse* myTMatrix = m_bigmatrix->makeTMatrix();
 
    ATH_MSG_DEBUG( "Created TMatrixDSparse" );
 
  
    double *val = new double[nDoF];
    for (int i=0;i<nDoF;i++) {
      val[i] = (*m_bigvector)[i];
    }
  
    TVectorD myTVector(nDoF, val);
    delete [] val;
  
    ATH_MSG_DEBUG( "Created TVectorD" );
 
  
    const AlignModuleList * moduleList = m_alignModuleTool->alignModules1D();
    int nModules = moduleList->size();
  
    double *hitmapA = new double[nModules];
    double *hitmapB = new double[nModules];
    AlignModuleList::const_iterator imod     = moduleList->begin();
    AlignModuleList::const_iterator imod_end = moduleList->end();
    int index(0);
    for(; imod != imod_end; ++imod) {
      AlignModule * module = *imod;
      hitmapA[index] = (double)module->nHits();
      hitmapB[index] = (double)module->nTracks();
      index++;
    }
  
    TVectorD hitmapHits(nModules, hitmapA); 
    TVectorD hitmapTracks(nModules, hitmapB);
    
    delete [] hitmapA;
    delete [] hitmapB;
    
 
 
    TFile myFile(filename,"recreate");
    hitmapHits.Write("Hits");
    hitmapTracks.Write("Tracks");
    myTMatrix->Write("Matrix");
    myTVector.Write("Vector");
    
    TVectorD scale(1, &m_scale) ;
    scale.Write("Scale");
    
   
    double *moduleInfoA = new double[nDoF];//unsigned long long 
    double *dofInfoA = new double[nDoF];//int 
  
    if (sizeof(unsigned long long) != sizeof(double))      
      ATH_MSG_ERROR("Module Identifiers will not be saved.  sizeof(double)!=sizeof(ulonglong)");
    else{      
      
      DataVector<AlignPar>* alignPars = m_alignModuleTool->alignParList1D();
      for (int i=0;i<(int)alignPars->size();i++) {
        //Do a direct memory copy to store unsigned long long in the momory of a double
        double target;
        uint64_t  id = (*alignPars)[i]->alignModule()->identify().get_compact();
        memcpy(&target, &id, sizeof(target));
        moduleInfoA[i]=target;
        //moduleInfoB[i]=(*alignPars)[i]->alignModule()->name();
        uint64_t dof = (*alignPars)[i]->paramType();
        memcpy(&target, &dof, sizeof(target));
        dofInfoA[i]=target;
      }
  
      TVectorD moduleIDs(nDoF, moduleInfoA) ;
      TVectorD moduleDoFs(nDoF,dofInfoA);
      delete [] moduleInfoA;
      moduleIDs.Write("ModuleID");
      moduleDoFs.Write("dof");
    }
    
    myFile.Write();
    myFile.Close();
    
    delete myTMatrix;
    ATH_MSG_DEBUG("Finshed writing TFILE");
 
  }

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 >, and AthCheckedComponent<::AthAlgTool >.

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.

                                               {
    // debug() << "updateVHKA for property " << p.name() << " " << p.toString() 
    //         << "  size: " << m_vhka.size() << endmsg;
    for (auto &a : m_vhka) {
      std::vector<SG::VarHandleKey*> keys = a->keys();
      for (auto k : keys) {
        k->setOwner(this);
      }
    }
  }

writeHitmap()

void MatrixTool::writeHitmap ( )

private

Definition at line 2283 of file MatrixTool.cxx.

  {
    ATH_MSG_INFO("writing the hitmap to file");
 
    const AlignModuleList * moduleList = m_alignModuleTool->alignModules1D();
    int nModules = moduleList->size();
 
    AlMat hitmap(nModules,2);
    AlignModuleList::const_iterator imod     = moduleList->begin();
    AlignModuleList::const_iterator imod_end = moduleList->end();
    int index(0);
    for(; imod != imod_end; ++imod) {
      AlignModule * module = *imod;
      hitmap[index][0] = module->nHits();
      hitmap[index][1] = module->nTracks();
      index++;
    }
 
    // Set Path for the hitmap matrix
    hitmap.SetPathBin(m_pathbin.value()+m_prefixName.value());
    hitmap.SetPathTxt(m_pathtxt.value()+m_prefixName.value());
 
    StatusCode sc = hitmap.Write("hitmap.bin",true); // write the hitmap matrix
 
    if (sc!=StatusCode::SUCCESS)
      ATH_MSG_ERROR("Problem writing hitmap matrix");
 
    if (m_writeHitmapTxt) {
      sc = hitmap.Write("hitmap.txt", false, 0);
      if (sc!=StatusCode::SUCCESS)
        ATH_MSG_ERROR("Problem writing hitmap matrix to text file");
    }
 
    ATH_MSG_DEBUG("hitmap written to: "<< m_pathbin.value()+m_prefixName.value() <<"hitmap.bin (.txt)");
  }

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.

{}; 

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.

                                        {this, "AlignCorrDBStep", 10.,
      "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.

{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.

                                                        {
      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.

{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.

      {this, "CalculateFullCovariance", true,
       "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.

                                {this, "MatrixDet", false,
      "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.

{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.

                                     {this, "Diagonalize", true,
      "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.

                                         {this, "EigenvalueStep", 1e3,
      "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.

{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.

      {this, "InputHitmapFiles", {"hitmap.bin"},
       "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.

      {this, "InputMatrixFiles", {"matrix.bin"},
       "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.

      {this, "InputTFiles", {"AlignmentTFile.root"},
       "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.

      {this, "InputVectorFiles", {"vector.bin"},
       "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.

                                      {this, "MaxReadErrors", 10,
      "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.

                                   {this, "MinNumHitsPerModule", 0,
      "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.

                                   {this, "MinNumTrksPerModule", 0,
      "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.

                               {this, "ModCut", 0,
      "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.

                                      {this, "PathBinName", "./",
      "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.

                                      {this, "PathTxtName", "./",
      "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.

                                         {this, "PrefixName", "",
      "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.

                                  {this, "PullCut", 1.0,
      "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.

                                     {this, "ReadHitmaps", false,
      "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.

                                    {this, "ReadTFile", false,
      "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.

{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.

{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.

{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.

{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.

{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.

{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.

{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.

{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.

{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.

{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.

{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.

{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.

                                          {this, "RemoveSpurious", false,
      "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.

                                  {this, "RunLocalMethod", true,
      "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.

{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.

{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.

                                     {this, "ScaleMatrix", false,
      "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.

                                            {this, "SoftEigenmodeCut", 0.,
      "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.

{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.

{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.

{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.

                                       {this, "WriteEigenMat", true,
      "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.

                                          {this, "WriteEigenMatTxt", true,
      "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.

                                     {this, "WriteHitmap", false,
      "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.

                                        {this, "WriteHitmapTxt", false,
      "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.

                                  {this, "WriteMat", true,
      "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.

                                     {this, "WriteMatTxt", true,
      "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.

                                          {this, "WriteModuleNames", false,
      "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.

                                    {this, "WriteTFile", false,
      "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.

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

---

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

• MatrixTool.h
• MatrixTool.cxx