基于快速k近邻的光子映射算法研究

王海波

摘要：光辐射强度估算是光子映射算法一个关键技术，传统使用简单、有效的k近邻（kNN）算法，但kNN具有计算复杂度高，内存需求量的缺点，新算法针对kNN的缺点，改进kNN搜索光子的方式，先将空间分割为多个固定长度的立方体，每个立方体体包含一定数量的光子数，通过测试各个立方体与光线接触点之间的位置搜索接触点周围的k个最近邻光子，进而估算光辐射强度，实验表明新算法搜索光子的速度更快，而且图形清晰度更高。

关键词：光子映射;光辐射强度估算;kNN;空间网格

中图分类号：TP3   文献标识码：A

文章编号：1009-3044（2019）28-0286-03

光子映射是较好的真实感图形渲染方式之一，与光线跟踪方法比较，光子映射能较好地渲染辉映、焦散效果[1]。光子映射分为两个阶段第一个阶段发射光子，跟踪光子，建立光子图;第二阶段利用光子图估算光照，从图形像素的角度发出光线，如果遇到反射或折射后，记录接触点，搜索接触点周围的光子，对接触点进行光辐射强度估算，渲染图形。因此光辐射强度估算是光子映射算法的第二阶段的一个关键技术。

为了准确估算光辐射强度，需要寻找到对接触点光辐射强度有影响的光子，传统的光子映射算法采用K-近邻算法对光辐射强度进行估算，k-近邻算法具有简单、有效的优点，但发射光子量較大时k-近邻算法的计算时间比较久。

为了改进计算精度，提高计算效率Garcia采用常核函数方法，并且放弃第K个光子辐射通量的计算，以避免偏差[2][3];PerChristensen主张对第K个光子计算1/2的光辐射通量[4];Reinhard Klein主张根据第k-1与第k个光子的平均距离来计算第k个光子对接触点的影响[5];对于核函数可选择圆锥核函数、高斯核函数、Epanechnikov核函数、Silverman核函数[6][7][8][9]。为实现光辐射强度估算无无偏差性，Hachisuka等提出渐进式光子映射算法，可从理论上对光辐射强度的计算达到无偏差计算，但需要循环发射大量光子[10][11]。

kNN算法具有简单、有效的特点，因此如何利用kNN算法的优点，提高光辐射强度估算也是新算法研究的方向之一。本算法先将空间分割为多个固定长度的立方体，每个立方体包含一定数量的光子数，通过测试各个球体与接触点之间的位置速搜索测试点周围的k个最近邻光子，进而估算光辐射强度。

1 光子映射

光子映射是一种全局光照，分两阶段。第一阶段从场景的光源发出大量的光子，跟踪光子建立光子图，建立密度低的全局光子图及密度高的焦散光子图。焦散光子图要求光子至少被镜面反射或折射一次，全局光子图存储各种路径的光子，当光子击中漫反射表面或略微光滑表面时，光子信息都会被存储到相应的光子图中。光子图存储每个光子，并记录光子能量或通量、输入方向和位置等信息，光子图采用平衡二叉树作为光子搜索的数据结构。

第二阶段是利用光子图计算光辐射强度，即从图形的每个像素发出光线，光线在场景中经过若干次镜面反射、折射后、漫反射后，记录接触点x既估算点，然后利用全局光子图及焦散光子图，搜索离估算点最近的若干光子如图1所示，并计算估算点的光辐射强度如式（1）所示

其中，x是估算点，fr是双向反射分布函数（BRDF），[Φi]是第i个光子的光通量，[rk]是x与 k个光子中距x最远的某个光子的距离。当光子数目越多，则光辐射强度的计算精度越高。因此光子的发射量影响到光辐射强度的计算，如何快速查找光子成为人们研究的方向。

2 快速K近邻模型

K近邻（KNN）算法是非常有效的基于距离的算法，被广泛应用于回归模型，主要思想是：待测样本的特征属性等于最近邻的k个样本特征属性的平均值。假设训练集合A包含s个样本，每个样本又含有t个属性既：{Ai=（Ai1，Ai2，…，Ait），i=1，2，…，s };待测样本为B，属性为（B1，B2，…，Bt）;则求得样本B的特征属的步骤如下：

3 快速kNN

3.1 构建新算法

