一种基于多GPU的高速计算全息图的方法

摘要

本发明涉及一种基于多GPU的高速计算全息图的方法，利用多个GPU同时计算，并且运用基于点元法的表格数据压缩方法计算全息图，在计算LUT表格数据时，将建立的三维空间分成多个虚拟的二维截面图层，由于每个二维平面所有点对应的传播方向的调制因子相同，所以在记录查找表数据时，将Z方向调制因子提取出来，只需记录二维平面内水平和竖直方向的调制因子，减少了层数倍的存储空间。结合并行实现压缩查找表的方法，不仅减少了查找表的存储空间，还加快了全息图的计算速度。同时，编程建立在主流GPU架构上，方便开发，利于推广。可广泛应用于各图形处理器当中。

说明书

技术领域

本发明涉一种高速计算全息图的方法，尤其涉及一种基于多GPU系统的全息图的快速计算方法。 

背景技术

国内外在全息图计算加速方面主要有两种方式，一是利用硬件进行快速运算，二是基于全息图计算原理，利用快速算法进行加速。 

新加坡的Yuechao Pan小组2009年在GPU上实现了该小组提出的加速算法S-LUT，通过GPU并行运算，达到了实时的计算速度。算法基于LUT方法，通过分析菲涅尔衍射的计算公式，将x方向与y方向分离开独立计算，记录x和y方向的调制因子，计算时，先将两个调制因子相乘，然后叠加得到整个全息图。 

现有的两种方法存在的问题是，首先记录的数据量仍然很大，较大的数据量会超出CPU缓存的容量；并且因为CPU的计算效率远远低于GPU，因此CPU从内存中读取数据的时间会增加数据的读取时间，从而降低计算速度；最后编程模块不通用，需要单独设计，且对硬件要求也高，因此不适合推广。 

发明内容

(一)要解决的技术问题 

本发明要解决的技术问题就是如何基于多GPU的高速计算全息图的方法，结合快速算法，提出一种将编程建立在主流GPU架构上的高速计算全息图的方法。 

(二)技术方案 

为了解决上述技术问题，本发明提供了一种基于多GPU的高速计 算全息图的方法，其特征在于，首先在Off-line下预先计算全息图，并存储于表格中；然后，在In-line下，首先从表格中读取构建目标物体点对应的X，Y方向调制因子，然后乘以对应的Z方向的调制因子，从而获得所需点的全息图； 

其中，Off-line的计算过程包括以下步骤： 

步骤一：在主机端申请LUT和设备端的内存空间； 

步骤二：建立LUT空间，并将数据传输到设备端； 

步骤三：调用设备端的kernel函数，该函数在设备端进行运算； 

步骤四：设备端保存主机端传来的物体数据； 

步骤五：设备端计算X，Y方向的调制因子，每个线程索引分别对应X，Y的值； 

步骤六：将计算结果传回主机端； 

步骤七：保存设备端传回的调制因子数据，释放申请的主机端和设备端空间； 

In-line的计算过程包括以下步骤： 

步骤一：在主机端申请LUT和设备端的内存空间； 

步骤二：读取物体数据和LUT中保存的数据； 

步骤三：将数据信息传输到设备端； 

步骤四：调用设备端的kernel函数，该函数在设备端进行运算； 

步骤五：设备端保存主机端传来的物体数据和LUT空间数据； 

步骤六：设备端进行运算，对于同一二维面上的物点，首先计算所有相同X坐标的点的Y方向的调制因子之和，然后乘以相应的相同Y坐标的点的X方向的调制因子之和，最后乘以对应平面的Z方向的调制因子便可以得出该平面物体的全息图； 

步骤七：将计算结果传回主机； 

步骤八：保存设备端传回的全息图数据，释放申请的主机端和设备端空间，保存并显示计算结果。 

优选地，所述Off-line或In-line的计算过程中，首先划分并行计算的任务，将全息图平面上的各个像素点分别交给GPU的各个thread进行计算。 

(三)有益效果 

