[Insight-users] Inverse Transform doesn't not catch singular matrices

[Francois Budin]

Hi all,

I was trying to compute the inverse of an affine transform (
itk::AffineTransform<double,3> ) using the following function:
bool GetInverse (Self *inverse) const

I couldn't succeed to detect when the matrix was non-invertible. I
tried to verify the value returned in the boolean but it always seems
to be 1. I also tried to detect if any exception were thrown but it
didn't catch anything.
My code is the following:


    bool noninvertible=0;
    itk::Transform<double,3,3>::Pointer inverseTransform;
    try
    {
     typename AffineTransformType::Pointer affine =
dynamic_cast<AffineTransformType* > ( transform.GetPointer() ) ;
    if(affine )
    {
      typename AffineTransformType::Pointer
affinetemp=AffineTransformType::New();
      std::cout<<affine->GetInverse(affinetemp)<<std::endl;//should
print 0 if non-invertible transform;
      inverseTransform=affinetemp;
    }
  else
    {
      typename RotationType::Pointer rigid =
dynamic_cast<RotationType* > ( transform.GetPointer() ) ;
      if(rigid)
      {
        typename RotationType::Pointer rigidtemp=RotationType::New();
        std::cout<<rigid->GetInverse(rigidtemp)<<std::endl;
        inverseTransform=rigidtemp;
      }
      else
      {
        noninvertible=1;//I don't need to handle the other kinds of transform
      }
    }
    }
    catch(...)
    {
      std::cout<<"Exception Detected"<<std::endl;
    noninvertible=1;
    }


Then I tried to write a small piece of code which copies what ITK is
doing to inverse an affine transform to test the inversion of the
matrix itself:

#include <itkMatrix.h>

int main(int argnb,char* arg[])
{
itk::Matrix<double,3,3> mat;
mat.Fill(0);
mat[0][0]=-1;
try
      {
     std::cout<<mat.GetInverse()<<std::endl;
      }
    catch(...)
      {
      std::cout<<"Singular"<<std::endl;
      }

return 0;
}

Nothing got caught, so I suspect that there is a problem somewhere and
that the vnl functions used to inverse the matrix does not throw any
exception when it is singular. If I am wrong, could anybody tell me
what I should change in my code?
Sincerely,

Francois