diff --git a/CMakeLists.txt b/CMakeLists.txt index 3a00be5..422b8f9 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -1,298 +1,297 @@ #----------------------------------------------------------------------------- # This is the root RTToolbox CMakeList file. #----------------------------------------------------------------------------- PROJECT(RTToolbox) CMAKE_MINIMUM_REQUIRED(VERSION 2.8.11.2) IF(CMAKE_BACKWARDS_COMPATIBILITY GREATER 2.8.11.2) SET(CMAKE_BACKWARDS_COMPATIBILITY 2.8.11.2 CACHE STRING "Latest version of CMake when this project was released." FORCE) ENDIF(CMAKE_BACKWARDS_COMPATIBILITY GREATER 2.8.11.2) # RTToolbox version number. SET(RTToolbox_VERSION_MAJOR "4") -SET(RTToolbox_VERSION_MINOR "0") +SET(RTToolbox_VERSION_MINOR "1") SET(RTToolbox_VERSION_PATCH "0") # Version string should not include patch level. The major.minor is # enough to distinguish available features of the toolbox. SET(RTToolbox_VERSION_STRING "${RTToolbox_VERSION_MAJOR}.${RTToolbox_VERSION_MINOR}") SET(RTToolbox_FULL_VERSION_STRING "${RTToolbox_VERSION_MAJOR}.${RTToolbox_VERSION_MINOR}.${RTToolbox_VERSION_PATCH}") # default build type SET(CMAKE_BUILD_TYPE Release) MARK_AS_ADVANCED(BUILD_SHARED_LIBS) IF (WIN32) IF (MSVC_VERSION LESS 1600) MESSAGE(FATAL_ERROR "RTToolbox requires at least Visual Studio 2010.") ENDIF(MSVC_VERSION LESS 1600) add_definitions(-D_SCL_SECURE_NO_WARNINGS) ELSE (WIN32) IF (CMAKE_COMPILER_IS_GNUCC) IF (CMAKE_CXX_COMPILER_VERSION VERSION_LESS 4.6 OR CMAKE_CXX_COMPILER_VERSION VERSION_EQUAL 5.0 OR CMAKE_CXX_COMPILER_VERSION VERSION_GREATER 5.0) MESSAGE(AUTHOR_WARNING "RTToolbox was only tested with GCC 4.6 and GCC 4.9. You are using GCC " ${CMAKE_CXX_COMPILER_VERSION} ". This compiler version might not work.") ENDIF() ENDIF() ENDIF(WIN32) IF (NOT (CMAKE_MAJOR_VERSION LESS 3)) IF(COMMAND CMAKE_POLICY) # Enable old CMake behaviour when dealing with export_library_dependencies(). # This is necessary to avoid warnings in CMake versions # greater than 3.0 # See http://www.cmake.org/cmake/help/v3.0/policy/CMP0033.html CMAKE_POLICY(SET CMP0033 OLD) ENDIF(COMMAND CMAKE_POLICY) ENDIF() -message(STATUS ${CMAKE_BUILD_TYPE}) #----------------------------------------------------------------------------- # CMake Function(s) and Macro(s) #----------------------------------------------------------------------------- include(cmake/MacroParseArguments.cmake) include(cmake/rttbMacroCreateModuleConf.cmake) include(cmake/rttbMacroCreateModule.cmake) include(cmake/rttbMacroCreateApplication.cmake) include(cmake/rttbMacroCheckModule.cmake) include(cmake/rttbMacroUseModule.cmake) include(cmake/rttbMacroCreateTestModule.cmake) include(cmake/rttbFunctionOrganizeSources.cmake) #----------------------------------------------------------------------------- # Basis config RTTB module infrastructure #----------------------------------------------------------------------------- set(RTTB_MODULES_CONF_DIR ${RTToolbox_BINARY_DIR}/modulesConf CACHE INTERNAL "Modules Conf") set(RTTB_MODULES_PACKAGE_DEPENDS_DIR ${RTToolbox_SOURCE_DIR}/cmake/PackageDepends) set(MODULES_PACKAGE_DEPENDS_DIRS ${RTTB_MODULES_PACKAGE_DEPENDS_DIR}) #----------------------------------------------------------------------------- # Testing setup # Configure Dart testing support. This should be done before any # MESSAGE(FATAL_ERROR ...) commands are invoked. #----------------------------------------------------------------------------- SET(CTEST_NEW_FORMAT 1) INCLUDE(CTest) ENABLE_TESTING() IF(BUILD_TESTING) CONFIGURE_FILE(${RTToolbox_SOURCE_DIR}/cmake/RemoveTemporaryFiles.cmake.in ${RTToolbox_BINARY_DIR}/cmake/RemoveTemporaryFiles.cmake @ONLY IMMEDIATE) CONFIGURE_FILE(${RTToolbox_SOURCE_DIR}/cmake/rttbSampleBuildTest.cmake.in ${RTToolbox_BINARY_DIR}/cmake/rttbSampleBuildTest.cmake @ONLY) CONFIGURE_FILE(${RTToolbox_SOURCE_DIR}/cmake/CTestCustom.ctest.in ${RTToolbox_BINARY_DIR}/cmake/CTestCustom.ctest @ONLY) FILE(WRITE ${RTToolbox_BINARY_DIR}/CTestCustom.cmake "INCLUDE(\"${RTToolbox_BINARY_DIR}/cmake/CTestCustom.ctest\")\n") SET(BUILDNAME "${BUILDNAME}" CACHE STRING "Name of build on the dashboard") MARK_AS_ADVANCED(BUILDNAME) ENDIF(BUILD_TESTING) #----------------------------------------------------------------------------- # Output directories. #----------------------------------------------------------------------------- IF(NOT LIBRARY_OUTPUT_PATH) SET (LIBRARY_OUTPUT_PATH ${RTToolbox_BINARY_DIR}/bin CACHE PATH "Single output directory for building all libraries.") ENDIF(NOT LIBRARY_OUTPUT_PATH) IF(NOT EXECUTABLE_OUTPUT_PATH) SET (EXECUTABLE_OUTPUT_PATH ${RTToolbox_BINARY_DIR}/bin CACHE PATH "Single output directory for building all executables.") ENDIF(NOT EXECUTABLE_OUTPUT_PATH) MARK_AS_ADVANCED(EXECUTABLE_OUTPUT_PATH LIBRARY_OUTPUT_PATH) MARK_AS_ADVANCED(LIBRARY_OUTPUT_PATH EXECUTABLE_OUTPUT_PATH) SET(RTToolbox_LIBRARY_PATH "${LIBRARY_OUTPUT_PATH}") SET(RTToolbox_EXECUTABLE_PATH "${EXECUTABLE_OUTPUT_PATH}") #----------------------------------------------------------------------------- # Find Doxygen. #----------------------------------------------------------------------------- FIND_PROGRAM(DOXYGEN_EXECUTABLE "doxygen") #----------------------------------------------------------------------------- # Installation vars. # RTToolbox_INSTALL_BIN_DIR - binary dir (executables) # RTToolbox_INSTALL_LIB_DIR - library dir (libs) # RTToolbox_INSTALL_INCLUDE_DIR - include dir (headers) # RTToolbox_INSTALL_NO_DEVELOPMENT - do not install development files # RTToolbox_INSTALL_NO_RUNTIME - do not install runtime files # RTToolbox_INSTALL_NO_DOCUMENTATION - do not install documentation files # Remark: needs directory are stored with no leading slash (CMake 2.4 and newer) #----------------------------------------------------------------------------- IF(NOT RTTOOLBOX_INSTALL_BIN_DIR) SET(RTTOOLBOX_INSTALL_BIN_DIR "bin") ENDIF(NOT RTTOOLBOX_INSTALL_BIN_DIR) IF(NOT RTTOOLBOX_INSTALL_LIB_DIR) SET(RTTOOLBOX_INSTALL_LIB_DIR "lib") ENDIF(NOT RTTOOLBOX_INSTALL_LIB_DIR) IF(NOT RTTOOLBOX_INSTALL_PACKAGE_DIR) SET(RTTOOLBOX_INSTALL_PACKAGE_DIR "lib") ENDIF(NOT RTTOOLBOX_INSTALL_PACKAGE_DIR) IF(NOT RTTOOLBOX_INSTALL_INCLUDE_DIR) SET(RTTOOLBOX_INSTALL_INCLUDE_DIR "include") ENDIF(NOT RTTOOLBOX_INSTALL_INCLUDE_DIR) IF(NOT RTTOOLBOX_INSTALL_NO_DEVELOPMENT) SET(RTTOOLBOX_INSTALL_NO_DEVELOPMENT 0) ENDIF(NOT RTTOOLBOX_INSTALL_NO_DEVELOPMENT) IF(NOT RTTOOLBOX_INSTALL_NO_RUNTIME) SET(RTTOOLBOX_INSTALL_NO_RUNTIME 0) ENDIF(NOT RTTOOLBOX_INSTALL_NO_RUNTIME) IF(NOT RTTOOLBOX_INSTALL_NO_DOCUMENTATION) SET(RTTOOLBOX_INSTALL_NO_DOCUMENTATION 0) ENDIF(NOT RTTOOLBOX_INSTALL_NO_DOCUMENTATION) SET(RTTOOLBOX_INSTALL_NO_LIBRARIES) IF(RTTOOLBOX_BUILD_SHARED_LIBS) IF(RTTOOLBOX_INSTALL_NO_RUNTIME AND RTTOOLBOX_INSTALL_NO_DEVELOPMENT) SET(RTTOOLBOX_INSTALL_NO_LIBRARIES 1) ENDIF(RTTOOLBOX_INSTALL_NO_RUNTIME AND RTTOOLBOX_INSTALL_NO_DEVELOPMENT) ELSE(RTTOOLBOX_BUILD_SHARED_LIBS) IF(RTTOOLBOX_INSTALL_NO_DEVELOPMENT) SET(RTTOOLBOX_INSTALL_NO_LIBRARIES 1) ENDIF(RTTOOLBOX_INSTALL_NO_DEVELOPMENT) ENDIF(RTTOOLBOX_BUILD_SHARED_LIBS) # set RTToolbox_DIR so it can be used by subprojects SET(RTToolbox_DIR "${CMAKE_BINARY_DIR}" CACHE INTERNAL "RTToolbox dir to be used by subprojects") #----------------------------------------------------------------------------- # DCMTK MT-Flag treat #----------------------------------------------------------------------------- option(RTTB_DCMTK_COMPLIANCE_ENFORCE_MT "This enforces the whole RTToolbox to be compiled with /MT,/MTd to be compliant with DCMTK" OFF) string(FIND ${CMAKE_GENERATOR} "Visual Studio" RTTB_VS_USED) if(RTTB_DCMTK_COMPLIANCE_ENFORCE_MT AND RTTB_VS_USED EQUAL 0) message(STATUS "Enforce DCMTK compliance: /MT and /MTd flags are used") string(REPLACE "/MD" "/MT" CMAKE_C_FLAGS_DEBUG ${CMAKE_C_FLAGS_DEBUG}) message(STATUS "CMAKE_C_FLAGS_DEBUG set to: ${CMAKE_C_FLAGS_DEBUG}") string(REPLACE "/MD" "/MT" CMAKE_C_FLAGS_RELEASE ${CMAKE_C_FLAGS_RELEASE}) message(STATUS "CMAKE_C_FLAGS_RELEASE set to: ${CMAKE_C_FLAGS_RELEASE}") string(REPLACE "/MD" "/MT" CMAKE_C_FLAGS_MINSIZEREL ${CMAKE_C_FLAGS_MINSIZEREL}) message(STATUS "CMAKE_C_FLAGS_MINSIZEREL set to: ${CMAKE_C_FLAGS_MINSIZEREL}") string(REPLACE "/MD" "/MT" CMAKE_C_FLAGS_RELWITHDEBINFO ${CMAKE_C_FLAGS_RELWITHDEBINFO}) message(STATUS "CMAKE_C_FLAGS_RELWITHDEBINFO set to: ${CMAKE_C_FLAGS_RELWITHDEBINFO}") string(REPLACE "/MD" "/MT" CMAKE_CXX_FLAGS_DEBUG ${CMAKE_CXX_FLAGS_DEBUG}) message(STATUS "CMAKE_CXX_FLAGS_DEBUG set to: ${CMAKE_CXX_FLAGS_DEBUG}") string(REPLACE "/MD" "/MT" CMAKE_CXX_FLAGS_RELEASE ${CMAKE_CXX_FLAGS_RELEASE}) message(STATUS "CMAKE_CXX_FLAGS_RELEASE set to: ${CMAKE_CXX_FLAGS_RELEASE}") string(REPLACE "/MD" "/MT" CMAKE_CXX_FLAGS_MINSIZEREL ${CMAKE_CXX_FLAGS_MINSIZEREL}) message(STATUS "CMAKE_CXX_FLAGS_MINSIZEREL set to: ${CMAKE_CXX_FLAGS_MINSIZEREL}") string(REPLACE "/MD" "/MT" CMAKE_CXX_FLAGS_RELWITHDEBINFO ${CMAKE_CXX_FLAGS_RELWITHDEBINFO}) message(STATUS "CMAKE_CXX_FLAGS_RELWITHDEBINFO set to: ${CMAKE_CXX_FLAGS_RELWITHDEBINFO}") endif() #----------------------------------------------------------------------------- # Advanced RTToolbox configuration #----------------------------------------------------------------------------- #----------------------------------------------------------------------------- # RTToolbox build configuration options. IF (WIN32) OPTION(BUILD_SHARED_LIBS "Build RTToolbox with shared libraries." OFF) ELSE (WIN32) OPTION(BUILD_SHARED_LIBS "Build RTToolbox with shared libraries." ON) ENDIF (WIN32) IF(BUILD_SHARED_LIBS) IF(WIN32) MESSAGE(FATAL_ERROR "RTToolbox currently does not support a dynamic build on Windows. We are working on that...") ENDIF(WIN32) ELSE(BUILD_SHARED_LIBS) IF(UNIX) MESSAGE(FATAL_ERROR "RTToolbox currently does not support a static build on unix like systems. We are working on that...") ENDIF(UNIX) ENDIF(BUILD_SHARED_LIBS) SET(RTToolbox_BUILD_SHARED_LIBS ${BUILD_SHARED_LIBS}) IF(NOT RTToolbox_NO_LIBRARY_VERSION) # This setting of SOVERSION assumes that any API change # will increment either the minor or major version number of RTToolbox. SET(RTToolbox_LIBRARY_PROPERTIES VERSION "${RTToolbox_VERSION_MAJOR}.${RTToolbox_VERSION_MINOR}.${RTToolbox_VERSION_PATCH}" SOVERSION "${RTToolbox_VERSION_MAJOR}.${RTToolbox_VERSION_MINOR}" ) ENDIF(NOT RTToolbox_NO_LIBRARY_VERSION) #----------------------------------------------------------------------------- # Configure files with settings for use by the build. # #----------------------------------------------------------------------------- CONFIGURE_FILE(${RTToolbox_SOURCE_DIR}/RTToolboxConfigure.h.in ${RTToolbox_BINARY_DIR}/RTToolboxConfigure.h) IF(NOT RTTOOLBOX_INSTALL_NO_DEVELOPMENT) INSTALL(FILES ${RTToolbox_BINARY_DIR}/RTToolboxConfigure.h DESTINATION ${RTTOOLBOX_INSTALL_INCLUDE_DIR} COMPONENT Development) ENDIF(NOT RTTOOLBOX_INSTALL_NO_DEVELOPMENT) #----------------------------------------------------------------------------- # The entire RTToolbox tree should use the same include path #----------------------------------------------------------------------------- #Default include dir. Others dirs will be defined by activated subprojects INCLUDE_DIRECTORIES(${RTToolbox_BINARY_DIR}) LINK_DIRECTORIES(${LIBARY_OUTPUT_PATH}) #Prepare the correct target information export by the subprojects SET(RTToolbox_TARGETS_FILE "${RTToolbox_BINARY_DIR}/RTToolboxTargets.cmake") FILE(WRITE ${RTToolbox_TARGETS_FILE} "# Generated by CMake, do not edit!") #----------------------------------------------------------------------------- # Dispatch the build into the proper subdirectories. #----------------------------------------------------------------------------- MESSAGE (STATUS "generating Project RTToolbox") ADD_SUBDIRECTORY (code) ADD_SUBDIRECTORY (demoapps) IF (BUILD_TESTING) ADD_SUBDIRECTORY (testing) ENDIF (BUILD_TESTING) ADD_SUBDIRECTORY (documentation) #----------------------------------------------------------------------------- # Help other projects use RTToolbox. #----------------------------------------------------------------------------- EXPORT(PACKAGE RTToolbox) # Copy the UseRTToolbox.cmake file to the binary tree for backward compatability. CONFIGURE_FILE(${RTToolbox_SOURCE_DIR}/UseRTToolbox.cmake.in ${RTToolbox_BINARY_DIR}/UseRTToolbox.cmake COPYONLY IMMEDIATE) # Save library dependencies. EXPORT_LIBRARY_DEPENDENCIES(${RTToolbox_BINARY_DIR}/RTToolboxLibraryDepends.cmake) # Create the RTToolboxConfig.cmake file containing the RTToolbox configuration. INCLUDE (${RTToolbox_SOURCE_DIR}/rttbGenerateRTToolboxConfig.cmake) IF(NOT RTToolbox_INSTALL_NO_DEVELOPMENT) INSTALL(FILES ${RTToolbox_BINARY_DIR}/RTToolboxConfig.cmake ${RTToolbox_BINARY_DIR}/RTToolboxTargets.cmake ${RTToolbox_BINARY_DIR}/RTToolboxLibraryDepends.cmake ${RTToolbox_BINARY_DIR}/UseRTToolbox.cmake DESTINATION ${RTTOOLBOX_INSTALL_PACKAGE_DIR} COMPONENT Development ) ENDIF(NOT RTToolbox_INSTALL_NO_DEVELOPMENT) diff --git a/code/models/files.cmake b/code/models/files.cmake index a6ac489..7174a52 100644 --- a/code/models/files.cmake +++ b/code/models/files.cmake @@ -1,25 +1,30 @@ SET(CPP_FILES rttbBioModel.cpp rttbBioModelCurve.cpp rttbBioModelScatterPlots.cpp rttbNTCPLKBModel.cpp rttbNTCPRSModel.cpp rttbTCPModel.cpp rttbTCPLQModel.cpp + rttbDoseBasedModels.cpp rttbDvhBasedModels.cpp rttbIntegration.cpp + rttbLQModelAccessor.cpp ) SET(H_FILES rttbBaseTypeModels.h rttbBioModel.h rttbBioModelCurve.h rttbBioModelScatterPlots.h rttbNTCPModel.h rttbNTCPLKBModel.h rttbNTCPRSModel.h rttbTCPModel.h rttbTCPLQModel.h + rttbDoseBasedModels.cpp rttbDvhBasedModels.h rttbIntegration.h + rttbBioModelAccessorInterface.h + rttbLQModelAccessor.h ) diff --git a/code/models/rttbBioModelAccessorInterface.h b/code/models/rttbBioModelAccessorInterface.h new file mode 100644 index 0000000..d81bb97 --- /dev/null +++ b/code/models/rttbBioModelAccessorInterface.h @@ -0,0 +1,82 @@ +// ----------------------------------------------------------------------- +// 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: 747 $ (last changed revision) +// @date $Date: 2014-09-17 12:01:00 +0200 (Mi, 17 Sep 2014) $ (last change date) +// @author $Author: hentsch $ (last changed by) +*/ +#ifndef __BIO_MODEL_ACCESSOR_INTERFACE_H +#define __BIO_MODEL_ACCESSOR_INTERFACE_H + +#include + +#include "rttbBaseType.h" +#include "rttbBaseTypeModels.h" +#include "rttbGeometricInfo.h" +#include "rttbIndexConversionInterface.h" + +namespace rttb +{ + namespace core + { + + /*! @class BioModelAccessorInterface + @brief Interface for all BioModelAccessor classes + */ + class BioModelAccessorInterface: public IndexConversionInterface + { + public: + typedef boost::shared_ptr BioModelAccessorPointer; + private: + BioModelAccessorInterface(const BioModelAccessorInterface&); //not implemented on purpose -> non-copyable + BioModelAccessorInterface& operator=(const + BioModelAccessorInterface&);//not implemented on purpose -> non-copyable + + public: + BioModelAccessorInterface() {}; + virtual ~BioModelAccessorInterface() {}; + + inline const core::GeometricInfo& getGeometricInfo() const + { + return _geoInfo; + }; + + inline GridSizeType getGridSize() const + { + return _geoInfo.getNumberOfVoxels(); + }; + + virtual models::BioModelValueType getBioModelValueAt(const VoxelGridID aID) const = 0; + + virtual models::BioModelValueType getBioModelValueAt(const VoxelGridIndex3D& aIndex) const = 0; + + /*! @brief is true if dose is on a homogeneous grid + @remarks Inhomogeneous grids are not supported at the moment, but if they will be supported in the future + the interface does not need to change. + */ + virtual bool isGridHomogeneous() const + { + return true; + } + + virtual const IDType getBioModelUID() const = 0; + protected: + core::GeometricInfo _geoInfo; + }; + } +} + +#endif diff --git a/code/models/rttbDoseBasedModels.cpp b/code/models/rttbDoseBasedModels.cpp new file mode 100644 index 0000000..90147ce --- /dev/null +++ b/code/models/rttbDoseBasedModels.cpp @@ -0,0 +1,41 @@ +// ----------------------------------------------------------------------- +// 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 "rttbDoseBasedModels.h" +#include "rttbInvalidParameterException.h" + +#include + +namespace rttb +{ + namespace models + { + rttb::models::BioModelValueType calcLQ(const DoseTypeGy dose, const DoseCalcType alpha, const DoseCalcType beta) + { + if (dose < 0 || alpha < 0 || beta < 0) + { + throw core::InvalidParameterException("Parameter invalid: dose, alpha, beta must be >=0!"); + } + + return exp(-((alpha * dose) + (beta * dose * dose))); + } + + } +} \ No newline at end of file diff --git a/code/models/rttbDoseBasedModels.h b/code/models/rttbDoseBasedModels.h new file mode 100644 index 0000000..ff8b6ce --- /dev/null +++ b/code/models/rttbDoseBasedModels.h @@ -0,0 +1,43 @@ +// ----------------------------------------------------------------------- +// 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: 798 $ (last changed revision) +// @date $Date: 2014-10-10 15:27:43 +0200 (Fr, 10 Okt 2014) $ (last change date) +// @author $Author: zhangl $ (last changed by) +*/ +#include "rttbBaseType.h" +#include "rttbBaseTypeModels.h" + + +namespace rttb +{ + namespace models + { + /*! @brief Calculate biological LinearQuadratic Model of a dose + @details LQ = exp(-(alpha*d+beta*d^2)) + @param dose + @param alpha + @param beta + @pre dose>=0 + @pre alpha>=0 + @pre beta>=0 + @return The LQ value + @exception InvalidParameterException Thrown if parameters were not set correctly. + */ + BioModelValueType calcLQ(const DoseTypeGy dose, const DoseCalcType alpha, + const DoseCalcType beta); + } +} \ No newline at end of file diff --git a/code/models/rttbDvhBasedModels.h b/code/models/rttbDvhBasedModels.h index c10b073..d7c9d17 100644 --- a/code/models/rttbDvhBasedModels.h +++ b/code/models/rttbDvhBasedModels.h @@ -1,50 +1,49 @@ #include #include "rttbDVH.h" #include "rttbBaseType.h" namespace rttb{ namespace models{ typedef core::DVH::DataDifferentialType DataDifferentialType; typedef core::DVH::DVHPointer DVHPointer; typedef std::map BEDDVHType; typedef std::map LQEDDVHType; - typedef double BioModelParameterType; /*! @brief Get Equivalent Uniform Dose (EUD) @pre dvh data differential is not empty, @pre aA is not zero, @return Return calculated EUD value, @exception InvalidParameterException Thrown if parameters were not set correctly. */ DoseStatisticType getEUD(const DVHPointer dvh, const DoseCalcType aA); /*! @brief Calculate Biological Effective/Equivalent Dose (BED) of dvh @param relativeVolume default false-> the corresponding volume value is the voxel number of the dose bin; if true-> the corresponding volume value is the relative volume % between 0 and 1, (the voxel number of this dose bin)/(number of voxels) @pre dvh should be an accumulated dvh of all fractions, not a single fraction dvh @pre dvh data differential is not empty @pre alpha_beta > 0 @pre numberOfFractions > 1 @return Return map: keys are BEDi in Gy, values are the volume of the dose bin @exception InvalidParameterException Thrown if parameters were not set correctly. */ BEDDVHType calcBEDDVH(const DVHPointer dvh, const int numberOfFractions, const DoseCalcType alpha_beta, const bool relativeVolume=false); /*! @brief Calculate Linear-quadratic equivalent dose for 2-Gy (LQED2) of dvh @param relativeVolume default false-> the corresponding volume value is the voxel number of the dose bin; if true-> the corresponding volume value is the relative volume % between 0 and 1, (the voxel number of this dose bin)/(number of voxels) @pre dvh should be an accumulated dvh of all fractions, not a single fraction dvh @pre dvh data differential is not empty @pre alpha_beta > 0 @pre numberOfFractions > 1 @return Return map: keys are LQED2 in Gy, values are the volume of the dose bin; return empty map if not initialized @exception InvalidParameterException Thrown if parameters were not set correctly. */ LQEDDVHType calcLQED2DVH(const DVHPointer dvh, const int numberOfFractions, const DoseCalcType alpha_beta, const bool relativeVolume=false); } } \ No newline at end of file diff --git a/code/models/rttbLQModelAccessor.cpp b/code/models/rttbLQModelAccessor.cpp new file mode 100644 index 0000000..7665d32 --- /dev/null +++ b/code/models/rttbLQModelAccessor.cpp @@ -0,0 +1,65 @@ +// ----------------------------------------------------------------------- +// 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: 793 $ (last changed revision) +// @date $Date: 2014-10-10 10:24:45 +0200 (Fr, 10 Okt 2014) $ (last change date) +// @author $Author: hentsch $ (last changed by) +*/ + +#include "rttbLQModelAccessor.h" +#include "rttbDoseBasedModels.h" +#include "rttbInvalidDoseException.h" + +namespace rttb +{ + namespace models + { + LQModelAccessor::~LQModelAccessor() + { + + } + + LQModelAccessor::LQModelAccessor(DoseAccessorPointer dose, BioModelParamType alpha, BioModelParamType beta) : + _dose(dose), _alpha(alpha), _beta(beta) + { + if (_dose == NULL) + { + throw core::InvalidDoseException("Dose is NULL"); + } + + assembleGeometricInfo(); + } + + BioModelValueType LQModelAccessor::getBioModelValueAt(const VoxelGridID aID) const + { + return calcLQ(_dose->getDoseAt(aID), _alpha, _beta); + } + + BioModelValueType LQModelAccessor::getBioModelValueAt(const VoxelGridIndex3D& aIndex) const + { + return calcLQ(_dose->getDoseAt(aIndex), _alpha, _beta); + } + + bool LQModelAccessor::assembleGeometricInfo() + { + _geoInfo = _dose->getGeometricInfo(); + + return true; + + } + } +} + diff --git a/code/models/rttbLQModelAccessor.h b/code/models/rttbLQModelAccessor.h new file mode 100644 index 0000000..f569d65 --- /dev/null +++ b/code/models/rttbLQModelAccessor.h @@ -0,0 +1,79 @@ +// ----------------------------------------------------------------------- +// 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: 793 $ (last changed revision) +// @date $Date: 2014-10-10 10:24:45 +0200 (Fr, 10 Okt 2014) $ (last change date) +// @author $Author: hentsch $ (last changed by) +*/ +#ifndef __LQ_MODEL_ACCESSOR_H +#define __LQ_MODEL_ACCESSOR_H + +#include "rttbBioModelAccessorInterface.h" +#include "rttbDoseAccessorInterface.h" +#include "rttbBaseTypeModels.h" + +namespace rttb +{ + namespace models + { + /*! @class LQModelAccessor + @brief This class gives access to the LQ Model information in an image + */ + class LQModelAccessor: public core::BioModelAccessorInterface + { + public: + typedef core::DoseAccessorInterface::DoseAccessorPointer DoseAccessorPointer; + private: + DoseAccessorPointer _dose; + BioModelParamType _alpha; + BioModelParamType _beta; + + IDType _bioModelUID; + + LQModelAccessor(); + + /*! @brief get all required data from the dose geometric info + */ + bool assembleGeometricInfo(); + + + public: + ~LQModelAccessor(); + + /*! @brief Constructor. + @pre dose must be a valid instance (and not null) + @exception InvalidDoseException if _dose is NULL + */ + LQModelAccessor(DoseAccessorPointer dose, BioModelParamType alpha, BioModelParamType beta); + + /*! @brief returns the LQ Model value for an id + */ + BioModelValueType getBioModelValueAt(const VoxelGridID aID) const; + + /*! @brief returns the LQ Model value for an index + */ + BioModelValueType getBioModelValueAt(const VoxelGridIndex3D& aIndex) const; + + const IDType getBioModelUID() const + { + return _bioModelUID; + }; + + }; + } +} + +#endif diff --git a/code/models/rttbTCPLQModel.cpp b/code/models/rttbTCPLQModel.cpp index 3c69723..71ba16e 100644 --- a/code/models/rttbTCPLQModel.cpp +++ b/code/models/rttbTCPLQModel.cpp @@ -1,273 +1,273 @@ // ----------------------------------------------------------------------- // 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 #define _USE_MATH_DEFINES #include #include #include #include #include "rttbTCPLQModel.h" #include "rttbDvhBasedModels.h" #include "rttbBaseTypeModels.h" #include "rttbIntegration.h" #include "rttbInvalidParameterException.h" #include "rttbExceptionMacros.h" namespace rttb { namespace models { TCPLQModel::TCPLQModel(): TCPModel(), _alphaMean(0), _alphaVariance(0), _alpha_beta(0), _rho(0) {} TCPLQModel::TCPLQModel(DVHPointer aDVH, BioModelParamType aAlphaMean, BioModelParamType aBeta, BioModelParamType aRho, int aNumberOfFractions): TCPModel(aDVH, aNumberOfFractions), _alphaMean(aAlphaMean), _alphaVariance(0), _alpha_beta(aAlphaMean / aBeta), _rho(aRho) {} TCPLQModel::TCPLQModel(DVHPointer aDVH, BioModelParamType aRho, int aNumberOfFractions, BioModelParamType aAlpha_Beta, BioModelParamType aAlphaMean, BioModelParamType aAlphaVariance): TCPModel(aDVH, aNumberOfFractions), _alphaMean(aAlphaMean), _alphaVariance(aAlphaVariance), _alpha_beta(aAlpha_Beta), _rho(aRho) {} void TCPLQModel::setParameters(const BioModelParamType aAlphaMean, const BioModelParamType aAlpha_Beta, const BioModelParamType aRho, const BioModelParamType aAlphaVariance) { _alphaMean = aAlphaMean; _alphaVariance = aAlphaVariance; _alpha_beta = aAlpha_Beta; _rho = aRho; //reset _value, because parameters have changed. _value = 0; } void TCPLQModel::setAlpha(const BioModelParamType aAlphaMean, const BioModelParamType aAlphaVariance) { _alphaVariance = aAlphaVariance; _alphaMean = aAlphaMean; } void TCPLQModel::setAlphaBeta(const BioModelParamType aAlpha_Beta) { _alpha_beta = aAlpha_Beta; } void TCPLQModel::setRho(const BioModelParamType aRho) { _rho = aRho; } - const BioModelParamType TCPLQModel::getAlpahBeta() + const BioModelParamType TCPLQModel::getAlphaBeta() { return _alpha_beta; } const BioModelParamType TCPLQModel::getAlphaMean() { return _alphaMean; } const BioModelParamType TCPLQModel::getAlphaVariance() { return _alphaVariance; } const BioModelParamType TCPLQModel::getRho() { return _rho; } long double TCPLQModel::calcTCPi(BioModelParamType aRho, BioModelParamType aAlphaMean, double vj, double bedj) { return exp(-aRho * vj * exp(-aAlphaMean * bedj)); } long double TCPLQModel::calcTCP(std::map aBEDDVH, BioModelParamType aRho, BioModelParamType aAlphaMean, double aDeltaV) { std::map::iterator it; long double tcp = 1; for (it = aBEDDVH.begin(); it != aBEDDVH.end(); ++it) { long double tcpi = this->calcTCPi(aRho, aAlphaMean, (*it).second * aDeltaV, (*it).first); tcp = tcp * tcpi; } return tcp; } long double TCPLQModel::calcTCPAlphaNormalDistribution(std::map aBEDDVH, BioModelParamType aRho, BioModelParamType aAlphaMean, BioModelParamType aAlphaVariance, double aDeltaV) { std::map::iterator it; std::vector volumeV2; std::vector bedV2; int i = 0; for (it = aBEDDVH.begin(); it != aBEDDVH.end(); ++it) { volumeV2.push_back((*it).second * aDeltaV); bedV2.push_back((*it).first); i++; } struct TcpParams params = {aAlphaMean, aAlphaVariance, aRho, volumeV2, bedV2}; double result = integrateTCP(0, params); if (result == -100) { std::cerr << "Integration error!\n"; return -1; } else { long double tcp = 1 / (pow(2 * M_PI, 0.5) * _alphaVariance) * result; return tcp; } } BioModelValueType TCPLQModel::calcModel(const double doseFactor) { core::DVH variantDVH = core::DVH(_dvh->getDataDifferential(), (DoseTypeGy)(_dvh->getDeltaD() * doseFactor), _dvh->getDeltaV(), "temporary", "temporary"); boost::shared_ptr spDVH = boost::make_shared(variantDVH); BioModelValueType value = 0; if (_alphaVariance == 0) { - if (_alphaMean <= 0 && _alpha_beta <= 0 && _rho <= 0) + if (_alphaMean <= 0 || _alpha_beta <= 0 || _rho <= 0) { throw core::InvalidParameterException("Parameter invalid: alpha, alpha/beta, rho and number of fractions must >0!"); } if (_numberOfFractions <= 1) { throw core::InvalidParameterException("Parameter invalid: numberOfFractions must be >1! The dvh should be an accumulated-dvh of all fractions, not a single fraction-dvh!"); } std::map dataBED = calcBEDDVH(spDVH, _numberOfFractions, _alpha_beta); value = (BioModelValueType)this->calcTCP(dataBED, _rho, _alphaMean, variantDVH.getDeltaV()); return value; } //if alpha normal distribution else { - if (this->_alpha_beta <= 0 && this->_alphaMean <= 0 && this->_alphaVariance < 0 && _rho <= 0 && _numberOfFractions <= 0) + if (this->_alpha_beta <= 0 || this->_alphaMean <= 0 || this->_alphaVariance < 0 || _rho <= 0) { throw core::InvalidParameterException("Parameter invalid: alpha/beta, alphaMean, rho and number of fractions must >0!"); } if (_numberOfFractions <= 1) { throw core::InvalidParameterException("Parameter invalid: numberOfFractions must be >1! The dvh should be an accumulated-dvh of all fractions, not a single fraction-dvh!"); } std::map dataBED = calcBEDDVH(spDVH, _numberOfFractions, _alpha_beta); value = (BioModelValueType)(this->calcTCPAlphaNormalDistribution(dataBED, _rho, _alphaMean, _alphaVariance, variantDVH.getDeltaV())); return value; } } void TCPLQModel::setParameterVector(const ParamVectorType& aParameterVector) { if (aParameterVector.size() != 4) { throw core::InvalidParameterException("Parameter invalid: aParameterVector.size must be 4! "); } else { _alphaMean = aParameterVector.at(0); _alphaVariance = aParameterVector.at(1); _alpha_beta = aParameterVector.at(2); _rho = aParameterVector.at(3); } } void TCPLQModel::setParameterByID(const int aParamId, const BioModelParamType aValue) { if (aParamId == 0) { _alphaMean = aValue; } else if (aParamId == 1) { _alphaVariance = aValue; } else if (aParamId == 2) { _alpha_beta = aValue; } else if (aParamId == 3) { _rho = aValue; } else { throw core::InvalidParameterException("Parameter invalid: aParamID must be 0(alphaMean) or 1(alphaVariance) or 2(alpha_beta) or 3(rho)! "); } } const int TCPLQModel::getParameterID(const std::string& aParamName) const { if (aParamName == "alphaMean") { return 0; } else if (aParamName == "alphaVariance") { return 1; } else if (aParamName == "alpha_beta") { return 2; } else if (aParamName == "rho") { return 3; } else { rttbExceptionMacro(core::InvalidParameterException, << "Parameter name " << aParamName << " invalid: it should be alphaMean or alphaVariance or alpha_beta or rho!"); } } }//end namespace models }//end namespace rttb diff --git a/code/models/rttbTCPLQModel.h b/code/models/rttbTCPLQModel.h index 78a2abe..560b33e 100644 --- a/code/models/rttbTCPLQModel.h +++ b/code/models/rttbTCPLQModel.h @@ -1,163 +1,163 @@ // ----------------------------------------------------------------------- // 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 __TCP_LQ_MODEL_H #define __TCP_LQ_MODEL_H #include #include #include #include "rttbBaseType.h" #include "rttbTCPModel.h" #include "rttbBaseTypeModels.h" namespace rttb { namespace models { /*! @class TCPLQModel @brief This class represents a TCP(Tumor Control Probability) LQ model (Nahum and Sanchez-Nieto 2001, Hall and Giaccia 2006) @see TCPModel */ class TCPLQModel: public TCPModel { public: typedef TCPModel::ParamVectorType ParamVectorType; typedef TCPModel::DVHPointer DVHPointer; private: /*! @brief Calculate intermediate tcp using alpha constant. This is a helper function for calcTCP() @see calcTCP */ long double calcTCPi(BioModelParamType aRho, BioModelParamType aAlphaMean, double vj, double bedj); /*! @brief Calculate tcp using alpha constant. */ long double calcTCP(std::map aBEDDVH, BioModelParamType aRho, BioModelParamType aAlphaMean, double aDeltaV); /*! @brief Calculate tcp using a normal distribution for alpha. */ long double calcTCPAlphaNormalDistribution(std::map aBEDDVH, BioModelParamType aRho, BioModelParamType aAlphaMean, BioModelParamType aAlphaVariance, double aDeltaV); protected: BioModelParamType _alphaMean; BioModelParamType _alphaVariance; BioModelParamType _alpha_beta; /*! Roh is the initial clonogenic cell density */ BioModelParamType _rho; /*! @brief Calculate the model value @param doseFactor scaling factor for prescribed dose. The model calculation will use the dvh with each di=old di*doseFactor. @pre _alphaMean >0 @pre _alphaVariance >= 0 @pre _alpha_beta > 0 @pre _rho > 0 @pre _numberOfFractions > 1 @exception InvalidParameterException Thrown if parameters were not set correctly. */ BioModelValueType calcModel(const double doseFactor = 1); public: TCPLQModel(); /*! @brief Constructor initializes member variables with given parameters. @pre aAlphaMean >0 @pre aBeta > 0 @pre aRho > 0 @pre aNumberOfFractions > 1 */ TCPLQModel(DVHPointer aDVH, BioModelParamType aAlphaMean, BioModelParamType aBeta, BioModelParamType aRho, int aNumberOfFractions); /*! @brief Constructor for alpha distribution initializes member variables with given parameters. @pre aAlphaMean >0 @pre aAlphaVariance >0 @pre aAlpha_Beta > 0 @pre aRho > 0 @pre aNumberOfFractions > 1 */ TCPLQModel(DVHPointer aDVH, BioModelParamType aRho, int aNumberOfFractions, BioModelParamType aAlpha_Beta, BioModelParamType aAlphaMean, BioModelParamType aAlphaVariance); const BioModelParamType getRho(); void setRho(const BioModelParamType aRho); const BioModelParamType getAlphaMean(); const BioModelParamType getAlphaVariance(); /*! @brief The distribution of the parameter alpha, which is characteristic for a population of cells, is described by the its mean and variance. If alpha is constant the variance is 0. @param aAlphaVariance The variance of alpha can be given, the default value is 0 resulting in constant alpha. */ void setAlpha(const BioModelParamType aAlphaMean, const BioModelParamType aAlphaVariance = 0); - const BioModelParamType getAlpahBeta(); + const BioModelParamType getAlphaBeta(); void setAlphaBeta(const BioModelParamType aAlpha_Beta); /*! @brief Set parameters for the TCP model. _value will be reset to 0. @param aAlpha_Beta alpha/beta constant . @param aAlphaMean mean of alpha distribution. @param aAlphaVariance variance of alpha distribution. */ void setParameters(const BioModelParamType aAlphaMean, const BioModelParamType aAlpha_Beta, const BioModelParamType aRho, const BioModelParamType aAlphaVariance = 0); /*! @brief Set parameter with ID. "alphaMean":0,"alphaVariance":1,"alpha_beta":2, "rho":3 @exception InvalidParameterException Thrown if aParamId is not 0 or 1 or 2 or 3. */ virtual void setParameterByID(const int aParamId, const BioModelParamType aValue); /*! @brief Set parameter vector, where index of vector is the parameter id. "alphaMean":0,"alphaVariance":1,"alpha_beta":2, "rho":3 @exception InvalidParameterException Thrown if aParamterVector.size()!=4. */ virtual void setParameterVector(const ParamVectorType& aParameterVector); /*! @brief Get parameter id. "alphaMean":0,"alphaVariance":1,"alpha_beta":2, "rho":3 @return 0 for "alphaMean", 1 for "alphaVariance", 2 for "alpha_beta", 3 for "rho" @exception InvalidParameterException Thrown if aParamName is not alphaMean or alphaVariance or alpha_beta or rho. */ virtual const int getParameterID(const std::string& aParamName) const; }; }//end algorithms }//end rttb #endif diff --git a/testing/examples/RTBioModelExampleTest.cpp b/testing/examples/RTBioModelExampleTest.cpp index 38acdf3..8b03060 100644 --- a/testing/examples/RTBioModelExampleTest.cpp +++ b/testing/examples/RTBioModelExampleTest.cpp @@ -1,374 +1,413 @@ // ----------------------------------------------------------------------- // 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) */ // this file defines the rttbCoreTests for the test driver // and all it expects is that you have a function called RegisterTests #include #include "litCheckMacros.h" #include "rttbBioModel.h" -#include "rttbDVHTxtFileReader.h" +#include "rttbDVHTxtFileReader.h" #include "rttbDVH.h" #include "rttbTCPLQModel.h" #include "rttbNTCPLKBModel.h" #include "rttbNTCPRSModel.h" #include "rttbBioModelScatterPlots.h" #include "rttbBioModelCurve.h" #include "rttbDvhBasedModels.h" #include "rttbDoseIteratorInterface.h" -namespace rttb{ - namespace testing{ +namespace rttb +{ + namespace testing + { - /*! @brief RTBioModelTest. - TCP calculated using a DVH PTV and LQ Model. + /*! @brief RTBioModelTest. + TCP calculated using a DVH PTV and LQ Model. NTCP tested using 3 Normal Tissue DVHs and LKB/RS Model. - Test if calculation in new architecture returns similar results to the + Test if calculation in new architecture returns similar results to the original implementation. WARNING: The values for comparison need to be adjusted if the input files are changed! */ - int RTBioModelExampleTest(int argc, char* argv[] ) + int RTBioModelExampleTest(int argc, char* argv[]) { PREPARE_DEFAULT_TEST_REPORTING; typedef rttb::models::CurveDataType CurveDataType; - typedef std::multimap > ScatterPlotType; + typedef std::multimap > ScatterPlotType; typedef core::DVH::DVHPointer DVHPointer; //increased accuracy requires double values in the calculation (rttbBaseType.h) double toleranceEUD = 1e-5; double tolerance = 1e-7; //ARGUMENTS: 1: ptv dvh file name // 2: normal tissue 1 dvh file name // 3: normal tissue 2 dvh file name // 4: normal tissue 3 dvh file name //...........5: Virtuos MPM_LR_ah dvh lung file name //...........6: Virtuos MPM_LR_ah dvh target file name std::string DVH_FILENAME_PTV; std::string DVH_FILENAME_NT1; std::string DVH_FILENAME_NT2; std::string DVH_FILENAME_NT3; std::string DVH_FILENAME_TV_TEST; std::string DVH_Virtuos_Target; std::string DVH_Virtuos_Lung; - if (argc>1) + if (argc > 1) { DVH_FILENAME_PTV = argv[1]; } - if (argc>2) + + if (argc > 2) { DVH_FILENAME_NT1 = argv[2]; } - if (argc>3) + + if (argc > 3) { DVH_FILENAME_NT2 = argv[3]; } - if (argc>4) + + if (argc > 4) { DVH_FILENAME_NT3 = argv[4]; } - if (argc>5) + + if (argc > 5) { DVH_FILENAME_TV_TEST = argv[5]; } - if(argc>6) + + if (argc > 6) { - DVH_Virtuos_Lung=argv[6]; + DVH_Virtuos_Lung = argv[6]; } - if(argc>7) + + if (argc > 7) { - DVH_Virtuos_Target=argv[7]; + DVH_Virtuos_Target = argv[7]; } - //DVH PTV - rttb::io::other::DVHTxtFileReader dvhReader=rttb::io::other::DVHTxtFileReader(DVH_FILENAME_PTV); - DVHPointer dvhPtr=dvhReader.generateDVH(); + //DVH PTV + rttb::io::other::DVHTxtFileReader dvhReader = rttb::io::other::DVHTxtFileReader(DVH_FILENAME_PTV); + DVHPointer dvhPtr = dvhReader.generateDVH(); - CHECK_CLOSE(6.04759613161786830000e+001,models::getEUD(dvhPtr,10),toleranceEUD); + CHECK_CLOSE(6.04759613161786830000e+001, models::getEUD(dvhPtr, 10), toleranceEUD); - rttb::io::other::DVHTxtFileReader dvhReader_test_tv=rttb::io::other::DVHTxtFileReader(DVH_FILENAME_TV_TEST); - DVHPointer dvh_test_tv=dvhReader_test_tv.generateDVH(); + rttb::io::other::DVHTxtFileReader dvhReader_test_tv = rttb::io::other::DVHTxtFileReader(DVH_FILENAME_TV_TEST); + DVHPointer dvh_test_tv = dvhReader_test_tv.generateDVH(); //test TCP LQ Model models::BioModelParamType alpha = 0.35; models::BioModelParamType beta = 0.023333333333333; models::BioModelParamType roh = 10000000; int numFractions = 2; DoseTypeGy normalizationDose = 68; - rttb::models::TCPLQModel tcplq=rttb::models::TCPLQModel(dvhPtr,alpha, beta, roh, numFractions); - CHECK_EQUAL(alpha,tcplq.getAlphaMean()); - CHECK_EQUAL(alpha/beta,tcplq.getAlpahBeta()); - CHECK_EQUAL(roh,tcplq.getRho()); + rttb::models::TCPLQModel tcplq = rttb::models::TCPLQModel(dvhPtr, alpha, beta, roh, numFractions); + CHECK_EQUAL(alpha, tcplq.getAlphaMean()); + CHECK_EQUAL(alpha / beta, tcplq.getAlphaBeta()); + CHECK_EQUAL(roh, tcplq.getRho()); CHECK_NO_THROW(tcplq.init()); - if(tcplq.init()){ - CHECK_CLOSE(1.00497232941856940000e-127,tcplq.getValue(),tolerance); + + if (tcplq.init()) + { + CHECK_CLOSE(1.00497232941856940000e-127, tcplq.getValue(), tolerance); } - CurveDataType curve=models::getCurveDoseVSBioModel(tcplq,normalizationDose); + CurveDataType curve = models::getCurveDoseVSBioModel(tcplq, normalizationDose); CurveDataType::iterator it; - for(it=curve.begin();it!=curve.end();it++){ - if ((*it).first < 72){ - CHECK_EQUAL(0,(*it).second); + + for (it = curve.begin(); it != curve.end(); it++) + { + if ((*it).first < 72) + { + CHECK_EQUAL(0, (*it).second); } - else if ((*it).first >150){ - CHECK((*it).second>0.9); + else if ((*it).first > 150) + { + CHECK((*it).second > 0.9); } } models::BioModelParamType alphaBeta = 10; - tcplq.setParameters(alpha,alphaBeta,roh,0.08); - CHECK_EQUAL(alpha,tcplq.getAlphaMean()); - CHECK_EQUAL(alphaBeta,tcplq.getAlpahBeta()); - CHECK_EQUAL(roh,tcplq.getRho()); - if(tcplq.init()){ - CHECK_CLOSE(1.84e-005,tcplq.getValue(),tolerance); + tcplq.setParameters(alpha, alphaBeta, roh, 0.08); + CHECK_EQUAL(alpha, tcplq.getAlphaMean()); + CHECK_EQUAL(alphaBeta, tcplq.getAlphaBeta()); + CHECK_EQUAL(roh, tcplq.getRho()); + + if (tcplq.init()) + { + CHECK_CLOSE(1.84e-005, tcplq.getValue(), tolerance); } normalizationDose = 40; - curve=models::getCurveDoseVSBioModel(tcplq,normalizationDose); + curve = models::getCurveDoseVSBioModel(tcplq, normalizationDose); alpha = 1; alphaBeta = 14.5; tcplq.setAlpha(alpha); tcplq.setAlphaBeta(alphaBeta); tcplq.setRho(roh); - CHECK_EQUAL(alpha,tcplq.getAlphaMean()); - CHECK_EQUAL(alphaBeta,tcplq.getAlpahBeta()); - CHECK_EQUAL(roh,tcplq.getRho()); - if(tcplq.init()){ + CHECK_EQUAL(alpha, tcplq.getAlphaMean()); + CHECK_EQUAL(alphaBeta, tcplq.getAlphaBeta()); + CHECK_EQUAL(roh, tcplq.getRho()); + + if (tcplq.init()) + { CHECK_CLOSE(0.954885, tcplq.getValue(), toleranceEUD); } alpha = 0.9; alphaBeta = 1; tcplq.setAlpha(alpha); tcplq.setAlphaBeta(alphaBeta); tcplq.setRho(roh); - CHECK_EQUAL(alpha,tcplq.getAlphaMean()); - CHECK_EQUAL(alphaBeta,tcplq.getAlpahBeta()); - CHECK_EQUAL(roh,tcplq.getRho()); - if(tcplq.init()){ - CHECK_EQUAL(1,tcplq.getValue()); + CHECK_EQUAL(alpha, tcplq.getAlphaMean()); + CHECK_EQUAL(alphaBeta, tcplq.getAlphaBeta()); + CHECK_EQUAL(roh, tcplq.getRho()); + + if (tcplq.init()) + { + CHECK_EQUAL(1, tcplq.getValue()); } //TCP LQ Test alpha = 0.3; beta = 0.03; roh = 10000000; numFractions = 20; - rttb::models::TCPLQModel tcplq_test=rttb::models::TCPLQModel(dvh_test_tv,alpha, beta, roh, numFractions); - CHECK_EQUAL(alpha,tcplq_test.getAlphaMean()); - CHECK_EQUAL(alpha/beta,tcplq_test.getAlpahBeta()); - CHECK_EQUAL(roh,tcplq_test.getRho()); + rttb::models::TCPLQModel tcplq_test = rttb::models::TCPLQModel(dvh_test_tv, alpha, beta, roh, numFractions); + CHECK_EQUAL(alpha, tcplq_test.getAlphaMean()); + CHECK_EQUAL(alpha / beta, tcplq_test.getAlphaBeta()); + CHECK_EQUAL(roh, tcplq_test.getRho()); CHECK_NO_THROW(tcplq_test.init()); - if(tcplq_test.init()){ - CHECK_CLOSE(9.79050278878883180000e-001,tcplq_test.getValue(),tolerance); + + if (tcplq_test.init()) + { + CHECK_CLOSE(9.79050278878883180000e-001, tcplq_test.getValue(), tolerance); } + normalizationDose = 60; - curve=models::getCurveDoseVSBioModel(tcplq_test,normalizationDose); + curve = models::getCurveDoseVSBioModel(tcplq_test, normalizationDose); //DVH HT 1 - rttb::io::other::DVHTxtFileReader dvhReader2=rttb::io::other::DVHTxtFileReader(DVH_FILENAME_NT1); - DVHPointer dvhPtr2=dvhReader2.generateDVH(); + rttb::io::other::DVHTxtFileReader dvhReader2 = rttb::io::other::DVHTxtFileReader(DVH_FILENAME_NT1); + DVHPointer dvhPtr2 = dvhReader2.generateDVH(); - CHECK_CLOSE(1.07920836034015810000e+001,models::getEUD(dvhPtr2,10),toleranceEUD); + CHECK_CLOSE(1.07920836034015810000e+001, models::getEUD(dvhPtr2, 10), toleranceEUD); //test RTNTCPLKBModel - rttb::models::NTCPLKBModel lkb=rttb::models::NTCPLKBModel(); + rttb::models::NTCPLKBModel lkb = rttb::models::NTCPLKBModel(); models::BioModelParamType aVal = 10; models::BioModelParamType mVal = 0.16; models::BioModelParamType d50Val = 55; - CHECK_EQUAL(0,lkb.getA()); - CHECK_EQUAL(0,lkb.getM()); - CHECK_EQUAL(0,lkb.getD50()); + CHECK_EQUAL(0, lkb.getA()); + CHECK_EQUAL(0, lkb.getM()); + CHECK_EQUAL(0, lkb.getD50()); lkb.setDVH(dvhPtr2); - CHECK_EQUAL(dvhPtr2,lkb.getDVH()); + CHECK_EQUAL(dvhPtr2, lkb.getDVH()); lkb.setA(aVal); - CHECK_EQUAL(aVal,lkb.getA()); + CHECK_EQUAL(aVal, lkb.getA()); lkb.setM(mVal); - CHECK_EQUAL(mVal,lkb.getM()); + CHECK_EQUAL(mVal, lkb.getM()); lkb.setD50(d50Val); - CHECK_EQUAL(d50Val,lkb.getD50()); + CHECK_EQUAL(d50Val, lkb.getD50()); CHECK_NO_THROW(lkb.init()); - if(lkb.init()){ - CHECK_CLOSE(2.53523522831366570000e-007,lkb.getValue(),tolerance); + + if (lkb.init()) + { + CHECK_CLOSE(2.53523522831366570000e-007, lkb.getValue(), tolerance); } //test RTNTCPRSModel - rttb::models::NTCPRSModel rs=rttb::models::NTCPRSModel(); + rttb::models::NTCPRSModel rs = rttb::models::NTCPRSModel(); models::BioModelParamType gammaVal = 1.7; models::BioModelParamType sVal = 1; - CHECK_EQUAL(0,rs.getGamma()); - CHECK_EQUAL(0,rs.getS()); - CHECK_EQUAL(0,rs.getD50()); + CHECK_EQUAL(0, rs.getGamma()); + CHECK_EQUAL(0, rs.getS()); + CHECK_EQUAL(0, rs.getD50()); rs.setDVH(dvhPtr2); - CHECK_EQUAL(dvhPtr2,rs.getDVH()); + CHECK_EQUAL(dvhPtr2, rs.getDVH()); rs.setD50(d50Val); - CHECK_EQUAL(d50Val,rs.getD50()); + CHECK_EQUAL(d50Val, rs.getD50()); rs.setGamma(gammaVal); - CHECK_EQUAL(gammaVal,rs.getGamma()); + CHECK_EQUAL(gammaVal, rs.getGamma()); rs.setS(sVal); - CHECK_EQUAL(sVal,rs.getS()); + CHECK_EQUAL(sVal, rs.getS()); CHECK_NO_THROW(rs.init()); - if(rs.init()){ - CHECK_CLOSE(3.70385888626145740000e-009,rs.getValue(),tolerance); + + if (rs.init()) + { + CHECK_CLOSE(3.70385888626145740000e-009, rs.getValue(), tolerance); } //DVH HT 2 - rttb::io::other::DVHTxtFileReader dvhReader3=rttb::io::other::DVHTxtFileReader(DVH_FILENAME_NT2); - DVHPointer dvhPtr3=dvhReader3.generateDVH(); - CHECK_CLOSE(1.26287047025885110000e+001,models::getEUD(dvhPtr3,10),toleranceEUD); + rttb::io::other::DVHTxtFileReader dvhReader3 = rttb::io::other::DVHTxtFileReader(DVH_FILENAME_NT2); + DVHPointer dvhPtr3 = dvhReader3.generateDVH(); + CHECK_CLOSE(1.26287047025885110000e+001, models::getEUD(dvhPtr3, 10), toleranceEUD); //test RTNTCPLKBModel aVal = 10; mVal = 0.16; d50Val = 55; lkb.setDVH(dvhPtr3); - CHECK_EQUAL(dvhPtr3,lkb.getDVH()); + CHECK_EQUAL(dvhPtr3, lkb.getDVH()); lkb.setA(aVal); - CHECK_EQUAL(aVal,lkb.getA()); + CHECK_EQUAL(aVal, lkb.getA()); lkb.setM(mVal); - CHECK_EQUAL(mVal,lkb.getM()); + CHECK_EQUAL(mVal, lkb.getM()); lkb.setD50(d50Val); - CHECK_EQUAL(d50Val,lkb.getD50()); - if(lkb.init()){ - CHECK_CLOSE(7.36294657754956700000e-007,lkb.getValue(),tolerance); + CHECK_EQUAL(d50Val, lkb.getD50()); + + if (lkb.init()) + { + CHECK_CLOSE(7.36294657754956700000e-007, lkb.getValue(), tolerance); } //test RTNTCPRSModel - rs=rttb::models::NTCPRSModel(); + rs = rttb::models::NTCPRSModel(); gammaVal = 1.7; sVal = 1; - CHECK_EQUAL(0,rs.getGamma()); - CHECK_EQUAL(0,rs.getS()); - CHECK_EQUAL(0,rs.getD50()); + CHECK_EQUAL(0, rs.getGamma()); + CHECK_EQUAL(0, rs.getS()); + CHECK_EQUAL(0, rs.getD50()); rs.setDVH(dvhPtr3); - CHECK_EQUAL(dvhPtr3,rs.getDVH()); + CHECK_EQUAL(dvhPtr3, rs.getDVH()); rs.setD50(d50Val); - CHECK_EQUAL(d50Val,rs.getD50()); + CHECK_EQUAL(d50Val, rs.getD50()); rs.setGamma(gammaVal); - CHECK_EQUAL(gammaVal,rs.getGamma()); + CHECK_EQUAL(gammaVal, rs.getGamma()); rs.setS(sVal); - CHECK_EQUAL(sVal,rs.getS()); - if(rs.init()){ - CHECK_CLOSE(1.76778795490939440000e-007,rs.getValue(),tolerance); + CHECK_EQUAL(sVal, rs.getS()); + + if (rs.init()) + { + CHECK_CLOSE(1.76778795490939440000e-007, rs.getValue(), tolerance); } //DVH HT 3 - rttb::io::other::DVHTxtFileReader dvhReader4=rttb::io::other::DVHTxtFileReader(DVH_FILENAME_NT3); - DVHPointer dvhPtr4=dvhReader4.generateDVH(); - CHECK_CLOSE(2.18212982041056310000e+001,models::getEUD(dvhPtr4,10),toleranceEUD); + rttb::io::other::DVHTxtFileReader dvhReader4 = rttb::io::other::DVHTxtFileReader(DVH_FILENAME_NT3); + DVHPointer dvhPtr4 = dvhReader4.generateDVH(); + CHECK_CLOSE(2.18212982041056310000e+001, models::getEUD(dvhPtr4, 10), toleranceEUD); //test RTNTCPLKBModel aVal = 10; mVal = 0.16; d50Val = 55; lkb.setDVH(dvhPtr4); - CHECK_EQUAL(dvhPtr4,lkb.getDVH()); + CHECK_EQUAL(dvhPtr4, lkb.getDVH()); lkb.setA(aVal); - CHECK_EQUAL(aVal,lkb.getA()); + CHECK_EQUAL(aVal, lkb.getA()); lkb.setM(mVal); - CHECK_EQUAL(mVal,lkb.getM()); + CHECK_EQUAL(mVal, lkb.getM()); lkb.setD50(d50Val); - CHECK_EQUAL(d50Val,lkb.getD50()); - if(lkb.init()){ - CHECK_CLOSE(8.15234192641929420000e-005,lkb.getValue(),tolerance); + CHECK_EQUAL(d50Val, lkb.getD50()); + + if (lkb.init()) + { + CHECK_CLOSE(8.15234192641929420000e-005, lkb.getValue(), tolerance); } //test RTNTCPRSModel - rs=rttb::models::NTCPRSModel(); + rs = rttb::models::NTCPRSModel(); gammaVal = 1.7; sVal = 1; - CHECK_EQUAL(0,rs.getGamma()); - CHECK_EQUAL(0,rs.getS()); - CHECK_EQUAL(0,rs.getD50()); + CHECK_EQUAL(0, rs.getGamma()); + CHECK_EQUAL(0, rs.getS()); + CHECK_EQUAL(0, rs.getD50()); rs.setDVH(dvhPtr4); - CHECK_EQUAL(dvhPtr4,rs.getDVH()); + CHECK_EQUAL(dvhPtr4, rs.getDVH()); rs.setD50(d50Val); - CHECK_EQUAL(d50Val,rs.getD50()); + CHECK_EQUAL(d50Val, rs.getD50()); rs.setGamma(gammaVal); - CHECK_EQUAL(gammaVal,rs.getGamma()); + CHECK_EQUAL(gammaVal, rs.getGamma()); rs.setS(sVal); - CHECK_EQUAL(sVal,rs.getS()); - if(rs.init()){ - CHECK_CLOSE(2.02607985020919480000e-004,rs.getValue(),tolerance); + CHECK_EQUAL(sVal, rs.getS()); + + if (rs.init()) + { + CHECK_CLOSE(2.02607985020919480000e-004, rs.getValue(), tolerance); } //test using Virtuos Pleuramesotheliom MPM_LR_ah - //DVH PTV + //DVH PTV - rttb::io::other::DVHTxtFileReader dR_Target=rttb::io::other::DVHTxtFileReader(DVH_Virtuos_Target); - DVHPointer dvhPtrTarget=dR_Target.generateDVH(); + rttb::io::other::DVHTxtFileReader dR_Target = rttb::io::other::DVHTxtFileReader(DVH_Virtuos_Target); + DVHPointer dvhPtrTarget = dR_Target.generateDVH(); - rttb::io::other::DVHTxtFileReader dR_Lung=rttb::io::other::DVHTxtFileReader(DVH_Virtuos_Lung); - DVHPointer dvhPtrLung=dR_Lung.generateDVH(); + rttb::io::other::DVHTxtFileReader dR_Lung = rttb::io::other::DVHTxtFileReader(DVH_Virtuos_Lung); + DVHPointer dvhPtrLung = dR_Lung.generateDVH(); //test TCP LQ Model models::BioModelParamType alphaMean = 0.34; - models::BioModelParamType alphaVarianz=0.02; + models::BioModelParamType alphaVarianz = 0.02; models::BioModelParamType alpha_beta = 28; models::BioModelParamType rho = 1200; int numFractionsVirtuos = 27; - rttb::models::TCPLQModel tcplqVirtuos=rttb::models::TCPLQModel(dvhPtrTarget,rho,numFractionsVirtuos,alpha_beta,alphaMean,alphaVarianz); - if(tcplqVirtuos.init()) + rttb::models::TCPLQModel tcplqVirtuos = rttb::models::TCPLQModel(dvhPtrTarget, rho, numFractionsVirtuos, alpha_beta, + alphaMean, alphaVarianz); + + if (tcplqVirtuos.init()) { - CHECK_CLOSE(0.8894,tcplqVirtuos.getValue(),1e-4); + CHECK_CLOSE(0.8894, tcplqVirtuos.getValue(), 1e-4); } - models::BioModelParamType d50Mean=20; - models::BioModelParamType d50Varianz=2; - models::BioModelParamType m=0.36; - models::BioModelParamType a=1.06; + models::BioModelParamType d50Mean = 20; + models::BioModelParamType d50Varianz = 2; + models::BioModelParamType m = 0.36; + models::BioModelParamType a = 1.06; - rttb::models::NTCPLKBModel lkbVirtuos=rttb::models::NTCPLKBModel(dvhPtrLung,d50Mean,m,a); - if(lkbVirtuos.init()){ - CHECK_CLOSE(0.0397,lkbVirtuos.getValue(),1e-4); + rttb::models::NTCPLKBModel lkbVirtuos = rttb::models::NTCPLKBModel(dvhPtrLung, d50Mean, m, a); + + if (lkbVirtuos.init()) + { + CHECK_CLOSE(0.0397, lkbVirtuos.getValue(), 1e-4); } RETURN_AND_REPORT_TEST_SUCCESS; } }//testing }//rttb \ No newline at end of file diff --git a/testing/examples/RTBioModelScatterPlotExampleTest.cpp b/testing/examples/RTBioModelScatterPlotExampleTest.cpp index ac78bd8..4cefe2f 100644 --- a/testing/examples/RTBioModelScatterPlotExampleTest.cpp +++ b/testing/examples/RTBioModelScatterPlotExampleTest.cpp @@ -1,480 +1,480 @@ // ----------------------------------------------------------------------- // 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) */ // this file defines the rttbCoreTests for the test driver // and all it expects is that you have a function called RegisterTests #include #include "litCheckMacros.h" #include "rttbBioModel.h" #include "rttbDVHTxtFileReader.h" #include "rttbDVH.h" #include "rttbTCPLQModel.h" #include "rttbNTCPLKBModel.h" #include "rttbNTCPRSModel.h" #include "rttbBioModelScatterPlots.h" #include "rttbBioModelCurve.h" #include "rttbDvhBasedModels.h" #include "../models/rttbScatterTester.h" #include "rttbInvalidParameterException.h" namespace rttb { namespace testing { /*! @brief RTBioModelScatterPlotExampleTest. calculating Curves and Scatterplots for TCP and NTCP models. The values on curve and scatterplot need to be similar for similar dose values. The range of difference is given by the variance used to generate the scatter. WARNING: The values for comparison need to be adjusted if the input files are changed! */ int RTBioModelScatterPlotExampleTest(int argc, char* argv[]) { PREPARE_DEFAULT_TEST_REPORTING; typedef rttb::models::CurveDataType CurveDataType; typedef rttb::models::ScatterPlotType ScatterPlotType; typedef core::DVH::DVHPointer DVHPointer; //increased accuracy requires double values in the calculation (rttbBaseType.h) double toleranceEUD = 1e-5; double tolerance = 1e-7; //ARGUMENTS: 1: ptv dvh file name // 2: normal tissue 1 dvh file name // 3: TV dvh file name std::string DVH_FILENAME_PTV; std::string DVH_FILENAME_NT1; std::string DVH_FILENAME_TV_TEST; if (argc > 1) { DVH_FILENAME_PTV = argv[1]; } if (argc > 2) { DVH_FILENAME_NT1 = argv[2]; } if (argc > 3) { DVH_FILENAME_TV_TEST = argv[3]; } //DVH PTV rttb::io::other::DVHTxtFileReader dvhReader = rttb::io::other::DVHTxtFileReader(DVH_FILENAME_PTV); DVHPointer dvhPtr = dvhReader.generateDVH(); rttb::io::other::DVHTxtFileReader dvhReader_test_tv = rttb::io::other::DVHTxtFileReader( DVH_FILENAME_TV_TEST); DVHPointer dvh_test_tv = dvhReader_test_tv.generateDVH(); //test TCP LQ Model models::BioModelParamType alpha = 0.35; models::BioModelParamType beta = 0.023333333333333; models::BioModelParamType roh = 10000000; int numFractions = 2; DoseTypeGy normalizationDose = 68; rttb::models::TCPLQModel tcplq = rttb::models::TCPLQModel(dvhPtr, alpha, beta, roh, numFractions); CHECK_NO_THROW(tcplq.init()); CurveDataType curve = models::getCurveDoseVSBioModel(tcplq, normalizationDose); CHECK_NO_THROW(models::getScatterPlotVary1Parameter(tcplq, 0, alpha, 0, normalizationDose, 100)); ScatterPlotType tcpScatter = models::getScatterPlotVary1Parameter(tcplq, 0, alpha, 0, normalizationDose, 100); CHECK_EQUAL(100, tcpScatter.size()); ScatterTester scatterCompare(curve, tcpScatter); CHECK_TESTER(scatterCompare); //test also with other parameter tcpScatter = models::getScatterPlotVary1Parameter(tcplq, 3, roh, 0, normalizationDose, 100); scatterCompare.setCompareScatter(tcpScatter); CHECK_TESTER(scatterCompare); std::vector paramIdVec; models::BioModel::ParamVectorType meanVec; models::BioModel::ParamVectorType meanVecTest; meanVecTest.push_back(alpha); models::BioModel::ParamVectorType varianceVec; //"alphaMean":0,"alphaVariance":1,"alpha_beta":2, "rho":3 paramIdVec.push_back(0); meanVec.push_back(tcplq.getAlphaMean()); varianceVec.push_back(0); //setting parameter 1 will change the resulting scatter plot dramatically - is it meant to? //this is unexpected since the value was taken from the original calculation //paramIdVec.push_back(1); meanVec.push_back(tcplq.getAlphaVariance()); varianceVec.push_back(0); paramIdVec.push_back(2); - meanVec.push_back(tcplq.getAlpahBeta()); + meanVec.push_back(tcplq.getAlphaBeta()); varianceVec.push_back(0); paramIdVec.push_back(3); meanVec.push_back(tcplq.getRho()); varianceVec.push_back(0); CHECK_THROW_EXPLICIT(models::getScatterPlotVaryParameters(tcplq, paramIdVec, meanVecTest, varianceVec, normalizationDose, 50), core::InvalidParameterException); ScatterPlotType scatterVary = models::getScatterPlotVaryParameters(tcplq, paramIdVec, meanVec, varianceVec, normalizationDose, 50); CHECK_EQUAL(50, scatterVary.size()); scatterCompare.setCompareScatter(scatterVary); CHECK_TESTER(scatterCompare); models::BioModelParamType variance = 0.00015; CHECK_NO_THROW(models::getScatterPlotVary1Parameter(tcplq, 0, alpha, variance, normalizationDose, 100)); tcpScatter = models::getScatterPlotVary1Parameter(tcplq, 0, alpha, variance, normalizationDose, 100); scatterCompare.setVariance(variance); scatterCompare.setCompareScatter(tcpScatter); //allow 5% of the points to deviate more scatterCompare.setAllowExceptions(true); CHECK_TESTER(scatterCompare); //test also with other parameter tcpScatter = models::getScatterPlotVary1Parameter(tcplq, 3, roh, variance, normalizationDose, 100); scatterCompare.setCompareScatter(tcpScatter); //allow 5% of the points to deviate more scatterCompare.setAllowExceptions(true); CHECK_TESTER(scatterCompare); paramIdVec.clear(); meanVec.clear(); varianceVec.clear(); paramIdVec.push_back(0); meanVec.push_back(tcplq.getAlphaMean()); varianceVec.push_back(variance); //paramIdVec.push_back(1); meanVec.push_back(tcplq.getAlphaVariance()); varianceVec.push_back(variance); paramIdVec.push_back(2); - meanVec.push_back(tcplq.getAlpahBeta()); + meanVec.push_back(tcplq.getAlphaBeta()); varianceVec.push_back(variance); paramIdVec.push_back(3); meanVec.push_back(tcplq.getRho()); varianceVec.push_back(variance); scatterVary = models::getScatterPlotVaryParameters(tcplq, paramIdVec, meanVec, varianceVec, normalizationDose, 50); scatterCompare.setCompareScatter(scatterVary); //allow 5% of the points to deviate more scatterCompare.setAllowExceptions(true); CHECK_TESTER(scatterCompare); models::BioModelParamType alphaBeta = 10; tcplq.setParameters(alpha, alphaBeta, roh, 0.08); tcplq.init(); normalizationDose = 40; curve = models::getCurveDoseVSBioModel(tcplq, normalizationDose); CHECK_NO_THROW(models::getScatterPlotVary1Parameter(tcplq, 0, alpha, 0, normalizationDose, 100)); tcpScatter = models::getScatterPlotVary1Parameter(tcplq, 0, alpha, 0, normalizationDose, 100); scatterCompare.setReferenceCurve(curve); scatterCompare.setVariance(0); //do not allow larger deviations scatterCompare.setAllowExceptions(false); scatterCompare.setCompareScatter(tcpScatter); CHECK_TESTER(scatterCompare); variance = 0.25; CHECK_NO_THROW(models::getScatterPlotVary1Parameter(tcplq, 0, alpha, variance, normalizationDose, 100)); tcpScatter = models::getScatterPlotVary1Parameter(tcplq, 0, alpha, variance, normalizationDose, 100); scatterCompare.setCompareScatter(tcpScatter); scatterCompare.setVariance(variance); //allow 5% of the points to deviate more scatterCompare.setAllowExceptions(true); CHECK_TESTER(scatterCompare); /*TCP LQ Test*/ alpha = 0.3; beta = 0.03; roh = 10000000; numFractions = 20; rttb::models::TCPLQModel tcplq_test = rttb::models::TCPLQModel(dvh_test_tv, alpha, beta, roh, numFractions); CHECK_NO_THROW(tcplq_test.init()); normalizationDose = 60; curve = models::getCurveDoseVSBioModel(tcplq_test, normalizationDose); CHECK_NO_THROW(models::getScatterPlotVary1Parameter(tcplq_test, 0, alpha, 0, normalizationDose, 100)); tcpScatter = models::getScatterPlotVary1Parameter(tcplq_test, 0, alpha, 0, normalizationDose, 100); scatterCompare.setReferenceCurve(curve); scatterCompare.setVariance(0); scatterCompare.setCompareScatter(tcpScatter); CHECK_TESTER(scatterCompare); //test also with other parameter tcpScatter = models::getScatterPlotVary1Parameter(tcplq_test, 3, roh, 0, normalizationDose, 100); scatterCompare.setCompareScatter(tcpScatter); CHECK_TESTER(scatterCompare); paramIdVec.clear(); meanVec.clear(); varianceVec.clear(); paramIdVec.push_back(0); meanVec.push_back(tcplq_test.getAlphaMean()); varianceVec.push_back(0); //paramIdVec.push_back(1); meanVec.push_back(tcplq_test.getAlphaVariance()); varianceVec.push_back(0); paramIdVec.push_back(2); - meanVec.push_back(tcplq_test.getAlpahBeta()); + meanVec.push_back(tcplq_test.getAlphaBeta()); varianceVec.push_back(0); paramIdVec.push_back(3); meanVec.push_back(tcplq_test.getRho()); varianceVec.push_back(0); scatterVary = models::getScatterPlotVaryParameters(tcplq_test, paramIdVec, meanVec, varianceVec, normalizationDose, 50); scatterCompare.setCompareScatter(scatterVary); CHECK_TESTER(scatterCompare); variance = 0.00025; CHECK_NO_THROW(models::getScatterPlotVary1Parameter(tcplq_test, 0, alpha, variance, normalizationDose, 100)); tcpScatter = models::getScatterPlotVary1Parameter(tcplq_test, 0, alpha, variance, normalizationDose, 100); scatterCompare.setCompareScatter(tcpScatter); scatterCompare.setVariance(variance); CHECK_TESTER(scatterCompare); //test also with other parameter tcpScatter = models::getScatterPlotVary1Parameter(tcplq_test, 3, roh, variance, normalizationDose, 100); scatterCompare.setCompareScatter(tcpScatter); scatterCompare.setAllowExceptions(true); CHECK_TESTER(scatterCompare); scatterCompare.setAllowExceptions(false); paramIdVec.clear(); meanVec.clear(); varianceVec.clear(); paramIdVec.push_back(0); meanVec.push_back(tcplq_test.getAlphaMean()); varianceVec.push_back(variance); //paramIdVec.push_back(1); meanVec.push_back(tcplq_test.getAlphaVariance()); varianceVec.push_back(variance); paramIdVec.push_back(2); - meanVec.push_back(tcplq_test.getAlpahBeta()); + meanVec.push_back(tcplq_test.getAlphaBeta()); varianceVec.push_back(variance); paramIdVec.push_back(3); meanVec.push_back(tcplq_test.getRho()); varianceVec.push_back(variance); scatterVary = models::getScatterPlotVaryParameters(tcplq_test, paramIdVec, meanVec, varianceVec, normalizationDose, 50); scatterCompare.setCompareScatter(scatterVary); //allow 5% of the points to deviate more scatterCompare.setAllowExceptions(true); CHECK_TESTER(scatterCompare); //DVH HT 1 rttb::io::other::DVHTxtFileReader dvhReader2 = rttb::io::other::DVHTxtFileReader(DVH_FILENAME_NT1); DVHPointer dvhPtr2 = dvhReader2.generateDVH(); CHECK_CLOSE(1.07920836034015810000e+001, models::getEUD(dvhPtr2, 10), toleranceEUD); //test RTNTCPLKBModel rttb::models::NTCPLKBModel lkb = rttb::models::NTCPLKBModel(); models::BioModelParamType aVal = 10; models::BioModelParamType mVal = 0.16; models::BioModelParamType d50Val = 55; lkb.setDVH(dvhPtr2); lkb.setA(aVal); lkb.setM(mVal); lkb.setD50(d50Val); CHECK_NO_THROW(lkb.init()); normalizationDose = 60; curve = models::getCurveDoseVSBioModel(lkb, normalizationDose); CHECK_NO_THROW(models::getScatterPlotVary1Parameter(lkb, 2, aVal, 0, normalizationDose, 100)); ScatterPlotType scatter = models::getScatterPlotVary1Parameter(lkb, 2, aVal, 0, normalizationDose, 100); scatterCompare.setReferenceCurve(curve); scatterCompare.setVariance(0); scatterCompare.setCompareScatter(scatter); CHECK_TESTER(scatterCompare); //"d50":0,"m":1,"a":2 //test also with other parameter scatter = models::getScatterPlotVary1Parameter(lkb, 0, d50Val, 0, normalizationDose, 100); scatterCompare.setCompareScatter(scatter); CHECK_TESTER(scatterCompare); paramIdVec.clear(); meanVec.clear(); varianceVec.clear(); paramIdVec.push_back(0); meanVec.push_back(lkb.getD50()); varianceVec.push_back(0); paramIdVec.push_back(1); meanVec.push_back(lkb.getM()); varianceVec.push_back(0); paramIdVec.push_back(2); meanVec.push_back(lkb.getA()); varianceVec.push_back(0); scatterVary = models::getScatterPlotVaryParameters(lkb, paramIdVec, meanVec, varianceVec, normalizationDose, 50); scatterCompare.setCompareScatter(scatterVary); CHECK_TESTER(scatterCompare); variance = 0.00025; CHECK_NO_THROW(models::getScatterPlotVary1Parameter(lkb, 2, aVal, variance, normalizationDose, 100)); scatter = models::getScatterPlotVary1Parameter(lkb, 2, aVal, variance, normalizationDose, 100); scatterCompare.setCompareScatter(scatter); scatterCompare.setVariance(variance); CHECK_TESTER(scatterCompare); //test also with other parameter scatter = models::getScatterPlotVary1Parameter(lkb, 0, d50Val, variance, normalizationDose, 100); scatterCompare.setCompareScatter(scatter); CHECK_TESTER(scatterCompare); paramIdVec.clear(); meanVec.clear(); varianceVec.clear(); paramIdVec.push_back(0); meanVec.push_back(lkb.getD50()); varianceVec.push_back(variance); paramIdVec.push_back(1); meanVec.push_back(lkb.getM()); varianceVec.push_back(variance); paramIdVec.push_back(2); meanVec.push_back(lkb.getA()); varianceVec.push_back(variance); scatterVary = models::getScatterPlotVaryParameters(lkb, paramIdVec, meanVec, varianceVec, normalizationDose, 50); scatterCompare.setCompareScatter(scatterVary); CHECK_TESTER(scatterCompare); //test RTNTCPRSModel rttb::models::NTCPRSModel rs = rttb::models::NTCPRSModel(); models::BioModelParamType gammaVal = 1.7; models::BioModelParamType sVal = 1; rs.setDVH(dvhPtr2); rs.setD50(d50Val); rs.setGamma(gammaVal); rs.setS(sVal); CHECK_NO_THROW(rs.init()); normalizationDose = 60; curve = models::getCurveDoseVSBioModel(rs, normalizationDose); CHECK_NO_THROW(models::getScatterPlotVary1Parameter(rs, 0, d50Val, 0, normalizationDose, 100)); scatter = models::getScatterPlotVary1Parameter(rs, 0, d50Val, 0, normalizationDose, 100); scatterCompare.setReferenceCurve(curve); scatterCompare.setVariance(0); scatterCompare.setCompareScatter(scatter); CHECK_TESTER(scatterCompare); //"d50":0,"gamma":1,"s":2 //test also with other parameter scatter = models::getScatterPlotVary1Parameter(rs, 1, gammaVal, 0, normalizationDose, 100); scatterCompare.setCompareScatter(scatter); CHECK_TESTER(scatterCompare); paramIdVec.clear(); meanVec.clear(); varianceVec.clear(); paramIdVec.push_back(0); meanVec.push_back(rs.getD50()); varianceVec.push_back(0); paramIdVec.push_back(1); meanVec.push_back(rs.getGamma()); varianceVec.push_back(0); paramIdVec.push_back(2); meanVec.push_back(rs.getS()); varianceVec.push_back(0); scatterVary = models::getScatterPlotVaryParameters(rs, paramIdVec, meanVec, varianceVec, normalizationDose, 50); scatterCompare.setCompareScatter(scatterVary); CHECK_TESTER(scatterCompare); variance = 0.0075; CHECK_NO_THROW(models::getScatterPlotVary1Parameter(rs, 0, d50Val, variance, normalizationDose, 100)); scatter = models::getScatterPlotVary1Parameter(rs, 0, d50Val, variance, normalizationDose, 100); scatterCompare.setCompareScatter(scatter); scatterCompare.setVariance(variance); CHECK_TESTER(scatterCompare); //test also with other parameter scatter = models::getScatterPlotVary1Parameter(rs, 2, sVal, variance, normalizationDose, 100); scatterCompare.setCompareScatter(scatter); CHECK_TESTER(scatterCompare); paramIdVec.clear(); meanVec.clear(); varianceVec.clear(); paramIdVec.push_back(0); meanVec.push_back(rs.getD50()); varianceVec.push_back(variance); paramIdVec.push_back(1); meanVec.push_back(rs.getGamma()); varianceVec.push_back(variance); paramIdVec.push_back(2); meanVec.push_back(rs.getS()); varianceVec.push_back(variance); scatterVary = models::getScatterPlotVaryParameters(rs, paramIdVec, meanVec, varianceVec, normalizationDose, 50); scatterCompare.setCompareScatter(scatterVary); CHECK_TESTER(scatterCompare); RETURN_AND_REPORT_TEST_SUCCESS; } }//testing }//rttb \ No newline at end of file diff --git a/testing/models/BioModelTest.cpp b/testing/models/BioModelTest.cpp index f79ef4f..c74d09a 100644 --- a/testing/models/BioModelTest.cpp +++ b/testing/models/BioModelTest.cpp @@ -1,288 +1,299 @@ // ----------------------------------------------------------------------- // RTToolbox - DKFZ radiotherapy quantitative evaluation library // // Copyright (c) German Cancer Research Center (DKFZ), // Software development for Integrated Diagnostics and Therapy (SIDT). // ALL RIGHTS RESERVED. // See rttbCopyright.txt or // http://www.dkfz.de/en/sidt/projects/rttb/copyright.html // // This software is distributed WITHOUT ANY WARRANTY; without even // the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR // PURPOSE. See the above copyright notices for more information. // //------------------------------------------------------------------------ /*! // @file // @version $Revision$ (last changed revision) // @date $Date$ (last change date) // @author $Author$ (last changed by) */ #include #include #include "litCheckMacros.h" #include "rttbBaseType.h" #include "rttbBaseTypeModels.h" #include "rttbBioModel.h" #include "rttbDVH.h" #include "rttbTCPLQModel.h" #include "rttbNTCPLKBModel.h" #include "rttbNTCPRSModel.h" #include "rttbBaseTypeModels.h" #include "rttbBioModelCurve.h" #include "rttbInvalidParameterException.h" #include "rttbBioModelScatterPlots.h" -namespace rttb{ - namespace testing{ +namespace rttb +{ + namespace testing + { typedef core::DVH::DataDifferentialType DataDifferentialType; /*! @brief RTBioModelTest. TCP calculated using a DVH PTV and LQ Model. NTCP tested using 3 Normal Tissue DVHs and LKB/RS Model. TCPLQ: 1) test constructors (values as expected?) 2) test init (calcTCPxxx) 3) test set/get NTCP (LKB): 1) test constructors (values as expected?) 2) test init (calcxxx) 3) test set/get NTCP (RS): 1) test constructors (values as expected?) 2) test init (calcxxx) 3) test set/get */ - int BioModelTest(int argc, char* argv[] ) + int BioModelTest(int argc, char* argv[]) { PREPARE_DEFAULT_TEST_REPORTING; //generate artificial DVH and corresponding statistical values DoseTypeGy binSize = DoseTypeGy(0.1); DoseVoxelVolumeType voxelVolume = 8; DataDifferentialType aDataDifferential; DoseCalcType value = 0; DVHVoxelNumber numberOfVoxels = 0; // creat default values - for(int i = 0; i < 98; i++){ + for (int i = 0; i < 98; i++) + { value = 0; - numberOfVoxels+=value; - aDataDifferential.push_back( value ); + numberOfVoxels += value; + aDataDifferential.push_back(value); } - aDataDifferential.push_back( 10 ); - aDataDifferential.push_back( 20 ); + + aDataDifferential.push_back(10); + aDataDifferential.push_back(20); const IDType structureID = "myStructure"; const IDType doseID = "myDose"; const IDType voxelizationID = "myVoxelization"; - core::DVH::DVHPointer dvhPtr = boost::make_shared(aDataDifferential, binSize, voxelVolume, structureID, doseID, voxelizationID); + core::DVH::DVHPointer dvhPtr = boost::make_shared(aDataDifferential, binSize, voxelVolume, structureID, + doseID, voxelizationID); //test TCP LQ Model models::BioModelParamType alpha = 0.35; models::BioModelParamType beta = 0.023333333333333; models::BioModelParamType roh = 10000000; int numFractions = 8; DoseTypeGy normalizationDose = 50; //1) test constructors (values as expected?) - rttb::models::TCPLQModel tcplq=rttb::models::TCPLQModel(); - CHECK_EQUAL(0,tcplq.getAlphaMean()); - CHECK_EQUAL(0,tcplq.getAlpahBeta()); - CHECK_EQUAL(0,tcplq.getRho()); - CHECK_EQUAL(0,tcplq.getValue()); - - tcplq=rttb::models::TCPLQModel(dvhPtr,roh, numFractions,alpha/beta, alpha,0.08); - CHECK_EQUAL(alpha,tcplq.getAlphaMean()); - CHECK_EQUAL(alpha/beta,tcplq.getAlpahBeta()); - CHECK_EQUAL(roh,tcplq.getRho()); - CHECK_EQUAL(0,tcplq.getValue()); - - tcplq=rttb::models::TCPLQModel(); - CHECK_EQUAL(0,tcplq.getAlphaMean()); - CHECK_EQUAL(0,tcplq.getAlpahBeta()); - CHECK_EQUAL(0,tcplq.getRho()); - CHECK_EQUAL(0,tcplq.getValue()); - - tcplq=rttb::models::TCPLQModel(dvhPtr,alpha, beta, roh, numFractions); - CHECK_EQUAL(alpha,tcplq.getAlphaMean()); - CHECK_EQUAL(alpha/beta,tcplq.getAlpahBeta()); - CHECK_EQUAL(roh,tcplq.getRho()); - CHECK_EQUAL(0,tcplq.getValue()); + rttb::models::TCPLQModel tcplq = rttb::models::TCPLQModel(); + CHECK_EQUAL(0, tcplq.getAlphaMean()); + CHECK_EQUAL(0, tcplq.getAlphaBeta()); + CHECK_EQUAL(0, tcplq.getRho()); + CHECK_EQUAL(0, tcplq.getValue()); + + tcplq = rttb::models::TCPLQModel(dvhPtr, roh, numFractions, alpha / beta, alpha, 0.08); + CHECK_EQUAL(alpha, tcplq.getAlphaMean()); + CHECK_EQUAL(alpha / beta, tcplq.getAlphaBeta()); + CHECK_EQUAL(roh, tcplq.getRho()); + CHECK_EQUAL(0, tcplq.getValue()); + + tcplq = rttb::models::TCPLQModel(); + CHECK_EQUAL(0, tcplq.getAlphaMean()); + CHECK_EQUAL(0, tcplq.getAlphaBeta()); + CHECK_EQUAL(0, tcplq.getRho()); + CHECK_EQUAL(0, tcplq.getValue()); + + tcplq = rttb::models::TCPLQModel(dvhPtr, alpha, beta, roh, numFractions); + CHECK_EQUAL(alpha, tcplq.getAlphaMean()); + CHECK_EQUAL(alpha / beta, tcplq.getAlphaBeta()); + CHECK_EQUAL(roh, tcplq.getRho()); + CHECK_EQUAL(0, tcplq.getValue()); //2) test init (calcTCPxxx) CHECK_NO_THROW(tcplq.init(1)); //3) test set/get - CHECK_EQUAL(0,tcplq.getValue()); - CHECK_NO_THROW(tcplq.setParameters(alpha,10,roh,0.08)); - CHECK_EQUAL(10,tcplq.getAlpahBeta()); - CHECK_EQUAL(0.08,tcplq.getAlphaVariance()); - CHECK_EQUAL(alpha,tcplq.getAlphaMean()); - CHECK_EQUAL(roh,tcplq.getRho()); + CHECK_EQUAL(0, tcplq.getValue()); + CHECK_NO_THROW(tcplq.setParameters(alpha, 10, roh, 0.08)); + CHECK_EQUAL(10, tcplq.getAlphaBeta()); + CHECK_EQUAL(0.08, tcplq.getAlphaVariance()); + CHECK_EQUAL(alpha, tcplq.getAlphaMean()); + CHECK_EQUAL(roh, tcplq.getRho()); - CHECK_NO_THROW(models::getCurveDoseVSBioModel(tcplq,normalizationDose)); + CHECK_NO_THROW(models::getCurveDoseVSBioModel(tcplq, normalizationDose)); std::vector aParameterVector; - aParameterVector.push_back(alpha+0.02); - CHECK_THROW_EXPLICIT(tcplq.setParameterVector(aParameterVector),core::InvalidParameterException); + aParameterVector.push_back(alpha + 0.02); + CHECK_THROW_EXPLICIT(tcplq.setParameterVector(aParameterVector), core::InvalidParameterException); aParameterVector.push_back(0.06); aParameterVector.push_back(8); - aParameterVector.push_back(roh/10); + aParameterVector.push_back(roh / 10); CHECK_NO_THROW(tcplq.setParameterVector(aParameterVector)); - CHECK_EQUAL(8,tcplq.getAlpahBeta()); - CHECK_EQUAL(0.06,tcplq.getAlphaVariance()); - CHECK_EQUAL(alpha+0.02,tcplq.getAlphaMean()); - CHECK_EQUAL(roh/10,tcplq.getRho()); - - for (int i = 0; i < 4; i++){ - CHECK_NO_THROW(tcplq.setParameterByID(i,models::BioModelParamType(i))); + CHECK_EQUAL(8, tcplq.getAlphaBeta()); + CHECK_EQUAL(0.06, tcplq.getAlphaVariance()); + CHECK_EQUAL(alpha + 0.02, tcplq.getAlphaMean()); + CHECK_EQUAL(roh / 10, tcplq.getRho()); + + for (int i = 0; i < 4; i++) + { + CHECK_NO_THROW(tcplq.setParameterByID(i, models::BioModelParamType(i))); } - CHECK_THROW_EXPLICIT(tcplq.setParameterByID(4,4.0),core::InvalidParameterException); + CHECK_THROW_EXPLICIT(tcplq.setParameterByID(4, 4.0), core::InvalidParameterException); - CHECK_EQUAL(0,tcplq.getParameterID("alphaMean")); - CHECK_EQUAL(0,tcplq.getAlphaMean()); - CHECK_EQUAL(1,tcplq.getParameterID("alphaVariance")); - CHECK_EQUAL(1,tcplq.getAlphaVariance()); - CHECK_EQUAL(2,tcplq.getParameterID("alpha_beta")); - CHECK_EQUAL(2,tcplq.getAlpahBeta()); - CHECK_EQUAL(3,tcplq.getParameterID("rho")); - CHECK_EQUAL(3,tcplq.getRho()); + CHECK_EQUAL(0, tcplq.getParameterID("alphaMean")); + CHECK_EQUAL(0, tcplq.getAlphaMean()); + CHECK_EQUAL(1, tcplq.getParameterID("alphaVariance")); + CHECK_EQUAL(1, tcplq.getAlphaVariance()); + CHECK_EQUAL(2, tcplq.getParameterID("alpha_beta")); + CHECK_EQUAL(2, tcplq.getAlphaBeta()); + CHECK_EQUAL(3, tcplq.getParameterID("rho")); + CHECK_EQUAL(3, tcplq.getRho()); //test NTCPLKBModel //1) test constructors (values as expected?) models::BioModelParamType aVal = 10; models::BioModelParamType mVal = 0.16; models::BioModelParamType d50Val = 35; CHECK_NO_THROW(rttb::models::NTCPLKBModel()); - rttb::models::NTCPLKBModel lkb=rttb::models::NTCPLKBModel(); - CHECK_EQUAL(0,lkb.getA()); - CHECK_EQUAL(0,lkb.getM()); - CHECK_EQUAL(0,lkb.getD50()); - CHECK_EQUAL(0,lkb.getValue()); + rttb::models::NTCPLKBModel lkb = rttb::models::NTCPLKBModel(); + CHECK_EQUAL(0, lkb.getA()); + CHECK_EQUAL(0, lkb.getM()); + CHECK_EQUAL(0, lkb.getD50()); + CHECK_EQUAL(0, lkb.getValue()); CHECK_NO_THROW(rttb::models::NTCPLKBModel(dvhPtr, d50Val, mVal, aVal)); - lkb=rttb::models::NTCPLKBModel(dvhPtr, d50Val, mVal, aVal); - CHECK_EQUAL(0,lkb.getValue()); - CHECK_EQUAL(dvhPtr,lkb.getDVH()); - CHECK_EQUAL(aVal,lkb.getA()); - CHECK_EQUAL(mVal,lkb.getM()); - CHECK_EQUAL(d50Val,lkb.getD50()); + lkb = rttb::models::NTCPLKBModel(dvhPtr, d50Val, mVal, aVal); + CHECK_EQUAL(0, lkb.getValue()); + CHECK_EQUAL(dvhPtr, lkb.getDVH()); + CHECK_EQUAL(aVal, lkb.getA()); + CHECK_EQUAL(mVal, lkb.getM()); + CHECK_EQUAL(d50Val, lkb.getD50()); //2) test init (calcxxx) CHECK_NO_THROW(lkb.init(1)); lkb.getValue(); //3) test set/get - lkb=rttb::models::NTCPLKBModel(); - CHECK_EQUAL(0,lkb.getA()); - CHECK_EQUAL(0,lkb.getM()); - CHECK_EQUAL(0,lkb.getD50()); + lkb = rttb::models::NTCPLKBModel(); + CHECK_EQUAL(0, lkb.getA()); + CHECK_EQUAL(0, lkb.getM()); + CHECK_EQUAL(0, lkb.getD50()); lkb.setDVH(dvhPtr); - CHECK_EQUAL(dvhPtr,lkb.getDVH()); + CHECK_EQUAL(dvhPtr, lkb.getDVH()); lkb.setA(aVal); - CHECK_EQUAL(aVal,lkb.getA()); + CHECK_EQUAL(aVal, lkb.getA()); lkb.setM(mVal); - CHECK_EQUAL(mVal,lkb.getM()); + CHECK_EQUAL(mVal, lkb.getM()); lkb.setD50(d50Val); - CHECK_EQUAL(d50Val,lkb.getD50()); + CHECK_EQUAL(d50Val, lkb.getD50()); CHECK_NO_THROW(models::getCurveEUDVSBioModel(lkb)); aParameterVector.clear(); - aParameterVector.push_back(d50Val+5); - CHECK_THROW_EXPLICIT(lkb.setParameterVector(aParameterVector),core::InvalidParameterException); - aParameterVector.push_back(mVal+0.2); - aParameterVector.push_back(aVal+0.5); + aParameterVector.push_back(d50Val + 5); + CHECK_THROW_EXPLICIT(lkb.setParameterVector(aParameterVector), core::InvalidParameterException); + aParameterVector.push_back(mVal + 0.2); + aParameterVector.push_back(aVal + 0.5); CHECK_NO_THROW(lkb.setParameterVector(aParameterVector)); - CHECK_EQUAL(aVal+0.5,lkb.getA()); - CHECK_EQUAL(mVal+0.2,lkb.getM()); - CHECK_EQUAL(d50Val+5,lkb.getD50()); + CHECK_EQUAL(aVal + 0.5, lkb.getA()); + CHECK_EQUAL(mVal + 0.2, lkb.getM()); + CHECK_EQUAL(d50Val + 5, lkb.getD50()); - for (int i = 0; i < 3; i++){ - CHECK_NO_THROW(lkb.setParameterByID(i,models::BioModelParamType(i))); + for (int i = 0; i < 3; i++) + { + CHECK_NO_THROW(lkb.setParameterByID(i, models::BioModelParamType(i))); } - CHECK_THROW_EXPLICIT(lkb.setParameterByID(3,4.0),core::InvalidParameterException); + CHECK_THROW_EXPLICIT(lkb.setParameterByID(3, 4.0), core::InvalidParameterException); - CHECK_EQUAL(0,lkb.getParameterID("d50")); - CHECK_EQUAL(0,lkb.getD50()); - CHECK_EQUAL(1,lkb.getParameterID("m")); - CHECK_EQUAL(1,lkb.getM()); - CHECK_EQUAL(2,lkb.getParameterID("a")); - CHECK_EQUAL(2,lkb.getA()); + CHECK_EQUAL(0, lkb.getParameterID("d50")); + CHECK_EQUAL(0, lkb.getD50()); + CHECK_EQUAL(1, lkb.getParameterID("m")); + CHECK_EQUAL(1, lkb.getM()); + CHECK_EQUAL(2, lkb.getParameterID("a")); + CHECK_EQUAL(2, lkb.getA()); //test NTCPRSModel //1) test constructors (values as expected?) CHECK_NO_THROW(rttb::models::NTCPRSModel()); models::BioModelParamType gammaVal = 1.7; models::BioModelParamType sVal = 1; - CHECK_NO_THROW(rttb::models::NTCPRSModel(dvhPtr,d50Val,gammaVal,sVal)); - rttb::models::NTCPRSModel rs=rttb::models::NTCPRSModel(dvhPtr,d50Val,gammaVal,sVal); - CHECK_EQUAL(dvhPtr,rs.getDVH()); - CHECK_EQUAL(d50Val,rs.getD50()); - CHECK_EQUAL(gammaVal,rs.getGamma()); - CHECK_EQUAL(sVal,rs.getS()); - - rs=rttb::models::NTCPRSModel(); - CHECK_EQUAL(0,rs.getGamma()); - CHECK_EQUAL(0,rs.getS()); - CHECK_EQUAL(0,rs.getD50()); + CHECK_NO_THROW(rttb::models::NTCPRSModel(dvhPtr, d50Val, gammaVal, sVal)); + rttb::models::NTCPRSModel rs = rttb::models::NTCPRSModel(dvhPtr, d50Val, gammaVal, sVal); + CHECK_EQUAL(dvhPtr, rs.getDVH()); + CHECK_EQUAL(d50Val, rs.getD50()); + CHECK_EQUAL(gammaVal, rs.getGamma()); + CHECK_EQUAL(sVal, rs.getS()); + + rs = rttb::models::NTCPRSModel(); + CHECK_EQUAL(0, rs.getGamma()); + CHECK_EQUAL(0, rs.getS()); + CHECK_EQUAL(0, rs.getD50()); //3) test set/get rs.setDVH(dvhPtr); - CHECK_EQUAL(dvhPtr,rs.getDVH()); + CHECK_EQUAL(dvhPtr, rs.getDVH()); rs.setD50(d50Val); - CHECK_EQUAL(d50Val,rs.getD50()); + CHECK_EQUAL(d50Val, rs.getD50()); rs.setGamma(gammaVal); - CHECK_EQUAL(gammaVal,rs.getGamma()); + CHECK_EQUAL(gammaVal, rs.getGamma()); rs.setS(sVal); - CHECK_EQUAL(sVal,rs.getS()); + CHECK_EQUAL(sVal, rs.getS()); //2) test init (calcxxx) CHECK_NO_THROW(rs.init(1)); //3) test set/get continued aParameterVector.clear(); - aParameterVector.push_back(d50Val+5); - CHECK_THROW_EXPLICIT(rs.setParameterVector(aParameterVector),core::InvalidParameterException); - aParameterVector.push_back(gammaVal+0.2); - aParameterVector.push_back(sVal+0.5); + aParameterVector.push_back(d50Val + 5); + CHECK_THROW_EXPLICIT(rs.setParameterVector(aParameterVector), core::InvalidParameterException); + aParameterVector.push_back(gammaVal + 0.2); + aParameterVector.push_back(sVal + 0.5); CHECK_NO_THROW(rs.setParameterVector(aParameterVector)); - CHECK_EQUAL(gammaVal+0.2,rs.getGamma()); - CHECK_EQUAL(sVal+0.5,rs.getS()); - CHECK_EQUAL(d50Val+5,rs.getD50()); + CHECK_EQUAL(gammaVal + 0.2, rs.getGamma()); + CHECK_EQUAL(sVal + 0.5, rs.getS()); + CHECK_EQUAL(d50Val + 5, rs.getD50()); - for (int i = 0; i < 3; i++){ - CHECK_NO_THROW(rs.setParameterByID(i,models::BioModelParamType(i))); + for (int i = 0; i < 3; i++) + { + CHECK_NO_THROW(rs.setParameterByID(i, models::BioModelParamType(i))); } - CHECK_THROW_EXPLICIT(rs.setParameterByID(3,4.0),core::InvalidParameterException); + CHECK_THROW_EXPLICIT(rs.setParameterByID(3, 4.0), core::InvalidParameterException); - CHECK_EQUAL(0,rs.getParameterID("d50")); - CHECK_EQUAL(0,rs.getD50()); - CHECK_EQUAL(1,rs.getParameterID("gamma")); - CHECK_EQUAL(1,rs.getGamma()); - CHECK_EQUAL(2,rs.getParameterID("s")); - CHECK_EQUAL(2,rs.getS()); + CHECK_EQUAL(0, rs.getParameterID("d50")); + CHECK_EQUAL(0, rs.getD50()); + CHECK_EQUAL(1, rs.getParameterID("gamma")); + CHECK_EQUAL(1, rs.getGamma()); + CHECK_EQUAL(2, rs.getParameterID("s")); + CHECK_EQUAL(2, rs.getS()); //Scatter plot tests - CHECK_NO_THROW(models::getScatterPlotVary1Parameter(tcplq,0,alpha,0,normalizationDose));//variance=0, will be set to 1e-30 - CHECK_THROW_EXPLICIT(models::getScatterPlotVary1Parameter(tcplq,0,alpha,alpha*0.1,0),core::InvalidParameterException);//normalisationdose=0 - CHECK_THROW_EXPLICIT(models::getScatterPlotVary1Parameter(tcplq,0,alpha,alpha*0.1,normalizationDose,10000,0,0),core::InvalidParameterException);//maxDose-minDose=0 + CHECK_NO_THROW(models::getScatterPlotVary1Parameter(tcplq, 0, alpha, 0, + normalizationDose)); //variance=0, will be set to 1e-30 + CHECK_THROW_EXPLICIT(models::getScatterPlotVary1Parameter(tcplq, 0, alpha, alpha * 0.1, 0), + core::InvalidParameterException);//normalisationdose=0 + CHECK_THROW_EXPLICIT(models::getScatterPlotVary1Parameter(tcplq, 0, alpha, alpha * 0.1, normalizationDose, 10000, 0, 0), + core::InvalidParameterException);//maxDose-minDose=0 RETURN_AND_REPORT_TEST_SUCCESS; } }//testing }//rttb \ No newline at end of file diff --git a/testing/models/CMakeLists.txt b/testing/models/CMakeLists.txt index 504fd5f..c81cb75 100644 --- a/testing/models/CMakeLists.txt +++ b/testing/models/CMakeLists.txt @@ -1,21 +1,22 @@ #----------------------------------------------------------------------------- # Setup the system information test. Write out some basic failsafe # information in case the test doesn't run. #----------------------------------------------------------------------------- SET(MODELS_TESTS ${EXECUTABLE_OUTPUT_PATH}/rttbModelTests) SET(MODELS_HEADER_TEST ${EXECUTABLE_OUTPUT_PATH}/rttbModelsHeaderTest) SET(TEST_DATA_ROOT ${RTTBTesting_SOURCE_DIR}/data) SET(TEMP ${RTTBTesting_BINARY_DIR}/temporary) #----------------------------------------------------------------------------- ADD_TEST(BioModelTest ${MODELS_TESTS} BioModelTest) ADD_TEST(BioModelScatterPlotTest ${MODELS_TESTS} BioModelScatterPlotTest) ADD_TEST(DvhBasedModelsTest ${MODELS_TESTS} DvhBasedModelsTest) +ADD_TEST(LQModelAccessorTest ${MODELS_TESTS} LQModelAccessorTest "${TEST_DATA_ROOT}/DICOM/TestDose/ConstantTwo.dcm" "${TEST_DATA_ROOT}/DICOM/TestDose/LinearIncreaseX.dcm") -RTTB_CREATE_TEST_MODULE(rttbModel DEPENDS RTTBModels RTTBOtherIO PACKAGE_DEPENDS Boost Litmus DCMTK) +RTTB_CREATE_TEST_MODULE(rttbModel DEPENDS RTTBModels RTTBOtherIO RTTBDicomIO PACKAGE_DEPENDS Boost Litmus DCMTK) diff --git a/testing/models/LQModelAccessorTest.cpp b/testing/models/LQModelAccessorTest.cpp new file mode 100644 index 0000000..b7c3c3c --- /dev/null +++ b/testing/models/LQModelAccessorTest.cpp @@ -0,0 +1,114 @@ +// ----------------------------------------------------------------------- +// 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: 4047 $ (last changed revision) +// @date $Date: 2012-10-29 16:19:15 +0100 (Mo, 29 Okt 2012) $ (last change date) +// @author $Author: mang $ (last changed by) +*/ +#include "litCheckMacros.h" + +#include + +#include "rttbLQModelAccessor.h" +#include "rttbGenericDoseIterator.h" +#include "rttbDicomFileDoseAccessorGenerator.h" +#include "rttbInvalidDoseException.h" + +namespace rttb +{ + namespace testing + { + + /*! @brief LQModelAccessorTest. + 1) Test constructor + 2) Test getGeometricInfo() + 3) Test getBioModelValueAt() + + */ + int LQModelAccessorTest(int argc, char* argv[]) + { + PREPARE_DEFAULT_TEST_REPORTING; + + std::string RTDOSE_FILENAME_CONSTANT_TWO; + std::string RTDOSE_FILENAME_INCREASE_X; + + if (argc > 2) + { + RTDOSE_FILENAME_CONSTANT_TWO = argv[1]; + RTDOSE_FILENAME_INCREASE_X = argv[2]; + } + else + { + + std::cout << "at least two parameters for LQModelAccessorTest are expected" << + std::endl; + return -1; + } + + typedef core::GenericDoseIterator::DoseAccessorPointer DoseAccessorPointer; + + rttb::io::dicom::DicomFileDoseAccessorGenerator doseAccessorGenerator1( + RTDOSE_FILENAME_CONSTANT_TWO.c_str()); + DoseAccessorPointer doseAccessor1(doseAccessorGenerator1.generateDoseAccessor()); + + rttb::io::dicom::DicomFileDoseAccessorGenerator doseAccessorGenerator2( + RTDOSE_FILENAME_INCREASE_X.c_str()); + DoseAccessorPointer doseAccessor2(doseAccessorGenerator2.generateDoseAccessor()); + + core::GeometricInfo doseAccessor1GeometricInfo = doseAccessor1->getGeometricInfo(); + core::GeometricInfo doseAccessor2GeometricInfo = doseAccessor2->getGeometricInfo(); + + DoseAccessorPointer doseAccessorNull; + + core::BioModelAccessorInterface::BioModelAccessorPointer LQWithConstantDose; + core::BioModelAccessorInterface::BioModelAccessorPointer LQWithIncreaseXDose; + + //1) test constructor + CHECK_THROW_EXPLICIT(models::LQModelAccessor(doseAccessorNull, 0, 0), core::InvalidDoseException); + + CHECK_NO_THROW(LQWithConstantDose = boost::make_shared(doseAccessor1, + 0.2, 0.02)); + CHECK_NO_THROW(LQWithIncreaseXDose = boost::make_shared(doseAccessor2, + 0.3, 0.01)); + + //2) Test getGeometricInfo() + CHECK_EQUAL(LQWithConstantDose->getGeometricInfo(), doseAccessor1GeometricInfo); + CHECK_EQUAL(LQWithIncreaseXDose->getGeometricInfo(), doseAccessor2GeometricInfo); + + //3) Test getBioModelValueAt() + + models::BioModelParamType expectedLQWithDoseTwo = exp(-(0.2 * 2 + (0.02 * 2 * 2))); + CHECK_EQUAL(LQWithConstantDose->getBioModelValueAt(0), expectedLQWithDoseTwo); + CHECK_EQUAL(LQWithConstantDose->getBioModelValueAt(LQWithConstantDose->getGridSize() - 1), expectedLQWithDoseTwo); + CHECK_EQUAL(LQWithConstantDose->getBioModelValueAt(VoxelGridIndex3D(1, 2, 6)), expectedLQWithDoseTwo); + CHECK_EQUAL(LQWithConstantDose->getBioModelValueAt(VoxelGridIndex3D(65, 40, 60)), expectedLQWithDoseTwo); + CHECK_EQUAL(LQWithIncreaseXDose->getBioModelValueAt(0), 1); + models::BioModelParamType expectedLQWithIncreaseX = exp(-(0.3 * 66 * 2.822386e-5 + (0.01 * 66 * 2.822386e-5 * 66 * + 2.822386e-5))); + CHECK_CLOSE(LQWithIncreaseXDose->getBioModelValueAt(LQWithIncreaseXDose->getGridSize() - 1), expectedLQWithIncreaseX, + errorConstant); + expectedLQWithIncreaseX = exp(-(0.3 * 1 * 2.822386e-5 + (0.01 * 1 * 2.822386e-5 * 1 * 2.822386e-5))); + CHECK_CLOSE(LQWithIncreaseXDose->getBioModelValueAt(VoxelGridIndex3D(1, 2, 6)), expectedLQWithIncreaseX, errorConstant); + expectedLQWithIncreaseX = exp(-(0.3 * 45 * 2.822386e-5 + (0.01 * 45 * 2.822386e-5 * 45 * 2.822386e-5))); + CHECK_CLOSE(LQWithIncreaseXDose->getBioModelValueAt(VoxelGridIndex3D(45, 40, 60)), expectedLQWithIncreaseX, + errorConstant); + + RETURN_AND_REPORT_TEST_SUCCESS; + + } + + }//testing +}//rttb \ No newline at end of file diff --git a/testing/models/files.cmake b/testing/models/files.cmake index 1e5c4fa..ca92f5f 100644 --- a/testing/models/files.cmake +++ b/testing/models/files.cmake @@ -1,15 +1,16 @@ SET(CPP_FILES rttbModelsTest.cpp BioModelTest.cpp rttbScatterTester.cpp BioModelScatterPlotTest.cpp DummyModel.cpp ../core/DummyDVHGenerator.cpp DvhBasedModelsTest.cpp + LQModelAccessorTest.cpp ) SET(H_FILES rttbScatterTester.h DummyModel.h ../core/DummyDVHGenerator.h ) diff --git a/testing/models/rttbModelsTest.cpp b/testing/models/rttbModelsTest.cpp index 59953f0..ee706c1 100644 --- a/testing/models/rttbModelsTest.cpp +++ b/testing/models/rttbModelsTest.cpp @@ -1,65 +1,66 @@ // ----------------------------------------------------------------------- // 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) */ // this file defines the rttbCoreTests for the test driver // and all it expects is that you have a function called RegisterTests #if defined(_MSC_VER) #pragma warning ( disable : 4786 ) #endif #include "litMultiTestsMain.h" namespace rttb { namespace testing { void registerTests() { LIT_REGISTER_TEST(BioModelTest); LIT_REGISTER_TEST(BioModelScatterPlotTest); LIT_REGISTER_TEST(DvhBasedModelsTest); + LIT_REGISTER_TEST(LQModelAccessorTest); } } } int main(int argc, char* argv[]) { int result = 0; rttb::testing::registerTests(); try { result = lit::multiTestsMain(argc, argv); } catch (const std::exception& /*e*/) { result = -1; } catch (...) { result = -1; } return result; }