当前位置： 首页 > news >正文

3D医学影像开发＜五＞：利用Hessian矩阵增强影像并设置固定阈值进行血管模型自动分割

news 2025/7/3 21:41:54

1、医学影像中对血管进行自动分割的意义

  在医学影像中对血管进行自动分割具有重要的临床和研究意义，主要体现在以下几个方面：

   1）自动分割血管可以帮助医生更准确地识别和定位血管病变，例如：动脉瘤、血管狭窄或阻塞、动静脉畸形、血管肿瘤，这样可以提高早期诊断的准确性和效率。

   2）在进行介入治疗或外科手术前，需要对病灶部位的血管结构有清晰了解，自动分割可以生成三维血管模型，帮助医生规划最佳的手术路径或介入方案。在术中导航系统中，分割结果还能用于实时定位。

   3）对于慢性血管疾病患者，如动脉粥样硬化、糖尿病视网膜病变等，自动分割可用于定期对血管变化进行量化分析（如血管直径、形态、密度变化），评估疾病进展和治疗效果。

   4）传统血管分割依赖人工标注，不仅耗时长，而且易受主观影响。自动分割可以显著提高工作效率。
在这里插入图片描述

2、Hessian矩阵算法

   ITK（Insight Segmentation and Registration Toolkit）中，Hessian 矩阵是一种重要的图像处理工具，广泛用于血管等细长结构的增强与分割，主要用于血管分割（CTA/MRA/眼底图像）、神经纤维束提取、细长结构（如气管、导管、神经鞘）增强等。以下是 ITK 中 Hessian 矩阵算法的详细介绍：

1.1 什么是 Hessian 矩阵？

  Hessian 矩阵是一个关于图像灰度的二阶偏导数矩阵，反映了图像在不同方向上的曲率信息。对于三维图像，Hessian 是一个 3×3 矩阵，其形式为：
在这里插入图片描述

1.2 ITK 中相关的类和流程

  在 ITK 中，常用如下几个类来计算并利用 Hessian 矩阵：

1） itk::HessianRecursiveGaussianImageFilter

  • 用途：用于计算图像的 Hessian 矩阵。
  • 原理：先对图像进行高斯平滑，然后计算二阶导数，避免噪声影响。
  • 输入：原始图像。
  • 输出：每个像素处的 Hessian 矩阵（封装为 symmetric tensor image）。

示例：

using HessianFilterType = itk::HessianRecursiveGaussianImageFilter<ImageType>;
HessianFilterType::Pointer hessianFilter = HessianFilterType::New();
hessianFilter->SetInput(inputImage);
hessianFilter->SetSigma(1.0);  // 控制平滑程度

2）itk::HessianToObjectnessMeasureImageFilter

  • 用途：基于 Hessian 矩阵的特征值来增强线状（如血管）、板状、点状结构。
  • 关键思想：使用 Frangi 等滤波器，通过 Hessian 的特征值组合判断像素属于何种结构。
  • 应用：常用于血管增强。

常见参数：

  • Alpha、Beta、Gamma：控制对不同形状的响应。
  • BrightObject：是否增强亮血管（True）或暗血管（False）。

示例：

using ObjectnessFilterType = itk::HessianToObjectnessMeasureImageFilter< HessianImageType, OutputImageType >;
ObjectnessFilterType::Pointer objectnessFilter = ObjectnessFilterType::New();
objectnessFilter->SetInput(hessianFilter->GetOutput());
objectnessFilter->SetBrightObject(true);

3）组合滤波器：itk::MultiScaleHessianBasedMeasureImageFilter

  • 用途：在多个尺度上自动搜索最佳血管结构响应。
  • 内部使用 Hessian 矩阵和对象性度量。
  • 支持血管增强和粗细结构检测。

示例：

using VesselnessFilterType = itk::MultiScaleHessianBasedMeasureImageFilter< ImageType, HessianImageType, OutputImageType >;
VesselnessFilterType::Pointer vesselness = VesselnessFilterType::New();
vesselness->SetInput(inputImage);
vesselness->SetSigmaMinimum(0.5);
vesselness->SetSigmaMaximum(4.0);
vesselness->SetNumberOfSigmaSteps(10);
vesselness->SetHessianToMeasureFilter(objectnessFilter);

