在路上●我的年青●悠忽两年

题记：悠忽两年，今天又到了我的生日，离开校园几近五年，大学毕业也快八年。偶然的机会帮着师弟分析一个问题，仔细整理后，发现用到的都是大学阶段所学。蓦然发觉，这些东西装入大脑十年有余，却不曾倒出来过。特发此文，回忆那些再也回不去的时光，想念我的老师、同学，缅怀青春！

MATLAB RADON及其相关函数分析

摘要：MATLAB Phantom函数产生的Shepp-Logan模型，和人脑结构相似，是医学影像学中用于CT图像仿真重建领域的经典模型。本文首先详细分析Shepp-Logan模型的结构，然后介绍对该模型数据求取Radon变换的两种方法。Radon变换后的数据求取Iradon变换，需要对数据预处理，预处理中用到滤波和卷积函数，第二节对其详述。模型数据的Fanbeam和IFanbeam变换，重要的是平行射线转扇形射线、以及扇形射线转平行射线的转换，而转换的过程就是插值，第三节详述常用的分段三次立方插值和三次样条插值方法。通篇都是分析，没有创新，涉及到的M文件有phantom，radon，iradon，fft，ifft，conv，pchip，spline。其中快速傅里叶变换实现频域滤波的物理意义，以及时域卷积的物理意义，并没有解释清楚，需要继续探讨。

关键字：Phantom函数C实现；Iradon函数C实现；离散傅里叶变换；离散卷积；埃尔米特插值时一阶导数求法；三次样条插值时非扭结条件；三对角矩阵的追赶法

1 投影变换

1.1 Radon变换

实际就是求取物体的截面(与接收器平行的那个截面)，在平行射线穿透下的投影值。本节先分析图像仿真重建用到的Shepp-Logan模型，然后介绍模型生成的数据，求取Randon变换的两种方法。

1.1.1 MATLAB phantom 函数

MATLAB phantom 函数用于产生头骨幻影数据，该数据是一个N*N的矩阵，N默认是256，可以用来验证randon和irandon、以及其他二维的图像重建算法的数值准确度。phantom函数产生的Shepp-Logan模型，和人脑结构相似，是医学影像学中最先用于CT图像重建仿真领域的经典模型，实际处理中为了视觉观察上的方便，一般选取修正的Shepp-Logan模型。phantom产生的图像，如图1-1所示，其中(a)图是模型，(b)图是二维灰度图像，(c)图是模拟的人脑结构模型。

(a)

(b)

(c)

图1-1 (a) 头骨幻影模型；(b) 头骨幻影灰度图像； (c) Shepp-Logan模拟的人脑模型

phantom函数产生的图像，有10个椭圆组成(圆是椭圆的一种特殊情况)，如图1-1(a)中的a，b，c，d … i，j，每个椭圆有6个特征数据，默认的数据如表1-1所示[1]。其中灰度值范围是[0 1]，图像显示时对应[0 255]，表示椭圆内部都是这个值，重复的地方要叠加，密度值的叠加就是灰度值；x轴和y轴的坐标范围都是从[-1 1]；旋转角度是度数，不是弧度。另外要注意的是，图1-1(a)倾斜的椭圆c和d，它们两个是由中心在坐标原点的椭圆，旋转一定角度，然后再将旋转后的椭圆，沿着x轴和y轴平移得来的，即先旋转后平移。

根据表1-1，结合图1-1(a)，则图1-1(b)中的灰度图像，灰度取值有6种情况，分别是0, 0.1, 0.2, 0.3, 0.4, 1。椭圆a外面的、以及椭圆c和d内部的都是0，椭圆d与椭圆e相交的地方、椭圆d与椭圆f相交的地方都是0.1，椭圆b内部其它地方是0.2，椭圆e、f、g、h、i、j内部的都是0.3，椭圆e和f相交的地方是0.4，椭圆a和b相交的地方是1。

表1-1 修正的Shepp-Logan模型中椭圆特征数据

(a)

(b)

图1-2 (a) 生成头骨幻影模型的MATLAB代码；(b) 生成头骨幻影数据的C代码