发明名称 一种心脏核磁共振图像的快速分割方法
摘要 本发明涉及一种心脏核磁共振图像的快速分割方法,包括以下步骤:一、对取得的心脏核磁共振图像进行高斯滤波预处理;二、在预处理后的图像上计算基于扩展邻域和噪声平滑的广义梯度矢量流的外力场;三、定义心脏左心室内膜初始化轮廓位置;四、对心脏左心室内膜进行分割;五、将心脏左心室内膜的最终分割轮廓结果定义为心脏左心室外膜的初始化轮廓位置;六、将原始边缘图中内膜轮廓所包围区域的边缘强度置为0,重新计算基于扩展邻域和噪声平滑的广义梯度矢量流的外力场;七、对心脏左心室外膜进行分割。本发明基于卷积运算,考虑了椭圆形状能量约束,具有运算速度快、捕捉范围大、抗噪能力强的优点,且在弱边界保护和深度凹陷区域的分割上性能卓越,能准确地分割心脏左心室内、外膜。
申请公布号 CN102509292A 申请公布日期 2012.06.20
申请号 CN201110342079.2 申请日期 2011.11.03
申请人 北京理工大学 发明人 刘利雄;赵恒博;魏军
分类号 G06T7/00(2006.01)I 主分类号 G06T7/00(2006.01)I
代理机构 代理人
主权项 1.一种心脏核磁共振图像的快速分割方法,其特征在于包括以下步骤:一、对取得的心脏核磁共振图像进行高斯滤波预处理根据方程<img file="FDA0000105033830000011.GIF" wi="261" he="50" />对取得的心脏核磁共振图像进行高斯滤波预处理,式中I<sub>0</sub>为输入的原始图像结构信息,G<sub>σ</sub>为标准差为σ的二维高斯函数,<img file="FDA0000105033830000012.GIF" wi="37" he="37" />表示卷积运算;二、在预处理后的图像上计算基于扩展邻域和噪声平滑的广义梯度矢量流的外力场,记为ENGGVF外力场,具体方法为:1)定义ENGGVF外力场F<sub>out</sub>的初始值定义图像I(x,y)的边缘映射为f(x,y),设f<sub>x</sub>和f<sub>y</sub>分别为边缘映射f沿x轴方向和y轴方向的一阶导数,边缘图的梯度向量<img file="FDA0000105033830000013.GIF" wi="399" he="55" />构成了一个向量场V(x,y)=[u(x,y),v(x,y)]=[f<sub>x</sub>,f<sub>y</sub>],作为外力场F<sub>out</sub>的初始值;2)根据外力场F<sub>out</sub>的初始值,计算ENGGVF外力场ENGGVF外力场的迭代公式如下:<maths num="0001"><![CDATA[<math><mrow><mfenced open='{' close=''><mtable><mtr><mtd><msub><mi>u</mi><mi>t</mi></msub><mo>=</mo><mi>g</mi><mrow><mo>(</mo><mo>|</mo><mo>&dtri;</mo><mi>f</mi><mo>|</mo><mo>)</mo></mrow><mrow><mo>(</mo><mi>u</mi><mo>&CircleTimes;</mo><msub><mi>NS</mi><mn>24</mn></msub><mo>-</mo><mi>u</mi><mo>&CircleTimes;</mo><msub><mi>RM</mi><mn>24</mn></msub><mo>)</mo></mrow><mo>-</mo><mi>h</mi><mrow><mo>(</mo><mo>|</mo><mo>&dtri;</mo><mi>f</mi><mo>|</mo><mo>)</mo></mrow><mrow><mo>(</mo><mi>u</mi><mo>-</mo><msub><mi>f</mi><mi>x</mi></msub><mo>)</mo></mrow></mtd></mtr><mtr><mtd><msub><mi>v</mi><mi>t</mi></msub><mo>=</mo><mi>g</mi><mrow><mo>(</mo><mo>|</mo><mo>&dtri;</mo><mi>f</mi><mo>|</mo><mo>)</mo></mrow><mrow><mo>(</mo><mi>v</mi><mo>&CircleTimes;</mo><msub><mi>NS</mi><mn>24</mn></msub><mo>-</mo><mi>v</mi><mo>&CircleTimes;</mo><msub><mi>RM</mi><mn>24</mn></msub><mo>)</mo></mrow><mo>-</mo><mi>h</mi><mrow><mo>(</mo><mo>|</mo><mo>&dtri;</mo><mi>f</mi><mo>|</mo><mo>)</mo></mrow><mrow><mo>(</mo><mi>v</mi><mo>-</mo><msub><mi>f</mi><mi>y</mi></msub><mo>)</mo></mrow></mtd></mtr></mtable></mfenced><mo>-</mo><mo>-</mo><mo>-</mo><mrow><mo>(</mo><mn>1</mn><mo>)</mo></mrow></mrow></math>]]></maths>其中,<img file="FDA0000105033830000015.GIF" wi="37" he="37" />表示卷积,<maths num="0002"><![CDATA[<math><mrow><mi>g</mi><mrow><mo>(</mo><mo>|</mo><mo>&dtri;</mo><mi>f</mi><mo>|</mo><mo>)</mo></mrow><mo>=</mo><mi>exp</mi><mrow><mo>(</mo><mo>-</mo><msup><mrow><mo>|</mo><mo>&dtri;</mo><mi>f</mi><mo>|</mo></mrow><mn>2</mn></msup><mo>/</mo><msup><mi>k</mi><mn>2</mn></msup><mo>)</mo></mrow><mo>,</mo></mrow></math>]]></maths>k为调整常量,<maths num="0003"><![CDATA[<math><mrow><mi>h</mi><mrow><mo>(</mo><mo>|</mo><mo>&dtri;</mo><mi>f</mi><mo>|</mo><mo>)</mo></mrow><mo>=</mo><mn>1</mn><mo>-</mo><mi>g</mi><mrow><mo>(</mo><mo>|</mo><mo>&dtri;</mo><mi>f</mi><mo>|</mo><mo>)</mo></mrow><mo>,</mo></mrow></math>]]></maths><maths num="0004"><![CDATA[<math><mrow><msub><mi>NS</mi><mn>24</mn></msub><mo>=</mo><mfenced open='[' close=']'><mtable><mtr><mtd><mn>0</mn></mtd><mtd><mn>0</mn></mtd><mtd><mn>1</mn><mo>/</mo><mn>12</mn></mtd><mtd><mn>0</mn></mtd><mtd><mn>0</mn></mtd></mtr><mtr><mtd><mn>0</mn></mtd><mtd><mn>0</mn></mtd><mtd><mn>1</mn><mo>/</mo><mn>12</mn></mtd><mtd><mn>0</mn></mtd><mtd><mn>0</mn></mtd></mtr><mtr><mtd><mn>1</mn><mo>/</mo><mn>12</mn></mtd><mtd><mn>1</mn><mo>/</mo><mn>12</mn></mtd><mtd><mn>1</mn><mo>/</mo><mn>3</mn></mtd><mtd><mn>1</mn><mo>/</mo><mn>12</mn></mtd><mtd><mn>1</mn><mo>/</mo><mn>12</mn></mtd></mtr><mtr><mtd><mn>0</mn></mtd><mtd><mn>0</mn></mtd><mtd><mn>1</mn><mo>/</mo><mn>12</mn></mtd><mtd><mn>0</mn></mtd><mtd><mn>0</mn></mtd></mtr><mtr><mtd><mn>0</mn></mtd><mtd><mn>0</mn></mtd><mtd><mn>1</mn><mo>/</mo><mn>12</mn></mtd><mtd><mn>0</mn></mtd><mtd><mn>0</mn></mtd></mtr></mtable></mfenced><mo>,</mo></mrow></math>]]></maths><maths num="0005"><![CDATA[<math><mrow><msub><mi>RM</mi><mn>24</mn></msub><mo>=</mo><mfenced open='[' close=']'><mtable><mtr><mtd><mn>1</mn><mo>/</mo><mn>24</mn></mtd><mtd><mn>1</mn><mo>/</mo><mn>24</mn></mtd><mtd><mn>1</mn><mo>/</mo><mn>24</mn></mtd><mtd><mn>1</mn><mo>/</mo><mn>24</mn></mtd><mtd><mn>1</mn><mo>/</mo><mn>24</mn></mtd></mtr><mtr><mtd><mn>1</mn><mo>/</mo><mn>24</mn></mtd><mtd><mn>1</mn><mo>/</mo><mn>24</mn></mtd><mtd><mn>1</mn><mo>/</mo><mn>24</mn></mtd><mtd><mn>1</mn><mo>/</mo><mn>24</mn></mtd><mtd><mn>1</mn><mo>/</mo><mn>24</mn></mtd></mtr><mtr><mtd><mn>1</mn><mo>/</mo><mn>24</mn></mtd><mtd><mn>1</mn><mo>/</mo><mn>24</mn></mtd><mtd><mn>0</mn></mtd><mtd><mn>1</mn><mo>/</mo><mn>24</mn></mtd><mtd><mn>1</mn><mo>/</mo><mn>24</mn></mtd></mtr><mtr><mtd><mn>1</mn><mo>/</mo><mn>24</mn></mtd><mtd><mn>1</mn><mo>/</mo><mn>24</mn></mtd><mtd><mn>1</mn><mo>/</mo><mn>24</mn></mtd><mtd><mn>1</mn><mo>/</mo><mn>24</mn></mtd><mtd><mn>1</mn><mo>/</mo><mn>24</mn></mtd></mtr><mtr><mtd><mn>1</mn><mo>/</mo><mn>24</mn></mtd><mtd><mn>1</mn><mo>/</mo><mn>24</mn></mtd><mtd><mn>1</mn><mo>/</mo><mn>24</mn></mtd><mtd><mn>1</mn><mo>/</mo><mn>24</mn></mtd><mtd><mn>1</mn><mo>/</mo><mn>24</mn></mtd></mtr></mtable></mfenced></mrow></math>]]></maths>为卷积模板,迭代次数由用户设定,迭代计算到ENGGVF外力场V(x,y)=[u(x,y),v(x,y)]稳定为止;三、在预处理后的图像上定义心脏左心室内膜初始化轮廓位置初始轮廓任意选取位于心脏左心室内膜范围内的一个圆;四、对心脏左心室内膜进行分割根据步骤二计算出的ENGGVF外力场,在初始化轮廓确定的情况下分割出内膜,在曲线演化过程中需添加椭圆形状能量约束,具体的分割过程如下:1)构造椭圆形状能量约束场定义轮廓线为曲线c(s)=(x(s),y(s)),其中s∈[0,1],构造椭圆形状能量约束场,引入椭圆形状能量约束项:<maths num="0006"><![CDATA[<math><mrow><msub><mi>E</mi><mi>ellipse</mi></msub><mo>=</mo><mfrac><mn>1</mn><mn>2</mn></mfrac><mo>&times;</mo><msubsup><mo>&Integral;</mo><mn>0</mn><mn>1</mn></msubsup><msup><mrow><mo>(</mo><mrow><mo>(</mo><mi>x</mi><mrow><mo>(</mo><mi>s</mi><mo>)</mo></mrow><mo>-</mo><msub><mi>x</mi><mi>c</mi></msub><mo>)</mo></mrow><mi>cos</mi><mrow><mo>(</mo><mi>&theta;</mi><mo>)</mo></mrow><mo>+</mo><mrow><mo>(</mo><mi>y</mi><mrow><mo>(</mo><mi>s</mi><mo>)</mo></mrow><mo>-</mo><msub><mi>y</mi><mi>c</mi></msub><mo>)</mo></mrow><mi>sin</mi><mrow><mo>(</mo><mi>&theta;</mi><mo>)</mo></mrow><mo>-</mo><msub><mi>r</mi><mn>1</mn></msub><mi>cos</mi><mrow><mo>(</mo><mn>2</mn><mi>&pi;s</mi><mo>-</mo><mi>&theta;</mi><mo>)</mo></mrow><mo>)</mo></mrow><mn>2</mn></msup><mi>ds</mi></mrow></math>]]></maths>(2)<maths num="0007"><![CDATA[<math><mrow><mo>+</mo><mfrac><mn>1</mn><mn>2</mn></mfrac><mo>&times;</mo><msubsup><mo>&Integral;</mo><mn>0</mn><mn>1</mn></msubsup><msup><mrow><mo>(</mo><mo>-</mo><mrow><mo>(</mo><mi>x</mi><mrow><mo>(</mo><mi>s</mi><mo>)</mo></mrow><mo>-</mo><msub><mi>x</mi><mi>c</mi></msub><mo>)</mo></mrow><mi>sin</mi><mrow><mo>(</mo><mi>&theta;</mi><mo>)</mo></mrow><mo>+</mo><mrow><mo>(</mo><mi>y</mi><mrow><mo>(</mo><mi>s</mi><mo>)</mo></mrow><mo>-</mo><msub><mi>y</mi><mi>c</mi></msub><mo>)</mo></mrow><mi>cos</mi><mrow><mo>(</mo><mi>&theta;</mi><mo>)</mo></mrow><mo>-</mo><msub><mi>r</mi><mn>2</mn></msub><mi>sin</mi><mrow><mo>(</mo><mn>2</mn><mi>&pi;s</mi><mo>-</mo><mi>&theta;</mi><mo>)</mo></mrow><mo>)</mo></mrow><mn>2</mn></msup><mi>ds</mi></mrow></math>]]></maths>其中,(x<sub>c</sub>,y<sub>c</sub>)为椭圆中心,θ为椭圆偏转角,r<sub>1</sub>,r<sub>2</sub>分别为椭圆的两半径,[c<sub>x</sub>,c<sub>y</sub>,θ,r<sub>1</sub>,r<sub>2</sub>]可以通过最小二乘法拟合;根据变分法原理,得到下列离散形式的欧拉方程<maths num="0008"><![CDATA[<math><mrow><mfenced open='{' close=''><mtable><mtr><mtd><msub><mi>x</mi><mi>i</mi></msub><mo>-</mo><msub><mi>x</mi><mi>c</mi></msub><mo>-</mo><msub><mi>r</mi><mn>1</mn></msub><mi>cos</mi><mrow><mo>(</mo><mn>2</mn><mi>&pi;i</mi><mo>/</mo><mi>n</mi><mo>-</mo><mi>&theta;</mi><mo>)</mo></mrow><mi>cos</mi><mrow><mo>(</mo><mi>&theta;</mi><mo>)</mo></mrow><msub><mrow><mo>+</mo><mi>r</mi></mrow><mn>2</mn></msub><mi>sin</mi><mrow><mo>(</mo><mn>2</mn><mi>&pi;i</mi><mo>/</mo><mi>n</mi><mo>-</mo><mi>&theta;</mi><mo>)</mo></mrow><mi>sin</mi><mrow><mo>(</mo><mi>&theta;</mi><mo>)</mo></mrow><mo>=</mo><mn>0</mn></mtd></mtr><mtr><mtd><msub><mi>y</mi><mi>i</mi></msub><mo>-</mo><msub><mi>y</mi><mi>c</mi></msub><mo>-</mo><msub><mi>r</mi><mn>1</mn></msub><mi>cos</mi><mrow><mo>(</mo><mn>2</mn><mi>&pi;i</mi><mo>/</mo><mi>n</mi><mo>-</mo><mi>&theta;</mi><mo>)</mo></mrow><mi>sin</mi><mrow><mo>(</mo><mi>&theta;</mi><mo>)</mo></mrow><mo>-</mo><msub><mi>r</mi><mn>2</mn></msub><mi>sin</mi><mrow><mo>(</mo><mn>2</mn><mi>&pi;i</mi><mo>/</mo><mi>n</mi><mo>-</mo><mi>&theta;</mi><mo>)</mo></mrow><mi>cos</mi><mrow><mo>(</mo><mi>&theta;</mi><mo>)</mo></mrow><mo>=</mo><mn>0</mn></mtd></mtr></mtable></mfenced><mo>-</mo><mo>-</mo><mo>-</mo><mrow><mo>(</mo><mn>3</mn><mo>)</mo></mrow></mrow></math>]]></maths>其中,i=0,1,…,n-1;从而得到椭圆形状能量约束场为:F<sub>ellipse</sub>=(x<sub>i</sub>,y<sub>i</sub>)=[x<sub>c</sub>+r<sub>1</sub>cos(2πi/n-θ)cos(θ)-r<sub>2</sub>sin(2πi/n-θ)sin(θ),                                                    (4)y<sub>c</sub>+r<sub>1</sub>cos(2πi/n-θ)sin(θ)+r<sub>2</sub>sin(2πi/n-θ)cos(θ)]2)构造曲线迭代公式为:c(s)=λ<sub>1</sub>F<sub>int</sub>+λ<sub>2</sub>F<sub>out</sub>+λ<sub>3</sub>F<sub>ellipse</sub>                            (5)其中,λ<sub>1</sub>,λ<sub>2</sub>和λ<sub>3</sub>分别为内力场、ENGGVF外力场和椭圆形状能量约束场的权重系数;内力场F<sub>int</sub>=αc<sub>ss</sub>(s)-βc<sub>ssss</sub>(s),其中,α和β为弹性和刚性系数,c<sub>ss</sub>(s)为曲线c(s)关于s的二阶导数,c<sub>ssss</sub>(s)为曲线c(S)关于s的四阶导数;而外力场为步骤二计算出的ENGGVF外力场;椭圆形状约束能量场采用公式(4);3)将初始化轮廓作为曲线迭代公式的初始值,迭代计算上述曲线迭代公式(5),得到一个稳定的解,即曲线收敛到左心室内膜的轮廓;五、将心脏左心室内膜的最终分割轮廓结果定义为心脏左心室外膜的初始化轮廓位置;六、将原始边缘图中内膜轮廓所包围区域的边缘强度置为0,这就抹平了左室壁内膜边缘及部分噪声,再采用这一改动的边缘图来重新计算ENGGVF外力场;七、对心脏左心室外膜进行分割在步骤六计算出的ENGGVF外力场作用下,根据步骤五定义的外膜初始化轮廓,对心脏左心室外膜进行分割,在曲线演化过程中需添加椭圆形状能量约束,具体的左心室外膜分割过程等同步骤四的左心室内膜分割过程。
地址 100081 北京市海淀区中关村南大街5号