<html xmlns:v="urn:schemas-microsoft-com:vml" xmlns:o="urn:schemas-microsoft-com:office:office" xmlns:w="urn:schemas-microsoft-com:office:word" xmlns:m="http://schemas.microsoft.com/office/2004/12/omml" xmlns="http://www.w3.org/TR/REC-html40"><head><META HTTP-EQUIV="Content-Type" CONTENT="text/html; charset=us-ascii"><meta name=Generator content="Microsoft Word 12 (filtered medium)"><style><!--
/* Font Definitions */
@font-face
        {font-family:"Cambria Math";
        panose-1:2 4 5 3 5 4 6 3 2 4;}
@font-face
        {font-family:Calibri;
        panose-1:2 15 5 2 2 2 4 3 2 4;}
@font-face
        {font-family:Tahoma;
        panose-1:2 11 6 4 3 5 4 4 2 4;}
/* Style Definitions */
p.MsoNormal, li.MsoNormal, div.MsoNormal
        {margin:0cm;
        margin-bottom:.0001pt;
        font-size:12.0pt;
        font-family:"Times New Roman","serif";}
a:link, span.MsoHyperlink
        {mso-style-priority:99;
        color:blue;
        text-decoration:underline;}
a:visited, span.MsoHyperlinkFollowed
        {mso-style-priority:99;
        color:purple;
        text-decoration:underline;}
