Radon变换图像重构

f ( x, y) F (u, v) exp[ j 2 (ux vy)]dudv


二维傅立叶反变换
f ( x, y) F (u, v) exp[ j 2 (ux vy)]dudv


作坐标变换，令： 可得出：
u cos

0
v sin
---我们从CT谈起
• 计算机层析成像（Computed Tomography, CT）是通过 对物体进行不同角度的射线投影测量来获取物体横截面信 息的成像技术。 • CT的核心技术是由投影数据来重建图像的理论，其实质 是由扫描所得到的的投影数据来求出成像平面上每个点的 衰减系数值。
二维投影与CT值
3 4
6
算法举例
5 • 根据反投影算法x1=p5 = 5 x6=p2+p3+p5=18
0
2
0
1
0 0
3 4 6
0
0
5
1
2
0
0
0
5 7 1 3
6 18 10 6
2 12 8 2
3 7 1 5
…
• 平均化处理，除以投影线数目 xi=xi/6
0
0
0
0
原像素值
0.83 1 0.33
反投影重建后
00.5
1.16
Radon逆变换
由测试结果可以看出：第一组 采用30个投影，效果较差；第二组 采用90个投影，效果较好；第三组 采用180个投影，效果很好，与原 始的图像非常接近。 这说明可以通过增加投影的数 目，来提高重建图像的质量。
Fra Baidu bibliotek
I ln 0 (y) dy I a
I (x) p(x) ln (x, y) dy I0 S
这里s表示射线经过的体内距离长度
b
CT任意角度扫描
而对于任意角度扫描，需要用旋转坐标来描述问题，建立置于扫描 ˆ 系统之上的旋转坐标系 x ˆ, y ˆ 的 y ˆ, y ˆ ，即让射线束与旋转坐标系 x 轴平行：
(b)反投影重建图像及星状伪迹
滤波反投影算法
• 滤波反投影法采用先修正、后反投影的做法，其基本方法 是：在某一投影角下取得了投影函数（一维函数）后，对 此一维投影函数作滤波处理，得到一个经过修正的投影函 数；然后再将此修正后的投影函数作反投影运算，得到所 需的密度函数。 • 滤波反投影法重建图像有以下几个步骤： • （1）对某一角度下的投影函数作一维傅立叶变换； • （2）对（1）的变换结果乘上一维权重因子； • （3）对（2）的加权结果作一维逆傅立叶变换； • （4）用（3）中得出的修正过的投影函数做直接反投影； • （5）改变投影角度，重复（1）~（4）的过程，直到完成 全部180度的反投影。


g ' (R) ( x cos y sin R) dR

g (R)
'


F ( , ) e j 2 R d
表示对投影函数的Fourier变换进行滤波变换，其 中 是滤波函数。
由傅立叶变换性质可知．频域中的滤波运算可等效地在空域中用卷积运算来完成 所以：
0.06 0.5
3
1.66 1
2
1.33 0.33
1.16
0.16 0.83
再除以投影线数，平均化
断层平面中某一点的密度值可看作这一平面内所有经 过该点的射线投影之和的平均值
伪迹
0 0 0 0
0 5 1 0
0 2 0 0
0 0 0 0
0.83 1.16 0.06 0.5
1 3 1.66 1
0.33 2 1.33 0.33
f ( x, y ) d
0

F ( , ) e j 2 ( x cos y sin ) d d




j 2 R F ( , ) e d ( x cos y sin R) dR
d
0
0.5 1.16 0.16 0.83
原像素值
再除以投影线数，平均化
• 反投影重建后，原来为0的点不再为0，形成伪迹
星状伪迹
• 我们考虑孤立点源反投影重建，中心点A 经n条投影线投影后，投影值均为1： p1=p2=...=pn=1 因此重建后
0 0 0 0 1 0 0 0 0
1 f A ( p1 p2 ... pn ) 1 n
（2）S-L滤波函数 与R-L滤波函数不同的是，S-L滤波函数它的关键是把 频域的陡峭截止改成缓慢截止。 用S-L滤波函数重建的图像中振荡相应较小，对含噪声 的数据重建出来的图像质量也较R-L滤波函数重建的图像 质量要好。但是，S-L滤波函数重建的图像在高频响应方 面不如R-L滤波函数好，这是因为S-L滤波函数在高频段偏 离了理想的滤波函数
Radon变换（续）
• 由此，可以沿任意角度 计算函数的投影，计算图像f(x,y)在任意角度 的Radon变换。
中心切片定理
• 密度函数在某一方向上的投影函数的一维傅立叶变换函数 是原密度函数的二维傅立叶变换函数在平面上沿同一方向 且过原点的直线上值。
滤波反投影算法的原理
1、在不同的角度下取得足够多的投影数据(Radon变换) 2、将这些投影数据做一维的Fourier变换，那么变换后的这些数据将充满 整个(u,v)平面。（许多过原点成不同夹角的直线） 3、也就是说，F(u,v)的全部值都为已知，那么我们将其做一次二维的 Fourier逆变换就可以得到原始的衰减系数函数f(x,y)
f ( x, y) d g (R) h(R) ( x cos y sin R) dR
0


