电解剖图EAM'即是粗配准后的CART0电解剖图EAM*。
[0036] 在本发明实施例中，CART0电解剖图EAM'到所述CT图像的距离可以定义为所有 CART0标测点到CT图像最近点欧式距离的均方根。因为仅需遍历48种主轴的对应关系，寻 找其中使得两图距离最小的一种作为最终的配准结果，所以粗配准算法中不涉及任何优化 搜索算法。
[0037]S13、在粗配准的基础上，以距离平方和为相似性度量函数，采用增强自由形变模 型和金字塔分层方法，精配准CART0电解剖图与CT图像，得到形变后的CT图像；精配准过 程中，CART0电解剖图为基准图像，CT图像为浮动图像：
[0038] 需要说明的是，以自由形变模型（FreeFormDeformation,FFD)配准CART0电解 剖图与CT图像的基本思路是：将由ΜΧΝΧ0个控制点P"^。构成的控制网格覆盖在CT图像 所在空间，调整形变参数，操控由控制点组成的控制网格，实现对CT图像形状的改变。
[0039]然而，FFD将每个控制点的最大位移约束在两个控制点间隔的一半以内，所以其形 变自由度较低，能描述的形变量有限。为了实现更为精确细致的CART0电解剖图与CT图像 配准，本发明采用增强自由形变模型（IncrementalFreeFormDeformation，IFFD)，增强 自由形变模型为每个控制点P"^。增加了形变量S进一步，本发明采用了基于B样条 的局部增强自由形变模型。
[0040] 如此，每个控制点的空间位置满足以下关系：
[0041] Ρη,η,ο=Ρη,η,ο〇+δΡη,η,0 (2)
[0042] 其中，Ρ&是控制点Ρ未发生形变时的空间位置， {d，。，spy一，δρζ…)}。
[0043] 依据上述定义，给定控制点从PminiC1°形变到Ρ。，形变量为δΡ^。， 控制网格覆盖下的CT图像各点Y= (x，y，z)亦随之发生形变，形变量为
点Y的最终空间位置为L(Y) =Y+δL(Y)。
[0044] 这里，
[0045]
[0046] U表示向下取整，Bk(u)表示B样条的第k个基函数，δPi+kij+1_((k，1 ，r)e[0, 3]X[0, 3]X[0, 3])是覆盖在CT图像点Y= (x，y，z)上的相邻64个控制点的形 变量。
[0047] 从上所述，CT图像中任意点Y的位置，完全由覆盖其上的相邻64个控制点的形变 量决定，改变即可实现CT图像和CART0电解剖图之间的配准。因此，基于 IFFD的CART0电解剖图与CT图像配准，即等同于寻找控制网格的形变量δΡ^。，使得形变 后的CT图像与CART0电解剖图空间位置达到一致。
[0048] 本发明实施例提供了步骤S13的具体实现方式。图3是本发明实施例中精配准 CART0电解剖图与CT图像的方法流程示意图，如图3所示，该精配准方法包括以下步骤：
[0049] S131、设定控制网格分辨率ΜΧΝΧ0和金字塔层数s的初始值，令Μ=Ν= 0 = 3， s= 1 ：
[0050] 为了减少弹性形变产生的大计算量，本发明实施例采用了金字塔分层技术，随着 金字塔层数的增加，其控制网格的分辨率逐渐增加。
[0051] S132、对CART0电解剖图每一个点X，运用迭代最近点方法寻找CT图像上 与之对应的点L( δ Ρηη,。;X);其中，CART0电解剖图位于空间Ω中，Ω={X}= {(x,y,ζ) 11 ^X^X, 1 ^y^Υ, 1 ^ζ^Ζ}；
[0052] S133、构建由ΜΧΝΧ0个控制点均匀分布构成的控制网格P，将控制网格P覆盖在 CT图像空间Γ之上：
[0053] 所述控制网格P具体描述如下：
[0054] p= {Pn,n,〇} = {(Pxn,n,〇,Pyn,n,〇,Pzn,n,〇)}((m,n,o)e[ 1,M]X[ 1,N]X[ 1, 0])
[0055] 本发明实施例采用增强自由形变模型，每个控制点的空间位置满足以下关系：
[0056] pn, n,。=pn,n,。〇+δPn, n,。 (4)
[0057]P"^。为控制网格的控制点发生形变后的空间位置；
[0058]P^。·"是控制点P"^。未发生形变时的空间位置；
[0059] 6?|^。是控制点？|^。的形变量，5?_={〇?:。，5?:,。，5?、。)};令控制 点的初始形变量δPm,niC]= 〇;
[0060]S134、以梯度下降法优化目标函数，获得每个控制点的形变量δΡ^。，然后 确定形变后控制点的位置：
[0061] 本发明的目标函数的表达式如下：
[0062] Ε(δΡηηο) =Esinilarity(δΡΠιn,0) +aEsnoothness (δΡηηο) (5)
[0063] Esinularity (δΡ^。）为相似性度量函数，本发明以距离平方和度量相似性时，该相似 性度量函数的表达式为：
[0065]X表示CART0电解剖图中的点，L(δΡη, η,。;X)表示CT图像中与X相对应的点。 [0066] 为平滑约束项，用于避免出现形状扭曲，保持配准网格的规整 性。 ^smoothness 的计算公式为：
[0068] α是平衡Esl_wty和Es_th_s两项的一个常量，α越小，收敛越快；α越大，则形 变区域越平滑。
[0069] 图4进一步显示了S134的具体过程如下：
[0070]S1341、对控制网格上的每一个控制点的形变量δΡ^。，根据公式（8)计算目标函 数的梯度▽￡;
[0074] S1342、判断是否成立，如果否，按照公式（10)更新该控制点形变量 3?|^。，并返回执行步骤31341;如果是，则5?|^。为该控制点？ |^。的形变量5?|1^。，输出δΡη,,η,ο；
[0076] 其中，μ表示每一步沿梯度向量迭代的步幅，ε表示阈值。
[0077]S135、由形变后的控制点位置，确定覆盖其下的CT图像中各点的位置。
[0078] 图5显示S135的具体过程如下：
[0079]S1351、根据所有控制点Ρ"^的形变量^，由Β样条的自由形变模型计算得出 CT图像上每一个点Y= (x，y，z)的形变量5L(Y)，计算公式如下
[0081]Bk(u)表示B样条的第k个基函数，SPi+kj+Wr((k，l，r)e[0,3]X[0,3]X[0,3]) 是覆盖在CT图像点Y= (x，y，z)上的相邻64个控制点的形变量；
[0082]S1352、根据CT图像上每一个点Y= (X，y，z)的形变量δL(Y)，确定形变后CT图 像各点Υ的最终空间位置L(Υ)，其中，L(Y) =Υ+δL(Υ)。
[0083]S136、更新控制网格的分辨率以及金字塔层数，令M=M+1，N=N+1，0 = 0+l，s= s+1〇
[0084]S137、循环执行步骤S132至步骤S136,直至s=s阈值：
[0085] 需要说明的是，在本发明实施例中，可以根据实际需要选取。选择需要考 虑计算量和配准精度之间的平衡，一般情况下，选择4或5即可。
[0086]S14、将精配准形变后的CT图像叠加在CART0电解剖图上，得到CART0电解剖图与 CT图像的最终配准图像。
[0087] 以上即为本发明实施例提供的基于增强弹性形变的CART0电解剖图与CT图像的 配准方法的【具体实施方式】。在该配准方法中，采用增强自由形变模型和金字塔分层技术，可 以对或大或小的弹性形变建模，构建的形变场平滑、连续，能够保持拓扑结构和一一对应关 系，能更为精确细致地描述CART0电解剖图与CT图像之间的形变差异，从而大幅提高CART0 电解剖图与CT图像的配准精度；同时，金字塔分层技术能减少因弹性形变产生的大的计算 量，在保证配准精度的同时，也能获得较快的配准速度。通过本发明提供的配准方法得到的