[Insight-users] Segmenation Error--after set advection image by GVF flow

[Baoyun Li]

Dear Amardeep:


After I called

 geodesicActiveContour->SetAdvectionImage(castflowfilter->GetOutput()); 

Do I need to call GenerateAdvectioniImage() since Atumatic Generate is set off.

After I add : geodesicActiveContour->GenerateAdvectionImage();

I still got error, please give me some suggestions.

Baoyun


________________________________
From: Baoyun Li <baoyun_li123 at yahoo.com>
To: Amardeep Singh <amar.singh at gmx.de>
Cc: insight-users at itk.org
Sent: Thursday, March 5, 2009 3:25:12 PM
Subject: Segmenation Error--after set advection image by GVF flow


Dear Amardeep and All:

Now I can build the program after seting advection image by GVF flow. However, when I run the program, the program shows segmeation error when geodesicActiveContour->Update();

The program stop at: itkSegmentationLevelSetFunction_txx  at  return ( m_AdvectionImage->GetPixel(idx) );


template <class TImageType, class TFeatureImageType>
typename SegmentationLevelSetFunction<TImageType, TFeatureImageType>::VectorType
SegmentationLevelSetFunction<TImageType, TFeatureImageType>
::AdvectionField(const NeighborhoodType &neighborhood,
                 const FloatOffsetType &offset, GlobalDataStruct *)  const
{
  IndexType idx = neighborhood.GetIndex();
  ContinuousIndexType cdx;
  for (unsigned i = 0; i < ImageDimension; ++i)
    {
    cdx[i] = static_cast<double>(idx[i]) - offset[i];
    }
  if ( m_VectorInterpolator->IsInsideBuffer(cdx) )
    {
    return ( m_VectorCast(m_VectorInterpolator->EvaluateAtContinuousIndex(cdx)));
    }
  //Just return the default else
    return ( m_AdvectionImage->GetPixel(idx) );
  
}

Clearly, something is still wrong with SetAdvectionImage. Below is related code.

Can you please help me figure out which setting is wrong? 

Thanks

Baoyun

 
  typedef itk::CovariantVector<double, Dim> myGradientType;
  typedef itk::Image<myGradientType, Dim>   myGradientImageType;
  typedef itk::GradientVectorFlowImageFilter<myGradientImageType, myGradientImageType>
                                              myGVFFilterType;
  typedef itk::GradientImageFilter<DoublelImageType, double,double>
                                              myGFilterType;
  typedef itk::GradientToMagnitudeImageFilter<myGradientImageType, DoublelImageType>
                                              myGToMFilterType;

  typedef itk::GradientRecursiveGaussianImageFilter<
                                            DoublelImageType,
                                            myGradientImageType
                                                  >  myFilterType;
  typedef typename GeodesicActiveContourFilterType::VectorImageType   AdvectionImageType; 
  typedef itk::CastImageFilter<myGradientImageType,AdvectionImageType> 
   CastFlowFilterType;
  typedef itk::LaplacianImageFilter<DoublelImageType, DoublelImageType> myLaplacianFilterType;
  typename myFilterType::Pointer filter = myFilterType::New();
  typename myGFilterType::Pointer gfilter = myGFilterType::New();
  typename myGToMFilterType::Pointer gtomfilter = myGToMFilterType::New();
  typename myLaplacianFilterType::Pointer m_LFilter = myLaplacianFilterType::New();
  typename myGVFFilterType::Pointer m_GVFFilter = myGVFFilterType::New();
  typename CastFlowFilterType::Pointer castflowfilter = CastFlowFilterType::New(); 
  typedef itk::ImageFileWriter< myGradientImageType >  WriterType2;
  typename WriterType2::Pointer writer21=WriterType2::New(); 
  

  typename WriterType1::Pointer writer4=WriterType1::New();
  writer4->SetInput(caster11->GetOutput());
  writer4->SetFileName("../data/intial_contour.hdr");
  writer4->Update();

  sigmoid->SetOutputMinimum(  0.0  );
  sigmoid->SetOutputMaximum(  1.0  );
  caster12->SetInput(inputfilter);
  smoothing->SetInput(caster12->GetOutput() );
  gradientMagnitude->SetInput( smoothing->GetOutput() );
  sigmoid->SetInput( gradientMagnitude->GetOutput() );
  caster31->SetInput(inputfilter);
  filter->SetInput( caster31->GetOutput() ); //caster short to double
  filter->SetSigma( 1.0);
  filter->Update();
  gtomfilter->SetInput(filter->GetOutput());
  gtomfilter->Update();
  gfilter->SetInput(gtomfilter->GetOutput());
  gfilter->Update();
  m_GVFFilter->SetInput(gfilter->GetOutput());
  m_GVFFilter->SetLaplacianFilter(m_LFilter);
  m_GVFFilter->SetNoiseLevel(500);
  m_GVFFilter->SetTimeStep(0.05);
  m_GVFFilter->SetIterationNum(5);
  m_GVFFilter->Update();
  writer21->SetInput(m_GVFFilter->GetOutput());
  writer21->SetFileName("../data/GVF_flow.hdr");
  writer21->Update();
//to get gradient manitude of GVF low
 
  geodesicActiveContour->SetAutoGenerateSpeedAdvection(false);
  //typename myGradientImageType::Pointer m_GVFField; 
  //m_GVFField = m_GVFFilter->GetOutput();
  castflowfilter->SetInput(m_GVFFilter->GetOutput()); 
  geodesicActiveContour->SetAdvectionImage(castflowfilter->GetOutput()); 
 // geodesicActiveContour->GenerateSpeedImage(); 
  //geodesicActiveContour->SetAdvectionImage(m_GVFFilter->GetOutput()); 
 
  
   
   
  
  geodesicActiveContour->SetInput( caster11->GetOutput() );
  geodesicActiveContour->SetFeatureImage( sigmoid->GetOutput() );
  geodesicActiveContour->GenerateSpeedImage(); 
  geodesicActiveContour->Update();
  thresholder->SetInput( geodesicActiveContour->GetOutput() );
  thresholder->Update();
  outputfilter=thresholder->GetOutput();


      
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://www.itk.org/pipermail/insight-users/attachments/20090305/2d08b06f/attachment-0001.htm>