式中h(R)为滤波函数纠的空域形式
• 所以要实现对投影数据实现图像重建，可以采取两步：首先将投影 数据和响应脉冲滤波器进行卷积，然后由式对不同旋转角θ求和， 就能实现图像重建。这就是卷积法进行图像重建的基本思路和方法。 • 卷积可看作一种滤波手段，卷积投影相当于对数据先滤波再将结果 逆投影回来，这样可以使模糊得到校正。
1、我们假设切片（物体横截面，断面）无限薄。 • 当强度为 I 0 的x-ray通过吸收 率为μ(x,y)的均匀吸收物体，由 2、我们认为,一幅图像在任意点(x,y) 上的灰度值 正比于那个点的相对线性衰减系数μ(x,y)。 于均匀吸收，则I必是指数下降， 则有
b I I 0 exp (y) dy a
图像重建概述
• 图像重建是图像处理中的一个重要分支，广泛地应用于物 体内部结构图像的检测和观察中，它是一种无损检测技术。 • 关于图像处理的一些基本内容，如对图像的几何处理，图 像的增强，还有复原等，均是从图像到图像，即输入的原 始数据是图像，处理后输出的仍是图像。而图像重建是从 数据到图像。 • 图像重建的三种常用检测模型：透射模型、发射模型、反 射模型
滤波函数
• 滤波函数的选取是滤波反投影法的关键问题 （1）R-L滤波函数 由于在频域中用矩形函数截断了滤波函数，在相应的 空域中造成振荡响应 ，重建的图像质量也不够满意
对应的频域形式为：
H ( ) rect ( / 20 )
理想的滤波函数 它是 在高频的权重很大，低频 的权重很小，所以高频噪 声就会很大，所以我们才 要对其进行修正
ˆ I (x) ˆ ) (x, ˆ y) ˆ dy ˆ ln( p(x, ) 是离散值，是测出值！ 经坐标系旋转变换后可得： I 0 s
ˆ ) 是测量值， 所以所谓投影 p (x, 是吸收系数沿着射线经过直线的积 分。实际上的问题是沿着若干条直 线的积分估算值来计算μ(x,y)值。
图像的细节对应的是高频部分， 轮廓对应的是图像的低频部分， 所以因为没有滤波，细节部分 恢复的不好，呈现很“模糊” 的情况
二、投影数据的多少对图像重建效果的影响
• 一个典型实例： • 在matlab图像处理工具箱中，有 一个phantom函数，可以用来创 建头部的剖视图，首先创建一个 头部的256×256剖视图，然后 分别计算3组不同的Radon变换， 第一组采用30个投影，第二组采 用90个投影，第三组采用180个 投影，用以比较采用不同组数的 投影参数重建的图像与原始图像 的差别。
所以θ角每旋转1度就可以取一组 投影数据，可得到180组不同的 投影。CT就是在收集各角度θ的 投影数据后，利用重建算法处理 得到物体的图像。
Radon变换
• Radon变换是计算图像在某一指定角度射线方向上的投影 的变换方法。二维函数f(x,y)的投影是其在确定方向上的线 积分，如下图所示，二维函数f(x,y)在水平方向的线积分就 是f(x,y)在y轴上的投影，二维函数f(x,y)在垂直方向的线性 积分就是f(x,y)在x轴上的投影。
对应的频域形式为：
H S L ( ) sin c( )rect ( ) 2 0 2 0
一、滤波反投影matlab实现
•
• • • • • • • • • • • •
%P=imread('lena.jpg'); P = phantom(256); %P= rgb2gray(O); R = radon(P,0:179); I0 = iradon(R,0:179,'linear','Ram-Lak'); I1 = iradon(R,0:179,'linear','Shepp-Logan'); I2=iradon(R,0:179,'linear','cosine'); I3 = iradon(R,0:179,'linear','none'); subplot(2,3,1), imshow(P), title('Original') subplot(2,3,2), imshow(I0,[]), title('FBP R-L') subplot(2,3,3), imshow(I1,[]), title('FBP S-L') subplot(2,3,4), imshow(I2,[]), title('FBP cosine') subplot(2,3,5), imshow(I3,[]), title('Unfiltered BP')
而其他点均为1/n 这类伪迹称为星状伪迹
1/n 1/n
1/n 1
1/n 1/n
1/n
1/n
1/n
星状伪迹
• 产生星状伪迹的原因在于：反投影重建的本质是把取 自有限物体空间的射线投影均匀地回抹(反投影)到射 线所及的无限空间的各点之上，包括原先像素值为零 的点（其实就是投影数据少产生的！！！）
(a)孤立点源
反投影算法举例
• 基本原理是将所测得的投影值按其原路径平均的分配到每 一点上，各个方向上投影值反投影后，在影像处进行叠加， 从而推体出原图像。 • 而滤波却是要投影函数的一维Fourier加上权重因子。
算法举例
5 2 1
p1 x3 x7 x11 x15 2 p2 x2 x6 x10 x14 6 p3 x5 x6 x7 x8 7 p4 x9 x10 x11 x12 1 p5 x1 x6 x11 x16 5 p6 x4 x7 x10 x13 3