一种计算锥束CT几何参数的方法与流程

本发明属于x射线计算机断层成像(computed tomography,ct)无损检测，具体涉及一种计算锥束ct几何参数的方法。

背景技术：

1、计算机断层成像技术(ct)作为无损检测领域最重要的成像技术之一，具有扫描速度快、成像质量好、分辨率高、通用性强等一系列优点，目前已被广泛应用于医学、安检、工业制造、考古学等各个领域。经过几十年的发展，ct不仅从二维扇束发展为目前主流的三维锥束扫描方式，其重建分辨率和扫描精度也得到了极大的提高。与此同时，ct设备在安装、制造和调试过程中的难度也随之增加。ct的几何结构如图1所示。

2、准确描述一个ct系统的几何参数共有7个，分别是射线源到旋转中心的距离(sod)；射线源经过旋转中心到探测器的距离(sdd)；探测器沿行和列方向的偏移(offsetu、offset v)；探测器沿行的旋转(θ)；探测器沿列的旋转(β)，探测器沿其平面法线的旋转(η)，如图2、图3、图4所示。

3、目前已有的锥束ct几何校正方法主要分为两种，一种是借助已知几何结构的标准模体，采集模体一个或多个角度的投影，通过模体与投影之间的关系来求解ct的几何参数；另一种是不借助任何标准模型，只根据投影数据直接求解几何参数。但标准模体为了保证精度通常加工比较困难、价格较高、通用性也会受到设备类型、分辨率等限制；而不需要标准模体的方法往往计算复杂或者难以实现。

技术实现思路

1、针对现有技术中存在的问题，本发明的目的在于提供一种计算锥束ct几何参数的方法。在最一般的情况下，锥束ct在每个角度下的几何结构可以通过4个向量来表示，如图5所示，其中表示向量os，表示向量od，其中o为坐标系中心点，s为射线源坐标，d为探测器边角坐标；向量os是以坐标中心o为起点、点s为终点的向量，向量od是以坐标中心o为起点、点d为终点的向量。

2、当中心坐标点o确定之后，射线源坐标s，探测器边角坐标d，探测器单位行向量和探测器单位列向量可以完全确定一个锥束ct的几何结构。

3、同时当被投影的体素点p的坐标已知时，可以通过公式(1)用来直接计算p投影在探测器上的坐标pd[u,v](其中为向量op)，

4、

5、我们对公式(1)进行变形，可得到公式(2)，

6、

7、将公式(2)写成矩阵形式，可得，

8、

9、进一步简化公式(3)，可得，

10、

11、其中a为3×4的投影矩阵，可用如下公式表示，

12、

13、本发明的主要内容是通过向量化表示ct的几何结构，从投影矩阵a中通过非线性方程组的形式求解出四个向量，然后再通过这四个向量求出ct的7个几何参数。

14、本发明的技术方案为：

15、一种计算锥束ct几何参数的方法，其步骤包括：

16、1)确定待扫描样品的中心坐标点o，然后根据该中心坐标点o确定锥束ct中射线源坐标s、探测器边角坐标d、探测器单位行向量和探测器单位列向量构造锥束ct的投影矩阵用于根据投影矩阵a计算该待扫描样品的任意体素点p投影在探测器上的坐标pd[u,v]，其中表示向量os，表示向量od；

17、2)对所述投影矩阵a进行求解，得到的具体值；

18、3)根据所得计算得到所述锥束ct的几何参数。

19、进一步的，对所述投影矩阵a进行求解，得到的具体值的方法为：

20、21)将所述投影矩阵a展开为其中

21、22)获取至少6组体素点及其对应的投影点；根据第i组中的体素点坐标pi＝[xi,yi,zi]及其对应的投影点坐标di＝[ui,vi]，得到

22、i＝1～6；

23、23)将步骤22)所得结果转为矩阵形式mia＝0并求解得到投影矩阵a的一维展开形式a＝(a11,a12,a13,a14,a21,a22,a23,a24,a31,a32,a33,a34)t；其中

24、

25、24)根据步骤23)所得a11,a12,a13,a14,a21,a22,a23,a24,a31,a32,a33,a34)和关系式

26、

27、得到sd；其中将投影矩阵a分为a1和a2两部分，

28、

29、

30、得到向量然后通过od＝sd+os求出向量od，即得到进一步的，对所述投影矩阵a进行求解，得到的具体值的方法为：

31、31)将所述投影矩阵a展开为其中

32、32)选取两个根据放置的小球作为校正模体，两小球的半径分别为r1、r2，相对高度为l，c为任意正实数；所述锥束ct沿椭圆轨迹对两小球进行扫描，得到每个小球在椭圆轨迹对应的四个顶点上的坐标其对应的投影点坐标；

33、33)对于步骤32)所得8组坐标，将第i组坐标中的体素点坐标pi＝[xi,yi,zi]及其对应的投影点坐标di＝[ui,vi]，得到

34、i＝1～8；

35、34)将步骤33)所得结果转为矩阵形式mia＝0并求解得到投影矩阵a的一维展开形式a＝(a11,a12,a13,a14,a21,a22,a23,a24,a31,a32,a33,a34)t；其中

36、

37、35)根据步骤34)所得a11,a12,a13,a14,a21,a22,a23,a24,a31,a32,a33,a34和关系式

38、

39、

40、

41、得到sd；其中将投影矩阵a分为a1和a2两部分，

42、

43、

44、得到然后通过od＝sd+os求出向量od，即得到进一步的，根据所得计算得到所述锥束ct的几何参数的方法为：

45、41)设探测器的法线向量为射线源垂直于旋转轴的连线与探测器平面的交点为

46、c，则

47、42)令使射线方向平行于旋转平台；则交点c的坐标为

48、

49、43)计算射线源s与原点o的距离sod、射线源点s和点c的距离sdd；探测器的行列探测单元数量分别为nu、nv，通过公式得到探测器偏移值

50、

51、44)设旋转矩阵为r3×3，当探测器没有偏转时标准向量为与探测器面的法线向量为根据关系式计算得到旋转矩根据旋转矩阵r得到对应的旋转角。

52、进一步的，所述旋转角分别为

53、进一步的，所述锥束ct的几何参数包括offsetu、offsetv、θ、β、η、sod、sdd；其中，sod为射线源到旋转中心的距离，sdd为射线源经过旋转中心到探测器的距离，offsetu为探测器沿行方向的偏移，offsetv为探测器沿列方向的偏移，θ为探测器沿行方向的旋转角度，β为探测器沿列方向的旋转角度，η为探测器沿其平面法线方向的旋转角度。

54、与现有技术相比，本发明的积极效果为：

55、本发明可以根据需要求解的参数灵活选择模体，以下分两种情况进行讨论。

56、(1)当需要的参数不包括sod(如图2所示)、探测器以列为轴的旋转角β(如图4(a)所示)、探测器以行为轴的旋转角θ(如图4(b)所示)时，模体小球的摆放和相对位置都可以忽视。然而如果探测器存在旋转角β和旋转角θ，这两个未知的旋转角一般会对其他参数的求解产生影响。相比于其他校正方法来说，本方法的参数求解精度受两个未知的旋转角的影响更小。

57、当两个随意摆放的小球相对高度已知时，可以进一步求得准确的sod值。

58、(2)当需要求解出ct的所有参数时，则需要7个以上小球的精确三维坐标和对应的投影坐标。在使用转台旋转360°的情况下，模体的条件可以进一步放宽为垂直转台放置的两个及以上小球，同时需要知道小球在转台的旋转半径，由上述信息也可以精确构建小球三维坐标和对应投影坐标关系，进而求出ct系统的所有几何参数。