本发明的一种基于多GPU的高速计算全息图的方法，利用多个GPU同时计算，并且运用基于点元法的表格数据压缩方法计算全息图，在计算LUT表格数据时，将建立的三维空间分成多个虚拟的二维截面图层，由于每个二维平面所有点对应的传播方向的调制因子相同，所以在记录查找表数据时，只需记录二维平面内水平和竖直方向的调制因子，减少了层数倍的存储空间。结合并行实现压缩查找表的方法，不仅减少了查找表的存储空间，还加快了全息图的计算速度。同时，编程建立在主流GPU编程架构上，方便开发，利于推广。 

附图说明

为了更清楚地说明本发明实施例或现有技术中的技术方案，下面将对实施例或现有技术描述中所需要使用的附图作简单地介绍，显而易见地，下面描述中的附图仅仅是本发明的一些实施例，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些附图获得其他的附图。 

图1：Off-line下的计算流程示意图； 

图2：in-line下的计算流程示意图； 

图3：计算全息图的任务划分结构示意图； 

具体实施方式

下面结合附图和实施例对本发明的实施方式作进一步详细描述。以下实施例用于说明本发明，但不能用来限制本发明的范围。点元法的基本思想是：将物体离散成独立的物点，将每个物点看成是点光源，利用球面波传播公式，分别计算每个物点传播到全息图平面上的光场分布，然后把每个物点所对应的全息图上的光场分布按物体的组合形式进行线性叠加，叠加的最后结果，即为物体传播到全息图上的光场分布。全息图平面的复振幅分布可以由式(1)求得： 

$H(x,y)=\underset{j=1}{\overset{N}{\Sigma }}{A}_j\exp \left[i({\mathrm{kr}}_j+{\phi }_j)\right]---\left(1\right)$

式中，H(x,y)是组成三维物体的各个离散点光源在全息图平面的复振幅分布的叠加。N为抽样点数，A_j为点光源的光强，φ_j为光源的初始相位(为了使物体漫反射化，通常设为随机位相)，r_j为三维物体上第j个点光源(x_j,y_j,z_j)与全息图平面上(x,y,0)点的距离，x_j,y_j,z_j分别是物体空间内，三维物体在x,y,z方向的坐标，x,y分别是全息图空间内，全息图上的x,y,z方向的像素坐标，d为三维物体所在的空间坐标系到全息图平面的垂直距离。： 

$r_j=\sqrt{{(x-x_j)}^2+{(y-y_j)}^2+{(d+z_j)}^2}---\left(2\right)$

表格数据压缩方法是基于夫琅禾费衍射。对于再现像，我们只关注它的强度信息，因此可以在菲涅尔衍射区观察到夫琅禾费衍射。全息平面上的复振幅可以近似表示为： 

$H(x,y)=\underset{j=0}{\overset{N-1}{\Sigma }}{A}_j\exp \left(\mathrm{ik}\frac{x^2+y^2}{2(d-z_j)}\right)\times \exp (-\mathrm{ik}\frac{x_{j}x+y_{j}y}{d-z_j}),\left(3\right)$

如果d＞＞z_j，上式可近似并进一步分解为： 

$\left(\begin{array}{c}H(x,y)=\underset{j=0}{\overset{N-1}{\Sigma }}{A}_j\exp \left(\mathrm{ik}\frac{x^2+y^2}{2(d-z_j)}\right)\exp (-\mathrm{ik}\frac{x_{i}x+y_{j}y}{d})\\ =\underset{j=0}{\overset{N-1}{\Sigma }}{A}_j\exp \left(\mathrm{ik}\frac{x^2+y^2}{2(d-z_j)}\right)\exp (-\mathrm{ik}\frac{x_{j}x}{d})\exp (-\mathrm{ik}\frac{y_{j}y}{d})\end{array}\right),\left(4\right)$

定义exp(ik(x²+y²)2(d-z_j))为Z方向的调制因子，H(x)＝exp(-ikxx_j/d)为x向的调制因子，V(y)＝exp(-ikyy_j/d)为y向的调制因子。则式(4)可以化简为： 

