ITK学习笔记(一)ITK的使用

  • 第一个ITK程序
    • 1、CMakeLists.txt
    • 2、 测试数据
    • 3、代码
    • 4、结果

ITK官方文档

第一个ITK程序

1、CMakeLists.txt

# This is the root ITK CMakeLists file.
cmake_minimum_required(VERSION 3.10)# This project is designed to be built outside the Insight source tree.
project(ITK_demo)# Find VTK
set(ITK_DIR D:/ProgramFiles/ITK-5.1.1/lib/cmake/ITK-5.1) # add wmz
find_package(ITK REQUIRED)
include_directories(${ITK_INCLUDE_DIRS})
message("ITK dir = ${ITK_INCLUDE_DIRS}")
message("ITK lib = ${ITK_LIBRARIES}") include(${ITK_USE_FILE})
#aux_source_directory(src SRC_LIST)
set(SRC_LIST./itk_demo.cpp)add_executable(itk_demo ${SRC_LIST} )
target_link_libraries(itk_demo ${ITK_LIBRARIES})

关于 include(${ITK_USE_FILE}) 的说明可以在 UseITK.cmake 中找到:

#  -------------
#
# This file is not part of the ITK API.  It exists purely as an
# implementation detail.  This CMake module may change from version to
# version without notice, or even be removed.
#
# We mean it.
## This file sets up include directories, link directories, IO settings and
# compiler settings for a project to use ITK.  It should not be
# included directly, but rather through the ITK_USE_FILE setting
# obtained from ITKConfig.cmake.

2、 测试数据

测试数据下载路径:https://github.com/InsightSoftwareConsortium/ITK/tree/master/Examples/Data
其实编译ITK时的目录下就有需要的测试数据,比如我的ITK-5.1.1目录下。

\ITK-5.1.1\Examples\Data

3、代码

作为第一个示例程序本来应该写一个很简单的像HelloWorld的程序,但是一些比较简单的官网的程序 要么依赖VTK,要么版本高于ITK5.1.1.
所以就找了一个比较长的程序,是一个配准的程序。
代码来自:https://github.com/InsightSoftwareConsortium/ITK/blob/master/Examples/RegistrationITKv4/MultiResImageRegistration1.cxx
我找了一个其他人做过的中文注释版

