diff --git a/Modules/DiffusionCmdApps/FiberQuantification/CMakeLists.txt b/Modules/DiffusionCmdApps/FiberQuantification/CMakeLists.txt index 235e276..3782c75 100644 --- a/Modules/DiffusionCmdApps/FiberQuantification/CMakeLists.txt +++ b/Modules/DiffusionCmdApps/FiberQuantification/CMakeLists.txt @@ -1,43 +1,40 @@ option(BUILD_DiffusionFiberQuantificationCmdApps "Build commandline tools for diffusion fiber processing" OFF) if(BUILD_DiffusionFiberQuantificationCmdApps OR MITK_BUILD_ALL_APPS) # needed include directories include_directories( ${CMAKE_CURRENT_SOURCE_DIR} ${CMAKE_CURRENT_BINARY_DIR} ) # list of diffusion cmdapps # if an app requires additional dependencies # they are added after a "^^" and separated by "_" set( diffusionFiberProcessingcmdapps TractDensity^^ PrintFiberStatistics^^ DistanceToSegmentation^^ - TractParcellation^^ - TractometryComparison^^ - EstimateNumSamplesForTractometry^^ ) foreach(diffusionFiberProcessingcmdapp ${diffusionFiberProcessingcmdapps}) # extract cmd app name and dependencies string(REPLACE "^^" "\\;" cmdapp_info ${diffusionFiberProcessingcmdapp}) set(cmdapp_info_list ${cmdapp_info}) list(GET cmdapp_info_list 0 appname) list(GET cmdapp_info_list 1 raw_dependencies) string(REPLACE "_" "\\;" dependencies "${raw_dependencies}") set(dependencies_list ${dependencies}) mitkFunctionCreateCommandLineApp( NAME ${appname} DEPENDS MitkDiffusionCmdApps MitkMriSimulation MitkFiberTracking ${dependencies_list} PACKAGE_DEPENDS ) endforeach() if(EXECUTABLE_IS_ENABLED) MITK_INSTALL_TARGETS(EXECUTABLES ${EXECUTABLE_TARGET}) endif() endif() diff --git a/Modules/DiffusionCmdApps/FiberQuantification/EstimateNumSamplesForTractometry.cpp b/Modules/DiffusionCmdApps/FiberQuantification/EstimateNumSamplesForTractometry.cpp deleted file mode 100644 index cbe48fc..0000000 --- a/Modules/DiffusionCmdApps/FiberQuantification/EstimateNumSamplesForTractometry.cpp +++ /dev/null @@ -1,111 +0,0 @@ -/*=================================================================== - -The Medical Imaging Interaction Toolkit (MITK) - -Copyright (c) German Cancer Research Center. - -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 "mitkDiffusionCommandLineParser.h" -#include -#include -#include -#include -#include - - - -/*! -\brief Modify input tractogram: fiber resampling, compression, pruning and transformation. -*/ -int main(int argc, char* argv[]) -{ - mitkDiffusionCommandLineParser parser; - - parser.setTitle("Tract Density"); - parser.setCategory("Fiber Quantification Methods"); - parser.setDescription("Estimate minimal number of sampling points for along tract-radiomics based on the input image spacing and the tracts. Each cell should span at least N voxels."); - parser.setContributor("MIC"); - - parser.setArgumentPrefix("--", "-"); - parser.addArgument("tracts", "", mitkDiffusionCommandLineParser::StringList, "Input tracts:", "input tracts", us::Any(), false); - parser.addArgument("images", "", mitkDiffusionCommandLineParser::StringList, "Input images:", "input images", us::Any(), false); - parser.addArgument("voxels", "", mitkDiffusionCommandLineParser::Int, "Number of voxels per parcel:", "", 3); - parser.addArgument("out_file", "", mitkDiffusionCommandLineParser::String, "Output File:", "output file", us::Any(), false); - - std::map parsedArgs = parser.parseArguments(argc, argv); - if (parsedArgs.size()==0) - return EXIT_FAILURE; - - mitkDiffusionCommandLineParser::StringContainerType inFibs = us::any_cast(parsedArgs["tracts"]); - mitkDiffusionCommandLineParser::StringContainerType inImages = us::any_cast(parsedArgs["images"]); - std::string outfile = us::any_cast(parsedArgs["out_file"]); - - int voxels = 3; - if (parsedArgs.count("voxels")) - voxels = us::any_cast(parsedArgs["voxels"]); - - if (inFibs.size()!=inImages.size()) - { - MITK_INFO << "Equal number of tracs and images required!"; - return EXIT_FAILURE; - } - - try - { - std::vector< std::string > fib_names, img_names; - auto input_tracts = mitk::DiffusionDataIOHelper::load_fibs(inFibs, &fib_names); - auto input_images = mitk::DiffusionDataIOHelper::load_mitk_images(inImages, &img_names); - - std::ofstream statistics_file; - statistics_file.open (outfile); - statistics_file << "tract_file;image_file;num_samples" << std::endl; - for (unsigned int i=0; i::Pointer itkImage = itk::Image::New(); - CastToItkImage(input_images.at(i), itkImage); - - auto v = mitk::Tractometry::EstimateNumSamplingPoints(itkImage, input_tracts.at(i), voxels); - statistics_file << fib_names.at(i) << ";" << img_names.at(i) << ";" << boost::lexical_cast(v) << std::endl; - } - - statistics_file.close(); - } - catch (const itk::ExceptionObject& e) - { - std::cout << e.what(); - return EXIT_FAILURE; - } - catch (std::exception& e) - { - std::cout << e.what(); - return EXIT_FAILURE; - } - catch (...) - { - std::cout << "ERROR!?!"; - return EXIT_FAILURE; - } - return EXIT_SUCCESS; -} diff --git a/Modules/DiffusionCmdApps/FiberQuantification/TractParcellation.cpp b/Modules/DiffusionCmdApps/FiberQuantification/TractParcellation.cpp deleted file mode 100644 index 20c7743..0000000 --- a/Modules/DiffusionCmdApps/FiberQuantification/TractParcellation.cpp +++ /dev/null @@ -1,162 +0,0 @@ -/*=================================================================== - -The Medical Imaging Interaction Toolkit (MITK) - -Copyright (c) German Cancer Research Center. - -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 "mitkDiffusionCommandLineParser.h" -#include -#include -#include -#include -#include - - -mitk::FiberBundle::Pointer LoadFib(std::string filename) -{ - std::vector fibInfile = mitk::IOUtil::Load(filename); - if( fibInfile.empty() ) - std::cout << "File " << filename << " could not be read!"; - mitk::BaseData::Pointer baseData = fibInfile.at(0); - return dynamic_cast(baseData.GetPointer()); -} - -/*! -\brief Modify input tractogram: fiber resampling, compression, pruning and transformation. -*/ -int main(int argc, char* argv[]) -{ - mitkDiffusionCommandLineParser parser; - - parser.setTitle("Tract Parcellation"); - parser.setCategory("Fiber Quantification Methods"); - parser.setDescription("Parcellate tract for along-tract radiomics."); - parser.setContributor("MIC"); - - parser.setArgumentPrefix("--", "-"); - parser.addArgument("", "i", mitkDiffusionCommandLineParser::String, "Input:", "input fiber bundle", us::Any(), false); - parser.addArgument("binary_output", "", mitkDiffusionCommandLineParser::String, "Binary Output:", "name of binary output parcels", us::Any()); - parser.addArgument("full_output", "", mitkDiffusionCommandLineParser::String, "Full Output:", "name of full parcellation image", us::Any()); - parser.addArgument("reference_image", "", mitkDiffusionCommandLineParser::String, "Reference image:", "output image will have geometry of this reference image", us::Any(), false); - - parser.addArgument("num_parcels", "", mitkDiffusionCommandLineParser::Int, "Number of parcels:", "", 15); - parser.addArgument("num_centroids", "", mitkDiffusionCommandLineParser::Int, "Number of centroids:", "", 0); - parser.addArgument("start_cluster_size", "", mitkDiffusionCommandLineParser::Float, "Cluster size (in mm):", "", 5.0); - - std::map parsedArgs = parser.parseArguments(argc, argv); - if (parsedArgs.size()==0) - return EXIT_FAILURE; - - std::string binary_output = ""; - if (parsedArgs.count("binary_output")) - binary_output = us::any_cast(parsedArgs["binary_output"]); - - std::string full_output = ""; - if (parsedArgs.count("full_output")) - full_output = us::any_cast(parsedArgs["full_output"]); - - int num_parcels = 15; - if (parsedArgs.count("num_parcels")) - num_parcels = us::any_cast(parsedArgs["num_parcels"]); - - int num_centroids = 0; - if (parsedArgs.count("num_centroids")) - num_centroids = us::any_cast(parsedArgs["num_centroids"]); - - float start_cluster_size = 5; - if (parsedArgs.count("start_cluster_size")) - start_cluster_size = us::any_cast(parsedArgs["start_cluster_size"]); - - std::string reference_image = us::any_cast(parsedArgs["reference_image"]); - std::string inFileName = us::any_cast(parsedArgs["i"]); - - try - { - mitk::FiberBundle::Pointer fib = LoadFib(inFileName); - mitk::Image::Pointer ref_img = mitk::IOUtil::Load(reference_image); - - typedef unsigned char OutPixType; - typedef itk::Image OutImageType; - - OutImageType::Pointer itkImage = OutImageType::New(); - CastToItkImage(ref_img, itkImage); - - itk::TractParcellationFilter< >::Pointer parcellator = itk::TractParcellationFilter< >::New(); - parcellator->SetInputImage(itkImage); - parcellator->SetNumParcels(num_parcels); - parcellator->SetInputTract(fib); - parcellator->SetNumCentroids(num_centroids); - parcellator->SetStartClusterSize(start_cluster_size); - parcellator->Update(); - OutImageType::Pointer out_image = parcellator->GetOutput(0); - OutImageType::Pointer out_image_pp = parcellator->GetOutput(1); - - if (!binary_output.empty()) - { - auto binary_segments = parcellator->GetBinarySplit(out_image_pp); - - int c=1; - for (auto itk_segment : binary_segments) - { - mitk::Image::Pointer mitk_segment = mitk::Image::New(); - mitk_segment->InitializeByItk(itk_segment.GetPointer()); - mitk_segment->SetVolume(itk_segment->GetBufferPointer()); - if (c<10) - mitk::IOUtil::Save(mitk_segment, binary_output + "_0" + boost::lexical_cast(c) + ".nrrd"); - else - mitk::IOUtil::Save(mitk_segment, binary_output + "_" + boost::lexical_cast(c) + ".nrrd"); - ++c; - } - } - if (binary_output.empty() && full_output.empty()) - { - full_output = "parcellation.nrrd"; - } - if (!full_output.empty()) - { - mitk::Image::Pointer mitk_img_pp = mitk::Image::New(); - mitk_img_pp->InitializeByItk(out_image_pp.GetPointer()); - mitk_img_pp->SetVolume(out_image_pp->GetBufferPointer()); - mitk::IOUtil::Save(mitk_img_pp, full_output ); - } - } - catch (const itk::ExceptionObject& e) - { - std::cout << e.what(); - return EXIT_FAILURE; - } - catch (std::exception& e) - { - std::cout << e.what(); - return EXIT_FAILURE; - } - catch (...) - { - std::cout << "ERROR!?!"; - return EXIT_FAILURE; - } - return EXIT_SUCCESS; -} diff --git a/Modules/DiffusionCmdApps/FiberQuantification/TractometryComparison.cpp b/Modules/DiffusionCmdApps/FiberQuantification/TractometryComparison.cpp deleted file mode 100644 index fe09690..0000000 --- a/Modules/DiffusionCmdApps/FiberQuantification/TractometryComparison.cpp +++ /dev/null @@ -1,137 +0,0 @@ -/*=================================================================== - -The Medical Imaging Interaction Toolkit (MITK) - -Copyright (c) German Cancer Research Center. - -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 "mitkDiffusionCommandLineParser.h" -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include - -#define _USE_MATH_DEFINES -#include - -/*! -\brief Join multiple tractograms -*/ -int main(int argc, char* argv[]) -{ - mitkDiffusionCommandLineParser parser; - - parser.setTitle("Fiber Statistics"); - parser.setCategory("Fiber Quantification Methods"); - parser.setContributor("MIC"); - parser.setDescription("Calculate mutliple tractometry methods."); - - parser.setArgumentPrefix("--", "-"); - parser.addArgument("tracts", "", mitkDiffusionCommandLineParser::StringList, "Input tracts:", "input tracts", us::Any(), false); - parser.addArgument("images", "", mitkDiffusionCommandLineParser::StringList, "Input images:", "input images", us::Any(), false); - parser.addArgument("out_file", "", mitkDiffusionCommandLineParser::String, "Output File:", "output file", us::Any(), false); - parser.addArgument("num_parcels", "", mitkDiffusionCommandLineParser::Int, "Number of parcels:", "", 15); - - std::map parsedArgs = parser.parseArguments(argc, argv); - if (parsedArgs.size()==0) - return EXIT_FAILURE; - - mitkDiffusionCommandLineParser::StringContainerType inFibs = us::any_cast(parsedArgs["tracts"]); - mitkDiffusionCommandLineParser::StringContainerType inImages = us::any_cast(parsedArgs["images"]); - std::string outfile = us::any_cast(parsedArgs["out_file"]); - - int num_parcels = 15; - if (parsedArgs.count("num_parcels")) - num_parcels = us::any_cast(parsedArgs["num_parcels"]); - - try - { - - std::vector< std::string > fib_names, img_names; - auto input_tracts = mitk::DiffusionDataIOHelper::load_fibs(inFibs, &fib_names); - auto input_images = mitk::DiffusionDataIOHelper::load_itk_images>(inImages, &img_names); - - bool add_header = !itksys::SystemTools::FileExists(outfile, true); - - - std::ofstream statistics_file; - statistics_file.open (outfile, std::ios_base::app); - if (add_header) - statistics_file << "method;tract_file;image_file;parcel;value" << std::endl; - - auto ref = input_tracts.at(0)->GetDeepCopy(); - for (unsigned int i=0; iResampleSpline(1.0); - - { - auto output = mitk::Tractometry::StaticResamplingTractometry(input_images.at(i), fib->GetDeepCopy(), num_parcels, ref); - for (unsigned int r=0; r(r) << ";" << boost::lexical_cast(output.get(r,c)) << std::endl; - } - - { - auto output = mitk::Tractometry::NearestCentroidPointTractometry(input_images.at(i), fib->GetDeepCopy(), num_parcels, 1, 99, ref); - for (unsigned int r=0; r(r) << ";" << boost::lexical_cast(output.at(r).get(c)) << std::endl; - } - { - itk::TractParcellationFilter< itk::Image, itk::Image >::Pointer parcellator = itk::TractParcellationFilter< itk::Image, itk::Image >::New(); - parcellator->SetInputImage(input_images.at(i)); - parcellator->SetNumParcels(num_parcels); - parcellator->SetInputTract(fib->GetDeepCopy()); - parcellator->SetReferenceTract(ref); - parcellator->Update(); - itk::Image::Pointer out_image_pp = parcellator->GetOutput(1); - - itk::ImageRegionConstIterator< itk::Image > it(input_images.at(i), input_images.at(i)->GetLargestPossibleRegion()); - itk::ImageRegionConstIterator< itk::Image > it2(out_image_pp, out_image_pp->GetLargestPossibleRegion()); - while (!it.IsAtEnd()) { - if (it2.Get()>0) - statistics_file << "TractParcellationFilter;" << fib_names.at(i) << ";" << img_names.at(i) << ";" << boost::lexical_cast(it2.Get()-1) << ";" << boost::lexical_cast(it.Get()) << std::endl; - ++it; - ++it2; - } - } - } - - statistics_file.close(); - } - catch (const itk::ExceptionObject& e) - { - std::cout << e.what(); - return EXIT_FAILURE; - } - catch (std::exception& e) - { - std::cout << e.what(); - return EXIT_FAILURE; - } - catch (...) - { - std::cout << "ERROR!?!"; - return EXIT_FAILURE; - } - return EXIT_SUCCESS; -} diff --git a/Modules/FiberTracking/Algorithms/itkTractDensityImageFilter.cpp b/Modules/FiberTracking/Algorithms/itkTractDensityImageFilter.cpp index 40ccff0..7c6094a 100644 --- a/Modules/FiberTracking/Algorithms/itkTractDensityImageFilter.cpp +++ b/Modules/FiberTracking/Algorithms/itkTractDensityImageFilter.cpp @@ -1,216 +1,209 @@ /*=================================================================== The Medical Imaging Interaction Toolkit (MITK) Coindex[1]right (c) German Cancer Research Center, 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 "itkTractDensityImageFilter.h" // VTK #include #include #include #include // misc #include #include #include #include namespace itk{ template< class OutputImageType, class RefImageType > TractDensityImageFilter< OutputImageType, RefImageType >::TractDensityImageFilter() : m_UpsamplingFactor(1) , m_InvertImage(false) , m_Mode(DENSITY) , m_UseImageGeometry(false) , m_OutputAbsoluteValues(false) , m_MaxDensity(0) , m_NumCoveredVoxels(0) { } template< class OutputImageType, class RefImageType > TractDensityImageFilter< OutputImageType, RefImageType >::~TractDensityImageFilter() { } template< class OutputImageType, class RefImageType > TDI_MODE TractDensityImageFilter< OutputImageType, RefImageType >::GetMode() const { return m_Mode; } template< class OutputImageType, class RefImageType > void TractDensityImageFilter< OutputImageType, RefImageType >::SetMode(const TDI_MODE &Mode) { m_Mode = Mode; } template< class OutputImageType, class RefImageType > void TractDensityImageFilter< OutputImageType, RefImageType >::GenerateData() { // generate upsampled image mitk::BaseGeometry::Pointer geometry = m_FiberBundle->GetGeometry(); typename OutputImageType::Pointer outImage = this->GetOutput(); // calculate new image parameters itk::Vector newSpacing; mitk::Point3D newOrigin; itk::Matrix newDirection; ImageRegion<3> upsampledRegion; if (m_UseImageGeometry && !m_InputImage.IsNull()) { MITK_INFO << "TractDensityImageFilter: using image geometry"; newSpacing = m_InputImage->GetSpacing()/m_UpsamplingFactor; upsampledRegion = m_InputImage->GetLargestPossibleRegion(); newOrigin = m_InputImage->GetOrigin(); typename OutputImageType::RegionType::SizeType size = upsampledRegion.GetSize(); size[0] *= m_UpsamplingFactor; size[1] *= m_UpsamplingFactor; size[2] *= m_UpsamplingFactor; upsampledRegion.SetSize(size); newDirection = m_InputImage->GetDirection(); } else { MITK_INFO << "TractDensityImageFilter: using fiber bundle geometry"; newSpacing = geometry->GetSpacing()/m_UpsamplingFactor; newOrigin = geometry->GetOrigin(); mitk::Geometry3D::BoundsArrayType bounds = geometry->GetBounds(); // we retrieve the origin from a vtk-polydata (corner-based) and hance have to translate it to an image geometry // i.e. center-based newOrigin[0] += bounds.GetElement(0) + 0.5 * newSpacing[0]; newOrigin[1] += bounds.GetElement(2) + 0.5 * newSpacing[1]; newOrigin[2] += bounds.GetElement(4) + 0.5 * newSpacing[2]; for (int i=0; i<3; i++) for (int j=0; j<3; j++) newDirection[j][i] = geometry->GetMatrixColumn(i)[j]; upsampledRegion.SetSize(0, ceil( geometry->GetExtent(0)*m_UpsamplingFactor ) ); upsampledRegion.SetSize(1, ceil( geometry->GetExtent(1)*m_UpsamplingFactor ) ); upsampledRegion.SetSize(2, ceil( geometry->GetExtent(2)*m_UpsamplingFactor ) ); } typename OutputImageType::RegionType::SizeType upsampledSize = upsampledRegion.GetSize(); // apply new image parameters outImage->SetSpacing( newSpacing ); outImage->SetOrigin( newOrigin ); outImage->SetDirection( newDirection ); outImage->SetLargestPossibleRegion( upsampledRegion ); outImage->SetBufferedRegion( upsampledRegion ); outImage->SetRequestedRegion( upsampledRegion ); outImage->Allocate(); outImage->FillBuffer(0.0); int w = upsampledSize[0]; int h = upsampledSize[1]; int d = upsampledSize[2]; // set/initialize output OutPixelType* outImageBufferPointer = (OutPixelType*)outImage->GetBufferPointer(); MITK_INFO << "TractDensityImageFilter: starting image generation"; vtkSmartPointer fiberPolyData = m_FiberBundle->GetFiberPolyData(); - m_AverageNumTraversedVoxels = 0; - m_AverageSegmentLength = 0; int numFibers = m_FiberBundle->GetNumFibers(); boost::progress_display disp(numFibers); for( int i=0; iGetCell(i); int numPoints = cell->GetNumberOfPoints(); vtkPoints* points = cell->GetPoints(); float weight = m_FiberBundle->GetFiberWeight(i); // fill output image for( int j=0; j startVertex = mitk::imv::GetItkPoint(points->GetPoint(j)); itk::Index<3> startIndex; itk::ContinuousIndex startIndexCont; outImage->TransformPhysicalPointToIndex(startVertex, startIndex); outImage->TransformPhysicalPointToContinuousIndex(startVertex, startIndexCont); itk::Point endVertex = mitk::imv::GetItkPoint(points->GetPoint(j + 1)); itk::Index<3> endIndex; itk::ContinuousIndex endIndexCont; outImage->TransformPhysicalPointToIndex(endVertex, endIndex); outImage->TransformPhysicalPointToContinuousIndex(endVertex, endIndexCont); std::vector< std::pair< itk::Index<3>, double > > segments = mitk::imv::IntersectImage(newSpacing, startIndex, endIndex, startIndexCont, endIndexCont); - m_AverageNumTraversedVoxels += segments.size(); for (std::pair< itk::Index<3>, double > segment : segments) { - m_AverageSegmentLength += segment.second; if (!outImage->GetLargestPossibleRegion().IsInside(segment.first)) continue; if (outImage->GetPixel(segment.first)==0) m_NumCoveredVoxels++; switch (m_Mode) { case BINARY: { outImage->SetPixel(segment.first, 1); break; } case VISITATION_COUNT: { outImage->SetPixel(segment.first, outImage->GetPixel(segment.first) + 1); break; } case DENSITY: { outImage->SetPixel(segment.first, outImage->GetPixel(segment.first)+segment.second * weight); break; } default: { outImage->SetPixel(segment.first, outImage->GetPixel(segment.first)+segment.second * weight); } } } } } - m_AverageSegmentLength /= m_AverageNumTraversedVoxels; - m_AverageNumTraversedVoxels = m_AverageNumTraversedVoxels/numFibers; - m_MaxDensity = 0; for (int i=0; i0) for (int i=0; i #include #include #include #include enum TDI_MODE : int { BINARY, VISITATION_COUNT, DENSITY }; namespace itk{ /** * \brief Generates tract density images from input fiberbundles (Calamante 2010). */ template< class OutputImageType, class RefImageType=OutputImageType > class TractDensityImageFilter : public ImageSource< OutputImageType > { public: typedef TractDensityImageFilter Self; typedef ProcessObject Superclass; typedef SmartPointer< Self > Pointer; typedef SmartPointer< const Self > ConstPointer; typedef typename OutputImageType::PixelType OutPixelType; itkFactorylessNewMacro(Self) itkCloneMacro(Self) itkTypeMacro( TractDensityImageFilter, ImageSource ) itkSetMacro( UpsamplingFactor, float) ///< use higher resolution for ouput image itkGetMacro( UpsamplingFactor, float) ///< use higher resolution for ouput image itkSetMacro( InvertImage, bool) ///< voxelvalue = 1-voxelvalue itkGetMacro( InvertImage, bool) ///< voxelvalue = 1-voxelvalue itkSetMacro( OutputAbsoluteValues, bool) ///< output absolute values of the number of fibers per voxel itkGetMacro( OutputAbsoluteValues, bool) ///< output absolute values of the number of fibers per voxel itkSetMacro( UseImageGeometry, bool) ///< use input image geometry to initialize output image itkGetMacro( UseImageGeometry, bool) ///< use input image geometry to initialize output image itkSetMacro( FiberBundle, mitk::FiberBundle::Pointer) ///< input fiber bundle itkSetMacro( InputImage, typename RefImageType::Pointer) ///< use input image geometry to initialize output image itkGetMacro( MaxDensity, OutPixelType) itkGetMacro( NumCoveredVoxels, unsigned int) - itkGetMacro( AverageNumTraversedVoxels, unsigned int) - itkGetMacro( AverageSegmentLength, float) void GenerateData() override; TDI_MODE GetMode() const; void SetMode(const TDI_MODE &Mode); protected: TractDensityImageFilter(); ~TractDensityImageFilter() override; typename RefImageType::Pointer m_InputImage; ///< use input image geometry to initialize output image mitk::FiberBundle::Pointer m_FiberBundle; ///< input fiber bundle float m_UpsamplingFactor; ///< use higher resolution for ouput image bool m_InvertImage; ///< voxelvalue = 1-voxelvalue TDI_MODE m_Mode; ///< what should the output look like bool m_UseImageGeometry; ///< use input image geometry to initialize output image bool m_OutputAbsoluteValues; ///< do not normalize image values to 0-1 bool m_UseTrilinearInterpolation; bool m_DoFiberResampling; bool m_WorkOnFiberCopy; OutPixelType m_MaxDensity; unsigned int m_NumCoveredVoxels; - unsigned int m_AverageNumTraversedVoxels; - float m_AverageSegmentLength; }; } #ifndef ITK_MANUAL_INSTANTIATION #include "itkTractDensityImageFilter.cpp" #endif #endif // __itkTractDensityImageFilter_h__ diff --git a/Modules/FiberTracking/Algorithms/itkTractParcellationFilter.cpp b/Modules/FiberTracking/Algorithms/itkTractParcellationFilter.cpp deleted file mode 100644 index 0513246..0000000 --- a/Modules/FiberTracking/Algorithms/itkTractParcellationFilter.cpp +++ /dev/null @@ -1,563 +0,0 @@ -/*=================================================================== - -The Medical Imaging Interaction Toolkit (MITK) - -Copyright (c) German Cancer Research Center. - -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 "itkTractParcellationFilter.h" - -// VTK -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include - -namespace itk{ - -template< class OutImageType, class InputImageType > -TractParcellationFilter< OutImageType, InputImageType >::TractParcellationFilter() - : m_UpsamplingFactor(1) - , m_NumParcels(0) - , m_NumCentroids(0) - , m_StartClusterSize(5) - , m_InputImage(nullptr) -{ - this->SetNumberOfRequiredOutputs(2); -} - -template< class OutImageType, class InputImageType > -TractParcellationFilter< OutImageType, InputImageType >::~TractParcellationFilter() -{ -} - -template< class OutImageType, class InputImageType > -bool TractParcellationFilter< OutImageType, InputImageType >::Flip(vtkSmartPointer< vtkPolyData > polydata1, int i, vtkSmartPointer< vtkPolyData > ref_poly, int ref_i) -{ - double d_direct = 0; - double d_flipped = 0; - - vtkCell* cell1 = polydata1->GetCell(ref_i); - if (ref_poly!=nullptr) - cell1 = ref_poly->GetCell(ref_i); - auto numPoints1 = cell1->GetNumberOfPoints(); - vtkPoints* points1 = cell1->GetPoints(); - - std::vector> ref_points; - for (int j=0; jGetPoint(j); - itk::Point itk_p; - itk_p[0] = p1[0]; - itk_p[1] = p1[1]; - itk_p[2] = p1[2]; - ref_points.push_back(itk_p); - } - - vtkCell* cell2 = polydata1->GetCell(i); - vtkPoints* points2 = cell2->GetPoints(); - - for (int j=0; jGetPoint(j); - d_direct += (p1[0]-p2[0])*(p1[0]-p2[0]) + (p1[1]-p2[1])*(p1[1]-p2[1]) + (p1[2]-p2[2])*(p1[2]-p2[2]); - - double* p3 = points2->GetPoint(numPoints1-j-1); - d_flipped += (p1[0]-p3[0])*(p1[0]-p3[0]) + (p1[1]-p3[1])*(p1[1]-p3[1]) + (p1[2]-p3[2])*(p1[2]-p3[2]); - } - - if (d_direct>d_flipped) - return true; - return false; -} - -template< class OutImageType, class InputImageType > -mitk::FiberBundle::Pointer TractParcellationFilter< OutImageType, InputImageType >::GetWorkingFib() -{ - mitk::FiberBundle::Pointer fib_static_resampled = m_InputTract->GetDeepCopy(); - mitk::Tractometry::ResampleIfNecessary(fib_static_resampled, m_NumParcels); - - // clustering - std::vector< mitk::ClusteringMetric* > metrics; - metrics.push_back({new mitk::ClusteringMetricEuclideanStd()}); -// metrics.push_back({new mitk::ClusteringMetricEuclideanMean()}); -// metrics.push_back({new mitk::ClusteringMetricLength()}); - - if (m_NumCentroids>0) - { - std::vector centroids; - std::shared_ptr< mitk::TractClusteringFilter > clusterer = std::make_shared(); - int c=0; - while (c<30 && (centroids.empty() || centroids.size()>static_cast(m_NumCentroids))) - { - float cluster_size = m_StartClusterSize + m_StartClusterSize*0.2*c; - float max_d = 0; - int i=1; - std::vector< float > distances; - while (max_d < m_InputTract->GetGeometry()->GetDiagonalLength()/2) - { - distances.push_back(cluster_size*i); - max_d = cluster_size*i; - ++i; - } - - clusterer->SetDistances(distances); - clusterer->SetTractogram(fib_static_resampled); - clusterer->SetMetrics(metrics); - clusterer->SetMergeDuplicateThreshold(cluster_size); - clusterer->SetDoResampling(false); - clusterer->SetNumPoints(m_NumParcels); - // clusterer->SetMaxClusters(m_Controls->m_MaxCentroids->value()); - clusterer->SetMinClusterSize(1); - clusterer->Update(); - centroids = clusterer->GetOutCentroids(); - ++c; - } - - mitk::FiberBundle::Pointer centroid_bundle = mitk::FiberBundle::New(); - centroid_bundle = centroid_bundle->AddBundles(centroids); - return centroid_bundle; - } - - return fib_static_resampled; -} - -template< class OutImageType, class InputImageType > -std::vector< typename itk::Image::Pointer > TractParcellationFilter< OutImageType, InputImageType >::GetBinarySplit(typename OutImageType::Pointer inImage) -{ - std::vector< typename itk::Image::Pointer > binary_maps; - - for (unsigned int i=0; i::Pointer parcel_image = itk::Image::New(); - parcel_image->SetSpacing( inImage->GetSpacing() ); - parcel_image->SetOrigin( inImage->GetOrigin() ); - parcel_image->SetDirection( inImage->GetDirection() ); - parcel_image->SetRegions( inImage->GetLargestPossibleRegion() ); - parcel_image->Allocate(); - parcel_image->FillBuffer(0); - - binary_maps.push_back(parcel_image); - } - - itk::ImageRegionIterator< itk::Image > p_it(inImage, inImage->GetLargestPossibleRegion()); - while(!p_it.IsAtEnd()) - { - if (p_it.Get()>0) - { - binary_maps.at(p_it.Get()-1)->SetPixel(p_it.GetIndex(), 1); - } - ++p_it; - } - - return binary_maps; -} - -template< class OutImageType, class InputImageType > -typename OutImageType::Pointer TractParcellationFilter< OutImageType, InputImageType >::PostprocessParcellation(typename OutImageType::Pointer inImage) -{ - itk::ImageRegionConstIterator< OutImageType > in_it(inImage, inImage->GetLargestPossibleRegion()); - - typename OutImageType::Pointer outImage = OutImageType::New(); - outImage->SetSpacing( inImage->GetSpacing() ); - outImage->SetOrigin( inImage->GetOrigin() ); - outImage->SetDirection( inImage->GetDirection() ); - outImage->SetRegions( inImage->GetLargestPossibleRegion() ); - outImage->Allocate(); - outImage->FillBuffer(0); - itk::ImageRegionIterator< OutImageType > out_it(outImage, outImage->GetLargestPossibleRegion()); - - MITK_INFO << "Postprocessing parcellation"; - while( !in_it.IsAtEnd() ) - { - if (in_it.Get()>0) - { - std::vector< OutPixelType > vote; vote.resize(m_NumParcels + 1, 0); - - typename OutImageType::SizeType regionSize; - regionSize[0] = 3; - regionSize[1] = 3; - regionSize[2] = 3; - - typename OutImageType::IndexType regionIndex = in_it.GetIndex(); - regionIndex[0] -= 1; - regionIndex[1] -= 1; - regionIndex[2] -= 1; - - typename OutImageType::RegionType region; - region.SetSize(regionSize); - region.SetIndex(regionIndex); - - itk::ImageRegionConstIterator rit(inImage, region); - while( !rit.IsAtEnd() ) - { - if (rit.GetIndex()!=regionIndex) - vote[rit.Get()]++; - ++rit; - } - - OutPixelType max = 0; - unsigned int max_parcel = -1; - for (unsigned int i=1; imax) - { - max = vote[i]; - max_parcel = i; - } - } - - out_it.Set(max_parcel); - } - ++out_it; - ++in_it; - } - MITK_INFO << "DONE"; - return outImage; -} - -template< class OutImageType, class InputImageType > -void TractParcellationFilter< OutImageType, InputImageType >::CentroidBasedParcellation(typename OutImageType::Pointer outImage) -{ - typename itk::TractDensityImageFilter< OutImageType >::Pointer generator = itk::TractDensityImageFilter< OutImageType >::New(); - generator->SetFiberBundle(m_InputTract); - generator->SetMode(TDI_MODE::BINARY); - generator->SetInputImage(outImage); - generator->SetUseImageGeometry(true); - generator->Update(); - auto tdi = generator->GetOutput(); - - if (m_NumParcels < 3) - m_NumParcels = mitk::Tractometry::EstimateNumSamplingPoints(tdi, m_InputTract, 5); - - m_WorkingTract = GetWorkingFib(); - vtkSmartPointer< vtkPolyData > polydata = m_WorkingTract->GetFiberPolyData(); - - itk::ImageRegionIterator< OutImageType > it(outImage, outImage->GetLargestPossibleRegion()); - itk::ImageRegionIterator< OutImageType > it_tdi(tdi, tdi->GetLargestPossibleRegion()); - - vtkSmartPointer< vtkPolyData > reference_polydata = nullptr; - if (m_ReferenceTract.IsNotNull()) - reference_polydata = m_ReferenceTract->GetFiberPolyData(); - - unsigned long num_vox = 0; - while( !it_tdi.IsAtEnd() ) - { - if (it_tdi.Get()>0) - ++num_vox; - ++it_tdi; - } - it_tdi.GoToBegin(); - - MITK_INFO << "Parcellating tract (centroid-based)"; - boost::progress_display disp(num_vox); - while( !it.IsAtEnd() ) - { - if (it_tdi.Get()>0) - { - int final_seg_id = -1; - float min_d = 99999999; - - itk::Image< float, 4 >::IndexType idx4; - idx4[0] = it_tdi.GetIndex()[0]; - idx4[1] = it_tdi.GetIndex()[1]; - idx4[2] = it_tdi.GetIndex()[2]; - - typename OutImageType::PointType image_point; - tdi->TransformIndexToPhysicalPoint(it.GetIndex(), image_point); - - for (unsigned int i=0; iGetNumFibers(); ++i) - { - vtkCell* cell = polydata->GetCell(i); - auto numPoints = cell->GetNumberOfPoints(); - vtkPoints* points = cell->GetPoints(); - - bool flip = Flip(polydata, i, reference_polydata); - - for (int j=0; j p; - int segment_id = -1; - if (flip) - { - segment_id = numPoints - j - 1; - p = mitk::imv::GetItkPoint(points->GetPoint(segment_id)); - } - else - { - p = mitk::imv::GetItkPoint(points->GetPoint(j)); - segment_id = j; - } - - float d = std::fabs( (p[0]-image_point[0]) ) + std::fabs( (p[1]-image_point[1]) ) + std::fabs( (p[2]-image_point[2]) ); - - if (d -void TractParcellationFilter< OutImageType, InputImageType >::VotingBasedParcellation(typename OutImageType::Pointer outImage) -{ - typename itk::TractDensityImageFilter< OutImageType >::Pointer generator = itk::TractDensityImageFilter< OutImageType >::New(); - generator->SetFiberBundle(m_InputTract); - generator->SetMode(TDI_MODE::BINARY); - generator->SetInputImage(outImage); - generator->SetUseImageGeometry(true); - generator->Update(); - auto tdi = generator->GetOutput(); - - if (m_NumParcels < 3) - m_NumParcels = mitk::Tractometry::EstimateNumSamplingPoints(tdi, m_InputTract, 5); - - itk::TractsToVectorImageFilter::Pointer fOdfFilter = itk::TractsToVectorImageFilter::New(); - fOdfFilter->SetMaskImage(tdi); - fOdfFilter->SetFiberBundle(m_InputTract); - fOdfFilter->SetNormalizationMethod(itk::TractsToVectorImageFilter::NormalizationMethods::SINGLE_VEC_NORM); - fOdfFilter->SetMaxNumDirections(1); - fOdfFilter->SetOnlyUseMaskGeometry(true); - fOdfFilter->Update(); - itk::Image< float, 4 >::Pointer dir_image = fOdfFilter->GetDirectionImage(); - - m_WorkingTract = GetWorkingFib(); - vtkSmartPointer< vtkPolyData > polydata = m_WorkingTract->GetFiberPolyData(); - - float maxd = m_InputTract->GetMeanFiberLength()/(0.5*m_NumParcels); - itk::ImageRegionIterator< OutImageType > it(outImage, outImage->GetLargestPossibleRegion()); - itk::ImageRegionIterator< OutImageType > it_tdi(tdi, tdi->GetLargestPossibleRegion()); - - vtkSmartPointer< vtkPolyData > reference_polydata = nullptr; - if (m_ReferenceTract.IsNotNull()) - reference_polydata = m_ReferenceTract->GetFiberPolyData(); - - // count voxels in mask - unsigned long num_vox = 0; - while( !it_tdi.IsAtEnd() ) - { - if (it_tdi.Get()>0) - ++num_vox; - ++it_tdi; - } - it_tdi.GoToBegin(); - - // - std::vector< std::vector< itk::Point > > streamlines; - for (unsigned int i=0; iGetNumFibers(); ++i) - { - vtkCell* cell = polydata->GetCell(i); - auto numPoints = cell->GetNumberOfPoints(); - vtkPoints* points = cell->GetPoints(); - - bool flip = Flip(polydata, i, reference_polydata); - - std::vector< itk::Point > streamline; - - for (int j=0; j p; - if (flip) - p = mitk::imv::GetItkPoint(points->GetPoint(numPoints - j - 1)); - else - p = mitk::imv::GetItkPoint(points->GetPoint(j)); - streamline.push_back(p); - } - streamlines.push_back(streamline); - } - - MITK_INFO << "Parcellating tract (parcel-voting)"; - boost::progress_display disp(num_vox); - while( !it.IsAtEnd() ) - { - if (it_tdi.Get()>0) - { - int final_seg_id = -1; - int mult = 2; - - itk::Image< float, 4 >::IndexType idx4; - idx4[0] = it_tdi.GetIndex()[0]; - idx4[1] = it_tdi.GetIndex()[1]; - idx4[2] = it_tdi.GetIndex()[2]; - - vnl_vector_fixed ref_dir; - idx4[3] = 0; - ref_dir[0] = dir_image->GetPixel(idx4); - idx4[3] = 1; - ref_dir[1] = dir_image->GetPixel(idx4); - idx4[3] = 2; - ref_dir[2] = dir_image->GetPixel(idx4); - - if (ref_dir.magnitude()>0.01) - { - ref_dir.normalize(); - - while(final_seg_id<0 && mult<5) - { - std::vector seg_vote; seg_vote.resize(m_NumParcels, 0); - typename OutImageType::PointType image_point; - tdi->TransformIndexToPhysicalPoint(it.GetIndex(), image_point); - -#pragma omp parallel for - for (int sid = 0; sid(streamlines.size()); ++sid) - { - - float local_d = 99999999; - int local_closest_seg = -1; - - auto streamline = streamlines[sid]; - for (unsigned int segment_id=0; segment_id p = streamline[segment_id]; - float d = std::fabs( (p[0]-image_point[0]) ) + std::fabs( (p[1]-image_point[1]) ) + std::fabs( (p[2]-image_point[2]) ); - - itk::Point p2; - if (segment_id dir; - dir[0] = p[0]-p2[0]; - dir[1] = p[1]-p2[1]; - dir[2] = p[2]-p2[2]; - - if (dir.magnitude()<0.0000001) - continue; - dir.normalize(); - - float a = std::fabs(dot_product(dir, ref_dir)); - if (a<0.0000001) - a += 0.0000001; - d += (1.0/a - 1.0) * maxd; - - if (dmax_count) - { - final_seg_id = i; - max_count = seg_vote.at(i); - } - } - - if (final_seg_id>=0) - it.Set(final_seg_id + 1); - - ++mult; - } - } - - ++disp; - } - ++it; - ++it_tdi; - } - MITK_INFO << "DONE"; -} - -template< class OutImageType, class InputImageType > -void TractParcellationFilter< OutImageType, InputImageType >::GenerateData() -{ - // generate upsampled image - mitk::BaseGeometry::Pointer geometry = m_InputTract->GetGeometry(); - - // calculate new image parameters - itk::Vector newSpacing; - mitk::Point3D newOrigin; - itk::Matrix newDirection; - ImageRegion<3> upsampledRegion; - if (!m_InputImage.IsNull()) - { - newSpacing = m_InputImage->GetSpacing()/m_UpsamplingFactor; - upsampledRegion = m_InputImage->GetLargestPossibleRegion(); - newOrigin = m_InputImage->GetOrigin(); - typename OutImageType::RegionType::SizeType size = upsampledRegion.GetSize(); - size[0] *= m_UpsamplingFactor; - size[1] *= m_UpsamplingFactor; - size[2] *= m_UpsamplingFactor; - upsampledRegion.SetSize(size); - newDirection = m_InputImage->GetDirection(); - } - else - { - newSpacing = geometry->GetSpacing()/m_UpsamplingFactor; - newOrigin = geometry->GetOrigin(); - mitk::Geometry3D::BoundsArrayType bounds = geometry->GetBounds(); - newOrigin[0] += bounds.GetElement(0); - newOrigin[1] += bounds.GetElement(2); - newOrigin[2] += bounds.GetElement(4); - - for (int i=0; i<3; i++) - for (int j=0; j<3; j++) - newDirection[j][i] = geometry->GetMatrixColumn(i)[j]; - upsampledRegion.SetSize(0, geometry->GetExtent(0)*m_UpsamplingFactor); - upsampledRegion.SetSize(1, geometry->GetExtent(1)*m_UpsamplingFactor); - upsampledRegion.SetSize(2, geometry->GetExtent(2)*m_UpsamplingFactor); - } - - // apply new image parameters - typename OutImageType::Pointer outImage = OutImageType::New(); - outImage->SetSpacing( newSpacing ); - outImage->SetOrigin( newOrigin ); - outImage->SetDirection( newDirection ); - outImage->SetRegions( upsampledRegion ); - outImage->Allocate(); - outImage->FillBuffer(0); - this->SetNthOutput(0, outImage); - - if (m_NumCentroids > 0) - CentroidBasedParcellation(outImage); - else - VotingBasedParcellation(outImage); - - auto outImage_pp = PostprocessParcellation(outImage); - outImage_pp = PostprocessParcellation(outImage_pp); - outImage_pp = PostprocessParcellation(outImage_pp); - this->SetNthOutput(1, outImage_pp); -} -} diff --git a/Modules/FiberTracking/Algorithms/itkTractParcellationFilter.h b/Modules/FiberTracking/Algorithms/itkTractParcellationFilter.h deleted file mode 100644 index 3b0d258..0000000 --- a/Modules/FiberTracking/Algorithms/itkTractParcellationFilter.h +++ /dev/null @@ -1,97 +0,0 @@ -/*=================================================================== - -The Medical Imaging Interaction Toolkit (MITK) - -Copyright (c) German Cancer Research Center. - -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 __itkTractParcellationFilter_h__ -#define __itkTractParcellationFilter_h__ - -#include -#include -#include - -namespace itk{ - -/** -* \brief Generates image where the pixel values are set according to the position along a fiber bundle. */ - -template< class OutImageType=itk::Image, class InputImageType=OutImageType > -class TractParcellationFilter : public ImageSource< OutImageType > -{ - -public: - typedef TractParcellationFilter Self; - typedef ProcessObject Superclass; - typedef SmartPointer< Self > Pointer; - typedef SmartPointer< const Self > ConstPointer; - - typedef typename OutImageType::PixelType OutPixelType; - - itkFactorylessNewMacro(Self) - itkCloneMacro(Self) - itkTypeMacro( TractParcellationFilter, ImageSource ) - - itkSetMacro( UpsamplingFactor, float) - itkGetMacro( UpsamplingFactor, float) - - itkSetMacro( NumParcels, unsigned int ) - itkGetMacro( NumParcels, unsigned int ) - - itkSetMacro( NumCentroids, unsigned int ) - itkGetMacro( NumCentroids, unsigned int ) - - itkSetMacro( StartClusterSize, float ) - itkGetMacro( StartClusterSize, float ) - - itkSetMacro( InputTract, mitk::FiberBundle::Pointer) - itkSetMacro( ReferenceTract, mitk::FiberBundle::Pointer) - itkGetMacro( WorkingTract, mitk::FiberBundle::Pointer) - - itkSetMacro( InputImage, typename InputImageType::Pointer) - - - std::vector< typename itk::Image::Pointer > GetBinarySplit(typename OutImageType::Pointer inImage); - - void GenerateData() override; - -protected: - - TractParcellationFilter(); - virtual ~TractParcellationFilter(); - - bool Flip(vtkSmartPointer< vtkPolyData > polydata1, int i, vtkSmartPointer< vtkPolyData > ref_poly=nullptr, int ref_i=0); - - mitk::FiberBundle::Pointer GetWorkingFib(); - - typename OutImageType::Pointer PostprocessParcellation(typename OutImageType::Pointer outImage); - - void VotingBasedParcellation(typename OutImageType::Pointer outImage); - void CentroidBasedParcellation(typename OutImageType::Pointer outImage); - - mitk::FiberBundle::Pointer m_ReferenceTract; - mitk::FiberBundle::Pointer m_InputTract; ///< input fiber bundle - mitk::FiberBundle::Pointer m_WorkingTract; - float m_UpsamplingFactor; ///< use higher resolution for ouput image - unsigned int m_NumParcels; - unsigned int m_NumCentroids; - float m_StartClusterSize; - typename InputImageType::Pointer m_InputImage; -}; - -} - -#ifndef ITK_MANUAL_INSTANTIATION -#include "itkTractParcellationFilter.cpp" -#endif - -#endif // __itkTractParcellationFilter_h__ diff --git a/Modules/FiberTracking/Algorithms/mitkTractometry.cpp b/Modules/FiberTracking/Algorithms/mitkTractometry.cpp index b433937..6f27c46 100644 --- a/Modules/FiberTracking/Algorithms/mitkTractometry.cpp +++ b/Modules/FiberTracking/Algorithms/mitkTractometry.cpp @@ -1,298 +1,289 @@ /*=================================================================== The Medical Imaging Interaction Toolkit (MITK) Copyright (c) German Cancer Research Center. 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 "mitkTractometry.h" #define _USE_MATH_DEFINES #include #include #include #include #include #include namespace mitk{ bool Tractometry::Flip(vtkSmartPointer< vtkPolyData > polydata1, int i, vtkSmartPointer< vtkPolyData > ref_poly) { double d_direct = 0; double d_flipped = 0; vtkCell* cell1 = polydata1->GetCell(0); if (ref_poly!=nullptr) cell1 = ref_poly->GetCell(0); auto numPoints1 = cell1->GetNumberOfPoints(); vtkPoints* points1 = cell1->GetPoints(); std::vector> ref_points; for (int j=0; jGetPoint(j); itk::Point itk_p; itk_p[0] = p1[0]; itk_p[1] = p1[1]; itk_p[2] = p1[2]; ref_points.push_back(itk_p); } vtkCell* cell2 = polydata1->GetCell(i); vtkPoints* points2 = cell2->GetPoints(); for (int j=0; jGetPoint(j); d_direct += (p1[0]-p2[0])*(p1[0]-p2[0]) + (p1[1]-p2[1])*(p1[1]-p2[1]) + (p1[2]-p2[2])*(p1[2]-p2[2]); double* p3 = points2->GetPoint(numPoints1-j-1); d_flipped += (p1[0]-p3[0])*(p1[0]-p3[0]) + (p1[1]-p3[1])*(p1[1]-p3[1]) + (p1[2]-p3[2])*(p1[2]-p3[2]); } if (d_direct>d_flipped) return true; return false; } void Tractometry::ResampleIfNecessary(mitk::FiberBundle::Pointer fib, unsigned int num_points) { auto poly = fib->GetFiberPolyData(); bool resample = false; for (int i=0; iGetNumberOfCells(); i++) { vtkCell* cell = poly->GetCell(i); if (cell->GetNumberOfPoints()!=num_points) { resample = true; MITK_INFO << "Resampling required!"; break; } } if (resample) fib->ResampleToNumPoints(num_points); } unsigned int Tractometry::EstimateNumSamplingPoints(itk::Image::Pointer ref_image, mitk::FiberBundle::Pointer fib, unsigned int voxels) { - typename itk::TractDensityImageFilter< itk::Image >::Pointer generator = itk::TractDensityImageFilter< itk::Image >::New(); - generator->SetFiberBundle(fib); - generator->SetMode(TDI_MODE::BINARY); - generator->SetInputImage(ref_image); - generator->SetUseImageGeometry(true); - generator->Update(); - auto spacing = ref_image->GetSpacing(); float f = (spacing[0] + spacing[1] + spacing[2])/3; - f /= generator->GetAverageSegmentLength(); - MITK_INFO << generator->GetAverageSegmentLength(); - MITK_INFO << generator->GetAverageNumTraversedVoxels(); - MITK_INFO << f; - unsigned int n = std::ceil(generator->GetAverageNumTraversedVoxels()/(voxels*f)); - if (n<3) - n = 3; + float num_voxels_passed = 0; + for (unsigned int i=0; iGetNumFibers(); ++i) + num_voxels_passed += fib->GetFiberLength(i)/f; + num_voxels_passed /= fib->GetNumFibers(); + unsigned int parcels = std::ceil(num_voxels_passed/voxels); - MITK_INFO << "Estimated number of sampling points " << n; + MITK_INFO << "Estimated number of sampling points " << parcels; - return n; + return parcels; } std::vector> Tractometry::NearestCentroidPointTractometry(itk::Image::Pointer itkImage, mitk::FiberBundle::Pointer working_fib, unsigned int num_points, unsigned int max_centroids, float cluster_size, mitk::FiberBundle::Pointer ref_fib) { vtkSmartPointer< vtkPolyData > polydata = working_fib->GetFiberPolyData(); vtkSmartPointer< vtkPolyData > ref_polydata = nullptr; if (ref_fib!=nullptr) { ResampleIfNecessary(ref_fib, num_points); ref_polydata = ref_fib->GetFiberPolyData(); } auto interpolator = itk::LinearInterpolateImageFunction< itk::Image, float >::New(); interpolator->SetInputImage(itkImage); mitk::LookupTable::Pointer lookupTable = mitk::LookupTable::New(); lookupTable->SetType(mitk::LookupTable::MULTILABEL); std::vector> output_temp; for(unsigned int i=0; i metrics; metrics.push_back({new mitk::ClusteringMetricEuclideanStd()}); mitk::FiberBundle::Pointer resampled = working_fib->GetDeepCopy(); ResampleIfNecessary(resampled, num_points); std::vector centroids; std::shared_ptr< mitk::TractClusteringFilter > clusterer = std::make_shared(); int c=0; while (c<30 && (centroids.empty() || centroids.size()>max_centroids)) { float cs = cluster_size + cluster_size*c*0.2; float max_d = 0; int i=1; std::vector< float > distances; while (max_d < resampled->GetGeometry()->GetDiagonalLength()/2) { distances.push_back(cs*i); max_d = cs*i; ++i; } clusterer->SetDistances(distances); clusterer->SetTractogram(resampled); clusterer->SetMetrics(metrics); clusterer->SetMergeDuplicateThreshold(cs); clusterer->SetDoResampling(false); clusterer->SetNumPoints(num_points); if (c==29) clusterer->SetMaxClusters(max_centroids); clusterer->SetMinClusterSize(1); clusterer->Update(); centroids = clusterer->GetOutCentroids(); ++c; } for (unsigned int i=0; iGetNumFibers(); ++i) { vtkCell* cell = polydata->GetCell(i); auto numPoints = cell->GetNumberOfPoints(); vtkPoints* points = cell->GetPoints(); vnl_vector values; values.set_size(numPoints); for (int j=0; jGetPoint(j); int min_bin = 0; float d=999999; for (auto centroid : centroids) { auto centroid_polydata = centroid->GetFiberPolyData(); vtkCell* centroid_cell = centroid_polydata->GetCell(0); auto centroid_numPoints = centroid_cell->GetNumberOfPoints(); vtkPoints* centroid_points = centroid_cell->GetPoints(); bool centroid_flip = Flip(centroid_polydata, 0, ref_polydata); for (int bin=0; binGetPoint(bin); float temp_d = std::sqrt((p[0]-centroid_p[0])*(p[0]-centroid_p[0]) + (p[1]-centroid_p[1])*(p[1]-centroid_p[1]) + (p[2]-centroid_p[2])*(p[2]-centroid_p[2])); if (temp_dGetColor(min_bin+1, rgb); working_fib->ColorSinglePoint(i, j, rgb); double pixelValue = mitk::imv::GetImageValue(mitk::imv::GetItkPoint(p), true, interpolator); output_temp.at(min_bin).push_back(pixelValue); } } std::vector> output; for (auto row_v : output_temp) { vnl_vector row; row.set_size(row_v.size()); int i = 0; for (auto v : row_v) { row.put(i, v); ++i; } output.push_back(row); } return output; } vnl_matrix Tractometry::StaticResamplingTractometry(itk::Image::Pointer itkImage, mitk::FiberBundle::Pointer working_fib, unsigned int num_points, mitk::FiberBundle::Pointer ref_fib) { ResampleIfNecessary(working_fib, num_points); vtkSmartPointer< vtkPolyData > polydata = working_fib->GetFiberPolyData(); vtkSmartPointer< vtkPolyData > ref_polydata = nullptr; if (ref_fib!=nullptr) { ResampleIfNecessary(ref_fib, num_points); ref_polydata = ref_fib->GetFiberPolyData(); } auto interpolator = itk::LinearInterpolateImageFunction< itk::Image, float >::New(); interpolator->SetInputImage(itkImage); mitk::LookupTable::Pointer lookupTable = mitk::LookupTable::New(); lookupTable->SetType(mitk::LookupTable::MULTILABEL); vnl_matrix output; output.set_size(num_points, working_fib->GetNumFibers()); output.fill(0.0); for (unsigned int i=0; iGetNumFibers(); ++i) { vtkCell* cell = polydata->GetCell(i); auto numPoints = cell->GetNumberOfPoints(); vtkPoints* points = cell->GetPoints(); bool flip = Flip(polydata, i, ref_polydata); for (int j=0; jGetColor(j+1, rgb); double* p; if (flip) { auto p_idx = numPoints - j - 1; p = points->GetPoint(p_idx); working_fib->ColorSinglePoint(i, p_idx, rgb); } else { p = points->GetPoint(j); working_fib->ColorSinglePoint(i, j, rgb); } double pixelValue = mitk::imv::GetImageValue(mitk::imv::GetItkPoint(p), true, interpolator); output.put(j, i, pixelValue); } } return output; } } diff --git a/Modules/FiberTracking/files.cmake b/Modules/FiberTracking/files.cmake index 81a3950..cc72a1f 100644 --- a/Modules/FiberTracking/files.cmake +++ b/Modules/FiberTracking/files.cmake @@ -1,75 +1,74 @@ set(CPP_FILES mitkStreamlineTractographyParameters.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 Algorithms/mitkTractClusteringFilter.cpp Algorithms/itkStreamlineTrackingFilter.cpp Algorithms/TrackingHandlers/mitkTrackingDataHandler.cpp Algorithms/TrackingHandlers/mitkTrackingHandlerTensor.cpp Algorithms/TrackingHandlers/mitkTrackingHandlerPeaks.cpp Algorithms/TrackingHandlers/mitkTrackingHandlerOdf.cpp Algorithms/mitkTractometry.cpp ) set(H_FILES mitkStreamlineTractographyParameters.h # Algorithms Algorithms/itkTractDensityImageFilter.h Algorithms/itkTractsToFiberEndingsImageFilter.h Algorithms/itkTractsToRgbaImageFilter.h Algorithms/itkTractsToVectorImageFilter.h Algorithms/itkEvaluateDirectionImagesFilter.h Algorithms/itkEvaluateTractogramDirectionsFilter.h Algorithms/itkFiberCurvatureFilter.h Algorithms/mitkTractClusteringFilter.h Algorithms/itkTractDistanceFilter.h Algorithms/itkFiberExtractionFilter.h Algorithms/itkTdiToVolumeFractionFilter.h Algorithms/itkDistanceFromSegmentationImageFilter.h - Algorithms/itkTractParcellationFilter.h Algorithms/mitkTractometry.h # Tractography Algorithms/TrackingHandlers/mitkTrackingDataHandler.h Algorithms/TrackingHandlers/mitkTrackingHandlerRandomForest.h Algorithms/TrackingHandlers/mitkTrackingHandlerTensor.h Algorithms/TrackingHandlers/mitkTrackingHandlerPeaks.h Algorithms/TrackingHandlers/mitkTrackingHandlerOdf.h Algorithms/itkGibbsTrackingFilter.h Algorithms/itkStochasticTractographyFilter.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 Algorithms/itkStreamlineTrackingFilter.h # Clustering Algorithms/ClusteringMetrics/mitkClusteringMetric.h Algorithms/ClusteringMetrics/mitkClusteringMetricEuclideanMean.h Algorithms/ClusteringMetrics/mitkClusteringMetricEuclideanMax.h Algorithms/ClusteringMetrics/mitkClusteringMetricEuclideanStd.h Algorithms/ClusteringMetrics/mitkClusteringMetricAnatomic.h Algorithms/ClusteringMetrics/mitkClusteringMetricScalarMap.h Algorithms/ClusteringMetrics/mitkClusteringMetricInnerAngles.h Algorithms/ClusteringMetrics/mitkClusteringMetricLength.h ) set(RESOURCE_FILES # Binary directory resources FiberTrackingLUTBaryCoords.bin FiberTrackingLUTIndices.bin ) diff --git a/Plugins/org.mitk.gui.qt.diffusionimaging.fiberprocessing/src/internal/QmitkTractometryView.cpp b/Plugins/org.mitk.gui.qt.diffusionimaging.fiberprocessing/src/internal/QmitkTractometryView.cpp index a6a01a0..70a1423 100644 --- a/Plugins/org.mitk.gui.qt.diffusionimaging.fiberprocessing/src/internal/QmitkTractometryView.cpp +++ b/Plugins/org.mitk.gui.qt.diffusionimaging.fiberprocessing/src/internal/QmitkTractometryView.cpp @@ -1,450 +1,354 @@ /*=================================================================== The Medical Imaging Interaction Toolkit (MITK) Copyright (c) German Cancer Research Center. 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 "QmitkTractometryView.h" #include #include #include #include #include #include #include #include #include #include #include #include #include #include #include #include #include #include #include #include #include #include -#include #include #include #include #include #include const std::string QmitkTractometryView::VIEW_ID = "org.mitk.views.tractometry"; using namespace mitk; QmitkTractometryView::QmitkTractometryView() : QmitkAbstractView() , m_Controls( nullptr ) , m_Visible(false) { } // Destructor QmitkTractometryView::~QmitkTractometryView() { } void QmitkTractometryView::CreateQtPartControl( QWidget *parent ) { // build up qt view, unless already done if ( !m_Controls ) { // create GUI widgets from the Qt Designer's .ui file m_Controls = new Ui::QmitkTractometryViewControls; m_Controls->setupUi( parent ); connect( m_Controls->m_MethodBox, SIGNAL(currentIndexChanged(int)), this, SLOT(UpdateGui()) ); connect( m_Controls->m_StartButton, SIGNAL(clicked()), this, SLOT(StartTractometry()) ); mitk::TNodePredicateDataType::Pointer imageP = mitk::TNodePredicateDataType::New(); mitk::NodePredicateDimension::Pointer dimP = mitk::NodePredicateDimension::New(3); m_Controls->m_ImageBox->SetDataStorage(this->GetDataStorage()); m_Controls->m_ImageBox->SetPredicate(mitk::NodePredicateAnd::New(imageP, dimP)); m_Controls->m_ChartWidget->SetXAxisLabel("Tract position"); m_Controls->m_ChartWidget->SetYAxisLabel("Image Value"); m_Controls->m_ChartWidget->SetTheme(QmitkChartWidget::ColorTheme::darkstyle); } } void QmitkTractometryView::SetFocus() { } void QmitkTractometryView::UpdateGui() { berry::IWorkbenchPart::Pointer nullPart; OnSelectionChanged(nullPart, QList(m_CurrentSelection)); } void QmitkTractometryView::StaticResamplingTractometry(mitk::Image::Pointer image, mitk::DataNode::Pointer node, std::vector > &data, std::string& clipboard_string) { itk::Image::Pointer itkImage = itk::Image::New(); CastToItkImage(image, itkImage); mitk::FiberBundle::Pointer fib = dynamic_cast(node->GetData()); unsigned int num_points = m_NumSamplingPoints; mitk::FiberBundle::Pointer working_fib = fib->GetDeepCopy(); vnl_matrix output = mitk::Tractometry::StaticResamplingTractometry(itkImage, working_fib, num_points, m_ReferenceFib); std::vector< double > std_values1; std::vector< double > std_values2; std::vector< double > mean_values; for (unsigned int i=0; i(mean); clipboard_string += " "; clipboard_string += boost::lexical_cast(stdev); clipboard_string += "\n"; mean_values.push_back(mean); std_values1.push_back(mean + stdev); std_values2.push_back(mean - stdev); } clipboard_string += "\n"; data.push_back(mean_values); data.push_back(std_values1); data.push_back(std_values2); if (m_Controls->m_ShowBinned->isChecked()) { mitk::DataNode::Pointer new_node = mitk::DataNode::New(); new_node->SetData(working_fib); new_node->SetName("binned_static"); new_node->SetVisibility(true); node->SetVisibility(false); GetDataStorage()->Add(new_node, node); } } void QmitkTractometryView::NearestCentroidPointTractometry(mitk::Image::Pointer image, mitk::DataNode::Pointer node, std::vector< std::vector< double > >& data, std::string& clipboard_string) { mitk::FiberBundle::Pointer fib = dynamic_cast(node->GetData()); unsigned int num_points = m_NumSamplingPoints; mitk::FiberBundle::Pointer working_fib = fib->GetDeepCopy(); working_fib->ResampleSpline(1.0); itk::Image::Pointer itkImage = itk::Image::New(); CastToItkImage(image, itkImage); auto output = mitk::Tractometry::NearestCentroidPointTractometry(itkImage, working_fib, num_points, m_Controls->m_MaxCentroids->value(), m_Controls->m_ClusterSize->value(), m_ReferenceFib); std::vector< double > std_values1; std::vector< double > std_values2; std::vector< double > mean_values; for (auto row : output) { float mean = row.mean(); double stdev = 0; for (unsigned int j=0; j(mean); clipboard_string += " "; clipboard_string += boost::lexical_cast(stdev); clipboard_string += "\n"; mean_values.push_back(mean); std_values1.push_back(mean + stdev); std_values2.push_back(mean - stdev); } clipboard_string += "\n"; data.push_back(mean_values); data.push_back(std_values1); data.push_back(std_values2); if (m_Controls->m_ShowBinned->isChecked()) { mitk::DataNode::Pointer new_node = mitk::DataNode::New(); new_node->SetData(working_fib); new_node->SetName("binned_centroid"); new_node->SetVisibility(true); node->SetVisibility(false); GetDataStorage()->Add(new_node, node); } } void QmitkTractometryView::Activated() { } void QmitkTractometryView::Deactivated() { } void QmitkTractometryView::Visible() { m_Visible = true; QList selection = GetDataManagerSelection(); berry::IWorkbenchPart::Pointer nullPart; OnSelectionChanged(nullPart, selection); } void QmitkTractometryView::Hidden() { m_Visible = false; } std::string QmitkTractometryView::RGBToHexString(double *rgb) { std::ostringstream os; for (int i = 0; i < 3; ++i) { os << std::setw(2) << std::setfill('0') << std::hex << static_cast(rgb[i] * 255); } return os.str(); } -void QmitkTractometryView::AlongTractRadiomicsPreprocessing(mitk::Image::Pointer image, mitk::DataNode::Pointer node, std::vector< std::vector< double > >& data, std::string& clipboard_string) -{ - mitk::FiberBundle::Pointer fib = dynamic_cast(node->GetData()); - - // calculate mask - typedef itk::Image ParcellationImageType; - ParcellationImageType::Pointer itkImage = ParcellationImageType::New(); - CastToItkImage(image, itkImage); - - itk::TractParcellationFilter< >::Pointer parcellator = itk::TractParcellationFilter< >::New(); - parcellator->SetInputImage(itkImage); - parcellator->SetNumParcels(m_NumSamplingPoints); - parcellator->SetInputTract(fib); - parcellator->SetNumCentroids(m_Controls->m_MaxCentroids->value()); - parcellator->SetStartClusterSize(m_Controls->m_ClusterSize->value()); - parcellator->Update(); - ParcellationImageType::Pointer out_image = parcellator->GetOutput(0); - ParcellationImageType::Pointer out_image_pp = parcellator->GetOutput(1); - auto binary_masks = parcellator->GetBinarySplit(out_image_pp); - - mitk::Image::Pointer seg_img = mitk::Image::New(); - seg_img->InitializeByItk(out_image.GetPointer()); - seg_img->SetVolume(out_image->GetBufferPointer()); - - mitk::Image::Pointer seg_img_pp = mitk::Image::New(); - seg_img_pp->InitializeByItk(out_image_pp.GetPointer()); - seg_img_pp->SetVolume(out_image_pp->GetBufferPointer()); - - std::vector< double > std_values1; - std::vector< double > std_values2; - std::vector< double > mean_values; - for (auto mask : binary_masks) - { - itk::Image::Pointer data_image = itk::Image::New(); - CastToItkImage(image, data_image); - itk::MaskedStatisticsImageFilter>::Pointer statisticsImageFilter = itk::MaskedStatisticsImageFilter>::New(); - statisticsImageFilter->SetInput(data_image); - statisticsImageFilter->SetMask(mask); - statisticsImageFilter->Update(); - double mean = statisticsImageFilter->GetMean(); - double stdev = std::sqrt(statisticsImageFilter->GetVariance()); - - std_values1.push_back(mean + stdev); - std_values2.push_back(mean - stdev); - mean_values.push_back(mean); - - clipboard_string += boost::lexical_cast(mean); - clipboard_string += " "; - clipboard_string += boost::lexical_cast(stdev); - clipboard_string += "\n"; - } - clipboard_string += "\n"; - - data.push_back(mean_values); - data.push_back(std_values1); - data.push_back(std_values2); - - mitk::LookupTable::Pointer lut = mitk::LookupTable::New(); - lut->SetType( mitk::LookupTable::MULTILABEL ); - mitk::LookupTableProperty::Pointer lut_prop = mitk::LookupTableProperty::New(); - lut_prop->SetLookupTable( lut ); - - mitk::LevelWindow lw; - lw.SetRangeMinMax(0, parcellator->GetNumParcels()); - -// mitk::DataNode::Pointer new_node = mitk::DataNode::New(); -// new_node->SetData(seg_img); -// new_node->SetName("tract parcellation"); -// new_node->SetVisibility(true); -// new_node->SetProperty("LookupTable", lut_prop ); -// new_node->SetProperty( "levelwindow", mitk::LevelWindowProperty::New( lw ) ); -// node->SetVisibility(false); -// GetDataStorage()->Add(new_node, node); - - mitk::DataNode::Pointer new_node2 = mitk::DataNode::New(); - new_node2->SetData(seg_img_pp); - new_node2->SetName("tract parcellation pp"); - new_node2->SetVisibility(false); - new_node2->SetProperty("LookupTable", lut_prop ); - new_node2->SetProperty( "levelwindow", mitk::LevelWindowProperty::New( lw ) ); - GetDataStorage()->Add(new_node2, node); - - mitk::DataNode::Pointer new_node3 = mitk::DataNode::New(); - auto working_fib = fib->GetDeepCopy(); - working_fib->ColorFibersByScalarMap(seg_img, false, false, mitk::LookupTable::LookupTableType::MULTILABEL, 0.9); - new_node3->SetData(working_fib); - new_node3->SetName("bin-colored fib"); - GetDataStorage()->Add(new_node3, node); -} - void QmitkTractometryView::StartTractometry() { m_ReferenceFib = dynamic_cast(m_CurrentSelection.at(0)->GetData())->GetDeepCopy(); mitk::Image::Pointer image = dynamic_cast(m_Controls->m_ImageBox->GetSelectedNode()->GetData()); MITK_INFO << "Resanmpling reference fibers"; if (m_Controls->m_SamplingPointsBox->value()<3) { typedef itk::Image ParcellationImageType; ParcellationImageType::Pointer itkImage = ParcellationImageType::New(); CastToItkImage(image, itkImage); m_NumSamplingPoints = mitk::Tractometry::EstimateNumSamplingPoints(itkImage, m_ReferenceFib, 5); } else m_NumSamplingPoints = m_Controls->m_SamplingPointsBox->value(); m_ReferenceFib->ResampleToNumPoints(m_NumSamplingPoints); double color[3] = {0,0,0}; mitk::LookupTable::Pointer lookupTable = mitk::LookupTable::New(); lookupTable->SetType(mitk::LookupTable::MULTILABEL); this->m_Controls->m_ChartWidget->Clear(); std::string clipboardString = ""; int c = 1; for (auto node : m_CurrentSelection) { clipboardString += node->GetName() + "\n"; clipboardString += "mean stdev\n"; std::vector< std::vector< double > > data; switch (m_Controls->m_MethodBox->currentIndex()) { case 0: { StaticResamplingTractometry( image, node, data, clipboardString ); break; } case 1: { NearestCentroidPointTractometry( image, node, data, clipboardString ); break; } - case 2: - { - AlongTractRadiomicsPreprocessing(image, node, data, clipboardString); - break; - } default: { StaticResamplingTractometry( image, node, data, clipboardString ); } } m_Controls->m_ChartWidget->AddData1D(data.at(0), node->GetName() + " Mean", QmitkChartWidget::ChartType::line); m_Controls->m_ChartWidget->SetLineStyle(node->GetName() + " Mean", QmitkChartWidget::LineStyle::solid); if (m_Controls->m_StDevBox->isChecked()) { m_Controls->m_ChartWidget->AddData1D(data.at(1), node->GetName() + " +STDEV", QmitkChartWidget::ChartType::line); m_Controls->m_ChartWidget->AddData1D(data.at(2), node->GetName() + " -STDEV", QmitkChartWidget::ChartType::line); m_Controls->m_ChartWidget->SetLineStyle(node->GetName() + " +STDEV", QmitkChartWidget::LineStyle::dashed); m_Controls->m_ChartWidget->SetLineStyle(node->GetName() + " -STDEV", QmitkChartWidget::LineStyle::dashed); } lookupTable->GetTableValue(c, color); this->m_Controls->m_ChartWidget->SetColor(node->GetName() + " Mean", RGBToHexString(color)); if (m_Controls->m_StDevBox->isChecked()) { color[0] *= 0.5; color[1] *= 0.5; color[2] *= 0.5; this->m_Controls->m_ChartWidget->SetColor(node->GetName() + " +STDEV", RGBToHexString(color)); this->m_Controls->m_ChartWidget->SetColor(node->GetName() + " -STDEV", RGBToHexString(color)); } this->m_Controls->m_ChartWidget->Show(true); this->m_Controls->m_ChartWidget->SetShowDataPoints(false); ++c; } QApplication::clipboard()->setText(clipboardString.c_str(), QClipboard::Clipboard); } void QmitkTractometryView::OnSelectionChanged(berry::IWorkbenchPart::Pointer /*part*/, const QList& nodes) { m_Controls->m_StartButton->setEnabled(false); if (!m_Visible) return; if (m_Controls->m_MethodBox->currentIndex()==0) m_Controls->m_ClusterFrame->setVisible(false); else m_Controls->m_ClusterFrame->setVisible(true); m_CurrentSelection.clear(); if(m_Controls->m_ImageBox->GetSelectedNode().IsNull()) return; for (auto node: nodes) if ( dynamic_cast(node->GetData()) ) m_CurrentSelection.push_back(node); if (!m_CurrentSelection.empty()) m_Controls->m_StartButton->setEnabled(true); } diff --git a/Plugins/org.mitk.gui.qt.diffusionimaging.fiberprocessing/src/internal/QmitkTractometryView.h b/Plugins/org.mitk.gui.qt.diffusionimaging.fiberprocessing/src/internal/QmitkTractometryView.h index ed7b3b5..b408767 100644 --- a/Plugins/org.mitk.gui.qt.diffusionimaging.fiberprocessing/src/internal/QmitkTractometryView.h +++ b/Plugins/org.mitk.gui.qt.diffusionimaging.fiberprocessing/src/internal/QmitkTractometryView.h @@ -1,84 +1,82 @@ /*=================================================================== The Medical Imaging Interaction Toolkit (MITK) Copyright (c) German Cancer Research Center. 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 QmitkTractometryView_h #define QmitkTractometryView_h #include #include "ui_QmitkTractometryViewControls.h" #include #include #include #include #include #include /*! \brief Weight fibers by linearly fitting them to the image data. */ class QmitkTractometryView : public QmitkAbstractView, public mitk::ILifecycleAwarePart { // 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; QmitkTractometryView(); virtual ~QmitkTractometryView(); virtual void CreateQtPartControl(QWidget *parent) override; void StaticResamplingTractometry(mitk::Image::Pointer image, mitk::DataNode::Pointer node, std::vector< std::vector< double > >& data, std::string& clipboard_string); void NearestCentroidPointTractometry(mitk::Image::Pointer image, mitk::DataNode::Pointer node, std::vector< std::vector< double > >& data, std::string& clipboard_string); - void AlongTractRadiomicsPreprocessing(mitk::Image::Pointer image, mitk::DataNode::Pointer node, std::vector > &data, std::string &clipboard_string); - virtual void SetFocus() override; virtual void Activated() override; virtual void Deactivated() override; virtual void Visible() override; virtual void Hidden() override; protected slots: void UpdateGui(); void StartTractometry(); protected: /// \brief called by QmitkAbstractView when DataManager's selection has changed virtual void OnSelectionChanged(berry::IWorkbenchPart::Pointer part, const QList& nodes) override; Ui::QmitkTractometryViewControls* m_Controls; std::string RGBToHexString(double *rgb); mitk::FiberBundle::Pointer m_ReferenceFib; QList m_CurrentSelection; unsigned int m_NumSamplingPoints; bool m_Visible; }; #endif // _QMITKFIBERTRACKINGVIEW_H_INCLUDED diff --git a/Plugins/org.mitk.gui.qt.diffusionimaging.fiberprocessing/src/internal/QmitkTractometryViewControls.ui b/Plugins/org.mitk.gui.qt.diffusionimaging.fiberprocessing/src/internal/QmitkTractometryViewControls.ui index a99b28f..effdb61 100644 --- a/Plugins/org.mitk.gui.qt.diffusionimaging.fiberprocessing/src/internal/QmitkTractometryViewControls.ui +++ b/Plugins/org.mitk.gui.qt.diffusionimaging.fiberprocessing/src/internal/QmitkTractometryViewControls.ui @@ -1,225 +1,223 @@ QmitkTractometryViewControls 0 0 484 951 Form QCommandLinkButton:disabled { border: none; } QGroupBox { background-color: transparent; } Qt::Vertical 20 40 Show STDEV true Show binned tract true 0 600 Centroid Options - 0 + 1 9999 - 0 + 1 Cluster Size: Max. Number of Centroids: + + 99999.000000000000000 + - 15.000000000000000 + 99999.000000000000000 QFrame::NoFrame QFrame::Raised 0 0 0 0 6 Binning Method: Input Image: - 2 + 0 Static Resampling Centroid Based - - - Segment voting - - Sampling Points: 0 99999 0 Start Tractometry QmitkDataStorageComboBox QComboBox
QmitkDataStorageComboBox.h
QmitkChartWidget QWidget
QmitkChartWidget.h