diff --git a/Modules/PhotoacousticsLib/files.cmake b/Modules/PhotoacousticsLib/files.cmake index 3f62323168..5bc7c84dfc 100644 --- a/Modules/PhotoacousticsLib/files.cmake +++ b/Modules/PhotoacousticsLib/files.cmake @@ -1,56 +1,65 @@ SET(H_FILES include/mitkPAPropertyCalculator.h include/mitkPAVector.h include/mitkPATissueGeneratorParameters.h include/mitkPAInSilicoTissueVolume.h - include/mitkPAPhantomTissueGenerator.h include/mitkPATissueGenerator.h include/mitkPAVesselTree.h include/mitkPAVessel.h - include/mitkPAVesselDrawer.h include/mitkPAVesselMeanderStrategy.h include/mitkPANoiseGenerator.h include/mitkPAVolume.h include/mitkPAComposedVolume.h include/mitkPASlicedVolumeGenerator.h include/mitkPAProbe.h include/mitkPALightSource.h include/mitkPAIOUtil.h include/mitkPAMonteCarloThreadHandler.h include/mitkPASimulationBatchGenerator.h include/mitkPAFluenceYOffsetPair.h include/mitkPAVolumeManipulator.h include/mitkPAVesselProperties.h include/mitkPASimulationBatchGeneratorParameters.h include/mitkPAExceptions.h + include/mitkPASpectralUnmixingFilterBase.h + include/mitkPALinearSpectralUnmixingFilter.h + include/mitkPASpectralUnmixingSO2.h + include/mitkPASpectralUnmixingFilterLagrange.h + include/mitkPASpectralUnmixingFilterSimplex.h + include/mitkPASpectralUnmixingFilterVigra.h ) set(CPP_FILES Domain/Vessel/mitkPAVesselTree.cpp Domain/Vessel/mitkPAVessel.cpp Domain/Vessel/mitkPAVesselMeanderStrategy.cpp Domain/Vessel/mitkPAVesselProperties.cpp Domain/Volume/mitkPAInSilicoTissueVolume.cpp Domain/Volume/mitkPAVolume.cpp Domain/Volume/mitkPAComposedVolume.cpp Domain/Volume/mitkPAFluenceYOffsetPair.cpp Generator/mitkPATissueGenerator.cpp - Generator/mitkPAPhantomTissueGenerator.cpp Generator/mitkPANoiseGenerator.cpp Generator/mitkPASlicedVolumeGenerator.cpp Generator/mitkPASimulationBatchGenerator.cpp Generator/mitkPASimulationBatchGeneratorParameters.cpp IO/mitkPAIOUtil.cpp Utils/mitkPAPropertyCalculator.cpp Utils/mitkPAVector.cpp Utils/mitkPATissueGeneratorParameters.cpp Utils/mitkPAVolumeManipulator.cpp Utils/ProbeDesign/mitkPAProbe.cpp Utils/ProbeDesign/mitkPALightSource.cpp Utils/Thread/mitkPAMonteCarloThreadHandler.cpp + SUFilter/mitkPASpectralUnmixingFilterBase.cpp + SUFilter/mitkPALinearSpectralUnmixingFilter.cpp + SUFilter/mitkPASpectralUnmixingSO2.cpp + SUFilter/mitkPASpectralUnmixingFilterSimplex.cpp + SUFilter/mitkPASpectralUnmixingFilterVigra.cpp + SUFilter/mitkPASpectralUnmixingFilterLagrange.cpp Utils/mitkPAVesselDrawer.cpp ) set(RESOURCE_FILES spectralLIB.dat ) diff --git a/Modules/PhotoacousticsLib/include/mitkPALinearSpectralUnmixingFilter.h b/Modules/PhotoacousticsLib/include/mitkPALinearSpectralUnmixingFilter.h new file mode 100644 index 0000000000..3f710b79a2 --- /dev/null +++ b/Modules/PhotoacousticsLib/include/mitkPALinearSpectralUnmixingFilter.h @@ -0,0 +1,103 @@ +/*=================================================================== + +The Medical Imaging Interaction Toolkit (MITK) + +Copyright (c) German Cancer Research Center, +Division of Medical and Biological Informatics. +All rights reserved. + +This software is distributed WITHOUT ANY WARRANTY; without +even the implied warranty of MERCHANTABILITY or FITNESS FOR +A PARTICULAR PURPOSE. + +See LICENSE.txt or http://www.mitk.org for details. + +===================================================================*/ + +#ifndef MITKLINEARPHOTOACOUSTICSPECTRALUNMIXINGFILTER_H +#define MITKLINEARPHOTOACOUSTICSPECTRALUNMIXINGFILTER_H + +#include +#include + +namespace mitk { + namespace pa { + /** + * \brief This filter is subclass of the spectral unmixing filter base. All algorithms in this class are + * based on the Eigen open source c++ library. It takes a multispectral pixel as input and returns a vector + * with the unmixing result. + * + * Input: + * - endmemberMatrix Eigen Matrix with number of chromophores colums and number of wavelengths rows so matrix element (i,j) contains + * the absorbtion of chromophore j @ wavelength i taken from the database by PropertyElement method. + * - inputVector Eigen Vector containing values of one pixel of XY-plane image with number of wavelength rows (z-dimension of a sequenece) + * so the pixelvalue of the first wavelength is stored in inputVector[0] and so on. + * + * Output: + * - Eigen vector with unmixing result of one multispectral pixel. The ith element of the vector corresponds to the ith entry of the + * m_Chromophore vector. + * + * Algorithms (see AlgortihmType enum): + * - HOUSEHOLDERQR computes the solution by QR decomposition + * - COLPIVHOUSEHOLDERQR computes the solution by QR decomposition + * - FULLPIVHOUSEHOLDERQR computes the solution by QR decomposition + * - LDLT computes the solution by Cholesky decomposition + * - LLT computes the solution by Cholesky decomposition + * - JACOBISVD computes the solution by singular value decomposition uses least square solver + * + * Possible exceptions: + * - "algorithm not working": doesn't work now (2018/06/19) + * - "404 VIGRA ALGORITHM NOT FOUND": Algorithm not implemented + */ + + class MITKPHOTOACOUSTICSLIB_EXPORT LinearSpectralUnmixingFilter : public SpectralUnmixingFilterBase + { + public: + + mitkClassMacro(LinearSpectralUnmixingFilter, SpectralUnmixingFilterBase) + itkFactorylessNewMacro(Self) + + /** + * \brief Contains all implemented Eigen algorithms for spectral unmixing. For detailed information of the algorithms look at the + * mitkPALinearSpectralUnmixingFilter.h documentation (section Algorithms). + */ + enum AlgortihmType + { + HOUSEHOLDERQR, + LDLT, + LLT, + COLPIVHOUSEHOLDERQR, + JACOBISVD, + FULLPIVLU, + FULLPIVHOUSEHOLDERQR + }; + + /** + * \brief Takes a mitk::pa::LinearSpectralUnmixingFilter::AlgortihmType and fix it for usage at the "SpectralUnmixingAlgorithm" method. + * @param algorithmName has to be a mitk::pa::LinearSpectralUnmixingFilter::AlgortihmType + */ + void SetAlgorithm(AlgortihmType inputAlgorithmName); + + protected: + LinearSpectralUnmixingFilter(); + virtual ~LinearSpectralUnmixingFilter(); + + /** + * \brief overrides the baseclass method with a mehtod to calculate the spectral unmixing result vector. Herain the class performs the + * algorithm set by the "SetAlgorithm" method and writes the result into a Eigen vector which is the return value. + * @param endmemberMatrix Matrix with number of chromophores colums and number of wavelengths rows so matrix element (i,j) contains + * the absorbtion of chromophore j @ wavelength i taken from the database by PropertyElement method. + * @param inputVector Vector containing values of one pixel of XY-plane image with number of wavelength rows (z-dimension of a sequenece) + * so the pixelvalue of the first wavelength is stored in inputVector[0] and so on. + * @throws if the algorithmName is not a member of the enum VigraAlgortihmType + * @throws if one chooses the ldlt/llt solver which doens't work yet + */ + virtual Eigen::VectorXf SpectralUnmixingAlgorithm(Eigen::Matrix endmemberMatrix, + Eigen::VectorXf inputVector) override; + + private: + AlgortihmType algorithmName; + }; + } +} +#endif // MITKLINEARPHOTOACOUSTICSPECTRALUNMIXINGFILTER_H diff --git a/Modules/PhotoacousticsLib/include/mitkPAPropertyCalculator.h b/Modules/PhotoacousticsLib/include/mitkPAPropertyCalculator.h index 232a4ef84c..89bf0be93d 100644 --- a/Modules/PhotoacousticsLib/include/mitkPAPropertyCalculator.h +++ b/Modules/PhotoacousticsLib/include/mitkPAPropertyCalculator.h @@ -1,74 +1,79 @@ /*=================================================================== The Medical Imaging Interaction Toolkit (MITK) Copyright (c) German Cancer Research Center, Division of Medical and Biological Informatics. All rights reserved. This software is distributed WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See LICENSE.txt or http://www.mitk.org for details. ===================================================================*/ #ifndef MITKPHOTOACOUSTICPROPERTYCALCULATOR_H #define MITKPHOTOACOUSTICPROPERTYCALCULATOR_H #include "MitkPhotoacousticsLibExports.h" //Includes for smart pointer usage #include "mitkCommon.h" #include "itkLightObject.h" namespace mitk { namespace pa { class MITKPHOTOACOUSTICSLIB_EXPORT PropertyCalculator : public itk::LightObject { public: mitkClassMacroItkParent(PropertyCalculator, itk::LightObject) itkFactorylessNewMacro(Self) - struct Properties + struct Properties { double mua; double mus; double g; }; enum TissueType { AIR = 1, BLOOD = 2, EPIDERMIS = 3, FAT = 4, STANDARD_TISSUE = 5 }; + enum ChromophoreType + { + OXYGENATED = 1, + DEOXYGENATED = 2, + WATER = 3, + FATTY = 4, + MELANIN = 5, + ONEENDMEMBER = 6 + }; + + double GetAbsorptionForWavelength( + ChromophoreType chromophoreType, int wavelength); + Properties CalculatePropertyForSpecificWavelength( TissueType tissueType, int wavelength, double oxygenSaturatonInFraction = 0); protected: PropertyCalculator(); ~PropertyCalculator() override; bool m_Valid = false; std::map> m_SpectralLibMap; - - enum MapType - { - OXYGENATED = 1, - DEOXYGENATED = 2, - WATER = 3, - FATTY = 4, - MELANIN = 5 - }; + }; } } #endif // MITKPHOTOACOUSTICPROPERTYCALCULATOR_H diff --git a/Modules/PhotoacousticsLib/include/mitkPASpectralUnmixingFilterBase.h b/Modules/PhotoacousticsLib/include/mitkPASpectralUnmixingFilterBase.h new file mode 100644 index 0000000000..a297b37b1a --- /dev/null +++ b/Modules/PhotoacousticsLib/include/mitkPASpectralUnmixingFilterBase.h @@ -0,0 +1,198 @@ +/*=================================================================== + +The Medical Imaging Interaction Toolkit (MITK) + +Copyright (c) German Cancer Research Center, +Division of Medical and Biological Informatics. +All rights reserved. + +This software is distributed WITHOUT ANY WARRANTY; without +even the implied warranty of MERCHANTABILITY or FITNESS FOR +A PARTICULAR PURPOSE. + +See LICENSE.txt or http://www.mitk.org for details. + +===================================================================*/ + +#ifndef MITKPHOTOACOUSTICSPECTRALUNMIXINGFILTERBASE_H +#define MITKPHOTOACOUSTICSPECTRALUNMIXINGFILTERBASE_H + +#include "mitkImageToImageFilter.h" +#include + +//Includes for smart pointer usage +#include "mitkCommon.h" +#include "itkLightObject.h" + +// Includes for AddEnmemberMatrix +#include "mitkPAPropertyCalculator.h" +#include + +namespace mitk { + namespace pa { + /** + * \brief The spectral unmixing filter base is designed as superclass for several spectral unmixing filter eg. Eigen- or Vigrabased ones. + * One can add wavelengths and chromophores and get a unmixed output images out of one MITK input image using algorithms from the subclasses. + * + * Input: + * The unmixing input has to be a 3D MITK image where the XY-plane is a image and the Z-direction represents recordings for different + * wavelengths. Herein a XY-plane for a specific Z-direction will be called "image". Every image has to be assigned to a certain wavelength. + * The "AddWavelength" uses "push_back" to write float values into a vector. The first wavelength will correspond to the first image!!! + * If there a more input images 'I' then added wavelengths 'w' the filter base interprets the next x images as repetition sequence of the same + * wavelengths. If I % w !=0 the surplus image(s) will be dropped. + * Addtionaly one has to add chromophores from the property calculator class enum "ChromophoreType" with the "AddChromophore" method. + * This method as well uses "push_back" but the chosen (arbitary) order will be the order of the outputs. + * + * Output: + * The output will be one MITK image per chosen chromophore. Where the XY-plane is a image and the Z-direction represents recordings for different + * sequences. Furthermore it is possible to creat an output image that contains the information about the relative error between unmixing result + * and the input image. + * + * Subclasses: + * - mitkPASpectralUnmixingFilterVigra + * - mitkPALinearSpectralUnmixingFilter (uses Eigen algorithms) + * - mitkPASpectralUnmixingFilterSimplex + * + * Possible exceptions: + * - "PIXELTYPE ERROR! FLOAT 32 REQUIRED": The MITK input image has to consist out of floats. + * - "ERROR! REMOVE WAVELENGTHS!": One needs at least the same amount of images (z-dimension) then added wavelengths. + * - "ADD MORE WAVELENGTHS!": One needs at least the same amount of wavelengths then added chromophores. + * - "WAVELENGTH XXX nm NOT SUPPORTED!": The wavelength is not part of the proptery calculater data base. The data base can be found @ + * [...]\mitk\Modules\PhotoacousticsLib\Resources\spectralLIB.dat + * - "ADD OUTPUTS HAS TO BE LARGER THEN ZERO!" + * - "NO WAVELENGHTS/CHROMOPHORES SELECZED! + * - "INDEX ERROR! NUMBER OF OUTPUTS DOESN'T FIT TO OTHER SETTIGNS!" + */ + + class MITKPHOTOACOUSTICSLIB_EXPORT SpectralUnmixingFilterBase : public mitk::ImageToImageFilter + { + public: + mitkClassMacro(SpectralUnmixingFilterBase, mitk::ImageToImageFilter); + + /** + * \brief AddChromophore takes mitk::pa::PropertyCalculator::ChromophoreType and writes them at the end of the m_Chromophore vector. + * The call of the method sets the order of the GetOutput method! + * + * @param chromophore has to be element of porperty calculater enum chromophore type + * @return for wavelength smaller then 300nm and larger then 1000nm the return will be 0, because not at the data base (2018/06/19) + */ + void AddChromophore(mitk::pa::PropertyCalculator::ChromophoreType chromophore); + + /** + * \brief AddWavelength takes integers and writes them at the end of the m_Wavelength vector. The first call of the method then + * corresponds to the first input image and so on. + * @param wavelength database supports integers between 300 and 1000 nm + */ + void AddWavelength(int wavelength); + + /* + * \brief Verbose gives more information to the console. Default value is false. + * @param verbose is the boolian to activate the MITK_INFO logged to the console + */ + virtual void Verbose(bool verbose); + + /** + * \brief AddOutputs takes an integer and sets indexed outputs + * @param outputs integer correponds to the number of output images + * @throws if outputs == 0 + */ + virtual void AddOutputs(unsigned int outputs); + + /* + * \brief RelativeError returns a image which compare the L2 norm of the input vector with the unmixing result + * @param relativeError is the boolian to activate this tool + */ + virtual void RelativeError(bool relativeError); + + /** + * \brief AddRelativeErrorSettings takes integers and writes them at the end of the m_RelativeErrorSettings vector. + * @param value has to be a integer + */ + virtual void AddRelativeErrorSettings(int value); + + ofstream myfile; // just for testing purposes; has to be removeed + + protected: + /** + * \brief Constructor creats proptery calculater smart pointer new() + */ + SpectralUnmixingFilterBase(); + virtual ~SpectralUnmixingFilterBase(); + + /** + * \brief The subclasses will override the mehtod to calculate the spectral unmixing result vector. + * @param endmemberMatrix Matrix with number of chromophores colums and number of wavelengths rows so matrix element (i,j) contains + * the absorbtion of chromophore j @ wavelength i taken from the database by PropertyElement method. + * @param inputVector Vector containing values of one pixel of XY-plane image with number of wavelength rows (z-dimension of a sequenece) + * so the pixelvalue of the first wavelength is stored in inputVector[0] and so on. + * @throws if algorithm implementiation fails (implemented for the algorithms with critical requirements) + */ + virtual Eigen::VectorXf SpectralUnmixingAlgorithm(Eigen::Matrix endmemberMatrix, + Eigen::VectorXf inputVector) = 0; + + bool m_Verbose = false; + bool m_RelativeError = false; + + std::vector m_Chromophore; + std::vector m_Wavelength; + std::vector m_RelativeErrorSettings; + + private: + /* + * \brief Initialized output images with the same XY-plane size like the input image and total size in z-direction equals number of sequences. + * The pixel type is set to float. + * @param totalNumberOfSequences = (unsigned int) (numberOfInputImages / numberOfWavelength) >= 1 + */ + virtual void InitializeOutputs(unsigned int totalNumberOfSequences); + + /* + * \brief Checks if there are a suitable amount of wavelengths and if the input image consists of floats. + * @param input pointer on the MITK input image + * @throws if there are more wavelength then images + * @throws if there are more chromophores then wavelengths + * @throws if the pixel type is not float 32 + * @throws if no chromophores or wavelengths selected as input + * @throws if the number of indexed outputs does'nt fit to the expected number + */ + virtual void CheckPreConditions(mitk::Image::Pointer input); + + /* + * \brief Inherit from the "ImageToImageFilter" Superclass. Herain it calls InitializeOutputs, CalculateEndmemberMatrix and + * CheckPreConditions methods and enables pixelwise access to do spectral unmixing with the "SpectralUnmixingAlgorithm" + * method. In the end the method writes the results into the new MITK output images. + */ + virtual void GenerateData() override; + + /* + * \brief Creats a Matrix with number of chromophores colums and number of wavelengths rows so matrix element (i,j) contains + * the absorbtion of chromophore j @ wavelength i. The absorbtion values are taken from the "PropertyElement" method. + * @param m_Chromophore is a vector of "PropertyCalculator::ChromophoreType" containing all selected chromophores for the unmixing + * @param m_Wavelength is a vector of integers containing all wavelengths of one sequence + */ + virtual Eigen::Matrix CalculateEndmemberMatrix( + std::vector m_Chromophore, std::vector m_Wavelength); + + /* + * \brief "PropertyElement" is the tool to access the absorbtion values out of the database using mitk::pa::PropertyCalculator::GetAbsorptionForWavelengt + * and checks if the requested wavelength is part of the database (not zero values). The "ONEENDMEMBER" is a pseudo absorber with static absorbtion 1 + * at every wavelength and is therefor not part of the database. If this one is the selected chromophore the return value is 1 for every wavelength. + * @param wavelength has to be integer between 300 and 1000 nm + * @param chromophore has to be a mitk::pa::PropertyCalculator::ChromophoreType + * @throws if mitk::pa::PropertyCalculator::GetAbsorptionForWavelengt returns 0, because this means that the delivered wavelength is not in the database. + */ + virtual float PropertyElement(mitk::pa::PropertyCalculator::ChromophoreType, int wavelength); + + /* + * \brief calculates the relative error between the input image and the unmixing result in the L2 norm + * @param endmemberMatrix is a Eigen matrix containing the endmember information + * @param inputVector is a Eigen vector containing the multispectral information of one pixel + * @param resultVector is a Eigen vector containing the spectral unmmixing result + */ + float CalculateRelativeError(Eigen::Matrix endmemberMatrix, + Eigen::VectorXf inputVector, Eigen::VectorXf resultVector); + + PropertyCalculator::Pointer m_PropertyCalculatorEigen; + }; + } +} +#endif // MITKPHOTOACOUSTICSPECTRALUNMIXINGFILTERBASE_ diff --git a/Modules/PhotoacousticsLib/include/mitkPASpectralUnmixingFilterLagrange.h b/Modules/PhotoacousticsLib/include/mitkPASpectralUnmixingFilterLagrange.h new file mode 100644 index 0000000000..def0212178 --- /dev/null +++ b/Modules/PhotoacousticsLib/include/mitkPASpectralUnmixingFilterLagrange.h @@ -0,0 +1,41 @@ +/*=================================================================== + +The Medical Imaging Interaction Toolkit (MITK) + +Copyright (c) German Cancer Research Center, +Division of Medical and Biological Informatics. +All rights reserved. + +This software is distributed WITHOUT ANY WARRANTY; without +even the implied warranty of MERCHANTABILITY or FITNESS FOR +A PARTICULAR PURPOSE. + +See LICENSE.txt or http://www.mitk.org for details. + +===================================================================*/ + +#ifndef MITKPHOTOACOUSTICSPECTRALUNMIXINGFILTERLAGRANGE_H +#define MITKPHOTOACOUSTICSPECTRALUNMIXINGFILTERLAGRANGE_H + +#include +#include + +namespace mitk { + namespace pa { + class MITKPHOTOACOUSTICSLIB_EXPORT SpectralUnmixingFilterLagrange : public SpectralUnmixingFilterBase + { + public: + + mitkClassMacro(SpectralUnmixingFilterLagrange, SpectralUnmixingFilterBase) + //itkFactorylessNewMacro(Self) + + protected: + SpectralUnmixingFilterLagrange(); + virtual ~SpectralUnmixingFilterLagrange(); + + private: + + }; + } +} +#endif // MITKPHOTOACOUSTICSPECTRALUNMIXINGFILTERLAGRANGE_H diff --git a/Modules/PhotoacousticsLib/include/mitkPASpectralUnmixingFilterSimplex.h b/Modules/PhotoacousticsLib/include/mitkPASpectralUnmixingFilterSimplex.h new file mode 100644 index 0000000000..a5920f4239 --- /dev/null +++ b/Modules/PhotoacousticsLib/include/mitkPASpectralUnmixingFilterSimplex.h @@ -0,0 +1,53 @@ +/*=================================================================== + +The Medical Imaging Interaction Toolkit (MITK) + +Copyright (c) German Cancer Research Center, +Division of Medical and Biological Informatics. +All rights reserved. + +This software is distributed WITHOUT ANY WARRANTY; without +even the implied warranty of MERCHANTABILITY or FITNESS FOR +A PARTICULAR PURPOSE. + +See LICENSE.txt or http://www.mitk.org for details. + +===================================================================*/ + +#ifndef MITKPHOTOACOUSTICSPECTRALUNMIXINGFILTERSIMPLEX_H +#define MITKPHOTOACOUSTICSPECTRALUNMIXINGFILTERSIMPLEX_H + +#include +#include + +namespace mitk { + namespace pa { + class MITKPHOTOACOUSTICSLIB_EXPORT SpectralUnmixingFilterSimplex : public SpectralUnmixingFilterBase + { + public: + + mitkClassMacro(SpectralUnmixingFilterSimplex, SpectralUnmixingFilterBase) + itkFactorylessNewMacro(Self) + + protected: + SpectralUnmixingFilterSimplex(); + virtual ~SpectralUnmixingFilterSimplex(); + + private: + + virtual Eigen::VectorXf SpectralUnmixingAlgorithm(Eigen::Matrix EndmemberMatrix, + Eigen::VectorXf inputVector) override; + + int factorial(int n); + virtual Eigen::Matrix GenerateA(Eigen::Matrix EndmemberMatrix, + Eigen::VectorXf inputVector, int i); + Eigen::Matrix GenerateD2(Eigen::Matrix A); + float simplexVolume(Eigen::Matrix Matrix); + + virtual Eigen::VectorXf Normalization(Eigen::Matrix EndmemberMatrix, + Eigen::VectorXf inputVector); + + }; + } +} +#endif // MITKPHOTOACOUSTICSPECTRALUNMIXINGFILTERSIMPLEX_H diff --git a/Modules/PhotoacousticsLib/include/mitkPASpectralUnmixingFilterVigra.h b/Modules/PhotoacousticsLib/include/mitkPASpectralUnmixingFilterVigra.h new file mode 100644 index 0000000000..a8b7926a5b --- /dev/null +++ b/Modules/PhotoacousticsLib/include/mitkPASpectralUnmixingFilterVigra.h @@ -0,0 +1,108 @@ +/*=================================================================== + +The Medical Imaging Interaction Toolkit (MITK) + +Copyright (c) German Cancer Research Center, +Division of Medical and Biological Informatics. +All rights reserved. + +This software is distributed WITHOUT ANY WARRANTY; without +even the implied warranty of MERCHANTABILITY or FITNESS FOR +A PARTICULAR PURPOSE. + +See LICENSE.txt or http://www.mitk.org for details. + +===================================================================*/ + +#ifndef MITKPHOTOACOUSTICSPECTRALUNMIXINGFILTERVIGRA_H +#define MITKPHOTOACOUSTICSPECTRALUNMIXINGFILTERVIGRA_H + +#include +#include +#include + + + +namespace mitk { + namespace pa { + /** + * \brief This filter is subclass of the spectral unmixing filter base. All algorithms in this class are + * based on the vigra open source c++ library. It takes a multispectral pixel as input and returns a vector + * with the unmixing result. + * + * Input: + * - endmemberMatrix Eigen Matrix with number of chromophores colums and number of wavelengths rows so matrix element (i,j) contains + * the absorbtion of chromophore j @ wavelength i taken from the database by PropertyElement method. + * - inputVector Eigen Vector containing values of one pixel of XY-plane image with number of wavelength rows (z-dimension of a sequenece) + * so the pixelvalue of the first wavelength is stored in inputVector[0] and so on. + * + * Output: + * - Eigen vector with unmixing result of one multispectral pixel. The ith element of the vector corresponds to the ith entry of the + * m_Chromophore vector. + * + * Algorithms (see VigraAlgortihmType enum): + * - LARS algorithm minimizes (A*x-b)^2 s.t. x>=0 using least angle regression + * - GOLDFARB minimizes (A*x-b)^2 s.t. x>=0 using the Goldfarb-Idnani algorithm + * - WEIGHTED minimizes transpose(A*x-b)*diag(weights)*(A*x-b) using QR decomposition + * - LS minimizes (A*x-b)^2 using QR decomposition + * + * Possible exceptions: + * - "404 VIGRA ALGORITHM NOT FOUND": Algorithm not implemented + */ + + class MITKPHOTOACOUSTICSLIB_EXPORT SpectralUnmixingFilterVigra : public SpectralUnmixingFilterBase + { + public: + mitkClassMacro(SpectralUnmixingFilterVigra, SpectralUnmixingFilterBase) + itkFactorylessNewMacro(Self) + + /** + * \brief Contains all implemented Vigra algorithms for spectral unmixing. For detailed information of the algorithms look at the + * mitkPASpectralUnmixingFilterVigra.h documentation (section Algorithms). + */ + enum VigraAlgortihmType + { + LARS, + GOLDFARB, + WEIGHTED, + LS + }; + + /** + * \brief AddWeight takes integers and writes them at the end of the weightsvec vector. The first call of the method then + * corresponds to the first input image/wavelength and so on. + * @param weight is a percentage (integer) between 1 and 100 + */ + void AddWeight(unsigned int weight); + + /** + * \brief Takes a mitk::pa::SpectralUnmixingFilterVigra::VigraAlgortihmType and fix it for usage at the "SpectralUnmixingAlgorithm" + * method. + * @param algorithmName has to be a mitk::pa::SpectralUnmixingFilterVigra::VigraAlgortihmType + */ + void SetAlgorithm(VigraAlgortihmType inputAlgorithmName); + + protected: + SpectralUnmixingFilterVigra(); + virtual ~SpectralUnmixingFilterVigra(); + + /** + * \brief overrides the baseclass method with a mehtod to calculate the spectral unmixing result vector. Herain it first converts the + * Eigen inputs to the Vigra class. Afterwards the class performs the algorithm set by the "SetAlgorithm" method and writes the result + * into a Eigen vector which is the return value. + * @param endmemberMatrix Matrix with number of chromophores colums and number of wavelengths rows so matrix element (i,j) contains + * the absorbtion of chromophore j @ wavelength i taken from the database by PropertyElement method. + * @param inputVector Vector containing values of one pixel of XY-plane image with number of wavelength rows (z-dimension of a sequenece) + * so the pixelvalue of the first wavelength is stored in inputVector[0] and so on. + * @throws if the algorithmName is not a member of the enum VigraAlgortihmType + */ + virtual Eigen::VectorXf SpectralUnmixingAlgorithm(Eigen::Matrix EndmemberMatrix, + Eigen::VectorXf inputVector) override; + + private: + std::vector weightsvec; + SpectralUnmixingFilterVigra::VigraAlgortihmType algorithmName; + }; + } +} +#endif // MITKPHOTOACOUSTICSPECTRALUNMIXINGFILTERVIGRA_H diff --git a/Modules/PhotoacousticsLib/include/mitkPASpectralUnmixingSO2.h b/Modules/PhotoacousticsLib/include/mitkPASpectralUnmixingSO2.h new file mode 100644 index 0000000000..c822e102f7 --- /dev/null +++ b/Modules/PhotoacousticsLib/include/mitkPASpectralUnmixingSO2.h @@ -0,0 +1,106 @@ +/*=================================================================== + +The Medical Imaging Interaction Toolkit (MITK) + +Copyright (c) German Cancer Research Center, +Division of Medical and Biological Informatics. +All rights reserved. + +This software is distributed WITHOUT ANY WARRANTY; without +even the implied warranty of MERCHANTABILITY or FITNESS FOR +A PARTICULAR PURPOSE. + +See LICENSE.txt or http://www.mitk.org for details. + +===================================================================*/ + +#ifndef MITKPHOTOACOUSTICSPECTRALUNMIXINGSO2_H +#define MITKPHOTOACOUSTICSPECTRALUNMIXINGSO2_H + +#include "mitkImageToImageFilter.h" +#include + +//Includes for smart pointer usage +#include "mitkCommon.h" +#include "itkLightObject.h" + + +namespace mitk { + namespace pa { + /** + * \brief derives out of two identical sized MITK images the oxygen saturation and return one MITK image as result. Furthermore + * it is possible to set settings that the result shows just SO2 values above a threshold, or above a input value for Hb, HbO2 to + * get just a oxygen saturation image of interessting structures. + * + * Input: + * The input has to be two 3D MITK images. The order of the inputs matter! The first input has to be the Hb image the second input + * has to be the HbO2 image. The settings are integer values. The SO2 threshold therefore is percentage value. + * + * Output: + * The output will be one MITK image. Where one can see the oxygen saturation of all pixels above the set threholds. If a pixel is + * below a threhold or NAN then the value will be set to zero. + */ + class MITKPHOTOACOUSTICSLIB_EXPORT SpectralUnmixingSO2 : public mitk::ImageToImageFilter + { + public: + mitkClassMacro(SpectralUnmixingSO2, mitk::ImageToImageFilter); + itkFactorylessNewMacro(Self); + + /** + * \brief AddSO2Settings takes integers and writes them at the end of the m_SO2Settings vector. + * @param value of the Setting + */ + virtual void AddSO2Settings(int value); + + /** + * \brief Verbose gives more information to the console. Default value is false. + * @param m_Verbose is the boolian to activate the MITK_INFO logged to the console + */ + virtual void Verbose(bool verbose); + + protected: + /** + * \brief Constructor sets number of input images to two and number of output images to one, respectively. + */ + SpectralUnmixingSO2(); + virtual ~SpectralUnmixingSO2(); + + std::vector m_SO2Settings; + bool m_Verbose = false; + + private: + /** + * \brief Inherit from the "ImageToImageFilter" Superclass. Herain it calls InitializeOutputs and the CheckPreConditions + * methods and enables pixelwise access to the inputs to calculate the oxygen saturation via the "calculate SO2" method. + */ + virtual void GenerateData() override; + + /** + * \brief Initialized output images with the same size like the input image. The pixel type is set to float. + */ + virtual void InitializeOutputs(); + + /** + * \brief Checks if the dimensions of the input images are equal and made out of floating poinhts. + * @throws if the inputs don't have the same size + * @throws if input images don't consist of floats + */ + virtual void CheckPreConditions(mitk::Image::Pointer inputHbO2, mitk::Image::Pointer inputHb); + + /** + * \brief calculates HbO2 / (Hb + HbO2) and afterwards checks if the result is significant (SO2ValueNotSiginificant method). + * If not the method returns zero otherwise it returns the calculated result. + * @param pixelHb is the pixel value of the Hb input. + * @param pixelHb is the pixel value of the Hb input. + * @warn if the HbO2 value is NAN (in patricular if Hb == -HbO2), but result will be set to zero + */ + float CalculateSO2(float pixelHb, float pixelHbO2); + + /** + * \brief return true if SO2 result is not significant by checking if the input values are above the threshold of the settings + */ + bool SO2ValueNotSiginificant(float Hb, float HbO2, float result); + }; + } +} +#endif // MITKPHOTOACOUSTICSPECTRALUNMIXINGSO2_ diff --git a/Modules/PhotoacousticsLib/src/SUFilter/mitkPALinearSpectralUnmixingFilter.cpp b/Modules/PhotoacousticsLib/src/SUFilter/mitkPALinearSpectralUnmixingFilter.cpp new file mode 100644 index 0000000000..ec6b845795 --- /dev/null +++ b/Modules/PhotoacousticsLib/src/SUFilter/mitkPALinearSpectralUnmixingFilter.cpp @@ -0,0 +1,84 @@ +/*=================================================================== + +The Medical Imaging Interaction Toolkit (MITK) + +Copyright (c) German Cancer Research Center, +Division of Medical and Biological Informatics. +All rights reserved. + +This software is distributed WITHOUT ANY WARRANTY; without +even the implied warranty of MERCHANTABILITY or FITNESS FOR +A PARTICULAR PURPOSE. + +See LICENSE.txt or http://www.mitk.org for details. + +===================================================================*/ + +#include "mitkPALinearSpectralUnmixingFilter.h" + +// Testing algorithms +#include + +// ImageAccessor +#include +#include + +mitk::pa::LinearSpectralUnmixingFilter::LinearSpectralUnmixingFilter() +{ +} + +mitk::pa::LinearSpectralUnmixingFilter::~LinearSpectralUnmixingFilter() +{ +} + +void mitk::pa::LinearSpectralUnmixingFilter::SetAlgorithm(mitk::pa::LinearSpectralUnmixingFilter::AlgortihmType inputAlgorithmName) +{ + algorithmName = inputAlgorithmName; +} + +Eigen::VectorXf mitk::pa::LinearSpectralUnmixingFilter::SpectralUnmixingAlgorithm( + Eigen::Matrix endmemberMatrix, Eigen::VectorXf inputVector) +{ + Eigen::VectorXf resultVector; + + if (mitk::pa::LinearSpectralUnmixingFilter::AlgortihmType::HOUSEHOLDERQR == algorithmName) + resultVector = endmemberMatrix.householderQr().solve(inputVector); + + else if (mitk::pa::LinearSpectralUnmixingFilter::AlgortihmType::LDLT == algorithmName) + { + Eigen::LLT lltOfA(endmemberMatrix); + if (lltOfA.info() == Eigen::NumericalIssue) + { + mitkThrow() << "Possibly non semi-positive definitie endmembermatrix!"; + } + else + resultVector = endmemberMatrix.ldlt().solve(inputVector); + } + + else if (mitk::pa::LinearSpectralUnmixingFilter::AlgortihmType::LLT == algorithmName) + { + Eigen::LLT lltOfA(endmemberMatrix); + if (lltOfA.info() == Eigen::NumericalIssue) + { + mitkThrow() << "Possibly non semi-positive definitie endmembermatrix!"; + } + else + resultVector = endmemberMatrix.llt().solve(inputVector); + } + + else if (mitk::pa::LinearSpectralUnmixingFilter::AlgortihmType::COLPIVHOUSEHOLDERQR == algorithmName) + resultVector = endmemberMatrix.colPivHouseholderQr().solve(inputVector); + + else if (mitk::pa::LinearSpectralUnmixingFilter::AlgortihmType::JACOBISVD == algorithmName) + resultVector = endmemberMatrix.jacobiSvd(Eigen::ComputeFullU | Eigen::ComputeFullV).solve(inputVector); + + else if (mitk::pa::LinearSpectralUnmixingFilter::AlgortihmType::FULLPIVLU == algorithmName) + resultVector = endmemberMatrix.fullPivLu().solve(inputVector); + + else if (mitk::pa::LinearSpectralUnmixingFilter::AlgortihmType::FULLPIVHOUSEHOLDERQR == algorithmName) + resultVector = endmemberMatrix.fullPivHouseholderQr().solve(inputVector); + else + mitkThrow() << "404 VIGRA ALGORITHM NOT FOUND"; + + return resultVector; +} diff --git a/Modules/PhotoacousticsLib/src/SUFilter/mitkPASpectralUnmixingFilterBase.cpp b/Modules/PhotoacousticsLib/src/SUFilter/mitkPASpectralUnmixingFilterBase.cpp new file mode 100644 index 0000000000..798354dc55 --- /dev/null +++ b/Modules/PhotoacousticsLib/src/SUFilter/mitkPASpectralUnmixingFilterBase.cpp @@ -0,0 +1,239 @@ +/*=================================================================== + +The Medical Imaging Interaction Toolkit (MITK) + +Copyright (c) German Cancer Research Center, +Division of Medical and Biological Informatics. +All rights reserved. + +This software is distributed WITHOUT ANY WARRANTY; without +even the implied warranty of MERCHANTABILITY or FITNESS FOR +A PARTICULAR PURPOSE. + +See LICENSE.txt or http://www.mitk.org for details. + +===================================================================*/ + +#include "mitkPASpectralUnmixingFilterBase.h" + +// Includes for AddEnmemberMatrix +#include "mitkPAPropertyCalculator.h" +#include + +// ImageAccessor +#include +#include + +mitk::pa::SpectralUnmixingFilterBase::SpectralUnmixingFilterBase() +{ + m_PropertyCalculatorEigen = mitk::pa::PropertyCalculator::New(); +} + +mitk::pa::SpectralUnmixingFilterBase::~SpectralUnmixingFilterBase() +{ + +} + +void mitk::pa::SpectralUnmixingFilterBase::AddOutputs(unsigned int outputs) +{ + if (outputs == 0) + mitkThrow() << "ADD OUTPUTS HAS TO BE LARGER THEN ZERO!"; + this->SetNumberOfIndexedOutputs(outputs); + for (unsigned int i = 0; iSetNthOutput(i, mitk::Image::New()); +} + +void mitk::pa::SpectralUnmixingFilterBase::AddWavelength(int wavelength) +{ + m_Wavelength.push_back(wavelength); +} + +void mitk::pa::SpectralUnmixingFilterBase::AddChromophore(mitk::pa::PropertyCalculator::ChromophoreType chromophore) +{ + m_Chromophore.push_back(chromophore); +} + +void mitk::pa::SpectralUnmixingFilterBase::Verbose(bool verbose) +{ + m_Verbose = verbose; +} + +void mitk::pa::SpectralUnmixingFilterBase::RelativeError(bool relativeError) +{ + m_RelativeError = relativeError; +} + +void mitk::pa::SpectralUnmixingFilterBase::AddRelativeErrorSettings(int value) +{ + m_RelativeErrorSettings.push_back(value); +} + +void mitk::pa::SpectralUnmixingFilterBase::GenerateData() +{ + MITK_INFO(m_Verbose) << "GENERATING DATA.."; + + mitk::Image::Pointer input = GetInput(0); + CheckPreConditions(input); + + unsigned int xDim = input->GetDimensions()[0]; + unsigned int yDim = input->GetDimensions()[1]; + unsigned int numberOfInputImages = input->GetDimensions()[2]; + + MITK_INFO(m_Verbose) << "x dimension: " << xDim; + MITK_INFO(m_Verbose) << "y dimension: " << yDim; + MITK_INFO(m_Verbose) << "z dimension: " << numberOfInputImages; + + // Copy input image into array + mitk::ImageReadAccessor readAccess(input); + const float* inputDataArray = ((const float*)readAccess.GetData()); + + unsigned int sequenceSize = m_Wavelength.size(); + unsigned int totalNumberOfSequences = numberOfInputImages / sequenceSize; + MITK_INFO(m_Verbose) << "TotalNumberOfSequences: " << totalNumberOfSequences; + + InitializeOutputs(totalNumberOfSequences); + + auto endmemberMatrix = CalculateEndmemberMatrix(m_Chromophore, m_Wavelength); + + // test to see pixel values @ txt file + myfile.open("SimplexNormalisation.txt"); + + unsigned int outputCounter = GetNumberOfIndexedOutputs(); + std::vector writteBufferVector; + for (unsigned int i = 0; i < outputCounter; ++i) + { + auto output = GetOutput(i); + mitk::ImageWriteAccessor writeOutput(output); + float* writeBuffer = (float *)writeOutput.GetData(); + writteBufferVector.push_back(writeBuffer); + } + + if (m_RelativeError == true) + { + // -1 because rel error is output[IndexedOutputs() - 1] and loop over chromophore outputs has to end at [IndexedOutputs() - 2] + outputCounter -= 1; + } + + for (unsigned int sequenceCounter = 0; sequenceCounter < totalNumberOfSequences; ++sequenceCounter) + { + MITK_INFO(m_Verbose) << "SequenceCounter: " << sequenceCounter; + //loop over every pixel in XY-plane + for (unsigned int x = 0; x < xDim; x++) + { + for (unsigned int y = 0; y < yDim; y++) + { + Eigen::VectorXf inputVector(sequenceSize); + for (unsigned int z = 0; z < sequenceSize; z++) + { + /** + * 'sequenceCounter*sequenceSize' has to be added to 'z' to ensure that one accesses the + * correct pixel, because the inputDataArray contains the information of all sequences and + * not just the one of the current sequence. + */ + unsigned int pixelNumber = (xDim*yDim*(z+sequenceCounter*sequenceSize)) + x * yDim + y; + auto pixel = inputDataArray[pixelNumber]; + + inputVector[z] = pixel; + } + Eigen::VectorXf resultVector = SpectralUnmixingAlgorithm(endmemberMatrix, inputVector); + + if (m_RelativeError == true) + { + float relativeError = CalculateRelativeError(endmemberMatrix, inputVector, resultVector); + writteBufferVector[outputCounter][(xDim*yDim * sequenceCounter) + x * yDim + y] = relativeError; + } + + for (unsigned int outputIdx = 0; outputIdx < outputCounter; ++outputIdx) + { + writteBufferVector[outputIdx][(xDim*yDim * sequenceCounter) + x * yDim + y] = resultVector[outputIdx]; + } + } + } + } + MITK_INFO(m_Verbose) << "GENERATING DATA...[DONE]"; + myfile.close(); +} + +void mitk::pa::SpectralUnmixingFilterBase::CheckPreConditions(mitk::Image::Pointer input) +{ + MITK_INFO(m_Verbose) << "CHECK PRECONDITIONS ..."; + + if (m_Chromophore.size() == 0 || m_Wavelength.size() == 0) + mitkThrow() << "NO WAVELENGHTS/CHROMOPHORES SELECTED!"; + + if (m_Wavelength.size() < input->GetDimensions()[2]) + MITK_WARN(m_Verbose) << "NUMBER OF WAVELENGTHS < NUMBER OF INPUT IMAGES"; + + if (m_Chromophore.size() > m_Wavelength.size()) + mitkThrow() << "ADD MORE WAVELENGTHS OR REMOVE ENDMEMBERS!"; + + if (input->GetPixelType() != mitk::MakeScalarPixelType()) + mitkThrow() << "PIXELTYPE ERROR! FLOAT REQUIRED"; + + if ((m_Chromophore.size()+ m_RelativeError )!= GetNumberOfIndexedOutputs() || input->GetDimensions()[2] < GetNumberOfIndexedOutputs()) + mitkThrow() << "INDEX ERROR! NUMBER OF OUTPUTS DOESN'T FIT TO OTHER SETTIGNS!"; + + MITK_INFO(m_Verbose) << "...[DONE]"; +} + +void mitk::pa::SpectralUnmixingFilterBase::InitializeOutputs(unsigned int totalNumberOfSequences) +{ + MITK_INFO(m_Verbose) << "Initialize Outputs ..."; + + unsigned int numberOfInputs = GetNumberOfIndexedInputs(); + unsigned int numberOfOutputs = GetNumberOfIndexedOutputs(); + MITK_INFO(m_Verbose) << "Inputs: " << numberOfInputs << " Outputs: " << numberOfOutputs; + + mitk::PixelType pixelType = mitk::MakeScalarPixelType(); + const int NUMBER_OF_SPATIAL_DIMENSIONS = 3; + auto* dimensions = new unsigned int[NUMBER_OF_SPATIAL_DIMENSIONS]; + for (unsigned int dimIdx = 0; dimIdx < 2; dimIdx++) + dimensions[dimIdx] = GetInput()->GetDimensions()[dimIdx]; + dimensions[2] = totalNumberOfSequences; + + for (unsigned int outputIdx = 0; outputIdx < numberOfOutputs; outputIdx++) + GetOutput(outputIdx)->Initialize(pixelType, NUMBER_OF_SPATIAL_DIMENSIONS, dimensions); + MITK_INFO(m_Verbose) << "...[DONE]"; +} + +Eigen::Matrix mitk::pa::SpectralUnmixingFilterBase::CalculateEndmemberMatrix( + std::vector m_Chromophore, std::vector m_Wavelength) +{ + unsigned int numberOfChromophores = m_Chromophore.size(); //columns + unsigned int numberOfWavelengths = m_Wavelength.size(); //rows + Eigen::Matrix endmemberMatrixEigen(numberOfWavelengths, numberOfChromophores); + + for (unsigned int j = 0; j < numberOfChromophores; ++j) + { + for (unsigned int i = 0; i < numberOfWavelengths; ++i) + endmemberMatrixEigen(i, j) = PropertyElement(m_Chromophore[j], m_Wavelength[i]); + } + MITK_INFO(m_Verbose) << "GENERATING ENMEMBERMATRIX [DONE]"; + return endmemberMatrixEigen; +} + +float mitk::pa::SpectralUnmixingFilterBase::PropertyElement(mitk::pa::PropertyCalculator::ChromophoreType chromophore, int wavelength) +{ + if (chromophore == mitk::pa::PropertyCalculator::ChromophoreType::ONEENDMEMBER) + return 1; + else + { + float value = m_PropertyCalculatorEigen->GetAbsorptionForWavelength(chromophore, wavelength); + if (value == 0) + mitkThrow() << "WAVELENGTH " << wavelength << "nm NOT SUPPORTED!"; + else + return value; + } +} + +float mitk::pa::SpectralUnmixingFilterBase::CalculateRelativeError(Eigen::Matrix endmemberMatrix, + Eigen::VectorXf inputVector, Eigen::VectorXf resultVector) +{ + float relativeError = (endmemberMatrix*resultVector - inputVector).norm() / inputVector.norm(); + for (int i = 0; i < 2; ++i) + { + if (resultVector[i] < m_RelativeErrorSettings[i]) + return 0; + } + return relativeError; +} diff --git a/Modules/PhotoacousticsLib/src/SUFilter/mitkPASpectralUnmixingFilterLagrange.cpp b/Modules/PhotoacousticsLib/src/SUFilter/mitkPASpectralUnmixingFilterLagrange.cpp new file mode 100644 index 0000000000..bfecfca51f --- /dev/null +++ b/Modules/PhotoacousticsLib/src/SUFilter/mitkPASpectralUnmixingFilterLagrange.cpp @@ -0,0 +1,32 @@ +/*=================================================================== + +The Medical Imaging Interaction Toolkit (MITK) + +Copyright (c) German Cancer Research Center, +Division of Medical and Biological Informatics. +All rights reserved. + +This software is distributed WITHOUT ANY WARRANTY; without +even the implied warranty of MERCHANTABILITY or FITNESS FOR +A PARTICULAR PURPOSE. + +See LICENSE.txt or http://www.mitk.org for details. + +===================================================================*/ + +#include "mitkPASpectralUnmixingFilterLagrange.h" + +// Includes for AddEnmemberMatrix +#include "mitkPAPropertyCalculator.h" +#include + + +mitk::pa::SpectralUnmixingFilterLagrange::SpectralUnmixingFilterLagrange() +{ + +} + +mitk::pa::SpectralUnmixingFilterLagrange::~SpectralUnmixingFilterLagrange() +{ + +} diff --git a/Modules/PhotoacousticsLib/src/SUFilter/mitkPASpectralUnmixingFilterSimplex.cpp b/Modules/PhotoacousticsLib/src/SUFilter/mitkPASpectralUnmixingFilterSimplex.cpp new file mode 100644 index 0000000000..7fb91f3159 --- /dev/null +++ b/Modules/PhotoacousticsLib/src/SUFilter/mitkPASpectralUnmixingFilterSimplex.cpp @@ -0,0 +1,176 @@ +/*=================================================================== + +The Medical Imaging Interaction Toolkit (MITK) + +Copyright (c) German Cancer Research Center, +Division of Medical and Biological Informatics. +All rights reserved. + +This software is distributed WITHOUT ANY WARRANTY; without +even the implied warranty of MERCHANTABILITY or FITNESS FOR +A PARTICULAR PURPOSE. + +See LICENSE.txt or http://www.mitk.org for details. + +===================================================================*/ + +#include + +// Includes for AddEnmemberMatrix +#include + +#include + +mitk::pa::SpectralUnmixingFilterSimplex::SpectralUnmixingFilterSimplex() +{ + +} + +mitk::pa::SpectralUnmixingFilterSimplex::~SpectralUnmixingFilterSimplex() +{ + +} + +Eigen::VectorXf mitk::pa::SpectralUnmixingFilterSimplex::SpectralUnmixingAlgorithm( + Eigen::Matrix EndmemberMatrix, Eigen::VectorXf inputVector) +{ + int numberOfChromophores = EndmemberMatrix.cols(); + + Eigen::VectorXf resultVector(numberOfChromophores); + Eigen::VectorXf normalizedInputVector(EndmemberMatrix.rows()); + normalizedInputVector = Normalization(EndmemberMatrix, inputVector); + //normalizedInputVector = inputVector; + + + float VolumeMax = simplexVolume(EndmemberMatrix); + for (int i = 0; i < numberOfChromophores; ++i) + { + Eigen::Matrix A = GenerateA(EndmemberMatrix, normalizedInputVector, i); + float Volume = simplexVolume(A); + + + resultVector[i] = Volume / VolumeMax; + + myfile << "resultVector["< mitk::pa::SpectralUnmixingFilterSimplex::GenerateA +(Eigen::Matrix EndmemberMatrix, Eigen::VectorXf inputVector, int i) +{ + Eigen::Matrix A = EndmemberMatrix; + int numberOfChromophores = EndmemberMatrix.cols(); + + for (int j = 0; j < numberOfChromophores; ++j) + { + A(i, j) = inputVector(j); + } + + return A; +} + +Eigen::Matrix mitk::pa::SpectralUnmixingFilterSimplex::GenerateD2 +(Eigen::Matrix A) +{ + int numberOfChromophores = A.cols(); + + Eigen::Matrix D2(numberOfChromophores, numberOfChromophores); + + for (int i = 0; i < numberOfChromophores; ++i) + { + for (int j = 0; j < numberOfChromophores; ++j) + { + Eigen::VectorXf x = A.col(i) - A.col(j); + //MITK_INFO << "a_col_i: " < Matrix) +{ + float Volume; + int numberOfChromophores = Matrix.cols(); + Eigen::Matrix C(numberOfChromophores + 1, numberOfChromophores + 1); + Eigen::Matrix D2 = GenerateD2(Matrix); + + for (int i = 0; i < numberOfChromophores; ++i) + { + for (int j = 0; j < numberOfChromophores; ++j) + { + C(i, j) = D2(i, j); + } + C(i, numberOfChromophores) = 1; + for (int k = 0; k < numberOfChromophores; ++k) + { + C(numberOfChromophores, k) = 1; + } + C(numberOfChromophores, numberOfChromophores) = 0; + } + + float detC = -C.determinant();// determinate von C + float denominator = (factorial(numberOfChromophores - 1)) ^ 2 * 2 ^ (numberOfChromophores - 1)*(-1) ^ numberOfChromophores; + Volume = std::sqrt(detC / denominator); + //MITK_INFO << "detC: " << detC; + + //MITK_INFO << "denominator: " << denominator; + + //MITK_INFO << "Volume: " << Volume; + + return Volume; +} + +int mitk::pa::SpectralUnmixingFilterSimplex::factorial(int n) +{ + if (n == 1) + return 1; + else + return n * factorial(n - 1); +} + +Eigen::VectorXf mitk::pa::SpectralUnmixingFilterSimplex::Normalization( + Eigen::Matrix EndmemberMatrix, Eigen::VectorXf inputVector) +{ + int numberOfWavelengths = inputVector.rows(); + Eigen::VectorXf result(numberOfWavelengths); + float NormalizationFactor = 1; + float foo; + float norm = 0; + + for (int i = 0; i < numberOfWavelengths; ++i) + { + foo = EndmemberMatrix(i, 0) - EndmemberMatrix(i, 1); + if (std::abs(foo) > norm) + norm = std::abs(foo); + } + +//ofstream myfile; +//myfile.open("SimplexNormalisation.txt"); + //NormalizationFactor = inputVector[0] * 2 / norm; + myfile << "Normalizationfactor " << NormalizationFactor << "\n"; + + for (int i = 0; i < numberOfWavelengths; ++i) + { + + result[i]=(inputVector[i]/NormalizationFactor); + } +//myfile.close(); + + return result; +} diff --git a/Modules/PhotoacousticsLib/src/SUFilter/mitkPASpectralUnmixingFilterVigra.cpp b/Modules/PhotoacousticsLib/src/SUFilter/mitkPASpectralUnmixingFilterVigra.cpp new file mode 100644 index 0000000000..f95d558033 --- /dev/null +++ b/Modules/PhotoacousticsLib/src/SUFilter/mitkPASpectralUnmixingFilterVigra.cpp @@ -0,0 +1,104 @@ +/*=================================================================== + +The Medical Imaging Interaction Toolkit (MITK) + +Copyright (c) German Cancer Research Center, +Division of Medical and Biological Informatics. +All rights reserved. + +This software is distributed WITHOUT ANY WARRANTY; without +even the implied warranty of MERCHANTABILITY or FITNESS FOR +A PARTICULAR PURPOSE. + +See LICENSE.txt or http://www.mitk.org for details. + +===================================================================*/ + +#include "mitkPASpectralUnmixingFilterVigra.h" + +// ImageAccessor +#include +#include + +//vigra +#include +#include +#include + +mitk::pa::SpectralUnmixingFilterVigra::SpectralUnmixingFilterVigra() +{ +} + +mitk::pa::SpectralUnmixingFilterVigra::~SpectralUnmixingFilterVigra() +{ +} + +void mitk::pa::SpectralUnmixingFilterVigra::SetAlgorithm(mitk::pa::SpectralUnmixingFilterVigra::VigraAlgortihmType inputAlgorithmName) +{ + algorithmName = inputAlgorithmName; +} + +void mitk::pa::SpectralUnmixingFilterVigra::AddWeight(unsigned int weight) +{ + double value = double(weight) / 100.0; + weightsvec.push_back(value); +} + +Eigen::VectorXf mitk::pa::SpectralUnmixingFilterVigra::SpectralUnmixingAlgorithm( + Eigen::Matrix endmemberMatrix, Eigen::VectorXf inputVector) +{ + unsigned int numberOfWavelengths = endmemberMatrix.rows(); + unsigned int numberOfChromophores = endmemberMatrix.cols(); + + // writes endmemberMatrix and inputVector into vigra::Matrix + std::vector aData; + std::vector bData; + for (unsigned int i = 0; i < numberOfWavelengths; ++i) + { + bData.push_back((double)inputVector(i)); + for (unsigned int j = 0; j < numberOfChromophores; ++j) + aData.push_back((double)endmemberMatrix(i, j)); + } + const double* aDat = aData.data(); + const double* bDat = bData.data(); + + vigra::Matrix A(vigra::Shape2(numberOfWavelengths, numberOfChromophores), aDat); + vigra::Matrix b(vigra::Shape2(numberOfWavelengths, 1), bDat); + vigra::Matrix x(vigra::Shape2(numberOfChromophores, 1)); + + if (mitk::pa::SpectralUnmixingFilterVigra::VigraAlgortihmType::LARS == algorithmName) + nonnegativeLeastSquares(A, b, x); + + else if (mitk::pa::SpectralUnmixingFilterVigra::VigraAlgortihmType::GOLDFARB == algorithmName) + { + vigra::linalg::Matrix eye(vigra::linalg::identityMatrix(numberOfChromophores)), + zeros(vigra::Shape2(numberOfChromophores, 1)), + empty, + U = vigra::linalg::transpose(A)*A, + // v= -transpose(A)*b replaced by -v used in "quadraticProgramming" + v = vigra::linalg::transpose(A)*b; + x = 0; + quadraticProgramming(U, -v, empty, empty, eye, zeros, x); + } + + else if (mitk::pa::SpectralUnmixingFilterVigra::VigraAlgortihmType::WEIGHTED == algorithmName) + { + if (weightsvec.size() != numberOfWavelengths) + mitkThrow() << "Number of weights and wavelengths doesn't match! OR Invalid weight!"; + const double* weightsdat = weightsvec.data(); + vigra::Matrix weigths(vigra::Shape2(numberOfWavelengths, 1), weightsdat); + vigra::linalg::weightedLeastSquares(A, b, weigths, x); + } + + else if (mitk::pa::SpectralUnmixingFilterVigra::VigraAlgortihmType::LS == algorithmName) + linearSolve(A, b, x); + + else + mitkThrow() << "404 VIGRA ALGORITHM NOT FOUND"; + + Eigen::VectorXf resultVector(numberOfChromophores); + for (unsigned int k = 0; k < numberOfChromophores; ++k) + resultVector[k] = (float)x(k, 0); + + return resultVector; +} diff --git a/Modules/PhotoacousticsLib/src/SUFilter/mitkPASpectralUnmixingSO2.cpp b/Modules/PhotoacousticsLib/src/SUFilter/mitkPASpectralUnmixingSO2.cpp new file mode 100644 index 0000000000..0c9bc04bac --- /dev/null +++ b/Modules/PhotoacousticsLib/src/SUFilter/mitkPASpectralUnmixingSO2.cpp @@ -0,0 +1,165 @@ +/*=================================================================== + +The Medical Imaging Interaction Toolkit (MITK) + +Copyright (c) German Cancer Research Center, +Division of Medical and Biological Informatics. +All rights reserved. + +This software is distributed WITHOUT ANY WARRANTY; without +even the implied warranty of MERCHANTABILITY or FITNESS FOR +A PARTICULAR PURPOSE. + +See LICENSE.txt or http://www.mitk.org for details. + +===================================================================*/ + +#include "mitkPASpectralUnmixingSO2.h" + +// ImageAccessor +#include +#include + + +mitk::pa::SpectralUnmixingSO2::SpectralUnmixingSO2() +{ + this->SetNumberOfIndexedInputs(2); + this->SetNumberOfIndexedOutputs(1); + this->SetNthOutput(0, mitk::Image::New()); +} + +mitk::pa::SpectralUnmixingSO2::~SpectralUnmixingSO2() +{ + +} + +void mitk::pa::SpectralUnmixingSO2::Verbose(bool verbose) +{ + m_Verbose = verbose; +} + +void mitk::pa::SpectralUnmixingSO2::GenerateData() +{ + MITK_INFO(m_Verbose) << "GENERATING DATA.."; + + // Get input image + mitk::Image::Pointer inputHbO2 = GetInput(0); + mitk::Image::Pointer inputHb = GetInput(1); + + CheckPreConditions(inputHbO2, inputHb); + + unsigned int xDim = inputHbO2->GetDimensions()[0]; + unsigned int yDim = inputHbO2->GetDimensions()[1]; + unsigned int zDim = inputHbO2->GetDimensions()[2]; + + InitializeOutputs(); + + mitk::ImageReadAccessor readAccessHbO2(inputHbO2); + mitk::ImageReadAccessor readAccessHb(inputHb); + + const float* inputDataArrayHbO2 = ((const float*)readAccessHbO2.GetData()); + const float* inputDataArrayHb = ((const float*)readAccessHb.GetData()); + + auto output = GetOutput(0); + mitk::ImageWriteAccessor writeOutput(output); + float* writeBuffer = (float *)writeOutput.GetData(); + + for (unsigned int x = 0; x < xDim; x++) + { + for (unsigned int y = 0; y < yDim; y++) + { + for (unsigned int z = 0;z < zDim; z++) + { + unsigned int pixelNumber = (xDim*yDim * z) + x * yDim + y; + float pixelHb = inputDataArrayHb[pixelNumber]; + float pixelHbO2 = inputDataArrayHbO2[pixelNumber]; + float resultSO2 = CalculateSO2(pixelHb, pixelHbO2); + writeBuffer[(xDim*yDim * z) + x * yDim + y] = resultSO2; + } + } + } + MITK_INFO(m_Verbose) << "GENERATING DATA...[DONE]"; +} + +void mitk::pa::SpectralUnmixingSO2::CheckPreConditions(mitk::Image::Pointer inputHbO2, mitk::Image::Pointer inputHb) +{ + unsigned int xDimHb = inputHb->GetDimensions()[0]; + unsigned int yDimHb = inputHb->GetDimensions()[1]; + unsigned int zDimHb = inputHb->GetDimensions()[2]; + + unsigned int xDimHbO2 = inputHbO2->GetDimensions()[0]; + unsigned int yDimHbO2 = inputHbO2->GetDimensions()[1]; + unsigned int zDimHbO2 = inputHbO2->GetDimensions()[2]; + + if (xDimHb != xDimHbO2 || yDimHb != yDimHbO2 || zDimHb != zDimHbO2) + mitkThrow() << "DIMENTIONALITY ERROR!"; + + if (inputHbO2->GetPixelType() != mitk::MakeScalarPixelType()) + mitkThrow() << "PIXELTYPE ERROR! FLOAT REQUIRED"; + + if (inputHb->GetPixelType() != mitk::MakeScalarPixelType()) + mitkThrow() << "PIXELTYPE ERROR! FLOAT REQUIRED"; + + MITK_INFO(m_Verbose) << "CHECK PRECONDITIONS ...[DONE]"; +} + +void mitk::pa::SpectralUnmixingSO2::InitializeOutputs() +{ + // UNUSED unsigned int numberOfInputs = GetNumberOfIndexedInputs(); + unsigned int numberOfOutputs = GetNumberOfIndexedOutputs(); + + mitk::PixelType pixelType = mitk::MakeScalarPixelType(); + const int NUMBER_OF_SPATIAL_DIMENSIONS = 3; + auto* dimensions = new unsigned int[NUMBER_OF_SPATIAL_DIMENSIONS]; + for(unsigned int dimIdx=0; dimIdxGetDimensions()[dimIdx]; + } + + for (unsigned int outputIdx = 0; outputIdx < numberOfOutputs; outputIdx++) + { + GetOutput(outputIdx)->Initialize(pixelType, NUMBER_OF_SPATIAL_DIMENSIONS, dimensions); + } +} + +float mitk::pa::SpectralUnmixingSO2::CalculateSO2(float Hb, float HbO2) +{ + float result = HbO2 / (Hb + HbO2); + + if (result != result) + { + MITK_WARN(m_Verbose) << "SO2 VALUE NAN! WILL BE SET TO ZERO!"; + return 0; + } + else + { + if (SO2ValueNotSiginificant(Hb, HbO2, result)) + { + return 0; + } + else return result; + } +} + +void mitk::pa::SpectralUnmixingSO2::AddSO2Settings(int value) +{ + m_SO2Settings.push_back(value); +} + +bool mitk::pa::SpectralUnmixingSO2::SO2ValueNotSiginificant(float Hb, float HbO2, float result) +{ + std::vector significant; + significant.push_back(HbO2); + significant.push_back(Hb); + significant.push_back(HbO2 + Hb); + significant.push_back(100*(result)); + + for (unsigned int i = 0; i < m_SO2Settings.size(); ++i) + { + if (m_SO2Settings[i] != 0 && m_SO2Settings[i] > significant[i] && (std::abs(m_SO2Settings[i] - significant[i]) > 1e-7)) + { + return true; + } + } + return false; +} diff --git a/Modules/PhotoacousticsLib/src/Utils/mitkPAPropertyCalculator.cpp b/Modules/PhotoacousticsLib/src/Utils/mitkPAPropertyCalculator.cpp index 03922d5835..1eeeb14258 100644 --- a/Modules/PhotoacousticsLib/src/Utils/mitkPAPropertyCalculator.cpp +++ b/Modules/PhotoacousticsLib/src/Utils/mitkPAPropertyCalculator.cpp @@ -1,155 +1,161 @@ /*=================================================================== The Medical Imaging Interaction Toolkit (MITK) Copyright (c) German Cancer Research Center, Division of Medical and Biological Informatics. All rights reserved. This software is distributed WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See LICENSE.txt or http://www.mitk.org for details. ===================================================================*/ #include "mitkPAPropertyCalculator.h" // us #include #include #include #include #include mitk::pa::PropertyCalculator::Properties mitk::pa::PropertyCalculator:: CalculatePropertyForSpecificWavelength( TissueType tissueType, int wavelength, double bloodOxygenInFraction) { Properties returnValue; if (!m_Valid) { mitkThrow() << "PhotoacousticPropertyGenerator was not loaded properly."; return returnValue; } double bloodOxygenation = bloodOxygenInFraction; double bloodVolumeFraction; double waterVolumeFraction; double fatVolumeFraction; double melanosomesVolumeFraction; double musp500; double fray; double bmie; double g; switch (tissueType) { case AIR: returnValue.mua = 1e-4; returnValue.mus = 1.0; returnValue.g = 1.0; return returnValue; case BLOOD: bloodVolumeFraction = 1.0; waterVolumeFraction = 0.95; fatVolumeFraction = 0.02; melanosomesVolumeFraction = 0; musp500 = 10; fray = 0; bmie = 1; g = 0.9; break; case EPIDERMIS: bloodVolumeFraction = 0; bloodOxygenation = 0; waterVolumeFraction = 0.75; fatVolumeFraction = 0.01; melanosomesVolumeFraction = 0.02; musp500 = 40; fray = 0; bmie = 1; g = 0.9; break; case FAT: bloodVolumeFraction = 0.01; bloodOxygenation = 0.75; waterVolumeFraction = 0.5; fatVolumeFraction = 0.6; melanosomesVolumeFraction = 0; musp500 = 17.47; fray = 0.2; bmie = 0.5; g = 0.9; break; case STANDARD_TISSUE: default: bloodVolumeFraction = 0.02; bloodOxygenation = 0.75; waterVolumeFraction = 0.8; fatVolumeFraction = 0.05; melanosomesVolumeFraction = 0; musp500 = 25; fray = 0.25; bmie = 1; g = 0.9; break; } // We want the reduced scattering coefficient directly. double musp = musp500 * (fray * pow(wavelength / 500.0, -4.0) + ((1 - fray) * pow(wavelength / 500.0, -bmie))); returnValue.mus = musp; returnValue.mus = 15;//musp; - double mua = bloodVolumeFraction*bloodOxygenation*m_SpectralLibMap[MapType::OXYGENATED][wavelength] + - bloodVolumeFraction*(1 - bloodOxygenation)*m_SpectralLibMap[MapType::DEOXYGENATED][wavelength] + - waterVolumeFraction*m_SpectralLibMap[MapType::WATER][wavelength] + - fatVolumeFraction*m_SpectralLibMap[MapType::FATTY][wavelength] + - melanosomesVolumeFraction*m_SpectralLibMap[MapType::MELANIN][wavelength]; + double mua = bloodVolumeFraction*bloodOxygenation*m_SpectralLibMap[ChromophoreType::OXYGENATED][wavelength] + + bloodVolumeFraction*(1 - bloodOxygenation)*m_SpectralLibMap[ChromophoreType::DEOXYGENATED][wavelength] + + waterVolumeFraction*m_SpectralLibMap[ChromophoreType::WATER][wavelength] + + fatVolumeFraction*m_SpectralLibMap[ChromophoreType::FATTY][wavelength] + + melanosomesVolumeFraction*m_SpectralLibMap[ChromophoreType::MELANIN][wavelength]; returnValue.mua = mua; returnValue.g = g; return returnValue; } mitk::pa::PropertyCalculator::PropertyCalculator() { us::ModuleResource resource = us::GetModuleContext()->GetModule()->GetResource("spectralLIB.dat"); if (resource.IsValid()) { us::ModuleResourceStream resourceStream(resource); std::string line; int wavelength = 300; int counter = 0; while (std::getline(resourceStream, line, '\t')) { int currentLineIdx = counter % 6; if (currentLineIdx == 0) wavelength = std::stoi(line); else { std::istringstream lineStream(line); double tempDouble; lineStream >> tempDouble; m_SpectralLibMap[currentLineIdx][wavelength] = tempDouble; } ++counter; } } else { m_Valid = false; } m_Valid = true; } mitk::pa::PropertyCalculator::~PropertyCalculator() { m_SpectralLibMap.clear(); m_Valid = false; } + +double mitk::pa::PropertyCalculator::GetAbsorptionForWavelength( + ChromophoreType ChromophoreType, int wavelength) +{ + return m_SpectralLibMap[ChromophoreType][wavelength]; +} diff --git a/Modules/PhotoacousticsLib/test/files.cmake b/Modules/PhotoacousticsLib/test/files.cmake index 0796a506a5..b3e91a9e1c 100644 --- a/Modules/PhotoacousticsLib/test/files.cmake +++ b/Modules/PhotoacousticsLib/test/files.cmake @@ -1,41 +1,42 @@ set(MODULE_TESTS # IMPORTANT: If you plan to deactivate / comment out a test please write a bug number to the commented out line of code. # # Example: #mitkMyTest #this test is commented out because of bug 12345 # # It is important that the bug is open and that the test will be activated again before the bug is closed. This assures that # no test is forgotten after it was commented out. If there is no bug for your current problem, please add a new one and # mark it as critical. ################## ON THE FENCE TESTS ################################################# # none ################## DISABLED TESTS ##################################################### ################# RUNNING TESTS ####################################################### mitkSlicedVolumeGeneratorTest.cpp mitkPhotoacousticTissueGeneratorTest.cpp mitkPhotoacousticVectorTest.cpp mitkPhotoacoustic3dVolumeTest.cpp mitkPhotoacousticVolumeTest.cpp mitkPhotoacousticVesselTest.cpp mitkPhotoacousticVesselTreeTest.cpp mitkPhotoacousticVesselMeanderStrategyTest.cpp mitkMcxyzXmlTest.cpp mitkPhotoacousticComposedVolumeTest.cpp mitkPhotoacousticNoiseGeneratorTest.cpp mitkPhotoacousticIOTest.cpp mitkMCThreadHandlerTest.cpp mitkSimulationBatchGeneratorTest.cpp mitkPropertyCalculatorTest.cpp + mitkSpectralUnmixingTest.cpp mitkPhotoacousticImageSavingTest.cpp ) set(RESOURCE_FILES pointsource.xml circlesource.xml rectanglesource.xml twopointsources.xml allsources.xml ) diff --git a/Modules/PhotoacousticsLib/test/mitkSpectralUnmixingTest.cpp b/Modules/PhotoacousticsLib/test/mitkSpectralUnmixingTest.cpp new file mode 100644 index 0000000000..ffdede367a --- /dev/null +++ b/Modules/PhotoacousticsLib/test/mitkSpectralUnmixingTest.cpp @@ -0,0 +1,684 @@ +//*=================================================================== + +//The Medical Imaging Interaction Toolkit (MITK) + +//Copyright (c) German Cancer Research Center, +//Division of Medical and Biological Informatics. +//All rights reserved. + +//This software is distributed WITHOUT ANY WARRANTY; without +//even the implied warranty of MERCHANTABILITY or FITNESS FOR +//A PARTICULAR PURPOSE. + +//See LICENSE.txt or http://www.mitk.org for details. + +//===================================================================*/ +#include +#include + +#include +#include +#include +#include +#include +#include + +class mitkSpectralUnmixingTestSuite : public mitk::TestFixture +{ + CPPUNIT_TEST_SUITE(mitkSpectralUnmixingTestSuite); + MITK_TEST(testEigenSUAlgorithm); + MITK_TEST(testVigraSUAlgorithm); + //MITK_TEST(testSimplexSUAlgorithm); + MITK_TEST(testSO2); + MITK_TEST(testExceptionSO2); + MITK_TEST(testWavelengthExceptions); + MITK_TEST(testNoChromophoresSelected); + MITK_TEST(testInputImage); + MITK_TEST(testAddOutput); + MITK_TEST(testWeightsError); + MITK_TEST(testOutputs); + CPPUNIT_TEST_SUITE_END(); + +private: + mitk::pa::SpectralUnmixingFilterBase::Pointer m_SpectralUnmixingFilter; + mitk::Image::Pointer inputImage; + std::vector m_inputWavelengths; + std::vector m_inputWeights; + std::vector m_CorrectResult; + float threshold; + +public: + + void setUp() override + { + //Set empty input image: + inputImage = mitk::Image::New(); + mitk::PixelType pixelType = mitk::MakeScalarPixelType(); + const int NUMBER_OF_SPATIAL_DIMENSIONS = 3; + auto* dimensions = new unsigned int[NUMBER_OF_SPATIAL_DIMENSIONS]; + + dimensions[0] = 1; + dimensions[1] = 1; + dimensions[2] = 5; + + inputImage->Initialize(pixelType, NUMBER_OF_SPATIAL_DIMENSIONS, dimensions); + + //Set wavelengths for unmixing: + m_inputWavelengths.push_back(750); + m_inputWavelengths.push_back(800); + + m_inputWeights.push_back(50); + m_inputWeights.push_back(100); + + //Set fraction of Hb and HbO2 to unmix: + float fracHb = 100; + float fracHbO2 = 300; + m_CorrectResult.push_back(fracHbO2); + m_CorrectResult.push_back(fracHb); + m_CorrectResult.push_back(fracHbO2 + 10); + m_CorrectResult.push_back(fracHb - 10); + threshold = 0.01; + + //Multiply values of wavelengths (750,800,850 nm) with fractions to get pixel values: + float px1 = fracHb * 7.52 + fracHbO2 * 2.77; + float px2 = fracHb * 4.08 + fracHbO2 * 4.37; + float px3 = (fracHb - 10) * 7.52 + (fracHbO2 + 10) * 2.77; + float px4 = (fracHb - 10) * 4.08 + (fracHbO2 + 10) * 4.37; + + float* data = new float[6]; + data[0] = px1; + data[1] = px2; + data[2] = px3; + data[3] = px4; + data[5] = 0; + + inputImage->SetImportVolume(data, mitk::Image::ImportMemoryManagementType::CopyMemory); + delete[] data; + } + + // Tests implemented EIGEN algortihms with correct inputs + void testEigenSUAlgorithm() + { + MITK_INFO << "START FILTER TEST ... "; + // Set input image + auto m_SpectralUnmixingFilter = mitk::pa::LinearSpectralUnmixingFilter::New(); + m_SpectralUnmixingFilter->Verbose(false); + m_SpectralUnmixingFilter->RelativeError(false); + m_SpectralUnmixingFilter->SetInput(inputImage); + m_SpectralUnmixingFilter->AddOutputs(2); + + //Set wavelengths to filter + for (unsigned int imageIndex = 0; imageIndex < m_inputWavelengths.size(); imageIndex++) + { + unsigned int wavelength = m_inputWavelengths[imageIndex]; + m_SpectralUnmixingFilter->AddWavelength(wavelength); + } + + //Set Chromophores to filter + m_SpectralUnmixingFilter->AddChromophore( + mitk::pa::PropertyCalculator::ChromophoreType::OXYGENATED); + m_SpectralUnmixingFilter->AddChromophore( + mitk::pa::PropertyCalculator::ChromophoreType::DEOXYGENATED); + + ofstream myfile; + myfile.open("EigenTestResult.txt"); + + std::vector m_Eigen = { + mitk::pa::LinearSpectralUnmixingFilter::AlgortihmType::HOUSEHOLDERQR, /* mitk::pa::LinearSpectralUnmixingFilter::AlgortihmType::LDLT, + mitk::pa::LinearSpectralUnmixingFilter::AlgortihmType::LLT,*/ mitk::pa::LinearSpectralUnmixingFilter::AlgortihmType::COLPIVHOUSEHOLDERQR, + mitk::pa::LinearSpectralUnmixingFilter::AlgortihmType::JACOBISVD, mitk::pa::LinearSpectralUnmixingFilter::AlgortihmType::FULLPIVLU, + mitk::pa::LinearSpectralUnmixingFilter::AlgortihmType::FULLPIVHOUSEHOLDERQR}; + + for (unsigned int Algorithmidx = 0; Algorithmidx < m_Eigen.size();++Algorithmidx) + { + m_SpectralUnmixingFilter->SetAlgorithm(m_Eigen[Algorithmidx]); + + m_SpectralUnmixingFilter->Update(); + + + /*For printed pixel values and results look at: [...]\mitk-superbuild\MITK-build\Modules\PhotoacousticsLib\test\*/ + + for (int i = 0; i < 2; ++i) + { + mitk::Image::Pointer output = m_SpectralUnmixingFilter->GetOutput(i); + mitk::ImageReadAccessor readAccess(output); + const float* inputDataArray = ((const float*)readAccess.GetData()); + auto pixel = inputDataArray[0]; + auto pixel2 = inputDataArray[1]; + + myfile << "Algorithmidx: " << Algorithmidx << "\n Output " << i << ": " << "\n" << inputDataArray[0] << "\n" << inputDataArray[1] << "\n"; + myfile << "Correct Result: " << "\n" << m_CorrectResult[i] << "\n" << m_CorrectResult[i + 2] << "\n"; + + CPPUNIT_ASSERT(std::abs(pixel - m_CorrectResult[i]) < threshold); + CPPUNIT_ASSERT(std::abs(pixel2 - m_CorrectResult[i + 2]) < threshold); + } + } + myfile.close(); + MITK_INFO << "EIGEN FILTER TEST SUCCESFULL :)"; + } + + // Tests implemented VIGRA algortihms with correct inputs + void testVigraSUAlgorithm() + { + MITK_INFO << "START FILTER TEST ... "; + // Set input image + auto m_SpectralUnmixingFilter = mitk::pa::SpectralUnmixingFilterVigra::New(); + m_SpectralUnmixingFilter->SetInput(inputImage); + m_SpectralUnmixingFilter->AddOutputs(2); + m_SpectralUnmixingFilter->Verbose(false); + m_SpectralUnmixingFilter->RelativeError(false); + + + //Set wavelengths to filter + for (unsigned int imageIndex = 0; imageIndex < m_inputWavelengths.size(); imageIndex++) + { + unsigned int wavelength = m_inputWavelengths[imageIndex]; + double Weight = m_inputWeights[imageIndex]; + m_SpectralUnmixingFilter->AddWavelength(wavelength); + m_SpectralUnmixingFilter->AddWeight(Weight); + } + + //Set Chromophores to filter + m_SpectralUnmixingFilter->AddChromophore( + mitk::pa::PropertyCalculator::ChromophoreType::OXYGENATED); + m_SpectralUnmixingFilter->AddChromophore( + mitk::pa::PropertyCalculator::ChromophoreType::DEOXYGENATED); + + + std::vector Vigra = { + mitk::pa::SpectralUnmixingFilterVigra::VigraAlgortihmType::LARS, mitk::pa::SpectralUnmixingFilterVigra::VigraAlgortihmType::GOLDFARB, + mitk::pa::SpectralUnmixingFilterVigra::VigraAlgortihmType::WEIGHTED, mitk::pa::SpectralUnmixingFilterVigra::VigraAlgortihmType::LS/*, + mitk::pa::SpectralUnmixingFilterVigra::VigraAlgortihmType::vigratest*/}; + + ofstream myfile; + myfile.open("VigraTestResult.txt"); + + for (unsigned int Algorithmidx = 0; Algorithmidx < Vigra.size();++Algorithmidx) + { + m_SpectralUnmixingFilter->SetAlgorithm(Vigra[0]); + + m_SpectralUnmixingFilter->Update(); + + + /*For printed pixel values and results look at: [...]\mitk-superbuild\MITK-build\Modules\PhotoacousticsLib\test\*/ + + for (int i = 0; i < 2; ++i) + { + mitk::Image::Pointer output = m_SpectralUnmixingFilter->GetOutput(i); + mitk::ImageReadAccessor readAccess(output); + const float* inputDataArray = ((const float*)readAccess.GetData()); + auto pixel = inputDataArray[0]; + auto pixel2 = inputDataArray[1]; + + myfile << "Algorithmidx: " << Algorithmidx << "\n Output " << i << ": " << "\n" << inputDataArray[0] << "\n" << inputDataArray[1] << "\n"; + myfile << "Correct Result: " << "\n" << m_CorrectResult[i] << "\n" << m_CorrectResult[i + 2] << "\n"; + + CPPUNIT_ASSERT(std::abs(pixel - m_CorrectResult[i]) < threshold); + CPPUNIT_ASSERT(std::abs(pixel2 - m_CorrectResult[i + 2]) < threshold); + } + } + myfile.close(); + MITK_INFO << "VIGRA FILTER TEST SUCCESFULL :)"; + } + + void testSimplexSUAlgorithm() + { + MITK_INFO << "START FILTER TEST ... "; + // Set input image + auto m_SpectralUnmixingFilter = mitk::pa::SpectralUnmixingFilterSimplex::New(); + m_SpectralUnmixingFilter->SetInput(inputImage); + m_SpectralUnmixingFilter->AddOutputs(2); + m_SpectralUnmixingFilter->Verbose(true); + m_SpectralUnmixingFilter->RelativeError(false); + + //Set wavelengths to filter + for (unsigned int imageIndex = 0; imageIndex < m_inputWavelengths.size(); imageIndex++) + { + unsigned int wavelength = m_inputWavelengths[imageIndex]; + m_SpectralUnmixingFilter->AddWavelength(wavelength); + } + + //Set Chromophores to filter + m_SpectralUnmixingFilter->AddChromophore( + mitk::pa::PropertyCalculator::ChromophoreType::OXYGENATED); + m_SpectralUnmixingFilter->AddChromophore( + mitk::pa::PropertyCalculator::ChromophoreType::DEOXYGENATED); + + m_SpectralUnmixingFilter->Update(); + + + /*For printed pixel values and results look at: [...]\mitk-superbuild\MITK-build\Modules\PhotoacousticsLib\test\*/ + ofstream myfile; + myfile.open("SimplexTestResult.txt"); + for (int i = 0; i < 2; ++i) + { + mitk::Image::Pointer output = m_SpectralUnmixingFilter->GetOutput(i); + mitk::ImageReadAccessor readAccess(output); + const float* inputDataArray = ((const float*)readAccess.GetData()); + auto pixel = inputDataArray[0]; + auto pixel2 = inputDataArray[1]; + + myfile << "Output " << i << ": " << "\n" << inputDataArray[0] << "\n" << inputDataArray[1] << "\n"; + myfile << "Correct Result: " << "\n" << m_CorrectResult[i] << "\n" << m_CorrectResult[i + 2] << "\n"; + + CPPUNIT_ASSERT(std::abs(pixel - m_CorrectResult[i])SetInput(0, inputImage); + + inputImage = nullptr; + inputImage = mitk::Image::New(); + mitk::PixelType pixelType = mitk::MakeScalarPixelType(); + const int NUMBER_OF_SPATIAL_DIMENSIONS = 3; + auto* dimensions = new unsigned int[NUMBER_OF_SPATIAL_DIMENSIONS]; + + dimensions[0] = 1; + dimensions[1] = 1; + dimensions[2] = 4; + + inputImage->Initialize(pixelType, NUMBER_OF_SPATIAL_DIMENSIONS, dimensions); + + float* data = new float[4]; + data[0] = 1; + data[1] = 2; + data[2] = 3; + data[3] = 4; + + inputImage->SetImportVolume(data, mitk::Image::ImportMemoryManagementType::CopyMemory); + delete[] data; + + m_sO2->SetInput(1, inputImage); + + MITK_TEST_FOR_EXCEPTION_BEGIN(itk::ExceptionObject) + m_sO2->Update(); + MITK_TEST_FOR_EXCEPTION_END(itk::ExceptionObject) + } + + // Tests SO2 Filter with correct inputs + void testSO2() + { + MITK_INFO << "START SO2 TEST ... "; + std::vector CorrectSO2Result1 = { 0, 0, 0, 0, 0 }; + std::vector Test1 = { 0,0,0,51 }; + std::vector CorrectSO2Result2 = { 0, 0.5, 0, 0.5, 0 }; + std::vector Test2 = { 1584, 0, 0, 0 }; + std::vector CorrectSO2Result3 = { 0.5, 0.5, 0, 0.5, 0 }; + std::vector Test3 = { 0, 1536, 0, 0 }; + std::vector CorrectSO2Result4 = { 0.5, 0.5, 0, 0.5, 0 }; + std::vector Test4 = { 0, 0, 3072, 49 }; + std::vector CorrectSO2Result5 = { 0.5, 0.5, 0.5, 0.5, 0 }; + std::vector Test5 = { 1, 1, 1, 49 }; + + std::vector> TestList; + std::vector> ResultList; + + TestList.push_back(Test1); + TestList.push_back(Test2); + TestList.push_back(Test3); + TestList.push_back(Test4); + TestList.push_back(Test5); + + ResultList.push_back(CorrectSO2Result1); + ResultList.push_back(CorrectSO2Result2); + ResultList.push_back(CorrectSO2Result3); + ResultList.push_back(CorrectSO2Result4); + ResultList.push_back(CorrectSO2Result5); + + /*For printed pixel values and results look at: [...]\mitk-superbuild\MITK-build\Modules\PhotoacousticsLib\test\*/ + ofstream myfile; + myfile.open("SO2TestResult.txt"); + + for (int k = 0; k < 5; ++k) + { + std::vector SO2Settings = TestList[k]; + std::vector m_CorrectSO2Result = ResultList[k]; + auto m_sO2 = mitk::pa::SpectralUnmixingSO2::New(); + m_sO2->SetInput(0, inputImage); + m_sO2->SetInput(1, inputImage); + for (unsigned int i = 0; i < SO2Settings.size(); ++i) + m_sO2->AddSO2Settings(SO2Settings[i]); + m_sO2->Update(); + + mitk::Image::Pointer output = m_sO2->GetOutput(0); + mitk::ImageReadAccessor readAccess(output); + const float* inputDataArray = ((const float*)readAccess.GetData()); + + for (unsigned int Pixel = 0; Pixel < inputImage->GetDimensions()[2]; ++Pixel) + { + auto Value = inputDataArray[Pixel]; + + myfile << "Output(Test " << k << ") " << Pixel << ": " << "\n" << Value << "\n"; + myfile << "Correct Result: " << "\n" << m_CorrectSO2Result[Pixel] << "\n"; + + CPPUNIT_ASSERT(std::abs(Value - m_CorrectSO2Result[Pixel]) < threshold); + } + } + myfile.close(); + MITK_INFO << "SO2 TEST SUCCESFULL :)"; + } + + // Test exceptions for wrong wavelength inputs + void testWavelengthExceptions() + { + MITK_INFO << "START WavelengthExceptions TEST ... "; + // Set input image + auto m_SpectralUnmixingFilter = mitk::pa::LinearSpectralUnmixingFilter::New(); + m_SpectralUnmixingFilter->Verbose(false); + m_SpectralUnmixingFilter->RelativeError(false); + m_SpectralUnmixingFilter->SetInput(inputImage); + m_SpectralUnmixingFilter->AddOutputs(2); + + m_SpectralUnmixingFilter->AddChromophore( + mitk::pa::PropertyCalculator::ChromophoreType::OXYGENATED); + m_SpectralUnmixingFilter->AddChromophore( + mitk::pa::PropertyCalculator::ChromophoreType::DEOXYGENATED); + + m_SpectralUnmixingFilter->SetAlgorithm(mitk::pa::LinearSpectralUnmixingFilter::AlgortihmType::HOUSEHOLDERQR); + + MITK_TEST_FOR_EXCEPTION_BEGIN(itk::ExceptionObject) + m_SpectralUnmixingFilter->Update(); + MITK_TEST_FOR_EXCEPTION_END(itk::ExceptionObject) + + m_SpectralUnmixingFilter->AddWavelength(300); + + MITK_TEST_FOR_EXCEPTION_BEGIN(itk::ExceptionObject) + m_SpectralUnmixingFilter->Update(); + MITK_TEST_FOR_EXCEPTION_END(itk::ExceptionObject) + + m_SpectralUnmixingFilter->AddWavelength(299); + + MITK_TEST_FOR_EXCEPTION_BEGIN(itk::ExceptionObject) + m_SpectralUnmixingFilter->Update(); + MITK_TEST_FOR_EXCEPTION_END(itk::ExceptionObject) + MITK_INFO << "DONE"; + } + + // Test exceptions for wrong chromophore inputs + void testNoChromophoresSelected() + { + MITK_INFO << "testNoChromophoresSelected"; + // Set input image + auto m_SpectralUnmixingFilter = mitk::pa::LinearSpectralUnmixingFilter::New(); + m_SpectralUnmixingFilter->Verbose(false); + m_SpectralUnmixingFilter->RelativeError(false); + m_SpectralUnmixingFilter->SetInput(inputImage); + m_SpectralUnmixingFilter->AddOutputs(2); + + //Set wavelengths to filter + for (unsigned int imageIndex = 0; imageIndex < m_inputWavelengths.size(); imageIndex++) + { + unsigned int wavelength = m_inputWavelengths[imageIndex]; + m_SpectralUnmixingFilter->AddWavelength(wavelength); + } + + m_SpectralUnmixingFilter->SetAlgorithm(mitk::pa::LinearSpectralUnmixingFilter::AlgortihmType::HOUSEHOLDERQR); + + MITK_TEST_FOR_EXCEPTION_BEGIN(itk::ExceptionObject) + m_SpectralUnmixingFilter->Update(); + MITK_TEST_FOR_EXCEPTION_END(itk::ExceptionObject) + } + + // Test exceptions for wrong input image + void testInputImage() + { + MITK_INFO << "INPUT IMAGE TEST"; + inputImage = nullptr; + // Set input image + auto m_SpectralUnmixingFilter = mitk::pa::LinearSpectralUnmixingFilter::New(); + m_SpectralUnmixingFilter->Verbose(false); + m_SpectralUnmixingFilter->RelativeError(false); + //m_SpectralUnmixingFilter->SetInput(inputImage); + m_SpectralUnmixingFilter->AddOutputs(2); + + //Set wavelengths to filter + for (unsigned int imageIndex = 0; imageIndex < m_inputWavelengths.size(); imageIndex++) + { + unsigned int wavelength = m_inputWavelengths[imageIndex]; + m_SpectralUnmixingFilter->AddWavelength(wavelength); + } + + m_SpectralUnmixingFilter->AddChromophore( + mitk::pa::PropertyCalculator::ChromophoreType::OXYGENATED); + m_SpectralUnmixingFilter->AddChromophore( + mitk::pa::PropertyCalculator::ChromophoreType::DEOXYGENATED); + + m_SpectralUnmixingFilter->SetAlgorithm(mitk::pa::LinearSpectralUnmixingFilter::AlgortihmType::HOUSEHOLDERQR); + + MITK_TEST_FOR_EXCEPTION_BEGIN(itk::ExceptionObject) + m_SpectralUnmixingFilter->Update(); + MITK_TEST_FOR_EXCEPTION_END(itk::ExceptionObject) + + inputImage = mitk::Image::New(); + mitk::PixelType pixelType = mitk::MakeScalarPixelType(); + const int NUMBER_OF_SPATIAL_DIMENSIONS = 3; + auto* dimensions = new unsigned int[NUMBER_OF_SPATIAL_DIMENSIONS]; + + dimensions[0] = 1; + dimensions[1] = 1; + dimensions[2] = 5; + + inputImage->Initialize(pixelType, NUMBER_OF_SPATIAL_DIMENSIONS, dimensions); + + double* data = new double[6]; + data[0] = 1; + data[1] = 2; + data[2] = 3; + data[3] = 4; + data[5] = 0; + + inputImage->SetImportVolume(data, mitk::Image::ImportMemoryManagementType::CopyMemory); + delete[] data; + + m_SpectralUnmixingFilter->SetInput(inputImage); + + MITK_TEST_FOR_EXCEPTION_BEGIN(itk::ExceptionObject) + m_SpectralUnmixingFilter->Update(); + MITK_TEST_FOR_EXCEPTION_END(itk::ExceptionObject) + } + + // Test exceptions for addOutputs method + void testAddOutput() + { + MITK_INFO << "addOutputs TEST"; + + // Set input image + auto m_SpectralUnmixingFilter = mitk::pa::LinearSpectralUnmixingFilter::New(); + m_SpectralUnmixingFilter->Verbose(false); + m_SpectralUnmixingFilter->RelativeError(false); + m_SpectralUnmixingFilter->SetInput(inputImage); + + //Set wavelengths to filter + for (unsigned int imageIndex = 0; imageIndex < m_inputWavelengths.size(); imageIndex++) + { + unsigned int wavelength = m_inputWavelengths[imageIndex]; + m_SpectralUnmixingFilter->AddWavelength(wavelength); + } + + m_SpectralUnmixingFilter->AddChromophore( + mitk::pa::PropertyCalculator::ChromophoreType::OXYGENATED); + m_SpectralUnmixingFilter->AddChromophore( + mitk::pa::PropertyCalculator::ChromophoreType::DEOXYGENATED); + + m_SpectralUnmixingFilter->SetAlgorithm(mitk::pa::LinearSpectralUnmixingFilter::AlgortihmType::HOUSEHOLDERQR); + + for (int i = 0; i < 4; ++i) + { + MITK_INFO << "i: " << i; + if (i != 2) + { + MITK_TEST_FOR_EXCEPTION_BEGIN(itk::ExceptionObject) + m_SpectralUnmixingFilter->AddOutputs(i); + m_SpectralUnmixingFilter->Update(); + MITK_TEST_FOR_EXCEPTION_END(itk::ExceptionObject) + } + else + { + m_SpectralUnmixingFilter->AddOutputs(2); + m_SpectralUnmixingFilter->Update(); + for (int i = 0; i < 2; ++i) + { + mitk::Image::Pointer output = m_SpectralUnmixingFilter->GetOutput(i); + mitk::ImageReadAccessor readAccess(output); + const float* inputDataArray = ((const float*)readAccess.GetData()); + auto pixel = inputDataArray[0]; + auto pixel2 = inputDataArray[1]; + + CPPUNIT_ASSERT(std::abs(pixel - m_CorrectResult[i]) < threshold); + CPPUNIT_ASSERT(std::abs(pixel2 - m_CorrectResult[i + 2]) < threshold); + } + } + } + } + + // Test exceptions for weights error + void testWeightsError() + { + MITK_INFO << "testWeightsError"; + + // Set input image + auto m_SpectralUnmixingFilter = mitk::pa::SpectralUnmixingFilterVigra::New(); + m_SpectralUnmixingFilter->Verbose(false); + m_SpectralUnmixingFilter->RelativeError(false); + m_SpectralUnmixingFilter->SetInput(inputImage); + m_SpectralUnmixingFilter->AddOutputs(2); + + //Set wavelengths to filter + for (unsigned int imageIndex = 0; imageIndex < m_inputWavelengths.size(); imageIndex++) + { + unsigned int wavelength = m_inputWavelengths[imageIndex]; + m_SpectralUnmixingFilter->AddWavelength(wavelength); + } + + m_SpectralUnmixingFilter->AddChromophore( + mitk::pa::PropertyCalculator::ChromophoreType::OXYGENATED); + m_SpectralUnmixingFilter->AddChromophore( + mitk::pa::PropertyCalculator::ChromophoreType::DEOXYGENATED); + + m_SpectralUnmixingFilter->SetAlgorithm(mitk::pa::SpectralUnmixingFilterVigra::VigraAlgortihmType::WEIGHTED); + + MITK_TEST_FOR_EXCEPTION_BEGIN(itk::ExceptionObject) + m_SpectralUnmixingFilter->Update(); + MITK_TEST_FOR_EXCEPTION_END(itk::ExceptionObject) + + m_SpectralUnmixingFilter->AddWeight(50); + MITK_TEST_FOR_EXCEPTION_BEGIN(itk::ExceptionObject) + m_SpectralUnmixingFilter->Update(); + MITK_TEST_FOR_EXCEPTION_END(itk::ExceptionObject) + + m_SpectralUnmixingFilter->AddWeight(50); + m_SpectralUnmixingFilter->Update(); + + for (int i = 0; i < 2; ++i) + { + mitk::Image::Pointer output = m_SpectralUnmixingFilter->GetOutput(i); + mitk::ImageReadAccessor readAccess(output); + const float* inputDataArray = ((const float*)readAccess.GetData()); + auto pixel = inputDataArray[0]; + auto pixel2 = inputDataArray[1]; + + CPPUNIT_ASSERT(std::abs(pixel - m_CorrectResult[i]) < threshold); + CPPUNIT_ASSERT(std::abs(pixel2 - m_CorrectResult[i + 2]) < threshold); + } + } + + // Test correct outputs + void testOutputs() + { + MITK_INFO << "TEST"; + + // Set input image + auto m_SpectralUnmixingFilter = mitk::pa::LinearSpectralUnmixingFilter::New(); + m_SpectralUnmixingFilter->Verbose(false); + m_SpectralUnmixingFilter->RelativeError(false); + m_SpectralUnmixingFilter->SetInput(inputImage); + m_SpectralUnmixingFilter->AddOutputs(2); + + //Set wavelengths to filter + for (unsigned int imageIndex = 0; imageIndex < m_inputWavelengths.size(); imageIndex++) + { + unsigned int wavelength = m_inputWavelengths[imageIndex]; + m_SpectralUnmixingFilter->AddWavelength(wavelength); + } + + m_SpectralUnmixingFilter->AddChromophore( + mitk::pa::PropertyCalculator::ChromophoreType::OXYGENATED); + m_SpectralUnmixingFilter->AddChromophore( + mitk::pa::PropertyCalculator::ChromophoreType::DEOXYGENATED); + + m_SpectralUnmixingFilter->SetAlgorithm(mitk::pa::LinearSpectralUnmixingFilter::AlgortihmType::HOUSEHOLDERQR); + + m_SpectralUnmixingFilter->Update(); + + for (int i = 0; i < 2; ++i) + { + mitk::Image::Pointer output = m_SpectralUnmixingFilter->GetOutput(i); + mitk::ImageReadAccessor readAccess(output); + const float* inputDataArray = ((const float*)readAccess.GetData()); + auto pixel = inputDataArray[0]; + auto pixel2 = inputDataArray[1]; + + CPPUNIT_ASSERT(std::abs(pixel - m_CorrectResult[i]) < threshold); + CPPUNIT_ASSERT(std::abs(pixel2 - m_CorrectResult[i + 2]) < threshold); + + // test correct output dimensions and pixel type + CPPUNIT_ASSERT(inputImage->GetDimensions()[0] == output->GetDimensions()[0]); + CPPUNIT_ASSERT(inputImage->GetDimensions()[0] == output->GetDimensions()[1]); + CPPUNIT_ASSERT(2 == output->GetDimensions()[2]); + CPPUNIT_ASSERT(output->GetPixelType() == mitk::MakeScalarPixelType()); + } + } + + // TEST TEMPLATE: + /* + // Test exceptions for + void test() + { + MITK_INFO << "TEST"; + + // Set input image + auto m_SpectralUnmixingFilter = mitk::pa::LinearSpectralUnmixingFilter::New(); + m_SpectralUnmixingFilter->Verbose(false); + m_SpectralUnmixingFilter->RelativeError(false); + m_SpectralUnmixingFilter->SetInput(inputImage); + m_SpectralUnmixingFilter->AddOutputs(2); + + //Set wavelengths to filter + for (unsigned int imageIndex = 0; imageIndex < m_inputWavelengths.size(); imageIndex++) + { + unsigned int wavelength = m_inputWavelengths[imageIndex]; + m_SpectralUnmixingFilter->AddWavelength(wavelength); + } + + m_SpectralUnmixingFilter->AddChromophore( + mitk::pa::PropertyCalculator::ChromophoreType::OXYGENATED); + m_SpectralUnmixingFilter->AddChromophore( + mitk::pa::PropertyCalculator::ChromophoreType::DEOXYGENATED); + + m_SpectralUnmixingFilter->SetAlgorithm(mitk::pa::LinearSpectralUnmixingFilter::AlgortihmType::HOUSEHOLDERQR); + + //MITK_TEST_FOR_EXCEPTION_BEGIN(itk::ExceptionObject) + m_SpectralUnmixingFilter->Update(); + //MITK_TEST_FOR_EXCEPTION_END(itk::ExceptionObject) + + }*/ + + void tearDown() override + { + m_SpectralUnmixingFilter = nullptr; + inputImage = nullptr; + m_inputWavelengths.clear(); + m_CorrectResult.clear(); + } +}; + +MITK_TEST_SUITE_REGISTRATION(mitkSpectralUnmixing) diff --git a/Plugins/PluginList.cmake b/Plugins/PluginList.cmake index 4ade7b895c..6857fe6ef0 100644 --- a/Plugins/PluginList.cmake +++ b/Plugins/PluginList.cmake @@ -1,112 +1,113 @@ # Plug-ins must be ordered according to their dependencies set(MITK_PLUGINS org.blueberry.core.runtime:ON org.blueberry.core.expressions:OFF org.blueberry.core.commands:OFF org.blueberry.core.jobs:OFF org.blueberry.ui.qt:OFF org.blueberry.ui.qt.help:ON org.blueberry.ui.qt.log:ON org.blueberry.ui.qt.objectinspector:OFF #org.blueberry.test:ON #org.blueberry.uitest:ON #Testing/org.blueberry.core.runtime.tests:ON #Testing/org.blueberry.osgi.tests:ON org.mitk.core.services:ON org.mitk.gui.common:ON org.mitk.planarfigure:ON org.mitk.core.ext:OFF org.mitk.core.jobs:OFF org.mitk.gui.qt.application:ON org.mitk.gui.qt.coreapplication:OFF org.mitk.gui.qt.ext:OFF org.mitk.gui.qt.extapplication:OFF org.mitk.gui.qt.common:ON org.mitk.gui.qt.stdmultiwidgeteditor:ON org.mitk.gui.qt.common.legacy:OFF org.mitk.gui.qt.cmdlinemodules:OFF org.mitk.gui.qt.diffusionimagingapp:OFF org.mitk.gui.qt.datamanager:ON org.mitk.gui.qt.datamanagerlight:OFF org.mitk.gui.qt.datastorageviewertest:OFF org.mitk.gui.qt.properties:ON org.mitk.gui.qt.basicimageprocessing:OFF org.mitk.gui.qt.dicom:OFF org.mitk.gui.qt.dicominspector:OFF org.mitk.gui.qt.diffusionimaging:OFF org.mitk.gui.qt.diffusionimaging.connectomics:OFF org.mitk.gui.qt.diffusionimaging.denoising:OFF org.mitk.gui.qt.diffusionimaging.fiberfox:OFF org.mitk.gui.qt.diffusionimaging.fiberprocessing:OFF org.mitk.gui.qt.diffusionimaging.ivim:OFF org.mitk.gui.qt.diffusionimaging.odfpeaks:OFF org.mitk.gui.qt.diffusionimaging.partialvolume:OFF org.mitk.gui.qt.diffusionimaging.preprocessing:OFF org.mitk.gui.qt.diffusionimaging.reconstruction:OFF org.mitk.gui.qt.diffusionimaging.registration:OFF org.mitk.gui.qt.diffusionimaging.tbss:OFF org.mitk.gui.qt.diffusionimaging.tractography:OFF org.mitk.gui.qt.diffusionimaging.python:OFF org.mitk.gui.qt.dosevisualization:OFF org.mitk.gui.qt.geometrytools:OFF org.mitk.gui.qt.igtexamples:OFF org.mitk.gui.qt.igttracking:OFF org.mitk.gui.qt.lasercontrol:OFF org.mitk.gui.qt.openigtlink:OFF org.mitk.gui.qt.imagecropper:OFF org.mitk.gui.qt.imagenavigator:ON org.mitk.gui.qt.viewnavigator:OFF org.mitk.gui.qt.materialeditor:OFF org.mitk.gui.qt.measurementtoolbox:OFF org.mitk.gui.qt.moviemaker:OFF org.mitk.gui.qt.pointsetinteraction:OFF org.mitk.gui.qt.pointsetinteractionmultispectrum:OFF org.mitk.gui.qt.python:OFF org.mitk.gui.qt.remeshing:OFF org.mitk.gui.qt.segmentation:OFF org.mitk.gui.qt.aicpregistration:OFF org.mitk.gui.qt.renderwindowmanager:OFF org.mitk.gui.qt.toftutorial:OFF org.mitk.gui.qt.tofutil:OFF org.mitk.gui.qt.tubegraph:OFF org.mitk.gui.qt.ugvisualization:OFF org.mitk.gui.qt.photoacoustics.pausviewer:OFF org.mitk.gui.qt.photoacoustics.pausmotioncompensation:OFF org.mitk.gui.qt.photoacoustics.imageprocessing:OFF org.mitk.gui.qt.photoacoustics.simulation:OFF + org.mitk.gui.qt.photoacoustics.spectralunmixing:OFF org.mitk.gui.qt.ultrasound:OFF org.mitk.gui.qt.volumevisualization:OFF org.mitk.gui.qt.eventrecorder:OFF org.mitk.gui.qt.xnat:OFF org.mitk.gui.qt.igt.app.echotrack:OFF org.mitk.gui.qt.spectrocamrecorder:OFF org.mitk.gui.qt.classificationsegmentation:OFF org.mitk.gui.qt.overlaymanager:OFF org.mitk.gui.qt.igt.app.hummelprotocolmeasurements:OFF org.mitk.gui.qt.multilabelsegmentation:OFF org.mitk.matchpoint.core.helper:OFF org.mitk.gui.qt.matchpoint.algorithm.browser:OFF org.mitk.gui.qt.matchpoint.algorithm.control:OFF org.mitk.gui.qt.matchpoint.algorithm.batch:OFF org.mitk.gui.qt.matchpoint.mapper:OFF org.mitk.gui.qt.matchpoint.framereg:OFF org.mitk.gui.qt.matchpoint.visualizer:OFF org.mitk.gui.qt.matchpoint.evaluator:OFF org.mitk.gui.qt.matchpoint.manipulator:OFF org.mitk.gui.qt.preprocessing.resampling:OFF org.mitk.gui.qt.cest:OFF org.mitk.gui.qt.fit.demo:OFF org.mitk.gui.qt.fit.inspector:OFF org.mitk.gui.qt.fit.genericfitting:OFF org.mitk.gui.qt.pharmacokinetics.mri:OFF org.mitk.gui.qt.pharmacokinetics.pet:OFF org.mitk.gui.qt.pharmacokinetics.simulation:OFF org.mitk.gui.qt.pharmacokinetics.curvedescriptor:OFF org.mitk.gui.qt.pharmacokinetics.concentration.mri:OFF ) diff --git a/Plugins/org.mitk.gui.qt.photoacoustics.spectralunmixing/CMakeLists.txt b/Plugins/org.mitk.gui.qt.photoacoustics.spectralunmixing/CMakeLists.txt new file mode 100644 index 0000000000..ac2049cb04 --- /dev/null +++ b/Plugins/org.mitk.gui.qt.photoacoustics.spectralunmixing/CMakeLists.txt @@ -0,0 +1,7 @@ +project(org_mitk_gui_qt_photoacoustics_spectralunmixing) + +mitk_create_plugin( + EXPORT_DIRECTIVE SPECTRALUNMIXING_EXPORT + EXPORTED_INCLUDE_SUFFIXES src + MODULE_DEPENDS MitkQtWidgetsExt MitkPhotoacousticsLib +) diff --git a/Plugins/org.mitk.gui.qt.photoacoustics.spectralunmixing/documentation/UserManual/GUI.PNG b/Plugins/org.mitk.gui.qt.photoacoustics.spectralunmixing/documentation/UserManual/GUI.PNG new file mode 100644 index 0000000000..23c2c5cdd0 Binary files /dev/null and b/Plugins/org.mitk.gui.qt.photoacoustics.spectralunmixing/documentation/UserManual/GUI.PNG differ diff --git a/Plugins/org.mitk.gui.qt.photoacoustics.spectralunmixing/documentation/UserManual/Manual.dox b/Plugins/org.mitk.gui.qt.photoacoustics.spectralunmixing/documentation/UserManual/Manual.dox new file mode 100644 index 0000000000..177ecacda1 --- /dev/null +++ b/Plugins/org.mitk.gui.qt.photoacoustics.spectralunmixing/documentation/UserManual/Manual.dox @@ -0,0 +1,37 @@ +/** +\page org_mitk_gui_qt_photoacoustics_spectralunmixing The spectral unmixing (SU) plugin + +\imageMacro{icon.png,"Icon of Spectralunmixing",2.00} + +\table of contents +
    +
  • Introduction +
  • Overview +
  • How to add an additional algorithm +
