图形学基础(5)——Sobol 采样
https://blog.csdn.net/qq_33274326/article/details/69952035
图形学基础(5)——Sobol 采样
2017年04月10日 17:24:24 OscarShen09 阅读数:2695 标签: 图形 采样器 随机数 光线追踪 更多
个人分类: 图形
一个光源会在三维空间以及角度空间中随机产生光线分布。比方说,一个点光源的发出位置 origin 不变,而光线方向按照余弦来均匀分布。当进行渲染时,必须发出足够多个光线,才能精确描述这个点光源。
光线的随机产生通常使用随机数产生器,随机数产生器的目标是制造一系列互相无关的随机数,然后追踪大量的随机光线,在运算量较小的情况下获得尽可能高的渲染结果。
所有基于现代 CPU 的随机数生成算法都是伪随机的(quasi-random)。它们受限于一个周期。当超过周期后就会重复出现,而不再是相互无关的随机数。这个周期的最终限定是由电脑的位数来决定的,因此,没有一个内建的随机数是“真正”随机的。
Sobol 采样使用了不同的方式来采样。比起随机数,Sobol 序列着重于在概率空间中产生均匀的分布。但这并不是单纯的使用网格来填满,而是使用一个本质上随机,但是巧妙的方法去“填满”概率空间,即之后产生的随机数会分布到之前没有采样到的区域。
Radical Inversion 运算和 Van der Corput 序列
本节来自低差异序列(一)- 常见序列的定义及性质,欢迎读者移步查看原文。
Radical Inversion 运算的定义:
i=∑l=0M−1al(i)bl
Φb,C(i)=(b−1...b−M)[C(a0(i)...aM−1(i))T]
上式中,如果将任意一个整数 i 表示成 b 进制的数字,然后把得到的数按位 al(i) 排成向量,再用这个向量和生成矩阵 C 相乘得到一个新向量,最后把新向量镜像到小数点右边得到另一个范围在 [0,1) 的数字,就称为 RadicalInversion 运算,记为 Φb,C(i)。
如果将上述过程简化一下,另 C 为单位矩阵,可以得到 Van der Corput 序列:
Φb(i)=(b−1...b−M)(a0(i)...aM−1(i))T=∑l=0M−1al(i)b−l−1
举个例子,正整数8以2为底数的radical inverse的计算过程如下。首先算出8的2进制表示,1000。此处假设 C 单位矩阵,所以直接将1000镜像到小数点右边,0.0001。这个二进制数的值就是最终结果,把它转换回10进制就得到1/16, 即 Φ2(8)=1/16=0.0625。下面的表给出了更多的以2为底数的Van der Corput序列的例子。
Van der Corput序列有几个属性:
- 每一个样本点都会落在当前已经有的点里“最没有被覆盖”的区域。例如 Φ2(3)=3/4 是刚好落在了 [0,1) 区间中被覆盖最少的区域(Φ2(0)=0,Φ2(1)=1/2,Φ2(2)=1/4)。
- 样本个数达到 bm 个点时对 [0,1) 会形成 uniform 划分。
- 很多时候并不能够代替伪随机数,因为点的位置和索引有很强的关系,例如在以2为底的Van der Corput序列中,索引为基数时候序列的值大于等于0.5,偶数时则小于0.5。
产生 Sobol 序列
本节来自低差异序列(二)- 高效实现以及应用,由于原文写很棒,我也就不献丑了,直接贴过来。
Sobol序列的每一个维度都是由底数为 2 的 radical inversion 组成,但每一个维度的radical inversion都有各自不同的矩阵 C
。
因为完全以2为底数,所以Sobol序列的生成可以直接使用bit位操作实现radical inversion,非常高效。Sobol序列的分布具有不仅均匀,而且当样本的个数为2的整数次幂时,在 [0,1)s
区间中以2为底的每个Elementary Interval中都有且只会有一个点,这意味着它可以生成和Stratified Sampling和Latin Hypercube同样高质量分布的样本(见下图),同时又不需要预先确定样本的数量或者将样本储存起来,并可以根据需要生成无限个样本,非常适合progressive的采样。这些性质也使得Sobol在需要一切对高维空间采样的应用中,例如图形,渲染以及金融领域,都非常流行。
Radical Inversion 实现
目前介绍过的所有序列都基于Radical Inversion这个操作,而这个操作尽管公式略微复杂,但它的实现非常直观和简单,下面贴上代码以及注。
double IntegerRadicalInverse(int Base, int i)
{
int numPoints, inverse;
numPoints = 1;
// 此循环将i在"Base"进制下的数字左右Flip
for(inverse = 0; i > 0; i /= Base)
{
inverse = inverse * Base + (i % Base);
numPoints = numPoints * Base;
}
// 除以Digit将这个数镜像到小数点右边
return inverse / (double) numPoints;
}
- 1
- 2
- 3
- 4
- 5
- 6
- 7
- 8
- 9
- 10
- 11
- 12
- 13
- 14
Sobol序列的所有维度都是基于2为底数的Radical Inversion。只不过每个维度有自己不同的生成矩阵(Generator Matrix)C
。因为是以2为底数,将数字从二进制中提取每一位的数字,以及和矩阵C
做运算,都可以用位操作(右移,异或等)来完成,非常高效。
double Sobol(uint64 i, uint Dimension)
{
double r;
// 将i依次右移,提取2进制里的每一位
for (uint k = 0; i; i >>= 1, k++)
if (i & 1) // 若当前位为1,则用异或和矩阵相乘
r ^= C[Dimension][k];
return r / (double) (1 << M); // 除以2^M,移到小数点右边
}
- 1
- 2
- 3
- 4
- 5
- 6
- 7
- 8
- 9
那么如何计算出能生成如此高质量分布的矩阵呢?Quasi Monte Carlo的学者们已经花了数10年的时间搜索这种矩阵,现在我们可以在这个网页(Sobol sequence generator )找到可以生成21201维度的Sobol序列的矩阵。