X射线造影图像序列中的心动周期恢复方法

摘要

本发明提出了一种X射线造影图像序列中的心动周期恢复方法，解决了临床获取覆盖多个心动周期或一个心动周期多个成像角度的造影图像序列后，需要结合心动周期获取同步的图像对的问题。步骤一、基于图像直方图能量相似性测度的特征点跟踪模块通过EHD相似性测度跟踪序列图像中的特征点；步骤二、基于特征点位置分布的矢量场分解模块根据所述特征点在图像中的位置分布，实现所述位移矢量场的分解；步骤三、心脏运动的三维简谐运动模型建立模块根据实际心脏运动形式的分析，完成三维简谐运动模型的构建，利用凸包算法实现简谐运动中心的估计；步骤四、心动周期的重建模块基于心脏运动矢量的估计，同步图像帧和心动周期所有相位，完成心动周期的重建。

说明书

技术领域

本发明涉及一种X射线造影图像序列中的心动周期恢复方法，可广泛应用 于临床研究中减影图像的配准和三维图像重建的过程中。

背景技术

心动周期对于同步时序造影图像具有重要的作用。为实现该目的，心电图 通常用于临床上心电信息的记录。然而，心电门控配置并不是传统单平面或双 平面造影系统的标准配置，为实现相关目标，额外的设备必带给临床医生一定 的麻烦，故重建心动周期为临床操作带来巨大的便利。

冠状动脉分布在心脏表面，其运动可间接反映心脏的基本运动，分析冠状 动脉的运动可作为心脏运动恢复、冠脉后续研究乃至相关疾病诊断的基础。冠 状动脉血管造影是临床冠脉疾病诊断和治疗常用的成像方式之一，其高分辨率 和高对比度使其成为心血管疾病诊断和治疗所用成像方式的“金标准”，在时 序造影图像中，包含大量血管的空间、结构及运动信息，通过分析不同时刻造 影图像中的冠脉树，对图像背景或冠脉上特征点进行提取及跟踪，获得相关运 动特性，建立精确的数学模型，可获得以运动场、形状变化等方式表达的整体 或局部信息，进而对冠脉运动进行定量分析。故在过去的二十年，基于造影图 像序列的冠脉运动分析成为国内外心脏病研究的焦点。

目前，基于二维造影图像的运动分析主要包括两个方向：一种是基于若干 幅不同时刻的造影图像进行运动补偿，插值获取连续时序的造影图像。在二维 造影图像中，因血管运动和变形幅度较大，临床应用中会丢失大量细节信息， 如狭窄部分，通过基于不同时刻血管的运动分析，在心脏舒缩方向平移像素单 位，实现对序列图像的运动补偿，获取连续的时序二维造影图像，进而完成造 影图像中血管运动特征的二维描述。另一种是在时序造影图像基础上，获取相 同时刻不同角度的造影图像，为三维重建提供数据基础。在临床应用中，单平 面造影系统及双平面造影系统均可根据心电触发信号(ECG)触发获取配准的 造影图像对，实现减影图像的配准或图像的三维重建。但前者存在精确度低的 缺点，后者临床费用较高，二者在实际应用中受到严重的限制。故通过对单平 面造影系统获取的二维造影图像的运动进行分析，恢复心电信号，逆序获取不 同角度相同时刻的造影图像，可综合解决上述两种问题。

现存心动周期重建存在的问题：

1.研究所用心动周期通过临床获取，设备复杂；

2.心动周期重建的数据源于超声图像序列，分辨率低；

3.心动周期的重建因未忽略呼吸运动的干扰，精度低。

发明内容

在临床获取覆盖多个心动周期或一个心动周期多个成像角度的造影图像 序列后，需要结合心动周期获取同步的图像对，但传统的造影图像系统中未合 并ECG触发获取影像的相关设备，导致图像对的获取过程的繁琐。本发明针对 此问题，提出了一种X射线造影图像序列中的心动周期恢复方法。

一种X射线造影图像序列中的心动周期恢复方法，包括如下步骤：

步骤一、基于图像直方图能量相似性测度的特征点跟踪模块通过EHD相似 性测度跟踪序列图像中的特征点，实现覆盖多个心动周期的序列图像中特征点 的位移矢量场的估计；

步骤二、基于特征点位置分布的矢量场分解模块根据所述特征点在图像中 的位置分布，实现所述位移矢量场的分解，实现冠状动脉和背景运动场的描述；

步骤三、心脏运动的三维简谐运动模型建立模块根据实际心脏运动形式的 分析，完成三维简谐运动模型的构建，利用凸包算法实现简谐运动中心的估计；

步骤四、心动周期的重建模块基于心脏运动矢量的估计，同步图像帧和心 动周期所有相位，完成心动周期的重建。

本发明的有益效果：

1、利用最大化EHD相似性测度获取造影序列中的特征点，全过程无需附加工 作指定图像中的相关结构，特征点的跟踪方法具有普适性；

2、根据特征点的位置，基于呼吸运动的描述可去除其对心脏运动分析的干扰， 精度高；

3、采用凸包算法获取运动模型的中心，运算量小。

附图说明

图1是本发明的模块图；

图2是本发明的方法流程图；

具体实施方式

下面结合附图详细说明本发明，具体包括以下几个步骤：

步骤1，对造影图像序列进行运动分析。首先，利用EHD相似性测度、Powell 参数优化算法和双三次图像插值，实现特征点的跟踪，完成特征点位移矢量场 的估计。

