图像处理第7章 图像投影重建

第七章图像投影重建图像投影重建，一般只从一个物体的多个（轴向）投影图重建目标图像的过程。输入是投影图，输出是重建图。投影重建与计算机层析成像（CT）关系密切。一个X射线源与接受胶片以相反方向运动。章节安排：7.1节介绍典型的投影重建方式，CT、MRI、电阻抗断层7.2节2-D、3-D投影原理分析讨论，给出由拉东变换得到的中心层定理7.3节介绍重建图像质量较差，但所需要计算量较小的傅里叶反变换重建法7.4节介绍几种容易用软件和硬件实现，切效果比较准确清晰的逆投影重建法7.5节介绍代数重建技术，可以通过迭代计算直接得到数值解7.6节讨论将变换法和级数展开法相结合的综合方法图象投影重建的历史1971年：第1台CT机器建成1972年：公布这台机器（1973年正式使用）1979年：两个发明人豪斯费尔德（G.H.Hounsfield）和柯马克（A.M.Cormack）获得了诺贝尔生理和医学奖，获得CT图象也被认为是第一次通过解决一个属于逆问题和病态问题（inverseandill-posedproblems）的数学问题来获得图象的成功实例图象投影重建根据一个物体的投影图重建目标图象的过程输入投影图

（图象处理）

输出重建图

一类特殊的图象恢复技术投影：退化过程（丢失了沿射线方向的分辨能力）重建：复原过程（恢复了2-D空间的分辨能力）7.1投影重建方式{如果传感器测量的数据具有物体某种感兴趣物理特性在空间分布的积分形式，就可以用投影重建的方法来获取物体内部、反应不同物理特性的图像}7.1.1透射断层成像TCT（简称CT），英文为TransmissionComputedTomography，CAT（Computer-AidedTomography，CT）发射源射出的射线穿透物体到达接收器，射线被物体吸收一部分，余下部分被接收器接收。接收器获得的射线强度实际上反映了物体各部分对射线的吸收情况若物体是均匀的，则：I：穿透物体的射线强度I0：射线源的强度K(s)沿射线方向物体点s的线性衰减系数/因子，L代表辐射的射线物体线性衰减系数7.1.1透射断层成像图（a）—（d）分别对应第1—4代CT系统的扫描成像结构示意图。圆代表拟成像的区域，经过发射源（X射线管）的虚线直线箭头表示发射源可沿箭头移动，从一个发射源到另一个发射源的虚线曲线箭头表示发射源可沿曲线转动第1代：发射源与接收器一对一，对向移动以覆盖整个拟成像区域第2代：每个发射源对应若干接收器，对向移动以覆盖整个拟成像区域第3代：每个发射源指对应一个接收器，发射源不移动，只需转动第4代：接收器构成完整的圆环，工作中没有运动，只需发射源转动第3、4代系统采用扇束投影方式，可尽量缩短投影时间，减少危害7.1.2发射断层成像ECT（emissioncomputedtomography）发射源在物体内部，接收器在物体外部（通常将放射性离子注入物体内部，从物体外接收其辐射），了解离子在物体内运动情况和分布，从而检测到与生理有关的状况/信息主要有两种：①PET，正电子发射（positronemissiontomography）正电子与负电子相撞湮灭而产生一对光子②SPECT，单光子发射（singlephotonemissionCT）使用在衰减中能产生γ射线的放射性离子。7.1.2发射断层成像①PET（positronemissiontomography）相对放置的两个接收器可以确定一条射线衰减时释放正电子的放射性粒子，放出的正电子很快与负电子相撞湮灭产生一对光子并以相反方向射出，所以相对放置的两个接收器可以确定一条射线P是投影数据7.1.2发射断层成像②SPECT（singlephotonemissionCT）将放射性物质注入物体内，不同的材料（如组织或器官）吸收后会发射γ射线，为确定射线方向，要用能阻止射线偏移的准直器来定向采集光子。一定方向的γ射线穿过准直器到达晶体，在那里γ射线光子转化为电信号。这些电信号提供了光子与晶体作用的位置，放射性物质的3-D分布就转化为2-D投影图像。7.1.3反射断层成像利用投影重建的原理工作。雷达系统获取的雷达图是由物体反射的回波所产生的雷达接收器在特定角度所接收到的回波强度是地面反射量在一个扫描段的积分（对比CT）投影重建就是要从这个积分获得地面（反射强度）的图象7.1.3反射断层成像(RCT)非聚焦合成孔径雷达(SAR)雷达运动，目标静止，相对运动增加分辨率目标A处回波信号的双程超前相位：

