diff --git a/code/masks/CMakeLists.txt b/code/masks/CMakeLists.txt
index 2ef7620..ea7fb93 100644
--- a/code/masks/CMakeLists.txt
+++ b/code/masks/CMakeLists.txt
@@ -1,5 +1,9 @@
-MESSAGE (STATUS "processing RTToolbox mask")
-
-RTTB_CREATE_MODULE(RTTBMasks DEPENDS RTTBCore PACKAGE_DEPENDS BoostBinaries)
-
-ADD_SUBDIRECTORY(boost)
+MESSAGE (STATUS "processing RTToolbox mask")
+
+RTTB_CREATE_MODULE(RTTBMasks DEPENDS RTTBCore PACKAGE_DEPENDS BoostBinaries)
+
+SET(RTTB_Boost_ADDITIONAL_COMPONENT thread)
+RTTB_CREATE_MODULE(RTTBBoostMask DEPENDS RTTBCore RTTBMasks PACKAGE_DEPENDS BoostBinaries)
+IF (CMAKE_COMPILER_IS_GNUCC)
+ SET(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -fext-numeric-literals")
+ENDIF()
diff --git a/code/masks/boost/CMakeLists.txt b/code/masks/boost/CMakeLists.txt
deleted file mode 100644
index 92bb848..0000000
--- a/code/masks/boost/CMakeLists.txt
+++ /dev/null
@@ -1,5 +0,0 @@
-SET(RTTB_Boost_ADDITIONAL_COMPONENT thread)
-RTTB_CREATE_MODULE(RTTBBoostMask DEPENDS RTTBCore RTTBMasks PACKAGE_DEPENDS BoostBinaries)
-IF (CMAKE_COMPILER_IS_GNUCC)
- SET(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -fext-numeric-literals")
-ENDIF()
diff --git a/code/masks/boost/files.cmake b/code/masks/boost/files.cmake
deleted file mode 100644
index faa905d..0000000
--- a/code/masks/boost/files.cmake
+++ /dev/null
@@ -1,13 +0,0 @@
-SET(CPP_FILES
- rttbBoostMask.cpp
- rttbBoostMaskAccessor.cpp
- rttbBoostMaskGenerateMaskVoxelListThread.cpp
- rttbBoostMaskVoxelizationThread.cpp
-)
-
-SET(H_FILES
- rttbBoostMask.h
- rttbBoostMaskAccessor.h
- rttbBoostMaskGenerateMaskVoxelListThread.h
- rttbBoostMaskVoxelizationThread.h
-)
diff --git a/code/masks/files.cmake b/code/masks/files.cmake
index 0cc43df..a3fdf77 100644
--- a/code/masks/files.cmake
+++ b/code/masks/files.cmake
@@ -1,7 +1,15 @@
SET(CPP_FILES
rttbGenericMutableMaskAccessor.cpp
+ rttbBoostMask.cpp
+ rttbBoostMaskAccessor.cpp
+ rttbBoostMaskGenerateMaskVoxelListThread.cpp
+ rttbBoostMaskVoxelizationThread.cpp
)
SET(H_FILES
rttbGenericMutableMaskAccessor.h
+ rttbBoostMask.h
+ rttbBoostMaskAccessor.h
+ rttbBoostMaskGenerateMaskVoxelListThread.h
+ rttbBoostMaskVoxelizationThread.h
)
diff --git a/code/masks/boost/rttbBoostMask.cpp b/code/masks/rttbBoostMask.cpp
similarity index 96%
rename from code/masks/boost/rttbBoostMask.cpp
rename to code/masks/rttbBoostMask.cpp
index bc8cf17..af20a7d 100644
--- a/code/masks/boost/rttbBoostMask.cpp
+++ b/code/masks/rttbBoostMask.cpp
@@ -1,565 +1,565 @@
-// -----------------------------------------------------------------------
-// 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 <limits>
-
-#include <boost/geometry/geometries/register/point.hpp>
-#include <boost/geometry/geometries/register/ring.hpp>
-#include <boost/make_shared.hpp>
-#include <boost/thread/thread.hpp>
-
-#include "rttbBoostMask.h"
-#include "rttbNullPointerException.h"
-#include "rttbInvalidParameterException.h"
-#include "rttbBoostMaskGenerateMaskVoxelListThread.h"
-#include "rttbBoostMaskVoxelizationThread.h"
-
-
-namespace rttb
-{
- namespace masks
- {
- namespace boost
- {
-
-
- BoostMask::BoostMask(core::GeometricInfo::Pointer aDoseGeoInfo,
- core::Structure::Pointer aStructure, bool strict, unsigned int numberOfThreads)
- : _geometricInfo(aDoseGeoInfo), _structure(aStructure),
- _strict(strict), _numberOfThreads(numberOfThreads), _voxelizationThickness(0.0),
- _voxelInStructure(::boost::make_shared<MaskVoxelList>())
- {
-
- _isUpToDate = false;
-
- if (_geometricInfo == nullptr)
- {
- throw rttb::core::NullPointerException("Error: Geometric info is nullptr!");
- }
- else if (_structure == nullptr)
- {
- throw rttb::core::NullPointerException("Error: Structure is nullptr!");
- }
-
- if (_numberOfThreads == 0)
- {
- _numberOfThreads = ::boost::thread::hardware_concurrency();
- if (_numberOfThreads == 0)
- {
- throw rttb::core::InvalidParameterException("Error: detection of the number of hardware threads is not possible. Please specify number of threads for voxelization explicitly as parameter in BoostMask.");
- }
- }
-
- }
-
- BoostMask::MaskVoxelListPointer BoostMask::getRelevantVoxelVector()
- {
- if (!_isUpToDate)
- {
- calcMask();
- }
-
- return _voxelInStructure;
- }
-
- void BoostMask::calcMask()
- {
- preprocessing();
- voxelization();
- generateMaskVoxelList();
- _isUpToDate = true;
- }
-
- void BoostMask::preprocessing()
- {
- rttb::PolygonSequenceType polygonSequence = _structure->getStructureVector();
-
- //Convert world coordinate polygons to the polygons with geometry coordinate
- rttb::PolygonSequenceType geometryCoordinatePolygonVector;
- rttb::PolygonSequenceType::iterator it;
- rttb::DoubleVoxelGridIndex3D globalMaxGridIndex(std::numeric_limits<double>::min(),
- std::numeric_limits<double>::min(), std::numeric_limits<double>::min());
- rttb::DoubleVoxelGridIndex3D globalMinGridIndex(_geometricInfo->getNumColumns(),
- _geometricInfo->getNumRows(), 0);
-
- for (it = polygonSequence.begin(); it != polygonSequence.end(); ++it)
- {
- PolygonType rttbPolygon = *it;
- PolygonType geometryCoordinatePolygon;
-
- //1. convert polygon to geometry coordinate polygons
- //2. calculate global min/max
- //3. check if polygon is planar
- if (!preprocessingPolygon(rttbPolygon, geometryCoordinatePolygon, globalMinGridIndex,
- globalMaxGridIndex, errorConstant))
- {
- throw rttb::core::Exception("TiltedMaskPlaneException");
- }
-
- geometryCoordinatePolygonVector.push_back(geometryCoordinatePolygon);
- }
-
- rttb::VoxelGridIndex3D minIndex = VoxelGridIndex3D(GridIndexType(globalMinGridIndex(0) ),
- GridIndexType(globalMinGridIndex(1) ), GridIndexType(globalMinGridIndex(2) ));
- rttb::VoxelGridIndex3D maxIndex = VoxelGridIndex3D(GridIndexType(globalMaxGridIndex(0) ),
- GridIndexType(globalMaxGridIndex(1) ), GridIndexType(globalMaxGridIndex(2) ));
-
- _globalBoundingBox.push_back(minIndex);
- _globalBoundingBox.push_back(maxIndex);
-
- //convert rttb polygon sequence to a map of z index and a vector of boost ring 2d (without holes)
- _ringMap = convertRTTBPolygonSequenceToBoostRingMap(geometryCoordinatePolygonVector);
-
- }
-
- void BoostMask::voxelization()
- {
-
- if (_globalBoundingBox.size() < 2)
- {
- throw rttb::core::InvalidParameterException("Bounding box calculation failed! ");
- }
-
- BoostRingMap::iterator itMap;
-
- size_t mapSizeInAThread = _ringMap.size() / _numberOfThreads;
- unsigned int count = 0;
- unsigned int countThread = 0;
- BoostPolygonMap polygonMap;
- std::vector<BoostPolygonMap> polygonMapVector;
-
- //check donut and convert to a map of z index and a vector of boost polygon 2d (with or without holes)
- for (itMap = _ringMap.begin(); itMap != _ringMap.end(); ++itMap)
- {
- //the vector of all boost 2d polygons with the same z grid index(donut polygon is accepted).
- BoostPolygonVector polygonVector = checkDonutAndConvert((*itMap).second);
-
- if (count == mapSizeInAThread && countThread < (_numberOfThreads - 1))
- {
- polygonMapVector.push_back(polygonMap);
- polygonMap.clear();
- count = 0;
- countThread++;
- }
-
- polygonMap.insert(std::pair<double, BoostPolygonVector>((*itMap).first,
- polygonVector));
- count++;
-
- }
-
- _voxelizationMap = ::boost::make_shared<std::map<double, BoostArray2DPointer> >();
-
- polygonMapVector.push_back(polygonMap); //insert the last one
-
- //generate voxelization map, multi-threading
- ::boost::thread_group threads;
- auto aMutex = ::boost::make_shared<::boost::shared_mutex>();
-
- for (const auto & i : polygonMapVector)
- {
- BoostMaskVoxelizationThread t(i, _globalBoundingBox,
- _voxelizationMap, aMutex, _strict);
- threads.create_thread(t);
- }
-
- threads.join_all();
- }
-
- void BoostMask::generateMaskVoxelList()
- {
- if (_globalBoundingBox.size() < 2)
- {
- throw rttb::core::InvalidParameterException("Bounding box calculation failed! ");
- }
-
- //check homogeneous of the voxelization plane (the contours plane)
- if (!calcVoxelizationThickness(_voxelizationThickness))
- {
- throw rttb::core::InvalidParameterException("Error: The contour plane should be homogeneous!");
- }
-
-
-
- ::boost::thread_group threads;
- auto aMutex = ::boost::make_shared<::boost::shared_mutex>();
-
- unsigned int sliceNumberInAThread = _geometricInfo->getNumSlices() / _numberOfThreads;
-
- //generate mask voxel list, multi-threading
- for (unsigned int i = 0; i < _numberOfThreads; ++i)
- {
- unsigned int beginSlice = i * sliceNumberInAThread;
- unsigned int endSlice;
-
- if (i < _numberOfThreads - 1)
- {
- endSlice = (i + 1) * sliceNumberInAThread;
- }
- else
- {
- endSlice = _geometricInfo->getNumSlices();
- }
-
-
- BoostMaskGenerateMaskVoxelListThread t(_globalBoundingBox, _geometricInfo, _voxelizationMap,
- _voxelizationThickness, beginSlice, endSlice,
- _voxelInStructure, aMutex);
-
- threads.create_thread(t);
-
- }
-
- threads.join_all();
-
- }
-
- bool BoostMask::preprocessingPolygon(const rttb::PolygonType& aRTTBPolygon,
- rttb::PolygonType& geometryCoordinatePolygon, rttb::DoubleVoxelGridIndex3D& minimum,
- rttb::DoubleVoxelGridIndex3D& maximum, double aErrorConstant) const
- {
-
- double minZ = _geometricInfo->getNumSlices();
- double maxZ = 0.0;
-
- for (auto worldCoordinatePoint : aRTTBPolygon)
- {
- //convert to geometry coordinate polygon
- rttb::DoubleVoxelGridIndex3D geometryCoordinatePoint;
- _geometricInfo->worldCoordinateToGeometryCoordinate(worldCoordinatePoint, geometryCoordinatePoint);
-
- geometryCoordinatePolygon.push_back(geometryCoordinatePoint);
-
- //calculate the current global min/max
- //min and max for x
- if (geometryCoordinatePoint(0) < minimum(0))
- {
- minimum(0) = geometryCoordinatePoint(0);
- }
-
- if (geometryCoordinatePoint(0) > maximum(0))
- {
- maximum(0) = geometryCoordinatePoint(0);
- }
-
- //min and max for y
- if (geometryCoordinatePoint(1) < minimum(1))
- {
- minimum(1) = geometryCoordinatePoint(1);
- }
-
- if (geometryCoordinatePoint(1) > maximum(1))
- {
- maximum(1) = geometryCoordinatePoint(1);
- }
-
- //min and max for z
- if (geometryCoordinatePoint(2) < minimum(2))
- {
- minimum(2) = geometryCoordinatePoint(2);
- }
-
- if (geometryCoordinatePoint(2) > maximum(2))
- {
- maximum(2) = geometryCoordinatePoint(2);
- }
-
- //check planar
- if (geometryCoordinatePoint(2) < minZ)
- {
- minZ = geometryCoordinatePoint(2);
- }
-
- if (geometryCoordinatePoint(2) > maxZ)
- {
- maxZ = geometryCoordinatePoint(2);
- }
-
- }
-
- return (std::abs(maxZ - minZ) <= aErrorConstant);
- }
-
-
- BoostMask::BoostRing2D BoostMask::convertRTTBPolygonToBoostRing(const rttb::PolygonType&
- aRTTBPolygon) const
- {
- BoostMask::BoostRing2D polygon2D;
- BoostPoint2D firstPoint;
-
- for (unsigned int i = 0; i < aRTTBPolygon.size(); i++)
- {
- rttb::WorldCoordinate3D rttbPoint = aRTTBPolygon.at(i);
- BoostPoint2D boostPoint(rttbPoint[0], rttbPoint[1]);
-
- if (i == 0)
- {
- firstPoint = boostPoint;
- }
-
- ::boost::geometry::append(polygon2D, boostPoint);
- }
-
- ::boost::geometry::append(polygon2D, firstPoint);
- return polygon2D;
- }
-
- BoostMask::BoostRingMap BoostMask::convertRTTBPolygonSequenceToBoostRingMap(
- const rttb::PolygonSequenceType& aRTTBPolygonVector) const
- {
- rttb::PolygonSequenceType::const_iterator it;
- BoostMask::BoostRingMap aRingMap;
-
- for (it = aRTTBPolygonVector.begin(); it != aRTTBPolygonVector.end(); ++it)
- {
- rttb::PolygonType rttbPolygon = *it;
- double zIndex = rttbPolygon.at(0)[2];//get the first z index of the polygon
- bool isFirstZ = true;
-
- if (!aRingMap.empty())
- {
- auto findIt = findNearestKey(aRingMap, zIndex, errorConstant);
-
- //if the z index is found (same slice), add the polygon to vector
- if (findIt != aRingMap.end())
- {
- //BoostRingVector ringVector = ;
- (*findIt).second.push_back(convertRTTBPolygonToBoostRing(rttbPolygon));
- isFirstZ = false;
- }
- }
-
- //if it is the first z index in the map, insert vector with the polygon
- if (isFirstZ)
- {
- BoostRingVector ringVector;
- ringVector.push_back(convertRTTBPolygonToBoostRing(rttbPolygon));
- aRingMap.insert(std::pair<double, BoostRingVector>(zIndex, ringVector));
- }
-
-
- }
-
- return aRingMap;
- }
-
- BoostMask::BoostRingMap::iterator BoostMask::findNearestKey(BoostMask::BoostRingMap&
- aBoostRingMap, double aIndex, double aErrorConstant) const
- {
- auto find = aBoostRingMap.find(aIndex);
-
- //if find a key equivalent to aIndex, found
- if (find != aBoostRingMap.end())
- {
- return find;
- }
- else
- {
- auto lowerBound = aBoostRingMap.lower_bound(aIndex);
-
- //if all keys go before aIndex, check the last key
- if (lowerBound == aBoostRingMap.end())
- {
- lowerBound = --aBoostRingMap.end();
- }
-
- //if the lower bound very close to aIndex, found
- if (std::abs((*lowerBound).first - aIndex) <= aErrorConstant)
- {
- return lowerBound;
- }
- else
- {
- //if the lower bound is the beginning, not found
- if (lowerBound == aBoostRingMap.begin())
- {
- return aBoostRingMap.end();
- }
- else
- {
- auto lowerBound1 = --lowerBound;//the key before the lower bound
-
- //if the key before the lower bound very close to a Index, found
- if (std::abs((*lowerBound1).first - aIndex) <= aErrorConstant)
- {
- return lowerBound1;
- }
- //else, not found
- else
- {
- return aBoostRingMap.end();
- }
- }
- }
-
- }
- }
-
- BoostMask::BoostPolygonVector BoostMask::checkDonutAndConvert(const BoostMask::BoostRingVector&
- aRingVector) const
- {
- //check donut
- BoostMask::BoostRingVector::const_iterator it1;
- BoostMask::BoostRingVector::const_iterator it2;
- BoostMask::BoostPolygonVector boostPolygonVector;
- std::vector<unsigned int> donutIndexVector;//store the outer and inner ring index
- BoostMask::BoostPolygonVector donutVector;//store new generated donut polygon
-
- //Get donut index and donut polygon
- unsigned int index1 = 0;
-
- for (it1 = aRingVector.begin(); it1 != aRingVector.end(); ++it1, index1++)
- {
- bool it1IsDonut = false;
-
- //check if the ring is already determined as a donut
- for (unsigned int i : donutIndexVector)
- {
- if (i == index1)
- {
- it1IsDonut = true;
- break;
- }
- }
-
- //if not jet, check now
- if (!it1IsDonut)
- {
- bool it2IsDonut = false;
- unsigned int index2 = 0;
-
- for (it2 = aRingVector.begin(); it2 != aRingVector.end(); ++it2, index2++)
- {
- if (it2 != it1)
- {
- BoostMask::BoostPolygon2D polygon2D;
-
- if (::boost::geometry::within(*it1, *it2))
- {
- ::boost::geometry::append(polygon2D, *it2);//append an outer ring to the polygon
- ::boost::geometry::interior_rings(polygon2D).resize(1);//create an interior ring
- ::boost::geometry::append(polygon2D, *it1, 0);//append a ring to the interior ring
- it2IsDonut = true;
- }
- //if donut
- else if (::boost::geometry::within(*it2, *it1))
- {
- ::boost::geometry::append(polygon2D, *it1);//append an outer ring to the polygon
- ::boost::geometry::interior_rings(polygon2D).resize(1);//create an interior ring
- ::boost::geometry::append(polygon2D, *it2, 0);//append a ring to the interior ring
- it2IsDonut = true;
- }
-
- if (it2IsDonut)
- {
- donutIndexVector.push_back(index1);
- donutIndexVector.push_back(index2);
- donutVector.push_back(polygon2D);//store donut polygon
- break;//Only store the first donut!
- }
- }
- }
- }
- }
-
- //Store no donut polygon to boostPolygonVector
- index1 = 0;
-
- for (it1 = aRingVector.begin(); it1 != aRingVector.end(); ++it1, index1++)
- {
- bool it1IsDonut = false;
-
- //check if the ring is the outer or inner of a donut
- for (unsigned int i : donutIndexVector)
- {
- if (i == index1)
- {
- it1IsDonut = true;
- break;
- }
- }
-
- if (!it1IsDonut)
- {
- BoostMask::BoostPolygon2D polygon2D;
- ::boost::geometry::append(polygon2D, *it1);
- boostPolygonVector.push_back(polygon2D);//insert the ring, which is not a part of donut
- }
- }
-
- //Append donut polygon to boostPolygonVector
- BoostMask::BoostPolygonVector::iterator itDonut;
-
- for (itDonut = donutVector.begin(); itDonut != donutVector.end(); ++itDonut)
- {
- boostPolygonVector.push_back(*itDonut);//append donuts
- }
-
- return boostPolygonVector;
- }
-
- bool BoostMask::calcVoxelizationThickness(double& aThickness) const
- {
-
- if (_voxelizationMap->size() <= 1)
- {
- aThickness = 1;
- return true;
- }
-
- double thickness = 0;
-
- auto it = _voxelizationMap->cbegin();
- auto it2 = ++_voxelizationMap->cbegin();
- for (;
- it != _voxelizationMap->cend() && it2 != _voxelizationMap->cend(); ++it, ++it2)
- {
- if (thickness == 0)
- {
- thickness = it2->first - it->first;
- }
- else
- {
- double curThickness = it2->first - it->first;
- //if not homogeneous (leave out double imprecisions), return false
- if (std::abs(thickness-curThickness)>errorConstant)
- {
- //return false;
- std::cout << "Two polygons are far from each other?" << std::endl;
- }
- }
-
- }
-
- if (thickness != 0)
- {
- aThickness = thickness;
- }
- else
- {
- aThickness = 1;
- }
-
- return true;
- }
- }
- }
-}
+// -----------------------------------------------------------------------
+// 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 <limits>
+
+#include <boost/geometry/geometries/register/point.hpp>
+#include <boost/geometry/geometries/register/ring.hpp>
+#include <boost/make_shared.hpp>
+#include <boost/thread/thread.hpp>
+
+#include "rttbBoostMask.h"
+#include "rttbNullPointerException.h"
+#include "rttbInvalidParameterException.h"
+#include "rttbBoostMaskGenerateMaskVoxelListThread.h"
+#include "rttbBoostMaskVoxelizationThread.h"
+
+
+namespace rttb
+{
+ namespace masks
+ {
+ namespace boost
+ {
+
+
+ BoostMask::BoostMask(core::GeometricInfo::Pointer aDoseGeoInfo,
+ core::Structure::Pointer aStructure, bool strict, unsigned int numberOfThreads)
+ : _geometricInfo(aDoseGeoInfo), _structure(aStructure),
+ _strict(strict), _numberOfThreads(numberOfThreads), _voxelizationThickness(0.0),
+ _voxelInStructure(::boost::make_shared<MaskVoxelList>())
+ {
+
+ _isUpToDate = false;
+
+ if (_geometricInfo == nullptr)
+ {
+ throw rttb::core::NullPointerException("Error: Geometric info is nullptr!");
+ }
+ else if (_structure == nullptr)
+ {
+ throw rttb::core::NullPointerException("Error: Structure is nullptr!");
+ }
+
+ if (_numberOfThreads == 0)
+ {
+ _numberOfThreads = ::boost::thread::hardware_concurrency();
+ if (_numberOfThreads == 0)
+ {
+ throw rttb::core::InvalidParameterException("Error: detection of the number of hardware threads is not possible. Please specify number of threads for voxelization explicitly as parameter in BoostMask.");
+ }
+ }
+
+ }
+
+ BoostMask::MaskVoxelListPointer BoostMask::getRelevantVoxelVector()
+ {
+ if (!_isUpToDate)
+ {
+ calcMask();
+ }
+
+ return _voxelInStructure;
+ }
+
+ void BoostMask::calcMask()
+ {
+ preprocessing();
+ voxelization();
+ generateMaskVoxelList();
+ _isUpToDate = true;
+ }
+
+ void BoostMask::preprocessing()
+ {
+ rttb::PolygonSequenceType polygonSequence = _structure->getStructureVector();
+
+ //Convert world coordinate polygons to the polygons with geometry coordinate
+ rttb::PolygonSequenceType geometryCoordinatePolygonVector;
+ rttb::PolygonSequenceType::iterator it;
+ rttb::DoubleVoxelGridIndex3D globalMaxGridIndex(std::numeric_limits<double>::min(),
+ std::numeric_limits<double>::min(), std::numeric_limits<double>::min());
+ rttb::DoubleVoxelGridIndex3D globalMinGridIndex(_geometricInfo->getNumColumns(),
+ _geometricInfo->getNumRows(), 0);
+
+ for (it = polygonSequence.begin(); it != polygonSequence.end(); ++it)
+ {
+ PolygonType rttbPolygon = *it;
+ PolygonType geometryCoordinatePolygon;
+
+ //1. convert polygon to geometry coordinate polygons
+ //2. calculate global min/max
+ //3. check if polygon is planar
+ if (!preprocessingPolygon(rttbPolygon, geometryCoordinatePolygon, globalMinGridIndex,
+ globalMaxGridIndex, errorConstant))
+ {
+ throw rttb::core::Exception("TiltedMaskPlaneException");
+ }
+
+ geometryCoordinatePolygonVector.push_back(geometryCoordinatePolygon);
+ }
+
+ rttb::VoxelGridIndex3D minIndex = VoxelGridIndex3D(GridIndexType(globalMinGridIndex(0) ),
+ GridIndexType(globalMinGridIndex(1) ), GridIndexType(globalMinGridIndex(2) ));
+ rttb::VoxelGridIndex3D maxIndex = VoxelGridIndex3D(GridIndexType(globalMaxGridIndex(0) ),
+ GridIndexType(globalMaxGridIndex(1) ), GridIndexType(globalMaxGridIndex(2) ));
+
+ _globalBoundingBox.push_back(minIndex);
+ _globalBoundingBox.push_back(maxIndex);
+
+ //convert rttb polygon sequence to a map of z index and a vector of boost ring 2d (without holes)
+ _ringMap = convertRTTBPolygonSequenceToBoostRingMap(geometryCoordinatePolygonVector);
+
+ }
+
+ void BoostMask::voxelization()
+ {
+
+ if (_globalBoundingBox.size() < 2)
+ {
+ throw rttb::core::InvalidParameterException("Bounding box calculation failed! ");
+ }
+
+ BoostRingMap::iterator itMap;
+
+ size_t mapSizeInAThread = _ringMap.size() / _numberOfThreads;
+ unsigned int count = 0;
+ unsigned int countThread = 0;
+ BoostPolygonMap polygonMap;
+ std::vector<BoostPolygonMap> polygonMapVector;
+
+ //check donut and convert to a map of z index and a vector of boost polygon 2d (with or without holes)
+ for (itMap = _ringMap.begin(); itMap != _ringMap.end(); ++itMap)
+ {
+ //the vector of all boost 2d polygons with the same z grid index(donut polygon is accepted).
+ BoostPolygonVector polygonVector = checkDonutAndConvert((*itMap).second);
+
+ if (count == mapSizeInAThread && countThread < (_numberOfThreads - 1))
+ {
+ polygonMapVector.push_back(polygonMap);
+ polygonMap.clear();
+ count = 0;
+ countThread++;
+ }
+
+ polygonMap.insert(std::pair<double, BoostPolygonVector>((*itMap).first,
+ polygonVector));
+ count++;
+
+ }
+
+ _voxelizationMap = ::boost::make_shared<std::map<double, BoostArray2DPointer> >();
+
+ polygonMapVector.push_back(polygonMap); //insert the last one
+
+ //generate voxelization map, multi-threading
+ ::boost::thread_group threads;
+ auto aMutex = ::boost::make_shared<::boost::shared_mutex>();
+
+ for (const auto & i : polygonMapVector)
+ {
+ BoostMaskVoxelizationThread t(i, _globalBoundingBox,
+ _voxelizationMap, aMutex, _strict);
+ threads.create_thread(t);
+ }
+
+ threads.join_all();
+ }
+
+ void BoostMask::generateMaskVoxelList()
+ {
+ if (_globalBoundingBox.size() < 2)
+ {
+ throw rttb::core::InvalidParameterException("Bounding box calculation failed! ");
+ }
+
+ //check homogeneous of the voxelization plane (the contours plane)
+ if (!calcVoxelizationThickness(_voxelizationThickness))
+ {
+ throw rttb::core::InvalidParameterException("Error: The contour plane should be homogeneous!");
+ }
+
+
+
+ ::boost::thread_group threads;
+ auto aMutex = ::boost::make_shared<::boost::shared_mutex>();
+
+ unsigned int sliceNumberInAThread = _geometricInfo->getNumSlices() / _numberOfThreads;
+
+ //generate mask voxel list, multi-threading
+ for (unsigned int i = 0; i < _numberOfThreads; ++i)
+ {
+ unsigned int beginSlice = i * sliceNumberInAThread;
+ unsigned int endSlice;
+
+ if (i < _numberOfThreads - 1)
+ {
+ endSlice = (i + 1) * sliceNumberInAThread;
+ }
+ else
+ {
+ endSlice = _geometricInfo->getNumSlices();
+ }
+
+
+ BoostMaskGenerateMaskVoxelListThread t(_globalBoundingBox, _geometricInfo, _voxelizationMap,
+ _voxelizationThickness, beginSlice, endSlice,
+ _voxelInStructure, aMutex);
+
+ threads.create_thread(t);
+
+ }
+
+ threads.join_all();
+
+ }
+
+ bool BoostMask::preprocessingPolygon(const rttb::PolygonType& aRTTBPolygon,
+ rttb::PolygonType& geometryCoordinatePolygon, rttb::DoubleVoxelGridIndex3D& minimum,
+ rttb::DoubleVoxelGridIndex3D& maximum, double aErrorConstant) const
+ {
+
+ double minZ = _geometricInfo->getNumSlices();
+ double maxZ = 0.0;
+
+ for (auto worldCoordinatePoint : aRTTBPolygon)
+ {
+ //convert to geometry coordinate polygon
+ rttb::DoubleVoxelGridIndex3D geometryCoordinatePoint;
+ _geometricInfo->worldCoordinateToGeometryCoordinate(worldCoordinatePoint, geometryCoordinatePoint);
+
+ geometryCoordinatePolygon.push_back(geometryCoordinatePoint);
+
+ //calculate the current global min/max
+ //min and max for x
+ if (geometryCoordinatePoint(0) < minimum(0))
+ {
+ minimum(0) = geometryCoordinatePoint(0);
+ }
+
+ if (geometryCoordinatePoint(0) > maximum(0))
+ {
+ maximum(0) = geometryCoordinatePoint(0);
+ }
+
+ //min and max for y
+ if (geometryCoordinatePoint(1) < minimum(1))
+ {
+ minimum(1) = geometryCoordinatePoint(1);
+ }
+
+ if (geometryCoordinatePoint(1) > maximum(1))
+ {
+ maximum(1) = geometryCoordinatePoint(1);
+ }
+
+ //min and max for z
+ if (geometryCoordinatePoint(2) < minimum(2))
+ {
+ minimum(2) = geometryCoordinatePoint(2);
+ }
+
+ if (geometryCoordinatePoint(2) > maximum(2))
+ {
+ maximum(2) = geometryCoordinatePoint(2);
+ }
+
+ //check planar
+ if (geometryCoordinatePoint(2) < minZ)
+ {
+ minZ = geometryCoordinatePoint(2);
+ }
+
+ if (geometryCoordinatePoint(2) > maxZ)
+ {
+ maxZ = geometryCoordinatePoint(2);
+ }
+
+ }
+
+ return (std::abs(maxZ - minZ) <= aErrorConstant);
+ }
+
+
+ BoostMask::BoostRing2D BoostMask::convertRTTBPolygonToBoostRing(const rttb::PolygonType&
+ aRTTBPolygon) const
+ {
+ BoostMask::BoostRing2D polygon2D;
+ BoostPoint2D firstPoint;
+
+ for (unsigned int i = 0; i < aRTTBPolygon.size(); i++)
+ {
+ rttb::WorldCoordinate3D rttbPoint = aRTTBPolygon.at(i);
+ BoostPoint2D boostPoint(rttbPoint[0], rttbPoint[1]);
+
+ if (i == 0)
+ {
+ firstPoint = boostPoint;
+ }
+
+ ::boost::geometry::append(polygon2D, boostPoint);
+ }
+
+ ::boost::geometry::append(polygon2D, firstPoint);
+ return polygon2D;
+ }
+
+ BoostMask::BoostRingMap BoostMask::convertRTTBPolygonSequenceToBoostRingMap(
+ const rttb::PolygonSequenceType& aRTTBPolygonVector) const
+ {
+ rttb::PolygonSequenceType::const_iterator it;
+ BoostMask::BoostRingMap aRingMap;
+
+ for (it = aRTTBPolygonVector.begin(); it != aRTTBPolygonVector.end(); ++it)
+ {
+ rttb::PolygonType rttbPolygon = *it;
+ double zIndex = rttbPolygon.at(0)[2];//get the first z index of the polygon
+ bool isFirstZ = true;
+
+ if (!aRingMap.empty())
+ {
+ auto findIt = findNearestKey(aRingMap, zIndex, errorConstant);
+
+ //if the z index is found (same slice), add the polygon to vector
+ if (findIt != aRingMap.end())
+ {
+ //BoostRingVector ringVector = ;
+ (*findIt).second.push_back(convertRTTBPolygonToBoostRing(rttbPolygon));
+ isFirstZ = false;
+ }
+ }
+
+ //if it is the first z index in the map, insert vector with the polygon
+ if (isFirstZ)
+ {
+ BoostRingVector ringVector;
+ ringVector.push_back(convertRTTBPolygonToBoostRing(rttbPolygon));
+ aRingMap.insert(std::pair<double, BoostRingVector>(zIndex, ringVector));
+ }
+
+
+ }
+
+ return aRingMap;
+ }
+
+ BoostMask::BoostRingMap::iterator BoostMask::findNearestKey(BoostMask::BoostRingMap&
+ aBoostRingMap, double aIndex, double aErrorConstant) const
+ {
+ auto find = aBoostRingMap.find(aIndex);
+
+ //if find a key equivalent to aIndex, found
+ if (find != aBoostRingMap.end())
+ {
+ return find;
+ }
+ else
+ {
+ auto lowerBound = aBoostRingMap.lower_bound(aIndex);
+
+ //if all keys go before aIndex, check the last key
+ if (lowerBound == aBoostRingMap.end())
+ {
+ lowerBound = --aBoostRingMap.end();
+ }
+
+ //if the lower bound very close to aIndex, found
+ if (std::abs((*lowerBound).first - aIndex) <= aErrorConstant)
+ {
+ return lowerBound;
+ }
+ else
+ {
+ //if the lower bound is the beginning, not found
+ if (lowerBound == aBoostRingMap.begin())
+ {
+ return aBoostRingMap.end();
+ }
+ else
+ {
+ auto lowerBound1 = --lowerBound;//the key before the lower bound
+
+ //if the key before the lower bound very close to a Index, found
+ if (std::abs((*lowerBound1).first - aIndex) <= aErrorConstant)
+ {
+ return lowerBound1;
+ }
+ //else, not found
+ else
+ {
+ return aBoostRingMap.end();
+ }
+ }
+ }
+
+ }
+ }
+
+ BoostMask::BoostPolygonVector BoostMask::checkDonutAndConvert(const BoostMask::BoostRingVector&
+ aRingVector) const
+ {
+ //check donut
+ BoostMask::BoostRingVector::const_iterator it1;
+ BoostMask::BoostRingVector::const_iterator it2;
+ BoostMask::BoostPolygonVector boostPolygonVector;
+ std::vector<unsigned int> donutIndexVector;//store the outer and inner ring index
+ BoostMask::BoostPolygonVector donutVector;//store new generated donut polygon
+
+ //Get donut index and donut polygon
+ unsigned int index1 = 0;
+
+ for (it1 = aRingVector.begin(); it1 != aRingVector.end(); ++it1, index1++)
+ {
+ bool it1IsDonut = false;
+
+ //check if the ring is already determined as a donut
+ for (unsigned int i : donutIndexVector)
+ {
+ if (i == index1)
+ {
+ it1IsDonut = true;
+ break;
+ }
+ }
+
+ //if not jet, check now
+ if (!it1IsDonut)
+ {
+ bool it2IsDonut = false;
+ unsigned int index2 = 0;
+
+ for (it2 = aRingVector.begin(); it2 != aRingVector.end(); ++it2, index2++)
+ {
+ if (it2 != it1)
+ {
+ BoostMask::BoostPolygon2D polygon2D;
+
+ if (::boost::geometry::within(*it1, *it2))
+ {
+ ::boost::geometry::append(polygon2D, *it2);//append an outer ring to the polygon
+ ::boost::geometry::interior_rings(polygon2D).resize(1);//create an interior ring
+ ::boost::geometry::append(polygon2D, *it1, 0);//append a ring to the interior ring
+ it2IsDonut = true;
+ }
+ //if donut
+ else if (::boost::geometry::within(*it2, *it1))
+ {
+ ::boost::geometry::append(polygon2D, *it1);//append an outer ring to the polygon
+ ::boost::geometry::interior_rings(polygon2D).resize(1);//create an interior ring
+ ::boost::geometry::append(polygon2D, *it2, 0);//append a ring to the interior ring
+ it2IsDonut = true;
+ }
+
+ if (it2IsDonut)
+ {
+ donutIndexVector.push_back(index1);
+ donutIndexVector.push_back(index2);
+ donutVector.push_back(polygon2D);//store donut polygon
+ break;//Only store the first donut!
+ }
+ }
+ }
+ }
+ }
+
+ //Store no donut polygon to boostPolygonVector
+ index1 = 0;
+
+ for (it1 = aRingVector.begin(); it1 != aRingVector.end(); ++it1, index1++)
+ {
+ bool it1IsDonut = false;
+
+ //check if the ring is the outer or inner of a donut
+ for (unsigned int i : donutIndexVector)
+ {
+ if (i == index1)
+ {
+ it1IsDonut = true;
+ break;
+ }
+ }
+
+ if (!it1IsDonut)
+ {
+ BoostMask::BoostPolygon2D polygon2D;
+ ::boost::geometry::append(polygon2D, *it1);
+ boostPolygonVector.push_back(polygon2D);//insert the ring, which is not a part of donut
+ }
+ }
+
+ //Append donut polygon to boostPolygonVector
+ BoostMask::BoostPolygonVector::iterator itDonut;
+
+ for (itDonut = donutVector.begin(); itDonut != donutVector.end(); ++itDonut)
+ {
+ boostPolygonVector.push_back(*itDonut);//append donuts
+ }
+
+ return boostPolygonVector;
+ }
+
+ bool BoostMask::calcVoxelizationThickness(double& aThickness) const
+ {
+
+ if (_voxelizationMap->size() <= 1)
+ {
+ aThickness = 1;
+ return true;
+ }
+
+ double thickness = 0;
+
+ auto it = _voxelizationMap->cbegin();
+ auto it2 = ++_voxelizationMap->cbegin();
+ for (;
+ it != _voxelizationMap->cend() && it2 != _voxelizationMap->cend(); ++it, ++it2)
+ {
+ if (thickness == 0)
+ {
+ thickness = it2->first - it->first;
+ }
+ else
+ {
+ double curThickness = it2->first - it->first;
+ //if not homogeneous (leave out double imprecisions), return false
+ if (std::abs(thickness-curThickness)>errorConstant)
+ {
+ //return false;
+ std::cout << "Two polygons are far from each other?" << std::endl;
+ }
+ }
+
+ }
+
+ if (thickness != 0)
+ {
+ aThickness = thickness;
+ }
+ else
+ {
+ aThickness = 1;
+ }
+
+ return true;
+ }
+ }
+ }
+}
diff --git a/code/masks/boost/rttbBoostMask.h b/code/masks/rttbBoostMask.h
similarity index 97%
rename from code/masks/boost/rttbBoostMask.h
rename to code/masks/rttbBoostMask.h
index 9481c17..cf994bb 100644
--- a/code/masks/boost/rttbBoostMask.h
+++ b/code/masks/rttbBoostMask.h
@@ -1,202 +1,202 @@
-// -----------------------------------------------------------------------
-// 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 __BOOST_MASK_R_H
-#define __BOOST_MASK_R_H
-
-#include "rttbBaseType.h"
-#include "rttbStructure.h"
-#include "rttbGeometricInfo.h"
-#include "rttbMaskAccessorInterface.h"
-
-#include <boost/geometry.hpp>
-#include <boost/geometry/geometries/point_xy.hpp>
-#include <boost/shared_ptr.hpp>
-#include <boost/multi_array.hpp>
-
-namespace rttb
-{
- namespace masks
- {
- namespace boost
- {
- /*! @class BoostMask
- * @brief Implementation of voxelization using boost::geometry.
- * @attention If "strict" is set to true, an exception will be thrown when the given structure has self intersection.
- * (A structure without self intersection means all contours of the structure have no self intersection, and
- * the polygons on the same slice have no intersection between each other, unless the case of a donut. A donut is accepted.)
- * If "strict" is set to false, debug information will be displayed when the given structure has self intersection. Self intersections will be ignored
- * and the mask will be calculated, however, it may cause errors in the mask results.
- */
- class BoostMask
- {
-
- public:
- using MaskVoxelList = core::MaskAccessorInterface::MaskVoxelList;
- using MaskVoxelListPointer = core::MaskAccessorInterface::MaskVoxelListPointer;
-
- /*! @brief Constructor
- * @exception rttb::core::NullPointerException thrown if aDoseGeoInfo or aStructure is nullptr
- * @param aDoseGeoInfo the GeometricInfo
- * @param aStructure the structure set
- * @param strict indicates whether to allow self intersection in the structure. If it is set to true, an exception will be thrown when the given structure has self intersection.
- * @param numberOfThreads number of threads used for voxelization. default value 0 means automatic detection, using the number of Hardware thread/cores
- * @exception InvalidParameterException thrown if strict is true and the structure has self intersections
- */
- BoostMask(core::GeometricInfo::Pointer aDoseGeoInfo, core::Structure::Pointer aStructure,
- bool strict = true, unsigned int numberOfThreads = 0);
-
- /*! @brief Generate mask and return the voxels in the mask
- * @exception rttb::core::InvalidParameterException thrown if the structure has self intersections
- */
- MaskVoxelListPointer getRelevantVoxelVector();
-
- private:
- using BoostPoint2D = ::boost::geometry::model::d2::point_xy<double>;
- using BoostPolygon2D = ::boost::geometry::model::polygon< ::boost::geometry::model::d2::point_xy<double> >;
- using BoostRing2D = ::boost::geometry::model::ring< ::boost::geometry::model::d2::point_xy<double> >;
- using BoostRingVector = std::vector<BoostRing2D>;//polygon without holes
- using BoostPolygonVector = std::vector<BoostPolygon2D>;//polygon with or without holes
- using VoxelIndexVector = std::vector<rttb::VoxelGridIndex3D>;
- typedef std::map<double, BoostPolygonVector>
- BoostPolygonMap;//map of the z index with the vector of boost 2d polygon
- typedef std::map<double, BoostRingVector>
- BoostRingMap;//map of the z index with the vector of boost 2d ring
- typedef ::boost::multi_array<double, 2> BoostArray2D;
- using BoostArray2DPointer = ::boost::shared_ptr<BoostArray2D>;
- typedef ::boost::shared_ptr<std::map<double, BoostArray2DPointer> > BoostArrayMapPointer;
-
- core::GeometricInfo::Pointer _geometricInfo;
-
- core::Structure::Pointer _structure;
-
- bool _strict;
-
- /*! @brief The number of threads
- */
- unsigned int _numberOfThreads;
-
- //@brief The thickness of the voxelization plane (the contour plane), in double dose grid index
- //@details for example, the first contour has the double grid index 0.1, the second 0.3, the third 0.5, then the thickness is 0.2
- double _voxelizationThickness;
-
- //@brief vector of the MaskVoxel inside the structure
- MaskVoxelListPointer _voxelInStructure;
-
- /*! @brief The map of z index and a vector of boost ring 2d (without holes)
- * @details Key: the double z grid index
- * Value: the vector of boost ring 2d (without holes)
- */
- BoostRingMap _ringMap;
-
- /*! @brief The min and max index of the global bounding box.
- * @details The first index has the minimum for x/y/z of the global bounding box.
- * The second index has the maximum for x/y/z of the global bounding index.
- */
- VoxelIndexVector _globalBoundingBox;
-
- /*! @brief The voxelization map
- * @details key: the converted double z grid index of a contour plane
- * value: the 2d mask, array[i][j] = the mask value of the position (i,j) in the global bounding box,
- * i: 0 - (_globalBoundingBoxSize0-1), j: 0 - (_globalBoundingBoxSize1-1)
- */
- BoostArrayMapPointer _voxelizationMap;
-
- /*! @brief If the mask is up to date
- */
- bool _isUpToDate;
-
- /*! @brief Voxelization and generate mask
- */
- void calcMask();
-
- /*! @brief The preprocessing step, wich consists of the following logic and Sub setps:
- * @details For all contours in a struct:
- * 1) Transfer the contour polygons into boost::geometry structures
- * 1a) Convert the contur points from world coordinates into geometry coordinates.
- * 1b) get min and max for x/y/z of a contour
- * 2) Tilt check: if difference of z_min and z_max is larger then a tolerance value -> there is a tilt. Throw rttb::TiltedMaskPlaneException.
- * 3) Get struct-bounding-box: get x_min_struct, y_min_struct, x_max_struct, y_max_struct to define the bounding box that containes all contours of a struct in x-y-dimensions.
- */
- void preprocessing();
-
- /*! @brief The voxelization step, which computes the voxelization planes (in x/y) for all contours of an struct.
-
- * @details For each contour (that is in the z-Range of the reference geometry) of the struct:
- * 1) Allocate result array (voxelization plane) based on the bounding box (see Preprocessing Step 3)
- * 2) Generate voxelization plane for the contour (based on the x-y-raster of the reference geometry).
- * 3) Add result Array (key is the z-Value of the contour)
- */
- void voxelization();
-
- /*! @brief mask voxel Generation step which transfers the voxelization planes into the (z-)geometry of the reference geometry.
- * @details It consists of following Sub steps :
- * For all "slices" in the reference geometry :
- * 1) generate weight vector for all voxelization planes for a given z - value of a slice
- * Iterate over the bounding box of a struct.For each voxel :
- * 2) Compute weighted sum of all voxelization planes(use weight vector, step 1)
- * 2a) If sum > 0 : Add mask voxel for the current x / y(inner Loop) and z value(outer Loop).
- * 3) return mask voxel list.
- */
- void generateMaskVoxelList();
-
- /*! @brief Convert the rttb polygon with world coordinate to the rttb polygon with double geometry coordinate, calculate the current min/max
- * and check if the polygon is planar
- * @param minimum the current global minimum
- * @param maximum the current global maximum
- * @return Return true if the polygon is planar, which means that the minimal and maximal z-coordinate of the polygon is not larger than a error constant
- */
- bool preprocessingPolygon(const rttb::PolygonType& aRTTBPolygon,
- rttb::PolygonType& geometryCoordinatePolygon, rttb::DoubleVoxelGridIndex3D& minimum,
- rttb::DoubleVoxelGridIndex3D& maximum, double aErrorConstant) const;
-
- /*! @brief Convert a rttb 3d polygon to a 2d boost ring*/
- BoostRing2D convertRTTBPolygonToBoostRing(const rttb::PolygonType& aRTTBPolygon) const;
-
- /*! @brief Convert a rttb 3d polygon to a map of z index with a vector of boost 2d ring, because of tilt check use the first z index of the polygon as the map key*/
- BoostRingMap convertRTTBPolygonSequenceToBoostRingMap(const rttb::PolygonSequenceType&
- aRTTBPolygonVector) const;
-
- /*! @brief Find the key with error constant to aIndex
- * @pre aBoostRingMap should not be empty
- * @return Return aBoostRingMap.end() if the key is not found
- */
- BoostMask::BoostRingMap::iterator findNearestKey(BoostMask::BoostRingMap& aBoostRingMap,
- double aIndex, double aErrorConstant) const;
-
- /*! @brief If 2 rings in the vector build a donut, convert the 2 rings to a donut polygon, other rings unchanged*/
- BoostPolygonVector checkDonutAndConvert(const BoostRingVector& aRingVector) const;
-
- /*! @brief Calculate the voxelization thickness.
- Return false, if the voxelization plane is not homogeneous
- */
- bool calcVoxelizationThickness(double& aThickness) const;
-
- };
-
- }
-
-
- }
-}
-
-#endif
+// -----------------------------------------------------------------------
+// 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 __BOOST_MASK_R_H
+#define __BOOST_MASK_R_H
+
+#include "rttbBaseType.h"
+#include "rttbStructure.h"
+#include "rttbGeometricInfo.h"
+#include "rttbMaskAccessorInterface.h"
+
+#include <boost/geometry.hpp>
+#include <boost/geometry/geometries/point_xy.hpp>
+#include <boost/shared_ptr.hpp>
+#include <boost/multi_array.hpp>
+
+namespace rttb
+{
+ namespace masks
+ {
+ namespace boost
+ {
+ /*! @class BoostMask
+ * @brief Implementation of voxelization using boost::geometry.
+ * @attention If "strict" is set to true, an exception will be thrown when the given structure has self intersection.
+ * (A structure without self intersection means all contours of the structure have no self intersection, and
+ * the polygons on the same slice have no intersection between each other, unless the case of a donut. A donut is accepted.)
+ * If "strict" is set to false, debug information will be displayed when the given structure has self intersection. Self intersections will be ignored
+ * and the mask will be calculated, however, it may cause errors in the mask results.
+ */
+ class BoostMask
+ {
+
+ public:
+ using MaskVoxelList = core::MaskAccessorInterface::MaskVoxelList;
+ using MaskVoxelListPointer = core::MaskAccessorInterface::MaskVoxelListPointer;
+
+ /*! @brief Constructor
+ * @exception rttb::core::NullPointerException thrown if aDoseGeoInfo or aStructure is nullptr
+ * @param aDoseGeoInfo the GeometricInfo
+ * @param aStructure the structure set
+ * @param strict indicates whether to allow self intersection in the structure. If it is set to true, an exception will be thrown when the given structure has self intersection.
+ * @param numberOfThreads number of threads used for voxelization. default value 0 means automatic detection, using the number of Hardware thread/cores
+ * @exception InvalidParameterException thrown if strict is true and the structure has self intersections
+ */
+ BoostMask(core::GeometricInfo::Pointer aDoseGeoInfo, core::Structure::Pointer aStructure,
+ bool strict = true, unsigned int numberOfThreads = 0);
+
+ /*! @brief Generate mask and return the voxels in the mask
+ * @exception rttb::core::InvalidParameterException thrown if the structure has self intersections
+ */
+ MaskVoxelListPointer getRelevantVoxelVector();
+
+ private:
+ using BoostPoint2D = ::boost::geometry::model::d2::point_xy<double>;
+ using BoostPolygon2D = ::boost::geometry::model::polygon< ::boost::geometry::model::d2::point_xy<double> >;
+ using BoostRing2D = ::boost::geometry::model::ring< ::boost::geometry::model::d2::point_xy<double> >;
+ using BoostRingVector = std::vector<BoostRing2D>;//polygon without holes
+ using BoostPolygonVector = std::vector<BoostPolygon2D>;//polygon with or without holes
+ using VoxelIndexVector = std::vector<rttb::VoxelGridIndex3D>;
+ typedef std::map<double, BoostPolygonVector>
+ BoostPolygonMap;//map of the z index with the vector of boost 2d polygon
+ typedef std::map<double, BoostRingVector>
+ BoostRingMap;//map of the z index with the vector of boost 2d ring
+ typedef ::boost::multi_array<double, 2> BoostArray2D;
+ using BoostArray2DPointer = ::boost::shared_ptr<BoostArray2D>;
+ typedef ::boost::shared_ptr<std::map<double, BoostArray2DPointer> > BoostArrayMapPointer;
+
+ core::GeometricInfo::Pointer _geometricInfo;
+
+ core::Structure::Pointer _structure;
+
+ bool _strict;
+
+ /*! @brief The number of threads
+ */
+ unsigned int _numberOfThreads;
+
+ //@brief The thickness of the voxelization plane (the contour plane), in double dose grid index
+ //@details for example, the first contour has the double grid index 0.1, the second 0.3, the third 0.5, then the thickness is 0.2
+ double _voxelizationThickness;
+
+ //@brief vector of the MaskVoxel inside the structure
+ MaskVoxelListPointer _voxelInStructure;
+
+ /*! @brief The map of z index and a vector of boost ring 2d (without holes)
+ * @details Key: the double z grid index
+ * Value: the vector of boost ring 2d (without holes)
+ */
+ BoostRingMap _ringMap;
+
+ /*! @brief The min and max index of the global bounding box.
+ * @details The first index has the minimum for x/y/z of the global bounding box.
+ * The second index has the maximum for x/y/z of the global bounding index.
+ */
+ VoxelIndexVector _globalBoundingBox;
+
+ /*! @brief The voxelization map
+ * @details key: the converted double z grid index of a contour plane
+ * value: the 2d mask, array[i][j] = the mask value of the position (i,j) in the global bounding box,
+ * i: 0 - (_globalBoundingBoxSize0-1), j: 0 - (_globalBoundingBoxSize1-1)
+ */
+ BoostArrayMapPointer _voxelizationMap;
+
+ /*! @brief If the mask is up to date
+ */
+ bool _isUpToDate;
+
+ /*! @brief Voxelization and generate mask
+ */
+ void calcMask();
+
+ /*! @brief The preprocessing step, wich consists of the following logic and Sub setps:
+ * @details For all contours in a struct:
+ * 1) Transfer the contour polygons into boost::geometry structures
+ * 1a) Convert the contur points from world coordinates into geometry coordinates.
+ * 1b) get min and max for x/y/z of a contour
+ * 2) Tilt check: if difference of z_min and z_max is larger then a tolerance value -> there is a tilt. Throw rttb::TiltedMaskPlaneException.
+ * 3) Get struct-bounding-box: get x_min_struct, y_min_struct, x_max_struct, y_max_struct to define the bounding box that containes all contours of a struct in x-y-dimensions.
+ */
+ void preprocessing();
+
+ /*! @brief The voxelization step, which computes the voxelization planes (in x/y) for all contours of an struct.
+
+ * @details For each contour (that is in the z-Range of the reference geometry) of the struct:
+ * 1) Allocate result array (voxelization plane) based on the bounding box (see Preprocessing Step 3)
+ * 2) Generate voxelization plane for the contour (based on the x-y-raster of the reference geometry).
+ * 3) Add result Array (key is the z-Value of the contour)
+ */
+ void voxelization();
+
+ /*! @brief mask voxel Generation step which transfers the voxelization planes into the (z-)geometry of the reference geometry.
+ * @details It consists of following Sub steps :
+ * For all "slices" in the reference geometry :
+ * 1) generate weight vector for all voxelization planes for a given z - value of a slice
+ * Iterate over the bounding box of a struct.For each voxel :
+ * 2) Compute weighted sum of all voxelization planes(use weight vector, step 1)
+ * 2a) If sum > 0 : Add mask voxel for the current x / y(inner Loop) and z value(outer Loop).
+ * 3) return mask voxel list.
+ */
+ void generateMaskVoxelList();
+
+ /*! @brief Convert the rttb polygon with world coordinate to the rttb polygon with double geometry coordinate, calculate the current min/max
+ * and check if the polygon is planar
+ * @param minimum the current global minimum
+ * @param maximum the current global maximum
+ * @return Return true if the polygon is planar, which means that the minimal and maximal z-coordinate of the polygon is not larger than a error constant
+ */
+ bool preprocessingPolygon(const rttb::PolygonType& aRTTBPolygon,
+ rttb::PolygonType& geometryCoordinatePolygon, rttb::DoubleVoxelGridIndex3D& minimum,
+ rttb::DoubleVoxelGridIndex3D& maximum, double aErrorConstant) const;
+
+ /*! @brief Convert a rttb 3d polygon to a 2d boost ring*/
+ BoostRing2D convertRTTBPolygonToBoostRing(const rttb::PolygonType& aRTTBPolygon) const;
+
+ /*! @brief Convert a rttb 3d polygon to a map of z index with a vector of boost 2d ring, because of tilt check use the first z index of the polygon as the map key*/
+ BoostRingMap convertRTTBPolygonSequenceToBoostRingMap(const rttb::PolygonSequenceType&
+ aRTTBPolygonVector) const;
+
+ /*! @brief Find the key with error constant to aIndex
+ * @pre aBoostRingMap should not be empty
+ * @return Return aBoostRingMap.end() if the key is not found
+ */
+ BoostMask::BoostRingMap::iterator findNearestKey(BoostMask::BoostRingMap& aBoostRingMap,
+ double aIndex, double aErrorConstant) const;
+
+ /*! @brief If 2 rings in the vector build a donut, convert the 2 rings to a donut polygon, other rings unchanged*/
+ BoostPolygonVector checkDonutAndConvert(const BoostRingVector& aRingVector) const;
+
+ /*! @brief Calculate the voxelization thickness.
+ Return false, if the voxelization plane is not homogeneous
+ */
+ bool calcVoxelizationThickness(double& aThickness) const;
+
+ };
+
+ }
+
+
+ }
+}
+
+#endif
diff --git a/code/masks/boost/rttbBoostMaskAccessor.cpp b/code/masks/rttbBoostMaskAccessor.cpp
similarity index 96%
rename from code/masks/boost/rttbBoostMaskAccessor.cpp
rename to code/masks/rttbBoostMaskAccessor.cpp
index 61d2525..0985637 100644
--- a/code/masks/boost/rttbBoostMaskAccessor.cpp
+++ b/code/masks/rttbBoostMaskAccessor.cpp
@@ -1,161 +1,161 @@
-// -----------------------------------------------------------------------
-// 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 "rttbBoostMaskAccessor.h"
-#include "rttbBoostMask.h"
-
-#include <boost/make_shared.hpp>
-
-#include <boost/uuid/uuid.hpp>
-#include <boost/uuid/uuid_generators.hpp>
-#include <boost/uuid/uuid_io.hpp>
-
-namespace rttb
-{
- namespace masks
- {
- namespace boost
- {
-
- BoostMaskAccessor::BoostMaskAccessor(StructTypePointer aStructurePointer,
- const core::GeometricInfo& aGeometricInfo, bool strict)
- : _spStructure(aStructurePointer), _geoInfo(aGeometricInfo), _strict(strict)
- {
- _spRelevantVoxelVector = MaskVoxelListPointer();
-
- //generate new structure set uid
- ::boost::uuids::uuid id;
- ::boost::uuids::random_generator generator;
- id = generator();
- std::stringstream ss;
- ss << id;
- _maskUID = "BoostMask_" + ss.str();
- }
-
-
- BoostMaskAccessor::~BoostMaskAccessor()
- = default;
-
- void BoostMaskAccessor::updateMask()
- {
- MaskVoxelList newRelevantVoxelVector;
-
- if (_spRelevantVoxelVector)
- {
- return; // already calculated
- }
-
- BoostMask mask(::boost::make_shared<core::GeometricInfo>(_geoInfo),
- _spStructure, _strict);
-
- _spRelevantVoxelVector = mask.getRelevantVoxelVector();
- }
-
- BoostMaskAccessor::MaskVoxelListPointer BoostMaskAccessor::getRelevantVoxelVector()
- {
- // if not already generated start voxelization here
- updateMask();
- return _spRelevantVoxelVector;
- }
-
- BoostMaskAccessor::MaskVoxelListPointer BoostMaskAccessor::getRelevantVoxelVector(
- float lowerThreshold)
- {
- auto filteredVoxelVectorPointer = ::boost::make_shared<MaskVoxelList>();
- updateMask();
- // filter relevant voxels
- auto it = _spRelevantVoxelVector->begin();
-
- while (it != _spRelevantVoxelVector->end())
- {
- if ((*it).getRelevantVolumeFraction() > lowerThreshold)
- {
- filteredVoxelVectorPointer->push_back(*it);
- }
-
- ++it;
- }
-
- // if mask calculation was not successful this is empty!
- return filteredVoxelVectorPointer;
- }
-
- bool BoostMaskAccessor::getMaskAt(VoxelGridID aID, core::MaskVoxel& voxel) const
- {
- //initialize return voxel
- voxel.setRelevantVolumeFraction(0);
-
- //check if ID is valid
- if (!_geoInfo.validID(aID))
- {
- return false;
- }
-
- //determine how a given voxel on the dose grid is masked
- if (_spRelevantVoxelVector)
- {
- auto it = _spRelevantVoxelVector->begin();
-
- while (it != _spRelevantVoxelVector->end())
- {
- if ((*it).getVoxelGridID() == aID)
- {
- voxel = (*it);
- return true;
- }
-
- ++it;
- }
-
- }
- // returns false if mask was not calculated without triggering calculation (otherwise not const!)
- else
- {
- return false;
- }
-
- return false;
-
- }
-
- bool BoostMaskAccessor::getMaskAt(const VoxelGridIndex3D& aIndex, core::MaskVoxel& voxel) const
- {
- //convert VoxelGridIndex3D to VoxelGridID
- VoxelGridID aVoxelGridID;
-
- if (_geoInfo.convert(aIndex, aVoxelGridID))
- {
- return getMaskAt(aVoxelGridID, voxel);
- }
- else
- {
- return false;
- }
- }
-
- const core::GeometricInfo& BoostMaskAccessor::getGeometricInfo() const
- {
- return _geoInfo;
- };
- }
-
- }
+// -----------------------------------------------------------------------
+// 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 "rttbBoostMaskAccessor.h"
+#include "rttbBoostMask.h"
+
+#include <boost/make_shared.hpp>
+
+#include <boost/uuid/uuid.hpp>
+#include <boost/uuid/uuid_generators.hpp>
+#include <boost/uuid/uuid_io.hpp>
+
+namespace rttb
+{
+ namespace masks
+ {
+ namespace boost
+ {
+
+ BoostMaskAccessor::BoostMaskAccessor(StructTypePointer aStructurePointer,
+ const core::GeometricInfo& aGeometricInfo, bool strict)
+ : _spStructure(aStructurePointer), _geoInfo(aGeometricInfo), _strict(strict)
+ {
+ _spRelevantVoxelVector = MaskVoxelListPointer();
+
+ //generate new structure set uid
+ ::boost::uuids::uuid id;
+ ::boost::uuids::random_generator generator;
+ id = generator();
+ std::stringstream ss;
+ ss << id;
+ _maskUID = "BoostMask_" + ss.str();
+ }
+
+
+ BoostMaskAccessor::~BoostMaskAccessor()
+ = default;
+
+ void BoostMaskAccessor::updateMask()
+ {
+ MaskVoxelList newRelevantVoxelVector;
+
+ if (_spRelevantVoxelVector)
+ {
+ return; // already calculated
+ }
+
+ BoostMask mask(::boost::make_shared<core::GeometricInfo>(_geoInfo),
+ _spStructure, _strict);
+
+ _spRelevantVoxelVector = mask.getRelevantVoxelVector();
+ }
+
+ BoostMaskAccessor::MaskVoxelListPointer BoostMaskAccessor::getRelevantVoxelVector()
+ {
+ // if not already generated start voxelization here
+ updateMask();
+ return _spRelevantVoxelVector;
+ }
+
+ BoostMaskAccessor::MaskVoxelListPointer BoostMaskAccessor::getRelevantVoxelVector(
+ float lowerThreshold)
+ {
+ auto filteredVoxelVectorPointer = ::boost::make_shared<MaskVoxelList>();
+ updateMask();
+ // filter relevant voxels
+ auto it = _spRelevantVoxelVector->begin();
+
+ while (it != _spRelevantVoxelVector->end())
+ {
+ if ((*it).getRelevantVolumeFraction() > lowerThreshold)
+ {
+ filteredVoxelVectorPointer->push_back(*it);
+ }
+
+ ++it;
+ }
+
+ // if mask calculation was not successful this is empty!
+ return filteredVoxelVectorPointer;
+ }
+
+ bool BoostMaskAccessor::getMaskAt(VoxelGridID aID, core::MaskVoxel& voxel) const
+ {
+ //initialize return voxel
+ voxel.setRelevantVolumeFraction(0);
+
+ //check if ID is valid
+ if (!_geoInfo.validID(aID))
+ {
+ return false;
+ }
+
+ //determine how a given voxel on the dose grid is masked
+ if (_spRelevantVoxelVector)
+ {
+ auto it = _spRelevantVoxelVector->begin();
+
+ while (it != _spRelevantVoxelVector->end())
+ {
+ if ((*it).getVoxelGridID() == aID)
+ {
+ voxel = (*it);
+ return true;
+ }
+
+ ++it;
+ }
+
+ }
+ // returns false if mask was not calculated without triggering calculation (otherwise not const!)
+ else
+ {
+ return false;
+ }
+
+ return false;
+
+ }
+
+ bool BoostMaskAccessor::getMaskAt(const VoxelGridIndex3D& aIndex, core::MaskVoxel& voxel) const
+ {
+ //convert VoxelGridIndex3D to VoxelGridID
+ VoxelGridID aVoxelGridID;
+
+ if (_geoInfo.convert(aIndex, aVoxelGridID))
+ {
+ return getMaskAt(aVoxelGridID, voxel);
+ }
+ else
+ {
+ return false;
+ }
+ }
+
+ const core::GeometricInfo& BoostMaskAccessor::getGeometricInfo() const
+ {
+ return _geoInfo;
+ };
+ }
+
+ }
}
\ No newline at end of file
diff --git a/code/masks/boost/rttbBoostMaskAccessor.h b/code/masks/rttbBoostMaskAccessor.h
similarity index 97%
rename from code/masks/boost/rttbBoostMaskAccessor.h
rename to code/masks/rttbBoostMaskAccessor.h
index 759b50c..83bf37c 100644
--- a/code/masks/boost/rttbBoostMaskAccessor.h
+++ b/code/masks/rttbBoostMaskAccessor.h
@@ -1,121 +1,121 @@
-// -----------------------------------------------------------------------
-// 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 __BOOST_MASK_R_ACCESSOR__H
-#define __BOOST_MASK_R_ACCESSOR__H
-
-#include "rttbBaseType.h"
-#include "rttbGeometricInfo.h"
-#include "rttbMaskAccessorInterface.h"
-#include "rttbStructure.h"
-
-#include "RTTBBoostMaskExports.h"
-
-namespace rttb
-{
- namespace masks
- {
- namespace boost
- {
- /*! @class BoostMaskAccessor
- * @brief Using the voxelization based on boost::geometry and generate the mask accessor.
- * @attention If "strict" is set to true, an exception will be thrown when the given structure has self intersection.
- * (A structure without self intersection means all contours of the structure have no self intersection, and
- * the polygons on the same slice have no intersection between each other, unless the case of a donut. A donut is accepted.)
- * If "strict" is set to false, debug information will be displayed when the given structure has self intersection. Self intersections will be ignored
- * and the mask will be calculated, however, it may cause errors in the mask results.
- */
- class RTTBBoostMask_EXPORT BoostMaskAccessor : public core::MaskAccessorInterface
- {
- public:
- using MaskVoxelList = core::MaskAccessorInterface::MaskVoxelList;
- using MaskVoxelListPointer = core::MaskAccessorInterface::MaskVoxelListPointer;
-
- using StructTypePointer = core::Structure::Pointer;
-
- private:
- StructTypePointer _spStructure;
- core::GeometricInfo _geoInfo;
- bool _strict;
-
- /*! vector containing list of mask voxels*/
- MaskVoxelListPointer _spRelevantVoxelVector;
-
- IDType _maskUID;
-
- public:
-
-
- /*! @brief Constructor with a structure pointer and a geometric info pointer
- * @param aStructurePointer smart pointer of the structure
- * @param aGeometricInfo smart pointer of the geometricInfo of the dose
- * @param strict indicates whether to allow self intersection in the structure. If it is set to true, an exception will be thrown when the given structure has self intersection.
- * @exception InvalidParameterException thrown if strict is true and the structure has self intersections
- */
- BoostMaskAccessor(StructTypePointer aStructurePointer, const core::GeometricInfo& aGeometricInfo,
- bool strict = true);
-
- /*! @brief destructor*/
- ~BoostMaskAccessor() override;
-
- /*! @brief voxelization of the given structures using boost algorithms*/
- void updateMask() override;
-
- /*! @brief get vector containing all relevant voxels that are inside the given structure*/
- MaskVoxelListPointer getRelevantVoxelVector() override;
-
- /*! @brief get vector containing all relevant voxels that have a relevant volume above the given threshold and are inside the given structure*/
- MaskVoxelListPointer getRelevantVoxelVector(float lowerThreshold) override;
-
- /*!@brief determine how a given voxel on the dose grid is masked
- * @param aID ID of the voxel in grid.
- * @param voxel Reference to the voxel.
- * @post after a valid call voxel containes the information of the specified grid voxel. If aID is not valid, voxel values are undefined.
- * The relevant volume fraction will be set to zero.
- * @return Indicates of the voxel exists and therefore if parameter voxel containes valid values.*/
- bool getMaskAt(const VoxelGridID aID, core::MaskVoxel& voxel) const override;
-
- /*!@brief determine how a given voxel on the dose grid is masked
- * @param aIndex 3d index of the voxel in grid.
- * @param voxel Reference to the voxel.
- * @return Indicates of the voxel exists and therefore if parameter voxel containes valid values.*/
- bool getMaskAt(const VoxelGridIndex3D& aIndex, core::MaskVoxel& voxel) const override;
-
- /*! @brief give access to GeometricInfo*/
- const core::GeometricInfo& getGeometricInfo() const override;
-
- /* @ brief is true if dose is on a homogeneous grid
- * @remark Inhomogeneous grids are not supported at the moment, but if they will be supported in the future the interface does not need to change.*/
- bool isGridHomogeneous() const override
- {
- return true;
- };
-
- IDType getMaskUID() const override
- {
- return _maskUID;
- };
-
- };
- }
- }
-}
-
-#endif
+// -----------------------------------------------------------------------
+// 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 __BOOST_MASK_R_ACCESSOR__H
+#define __BOOST_MASK_R_ACCESSOR__H
+
+#include "rttbBaseType.h"
+#include "rttbGeometricInfo.h"
+#include "rttbMaskAccessorInterface.h"
+#include "rttbStructure.h"
+
+#include "RTTBBoostMaskExports.h"
+
+namespace rttb
+{
+ namespace masks
+ {
+ namespace boost
+ {
+ /*! @class BoostMaskAccessor
+ * @brief Using the voxelization based on boost::geometry and generate the mask accessor.
+ * @attention If "strict" is set to true, an exception will be thrown when the given structure has self intersection.
+ * (A structure without self intersection means all contours of the structure have no self intersection, and
+ * the polygons on the same slice have no intersection between each other, unless the case of a donut. A donut is accepted.)
+ * If "strict" is set to false, debug information will be displayed when the given structure has self intersection. Self intersections will be ignored
+ * and the mask will be calculated, however, it may cause errors in the mask results.
+ */
+ class RTTBBoostMask_EXPORT BoostMaskAccessor : public core::MaskAccessorInterface
+ {
+ public:
+ using MaskVoxelList = core::MaskAccessorInterface::MaskVoxelList;
+ using MaskVoxelListPointer = core::MaskAccessorInterface::MaskVoxelListPointer;
+
+ using StructTypePointer = core::Structure::Pointer;
+
+ private:
+ StructTypePointer _spStructure;
+ core::GeometricInfo _geoInfo;
+ bool _strict;
+
+ /*! vector containing list of mask voxels*/
+ MaskVoxelListPointer _spRelevantVoxelVector;
+
+ IDType _maskUID;
+
+ public:
+
+
+ /*! @brief Constructor with a structure pointer and a geometric info pointer
+ * @param aStructurePointer smart pointer of the structure
+ * @param aGeometricInfo smart pointer of the geometricInfo of the dose
+ * @param strict indicates whether to allow self intersection in the structure. If it is set to true, an exception will be thrown when the given structure has self intersection.
+ * @exception InvalidParameterException thrown if strict is true and the structure has self intersections
+ */
+ BoostMaskAccessor(StructTypePointer aStructurePointer, const core::GeometricInfo& aGeometricInfo,
+ bool strict = true);
+
+ /*! @brief destructor*/
+ ~BoostMaskAccessor() override;
+
+ /*! @brief voxelization of the given structures using boost algorithms*/
+ void updateMask() override;
+
+ /*! @brief get vector containing all relevant voxels that are inside the given structure*/
+ MaskVoxelListPointer getRelevantVoxelVector() override;
+
+ /*! @brief get vector containing all relevant voxels that have a relevant volume above the given threshold and are inside the given structure*/
+ MaskVoxelListPointer getRelevantVoxelVector(float lowerThreshold) override;
+
+ /*!@brief determine how a given voxel on the dose grid is masked
+ * @param aID ID of the voxel in grid.
+ * @param voxel Reference to the voxel.
+ * @post after a valid call voxel containes the information of the specified grid voxel. If aID is not valid, voxel values are undefined.
+ * The relevant volume fraction will be set to zero.
+ * @return Indicates of the voxel exists and therefore if parameter voxel containes valid values.*/
+ bool getMaskAt(const VoxelGridID aID, core::MaskVoxel& voxel) const override;
+
+ /*!@brief determine how a given voxel on the dose grid is masked
+ * @param aIndex 3d index of the voxel in grid.
+ * @param voxel Reference to the voxel.
+ * @return Indicates of the voxel exists and therefore if parameter voxel containes valid values.*/
+ bool getMaskAt(const VoxelGridIndex3D& aIndex, core::MaskVoxel& voxel) const override;
+
+ /*! @brief give access to GeometricInfo*/
+ const core::GeometricInfo& getGeometricInfo() const override;
+
+ /* @ brief is true if dose is on a homogeneous grid
+ * @remark Inhomogeneous grids are not supported at the moment, but if they will be supported in the future the interface does not need to change.*/
+ bool isGridHomogeneous() const override
+ {
+ return true;
+ };
+
+ IDType getMaskUID() const override
+ {
+ return _maskUID;
+ };
+
+ };
+ }
+ }
+}
+
+#endif
diff --git a/code/masks/boost/rttbBoostMaskGenerateMaskVoxelListThread.cpp b/code/masks/rttbBoostMaskGenerateMaskVoxelListThread.cpp
similarity index 97%
rename from code/masks/boost/rttbBoostMaskGenerateMaskVoxelListThread.cpp
rename to code/masks/rttbBoostMaskGenerateMaskVoxelListThread.cpp
index cf0f169..02c39ce 100644
--- a/code/masks/boost/rttbBoostMaskGenerateMaskVoxelListThread.cpp
+++ b/code/masks/rttbBoostMaskGenerateMaskVoxelListThread.cpp
@@ -1,156 +1,156 @@
-// -----------------------------------------------------------------------
-// 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: 1127 $ (last changed revision)
-// @date $Date: 2015-10-01 13:33:33 +0200 (Do, 01 Okt 2015) $ (last change date)
-// @author $Author: hentsch $ (last changed by)
-*/
-
-#include "rttbBoostMaskGenerateMaskVoxelListThread.h"
-
-#include "rttbInvalidParameterException.h"
-#include <boost/thread.hpp>
-
-
-namespace rttb
-{
- namespace masks
- {
- namespace boost
- {
- BoostMaskGenerateMaskVoxelListThread::BoostMaskGenerateMaskVoxelListThread(
- const VoxelIndexVector& aGlobalBoundingBox,
- core::GeometricInfo::Pointer aGeometricInfo,
- BoostArrayMapPointer aVoxelizationMap,
- double aVoxelizationThickness,
- unsigned int aBeginSlice,
- unsigned int aEndSlice,
- MaskVoxelListPointer aMaskVoxelList,
- ::boost::shared_ptr<::boost::shared_mutex> aMutex) :
- _globalBoundingBox(aGlobalBoundingBox), _geometricInfo(aGeometricInfo),
- _voxelizationMap(aVoxelizationMap), _voxelizationThickness(aVoxelizationThickness),
- _beginSlice(aBeginSlice), _endSlice(aEndSlice),
- _resultMaskVoxelList(aMaskVoxelList), _mutex(aMutex)
- {}
-
- void BoostMaskGenerateMaskVoxelListThread::operator()()
- {
- rttb::VoxelGridIndex3D minIndex = _globalBoundingBox.at(0);
- rttb::VoxelGridIndex3D maxIndex = _globalBoundingBox.at(1);
- unsigned int globalBoundingBoxSize0 = maxIndex[0] - minIndex[0] + 1;
- unsigned int globalBoundingBoxSize1 = maxIndex[1] - minIndex[1] + 1;
-
- std::vector<core::MaskVoxel> maskVoxelsInThread;
-
- for (unsigned int indexZ = _beginSlice; indexZ < _endSlice; ++indexZ)
- {
- //calculate weight vector
- std::map<double, double> weightVectorForZ;
- calcWeightVector(indexZ, weightVectorForZ);
-
- //For each x,y, calc sum of all voxelization plane, use weight vector
- for (unsigned int x = 0; x < globalBoundingBoxSize0; ++x)
- {
- for (unsigned int y = 0; y < globalBoundingBoxSize1; ++y)
- {
- rttb::VoxelGridIndex3D currentIndex;
- currentIndex[0] = x + minIndex[0];
- currentIndex[1] = y + minIndex[1];
- currentIndex[2] = indexZ;
- rttb::VoxelGridID gridID;
- _geometricInfo->convert(currentIndex, gridID);
- double volumeFraction = 0;
-
- auto it = weightVectorForZ.cbegin();
- auto itMap = _voxelizationMap->cbegin();
-
- for (; it != weightVectorForZ.cend()
- && itMap != _voxelizationMap->cend(); ++it, ++itMap)
- {
- double weight = it->second;
- if (weight > 0){
- BoostArray2DPointer voxelizationArray = itMap->second;
-
- //calc sum of all voxelization plane, use weight
- volumeFraction += (*voxelizationArray)[x][y] * weight;
- }
- }
-
- if (volumeFraction > 1 && (volumeFraction - 1) <= errorConstant)
- {
- volumeFraction = 1;
- }
- else if (volumeFraction < 0 || (volumeFraction - 1) > errorConstant)
- {
- throw rttb::core::InvalidParameterException("Mask calculation failed! The volume fraction should >= 0 and <= 1!");
- }
-
- //insert mask voxel if volumeFraction > 0
- if (volumeFraction > 0)
- {
- core::MaskVoxel maskVoxelPtr = core::MaskVoxel(gridID, volumeFraction);
- maskVoxelsInThread.push_back(maskVoxelPtr);
- }
- }
-
- }
- }
-
- ::boost::unique_lock<::boost::shared_mutex> lock(*_mutex);
- _resultMaskVoxelList->insert(_resultMaskVoxelList->end(), maskVoxelsInThread.begin(), maskVoxelsInThread.end());
- }
-
- void BoostMaskGenerateMaskVoxelListThread::calcWeightVector(const rttb::VoxelGridID& aIndexZ,
- std::map<double, double>& weightVector) const
- {
- double indexZMin = aIndexZ - 0.5;
- double indexZMax = aIndexZ + 0.5;
-
- for (auto & it : *_voxelizationMap)
- {
- double voxelizationPlaneIndexMin = it.first - 0.5 * _voxelizationThickness;
- double voxelizationPlaneIndexMax = it.first + 0.5 * _voxelizationThickness;
- double weight = 0;
-
- if ((voxelizationPlaneIndexMin < indexZMin) && (voxelizationPlaneIndexMax > indexZMin))
- {
- if (voxelizationPlaneIndexMax < indexZMax)
- {
- weight = voxelizationPlaneIndexMax - indexZMin;
- }
- else
- {
- weight = 1;
- }
- }
- else if ((voxelizationPlaneIndexMin >= indexZMin) && (voxelizationPlaneIndexMin < indexZMax))
- {
- if (voxelizationPlaneIndexMax < indexZMax)
- {
- weight = _voxelizationThickness;
- }
- else
- {
- weight = indexZMax - voxelizationPlaneIndexMin;
- }
- }
-
- weightVector.insert(std::pair<double, double>(it.first, weight));
- }
- }
- }
- }
-}
+// -----------------------------------------------------------------------
+// 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: 1127 $ (last changed revision)
+// @date $Date: 2015-10-01 13:33:33 +0200 (Do, 01 Okt 2015) $ (last change date)
+// @author $Author: hentsch $ (last changed by)
+*/
+
+#include "rttbBoostMaskGenerateMaskVoxelListThread.h"
+
+#include "rttbInvalidParameterException.h"
+#include <boost/thread.hpp>
+
+
+namespace rttb
+{
+ namespace masks
+ {
+ namespace boost
+ {
+ BoostMaskGenerateMaskVoxelListThread::BoostMaskGenerateMaskVoxelListThread(
+ const VoxelIndexVector& aGlobalBoundingBox,
+ core::GeometricInfo::Pointer aGeometricInfo,
+ BoostArrayMapPointer aVoxelizationMap,
+ double aVoxelizationThickness,
+ unsigned int aBeginSlice,
+ unsigned int aEndSlice,
+ MaskVoxelListPointer aMaskVoxelList,
+ ::boost::shared_ptr<::boost::shared_mutex> aMutex) :
+ _globalBoundingBox(aGlobalBoundingBox), _geometricInfo(aGeometricInfo),
+ _voxelizationMap(aVoxelizationMap), _voxelizationThickness(aVoxelizationThickness),
+ _beginSlice(aBeginSlice), _endSlice(aEndSlice),
+ _resultMaskVoxelList(aMaskVoxelList), _mutex(aMutex)
+ {}
+
+ void BoostMaskGenerateMaskVoxelListThread::operator()()
+ {
+ rttb::VoxelGridIndex3D minIndex = _globalBoundingBox.at(0);
+ rttb::VoxelGridIndex3D maxIndex = _globalBoundingBox.at(1);
+ unsigned int globalBoundingBoxSize0 = maxIndex[0] - minIndex[0] + 1;
+ unsigned int globalBoundingBoxSize1 = maxIndex[1] - minIndex[1] + 1;
+
+ std::vector<core::MaskVoxel> maskVoxelsInThread;
+
+ for (unsigned int indexZ = _beginSlice; indexZ < _endSlice; ++indexZ)
+ {
+ //calculate weight vector
+ std::map<double, double> weightVectorForZ;
+ calcWeightVector(indexZ, weightVectorForZ);
+
+ //For each x,y, calc sum of all voxelization plane, use weight vector
+ for (unsigned int x = 0; x < globalBoundingBoxSize0; ++x)
+ {
+ for (unsigned int y = 0; y < globalBoundingBoxSize1; ++y)
+ {
+ rttb::VoxelGridIndex3D currentIndex;
+ currentIndex[0] = x + minIndex[0];
+ currentIndex[1] = y + minIndex[1];
+ currentIndex[2] = indexZ;
+ rttb::VoxelGridID gridID;
+ _geometricInfo->convert(currentIndex, gridID);
+ double volumeFraction = 0;
+
+ auto it = weightVectorForZ.cbegin();
+ auto itMap = _voxelizationMap->cbegin();
+
+ for (; it != weightVectorForZ.cend()
+ && itMap != _voxelizationMap->cend(); ++it, ++itMap)
+ {
+ double weight = it->second;
+ if (weight > 0){
+ BoostArray2DPointer voxelizationArray = itMap->second;
+
+ //calc sum of all voxelization plane, use weight
+ volumeFraction += (*voxelizationArray)[x][y] * weight;
+ }
+ }
+
+ if (volumeFraction > 1 && (volumeFraction - 1) <= errorConstant)
+ {
+ volumeFraction = 1;
+ }
+ else if (volumeFraction < 0 || (volumeFraction - 1) > errorConstant)
+ {
+ throw rttb::core::InvalidParameterException("Mask calculation failed! The volume fraction should >= 0 and <= 1!");
+ }
+
+ //insert mask voxel if volumeFraction > 0
+ if (volumeFraction > 0)
+ {
+ core::MaskVoxel maskVoxelPtr = core::MaskVoxel(gridID, volumeFraction);
+ maskVoxelsInThread.push_back(maskVoxelPtr);
+ }
+ }
+
+ }
+ }
+
+ ::boost::unique_lock<::boost::shared_mutex> lock(*_mutex);
+ _resultMaskVoxelList->insert(_resultMaskVoxelList->end(), maskVoxelsInThread.begin(), maskVoxelsInThread.end());
+ }
+
+ void BoostMaskGenerateMaskVoxelListThread::calcWeightVector(const rttb::VoxelGridID& aIndexZ,
+ std::map<double, double>& weightVector) const
+ {
+ double indexZMin = aIndexZ - 0.5;
+ double indexZMax = aIndexZ + 0.5;
+
+ for (auto & it : *_voxelizationMap)
+ {
+ double voxelizationPlaneIndexMin = it.first - 0.5 * _voxelizationThickness;
+ double voxelizationPlaneIndexMax = it.first + 0.5 * _voxelizationThickness;
+ double weight = 0;
+
+ if ((voxelizationPlaneIndexMin < indexZMin) && (voxelizationPlaneIndexMax > indexZMin))
+ {
+ if (voxelizationPlaneIndexMax < indexZMax)
+ {
+ weight = voxelizationPlaneIndexMax - indexZMin;
+ }
+ else
+ {
+ weight = 1;
+ }
+ }
+ else if ((voxelizationPlaneIndexMin >= indexZMin) && (voxelizationPlaneIndexMin < indexZMax))
+ {
+ if (voxelizationPlaneIndexMax < indexZMax)
+ {
+ weight = _voxelizationThickness;
+ }
+ else
+ {
+ weight = indexZMax - voxelizationPlaneIndexMin;
+ }
+ }
+
+ weightVector.insert(std::pair<double, double>(it.first, weight));
+ }
+ }
+ }
+ }
+}
diff --git a/code/masks/boost/rttbBoostMaskGenerateMaskVoxelListThread.h b/code/masks/rttbBoostMaskGenerateMaskVoxelListThread.h
similarity index 97%
rename from code/masks/boost/rttbBoostMaskGenerateMaskVoxelListThread.h
rename to code/masks/rttbBoostMaskGenerateMaskVoxelListThread.h
index b77f70e..76dc429 100644
--- a/code/masks/boost/rttbBoostMaskGenerateMaskVoxelListThread.h
+++ b/code/masks/rttbBoostMaskGenerateMaskVoxelListThread.h
@@ -1,89 +1,89 @@
-// -----------------------------------------------------------------------
-// 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: 1127 $ (last changed revision)
-// @date $Date: 2015-10-01 13:33:33 +0200 (Do, 01 Okt 2015) $ (last change date)
-// @author $Author: hentsch $ (last changed by)
-*/
-
-
-#ifndef __BOOST_MASK_GENERATE_MASK_VOXEL_LIST_H
-#define __BOOST_MASK_GENERATE_MASK_VOXEL_LIST_H
-
-#include "rttbBaseType.h"
-#include "rttbGeometricInfo.h"
-#include "rttbMaskAccessorInterface.h"
-
-#include <boost/multi_array.hpp>
-#include <boost/shared_ptr.hpp>
-#include <boost/thread/shared_mutex.hpp>
-
-namespace rttb
-{
- namespace masks
- {
- namespace boost
- {
- /*! @class BoostMaskGenerateMaskVoxelListThread
- *
- */
- class BoostMaskGenerateMaskVoxelListThread
- {
-
- public:
- typedef ::boost::multi_array<double, 2> BoostArray2D;
- using BoostArray2DPointer = ::boost::shared_ptr<BoostArray2D>;
- typedef ::boost::shared_ptr<std::map<double, BoostArray2DPointer> > BoostArrayMapPointer;
- using VoxelIndexVector = std::vector<rttb::VoxelGridIndex3D>;
- using MaskVoxelListPointer = core::MaskAccessorInterface::MaskVoxelListPointer;
-
- BoostMaskGenerateMaskVoxelListThread(const VoxelIndexVector& aGlobalBoundingBox,
- core::GeometricInfo::Pointer aGeometricInfo,
- BoostArrayMapPointer aVoxelizationMap,
- double aVoxelizationThickness,
- unsigned int aBeginSlice,
- unsigned int aEndSlice,
- MaskVoxelListPointer aMaskVoxelList, ::boost::shared_ptr<::boost::shared_mutex> aMutex);
- void operator()();
-
- private:
- VoxelIndexVector _globalBoundingBox;
- core::GeometricInfo::Pointer _geometricInfo;
- BoostArrayMapPointer _voxelizationMap;
- //(for example, the first contour has the double grid index 0.1, the second 0.3, the third 0.5, then the thickness is 0.2)
- double _voxelizationThickness;
-
- unsigned int _beginSlice;
- /*! @brief _endSlice is the first index not to be processed (like a end iterator)
- */
- unsigned int _endSlice;
-
- MaskVoxelListPointer _resultMaskVoxelList;
- ::boost::shared_ptr<::boost::shared_mutex> _mutex;
-
- /*! @brief For each dose grid index z, calculate the weight vector for each structure contour
- */
- void calcWeightVector(const rttb::VoxelGridID& aIndexZ,
- std::map<double, double>& weightVector) const;
- };
-
- }
-
-
- }
-}
-
+// -----------------------------------------------------------------------
+// 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: 1127 $ (last changed revision)
+// @date $Date: 2015-10-01 13:33:33 +0200 (Do, 01 Okt 2015) $ (last change date)
+// @author $Author: hentsch $ (last changed by)
+*/
+
+
+#ifndef __BOOST_MASK_GENERATE_MASK_VOXEL_LIST_H
+#define __BOOST_MASK_GENERATE_MASK_VOXEL_LIST_H
+
+#include "rttbBaseType.h"
+#include "rttbGeometricInfo.h"
+#include "rttbMaskAccessorInterface.h"
+
+#include <boost/multi_array.hpp>
+#include <boost/shared_ptr.hpp>
+#include <boost/thread/shared_mutex.hpp>
+
+namespace rttb
+{
+ namespace masks
+ {
+ namespace boost
+ {
+ /*! @class BoostMaskGenerateMaskVoxelListThread
+ *
+ */
+ class BoostMaskGenerateMaskVoxelListThread
+ {
+
+ public:
+ typedef ::boost::multi_array<double, 2> BoostArray2D;
+ using BoostArray2DPointer = ::boost::shared_ptr<BoostArray2D>;
+ typedef ::boost::shared_ptr<std::map<double, BoostArray2DPointer> > BoostArrayMapPointer;
+ using VoxelIndexVector = std::vector<rttb::VoxelGridIndex3D>;
+ using MaskVoxelListPointer = core::MaskAccessorInterface::MaskVoxelListPointer;
+
+ BoostMaskGenerateMaskVoxelListThread(const VoxelIndexVector& aGlobalBoundingBox,
+ core::GeometricInfo::Pointer aGeometricInfo,
+ BoostArrayMapPointer aVoxelizationMap,
+ double aVoxelizationThickness,
+ unsigned int aBeginSlice,
+ unsigned int aEndSlice,
+ MaskVoxelListPointer aMaskVoxelList, ::boost::shared_ptr<::boost::shared_mutex> aMutex);
+ void operator()();
+
+ private:
+ VoxelIndexVector _globalBoundingBox;
+ core::GeometricInfo::Pointer _geometricInfo;
+ BoostArrayMapPointer _voxelizationMap;
+ //(for example, the first contour has the double grid index 0.1, the second 0.3, the third 0.5, then the thickness is 0.2)
+ double _voxelizationThickness;
+
+ unsigned int _beginSlice;
+ /*! @brief _endSlice is the first index not to be processed (like a end iterator)
+ */
+ unsigned int _endSlice;
+
+ MaskVoxelListPointer _resultMaskVoxelList;
+ ::boost::shared_ptr<::boost::shared_mutex> _mutex;
+
+ /*! @brief For each dose grid index z, calculate the weight vector for each structure contour
+ */
+ void calcWeightVector(const rttb::VoxelGridID& aIndexZ,
+ std::map<double, double>& weightVector) const;
+ };
+
+ }
+
+
+ }
+}
+
#endif
\ No newline at end of file
diff --git a/code/masks/boost/rttbBoostMaskVoxelizationThread.cpp b/code/masks/rttbBoostMaskVoxelizationThread.cpp
similarity index 97%
rename from code/masks/boost/rttbBoostMaskVoxelizationThread.cpp
rename to code/masks/rttbBoostMaskVoxelizationThread.cpp
index 9b4deec..d184883 100644
--- a/code/masks/boost/rttbBoostMaskVoxelizationThread.cpp
+++ b/code/masks/rttbBoostMaskVoxelizationThread.cpp
@@ -1,174 +1,174 @@
-// -----------------------------------------------------------------------
-// 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: 1127 $ (last changed revision)
-// @date $Date: 2015-10-01 13:33:33 +0200 (Do, 01 Okt 2015) $ (last change date)
-// @author $Author: hentsch $ (last changed by)
-*/
-
-#include "rttbBoostMaskVoxelizationThread.h"
-
-#include "rttbInvalidParameterException.h"
-
-#include <boost/geometry.hpp>
-#include <boost/thread.hpp>
-#include <boost/make_shared.hpp>
-
-namespace rttb
-{
- namespace masks
- {
- namespace boost
- {
- BoostMaskVoxelizationThread::BoostMaskVoxelizationThread(const BoostPolygonMap& APolygonMap,
- const VoxelIndexVector& aGlobalBoundingBox, BoostArrayMapPointer anArrayMap, ::boost::shared_ptr<::boost::shared_mutex> aMutex, bool strict) : _geometryCoordinateBoostPolygonMap(APolygonMap),
- _globalBoundingBox(aGlobalBoundingBox), _resultVoxelization(anArrayMap), _mutex(aMutex), _strict(strict)
- {
- }
-
- void BoostMaskVoxelizationThread::operator()()
- {
- rttb::VoxelGridIndex3D minIndex = _globalBoundingBox.at(0);
- rttb::VoxelGridIndex3D maxIndex = _globalBoundingBox.at(1);
- const unsigned int globalBoundingBoxSize0 = maxIndex[0] - minIndex[0] + 1;
- const unsigned int globalBoundingBoxSize1 = maxIndex[1] - minIndex[1] + 1;
-
- std::map<double, ::boost::shared_ptr<BoostArray2D> > voxelizationMapInThread;
-
- for (auto & it : _geometryCoordinateBoostPolygonMap)
- {
- BoostArray2D maskArray(::boost::extents[globalBoundingBoxSize0][globalBoundingBoxSize1]);
-
- BoostPolygonVector boostPolygonVec = it.second;
-
- for (unsigned int x = 0; x < globalBoundingBoxSize0; ++x)
- {
- for (unsigned int y = 0; y < globalBoundingBoxSize1; ++y)
- {
- rttb::VoxelGridIndex3D currentIndex;
- currentIndex[0] = x + minIndex[0];
- currentIndex[1] = y + minIndex[1];
- currentIndex[2] = 0;
-
- //Get intersection polygons of the dose voxel and the structure
- BoostPolygonDeque polygons = getIntersections(currentIndex, boostPolygonVec);
-
- //Calc areas of all intersection polygons
- double volumeFraction = calcArea(polygons);
- volumeFraction = correctForErrorAndStrictness(volumeFraction, _strict);
- if (volumeFraction < 0 || volumeFraction > 1 )
- {
- throw rttb::core::InvalidParameterException("Mask calculation failed! The volume fraction should >= 0 and <= 1!");
- }
-
- maskArray[x][y] = volumeFraction;
- }
- }
-
- voxelizationMapInThread.insert(std::pair<double, BoostArray2DPointer>(it.first, ::boost::make_shared<BoostArray2D>(maskArray)));
- }
- //insert gathered values into voxelization map
- ::boost::unique_lock<::boost::shared_mutex> lock(*_mutex);
- _resultVoxelization->insert(voxelizationMapInThread.begin(), voxelizationMapInThread.end());
-
- }
-
- /*Get intersection polygons of the contour and a voxel polygon*/
- BoostMaskVoxelizationThread::BoostPolygonDeque BoostMaskVoxelizationThread::getIntersections(
- const rttb::VoxelGridIndex3D&
- aVoxelIndex3D, const BoostPolygonVector& intersectionSlicePolygons)
- {
- BoostMaskVoxelizationThread::BoostPolygonDeque polygonDeque;
-
- BoostRing2D voxelPolygon = get2DContour(aVoxelIndex3D);
- ::boost::geometry::correct(voxelPolygon);
-
- BoostPolygonVector::const_iterator it;
-
- for (it = intersectionSlicePolygons.begin(); it != intersectionSlicePolygons.end(); ++it)
- {
- BoostPolygon2D contour = *it;
- ::boost::geometry::correct(contour);
-
- BoostPolygonDeque intersection;
- ::boost::geometry::intersection(voxelPolygon, contour, intersection);
- polygonDeque.insert(polygonDeque.end(), intersection.begin(), intersection.end());
- }
-
- return polygonDeque;
- }
-
- BoostMaskVoxelizationThread::BoostRing2D BoostMaskVoxelizationThread::get2DContour(
- const rttb::VoxelGridIndex3D& aVoxelGrid3D)
- {
- BoostRing2D polygon;
-
-
- BoostPoint2D point1(aVoxelGrid3D[0] - 0.5, aVoxelGrid3D[1] - 0.5);
- ::boost::geometry::append(polygon, point1);
-
- BoostPoint2D point2(aVoxelGrid3D[0] + 0.5, aVoxelGrid3D[1] - 0.5);
- ::boost::geometry::append(polygon, point2);
-
- BoostPoint2D point3(aVoxelGrid3D[0] + 0.5, aVoxelGrid3D[1] + 0.5);
- ::boost::geometry::append(polygon, point3);
-
- BoostPoint2D point4(aVoxelGrid3D[0] - 0.5, aVoxelGrid3D[1] + 0.5);
- ::boost::geometry::append(polygon, point4);
-
- ::boost::geometry::append(polygon, point1);
-
- return polygon;
-
- }
-
- /*Calculate the intersection area*/
- double BoostMaskVoxelizationThread::calcArea(const BoostPolygonDeque& aPolygonDeque)
- {
- double area = 0;
-
- BoostPolygonDeque::const_iterator it;
-
- for (it = aPolygonDeque.begin(); it != aPolygonDeque.end(); ++it)
- {
- area += ::boost::geometry::area(*it);
- }
-
- return area;
- }
-
- double BoostMaskVoxelizationThread::correctForErrorAndStrictness(double volumeFraction, bool strict) const
- {
- if (strict){
- if (volumeFraction > 1 && (volumeFraction - 1) <= errorConstant)
- {
- volumeFraction = 1;
- }
- }
- else {
- if (volumeFraction > 1){
- volumeFraction = 1;
- }
- else if (volumeFraction < 0){
- volumeFraction = 0;
- }
- }
- return volumeFraction;
- }
-
- }
- }
-}
+// -----------------------------------------------------------------------
+// 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: 1127 $ (last changed revision)
+// @date $Date: 2015-10-01 13:33:33 +0200 (Do, 01 Okt 2015) $ (last change date)
+// @author $Author: hentsch $ (last changed by)
+*/
+
+#include "rttbBoostMaskVoxelizationThread.h"
+
+#include "rttbInvalidParameterException.h"
+
+#include <boost/geometry.hpp>
+#include <boost/thread.hpp>
+#include <boost/make_shared.hpp>
+
+namespace rttb
+{
+ namespace masks
+ {
+ namespace boost
+ {
+ BoostMaskVoxelizationThread::BoostMaskVoxelizationThread(const BoostPolygonMap& APolygonMap,
+ const VoxelIndexVector& aGlobalBoundingBox, BoostArrayMapPointer anArrayMap, ::boost::shared_ptr<::boost::shared_mutex> aMutex, bool strict) : _geometryCoordinateBoostPolygonMap(APolygonMap),
+ _globalBoundingBox(aGlobalBoundingBox), _resultVoxelization(anArrayMap), _mutex(aMutex), _strict(strict)
+ {
+ }
+
+ void BoostMaskVoxelizationThread::operator()()
+ {
+ rttb::VoxelGridIndex3D minIndex = _globalBoundingBox.at(0);
+ rttb::VoxelGridIndex3D maxIndex = _globalBoundingBox.at(1);
+ const unsigned int globalBoundingBoxSize0 = maxIndex[0] - minIndex[0] + 1;
+ const unsigned int globalBoundingBoxSize1 = maxIndex[1] - minIndex[1] + 1;
+
+ std::map<double, ::boost::shared_ptr<BoostArray2D> > voxelizationMapInThread;
+
+ for (auto & it : _geometryCoordinateBoostPolygonMap)
+ {
+ BoostArray2D maskArray(::boost::extents[globalBoundingBoxSize0][globalBoundingBoxSize1]);
+
+ BoostPolygonVector boostPolygonVec = it.second;
+
+ for (unsigned int x = 0; x < globalBoundingBoxSize0; ++x)
+ {
+ for (unsigned int y = 0; y < globalBoundingBoxSize1; ++y)
+ {
+ rttb::VoxelGridIndex3D currentIndex;
+ currentIndex[0] = x + minIndex[0];
+ currentIndex[1] = y + minIndex[1];
+ currentIndex[2] = 0;
+
+ //Get intersection polygons of the dose voxel and the structure
+ BoostPolygonDeque polygons = getIntersections(currentIndex, boostPolygonVec);
+
+ //Calc areas of all intersection polygons
+ double volumeFraction = calcArea(polygons);
+ volumeFraction = correctForErrorAndStrictness(volumeFraction, _strict);
+ if (volumeFraction < 0 || volumeFraction > 1 )
+ {
+ throw rttb::core::InvalidParameterException("Mask calculation failed! The volume fraction should >= 0 and <= 1!");
+ }
+
+ maskArray[x][y] = volumeFraction;
+ }
+ }
+
+ voxelizationMapInThread.insert(std::pair<double, BoostArray2DPointer>(it.first, ::boost::make_shared<BoostArray2D>(maskArray)));
+ }
+ //insert gathered values into voxelization map
+ ::boost::unique_lock<::boost::shared_mutex> lock(*_mutex);
+ _resultVoxelization->insert(voxelizationMapInThread.begin(), voxelizationMapInThread.end());
+
+ }
+
+ /*Get intersection polygons of the contour and a voxel polygon*/
+ BoostMaskVoxelizationThread::BoostPolygonDeque BoostMaskVoxelizationThread::getIntersections(
+ const rttb::VoxelGridIndex3D&
+ aVoxelIndex3D, const BoostPolygonVector& intersectionSlicePolygons)
+ {
+ BoostMaskVoxelizationThread::BoostPolygonDeque polygonDeque;
+
+ BoostRing2D voxelPolygon = get2DContour(aVoxelIndex3D);
+ ::boost::geometry::correct(voxelPolygon);
+
+ BoostPolygonVector::const_iterator it;
+
+ for (it = intersectionSlicePolygons.begin(); it != intersectionSlicePolygons.end(); ++it)
+ {
+ BoostPolygon2D contour = *it;
+ ::boost::geometry::correct(contour);
+
+ BoostPolygonDeque intersection;
+ ::boost::geometry::intersection(voxelPolygon, contour, intersection);
+ polygonDeque.insert(polygonDeque.end(), intersection.begin(), intersection.end());
+ }
+
+ return polygonDeque;
+ }
+
+ BoostMaskVoxelizationThread::BoostRing2D BoostMaskVoxelizationThread::get2DContour(
+ const rttb::VoxelGridIndex3D& aVoxelGrid3D)
+ {
+ BoostRing2D polygon;
+
+
+ BoostPoint2D point1(aVoxelGrid3D[0] - 0.5, aVoxelGrid3D[1] - 0.5);
+ ::boost::geometry::append(polygon, point1);
+
+ BoostPoint2D point2(aVoxelGrid3D[0] + 0.5, aVoxelGrid3D[1] - 0.5);
+ ::boost::geometry::append(polygon, point2);
+
+ BoostPoint2D point3(aVoxelGrid3D[0] + 0.5, aVoxelGrid3D[1] + 0.5);
+ ::boost::geometry::append(polygon, point3);
+
+ BoostPoint2D point4(aVoxelGrid3D[0] - 0.5, aVoxelGrid3D[1] + 0.5);
+ ::boost::geometry::append(polygon, point4);
+
+ ::boost::geometry::append(polygon, point1);
+
+ return polygon;
+
+ }
+
+ /*Calculate the intersection area*/
+ double BoostMaskVoxelizationThread::calcArea(const BoostPolygonDeque& aPolygonDeque)
+ {
+ double area = 0;
+
+ BoostPolygonDeque::const_iterator it;
+
+ for (it = aPolygonDeque.begin(); it != aPolygonDeque.end(); ++it)
+ {
+ area += ::boost::geometry::area(*it);
+ }
+
+ return area;
+ }
+
+ double BoostMaskVoxelizationThread::correctForErrorAndStrictness(double volumeFraction, bool strict) const
+ {
+ if (strict){
+ if (volumeFraction > 1 && (volumeFraction - 1) <= errorConstant)
+ {
+ volumeFraction = 1;
+ }
+ }
+ else {
+ if (volumeFraction > 1){
+ volumeFraction = 1;
+ }
+ else if (volumeFraction < 0){
+ volumeFraction = 0;
+ }
+ }
+ return volumeFraction;
+ }
+
+ }
+ }
+}
diff --git a/code/masks/boost/rttbBoostMaskVoxelizationThread.h b/code/masks/rttbBoostMaskVoxelizationThread.h
similarity index 97%
rename from code/masks/boost/rttbBoostMaskVoxelizationThread.h
rename to code/masks/rttbBoostMaskVoxelizationThread.h
index 9743c4e..3ffd14c 100644
--- a/code/masks/boost/rttbBoostMaskVoxelizationThread.h
+++ b/code/masks/rttbBoostMaskVoxelizationThread.h
@@ -1,109 +1,109 @@
-// -----------------------------------------------------------------------
-// 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: 1127 $ (last changed revision)
-// @date $Date: 2015-10-01 13:33:33 +0200 (Do, 01 Okt 2015) $ (last change date)
-// @author $Author: hentsch $ (last changed by)
-*/
-
-
-#ifndef __BOOST_MASK_VOXELIZATION_THREAD_H
-#define __BOOST_MASK_VOXELIZATION_THREAD_H
-
-#include <deque>
-
-#include "rttbBaseType.h"
-
-#include <boost/multi_array.hpp>
-#include <boost/shared_ptr.hpp>
-#include <boost/thread/shared_mutex.hpp>
-#include <boost/geometry/geometries/point_xy.hpp>
-#include <boost/geometry/geometries/polygon.hpp>
-
-namespace rttb
-{
- namespace masks
- {
- namespace boost
- {
- /*! @class BoostMaskGenerateMaskVoxelListThread
- *
- */
- class BoostMaskVoxelizationThread
- {
-
- public:
- using BoostPolygon2D = ::boost::geometry::model::polygon< ::boost::geometry::model::d2::point_xy<double> >;
- using BoostPolygonVector = std::vector<BoostPolygon2D>;//polygon with or without holes
- typedef std::map<double, BoostPolygonVector>
- BoostPolygonMap;//map of the z index with the vector of boost 2d polygon
- using VoxelIndexVector = std::vector<rttb::VoxelGridIndex3D>;
- typedef ::boost::multi_array<double, 2> BoostArray2D;
- using BoostArray2DPointer = ::boost::shared_ptr<BoostArray2D>;
- typedef ::boost::shared_ptr<std::map<double, BoostArray2DPointer > > BoostArrayMapPointer;
-
- /*! @brief Constructor
- * @param aMutex a mutex for thread-safe handling of the _resultVoxelization
- * @param strict true means that volumeFractions of <0 and >1 are NOT corrected. Otherwise, they are automatically corrected to 0 or 1, respectively.
- */
- BoostMaskVoxelizationThread(const BoostPolygonMap& APolygonMap,
- const VoxelIndexVector& aGlobalBoundingBox, BoostArrayMapPointer anArrayMap, ::boost::shared_ptr<::boost::shared_mutex> aMutex, bool strict);
-
- void operator()();
-
-
- private:
- using BoostPolygonDeque = std::deque<BoostPolygon2D>;
- using BoostRing2D = ::boost::geometry::model::ring< ::boost::geometry::model::d2::point_xy<double> >;
- using BoostPoint2D = ::boost::geometry::model::d2::point_xy<double>;
-
-
- BoostPolygonMap _geometryCoordinateBoostPolygonMap;
- VoxelIndexVector _globalBoundingBox;
- BoostArrayMapPointer _resultVoxelization;
- ::boost::shared_ptr<::boost::shared_mutex> _mutex;
- bool _strict;
-
- /*! @brief Get intersection polygons of the contour and a voxel polygon
- * @param aVoxelIndex3D The 3d grid index of the voxel
- * @param intersectionSlicePolygons The polygons of the slice intersecting the voxel
- * @return Return all intersection polygons of the structure and the voxel
- */
- static BoostPolygonDeque getIntersections(const rttb::VoxelGridIndex3D& aVoxelIndex3D,
- const BoostPolygonVector& intersectionSlicePolygons);
-
- /*! @brief Get the voxel 2d contour polygon in geometry coordinate*/
- static BoostRing2D get2DContour(const rttb::VoxelGridIndex3D& aVoxelGrid3D);
-
- /*! @brief Calculate the area of all polygons
- * @param aPolygonDeque The deque of polygons
- * @return Return the area of all polygons
- */
- static double calcArea(const BoostPolygonDeque& aPolygonDeque);
- /*! @brief Corrects the volumeFraction
- * @details the volume fraction is corrected in case of strict=true. Otherwise, it's only corrected for double imprecision
- * @return The corrected volumeFraction
- */
- double correctForErrorAndStrictness(double volumeFraction, bool strict) const;
- };
-
- }
-
-
- }
-}
-
+// -----------------------------------------------------------------------
+// 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: 1127 $ (last changed revision)
+// @date $Date: 2015-10-01 13:33:33 +0200 (Do, 01 Okt 2015) $ (last change date)
+// @author $Author: hentsch $ (last changed by)
+*/
+
+
+#ifndef __BOOST_MASK_VOXELIZATION_THREAD_H
+#define __BOOST_MASK_VOXELIZATION_THREAD_H
+
+#include <deque>
+
+#include "rttbBaseType.h"
+
+#include <boost/multi_array.hpp>
+#include <boost/shared_ptr.hpp>
+#include <boost/thread/shared_mutex.hpp>
+#include <boost/geometry/geometries/point_xy.hpp>
+#include <boost/geometry/geometries/polygon.hpp>
+
+namespace rttb
+{
+ namespace masks
+ {
+ namespace boost
+ {
+ /*! @class BoostMaskGenerateMaskVoxelListThread
+ *
+ */
+ class BoostMaskVoxelizationThread
+ {
+
+ public:
+ using BoostPolygon2D = ::boost::geometry::model::polygon< ::boost::geometry::model::d2::point_xy<double> >;
+ using BoostPolygonVector = std::vector<BoostPolygon2D>;//polygon with or without holes
+ typedef std::map<double, BoostPolygonVector>
+ BoostPolygonMap;//map of the z index with the vector of boost 2d polygon
+ using VoxelIndexVector = std::vector<rttb::VoxelGridIndex3D>;
+ typedef ::boost::multi_array<double, 2> BoostArray2D;
+ using BoostArray2DPointer = ::boost::shared_ptr<BoostArray2D>;
+ typedef ::boost::shared_ptr<std::map<double, BoostArray2DPointer > > BoostArrayMapPointer;
+
+ /*! @brief Constructor
+ * @param aMutex a mutex for thread-safe handling of the _resultVoxelization
+ * @param strict true means that volumeFractions of <0 and >1 are NOT corrected. Otherwise, they are automatically corrected to 0 or 1, respectively.
+ */
+ BoostMaskVoxelizationThread(const BoostPolygonMap& APolygonMap,
+ const VoxelIndexVector& aGlobalBoundingBox, BoostArrayMapPointer anArrayMap, ::boost::shared_ptr<::boost::shared_mutex> aMutex, bool strict);
+
+ void operator()();
+
+
+ private:
+ using BoostPolygonDeque = std::deque<BoostPolygon2D>;
+ using BoostRing2D = ::boost::geometry::model::ring< ::boost::geometry::model::d2::point_xy<double> >;
+ using BoostPoint2D = ::boost::geometry::model::d2::point_xy<double>;
+
+
+ BoostPolygonMap _geometryCoordinateBoostPolygonMap;
+ VoxelIndexVector _globalBoundingBox;
+ BoostArrayMapPointer _resultVoxelization;
+ ::boost::shared_ptr<::boost::shared_mutex> _mutex;
+ bool _strict;
+
+ /*! @brief Get intersection polygons of the contour and a voxel polygon
+ * @param aVoxelIndex3D The 3d grid index of the voxel
+ * @param intersectionSlicePolygons The polygons of the slice intersecting the voxel
+ * @return Return all intersection polygons of the structure and the voxel
+ */
+ static BoostPolygonDeque getIntersections(const rttb::VoxelGridIndex3D& aVoxelIndex3D,
+ const BoostPolygonVector& intersectionSlicePolygons);
+
+ /*! @brief Get the voxel 2d contour polygon in geometry coordinate*/
+ static BoostRing2D get2DContour(const rttb::VoxelGridIndex3D& aVoxelGrid3D);
+
+ /*! @brief Calculate the area of all polygons
+ * @param aPolygonDeque The deque of polygons
+ * @return Return the area of all polygons
+ */
+ static double calcArea(const BoostPolygonDeque& aPolygonDeque);
+ /*! @brief Corrects the volumeFraction
+ * @details the volume fraction is corrected in case of strict=true. Otherwise, it's only corrected for double imprecision
+ * @return The corrected volumeFraction
+ */
+ double correctForErrorAndStrictness(double volumeFraction, bool strict) const;
+ };
+
+ }
+
+
+ }
+}
+
#endif
\ No newline at end of file