X-Git-Url: https://git.creatis.insa-lyon.fr/pubgit/?a=blobdiff_plain;f=itk%2FclitkAddRelativePositionConstraintToLabelImageFilter.txx;h=a3759f49fa2f69606401d0d5289faaf96c70f3d3;hb=3a823a03d4a2217e41e4dd3b05a8ea825bb762a5;hp=4078b8463e0ca5782a5765d6eb450c884151b4be;hpb=de2db2051ac9b17a5e3e2b7dc34f42deea3f1526;p=clitk.git diff --git a/itk/clitkAddRelativePositionConstraintToLabelImageFilter.txx b/itk/clitkAddRelativePositionConstraintToLabelImageFilter.txx index 4078b84..a3759f4 100644 --- a/itk/clitkAddRelativePositionConstraintToLabelImageFilter.txx +++ b/itk/clitkAddRelativePositionConstraintToLabelImageFilter.txx @@ -3,7 +3,7 @@ Authors belong to: - University of LYON http://www.universite-lyon.fr/ - - Léon Bérard cancer center http://oncora1.lyon.fnclcc.fr + - Léon Bérard cancer center http://www.centreleonberard.fr - CREATIS CNRS laboratory http://www.creatis.insa-lyon.fr This software is distributed WITHOUT ANY WARRANTY; without even @@ -14,7 +14,7 @@ - BSD See included LICENSE.txt file - CeCILL-B http://www.cecill.info/licences/Licence_CeCILL-B_V1-en.html - ======================================================================-====*/ + ===========================================================================**/ // clitk #include "clitkCommon.h" @@ -22,66 +22,60 @@ #include "clitkAutoCropFilter.h" #include "clitkResampleImageWithOptionsFilter.h" #include "clitkBooleanOperatorLabelImageFilter.h" +#include "clitkCropLikeImageFilter.h" // itk #include -#include "itkStatisticsLabelMapFilter.h" -#include "itkLabelImageToStatisticsLabelMapFilter.h" -#include "itkRegionOfInterestImageFilter.h" -#include "itkBinaryThresholdImageFilter.h" -#include "itkBinaryErodeImageFilter.h" -#include "itkBinaryBallStructuringElement.h" +#include +#include +#include +#include +#include +#include +#include +#if ITK_VERSION_MAJOR >= 4 + #include +#else + #include +#endif // itk [Bloch et al] #include "RelativePositionPropImageFilter.h" //-------------------------------------------------------------------- -template -clitk::AddRelativePositionConstraintToLabelImageFilter:: +template +clitk::AddRelativePositionConstraintToLabelImageFilter:: AddRelativePositionConstraintToLabelImageFilter(): clitk::FilterBase(), - itk::ImageToImageFilter() + itk::ImageToImageFilter() { this->SetNumberOfRequiredInputs(2); SetFuzzyThreshold(0.6); SetBackgroundValue(0); SetObjectBackgroundValue(0); - SetOrientationType(LeftTo); - ResampleBeforeRelativePositionFilterOn(); + ClearOrientationType(); + IntermediateSpacingFlagOn(); SetIntermediateSpacing(10); - - // Step 1 : resample (option=sampling) - // Step 2 : Pad (no) - // Step 3 : relative position (option = angle) - // Step 4 : Threshold - // Step 5 : Erode for boundary - // Step 6 : resample to initial spacing - // Step 7 : pad if not the same size : it can occur when downsample and upsample - // Step 6: combine input+thresholded relpos - // Step 7: autocrop - - // Step 1 : resample - ResampleBeforeRelativePositionFilterOn(); - SetIntermediateSpacing(10); - SetOrientationType(LeftTo); - // SegmentationStep * step = new SegmentationStep; - // step->SetName("Initial resampling and relative position map"); - // SetStep(step, &Self::ResampleAndRelativePositionMap); - - // Step 3 : threshold + postprocess - // SetFuzzyThreshold(0.4); - // step = new SegmentationStep; - // step->SetName("Map threshold and post-processing"); - // SetStep(step, &Self::MapThresholdAndPostProcessing); + AutoCropFlagOn(); + InverseOrientationFlagOff(); + RemoveObjectFlagOn(); + CombineWithOrFlagOff(); + VerboseStepFlagOff(); + WriteStepFlagOff(); + FuzzyMapOnlyFlagOff(); + FastFlagOff(); + SetRadius(2.0); + SetK1(vcl_acos(-1.0)/2); } //-------------------------------------------------------------------- //-------------------------------------------------------------------- -template +template void -clitk::AddRelativePositionConstraintToLabelImageFilter:: -SetInput(const ImageType * image) { +clitk::AddRelativePositionConstraintToLabelImageFilter:: +SetInput(const ImageType * image) +{ // Process object is not const-correct so the const casting is required. this->SetNthInput(0, const_cast(image)); } @@ -89,10 +83,11 @@ SetInput(const ImageType * image) { //-------------------------------------------------------------------- -template +template void -clitk::AddRelativePositionConstraintToLabelImageFilter:: -SetInputObject(const ImageType * image) { +clitk::AddRelativePositionConstraintToLabelImageFilter:: +SetInputObject(const ImageType * image) +{ // Process object is not const-correct so the const casting is required. this->SetNthInput(1, const_cast(image)); } @@ -100,11 +95,66 @@ SetInputObject(const ImageType * image) { //-------------------------------------------------------------------- -template +template +void +clitk::AddRelativePositionConstraintToLabelImageFilter:: +ClearOrientationType() +{ + m_OrientationTypeString.clear(); + m_OrientationType.clear(); + m_Angle1.clear(); + m_Angle2.clear(); +} +//-------------------------------------------------------------------- + + +//-------------------------------------------------------------------- +template +int +clitk::AddRelativePositionConstraintToLabelImageFilter:: +GetNumberOfAngles() +{ + return m_OrientationType.size(); +} +//-------------------------------------------------------------------- + + +//-------------------------------------------------------------------- +template +void +clitk::AddRelativePositionConstraintToLabelImageFilter:: +AddOrientationTypeString(std::string t) +{ + m_OrientationTypeString.push_back(t); + + if (t == "LeftTo") { AddOrientationType(LeftTo); return; } + if (t == "RightTo") { AddOrientationType(RightTo); return; } + if (t == "AntTo") { AddOrientationType(AntTo); return; } + if (t == "PostTo") { AddOrientationType(PostTo); return; } + if (t == "SupTo") { AddOrientationType(SupTo); return; } + if (t == "InfTo") { AddOrientationType(InfTo); return; } + + if (t == "NotLeftTo") { AddOrientationType(LeftTo); InverseOrientationFlagOn(); return; } + if (t == "NotRightTo") { AddOrientationType(RightTo); InverseOrientationFlagOn(); return; } + if (t == "NotAntTo") { AddOrientationType(AntTo); InverseOrientationFlagOn(); return; } + if (t == "NotPostTo") { AddOrientationType(PostTo); InverseOrientationFlagOn(); return; } + if (t == "NotSupTo") { AddOrientationType(SupTo); InverseOrientationFlagOn(); return; } + if (t == "NotInfTo") { AddOrientationType(InfTo); InverseOrientationFlagOn(); return; } + + if (t == "Angle") return; + + clitkExceptionMacro("Error, you must provide LeftTo,RightTo or AntTo,PostTo or SupTo,InfTo (or NotLeftTo, NotRightTo etc) but you give " << t); +} +//-------------------------------------------------------------------- + + +//-------------------------------------------------------------------- +template void -clitk::AddRelativePositionConstraintToLabelImageFilter:: -GenerateOutputInformation() { - ImagePointer input = dynamic_cast(itk::ProcessObject::GetInput(0)); +clitk::AddRelativePositionConstraintToLabelImageFilter:: +GenerateOutputInformation() +{ + ImagePointer input = dynamic_cast(itk::ProcessObject::GetInput(0)); ImagePointer outputImage = this->GetOutput(0); outputImage->SetRegions(outputImage->GetLargestPossibleRegion()); } @@ -112,15 +162,16 @@ GenerateOutputInformation() { //-------------------------------------------------------------------- -template +template void -clitk::AddRelativePositionConstraintToLabelImageFilter:: -GenerateInputRequestedRegion() { +clitk::AddRelativePositionConstraintToLabelImageFilter:: +GenerateInputRequestedRegion() +{ // Call default - itk::ImageToImageFilter::GenerateInputRequestedRegion(); + itk::ImageToImageFilter::GenerateInputRequestedRegion(); // Get input pointers and set requested region to common region - ImagePointer input1 = dynamic_cast(itk::ProcessObject::GetInput(0)); - ImagePointer input2 = dynamic_cast(itk::ProcessObject::GetInput(1)); + ImagePointer input1 = dynamic_cast(itk::ProcessObject::GetInput(0)); + ImagePointer input2 = dynamic_cast(itk::ProcessObject::GetInput(1)); input1->SetRequestedRegion(input1->GetLargestPossibleRegion()); input2->SetRequestedRegion(input2->GetLargestPossibleRegion()); } @@ -128,41 +179,63 @@ GenerateInputRequestedRegion() { //-------------------------------------------------------------------- -template +template void -clitk::AddRelativePositionConstraintToLabelImageFilter:: -SetAngle1(double a) { - SetOrientationType(Angle); - m_Angle1 = a; +clitk::AddRelativePositionConstraintToLabelImageFilter:: +AddAnglesInRad(double a, double b) +{ + m_OrientationTypeString.push_back("Angle"); + m_OrientationType.push_back(Angle); + m_Angle1.push_back(a); + m_Angle2.push_back(b); } //-------------------------------------------------------------------- //-------------------------------------------------------------------- -template +template void -clitk::AddRelativePositionConstraintToLabelImageFilter:: -SetAngle2(double a) { - SetOrientationType(Angle); - m_Angle2 = a; +clitk::AddRelativePositionConstraintToLabelImageFilter:: +AddAnglesInDeg(double a, double b) +{ + AddAnglesInRad(clitk::deg2rad(a), clitk::deg2rad(b)); } //-------------------------------------------------------------------- //-------------------------------------------------------------------- -template +template void -clitk::AddRelativePositionConstraintToLabelImageFilter:: -SetOrientationType(OrientationTypeEnumeration orientation) { - m_OrientationType = orientation; - switch (m_OrientationType) { - case LeftTo: m_Angle1 = clitk::deg2rad(0); m_Angle2 = clitk::deg2rad(0); break; - case RightTo: m_Angle1 = clitk::deg2rad(180); m_Angle2 = clitk::deg2rad(0); break; - case AntTo: m_Angle1 = clitk::deg2rad(90); m_Angle2 = clitk::deg2rad(0); break; - case PostTo: m_Angle1 = clitk::deg2rad(-90); m_Angle2 = clitk::deg2rad(0); break; - case InfTo: m_Angle1 = clitk::deg2rad(0); m_Angle2 = clitk::deg2rad(90); break; - case SupTo: m_Angle1 = clitk::deg2rad(0); m_Angle2 = clitk::deg2rad(-90); break; - case Angle: break; +clitk::AddRelativePositionConstraintToLabelImageFilter:: +AddOrientationType(OrientationTypeEnumeration orientation) +{ + m_OrientationType.push_back(orientation); + switch (orientation) { + case RightTo: + m_Angle1.push_back(clitk::deg2rad(0)); + m_Angle2.push_back(clitk::deg2rad(0)); + break; + case LeftTo: + m_Angle1.push_back(clitk::deg2rad(180)); + m_Angle2.push_back(clitk::deg2rad(0)); + break; + case AntTo: + m_Angle1.push_back(clitk::deg2rad(90)); + m_Angle2.push_back(clitk::deg2rad(0)); + break; + case PostTo: + m_Angle1.push_back(clitk::deg2rad(-90)); + m_Angle2.push_back(clitk::deg2rad(0)); + break; + case InfTo: + m_Angle1.push_back(clitk::deg2rad(0)); + m_Angle2.push_back(clitk::deg2rad(90)); + break; + case SupTo: + m_Angle1.push_back(clitk::deg2rad(0)); + m_Angle2.push_back(clitk::deg2rad(-90)); + break; + case Angle: break; } /* A1 A2 Left 0 0 @@ -177,86 +250,228 @@ SetOrientationType(OrientationTypeEnumeration orientation) { //-------------------------------------------------------------------- -template +template void -clitk::AddRelativePositionConstraintToLabelImageFilter:: -GenerateData() { +clitk::AddRelativePositionConstraintToLabelImageFilter:: +PrintOptions() +{ + DD((int)this->GetBackgroundValue()); + DD((int)this->GetObjectBackgroundValue()); + DDV(this->GetOrientationTypeString(), (uint)this->GetNumberOfAngles()); + DD(this->GetIntermediateSpacingFlag()); + DD(this->GetIntermediateSpacing()); + DD(this->GetFuzzyThreshold()); + DD(this->GetAutoCropFlag()); + DD(this->GetInverseOrientationFlag()); + DD(this->GetRemoveObjectFlag()); + DD(this->GetCombineWithOrFlag()); +} +//-------------------------------------------------------------------- + + +//-------------------------------------------------------------------- +template +void +clitk::AddRelativePositionConstraintToLabelImageFilter:: +GenerateData() +{ + if (GetNumberOfAngles() <1) { + clitkExceptionMacro("Add at least one orientation type"); + } + + if (GetVerboseOptionFlag()) { + for(int i=0; i(itk::ProcessObject::GetInput(0)); object = dynamic_cast(itk::ProcessObject::GetInput(1)); + static const unsigned int dim = ImageType::ImageDimension; + + // Step 2: object pad to input image -> we want to compute the + // relative position for each point belonging to the input image + // domain, so we have to extend (pad) the object image to fit the + // domain size + working_image = object; + if (!clitk::HaveSameSizeAndSpacing(input, working_image)) { + StartNewStep("Pad (resize) object to input size"); + + if (0) { // OLD VERSION (TO REMOVE) + StartNewStep("Pad object to image size"); + typename ImageType::Pointer output = ImageType::New(); + SizeType size; + for(unsigned int i=0; iGetLargestPossibleRegion().GetSize()[i]* + input->GetSpacing()[i])/(double)working_image->GetSpacing()[i]); + } + + // The index of the input is not necessarily zero, so we have to + // take it into account (not done) + RegionType region; + IndexType index = input->GetLargestPossibleRegion().GetIndex(); + region.SetSize(size); + for(unsigned int i=0; iSetLargestPossibleRegion(region); + output->SetRegions(region); + output->SetSpacing(working_image->GetSpacing()); + PointType origin = input->GetOrigin(); + for(unsigned int i=0; iGetSpacing()[i] + input->GetOrigin()[i]; + } + output->SetOrigin(origin); + // output->SetOrigin(input->GetOrigin()); + + output->Allocate(); + output->FillBuffer(m_BackgroundValue); + typename PasteFilterType::Pointer padFilter = PasteFilterType::New(); + // typename PasteFilterType::InputImageIndexType index; + for(unsigned int i=0; iGetSpacing()[i]/(double)working_image->GetSpacing()[i] + + lrint((working_image->GetOrigin()[i] - input->GetOrigin()[i])/working_image->GetSpacing()[i]); + } + padFilter->SetSourceImage(working_image); + padFilter->SetDestinationImage(output); + padFilter->SetDestinationIndex(index); + padFilter->SetSourceRegion(working_image->GetLargestPossibleRegion()); + padFilter->Update(); + working_image = padFilter->GetOutput(); + } + + // Resize object like input + working_image = clitk::ResizeImageLike(working_image, input, GetBackgroundValue()); + StopCurrentStep(working_image); + } //-------------------------------------------------------------------- //-------------------------------------------------------------------- - StartNewStep("Resample And Relative Position Map"); - - static const unsigned int dim = ImageType::ImageDimension; // Step 1 : resample - if (m_ResampleBeforeRelativePositionFilter) { + if (m_IntermediateSpacingFlag) { + StartNewStep("Resample object to intermediate spacing (" + toString(m_IntermediateSpacing) + ")"); typedef clitk::ResampleImageWithOptionsFilter ResampleFilterType; typename ResampleFilterType::Pointer resampleFilter = ResampleFilterType::New(); - resampleFilter->SetInput(object); + resampleFilter->SetInput(working_image); + resampleFilter->SetDefaultPixelValue(0); resampleFilter->SetOutputIsoSpacing(m_IntermediateSpacing); resampleFilter->SetGaussianFilteringEnabled(false); - // resampleFilter->SetVerboseOptions(true); + //resampleFilter->SetVerboseOptions(true); resampleFilter->Update(); working_image = resampleFilter->GetOutput(); + StopCurrentStep(working_image); } - else { - working_image = object; - } - // writeImage(working_image, "res.mhd"); - - // Step 2: object pad to input image -> we want to compute the - // relative position for each point belonging to the input image - // domain, so we have to extend (pad) the object image to fit the - // domain size - typename ImageType::Pointer output = ImageType::New(); - SizeType size; - for(unsigned int i=0; iGetLargestPossibleRegion().GetSize()[i]*input->GetSpacing()[i])/(double)working_image->GetSpacing()[i]); - } - RegionType region; - region.SetSize(size); - // output->SetLargestPossibleRegion(region); - output->SetRegions(region); - output->SetSpacing(working_image->GetSpacing()); - output->SetOrigin(input->GetOrigin()); - output->Allocate(); - output->FillBuffer(m_BackgroundValue); - typename PadFilterType::Pointer padFilter = PadFilterType::New(); - typename PadFilterType::InputImageIndexType index; - for(unsigned int i=0; iGetOrigin()[i] - input->GetOrigin()[i])/(double)m_IntermediateSpacing); - } - padFilter->SetSourceImage(working_image); - padFilter->SetDestinationImage(output); - padFilter->SetDestinationIndex(index); - padFilter->SetSourceRegion(working_image->GetLargestPossibleRegion()); - padFilter->Update(); - working_image = padFilter->GetOutput(); // Keep object image (with resampline and pad) object_resampled = working_image; - // writeImage(working_image, "pad.mhd"); // Step 3: compute rel pos in object + StartNewStep("Relative Position Map"); typedef itk::RelativePositionPropImageFilter RelPosFilterType; - typename RelPosFilterType::Pointer relPosFilter = RelPosFilterType::New(); - relPosFilter->SetInput(working_image); - relPosFilter->SetAlpha1(m_Angle1); // xy plane - relPosFilter->SetAlpha2(m_Angle2); - relPosFilter->SetK1(M_PI/2.0); // Opening parameter, default = pi/2 - relPosFilter->SetFast(true); - relPosFilter->SetRadius(1); // seems sufficient in this cas - // relPosFilter->SetVerboseProgress(true); - relPosFilter->Update(); - relPos = relPosFilter->GetOutput(); - this->template StopCurrentStep(relPos); + typename RelPosFilterType::Pointer relPosFilter; + + for(int i=0; iSetFast(GetFastFlag()); + relPosFilter->SetRadius(GetRadius()); + relPosFilter->SetInput(working_image); + relPosFilter->SetAlpha1(m_Angle1[i]); // xy plane + relPosFilter->SetAlpha2(m_Angle2[i]); + relPosFilter->SetK1(GetK1());// M_PI/2.0); // Opening parameter, default = pi/2 + // relPosFilter->SetVerboseProgress(true); + + relPosFilter->Update(); + relPos = relPosFilter->GetOutput(); + + if (GetNumberOfAngles() != 1) { + // Creation of the first m_FuzzyMap + if (i==0) { + m_FuzzyMap = clitk::NewImageLike(relPos, true); + m_FuzzyMap->FillBuffer(0.0); + } + + // Add to current fuzzy map + typedef itk::AddImageFilter AddImageFilter; + typename AddImageFilter::Pointer addFilter = AddImageFilter::New(); + addFilter->SetInput1(m_FuzzyMap); + addFilter->SetInput2(relPos); + addFilter->Update(); + m_FuzzyMap = addFilter->GetOutput(); + } + else m_FuzzyMap = relPos; + } + + // Divide by the number of relpos + if (GetNumberOfAngles() != 1) { +#if ITK_VERSION_MAJOR >= 4 + typedef itk::DivideImageFilter DivideFilter; + typename DivideFilter::Pointer divideFilter = DivideFilter::New(); + divideFilter->SetConstant2(GetNumberOfAngles()); +#else + typedef itk::DivideByConstantImageFilter DivideFilter; + typename DivideFilter::Pointer divideFilter = DivideFilter::New(); + divideFilter->SetConstant(GetNumberOfAngles()); +#endif + divideFilter->SetInput(m_FuzzyMap); + divideFilter->Update(); + m_FuzzyMap = divideFilter->GetOutput(); + } + + relPos = m_FuzzyMap; + StopCurrentStep(relPos); + if (GetFuzzyMapOnlyFlag()) { + // If the spacing is used, recompute fuzzy map with the input size/spacing + if (m_IntermediateSpacingFlag) { + StartNewStep("Resample FuzzyMap to come back to the same sampling than input"); + typedef clitk::ResampleImageWithOptionsFilter ResampleFilterType; + typename ResampleFilterType::Pointer resampleFilter = ResampleFilterType::New(); + resampleFilter->SetDefaultPixelValue(0.0); //Default fuzzymap value FIXME + resampleFilter->SetInput(relPos); + resampleFilter->SetOutputSpacing(input->GetSpacing()); + resampleFilter->SetGaussianFilteringEnabled(false); + resampleFilter->Update(); + relPos = m_FuzzyMap = resampleFilter->GetOutput(); + StopCurrentStep(relPos); + + // Need to put exactly the same size + if (relPos->GetLargestPossibleRegion() != input->GetLargestPossibleRegion()) { + StartNewStep("Pad to get the same size than input"); + typename FloatImageType::Pointer temp = FloatImageType::New(); + temp->CopyInformation(input); + temp->SetRegions(input->GetLargestPossibleRegion()); // Do not forget !! + temp->Allocate(); + temp->FillBuffer(0.0); // Default fuzzymap value FIXME + typename PasteFloatFilterType::Pointer padFilter2 = PasteFloatFilterType::New(); + padFilter2->SetSourceImage(relPos); + padFilter2->SetDestinationImage(temp); + padFilter2->SetDestinationIndex(input->GetLargestPossibleRegion().GetIndex()); + padFilter2->SetSourceRegion(relPos->GetLargestPossibleRegion()); + padFilter2->Update(); + relPos = padFilter2->GetOutput(); + StopCurrentStep(relPos); + m_FuzzyMap = relPos; + } + } + else { + StopCurrentStep(relPos); + } + return; + } //-------------------------------------------------------------------- //-------------------------------------------------------------------- - StartNewStep("Map Threshold And Post Processing"); - + StartNewStep("Map Threshold"); // Step 1: threshold typedef itk::BinaryThresholdImageFilter BinaryThresholdImageFilterType; typename BinaryThresholdImageFilterType::Pointer thresholdFilter = BinaryThresholdImageFilterType::New(); @@ -266,10 +481,14 @@ GenerateData() { thresholdFilter->SetLowerThreshold(m_FuzzyThreshold); thresholdFilter->Update(); working_image = thresholdFilter->GetOutput(); + StopCurrentStep(working_image); + //-------------------------------------------------------------------- + //-------------------------------------------------------------------- + StartNewStep("Post Processing: erosion with initial mask"); // Step 2 : erosion with initial mask to exclude pixels that were // inside the resampled version and outside the original mask - typedef itk::BinaryBallStructuringElement StructuringElementType; + typedef itk::BinaryBallStructuringElement StructuringElementType; StructuringElementType kernel; kernel.SetRadius(1); kernel.CreateStructuringElement(); @@ -281,9 +500,13 @@ GenerateData() { erodeFilter->SetErodeValue(m_BackgroundValue+1); erodeFilter->Update(); working_image = erodeFilter->GetOutput(); + StopCurrentStep(working_image); + //-------------------------------------------------------------------- + //-------------------------------------------------------------------- // Step 5: resample to initial spacing - if (m_ResampleBeforeRelativePositionFilter) { + if (m_IntermediateSpacingFlag) { + StartNewStep("Resample to come back to the same sampling than input"); typedef clitk::ResampleImageWithOptionsFilter ResampleFilterType; typename ResampleFilterType::Pointer resampleFilter = ResampleFilterType::New(); resampleFilter->SetDefaultPixelValue(m_BackgroundValue); @@ -293,28 +516,37 @@ GenerateData() { // resampleFilter->SetVerboseOptions(true); resampleFilter->Update(); working_image = resampleFilter->GetOutput(); + StopCurrentStep(working_image); } - // writeImage(working_image, "resinitial.mhd"); - + //-------------------------------------------------------------------- + //-------------------------------------------------------------------- // Pre Step 6: pad if not the same size : it can occur when downsample and upsample + //if (!HaveSameSizeAndSpacing(working_image, input)) { if (working_image->GetLargestPossibleRegion() != input->GetLargestPossibleRegion()) { + StartNewStep("Pad to get the same size than input"); typename ImageType::Pointer temp = ImageType::New(); temp->CopyInformation(input); temp->SetRegions(input->GetLargestPossibleRegion()); // Do not forget !! temp->Allocate(); - temp->FillBuffer(m_BackgroundValue); - typename PadFilterType::Pointer padFilter2 = PadFilterType::New(); // if yes : redo relpos + temp->FillBuffer(m_BackgroundValue); + typename PasteFilterType::Pointer padFilter2 = PasteFilterType::New(); padFilter2->SetSourceImage(working_image); padFilter2->SetDestinationImage(temp); padFilter2->SetDestinationIndex(input->GetLargestPossibleRegion().GetIndex()); padFilter2->SetSourceRegion(working_image->GetLargestPossibleRegion()); padFilter2->Update(); working_image = padFilter2->GetOutput(); + StopCurrentStep(working_image); + } + else { + //DD("[debug] Rel Pos : no padding after"); } - // writeImage(working_image, "pad2.mhd"); + //-------------------------------------------------------------------- + //-------------------------------------------------------------------- // Step 6: combine input+thresholded relpos + StartNewStep("Combine with initial input (boolean And)"); typedef clitk::BooleanOperatorLabelImageFilter BoolFilterType; typename BoolFilterType::Pointer combineFilter = BoolFilterType::New(); combineFilter->SetBackgroundValue(m_BackgroundValue); @@ -323,36 +555,69 @@ GenerateData() { combineFilter->SetForegroundValue(m_BackgroundValue+1); combineFilter->SetInput1(input); combineFilter->SetInput2(working_image); - combineFilter->SetOperationType(BoolFilterType::And); - combineFilter->InPlaceOn(); + if (GetInverseOrientationFlag()) + combineFilter->SetOperationType(BoolFilterType::AndNot); + else { + if (GetCombineWithOrFlag()) + combineFilter->SetOperationType(BoolFilterType::Or); + else + combineFilter->SetOperationType(BoolFilterType::And); + } + combineFilter->InPlaceOff(); // Do not modify initial input (!) combineFilter->Update(); working_image = combineFilter->GetOutput(); - - combineFilter = BoolFilterType::New(); - combineFilter->SetInput1(working_image); - combineFilter->SetInput2(object); - combineFilter->SetOperationType(BoolFilterType::AndNot); - combineFilter->InPlaceOn(); - combineFilter->Update(); - working_image = combineFilter->GetOutput(); - // writeImage(working_image, "combine.mhd"); + // Remove (if needed the object from the support) + if (GetRemoveObjectFlag()) { + combineFilter = BoolFilterType::New(); + combineFilter->SetInput1(working_image); + combineFilter->SetInput2(object); + combineFilter->SetOperationType(BoolFilterType::AndNot); + combineFilter->InPlaceOn(); + combineFilter->Update(); + working_image = combineFilter->GetOutput(); + } + // In the other case, we must *add* the initial object to keep it + // but not more than the initial support + else { + combineFilter = BoolFilterType::New(); + combineFilter->SetInput1(working_image); + combineFilter->SetInput2(object); + combineFilter->SetOperationType(BoolFilterType::Or); + combineFilter->InPlaceOn(); + combineFilter->Update(); + working_image = combineFilter->GetOutput(); // not needed because InPlaceOn ? + combineFilter = BoolFilterType::New(); + combineFilter->SetInput1(working_image); + combineFilter->SetInput2(input); + combineFilter->SetOperationType(BoolFilterType::And); + combineFilter->InPlaceOn(); + combineFilter->Update(); + working_image = combineFilter->GetOutput(); + } + + StopCurrentStep(working_image); + //-------------------------------------------------------------------- + //-------------------------------------------------------------------- // Step 7: autocrop - typedef clitk::AutoCropFilter CropFilterType; - typename CropFilterType::Pointer cropFilter = CropFilterType::New(); - cropFilter->SetInput(working_image); - cropFilter->ReleaseDataFlagOff(); - cropFilter->Update(); - working_image = cropFilter->GetOutput(); - this->template StopCurrentStep(working_image); + if (GetAutoCropFlag()) { + StartNewStep("Final AutoCrop"); + typedef clitk::AutoCropFilter CropFilterType; + typename CropFilterType::Pointer cropFilter = CropFilterType::New(); + cropFilter->SetInput(working_image); + cropFilter->ReleaseDataFlagOff(); + cropFilter->Update(); + working_image = cropFilter->GetOutput(); + StopCurrentStep(working_image); + } //-------------------------------------------------------------------- //-------------------------------------------------------------------- // Final Step -> set output this->SetNthOutput(0, working_image); - return; + // this->GraftOutput(working_image); } //--------------------------------------------------------------------