目标B处回波信号双程超前相位B处回波响应R为目标A与雷达间最近距离，d为B与A间的位移量，v是雷达沿Y轴速度，T是有效积累时间7.1.4电阻抗断层成像(EIT)（EIT）采用交流电场对物体进行激励。对电导或电抗的改变比较敏感。将低频率的电流注入物体内部并测量在物体外表处的电势场；采用图象重建算法重建出物体内部区域的电导和电抗的分布或变化的图象(基于边界测量值估计场域电导率分布)EIT图象能反映组织或器官携带的病理和生理信息：不同的生物组织或器官在不同的生理、病理条件下其电阻抗特性不同。（图7.1.5）EIT技术分析数学角度上，EIT和各类CT有类似之处，都需要处理外部信息来获得反应内部结构的图像，成像常针对穿过物体的一个2-D截面进行。EIT借助电流的扩散来获得电导的分布，与各CT不同；成像技术安全、简便，分辨率较差，分辨率依赖于电极数量，可以同时接触到物体的电极数量常受到很多限制。缺点是由于非线性的病态问题，如果测量又很小的误差就有可能对电导的计算产生很大影响。7.1.5磁共振成像（MRI）早期称核磁共振。质子在磁场中会进动.当一定强度和频率的共振场信号作用于物体时，质子吸收能量并转向与磁场相交的朝向质子吸收的能量释放并被接收器检测到。根据检测到的信号就可以确定质子的密度检测到的信号是MRI信号沿直线的积分。所以检测目标的工作成了由投影重建的问题每个MRI系统都包括磁场子系统、发射/接收子系统、计算机图像重建和显示子系统。磁共振成像根据时变非均匀磁场和射频磁场及其激励而产生的磁共振信号来重建物体的自旋密度分布函数磁共振成像（MRI）7.2投影重建原理{各种成象方式在原理不完全相同但有一些共性}7.2.1基本模型将需要投影重建的物质材料限制在一个无限薄的平面上，使得重建图像在任意点的灰度值正比于射线投影到的那个点所固有的相对线性衰减系数，（f(x,y)代表某种物理量在2-D平面上的分布）7.2.1投影重建图像示意图设f(x,y)在以原点为圆心的单位圆Q外为0，现考虑有一条由发射源到接收器的直线在平面上与f(x,y)在Q内相交，这条直线用两个参数来确定：1，它与原点的距离s；2，它与Y轴的夹角θ。7.2.1基本模型设Q为单位圆，积分上下限分别为t和-t沿直线（s,θ）对f（x,y)的积分(7.2.1)7.2.2拉东变换拉东变换主要用来接方程式（7.2.1）线积分直线l

的方程借助δ函数P和θ定义的图中直线l的线积分（7.2.4）7.2.2拉东变换中心层定理可以证明，对f(x,y)的2-D傅里叶变换与对f(x,y)先进行拉东变换后再进行1-D傅里叶变换得到的结果相等（7.2.5）7.2.2拉东变换对f(x,y)沿一个固定角度投影结果的1-D傅里叶变换对应f(x,y)的2-D傅里叶变换中沿相同角度的一个剖面/层，如图7.2.3.7.3傅里叶反变换重建基于变换的重建方法，它是首先在投影中得到应用的方法1.基本步骤和定义(1)建立数学模型，其中已知量和未知量都是连续实数的函数(2)利用反变换公式（可有多个等价的）解未知量(3)调节反变换公式以适应离散、有噪声应用的需求重建算法：设图象区被1个直角网格所覆盖，K为X方向上的点数，L为Y方向上的点数，要通过M×N个测量值g(m△s,n△θ)估计出在K×L个采样点的f(k△x,l△y)7.3傅里叶反变换重建考虑在s和θ上都均匀采样的情况M个间距为△

s的射线，N个相差△θ

的角度进行投影保证一系列射线覆盖单位圆，选取和

此时g(m△s,n△θ)为平行投影的射线数据。7.3傅里叶反变换重建2傅里叶变换投影定理G(R，θ)是

g（s,θ）对应第一个变量s的1-D傅里叶变换F(X,Y)是f(x,y)的2-D傅里叶变换f(x,y)以θ角进行投影的傅里叶变换等于f(x,y)的傅里叶变换在傅里叶空间(R,θ)处的值

f(x,y)在与X轴成θ

