From ca6fa6af7cd429b8cb9b7ad96228957f373bb71b Mon Sep 17 00:00:00 2001 From: Matt McCormick Date: Sun, 18 Aug 2019 18:50:48 -0400 Subject: [PATCH] ENH: GradientBasedAngleOfIncidentImageFilter and AcousticImpulseResponseImageFilter. Import from ITKTubeTK/TubeTKLib/USTK. --- ITKKWStyleOverwrite.txt | 1 + .../itkAcousticImpulseResponseImageFilter.h | 135 ++++++++++++ .../itkAcousticImpulseResponseImageFilter.hxx | 134 ++++++++++++ ...GradientBasedAngleOfIncidenceImageFilter.h | 168 +++++++++++++++ ...adientBasedAngleOfIncidenceImageFilter.hxx | 195 ++++++++++++++++++ ...eFilterTest1GradientImageFilter.mha.sha512 | 1 + ...entRecursiveGaussianImageFilter.mha.sha512 | 1 + ...eFilterTest2GradientImageFilter.mha.sha512 | 1 + ...entRecursiveGaussianImageFilter.mha.sha512 | 1 + ...eFilterTest3GradientImageFilter.mha.sha512 | 1 + ...entRecursiveGaussianImageFilter.mha.sha512 | 1 + ...eFilterTest4GradientImageFilter.mha.sha512 | 1 + ...entRecursiveGaussianImageFilter.mha.sha512 | 1 + ...ilterTest2GradientImageFilter.1.mha.sha512 | 1 + ...eFilterTest2GradientImageFilter.mha.sha512 | 1 + ...tRecursiveGaussianImageFilter.1.mha.sha512 | 1 + ...entRecursiveGaussianImageFilter.mha.sha512 | 1 + ...eFilterTest3GradientImageFilter.mha.sha512 | 1 + ...entRecursiveGaussianImageFilter.mha.sha512 | 1 + ...ilterTest4GradientImageFilter.1.mha.sha512 | 1 + ...eFilterTest4GradientImageFilter.mha.sha512 | 1 + ...tRecursiveGaussianImageFilter.1.mha.sha512 | 1 + ...entRecursiveGaussianImageFilter.mha.sha512 | 1 + ...FilterTestGradientImageFilter.1.mha.sha512 | 1 + ...geFilterTestGradientImageFilter.mha.sha512 | 1 + ...tRecursiveGaussianImageFilter.1.mha.sha512 | 1 + ...entRecursiveGaussianImageFilter.mha.sha512 | 1 + test/CMakeLists.txt | 138 +++++++++++++ ...cousticImpedanceImageFilterTest.mha.sha512 | 1 + ...AcousticImpulseResponseImageFilterTest.cxx | 180 ++++++++++++++++ ...ntBasedAngleOfIncidenceImageFilterTest.cxx | 138 +++++++++++++ 31 files changed, 1112 insertions(+) create mode 100644 include/itkAcousticImpulseResponseImageFilter.h create mode 100644 include/itkAcousticImpulseResponseImageFilter.hxx create mode 100644 include/itkGradientBasedAngleOfIncidenceImageFilter.h create mode 100644 include/itkGradientBasedAngleOfIncidenceImageFilter.hxx create mode 100644 test/Baseline/itkAcousticImpulseResponseImageFilterTest1GradientImageFilter.mha.sha512 create mode 100644 test/Baseline/itkAcousticImpulseResponseImageFilterTest1GradientRecursiveGaussianImageFilter.mha.sha512 create mode 100644 test/Baseline/itkAcousticImpulseResponseImageFilterTest2GradientImageFilter.mha.sha512 create mode 100644 test/Baseline/itkAcousticImpulseResponseImageFilterTest2GradientRecursiveGaussianImageFilter.mha.sha512 create mode 100644 test/Baseline/itkAcousticImpulseResponseImageFilterTest3GradientImageFilter.mha.sha512 create mode 100644 test/Baseline/itkAcousticImpulseResponseImageFilterTest3GradientRecursiveGaussianImageFilter.mha.sha512 create mode 100644 test/Baseline/itkAcousticImpulseResponseImageFilterTest4GradientImageFilter.mha.sha512 create mode 100644 test/Baseline/itkAcousticImpulseResponseImageFilterTest4GradientRecursiveGaussianImageFilter.mha.sha512 create mode 100644 test/Baseline/itkGradientBasedAngleOfIncidenceImageFilterTest2GradientImageFilter.1.mha.sha512 create mode 100644 test/Baseline/itkGradientBasedAngleOfIncidenceImageFilterTest2GradientImageFilter.mha.sha512 create mode 100644 test/Baseline/itkGradientBasedAngleOfIncidenceImageFilterTest2GradientRecursiveGaussianImageFilter.1.mha.sha512 create mode 100644 test/Baseline/itkGradientBasedAngleOfIncidenceImageFilterTest2GradientRecursiveGaussianImageFilter.mha.sha512 create mode 100644 test/Baseline/itkGradientBasedAngleOfIncidenceImageFilterTest3GradientImageFilter.mha.sha512 create mode 100644 test/Baseline/itkGradientBasedAngleOfIncidenceImageFilterTest3GradientRecursiveGaussianImageFilter.mha.sha512 create mode 100644 test/Baseline/itkGradientBasedAngleOfIncidenceImageFilterTest4GradientImageFilter.1.mha.sha512 create mode 100644 test/Baseline/itkGradientBasedAngleOfIncidenceImageFilterTest4GradientImageFilter.mha.sha512 create mode 100644 test/Baseline/itkGradientBasedAngleOfIncidenceImageFilterTest4GradientRecursiveGaussianImageFilter.1.mha.sha512 create mode 100644 test/Baseline/itkGradientBasedAngleOfIncidenceImageFilterTest4GradientRecursiveGaussianImageFilter.mha.sha512 create mode 100644 test/Baseline/itkGradientBasedAngleOfIncidenceImageFilterTestGradientImageFilter.1.mha.sha512 create mode 100644 test/Baseline/itkGradientBasedAngleOfIncidenceImageFilterTestGradientImageFilter.mha.sha512 create mode 100644 test/Baseline/itkGradientBasedAngleOfIncidenceImageFilterTestGradientRecursiveGaussianImageFilter.1.mha.sha512 create mode 100644 test/Baseline/itkGradientBasedAngleOfIncidenceImageFilterTestGradientRecursiveGaussianImageFilter.mha.sha512 create mode 100644 test/Input/itkLabelMapToAcousticImpedanceImageFilterTest.mha.sha512 create mode 100644 test/itkAcousticImpulseResponseImageFilterTest.cxx create mode 100644 test/itkGradientBasedAngleOfIncidenceImageFilterTest.cxx diff --git a/ITKKWStyleOverwrite.txt b/ITKKWStyleOverwrite.txt index 8d04f4aa..cc9356b2 100644 --- a/ITKKWStyleOverwrite.txt +++ b/ITKKWStyleOverwrite.txt @@ -3,3 +3,4 @@ include/itkBlockMatchingBayesianRegularizationDisplacementCalculator\.h Comments include/itkBlockMatchingMetricImageToDisplacementCalculator\.h Comments Disable include/itkBlockMatchingOptimizingInterpolationDisplacementCalculator\.h Comments Disable include/itkBlockMatchingParabolicInterpolationDisplacementCalculator\.h Comments Disable +test/itkAcousticImpulseResponseImageFilterTest\.cxx Namespace Disable diff --git a/include/itkAcousticImpulseResponseImageFilter.h b/include/itkAcousticImpulseResponseImageFilter.h new file mode 100644 index 00000000..eacafcd5 --- /dev/null +++ b/include/itkAcousticImpulseResponseImageFilter.h @@ -0,0 +1,135 @@ +/*========================================================================= + * + * Copyright Insight Software Consortium + * + * Licensed under the Apache License, Version 2.0 (the "License"); + * you may not use this file except in compliance with the License. + * You may obtain a copy of the License at + * + * http://www.apache.org/licenses/LICENSE-2.0.txt + * + * Unless required by applicable law or agreed to in writing, software + * distributed under the License is distributed on an "AS IS" BASIS, + * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. + * See the License for the specific language governing permissions and + * limitations under the License. + * + *=========================================================================*/ +#ifndef itkAcousticImpulseResponseImageFilter_h +#define itkAcousticImpulseResponseImageFilter_h + +#include "itkCastImageFilter.h" +#include "itkCovariantVector.h" +#include "itkImageToImageFilter.h" + +namespace itk +{ +/** \class AcousticImpulseResponseImageFilter + * + * \brief Compute the acoustic impulse response along the beam direction. + * + * This filter computes the acoustic pressure impulse response along a beam + * direction. + * + * \f[ + * T( \mathbf{x} ) = + * \cos^n \theta \left( \frac{|\nabla Z( \mathbf{x} )|} + * {2 * Z( \mathbf{x} )} \right ) + * \f] + * + * where: + * + * \f{eqnarray*} + * T( \mathbf{x} ) &=& \mbox{acoustic pressure impulse response} + * n &=& \mbox{specifies angle dependence} + * Z &=& \mbox{acoustic impedance} + * \cos \theta &=& \mbox{angle of incidence} + * \f} + * + * This filter requires two inputs. The first input is the input acoustic + * impedance image. The second input is an angle of incidence image. + * \f$n\f$ is specified with the \c AngleDependence parameter, and + * defaults to one. + * + * It is possible to specify the filter used to calculate the gradient + * magnitude with \c SetGradientMagnitudeFilter. + * + * \ingroup ImageToImageFilter + * \ingroup Ultrasound + */ +template< typename TInputImage, typename TOutputImage, + typename TOperatorValue = float > +class AcousticImpulseResponseImageFilter + : public ImageToImageFilter< TInputImage, TOutputImage > +{ +public: + /** Standard class typedefs. */ + typedef AcousticImpulseResponseImageFilter Self; + typedef ImageToImageFilter< TInputImage, TOutputImage > Superclass; + typedef SmartPointer< Self > Pointer; + typedef SmartPointer< const Self > ConstPointer; + + /** Method for creation through the object factory. */ + itkNewMacro( Self ); + + /** Run-time type information ( and related methods ). */ + itkTypeMacro( AcousticImpulseResponseImageFilter, ImageToImageFilter ); + + /** Some convenient typedefs. */ + typedef TInputImage InputImageType; + typedef typename InputImageType::PointType OriginType; + typedef TOutputImage OutputImageType; + typedef typename OutputImageType::RegionType OutputImageRegionType; + + itkStaticConstMacro( ImageDimension, unsigned int, + InputImageType::ImageDimension ); + + typedef TOperatorValue + OperatorValueType; + typedef Image< OperatorValueType, ImageDimension > + OperatorImageType; + typedef ImageToImageFilter< OperatorImageType, OperatorImageType > + GradientMagnitudeFilterType; + + /** Set/Get the angle dependence term \c n in the class documentation. */ + itkSetMacro( AngleDependence, double ); + itkGetConstMacro( AngleDependence, double ); + + /** Set/Get the filter used to compute the gradient magnitude. */ + itkSetObjectMacro( GradientMagnitudeFilter, GradientMagnitudeFilterType ); + itkGetObjectMacro( GradientMagnitudeFilter, GradientMagnitudeFilterType ); + +protected: + AcousticImpulseResponseImageFilter( void ); + virtual ~AcousticImpulseResponseImageFilter( void ) {} + + void PrintSelf( std::ostream & os, Indent indent ) const override; + + void BeforeThreadedGenerateData( void ) override; + + void ThreadedGenerateData( const OutputImageRegionType & + outputRegionForThread, ThreadIdType threadId ) override; + +private: + // purposely not implemented + AcousticImpulseResponseImageFilter( const Self & ); + // purposely not implemented + void operator=( const Self & ); + + typename GradientMagnitudeFilterType::Pointer m_GradientMagnitudeFilter; + + double m_AngleDependence; + + typedef CastImageFilter< InputImageType, OperatorImageType > + CastImageFilterType; + typename CastImageFilterType::Pointer m_CastImageFilter; + +}; // End class AcousticImpulseResponseImageFilter + +} // End namespace itk + +#ifndef ITK_MANUAL_INSTANTIATION +#include "itkAcousticImpulseResponseImageFilter.hxx" +#endif + +#endif // End !defined( __itkAcousticImpulseResponseImageFilter_h ) diff --git a/include/itkAcousticImpulseResponseImageFilter.hxx b/include/itkAcousticImpulseResponseImageFilter.hxx new file mode 100644 index 00000000..e1651844 --- /dev/null +++ b/include/itkAcousticImpulseResponseImageFilter.hxx @@ -0,0 +1,134 @@ +/*========================================================================= + * + * Copyright Insight Software Consortium + * + * Licensed under the Apache License, Version 2.0 (the "License"); + * you may not use this file except in compliance with the License. + * You may obtain a copy of the License at + * + * http://www.apache.org/licenses/LICENSE-2.0.txt + * + * Unless required by applicable law or agreed to in writing, software + * distributed under the License is distributed on an "AS IS" BASIS, + * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. + * See the License for the specific language governing permissions and + * limitations under the License. + * + *=========================================================================*/ +#ifndef itkAcousticImpulseResponseImageFilter_hxx +#define itkAcousticImpulseResponseImageFilter_hxx + +#include "itkAcousticImpulseResponseImageFilter.h" +#include "itkGradientMagnitudeImageFilter.h" + +namespace itk +{ + +template< typename TInputImage, class TOutputImage, class TOperatorValue > +AcousticImpulseResponseImageFilter< TInputImage, TOutputImage, + TOperatorValue > +::AcousticImpulseResponseImageFilter( void ) + : m_AngleDependence( 1.0 ) +{ + this->SetNumberOfRequiredInputs( 2 ); + + typedef GradientMagnitudeImageFilter< OperatorImageType, + OperatorImageType > DefaultGradientMagnitudeFilterType; + this->m_GradientMagnitudeFilter = DefaultGradientMagnitudeFilterType:: + New(); + + this->m_CastImageFilter = CastImageFilterType::New(); + //Use the ITKv4 Threading Model (call ThreadedGenerateData instead of DynamicThreadedGenerateData) + this->DynamicMultiThreadingOff(); +} + + +template< typename TInputImage, class TOutputImage, class TOperatorValue > +void +AcousticImpulseResponseImageFilter< TInputImage, TOutputImage, + TOperatorValue > +::BeforeThreadedGenerateData( void ) +{ + this->m_CastImageFilter->SetInput( this->GetInput() ); + this->m_GradientMagnitudeFilter->SetNumberOfWorkUnits( + this->GetNumberOfWorkUnits() ); + this->m_GradientMagnitudeFilter->SetInput( + this->m_CastImageFilter->GetOutput() ); + this->m_GradientMagnitudeFilter->Update(); +} + + +template< typename TInputImage, class TOutputImage, class TOperatorValue > +void +AcousticImpulseResponseImageFilter< TInputImage, TOutputImage, + TOperatorValue > +::ThreadedGenerateData( const OutputImageRegionType & outputRegionForThread, + ThreadIdType threadId ) +{ + const InputImageType * input = this->GetInput( 0 ); + const OperatorImageType * angleOfIncidence = this->GetInput( 1 ); + const OperatorImageType * gradientMagnitude = + this->m_GradientMagnitudeFilter->GetOutput(); + OutputImageType * output = this->GetOutput(); + + typedef ImageRegionConstIterator< InputImageType > InputIteratorType; + typedef ImageRegionIterator< OutputImageType > OutputIteratorType; + + typedef ImageRegionConstIterator< OperatorImageType > + AngleOfIncidenceIteratorType; + typedef ImageRegionConstIterator< OperatorImageType > + GradientMagnitudeIteratorType; + + InputIteratorType inputIt( input, outputRegionForThread ); + AngleOfIncidenceIteratorType angleOfIncidenceIt( angleOfIncidence, + outputRegionForThread ); + GradientMagnitudeIteratorType gradientMagnitudeIt( gradientMagnitude, + outputRegionForThread ); + OutputIteratorType outputIt( output, outputRegionForThread ); + + ProgressReporter progress( this, threadId, + outputRegionForThread.GetNumberOfPixels() ); + + if( this->m_AngleDependence == 1.0 ) // avoid the pow call + { + for( inputIt.GoToBegin(), angleOfIncidenceIt.GoToBegin(), + gradientMagnitudeIt.GoToBegin(), outputIt.GoToBegin(); + !outputIt.IsAtEnd(); + ++inputIt, ++angleOfIncidenceIt, ++gradientMagnitudeIt, ++outputIt ) + { + outputIt.Set( angleOfIncidenceIt.Get() * + gradientMagnitudeIt.Get() / ( 2.0 * inputIt.Get() ) ); + } + } + else + { + for( inputIt.GoToBegin(), angleOfIncidenceIt.GoToBegin(), + gradientMagnitudeIt.GoToBegin(), outputIt.GoToBegin(); + !outputIt.IsAtEnd(); + ++inputIt, ++angleOfIncidenceIt, ++gradientMagnitudeIt, ++outputIt ) + { + outputIt.Set( static_cast< typename OutputImageType::PixelType >( + std::pow( static_cast< OperatorValueType >( + angleOfIncidenceIt.Get() ), + static_cast< OperatorValueType >( this->m_AngleDependence * + gradientMagnitudeIt.Get() / ( 2.0 * inputIt.Get() ) ) ) ) ); + } + } +} + + +template< typename TInputImage, typename TOutputImage, typename TOperatorValue > +void +AcousticImpulseResponseImageFilter< TInputImage, TOutputImage, + TOperatorValue > +::PrintSelf( std::ostream & os, Indent indent ) const +{ + Superclass::PrintSelf( os, indent ); + os << indent << "AngleDependence: " + << this->m_AngleDependence + << std::endl; +} + +} // End namespace itk + +#endif // End !defined( __itkAcousticImpulseResponseImageFilter_hxx ) diff --git a/include/itkGradientBasedAngleOfIncidenceImageFilter.h b/include/itkGradientBasedAngleOfIncidenceImageFilter.h new file mode 100644 index 00000000..fd9f4517 --- /dev/null +++ b/include/itkGradientBasedAngleOfIncidenceImageFilter.h @@ -0,0 +1,168 @@ +/*========================================================================= + * + * Copyright Insight Software Consortium + * + * Licensed under the Apache License, Version 2.0 (the "License"); + * you may not use this file except in compliance with the License. + * You may obtain a copy of the License at + * + * http://www.apache.org/licenses/LICENSE-2.0.txt + * + * Unless required by applicable law or agreed to in writing, software + * distributed under the License is distributed on an "AS IS" BASIS, + * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. + * See the License for the specific language governing permissions and + * limitations under the License. + * + *=========================================================================*/ +#ifndef itkGradientBasedAngleOfIncidenceImageFilter_h +#define itkGradientBasedAngleOfIncidenceImageFilter_h + +#include "itkCastImageFilter.h" +#include "itkCovariantVector.h" +#include "itkImageToImageFilter.h" + +namespace itk +{ +/** \class GradientBasedAngleOfIncidenceImageFilter + * \brief Computes cosine of the angle of incidence. + * + * The cosine of the angle of incidence is computed as the angle between the + * ultrasound beam direction and the normal of the "local surface", which is + * computed as the local gradient. + * + * For every input pixel location, the beam direction is computed by + * normalizing + * the vector from that location to the center of rotation of a phased + * array or + * curvilinear array probe, specified with the \c UltrasoundProbeOrigin. + * The + * gradient is computed with a gradient filter of the user's choice -- the + * default is a GradientImageFilter, but a difference filter, e.g. a + * GradientRecursiveGaussianImageFilter could be used instead. + * + * The cosine of the angle of incidence is computed as the dot product of + * the + * two normalized vectors. + * + * \ingroup ImageToImageFilter + * \ingroup Ultrasound + */ +template< typename TInputImage, typename TOutputImage, + typename TOperatorValue = float > +class GradientBasedAngleOfIncidenceImageFilter + : public ImageToImageFilter< TInputImage, TOutputImage > +{ +public: + /** Standard class typedefs. */ + typedef GradientBasedAngleOfIncidenceImageFilter Self; + typedef ImageToImageFilter< TInputImage, TOutputImage > Superclass; + typedef SmartPointer< Self > Pointer; + typedef SmartPointer< const Self > ConstPointer; + + /** Method for creation through the object factory. */ + itkNewMacro( Self ); + + /** Run-time type information ( and related methods ). */ + itkTypeMacro( GradientBasedAngleOfIncidenceImageFilter, + ImageToImageFilter ); + + /** Some convenient typedefs. */ + typedef TInputImage InputImageType; + typedef typename InputImageType::PointType OriginType; + typedef TOutputImage OutputImageType; + typedef typename OutputImageType::RegionType OutputImageRegionType; + + itkStaticConstMacro( ImageDimension, unsigned int, + InputImageType::ImageDimension ); + + typedef TOperatorValue + OperatorValueType; + typedef Image< OperatorValueType, ImageDimension > + OperatorImageType; + typedef CovariantVector< OperatorValueType, ImageDimension > + GradientOutputPixelType; + typedef Image< GradientOutputPixelType, ImageDimension > + GradientOutputImageType; + typedef ImageToImageFilter< OperatorImageType, GradientOutputImageType > + GradientFilterType; + typedef Vector< TOperatorValue, InputImageType::ImageDimension > + BeamDirectionType; + + /** Probe type. Determines how the beam angle is calculated. For + * CURVILINEAR or PHASED, the UltrasoundProbeOrigin must be set. For + * a LINEAR probe, the UltrasoundProbeDirection must be set. */ + typedef enum + { + CURVILINEAR, + PHASED, + LINEAR + } + ProbeType; + + /** Set/Get the probe type. This determines how the beam direction is + * computed. */ + itkSetMacro( UltrasoundProbeType, ProbeType ); + itkGetConstMacro( UltrasoundProbeType, ProbeType ); + + /** Set/Get the location of the ultrasound beam probe center of rotation. + * This is only valid when the UltrasoundProbeType is CURVILINEAR or + * PHASED. */ + itkSetMacro( UltrasoundProbeOrigin, OriginType ); + itkGetConstMacro( UltrasoundProbeOrigin, OriginType ); + + /** Set/Get the direction of the ultrasound beam. This is only valid + * when the UltrasoundProbeType is LINEAR. */ + void SetUltrasoundProbeBeamDirection( const BeamDirectionType & + beamDirection ); + itkGetConstMacro( UltrasoundProbeBeamDirection, BeamDirectionType ); + + /** Set/Get the filter used to calculate the gradients of the input image. + * The default is a simple GradientImageFilter. */ + itkSetObjectMacro( GradientFilter, GradientFilterType ); + itkGetObjectMacro( GradientFilter, GradientFilterType ); + + /** Set/Get the tolerance for the gradient magnitude. If the gradient + * magnitude is below this value, the output is set to zero. */ + itkSetMacro( GradientMagnitudeTolerance, double ); + itkGetConstMacro( GradientMagnitudeTolerance, double ); + +protected: + GradientBasedAngleOfIncidenceImageFilter( void ); + virtual ~GradientBasedAngleOfIncidenceImageFilter( void ) {} + + void PrintSelf( std::ostream & os, Indent indent ) const override; + + void BeforeThreadedGenerateData( void ) override; + + void ThreadedGenerateData( + const OutputImageRegionType & outputRegionForThread, + ThreadIdType threadId ) override; + +private: + //purposely not implemented + GradientBasedAngleOfIncidenceImageFilter( const Self & ); + //purposely not implemented + void operator=( const Self & ); + + typedef CastImageFilter< InputImageType, OperatorImageType > + CastImageFilterType; + typename CastImageFilterType::Pointer m_CastImageFilter; + + + typename GradientFilterType::Pointer m_GradientFilter; + + double m_GradientMagnitudeTolerance; + ProbeType m_UltrasoundProbeType; + OriginType m_UltrasoundProbeOrigin; + BeamDirectionType m_UltrasoundProbeBeamDirection; + +}; // End class GradientBasedAngleOfIncidenceImageFilter + +} // End namespace itk + +#ifndef ITK_MANUAL_INSTANTIATION +#include "itkGradientBasedAngleOfIncidenceImageFilter.hxx" +#endif + +#endif // End !defined( __itkGradientBasedAngleOfIncidenceImageFilter_h ) diff --git a/include/itkGradientBasedAngleOfIncidenceImageFilter.hxx b/include/itkGradientBasedAngleOfIncidenceImageFilter.hxx new file mode 100644 index 00000000..cab7ad17 --- /dev/null +++ b/include/itkGradientBasedAngleOfIncidenceImageFilter.hxx @@ -0,0 +1,195 @@ +/*========================================================================= + * + * Copyright Insight Software Consortium + * + * Licensed under the Apache License, Version 2.0 (the "License"); + * you may not use this file except in compliance with the License. + * You may obtain a copy of the License at + * + * http://www.apache.org/licenses/LICENSE-2.0.txt + * + * Unless required by applicable law or agreed to in writing, software + * distributed under the License is distributed on an "AS IS" BASIS, + * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. + * See the License for the specific language governing permissions and + * limitations under the License. + * + *=========================================================================*/ +#ifndef itkGradientBasedAngleOfIncidenceImageFilter_hxx +#define itkGradientBasedAngleOfIncidenceImageFilter_hxx + +#include "itkGradientBasedAngleOfIncidenceImageFilter.h" + +#include "itkGradientImageFilter.h" +#include "itkImageRegionConstIteratorWithIndex.h" +#include "itkProgressReporter.h" + +namespace itk +{ + +template< typename TInputImage, typename TOutputImage, typename TOperatorValue > +GradientBasedAngleOfIncidenceImageFilter< TInputImage, + TOutputImage, + TOperatorValue > +::GradientBasedAngleOfIncidenceImageFilter() +{ + typedef GradientImageFilter< + OperatorImageType, OperatorValueType, OperatorValueType > + DefaultGradientFilterType; + this->m_GradientFilter = DefaultGradientFilterType::New().GetPointer(); + + this->m_CastImageFilter = CastImageFilterType::New(); + + this->m_UltrasoundProbeType = CURVILINEAR; + this->m_UltrasoundProbeOrigin.Fill( 0.0 ); + this->m_UltrasoundProbeBeamDirection.Fill( 0.0 ); + + this->m_GradientMagnitudeTolerance = 1.0e-7; + //Use the ITKv4 Threading Model (call ThreadedGenerateData instead of DynamicThreadedGenerateData) + this->DynamicMultiThreadingOff(); +} + + +template< typename TInputImage, typename TOutputImage, typename TOperatorValue > +void +GradientBasedAngleOfIncidenceImageFilter< TInputImage, + TOutputImage, + TOperatorValue > +::SetUltrasoundProbeBeamDirection( const BeamDirectionType & beamDirection ) +{ + if( beamDirection != this->m_UltrasoundProbeBeamDirection ) + { + this->m_UltrasoundProbeBeamDirection = beamDirection; + this->m_UltrasoundProbeBeamDirection.Normalize(); + this->Modified(); + } +} + + +template< typename TInputImage, typename TOutputImage, typename TOperatorValue > +void +GradientBasedAngleOfIncidenceImageFilter< TInputImage, + TOutputImage, + TOperatorValue > +::BeforeThreadedGenerateData() +{ + this->m_CastImageFilter->SetInput( this->GetInput() ); + this->m_GradientFilter->SetNumberOfWorkUnits( this->GetNumberOfWorkUnits() ); + this->m_GradientFilter->SetInput( this->m_CastImageFilter->GetOutput() ); + this->m_GradientFilter->Update(); + + if( this->m_UltrasoundProbeType == LINEAR && + this->m_UltrasoundProbeBeamDirection.GetNorm() + == static_cast< OperatorValueType >( 0.0 ) ) + { + itkExceptionMacro( + << "The BeamDirection must be specified with a linear probe." ); + } +} + + +template< typename TInputImage, typename TOutputImage, typename TOperatorValue > +void +GradientBasedAngleOfIncidenceImageFilter< TInputImage, + TOutputImage, + TOperatorValue > +::ThreadedGenerateData( const OutputImageRegionType & outputRegionForThread, + ThreadIdType threadId ) +{ + const InputImageType * input = this->GetInput(); + const GradientOutputImageType * gradient = + this->m_GradientFilter->GetOutput(); + OutputImageType * output = this->GetOutput(); + + const OriginType origin = this->m_UltrasoundProbeOrigin; + const double gradientMagnitudeTolerance = + this->m_GradientMagnitudeTolerance; + + typedef ImageRegionConstIteratorWithIndex< InputImageType > + InputIteratorType; + typedef ImageRegionConstIterator< GradientOutputImageType > + GradientIteratorType; + typedef ImageRegionIterator< OutputImageType > + OutputIteratorType; + + InputIteratorType inputIt( input, outputRegionForThread ); + GradientIteratorType gradientIt( gradient, outputRegionForThread ); + OutputIteratorType outputIt( output, outputRegionForThread ); + + ProgressReporter progress( this, threadId, + outputRegionForThread.GetNumberOfPixels() ); + + BeamDirectionType beamDirection = this->m_UltrasoundProbeBeamDirection; + + for( inputIt.GoToBegin(), gradientIt.GoToBegin(), outputIt.GoToBegin(); + !outputIt.IsAtEnd(); + ++inputIt, ++gradientIt, ++outputIt ) + { + typename InputImageType::IndexType index = inputIt.GetIndex(); + typename InputImageType::PointType point; + input->TransformIndexToPhysicalPoint( index, point ); + + if( this->m_UltrasoundProbeType == CURVILINEAR || + this->m_UltrasoundProbeType == PHASED ) + { + beamDirection = point - origin; + beamDirection.Normalize(); + } + + GradientOutputPixelType gradientPixel = gradientIt.Get(); + const typename GradientOutputPixelType::RealValueType gradientNorm = + gradientPixel.GetNorm(); + gradientPixel /= gradientNorm; + + // output scalar product of the two normalized vectors + typedef typename OutputImageType::PixelType OutputPixelType; + const OutputPixelType outputPixel = gradientPixel * beamDirection; + if( vnl_math_isnan( outputPixel ) + || gradientNorm < gradientMagnitudeTolerance ) + { + outputIt.Set( NumericTraits< OutputPixelType >::Zero ); + } + else + { + outputIt.Set( gradientPixel * beamDirection ); + } + progress.CompletedPixel(); + } +} + + +template< typename TInputImage, typename TOutputImage, typename TOperatorValue > +void +GradientBasedAngleOfIncidenceImageFilter< TInputImage, + TOutputImage, + TOperatorValue > +::PrintSelf( std::ostream & os, Indent indent ) const +{ + Superclass::PrintSelf( os, indent ); + os << indent << "GradientMagnitudeTolerance: " + << this->m_GradientMagnitudeTolerance + << indent << "UltrasoundProbeType: "; + switch( this->m_UltrasoundProbeType ) + { + case CURVILINEAR: + os << "CURVILINEAR"; + break; + case PHASED: + os << "PHASED"; + break; + case LINEAR: + os << "LINEAR"; + break; + default: + os << "INVALID"; + } + os << indent << "UltrasoundProbeOrigin: " + << this->m_UltrasoundProbeOrigin + << indent << "BeamDirection: " + << this->m_UltrasoundProbeBeamDirection + << std::endl; +} + +} // End namespace itk + +#endif // End !defined( __itkGradientBasedAngleOfIncidenceImageFilter_hxx ) diff --git a/test/Baseline/itkAcousticImpulseResponseImageFilterTest1GradientImageFilter.mha.sha512 b/test/Baseline/itkAcousticImpulseResponseImageFilterTest1GradientImageFilter.mha.sha512 new file mode 100644 index 00000000..92c129a9 --- /dev/null +++ b/test/Baseline/itkAcousticImpulseResponseImageFilterTest1GradientImageFilter.mha.sha512 @@ -0,0 +1 @@ +fe5cdaeb057eae3a903df479fccd522f270095d92ed9345be04f72a137abde21c752252b29ee5aeef0d4c154df922240adc2cbed2328d96b6454819515c78bf7 diff --git a/test/Baseline/itkAcousticImpulseResponseImageFilterTest1GradientRecursiveGaussianImageFilter.mha.sha512 b/test/Baseline/itkAcousticImpulseResponseImageFilterTest1GradientRecursiveGaussianImageFilter.mha.sha512 new file mode 100644 index 00000000..81e0f157 --- /dev/null +++ b/test/Baseline/itkAcousticImpulseResponseImageFilterTest1GradientRecursiveGaussianImageFilter.mha.sha512 @@ -0,0 +1 @@ +9e02fa8eb78fd541e4773cca77a7a19fce60248f1ee65dfce23347628c59a3f5a8807d2d20dde1ab17932e5d5dbe187d1b5a5c33f641b3da08cd77473f0cc711 diff --git a/test/Baseline/itkAcousticImpulseResponseImageFilterTest2GradientImageFilter.mha.sha512 b/test/Baseline/itkAcousticImpulseResponseImageFilterTest2GradientImageFilter.mha.sha512 new file mode 100644 index 00000000..e303a72b --- /dev/null +++ b/test/Baseline/itkAcousticImpulseResponseImageFilterTest2GradientImageFilter.mha.sha512 @@ -0,0 +1 @@ +0a249eb2d5a088bc87069d34c9a64f83c03bb3cd3a0ee896f2052935f0607913376efa3b27c2e565730d01e11814d0d4ea129d0e69a29e3844748cf0c73b4be4 diff --git a/test/Baseline/itkAcousticImpulseResponseImageFilterTest2GradientRecursiveGaussianImageFilter.mha.sha512 b/test/Baseline/itkAcousticImpulseResponseImageFilterTest2GradientRecursiveGaussianImageFilter.mha.sha512 new file mode 100644 index 00000000..a64befaa --- /dev/null +++ b/test/Baseline/itkAcousticImpulseResponseImageFilterTest2GradientRecursiveGaussianImageFilter.mha.sha512 @@ -0,0 +1 @@ +3e36147ec2a8f287f8ebc9a20e864794673c2d48d3310dd9c16f26a3773d75b315e650c5f66c732d3a9b2a818c53a5d5b6c105a25f2a917fde315f6168ed5a76 diff --git a/test/Baseline/itkAcousticImpulseResponseImageFilterTest3GradientImageFilter.mha.sha512 b/test/Baseline/itkAcousticImpulseResponseImageFilterTest3GradientImageFilter.mha.sha512 new file mode 100644 index 00000000..b62df819 --- /dev/null +++ b/test/Baseline/itkAcousticImpulseResponseImageFilterTest3GradientImageFilter.mha.sha512 @@ -0,0 +1 @@ +04866b68af3bf907a13bfc2eec326da151da84a22e19a8f12e103fb59ec3b1a64038bf5c6374c53f6cbcf8b4f4bd03aa08b5573dde3a471a0b6916022cde6e8e diff --git a/test/Baseline/itkAcousticImpulseResponseImageFilterTest3GradientRecursiveGaussianImageFilter.mha.sha512 b/test/Baseline/itkAcousticImpulseResponseImageFilterTest3GradientRecursiveGaussianImageFilter.mha.sha512 new file mode 100644 index 00000000..d2953e2e --- /dev/null +++ b/test/Baseline/itkAcousticImpulseResponseImageFilterTest3GradientRecursiveGaussianImageFilter.mha.sha512 @@ -0,0 +1 @@ +ccc87a349696269bc93fac23b3563e6876673e248271225abccd50720215b8b10e66577fb87bc7cc859df9f3db77ac444dcc4da0294983ccdfa9e48f85603fa3 diff --git a/test/Baseline/itkAcousticImpulseResponseImageFilterTest4GradientImageFilter.mha.sha512 b/test/Baseline/itkAcousticImpulseResponseImageFilterTest4GradientImageFilter.mha.sha512 new file mode 100644 index 00000000..0b6e486f --- /dev/null +++ b/test/Baseline/itkAcousticImpulseResponseImageFilterTest4GradientImageFilter.mha.sha512 @@ -0,0 +1 @@ +3cbfeed77ad77974291b509c4bc4fe0bcb2f728f865f05681402c1a739d7b0a6130537e1b4666537457ce69b880fe5e9398a0a7f1866284587f53e949569693a diff --git a/test/Baseline/itkAcousticImpulseResponseImageFilterTest4GradientRecursiveGaussianImageFilter.mha.sha512 b/test/Baseline/itkAcousticImpulseResponseImageFilterTest4GradientRecursiveGaussianImageFilter.mha.sha512 new file mode 100644 index 00000000..856d6c32 --- /dev/null +++ b/test/Baseline/itkAcousticImpulseResponseImageFilterTest4GradientRecursiveGaussianImageFilter.mha.sha512 @@ -0,0 +1 @@ +dd133bf5f1e23436d7852c1326db9d6fde4228513d84c1135993afdaaf248be520b453a24d3dcefc3016908da78e8c7763d791e5b3dbfb66082a6e480866213e diff --git a/test/Baseline/itkGradientBasedAngleOfIncidenceImageFilterTest2GradientImageFilter.1.mha.sha512 b/test/Baseline/itkGradientBasedAngleOfIncidenceImageFilterTest2GradientImageFilter.1.mha.sha512 new file mode 100644 index 00000000..965ae31f --- /dev/null +++ b/test/Baseline/itkGradientBasedAngleOfIncidenceImageFilterTest2GradientImageFilter.1.mha.sha512 @@ -0,0 +1 @@ +2b5a4ff8f3e6f25e2fb570c7b10bdaff70308b0ca7bc5e4e7bdddef58e497c33c5c06dbcfdab9fe577999704ca832fe741518ac4632aef299bb514e44cb74dde diff --git a/test/Baseline/itkGradientBasedAngleOfIncidenceImageFilterTest2GradientImageFilter.mha.sha512 b/test/Baseline/itkGradientBasedAngleOfIncidenceImageFilterTest2GradientImageFilter.mha.sha512 new file mode 100644 index 00000000..50197dc8 --- /dev/null +++ b/test/Baseline/itkGradientBasedAngleOfIncidenceImageFilterTest2GradientImageFilter.mha.sha512 @@ -0,0 +1 @@ +228632bef5ccc9ef33587260cb09a910f6c9c0ace6484df32016902b4cd768549d543f80213a9ecef8759fdb97d93bfe6b965c740864fbc45606c749962e9308 diff --git a/test/Baseline/itkGradientBasedAngleOfIncidenceImageFilterTest2GradientRecursiveGaussianImageFilter.1.mha.sha512 b/test/Baseline/itkGradientBasedAngleOfIncidenceImageFilterTest2GradientRecursiveGaussianImageFilter.1.mha.sha512 new file mode 100644 index 00000000..bb5125e9 --- /dev/null +++ b/test/Baseline/itkGradientBasedAngleOfIncidenceImageFilterTest2GradientRecursiveGaussianImageFilter.1.mha.sha512 @@ -0,0 +1 @@ +7bd4d9cb56120165d2658362124c33d6d9ccb3fc96472f625efbd3ef080def59df42dfd38c59696cfef1297920c3dc2ab37aba28f79466cf5e45c1e389966939 diff --git a/test/Baseline/itkGradientBasedAngleOfIncidenceImageFilterTest2GradientRecursiveGaussianImageFilter.mha.sha512 b/test/Baseline/itkGradientBasedAngleOfIncidenceImageFilterTest2GradientRecursiveGaussianImageFilter.mha.sha512 new file mode 100644 index 00000000..50d07a51 --- /dev/null +++ b/test/Baseline/itkGradientBasedAngleOfIncidenceImageFilterTest2GradientRecursiveGaussianImageFilter.mha.sha512 @@ -0,0 +1 @@ +41397f76f1eb96121c67ab770e627163856666f918a6823c6358fe39ccec6a1bbb53315bb461e274a49c0e5c1a915eed445391928c16c66a9a0d6385e611ab62 diff --git a/test/Baseline/itkGradientBasedAngleOfIncidenceImageFilterTest3GradientImageFilter.mha.sha512 b/test/Baseline/itkGradientBasedAngleOfIncidenceImageFilterTest3GradientImageFilter.mha.sha512 new file mode 100644 index 00000000..c6671d8f --- /dev/null +++ b/test/Baseline/itkGradientBasedAngleOfIncidenceImageFilterTest3GradientImageFilter.mha.sha512 @@ -0,0 +1 @@ +65241d839d24720638e42ea7039c94bb70a21ab9e0358ea88f722e6bf4ecc29e97737c9c2da52db3e06560c9b7257ad0e155ab0a4a2837561af4559970698b8f diff --git a/test/Baseline/itkGradientBasedAngleOfIncidenceImageFilterTest3GradientRecursiveGaussianImageFilter.mha.sha512 b/test/Baseline/itkGradientBasedAngleOfIncidenceImageFilterTest3GradientRecursiveGaussianImageFilter.mha.sha512 new file mode 100644 index 00000000..e6eca043 --- /dev/null +++ b/test/Baseline/itkGradientBasedAngleOfIncidenceImageFilterTest3GradientRecursiveGaussianImageFilter.mha.sha512 @@ -0,0 +1 @@ +800dbe7ef151c1ed8f1f9f3fa76360018cbe1cee375acd4f4750913af797a299c5cd3b676be211b19fb34641bc2d7887f2e089142e45973d796667e054c8b76f diff --git a/test/Baseline/itkGradientBasedAngleOfIncidenceImageFilterTest4GradientImageFilter.1.mha.sha512 b/test/Baseline/itkGradientBasedAngleOfIncidenceImageFilterTest4GradientImageFilter.1.mha.sha512 new file mode 100644 index 00000000..6ce8f222 --- /dev/null +++ b/test/Baseline/itkGradientBasedAngleOfIncidenceImageFilterTest4GradientImageFilter.1.mha.sha512 @@ -0,0 +1 @@ +7efcfd2db2e0bba4615c3c048a5a2553e64ee69847fd093cf49ee7470488caa1afa1d4c38825f874839ccb0ad85b25f5bd197a1b978342b0afc8a7c6e491c654 diff --git a/test/Baseline/itkGradientBasedAngleOfIncidenceImageFilterTest4GradientImageFilter.mha.sha512 b/test/Baseline/itkGradientBasedAngleOfIncidenceImageFilterTest4GradientImageFilter.mha.sha512 new file mode 100644 index 00000000..69f3e6bc --- /dev/null +++ b/test/Baseline/itkGradientBasedAngleOfIncidenceImageFilterTest4GradientImageFilter.mha.sha512 @@ -0,0 +1 @@ +338cb0779270d2d57544f1f40450cb5f2931cbf5841b568d290672a951461744f5e75a47610e8f68df128822d598c350372c88de4eba071d46bd02550d32436d diff --git a/test/Baseline/itkGradientBasedAngleOfIncidenceImageFilterTest4GradientRecursiveGaussianImageFilter.1.mha.sha512 b/test/Baseline/itkGradientBasedAngleOfIncidenceImageFilterTest4GradientRecursiveGaussianImageFilter.1.mha.sha512 new file mode 100644 index 00000000..13933dfb --- /dev/null +++ b/test/Baseline/itkGradientBasedAngleOfIncidenceImageFilterTest4GradientRecursiveGaussianImageFilter.1.mha.sha512 @@ -0,0 +1 @@ +35d22ab161a0b971be207a6390b03213c39e9bbf6c4fbf2c2e456c6dd3cd6737320ff3956b3633bfd140bdb629bc94465dfaa9ffcd27e6b5070666fceb967059 diff --git a/test/Baseline/itkGradientBasedAngleOfIncidenceImageFilterTest4GradientRecursiveGaussianImageFilter.mha.sha512 b/test/Baseline/itkGradientBasedAngleOfIncidenceImageFilterTest4GradientRecursiveGaussianImageFilter.mha.sha512 new file mode 100644 index 00000000..218be9b4 --- /dev/null +++ b/test/Baseline/itkGradientBasedAngleOfIncidenceImageFilterTest4GradientRecursiveGaussianImageFilter.mha.sha512 @@ -0,0 +1 @@ +1627d811cbc539a0514cd09752e73de0fcbf5dfd8a07de6cc759726bf3d270086e7989653e9035e88260908031766780c0a489a9d5cb0b8924c41f2f1ed9454f diff --git a/test/Baseline/itkGradientBasedAngleOfIncidenceImageFilterTestGradientImageFilter.1.mha.sha512 b/test/Baseline/itkGradientBasedAngleOfIncidenceImageFilterTestGradientImageFilter.1.mha.sha512 new file mode 100644 index 00000000..4808086d --- /dev/null +++ b/test/Baseline/itkGradientBasedAngleOfIncidenceImageFilterTestGradientImageFilter.1.mha.sha512 @@ -0,0 +1 @@ +0b29f8057e8a23be31e0468097bc4a9da5ab26a836582bc12c422fcf8ac01eaaa601f6013b733baee3232fa1f76142ed7355d25774a182121f4eb5cec4b02ec9 diff --git a/test/Baseline/itkGradientBasedAngleOfIncidenceImageFilterTestGradientImageFilter.mha.sha512 b/test/Baseline/itkGradientBasedAngleOfIncidenceImageFilterTestGradientImageFilter.mha.sha512 new file mode 100644 index 00000000..bcdbca9d --- /dev/null +++ b/test/Baseline/itkGradientBasedAngleOfIncidenceImageFilterTestGradientImageFilter.mha.sha512 @@ -0,0 +1 @@ +a737d425c13be96e9f6d91ff4592198f241aa1409d24e248019fe772cf04afca79087df5d7b21f79c334f0f742f0a0e2613f1a26eaf0c2d9efec983b86779114 diff --git a/test/Baseline/itkGradientBasedAngleOfIncidenceImageFilterTestGradientRecursiveGaussianImageFilter.1.mha.sha512 b/test/Baseline/itkGradientBasedAngleOfIncidenceImageFilterTestGradientRecursiveGaussianImageFilter.1.mha.sha512 new file mode 100644 index 00000000..5fa86c70 --- /dev/null +++ b/test/Baseline/itkGradientBasedAngleOfIncidenceImageFilterTestGradientRecursiveGaussianImageFilter.1.mha.sha512 @@ -0,0 +1 @@ +7e58ff93e8ca99f53c1c10475f9c44190078af477f42995c30086f550781cf712273ad0aa79b3377554179452027a1ae0a6d0c8f8d8563b4438d058b168abb2b diff --git a/test/Baseline/itkGradientBasedAngleOfIncidenceImageFilterTestGradientRecursiveGaussianImageFilter.mha.sha512 b/test/Baseline/itkGradientBasedAngleOfIncidenceImageFilterTestGradientRecursiveGaussianImageFilter.mha.sha512 new file mode 100644 index 00000000..7c27bc14 --- /dev/null +++ b/test/Baseline/itkGradientBasedAngleOfIncidenceImageFilterTestGradientRecursiveGaussianImageFilter.mha.sha512 @@ -0,0 +1 @@ +64094f9e713d179037dc5cef886f6483e29b4bc705f919f4cc6a1e46b4b66c2dcb9dbe96ec0b5c5f0d5893c9983b835d0c17b1cdf01098b665440ed1ab10c40a diff --git a/test/CMakeLists.txt b/test/CMakeLists.txt index 4bbe0672..c8134356 100644 --- a/test/CMakeLists.txt +++ b/test/CMakeLists.txt @@ -7,12 +7,14 @@ else() endif() set(UltrasoundTests + itkAcousticImpulseResponseImageFilterTest.cxx itkAnalyticSignalImageFilterTest.cxx itkBModeImageFilterTestTiming.cxx itkCurvilinearArraySpecialCoordinatesImageTest.cxx itkCurvilinearArrayUltrasoundImageFileReaderTest.cxx itkComplexToComplex1DFFTImageFilterTest.cxx itkFFT1DImageFilterTest.cxx + itkGradientBasedAngleOfIncidenceImageFilterTest.cxx itkHDF5BModeUltrasoundImageFileReaderTest.cxx itkHDF5UltrasoundImageIOTest.cxx itkHDF5UltrasoundImageIOCanReadITKImageTest.cxx @@ -52,6 +54,142 @@ endif() CreateTestDriver(Ultrasound "${Ultrasound-Test_LIBRARIES}" "${UltrasoundTests}") +itk_add_test(NAME itkAcousticImpulseResponseImageFilterTest1 + COMMAND UltrasoundTestDriver + --compareNumberOfPixelsTolerance 50 + --compare + DATA{Baseline/itkAcousticImpulseResponseImageFilterTest1GradientImageFilter.mha} + ${ITK_TEST_OUTPUT_DIR}/itkAcousticImpulseResponseImageFilterTest1GradientImageFilter.mha + --compare + DATA{Baseline/itkAcousticImpulseResponseImageFilterTest1GradientRecursiveGaussianImageFilter.mha} + ${ITK_TEST_OUTPUT_DIR}/itkAcousticImpulseResponseImageFilterTest1GradientRecursiveGaussianImageFilter.mha + itkAcousticImpulseResponseImageFilterTest + DATA{Input/itkLabelMapToAcousticImpedanceImageFilterTest.mha} + ${ITK_TEST_OUTPUT_DIR}/itkAcousticImpulseResponseImageFilterTest1GradientImageFilter.mha + ${ITK_TEST_OUTPUT_DIR}/itkAcousticImpulseResponseImageFilterTest1GradientRecursiveGaussianImageFilter.mha + 31.8 + -24.8 + 1.0 + ) +itk_add_test(NAME itkAcousticImpulseResponseImageFilterTest2 + COMMAND UltrasoundTestDriver + --compareNumberOfPixelsTolerance 50 + --compare + DATA{Baseline/itkAcousticImpulseResponseImageFilterTest2GradientImageFilter.mha} + ${ITK_TEST_OUTPUT_DIR}/itkAcousticImpulseResponseImageFilterTest2GradientImageFilter.mha + --compare + DATA{Baseline/itkAcousticImpulseResponseImageFilterTest2GradientRecursiveGaussianImageFilter.mha} + ${ITK_TEST_OUTPUT_DIR}/itkAcousticImpulseResponseImageFilterTest2GradientRecursiveGaussianImageFilter.mha + itkAcousticImpulseResponseImageFilterTest + DATA{Input/itkLabelMapToAcousticImpedanceImageFilterTest.mha} + ${ITK_TEST_OUTPUT_DIR}/itkAcousticImpulseResponseImageFilterTest2GradientImageFilter.mha + ${ITK_TEST_OUTPUT_DIR}/itkAcousticImpulseResponseImageFilterTest2GradientRecursiveGaussianImageFilter.mha + 136.5 + -90.8 + 1.0 + ) +itk_add_test(NAME itkAcousticImpulseResponseImageFilterTest3 + COMMAND UltrasoundTestDriver + --compareNumberOfPixelsTolerance 50 + --compare + DATA{Baseline/itkAcousticImpulseResponseImageFilterTest3GradientImageFilter.mha} + ${ITK_TEST_OUTPUT_DIR}/itkAcousticImpulseResponseImageFilterTest3GradientImageFilter.mha + --compare + DATA{Baseline/itkAcousticImpulseResponseImageFilterTest3GradientRecursiveGaussianImageFilter.mha} + ${ITK_TEST_OUTPUT_DIR}/itkAcousticImpulseResponseImageFilterTest3GradientRecursiveGaussianImageFilter.mha + itkAcousticImpulseResponseImageFilterTest + DATA{Input/itkLabelMapToAcousticImpedanceImageFilterTest.mha} + ${ITK_TEST_OUTPUT_DIR}/itkAcousticImpulseResponseImageFilterTest3GradientImageFilter.mha + ${ITK_TEST_OUTPUT_DIR}/itkAcousticImpulseResponseImageFilterTest3GradientRecursiveGaussianImageFilter.mha + 136.5 + -90.8 + 0.5 + ) +itk_add_test(NAME itkAcousticImpulseResponseImageFilterTest4 + COMMAND UltrasoundTestDriver + --compareNumberOfPixelsTolerance 50 + --compare + DATA{Baseline/itkAcousticImpulseResponseImageFilterTest4GradientImageFilter.mha} + ${ITK_TEST_OUTPUT_DIR}/itkAcousticImpulseResponseImageFilterTest4GradientImageFilter.mha + --compare + DATA{Baseline/itkAcousticImpulseResponseImageFilterTest4GradientRecursiveGaussianImageFilter.mha} + ${ITK_TEST_OUTPUT_DIR}/itkAcousticImpulseResponseImageFilterTest4GradientRecursiveGaussianImageFilter.mha + itkAcousticImpulseResponseImageFilterTest + DATA{Input/itkLabelMapToAcousticImpedanceImageFilterTest.mha} + ${ITK_TEST_OUTPUT_DIR}/itkAcousticImpulseResponseImageFilterTest4GradientImageFilter.mha + ${ITK_TEST_OUTPUT_DIR}/itkAcousticImpulseResponseImageFilterTest4GradientRecursiveGaussianImageFilter.mha + 136.5 + -90.8 + 2.0 + ) +itk_add_test(NAME itkGradientBasedAngleOfIncidenceImageFilterTest + COMMAND UltrasoundTestDriver + --compareNumberOfPixelsTolerance 50 + --compare + DATA{Baseline/itkGradientBasedAngleOfIncidenceImageFilterTestGradientImageFilter.mha} + ${ITK_TEST_OUTPUT_DIR}/itkGradientBasedAngleOfIncidenceImageFilterTestGradientImageFilter.mha + --compare + DATA{Baseline/itkGradientBasedAngleOfIncidenceImageFilterTestGradientRecursiveGaussianImageFilter.mha} + ${ITK_TEST_OUTPUT_DIR}/itkGradientBasedAngleOfIncidenceImageFilterTestGradientRecursiveGaussianImageFilter.mha + itkGradientBasedAngleOfIncidenceImageFilterTest + DATA{Input/itkLabelMapToAcousticImpedanceImageFilterTest.mha} + ${ITK_TEST_OUTPUT_DIR}/itkGradientBasedAngleOfIncidenceImageFilterTestGradientImageFilter.mha + ${ITK_TEST_OUTPUT_DIR}/itkGradientBasedAngleOfIncidenceImageFilterTestGradientRecursiveGaussianImageFilter.mha + CURVILINEAR + 31.8 + -24.8 + ) +itk_add_test(NAME itkGradientBasedAngleOfIncidenceImageFilterTest2 + COMMAND UltrasoundTestDriver + --compareNumberOfPixelsTolerance 50 + --compare + DATA{Baseline/itkGradientBasedAngleOfIncidenceImageFilterTest2GradientImageFilter.mha} + ${ITK_TEST_OUTPUT_DIR}/itkGradientBasedAngleOfIncidenceImageFilterTest2GradientImageFilter.mha + --compare + DATA{Baseline/itkGradientBasedAngleOfIncidenceImageFilterTest2GradientRecursiveGaussianImageFilter.mha} + ${ITK_TEST_OUTPUT_DIR}/itkGradientBasedAngleOfIncidenceImageFilterTest2GradientRecursiveGaussianImageFilter.mha + itkGradientBasedAngleOfIncidenceImageFilterTest + DATA{Input/itkLabelMapToAcousticImpedanceImageFilterTest.mha} + ${ITK_TEST_OUTPUT_DIR}/itkGradientBasedAngleOfIncidenceImageFilterTest2GradientImageFilter.mha + ${ITK_TEST_OUTPUT_DIR}/itkGradientBasedAngleOfIncidenceImageFilterTest2GradientRecursiveGaussianImageFilter.mha + CURVILINEAR + 136.5 + -90.8 + ) +itk_add_test(NAME itkGradientBasedAngleOfIncidenceImageFilterTest3 + COMMAND UltrasoundTestDriver + --compareNumberOfPixelsTolerance 50 + --compare + DATA{Baseline/itkGradientBasedAngleOfIncidenceImageFilterTest3GradientImageFilter.mha} + ${ITK_TEST_OUTPUT_DIR}/itkGradientBasedAngleOfIncidenceImageFilterTest3GradientImageFilter.mha + --compare + DATA{Baseline/itkGradientBasedAngleOfIncidenceImageFilterTest3GradientRecursiveGaussianImageFilter.mha} + ${ITK_TEST_OUTPUT_DIR}/itkGradientBasedAngleOfIncidenceImageFilterTest3GradientRecursiveGaussianImageFilter.mha + itkGradientBasedAngleOfIncidenceImageFilterTest + DATA{Input/itkLabelMapToAcousticImpedanceImageFilterTest.mha} + ${ITK_TEST_OUTPUT_DIR}/itkGradientBasedAngleOfIncidenceImageFilterTest3GradientImageFilter.mha + ${ITK_TEST_OUTPUT_DIR}/itkGradientBasedAngleOfIncidenceImageFilterTest3GradientRecursiveGaussianImageFilter.mha + LINEAR + 0.0 + 1.0 + ) +itk_add_test(NAME itkGradientBasedAngleOfIncidenceImageFilterTest4 + COMMAND UltrasoundTestDriver + --compareNumberOfPixelsTolerance 50 + --compare + DATA{Baseline/itkGradientBasedAngleOfIncidenceImageFilterTest4GradientImageFilter.mha} + ${ITK_TEST_OUTPUT_DIR}/itkGradientBasedAngleOfIncidenceImageFilterTest4GradientImageFilter.mha + --compare + DATA{Baseline/itkGradientBasedAngleOfIncidenceImageFilterTest4GradientRecursiveGaussianImageFilter.mha} + ${ITK_TEST_OUTPUT_DIR}/itkGradientBasedAngleOfIncidenceImageFilterTest4GradientRecursiveGaussianImageFilter.mha + itkGradientBasedAngleOfIncidenceImageFilterTest + DATA{Input/itkLabelMapToAcousticImpedanceImageFilterTest.mha} + ${ITK_TEST_OUTPUT_DIR}/itkGradientBasedAngleOfIncidenceImageFilterTest4GradientImageFilter.mha + ${ITK_TEST_OUTPUT_DIR}/itkGradientBasedAngleOfIncidenceImageFilterTest4GradientRecursiveGaussianImageFilter.mha + LINEAR + 0.5 + 0.5 + ) itk_add_test(NAME itkAnalyticSignalImageFilterTest COMMAND UltrasoundTestDriver --compare diff --git a/test/Input/itkLabelMapToAcousticImpedanceImageFilterTest.mha.sha512 b/test/Input/itkLabelMapToAcousticImpedanceImageFilterTest.mha.sha512 new file mode 100644 index 00000000..9e84aa92 --- /dev/null +++ b/test/Input/itkLabelMapToAcousticImpedanceImageFilterTest.mha.sha512 @@ -0,0 +1 @@ +0876c9975d022950c73e617a6117bc3d1614a42754ecdf17a36efdc1dbcb567586fa8846052acdbd5c759126e7a6e758bcef39058c50f6f07c76a946007e9a20 diff --git a/test/itkAcousticImpulseResponseImageFilterTest.cxx b/test/itkAcousticImpulseResponseImageFilterTest.cxx new file mode 100644 index 00000000..296c3768 --- /dev/null +++ b/test/itkAcousticImpulseResponseImageFilterTest.cxx @@ -0,0 +1,180 @@ +/*========================================================================= + * + * Copyright Insight Software Consortium + * + * Licensed under the Apache License, Version 2.0 (the "License"); + * you may not use this file except in compliance with the License. + * You may obtain a copy of the License at + * + * http://www.apache.org/licenses/LICENSE-2.0.txt + * + * Unless required by applicable law or agreed to in writing, software + * distributed under the License is distributed on an "AS IS" BASIS, + * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. + * See the License for the specific language governing permissions and + * limitations under the License. + * + *=========================================================================*/ +#include "itkAcousticImpulseResponseImageFilter.h" +#include "itkGradientBasedAngleOfIncidenceImageFilter.h" + +#include +#include +#include +#include +#include +#include +#include + +// Very rough B-Mode that does not take into account system properties or a +// a realistic image formation process. +template< class TInputImage > +void +RoughBMode( TInputImage * inputImage, const char * fileName ) +{ + using InputImageType = TInputImage; + using AbsAdaptorType = itk::AbsImageAdaptor< InputImageType, + typename InputImageType::PixelType >; + typename AbsAdaptorType::Pointer absAdaptor = AbsAdaptorType::New(); + absAdaptor->SetImage( inputImage ); + + typedef itk::AddImageFilter< + AbsAdaptorType, + InputImageType, + InputImageType > + AddConstantFilterType; + typename AddConstantFilterType::Pointer addConstantFilter = + AddConstantFilterType::New(); + addConstantFilter->SetInput1( absAdaptor ); + addConstantFilter->SetConstant2( 1.0e-12 ); + addConstantFilter->Update(); + + typedef itk::Log10ImageAdaptor< InputImageType, + typename InputImageType::PixelType > + Log10AdaptorType; + typename Log10AdaptorType::Pointer log10Adaptor = + Log10AdaptorType::New(); + log10Adaptor->SetImage( addConstantFilter->GetOutput() ); + + typedef unsigned char OutputPixelType; + typedef itk::Image< OutputPixelType, InputImageType::ImageDimension > + OutputImageType; + typedef itk::IntensityWindowingImageFilter< Log10AdaptorType, OutputImageType > + IntensityWindowingFilterType; + typename IntensityWindowingFilterType::Pointer + intensityWindowingFilter = IntensityWindowingFilterType::New(); + intensityWindowingFilter->SetWindowMinimum( -5.0 ); + intensityWindowingFilter->SetWindowMaximum( 2.0 ); + intensityWindowingFilter->SetInput( log10Adaptor ); + + typedef itk::ImageFileWriter< OutputImageType > WriterType; + typename WriterType::Pointer writer = WriterType::New(); + writer->SetInput( intensityWindowingFilter->GetOutput() ); + writer->SetFileName( fileName ); + writer->Update(); +} + +int itkAcousticImpulseResponseImageFilterTest( int argc, char * argv[] ) +{ + // Argument parsing. + if( argc < 7 ) + { + std::cerr << "Missing arguments." << std::endl; + std::cerr << "Usage: " + << argv[0] + << " inputImage" + << " outputImage" + << " outputImageWithGradientRecursiveGaussian" + << " originX" + << " originY" + << " angleDependence" + << std::endl; + return EXIT_FAILURE; + } + const char * inputImage = argv[1]; + const char * outputImage = argv[2]; + const char * outputImageWithGradientRecursiveGaussian = argv[3]; + const char * originX = argv[4]; + const char * originY = argv[5]; + const char * angleDependence = argv[6]; + + // Types + enum { Dimension = 2 }; + + typedef float PixelType; + typedef itk::Image< PixelType, Dimension > ImageType; + + // Reader + typedef itk::ImageFileReader< ImageType > ReaderType; + ReaderType::Pointer reader = ReaderType::New(); + reader->SetFileName( inputImage ); + + // Calculate the angle of incidence + typedef itk::GradientBasedAngleOfIncidenceImageFilter< ImageType, ImageType > + AngleOfIncidenceFilterType; + AngleOfIncidenceFilterType::Pointer angleOfIncidenceFilter = + AngleOfIncidenceFilterType::New(); + angleOfIncidenceFilter->SetInput( reader->GetOutput() ); + AngleOfIncidenceFilterType::OriginType probeOrigin; + std::istringstream istrm; + istrm.str( originX ); + istrm >> probeOrigin[0]; + istrm.clear(); + istrm.str( originY ); + istrm >> probeOrigin[1]; + angleOfIncidenceFilter->SetUltrasoundProbeOrigin( probeOrigin ); + + // Calculate the acoustic impulse response + typedef itk::AcousticImpulseResponseImageFilter< ImageType, ImageType > + AcousticImpulseResponseFilterType; + AcousticImpulseResponseFilterType::Pointer acousticImpulseResponseFilter = + AcousticImpulseResponseFilterType::New(); + acousticImpulseResponseFilter->SetInput( 0, + reader->GetOutput() ); + acousticImpulseResponseFilter->SetInput( 1, + angleOfIncidenceFilter->GetOutput() ); + istrm.clear(); + istrm.str( angleDependence ); + double angleDependenceDouble; + istrm >> angleDependenceDouble; + acousticImpulseResponseFilter->SetAngleDependence( angleDependenceDouble ); + + // Write image. + try + { + acousticImpulseResponseFilter->Update(); + RoughBMode( acousticImpulseResponseFilter->GetOutput(), outputImage ); + } + catch( itk::ExceptionObject & error ) + { + std::cerr << "Error: " << error << std::endl; + return EXIT_FAILURE; + } + + // Use a different gradient filter. + typedef itk::GradientMagnitudeRecursiveGaussianImageFilter< + AcousticImpulseResponseFilterType::OperatorImageType, + AcousticImpulseResponseFilterType::OperatorImageType + > + GradientMagnitudeRecursiveGaussianFilterType; + GradientMagnitudeRecursiveGaussianFilterType::Pointer + gradientMagnitudeRecursiveGaussianFilter = + GradientMagnitudeRecursiveGaussianFilterType::New(); + gradientMagnitudeRecursiveGaussianFilter->SetSigma( 1.0 ); + acousticImpulseResponseFilter->SetGradientMagnitudeFilter( + gradientMagnitudeRecursiveGaussianFilter ); + + try + { + acousticImpulseResponseFilter->Update(); + RoughBMode( acousticImpulseResponseFilter->GetOutput(), + outputImageWithGradientRecursiveGaussian ); + } + catch( itk::ExceptionObject & error ) + { + std::cerr << "Error: " << error << std::endl; + return EXIT_FAILURE; + } + + return EXIT_SUCCESS; +} diff --git a/test/itkGradientBasedAngleOfIncidenceImageFilterTest.cxx b/test/itkGradientBasedAngleOfIncidenceImageFilterTest.cxx new file mode 100644 index 00000000..371a93fe --- /dev/null +++ b/test/itkGradientBasedAngleOfIncidenceImageFilterTest.cxx @@ -0,0 +1,138 @@ +/*========================================================================= + * + * Copyright Insight Software Consortium + * + * Licensed under the Apache License, Version 2.0 ( the "License" ); + * you may not use this file except in compliance with the License. + * You may obtain a copy of the License at + * + * http://www.apache.org/licenses/LICENSE-2.0.txt + * + * Unless required by applicable law or agreed to in writing, software + * distributed under the License is distributed on an "AS IS" BASIS, + * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. + * See the License for the specific language governing permissions and + * limitations under the License. + * + *=========================================================================*/ + +#include "itkGradientBasedAngleOfIncidenceImageFilter.h" + +#include +#include +#include + +int itkGradientBasedAngleOfIncidenceImageFilterTest( int argc, char * argv[] ) +{ + // Argument parsing. + if( argc < 7 ) + { + std::cerr << "Missing arguments." << std::endl; + std::cerr << "Usage: " + << argv[0] + << " inputImage" + << " outputImage" + << " outputImageWithGradientRecursiveGaussian" + << " probeType( CURVILINEAR|LINEAR )" + << " ( originX|beamDirectionX )" + << " ( originY|beamDirectionY )" + << std::endl; + return EXIT_FAILURE; + } + const char * inputImage = argv[1]; + const char * outputImage = argv[2]; + const char * outputImageWithGradientRecursiveGaussian = argv[3]; + const char * probeType = argv[4]; + const char * originX = argv[5]; + const char * originY = argv[6]; + const char * beamDirectionX = argv[5]; + const char * beamDirectionY = argv[6]; + + // Types + enum { Dimension = 2 }; + + typedef float PixelType; + typedef itk::Image< PixelType, Dimension > ImageType; + + // Reader + typedef itk::ImageFileReader< ImageType > ReaderType; + ReaderType::Pointer reader = ReaderType::New(); + reader->SetFileName( inputImage ); + + // Calculate the angle of incidence + typedef itk::GradientBasedAngleOfIncidenceImageFilter< ImageType, ImageType > + AngleOfIncidenceFilterType; + AngleOfIncidenceFilterType::Pointer angleOfIncidenceFilter = + AngleOfIncidenceFilterType::New(); + angleOfIncidenceFilter->SetInput( reader->GetOutput() ); + AngleOfIncidenceFilterType::OriginType probeOrigin; + std::istringstream istrm; + istrm.str( originX ); + istrm >> probeOrigin[0]; + istrm.clear(); + istrm.str( originY ); + istrm >> probeOrigin[1]; + istrm.clear(); + AngleOfIncidenceFilterType::BeamDirectionType beamDirection; + istrm.str( beamDirectionX ); + istrm >> beamDirection[0]; + istrm.clear(); + istrm.str( beamDirectionY ); + istrm >> beamDirection[1]; + const std::string probeTypeStr( probeType ); + if( probeTypeStr.compare( "CURVILINEAR" ) == 0 ) + { + angleOfIncidenceFilter->SetUltrasoundProbeType( AngleOfIncidenceFilterType::CURVILINEAR ); + angleOfIncidenceFilter->SetUltrasoundProbeOrigin( probeOrigin ); + } + else if( probeTypeStr.compare( "LINEAR" ) == 0 ) + { + angleOfIncidenceFilter->SetUltrasoundProbeType( AngleOfIncidenceFilterType::LINEAR ); + angleOfIncidenceFilter->SetUltrasoundProbeBeamDirection( beamDirection ); + } + else + { + std::cerr << "Unknown probe type" << std::endl; + return EXIT_FAILURE; + } + + // Writer + typedef itk::ImageFileWriter< ImageType > ImageWriterType; + ImageWriterType::Pointer writer = ImageWriterType::New(); + writer->SetFileName( outputImage ); + writer->SetInput( angleOfIncidenceFilter->GetOutput() ); + try + { + writer->Update(); + } + catch( itk::ExceptionObject & error ) + { + std::cerr << "Error: " << error << std::endl; + return EXIT_FAILURE; + } + + // Use a different gradient filter. + typedef itk::GradientRecursiveGaussianImageFilter< + AngleOfIncidenceFilterType::OperatorImageType, + AngleOfIncidenceFilterType::GradientOutputImageType + > + GradientRecursiveGaussianFilterType; + GradientRecursiveGaussianFilterType::Pointer gradientRecursiveGaussianFilter = + GradientRecursiveGaussianFilterType::New(); + gradientRecursiveGaussianFilter->SetSigma( 0.5 ); + angleOfIncidenceFilter->SetGradientFilter( gradientRecursiveGaussianFilter ); + angleOfIncidenceFilter->SetGradientMagnitudeTolerance( 0.5e-3 ); + + writer->SetFileName( outputImageWithGradientRecursiveGaussian ); + try + { + writer->Update(); + } + catch( itk::ExceptionObject & error ) + { + std::cerr << "Error: " << error << std::endl; + return EXIT_FAILURE; + } + + return EXIT_SUCCESS; +}