#include "itkImageRegistrationMethodv4.h"
#include "itkTranslationTransform.h"
#include "itkMeanSquaresImageToImageMetricv4.h"
#include "itkRegularStepGradientDescentOptimizerv4.h"
#include "itkImageFileReader.h"
#include "itkImageFileWriter.h"
#include "itkPNGImageIOFactory.h"
#include "itkResampleImageFilter.h"
#include "itkCastImageFilter.h"
#include "itkRescaleIntensityImageFilter.h"
#include "itkSubtractImageFilter.h"
/*****************************************************************************************************************
* 本例子是一个图像配准的Demo
*   0、创建了一个Command对象,用于监控配准的过程,被后面的对象调用
*   1、首先要定义像素的维度以及像素类型:进进而链接参考图像以及浮动图像
*   2、定义框架的基本组件:
*               确定变换种类:TransformType:二维变换
*               确定优化方法:OptimizerType:梯度下降
*               确定相似度度量:MetricType:链接两个图像:浮动图像以及参考图像
*   3、创建图像组件,并且通过创建上述框架,进而进行设置(链接)
*   4、设置插值方法:LinearInterpolateImageFunction并且链接在一起
*   5、6:通过ImageFileReader方法进行读取,链接到 registration并更新
*   7、针对前面的TransformType进行实例化:平移变换用于配准SetInitialTransformParameters:用于设置初始值
*   8、针对优化方法的设置:OptimizerType:前面在创建的时候已经设置了其梯度下降方法,此步骤用于对其微调:初始步长,收敛公差,最大迭代次数
*   9、通过RegistrationParameterScalesFromPhysicalShift:将每一个配准要素链接到配准方法中执行,
*   10、实例化Common对象,监控配准过程的执行,触发配准过程--迭代
*   11、通过update函数触发配准的执行
*   12、配准结果定义空间变换的参数序列:其结果由GetLastTransformParameters( )获得并且输出
*           X、Y的变换:TranslationAlongX;TranslationAlongY
*           迭代次数:numberOfIterations
*           最后的结果:bestValue
*       通过CompositeTransform:AddTransform将转换添加到堆栈的背面,并且拥有可优化的参数。
*           也就是说:添加堆栈,副本??
*   13、14、15、16、ResampleFilterType方法:
*           用变换参数将两幅图像进行叠加比较,并设置重采样滤波器:输入两幅图像
*           输出的是一个变换
*           对滤波器进行相关参数的设置:大小、原点、间距、位置
*           并通过CastFilterType:setInput:weiter进行相关的输出
*           此时:这个图象就是配准结束后的图像
*   17、通过itk::SubtractImageFilter对两幅图像进行比较:
*           fixedImageReader;resampler
*   18、对图像进行处理 itk::RescaleIntensityImageFilter:调节一下亮度;并进行输出
*   19、一致性转发计算参考图像与正在移动图像之间的不同,输出图片5
******************************************************************************************************************/
/*CommandIterationUpdate 类:继承Command,监视配准过程的执行。每调用一次,输出相应参数object类指向事件的观察者Execute方法,类似cellbake,回转observer方法:*/
class CommandIterationUpdate : public itk::Command
{public:typedef CommandIterationUpdate   Self;typedef itk::Command             Superclass;typedef itk::SmartPointer<Self>  Pointer;itkNewMacro(Self);//宏,包装了所有的new()所有代码
protected:CommandIterationUpdate() {};public:typedef itk::RegularStepGradientDescentOptimizerv4<double> OptimizerType;typedef const OptimizerType* OptimizerPointer;void Execute(itk::Object* caller, const itk::EventObject& event) ITK_OVERRIDE{Execute((const itk::Object*)caller, event);}//Object表示激活事件的对象,event表示被激活的事件void Execute(const itk::Object* object, const itk::EventObject& event) ITK_OVERRIDE{OptimizerPointer optimizer = static_cast<OptimizerPointer>(object);//checkEvent表示是否观察的对象if (!itk::IterationEvent().CheckEvent(&event)){return;}std::cout << optimizer->GetCurrentIteration() << " = ";std::cout << optimizer->GetValue() << " : ";std::cout << optimizer->GetCurrentPosition() << std::endl;}};int main()
{//1、定义图像的维度以及像素执行const    unsigned int    Dimension = 2;//定义维度typedef  float           PixelType;//图像像素类型typedef itk::Image< PixelType, Dimension >  FixedImageType;//输入数据的类型:参考图像 typedef itk::Image< PixelType, Dimension >  MovingImageType;//浮动图像//2、定义配准框架的基本组件:变换、优化、测度配准组件//用于2D图像的一个刚性配准,变换的唯一参数是:空间坐标类型//配准typedef itk::TranslationTransform< double, Dimension > TransformType;//把参考图像的空间映射到待配准图像的映射//优化typedef itk::RegularStepGradientDescentOptimizerv4<double> OptimizerType;//优化算法:牛顿梯度下降法//度量typedef itk::MeanSquaresImageToImageMetricv4<//相似度测量:均方根FixedImageType,MovingImageType >    MetricType;//3、该组件用用于连接其他组件typedef itk::ImageRegistrationMethodv4<FixedImageType,MovingImageType,TransformType   >    RegistrationType;MetricType::Pointer         metric = MetricType::New();OptimizerType::Pointer      optimizer = OptimizerType::New();RegistrationType::Pointer   registration = RegistrationType::New();//连接组件:变换、优化组件registration->SetMetric(metric);registration->SetOptimizer(optimizer);//4、插值方法typedef itk::LinearInterpolateImageFunction<//选择校对机类型,校对机会对配准图像在非网格位置的程度进行评估FixedImageType,double > FixedLinearInterpolatorType;typedef itk::LinearInterpolateImageFunction<MovingImageType,double > MovingLinearInterpolatorType;FixedLinearInterpolatorType::Pointer fixedInterpolator =//每一个配准要素需要其new创建FixedLinearInterpolatorType::New();MovingLinearInterpolatorType::Pointer movingInterpolator =MovingLinearInterpolatorType::New();metric->SetFixedInterpolator(fixedInterpolator);metric->SetMovingInterpolator(movingInterpolator);//5、设置待配准图像以及变换区域typedef itk::ImageFileReader< FixedImageType  >   FixedImageReaderType;typedef itk::ImageFileReader< MovingImageType >   MovingImageReaderType;FixedImageReaderType::Pointer   fixedImageReader = FixedImageReaderType::New();MovingImageReaderType::Pointer  movingImageReader = MovingImageReaderType::New();//6、读图像!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!itk::PNGImageIOFactory::RegisterOneFactory();fixedImageReader->SetFileName("E:\\documents\\vs2019\\itk_demo\\build\\RelWithDebInfo\\data\\BrainProtonDensitySliceBorder20.png");//输入图像文件movingImageReader->SetFileName("E:\\documents\\vs2019\\itk_demo\\build\\RelWithDebInfo\\data\\BrainProtonDensitySliceShifted13x17y.png");//因为图像是从文件读取的,所以下面方法是用于获取图像数据的//需要 itk::ImageRegistrationMethod 从 file readers 的输出获取输入registration->SetFixedImage(fixedImageReader->GetOutput());registration->SetMovingImage(movingImageReader->GetOutput());//更新reader,确保其有效// fixedImageReader->Update();// movingImageReader->Update();//7、平移变换用于配准TransformType::Pointer movingInitialTransform = TransformType::New();TransformType::ParametersType initialParameters(movingInitialTransform->GetNumberOfParameters());initialParameters[0] = 0.0;  // Initial offset in mm along XinitialParameters[1] = 0.0;  // Initial offset in mm along YmovingInitialTransform->SetParameters(initialParameters);registration->SetMovingInitialTransform(movingInitialTransform);//8、准备执行配准方法:对优化器参数进行微调TransformType::Pointer   identityTransform = TransformType::New();identityTransform->SetIdentity();registration->SetFixedInitialTransform(identityTransform);//初始振幅的长度用SetMaximumStepLength( ) 定义//建立迭代的次数需要谨慎。最大数用SetNumberOfIterations()定义:optimizer->SetLearningRate(4);optimizer->SetMinimumStepLength(0.001);//优化器的收敛公差optimizer->SetRelaxationFactor(0.5);//9、将每一个配准要素连接到配准方法执行中bool useEstimator = false;//useEstimator = atoi(argv[6]) != 0;if (useEstimator){typedef itk::RegistrationParameterScalesFromPhysicalShift<MetricType> ScalesEstimatorType;ScalesEstimatorType::Pointer scalesEstimator = ScalesEstimatorType::New();scalesEstimator->SetMetric(metric);scalesEstimator->SetTransformForward(true);optimizer->SetScalesEstimator(scalesEstimator);optimizer->SetDoEstimateLearningRateOnce(true);}optimizer->SetNumberOfIterations(200);//最大迭代次数//10、实例化commend对象,监视配准过程的执行,并处触发配准过程CommandIterationUpdate::Pointer observer = CommandIterationUpdate::New();optimizer->AddObserver(itk::IterationEvent(), observer);const unsigned int numberOfLevels = 1;RegistrationType::ShrinkFactorsArrayType shrinkFactorsPerLevel;shrinkFactorsPerLevel.SetSize(1);shrinkFactorsPerLevel[0] = 1;RegistrationType::SmoothingSigmasArrayType smoothingSigmasPerLevel;smoothingSigmasPerLevel.SetSize(1);smoothingSigmasPerLevel[0] = 0;registration->SetNumberOfLevels(numberOfLevels);registration->SetSmoothingSigmasPerLevel(smoothingSigmasPerLevel);registration->SetShrinkFactorsPerLevel(shrinkFactorsPerLevel);//11、通过调用Update函数触发配准执行try{registration->Update();std::cout << "Optimizer stop condition: "<< registration->GetOptimizer()->GetStopConditionDescription()<< std::endl;}catch (itk::ExceptionObject& err){std::cerr << "ExceptionObject caught !" << std::endl;std::cerr << err << std::endl;return EXIT_FAILURE;}//12、配准结果是一系列定义空间变换的参数序列,结果由get获得TransformType::ConstPointer transform = registration->GetTransform();TransformType::ParametersType finalParameters = transform->GetParameters();const double TranslationAlongX = finalParameters[0];//队列中每个元素对应着沿着一个空间维度的平移const double TranslationAlongY = finalParameters[1];//优化器能够询问抵达收敛的迭代的实际次数并通过GetCurrentIteration()返回出来const unsigned int numberOfIterations = optimizer->GetCurrentIteration();//迭代次数//最终参数集合的图像量规值通过优化器的GetValue();const double bestValue = optimizer->GetValue();//最优化的度量//将上述输出std::cout << "Result = " << std::endl;std::cout << " Translation X = " << TranslationAlongX << std::endl;//输出移动X的值std::cout << " Translation Y = " << TranslationAlongY << std::endl;//输出移动Y的值std::cout << " Iterations    = " << numberOfIterations << std::endl;//输出迭代次数std::cout << " Metric value  = " << bestValue << std::endl;//输出优化的度量typedef itk::CompositeTransform<double,Dimension > CompositeTransformType;CompositeTransformType::Pointer outputCompositeTransform =CompositeTransformType::New();outputCompositeTransform->AddTransform(movingInitialTransform);outputCompositeTransform->AddTransform(registration->GetModifiableTransform());//13、用变换结果将待配准图映射到参考图像中typedef itk::ResampleImageFilter<MovingImageType,FixedImageType >    ResampleFilterType;//14、创建一个重采样滤波器,输入待配准图像ResampleFilterType::Pointer resampler = ResampleFilterType::New();resampler->SetInput(movingImageReader->GetOutput());//配准函数生成的变换也作为重采样滤波器的输入被传递resampler->SetTransform(outputCompositeTransform);//15、ResampleImageFilter要求指定额外的参数,特别是输出图像的间 距、原点和大小FixedImageType::Pointer fixedImage = fixedImageReader->GetOutput();resampler->SetSize(fixedImage->GetLargestPossibleRegion().GetSize());//尺寸resampler->SetOutputOrigin(fixedImage->GetOrigin());//原点resampler->SetOutputSpacing(fixedImage->GetSpacing());//间距resampler->SetOutputDirection(fixedImage->GetDirection());//位置resampler->SetDefaultPixelValue(100);//16、滤波器的输出被传递给一个在文件中存储图像的writertypedef unsigned char OutputPixelType;typedef itk::Image< OutputPixelType, Dimension > OutputImageType;typedef itk::CastImageFilter<//转化重采样的像素类型到最终的writer类型FixedImageType,OutputImageType >          CastFilterType;typedef itk::ImageFileWriter< OutputImageType >  WriterType;//调用new函数创建新的滤波器WriterType::Pointer      writer = WriterType::New();CastFilterType::Pointer  caster = CastFilterType::New();writer->SetFileName("E:\\documents\\vs2019\\itk_demo\\build\\RelWithDebInfo\\output\\RegistrationITKv4Moving13x17yInputType.png");//写到文件夹位置caster->SetInput(resampler->GetOutput());writer->SetInput(caster->GetOutput());writer->Update();//触发更新//17、参照图像和被变换的待配准图像很容易用itk::SubtractImageFilter比较//pixel-wise滤波器 计算两幅输入的同源像素的不同:typedef itk::SubtractImageFilter<FixedImageType,FixedImageType,FixedImageType > DifferenceFilterType;DifferenceFilterType::Pointer difference = DifferenceFilterType::New();difference->SetInput1(fixedImageReader->GetOutput());//不同difference->SetInput2(resampler->GetOutput());//18、两幅图像的不同也许比较暗,我们用下面方法对其进行调节亮度,使之更加的明显typedef itk::RescaleIntensityImageFilter<FixedImageType,OutputImageType >   RescalerType;RescalerType::Pointer intensityRescaler = RescalerType::New();intensityRescaler->SetInput(difference->GetOutput());intensityRescaler->SetOutputMinimum(0);intensityRescaler->SetOutputMaximum(255);resampler->SetDefaultPixelValue(1);//输出到另外一个位置(调亮)WriterType::Pointer writer2 = WriterType::New();writer2->SetInput(intensityRescaler->GetOutput());writer2->SetFileName("E:\\documents\\vs2019\\itk_demo\\build\\RelWithDebInfo\\output\\Moving13x17yInputType.png");writer2->Update();//设置了一致性转换,计算参考图像的不同resampler->SetTransform(identityTransform);writer2->SetFileName("E:\\documents\\vs2019\\itk_demo\\build\\RelWithDebInfo\\output\\DifferenceBeforeRegistration.png");writer2->Update();return EXIT_SUCCESS;}

4、结果

输入图像

输出图像

ITK学习笔记(一)ITK的使用相关推荐

  1. ITK学习笔记(九) simpleITK重采样数据到固定尺寸

    ITK学习笔记(九) simpleITK重采样数据到固定尺寸 1.resampleSize 2.CT图像预处理之重采样 3.SKIMAGE.TRANSFORM-两种缩放皆可 4.Python Simp ...

  2. ITK学习笔记(二)最大圆度连通域提取

    ITK学习笔记(二)最大连通域提取 1.代码 1.1 code 1.2.生成的mask图片 1.3.打印信息 关于编译部分可以参考我的:VTK学习笔记(十)基于VTK和ITK程序 最大圆度概念: 圆度 ...

  3. ITK学习笔记(八) ITK高斯混合模型 GMM EM

    ITK学习笔记(八) ITK高斯混合模型 GMM EM 1.高斯混合模型 2.变分贝叶斯高斯混合 3.ITK中的GMM.EM 1.高斯混合模型 sklearn.mixture是一个能够学习高斯混合模型 ...

  4. ITK学习笔记(四)基于C++和ITK获取3D医学图像的最大连通区域

    ITK学习笔记(四)基于C++和ITK获取3D医学图像的最大连通区域 1.代码 2.效果 3.主要函数介绍 1.代码 #include "itkConnectedComponentImage ...

  5. ITK学习笔记(七) ITK旋转方向位置不变

    ITK学习笔记(七) ITK旋转方向位置不变 sitk使用的使LPS坐标系,与slicer中使用的RAS坐标系不同. 在sitk中的方向(1,1,1)(1,1,1)(1,1,1),在slicer中的方 ...

  6. ITK学习笔记(五) ITK计算质心

    ITK学习笔记(五) ITK计算质心 通过LabelGeometryImageFilter可以求Label Image各个区域的质心和BoundingBox.然而,如果图像的spacing不是1,它并 ...

  7. ITK学习笔记(十一) SimpleITK相关使用

    ITK学习笔记(十一) SimpleITK相关使用 参考:[SimpleITK]医疗影像分割结果评价指标计算 参考:Python SimpleITK.GetArrayViewFromImage方法代码 ...

  8. ITK学习笔记(十二) SimpleITK获取二值图像bbox

    ITK学习笔记(十二) SimpleITK获取二值图像bbox sitk_ seg_ img = sitk . ReadImage( m ) bbox =get_ bbox_ from mask(si ...

  9. ITK学习笔记——OSTU自动选取最佳阈值进行二值化

    原理介绍: Ostu大律法方法又名最大类间差方法,通过统计整个图像的直方图特性来实现全局阈值T的自动选取,其算法步骤为: 先计算图像的直方图,将图像所有的像素点按照0~255共256个bin,统计落在 ...

最新文章

  1. python ** 运算符_Python从入门到放弃运算符(2)
  2. rocketmq docker集群_【解构云原生】RocketMQ高可用方案调研及On K8S设计展望(下)...
  3. Activity 在横竖屏切换情况下的生命周期变化
  4. Kerberos认证过程学习理解
  5. hypertrm超级终端_win7电脑怎么添加windows超级终端
  6. CentOS二进制安装Kubernetes
  7. js 父子标签同时设置onclick,子标签触发父标签onclick解决办法
  8. mysql练习_MySQL基础知识—习题练习
  9. 海康设备添加到萤石云
  10. clockworkmod CWM简单介绍
  11. Java的图标和由来
  12. 在AndroidMenifest动态引用ApplicationId,解决INSTALL_FAILURE_CONFLICTION_PROVIDER
  13. better-scroll 在vue中使用on监听事件,事件无效问题解决
  14. mac 升级java_mac版java更新升级方法
  15. 爬虫实战【6】Ajax内容解析-今日头条图集
  16. 全量发布与灰度发布(APP版本发布)
  17. Java对base64编解码总结
  18. 简单概述PCIe总线协议
  19. 一键分析Android的BugReport
  20. IP地址更改小工具(bat命令)

热门文章

  1. 数据库自动收缩造成的阻塞
  2. 日本金融巨头SBI的加密货币子公司支持移动钱包BRD
  3. 树莓派开发之旅 (三)【登录树莓派的4种方式】
  4. 手把手教你写一个获取短信验证码业务(基于Spring Boot+Redis+阿里云短信服务)
  5. 编译错误jni/Android.mk:12: *** commands commence before first targe
  6. 安装rocon、kobuki、turtlebot
  7. 直播预告 | 腾讯云工业AI系列直播
  8. ceph源码分析--monitor的lease机制
  9. python画混淆矩阵图
  10. 新版TCGA的甲基化数据分析