#ifndef __itkSignedMaurerParallelDistanceMapImageFilter_h #define __itkSignedMaurerParallelDistanceMapImageFilter_h #include "itkImage.h" #include "itkImageToImageFilter.h" #include "itkMultiThreader.h" #include "itkBarrier.h" /** \class SignedMaurerParallelDistanceMapImageFilter * * \brief This filter computes the exact Euclidean distance transform of a * binary image in three dimensions in linear time and is parallelized to * compute it in O(n/p) time for image size n and number of threads p. * * \par Inputs and Outputs * This is an image-to-image filter. It is three-dimensional. The input and * output images must be of the same size. Input should be of an integer * type so as to preserve integer arithmetic as much as possible. The output * type should be "int" if squared signed distance is used and "float" or * "double" if the unsquared distance or distance using spacing is desired. * * The default behavior is to generate the squared signed distance transform, * disregarding image spacing and treating interior distances as negative. * The filter also assumes that the input is not already just the values 1 * and 0. These conventions can be changed. * * \par Parameters * The Set/Get macros for ForegroundValue specify the foreground of input * binary image. This is normally 1 and this is the default. Whether a * conversion to 1 and 0 is needed can be set with SetBinaryConvert, with the * default having it required (true). This is a multithreaded filter and, as * such, supports SetNumberOfThreads. The other parameters are the same as * SignedMaurerDistanceMapImageFilter and DanielssonDistanceImageFilter (this * filter does not return the Voronoi map, however).\ * * \cite C. R. Maurer, Jr., R. Qi, and V. Raghavan, "A Linear Time Algorithm * for Computing Exact Euclidean Distance Transforms of Binary Images in * Arbitrary Dimensions", IEEE - Transactions on Pattern Analysis and Machine * Intelligence, 25(2): 265-270, 2003. * * \ingroup ImageFeatureExtraction * */ namespace itk { template class SignedMaurerParallelDistanceMapImageFilter : public ImageToImageFilter< TInputImage, TOutputImage > { public: /** Extract dimension from input and output image. */ itkStaticConstMacro(InputImageDimension, unsigned int, TInputImage::ImageDimension); itkStaticConstMacro(OutputImageDimension, unsigned int, TOutputImage::ImageDimension); /** Convenient typedefs for simplifying declarations. */ typedef TInputImage InputImageType; typedef TOutputImage OutputImageType; /** Standard class typedefs. */ typedef SignedMaurerParallelDistanceMapImageFilter Self; typedef ImageToImageFilter Superclass; typedef SmartPointer Pointer; typedef SmartPointer ConstPointer; /** Method for creation through the object factory. */ itkNewMacro(Self); /** Image typedef support. */ typedef typename InputImageType::PixelType InputPixelType; typedef typename OutputImageType::PixelType OutputPixelType; typedef typename InputImageType::SizeType InputSizeType; typedef typename OutputImageType::SizeType OutputSizeType; typedef typename InputImageType::IndexType InputIndexType; typedef typename OutputImageType::IndexType OutputIndexType; typedef typename InputImageType::SpacingType InputSpacingType; typedef typename OutputImageType::SpacingType OutputSpacingType; typedef typename InputImageType::RegionType InputImageRegionType; typedef typename OutputImageType::RegionType OutputImageRegionType; /** Set if the distance should be squared. */ itkSetMacro(SquaredDistance, bool); /** Get the distance squared. */ itkGetConstReferenceMacro(SquaredDistance, bool); /** Set On/Off if the distance is squared. */ itkBooleanMacro(SquaredDistance); /** Set if the inside represents positive values in the signed distance * map. By convention ON pixels are treated as inside pixels. */ itkSetMacro(InsideIsPositive, bool); /** Get if the inside represents positive values in the signed distance map. * See GetInsideIsPositive() */ itkGetConstReferenceMacro(InsideIsPositive, bool); /** Set if the inside represents positive values in the signed distance * map. By convention ON pixels are treated as inside pixels. Default is * true. */ itkBooleanMacro(InsideIsPositive); /** Set if image spacing should be used in computing distances. */ itkSetMacro(UseImageSpacing, bool); /** Get whether spacing is used. */ itkGetConstReferenceMacro(UseImageSpacing, bool); /** Set On/Off whether spacing is used. */ itkBooleanMacro(UseImageSpacing); /** Set if distance maps should be signed. */ itkSetMacro(Signed, bool); /** Get whether distance map is signed. */ itkGetConstReferenceMacro(Signed, bool); /** Set On/Off whether distance map is signed. */ itkBooleanMacro(Signed); /** Set if image is not yet values of 1 and 0. */ itkSetMacro(BinaryConvert, bool); /** Get whether image is not yet values of 1 and 0. */ itkGetConstReferenceMacro(BinaryConvert, bool); /** Set On/Off whether image is not yet values of 1 and 0. */ itkBooleanMacro(BinaryConvert); /** * Set the foreground value which defines the object. Usually this * value is = 1. */ itkSetMacro(ForegroundValue, InputPixelType); itkGetConstReferenceMacro(ForegroundValue, InputPixelType); protected: SignedMaurerParallelDistanceMapImageFilter() : m_ForegroundValue(1), m_InsideIsPositive(false), m_SquaredDistance(false), m_UseImageSpacing(false), m_Signed(true), m_BinaryConvert(true) { m_Barrier = Barrier::New(); }; virtual ~SignedMaurerParallelDistanceMapImageFilter() {} void PrintSelf(std::ostream& os, Indent indent) const; void GenerateData(); void WaitForAll(); private: struct ClosestPoint { int x; int y; int z; }; struct MaurerThreadStruct { SignedMaurerParallelDistanceMapImageFilter *Filter; InputPixelType Feature; bool Write; ClosestPoint *Closest; }; void ComputeDistanceTransform(ClosestPoint *closest, InputPixelType featureVal, bool writeOutput); static ITK_THREAD_RETURN_TYPE ComputeDistanceTransformThreaderCallback( void *arg ); virtual void ThreadedComputeDistanceTransform(const OutputImageRegionType& outputZPlanesRegionForThread, const OutputImageRegionType& outputXPlanesRegionForThread, unsigned int threadID, unsigned int totalZ, unsigned int totalX, ClosestPoint *closest, InputPixelType featureVal, bool writeOutput); unsigned int SplitRequestedRegionOnNPlane(unsigned int i, unsigned int num, OutputImageRegionType& splitRegion, unsigned int dimension); void distanceForXPlane(ClosestPoint* closest, InputPixelType featureVal, const OutputImageRegionType& threadOutputRegion, unsigned int threadID); void distanceForYPlane(ClosestPoint* closest, const OutputImageRegionType& threadOutputRegion, unsigned int threadID); void distanceForZPlane(ClosestPoint* closest, const OutputImageRegionType& threadOutputRegion, unsigned int threadID, bool writeOutput); ClosestPoint *m_ClosestExterior; ClosestPoint *m_ClosestInterior; typename Barrier::Pointer m_Barrier; typename InputImageType::ConstPointer m_BaseImage; InputSpacingType m_SpacingSquared; InputSizeType m_Size; InputPixelType m_ForegroundValue; bool m_InsideIsPositive; bool m_UseImageSpacing; bool m_SquaredDistance; bool m_Signed; bool m_BinaryConvert; }; } // end namespace itk #ifndef ITK_MANUAL_INSTANTIATION #include "itkSignedMaurerParallelDistanceMapImageFilter.txx" #endif #endif