【图像重建】基于先验和运动的PRIMOR算法实现呼吸门控 CT图像重建附matlab代码

发布于:2023-02-02 ⋅ 阅读:(337) ⋅ 点赞:(0)

 1 内容介绍

心胸小动物成像中呼吸门控的低剂量方案导致使用 Feldkamp-Davis-Kress (FDK) 方法重建的图像中的条纹伪影。我们提出了一种新颖的基于先验和基于运动的重建(PRIMOR)方法,该方法通过添加一个惩罚函数来改进基于先验的重建 (PBR)运动模型。先验图像是作为所有呼吸门的平均值生成的,用FDK重建。使用非刚性估计呼吸门之间的运动基于层次B样条的配准方法我们将 PRIMOR 与没有运动估计的等效 PBR 方法进行比较,使用高剂量数据。从这些用微型 CT 扫描仪获得的数据中,我们发现了不同的场景通过改变光子通量和投影数量来模拟。方法在对比度噪声比 (CNR)、均方误差 (MSE)、条纹伪影指标(SAI)、解误差范数 (SEN) 和呼吸运动校正。此外,要评估每种方法对肺研究量化的影响,我们计算了从分割每个图像获得的掩码的 Jaccard 相似性指数,与获得的相比较从高剂量重建。在所有情况下,这两种迭代方法都极大地改进了 FDK 重建。PBR 容易出现条纹伪影并在骨骼中呈现模糊效果和肺组织同时使用少量投射和低剂量。采用 PBR作为参考,PRIMOR 将 CNR 提高了 33%,并降低了 MSE、SAI 和 SEN分别为 20%、4% 和 13%。PRIMOR 还为呼吸运动提供了更好的补偿和更高的 Jaccard 相似指数。总之,提出的新方法小动物扫描仪中的低剂量呼吸门控显示图像改善质量,并允许减少剂量或减少之间的投影数量是以前 PBR 方法的两倍和三倍。

2 仿真代码

function [GridAll,SpacingAll,GridDAll,SpacingDAll] = ComputeSplineRegTwoDirectionsPoolDouble(u0,Options)
% [GridAll,SpacingAll,GridDAll,SpacingDAll] = ComputeSplineRegTwoDirectionsPoolDouble(u0,Options)
%
% Computes the registration between consecutive gates, forward and
% backwards to be used by PRIMOR method. 
%
% Inputs:
%
% u0        = image of all gates, size n_x x n_y x number gates
% Options   = parametes for the FFD-based registration package
% Options.Penalty     = 1e-4;
% Options.Registration= 'NonRigid';
% Options.MaxRef      = 3;
% Options.Similarity  = 'sd';
%
% Outputs: 
%
% GridAll,SpacingAll,GridDAll,SpacingDAll = mesh of grid points and
% spacing for forward and backward registration required by PRIMOR method

% This version uses the double of processors for registration forward and
% backwards
%
% Requirements: 
%
% FFD-based registration package: For the motion estimation step we used
% the FFD-based registration method available in MATLAB Central (Dirk-Jan
% Kroon; B-spline grid,  image and point based registration; 2012,
% retrieved from 
% http://www.mathworks.com/matlabcentral/fileexchange/20057-b-spline-grid-image-and-point-based-registration), 
%
% Juan FPJ Abascal, Monica Abella
% Departamento de Bioingenieria e Ingenieria Aeroespacial
% Universidad Carlos III de Madrid, Madrid, Spain
% mabella@hggm.es, juanabascal78@gmail.com, desco@hggm.es

dimIm       = size(u0);
numTime     = dimIm(3);

parfor ih = 1:2*numTime
    if ih <= numTime
        % Backward
        Imoving             = abs(u0(:,:,ih));
    
        % Up-transformation and down-transformations
        if ih == numTime
            [Ireg,Grid,Spacing,M,B,F]=image_registration(Imoving,abs(u0(:,:,1)),Options);
        elseif ih == 1
            [Ireg,Grid,Spacing,M,B,F]=image_registration(Imoving,abs(u0(:,:,ih+1)),Options);
        else
            [Ireg,Grid,Spacing,M,B,F]=image_registration(Imoving,abs(u0(:,:,ih+1)),Options);
        end % numtime

        GridAll(:,:,:,ih)   = Grid;
        SpacingAll(:,:,ih)  = Spacing;       
    else 
        % Forward transformation
        ih2                 = ih-numTime;
        Imoving             = abs(u0(:,:,ih2));
        
        if ih2 == 2*numTime
            [IregD,GridD,SpacingD,MD,BD,FD]=image_registration(Imoving,abs(u0(:,:,ih2-1)),Options);
        elseif ih2 == 1
            [IregD,GridD,SpacingD,MD,BD,FD]=image_registration(Imoving,abs(u0(:,:,end)),Options);
        else
            [IregD,GridD,SpacingD,MD,BD,FD]=image_registration(Imoving,abs(u0(:,:,ih2-1)),Options);        
        end % numtime       
        GridDAll(:,:,:,ih-numTime)  = GridD;
        SpacingDAll(:,:,ih-numTime) = SpacingD;        
    end % ih2
end % ih

3 运行结果

4 参考文献

博主简介:擅长智能优化算法、神经网络预测、信号处理、元胞自动机、图像处理、路径规划、无人机等多种领域的Matlab仿真,相关matlab代码问题可私信交流。

部分理论引用网络文献,若有侵权联系博主删除。