diff --git a/code/core/rttbGeometricInfo.cpp b/code/core/rttbGeometricInfo.cpp index 685812a..6cde4fc 100644 --- a/code/core/rttbGeometricInfo.cpp +++ b/code/core/rttbGeometricInfo.cpp @@ -1,304 +1,304 @@ // ----------------------------------------------------------------------- // RTToolbox - DKFZ radiotherapy quantitative evaluation library // // Copyright (c) German Cancer Research Center (DKFZ), // Software development for Integrated Diagnostics and Therapy (SIDT). // ALL RIGHTS RESERVED. // See rttbCopyright.txt or // http://www.dkfz.de/en/sidt/projects/rttb/copyright.html // // This software is distributed WITHOUT ANY WARRANTY; without even // the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR // PURPOSE. See the above copyright notices for more information. // //------------------------------------------------------------------------ #include "rttbGeometricInfo.h" #include #include #include namespace rttb { namespace core { void GeometricInfo::setSpacing(const SpacingVectorType3D& aSpacingVector) { _spacing = aSpacingVector; } const SpacingVectorType3D& GeometricInfo::getSpacing() const { return _spacing; } void GeometricInfo::setImagePositionPatient(const WorldCoordinate3D& aImagePositionPatient) { _imagePositionPatient = aImagePositionPatient; } const WorldCoordinate3D& GeometricInfo::getImagePositionPatient() const { return _imagePositionPatient; } void GeometricInfo::setOrientationMatrix(const OrientationMatrix& anOrientationMatrix) { _orientationMatrix = anOrientationMatrix; computeInvertOrientation(); } bool GeometricInfo::computeInvertOrientation() { using pmatrix = boost::numeric::ublas::permutation_matrix; boost::numeric::ublas::matrix A(_orientationMatrix); // create a permutation matrix for the LU-factorization pmatrix pm(A.size1()); size_t res = boost::numeric::ublas::lu_factorize(A, pm); if (res != 0) { return false; } _invertedOrientationMatrix.assign(boost::numeric::ublas::identity_matrix (A.size1())); // backsubstitute to get the inverse boost::numeric::ublas::lu_substitute(A, pm, _invertedOrientationMatrix); return true; } void GeometricInfo::setImageSize(const ImageSize& aSize) { setNumColumns(aSize(0)); setNumRows(aSize(1)); setNumSlices(aSize(2)); } const ImageSize GeometricInfo::getImageSize() const { return ImageSize(static_cast(getNumColumns()), static_cast(getNumRows()), static_cast(getNumSlices())); } void GeometricInfo::setNumColumns(const VoxelGridDimensionType aValue) { _numberOfColumns = aValue; } const VoxelGridDimensionType GeometricInfo::getNumColumns() const { return _numberOfColumns; } void GeometricInfo::setNumRows(const VoxelGridDimensionType aValue) { _numberOfRows = aValue; } const VoxelGridDimensionType GeometricInfo::getNumRows() const { return _numberOfRows; } void GeometricInfo::setNumSlices(const VoxelGridDimensionType aValue) { _numberOfFrames = aValue; } const VoxelGridDimensionType GeometricInfo::getNumSlices() const { return _numberOfFrames; } bool operator==(const GeometricInfo& gInfo, const GeometricInfo& gInfo1) { return (gInfo.getImagePositionPatient() == gInfo1.getImagePositionPatient() && gInfo.getOrientationMatrix() == gInfo1.getOrientationMatrix() && gInfo.getSpacing() == gInfo1.getSpacing() && gInfo.getNumColumns() == gInfo1.getNumColumns() && gInfo.getNumRows() == gInfo1.getNumRows() && gInfo.getNumSlices() == gInfo1.getNumSlices()); } bool GeometricInfo::equalsAlmost(const GeometricInfo& another, double errorConstantGI /*= 1e-5*/) const { return (getImagePositionPatient().equalsAlmost(another.getImagePositionPatient(), errorConstantGI) && getOrientationMatrix().equalsAlmost(another.getOrientationMatrix(), errorConstantGI) && getSpacing().equalsAlmost(another.getSpacing(), errorConstantGI) && getNumColumns() == another.getNumColumns() && getNumRows() == another.getNumRows() && getNumSlices() == another.getNumSlices()); } bool GeometricInfo::worldCoordinateToGeometryCoordinate(const WorldCoordinate3D& aWorldCoordinate, ContinousVoxelGridIndex3D& aIndex) const { WorldCoordinate3D distanceToIP; distanceToIP = aWorldCoordinate - _imagePositionPatient; boost::numeric::ublas::vector result = boost::numeric::ublas::prod( _invertedOrientationMatrix, distanceToIP); boost::numeric::ublas::vector resultS = boost::numeric::ublas::element_div(result, _spacing); aIndex = ContinousVoxelGridIndex3D(resultS(0), resultS(1), resultS(2)); //if we convert ContinousVoxelGridIndex3D (double) to VoxelGridIndex3D (unsigned int), we can't find out if it's negative. //So we have to check before. if (aIndex(0) < -0.5 || aIndex(1) < -0.5 || aIndex(2) < -0.5){ return false; } else { //check if it is inside VoxelGridIndex3D indexInt = VoxelGridIndex3D(GridIndexType(aIndex(0) + 0.5), GridIndexType(aIndex(1) + 0.5), GridIndexType(aIndex(2) + 0.5)); return isInside(indexInt); } } bool GeometricInfo::worldCoordinateToIndex(const WorldCoordinate3D& aWorldCoordinate, VoxelGridIndex3D& aIndex) const { ContinousVoxelGridIndex3D doubleIndex; bool inside = worldCoordinateToGeometryCoordinate(aWorldCoordinate, doubleIndex); aIndex = VoxelGridIndex3D(GridIndexType(doubleIndex(0)+0.5), GridIndexType(doubleIndex(1)+0.5), GridIndexType(doubleIndex(2)+0.5)); return inside; } - bool GeometricInfo::geometryCoordinateToWorldCoordinate(const ContinousVoxelGridIndex3D& aIndex, + bool GeometricInfo::continousIndexToWorldCoordinate(const ContinousVoxelGridIndex3D& aIndex, WorldCoordinate3D& aWorldCoordinate) const { boost::numeric::ublas::vector resultS = boost::numeric::ublas::element_prod( aIndex, _spacing); boost::numeric::ublas::vector result = boost::numeric::ublas::prod( _orientationMatrix, resultS); aWorldCoordinate = result + _imagePositionPatient; //if we convert ContinousVoxelGridIndex3D (double) to VoxelGridIndex3D (unsigned int), we can't find out if it's negative. //So we have to check before. if (aIndex(0) < -0.5 || aIndex(1) < -0.5 || aIndex(2) < -0.5){ return false; } else { VoxelGridIndex3D indexInt = VoxelGridIndex3D(GridIndexType(aIndex(0) + 0.5), GridIndexType(aIndex(1) + 0.5), GridIndexType(aIndex(2) + 0.5)); return isInside(indexInt); } } bool GeometricInfo::indexToWorldCoordinate(const VoxelGridIndex3D& aIndex, WorldCoordinate3D& aWorldCoordinate) const { ContinousVoxelGridIndex3D indexDouble = ContinousVoxelGridIndex3D(aIndex(0), aIndex(1), aIndex(2)); - return geometryCoordinateToWorldCoordinate(indexDouble, aWorldCoordinate); + return continousIndexToWorldCoordinate(indexDouble, aWorldCoordinate); } bool GeometricInfo::isInside(const VoxelGridIndex3D& aIndex) const { return (aIndex(0) >= 0 && aIndex(1) >= 0 && aIndex(2) >= 0 && aIndex(0) < static_cast(_numberOfColumns) && aIndex(1) < static_cast(_numberOfRows) && aIndex(2) < static_cast(_numberOfFrames)); } bool GeometricInfo::isInside(const WorldCoordinate3D& aWorldCoordinate) const { VoxelGridIndex3D currentIndex; return (worldCoordinateToIndex(aWorldCoordinate, currentIndex)); } const GridSizeType GeometricInfo::getNumberOfVoxels() const { auto nVoxels = static_cast(_numberOfRows * _numberOfColumns * _numberOfFrames); return nVoxels; } bool GeometricInfo::convert(const VoxelGridID& gridID, VoxelGridIndex3D& gridIndex) const { if (validID(gridID)) { gridIndex(0) = gridID % getNumColumns(); VoxelGridID tempID = (gridID - gridIndex.x()) / getNumColumns(); gridIndex(1) = tempID % getNumRows(); gridIndex(2) = (tempID - gridIndex.y()) / getNumRows(); return true; } return false; } bool GeometricInfo::convert(const VoxelGridIndex3D& gridIndex, VoxelGridID& gridID) const { if ((gridIndex.x() >= static_cast(getNumColumns())) || (gridIndex.y() >= static_cast(getNumRows())) || (gridIndex.z() >= static_cast(getNumSlices()))) { return false; } else { gridID = gridIndex.z() * getNumColumns() * getNumRows() + gridIndex.y() * getNumColumns() + gridIndex.x(); return validID(gridID); } } bool GeometricInfo::validID(const VoxelGridID aID) const { return (aID >= 0 && aID < getNumberOfVoxels()); } bool GeometricInfo::validIndex(const VoxelGridIndex3D& aIndex) const { VoxelGridID aID; if (!convert(aIndex, aID)) { return false; } else { return validID(aID); } } std::ostream& operator<<(std::ostream& s, const GeometricInfo& anGeometricInfo) { s << "[ " << anGeometricInfo.getImagePositionPatient() << "; " << anGeometricInfo.getOrientationMatrix() << "; " << anGeometricInfo.getSpacing() << "; " << "; " << anGeometricInfo.getNumColumns() << "; " << anGeometricInfo.getNumRows() << "; " << anGeometricInfo.getNumSlices() << " ]"; return s; } }//end namespace core }//end namespace rttb diff --git a/code/core/rttbGeometricInfo.h b/code/core/rttbGeometricInfo.h index db45406..e694bd9 100644 --- a/code/core/rttbGeometricInfo.h +++ b/code/core/rttbGeometricInfo.h @@ -1,187 +1,187 @@ // ----------------------------------------------------------------------- // RTToolbox - DKFZ radiotherapy quantitative evaluation library // // Copyright (c) German Cancer Research Center (DKFZ), // Software development for Integrated Diagnostics and Therapy (SIDT). // ALL RIGHTS RESERVED. // See rttbCopyright.txt or // http://www.dkfz.de/en/sidt/projects/rttb/copyright.html // // This software is distributed WITHOUT ANY WARRANTY; without even // the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR // PURPOSE. See the above copyright notices for more information. // //------------------------------------------------------------------------ #ifndef __GEOMETRIC_INFO_NEW_H #define __GEOMETRIC_INFO_NEW_H #include #include #include "rttbBaseType.h" #include #include "rttbCommon.h" #include "RTTBCoreExports.h" #ifdef _MSC_VER #pragma warning(push) #pragma warning(disable: 4251) #endif namespace rttb { namespace core { /*! @brief GeometricInfo objects contain all the information required for transformations between voxel grid coordinates and world coordinates. Corresponding converter functions are also available. @note ITK Pixel Indexing used (http://www.itk.org/Doxygen45/html/classitk_1_1Image.html): The Index type reverses the order so that with Index[0] = col, Index[1] = row, Index[2] = slice. */ class RTTBCore_EXPORT GeometricInfo { public: rttbClassMacroNoParent(GeometricInfo) private: WorldCoordinate3D _imagePositionPatient{ 0 }; OrientationMatrix _orientationMatrix{ 0 }; OrientationMatrix _invertedOrientationMatrix{ 0 }; SpacingVectorType3D _spacing{ 0 }; VoxelGridDimensionType _numberOfColumns{0}; VoxelGridDimensionType _numberOfRows{0}; VoxelGridDimensionType _numberOfFrames{0}; /* @brief Matrix inversion routine. Uses lu_factorize and lu_substitute in uBLAS to invert a matrix http://savingyoutime.wordpress.com/2009/09/21/c-matrix-inversion-boostublas/ */ bool computeInvertOrientation(); public: /*! @brief Constructor, initializes orientation matrix, spacing vector and patient position with zeros. */ GeometricInfo() = default; void setSpacing(const SpacingVectorType3D& aSpacingVector); const SpacingVectorType3D& getSpacing() const; void setImagePositionPatient(const WorldCoordinate3D& aImagePositionPatient); const WorldCoordinate3D& getImagePositionPatient() const; void setOrientationMatrix(const OrientationMatrix& anOrientationMatrix); const OrientationMatrix getOrientationMatrix() const { return _orientationMatrix; }; void setImageSize(const ImageSize& aSize); const ImageSize getImageSize() const; void setNumColumns(const VoxelGridDimensionType aValue); const VoxelGridDimensionType getNumColumns() const; void setNumRows(const VoxelGridDimensionType aValue); const VoxelGridDimensionType getNumRows() const; void setNumSlices(const VoxelGridDimensionType aValue); const VoxelGridDimensionType getNumSlices() const; /*! @brief determines equality of two GeometricInfo objects. */ friend bool RTTBCore_EXPORT operator == (const GeometricInfo& gInfo, const GeometricInfo& gInfo1); bool equalsAlmost(const GeometricInfo& another, double errorConstantGI = 1e-5) const; /*! @brief converts world coordinates to voxel grid index. @details the voxels world coordinates are defined by spacing, orientation and imagePositionPatient. (-0.5/-0.5/-0.5) --> (0/0/0) and (0.4999/0.4999/0.4999) --> (0/0/0) define the outer coordinates of a voxel with spacing=1, orientation= x y z (identity matrix) and imagePositionPatient=(0/0/0). @sa WorldCoordinate3D VoxelGridIndex3D @note The conversion of values is done even if the target index is not inside the given voxel grid. @returns false if aWorldCoordinate is outside the voxel grid, true otherwise. */ bool worldCoordinateToIndex(const WorldCoordinate3D& aWorldCoordinate, VoxelGridIndex3D& aIndex) const; /*! @brief converts world coordinates to double geometry coordinate. @details This is needed because of a double precision voxel coordinate system for voxelization. The world coordinate of the image position patient is the center of the first voxel (0.0/0.0/0.0). (-0.5/-0.5/-0.5) --> (-0.5/-0.5/-0.5) and (0.4999/0.4999/0.4999) --> (0.4999/0.4999/0.4999) with spacing=1, orientation= x y z (identity matrix) and imagePositionPatient=(0/0/0). @sa WorldCoordinate3D, ContinousVoxelGridIndex3D @note The conversion of values is done even if the target index is not inside the given voxel grid. @returns false if aWorldCoordinate is outside the voxel grid, true otherwise. */ bool worldCoordinateToGeometryCoordinate(const WorldCoordinate3D& aWorldCoordinate, ContinousVoxelGridIndex3D& aIndex) const; /*! @brief converts double geometry coordinate to world coordinates. @details This is needed because of a double precision voxel coordinate system for voxelization. The world coordinate of the image position patient is the center of the first voxel (0.0/0.0/0.0). (-0.5/-0.5/-0.5) --> (-0.5/-0.5/-0.5) and (5.5/3.2/1.0) --> (5.5/3.2/1.0) with spacing=1, orientation= x y z (identity matrix) and imagePositionPatient=(0/0/0). @sa ContinousVoxelGridIndex3D, WorldCoordinate3D @note The conversion of values is done even if the target index is not inside the given voxel grid. @returns false if aWorldCoordinate is outside the voxel grid, true otherwise. */ - bool geometryCoordinateToWorldCoordinate(const ContinousVoxelGridIndex3D& aIndex, + bool continousIndexToWorldCoordinate(const ContinousVoxelGridIndex3D& aIndex, WorldCoordinate3D& aWorldCoordinate) const; /*! @brief convert voxel grid index to world coordinates @details The world coordinate of the image position patient (center of the first voxel) is the center of the first voxel (0.0/0.0/0.0) (0/0/0) --> (0.0/0.0/0.0) and (1/1/2) --> (1.0/1.0/2.0) with spacing=1, orientation= x y z (identity matrix) and imagePositionPatient=(0/0/0). Thus, the center of the voxel is taken and converted. @sa VoxelGridIndex3D, WorldCoordinate3D @note The conversion of values is done even if the target index is not inside the given voxel grid. @returns false if aWorldCoordinate is outside the voxel grid, true otherwise. */ bool indexToWorldCoordinate(const VoxelGridIndex3D& aIndex, WorldCoordinate3D& aWorldCoordinate) const; /*! @brief check if a given voxel grid index is inside the given voxel grid.*/ bool isInside(const VoxelGridIndex3D& aIndex) const; /*! @brief check if a given world coordinate is inside the given voxel grid.*/ bool isInside(const WorldCoordinate3D& aWorldCoordinate) const; const GridSizeType getNumberOfVoxels() const; bool convert(const VoxelGridID& gridID, VoxelGridIndex3D& gridIndex) const; bool convert(const VoxelGridIndex3D& gridIndex, VoxelGridID& gridID) const; /*! @brief test if given ID is inside current dose grid */ bool validID(const VoxelGridID aID) const; /*! @brief test if given index is inside current dose grid */ bool validIndex(const VoxelGridIndex3D& aIndex) const; /*!@ brief generates string stream representation of the GeometricInfo object. */ friend std::ostream& operator << (std::ostream& s, const GeometricInfo& anGeometricInfo); }; } } #ifdef _MSC_VER #pragma warning(pop) #endif #endif diff --git a/testing/core/CreateTestStructure.cpp b/testing/core/CreateTestStructure.cpp index 110dee9..11d524e 100644 --- a/testing/core/CreateTestStructure.cpp +++ b/testing/core/CreateTestStructure.cpp @@ -1,681 +1,681 @@ // ----------------------------------------------------------------------- // RTToolbox - DKFZ radiotherapy quantitative evaluation library // // Copyright (c) German Cancer Research Center (DKFZ), // Software development for Integrated Diagnostics and Therapy (SIDT). // ALL RIGHTS RESERVED. // See rttbCopyright.txt or // http://www.dkfz.de/en/sidt/projects/rttb/copyright.html // // This software is distributed WITHOUT ANY WARRANTY; without even // the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR // PURPOSE. See the above copyright notices for more information. // //------------------------------------------------------------------------ // this file defines the rttbCoreTests for the test driver // and all it expects is that you have a function called RegisterTests #include #include "CreateTestStructure.h" #include "rttbNullPointerException.h" #include "rttbInvalidParameterException.h" #include "rttbInvalidDoseException.h" namespace rttb { namespace testing { CreateTestStructure::~CreateTestStructure() {} CreateTestStructure::CreateTestStructure(const core::GeometricInfo& aGeoInfo) { _geoInfo = aGeoInfo; } PolygonType CreateTestStructure::createPolygonLeftUpper(std::vector aVoxelVector, GridIndexType sliceNumber) { PolygonType polygon; for (size_t i = 0; i < aVoxelVector.size(); i++) { VoxelGridIndex3D voxelIndex; voxelIndex(0) = (aVoxelVector.at(i)).x(); voxelIndex(1) = (aVoxelVector.at(i)).y(); voxelIndex(2) = sliceNumber; WorldCoordinate3D p1; _geoInfo.indexToWorldCoordinate(voxelIndex, p1); polygon.push_back(p1); std::cout << "(" << p1.x() << "," << p1.y() << "," << p1.z() << ")" << "; "; } std::cout << std::endl; return polygon; } PolygonType CreateTestStructure::createPolygonCenter(std::vector aVoxelVector, GridIndexType sliceNumber) { PolygonType polygon; for (size_t i = 0; i < aVoxelVector.size(); i++) { VoxelGridIndex3D voxelIndex; voxelIndex(0) = (aVoxelVector.at(i)).x(); voxelIndex(1) = (aVoxelVector.at(i)).y(); voxelIndex(2) = sliceNumber; WorldCoordinate3D p1; _geoInfo.indexToWorldCoordinate(voxelIndex, p1); WorldCoordinate3D p; p(0) = p1.x() ; p(1) = p1.y() ; //This can be done directly for x/y coordinates, but not for z. Thus, we do it in this function. p(2) = p1.z() - _geoInfo.getSpacing().z() / 2; polygon.push_back(p); } return polygon; } PolygonType CreateTestStructure::createPolygonCenterOnPlaneCenter(std::vector aVoxelVector, GridIndexType sliceNumber) { PolygonType polygon; for (size_t i = 0; i < aVoxelVector.size(); i++) { VoxelGridIndex3D voxelIndex; ContinousVoxelGridIndex3D indexDouble = ContinousVoxelGridIndex3D((aVoxelVector.at(i)).x(), (aVoxelVector.at(i)).y(), sliceNumber); WorldCoordinate3D p1; - _geoInfo.geometryCoordinateToWorldCoordinate(indexDouble, p1); + _geoInfo.continousIndexToWorldCoordinate(indexDouble, p1); polygon.push_back(p1); } return polygon; } PolygonType CreateTestStructure::createPolygonBetweenUpperLeftAndCenter( std::vector aVoxelVector, GridIndexType sliceNumber) { PolygonType polygon; for (size_t i = 0; i < aVoxelVector.size(); i++) { VoxelGridIndex3D voxelIndex; voxelIndex(0) = (aVoxelVector.at(i)).x(); voxelIndex(1) = (aVoxelVector.at(i)).y(); voxelIndex(2) = sliceNumber; WorldCoordinate3D p1; _geoInfo.indexToWorldCoordinate(voxelIndex, p1); SpacingVectorType3D delta = _geoInfo.getSpacing(); WorldCoordinate3D p; p(0) = p1.x() + delta.x() / 4; p(1) = p1.y() + delta.y() / 4; p(2) = p1.z(); polygon.push_back(p); std::cout << "(" << p.x() << "," << p.y() << "," << p.z() << ")" << "; "; } std::cout << std::endl; return polygon; } PolygonType CreateTestStructure::createPolygonBetweenLowerRightAndCenter( std::vector aVoxelVector, GridIndexType sliceNumber) { PolygonType polygon; for (size_t i = 0; i < aVoxelVector.size(); i++) { VoxelGridIndex3D voxelIndex; voxelIndex(0) = (aVoxelVector.at(i)).x(); voxelIndex(1) = (aVoxelVector.at(i)).y(); voxelIndex(2) = sliceNumber; WorldCoordinate3D p1; _geoInfo.indexToWorldCoordinate(voxelIndex, p1); SpacingVectorType3D delta = _geoInfo.getSpacing(); WorldCoordinate3D p; p(0) = p1.x() + delta.x() * 0.75; p(1) = p1.y() + delta.y() * 0.75; p(2) = p1.z(); polygon.push_back(p); std::cout << "(" << p.x() << "," << p.y() << "," << p.z() << ")" << "; "; } std::cout << std::endl; return polygon; } PolygonType CreateTestStructure::createPolygonLeftEdgeMiddle(std::vector aVoxelVector, GridIndexType sliceNumber) { PolygonType polygon; for (size_t i = 0; i < aVoxelVector.size(); i++) { VoxelGridIndex3D voxelIndex; voxelIndex(0) = (aVoxelVector.at(i)).x(); voxelIndex(1) = (aVoxelVector.at(i)).y(); voxelIndex(2) = sliceNumber; WorldCoordinate3D p1; _geoInfo.indexToWorldCoordinate(voxelIndex, p1); SpacingVectorType3D delta = _geoInfo.getSpacing(); WorldCoordinate3D p; p(0) = p1.x(); p(1) = p1.y() + delta.y() * 0.5; p(2) = p1.z(); polygon.push_back(p); std::cout << "(" << p.x() << "," << p.y() << "," << p.z() << ")" << "; "; } std::cout << std::endl; return polygon; } PolygonType CreateTestStructure::createPolygonUpperCenter(std::vector aVoxelVector, GridIndexType sliceNumber) { PolygonType polygon; for (size_t i = 0; i < aVoxelVector.size(); i++) { VoxelGridIndex3D voxelIndex; voxelIndex(0) = (aVoxelVector.at(i)).x(); voxelIndex(1) = (aVoxelVector.at(i)).y(); voxelIndex(2) = sliceNumber; WorldCoordinate3D p1; _geoInfo.indexToWorldCoordinate(voxelIndex, p1); SpacingVectorType3D delta = _geoInfo.getSpacing(); WorldCoordinate3D p; p(0) = p1.x() + delta.x() * 0.5; p(1) = p1.y(); p(2) = p1.z(); polygon.push_back(p); std::cout << "(" << p.x() << "," << p.y() << "," << p.z() << ")" << "; "; } std::cout << std::endl; return polygon; } PolygonType CreateTestStructure::createPolygonIntermediatePoints(std::vector aVoxelVector, GridIndexType sliceNumber) { PolygonType polygon; for (size_t i = 0; i < aVoxelVector.size(); i++) { VoxelGridIndex3D voxelIndex; VoxelGridIndex3D voxelIndex2; voxelIndex(0) = (aVoxelVector.at(i)).x(); voxelIndex(1) = (aVoxelVector.at(i)).y(); voxelIndex(2) = sliceNumber; if (i < (aVoxelVector.size() - 1)) { voxelIndex2(0) = (aVoxelVector.at(i + 1)).x(); voxelIndex2(1) = (aVoxelVector.at(i + 1)).y(); voxelIndex2(2) = sliceNumber; } else { voxelIndex2(0) = (aVoxelVector.at(0)).x(); voxelIndex2(1) = (aVoxelVector.at(0)).y(); voxelIndex2(2) = sliceNumber; } WorldCoordinate3D p1; _geoInfo.indexToWorldCoordinate(voxelIndex, p1); SpacingVectorType3D delta = _geoInfo.getSpacing(); WorldCoordinate3D p2; _geoInfo.indexToWorldCoordinate(voxelIndex2, p2); SpacingVectorType3D delta2 = _geoInfo.getSpacing(); WorldCoordinate3D wp1; wp1(0) = p1.x() + delta.x() * 0.75; wp1(1) = p1.y() + delta.y() * 0.75; wp1(2) = p1.z(); WorldCoordinate3D wp2; wp2(0) = p2.x() + delta.x() * 0.75; wp2(1) = p2.y() + delta.y() * 0.75; wp2(2) = p2.z(); polygon.push_back(wp1); double diffX = (wp2.x() - wp1.x()) / 1000.0; double diffY = (wp2.y() - wp1.y()) / 1000.0; WorldCoordinate3D wp_intermediate; wp_intermediate(0) = 0; wp_intermediate(1) = 0; for (int j = 0; j < 1000; j++) { wp_intermediate(0) = wp1.x() + j * diffX; wp_intermediate(1) = wp1.y() + j * diffY; polygon.push_back(wp_intermediate); } } std::cout << std::endl; return polygon; } PolygonType CreateTestStructure::createPolygonCircle(std::vector aVoxelVector, GridIndexType sliceNumber) { PolygonType polygon; if (aVoxelVector.size() > 0) { unsigned int i = 0; VoxelGridIndex3D voxelIndex; voxelIndex(0) = (aVoxelVector.at(i)).x(); voxelIndex(1) = (aVoxelVector.at(i)).y(); voxelIndex(2) = sliceNumber; WorldCoordinate3D p1; _geoInfo.indexToWorldCoordinate(voxelIndex, p1); SpacingVectorType3D delta = _geoInfo.getSpacing(); WorldCoordinate3D wp1; wp1(0) = p1.x(); wp1(1) = p1.y(); wp1(2) = p1.z(); WorldCoordinate3D wp_intermediate; wp_intermediate(0) = 0; wp_intermediate(1) = 0; wp_intermediate(2) = p1.z(); double radius = 2 * delta.x(); double frac_radius = (radius * 0.001); double correct_y = (delta.x() / delta.y()); for (unsigned int j = 0; j <= 1000; j++) { double y_offset = sqrt((radius * radius) - ((frac_radius * j) * (frac_radius * j))); wp_intermediate(0) = wp1.x() + frac_radius * j; wp_intermediate(1) = wp1.y() - (y_offset * correct_y); polygon.push_back(wp_intermediate); } for (unsigned int j = 1000; j <= 1000; j--) { double y_offset = sqrt((radius * radius) - ((frac_radius * j) * (frac_radius * j))); wp_intermediate(0) = wp1.x() + frac_radius * j; wp_intermediate(1) = wp1.y() + (y_offset * correct_y); polygon.push_back(wp_intermediate); } for (unsigned int j = 0; j <= 1000; j++) { double y_offset = sqrt((radius * radius) - ((frac_radius * j) * (frac_radius * j))); wp_intermediate(0) = wp1.x() - frac_radius * j; wp_intermediate(1) = wp1.y() + y_offset * correct_y; polygon.push_back(wp_intermediate); } for (unsigned int j = 1000; j <= 1000; j--) { double y_offset = sqrt((radius * radius) - ((frac_radius * j) * (frac_radius * j))); wp_intermediate(0) = wp1.x() + frac_radius * j; wp_intermediate(1) = wp1.y() + (y_offset * correct_y); polygon.push_back(wp_intermediate); } } std::cout << std::endl; return polygon; } PolygonType CreateTestStructure::createStructureSeveralSectionsInsideOneVoxelA( std::vector aVoxelVector, GridIndexType sliceNumber) { PolygonType polygon; if (aVoxelVector.size() > 0) { int i = 0; VoxelGridIndex3D voxelIndex; voxelIndex(0) = (aVoxelVector.at(i)).x(); voxelIndex(1) = (aVoxelVector.at(i)).y(); voxelIndex(2) = sliceNumber; WorldCoordinate3D p1; _geoInfo.indexToWorldCoordinate(voxelIndex, p1); SpacingVectorType3D delta = _geoInfo.getSpacing(); WorldCoordinate3D wp1; wp1(0) = p1.x(); wp1(1) = p1.y(); wp1(2) = p1.z(); WorldCoordinate3D wp_intermediate; wp_intermediate(0) = 0; wp_intermediate(1) = 0; wp_intermediate(2) = p1.z(); wp_intermediate(0) = wp1.x() + (delta.x() * 0.25); wp_intermediate(1) = wp1.y() + (delta.y() * 0.75); polygon.push_back(wp_intermediate); wp_intermediate(0) = wp1.x() + (delta.x() * 0.25); wp_intermediate(1) = wp1.y() + (delta.y() * 2.75); polygon.push_back(wp_intermediate); wp_intermediate(0) = wp1.x() + (delta.x() * 0.5); wp_intermediate(1) = wp1.y() + (delta.y() * 2.75); polygon.push_back(wp_intermediate); wp_intermediate(0) = wp1.x() + (delta.x() * 0.5); wp_intermediate(1) = wp1.y() + (delta.y() * 0.75); polygon.push_back(wp_intermediate); wp_intermediate(0) = wp1.x() + (delta.x() * 0.75); wp_intermediate(1) = wp1.y() + (delta.y() * 0.75); polygon.push_back(wp_intermediate); wp_intermediate(0) = wp1.x() + (delta.x() * 0.75); wp_intermediate(1) = wp1.y() + (delta.y() * 2.75); polygon.push_back(wp_intermediate); wp_intermediate(0) = wp1.x() + (delta.x() * 1.0); wp_intermediate(1) = wp1.y() + (delta.y() * 2.75); polygon.push_back(wp_intermediate); wp_intermediate(0) = wp1.x() + (delta.x() * 1.0); wp_intermediate(1) = wp1.y() + (delta.y() * 0.75); polygon.push_back(wp_intermediate); wp_intermediate(0) = wp1.x() + (delta.x() * 1.25); wp_intermediate(1) = wp1.y() + (delta.y() * 0.75); polygon.push_back(wp_intermediate); wp_intermediate(0) = wp1.x() + (delta.x() * 1.25); wp_intermediate(1) = wp1.y() + (delta.y() * 2.75); polygon.push_back(wp_intermediate); wp_intermediate(0) = wp1.x() + (delta.x() * 1.5); wp_intermediate(1) = wp1.y() + (delta.y() * 2.75); polygon.push_back(wp_intermediate); wp_intermediate(0) = wp1.x() + (delta.x() * 1.5); wp_intermediate(1) = wp1.y() + (delta.y() * 0.75); polygon.push_back(wp_intermediate); wp_intermediate(0) = wp1.x() + (delta.x() * 1.75); wp_intermediate(1) = wp1.y() + (delta.y() * 0.75); polygon.push_back(wp_intermediate); wp_intermediate(0) = wp1.x() + (delta.x() * 1.75); wp_intermediate(1) = wp1.y() + (delta.y() * 3.0); polygon.push_back(wp_intermediate); wp_intermediate(0) = wp1.x(); wp_intermediate(1) = wp1.y() + (delta.y() * 3.0); polygon.push_back(wp_intermediate); wp_intermediate(0) = wp1.x(); wp_intermediate(1) = wp1.y() + (delta.y() * 3.0); polygon.push_back(wp_intermediate); wp_intermediate(0) = wp1.x(); wp_intermediate(1) = wp1.y() + (delta.y() * 0.75); polygon.push_back(wp_intermediate); } std::cout << std::endl; return polygon; } PolygonType CreateTestStructure::createStructureSelfTouchingA(std::vector aVoxelVector, GridIndexType sliceNumber) { PolygonType polygon; if (aVoxelVector.size() > 0) { int i = 0; VoxelGridIndex3D voxelIndex; voxelIndex(0) = (aVoxelVector.at(i)).x(); voxelIndex(1) = (aVoxelVector.at(i)).y(); voxelIndex(2) = sliceNumber; WorldCoordinate3D p1; _geoInfo.indexToWorldCoordinate(voxelIndex, p1); SpacingVectorType3D delta = _geoInfo.getSpacing(); WorldCoordinate3D wp1; wp1(0) = p1.x(); wp1(1) = p1.y(); wp1(2) = p1.z(); WorldCoordinate3D wp_intermediate; wp_intermediate(0) = 0; wp_intermediate(1) = 0; wp_intermediate(2) = p1.z(); wp_intermediate(0) = wp1.x(); wp_intermediate(1) = wp1.y() + (delta.y() * 0.5); polygon.push_back(wp_intermediate); wp_intermediate(0) = wp1.x() + (delta.x() * 0.5); wp_intermediate(1) = wp1.y(); polygon.push_back(wp_intermediate); wp_intermediate(0) = wp1.x() + (delta.x() * 0.5); wp_intermediate(1) = wp1.y() + (delta.y() * 0.5); polygon.push_back(wp_intermediate); wp_intermediate(0) = wp1.x() + (delta.x() * 1.0); wp_intermediate(1) = wp1.y() + (delta.y() * 0.5); polygon.push_back(wp_intermediate); wp_intermediate(0) = wp1.x() + (delta.x() * 1.0); wp_intermediate(1) = wp1.y() + (delta.y() * 1.0); polygon.push_back(wp_intermediate); wp_intermediate(0) = wp1.x() + (delta.x() * 0.5); wp_intermediate(1) = wp1.y() + (delta.y() * 1.0); polygon.push_back(wp_intermediate); wp_intermediate(0) = wp1.x() + (delta.x() * 0.5); wp_intermediate(1) = wp1.y() + (delta.y() * 0.5); polygon.push_back(wp_intermediate); wp_intermediate(0) = wp1.x() + (delta.x() * 0.5); wp_intermediate(1) = wp1.y() + (delta.y() * 1.0); polygon.push_back(wp_intermediate); wp_intermediate(0) = wp1.x(); wp_intermediate(1) = wp1.y() + (delta.y() * 0.5); polygon.push_back(wp_intermediate); } std::cout << std::endl; return polygon; } PolygonType CreateTestStructure::createStructureSelfTouchingB(std::vector aVoxelVector, GridIndexType sliceNumber) { PolygonType polygon; if (aVoxelVector.size() > 0) { int i = 0; VoxelGridIndex3D voxelIndex; voxelIndex(0) = (aVoxelVector.at(i)).x(); voxelIndex(1) = (aVoxelVector.at(i)).y(); voxelIndex(2) = sliceNumber; WorldCoordinate3D p1; _geoInfo.indexToWorldCoordinate(voxelIndex, p1); SpacingVectorType3D delta = _geoInfo.getSpacing(); WorldCoordinate3D wp1; wp1(0) = p1.x(); wp1(1) = p1.y(); wp1(2) = p1.z(); WorldCoordinate3D wp_intermediate; wp_intermediate(0) = 0; wp_intermediate(1) = 0; wp_intermediate(2) = p1.z(); wp_intermediate(0) = wp1.x(); wp_intermediate(1) = wp1.y() + (delta.y() * 0.5); polygon.push_back(wp_intermediate); wp_intermediate(0) = wp1.x() + (delta.x() * 0.5); wp_intermediate(1) = wp1.y(); polygon.push_back(wp_intermediate); wp_intermediate(0) = wp1.x() + (delta.x() * 0.5); wp_intermediate(1) = wp1.y() + (delta.y() * 0.5); polygon.push_back(wp_intermediate); wp_intermediate(0) = wp1.x() + (delta.x() * 1.0); wp_intermediate(1) = wp1.y() + (delta.y() * 0.5); polygon.push_back(wp_intermediate); wp_intermediate(0) = wp1.x() + (delta.x() * 1.0); wp_intermediate(1) = wp1.y() + (delta.y() * 1.0); polygon.push_back(wp_intermediate); wp_intermediate(0) = wp1.x() + (delta.x() * 0.5); wp_intermediate(1) = wp1.y() + (delta.y() * 1.0); polygon.push_back(wp_intermediate); wp_intermediate(0) = wp1.x() + (delta.x() * 0.5); wp_intermediate(1) = wp1.y() + (delta.y() * 0.5); polygon.push_back(wp_intermediate); wp_intermediate(0) = wp1.x(); wp_intermediate(1) = wp1.y() + (delta.y() * 0.5); polygon.push_back(wp_intermediate); wp_intermediate(0) = wp1.x() + (delta.x() * 0.5); wp_intermediate(1) = wp1.y() + (delta.y() * 0.5); polygon.push_back(wp_intermediate); wp_intermediate(0) = wp1.x() + (delta.x() * 0.5); wp_intermediate(1) = wp1.y() + (delta.y() * 1.0); polygon.push_back(wp_intermediate); wp_intermediate(0) = wp1.x(); wp_intermediate(1) = wp1.y() + (delta.y() * 0.5); polygon.push_back(wp_intermediate); } std::cout << std::endl; return polygon; } }//testing }//rttb diff --git a/testing/core/GeometricInfoTest.cpp b/testing/core/GeometricInfoTest.cpp index 6b51e8d..df99308 100644 --- a/testing/core/GeometricInfoTest.cpp +++ b/testing/core/GeometricInfoTest.cpp @@ -1,555 +1,555 @@ // ----------------------------------------------------------------------- // RTToolbox - DKFZ radiotherapy quantitative evaluation library // // Copyright (c) German Cancer Research Center (DKFZ), // Software development for Integrated Diagnostics and Therapy (SIDT). // ALL RIGHTS RESERVED. // See rttbCopyright.txt or // http://www.dkfz.de/en/sidt/projects/rttb/copyright.html // // This software is distributed WITHOUT ANY WARRANTY; without even // the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR // PURPOSE. See the above copyright notices for more information. // //------------------------------------------------------------------------ #include "litCheckMacros.h" #include "rttbBaseType.h" #include "rttbGeometricInfo.h" namespace rttb { namespace testing { /*!@brief GeometricInfoTest - test the API of GeometricInfo @note ITK pixel indexing: Index[0] = col, Index[1] = row, Index[2] = slice. 1) test default constructor (values as expected?) 2) test set/getImagePositionPatient 4) test set/getSpacing 5) test set/getNumColumns/Rows/Slices 6) test get/setOrientationMatrix 8) test operators "==" 9) test equalsAlmost 10) test world to index coordinate conversion 11) test isInside and index to world coordinate conversion - 12) test with simple Geometry: isInside, geometryCoordinateToWorldCoordinate(), worldCoordinateToGeometryCoordinate(), indexToWorldCoordinate() + 12) test with simple Geometry: isInside, continousIndexToWorldCoordinate(), worldCoordinateToGeometryCoordinate(), indexToWorldCoordinate() 13) test getNumberOfVoxels 14) Test convert, validID and validIndex */ int GeometricInfoTest(int argc, char* argv[]) { PREPARE_DEFAULT_TEST_REPORTING; //1) test default constructor (values as expected?) CHECK_NO_THROW(core::GeometricInfo()); core::GeometricInfo geoInfo; SpacingVectorType3D testNullSV(0); WorldCoordinate3D testNullWC(0); OrientationMatrix testNullOM(0); CHECK_EQUAL(testNullSV, geoInfo.getSpacing()); CHECK_EQUAL(testNullWC, geoInfo.getImagePositionPatient()); CHECK_EQUAL(testNullOM, geoInfo.getOrientationMatrix()); //2) test set/getImagePositionPatient WorldCoordinate3D testIPP(1.2, 3.4, 5.6); CHECK_NO_THROW(geoInfo.setImagePositionPatient(testIPP)); geoInfo.setImagePositionPatient(testIPP); CHECK_EQUAL(testIPP, geoInfo.getImagePositionPatient()); //4) test set/getSpacing //negative spacing does not make sense! /*!@is related to #2028 Should SpacingTypeVector/GridVolumeType/OrientationMatrix be forced to be non-negative?*/ SpacingVectorType3D expectedSpacing(4.15, 2.35, 100); expectedSpacing(0) = 4.15; expectedSpacing(1) = 2.35; expectedSpacing(2) = 100; CHECK_NO_THROW(geoInfo.setSpacing(expectedSpacing)); geoInfo.setSpacing(expectedSpacing); CHECK_EQUAL(expectedSpacing, geoInfo.getSpacing()); //5) test set/getNumColumns/Rows/Slices const VoxelGridIndex3D expectedVoxelDims(10, 5, 3); //CHECK_THROW(geoInfo.setNumColumns(1.2)); -> implicit conversion will prevent exception CHECK_NO_THROW(geoInfo.setNumColumns(expectedVoxelDims(0))); geoInfo.setNumColumns(expectedVoxelDims(0)); CHECK_NO_THROW(geoInfo.setNumRows(expectedVoxelDims(1))); geoInfo.setNumRows(expectedVoxelDims(1)); //CHECK_THROW(geoInfo.setNumSlices(4.2)); -> implicit conversion will prevent exception CHECK_NO_THROW(geoInfo.setNumSlices(expectedVoxelDims(2))); geoInfo.setNumSlices(expectedVoxelDims(2)); ImageSize rttbSize = geoInfo.getImageSize(); CHECK_EQUAL(rttbSize(0), geoInfo.getNumColumns()); CHECK_EQUAL(rttbSize(1), geoInfo.getNumRows()); CHECK_EQUAL(rttbSize(2), geoInfo.getNumSlices()); rttbSize = ImageSize(11, 99, 6); core::GeometricInfo geoInfo3; geoInfo3.setImageSize(rttbSize); CHECK_EQUAL(rttbSize(0), geoInfo3.getNumColumns()); CHECK_EQUAL(rttbSize(1), geoInfo3.getNumRows()); CHECK_EQUAL(rttbSize(2), geoInfo3.getNumSlices()); //6) test get/setOrientationMatrix CHECK_EQUAL(testNullOM, geoInfo.getOrientationMatrix()); OrientationMatrix testOM(0); const WorldCoordinate3D testIORow(5.5, 4.7, 3.2); const WorldCoordinate3D testIOColumn(2.5, 1.8, 9.1); WorldCoordinate3D ortho = testIORow.cross(testIOColumn); testOM(0, 0) = testIORow(0); testOM(1, 0) = testIORow(1); testOM(2, 0) = testIORow(2); CHECK_NO_THROW(geoInfo.setOrientationMatrix(testOM)); geoInfo.setOrientationMatrix(testOM); CHECK_EQUAL(testOM, geoInfo.getOrientationMatrix()); testOM(0, 1) = testIOColumn(0); testOM(1, 1) = testIOColumn(1); testOM(2, 1) = testIOColumn(2); CHECK_NO_THROW(geoInfo.setOrientationMatrix(testOM)); geoInfo.setOrientationMatrix(testOM); CHECK_EQUAL(testOM, geoInfo.getOrientationMatrix()); testOM(0, 2) = ortho(0); testOM(1, 2) = ortho(1); testOM(2, 2) = ortho(2); CHECK_NO_THROW(geoInfo.setOrientationMatrix(testOM)); geoInfo.setOrientationMatrix(testOM); CHECK_EQUAL(testOM, geoInfo.getOrientationMatrix()); //8) test operators "==" core::GeometricInfo geoInfo2; CHECK_EQUAL(geoInfo, geoInfo); CHECK(!(geoInfo == geoInfo2)); CHECK_EQUAL(geoInfo.getOrientationMatrix(), testOM); CHECK(!(geoInfo.getOrientationMatrix() == testNullOM)); //9) test equalsALmost OrientationMatrix testOM2 = testOM; SpacingVectorType3D testSPV2 = expectedSpacing; WorldCoordinate3D testIPP2 = testIPP; core::GeometricInfo testGI2, testGIEmpty; testGI2.setImagePositionPatient(testIPP2); testGI2.setOrientationMatrix(testOM2); testGI2.setSpacing(testSPV2); double smallValue = 0.000000001; testOM(0, 0) += smallValue; testSPV2(2) += smallValue; testIPP2(1) += smallValue; core::GeometricInfo testGI2similar; testGI2similar.setImagePositionPatient(testIPP2); testGI2similar.setOrientationMatrix(testOM2); testGI2similar.setSpacing(testSPV2); CHECK_EQUAL(testGI2.equalsAlmost(testGI2similar), true); CHECK_EQUAL(testGI2similar.equalsAlmost(testGI2), true); CHECK_EQUAL(testGI2.equalsAlmost(testGI2similar, smallValue * 0.001), false); CHECK_EQUAL(testGIEmpty.equalsAlmost(testGI2), false); CHECK_EQUAL(testGI2.equalsAlmost(testGIEmpty), false); //10) test world to index coordinate conversion //use unit matrix as orientation matrix CHECK_NO_THROW(geoInfo.setOrientationMatrix(OrientationMatrix())); //origin (inside) WorldCoordinate3D insideTestWC1 = geoInfo.getImagePositionPatient(); //inside const VoxelGridIndex3D expectedIndex(8, 3, 2); WorldCoordinate3D insideTestWC2(expectedIndex(0)*expectedSpacing(0) + testIPP(0), expectedIndex(1)*expectedSpacing(1) + testIPP(1), expectedIndex(2)*expectedSpacing(2) + testIPP(2)); //outside WorldCoordinate3D insideTestWC3(-33.12, 0, 14); // outside (dimension of grid) WorldCoordinate3D insideTestWC4(expectedVoxelDims(0)*expectedSpacing(0) + testIPP(0), expectedVoxelDims(1)*expectedSpacing(1) + testIPP(1), expectedVoxelDims(2)*expectedSpacing(2) + testIPP(2)); CHECK(geoInfo.isInside(insideTestWC1)); CHECK(geoInfo.isInside(insideTestWC2)); CHECK(!(geoInfo.isInside(insideTestWC3))); CHECK(!(geoInfo.isInside(insideTestWC4))); VoxelGridIndex3D testConvert(0); CHECK(geoInfo.worldCoordinateToIndex(insideTestWC1, testConvert)); CHECK(geoInfo.isInside(testConvert)); CHECK_EQUAL(VoxelGridIndex3D(0), testConvert); CHECK(geoInfo.worldCoordinateToIndex(insideTestWC2, testConvert)); CHECK(geoInfo.isInside(testConvert)); CHECK_EQUAL(expectedIndex, testConvert); CHECK(!(geoInfo.worldCoordinateToIndex(insideTestWC3, testConvert))); //CHECK_EQUAL(VoxelGridIndex3D(0),testConvert); //if value is in a negative grid position it will be converted //to a very large unrelated number. CHECK(!(geoInfo.isInside(testConvert))); CHECK(!(geoInfo.worldCoordinateToIndex(insideTestWC4, testConvert))); CHECK_EQUAL(expectedVoxelDims, testConvert); CHECK(!(geoInfo.isInside(testConvert))); //use a more complicated orientation matrix OrientationMatrix newOrientation(0); newOrientation(0, 0) = 0.5; newOrientation(1, 2) = -3; newOrientation(2, 1) = 1; CHECK_NO_THROW(geoInfo.setOrientationMatrix(newOrientation)); testIPP = WorldCoordinate3D(20, 100, -1000); CHECK_NO_THROW(geoInfo.setImagePositionPatient(testIPP)); CHECK_NO_THROW(geoInfo.setSpacing(SpacingVectorType3D(1))); //values for testing were generated with a dedicated MeVisLab routine insideTestWC1 = geoInfo.getImagePositionPatient(); //origin (inside) const VoxelGridIndex3D expectedIndexWC1(0, 0, 0); insideTestWC2 = WorldCoordinate3D(22.5, 97, -998); //inside const VoxelGridIndex3D expectedIndexWC2(5, 2, 1); insideTestWC3 = WorldCoordinate3D(26, 88, -996); //outside const VoxelGridIndex3D expectedIndexWC3(12, 4, 4); insideTestWC4 = WorldCoordinate3D(25, 91, -995); // outside: Grid dimension = [10,5,3] const VoxelGridIndex3D expectedIndexWC4(10, 5, 3); CHECK(geoInfo.isInside(insideTestWC1)); CHECK_EQUAL(geoInfo.isInside(insideTestWC1), geoInfo.isInside(expectedIndexWC1)); CHECK(geoInfo.isInside(insideTestWC2)); CHECK_EQUAL(geoInfo.isInside(insideTestWC2), geoInfo.isInside(expectedIndexWC2)); CHECK(!(geoInfo.isInside(insideTestWC3))); CHECK_EQUAL(geoInfo.isInside(insideTestWC3), geoInfo.isInside(expectedIndexWC3)); CHECK(!(geoInfo.isInside(insideTestWC4))); CHECK_EQUAL(geoInfo.isInside(insideTestWC4), geoInfo.isInside(expectedIndexWC4)); CHECK(geoInfo.worldCoordinateToIndex(insideTestWC1, testConvert)); CHECK(geoInfo.isInside(testConvert)); CHECK_EQUAL(expectedIndexWC1, testConvert); CHECK(geoInfo.worldCoordinateToIndex(insideTestWC2, testConvert)); CHECK(geoInfo.isInside(testConvert)); CHECK_EQUAL(expectedIndexWC2, testConvert); CHECK(!(geoInfo.worldCoordinateToIndex(insideTestWC3, testConvert))); CHECK(!(geoInfo.isInside(testConvert))); CHECK_EQUAL(expectedIndexWC3, testConvert); CHECK(!(geoInfo.worldCoordinateToIndex(insideTestWC4, testConvert))); CHECK(!(geoInfo.isInside(testConvert))); CHECK_EQUAL(expectedIndexWC4, testConvert); //11) test isInside and index to world coordinate conversion //use unit matrix as orientation matrix CHECK_NO_THROW(geoInfo.setOrientationMatrix(OrientationMatrix())); VoxelGridIndex3D insideTest1(0, 0, 0); //origin (inside) VoxelGridIndex3D insideTest2(2, 3, 1); //inside VoxelGridIndex3D insideTest3(0, 6, 14); //outside VoxelGridIndex3D insideTest4 = expectedVoxelDims; // outside CHECK(geoInfo.isInside(insideTest1)); CHECK(geoInfo.isInside(insideTest2)); CHECK(!(geoInfo.isInside(insideTest3))); CHECK(!(geoInfo.isInside(insideTest4))); WorldCoordinate3D testInside(0); CHECK(geoInfo.indexToWorldCoordinate(insideTest1, testInside)); CHECK(geoInfo.isInside(testInside)); //CHECK_EQUAL(geoInfo.getImagePositionPatient(),testInside); //half voxel shift prevents equality! CHECK(geoInfo.indexToWorldCoordinate(insideTest2, testInside)); CHECK(geoInfo.isInside(testInside)); CHECK(!(geoInfo.indexToWorldCoordinate(insideTest3, testInside))); CHECK(!(geoInfo.isInside(testInside))); CHECK(!(geoInfo.indexToWorldCoordinate(insideTest4, testInside))); CHECK(!(geoInfo.isInside(testInside))); WorldCoordinate3D testWorldCoordinate; ContinousVoxelGridIndex3D testDoubleIndex; ContinousVoxelGridIndex3D doubleIndex1 = ContinousVoxelGridIndex3D(0.1, 0, -0.3); const WorldCoordinate3D expectedDoubleIndex1(20.1, 100, -1000.3); ContinousVoxelGridIndex3D doubleIndex2 = ContinousVoxelGridIndex3D(11, 6, 15); //outside const WorldCoordinate3D expectedDoubleIndex2(31, 106, -985); ContinousVoxelGridIndex3D doubleIndex3 = ContinousVoxelGridIndex3D(10.1, 5.0, 3.0); // outside: Grid dimension = [10,5,3] const WorldCoordinate3D expectedDoubleIndex3(30.1, 105, -997); ContinousVoxelGridIndex3D doubleIndex4 = ContinousVoxelGridIndex3D(0.0, 0.0, 0.0); const WorldCoordinate3D expectedDoubleIndex4 = geoInfo.getImagePositionPatient(); //test double index to world coordinate - geoInfo.geometryCoordinateToWorldCoordinate(doubleIndex1, testWorldCoordinate); + geoInfo.continousIndexToWorldCoordinate(doubleIndex1, testWorldCoordinate); CHECK_EQUAL(testWorldCoordinate, expectedDoubleIndex1); - geoInfo.geometryCoordinateToWorldCoordinate(doubleIndex2, testWorldCoordinate); + geoInfo.continousIndexToWorldCoordinate(doubleIndex2, testWorldCoordinate); CHECK_EQUAL(testWorldCoordinate, expectedDoubleIndex2); - geoInfo.geometryCoordinateToWorldCoordinate(doubleIndex3, testWorldCoordinate); + geoInfo.continousIndexToWorldCoordinate(doubleIndex3, testWorldCoordinate); CHECK_EQUAL(testWorldCoordinate, expectedDoubleIndex3); - geoInfo.geometryCoordinateToWorldCoordinate(doubleIndex4, testWorldCoordinate); + geoInfo.continousIndexToWorldCoordinate(doubleIndex4, testWorldCoordinate); CHECK_EQUAL(testWorldCoordinate, expectedDoubleIndex4); geoInfo.worldCoordinateToGeometryCoordinate(expectedDoubleIndex4, testDoubleIndex); CHECK_EQUAL(testDoubleIndex, doubleIndex4); geoInfo.worldCoordinateToGeometryCoordinate(expectedDoubleIndex3, testDoubleIndex); CHECK_CLOSE(testDoubleIndex(0), doubleIndex3(0), errorConstant); CHECK_CLOSE(testDoubleIndex(1), doubleIndex3(1), errorConstant); CHECK_CLOSE(testDoubleIndex(2), doubleIndex3(2), errorConstant); geoInfo.worldCoordinateToGeometryCoordinate(expectedDoubleIndex2, testDoubleIndex); CHECK_EQUAL(testDoubleIndex, doubleIndex2); geoInfo.worldCoordinateToGeometryCoordinate(expectedDoubleIndex1, testDoubleIndex); CHECK_CLOSE(testDoubleIndex(0), doubleIndex1(0), errorConstant); CHECK_CLOSE(testDoubleIndex(1), doubleIndex1(1), errorConstant); CHECK_CLOSE(testDoubleIndex(2), doubleIndex1(2), errorConstant); VoxelGridIndex3D testIntIndex; geoInfo.worldCoordinateToIndex(expectedDoubleIndex4, testIntIndex); CHECK_EQUAL(testIntIndex, insideTest1); geoInfo.worldCoordinateToIndex(expectedDoubleIndex1, testIntIndex); CHECK_EQUAL(testIntIndex, insideTest1); geoInfo.worldCoordinateToIndex(expectedDoubleIndex3, testIntIndex); CHECK_EQUAL(testIntIndex, expectedVoxelDims); //use a more complicated orientation matrix newOrientation = OrientationMatrix(0); newOrientation(0, 0) = 0.5; newOrientation(1, 2) = -3; newOrientation(2, 1) = 1; CHECK_NO_THROW(geoInfo.setOrientationMatrix(newOrientation)); testIPP = WorldCoordinate3D(20, 100, -1000); CHECK_NO_THROW(geoInfo.setImagePositionPatient(testIPP)); CHECK_NO_THROW(geoInfo.setSpacing(SpacingVectorType3D(1))); //values for testing were generated with a dedicated MeVisLab routine //no half voxel shift anymore because we changed indexToWorldCoordinate/worldCoordinateToIndex insideTest1 = VoxelGridIndex3D(0, 0, 0); //origin (inside) const WorldCoordinate3D expectedIndex1(20, 100, -1000); insideTest2 = VoxelGridIndex3D(6, 0, 2); //inside const WorldCoordinate3D expectedIndex2(23, 94, -1000); insideTest3 = VoxelGridIndex3D(11, 6, 15); //outside const WorldCoordinate3D expectedIndex3(25.5, 55, -994); insideTest4 = VoxelGridIndex3D(10, 5, 3); // outside: Grid dimension = [10,5,3] const WorldCoordinate3D expectedIndex4(25, 91, -995); CHECK(geoInfo.isInside(insideTest1)); CHECK_EQUAL(geoInfo.isInside(insideTest1), geoInfo.isInside(expectedIndex1)); CHECK(geoInfo.isInside(insideTest2)); CHECK_EQUAL(geoInfo.isInside(insideTest2), geoInfo.isInside(expectedIndex2)); CHECK(!(geoInfo.isInside(insideTest3))); CHECK_EQUAL(geoInfo.isInside(insideTest3), geoInfo.isInside(expectedIndex3)); CHECK(!(geoInfo.isInside(insideTest4))); CHECK_EQUAL(geoInfo.isInside(insideTest4), geoInfo.isInside(expectedIndex4)); CHECK(geoInfo.indexToWorldCoordinate(insideTest1, testInside)); CHECK(geoInfo.isInside(testInside)); CHECK_EQUAL(expectedIndex1, testInside); CHECK(geoInfo.indexToWorldCoordinate(insideTest2, testInside)); CHECK(geoInfo.isInside(testInside)); CHECK_EQUAL(expectedIndex2, testInside); CHECK(!(geoInfo.indexToWorldCoordinate(insideTest3, testInside))); CHECK(!(geoInfo.isInside(testInside))); CHECK_EQUAL(expectedIndex3, testInside); CHECK(!(geoInfo.indexToWorldCoordinate(insideTest4, testInside))); CHECK(!(geoInfo.isInside(testInside))); CHECK_EQUAL(expectedIndex4, testInside); - //12) test with simple Geometry: isInside, geometryCoordinateToWorldCoordinate(), worldCoordinateToGeometryCoordinate(), indexToWorldCoordinate() + //12) test with simple Geometry: isInside, continousIndexToWorldCoordinate(), worldCoordinateToGeometryCoordinate(), indexToWorldCoordinate() core::GeometricInfo geoInfoSimple; ImageSize rttbSimpleSize = ImageSize(10, 10, 10); geoInfoSimple.setImageSize(rttbSimpleSize); SpacingVectorType3D spacingSimple(1, 1, 1); geoInfoSimple.setSpacing(spacingSimple); OrientationMatrix OMOnes; geoInfoSimple.setOrientationMatrix(OMOnes); const ContinousVoxelGridIndex3D doubleIndexPixelOutside1 = ContinousVoxelGridIndex3D(-0.501, 0.0, 0.0); const ContinousVoxelGridIndex3D doubleIndexPixelOutside2 = ContinousVoxelGridIndex3D(0.0, 9.501, 0.0); const ContinousVoxelGridIndex3D doubleIndexPixelZero1 = ContinousVoxelGridIndex3D(0.0, 0.0, 0.0); const ContinousVoxelGridIndex3D doubleIndexPixelZero2 = ContinousVoxelGridIndex3D(-0.5, -0.5, -0.5); const ContinousVoxelGridIndex3D doubleIndexPixelZero3 = ContinousVoxelGridIndex3D(0.499999, 0.499999, 0.499999); const ContinousVoxelGridIndex3D doubleIndexPixelOne1 = ContinousVoxelGridIndex3D(1.0, 0.0, 0.0); const ContinousVoxelGridIndex3D doubleIndexPixelOne2 = ContinousVoxelGridIndex3D(0.5, 0.499999, 0.499999); const ContinousVoxelGridIndex3D doubleIndexPixelOne3 = ContinousVoxelGridIndex3D(1.49, -0.5, -0.5); const ContinousVoxelGridIndex3D doubleIndexPixelLast1 = ContinousVoxelGridIndex3D(9.0, 9.0, 9.0); const ContinousVoxelGridIndex3D doubleIndexPixelLast2 = ContinousVoxelGridIndex3D(9.4999, 9.4999, 9.4999); const ContinousVoxelGridIndex3D doubleIndexPixelLast3 = ContinousVoxelGridIndex3D(8.501, 8.501, 8.501); const VoxelGridIndex3D indexPixelOutside = VoxelGridIndex3D(11, 0, 0); const VoxelGridIndex3D indexPixelZero = VoxelGridIndex3D(0, 0, 0); const VoxelGridIndex3D indexPixelOne = VoxelGridIndex3D(1, 0, 0); const VoxelGridIndex3D indexPixelLast = VoxelGridIndex3D(9, 9, 9); const WorldCoordinate3D worldCoordinateOutside1(doubleIndexPixelOutside1(0), doubleIndexPixelOutside1(1), doubleIndexPixelOutside1(2)); const WorldCoordinate3D worldCoordinateOutside2(doubleIndexPixelOutside2(0), doubleIndexPixelOutside2(1), doubleIndexPixelOutside2(2)); const WorldCoordinate3D worldCoordinatePixelZero1(doubleIndexPixelZero1(0), doubleIndexPixelZero1(1), doubleIndexPixelZero1(2)); const WorldCoordinate3D worldCoordinatePixelZero2(doubleIndexPixelZero2(0), doubleIndexPixelZero2(1), doubleIndexPixelZero2(2)); const WorldCoordinate3D worldCoordinatePixelZero3(doubleIndexPixelZero3(0), doubleIndexPixelZero3(1), doubleIndexPixelZero3(2)); const WorldCoordinate3D worldCoordinatePixelOne1(doubleIndexPixelOne1(0), doubleIndexPixelOne1(1), doubleIndexPixelOne1(2)); const WorldCoordinate3D worldCoordinatePixelOne2(doubleIndexPixelOne2(0), doubleIndexPixelOne2(1), doubleIndexPixelOne2(2)); const WorldCoordinate3D worldCoordinatePixelOne3(doubleIndexPixelOne3(0), doubleIndexPixelOne3(1), doubleIndexPixelOne3(2)); const WorldCoordinate3D worldCoordinatePixelLast1(doubleIndexPixelLast1(0), doubleIndexPixelLast1(1), doubleIndexPixelLast1(2)); const WorldCoordinate3D worldCoordinatePixelLast2(doubleIndexPixelLast2(0), doubleIndexPixelLast2(1), doubleIndexPixelLast2(2)); const WorldCoordinate3D worldCoordinatePixelLast3(doubleIndexPixelLast3(0), doubleIndexPixelLast3(1), doubleIndexPixelLast3(2)); bool isInside; - isInside = geoInfoSimple.geometryCoordinateToWorldCoordinate(doubleIndexPixelOutside1, testWorldCoordinate); + isInside = geoInfoSimple.continousIndexToWorldCoordinate(doubleIndexPixelOutside1, testWorldCoordinate); CHECK(!isInside); - isInside = geoInfoSimple.geometryCoordinateToWorldCoordinate(doubleIndexPixelOutside2, testWorldCoordinate); + isInside = geoInfoSimple.continousIndexToWorldCoordinate(doubleIndexPixelOutside2, testWorldCoordinate); CHECK(!isInside); - isInside = geoInfoSimple.geometryCoordinateToWorldCoordinate(doubleIndexPixelZero1, testWorldCoordinate); + isInside = geoInfoSimple.continousIndexToWorldCoordinate(doubleIndexPixelZero1, testWorldCoordinate); CHECK_EQUAL(testWorldCoordinate, worldCoordinatePixelZero1); CHECK(isInside); - isInside = geoInfoSimple.geometryCoordinateToWorldCoordinate(doubleIndexPixelZero2, testWorldCoordinate); + isInside = geoInfoSimple.continousIndexToWorldCoordinate(doubleIndexPixelZero2, testWorldCoordinate); CHECK_EQUAL(testWorldCoordinate, worldCoordinatePixelZero2); CHECK(isInside); - isInside = geoInfoSimple.geometryCoordinateToWorldCoordinate(doubleIndexPixelZero3, testWorldCoordinate); + isInside = geoInfoSimple.continousIndexToWorldCoordinate(doubleIndexPixelZero3, testWorldCoordinate); CHECK_EQUAL(testWorldCoordinate, worldCoordinatePixelZero3); CHECK(isInside); - isInside = geoInfoSimple.geometryCoordinateToWorldCoordinate(doubleIndexPixelOne1, testWorldCoordinate); + isInside = geoInfoSimple.continousIndexToWorldCoordinate(doubleIndexPixelOne1, testWorldCoordinate); CHECK_EQUAL(testWorldCoordinate, worldCoordinatePixelOne1); CHECK(isInside); - isInside = geoInfoSimple.geometryCoordinateToWorldCoordinate(doubleIndexPixelOne2, testWorldCoordinate); + isInside = geoInfoSimple.continousIndexToWorldCoordinate(doubleIndexPixelOne2, testWorldCoordinate); CHECK_EQUAL(testWorldCoordinate, worldCoordinatePixelOne2); CHECK(isInside); - isInside = geoInfoSimple.geometryCoordinateToWorldCoordinate(doubleIndexPixelOne3, testWorldCoordinate); + isInside = geoInfoSimple.continousIndexToWorldCoordinate(doubleIndexPixelOne3, testWorldCoordinate); CHECK_EQUAL(testWorldCoordinate, worldCoordinatePixelOne3); CHECK(isInside); - isInside = geoInfoSimple.geometryCoordinateToWorldCoordinate(doubleIndexPixelLast1, testWorldCoordinate); + isInside = geoInfoSimple.continousIndexToWorldCoordinate(doubleIndexPixelLast1, testWorldCoordinate); CHECK_EQUAL(testWorldCoordinate, worldCoordinatePixelLast1); CHECK(isInside); - isInside = geoInfoSimple.geometryCoordinateToWorldCoordinate(doubleIndexPixelLast2, testWorldCoordinate); + isInside = geoInfoSimple.continousIndexToWorldCoordinate(doubleIndexPixelLast2, testWorldCoordinate); CHECK_EQUAL(testWorldCoordinate, worldCoordinatePixelLast2); CHECK(isInside); - isInside = geoInfoSimple.geometryCoordinateToWorldCoordinate(doubleIndexPixelLast3, testWorldCoordinate); + isInside = geoInfoSimple.continousIndexToWorldCoordinate(doubleIndexPixelLast3, testWorldCoordinate); CHECK_EQUAL(testWorldCoordinate, worldCoordinatePixelLast3); CHECK(isInside); isInside = geoInfoSimple.worldCoordinateToGeometryCoordinate(worldCoordinateOutside1, testDoubleIndex); CHECK(!isInside); isInside = geoInfoSimple.worldCoordinateToGeometryCoordinate(worldCoordinateOutside2, testDoubleIndex); CHECK(!isInside); isInside = geoInfoSimple.worldCoordinateToGeometryCoordinate(worldCoordinatePixelZero1, testDoubleIndex); CHECK_EQUAL(testDoubleIndex, doubleIndexPixelZero1); CHECK(isInside); isInside = geoInfoSimple.worldCoordinateToGeometryCoordinate(worldCoordinatePixelZero2, testDoubleIndex); CHECK_EQUAL(testDoubleIndex, doubleIndexPixelZero2); CHECK(isInside); isInside = geoInfoSimple.worldCoordinateToGeometryCoordinate(worldCoordinatePixelZero3, testDoubleIndex); CHECK_EQUAL(testDoubleIndex, doubleIndexPixelZero3); CHECK(isInside); isInside = geoInfoSimple.worldCoordinateToGeometryCoordinate(worldCoordinatePixelOne1, testDoubleIndex); CHECK_EQUAL(testDoubleIndex, doubleIndexPixelOne1); CHECK(isInside); isInside = geoInfoSimple.worldCoordinateToGeometryCoordinate(worldCoordinatePixelOne2, testDoubleIndex); CHECK_EQUAL(testDoubleIndex, doubleIndexPixelOne2); CHECK(isInside); isInside = geoInfoSimple.worldCoordinateToGeometryCoordinate(worldCoordinatePixelOne3, testDoubleIndex); CHECK_EQUAL(testDoubleIndex, doubleIndexPixelOne3); CHECK(isInside); isInside = geoInfoSimple.worldCoordinateToGeometryCoordinate(worldCoordinatePixelLast1, testDoubleIndex); CHECK_EQUAL(testDoubleIndex, doubleIndexPixelLast1); CHECK(isInside); isInside = geoInfoSimple.worldCoordinateToGeometryCoordinate(worldCoordinatePixelLast2, testDoubleIndex); CHECK_EQUAL(testDoubleIndex, doubleIndexPixelLast2); CHECK(isInside); isInside = geoInfoSimple.worldCoordinateToGeometryCoordinate(worldCoordinatePixelLast3, testDoubleIndex); CHECK_EQUAL(testDoubleIndex, doubleIndexPixelLast3); CHECK(isInside); isInside = geoInfoSimple.indexToWorldCoordinate(indexPixelOutside, testWorldCoordinate); CHECK(!isInside); isInside = geoInfoSimple.indexToWorldCoordinate(indexPixelZero, testWorldCoordinate); CHECK(isInside); CHECK_EQUAL(testWorldCoordinate, worldCoordinatePixelZero1); isInside = geoInfoSimple.indexToWorldCoordinate(indexPixelOne, testWorldCoordinate); CHECK(isInside); CHECK_EQUAL(testWorldCoordinate, worldCoordinatePixelOne1); isInside = geoInfoSimple.indexToWorldCoordinate(indexPixelLast, testWorldCoordinate); CHECK(isInside); CHECK_EQUAL(testWorldCoordinate, worldCoordinatePixelLast1); isInside = geoInfoSimple.worldCoordinateToIndex(worldCoordinateOutside1, testIntIndex); CHECK(!isInside); isInside = geoInfoSimple.worldCoordinateToIndex(worldCoordinateOutside2, testIntIndex); CHECK(!isInside); isInside = geoInfoSimple.worldCoordinateToIndex(worldCoordinatePixelZero1, testIntIndex); CHECK(isInside); CHECK_EQUAL(testIntIndex, indexPixelZero); isInside = geoInfoSimple.worldCoordinateToIndex(worldCoordinatePixelZero2, testIntIndex); CHECK(isInside); CHECK_EQUAL(testIntIndex, indexPixelZero); isInside = geoInfoSimple.worldCoordinateToIndex(worldCoordinatePixelZero3, testIntIndex); CHECK(isInside); CHECK_EQUAL(testIntIndex, indexPixelZero); isInside = geoInfoSimple.worldCoordinateToIndex(worldCoordinatePixelOne1, testIntIndex); CHECK(isInside); CHECK_EQUAL(testIntIndex, indexPixelOne); isInside = geoInfoSimple.worldCoordinateToIndex(worldCoordinatePixelOne2, testIntIndex); CHECK(isInside); CHECK_EQUAL(testIntIndex, indexPixelOne); isInside = geoInfoSimple.worldCoordinateToIndex(worldCoordinatePixelOne3, testIntIndex); CHECK(isInside); CHECK_EQUAL(testIntIndex, indexPixelOne); isInside = geoInfoSimple.worldCoordinateToIndex(worldCoordinatePixelLast1, testIntIndex); CHECK(isInside); CHECK_EQUAL(testIntIndex, indexPixelLast); isInside = geoInfoSimple.worldCoordinateToIndex(worldCoordinatePixelLast2, testIntIndex); CHECK(isInside); CHECK_EQUAL(testIntIndex, indexPixelLast); isInside = geoInfoSimple.worldCoordinateToIndex(worldCoordinatePixelLast3, testIntIndex); CHECK(isInside); CHECK_EQUAL(testIntIndex, indexPixelLast); //13) test getNumberOfVoxels CHECK_EQUAL(expectedVoxelDims(0)*expectedVoxelDims(1)*expectedVoxelDims(2), geoInfo.getNumberOfVoxels()); //14) Test convert, validID and validIndex geoInfo.setNumColumns(50); geoInfo.setNumRows(30); geoInfo.setNumSlices(40); VoxelGridIndex3D startIndex(0, 0, 0); VoxelGridID startId(0); VoxelGridIndex3D endIndex(geoInfo.getNumColumns() - 1, geoInfo.getNumRows() - 1, geoInfo.getNumSlices() - 1); VoxelGridID endId((geoInfo.getNumColumns()*geoInfo.getNumRows()*geoInfo.getNumSlices()) - 1); VoxelGridIndex3D indexInvalid(geoInfo.getNumColumns(), geoInfo.getNumRows(), geoInfo.getNumSlices()); VoxelGridID idInvalid(geoInfo.getNumColumns()*geoInfo.getNumRows()*geoInfo.getNumSlices()); CHECK(geoInfo.validID(startId)); CHECK(geoInfo.validIndex(startIndex)); VoxelGridIndex3D aIndex; VoxelGridID aId; CHECK(geoInfo.convert(startIndex, aId)); CHECK(geoInfo.convert(startId, aIndex)); CHECK_EQUAL(aId, startId); CHECK_EQUAL(aIndex, startIndex); CHECK(geoInfo.validID(endId)); CHECK(geoInfo.validIndex(endIndex)); CHECK(geoInfo.convert(endIndex, aId)); CHECK(geoInfo.convert(endId, aIndex)); CHECK_EQUAL(aId, endId); CHECK_EQUAL(aIndex, endIndex); CHECK(!geoInfo.validID(idInvalid)); CHECK(!geoInfo.validIndex(indexInvalid)); CHECK(!geoInfo.convert(idInvalid, aIndex)); CHECK(!geoInfo.convert(indexInvalid, aId)); RETURN_AND_REPORT_TEST_SUCCESS; } }//end namespace testing }//end namespace rttb