span.EmailStyle17
        {mso-style-type:personal-reply;
        font-family:"Calibri","sans-serif";
        color:#1F497D;}
.MsoChpDefault
        {mso-style-type:export-only;}
.MsoPapDefault
        {mso-style-type:export-only;
        mso-margin-top-alt:auto;
        mso-margin-bottom-alt:auto;}
@page WordSection1
        {size:612.0pt 792.0pt;
        margin:70.85pt 70.85pt 2.0cm 70.85pt;}
div.WordSection1
        {page:WordSection1;}
--></style><!--[if gte mso 9]><xml>
<o:shapedefaults v:ext="edit" spidmax="1026" />
</xml><![endif]--><!--[if gte mso 9]><xml>
<o:shapelayout v:ext="edit">
<o:idmap v:ext="edit" data="1" />
</o:shapelayout></xml><![endif]--></head><body lang=DE-AT link=blue vlink=purple><div class=WordSection1><p class=MsoNormal><span style='font-size:11.0pt;font-family:"Calibri","sans-serif";color:#1F497D'>Hi Sundar,<o:p></o:p></span></p><p class=MsoNormal><span style='font-size:11.0pt;font-family:"Calibri","sans-serif";color:#1F497D'><o:p> </o:p></span></p><p class=MsoNormal><span lang=EN-GB style='font-size:11.0pt;font-family:"Calibri","sans-serif";color:#1F497D'>I use the eigenvalue analysis too. So i cross verified with my code whether you were missing something. My code is:<o:p></o:p></span></p><p class=MsoNormal><span lang=EN-GB style='font-size:11.0pt;font-family:"Calibri","sans-serif";color:#1F497D'><o:p> </o:p></span></p><p class=MsoNormal style='text-autospace:none'><span lang=EN-GB style='font-size:11.0pt;font-family:"Calibri","sans-serif";color:#1F497D'> typedef itk::SymmetricEigenAnalysis < HessianRecursiveGaussianFilterType::OutputImageType::PixelType, VectorPixelType, EV_PixelType > EigAnalysisType;<o:p></o:p></span></p><p class=MsoNormal style='text-autospace:none'><span lang=EN-GB style='font-size:11.0pt;font-family:"Calibri","sans-serif";color:#1F497D'> //set eigen analysis<o:p></o:p></span></p><p class=MsoNormal style='text-autospace:none'><span lang=EN-GB style='font-size:11.0pt;font-family:"Calibri","sans-serif";color:#1F497D'> EigAnalysisType eig;<o:p></o:p></span></p><p class=MsoNormal style='text-autospace:none'><span lang=EN-GB style='font-size:11.0pt;font-family:"Calibri","sans-serif";color:#1F497D'> <span style='background:red;mso-highlight:red'>eig.SetDimension( 3 );</span><o:p></o:p></span></p><p class=MsoNormal><span lang=EN-GB style='font-size:11.0pt;font-family:"Calibri","sans-serif";color:#1F497D'> eig.SetOrderEigenMagnitudes( true );<o:p></o:p></span></p><p class=MsoNormal><span lang=EN-GB style='font-size:11.0pt;font-family:"Calibri","sans-serif";color:#1F497D'><o:p> </o:p></span></p><p class=MsoNormal><span lang=EN-GB style='font-size:11.0pt;font-family:"Calibri","sans-serif";color:#1F497D'>I found that the input matrix dimension was not set. Please check whether this is causing the error.<o:p></o:p></span></p><p class=MsoNormal><span lang=EN-GB style='font-size:11.0pt;font-family:"Calibri","sans-serif";color:#1F497D'>All the best.<o:p></o:p></span></p><p class=MsoNormal><span lang=EN-GB style='font-size:11.0pt;font-family:"Calibri","sans-serif";color:#1F497D'>Kana<o:p></o:p></span></p><p class=MsoNormal><span lang=EN-GB style='font-size:11.0pt;font-family:"Calibri","sans-serif";color:#1F497D'><o:p> </o:p></span></p><div style='border:none;border-top:solid #B5C4DF 1.0pt;padding:3.0pt 0cm 0cm 0cm'><p class=MsoNormal><b><span lang=EN-US style='font-size:10.0pt;font-family:"Tahoma","sans-serif"'>From:</span></b><span lang=EN-US style='font-size:10.0pt;font-family:"Tahoma","sans-serif"'> insight-users-bounces@itk.org [mailto:insight-users-bounces@itk.org] <b>On Behalf Of </b>Sundaresan R<br><b>Sent:</b> 22 February 2011 22:58<br><b>To:</b> insight-users@itk.org<br><b>Subject:</b> [Insight-users] EigenVector Computation<o:p></o:p></span></p></div><p class=MsoNormal><span lang=EN-GB><o:p> </o:p></span></p><p class=MsoNormal style='mso-margin-top-alt:auto;mso-margin-bottom-alt:auto'>Hi<br><br>I'm computing the Hessians and then use SymmetricEigenAnalysis to compute the eigenVectors as follows:<br>I get a segmentation fault at the call to ComputeEigenValuesAndVectors method. Can anyone help me fix this issue?<br><br>Thanks,<br>Sundar<br><br>**************************************************************************************<br> <br>#include "itkHessianRecursiveGaussianImageFilter.h" <br>#include "itkSymmetricEigenAnalysis.h"<br>#include "itkSymmetricEigenAnalysisImageFilter.h" <br>#include "itkSymmetricSecondRankTensor.h" <br>#include "itkCastImageFilter.h" <br>#include "itkImage.h" <br>#include "itkImageFileReader.h" <br>#include "itkImageFileWriter.h" <br>#include <iostream><br><br> <br> const unsigned int Dimension = 3; <br> typedef double InputPixelType; <br> typedef double OutputPixelType; <br> typedef itk::Image< InputPixelType, Dimension > InputImageType; <br> typedef itk::Image< OutputPixelType, Dimension > OutputImageType; <br> typedef itk::ImageFileReader< InputImageType > ReaderType; <br> typedef itk::ImageFileWriter< OutputImageType > WriterType; <br> typedef itk::HessianRecursiveGaussianImageFilter< InputImageType > HessianFilterType; <br> typedef itk::Vector< InputPixelType, Dimension > VectorType;<br> typedef itk::Image< VectorType, Dimension > VectorImageType; <br> typedef itk::Matrix< InputPixelType, Dimension, Dimension > MatrixType;<br> typedef itk::Image< MatrixType, Dimension > MatrixImageType;<br> typedef itk::ImageRegionConstIterator< MatrixImageType > MatrixIteratorType;<br> typedef itk::ImageRegionIterator< MatrixImageType > OutputIteratorType;<br> typedef itk::SymmetricEigenAnalysis< MatrixType, VectorType, MatrixType > EigenVectorAnalysisType;<br> typedef MatrixImageType::Pointer MatrixImagePointer;<br> typedef MatrixImageType::RegionType RegionType;<br> typedef MatrixImageType::PointType PointType;<br> typedef MatrixImageType::SpacingType SpacingType;<br> typedef itk::ImageRegionIteratorWithIndex<HessianFilterType::OutputImageType> HessianIteratorType;<br> typedef itk::ImageFileWriter< MatrixImageType > MatrixWriterType;<br> <br>int main( int argc, char *argv[] ) <br>{ <br> if( argc < 3 ) <br> { <br> std::cerr << "Usage: ./" << argv[0] << " inputImage sigma" << std::endl; <br> } <br> <br> ReaderType::Pointer reader = ReaderType::New(); <br> reader->SetFileName( argv[1] ); <br> <br> HessianFilterType::Pointer hessianFilter = HessianFilterType::New(); <br> hessianFilter->SetInput( reader->GetOutput() ); <br> if( argc >= 3 ) <br> hessianFilter->SetSigma( static_cast< double >( atof(argv[2])) ); <br> std::cout << "Calculating Hessian matrix" << std::endl; <br> hessianFilter->Update(); <br> std::cout << "Done" << std::endl; <br> <br> EigenVectorAnalysisType evap; <br> HessianFilterType::OutputImageType::Pointer input = hessianFilter->GetOutput();<br> RegionType region = input->GetLargestPossibleRegion();<br> SpacingType spacing = input->GetSpacing();<br> PointType origin = input->GetOrigin();<br><br> MatrixType eigenMatrix, myMatrix;<br> eigenMatrix[0][0] = 0; eigenMatrix[0][1]=0;eigenMatrix[0][2]=0;<br> eigenMatrix[1][0] = 0; eigenMatrix[1][1]=0;eigenMatrix[1][2]=0;<br> eigenMatrix[2][0] = 0; eigenMatrix[2][1]=0;eigenMatrix[2][2]=0;<br><br> // Allocate the output image<br> MatrixImagePointer m_Output = MatrixImageType::New();<br> m_Output->SetRegions( region );<br> m_Output->SetOrigin( origin );<br> m_Output->SetSpacing( spacing );<br> m_Output->Allocate();<br> m_Output->FillBuffer( eigenMatrix );<br><br> VectorType eigenVal;<br> eigenVal[0]=eigenVal[1]=eigenVal[2]=0;<br> HessianIteratorType mIt( input, region );<br> OutputIteratorType oIt( m_Output, region );<br> itk::SymmetricSecondRankTensor<double, 3> myTensor;<br> mIt.GoToBegin();<br> oIt.GoToBegin();<br> while( !mIt.IsAtEnd() )<br> {<br> // Compute eigen values and eigen matrix<br> evap.SetOrderEigenValues( 1 );<br> myTensor = mIt.Get();<br> myMatrix[0][0] = myTensor[0];<br> myMatrix[0][1] = myTensor[1];<br> myMatrix[0][2] = myTensor[2];<br> myMatrix[1][0] = myTensor[1];<br> myMatrix[1][1] = myTensor[3];<br> myMatrix[1][2] = myTensor[4];<br> myMatrix[2][0] = myTensor[2];<br> myMatrix[2][1] = myTensor[4];<br> myMatrix[2][2] = myTensor[5];<br> std::cout << "myTensor " << myTensor << std::endl;<br> std::cout << "myMatrix " << myMatrix << std::endl; <br> std::cout << "eigenVal " << eigenVal << std::endl;<br> std::cout << "eigenMatrix " << eigenMatrix << std::endl; <br> evap.ComputeEigenValuesAndVectors( myMatrix, eigenVal, eigenMatrix );<br> std::cout << " After ComputeEigenValuesAndVectors" << std::endl; <br> oIt.Set( eigenMatrix );<br> ++oIt;<br> ++mIt;<br> }<br> MatrixWriterType::Pointer writer = MatrixWriterType::New(); <br> writer->SetInput( m_Output );<br> writer->SetFileName( "EigenMatrix.mhd" ); <br> try<br> {<br> writer->Update();<br> }<br> catch ( itk::ExceptionObject e )<br> {<br> std::cerr << "Error: " << e << std::endl;<br> }<br> return EXIT_SUCCESS; <br>} <br>***************************************************************************************************<o:p></o:p></p></div></body></html>