股骨三维模型可视化方法技术领域
本发明涉及图像处理领域。
背景技术
股骨是一种典型的人体长骨,具有不规则的空间结构,由两端的松质骨和中间的管状
皮质骨构成。股骨三维模型的可视化,对于股骨骨折的术前诊断和手术计划的制定都具有
极其重要的意义。
目前国内骨折的影像显示主要采用直接数字X光图像DR(Digital Radiography),
但该方法最大的缺陷在于二维投影成像,导致影像重叠,丢失了大量三维的空间信息,而
这些三维信息对于医生的临床诊断恰恰是十分重要的;另一方面,现有的医学三维可视化
技术如计算机断层扫描成像CT(Computer Tomography)虽然可以获得研究对象的精确三
维结构形态,但其体积庞大,不适于在突发情况下的骨折影像显示。
发明内容
本发明目的是为了解决现有采用直接数字X光图像DR无法获取三维空间信息,而
采用CT由于体积庞大不适于在突发情况下的股骨影像显示的问题,提供了一种股骨三维
模型可视化方法。
本发明所述股骨三维模型可视化方法,它包括以下步骤:
步骤一、建立通用股骨三维模型,通过通用股骨三维模型的正位投影获取通用股骨的
正位股骨三维模型投影轮廓,通过通用股骨三维模型的侧位投影获取通用股骨的侧位股骨
三维模型投影轮廓;
步骤二、对待建模的股骨正位DR图像和股骨侧位DR图像进行预处理,提取股骨正
位DR图像中的股骨边缘轮廓获得正位股骨边缘轮廓,提取股骨侧位DR图像中的股骨边
缘轮廓获得侧位股骨边缘轮廓;
步骤三、对步骤二提取的待建模的股骨正位DR图像的边缘轮廓与步骤一所述的通用
正位股骨三维模型投影轮廓进行仿射配准,确定正位图像的配准变换关系;
对步骤二提取的待建模的股骨侧位DR图像的边缘轮廓与步骤一所述的通用侧位股
骨三维模型投影轮廓进行仿射配准,确定侧位图像的配准变换关系;
步骤四、根据步骤三获得的正位图像的配准变换关系和侧位图像的配准变换关系确定
股骨三维模型的姿态,实现股骨的三维模型的可视化。
本发明的优点:相对于已有的医学三维可视化方法,本发明利用双平面DR图像实现
股骨的三维模型可视化,成本低,操作简单,有害辐射小,重建图像比较迅速,能为股骨
粉碎性骨折手术的诊断和手术计划的制定提供可靠的三维影像学信息。
采用本发明所述方法使骨折断端移位情况显示得更加立体直观,医生不必再凭借经验
估计骨折的类型和损伤的情况,减少了漏检和误诊的发生率,对术前手术计划的制定和手
术方式的选择有重要意义。
附图说明
图1为本发明的框图示意图,图2至图5为具体实施方式四中能量泛函与边界曲线关
系的示意图,其中图2是边界曲线包含边缘轮廓的示意图,图3是边界曲线在边缘轮廓内
部的示意图,图4是边界曲线穿过边缘轮廓的示意图,图5是边界曲线对应边缘轮廓的示
意图,图6为具体实施方式二中三维点与其在二维图像上的投影点之间变换关系的示意
图,图7为具体实施方式二中通用股骨的三维模型示意图,图8为具体实施方式五中仿射
配准前待建模的股骨正位DR图像轮廓和通用股骨正位投影轮廓的示意图,图中粗线表示
待建模的股骨正位DR图像轮廓,细线表示通用股骨正位投影轮廓;图9为具体实施方式
五中仿射配准后待建模的股骨正位DR图像轮廓和通用股骨正位投影轮廓的示意图,图
10是三维姿态估计后的股骨三维模型示意图。
具体实施方式
具体实施方式一:结合图1说明本实施方式,本实施方式所述股骨三维模型可视化
方法,它包括以下步骤:
步骤一、建立通用股骨三维模型,通过通用股骨三维模型的正位投影获取通用股骨的
正位股骨三维模型投影轮廓,通过通用股骨三维模型的侧位投影获取通用股骨的侧位股骨
三维模型投影轮廓;
步骤二、对待建模的股骨正位DR图像和股骨侧位DR图像进行预处理,提取股骨正
位DR图像中的股骨边缘轮廓获得正位股骨边缘轮廓,提取股骨侧位DR图像中的股骨边
缘轮廓获得侧位股骨边缘轮廓;
步骤三、对步骤二提取的待建模的股骨正位DR图像的边缘轮廓与步骤一所述的通用
正位股骨三维模型投影轮廓进行仿射配准,确定正位图像的配准变换关系;
对步骤二提取的待建模的股骨侧位DR图像的边缘轮廓与步骤一所述的通用侧位股
骨三维模型投影轮廓进行仿射配准,确定侧位图像的配准变换关系;
步骤四、根据步骤三获得的正位图像的配准变换关系和侧位图像的配准变换关系确定
股骨三维模型的姿态,实现股骨的三维模型的可视化。
具体实施方式二:结合图6和图7说明本实施方式,本实施方式是对实施方式一所
述股骨三维模型可视化方法的进一步限定,步骤一中建立通用股骨三维模型,通过通用股
骨三维模型进行正侧位投影获取通用股骨三维模型投影轮廓的过程为:
步骤一一、选取一个通用股骨对其进行CT扫描,获取通用股骨的CT片,根据通用
股骨的CT片构建出通用股骨三维模型;
步骤一二、对通用股骨三维模型进行正位投影获取正位通用股骨三维模型的投影轮
廓,对通用股骨三维模型进行侧位投影获取侧位通用股骨三维模型的投影轮廓,具体为:
三维点与其在二维图像上的投影点之间有如下变换关系:
u
v
1
=
A
1
A
2
X
w
Y
w
Z
w
1
-
-
-
(
1
)
]]>
其中(Xw,Yw,Zw)是点P在世界坐标系下的三维坐标;(u,v)是该点二维像素坐标系下的
二维投影坐标;A1为内参数矩阵,描述了相机几何及光学特性参数;A2为外参数矩阵,
描述了相机在世界坐标系的位置和姿态。
这里设模板平面在世界坐标系的Zw=0的平面上,不失一般性。则(1)式能改写为:
u
v
1
=
A
1
r
1
r
2
r
3
t
X
w
Y
w
0
1
=
A
1
r
1
r
2
t
X
w
Y
w
1
-
-
-
(
2
)
]]>
其中r1,r2,r3,t分别是矩阵A2的列向量,即A2=[r1 r2 r3 t]。由于每幅图像都能确定
一个3*3的单应性矩阵H,该单应性矩阵能通过棋盘标定法计算来确定,并根据单应性矩
阵的性质,有如下关系式:
H
=
h
1
h
2
h
3
=
η
A
1
r
1
r
2
t
-
-
-
(
3
)
]]>
其中η为任意标量。旋转矩阵的性质具有如下两个性质:
r
1
T
r
2
=
0
|
|
r
1
|
|
=
|
|
r
2
|
|
=
1
-
-
-
(
4
)
]]>
由(3)式可得:
r
1
=
η
-
1
A
1
-
1
h
1
r
2
=
η
-
1
A
1
-
1
h
2
-
-
-
(
5
)
]]>
将(5)式代入(4)式中,每个图像就能获得以下两个对内参数矩阵的约束:
h
1
T
A
1
-
T
A
1
-
1
h
2
=
0
h
1
T
A
1
-
T
A
1
-
1
h
1
=
h
2
T
A
1
-
T
A
1
-
1
h
2
-
-
-
(
6
)
]]>
由公式(6)即可求解出内参数矩阵A1。再根据A1和单应性矩阵H能进一步计算出
外参数r1,r2,r3,t:
r
1
=
η
-
1
A
1
-
1
h
1
r
2
=
η
-
1
A
1
-
1
h
2
r
3
=
r
1
×
r
2
t
=
η
-
1
A
1
-
1
h
3
-
-
-
(
7
)
;
]]>
将公式(6)和(7)得到的参数代入公式(2),得到三维点在二维图像上对应的投影点;
对通用股骨三维模型进行正位投影,所述正位投影的外边缘点构成正位股骨三维模型
的投影轮廓;
对通用股骨三维模型进行侧位投影,所述侧位投影的外边缘点构成侧位股骨三维模型
的投影轮廓。
本发明对棋盘标定法进行了改进,采用覆细铜线网格构成的PCB板取代经典的黑白
格棋盘作为标定板,以适用于X射线下的拍摄环境,利用以上标定的参数,分别对通用
模型进行正位投影和侧位投影,获得投影轮廓,完成对通用股骨三维模型投影轮廓的获取。
具体实施方式三:本实施方式是对实施方式一所述股骨三维模型可视化方法的进一步
限定,步骤二中对待建模的股骨正位DR图像和股骨侧位DR图像进行预处理的过程为:
对待建模的股骨正位DR图像和股骨侧位DR图像进行中值滤波。
利用中值滤波方法对股骨正位DR图像和股骨侧位DR进行图像增强预处理,达到消
减噪声、突出解剖结构特征的目的。
具体实施方式四:结合图2至图5说明本实施方式,本实施方式是对实施方式一所
述股骨三维模型可视化方法的进一步限定,步骤二中提取待建模的正位股骨边缘轮廓和待
建模的侧位股骨边缘轮廓的过程为:
步骤二一、C为待提取的正位股骨边缘轮廓曲线,初始化该闭合曲线C=C0,所述
闭合曲线置于股骨正位DR图像或股骨侧位DR图像上,C0为任意形状的闭合曲线;
采用无边缘活动轮廓模型在股骨正位DR图像中进行全局边缘检测,无边缘活动轮廓
模型定义最小化能量泛函为:
其中,c1是曲线内部的平均灰度,c2是曲线外部的平均灰度,Fin表示曲线内部区域
的灰度值与c1的平方误差,Fout表示曲线外部区域的灰度值与c2的平方误差,I表示DR
图像,Length(C)表示边界曲线的长度,μ、λ1和λ2是权重系数,其中,μ≥0;
能量泛函与边界曲线具有如下关系:
(一)当Fin>0,Fout=0,F(C,c1,c2)=Fin+Fout>0时,边界曲线包含边缘轮廓;见图2;
(二)当Fin=0,Fout>0,F(C,c1,c2)=Fin+Fout>0时,边界曲线在边缘轮廓内部;见
图3;
(三)当Fin>0,Fout>0,F(C,c1,c2)=Fin+Fout>0时,边界曲线穿过边缘轮廓;见图4;
(四)当Fin=0,Fout=0,F(C,c1,c2)=Fin+Fout=0时,边界曲线对应边缘轮廓;见图5;
仅当能量函数达到最小值时,对应的边界曲线为图像的边缘轮廓;
步骤二二、采用变分水平集的方法,将公式(8)中引入Heaviside函数,将它修改为
水平集函数u的泛函得到公式(9):
F
(
u
,
c
1
,
c
2
)
=
λ
1
∫
∫
|
I
-
c
1
|
2
H
(
u
)
dxdy
+
λ
2
∫
∫
|
I
-
c
2
|
2
(
1
-
H
(
u
)
)
dxdy
+
μ
∫
∫
δ
(
u
)
|
▿
u
|
dxdy
-
-
-
(
9
)
]]>
其中,u(x,y)为水平集函数,与曲线存在如下关系:
u
(
x
,
y
)
<
0
(
x
,
y
)
inside
(
C
)
u
(
x
,
y
)
=
0
(
x
,
y
)
on
(
C
)
u
(
x
,
y
)
>
0
(
x
,
y
)
outside
(
C
)
-
-
-
(
10
)
]]>
Heaviside函数表达式为:
H
(
z
)
=
1
,
z
≥
0
0
,
z
<
0
-
-
-
(
11
)
]]>
δ
(
z
)
=
dH
(
z
)
dz
-
-
-
(
12
)
]]>
将公式(8)关于C曲线环路积分项改写为面积分的形式:
将公式(10)、(11)、(12)代入公式(8)能得到公式(9);
步骤二三、利用梯度下降流方法求解公式(9)的最小值,得公式:
∂
u
∂
t
=
δ
(
u
)
[
μ
·
div
(
▿
u
|
▿
u
|
)
-
λ
1
(
1
-
c
1
)
2
+
λ
2
(
1
-
c
2
)
2
]
-
-
-
(
14
)
]]>
其中,按照差分的形式能表示为
其中,τ表示时间步长,设为0.1;
步骤二四、将初始轮廓值C0=u0(x,y)代入公式(14),设定迭代次数,求解出此时刻
的水平集函数u,根据公式(10),令u(x,y)=0,得到边缘轮廓曲线C。
具体实施方式五:本实施方式是对实施方式一所述股骨三维模型可视化方法的进一
步限定,步骤三中的变换关系的获取过程为:
步骤三一一、由患者股骨边缘轮廓上的点构成参考图像点集P,P={p1,p2,...,pi},
i=1,2,...n,n为正整数;由通用股骨三维模型投影轮廓上的点构成浮动图像点集Q,
Q={q1,q2,...,qi},并对两个点集进行采样,保证数据量相等。仿射配准前待建模的股骨正
位DR图像轮廓和通用股骨正位投影轮廓见图8;
步骤三一二、将参考图像点集P与浮动图像点集Q代入公式(15)中,依据最小二
乘准则进行参数优化,获得α,β,tx,ty,tz,sx,sy,sz共八个自由度的配准变换参数;配准算
法迭代终止条件为:超过算法的最大迭代次数(一般取值范围为100次~1000次),或浮
动图像点集近两次的迭代误差小于设定阈值(一般取值范围为10-4~10-10),即终止配准算
法的迭代过程。仿射配准后待建模的股骨正位DR图像轮廓和通用股骨正位投影轮廓见图
9。
(
α
,
β
,
t
x
,
t
y
,
t
z
,
s
x
,
s
y
,
s
z
)
=
arg
min
(
Σ
i
=
1
n
|
|
p
i
-
(
S
2
×
2
R
2
×
2
q
i
+
T
2
×
1
)
|
|
2
2
)
-
-
-
(
15
)
]]>
其中α,β分别为绕X轴、Y轴的旋转角度,tx,ty,tz分别为沿X轴、Y轴、Z轴的平
移量,sx,sy,sz分别为沿X轴、Y轴、Z轴的比例系数,S2×2为二维比例变换矩阵,R2×2
为二维旋转变换矩阵,T2×1为二维平移向量,
q
i
=
x
i
y
i
]]>为浮动点集中某点的二维坐标向
量,
p
i
=
a
i
b
i
]]>为参考点集中某点的二维坐标向量。
对于股骨的正位图像,配准变换的具体表达式为:
S
AP
R
AP
q
i
+
T
AP
=
s
y
0
0
s
z
cos
α
-
sin
α
sin
α
cos
α
q
i
+
t
y
t
z
-
-
-
(
16
)
]]>
对于股骨的侧位图像,配准变换关系的具体表达式为:
S
LAT
R
LAT
q
i
+
T
LAT
=
s
x
0
0
s
z
cos
β
sin
β
-
sin
β
cos
β
q
i
+
t
x
t
z
-
-
-
(
17
)
.
]]>
具体实施方式六:本实施方式是对实施方式一所述股骨三维模型可视化方法的进一
步限定,步骤四所述三维变换关系矩阵G满足公式:
X
w
′
Y
w
′
Z
w
′
1
=
G
X
w
Y
w
Z
w
1
=
R
3
×
3
R
3
×
3
T
3
×
1
0
1
X
w
Y
w
Z
w
1
-
-
-
(
18
)
]]>
其中,
X
w
Y
w
Z
w
1
]]>为原始通用模型的三维齐次坐标;
X
w
′
Y
w
′
Z
w
′
1
]]>为姿态调整后股骨模型的三维齐次坐标;
S
3
×
3
=
s
x
0
0
0
s
y
0
0
0
s
z
,
]]>S3×3为三维比例变换矩阵;
R
3
×
3
=
R
z
R
y
R
x
=
cos
β
cos
γ
sin
α
sin
β
cos
γ
-
cos
α
sin
γ
cos
α
sin
β
cos
γ
+
sin
α
sin
γ
cos
β
sin
γ
cos
α
cos
γ
+
sin
α
sin
β
sin
γ
cos
α
sin
β
sin
γ
-
sin
α
cos
γ
-
sin
β
sin
α
cos
β
cos
α
cos
β
]]>
,R3×3为三维旋转变换矩阵;
R
x
=
1
0
0
0
cos
α
-
sin
α
0
sin
α
cos
α
]]>为绕X轴的矩阵分量,
R
y
=
cos
β
0
sin
β
0
1
0
-
sin
β
0
cos
β
]]>为绕Y
轴的矩阵分量,
R
z
=
cos
γ
-
sin
γ
0
sin
γ
cos
γ
0
0
0
1
]]>为Z轴旋转的矩阵分量;
T
3
×
1
=
t
x
t
t
t
z
,
]]>T3×1为三维平移变换矩阵;
G
=
R
3
×
3
R
3
×
3
T
3
×
1
0
1
,
]]>G为总的三维齐次坐标变换矩阵。
这里,除了绕Z轴的旋转角度γ未知以外,{α,β,tx,ty,tz,sx,sy,sz}等八个配准参数均已
求得。参数γ的估计方法如下:
步骤四一、医生根据股骨骨折手术的临床经验,确定出γ的变化范围(一般范围为:
-10°~+10°)。设定γ的初始值为γ0,并与其他八个配准参数构成三维变换矩阵G,根据
公式(18)将G作用到三维通用模型上,对股骨三维模型进行更新;
步骤四二、根据公式(2),将更新后的三维模型重新做正位投影,获取投影轮廓数据。
按公式(19)计算正位投影轮廓与患者DR图像轮廓间的平均误差:
E
=
1
n
Σ
i
=
1
n
min
|
|
p
i
-
q
j
′
|
|
2
2
-
-
-
(
19
)
]]>
其中,Q′={q′1,q′2,...,q′j},j=1,2,...m,m为正整数,为更新后三维模型的正位投影轮
廓点集,P={p1,p2,...,pi},i=1,2,...n,n为正整数,为患者股骨的正位DR图像轮廓点
集。
步骤四三、以0.5°为步长改变γ的角度值,重复步骤四一和步骤四二,在γ的变化
范围获得20组误差E,其中最小的误差E对应的γ值为最优的Z轴旋转角度,至此完成
全部的股骨三维姿态估计。见图10。
此时更新后的股骨三维模型即为经过姿态估计后的患者个体化股骨三维模型。