3、血管模型自动分割用法详解

3.1 读取原始图像

在这里插入图片描述

	QString sInputPath = QDir::currentPath() + "vessel.nii.gz";mitk::Image::Pointer mitkBrainImage = nullptr;try{// 自动返回一个 BaseData 的向量（实际可为 Image、Surface 等）std::vector<mitk::BaseData::Pointer> baseData = mitk::IOUtil::Load(sInputPath .toStdString());if (baseData.empty()){std::cerr << "Failed to load the .nii.gz file." << std::endl;return;}mitkBrainImage = dynamic_cast<mitk::Image*>(baseData[0].GetPointer());if (!mitkBrainImage.IsNotNull()){std::cerr << "Failed to load image!" << std::endl;}}catch (const std::exception& e){std::cerr << "An error occurred: " << e.what() << std::endl;return;}// 定义图像类型const unsigned int Dimension = 3;using PixelType = double;using ImageType = itk::Image<PixelType, Dimension>;ImageType::Pointer resultIamge;try{mitk::CastToItkImage(mitkBrainImage, resultIamge);}catch (itk::ExceptionObject& e){std::cerr << "Error reading input image: " << e << std::endl;return;}
3.2 利用bet2算法提取脑组织（可以跳过）

  可以参考另一篇博文进行脑组织提取。当然也可以不去除颅骨，获取整个脑袋的血管。
在这里插入图片描述

3.3 采用Hessian矩阵进行血管增强

在这里插入图片描述

using ObjectnessFilterType = itk::HessianToObjectnessMeasureImageFilter<HessianImageType, ImageType>;ObjectnessFilterType::Pointer objectness_filter = ObjectnessFilterType::New();objectness_filter->SetBrightObject(true);objectness_filter->SetScaleObjectnessMeasure(true);objectness_filter->SetAlpha(0.5);objectness_filter->SetBeta(1.0);objectness_filter->SetGamma(5.0);using MultiScaleFilterType = itk::MultiScaleHessianBasedMeasureImageFilter<ImageType, HessianImageType, ImageType>;MultiScaleFilterType::Pointer multi_scale_filter = MultiScaleFilterType::New();multi_scale_filter->SetInput(resultIamge);multi_scale_filter->SetHessianToMeasureFilter(objectness_filter);multi_scale_filter->SetSigmaStepMethodToLogarithmic();multi_scale_filter->SetSigmaMinimum(sigma_minimum);multi_scale_filter->SetSigmaMaximum(sigma_maximum);multi_scale_filter->SetNumberOfSigmaSteps(number_of_sigma_steps);// 保存血管增强后的图像using WriterType1 = itk::ImageFileWriter<ImageType>;WriterType1::Pointer writer1 = WriterType1::New();writer1->SetInput(multi_scale_filter->GetOutput());QString sPathName = "Vessel_step1.mha";writer1->SetFileName(sPathName.toStdString());try{writer1->Update();}catch (itk::ExceptionObject& e){std::cerr << "Error writing step1.mha: " << e << std::endl;return;}
3.4 将增强后的血管归一化到0-255范围内

在这里插入图片描述

		using RescaleFilterType = itk::RescaleIntensityImageFilter<ImageType, OutputImageType>;RescaleFilterType::Pointer rescale_filter = RescaleFilterType::New();rescale_filter->SetInput(multi_scale_filter->GetOutput());// 保存归一化后的图像using WriterType2 = itk::ImageFileWriter<OutputImageType>;WriterType2::Pointer writer2 = WriterType2::New();writer2->SetInput(rescale_filter->GetOutput());sPathName = "Vessel_step2.mha";writer2->SetFileName(sPathName.toStdString());try{writer2->Update();}catch (itk::ExceptionObject& e){std::cerr << "Error writing step2.mha: " << e << std::endl;return;}
3.5 对归一化的血管进行固定阈值分割

在这里插入图片描述

		using ThresholdFilterType = itk::BinaryThresholdImageFilter<OutputImageType, OutputImageType>;ThresholdFilterType::Pointer thresholdFilter = ThresholdFilterType::New();thresholdFilter->SetInput(rescale_filter->GetOutput());thresholdFilter->SetLowerThreshold(lowerThreshold);thresholdFilter->SetUpperThreshold(255);thresholdFilter->SetOutsideValue(0);thresholdFilter->SetInsideValue(255);// 保存阈值分割后的图像using WriterType3 = itk::ImageFileWriter<OutputImageType>;WriterType3::Pointer writer3 = WriterType3::New();writer3->SetInput(thresholdFilter->GetOutput());sPathName = "Vessel_step3.mha";writer3->SetFileName(sPathName.toStdString());try{writer3->Update();}catch (itk::ExceptionObject& e){std::cerr << "Error writing step3.mha: " << e << std::endl;return;}
3.6 对分割的血管进行平滑处理去除噪点

在这里插入图片描述

mitk::Image::Pointer mitkVesselImage;try{std::vector<mitk::BaseData::Pointer> baseData = mitk::IOUtil::Load(sPathName.toStdString());if (baseData.empty()){std::cerr << "Failed to load the .nrrd file." << std::endl;return;}// 假设读取的是图像数据，将其转换为 mitk::Image 类型mitkVesselImage = dynamic_cast<mitk::Image*>(baseData[0].GetPointer());if (!mitkVesselImage.IsNotNull()){std::cerr << "The loaded data is not a valid image." << std::endl;return;}}catch (const std::exception& e){std::cerr << "An error occurred: " << e.what() << std::endl;return;}mitk::ImageToSurfaceFilter::Pointer imageToSurfaceFilter = mitk::ImageToSurfaceFilter::New();imageToSurfaceFilter->SetInput(mitkVesselImage);imageToSurfaceFilter->SetThreshold(0.5);imageToSurfaceFilter->SetSmooth(true);imageToSurfaceFilter->SetDecimate(mitk::ImageToSurfaceFilter::NoDecimation); // decimate laterimageToSurfaceFilter->Update();vtkSmartPointer<vtkPolyData> polyData = imageToSurfaceFilter->GetOutput()->GetVtkPolyData();
http://www.xdnf.cn/news/337645.html

相关文章：

  • 认识不同格式的点云数据 -PCD点云数据 文本点云数据
  • 【前端】webstorm运行程序浏览器报network error
  • PyTorch 版本、torchvision 版本和 Python 版本的对应关系
  • iptables 访问控制列表使用记录
  • 在UI原型设计中，低、高保真原型图有什么区别？
  • 用NVivo革新企业创新：洞悉市场情绪，引领金融未来
  • 使用thymeleaf模版导出swagger3的word格式接口文档
  • 【每天学习一点点】使用Python的pathlib模块分割文件路径
  • HBuilderX安卓真机运行安装失败解决汇总
  • Git实战经验分享：深入掌握git commit --amend的进阶技巧
  • zookeeper实现分布式获取全局唯一自增ID的案例。
  • 论文速读《DARE：基于扩散模型的自主机器人探索新范式》
  • 【Linux网络】网络命令
  • 基于LSTM与SHAP可解释性分析的神经网络回归预测模型【MATLAB】
  • 基于vueflow可拖拽元素的示例(基于官网示例的单文件示例)
  • 深入解析 C# 常用数据结构：特点、区别与优缺点分析
  • C/C++内存分布
  • JVM——Java虚拟机是怎么实现synchronized的？
  • 力扣刷题Day 43：矩阵置零（73）
  • 【随笔】Google学术：but your computer or network may be sending automated queries.
  • 红黑树的应用场景 —— 进程调度 CFS 与内存管理
  • 基于SpringBoot的校园周边美食探索及分享平台的设计与实现
  • Linux系统下使用Kafka和Zookeeper
  • C++ | 常用语法笔记
  • 宝塔面板部署 springboot + mysql 项目
  • CMake笔记（简易教程）
  • 【探寻C++之旅】第十三章：红黑树
  • 第8章-3 查询性能优化1
  • kotlin @JvmStatic注解的作用和使用场景
  • 《信息论与编码课程笔记》——信源编码（1）