CN107221013A - 一种基于变分光流估计肺部4d‑ct图像超分辨率处理方法 - Google Patents

一种基于变分光流估计肺部4d‑ct图像超分辨率处理方法 Download PDF

Info

Publication number
CN107221013A
CN107221013A CN201710344897.3A CN201710344897A CN107221013A CN 107221013 A CN107221013 A CN 107221013A CN 201710344897 A CN201710344897 A CN 201710344897A CN 107221013 A CN107221013 A CN 107221013A
Authority
CN
China
Prior art keywords
mrow
msub
msubsup
image
munderover
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Granted
Application number
CN201710344897.3A
Other languages
English (en)
Other versions
CN107221013B (zh
Inventor
刘慧�
耿凤欢
郭强
张小峰
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Ludong University
Shandong University of Finance and Economics
Original Assignee
Ludong University
Shandong University of Finance and Economics
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Ludong University, Shandong University of Finance and Economics filed Critical Ludong University
Priority to CN201710344897.3A priority Critical patent/CN107221013B/zh
Publication of CN107221013A publication Critical patent/CN107221013A/zh
Application granted granted Critical
Publication of CN107221013B publication Critical patent/CN107221013B/zh
Expired - Fee Related legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Classifications

    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T11/00—2D [Two Dimensional] image generation
    • G06T11/003—Reconstruction from projections, e.g. tomography
    • G06T11/006—Inverse problem, transformation from projection-space into object-space, e.g. transform methods, back-projection, algebraic methods
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00—Image analysis
    • G06T7/0002—Inspection of images, e.g. flaw detection
    • G06T7/0012—Biomedical image inspection
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00—Image analysis
    • G06T7/20—Analysis of motion
    • G06T7/269—Analysis of motion using gradient-based methods
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00—Indexing scheme for image analysis or image enhancement
    • G06T2207/10—Image acquisition modality
    • G06T2207/10072—Tomographic images
    • G06T2207/10081—Computed x-ray tomography [CT]
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00—Indexing scheme for image analysis or image enhancement
    • G06T2207/20—Special algorithmic details
    • G06T2207/20024—Filtering details
    • G06T2207/20028—Bilateral filtering
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00—Indexing scheme for image analysis or image enhancement
    • G06T2207/20—Special algorithmic details
    • G06T2207/20172—Image enhancement details
    • G06T2207/20192—Edge enhancement; Edge preservation
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2211/00—Image generation
    • G06T2211/40—Computed tomography
    • G06T2211/421—Filtered back projection [FBP]
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2211/00—Image generation
    • G06T2211/40—Computed tomography
    • G06T2211/424—Iterative

Landscapes

  • Engineering & Computer Science (AREA)
  • Theoretical Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Computer Vision & Pattern Recognition (AREA)
  • Radiology & Medical Imaging (AREA)
  • Algebra (AREA)
  • General Health & Medical Sciences (AREA)
  • Medical Informatics (AREA)
  • Nuclear Medicine, Radiotherapy & Molecular Imaging (AREA)
  • Multimedia (AREA)
  • Quality & Reliability (AREA)
  • Health & Medical Sciences (AREA)
  • Mathematical Analysis (AREA)
  • Mathematical Optimization (AREA)
  • Mathematical Physics (AREA)
  • Pure & Applied Mathematics (AREA)
  • Other Investigation Or Analysis Of Materials By Electrical Means (AREA)
  • Image Processing (AREA)
  • Image Analysis (AREA)

Abstract

本发明提供一种基于变分光流估计肺部4D‑CT图像超分辨率处理方法，构建了一个用于求解肺部4D‑CT不同相位图像之间的光流场的变分光流模型；利用快速交替方向乘子法求解该模型，得到不同相位图像之间的光流场；基于光流场，并利用非局部迭代反投影超分辨率重建算法，实现了高分辨率肺部图像的重建。本方法在增强图像纹理结构的同时更好地保留了图像的轮廓。影响超分辨率重建算法效果的主要因素是图像配准的精度，采用了局部和全局相结合CLG的变分光流估计模型获得更加精准稠密的光流场，并采用快速交替方向乘子法实现光流模型的快速求解；其二，采用基于边缘的非局部迭代反投影算法，有效增强超分辨率图像中的边缘和纹理细。

Description

一种基于变分光流估计肺部4D-CT图像超分辨率处理方法
技术领域
本发明涉及肺部图像处理领域，尤其涉及一种基于变分光流估计肺部4D-CT图像超分辨率处理方法。
背景技术
放射治疗是医学上治疗肺癌最常用的手段。肺部4D-CT图像因其能够提供必要的呼吸运动信息来引导医生进行精确地放射治疗，从而得到越来越多的关注。肺部4D-CT图像是由呼吸运动周期中的多个相位各自对应的3D-CT图像按照获取时间的先后顺序排序得到的，能够包含呼吸周期内肺部器官以及肿瘤区域的空间结构，以及运动信息。这些信息可以辅助医生精确定位靶区，有助于病人的放射治疗。但是，众所周知CT图像的获取伴随着高剂量照射，为了降低对病人的辐射量，往往只能降低沿纵向的采样，从而导致肺部4D-CT图像的层间分辨率远低于层内分辨率。在观察3D肺部图像的冠状图时，为了还原真实肺部形状，通常会将3D图像沿Z轴方向进行插值放大。传统插值算法会使图像变模糊，因此，如何提高4D-CT图像超分辨率重建图像质量是亟待解决的问题。
发明内容
为了克服上述现有技术中的不足，本发明提供一种基于变分光流估计肺部4D-CT图像超分辨率处理方法，处理方法包括：
S1：采用CLG变分光流估计算法得到浮动图像，浮动图像与其余相位图像形成光流场；
S2：结合光流场以及图像放大倍数，得到图像间的仿射变换；
S3：采用改进的非局部迭代反投影算法重建出高分辨率图像。
优选地，步骤S1包括：
基于光流场基本理论，将肺部4D-CT图像看作肺部图像I随时间发生形变而生成的图像序列；即呼吸运动过程中肺部图像各个像素位置的亮度沿着运动轨迹保持不变，即：
I(x(t),y(t),t)＝C (1)
对光流的l1范数约束和各项异性扩散约束，以及对数据项的双边滤波操作，得到肺部4D-CT图像的CLG变分光流估计模型，如下式所示：
其中λ为数据残差项与数据平滑项之间的权重系数，bfw表示对数据项进行双边滤波，U＝(u,v)是待求解的光流场；模型引入了各向异性扩散因子以实现光流的各向异性扩散；F(U)为配准后图像I1与参考帧图像I0的差，定义如下：
其中X＝(x,y)为像素点坐标，U0是光流场的初始估计。
优选地，步骤S2结合光流场以及图像放大倍数，得到图像间的仿射变换具体包括：
定义Ψi(X)为张量X在第i维度上的融合算子，光流模型可转化成如下等价：
其中，Vi,1≤i≤2是松弛变量。
优选地，步骤S3采用改进的非局部迭代反投影算法重建出高分辨率图像具体包括：
设置高分辨率图像的初始估计对初始估计按照图像退化模型生成与集合相对应的低分辨率图像估计集通过反向投影重建误差{I-I(0)}到高分辨率估计来提高超分辨率重建效果；迭代此过程以最小化误差函数；
第n次迭代中，低分辨率图像估计的生成可通过下面的模型模拟实现：
其中，是指图像的仿射变换；h是点扩散函数；↓s为下采样算子；反投影重建过程可表示为：
式中，是Tk的逆；p是反投影算子，p的取值将影响算法迭代时间；↑s表示上采样算子；
对IBP算法中的投影重建误差进行非局部均值滤波，保留图像中的高频细节；
利用数据加权平均法进行融合，根据低分辨率图像光流估计误差确定该图像的重建误差在反投影过程中的权重，抑制光流估计误差对重建结果的影响，反投影过程可表示为：
将ωk设置为Ik与间欧式距离dk的指数函数，即：
其中t是控制核衰变速度的参数；
将图像所包含的所有相位图像作为已知的低分辨率图像集合，从中选出某一相位图像作为高分辨率图像。
优选地，为了求解式(4)中的条件极值问题，将该约束问题转换成无约束的增广拉格朗日函数，如下式所示：
应用ADMM算法得到下面的迭代求解步骤：
Step1：固定Vi和Wi，更新Uk+1：
转换为：
Step2：固定U和Wi，求解Vi,i＝1,2
其中，{Vi}表示集合{V1,V2}；计算Vi,i＝1,2，且等价写成下式
式中的式(9)中分解成规模更小的、并行求解的Lasso问题，即
ν是张量Vi所包含的一维张量，t是Γi中与ν相对应的一维张量，k是ν的维数，di是D中与ν相对应一维张量的第i个元素，式(10)是一维全变分问题求解，并行计算的引入使得Vi的计算可以在毫秒时间内完成；
Step3：计算Wi,i＝1,2
该迭代求解过程在原始残差和对偶残差满足特定条件时结束。
从以上技术方案可以看出，本发明具有以下优点：
基于变分光流估计肺部4D-CT图像超分辨率处理方法构建了一个用于求解肺部4D-CT不同相位图像之间的光流场的变分光流模型；然后，利用快速交替方向乘子法求解该模型，得到不同相位图像之间的光流场；最后，基于光流场，并利用非局部迭代反投影超分辨率重建算法，实现了高分辨率肺部图像的重建.实验结果表明，与已有算法相比，本方法在增强图像纹理结构的同时更好地保留了图像的轮廓。
为了提高重建图像的质量和速度，采用迭代反投影方法，充分利用4D-CT不同相位间的互补信息实现快速的CT图像超分辨率重建，并在以下两个方面进行了创新：其一，鉴于影响超分辨率重建算法效果的主要因素是图像配准的精度，采用了局部和全局相结合CLG的变分光流估计模型获得更加精准稠密的光流场，并采用快速交替方向乘子法实现光流模型的快速求解；其二，采用基于边缘的非局部迭代反投影算法，有效增强超分辨率图像中的边缘和纹理细。
附图说明
为了更清楚地说明本发明的技术方案，下面将对描述中所需要使用的附图作简单地介绍，显而易见地，下面描述中的附图仅仅是本发明的一些实施例，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些附图获得其他的附图。
图1为滤波器大小对光流精度的影响示意图；
图2为CLG变分光流估计算法得到的光流场和根据光流场对基准图像运动补偿后得到补偿图像以及补偿图像与其对应的浮动图像间的差图；
图3为不同算法得到的冠状图超分辨率图像视觉效果对比.前三行图片分别对应以第二份4D-CT数据中第一相位、第五相位以及第八相位为基准图像得到的超分辨率结果，第四行图像依次为第三行图像中红框区域放大后的图像；
图4为不同算法得到的矢状面超分辨率图像视觉效果对比，图片排列方式与图3相一致。
具体实施方式
为使得本发明的发明目的、特征、优点能够更加的明显和易懂，下面将运用具体的实施例及附图，对本发明保护的技术方案进行清楚、完整地描述，显然，下面所描述的实施例仅仅是本发明一部分实施例，而非全部的实施例。基于本专利中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其它实施例，都属于本专利保护的范围。
在本发明中，CLG变分光流估计模型是用于计算视频中相邻帧间光流场的模型，4D-CT数据中不同相位间光流场亦可用此模型求解。本发明针对肺部4D-CT图像的CLG变分光流估计模型，以及利用快速交替方向乘子算法求解该光流模型的过程。
光流估计是利用图像序列中参考帧与当前帧之间像素值数据的时域相关性来确定当前帧中各像素位置相对于参考帧的“运动”。光流计算模型都是基于光流基本等式，通过引入光流约束条件来实现光流求解的。CLG变分光流估计是一种对噪声具有鲁棒性且更加准确的光流计算模型。该模型中引入了对光流的l1范数约束，从而保证了模型对噪声的鲁棒性；模型引入了对数据项的双边滤波操作，以减小相邻像素之间的相互影响；该模型基于图像局部特性实现了光流的各向异性传播，有效地抑制光流过传播带来的负影响。本发明构建适用于肺部4D-CT数据不同相位图像间光流估计的计算模型。
本发明基于光流场基本理论，将肺部4D-CT图像看作肺部图像I随时间发生形变而生成的图像序列，也称为图像一致性假设，即呼吸运动过程中肺部图像各个像素位置的亮度沿着运动轨迹保持不变，即：
I(x(t),y(t),t)＝C (1)
引入对光流的l1范数约束和各项异性扩散约束，以及对数据项的双边滤波操作，得到肺部4D-CT图像的CLG变分光流估计模型，如下式所示：
其中λ为数据残差项与数据平滑项之间的权重系数，bfw表示对数据项进行双边滤波，U＝(u,v)是待求解的光流场；模型引入了各向异性扩散因子以实现光流的各向异性扩散；F(U)为配准后图像I1(即除参考帧外的所有帧图像)与参考帧图像I0的差，定义如下：
其中X＝(x,y)为像素点坐标，U0是光流场的初始估计。数据残差项的局部约束使得邻域像素的光流会影响中心像素点光流的计算，为了确保在得到稠密光流场的同时弱化邻域像素的影响，引入对数据残差项的滤波操作.滤波算法有双边滤波以及基于奇异值分解的滤波算法，其中双边滤波算法实现简单且能满足需求，故本发明采用双边滤波对数据残差项进行约束。此外，在CLG变分光流估计模型的求解过程中，光流场的各向同性传播会造成光流的过传播，从而降低光流场的精度。本发明中CLG变分光流模型通过加入光流的各向异性滤波操作，有效地抑制了光流的过传播。
在本发明中，为了实现3D-CT序列图像间光流场的快速求解，采用快速交替方向乘子法(Alternating Direction Method of Multipliers,ADMM)算法求解本文中的变分光流模型。ADMM是一种求解优化问题的计算框架，通过将全局问题分解成多个相互独立、较小且易于求解的局部子问题，之后协调子问题的解得到全局问题的解.ADMM算法结合了增广拉格朗日方法以及求解凸优化问题的对偶分解策略的优势，实现了对凸优化问题的快速求解.定义Ψi(X)为张量X在第i维度上的融合算子，光流模型可转化成如下等价问题：
其中，Vi,1≤i≤2是松弛变量.为了求解式(4)中的条件极值问题，将该约束问题转换成无约束的增广拉格朗日函数，如下式所示：
应用ADMM算法得到下面的迭代求解步骤：
Step1：固定Vi和Wi，更新Uk+1：
最优解为：
Step2：固定U和Wi，求解Vi,i＝1,2
其中，{Vi}表示集合{V1,V2}。式(8)是可分解的，独立地计算Vi,i＝1,2，且可等价写成下式
式中的式(9)中的问题可以分解成规模更小的、可并行求解的Lasso问题，即
ν是张量Vi所包含的一维张量，t是Γi中与ν相对应的一维张量，k是ν的维数，di是D中与ν相对应一维张量的第i个元素。式(10)是一维全变分问题有很多算法可以快速求解，并行计算的引入使得Vi的计算可以在毫秒时间内完成。
Step3：计算Wi,i＝1,2
该迭代求解过程在原始残差和对偶残差满足特定条件时结束。
由于迭代反投影重建过程中并未考虑边缘的方向和强度信息，导致最终收敛于边缘模糊的高分辨率图像。为了提高图像中边缘区域的质量，非局部迭代反投影算法提出借助非局部冗余信息减小反投影重建误差的单幅图像超分辨率重建算法。
在本发明中，通过迭代地最小化重建误差实现多幅低分辨率图像序列重建得到高分辨率图像的迭代反投影算法。对于输入的低分辨率图像集合首先给出高分辨率图像的初始估计然后对初始估计按照图像退化模型生成与集合相对应的低分辨率图像估计集通过反向投影重建误差{I-I(0)}到高分辨率估计来提高超分辨率重建效果.迭代此过程以最小化误差函数。
第n次迭代中，低分辨率图像估计的生成可通过下面的模型模拟实现:
其中，是指图像的仿射变换；h是点扩散函数；↓s为下采样算子，反投影重建过程可表示为：
式中，是Tk的逆；p是反投影算子，p的取值将影响算法迭代时间；↑s表示上采样算子。
在本发明中，为了更好地利用不同帧图像间的互补信息，结合迭代反投影算法和非局部迭代反投影算法提出了改进的非局部迭代反投影，实现4D-CT图像超分辨率重建。本发明算法在以下两个方面具有创新性：一方面，对IBP算法中的投影重建误差进行非局部均值滤波，来保留图像中的高频细节；另一方面，对反投影重建过程中的图像融合策略进行了改进。
IBP算法在反投影重建过程中采用了灰度平均的像素级图像数据融合策略，该策略简单且便于实现，但在输入低分辨率图像间光流场估计误差较大的情况下会严重影响超分辨率重建结果的质量。基于此，本发明利用数据加权平均法进行融合，根据低分辨率图像光流估计误差确定该图像的重建误差在反投影过程中的权重。这样便可以很好地抑制光流估计误差对重建结果的影响，从而得到高质量的高分辨率图像.修改后的反投影过程可表示为：
将ωk设置为Ik与间欧式距离dk的指数函数，即：
其中t是控制核衰变速度的参数。
本发明将4D-CT图像所包含的所有相位图像作为已知的低分辨率图像集合，从中选出某一相位图像作为高分辨率图像。
在本发明中，通过将多维变分问题转换成标准的l1范数问题，并证明了转换后的问题满足收敛的约束条件，且时间复杂度为O(1/k)，而光流求解算法仅是该算法中的特例(二维变分问题)，所以收敛且时间复杂度为O(1/2)。
当待配准图像与基准间存在较大位移时，式(1)图像灰度一致性假设会造成估计结果误差较大，为此，在求解过程中，利用由粗到细的图像金字塔方法提高光流估计精度。式(2)中各项异性扩散因子D的参数设置为典型值α＝5,β＝1/2，权重系数λ＝1000。增广拉格朗日函数中的惩罚项参数设置为ρ＝10，由光流求解过程可知，双边滤波器的邻域大小的选取对光流估计结果影响最大。因此，采用几种不同大小滤波器并分别对两幅图像间存在大位移和小位移的情况进行了对比实验，如图1所示。
图1中，横轴表示滤波器邻域大小，纵轴表示图像插值误差.插值误差是一种光流估计效果评价标准，其值为真实图像与基于光流补偿得到估计图像间的均方差。由图1可知，当滤波器大小为3×3，5×5时，对于大位移图像间的光流估计效果较差；当滤波器大小为9×9,11×11,13×13时，虽然对图像间存在大位移的情况下光流估计结果较精确，但计算时间较长且对图像间存在小位移情况估计偏差变大了，这是由于滤波器的邻域变大，使得影响当前中心像素的邻域像素变多，从而影响估计结果；当滤波器大小为7×7时，能够很好地处理大位移光流估计问题，且处理时间也能达到要求。因此，实验中将双边滤波器的大小设置为7×7。
图2显示了使用基于CLG变分光流估计算法对不同相位中的三幅冠状面图像进行光流估计的两个实例。这两幅冠状面图像是从同一幅4D-CT图像中选取不同相位上的冠状面图像(a)、(b)和(f)，其中(a)作为参考帧图像，(b)和(f)分别作为图像，使用CLG变分光流估计算法分别得到(a)和(b)之间的光流场图(c)以及(a)和(f)之间的光流场图(g)。利用得到的光流场(c)对图(a)进行运动补偿后得到图(d)；利用得到的光流场(g)对图(a)进行运动补偿后得到图(h)。图(e)是图(b)和图(d)之间的差图像，图(i)是图(f)和图(h)之间的差图像。观察可知，图(b)和图(f)中肺部相对于(a)整体向下运动且(a)和(b)之间的位移较小、(a)和(f)之间的位移较大，这与得到光流场一致.结果可见，(b)和(d)以及(f)和(h)的差异很小，说明CLG变分光流估计能准确估计图像间的位移场，有助于更好地进行超分辨率重建。
实验所用数据是从可形变图像配准库获取的肺部4D-CT数据集，从中选取了5组4D-CT数据(数据1-数据5).每组数据含有呼吸运动过程中的10个不同相位，包括呼气末端相位和吸气末端相位.数据的层内像素尺寸范围从(0.97×0.97)mm2到(1.116×1.16)mm2，层厚为2.5mm.下面将分别展示CLG变分光流估计和最终的超分辨率重建结果。
非局部迭代反投影算法可以很好地抑制高对比度边缘的“振铃”现象。本发明将算法与非局部反投影算法进行了对比，以评估算法在边界保持方面的效果。同时，将基于全搜索的运动估计超分辨算法加入对比，论证光流估计对超分辨重建效果的影响。
图3展示了采用不同算法对数据2进行冠状面图像重建的结果图像.由于数据2中层内像素尺寸为(1.16×1.16)mm2且层厚2.5mm，为了使轴向采样率和层内采样率相同，本发明将数据2冠状图上采样倍数设置为2.15，图3中列图像分别是利用双三次插值(Bicubic)算法、基于全搜索运动估计(Full Search)的迭代反投影重建算法、非局部迭代(NLIBP)反投影重建算法以及算法得到的冠状面图像超分辨率结果。前三行图像是分别对应数据2中不同相位(相位5、相位8和相位1)的超分辨率重建结果。为了更好地对比不同算法重建效果，在第四行给出了第三行图像中红框区域细节的放大图。对比前三行图像可以发现，当以不同相位图像作为基准图像时得到的超分辨率结果图像间存在差异，但在图像质量方面相差无几。对比放大后的肺实质周围的血管与组织边缘处的细节可以看出，本发明算法重建出的图像具有更清晰的结构，同时边缘和细节信息也得到加强。
图4则给出采用不同算法对数据2进行矢状面图像重建得到的结果图像。图4中图片排列方式与图3一致.结果表明算法的超分辨率重建效果最好，与其他算法相比，算法在重建图像的边缘处理上做得最好。
为了更好地评价重建结果图像的清晰度，采用平均梯度作为量化标准。图像边缘附近的灰度值变化率较大，这种变化率大小可以用来表示图像的清晰度。图像平均梯度就是用图像多维方向上灰度变化率之和来表征图像的相对清晰度，值越大图像越清晰，可用公式表示：
其中，I(x,y)是图像在该位置上的灰度值；和分别表示图像在x方向和y方向上的梯度；M和N分别是图像的行数以及列数。
表(1)列出了分别利用双三次插值算法、全搜索运动估计、CLG光流估计以及算法对5幅4D-CT图像选择不同相位图像为浮动图像，对不同平面进行超分辨率重建，然后计算得到所有重建结果图像平均梯度值的平均值。对比表(1)中数据可知，较前三个方法，方法得到的图像平均梯度更高，图像清晰度明显增强。
表1对以下五组数据，用四种不同算法得到的平均梯度
从图3和图4中的细节放大图中可以看出，全搜索算法由于得到的光流场不够精确，导致重建图像中出现一些噪声。本发明算法利用基于CLG变分光流估计算法计算得到比全搜索估计更加精确和稠密的光流场，所以超分辨结果要优于全搜索算法。由于非局部迭代反投影只利用低分辨率图像自身的冗余信息，导致超分辨率结果图像在具有较少信息冗余的器官组织边缘出现了模糊。而本发明算法则是很好地利用不同帧图像间的互补信息，在保留细节信息的同时提高了重建图像的清晰度。对比表1中5份4D-CT数据重建结果图像的平均梯度值可知，本发明算法得到的平均梯度值都是最大的，这表明本发明算法更好地保留了图像的细节信息。
实验用C++语言实现了本发明中光流求解算法，并在配备英特尔I5处理器(主频为3.2GHz)的主机上进行了测试。对于分辨率为256×99的浮动图像和基准图像，光流计算时间大约为300ms，极大地缩短了光流计算时间，得益于更短的光流计算时间，对于10幅图像组成的实验图像序列，使用本发明超分辨率算法可以在30s内得到超分辨结果。
本发明提出了一种基于CLG变分光流估计的超分辨率重建算法，用于提高4D-CT冠状面图像分辨率。将不同相位的3D-CT图像看作不同“帧”图像，基于CLG变分光流模型，并利用快速ADMM算法迅速收敛的性质，实现了浮动“帧”图像与基准“帧”图像间光流场的快速计算；然后，基于得到的光流场，改进了非局部迭代反投影算法，并利用不同相位间的互补信息增强重建了图像的结构细节信息.实验结果表明，本发明算法不仅有效增强了图像的纹理结构，而且能够更好地保留图像中的轮廓信息。
对所公开的实施例的上述说明，使本领域专业技术人员能够实现或使用本发明。对这些实施例的多种修改对本领域的专业技术人员来说将是显而易见的，本文中所定义的一般原理可以在不脱离本发明的精神或范围的情况下，在其它实施例中实现。因此，本发明将不会被限制于本文所示的这些实施例，而是要符合与本文所公开的原理和新颖特点相一致的最宽的范围。

Claims (5)

1.一种基于变分光流估计肺部4D-CT图像超分辨率处理方法，其特征在于，处理方法包括：
S1：采用CLG变分光流估计算法得到浮动图像，浮动图像与其余相位图像形成光流场；
S2：结合光流场以及图像放大倍数，得到图像间的仿射变换；
S3：采用改进的非局部迭代反投影算法重建出高分辨率图像。
2.根据权利要求1所述的基于变分光流估计肺部4D-CT图像超分辨率处理方法，其特征在于，
步骤S1包括：
基于光流场基本理论，将肺部4D-CT图像看作肺部图像I随时间发生形变而生成的图像序列；即呼吸运动过程中肺部图像各个像素位置的亮度沿着运动轨迹保持不变，即：
I(x(t),y(t),t)＝C (1)
对光流的l1范数约束和各项异性扩散约束，以及对数据项的双边滤波操作，得到肺部4D-CT图像的CLG变分光流估计模型，如下式所示：
<mrow> <msub> <mi>E</mi> <mrow> <mi>C</mi> <mi>L</mi> <mi>G</mi> <mo>-</mo> <mi>T</mi> <mi>V</mi> </mrow> </msub> <mo>=</mo> <munder> <mo>&amp;Integral;</mo> <mi>&amp;Omega;</mi> </munder> <mo>(</mo> <mrow> <mi>&amp;lambda;</mi> <munder> <mi>&amp;Sigma;</mi> <mrow> <mi>r</mi> <mi>e</mi> <mi>g</mi> <mi>i</mi> <mi>o</mi> <mi>n</mi> </mrow> </munder> <mi>b</mi> <mi>f</mi> <mi>w</mi> <mo>&amp;CenterDot;</mo> <mi>F</mi> <msup> <mrow> <mo>(</mo> <mi>U</mi> <mo>)</mo> </mrow> <mn>2</mn> </msup> <mo>+</mo> <mi>D</mi> <mo>&amp;CenterDot;</mo> <mo>|</mo> <mo>|</mo> <mi>U</mi> <mo>|</mo> <msub> <mo>|</mo> <mn>1</mn> </msub> </mrow> <mo>)</mo> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>2</mn> <mo>)</mo> </mrow> </mrow>
其中λ为数据残差项与数据平滑项之间的权重系数，bfw表示对数据项进行双边滤波，U＝(u,v)是待求解的光流场；模型引入了各向异性扩散因子以实现光流的各向异性扩散；F(U)为配准后图像I1与参考帧图像I0的差，定义如下：
其中X＝(x,y)为像素点坐标，U0是光流场的初始估计。
3.根据权利要求1所述的基于变分光流估计肺部4D-CT图像超分辨率处理方法，其特征在于，
步骤S2结合光流场以及图像放大倍数，得到图像间的仿射变换具体包括：定义Ψi(X)为张量X在第i维度上的融合算子，光流模型可转化成如下等价：
<mrow> <mtable> <mtr> <mtd> <mrow> <munder> <mrow> <mi>m</mi> <mi>i</mi> <mi>n</mi> </mrow> <mrow> <mi>U</mi> <mo>,</mo> <mi>V</mi> </mrow> </munder> <mrow> <mo>(</mo> <mi>&amp;lambda;</mi> <mo>&amp;CenterDot;</mo> <munder> <mo>&amp;Sigma;</mo> <mrow> <mi>r</mi> <mi>e</mi> <mi>g</mi> <mi>i</mi> <mi>o</mi> <mi>n</mi> </mrow> </munder> <mi>b</mi> <mi>f</mi> <mi>w</mi> <mo>&amp;CenterDot;</mo> <mi>F</mi> <msup> <mrow> <mo>(</mo> <mi>U</mi> <mo>)</mo> </mrow> <mn>2</mn> </msup> <munderover> <mo>&amp;Sigma;</mo> <mrow> <mi>i</mi> <mo>=</mo> <mn>1</mn> </mrow> <mn>2</mn> </munderover> <mi>D</mi> <mo>&amp;CenterDot;</mo> <mi>&amp;Psi;</mi> <mo>(</mo> <msub> <mi>V</mi> <mi>i</mi> </msub> <mo>)</mo> <mo>)</mo> </mrow> </mrow> </mtd> </mtr> <mtr> <mtd> <mtable> <mtr> <mtd> <mrow> <mi>s</mi> <mo>.</mo> <mi>t</mi> <mo>.</mo> </mrow> </mtd> <mtd> <mrow> <mi>U</mi> <mo>=</mo> <msub> <mi>V</mi> <mi>i</mi> </msub> <mo>,</mo> <mi>f</mi> <mi>o</mi> <mi>r</mi> </mrow> </mtd> <mtd> <mrow> <mn>1</mn> <mo>&amp;le;</mo> <mi>i</mi> <mo>&amp;le;</mo> <mn>2</mn> </mrow> </mtd> </mtr> </mtable> </mtd> </mtr> </mtable> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>4</mn> <mo>)</mo> </mrow> </mrow>
其中，Vi,1≤i≤2是松弛变量。
4.根据权利要求1所述的基于变分光流估计肺部4D-CT图像超分辨率处理方法，其特征在于，
步骤S3采用改进的非局部迭代反投影算法重建出高分辨率图像具体包括：
设置高分辨率图像的初始估计对初始估计按照图像退化模型生成与集合相对应的低分辨率图像估计集通过反向投影重建误差{I-I(0)}到高分辨率估计来提高超分辨率重建效果；迭代此过程以最小化误差函数；
<mrow> <msup> <mi>e</mi> <mrow> <mo>(</mo> <mi>n</mi> <mo>)</mo> </mrow> </msup> <mo>=</mo> <msqrt> <mrow> <mfrac> <mn>1</mn> <mi>k</mi> </mfrac> <munderover> <mo>&amp;Sigma;</mo> <mrow> <mi>k</mi> <mo>=</mo> <mn>1</mn> </mrow> <mi>K</mi> </munderover> <mo>|</mo> <mo>|</mo> <msub> <mi>I</mi> <mi>k</mi> </msub> <mo>-</mo> <msubsup> <mi>I</mi> <mi>k</mi> <mrow> <mo>(</mo> <mn>0</mn> <mo>)</mo> </mrow> </msubsup> <mo>|</mo> <msubsup> <mo>|</mo> <mn>2</mn> <mn>2</mn> </msubsup> </mrow> </msqrt> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>12</mn> <mo>)</mo> </mrow> </mrow>
第n次迭代中，低分辨率图像估计的生成可通过下面的模型模拟实现：
其中，是指图像的仿射变换；h是点扩散函数；↓s为下采样算子；反投影重建过程可表示为：
<mrow> <msubsup> <mi>I</mi> <mi>h</mi> <mrow> <mo>(</mo> <mi>n</mi> <mo>+</mo> <mn>1</mn> <mo>)</mo> </mrow> </msubsup> <mo>=</mo> <msubsup> <mi>I</mi> <mi>h</mi> <mrow> <mo>(</mo> <mi>n</mi> <mo>)</mo> </mrow> </msubsup> <mo>+</mo> <mfrac> <mn>1</mn> <mi>K</mi> </mfrac> <munderover> <mo>&amp;Sigma;</mo> <mrow> <mi>k</mi> <mo>=</mo> <mn>1</mn> </mrow> <mi>K</mi> </munderover> <msubsup> <mi>T</mi> <mi>k</mi> <mrow> <mo>-</mo> <mn>1</mn> </mrow> </msubsup> <mrow> <mo>(</mo> <mo>(</mo> <mrow> <mrow> <mo>(</mo> <mrow> <msub> <mi>I</mi> <mi>k</mi> </msub> <mo>-</mo> <msubsup> <mi>I</mi> <mi>k</mi> <mrow> <mo>(</mo> <mi>n</mi> <mo>)</mo> </mrow> </msubsup> </mrow> <mo>)</mo> </mrow> <mo>&amp;UpArrow;</mo> <mi>s</mi> </mrow> <mo>)</mo> <mi>p</mi> <mo>)</mo> </mrow> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>14</mn> <mo>)</mo> </mrow> </mrow>
式中，是Tk的逆；p是反投影算子，p的取值将影响算法迭代时间；↑s表示上采样算子；
对IBP算法中的投影重建误差进行非局部均值滤波，保留图像中的高频细节；
利用数据加权平均法进行融合，根据低分辨率图像光流估计误差确定该图像的重建误差在反投影过程中的权重，抑制光流估计误差对重建结果的影响，反投影过程可表示为：
<mrow> <msubsup> <mi>I</mi> <mi>h</mi> <mrow> <mo>(</mo> <mi>n</mi> <mo>+</mo> <mn>1</mn> <mo>)</mo> </mrow> </msubsup> <mo>=</mo> <msubsup> <mi>I</mi> <mi>h</mi> <mrow> <mo>(</mo> <mi>n</mi> <mo>)</mo> </mrow> </msubsup> <mo>+</mo> <msub> <mi>&amp;omega;</mi> <mi>k</mi> </msub> <munderover> <mo>&amp;Sigma;</mo> <mrow> <mi>k</mi> <mo>=</mo> <mn>1</mn> </mrow> <mi>K</mi> </munderover> <msubsup> <mi>T</mi> <mi>k</mi> <mrow> <mo>-</mo> <mn>1</mn> </mrow> </msubsup> <mrow> <mo>(</mo> <mo>(</mo> <mrow> <mrow> <mo>(</mo> <mrow> <msub> <mi>I</mi> <mi>k</mi> </msub> <mo>-</mo> <msubsup> <mi>I</mi> <mi>k</mi> <mrow> <mo>(</mo> <mi>n</mi> <mo>)</mo> </mrow> </msubsup> </mrow> <mo>)</mo> </mrow> <mo>&amp;UpArrow;</mo> <mi>s</mi> </mrow> <mo>)</mo> <mi>p</mi> <mo>)</mo> </mrow> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>15</mn> <mo>)</mo> </mrow> </mrow>
将ωk设置为Ik与间欧式距离dk的指数函数，即：
<mrow> <msub> <mi>&amp;omega;</mi> <mi>k</mi> </msub> <mo>=</mo> <mfrac> <mn>1</mn> <mrow> <munderover> <mo>&amp;Sigma;</mo> <mrow> <mi>k</mi> <mo>=</mo> <mn>1</mn> </mrow> <mi>K</mi> </munderover> <msub> <mi>&amp;omega;</mi> <mi>k</mi> </msub> </mrow> </mfrac> <msup> <mi>e</mi> <mrow> <mo>-</mo> <mfrac> <msub> <mi>d</mi> <mi>k</mi> </msub> <mi>t</mi> </mfrac> </mrow> </msup> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>16</mn> <mo>)</mo> </mrow> </mrow>
其中t是控制核衰变速度的参数；
将图像所包含的所有相位图像作为已知的低分辨率图像集合，从中选出某一相位图像作为高分辨率图像。
5.根据权利要求3所述的基于变分光流估计肺部4D-CT图像超分辨率处理方法，其特征在于，
为了求解式(4)中的条件极值问题，将该约束问题转换成无约束的增广拉格朗日函数，如下式所示：
<mrow> <mtable> <mtr> <mtd> <mrow> <mi>&amp;Gamma;</mi> <mrow> <mo>(</mo> <mi>U</mi> <mo>,</mo> <msub> <mi>V</mi> <mi>i</mi> </msub> <mo>,</mo> <msub> <mi>W</mi> <mi>i</mi> </msub> <mo>)</mo> </mrow> <mo>=</mo> <mi>&amp;lambda;</mi> <mo>&amp;CenterDot;</mo> <munder> <mo>&amp;Sigma;</mo> <mrow> <mi>r</mi> <mi>e</mi> <mi>g</mi> <mi>i</mi> <mi>o</mi> <mi>n</mi> </mrow> </munder> <mi>b</mi> <mi>f</mi> <mi>w</mi> <mo>&amp;CenterDot;</mo> <mi>F</mi> <msup> <mrow> <mo>(</mo> <mi>U</mi> <mo>)</mo> </mrow> <mn>2</mn> </msup> <mo>+</mo> <munderover> <mo>&amp;Sigma;</mo> <mrow> <mi>i</mi> <mo>=</mo> <mn>1</mn> </mrow> <mn>2</mn> </munderover> <mi>D</mi> <mo>&amp;CenterDot;</mo> <mi>&amp;Psi;</mi> <mrow> <mo>(</mo> <msub> <mi>V</mi> <mi>i</mi> </msub> <mo>)</mo> </mrow> </mrow> </mtd> </mtr> <mtr> <mtd> <mrow> <mo>+</mo> <munderover> <mi>&amp;Sigma;</mi> <mrow> <mi>i</mi> <mo>=</mo> <mn>1</mn> </mrow> <mn>2</mn> </munderover> <mrow> <mo>(</mo> <mrow> <msub> <mi>W</mi> <mi>i</mi> </msub> <mo>&amp;CenterDot;</mo> <mrow> <mo>(</mo> <mrow> <msub> <mi>V</mi> <mi>i</mi> </msub> <mo>-</mo> <mi>U</mi> </mrow> <mo>)</mo> </mrow> </mrow> <mo>)</mo> </mrow> <mo>+</mo> <mfrac> <mi>&amp;rho;</mi> <mn>2</mn> </mfrac> <munderover> <mo>&amp;Sigma;</mo> <mrow> <mi>i</mi> <mo>=</mo> <mn>1</mn> </mrow> <mn>2</mn> </munderover> <mo>|</mo> <mo>|</mo> <msub> <mi>V</mi> <mi>i</mi> </msub> <mo>-</mo> <mi>U</mi> <mo>|</mo> <msubsup> <mo>|</mo> <mi>F</mi> <mn>2</mn> </msubsup> </mrow> </mtd> </mtr> </mtable> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>5</mn> <mo>)</mo> </mrow> </mrow>
应用ADMM算法得到下面的迭代求解步骤：
Step1：固定Vi和Wi，更新Uk+1：
<mrow> <mtable> <mtr> <mtd> <mrow> <msup> <mi>U</mi> <mrow> <mi>k</mi> <mo>+</mo> <mn>1</mn> </mrow> </msup> <mo>=</mo> <munder> <mi>argmin</mi> <mi>U</mi> </munder> <mi>&amp;lambda;</mi> <mo>&amp;CenterDot;</mo> <munder> <mi>&amp;Sigma;</mi> <mrow> <mi>r</mi> <mi>e</mi> <mi>g</mi> <mi>i</mi> <mi>o</mi> <mi>n</mi> </mrow> </munder> <mi>b</mi> <mi>f</mi> <mi>w</mi> <mo>&amp;CenterDot;</mo> <mi>F</mi> <msup> <mrow> <mo>(</mo> <mi>U</mi> <mo>)</mo> </mrow> <mn>2</mn> </msup> <mo>-</mo> <munderover> <mo>&amp;Sigma;</mo> <mrow> <mi>i</mi> <mo>=</mo> <mn>1</mn> </mrow> <mn>2</mn> </munderover> <msub> <mi>W</mi> <mi>i</mi> </msub> <mo>&amp;CenterDot;</mo> <mi>U</mi> </mrow> </mtd> </mtr> <mtr> <mtd> <mrow> <mo>+</mo> <mfrac> <mi>&amp;rho;</mi> <mn>2</mn> </mfrac> <munderover> <mo>&amp;Sigma;</mo> <mrow> <mi>i</mi> <mo>=</mo> <mn>1</mn> </mrow> <mn>2</mn> </munderover> <mo>|</mo> <mo>|</mo> <msub> <mi>V</mi> <mi>i</mi> </msub> <mo>-</mo> <mi>U</mi> <mo>|</mo> <msubsup> <mo>|</mo> <mi>F</mi> <mn>2</mn> </msubsup> </mrow> </mtd> </mtr> </mtable> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>6</mn> <mo>)</mo> </mrow> </mrow>
转换为：
<mrow> <msup> <mi>U</mi> <mrow> <mi>k</mi> <mo>+</mo> <mn>1</mn> </mrow> </msup> <mo>=</mo> <mfrac> <mrow> <mo>-</mo> <mn>2</mn> <mi>&amp;lambda;</mi> <munder> <mi>&amp;Sigma;</mi> <mrow> <mi>r</mi> <mi>e</mi> <mi>g</mi> <mi>i</mi> <mi>o</mi> <mi>n</mi> </mrow> </munder> <mi>b</mi> <mi>f</mi> <mi>w</mi> <mo>&amp;CenterDot;</mo> <mi>F</mi> <mrow> <mo>(</mo> <mi>U</mi> <mo>)</mo> </mrow> <mo>&amp;CenterDot;</mo> <mo>&amp;dtri;</mo> <msub> <mi>I</mi> <mrow> <mn>1</mn> <mi>w</mi> </mrow> </msub> <mo>+</mo> <munderover> <mi>&amp;Sigma;</mi> <mrow> <mi>i</mi> <mo>=</mo> <mn>1</mn> </mrow> <mn>2</mn> </munderover> <mrow> <mo>(</mo> <mrow> <msubsup> <mi>W</mi> <mi>i</mi> <mi>k</mi> </msubsup> <mo>+</mo> <msub> <mi>&amp;rho;V</mi> <mi>i</mi> </msub> </mrow> <mo>)</mo> </mrow> </mrow> <mrow> <mn>2</mn> <mi>&amp;rho;</mi> <mo>+</mo> <mn>2</mn> <mi>&amp;lambda;</mi> <munder> <mi>&amp;Sigma;</mi> <mrow> <mi>r</mi> <mi>e</mi> <mi>g</mi> <mi>i</mi> <mi>o</mi> <mi>n</mi> </mrow> </munder> <mi>b</mi> <mi>f</mi> <mi>w</mi> <mo>&amp;CenterDot;</mo> <mo>&amp;dtri;</mo> <msup> <msub> <mi>I</mi> <mrow> <mn>1</mn> <mi>w</mi> </mrow> </msub> <mn>2</mn> </msup> </mrow> </mfrac> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>7</mn> <mo>)</mo> </mrow> </mrow>
Step2：固定U和Wi，求解Vi,i＝1,2
<mrow> <mtable> <mtr> <mtd> <mrow> <mo>{</mo> <msub> <mi>V</mi> <mi>i</mi> </msub> <mo>}</mo> <mo>=</mo> <munder> <mi>argmin</mi> <mrow> <mo>{</mo> <msub> <mi>V</mi> <mi>i</mi> </msub> <mo>}</mo> </mrow> </munder> <mfrac> <mi>&amp;rho;</mi> <mn>2</mn> </mfrac> <munderover> <mo>&amp;Sigma;</mo> <mrow> <mi>i</mi> <mo>=</mo> <mn>1</mn> </mrow> <mn>2</mn> </munderover> <mo>|</mo> <mo>|</mo> <msub> <mi>V</mi> <mi>i</mi> </msub> <mo>-</mo> <msup> <mi>U</mi> <mrow> <mi>k</mi> <mo>+</mo> <mn>1</mn> </mrow> </msup> <mo>|</mo> <msubsup> <mo>|</mo> <mi>F</mi> <mn>2</mn> </msubsup> <mo>+</mo> <munderover> <mo>&amp;Sigma;</mo> <mrow> <mi>i</mi> <mo>=</mo> <mn>1</mn> </mrow> <mn>2</mn> </munderover> <mrow> <mo>(</mo> <msub> <mi>W</mi> <mi>i</mi> </msub> <mo>&amp;CenterDot;</mo> <msub> <mi>V</mi> <mi>i</mi> </msub> <mo>)</mo> </mrow> </mrow> </mtd> </mtr> <mtr> <mtd> <mrow> <mo>+</mo> <mi>&amp;lambda;</mi> <munderover> <mo>&amp;Sigma;</mo> <mrow> <mi>i</mi> <mo>=</mo> <mn>1</mn> </mrow> <mn>2</mn> </munderover> <mi>D</mi> <mo>&amp;CenterDot;</mo> <msub> <mi>&amp;Psi;</mi> <mi>i</mi> </msub> <mrow> <mo>(</mo> <msub> <mi>V</mi> <mi>i</mi> </msub> <mo>)</mo> </mrow> </mrow> </mtd> </mtr> </mtable> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>8</mn> <mo>)</mo> </mrow> </mrow>
其中，{Vi}表示集合{V1,V2}；计算Vi,i＝1,2，且等价写成下式
<mrow> <msubsup> <mi>V</mi> <mi>i</mi> <mrow> <mi>k</mi> <mo>+</mo> <mn>1</mn> </mrow> </msubsup> <mo>=</mo> <munder> <mi>argmin</mi> <msub> <mi>V</mi> <mi>i</mi> </msub> </munder> <mfrac> <mn>1</mn> <mn>2</mn> </mfrac> <mo>|</mo> <mo>|</mo> <msub> <mi>V</mi> <mi>i</mi> </msub> <mo>-</mo> <msub> <mi>&amp;Gamma;</mi> <mi>i</mi> </msub> <mo>|</mo> <msubsup> <mo>|</mo> <mi>F</mi> <mn>2</mn> </msubsup> <mo>+</mo> <mfrac> <mi>&amp;lambda;</mi> <mi>&amp;rho;</mi> </mfrac> <mi>D</mi> <mo>&amp;CenterDot;</mo> <msub> <mi>&amp;Psi;</mi> <mi>i</mi> </msub> <mrow> <mo>(</mo> <msub> <mi>V</mi> <mi>i</mi> </msub> <mo>)</mo> </mrow> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>9</mn> <mo>)</mo> </mrow> </mrow>
式中的式(9)中分解成规模更小的、并行求解的Lasso问题，即
<mrow> <msup> <mi>v</mi> <mrow> <mi>k</mi> <mo>+</mo> <mn>1</mn> </mrow> </msup> <mo>=</mo> <munder> <mrow> <mi>arg</mi> <mi>min</mi> </mrow> <mi>v</mi> </munder> <mfrac> <mn>1</mn> <mn>2</mn> </mfrac> <mo>|</mo> <mo>|</mo> <mi>v</mi> <mo>-</mo> <mi>t</mi> <mo>|</mo> <msup> <mo>|</mo> <mn>2</mn> </msup> <mo>+</mo> <mfrac> <mi>&amp;lambda;</mi> <mi>&amp;rho;</mi> </mfrac> <munderover> <mo>&amp;Sigma;</mo> <mrow> <mi>i</mi> <mo>=</mo> <mn>1</mn> </mrow> <mrow> <mi>k</mi> <mo>-</mo> <mn>1</mn> </mrow> </munderover> <msub> <mi>d</mi> <mi>i</mi> </msub> <mo>|</mo> <msub> <mi>v</mi> <mi>i</mi> </msub> <mo>-</mo> <msub> <mi>v</mi> <mrow> <mi>i</mi> <mo>+</mo> <mn>1</mn> </mrow> </msub> <mo>|</mo> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>10</mn> <mo>)</mo> </mrow> </mrow>
ν是张量Vi所包含的一维张量，t是Γi中与ν相对应的一维张量，k是ν的维数，di是D中与ν相对应一维张量的第i个元素，式(10)是一维全变分问题求解，并行计算的引入使得Vi的计算可以在毫秒时间内完成；
Step3：计算Wi,i＝1,2
<mrow> <msubsup> <mi>W</mi> <mi>i</mi> <mrow> <mi>k</mi> <mo>+</mo> <mn>1</mn> </mrow> </msubsup> <mo>=</mo> <msubsup> <mi>W</mi> <mi>i</mi> <mi>k</mi> </msubsup> <mo>+</mo> <mi>&amp;rho;</mi> <mrow> <mo>(</mo> <msubsup> <mi>V</mi> <mi>i</mi> <mrow> <mi>k</mi> <mo>+</mo> <mn>1</mn> </mrow> </msubsup> <mo>-</mo> <msubsup> <mi>U</mi> <mi>i</mi> <mrow> <mi>k</mi> <mo>+</mo> <mn>1</mn> </mrow> </msubsup> <mo>)</mo> </mrow> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>11</mn> <mo>)</mo> </mrow> </mrow>
该迭代求解过程在原始残差和对偶残差满足特定条件时结束。
CN201710344897.3A 2017-05-16 2017-05-16 一种基于变分光流估计肺部4d-ct图像超分辨率处理方法 Expired - Fee Related CN107221013B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201710344897.3A CN107221013B (zh) 2017-05-16 2017-05-16 一种基于变分光流估计肺部4d-ct图像超分辨率处理方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201710344897.3A CN107221013B (zh) 2017-05-16 2017-05-16 一种基于变分光流估计肺部4d-ct图像超分辨率处理方法

