c++image-processingitk

ITK Fast Marching output


I'm using ITK to do some preprocessing and I wanted to test something with the Fast Marching filter and the Geodesic Active Contour filter. I'm following the algorithm described in the ITK software guide, section 9.3.3. However, I'm not getting the expected results. I'm working with a 3D image.

Here is my code:

AnisotropicDiffusionFilter::Pointer anisotropic_filter = AnisotropicDiffusionFilter::New();
anisotropic_filter->SetInput(itk_image_in);
anisotropic_filter->SetTimeStep(0.0625);
anisotropic_filter->SetNumberOfIterations(5);
anisotropic_filter->SetConductanceParameter(3.0); 
anisotropic_filter->Update();
 
GradientFilter::Pointer gradient_filter = GradientFilter::New();
gradient_filter->SetInput(anisotropic_filter->GetOutput());
gradient_filter->SetSigma(0.5);
gradient_filter->Update();
 
SigmoidFilter::Pointer sigmoid_filter = SigmoidFilter::New();
sigmoid_filter->SetInput(gradient_filter->GetOutput());
sigmoid_filter->SetOutputMinimum(0.0);
sigmoid_filter->SetOutputMaximum(1.0);
sigmoid_filter->SetAlpha(-1.5);
sigmoid_filter->SetBeta(4.0);
sigmoid_filter->Update();
 
FastMarchingFilter::Pointer fast_marching = FastMarchingFilter::New();
NodeContainer::Pointer seeds = NodeContainer::New();
Node node;
const double seedValue = -50.0;
node.SetValue(seedValue);
seeds->Initialize();
 
vector<GeoVec3s>::iterator it = m_clicks_.begin();
int i=0;
for(; it != m_clicks_.end(); it++)
{
    itkIndex index;
    index[0] = (*it)[0];
    index[1] = (*it)[1];
    index[2] = (*it)[2];
    node.SetIndex(index);
    seeds->InsertElement(i++, node);
}
fast_marching->SetTrialPoints(seeds);
fast_marching->SetSpeedConstant(1.0);
fast_marching->SetStoppingValue(100);
//fast_marching->SetInput(sigmoid_filter->GetOutput());
fast_marching->SetOutputSize(sigmoid_filter->GetOutput()->GetBufferedRegion().GetSize());
fast_marching->Update();
 
GeodesicFilter::Pointer geodesic_filter = GeodesicFilter::New();
geodesic_filter->SetInput(fast_marching->GetOutput());
geodesic_filter->SetFeatureImage(sigmoid_filter->GetOutput());
 
geodesic_filter->SetPropagationScaling(0.5);
geodesic_filter->SetCurvatureScaling(5.0);
geodesic_filter->SetAdvectionScaling(1.0);
geodesic_filter->SetMaximumRMSError( 0.02 );
geodesic_filter->Update();
 
BinaryThresholdFilter::Pointer thresholder = BinaryThresholdFilter::New();
thresholder->SetLowerThreshold(-1000);
thresholder->SetUpperThreshold(0);
thresholder->SetOutsideValue(0);
thresholder->SetInsideValue(255);
thresholder->SetInput( geodesic_filter->GetOutput() );

I'm using metrics described in this paper which goal is the same as mine.

I have a few questions:

  1. The fast marching filter should output a distance map right? Instead, when I output my volume to a series of png (between values 0 and 4095) I have a binary image (pixels are either 0 or 4095). I think I should get a greyscale volume indicating the time needed for each pixel to be attained from the seeds.
  2. Following the procedure described by Suzuki I succeeded to make the algorithm work more or less however I changed the values of the parameters of the geodesic filter. I don't remember the exact values but it wasn't close to those described in the paper. As we are working with the sigmoid input which is normalized between 0 and 1, what is happening?
  3. Should I rather use a constant speed function for the fast marching filter or the sigmoid image? When should either method be preferred?
  4. I'm using a re-scaler to output my float images (output from the filters). Could this be the reason for the inconsistencies I'm seeing?
  5. Any advice on what I could be doing wrong?

Thanks.


Solution

  • Ok so I found my problem. The Fast Marching filter does output a time crossing map (distance map) but as I specified a stopping value in the algorithm all the pixels that weren't visited had a high value (1.7e+38 as it is half the max value of the type used for the output image which were float in my case 3.4e+38). So it squeezed all my image dynamic when I used the rescale filter and the result was an binary image. I think better results are achieved with a sigmoid image as input for the fast marching filter. Thanks @nav for the advice.