医学图像重建 1 | Radon 变换、过滤反投影算法、中心切片定理
- 转自微信公众号:机器学习炼丹术
- 笔记:陈亦新
概述
医学图像重建的目的就是得到上图的f(x,y)的图像。我们只能获取到投影的数据,也就是右边的sensor检测到的强度信息。当然上图来看,是把一个2D的图像投影成了1D的数据,那么这样肯定是无法复原的。
在投影的过程中,并不是上述这一个角度。上述的投影角度为0,是水平从左到右的。假设我们对角度旋转180度,每旋转1度都进行一次投影。那么最终我们会有180个1D的投影数据,然后如何从这些1D投影数据还原2D的原始图像就是我们所说的重建算法。
Radon变换
这个变换讲述的就是将2D物体投影成1D的过程。2D的两个维度记作x和y,1D的数据只有1个维度,我们记作s。但是我们还需要考虑这个radon变成的1D其实是在某一个特定投影角度下的1D数据,所以其实上是还要加上角度的变量.
用作数学的表示,那就是radon变换就是将f(x,y)变成的过程。这个公式的推导我是这样理解的,在某一个特定角度下*(水平方向投影)的s和f(x,y)存在这个关系:
如果考虑到角度,那么就变成:
现在我们来理解什么是垂直sensor的直线上。
这个图展示了一个投影角度为的通用情况。的s的物理含义是直线与2D坐标原点的垂直距离,也是1D投影距离1D投影坐标原点的距离,就是上图中的p。
现在我们需要吧垂直sensor的直线上约束转换成数学的表示,直线的表示常见就是y-kx-b=0的形式,但是我们需要用p和来表示直线,因此直线可以表示为:
那么上面公式就变成
但是这还不够好,你在积分下面加个约束,那后续做什么运算都不行。就像是带约束的优化问题有拉格朗日算子法一样,这里我们引入狄拉克分布函数。看下百度的结果:
再看下书中的定义:
再看下我的理解:
- 狄拉克函数是一个概念,理想情况下就是一个类似在0的地方无穷大,非0地方等于0的脉冲函数。当然这种函数可能不可导或者没有什么很好的性质。因此我们可以用一些性质比较好的函数来模拟这种效应。也就是高斯函数。公式为:
其中n越大,曲线越来越窄,图像如下:
为什么这样定义狄拉克函数呢?首先,这样的定义是符合狄拉克函数的广义概念的,所以他是狄拉克函数的一种定义。然后就是这样定义会有很多很方便的性质,后续用到的时候我再做介绍。
我们积分约束那里,现在我们可以写成这样的形式:
因为当xcos\theta+ysin\theta-p不等于0的时候,其实可以看作是0(我的理解哈).
至此,我们理解了什么是radon变换,是一个多角度投影的正向过程。
中心切片定理
中心切片定理是断层扫描成像的理论基础。这个定理还可以叫做:投影切片定理和傅里叶中心切片定理。二维图像的中心切片定义指出:二维图像f(x,y)的角度的投影的傅里叶变换等于函数f(x,y)的傅里叶变换同样沿着\theta角度过原点的片段
看下图:
-
右边黑乎乎的是f(x,y)的傅里叶变化图像。
-
f(x,y)沿着某一个方向投影得到绿色的1D分布,这个是radon过程。
-
然后把1D投影分布做傅里叶变换得到红色1D频域分布。
-
这个中心切片定理关键就是说,这个红色的1D分布,其实是等于右图当中红线上的数据。
-
这样,我们就建立起来了,投影数据和f(x,y)的傅里叶变换图像的关系,之后通过2D反傅里叶变换就可以得到f(x,y)的图像了。这就是重建。
关键在于,中心切片定理是如何证明的。 首先定义为投影的的傅里叶变换,为f(x,y)的傅里叶变换。
我们需要证明:
这里我们需要先知道狄拉克函数的两个性质:
我是这样理解的,狄拉克函数因为积分和为1,所以可以看成对f(x)求期望的过程。如果x变成了ax,倍数变化。那么意味着期望概率的稀释。所以稀释的倍数越大,最终要乘上稀释倍数的倒数。
这个也简单,我理解为是期望的平移,从f(0)移动到了f(a).
下面推导过程我手写了:
傅里叶变换与反傅里叶变换公式
- 1维傅里叶变换
- 1维傅里叶反变换
- 2维傅里叶反变换
FBP filtered backprojection 滤波反投影算法
其实通过上面的过程进行重建,也就是反投影,是会得到比较模糊的图像的。因为在构建的时候,是会出现密度不均匀的情况。因为每一个角度都是过原点的,所以越靠近原点的密度就越高,约远离原点的区域的密度越低。或者可以说,在傅里叶空间,低频区域的密度高,权重就会过大。因此为了消除这种模糊的效果,要对傅里叶空间进行加权矫正,使其密度均匀。
做法就是,在得到的投影傅里叶变换空间中,乘上。
现在我们来推导FBP算法,也就是反投影算法。
我们需要得到的f(x,y)就是通过反傅里叶变化的方法:
根据中心切片定理,原始图像的傅里叶变换等于投影的傅里叶变换,投影的因为密度不均匀需要通过进行矫正,这个w的绝对值其实就是。矫正后的投影的傅里叶变换我们写作,公式如下:
我们对w变量做1维的反傅里叶变换,
推荐阅读
-
构建下一代去中心化应用程序:基于 BASE 链的 DApp 开发
-
详细解释 Python 中的 __getitem__ 方法和切片对象
-
每周播客 | 1-18,阿里的测试太好还是开发者太弱?
-
[校园生活] 关于监考的思考 (1)
-
Webpack 源代码分析 (1) ----- webpack.cmd
-
goip技术原理图解_12式木桩模块化技术定型训练 1 小年
-
Java 实现 蓝桥杯 VIP 算法改进 分分钟碎纸
-
用 C 语言读取控制台的上下左右箭头键命令--方法 1
-
m1 MacBook Air 3 个月使用总结及原生运行于苹果芯片架构的软件推荐
-
2021 苹果电脑选购建议,资深用户带你了解 MacBook Air MacBook Pro(带 M1)