diff --git a/Modules/DiffusionImaging/FiberTracking/IODataStructures/mitkFiberfoxNoTemplateParameters.cpp b/Modules/DiffusionImaging/FiberTracking/IODataStructures/mitkFiberfoxNoTemplateParameters.cpp new file mode 100644 index 0000000000..e3c89ef9ed --- /dev/null +++ b/Modules/DiffusionImaging/FiberTracking/IODataStructures/mitkFiberfoxNoTemplateParameters.cpp @@ -0,0 +1,150 @@ +/*=================================================================== + +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 "mitkFiberfoxNoTemplateParameters.h" + +void mitk::SignalGenerationParameters::GenerateGradientHalfShell() +{ + int NPoints = 2*m_NumGradients; + m_GradientDirections.clear(); + + m_NumBaseline = NPoints/20; + if (m_NumBaseline==0) + m_NumBaseline=1; + + GradientType g; + g.Fill(0.0); + for (unsigned int i=0; i theta; theta.set_size(NPoints); + vnl_vector phi; phi.set_size(NPoints); + double C = sqrt(4*M_PI); + phi(0) = 0.0; + phi(NPoints-1) = 0.0; + for(int i=0; i0 && i mitk::SignalGenerationParameters::GetBaselineIndices() +{ + std::vector< int > result; + for( unsigned int i=0; im_GradientDirections.size(); i++) + if (m_GradientDirections.at(i).GetNorm()<0.0001) + result.push_back(i); + return result; +} + +unsigned int mitk::SignalGenerationParameters::GetFirstBaselineIndex() +{ + for( unsigned int i=0; im_GradientDirections.size(); i++) + if (m_GradientDirections.at(i).GetNorm()<0.0001) + return i; + return -1; +} + +bool mitk::SignalGenerationParameters::IsBaselineIndex(unsigned int idx) +{ + if (m_GradientDirections.size()>idx && m_GradientDirections.at(idx).GetNorm()<0.0001) + return true; + return false; +} + +unsigned int mitk::SignalGenerationParameters::GetNumWeightedVolumes() +{ + return m_NumGradients; +} + +unsigned int mitk::SignalGenerationParameters::GetNumBaselineVolumes() +{ + return m_NumBaseline; +} + +unsigned int mitk::SignalGenerationParameters::GetNumVolumes() +{ + return m_GradientDirections.size(); +} + +mitk::SignalGenerationParameters::GradientListType mitk::SignalGenerationParameters::GetGradientDirections() +{ + return m_GradientDirections; +} + +mitk::SignalGenerationParameters::GradientType mitk::SignalGenerationParameters::GetGradientDirection(unsigned int i) +{ + return m_GradientDirections.at(i); +} + +void mitk::SignalGenerationParameters::SetNumWeightedVolumes(int numGradients) +{ + m_NumGradients = numGradients; + GenerateGradientHalfShell(); +} + +void mitk::SignalGenerationParameters::SetGradienDirections(GradientListType gradientList) +{ + m_GradientDirections = gradientList; + m_NumGradients = 0; + m_NumBaseline = 0; + for( unsigned int i=0; im_GradientDirections.size(); i++) + { + if (m_GradientDirections.at(i).GetNorm()>0.0001) + m_NumGradients++; + else + m_NumBaseline++; + } +} + +void mitk::SignalGenerationParameters::SetGradienDirections(mitk::DiffusionImage::GradientDirectionContainerType::Pointer gradientList) +{ + m_NumGradients = 0; + m_NumBaseline = 0; + m_GradientDirections.clear(); + + for( unsigned int i=0; iSize(); i++) + { + GradientType g; + g[0] = gradientList->at(i)[0]; + g[1] = gradientList->at(i)[1]; + g[2] = gradientList->at(i)[2]; + m_GradientDirections.push_back(g); + + if (m_GradientDirections.at(i).GetNorm()>0.0001) + m_NumGradients++; + else + m_NumBaseline++; + } +} diff --git a/Modules/DiffusionImaging/FiberTracking/IODataStructures/mitkFiberfoxParameters.h b/Modules/DiffusionImaging/FiberTracking/IODataStructures/mitkFiberfoxNoTemplateParameters.h similarity index 57% copy from Modules/DiffusionImaging/FiberTracking/IODataStructures/mitkFiberfoxParameters.h copy to Modules/DiffusionImaging/FiberTracking/IODataStructures/mitkFiberfoxNoTemplateParameters.h index 4e43573b4b..085e571b7f 100644 --- a/Modules/DiffusionImaging/FiberTracking/IODataStructures/mitkFiberfoxParameters.h +++ b/Modules/DiffusionImaging/FiberTracking/IODataStructures/mitkFiberfoxNoTemplateParameters.h @@ -1,321 +1,228 @@ /*=================================================================== 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 _MITK_FiberfoxParameters_H -#define _MITK_FiberfoxParameters_H +#ifndef _MITK_FiberfoxNoTemplateParameters_H +#define _MITK_FiberfoxNoTemplateParameters_H #include -#include -#include -#include +#include #include -#include -#include #include #include -#include -#include -#include -#include -#include -#include - -using namespace std; - namespace mitk { /** Signal generation */ -class MitkFiberTracking_EXPORT SignalGenerationParameters +class SignalGenerationParameters { public: typedef itk::Image ItkDoubleImgType; typedef itk::Image ItkUcharImgType; typedef itk::Vector GradientType; typedef std::vector GradientListType; enum DiffusionDirectionMode { FIBER_TANGENT_DIRECTIONS, MAIN_FIBER_DIRECTIONS, RANDOM_DIRECTIONS }; SignalGenerationParameters() - : m_FrequencyMap(NULL) - , m_MaskImage(NULL) - , m_SignalScale(100) + : m_SignalScale(100) , m_tEcho(100) , m_tLine(1) , m_tInhom(50) , m_Bvalue(1000) , m_SimulateKspaceAcquisition(false) , m_AxonRadius(0) , m_DiffusionDirectionMode(SignalGenerationParameters::FIBER_TANGENT_DIRECTIONS) , m_FiberSeparationThreshold(30) , m_Spikes(0) , m_SpikeAmplitude(1) , m_KspaceLineOffset(0) , m_EddyStrength(0) , m_Tau(70) , m_CroppingFactor(1) , m_DoAddGibbsRinging(false) , m_DoSimulateRelaxation(true) , m_DoDisablePartialVolume(false) , m_DoAddMotion(false) , m_DoRandomizeMotion(true) + , m_FrequencyMap(NULL) + , m_MaskImage(NULL) { m_ImageRegion.SetSize(0, 12); m_ImageRegion.SetSize(1, 12); m_ImageRegion.SetSize(2, 3); m_ImageSpacing.Fill(2.0); m_ImageOrigin.Fill(0.0); m_ImageDirection.SetIdentity(); m_Translation.Fill(0.0); m_Rotation.Fill(0.0); SetNumWeightedVolumes(6); } + + ~SignalGenerationParameters() + { + } + /** input/output image specifications */ itk::ImageRegion<3> m_ImageRegion; ///< Image size. itk::Vector m_ImageSpacing; ///< Image voxel size. itk::Point m_ImageOrigin; ///< Image origin. itk::Matrix m_ImageDirection; ///< Image rotation matrix. /** Other acquisitions parameters */ double m_SignalScale; ///< Scaling factor for output signal (before noise is added). double m_tEcho; ///< Echo time TE. double m_tLine; ///< k-space line readout time. double m_tInhom; ///< T2' double m_Bvalue; ///< Acquisition b-value bool m_SimulateKspaceAcquisition;///< Flag to enable/disable k-space acquisition simulation double m_AxonRadius; ///< Determines compartment volume fractions (0 == automatic axon radius estimation) bool m_DoDisablePartialVolume; ///< Disable partial volume effects. Each voxel is either all fiber or all non-fiber. DiffusionDirectionMode m_DiffusionDirectionMode; ///< Determines how the main diffusion direction of the signal models is selected double m_FiberSeparationThreshold; ///< Used for random and and mein fiber deriction DiffusionDirectionMode /** Artifacts and other effects */ unsigned int m_Spikes; ///< Number of spikes randomly appearing in the image double m_SpikeAmplitude; ///< amplitude of spikes relative to the largest signal intensity (magnitude of complex) double m_KspaceLineOffset; ///< Causes N/2 ghosts. Larger offset means stronger ghost. double m_EddyStrength; ///< Strength of eddy current induced gradients in mT/m. double m_Tau; ///< Eddy current decay constant (in ms) double m_CroppingFactor; ///< FOV size in y-direction is multiplied by this factor. Causes aliasing artifacts. bool m_DoAddGibbsRinging; ///< Add Gibbs ringing artifact bool m_DoSimulateRelaxation; ///< Add T2 relaxation effects bool m_DoAddMotion; ///< Enable motion artifacts. bool m_DoRandomizeMotion; ///< Toggles between random and linear motion. itk::Vector m_Translation; ///< Maximum translational motion. itk::Vector m_Rotation; ///< Maximum rotational motion. ItkDoubleImgType::Pointer m_FrequencyMap; ///< If != NULL, distortions are added to the image using this frequency map. ItkUcharImgType::Pointer m_MaskImage; ///< Signal is only genrated inside of the mask image. - inline void GenerateGradientHalfShell(); ///< Generates half shell of gradient directions (with m_NumGradients non-zero directions) - inline std::vector< int > GetBaselineIndices(); ///< Returns list of nun-diffusion-weighted image volume indices - inline unsigned int GetFirstBaselineIndex(); ///< Returns index of first non-diffusion-weighted image volume - inline bool IsBaselineIndex(unsigned int idx); ///< Checks if image volume with given index is non-diffusion-weighted volume or not. - inline unsigned int GetNumWeightedVolumes(); ///< Get number of diffusion-weighted image volumes - inline unsigned int GetNumBaselineVolumes(); ///< Get number of non-diffusion-weighted image volumes - inline unsigned int GetNumVolumes(); ///< Get number of baseline and diffusion-weighted image volumes - inline GradientListType GetGradientDirections(); ///< Return gradient direction container - inline GradientType GetGradientDirection(unsigned int i); + void GenerateGradientHalfShell(); ///< Generates half shell of gradient directions (with m_NumGradients non-zero directions) + std::vector< int > GetBaselineIndices(); ///< Returns list of nun-diffusion-weighted image volume indices + unsigned int GetFirstBaselineIndex(); ///< Returns index of first non-diffusion-weighted image volume + bool IsBaselineIndex(unsigned int idx); ///< Checks if image volume with given index is non-diffusion-weighted volume or not. + unsigned int GetNumWeightedVolumes(); ///< Get number of diffusion-weighted image volumes + unsigned int GetNumBaselineVolumes(); ///< Get number of non-diffusion-weighted image volumes + unsigned int GetNumVolumes(); ///< Get number of baseline and diffusion-weighted image volumes + GradientListType GetGradientDirections(); ///< Return gradient direction container + GradientType GetGradientDirection(unsigned int i); - inline void SetNumWeightedVolumes(int numGradients); ///< Automaticall calls GenerateGradientHalfShell() afterwards. - inline void SetGradienDirections(GradientListType gradientList); - inline void SetGradienDirections(mitk::DiffusionImage::GradientDirectionContainerType::Pointer gradientList); + void SetNumWeightedVolumes(int numGradients); ///< Automaticall calls GenerateGradientHalfShell() afterwards. + void SetGradienDirections(GradientListType gradientList); + void SetGradienDirections(mitk::DiffusionImage::GradientDirectionContainerType::Pointer gradientList); protected: unsigned int m_NumGradients; ///< Number of diffusion-weighted image volumes. unsigned int m_NumBaseline; ///< Number of non-diffusion-weighted image volumes. GradientListType m_GradientDirections; ///< Total number of image volumes. }; /** Fiber generation */ -class MitkFiberTracking_EXPORT FiberGenerationParameters +class FiberGenerationParameters { public: enum FiberDistribution{ DISTRIBUTE_UNIFORM, // distribute fibers uniformly in the ROIs DISTRIBUTE_GAUSSIAN // distribute fibers using a 2D gaussian }; - typedef vector< vector< mitk::PlanarEllipse::Pointer > > FiducialListType; - typedef vector< vector< unsigned int > > FlipListType; + typedef std::vector< std::vector< mitk::PlanarEllipse::Pointer > > FiducialListType; + typedef std::vector< std::vector< unsigned int > > FlipListType; FiberGenerationParameters() : m_Distribution(DISTRIBUTE_UNIFORM) , m_Density(100) , m_Variance(100) , m_Sampling(1) , m_Tension(0) , m_Continuity(0) , m_Bias(0) { m_Rotation.Fill(0.0); m_Translation.Fill(0.0); m_Scale.Fill(1.0); } FiberDistribution m_Distribution; unsigned int m_Density; double m_Variance; double m_Sampling; double m_Tension; double m_Continuity; double m_Bias; mitk::Vector3D m_Rotation; mitk::Vector3D m_Translation; mitk::Vector3D m_Scale; FlipListType m_FlipList; ///< contains flags indicating a flip of the 2D fiber x-coordinates (needed to resolve some unwanted fiber twisting) FiducialListType m_Fiducials; ///< container of the planar ellipses used as fiducials for the fiber generation process }; /** GUI persistence, input, output, ... */ -class MitkFiberTracking_EXPORT MiscFiberfoxParameters +class MiscFiberfoxParameters { public: MiscFiberfoxParameters() : m_ResultNode(DataNode::New()) , m_ParentNode(NULL) , m_SignalModelString("") , m_ArtifactModelString("") , m_OutputPath("") , m_CheckOutputVolumeFractionsBox(false) , m_CheckAdvancedSignalOptionsBox(false) , m_CheckAddNoiseBox(false) , m_CheckAddGhostsBox(false) , m_CheckAddAliasingBox(false) , m_CheckAddSpikesBox(false) , m_CheckAddEddyCurrentsBox(false) , m_CheckAddDistortionsBox(false) , m_CheckRealTimeFibersBox(true) , m_CheckAdvancedFiberOptionsBox(false) , m_CheckConstantRadiusBox(false) , m_CheckIncludeFiducialsBox(true) {} DataNode::Pointer m_ResultNode; ///< Stores resulting image. DataNode::Pointer m_ParentNode; ///< Parent node of result node. - string m_SignalModelString; ///< Appendet to the name of the result node - string m_ArtifactModelString; ///< Appendet to the name of the result node - string m_OutputPath; ///< Image is automatically saved to the specified folder after simulation is finished. + std::string m_SignalModelString; ///< Appendet to the name of the result node + std::string m_ArtifactModelString; ///< Appendet to the name of the result node + std::string m_OutputPath; ///< Image is automatically saved to the specified folder after simulation is finished. /** member variables that store the check-state of GUI checkboxes */ // image generation bool m_CheckOutputVolumeFractionsBox; bool m_CheckAdvancedSignalOptionsBox; bool m_CheckAddNoiseBox; bool m_CheckAddGhostsBox; bool m_CheckAddAliasingBox; bool m_CheckAddSpikesBox; bool m_CheckAddEddyCurrentsBox; bool m_CheckAddDistortionsBox; // fiber generation bool m_CheckRealTimeFibersBox; bool m_CheckAdvancedFiberOptionsBox; bool m_CheckConstantRadiusBox; bool m_CheckIncludeFiducialsBox; }; - -/** - * \brief Datastructure to manage the Fiberfox signal generation parameters. - * - */ -template< class ScalarType = double > -class MitkFiberTracking_EXPORT FiberfoxParameters -{ -public: - - typedef itk::Image ItkDoubleImgType; - typedef itk::Image ItkUcharImgType; - typedef DiffusionSignalModel DiffusionModelType; - typedef std::vector< DiffusionModelType* > DiffusionModelListType; - typedef DiffusionNoiseModel NoiseModelType; - - FiberfoxParameters(); - ~FiberfoxParameters(); - - /** Get same parameter object with different template parameter */ - template< class OutType > - FiberfoxParameters< OutType > CopyParameters() - { - FiberfoxParameters< OutType > out; - - out.m_FiberGen = m_FiberGen; - out.m_SignalGen = m_SignalGen; - out.m_Misc = m_Misc; - - if (m_NoiseModel!=NULL) - { - if (dynamic_cast*>(m_NoiseModel)) - out.m_NoiseModel = new mitk::RicianNoiseModel(); - else if (dynamic_cast*>(m_NoiseModel)) - out.m_NoiseModel = new mitk::ChiSquareNoiseModel(); - out.m_NoiseModel->SetNoiseVariance(m_NoiseModel->GetNoiseVariance()); - } - - for (unsigned int i=0; i* outModel = NULL; - mitk::DiffusionSignalModel* signalModel = NULL; - if (i*>(signalModel)) - outModel = new mitk::StickModel(dynamic_cast*>(signalModel)); - else if (dynamic_cast*>(signalModel)) - outModel = new mitk::TensorModel(dynamic_cast*>(signalModel)); - else if (dynamic_cast*>(signalModel)) - outModel = new mitk::RawShModel(dynamic_cast*>(signalModel)); - else if (dynamic_cast*>(signalModel)) - outModel = new mitk::BallModel(dynamic_cast*>(signalModel)); - else if (dynamic_cast*>(signalModel)) - outModel = new mitk::AstroStickModel(dynamic_cast*>(signalModel)); - else if (dynamic_cast*>(signalModel)) - outModel = new mitk::DotModel(dynamic_cast*>(signalModel)); - - if (i #include #include #include #include #include template< class ScalarType > mitk::FiberfoxParameters< ScalarType >::FiberfoxParameters() : m_NoiseModel(NULL) { } template< class ScalarType > mitk::FiberfoxParameters< ScalarType >::~FiberfoxParameters() { // if (m_NoiseModel!=NULL) // delete m_NoiseModel; } -void mitk::SignalGenerationParameters::GenerateGradientHalfShell() -{ - int NPoints = 2*m_NumGradients; - m_GradientDirections.clear(); - - m_NumBaseline = NPoints/20; - if (m_NumBaseline==0) - m_NumBaseline=1; - - GradientType g; - g.Fill(0.0); - for (unsigned int i=0; i theta; theta.set_size(NPoints); - vnl_vector phi; phi.set_size(NPoints); - double C = sqrt(4*M_PI); - phi(0) = 0.0; - phi(NPoints-1) = 0.0; - for(int i=0; i0 && i mitk::SignalGenerationParameters::GetBaselineIndices() -{ - std::vector< int > result; - for( unsigned int i=0; im_GradientDirections.size(); i++) - if (m_GradientDirections.at(i).GetNorm()<0.0001) - result.push_back(i); - return result; -} - -unsigned int mitk::SignalGenerationParameters::GetFirstBaselineIndex() -{ - for( unsigned int i=0; im_GradientDirections.size(); i++) - if (m_GradientDirections.at(i).GetNorm()<0.0001) - return i; - return -1; -} - -bool mitk::SignalGenerationParameters::IsBaselineIndex(unsigned int idx) -{ - if (m_GradientDirections.size()>idx && m_GradientDirections.at(idx).GetNorm()<0.0001) - return true; - return false; -} - -unsigned int mitk::SignalGenerationParameters::GetNumWeightedVolumes() -{ - return m_NumGradients; -} - -unsigned int mitk::SignalGenerationParameters::GetNumBaselineVolumes() -{ - return m_NumBaseline; -} - -unsigned int mitk::SignalGenerationParameters::GetNumVolumes() -{ - return m_GradientDirections.size(); -} - -mitk::SignalGenerationParameters::GradientListType mitk::SignalGenerationParameters::GetGradientDirections() -{ - return m_GradientDirections; -} - -mitk::SignalGenerationParameters::GradientType mitk::SignalGenerationParameters::GetGradientDirection(unsigned int i) -{ - return m_GradientDirections.at(i); -} - -void mitk::SignalGenerationParameters::SetNumWeightedVolumes(int numGradients) -{ - m_NumGradients = numGradients; - GenerateGradientHalfShell(); -} - -void mitk::SignalGenerationParameters::SetGradienDirections(GradientListType gradientList) -{ - m_GradientDirections = gradientList; - m_NumGradients = 0; - m_NumBaseline = 0; - for( unsigned int i=0; im_GradientDirections.size(); i++) - { - if (m_GradientDirections.at(i).GetNorm()>0.0001) - m_NumGradients++; - else - m_NumBaseline++; - } -} - -void mitk::SignalGenerationParameters::SetGradienDirections(mitk::DiffusionImage::GradientDirectionContainerType::Pointer gradientList) -{ - m_NumGradients = 0; - m_NumBaseline = 0; - m_GradientDirections.clear(); - - for( unsigned int i=0; iSize(); i++) - { - GradientType g; - g[0] = gradientList->at(i)[0]; - g[1] = gradientList->at(i)[1]; - g[2] = gradientList->at(i)[2]; - m_GradientDirections.push_back(g); - - if (m_GradientDirections.at(i).GetNorm()>0.0001) - m_NumGradients++; - else - m_NumBaseline++; - } -} - template< class ScalarType > -void mitk::FiberfoxParameters< ScalarType >::SaveParameters(string filename) +void mitk::FiberfoxParameters< ScalarType >::SaveParameters(std::string filename) { if(filename.empty()) return; if(".ffp"!=filename.substr(filename.size()-4, 4)) filename += ".ffp"; boost::property_tree::ptree parameters; // fiber generation parameters parameters.put("fiberfox.fibers.distribution", m_FiberGen.m_Distribution); parameters.put("fiberfox.fibers.variance", m_FiberGen.m_Variance); parameters.put("fiberfox.fibers.density", m_FiberGen.m_Density); parameters.put("fiberfox.fibers.spline.sampling", m_FiberGen.m_Sampling); parameters.put("fiberfox.fibers.spline.tension", m_FiberGen.m_Tension); parameters.put("fiberfox.fibers.spline.continuity", m_FiberGen.m_Continuity); parameters.put("fiberfox.fibers.spline.bias", m_FiberGen.m_Bias); parameters.put("fiberfox.fibers.rotation.x", m_FiberGen.m_Rotation[0]); parameters.put("fiberfox.fibers.rotation.y", m_FiberGen.m_Rotation[1]); parameters.put("fiberfox.fibers.rotation.z", m_FiberGen.m_Rotation[2]); parameters.put("fiberfox.fibers.translation.x", m_FiberGen.m_Translation[0]); parameters.put("fiberfox.fibers.translation.y", m_FiberGen.m_Translation[1]); parameters.put("fiberfox.fibers.translation.z", m_FiberGen.m_Translation[2]); parameters.put("fiberfox.fibers.scale.x", m_FiberGen.m_Scale[0]); parameters.put("fiberfox.fibers.scale.y", m_FiberGen.m_Scale[1]); parameters.put("fiberfox.fibers.scale.z", m_FiberGen.m_Scale[2]); // image generation parameters parameters.put("fiberfox.image.basic.size.x", m_SignalGen.m_ImageRegion.GetSize(0)); parameters.put("fiberfox.image.basic.size.y", m_SignalGen.m_ImageRegion.GetSize(1)); parameters.put("fiberfox.image.basic.size.z", m_SignalGen.m_ImageRegion.GetSize(2)); parameters.put("fiberfox.image.basic.spacing.x", m_SignalGen.m_ImageSpacing[0]); parameters.put("fiberfox.image.basic.spacing.y", m_SignalGen.m_ImageSpacing[1]); parameters.put("fiberfox.image.basic.spacing.z", m_SignalGen.m_ImageSpacing[2]); parameters.put("fiberfox.image.basic.origin.x", m_SignalGen.m_ImageOrigin[0]); parameters.put("fiberfox.image.basic.origin.y", m_SignalGen.m_ImageOrigin[1]); parameters.put("fiberfox.image.basic.origin.z", m_SignalGen.m_ImageOrigin[2]); parameters.put("fiberfox.image.basic.direction.1", m_SignalGen.m_ImageDirection[0][0]); parameters.put("fiberfox.image.basic.direction.2", m_SignalGen.m_ImageDirection[0][1]); parameters.put("fiberfox.image.basic.direction.3", m_SignalGen.m_ImageDirection[0][2]); parameters.put("fiberfox.image.basic.direction.4", m_SignalGen.m_ImageDirection[1][0]); parameters.put("fiberfox.image.basic.direction.5", m_SignalGen.m_ImageDirection[1][1]); parameters.put("fiberfox.image.basic.direction.6", m_SignalGen.m_ImageDirection[1][2]); parameters.put("fiberfox.image.basic.direction.7", m_SignalGen.m_ImageDirection[2][0]); parameters.put("fiberfox.image.basic.direction.8", m_SignalGen.m_ImageDirection[2][1]); parameters.put("fiberfox.image.basic.direction.9", m_SignalGen.m_ImageDirection[2][2]); parameters.put("fiberfox.image.basic.numgradients", m_SignalGen.GetNumWeightedVolumes()); for( unsigned int i=0; im_SignalGen.GetNumVolumes(); i++) { - parameters.put("fiberfox.image.gradients."+boost::lexical_cast(i)+".x", m_SignalGen.GetGradientDirection(i)[0]); - parameters.put("fiberfox.image.gradients."+boost::lexical_cast(i)+".y", m_SignalGen.GetGradientDirection(i)[1]); - parameters.put("fiberfox.image.gradients."+boost::lexical_cast(i)+".z", m_SignalGen.GetGradientDirection(i)[2]); + parameters.put("fiberfox.image.gradients."+boost::lexical_cast(i)+".x", m_SignalGen.GetGradientDirection(i)[0]); + parameters.put("fiberfox.image.gradients."+boost::lexical_cast(i)+".y", m_SignalGen.GetGradientDirection(i)[1]); + parameters.put("fiberfox.image.gradients."+boost::lexical_cast(i)+".z", m_SignalGen.GetGradientDirection(i)[2]); } parameters.put("fiberfox.image.signalScale", m_SignalGen.m_SignalScale); parameters.put("fiberfox.image.tEcho", m_SignalGen.m_tEcho); parameters.put("fiberfox.image.tLine", m_SignalGen.m_tLine); parameters.put("fiberfox.image.tInhom", m_SignalGen.m_tInhom); parameters.put("fiberfox.image.bvalue", m_SignalGen.m_Bvalue); parameters.put("fiberfox.image.simulatekspace", m_SignalGen.m_SimulateKspaceAcquisition); parameters.put("fiberfox.image.axonRadius", m_SignalGen.m_AxonRadius); parameters.put("fiberfox.image.diffusiondirectionmode", m_SignalGen.m_DiffusionDirectionMode); parameters.put("fiberfox.image.fiberseparationthreshold", m_SignalGen.m_FiberSeparationThreshold); parameters.put("fiberfox.image.doSimulateRelaxation", m_SignalGen.m_DoSimulateRelaxation); parameters.put("fiberfox.image.doDisablePartialVolume", m_SignalGen.m_DoDisablePartialVolume); parameters.put("fiberfox.image.artifacts.spikesnum", m_SignalGen.m_Spikes); parameters.put("fiberfox.image.artifacts.spikesscale", m_SignalGen.m_SpikeAmplitude); parameters.put("fiberfox.image.artifacts.kspaceLineOffset", m_SignalGen.m_KspaceLineOffset); parameters.put("fiberfox.image.artifacts.eddyStrength", m_SignalGen.m_EddyStrength); parameters.put("fiberfox.image.artifacts.eddyTau", m_SignalGen.m_Tau); parameters.put("fiberfox.image.artifacts.aliasingfactor", m_SignalGen.m_CroppingFactor); parameters.put("fiberfox.image.artifacts.addringing", m_SignalGen.m_DoAddGibbsRinging); parameters.put("fiberfox.image.artifacts.doAddMotion", m_SignalGen.m_DoAddMotion); parameters.put("fiberfox.image.artifacts.randomMotion", m_SignalGen.m_DoRandomizeMotion); parameters.put("fiberfox.image.artifacts.translation0", m_SignalGen.m_Translation[0]); parameters.put("fiberfox.image.artifacts.translation1", m_SignalGen.m_Translation[1]); parameters.put("fiberfox.image.artifacts.translation2", m_SignalGen.m_Translation[2]); parameters.put("fiberfox.image.artifacts.rotation0", m_SignalGen.m_Rotation[0]); parameters.put("fiberfox.image.artifacts.rotation1", m_SignalGen.m_Rotation[1]); parameters.put("fiberfox.image.artifacts.rotation2", m_SignalGen.m_Rotation[2]); parameters.put("fiberfox.image.artifacts.addnoise", m_Misc.m_CheckAddNoiseBox); parameters.put("fiberfox.image.artifacts.addghosts", m_Misc.m_CheckAddGhostsBox); parameters.put("fiberfox.image.artifacts.addaliasing", m_Misc.m_CheckAddAliasingBox); parameters.put("fiberfox.image.artifacts.addspikes", m_Misc.m_CheckAddSpikesBox); parameters.put("fiberfox.image.artifacts.addeddycurrents", m_Misc.m_CheckAddEddyCurrentsBox); parameters.put("fiberfox.image.artifacts.doAddDistortions", m_Misc.m_CheckAddDistortionsBox); parameters.put("fiberfox.image.outputvolumefractions", m_Misc.m_CheckOutputVolumeFractionsBox); parameters.put("fiberfox.image.showadvanced", m_Misc.m_CheckAdvancedSignalOptionsBox); parameters.put("fiberfox.image.signalmodelstring", m_Misc.m_SignalModelString); parameters.put("fiberfox.image.artifactmodelstring", m_Misc.m_ArtifactModelString); parameters.put("fiberfox.image.outpath", m_Misc.m_OutputPath); parameters.put("fiberfox.fibers.realtime", m_Misc.m_CheckRealTimeFibersBox); parameters.put("fiberfox.fibers.showadvanced", m_Misc.m_CheckAdvancedFiberOptionsBox); parameters.put("fiberfox.fibers.constantradius", m_Misc.m_CheckConstantRadiusBox); parameters.put("fiberfox.fibers.includeFiducials", m_Misc.m_CheckIncludeFiducialsBox); if (m_NoiseModel!=NULL) { parameters.put("fiberfox.image.artifacts.noisevariance", m_NoiseModel->GetNoiseVariance()); if (dynamic_cast*>(m_NoiseModel)) parameters.put("fiberfox.image.artifacts.noisetype", "rice"); else if (dynamic_cast*>(m_NoiseModel)) parameters.put("fiberfox.image.artifacts.noisetype", "chisquare"); } for (int i=0; i* signalModel = NULL; if (i(i)+".type", "fiber"); + parameters.put("fiberfox.image.compartments."+boost::lexical_cast(i)+".type", "fiber"); } else { signalModel = m_NonFiberModelList.at(i-m_FiberModelList.size()); - parameters.put("fiberfox.image.compartments."+boost::lexical_cast(i)+".type", "non-fiber"); + parameters.put("fiberfox.image.compartments."+boost::lexical_cast(i)+".type", "non-fiber"); } if (dynamic_cast*>(signalModel)) { mitk::StickModel* model = dynamic_cast*>(signalModel); - parameters.put("fiberfox.image.compartments."+boost::lexical_cast(i)+".stick.d", model->GetDiffusivity()); - parameters.put("fiberfox.image.compartments."+boost::lexical_cast(i)+".stick.t2", model->GetT2()); + parameters.put("fiberfox.image.compartments."+boost::lexical_cast(i)+".stick.d", model->GetDiffusivity()); + parameters.put("fiberfox.image.compartments."+boost::lexical_cast(i)+".stick.t2", model->GetT2()); } else if (dynamic_cast*>(signalModel)) { mitk::TensorModel* model = dynamic_cast*>(signalModel); - parameters.put("fiberfox.image.compartments."+boost::lexical_cast(i)+".tensor.d1", model->GetDiffusivity1()); - parameters.put("fiberfox.image.compartments."+boost::lexical_cast(i)+".tensor.d2", model->GetDiffusivity2()); - parameters.put("fiberfox.image.compartments."+boost::lexical_cast(i)+".tensor.d3", model->GetDiffusivity3()); - parameters.put("fiberfox.image.compartments."+boost::lexical_cast(i)+".tensor.t2", model->GetT2()); + parameters.put("fiberfox.image.compartments."+boost::lexical_cast(i)+".tensor.d1", model->GetDiffusivity1()); + parameters.put("fiberfox.image.compartments."+boost::lexical_cast(i)+".tensor.d2", model->GetDiffusivity2()); + parameters.put("fiberfox.image.compartments."+boost::lexical_cast(i)+".tensor.d3", model->GetDiffusivity3()); + parameters.put("fiberfox.image.compartments."+boost::lexical_cast(i)+".tensor.t2", model->GetT2()); } else if (dynamic_cast*>(signalModel)) { mitk::RawShModel* model = dynamic_cast*>(signalModel); - parameters.put("fiberfox.image.compartments."+boost::lexical_cast(i)+".prototype.minFA", model->GetFaRange().first); - parameters.put("fiberfox.image.compartments."+boost::lexical_cast(i)+".prototype.maxFA", model->GetFaRange().second); - parameters.put("fiberfox.image.compartments."+boost::lexical_cast(i)+".prototype.minADC", model->GetAdcRange().first); - parameters.put("fiberfox.image.compartments."+boost::lexical_cast(i)+".prototype.maxADC", model->GetAdcRange().second); - parameters.put("fiberfox.image.compartments."+boost::lexical_cast(i)+".prototype.numSamples", model->GetMaxNumKernels()); + parameters.put("fiberfox.image.compartments."+boost::lexical_cast(i)+".prototype.minFA", model->GetFaRange().first); + parameters.put("fiberfox.image.compartments."+boost::lexical_cast(i)+".prototype.maxFA", model->GetFaRange().second); + parameters.put("fiberfox.image.compartments."+boost::lexical_cast(i)+".prototype.minADC", model->GetAdcRange().first); + parameters.put("fiberfox.image.compartments."+boost::lexical_cast(i)+".prototype.maxADC", model->GetAdcRange().second); + parameters.put("fiberfox.image.compartments."+boost::lexical_cast(i)+".prototype.numSamples", model->GetMaxNumKernels()); } else if (dynamic_cast*>(signalModel)) { mitk::BallModel* model = dynamic_cast*>(signalModel); - parameters.put("fiberfox.image.compartments."+boost::lexical_cast(i)+".ball.d", model->GetDiffusivity()); - parameters.put("fiberfox.image.compartments."+boost::lexical_cast(i)+".ball.t2", model->GetT2()); + parameters.put("fiberfox.image.compartments."+boost::lexical_cast(i)+".ball.d", model->GetDiffusivity()); + parameters.put("fiberfox.image.compartments."+boost::lexical_cast(i)+".ball.t2", model->GetT2()); } else if (dynamic_cast*>(signalModel)) { mitk::AstroStickModel* model = dynamic_cast*>(signalModel); - parameters.put("fiberfox.image.compartments."+boost::lexical_cast(i)+".astrosticks.d", model->GetDiffusivity()); - parameters.put("fiberfox.image.compartments."+boost::lexical_cast(i)+".astrosticks.t2", model->GetT2()); - parameters.put("fiberfox.image.compartments."+boost::lexical_cast(i)+".astrosticks.randomize", model->GetRandomizeSticks()); + parameters.put("fiberfox.image.compartments."+boost::lexical_cast(i)+".astrosticks.d", model->GetDiffusivity()); + parameters.put("fiberfox.image.compartments."+boost::lexical_cast(i)+".astrosticks.t2", model->GetT2()); + parameters.put("fiberfox.image.compartments."+boost::lexical_cast(i)+".astrosticks.randomize", model->GetRandomizeSticks()); } else if (dynamic_cast*>(signalModel)) { mitk::DotModel* model = dynamic_cast*>(signalModel); - parameters.put("fiberfox.image.compartments."+boost::lexical_cast(i)+".dot.t2", model->GetT2()); + parameters.put("fiberfox.image.compartments."+boost::lexical_cast(i)+".dot.t2", model->GetT2()); } if (signalModel!=NULL) { - parameters.put("fiberfox.image.compartments."+boost::lexical_cast(i)+".ID", signalModel->m_CompartmentId); + parameters.put("fiberfox.image.compartments."+boost::lexical_cast(i)+".ID", signalModel->m_CompartmentId); if (signalModel->GetVolumeFractionImage().IsNotNull()) { try{ itk::ImageFileWriter::Pointer writer = itk::ImageFileWriter::New(); - writer->SetFileName(filename+"_VOLUME"+boost::lexical_cast(signalModel->m_CompartmentId)+".nrrd"); + writer->SetFileName(filename+"_VOLUME"+boost::lexical_cast(signalModel->m_CompartmentId)+".nrrd"); writer->SetInput(signalModel->GetVolumeFractionImage()); writer->Update(); - MITK_INFO << "Volume fraction image for compartment "+boost::lexical_cast(signalModel->m_CompartmentId)+" saved."; + MITK_INFO << "Volume fraction image for compartment "+boost::lexical_cast(signalModel->m_CompartmentId)+" saved."; } catch(...) { } } } } boost::property_tree::xml_writer_settings writerSettings(' ', 2); boost::property_tree::xml_parser::write_xml(filename, parameters, std::locale(), writerSettings); try{ itk::ImageFileWriter::Pointer writer = itk::ImageFileWriter::New(); writer->SetFileName(filename+"_FMAP.nrrd"); writer->SetInput(m_SignalGen.m_FrequencyMap); writer->Update(); } catch(...) { MITK_INFO << "No frequency map saved."; } try{ itk::ImageFileWriter::Pointer writer = itk::ImageFileWriter::New(); writer->SetFileName(filename+"_MASK.nrrd"); writer->SetInput(m_SignalGen.m_MaskImage); writer->Update(); } catch(...) { MITK_INFO << "No mask image saved."; } } template< class ScalarType > -void mitk::FiberfoxParameters< ScalarType >::LoadParameters(string filename) +void mitk::FiberfoxParameters< ScalarType >::LoadParameters(std::string filename) { if(filename.empty()) return; boost::property_tree::ptree parameterTree; boost::property_tree::xml_parser::read_xml(filename, parameterTree); m_FiberModelList.clear(); m_NonFiberModelList.clear(); if (m_NoiseModel!=NULL) delete m_NoiseModel; BOOST_FOREACH( boost::property_tree::ptree::value_type const& v1, parameterTree.get_child("fiberfox") ) { if( v1.first == "fibers" ) { m_Misc.m_CheckRealTimeFibersBox = v1.second.get("realtime", m_Misc.m_CheckRealTimeFibersBox); m_Misc.m_CheckAdvancedFiberOptionsBox = v1.second.get("showadvanced", m_Misc.m_CheckAdvancedFiberOptionsBox); m_Misc.m_CheckConstantRadiusBox = v1.second.get("constantradius", m_Misc.m_CheckConstantRadiusBox); m_Misc.m_CheckIncludeFiducialsBox = v1.second.get("includeFiducials", m_Misc.m_CheckIncludeFiducialsBox); switch (v1.second.get("distribution", 0)) { case 0: m_FiberGen.m_Distribution = FiberGenerationParameters::DISTRIBUTE_UNIFORM; break; case 1: m_FiberGen.m_Distribution = FiberGenerationParameters::DISTRIBUTE_GAUSSIAN; break; default: m_FiberGen.m_Distribution = FiberGenerationParameters::DISTRIBUTE_UNIFORM; } m_FiberGen.m_Variance = v1.second.get("variance", m_FiberGen.m_Variance); m_FiberGen.m_Density = v1.second.get("density", m_FiberGen.m_Density); m_FiberGen.m_Sampling = v1.second.get("spline.sampling", m_FiberGen.m_Sampling); m_FiberGen.m_Tension = v1.second.get("spline.tension", m_FiberGen.m_Tension); m_FiberGen.m_Continuity = v1.second.get("spline.continuity", m_FiberGen.m_Continuity); m_FiberGen.m_Bias = v1.second.get("spline.bias", m_FiberGen.m_Bias); m_FiberGen.m_Rotation[0] = v1.second.get("rotation.x", m_FiberGen.m_Rotation[0]); m_FiberGen.m_Rotation[1] = v1.second.get("rotation.y", m_FiberGen.m_Rotation[1]); m_FiberGen.m_Rotation[2] = v1.second.get("rotation.z", m_FiberGen.m_Rotation[2]); m_FiberGen.m_Translation[0] = v1.second.get("translation.x", m_FiberGen.m_Translation[0]); m_FiberGen.m_Translation[1] = v1.second.get("translation.y", m_FiberGen.m_Translation[1]); m_FiberGen.m_Translation[2] = v1.second.get("translation.z", m_FiberGen.m_Translation[2]); m_FiberGen.m_Scale[0] = v1.second.get("scale.x", m_FiberGen.m_Scale[0]); m_FiberGen.m_Scale[1] = v1.second.get("scale.y", m_FiberGen.m_Scale[1]); m_FiberGen.m_Scale[2] = v1.second.get("scale.z", m_FiberGen.m_Scale[2]); } else if ( v1.first == "image" ) { - m_Misc.m_SignalModelString = v1.second.get("signalmodelstring", m_Misc.m_SignalModelString); - m_Misc.m_ArtifactModelString = v1.second.get("artifactmodelstring", m_Misc.m_ArtifactModelString); - m_Misc.m_OutputPath = v1.second.get("outpath", m_Misc.m_OutputPath); + m_Misc.m_SignalModelString = v1.second.get("signalmodelstring", m_Misc.m_SignalModelString); + m_Misc.m_ArtifactModelString = v1.second.get("artifactmodelstring", m_Misc.m_ArtifactModelString); + m_Misc.m_OutputPath = v1.second.get("outpath", m_Misc.m_OutputPath); m_Misc.m_CheckOutputVolumeFractionsBox = v1.second.get("outputvolumefractions", m_Misc.m_CheckOutputVolumeFractionsBox); m_Misc.m_CheckAdvancedSignalOptionsBox = v1.second.get("showadvanced", m_Misc.m_CheckAdvancedSignalOptionsBox); m_Misc.m_CheckAddDistortionsBox = v1.second.get("artifacts.doAddDistortions", m_Misc.m_CheckAddDistortionsBox); m_Misc.m_CheckAddNoiseBox = v1.second.get("artifacts.addnoise", m_Misc.m_CheckAddNoiseBox); m_Misc.m_CheckAddGhostsBox = v1.second.get("artifacts.addghosts", m_Misc.m_CheckAddGhostsBox); m_Misc.m_CheckAddAliasingBox = v1.second.get("artifacts.addaliasing", m_Misc.m_CheckAddAliasingBox); m_Misc.m_CheckAddSpikesBox = v1.second.get("artifacts.addspikes", m_Misc.m_CheckAddSpikesBox); m_Misc.m_CheckAddEddyCurrentsBox = v1.second.get("artifacts.addeddycurrents", m_Misc.m_CheckAddEddyCurrentsBox); m_SignalGen.m_ImageRegion.SetSize(0, v1.second.get("basic.size.x",m_SignalGen.m_ImageRegion.GetSize(0))); m_SignalGen.m_ImageRegion.SetSize(1, v1.second.get("basic.size.y",m_SignalGen.m_ImageRegion.GetSize(1))); m_SignalGen.m_ImageRegion.SetSize(2, v1.second.get("basic.size.z",m_SignalGen.m_ImageRegion.GetSize(2))); m_SignalGen.m_ImageSpacing[0] = v1.second.get("basic.spacing.x",m_SignalGen.m_ImageSpacing[0]); m_SignalGen.m_ImageSpacing[1] = v1.second.get("basic.spacing.y",m_SignalGen.m_ImageSpacing[1]); m_SignalGen.m_ImageSpacing[2] = v1.second.get("basic.spacing.z",m_SignalGen.m_ImageSpacing[2]); m_SignalGen.m_ImageOrigin[0] = v1.second.get("basic.origin.x",m_SignalGen.m_ImageOrigin[0]); m_SignalGen.m_ImageOrigin[1] = v1.second.get("basic.origin.y",m_SignalGen.m_ImageOrigin[1]); m_SignalGen.m_ImageOrigin[2] = v1.second.get("basic.origin.z",m_SignalGen.m_ImageOrigin[2]); m_SignalGen.m_ImageDirection[0][0] = v1.second.get("basic.direction.1",m_SignalGen.m_ImageDirection[0][0]); m_SignalGen.m_ImageDirection[0][1] = v1.second.get("basic.direction.2",m_SignalGen.m_ImageDirection[0][1]); m_SignalGen.m_ImageDirection[0][2] = v1.second.get("basic.direction.3",m_SignalGen.m_ImageDirection[0][2]); m_SignalGen.m_ImageDirection[1][0] = v1.second.get("basic.direction.4",m_SignalGen.m_ImageDirection[1][0]); m_SignalGen.m_ImageDirection[1][1] = v1.second.get("basic.direction.5",m_SignalGen.m_ImageDirection[1][1]); m_SignalGen.m_ImageDirection[1][2] = v1.second.get("basic.direction.6",m_SignalGen.m_ImageDirection[1][2]); m_SignalGen.m_ImageDirection[2][0] = v1.second.get("basic.direction.7",m_SignalGen.m_ImageDirection[2][0]); m_SignalGen.m_ImageDirection[2][1] = v1.second.get("basic.direction.8",m_SignalGen.m_ImageDirection[2][1]); m_SignalGen.m_ImageDirection[2][2] = v1.second.get("basic.direction.9",m_SignalGen.m_ImageDirection[2][2]); m_SignalGen.m_SignalScale = v1.second.get("signalScale", m_SignalGen.m_SignalScale); m_SignalGen.m_tEcho = v1.second.get("tEcho", m_SignalGen.m_tEcho); m_SignalGen.m_tLine = v1.second.get("tLine", m_SignalGen.m_tLine); m_SignalGen.m_tInhom = v1.second.get("tInhom", m_SignalGen.m_tInhom); m_SignalGen.m_Bvalue = v1.second.get("bvalue", m_SignalGen.m_Bvalue); m_SignalGen.m_SimulateKspaceAcquisition = v1.second.get("simulatekspace", m_SignalGen.m_SimulateKspaceAcquisition); m_SignalGen.m_AxonRadius = v1.second.get("axonRadius", m_SignalGen.m_AxonRadius); switch (v1.second.get("diffusiondirectionmode", 0)) { case 0: m_SignalGen.m_DiffusionDirectionMode = SignalGenerationParameters::FIBER_TANGENT_DIRECTIONS; break; case 1: m_SignalGen.m_DiffusionDirectionMode = SignalGenerationParameters::MAIN_FIBER_DIRECTIONS; break; case 2: m_SignalGen.m_DiffusionDirectionMode = SignalGenerationParameters::RANDOM_DIRECTIONS; break; default: m_SignalGen.m_DiffusionDirectionMode = SignalGenerationParameters::FIBER_TANGENT_DIRECTIONS; } m_SignalGen.m_FiberSeparationThreshold = v1.second.get("fiberseparationthreshold", m_SignalGen.m_FiberSeparationThreshold); m_SignalGen.m_Spikes = v1.second.get("artifacts.spikesnum", m_SignalGen.m_Spikes); m_SignalGen.m_SpikeAmplitude = v1.second.get("artifacts.spikesscale", m_SignalGen.m_SpikeAmplitude); m_SignalGen.m_KspaceLineOffset = v1.second.get("artifacts.kspaceLineOffset", m_SignalGen.m_KspaceLineOffset); m_SignalGen.m_EddyStrength = v1.second.get("artifacts.eddyStrength", m_SignalGen.m_EddyStrength); m_SignalGen.m_Tau = v1.second.get("artifacts.eddyTau", m_SignalGen.m_Tau); m_SignalGen.m_CroppingFactor = v1.second.get("artifacts.aliasingfactor", m_SignalGen.m_CroppingFactor); m_SignalGen.m_DoAddGibbsRinging = v1.second.get("artifacts.addringing", m_SignalGen.m_DoAddGibbsRinging); m_SignalGen.m_DoSimulateRelaxation = v1.second.get("doSimulateRelaxation", m_SignalGen.m_DoSimulateRelaxation); m_SignalGen.m_DoDisablePartialVolume = v1.second.get("doDisablePartialVolume", m_SignalGen.m_DoDisablePartialVolume); m_SignalGen.m_DoAddMotion = v1.second.get("artifacts.doAddMotion", m_SignalGen.m_DoAddMotion); m_SignalGen.m_DoRandomizeMotion = v1.second.get("artifacts.randomMotion", m_SignalGen.m_DoRandomizeMotion); m_SignalGen.m_Translation[0] = v1.second.get("artifacts.translation0", m_SignalGen.m_Translation[0]); m_SignalGen.m_Translation[1] = v1.second.get("artifacts.translation1", m_SignalGen.m_Translation[1]); m_SignalGen.m_Translation[2] = v1.second.get("artifacts.translation2", m_SignalGen.m_Translation[2]); m_SignalGen.m_Rotation[0] = v1.second.get("artifacts.rotation0", m_SignalGen.m_Rotation[0]); m_SignalGen.m_Rotation[1] = v1.second.get("artifacts.rotation1", m_SignalGen.m_Rotation[1]); m_SignalGen.m_Rotation[2] = v1.second.get("artifacts.rotation2", m_SignalGen.m_Rotation[2]); // m_SignalGen.SetNumWeightedVolumes(v1.second.get("numgradients", m_SignalGen.GetNumWeightedVolumes())); SignalGenerationParameters::GradientListType gradients; BOOST_FOREACH( boost::property_tree::ptree::value_type const& v2, v1.second.get_child("gradients") ) { SignalGenerationParameters::GradientType g; g[0] = v2.second.get("x"); g[1] = v2.second.get("y"); g[2] = v2.second.get("z"); gradients.push_back(g); } m_SignalGen.SetGradienDirections(gradients); try { - if (v1.second.get("artifacts.noisetype")=="rice") + if (v1.second.get("artifacts.noisetype")=="rice") { m_NoiseModel = new mitk::RicianNoiseModel(); m_NoiseModel->SetNoiseVariance(v1.second.get("artifacts.noisevariance")); } } catch(...) {} try { - if (v1.second.get("artifacts.noisetype")=="chisquare") + if (v1.second.get("artifacts.noisetype")=="chisquare") { m_NoiseModel = new mitk::ChiSquareNoiseModel(); m_NoiseModel->SetNoiseVariance(v1.second.get("artifacts.noisevariance")); } } catch(...){ } BOOST_FOREACH( boost::property_tree::ptree::value_type const& v2, v1.second.get_child("compartments") ) { mitk::DiffusionSignalModel* signalModel = NULL; try // stick model { mitk::StickModel* model = new mitk::StickModel(); model->SetDiffusivity(v2.second.get("stick.d")); model->SetT2(v2.second.get("stick.t2")); model->m_CompartmentId = v2.second.get("ID"); - if (v2.second.get("type")=="fiber") + if (v2.second.get("type")=="fiber") m_FiberModelList.push_back(model); - else if (v2.second.get("type")=="non-fiber") + else if (v2.second.get("type")=="non-fiber") m_NonFiberModelList.push_back(model); signalModel = model; } catch (...) { } try // tensor model { mitk::TensorModel* model = new mitk::TensorModel(); model->SetDiffusivity1(v2.second.get("tensor.d1")); model->SetDiffusivity2(v2.second.get("tensor.d2")); model->SetDiffusivity3(v2.second.get("tensor.d3")); model->SetT2(v2.second.get("tensor.t2")); model->m_CompartmentId = v2.second.get("ID"); - if (v2.second.get("type")=="fiber") + if (v2.second.get("type")=="fiber") m_FiberModelList.push_back(model); - else if (v2.second.get("type")=="non-fiber") + else if (v2.second.get("type")=="non-fiber") m_NonFiberModelList.push_back(model); signalModel = model; } catch (...) { } try // ball model { mitk::BallModel* model = new mitk::BallModel(); model->SetDiffusivity(v2.second.get("ball.d")); model->SetT2(v2.second.get("ball.t2")); model->m_CompartmentId = v2.second.get("ID"); - if (v2.second.get("type")=="fiber") + if (v2.second.get("type")=="fiber") m_FiberModelList.push_back(model); - else if (v2.second.get("type")=="non-fiber") + else if (v2.second.get("type")=="non-fiber") m_NonFiberModelList.push_back(model); signalModel = model; } catch (...) { } try // astrosticks model { mitk::AstroStickModel* model = new AstroStickModel(); model->SetDiffusivity(v2.second.get("astrosticks.d")); model->SetT2(v2.second.get("astrosticks.t2")); model->SetRandomizeSticks(v2.second.get("astrosticks.randomize")); model->m_CompartmentId = v2.second.get("ID"); - if (v2.second.get("type")=="fiber") + if (v2.second.get("type")=="fiber") m_FiberModelList.push_back(model); - else if (v2.second.get("type")=="non-fiber") + else if (v2.second.get("type")=="non-fiber") m_NonFiberModelList.push_back(model); signalModel = model; } catch (...) { } try // dot model { mitk::DotModel* model = new mitk::DotModel(); model->SetT2(v2.second.get("dot.t2")); model->m_CompartmentId = v2.second.get("ID"); - if (v2.second.get("type")=="fiber") + if (v2.second.get("type")=="fiber") m_FiberModelList.push_back(model); - else if (v2.second.get("type")=="non-fiber") + else if (v2.second.get("type")=="non-fiber") m_NonFiberModelList.push_back(model); signalModel = model; } catch (...) { } try // prototype model { mitk::RawShModel* model = new mitk::RawShModel(); model->SetMaxNumKernels(v2.second.get("prototype.numSamples")); model->SetFaRange(v2.second.get("prototype.minFA"), v2.second.get("prototype.maxFA")); model->SetAdcRange(v2.second.get("prototype.minADC"), v2.second.get("prototype.maxADC")); model->m_CompartmentId = v2.second.get("ID"); - if (v2.second.get("type")=="fiber") + if (v2.second.get("type")=="fiber") m_FiberModelList.push_back(model); - else if (v2.second.get("type")=="non-fiber") + else if (v2.second.get("type")=="non-fiber") m_NonFiberModelList.push_back(model); signalModel = model; } catch (...) { } if (signalModel!=NULL) { signalModel->SetGradientList(gradients); try{ itk::ImageFileReader::Pointer reader = itk::ImageFileReader::New(); - reader->SetFileName(filename+"_VOLUME"+v2.second.get("ID")+".nrrd"); + reader->SetFileName(filename+"_VOLUME"+v2.second.get("ID")+".nrrd"); reader->Update(); signalModel->SetVolumeFractionImage(reader->GetOutput()); } catch(...) { } } } } } try{ itk::ImageFileReader::Pointer reader = itk::ImageFileReader::New(); reader->SetFileName(filename+"_FMAP.nrrd"); reader->Update(); m_SignalGen.m_FrequencyMap = reader->GetOutput(); } catch(...) { MITK_INFO << "No frequency map saved."; } try{ itk::ImageFileReader::Pointer reader = itk::ImageFileReader::New(); reader->SetFileName(filename+"_MASK.nrrd"); reader->Update(); m_SignalGen.m_MaskImage = reader->GetOutput(); } catch(...) { MITK_INFO << "No mask image saved."; } } template< class ScalarType > void mitk::FiberfoxParameters< ScalarType >::PrintSelf() { MITK_INFO << "Not implemented :("; } diff --git a/Modules/DiffusionImaging/FiberTracking/IODataStructures/mitkFiberfoxParameters.h b/Modules/DiffusionImaging/FiberTracking/IODataStructures/mitkFiberfoxParameters.h index 4e43573b4b..7d110c4b8f 100644 --- a/Modules/DiffusionImaging/FiberTracking/IODataStructures/mitkFiberfoxParameters.h +++ b/Modules/DiffusionImaging/FiberTracking/IODataStructures/mitkFiberfoxParameters.h @@ -1,321 +1,122 @@ /*=================================================================== 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 _MITK_FiberfoxParameters_H #define _MITK_FiberfoxParameters_H -#include -#include +#include #include #include -#include #include #include -#include -#include #include #include #include #include #include #include -using namespace std; - namespace mitk { -/** Signal generation */ -class MitkFiberTracking_EXPORT SignalGenerationParameters -{ -public: - typedef itk::Image ItkDoubleImgType; - typedef itk::Image ItkUcharImgType; - typedef itk::Vector GradientType; - typedef std::vector GradientListType; - - enum DiffusionDirectionMode { - FIBER_TANGENT_DIRECTIONS, - MAIN_FIBER_DIRECTIONS, - RANDOM_DIRECTIONS - }; - - SignalGenerationParameters() - : m_FrequencyMap(NULL) - , m_MaskImage(NULL) - , m_SignalScale(100) - , m_tEcho(100) - , m_tLine(1) - , m_tInhom(50) - , m_Bvalue(1000) - , m_SimulateKspaceAcquisition(false) - , m_AxonRadius(0) - , m_DiffusionDirectionMode(SignalGenerationParameters::FIBER_TANGENT_DIRECTIONS) - , m_FiberSeparationThreshold(30) - , m_Spikes(0) - , m_SpikeAmplitude(1) - , m_KspaceLineOffset(0) - , m_EddyStrength(0) - , m_Tau(70) - , m_CroppingFactor(1) - , m_DoAddGibbsRinging(false) - , m_DoSimulateRelaxation(true) - , m_DoDisablePartialVolume(false) - , m_DoAddMotion(false) - , m_DoRandomizeMotion(true) - { - m_ImageRegion.SetSize(0, 12); - m_ImageRegion.SetSize(1, 12); - m_ImageRegion.SetSize(2, 3); - m_ImageSpacing.Fill(2.0); - m_ImageOrigin.Fill(0.0); - m_ImageDirection.SetIdentity(); - m_Translation.Fill(0.0); - m_Rotation.Fill(0.0); - SetNumWeightedVolumes(6); - } - - /** input/output image specifications */ - itk::ImageRegion<3> m_ImageRegion; ///< Image size. - itk::Vector m_ImageSpacing; ///< Image voxel size. - itk::Point m_ImageOrigin; ///< Image origin. - itk::Matrix m_ImageDirection; ///< Image rotation matrix. - - /** Other acquisitions parameters */ - double m_SignalScale; ///< Scaling factor for output signal (before noise is added). - double m_tEcho; ///< Echo time TE. - double m_tLine; ///< k-space line readout time. - double m_tInhom; ///< T2' - double m_Bvalue; ///< Acquisition b-value - bool m_SimulateKspaceAcquisition;///< Flag to enable/disable k-space acquisition simulation - double m_AxonRadius; ///< Determines compartment volume fractions (0 == automatic axon radius estimation) - bool m_DoDisablePartialVolume; ///< Disable partial volume effects. Each voxel is either all fiber or all non-fiber. - DiffusionDirectionMode m_DiffusionDirectionMode; ///< Determines how the main diffusion direction of the signal models is selected - double m_FiberSeparationThreshold; ///< Used for random and and mein fiber deriction DiffusionDirectionMode - - /** Artifacts and other effects */ - unsigned int m_Spikes; ///< Number of spikes randomly appearing in the image - double m_SpikeAmplitude; ///< amplitude of spikes relative to the largest signal intensity (magnitude of complex) - double m_KspaceLineOffset; ///< Causes N/2 ghosts. Larger offset means stronger ghost. - double m_EddyStrength; ///< Strength of eddy current induced gradients in mT/m. - double m_Tau; ///< Eddy current decay constant (in ms) - double m_CroppingFactor; ///< FOV size in y-direction is multiplied by this factor. Causes aliasing artifacts. - bool m_DoAddGibbsRinging; ///< Add Gibbs ringing artifact - bool m_DoSimulateRelaxation; ///< Add T2 relaxation effects - bool m_DoAddMotion; ///< Enable motion artifacts. - bool m_DoRandomizeMotion; ///< Toggles between random and linear motion. - itk::Vector m_Translation; ///< Maximum translational motion. - itk::Vector m_Rotation; ///< Maximum rotational motion. - ItkDoubleImgType::Pointer m_FrequencyMap; ///< If != NULL, distortions are added to the image using this frequency map. - ItkUcharImgType::Pointer m_MaskImage; ///< Signal is only genrated inside of the mask image. - - inline void GenerateGradientHalfShell(); ///< Generates half shell of gradient directions (with m_NumGradients non-zero directions) - inline std::vector< int > GetBaselineIndices(); ///< Returns list of nun-diffusion-weighted image volume indices - inline unsigned int GetFirstBaselineIndex(); ///< Returns index of first non-diffusion-weighted image volume - inline bool IsBaselineIndex(unsigned int idx); ///< Checks if image volume with given index is non-diffusion-weighted volume or not. - inline unsigned int GetNumWeightedVolumes(); ///< Get number of diffusion-weighted image volumes - inline unsigned int GetNumBaselineVolumes(); ///< Get number of non-diffusion-weighted image volumes - inline unsigned int GetNumVolumes(); ///< Get number of baseline and diffusion-weighted image volumes - inline GradientListType GetGradientDirections(); ///< Return gradient direction container - inline GradientType GetGradientDirection(unsigned int i); - - inline void SetNumWeightedVolumes(int numGradients); ///< Automaticall calls GenerateGradientHalfShell() afterwards. - inline void SetGradienDirections(GradientListType gradientList); - inline void SetGradienDirections(mitk::DiffusionImage::GradientDirectionContainerType::Pointer gradientList); - -protected: - - unsigned int m_NumGradients; ///< Number of diffusion-weighted image volumes. - unsigned int m_NumBaseline; ///< Number of non-diffusion-weighted image volumes. - GradientListType m_GradientDirections; ///< Total number of image volumes. -}; - -/** Fiber generation */ -class MitkFiberTracking_EXPORT FiberGenerationParameters -{ -public: - - enum FiberDistribution{ - DISTRIBUTE_UNIFORM, // distribute fibers uniformly in the ROIs - DISTRIBUTE_GAUSSIAN // distribute fibers using a 2D gaussian - }; - - typedef vector< vector< mitk::PlanarEllipse::Pointer > > FiducialListType; - typedef vector< vector< unsigned int > > FlipListType; - - FiberGenerationParameters() - : m_Distribution(DISTRIBUTE_UNIFORM) - , m_Density(100) - , m_Variance(100) - , m_Sampling(1) - , m_Tension(0) - , m_Continuity(0) - , m_Bias(0) - { - m_Rotation.Fill(0.0); - m_Translation.Fill(0.0); - m_Scale.Fill(1.0); - } - - FiberDistribution m_Distribution; - unsigned int m_Density; - double m_Variance; - double m_Sampling; - double m_Tension; - double m_Continuity; - double m_Bias; - mitk::Vector3D m_Rotation; - mitk::Vector3D m_Translation; - mitk::Vector3D m_Scale; - FlipListType m_FlipList; ///< contains flags indicating a flip of the 2D fiber x-coordinates (needed to resolve some unwanted fiber twisting) - FiducialListType m_Fiducials; ///< container of the planar ellipses used as fiducials for the fiber generation process -}; - -/** GUI persistence, input, output, ... */ -class MitkFiberTracking_EXPORT MiscFiberfoxParameters -{ -public: - MiscFiberfoxParameters() - : m_ResultNode(DataNode::New()) - , m_ParentNode(NULL) - , m_SignalModelString("") - , m_ArtifactModelString("") - , m_OutputPath("") - , m_CheckOutputVolumeFractionsBox(false) - , m_CheckAdvancedSignalOptionsBox(false) - , m_CheckAddNoiseBox(false) - , m_CheckAddGhostsBox(false) - , m_CheckAddAliasingBox(false) - , m_CheckAddSpikesBox(false) - , m_CheckAddEddyCurrentsBox(false) - , m_CheckAddDistortionsBox(false) - , m_CheckRealTimeFibersBox(true) - , m_CheckAdvancedFiberOptionsBox(false) - , m_CheckConstantRadiusBox(false) - , m_CheckIncludeFiducialsBox(true) - {} - - DataNode::Pointer m_ResultNode; ///< Stores resulting image. - DataNode::Pointer m_ParentNode; ///< Parent node of result node. - string m_SignalModelString; ///< Appendet to the name of the result node - string m_ArtifactModelString; ///< Appendet to the name of the result node - string m_OutputPath; ///< Image is automatically saved to the specified folder after simulation is finished. - - /** member variables that store the check-state of GUI checkboxes */ - // image generation - bool m_CheckOutputVolumeFractionsBox; - bool m_CheckAdvancedSignalOptionsBox; - bool m_CheckAddNoiseBox; - bool m_CheckAddGhostsBox; - bool m_CheckAddAliasingBox; - bool m_CheckAddSpikesBox; - bool m_CheckAddEddyCurrentsBox; - bool m_CheckAddDistortionsBox; - // fiber generation - bool m_CheckRealTimeFibersBox; - bool m_CheckAdvancedFiberOptionsBox; - bool m_CheckConstantRadiusBox; - bool m_CheckIncludeFiducialsBox; -}; - /** * \brief Datastructure to manage the Fiberfox signal generation parameters. * */ template< class ScalarType = double > -class MitkFiberTracking_EXPORT FiberfoxParameters +class FiberfoxParameters { public: typedef itk::Image ItkDoubleImgType; typedef itk::Image ItkUcharImgType; typedef DiffusionSignalModel DiffusionModelType; typedef std::vector< DiffusionModelType* > DiffusionModelListType; typedef DiffusionNoiseModel NoiseModelType; FiberfoxParameters(); ~FiberfoxParameters(); /** Get same parameter object with different template parameter */ template< class OutType > FiberfoxParameters< OutType > CopyParameters() { FiberfoxParameters< OutType > out; out.m_FiberGen = m_FiberGen; out.m_SignalGen = m_SignalGen; out.m_Misc = m_Misc; if (m_NoiseModel!=NULL) { if (dynamic_cast*>(m_NoiseModel)) out.m_NoiseModel = new mitk::RicianNoiseModel(); else if (dynamic_cast*>(m_NoiseModel)) out.m_NoiseModel = new mitk::ChiSquareNoiseModel(); out.m_NoiseModel->SetNoiseVariance(m_NoiseModel->GetNoiseVariance()); } for (unsigned int i=0; i* outModel = NULL; mitk::DiffusionSignalModel* signalModel = NULL; if (i*>(signalModel)) outModel = new mitk::StickModel(dynamic_cast*>(signalModel)); else if (dynamic_cast*>(signalModel)) outModel = new mitk::TensorModel(dynamic_cast*>(signalModel)); else if (dynamic_cast*>(signalModel)) outModel = new mitk::RawShModel(dynamic_cast*>(signalModel)); else if (dynamic_cast*>(signalModel)) outModel = new mitk::BallModel(dynamic_cast*>(signalModel)); else if (dynamic_cast*>(signalModel)) outModel = new mitk::AstroStickModel(dynamic_cast*>(signalModel)); else if (dynamic_cast*>(signalModel)) outModel = new mitk::DotModel(dynamic_cast*>(signalModel)); if (i #include #include #include #include #include +#include +#include +#include using namespace mitk; using namespace boost::math; template< class ScalarType > RawShModel< ScalarType >::RawShModel() : m_ShOrder(0) , m_ModelIndex(-1) , m_MaxNumKernels(1000) { this->m_RandGen = itk::Statistics::MersenneTwisterRandomVariateGenerator::New(); this->m_RandGen->SetSeed(); m_AdcRange.first = 0; m_AdcRange.second = 0.004; m_FaRange.first = 0; m_FaRange.second = 1; } template< class ScalarType > RawShModel< ScalarType >::~RawShModel() { } template< class ScalarType > void RawShModel< ScalarType >::Clear() { m_ShCoefficients.clear(); m_PrototypeMaxDirection.clear(); m_B0Signal.clear(); } template< class ScalarType > void RawShModel< ScalarType >::RandomModel() { m_ModelIndex = this->m_RandGen->GetIntegerVariate(m_B0Signal.size()-1); } template< class ScalarType > unsigned int RawShModel< ScalarType >::GetNumberOfKernels() { return m_B0Signal.size(); } +template< class ScalarType > +bool RawShModel< ScalarType >::SampleKernels(DiffusionImage::Pointer diffImg, ItkUcharImageType::Pointer maskImage, TensorImageType::Pointer tensorImage, QballFilterType::CoefficientImageType::Pointer itkFeatureImage, ItkDoubleImageType::Pointer adcImage) +{ + if (diffImg.IsNull()) + return false; + + const int shOrder = 2; + + if (tensorImage.IsNull()) + { + typedef itk::DiffusionTensor3DReconstructionImageFilter< short, short, double > TensorReconstructionImageFilterType; + TensorReconstructionImageFilterType::Pointer filter = TensorReconstructionImageFilterType::New(); + filter->SetGradientImage( diffImg->GetDirections(), diffImg->GetVectorImage() ); + filter->SetBValue(diffImg->GetReferenceBValue()); + filter->Update(); + tensorImage = filter->GetOutput(); + } + + const int NumCoeffs = (shOrder*shOrder + shOrder + 2)/2 + shOrder; + if (itkFeatureImage.IsNull()) + { + QballFilterType::Pointer qballfilter = QballFilterType::New(); + qballfilter->SetGradientImage( diffImg->GetDirections(), diffImg->GetVectorImage() ); + qballfilter->SetBValue(diffImg->GetReferenceBValue()); + qballfilter->SetLambda(0.006); + qballfilter->SetNormalizationMethod(QballFilterType::QBAR_RAW_SIGNAL); + qballfilter->Update(); + itkFeatureImage = qballfilter->GetCoefficientImage(); + } + + if (adcImage.IsNull()) + { + itk::AdcImageFilter< short, double >::Pointer adcFilter = itk::AdcImageFilter< short, double >::New(); + adcFilter->SetInput(diffImg->GetVectorImage()); + adcFilter->SetGradientDirections(diffImg->GetDirections()); + adcFilter->SetB_value(diffImg->GetReferenceBValue()); + adcFilter->Update(); + adcImage = adcFilter->GetOutput(); + } + + int b0Index; + for (unsigned int i=0; iGetDirectionsWithoutMeasurementFrame()->Size(); i++) + if ( diffImg->GetDirectionsWithoutMeasurementFrame()->GetElement(i).magnitude()<0.001 ) + { + b0Index = i; + break; + } + + double max = 0; + { + itk::ImageRegionIterator< itk::VectorImage< short, 3 > > it(diffImg->GetVectorImage(), diffImg->GetVectorImage()->GetLargestPossibleRegion()); + while(!it.IsAtEnd()) + { + if (maskImage.IsNotNull() && maskImage->GetPixel(it.GetIndex())<=0) + { + ++it; + continue; + } + if (it.Get()[b0Index]>max) + max = it.Get()[b0Index]; + ++it; + } + } + + MITK_INFO << "Sampling signal kernels."; + unsigned int count = 0; + itk::ImageRegionIterator< itk::Image< itk::DiffusionTensor3D< double >, 3 > > it(tensorImage, tensorImage->GetLargestPossibleRegion()); + while(!it.IsAtEnd()) + { + bool skipPixel = false; + if (maskImage.IsNotNull() && maskImage->GetPixel(it.GetIndex())<=0) + { + ++it; + continue; + } + for (unsigned int i=0; iGetDirections()->Size(); i++) + { + if (diffImg->GetVectorImage()->GetPixel(it.GetIndex())[i]!=diffImg->GetVectorImage()->GetPixel(it.GetIndex())[i] || diffImg->GetVectorImage()->GetPixel(it.GetIndex())[i]<=0 || diffImg->GetVectorImage()->GetPixel(it.GetIndex())[i]>diffImg->GetVectorImage()->GetPixel(it.GetIndex())[b0Index]) + { + skipPixel = true; + break; + } + } + if (skipPixel) + { + ++it; + continue; + } + + typedef itk::DiffusionTensor3D TensorType; + TensorType::EigenValuesArrayType eigenvalues; + TensorType::EigenVectorsMatrixType eigenvectors; + TensorType tensor = it.Get(); + double FA = tensor.GetFractionalAnisotropy(); + double ADC = adcImage->GetPixel(it.GetIndex()); + QballFilterType::CoefficientImageType::PixelType itkv = itkFeatureImage->GetPixel(it.GetIndex()); + vnl_vector_fixed< double, NumCoeffs > coeffs; + for (unsigned int c=0; cGetMaxNumKernels()>this->GetNumberOfKernels() && FA>m_FaRange.first && FAm_AdcRange.first && ADCSetShCoefficients( coeffs, (double)diffImg->GetVectorImage()->GetPixel(it.GetIndex())[b0Index]/max )) + { + tensor.ComputeEigenAnalysis(eigenvalues, eigenvectors); + itk::Vector dir; + dir[0] = eigenvectors(2, 0); + dir[1] = eigenvectors(2, 1); + dir[2] = eigenvectors(2, 2); + m_PrototypeMaxDirection.push_back(dir); + count++; + MITK_INFO << "KERNEL: " << it.GetIndex() << " (" << count << ")"; + } + } + ++it; + } + + if (m_ShCoefficients.size()>0) + return true; + return false; +} + // convert cartesian to spherical coordinates template< class ScalarType > void RawShModel< ScalarType >::Cart2Sph( GradientListType gradients ) { m_SphCoords.set_size(gradients.size(), 2); m_SphCoords.fill(0.0); for (unsigned int i=0; i 0.0001 ) { gradients[i].Normalize(); m_SphCoords(i,0) = acos(gradients[i][2]); // theta m_SphCoords(i,1) = atan2(gradients[i][1], gradients[i][0]); // phi } } } template< class ScalarType > void RawShModel< ScalarType >::SetFiberDirection(GradientType fiberDirection) { this->m_FiberDirection = fiberDirection; this->m_FiberDirection.Normalize(); RandomModel(); GradientType axis = itk::CrossProduct(this->m_FiberDirection, m_PrototypeMaxDirection.at(m_ModelIndex)); axis.Normalize(); vnl_quaternion rotation(axis.GetVnlVector(), acos(dot_product(this->m_FiberDirection.GetVnlVector(), m_PrototypeMaxDirection.at(m_ModelIndex).GetVnlVector()))); rotation.normalize(); GradientListType gradients; for (unsigned int i=0; im_GradientList.size(); i++) { GradientType dir = this->m_GradientList.at(i); if( dir.GetNorm() > 0.0001 ) { dir.Normalize(); vnl_vector_fixed< double, 3 > vnlDir = rotation.rotate(dir.GetVnlVector()); dir[0] = vnlDir[0]; dir[1] = vnlDir[1]; dir[2] = vnlDir[2]; dir.Normalize(); } gradients.push_back(dir); } Cart2Sph( gradients ); } template< class ScalarType > bool RawShModel< ScalarType >::SetShCoefficients(vnl_vector< double > shCoefficients, double b0 ) { m_ShOrder = 2; while ( (m_ShOrder*m_ShOrder + m_ShOrder + 2)/2 + m_ShOrder <= shCoefficients.size() ) m_ShOrder += 2; m_ShOrder -= 2; m_ModelIndex = m_B0Signal.size(); m_B0Signal.push_back(b0); m_ShCoefficients.push_back(shCoefficients); -// itk::OrientationDistributionFunction odf; -// GradientListType gradients; -// for (unsigned int i=0; im_GradientList ); -// PixelType signal = SimulateMeasurement(); - -// int minDirIdx = 0; -// double min = itk::NumericTraits::max(); -// for (unsigned int i=0; ib0 || signal[i]<0) -// { -// MITK_INFO << "Corrupted signal value detected. Kernel rejected."; -// m_B0Signal.pop_back(); -// m_ShCoefficients.pop_back(); -// return false; -// } -// if (signal[i]m_GradientList.at(minDirIdx); -// maxDir.Normalize(); -// m_PrototypeMaxDirection.push_back(maxDir); + // itk::OrientationDistributionFunction odf; + // GradientListType gradients; + // for (unsigned int i=0; im_GradientList ); + // PixelType signal = SimulateMeasurement(); + + // int minDirIdx = 0; + // double min = itk::NumericTraits::max(); + // for (unsigned int i=0; ib0 || signal[i]<0) + // { + // MITK_INFO << "Corrupted signal value detected. Kernel rejected."; + // m_B0Signal.pop_back(); + // m_ShCoefficients.pop_back(); + // return false; + // } + // if (signal[i]m_GradientList.at(minDirIdx); + // maxDir.Normalize(); + // m_PrototypeMaxDirection.push_back(maxDir); Cart2Sph( this->m_GradientList ); m_ModelIndex = -1; return true; } template< class ScalarType > ScalarType RawShModel< ScalarType >::SimulateMeasurement(unsigned int dir) { ScalarType signal = 0; if (m_ModelIndex==-1) RandomModel(); if (dir>=this->m_GradientList.size()) return signal; int j, m; double mag, plm; if (this->m_GradientList[dir].GetNorm()>0.001) { j=0; for (int l=0; l<=m_ShOrder; l=l+2) for (m=-l; m<=l; m++) { plm = legendre_p(l,abs(m),cos(m_SphCoords(dir,0))); mag = sqrt((double)(2*l+1)/(4.0*M_PI)*factorial(l-abs(m))/factorial(l+abs(m)))*plm; double basis; if (m<0) basis = sqrt(2.0)*mag*cos(fabs((double)m)*m_SphCoords(dir,1)); else if (m==0) basis = mag; else basis = pow(-1.0, m)*sqrt(2.0)*mag*sin(m*m_SphCoords(dir,1)); signal += m_ShCoefficients.at(m_ModelIndex)[j]*basis; j++; } } else signal = m_B0Signal.at(m_ModelIndex); m_ModelIndex = -1; return signal; } template< class ScalarType > typename RawShModel< ScalarType >::PixelType RawShModel< ScalarType >::SimulateMeasurement() { if (m_ModelIndex==-1) RandomModel(); PixelType signal; signal.SetSize(this->m_GradientList.size()); int M = this->m_GradientList.size(); int j, m; double mag, plm; vnl_matrix< double > shBasis; shBasis.set_size(M, m_ShCoefficients.at(m_ModelIndex).size()); shBasis.fill(0.0); for (int p=0; pm_GradientList[p].GetNorm()>0.001) { j=0; for (int l=0; l<=m_ShOrder; l=l+2) for (m=-l; m<=l; m++) { plm = legendre_p(l,abs(m),cos(m_SphCoords(p,0))); mag = sqrt((double)(2*l+1)/(4.0*M_PI)*factorial(l-abs(m))/factorial(l+abs(m)))*plm; if (m<0) shBasis(p,j) = sqrt(2.0)*mag*cos(fabs((double)m)*m_SphCoords(p,1)); else if (m==0) shBasis(p,j) = mag; else shBasis(p,j) = pow(-1.0, m)*sqrt(2.0)*mag*sin(m*m_SphCoords(p,1)); j++; } double val = 0; for (unsigned int cidx=0; cidx +#include +#include namespace mitk { /** * \brief The spherical harmonic representation of a prototype diffusion weighted MR signal is used to obtain the direction dependent signal. * */ template< class ScalarType = double > class RawShModel : public DiffusionSignalModel< ScalarType > { public: RawShModel(); template< class OtherType >RawShModel(RawShModel* model) { this->m_CompartmentId = model->m_CompartmentId; this->m_T2 = model->GetT2(); this->m_FiberDirection = model->GetFiberDirection(); this->m_GradientList = model->GetGradientList(); this->m_VolumeFractionImage = model->GetVolumeFractionImage(); this->m_RandGen = model->GetRandomGenerator(); this->m_AdcRange = model->GetAdcRange(); this->m_FaRange = model->GetFaRange(); this->m_ShCoefficients = model->GetShCoefficients(); this->m_B0Signal = model->GetB0Signals(); this->m_SphCoords = model->GetSphericalCoordinates(); this->m_ShOrder = model->GetShOrder(); this->m_ModelIndex = model->GetModelIndex(); this->m_MaxNumKernels = model->GetMaxNumKernels(); } ~RawShModel(); + typedef itk::Image< double, 3 > ItkDoubleImageType; + typedef itk::Image< unsigned char, 3 > ItkUcharImageType; + typedef itk::Image< itk::DiffusionTensor3D< double >, 3 > TensorImageType; + typedef itk::AnalyticalDiffusionQballReconstructionImageFilter QballFilterType; typedef typename DiffusionSignalModel< ScalarType >::PixelType PixelType; typedef typename DiffusionSignalModel< ScalarType >::GradientType GradientType; typedef typename DiffusionSignalModel< ScalarType >::GradientListType GradientListType; typedef itk::Matrix< double, 3, 3 > MatrixType; /** Actual signal generation **/ PixelType SimulateMeasurement(); ScalarType SimulateMeasurement(unsigned int dir); bool SetShCoefficients(vnl_vector< double > shCoefficients, double b0); void SetFiberDirection(GradientType fiberDirection); void SetGradientList(GradientListType gradientList) { this->m_GradientList = gradientList; } void SetFaRange(double min, double max){ m_FaRange.first = min; m_FaRange.second = max; } void SetAdcRange(double min, double max){ m_AdcRange.first = min; m_AdcRange.second = max; } void SetMaxNumKernels(unsigned int max){ m_MaxNumKernels = max; } unsigned int GetNumberOfKernels(); std::pair< double, double > GetFaRange(){ return m_FaRange; } std::pair< double, double > GetAdcRange(){ return m_AdcRange; } unsigned int GetMaxNumKernels(){ return m_MaxNumKernels; } void Clear(); - std::vector< GradientType > m_PrototypeMaxDirection; std::vector< vnl_vector< double > > GetShCoefficients(){ return m_ShCoefficients; } std::vector< double > GetB0Signals(){ return m_B0Signal; } vnl_matrix GetSphericalCoordinates(){ return m_SphCoords; } unsigned int GetShOrder(){ return m_ShOrder; } int GetModelIndex(){ return m_ModelIndex; } double GetBaselineSignal(int index){ return m_B0Signal.at(index); } vnl_vector< double > GetCoefficients(int index){ return m_ShCoefficients.at(index); } + bool SampleKernels(DiffusionImage::Pointer diffImg, ItkUcharImageType::Pointer maskImage, TensorImageType::Pointer tensorImage=NULL, QballFilterType::CoefficientImageType::Pointer coeffImage=NULL, ItkDoubleImageType::Pointer adcImage=NULL); + protected: void Cart2Sph( GradientListType gradients ); void RandomModel(); - std::pair< double, double > m_AdcRange; - std::pair< double, double > m_FaRange; std::vector< vnl_vector< double > > m_ShCoefficients; std::vector< double > m_B0Signal; + std::vector< GradientType > m_PrototypeMaxDirection; vnl_matrix m_SphCoords; + std::pair< double, double > m_AdcRange; + std::pair< double, double > m_FaRange; unsigned int m_ShOrder; int m_ModelIndex; unsigned int m_MaxNumKernels; }; } #include "mitkRawShModel.cpp" #endif diff --git a/Modules/DiffusionImaging/FiberTracking/Testing/mitkFiberfoxSignalGenerationTest.cpp b/Modules/DiffusionImaging/FiberTracking/Testing/mitkFiberfoxSignalGenerationTest.cpp index 1f28408a0a..6e45aa8009 100644 --- a/Modules/DiffusionImaging/FiberTracking/Testing/mitkFiberfoxSignalGenerationTest.cpp +++ b/Modules/DiffusionImaging/FiberTracking/Testing/mitkFiberfoxSignalGenerationTest.cpp @@ -1,276 +1,277 @@ /*=================================================================== 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 #include #include #include #include #include #include #include #include #include /**Documentation * Test the Fiberfox simulation functions (fiberBundle -> diffusion weighted image) */ bool CompareDwi(itk::VectorImage< short, 3 >* dwi1, itk::VectorImage< short, 3 >* dwi2) { typedef itk::VectorImage< short, 3 > DwiImageType; try{ itk::ImageRegionIterator< DwiImageType > it1(dwi1, dwi1->GetLargestPossibleRegion()); itk::ImageRegionIterator< DwiImageType > it2(dwi2, dwi2->GetLargestPossibleRegion()); while(!it1.IsAtEnd()) { if (it1.Get()!=it2.Get()) return false; ++it1; ++it2; } } catch(...) { return false; } return true; } void StartSimulation(FiberfoxParameters parameters, FiberBundleX::Pointer fiberBundle, mitk::DiffusionImage::Pointer refImage, string message) { itk::TractsToDWIImageFilter< short >::Pointer tractsToDwiFilter = itk::TractsToDWIImageFilter< short >::New(); tractsToDwiFilter->SetUseConstantRandSeed(true); tractsToDwiFilter->SetParameters(parameters); tractsToDwiFilter->SetFiberBundle(fiberBundle); tractsToDwiFilter->Update(); mitk::DiffusionImage::Pointer testImage = mitk::DiffusionImage::New(); testImage->SetVectorImage( tractsToDwiFilter->GetOutput() ); testImage->SetReferenceBValue(parameters.m_SignalGen.m_Bvalue); testImage->SetDirections(parameters.m_SignalGen.GetGradientDirections()); testImage->InitializeFromVectorImage(); if (refImage.IsNotNull()) { bool cond = CompareDwi(testImage->GetVectorImage(), refImage->GetVectorImage()); if (!cond) { MITK_INFO << "Saving test and rference image to " << mitk::IOUtil::GetTempPath(); mitk::IOUtil::SaveBaseData(testImage, mitk::IOUtil::GetTempPath()+"testImage.dwi"); mitk::IOUtil::SaveBaseData(refImage, mitk::IOUtil::GetTempPath()+"refImage.dwi"); } MITK_TEST_CONDITION_REQUIRED(cond, message); } } int mitkFiberfoxSignalGenerationTest(int argc, char* argv[]) { MITK_TEST_BEGIN("mitkFiberfoxSignalGenerationTest"); MITK_TEST_CONDITION_REQUIRED(argc>=19,"check for input data"); // input fiber bundle FiberBundleXReader::Pointer fibReader = FiberBundleXReader::New(); fibReader->SetFileName(argv[1]); fibReader->Update(); FiberBundleX::Pointer fiberBundle = dynamic_cast(fibReader->GetOutput()); // reference diffusion weighted images mitk::DiffusionImage::Pointer stickBall = dynamic_cast*>(mitk::IOUtil::LoadDataNode(argv[2])->GetData()); mitk::DiffusionImage::Pointer stickAstrosticks = dynamic_cast*>(mitk::IOUtil::LoadDataNode(argv[3])->GetData()); mitk::DiffusionImage::Pointer stickDot = dynamic_cast*>(mitk::IOUtil::LoadDataNode(argv[4])->GetData()); mitk::DiffusionImage::Pointer tensorBall = dynamic_cast*>(mitk::IOUtil::LoadDataNode(argv[5])->GetData()); mitk::DiffusionImage::Pointer stickTensorBall = dynamic_cast*>(mitk::IOUtil::LoadDataNode(argv[6])->GetData()); mitk::DiffusionImage::Pointer stickTensorBallAstrosticks = dynamic_cast*>(mitk::IOUtil::LoadDataNode(argv[7])->GetData()); mitk::DiffusionImage::Pointer gibbsringing = dynamic_cast*>(mitk::IOUtil::LoadDataNode(argv[8])->GetData()); mitk::DiffusionImage::Pointer ghost = dynamic_cast*>(mitk::IOUtil::LoadDataNode(argv[9])->GetData()); mitk::DiffusionImage::Pointer aliasing = dynamic_cast*>(mitk::IOUtil::LoadDataNode(argv[10])->GetData()); mitk::DiffusionImage::Pointer eddy = dynamic_cast*>(mitk::IOUtil::LoadDataNode(argv[11])->GetData()); mitk::DiffusionImage::Pointer linearmotion = dynamic_cast*>(mitk::IOUtil::LoadDataNode(argv[12])->GetData()); mitk::DiffusionImage::Pointer randommotion = dynamic_cast*>(mitk::IOUtil::LoadDataNode(argv[13])->GetData()); mitk::DiffusionImage::Pointer spikes = dynamic_cast*>(mitk::IOUtil::LoadDataNode(argv[14])->GetData()); mitk::DiffusionImage::Pointer riciannoise = dynamic_cast*>(mitk::IOUtil::LoadDataNode(argv[15])->GetData()); mitk::DiffusionImage::Pointer chisquarenoise = dynamic_cast*>(mitk::IOUtil::LoadDataNode(argv[16])->GetData()); mitk::DiffusionImage::Pointer distortions = dynamic_cast*>(mitk::IOUtil::LoadDataNode(argv[17])->GetData()); mitk::Image::Pointer mitkFMap = dynamic_cast(mitk::IOUtil::LoadDataNode(argv[18])->GetData()); typedef itk::Image ItkDoubleImgType; ItkDoubleImgType::Pointer fMap = ItkDoubleImgType::New(); mitk::CastToItkImage(mitkFMap, fMap); FiberfoxParameters parameters; parameters.m_SignalGen.m_SimulateKspaceAcquisition = true; parameters.m_SignalGen.m_SignalScale = 10000; parameters.m_SignalGen.m_ImageRegion = stickBall->GetVectorImage()->GetLargestPossibleRegion(); parameters.m_SignalGen.m_ImageSpacing = stickBall->GetVectorImage()->GetSpacing(); parameters.m_SignalGen.m_ImageOrigin = stickBall->GetVectorImage()->GetOrigin(); parameters.m_SignalGen.m_ImageDirection = stickBall->GetVectorImage()->GetDirection(); parameters.m_SignalGen.m_Bvalue = stickBall->GetReferenceBValue(); parameters.m_SignalGen.SetGradienDirections(stickBall->GetDirections()); // intra and inter axonal compartments mitk::StickModel stickModel; stickModel.SetBvalue(parameters.m_SignalGen.m_Bvalue); stickModel.SetT2(110); stickModel.SetDiffusivity(0.001); stickModel.SetGradientList(parameters.m_SignalGen.GetGradientDirections()); mitk::TensorModel tensorModel; tensorModel.SetT2(110); stickModel.SetBvalue(parameters.m_SignalGen.m_Bvalue); tensorModel.SetDiffusivity1(0.001); tensorModel.SetDiffusivity2(0.00025); tensorModel.SetDiffusivity3(0.00025); tensorModel.SetGradientList(parameters.m_SignalGen.GetGradientDirections()); // extra axonal compartment models mitk::BallModel ballModel; ballModel.SetT2(80); ballModel.SetBvalue(parameters.m_SignalGen.m_Bvalue); ballModel.SetDiffusivity(0.001); ballModel.SetGradientList(parameters.m_SignalGen.GetGradientDirections()); mitk::AstroStickModel astrosticksModel; astrosticksModel.SetT2(80); astrosticksModel.SetBvalue(parameters.m_SignalGen.m_Bvalue); astrosticksModel.SetDiffusivity(0.001); astrosticksModel.SetRandomizeSticks(true); astrosticksModel.SetSeed(0); astrosticksModel.SetGradientList(parameters.m_SignalGen.GetGradientDirections()); mitk::DotModel dotModel; dotModel.SetT2(80); dotModel.SetGradientList(parameters.m_SignalGen.GetGradientDirections()); // noise models mitk::RicianNoiseModel* ricianNoiseModel = new mitk::RicianNoiseModel(); ricianNoiseModel->SetNoiseVariance(1000000); ricianNoiseModel->SetSeed(0); // Rician noise mitk::ChiSquareNoiseModel* chiSquareNoiseModel = new mitk::ChiSquareNoiseModel(); chiSquareNoiseModel->SetNoiseVariance(1000000); chiSquareNoiseModel->SetSeed(0); try{ // Stick-Ball parameters.m_FiberModelList.push_back(&stickModel); parameters.m_NonFiberModelList.push_back(&ballModel); StartSimulation(parameters, fiberBundle, stickBall, argv[2]); // Srick-Astrosticks parameters.m_NonFiberModelList.clear(); parameters.m_NonFiberModelList.push_back(&astrosticksModel); StartSimulation(parameters, fiberBundle, stickAstrosticks, argv[3]); // Stick-Dot parameters.m_NonFiberModelList.clear(); parameters.m_NonFiberModelList.push_back(&dotModel); StartSimulation(parameters, fiberBundle, stickDot, argv[4]); // Tensor-Ball parameters.m_FiberModelList.clear(); parameters.m_FiberModelList.push_back(&tensorModel); parameters.m_NonFiberModelList.clear(); parameters.m_NonFiberModelList.push_back(&ballModel); StartSimulation(parameters, fiberBundle, tensorBall, argv[5]); // Stick-Tensor-Ball parameters.m_FiberModelList.clear(); parameters.m_FiberModelList.push_back(&stickModel); parameters.m_FiberModelList.push_back(&tensorModel); parameters.m_NonFiberModelList.clear(); parameters.m_NonFiberModelList.push_back(&ballModel); StartSimulation(parameters, fiberBundle, stickTensorBall, argv[6]); // Stick-Tensor-Ball-Astrosticks parameters.m_NonFiberModelList.push_back(&astrosticksModel); StartSimulation(parameters, fiberBundle, stickTensorBallAstrosticks, argv[7]); // Gibbs ringing parameters.m_FiberModelList.clear(); parameters.m_FiberModelList.push_back(&stickModel); parameters.m_NonFiberModelList.clear(); parameters.m_NonFiberModelList.push_back(&ballModel); parameters.m_SignalGen.m_DoAddGibbsRinging = true; StartSimulation(parameters, fiberBundle, gibbsringing, argv[8]); // Ghost parameters.m_SignalGen.m_DoAddGibbsRinging = false; parameters.m_SignalGen.m_KspaceLineOffset = 0.25; StartSimulation(parameters, fiberBundle, ghost, argv[9]); // Aliasing parameters.m_SignalGen.m_KspaceLineOffset = 0; parameters.m_SignalGen.m_CroppingFactor = 0.4; parameters.m_SignalGen.m_SignalScale = 1000; StartSimulation(parameters, fiberBundle, aliasing, argv[10]); // Eddy currents parameters.m_SignalGen.m_CroppingFactor = 1; parameters.m_SignalGen.m_SignalScale = 10000; parameters.m_SignalGen.m_EddyStrength = 0.05; StartSimulation(parameters, fiberBundle, eddy, argv[11]); // Motion (linear) parameters.m_SignalGen.m_EddyStrength = 0.0; parameters.m_SignalGen.m_DoAddMotion = true; parameters.m_SignalGen.m_DoRandomizeMotion = false; parameters.m_SignalGen.m_Translation[1] = 10; parameters.m_SignalGen.m_Rotation[2] = 90; StartSimulation(parameters, fiberBundle, linearmotion, argv[12]); // Motion (random) parameters.m_SignalGen.m_DoRandomizeMotion = true; parameters.m_SignalGen.m_Translation[1] = 5; parameters.m_SignalGen.m_Rotation[2] = 45; StartSimulation(parameters, fiberBundle, randommotion, argv[13]); // Spikes parameters.m_SignalGen.m_DoAddMotion = false; parameters.m_SignalGen.m_Spikes = 5; parameters.m_SignalGen.m_SpikeAmplitude = 1; StartSimulation(parameters, fiberBundle, spikes, argv[14]); // Rician noise parameters.m_SignalGen.m_Spikes = 0; parameters.m_NoiseModel = ricianNoiseModel; StartSimulation(parameters, fiberBundle, riciannoise, argv[15]); delete parameters.m_NoiseModel; // Chi-square noise parameters.m_NoiseModel = chiSquareNoiseModel; StartSimulation(parameters, fiberBundle, chisquarenoise, argv[16]); delete parameters.m_NoiseModel; // Distortions parameters.m_NoiseModel = NULL; parameters.m_SignalGen.m_FrequencyMap = fMap; StartSimulation(parameters, fiberBundle, distortions, argv[17]); } catch (std::exception &e) { MITK_TEST_CONDITION_REQUIRED(false, e.what()); } // always end with this! MITK_TEST_END(); } diff --git a/Modules/DiffusionImaging/FiberTracking/files.cmake b/Modules/DiffusionImaging/FiberTracking/files.cmake index c0f05f474b..456eca1f99 100644 --- a/Modules/DiffusionImaging/FiberTracking/files.cmake +++ b/Modules/DiffusionImaging/FiberTracking/files.cmake @@ -1,79 +1,81 @@ set(CPP_FILES ## IO datastructures IODataStructures/FiberBundleX/mitkFiberBundleX.cpp IODataStructures/FiberBundleX/mitkFiberBundleXWriter.cpp IODataStructures/FiberBundleX/mitkFiberBundleXReader.cpp IODataStructures/FiberBundleX/mitkFiberBundleXSerializer.cpp IODataStructures/FiberBundleX/mitkTrackvis.cpp IODataStructures/PlanarFigureComposite/mitkPlanarFigureComposite.cpp + IODataStructures/mitkFiberfoxNoTemplateParameters.cpp # Interactions Interactions/mitkFiberBundleInteractor.cpp # Tractography Algorithms/GibbsTracking/mitkParticleGrid.cpp Algorithms/GibbsTracking/mitkMetropolisHastingsSampler.cpp Algorithms/GibbsTracking/mitkEnergyComputer.cpp Algorithms/GibbsTracking/mitkGibbsEnergyComputer.cpp Algorithms/GibbsTracking/mitkFiberBuilder.cpp Algorithms/GibbsTracking/mitkSphereInterpolator.cpp ) set(H_FILES # DataStructures -> FiberBundleX IODataStructures/FiberBundleX/mitkFiberBundleX.h IODataStructures/FiberBundleX/mitkFiberBundleXWriter.h IODataStructures/FiberBundleX/mitkFiberBundleXReader.h IODataStructures/FiberBundleX/mitkFiberBundleXSerializer.h IODataStructures/FiberBundleX/mitkTrackvis.h IODataStructures/mitkFiberfoxParameters.h + IODataStructures/mitkFiberfoxNoTemplateParameters.h # Algorithms Algorithms/itkTractDensityImageFilter.h Algorithms/itkTractsToFiberEndingsImageFilter.h Algorithms/itkTractsToRgbaImageFilter.h Algorithms/itkElectrostaticRepulsionDiffusionGradientReductionFilter.h Algorithms/itkFibersFromPlanarFiguresFilter.h Algorithms/itkTractsToDWIImageFilter.h Algorithms/itkTractsToVectorImageFilter.h Algorithms/itkKspaceImageFilter.h Algorithms/itkDftImageFilter.h Algorithms/itkAddArtifactsToDwiImageFilter.h Algorithms/itkFieldmapGeneratorFilter.h Algorithms/itkEvaluateDirectionImagesFilter.h Algorithms/itkEvaluateTractogramDirectionsFilter.h # (old) Tractography Algorithms/itkGibbsTrackingFilter.h Algorithms/itkStochasticTractographyFilter.h Algorithms/itkStreamlineTrackingFilter.h Algorithms/GibbsTracking/mitkParticle.h Algorithms/GibbsTracking/mitkParticleGrid.h Algorithms/GibbsTracking/mitkMetropolisHastingsSampler.h Algorithms/GibbsTracking/mitkSimpSamp.h Algorithms/GibbsTracking/mitkEnergyComputer.h Algorithms/GibbsTracking/mitkGibbsEnergyComputer.h Algorithms/GibbsTracking/mitkSphereInterpolator.h Algorithms/GibbsTracking/mitkFiberBuilder.h # Signal Models SignalModels/mitkDiffusionSignalModel.h SignalModels/mitkTensorModel.h SignalModels/mitkBallModel.h SignalModels/mitkDotModel.h SignalModels/mitkAstroStickModel.h SignalModels/mitkStickModel.h SignalModels/mitkRawShModel.h SignalModels/mitkDiffusionNoiseModel.h SignalModels/mitkRicianNoiseModel.h SignalModels/mitkChiSquareNoiseModel.h ) set(RESOURCE_FILES # Binary directory resources FiberTrackingLUTBaryCoords.bin FiberTrackingLUTIndices.bin # Shaders Shaders/mitkShaderFiberClipping.xml )