CT系统参数标定数学模型.docx
- 文档编号:3294835
- 上传时间:2022-11-21
- 格式:DOCX
- 页数:25
- 大小:326.52KB
CT系统参数标定数学模型.docx
《CT系统参数标定数学模型.docx》由会员分享,可在线阅读,更多相关《CT系统参数标定数学模型.docx(25页珍藏版)》请在冰豆网上搜索。
CT系统参数标定数学模型
CT系统参数标定数学模型
摘要
本文通过对体模系统和投影点坐标的分析,提出了较为充分合理的假设,并对数据进行了拟合处理。
先后建立了两个模型(理想状态下的模型,只考虑载物台和探测器等一些因素误差的模型)来求解CT系统参数标定的问题。
在对问题一的分析中,只考虑在理想情况下:
即载物台水平和探测器无偏转,并不考虑机械误差的情况下各个参数之间的关系,我们以探测器左下角为(0,0),右方向为x轴正方向,上方向为y轴正方向,垂直于xoy平面的走直线为z轴建立如图4所示的空间直角坐标系。
运用几何知识作图,通过建立相似图形的等比例方程来研究确定几个系统参数之间的关系(X射线源的位置参数、载物中心线的位置参数和探测器的位置参数)。
最后讨论了引入载物台中心轴倾斜、载物转台不均匀转动和数据不确定度对误差的影响。
在分析问题二时,考虑到所给椭圆点集可能有偏差,于是首先运用MATLAB软件对数据进行处理拟合,得到更为确信的数据并由此求出修正后的“椭圆”方程。
从拟合出来的椭圆方程,我们可以得知椭圆的中心坐标还有其大小,并因此看出我们拟合出来的椭圆并不和我们的理想模型相符合。
所以我们想到了有机械误差存在。
由于同时考虑多个误差因素将使模型计算变得极为复杂,对此,我们假设探测器无偏转,只有载物台倾斜,建立模型二,并对探测器上的两个椭圆进行旋转修正。
易知修正前后两个椭圆的大小并没有改变,只是相对位置发生了变化,发射源A在探测器上的投影坐标也没有改变。
对此,我们以A点为旋转中心对其进行旋转修正并根据拟合后的数据求出两个类椭圆中心连线的斜率。
从而求出载物台倾斜角
0.00210563,并依据此进一步修正数据。
最后在问题三,对模型三的机械误差进行了分析计算。
机械误差主要来源为探测器上数据像素误差、载物台的转轴偏斜,圆心距离的误差,并对其分析了原因。
关键词
CT系统载物台转轴偏斜数据修正拟合机械误差参数标定
问题重述
CT自发明以来,被公认为自伦琴发现X射线后是放射领域的最重要的发明之一。
工业CT一般采用扫描物件旋转和面探测器的结构,工作示意图如图1示:
图1工业CT系统工作原理示意图
扫描物件围绕某一固定转轴旋转,每隔一定角度采集一张图像,然后根据采集的图像采用3D图像重建算法即可将原始3D物件重建出来。
在实际中,由X射线源发出X射线,经过扫描物体衰减后照射在探测器上,探测器根据接收到的光子数的计数实现光电转换,从而形成灰度图像。
一般平板探测器大小为3000×2000像素,每个像素为0.127mm。
在CT系统安装过程中往往存在机械误差,而这些误差对于物件重建的准确性往往是致关重要的,实际中就需要对安装好的CT系统进行参数标定,系统参数主要包括X射线源的位置、载物中心线的位置和探测器的位置参数等。
传统的方法是采用尺子测量法,对于CT系统参数标定精度要求是不够的(通常CT重建过程中系统参数的机械误差不允许超过一个像素)。
目前,一般通过实验的方法实现参数标定。
所谓实验的方法标定参数就是通过扫描已知参数的体模,分析投影数据,估计系统参数的方法。
针对工业CT要求解决如下问题:
(1)建立合适的坐标系,正确描述CT系统的各种参数和机械误差,并建模分析这些参数的关系和可能的机械误差。
(2)通常采用轴承钢球作为CT参数标定体模。
因为它具有各向投影一致,边缘清晰等优点,在里采用两个钢球实现CT的参数标定实验(如图3所示)。
两个钢球置于有机玻璃管中(X射线容易透过有机玻璃,容易后期球心获取的图像处理),钢球直径为8mm±0.0008mm,两球心距离约为100mm±1mm。
将有机玻璃管固定在旋转载物台上,载物台携带有机玻璃管以均匀速度旋转,与此同时X射线源发出射线,探测器采集数据(旋转一周,等间距采集180张)。
附件1给出了根据180张采集图像提取出的球心投影坐标。
两个钢球球心投影的轨迹是两个椭圆(如图4所示)。
试根据实验数据估计该CT系统的参数值,即给出CT系统的标定。
(3)在参数标定过程中可能存在多种可能的机械误差,试就你的参数标定可能的误差进行分析。
图2双球体模型示意图图3双球投影轨迹
(双球轨迹坐标见附录1)
问题分析
题目要求我们求出工业CT系统的各项参数,这些参数包括以下六个:
放射源到转盘的转轴的水平距离
;
转盘的转轴到探测器平面的水平距离
;
放射源到探测器平面的水平距离
;
放射源到xoz平面的垂直距离
;(平面xoz具体位置见问题
(1)解答中的坐标系)
放射源与yoz平面的水平距离
;
转轴到yoz平面的水平距离d。
如图所示:
为了使参数更容易确定,我们选取探测器作为参照物,这样各个参数就有了参考的标准。
问题一要求我们求出系统各个参数之间的关系,还有根据我们所列的关系式进行误差分析。
本题误差主要来源于工业CT系统安装过程。
除了装置本身的由于在装系统时,可能会导致物台中心轴倾斜、载物转台不均匀转动,探测器前后俯仰等,都会导致机械误差的存在。
我们在计算时会考虑这些因素,并建立相应的模型。
所以,在标定时应该给出其机械误差大小。
可以根据体模的投影数据,利用几何知识来求解主要包括X射线源的位置、载物中心线的位置和探测器的位置并求出他们的关系。
问题二则以题目给出的360个坐标为基础,运用matlab软件来拟合两个椭圆,考虑到这些点的坐标本身就有一点的误差,经过一定的修改和筛选后,进行多次拟合进而求出椭圆的方程。
求出椭圆的方程后,再利用在问题
(1)中各个参数的关系,就可以求出各个参数的值。
根据分析,我们知道了拟合出来的两个椭圆的中心连线并不垂直于X轴,我们可以利用绕点旋转的方法对其进行修正。
问题三中的误差分析,本题中考虑的机械误差来源有两圆心距离误差,载物台的中心转轴倾斜,这些根据问题二中的数据,我们已经对它进行了修正。
模型的假设与符号说明
模型假设
1、假设光子波动不明显,即忽略光的衍射现象;
2、将两个钢球的球心看作一个质点;
3、假设实验所收集的数据能客观反映实际情况;
4、玻璃管在随着转盘转动的过程中,始终垂直固定于转盘;
5、载物台在转动过程中,玻璃管与两个钢球始终保持相对静止。
符号说明及名词解释
符号
符号含义
数据类型
放射源到转轴的水平距离
待求参数
转轴到探测器平面的水平距离
待求参数
放射源到探测器平面的水平距离
待求参数
放射源到yoz平面的距离
待求参数
放射源到xoz平面的垂直距离
待求参数
d
转轴到yoz平面的距离
待求参数
r
有机玻璃管和小球的轨迹圆半径
中间值
上椭圆的短轴
中间值
下椭圆的短轴
中间值
b
两椭圆长轴的平均值
中间值
m
两椭圆最高点间的距离
中间值
k
两椭圆最低点间的距离
中间值
t
放射源在探测器上的投影与上椭圆最低点的距离
中间值
h
两钢球球心间的距离
已知量
一、模型建立及问题求解
问题一
建立合适的坐标系,正确描述CT系统的各种参数和机械误差,并建模分析这些参数的关系和可能的机械误差。
将探测器所在的平面定为参照物,结合题目已给出的二维坐标系xoy,在探测器平面左下方原点(0,0)引一条Z轴,形成三维坐标系,如下图示:
依照题意,我们首先得求出载物台的转盘半径r,为了求出转盘的半径,我们先作了这样的一个实验,我们以探测器左下角为(0,0),右方向为x轴正方向,上方向为y轴正方向,垂直于xoy平面的走直线为z轴建立如图1所示的空间直角坐标系。
,再让发射源与玻璃管中的下边小球处于同一个小平面上进行实验。
在实验过程中,我们能在探测器上得到玻璃管中上面小球所形成的一条直线为2b,而下面小球的投影则形成一个类椭圆,我们可以通过探测器上的数据拟合得到椭圆的方程,从而得到我们想要椭圆的长轴长和短轴长还有中心坐标,进而可以计算出有机玻璃管和小球的轨迹圆半径r。
图4
为了简单分析,我们分别画出图1的正视图和俯视图如下:
俯视图
正视图
根据正视图和俯视图,我们利用相似三角形和勾股定理得到了以下关系式:
(1)
(2)
(3)
联立
(1)
(2)(3)可解得:
r=
(4)
根据以上的式子,我们可以根据实验得出的数据求出r的值,然而对于相同的这个装置来说,以下的实验中有机玻璃管和小球的轨迹圆半径r是不会改变的,还有因为其俯视图中的2b的值是发射源垂直于xoy平面发射出来的投影的长度,而在其他情况下是求不出那个值,为了计算方便,在以下的运算中,我们将有机玻璃管和小球的轨迹圆半径r作为一个常数a进行计算,即令r=a(a则为由上面(4)求出的数值)
接着,我们按照题目的要求保持发射源A与载物台和探测器的水平距离不变,即相对位置不变,而在竖直方向上移动发射源A,使其与两个小球都不在同一个水平面进行以下实验,其三维图如下:
为了方便分析,我们作出它的正视图如下:
根据其正视图,利用相似三角形定理和勾股定理可以得到以下关系式:
(5)
(6)
(7)
r=a(其中a为一个常数) (8)
根据CT系统的工作原理,即令扫描物件围绕某一固定转轴旋转,并每隔一定角度采集一张图像,然后根据采集的图像采用3D图像重建算法将原始3D物件重建出来。
再根据重建出来的图像来标定我们的系统参数,由此可知,凡是在重建图像所需要的数据误差都会影响我们的参数标定。
因此,根据上方我们所列出来的关系式可以看出此过程中存在的误差有以下这几种:
(1)载物台没有正常匀速转动,引起重建图像不连续
(2)扫描物件的放置发生倾斜或在运动过程产生抖动,引起重建图像不规则
(3)载物台转动角度发生偏离或载物台的转轴发生倾斜,引起重建图像的数据不好处理
(4)实验中两个钢球的半径大小有一定的误差,使得有机玻璃管和小球的轨迹圆半径r存在误差从而影响参数标定。
问题二
根据照题目给出的数据,我们用MATLAB软件进行拟合,绘制出上下两个椭圆的图形(算法详见附录2),由图形反映拟合情况,并解出椭圆一般方程的系数
、
,进而求解椭圆方程,得到如下结果:
数据图形拟合
上方椭圆的中心坐标为(
,
,0)长半轴长为
,短半轴长为
(单位为像素)
下方椭圆的中心坐标为(
,0)长半轴长为
,短半轴长为
(单位为像素)
上方椭圆的方程为:
下方椭圆的方程为:
由此,我们可以得出B,D,C,E各点的纵坐标分别为
将其坐标代入(5)(6)(7)(8)式可得:
,
,
即我们可以得出发射源A的坐标为(1423.832,811.83225,11.58469a)
综上所述,我们可以对给定CT系统参数进行标定(单位为像素):
X射线源的位置与载物台中心线的距离:
X射线源的位置与探测器的距离:
若以坐标和方程的形式给出,则表示为:
发射源A的坐标为(1423.832,811.83225,11.58469a)
载物台中心线的方程
探测器的位置:
x-o-y平面
模型的进一步分析与修正
在上述的分析过程中,我们发现经过拟合出来的两个椭圆的中心连线并不垂直于X轴,与实际中我们预测分析的有所不同。
由于影响系统的机械误差有很多,在此,我们只认为这是因为载物台在转动过程中所绕的转轴倾斜产生偏角所导致的,与其它因素无关。
所以,我们对其作了修正如下(左为修正前,右为修正后):
由问题二我们已经求出了发射源A的坐标为(1423.832,811.83225,11.58469),由于此过程中只在xoy平面上修正,所以修正后发射源的坐标是不会改变的。
对此,我们以A点在xoy平面上的投影A’为中心,线段
绕着中心A’旋转到线段
垂直于X轴,此时,上下两个椭圆的中心坐标分别为
(1423.832,1812.56265,0),
(1423.832,188.49018,0),由于修正后的椭圆大小是不会改变的,只改变了其位置,所以修正后的椭圆方程分别为:
上方椭圆的方程:
下方椭圆的方程:
此时我们还可以求出修正的角度a,其修正的图形在xoy平面的角度a如下图所示:
以上方椭圆的中心坐标为标准,则修正前的角度a有如下关系:
代入数据可以求出:
=0.00210563
即修正的角度的大小为:
0.00210563
经过了这样的修正,我们可以更加准确的利用探测器上的数据进行拟合图像,使得重建的图像更加具有准确性,我们能够更好的进行标定CT系统的参数。
问题三
除了问题二中所提到的载物台的转轴偏斜引起误差之外,还存在以下两种误差可能引起参数标定不准确。
首先,如果载物台没有正常匀速转动,那么扫描物件的旋转速度会受到影响,探测器所采集到的具体数据可能会发生断点或重点,如下图的断点情况示意图:
其次,扫描物件的放置发生倾斜或在运动过程产生抖动会导致探测器所采集到的图像数据发生偏差,影响了Matlab的拟合效果以及相关数据的准确性。
如图12所示,很明显地,投影得到的图形发生了扭曲,形成机械误差,影响到了参数的标定。
模型评价与改进
优点:
1、对体模的投影点坐标进行数据的修正,使得计算的结果更加准确;
2、考虑到的误差比较合理,有较大的实际意义
3、运用MATLAB软件进行数据运算和图表拟合,使得计算更加简单准确。
缺点:
1、机械误差存在的因素考虑还不够全面,且仅限于理论上的分析,有待实践过程中的进一步检验与改进
2、在数据误差处理时仅用理论知识进行处理,实际操作有些难度。
参考文献
[1]XX文库
[2]中国知网
[3]谢照鸿,范正森等人编著,《数学建模技术》,中国水利水电出版社,2003年9月
[4]周永正,詹棠森,方成鸿,邱望仁,《数学建模》,同济大学出版社,2011年1月
附录
附录1
双球模型投影在探测器上的坐标
上方球心坐标
下方球心坐标
x轴
y轴
x轴
y轴
2762
1792
2737
214
2754
1786
2727
217
2745
1781
2715
220
2735
1775
2701
224
2723
1769
2687
227
2710
1764
2671
231
2695
1758
2653
234
2679
1753
2635
237
2661
1747
2615
240
2643
1742
2594
243
2623
1737
2572
246
2602
1732
2549
249
2580
1727
2525
252
2557
1722
2500
255
2532
1717
2474
257
2507
1713
2447
260
2481
1708
2419
262
2454
1704
2390
265
2425
1700
2360
267
2396
1696
2330
269
2367
1692
2299
271
2336
1688
2267
274
2305
1684
2235
275
2273
1681
2202
277
2240
1677
2168
279
2207
1674
2134
281
2173
1671
2099
283
2138
1668
2064
284
2103
1665
2028
285
2068
1663
1992
287
2032
1660
1955
288
1995
1658
1918
289
1959
1655
1881
290
1921
1653
1843
291
1884
1651
1805
292
1846
1650
1767
293
1807
1648
1728
294
1769
1646
1690
294
1730
1645
1651
295
1691
1644
1612
295
1652
1643
1573
296
1613
1642
1533
296
1573
1641
1494
296
1534
1640
1455
296
1494
1640
1415
296
1454
1640
1376
296
1415
1639
1337
296
1375
1639
1297
296
1336
1640
1258
295
1296
1640
1219
295
1256
1640
1180
294
1217
1641
1141
293
1178
1641
1103
293
1139
1642
1064
292
1100
1643
1026
291
1061
1644
988
290
1023
1646
951
289
985
1647
914
288
947
1649
877
286
910
1651
841
285
873
1652
805
283
836
1654
769
282
800
1657
734
280
764
1659
700
278
729
1661
666
277
694
1664
632
275
660
1667
600
273
626
1670
568
271
593
1673
536
268
561
1676
506
266
529
1679
476
264
499
1683
447
261
468
1686
418
259
439
1690
391
256
410
1694
364
254
383
1698
339
251
356
1702
314
248
330
1706
291
245
305
1711
268
242
281
1715
247
239
259
1720
227
236
237
1725
207
233
216
1730
190
230
197
1735
173
226
179
1740
158
223
162
1745
144
219
146
1750
131
216
132
1756
120
212
119
1761
110
209
107
1767
102
205
97
1773
95
202
89
1778
90
198
82
1784
86
194
76
1790
84
190
72
1796
83
187
70
1802
85
183
69
1808
88
179
70
1814
92
175
73
1821
99
171
77
1827
107
167
84
1833
117
164
92
1839
129
160
102
1845
143
156
113
1852
158
152
127
1858
176
148
143
1864
195
145
160
1870
217
141
179
1877
240
137
201
1883
265
134
224
1889
291
130
249
1895
320
127
276
1900
351
123
304
1906
383
120
335
1912
417
117
367
1917
453
114
401
1923
490
110
438
1928
529
108
475
1933
570
105
514
1938
612
102
555
1942
656
100
598
1947
701
97
642
1951
747
95
687
1955
795
93
734
1959
844
91
782
1963
893
89
832
1966
944
87
882
1970
996
86
933
1972
1049
85
986
1975
1103
84
1039
1977
1157
83
1093
1979
1211
82
1148
1981
1267
81
1203
1983
1322
81
1259
1984
1378
81
1315
1985
1434
81
1371
1985
1489
81
1427
1986
1545
81
1484
1986
1600
82
1540
1985
1655
83
1596
1985
1710
83
1652
1984
1764
84
1707
1983
1817
86
1761
1981
1870
87
1815
1980
1922
89
1868
1978
1973
91
1921
1975
2023
92
1972
1973
2071
95
2022
1970
2119
97
2072
1967
2165
99
2120
1963
2210
102
2166
1960
2253
104
2212
1956
2295
107
2255
1952
2336
110
2297
1948
2374
113
2339
1943
2412
116
2378
1938
- 配套讲稿:
如PPT文件的首页显示word图标,表示该PPT已包含配套word讲稿。双击word图标可打开word文档。
- 特殊限制:
部分文档作品中含有的国旗、国徽等图片,仅作为作品整体效果示例展示,禁止商用。设计者仅对作品中独创性部分享有著作权。
- 关 键 词:
- CT 系统 参数 标定 数学模型