Improved Elastix conversion

diff --git a/common/clitkElastix.h b/common/clitkElastix.h
index 738181fd6b9595f0951f6ff7bdf213436701b5c2..bc333caef5d48d4c3f681a1cc0788ef5e394579b 100644 (file)
--- a/common/clitkElastix.h
+++ b/common/clitkElastix.h
@@ -66,88 +66,80 @@ GetValuesFromValue(const std::string & s,
 //-------------------------------------------------------------------
 template<unsigned int Dimension>
 typename itk::Matrix<double, Dimension+1, Dimension+1>
-createMatrixFromElastixFile(std::vector<std::string> & filename, bool verbose=true) {
+createMatrixFromElastixFile(std::string& filename, bool verbose=true) {
   if (Dimension != 3) {
     FATAL("Only 3D yet" << std::endl);
   }
-  typename itk::Matrix<double, Dimension+1, Dimension+1> matrix;
+  typename itk::Matrix<double, Dimension+1, Dimension+1> matrix, init;
 
   itk::Euler3DTransform<double>::Pointer mat = itk::Euler3DTransform<double>::New();
   itk::Euler3DTransform<double>::Pointer previous;
-  for(uint j=0; j<filename.size(); j++) {
-
-    // Open file
-    if (verbose) std::cout << "Read elastix parameters in " << filename[j] << std::endl;
-    std::ifstream is;
-    clitk::openFileForReading(is, filename[j]);
-
-    // Check Transform
-    std::string s;
-    bool b = GetElastixValueFromTag(is, "Transform ", s);
-    if (!b) {
-      FATAL("Error must read 'Transform' in " << filename[j] << std::endl);
-    }
-    if (s != "EulerTransform") {
-      FATAL("Sorry only 'EulerTransform'" << std::endl);
-    }
 
-    // FIXME check
-    //    (InitialTransformParametersFilename[j] "NoInitialTransform")
+  // Open file
+  if (verbose) std::cout << "Read elastix parameters in " << filename << std::endl;
+  std::ifstream is;
+  clitk::openFileForReading(is, filename);
 
-    // Get CenterOfRotationPoint
-    GetElastixValueFromTag(is, "CenterOfRotationPoint ", s); // space is needed
-    if (!b) {
-      FATAL("Error must read 'CenterOfRotationPoint' in " << filename[j] << std::endl);
-    }
-    std::vector<std::string> cor;
-    GetValuesFromValue(s, cor);
-    itk::Euler3DTransform<double>::CenterType c;
-    for(uint i=0; i<3; i++)
-      c[i] = atof(cor[i].c_str());
-    mat->SetCenter(c);
-
-    // Get Transformparameters
-    GetElastixValueFromTag(is, "ComputeZYX ", s); // space is needed
-    mat->SetComputeZYX( s==std::string("true") );
-
-    // Get Transformparameters
-    GetElastixValueFromTag(is, "TransformParameters ", s); // space is needed
-    if (!b) {
-      FATAL("Error must read 'TransformParameters' in " << filename[j] << std::endl);
-    }
-    std::vector<std::string> results;
-    GetValuesFromValue(s, results);
-
-    // construct a stream from the string
-    itk::Euler3DTransform<double>::ParametersType p;
-    p.SetSize(6);
-    for(uint i=0; i<3; i++)
-      p[i] = atof(results[i].c_str()); // Rotation
-    for(uint i=0; i<3; i++)
-      p[i+3] = atof(results[i+3].c_str()); // Translation
-    mat->SetParameters(p);
-
-    if (verbose) {
-      std::cout << "Rotation      (deg) : " << rad2deg(p[0]) << " " << rad2deg(p[1]) << " " << rad2deg(p[2]) << std::endl;
-      std::cout << "Center of rot (phy) : " << c[0] << " " << c[1] << " " << c[2] << std::endl;
-      std::cout << "Translation   (phy) : " << p[3] << " " << p[4] << " " << p[5] << std::endl;
-    }
+  // Check Transform
+  std::string s;
+  bool b = GetElastixValueFromTag(is, "Transform ", s);
+  if (!b) {
+    FATAL("Error must read 'Transform' in " << filename << std::endl);
+  }
+  if (s != "EulerTransform") {
+    FATAL("Sorry only 'EulerTransform'" << std::endl);
+  }
 
-    // Compose with previous if needed
-    if (j!=0) {
-      mat->Compose(previous);
-      if (verbose) {
-        std::cout << "Composed rotation      (deg) : " << rad2deg(mat->GetAngleX()) << " " << rad2deg(mat->GetAngleY()) << " " << rad2deg(mat->GetAngleZ()) << std::endl;
-        std::cout << "Composed center of rot (phy) : " << mat->GetCenter() << std::endl;
-        std::cout << "Compsoed translation   (phy) : " << mat->GetTranslation() << std::endl;
-      }
-    }
-    // previous = mat->Clone(); // ITK4
-    previous = itk::Euler3DTransform<double>::New();
-    previous->SetParameters(mat->GetParameters());
-    previous->SetCenter(c);
-    previous->SetComputeZYX(mat->GetComputeZYX());
+  // Get previous
+  b = GetElastixValueFromTag(is, "InitialTransformParametersFileName ", s);
+  if(s == "NoInitialTransform")
+    init.SetIdentity();
+  else
+    init = createMatrixFromElastixFile<Dimension>(s, verbose);
+
+  // Get CenterOfRotationPoint
+  GetElastixValueFromTag(is, "CenterOfRotationPoint ", s); // space is needed
+  if (!b) {
+    FATAL("Error must read 'CenterOfRotationPoint' in " << filename << std::endl);
+  }
+  std::vector<std::string> cor;
+  GetValuesFromValue(s, cor);
+  itk::Euler3DTransform<double>::CenterType c;
+  for(uint i=0; i<3; i++)
+    c[i] = atof(cor[i].c_str());
+  mat->SetCenter(c);
+
+  // Get Transformparameters
+  GetElastixValueFromTag(is, "ComputeZYX ", s); // space is needed
+  mat->SetComputeZYX( s==std::string("true") );
+
+  // Get Transformparameters
+  GetElastixValueFromTag(is, "TransformParameters ", s); // space is needed
+  if (!b) {
+    FATAL("Error must read 'TransformParameters' in " << filename << std::endl);
   }
+  std::vector<std::string> results;
+  GetValuesFromValue(s, results);
+
+  // construct a stream from the string
+  itk::Euler3DTransform<double>::ParametersType p;
+  p.SetSize(6);
+  for(uint i=0; i<3; i++)
+    p[i] = atof(results[i].c_str()); // Rotation
+  for(uint i=0; i<3; i++)
+    p[i+3] = atof(results[i+3].c_str()); // Translation
+  mat->SetParameters(p);
+
+  if (verbose) {
+    std::cout << "Rotation      (deg) : " << rad2deg(p[0]) << " " << rad2deg(p[1]) << " " << rad2deg(p[2]) << std::endl;
+    std::cout << "Center of rot (phy) : " << c[0] << " " << c[1] << " " << c[2] << std::endl;
+    std::cout << "Translation   (phy) : " << p[3] << " " << p[4] << " " << p[5] << std::endl;
+  }
+
+  previous = itk::Euler3DTransform<double>::New();
+  previous->SetParameters(mat->GetParameters());
+  previous->SetCenter(c);
+  previous->SetComputeZYX(mat->GetComputeZYX());
 
   mat = previous;
   for(uint i=0; i<3; i++)
@@ -159,7 +151,7 @@ createMatrixFromElastixFile(std::vector<std::string> & filename, bool verbose=tr
   matrix[2][3] = mat->GetOffset()[2];
   matrix[3][3] = 1;
 
-  return matrix;
+  return matrix*init;
 }
 }
 //-------------------------------------------------------------------

diff --git a/common/vvImageReader.cxx b/common/vvImageReader.cxx
index 23c78f155dd9e4c51ddad251cf87c549e6d7b70c..2c094a8c3eda3cfb67e85acddafb9cceb1bf562c 100644 (file)
--- a/common/vvImageReader.cxx
+++ b/common/vvImageReader.cxx
@@ -188,9 +188,7 @@ void vvImageReader::ReadMatImageTransform()
   f.open(filename.c_str());
   if(!itkMatRead && f.is_open()) {
     itkMatRead = true;
-    std::vector<std::string> l;
-    l.push_back(filename);
-    itkMat = clitk::createMatrixFromElastixFile<3>(l, true);
+    itkMat = clitk::createMatrixFromElastixFile<3>(filename, true);
   }
   f.close();
 

diff --git a/tools/clitkAffineTransform.ggo b/tools/clitkAffineTransform.ggo
index 3cc4eae503b299a8d23e820340870928159a76c2..ead649ff4eb809b5e5789023921685b72257fb65 100644 (file)
--- a/tools/clitkAffineTransform.ggo
+++ b/tools/clitkAffineTransform.ggo
@@ -18,7 +18,7 @@ option "spacing"      -       "New output spacing if different from input"    double  no      mul
 option "spacinglike"   -       "New output spacing like this image"            string  no
 option "origin"                -       "New output origin if different from input"     double  no      multiple
 option "matrix"                m       "Affine matrix (homogene) filename"             string  no
-option "elastix"       e       "Read EulerTransform from elastix output file (combine if multiple)"    string  no      multiple
+option "elastix"       e       "Read EulerTransform from elastix output file (combine if multiple)"    string  no
 option "rotate"                r       "Rotation to apply (radians)"                   double  no      multiple
 option "translate"     t       "Translation to apply (mm)"                     double  no      multiple
 option "pad"           -       "Edge padding value"                            double  no      default="0.0"

diff --git a/tools/clitkAffineTransformGenericFilter.txx b/tools/clitkAffineTransformGenericFilter.txx
index 150c8c172187c9455af05d8a7861e05d34c46855..1fa58b582919d4d7b03447f4be7aecb14a8738ac 100644 (file)
--- a/tools/clitkAffineTransformGenericFilter.txx
+++ b/tools/clitkAffineTransformGenericFilter.txx
@@ -182,9 +182,8 @@ namespace clitk
           }
         else {
           if (m_ArgsInfo.elastix_given) {
-            std::vector<std::string> s;
-            for(uint i=0; i<m_ArgsInfo.elastix_given; i++) s.push_back(m_ArgsInfo.elastix_arg[i]);
-            matrix = createMatrixFromElastixFile<Dimension>(s, m_Verbose);
+            std::string filename(m_ArgsInfo.elastix_arg);
+            matrix = createMatrixFromElastixFile<Dimension>(filename, m_Verbose);
           }
           else 
             matrix.SetIdentity();

diff --git a/tools/clitkElastixTransformToMatrix.cxx b/tools/clitkElastixTransformToMatrix.cxx
index 7efc4e965161c8171f80f1d4157431fcdb02abd2..08113b2441481f38b7b3a3ca04bcab86d1cad160 100644 (file)
--- a/tools/clitkElastixTransformToMatrix.cxx
+++ b/tools/clitkElastixTransformToMatrix.cxx
@@ -31,9 +31,9 @@ int main(int argc, char * argv[])
   CLITK_INIT;
 
   // Use static fct of AffineTransformGenericFilter
-  std::vector<std::string> l;
-  l.push_back(args_info.input_arg);
-  itk::Matrix<double, 4, 4> m = clitk::createMatrixFromElastixFile<3>(l, args_info.verbose_flag);
+  std::string filename(args_info.input_arg);
+  itk::Matrix<double, 4, 4> m = clitk::createMatrixFromElastixFile<3>(filename,
+                                                                      args_info.verbose_flag);
 
   // Print matrix
   std::ofstream os;