数学中的图像重构 -- CT中的 Radon 变换 图解
这是我之前在剑桥大学上的一节研究生应数选修课 Image Reconstruction,之前没怎么听懂,所以这段时间想把它补起来。这节课老师没有明确的讲义,所以我就照记着的一些书的顺序,把它复习了。整堂课只有我一个人上 QAQ,所以应该算是在数学系里比较小众的方向吧。因此这篇笔记 基本上是为了 我自己以后查资料或公式好找一点。笔记部分摘自 Mathematical Methods in ...
这是我之前在剑桥大学上的一节研究生应数选修课 Image Reconstruction,之前没怎么听懂,所以这段时间想把它补起来。
这节课老师没有明确的讲义,所以我就照记着的一些书的顺序,把它复习了。
整堂课只有我一个人上 QAQ,所以应该算是在数学系里比较小众的方向吧。
因此这篇笔记 基本上是为了 我自己以后查资料或公式好找一点。
笔记部分摘自 Mathematical Methods in Image Reconstruction. F.Natterer and F.Wuebbeling.
文章目录
引言
图像重构 Intro
Radon Transform 拉东变换 R f {\color{red}\mathcal{R}f} Rf
-
拉东变换是一个积分变换
-
先它将定义在二维平面上的一个函数 f ( x 1 , x 2 ) f(x_1,x_2) f(x1,x2) 沿着平面上的任意一条直线做线积分,相当于对函数 f ( x 1 , x 2 ) f(x_1,x_2) f(x1,x2) 做 CT扫描
-
令 X射线 在 x x x 点 线性衰减后 为 f ( x ) f(x) f(x)。 那么 X射线光束 L L L 经过的部分成像就是 g ( L ) = ∫ L f ( x ) d x g(L)=\int_L f(x)dx g(L)=∫Lf(x)dx 。
-
如果我们用 s s s (原点到射线的距离)与 θ \theta θ (距离的夹角)表示射线 L L L, 如图:
-
可得这条直线为 L ( θ , s ) = { x ∈ R 2 : x ⋅ θ = s } , θ ∈ S 1 , s ∈ R 1 L(\theta,s) = \{x\in \R^2: x\cdot \theta = s\},\; \theta \in S^1, \; s\in \R^1 L(θ,s)={x∈R2:x⋅θ=s},θ∈S1,s∈R1
-
可得如下方程:
g ( θ , s ) = ∫ x ⋅ θ = s f ( x ) d x = ( R f ) ( θ , s ) g(\theta,s) = \int_{x\cdot \theta = s} f(x)dx = (\mathcal{R}f)(\theta,s) g(θ,s)=∫x⋅θ=sf(x)dx=(Rf)(θ,s)
- 用python把它表示出来:
import matplotlib.pyplot as plt
import numpy as np
# 为了之后扩展到 n 维 这里用 x_1 x_2 表示
x1 = np.arange(-10,10.01,.01) # x_1 可以理解为 二维(x,y) 坐标轴的x
x2 = np.arange(-10,10.01,.01) # x_2 可以理解为 二维(x,y) 坐标轴的y
xv = np.meshgrid(x1,x2) # 网格化
degree = 45 # theta 角度
theta = degree/180 *np.pi
s = 2. # s 长度
Lv = abs(xv[0]*np.cos(theta) + xv[1]*np.sin(theta)-s) # x1*cos(theta) + x2*cos(theta) - s = 0
L = np.max((Lv<0.01)*xv[1],axis=0)+np.min((Lv<0.01)*xv[1],axis=0) # 筛选出 满足条件的值
# L.nonzero() L 的非零值坐标
plt.plot(x1[L.nonzero()],L[L.nonzero()],'blue')
plt.show()
- 我还做了变量为 s s s 和 θ \theta θ 的动图:
- 接下来,把它拓展到更高维度 n 维。
- 构筑 n 维超平面 H ( θ , s ) = { x ∈ R n : x ⋅ θ = s } , θ ∈ S n − 1 (unit sphere) , s ∈ R 1 H(\theta,s) = \{x\in \R^n: x\cdot \theta = s\},\; \theta \in S^{n-1} \text{ (unit sphere) }, \; s\in \R^1 H(θ,s)={x∈Rn:x⋅θ=s},θ∈Sn−1 (unit sphere) ,s∈R1
- 超平面垂直于 θ \theta θ 距离为 s s s
- 同时 H ( − θ , − s ) = H ( θ , s ) H(-\theta,-s) = H(\theta,s) H(−θ,−s)=H(θ,s) 为 偶函数
- 可定义
( R f ) ( θ , s ) = ∫ x ⋅ θ = s f ( x ) d x = ∫ H ( θ , s ) f ( x ) d x (\mathcal{R}f) (\theta, s)= \int_{x\cdot \theta = s} f(x)dx = \int_{H(\theta,s)}f(x) dx (Rf)(θ,s)=∫x⋅θ=sf(x)dx=∫H(θ,s)f(x)dx
- R f \mathcal{R}f Rf 函数在 单位圆柱 unit sylinder 上
C n = { ( θ , s ) : θ ∈ S n − 1 , s ∈ R 1 } {\color{red}C^n} = \{(\theta,s): \theta \in S^{n-1} ,\; s\in \R^1\} Cn={(θ,s):θ∈Sn−1,s∈R1}
-
因为 H ( − θ , − s ) = H ( θ , s ) H(-\theta,-s) = H(\theta,s) H(−θ,−s)=H(θ,s) , 所以 ( R f ) ( − θ , − s ) = R f ( θ , s ) (\mathcal{R}f)(-\theta,-s) = \mathcal{R}f(\theta,s) (Rf)(−θ,−s)=Rf(θ,s) 为偶函数
-
由于 s − x ⋅ θ = 0 s - x\cdot \theta = 0 s−x⋅θ=0 变换可得:
( R f ) ( θ , s ) = ∫ R n δ ( s − x ⋅ θ ) f ( x ) d x (\mathcal{R}f)(\theta,s) = \int_{\R^n} \delta(s - x \cdot \theta) f(x) dx (Rf)(θ,s)=∫Rnδ(s−x⋅θ)f(x)dx
- θ ⊥ = { x ∈ R n : x ⋅ θ = 0 } \theta^\perp = \{x\in\R^n: x\cdot \theta = 0\} θ⊥={x∈Rn:x⋅θ=0} 正交 θ \theta θ 的子空间
( R f ) ( θ , s ) = ∫ θ ⊥ f ( s θ + y ) d y (\mathcal{R}f)(\theta,s) = \int_{\theta^\perp}f(s\theta + y) dy (Rf)(θ,s)=∫θ⊥f(sθ+y)dy
- 由于 f ∈ S ( R n ) f \in {\color{red}\mathcal{S}(\R^n)} f∈S(Rn) 在速降函数空间内(Schwartz space)
- 那么 R f ∈ S ( C n ) \mathcal{R}f \in \mathcal{S}(C^n) Rf∈S(Cn)
- S ( C n ) {\color{red}\mathcal{S}(C^n)} S(Cn) 为
S ( C n ) = { g ∈ C ∞ ( C n ) : s l ∂ k ∂ s k g ( θ , s ) bounded , l , k = 0 , 1 , ⋯ } \mathcal{S}(C^n) = \Big \{ g \in \mathcal{C}^\infty(C^n): s^l \frac{\partial^k}{\partial s^k} g(\theta,s) \text{ bounded }, \; l,k = 0,1,\cdots \Big \} S(Cn)={g∈C∞(Cn):sl∂sk∂kg(θ,s) bounded ,l,k=0,1,⋯}
-
C ∞ ( C n ) {\color{red}\mathcal{C}^\infty}(C^n) C∞(Cn) 是所有从 C n C^n Cn 射到 C \mathcal{C} C 的光滑函数 (无穷阶可导的函数)
-
卷积运算
( g ⋆ h ) ( θ , s ) = ∫ R 1 g ( θ , s − t ) h ( θ , t ) d t (g \star h)(\theta, s) = \int_{\R^1} g(\theta, s-t)h(\theta, t)dt (g⋆h)(θ,s)=∫R1g(θ,s−t)h(θ,t)dt
傅里叶变换 F {\color{red}\mathcal{F}} F
F ( g ( θ , s ) ) = ( 2 π ) − 1 / 2 ∫ R 1 g ( θ , s ) e − i s σ d s \mathcal{F}(g(\theta, s)) = (2\pi)^{-1/2} \int_{\R^1}g(\theta, s) e^{-is\sigma} ds F(g(θ,s))=(2π)−1/2∫R1g(θ,s)e−isσds
- 当 f ∈ S ( R n ) , θ ∈ S n − 1 (unit sphere) , σ ∈ R 1 f\in \mathcal{S}(\R^n), \; \theta \in S^{n-1} \text{ (unit sphere) }, \; \sigma \in \R^1 f∈S(Rn),θ∈Sn−1 (unit sphere) ,σ∈R1
F [ ( R f ) ( θ , σ ) ] C n = ( 2 π ) n − 1 2 F [ f ( σ θ ) ] R n \mathcal{F}\big[(\mathcal{R}f)(\theta, \sigma)\big]_{C^n} = (2\pi)^{\frac{n-1}{2}}\mathcal{F}[f(\sigma \theta)]_{\R^n} F[(Rf)(θ,σ)]Cn=(2π)2n−1F[f(σθ)]Rn
- 当 f , g ∈ S ( R n ) f, g \in \mathcal{S}(\R^n) f,g∈S(Rn)
R f ⋆ R g C n = R ( f ⋆ g ) R n \underset{C^n}{\mathcal{R}f \star \mathcal{R}g} = \underset{\R^n}{\mathcal{R}(f\star g)} CnRf⋆Rg=RnR(f⋆g)
BackProjection operator 后向投影 算子 R ∗ {\color{red}\mathcal{R}^*} R∗
( R ∗ g ) ( x ) = ∫ S n − 1 g ( θ , x ⋅ θ ) d θ , g ∈ S ( C n ) (\mathcal{R}^* g)(x) = \int_{S^{n-1}} g(\theta, x\cdot \theta) d\theta, \; g\in \mathcal{S}(C^n) (R∗g)(x)=∫Sn−1g(θ,x⋅θ)dθ,g∈S(Cn)
-
因此 对于 g = R f , ( R ∗ g ) ( x ) g = \mathcal{R}f, \; (\mathcal{R}^*g)(x) g=Rf,(R∗g)(x) 是 所有超平面 经过 x x x, f f f 的积分 的均值。
-
在数学中 可以说 R ∗ \mathcal{R}^* R∗ 是 R \mathcal{R} R 的共轭
-
对于 g ∈ S ( R 1 ) , f ∈ S ( R n ) g \in \mathcal{S}(\R^1), \; f\in \mathcal{S}(\R^n) g∈S(R1),f∈S(Rn)
∫ R 1 g ( s ) ( R f ) ( θ , s ) d s = ∫ R n g ( x ⋅ θ ) f ( x ) d x \int_{\R^1} g(s) (\mathcal{R}f)(\theta,s) ds = \int_{\R^n} g(x \cdot \theta)f(x) dx ∫R1g(s)(Rf)(θ,s)ds=∫Rng(x⋅θ)f(x)dx
- 对于 g ∈ S ( C n ) , f ∈ S ( R n ) g \in \mathcal{S}(C^n), \; f\in \mathcal{S}(\R^n) g∈S(Cn),f∈S(Rn)
∫ S n − 1 ∫ R 1 g R f d θ d s = ∫ R n ( R ∗ g ) f ( x ) d ( x ) d x \int_{S^{n-1}}\int_{\R^1} g \mathcal{R}f d\theta ds = \int_{\R^n} (\mathcal{R}^*g)f(x)d(x) dx ∫Sn−1∫R1gRfdθds=∫Rn(R∗g)f(x)d(x)dx
- 当 f ∈ S ( R n ) f\in \mathcal{S}(\R^n) f∈S(Rn) 和 g ∈ S ( C n ) g\in \mathcal{S}(C^n) g∈S(Cn)
( R ∗ g ) ⋆ f = R ∗ ( g ⋆ R f ) (\mathcal{R}^*g ) \star f = \mathcal{R}^* (g \star \mathcal{R}f) (R∗g)⋆f=R∗(g⋆Rf)
- 当 g ∈ S ( C n ) g\in\mathcal{S}(C^n) g∈S(Cn) 为偶 (例 g ( θ , s ) = g ( − θ , − s ) g(\theta,s)= g(-\theta,-s) g(θ,s)=g(−θ,−s) )
F [ ( R ∗ g ) ( ξ ) ] = 2 ( 2 π ) n − 1 2 ∣ ξ ∣ 1 − n F [ g ( ξ ∣ ξ ∣ , ∣ ξ ∣ ) ] \mathcal{F}\big[(\mathcal{R}^* g)(\xi)\big] = 2 (2\pi)^{\frac{n-1}{2}}\lvert \xi \rvert^{1-n} \mathcal{F}\Big[g(\frac{\xi}{\lvert \xi \rvert},\lvert \xi \rvert)\Big] F[(R∗g)(ξ)]=2(2π)2n−1∣ξ∣1−nF[g(∣ξ∣ξ,∣ξ∣)]
- 通过 Riesz potential 里斯位势 在 C n C^n Cn
F [ ( I α g ) ( θ , σ ) ] = ∣ σ ∣ − α F [ g ( θ , σ ) ] , α < 1 \mathcal{F}\big[(\mathcal{I}^\alpha g)(\theta, \sigma)\big] = \lvert \sigma \rvert^{-\alpha} \mathcal{F}\Big[g(\theta,\sigma)\Big], \; \alpha <1 F[(Iαg)(θ,σ)]=∣σ∣−αF[g(θ,σ)],α<1
更多推荐
所有评论(0)