基于Operator Discretization Library (ODL)的CT图像正反投

持续更新…

Operator Discretization Library (ODL)是一个可以快速求解逆问题的python库。而计算机断层成像的解析或迭代重建方法本质上也是对目标函数进行求逆的过程，odl库中的tomo模块可以很好的解决2D/3D CT图像重建问题且可以利用GPU加速。有关odl库的详细文档描述见：

https://odl.readthedocs.io/

调用odl.tomo进行正反投操作会调用到astra工具箱，有关astra的详细文档描述见：

http://www.astra-toolbox.com/docs/

有关astra工具箱调用代码示例的ppt可见：

https://www.astra-toolbox.com/files/misc/ICTMS2019/20190722_ICTMS_ASTRA_workshop_2d.pdf

1. FBP（2D）

生成平行束投影的几何：
若重建空间recon_space为2D数据则parallel_beam_geometry调用Parallel2dGeometry，若recon_space为3D数据则parallel_beam_geometry调用Parallel3dAxisGeometry。

geometry = odl.tomo.parallel_beam_geometry(reco_space)

使用anaconda3安装odl，parallel_beam_geometry函数在odl库中的位置：

C:\ProgramData\Anaconda3\Lib\site-packages\odl\tomo\geometry\parallel.py

生成扇形束投影的几何：
若重建空间recon_space为2D数据则cone_beam_geometry调用FanFlatGeometry，若recon_space为3D数据则cone_beam_geometry调用ConeFlatGeometry。
若不设置num_angle的值则会使用默认的足够对数据进行完全采样的数值。

geometry = odl.tomo.cone_beam_geometry(reco_space,src_radius=40, det_radius=80,num_angles=None)

cone_beam_geometry函数在odl库中的位置：

C:\ProgramData\Anaconda3\Lib\site-packages\odl\tomo\geometry\conebeam.py

import odl
import matplotlib.pyplot as plt
import numpy as np
reco_space = odl.uniform_discr(min_pt=[-20, -20], max_pt=[20, 20], shape=[512, 512])
# parallel beam
geometry = odl.tomo.parallel_beam_geometry(reco_space)
# fan beam
#geometry = odl.tomo.cone_beam_geometry(reco_space,src_radius=40, det_radius=80,num_angles=None)
projector = odl.tomo.RayTransform(reco_space,geometry, impl='astra_cuda')
vol = odl.phantom.shepp_logan(projector.domain, modified=True)
proj = projector(vol)
backproj = projector.adjoint(proj) # 直接反投影
FBP = odl.tomo.fbp_op(projector) 
img = FBP(proj) # 滤波反投影
plt.imshow(np.flipud(proj.data.T),'gray')
plt.figure()
plt.imshow(np.flipud(img.data.T),'gray')
plt.figure()
plt.imshow(np.flipud(backproj.data.T),'gray')
plt.show()

图1 平行束：正投影（左），直接反投影（中），滤波反投影（右）

图2 扇形束：正投影（左），直接反投影（中），滤波反投影（右）