$H(x,y)=\underset{j=0}{\overset{N-1}{\Sigma }}{A}_{j}H\left(x\right)V\left(y\right)L\left(z\right),\left(5\right)$

可以认为三维物体是由N_z个二维点阵组成，每一个二维点阵共有物点N_xy个。对于每一层上的采样点，它们具有相同的轴向调制因子。对于每一层，若有N_y个点分布在同一垂直方向，那么这些点具有相同的水平方向的调制因子，如果每层有N_x个不同的垂直方向，最后可得全息图计算公式： 

$H(x,y)=\underset{j_z=0}{\overset{N_z-1}{\Sigma }}\left\{\underset{j_x=0}{\overset{N_x-1}{\Sigma }}\right[\underset{j_y=0}{\overset{N_y-1}{\Sigma }}{A}_{\mathrm{jy}}V(y,y_{\mathrm{jy}})\left]H(x,x_{j}x)\right\}L(z,z_{\mathrm{jz}})---\left(6\right)$

计算主要分为两步，第一步，将Z方向的调制因子提取出来，在Off-line下预先计算X和Y方向调制因子，并存储于表格中。第二步，在In-line下，首先从表格中读取构建目标物体点对应的X，Y方向调制因子，然后乘以对应的Z方向的调制因子，从而获得所需点的全息图。 

Off-line的计算过程，如图1所示： 

步骤一：在主机端申请LUT和设备端的内存空间； 

步骤二：建立LUT空间，并将数据传输到GPU端； 

步骤三：调用设备端的kernel函数，该函数在GPU设备端进行运算； 

步骤四：GPU端保存CPU端传来的物体数据； 

步骤五：GPU端并行计算x和y方向的调制因子，每个线程索引分别对应x和y的值； 

步骤六：将计算结果传回主机。 

步骤七：保存设备端传回的调制因子数据，释放申请的主机端和设备端空间。 

In-line的计算过程，如图2所示： 

步骤一：在主机端申请LUT和设备端的内存空间； 

步骤二：读取物体数据和LUT中保存的数据； 

步骤三：将数据信息传输到GPU端； 

步骤四：调用设备端的kernel函数，该函数在GPU设备端进行运算； 

步骤五：GPU端保存CPU端传来的物体数据和LUT空间数据； 

步骤六：GPU端进行运算，对于同一二维面上的物点，首先计算所有相同x坐标的点的y方向的调制因子之和，然后乘以相应的相同y坐标的点的x方向的调制因子之和，最后乘以对应平面的z方向的调制因子便可以得出该平面物体的全息图。 

步骤七：将计算结果传回主机。 

步骤八：保存设备端传回的全息图数据，释放申请的主机端和设备端空间，保存并显示计算结果。 

GPU实现该算法时，在Off-line和In-line下，都可以用并行思想进行优化。计算时，首先需要划分并行计算的任务。针对全息图平面上的M×N个抽样点，可以将其交给M×N个threads进行计算，每个thread负责一个全息图抽样点的计算，这些threads分布于不同的线程块block和线程格grid中。其划分方式如图3所示，A是线程格grid组块，B是线程块block组块，C表示全息图。 

本实施例用并行的思想计算全息图，比传统的串行算法更快；高速的硬件GPU系统具有比CPU更高的浮点运算能力，对于计算性的工作，GPU效率比CPU高几百倍；多个GPU并行计算可以计算大尺寸的全息图及较大数据量的三维场景，提高显示的质量；提出的快速算法有效的节省了存储空间，降低了数据传输造成的时间 损耗，也加快了记录的速度。 

以上实施方式仅用于说明本发明，而非对本发明的限制。尽管参照实施例对本发明进行了详细说明，本领域的普通技术人员应当理解，对本发明的技术方案进行各种组合、修改或者等同替换，都不脱离本发明技术方案的精神和范围，均应涵盖在本发明的权利要求范围当中。