+ +\section org_mitk_gui_qt_photoacoustics_spectralunmixingIntroduction +The spectral unmixing plugin provides a GUI tool to perform spectral unmixing of multispectral MITK images. It was designed to unmix beamformed photoacoustic images. The outputs are MITK images for every chosen absorber (endmember). Furthermore it is possible to calculate the oxygen saturation of the multispectral input if the endmembers oxy- and deoxyhemoglobin are selected in the GUI as well as an output image that contains the information about the relative error between unmixing result and the input image. +Detailed information about the Plugin, the baseclass and its subclasses can be found in their header files. If you want to call the SU filter from your own class have a look at the “mitkSpectralUnmixingTest.cpp”. There you find information about which functions are callable or have to be called from your class to guarantee the promised functionality of the SU filter. + +\section org_mitk_gui_qt_photoacoustics_spectralunmixingOverview +
    +
  • supports several linear Eigen solvers collected in the “mitkPALinearSpectralUnmixingFilter” +
  • supports several linear Vigra solvers collected in the “mitkPASpectralUnmixingVigraFilter” +
  • calculates oxygen saturation with the “mitkPASpectralUnmixingSO2” +
+ +\section How to add an additional algorithm: +If the algorithm fits in one of the existing classes you can ignore steps 0. – 3. +0. Have a look at the commit rMITK36cfd1731089: implement three empty classes for Simplex, Lagrange and Vigra SU algorithms. Which actually are exactly the first (not all!) steps to implement a new algorithm. +1. Add your future header and cpp file to files.cmake +2. Create a header file which needs at least the methods shown in header.png +\imageMacro{header.png,"empty header for a new SU algorithm",2.00} +3. Create a cpp file which takes an Eigen endmember matrix and an Eigen input vector as inputs and returns an Eigen vector as result. A structure like in the cpp.png is recommended. If your class will consist of more than one algorithm you should have an if/else decision between them with an enum like in the cpp.png otherwise you can directly return your result. +\imageMacro{cpp.png,"example cpp file for a new SU algorithm",2.00} +4. In the Plugin you just have to add another “else if” like in the plugin.png. The string in the else if has to be the same then selectable in the GUI(step 5). +\imageMacro{plugin.png,"changes of Plugin for a new SU algorithm",2.00} +5. To make you algorithm selectable you have to add to the GUI Combobox. Klick at 1. (GUI.png), then at 2. and then name your algorithm 3. like in step 4. +\imageMacro{GUI.png,"changes of GUI for a new SU algorithm",2.00} + +*/ diff --git a/Plugins/org.mitk.gui.qt.photoacoustics.spectralunmixing/documentation/UserManual/cpp.PNG b/Plugins/org.mitk.gui.qt.photoacoustics.spectralunmixing/documentation/UserManual/cpp.PNG new file mode 100644 index 0000000000..b47993b6cc Binary files /dev/null and b/Plugins/org.mitk.gui.qt.photoacoustics.spectralunmixing/documentation/UserManual/cpp.PNG differ diff --git a/Plugins/org.mitk.gui.qt.photoacoustics.spectralunmixing/documentation/UserManual/header.PNG b/Plugins/org.mitk.gui.qt.photoacoustics.spectralunmixing/documentation/UserManual/header.PNG new file mode 100644 index 0000000000..1f02be12c8 Binary files /dev/null and b/Plugins/org.mitk.gui.qt.photoacoustics.spectralunmixing/documentation/UserManual/header.PNG differ diff --git a/Plugins/org.mitk.gui.qt.photoacoustics.spectralunmixing/documentation/UserManual/icon.png b/Plugins/org.mitk.gui.qt.photoacoustics.spectralunmixing/documentation/UserManual/icon.png new file mode 100644 index 0000000000..3384467830 Binary files /dev/null and b/Plugins/org.mitk.gui.qt.photoacoustics.spectralunmixing/documentation/UserManual/icon.png differ diff --git a/Plugins/org.mitk.gui.qt.photoacoustics.spectralunmixing/documentation/UserManual/plugin.PNG b/Plugins/org.mitk.gui.qt.photoacoustics.spectralunmixing/documentation/UserManual/plugin.PNG new file mode 100644 index 0000000000..89d7559d66 Binary files /dev/null and b/Plugins/org.mitk.gui.qt.photoacoustics.spectralunmixing/documentation/UserManual/plugin.PNG differ diff --git a/Plugins/org.mitk.gui.qt.photoacoustics.spectralunmixing/documentation/doxygen/modules.dox b/Plugins/org.mitk.gui.qt.photoacoustics.spectralunmixing/documentation/doxygen/modules.dox new file mode 100644 index 0000000000..de1394cc0e --- /dev/null +++ b/Plugins/org.mitk.gui.qt.photoacoustics.spectralunmixing/documentation/doxygen/modules.dox @@ -0,0 +1,16 @@ +/** + \defgroup org_mitk_gui_qt_photoacoustics_spectralunmixing org.mitk.gui.qt.photoacoustics.spectralunmixing + \ingroup MITKPlugins + + \brief Describe your plugin here. + +*/ + +/** + \defgroup org_mitk_gui_qt_photoacoustics_spectralunmixing_internal Internal + \ingroup org_mitk_gui_qt_photoacoustics_spectralunmixing + + \brief This subcategory includes the internal classes of the org.mitk.gui.qt.photoacoustics.spectralunmixing plugin. Other + plugins must not rely on these classes. They contain implementation details and their interface + may change at any time. We mean it. +*/ diff --git a/Plugins/org.mitk.gui.qt.photoacoustics.spectralunmixing/files.cmake b/Plugins/org.mitk.gui.qt.photoacoustics.spectralunmixing/files.cmake new file mode 100644 index 0000000000..2940e51e6b --- /dev/null +++ b/Plugins/org.mitk.gui.qt.photoacoustics.spectralunmixing/files.cmake @@ -0,0 +1,42 @@ +set(SRC_CPP_FILES + +) + +set(INTERNAL_CPP_FILES + org_mitk_gui_qt_photoacoustics_spectralunmixing_Activator.cpp + SpectralUnmixing.cpp +) + +set(UI_FILES + src/internal/SpectralUnmixingControls.ui +) + +set(MOC_H_FILES + src/internal/org_mitk_gui_qt_photoacoustics_spectralunmixing_Activator.h + src/internal/SpectralUnmixing.h +) + +# list of resource files which can be used by the plug-in +# system without loading the plug-ins shared library, +# for example the icon used in the menu and tabs for the +# plug-in views in the workbench +set(CACHED_RESOURCE_FILES + resources/spectralunmixing.svg + plugin.xml +) + +# list of Qt .qrc files which contain additional resources +# specific to this plugin +set(QRC_FILES + +) + +set(CPP_FILES ) + +foreach(file ${SRC_CPP_FILES}) + set(CPP_FILES ${CPP_FILES} src/${file}) +endforeach(file ${SRC_CPP_FILES}) + +foreach(file ${INTERNAL_CPP_FILES}) + set(CPP_FILES ${CPP_FILES} src/internal/${file}) +endforeach(file ${INTERNAL_CPP_FILES}) diff --git a/Plugins/org.mitk.gui.qt.photoacoustics.spectralunmixing/manifest_headers.cmake b/Plugins/org.mitk.gui.qt.photoacoustics.spectralunmixing/manifest_headers.cmake new file mode 100644 index 0000000000..a537cfdfe3 --- /dev/null +++ b/Plugins/org.mitk.gui.qt.photoacoustics.spectralunmixing/manifest_headers.cmake @@ -0,0 +1,5 @@ +set(Plugin-Name "Spectralunmixing") +set(Plugin-Version "0.1") +set(Plugin-Vendor "DKFZ") +set(Plugin-ContactAddress "") +set(Require-Plugin org.mitk.gui.qt.common) diff --git a/Plugins/org.mitk.gui.qt.photoacoustics.spectralunmixing/plugin.xml b/Plugins/org.mitk.gui.qt.photoacoustics.spectralunmixing/plugin.xml new file mode 100644 index 0000000000..0e3fe22367 --- /dev/null +++ b/Plugins/org.mitk.gui.qt.photoacoustics.spectralunmixing/plugin.xml @@ -0,0 +1,12 @@ + + + + + + + + diff --git a/Plugins/org.mitk.gui.qt.photoacoustics.spectralunmixing/resources/spectralunmixing.svg b/Plugins/org.mitk.gui.qt.photoacoustics.spectralunmixing/resources/spectralunmixing.svg new file mode 100644 index 0000000000..d021376ed0 --- /dev/null +++ b/Plugins/org.mitk.gui.qt.photoacoustics.spectralunmixing/resources/spectralunmixing.svg @@ -0,0 +1,103 @@ + + + + + + + + + + image/svg+xml + + + + + + + + + + + + + + + + + + diff --git a/Plugins/org.mitk.gui.qt.photoacoustics.spectralunmixing/src/internal/SpectralUnmixing.cpp b/Plugins/org.mitk.gui.qt.photoacoustics.spectralunmixing/src/internal/SpectralUnmixing.cpp new file mode 100644 index 0000000000..acdfcdd956 --- /dev/null +++ b/Plugins/org.mitk.gui.qt.photoacoustics.spectralunmixing/src/internal/SpectralUnmixing.cpp @@ -0,0 +1,482 @@ +/*=================================================================== + +The Medical Imaging Interaction Toolkit (MITK) + +Copyright (c) German Cancer Research Center, +Division of Medical and Biological Informatics. +All rights reserved. + +This software is distributed WITHOUT ANY WARRANTY; without +even the implied warranty of MERCHANTABILITY or FITNESS FOR +A PARTICULAR PURPOSE. + +See LICENSE.txt or http://www.mitk.org for details. + +===================================================================*/ + +// Blueberry +#include +#include + +// Qmitk +#include "SpectralUnmixing.h" + +// Qt +#include + +// mitk image +#include + +// Include to perform Spectral Unmixing +#include "mitkPASpectralUnmixingFilterBase.h" +#include "mitkPALinearSpectralUnmixingFilter.h" +#include "mitkPASpectralUnmixingSO2.h" +#include "mitkPASpectralUnmixingFilterVigra.h" +#include "mitkPASpectralUnmixingFilterLagrange.h" +#include "mitkPASpectralUnmixingFilterSimplex.h" + +#include +#include + +const std::string SpectralUnmixing::VIEW_ID = "org.mitk.views.spectralunmixing"; + +void SpectralUnmixing::SetFocus() +{ + m_Controls.buttonPerformImageProcessing->setFocus(); +} + +void SpectralUnmixing::CreateQtPartControl(QWidget *parent) +{ + // create GUI widgets from the Qt Designer's .ui file + m_Controls.setupUi(parent); + connect(m_Controls.buttonPerformImageProcessing, &QPushButton::clicked, this, &SpectralUnmixing::DoImageProcessing); + m_Controls.tableWeight->hide(); + m_Controls.tableSO2->hide(); + m_Controls.tableError->hide(); + connect((QObject*)(m_Controls.QComboBoxAlgorithm), SIGNAL(currentIndexChanged(int)), this, SLOT(EnableGUIWeight())); + connect((QObject*)(m_Controls.checkBoxsO2), SIGNAL(clicked()), this, SLOT(EnableGUISO2())); + connect((QObject*)(m_Controls.checkBoxError), SIGNAL(clicked()), this, SLOT(EnableGUIError())); + this->connect(this, SIGNAL(finishSignal()), this, SLOT(storeOutputs())); + this->connect(this, SIGNAL(crashSignal()), this, SLOT(crashInfo())); + this->connect(this, SIGNAL(enableSignal()), this, SLOT(EnableGUIControl())); +} + +void SpectralUnmixing::EnableGUIControl() +{ + SwitchGUIControl(false); +} + +void SpectralUnmixing::SwitchGUIControl(bool change) +{ + m_Controls.inputtable->setEnabled(change); + m_Controls.checkBoxOx->setEnabled(change); + m_Controls.checkBoxDeOx->setEnabled(change); + m_Controls.checkBoxMelanin->setEnabled(change); + m_Controls.checkBoxAdd->setEnabled(change); + m_Controls.QComboBoxAlgorithm->setEnabled(change); + m_Controls.tableWeight->setEnabled(change); + m_Controls.checkBoxsO2->setEnabled(change); + m_Controls.tableSO2->setEnabled(change); + m_Controls.checkBoxVerbose->setEnabled(change); + m_Controls.checkBoxChrono->setEnabled(change); + m_Controls.buttonPerformImageProcessing->setEnabled(change); + m_Controls.checkBoxError->setEnabled(change); +} + +void SpectralUnmixing::EnableGUIWeight() +{ + auto qs = m_Controls.QComboBoxAlgorithm->currentText(); + std::string Algorithm = qs.toUtf8().constData(); + if (Algorithm == "weighted") + m_Controls.tableWeight->show(); + else + m_Controls.tableWeight->hide(); +} + +void SpectralUnmixing::EnableGUISO2() +{ + if (m_Controls.checkBoxsO2->isChecked()) + m_Controls.tableSO2->show(); + else + m_Controls.tableSO2->hide(); +} + +void SpectralUnmixing::EnableGUIError() +{ + if (m_Controls.checkBoxError->isChecked()) + m_Controls.tableError->show(); + else + m_Controls.tableError->hide(); +} + +void SpectralUnmixing::SetVerboseMode(mitk::pa::SpectralUnmixingFilterBase::Pointer m_SpectralUnmixingFilter, bool PluginVerbose) +{ + m_SpectralUnmixingFilter->Verbose(PluginVerbose); +} + +void SpectralUnmixing::SetWavlength(mitk::pa::SpectralUnmixingFilterBase::Pointer m_SpectralUnmixingFilter) +{ + int col = 0; + int Wavelength = 1; + while (m_Controls.inputtable->item(0, col) && Wavelength > 0) + { + QString Text = m_Controls.inputtable->item(0, col)->text(); + Wavelength = Text.toInt(); + if (Wavelength > 0) + { + m_SpectralUnmixingFilter->AddWavelength(Wavelength); + MITK_INFO(PluginVerbose) << "Wavelength: " << Wavelength << "nm \n"; + } + ++col; + } +} + +void SpectralUnmixing::SetChromophore(mitk::pa::SpectralUnmixingFilterBase::Pointer m_SpectralUnmixingFilter, std::vector boolVec, + std::vector chromophoreNameVec) +{ + unsigned int numberofChromophores = 0; + + std::vector m_ChromoType = { mitk::pa::PropertyCalculator::ChromophoreType::OXYGENATED, + mitk::pa::PropertyCalculator::ChromophoreType::DEOXYGENATED, mitk::pa::PropertyCalculator::ChromophoreType::MELANIN, + mitk::pa::PropertyCalculator::ChromophoreType::ONEENDMEMBER}; + + for (unsigned int chromo = 0; chromo < m_ChromoType.size(); ++chromo) + { + if (boolVec[chromo] == true) + { + MITK_INFO(PluginVerbose) << "Chromophore: " << chromophoreNameVec[chromo]; + m_SpectralUnmixingFilter->AddChromophore(m_ChromoType[chromo]); + numberofChromophores += 1; + } + } + + if (numberofChromophores == 0) + mitkThrow() << "PRESS 'IGNORE' AND CHOOSE A CHROMOPHORE!"; +} + +void SpectralUnmixing::OnSelectionChanged(berry::IWorkbenchPart::Pointer /*source*/, + const QList &nodes) +{ + // iterate all selected objects, adjust warning visibility + foreach (mitk::DataNode::Pointer node, nodes) + { + if (node.IsNotNull() && dynamic_cast(node->GetData())) + { + m_Controls.labelWarning->setVisible(false); + m_Controls.buttonPerformImageProcessing->setEnabled(true); + return; + } + } + + m_Controls.labelWarning->setVisible(true); + m_Controls.buttonPerformImageProcessing->setEnabled(false); +} + +mitk::pa::SpectralUnmixingFilterBase::Pointer SpectralUnmixing::GetFilterInstance(std::string algorithm) +{ + mitk::pa::SpectralUnmixingFilterBase::Pointer spectralUnmixingFilter; + + if (algorithm == "householderQr") + { + spectralUnmixingFilter = mitk::pa::LinearSpectralUnmixingFilter::New(); + dynamic_cast(spectralUnmixingFilter.GetPointer()) + ->SetAlgorithm(mitk::pa::LinearSpectralUnmixingFilter::AlgortihmType::HOUSEHOLDERQR); + } + + else if (algorithm == "ldlt") + { + MITK_WARN << "Unfortunaly algorithm is likley to fail!"; + spectralUnmixingFilter = mitk::pa::LinearSpectralUnmixingFilter::New(); + dynamic_cast(spectralUnmixingFilter.GetPointer()) + ->SetAlgorithm(mitk::pa::LinearSpectralUnmixingFilter::AlgortihmType::LDLT); + } + + else if (algorithm == "llt") + { + MITK_WARN << "Unfortunaly algorithm is likley to fail!"; + spectralUnmixingFilter = mitk::pa::LinearSpectralUnmixingFilter::New(); + dynamic_cast(spectralUnmixingFilter.GetPointer()) + ->SetAlgorithm(mitk::pa::LinearSpectralUnmixingFilter::AlgortihmType::LLT); + } + + else if (algorithm == "colPivHouseholderQr") + { + spectralUnmixingFilter = mitk::pa::LinearSpectralUnmixingFilter::New(); + dynamic_cast(spectralUnmixingFilter.GetPointer()) + ->SetAlgorithm(mitk::pa::LinearSpectralUnmixingFilter::AlgortihmType::COLPIVHOUSEHOLDERQR); + } + + else if (algorithm == "jacobiSvd") + { + spectralUnmixingFilter = mitk::pa::LinearSpectralUnmixingFilter::New(); + dynamic_cast(spectralUnmixingFilter.GetPointer()) + ->SetAlgorithm(mitk::pa::LinearSpectralUnmixingFilter::AlgortihmType::JACOBISVD); + } + + else if (algorithm == "fullPivLu") + { + spectralUnmixingFilter = mitk::pa::LinearSpectralUnmixingFilter::New(); + dynamic_cast(spectralUnmixingFilter.GetPointer()) + ->SetAlgorithm(mitk::pa::LinearSpectralUnmixingFilter::AlgortihmType::FULLPIVLU); + } + + else if (algorithm == "fullPivHouseholderQr") + { + spectralUnmixingFilter = mitk::pa::LinearSpectralUnmixingFilter::New(); + dynamic_cast(spectralUnmixingFilter.GetPointer()) + ->SetAlgorithm(mitk::pa::LinearSpectralUnmixingFilter::AlgortihmType::FULLPIVHOUSEHOLDERQR); + } + + else if (algorithm == "NNLARS") + { + spectralUnmixingFilter = mitk::pa::SpectralUnmixingFilterVigra::New(); + dynamic_cast(spectralUnmixingFilter.GetPointer()) + ->SetAlgorithm(mitk::pa::SpectralUnmixingFilterVigra::VigraAlgortihmType::LARS); + } + + else if (algorithm == "NNGoldfarb") + { + spectralUnmixingFilter = mitk::pa::SpectralUnmixingFilterVigra::New(); + dynamic_cast(spectralUnmixingFilter.GetPointer()) + ->SetAlgorithm(mitk::pa::SpectralUnmixingFilterVigra::VigraAlgortihmType::GOLDFARB); + } + + else if (algorithm == "weighted") + { + spectralUnmixingFilter = mitk::pa::SpectralUnmixingFilterVigra::New(); + dynamic_cast(spectralUnmixingFilter.GetPointer()) + ->SetAlgorithm(mitk::pa::SpectralUnmixingFilterVigra::VigraAlgortihmType::WEIGHTED); + //Tranfer GUI information(Weights) to filter + unsigned int colunm = 0; + int Weight = 1; + while (m_Controls.tableWeight->item(0, colunm) && Weight > 0) + { + QString Text = m_Controls.tableWeight->item(0, colunm)->text(); + Weight = Text.toInt(); + if (Weight > 0) + { + MITK_INFO(PluginVerbose) << "Weight: " << Weight; + dynamic_cast(spectralUnmixingFilter.GetPointer()) + ->AddWeight(Weight); + } + ++colunm; + } + } + + else if (algorithm == "LS") + { + spectralUnmixingFilter = mitk::pa::SpectralUnmixingFilterVigra::New(); + dynamic_cast(spectralUnmixingFilter.GetPointer()) + ->SetAlgorithm(mitk::pa::SpectralUnmixingFilterVigra::VigraAlgortihmType::LS); + } + + else if (algorithm == "SimplexMax") + { + spectralUnmixingFilter = mitk::pa::SpectralUnmixingFilterSimplex::New(); + } + + else + mitkThrow() << "404 ALGORITHM NOT FOUND!"; + + return spectralUnmixingFilter; +} + + +void SpectralUnmixing::SetSO2Settings(mitk::pa::SpectralUnmixingSO2::Pointer m_sO2) +{ + for (unsigned int i = 0; i < 4; ++i) + { + if (m_Controls.tableSO2->item(0, i)) + { + QString Text = m_Controls.tableSO2->item(0, i)->text(); + int value = Text.toInt(); + MITK_INFO(PluginVerbose) << "SO2 setting value: " << value; + m_sO2->AddSO2Settings(value); + } + else + m_sO2->AddSO2Settings(0); + } +} + +void SpectralUnmixing::SetRelativeErrorSettings(mitk::pa::SpectralUnmixingFilterBase::Pointer m_SpectralUnmixingFilter) +{ + for (unsigned int i = 0; i < 2; ++i) + { + if (m_Controls.tableError->item(0, i)) + { + QString Text = m_Controls.tableError->item(0, i)->text(); + int value = Text.toInt(); + MITK_INFO(PluginVerbose) << "Relative error setting value: " << value; + m_SpectralUnmixingFilter->AddRelativeErrorSettings(value); + } + else + m_SpectralUnmixingFilter->AddRelativeErrorSettings(0); + } +} + +void SpectralUnmixing::CalculateSO2(mitk::pa::SpectralUnmixingFilterBase::Pointer m_SpectralUnmixingFilter, std::vector boolVec) +{ + MITK_INFO(PluginVerbose) << "CALCULATE OXYGEN SATURATION ..."; + + if (!boolVec[0]) + mitkThrow() << "SELECT CHROMOPHORE DEOXYHEMOGLOBIN!"; + if (!boolVec[1]) + mitkThrow() << "SELECT CHROMOPHORE OXYHEMOGLOBIN!"; + auto m_sO2 = mitk::pa::SpectralUnmixingSO2::New(); + m_sO2->Verbose(PluginVerbose); + SetSO2Settings(m_sO2); + + // Initialize pipeline from SU filter class to SO2 class + auto output1 = m_SpectralUnmixingFilter->GetOutput(0); + auto output2 = m_SpectralUnmixingFilter->GetOutput(1); + m_sO2->SetInput(0, output1); + m_sO2->SetInput(1, output2); + + m_sO2->Update(); + + mitk::Image::Pointer sO2 = m_sO2->GetOutput(0); + sO2->SetSpacing(output1->GetGeometry()->GetSpacing()); + WriteOutputToDataStorage(sO2, "sO2"); + MITK_INFO(PluginVerbose) << "[DONE]"; +} + +void SpectralUnmixing::WriteOutputToDataStorage(mitk::Image::Pointer m_Image, std::string name) +{ + mitk::DataNode::Pointer dataNodeOutput = mitk::DataNode::New(); + dataNodeOutput->SetData(m_Image); + dataNodeOutput->SetName(name); + this->GetDataStorage()->Add(dataNodeOutput); +} + +void SpectralUnmixing::Settings(mitk::Image::Pointer image) +{ + boolVec = { m_Controls.checkBoxOx->isChecked(), m_Controls.checkBoxDeOx->isChecked(), + m_Controls.checkBoxMelanin->isChecked(), m_Controls.checkBoxAdd->isChecked() }; + outputNameVec = { "HbO2", "Hb", "Melanin", "Static Endmember" }; + sO2Bool = (m_Controls.checkBoxsO2->isChecked()); + + //Read GUI information(algorithm) + auto qs = m_Controls.QComboBoxAlgorithm->currentText(); + Algorithm = qs.toUtf8().constData(); + + m_SpectralUnmixingFilter = GetFilterInstance(Algorithm); + SetVerboseMode(m_SpectralUnmixingFilter, PluginVerbose); + m_SpectralUnmixingFilter->RelativeError(m_Controls.checkBoxError->isChecked()); + m_SpectralUnmixingFilter->SetInput(image); + SetWavlength(m_SpectralUnmixingFilter); + SetChromophore(m_SpectralUnmixingFilter, boolVec, outputNameVec); + + boolVec.push_back(m_Controls.checkBoxError->isChecked()); + outputNameVec.push_back("Relative Error"); + if (m_Controls.checkBoxError->isChecked()) + SetRelativeErrorSettings(m_SpectralUnmixingFilter); + + m_SpectralUnmixingFilter->AddOutputs(std::accumulate(boolVec.begin(), boolVec.end(), 0)); + MITK_INFO(PluginVerbose) << "Number of indexed outputs: " << std::accumulate(boolVec.begin(), boolVec.end(), 0); +} + +void SpectralUnmixing::storeOutputs() +{ + int outputCounter = 0; + mitk::Image::Pointer m_Output; + for (unsigned int chromophore = 0; chromophore < outputNameVec.size(); ++chromophore) + { + if (boolVec[chromophore] != false) + { + m_Output = m_SpectralUnmixingFilter->GetOutput(outputCounter++); + m_Output->SetSpacing(image->GetGeometry()->GetSpacing()); + WriteOutputToDataStorage(m_Output, outputNameVec[chromophore] + Algorithm); + } + } + + if (sO2Bool) + CalculateSO2(m_SpectralUnmixingFilter, boolVec); + + mitk::RenderingManager::GetInstance()->InitializeViewsByBoundingObjects(this->GetDataStorage()); + MITK_INFO(PluginVerbose) << "Adding images to DataStorage...[DONE]"; + + std::chrono::steady_clock::time_point _end(std::chrono::steady_clock::now()); + MITK_INFO(m_Controls.checkBoxChrono->isChecked()) << "Time for image Processing: " + << std::chrono::duration_cast>(_end - _start).count(); + QApplication::setOverrideCursor(Qt::ArrowCursor); + SwitchGUIControl(true); +} + +void SpectralUnmixing::WorkingThreadUpdateFilter(mitk::pa::SpectralUnmixingFilterBase::Pointer m_SpectralUnmixingFilter) +{ + //SwitchGUIControl(false); + emit enableSignal(); + try + { + m_SpectralUnmixingFilter->Update(); + emit finishSignal(); + } + catch (const mitk::Exception& e) + { + QApplication::setOverrideCursor(Qt::ArrowCursor); + SwitchGUIControl(true); + errorMessage = e.GetDescription(); + emit crashSignal(); + } +} + +void SpectralUnmixing::crashInfo() +{ + const char *error = errorMessage.c_str(); + QMessageBox::information(nullptr, "Template", error); +} + +void SpectralUnmixing::DoImageProcessing() +{ + QList nodes = this->GetDataManagerSelection(); + if (nodes.empty()) + return; + + mitk::DataNode *node = nodes.front(); + + if (!node) + { + // Nothing selected. Inform the user and return + QMessageBox::information(nullptr, "Template", "Please load and select an image before starting image processing."); + return; + } + + // here we have a valid mitk::DataNode + + // a node itself is not very useful, we need its data item (the image) + mitk::BaseData *data = node->GetData(); + if (data) + { + // test if this data item is an image or not (could also be a surface or something totally different) + image = dynamic_cast(data); + if (image) + { + std::stringstream message; + std::string name; + message << "PERFORMING SPECTRAL UNMIXING "; + if (node->GetName(name)) + { + // a property called "name" was found for this DataNode + message << "'" << name << "'"; + } + message << "."; + + _start = std::chrono::steady_clock::now(); + PluginVerbose = m_Controls.checkBoxVerbose->isChecked(); + MITK_INFO(PluginVerbose) << message.str(); + + try + { + Settings(image); + MITK_INFO(PluginVerbose) << "Updating Filter..."; + QApplication::setOverrideCursor(Qt::WaitCursor); + QtConcurrent::run(this, &SpectralUnmixing::WorkingThreadUpdateFilter, m_SpectralUnmixingFilter); + } + catch (const mitk::Exception& e) + { + QApplication::setOverrideCursor(Qt::ArrowCursor); + QMessageBox::information(nullptr, "Template", e.GetDescription()); + } + } + } +} diff --git a/Plugins/org.mitk.gui.qt.photoacoustics.spectralunmixing/src/internal/SpectralUnmixing.h b/Plugins/org.mitk.gui.qt.photoacoustics.spectralunmixing/src/internal/SpectralUnmixing.h new file mode 100644 index 0000000000..175d3e6bdf --- /dev/null +++ b/Plugins/org.mitk.gui.qt.photoacoustics.spectralunmixing/src/internal/SpectralUnmixing.h @@ -0,0 +1,200 @@ +/*=================================================================== + +The Medical Imaging Interaction Toolkit (MITK) + +Copyright (c) German Cancer Research Center, +Division of Medical and Biological Informatics. +All rights reserved. + +This software is distributed WITHOUT ANY WARRANTY; without +even the implied warranty of MERCHANTABILITY or FITNESS FOR +A PARTICULAR PURPOSE. + +See LICENSE.txt or http://www.mitk.org for details. + +===================================================================*/ + +#ifndef SpectralUnmixing_h +#define SpectralUnmixing_h + +#include + +#include +#include +#include + +#include "ui_SpectralUnmixingControls.h" + +/** +* \brief The spectral unmixing plugin provides a GUI tool to perform spectral unmixing of multispectral MITK images. +* It was designed to unmix beamformed photoacoustic imgaes. The outputs are as well MITK images for every chosen absorber +* (endmember). Furthermore it is possible to calculate the oxygen saturation of the multispectral input if the endmembers +* oxy- and deoxyhemoglobin are selected in the GUI. +* +* For further information look at the documentation of the mitkPASpectralUnmixingFilterBase.h +* +* @exeption if the GenerateOutput method throws a exception the plugin will show a QMessageBox with the exception +* message at the GUI +*/ +class SpectralUnmixing : public QmitkAbstractView +{ + // this is needed for all Qt objects that should have a Qt meta-object + // (everything that derives from QObject and wants to have signal/slots) + Q_OBJECT + + public: + static const std::string VIEW_ID; + + protected: + virtual void CreateQtPartControl(QWidget *parent) override; + virtual void SetFocus() override; + + /// \brief called by QmitkFunctionality when DataManager's selection has changed + virtual void OnSelectionChanged(berry::IWorkbenchPart::Pointer source, + const QList &nodes) override; + + + /** + * \brief Called when the user clicks the GUI button. Checks if the selected data is an image. Then passen on the GUI + * information using the Settings method. Afterwards it performs spectral unmixing via the WorkingThreadUpdateFilter + * method in a own thread. The spectral unmixing is based on the spectral unmixing filter base and its subclasses. + * @exception if nothing is selected. Inform the user and return + * @exception if settings fails. Informs with the mitkthorw information of the filter as QMessageBox + */ + void DoImageProcessing(); + + /** + * \brief slots are there to show/hide input tables for weights-, relative error and SO2 settings ig they are not needed + */ + public slots: + void EnableGUIWeight(); + void EnableGUISO2(); + void EnableGUIError(); + + /** + * \brief slot waits for finishSignal of the working thread and starts storeOutputs + */ + public slots: + /** + * \brief slot does the image post processing + * - GetOutput from m_SpectralUnmixingFilter + * - calles WriteOutputToDataStorage + * - if (true) calls CalculateSO2 + * - does the rendering + * - if (true) shows the chrono result + * - switches the GUI back on + */ + void storeOutputs(); + signals: + void finishSignal(); + + public slots: + void EnableGUIControl(); + signals: + void enableSignal(); + + /** + * \brief slot waits for crashSignal and if neccessary ends working thread and shows QMessageBox with the error message + */ + public slots: + void crashInfo(); + signals: + void crashSignal(); + + protected: + Ui::SpectralUnmixingControls m_Controls; + + /** + * \brief passes the algorithm information from the GUI on to the spectral unmixing filter base subclass method + * "SetAlgortihm" and initializes the subclassFilter::Pointer. + * @param algorithm has to be a string which can be assigned to the mitk::pa::LinearSpectralUnmixingFilter::AlgortihmType + * @throws if the algorithm string doesn't match to an implemented algorithm + */ + mitk::pa::SpectralUnmixingFilterBase::Pointer GetFilterInstance(std::string algorithm); + + bool PluginVerbose = true; + mitk::pa::SpectralUnmixingFilterBase::Pointer m_SpectralUnmixingFilter; + std::vector outputNameVec; + std::vector boolVec; + std::string Algorithm; + bool sO2Bool; + mitk::Image *image; + std::chrono::steady_clock::time_point _start; + std::string errorMessage; + + private: + /** + * \brief thread + * - disables GUI + * - tries Filter->Update() method + * - gives finishSignal which calls storeOutputs + * - cathes by enables GUI and gives crashSignal + */ + void WorkingThreadUpdateFilter(mitk::pa::SpectralUnmixingFilterBase::Pointer m_SpectralUnmixingFilter); + + /** + * \brief takes an MITK image as input and performs spectral unmixing based on the spectral unmixing filter base and its subclasses. + * Therefor it first passes all information from the GUI into the filter by using the "set"methods of the plugin, which then are calling + * the "add" methods of the filter(base). + * @param image has to be an MITK image (pointer). For the request on the image look at the docu of the mitkPASpectralUnmixngFilterBase.h + */ + virtual void Settings(mitk::Image::Pointer image); + + /** + * \brief The method takes a image pointer and a file name which then will get to the data storage. + * @param m_Image is a mitk_::Image::Pointer pointing at the output which one wants to get stored + * @param name has to be a string and will be the file name + */ + virtual void WriteOutputToDataStorage(mitk::Image::Pointer m_Image, std::string name); + + /** + * \brief passes the algorithm information if verbose mode is requested from the GUI on to the spectral unmixing filter + * @param m_SpectralUnmixingFilter is a pointer of the spectral unmixing filter base + * @param PluginVerbose is the GUI information bool + */ + virtual void SetVerboseMode(mitk::pa::SpectralUnmixingFilterBase::Pointer m_SpectralUnmixingFilter, bool PluginVerbose); + + /** + * \brief passes the wavelength information from the GUI on to the spectral unmixing filter base method "AddWavelength". + * @param m_SpectralUnmixingFilter is a pointer of the spectral unmixing filter base + */ + virtual void SetWavlength(mitk::pa::SpectralUnmixingFilterBase::Pointer m_SpectralUnmixingFilter); + + /** + * \brief passes the chromophore information from the GUI on to the spectral unmixing filter base method "AddChromophore". + * @param m_SpectralUnmixingFilter is a pointer of the spectral unmixing filter base + * @param boolVec is a vector which contains the information which chromophore was checked in the GUI + * @param chromophoreNameVec contains the names of all chromophores as strings + * @throws "PRESS 'IGNORE' AND CHOOSE A CHROMOPHORE!" if no chromophore was chosen + */ + virtual void SetChromophore(mitk::pa::SpectralUnmixingFilterBase::Pointer m_SpectralUnmixingFilter, std::vector boolVec, std::vector chromophoreNameVec); + + /** + * \brief passes the SetRelativeErrorSettings information from the GUI on to the spectral unmixing filter base method "AddRelativeErrorSettings". + * @param m_SpectralUnmixingFilter is a pointer of the spectral unmixing filter base# + */ + virtual void SetRelativeErrorSettings(mitk::pa::SpectralUnmixingFilterBase::Pointer m_SpectralUnmixingFilter); + + /** + * \brief passes the SetSO2Settings information from the GUI on to the spectral unmixing SO2 filter method "AddSO2Settings". + * @param m_sO2 is a pointer of the spectral unmixing SO2 filter + */ + virtual void SetSO2Settings(mitk::pa::SpectralUnmixingSO2::Pointer m_sO2); + + /** + * \brief calcultes out of two identical sized MITK images the oxygen saturation and stores the result in an image. Herein the two + * input images are the output for oxy- and deoxyhemoglobin from the GenerateOutput method (spectral unmixing filter results). + * @param m_SpectralUnmixingFilter is a pointer of the spectral unmixing filter base to get the filter output images as sO2 input + * @param boolVec is a vector which contains the information which chromophore was checked in the GUI + * @throws if oxy- or deoxyhemoglobin was not selected in the GUI + */ + virtual void CalculateSO2(mitk::pa::SpectralUnmixingFilterBase::Pointer m_SpectralUnmixingFilter, std::vector boolVec); + + /** + * \brief enables/disables GUI + * @param change true means GUI buttons enabled, false disabled respectively + */ + virtual void SwitchGUIControl(bool change); +}; + +#endif // SpectralUnmixing_h diff --git a/Plugins/org.mitk.gui.qt.photoacoustics.spectralunmixing/src/internal/SpectralUnmixingControls.ui b/Plugins/org.mitk.gui.qt.photoacoustics.spectralunmixing/src/internal/SpectralUnmixingControls.ui new file mode 100644 index 0000000000..16095add30 --- /dev/null +++ b/Plugins/org.mitk.gui.qt.photoacoustics.spectralunmixing/src/internal/SpectralUnmixingControls.ui @@ -0,0 +1,1232 @@ + + + SpectralUnmixingControls + + + + 0 + 0 + 325 + 742 + + + + + 0 + 0 + + + + QmitkTemplate + + + + + + QLabel { color: rgb(255, 0, 0) } + + + Please select an image! + + + + + + + Qt::Horizontal + + + + + + + + 75 + true + + + + Wavelengths settings + + + + + + + + 0 + 75 + + + + + 16777215 + 75 + + + + + 10 + + + + <html><head/><body><p>* \brief All values have to be intergers. They need to have the same order than the wavelength at the image.</p></body></html> + + + true + + + Qt::SolidLine + + + 42 + + + 50 + + + 30 + + + 30 + + + + λ [nm] + + + + + 1 + + + + + 2 + + + + + 3 + + + + + 4 + + + + + 5 + + + + + 6 + + + + + 7 + + + + + 8 + + + + + 9 + + + + + 10 + + + + + 11 + + + + + 12 + + + + + 13 + + + + + 14 + + + + + 15 + + + + + 16 + + + + + 17 + + + + + 18 + + + + + 19 + + + + + 20 + + + + + 21 + + + + + 22 + + + + + 23 + + + + + 24 + + + + + 25 + + + + + 26 + + + + + 27 + + + + + 28 + + + + + 29 + + + + + 30 + + + + + 31 + + + + + 32 + + + + + 33 + + + + + 34 + + + + + 35 + + + + + 36 + + + + + 37 + + + + + 38 + + + + + 39 + + + + + 40 + + + + + 700 + + + + + 710 + + + + + 720 + + + + + 730 + + + + + 740 + + + + + 750 + + + + + 760 + + + + + 770 + + + + + 780 + + + + + 790 + + + + + 800 + + + + + 810 + + + + + 820 + + + + + 830 + + + + + 840 + + + + + 850 + + + + + 860 + + + + + 870 + + + + + 880 + + + + + 890 + + + + + 900 + + + + + 910 + + + + + 920 + + + + + 930 + + + + + 940 + + + + + 950 + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + Qt::Horizontal + + + + + + + + 75 + true + + + + <html><head/><body><p><span style=" font-weight:400;">* \brief Select as least one Absorber. It's not possible to select more absorbers then wavelengths.</span></p></body></html> + + + Chromophore selection + + + + + + + <html><head/><body><p>* \brief One of the main absorbers in near infrared spectrum.</p></body></html> + + + Oxyhemoglobin + + + true + + + + + + + true + + + <html><head/><body><p>* \brief One of the main absorbers in near infrared spectrum.</p></body></html> + + + Deoxygenated hemoglobin + + + true + + + true + + + false + + + + + + + <html><head/><body><p>* \brief One of the main absorbers in near infrared spectrum.</p></body></html> + + + Melanin + + + + + + + <html><head/><body><p>* \brief This endmember will be unmixed with 1 at all wavelgnths.</p></body></html> + + + Static Endmember + + + + + + + Qt::Horizontal + + + + + + + + 75 + true + + + + <html><head/><body><p><span style=" font-weight:400;">* \brief One needs to choose an spectral unmixing algorithm.</span></p></body></html> + + + Unmixing algorithm + + + + + + + true + + + + MS Shell Dlg 2 + + + + <html><head/><body><p>* \brief For detailed information about the algorithms please have a look at the documentation.</p></body></html> + + + false + + + 21 + + + 2147483647 + + + + ==CHOSE ALGORITHM== + + + + + ==QR decomposition== + + + + + householderQr + + + + + colPivHouseholderQr + + + + + fullPivHouseholderQr + + + + + ==LU decompositon + + + + + fullPivLu + + + + + ==Cholesky decompostion== + + + + + ldlt + + + + + llt + + + + .. + + + + + ==Least squares== + + + + + LS + + + + + jacobiSvd + + + + + NNLARS + + + + + NNGoldfarb + + + + + weighted + + + + + ==Others== + + + + + SimplexMax + + + + + + + + true + + + + 0 + 75 + + + + + 380 + 81 + + + + <html><head/><body><p>* \brief the weights are at the same order as the wavelength</p></body></html> + + + false + + + false + + + + Weights [%] + + + + + 1 + + + + + 2 + + + + + 3 + + + + + 4 + + + + + 5 + + + + + 6 + + + + + 7 + + + + + 8 + + + + + 9 + + + + + 10 + + + + + 11 + + + + + 12 + + + + + 13 + + + + + 14 + + + + + 15 + + + + + 16 + + + + + 17 + + + + + 18 + + + + + 19 + + + + + 20 + + + + + 21 + + + + + 22 + + + + + 23 + + + + + 24 + + + + + 25 + + + + + 26 + + + + + 11 + + + + + 10 + + + + + 9 + + + + + 8 + + + + + 7 + + + + + 7 + + + + + 6 + + + + + 6 + + + + + 5 + + + + + 5 + + + + + 5 + + + + + 5 + + + + + 5 + + + + + 5 + + + + + 5 + + + + + 5 + + + + + 5 + + + + + 5 + + + + + 5 + + + + + 5 + + + + + 5 + + + + + 5 + + + + + 5 + + + + + 6 + + + + + 6 + + + + + 6 + + + + + + + + Qt::Horizontal + + + + + + + + 75 + true + + + + Oxygen saturation + + + + + + + <html><head/><body><p>* \brief calculates HbO2/(HbO2+Hb) if De- and oxyhemoglobin are selected.</p></body></html> + + + calculate sO2 + + + false + + + + + + + + 380 + 82 + + + + <html><head/><body><p>* \brief below threshold calculated sO2 value will set to zero</p></body></html> + + + + Threshold + + + + + HbO2 + + + + + Hb + + + + + + + + Sum + + + + + SO2 % + + + + + + + + + + + + + + + + + + + + + + + + + + + + Qt::Horizontal + + + + + + + + 75 + true + + + + Additional Settings + + + + + + + <html><head/><body><p>* \brief This mode will give additional console outputs for debugging.</p></body></html> + + + Verbose Mode (additional console outputs) + + + false + + + + + + + <html><head/><body><p>* \brief This checkbox will start Chrono and takes the time between clicking of the &quot;Perform spectral unmixing&quot; button until the GUI enables again.</p></body></html> + + + Chrono + + + + + + + <html><head/><body><p>* \brief Calculates the realtive error between unmixing result and the input image in the L2 norm.</p></body></html> + + + Relative error image + + + + + + + + 307 + 70 + + + + <html><head/><body><p>* \brief below the threshold calculated relative error will set to zero</p></body></html> + + + + Threshold + + + + + HbO2 + + + + + Hb + + + + + + + + + + + + + + + + + + Do image processing + + + Perform spectral unmixing + + + + + + + Qt::Horizontal + + + + + + + Qt::Vertical + + + + 17 + 54 + + + + + + + + + + diff --git a/Plugins/org.mitk.gui.qt.photoacoustics.spectralunmixing/src/internal/org_mitk_gui_qt_photoacoustics_spectralunmixing_Activator.cpp b/Plugins/org.mitk.gui.qt.photoacoustics.spectralunmixing/src/internal/org_mitk_gui_qt_photoacoustics_spectralunmixing_Activator.cpp new file mode 100644 index 0000000000..000fc8d83d --- /dev/null +++ b/Plugins/org.mitk.gui.qt.photoacoustics.spectralunmixing/src/internal/org_mitk_gui_qt_photoacoustics_spectralunmixing_Activator.cpp @@ -0,0 +1,29 @@ +/*=================================================================== + +The Medical Imaging Interaction Toolkit (MITK) + +Copyright (c) German Cancer Research Center, +Division of Medical and Biological Informatics. +All rights reserved. + +This software is distributed WITHOUT ANY WARRANTY; without +even the implied warranty of MERCHANTABILITY or FITNESS FOR +A PARTICULAR PURPOSE. + +See LICENSE.txt or http://www.mitk.org for details. + +===================================================================*/ + + +#include "org_mitk_gui_qt_photoacoustics_spectralunmixing_Activator.h" +#include "SpectralUnmixing.h" + +namespace mitk +{ + void org_mitk_gui_qt_photoacoustics_spectralunmixing_Activator::start(ctkPluginContext *context) + { + BERRY_REGISTER_EXTENSION_CLASS(SpectralUnmixing, context) + } + + void org_mitk_gui_qt_photoacoustics_spectralunmixing_Activator::stop(ctkPluginContext *context) { Q_UNUSED(context) } +} diff --git a/Plugins/org.mitk.gui.qt.photoacoustics.spectralunmixing/src/internal/org_mitk_gui_qt_photoacoustics_spectralunmixing_Activator.h b/Plugins/org.mitk.gui.qt.photoacoustics.spectralunmixing/src/internal/org_mitk_gui_qt_photoacoustics_spectralunmixing_Activator.h new file mode 100644 index 0000000000..08d0614225 --- /dev/null +++ b/Plugins/org.mitk.gui.qt.photoacoustics.spectralunmixing/src/internal/org_mitk_gui_qt_photoacoustics_spectralunmixing_Activator.h @@ -0,0 +1,38 @@ +/*=================================================================== + +The Medical Imaging Interaction Toolkit (MITK) + +Copyright (c) German Cancer Research Center, +Division of Medical and Biological Informatics. +All rights reserved. + +This software is distributed WITHOUT ANY WARRANTY; without +even the implied warranty of MERCHANTABILITY or FITNESS FOR +A PARTICULAR PURPOSE. + +See LICENSE.txt or http://www.mitk.org for details. + +===================================================================*/ + + +#ifndef org_mitk_gui_qt_photoacoustics_spectralunmixing_Activator_h +#define org_mitk_gui_qt_photoacoustics_spectralunmixing_Activator_h + +#include + +namespace mitk +{ + class org_mitk_gui_qt_photoacoustics_spectralunmixing_Activator : public QObject, public ctkPluginActivator + { + Q_OBJECT + Q_PLUGIN_METADATA(IID "org_mitk_gui_qt_photoacoustics_spectralunmixing") + Q_INTERFACES(ctkPluginActivator) + + public: + void start(ctkPluginContext *context); + void stop(ctkPluginContext *context); + + }; // org_mitk_gui_qt_photoacoustics_spectralunmixing_Activator +} + +#endif // org_mitk_gui_qt_photoacoustics_spectralunmixing_Activator_h