From db358353c7af7e999daac01b0dce3ed72be10a65 Mon Sep 17 00:00:00 2001 From: David Sarrut Date: Fri, 3 Feb 2012 07:57:12 +0100 Subject: [PATCH] Allow to combine multiple elastix transfo --- tools/clitkAffineTransform.ggo | 2 +- tools/clitkAffineTransformGenericFilter.h | 2 +- tools/clitkAffineTransformGenericFilter.txx | 139 +++++++++++--------- 3 files changed, 81 insertions(+), 62 deletions(-) diff --git a/tools/clitkAffineTransform.ggo b/tools/clitkAffineTransform.ggo index dc64bff..3cc4eae 100644 --- 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" string no +option "elastix" e "Read EulerTransform from elastix output file (combine if multiple)" string no multiple 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.h b/tools/clitkAffineTransformGenericFilter.h index f990fa5..6cd2a91 100644 --- a/tools/clitkAffineTransformGenericFilter.h +++ b/tools/clitkAffineTransformGenericFilter.h @@ -105,7 +105,7 @@ namespace clitk template typename itk::Matrix - createMatrixFromElastixFile(std::string filename); + createMatrixFromElastixFile(std::vector & filename); bool GetElastixValueFromTag(std::ifstream & is, std::string tag, std::string & value); void GetValuesFromValue(const std::string & s, diff --git a/tools/clitkAffineTransformGenericFilter.txx b/tools/clitkAffineTransformGenericFilter.txx index b2d7a53..68532a6 100644 --- a/tools/clitkAffineTransformGenericFilter.txx +++ b/tools/clitkAffineTransformGenericFilter.txx @@ -52,12 +52,12 @@ namespace clitk // Call UpdateWithDim if(Dimension==2) UpdateWithDim<2>(PixelType, Components); else - if(Dimension==3) UpdateWithDim<3>(PixelType, Components); - else if (Dimension==4)UpdateWithDim<4>(PixelType, Components); - else { - std::cout<<"Error, Only for 2, 3 or 4 Dimensions!!!"<(PixelType, Components); + else if (Dimension==4)UpdateWithDim<4>(PixelType, Components); + else { + std::cout<<"Error, Only for 2, 3 or 4 Dimensions!!!"<(); } - // else if(PixelType == "unsigned_short"){ - // if (m_Verbose) std::cout << "Launching filter in "<< Dimension <<"D and unsigned_short..." << std::endl; - // UpdateWithDimAndPixelType(); - // } + else if(PixelType == "unsigned_short"){ + if (m_Verbose) std::cout << "Launching filter in "<< Dimension <<"D and unsigned_short..." << std::endl; + UpdateWithDimAndPixelType(); + } else if (PixelType == "unsigned_char") { if (m_Verbose) std::cout << "Launching filter in "<< Dimension <<"D and unsigned_char..." << std::endl; @@ -181,7 +181,9 @@ namespace clitk } else { if (m_ArgsInfo.elastix_given) { - matrix = createMatrixFromElastixFile(m_ArgsInfo.elastix_arg); + std::vector s; + for(uint i=0; i(s); } else matrix.SetIdentity(); @@ -495,64 +497,81 @@ namespace clitk template template typename itk::Matrix - AffineTransformGenericFilter::createMatrixFromElastixFile(std::string filename) + AffineTransformGenericFilter::createMatrixFromElastixFile(std::vector & filename) { if (Dimension != 3) { FATAL("Only 3D yet" << std::endl); } typename itk::Matrix matrix; - // Open file - std::ifstream is; - clitk::openFileForReading(is, filename); - - // 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); - } - - // FIXME check - // (InitialTransformParametersFileName "NoInitialTransform") - - // Get CenterOfRotationPoint - GetElastixValueFromTag(is, "CenterOfRotationPoint ", s); // space is needed - if (!b) { - FATAL("Error must read 'CenterOfRotationPoint' in " << filename << std::endl); - } - std::vector cor; - GetValuesFromValue(s, cor); - - // Get Transformparameters - GetElastixValueFromTag(is, "TransformParameters ", s); // space is needed - if (!b) { - FATAL("Error must read 'TransformParameters' in " << filename << std::endl); - } - std::vector results; - GetValuesFromValue(s, results); - - // construct a stream from the string itk::CenteredEuler3DTransform::Pointer mat = itk::CenteredEuler3DTransform::New(); - itk::CenteredEuler3DTransform::ParametersType p; - p.SetSize(9); - 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(cor[i].c_str()); // Centre of rotation - for(uint i=0; i<3; i++) - p[i+6] = atof(results[i+3].c_str()); // Translation - mat->SetParameters(p); + itk::CenteredEuler3DTransform::Pointer previous; + for(uint j=0; j cor; + GetValuesFromValue(s, cor); + + // Get Transformparameters + GetElastixValueFromTag(is, "TransformParameters ", s); // space is needed + if (!b) { + FATAL("Error must read 'TransformParameters' in " << filename[j] << std::endl); + } + std::vector results; + GetValuesFromValue(s, results); + + // construct a stream from the string + itk::CenteredEuler3DTransform::ParametersType p; + p.SetSize(9); + 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(cor[i].c_str()); // Centre of rotation + for(uint i=0; i<3; i++) + p[i+6] = atof(results[i+3].c_str()); // Translation + mat->SetParameters(p); - if (m_Verbose) { - std::cout << "Rotation (deg) : " << rad2deg(p[0]) << " " << rad2deg(p[1]) << " " << rad2deg(p[2]) << std::endl; - std::cout << "Translation (phy) : " << p[3] << " " << p[4] << " " << p[5] << std::endl; - std::cout << "Center of rot (phy) : " << p[6] << " " << p[7] << " " << p[8] << std::endl; + if (m_Verbose) { + std::cout << "Rotation (deg) : " << rad2deg(p[0]) << " " << rad2deg(p[1]) << " " << rad2deg(p[2]) << std::endl; + std::cout << "Center of rot (phy) : " << p[3] << " " << p[4] << " " << p[5] << std::endl; + std::cout << "Translation (phy) : " << p[6] << " " << p[7] << " " << p[8] << std::endl; + } + + // Compose with previous if needed + if (j!=0) { + mat->Compose(previous); + if (m_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(); } + mat = previous; for(uint i=0; i<3; i++) for(uint j=0; j<3; j++) matrix[i][j] = mat->GetMatrix()[i][j]; @@ -582,7 +601,7 @@ namespace clitk value.erase (std::remove (value.begin(), value.end(), ')'), value.end()); return true; } - } + } return false; } //------------------------------------------------------------------- -- 2.46.0