Publications (2)

Publication Number Publication Date
CN107221013A true CN107221013A (zh) 2017-09-29
CN107221013B CN107221013B (zh) 2020-07-17

Family

ID=59944142

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201710344897.3A Expired - Fee Related CN107221013B (zh) 2017-05-16 2017-05-16 一种基于变分光流估计肺部4d-ct图像超分辨率处理方法

Country Status (1)

Country Link
CN (1) CN107221013B (zh)

Cited By (9)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN108765290A (zh) * 2018-05-29 2018-11-06 天津大学 一种基于改进的稠密卷积神经网络的超分辨率重建方法
CN109448019A (zh) * 2018-09-27 2019-03-08 华南理工大学 一种变分光流模型平滑参数的自适应方法
CN109453473A (zh) * 2018-12-30 2019-03-12 上海联影医疗科技有限公司 放射治疗计划系统、确定装置及存储介质
CN109509149A (zh) * 2018-10-15 2019-03-22 天津大学 一种基于双通道卷积网络特征融合的超分辨率重建方法
CN110211193A (zh) * 2019-05-17 2019-09-06 山东财经大学 三维ct层间图像插值修复与超分辨处理方法及装置
CN110400357A (zh) * 2019-07-05 2019-11-01 北京航空航天大学 一种基于运动感知图像约束的4d-cbct重建方法
CN111724359A (zh) * 2020-06-12 2020-09-29 深圳技术大学 一种确定肺叶运动轨迹的方法、装置和存储介质
CN112016456A (zh) * 2020-08-27 2020-12-01 南京理工大学 基于自适应反向投影深度学习的视频超分辨率方法及系统
CN112634126A (zh) * 2020-12-22 2021-04-09 厦门美图之家科技有限公司 人像减龄处理方法、训练方法、装置、设备及存储介质

Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN103886568A (zh) * 2014-03-18 2014-06-25 南方医科大学 基于配准的肺4d-ct图像超分辨率重建方法
CN104299209A (zh) * 2014-09-18 2015-01-21 南方医科大学 基于快速亚像素运动估计的肺4d-ct图像超分辨率重建方法
CN106570862A (zh) * 2016-10-25 2017-04-19 中国人民解放军信息工程大学 一种超分辨率重建质量的评价方法及装置

