X-Git-Url: https://git.creatis.insa-lyon.fr/pubgit/?a=blobdiff_plain;f=itk%2FclitkAddRelativePositionConstraintToLabelImageFilter.txx;h=2b159bf9e9a98874f2a5f167ad41cdd5a27850e9;hb=302ee878f46ba0f96a4913b1715ae768bdb665b0;hp=66624a07a571df23564e758e09215d332320ddc0;hpb=0e9612d2827fc3e3ea10a9e8c2e365dd992692ae;p=clitk.git diff --git a/itk/clitkAddRelativePositionConstraintToLabelImageFilter.txx b/itk/clitkAddRelativePositionConstraintToLabelImageFilter.txx index 66624a0..2b159bf 100644 --- a/itk/clitkAddRelativePositionConstraintToLabelImageFilter.txx +++ b/itk/clitkAddRelativePositionConstraintToLabelImageFilter.txx @@ -36,11 +36,11 @@ #include "RelativePositionPropImageFilter.h" //-------------------------------------------------------------------- -template -clitk::AddRelativePositionConstraintToLabelImageFilter:: +template +clitk::AddRelativePositionConstraintToLabelImageFilter:: AddRelativePositionConstraintToLabelImageFilter(): clitk::FilterBase(), - itk::ImageToImageFilter() + itk::ImageToImageFilter() { this->SetNumberOfRequiredInputs(2); SetFuzzyThreshold(0.6); @@ -49,39 +49,17 @@ AddRelativePositionConstraintToLabelImageFilter(): SetOrientationType(LeftTo); ResampleBeforeRelativePositionFilterOn(); 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); + AutoCropOn(); } //-------------------------------------------------------------------- //-------------------------------------------------------------------- -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 +67,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 +79,12 @@ SetInputObject(const ImageType * image) { //-------------------------------------------------------------------- -template +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 +92,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,10 +109,11 @@ GenerateInputRequestedRegion() { //-------------------------------------------------------------------- -template +template void -clitk::AddRelativePositionConstraintToLabelImageFilter:: -SetAngle1(double a) { +clitk::AddRelativePositionConstraintToLabelImageFilter:: +SetAngle1(double a) +{ SetOrientationType(Angle); m_Angle1 = a; } @@ -139,10 +121,11 @@ SetAngle1(double a) { //-------------------------------------------------------------------- -template +template void -clitk::AddRelativePositionConstraintToLabelImageFilter:: -SetAngle2(double a) { +clitk::AddRelativePositionConstraintToLabelImageFilter:: +SetAngle2(double a) +{ SetOrientationType(Angle); m_Angle2 = a; } @@ -150,10 +133,11 @@ SetAngle2(double a) { //-------------------------------------------------------------------- -template +template void -clitk::AddRelativePositionConstraintToLabelImageFilter:: -SetOrientationType(OrientationTypeEnumeration orientation) { +clitk::AddRelativePositionConstraintToLabelImageFilter:: +SetOrientationType(OrientationTypeEnumeration orientation) +{ m_OrientationType = orientation; switch (m_OrientationType) { case LeftTo: m_Angle1 = clitk::deg2rad(0); m_Angle2 = clitk::deg2rad(0); break; @@ -177,19 +161,19 @@ SetOrientationType(OrientationTypeEnumeration orientation) { //-------------------------------------------------------------------- -template +template void -clitk::AddRelativePositionConstraintToLabelImageFilter:: -GenerateData() { +clitk::AddRelativePositionConstraintToLabelImageFilter:: +GenerateData() +{ // Get input pointer input = dynamic_cast(itk::ProcessObject::GetInput(0)); object = dynamic_cast(itk::ProcessObject::GetInput(1)); //-------------------------------------------------------------------- //-------------------------------------------------------------------- - StartNewStep("Resample And Relative Position Map"); - static const unsigned int dim = ImageType::ImageDimension; + StartNewStep("Initial resample"); // Step 1 : resample if (m_ResampleBeforeRelativePositionFilter) { typedef clitk::ResampleImageWithOptionsFilter ResampleFilterType; @@ -197,37 +181,55 @@ GenerateData() { resampleFilter->SetInput(object); resampleFilter->SetOutputIsoSpacing(m_IntermediateSpacing); resampleFilter->SetGaussianFilteringEnabled(false); - // resampleFilter->SetVerboseOptions(true); + // resampleFilter->SetVerboseOptions(true); resampleFilter->Update(); working_image = resampleFilter->GetOutput(); } else { working_image = object; - } - // writeImage(working_image, "res.mhd"); + } + StopCurrentStep(working_image); // 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 if (!clitk::HaveSameSizeAndSpacing(input, working_image)) { + 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()); - output->SetOrigin(input->GetOrigin()); + 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 PadFilterType::Pointer padFilter = PadFilterType::New(); - typename PadFilterType::InputImageIndexType index; + // typename PadFilterType::InputImageIndexType index; for(unsigned int i=0; iGetOrigin()[i] - input->GetOrigin()[i])/(double)m_IntermediateSpacing); + index[i] = -index[i]*input->GetSpacing()[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); @@ -235,16 +237,17 @@ GenerateData() { padFilter->SetSourceRegion(working_image->GetLargestPossibleRegion()); padFilter->Update(); working_image = padFilter->GetOutput(); + StopCurrentStep(working_image); } else { - DD("[debug] RelPos : same size and spacing : no padding"); + // DD("[debug] RelPos : same size and spacing : no padding"); } - // Keep object image (with resampline and pad) object_resampled = working_image; - // writeImage(working_image, "pad.mhd"); + // StopCurrentStep(working_image); // 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); @@ -256,12 +259,11 @@ GenerateData() { // relPosFilter->SetVerboseProgress(true); relPosFilter->Update(); relPos = relPosFilter->GetOutput(); - this->template StopCurrentStep(relPos); + StopCurrentStep(relPos); //-------------------------------------------------------------------- //-------------------------------------------------------------------- - StartNewStep("Map Threshold And Post Processing"); - + StartNewStep("Map Threshold"); // Step 1: threshold typedef itk::BinaryThresholdImageFilter BinaryThresholdImageFilterType; typename BinaryThresholdImageFilterType::Pointer thresholdFilter = BinaryThresholdImageFilterType::New(); @@ -271,10 +273,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(); @@ -286,9 +292,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) { + StartNewStep("Resample to get the same sampling than input"); typedef clitk::ResampleImageWithOptionsFilter ResampleFilterType; typename ResampleFilterType::Pointer resampleFilter = ResampleFilterType::New(); resampleFilter->SetDefaultPixelValue(m_BackgroundValue); @@ -298,12 +308,14 @@ 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 (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 !! @@ -316,15 +328,20 @@ GenerateData() { padFilter2->SetSourceRegion(working_image->GetLargestPossibleRegion()); padFilter2->Update(); working_image = padFilter2->GetOutput(); + StopCurrentStep(working_image); } else { - DD("[debug] Rel Pos : no padding after"); + //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(); + writeImage(input, "i.mhd"); + writeImage(working_image, "w.mhd"); combineFilter->SetBackgroundValue(m_BackgroundValue); combineFilter->SetBackgroundValue1(m_BackgroundValue); combineFilter->SetBackgroundValue2(m_BackgroundValue); @@ -335,6 +352,7 @@ GenerateData() { combineFilter->InPlaceOn(); combineFilter->Update(); working_image = combineFilter->GetOutput(); + // writeImage(working_image, "res.mhd"); combineFilter = BoolFilterType::New(); combineFilter->SetInput1(working_image); @@ -344,16 +362,21 @@ GenerateData() { combineFilter->Update(); working_image = combineFilter->GetOutput(); - // writeImage(working_image, "combine.mhd"); + 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 (GetAutoCrop()) { + 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); + } //-------------------------------------------------------------------- //--------------------------------------------------------------------