2017数学建模A题二等奖论文
相关主题
- 1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。
2.利用(1)中得到的标定参数，确定附件 3 中的未知介质在正方形托盘中的 位置、几何形状和吸收率等信息。另外，请具体给出图 3 所给的 10 个位置处的 吸收率，相应的数据文件见附件 4。
3. 附件 5 是利用上述 CT 系统得到的另一个未知介质的接收信息。利用(1) 中得到的标定参数，给出该未知介质的相关信息。另外，请具体给出图 3 所给的 10 个位置处的吸收率。
到下表（见附件“180 个 ”.xls）：
n 角度θ
n 角度θ
n 角度θ
n 角度θ
n 角度θ
n 角度θ
n 角度θ
n 角度θ
n 角度θ
n
角度θ n
角度θ n
角度θ n
角度θ n
角度θ n
角度θ n
角度θ n
角度θ n
角度θ
1 29.4984
11 39.9991
21 53.9763
31 61.8652
41 70.7043
2
里叶级数分析求出 k 值。联立两次 K 探测器接收的光线的直线方程，交点就是 旋转中心，投影数据分析求探测器间距、180 个方向。 2.2 问题二分析
问题二是一个已知角度下的投影，还原原函数密度的问题。根据中心切片 定理：函数密度在某一方向上的投影的一维傅里叶变换函数是原函数密度的二 维傅里叶变换函数在同一角度且过原点直线上的值。据此定理结合附件数据求 解变换后的傅里叶函数在每一个角度下的频域函数，再取傅立叶逆变换。 2.3 问题三分析
三、模型假设
1.忽略空气对 X 射线的吸收率。 2.不考虑光的折射和反射。 3.假设所测量的数据均为正常测量的情况 ，不考虑测量过程的抖动因素。 4.假设 CT 的滤波器是理想的。
四、符号说明
θ :光线的方向角，光线矢量与 y 轴正方向夹角,逆时针旋转为增大的正方向 β :探测器和 x 轴正方向夹角 i:探测器编（当探测器旋转到平行于正方形上方时，从左向右以 0―512 编号） K:距离旋转中心最近的探测器编号 Q:吸收率 d:椭圆中心和圆心之间对应多个探测器的距离 l :相邻探测器间距离
51 80.1570
61 89.9579
71 99.4865
81 108.5526
91 117.3405
101 125.5902
111 139.4572
121 149.1036
131 158.7065
141 167.0810
151 182.2443
161 191.2388
171 199.2646
2 31.5744ຫໍສະໝຸດ Baidu
β 180o sinβ = d l
=arcsin d 180o l
d | ymax1 ymax2 |gl
经过 K 探测器的光路直线方程： yk tan( 180o)gxk [b a tan( 180o)]
求K 由几何关系可知，对于单一的几何图形来说，在任一确定角度上，通过椭圆中 心的弦和圆的圆心的弦是最长的，即对应的探测器接收数据最大。将附件 2 的数
由于两条包络是围绕同一旋转中心得到的，所以它们的相位相等，将两条包 络线进行取平均得到 y。
图 5-5 椭圆和圆连线中点在接收器上的投影位置图
7
假使探测器不动，模板相对探测器旋转，圆心绕旋转中心的轨迹，分时投影到 探测器上的曲线是余弦函数曲线形状，旋转后的结果图如下：（旋转动图见附件）。
5-6 几何中心点围绕旋转中心旋转产生的余弦图像
CT 系统安装时往往存在误差，从而影响成像质量，因此需要对安装好的 CT 系统进行参数标定，即借助于已知结构的样品（称为模板）标定 CT 系统的参数， 并据此对未知结构的样品进行成像。
1.2 问题的提出
1.在正方形托盘上放置两个均匀固体介质组成的标定模板，其中每一点的数 值反映了该点的吸收强度，这里称为“吸收率”。请根据这一模板及其接收信息， 确定 CT 系统旋转中心在正方形托盘中的位置、探测器单元之间的距离以及该 CT 系统使用的 X 射线的 180 个方向。
3
y1,y2:包络线位置函数 ymax1 、 ymax 2 ：椭圆中心和圆心对应探测器的编号 n：测量次数（1≤n≤180） w：傅里叶展开函数的基频率 （a,b）：旋转中心坐标 f1 、 f2 ：圆的 8mm 直径对应 28 个间距、29 个间距的探测器，出现频率 f (x, y) ：未知物体的密度函数
同理可得：
图 5-2（ 90o 180o 旋转过程图） β 90o
5
sinβ = d l
=arcsin d 90o l
经过 K 探测器的光路直线方程： yk tan( 90o)gxk [b a tan( 90o)] （3）180o 270o
图 5-3（180o 270o旋转过程图） 同理可得
6
据，对每一列数据使用牛顿迭代法找到极大值 ymax1 、 ymax2 极大值对应的探测器 编号即为穿过圆心的光线对应的探测器。横轴代表 180 次测量的转次（暂未知角 度），纵轴代表一次测量时，每一个的传感器测量状态（状态：传感器有无吸收 率，1≤i≤512）。
图 5-4 圆与椭圆双旋转图 曲线拟合方法：
12 40.5428
22 54.8409
32 62.2631
42 71.4474
52 81.2255
62 89.9579
72 100.1999
82 109.6686
92 118.5342
102 126.4595
112 140.0010
122 150.5016
对于问题二和问题三，本文根据中心切片定理，建立起投影物体的密度与其 二维傅里叶变换后的函数的关系模型。在每一个确定角度下，分别使用直接反投 影法和滤波反投影法
f (x, y) | f (R cos, Rsin) | G (R)
G () F(, ) F [ f (x, y) | ]
将变换后的函数作傅立叶逆变换就可以得到原物体的密度函数，并进行插值 填充由 MATLAB 生成吸收率矩阵。问题二的介质几何中心为（43.1641,47.4609）， 形状大致是椭圆形的嵌套（见图 5-9）。问题二的 10 个点吸收率（附件 4 坐标 点自上到下）分别为：0、0、2.0053、1.9714、1.9613、1.9644、0、0、0、0。 图片见本文 4.2.2 模型求解部分。问题三的 10 个点吸收率（附件 4 坐标点自上 到下）分别为：0、0、0.6479、0、0、4.7242、、0、10.5331、0、0。图片轮廓 形状类似于组织细胞图，图片见本文 4.2.3 模型求解部分。
1
1
yn a0 a1 cos nw b1 sin nw
其中 y1、y2、y3 已知，a0、a1、b1 是未知。
用上面求出的θ 值对应的每一次位置，计算出每一次经过 K 探测器的光路直线方程， 求两两直线的焦点，求出旋转中心。由平面几何分析可得，每一次经过 K 探测器的光路 直线方程为：
yk tan gxk (b a tan )
对它利用傅里叶级数展开可以得到直流分量和余弦分量，直流分量的大小就是 K。
yn a0 a1 cos nw b1 sin nw
进行一级傅立叶级数展开：
y1 a0 a1 cos w b1 sin w
……
y2 a0
a cos 2w b sin 2w
1
1
y3 a0
a cos 3w b sin 3w
4. 分析(1)中参数标定的精度和稳定性。在此基础上自行设计新模板、建立 对应的标定模型，以改进标定精度和稳定性，并说明理由。
二、问题分析
2.1 问题一分析
以正方形托盘左下角为坐标原点建立平面直角坐标系。讨论光线方向角和 探测器与 x 轴夹角的关系，根据投影方向角、旋转中心的坐标、旋转中心到最 近的探测器的标号之间的关系，建立光路方程。光路方程被模板的几何图形所 截长度正比于所给的吸收率。经过旋转，圆心的运动轨迹随角度的变化是余弦 图并且其直流分量 k 即为旋转中心所对应的探测器的编号，对其曲线拟合和傅
yn a0 a1 cos nw b1 sin nw ， yk tan gxk (b a tan )
求解模型时，将附件的探测器接收数据和带入，利用光路方向角与椭圆中心 和圆心之间的距离的几何关系求出 180 个方向角度，旋转中心坐标 （50.2775,56.1050）,相邻探测器间距离为 0.2775mm。
问题三与问题二分析和求解类似，并且都是建立在第一问的已知 180 次测量 角度上的，由于上述处理方法非特殊取值的方法，因此该方法具有普遍适用型， 此处不再赘述。 2.4 问题四的分析
模板对稳定性的影响在于模板在正方形托盘位置的分布，由几何概型可知， 模板分布越集中于托盘中部，被探测器检测到的概率越大，稳定性越高。因此 本文设计出的模板模型是圆心位于正方形托盘重心处的圆。
G (R) ：密度函数在 方向上的投影函数 G () ：密度函数在 方向上的投影函数一维傅里叶变换函数
F(, ) ： f (x, y) 的二维傅里叶变换函数 fW (x, y) ：重建后的密度函数 H(w) ：滤波器的频域表达式 FWH (, ) ： F(, ) 滤波后的结果 (n) ：每一次测量角度的弧度制 t1、t2 :存储相邻两次旋转的稳定因子矩阵 lt ：程序中的矩阵的长度 N：值为 512，代表 512 个探测器
联立任意两个角度下的经过 K 探测器的光路直线方程，计算出多组两直线的交点， 利用统计求取均值，所得就是旋转中心点。
由于圆在探测器上投影的永远是直径，为定值，因此在每一个角度确定的情况下，通
过圆直径所对应的探测器的间距个数，求出 l
l
8 29
f1
8 28
f2
8
5.1.2 问题 1 模型的求解
根据数据按列迭代求取极大值，从而确定过椭圆中心和圆心每一个θ 下对应的 d 的大小，最后利用 4.1.2 中三个情况下的θ 数学表达式，求出每一个角度，得
关键词：余弦轨迹 滤波反投影法 中心切片定理,傅里叶变换 傅里叶级数分 析
1
一、问题重述
1.1
CT(Computed Tomography)可以在不破坏样品的情况下，利用样品对射线能 量的吸收特性对生物组织和工程材料的样品进行断层成像，由此获取样品内部的 结构信息。一种典型的二维 CT 系统如图 1 所示，平行入射的 X 射线垂直于探测 器平面，每个探测器单元看成一个接收点，且等距排列。X 射线的发射器和探测 器相对位置固定不变，整个发射-接收系统绕某固定的旋转中心逆时针旋转 180 次。对每一个 X 射线方向，在具有 512 个等距单元的探测器上测量经位置固定不 动的二维待检测介质吸收衰减后的射线能量，并经过增益等处理后得到 180 组接 收信息。
五、模型的建立与求解
5.1.1 问题 1 的模型建立 建立平面直角坐标系求解θ 和经过 K 探测器的光路直线方程（分三类讨论） （1） 0 90o
4
图 5-1（ 0 90o 旋转过程图）
由几何关系可知：β =θ
cosβ d l
arccos
d l
经过 K 探测器的光路直线方程： yk tan gxk (b a tan ) （2） 90o 180o
对于问题四，问题 1 中的相邻探测器间距离参数，是影响成像精度的重要参 数，距离越小精确度越高。旋转中心在托盘上的位置和探测器的个数，是影响成 像稳定性的重要参数，增加单位长度的探测器个数，并且调整旋转中心和探测器 中心相对位置，使得探测器每转一个角度，接收到信息探测器尽可能在中部（避 免下次转过一个角度时，介质的投影信息溢出光屏）。探测器绕旋转中心扫过的 方向角度尽量稍微大于 180°，使信息采集的更加完备，从而增加了 CT 工作的 稳定性。
CT 系统参数标定及成像
摘要
目前，如何通过 CT 机的探测器已知的投影信息还原出物体几何形状以及内 部结构，探测已知模块得到的数据推导 CT 机的参数以及校正，是研究精密仪器 的重要方面。本文通过建立数学模型，给出了通过分析表格数据来解决上述所有 问题的方法。
对于问题一，以正方形左下角为原点建立平面直角坐标系对方向角进行几何 抽象，通过介质上任意一点在转动过程中在光屏上的投影点的运动轨迹都是在旋 转中心的投影点的上下波动这一几何特性，引入三角函数、直流分量和交流分量 并与几何方法结合，建立了已知模板接收信息、吸收率、CT 参数、余弦轨迹的 傅里叶展开之间关系的综合模型。
3. 附件 5 是利用上述 CT 系统得到的另一个未知介质的接收信息。利用(1) 中得到的标定参数，给出该未知介质的相关信息。另外，请具体给出图 3 所给的 10 个位置处的吸收率。
到下表（见附件“180 个 ”.xls）：
n 角度θ
n 角度θ
n 角度θ
n 角度θ
n 角度θ
n 角度θ
n 角度θ
n 角度θ
n 角度θ
n
角度θ n
角度θ n
角度θ n
角度θ n
角度θ n
角度θ n
角度θ n
角度θ n
角度θ
1 29.4984
11 39.9991
21 53.9763
31 61.8652
41 70.7043
2
里叶级数分析求出 k 值。联立两次 K 探测器接收的光线的直线方程，交点就是 旋转中心，投影数据分析求探测器间距、180 个方向。 2.2 问题二分析
问题二是一个已知角度下的投影，还原原函数密度的问题。根据中心切片 定理：函数密度在某一方向上的投影的一维傅里叶变换函数是原函数密度的二 维傅里叶变换函数在同一角度且过原点直线上的值。据此定理结合附件数据求 解变换后的傅里叶函数在每一个角度下的频域函数，再取傅立叶逆变换。 2.3 问题三分析
三、模型假设
1.忽略空气对 X 射线的吸收率。 2.不考虑光的折射和反射。 3.假设所测量的数据均为正常测量的情况 ，不考虑测量过程的抖动因素。 4.假设 CT 的滤波器是理想的。
四、符号说明
θ :光线的方向角，光线矢量与 y 轴正方向夹角,逆时针旋转为增大的正方向 β :探测器和 x 轴正方向夹角 i:探测器编（当探测器旋转到平行于正方形上方时，从左向右以 0―512 编号） K:距离旋转中心最近的探测器编号 Q:吸收率 d:椭圆中心和圆心之间对应多个探测器的距离 l :相邻探测器间距离
51 80.1570
61 89.9579
71 99.4865
81 108.5526
91 117.3405
101 125.5902
111 139.4572
121 149.1036
131 158.7065
141 167.0810
151 182.2443
161 191.2388
171 199.2646
2 31.5744ຫໍສະໝຸດ Baidu
β 180o sinβ = d l
=arcsin d 180o l
d | ymax1 ymax2 |gl
经过 K 探测器的光路直线方程： yk tan( 180o)gxk [b a tan( 180o)]
求K 由几何关系可知，对于单一的几何图形来说，在任一确定角度上，通过椭圆中 心的弦和圆的圆心的弦是最长的，即对应的探测器接收数据最大。将附件 2 的数
由于两条包络是围绕同一旋转中心得到的，所以它们的相位相等，将两条包 络线进行取平均得到 y。
图 5-5 椭圆和圆连线中点在接收器上的投影位置图
7
假使探测器不动，模板相对探测器旋转，圆心绕旋转中心的轨迹，分时投影到 探测器上的曲线是余弦函数曲线形状，旋转后的结果图如下：（旋转动图见附件）。
5-6 几何中心点围绕旋转中心旋转产生的余弦图像
CT 系统安装时往往存在误差，从而影响成像质量，因此需要对安装好的 CT 系统进行参数标定，即借助于已知结构的样品（称为模板）标定 CT 系统的参数， 并据此对未知结构的样品进行成像。
1.2 问题的提出
1.在正方形托盘上放置两个均匀固体介质组成的标定模板，其中每一点的数 值反映了该点的吸收强度，这里称为“吸收率”。请根据这一模板及其接收信息， 确定 CT 系统旋转中心在正方形托盘中的位置、探测器单元之间的距离以及该 CT 系统使用的 X 射线的 180 个方向。
3
y1,y2:包络线位置函数 ymax1 、 ymax 2 ：椭圆中心和圆心对应探测器的编号 n：测量次数（1≤n≤180） w：傅里叶展开函数的基频率 （a,b）：旋转中心坐标 f1 、 f2 ：圆的 8mm 直径对应 28 个间距、29 个间距的探测器，出现频率 f (x, y) ：未知物体的密度函数
同理可得：
图 5-2（ 90o 180o 旋转过程图） β 90o
5
sinβ = d l
=arcsin d 90o l
经过 K 探测器的光路直线方程： yk tan( 90o)gxk [b a tan( 90o)] （3）180o 270o
图 5-3（180o 270o旋转过程图） 同理可得
6
据，对每一列数据使用牛顿迭代法找到极大值 ymax1 、 ymax2 极大值对应的探测器 编号即为穿过圆心的光线对应的探测器。横轴代表 180 次测量的转次（暂未知角 度），纵轴代表一次测量时，每一个的传感器测量状态（状态：传感器有无吸收 率，1≤i≤512）。
图 5-4 圆与椭圆双旋转图 曲线拟合方法：
12 40.5428
22 54.8409
32 62.2631
42 71.4474
52 81.2255
62 89.9579
72 100.1999
82 109.6686
92 118.5342
102 126.4595
112 140.0010
122 150.5016
对于问题二和问题三，本文根据中心切片定理，建立起投影物体的密度与其 二维傅里叶变换后的函数的关系模型。在每一个确定角度下，分别使用直接反投 影法和滤波反投影法
f (x, y) | f (R cos, Rsin) | G (R)
G () F(, ) F [ f (x, y) | ]
将变换后的函数作傅立叶逆变换就可以得到原物体的密度函数，并进行插值 填充由 MATLAB 生成吸收率矩阵。问题二的介质几何中心为（43.1641,47.4609）， 形状大致是椭圆形的嵌套（见图 5-9）。问题二的 10 个点吸收率（附件 4 坐标 点自上到下）分别为：0、0、2.0053、1.9714、1.9613、1.9644、0、0、0、0。 图片见本文 4.2.2 模型求解部分。问题三的 10 个点吸收率（附件 4 坐标点自上 到下）分别为：0、0、0.6479、0、0、4.7242、、0、10.5331、0、0。图片轮廓 形状类似于组织细胞图，图片见本文 4.2.3 模型求解部分。
1
1
yn a0 a1 cos nw b1 sin nw
其中 y1、y2、y3 已知，a0、a1、b1 是未知。
用上面求出的θ 值对应的每一次位置，计算出每一次经过 K 探测器的光路直线方程， 求两两直线的焦点，求出旋转中心。由平面几何分析可得，每一次经过 K 探测器的光路 直线方程为：
yk tan gxk (b a tan )
对它利用傅里叶级数展开可以得到直流分量和余弦分量，直流分量的大小就是 K。
yn a0 a1 cos nw b1 sin nw
进行一级傅立叶级数展开：
y1 a0 a1 cos w b1 sin w
……
y2 a0
a cos 2w b sin 2w
1
1
y3 a0
a cos 3w b sin 3w
4. 分析(1)中参数标定的精度和稳定性。在此基础上自行设计新模板、建立 对应的标定模型，以改进标定精度和稳定性，并说明理由。
二、问题分析
2.1 问题一分析
以正方形托盘左下角为坐标原点建立平面直角坐标系。讨论光线方向角和 探测器与 x 轴夹角的关系，根据投影方向角、旋转中心的坐标、旋转中心到最 近的探测器的标号之间的关系，建立光路方程。光路方程被模板的几何图形所 截长度正比于所给的吸收率。经过旋转，圆心的运动轨迹随角度的变化是余弦 图并且其直流分量 k 即为旋转中心所对应的探测器的编号，对其曲线拟合和傅
yn a0 a1 cos nw b1 sin nw ， yk tan gxk (b a tan )
求解模型时，将附件的探测器接收数据和带入，利用光路方向角与椭圆中心 和圆心之间的距离的几何关系求出 180 个方向角度，旋转中心坐标 （50.2775,56.1050）,相邻探测器间距离为 0.2775mm。
问题三与问题二分析和求解类似，并且都是建立在第一问的已知 180 次测量 角度上的，由于上述处理方法非特殊取值的方法，因此该方法具有普遍适用型， 此处不再赘述。 2.4 问题四的分析
模板对稳定性的影响在于模板在正方形托盘位置的分布，由几何概型可知， 模板分布越集中于托盘中部，被探测器检测到的概率越大，稳定性越高。因此 本文设计出的模板模型是圆心位于正方形托盘重心处的圆。
G (R) ：密度函数在 方向上的投影函数 G () ：密度函数在 方向上的投影函数一维傅里叶变换函数
F(, ) ： f (x, y) 的二维傅里叶变换函数 fW (x, y) ：重建后的密度函数 H(w) ：滤波器的频域表达式 FWH (, ) ： F(, ) 滤波后的结果 (n) ：每一次测量角度的弧度制 t1、t2 :存储相邻两次旋转的稳定因子矩阵 lt ：程序中的矩阵的长度 N：值为 512，代表 512 个探测器
联立任意两个角度下的经过 K 探测器的光路直线方程，计算出多组两直线的交点， 利用统计求取均值，所得就是旋转中心点。
由于圆在探测器上投影的永远是直径，为定值，因此在每一个角度确定的情况下，通
过圆直径所对应的探测器的间距个数，求出 l
l
8 29
f1
8 28
f2
8
5.1.2 问题 1 模型的求解
根据数据按列迭代求取极大值，从而确定过椭圆中心和圆心每一个θ 下对应的 d 的大小，最后利用 4.1.2 中三个情况下的θ 数学表达式，求出每一个角度，得
关键词：余弦轨迹 滤波反投影法 中心切片定理,傅里叶变换 傅里叶级数分 析
1
一、问题重述
1.1
CT(Computed Tomography)可以在不破坏样品的情况下，利用样品对射线能 量的吸收特性对生物组织和工程材料的样品进行断层成像，由此获取样品内部的 结构信息。一种典型的二维 CT 系统如图 1 所示，平行入射的 X 射线垂直于探测 器平面，每个探测器单元看成一个接收点，且等距排列。X 射线的发射器和探测 器相对位置固定不变，整个发射-接收系统绕某固定的旋转中心逆时针旋转 180 次。对每一个 X 射线方向，在具有 512 个等距单元的探测器上测量经位置固定不 动的二维待检测介质吸收衰减后的射线能量，并经过增益等处理后得到 180 组接 收信息。
五、模型的建立与求解
5.1.1 问题 1 的模型建立 建立平面直角坐标系求解θ 和经过 K 探测器的光路直线方程（分三类讨论） （1） 0 90o
4
图 5-1（ 0 90o 旋转过程图）
由几何关系可知：β =θ
cosβ d l
arccos
d l
经过 K 探测器的光路直线方程： yk tan gxk (b a tan ) （2） 90o 180o
对于问题四，问题 1 中的相邻探测器间距离参数，是影响成像精度的重要参 数，距离越小精确度越高。旋转中心在托盘上的位置和探测器的个数，是影响成 像稳定性的重要参数，增加单位长度的探测器个数，并且调整旋转中心和探测器 中心相对位置，使得探测器每转一个角度，接收到信息探测器尽可能在中部（避 免下次转过一个角度时，介质的投影信息溢出光屏）。探测器绕旋转中心扫过的 方向角度尽量稍微大于 180°，使信息采集的更加完备，从而增加了 CT 工作的 稳定性。
CT 系统参数标定及成像
摘要
目前，如何通过 CT 机的探测器已知的投影信息还原出物体几何形状以及内 部结构，探测已知模块得到的数据推导 CT 机的参数以及校正，是研究精密仪器 的重要方面。本文通过建立数学模型，给出了通过分析表格数据来解决上述所有 问题的方法。
对于问题一，以正方形左下角为原点建立平面直角坐标系对方向角进行几何 抽象，通过介质上任意一点在转动过程中在光屏上的投影点的运动轨迹都是在旋 转中心的投影点的上下波动这一几何特性，引入三角函数、直流分量和交流分量 并与几何方法结合，建立了已知模板接收信息、吸收率、CT 参数、余弦轨迹的 傅里叶展开之间关系的综合模型。