王航 崔尚 郑维连 张振荣
摘要:将数据进行解析变换,引入 Radon 逆变换,推导出各种不同形式的逼近表达式,再离散计算,最后运用MATLAB进行编程得到关于未知介质形状的几何图形与图形的吸收率数据。
关键词:MATLAB;解析变换;Radon 变换
1 未知介质几何形状与吸收率的确定
首先将数据导入到MATLAB,画图。
根据已有结果知该系统的初始扫描位置为,铅直方向直线逆时针旋转29.6463°处(铅直线向下正方向)。需要对图形数据进行修正。经修正后,才能真正反映未知介质在正方形托盘中的位置关系。修正后附件3的数据画出图形:
由图可知,未知介质在xy 平面内分布不均匀,即衰减系数 μ =μ(x, y) ,则可得在某一方向l ,沿某一路径 L 的总衰减为:
[JZ(][XCimage310.tif;%50%50,JZ][JZ)][JY](1)
衰减的集合[XCimage311.tif;%50%50,JZ]称为投影,即接收信息,反投影重建图像的原理就是根据一系列的投影[XCimage311.tif;%50%50,JZ],求出被积函数[XCimage312.tif;%50%50,JZ],从而得出[XCimage313.tif;%50%50,JZ]分布(密度分布)的
图像,即我们所要求的图形经过CT扫描之后的最终图像以及个点的吸收率。所以我们引入Radon 变换。
Radon 变换的定义为:若[XCimage314.tif;%50%50,JZ]基本空间[XCimage315.tif;%50%50,JZ],[XCimage316.tif;%50%50,JZ]对于所有的[XCimage317.tif;%50%50,JZ]是已知的,则称[XCimage316.tif;%50%50,JZ]是二维函数[XCimage318.tif;%50%50,JZ]的Radon变换,记为Rf。
由于研究的图形为二维图形所以,在Radon 变换中考虑直线表达为[XCimage319.tif;%50%50,JZ]
其中[XCimage320.tif;%50%50,JZ]为法线长(即原点到直线的垂线长),[XCimage321.tif;%50%50,JZ]为法线与[XCimage322.tif;%50%50,JZ]轴的交点。[XCimage323.tif;%50%50,JZ]为直
线[XCimage324.tif;%50%50,JZ]的位置参数。
积分表达为:
[JZ][XCimage325.tif;%50%50,JZ]
现令[XCimage326.tif;%50%50,JZ]为直线[XCimage327.tif;%50%50,JZ]的方向,[XCimage328.tif;%50%50,JZ]为[XCimage329.tif;%50%50,JZ]的法线方向,则可得到[XCimage330.tif;%50%50,JZ]得新的坐标系,
[JZ(][XCimage331.tif;%50%50,JZ][JZ)][JY](2)
则有积分[XCimage332.tif;%50%50,JZ][JY](3)
记[XCimage333.tif;%50%50,JZ]
[JZ][XCimage334.tif;%50%50,JZ]
其中[XCimage335.tif;%50%50,JZ]表示[XCimage336.tif;%50%50,JZ]方向,[XCimage337.tif;%50%50,JZ]表示[XCimage338.tif;%50%50,JZ]的法线上的单位向量。
最后得到直线[XCimage339.tif;%50%50,JZ]的方程为
[JZ(][XCimage340.tif;%50%50,JZ][JZ)][JY](4)
即[XCimage341.tif;%50%50,JZ]
所以经过Radon变换的Rf函数为[XCimage342.tif;%50%50,JZ]
[JZ(][XCimage343.tif;%50%50,JZ][JZ)][JY](5)
Radon 变换的求逆公式:
[JZ(][XCimage344.tif;%65%65,JZ][JZ)][JY](6)
其中被积函数[XCimage345.tif;%50%50,JZ]为的集合矩阵为我们要求的未知介质的几何形状的图形;[XCimage346.tif;%50%50,JZ]为[XCimage347.tif;%50%50,JZ]的投影,即我们已知的数据;图形[XCimage348.tif;%50%50,JZ]在点[XCimage349.tif;%50%50,JZ]的值称为[XCimage349.tif;%50%50,JZ]的密度,即未知介质每个点的吸收率。[1]
通过MATLAB对Radon变换模型变成进行运算,可以得到结果如下:
未知介质的形状为一椭圆形,椭圆内部又分布着5个不同的小椭圆。
2 确定介质在正方形托盘上的位置
由已有的结论可知:附件提供的是CT系统512个射线与探测单元以第256与257个单元中点为中心旋转180度得到的接收数据,边界数据是CT系统探测器的探测边缘,而不是正方形托盘的的边缘,所以,需要先确定正方形托盘在图形中的具体位置。
由已知数据可知正方形托盘的边长,中心点,和CT系统旋转中心在托盘上的位置。所以旋转中心至中心点的水平距离与竖直距离为:9.2663mm、6.2729mm。
将两个距离分别除以问题一中求解的两个单元点的距离得到两点相差33、19.5个单元点。最终求出关于正方形托盘中心在图形上的位置(289.5,276)。
然后求出正方形边长为[XCimage350.tif;%50%50,JZ],則正方形边界到托盘中心的距离为180.6358。由此可列出求解正方形边界单元点的式子:
[JZ][XCimage351.tif;%50%50,JZ];[XCimage352.tif;%50%50,JZ]。
3 确定给定的托盘上的10个点的吸收率
根据上面求出的正方形托盘在图形中的位置,将10个点在图中的位置用单元点表示出来并找到对应的吸收率。例如,在正方形托盘中坐标为(10,18)的一个点,它在由附件3重建出来的图形中的位置可由以下公式表示:
得到它的吸收率为0。
同理得到十个数据点在图中的位置,由于接受信息与图形形状吸收率已知,所以根据附件的数据用以上方法进行推导,经过求解后,得到图形与成像图基本一致。所以证明此方法有效,模型通过验证。
参考文献:
[1]孔慧华.加速图像重建的迭代算法研究[D].中北大学博士论文,2006:511.
[2]刘成龙.精通MATLAB图像处理[M].北京:清华大学出版社,2015.
基金项目:天津农学院高校教师教育改革创新引导发展项目“高校数学教师创新创业教育能力的提升”(20170301)
*通讯作者:张振荣(1978),女,河北正定人,硕士,讲师,研究方向为数理统计应用及高等数学教学。