From ff4abcacee4dbd9d93dfd6ee1ebb15918ebf1f95 Mon Sep 17 00:00:00 2001 From: dsarrut Date: Wed, 30 Jun 2010 06:00:50 +0000 Subject: [PATCH] boolean operator between 2 images (take into account overlapping region) --- itk/clitkBooleanOperatorLabelImageFilter.h | 138 ++++++++++ itk/clitkBooleanOperatorLabelImageFilter.txx | 252 +++++++++++++++++++ 2 files changed, 390 insertions(+) create mode 100644 itk/clitkBooleanOperatorLabelImageFilter.h create mode 100644 itk/clitkBooleanOperatorLabelImageFilter.txx diff --git a/itk/clitkBooleanOperatorLabelImageFilter.h b/itk/clitkBooleanOperatorLabelImageFilter.h new file mode 100644 index 0000000..f81b440 --- /dev/null +++ b/itk/clitkBooleanOperatorLabelImageFilter.h @@ -0,0 +1,138 @@ +/*========================================================================= + Program: vv http://www.creatis.insa-lyon.fr/rio/vv + + Authors belong to: + - University of LYON http://www.universite-lyon.fr/ + - Léon Bérard cancer center http://oncora1.lyon.fnclcc.fr + - CREATIS CNRS laboratory http://www.creatis.insa-lyon.fr + + This software is distributed WITHOUT ANY WARRANTY; without even + the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR + PURPOSE. See the copyright notices for more information. + + It is distributed under dual licence + + - BSD See included LICENSE.txt file + - CeCILL-B http://www.cecill.info/licences/Licence_CeCILL-B_V1-en.html + ======================================================================-====*/ + +#ifndef CLITKBOOLEANOPERATORLABELIMAGEFILTER_H +#define CLITKBOOLEANOPERATORLABELIMAGEFILTER_H + +#include "itkInPlaceImageFilter.h" +#include "itkImageRegionIteratorWithIndex.h" + +namespace clitk { + + //-------------------------------------------------------------------- + /* + Perform boolan operation between two mask images. Like itkAnd and + others, but take care of: + - origin of the images (spacing must be equal) + - in place or not output + - Binary or Label images as inputs. Label i BG only ; Binary is FG only + + - NO THREAD -> I dont know how (yet) to manage two different inputRegionForThread + + */ + //-------------------------------------------------------------------- + + template + class ITK_EXPORT BooleanOperatorLabelImageFilter: + public itk::InPlaceImageFilter { + + public: + /** Standard class typedefs. */ + typedef BooleanOperatorLabelImageFilter Self; + typedef itk::InPlaceImageFilter Superclass; + typedef itk::SmartPointer Pointer; + typedef itk::SmartPointer ConstPointer; + + /** Method for creation through the object factory. */ + itkNewMacro(Self); + + /** Run-time type information (and related methods). */ + itkTypeMacro(BooleanOperatorLabelImageFilter, InPlaceImageFilter); + + /** Some convenient typedefs. */ + typedef TInputImage1 Input1ImageType; + typedef typename Input1ImageType::ConstPointer Input1ImageConstPointer; + typedef typename Input1ImageType::Pointer Input1ImagePointer; + typedef typename Input1ImageType::RegionType Input1ImageRegionType; + typedef typename Input1ImageType::PixelType Input1ImagePixelType; + + typedef TInputImage2 Input2ImageType; + typedef typename Input2ImageType::ConstPointer Input2ImageConstPointer; + typedef typename Input2ImageType::Pointer Input2ImagePointer; + typedef typename Input2ImageType::RegionType Input2ImageRegionType; + typedef typename Input2ImageType::PixelType Input2ImagePixelType; + + typedef TOutputImage OutputImageType; + typedef typename OutputImageType::Pointer OutputImagePointer; + typedef typename OutputImageType::RegionType OutputImageRegionType; + typedef typename OutputImageType::PixelType OutputImagePixelType; + + /** Connect one of the operands for pixel-wise addition */ + void SetInput1( const TInputImage1 * image1); + + /** Connect one of the operands for pixel-wise addition */ + void SetInput2( const TInputImage2 * image2); + + // Set type of operation + typedef enum { + And = 0, + AndNot = 1 + } OperationTypeEnumeration; + itkGetMacro(OperationType, OperationTypeEnumeration); + itkSetMacro(OperationType, OperationTypeEnumeration); + + // LabelImage information (BG and FG) + void SetBackgroundValue1(Input1ImagePixelType p); + void SetBackgroundValue2(Input2ImagePixelType p); + void SetBackgroundValue(OutputImagePixelType p); + void SetForegroundValue(OutputImagePixelType p); + + /** ImageDimension constants */ + itkStaticConstMacro(InputImage1Dimension, unsigned int, TInputImage1::ImageDimension); + itkStaticConstMacro(InputImage2Dimension, unsigned int, TInputImage2::ImageDimension); + itkStaticConstMacro(OutputImageDimension, unsigned int, TOutputImage::ImageDimension); + + protected: + BooleanOperatorLabelImageFilter(); + virtual ~BooleanOperatorLabelImageFilter() {} + + virtual void GenerateOutputInformation(); + virtual void GenerateInputRequestedRegion(); + virtual void GenerateData(); + + virtual void ReleaseInputs() { } // Do not release date to keep input in memory and continue ... + + Input1ImagePixelType mBackgroundValue1; + Input2ImagePixelType mBackgroundValue2; + OutputImagePixelType mBackgroundValue; + OutputImagePixelType mForegroundValue; + + Input1ImageRegionType input1Region; + Input2ImageRegionType input2Region; + OutputImageRegionType outputRegion; + + OperationTypeEnumeration m_OperationType; + + template void LoopAndNot(Iter1 it1, Iter1 it2, Iter2 ot); + template void LoopAnd(Iter1 it1, Iter1 it2, Iter2 ot); + + private: + BooleanOperatorLabelImageFilter(const Self&); //purposely not implemented + void operator=(const Self&); //purposely not implemented + + }; // end class + //-------------------------------------------------------------------- + +} // end namespace clitk +//-------------------------------------------------------------------- + +#ifndef ITK_MANUAL_INSTANTIATION +#include "clitkBooleanOperatorLabelImageFilter.txx" +#endif + +#endif diff --git a/itk/clitkBooleanOperatorLabelImageFilter.txx b/itk/clitkBooleanOperatorLabelImageFilter.txx new file mode 100644 index 0000000..254085b --- /dev/null +++ b/itk/clitkBooleanOperatorLabelImageFilter.txx @@ -0,0 +1,252 @@ +/*========================================================================= + Program: vv http://www.creatis.insa-lyon.fr/rio/vv + + Authors belong to: + - University of LYON http://www.universite-lyon.fr/ + - Léon Bérard cancer center http://oncora1.lyon.fnclcc.fr + - CREATIS CNRS laboratory http://www.creatis.insa-lyon.fr + + This software is distributed WITHOUT ANY WARRANTY; without even + the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR + PURPOSE. See the copyright notices for more information. + + It is distributed under dual licence + + - BSD See included LICENSE.txt file + - CeCILL-B http://www.cecill.info/licences/Licence_CeCILL-B_V1-en.html + ======================================================================-====*/ + +#ifndef CLITKBOOLEANOPERATORLABELIMAGEFILTER_TXX +#define CLITKBOOLEANOPERATORLABELIMAGEFILTER_TXX + +#include "clitkCommon.h" +#include "clitkBooleanOperatorLabelImageFilter.h" +#include "clitkSegmentationUtils.h" + +namespace clitk { + + //-------------------------------------------------------------------- + template + BooleanOperatorLabelImageFilter:: + BooleanOperatorLabelImageFilter():itk::InPlaceImageFilter() { + this->SetNumberOfRequiredInputs( 2 ); + this->InPlaceOn(); + mBackgroundValue1 = 0; + mBackgroundValue2 = 0; + mBackgroundValue = 0; + mForegroundValue = 1; + m_OperationType = AndNot; + } + //-------------------------------------------------------------------- + + + //-------------------------------------------------------------------- + template + void + BooleanOperatorLabelImageFilter:: + SetInput1(const TInputImage1 * image1) { + // Process object is not const-correct so the const casting is required. + this->SetNthInput(0, const_cast( image1 )); + } + //-------------------------------------------------------------------- + + + //-------------------------------------------------------------------- + template + void + BooleanOperatorLabelImageFilter:: + SetBackgroundValue1(Input1ImagePixelType p) { + mBackgroundValue1 = p; + } + //-------------------------------------------------------------------- + + + //-------------------------------------------------------------------- + template + void + BooleanOperatorLabelImageFilter:: + SetBackgroundValue2(Input2ImagePixelType p) { + mBackgroundValue2 = p; + } + //-------------------------------------------------------------------- + + + //-------------------------------------------------------------------- + template + void + BooleanOperatorLabelImageFilter:: + SetBackgroundValue(OutputImagePixelType p) { + mBackgroundValue = p; + } + //-------------------------------------------------------------------- + + + //-------------------------------------------------------------------- + template + void + BooleanOperatorLabelImageFilter:: + SetForegroundValue(OutputImagePixelType p) { + mForegroundValue = p; + } + //-------------------------------------------------------------------- + + + //-------------------------------------------------------------------- + template + void + BooleanOperatorLabelImageFilter:: + SetInput2(const TInputImage2 * image2) { + // Process object is not const-correct so the const casting is required. + this->SetNthInput(1, const_cast( image2 )); + } + //-------------------------------------------------------------------- + + + //-------------------------------------------------------------------- + template + void + BooleanOperatorLabelImageFilter:: + GenerateOutputInformation() { + + // Get input pointers + Input1ImagePointer input1 = dynamic_cast(itk::ProcessObject::GetInput(0)); + Input2ImagePointer input2 = dynamic_cast(itk::ProcessObject::GetInput(1)); + + // Check spacing + static const unsigned int Dim = Input1ImageType::ImageDimension; + for(unsigned int i=0; iGetSpacing()[i] != input2->GetSpacing()[i]) { + itkExceptionMacro(<< "Input 1&2 must have the same spacing. " << std::endl + << "\t input1 = " << input1->GetSpacing() << std::endl + << "\t input2 = " << input2->GetSpacing() << std::endl); + } + // if (input1->GetLargestPossibleRegion().GetSize()[i] != input2->GetLargestPossibleRegion().GetSize()[i]) { +// itkExceptionMacro(<< "Input 1&2 must have the same size. " << std::endl +// << "\t input1 = " << input1->GetLargestPossibleRegion().GetSize() << std::endl +// << "\t input2 = " << input2->GetLargestPossibleRegion().GetSize() << std::endl); +// } + } + + // Perform default implementation + Superclass::GenerateOutputInformation(); + + // Get output pointer + OutputImagePointer outputImage = this->GetOutput(0); + + // If InPlace, do not create output + // DD(this->GetInPlace()); + if (this->GetInPlace() && this->CanRunInPlace()) { + OutputImagePointer inputAsOutput + = dynamic_cast(const_cast(this->GetInput())); + inputAsOutput->SetRequestedRegion(outputImage->GetLargestPossibleRegion()); + inputAsOutput->SetBufferedRegion(outputImage->GetLargestPossibleRegion()); + // inputAsOutput->SetRegions(outputImage->GetLargestPossibleRegion()); + this->GraftOutput( inputAsOutput ); + } + else { + outputImage->SetRequestedRegion(outputImage->GetLargestPossibleRegion()); + outputImage->SetBufferedRegion(outputImage->GetLargestPossibleRegion()); + outputImage->SetRegions(outputImage->GetLargestPossibleRegion()); + outputImage->Allocate(); + OutputImagePointer inputAsOutput + = dynamic_cast(const_cast(this->GetInput())); + CopyValues(inputAsOutput, outputImage); + } + + // Compute intersection bounding box (in physical coordinate) and regions (in pixel coordinate) + typedef itk::BoundingBox BBType; + typename BBType::Pointer bbInput1 = BBType::New(); + ComputeBBFromImageRegion(input1, input1->GetLargestPossibleRegion(), bbInput1); + typename BBType::Pointer bbInput2 = BBType::New(); + ComputeBBFromImageRegion(input2, input2->GetLargestPossibleRegion(), bbInput2); + typename BBType::Pointer bbOutput = BBType::New(); + ComputeBBFromImageRegion(outputImage, outputImage->GetLargestPossibleRegion(), bbOutput); + + typename BBType::Pointer bb = BBType::New(); + ComputeBBIntersection(bb, bbInput1, bbInput2); + ComputeBBIntersection(bb, bb, bbOutput); + + ComputeRegionFromBB(input1, bb, input1Region); + ComputeRegionFromBB(input2, bb, input2Region); + ComputeRegionFromBB(outputImage, bb, outputRegion); + + } + //-------------------------------------------------------------------- + + + //-------------------------------------------------------------------- + template + void + BooleanOperatorLabelImageFilter:: + GenerateInputRequestedRegion() { + // Call default + itk::InPlaceImageFilter::GenerateInputRequestedRegion(); + // Get input pointers and set requested region to common region + Input1ImagePointer input1 = dynamic_cast(itk::ProcessObject::GetInput(0)); + Input2ImagePointer input2 = dynamic_cast(itk::ProcessObject::GetInput(1)); + input1->SetRequestedRegion(input1Region); + input2->SetRequestedRegion(input2Region); + } + //-------------------------------------------------------------------- + + + //-------------------------------------------------------------------- + template + void + BooleanOperatorLabelImageFilter:: + GenerateData() { + // Get input pointers + Input1ImageConstPointer input1 = dynamic_cast(itk::ProcessObject::GetInput(0)); + Input2ImageConstPointer input2 = dynamic_cast(itk::ProcessObject::GetInput(1)); + + // Get output pointer + OutputImagePointer output = this->GetOutput(0); + + // Get Region iterators + itk::ImageRegionConstIterator it1(input1, input1Region); + itk::ImageRegionConstIterator it2(input2, input2Region); + itk::ImageRegionIterator ot (output, outputRegion); + it1.GoToBegin(); + it2.GoToBegin(); + ot.GoToBegin(); + + switch (m_OperationType) { + case AndNot: LoopAndNot(it1, it2, ot); break; + case And: LoopAnd(it1, it2, ot); break; + } + } + //-------------------------------------------------------------------- + + + //-------------------------------------------------------------------- +#define LOOP_BEGIN(FUNCTION_NAME) \ + template \ + template \ + void \ + BooleanOperatorLabelImageFilter:: \ + FUNCTION_NAME(Iter1 it1, Iter1 it2, Iter2 ot) { \ + while (!ot.IsAtEnd()) { + +#define LOOP_END ++it1; ++it2; ++ot; }} + //-------------------------------------------------------------------- + + + //-------------------------------------------------------------------- + LOOP_BEGIN(LoopAndNot) + if ((it1.Get() != mBackgroundValue1) && (it2.Get() == mBackgroundValue2)) { ot.Set(mForegroundValue); } + else { ot.Set(mBackgroundValue); } + LOOP_END + //-------------------------------------------------------------------- + + + //-------------------------------------------------------------------- + LOOP_BEGIN(LoopAnd) + if ((it1.Get() != mBackgroundValue1) && (it2.Get() != mBackgroundValue2)) { ot.Set(mForegroundValue); } + else { ot.Set(mBackgroundValue); } + LOOP_END + //-------------------------------------------------------------------- + + +}//end clitk + +#endif //#define CLITKBOOLEANOPERATORLABELIMAGEFILTER_TXX -- 2.47.1