特征点的跟踪可转化为寻找最优变换T_p。如果用F_i,k(x,y)和L_i,k+1(x,y)表示 第k帧和第k+1帧图像中第i个特征点所对应子块图像中像素点的灰度值，则 最优变换可认为使目标函数f(T_p(F_i,k),L_i,k+1)达到其最大值或最小值，其中p表 示一系列变换参数。

则跟踪算法可用下列式子表示：

I＝function_max/min(T_p(F_i,k),L_i,k+1)

作为初期研究，搜寻空间被局限于二维旋转和平移。因此，最优变换Tp可 表示如下：

$T_P(x,y)=\left(\begin{array}{ccc}\cos \left(\theta \right)& \sin \left(\theta \right)& t_x\\ -\sin \left(\theta \right)& \cos \left(\theta \right)& t_y\\ 0& 0& 1\end{array}\right)\left(\begin{array}{c}x\\ y\\ 1\end{array}\right)$

在寻求最优变化的过程中，跟踪过程可分为以下三个重要步骤：(1)相似 性测度，(2)最优化方法，(3)插值方法。

EHD定义如下：

$\mathrm{EHD}(F,L)=\underset{g=-G}{\overset{G}{\Sigma }}{H}^2\left(g\right)$

其中H(g)表示相邻时序帧图像灰度差的归一化直方图分布，表示如下：

$H\left(g\right)=\frac{1}{N}\underset{x}{\Sigma }\underset{y}{\Sigma }\delta (L_{i,k+1}(x,y)-F_{i,k}(x,y),g)$

其中表示图像灰度差，δ(x,y)是二维的狄拉克函数。

考虑到子图像变形幅度，本发明采用双三次插值保证图像的平滑性和连通 性。考虑到优化算法的计算复杂度，本发明采用Powell参数优化算法优化最 优变换中的变换参数。

在本发明中，有关图像中特征点的提取分为以下两类：第一帧图像中通过 手动获取特征点，在后续的图像中，利用上面介绍的方法进行特征点的跟踪。

其次，完成所有帧图像中特征点的跟踪之后，根据特征点在帧图像中的像 素坐标，重建图像中相关特征点的位移场。公式表达如下：

${\overrightarrow{S}}_{i,k}^{\mathrm{feature}}=X_{i,k+1}^{\mathrm{feature}}-X_{i,k}^{\mathrm{feature}}$

步骤2，简谐运动模型的建立。首先，按照图像中特征点的位置分布，将 步骤1中所获取的位移场分解为可描述冠状动脉和背景运动的位移场。公式描 述如下：

${\overrightarrow{S}}_{\mathrm{coronaryartery}}={\overrightarrow{S}}_{\mathrm{heart}}+{\overrightarrow{S}}_{\mathrm{respiratory}}$

${\overrightarrow{S}}_{\mathrm{background}}={\overrightarrow{S}}_{\mathrm{respiratory}}$

${\overrightarrow{S}}_{\mathrm{coronrayartery}}-{\overrightarrow{S}}_{\mathrm{background}}={\overrightarrow{S}}_{\mathrm{heart}}$

则有关心脏运动的速度矢量场可描述如下：

${\overrightarrow{V}}_{\mathrm{heart}}={\overrightarrow{S}}_{\mathrm{heart}}/\mathrm{\Delta t}$

$\left|{\overrightarrow{V}}_{\mathrm{heart}}\right|={[{\left(V_{\mathrm{heart}}^x\right)}^2+{\left(V_{\mathrm{heart}}^y\right)}^2]}^{1/2}$

然后，综合实际心脏运动的详细描述和估计的心脏速度矢量场，构建描述 心脏运动的简谐运动模型，并利用凸包算法，结合心脏的速度矢量，即可获得 运动模型的中心坐标。

由于冠状动脉分布在心脏的表面，通过基于背景特征点的分析，即可去除 冠脉运动中因呼吸引起的部分，剩余的运动矢量即为心脏运动所引起的，此时 冠脉的运动矢量即可反映心脏的运动。心脏运动速度描述如下：

本发明假设心脏在运动时只存在收缩和舒张两种运动模式，因此心脏表面 的特征点的运动即可被认为是以心脏的运动中心为中心的简谐运动，此种描述 即为本专利所构建的心脏的简谐运动模型。

$V_k^{\mathrm{systole}/\mathrm{diastole}}=\underset{i-1}{\overset{n}{\Sigma }}\frac{{\overrightarrow{V}}_{i,k}^{\mathrm{heart}}·\overrightarrow{{\mathrm{OA}}_{i,k}}}{\left|\overrightarrow{{\mathrm{OA}}_{i,k}}\right|}$

完成模型的构建后，求取冠脉因心脏运动引起的速度矢量的交点，利用凸 包算法获取交点所对应的最大凸多边形，此时凸多边形的重心即可认为是心脏 运动的中心。

步骤3，心动周期的重建。通过对实际心动周期每个相位心肌运动状况的 分析，依据表1的分析结果，完成心电相位和图像帧之间的对应，最终实现心 动周期的重建。

表1心动相位运动描述

心动阶段 图像中特征点运动规律 等容收缩期 图像中特征点收缩运动，但速度较小。 快速射血期 图像中特征点快速收缩运动。 减慢射血期 图像中特征点继续收缩运动，但速度放缓。 等容舒张期 图像中特征点由收缩运动转变为舒张运动。 快速充盈期 图像中特征点快速舒张运动。 减慢充盈期 图像中特征点继续舒张运动，但速度放缓 心房收缩期 图像中特征点由舒张运动转变为收缩运动。