diff --git a/Modules/DiffusionImaging/Reconstruction/itkDiffusionMultiShellQballReconstructionImageFilter.cpp b/Modules/DiffusionImaging/Reconstruction/itkDiffusionMultiShellQballReconstructionImageFilter.cpp index 4b8d2f2ca3..72064fc5d3 100644 --- a/Modules/DiffusionImaging/Reconstruction/itkDiffusionMultiShellQballReconstructionImageFilter.cpp +++ b/Modules/DiffusionImaging/Reconstruction/itkDiffusionMultiShellQballReconstructionImageFilter.cpp @@ -1,583 +1,785 @@ /*========================================================================= Program: Medical Imaging & Interaction Toolkit Language: C++ Date: $Date: 2009-07-14 19:11:20 +0200 (Tue, 14 Jul 2009) $ Version: $Revision: 18127 $ Copyright (c) German Cancer Research Center, Division of Medical and Biological Informatics. All rights reserved. See MITKCopyright.txt or http://www.mitk.org/copyright.html for details. This software is distributed WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the above copyright notices for more information. =========================================================================*/ #ifndef __itkDiffusionMultiShellQballReconstructionImageFilter_cpp #define __itkDiffusionMultiShellQballReconstructionImageFilter_cpp #include #include #include #include #include #include #include #include #include #include #define _USE_MATH_DEFINES #include #include "mitkSphericalHarmonicsFunctions.h" #include "itkPointShell.h" #include namespace itk { #define QBALL_ANAL_RECON_PI M_PI template< class T, class TG, class TO, int L, int NODF> DiffusionMultiShellQballReconstructionImageFilter ::DiffusionMultiShellQballReconstructionImageFilter() : m_GradientDirectionContainer(NULL), m_GradientIndexMap(NULL), m_NumberOfGradientDirections(0), m_NumberOfBaselineImages(1), m_Threshold(NumericTraits< ReferencePixelType >::NonpositiveMin()), m_BValue(1.0), m_Lambda(0.0), - m_IsHemisphericalArrangementOfGradientDirections(false) + m_IsHemisphericalArrangementOfGradientDirections(false), + m_IsArithmeticProgession(false), + m_ReconstructionType(Standard1Shell) { // At least 1 inputs is necessary for a vector image. // For images added one at a time we need at least six this->SetNumberOfRequiredInputs( 1 ); } template< class TReferenceImagePixelType, class TGradientImagePixelType, class TOdfPixelType, int NOrderL, int NrOdfDirections> typename itk::DiffusionMultiShellQballReconstructionImageFilter< TReferenceImagePixelType,TGradientImagePixelType,TOdfPixelType, NOrderL,NrOdfDirections>::OdfPixelType itk::DiffusionMultiShellQballReconstructionImageFilter ::Normalize( OdfPixelType odf, typename NumericTraits::AccumulateType b0 ) { for(int i=0; i0) odf /= sum; return odf; } template vnl_vectoritk::DiffusionMultiShellQballReconstructionImageFilter ::PreNormalize( vnl_vector vec, typename NumericTraits::AccumulateType b0 ) { // Threshold double sigma = 0.001; double temp; - int n = vec.size(); + double b0f = (double)b0; - for(int i=0; i " << vec[i]; vec[i] = log(-log(vec[i])); } return vec; } template< class T, class TG, class TO, int L, int NODF> void DiffusionMultiShellQballReconstructionImageFilter ::SetGradientImage( GradientDirectionContainerType *gradientDirection, const GradientImagesType *gradientImage ) { this->m_GradientDirectionContainer = gradientDirection; this->m_NumberOfBaselineImages = 0; for(GradientDirectionContainerType::Iterator it = gradientDirection->Begin(); it != gradientDirection->End(); it++) { if( it.Value().one_norm() <= 0.0 ) this->m_NumberOfBaselineImages ++; else it.Value() = it.Value() / it.Value().two_norm(); } this->m_NumberOfGradientDirections = gradientDirection->Size() - this->m_NumberOfBaselineImages; // ensure that the gradient image we received has as many components as // the number of gradient directions if( gradientImage->GetVectorLength() != this->m_NumberOfBaselineImages + m_NumberOfGradientDirections ) { itkExceptionMacro( << m_NumberOfGradientDirections << " gradients + " << this->m_NumberOfBaselineImages << "baselines = " << m_NumberOfGradientDirections + this->m_NumberOfBaselineImages << " directions specified but image has " << gradientImage->GetVectorLength() << " components."); } this->ProcessObject::SetNthInput( 0, const_cast< GradientImagesType* >(gradientImage) ); } template< class T, class TG, class TO, int L, int NODF> void DiffusionMultiShellQballReconstructionImageFilter ::BeforeThreadedGenerateData() { if( m_NumberOfGradientDirections < (((L+1)*(L+2))/2) /* && m_NumberOfGradientDirections < 6 */ ) { itkExceptionMacro( << "At least " << ((L+1)*(L+2))/2 << " gradient directions are required" ); } // Input must be an itk::VectorImage. std::string gradientImageClassName(this->ProcessObject::GetInput(0)->GetNameOfClass()); if ( strcmp(gradientImageClassName.c_str(),"VectorImage") != 0 ) itkExceptionMacro( << "There is only one Gradient image. I expect that to be a VectorImage. But its of type: " << gradientImageClassName ); - this->ComputeReconstructionMatrix(); - m_BZeroImage = BZeroImageType::New(); typename GradientImagesType::Pointer img = static_cast< GradientImagesType * >( this->ProcessObject::GetInput(0) ); m_BZeroImage->SetSpacing( img->GetSpacing() ); // Set the image spacing m_BZeroImage->SetOrigin( img->GetOrigin() ); // Set the image origin m_BZeroImage->SetDirection( img->GetDirection() ); // Set the image direction m_BZeroImage->SetLargestPossibleRegion( img->GetLargestPossibleRegion()); m_BZeroImage->SetBufferedRegion( img->GetLargestPossibleRegion() ); m_BZeroImage->Allocate(); -} - -template< class T, class TG, class TO, int L, int NODF> -void DiffusionMultiShellQballReconstructionImageFilter -::ThreadedGenerateData(const OutputImageRegionType& outputRegionForThread, int NumberOfThreads) -{ - // Get output image pointer - typename OutputImageType::Pointer outputImage = static_cast< OutputImageType * >(this->ProcessObject::GetOutput(0)); - - // ImageRegionIterator for the output image - ImageRegionIterator< OutputImageType > oit(outputImage, outputRegionForThread); - oit.GoToBegin(); - - // ImageRegionIterator for the BZero (output) image - ImageRegionIterator< BZeroImageType > oit2(m_BZeroImage, outputRegionForThread); - oit2.GoToBegin(); - - // if no GradientIndexMap is set take all Gradients to GradientIndicies-Vector and add this to GradientIndexMap[ALL] // add All Bzero if(m_GradientIndexMap == 0){ m_GradientIndexMap = new GradientIndexMap(); // split for all gradient directions the image in baseline-indicies and gradient-indicies for a fast access GradientDirectionContainerType::ConstIterator gdcit; for( gdcit = this->m_GradientDirectionContainer->Begin(); gdcit != this->m_GradientDirectionContainer->End(); ++gdcit) { if(gdcit.Value().one_norm() <= 0.0) // if vector length at current position is 0, it is a basline signal { - (*m_GradientIndexMap)[BZero].push_back(gdcit.Index()); + (*m_GradientIndexMap)[0].push_back(gdcit.Index()); }else{ // it is a gradient signal of length != 0 - (*m_GradientIndexMap)[AllGradients].push_back(gdcit.Index()); + (*m_GradientIndexMap)[1].push_back(gdcit.Index()); } } } + if(m_GradientIndexMap->size() == 4){ + + GradientIndexMapIteraotr it = m_GradientIndexMap->begin(); + it++; + int b1 = (*it).first; + it++; + int b2 = (*it).first; + it++; + int b3 = (*it).first; + + + if(b2 - b1 == b1 && b3 - b2 == b1 ) + { + m_ReconstructionType = Analytical3Shells; + }else + { + m_ReconstructionType = NumericalNShells; + } + }else if(m_GradientIndexMap->size() > 2) + { + m_ReconstructionType = NumericalNShells; + }else{ + m_ReconstructionType = Standard1Shell; + } + + + + switch(m_ReconstructionType) + { + case Analytical3Shells: + this->ComputeReconstructionMatrix((*(m_GradientIndexMap->end())).second); + break; + case NumericalNShells: + this->ComputeReconstructionMatrix((*m_GradientIndexMap)[1]); + break; + case Standard1Shell: + this->ComputeReconstructionMatrix((*m_GradientIndexMap)[1]); + break; + } + + +} + + + +template< class T, class TG, class TO, int L, int NODF> +void DiffusionMultiShellQballReconstructionImageFilter +::StandardOneShellReconstruction(const OutputImageRegionType& outputRegionForThread) +{ + // Get output image pointer + typename OutputImageType::Pointer outputImage = static_cast< OutputImageType * >(this->ProcessObject::GetOutput(0)); + + // ImageRegionIterator for the output image + ImageRegionIterator< OutputImageType > oit(outputImage, outputRegionForThread); + oit.GoToBegin(); + + // ImageRegionIterator for the BZero (output) image + ImageRegionIterator< BZeroImageType > oit2(m_BZeroImage, outputRegionForThread); + oit2.GoToBegin(); + + IndiciesVector BZeroIndicies = (*m_GradientIndexMap)[0]; + IndiciesVector SignalIndicies = (*m_GradientIndexMap)[1]; + // if the gradient directiosn aragement is hemispherical, duplicate all gradient directions // alone, interested in the value, the direction can be neglected if(m_IsHemisphericalArrangementOfGradientDirections){ - int NumbersOfGradientIndicies = (*m_GradientIndexMap)[AllGradients].size(); + int NumbersOfGradientIndicies = SignalIndicies.size(); for (int i = 0 ; i < NumbersOfGradientIndicies; i++) - (*m_GradientIndexMap)[AllGradients].push_back((*m_GradientIndexMap)[AllGradients][i]); + SignalIndicies.push_back(SignalIndicies[i]); } // Get input gradient image pointer typename GradientImagesType::Pointer gradientImagePointer = static_cast< GradientImagesType * >( this->ProcessObject::GetInput(0) ); // Const ImageRegionIterator for input gradient image typedef ImageRegionConstIterator< GradientImagesType > GradientIteratorType; GradientIteratorType git(gradientImagePointer, outputRegionForThread ); git.GoToBegin(); typedef typename GradientImagesType::PixelType GradientVectorType; // iterate overall voxels of the gradient image region while( ! git.IsAtEnd() ) { GradientVectorType b = git.Get(); // compute the average bzero signal typename NumericTraits::AccumulateType b0 = NumericTraits::Zero; - for(unsigned int i = 0; i < (*m_GradientIndexMap)[BZero].size(); ++i) + for(unsigned int i = 0; i < BZeroIndicies.size(); ++i) { - b0 += b[(*m_GradientIndexMap)[BZero][i]]; + b0 += b[BZeroIndicies[i]]; } b0 /= this->m_NumberOfBaselineImages; // ODF Vector OdfPixelType odf(0.0); // Create the Signal Vector vnl_vector SignalVector(m_NumberOfGradientDirections); if( (b0 != 0) && (b0 >= m_Threshold) ) { - IndiciesVector & gradientIndicies = (*m_GradientIndexMap)[AllGradients]; - for( unsigned int i = 0; i< m_NumberOfGradientDirections; i++ ) + + for( unsigned int i = 0; i< SignalIndicies.size(); i++ ) { - unsigned int index = gradientIndicies[i]; - SignalVector[i] = static_cast(b[index]); + SignalVector[i] = static_cast(b[SignalIndicies[i]]); } // apply threashold an generate ln(-ln(E)) signal // Replace SignalVector with PreNormalized SignalVector SignalVector = PreNormalize(SignalVector, b0); // ODF coeffs-vector vnl_vector coeffs(m_NumberCoefficients); // approximate ODF coeffs coeffs = ( (*m_CoeffReconstructionMatrix) * SignalVector ); coeffs[0] += 1.0/(2.0*sqrt(QBALL_ANAL_RECON_PI)); odf = ( (*m_SphericalHarmonicBasisMatrix) * coeffs ).data_block(); } // set ODF to ODF-Image oit.Set( odf ); ++oit; // set b0(average) to b0-Image oit2.Set( b0 ); ++oit2; ++git; } MITK_INFO << "One Thread finished reconstruction"; +} + +template< class T, class TG, class TO, int L, int NODF> +void DiffusionMultiShellQballReconstructionImageFilter +::AnalyticalThreeShellReconstruction(const OutputImageRegionType& outputRegionForThread) +{ + + typename OutputImageType::Pointer outputImage = static_cast< OutputImageType * >(this->ProcessObject::GetOutput(0)); + + ImageRegionIterator< OutputImageType > oit(outputImage, outputRegionForThread); + oit.GoToBegin(); + + ImageRegionIterator< BZeroImageType > oit2(m_BZeroImage, outputRegionForThread); + oit2.GoToBegin(); + + IndiciesVector BZeroIndicies = (*m_GradientIndexMap)[0]; + + GradientIndexMapIteraotr it = m_GradientIndexMap->begin(); + it++; + IndiciesVector Shell1Indiecies = (*it).second; + it++; + IndiciesVector Shell2Indiecies = (*it).second; + it++; + IndiciesVector Shell3Indiecies = (*it).second; + + assert(Shell1Indiecies.size() == Shell2Indiecies.size() && Shell1Indiecies.size() == Shell3Indiecies.size()); + + if(m_IsHemisphericalArrangementOfGradientDirections){ + int NumbersOfGradientIndicies = Shell1Indiecies.size(); + for (int i = 0 ; i < NumbersOfGradientIndicies; i++){ + Shell1Indiecies.push_back(Shell1Indiecies[i]); + Shell2Indiecies.push_back(Shell2Indiecies[i]); + Shell3Indiecies.push_back(Shell3Indiecies[i]); + } + } + + // Get input gradient image pointer + typename GradientImagesType::Pointer gradientImagePointer = static_cast< GradientImagesType * >( this->ProcessObject::GetInput(0) ); + + // Const ImageRegionIterator for input gradient image + typedef ImageRegionConstIterator< GradientImagesType > GradientIteratorType; + GradientIteratorType git(gradientImagePointer, outputRegionForThread ); + git.GoToBegin(); + + typedef typename GradientImagesType::PixelType GradientVectorType; + + // iterate overall voxels of the gradient image region + while( ! git.IsAtEnd() ) + { + GradientVectorType b = git.Get(); + + // compute the average bzero signal + typename NumericTraits::AccumulateType b0 = NumericTraits::Zero; + for(unsigned int i = 0; i < BZeroIndicies.size(); ++i) + { + b0 += b[BZeroIndicies[i]]; + } + b0 /= this->m_NumberOfBaselineImages; + + // ODF Vector + OdfPixelType odf(0.0); + + // Create the Signal Vector + //vnl_vector SignalVector(m_NumberOfGradientDirections); + if( (b0 != 0) && (b0 >= m_Threshold) ) + { + + + + vnl_vector AlphaValues(Shell1Indiecies.size()); + vnl_vector BetaValues(Shell1Indiecies.size()); + vnl_vector LAValues(Shell1Indiecies.size()); + + for( unsigned int i = 0; i< Shell1Indiecies.size(); i++ ) + { + double E1 = b[Shell1Indiecies[i]]; + double E2 = b[Shell1Indiecies[i]]; + double E3 = b[Shell1Indiecies[i]]; + // CONSTRAINS INSERT HERE + double A = (E3 -E1*E2) / (2*(E2-E1*E1)); + double B = sqrt( A * A - ((E1*E3 - E2 * E2) / (E2-E1*E1)) ); + + + AlphaValues[i] = A + B; + BetaValues[i] = A - B; + LAValues[i] = 0.5 + ((E1 - A)/(2*B)); + } + + + // for( unsigned int i = 0; i< AlphaValues.size(); i++ ) + // { + // AlphaValues[i] = static_cast(AlphaValues[i]); + // } + + // ln(-ln(Alphak)) + AlphaValues = PreNormalize(AlphaValues, b0); + + BetaValues = PreNormalize(BetaValues, b0); + + + + vnl_vector SignalVector(Shell1Indiecies.size()); + for( unsigned int i = 0; i< AlphaValues.size(); i++ ) + { + SignalVector = LAValues[i] * AlphaValues + LAValues[i] * BetaValues; + } + + // ODF coeffs-vector + vnl_vector coeffs(m_NumberCoefficients); + + // approximate ODF coeffs + coeffs = ( (*m_CoeffReconstructionMatrix) * SignalVector ); + coeffs[0] += 1.0/(2.0*sqrt(QBALL_ANAL_RECON_PI)); + + odf = ( (*m_SphericalHarmonicBasisMatrix) * coeffs ).data_block(); + } + // set ODF to ODF-Image + oit.Set( odf ); + ++oit; + + // set b0(average) to b0-Image + oit2.Set( b0 ); + ++oit2; + + ++git; + + } + + +} + + +template< class T, class TG, class TO, int L, int NODF> +vnl_vector DiffusionMultiShellQballReconstructionImageFilter +::AnalyticalThreeShellParameterEstimation(const IndiciesVector * shell1Indicies,const IndiciesVector * shell2Indicies ,const IndiciesVector * shell3Indicies, vnl_vector) +{ + + + +} + +template< class T, class TG, class TO, int L, int NODF> +void DiffusionMultiShellQballReconstructionImageFilter +::ThreadedGenerateData(const OutputImageRegionType& outputRegionForThread, int NumberOfThreads) +{ + switch(m_ReconstructionType) + { + case Standard1Shell: + StandardOneShellReconstruction(outputRegionForThread); + break; + case Analytical3Shells: + AnalyticalThreeShellReconstruction(outputRegionForThread); + break; + case NumericalNShells: + break; + } + + } template< class T, class TG, class TO, int L, int NODF> void DiffusionMultiShellQballReconstructionImageFilter:: ComputeSphericalHarmonicsBasis(vnl_matrix * QBallReference, vnl_matrix *SHBasisOutput, vnl_matrix* LaplaciaBaltramiOutput, vnl_vector* SHOrderAssociation, vnl_matrix* SHEigenvalues ) { for(unsigned int i=0; i< (*SHBasisOutput).rows(); i++) { for(int k = 0; k <= L; k += 2) { for(int m =- k; m <= k; m++) { int j = ( k * k + k + 2 ) / 2 + m - 1; // Compute SHBasisFunctions double phi = (*QBallReference)(0,i); double th = (*QBallReference)(1,i); (*SHBasisOutput)(i,j) = mitk::ShericalHarmonicsFunctions::Yj(m,k,th,phi); // Laplacian Baltrami Order Association if(LaplaciaBaltramiOutput) (*LaplaciaBaltramiOutput)(j,j) = k*k*(k + 1)*(k+1); // SHEigenvalues with order Accosiation kj if(SHEigenvalues) (*SHEigenvalues)(j,j) = -k* (k+1); // Order Association if(SHOrderAssociation) (*SHOrderAssociation)[j] = k; } } } } template< class T, class TG, class TO, int L, int NODF> void DiffusionMultiShellQballReconstructionImageFilter ::ComputeFunkRadonTransformationMatrix(vnl_vector* SHOrderAssociationReference, vnl_matrix* FRTMatrixOutput ) { for(int i=0; i bool DiffusionMultiShellQballReconstructionImageFilter ::CheckHemisphericalArrangementOfGradientDirections() { // handle acquisition schemes where only half of the spherical // shell is sampled by the gradient directions. In this case, // each gradient direction is duplicated in negative direction. vnl_vector centerMass(3); centerMass.fill(0.0); int count = 0; GradientDirectionContainerType::ConstIterator gdcit1; for( gdcit1 = this->m_GradientDirectionContainer->Begin(); gdcit1 != this->m_GradientDirectionContainer->End(); ++gdcit1) { if(gdcit1.Value().one_norm() > 0.0) { centerMass += gdcit1.Value(); count ++; } } centerMass /= count; if(centerMass.two_norm() > 0.1) { return false; } return true; } template< class T, class TG, class TO, int L, int NOdfDirections> void DiffusionMultiShellQballReconstructionImageFilter -::ComputeReconstructionMatrix() +::ComputeReconstructionMatrix(std::vector< unsigned int > gradientIndiciesVector) { typedef std::auto_ptr< vnl_matrix< double> > MatrixDoublePtr; typedef std::auto_ptr< vnl_vector< int > > VectorIntPtr; typedef std::auto_ptr< vnl_matrix_inverse< double > > InverseMatrixDoublePtr; - if( m_NumberOfGradientDirections < (((L+1)*(L+2))/2) || m_NumberOfGradientDirections < 6 ) + int numberOfGradientDirections = gradientIndiciesVector.size(); + + if( numberOfGradientDirections < (((L+1)*(L+2))/2) || numberOfGradientDirections < 6 ) { itkExceptionMacro( << "At least (L+1)(L+2)/2 gradient directions are required" ); } CheckDuplicateDiffusionGradients(); // check if gradient directions are arrangement as a hemisphere(true) or sphere(false) m_IsHemisphericalArrangementOfGradientDirections = CheckHemisphericalArrangementOfGradientDirections(); - if(m_IsHemisphericalArrangementOfGradientDirections) m_NumberOfGradientDirections *= 2; + if(m_IsHemisphericalArrangementOfGradientDirections) numberOfGradientDirections *= 2; - MatrixDoublePtr Q(new vnl_matrix(3, m_NumberOfGradientDirections)); + MatrixDoublePtr Q(new vnl_matrix(3, numberOfGradientDirections)); Q->fill(0.0); // Cartesian to spherical coordinates { - int i = 0; - GradientDirectionContainerType::ConstIterator gdcit; - for( gdcit = this->m_GradientDirectionContainer->Begin(); gdcit != this->m_GradientDirectionContainer->End(); ++gdcit) + int j = 0; + + for(int i = 0; i < gradientIndiciesVector.size(); i++) { - if(gdcit.Value().one_norm() > 0.0) - { - double x = gdcit.Value().get(0); - double y = gdcit.Value().get(1); - double z = gdcit.Value().get(2); - double cart[3]; - mitk::ShericalHarmonicsFunctions::Cart2Sph(x,y,z,cart); - (*Q)(0,i) = cart[0]; - (*Q)(1,i) = cart[1]; - (*Q)(2,i++) = cart[2]; - } + + double x = this->m_GradientDirectionContainer->ElementAt(gradientIndiciesVector[i]).get(0); + double y = this->m_GradientDirectionContainer->ElementAt(gradientIndiciesVector[i]).get(1); + double z = this->m_GradientDirectionContainer->ElementAt(gradientIndiciesVector[i]).get(2); + double cart[3]; + mitk::ShericalHarmonicsFunctions::Cart2Sph(x,y,z,cart); + (*Q)(0,j) = cart[0]; + (*Q)(1,j) = cart[1]; + (*Q)(2,j++) = cart[2]; + } if(m_IsHemisphericalArrangementOfGradientDirections) { - for( gdcit = this->m_GradientDirectionContainer->Begin(); gdcit != this->m_GradientDirectionContainer->End(); ++gdcit) + for(int i = 0; i < gradientIndiciesVector.size(); i++) { - if(gdcit.Value().one_norm() > 0.0) - { - double x = gdcit.Value().get(0); - double y = gdcit.Value().get(1); - double z = gdcit.Value().get(2); - double cart[3]; - mitk::ShericalHarmonicsFunctions::Cart2Sph(x,y,z,cart); - (*Q)(0,i) = cart[0]; - (*Q)(1,i) = cart[1]; - (*Q)(2,i++) = cart[2]; - } + double x = this->m_GradientDirectionContainer->ElementAt(gradientIndiciesVector[i]).get(0); + double y = this->m_GradientDirectionContainer->ElementAt(gradientIndiciesVector[i]).get(1); + double z = this->m_GradientDirectionContainer->ElementAt(gradientIndiciesVector[i]).get(2); + double cart[3]; + mitk::ShericalHarmonicsFunctions::Cart2Sph(x,y,z,cart); + (*Q)(0,j) = cart[0]; + (*Q)(1,j) = cart[1]; + (*Q)(2,j++) = cart[2]; } } } int LOrder = L; m_NumberCoefficients = (int)(LOrder*LOrder + LOrder + 2.0)/2.0 + LOrder; - MatrixDoublePtr SHBasis(new vnl_matrix(m_NumberOfGradientDirections,m_NumberCoefficients)); + MatrixDoublePtr SHBasis(new vnl_matrix(numberOfGradientDirections,m_NumberCoefficients)); SHBasis->fill(0.0); VectorIntPtr SHOrderAssociation(new vnl_vector(m_NumberCoefficients)); SHOrderAssociation->fill(0.0); MatrixDoublePtr LaplacianBaltrami(new vnl_matrix(m_NumberCoefficients,m_NumberCoefficients)); LaplacianBaltrami->fill(0.0); MatrixDoublePtr FRTMatrix(new vnl_matrix(m_NumberCoefficients,m_NumberCoefficients)); FRTMatrix->fill(0.0); MatrixDoublePtr SHEigenvalues(new vnl_matrix(m_NumberCoefficients,m_NumberCoefficients)); SHEigenvalues->fill(0.0); // SHBasis-Matrix + LaplacianBaltrami-Matrix + SHOrderAssociationVector ComputeSphericalHarmonicsBasis(Q.get() ,SHBasis.get(), LaplacianBaltrami.get(), SHOrderAssociation.get(), SHEigenvalues.get()); // Compute FunkRadon Transformation Matrix Associated to SHBasis Order lj ComputeFunkRadonTransformationMatrix(SHOrderAssociation.get() ,FRTMatrix.get()); MatrixDoublePtr temp(new vnl_matrix(((SHBasis->transpose()) * (*SHBasis)) + (m_Lambda * (*LaplacianBaltrami)))); InverseMatrixDoublePtr pseudo_inv(new vnl_matrix_inverse((*temp))); MatrixDoublePtr inverse(new vnl_matrix(m_NumberCoefficients,m_NumberCoefficients)); inverse->fill(0.0); (*inverse) = pseudo_inv->inverse(); // ODF Factor ( missing 1/4PI ?? ) double factor = (1.0/(16.0*QBALL_ANAL_RECON_PI*QBALL_ANAL_RECON_PI)); - MatrixDoublePtr TransformationMatrix(new vnl_matrix( factor * ((*FRTMatrix) * ((*SHEigenvalues) * ((*inverse) * (SHBasis->transpose()) ) ) ))); + m_SignalInterpolationReonstructionMatrix = new vnl_matrix( (*inverse) * (SHBasis->transpose()) ); + + MatrixDoublePtr TransformationMatrix(new vnl_matrix( factor * ((*FRTMatrix) * ((*SHEigenvalues) * (*m_SignalInterpolationReonstructionMatrix))) ) ); + + m_CoeffReconstructionMatrix = new vnl_matrix(m_NumberCoefficients,numberOfGradientDirections); - m_CoeffReconstructionMatrix = new vnl_matrix(m_NumberCoefficients,m_NumberOfGradientDirections); // Cast double to float for(int i=0; iodfs later vnl_matrix_fixed* U = itk::PointShell >::DistributePointShell(); m_SphericalHarmonicBasisMatrix = new vnl_matrix(NOdfDirections,m_NumberCoefficients); for(int i=0; i(NOdfDirections,m_NumberOfGradientDirections); - *m_ReconstructionMatrix = (*m_SphericalHarmonicBasisMatrix) * (*m_CoeffReconstructionMatrix); + m_ReconstructionMatrix = new vnl_matrix((*m_SphericalHarmonicBasisMatrix) * (*m_CoeffReconstructionMatrix)); } template< class T, class TG, class TO, int L, int NODF> template< class VNLType > void DiffusionMultiShellQballReconstructionImageFilter ::printMatrix( VNLType * mat ) { std::stringstream stream; for(int i = 0 ; i < mat->rows(); i++) { stream.str(""); for(int j = 0; j < mat->cols(); j++) { stream << (*mat)(i,j) << " "; } MITK_INFO << stream.str(); } } template< class T, class TG, class TO, int L, int NODF> bool DiffusionMultiShellQballReconstructionImageFilter ::CheckDuplicateDiffusionGradients() { GradientDirectionContainerType::ConstIterator gdcit1; GradientDirectionContainerType::ConstIterator gdcit2; for(gdcit1 = this->m_GradientDirectionContainer->Begin(); gdcit1 != this->m_GradientDirectionContainer->End(); ++gdcit1) { for(gdcit2 = this->m_GradientDirectionContainer->Begin(); gdcit2 != this->m_GradientDirectionContainer->End(); ++gdcit2) { if(gdcit1.Value() == gdcit2.Value() && gdcit1.Index() != gdcit2.Index()) { itkWarningMacro( << "Some of the Diffusion Gradients equal each other. Corresponding image data should be averaged before calling this filter." ); return true; } } } return false; } template< class T, class TG, class TO, int L, int NODF> void DiffusionMultiShellQballReconstructionImageFilter ::PrintSelf(std::ostream& os, Indent indent) const { std::locale C("C"); std::locale originalLocale = os.getloc(); os.imbue(C); Superclass::PrintSelf(os,indent); os << indent << "OdfReconstructionMatrix: " << m_ReconstructionMatrix << std::endl; if ( m_GradientDirectionContainer ) { os << indent << "GradientDirectionContainer: " << m_GradientDirectionContainer << std::endl; } else { os << indent << "GradientDirectionContainer: (Gradient directions not set)" << std::endl; } os << indent << "NumberOfGradientDirections: " << m_NumberOfGradientDirections << std::endl; os << indent << "NumberOfBaselineImages: " << m_NumberOfBaselineImages << std::endl; os << indent << "Threshold for reference B0 image: " << m_Threshold << std::endl; os << indent << "BValue: " << m_BValue << std::endl; os.imbue( originalLocale ); } } #endif // __itkDiffusionMultiShellQballReconstructionImageFilter_cpp diff --git a/Modules/DiffusionImaging/Reconstruction/itkDiffusionMultiShellQballReconstructionImageFilter.h b/Modules/DiffusionImaging/Reconstruction/itkDiffusionMultiShellQballReconstructionImageFilter.h index 06a65027f1..62babc6352 100644 --- a/Modules/DiffusionImaging/Reconstruction/itkDiffusionMultiShellQballReconstructionImageFilter.h +++ b/Modules/DiffusionImaging/Reconstruction/itkDiffusionMultiShellQballReconstructionImageFilter.h @@ -1,194 +1,209 @@ /*========================================================================= Program: Medical Imaging & Interaction Toolkit Language: C++ Date: $Date: 2009-07-14 19:11:20 +0200 (Tue, 14 Jul 2009) $ Version: $Revision: 18127 $ Copyright (c) German Cancer Research Center, Division of Medical and Biological Informatics. All rights reserved. See MITKCopyright.txt or http://www.mitk.org/copyright.html for details. This software is distributed WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the above copyright notices for more information. =========================================================================*/ #ifndef __itkDiffusionMultiShellQballReconstructionImageFilter_h_ #define __itkDiffusionMultiShellQballReconstructionImageFilter_h_ #include "itkImageToImageFilter.h" #include "vnl/vnl_vector_fixed.h" #include "vnl/vnl_matrix.h" #include "vnl/algo/vnl_svd.h" #include "itkVectorContainer.h" #include "itkVectorImage.h" #include namespace itk{ /** \class DiffusionMultiShellQballReconstructionImageFilter Aganj_2010 */ template< class TReferenceImagePixelType, class TGradientImagePixelType, class TOdfPixelType, int NOrderL, int NrOdfDirections> class DiffusionMultiShellQballReconstructionImageFilter : public ImageToImageFilter< Image< TReferenceImagePixelType, 3 >, Image< Vector< TOdfPixelType, NrOdfDirections >, 3 > > { public: typedef DiffusionMultiShellQballReconstructionImageFilter Self; typedef SmartPointer Pointer; typedef SmartPointer ConstPointer; typedef ImageToImageFilter< Image< TReferenceImagePixelType, 3>, Image< Vector< TOdfPixelType, NrOdfDirections >, 3 > > Superclass; typedef TReferenceImagePixelType ReferencePixelType; typedef TGradientImagePixelType GradientPixelType; typedef Vector< TOdfPixelType, NrOdfDirections > OdfPixelType; typedef typename Superclass::InputImageType ReferenceImageType; typedef Image< OdfPixelType, 3 > OdfImageType; typedef OdfImageType OutputImageType; typedef TOdfPixelType BZeroPixelType; typedef Image< BZeroPixelType, 3 > BZeroImageType; typedef typename Superclass::OutputImageRegionType OutputImageRegionType; /** Typedef defining one (of the many) gradient images. */ typedef Image< GradientPixelType, 3 > GradientImageType; /** An alternative typedef defining one (of the many) gradient images. * It will be assumed that the vectorImage has the same dimension as the * Reference image and a vector length parameter of \c n (number of * gradient directions)*/ typedef VectorImage< GradientPixelType, 3 > GradientImagesType; /** Holds the ODF reconstruction matrix */ typedef vnl_matrix< TOdfPixelType >* OdfReconstructionMatrixType; - typedef vnl_matrix< double > CoefficientMatrixType; + typedef vnl_matrix< double > * CoefficientMatrixType; /** Holds each magnetic field gradient used to acquire one DWImage */ typedef vnl_vector_fixed< double, 3 > GradientDirectionType; /** Container to hold gradient directions of the 'n' DW measurements */ typedef VectorContainer< unsigned int, GradientDirectionType > GradientDirectionContainerType; typedef std::map > GradientIndexMap; typedef std::map >::iterator GradientIndexMapIteraotr; typedef std::vector IndiciesVector; // --------------------------------------------------------------------------------------------// /** Method for creation through the object factory. */ itkNewMacro(Self); /** Runtime information support. */ itkTypeMacro(DiffusionMultiShellQballReconstructionImageFilter, ImageToImageFilter); /** set method to add gradient directions and its corresponding * image. The image here is a VectorImage. The user is expected to pass the * gradient directions in a container. The ith element of the container * corresponds to the gradient direction of the ith component image the * VectorImage. For the baseline image, a vector of all zeros * should be set.*/ void SetGradientImage( GradientDirectionContainerType *, const GradientImagesType *image); - void SetGradientIndexMap(const GradientIndexMap * gradientIdexMap) + void SetGradientIndexMap(GradientIndexMap * gradientIdexMap) { m_GradientIndexMap = gradientIdexMap; } /** Get reference image */ virtual ReferenceImageType * GetReferenceImage() { return ( static_cast< ReferenceImageType *>(this->ProcessObject::GetInput(0)) ); } /** Return the gradient direction. idx is 0 based */ virtual GradientDirectionType GetGradientDirection( unsigned int idx) const { - if( idx >= m_NumberOfGradientDirections ) + if( idx >= m_GradientDirectionContainer->Size() ) { itkExceptionMacro( << "Gradient direction " << idx << "does not exist" ); } return m_GradientDirectionContainer->ElementAt( idx+1 ); } OdfPixelType Normalize(OdfPixelType odf, typename NumericTraits::AccumulateType b0 ); vnl_vector PreNormalize( vnl_vector vec, typename NumericTraits::AccumulateType b0 ); /** Threshold on the reference image data. The output ODF will be a null * pdf for pixels in the reference image that have a value less than this * threshold. */ itkSetMacro( Threshold, ReferencePixelType ); itkGetMacro( Threshold, ReferencePixelType ); itkGetMacro( BZeroImage, typename BZeroImageType::Pointer); //itkGetMacro( ODFSumImage, typename BlaImage::Pointer); itkSetMacro( BValue, TOdfPixelType); itkSetMacro( Lambda, double ); itkGetMacro( Lambda, double ); itkGetConstReferenceMacro( BValue, TOdfPixelType); protected: DiffusionMultiShellQballReconstructionImageFilter(); ~DiffusionMultiShellQballReconstructionImageFilter() {}; void PrintSelf(std::ostream& os, Indent indent) const; - void ComputeReconstructionMatrix(); + void ComputeReconstructionMatrix(std::vector< unsigned int > gradientIndiciesVector); bool CheckDuplicateDiffusionGradients(); void ComputeSphericalHarmonicsBasis(vnl_matrix* QBallReference, vnl_matrix* SHBasisOutput, vnl_matrix* LaplaciaBaltramiOutput, vnl_vector* SHOrderAssociation , vnl_matrix * SHEigenvalues); void ComputeFunkRadonTransformationMatrix(vnl_vector* SHOrderAssociationReference, vnl_matrix* FRTMatrixOutput ); bool CheckHemisphericalArrangementOfGradientDirections(); void BeforeThreadedGenerateData(); void ThreadedGenerateData( const OutputImageRegionType &outputRegionForThread, int NumberOfThreads ); - + vnl_vector AnalyticalThreeShellParameterEstimation(const IndiciesVector * shell1, const IndiciesVector * shell2, const IndiciesVector * shell3, vnl_vector b); + void StandardOneShellReconstruction(const OutputImageRegionType& outputRegionForThread); + void AnalyticalThreeShellReconstruction(const OutputImageRegionType& outputRegionForThread); private: + enum ReconstructionType + { + Analytical3Shells, + NumericalNShells, + Standard1Shell + }; + OdfReconstructionMatrixType m_ReconstructionMatrix; OdfReconstructionMatrixType m_CoeffReconstructionMatrix; OdfReconstructionMatrixType m_SphericalHarmonicBasisMatrix; + + CoefficientMatrixType m_SignalInterpolationReonstructionMatrix; + + /** container to hold gradient directions */ GradientDirectionContainerType::Pointer m_GradientDirectionContainer; /** Number of gradient measurements */ unsigned int m_NumberOfGradientDirections; /** Number of baseline images */ unsigned int m_NumberOfBaselineImages; /** Threshold on the reference image data */ ReferencePixelType m_Threshold; /** LeBihan's b-value for normalizing tensors */ TOdfPixelType m_BValue; typename BZeroImageType::Pointer m_BZeroImage; GradientIndexMap * m_GradientIndexMap; double m_Lambda; bool m_IsHemisphericalArrangementOfGradientDirections; + bool m_IsArithmeticProgession; + int m_NumberCoefficients; + ReconstructionType m_ReconstructionType; + + + template< class VNLType > void printMatrix( VNLType * mat ); - enum GradientDirections - { - BZero = 0, - AllGradients = -1 - }; + }; } #ifndef ITK_MANUAL_INSTANTIATION #include "itkDiffusionMultiShellQballReconstructionImageFilter.cpp" #endif #endif //__itkDiffusionMultiShellQballReconstructionImageFilter_h_ diff --git a/Plugins/org.mitk.gui.qt.diffusionimaging/src/internal/QmitkQBallReconstructionView.cpp b/Plugins/org.mitk.gui.qt.diffusionimaging/src/internal/QmitkQBallReconstructionView.cpp index 30b85538e8..0962b87199 100644 --- a/Plugins/org.mitk.gui.qt.diffusionimaging/src/internal/QmitkQBallReconstructionView.cpp +++ b/Plugins/org.mitk.gui.qt.diffusionimaging/src/internal/QmitkQBallReconstructionView.cpp @@ -1,999 +1,1004 @@ /*========================================================================= Program: Medical Imaging & Interaction Toolkit Module: $RCSfile$ Language: C++ Date: $Date: 2009-05-28 17:19:30 +0200 (Do, 28 Mai 2009) $ Version: $Revision: 17495 $ Copyright (c) German Cancer Research Center, Division of Medical and Biological Informatics. All rights reserved. See MITKCopyright.txt or http://www.mitk.org/copyright.html for details. This software is distributed WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the above copyright notices for more information. =========================================================================*/ //#define MBILOG_ENABLE_DEBUG #include "QmitkQBallReconstructionView.h" #include "mitkDiffusionImagingConfigure.h" // qt includes #include // itk includes #include "itkTimeProbe.h" // mitk includes #include "mitkProgressBar.h" #include "mitkStatusBar.h" #include "mitkNodePredicateDataType.h" #include "QmitkDataStorageComboBox.h" #include "QmitkStdMultiWidget.h" #include "itkDiffusionQballReconstructionImageFilter.h" #include "itkAnalyticalDiffusionQballReconstructionImageFilter.h" #include "itkDiffusionMultiShellQballReconstructionImageFilter.h" #include "itkVectorContainer.h" #include "mitkQBallImage.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 "berryIStructuredSelection.h" #include "berryIWorkbenchWindow.h" #include "berryISelectionService.h" #include const std::string QmitkQBallReconstructionView::VIEW_ID = "org.mitk.views.qballreconstruction"; #define DI_INFO MITK_INFO("DiffusionImaging") typedef float TTensorPixelType; const int QmitkQBallReconstructionView::nrconvkernels = 252; using namespace berry; struct QbrSelListener : ISelectionListener { berryObjectMacro(QbrSelListener); QbrSelListener(QmitkQBallReconstructionView* view) { m_View = view; } void DoSelectionChanged(ISelection::ConstPointer selection) { // save current selection in member variable m_View->m_CurrentSelection = selection.Cast(); // do something with the selected items if(m_View->m_CurrentSelection) { bool foundDwiVolume = false; // iterate selection for (IStructuredSelection::iterator i = m_View->m_CurrentSelection->Begin(); i != m_View->m_CurrentSelection->End(); ++i) { // extract datatree node if (mitk::DataNodeObject::Pointer nodeObj = i->Cast()) { mitk::DataNode::Pointer node = nodeObj->GetDataNode(); // only look at interesting types if(QString("DiffusionImage").compare(node->GetData()->GetNameOfClass())==0) { foundDwiVolume = true; } } } m_View->m_Controls->m_ButtonStandard->setEnabled(foundDwiVolume); m_View->m_Controls->m_groupQShellSelecton->setEnabled(foundDwiVolume); } } void SelectionChanged(IWorkbenchPart::Pointer part, ISelection::ConstPointer selection) { // check, if selection comes from datamanager if (part) { QString partname(part->GetPartName().c_str()); if(partname.compare("Datamanager")==0) { // apply selection DoSelectionChanged(selection); } } } QmitkQBallReconstructionView* m_View; }; QmitkQBallReconstructionView::QmitkQBallReconstructionView() : QmitkFunctionality(), m_Controls(NULL), m_MultiWidget(NULL) { } QmitkQBallReconstructionView::QmitkQBallReconstructionView(const QmitkQBallReconstructionView& other) { Q_UNUSED(other); throw std::runtime_error("Copy constructor not implemented"); } /*//void QmitkQBallReconstructionView::OpactiyChanged(int value) //{ // if (m_CurrentSelection) // { // if (mitk::DataNodeObject::Pointer nodeObj = m_CurrentSelection->Begin()->Cast()) // { // mitk::DataNode::Pointer node = nodeObj->GetDataNode(); // if(QString("DiffusionImage").compare(node->GetData()->GetNameOfClass())==0) // { // node->SetIntProperty("DisplayChannel", value); // mitk::RenderingManager::GetInstance()->RequestUpdateAll(); // } // } // } //} // //void QmitkQBallReconstructionView::OpactiyActionChanged() //{ // if (m_CurrentSelection) // { // if (mitk::DataNodeObject::Pointer nodeObj = m_CurrentSelection->Begin()->Cast()) // { // mitk::DataNode::Pointer node = nodeObj->GetDataNode(); // if(QString("DiffusionImage").compare(node->GetData()->GetNameOfClass())==0) // { // int displayChannel = 0.0; // if(node->GetIntProperty("DisplayChannel", displayChannel)) // { // m_OpacitySlider->setValue(displayChannel); // } // } // } // } // // MITK_INFO << "changed"; //} */ QmitkQBallReconstructionView::~QmitkQBallReconstructionView() { this->GetSite()->GetWorkbenchWindow()->GetSelectionService()->RemovePostSelectionListener(/*"org.mitk.views.datamanager",*/ m_SelListener); } void QmitkQBallReconstructionView::CreateQtPartControl(QWidget *parent) { if (!m_Controls) { // create GUI widgets m_Controls = new Ui::QmitkQBallReconstructionViewControls; m_Controls->setupUi(parent); this->CreateConnections(); QStringList items; items << "2" << "4" << "6" << "8" << "10" << "12"; m_Controls->m_QBallReconstructionMaxLLevelComboBox->addItems(items); m_Controls->m_QBallReconstructionMaxLLevelComboBox->setCurrentIndex(1); MethodChoosen(m_Controls->m_QBallReconstructionMethodComboBox->currentIndex()); #ifndef DIFFUSION_IMAGING_EXTENDED m_Controls->m_QBallReconstructionMethodComboBox->removeItem(3); #endif AdvancedCheckboxClicked(); // define data type for combobox //m_Controls->m_ImageSelector->SetDataStorage( this->GetDefaultDataStorage() ); //m_Controls->m_ImageSelector->SetPredicate( mitk::NodePredicateDataType::New("DiffusionImage") ); } m_SelListener = berry::ISelectionListener::Pointer(new QbrSelListener(this)); this->GetSite()->GetWorkbenchWindow()->GetSelectionService()->AddPostSelectionListener(/*"org.mitk.views.datamanager",*/ m_SelListener); berry::ISelection::ConstPointer sel( this->GetSite()->GetWorkbenchWindow()->GetSelectionService()->GetSelection("org.mitk.views.datamanager")); m_CurrentSelection = sel.Cast(); m_SelListener.Cast()->DoSelectionChanged(sel); } void QmitkQBallReconstructionView::StdMultiWidgetAvailable (QmitkStdMultiWidget &stdMultiWidget) { m_MultiWidget = &stdMultiWidget; } void QmitkQBallReconstructionView::StdMultiWidgetNotAvailable() { m_MultiWidget = NULL; } void QmitkQBallReconstructionView::CreateConnections() { if ( m_Controls ) { connect( (QObject*)(m_Controls->m_ButtonStandard), SIGNAL(clicked()), this, SLOT(ReconstructStandard()) ); connect( (QObject*)(m_Controls->m_AdvancedCheckbox), SIGNAL(clicked()), this, SLOT(AdvancedCheckboxClicked()) ); connect( (QObject*)(m_Controls->m_QBallReconstructionMethodComboBox), SIGNAL(currentIndexChanged(int)), this, SLOT(MethodChoosen(int)) ); } } void QmitkQBallReconstructionView::Activated() { QmitkFunctionality::Activated(); berry::ISelection::ConstPointer sel( this->GetSite()->GetWorkbenchWindow()->GetSelectionService()->GetSelection("org.mitk.views.datamanager")); m_CurrentSelection = sel.Cast(); m_SelListener.Cast()->DoSelectionChanged(sel); } void QmitkQBallReconstructionView::Deactivated() { QmitkFunctionality::Deactivated(); } void QmitkQBallReconstructionView::ReconstructStandard() { int index = m_Controls->m_QBallReconstructionMethodComboBox->currentIndex(); #ifndef DIFFUSION_IMAGING_EXTENDED if(index>=3) { index = index + 1; } #endif switch(index) { case 0: { // Numerical Reconstruct(0,0); break; } case 1: { // Standard Reconstruct(1,0); break; } case 2: { // Solid Angle Reconstruct(1,6); break; } case 3: { // Constrained Solid Angle Reconstruct(1,7); break; } case 4: { // ADC Reconstruct(1,4); break; } case 5: { // Raw Signal Reconstruct(1,5); break; } case 6: { // Multi Shell Reconstruct(2,7); break; } } } void QmitkQBallReconstructionView::MethodChoosen(int method) { #ifndef DIFFUSION_IMAGING_EXTENDED if(method>=3) { method = method + 1; } #endif m_Controls->m_groupQShellSelecton->setVisible(false); switch(method) { case 0: m_Controls->m_Description->setText("Numerical recon. (Tuch2004)"); break; case 1: m_Controls->m_Description->setText("Spherical harmonics recon. (Descoteaux2007)"); break; case 2: m_Controls->m_Description->setText("SH recon. with solid angle consideration (Aganj2009)"); break; case 3: m_Controls->m_Description->setText("SH solid angle with non-neg. constraint (Goh2009)"); break; case 4: m_Controls->m_Description->setText("SH recon. of the plain ADC-profiles"); break; case 5: m_Controls->m_Description->setText("SH recon. of the raw diffusion signal"); break; case 6: m_Controls->m_Description->setText("SH solid angle for multiple Q-Shell acquisitions"); if (m_CurrentSelection) { QVBoxLayout * layout = dynamic_cast(m_Controls->m_groupQShellSelecton->layout()); for(int i =0; i < layout->count(); i++){ layout->removeItem(layout->itemAt(i)); } for (IStructuredSelection::iterator i = m_CurrentSelection->Begin(); i != m_CurrentSelection->End(); ++i) { if (mitk::DataNodeObject::Pointer nodeObj = i->Cast()) { mitk::DataNode::Pointer node = nodeObj->GetDataNode(); typedef itk::DiffusionMultiShellQballReconstructionImageFilter::GradientIndexMap GradientIndexMap; - GradientIndexMap * bValueMap = new GradientIndexMap(); + m_GradientIndexMap = new GradientIndexMap(); GradientIndexMap::iterator it; if(QString("DiffusionImage").compare(node->GetData()->GetNameOfClass())==0) { mitk::DiffusionImage* vols = static_cast*>(node->GetData()); - vols->GetBvalueList(bValueMap); + vols->GetBvalueList(m_GradientIndexMap); QString str = QString((dynamic_cast(node->GetPropertyList()->GetProperty("name")))->GetValue()); QLabel *label = new QLabel("Image: " + str ); layout->addWidget(label); - for(it = bValueMap->begin(); it != bValueMap->end(); it ++){ + for(it = m_GradientIndexMap->begin(); it != m_GradientIndexMap->end(); it ++){ if(it->first == 0) continue; QString checkboxName ; checkboxName.setNum(it->first); QCheckBox *bValueCheckBox = new QCheckBox( "B-Value " + checkboxName ); bValueCheckBox->setChecked(true); layout->addWidget(bValueCheckBox); } layout->addStretch(1); } } } } m_Controls->m_groupQShellSelecton->setVisible(true); break; } } void QmitkQBallReconstructionView::AdvancedCheckboxClicked() { bool check = m_Controls-> m_AdvancedCheckbox->isChecked(); m_Controls->m_QBallReconstructionMaxLLevelTextLabel_2->setVisible(check); m_Controls->m_QBallReconstructionMaxLLevelComboBox->setVisible(check); m_Controls->m_QBallReconstructionLambdaTextLabel_2->setVisible(check); m_Controls->m_QBallReconstructionLambdaLineEdit->setVisible(check); m_Controls->m_QBallReconstructionThresholdLabel_2->setVisible(check); m_Controls->m_QBallReconstructionThreasholdEdit->setVisible(check); m_Controls->m_OutputB0Image->setVisible(check); m_Controls->label_2->setVisible(check); //m_Controls->textLabel1_2->setVisible(check); //m_Controls->m_QBallReconstructionLambdaStepLineEdit->setVisible(check); //m_Controls->textLabel1_3->setVisible(check); m_Controls->frame_2->setVisible(check); } void QmitkQBallReconstructionView::Reconstruct(int method, int normalization) { if (m_CurrentSelection) { mitk::DataStorage::SetOfObjects::Pointer set = mitk::DataStorage::SetOfObjects::New(); int at = 0; for (IStructuredSelection::iterator i = m_CurrentSelection->Begin(); i != m_CurrentSelection->End(); ++i) { if (mitk::DataNodeObject::Pointer nodeObj = i->Cast()) { mitk::DataNode::Pointer node = nodeObj->GetDataNode(); if(QString("DiffusionImage").compare(node->GetData()->GetNameOfClass())==0) { set->InsertElement(at++, node); } } } if(method == 0) { NumericalQBallReconstruction(set, normalization); } else { #if BOOST_VERSION / 100000 > 0 #if BOOST_VERSION / 100 % 1000 > 34 if(method == 1) { AnalyticalQBallReconstruction(set, normalization); } if(method == 2) { MultiShellQBallReconstruction(set, normalization); } #else std::cout << "ERROR: Boost 1.35 minimum required" << std::endl; QMessageBox::warning(NULL,"ERROR","Boost 1.35 minimum required"); #endif #else std::cout << "ERROR: Boost 1.35 minimum required" << std::endl; QMessageBox::warning(NULL,"ERROR","Boost 1.35 minimum required"); #endif } } } void QmitkQBallReconstructionView::NumericalQBallReconstruction (mitk::DataStorage::SetOfObjects::Pointer inImages, int normalization) { try { itk::TimeProbe clock; int nrFiles = inImages->size(); if (!nrFiles) return; QString status; mitk::ProgressBar::GetInstance()->AddStepsToDo(nrFiles); mitk::DataStorage::SetOfObjects::const_iterator itemiter( inImages->begin() ); mitk::DataStorage::SetOfObjects::const_iterator itemiterend( inImages->end() ); std::vector nodes; while ( itemiter != itemiterend ) // for all items { mitk::DiffusionImage* vols = static_cast*>( (*itemiter)->GetData()); std::string nodename; (*itemiter)->GetStringProperty("name", nodename); ++itemiter; // QBALL RECONSTRUCTION clock.Start(); MBI_INFO << "QBall reconstruction "; mitk::StatusBar::GetInstance()->DisplayText(status.sprintf( "QBall reconstruction for %s", nodename.c_str()).toAscii()); typedef itk::DiffusionQballReconstructionImageFilter QballReconstructionImageFilterType; QballReconstructionImageFilterType::Pointer filter = QballReconstructionImageFilterType::New(); filter->SetGradientImage( vols->GetDirections(), vols->GetVectorImage() ); filter->SetBValue(vols->GetB_Value()); filter->SetThreshold( m_Controls->m_QBallReconstructionThreasholdEdit->text().toFloat() ); switch(normalization) { case 0: { filter->SetNormalizationMethod(QballReconstructionImageFilterType::QBR_STANDARD); break; } case 1: { filter->SetNormalizationMethod(QballReconstructionImageFilterType::QBR_B_ZERO_B_VALUE); break; } case 2: { filter->SetNormalizationMethod(QballReconstructionImageFilterType::QBR_B_ZERO); break; } case 3: { filter->SetNormalizationMethod(QballReconstructionImageFilterType::QBR_NONE); break; } default: { filter->SetNormalizationMethod(QballReconstructionImageFilterType::QBR_STANDARD); } } filter->Update(); clock.Stop(); MBI_DEBUG << "took " << clock.GetMeanTime() << "s." ; // ODFs TO DATATREE mitk::QBallImage::Pointer image = mitk::QBallImage::New(); image->InitializeByItk( filter->GetOutput() ); //image->SetImportVolume( filter->GetOutput()->GetBufferPointer(), 0, 0, mitk::Image::ImportMemoryManagementType::ManageMemory ); image->SetVolume( filter->GetOutput()->GetBufferPointer() ); mitk::DataNode::Pointer node=mitk::DataNode::New(); node->SetData( image ); QString newname; newname = newname.append(nodename.c_str()); newname = newname.append("_QN%1").arg(normalization); SetDefaultNodeProperties(node, newname.toStdString()); nodes.push_back(node); // B-Zero TO DATATREE if(m_Controls->m_OutputB0Image->isChecked()) { mitk::Image::Pointer image4 = mitk::Image::New(); image4->InitializeByItk( filter->GetBZeroImage().GetPointer() ); image4->SetVolume( filter->GetBZeroImage()->GetBufferPointer() ); mitk::DataNode::Pointer node4=mitk::DataNode::New(); node4->SetData( image4 ); node4->SetProperty( "name", mitk::StringProperty::New( QString(nodename.c_str()).append("_b0").toStdString()) ); nodes.push_back(node4); } mitk::ProgressBar::GetInstance()->Progress(); } std::vector::iterator nodeIt; for(nodeIt = nodes.begin(); nodeIt != nodes.end(); ++nodeIt) GetDefaultDataStorage()->Add(*nodeIt); mitk::StatusBar::GetInstance()->DisplayText(status.sprintf("Finished Processing %d Files", nrFiles).toAscii()); m_MultiWidget->RequestUpdate(); } catch (itk::ExceptionObject &ex) { MBI_INFO << ex ; return ; } } void QmitkQBallReconstructionView::AnalyticalQBallReconstruction( mitk::DataStorage::SetOfObjects::Pointer inImages, int normalization) { try { itk::TimeProbe clock; int nrFiles = inImages->size(); if (!nrFiles) return; std::vector lambdas; float minLambda = m_Controls->m_QBallReconstructionLambdaLineEdit->text().toFloat(); lambdas.push_back(minLambda); int nLambdas = lambdas.size(); QString status; mitk::ProgressBar::GetInstance()->AddStepsToDo(nrFiles*nLambdas); mitk::DataStorage::SetOfObjects::const_iterator itemiter( inImages->begin() ); mitk::DataStorage::SetOfObjects::const_iterator itemiterend( inImages->end() ); std::vector* nodes = new std::vector(); while ( itemiter != itemiterend ) // for all items { mitk::DiffusionImage* vols = static_cast*>( (*itemiter)->GetData()); std::string nodename; (*itemiter)->GetStringProperty("name",nodename); itemiter++; // QBALL RECONSTRUCTION clock.Start(); MBI_INFO << "QBall reconstruction "; mitk::StatusBar::GetInstance()->DisplayText(status.sprintf( "QBall reconstruction for %s", nodename.c_str()).toAscii()); for(int i=0; im_QBallReconstructionMaxLLevelComboBox->currentIndex()) { case 0: { TemplatedAnalyticalQBallReconstruction<2>(vols, currentLambda, nodename, nodes, normalization); break; } case 1: { TemplatedAnalyticalQBallReconstruction<4>(vols, currentLambda, nodename, nodes, normalization); break; } case 2: { TemplatedAnalyticalQBallReconstruction<6>(vols, currentLambda, nodename, nodes, normalization); break; } case 3: { TemplatedAnalyticalQBallReconstruction<8>(vols, currentLambda, nodename, nodes, normalization); break; } case 4: { TemplatedAnalyticalQBallReconstruction<10>(vols, currentLambda, nodename, nodes, normalization); break; } case 5: { TemplatedAnalyticalQBallReconstruction<12>(vols, currentLambda, nodename, nodes, normalization); break; } } clock.Stop(); MBI_DEBUG << "took " << clock.GetMeanTime() << "s." ; mitk::ProgressBar::GetInstance()->Progress(); } } std::vector::iterator nodeIt; for(nodeIt = nodes->begin(); nodeIt != nodes->end(); ++nodeIt) GetDefaultDataStorage()->Add(*nodeIt); m_MultiWidget->RequestUpdate(); mitk::StatusBar::GetInstance()->DisplayText(status.sprintf("Finished Processing %d Files", nrFiles).toAscii()); } catch (itk::ExceptionObject &ex) { MBI_INFO << ex ; return ; } } template void QmitkQBallReconstructionView::TemplatedAnalyticalQBallReconstruction( mitk::DiffusionImage* vols, float lambda, std::string nodename, std::vector* nodes, int normalization) { typedef itk::AnalyticalDiffusionQballReconstructionImageFilter FilterType; typename FilterType::Pointer filter = FilterType::New(); filter->SetGradientImage( vols->GetDirections(), vols->GetVectorImage() ); filter->SetBValue(vols->GetB_Value()); filter->SetThreshold( m_Controls->m_QBallReconstructionThreasholdEdit->text().toFloat() ); filter->SetLambda(lambda); switch(normalization) { case 0: { filter->SetNormalizationMethod(FilterType::QBAR_STANDARD); break; } case 1: { filter->SetNormalizationMethod(FilterType::QBAR_B_ZERO_B_VALUE); break; } case 2: { filter->SetNormalizationMethod(FilterType::QBAR_B_ZERO); break; } case 3: { filter->SetNormalizationMethod(FilterType::QBAR_NONE); break; } case 4: { filter->SetNormalizationMethod(FilterType::QBAR_ADC_ONLY); break; } case 5: { filter->SetNormalizationMethod(FilterType::QBAR_RAW_SIGNAL); break; } case 6: { filter->SetNormalizationMethod(FilterType::QBAR_SOLID_ANGLE); break; } case 7: { filter->SetNormalizationMethod(FilterType::QBAR_NONNEG_SOLID_ANGLE); break; } default: { filter->SetNormalizationMethod(FilterType::QBAR_STANDARD); } } filter->Update(); // ODFs TO DATATREE mitk::QBallImage::Pointer image = mitk::QBallImage::New(); image->InitializeByItk( filter->GetOutput() ); image->SetVolume( filter->GetOutput()->GetBufferPointer() ); mitk::DataNode::Pointer node=mitk::DataNode::New(); node->SetData( image ); QString newname; newname = newname.append(nodename.c_str()); newname = newname.append("_QA%1").arg(normalization); SetDefaultNodeProperties(node, newname.toStdString()); nodes->push_back(node); // mitk::Image::Pointer image5 = mitk::Image::New(); // image5->InitializeByItk( filter->GetODFSumImage().GetPointer() ); // image5->SetVolume( filter->GetODFSumImage()->GetBufferPointer() ); // mitk::DataNode::Pointer node5=mitk::DataNode::New(); // node5->SetData( image5 ); // node5->SetProperty( "name", mitk::StringProperty::New( // QString(nodename.c_str()).append("_ODF").toStdString()) ); // nodes->push_back(node5); // B-Zero TO DATATREE if(m_Controls->m_OutputB0Image->isChecked()) { mitk::Image::Pointer image4 = mitk::Image::New(); image4->InitializeByItk( filter->GetBZeroImage().GetPointer() ); image4->SetVolume( filter->GetBZeroImage()->GetBufferPointer() ); mitk::DataNode::Pointer node4=mitk::DataNode::New(); node4->SetData( image4 ); node4->SetProperty( "name", mitk::StringProperty::New( QString(nodename.c_str()).append("_b0").toStdString()) ); nodes->push_back(node4); } } void QmitkQBallReconstructionView::SetDefaultNodeProperties(mitk::DataNode::Pointer node, std::string name) { node->SetProperty( "ShowMaxNumber", mitk::IntProperty::New( 500 ) ); node->SetProperty( "Scaling", mitk::FloatProperty::New( 1.0 ) ); node->SetProperty( "Normalization", mitk::OdfNormalizationMethodProperty::New()); node->SetProperty( "ScaleBy", mitk::OdfScaleByProperty::New()); node->SetProperty( "IndexParam1", mitk::FloatProperty::New(2)); node->SetProperty( "IndexParam2", mitk::FloatProperty::New(1)); node->SetProperty( "visible", mitk::BoolProperty::New( true ) ); node->SetProperty( "VisibleOdfs", mitk::BoolProperty::New( false ) ); node->SetProperty ("layer", mitk::IntProperty::New(100)); node->SetProperty( "DoRefresh", mitk::BoolProperty::New( true ) ); //node->SetProperty( "opacity", mitk::FloatProperty::New(1.0f) ); node->SetProperty( "name", mitk::StringProperty::New(name) ); } /*//node->SetProperty( "volumerendering", mitk::BoolProperty::New( false ) ); //node->SetProperty( "use color", mitk::BoolProperty::New( true ) ); //node->SetProperty( "texture interpolation", mitk::BoolProperty::New( true ) ); //node->SetProperty( "reslice interpolation", mitk::VtkResliceInterpolationProperty::New() ); //node->SetProperty( "layer", mitk::IntProperty::New(0)); //node->SetProperty( "in plane resample extent by geometry", mitk::BoolProperty::New( false ) ); //node->SetOpacity(1.0f); //node->SetColor(1.0,1.0,1.0); //node->SetVisibility(true); //node->SetProperty( "IsQBallVolume", mitk::BoolProperty::New( true ) ); //mitk::LevelWindowProperty::Pointer levWinProp = mitk::LevelWindowProperty::New(); //mitk::LevelWindow levelwindow; //// levelwindow.SetAuto( image ); //levWinProp->SetLevelWindow( levelwindow ); //node->GetPropertyList()->SetProperty( "levelwindow", levWinProp ); //// add a default rainbow lookup table for color mapping //if(!node->GetProperty("LookupTable")) //{ // mitk::LookupTable::Pointer mitkLut = mitk::LookupTable::New(); // vtkLookupTable* vtkLut = mitkLut->GetVtkLookupTable(); // vtkLut->SetHueRange(0.6667, 0.0); // vtkLut->SetTableRange(0.0, 20.0); // vtkLut->Build(); // mitk::LookupTableProperty::Pointer mitkLutProp = mitk::LookupTableProperty::New(); // mitkLutProp->SetLookupTable(mitkLut); // node->SetProperty( "LookupTable", mitkLutProp ); //} //if(!node->GetProperty("binary")) // node->SetProperty( "binary", mitk::BoolProperty::New( false ) ); //// add a default transfer function //mitk::TransferFunction::Pointer tf = mitk::TransferFunction::New(); //node->SetProperty ( "TransferFunction", mitk::TransferFunctionProperty::New ( tf.GetPointer() ) ); //// set foldername as string property //mitk::StringProperty::Pointer nameProp = mitk::StringProperty::New( name ); //node->SetProperty( "name", nameProp );*/ void QmitkQBallReconstructionView::MultiShellQBallReconstruction(mitk::DataStorage::SetOfObjects::Pointer inImages, int normalization) { try { itk::TimeProbe clock; int nrFiles = inImages->size(); if (!nrFiles) return; std::vector lambdas; float minLambda = m_Controls->m_QBallReconstructionLambdaLineEdit->text().toFloat(); lambdas.push_back(minLambda); int nLambdas = lambdas.size(); QString status; mitk::ProgressBar::GetInstance()->AddStepsToDo(nrFiles*nLambdas); mitk::DataStorage::SetOfObjects::const_iterator itemiter( inImages->begin() ); mitk::DataStorage::SetOfObjects::const_iterator itemiterend( inImages->end() ); std::vector* nodes = new std::vector(); while ( itemiter != itemiterend ) // for all items { mitk::DiffusionImage* vols = static_cast*>( (*itemiter)->GetData()); std::string nodename; (*itemiter)->GetStringProperty("name",nodename); itemiter++; // QBALL RECONSTRUCTION clock.Start(); MBI_INFO << "QBall reconstruction "; mitk::StatusBar::GetInstance()->DisplayText(status.sprintf( "QBall reconstruction for %s", nodename.c_str()).toAscii()); for(int i=0; im_QBallReconstructionMaxLLevelComboBox->currentIndex()) { case 0: { TemplatedMultiShellQBallReconstruction<2>(vols, currentLambda, nodename, nodes, normalization); break; } case 1: { TemplatedMultiShellQBallReconstruction<4>(vols, currentLambda, nodename, nodes, normalization); break; } case 2: { TemplatedMultiShellQBallReconstruction<6>(vols, currentLambda, nodename, nodes, normalization); break; } case 3: { TemplatedMultiShellQBallReconstruction<8>(vols, currentLambda, nodename, nodes, normalization); break; } case 4: { TemplatedMultiShellQBallReconstruction<10>(vols, currentLambda, nodename, nodes, normalization); break; } case 5: { TemplatedMultiShellQBallReconstruction<12>(vols, currentLambda, nodename, nodes, normalization); break; } } clock.Stop(); MBI_DEBUG << "took " << clock.GetMeanTime() << "s." ; mitk::ProgressBar::GetInstance()->Progress(); } } std::vector::iterator nodeIt; for(nodeIt = nodes->begin(); nodeIt != nodes->end(); ++nodeIt) GetDefaultDataStorage()->Add(*nodeIt); m_MultiWidget->RequestUpdate(); mitk::StatusBar::GetInstance()->DisplayText(status.sprintf("Finished Processing %d Files", nrFiles).toAscii()); } catch (itk::ExceptionObject &ex) { MBI_INFO << ex ; return ; } } template void QmitkQBallReconstructionView::TemplatedMultiShellQBallReconstruction(mitk::DiffusionImage* vols, float lambda, std::string nodename, std::vector* nodes, int normalization) { typedef itk::DiffusionMultiShellQballReconstructionImageFilter FilterType; typename FilterType::Pointer filter = FilterType::New(); filter->SetGradientImage( vols->GetDirections(), vols->GetVectorImage() ); filter->SetBValue(vols->GetB_Value()); filter->SetThreshold( m_Controls->m_QBallReconstructionThreasholdEdit->text().toFloat() ); filter->SetLambda(lambda); + //filter->SetGradientIndexMap(m_GradientIndexMap); + + for(std::map >::iterator it = m_GradientIndexMap->begin() ; it != m_GradientIndexMap->end(); it++) + MITK_INFO << (*it).first << " Count " << (*it).second.size(); + filter->Update(); // ODFs TO DATATREE mitk::QBallImage::Pointer image = mitk::QBallImage::New(); image->InitializeByItk( filter->GetOutput() ); image->SetVolume( filter->GetOutput()->GetBufferPointer() ); mitk::DataNode::Pointer node=mitk::DataNode::New(); node->SetData( image ); QString newname; newname = newname.append(nodename.c_str()); newname = newname.append("_QA%1").arg(normalization); SetDefaultNodeProperties(node, newname.toStdString()); nodes->push_back(node); // B-Zero TO DATATREE if(m_Controls->m_OutputB0Image->isChecked()) { mitk::Image::Pointer image4 = mitk::Image::New(); image4->InitializeByItk( filter->GetBZeroImage().GetPointer() ); image4->SetVolume( filter->GetBZeroImage()->GetBufferPointer() ); mitk::DataNode::Pointer node4=mitk::DataNode::New(); node4->SetData( image4 ); node4->SetProperty( "name", mitk::StringProperty::New( QString(nodename.c_str()).append("_b0").toStdString()) ); nodes->push_back(node4); } } diff --git a/Plugins/org.mitk.gui.qt.diffusionimaging/src/internal/QmitkQBallReconstructionView.h b/Plugins/org.mitk.gui.qt.diffusionimaging/src/internal/QmitkQBallReconstructionView.h index e3e3981a01..f8efa140cd 100644 --- a/Plugins/org.mitk.gui.qt.diffusionimaging/src/internal/QmitkQBallReconstructionView.h +++ b/Plugins/org.mitk.gui.qt.diffusionimaging/src/internal/QmitkQBallReconstructionView.h @@ -1,124 +1,128 @@ /*========================================================================= Program: Medical Imaging & Interaction Toolkit Module: $RCSfile$ Language: C++ Date: $Date: 2009-05-28 17:19:30 +0200 (Do, 28 Mai 2009) $ Version: $Revision: 17495 $ Copyright (c) German Cancer Research Center, Division of Medical and Biological Informatics. All rights reserved. See MITKCopyright.txt or http://www.mitk.org/copyright.html for details. This software is distributed WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the above copyright notices for more information. =========================================================================*/ #ifndef _QMITKQBALLRECONSTRUCTIONVIEW_H_INCLUDED #define _QMITKQBALLRECONSTRUCTIONVIEW_H_INCLUDED #include #include #include "ui_QmitkQBallReconstructionViewControls.h" #include "mitkDiffusionImage.h" #include #include #include typedef short DiffusionPixelType; struct QbrSelListener; /*! * \ingroup org_mitk_gui_qt_qballreconstruction_internal * * \brief QmitkQBallReconstructionView * * Document your class here. * * \sa QmitkFunctionality */ class QmitkQBallReconstructionView : public QmitkFunctionality { friend struct QbrSelListener; // this is needed for all Qt objects that should have a MOC object (everything that derives from QObject) Q_OBJECT public: static const std::string VIEW_ID; QmitkQBallReconstructionView(); QmitkQBallReconstructionView(const QmitkQBallReconstructionView& other); virtual ~QmitkQBallReconstructionView(); virtual void CreateQtPartControl(QWidget *parent); /// \brief Creation of the connections of main and control widget virtual void CreateConnections(); /// \brief Called when the functionality is activated virtual void Activated(); virtual void Deactivated(); virtual void StdMultiWidgetAvailable (QmitkStdMultiWidget &stdMultiWidget); virtual void StdMultiWidgetNotAvailable(); static const int nrconvkernels; + + protected slots: void ReconstructStandard(); //void ReconstructNormalized1(); //void ReconstructNormalized2(); //void ReconstructNonNormalized(); //void AnalyticallyReconstructStandard(); //void AnalyticallyReconstructSolidAngle(); //void AnalyticallyReconstructNonNegSolidAngle(); //void AnalyticallyReconstructAdc(); //void AnalyticallyReconstructRaw(); void AdvancedCheckboxClicked(); void MethodChoosen(int method); void Reconstruct(int method, int normalization); void NumericalQBallReconstruction(mitk::DataStorage::SetOfObjects::Pointer inImages, int normalization); void AnalyticalQBallReconstruction(mitk::DataStorage::SetOfObjects::Pointer inImages, int normalization); void MultiShellQBallReconstruction(mitk::DataStorage::SetOfObjects::Pointer inImages, int normalization); protected: Ui::QmitkQBallReconstructionViewControls* m_Controls; QmitkStdMultiWidget* m_MultiWidget; template void TemplatedAnalyticalQBallReconstruction(mitk::DiffusionImage* vols, float lambda, std::string nodename, std::vector* nodes, int normalization); template void TemplatedMultiShellQBallReconstruction(mitk::DiffusionImage* vols, float lambda, std::string nodename, std::vector* nodes, int normalization); void SetDefaultNodeProperties(mitk::DataNode::Pointer node, std::string name); berry::ISelectionListener::Pointer m_SelListener; berry::IStructuredSelection::ConstPointer m_CurrentSelection; + std::map > * m_GradientIndexMap; + }; #endif // _QMITKQBALLRECONSTRUCTIONVIEW_H_INCLUDED