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

[Amardeep Singh]

Dear Baoyun

I am not sure, I understand what you mean by saying "Since you just 
recently fixed the GVF as advection term".
I had hacked the code, because I did not know of the flag 
"SetAutoGenerateSpeedAdvection(...)". In the end, I ended up using the 
same code as you.
You might want to have a look at the full posting:
http://www.nabble.com/Use-of-advection-image-to-guide-geodesic-level-set-td22083962.html#a22159804

Best regards
Amardeep

Baoyun Li wrote:
> Dear Amardeep:
>  
> Thank you for your help. Since you just recently fixed the GVF as 
> advection term, can you please send me a copy of your test code.
>  
> I saw your original code (GVFtest.cc) pasted on the web. But that is 
> not the working version.
>  
> Thanks.
>  
> Baoyun
>
> ------------------------------------------------------------------------
> *From:* Amardeep Singh <amar.singh at gmx.de>
> *To:* Baoyun Li <baoyun_li123 at yahoo.com>
> *Cc:* insight-users at itk.org
> *Sent:* Friday, March 6, 2009 4:36:35 AM
> *Subject:* Re: Segmenation Error--after set advection image by GVF flow
>
> Dear Baoyun
>
> Unfortunately, I have not yet encountered your error message, so I am 
> not sure, if I can help you.
> I would try to save the gradient vector field as a *.vtk file, not in 
> Analyze format, as the gradient vector flow output is a vector field. 
> Can Analyze deal with that? (I don't know.)
> Furthermore, I would suggest that you encapsulate every Update() call 
> in a try-catch block:
> try
>       {
>           ...->Update();
>       }
>       catch( itk::ExceptionObject & excep )
>       {
>         std::cerr << "Exception caught !" << std::endl;
>         std::cerr << excep << std::endl;
>       }
>
> This might give you more meaningful error messages.
>
> Best regards
> Amardeep
>
> Baoyun Li wrote:
> > Dear All:
> >  The error message is: Program *received* signal *SIGFPE*, 
> *Arithmetic exception*.
> >  It happens in itkImageHelp.h
> >  index[NLoop] = static_cast<IndexValueType>(offset / 
> offsetTable[NLoop]);
> >  inline static void ComputeIndexInner(const IndexType 
> &bufferedRegionIndex,
> >                                        OffsetValueType &offset,
> >                                        const OffsetValueType 
> offsetTable[],
> >                                        IndexType &index,
> >                                        const UniqueTypeBoolFalse& )
> >    {
> >    index[NLoop] = static_cast<IndexValueType>(offset / 
> offsetTable[NLoop]);
> >    offset = offset - (index[NLoop] * offsetTable[NLoop]);
> >    index[NLoop] = index[NLoop] + bufferedRegionIndex[NLoop];
> >    ImageHelper<NImageDimension, NLoop-1>::
> >      ComputeIndexInner(bufferedRegionIndex,
> >                        offset,
> >                        offsetTable,
> >                        index,
> >                        Concept::Detail::UniqueType_bool<(NLoop==1)>());
> >    }
> >  Does this mean some overflow, or my datatype has some problme.
> >  Thanks
> >  Baoyun
> >
> > 
> > ------------------------------------------------------------------------
> > *From:* Baoyun Li <baoyun_li123 at yahoo.com 
> <mailto:baoyun_li123 at yahoo.com>>
> > *To:* Amardeep Singh <amar.singh at gmx.de <mailto:amar.singh at gmx.de>>
> > *Cc:* insight-users at itk.org <mailto:insight-users at itk.org>
> > *Sent:* Thursday, March 5, 2009 4:51:41 PM
> > *Subject:* Re: Segmenation Error--after set advection image by GVF flow
> >
> > 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 
> <mailto:baoyun_li123 at yahoo.com>>
> > *To:* Amardeep Singh <amar.singh at gmx.de <mailto:amar.singh at gmx.de>>
> > *Cc:* insight-users at itk.org <mailto: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();
> >
> >
> >
> >
>