mitk: mitkGeometry3DTest.cpp Source File

Go to the documentation of this file.
00001 /*=========================================================================
00002 
00003 Program:   Medical Imaging & Interaction Toolkit
00004 Language:  C++
00005 Date:      $Date$
00006 Version:   $Revision$
00007 
00008 Copyright (c) German Cancer Research Center, Division of Medical and
00009 Biological Informatics. All rights reserved.
00010 See MITKCopyright.txt or https://www.mitk.org/copyright.html for details.
00011 
00012 This software is distributed WITHOUT ANY WARRANTY; without even
00013 the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
00014 PURPOSE.  See the above copyright notices for more information.
00015 
00016 =========================================================================*/
00017 
00018 #include "mitkGeometry3D.h"
00019 
00020 #include <vnl/vnl_quaternion.h>
00021 #include <vnl/vnl_quaternion.txx>
00022 
00023 #include "mitkRotationOperation.h"
00024 #include "mitkInteractionConst.h"
00025 #include <mitkMatrixConvert.h>
00026 
00027 #include "mitkTestingMacros.h"
00028 #include <fstream>
00029 
00030 bool testGetAxisVectorVariants(mitk::Geometry3D* geometry)
00031 {
00032   int direction;
00033   for(direction=0; direction<3; ++direction)
00034   {
00035     mitk::Vector3D frontToBack;
00036     switch(direction)
00037     {
00038       case 0: frontToBack = geometry->GetCornerPoint(false, false, false)-geometry->GetCornerPoint(true , false, false); break; //7-3
00039       case 1: frontToBack = geometry->GetCornerPoint(false, false, false)-geometry->GetCornerPoint(false, true , false); break; //7-5
00040       case 2: frontToBack = geometry->GetCornerPoint(false, false, false)-geometry->GetCornerPoint(false , false, true);  break; //7-2
00041     }
00042     std::cout << "Testing GetAxisVector(int) vs GetAxisVector(bool, bool, bool): ";
00043     if(mitk::Equal(geometry->GetAxisVector(direction), frontToBack) == false)
00044     {
00045       std::cout<<"[FAILED]"<<std::endl;
00046       return false;
00047     }
00048     std::cout<<"[PASSED]"<<std::endl;
00049   }
00050   return true;
00051 }
00052 
00053 bool testGetAxisVectorExtent(mitk::Geometry3D* geometry)
00054 {
00055   int direction;
00056   for(direction=0; direction<3; ++direction)
00057   {
00058     if(mitk::Equal(geometry->GetAxisVector(direction).GetNorm(), geometry->GetExtentInMM(direction)) == false)
00059     {
00060       std::cout<<"[FAILED]"<<std::endl;
00061       return false;
00062     }
00063     std::cout<<"[PASSED]"<<std::endl;
00064   }
00065   return true;
00066 }
00067 
00068 // a part of the test requires axis-parallel coordinates
00069 int testIndexAndWorldConsistency(mitk::Geometry3D* geometry3d)
00070 {
00071   MITK_TEST_OUTPUT( << "Testing consistency of index and world coordinate systems: ");
00072   mitk::Point3D origin = geometry3d->GetOrigin();
00073   mitk::Point3D dummy;
00074 
00075   MITK_TEST_OUTPUT( << " Testing index->world->index conversion consistency");
00076   geometry3d->WorldToIndex(origin, dummy);
00077   geometry3d->IndexToWorld(dummy, dummy);
00078   MITK_TEST_CONDITION_REQUIRED(dummy == origin, "");
00079 
00080   MITK_TEST_OUTPUT( << " Testing WorldToIndex(origin, mitk::Point3D)==(0,0,0)");
00081   mitk::Point3D globalOrigin;
00082   mitk::FillVector3D(globalOrigin, 0,0,0);
00083   
00084   mitk::Point3D originContinuousIndex;
00085   geometry3d->WorldToIndex(origin, originContinuousIndex);
00086   MITK_TEST_CONDITION_REQUIRED(originContinuousIndex == globalOrigin, "");
00087 
00088   MITK_TEST_OUTPUT( << " Testing WorldToIndex(origin, itk::Index)==(0,0,0)");
00089   itk::Index<3> itkindex;
00090   geometry3d->WorldToIndex(origin, itkindex);
00091   itk::Index<3> globalOriginIndex;
00092   mitk::vtk2itk(globalOrigin, globalOriginIndex);
00093   MITK_TEST_CONDITION_REQUIRED(itkindex == globalOriginIndex, "");
00094 
00095   MITK_TEST_OUTPUT( << " Testing WorldToIndex(origin-0.5*spacing, itk::Index)==(0,0,0)");
00096   mitk::Vector3D halfSpacingStep = geometry3d->GetSpacing()*0.5;
00097   mitk::Matrix3D rotation;
00098   mitk::Point3D originOffCenter = origin-halfSpacingStep;
00099   geometry3d->WorldToIndex(originOffCenter, itkindex);
00100   MITK_TEST_CONDITION_REQUIRED(itkindex == globalOriginIndex, "");
00101 
00102   MITK_TEST_OUTPUT( << " Testing WorldToIndex(origin+0.5*spacing-eps, itk::Index)==(0,0,0)");  
00103   originOffCenter = origin+halfSpacingStep;
00104   originOffCenter -= 0.0001;
00105   geometry3d->WorldToIndex( originOffCenter, itkindex);
00106   MITK_TEST_CONDITION_REQUIRED(itkindex == globalOriginIndex, "");
00107 
00108   MITK_TEST_OUTPUT( << " Testing WorldToIndex(origin+0.5*spacing, itk::Index)==(1,1,1)");
00109   originOffCenter = origin+halfSpacingStep;
00110   itk::Index<3> global111;
00111   mitk::FillVector3D(global111, 1,1,1);
00112   geometry3d->WorldToIndex( originOffCenter, itkindex);
00113   MITK_TEST_CONDITION_REQUIRED(itkindex == global111, "");
00114 
00115   MITK_TEST_OUTPUT( << " Testing WorldToIndex(GetCenter())==BoundingBox.GetCenter: ");
00116   mitk::Point3D center = geometry3d->GetCenter();
00117   mitk::Point3D centerContIndex;
00118   geometry3d->WorldToIndex(center, centerContIndex);
00119   mitk::BoundingBox::ConstPointer boundingBox = geometry3d->GetBoundingBox();
00120   mitk::BoundingBox::PointType centerBounds = boundingBox->GetCenter();
00121   MITK_TEST_CONDITION_REQUIRED(mitk::Equal(centerContIndex,centerBounds), "");
00122 
00123   MITK_TEST_OUTPUT( << " Testing GetCenter()==IndexToWorld(BoundingBox.GetCenter): ");
00124   center = geometry3d->GetCenter();
00125   mitk::Point3D centerBoundsInWorldCoords;
00126   geometry3d->IndexToWorld(centerBounds, centerBoundsInWorldCoords);
00127   MITK_TEST_CONDITION_REQUIRED(mitk::Equal(center,centerBoundsInWorldCoords), "");
00128 
00129   return EXIT_SUCCESS;
00130 }
00131 
00132 #include <itkImage.h>
00133 
00134 int testItkImageIsCenterBased()
00135 {
00136   MITK_TEST_OUTPUT(<< "Testing whether itk::Image coordinates are center-based.");
00137   typedef itk::Image<int,3> ItkIntImage3D;
00138   ItkIntImage3D::Pointer itkintimage = ItkIntImage3D::New();
00139   ItkIntImage3D::SizeType size;
00140   size.Fill(10);
00141   mitk::Point3D origin;
00142   mitk::FillVector3D(origin, 2,3,7);
00143   itkintimage->Initialize();
00144   itkintimage->SetRegions(size);
00145   itkintimage->SetOrigin(origin);
00146   std::cout<<"[PASSED]"<<std::endl;
00147 
00148   MITK_TEST_OUTPUT( << " Testing itk::Image::TransformPhysicalPointToContinuousIndex(origin)==(0,0,0)");
00149   mitk::Point3D globalOrigin;
00150   mitk::FillVector3D(globalOrigin, 0,0,0);
00151 
00152   itk::ContinuousIndex<mitk::ScalarType, 3> originContinuousIndex;
00153   itkintimage->TransformPhysicalPointToContinuousIndex(origin, originContinuousIndex);
00154   MITK_TEST_CONDITION_REQUIRED(originContinuousIndex == globalOrigin, "");
00155 
00156   MITK_TEST_OUTPUT( << " Testing itk::Image::TransformPhysicalPointToIndex(origin)==(0,0,0)");
00157   itk::Index<3> itkindex;
00158   itkintimage->TransformPhysicalPointToIndex(origin, itkindex);
00159   itk::Index<3> globalOriginIndex;
00160   mitk::vtk2itk(globalOrigin, globalOriginIndex);
00161   MITK_TEST_CONDITION_REQUIRED(itkindex == globalOriginIndex, "");
00162 
00163   MITK_TEST_OUTPUT( << " Testing itk::Image::TransformPhysicalPointToIndex(origin-0.5*spacing)==(0,0,0)");
00164   mitk::Vector3D halfSpacingStep = itkintimage->GetSpacing()*0.5;
00165   mitk::Matrix3D rotation;
00166   mitk::Point3D originOffCenter = origin-halfSpacingStep;
00167   itkintimage->TransformPhysicalPointToIndex(originOffCenter, itkindex);
00168   MITK_TEST_CONDITION_REQUIRED(itkindex == globalOriginIndex, "");
00169 
00170   MITK_TEST_OUTPUT( << " Testing itk::Image::TransformPhysicalPointToIndex(origin+0.5*spacing-eps, itk::Index)==(0,0,0)");  
00171   originOffCenter = origin+halfSpacingStep;
00172   originOffCenter -= 0.0001;
00173   itkintimage->TransformPhysicalPointToIndex( originOffCenter, itkindex);
00174   MITK_TEST_CONDITION_REQUIRED(itkindex == globalOriginIndex, "");
00175 
00176   MITK_TEST_OUTPUT( << " Testing itk::Image::TransformPhysicalPointToIndex(origin+0.5*spacing, itk::Index)==(1,1,1)");
00177   originOffCenter = origin+halfSpacingStep;
00178   itk::Index<3> global111;
00179   mitk::FillVector3D(global111, 1,1,1);
00180   itkintimage->TransformPhysicalPointToIndex( originOffCenter, itkindex);
00181   MITK_TEST_CONDITION_REQUIRED(itkindex == global111, "");
00182 
00183   MITK_TEST_OUTPUT( << "=> Yes, itk::Image coordinates are center-based.");
00184 
00185   return EXIT_SUCCESS;
00186 }
00187 
00188 int testGeometry3D(bool imageGeometry)
00189 {
00190   float bounds[ ] = {-10.0, 17.0, -12.0, 188.0, 13.0, 211.0};
00191 
00192   mitk::Geometry3D::Pointer geometry3d = mitk::Geometry3D::New();
00193 
00194   std::cout << "Initializing: ";
00195   geometry3d->Initialize();
00196   std::cout<<"[PASSED]"<<std::endl;
00197 
00198   MITK_TEST_OUTPUT(<< "Setting ImageGeometry to " << imageGeometry);
00199   geometry3d->SetImageGeometry(imageGeometry);
00200   std::cout<<"[PASSED]"<<std::endl;
00201 
00202   std::cout << "Setting bounds by SetFloatBounds(): ";
00203   geometry3d->SetFloatBounds(bounds);
00204   std::cout<<"[PASSED]"<<std::endl;
00205  
00206   if(testGetAxisVectorVariants(geometry3d) == false)
00207     return EXIT_FAILURE;
00208 
00209   if(testGetAxisVectorExtent(geometry3d) == false)
00210     return EXIT_FAILURE;
00211 
00212   std::cout << "Creating an AffineTransform3D transform: ";
00213   mitk::AffineTransform3D::MatrixType matrix;
00214   matrix.SetIdentity();
00215   matrix(1,1) = 2;
00216   mitk::AffineTransform3D::Pointer transform;
00217   transform = mitk::AffineTransform3D::New();
00218   transform->SetMatrix(matrix);
00219   std::cout<<"[PASSED]"<<std::endl;
00220 
00221   std::cout << "Testing a SetIndexToWorldTransform: ";
00222   geometry3d->SetIndexToWorldTransform(transform);
00223   std::cout<<"[PASSED]"<<std::endl;
00224 
00225   std::cout << "Testing correctness of value returned by GetSpacing: ";
00226   const mitk::Vector3D& spacing1 = geometry3d->GetSpacing();
00227   mitk::Vector3D expectedSpacing;
00228   expectedSpacing.Fill(1.0);
00229   expectedSpacing[1] = 2;
00230   if( mitk::Equal(spacing1, expectedSpacing) == false )
00231   {
00232       std::cout<<"[FAILED]"<<std::endl;
00233       return EXIT_FAILURE;
00234   }
00235   std::cout<<"[PASSED]"<<std::endl;
00236 
00237   std::cout << "Testing a Compose(transform): ";
00238   geometry3d->Compose(transform);
00239   std::cout<<"[PASSED]"<<std::endl;
00240 
00241   std::cout << "Testing correctness of value returned by GetSpacing: ";
00242   const mitk::Vector3D& spacing2 = geometry3d->GetSpacing();
00243   expectedSpacing[1] = 4;
00244   if( mitk::Equal(spacing2, expectedSpacing) == false )
00245   {
00246       std::cout<<"[FAILED]"<<std::endl;
00247       return EXIT_FAILURE;
00248   }
00249   std::cout<<"[PASSED]"<<std::endl;
00250 
00251   std::cout << "Testing correctness of SetSpacing: ";
00252   mitk::Vector3D newspacing;
00253   mitk::FillVector3D(newspacing, 1.5, 2.5, 3.5);
00254   geometry3d->SetSpacing(newspacing);
00255   const mitk::Vector3D& spacing3 = geometry3d->GetSpacing();
00256   if( mitk::Equal(spacing3, newspacing) == false )
00257   {
00258       std::cout<<"[FAILED]"<<std::endl;
00259       return EXIT_FAILURE;
00260   }
00261   std::cout<<"[PASSED]"<<std::endl;
00262 
00263   testIndexAndWorldConsistency(geometry3d);
00264 
00265   std::cout << "Testing a rotation of the geometry: ";
00266   double angle = 35.0;
00267   mitk::Vector3D rotationVector; mitk::FillVector3D( rotationVector, 1, 0, 0 );
00268   mitk::Point3D center = geometry3d->GetCenter();
00269   mitk::RotationOperation* op = new mitk::RotationOperation( mitk::OpROTATE, center, rotationVector, angle );
00270   geometry3d->ExecuteOperation(op);
00271 
00272   MITK_TEST_OUTPUT( << " Testing mitk::GetRotation() and success of rotation");
00273   mitk::Matrix3D rotation;
00274   mitk::GetRotation(geometry3d, rotation);
00275   mitk::Vector3D voxelStep=rotation*newspacing;  
00276   mitk::Vector3D voxelStepIndex;
00277   geometry3d->WorldToIndex(geometry3d->GetOrigin(), voxelStep, voxelStepIndex);
00278   mitk::Vector3D expectedVoxelStepIndex;
00279   expectedVoxelStepIndex.Fill(1);
00280   MITK_TEST_CONDITION_REQUIRED(mitk::Equal(voxelStepIndex,expectedVoxelStepIndex), "");
00281   delete op;
00282   std::cout<<"[PASSED]"<<std::endl;
00283 
00284   MITK_TEST_OUTPUT( << "Testing that ImageGeometry is still " << imageGeometry);
00285   MITK_TEST_CONDITION_REQUIRED(geometry3d->GetImageGeometry() == imageGeometry, "");
00286 
00287   return EXIT_SUCCESS;
00288 }
00289 
00290 int mitkGeometry3DTest(int /*argc*/, char* /*argv*/[])
00291 {
00292   MITK_TEST_BEGIN(mitkGeometry3DTest);
00293 
00294   int result;
00295 
00296   MITK_TEST_CONDITION_REQUIRED( (result = testItkImageIsCenterBased()) == EXIT_SUCCESS, "");
00297 
00298   MITK_TEST_OUTPUT(<< "Running main part of test with ImageGeometry = false");
00299   MITK_TEST_CONDITION_REQUIRED( (result = testGeometry3D(false)) == EXIT_SUCCESS, "");
00300 
00301   MITK_TEST_OUTPUT(<< "Running main part of test with ImageGeometry = true");
00302   MITK_TEST_CONDITION_REQUIRED( (result = testGeometry3D(true)) == EXIT_SUCCESS, "");
00303 
00304   MITK_TEST_END();
00305 
00306   return EXIT_SUCCESS;
00307 }