哪一个军事网站做的比较好,河南建筑公司实力排名,重庆发布微信公众号,96个html静态网站模板打包专栏地址#xff1a;『youcans 的 OpenCV 例程 300篇 - 总目录』 【第 7 章#xff1a;图像复原与重建】 110. 投影和雷登变换 111. 雷登变换反投影重建图像 112. 滤波反投影重建图像 【youcans 的 OpenCV 例程 300 篇】112. 滤波反投影重建图像 7. 投影重建图像
图像重建的…专栏地址『youcans 的 OpenCV 例程 300篇 - 总目录』 【第 7 章图像复原与重建】 110. 投影和雷登变换 111. 雷登变换反投影重建图像 112. 滤波反投影重建图像 【youcans 的 OpenCV 例程 300 篇】112. 滤波反投影重建图像 7. 投影重建图像
图像重建的基本思想就是通过探测物体的投影数据重建物体的实际内部构造。
7.3 滤波反投影重建图像 Filtered back projection
直接由正弦图得到反投影图像会存在严重的模糊这是早期 CT 系统所存在的问题。
傅立叶中心切片定理表明投影的一维傅立叶变换是得到投影区域的二维傅立叶变换的切片。滤波反投影重建算法在反投影前将每一个采集投影角度下的投影进行卷积处理从而改善点扩散函数引起的形状伪影有效地改善了重建的图像质量。
以 f(x,y)f(x,y)f(x,y) 表示需要重建的图像F(u,v)F(u,v)F(u,v) 的二维傅里叶反变换是 f(x,y)∫−∞∞∫−∞∞F(u,v)ej2π(uxvy)dudvf(x,y) \int^{\infty}_{-\infty} \int^{\infty}_{-\infty} F(u,v) e^{j 2 \pi (uxvy)} dudv f(x,y)∫−∞∞∫−∞∞F(u,v)ej2π(uxvy)dudv 利用傅里叶切片定理可以得到 f(x,y)∫0π[∫−∞∞∣ω∣G(ω,θ)ej2πωρdω]dθf(x,y) \int^{\pi}_{0} \big[ \int^{\infty}_{-\infty} |\omega|G(\omega,\theta) e^{j 2 \pi \omega \rho} d\omega \big] d\theta f(x,y)∫0π[∫−∞∞∣ω∣G(ω,θ)ej2πωρdω]dθ 括号 [] 内部是一个一维傅里叶反变换可以认为这是一个一维滤波器的传递函数。由于 ∣ω∣|\omega|∣ω∣ 是一个不可积的斜坡函数Slope function可以通过对斜坡加窗进行限制典型地如汉明窗Hamming window、韩窗Hann window。
该式也可以使用空间卷积来实现 f(x,y)∫0π[∫−∞∞∣ω∣G(ω,θ)ej2πωρdω]dθ∫0π[∫−∞∞g(ρ,θ)s(xcosθysinθ−ρ)dρ]dθf(x,y) \int^{\pi}_{0} \big[ \int^{\infty}_{-\infty} |\omega|G(\omega,\theta) e^{j 2 \pi \omega \rho} d\omega \big] d\theta \\ \int^{\pi}_{0} \big[ \int^{\infty}_{-\infty} g(\rho,\theta) s(xcos \theta ysin \theta - \rho) d\rho \big] d\theta \\ f(x,y)∫0π[∫−∞∞∣ω∣G(ω,θ)ej2πωρdω]dθ∫0π[∫−∞∞g(ρ,θ)s(xcosθysinθ−ρ)dρ]dθ 这表明将对应的投影 g(ρ,θ)g(\rho, \theta)g(ρ,θ) 与斜坡滤波器传递函数 s(ρ)s(\rho)s(ρ) 的傅里叶反变换进行卷积可以得到角度 θ\thetaθ 的各个反投影整个反投影图像可以通过对所有反投影图像积分得到。 例程 9.25滤波反投影重建图像
SL 滤波器是使用 Sinc 函数对斜坡滤波器进行截断产生其离散形式为 hSL(nδ)−2π2δ2(4∗n2−1),n0,±1,±2,...h_{SL}(n \delta) - \frac{2}{\pi ^2 \delta ^2 (4*n^2 -1)}, \ n0,\pm 1,\pm 2,... hSL(nδ)−π2δ2(4∗n2−1)2, n0,±1,±2,... 利用投影值和 SL 滤波器进行卷积然后再进行反投影就可以实现图像重建。 # 9.25: 滤波反投影重建图像from scipy import ndimagedef discreteRadonTransform(image, steps): # 离散雷登变换channels image.shape[0]resRadon np.zeros((channels, channels), dtypenp.float32)for s in range(steps):rotation ndimage.rotate(image, -s * 180/steps, reshapeFalse).astype(np.float32)resRadon[:, s] sum(rotation)return resRadondef inverseRadonTransform(image, steps): # 雷登变换反投影重建图像channels image.shape[0]backproject np.zeros((steps, channels, channels))for s in range(steps):expandDims np.expand_dims(image[:, s], axis0)repeat expandDims.repeat(channels, axis0)backproject[s] ndimage.rotate(repeat, s * 180/steps, reshapeFalse).astype(np.float32)invRadon np.sum(backproject, axis0)return invRadondef SLFilter(N, d): # SL 滤波器, Sinc 函数对斜坡滤波器进行截断filterSL np.zeros((N,))for i in range(N):filterSL[i] - 2 / (np.pi**2 * d**2 * (4 * (i-N/2)**2 - 1))return filterSLdef filterInvRadonTransform(image, steps): # 滤波反投影重建图像channels len(image[0])backproject np.zeros((steps, channels, channels)) # 反投影filter SLFilter(channels, 1) # SL 滤波器for s in range(steps):value image[:, s] # 投影值valueFiltered np.convolve(filter, value, same) # 投影值和 SL 滤波器进行卷积filterExpandDims np.expand_dims(valueFiltered, axis0)filterRepeat filterExpandDims.repeat(channels, axis0)backproject[s] ndimage.rotate(filterRepeat, s * 180/steps, reshapeFalse).astype(np.float32)filterInvRadon np.sum(backproject, axis0)return filterInvRadon# 读取原始图像img1 cv2.imread(../images/Fig0534a.tif, 0) # flags0 读取为灰度图像img2 cv2.imread(../images/Fig0534c.tif, 0)# 雷登变换imgRadon1 discreteRadonTransform(img1, img1.shape[0]) # Radon 变换imgRadon2 discreteRadonTransform(img2, img2.shape[0])# 雷登变换反投影重建图像imgInvRadon1 inverseRadonTransform(imgRadon1, imgRadon1.shape[0])imgInvRadon2 inverseRadonTransform(imgRadon2, imgRadon2.shape[0])# 滤波反投影重建图像imgFilterInvRadon1 filterInvRadonTransform(imgRadon1, imgRadon1.shape[0])imgFilterInvRadon2 filterInvRadonTransform(imgRadon2, imgRadon2.shape[0])plt.figure(figsize(9, 7))plt.subplot(231), plt.axis(off), plt.title(origin image), plt.imshow(img1, gray) # 绘制原始图像plt.subplot(232), plt.axis(off), plt.title(inverse Radon trans), plt.imshow(imgInvRadon1, gray)plt.subplot(233), plt.axis(off), plt.title(filter inv-Radon trans), plt.imshow(imgFilterInvRadon1, gray)plt.subplot(234), plt.axis(off), plt.title(origin image), plt.imshow(img2, gray)plt.subplot(235), plt.axis(off), plt.title(inverse Radon trans), plt.imshow(imgInvRadon2, gray)plt.subplot(236), plt.axis(off), plt.title(filter inv-Radon trans), plt.imshow(imgFilterInvRadon2, gray)plt.tight_layout()plt.show()本节完 版权声明 youcansxupt 原创作品转载必须标注原文链接(https://blog.csdn.net/youcans/article/details/123088438) Copyright 2022 youcans, XUPT Crated2022-2-28