Patent Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN103886568A (zh) * 2014-03-18 2014-06-25 南方医科大学 基于配准的肺4d-ct图像超分辨率重建方法
CN104299209A (zh) * 2014-09-18 2015-01-21 南方医科大学 基于快速亚像素运动估计的肺4d-ct图像超分辨率重建方法
CN106570862A (zh) * 2016-10-25 2017-04-19 中国人民解放军信息工程大学 一种超分辨率重建质量的评价方法及装置

Non-Patent Citations (4)

* Cited by examiner, † Cited by third party
Title
YU ZHANG 等: "Super-resolution reconstruction for 4D computed tomography of the lung via the projections onto convex sets approach", 《MEDICAL PHYSICS》 *
吴秀秀 等: "基于配准的肺4D-CT图像超分辨率重建", 《电子学报》 *
李秋萍 等: "基于改进的半监督FCM聚类算法的肺结节分类与识别", 《图学学报》 *
赵永洁: "基于光流场模型的非刚性医学图像配准算法研究", 《中国优秀硕士学位论文全文数据库 信息科技辑》 *

Cited By (12)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN108765290A (zh) * 2018-05-29 2018-11-06 天津大学 一种基于改进的稠密卷积神经网络的超分辨率重建方法
CN109448019A (zh) * 2018-09-27 2019-03-08 华南理工大学 一种变分光流模型平滑参数的自适应方法
CN109509149A (zh) * 2018-10-15 2019-03-22 天津大学 一种基于双通道卷积网络特征融合的超分辨率重建方法
CN109453473A (zh) * 2018-12-30 2019-03-12 上海联影医疗科技有限公司 放射治疗计划系统、确定装置及存储介质
CN110211193A (zh) * 2019-05-17 2019-09-06 山东财经大学 三维ct层间图像插值修复与超分辨处理方法及装置
CN110211193B (zh) * 2019-05-17 2020-08-04 山东财经大学 三维ct层间图像插值修复与超分辨处理方法及装置
CN110400357A (zh) * 2019-07-05 2019-11-01 北京航空航天大学 一种基于运动感知图像约束的4d-cbct重建方法
CN111724359A (zh) * 2020-06-12 2020-09-29 深圳技术大学 一种确定肺叶运动轨迹的方法、装置和存储介质
CN111724359B (zh) * 2020-06-12 2023-06-02 深圳技术大学 一种确定肺叶运动轨迹的方法、装置和存储介质
CN112016456A (zh) * 2020-08-27 2020-12-01 南京理工大学 基于自适应反向投影深度学习的视频超分辨率方法及系统
CN112016456B (zh) * 2020-08-27 2021-06-22 南京理工大学 基于自适应反向投影深度学习的视频超分辨率方法及系统
CN112634126A (zh) * 2020-12-22 2021-04-09 厦门美图之家科技有限公司 人像减龄处理方法、训练方法、装置、设备及存储介质