角的直线上投影的傅里叶变换是f(x,y)的傅里叶变换在朝向角θ上的1个截面。7.3傅里叶反变换重建3.模型重建

为检验重建公式的正确性和把握重建算法中各个参数对重建效果的影响，利用幻影模型进行试验7.3傅里叶反变换重建图7.3.2给出的一幅改进的Sheep-Logan头部模型图(尺寸115×115,256级灰度)。10个圆/椭圆(见表7.3.1数据)原Sheep-Logan的对比度小，图中小的椭圆看不清楚，改进图调整了各个椭圆的密度，使得各个椭圆的对比度有所增强，视觉效果更好一些。7.4逆投影重建7.4.1逆投影重建原理将从各个方向得到的投影逆向返回到该方向的各个位置，如果对多个投影方向都进行这样的逆投影并叠加结果，就有可能建立平面上的一个分布。（a）分别给出水平投影和逆投影的示意图，发射源发出均匀射线，由于所穿透物体各处密度不同，各接收器得到的响应不同。（b）给出垂直投影和逆投影的示意图，与水平方向的效果类似7.4.1逆投影重建原理利用逆投影可以重建原始空间密度分布不同方向逆投影结果叠加逆投影重建的效果（a）表示物体密度分布；（b）、（c）分别表示水平和垂直逆投影结果；（d）表示将水平和垂直逆投影结果叠加的效果；(e)图表示将更多逆投影结果叠加的效果，随着逆投影结果不断叠加，将越来越反映原始物体的密度分布。7.4.2

卷积逆投影重建1.连续公式推导投影定理将每个投影得到的数据像图像采集那样扩散回图像，并不需要像傅里叶方法那样存储复频率空间图，在极坐标系中取傅里叶反变换：只在有限带宽|R|<1/(2△s)的情况下对G(R,θ)进行估计（7.4.1）（7.4.2）7.4.2卷积逆投影重建将式(7.4.2)带入(7.4.1)并交换对s和R的积分次序

（7.4.5）g'(s',θ)是f(x,y)在θ角方向的投影与h(s)的卷积，可称为在θ方向上卷积了的投影；h(s)称为卷积函数。因为g'(·)中的参数是一条以θ角通过(x,y)点的射线的参数，所以fw(x,y)是所有与过(x,y)的射线所对应的卷积后投影的积分。（7.4.3）7.4.2卷积逆投影重建2.离散计算对（7.4.5）代表的逆投影过程近似对每个θn需要对K×L个s'计算g'(s',θn').由于K和L一般都很大，直接计算工作量很大，所以卷积过程的近似（分两步以减少计算量）

离散卷积

一次插值(7.4.6)7.4.2

卷积逆投影重建3.扇束投影重建实际应用中需要尽量缩短投影时间，减少由于物体在投影期间的运动而造成的图像失真以及对患者的伤害。扇束投影主要有两种集合测量类型，分别对应第二代和第三代CT系统，使用1个发射器和1组接收器

要在扇形投影的情况下进行重建，可通过将中心投影转化为平行投影，在用平行投影重建技术进行重建。讨论接收器在一段弧上等角度间隔排列的情况，用（s,θ）所指定的一条射线可看做是一组用（α，β）指定的射线中的一条，其中α是该射线与中心射线的离散角，β是源与原点连线和Y轴夹角，它确定了源的方向。线积分g(s,θ)现在记为p(α，β)(对|s|<D,D是源到原点距离)。假设源处于物体外部，所以对所有的|α|<δ,有p(α，β)=0。δ为与目标相切的射线与中心射线的夹角。s=Dsinαθ=α+βg(s,θ)=p(α，β)(7.4.9)设U是从源到需要重建点P的距离（P的位置用r和ϕ表示），ψ是从源到P的连线与中心发射线的夹角，则得到（7.4.10）（7.4.11）（7.4.12）综合以上各式（7.4.3）改写为：（7.4.13）将(s,θ)换成(α,β):将式(7.4.2)带入，并把对β的积分上下限换为2π和0，则上式可分解转化成两步：（7.4.15）（7.4.16）（7.4.14）7.4.2卷积逆投影重建傅里叶反变换重建法和卷积逆投影法的比较傅里叶反变换重建法卷积逆投影法2-D傅里叶反变换用直角坐标表示2-D傅里叶反变换用极坐标表示所需计算量比较小，较适合图像尺寸或数据量较大基本算法容易用软件和硬件实现，可重建出准确清晰的图像平行投影时导出的不能在保持原有效率的条件下进行修改以适用于扇形扫描投影平行投影时可用不同的方法修改以适用于扇形扫描投影7.4.3其他逆投影重建方法1.逆投影滤波逆投影滤波器（filteroftheback-projections）先逆投影，再滤波（卷积）对投影进行逆投影操作得到的结果是模糊了的图像，它是真实图像与1/(x2+y2)1/2卷积的结果.令对投影结果进行逆投影操作得到的模糊图像为（7.4.17）真实图像与模糊图像由下式联系（7.4.18）2-D卷积具体推导过程见书本165页实现逆投影滤波在空域和在频域的流程分别如图7.4.7左右部分所示拉东变换傅里叶反变换傅里叶变换逆投影7.4.3其他逆投影重建方法2.滤波逆投影（back-projectionofthefilteredprojections）先对投影滤波，再逆投影

