【关键词】 体素; 插值; FDK算法
A new interpolation method in the FDK algorithmZHAI Jing, PAN Jinxiao
(National Key Lab for Electronic Measurement and Technology, North University of China,
Taiyuan 030051, China)
Abstract: In the backprojection process of the FDK algorithm, every voxels project distribution in the detector have itself characteristics, This paper presents a new interpolation method. The voxel projection based on the characteristics used in the reconstruction process, According to its different scanning angle detector modules in all the occupied area of the projection of the size as the voxel projection. Actual experimental results show that the FDK algorithm of this new interpolation method give an better result in the reconstructed images verge than the traditional interpolation methods (such as replication interpolation, bilinear interpolation), and this new interpolation method can restrain noise effect.
Key words: voxelinterpolationFDK algorithm
引言考试大
在FDK算法的反投影过程中,由于数据的离散性,会出现象素的投影地址“对不准”现象,一般需要插值运算。插值是指在己知的坐标范围内,一种基于模型的从离散数据估计连续数据的方法。经典的线性插值技术包括最邻近插值(replication)[1],双线性插值(bilinear)[2],双三次(Bicubic)插值[1,3~6]等。本文考虑了三维重建图像的每个像素在不同角度时在探测器上的近似面积,(此近似面积大于1)并将在此面积的不同探测器上的投影值的加权和作为该像素在此扫描角度的投影值。
1 像素投影图形的取法
通常情状下,我们认为物体的像素在探测器上的投影是一个点,但在实际中,在某一个扫描角度下,经过每个体素的射线在探测器上形成一个几何图形。在重建过程中,这个几何面的形状不规则而且其面积很不易求出,因此,我们可近似考虑每个像素的一部分点在探测器上所形成的图形。具体步骤如下:
如图1 所示,在锥束圆轨迹扫描结构中,定义锥束的投影角为β,扇角为γ,锥角为κ。设ABCD—EFGH是要重建图像的某一个体素(i,j,k),如图2所示,A: i-12,j+12,k+12,
B: i-12,j-12,k+12,
C: i+12,j-12,k+12,
D: i+12,j+12,k+12,
E: i-12,j+12,k-12,
F: i-12,j-12,k-12,
G: i+12,j-12,k-12,
H: i+12,j+12,k-12.当射线源介于x正半轴到y正半轴之间这12π弧度(即0°≤β<90°)时,考虑ACGE这个对角面在探测器上的投影图形。当射线源介于y正半轴到x负半轴之间这12π弧度(即90°≤β<180°)时,考虑BDHF这个对角面在探测器上的投影图形。当射线源介于x负半轴到y负半轴之间这12π弧度(即180°≤β<270°)时,考虑ACGE这个对角面在探测器上的投影图形。当射线源介于y负半轴到x正半轴之间这12π弧度(即270°≤β<360°)时,考虑BDHF这个对角面在探测器上的投影图形。分别计算在上述条件下,其对角面的投影在探测器上分布的情状。
2 像素反投影值的计算考试大考试大http:/
一般认为,像素的反投影值是由点的插值取得。有近邻插值、双线性插值等。在本文中,我们考虑上诉投影图形的加权值。任取某一待重建体素上的顶点(x,y,z),β为扫描角度,lso为探源到物体中心距离, lso2为探源到探测器距离,那么它在探测器上的落点p的坐标[7]: x′=(lso2/(lso-x×cos(β)+y×sin(β)))
·(x×sin(β)-y×cos(β)) (1)
y′=(lso2/(lso-x×cos(β)+y×sin(β)))×z (2)因为物体离射线源距离比较远而离探测器又很近,又由经验可知,当射线源的扫描角度是14π、34π、54π、74π时,投影的几何面达到,大约就是一个宽是1、长是2的近似矩形,当射线源的扫描角度是12π、π、32π、2π时,投影的几何面投影的几何图形的面达到最小,大约就是一个宽是1、长是1的近似矩形。如图3所示,其在探测器上的分布共有六种情况,根据具体情况计算出该投影在每个探测器的探元中的面积,记为Si,设p′(xi,yi,β)是在扫描角度β时该面积上的投影值,(i,j,k,β)是在β时体素(i,j,k)要取的反投影值,N是投影面占据探元的个数,其值是6。我们得到计算像素反投影值的公式:(i,j,k,β)=Ni=1p′(xi,yi,β)*Si.(3)图3 像素投影在探测器上的近似分布情况示意图
3 实验结果
实验采用220 kv,10 mA的X射线源。探元的大小为0.127 mm。采用的探测器为 PAXSCAN 2520。工作模式:数字视频。数据类型:unsigned short。A/D:12bit。射线源—标准件—探测器间距:850 mm~200 mm,旋转一周采样间隔为1度,某一角度下的投影如图3示,大小为256*256,分别采用最近邻插值,双线性插值以及本文提到的新插值方法重建图像,图像大小为256*256。
由实验结果表明:新的插值方法比我们通常用的最近邻和双线性插值法对图像有明显的改良,由图像的灰度曲线也可表明新方法得到的图像边缘要好于前两种方法。也能看出,由新插值得到图像的边缘比较清晰,而且还有抑制噪声的效果。让反投影的信息量的值由与周围差别比较大的某一探元上的获取改为在周围的信息量更接近的相邻几个探元上综合获取,这样能够有效地抑制噪声和孤立点,但是同时图像的对比度也会有所下降。
4 总结
通过以上的分析和实验表明,在基于圆轨迹的锥束CT扫描和重建过程中,如果考虑待重建体素落在探测器上一定的投影面积,并且在反投影重建这个体素点的时候,考虑那些探元上的信息量,将会对重建图像的质量有很大的提高,不仅图像的边缘更清晰,而且还能达到抑制噪声的目的,但是图像的对比度有所下降。当然,这种新的插值方法不仅适用于圆轨迹的扫描方,而且对所有的锥束扫描方式,并采用平板探测器采集数据的各种反投影的锥束重建都是适用的。
【参考文献】
[1]Parker J A, Kenyon R V, Troxe L D E. Comparison of interpolating methods for image resampling. IEEE Transaction on Medical Imaging. 1983, 2(1), 31-390.
[2]Jain A K. Fundaments of Digital Image Processing. Englewood Cliffs. NJ: PrenticeHall, 1989.
[3]Chen T C, deFigueiredo R J P. Twodimensional interpolation by generalized spline filters based on partial differential equation image models. IEEE Transaction on Acoustics, Speech,Signal Processing (ASSP). 1985, 33(3),31-642.
[4]Hou H S, Andrews H C. Cubic splines for image interpolation and digital filtering. IEEE Transaction on Acoustics, Speech, Signal Processing (ASSP). 1978, 26(6), 508-517.
[5]G Keys R. Cubic convolution interpolation for digital image processing. IEEE Transaction on Acoustics, Speech, Signal Processing (ASSP). 1981,29(6),1153-11600.
[6]Unser M, Aldroubi A, Eden M. Fast bspline transforms for continuous image representation and1991 interpolation. IEEE Transaction on Pattern Analysis and Machine Intelligence (TPAMI), 13(3),277-285.[7]孙怡,侯颖,胡家升.体积CT投影数据的模拟方法 [J]. CT理论与应用研究,2005,14(1-6).
java产生随机数的几种方式n 在j2se里我们可以使用Math.random()方法来产生一个随机数,这个产生的随机数是0-1之间的一个double,我们可以把他乘以一定的数,比如说乘以100,他就是个100以内的随机,这个在j2me中没有。
n 在java.util这个包里面提供了一个Random的类,我们可以新建一个Random的对象来产生随机数,他可以产生随机整数、随机float、随机double,随机long,这个也是我们在j2me的程序里经常用的一个取随机数的方法。
n 在我们的System类中有一个currentTimeMillis()方法,这个方法返回一个从1970年1月1号0点0分0秒到目前的一个毫秒数,返回类型是long,我们可以拿他作为一个随机数,我们可以拿他对一些数取模,就可以把他限制在一个范围之内啦
其实在Random的默认构造方法里也是使用上面第三种方法进行随机数的产生的
对于方法二中的Random类有以下说明:
java.util.Random类有两种方式构建方式:带种子和不带种子
不带种子:
此种方式将会返回随机的数字,每次运行结果不一样
public class RandomTest {
public static void main(String[] args) {
java.util.Random r=new java.util.Random()
for(int i=0i<10i++){
System.out.println(r.nextInt())
}
}
带种子:
此种方式,无论程序运行多少次,返回结果都是一样的
public static void main(String[] args) {
java.util.Random r=new java.util.Random(10)
for(int i=0i<10i++){
System.out.println(r.nextInt())
}
}
两种方式的差别在于
(1) 首先请打开Java Doc,我们会看到Random类的说明:
此类的实例用于生成伪随机数流,此类使用 48 位的种子,该种子可以使用线性同余公式对其进行修改(请参阅 Donald Knuth 的《The Art of Computer Programming, Volume 2》,第 3.2.1 节)。
如果用相同的种子创建两个 Random 实例,则对每个实例进行相同的方法调用序列,它们将生成并返回相同的数字序列。为了保证实现这种特性,我们为类Random指定了特定的算法。为了 Java 代码的完全可移植性,Java 实现必须让类 Random 使用此处所示的所有算法。但是允许 Random 类的子类使用其他算法,只要其符合所有方法的常规协定即可。
Java Doc对Random类已经解释得非常明白,我们的测试也验证了这一点。
(2) 如果没有提供种子数,Random实例的种子数将是当前时间的毫秒数,可以通过System.currentTimeMillis()来获得当前时间的毫秒数。打开JDK的源代码,我们可以非常明确地看到这一点。
/**
* Creates a new random number generator. Its seed is initialized to
* a value based on the current time:
* Random() { this(System.currentTimeMillis())}java.lang.System#currentTimeMillis()
*/
public Random() { this(System.currentTimeMillis())}
另外:
random对象的nextInt(),nextInt(int n)方法的说明:
int nextInt()
返回下一个伪随机数,它是此随机数生成器的序列中均匀分布的 int 值。
int nextInt(int n)
返回一个伪随机数,它是从此随机数生成器的序列中取出的、在 0(包括)和指定值(不包括)之间均匀分布的 int值。
欢迎分享,转载请注明来源:内存溢出
评论列表(0条)