diff --git a/Modules/DiffusionImaging/DiffusionCore/IODataStructures/DiffusionWeightedImages/mitkDiffusionImage.h b/Modules/DiffusionImaging/DiffusionCore/IODataStructures/DiffusionWeightedImages/mitkDiffusionImage.h index ab89e3ac19..0c3092e558 100644 --- a/Modules/DiffusionImaging/DiffusionCore/IODataStructures/DiffusionWeightedImages/mitkDiffusionImage.h +++ b/Modules/DiffusionImaging/DiffusionCore/IODataStructures/DiffusionWeightedImages/mitkDiffusionImage.h @@ -1,224 +1,218 @@ /*=================================================================== 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 __mitkDiffusionImage__h #define __mitkDiffusionImage__h #include "mitkImage.h" #include "itkVectorImage.h" #include "itkVectorImageToImageAdaptor.h" #include namespace mitk { /** * \brief this class encapsulates diffusion volumes (vectorimages not * yet supported by mitkImage) */ template< class TPixelType > class DiffusionImage : public Image { public: typedef typename itk::VectorImage ImageType; typedef vnl_vector_fixed< double, 3 > GradientDirectionType; typedef itk::VectorContainer< unsigned int,GradientDirectionType > GradientDirectionContainerType; typedef GradientDirectionContainerType::Pointer GradientDirectionContainerTypePointer; typedef vnl_matrix_fixed< double, 3, 3 > MeasurementFrameType; typedef std::vector< unsigned int > IndicesVector; typedef itk::VectorImage< TPixelType, 3 > ItkDwiType; typedef itk::ImageDuplicator< ItkDwiType > DwiDuplicatorType; /** * \brief The BValueMap contains seperated IndicesVectors for each b value (index for GradientDirectionContainer) * key := b value * value := indicesVector */ typedef std::map< unsigned int , IndicesVector > BValueMap; mitkClassMacro( DiffusionImage, Image ) itkFactorylessNewMacro(Self) itkCloneMacro(Self) /** * \brief Return the itkVectorImage as pointer */ typename ImageType::Pointer GetVectorImage(); /** * \brief Return the itkVectorImage as const pointer */ const typename ImageType::Pointer GetVectorImage() const; /** * \brief Set the itkVectorImage * * \param SmartPointer ofitk::VectorImage */ void SetVectorImage(typename ImageType::Pointer image ); /** * \brief Initialize the mitkImage (base-class) using the first b-zero value of the itkVectorImage * * \warning First set itkVectorImage, GradientDirectionContainer and the b value of the image. * Has to be called. */ void InitializeFromVectorImage(); /** * \brief Return the directions as GradientDirectionContainer const pointer * * \warning MeasurementFrame already applied */ GradientDirectionContainerType::Pointer GetDirections() const; /** * \brief Return the directions as GradientDirectionContainer const pointer * * \warning no MeasurmentFrame applied */ GradientDirectionContainerType::Pointer GetDirectionsWithoutMeasurementFrame() const; /** * \brief Set the original and current GradientDirectionContainer * \param GradientdirectionContainer as pointer * \warning Replace the original (m_OriginalDirections) and current (m_Directions) GradientDirectionContainer. * For the current directions the MeasurmentFrame is applied. Remove existing GradientDirectionContainer observer * of m_Directions and add a new observer to it (observer update the m_B_ValueMap if modifications occours on the container) */ void SetDirections( GradientDirectionContainerType::Pointer directions ); /** * \brief Set the original and current GradientDirectionContainer * \param std::vector > * \warning Replace the original (m_OriginalDirections) and current (m_Directions) GradientDirectionContainer. * For the current directions the MeasurmentFrame is applied. * Remove existing GradientDirectionContainer observer of m_Directions and add a new observer to it * (observer update the m_B_ValueMap if modifications occours on the container). * Calls the overloaded GradientDirectionContainer pointer version of the method. */ void SetDirections( const std::vector > & directions ); /** * \brief Return a copy of the MeasurmentFrame (m_MeasurementFrame) */ MeasurementFrameType GetMeasurementFrame() const; /** * \brief Set the MeasurementFrame * * \param Const reference of MeasurementFrameType */ void SetMeasurementFrame( const MeasurementFrameType & mFrame ); /** * \brief Return true if gradients of with diffrent b-values exist */ bool GetNumberOfBValues() const; /** * \brief Return a BValueMap (key: b value, Value: IndicesVector of the GradientDirectionContainer) */ const BValueMap & GetBValueMap() const; /** * \brief Return the reference b value * GradientsLength * referenceBValue encoding the true b value of the corresponding GradientDirection */ float GetReferenceBValue() const; /** * \brief Set the reference b value s */ void SetReferenceBValue( float val ); GradientDirectionContainerType::Pointer CalcAveragedDirectionSet(double precision, GradientDirectionContainerType::Pointer directions); void AverageRedundantGradients(double precision); void SetDisplayIndexForRendering(int displayIndex); - - /** - * \brief Return deep copy of this image - */ - mitk::DiffusionImage< TPixelType >* GetDeepCopy(); - protected: mitkCloneMacro(Self); DiffusionImage(); DiffusionImage(const DiffusionImage &); virtual ~DiffusionImage(); // Helper Methods bool AreAlike(GradientDirectionType g1, GradientDirectionType g2, double precision); float GetB_Value(unsigned int i); /** * \brief Apply the previouse set MeasurmentFrame to all gradients in the GradientsDirectionContainer (m_Directions) * * \warning first set the MeasurmentFrame */ void ApplyMeasurementFrame(); /** * \brief Update the BValueMap (m_B_ValueMap) using the current gradient directions (m_Directions) * * \warning Have to be called after each manipulation on the GradientDirectionContainer * !especially after manipulation of the m_Directions (GetDirections()) container via pointer access! */ void UpdateBValueMap(); typename ImageType::Pointer m_VectorImage; float m_B_Value; MeasurementFrameType m_MeasurementFrame; GradientDirectionContainerType::Pointer m_OriginalDirections; GradientDirectionContainerType::Pointer m_Directions; int m_DisplayIndex; BValueMap m_B_ValueMap; }; /** * @brief Equal A function comparing two images for beeing equal in meta- and imagedata * * @ingroup MITKTestingAPI * * Following aspects are tested for equality: * - mitk image equal test * - GradientDirectionContainer * - MeasurementFrame * - reference BValue * - BValueMap * - itkVectorImage * * @param rightHandSide An image to be compared * @param leftHandSide An image to be compared * @param eps Tolarence for comparison. You can use mitk::eps in most cases. * @param verbose Flag indicating if the user wants detailed console output or not. * @return true, if all subsequent comparisons are true, false otherwise */ template inline bool Equal(const mitk::DiffusionImage &leftHandSide, const mitk::DiffusionImage &rightHandSide, ScalarType eps, bool verbose ); } // namespace mitk #include "mitkDiffusionImage.txx" #endif /* __mitkDiffusionImage__h */ diff --git a/Modules/DiffusionImaging/DiffusionCore/IODataStructures/DiffusionWeightedImages/mitkDiffusionImage.txx b/Modules/DiffusionImaging/DiffusionCore/IODataStructures/DiffusionWeightedImages/mitkDiffusionImage.txx index f0e190a48f..9bfe0a5d0e 100644 --- a/Modules/DiffusionImaging/DiffusionCore/IODataStructures/DiffusionWeightedImages/mitkDiffusionImage.txx +++ b/Modules/DiffusionImaging/DiffusionCore/IODataStructures/DiffusionWeightedImages/mitkDiffusionImage.txx @@ -1,655 +1,631 @@ /*=================================================================== 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 "itkImageRegionIterator.h" #include "itkImageRegionConstIterator.h" #include "mitkImageCast.h" #include "mitkDiffusionImage.h" #include "itkImageDuplicator.h" #include "itkTestingComparisonImageFilter.h" template mitk::DiffusionImage::DiffusionImage() : m_VectorImage(0), m_B_Value(-1.0), m_OriginalDirections(0), m_Directions(0) { MeasurementFrameType mf; for(int i=0; i<3; i++) for(int j=0; j<3; j++) mf[i][j] = 0; for(int i=0; i<3; i++) mf[i][i] = 1; m_MeasurementFrame = mf; } template mitk::DiffusionImage::DiffusionImage(const DiffusionImage & orig) : mitk::Image(orig), m_VectorImage( 0 ), m_B_Value(0), m_OriginalDirections(0), m_Directions(0) { // Deep copy VectorImage typename itk::ImageDuplicator::Pointer duplicator = itk::ImageDuplicator::New(); duplicator->SetInputImage( orig.m_VectorImage ); duplicator->Update(); GradientDirectionContainerType::ConstIterator origIt; GradientDirectionContainerType::ConstIterator origItEnd; // Deep Copy OrignalDirectioncontainer GradientDirectionContainerType::Pointer OriginalDirectionscontainer = GradientDirectionContainerType::New(); origIt = orig.GetDirectionsWithoutMeasurementFrame()->Begin(); origItEnd = orig.GetDirectionsWithoutMeasurementFrame()->End(); for(;origIt != origItEnd; ++origIt) OriginalDirectionscontainer->push_back(origIt.Value()); // Deep Copy Directioncontainer GradientDirectionContainerType::Pointer DirectionsContainer = GradientDirectionContainerType::New(); origIt = orig.GetDirections()->Begin(); origItEnd = orig.GetDirections()->End(); for(;origIt != origItEnd; ++origIt) DirectionsContainer->push_back(origIt.Value()); // Set member of the new clone m_VectorImage = duplicator->GetOutput(); m_B_Value = orig.GetReferenceBValue(); m_OriginalDirections = OriginalDirectionscontainer; m_Directions = DirectionsContainer; m_MeasurementFrame = orig.GetMeasurementFrame(); ApplyMeasurementFrame(); UpdateBValueMap(); InitializeFromVectorImage(); // m_VectorImage = duplicator->GetOutput(); // m_B_Value = orig.GetB_Value(); // m_MeasurementFrame = orig.GetMeasurementFrame(); // m_Directions = copyInContainer; // m_OriginalDirections = copyInContainer; } template mitk::DiffusionImage::~DiffusionImage() { // Remove Observer for m_Directions //RemoveObserver(); } template void mitk::DiffusionImage ::InitializeFromVectorImage() { if(!m_VectorImage || !m_Directions || m_B_Value==-1.0) { MITK_INFO << "DiffusionImage could not be initialized. Set all members first!" << std::endl; // Fehjlermedlung updaten return; } // find bzero index int firstZeroIndex = -1; for(GradientDirectionContainerType::ConstIterator it = m_Directions->Begin(); it != m_Directions->End(); ++it) { firstZeroIndex++; GradientDirectionType g = it.Value(); if(g[0] == 0 && g[1] == 0 && g[2] == 0 ) break; } typedef itk::Image ImgType; typename ImgType::Pointer img = ImgType::New(); img->SetSpacing( m_VectorImage->GetSpacing() ); // Set the image spacing img->SetOrigin( m_VectorImage->GetOrigin() ); // Set the image origin img->SetDirection( m_VectorImage->GetDirection() ); // Set the image direction img->SetLargestPossibleRegion( m_VectorImage->GetLargestPossibleRegion()); img->SetBufferedRegion( m_VectorImage->GetLargestPossibleRegion() ); img->Allocate(); int vecLength = m_VectorImage->GetVectorLength(); InitializeByItk( img.GetPointer(), 1, vecLength ); itk::ImageRegionIterator itw (img, img->GetLargestPossibleRegion() ); itw.GoToBegin(); itk::ImageRegionConstIterator itr (m_VectorImage, m_VectorImage->GetLargestPossibleRegion() ); itr.GoToBegin(); while(!itr.IsAtEnd()) { itw.Set(itr.Get().GetElement(firstZeroIndex)); ++itr; ++itw; } // init SetImportVolume(img->GetBufferPointer()); m_DisplayIndex = firstZeroIndex; MITK_INFO << "Diffusion-Image successfully initialized."; } - -template -typename mitk::DiffusionImage< TPixelType >* mitk::DiffusionImage::GetDeepCopy() -{ - typename DwiDuplicatorType::Pointer duplicator = DwiDuplicatorType::New(); - duplicator->SetInputImage(this->GetVectorImage()); - duplicator->Update(); - - mitk::DiffusionImage< TPixelType >::Pointer newDwi = mitk::DiffusionImage< TPixelType >::New(); - newDwi->SetVectorImage( duplicator->GetOutput() ); - newDwi->SetReferenceBValue(this->GetReferenceBValue()); - - GradientDirectionContainerType::Pointer newDirectionContainer = GradientDirectionContainerType::New(); - for (unsigned int i=0; iGetDirectionsWithoutMeasurementFrame()->Size(); i++) - { - newDirectionContainer->InsertElement(i, this->GetDirectionsWithoutMeasurementFrame()->GetElement(i)); - } - newDwi->SetDirections(newDirectionContainer); - newDwi->SetMeasurementFrame(this->GetMeasurementFrame()); - newDwi->InitializeFromVectorImage(); - - return newDwi; -} - template void mitk::DiffusionImage ::SetDisplayIndexForRendering(int displayIndex) { int index = displayIndex; int vecLength = m_VectorImage->GetVectorLength(); index = index > vecLength-1 ? vecLength-1 : index; if( m_DisplayIndex != index ) { typedef itk::Image ImgType; typename ImgType::Pointer img = ImgType::New(); CastToItkImage(this, img); itk::ImageRegionIterator itw (img, img->GetLargestPossibleRegion() ); itw.GoToBegin(); itk::ImageRegionConstIterator itr (m_VectorImage, m_VectorImage->GetLargestPossibleRegion() ); itr.GoToBegin(); while(!itr.IsAtEnd()) { itw.Set(itr.Get().GetElement(index)); ++itr; ++itw; } } m_DisplayIndex = index; } template bool mitk::DiffusionImage::AreAlike(GradientDirectionType g1, GradientDirectionType g2, double precision) { GradientDirectionType diff = g1 - g2; GradientDirectionType diff2 = g1 + g2; return diff.two_norm() < precision || diff2.two_norm() < precision; } template mitk::DiffusionImage::GradientDirectionContainerType::Pointer mitk::DiffusionImage::CalcAveragedDirectionSet(double precision, GradientDirectionContainerType::Pointer directions) { // save old and construct new direction container GradientDirectionContainerType::Pointer newDirections = GradientDirectionContainerType::New(); // fill new direction container for(GradientDirectionContainerType::ConstIterator gdcitOld = directions->Begin(); gdcitOld != directions->End(); ++gdcitOld) { // already exists? bool found = false; for(GradientDirectionContainerType::ConstIterator gdcitNew = newDirections->Begin(); gdcitNew != newDirections->End(); ++gdcitNew) { if(AreAlike(gdcitNew.Value(), gdcitOld.Value(), precision)) { found = true; break; } } // if not found, add it to new container if(!found) { newDirections->push_back(gdcitOld.Value()); } } return newDirections; } template void mitk::DiffusionImage::AverageRedundantGradients(double precision) { GradientDirectionContainerType::Pointer newDirs = CalcAveragedDirectionSet(precision, m_Directions); GradientDirectionContainerType::Pointer newOriginalDirs = CalcAveragedDirectionSet(precision, m_OriginalDirections); // if sizes equal, we do not need to do anything in this function if(m_Directions->size() == newDirs->size() || m_OriginalDirections->size() == newOriginalDirs->size()) return; GradientDirectionContainerType::Pointer oldDirections = m_OriginalDirections; m_Directions = newDirs; m_OriginalDirections = newOriginalDirs; // new image typename ImageType::Pointer oldImage = m_VectorImage; m_VectorImage = ImageType::New(); m_VectorImage->SetSpacing( oldImage->GetSpacing() ); // Set the image spacing m_VectorImage->SetOrigin( oldImage->GetOrigin() ); // Set the image origin m_VectorImage->SetDirection( oldImage->GetDirection() ); // Set the image direction m_VectorImage->SetLargestPossibleRegion( oldImage->GetLargestPossibleRegion() ); m_VectorImage->SetVectorLength( m_Directions->size() ); m_VectorImage->SetBufferedRegion( oldImage->GetLargestPossibleRegion() ); m_VectorImage->Allocate(); // average image data that corresponds to identical directions itk::ImageRegionIterator< ImageType > newIt(m_VectorImage, m_VectorImage->GetLargestPossibleRegion()); newIt.GoToBegin(); itk::ImageRegionIterator< ImageType > oldIt(oldImage, oldImage->GetLargestPossibleRegion()); oldIt.GoToBegin(); // initial new value of voxel typename ImageType::PixelType newVec; newVec.SetSize(m_Directions->size()); newVec.AllocateElements(m_Directions->size()); std::vector > dirIndices; for(GradientDirectionContainerType::ConstIterator gdcitNew = m_Directions->Begin(); gdcitNew != m_Directions->End(); ++gdcitNew) { dirIndices.push_back(std::vector(0)); for(GradientDirectionContainerType::ConstIterator gdcitOld = oldDirections->Begin(); gdcitOld != oldDirections->End(); ++gdcitOld) { if(AreAlike(gdcitNew.Value(), gdcitOld.Value(), precision)) { //MITK_INFO << gdcitNew.Value() << " " << gdcitOld.Value(); dirIndices[gdcitNew.Index()].push_back(gdcitOld.Index()); } } } //int ind1 = -1; while(!newIt.IsAtEnd()) { // progress //typename ImageType::IndexType ind = newIt.GetIndex(); //ind1 = ind.m_Index[2]; // init new vector with zeros newVec.Fill(0.0); // the old voxel value with duplicates typename ImageType::PixelType oldVec = oldIt.Get(); for(unsigned int i=0; i void mitk::DiffusionImage::ApplyMeasurementFrame() { if(m_OriginalDirections.IsNull()) { // original direction container was not set return; } m_Directions = GradientDirectionContainerType::New(); int c = 0; for(GradientDirectionContainerType::ConstIterator gdcit = m_OriginalDirections->Begin(); gdcit != m_OriginalDirections->End(); ++gdcit) { vnl_vector vec = gdcit.Value(); vec = vec.pre_multiply(m_MeasurementFrame); m_Directions->InsertElement(c, vec); c++; } } template void mitk::DiffusionImage::UpdateBValueMap() { if(!m_B_ValueMap.empty()) m_B_ValueMap.clear(); GradientDirectionContainerType::ConstIterator gdcit; for( gdcit = this->m_Directions->Begin(); gdcit != this->m_Directions->End(); ++gdcit) m_B_ValueMap[GetB_Value(gdcit.Index())].push_back(gdcit.Index()); } template float mitk::DiffusionImage::GetB_Value(unsigned int i) { if(i > m_Directions->Size()-1) return -1; if(m_Directions->ElementAt(i).one_norm() <= 0.0) { return 0; } else { double twonorm = m_Directions->ElementAt(i).two_norm(); double bval = m_B_Value*twonorm*twonorm; if (bval<0) bval = ceil(bval - 0.5); else bval = floor(bval + 0.5); return bval; } } template void mitk::DiffusionImage::SetDirections( GradientDirectionContainerType::Pointer directions ) { this->m_OriginalDirections = directions; ApplyMeasurementFrame(); UpdateBValueMap(); } template void mitk::DiffusionImage::SetDirections(const std::vector > &directions) { GradientDirectionContainerType::Pointer tempContainer = GradientDirectionContainerType::New(); for(unsigned int i=0; iInsertElement( i, directions[i].GetVnlVector() ); SetDirections(tempContainer); } template typename mitk::DiffusionImage::MeasurementFrameType mitk::DiffusionImage::GetMeasurementFrame() const { return m_MeasurementFrame; } template void mitk::DiffusionImage::SetMeasurementFrame( const MeasurementFrameType & mFrame ) { m_MeasurementFrame = mFrame; ApplyMeasurementFrame(); } template bool mitk::DiffusionImage::GetNumberOfBValues() const { return m_B_ValueMap.size(); } template const typename mitk::DiffusionImage::BValueMap & mitk::DiffusionImage::GetBValueMap() const { return m_B_ValueMap; } template float mitk::DiffusionImage::GetReferenceBValue() const { return m_B_Value; } template void mitk::DiffusionImage::SetReferenceBValue( float val ) { m_B_Value = val; } template typename mitk::DiffusionImage::GradientDirectionContainerTypePointer mitk::DiffusionImage::GetDirections() const { return m_Directions; } template typename mitk::DiffusionImage::GradientDirectionContainerTypePointer mitk::DiffusionImage::GetDirectionsWithoutMeasurementFrame() const { return m_OriginalDirections; } template typename itk::VectorImage::Pointer mitk::DiffusionImage::GetVectorImage() { return m_VectorImage; } template const typename itk::VectorImage::Pointer mitk::DiffusionImage::GetVectorImage() const { return m_VectorImage; } template void mitk::DiffusionImage::SetVectorImage(typename ImageType::Pointer image ) { m_VectorImage = image; } template inline bool mitk::Equal(const mitk::DiffusionImage& leftHandSide, const mitk::DiffusionImage &rightHandSide, ScalarType eps, bool verbose ) { bool returnValue = true; if(leftHandSide.GetReferenceBValue() != rightHandSide.GetReferenceBValue()) { if(verbose) MITK_INFO << "[( DiffusionImage )] Reference BValue is not Equal."; returnValue = false; } if(leftHandSide.GetMeasurementFrame() != rightHandSide.GetMeasurementFrame()) { if(verbose) MITK_INFO << "[( DiffusionImage )] MeasurementFrame is not Equal."; returnValue = false; } if(leftHandSide.GetBValueMap().size() != rightHandSide.GetBValueMap().size()) { if(verbose) MITK_INFO << "[( DiffusionImage )] BValue Map: Size is not Equal."; returnValue = false; } if(leftHandSide.GetNumberOfBValues() != rightHandSide.GetNumberOfBValues()) { if(verbose) MITK_INFO << "[( DiffusionImage )] BValue Map (GetNumberOfBValues): Size is not Equal."; returnValue = false; } // Slow testing area const mitk::Image* img1 = dynamic_cast(&leftHandSide); const mitk::Image* img2 = dynamic_cast(&rightHandSide); if(mitk::Equal(*img1,*img2,eps,verbose) == false) { if(verbose) MITK_INFO << "[( DiffusionImage )] Base-Class (mitk::Image) is not Equal."; returnValue = false; } { typename mitk::DiffusionImage::GradientDirectionContainerType::Iterator lhsIt = leftHandSide.GetDirectionsWithoutMeasurementFrame()->Begin(); typename mitk::DiffusionImage::GradientDirectionContainerType::Iterator lhsItEnd = leftHandSide.GetDirectionsWithoutMeasurementFrame()->End(); typename mitk::DiffusionImage::GradientDirectionContainerType::Iterator rhsIt = rightHandSide.GetDirectionsWithoutMeasurementFrame()->Begin(); for(;lhsIt != lhsItEnd;) { bool vectorNotEqual = false; for(unsigned int i = 0 ; i < lhsIt.Value().size(); i++) vectorNotEqual |= !mitk::Equal(lhsIt.Value().get(i),rhsIt.Value().get(i),0.0001); if(vectorNotEqual) { if(verbose) MITK_INFO << "[( DiffusionImage )] Original GradientDirections are not Equal."; returnValue = false; break; } ++rhsIt; ++lhsIt; } } { typename mitk::DiffusionImage::GradientDirectionContainerType::Iterator lhsIt = leftHandSide.GetDirections()->Begin(); typename mitk::DiffusionImage::GradientDirectionContainerType::Iterator lhsItEnd = leftHandSide.GetDirections()->End(); typename mitk::DiffusionImage::GradientDirectionContainerType::Iterator rhsIt = rightHandSide.GetDirections()->Begin(); for(;lhsIt != lhsItEnd;) { bool vectorNotEqual = false; for(unsigned int i = 0 ; i < lhsIt.Value().size(); i++) vectorNotEqual |= !mitk::Equal(lhsIt.Value().get(i),rhsIt.Value().get(i),0.0001); if(vectorNotEqual) { if(verbose) MITK_INFO << "[( DiffusionImage )] GradientDirections are not Equal."; returnValue = false; break; } ++rhsIt; ++lhsIt; } } { typename mitk::DiffusionImage::BValueMap rhsMap = rightHandSide.GetBValueMap(); typename mitk::DiffusionImage::BValueMap lhsMap = leftHandSide.GetBValueMap(); typename mitk::DiffusionImage::BValueMap::const_iterator lhsIt = lhsMap.begin(); typename mitk::DiffusionImage::BValueMap::const_iterator lhsItEnd = lhsMap.end(); typename mitk::DiffusionImage::BValueMap::const_iterator rhsIt = rhsMap.begin(); for(;lhsIt != lhsItEnd;) { if(lhsIt->first != rhsIt->first) { if(verbose) MITK_INFO << "[( DiffusionImage )] BValue Map: lhsKey " << lhsIt->first << " != rhsKey " << rhsIt->first << " is not Equal."; returnValue = false; break; } if(lhsIt->second.size() != rhsIt->second.size()) { if(verbose) MITK_INFO << "[( DiffusionImage )] BValue Map: Indices vector size is not Equal. (Key: " << lhsIt->first <<")"; returnValue = false; break; } typename mitk::DiffusionImage::IndicesVector::const_iterator lhsIndVecIt = lhsIt->second.begin(); typename mitk::DiffusionImage::IndicesVector::const_iterator lhsIndVecItEnd = lhsIt->second.end(); typename mitk::DiffusionImage::IndicesVector::const_iterator rhsIndVecIt = rhsIt->second.begin(); for(;lhsIndVecIt != lhsIndVecItEnd;) { if(*lhsIndVecIt != *rhsIndVecIt) { if(verbose) MITK_INFO << "[( DiffusionImage )] BValue Map: Indices are not Equal. (Key: " << lhsIt->first <<")"; returnValue = false; break; } ++rhsIndVecIt; ++lhsIndVecIt; } lhsIt++; rhsIt++; } } try{ itk::ImageRegionIterator< itk::VectorImage< TPixelType, 3 > > it1(leftHandSide.GetVectorImage(), leftHandSide.GetVectorImage()->GetLargestPossibleRegion()); itk::ImageRegionIterator< itk::VectorImage< TPixelType, 3 > > it2(rightHandSide.GetVectorImage(), rightHandSide.GetVectorImage()->GetLargestPossibleRegion()); while(!it1.IsAtEnd()) { if (it1.Get()!=it2.Get()){ if(verbose) MITK_INFO << "[( DiffusionImage )] Capsulated itk::VectorImage is not Equal."; returnValue = false; break; } ++it1; ++it2; } } catch(...) { if(verbose) MITK_INFO << "[( DiffusionImage )] Comparision of itk Vector image abort by exception."; returnValue = false; } return returnValue; } diff --git a/Plugins/org.mitk.gui.qt.diffusionimaging/src/internal/QmitkPreprocessingView.cpp b/Plugins/org.mitk.gui.qt.diffusionimaging/src/internal/QmitkPreprocessingView.cpp index b5d2fb1307..9f6843afac 100644 --- a/Plugins/org.mitk.gui.qt.diffusionimaging/src/internal/QmitkPreprocessingView.cpp +++ b/Plugins/org.mitk.gui.qt.diffusionimaging/src/internal/QmitkPreprocessingView.cpp @@ -1,1389 +1,1390 @@ /*=================================================================== 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. ===================================================================*/ //#define MBILOG_ENABLE_DEBUG #include "QmitkPreprocessingView.h" #include "mitkDiffusionImagingConfigure.h" // qt includes #include // itk includes #include "itkTimeProbe.h" #include "itkB0ImageExtractionImageFilter.h" #include "itkB0ImageExtractionToSeparateImageFilter.h" #include "itkBrainMaskExtractionImageFilter.h" #include "itkCastImageFilter.h" #include "itkVectorContainer.h" #include #include #include #include // Multishell includes #include // Multishell Functors #include #include #include #include // mitk includes #include "QmitkDataStorageComboBox.h" #include "QmitkStdMultiWidget.h" #include "mitkProgressBar.h" #include "mitkStatusBar.h" #include "mitkNodePredicateDataType.h" #include "mitkProperties.h" #include "mitkVtkResliceInterpolationProperty.h" #include "mitkLookupTable.h" #include "mitkLookupTableProperty.h" #include "mitkTransferFunction.h" #include "mitkTransferFunctionProperty.h" #include "mitkDataNodeObject.h" #include "mitkOdfNormalizationMethodProperty.h" #include "mitkOdfScaleByProperty.h" #include #include #include #include #include #include #include #include #include #include #include const std::string QmitkPreprocessingView::VIEW_ID = "org.mitk.views.preprocessing"; #define DI_INFO MITK_INFO("DiffusionImaging") typedef float TTensorPixelType; QmitkPreprocessingView::QmitkPreprocessingView() : QmitkFunctionality(), m_Controls(NULL), m_MultiWidget(NULL), m_DiffusionImage(NULL) { } QmitkPreprocessingView::~QmitkPreprocessingView() { } void QmitkPreprocessingView::CreateQtPartControl(QWidget *parent) { if (!m_Controls) { // create GUI widgets m_Controls = new Ui::QmitkPreprocessingViewControls; m_Controls->setupUi(parent); this->CreateConnections(); m_Controls->m_MeasurementFrameTable->horizontalHeader()->setResizeMode(QHeaderView::Stretch); m_Controls->m_MeasurementFrameTable->verticalHeader()->setResizeMode(QHeaderView::Stretch); m_Controls->m_DirectionMatrixTable->horizontalHeader()->setResizeMode(QHeaderView::Stretch); m_Controls->m_DirectionMatrixTable->verticalHeader()->setResizeMode(QHeaderView::Stretch); } } void QmitkPreprocessingView::StdMultiWidgetAvailable (QmitkStdMultiWidget &stdMultiWidget) { m_MultiWidget = &stdMultiWidget; } void QmitkPreprocessingView::StdMultiWidgetNotAvailable() { m_MultiWidget = NULL; } void QmitkPreprocessingView::CreateConnections() { if ( m_Controls ) { m_Controls->m_ExtractBrainMask->setVisible(false); m_Controls->m_BrainMaskIterationsBox->setVisible(false); m_Controls->m_ResampleIntFrame->setVisible(false); connect( (QObject*)(m_Controls->m_ButtonAverageGradients), SIGNAL(clicked()), this, SLOT(AverageGradients()) ); connect( (QObject*)(m_Controls->m_ButtonExtractB0), SIGNAL(clicked()), this, SLOT(ExtractB0()) ); connect( (QObject*)(m_Controls->m_ModifyMeasurementFrame), SIGNAL(clicked()), this, SLOT(DoApplyMesurementFrame()) ); connect( (QObject*)(m_Controls->m_ReduceGradientsButton), SIGNAL(clicked()), this, SLOT(DoReduceGradientDirections()) ); connect( (QObject*)(m_Controls->m_ShowGradientsButton), SIGNAL(clicked()), this, SLOT(DoShowGradientDirections()) ); connect( (QObject*)(m_Controls->m_MirrorGradientToHalfSphereButton), SIGNAL(clicked()), this, SLOT(DoHalfSphereGradientDirections()) ); connect( (QObject*)(m_Controls->m_MergeDwisButton), SIGNAL(clicked()), this, SLOT(MergeDwis()) ); connect( (QObject*)(m_Controls->m_ProjectSignalButton), SIGNAL(clicked()), this, SLOT(DoProjectSignal()) ); connect( (QObject*)(m_Controls->m_B_ValueMap_Rounder_SpinBox), SIGNAL(valueChanged(int)), this, SLOT(UpdateDwiBValueMapRounder(int))); connect( (QObject*)(m_Controls->m_CreateLengthCorrectedDwi), SIGNAL(clicked()), this, SLOT(DoLengthCorrection()) ); connect( (QObject*)(m_Controls->m_CalcAdcButton), SIGNAL(clicked()), this, SLOT(DoAdcCalculation()) ); connect( (QObject*)(m_Controls->m_NormalizeImageValuesButton), SIGNAL(clicked()), this, SLOT(DoDwiNormalization()) ); connect( (QObject*)(m_Controls->m_ModifyDirection), SIGNAL(clicked()), this, SLOT(DoApplyDirectionMatrix()) ); connect( (QObject*)(m_Controls->m_ResampleImageButton), SIGNAL(clicked()), this, SLOT(DoResampleImage()) ); connect( (QObject*)(m_Controls->m_ResampleTypeBox), SIGNAL(currentIndexChanged(int)), this, SLOT(DoUpdateInterpolationGui(int)) ); // connect( (QObject*)(m_Controls->m_ExtractBrainMask), SIGNAL(clicked()), this, SLOT(DoExtractBrainMask()) ); } } void QmitkPreprocessingView::DoUpdateInterpolationGui(int i) { switch (i) { case 0: { m_Controls->m_ResampleIntFrame->setVisible(false); m_Controls->m_ResampleDoubleFrame->setVisible(true); m_Controls->m_ResampleDoubleX->setValue(2); m_Controls->m_ResampleDoubleY->setValue(2); m_Controls->m_ResampleDoubleZ->setValue(2); break; } case 1: { m_Controls->m_ResampleIntFrame->setVisible(false); m_Controls->m_ResampleDoubleFrame->setVisible(true); if (m_DiffusionImage.IsNotNull()) { ItkDwiType::Pointer itkDwi = m_DiffusionImage->GetVectorImage(); m_Controls->m_ResampleDoubleX->setValue(itkDwi->GetSpacing()[0]); m_Controls->m_ResampleDoubleY->setValue(itkDwi->GetSpacing()[1]); m_Controls->m_ResampleDoubleZ->setValue(itkDwi->GetSpacing()[2]); } else if (m_SelectedImage.IsNotNull()) { mitk::BaseGeometry* geom = m_SelectedImage->GetGeometry(); m_Controls->m_ResampleDoubleX->setValue(geom->GetSpacing()[0]); m_Controls->m_ResampleDoubleY->setValue(geom->GetSpacing()[1]); m_Controls->m_ResampleDoubleZ->setValue(geom->GetSpacing()[2]); } break; } case 2: { m_Controls->m_ResampleIntFrame->setVisible(true); m_Controls->m_ResampleDoubleFrame->setVisible(false); if (m_DiffusionImage.IsNotNull()) { ItkDwiType::Pointer itkDwi = m_DiffusionImage->GetVectorImage(); m_Controls->m_ResampleIntX->setValue(itkDwi->GetLargestPossibleRegion().GetSize(0)); m_Controls->m_ResampleIntY->setValue(itkDwi->GetLargestPossibleRegion().GetSize(1)); m_Controls->m_ResampleIntZ->setValue(itkDwi->GetLargestPossibleRegion().GetSize(2)); } else if (m_SelectedImage.IsNotNull()) { mitk::BaseGeometry* geom = m_SelectedImage->GetGeometry(); m_Controls->m_ResampleIntX->setValue(geom->GetExtent(0)); m_Controls->m_ResampleIntY->setValue(geom->GetExtent(1)); m_Controls->m_ResampleIntZ->setValue(geom->GetExtent(2)); } break; } default: { m_Controls->m_ResampleDoubleX->setValue(2); m_Controls->m_ResampleDoubleY->setValue(2); m_Controls->m_ResampleDoubleZ->setValue(2); m_Controls->m_ResampleIntFrame->setVisible(false); m_Controls->m_ResampleDoubleFrame->setVisible(true); } } } void QmitkPreprocessingView::DoExtractBrainMask() { // if (m_SelectedImage.IsNull()) // return; // typedef itk::Image ShortImageType; // ShortImageType::Pointer itkImage = ShortImageType::New(); // mitk::CastToItkImage(m_SelectedImage, itkImage); // typedef itk::BrainMaskExtractionImageFilter< unsigned char > FilterType; // FilterType::Pointer filter = FilterType::New(); // filter->SetInput(itkImage); // filter->SetMaxNumIterations(m_Controls->m_BrainMaskIterationsBox->value()); // filter->Update(); // mitk::Image::Pointer image = mitk::Image::New(); // image->InitializeByItk( filter->GetOutput() ); // image->SetVolume( filter->GetOutput()->GetBufferPointer() ); // mitk::DataNode::Pointer imageNode = mitk::DataNode::New(); // imageNode->SetData( image ); // imageNode->SetName("BRAINMASK"); // GetDefaultDataStorage()->Add(imageNode); } void QmitkPreprocessingView::DoResampleImage() { if (m_DiffusionImage.IsNotNull()) { typedef itk::ResampleDwiImageFilter< short > ResampleFilter; typename ResampleFilter::Pointer resampler = ResampleFilter::New(); resampler->SetInput(m_DiffusionImage->GetVectorImage()); switch (m_Controls->m_ResampleTypeBox->currentIndex()) { case 0: { itk::Vector< double, 3 > samplingFactor; samplingFactor[0] = m_Controls->m_ResampleDoubleX->value(); samplingFactor[1] = m_Controls->m_ResampleDoubleY->value(); samplingFactor[2] = m_Controls->m_ResampleDoubleZ->value(); resampler->SetSamplingFactor(samplingFactor); break; } case 1: { itk::Vector< double, 3 > newSpacing; newSpacing[0] = m_Controls->m_ResampleDoubleX->value(); newSpacing[1] = m_Controls->m_ResampleDoubleY->value(); newSpacing[2] = m_Controls->m_ResampleDoubleZ->value(); resampler->SetNewSpacing(newSpacing); break; } case 2: { itk::ImageRegion<3> newRegion; newRegion.SetSize(0, m_Controls->m_ResampleIntX->value()); newRegion.SetSize(1, m_Controls->m_ResampleIntY->value()); newRegion.SetSize(2, m_Controls->m_ResampleIntZ->value()); resampler->SetNewImageSize(newRegion); break; } default: { MITK_WARN << "Unknown resampling parameters!"; return; } } QString outAdd; switch (m_Controls->m_InterpolatorBox->currentIndex()) { case 0: { resampler->SetInterpolation(ResampleFilter::Interpolate_NearestNeighbour); outAdd = "NearestNeighbour"; break; } case 1: { resampler->SetInterpolation(ResampleFilter::Interpolate_Linear); outAdd = "Linear"; break; } case 2: { resampler->SetInterpolation(ResampleFilter::Interpolate_BSpline); outAdd = "BSpline"; break; } case 3: { resampler->SetInterpolation(ResampleFilter::Interpolate_WindowedSinc); outAdd = "WindowedSinc"; break; } default: { resampler->SetInterpolation(ResampleFilter::Interpolate_NearestNeighbour); outAdd = "NearestNeighbour"; } } resampler->Update(); typedef mitk::DiffusionImage DiffusionImageType; DiffusionImageType::Pointer image = DiffusionImageType::New(); image->SetVectorImage( resampler->GetOutput() ); image->SetReferenceBValue( m_DiffusionImage->GetReferenceBValue() ); image->SetDirections( m_DiffusionImage->GetDirections() ); image->InitializeFromVectorImage(); mitk::DataNode::Pointer imageNode = mitk::DataNode::New(); imageNode->SetData( image ); QString name = m_SelectedDiffusionNodes.back()->GetName().c_str(); imageNode->SetName((name+"_resampled_"+outAdd).toStdString().c_str()); GetDefaultDataStorage()->Add(imageNode, m_SelectedDiffusionNodes.back()); } else if (m_SelectedImage.IsNotNull()) { AccessFixedDimensionByItk(m_SelectedImage, TemplatedResampleImage,3); } } template < typename TPixel, unsigned int VImageDimension > void QmitkPreprocessingView::TemplatedResampleImage( itk::Image* itkImage) { itk::Vector< double, 3 > newSpacing; itk::ImageRegion<3> newRegion; switch (m_Controls->m_ResampleTypeBox->currentIndex()) { case 0: { itk::Vector< double, 3 > sampling; sampling[0] = m_Controls->m_ResampleDoubleX->value(); sampling[1] = m_Controls->m_ResampleDoubleY->value(); sampling[2] = m_Controls->m_ResampleDoubleZ->value(); newSpacing = itkImage->GetSpacing(); newSpacing[0] /= sampling[0]; newSpacing[1] /= sampling[1]; newSpacing[2] /= sampling[2]; newRegion = itkImage->GetLargestPossibleRegion(); newRegion.SetSize(0, newRegion.GetSize(0)*sampling[0]); newRegion.SetSize(1, newRegion.GetSize(1)*sampling[1]); newRegion.SetSize(2, newRegion.GetSize(2)*sampling[2]); break; } case 1: { newSpacing[0] = m_Controls->m_ResampleDoubleX->value(); newSpacing[1] = m_Controls->m_ResampleDoubleY->value(); newSpacing[2] = m_Controls->m_ResampleDoubleZ->value(); itk::Vector< double, 3 > oldSpacing = itkImage->GetSpacing(); itk::Vector< double, 3 > sampling; sampling[0] = oldSpacing[0]/newSpacing[0]; sampling[1] = oldSpacing[1]/newSpacing[1]; sampling[2] = oldSpacing[2]/newSpacing[2]; newRegion = itkImage->GetLargestPossibleRegion(); newRegion.SetSize(0, newRegion.GetSize(0)*sampling[0]); newRegion.SetSize(1, newRegion.GetSize(1)*sampling[1]); newRegion.SetSize(2, newRegion.GetSize(2)*sampling[2]); break; } case 2: { newRegion.SetSize(0, m_Controls->m_ResampleIntX->value()); newRegion.SetSize(1, m_Controls->m_ResampleIntY->value()); newRegion.SetSize(2, m_Controls->m_ResampleIntZ->value()); itk::ImageRegion<3> oldRegion = itkImage->GetLargestPossibleRegion(); itk::Vector< double, 3 > sampling; sampling[0] = (double)newRegion.GetSize(0)/oldRegion.GetSize(0); sampling[1] = (double)newRegion.GetSize(1)/oldRegion.GetSize(1); sampling[2] = (double)newRegion.GetSize(2)/oldRegion.GetSize(2); newSpacing = itkImage->GetSpacing(); newSpacing[0] /= sampling[0]; newSpacing[1] /= sampling[1]; newSpacing[2] /= sampling[2]; break; } default: { MITK_WARN << "Unknown resampling parameters!"; return; } } itk::Point origin = itkImage->GetOrigin(); origin[0] -= itkImage->GetSpacing()[0]/2; origin[1] -= itkImage->GetSpacing()[1]/2; origin[2] -= itkImage->GetSpacing()[2]/2; origin[0] += newSpacing[0]/2; origin[1] += newSpacing[1]/2; origin[2] += newSpacing[2]/2; typedef itk::Image ImageType; typename ImageType::Pointer outImage = ImageType::New(); outImage->SetSpacing( newSpacing ); outImage->SetOrigin( origin ); outImage->SetDirection( itkImage->GetDirection() ); outImage->SetLargestPossibleRegion( newRegion ); outImage->SetBufferedRegion( newRegion ); outImage->SetRequestedRegion( newRegion ); outImage->Allocate(); typedef itk::ResampleImageFilter ResampleFilter; typename ResampleFilter::Pointer resampler = ResampleFilter::New(); resampler->SetInput(itkImage); resampler->SetOutputParametersFromImage(outImage); QString outAdd; switch (m_Controls->m_InterpolatorBox->currentIndex()) { case 0: { typename itk::NearestNeighborInterpolateImageFunction::Pointer interp = itk::NearestNeighborInterpolateImageFunction::New(); resampler->SetInterpolator(interp); outAdd = "NearestNeighbour"; break; } case 1: { typename itk::LinearInterpolateImageFunction::Pointer interp = itk::LinearInterpolateImageFunction::New(); resampler->SetInterpolator(interp); outAdd = "Linear"; break; } case 2: { typename itk::BSplineInterpolateImageFunction::Pointer interp = itk::BSplineInterpolateImageFunction::New(); resampler->SetInterpolator(interp); outAdd = "BSpline"; break; } case 3: { typename itk::WindowedSincInterpolateImageFunction::Pointer interp = itk::WindowedSincInterpolateImageFunction::New(); resampler->SetInterpolator(interp); outAdd = "WindowedSinc"; break; } default: { typename itk::NearestNeighborInterpolateImageFunction::Pointer interp = itk::NearestNeighborInterpolateImageFunction::New(); resampler->SetInterpolator(interp); outAdd = "NearestNeighbour"; } } resampler->Update(); mitk::Image::Pointer image = mitk::Image::New(); image->InitializeByItk( resampler->GetOutput() ); image->SetVolume( resampler->GetOutput()->GetBufferPointer() ); mitk::DataNode::Pointer imageNode = mitk::DataNode::New(); imageNode->SetData( image ); QString name = m_SelectedImageNode->GetName().c_str(); imageNode->SetName((name+"_resampled_"+outAdd).toStdString().c_str()); GetDefaultDataStorage()->Add(imageNode, m_SelectedImageNode); } void QmitkPreprocessingView::DoApplyDirectionMatrix() { if (m_DiffusionImage.IsNotNull()) { - MitkDwiType::Pointer newDwi = m_DiffusionImage->GetDeepCopy(); + MitkDwiType::Pointer newDwi = m_DiffusionImage->Clone(); ItkDwiType::DirectionType newDirection; for (int r=0; r<3; r++) for (int c=0; c<3; c++) { QTableWidgetItem* item = m_Controls->m_DirectionMatrixTable->item(r,c); if (!item) return; newDirection[r][c] = item->text().toDouble(); } ItkDwiType::Pointer itkDwi = newDwi->GetVectorImage(); typedef mitk::DiffusionImage MitkDwiType; vnl_matrix_fixed< double,3,3 > oldInverseDirection = itkDwi->GetDirection().GetInverse(); MitkDwiType::GradientDirectionContainerType::Pointer oldGradients = m_DiffusionImage->GetDirectionsWithoutMeasurementFrame(); MitkDwiType::GradientDirectionContainerType::Pointer newGradients = MitkDwiType::GradientDirectionContainerType::New(); for (unsigned int i=0; iSize(); i++) { MitkDwiType::GradientDirectionType g = oldGradients->GetElement(i); double mag = g.magnitude(); MitkDwiType::GradientDirectionType newG = oldInverseDirection*g; newG = newDirection.GetVnlMatrix()*newG; newG.normalize(); newGradients->InsertElement(i, newG*mag); } newDwi->SetDirections(newGradients); itkDwi->SetDirection(newDirection); newDwi->InitializeFromVectorImage(); mitk::DataNode::Pointer imageNode = mitk::DataNode::New(); imageNode->SetData( newDwi ); QString name = m_SelectedDiffusionNodes.back()->GetName().c_str(); imageNode->SetName((name+"_newdirection").toStdString().c_str()); GetDefaultDataStorage()->Add(imageNode, m_SelectedDiffusionNodes.back()); mitk::RenderingManager::GetInstance()->InitializeViews( imageNode->GetData()->GetTimeGeometry(), mitk::RenderingManager::REQUEST_UPDATE_ALL, true ); mitk::RenderingManager::GetInstance()->RequestUpdateAll(); } else if (m_SelectedImage.IsNotNull()) { AccessFixedDimensionByItk(m_SelectedImage, TemplatedApplyRotation,3); } } template < typename TPixel, unsigned int VImageDimension > void QmitkPreprocessingView::TemplatedApplyRotation( itk::Image* itkImage) { ItkDwiType::DirectionType newDirection; for (int r=0; r<3; r++) for (int c=0; c<3; c++) { QTableWidgetItem* item = m_Controls->m_DirectionMatrixTable->item(r,c); if (!item) return; newDirection[r][c] = item->text().toDouble(); } typedef itk::Image ImageType; typename ImageType::Pointer newImage = ImageType::New(); newImage->SetSpacing( itkImage->GetSpacing() ); newImage->SetOrigin( itkImage->GetOrigin() ); newImage->SetDirection( newDirection ); newImage->SetLargestPossibleRegion( itkImage->GetLargestPossibleRegion() ); newImage->SetBufferedRegion( itkImage->GetLargestPossibleRegion() ); newImage->SetRequestedRegion( itkImage->GetLargestPossibleRegion() ); newImage->Allocate(); newImage->FillBuffer(0); itk::ImageRegionIterator< itk::Image > it(itkImage, itkImage->GetLargestPossibleRegion()); while(!it.IsAtEnd()) { newImage->SetPixel(it.GetIndex(), it.Get()); ++it; } mitk::Image::Pointer newMitkImage = mitk::Image::New(); newMitkImage->InitializeByItk(newImage.GetPointer()); newMitkImage->SetVolume(newImage->GetBufferPointer()); mitk::DataNode::Pointer imageNode = mitk::DataNode::New(); imageNode->SetData( newMitkImage ); QString name = m_SelectedImageNode->GetName().c_str(); imageNode->SetName((name+"_newdirection").toStdString().c_str()); GetDefaultDataStorage()->Add(imageNode, m_SelectedImageNode); mitk::RenderingManager::GetInstance()->InitializeViews( m_SelectedImageNode->GetData()->GetTimeGeometry(), mitk::RenderingManager::REQUEST_UPDATE_ALL, true ); mitk::RenderingManager::GetInstance()->RequestUpdateAll(); } void QmitkPreprocessingView::DoProjectSignal() { switch(m_Controls->m_ProjectionMethodBox->currentIndex()) { case 0: DoADCAverage(); break; case 1: DoAKCFit(); break; case 2: DoBiExpFit(); break; default: DoADCAverage(); } } void QmitkPreprocessingView::DoDwiNormalization() { if (m_DiffusionImage.IsNull()) return; typedef mitk::DiffusionImage DiffusionImageType; typedef itk::DwiNormilzationFilter FilterType; FilterType::Pointer filter = FilterType::New(); filter->SetInput(m_DiffusionImage->GetVectorImage()); filter->SetGradientDirections(m_DiffusionImage->GetDirections()); if (m_Controls->m_NormalizationMethodBox->currentIndex()==1) filter->SetUseGlobalMax(true); filter->Update(); DiffusionImageType::Pointer image = DiffusionImageType::New(); image->SetVectorImage( filter->GetOutput() ); image->SetReferenceBValue( m_DiffusionImage->GetReferenceBValue() ); image->SetDirections( m_DiffusionImage->GetDirections() ); image->InitializeFromVectorImage(); mitk::DataNode::Pointer imageNode = mitk::DataNode::New(); imageNode->SetData( image ); QString name = m_SelectedDiffusionNodes.back()->GetName().c_str(); imageNode->SetName((name+"_normalized").toStdString().c_str()); GetDefaultDataStorage()->Add(imageNode, m_SelectedDiffusionNodes.back()); } void QmitkPreprocessingView::DoLengthCorrection() { if (m_DiffusionImage.IsNull()) return; typedef mitk::DiffusionImage DiffusionImageType; typedef itk::DwiGradientLengthCorrectionFilter FilterType; FilterType::Pointer filter = FilterType::New(); filter->SetRoundingValue( m_Controls->m_B_ValueMap_Rounder_SpinBox->value()); filter->SetReferenceBValue(m_DiffusionImage->GetReferenceBValue()); filter->SetReferenceGradientDirectionContainer(m_DiffusionImage->GetDirections()); filter->Update(); DiffusionImageType::Pointer image = DiffusionImageType::New(); image->SetVectorImage( m_DiffusionImage->GetVectorImage()); image->SetReferenceBValue( filter->GetNewBValue() ); image->SetDirections( filter->GetOutputGradientDirectionContainer()); image->InitializeFromVectorImage(); mitk::DataNode::Pointer imageNode = mitk::DataNode::New(); imageNode->SetData( image ); QString name = m_SelectedDiffusionNodes.back()->GetName().c_str(); imageNode->SetName((name+"_rounded").toStdString().c_str()); GetDefaultDataStorage()->Add(imageNode, m_SelectedDiffusionNodes.back()); } void QmitkPreprocessingView::UpdateDwiBValueMapRounder(int i) { if (m_DiffusionImage.IsNull()) return; //m_DiffusionImage->UpdateBValueMap(); UpdateBValueTableWidget(i); } void QmitkPreprocessingView::CallMultishellToSingleShellFilter(itk::DWIVoxelFunctor * functor, mitk::DiffusionImage::Pointer ImPtr, QString imageName, mitk::DataNode* parent) { typedef itk::RadialMultishellToSingleshellImageFilter FilterType; // filter input parameter const mitk::DiffusionImage::BValueMap &originalShellMap = ImPtr->GetBValueMap(); const mitk::DiffusionImage::ImageType *vectorImage = ImPtr->GetVectorImage(); const mitk::DiffusionImage::GradientDirectionContainerType::Pointer gradientContainer = ImPtr->GetDirections(); const unsigned int &bValue = ImPtr->GetReferenceBValue(); mitk::DataNode::Pointer imageNode = 0; // filter call FilterType::Pointer filter = FilterType::New(); filter->SetInput(vectorImage); filter->SetOriginalGradientDirections(gradientContainer); filter->SetOriginalBValueMap(originalShellMap); filter->SetOriginalBValue(bValue); filter->SetFunctor(functor); filter->Update(); // create new DWI image mitk::DiffusionImage::Pointer outImage = mitk::DiffusionImage::New(); outImage->SetVectorImage( filter->GetOutput() ); outImage->SetReferenceBValue( m_Controls->m_targetBValueSpinBox->value() ); outImage->SetDirections( filter->GetTargetGradientDirections() ); outImage->InitializeFromVectorImage(); imageNode = mitk::DataNode::New(); imageNode->SetData( outImage ); imageNode->SetName(imageName.toStdString().c_str()); GetDefaultDataStorage()->Add(imageNode, parent); // if(m_Controls->m_OutputRMSErrorImage->isChecked()){ // // create new Error image // FilterType::ErrorImageType::Pointer errImage = filter->GetErrorImage(); // mitk::Image::Pointer mitkErrImage = mitk::Image::New(); // mitkErrImage->InitializeByItk(errImage); // mitkErrImage->SetVolume(errImage->GetBufferPointer()); // imageNode = mitk::DataNode::New(); // imageNode->SetData( mitkErrImage ); // imageNode->SetName((imageName+"_Error").toStdString().c_str()); // GetDefaultDataStorage()->Add(imageNode); // } } void QmitkPreprocessingView::DoBiExpFit() { itk::BiExpFitFunctor::Pointer functor = itk::BiExpFitFunctor::New(); for (unsigned int i=0; i::Pointer inImage = dynamic_cast< mitk::DiffusionImage* >(m_SelectedDiffusionNodes.at(i)->GetData()); QString name(m_SelectedDiffusionNodes.at(i)->GetName().c_str()); const mitk::DiffusionImage::BValueMap & originalShellMap = inImage->GetBValueMap(); mitk::DiffusionImage::BValueMap::const_iterator it = originalShellMap.begin(); ++it;/* skip b=0*/ unsigned int s = 0; /*shell index */ vnl_vector bValueList(originalShellMap.size()-1); while(it != originalShellMap.end()) bValueList.put(s++,(it++)->first); const double targetBValue = m_Controls->m_targetBValueSpinBox->value(); functor->setListOfBValues(bValueList); functor->setTargetBValue(targetBValue); CallMultishellToSingleShellFilter(functor,inImage,name + "_BiExp", m_SelectedDiffusionNodes.at(i)); } } void QmitkPreprocessingView::DoAKCFit() { itk::KurtosisFitFunctor::Pointer functor = itk::KurtosisFitFunctor::New(); for (unsigned int i=0; i::Pointer inImage = dynamic_cast< mitk::DiffusionImage* >(m_SelectedDiffusionNodes.at(i)->GetData()); QString name(m_SelectedDiffusionNodes.at(i)->GetName().c_str()); const mitk::DiffusionImage::BValueMap & originalShellMap = inImage->GetBValueMap(); mitk::DiffusionImage::BValueMap::const_iterator it = originalShellMap.begin(); ++it;/* skip b=0*/ unsigned int s = 0; /*shell index */ vnl_vector bValueList(originalShellMap.size()-1); while(it != originalShellMap.end()) bValueList.put(s++,(it++)->first); const double targetBValue = m_Controls->m_targetBValueSpinBox->value(); functor->setListOfBValues(bValueList); functor->setTargetBValue(targetBValue); CallMultishellToSingleShellFilter(functor,inImage,name + "_AKC", m_SelectedDiffusionNodes.at(i)); } } void QmitkPreprocessingView::DoADCFit() { // later } void QmitkPreprocessingView::DoADCAverage() { itk::ADCAverageFunctor::Pointer functor = itk::ADCAverageFunctor::New(); for (unsigned int i=0; i::Pointer inImage = dynamic_cast< mitk::DiffusionImage* >(m_SelectedDiffusionNodes.at(i)->GetData()); QString name(m_SelectedDiffusionNodes.at(i)->GetName().c_str()); const mitk::DiffusionImage::BValueMap & originalShellMap = inImage->GetBValueMap(); mitk::DiffusionImage::BValueMap::const_iterator it = originalShellMap.begin(); ++it;/* skip b=0*/ unsigned int s = 0; /*shell index */ vnl_vector bValueList(originalShellMap.size()-1); while(it != originalShellMap.end()) bValueList.put(s++,(it++)->first); const double targetBValue = m_Controls->m_targetBValueSpinBox->value(); functor->setListOfBValues(bValueList); functor->setTargetBValue(targetBValue); CallMultishellToSingleShellFilter(functor,inImage,name + "_ADC", m_SelectedDiffusionNodes.at(i)); } } void QmitkPreprocessingView::DoAdcCalculation() { if (m_DiffusionImage.IsNull()) return; typedef mitk::DiffusionImage< DiffusionPixelType > DiffusionImageType; typedef itk::AdcImageFilter< DiffusionPixelType, double > FilterType; for (unsigned int i=0; i(m_SelectedDiffusionNodes.at(i)->GetData()); FilterType::Pointer filter = FilterType::New(); filter->SetInput(inImage->GetVectorImage()); filter->SetGradientDirections(inImage->GetDirections()); filter->SetB_value(inImage->GetReferenceBValue()); filter->Update(); mitk::Image::Pointer image = mitk::Image::New(); image->InitializeByItk( filter->GetOutput() ); image->SetVolume( filter->GetOutput()->GetBufferPointer() ); mitk::DataNode::Pointer imageNode = mitk::DataNode::New(); imageNode->SetData( image ); QString name = m_SelectedDiffusionNodes.at(i)->GetName().c_str(); imageNode->SetName((name+"_ADC").toStdString().c_str()); GetDefaultDataStorage()->Add(imageNode, m_SelectedDiffusionNodes.at(i)); } } void QmitkPreprocessingView::UpdateBValueTableWidget(int i) { if (m_DiffusionImage.IsNull()) { m_Controls->m_B_ValueMap_TableWidget->clear(); m_Controls->m_B_ValueMap_TableWidget->setRowCount(1); QStringList headerList; headerList << "b-Value" << "Number of gradients"; m_Controls->m_B_ValueMap_TableWidget->setHorizontalHeaderLabels(headerList); m_Controls->m_B_ValueMap_TableWidget->setItem(0,0,new QTableWidgetItem("-")); m_Controls->m_B_ValueMap_TableWidget->setItem(0,1,new QTableWidgetItem("-")); }else{ typedef mitk::DiffusionImage::BValueMap BValueMap; typedef mitk::DiffusionImage::BValueMap::iterator BValueMapIterator; BValueMapIterator it; BValueMap roundedBValueMap = m_DiffusionImage->GetBValueMap(); m_Controls->m_B_ValueMap_TableWidget->clear(); m_Controls->m_B_ValueMap_TableWidget->setRowCount(roundedBValueMap.size() ); QStringList headerList; headerList << "b-Value" << "Number of gradients"; m_Controls->m_B_ValueMap_TableWidget->setHorizontalHeaderLabels(headerList); int i = 0 ; for(it = roundedBValueMap.begin() ;it != roundedBValueMap.end(); it++) { m_Controls->m_B_ValueMap_TableWidget->setItem(i,0,new QTableWidgetItem(QString::number(it->first))); QTableWidgetItem* item = m_Controls->m_B_ValueMap_TableWidget->item(i,0); item->setFlags(item->flags() & ~Qt::ItemIsEditable); m_Controls->m_B_ValueMap_TableWidget->setItem(i,1,new QTableWidgetItem(QString::number(it->second.size()))); i++; } } } template < typename TPixel, unsigned int VImageDimension > void QmitkPreprocessingView::TemplatedUpdateGui( itk::Image* itkImage) { for (int r=0; r<3; r++) for (int c=0; c<3; c++) { QTableWidgetItem* item = m_Controls->m_DirectionMatrixTable->item(r,c); delete item; item = new QTableWidgetItem(); item->setTextAlignment(Qt::AlignCenter | Qt::AlignVCenter); item->setText(QString::number(itkImage->GetDirection()[r][c])); m_Controls->m_DirectionMatrixTable->setItem(r,c,item); } } void QmitkPreprocessingView::OnSelectionChanged( std::vector nodes ) { bool foundDwiVolume = false; bool foundImageVolume = false; m_DiffusionImage = NULL; m_SelectedImage = NULL; m_SelectedImageNode = NULL; m_SelectedDiffusionNodes.clear(); // iterate selection for( std::vector::iterator it = nodes.begin(); it != nodes.end(); ++it ) { mitk::DataNode::Pointer node = *it; if( node.IsNotNull() && dynamic_cast*>(node->GetData()) ) { foundDwiVolume = true; foundImageVolume = true; m_DiffusionImage = dynamic_cast*>(node->GetData()); m_Controls->m_DiffusionImageLabel->setText(node->GetName().c_str()); m_SelectedDiffusionNodes.push_back(node); } else if( node.IsNotNull() && dynamic_cast(node->GetData()) ) { foundImageVolume = true; m_SelectedImage = dynamic_cast(node->GetData()); m_Controls->m_DiffusionImageLabel->setText(node->GetName().c_str()); m_SelectedImageNode = node; } } m_Controls->m_ButtonAverageGradients->setEnabled(foundDwiVolume); m_Controls->m_ButtonExtractB0->setEnabled(foundDwiVolume); m_Controls->m_CheckExtractAll->setEnabled(foundDwiVolume); m_Controls->m_ModifyMeasurementFrame->setEnabled(foundDwiVolume); m_Controls->m_MeasurementFrameTable->setEnabled(foundDwiVolume); m_Controls->m_ReduceGradientsButton->setEnabled(foundDwiVolume); m_Controls->m_ShowGradientsButton->setEnabled(foundDwiVolume); m_Controls->m_MirrorGradientToHalfSphereButton->setEnabled(foundDwiVolume); m_Controls->m_MergeDwisButton->setEnabled(foundDwiVolume); m_Controls->m_B_ValueMap_Rounder_SpinBox->setEnabled(foundDwiVolume); m_Controls->m_ProjectSignalButton->setEnabled(foundDwiVolume); m_Controls->m_CreateLengthCorrectedDwi->setEnabled(foundDwiVolume); m_Controls->m_CalcAdcButton->setEnabled(foundDwiVolume); m_Controls->m_targetBValueSpinBox->setEnabled(foundDwiVolume); m_Controls->m_NormalizeImageValuesButton->setEnabled(foundDwiVolume); m_Controls->m_DirectionMatrixTable->setEnabled(foundImageVolume); m_Controls->m_ModifyDirection->setEnabled(foundImageVolume); m_Controls->m_ExtractBrainMask->setEnabled(foundImageVolume); m_Controls->m_ResampleImageButton->setEnabled(foundImageVolume); // reset sampling frame to 1 and update all ealted components m_Controls->m_B_ValueMap_Rounder_SpinBox->setValue(1); UpdateBValueTableWidget(m_Controls->m_B_ValueMap_Rounder_SpinBox->value()); DoUpdateInterpolationGui(m_Controls->m_ResampleTypeBox->currentIndex()); if (foundDwiVolume) { m_Controls->m_InputData->setTitle("Input Data"); vnl_matrix_fixed< double, 3, 3 > mf = m_DiffusionImage->GetMeasurementFrame(); for (int r=0; r<3; r++) for (int c=0; c<3; c++) { // Measurement frame { QTableWidgetItem* item = m_Controls->m_MeasurementFrameTable->item(r,c); delete item; item = new QTableWidgetItem(); item->setTextAlignment(Qt::AlignCenter | Qt::AlignVCenter); item->setText(QString::number(mf.get(r,c))); m_Controls->m_MeasurementFrameTable->setItem(r,c,item); } // Direction matrix { QTableWidgetItem* item = m_Controls->m_DirectionMatrixTable->item(r,c); delete item; item = new QTableWidgetItem(); item->setTextAlignment(Qt::AlignCenter | Qt::AlignVCenter); item->setText(QString::number(m_DiffusionImage->GetVectorImage()->GetDirection()[r][c])); m_Controls->m_DirectionMatrixTable->setItem(r,c,item); } } //calculate target bValue for MultishellToSingleShellfilter const mitk::DiffusionImage::BValueMap & bValMap = m_DiffusionImage->GetBValueMap(); mitk::DiffusionImage::BValueMap::const_iterator it = bValMap.begin(); unsigned int targetBVal = 0; while(it != bValMap.end()) targetBVal += (it++)->first; targetBVal /= (float)bValMap.size()-1; m_Controls->m_targetBValueSpinBox->setValue(targetBVal); } else if (foundImageVolume) { for (int r=0; r<3; r++) for (int c=0; c<3; c++) { QTableWidgetItem* item = m_Controls->m_MeasurementFrameTable->item(r,c); delete item; item = new QTableWidgetItem(); m_Controls->m_MeasurementFrameTable->setItem(r,c,item); } AccessFixedDimensionByItk(m_SelectedImage, TemplatedUpdateGui,3); } else { for (int r=0; r<3; r++) for (int c=0; c<3; c++) { { QTableWidgetItem* item = m_Controls->m_MeasurementFrameTable->item(r,c); delete item; item = new QTableWidgetItem(); m_Controls->m_MeasurementFrameTable->setItem(r,c,item); } { QTableWidgetItem* item = m_Controls->m_DirectionMatrixTable->item(r,c); delete item; item = new QTableWidgetItem(); m_Controls->m_DirectionMatrixTable->setItem(r,c,item); } } m_Controls->m_DiffusionImageLabel->setText("mandatory"); m_Controls->m_InputData->setTitle("Please Select Input Data"); } } void QmitkPreprocessingView::Activated() { QmitkFunctionality::Activated(); } void QmitkPreprocessingView::Deactivated() { QmitkFunctionality::Deactivated(); } void QmitkPreprocessingView::DoHalfSphereGradientDirections() { - MitkDwiType::Pointer newDwi = m_DiffusionImage->GetDeepCopy(); + + MitkDwiType::Pointer newDwi = m_DiffusionImage->Clone(); GradientDirectionContainerType::Pointer gradientContainer = newDwi->GetDirections(); for (unsigned int j=0; jSize(); j++) if (gradientContainer->at(j)[0]<0) gradientContainer->at(j) = -gradientContainer->at(j); newDwi->SetDirections(gradientContainer); mitk::DataNode::Pointer imageNode = mitk::DataNode::New(); imageNode->SetData( newDwi ); QString name = m_SelectedDiffusionNodes.back()->GetName().c_str(); imageNode->SetName((name+"_halfsphere").toStdString().c_str()); GetDefaultDataStorage()->Add(imageNode, m_SelectedDiffusionNodes.back()); } void QmitkPreprocessingView::DoApplyMesurementFrame() { if (m_DiffusionImage.IsNull()) return; vnl_matrix_fixed< double, 3, 3 > mf; for (int r=0; r<3; r++) for (int c=0; c<3; c++) { QTableWidgetItem* item = m_Controls->m_MeasurementFrameTable->item(r,c); if (!item) return; mf[r][c] = item->text().toDouble(); } - MitkDwiType::Pointer newDwi = m_DiffusionImage->GetDeepCopy(); + MitkDwiType::Pointer newDwi = m_DiffusionImage->Clone(); newDwi->SetMeasurementFrame(mf); mitk::DataNode::Pointer imageNode = mitk::DataNode::New(); imageNode->SetData( newDwi ); QString name = m_SelectedDiffusionNodes.back()->GetName().c_str(); imageNode->SetName((name+"_new-MF").toStdString().c_str()); GetDefaultDataStorage()->Add(imageNode, m_SelectedDiffusionNodes.back()); } void QmitkPreprocessingView::DoShowGradientDirections() { if (m_DiffusionImage.IsNull()) return; int maxIndex = 0; unsigned int maxSize = m_DiffusionImage->GetDimension(0); if (maxSizeGetDimension(1)) { maxSize = m_DiffusionImage->GetDimension(1); maxIndex = 1; } if (maxSizeGetDimension(2)) { maxSize = m_DiffusionImage->GetDimension(2); maxIndex = 2; } mitk::Point3D origin = m_DiffusionImage->GetGeometry()->GetOrigin(); mitk::PointSet::Pointer originSet = mitk::PointSet::New(); typedef mitk::DiffusionImage::BValueMap BValueMap; typedef mitk::DiffusionImage::BValueMap::iterator BValueMapIterator; BValueMap bValMap = m_DiffusionImage->GetBValueMap(); GradientDirectionContainerType::Pointer gradientContainer = m_DiffusionImage->GetDirections(); mitk::BaseGeometry::Pointer geometry = m_DiffusionImage->GetGeometry(); int shellCount = 1; for(BValueMapIterator it = bValMap.begin(); it!=bValMap.end(); ++it) { mitk::PointSet::Pointer pointset = mitk::PointSet::New(); for (unsigned int j=0; jsecond.size(); j++) { mitk::Point3D ip; vnl_vector_fixed< double, 3 > v = gradientContainer->at(it->second[j]); if (v.magnitude()>mitk::eps) { ip[0] = v[0]*maxSize*geometry->GetSpacing()[maxIndex]/2 + origin[0]-0.5*geometry->GetSpacing()[0] + geometry->GetSpacing()[0]*m_DiffusionImage->GetDimension(0)/2; ip[1] = v[1]*maxSize*geometry->GetSpacing()[maxIndex]/2 + origin[1]-0.5*geometry->GetSpacing()[1] + geometry->GetSpacing()[1]*m_DiffusionImage->GetDimension(1)/2; ip[2] = v[2]*maxSize*geometry->GetSpacing()[maxIndex]/2 + origin[2]-0.5*geometry->GetSpacing()[2] + geometry->GetSpacing()[2]*m_DiffusionImage->GetDimension(2)/2; pointset->InsertPoint(j, ip); } else if (originSet->IsEmpty()) { ip[0] = v[0]*maxSize*geometry->GetSpacing()[maxIndex]/2 + origin[0]-0.5*geometry->GetSpacing()[0] + geometry->GetSpacing()[0]*m_DiffusionImage->GetDimension(0)/2; ip[1] = v[1]*maxSize*geometry->GetSpacing()[maxIndex]/2 + origin[1]-0.5*geometry->GetSpacing()[1] + geometry->GetSpacing()[1]*m_DiffusionImage->GetDimension(1)/2; ip[2] = v[2]*maxSize*geometry->GetSpacing()[maxIndex]/2 + origin[2]-0.5*geometry->GetSpacing()[2] + geometry->GetSpacing()[2]*m_DiffusionImage->GetDimension(2)/2; originSet->InsertPoint(j, ip); } } if (it->firstSetData(pointset); QString name = m_SelectedDiffusionNodes.front()->GetName().c_str(); name += "_Shell_"; name += QString::number(it->first); node->SetName(name.toStdString().c_str()); node->SetProperty("pointsize", mitk::FloatProperty::New((float)maxSize/50)); int b0 = shellCount%2; int b1 = 0; int b2 = 0; if (shellCount>4) b2 = 1; if (shellCount%4 >= 2) b1 = 1; node->SetProperty("color", mitk::ColorProperty::New(b2, b1, b0)); GetDefaultDataStorage()->Add(node, m_SelectedDiffusionNodes.front()); shellCount++; } // add origin to datastorage mitk::DataNode::Pointer node = mitk::DataNode::New(); node->SetData(originSet); QString name = m_SelectedDiffusionNodes.front()->GetName().c_str(); name += "_Origin"; node->SetName(name.toStdString().c_str()); node->SetProperty("pointsize", mitk::FloatProperty::New((float)maxSize/50)); node->SetProperty("color", mitk::ColorProperty::New(1,1,1)); GetDefaultDataStorage()->Add(node, m_SelectedDiffusionNodes.front()); } void QmitkPreprocessingView::DoReduceGradientDirections() { if (m_DiffusionImage.IsNull()) return; typedef mitk::DiffusionImage DiffusionImageType; typedef itk::ElectrostaticRepulsionDiffusionGradientReductionFilter FilterType; typedef DiffusionImageType::BValueMap BValueMap; // GetShellSelection from GUI BValueMap shellSlectionMap; BValueMap originalShellMap = m_DiffusionImage->GetBValueMap(); std::vector newNumGradientDirections; int shellCounter = 0; QString name = m_SelectedDiffusionNodes.front()->GetName().c_str(); for (int i=0; im_B_ValueMap_TableWidget->rowCount(); i++) { double BValue = m_Controls->m_B_ValueMap_TableWidget->item(i,0)->text().toDouble(); shellSlectionMap[BValue] = originalShellMap[BValue]; unsigned int num = m_Controls->m_B_ValueMap_TableWidget->item(i,1)->text().toUInt(); newNumGradientDirections.push_back(num); name += "_"; name += QString::number(num); shellCounter++; } if (newNumGradientDirections.empty()) return; GradientDirectionContainerType::Pointer gradientContainer = m_DiffusionImage->GetDirections(); FilterType::Pointer filter = FilterType::New(); filter->SetInput(m_DiffusionImage->GetVectorImage()); filter->SetOriginalGradientDirections(gradientContainer); filter->SetNumGradientDirections(newNumGradientDirections); filter->SetOriginalBValueMap(originalShellMap); filter->SetShellSelectionBValueMap(shellSlectionMap); filter->Update(); DiffusionImageType::Pointer image = DiffusionImageType::New(); image->SetVectorImage( filter->GetOutput() ); image->SetReferenceBValue(m_DiffusionImage->GetReferenceBValue()); image->SetDirections(filter->GetGradientDirections()); image->SetMeasurementFrame(m_DiffusionImage->GetMeasurementFrame()); image->InitializeFromVectorImage(); mitk::DataNode::Pointer imageNode = mitk::DataNode::New(); imageNode->SetData( image ); imageNode->SetName(name.toStdString().c_str()); GetDefaultDataStorage()->Add(imageNode, m_SelectedDiffusionNodes.back()); } void QmitkPreprocessingView::MergeDwis() { typedef mitk::DiffusionImage DiffusionImageType; typedef DiffusionImageType::GradientDirectionContainerType GradientContainerType; if (m_SelectedDiffusionNodes.size()<2) return; typedef itk::VectorImage DwiImageType; typedef DwiImageType::PixelType DwiPixelType; typedef DwiImageType::RegionType DwiRegionType; typedef std::vector< DwiImageType::Pointer > DwiImageContainerType; typedef std::vector< GradientContainerType::Pointer > GradientListContainerType; DwiImageContainerType imageContainer; GradientListContainerType gradientListContainer; std::vector< double > bValueContainer; QString name = m_SelectedDiffusionNodes.front()->GetName().c_str(); for (unsigned int i=0; i* >( m_SelectedDiffusionNodes.at(i)->GetData() ); if ( dwi.IsNotNull() ) { imageContainer.push_back(dwi->GetVectorImage()); gradientListContainer.push_back(dwi->GetDirections()); bValueContainer.push_back(dwi->GetReferenceBValue()); if (i>0) { name += "+"; name += m_SelectedDiffusionNodes.at(i)->GetName().c_str(); } } } typedef itk::MergeDiffusionImagesFilter FilterType; FilterType::Pointer filter = FilterType::New(); filter->SetImageVolumes(imageContainer); filter->SetGradientLists(gradientListContainer); filter->SetBValues(bValueContainer); filter->Update(); vnl_matrix_fixed< double, 3, 3 > mf; mf.set_identity(); DiffusionImageType::Pointer image = DiffusionImageType::New(); image->SetVectorImage( filter->GetOutput() ); image->SetReferenceBValue(filter->GetB_Value()); image->SetDirections(filter->GetOutputGradients()); image->SetMeasurementFrame(mf); image->InitializeFromVectorImage(); mitk::DataNode::Pointer imageNode = mitk::DataNode::New(); imageNode->SetData( image ); imageNode->SetName(name.toStdString().c_str()); GetDefaultDataStorage()->Add(imageNode); } void QmitkPreprocessingView::ExtractB0() { typedef mitk::DiffusionImage DiffusionImageType; typedef DiffusionImageType::GradientDirectionContainerType GradientContainerType; int nrFiles = m_SelectedDiffusionNodes.size(); if (!nrFiles) return; // call the extraction withou averaging if the check-box is checked if( this->m_Controls->m_CheckExtractAll->isChecked() ) { DoExtractBOWithoutAveraging(); return; } mitk::DataStorage::SetOfObjects::const_iterator itemiter( m_SelectedDiffusionNodes.begin() ); mitk::DataStorage::SetOfObjects::const_iterator itemiterend( m_SelectedDiffusionNodes.end() ); std::vector nodes; while ( itemiter != itemiterend ) // for all items { DiffusionImageType* vols = static_cast( (*itemiter)->GetData()); std::string nodename; (*itemiter)->GetStringProperty("name", nodename); // Extract image using found index typedef itk::B0ImageExtractionImageFilter FilterType; FilterType::Pointer filter = FilterType::New(); filter->SetInput(vols->GetVectorImage()); filter->SetDirections(vols->GetDirections()); filter->Update(); mitk::Image::Pointer mitkImage = mitk::Image::New(); mitkImage->InitializeByItk( filter->GetOutput() ); mitkImage->SetVolume( filter->GetOutput()->GetBufferPointer() ); mitk::DataNode::Pointer node=mitk::DataNode::New(); node->SetData( mitkImage ); node->SetProperty( "name", mitk::StringProperty::New(nodename + "_B0")); GetDefaultDataStorage()->Add(node, (*itemiter)); ++itemiter; } } void QmitkPreprocessingView::DoExtractBOWithoutAveraging() { // typedefs typedef mitk::DiffusionImage DiffusionImageType; typedef DiffusionImageType::GradientDirectionContainerType GradientContainerType; typedef itk::B0ImageExtractionToSeparateImageFilter< short, short> FilterType; // check number of selected objects, return if empty int nrFiles = m_SelectedDiffusionNodes.size(); if (!nrFiles) return; mitk::DataStorage::SetOfObjects::const_iterator itemiter( m_SelectedDiffusionNodes.begin() ); mitk::DataStorage::SetOfObjects::const_iterator itemiterend( m_SelectedDiffusionNodes.end() ); while ( itemiter != itemiterend ) // for all items { DiffusionImageType* vols = static_cast( (*itemiter)->GetData()); std::string nodename; (*itemiter)->GetStringProperty("name", nodename); // Extract image using found index FilterType::Pointer filter = FilterType::New(); filter->SetInput(vols->GetVectorImage()); filter->SetDirections(vols->GetDirections()); filter->Update(); mitk::Image::Pointer mitkImage = mitk::Image::New(); mitkImage->InitializeByItk( filter->GetOutput() ); mitkImage->SetImportChannel( filter->GetOutput()->GetBufferPointer() ); mitk::DataNode::Pointer node=mitk::DataNode::New(); node->SetData( mitkImage ); node->SetProperty( "name", mitk::StringProperty::New(nodename + "_B0_ALL")); GetDefaultDataStorage()->Add(node, (*itemiter)); /*A reinitialization is needed to access the time channels via the ImageNavigationController The Global-Geometry can not recognize the time channel without a re-init. (for a new selection in datamanger a automatically updated of the Global-Geometry should be done - if it contains the time channel)*/ mitk::RenderingManager::GetInstance()->InitializeViews(node->GetData()->GetTimeGeometry(),mitk::RenderingManager::REQUEST_UPDATE_ALL, true); ++itemiter; } } void QmitkPreprocessingView::AverageGradients() { int nrFiles = m_SelectedDiffusionNodes.size(); if (!nrFiles) return; mitk::DataStorage::SetOfObjects::const_iterator itemiter( m_SelectedDiffusionNodes.begin() ); mitk::DataStorage::SetOfObjects::const_iterator itemiterend( m_SelectedDiffusionNodes.end() ); while ( itemiter != itemiterend ) // for all items { mitk::DiffusionImage* mitkDwi = static_cast*>( (*itemiter)->GetData()); - MitkDwiType::Pointer newDwi = mitkDwi->GetDeepCopy(); + MitkDwiType::Pointer newDwi = mitkDwi->Clone(); newDwi->AverageRedundantGradients(m_Controls->m_Blur->value()); mitk::DataNode::Pointer imageNode = mitk::DataNode::New(); imageNode->SetData( newDwi ); QString name = (*itemiter)->GetName().c_str(); imageNode->SetName((name+"_averaged").toStdString().c_str()); GetDefaultDataStorage()->Add(imageNode, (*itemiter)); ++itemiter; } }