基本思路：每个投影中的衰减（由目标吸收造成）与沿每个投影的目标结构有关，如果仅从一个投影是无法获得沿该方向各个位置的吸收的，但可以把对吸收值的测量平均分布到该方向上，要是能对多个方向中的每个都进行这样的分派，则将吸收值叠加起来就可得到反映目标结构的特征值。这种重建方法与收集足够多的投影并进行傅里叶空间重建是等价的，但计算量要少得多。滤波逆投影方法可看做拉东反变换的一种近似的计算机实现方法。它利用傅里叶反变换将对q函数的计算转化为对径向量s的计算，在空域和频域的实现流程分别如图所示。拉东变换傅里叶变换傅里叶反变换逆投影7.5代数重建技术1.重建模型

代数重建技术也称迭代算法或优化技术，它从一开始就在离散域中建模，且求解是迭代进行的。y：测量矢量x：图像矢量1.重建模型将要重建的目标放在一个直角坐标网格中，发射源和接收器都是点状的，他们之间连线对应一条射线，将每个像素按扫描次序排为1到N（网格总数）第i条射线与第j个象素相交段的长度aij，代表对第j个象素沿第i条射线的贡献的权值.（如图）重建的输入数据是沿每条射线的投影积分，在沿射线方向上每个像素位置对线性衰减系数贡献的求和值就等于测量到的吸收数值，即投影结果。在投影已知的情况下，对每条射线的积分都能提供一个方程，合起来构成一组齐次方程，未知数的数量就是图像平面中像素的数量，方程的数量就是线积分的数量。重建看做是解一组齐次方程组。如果用yi表示沿射线方向上的总吸收测量值，则（7.5.1）写成矩阵形式：y=Ax(7.5.2)y：测量矢量x：图像矢量A是M×N投影矩阵，为获得高质量的图像，M和N都需要在105量级，，所以A是一个非常大的矩阵，但由于对每条射线来说，它只与很少的像素相交，因此A中长只有少于1%的元素不为0.7.5代数重建技术2.无松弛的代数重建技术首先初始化一个图像矢量x作为迭代起点，然后如下迭代：

思路：每次取1条射线，改变图象中与该线相交象素的灰度值，从而把当前的图象矢量x

(k)更新为x

(k+1)。具体运算中就是将测量值与由当前算得的投影数据的差正比于aij重新分配到（该射线的）各个象素上去。（7.5.3）7.5代数重建技术3.松弛的代数重建技术将无松弛的代数重建技术推广，迭代表达式：

r是松弛系数，它控制了收敛的速度r的取值一般在0~2之间，当r=1时，式(7.5.4)变成式(7.5.3)；当r很小时，式(7.5.4)变成与传统的最小平方解等价。（7.5.4）7.5代数重建技术4.级数法的一些特点前面介绍的代数重建技术都利用了级数展开的方法，一般变化法比级数法所需要的计算量小得多，因而大多数实用系统都采用变换法。但级数法有一些独特的优点：(1)由于在空域中比较灵活，所以容易调整以适应新的应用(2)能重建出较高对比度的图象（对密度突变的材料）(3)借助多次迭代可从较少数量的投影（<10）重建图象(4)比变换法更适合于ECT(5)比变换法更适合于3-D重建问题(6)比变换法更适合于不完整投影情况7.6综合重建方法“综合”

有时体现在公式的推导里，有时体现在实现的方法上，还有时体现在实际的应用中综合重建法一些例子：

(1)迭代变换法（iterativetransform）(2)迭代重建重投影（iterativer