发明名称 基于通用图形显示卡的被测体正投影与反投影方法
摘要 本发明属于投影探测与成像技术领域,其特点在于该方法含有以下步骤:把被测物体的重建体数据用三维纹理存储表示;预先计算射线出射点至平面探测器上个探测器之间的方向向量数组,并存储到通用图形显示卡中形成方向向量纹理;利用方向向量纹理与方向向量绕旋转轴旋转来实现任一投影角度的正投影过程的一次完成,无中间结果转存;在反投影过程中,利用三维纹理索引技术使虚拟切片组的放置位置不再受到限制。本发明提高了正投影过程的速度与精度,也增加了反投影过程的灵活性。
申请公布号 CN100386779C 申请公布日期 2008.05.07
申请号 CN200610012090.1 申请日期 2006.06.02
申请人 清华大学 发明人 陆文凯;张雷
分类号 G06T11/00(2006.01) 主分类号 G06T11/00(2006.01)
代理机构 代理人
主权项 1.基于通用图形显示卡的被测体正投影与反投影方法,其特征在于:所述方法分为正投影与反投影两个阶段,其中:所述的正投影方法依次含有以下步骤:步骤(1),把所述被测体置于一个立方体中,并将此立方体在计算机中存储为三维数组,称为重建数据体,在以该立方体的几何中心为原点,以垂直于该立方体表面方向建立x、y、z三维直角坐标系;步骤(2),以y轴为旋转轴,以位于所述z轴延长线上一点A(0,0,zA 0)为x射线出射点,以呈面阵排列的x射线平板探测器为接收平面,连接x射线出射点A与所述平板探测器上的一个探测器D(xD 0,yD 0,zD 0),得到该探测器在投影角度为0度时的方向向量(sD,tD,pD):<math><mrow><mfenced open='[' close=']'><mtable><mtr><mtd><msub><mi>s</mi><mi>D</mi></msub></mtd></mtr><mtr><mtd><msub><mi>t</mi><mi>D</mi></msub></mtd></mtr><mtr><mtd><msub><mi>p</mi><mi>D</mi></msub></mtd></mtr></mtable></mfenced><mo>=</mo><mfenced open='[' close=']'><mtable><mtr><mtd><msubsup><mi>x</mi><mi>D</mi><mn>0</mn></msubsup><mo>/</mo><mi>M</mi></mtd></mtr><mtr><mtd><msubsup><mi>y</mi><mi>D</mi><mn>0</mn></msubsup><mo>/</mo><mi>M</mi></mtd></mtr><mtr><mtd><mrow><mo>(</mo><msubsup><mi>z</mi><mi>D</mi><mn>0</mn></msubsup><mo>-</mo><msubsup><mi>z</mi><mi>A</mi><mn>0</mn></msubsup><mo>)</mo></mrow><mo>/</mo><mi>M</mi></mtd></mtr></mtable></mfenced></mrow></math> 其中<math><mrow><mi>M</mi><mo>=</mo><msup><mrow><mo>(</mo><msup><mrow><mo>(</mo><msubsup><mi>x</mi><mi>D</mi><mn>0</mn></msubsup><mo>)</mo></mrow><mn>2</mn></msup><mo>+</mo><msup><mrow><mo>(</mo><msubsup><mi>y</mi><mi>D</mi><mn>0</mn></msubsup><mo>)</mo></mrow><mn>2</mn></msup><mo>+</mo><msup><mrow><mo>(</mo><msubsup><mi>z</mi><mi>D</mi><mn>0</mn></msubsup><mo>-</mo><msubsup><mi>z</mi><mi>A</mi><mn>0</mn></msubsup><mo>)</mo></mrow><mn>2</mn></msup><mo>)</mo></mrow><mfrac><mn>1</mn><mn>2</mn></mfrac></msup></mrow></math> 对于所述平板探测器内的所有探测器,则得到一个与该探测器大小相同的方向向量数组;步骤(3),利用所述通用图形显示卡将步骤(2)得到的方向向量数组存储为一个二维纹理,其中每个方向向量的三个分量sD,tD,pD分别存储在纹理的R、G、B三个颜色通道中,形成二维方向向量纹理;步骤(4),在所述x射线出射点与平板探测器相对位置保持不变的条件下,围绕y轴旋转θ度,该计算机计算各个探测器的方向向量(s,t,p):<math><mrow><mfenced open='[' close=']'><mtable><mtr><mtd><mi>s</mi></mtd></mtr><mtr><mtd><mi>t</mi></mtd></mtr><mtr><mtd><mi>p</mi></mtd></mtr></mtable></mfenced><mo>=</mo><mfenced open='[' close=']'><mtable><mtr><mtd><mi>cos</mi><mi>&theta;</mi></mtd><mtd><mn>0</mn></mtd><mtd><mi>sin</mi><mi>&theta;</mi></mtd></mtr><mtr><mtd><mn>0</mn></mtd><mtd><mn>1</mn></mtd><mtd><mn>0</mn></mtd></mtr><mtr><mtd><mo>-</mo><mi>sin</mi><mi>&theta;</mi></mtd><mtd><mn>0</mn></mtd><mtd><mi>cos</mi><mi>&theta;</mi></mtd></mtr></mtable></mfenced><mfenced open='[' close=']'><mtable><mtr><mtd><msub><mi>s</mi><mi>D</mi></msub></mtd></mtr><mtr><mtd><msub><mi>t</mi><mi>D</mi></msub></mtd></mtr><mtr><mtd><msub><mi>p</mi><mi>D</mi></msub></mtd></mtr></mtable></mfenced></mrow></math> 此时射线出射点A的坐标为(xA θ,yA θ,zA θ)为<math><mrow><mfenced open='[' close=']'><mtable><mtr><mtd><msubsup><mi>x</mi><mi>A</mi><mi>&theta;</mi></msubsup></mtd></mtr><mtr><mtd><msubsup><mi>y</mi><mi>A</mi><mi>&theta;</mi></msubsup></mtd></mtr><mtr><mtd><msubsup><mi>z</mi><mi>A</mi><mi>&theta;</mi></msubsup></mtd></mtr></mtable></mfenced><mo>=</mo><mfenced open='[' close=']'><mtable><mtr><mtd><mi>cos</mi><mi>&theta;</mi></mtd><mtd><mn>0</mn></mtd><mtd><mi>sin</mi><mi>&theta;</mi></mtd></mtr><mtr><mtd><mn>0</mn></mtd><mtd><mn>1</mn></mtd><mtd><mn>0</mn></mtd></mtr><mtr><mtd><mo>-</mo><mi>sin</mi><mi>&theta;</mi></mtd><mtd><mn>0</mn></mtd><mtd><mi>cos</mi><mi>&theta;</mi></mtd></mtr></mtable></mfenced><mfenced open='[' close=']'><mtable><mtr><mtd><mn>0</mn></mtd></mtr><mtr><mtd><mn>0</mn></mtd></mtr><mtr><mtd><mo>-</mo><mi>OE</mi></mtd></mtr></mtable></mfenced></mrow></math> 其中OE为射线出射点A到旋转中心即坐标原点的距离;步骤(5),该计算机在X射线出射点A到探测器平面之间按设定的对被测重建数据体的采样密度划分采样切片,再计算X射线出射点到采样点的距离即采样深度为R1时的纹理采样坐标(x,y,z):<math><mrow><mfenced open='[' close=']'><mtable><mtr><mtd><mi>x</mi></mtd></mtr><mtr><mtd><mi>y</mi></mtd></mtr><mtr><mtd><mi>z</mi></mtd></mtr></mtable></mfenced><mo>=</mo><msub><mi>R</mi><mn>1</mn></msub><mfenced open='[' close=']'><mtable><mtr><mtd><mi>s</mi></mtd></mtr><mtr><mtd><mi>t</mi></mtd></mtr><mtr><mtd><mi>p</mi></mtd></mtr></mtable></mfenced><mo>+</mo><mfenced open='[' close=']'><mtable><mtr><mtd><msubsup><mi>x</mi><mi>A</mi><mi>&theta;</mi></msubsup></mtd></mtr><mtr><mtd><msubsup><mi>y</mi><mi>A</mi><mi>&theta;</mi></msubsup></mtd></mtr><mtr><mtd><msubsup><mi>z</mi><mi>A</mi><mi>&theta;</mi></msubsup></mtd></mtr></mtable></mfenced></mrow></math> 同理得到各个探测器在不同采样深度下的纹理采样坐标;步骤(6),将存储内容为该重建数据体的三维数组及步骤(5)计算得到的同一探测器不同采样深度的纹理采样坐标分别代入到OpenGL图形开发语言提供的三维纹理索引函数中,得到同一探测器对应的不同采样深度的纹理值,将这些纹理值相加就得到该探测器的投影值,同理,该计算机可以计算得到在此投影角度θ下其他探测器的投影值;所述的反投影方法依次含有以下步骤:步骤(I),构筑一个与正投影过程相同的立方体,在以该立方体的几何中心为原点,以垂直于该立方体表面方向建立x、y、z三维直角坐标系,把该立方体在计算机中存储为三维数组,称为重建数据体,并将此重建数据体分割为多个二维切片,以下把这些二维切片称为虚拟切片;步骤(II),当虚拟切片位置与原重建数据体正切片位置相同时,按以下步骤进行,所述重建数据体正切片是指所述立方体的一个二维切片,且该切片垂直于直角坐标系主轴即旋转轴y轴,而切片的各边又分别平行于其余两个坐标系主轴;步骤(II.a),对虚拟切片上的任一体素V(x,y,z),计算其对应的投影纹理坐标Tpro(sT,tT,pT,qT):Tpro(sT,tT,pT,qT)=V(x,y,z,1)·TM其中TM为纹理矩阵TM=Tans1·Rot·Tans2·Proj·Sc·Tans3 Tans1、Rot、Tans2、Proj、Sc、Tans3为OpenGL图形开发语言提供的基本操作矩阵,其中Tans1代表平移DAO个单位距离,DAO为旋转轨道半径,Rot代表旋转θ角,θ为投影角度,Tans2代表平移-DAO个单位距离,Proj代表透视关系矩阵,由OpenGL图形开发语言根据射线锥角决定,Sc代表放大0.5倍,Tans3代表平移0.5个单位距离;步骤(II.b),将二维数组和步骤(II.a)所述的投影纹理坐标代入到OpenGL图形开发语言提供的投影纹理索引函数,即输出所述体素V(x,y,z)的更新值dv,同理,其他体素的更新值也可得到;步骤(IIi),当虚拟切片位置与原重建数据体正切片位置不一致时,按以下步骤进行:步骤(III.a),按照步骤(II.a)和步骤(II.b)的方法计算此时虚拟切片上各个体素的更新值;步骤(III.b),将步骤(III.a)得到的所有虚拟切片上的更新值依照虚拟切片的空间位置组合为一个虚拟更新数据体,若在原重建数据体坐标系下的坐标为V(x,y,z),在所述虚拟更新数据体坐标系下的坐标为V′(x′,y′,z′),则V′(x′,y′,z′)=TN·V(x,y,z)其中TN为纹理矩阵,根据原重建数据体坐标系和虚拟更新数据体坐标系的位置关系由OpenGL图形开发语言提供的基本操作矩阵依先后顺序相乘得到,所述先后顺序是指该原重建数据体坐标系通过平移或/和旋转坐标的方式重合到虚拟更新数据体坐标系的先后顺序;步骤(III.c),将所述的虚拟更新数据体和步骤(III.b)的得到的虚拟更新数据体坐标系下的坐标V′(x′,y′,z′)代入OpenGL图形开发语言提供的三维纹理索引函数,即可得到步骤(II.b)所述的体素V(x,y,z)的更新值dv;步骤(N),用步骤(II.b)或步骤(III.b)所述的更新值dv更新原重建数据体的值。
地址 100084北京市100084-82信箱