Also Published As

Publication number Publication date
CN107221013B (zh) 2020-07-17

Similar Documents

Publication Publication Date Title
CN107221013B (zh) 一种基于变分光流估计肺部4d-ct图像超分辨率处理方法
Chen et al. LEARN: Learned experts’ assessment-based reconstruction network for sparse-data CT
Trinh et al. Novel example-based method for super-resolution and denoising of medical images
Kudo et al. Virtual thin slice: 3D conditional GAN-based super-resolution for CT slice interval
JP4885138B2 (ja) 一連の画像における動き修正のための方法およびシステム
Woo et al. Reconstruction of high-resolution tongue volumes from MRI
Isola et al. Fully automatic nonrigid registration‐based local motion estimation for motion‐corrected iterative cardiac CT reconstruction
Wang et al. Fast image upsampling via the displacement field
Gu et al. A contour-guided deformable image registration algorithm for adaptive radiotherapy
EP3486864B1 (en) Isotropic 3d image reconstruction using 3d patches-based self-similarity learning
CN110211193B (zh) 三维ct层间图像插值修复与超分辨处理方法及装置
Bruder et al. Adaptive iterative reconstruction
CN111724307B (zh) 一种基于最大后验概率和非局部低秩先验的图像超分辨重建方法，终端及可读存储介质
Zhang et al. MR image super-resolution reconstruction using sparse representation, nonlocal similarity and sparse derivative prior
Müller et al. Image artefact propagation in motion estimation and reconstruction in interventional cardiac C-arm CT
CN109544478A (zh) 一种基于奇异值分解的非局部均值ct图像降噪方法
Dillon et al. Evaluating reconstruction algorithms for respiratory motion guided acquisition
Guan et al. Generative modeling in sinogram domain for sparse-view CT reconstruction
Li et al. 3D coronary artery reconstruction by 2D motion compensation based on mutual information
CN104299209A (zh) 基于快速亚像素运动估计的肺4d-ct图像超分辨率重建方法
Isola et al. Cardiac motion-corrected iterative cone-beam CT reconstruction using a semi-automatic minimum cost path-based coronary centerline extraction
Huang et al. Mambamir: An arbitrary-masked mamba for joint medical image reconstruction and uncertainty estimation
Cierniak et al. A practical statistical approach to the reconstruction problem using a single slice rebinning method
Nagare et al. A bias-reducing loss function for CT image denoising
Barbano et al. Steerable conditional diffusion for out-of-distribution adaptation in imaging inverse problems

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant
CF01 Termination of patent right due to non-payment of annual fee
CF01 Termination of patent right due to non-payment of annual fee

Granted publication date: 20200717