diff --git a/code/interpolation/rttbInterpolationBase.cpp b/code/interpolation/rttbInterpolationBase.cpp index 65f134a..3d95507 100644 --- a/code/interpolation/rttbInterpolationBase.cpp +++ b/code/interpolation/rttbInterpolationBase.cpp @@ -1,191 +1,196 @@ // ----------------------------------------------------------------------- // 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. // //------------------------------------------------------------------------ /*! // @file // @version $Revision$ (last changed revision) // @date $Date$ (last change date) // @author $Author$ (last changed by) */ +#include + #include "rttbInterpolationBase.h" #include "rttbInvalidParameterException.h" #include "rttbNullPointerException.h" #include "rttbMappingOutsideOfImageException.h" namespace rttb { namespace interpolation { void InterpolationBase::setDoseAccessorPointer(const DoseAccessorPointer originalDose) { if (originalDose != NULL) { _spOriginalDose = originalDose; } else { throw core::NullPointerException("originalDose is NULL!"); } }; void InterpolationBase::getNeighborhoodVoxelValues( const WorldCoordinate3D& aWorldCoordinate, unsigned int neighborhood, boost::array& target, - boost::shared_array values) const + boost::shared_ptr values) const { if (_spOriginalDose == NULL) { throw core::NullPointerException("originalDose is NULL!"); } //Determine target (abs(desired worldCoordinate- corner pixel world coordinate/pixel spacing) and values of corner pixels (from originalDose) VoxelGridIndex3D aIndex; if (_spOriginalDose->getGeometricInfo().worldCoordinateToIndex(aWorldCoordinate, aIndex)) { //determine just the nearest voxel to the world coordinate if (neighborhood == 0) { values[0] = _spOriginalDose->getDoseAt(aIndex); } //determine the 8 voxels around the world coordinate else if (neighborhood == 8) { - std::vector cornerPoints; + std::list cornerPoints; WorldCoordinate3D theNextVoxel; _spOriginalDose->getGeometricInfo().indexToWorldCoordinate(aIndex, theNextVoxel); SpacingVectorType3D pixelSpacing = (_spOriginalDose->getGeometricInfo()).getSpacing(); VoxelGridIndex3D leftTopFrontCoordinate; //find the voxel with the smallest coordinate values in each dimension. This defines the standard cube for (int i = 0; i < 3; i++) { if (aWorldCoordinate[i] < theNextVoxel[i]) { if (aIndex[i] > 0) { leftTopFrontCoordinate[i] = aIndex[i] - 1; target[i] = (aWorldCoordinate[i] - (theNextVoxel[i] - pixelSpacing[i])) / pixelSpacing[i]; } //@todo: this is a workaround, not a good solution else { leftTopFrontCoordinate[i] = aIndex[i]; target[i] = (aWorldCoordinate[i] - (theNextVoxel[i] - pixelSpacing[i])) / pixelSpacing[i]; } } else { leftTopFrontCoordinate[i] = aIndex[i]; target[i] = (aWorldCoordinate[i] - theNextVoxel[i]) / pixelSpacing[i]; } } for (int zIncr = 0; zIncr < 2; zIncr++) { for (int yIncr = 0; yIncr < 2; yIncr++) { for (int xIncr = 0; xIncr < 2; xIncr++) { cornerPoints.push_back(VoxelGridIndex3D(leftTopFrontCoordinate[0] + xIncr, leftTopFrontCoordinate[1] + yIncr, leftTopFrontCoordinate[2] + zIncr)); } } } //target range has to be always [0,1] for (int i = 0; i < 3; i++) { assert(target[i] >= 0.0 && target[i] <= 1.0); } + unsigned int count = 0; + //now just get the values of all (dose) voxels and store them in values - for (int i = 0; i < cornerPoints.size(); ++i) + for (auto cornerPointsIterator = std::begin(cornerPoints); cornerPointsIterator != std::end(cornerPoints); + ++cornerPointsIterator, ++count) { - if (_spOriginalDose->getGeometricInfo().isInside(cornerPoints.at(i))) + if (_spOriginalDose->getGeometricInfo().isInside(*cornerPointsIterator)) { - values[i] = _spOriginalDose->getDoseAt(cornerPoints.at(i)); + values[count] = _spOriginalDose->getDoseAt(*cornerPointsIterator); } else { //outside value! boundary treatment - values[i] = getNearestInsideVoxelValue(cornerPoints.at(i)); + values[count] = getNearestInsideVoxelValue(*cornerPointsIterator); } - assert(values[i] != -1); + assert(values[count] != -1); } } else { throw core::InvalidParameterException("neighborhoods other than 0 and 8 not yet supported in Interpolation"); } } else { throw core::MappingOutsideOfImageException("Error in conversion from world coordinates to index"); } } DoseTypeGy InterpolationBase::getNearestInsideVoxelValue(const VoxelGridIndex3D& currentVoxelIndex) const { VoxelGridIndex3D voxelChangedXYZ[] = {currentVoxelIndex, currentVoxelIndex, currentVoxelIndex, currentVoxelIndex, currentVoxelIndex, currentVoxelIndex, currentVoxelIndex}; int runningIndex; //x,y,z for (runningIndex = 0; runningIndex < 3; ++runningIndex) { voxelChangedXYZ[runningIndex][runningIndex] -= 1; } //xy voxelChangedXYZ[runningIndex][0] -= 1; voxelChangedXYZ[runningIndex][1] -= 1; ++runningIndex; //xz voxelChangedXYZ[runningIndex][0] -= 1; voxelChangedXYZ[runningIndex][2] -= 1; ++runningIndex; //yz voxelChangedXYZ[runningIndex][1] -= 1; voxelChangedXYZ[runningIndex][2] -= 1; ++runningIndex; //xyz voxelChangedXYZ[runningIndex][0] -= 1; voxelChangedXYZ[runningIndex][1] -= 1; voxelChangedXYZ[runningIndex][2] -= 1; ++runningIndex; int replacementVoxelIndex = 0; while (replacementVoxelIndex < runningIndex) { if (_spOriginalDose->getGeometricInfo().validIndex(voxelChangedXYZ[replacementVoxelIndex])) { return _spOriginalDose->getDoseAt(voxelChangedXYZ[replacementVoxelIndex]); } ++replacementVoxelIndex; } return -1; } }//end namespace core }//end namespace rttb diff --git a/code/interpolation/rttbInterpolationBase.h b/code/interpolation/rttbInterpolationBase.h index b0a28a8..2f0d46f 100644 --- a/code/interpolation/rttbInterpolationBase.h +++ b/code/interpolation/rttbInterpolationBase.h @@ -1,95 +1,92 @@ // ----------------------------------------------------------------------- // 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. // //------------------------------------------------------------------------ /*! // @file // @version $Revision$ (last changed revision) // @date $Date$ (last change date) // @author $Author$ (last changed by) */ #ifndef __INTERPOLATION_BASE_H #define __INTERPOLATION_BASE_H -#include #include #include -#include #include #include "rttbDoseAccessorInterface.h" #include "rttbBaseType.h" namespace rttb { namespace interpolation { /*! @class InterpolationBase @brief Base class for interpolation. @ingroup interpolation */ class InterpolationBase { public: typedef boost::shared_ptr Pointer; typedef rttb::core::DoseAccessorInterface::DoseAccessorPointer DoseAccessorPointer; /*! @brief Constructor */ InterpolationBase() {}; /*! @brief Virtual destructor of base class */ virtual ~InterpolationBase() {}; /*! @brief Sets the DoseAccessPointer @pre originalDose initialized @exception core::NullPointerException if originalDose==NULL */ void setDoseAccessorPointer(const DoseAccessorPointer originalDose); /*! @brief Returns the interpolated value for the given world coordinate */ virtual DoseTypeGy getValue(const WorldCoordinate3D& aWorldCoordinate) const = 0; protected: DoseAccessorPointer _spOriginalDose; /*! @brief determines voxels in a certain neighborhood of a physical based coordinate and converts in a standard cube with corner points [0 0 0], [1 0 0], [0 1 0], [1 1 0], [0 0 1], [1 0 1], [0 1 1], [1 1 1]. @param aWorldCoordinate the coordinate where to start @param neighborhood voxel around coordinate (currently only 0 and 8 implemented) @param target coordinates inside the standard cube with values [0 1] in each dimension. - @param values dose values at all corner points of the standard cube. Is of type boost:shared_array[neighborhood] - @pre target and values have to be correctly initialized (e.g. boost::array target = {0.0, 0.0, 0.0}; boost::shared_array values(new DoseTypeGy[8]()); where 8 is neighborhood) + @param values dose values at all corner points of the standard cube. Is of type boost:shared_ptr[neighborhood] + @pre target and values have to be correctly initialized (e.g. boost::array target = {0.0, 0.0, 0.0}; boost::shared_ptr values(new DoseTypeGy[8]()); where 8 is neighborhood) @exception core::InvalidParameterException if neighborhood =! 0 && !=8 @exception core::MappingOutsideOfImageException if initial mapping of aWorldCoordinate is outside image @exception core::NullPointerException if dose is NULL - @todo use shared_ptr[] with boost >1.52 instead of shared_array */ void getNeighborhoodVoxelValues(const WorldCoordinate3D& aWorldCoordinate, - unsigned int neighborhood, boost::array& target, - boost::shared_array values) const; + unsigned int neighborhood, boost::array& target, + boost::shared_ptr values) const; /*! @brief returns the nearest inside voxel value @pre the voxelGridIndex is outside the image and voxelGridIndex>image.size() for all dimensions. Also voxelGridIndex[]>=0 for all dimensions @note used for virtually expanding the image by one voxel as edge handling */ DoseTypeGy getNearestInsideVoxelValue(const VoxelGridIndex3D& currentVoxelIndex) const; }; } } #endif diff --git a/code/interpolation/rttbLinearInterpolation.cpp b/code/interpolation/rttbLinearInterpolation.cpp index e54ab53..60a1ae1 100644 --- a/code/interpolation/rttbLinearInterpolation.cpp +++ b/code/interpolation/rttbLinearInterpolation.cpp @@ -1,57 +1,57 @@ // ----------------------------------------------------------------------- // 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. // //------------------------------------------------------------------------ /*! // @file // @version $Revision$ (last changed revision) // @date $Date$ (last change date) // @author $Author$ (last changed by) */ #include "rttbLinearInterpolation.h" namespace rttb { namespace interpolation { DoseTypeGy LinearInterpolation::trilinear(boost::array target, - boost::shared_array values) const + boost::shared_ptr values) const { //4 linear interpolation in x direction DoseTypeGy c_00 = values[0] * (1.0 - target[0]) + values[1] * target[0]; DoseTypeGy c_10 = values[2] * (1.0 - target[0]) + values[3] * target[0]; DoseTypeGy c_01 = values[4] * (1.0 - target[0]) + values[5] * target[0]; DoseTypeGy c_11 = values[6] * (1.0 - target[0]) + values[7] * target[0]; //combine result in y direction DoseTypeGy c_0 = c_00 * (1.0 - target[1]) + c_10 * target[1]; DoseTypeGy c_1 = c_01 * (1.0 - target[1]) + c_11 * target[1]; //finally incorporate z direction return (c_0 * (1.0 - target[2]) + c_1 * target[2]); } DoseTypeGy LinearInterpolation::getValue(const WorldCoordinate3D& aWorldCoordinate) const { //proper initialization of target and values boost::array target = {0.0, 0.0, 0.0}; - boost::shared_array values(new DoseTypeGy[8]()); + boost::shared_ptr values(new DoseTypeGy[8]()); getNeighborhoodVoxelValues(aWorldCoordinate, 8, target, values); return trilinear(target, values); } } } diff --git a/code/interpolation/rttbLinearInterpolation.h b/code/interpolation/rttbLinearInterpolation.h index 0556017..ea3f3c3 100644 --- a/code/interpolation/rttbLinearInterpolation.h +++ b/code/interpolation/rttbLinearInterpolation.h @@ -1,62 +1,61 @@ // ----------------------------------------------------------------------- // 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. // //------------------------------------------------------------------------ /*! // @file // @version $Revision$ (last changed revision) // @date $Date$ (last change date) // @author $Author$ (last changed by) */ #ifndef __LINEAR_INTERPOLATION_H #define __LINEAR_INTERPOLATION_H -#include #include #include "rttbInterpolationBase.h" namespace rttb { namespace interpolation { /*! @class LinearInterpolation @brief Linear interpolation. @ingroup interpolation */ class LinearInterpolation : public InterpolationBase { public: typedef rttb::core::DoseAccessorInterface::DoseAccessorPointer DoseAccessorPointer; /*! @brief Constructor */ LinearInterpolation() {}; /*! @brief Returns the interpolated value */ DoseTypeGy getValue(const WorldCoordinate3D& aWorldCoordinate) const; private: /*! @brief Trilinar interpolation @sa InterpolationBase for details about target and values @note Source: http://en.wikipedia.org/wiki/Trilinear_interpolation */ - DoseTypeGy trilinear(boost::array target, boost::shared_array values) const; + DoseTypeGy trilinear(boost::array target, boost::shared_ptr values) const; }; } } #endif diff --git a/code/interpolation/rttbNearestNeighborInterpolation.cpp b/code/interpolation/rttbNearestNeighborInterpolation.cpp index e3e0484..69edbf0 100644 --- a/code/interpolation/rttbNearestNeighborInterpolation.cpp +++ b/code/interpolation/rttbNearestNeighborInterpolation.cpp @@ -1,38 +1,40 @@ // ----------------------------------------------------------------------- // 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. // //------------------------------------------------------------------------ /*! // @file // @version $Revision$ (last changed revision) // @date $Date$ (last change date) // @author $Author$ (last changed by) */ #include "rttbNearestNeighborInterpolation.h" +#include + namespace rttb { namespace interpolation { DoseTypeGy NearestNeighborInterpolation::getValue(const WorldCoordinate3D& aWorldCoordinate) const { //proper initialization of target and values (although target is irrelevant in nearest neighbor case) boost::array target = {0.0, 0.0, 0.0}; - boost::shared_array values(new DoseTypeGy[8]()); + boost::shared_ptr values(new DoseTypeGy[8]()); getNeighborhoodVoxelValues(aWorldCoordinate, 0, target, values); return values[0]; } } } diff --git a/code/interpolation/rttbNearestNeighborInterpolation.h b/code/interpolation/rttbNearestNeighborInterpolation.h index 1baf3fa..ffc8079 100644 --- a/code/interpolation/rttbNearestNeighborInterpolation.h +++ b/code/interpolation/rttbNearestNeighborInterpolation.h @@ -1,54 +1,51 @@ // ----------------------------------------------------------------------- // 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. // //------------------------------------------------------------------------ /*! // @file // @version $Revision$ (last changed revision) // @date $Date$ (last change date) // @author $Author$ (last changed by) */ #ifndef __NEAREST_NEIGHBOR_INTERPOLATION_H #define __NEAREST_NEIGHBOR_INTERPOLATION_H -#include -#include - #include "rttbInterpolationBase.h" namespace rttb { namespace interpolation { /*! @class NearestNeighborInterpolation @brief Nearest Neighbor interpolation @ingroup interpolation */ class NearestNeighborInterpolation : public InterpolationBase { public: /*! @brief Constructor */ NearestNeighborInterpolation() {}; /*! @brief Returns the interpolated value (the nearest voxel value given by _spOriginalDose->getGeometricInfo().worldCoordinateToIndex()) */ DoseTypeGy getValue(const WorldCoordinate3D& aWorldCoordinate) const; }; } } #endif diff --git a/testing/interpolation/MatchPointBinding/simpleRegistrationWorkflow.h b/testing/interpolation/MatchPointBinding/simpleRegistrationWorkflow.h index 73db052..71912da 100644 --- a/testing/interpolation/MatchPointBinding/simpleRegistrationWorkflow.h +++ b/testing/interpolation/MatchPointBinding/simpleRegistrationWorkflow.h @@ -1,64 +1,63 @@ // ----------------------------------------------------------------------- // 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. // //------------------------------------------------------------------------ /*! // @file // @version $Revision: 741 $ (last changed revision) // @date $Date: 2014-09-16 16:34:22 +0200 (Di, 16 Sep 2014) $ (last change date) // @author $Author: hentsch $ (last changed by) */ #ifndef __SIMPLE_REGISTRATION_WORKFLOW_H #define __SIMPLE_REGISTRATION_WORKFLOW_H #include -#include #include "registrationHelper.h" #include "mapImageMappingTask.h" #include "mapITKEuler3DMattesMIRegistrationAlgorithmTemplate.h" #include "mapExceptionObject.h" /*! @class simpleRegistrationWorkflow @brief implements a concrete registration algorithm of MatchPoint */ class simpleRegistrationWorkflow { public: typedef map::core::Registration<3, 3> Registration3D3DType; typedef Registration3D3DType::Pointer Registration3D3DPointer; typedef map::algorithm::boxed::ITKEuler3DMattesMIRegistrationAlgorithm AlgorithmTypeEuler; private: std::string _targetFilename; std::string _movingFilename; std::string _targetDirectory; std::string _movingDirectory; AppGlobals globals; AlgorithmTypeEuler::Pointer _spAlgorithmEuler; public: /*! @brief Constructor */ simpleRegistrationWorkflow(std::string targetFilename, std::string movingFilename, bool isDirectory = false); map::core::Registration<3, 3>::Pointer getRegistration(); const itk::Image* getTargetImage(); vnl_vector getRegistrationParameters(Registration3D3DPointer reg); protected: void initializeAndPerformRegistration(); }; #endif \ No newline at end of file