void test_TMatrixDSymEigen() {

  TMatrixDSym m(3);
  m(0,0) = m(1,1) = m(2,2) = 1;
  m(0,1) = m(1,0) = m(1,2) = m(2,1) = 0.01;
  m(0,2) = m(2,0) = 0.01;
  cout << "Initial matrix : " << endl;
  m.Print();

  cout << "Determinant : " << endl;
  cout << m.Determinant() << endl;

  TMatrixDSymEigen mtm(m);
  TVectorD eigenVal = mtm.GetEigenValues();
  TMatrixD eigenValMat(3,3);
  for (int i = 0; i < 3; ++i) eigenValMat(i,i) = eigenVal(i); 

  cout << "Eigenvalues (D): " << endl;
  eigenVal.Print();
  
  TMatrixD eigenVec  = mtm.GetEigenVectors();
  TMatrixD eigenVecT = mtm.GetEigenVectors();
  eigenVecT.T();
  cout << "EigenVectors (Q): " << endl;
  eigenVec.Print();
  
  
  TMatrixD mCheck = eigenVec*eigenValMat*eigenVecT;
  cout << "Q*D*QT : " << endl;
  mCheck.Print();
  
}