失眠网 > ITK—hessian矩阵计算实现及应用

ITK—hessian矩阵计算实现及应用

时间：2019-04-04 09:56:03

文章目录

ITK—hessian矩阵计算实现及应用一、hessian矩阵定义二、hessian矩阵特征值1.特征值2.特征值的应用3.基于itk实现hessian矩阵计算总结

一、hessian矩阵定义

表示一个点各个方向的曲率：

三维空间中的每个点的hessian矩阵为：

二、hessian矩阵特征值

1.特征值

hessian矩阵的三个特征值为lamda1，lamda2，lamda3

特征值越大，在其对应的特征向量方向上，该点的空间曲率越大，将特征值默认按绝对值从小到大排列。

2.特征值的应用

根据特征值间的差异可以判断当前点的空间结构

1、lamda1和lamda2趋近于0，lamda3远小于lamda1和lamda2时，空间结构为薄面体

2、lamda1和lamda2渐渐减小趋近于lamda3，且三个特征值都远小于0时，薄面体趋近于为节点结构

3、lamda1保持不变，lamda2渐渐趋近于lamda3，薄面体结构趋近于管状结构

基于上述特征值的特性，可以用作医学图像上血管增强。

3.基于itk实现hessian矩阵计算

/计算hessian矩阵typedef itk::HessianRecursiveGaussianImageFilter<ITKImageType > HessianFilterType;typename HessianFilterType::Pointer hessianFilter = HessianFilterType::New();hessianFilter->SetInput(input);hessianFilter->SetSigma(sigma);hessianFilter->Update();typedef itk::Image< itk::SymmetricSecondRankTensor< double, 3 >, 3 > HessianImageType;typedef itk::FixedArray< double, 3 >EigenValueArrayType;typedef itk::Image< EigenValueArrayType, 3 >EigenValueImageType;//hesian矩阵的特征值分析typedef itk::SymmetricEigenAnalysisImageFilter<HessianImageType, EigenValueImageType >EigenAnalysisFilterType;typename EigenAnalysisFilterType::Pointer m_SymmetricEigenValueFilter = EigenAnalysisFilterType::New();m_SymmetricEigenValueFilter->SetDimension(3);m_SymmetricEigenValueFilter->OrderEigenValuesBy(EigenAnalysisFilterType::FunctorType::OrderByMagnitude);m_SymmetricEigenValueFilter->SetInput(hessianFilter->GetOutput());m_SymmetricEigenValueFilter->Update();typedef typename EigenAnalysisFilterType::OutputImageType EigenValueOutputImageType;const typename EigenValueOutputImageType::ConstPointer eigenImage = m_SymmetricEigenValueFilter->GetOutput();typedef itk::ImageRegionConstIterator< EigenValueOutputImageType > TraceItType;TraceItType traceIt(eigenImage, eigenImage->GetRequestedRegion());traceIt.GoToBegin();//遍历每个体素，可以依据特征值的情况做具体应用，比如上述各种特定结构增强，具体可以参考血管增强实例，或者其他论文理论。vector<double> Degree;vector<IndexType> Points;while (!traceIt.IsAtEnd()){EigenValueArrayType eigenValue;IndexType index;eigenValue = traceIt.Get();index = traceIt.GetIndex();double s = 0;按绝对值从小到大排列特征值for (int i = 0; i < 3; i++){for (int j = 0; j < 2 - i; j++){if (abs(eigenValue[j]) > abs(eigenValue[j + 1])){double temp = eigenValue[j];eigenValue[j] = eigenValue[j + 1];eigenValue[j + 1] = temp;}}}//具体应用..........}