本发明涉及医学图像处理技术领域，特别涉及一种基于b样条的医学图像倾斜校正方法以及一种非临时性计算机可读存储介质。
背景技术：
在医学图像成像过程中，尤其在ct(computedtomography，电子计算机断层扫描)和mri(magneticresonance，磁共振成像)成像过程中，由于存在相关设备欠稳定、患者随意欲动或情绪紧张等问题，使生成的图像常常存在较为严重的倾斜现象，这给后续的图像匹配和图像融合带来极大的困难。
为了解决图像倾斜的问题，国内外专家学者进行了各种有益的探索，提出一些较为实用的方法，例如投影法(projectionprofile)、霍夫变换法(houghtransform，霍夫变换)和近邻法(componentnearestneighborclustering)。其中，投影法是常用的倾斜校正方法之一，它是基于对投影图形形状进行分析的方法，但是由于需要计算每个倾斜角的投影形状，因此这种方法的计算量非常大；hough变换法是借助于hough变换，用图像空间目标像素的坐标去计算参数空间中参考点的可能轨迹，该方法对于直线图形效果较好，校正精度较高；近邻法是通过找出所有连通区中心点的k个最近邻点，计算每对最近邻点的矢量方向并统计生成直方图，直方图的峰值对应于整个页面的倾斜角，由于在图像中连通成分较多，因此也相当费时。
但是，上述方法存在的问题是，无法对灰度不均匀的倾斜图像进行倾斜面校正。
技术实现要素：
本发明的目的是提供一种基于b样条的医学图像倾斜校正方法，能够实现对灰度不均匀的倾斜图像的校正。
为解决上述问题，本发明的第一方面提供了一种基于b样条的医学图像倾斜校正方法，包括以下步骤：获取待校正的医学图像，并对所述待校正的医学图像进行二值化处理以获取二值化的医学图像；根据b样条梯度算子对所述二值化的医学图像进行边缘检测，并获得二值化的边缘图像；计算所述边缘图像的包围盒，以得到所述边缘图像的实际大小，进而根据所述边缘图像的实际大小截取所述待校正的医学图像的子图像；计算所述子图像的质心和归一化的中心矩，并得到旋转角；将待校正的医学图像按照所述旋转角反向旋转，以获取校正后的医学图像。
进一步地，所述根据b样条梯度算子对所述二值化的医学图像进行边缘检测采用如下公式：
其中，▽f为像素点(x,y)的梯度大小，bx和by为b样条梯度算子。
进一步地，所述计算所述二值化的边缘图像的包围盒为计算所述二值化的边缘图像的上、下、左、右的边界。
进一步地，所述子图像的质心通过所述子图像的零阶矩和一阶矩计算。
进一步地，所述旋转角采用如下公式获得：
其中，m′1,1、m′2,0和m'0,2为中心矩。
根据本发明的另一个方面提供一种非临时性计算机可读存储介质，其上存储有计算机程序，该程序被处理器执行时实现所述的医学图像倾斜校正方法。
根据本发明提出的基于b样条的医学图像倾斜校正方法，对待校正的医学图像进行二值化处理并获取二值化的医学图像，根据b样条梯度算子对二值化的医学图像进行边缘检测，并获得二值化的边缘图像，计算二值化的边缘图像的包围盒，得到边缘图像的实际大小，根据边缘图像的事迹大小截取待校正的医学图像的子图像，计算子图像的质心和归一化的中心矩，并得到旋转角，对待校正的医学图像按照旋转角反向旋转，获取校正后的医学图像。
本发明的上述技术方案具有如下有益的技术效果：能够快速准确的检测出图像的倾斜角，还适用于灰度不均的倾斜图像。
附图说明
图1为根据本发明实施例的基于b样条的医学图像倾斜校正方法的流程图；
图2a为根据本发明一个实施例的医学图像正倾斜模型；
图2b为根据本发明另一个实施例的医学图像负倾斜模型；
图3a为根据本发明一个具体实施例的待校正的医学图像；
图3b为根据本发明另一个具体实施例的待校正的医学图像；
图4a为根据本发明一个具体实施例的待校正的医学图像的边缘图像；
图4b为根据本发明另一个具体实施例的待校正的医学图像的边缘图像；
图5a为根据本发明一个具体实施例的边缘图像的包围盒图像；
图5b为根据本发明另一个具体实施例的边缘图像的包围盒图像；
图6a为根据本发明一个具体实施例的待校正图像的子图像；
图6b为根据本发明另一个具体实施例的待校正图像的子图像；
图7a为根据本发明一个具体实施例的校正后的医学图像；
图7b为根据本发明另一个具体实施例的校正后的医学图像；
其中，一个具体实例为对ct图像进行校正，另一个具体实施例为对mri图像进行校正。
具体实施方式
为使本发明的目的、技术方案和优点更加清楚明了，下面结合具体实施方式并参照附图，对本发明进一步详细说明。应该理解，这些描述只是示例性的，而并非要限制本发明的范围。此外，在以下说明中，省略了对公知结构和技术的描述，以避免不必要地混淆本发明的概念。
在附图中示出了根据本发明实施例的层结构示意图。这些图并非是按比例绘制的，其中为了清楚的目的，放大了某些细节，并且可能省略了某些细节。图中所示出的各种区域、层的形状以及它们之间的相对大小、位置关系仅是示例性的，实际中可能由于制造公差或技术限制而有所偏差，并且本领域技术人员根据实际所需可以另外设计具有不同形状、大小、相对位置的区域/层。
显然，所描述的实施例是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。
在本发明的描述中，需要说明的是，术语“第一”、“第二”、“第三”仅用于描述目的，而不能理解为指示或暗示相对重要性。
下面参考附图来描述本发明实施例提出的基于b样条的图像倾斜校正方法。
图1为根据本发明实施例的基于b样条的医学图像倾斜校正方法的流程图。根据图像区域倾斜方向，医学图像有两种可能的倾斜方式：正倾斜和负倾斜，如图2a和2b所示，正倾斜表示相对于x轴，图像逆时针旋转了α度，负倾斜表示相对于x轴，图像顺时针旋转了α度，由此可知，图像倾斜的实质是图像主轴x‘与坐标轴x之间有一个倾斜角α，只要求得该倾斜角α，即可把图像绕质心o旋转-α度即可。
如图1所示，本发明实施例的基于b样条的医学图像倾斜校正方法包括以下步骤：
s1：获取待校正的医学图像，并对待校正的医学图像进行二值化处理以获取二值化的医学图像。
需要说明的是，对医学图像进行二值化处理即使待校正的医学图像的灰度值为0或1。其中，待处理的医学图像f的大写可为m×n，左上角像素位置为(1，1)，像素点(x,y)的灰度值为f(x,y)。如图3a和3b所示，分别是本发明一个具体实施例中待处理的ct图像和另一个具体实施例中待处理的mri图像。
s2：根据b样条梯度算子对二值化的医学图像进行边缘检测，并获得二值化的边缘图像。
需要说明的是，n次b样条函数的定义为：
其中，
由式(1)可知b3(x)的表达式为：
对于3次b样条函数的一阶导数为：
其中，sgn(x)为符号函数，并且有
bk×l(x-i,y-j)＝bk(x-i)bl(y-j)5)
由式(4)和(5)可知，曲面是待处理的医学图像f和k次样条及l次样条函数离散卷积的结果，在本发明实施例中，k＝l＝3，当x,y取整数时，nk,l(x,y)为像素点(x,y)的3×3领域。由式(4)可知曲面是二次连续可微的，即该曲面是一个光滑曲面。根据b样条函数的局部性质，利用式(4)和(5)可以进一步得到：
根据式(6)，可以得到沿和方向倒数：
由式(2)、(3)和(7)，分别计算3次b样条拟合曲面沿和方向倒数卷积模版：
其中，从边缘检测的角度来看，bx和by分别为水平边缘和垂直边缘卷积核，即b样条的梯度算子。
进一步地，根据b样条梯度算子对二值化的医学图像进行边缘检测采用如下公式：
其中，▽f为像素点(x,y)的梯度大小，bx和by为b样条梯度算子。如图4a和4b所示，分别根据b样条梯度算子对二值化的医学图像进行边缘检测获得ct图像和mri图像的边缘图像。
s3：计算边缘图像的包围盒，以得到边缘图像的实际大小，进而根据边缘图像的实际大小截取待校正的医学图像的子图像。
其中，计算二值化的边缘图像的包围盒为计算二值化的边缘图像的上、下、左、右的边界。
具体地，b_flag可设定为边界标记，b_flag∈{1,2,3,4}分别表示寻找包围盒的上、下、左、右的边界，value为边界的初始行号或列号，b_value为边界值；seek_f表示当前边界标志，seek_f＝1表示还没有找到当前边界；seek_f＝0表示找到当前边界。
如果b_flag＝1，则在当前行所有像素灰度值累加和大于等于1时，即在当前行至少有1个白色像素点时，判断找到包围盒的上边界，令seek_f＝0；否则b_value＝b_value+1，继续逐行寻找。
如果b_flag＝2，则在当前行所有像素灰度值累加和大于等于1时，即在当前行至少有1个白色像素点时，判断找到包围盒的下边界，令seek_f＝0；否则b_value＝b_value-1，继续逐行寻找。
如果b_flag＝3，则在当前行所有像素灰度值累加和大于等于1时，即在当前列至少有1个白色像素点时，判断找到包围盒的左边界，令seek_f＝0；否则b_value＝b_value+1，继续逐列寻找。
如果b_flag＝4，则在当前行所有像素灰度值累加和大于等于1时，即在当前列至少有1个白色像素点时，判断找到包围盒的右边界，令seek_f＝0；否则b_value＝b_value-1，继续逐列寻找。
其中，得到边缘图像的包围盒即为得到边缘图像的实际大小，可根据边缘图像的实际大小截取待校正的医学图像的子图像f_sub。如图5a和5b分别为ct图像和mri图像对应的包围盒图像，图6a和6b为根据包围盒图像截取的ct图像和mri图像的子图像。
s4：计算子图像的质心和归一化的中心矩，并得到旋转角。
需要说明的是，由于图像区域的某些矩对于平移、旋转、尺度等几何变换具有一些不变的特性，因此矩的表示方法在物体分类和识别方面具有重要意义。
对于二元有界函数g(x,y)，它的(p+q)阶矩定义为：
对二维离散函数m(x,y)，(9)式可表示为：
其中，(p+q)为矩的阶。
由此，零阶矩是物体的面积，二元有界函数g(x,y)的零阶矩为：
对于一幅数字图像即待校正的医学图像，可以认为是一个二维离散函数m(x,y)，因此，其零阶矩为：
由此可知，二维离散函数m(x,y)的一阶矩或高阶矩，除以m0,0后，均与物体的大小无关，即待校正的图像一阶矩或高阶矩除以m0,0后，与待校正的医学图像的子图像实际大小无关。
当p＝1，q＝0和p＝0，q＝1时，
则
具体地，由于
因为中心矩以质心作为原点进行计算，所以中心矩具有位置无关性。
进一步地，使二阶中心矩从m′1,1变得最小的旋转角α可以由下式得出：
进一步地，可以获取旋转角α为：
其中，x，y为物体的主轴，m′1,1、m'2,0和m'0,2为中心矩。
基于上述分析，可以认为，医学图像倾斜的原因是由于图像的主轴绕质心相对于水平轴旋转了一定角度α，上述物体为待校正的医学图像的子图像，其中心点o即为质心。
s5：将待校正的医学图像按照旋转角反向旋转，以获取校正后的医学图像。
由此，在计算获得旋转角α后，将待校正的医学图像以中心点为旋转轴按照旋转角进行反向旋转，即可获取校正后的医学图像。如图7a和7b所示，分别为经过校正后的ct图像和mri图像。
综上所述，根据本发明提出的基于b样条的医学图像倾斜校正方法，对待校正的医学图像进行二值化处理并获取二值化的医学图像，根据b样条梯度算子对二值化的医学图像进行边缘检测，并获得二值化的边缘图像，计算二值化的边缘图像的包围盒，得到边缘图像的实际大小，根据边缘图像的事迹大小截取待校正的医学图像的子图像，计算子图像的质心和归一化的中心矩，并得到旋转角，将待校正的医学图像按照旋转角反向旋转，获取校正后的医学图像。由此，本发明实施例的医学图像校正方法能够快速准确的检测出图像的倾斜角，还适用于灰度不均的倾斜图像。
本发明实施例还提出了一种非临时性计算机可读存储介质，其上存储有计算机程序，该程序被处理器执行时实现的医学图像倾斜校正方法。
根据本发明实施例的非临时性计算机可读存储介质，通过基于b样条的医学图像倾斜校正方法能够快速准确的检测出图像的倾斜角，还适用于灰度不均的倾斜图像。
应当理解的是，本发明的上述具体实施方式仅仅用于示例性说明或解释本发明的原理，而不构成对本发明的限制。因此，在不偏离本发明的精神和范围的情况下所做的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。此外，本发明所附权利要求旨在涵盖落入所附权利要求范围和边界、或者这种范围和边界的等同形式内的全部变化和修改例。
以上参照本发明的实施例对本发明予以了说明。但是，这些实施例仅仅是为了说明的目的，而并非为了限制本发明的范围。本发明的范围由所附权利要求及其等价物限定。不脱离本发明的范围，本领域技术人员可以做出多种替换和修改，这些替换和修改都应落在本发明的范围之内。
尽管已经详细描述了本发明的实施方式，但是应该理解的是，在不偏离本发明的精神和范围的情况下，可以对本发明的实施方式做出各种改变、替换和变更。
显然，上述实施例仅仅是为清楚地说明所作的举例，而并非对实施方式的限定。对于所属领域的普通技术人员来说，在上述说明的基础上还可以做出其它不同形式的变化或变动。这里无需也无法对所有的实施方式予以穷举。而由此所引伸出的显而易见的变化或变动仍处于本发明创造的保护范围之中。
本领域内的技术人员应明白，本发明的实施例可提供为方法、系统、或计算机程序产品。因此，本发明可采用完全硬件实施例、完全软件实施例、或结合软件和硬件方面的实施例的形式。而且，本发明可采用在一个或多个其中包含有计算机可用程序代码的计算机可用存储介质(包括但不限于磁盘存储器、cd-rom、光学存储器等)上实施的计算机程序产品的形式。
本发明是参照根据本发明实施例的方法、设备(系统)、和计算机程序产品的流程图和/或方框图来描述的。应理解可由计算机程序指令实现流程图和/或方框图中的每一流程和/或方框、以及流程图和/或方框图中的流程和/或方框的结合。可提供这些计算机程序指令到通用计算机、专用计算机、嵌入式处理机或其他可编程数据处理设备的处理器以产生一个机器，使得通过计算机或其他可编程数据处理设备的处理器执行的指令产生用于实现在流程图一个流程或多个流程和/或方框图一个方框或多个方框中指定的功能的装置。
这些计算机程序指令也可存储在能引导计算机或其他可编程数据处理设备以特定方式工作的计算机可读存储器中，使得存储在该计算机可读存储器中的指令产生包括指令装置的制造品，该指令装置实现在流程图一个流程或多个流程和/或方框图一个方框或多个方框中指定的功能。
这些计算机程序指令也可装载到计算机或其他可编程数据处理设备上，使得在计算机或其他可编程设备上执行一系列操作步骤以产生计算机实现的处理，从而在计算机或其他可编程设备上执行的指令提供用于实现在流程图一个流程或多个流程和/或方框图一个方框或多个方框中指定的功能的步骤。