新算法使用空间网格的方法先将空间划分为若干立方体，将所有光子置入到这些立方体中。立方体网格单元太大太小，都会对整个查找过程产生不良影响，若立方体网格单元太小，会增加存储量，降低效率;若立方体网格单元太大，则每个立方体单元会包含过多面片，对求交造成困难，因此新算法划分立方体的个数为p=n/t其中n是光子总数，t是未知数，一般取值为20，立方体的边长为L，满足[L3=n/t]。如果有多个不包含光子的立方体连续在一起，则合并为一个立方体，保证生成的空间单元数不超过O（n），如算法1所示。接着搜索估算点周围的k个立方体，并从k个立方体中搜索最近的k个光子，完成对估算点光辐射强度的估算如算法2，图2所示。

3.2 算法分析

设发射光子数n，k为kNN算法参数，m为估算点的个数。传统KNN的时间复杂度为O（mnk），表示每个估算点要计算同n个光子的距离，同时为求出k个最近邻的光子，内存中要维持一个k长度的插入排序表，在排序表中，每插入一个新值，对表的最大操作次数为k。新算法中设立方体的总数为p，时间复杂性包含求k个最小立方体及其所包含样本中的k个最近邻样本。每个立方体平均包含的光子数为n/p，因此时间复杂性为O（mpk+mnk2/p）=O（mk（p+nk/p）），故当p+nk/pp2/（p-k）时，新算法的时间复杂性低于传统kNN算法，由于k取值同立方体数相比要小得多，故当n>p>k时，实际取值p<=n/20，因此n>p>k条件是成立的，新算法的运行效率优于传统的kNN算法。

4 算法实现

采用vs2013和OpenGL的编程环境，在一台配置为Intel（R） Core（TM） i5，8GB内存，NVIDIA GeForce 610M显卡，win7下进行实验。

实验所用的测试场景如图1，图2，图3 所示，表1给出了各图的参数及渲染时间，其中图1的发射的光子数为10M，图2发射的光子数为15M，图3发射的光子数为12M。由表中可以看出：新算法与传统kNN的算法发射的光子数及需搜索的光子数是一样的，新算法的每个立方体包含的光子数为20，从表中可以看出新算法的渲染时间比传统kNN算法要快，图1快24.7%，图2快33.0%，图3快23.1%，图2增快的幅度最大是因为空间分割最紧密。从图1、图2、图3的图形渲染质量方面来看新算法也比传统的kNN算法更清晰。

5 总结

新算法针对kNN搜索光子的计算量大的缺点，改进了kNN搜索光子的方式，先将空间分割为多个固定长度的立方体，每个立方体体包含一定数量的光子数，通过测试各个球体与接触点之间的位置速搜索测试点周围的k个最近邻光子，进而估算光辐射强度，实验表明新算法搜索光子的速度更快，而且图形清晰度更高。

参考文献：

[1] Per H Christensen， Henrik W Jensen， Toshi Kato， Frank Suykens. A Practical Guide to Global Illumination Using Photon Mapping [R].USA： Siggraph， 2002

[2]Garcia， R.， Urena， C.， Sbert， M.. Description and solution of an unreported intrinsic bias in photon mapping density estimation with constant kernel[J]. Comput. Graph. Forum， 2012，31（1）：33-41.

[3]Garcia， R.， Urena， C.， Poch， J.， et al.， 2014. Overestimation and underestimation biases in photon mapping with non-constant kernels[J]. IEEE Trans. Visual. Comput. Graph.， 20（10）：1441-1450.

[4]Jensen H W， Christensen P. High quality rendering using ray tracing and photon mapping[C]// Siggraph 07： Acm Siggraph Courses， Acm. 2007.

[5] R. Garc??a， C. Ure ? na， and M. Sbert， “Description and solution of an unreported intrinsic bias in photon mapping density estimation

with constant kernel.” Presented at Eurographics 24th Symposium on Rendering. Invited paper from Computer GraphicsForum， 2013.

[6] Roland， S. Bias compensation for photon maps Comput[J]. Graph. Forum， 2003，22（4）：729-742.

[7]Myszkowski， K.， 1997. Lighting reconstruction using fast and adaptive density estimation techniques. Proc. Eurographics Workshop on Rendering Techniques， p.251-262.

[8] Schjoth， L.， 2009. Anisotropic Density Estimation in Global Illumination. PhD Thesis， University of Copenhagen， Denmark

[9]Schjoth， L.， Sporring， J.， Olsen， O.F.， 2008. Diffusion based photon mapping.

[10] T. Hachisuka， S. Ogaki， and H. W. Jensen， “Progressivephoton mapping，” ACM Trans. Graph.， vol. 27， pp. 130：1–130：8， December 2008.

[11] C. Knaus and M. Zwicker， “Progressive photon mapping： A probabilistic approach，” ACM Trans. Graph.， vol. 30，no. 3， pp. 25：1–25：13， May 2011.

【通聯编辑：梁书】